- Timestamp:
- 2020-09-15T12:56:56+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/temporary_r4_trunk/src/OCE/DYN/dynspg_ts.F90
r13469 r13470 253 253 IF( ln_wd_il ) THEN ! W/D : limiter applied to spgspg 254 254 CALL wad_spg( sshn, zcpx, zcpy ) ! Calculating W/D gravity filters, zcpx and zcpy 255 DO_2D _00_00255 DO_2D( 0, 0, 0, 0 ) 256 256 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 257 257 & * r1_e1u(ji,jj) * zcpx(ji,jj) * wdrampu(ji,jj) !jth … … 260 260 END_2D 261 261 ELSE ! now suface pressure gradient 262 DO_2D _00_00262 DO_2D( 0, 0, 0, 0 ) 263 263 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) * r1_e1u(ji,jj) 264 264 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) * r1_e2v(ji,jj) … … 268 268 ENDIF 269 269 ! 270 DO_2D _00_00270 DO_2D( 0, 0, 0, 0 ) 271 271 zu_frc(ji,jj) = zu_frc(ji,jj) - zu_trd(ji,jj) * ssumask(ji,jj) 272 272 zv_frc(ji,jj) = zv_frc(ji,jj) - zv_trd(ji,jj) * ssvmask(ji,jj) … … 281 281 IF( ln_apr_dyn ) THEN 282 282 IF( ln_bt_fw ) THEN ! FORWARD integration: use kt+1/2 pressure (NOW+1/2) 283 DO_2D _00_00283 DO_2D( 0, 0, 0, 0 ) 284 284 zu_frc(ji,jj) = zu_frc(ji,jj) + grav * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) ) * r1_e1u(ji,jj) 285 285 zv_frc(ji,jj) = zv_frc(ji,jj) + grav * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) ) * r1_e2v(ji,jj) … … 287 287 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 pressure (NOW) 288 288 zztmp = grav * r1_2 289 DO_2D _00_00289 DO_2D( 0, 0, 0, 0 ) 290 290 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) & 291 291 & + ssh_ibb(ji+1,jj ) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) … … 299 299 ! ! ---------------------------------- ! 300 300 IF( ln_bt_fw ) THEN ! Add wind forcing 301 DO_2D _00_00301 DO_2D( 0, 0, 0, 0 ) 302 302 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 303 303 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) … … 305 305 ELSE 306 306 zztmp = r1_rau0 * r1_2 307 DO_2D _00_00307 DO_2D( 0, 0, 0, 0 ) 308 308 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 309 309 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) … … 443 443 ! 444 444 ! ! ocean u- and v-depth at mid-step (separate DO-loops remove the need of a lbc_lnk) 445 DO_2D _11_10445 DO_2D( 1, 1, 1, 0 ) 446 446 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & 447 447 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 448 448 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 449 449 END_2D 450 DO_2D _10_11450 DO_2D( 1, 0, 1, 1 ) 451 451 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj) & 452 452 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & … … 508 508 !-- ssh = ssh - delta_t' * [ frc + div( flux ) ] --! 509 509 !-------------------------------------------------------------------------! 510 DO_2D _00_00510 DO_2D( 0, 0, 0, 0 ) 511 511 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 512 512 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rdtbt * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) … … 533 533 ! Sea Surface Height at u-,v-points (vvl case only) 534 534 IF( .NOT.ln_linssh ) THEN 535 DO_2D _00_00535 DO_2D( 0, 0, 0, 0 ) 536 536 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 537 537 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & … … 553 553 ! ! Surface pressure gradient 554 554 zldg = ( 1._wp - rn_scal_load ) * grav ! local factor 555 DO_2D _00_00555 DO_2D( 0, 0, 0, 0 ) 556 556 zu_spg(ji,jj) = - zldg * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 557 557 zv_spg(ji,jj) = - zldg * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) … … 571 571 ! Add tidal astronomical forcing if defined 572 572 IF ( ln_tide .AND. ln_tide_pot ) THEN 573 DO_2D _00_00573 DO_2D( 0, 0, 0, 0 ) 574 574 zu_trd(ji,jj) = zu_trd(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 575 575 zv_trd(ji,jj) = zv_trd(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) … … 580 580 !jth do implicitly instead 581 581 IF ( .NOT. ll_wd ) THEN ! Revert to explicit for bit comparison tests in non wad runs 582 DO_2D _00_00582 DO_2D( 0, 0, 0, 0 ) 583 583 zu_trd(ji,jj) = zu_trd(ji,jj) + zCdU_u(ji,jj) * un_e(ji,jj) * hur_e(ji,jj) 584 584 zv_trd(ji,jj) = zv_trd(ji,jj) + zCdU_v(ji,jj) * vn_e(ji,jj) * hvr_e(ji,jj) … … 598 598 !------------------------------------------------------------------------------------------------------------------------! 599 599 IF( ln_dynadv_vec .OR. ln_linssh ) THEN !* Vector form 600 DO_2D _00_00600 DO_2D( 0, 0, 0, 0 ) 601 601 ua_e(ji,jj) = ( un_e(ji,jj) & 602 602 & + rdtbt * ( zu_spg(ji,jj) & … … 613 613 ! 614 614 ELSE !* Flux form 615 DO_2D _00_00615 DO_2D( 0, 0, 0, 0 ) 616 616 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 617 617 ! ! backward interpolated depth used in spg terms at jn+1/2 … … 637 637 !jth implicit bottom friction: 638 638 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 639 DO_2D _00_00639 DO_2D( 0, 0, 0, 0 ) 640 640 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rdtbt * zCdU_u(ji,jj) * hur_e(ji,jj)) 641 641 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rdtbt * zCdU_v(ji,jj) * hvr_e(ji,jj)) … … 703 703 IF (ln_bt_fw) THEN 704 704 IF( .NOT.( kt == nit000 .AND. neuler==0 ) ) THEN 705 DO_2D _11_11705 DO_2D( 1, 1, 1, 1 ) 706 706 zun_save = un_adv(ji,jj) 707 707 zvn_save = vn_adv(ji,jj) … … 734 734 ELSE 735 735 ! At this stage, ssha has been corrected: compute new depths at velocity points 736 DO_2D _10_10736 DO_2D( 1, 0, 1, 0 ) 737 737 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 738 738 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & … … 969 969 ! Max courant number for ext. grav. waves 970 970 ! 971 DO_2D _11_11971 DO_2D( 1, 1, 1, 1 ) 972 972 zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 973 973 zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) … … 1093 1093 SELECT CASE( nn_een_e3f ) !* ff_f/e3 at F-point 1094 1094 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 1095 DO_2D _10_101095 DO_2D( 1, 0, 1, 0 ) 1096 1096 zwz(ji,jj) = ( ht_n(ji ,jj+1) + ht_n(ji+1,jj+1) + & 1097 1097 & ht_n(ji ,jj ) + ht_n(ji+1,jj ) ) * 0.25_wp … … 1099 1099 END_2D 1100 1100 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 1101 DO_2D _10_101101 DO_2D( 1, 0, 1, 0 ) 1102 1102 zwz(ji,jj) = ( ht_n (ji ,jj+1) + ht_n (ji+1,jj+1) & 1103 1103 & + ht_n (ji ,jj ) + ht_n (ji+1,jj ) ) & … … 1110 1110 ! 1111 1111 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1112 DO_2D _01_011112 DO_2D( 0, 1, 0, 1 ) 1113 1113 ftne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) 1114 1114 ftnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj ) … … 1119 1119 CASE( np_EET ) != EEN scheme using e3t (energy conserving scheme) 1120 1120 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1121 DO_2D _01_011121 DO_2D( 0, 1, 0, 1 ) 1122 1122 z1_ht = ssmask(ji,jj) / ( ht_n(ji,jj) + 1._wp - ssmask(ji,jj) ) 1123 1123 ftne(ji,jj) = ( ff_f(ji-1,jj ) + ff_f(ji ,jj ) + ff_f(ji ,jj-1) ) * z1_ht … … 1152 1152 ! 1153 1153 !zhf(:,:) = hbatf(:,:) 1154 DO_2D _10_101154 DO_2D( 1, 0, 1, 0 ) 1155 1155 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) & 1156 1156 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) & … … 1171 1171 CALL lbc_lnk( 'dynspg_ts', zhf, 'F', 1._wp ) 1172 1172 ! JC: TBC. hf should be greater than 0 1173 DO_2D _11_111173 DO_2D( 1, 1, 1, 1 ) 1174 1174 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj) 1175 1175 END_2D … … 1194 1194 SELECT CASE( nvor_scheme ) 1195 1195 CASE( np_ENT ) ! enstrophy conserving scheme (f-point) 1196 DO_2D _00_001196 DO_2D( 0, 0, 0, 0 ) 1197 1197 z1_hu = ssumask(ji,jj) / ( hu_n(ji,jj) + 1._wp - ssumask(ji,jj) ) 1198 1198 z1_hv = ssvmask(ji,jj) / ( hv_n(ji,jj) + 1._wp - ssvmask(ji,jj) ) … … 1207 1207 ! 1208 1208 CASE( np_ENE , np_MIX ) ! energy conserving scheme (t-point) ENE or MIX 1209 DO_2D _00_001209 DO_2D( 0, 0, 0, 0 ) 1210 1210 zy1 = ( zhV(ji,jj-1) + zhV(ji+1,jj-1) ) * r1_e1u(ji,jj) 1211 1211 zy2 = ( zhV(ji,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) … … 1218 1218 ! 1219 1219 CASE( np_ENS ) ! enstrophy conserving scheme (f-point) 1220 DO_2D _00_001220 DO_2D( 0, 0, 0, 0 ) 1221 1221 zy1 = r1_8 * ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) & 1222 1222 & + zhV(ji ,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) … … 1228 1228 ! 1229 1229 CASE( np_EET , np_EEN ) ! energy & enstrophy scheme (using e3t or e3f) 1230 DO_2D _00_001230 DO_2D( 0, 0, 0, 0 ) 1231 1231 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zhV(ji ,jj ) & 1232 1232 & + ftnw(ji+1,jj) * zhV(ji+1,jj ) & … … 1262 1262 ! 1263 1263 IF( ln_wd_dl_rmp ) THEN 1264 DO_2D _11_111264 DO_2D( 1, 1, 1, 1 ) 1265 1265 IF ( pssh(ji,jj) + ht_0(ji,jj) > 2._wp * rn_wdmin1 ) THEN 1266 1266 ! IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin2 ) THEN … … 1273 1273 END_2D 1274 1274 ELSE 1275 DO_2D _11_111275 DO_2D( 1, 1, 1, 1 ) 1276 1276 IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN ; ptmsk(ji,jj) = 1._wp 1277 1277 ELSE ; ptmsk(ji,jj) = 0._wp … … 1301 1301 !!---------------------------------------------------------------------- 1302 1302 ! 1303 DO_2D _11_101303 DO_2D( 1, 1, 1, 0 ) 1304 1304 IF ( phU(ji,jj) > 0._wp ) THEN ; pUmsk(ji,jj) = pTmsk(ji ,jj) 1305 1305 ELSE ; pUmsk(ji,jj) = pTmsk(ji+1,jj) … … 1309 1309 END_2D 1310 1310 ! 1311 DO_2D _10_111311 DO_2D( 1, 0, 1, 1 ) 1312 1312 IF ( phV(ji,jj) > 0._wp ) THEN ; pVmsk(ji,jj) = pTmsk(ji,jj ) 1313 1313 ELSE ; pVmsk(ji,jj) = pTmsk(ji,jj+1) … … 1331 1331 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: zcpx, zcpy 1332 1332 !!---------------------------------------------------------------------- 1333 DO_2D _00_001333 DO_2D( 0, 0, 0, 0 ) 1334 1334 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1335 1335 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) .AND. & … … 1397 1397 IF( ln_isfcav.OR.ln_drgice_imp ) THEN ! top+bottom friction (ocean cavities) 1398 1398 1399 DO_2D _00_001399 DO_2D( 0, 0, 0, 0 ) 1400 1400 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 1401 1401 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) 1402 1402 END_2D 1403 1403 ELSE ! bottom friction only 1404 DO_2D _00_001404 DO_2D( 0, 0, 0, 0 ) 1405 1405 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 1406 1406 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) … … 1412 1412 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW bottom baroclinic velocities 1413 1413 1414 DO_2D _00_001414 DO_2D( 0, 0, 0, 0 ) 1415 1415 ikbu = mbku(ji,jj) 1416 1416 ikbv = mbkv(ji,jj) … … 1420 1420 ELSE ! CENTRED integration: use BEFORE bottom baroclinic velocities 1421 1421 1422 DO_2D _00_001422 DO_2D( 0, 0, 0, 0 ) 1423 1423 ikbu = mbku(ji,jj) 1424 1424 ikbv = mbkv(ji,jj) … … 1430 1430 IF( ln_wd_il ) THEN ! W/D : use the "clipped" bottom friction !!gm explain WHY, please ! 1431 1431 zztmp = -1._wp / rdtbt 1432 DO_2D _00_001432 DO_2D( 0, 0, 0, 0 ) 1433 1433 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( & 1434 1434 & r1_hu_n(ji,jj) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) … … 1438 1438 ELSE ! use "unclipped" drag (even if explicit friction is used in 3D calculation) 1439 1439 1440 DO_2D _00_001440 DO_2D( 0, 0, 0, 0 ) 1441 1441 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + r1_hu_n(ji,jj) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zu_i(ji,jj) 1442 1442 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + r1_hv_n(ji,jj) * r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zv_i(ji,jj) … … 1450 1450 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW top baroclinic velocity 1451 1451 1452 DO_2D _00_001452 DO_2D( 0, 0, 0, 0 ) 1453 1453 iktu = miku(ji,jj) 1454 1454 iktv = mikv(ji,jj) … … 1458 1458 ELSE ! CENTRED integration: use BEFORE top baroclinic velocity 1459 1459 1460 DO_2D _00_001460 DO_2D( 0, 0, 0, 0 ) 1461 1461 iktu = miku(ji,jj) 1462 1462 iktv = mikv(ji,jj) … … 1468 1468 ! ! use "unclipped" top drag (even if explicit friction is used in 3D calculation) 1469 1469 1470 DO_2D _00_001470 DO_2D( 0, 0, 0, 0 ) 1471 1471 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + r1_hu_n(ji,jj) * r1_2*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zu_i(ji,jj) 1472 1472 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + r1_hv_n(ji,jj) * r1_2*( rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) * zv_i(ji,jj)
Note: See TracChangeset
for help on using the changeset viewer.