Changeset 5224 for LMDZ6/branches/Amaury_dev/libf/phylmdiso

Timestamp:

Sep 24, 2024, 10:47:17 AM (4 months ago)

Author:

abarral

Message:

Merge r5204 r5205
Light lint
Correct missing IOIPSL includes

Location:

LMDZ6/branches/Amaury_dev/libf/phylmdiso

Files:

• ice_sursat_mod.F90 (deleted)
• lmdz_lscp_condensation.f90 (copied) (copied from LMDZ6/branches/cirrus/libf/phylmdiso/lmdz_lscp_condensation.F90) (1 diff, 1 prop)
• phyetat0_mod.F90 (modified) (6 diffs)
• phyredem.F90 (modified) (2 diffs)
• physiq_mod.F90 (modified) (21 diffs)

LMDZ6/branches/Amaury_dev/libf/phylmdiso/lmdz_lscp_condensation.f90

@@ -1 @@
-link ../phylmd/lmdz_lscp_condensation.F90
+link ../phylmd/lmdz_lscp_condensation.f90

LMDZ6/branches/Amaury_dev/libf/phylmdiso/phyetat0_mod.F90

@@ -1 @@
-! $Id: phyetat0.F90 3890 2021-05-05 15:15:06Z jyg $
-
 MODULE phyetat0_mod
   USE lmdz_abort_physic, ONLY: abort_physic
@@ -26 +24 @@
        du_gwd_rando, du_gwd_front, entr_therm, f0, fm_therm, &
        falb_dir, falb_dif, prw_ancien, prlw_ancien, prsw_ancien, prbsw_ancien, &
-       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, qbs_ancien, rneb_ancien, radpas, radsol, rain_fall, ratqs, &
+       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, qbs_ancien, &
+       cf_ancien, rvc_ancien, radpas, radsol, rain_fall, ratqs, &
        rnebcon, rugoro, sig1, snow_fall, bs_fall, solaire_etat0, sollw, sollwdown, &
        solsw, solswfdiff, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
@@ -421 +420 @@
   ancien_ok=ancien_ok.AND.phyetat0_get(ql_ancien,"QLANCIEN","QLANCIEN",0.)
   ancien_ok=ancien_ok.AND.phyetat0_get(qs_ancien,"QSANCIEN","QSANCIEN",0.)
-  ancien_ok=ancien_ok.AND.phyetat0_get(rneb_ancien,"RNEBANCIEN","RNEBANCIEN",0.)
   ancien_ok=ancien_ok.AND.phyetat0_get(u_ancien,"UANCIEN","UANCIEN",0.)
   ancien_ok=ancien_ok.AND.phyetat0_get(v_ancien,"VANCIEN","VANCIEN",0.)
@@ -437 +435 @@
   ENDIF
 
+  ! cas specifique des variables de la sursaturation par rapport a la glace
+  IF ( ok_ice_supersat ) THEN
+    ancien_ok=ancien_ok.AND.phyetat0_get(cf_ancien,"CFANCIEN","CFANCIEN",0.)
+    ancien_ok=ancien_ok.AND.phyetat0_get(rvc_ancien,"RVCANCIEN","RVCANCIEN",0.)
+  ELSE
+    cf_ancien(:,:)=0.
+    rvc_ancien(:,:)=0.
+  ENDIF
+
   ! Ehouarn: addtional tests to check if t_ancien, q_ancien contain
   !          dummy values (as is the case when generated by ce0l,
@@ -443 +450 @@
        (maxval(ql_ancien)==minval(ql_ancien))     .OR. &
        (maxval(qs_ancien)==minval(qs_ancien))     .OR. &
-       (maxval(rneb_ancien)==minval(rneb_ancien)) .OR. &
        (maxval(prw_ancien)==minval(prw_ancien))   .OR. &
        (maxval(prlw_ancien)==minval(prlw_ancien)) .OR. &
@@ -456 +462 @@
        ancien_ok=.FALSE.
     ENDIF
+  ENDIF
+
+  IF ( ok_ice_supersat ) THEN
+    IF ( (maxval(cf_ancien)==minval(cf_ancien))     .OR. &
+         (maxval(rvc_ancien)==minval(rvc_ancien)) ) THEN
+       ancien_ok=.false.
+     ENDIF
   ENDIF
 

LMDZ6/branches/Amaury_dev/libf/phylmdiso/phyredem.F90

@@ -19 +19 @@
                                 zval, rugoro, t_ancien, q_ancien,            &
                                 prw_ancien, prlw_ancien, prsw_ancien, prbsw_ancien,      &
-                                ql_ancien, qs_ancien, qbs_ancien,  u_ancien, &
-                                v_ancien, clwcon, rnebcon, ratqs, pbl_tke,   &
+                                ql_ancien, qs_ancien, qbs_ancien, cf_ancien, &
+                                rvc_ancien, u_ancien, v_ancien,              &
+                                clwcon, rnebcon, ratqs, pbl_tke,             &
                                 wake_delta_pbl_tke, zmax0, f0, sig1, w01,    &
                                 wake_deltat, wake_deltaq, wake_s, wake_dens, &
@@ -266 +267 @@
        CALL put_field(pass,"QBSANCIEN", "QBSANCIEN", qbs_ancien)
        CALL put_field(pass,"PRBSWANCIEN", "PRBSWANCIEN", prbsw_ancien)
+    ENDIF
+
+    IF ( ok_ice_supersat ) THEN
+      CALL put_field(pass,"CFANCIEN", "CFANCIEN", cf_ancien)
+      CALL put_field(pass,"RVCANCIEN", "RVCANCIEN", rvc_ancien)
     ENDIF
 

LMDZ6/branches/Amaury_dev/libf/phylmdiso/physiq_mod.F90

@@ -1 @@
-
-! $Id: physiq_mod.F90 3908 2021-05-20 07:11:13Z idelkadi $
-
-!#define IO_DEBUG
 MODULE physiq_mod
 
@@ -72 +68 @@
     USE tracinca_mod, ONLY: config_inca
     USE tropopause_m,     ONLY: dyn_tropopause
-    USE ice_sursat_mod,  ONLY: flight_init, airplane
     USE lmdz_vampir
     USE lmdz_writefield_phy
@@ -181 +176 @@
        ! [Variables internes non sauvegardees de la physique]
        ! Variables locales pour effectuer les appels en serie
-       t_seri,q_seri,ql_seri,qs_seri,qbs_seri,u_seri,v_seri,tr_seri,rneb_seri, &
+       t_seri,q_seri,ql_seri,qs_seri,qbs_seri, &
+       u_seri,v_seri,cf_seri,rvc_seri,tr_seri, &
+       rhcl, &
        qx_seri, & ! CR
        rhcl, &        
        ! Dynamic tendencies (diagnostics)
-       d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_qbs_dyn,d_u_dyn,d_v_dyn,d_tr_dyn,d_rneb_dyn, &
+       d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_qbs_dyn, &
+       d_u_dyn,d_v_dyn,d_cf_dyn,d_rvc_dyn,d_tr_dyn, &
        d_q_dyn2d,d_ql_dyn2d,d_qs_dyn2d,d_qbs_dyn2d, &
        ! Physic tendencies
@@ -360 +358 @@
        pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb,  &
        distcltop, temp_cltop,  &
-       zqsatl, zqsats, &
-       qclr, qcld, qss, qvc, rnebclr, rnebss, gamma_ss, &
+       !-- LSCP - condensation and ice supersaturation variables
+       qsub, qissr, qcld, subfra, issrfra, gamma_cond, ratio_qi_qtot, &
+       dcf_sub, dcf_con, dcf_mix, dqi_adj, dqi_sub, dqi_con, dqi_mix, &
+       dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, qsatliq, qsatice, &
+       !-- LSCP - aviation and contrails variables
        Tcontr, qcontr, qcontr2, fcontrN, fcontrP, &
+       dcf_avi, dqi_avi, dqvc_avi, flight_dist, flight_h2o, &
+       !
        cldemi,  &
        cldfra, cldtau, fiwc,  &
@@ -578 +581 @@
 ! reevap -> je commente les 2 lignes au dessus et je laisse la definition
 ! plutot dans infotrac_phy
-    INTEGER,SAVE :: irneb, ibs
-!$OMP THREADPRIVATE(irneb, ibs)
+    INTEGER,SAVE :: irneb, ibs, icf,irvc
+!$OMP THREADPRIVATE(irneb, ibs, icf,irvc)
 
 
@@ -1450 +1453 @@
        iliq = strIdx(tracers(:)%name, addPhase('H2O', 'l'))
        isol = strIdx(tracers(:)%name, addPhase('H2O', 's'))
-       irneb= strIdx(tracers(:)%name, addPhase('H2O', 'r'))
        ibs  = strIdx(tracers(:)%name, addPhase('H2O', 'b'))
+       icf  = strIdx(tracers(:)%name, addPhase('H2O', 'f'))
+       irvc = strIdx(tracers(:)%name, addPhase('H2O', 'c'))
 !       CALL init_etat0_limit_unstruct
-!       IF (.NOT. create_etat0_limit) CALL init_limit_read(days_elapsed)
+       IF (.NOT. create_etat0_limit) CALL init_limit_read(days_elapsed)
        !CR:nvelles variables convection/poches froides
 
@@ -1508 +1512 @@
        ENDIF
 
-       IF (ok_ice_sursat.AND.(iflag_ice_thermo==0)) THEN
-          WRITE (lunout, *) ' ok_ice_sursat=y requires iflag_ice_thermo=1 as well'
+       IF (ok_ice_supersat.AND.(iflag_ice_thermo==0)) THEN
+          WRITE (lunout, *) ' ok_ice_supersat=y requires iflag_ice_thermo=1 as well'
           abort_message='see above'
           CALL abort_physic(modname,abort_message,1)
        ENDIF
 
-       IF (ok_ice_sursat.AND.(nqo<4)) THEN
-          WRITE (lunout, *) ' ok_ice_sursat=y requires 4 H2O tracers ', &
-               '(H2O_g, H2O_l, H2O_s, H2O_r) but nqo=', nqo, '. Might as well stop here.'
+       IF (ok_ice_supersat.AND.(nqo<5)) THEN
+          WRITE (lunout, *) ' ok_ice_supersat=y requires 4 H2O tracers ', &
+               '(H2O_g, H2O_l, H2O_s, H2O_f, H2O_c) but nqo=', nqo, '. Might as well stop here.'
           abort_message='see above'
           CALL abort_physic(modname,abort_message,1)
        ENDIF
 
-       IF (ok_plane_h2o.AND..NOT.ok_ice_sursat) THEN
-          WRITE (lunout, *) ' ok_plane_h2o=y requires ok_ice_sursat=y '
+       IF (ok_plane_h2o.AND..NOT.ok_ice_supersat) THEN
+          WRITE (lunout, *) ' ok_plane_h2o=y requires ok_ice_supersat=y '
           abort_message='see above'
           CALL abort_physic(modname,abort_message,1)
        ENDIF
 
-       IF (ok_plane_contrail.AND..NOT.ok_ice_sursat) THEN
-          WRITE (lunout, *) ' ok_plane_contrail=y requires ok_ice_sursat=y '
+       IF (ok_plane_contrail.AND..NOT.ok_ice_supersat) THEN
+          WRITE (lunout, *) ' ok_plane_contrail=y requires ok_ice_supersat=y '
           abort_message='see above'
           CALL abort_physic(modname,abort_message,1)
@@ -1538 +1542 @@
           CALL abort_physic(modname,abort_message, 1)
 #endif
-         IF ((ok_ice_sursat.AND.nqo <5).OR.(.NOT.ok_ice_sursat.AND.nqo<4)) THEN
+         IF ((ok_ice_supersat.AND.nqo <6).OR.(.NOT.ok_ice_supersat.AND.nqo<4)) THEN
              WRITE (lunout, *) 'activation of blowing snow needs a specific H2O tracer', &
                                'but nqo=', nqo
@@ -1981 +1985 @@
         RG,RD,RCPD,RKAPPA,RLVTT,RETV)
        CALL ratqs_ini(klon,klev,iflag_thermals,lunout,nbsrf,is_lic,is_ter,RG,RV,RD,RCPD,RLSTT,RLVTT,RTT)
-       CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_sursat,iflag_ratqs,fl_cor_ebil,RCPD,RLSTT,RLVTT,RLMLT,RVTMP2,RTT,RD,RG,RV,RPI)
+       CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_supersat,iflag_ratqs,fl_cor_ebil, &
+                     RCPD,RLSTT,RLVTT,RLMLT,RVTMP2,RTT,RD,RV,RG,RPI,EPS_W)
        CALL blowing_snow_ini(RCPD, RLSTT, RLVTT, RLMLT, &
                              RVTMP2, RTT,RD,RG, RV, RPI)
@@ -2458 +2463 @@
                     sollwdown(:))
 
+      !--Init for LSCP - condensation
+      ratio_qi_qtot(:,:) = 0.
 
 
@@ -2565 +2572 @@
           q_seri(i,k)  = qx(i,k,ivap)
           ql_seri(i,k) = qx(i,k,iliq)
-          qbs_seri(i,k) = 0.
+          qbs_seri(i,k)= 0.
+          cf_seri(i,k) = 0.
+          rvc_seri(i,k)= 0.
           !CR: ATTENTION, on rajoute la variable glace
           IF (nqo==2) THEN             !--vapour and liquid only
              qs_seri(i,k) = 0.
-             rneb_seri(i,k) = 0.
           ELSE IF (nqo==3) THEN        !--vapour, liquid and ice
              qs_seri(i,k) = qx(i,k,isol)
-             rneb_seri(i,k) = 0.
-          ELSE IF (nqo>=4) THEN        !--vapour, liquid, ice and rneb and blowing snow
+          ELSE IF (nqo>=4) THEN        !--vapour, liquid, ice, blowing snow, cloud fraction and cloudy water vapor to total water vapor ratio
              qs_seri(i,k) = qx(i,k,isol)
-             IF (ok_ice_sursat) THEN
-               rneb_seri(i,k) = qx(i,k,irneb)
+             IF (ok_ice_supersat) THEN
+               cf_seri(i,k) = qx(i,k,icf)
+               rvc_seri(i,k) = qx(i,k,irvc)
              ENDIF
              IF (ok_bs) THEN
@@ -2748 +2756 @@
        d_ql_dyn(:,:) = (ql_seri(:,:)-ql_ancien(:,:))/phys_tstep
        d_qs_dyn(:,:) = (qs_seri(:,:)-qs_ancien(:,:))/phys_tstep
-       d_qbs_dyn(:,:) = (qbs_seri(:,:)-qbs_ancien(:,:))/phys_tstep
+       d_qbs_dyn(:,:)= (qbs_seri(:,:)-qbs_ancien(:,:))/phys_tstep
+       d_cf_dyn(:,:) = (cf_seri(:,:)-cf_ancien(:,:))/phys_tstep
+       d_rvc_dyn(:,:)= (rvc_seri(:,:)-rvc_ancien(:,:))/phys_tstep
        CALL water_int(klon,klev,q_seri,zmasse,zx_tmp_fi2d)
        d_q_dyn2d(:)=(zx_tmp_fi2d(:)-prw_ancien(:))/phys_tstep
@@ -2760 +2770 @@
        IF (nqtot > nqo) d_tr_dyn(:,:,:)=(tr_seri(:,:,:)-tr_ancien(:,:,:))/phys_tstep
        ! !! RomP <<<
-       !!d_rneb_dyn(:,:)=(rneb_seri(:,:)-rneb_ancien(:,:))/phys_tstep
-       d_rneb_dyn(:,:)=0.0
 
 #ifdef ISO
@@ -2843 +2851 @@
        d_ql_dyn(:,:) = 0.0
        d_qs_dyn(:,:) = 0.0
+       d_qbs_dyn(:,:)= 0.0
+       d_cf_dyn(:,:) = 0.0
+       d_rvc_dyn(:,:)= 0.0
        d_q_dyn2d(:)  = 0.0
        d_ql_dyn2d(:) = 0.0
@@ -2869 +2880 @@
        IF (nqtot > nqo) d_tr_dyn(:,:,:)= 0.0
        ! !! RomP <<<
-       d_rneb_dyn(:,:)=0.0
-       d_qbs_dyn(:,:)=0.0
        ancien_ok = .TRUE.
 #ifdef ISO
@@ -2980 +2989 @@
           ! "zmasse" changes a little.)
        ENDIF
+    ENDIF
+
+    !-- Needed for LSCP - condensation and ice supersaturation
+    IF (ok_ice_supersat) THEN
+      DO k = 1, klev
+        DO i = 1, klon
+          IF ( ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) ) > 0. ) THEN
+            ratio_qi_qtot(i,k) = qs_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
+            rvc_seri(i,k) = rvc_seri(i,k) * q_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
+          ELSE
+            ratio_qi_qtot(i,k) = 0.
+            rvc_seri(i,k) = 0.
+          ENDIF
+        ENDDO
+      ENDDO
     ENDIF
 
@@ -5013 +5037 @@
 
     !--mise à jour de flight_m et flight_h2o dans leur module
-    IF (ok_plane_h2o .OR. ok_plane_contrail) THEN
-      CALL airplane(debut,pphis,pplay,paprs,t_seri)
-    ENDIF
+    !IF (ok_plane_h2o .OR. ok_plane_contrail) THEN
+    !  CALL airplane(debut,pphis,pplay,paprs,t_seri)
+    !ENDIF
 
     CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
          t_seri, q_seri,qs_ini,ptconv,ratqs, &
-         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, rneb_seri, & 
+         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, &
          pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb,  &
          radocond, picefra, rain_lsc, snow_lsc, &
@@ -5025 +5049 @@
          prfl, psfl, rhcl,  &
          zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cld_th, &
-         iflag_ice_thermo, ok_ice_sursat, zqsatl, zqsats, distcltop, temp_cltop,  &
-         pbl_tke(:,:,is_ave), pbl_eps(:,:,is_ave), qclr, qcld, qss, qvc, rnebclr, rnebss, gamma_ss, &
-         Tcontr, qcontr, qcontr2, fcontrN, fcontrP , &
+         iflag_ice_thermo, distcltop, temp_cltop, &
+         pbl_tke(:,:,is_ave), pbl_eps(:,:,is_ave), &
+         cell_area, &
+         cf_seri, rvc_seri, u_seri, v_seri, &
+         qsub, qissr, qcld, subfra, issrfra, gamma_cond, ratio_qi_qtot, &
+         dcf_sub, dcf_con, dcf_mix, dqi_adj, dqi_sub, dqi_con, dqi_mix, &
+         dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, qsatliq, qsatice, &
+         Tcontr, qcontr, qcontr2, fcontrN, fcontrP, &
+         dcf_avi, dqi_avi, dqvc_avi, flight_dist, flight_h2o, &
          cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv, &
          qraindiag, qsnowdiag, dqreva, dqssub, dqrauto, dqrcol, dqrmelt, &
@@ -7072 +7102 @@
              d_qx(i,k,isol) = ( qs_seri(i,k) - qx(i,k,isol) ) / phys_tstep
           ENDIF
-          !--ice_sursat: nqo=4, on ajoute rneb
-          IF (nqo>=4 .AND. ok_ice_sursat) THEN
-             d_qx(i,k,irneb) = ( rneb_seri(i,k) - qx(i,k,irneb) ) / phys_tstep
+          !--ice_supersat: nqo=5, we add cloud fraction and cloudy water vapor to total water vapor ratio
+          IF (nqo>=5 .and. ok_ice_supersat) THEN
+            d_qx(i, k, icf) = (cf_seri(i, k) - qx(i, k, icf)) / phys_tstep
+            d_qx(i, k, irvc) = (rvc_seri(i, k) - qx(i, k, irvc)) / phys_tstep
           ENDIF
 
@@ -7080 +7111 @@
              d_qx(i,k,ibs) = ( qbs_seri(i,k) - qx(i,k,ibs) ) / phys_tstep
           ENDIF
-
 
        ENDDO
@@ -7126 +7156 @@
     ql_ancien(:,:) = ql_seri(:,:)
     qs_ancien(:,:) = qs_seri(:,:)
-    qbs_ancien(:,:) = qbs_seri(:,:)
-    rneb_ancien(:,:) = rneb_seri(:,:)
+    qbs_ancien(:,:)= qbs_seri(:,:)
+    cf_ancien(:,:) = cf_seri(:,:)
+    rvc_ancien(:,:)= rvc_seri(:,:)
 #ifdef ISO
     xt_ancien(:,:,:)=xt_seri(:,:,:)