Changeset 13295 for NEMO/trunk/src/OCE/DYN/dynspg_ts.F90
- Timestamp:
- 2020-07-10T20:24:21+02:00 (4 years ago)
- File:
-
- 1 edited
-
NEMO/trunk/src/OCE/DYN/dynspg_ts.F90 (modified) (42 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/DYN/dynspg_ts.F90
r13289 r13295 264 264 IF( ln_wd_il ) THEN ! W/D : limiter applied to spgspg 265 265 CALL wad_spg( pssh(:,:,Kmm), zcpx, zcpy ) ! Calculating W/D gravity filters, zcpx and zcpy 266 DO_2D _00_00266 DO_2D( 0, 0, 0, 0 ) 267 267 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( pssh(ji+1,jj ,Kmm) - pssh(ji ,jj ,Kmm) ) & 268 268 & * r1_e1u(ji,jj) * zcpx(ji,jj) * wdrampu(ji,jj) !jth … … 271 271 END_2D 272 272 ELSE ! now suface pressure gradient 273 DO_2D _00_00273 DO_2D( 0, 0, 0, 0 ) 274 274 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( pssh(ji+1,jj ,Kmm) - pssh(ji ,jj ,Kmm) ) * r1_e1u(ji,jj) 275 275 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( pssh(ji ,jj+1,Kmm) - pssh(ji ,jj ,Kmm) ) * r1_e2v(ji,jj) … … 279 279 ENDIF 280 280 ! 281 DO_2D _00_00281 DO_2D( 0, 0, 0, 0 ) 282 282 zu_frc(ji,jj) = zu_frc(ji,jj) - zu_trd(ji,jj) * ssumask(ji,jj) 283 283 zv_frc(ji,jj) = zv_frc(ji,jj) - zv_trd(ji,jj) * ssvmask(ji,jj) … … 291 291 IF( ln_apr_dyn ) THEN 292 292 IF( ln_bt_fw ) THEN ! FORWARD integration: use kt+1/2 pressure (NOW+1/2) 293 DO_2D _00_00293 DO_2D( 0, 0, 0, 0 ) 294 294 zu_frc(ji,jj) = zu_frc(ji,jj) + grav * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) ) * r1_e1u(ji,jj) 295 295 zv_frc(ji,jj) = zv_frc(ji,jj) + grav * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) ) * r1_e2v(ji,jj) … … 297 297 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 pressure (NOW) 298 298 zztmp = grav * r1_2 299 DO_2D _00_00299 DO_2D( 0, 0, 0, 0 ) 300 300 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) & 301 301 & + ssh_ibb(ji+1,jj ) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) … … 309 309 ! ! ---------------------------------- ! 310 310 IF( ln_bt_fw ) THEN ! Add wind forcing 311 DO_2D _00_00311 DO_2D( 0, 0, 0, 0 ) 312 312 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rho0 * utau(ji,jj) * r1_hu(ji,jj,Kmm) 313 313 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rho0 * vtau(ji,jj) * r1_hv(ji,jj,Kmm) … … 315 315 ELSE 316 316 zztmp = r1_rho0 * r1_2 317 DO_2D _00_00317 DO_2D( 0, 0, 0, 0 ) 318 318 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu(ji,jj,Kmm) 319 319 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv(ji,jj,Kmm) … … 475 475 ! 476 476 ! ! ocean u- and v-depth at mid-step (separate DO-loops remove the need of a lbc_lnk) 477 DO_2D _11_10477 DO_2D( 1, 1, 1, 0 ) 478 478 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & 479 479 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 480 480 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 481 481 END_2D 482 DO_2D _10_11482 DO_2D( 1, 0, 1, 1 ) 483 483 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj) & 484 484 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & … … 515 515 !-- ssh = ssh - delta_t' * [ frc + div( flux ) ] --! 516 516 !-------------------------------------------------------------------------! 517 DO_2D _00_00517 DO_2D( 0, 0, 0, 0 ) 518 518 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 519 519 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rDt_e * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) … … 541 541 ! Sea Surface Height at u-,v-points (vvl case only) 542 542 IF( .NOT.ln_linssh ) THEN 543 DO_2D _00_00543 DO_2D( 0, 0, 0, 0 ) 544 544 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 545 545 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & … … 561 561 ! ! Surface pressure gradient 562 562 zldg = ( 1._wp - rn_scal_load ) * grav ! local factor 563 DO_2D _00_00563 DO_2D( 0, 0, 0, 0 ) 564 564 zu_spg(ji,jj) = - zldg * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 565 565 zv_spg(ji,jj) = - zldg * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) … … 579 579 ! Add tidal astronomical forcing if defined 580 580 IF ( ln_tide .AND. ln_tide_pot ) THEN 581 DO_2D _00_00581 DO_2D( 0, 0, 0, 0 ) 582 582 zu_trd(ji,jj) = zu_trd(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 583 583 zv_trd(ji,jj) = zv_trd(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) … … 588 588 !jth do implicitly instead 589 589 IF ( .NOT. ll_wd ) THEN ! Revert to explicit for bit comparison tests in non wad runs 590 DO_2D _00_00590 DO_2D( 0, 0, 0, 0 ) 591 591 zu_trd(ji,jj) = zu_trd(ji,jj) + zCdU_u(ji,jj) * un_e(ji,jj) * hur_e(ji,jj) 592 592 zv_trd(ji,jj) = zv_trd(ji,jj) + zCdU_v(ji,jj) * vn_e(ji,jj) * hvr_e(ji,jj) … … 606 606 !------------------------------------------------------------------------------------------------------------------------! 607 607 IF( ln_dynadv_vec .OR. ln_linssh ) THEN !* Vector form 608 DO_2D _00_00608 DO_2D( 0, 0, 0, 0 ) 609 609 ua_e(ji,jj) = ( un_e(ji,jj) & 610 610 & + rDt_e * ( zu_spg(ji,jj) & … … 621 621 ! 622 622 ELSE !* Flux form 623 DO_2D _00_00623 DO_2D( 0, 0, 0, 0 ) 624 624 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 625 625 ! ! backward interpolated depth used in spg terms at jn+1/2 … … 645 645 !jth implicit bottom friction: 646 646 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 647 DO_2D _00_00647 DO_2D( 0, 0, 0, 0 ) 648 648 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rDt_e * zCdU_u(ji,jj) * hur_e(ji,jj)) 649 649 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rDt_e * zCdU_v(ji,jj) * hvr_e(ji,jj)) … … 712 712 IF (ln_bt_fw) THEN 713 713 IF( .NOT.( kt == nit000 .AND. l_1st_euler ) ) THEN 714 DO_2D _11_11714 DO_2D( 1, 1, 1, 1 ) 715 715 zun_save = un_adv(ji,jj) 716 716 zvn_save = vn_adv(ji,jj) … … 743 743 ELSE 744 744 ! At this stage, pssh(:,:,:,Krhs) has been corrected: compute new depths at velocity points 745 DO_2D _10_10745 DO_2D( 1, 0, 1, 0 ) 746 746 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 747 747 & * ( e1e2t(ji ,jj) * pssh(ji ,jj,Kaa) & … … 975 975 ! Max courant number for ext. grav. waves 976 976 ! 977 DO_2D _00_00977 DO_2D( 0, 0, 0, 0 ) 978 978 zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 979 979 zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) … … 1100 1100 SELECT CASE( nn_een_e3f ) !* ff_f/e3 at F-point 1101 1101 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 1102 DO_2D _10_101102 DO_2D( 1, 0, 1, 0 ) 1103 1103 zwz(ji,jj) = ( ht(ji ,jj+1) + ht(ji+1,jj+1) + & 1104 1104 & ht(ji ,jj ) + ht(ji+1,jj ) ) * 0.25_wp … … 1106 1106 END_2D 1107 1107 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 1108 DO_2D _10_101108 DO_2D( 1, 0, 1, 0 ) 1109 1109 zwz(ji,jj) = ( ht (ji ,jj+1) + ht (ji+1,jj+1) & 1110 1110 & + ht (ji ,jj ) + ht (ji+1,jj ) ) & … … 1117 1117 ! 1118 1118 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1119 DO_2D _01_011119 DO_2D( 0, 1, 0, 1 ) 1120 1120 ftne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) 1121 1121 ftnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj ) … … 1126 1126 CASE( np_EET ) != EEN scheme using e3t energy conserving scheme 1127 1127 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1128 DO_2D _01_011128 DO_2D( 0, 1, 0, 1 ) 1129 1129 z1_ht = ssmask(ji,jj) / ( ht(ji,jj) + 1._wp - ssmask(ji,jj) ) 1130 1130 ftne(ji,jj) = ( ff_f(ji-1,jj ) + ff_f(ji ,jj ) + ff_f(ji ,jj-1) ) * z1_ht … … 1159 1159 ! 1160 1160 !zhf(:,:) = hbatf(:,:) 1161 DO_2D _10_101161 DO_2D( 1, 0, 1, 0 ) 1162 1162 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) & 1163 1163 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) & … … 1178 1178 CALL lbc_lnk( 'dynspg_ts', zhf, 'F', 1._wp ) 1179 1179 ! JC: TBC. hf should be greater than 0 1180 DO_2D _11_111180 DO_2D( 1, 1, 1, 1 ) 1181 1181 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj) 1182 1182 END_2D … … 1201 1201 SELECT CASE( nvor_scheme ) 1202 1202 CASE( np_ENT ) ! enstrophy conserving scheme (f-point) 1203 DO_2D _00_001203 DO_2D( 0, 0, 0, 0 ) 1204 1204 z1_hu = ssumask(ji,jj) / ( phu(ji,jj) + 1._wp - ssumask(ji,jj) ) 1205 1205 z1_hv = ssvmask(ji,jj) / ( phv(ji,jj) + 1._wp - ssvmask(ji,jj) ) … … 1214 1214 ! 1215 1215 CASE( np_ENE , np_MIX ) ! energy conserving scheme (t-point) ENE or MIX 1216 DO_2D _00_001216 DO_2D( 0, 0, 0, 0 ) 1217 1217 zy1 = ( zhV(ji,jj-1) + zhV(ji+1,jj-1) ) * r1_e1u(ji,jj) 1218 1218 zy2 = ( zhV(ji,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) … … 1225 1225 ! 1226 1226 CASE( np_ENS ) ! enstrophy conserving scheme (f-point) 1227 DO_2D _00_001227 DO_2D( 0, 0, 0, 0 ) 1228 1228 zy1 = r1_8 * ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) & 1229 1229 & + zhV(ji ,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) … … 1235 1235 ! 1236 1236 CASE( np_EET , np_EEN ) ! energy & enstrophy scheme (using e3t or e3f) 1237 DO_2D _00_001237 DO_2D( 0, 0, 0, 0 ) 1238 1238 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zhV(ji ,jj ) & 1239 1239 & + ftnw(ji+1,jj) * zhV(ji+1,jj ) & … … 1269 1269 ! 1270 1270 IF( ln_wd_dl_rmp ) THEN 1271 DO_2D _11_111271 DO_2D( 1, 1, 1, 1 ) 1272 1272 IF ( pssh(ji,jj) + ht_0(ji,jj) > 2._wp * rn_wdmin1 ) THEN 1273 1273 ! IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin2 ) THEN … … 1280 1280 END_2D 1281 1281 ELSE 1282 DO_2D _11_111282 DO_2D( 1, 1, 1, 1 ) 1283 1283 IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN ; ptmsk(ji,jj) = 1._wp 1284 1284 ELSE ; ptmsk(ji,jj) = 0._wp … … 1308 1308 !!---------------------------------------------------------------------- 1309 1309 ! 1310 DO_2D _11_101310 DO_2D( 1, 1, 1, 0 ) 1311 1311 IF ( phU(ji,jj) > 0._wp ) THEN ; pUmsk(ji,jj) = pTmsk(ji ,jj) 1312 1312 ELSE ; pUmsk(ji,jj) = pTmsk(ji+1,jj) … … 1316 1316 END_2D 1317 1317 ! 1318 DO_2D _10_111318 DO_2D( 1, 0, 1, 1 ) 1319 1319 IF ( phV(ji,jj) > 0._wp ) THEN ; pVmsk(ji,jj) = pTmsk(ji,jj ) 1320 1320 ELSE ; pVmsk(ji,jj) = pTmsk(ji,jj+1) … … 1338 1338 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: zcpx, zcpy 1339 1339 !!---------------------------------------------------------------------- 1340 DO_2D _00_001340 DO_2D( 0, 0, 0, 0 ) 1341 1341 ll_tmp1 = MIN( pshn(ji,jj) , pshn(ji+1,jj) ) > & 1342 1342 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) .AND. & … … 1407 1407 IF( ln_isfcav ) THEN ! top+bottom friction (ocean cavities) 1408 1408 1409 DO_2D _00_001409 DO_2D( 0, 0, 0, 0 ) 1410 1410 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 1411 1411 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) 1412 1412 END_2D 1413 1413 ELSE ! bottom friction only 1414 DO_2D _00_001414 DO_2D( 0, 0, 0, 0 ) 1415 1415 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 1416 1416 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) … … 1422 1422 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW bottom baroclinic velocities 1423 1423 1424 DO_2D _00_001424 DO_2D( 0, 0, 0, 0 ) 1425 1425 ikbu = mbku(ji,jj) 1426 1426 ikbv = mbkv(ji,jj) … … 1430 1430 ELSE ! CENTRED integration: use BEFORE bottom baroclinic velocities 1431 1431 1432 DO_2D _00_001432 DO_2D( 0, 0, 0, 0 ) 1433 1433 ikbu = mbku(ji,jj) 1434 1434 ikbv = mbkv(ji,jj) … … 1440 1440 IF( ln_wd_il ) THEN ! W/D : use the "clipped" bottom friction !!gm explain WHY, please ! 1441 1441 zztmp = -1._wp / rDt_e 1442 DO_2D _00_001442 DO_2D( 0, 0, 0, 0 ) 1443 1443 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( & 1444 1444 & r1_hu(ji,jj,Kmm) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) … … 1448 1448 ELSE ! use "unclipped" drag (even if explicit friction is used in 3D calculation) 1449 1449 1450 DO_2D _00_001450 DO_2D( 0, 0, 0, 0 ) 1451 1451 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + r1_hu(ji,jj,Kmm) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) * zu_i(ji,jj) 1452 1452 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + r1_hv(ji,jj,Kmm) * r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) * zv_i(ji,jj) … … 1460 1460 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW top baroclinic velocity 1461 1461 1462 DO_2D _00_001462 DO_2D( 0, 0, 0, 0 ) 1463 1463 iktu = miku(ji,jj) 1464 1464 iktv = mikv(ji,jj) … … 1468 1468 ELSE ! CENTRED integration: use BEFORE top baroclinic velocity 1469 1469 1470 DO_2D _00_001470 DO_2D( 0, 0, 0, 0 ) 1471 1471 iktu = miku(ji,jj) 1472 1472 iktv = mikv(ji,jj) … … 1478 1478 ! ! use "unclipped" top drag (even if explicit friction is used in 3D calculation) 1479 1479 1480 DO_2D _00_001480 DO_2D( 0, 0, 0, 0 ) 1481 1481 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + r1_hu(ji,jj,Kmm) * r1_2*( rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) * zu_i(ji,jj) 1482 1482 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + r1_hv(ji,jj,Kmm) * 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.
