Changeset 14555

Timestamp:

2021-02-26T21:49:11+01:00 (3 years ago)

Author:

smueller

Message:

Synchronisation of the OSMOSIS boundary layer scheme with the version developed in branch ^{/NEMO/branches/NERC/dev_r11078_OSMOSIS_IMMERSE_Nurser_4.0: adoption of changeset [14441] (ticket #2353)}

File:

• NEMO/branches/2021/dev_r14122_HPC-08_Mueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90 (modified) (4 diffs)

NEMO/branches/2021/dev_r14122_HPC-08_Mueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90

@@ -1462 +1462 @@
                 zwb_fk_b(ji,jj) = 0.5 * zwb_fk(ji,jj) * hmle(ji,jj) / hbl(ji,jj)
              ENDIF
-             zvel_max = ( zwstrl(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj)
+             zvel_max = ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj)
              IF ( ( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) < 0.0 ) THEN
                 ! OSBL is deepening, entrainment > restratification
-                IF ( zdb_bl(ji,jj) > 0.0 .and. zdbdz_bl_ext(ji,jj) > 0.0 ) THEN
-                   zgamma_b_nd = zdbdz_bl_ext(ji,jj) * zdh(ji,jj) / zdb_ml(ji,jj)
-                   zpsi = -1.0_wp * zalpha_pyc(ji,jj) * ( zwb0(ji,jj) - ( zwb_min(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ) ) * zdh(ji,jj) / zhbl(ji,jj)
-                   zpsi = zpsi - 4.0_wp * ( 1.0_wp + zdh(ji,jj) / zhbl(ji,jj) ) * zgamma_b_nd * ( zwb_min(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) )
+                IF ( zdb_bl(ji,jj) > 1e-15 ) THEN
+                   zgamma_b_nd = MAX( zdbdz_bl_ext(ji,jj), 0.0_wp ) * zdh(ji,jj) / zdb_ml(ji,jj)
+                   zpsi = ( 1.0_wp - 0.5_wp * zdh(ji,jj) / zhbl(ji,jj) ) *   &
+                      &   ( zwb0(ji,jj) - MIN( ( zwb_min(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ), 0.0_wp ) ) * zdh(ji,jj) / zhbl(ji,jj)
+                   zpsi = zpsi + 1.75_wp * ( 1.0_wp - 0.5_wp * zdh(ji,jj) / zhbl(ji,jj) ) *   &
+                      &   ( zdh(ji,jj) / zhbl(ji,jj) + zgamma_b_nd ) * MIN( ( zwb_min(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ), 0.0_wp )
                    zpsi = zalpha_b * MAX( zpsi, 0.0_wp )
-                   zdhdt(ji,jj) = -1.0_wp * ( zwb_ent(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ) / ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15_wp ) )   &
-                      &           + zpsi / MAX( zdb_bl(ji,jj), 1e-15 )
+                   zdhdt(ji,jj) = -1.0_wp * ( zwb_ent(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ) /   &
+                      &                     ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15_wp ) ) +   &
+                      &           zpsi / ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15 ) )
                    IF ( j_ddh(ji,jj) == 1 ) THEN
                      IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN
@@ -1488 +1491 @@
                      zddhdt = 0.0_wp
                   ENDIF ! j_ddh
-                  zdhdt(ji,jj) = zdhdt(ji,jj) + zalpha_b * ( zalpha_pyc(ji,jj) * zdb_ml(ji,jj) + 2.0_wp * zdbdz_bl_ext(ji,jj) * zdh(ji,jj) )   &
-                     &                          * zddhdt / zdb_bl(ji,jj)
+                  zdhdt(ji,jj) = zdhdt(ji,jj) + zalpha_b * ( 1.0_wp - 0.5_wp * zdh(ji,jj) / zhbl(ji,jj) ) * zddhdt / ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15 ) )
                 ELSE    ! zdb_bl >0
                    zdhdt(ji,jj) = -( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) /  MAX( zvel_max, 1.0e-15)
@@ -1695 +1697 @@
       INTEGER :: jj, ji
       INTEGER :: inhml
-      REAL(wp) :: zari, ztau, zdh_ref, zddhdt
+      REAL(wp) :: zari, ztau, zdh_ref, zddhdt, zvel_max
       REAL, PARAMETER :: a_ddh = 2.5, a_ddh_2 = 3.5 ! also in pycnocline_depth
 
@@ -1703 +1705 @@
       IF ( lshear(ji,jj) ) THEN
          IF ( lconv(ji,jj) ) THEN
-          IF ( zdb_bl(ji,jj) > 0.0_wp ) THEN
+          IF ( zdb_bl(ji,jj) > 1e-15_wp ) THEN
            IF ( j_ddh(ji,jj) == 0 ) THEN
+              zvel_max = ( zvstr(ji,jj)**3 + 0.5_wp * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj)
 ! ddhdt for pycnocline determined in osm_calculate_dhdt
-!             zddhdt = -a_ddh * ( 1.0 - 1.6 * zdh(ji,jj) / zhbl(ji,jj) ) * zwb_ent(ji,jj) / zdb_bl(ji,jj)
-!             zddhdt = EXP( - 4.0 * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX(zustar(ji,jj), 1.e-8 ) ) * zddhdt
+              zddhdt = -a_ddh * ( 1.0 - 1.6 * zdh(ji,jj) / zhbl(ji,jj) ) * zwb_ent(ji,jj) / ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15 ) )
+              zddhdt = EXP( -4.0_wp * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX( zustar(ji,jj), 1e-8 ) ) * zddhdt
 ! maximum limit for how thick the shear layer can grow relative to the thickness of the boundary kayer
-!             dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_rdt, 0.625 * hbl(ji,jj) )
-             ztau = MAX( zdb_bl(ji,jj) * ( 0.625_wp * hbl(ji,jj) ) / ( a_ddh * MAX( -1.0_wp * zwb_ent(ji,jj), 1e-12_wp ) ), 2.0_wp * rn_Dt )
-             dh(ji,jj) = dh(ji,jj) * EXP( -1.0_wp * rn_Dt / ztau ) + 0.625_wp * zhbl(ji,jj) * ( 1.0_wp - EXP( -1.0_wp * rn_Dt / ztau ) )
-             IF ( dh(ji,jj) >= hbl(ji,jj) ) dh(ji,jj) = 0.8_wp * zhbl(ji,jj)
+             dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_Dt, 0.625_wp * hbl(ji,jj) )
            ELSE
-! Temporary (probably) Recalculate dh_ref to ensure dh doesn't go negative. Can't do this using zddhdt from calculate_dhdt
+! Need to recalculate because hbl has been updated.
              IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN
                zari = MIN( 1.5 * zdb_bl(ji,jj) / ( zhbl(ji,jj) * ( MAX(zdbdz_bl_ext(ji,jj),0._wp) + zdb_bl(ji,jj)**2 / MAX(4.5 * zvstr(ji,jj)**2 , 1.e-12 )) ), 0.2d0 )