Changeset 1926 for trunk/LMDZ.TITAN/libf

Timestamp:

May 28, 2018, 11:47:03 PM (7 years ago)

Author:

jvatant

Message:

1) Microphysics diags / outputs :

---

+ Add supplementary diagnostics outputs for microphysics ( precip, flux, rc ... ) ( new muphy_diag.F90 module )
+ Correct the outputs of microphys tracers to be in X/m-3 to be comparable to "standard values"

+ Also update the deftank callphys.def with latest revs modifs for microphysics

2) Condensation / chemistry updates :

---

+ Moved chemistry AFTER microphysics

• To have mufi condensation before photochem
• Chemistry called last coherent with the fact that it brings back fields to equilibrium

+ If 2D chemistry, make zonally averaged fields go through mufi and chem condensation

to have non saturated profiles in input of photochemistry
( other 'short' processes neglected as 2D -> no diurnal cycle, just seasonal evolution )

+ Also corrected the positivity check ( took Mars GCM syntax ) after chemistry ( could previously lead to negs )

3) Noticed a weird behaviour ( bug? ) :

---

+ If generalize the use of arrays *_indx for tracers, to get rid of ugly "iq+nmicro",

it ends up with weird results / crash in optim mode ( ok in debug ) but didn't find out why ...

--JVO

Location:

trunk/LMDZ.TITAN/libf

Files:

• muphytitan/mm_microphysic.f90 (modified) (1 diff)
• muphytitan/mmp_gcm.f90 (modified) (2 diffs)
• muphytitan/mmp_globals.f90 (modified) (1 diff)
• muphytitan/mmp_optics.f90 (modified) (1 diff)
• phytitan/calmufi.F90 (modified) (2 diffs)
• phytitan/inimufi.F90 (modified) (3 diffs)
• phytitan/muphy_diag.F90 (added)
• phytitan/physiq_mod.F90 (modified) (14 diffs)
• phytitan/tracer_h.F90 (modified) (2 diffs)

trunk/LMDZ.TITAN/libf/muphytitan/mm_microphysic.f90 ¶

@@ -202 +202 @@
     !! microphysics is disabled (see [[mm_globals(module):mm_w_clouds(variable)]] documentation).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: aer_prec   !! Aerosols precipitations (both modes) (m).
-    REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: ccn_prec   !! CCN precipitations (m).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: aer_s_flux !! Spherical aerosol mass flux (\(kg.m^{-2}.s^{-1}\)).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: aer_f_flux !! Fractal aerosol mass flux (\(kg.m^{-2}.s^{-1}\)).
+    REAL(kind=mm_wp), INTENT(out), OPTIONAL                 :: ccn_prec   !! CCN precipitations (m).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: ccn_flux   !! CCN mass flux (\(kg.m^{-2}.s^{-1}\)).
+    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: ice_prec   !! Ice precipitations (m).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: ice_fluxes !! Ice sedimentation fluxes (\(kg.m^{-2}.s^{-1}\)).
     REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:,:) :: gazs_sat   !! Condensible gaz saturation ratios (--).
-    REAL(kind=mm_wp), INTENT(out), OPTIONAL, DIMENSION(:)   :: ice_prec   !! Ice precipitations (m).
 
     IF (PRESENT(aer_prec))   aer_prec   = ABS(mm_aer_prec) 

trunk/LMDZ.TITAN/libf/muphytitan/mmp_gcm.f90 ¶

@@ -93 +93 @@
     INTEGER                                           :: i
     TYPE(cfgparser)                                   :: cparser
-    CHARACTER(len=st_slen)                            :: spcpath,pssfile,mqfile
+    CHARACTER(len=st_slen)                            :: spcpath,pssfile,mqfile,opt_file
     CHARACTER(len=st_slen), DIMENSION(:), ALLOCATABLE :: species
     REAL(kind=mm_wp), DIMENSION(:), ALLOCATABLE       :: tmp
@@ -138 +138 @@
     err = mm_check_opt(cfg_get_value(cparser,"haze_coag_interactions",coag_choice),coag_choice,7,mm_log)
 
+    ! optic look-up table file path.
+    mmp_optic_file = ''
+    opt_file = ''
+    err = mm_check_opt(cfg_get_value(cparser,"optics/optic_file",opt_file),opt_file,'',mm_log)
+    IF (err /= 0) THEN
+      WRITE(*,'(a)') "Warning: I was unable to retrieve the path of the optic look-up table file:"
+      WRITE(*,'(a)') "  The GCM may abort if it uses YAMMS optical properties calculation module !"
+    ELSE
+      mmp_optic_file = TRIM(opt_file)
+    ENDIF
+
     ! Retrieve clouds species configuration file
     spcpath = ''

trunk/LMDZ.TITAN/libf/muphytitan/mmp_globals.f90 ¶

@@ -92 +92 @@
   !> Aerosol electric charge correction control flag.
   LOGICAL, SAVE :: mmp_w_qe   = .true.
-
+  !> Optic look-up table file path.
+  CHARACTER(len=:), ALLOCATABLE, SAVE :: mmp_optic_file
 
   CONTAINS 

trunk/LMDZ.TITAN/libf/muphytitan/mmp_optics.f90 ¶

@@ -62 +62 @@
 
   PRIVATE
+  PUBLIC :: mmp_optic_file ! from mmp_globals :)
   PUBLIC :: mmp_initialize_optics
   PUBLIC :: mmp_sph_optics_vis,mmp_sph_optics_ir

trunk/LMDZ.TITAN/libf/phytitan/calmufi.F90 ¶

@@ -20 +20 @@
   USE comcstfi_mod, only : g
   USE callkeys_mod, only : callclouds
+  USE muphy_diag
   IMPLICIT NONE
 
@@ -141 +142 @@
         call abort_program(error("mm_muphys aborted -> initialization not done !",-1))
     ENDIF
+    ! save diags (if no clouds, relevant arrays will be set to 0 !)
+    call mm_diagnostics(mmd_aer_prec(ilon),mmd_aer_s_flux(ilon,:),mmd_aer_f_flux(ilon,:),  &
+                        mmd_ccn_prec(ilon),mmd_ccn_flux(ilon,:), mmd_ice_prec(ilon,:),   &
+                        mmd_ice_fluxes(ilon,:,:),mmd_gazs_sat(ilon,:,:))
+    call mm_get_radii(mmd_rc_sph(ilon,:),mmd_rc_fra(ilon,:),mmd_rc_cld(ilon,:))
 
     ! Convert tracers back to intensives ( except for gazs where we work with molar mass ratio )

trunk/LMDZ.TITAN/libf/phytitan/inimufi.F90 ¶

@@ -1 @@
-subroutine inimufi(nq,ptimestep)
+subroutine inimufi(ptimestep)
 
   use mmp_gcm
@@ -31 +31 @@
   real, intent(in)                   :: ptimestep    ! Timestep (s)
 
-  integer, intent(in)                :: nq           ! Total number of tracers 
-
   integer :: i,idx
   character(len=20), dimension(4), parameter :: aernames = &
@@ -49 +47 @@
   ! ---------------------------------------------------
 
-   call mmp_initialize(ptimestep,p_prod,tx_prod,rc_prod, &
+  ! enable log for what it's worth...
+  ! mm_log = .true.
+  
+  call mmp_initialize(ptimestep,p_prod,tx_prod,rc_prod, &
         rad,g,air_rad,mugaz,callclouds,config_mufi)
 

trunk/LMDZ.TITAN/libf/phytitan/physiq_mod.F90 ¶

@@ -1 @@
-#define USE_QTEST
-
       module physiq_mod
       
@@ -19 +17 @@
       use radcommon_h, only: sigma, glat, grav, BWNV
       use surfdat_h, only: phisfi, zmea, zstd, zsig, zgam, zthe
-      use comchem_h, only: nkim
+      use comchem_h, only: nkim, cnames
       use comdiurn_h, only: coslat, sinlat, coslon, sinlon
       use comsaison_h, only: mu0, fract, dist_star, declin, right_ascen
@@ -50 +48 @@
 #endif
       use MMP_OPTICS
+      use muphy_diag
       implicit none
 
@@ -80 +79 @@
 !
 !      V. Tracers
-!         V.1. Chemistry
-!         V.2. Microphysics
+!         V.1. Microphysics
+!         V.2. Chemistry
 !         V.3. Updates (pressure variations, surface budget).
 !         V.4. Surface Tracer Update.
@@ -291 +290 @@
       
       real zdqmufi(ngrid,nlayer,nq)   ! Microphysical tendency.
+      
+      real zdqfibar(ngrid,nlayer,nq)   ! For 2D chemistry
+      real zdqmufibar(ngrid,nlayer,nq) ! For 2D chemistry
                         
       ! For Winds : (m/s/s)
@@ -372 +374 @@
       ! Chemical tracers in molar fraction 
       real, dimension(ngrid,nlayer,nkim)          :: ychim ! (mol/mol)
+      real, dimension(ngrid,nlayer,nkim)          :: ychimbar ! For 2D chemistry
 
       ! Molar fraction tendencies ( chemistry and condensation ) for tracers (mol/mol/s)
       real, dimension(ngrid,nlayer,nq)            :: dyccond ! All tracers, we want to use indx on it.
+      real, dimension(ngrid,nlayer,nq)            :: dyccondbar ! For 2D chemistry
       real, dimension(:,:,:), allocatable, save   :: dycchi  ! Only for chem tracers. Saved since chemistry is not called every step.
 !$OMP THREADPRIVATE(dycchi)
@@ -393 +397 @@
       ! -----******----- END FOR MUPHYS OPTICS -----******----- 
 
+      real,dimension(:,:), allocatable :: i2e2d ! factor to convert X.kg-1 in X.m-3 (for microphysics diags) 
+      
       real,save,dimension(:,:,:), allocatable :: tpq ! Tracers for decoupled microphysical tests ( temporary in 01/18 )
 !$OMP THREADPRIVATE(tpq)
@@ -432 +438 @@
         ! Initialisation of nmicro as well as tracers names, indexes ...
         if (ngrid.ne.1) then ! Already done in rcm1d
-           call initracer2(nq,nametrac)
+           call initracer2(nq,nametrac) ! WARNING JB (27/03/2018): should be wrapped in an OMP SINGLE statement (see module notes)
         endif
 
@@ -539 +545 @@
 
          IF ( callmufi ) THEN
-
-           call inimufi(nq,ptimestep)
-           
+           ! WARNING JB (27/03/2018): inimufi AND mmp_initialize_optics should be wrapped in an OMP SINGLE statement.
+           call inimufi(ptimestep)
            ! Optical coupling of YAMMS is plugged but inactivated for now
            ! as long as the microphysics only isn't fully debugged -- JVO 01/18
-           IF (.NOT.uncoupl_optic_haze) call mmp_initialize_optics("/path/to/mmp_optic_table.nc")
+           ! NOTE JB (03/18): mmp_optic_file is initialized in inimufi => mmp_initialize_optics must be called AFTER inimufi
+           IF (.NOT.uncoupl_optic_haze) call mmp_initialize_optics(mmp_optic_file)
+
+           ! initialize microphysics diagnostics arrays.
+           call ini_diag_arrays(ngrid,nlayer,nice)
 
          ENDIF
@@ -1068 +1077 @@
       if (tracer) then
 
-  ! -------------------------
-  !   V.1. Chemistry
-  ! -------------------------
-
-         if (callchim) then
-
-            ! o. Convert updated tracers to molar fraction
-            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-            do iq = 1,nkim
-               ychim(:,:,iq) = ( pq(:,:,iq+nmicro) + pdq(:,:,iq+nmicro) ) / rat_mmol(iq+nmicro)
-            enddo
-
-            ! i. Condensation after the transport
-            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-            do iq=1,nkim
-               do l=1,nlayer
-                  do ig=1,ngrid
-                     if ( ychim(ig,l,iq).gt.qysat(l,iq) ) then
-                        dyccond(ig,l,iq+nmicro) = ( -ychim(ig,l,iq)+qysat(l,iq) ) / ptimestep
-                     else
-                        dyccond(ig,l,iq+nmicro) = 0.0 ! since array not saved
-                     endif
-                  enddo
-               enddo
-            enddo
-
-            if ( callclouds ) dyccond(:,:,ices_indx) = 0.0 ! Condensation will be calculated in the cloud microphysics 
-
-            do iq=1,nkim
-              ychim(:,:,iq) = ychim(:,:,iq) + dyccond(:,:,iq+nmicro) ! update molar ychim for following calchim
-              zdqcond(:,:,iq+nmicro) = dyccond(:,:,iq+nmicro)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio
-            enddo
-            
-            pdq(:,:,:) = pdq(:,:,:) + zdqcond(:,:,:)
-
-            ! ii. Photochemistry
-            ! ~~~~~~~~~~~~~~~~~~
-            if( mod(icount-1,ichim).eq.0. ) then
-               
-               print *, "We enter in the chemistry ..."
-
-               if (moyzon_ch) then ! 2D zonally averaged chemistry
-
-                 do iq = 1,nkim ! In this case we send zonal average from dynamics to chem. module
-                   ychim(:,:,iq) = zqfibar(:,:,iq+nmicro) / rat_mmol(iq+nmicro)
-                 enddo
-
-                 call calchim(ngrid,ychim,declin,ctimestep,ztfibar,zphibar,  &
-                              zplaybar,zplevbar,zzlaybar,zzlevbar,dycchi)
-
-               else ! 3D chemistry (or 1D run)
-                 call calchim(ngrid,ychim,declin,ctimestep,pt,pphi,  &
-                              pplay,pplev,zzlay,zzlev,dycchi)
-               endif ! if moyzon
-
-            endif
-            
-            do iq=1,nkim
-               zdqchi(:,:,iq+nmicro) = dycchi(:,:,iq)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio
-  
-               where( (pq(:,:,iq+nmicro)+zdqchi(:,:,iq+nmicro) ).LT.1e-40)  & ! When using zonal means we set the same tendency
-                       zdqchi(:,:,iq+nmicro) = 1e-40 - ( pq(:,:,iq+nmicro) )  ! everywhere in longitude -> can lead to negs !
-             enddo
-
-            pdq(:,:,:) = pdq(:,:,:) + zdqchi(:,:,:)
-            
-         endif ! end of 'callchim'
-
   ! -------------------
-  !   V.2 Microphysics
+  !   V.1 Microphysics
   ! -------------------
+
+         ! JVO 05/18 : We must call microphysics before chemistry, for condensation !
  
          if (callmufi) then
@@ -1149 +1092 @@
             pdq(:,:,:) = pdq(:,:,:) + zdqmufi(:,:,:)
 #endif
-         endif
-      
+
+            ! Microphysics condensation for 2D fields to sent non-saturated fields to photochem
+            if ( callclouds .and. moyzon_ch .and. mod(icount-1,ichim).eq.0 ) then
+              zdqfibar(:,:,:) = 0.0 ! We work in zonal average -> forget processes other than condensation
+              call calmufi(ptimestep,zplevbar,zzlevbar,zplaybar,zzlaybar, &
+                           ztfibar,zqfibar,zdqfibar,zdqmufibar)
+            endif
+
+         endif
+      
+  ! -----------------
+  !   V.2. Chemistry
+  ! -----------------
+  !   NB : Must be call last ( brings fields back to an equilibrium )
+
+  ! Known bug ? ( JVO 18 ) : If you try to use chimi_indx instead of iq+nmicro
+  ! it leads to weird results / crash on dev mod ( ok in debug ) ... Why ? Idk ...
+
+         if (callchim) then
+
+            ! o. Convert updated tracers to molar fraction
+            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+            do iq = 1,nkim
+               ychim(:,:,iq) = ( pq(:,:,iq+nmicro) + pdq(:,:,iq+nmicro) ) / rat_mmol(iq+nmicro)
+            enddo
+
+            ! JVO 05/18 : We update zonally averaged fields with condensation 
+            ! as it is compulsory to have correct photochem production. But for other 
+            ! processes ( convadj ... ) we miss them in any case as we work in zonally/diurnal 
+            ! mean -> no fine diurnal/short time evolution, only seasonal evolution only.
+            if ( moyzon_ch .and. mod(icount-1,ichim).eq. 0 ) then
+              do iq = 1,nkim
+                ychimbar(:,:,iq) =  zqfibar(:,:,iq+nmicro) / rat_mmol(iq+nmicro)
+                  if ( callclouds ) then
+                    ychimbar(:,:,iq) =  ychimbar(:,:,iq) + ( zdqmufibar(:,:,iq+nmicro) / rat_mmol(iq+nmicro) )
+                  endif
+              enddo
+            endif
+
+            ! i. Condensation after the transport
+            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+            do iq=1,nkim
+               do l=1,nlayer
+                  do ig=1,ngrid
+                     if ( ychim(ig,l,iq).gt.qysat(l,iq) ) then
+                        dyccond(ig,l,iq+nmicro) = ( -ychim(ig,l,iq)+qysat(l,iq) ) / ptimestep
+                     else
+                        dyccond(ig,l,iq+nmicro) = 0.0
+                     endif
+                  enddo
+               enddo
+            enddo
+
+            if ( callclouds ) dyccond(:,:,ices_indx) = 0.0 ! Condensation have been calculated in the cloud microphysics 
+
+            do iq=1,nkim
+              ychim(:,:,iq) = ychim(:,:,iq) + dyccond(:,:,iq+nmicro) ! update molar ychim for following calchim
+              zdqcond(:,:,iq+nmicro) = dyccond(:,:,iq+nmicro)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio
+            enddo
+            
+            pdq(:,:,:) = pdq(:,:,:) + zdqcond(:,:,:)
+
+            ! 2D zonally averaged fields needed to condense before photochem
+            if ( moyzon_ch .and. mod(icount-1,ichim).eq. 0 ) then
+              do iq = 1,nkim
+                do l=1,nlayer
+                  do ig=1,ngrid
+                    if ( ychimbar(ig,l,iq).gt.qysat(l,iq) ) then
+                       dyccondbar(ig,l,iq+nmicro) = ( -ychimbar(ig,l,iq)+qysat(l,iq) ) / ptimestep
+                    else
+                       dyccondbar(ig,l,iq+nmicro) = 0.0
+                    endif
+                  enddo
+                enddo
+              enddo
+
+              if ( callclouds ) dyccondbar(:,:,ices_indx) = 0.0 ! Condensation have been calculated in the cloud microphysics 
+
+              do iq=1,nkim
+                ychimbar(:,:,iq) = ychimbar(:,:,iq) + dyccondbar(:,:,iq+nmicro)
+              enddo
+
+            endif ! if ( moyzon_ch .and. mod(icount-1,ichim).eq. 0 )
+
+            ! ii. Photochemistry ( must be call after condensation)
+            ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+            if( mod(icount-1,ichim).eq.0. ) then
+               
+               print *, "We enter in the chemistry ..."
+
+               if (moyzon_ch) then ! 2D zonally averaged chemistry
+
+                 ! Here we send zonal average fields ( corrected with cond ) from dynamics to chem. module
+                 call calchim(ngrid,ychimbar,declin,ctimestep,ztfibar,zphibar,  &
+                              zplaybar,zplevbar,zzlaybar,zzlevbar,dycchi)
+               else ! 3D chemistry (or 1D run)
+                 call calchim(ngrid,ychim,declin,ctimestep,pt,pphi,  &
+                              pplay,pplev,zzlay,zzlev,dycchi)
+               endif ! if moyzon
+
+            endif
+            
+            do iq=1,nkim
+              zdqchi(:,:,iq+nmicro) = dycchi(:,:,iq)*rat_mmol(iq+nmicro) ! convert tendencies to mass mixing ratio
+
+              where( (pq(:,:,iq+nmicro) + ( pdq(:,:,iq+nmicro)+zdqchi(:,:,iq+nmicro) )*ptimestep ) .LT. 0.)  & ! When using zonal means we set the same tendency
+                      zdqchi(:,:,iq+nmicro) = 1.e-30 - pdq(:,:,iq+nmicro) - pq(:,:,iq+nmicro)/ptimestep        ! everywhere in longitude -> could lead to negs !
+            enddo
+
+            pdq(:,:,:) = pdq(:,:,:) + zdqchi(:,:,:)
+            
+         endif ! end of 'callchim'
+
   ! ---------------
   !   V.3 Updates
@@ -1224 +1278 @@
 
       ! Note : For output only: the actual model integration is performed in the dynamics.
-
-
  
       ! Temperature, zonal and meridional winds.
@@ -1457 +1509 @@
       if (tracer) then
 
-         if (callmufi) then ! For now we assume an given order for tracers !
+         if (callmufi) then
+           ! Microphysical tracers are expressed in unit/m3.
+           ! convert X.kg-1 --> X.m-3
+           ALLOCATE(i2e2d(ngrid,nlayer))
+           i2e2d(:,:) = ( pplev(:,1:nlayer)-pplev(:,2:nlayer+1) ) / g /(zzlev(:,2:nlayer+1)-zzlev(:,1:nlayer))
+
 #ifdef USE_QTEST
             ! Microphysical tracers passed through dyn+phys(except mufi)
-            call writediagfi(ngrid,"mu_m0as_dp","Dynphys only spherical mode 0th order moment",'kg/kg',3,zq(:,:,1))
-            call writediagfi(ngrid,"mu_m3as_dp","Dynphys only spherical mode 3rd order moment",'kg/kg',3,zq(:,:,2))
-            call writediagfi(ngrid,"mu_m0af_dp","Dynphys only fractal mode 0th order moment",'kg/kg',3,zq(:,:,3))
-            call writediagfi(ngrid,"mu_m3af_dp","Dynphys only fractal mode 3rd order moment",'kg/kg',3,zq(:,:,4))
+            call writediagfi(ngrid,"mu_m0as_dp","Dynphys only spherical mode 0th order moment",'m-3',3,zq(:,:,micro_indx(1))*i2e2d)
+            call writediagfi(ngrid,"mu_m3as_dp","Dynphys only spherical mode 3rd order moment",'m3/m3',3,zq(:,:,micro_indx(2))*i2e2d)
+            call writediagfi(ngrid,"mu_m0af_dp","Dynphys only fractal mode 0th order moment",'m-3',3,zq(:,:,micro_indx(3))*i2e2d)
+            call writediagfi(ngrid,"mu_m3af_dp","Dynphys only fractal mode 3rd order moment",'m3/m3',3,zq(:,:,micro_indx(4))*i2e2d)
             ! Microphysical tracers passed through mufi only
-            call writediagfi(ngrid,"mu_m0as_mo","Mufi only spherical mode 0th order moment",'kg/kg',3,tpq(:,:,1))
-            call writediagfi(ngrid,"mu_m3as_mo","Mufi only spherical mode 3rd order moment",'kg/kg',3,tpq(:,:,2))
-            call writediagfi(ngrid,"mu_m0af_mo","Mufi only fractal mode 0th order moment",'kg/kg',3,tpq(:,:,3))
-            call writediagfi(ngrid,"mu_m3af_mo","Mufi only fractal mode 3rd order moment",'kg/kg',3,tpq(:,:,4))
+            call writediagfi(ngrid,"mu_m0as_mo","Mufi only spherical mode 0th order moment",'m-3',3,tpq(:,:,micro_indx(1))*i2e2d)
+            call writediagfi(ngrid,"mu_m3as_mo","Mufi only spherical mode 3rd order moment",'m3/m3',3,tpq(:,:,micro_indx(2))*i2e2d)
+            call writediagfi(ngrid,"mu_m0af_mo","Mufi only fractal mode 0th order moment",'m-3',3,tpq(:,:,micro_indx(3))*i2e2d)
+            call writediagfi(ngrid,"mu_m3af_mo","Mufi only fractal mode 3rd order moment",'m3/m3',3,tpq(:,:,micro_indx(4))*i2e2d)
 #else
-            call writediagfi(ngrid,"mu_m0as","Spherical mode 0th order moment",'kg/kg',3,zq(:,:,1))
-            call writediagfi(ngrid,"mu_m3as","Spherical mode 3rd order moment",'kg/kg',3,zq(:,:,2))
-            call writediagfi(ngrid,"mu_m0af","Fractal mode 0th order moment",'kg/kg',3,zq(:,:,3))
-            call writediagfi(ngrid,"mu_m3af","Fractal mode 3rd order moment",'kg/kg',3,zq(:,:,4))
+            call writediagfi(ngrid,"mu_m0as","Spherical mode 0th order moment",'m-3',3,zq(:,:,micro_indx(1))*i2e2d)
+            call writediagfi(ngrid,"mu_m3as","Spherical mode 3rd order moment",'m3/m3',3,zq(:,:,micro_indx(2))*i2e2d)
+            call writediagfi(ngrid,"mu_m0af","Fractal mode 0th order moment",'m-3',3,zq(:,:,micro_indx(3))*i2e2d)
+            call writediagfi(ngrid,"mu_m3af","Fractal mode 3rd order moment",'m3/m3',3,zq(:,:,micro_indx(4))*i2e2d)
 #endif
+            
+            DEALLOCATE(i2e2d)
+
+            ! Microphysical diagnostics
+            call writediagfi(ngrid,"mmd_aer_prec","Total aerosols precipitations",'m',2,mmd_aer_prec)
+            call writediagfi(ngrid,"mmd_aer_s_flux","Spherical aerosols sedimentation flux",'kg.m-2.s-1',3,mmd_aer_s_flux)
+            call writediagfi(ngrid,"mmd_aer_f_flux","Fractal aerosols sedimentation flux",'kg.m-2.s-1',3,mmd_aer_f_flux)
+            call writediagfi(ngrid,"mmd_rc_sph","Spherical mode caracteristic radius",'m',3,mmd_rc_sph)
+            call writediagfi(ngrid,"mmd_rc_fra","Fractal mode caracteristic radius",'m',3,mmd_rc_fra)
+
          endif ! end of 'callmufi'
 
-        if ( callchim .and. (.not. callmufi) ) then
-           call writediagfi(ngrid,"C2H2","C2H2",'kg/kg',3,zq(:,:,10))
-           call writediagfi(ngrid,"C2H4","C2H4",'kg/kg',3,zq(:,:,12))
-           call writediagfi(ngrid,"C2H6","C2H6",'kg/kg',3,zq(:,:,14))
-           call writediagfi(ngrid,"C4H2","C4H2",'kg/kg',3,zq(:,:,26))
-           call writediagfi(ngrid,"HCN","HCN",'kg/kg',3,zq(:,:,36))
-           call writediagfi(ngrid,"HC3N","HC3N",'kg/kg',3,zq(:,:,42))
-        else if ( callchim .and. callmufi ) then
-           call writediagfi(ngrid,"C2H2","C2H2",'kg/kg',3,zq(:,:,14))
-           call writediagfi(ngrid,"C2H4","C2H4",'kg/kg',3,zq(:,:,16))
-           call writediagfi(ngrid,"C2H6","C2H6",'kg/kg',3,zq(:,:,18))
-           call writediagfi(ngrid,"C4H2","C4H2",'kg/kg',3,zq(:,:,30))
-           call writediagfi(ngrid,"HCN","HCN",'kg/kg',3,zq(:,:,40))
-           call writediagfi(ngrid,"HC3N","HC3N",'kg/kg',3,zq(:,:,46))
-        endif
+         ! Chemical tracers
+         if (callchim) then
+           do iq=1,nkim
+             call writediagfi(ngrid,cnames(iq),cnames(iq),'mol/mol',3,zq(:,:,iq+nmicro)/rat_mmol(iq+nmicro))
+           enddo
+         endif
 
        endif ! end of 'tracer'
@@ -1532 +1590 @@
       CALL send_xios_field("dtrad",dtrad)
       CALL send_xios_field("dtdyn",zdtdyn)
-        
+
+      ! Chemical tracers
+      if (callchim) then
+        do iq=1,nkim
+          CALL send_xios_field(cnames(iq),zq(:,:,iq+nmicro)/rat_mmol(iq+nmicro))
+        enddo
+      endif
+
+      ! Microphysical tracers
+      if (callmufi) then
+        CALL send_xios_field("mu_m0as",zq(:,:,micro_indx(1))*i2e2d)
+        CALL send_xios_field("mu_m3as",zq(:,:,micro_indx(2))*i2e2d)
+        CALL send_xios_field("mu_m0af",zq(:,:,micro_indx(3))*i2e2d)
+        CALL send_xios_field("mu_m3af",zq(:,:,micro_indx(4))*i2e2d)
+      endif        
+
       if (lastcall.and.is_omp_master) then
         write(*,*) "physiq: call xios_context_finalize"

trunk/LMDZ.TITAN/libf/phytitan/tracer_h.F90 ¶

@@ -1 @@
-
+! WARNING (JB 27/03/2018)
+!
+! OpenMP directives in this module are inconsistent:
+!   - sizes (nmicro...) are thread private \___ BUT BOTH ARE INITIALIZED IN THE SAME ROUTINE 
+!   - indexes array are common             /
+!
+!  Tracers sizes do not need to be private.
+!  In such case, OMP THREADPRIVATE should be removed and initracer2 should be called within an OMP SINGLE statement.
+!
 MODULE tracer_h
   !! Stores data related to physics tracers.
@@ -51 +59 @@
     !! dynamics counterpart.
     !!
-    !! In addition, it initializes arrays of indexes for the different sub-processs of the physics:
+    !! In addition, it initializes arrays of indexes for the different sub-processes of the physics:
     !!
     !!   - tracers_h:micro_indxs, the array of tracers indexes used for the microphysics.