- Timestamp:
- 2020-09-15T12:49:18+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/temporary_r4_trunk/src/OCE/DYN/dynspg_ts.F90
r13466 r13469 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 jj = 2, jpjm1 256 DO ji = 2, jpim1 ! SPG with the application of W/D gravity filters 257 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 258 & * r1_e1u(ji,jj) * zcpx(ji,jj) * wdrampu(ji,jj) !jth 259 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) & 260 & * r1_e2v(ji,jj) * zcpy(ji,jj) * wdrampv(ji,jj) !jth 261 END DO 262 END DO 255 DO_2D_00_00 256 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 257 & * r1_e1u(ji,jj) * zcpx(ji,jj) * wdrampu(ji,jj) !jth 258 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) & 259 & * r1_e2v(ji,jj) * zcpy(ji,jj) * wdrampv(ji,jj) !jth 260 END_2D 263 261 ELSE ! now suface pressure gradient 264 DO jj = 2, jpjm1 265 DO ji = fs_2, fs_jpim1 ! vector opt. 266 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) * r1_e1u(ji,jj) 267 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) * r1_e2v(ji,jj) 268 END DO 269 END DO 270 ENDIF 271 ! 272 ENDIF 273 ! 274 DO jj = 2, jpjm1 ! Remove coriolis term (and possibly spg) from barotropic trend 275 DO ji = fs_2, fs_jpim1 276 zu_frc(ji,jj) = zu_frc(ji,jj) - zu_trd(ji,jj) * ssumask(ji,jj) 277 zv_frc(ji,jj) = zv_frc(ji,jj) - zv_trd(ji,jj) * ssvmask(ji,jj) 278 END DO 279 END DO 262 DO_2D_00_00 263 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) * r1_e1u(ji,jj) 264 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) * r1_e2v(ji,jj) 265 END_2D 266 ENDIF 267 ! 268 ENDIF 269 ! 270 DO_2D_00_00 271 zu_frc(ji,jj) = zu_frc(ji,jj) - zu_trd(ji,jj) * ssumask(ji,jj) 272 zv_frc(ji,jj) = zv_frc(ji,jj) - zv_trd(ji,jj) * ssvmask(ji,jj) 273 END_2D 280 274 ! 281 275 ! != Add bottom stress contribution from baroclinic velocities =! … … 287 281 IF( ln_apr_dyn ) THEN 288 282 IF( ln_bt_fw ) THEN ! FORWARD integration: use kt+1/2 pressure (NOW+1/2) 289 DO jj = 2, jpjm1 290 DO ji = fs_2, fs_jpim1 ! vector opt. 291 zu_frc(ji,jj) = zu_frc(ji,jj) + grav * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) ) * r1_e1u(ji,jj) 292 zv_frc(ji,jj) = zv_frc(ji,jj) + grav * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) ) * r1_e2v(ji,jj) 293 END DO 294 END DO 283 DO_2D_00_00 284 zu_frc(ji,jj) = zu_frc(ji,jj) + grav * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) ) * r1_e1u(ji,jj) 285 zv_frc(ji,jj) = zv_frc(ji,jj) + grav * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) ) * r1_e2v(ji,jj) 286 END_2D 295 287 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 pressure (NOW) 296 288 zztmp = grav * r1_2 297 DO jj = 2, jpjm1 298 DO ji = fs_2, fs_jpim1 ! vector opt. 299 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) & 300 & + ssh_ibb(ji+1,jj ) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) 301 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) & 302 & + ssh_ibb(ji ,jj+1) - ssh_ibb(ji,jj) ) * r1_e2v(ji,jj) 303 END DO 304 END DO 289 DO_2D_00_00 290 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( ssh_ib (ji+1,jj ) - ssh_ib (ji,jj) & 291 & + ssh_ibb(ji+1,jj ) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) 292 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( ssh_ib (ji ,jj+1) - ssh_ib (ji,jj) & 293 & + ssh_ibb(ji ,jj+1) - ssh_ibb(ji,jj) ) * r1_e2v(ji,jj) 294 END_2D 305 295 ENDIF 306 296 ENDIF … … 309 299 ! ! ---------------------------------- ! 310 300 IF( ln_bt_fw ) THEN ! Add wind forcing 311 DO jj = 2, jpjm1 312 DO ji = fs_2, fs_jpim1 ! vector opt. 313 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 314 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) 315 END DO 316 END DO 301 DO_2D_00_00 302 zu_frc(ji,jj) = zu_frc(ji,jj) + r1_rau0 * utau(ji,jj) * r1_hu_n(ji,jj) 303 zv_frc(ji,jj) = zv_frc(ji,jj) + r1_rau0 * vtau(ji,jj) * r1_hv_n(ji,jj) 304 END_2D 317 305 ELSE 318 306 zztmp = r1_rau0 * r1_2 319 DO jj = 2, jpjm1 320 DO ji = fs_2, fs_jpim1 ! vector opt. 321 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 322 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 323 END DO 324 END DO 307 DO_2D_00_00 308 zu_frc(ji,jj) = zu_frc(ji,jj) + zztmp * ( utau_b(ji,jj) + utau(ji,jj) ) * r1_hu_n(ji,jj) 309 zv_frc(ji,jj) = zv_frc(ji,jj) + zztmp * ( vtau_b(ji,jj) + vtau(ji,jj) ) * r1_hv_n(ji,jj) 310 END_2D 325 311 ENDIF 326 312 ! … … 457 443 ! 458 444 ! ! ocean u- and v-depth at mid-step (separate DO-loops remove the need of a lbc_lnk) 459 DO jj = 1, jpj 460 DO ji = 1, jpim1 ! not jpi-column 461 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & 462 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 463 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 464 END DO 465 END DO 466 DO jj = 1, jpjm1 ! not jpj-row 467 DO ji = 1, jpi 468 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj) & 469 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 470 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 471 END DO 472 END DO 445 DO_2D_11_10 446 zhup2_e(ji,jj) = hu_0(ji,jj) + r1_2 * r1_e1e2u(ji,jj) & 447 & * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 448 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 449 END_2D 450 DO_2D_10_11 451 zhvp2_e(ji,jj) = hv_0(ji,jj) + r1_2 * r1_e1e2v(ji,jj) & 452 & * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 453 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 454 END_2D 473 455 ! 474 456 ENDIF … … 526 508 !-- ssh = ssh - delta_t' * [ frc + div( flux ) ] --! 527 509 !-------------------------------------------------------------------------! 528 DO jj = 2, jpjm1 ! INNER domain 529 DO ji = 2, jpim1 530 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 531 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rdtbt * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) 532 END DO 533 END DO 510 DO_2D_00_00 511 zhdiv = ( zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1) ) * r1_e1e2t(ji,jj) 512 ssha_e(ji,jj) = ( sshn_e(ji,jj) - rdtbt * ( zssh_frc(ji,jj) + zhdiv ) ) * ssmask(ji,jj) 513 END_2D 534 514 ! 535 515 CALL lbc_lnk_multi( 'dynspg_ts', ssha_e, 'T', 1._wp, zhU, 'U', -1._wp, zhV, 'V', -1._wp ) … … 553 533 ! Sea Surface Height at u-,v-points (vvl case only) 554 534 IF( .NOT.ln_linssh ) THEN 555 DO jj = 2, jpjm1 ! INNER domain, will be extended to whole domain later 556 DO ji = 2, jpim1 ! NO Vector Opt. 557 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 558 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 559 & + e1e2t(ji+1,jj ) * ssha_e(ji+1,jj ) ) 560 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 561 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 562 & + e1e2t(ji ,jj+1) * ssha_e(ji ,jj+1) ) 563 END DO 564 END DO 535 DO_2D_00_00 536 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 537 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 538 & + e1e2t(ji+1,jj ) * ssha_e(ji+1,jj ) ) 539 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 540 & * ( e1e2t(ji ,jj ) * ssha_e(ji ,jj ) & 541 & + e1e2t(ji ,jj+1) * ssha_e(ji ,jj+1) ) 542 END_2D 565 543 ENDIF 566 544 ! … … 575 553 ! ! Surface pressure gradient 576 554 zldg = ( 1._wp - rn_scal_load ) * grav ! local factor 577 DO jj = 2, jpjm1 578 DO ji = 2, jpim1 579 zu_spg(ji,jj) = - zldg * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 580 zv_spg(ji,jj) = - zldg * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) 581 END DO 582 END DO 555 DO_2D_00_00 556 zu_spg(ji,jj) = - zldg * ( zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj) ) * r1_e1u(ji,jj) 557 zv_spg(ji,jj) = - zldg * ( zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj) ) * r1_e2v(ji,jj) 558 END_2D 583 559 IF( ln_wd_il ) THEN ! W/D : gravity filters applied on pressure gradient 584 560 CALL wad_spg( zsshp2_e, zcpx, zcpy ) ! Calculating W/D gravity filters … … 595 571 ! Add tidal astronomical forcing if defined 596 572 IF ( ln_tide .AND. ln_tide_pot ) THEN 597 DO jj = 2, jpjm1 598 DO ji = fs_2, fs_jpim1 ! vector opt. 599 zu_trd(ji,jj) = zu_trd(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 600 zv_trd(ji,jj) = zv_trd(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) 601 END DO 602 END DO 573 DO_2D_00_00 574 zu_trd(ji,jj) = zu_trd(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 575 zv_trd(ji,jj) = zv_trd(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) 576 END_2D 603 577 ENDIF 604 578 ! … … 606 580 !jth do implicitly instead 607 581 IF ( .NOT. ll_wd ) THEN ! Revert to explicit for bit comparison tests in non wad runs 608 DO jj = 2, jpjm1 609 DO ji = fs_2, fs_jpim1 ! vector opt. 610 zu_trd(ji,jj) = zu_trd(ji,jj) + zCdU_u(ji,jj) * un_e(ji,jj) * hur_e(ji,jj) 611 zv_trd(ji,jj) = zv_trd(ji,jj) + zCdU_v(ji,jj) * vn_e(ji,jj) * hvr_e(ji,jj) 612 END DO 613 END DO 582 DO_2D_00_00 583 zu_trd(ji,jj) = zu_trd(ji,jj) + zCdU_u(ji,jj) * un_e(ji,jj) * hur_e(ji,jj) 584 zv_trd(ji,jj) = zv_trd(ji,jj) + zCdU_v(ji,jj) * vn_e(ji,jj) * hvr_e(ji,jj) 585 END_2D 614 586 ENDIF 615 587 ! … … 626 598 !------------------------------------------------------------------------------------------------------------------------! 627 599 IF( ln_dynadv_vec .OR. ln_linssh ) THEN !* Vector form 628 DO jj = 2, jpjm1 629 DO ji = fs_2, fs_jpim1 ! vector opt. 630 ua_e(ji,jj) = ( un_e(ji,jj) & 631 & + rdtbt * ( zu_spg(ji,jj) & 632 & + zu_trd(ji,jj) & 633 & + zu_frc(ji,jj) ) & 634 & ) * ssumask(ji,jj) 635 636 va_e(ji,jj) = ( vn_e(ji,jj) & 637 & + rdtbt * ( zv_spg(ji,jj) & 638 & + zv_trd(ji,jj) & 639 & + zv_frc(ji,jj) ) & 640 & ) * ssvmask(ji,jj) 641 END DO 642 END DO 600 DO_2D_00_00 601 ua_e(ji,jj) = ( un_e(ji,jj) & 602 & + rdtbt * ( zu_spg(ji,jj) & 603 & + zu_trd(ji,jj) & 604 & + zu_frc(ji,jj) ) & 605 & ) * ssumask(ji,jj) 606 607 va_e(ji,jj) = ( vn_e(ji,jj) & 608 & + rdtbt * ( zv_spg(ji,jj) & 609 & + zv_trd(ji,jj) & 610 & + zv_frc(ji,jj) ) & 611 & ) * ssvmask(ji,jj) 612 END_2D 643 613 ! 644 614 ELSE !* Flux form 645 DO jj = 2, jpjm1 646 DO ji = 2, jpim1 647 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 648 ! ! backward interpolated depth used in spg terms at jn+1/2 649 zhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 650 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 651 zhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 652 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 653 ! ! inverse depth at jn+1 654 z1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) 655 z1_hv = ssvmask(ji,jj) / ( hv_0(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj) ) 656 ! 657 ua_e(ji,jj) = ( hu_e (ji,jj) * un_e (ji,jj) & 658 & + rdtbt * ( zhu_bck * zu_spg (ji,jj) & ! 659 & + zhup2_e(ji,jj) * zu_trd (ji,jj) & ! 660 & + hu_n (ji,jj) * zu_frc (ji,jj) ) ) * z1_hu 661 ! 662 va_e(ji,jj) = ( hv_e (ji,jj) * vn_e (ji,jj) & 663 & + rdtbt * ( zhv_bck * zv_spg (ji,jj) & ! 664 & + zhvp2_e(ji,jj) * zv_trd (ji,jj) & ! 665 & + hv_n (ji,jj) * zv_frc (ji,jj) ) ) * z1_hv 666 END DO 667 END DO 615 DO_2D_00_00 616 ! ! hu_e, hv_e hold depth at jn, zhup2_e, zhvp2_e hold extrapolated depth at jn+1/2 617 ! ! backward interpolated depth used in spg terms at jn+1/2 618 zhu_bck = hu_0(ji,jj) + r1_2*r1_e1e2u(ji,jj) * ( e1e2t(ji ,jj) * zsshp2_e(ji ,jj) & 619 & + e1e2t(ji+1,jj) * zsshp2_e(ji+1,jj) ) * ssumask(ji,jj) 620 zhv_bck = hv_0(ji,jj) + r1_2*r1_e1e2v(ji,jj) * ( e1e2t(ji,jj ) * zsshp2_e(ji,jj ) & 621 & + e1e2t(ji,jj+1) * zsshp2_e(ji,jj+1) ) * ssvmask(ji,jj) 622 ! ! inverse depth at jn+1 623 z1_hu = ssumask(ji,jj) / ( hu_0(ji,jj) + zsshu_a(ji,jj) + 1._wp - ssumask(ji,jj) ) 624 z1_hv = ssvmask(ji,jj) / ( hv_0(ji,jj) + zsshv_a(ji,jj) + 1._wp - ssvmask(ji,jj) ) 625 ! 626 ua_e(ji,jj) = ( hu_e (ji,jj) * un_e (ji,jj) & 627 & + rdtbt * ( zhu_bck * zu_spg (ji,jj) & ! 628 & + zhup2_e(ji,jj) * zu_trd (ji,jj) & ! 629 & + hu_n (ji,jj) * zu_frc (ji,jj) ) ) * z1_hu 630 ! 631 va_e(ji,jj) = ( hv_e (ji,jj) * vn_e (ji,jj) & 632 & + rdtbt * ( zhv_bck * zv_spg (ji,jj) & ! 633 & + zhvp2_e(ji,jj) * zv_trd (ji,jj) & ! 634 & + hv_n (ji,jj) * zv_frc (ji,jj) ) ) * z1_hv 635 END_2D 668 636 ENDIF 669 637 !jth implicit bottom friction: 670 638 IF ( ll_wd ) THEN ! revert to explicit for bit comparison tests in non wad runs 671 DO jj = 2, jpjm1 672 DO ji = fs_2, fs_jpim1 ! vector opt. 673 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rdtbt * zCdU_u(ji,jj) * hur_e(ji,jj)) 674 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rdtbt * zCdU_v(ji,jj) * hvr_e(ji,jj)) 675 END DO 676 END DO 639 DO_2D_00_00 640 ua_e(ji,jj) = ua_e(ji,jj) /(1.0 - rdtbt * zCdU_u(ji,jj) * hur_e(ji,jj)) 641 va_e(ji,jj) = va_e(ji,jj) /(1.0 - rdtbt * zCdU_v(ji,jj) * hvr_e(ji,jj)) 642 END_2D 677 643 ENDIF 678 644 … … 737 703 IF (ln_bt_fw) THEN 738 704 IF( .NOT.( kt == nit000 .AND. neuler==0 ) ) THEN 739 DO jj = 1, jpj 740 DO ji = 1, jpi 741 zun_save = un_adv(ji,jj) 742 zvn_save = vn_adv(ji,jj) 743 ! ! apply the previously computed correction 744 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - atfp * un_bf(ji,jj) ) 745 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - atfp * vn_bf(ji,jj) ) 746 ! ! Update corrective fluxes for next time step 747 un_bf(ji,jj) = atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) ) 748 vn_bf(ji,jj) = atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) ) 749 ! ! Save integrated transport for next computation 750 ub2_b(ji,jj) = zun_save 751 vb2_b(ji,jj) = zvn_save 752 END DO 753 END DO 705 DO_2D_11_11 706 zun_save = un_adv(ji,jj) 707 zvn_save = vn_adv(ji,jj) 708 ! ! apply the previously computed correction 709 un_adv(ji,jj) = r1_2 * ( ub2_b(ji,jj) + zun_save - atfp * un_bf(ji,jj) ) 710 vn_adv(ji,jj) = r1_2 * ( vb2_b(ji,jj) + zvn_save - atfp * vn_bf(ji,jj) ) 711 ! ! Update corrective fluxes for next time step 712 un_bf(ji,jj) = atfp * un_bf(ji,jj) + ( zun_save - ub2_b(ji,jj) ) 713 vn_bf(ji,jj) = atfp * vn_bf(ji,jj) + ( zvn_save - vb2_b(ji,jj) ) 714 ! ! Save integrated transport for next computation 715 ub2_b(ji,jj) = zun_save 716 vb2_b(ji,jj) = zvn_save 717 END_2D 754 718 ELSE 755 719 un_bf(:,:) = 0._wp ! corrective fluxes for next time step set to zero … … 770 734 ELSE 771 735 ! At this stage, ssha has been corrected: compute new depths at velocity points 772 DO jj = 1, jpjm1 773 DO ji = 1, jpim1 ! NO Vector Opt. 774 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 775 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 776 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 777 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 778 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 779 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 780 END DO 781 END DO 736 DO_2D_10_10 737 zsshu_a(ji,jj) = r1_2 * ssumask(ji,jj) * r1_e1e2u(ji,jj) & 738 & * ( e1e2t(ji ,jj) * ssha(ji ,jj) & 739 & + e1e2t(ji+1,jj) * ssha(ji+1,jj) ) 740 zsshv_a(ji,jj) = r1_2 * ssvmask(ji,jj) * r1_e1e2v(ji,jj) & 741 & * ( e1e2t(ji,jj ) * ssha(ji,jj ) & 742 & + e1e2t(ji,jj+1) * ssha(ji,jj+1) ) 743 END_2D 782 744 CALL lbc_lnk_multi( 'dynspg_ts', zsshu_a, 'U', 1._wp, zsshv_a, 'V', 1._wp ) ! Boundary conditions 783 745 ! … … 794 756 ! Correct velocities so that the barotropic velocity equals (un_adv, vn_adv) (in all cases) 795 757 DO jk = 1, jpkm1 796 u n(:,:,jk) = ( un(:,:,jk) + un_adv(:,:)*r1_hu_n(:,:) - un_b(:,:) ) * umask(:,:,jk)797 v n(:,:,jk) = ( vn(:,:,jk) + vn_adv(:,:)*r1_hv_n(:,:) - vn_b(:,:) ) * vmask(:,:,jk)758 uu(:,:,jk,Nii) = ( uu(:,:,jk,Nii) + un_adv(:,:)*r1_hu_n(:,:) - un_b(:,:) ) * umask(:,:,jk) 759 vv(:,:,jk,Nii) = ( vv(:,:,jk,Nii) + vn_adv(:,:)*r1_hv_n(:,:) - vn_b(:,:) ) * vmask(:,:,jk) 798 760 END DO 799 761 … … 802 764 CALL lbc_lnk_multi( 'dynspg_ts', zuwdav2, 'U', 1._wp, zvwdav2, 'V', 1._wp) 803 765 DO jk = 1, jpkm1 804 u n(:,:,jk) = ( un_adv(:,:)*r1_hu_n(:,:) &805 & + zuwdav2(:,:)*(u n(:,:,jk) - un_adv(:,:)*r1_hu_n(:,:)) ) * umask(:,:,jk)806 v n(:,:,jk) = ( vn_adv(:,:)*r1_hv_n(:,:) &807 & + zvwdav2(:,:)*(v n(:,:,jk) - vn_adv(:,:)*r1_hv_n(:,:)) ) * vmask(:,:,jk)766 uu(:,:,jk,Nii) = ( un_adv(:,:)*r1_hu_n(:,:) & 767 & + zuwdav2(:,:)*(uu(:,:,jk,Nii) - un_adv(:,:)*r1_hu_n(:,:)) ) * umask(:,:,jk) 768 vv(:,:,jk,Nii) = ( vn_adv(:,:)*r1_hv_n(:,:) & 769 & + zvwdav2(:,:)*(vv(:,:,jk,Nii) - vn_adv(:,:)*r1_hv_n(:,:)) ) * vmask(:,:,jk) 808 770 END DO 809 771 END IF … … 1007 969 ! Max courant number for ext. grav. waves 1008 970 ! 1009 DO jj = 1, jpj 1010 DO ji =1, jpi 1011 zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 1012 zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) 1013 zcu(ji,jj) = SQRT( grav * MAX(ht_0(ji,jj),0._wp) * (zxr2 + zyr2) ) 1014 END DO 1015 END DO 971 DO_2D_11_11 972 zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 973 zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) 974 zcu(ji,jj) = SQRT( grav * MAX(ht_0(ji,jj),0._wp) * (zxr2 + zyr2) ) 975 END_2D 1016 976 ! 1017 977 zcmax = MAXVAL( zcu(:,:) ) … … 1133 1093 SELECT CASE( nn_een_e3f ) !* ff_f/e3 at F-point 1134 1094 CASE ( 0 ) ! original formulation (masked averaging of e3t divided by 4) 1135 DO jj = 1, jpjm1 1136 DO ji = 1, jpim1 1137 zwz(ji,jj) = ( ht_n(ji ,jj+1) + ht_n(ji+1,jj+1) + & 1138 & ht_n(ji ,jj ) + ht_n(ji+1,jj ) ) * 0.25_wp 1139 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1140 END DO 1141 END DO 1095 DO_2D_10_10 1096 zwz(ji,jj) = ( ht_n(ji ,jj+1) + ht_n(ji+1,jj+1) + & 1097 & ht_n(ji ,jj ) + ht_n(ji+1,jj ) ) * 0.25_wp 1098 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1099 END_2D 1142 1100 CASE ( 1 ) ! new formulation (masked averaging of e3t divided by the sum of mask) 1143 DO jj = 1, jpjm1 1144 DO ji = 1, jpim1 1145 zwz(ji,jj) = ( ht_n (ji ,jj+1) + ht_n (ji+1,jj+1) & 1146 & + ht_n (ji ,jj ) + ht_n (ji+1,jj ) ) & 1147 & / ( MAX( 1._wp, ssmask(ji ,jj+1) + ssmask(ji+1,jj+1) & 1148 & + ssmask(ji ,jj ) + ssmask(ji+1,jj ) ) ) 1149 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1150 END DO 1151 END DO 1101 DO_2D_10_10 1102 zwz(ji,jj) = ( ht_n (ji ,jj+1) + ht_n (ji+1,jj+1) & 1103 & + ht_n (ji ,jj ) + ht_n (ji+1,jj ) ) & 1104 & / ( MAX( 1._wp, ssmask(ji ,jj+1) + ssmask(ji+1,jj+1) & 1105 & + ssmask(ji ,jj ) + ssmask(ji+1,jj ) ) ) 1106 IF( zwz(ji,jj) /= 0._wp ) zwz(ji,jj) = ff_f(ji,jj) / zwz(ji,jj) 1107 END_2D 1152 1108 END SELECT 1153 1109 CALL lbc_lnk( 'dynspg_ts', zwz, 'F', 1._wp ) 1154 1110 ! 1155 1111 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1156 DO jj = 2, jpj 1157 DO ji = 2, jpi 1158 ftne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) 1159 ftnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj ) 1160 ftse(ji,jj) = zwz(ji ,jj ) + zwz(ji ,jj-1) + zwz(ji-1,jj-1) 1161 ftsw(ji,jj) = zwz(ji ,jj-1) + zwz(ji-1,jj-1) + zwz(ji-1,jj ) 1162 END DO 1163 END DO 1112 DO_2D_01_01 1113 ftne(ji,jj) = zwz(ji-1,jj ) + zwz(ji ,jj ) + zwz(ji ,jj-1) 1114 ftnw(ji,jj) = zwz(ji-1,jj-1) + zwz(ji-1,jj ) + zwz(ji ,jj ) 1115 ftse(ji,jj) = zwz(ji ,jj ) + zwz(ji ,jj-1) + zwz(ji-1,jj-1) 1116 ftsw(ji,jj) = zwz(ji ,jj-1) + zwz(ji-1,jj-1) + zwz(ji-1,jj ) 1117 END_2D 1164 1118 ! 1165 1119 CASE( np_EET ) != EEN scheme using e3t (energy conserving scheme) 1166 1120 ftne(1,:) = 0._wp ; ftnw(1,:) = 0._wp ; ftse(1,:) = 0._wp ; ftsw(1,:) = 0._wp 1167 DO jj = 2, jpj 1168 DO ji = 2, jpi 1169 z1_ht = ssmask(ji,jj) / ( ht_n(ji,jj) + 1._wp - ssmask(ji,jj) ) 1170 ftne(ji,jj) = ( ff_f(ji-1,jj ) + ff_f(ji ,jj ) + ff_f(ji ,jj-1) ) * z1_ht 1171 ftnw(ji,jj) = ( ff_f(ji-1,jj-1) + ff_f(ji-1,jj ) + ff_f(ji ,jj ) ) * z1_ht 1172 ftse(ji,jj) = ( ff_f(ji ,jj ) + ff_f(ji ,jj-1) + ff_f(ji-1,jj-1) ) * z1_ht 1173 ftsw(ji,jj) = ( ff_f(ji ,jj-1) + ff_f(ji-1,jj-1) + ff_f(ji-1,jj ) ) * z1_ht 1174 END DO 1175 END DO 1121 DO_2D_01_01 1122 z1_ht = ssmask(ji,jj) / ( ht_n(ji,jj) + 1._wp - ssmask(ji,jj) ) 1123 ftne(ji,jj) = ( ff_f(ji-1,jj ) + ff_f(ji ,jj ) + ff_f(ji ,jj-1) ) * z1_ht 1124 ftnw(ji,jj) = ( ff_f(ji-1,jj-1) + ff_f(ji-1,jj ) + ff_f(ji ,jj ) ) * z1_ht 1125 ftse(ji,jj) = ( ff_f(ji ,jj ) + ff_f(ji ,jj-1) + ff_f(ji-1,jj-1) ) * z1_ht 1126 ftsw(ji,jj) = ( ff_f(ji ,jj-1) + ff_f(ji-1,jj-1) + ff_f(ji-1,jj ) ) * z1_ht 1127 END_2D 1176 1128 ! 1177 1129 CASE( np_ENE, np_ENS , np_MIX ) != all other schemes (ENE, ENS, MIX) except ENT ! … … 1200 1152 ! 1201 1153 !zhf(:,:) = hbatf(:,:) 1202 DO jj = 1, jpjm1 1203 DO ji = 1, jpim1 1204 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) & 1205 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) & 1206 & / MAX( ssmask(ji,jj ) + ssmask(ji+1,jj ) & 1207 & + ssmask(ji,jj+1) + ssmask(ji+1,jj+1) , 1._wp ) 1208 END DO 1209 END DO 1154 DO_2D_10_10 1155 zhf(ji,jj) = ( ht_0 (ji,jj ) + ht_0 (ji+1,jj ) & 1156 & + ht_0 (ji,jj+1) + ht_0 (ji+1,jj+1) ) & 1157 & / MAX( ssmask(ji,jj ) + ssmask(ji+1,jj ) & 1158 & + ssmask(ji,jj+1) + ssmask(ji+1,jj+1) , 1._wp ) 1159 END_2D 1210 1160 ENDIF 1211 1161 ! … … 1221 1171 CALL lbc_lnk( 'dynspg_ts', zhf, 'F', 1._wp ) 1222 1172 ! JC: TBC. hf should be greater than 0 1223 DO jj = 1, jpj 1224 DO ji = 1, jpi 1225 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj) 1226 END DO 1227 END DO 1173 DO_2D_11_11 1174 IF( zhf(ji,jj) /= 0._wp ) zwz(ji,jj) = 1._wp / zhf(ji,jj) 1175 END_2D 1228 1176 zwz(:,:) = ff_f(:,:) * zwz(:,:) 1229 1177 END SELECT … … 1246 1194 SELECT CASE( nvor_scheme ) 1247 1195 CASE( np_ENT ) ! enstrophy conserving scheme (f-point) 1248 DO jj = 2, jpjm1 1249 DO ji = 2, jpim1 1250 z1_hu = ssumask(ji,jj) / ( hu_n(ji,jj) + 1._wp - ssumask(ji,jj) ) 1251 z1_hv = ssvmask(ji,jj) / ( hv_n(ji,jj) + 1._wp - ssvmask(ji,jj) ) 1252 zu_trd(ji,jj) = + r1_4 * r1_e1e2u(ji,jj) * z1_hu & 1253 & * ( e1e2t(ji+1,jj)*ht_n(ji+1,jj)*ff_t(ji+1,jj) * ( vn_b(ji+1,jj) + vn_b(ji+1,jj-1) ) & 1254 & + e1e2t(ji ,jj)*ht_n(ji ,jj)*ff_t(ji ,jj) * ( vn_b(ji ,jj) + vn_b(ji ,jj-1) ) ) 1255 ! 1256 zv_trd(ji,jj) = - r1_4 * r1_e1e2v(ji,jj) * z1_hv & 1257 & * ( e1e2t(ji,jj+1)*ht_n(ji,jj+1)*ff_t(ji,jj+1) * ( un_b(ji,jj+1) + un_b(ji-1,jj+1) ) & 1258 & + e1e2t(ji,jj )*ht_n(ji,jj )*ff_t(ji,jj ) * ( un_b(ji,jj ) + un_b(ji-1,jj ) ) ) 1259 END DO 1260 END DO 1196 DO_2D_00_00 1197 z1_hu = ssumask(ji,jj) / ( hu_n(ji,jj) + 1._wp - ssumask(ji,jj) ) 1198 z1_hv = ssvmask(ji,jj) / ( hv_n(ji,jj) + 1._wp - ssvmask(ji,jj) ) 1199 zu_trd(ji,jj) = + r1_4 * r1_e1e2u(ji,jj) * z1_hu & 1200 & * ( e1e2t(ji+1,jj)*ht_n(ji+1,jj)*ff_t(ji+1,jj) * ( vn_b(ji+1,jj) + vn_b(ji+1,jj-1) ) & 1201 & + e1e2t(ji ,jj)*ht_n(ji ,jj)*ff_t(ji ,jj) * ( vn_b(ji ,jj) + vn_b(ji ,jj-1) ) ) 1202 ! 1203 zv_trd(ji,jj) = - r1_4 * r1_e1e2v(ji,jj) * z1_hv & 1204 & * ( e1e2t(ji,jj+1)*ht_n(ji,jj+1)*ff_t(ji,jj+1) * ( un_b(ji,jj+1) + un_b(ji-1,jj+1) ) & 1205 & + e1e2t(ji,jj )*ht_n(ji,jj )*ff_t(ji,jj ) * ( un_b(ji,jj ) + un_b(ji-1,jj ) ) ) 1206 END_2D 1261 1207 ! 1262 1208 CASE( np_ENE , np_MIX ) ! energy conserving scheme (t-point) ENE or MIX 1263 DO jj = 2, jpjm1 1264 DO ji = fs_2, fs_jpim1 ! vector opt. 1265 zy1 = ( zhV(ji,jj-1) + zhV(ji+1,jj-1) ) * r1_e1u(ji,jj) 1266 zy2 = ( zhV(ji,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) 1267 zx1 = ( zhU(ji-1,jj) + zhU(ji-1,jj+1) ) * r1_e2v(ji,jj) 1268 zx2 = ( zhU(ji ,jj) + zhU(ji ,jj+1) ) * r1_e2v(ji,jj) 1269 ! energy conserving formulation for planetary vorticity term 1270 zu_trd(ji,jj) = r1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 1271 zv_trd(ji,jj) = - r1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 ) 1272 END DO 1273 END DO 1209 DO_2D_00_00 1210 zy1 = ( zhV(ji,jj-1) + zhV(ji+1,jj-1) ) * r1_e1u(ji,jj) 1211 zy2 = ( zhV(ji,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) 1212 zx1 = ( zhU(ji-1,jj) + zhU(ji-1,jj+1) ) * r1_e2v(ji,jj) 1213 zx2 = ( zhU(ji ,jj) + zhU(ji ,jj+1) ) * r1_e2v(ji,jj) 1214 ! energy conserving formulation for planetary vorticity term 1215 zu_trd(ji,jj) = r1_4 * ( zwz(ji ,jj-1) * zy1 + zwz(ji,jj) * zy2 ) 1216 zv_trd(ji,jj) = - r1_4 * ( zwz(ji-1,jj ) * zx1 + zwz(ji,jj) * zx2 ) 1217 END_2D 1274 1218 ! 1275 1219 CASE( np_ENS ) ! enstrophy conserving scheme (f-point) 1276 DO jj = 2, jpjm1 1277 DO ji = fs_2, fs_jpim1 ! vector opt. 1278 zy1 = r1_8 * ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) & 1279 & + zhV(ji ,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) 1280 zx1 = - r1_8 * ( zhU(ji-1,jj ) + zhU(ji-1,jj+1) & 1281 & + zhU(ji ,jj ) + zhU(ji ,jj+1) ) * r1_e2v(ji,jj) 1282 zu_trd(ji,jj) = zy1 * ( zwz(ji ,jj-1) + zwz(ji,jj) ) 1283 zv_trd(ji,jj) = zx1 * ( zwz(ji-1,jj ) + zwz(ji,jj) ) 1284 END DO 1285 END DO 1220 DO_2D_00_00 1221 zy1 = r1_8 * ( zhV(ji ,jj-1) + zhV(ji+1,jj-1) & 1222 & + zhV(ji ,jj ) + zhV(ji+1,jj ) ) * r1_e1u(ji,jj) 1223 zx1 = - r1_8 * ( zhU(ji-1,jj ) + zhU(ji-1,jj+1) & 1224 & + zhU(ji ,jj ) + zhU(ji ,jj+1) ) * r1_e2v(ji,jj) 1225 zu_trd(ji,jj) = zy1 * ( zwz(ji ,jj-1) + zwz(ji,jj) ) 1226 zv_trd(ji,jj) = zx1 * ( zwz(ji-1,jj ) + zwz(ji,jj) ) 1227 END_2D 1286 1228 ! 1287 1229 CASE( np_EET , np_EEN ) ! energy & enstrophy scheme (using e3t or e3f) 1288 DO jj = 2, jpjm1 1289 DO ji = fs_2, fs_jpim1 ! vector opt. 1290 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zhV(ji ,jj ) & 1291 & + ftnw(ji+1,jj) * zhV(ji+1,jj ) & 1292 & + ftse(ji,jj ) * zhV(ji ,jj-1) & 1293 & + ftsw(ji+1,jj) * zhV(ji+1,jj-1) ) 1294 zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zhU(ji-1,jj+1) & 1295 & + ftse(ji,jj+1) * zhU(ji ,jj+1) & 1296 & + ftnw(ji,jj ) * zhU(ji-1,jj ) & 1297 & + ftne(ji,jj ) * zhU(ji ,jj ) ) 1298 END DO 1299 END DO 1230 DO_2D_00_00 1231 zu_trd(ji,jj) = + r1_12 * r1_e1u(ji,jj) * ( ftne(ji,jj ) * zhV(ji ,jj ) & 1232 & + ftnw(ji+1,jj) * zhV(ji+1,jj ) & 1233 & + ftse(ji,jj ) * zhV(ji ,jj-1) & 1234 & + ftsw(ji+1,jj) * zhV(ji+1,jj-1) ) 1235 zv_trd(ji,jj) = - r1_12 * r1_e2v(ji,jj) * ( ftsw(ji,jj+1) * zhU(ji-1,jj+1) & 1236 & + ftse(ji,jj+1) * zhU(ji ,jj+1) & 1237 & + ftnw(ji,jj ) * zhU(ji-1,jj ) & 1238 & + ftne(ji,jj ) * zhU(ji ,jj ) ) 1239 END_2D 1300 1240 ! 1301 1241 END SELECT … … 1322 1262 ! 1323 1263 IF( ln_wd_dl_rmp ) THEN 1324 DO jj = 1, jpj 1325 DO ji = 1, jpi 1326 IF ( pssh(ji,jj) + ht_0(ji,jj) > 2._wp * rn_wdmin1 ) THEN 1327 ! IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin2 ) THEN 1328 ptmsk(ji,jj) = 1._wp 1329 ELSEIF( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN 1330 ptmsk(ji,jj) = TANH( 50._wp*( ( pssh(ji,jj) + ht_0(ji,jj) - rn_wdmin1 )*r_rn_wdmin1) ) 1331 ELSE 1332 ptmsk(ji,jj) = 0._wp 1333 ENDIF 1334 END DO 1335 END DO 1264 DO_2D_11_11 1265 IF ( pssh(ji,jj) + ht_0(ji,jj) > 2._wp * rn_wdmin1 ) THEN 1266 ! IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin2 ) THEN 1267 ptmsk(ji,jj) = 1._wp 1268 ELSEIF( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN 1269 ptmsk(ji,jj) = TANH( 50._wp*( ( pssh(ji,jj) + ht_0(ji,jj) - rn_wdmin1 )*r_rn_wdmin1) ) 1270 ELSE 1271 ptmsk(ji,jj) = 0._wp 1272 ENDIF 1273 END_2D 1336 1274 ELSE 1337 DO jj = 1, jpj 1338 DO ji = 1, jpi 1339 IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN ; ptmsk(ji,jj) = 1._wp 1340 ELSE ; ptmsk(ji,jj) = 0._wp 1341 ENDIF 1342 END DO 1343 END DO 1275 DO_2D_11_11 1276 IF ( pssh(ji,jj) + ht_0(ji,jj) > rn_wdmin1 ) THEN ; ptmsk(ji,jj) = 1._wp 1277 ELSE ; ptmsk(ji,jj) = 0._wp 1278 ENDIF 1279 END_2D 1344 1280 ENDIF 1345 1281 ! … … 1365 1301 !!---------------------------------------------------------------------- 1366 1302 ! 1367 DO jj = 1, jpj 1368 DO ji = 1, jpim1 ! not jpi-column 1369 IF ( phU(ji,jj) > 0._wp ) THEN ; pUmsk(ji,jj) = pTmsk(ji ,jj) 1370 ELSE ; pUmsk(ji,jj) = pTmsk(ji+1,jj) 1371 ENDIF 1372 phU(ji,jj) = pUmsk(ji,jj)*phU(ji,jj) 1373 pu (ji,jj) = pUmsk(ji,jj)*pu (ji,jj) 1374 END DO 1375 END DO 1376 ! 1377 DO jj = 1, jpjm1 ! not jpj-row 1378 DO ji = 1, jpi 1379 IF ( phV(ji,jj) > 0._wp ) THEN ; pVmsk(ji,jj) = pTmsk(ji,jj ) 1380 ELSE ; pVmsk(ji,jj) = pTmsk(ji,jj+1) 1381 ENDIF 1382 phV(ji,jj) = pVmsk(ji,jj)*phV(ji,jj) 1383 pv (ji,jj) = pVmsk(ji,jj)*pv (ji,jj) 1384 END DO 1385 END DO 1303 DO_2D_11_10 1304 IF ( phU(ji,jj) > 0._wp ) THEN ; pUmsk(ji,jj) = pTmsk(ji ,jj) 1305 ELSE ; pUmsk(ji,jj) = pTmsk(ji+1,jj) 1306 ENDIF 1307 phU(ji,jj) = pUmsk(ji,jj)*phU(ji,jj) 1308 pu (ji,jj) = pUmsk(ji,jj)*pu (ji,jj) 1309 END_2D 1310 ! 1311 DO_2D_10_11 1312 IF ( phV(ji,jj) > 0._wp ) THEN ; pVmsk(ji,jj) = pTmsk(ji,jj ) 1313 ELSE ; pVmsk(ji,jj) = pTmsk(ji,jj+1) 1314 ENDIF 1315 phV(ji,jj) = pVmsk(ji,jj)*phV(ji,jj) 1316 pv (ji,jj) = pVmsk(ji,jj)*pv (ji,jj) 1317 END_2D 1386 1318 ! 1387 1319 END SUBROUTINE wad_Umsk … … 1399 1331 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: zcpx, zcpy 1400 1332 !!---------------------------------------------------------------------- 1401 DO jj = 2, jpjm1 1402 DO ji = 2, jpim1 1403 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1404 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) .AND. & 1405 & MAX( sshn(ji,jj) + ht_0(ji,jj) , sshn(ji+1,jj) + ht_0(ji+1,jj) ) & 1406 & > rn_wdmin1 + rn_wdmin2 1407 ll_tmp2 = ( ABS( sshn(ji+1,jj) - sshn(ji ,jj)) > 1.E-12 ).AND.( & 1408 & MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1409 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 1410 IF(ll_tmp1) THEN 1411 zcpx(ji,jj) = 1.0_wp 1412 ELSEIF(ll_tmp2) THEN 1413 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 1414 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) & 1415 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 1416 zcpx(ji,jj) = max(min( zcpx(ji,jj) , 1.0_wp),0.0_wp) 1417 ELSE 1418 zcpx(ji,jj) = 0._wp 1419 ENDIF 1420 ! 1421 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1422 & MAX( -ht_0(ji,jj) , -ht_0(ji,jj+1) ) .AND. & 1423 & MAX( sshn(ji,jj) + ht_0(ji,jj) , sshn(ji,jj+1) + ht_0(ji,jj+1) ) & 1424 & > rn_wdmin1 + rn_wdmin2 1425 ll_tmp2 = ( ABS( sshn(ji,jj) - sshn(ji,jj+1)) > 1.E-12 ).AND.( & 1426 & MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1427 & MAX( -ht_0(ji,jj) , -ht_0(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 ) 1428 1429 IF(ll_tmp1) THEN 1430 zcpy(ji,jj) = 1.0_wp 1431 ELSE IF(ll_tmp2) THEN 1432 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 1433 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) & 1434 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 1435 zcpy(ji,jj) = MAX( 0._wp , MIN( zcpy(ji,jj) , 1.0_wp ) ) 1436 ELSE 1437 zcpy(ji,jj) = 0._wp 1438 ENDIF 1439 END DO 1440 END DO 1333 DO_2D_00_00 1334 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1335 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) .AND. & 1336 & MAX( sshn(ji,jj) + ht_0(ji,jj) , sshn(ji+1,jj) + ht_0(ji+1,jj) ) & 1337 & > rn_wdmin1 + rn_wdmin2 1338 ll_tmp2 = ( ABS( sshn(ji+1,jj) - sshn(ji ,jj)) > 1.E-12 ).AND.( & 1339 & MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1340 & MAX( -ht_0(ji,jj) , -ht_0(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 ) 1341 IF(ll_tmp1) THEN 1342 zcpx(ji,jj) = 1.0_wp 1343 ELSEIF(ll_tmp2) THEN 1344 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 1345 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) & 1346 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 1347 zcpx(ji,jj) = max(min( zcpx(ji,jj) , 1.0_wp),0.0_wp) 1348 ELSE 1349 zcpx(ji,jj) = 0._wp 1350 ENDIF 1351 ! 1352 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1353 & MAX( -ht_0(ji,jj) , -ht_0(ji,jj+1) ) .AND. & 1354 & MAX( sshn(ji,jj) + ht_0(ji,jj) , sshn(ji,jj+1) + ht_0(ji,jj+1) ) & 1355 & > rn_wdmin1 + rn_wdmin2 1356 ll_tmp2 = ( ABS( sshn(ji,jj) - sshn(ji,jj+1)) > 1.E-12 ).AND.( & 1357 & MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1358 & MAX( -ht_0(ji,jj) , -ht_0(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 ) 1359 1360 IF(ll_tmp1) THEN 1361 zcpy(ji,jj) = 1.0_wp 1362 ELSE IF(ll_tmp2) THEN 1363 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 1364 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) & 1365 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 1366 zcpy(ji,jj) = MAX( 0._wp , MIN( zcpy(ji,jj) , 1.0_wp ) ) 1367 ELSE 1368 zcpy(ji,jj) = 0._wp 1369 ENDIF 1370 END_2D 1441 1371 1442 1372 END SUBROUTINE wad_spg … … 1467 1397 IF( ln_isfcav.OR.ln_drgice_imp ) THEN ! top+bottom friction (ocean cavities) 1468 1398 1469 DO jj = 2, jpjm1 1470 DO ji = 2, jpim1 ! INNER domain 1471 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) + rCdU_top(ji+1,jj)+rCdU_top(ji,jj) ) 1472 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) + rCdU_top(ji,jj+1)+rCdU_top(ji,jj) ) 1473 END DO 1474 END DO 1399 DO_2D_00_00 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 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 END_2D 1475 1403 ELSE ! bottom friction only 1476 DO jj = 2, jpjm1 1477 DO ji = 2, jpim1 ! INNER domain 1478 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 1479 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 1480 END DO 1481 END DO 1404 DO_2D_00_00 1405 pCdU_u(ji,jj) = r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) 1406 pCdU_v(ji,jj) = r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) 1407 END_2D 1482 1408 ENDIF 1483 1409 ! … … 1486 1412 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW bottom baroclinic velocities 1487 1413 1488 DO jj = 2, jpjm1 1489 DO ji = 2, jpim1 ! INNER domain 1490 ikbu = mbku(ji,jj) 1491 ikbv = mbkv(ji,jj) 1492 zu_i(ji,jj) = un(ji,jj,ikbu) - un_b(ji,jj) 1493 zv_i(ji,jj) = vn(ji,jj,ikbv) - vn_b(ji,jj) 1494 END DO 1495 END DO 1414 DO_2D_00_00 1415 ikbu = mbku(ji,jj) 1416 ikbv = mbkv(ji,jj) 1417 zu_i(ji,jj) = uu(ji,jj,ikbu,Nii) - un_b(ji,jj) 1418 zv_i(ji,jj) = vv(ji,jj,ikbv,Nii) - vn_b(ji,jj) 1419 END_2D 1496 1420 ELSE ! CENTRED integration: use BEFORE bottom baroclinic velocities 1497 1421 1498 DO jj = 2, jpjm1 1499 DO ji = 2, jpim1 ! INNER domain 1500 ikbu = mbku(ji,jj) 1501 ikbv = mbkv(ji,jj) 1502 zu_i(ji,jj) = ub(ji,jj,ikbu) - ub_b(ji,jj) 1503 zv_i(ji,jj) = vb(ji,jj,ikbv) - vb_b(ji,jj) 1504 END DO 1505 END DO 1422 DO_2D_00_00 1423 ikbu = mbku(ji,jj) 1424 ikbv = mbkv(ji,jj) 1425 zu_i(ji,jj) = uu(ji,jj,ikbu,Nnn) - ub_b(ji,jj) 1426 zv_i(ji,jj) = vv(ji,jj,ikbv,Nnn) - vb_b(ji,jj) 1427 END_2D 1506 1428 ENDIF 1507 1429 ! 1508 1430 IF( ln_wd_il ) THEN ! W/D : use the "clipped" bottom friction !!gm explain WHY, please ! 1509 1431 zztmp = -1._wp / rdtbt 1510 DO jj = 2, jpjm1 1511 DO ji = 2, jpim1 ! INNER domain 1512 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( & 1513 & r1_hu_n(ji,jj) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) 1514 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + zv_i(ji,jj) * wdrampv(ji,jj) * MAX( & 1515 & r1_hv_n(ji,jj) * r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) 1516 END DO 1517 END DO 1432 DO_2D_00_00 1433 pu_RHSi(ji,jj) = pu_RHSi(ji,jj) + zu_i(ji,jj) * wdrampu(ji,jj) * MAX( & 1434 & r1_hu_n(ji,jj) * r1_2*( rCdU_bot(ji+1,jj)+rCdU_bot(ji,jj) ) , zztmp ) 1435 pv_RHSi(ji,jj) = pv_RHSi(ji,jj) + zv_i(ji,jj) * wdrampv(ji,jj) * MAX( & 1436 & r1_hv_n(ji,jj) * r1_2*( rCdU_bot(ji,jj+1)+rCdU_bot(ji,jj) ) , zztmp ) 1437 END_2D 1518 1438 ELSE ! use "unclipped" drag (even if explicit friction is used in 3D calculation) 1519 1439 1520 DO jj = 2, jpjm1 1521 DO ji = 2, jpim1 ! INNER domain 1522 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) 1523 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) 1524 END DO 1525 END DO 1440 DO_2D_00_00 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 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) 1443 END_2D 1526 1444 END IF 1527 1445 ! … … 1532 1450 IF( ln_bt_fw ) THEN ! FORWARD integration: use NOW top baroclinic velocity 1533 1451 1534 DO jj = 2, jpjm1 1535 DO ji = 2, jpim1 ! INNER domain 1536 iktu = miku(ji,jj) 1537 iktv = mikv(ji,jj) 1538 zu_i(ji,jj) = un(ji,jj,iktu) - un_b(ji,jj) 1539 zv_i(ji,jj) = vn(ji,jj,iktv) - vn_b(ji,jj) 1540 END DO 1541 END DO 1452 DO_2D_00_00 1453 iktu = miku(ji,jj) 1454 iktv = mikv(ji,jj) 1455 zu_i(ji,jj) = uu(ji,jj,iktu,Nii) - un_b(ji,jj) 1456 zv_i(ji,jj) = vv(ji,jj,iktv,Nii) - vn_b(ji,jj) 1457 END_2D 1542 1458 ELSE ! CENTRED integration: use BEFORE top baroclinic velocity 1543 1459 1544 DO jj = 2, jpjm1 1545 DO ji = 2, jpim1 ! INNER domain 1546 iktu = miku(ji,jj) 1547 iktv = mikv(ji,jj) 1548 zu_i(ji,jj) = ub(ji,jj,iktu) - ub_b(ji,jj) 1549 zv_i(ji,jj) = vb(ji,jj,iktv) - vb_b(ji,jj) 1550 END DO 1551 END DO 1460 DO_2D_00_00 1461 iktu = miku(ji,jj) 1462 iktv = mikv(ji,jj) 1463 zu_i(ji,jj) = uu(ji,jj,iktu,Nnn) - ub_b(ji,jj) 1464 zv_i(ji,jj) = vv(ji,jj,iktv,Nnn) - vb_b(ji,jj) 1465 END_2D 1552 1466 ENDIF 1553 1467 ! 1554 1468 ! ! use "unclipped" top drag (even if explicit friction is used in 3D calculation) 1555 1469 1556 DO jj = 2, jpjm1 1557 DO ji = 2, jpim1 ! INNER domain 1558 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) 1559 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) 1560 END DO 1561 END DO 1470 DO_2D_00_00 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 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) 1473 END_2D 1562 1474 ! 1563 1475 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
