Changeset 14441
- Timestamp:
- 2021-02-11T17:34:32+01:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser_4.0/src/OCE/ZDF/zdfosm.F90
r14412 r14441 2103 2103 zwb_fk_b(ji,jj) = 0.5 * zwb_fk(ji,jj) * hmle(ji,jj) / hbl(ji,jj) 2104 2104 ENDIF 2105 zvel_max = ( z wstrl(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj)2105 zvel_max = ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj) 2106 2106 IF ( ( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) < 0.0 ) THEN 2107 2107 ! OSBL is deepening, entrainment > restratification 2108 IF ( zdb_bl(ji,jj) > 0.0 .and. zdbdz_bl_ext(ji,jj) > 0.0) THEN2109 zgamma_b_nd = zdbdz_bl_ext(ji,jj) * zdh(ji,jj) / zdb_ml(ji,jj)2110 zpsi = -zalpha_pyc(ji,jj) * ( zwb0tot(ji,jj) - ( zwb_min(ji,jj) + 2.0 * zwb_fk_b(ji,jj)) ) * zdh(ji,jj) / zhbl(ji,jj)2111 zpsi = zpsi - 4.0 * ( 1.0 + zdh(ji,jj) /zhbl(ji,jj) ) * zgamma_b_nd * ( zwb_min(ji,jj) + 2.0 * zwb_fk_b(ji,jj))2108 IF ( zdb_bl(ji,jj) > 1.0e-15 ) THEN 2109 zgamma_b_nd = MAX( zdbdz_bl_ext(ji,jj), 0._wp ) * zdh(ji,jj) / zdb_ml(ji,jj) 2110 zpsi = ( 1.0 - 0.5 * zdh(ji,jj) / zhbl(ji,jj) ) * ( zwb0(ji,jj) - MIN( ( zwb_min(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ), 0._wp ) ) * zdh(ji,jj) / zhbl(ji,jj) 2111 zpsi = zpsi + 1.75 * ( 1.0 - 0.5 * zdh(ji,jj) / zhbl(ji,jj) )*( zdh(ji,jj) / zhbl(ji,jj) + zgamma_b_nd ) * MIN( ( zwb_min(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ), 0._wp ) 2112 2112 zpsi = zalpha_b * MAX ( zpsi, 0._wp ) 2113 zdhdt(ji,jj) = -( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) / ( zvel_max + MAX(zdb_bl(ji,jj), 1.0e-15) ) + zpsi / MAX( zdb_bl(ji,jj), 1.e-15)2113 zdhdt(ji,jj) = -( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) / ( zvel_max + MAX(zdb_bl(ji,jj), 1.0e-15) ) + zpsi / ( zvel_max + MAX( zdb_bl(ji,jj), 1.e-15 ) ) 2114 2114 IF ( j_ddh(ji,jj) == 1 ) THEN 2115 2115 IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN … … 2128 2128 zddhdt = 0._wp 2129 2129 ENDIF ! j_ddh 2130 zdhdt(ji,jj) = zdhdt(ji,jj) + zalpha_b * ( zalpha_pyc(ji,jj) * zdb_ml(ji,jj) + 2.0 * zdbdz_bl_ext(ji,jj) * zdh(ji,jj) ) * zddhdt / zdb_bl(ji,jj)2130 zdhdt(ji,jj) = zdhdt(ji,jj) + zalpha_b * ( 1.0 -0.5 * zdh(ji,jj) / zhbl(ji,jj) ) * zddhdt / ( zvel_max + MAX( zdb_bl(ji,jj), 1.0e-15 ) ) 2131 2131 ELSE ! zdb_bl >0 2132 2132 zdhdt(ji,jj) = -( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) / MAX( zvel_max, 1.0e-15) … … 2334 2334 INTEGER :: jj, ji 2335 2335 INTEGER :: inhml 2336 REAL(wp) :: zari, ztau, zdh_ref, zddhdt 2336 REAL(wp) :: zari, ztau, zdh_ref, zddhdt, zvel_max 2337 2337 REAL, PARAMETER :: a_ddh = 2.5, a_ddh_2 = 3.5 ! also in pycnocline_depth 2338 2338 … … 2342 2342 IF ( lshear(ji,jj) ) THEN 2343 2343 IF ( lconv(ji,jj) ) THEN 2344 IF ( zdb_bl(ji,jj) > 0._wp) THEN2344 IF ( zdb_bl(ji,jj) > 1.0e-15) THEN 2345 2345 IF ( j_ddh(ji,jj) == 0 ) THEN 2346 zvel_max = ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj) 2346 2347 ! ddhdt for pycnocline determined in osm_calculate_dhdt 2347 ! zddhdt = -a_ddh * ( 1.0 - 1.6 * zdh(ji,jj) / zhbl(ji,jj) ) * zwb_ent(ji,jj) / zdb_bl(ji,jj)2348 !zddhdt = EXP( - 4.0 * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX(zustar(ji,jj), 1.e-8 ) ) * zddhdt2348 zddhdt = -a_ddh * ( 1.0 - 1.6 * zdh(ji,jj) / zhbl(ji,jj) ) * zwb_ent(ji,jj) / ( zvel_max + MAX( zdb_bl(ji,jj), 1.0e-15 ) ) 2349 zddhdt = EXP( - 4.0 * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX(zustar(ji,jj), 1.e-8 ) ) * zddhdt 2349 2350 ! maximum limit for how thick the shear layer can grow relative to the thickness of the boundary kayer 2350 ! dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_rdt, 0.625 * hbl(ji,jj) ) 2351 ztau = MAX( zdb_bl(ji,jj) * ( 0.625 * hbl(ji,jj) ) / ( a_ddh * MAX(-zwb_ent(ji,jj), 1.e-12) ), 2.0 * rn_rdt ) 2352 dh(ji,jj) = dh(ji,jj) * EXP( -rn_rdt / ztau ) + 0.625 * zhbl(ji,jj) * ( 1.0 - EXP( -rn_rdt / ztau ) ) 2353 IF ( dh(ji,jj) >= hbl(ji,jj) ) dh(ji,jj) = 0.8 * zhbl(ji,jj) 2351 dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_rdt, 0.625 * hbl(ji,jj) ) 2354 2352 ELSE 2355 ! Temporary (probably) Recalculate dh_ref to ensure dh doesn't go negative. Can't do this using zddhdt from calculate_dhdt2353 ! Need to recalculate because hbl has been updated. 2356 2354 IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN 2357 2355 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 )
Note: See TracChangeset
for help on using the changeset viewer.