Changeset 7753 for trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90
- Timestamp:
- 2017-03-03T12:46:59+01:00 (7 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90 (modified) (36 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftmx.F90
r7698 r7753 121 121 ! ! ----------------------- ! 122 122 ! !* First estimation (with n2 bound by rn_n2min) bounded by 60 cm2/s 123 !$OMP PARALLEL 124 !$OMP DO schedule(static) private(jk,jj,ji) 125 DO jk = 1, jpk 126 DO jj = 1, jpj 127 DO ji = 1, jpi 128 zav_tide(ji,jj,jk) = MIN( 60.e-4, az_tmx(ji,jj,jk) / MAX( rn_n2min, rn2(ji,jj,jk) ) ) 129 END DO 130 END DO 131 END DO 132 !$OMP END DO NOWAIT 133 134 !$OMP DO schedule(static) private(jj, ji) 135 DO jj = 1, jpj 136 DO ji = 1, jpi 137 zkz(ji,jj) = 0.e0 !* Associated potential energy consummed over the whole water column 138 END DO 139 END DO 123 zav_tide(:,:,:) = MIN( 60.e-4, az_tmx(:,:,:) / MAX( rn_n2min, rn2(:,:,:) ) ) 124 125 zkz(:,:) = 0.e0 !* Associated potential energy consummed over the whole water column 140 126 DO jk = 2, jpkm1 141 !$OMP DO schedule(static) private(jj, ji) 142 DO jj = 1, jpj 143 DO ji = 1, jpi 144 zkz(ji,jj) = zkz(ji,jj) + e3w_n(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) ) * rau0 * zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 145 END DO 146 END DO 147 END DO 148 149 !$OMP DO schedule(static) private(jj, ji) 127 zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zav_tide(:,:,jk) * wmask(:,:,jk) 128 END DO 129 150 130 DO jj = 1, jpj !* Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 151 131 DO ji = 1, jpi … … 155 135 156 136 DO jk = 2, jpkm1 !* Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zav_tide bound by 300 cm2/s 157 !$OMP DO schedule(static) private(jj, ji) 158 DO jj = 1, jpj 159 DO ji = 1, jpi 160 zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * MIN( zkz(ji,jj), 30./6. ) * wmask(ji,jj,jk) !kz max = 300 cm2/s 161 END DO 162 END DO 163 END DO 164 !$OMP END PARALLEL 137 zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk) !kz max = 300 cm2/s 138 END DO 165 139 166 140 IF( kt == nit000 ) THEN !* check at first time-step: diagnose the energy consumed by zav_tide 167 141 ztpc = 0._wp 168 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji,ztpc)169 142 DO jk= 1, jpk 170 143 DO jj= 1, jpj … … 189 162 ! ! Update mixing coefs ! 190 163 ! ! ----------------------- ! 191 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)192 164 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with tidal mixing 193 DO jj = 1, jpj 194 DO ji = 1, jpi ! vector opt. 195 avt(ji,jj,jk) = avt(ji,jj,jk) + zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 196 avm(ji,jj,jk) = avm(ji,jj,jk) + zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 197 END DO 198 END DO 165 avt(:,:,jk) = avt(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 166 avm(:,:,jk) = avm(:,:,jk) + zav_tide(:,:,jk) * wmask(:,:,jk) 199 167 DO jj = 2, jpjm1 200 168 DO ji = fs_2, fs_jpim1 ! vector opt. … … 257 225 258 226 ! ! compute the form function using N2 at each time step 259 !$OMP PARALLEL 260 !$OMP DO schedule(static) private(jj, ji) 261 DO jj = 1, jpj 262 DO ji = 1, jpi 263 zempba_3d_1(ji,jj,jpk) = 0.e0 264 zempba_3d_2(ji,jj,jpk) = 0.e0 265 END DO 266 END DO 267 !$OMP DO schedule(static) private(jk,jj,ji) 227 zempba_3d_1(:,:,jpk) = 0.e0 228 zempba_3d_2(:,:,jpk) = 0.e0 268 229 DO jk = 1, jpkm1 269 DO jj = 1, jpj 270 DO ji = 1, jpi 271 zdn2dz (ji,jj,jk) = rn2(ji,jj,jk) - rn2(ji,jj,jk+1) ! Vertical profile of dN2/dz 272 zempba_3d_1(ji,jj,jk) = SQRT( MAX( 0.e0, rn2(ji,jj,jk) ) ) ! - - of N 273 zempba_3d_2(ji,jj,jk) = MAX( 0.e0, rn2(ji,jj,jk) ) ! - - of N^2 274 END DO 275 END DO 276 END DO 277 !$OMP END DO NOWAIT 278 ! 279 !$OMP DO schedule(static) private(jj, ji) 280 DO jj = 1, jpj 281 DO ji = 1, jpi 282 zsum (ji,jj) = 0.e0 283 zsum1(ji,jj) = 0.e0 284 zsum2(ji,jj) = 0.e0 285 END DO 286 END DO 230 zdn2dz (:,:,jk) = rn2(:,:,jk) - rn2(:,:,jk+1) ! Vertical profile of dN2/dz 231 zempba_3d_1(:,:,jk) = SQRT( MAX( 0.e0, rn2(:,:,jk) ) ) ! - - of N 232 zempba_3d_2(:,:,jk) = MAX( 0.e0, rn2(:,:,jk) ) ! - - of N^2 233 END DO 234 ! 235 zsum (:,:) = 0.e0 236 zsum1(:,:) = 0.e0 237 zsum2(:,:) = 0.e0 287 238 DO jk= 2, jpk 288 !$OMP DO schedule(static) private(jj,ji) 289 DO jj= 1, jpj 290 DO ji= 1, jpi 291 zsum1(ji,jj) = zsum1(ji,jj) + zempba_3d_1(ji,jj,jk) * e3w_n(ji,jj,jk) * wmask(ji,jj,jk) 292 zsum2(ji,jj) = zsum2(ji,jj) + zempba_3d_2(ji,jj,jk) * e3w_n(ji,jj,jk) * wmask(ji,jj,jk) 293 END DO 294 END DO 295 END DO 296 !$OMP DO schedule(static) private(jj,ji) 239 zsum1(:,:) = zsum1(:,:) + zempba_3d_1(:,:,jk) * e3w_n(:,:,jk) * wmask(:,:,jk) 240 zsum2(:,:) = zsum2(:,:) + zempba_3d_2(:,:,jk) * e3w_n(:,:,jk) * wmask(:,:,jk) 241 END DO 297 242 DO jj = 1, jpj 298 243 DO ji = 1, jpi … … 303 248 304 249 DO jk= 1, jpk 305 !$OMP DO schedule(static) private(jj,ji,zcoef,ztpc)306 250 DO jj = 1, jpj 307 251 DO ji = 1, jpi … … 315 259 END DO 316 260 END DO 317 !$OMP DO schedule(static) private(jj,ji)318 261 DO jj = 1, jpj 319 262 DO ji = 1, jpi … … 324 267 ! ! first estimation bounded by 10 cm2/s (with n2 bounded by rn_n2min) 325 268 zcoef = rn_tfe_itf / ( rn_tfe * rau0 ) 326 !$OMP DO schedule(static) private(jk,jj,ji)327 269 DO jk = 1, jpk 328 DO jj = 1, jpj 329 DO ji = 1, jpi 330 zavt_itf(ji,jj,jk) = MIN( 10.e-4, zcoef * en_tmx(ji,jj) * zsum(ji,jj) * zempba_3d(ji,jj,jk) & 331 & / MAX( rn_n2min, rn2(ji,jj,jk) ) * tmask(ji,jj,jk) ) 332 END DO 333 END DO 334 END DO 335 336 !$OMP DO schedule(static) private(jj, ji) 337 DO jj = 1, jpj 338 DO ji = 1, jpi 339 zkz(ji,jj) = 0.e0 ! Associated potential energy consummed over the whole water column 340 END DO 341 END DO 270 zavt_itf(:,:,jk) = MIN( 10.e-4, zcoef * en_tmx(:,:) * zsum(:,:) * zempba_3d(:,:,jk) & 271 & / MAX( rn_n2min, rn2(:,:,jk) ) * tmask(:,:,jk) ) 272 END DO 273 274 zkz(:,:) = 0.e0 ! Associated potential energy consummed over the whole water column 342 275 DO jk = 2, jpkm1 343 !$OMP DO schedule(static) private(jj,ji) 344 DO jj = 1, jpj 345 DO ji = 1, jpi 346 zkz(ji,jj) = zkz(ji,jj) + e3w_n(ji,jj,jk) * MAX( 0.e0, rn2(ji,jj,jk) ) * rau0 * zavt_itf(ji,jj,jk) * wmask(ji,jj,jk) 347 END DO 348 END DO 349 END DO 350 351 !$OMP DO schedule(static) private(jj,ji) 276 zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX( 0.e0, rn2(:,:,jk) ) * rau0 * zavt_itf(:,:,jk) * wmask(:,:,jk) 277 END DO 278 352 279 DO jj = 1, jpj ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz to recover en_tmx 353 280 DO ji = 1, jpi … … 356 283 END DO 357 284 358 !$OMP DO schedule(static) private(jk,jj,ji)359 285 DO jk = 2, jpkm1 ! Mutiply by zkz to recover en_tmx, BUT bound by 30/6 ==> zavt_itf bound by 300 cm2/s 360 DO jj = 1, jpj 361 DO ji = 1, jpi 362 zavt_itf(ji,jj,jk) = zavt_itf(ji,jj,jk) * MIN( zkz(ji,jj), 120./10. ) * wmask(ji,jj,jk) ! kz max = 120 cm2/s 363 END DO 364 END DO 365 END DO 366 !$OMP END PARALLEL 286 zavt_itf(:,:,jk) = zavt_itf(:,:,jk) * MIN( zkz(:,:), 120./10. ) * wmask(:,:,jk) ! kz max = 120 cm2/s 287 END DO 367 288 368 289 IF( kt == nit000 ) THEN ! diagnose the nergy consumed by zavt_itf 369 290 ztpc = 0.e0 370 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji,ztpc)371 291 DO jk= 1, jpk 372 292 DO jj= 1, jpj … … 383 303 384 304 ! ! Update pav with the ITF mixing coefficient 385 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)386 305 DO jk = 2, jpkm1 387 DO jj= 1, jpj 388 DO ji= 1, jpi 389 pav(ji,jj,jk) = pav (ji,jj,jk) * ( 1.e0 - mask_itf(ji,jj) ) & 390 & + zavt_itf(ji,jj,jk) * mask_itf(ji,jj) 391 END DO 392 END DO 306 pav(:,:,jk) = pav (:,:,jk) * ( 1.e0 - mask_itf(:,:) ) & 307 & + zavt_itf(:,:,jk) * mask_itf(:,:) 393 308 END DO 394 309 ! … … 494 409 ! ! only the energy available for mixing is taken into account, 495 410 ! ! (mixing efficiency tidal dissipation efficiency) 496 !$OMP PARALLEL 497 498 !$OMP DO schedule(static) private(jj, ji) 499 DO jj = 1, jpj 500 DO ji = 1, jpi 501 en_tmx(ji,jj) = - rn_tfe * rn_me * ( zem2(ji,jj) * 1.25 + zek1(ji,jj) ) * ssmask(ji,jj) 502 END DO 503 END DO 411 en_tmx(:,:) = - rn_tfe * rn_me * ( zem2(:,:) * 1.25 + zek1(:,:) ) * ssmask(:,:) 504 412 505 413 !============ … … 508 416 !! the error is thus ~1% which I feel comfortable with, compared to uncertainties in tidal energy dissipation. 509 417 ! ! Vertical structure (az_tmx) 510 !$OMP DO schedule(static) private(jj, ji)511 418 DO jj = 1, jpj ! part independent of the level 512 419 DO ji = 1, jpi … … 516 423 END DO 517 424 END DO 518 !$OMP DO schedule(static) private(jk, jj, ji)519 425 DO jk= 1, jpk ! complete with the level-dependent part 520 426 DO jj = 1, jpj … … 524 430 END DO 525 431 END DO 526 !$OMP END PARALLEL527 432 !=========== 528 433 ! … … 531 436 ! Total power consumption due to vertical mixing 532 437 ! zpc = rau0 * 1/rn_me * rn2 * zav_tide 438 zav_tide(:,:,:) = 0.e0 439 DO jk = 2, jpkm1 440 zav_tide(:,:,jk) = az_tmx(:,:,jk) / MAX( rn_n2min, rn2(:,:,jk) ) 441 END DO 442 ! 533 443 ztpc = 0._wp 534 !$OMP PARALLEL 535 !$OMP DO schedule(static) private(jk, jj, ji) 536 DO jk = 1, jpk 537 DO jj = 1, jpj 538 DO ji = 1, jpi 539 zav_tide(ji,jj,jk) = 0.e0 540 END DO 541 END DO 542 END DO 543 !$OMP DO schedule(static) private(jk,jj,ji) 544 DO jk = 2, jpkm1 545 DO jj = 1, jpj 546 DO ji = 1, jpi 547 zav_tide(ji,jj,jk) = az_tmx(ji,jj,jk) / MAX( rn_n2min, rn2(ji,jj,jk) ) 548 END DO 549 END DO 550 END DO 551 ! 552 !$OMP DO schedule(static) private(jk, jj, ji) 553 DO jk= 1, jpk 554 DO jj = 1, jpj 555 DO ji = 1, jpi 556 zpc(ji,jj,jk) = MAX(rn_n2min,rn2(ji,jj,jk)) * zav_tide(ji,jj,jk) 557 END DO 558 END DO 559 END DO 560 !$OMP DO schedule(static) private(jk, jj, ji, ztpc) 444 zpc(:,:,:) = MAX(rn_n2min,rn2(:,:,:)) * zav_tide(:,:,:) 561 445 DO jk= 2, jpkm1 562 446 DO jj = 1, jpj … … 566 450 END DO 567 451 END DO 568 !$OMP END PARALLEL569 452 IF( lk_mpp ) CALL mpp_sum( ztpc ) 570 453 ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc … … 574 457 ! 575 458 ! control print 2 576 !$OMP PARALLEL 577 !$OMP DO schedule(static) private(jk, jj, ji) 578 DO jk= 1, jpk 579 DO jj = 1, jpj 580 DO ji = 1, jpi 581 zav_tide(ji,jj,jk) = MIN( zav_tide(ji,jj,jk), 60.e-4 ) 582 zkz(ji,jj) = 0._wp 583 END DO 584 END DO 585 END DO 586 459 zav_tide(:,:,:) = MIN( zav_tide(:,:,:), 60.e-4 ) 460 zkz(:,:) = 0._wp 587 461 DO jk = 2, jpkm1 588 !$OMP DO schedule(static) private(jj, ji) 589 DO jj = 1, jpj 590 DO ji = 1, jpi 591 zkz(ji,jj) = zkz(ji,jj) + e3w_n(ji,jj,jk) * MAX(0.e0, rn2(ji,jj,jk)) * rau0 * zav_tide(ji,jj,jk) * wmask(ji,jj,jk) 592 END DO 593 END DO 462 zkz(:,:) = zkz(:,:) + e3w_n(:,:,jk) * MAX(0.e0, rn2(:,:,jk)) * rau0 * zav_tide(:,:,jk) * wmask(:,:,jk) 594 463 END DO 595 464 ! Here zkz should be equal to en_tmx ==> multiply by en_tmx/zkz 596 !$OMP DO schedule(static) private(jj, ji)597 465 DO jj = 1, jpj 598 466 DO ji = 1, jpi … … 603 471 END DO 604 472 ztpc = 1.e50 605 !$OMP DO schedule(static) private(jj, ji, ztpc)606 473 DO jj = 1, jpj 607 474 DO ji = 1, jpi … … 611 478 END DO 612 479 END DO 613 !$OMP END PARALLEL614 480 WRITE(numout,*) ' Min de zkz ', ztpc, ' Max = ', maxval(zkz(:,:) ) 615 !$OMP PARALLEL616 481 ! 617 !$OMP DO schedule(static) private(jk,jj,ji)618 482 DO jk = 2, jpkm1 619 DO jj = 1, jpj 620 DO ji = 1, jpi 621 zav_tide(ji,jj,jk) = zav_tide(ji,jj,jk) * MIN( zkz(ji,jj), 30./6. ) * wmask(ji,jj,jk) !kz max = 300 cm2/s 622 END DO 623 END DO 483 zav_tide(:,:,jk) = zav_tide(:,:,jk) * MIN( zkz(:,:), 30./6. ) * wmask(:,:,jk) !kz max = 300 cm2/s 624 484 END DO 625 485 ztpc = 0._wp 626 !$OMP DO schedule(static) private(jk, jj, ji) 627 DO jk= 1, jpk 628 DO jj = 1, jpj 629 DO ji = 1, jpi 630 zpc(ji,jj,jk) = Max(0.e0,rn2(ji,jj,jk)) * zav_tide(ji,jj,jk) 631 END DO 632 END DO 633 END DO 634 !$OMP DO schedule(static) private(jk, jj, ji, ztpc) 486 zpc(:,:,:) = Max(0.e0,rn2(:,:,:)) * zav_tide(:,:,:) 635 487 DO jk= 1, jpk 636 488 DO jj = 1, jpj … … 640 492 END DO 641 493 END DO 642 !$OMP END PARALLEL643 494 IF( lk_mpp ) CALL mpp_sum( ztpc ) 644 495 ztpc= rau0 * 1/(rn_tfe * rn_me) * ztpc … … 649 500 & / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask (:,:,jk) * tmask_i(:,:) ) ) 650 501 ztpc = 1.e50 651 !$OMP PARALLEL DO schedule(static) private(ztpc, jj, ji)652 502 DO jj = 1, jpj 653 503 DO ji = 1, jpi … … 663 513 WRITE(numout,*) ' Initial profile of tidal vertical mixing' 664 514 DO jk = 1, jpk 665 !$OMP PARALLEL DO schedule(static) private(jj, ji)666 515 DO jj = 1,jpj 667 516 DO ji = 1,jpi … … 674 523 END DO 675 524 DO jk = 1, jpk 676 !$OMP PARALLEL DO schedule(static) private(jj, ji) 677 DO jj = 1,jpj 678 DO ji = 1,jpi 679 zkz(ji,jj) = az_tmx(ji,jj,jk) /rn_n2min 680 END DO 681 END DO 525 zkz(:,:) = az_tmx(:,:,jk) /rn_n2min 682 526 ze_z = SUM( e1e2t(:,:) * zkz (:,:) * tmask_i(:,:) ) & 683 527 & / MAX( 1.e-20, SUM( e1e2t(:,:) * wmask(:,:,jk) * tmask_i(:,:) ) ) … … 845 689 ! !* Critical slope mixing: distribute energy over the time-varying ocean depth, 846 690 ! using an exponential decay from the seafloor. 847 !$OMP PARALLEL848 !$OMP DO schedule(static) private(jj,ji)849 691 DO jj = 1, jpj ! part independent of the level 850 692 DO ji = 1, jpi … … 855 697 END DO 856 698 857 !$OMP DO schedule(static) private(jk,jj,ji)858 699 DO jk = 2, jpkm1 ! complete with the level-dependent part 859 DO jj = 1, jpj 860 DO ji = 1, jpi 861 emix_tmx(ji,jj,jk) = zfact(ji,jj) * ( EXP( ( gde3w_n(ji,jj,jk ) - zhdep(ji,jj) ) / hcri_tmx(:,:) ) & 862 & - EXP( ( gde3w_n(ji,jj,jk-1) - zhdep(ji,jj) ) / hcri_tmx(ji,jj) ) ) * wmask(ji,jj,jk) & 863 & / ( gde3w_n(ji,jj,jk) - gde3w_n(ji,jj,jk-1) ) 864 END DO 865 END DO 866 END DO 867 !$OMP END PARALLEL 700 emix_tmx(:,:,jk) = zfact(:,:) * ( EXP( ( gde3w_n(:,:,jk ) - zhdep(:,:) ) / hcri_tmx(:,:) ) & 701 & - EXP( ( gde3w_n(:,:,jk-1) - zhdep(:,:) ) / hcri_tmx(:,:) ) ) * wmask(:,:,jk) & 702 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 703 END DO 868 704 869 705 ! !* Pycnocline-intensified mixing: distribute energy over the time-varying … … 874 710 CASE ( 1 ) ! Dissipation scales as N (recommended) 875 711 876 !$OMP PARALLEL 877 !$OMP DO schedule(static) private(jj, ji) 878 DO jj = 1, jpj 879 DO ji = 1, jpi 880 zfact(ji,jj) = 0._wp 881 END DO 882 END DO 712 zfact(:,:) = 0._wp 883 713 DO jk = 2, jpkm1 ! part independent of the level 884 !$OMP DO schedule(static) private(jj,ji) 885 DO jj = 1, jpj ! part independent of the level 886 DO ji = 1, jpi 887 zfact(ji,jj) = zfact(ji,jj) + e3w_n(ji,jj,jk) * SQRT( MAX( 0._wp, rn2(ji,jj,jk) ) ) * wmask(ji,jj,jk) 888 END DO 889 END DO 890 END DO 891 892 !$OMP DO schedule(static) private(jj,ji) 714 zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk) 715 END DO 716 893 717 DO jj = 1, jpj 894 718 DO ji = 1, jpi … … 897 721 END DO 898 722 899 !$OMP DO schedule(static) private(jk,jj,ji)900 723 DO jk = 2, jpkm1 ! complete with the level-dependent part 901 DO jj = 1, jpj 902 DO ji = 1, jpi 903 emix_tmx(ji,jj,jk) = emix_tmx(ji,jj,jk) + zfact(ji,jj) * SQRT( MAX( 0._wp, rn2(ji,jj,jk) ) ) * wmask(ji,ji,jk) 904 END DO 905 END DO 906 END DO 907 !$OMP END PARALLEL 724 emix_tmx(:,:,jk) = emix_tmx(:,:,jk) + zfact(:,:) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk) 725 END DO 908 726 909 727 CASE ( 2 ) ! Dissipation scales as N^2 910 728 911 !$OMP PARALLEL 912 !$OMP DO schedule(static) private(jj, ji) 913 DO jj = 1, jpj 914 DO ji = 1, jpi 915 zfact(ji,jj) = 0._wp 916 END DO 917 END DO 918 919 DO jk = 2, jpkm1 920 !$OMP DO schedule(static) private(jj,ji) 921 DO jj = 1, jpj 922 DO ji = 1, jpi 923 zfact(ji,jj) = zfact(ji,jj) + e3w_n(ji,jj,jk) * MAX( 0._wp, rn2(ji,jj,jk) ) * wmask(ji,jj,jk) 924 END DO 925 END DO 926 END DO 927 928 !$OMP DO schedule(static) private(jj,ji) 729 zfact(:,:) = 0._wp 730 DO jk = 2, jpkm1 ! part independent of the level 731 zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 732 END DO 733 929 734 DO jj= 1, jpj 930 735 DO ji = 1, jpi … … 933 738 END DO 934 739 935 !$OMP DO schedule(static) private(jk,jj,ji)936 740 DO jk = 2, jpkm1 ! complete with the level-dependent part 937 DO jj = 1, jpj 938 DO ji = 1, jpi 939 emix_tmx(ji,jj,jk) = emix_tmx(ji,jj,jk) + zfact(ji,jj) * MAX( 0._wp, rn2(ji,jj,jk) ) * wmask(ji,ji,jk) 940 END DO 941 END DO 942 END DO 943 !$OMP END PARALLEL 741 emix_tmx(:,:,jk) = emix_tmx(:,:,jk) + zfact(:,:) * MAX( 0._wp, rn2(:,:,jk) ) * wmask(:,:,jk) 742 END DO 944 743 945 744 END SELECT … … 948 747 ! !* ocean depth as proportional to rn2 * exp(-z_wkb/rn_hbot) 949 748 950 !$OMP PARALLEL 951 !$OMP DO schedule(static) private(jk,jj,ji) 952 DO jk = 1, jpk 953 DO jj = 1, jpj 954 DO ji = 1, jpi 955 zwkb(ji,jj,jk) = 0._wp 956 END DO 957 END DO 958 END DO 959 !$OMP DO schedule(static) private(jj,ji) 960 DO jj = 1, jpj 961 DO ji = 1, jpi 962 zfact(ji,jj) = 0._wp 963 END DO 964 END DO 749 zwkb(:,:,:) = 0._wp 750 zfact(:,:) = 0._wp 965 751 DO jk = 2, jpkm1 966 !$OMP DO schedule(static) private(jj,ji) 967 DO jj = 1, jpj 968 DO ji = 1, jpi 969 zfact(ji,jj) = zfact(ji,jj) + e3w_n(ji,jj,jk) * SQRT( MAX( 0._wp, rn2(ji,jj,jk) ) ) * wmask(ji,jj,jk) 970 zwkb(ji,jj,jk) = zfact(ji,jj) 971 END DO 972 END DO 973 END DO 974 975 !$OMP DO schedule(static) private(jk,jj,ji) 752 zfact(:,:) = zfact(:,:) + e3w_n(:,:,jk) * SQRT( MAX( 0._wp, rn2(:,:,jk) ) ) * wmask(:,:,jk) 753 zwkb(:,:,jk) = zfact(:,:) 754 END DO 755 976 756 DO jk = 2, jpkm1 977 757 DO jj = 1, jpj … … 982 762 END DO 983 763 END DO 984 985 !$OMP DO schedule(static) private(jj, ji) 986 DO jj = 1, jpj 987 DO ji = 1, jpi 988 zwkb(ji,jj,1) = zhdep(ji,jj) * tmask(ji,jj,1) 989 END DO 990 END DO 991 !$OMP END DO NOWAIT 992 !$OMP DO schedule(static) private(jk,jj,ji) 993 DO jk = 1, jpk 994 DO jj = 1, jpj 995 DO ji = 1, jpi 996 zweight(ji,jj,jk) = 0._wp 997 END DO 998 END DO 999 END DO 1000 1001 !$OMP DO schedule(static) private(jk,jj,ji) 764 zwkb(:,:,1) = zhdep(:,:) * tmask(:,:,1) 765 766 zweight(:,:,:) = 0._wp 1002 767 DO jk = 2, jpkm1 1003 DO jj = 1, jpj 1004 DO ji = 1, jpi 1005 zweight(ji,jj,jk) = MAX( 0._wp, rn2(ji,jj,jk) ) * hbot_tmx(ji,jj) * wmask(ji,jj,jk) & 1006 & * ( EXP( -zwkb(ji,jj,jk) / hbot_tmx(ji,jj) ) - EXP( -zwkb(ji,jj,jk-1) / hbot_tmx(ji,jj) ) ) 1007 END DO 1008 END DO 1009 END DO 1010 1011 !$OMP DO schedule(static) private(jj, ji) 1012 DO jj = 1, jpj 1013 DO ji = 1, jpi 1014 zfact(ji,jj) = 0._wp 1015 END DO 1016 END DO 1017 768 zweight(:,:,jk) = MAX( 0._wp, rn2(:,:,jk) ) * hbot_tmx(:,:) * wmask(:,:,jk) & 769 & * ( EXP( -zwkb(:,:,jk) / hbot_tmx(:,:) ) - EXP( -zwkb(:,:,jk-1) / hbot_tmx(:,:) ) ) 770 END DO 771 772 zfact(:,:) = 0._wp 1018 773 DO jk = 2, jpkm1 ! part independent of the level 1019 !$OMP DO schedule(static) private(jj,ji) 1020 DO jj = 1, jpj 1021 DO ji = 1, jpi 1022 zfact(ji,jj) = zfact(ji,jj) + zweight(ji,jj,jk) 1023 END DO 1024 END DO 1025 END DO 1026 1027 !$OMP DO schedule(static) private(jj,ji) 774 zfact(:,:) = zfact(:,:) + zweight(:,:,jk) 775 END DO 776 1028 777 DO jj = 1, jpj 1029 778 DO ji = 1, jpi … … 1032 781 END DO 1033 782 1034 !$OMP DO schedule(static) private(jk,jj,ji)1035 783 DO jk = 2, jpkm1 ! complete with the level-dependent part 1036 DO jj = 1, jpj 1037 DO ji = 1, jpi 1038 emix_tmx(ji,jj,jk) = emix_tmx(ji,jj,jk) + zweight(ji,jj,jk) * zfact(ji,jj) * wmask(ji,ji,jk) 1039 & / ( gde3w_n(ji,jj,jk) - gde3w_n(ji,jj,jk-1) ) 1040 END DO 1041 END DO 1042 END DO 1043 !$OMP END DO NOWAIT 784 emix_tmx(:,:,jk) = emix_tmx(:,:,jk) + zweight(:,:,jk) * zfact(:,:) * wmask(:,:,jk) & 785 & / ( gde3w_n(:,:,jk) - gde3w_n(:,:,jk-1) ) 786 END DO 1044 787 1045 788 1046 789 ! Calculate molecular kinematic viscosity 1047 !$OMP DO schedule(static) private(jj, ji) 1048 DO jj = 1, jpj 1049 DO ji = 1, jpi 1050 znu_t(ji,jj,jk) = 1.e-4_wp * ( 17.91_wp - 0.53810_wp * tsn(ji,jj,jk,jp_tem) & 1051 & + 0.00694_wp * tsn(ji,jj,jk,jp_tem) * tsn(ji,jj,jk,jp_tem) & 1052 & + 0.02305_wp * tsn(ji,jj,jk,jp_sal) ) * tmask(ji,jj,jk) * r1_rau0 1053 END DO 1054 END DO 1055 !$OMP DO schedule(static) private(jk,jj,ji) 790 znu_t(:,:,:) = 1.e-4_wp * ( 17.91_wp - 0.53810_wp * tsn(:,:,:,jp_tem) + 0.00694_wp * tsn(:,:,:,jp_tem) * tsn(:,:,:,jp_tem) & 791 & + 0.02305_wp * tsn(:,:,:,jp_sal) ) * tmask(:,:,:) * r1_rau0 1056 792 DO jk = 2, jpkm1 1057 DO jj = 1, jpj 1058 DO ji = 1, jpi 1059 znu_w(ji,jj,jk) = 0.5_wp * ( znu_t(ji,jj,jk-1) + znu_t(ji,jj,jk) ) * wmask(ji,jj,jk) 1060 END DO 1061 END DO 793 znu_w(:,:,jk) = 0.5_wp * ( znu_t(:,:,jk-1) + znu_t(:,:,jk) ) * wmask(:,:,jk) 1062 794 END DO 1063 795 1064 796 ! Calculate turbulence intensity parameter Reb 1065 !$OMP DO schedule(static) private(jk,jj,ji)1066 797 DO jk = 2, jpkm1 1067 DO jj = 1, jpj 1068 DO ji = 1, jpi 1069 zReb(ji,jj,jk) = emix_tmx(ji,jj,jk) / MAX( 1.e-20_wp, znu_w(ji,jj,jk) * rn2(ji,jj,jk) ) 1070 END DO 1071 END DO 798 zReb(:,:,jk) = emix_tmx(:,:,jk) / MAX( 1.e-20_wp, znu_w(:,:,jk) * rn2(:,:,jk) ) 1072 799 END DO 1073 800 1074 801 ! Define internal wave-induced diffusivity 1075 !$OMP DO schedule(static) private(jk,jj,ji)1076 802 DO jk = 2, jpkm1 1077 DO jj = 1, jpj 1078 DO ji = 1, jpi 1079 zav_wave(ji,jj,jk) = znu_w(ji,jj,jk) * zReb(ji,jj,jk) * r1_6 ! This corresponds to a constant mixing efficiency of 1/6 1080 END DO 1081 END DO 1082 END DO 1083 !$OMP END PARALLEL 803 zav_wave(:,:,jk) = znu_w(:,:,jk) * zReb(:,:,jk) * r1_6 ! This corresponds to a constant mixing efficiency of 1/6 804 END DO 1084 805 1085 806 IF( ln_mevar ) THEN ! Variable mixing efficiency case : modify zav_wave in the 1086 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)1087 807 DO jk = 2, jpkm1 ! energetic (Reb > 480) and buoyancy-controlled (Reb <10.224 ) regimes 1088 808 DO jj = 1, jpj … … 1098 818 ENDIF 1099 819 1100 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)1101 820 DO jk = 2, jpkm1 ! Bound diffusivity by molecular value and 100 cm2/s 1102 DO jj = 1, jpj 1103 DO ji = 1, jpi 1104 zav_wave(ji,jj,jk) = MIN( MAX( 1.4e-7_wp, zav_wave(ji,jj,jk) ), 1.e-2_wp ) * wmask(ji,jj,jk) 1105 END DO 1106 END DO 821 zav_wave(:,:,jk) = MIN( MAX( 1.4e-7_wp, zav_wave(:,:,jk) ), 1.e-2_wp ) * wmask(:,:,jk) 1107 822 END DO 1108 823 1109 824 IF( kt == nit000 ) THEN !* Control print at first time-step: diagnose the energy consumed by zav_wave 1110 825 ztpc = 0._wp 1111 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji,ztpc)1112 826 DO jk = 2, jpkm1 1113 827 DO jj = 1, jpj … … 1135 849 ! 1136 850 IF( ln_tsdiff ) THEN !* Option for differential mixing of salinity and temperature 1137 !$OMP PARALLEL1138 !$OMP DO schedule(static) private(jk,jj,ji)1139 851 DO jk = 2, jpkm1 ! Calculate S/T diffusivity ratio as a function of Reb 1140 852 DO jj = 1, jpj … … 1146 858 END DO 1147 859 END DO 1148 !$OMP DO schedule(static) private(jk,jj,ji)860 CALL iom_put( "av_ratio", zav_ratio ) 1149 861 DO jk = 2, jpkm1 !* update momentum & tracer diffusivity with wave-driven mixing 1150 DO jj = 1, jpj 1151 DO ji = 1, jpi 1152 fsavs(ji,jj,jk) = avt(ji,jj,jk) + zav_wave(ji,jj,jk) * zav_ratio(ji,jj,jk) 1153 avt (ji,jj,jk) = avt(ji,jj,jk) + zav_wave(ji,jj,jk) 1154 avm (ji,jj,jk) = avm(ji,jj,jk) + zav_wave(ji,jj,jk) 1155 END DO 1156 END DO 1157 END DO 1158 !$OMP END PARALLEL 1159 CALL iom_put( "av_ratio", zav_ratio ) 862 fsavs(:,:,jk) = avt(:,:,jk) + zav_wave(:,:,jk) * zav_ratio(:,:,jk) 863 avt (:,:,jk) = avt(:,:,jk) + zav_wave(:,:,jk) 864 avm (:,:,jk) = avm(:,:,jk) + zav_wave(:,:,jk) 865 END DO 1160 866 ! 1161 867 ELSE !* update momentum & tracer diffusivity with wave-driven mixing 1162 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji)1163 868 DO jk = 2, jpkm1 1164 DO jj = 1, jpj 1165 DO ji = 1, jpi 1166 fsavs(ji,jj,jk) = avt(ji,jj,jk) + zav_wave(ji,jj,jk) 1167 avt (ji,jj,jk) = avt(ji,jj,jk) + zav_wave(ji,jj,jk) 1168 avm (ji,jj,jk) = avm(ji,jj,jk) + zav_wave(ji,jj,jk) 1169 END DO 1170 END DO 1171 END DO 1172 ENDIF 1173 1174 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 869 fsavs(:,:,jk) = avt(:,:,jk) + zav_wave(:,:,jk) 870 avt (:,:,jk) = avt(:,:,jk) + zav_wave(:,:,jk) 871 avm (:,:,jk) = avm(:,:,jk) + zav_wave(:,:,jk) 872 END DO 873 ENDIF 874 1175 875 DO jk = 2, jpkm1 !* update momentum diffusivity at wu and wv points 1176 876 DO jj = 2, jpjm1 … … 1188 888 ! vertical integral of rau0 * Kz * N^2 (pcmap_tmx), energy density (emix_tmx) 1189 889 IF( iom_use("bflx_tmx") .OR. iom_use("pcmap_tmx") ) THEN 1190 !$OMP PARALLEL 1191 !$OMP DO schedule(static) private(jk,jj,ji) 1192 DO jk = 1, jpk 1193 DO jj = 1, jpj 1194 DO ji = 1, jpi 1195 bflx_tmx(ji,jj,jk) = MAX( 0._wp, rn2(ji,jj,jk) ) * zav_wave(ji,jj,jk) 1196 END DO 1197 END DO 1198 END DO 1199 !$OMP END DO NOWAIT 1200 !$OMP DO schedule(static) private(jj, ji) 1201 DO jj = 1, jpj 1202 DO ji = 1, jpi 1203 pcmap_tmx(ji,jj) = 0._wp 1204 END DO 1205 END DO 1206 DO jk = 2, jpkm1 1207 !$OMP DO schedule(static) private(jj, ji) 1208 DO jj = 1, jpj 1209 DO ji = 1, jpi 1210 pcmap_tmx(ji,jj) = pcmap_tmx(ji,jj) + e3w_n(ji,jj,jk) * bflx_tmx(ji,jj,jk) * wmask(ji,jj,jk) 1211 END DO 1212 END DO 1213 END DO 1214 !$OMP DO schedule(static) private(jj, ji) 1215 DO jj = 1, jpj 1216 DO ji = 1, jpi 1217 pcmap_tmx(ji,jj) = rau0 * pcmap_tmx(ji,jj) 1218 END DO 1219 END DO 1220 !$OMP END PARALLEL 890 bflx_tmx(:,:,:) = MAX( 0._wp, rn2(:,:,:) ) * zav_wave(:,:,:) 891 pcmap_tmx(:,:) = 0._wp 892 DO jk = 2, jpkm1 893 pcmap_tmx(:,:) = pcmap_tmx(:,:) + e3w_n(:,:,jk) * bflx_tmx(:,:,jk) * wmask(:,:,jk) 894 END DO 895 pcmap_tmx(:,:) = rau0 * pcmap_tmx(:,:) 1221 896 CALL iom_put( "bflx_tmx", bflx_tmx ) 1222 897 CALL iom_put( "pcmap_tmx", pcmap_tmx ) … … 1295 970 avmb(:) = 1.4e-6_wp ! viscous molecular value 1296 971 avtb(:) = 1.e-10_wp ! very small diffusive minimum (background avt is specified in zdf_tmx) 1297 !$OMP PARALLEL DO schedule(static) private(jj, ji) 1298 DO jj = 1, jpj 1299 DO ji = 1, jpi 1300 avtb_2d(ji,jj) = 1.e0_wp ! uniform 1301 END DO 1302 END DO 972 avtb_2d(:,:) = 1.e0_wp ! uniform 1303 973 IF(lwp) THEN ! Control print 1304 974 WRITE(numout,*) … … 1333 1003 CALL iom_close(inum) 1334 1004 1335 !$OMP PARALLEL DO schedule(static) private(jj, ji) 1336 DO jj = 1, jpj 1337 DO ji = 1, jpi 1338 ebot_tmx(ji,jj) = ebot_tmx(ji,jj) * ssmask(ji,jj) 1339 epyc_tmx(ji,jj) = epyc_tmx(ji,jj) * ssmask(ji,jj) 1340 ecri_tmx(ji,jj) = ecri_tmx(ji,jj) * ssmask(ji,jj) 1341 1342 ! Set once for all to zero the first and last vertical levels of appropriate variables 1343 emix_tmx (ji,jj, 1 ) = 0._wp 1344 emix_tmx (ji,jj,jpk) = 0._wp 1345 zav_ratio(ji,jj, 1 ) = 0._wp 1346 zav_ratio(ji,jj,jpk) = 0._wp 1347 zav_wave (ji,jj, 1 ) = 0._wp 1348 zav_wave (ji,jj,jpk) = 0._wp 1349 END DO 1350 END DO 1005 ebot_tmx(:,:) = ebot_tmx(:,:) * ssmask(:,:) 1006 epyc_tmx(:,:) = epyc_tmx(:,:) * ssmask(:,:) 1007 ecri_tmx(:,:) = ecri_tmx(:,:) * ssmask(:,:) 1008 1009 ! Set once for all to zero the first and last vertical levels of appropriate variables 1010 emix_tmx (:,:, 1 ) = 0._wp 1011 emix_tmx (:,:,jpk) = 0._wp 1012 zav_ratio(:,:, 1 ) = 0._wp 1013 zav_ratio(:,:,jpk) = 0._wp 1014 zav_wave (:,:, 1 ) = 0._wp 1015 zav_wave (:,:,jpk) = 0._wp 1351 1016 1352 1017 zbot = glob_sum( e1e2t(:,:) * ebot_tmx(:,:) )
Note: See TracChangeset
for help on using the changeset viewer.
