Changeset 1519 for trunk/LMDZ.VENUS/libf

Timestamp:

Mar 23, 2016, 12:09:36 PM (9 years ago)

Author:

slebonnois

Message:

SL: volcanisme et traceurs sur Venus

Location:

trunk/LMDZ.VENUS/libf/phyvenus

Files:

• physiq.F (modified) (1 diff)
• phytrac_emiss.F (modified) (4 diffs)
• phytrac_relax.F (modified) (5 diffs)

trunk/LMDZ.VENUS/libf/phyvenus/physiq.F

@@ -227 +227 @@
       EXTERNAL diagphy
       EXTERNAL mucorr
-      EXTERNAL phytrac
       EXTERNAL nirco2abs
       EXTERNAL nir_leedat

trunk/LMDZ.VENUS/libf/phyvenus/phytrac_emiss.F

@@ -33 +33 @@
       use dimphy
       USE comgeomphy
+      USE chemparam_mod,only:M_tr
       IMPLICIT none
 #include "YOMCST.h"
@@ -65 +66 @@
 cAA ----------------------------
 
-C les traceurs
-
-c===================
-c it--------indice de traceur
-c k,i---------indices long, vert
-c===================
 c pour emission volcan
       real :: deltatr(klon,klev,nqtot)
 
-      integer,parameter :: nbsrc=2,nblat=5,nblon=4
-!     integer,parameter :: Nemiss=1   ! duree emission (Ed)
-      integer,save :: Nemiss(nbsrc)      ! duration emission (Ed)
-      real,save :: source_volcan(nbsrc)  ! flux emission (kg/s)
+      integer,parameter :: nblat=5,nblon=4,nbz=3
+      integer,parameter :: Nemiss=1     ! duree emission (Ed)
+      integer,save :: Nheight(nbz)      ! layer emission
+      real,save :: so2_quantity         ! quantity so2 (kg)
       real,save :: lat_volcan(nblat),lon_volcan(nblon)
       real,save :: area_emiss(nblat,nblon)
       integer,save :: ig_volcan(nblat,nblon)
 
-c======================================================================
-
       INTEGER i, k, it
-      integer ilat,ilon,iemiss
+      integer ilat,ilon,iz
       real    deltalat,deltalon
-
-c Variables liees a l'ecriture de la bande histoire physique
-
-c Variables locales pour effectuer les appels en serie
-c----------------------------------------------------
-
-      REAL d_tr(nlon,nlev) ! tendances de traceurs 
-
-      character*20 modname
-      character*80 abort_message
-
-c======================================================================
-
-      modname = 'phytrac_emiss'
+c======================================================================
+
 c EMISSION TRACEURS 
 
 c---------
 c debutphy
-c---------
+c---------      
       if (debutphy) then
-         print*,"DEBUT PHYTRAC"
-         print*,"PHYTRAC: EMISSION"
-
-c=============================================================
-c=============================================================
-c=============================================================
-c   Initialisation des traceurs
-c=============================================================
-c=============================================================
-c=============================================================
-
-C=========================================================================
+
+        print*,"DEBUT PHYTRAC"
+        print*,"PHYTRAC: EMISSION"
+
+        ALLOCATE(M_tr(nqtot))
+        M_tr(:)=64.                 ! SO2
+        
 C=========================================================================
 c Caracteristiques des traceurs emis:
+C=========================================================================
 
 c nombre total de traceur
-         if (nbsrc*nblat*nblon .gt. nqtot) then
+         if (nbz*nblat*nblon .gt. nqtot) then
+            print*, nbz*nblat*nblon, nqtot
             write(*,*) "Attention, pas assez de traceurs"
             write(*,*) "le dernier sera bien le dernier"
          endif
 
-c source en kg/s
-         source_volcan(1) = 1.
-         source_volcan(2) = 1000.
-c duration in Ed
-         Nemiss(1) = 1
-         Nemiss(2) = 10
+c quantite en kg
+         so2_quantity = 20.*10.**9.
+
+c height (in layer index)
+         Nheight(1) =  6  ! ~ 1 km
+         Nheight(2) = 16  ! ~ 25 km
+         Nheight(3) = 24  ! ~ 50 km
+
 c localisation volcan
          lat_volcan(1) =  70.
@@ -140 +120 @@
          lat_volcan(4) = -35.
          lat_volcan(5) = -70.
-         lon_volcan(1) = -120.
-         lon_volcan(2) =  -30.
-         lon_volcan(3) =   60.
-         lon_volcan(4) =  150.
-
+         lon_volcan(1) = -125.
+         lon_volcan(2) =  -35.
+         lon_volcan(3) =   55.
+         lon_volcan(4) =  145.
+         
+         ig_volcan(ilat,ilon)= 0
          deltalat = 180./jjm
-         deltalon = 360./jjm
+         deltalon = 360./iim
+
          do i=1,nlon
           do ilat=1,nblat
@@ -156 +138 @@
              ig_volcan(ilat,ilon)= i
              area_emiss(ilat,ilon) = airephy(i)
-            endif
-           enddo
-          enddo
-         enddo
-
-C=========================================================================
-C=========================================================================
-
+             print*,"Lat,lon=",ilat,ilon," OK"
+            end if
+           end do
+          end do
+         end do
+
+c Reinit des traceurs si necessaire
+         if (reinit_trac) then
+           tr_seri(:,:,:)=0.
+         endif
+         
+C=========================================================================
+C=========================================================================
+      ENDIF  ! fin debutphy 
 c-------------
 c fin debutphy
 c-------------
-      ENDIF  ! fin debutphy 
 
 c======================================================================
 c Emission d'un traceur pendant un certain temps
-c======================================================================
-         do i = 1,nlon
-          do iemiss = 1,nbsrc
-           do ilat  = 1,nblat
-            do ilon  = 1,nblon
-             it=(iemiss-1)*nblat*nblon+(ilat-1)*nblon+ilon
-             it=min(it,nqtot)
-             deltatr(i,1,it) = 0.
-
-             if (i .eq. ig_volcan(ilat,ilon)) then
+c necessite raz_date=1 dans run.def
+c et reinit_trac=y
+c======================================================================
+       deltatr(:,:,:) = 0.
 
 c source appliquee pendant Nemiss Ed
-               if (timesimu .lt. 86400.*Nemiss(iemiss)) then
-
+       if (timesimu .lt. 86400*Nemiss) then
+
+c emet les traceurs qui sont presents sur la grille
+        do ilat  = 1,nblat
+        do ilon  = 1,nblon
+         if (ig_volcan(ilat,ilon).ne.0) then
+         
+          do iz = 1,nbz
+           it=min( (iz-1)*nblat*nblon+(ilat-1)*nblon+ilon , nqtot )         
+           i=ig_volcan(ilat,ilon)
+           
+c injection dans une seule cellule:
 c source en kg/kg/s
-           deltatr(i,1,it) = source_volcan(iemiss)*RG
-     $     /(area_emiss(ilat,ilon)*(paprs(i,1) - paprs(i,2)))
-           tr_seri(i,1,it) = tr_seri(i,1,it) + deltatr(i,1,it)*pdtphys
-
-               end if  ! duree emission
-             end if ! i localisation
-            end do
-           end do
+c            deltatr(i,Nheight(iz),it) = so2_quantity/(86400.*Nemiss) ! kg/s
+c     $ *RG/( area_emiss(ilat,ilon)
+c     $      *(paprs(i,Nheight(iz))-paprs(i,Nheight(iz)+1)) )    ! /kg (masse cellule)
+     
+c            tr_seri(i,Nheight(iz),it) = tr_seri(i,Nheight(iz),it)
+c     $      + deltatr(i,Nheight(iz),it)*pdtphys
+
+c injection dans toute la colonne (a faire):
+            do k=1,Nheight(iz)
+             deltatr(i,k,it) = so2_quantity/(86400.*Nemiss) ! kg/s
+     $  *RG/( area_emiss(ilat,ilon)
+     $       *(paprs(i,1)-paprs(i,Nheight(iz)+1)) )    ! /kg (masse colonne)
+     
+             tr_seri(i,k,it) = tr_seri(i,k,it)+deltatr(i,k,it)*pdtphys
+            end do
+            
           end do
-         end do
+          
+         endif  ! ig_volcan!=0
+        end do
+        end do
+
+       end if  ! duree emission
+       
 c======================================================================
 c======================================================================

trunk/LMDZ.VENUS/libf/phyvenus/phytrac_relax.F

@@ -30 +30 @@
       use dimphy
       USE comgeomphy
+      USE chemparam_mod,only:M_tr
       IMPLICIT none
 #include "YOMCST.h"
@@ -95 +96 @@
       nqCO_OCS   = 6
 
+c !!!!!! ATTENTION: A REVOIR / A VERIFIER
+c   les traceurs sont en mass mixing ratio dans la dyn
+c   et convertis en frac mol avec M_tr dans les sorties...
+c   Mettre le profil de rappel en mass mixing ratio !!
+
+      print*,"METTRE A JOUR phytrac_relax"
+      stop
+
 c---------
 c debutphy
@@ -103 +112 @@
          allocate(zprof(nlev,nqtot),tau(nlev,nqtot))
 
+         ALLOCATE(M_tr(nqtot))
+      
 c=============================================================
 c=============================================================
@@ -144 +155 @@
       do it=1,3 
        print*,"INIT ZPROF ",tname(it)
+       M_tr(it)=28.       ! CO
        do k=1,nlev
          zprof(k,it)=0.
@@ -183 +195 @@
 c OCS
        print*,"INIT ZPROF ",tname(it+3)
+       M_tr(it+3)=60.       ! OCS
        do k=1,nlev
          zprof(k,it+3)=0.