Changeset 3974

Timestamp:

Jul 30, 2021, 4:27:30 PM (3 years ago)

Author:

asima

Message:

In the "bucket" soil model, soil thermal inertia (formerly a constant value) can also be now a function of soil moisture "qsol"
(F Cheruy's idea; results from "Stage L3" Eve Rebouillat, with E Vignon, A Sima, F Cheruy).

The function depends on iflag_inertie, read from physiq.def (via conf_phys_m.F90), as follows :
(only in physiq.def_NPv6.1 for the time being)

---

# Controle de la definition de l'inertie thermique du sol dans le modele "bucket":
# 0 (default) : constant
# 1 : fonction lineaire de qsol
# 2 : 2 fonctions lineaires de qsol :
# une pour "sable" (boite lat/lon pour Sahara), une pour le reste (limons et argile)
# 3 : fonction lineaire de qsol, valeurs entre 900 pour qsol=0 et 2000 pour qsol=150

iflag_inertie=0

---

Relevant changes are in soil.F90 and physiq.def_NPv6.1 ; the other files are adjusted accordingly.

Location:

LMDZ6/trunk

Files:

• DefLists/physiq.def_NPv6.1 (modified) (1 diff)
• libf/phylmd/comsoil.h (modified) (1 diff)
• libf/phylmd/conf_phys_m.F90 (modified) (4 diffs)
• libf/phylmd/ocean_forced_mod.F90 (modified) (2 diffs)
• libf/phylmd/soil.F90 (modified) (5 diffs)
• libf/phylmd/surf_land_bucket_mod.F90 (modified) (2 diffs)
• libf/phylmd/surf_landice_mod.F90 (modified) (2 diffs)

LMDZ6/trunk/DefLists/physiq.def_NPv6.1

@@ -169 +169 @@
 iflag_albedo=1
 
+    # Controle de la definition de l'inertie thermique du sol dans le modele "bucket":
+    #   0 (default) : constant 
+    #   1 : fonction lineaire de qsol
+    #   2 : 2 fonctions lineaires de qsol :
+    #       une pour "sable" (boite lat/lon pour Sahara), une pour le reste (limons et argile)
+    #   3 : fonction lineaire de qsol, valeurs entre 900 pour qsol=0 et 2000 pour qsol=150
+iflag_inertie=0
+
           # Frequence appel convection. Nombre appels par jour
 nbapp_cv=48

LMDZ6/trunk/libf/phylmd/comsoil.h

@@ -4 +4 @@
 
       common /comsoil/inertie_sol,inertie_sno,inertie_sic,inertie_lic,  &
-     &                iflag_sic
+     &                iflag_sic,iflag_inertie
       real inertie_sol,inertie_sno,inertie_sic,inertie_lic
-      integer iflag_sic
+      integer iflag_sic,iflag_inertie
 !$OMP THREADPRIVATE(/comsoil/)

LMDZ6/trunk/libf/phylmd/conf_phys_m.F90

@@ -183 +183 @@
     REAL,SAVE :: exposant_glace_omp
     REAL,SAVE :: rei_min_omp, rei_max_omp
-    INTEGER,SAVE :: iflag_sic_omp
+    INTEGER,SAVE :: iflag_sic_omp, iflag_inertie_omp
     REAL,SAVE :: inertie_sol_omp,inertie_sno_omp,inertie_sic_omp
     REAL,SAVE :: inertie_lic_omp
@@ -1325 +1325 @@
     CALL getin('iflag_sic',iflag_sic_omp)
     !
+    !Config Key  = iflag_inertie
+    !Config Desc =
+    !Config Def  = 0
+    !Config Help =
+    !
+    iflag_inertie_omp = 0
+    CALL getin('iflag_inertie',iflag_inertie_omp)
+    !
     !Config Key  = inertie_sic
     !Config Desc =  
@@ -2395 +2403 @@
     albsno0 = albsno0_omp
     iflag_sic = iflag_sic_omp
+    iflag_inertie = iflag_inertie_omp
     inertie_sol = inertie_sol_omp
     inertie_sic = inertie_sic_omp
@@ -2881 +2890 @@
     WRITE(lunout,*) ' albsno0 = ', albsno0
     WRITE(lunout,*) ' iflag_sic = ', iflag_sic
+    WRITE(lunout,*) ' iflag_inertie = ', iflag_inertie
     WRITE(lunout,*) ' inertie_sol = ', inertie_sol
     WRITE(lunout,*) ' inertie_sic = ', inertie_sic

LMDZ6/trunk/libf/phylmd/ocean_forced_mod.F90

@@ -180 +180 @@
 !
     USE dimphy
+    USE geometry_mod, ONLY: longitude,latitude
     USE calcul_fluxs_mod
     USE surface_data,     ONLY : calice, calsno
@@ -260 +261 @@
     IF (soil_model) THEN 
 ! update tsoil and calculate soilcap and soilflux
-       CALL soil(dtime, is_sic, knon, snow, tsurf_tmp, tsoil,soilcap, soilflux)
+       CALL soil(dtime, is_sic, knon, snow, tsurf_tmp, qsol, &
+        & longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil,soilcap, soilflux)
        cal(1:knon) = RCPD / soilcap(1:knon)
        radsol(1:knon) = radsol(1:knon)  + soilflux(1:knon)

LMDZ6/trunk/libf/phylmd/soil.F90

@@ -2 +2 @@
 ! $Header$
 !
-SUBROUTINE soil(ptimestep, indice, knon, snow, ptsrf, &
-     ptsoil, pcapcal, pfluxgrd)
+SUBROUTINE soil(ptimestep, indice, knon, snow, ptsrf, qsol, &
+     lon, lat, ptsoil, pcapcal, pfluxgrd)
   
   USE dimphy
@@ -21 +21 @@
 !                            the surface conduction flux pcapcal
 !
+!   Update: 2021/07 : soil thermal inertia, formerly a constant value,
+!   ------   can also be now a function of soil moisture (F Cheruy's idea)
+!            depending on iflag_inertie, read from physiq.def via conf_phys_m.F90
+!            ("Stage L3" Eve Rebouillat, with E Vignon, A Sima, F Cheruy)
 !
 !   Method: Implicit time integration
@@ -48 +52 @@
 !   snow(klon)           snow
 !   ptsrf(klon)          surface temperature at time-step t (K)
+!   qsol(klon)           soil moisture (kg/m2 or mm)
+!   lon(klon)            longitude in radian
+!   lat(klon)            latitude in radian
 !   ptsoil(klon,nsoilmx) temperature inside the ground (K)
 !   pcapcal(klon)        surfacic specific heat (W*m-2*s*K-1)
@@ -60 +67 @@
 ! ---------
   REAL, INTENT(IN)                     :: ptimestep
-  INTEGER, INTENT(IN)                  :: indice, knon
+  INTEGER, INTENT(IN)                  :: indice, knon !, knindex
   REAL, DIMENSION(klon), INTENT(IN)    :: snow
   REAL, DIMENSION(klon), INTENT(IN)    :: ptsrf
-  
+  REAL, DIMENSION(klon), INTENT(IN)    :: qsol
+  REAL, DIMENSION(klon), INTENT(IN)    :: lon
+  REAL, DIMENSION(klon), INTENT(IN)    :: lat
+
   REAL, DIMENSION(klon,nsoilmx), INTENT(INOUT) :: ptsoil
   REAL, DIMENSION(klon), INTENT(OUT)           :: pcapcal
@@ -182 +192 @@
 !      knon, knindex, ztherm_i)
   ELSE IF (indice == is_ter) THEN
+     ! 
+     ! La relation entre l'inertie thermique du sol et qsol change d'apres 
+     !   iflag_inertie, defini dans physiq.def, et appele via comsoil.h
+     ! 
      DO ig = 1, knon
-        ztherm_i(ig)   = inertie_sol
+        ! iflag_inertie=0 correspond au cas inertie=constant, comme avant
+        IF (iflag_inertie==0) THEN         
+           ztherm_i(ig)   = inertie_sol
+        ELSE IF (iflag_inertie == 1) THEN
+          ! I = a_qsol * qsol + b  modele lineaire deduit d'une
+          ! regression lineaire I = a_mrsos * mrsos + b obtenue sur 
+          ! sorties MO d'une simulation LMDZOR(CMIP6) sur l'annee 2000 
+          ! sur tous les points avec frac_snow=0 
+          ! Difference entre qsol et mrsos prise en compte par un
+          ! facteur d'echelle sur le coefficient directeur de regression: 
+          ! fact = 35./150. = mrsos_max/qsol_max
+          ! et a_qsol = a_mrsos * fact (car a = dI/dHumidite)
+            ztherm_i(ig) = 30.0 *35.0/150.0 *qsol(ig) +770.0
+          ! AS : pour qsol entre 0 - 150, on a I entre 770 - 1820
+        ELSE IF (iflag_inertie == 2) THEN
+          ! deux regressions lineaires, sur les memes sorties,  
+          ! distinguant le type de sol : sable ou autre (limons/argile)
+          ! Implementation simple : regression type "sable" seulement pour 
+          ! Sahara, defini par une "boite" lat/lon (NB : en radians !! )
+          IF (lon(ig)>-0.35 .AND. lon(ig)<0.70 .AND. lat(ig)>0.17 .AND. lat(ig)<0.52) THEN 
+              ! Valeurs theoriquement entre 728 et 2373 ; qsol valeurs basses
+              ztherm_i(ig) = 47. *35.0/150.0 *qsol(ig) +728.  ! boite type "sable" pour Sahara
+          ELSE 
+              ! Valeurs theoriquement entre 550 et 1940 ; qsol valeurs moyennes et hautes
+              ztherm_i(ig) = 41. *35.0/150.0 *qsol(ig) +505.
+          ENDIF
+        ELSE IF (iflag_inertie == 3) THEN
+          ! AS : idee a tester : 
+          ! si la relation doit etre une droite, 
+          ! definissons-la en fonction des valeurs min et max de qsol (0:150),
+          ! et de l'inertie (900 : 2000 ou 2400 ; choix ici: 2000)
+          ! I = I_min + qsol * (I_max - I_min)/(qsol_max - qsol_min)
+              ztherm_i(ig) = 900. + qsol(ig) * (2000. - 900.)/150.
+        ELSE          
+          WRITE (lunout,*) "Le choix iflag_inertie = ",iflag_inertie," n'est pas defini. Veuillez choisir un entier entre 0 et 3" 
+        ENDIF
+     !
+     ! Fin de l'introduction de la relation entre l'inertie thermique du sol et qsol
+     !-------------------------------------------
+        !AS : donc le moindre flocon de neige sur un point de grid
+        ! fait que l'inertie du point passe a la valeur pour neige ! 
         IF (snow(ig) > 0.0) ztherm_i(ig)   = inertie_sno
+       
      ENDDO
 !    CALL iophys_ecrit_index('ztherm_ter', 1, 'ztherm_ter', 'USI', &

LMDZ6/trunk/libf/phylmd/surf_land_bucket_mod.F90

@@ -24 +24 @@
     USE cpl_mod
     USE dimphy
-    USE geometry_mod, ONLY: latitude 
+    USE geometry_mod, ONLY: longitude,latitude 
     USE mod_grid_phy_lmdz
     USE mod_phys_lmdz_para
@@ -103 +103 @@
        
 ! calculate temperature, heat capacity and conduction flux in soil
-    IF (soil_model) THEN 
-       CALL soil(dtime, is_ter, knon, snow, tsurf, tsoil, soilcap, soilflux)
+    IF (soil_model) THEN
+       CALL soil(dtime, is_ter, knon, snow, tsurf, qsol,  & 
+        & longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil, soilcap, soilflux)
+
        DO i=1, knon
           cal(i) = RCPD / soilcap(i)

LMDZ6/trunk/libf/phylmd/surf_landice_mod.F90

@@ -25 +25 @@
 
     USE dimphy
+    USE geometry_mod,     ONLY : longitude,latitude
     USE surface_data,     ONLY : type_ocean, calice, calsno, landice_opt, iflag_albcalc
     USE fonte_neige_mod,  ONLY : fonte_neige,run_off_lic,fqcalving_global,ffonte_global,fqfonte_global,runofflic_global
@@ -278 +279 @@
     ! use soil model and recalculate properly cal
     IF (soil_model) THEN 
-       CALL soil(dtime, is_lic, knon, snow, tsurf, tsoil, soilcap, soilflux)
+       CALL soil(dtime, is_lic, knon, snow, tsurf, qsol, &
+        & longitude(knindex(1:knon)), latitude(knindex(1:knon)), tsoil, soilcap, soilflux)
        cal(1:knon) = RCPD / soilcap(1:knon)
        radsol(1:knon)  = radsol(1:knon) + soilflux(1:knon)