Changeset 2635 for trunk/LMDZ.GENERIC

Timestamp:

Mar 4, 2022, 6:29:05 PM (3 years ago)

Author:

aslmd

Message:

Introduction of cpp_mugaz_mode to decide what source should be used for cpp and mugaz. Force_cpp and check_cpp_match are deprecated (see issue la-communaut-des-mod-les-atmosph-riques-plan-taires/git-trunk#27).

Location:

trunk/LMDZ.GENERIC/libf/phystd

Files:

• calc_cpp_mugaz.F90 (modified) (4 diffs)
• callkeys_mod.F90 (modified) (1 diff)
• dyn1d/kcm1d.F90 (modified) (2 diffs)
• dyn1d/rcm1d.F (modified) (2 diffs)
• inifis_mod.F90 (modified) (4 diffs)
• tabfi_mod.F90 (modified) (2 diffs)

trunk/LMDZ.GENERIC/libf/phystd/calc_cpp_mugaz.F90

@@ -4 +4 @@
 !     Purpose
 !     -------
-!     Check to see if the atmospheric specific heat capacity and
+!     Computes atmospheric specific heat capacity and
 !     mean molar mass for the gas mixture defined in gases.def
-!     corresponds to what we're using. If it doesn't, abort run
-!     unless option 'check_cpp_match' is set to false in 
-!     callphys.def.
+!     with standard values.
+!     In 1d, one can force the model to adopt these standard values
+!     by setting cpp_mugaz_mode=2 in callfis.def.
 !
 !     Authors
@@ -14 +14 @@
 !     Robin Wordsworth (2009)
 !     A. Spiga: make the routine OK with latest changes in rcm1d
+!     Jeremy Leconte (2022)
 !
 !==================================================================
@@ -19 +20 @@
       use gases_h
       use comcstfi_mod, only: cpp, mugaz
-      use callkeys_mod, only: check_cpp_match,force_cpp
+      use callkeys_mod, only: cpp_mugaz_mode
+      use mod_phys_lmdz_para, only : is_master
       implicit none
 
@@ -122 +124 @@
       print*,'Cp in calc_cpp_mugaz is ',cpp_c,'J kg^-1 K^-1'
       print*,'Mg in calc_cpp_mugaz is ',mugaz_c,'amu'
-      print*,'Predefined Cp in physics is ',cpp,'J kg^-1 K^-1'
-      print*,'Predefined Mg in physics is ',mugaz,'amu'
 
-      if (check_cpp_match) then
-         print*,'REQUEST TO CHECK cpp_match :'
-         if((abs(1.-cpp/cpp_c).gt.1.e-6) .or.  &
-              (abs(1.-mugaz/mugaz_c).gt.1.e-6)) then
-            ! Ehouarn: tolerate a small mismatch between computed/stored values
-            print*,'--> Values do not match!'
-            print*,'    Either adjust cpp / mugaz via newstart to calculated values,'
-            print*,'    or set check_cpp_match to .false. in callphys.def.'
-            stop
-         else
-            print*,'--> OK. Settings match composition.'
-         endif
+      if(((abs(1.-cpp/cpp_c).gt.1.e-6) .or.  &
+           (abs(1.-mugaz/mugaz_c).gt.1.e-6)).and. is_master ) then
+         ! Ehouarn: tolerate a small mismatch between computed/stored values
+         print*,'--> Values do not match with the predefined one ! (',cpp,',',mugaz,')'
+         print*,'    Because cp varies with temperature and that some gases may not appear in gases.def,'
+         print*,'    a small discrepancy might be completely normal.'
+         print*,'    But you might want to check that!'
+         print*,'    If you want to use the values calculated here, adjust cpp / mugaz in the dynamics via newstart (3d)'
+         print*,'    or use cpp_mugaz_mode=2 (if you are in 1d).'
       endif
 
-      if (.not.force_cpp) then
-          print*,'*** Setting cpp & mugaz to computations in calc_cpp_mugaz.'
+      if (cpp_mugaz_mode==2) then
+          if (is_master) print*,'*** cpp_mugaz_mode==2, so setting cpp & mugaz to computations in calc_cpp_mugaz.'
           mugaz = mugaz_c
           cpp = cpp_c
       else
-          print*,'*** Setting cpp & mugaz to predefined values.'
+          if (is_master) print*,'*** Leaving cpp & mugaz equal to predefined values'
+          if (is_master) print*,'(either from dynamics (cpp_mugaz_mode=0) or callfis (cpp_mugaz_mode=1)).'
       endif
 

trunk/LMDZ.GENERIC/libf/phystd/callkeys_mod.F90

@@ -25 +25 @@
 !$OMP THREADPRIVATE(enertest,nonideal,meanOLR,kastprof,diagdtau)
       logical,save :: newtonian
-      logical,save :: check_cpp_match
       logical,save :: force_cpp
+      integer,save :: cpp_mugaz_mode
       logical,save :: testradtimes
       logical,save :: rayleigh
-!$OMP THREADPRIVATE(newtonian,check_cpp_match,force_cpp,testradtimes,rayleigh)
+!$OMP THREADPRIVATE(newtonian,force_cpp,cpp_mugaz_mode,testradtimes,rayleigh)
       logical,save :: stelbbody
       logical,save :: nearco2cond

trunk/LMDZ.GENERIC/libf/phystd/dyn1d/kcm1d.F90

@@ -8 +8 @@
   use comsaison_h, only: mu0, fract, dist_star
   use planete_mod
-  use callkeys_mod, only: check_cpp_match, pceil, tstrat, tracer, global1d
+  use callkeys_mod, only: pceil, tstrat, tracer, global1d
   use inifis_mod, only: inifis
   use comcstfi_mod
@@ -267 +267 @@
   write(*,*) " pceil = ", pceil
 
-  check_cpp_match = .false.
-  call getin("check_cpp_match",check_cpp_match)
-  if (check_cpp_match) then
-     print*,"In 1D modeling, check_cpp_match is supposed to be F"
-     print*,"Please correct callphys.def"
-     stop
-  endif
-
-
 !!! GEOGRAPHICAL INITIALIZATIONS
      !!! AREA. useless in 1D

trunk/LMDZ.GENERIC/libf/phystd/dyn1d/rcm1d.F

@@ -23 +23 @@
       use time_phylmdz_mod, only: daysec, dtphys, day_step, ecritphy,
      &                            nday, iphysiq
-      use callkeys_mod, only: tracer,check_cpp_match,rings_shadow,
+      use callkeys_mod, only: tracer,rings_shadow,
      &                        specOLR,water,pceil,ok_slab_ocean,photochem
       USE comvert_mod, ONLY: ap,bp,aps,bps,pa,preff, sig,
@@ -360 +360 @@
       endif ! of if (tracer)
 
-!!! We have to check that check_cpp_match is F for 1D computations
-!!! We think this check is better than make a particular case for 1D in inifis or calc_cpp_mugaz
-      check_cpp_match = .false.
-      call getin("check_cpp_match",check_cpp_match)
-      if (check_cpp_match) then
-          print*,"In 1D modeling, check_cpp_match is supposed to be F"
-          print*,"Please correct callphys.def"
-          stop
-      endif
-
 !!! GEOGRAPHICAL INITIALIZATIONS
      !!! AREA. useless in 1D

trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90

@@ -240 +240 @@
      if (is_master) write(*,*) trim(rname)//": enertest = ",enertest
 
-     if (is_master) write(*,*) trim(rname)//&
-       ": Check to see if cpp values used match gases.def ?"
-     check_cpp_match=.true. ! default value
-     call getin_p("check_cpp_match",check_cpp_match)
-     if (is_master) write(*,*) trim(rname)//&
-       ": check_cpp_match = ",check_cpp_match
-
      if (is_master) write(*,*) trim(rname)//": call radiative transfer ?"
      callrad=.true. ! default value
@@ -1008 +1001 @@
      if (is_master) write(*,*)trim(rname)//": kmixmin = ",kmixmin
 
-     if (is_master) write(*,*)trim(rname)//&
-       ": Does user want to force cpp and mugaz?"
+     if (is_master) write(*,*)'Predefined Cp from dynamics is ',cpp,'J kg^-1 K^-1'
+     if (is_master) write(*,*)'Predefined Mg from dynamics is ',mugaz,'amu'
+
      force_cpp=.false. ! default value
      call getin_p("force_cpp",force_cpp)
-     if (is_master) write(*,*)trim(rname)//": force_cpp = ",force_cpp
-
      if (force_cpp) then
+      if (is_master) write(*,*)trim(rname)//": force_cpp = ",force_cpp
+      if (is_master) write(*,*)trim(rname)//": force_cpp is deprecated.",&
+      "Set cpp_mugaz_mode=1 in callfis to emulate force_cpp=.true."
+      call abort_physic(rname,"Anyway, you need to set force_cpp=.false. to continue.",1)
+     endif
+
+     if (is_master) write(*,*)trim(rname)//&
+     ": where do you want your cpp/mugaz value to come from?",&
+     "=> 0: dynamics (3d), 1: forced in callfis (1d), 2: computed from gases.def (1d)?"
+     cpp_mugaz_mode = 0 ! default value
+     call getin_p("cpp_mugaz_mode",cpp_mugaz_mode)
+     if (is_master) write(*,*)trim(rname)//": cpp_mugaz_mode = ",cpp_mugaz_mode
+
+     if ((cpp_mugaz_mode >= 1).and.(is_master).and.(ngrid>1)) then
+        write(*,*)'    !!!!! Be aware that having different values for cpp and mugaz in the dynamics and physics'
+        write(*,*)'    in 3D can result in very pathological behavior. You have been warned !!!!!'
+     endif
+
+     if (cpp_mugaz_mode == 1) then
        mugaz = -99999.
        if (is_master) write(*,*)trim(rname)//&
@@ -1020 +1031 @@
        call getin_p("mugaz",mugaz)
        IF (mugaz.eq.-99999.) THEN
-         call abort_physic(rname,"mugaz must be set if force_cpp = T",1)
-       ELSE
-         if (is_master) write(*,*)trim(rname)//": mugaz=",mugaz
+         call abort_physic(rname,"mugaz must be set if cpp_mugaz_mode = 1",1)
        ENDIF
        cpp = -99999.
@@ -1029 +1038 @@
        call getin_p("cpp",cpp)
        IF (cpp.eq.-99999.) THEN
-           PRINT *, "cpp must be set if force_cpp = T"
+           PRINT *, "cpp must be set if cpp_mugaz_mode = 1"
            STOP
-       ELSE
-           if (is_master) write(*,*)trim(rname)//": cpp=",cpp
        ENDIF
-     endif ! of if (force_cpp)
+       if (is_master) write(*,*)'New Cp from callfis is ',cpp,'J kg^-1 K^-1'
+       if (is_master) write(*,*)'New Mg from callfis is ',mugaz,'amu'
+  
+     endif ! of if (cpp_mugaz_mode == 1)
      call su_gases
      call calc_cpp_mugaz

trunk/LMDZ.GENERIC/libf/phystd/tabfi_mod.F90

@@ -61 +61 @@
       use comcstfi_mod, only: rad, omeg, g, mugaz, rcp, cpp, r
       use time_phylmdz_mod, only: dtphys, daysec
-      use callkeys_mod, only: check_cpp_match,force_cpp
+      use callkeys_mod, only: cpp_mugaz_mode
       implicit none
  
@@ -494 +494 @@
         else if (modif(1:len_trim(modif)).eq.'calc_cpp_mugaz') then
           write(*,*) 'current value rcp, mugaz:',rcp,mugaz
-          check_cpp_match=.false.
-          force_cpp=.false.
-          call su_gases
-          call calc_cpp_mugaz
+          cpp_mugaz_mode = 2
+          call su_gases
+          call calc_cpp_mugaz
           write(*,*) 
           write(*,*) ' cpp (new value):',cpp