- Timestamp:
- 2021-04-27T17:33:44+02:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/ZDF/zdfgls.F90
r14601 r14757 179 179 180 180 ! Compute surface, top and bottom friction at T-points 181 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) !== surface ocean friction182 181 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) !== surface ocean friction 183 182 ustar2_surf(ji,jj) = r1_rho0 * taum(ji,jj) * tmask(ji,jj,1) ! surface friction … … 187 186 ! 188 187 IF( .NOT.ln_drg_OFF ) THEN !== top/bottom friction (explicit before friction) 189 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) ! bottom friction (explicit before friction)190 188 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! bottom friction (explicit before friction) 191 189 zmsku = 0.5_wp * ( 2._wp - umask(ji-1,jj,mbkt(ji,jj)) * umask(ji,jj,mbkt(ji,jj)) ) … … 195 193 END_2D 196 194 IF( ln_isfcav ) THEN 197 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) ! top friction198 195 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! top friction 199 196 zmsku = 0.5_wp * ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) … … 223 220 zhsro(:,:) = ( (1._wp-zice_fra(:,:)) * zhsro(:,:) + zice_fra(:,:) * rn_hsri )*tmask(:,:,1) + (1._wp - tmask(:,:,1))*rn_hsro 224 221 ! 225 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) !== Compute dissipation rate ==!226 222 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) !== Compute dissipation rate ==! 227 223 eps(ji,jj,jk) = rc03 * en(ji,jj,jk) * SQRT( en(ji,jj,jk) ) / hmxl_n(ji,jj,jk) … … 233 229 234 230 IF( nn_clos == 0 ) THEN ! Mellor-Yamada 235 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )236 231 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 237 232 zup = hmxl_n(ji,jj,jk) * gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) … … 255 250 ! Warning : after this step, en : right hand side of the matrix 256 251 257 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )258 252 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 259 253 ! … … 333 327 ! at k=2, set de/dz=Fw 334 328 !cbr 335 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) ! zdiag zd_lw not defined/used on the halo336 329 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! zdiag zd_lw not defined/used on the halo 337 330 zdiag(ji,jj,2) = zdiag(ji,jj,2) + zd_lw(ji,jj,2) ! Remove zd_lw from zdiag … … 355 348 ! ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = rn_lmin 356 349 ! ! Balance between the production and the dissipation terms 357 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )358 350 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 359 351 !!gm This means that bottom and ocean w-level above have a specified "en" value. Sure ???? … … 374 366 ! 375 367 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 376 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )377 368 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 378 369 itop = mikt(ji,jj) ! k top w-point … … 393 384 CASE ( 1 ) ! Neumman boundary condition 394 385 ! 395 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )396 386 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 397 387 ibot = mbkt(ji,jj) + 1 ! k bottom level of w-point … … 409 399 END_2D 410 400 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 411 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )412 401 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 413 402 itop = mikt(ji,jj) ! k top w-point … … 431 420 ! ---------------------------------------------------------- 432 421 ! 433 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1434 422 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 435 423 zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 436 424 END_3D 437 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1438 425 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1 439 426 zd_lw(ji,jj,jk) = en(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 440 427 END_3D 441 ! [comm_cleanup] ! DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! Third recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk442 428 DO_3DS( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, jpkm1, 2, -1 ) ! Third recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk 443 429 en(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * en(ji,jj,jk+1) ) / zdiag(ji,jj,jk) … … 455 441 ! 456 442 CASE( 0 ) ! k-kl (Mellor-Yamada) 457 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )458 443 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 459 444 psi(ji,jj,jk) = eb(ji,jj,jk) * hmxl_b(ji,jj,jk) … … 461 446 ! 462 447 CASE( 1 ) ! k-eps 463 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )464 448 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 465 449 psi(ji,jj,jk) = eps(ji,jj,jk) … … 467 451 ! 468 452 CASE( 2 ) ! k-w 469 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )470 453 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 471 454 psi(ji,jj,jk) = SQRT( eb(ji,jj,jk) ) / ( rc0 * hmxl_b(ji,jj,jk) ) … … 473 456 ! 474 457 CASE( 3 ) ! generic 475 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )476 458 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 477 459 psi(ji,jj,jk) = rc02 * eb(ji,jj,jk) * hmxl_b(ji,jj,jk)**rnn … … 487 469 ! Warning : after this step, en : right hand side of the matrix 488 470 489 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )490 471 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 491 472 ! … … 560 541 ! 561 542 ! Neumann condition at k=2 562 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) ! zdiag zd_lw not defined/used on the halo563 543 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! zdiag zd_lw not defined/used on the halo 564 544 zdiag(ji,jj,2) = zdiag(ji,jj,2) + zd_lw(ji,jj,2) ! Remove zd_lw from zdiag … … 589 569 ! ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = vkarmn * r_z0_bot 590 570 ! ! Balance between the production and the dissipation terms 591 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )592 571 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 593 572 ibot = mbkt(ji,jj) + 1 ! k bottom level of w-point … … 609 588 CASE ( 1 ) ! Neumman boundary condition 610 589 ! 611 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )612 590 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 613 591 ibot = mbkt(ji,jj) + 1 ! k bottom level of w-point … … 638 616 ! ---------------- 639 617 ! 640 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1641 618 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 642 619 zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 643 620 END_3D 644 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1645 621 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1 646 622 zd_lw(ji,jj,jk) = psi(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 647 623 END_3D 648 ! [comm_cleanup] ! DO_3DS( 0, 0, 0, 0, jpkm1, 2, -1 ) ! Third recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk649 624 DO_3DS( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, jpkm1, 2, -1 ) ! Third recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk 650 625 psi(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * psi(ji,jj,jk+1) ) / zdiag(ji,jj,jk) … … 657 632 ! 658 633 CASE( 0 ) ! k-kl (Mellor-Yamada) 659 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )660 634 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 661 635 eps(ji,jj,jk) = rc03 * en(ji,jj,jk) * en(ji,jj,jk) * SQRT( en(ji,jj,jk) ) / MAX( psi(ji,jj,jk), rn_epsmin) … … 663 637 ! 664 638 CASE( 1 ) ! k-eps 665 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )666 639 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 667 640 eps(ji,jj,jk) = psi(ji,jj,jk) … … 669 642 ! 670 643 CASE( 2 ) ! k-w 671 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )672 644 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 673 645 eps(ji,jj,jk) = rc04 * en(ji,jj,jk) * psi(ji,jj,jk) … … 678 650 zex1 = ( 1.5_wp + rmm/rnn ) 679 651 zex2 = -1._wp / rnn 680 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )681 652 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 682 653 eps(ji,jj,jk) = zcoef * en(ji,jj,jk)**zex1 * psi(ji,jj,jk)**zex2 … … 687 658 ! Limit dissipation rate under stable stratification 688 659 ! -------------------------------------------------- 689 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! Note that this set boundary conditions on hmxl_n at the same time690 660 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) ! Note that this set boundary conditions on hmxl_n at the same time 691 661 ! limitation … … 704 674 ! 705 675 CASE ( 0 , 1 ) ! Galperin or Kantha-Clayson stability functions 706 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )707 676 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 708 677 ! zcof = l²/q² … … 722 691 ! 723 692 CASE ( 2, 3 ) ! Canuto stability functions 724 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 )725 693 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 726 694 ! zcof = l²/q² … … 755 723 ! default value, in case jpk > mbkt(ji,jj)+1. Not needed but avoid a bug when looking for undefined values (-fpe0) 756 724 zstm(:,:,jpk) = 0. 757 ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) ! update bottom with good values758 725 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) ! update bottom with good values 759 726 zstm(ji,jj,mbkt(ji,jj)+1) = zstm(ji,jj,mbkt(ji,jj)) … … 771 738 ! later overwritten by surface/bottom boundaries conditions, so we don't really care of p_avm(:,:1) and p_avm(:,:jpk) 772 739 ! for zd_lw and zd_up but they have to be defined to avoid a bug when looking for undefined values (-fpe0) 773 ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpk )774 740 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk ) 775 741 zsqen = SQRT( 2._wp * en(ji,jj,jk) ) * hmxl_n(ji,jj,jk)
Note: See TracChangeset
for help on using the changeset viewer.
