Changeset 7753 for trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
- Timestamp:
- 2017-03-03T12:46:59+01:00 (7 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90 (modified) (39 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
r7698 r7753 171 171 !!---------------------------------------------------------------------- 172 172 INTEGER, INTENT(in) :: kt ! ocean time step 173 INTEGER :: jk, jj, ji174 173 !!---------------------------------------------------------------------- 175 174 ! … … 180 179 ! 181 180 IF( kt /= nit000 ) THEN ! restore before value to compute tke 182 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 183 DO jk = 1, jpk 184 DO jj = 1, jpj 185 DO ji = 1, jpi 186 avt (ji,jj,jk) = avt_k (ji,jj,jk) 187 avm (ji,jj,jk) = avm_k (ji,jj,jk) 188 avmu(ji,jj,jk) = avmu_k(ji,jj,jk) 189 avmv(ji,jj,jk) = avmv_k(ji,jj,jk) 190 END DO 191 END DO 192 END DO 181 avt (:,:,:) = avt_k (:,:,:) 182 avm (:,:,:) = avm_k (:,:,:) 183 avmu(:,:,:) = avmu_k(:,:,:) 184 avmv(:,:,:) = avmv_k(:,:,:) 193 185 ENDIF 194 186 ! … … 197 189 CALL tke_avn ! now avt, avm, avmu, avmv 198 190 ! 199 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 200 DO jk = 1, jpk 201 DO jj = 1, jpj 202 DO ji = 1, jpi 203 avt_k (ji,jj,jk) = avt (ji,jj,jk) 204 avm_k (ji,jj,jk) = avm (ji,jj,jk) 205 avmu_k(ji,jj,jk) = avmu(ji,jj,jk) 206 avmv_k(ji,jj,jk) = avmv(ji,jj,jk) 207 END DO 208 END DO 209 END DO 191 avt_k (:,:,:) = avt (:,:,:) 192 avm_k (:,:,:) = avm (:,:,:) 193 avmu_k(:,:,:) = avmu(:,:,:) 194 avmv_k(:,:,:) = avmv(:,:,:) 210 195 ! 211 196 #if defined key_agrif … … 268 253 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 269 254 IF ( ln_isfcav ) THEN 270 !$OMP PARALLEL DO schedule(static) private(jj, ji)271 255 DO jj = 2, jpjm1 ! en(mikt(ji,jj)) = rn_emin 272 256 DO ji = fs_2, fs_jpim1 ! vector opt. … … 275 259 END DO 276 260 END IF 277 !$OMP PARALLEL DO schedule(static) private(jj, ji)278 261 DO jj = 2, jpjm1 ! en(1) = rn_ebb taum / rau0 (min value rn_emin0) 279 262 DO ji = fs_2, fs_jpim1 ! vector opt. … … 310 293 ! 311 294 ! !* total energy produce by LC : cumulative sum over jk 312 !$OMP PARALLEL 313 !$OMP DO schedule(static) private(jj, ji) 314 DO jj =1, jpj 315 DO ji=1, jpi 316 zpelc(ji,jj,1) = MAX( rn2b(ji,jj,1), 0._wp ) * gdepw_n(ji,jj,1) * e3w_n(ji,jj,1) 317 END DO 318 END DO 295 zpelc(:,:,1) = MAX( rn2b(:,:,1), 0._wp ) * gdepw_n(:,:,1) * e3w_n(:,:,1) 319 296 DO jk = 2, jpk 320 !$OMP DO schedule(static) private(jj, ji) 321 DO jj =1, jpj 322 DO ji=1, jpi 323 zpelc(ji,jj,jk) = zpelc(ji,jj,jk-1) + MAX( rn2b(ji,jj,jk), 0._wp ) * gdepw_n(ji,jj,jk) * e3w_n(ji,jj,jk) 324 END DO 325 END DO 297 zpelc(:,:,jk) = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0._wp ) * gdepw_n(:,:,jk) * e3w_n(:,:,jk) 326 298 END DO 327 299 ! !* finite Langmuir Circulation depth 328 300 zcof = 0.5 * 0.016 * 0.016 / ( zrhoa * zcdrag ) 329 !$OMP DO schedule(static) private(jj,ji) 330 DO jj = 1, jpj 331 DO ji = 1, jpi 332 imlc(ji,jj) = mbkt(ji,jj) + 1 ! Initialization to the number of w ocean point (=2 over land) 333 END DO 334 END DO 301 imlc(:,:) = mbkt(:,:) + 1 ! Initialization to the number of w ocean point (=2 over land) 335 302 DO jk = jpkm1, 2, -1 336 !$OMP DO schedule(static) private(jj, ji, zus)337 303 DO jj = 1, jpj ! Last w-level at which zpelc>=0.5*us*us 338 304 DO ji = 1, jpi ! with us=0.016*wind(starting from jpk-1) … … 343 309 END DO 344 310 ! ! finite LC depth 345 !$OMP DO schedule(static) private(jj, ji)346 311 DO jj = 1, jpj 347 312 DO ji = 1, jpi … … 350 315 END DO 351 316 zcof = 0.016 / SQRT( zrhoa * zcdrag ) 352 !$OMP DO schedule(static) private(jk, jj, ji, zus, zind, zwlc)353 317 DO jk = 2, jpkm1 !* TKE Langmuir circulation source term added to en 354 318 DO jj = 2, jpjm1 … … 364 328 END DO 365 329 END DO 366 !$OMP END PARALLEL367 330 ! 368 331 ENDIF … … 375 338 ! ! zdiag : diagonal zd_up : upper diagonal zd_lw : lower diagonal 376 339 ! 377 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji)378 340 DO jk = 2, jpkm1 !* Shear production at uw- and vw-points (energy conserving form) 379 341 DO jj = 1, jpjm1 … … 394 356 ! Note that zesh2 is also computed in the next loop. 395 357 ! We decided to compute it twice to keep code readability and avoid an IF case in the DO loops 396 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zesh2, zri)397 358 DO jk = 2, jpkm1 398 359 DO jj = 2, jpjm1 … … 411 372 ENDIF 412 373 ! 413 !$OMP PARALLEL414 !$OMP DO schedule(static) private(jk, jj, ji, zcof, zzd_up, zzd_lw, zesh2)415 374 DO jk = 2, jpkm1 !* Matrix and right hand side in en 416 375 DO jj = 2, jpjm1 … … 446 405 ! !* Matrix inversion from level 2 (tke prescribed at level 1) 447 406 DO jk = 3, jpkm1 ! First recurrence : Dk = Dk - Lk * Uk-1 / Dk-1 448 !$OMP DO schedule(static) private(jj, ji)449 407 DO jj = 2, jpjm1 450 408 DO ji = fs_2, fs_jpim1 ! vector opt. … … 453 411 END DO 454 412 END DO 455 !$OMP DO schedule(static) private(jj, ji)456 413 DO jj = 2, jpjm1 ! Second recurrence : Lk = RHSk - Lk / Dk-1 * Lk-1 457 414 DO ji = fs_2, fs_jpim1 ! vector opt. … … 460 417 END DO 461 418 DO jk = 3, jpkm1 462 !$OMP DO schedule(static) private(jj, ji)463 419 DO jj = 2, jpjm1 464 420 DO ji = fs_2, fs_jpim1 ! vector opt. … … 467 423 END DO 468 424 END DO 469 !$OMP DO schedule(static) private(jj, ji)470 425 DO jj = 2, jpjm1 ! thrid recurrence : Ek = ( Lk - Uk * Ek+1 ) / Dk 471 426 DO ji = fs_2, fs_jpim1 ! vector opt. … … 474 429 END DO 475 430 DO jk = jpk-2, 2, -1 476 !$OMP DO schedule(static) private(jj, ji)477 431 DO jj = 2, jpjm1 478 432 DO ji = fs_2, fs_jpim1 ! vector opt. … … 481 435 END DO 482 436 END DO 483 !$OMP DO schedule(static) private(jk,jj, ji)484 437 DO jk = 2, jpkm1 ! set the minimum value of tke 485 438 DO jj = 2, jpjm1 … … 489 442 END DO 490 443 END DO 491 !$OMP END PARALLEL492 444 493 445 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< … … 498 450 499 451 IF( nn_etau == 1 ) THEN !* penetration below the mixed layer (rn_efr fraction) 500 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji)501 452 DO jk = 2, jpkm1 502 453 DO jj = 2, jpjm1 … … 508 459 END DO 509 460 ELSEIF( nn_etau == 2 ) THEN !* act only at the base of the mixed layer (jk=nmln) (rn_efr fraction) 510 !$OMP PARALLEL DO schedule(static) private(jj, ji, jk)511 461 DO jj = 2, jpjm1 512 462 DO ji = fs_2, fs_jpim1 ! vector opt. … … 517 467 END DO 518 468 ELSEIF( nn_etau == 3 ) THEN !* penetration belox the mixed layer (HF variability) 519 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, ztx2, zty2, ztau, zdif)520 469 DO jk = 2, jpkm1 521 470 DO jj = 2, jpjm1 … … 596 545 ! 597 546 ! initialisation of interior minimum value (avoid a 2d loop with mikt) 598 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 599 DO jk = 1, jpk 600 DO jj = 1, jpj 601 DO ji = 1, jpi 602 zmxlm(ji,jj,jk) = rmxl_min 603 zmxld(ji,jj,jk) = rmxl_min 604 END DO 605 END DO 606 END DO 547 zmxlm(:,:,:) = rmxl_min 548 zmxld(:,:,:) = rmxl_min 607 549 ! 608 550 IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rau0*g) 609 !$OMP PARALLEL DO schedule(static) private(jj, ji, zraug)610 551 DO jj = 2, jpjm1 611 552 DO ji = fs_2, fs_jpim1 … … 615 556 END DO 616 557 ELSE 617 !$OMP PARALLEL DO schedule(static) private(jj,ji) 618 DO jj = 1, jpj 619 DO ji = 1, jpi 620 zmxlm(ji,jj,1) = rn_mxl0 621 END DO 622 END DO 558 zmxlm(:,:,1) = rn_mxl0 623 559 ENDIF 624 560 ! 625 !$OMP PARALLEL626 !$OMP DO schedule(static) private(jk, jj, ji, zrn2)627 561 DO jk = 2, jpkm1 ! interior value : l=sqrt(2*e/n^2) 628 562 DO jj = 2, jpjm1 … … 636 570 ! !* Physical limits for the mixing length 637 571 ! 638 !$OMP DO schedule(static) private(jj,ji) 639 DO jj = 1, jpj 640 DO ji = 1, jpi 641 zmxld(ji,jj, 1 ) = zmxlm(ji,jj,1) ! surface set to the minimum value 642 zmxld(ji,jj,jpk) = rmxl_min ! last level set to the minimum value 643 END DO 644 END DO 645 !$OMP END PARALLEL 572 zmxld(:,:, 1 ) = zmxlm(:,:,1) ! surface set to the minimum value 573 zmxld(:,:,jpk) = rmxl_min ! last level set to the minimum value 646 574 ! 647 575 SELECT CASE ( nn_mxl ) … … 650 578 ! where wmask = 0 set zmxlm == e3w_n 651 579 CASE ( 0 ) ! bounded by the distance to surface and bottom 652 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zemxl)653 580 DO jk = 2, jpkm1 654 581 DO jj = 2, jpjm1 … … 664 591 ! 665 592 CASE ( 1 ) ! bounded by the vertical scale factor 666 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zemxl)667 593 DO jk = 2, jpkm1 668 594 DO jj = 2, jpjm1 … … 676 602 ! 677 603 CASE ( 2 ) ! |dk[xml]| bounded by e3t : 678 !$OMP PARALLEL679 604 DO jk = 2, jpkm1 ! from the surface to the bottom : 680 !$OMP DO schedule(static) private(jj, ji)681 605 DO jj = 2, jpjm1 682 606 DO ji = fs_2, fs_jpim1 ! vector opt. … … 686 610 END DO 687 611 DO jk = jpkm1, 2, -1 ! from the bottom to the surface : 688 !$OMP DO schedule(static) private(jj, ji, zemxl)689 612 DO jj = 2, jpjm1 690 613 DO ji = fs_2, fs_jpim1 ! vector opt. … … 695 618 END DO 696 619 END DO 697 !$OMP END PARALLEL698 620 ! 699 621 CASE ( 3 ) ! lup and ldown, |dk[xml]| bounded by e3t : 700 !$OMP PARALLEL701 622 DO jk = 2, jpkm1 ! from the surface to the bottom : lup 702 !$OMP DO schedule(static) private(jj, ji)703 623 DO jj = 2, jpjm1 704 624 DO ji = fs_2, fs_jpim1 ! vector opt. … … 708 628 END DO 709 629 DO jk = jpkm1, 2, -1 ! from the bottom to the surface : ldown 710 !$OMP DO schedule(static) private(jj, ji)711 630 DO jj = 2, jpjm1 712 631 DO ji = fs_2, fs_jpim1 ! vector opt. … … 715 634 END DO 716 635 END DO 717 !$OMP DO schedule(static) private(jk, jj, ji, zemlm, zemlp)718 636 DO jk = 2, jpkm1 719 637 DO jj = 2, jpjm1 … … 726 644 END DO 727 645 END DO 728 !$OMP END PARALLEL729 646 ! 730 647 END SELECT 731 648 ! 732 649 # if defined key_c1d 733 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 734 DO jk = 1, jpk 735 DO jj = 1, jpj 736 DO ji = 1, jpi 737 e_dis(ji,jj,jk) = zmxld(ji,jj,jk) ! c1d configuration : save mixing and dissipation turbulent length scales 738 e_mix(ji,jj,jk) = zmxlm(ji,jj,jk) 739 END DO 740 END DO 741 END DO 650 e_dis(:,:,:) = zmxld(:,:,:) ! c1d configuration : save mixing and dissipation turbulent length scales 651 e_mix(:,:,:) = zmxlm(:,:,:) 742 652 # endif 743 653 … … 745 655 ! ! Vertical eddy viscosity and diffusivity (avmu, avmv, avt) 746 656 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 747 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji, zsqen, zav)748 657 DO jk = 1, jpkm1 !* vertical eddy viscosity & diffivity at w-points 749 658 DO jj = 2, jpjm1 … … 759 668 CALL lbc_lnk( avm, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) 760 669 ! 761 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji)762 670 DO jk = 2, jpkm1 !* vertical eddy viscosity at wu- and wv-points 763 671 DO jj = 2, jpjm1 … … 771 679 ! 772 680 IF( nn_pdl == 1 ) THEN !* Prandtl number case: update avt 773 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji)774 681 DO jk = 2, jpkm1 775 682 DO jj = 2, jpjm1 … … 891 798 SELECT CASE( nn_htau ) ! Choice of the depth of penetration 892 799 CASE( 0 ) ! constant depth penetration (here 10 meters) 893 !$OMP PARALLEL DO schedule(static) private(jj,ji) 894 DO jj = 1, jpj 895 DO ji = 1, jpi 896 htau(ji,jj) = 10._wp 897 END DO 898 END DO 800 htau(:,:) = 10._wp 899 801 CASE( 1 ) ! F(latitude) : 0.5m to 30m poleward of 40 degrees 900 !$OMP PARALLEL DO schedule(static) private(jj,ji) 901 DO jj = 1, jpj 902 DO ji = 1, jpi 903 htau(ji,jj) = MAX( 0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(ji,jj) ) ) ) ) 904 END DO 905 END DO 802 htau(:,:) = MAX( 0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) ) ) 906 803 END SELECT 907 804 ENDIF 908 805 ! !* set vertical eddy coef. to the background value 909 !$OMP PARALLEL910 !$OMP DO schedule(static) private(jk,jj,ji)911 806 DO jk = 1, jpk 912 DO jj = 1, jpj 913 DO ji = 1, jpi 914 avt (ji,jj,jk) = avtb(jk) * wmask (ji,jj,jk) 915 avm (ji,jj,jk) = avmb(jk) * wmask (ji,jj,jk) 916 avmu(ji,jj,jk) = avmb(jk) * wumask(ji,jj,jk) 917 avmv(ji,jj,jk) = avmb(jk) * wvmask(ji,jj,jk) 918 END DO 919 END DO 920 END DO 921 !$OMP END DO NOWAIT 922 !$OMP DO schedule(static) private(jk,jj,ji) 923 DO jk = 1, jpk 924 DO jj = 1, jpj 925 DO ji = 1, jpi 926 dissl(ji,jj,jk) = 1.e-12_wp 927 END DO 928 END DO 929 END DO 930 !$OMP END PARALLEL 807 avt (:,:,jk) = avtb(jk) * wmask (:,:,jk) 808 avm (:,:,jk) = avmb(jk) * wmask (:,:,jk) 809 avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk) 810 avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk) 811 END DO 812 dissl(:,:,:) = 1.e-12_wp 931 813 ! 932 814 CALL tke_rst( nit000, 'READ' ) !* read or initialize all required files … … 948 830 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 949 831 ! 950 INTEGER :: jit, jk , jj, ji! dummy loop indices832 INTEGER :: jit, jk ! dummy loop indices 951 833 INTEGER :: id1, id2, id3, id4, id5, id6 ! local integers 952 834 !!---------------------------------------------------------------------- … … 975 857 ELSE ! No TKE array found: initialisation 976 858 IF(lwp) WRITE(numout,*) ' ===>>>> : previous run without tke scheme, en computed by iterative loop' 977 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 978 DO jk = 1, jpk 979 DO jj = 1, jpj 980 DO ji = 1, jpi 981 en (ji,jj,jk) = rn_emin * tmask(ji,jj,jk) 982 END DO 983 END DO 984 END DO 859 en (:,:,:) = rn_emin * tmask(:,:,:) 985 860 CALL tke_avn ! recompute avt, avm, avmu, avmv and dissl (approximation) 986 861 ! 987 !$OMP PARALLEL DO schedule(static) private(jk,jj,ji) 988 DO jk = 1, jpk 989 DO jj = 1, jpj 990 DO ji = 1, jpi 991 avt_k (ji,jj,jk) = avt (ji,jj,jk) 992 avm_k (ji,jj,jk) = avm (ji,jj,jk) 993 avmu_k(ji,jj,jk) = avmu(ji,jj,jk) 994 avmv_k(ji,jj,jk) = avmv(ji,jj,jk) 995 END DO 996 END DO 997 END DO 862 avt_k (:,:,:) = avt (:,:,:) 863 avm_k (:,:,:) = avm (:,:,:) 864 avmu_k(:,:,:) = avmu(:,:,:) 865 avmv_k(:,:,:) = avmv(:,:,:) 998 866 ! 999 867 DO jit = nit000 + 1, nit000 + 10 ; CALL zdf_tke( jit ) ; END DO 1000 868 ENDIF 1001 869 ELSE !* Start from rest 1002 !$OMP PARALLEL 1003 !$OMP DO schedule(static) private(jk,jj,ji) 1004 DO jk = 1, jpk 1005 DO jj = 1, jpj 1006 DO ji = 1, jpi 1007 en(ji,jj,jk) = rn_emin * tmask(ji,jj,jk) 1008 END DO 1009 END DO 870 en(:,:,:) = rn_emin * tmask(:,:,:) 871 DO jk = 1, jpk ! set the Kz to the background value 872 avt (:,:,jk) = avtb(jk) * wmask (:,:,jk) 873 avm (:,:,jk) = avmb(jk) * wmask (:,:,jk) 874 avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk) 875 avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk) 1010 876 END DO 1011 !$OMP END DO NOWAIT1012 !$OMP DO schedule(static) private(jk)1013 DO jk = 1, jpk ! set the Kz to the background value1014 DO jj = 1, jpj1015 DO ji = 1, jpi1016 avt (ji,jj,jk) = avtb(jk) * wmask (ji,jj,jk)1017 avm (ji,jj,jk) = avmb(jk) * wmask (ji,jj,jk)1018 avmu(ji,jj,jk) = avmb(jk) * wumask(ji,jj,jk)1019 avmv(ji,jj,jk) = avmb(jk) * wvmask(ji,jj,jk)1020 END DO1021 END DO1022 END DO1023 !$OMP END PARALLEL1024 877 ENDIF 1025 878 !
Note: See TracChangeset
for help on using the changeset viewer.
