Changeset 4482 for LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/pbl_surface_mod.F90

Timestamp:

Mar 29, 2023, 3:14:27 PM (15 months ago)

Author:

lguez

Message:

Sync latest trunk changes to branch LMDZ_ECRad

Location:

LMDZ6/branches/LMDZ_ECRad

Files:

• . (modified) (1 prop)
• libf/phylmdiso/pbl_surface_mod.F90 (modified) (20 diffs)

LMDZ6/branches/LMDZ_ECRad/libf/phylmdiso/pbl_surface_mod.F90

@@ -14 +14 @@
   USE mod_grid_phy_lmdz,   ONLY : klon_glo
   USE ioipsl
-  USE surface_data,        ONLY : type_ocean, ok_veget
+  USE surface_data,        ONLY : type_ocean, ok_veget, landice_opt
   USE surf_land_mod,       ONLY : surf_land
   USE surf_landice_mod,    ONLY : surf_landice
@@ -23 +23 @@
   USE climb_wind_mod,      ONLY : climb_wind_down, climb_wind_up
   USE coef_diff_turb_mod,  ONLY : coef_diff_turb
+  USE atke_exchange_coeff_mod, ONLY :  atke_compute_km_kh
   USE ioipsl_getin_p_mod,  ONLY : getin_p
   USE cdrag_mod
@@ -406 +407 @@
 #endif
     USE ioipsl_getin_p_mod, ONLY : getin_p
-    use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, zsig, zmea
-    use phys_output_var_mod, only: dter, dser, tkt, tks, taur, sss
+    use phys_state_var_mod, only: ds_ns, dt_ns, delta_sst, delta_sal, dter, &
+         dser, dt_ds, zsig, zmea
+    use phys_output_var_mod, only: tkt, tks, taur, sss
 #ifdef CPP_XIOS
     USE wxios, ONLY: missing_val
@@ -1028 +1030 @@
     REAL, DIMENSION(klon)              :: yrmu0
     ! Martin
-
-    REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, ydser, &
-         ytkt, ytks, ytaur, ysss
-    ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser, tkt, tks,
-    ! taur, sss on ocean points
+    REAL, DIMENSION(klon)              :: yri0
+    REAL, DIMENSION(klon):: ydelta_sst, ydelta_sal, yds_ns, ydt_ns, ydter, &
+         ydser, ydt_ds, ytkt, ytks, ytaur, ysss
+    ! compression of delta_sst, delta_sal, ds_ns, dt_ns, dter, dser,
+    ! dt_ds, tkt, tks, taur, sss on ocean points
 
 #ifdef ISO
@@ -1277 +1279 @@
 
     ytke=0.
+    yri0(:)=0.
 !FC 
     y_treedrg=0.
@@ -1813 +1816 @@
              ydelta_sal(:knon) = delta_sal(ni(:knon))
              ydelta_sst(:knon) = delta_sst(ni(:knon))
+             ydter(:knon) = dter(ni(:knon))
+             ydser(:knon) = dser(ni(:knon))
+             ydt_ds(:knon) = dt_ds(ni(:knon))
           end if
           
@@ -1844 +1850 @@
         CALL cdrag(knon, nsrf, &
             speed, yt(:,1), yq(:,1), zgeo1, ypaprs(:,1),&
-            yts, yqsurf, yz0m, yz0h, &
+            yts, yqsurf, yz0m, yz0h, yri0, 0, &
             ycdragm, ycdragh, zri1, pref )
 
@@ -1878 +1884 @@
             CALL cdrag(knon, nsrf, &
             speed_x, yt_x(:,1), yq_x(:,1), zgeo1_x, ypaprs(:,1),&
-            yts_x, yqsurf_x, yz0m, yz0h, &
+            yts_x, yqsurf_x, yz0m, yz0h, yri0, 0, &
             ycdragm_x, ycdragh_x, zri1_x, pref_x )
 
@@ -1905 +1911 @@
         CALL cdrag(knon, nsrf, &
             speed_w, yt_w(:,1), yq_w(:,1), zgeo1_w, ypaprs(:,1),&
-            yts_w, yqsurf_w, yz0m, yz0h, &
+            yts_w, yqsurf_w, yz0m, yz0h, yri0, 0, &
             ycdragm_w, ycdragh_w, zri1_w, pref_w )
 !
@@ -1947 +1953 @@
 
        ENDIF
+
+        IF (iflag_pbl>=50) THEN
+
+        CALL atke_compute_km_kh(knon,klev,yu,yv,yt, &
+             ypplay,ypaprs,ytke,ycoefm, ycoefh)
+
+        ELSE
+
+
         CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
             ypaprs, ypplay, yu, yv, yq, yt, yts, yqsurf, ycdragm, &
@@ -1962 +1977 @@
        ENDDO
        ENDIF
+
+       ENDIF ! iflag_pbl >= 50
+
         IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh ',ycoefh
 !
@@ -1977 +1995 @@
 
        ENDIF
+
+        IF (iflag_pbl>=50) THEN
+
+        CALL atke_compute_km_kh(knon,klev,yu_x,yv_x,yt_x, &
+             ypplay,ypaprs,ytke_x,ycoefm_x, ycoefh_x)
+
+        ELSE
+
+
         CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
             ypaprs, ypplay, yu_x, yv_x, yq_x, yt_x, yts_x, yqsurf_x, ycdragm_x, &
@@ -1992 +2019 @@
        ENDDO
        ENDIF
+
+       ENDIF ! iflag_pbl >= 50
+
         IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh_x ',ycoefh_x
 !
@@ -2005 +2035 @@
       print *,' args coef_diff_turb: ytke_w ', ytke_w
        ENDIF
+
+        IF (iflag_pbl>=50) THEN
+
+        CALL atke_compute_km_kh(knon,klev,yu_w,yv_w,yt_w, &
+             ypplay,ypaprs,ytke_w,ycoefm_w, ycoefh_w)
+
+        ELSE
+
+
         CALL coef_diff_turb(dtime, nsrf, knon, ni,  &
             ypaprs, ypplay, yu_w, yv_w, yq_w, yt_w, yts_w, yqsurf_w, ycdragm_w, &
@@ -2019 +2058 @@
        ENDDO
        ENDIF
+
+       ENDIF ! iflag_pbl >= 50
+
         IF (prt_level >=10) print *,'coef_diff_turb -> ycoefh_w ',ycoefh_w
 !
@@ -2382 +2424 @@
        CASE(is_lic)
           ! Martin
-          CALL surf_landice(itap, dtime, knon, ni, &
-               rlon, rlat, debut, lafin, &
-               yrmu0, ylwdown, yalb, zgeo1, &
-               ysolsw, ysollw, yts, ypplay(:,1), &
-!!jyg               ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),&
-               ycdragh, ycdragm, yrain_f, ysnow_f, yt1, yq1,&
-               AcoefH, AcoefQ, BcoefH, BcoefQ, &
-               AcoefU, AcoefV, BcoefU, BcoefV, &
-               ypsref, yu1, yv1, ygustiness, yrugoro, pctsrf, &
-               ysnow, yqsurf, yqsol, yagesno, &
-               ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, &
-               ytsurf_new, y_dflux_t, y_dflux_q, &
-               yzmea, yzsig, ycldt, &
-               ysnowhgt, yqsnow, ytoice, ysissnow, &
-               yalb3_new, yrunoff, &
-               y_flux_u1, y_flux_v1 &
-#ifdef ISO
-           &    ,yxtrain_f, yxtsnow_f,yxt1,yRland_ice &
-           &    ,yxtsnow,yxtsol,yxtevap &
+          IF (landice_opt .LT. 2) THEN
+             ! Land ice is treated by LMDZ and not by ORCHIDEE
+             
+             CALL surf_landice(itap, dtime, knon, ni, &
+                  rlon, rlat, debut, lafin, &
+                  yrmu0, ylwdown, yalb, zgeo1, &
+                  ysolsw, ysollw, yts, ypplay(:,1), &
+                  !!jyg               ycdragh, ycdragm, yrain_f, ysnow_f, yt(:,1), yq(:,1),&
+                  ycdragh, ycdragm, yrain_f, ysnow_f, yt1, yq1,&
+                  AcoefH, AcoefQ, BcoefH, BcoefQ, &
+                  AcoefU, AcoefV, BcoefU, BcoefV, &
+                  ypsref, yu1, yv1, ygustiness, yrugoro, pctsrf, &
+                  ysnow, yqsurf, yqsol, yagesno, &
+                  ytsoil, yz0m, yz0h, SFRWL, yalb_dir_new, yalb_dif_new, yevap,yfluxsens,yfluxlat, &
+                  ytsurf_new, y_dflux_t, y_dflux_q, &
+                  yzmea, yzsig, ycldt, &
+                  ysnowhgt, yqsnow, ytoice, ysissnow, &
+                  yalb3_new, yrunoff, &
+                  y_flux_u1, y_flux_v1 &
+#ifdef ISO
+                  &    ,yxtrain_f, yxtsnow_f,yxt1,yRland_ice &
+                  &    ,yxtsnow,yxtsol,yxtevap &
 #endif              
-           &    )
-
-!jyg<
-!!          alb3_lic(:)=0.
-!>jyg
-          DO j = 1, knon
-             i = ni(j)
-             alb3_lic(i) = yalb3_new(j)
-             snowhgt(i)   = ysnowhgt(j)
-             qsnow(i)     = yqsnow(j)
-             to_ice(i)    = ytoice(j)
-             sissnow(i)   = ysissnow(j)
-             runoff(i)    = yrunoff(j)
-          ENDDO
-          ! Martin
-! Special DICE MPL 05082013 puis BOMEX MPL 20150410
-       IF (ok_prescr_ust) THEN
-          DO j=1,knon
-          y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
-          y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
-          ENDDO
-      ENDIF
-
+                  &    )
+             
+             !jyg<
+             !!          alb3_lic(:)=0.
+             !>jyg
+             DO j = 1, knon
+                i = ni(j)
+                alb3_lic(i) = yalb3_new(j)
+                snowhgt(i)   = ysnowhgt(j)
+                qsnow(i)     = yqsnow(j)
+                to_ice(i)    = ytoice(j)
+                sissnow(i)   = ysissnow(j)
+                runoff(i)    = yrunoff(j)
+             ENDDO
+             ! Martin
+             ! Special DICE MPL 05082013 puis BOMEX MPL 20150410
+             IF (ok_prescr_ust) THEN
+                DO j=1,knon
+                   y_flux_u1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yu(j,1)*ypplay(j,1)/RD/yt(j,1)
+                   y_flux_v1(j)=ycdragm(j)*(1.+sqrt(yu(j,1)*yu(j,1)+yv(j,1)*yv(j,1)))*yv(j,1)*ypplay(j,1)/RD/yt(j,1)
+                ENDDO
+             ENDIF
+             
 #ifdef ISOVERIF
-        do j=1,knon
-          do ixt=1,ntraciso
-            call iso_verif_noNaN(yxtevap(ixt,j), &
-         &      'pbl_surface 1095a: apres surf_landice')
-            call iso_verif_noNaN(yxtsol(ixt,j), &
-         &      'pbl_surface 1095b: apres surf_landice')
-          enddo
-        enddo
+             do j=1,knon
+                do ixt=1,ntraciso
+                   call iso_verif_noNaN(yxtevap(ixt,j), &
+                        &      'pbl_surface 1095a: apres surf_landice')
+                   call iso_verif_noNaN(yxtsol(ixt,j), &
+                        &      'pbl_surface 1095b: apres surf_landice')
+                enddo
+             enddo
 #endif
 #ifdef ISOVERIF
-        !write(*,*) 'pbl_surface_mod 1060: sortie surf_landice'
-        do j=1,knon
-          if (iso_eau.gt.0) then     
-                 call iso_verif_egalite(yxtsnow(iso_eau,j), &
-     &                  ysnow(j),'pbl_surf_mod 1064')
-           endif !if (iso_eau.gt.0) then
-        enddo !do i=1,klon
-#endif
-          
+             !write(*,*) 'pbl_surface_mod 1060: sortie surf_landice'
+             do j=1,knon
+                if (iso_eau.gt.0) then     
+                   call iso_verif_egalite(yxtsnow(iso_eau,j), &
+                        &                  ysnow(j),'pbl_surf_mod 1064')
+                endif !if (iso_eau.gt.0) then
+             enddo !do i=1,klon
+#endif
+          END IF
        CASE(is_oce)
            CALL surf_ocean(rlon, rlat, ysolsw, ysollw, yalb_vis, &
@@ -2459 +2504 @@
                y_flux_u1, y_flux_v1, ydelta_sst(:knon), ydelta_sal(:knon), &
                yds_ns(:knon), ydt_ns(:knon), ydter(:knon), ydser(:knon), &
-               ytkt(:knon), ytks(:knon), ytaur(:knon), ysss &
+               ydt_ds(:knon), ytkt(:knon), ytks(:knon), ytaur(:knon), ysss &
 #ifdef ISO
          &      ,yxtrain_f, yxtsnow_f,yxt1,Roce, &
@@ -3391 +3436 @@
           taur(ni(:knon)) = ytaur(:knon)
           sss(ni(:knon)) = ysss(:knon)
+
+          if (activate_ocean_skin == 2 .and. type_ocean == "couple") then
+             dt_ds = missing_val
+             dt_ds(ni(:knon)) = ydt_ds(:knon)
+          end if
        end if
 
@@ -4128 +4178 @@
 
     USE indice_sol_mod
-    use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst
+    use phys_state_var_mod, only: delta_sal, ds_ns, dt_ns, delta_sst, dter, &
+         dser, dt_ds
     use config_ocean_skin_m, only: activate_ocean_skin
 
@@ -4234 +4285 @@
                          delta_sal(i) = 0.
                          delta_sst(i) = 0.
+                         dter(i) = 0.
+                         dser(i) = 0.
+                         dt_ds(i) = 0.
                       end if