Changeset 2126 for trunk/LMDZ.COMMON/libf/dyn3dpar

Timestamp:

Apr 26, 2019, 11:18:52 AM (6 years ago)

Author:

emillour

Message:

Common dynamics:
Some work to enforce total tracer mass conservation in the dynamics.
Still to be further studied and validated.
For now these changes are triggered by setting a "force_conserv_tracer"
flag to ".true." in run.def (default is ".false." to not change anything
with respect to previous versions).
When force_conserv_tracer=.true. then:

1. Rescale tracer mass in caladvtrac after tracer advection computations
2. Recompute q ratios once atmospheric mass has been updated in integrd

These steps technically ensure total tracer mass conservation but it
might be the tracer advection scheme and/or time-stepping updating
sequence of fields that should be rethought or fixed.
EM

Location:

trunk/LMDZ.COMMON/libf/dyn3dpar

Files:

• advtrac_p.F90 (modified) (3 diffs)
• conf_gcm.F90 (modified) (2 diffs)
• integrd_p.F (modified) (3 diffs)
• planetary_operations_p.F90 (added)

trunk/LMDZ.COMMON/libf/dyn3dpar/advtrac_p.F90

@@ -18 +18 @@
   USE times
   USE infotrac, ONLY: nqtot, iadv
-  USE control_mod, ONLY: iapp_tracvl, day_step, planet_type
+  USE control_mod, ONLY: iapp_tracvl, day_step, planet_type, &
+                         force_conserv_tracer
   USE comconst_mod, ONLY: dtvr
+  USE planetary_operations_p, ONLY: planetary_tracer_amount_from_mass_p
   IMPLICIT NONE
   !
@@ -77 +79 @@
   REAL,SAVE :: teta_tmp(ip1jmp1,llm)
   REAL,SAVE :: pk_tmp(ip1jmp1,llm)
+  REAL :: totaltracer_old(nqtot),totaltracer_new(nqtot)
+  REAL :: ratio
+
+! Ehouarn : try to fix tracer conservation issues:
+  if (force_conserv_tracer) then
+    do iq=1,nqtot
+       call planetary_tracer_amount_from_mass_p(masse,q(:,:,iq), &
+                                          totaltracer_old(iq))
+    enddo
+  endif
 
   ijb_u=ij_begin
@@ -467 +479 @@
      endif ! of if (planet_type=="earth")
 
+     ! Ehouarn : try to fix tracer conservation after tracer advection
+     if (force_conserv_tracer) then
+       ijb=ij_begin
+       ije=ij_end
+       do iq=1,nqtot
+         call planetary_tracer_amount_from_mass_p(masse,q(:,:,iq), &
+                                     totaltracer_new(iq))
+         ratio=totaltracer_old(iq)/totaltracer_new(iq)
+         q(ijb:ije,1:llm,iq)=q(ijb:ije,1:llm,iq)*ratio
+       enddo
+     endif !of if (force_conserv_tracer)
+
         !------------------------------------------------------------------
         !   on reinitialise a zero les flux de masse cumules

trunk/LMDZ.COMMON/libf/dyn3dpar/conf_gcm.F90

@@ -255 +255 @@
   iapp_tracvl = iperiod
   CALL getin('iapp_tracvl',iapp_tracvl)
+
+! Enforce tracer conservation in integrd & caladvtrac ?
+  force_conserv_tracer = .false. 
+  CALL getin('force_conserv_tracer',force_conserv_tracer)
 
 !Config  Key  = iconser
@@ -987 +991 @@
       write(lunout,*)' iphysiq = ', iphysiq
       write(lunout,*)' iflag_trac = ', iflag_trac
+      write(lunout,*)' iapp_tracvl = ', iapp_tracvl
+      write(lunout,*)' force_conserv_tracer = ', force_conserv_tracer
       write(lunout,*)' clon = ', clon
       write(lunout,*)' clat = ', clat

trunk/LMDZ.COMMON/libf/dyn3dpar/integrd_p.F

@@ -7 +7 @@
       USE parallel_lmdz, ONLY: ij_begin, ij_end, pole_nord, pole_sud,
      &                         omp_chunk
-      USE control_mod, only : planet_type
+      USE control_mod, only : planet_type,force_conserv_tracer
       USE comvert_mod, ONLY: ap,bp
       USE comconst_mod, ONLY: pi
@@ -65 +65 @@
       REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
       REAL massescr( ip1jmp1,llm )
+      REAL :: massratio(ip1jmp1,llm)
 !      REAL finvmasse(ip1jmp1,llm)
       REAL,SAVE :: p(ip1jmp1,llmp1)
@@ -361 +362 @@
       endif ! of if (planet_type.eq."earth")
 
+
+      if (force_conserv_tracer) then
+        ! Ehouarn: try to keep total amont of tracers fixed 
+        ! by acounting for mass change in each cell
+        massratio(ijb:ije,1:llm)=massescr(ijb:ije,1:llm)
+     &                             /masse(ijb:ije,1:llm)
+        do iq=1,nq
+        q(ijb:ije,1:llm,iq)=q(ijb:ije,1:llm,iq)
+     &                        *massratio(ijb:ije,1:llm)
+        enddo
+      endif ! of if (force_conserv_tracer)
 c
 c