Changeset 14555 for NEMO/branches/2021/dev_r14122_HPC-08_Mueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90
- Timestamp:
- 2021-02-26T21:49:11+01:00 (2 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14122_HPC-08_Mueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90
r14554 r14555 1462 1462 zwb_fk_b(ji,jj) = 0.5 * zwb_fk(ji,jj) * hmle(ji,jj) / hbl(ji,jj) 1463 1463 ENDIF 1464 zvel_max = ( z wstrl(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj)1464 zvel_max = ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj) 1465 1465 IF ( ( zwb_ent(ji,jj) + 2.0 * zwb_fk_b(ji,jj) ) < 0.0 ) THEN 1466 1466 ! OSBL is deepening, entrainment > restratification 1467 IF ( zdb_bl(ji,jj) > 0.0 .and. zdbdz_bl_ext(ji,jj) > 0.0 ) THEN 1468 zgamma_b_nd = zdbdz_bl_ext(ji,jj) * zdh(ji,jj) / zdb_ml(ji,jj) 1469 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) 1470 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) ) 1467 IF ( zdb_bl(ji,jj) > 1e-15 ) THEN 1468 zgamma_b_nd = MAX( zdbdz_bl_ext(ji,jj), 0.0_wp ) * zdh(ji,jj) / zdb_ml(ji,jj) 1469 zpsi = ( 1.0_wp - 0.5_wp * zdh(ji,jj) / zhbl(ji,jj) ) * & 1470 & ( 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) 1471 zpsi = zpsi + 1.75_wp * ( 1.0_wp - 0.5_wp * zdh(ji,jj) / zhbl(ji,jj) ) * & 1472 & ( 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 ) 1471 1473 zpsi = zalpha_b * MAX( zpsi, 0.0_wp ) 1472 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 ) ) & 1473 & + zpsi / MAX( zdb_bl(ji,jj), 1e-15 ) 1474 zdhdt(ji,jj) = -1.0_wp * ( zwb_ent(ji,jj) + 2.0_wp * zwb_fk_b(ji,jj) ) / & 1475 & ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15_wp ) ) + & 1476 & zpsi / ( zvel_max + MAX( zdb_bl(ji,jj), 1e-15 ) ) 1474 1477 IF ( j_ddh(ji,jj) == 1 ) THEN 1475 1478 IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN … … 1488 1491 zddhdt = 0.0_wp 1489 1492 ENDIF ! j_ddh 1490 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) ) & 1491 & * zddhdt / zdb_bl(ji,jj) 1493 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 ) ) 1492 1494 ELSE ! zdb_bl >0 1493 1495 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 1696 1698 INTEGER :: inhml 1697 REAL(wp) :: zari, ztau, zdh_ref, zddhdt 1699 REAL(wp) :: zari, ztau, zdh_ref, zddhdt, zvel_max 1698 1700 REAL, PARAMETER :: a_ddh = 2.5, a_ddh_2 = 3.5 ! also in pycnocline_depth 1699 1701 … … 1703 1705 IF ( lshear(ji,jj) ) THEN 1704 1706 IF ( lconv(ji,jj) ) THEN 1705 IF ( zdb_bl(ji,jj) > 0.0_wp ) THEN1707 IF ( zdb_bl(ji,jj) > 1e-15_wp ) THEN 1706 1708 IF ( j_ddh(ji,jj) == 0 ) THEN 1709 zvel_max = ( zvstr(ji,jj)**3 + 0.5_wp * zwstrc(ji,jj)**3 )**p2third / hbl(ji,jj) 1707 1710 ! ddhdt for pycnocline determined in osm_calculate_dhdt 1708 ! zddhdt = -a_ddh * ( 1.0 - 1.6 * zdh(ji,jj) / zhbl(ji,jj) ) * zwb_ent(ji,jj) / zdb_bl(ji,jj)1709 ! zddhdt = EXP( - 4.0 * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX(zustar(ji,jj), 1.e-8 ) ) * zddhdt1711 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 ) ) 1712 zddhdt = EXP( -4.0_wp * ABS( ff_t(ji,jj) ) * zhbl(ji,jj) / MAX( zustar(ji,jj), 1e-8 ) ) * zddhdt 1710 1713 ! maximum limit for how thick the shear layer can grow relative to the thickness of the boundary kayer 1711 ! dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_rdt, 0.625 * hbl(ji,jj) ) 1712 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 ) 1713 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 ) ) 1714 IF ( dh(ji,jj) >= hbl(ji,jj) ) dh(ji,jj) = 0.8_wp * zhbl(ji,jj) 1714 dh(ji,jj) = MIN( dh(ji,jj) + zddhdt * rn_Dt, 0.625_wp * hbl(ji,jj) ) 1715 1715 ELSE 1716 ! Temporary (probably) Recalculate dh_ref to ensure dh doesn't go negative. Can't do this using zddhdt from calculate_dhdt1716 ! Need to recalculate because hbl has been updated. 1717 1717 IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN 1718 1718 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.