- Timestamp:
- 2017-04-29T17:24:54+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
r7953 r7990 27 27 !! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase 28 28 !! 3.6 ! 2014-11 (P. Mathiot) add ice shelf capability 29 !! 4.0 ! 2017-04 (G. Madec) Remove CPP keys29 !! 4.0 ! 2017-04 (G. Madec) remove CPP ddm key & avm at t-point only 30 30 !!---------------------------------------------------------------------- 31 31 … … 159 159 !! 160 160 !! ** Action : compute en (now turbulent kinetic energy) 161 !! update avt, avm u, avmv(before vertical eddy coef.)161 !! update avt, avm (before vertical eddy coef.) 162 162 !! 163 163 !! References : Gaspar et al., JGR, 1990, … … 175 175 #endif 176 176 ! 177 IF( kt /= nit000 ) THEN ! restore before value to compute tke 178 avt (:,:,:) = avt_k (:,:,:) 179 avm (:,:,:) = avm_k (:,:,:) 180 avmu(:,:,:) = avmu_k(:,:,:) 181 avmv(:,:,:) = avmv_k(:,:,:) 182 ENDIF 183 ! 184 CALL tke_tke ! now tke (en) 185 ! 186 CALL tke_avn ! now avt, avm, avmu, avmv 187 ! 188 avt_k (:,:,:) = avt (:,:,:) 189 avm_k (:,:,:) = avm (:,:,:) 190 avmu_k(:,:,:) = avmu(:,:,:) 191 avmv_k(:,:,:) = avmv(:,:,:) 177 avt(:,:,:) = avt_k(:,:,:) ! restore before Kz computed with TKE alone 178 avm(:,:,:) = avm_k(:,:,:) 179 ! 180 CALL tke_tke ! now tke (en) 181 ! 182 CALL tke_avn ! now avt, avm, dissl 183 ! 184 avt_k(:,:,:) = avt(:,:,:) ! store Kz computed with TKE alone 185 avm_k(:,:,:) = avm(:,:,:) 192 186 ! 193 187 #if defined key_agrif … … 195 189 IF( .NOT.Agrif_Root() ) CALL Agrif_Update_Tke( kt ) ! children only 196 190 #endif 197 ! 191 ! 192 IF( lrst_oce ) CALL tke_rst( kt, 'WRITE' ) 193 ! 198 194 END SUBROUTINE zdf_tke 199 195 … … 207 203 !! ** Method : - TKE surface boundary condition 208 204 !! - source term due to Langmuir cells (Axell JGR 2002) (ln_lc=T) 209 !! - source term due to shear ( saved in avmu, avmv arrays)205 !! - source term due to shear (= Kz dz[Ub] * dz[Un] ) 210 206 !! - Now TKE : resolution of the TKE equation by inverting 211 207 !! a tridiagonal linear system by a "methode de chasse" … … 213 209 !! 214 210 !! ** Action : - en : now turbulent kinetic energy) 215 !! - avmu, avmv : production of TKE by shear at u and v-points216 !! (= Kz dz[Ub] * dz[Un] )217 211 !! --------------------------------------------------------------------- 218 212 INTEGER :: ji, jj, jk ! dummy loop arguments … … 228 222 REAL(wp) :: zzd_up, zzd_lw ! - - 229 223 !!bfr REAL(wp) :: zebot ! - - 230 INTEGER , POINTER, DIMENSION(:,:) :: imlc224 INTEGER , DIMENSION(jpi,jpj) :: imlc 231 225 REAL(wp), POINTER, DIMENSION(:,: ) :: zhlc 232 226 REAL(wp), POINTER, DIMENSION(:,:,:) :: zpelc, zdiag, zd_up, zd_lw, z3du, z3dv … … 236 230 IF( nn_timing == 1 ) CALL timing_start('tke_tke') 237 231 ! 238 CALL wrk_alloc( jpi,jpj, imlc ) ! integer239 232 CALL wrk_alloc( jpi,jpj, zhlc ) 240 233 CALL wrk_alloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw, z3du, z3dv ) … … 473 466 END DO 474 467 ENDIF 468 !!gm not sure we need this lbc .... 475 469 CALL lbc_lnk( en, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) 476 470 ! 477 CALL wrk_dealloc( jpi,jpj, imlc ) ! integer478 471 CALL wrk_dealloc( jpi,jpj, zhlc ) 479 472 CALL wrk_dealloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw, z3du, z3dv ) … … 516 509 !! with pdlr=1 if nn_pdl=0, pdlr=1/pdl=F(Ri) otherwise. 517 510 !! 518 !! ** Action : - avt : now vertical eddy diffusivity (w-point) 519 !! - avmu, avmv : now vertical eddy viscosity at uw- and vw-points 511 !! ** Action : - avt, avm : now vertical eddy diffusivity and viscosity (w-point) 520 512 !!---------------------------------------------------------------------- 521 513 INTEGER :: ji, jj, jk ! dummy loop indices 522 REAL(wp) :: zrn2, zraug, zcoef, zav 523 REAL(wp) :: zdku, zri, zsqen! - -524 REAL(wp) :: z dkv, zemxl, zemlm, zemlp! - -514 REAL(wp) :: zrn2, zraug, zcoef, zav ! local scalars 515 REAL(wp) :: zdku, zdkv, zsqen ! - - 516 REAL(wp) :: zemxl, zemlm, zemlp ! - - 525 517 REAL(wp), POINTER, DIMENSION(:,:,:) :: zmpdl, zmxlm, zmxld 526 518 !!-------------------------------------------------------------------- … … 645 637 646 638 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 647 ! ! Vertical eddy viscosity and diffusivity (avm u, avmv,avt)639 ! ! Vertical eddy viscosity and diffusivity (avm and avt) 648 640 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 649 641 DO jk = 1, jpkm1 !* vertical eddy viscosity & diffivity at w-points … … 658 650 END DO 659 651 END DO 660 CALL lbc_lnk( avm, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) 661 ! 662 DO jk = 2, jpkm1 !* vertical eddy viscosity at wu- and wv-points 663 DO jj = 2, jpjm1 664 DO ji = fs_2, fs_jpim1 ! vector opt. 665 avmu(ji,jj,jk) = 0.5 * ( avm(ji,jj,jk) + avm(ji+1,jj ,jk) ) * wumask(ji,jj,jk) 666 avmv(ji,jj,jk) = 0.5 * ( avm(ji,jj,jk) + avm(ji ,jj+1,jk) ) * wvmask(ji,jj,jk) 667 END DO 668 END DO 669 END DO 670 CALL lbc_lnk( avmu, 'U', 1. ) ; CALL lbc_lnk( avmv, 'V', 1. ) ! Lateral boundary conditions 652 ! 671 653 ! 672 654 IF( nn_pdl == 1 ) THEN !* Prandtl number case: update avt … … 677 659 # if defined key_c1d 678 660 e_pdl(ji,jj,jk) = apdlr(ji,jj,jk) * wmask(ji,jj,jk) ! c1d configuration : save masked Prandlt number 679 !!gm bug NO zri here....680 !!gm remove the specific diag for c1d !681 e_ric(ji,jj,jk) = zri * wmask(ji,jj,jk) ! c1d config. : save Ri682 661 # endif 683 662 END DO … … 685 664 END DO 686 665 ENDIF 687 CALL lbc_lnk( avt, 'W', 1. ) ! Lateral boundary conditions on avt (sign unchanged) 666 !!gm I'm sure this is needed to compute the shear term at next time-step 667 CALL lbc_lnk( avm, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) 668 CALL lbc_lnk( avt, 'W', 1. ) ! Lateral boundary conditions on avt (sign unchanged) 688 669 689 670 IF(ln_ctl) THEN 690 CALL prt_ctl( tab3d_1=en , clinfo1=' tke - e: ', tab3d_2=avt, clinfo2=' t: ', ovlap=1, kdim=jpk) 691 CALL prt_ctl( tab3d_1=avmu, clinfo1=' tke - u: ', mask1=umask, & 692 & tab3d_2=avmv, clinfo2= ' v: ', mask2=vmask, ovlap=1, kdim=jpk ) 671 CALL prt_ctl( tab3d_1=en , clinfo1=' tke - e: ', tab3d_2=avt, clinfo2=' t: ', ovlap=1, kdim=jpk) 672 CALL prt_ctl( tab3d_1=avm, clinfo1=' tke - m: ', ovlap=1, kdim=jpk ) 693 673 ENDIF 694 674 ! … … 759 739 ENDIF 760 740 ! 761 IF( ln_zdftmx ) THEN ! Internal wave driven mixing 762 ! ! specific values of rn_emin & rmxl_min are used 763 rn_emin = 1.e-10_wp 764 rmxl_min = 1.e-03_wp 741 IF( ln_zdfiwm ) THEN ! Internal wave-driven mixing 742 rn_emin = 1.e-10_wp ! specific values of rn_emin & rmxl_min are used 743 rmxl_min = 1.e-03_wp ! associated avt minimum = molecular salt diffusivity (10^-9 m2/s) 765 744 IF(lwp) WRITE(numout,*) ' Internal wave-driven mixing case: force rn_emin = 1.e-10 and rmxl_min = 1.e-3 ' 766 ELSE 745 ELSE ! standard case : associated avt minimum = molecular viscosity (10^-6 m2/s) 767 746 rmxl_min = 1.e-6_wp / ( rn_ediff * SQRT( rn_emin ) ) ! resulting minimum length to recover molecular viscosity 768 747 IF(lwp) WRITE(numout,*) ' minimum mixing length with your parameters rmxl_min = ', rmxl_min … … 798 777 avt (:,:,jk) = avtb(jk) * wmask (:,:,jk) 799 778 avm (:,:,jk) = avmb(jk) * wmask (:,:,jk) 800 avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk)801 avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk)802 779 END DO 803 780 dissl(:,:,:) = 1.e-12_wp … … 821 798 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 822 799 ! 823 INTEGER :: jit, jk ! dummy loop indices824 INTEGER :: id1, id2, id3, id4 , id5, id6! local integers800 INTEGER :: jit, jk ! dummy loop indices 801 INTEGER :: id1, id2, id3, id4 ! local integers 825 802 !!---------------------------------------------------------------------- 826 803 ! … … 829 806 IF( ln_rstart ) THEN !* Read the restart file 830 807 id1 = iom_varid( numror, 'en' , ldstop = .FALSE. ) 831 id2 = iom_varid( numror, 'avt' , ldstop = .FALSE. ) 832 id3 = iom_varid( numror, 'avm' , ldstop = .FALSE. ) 833 id4 = iom_varid( numror, 'avmu' , ldstop = .FALSE. ) 834 id5 = iom_varid( numror, 'avmv' , ldstop = .FALSE. ) 835 id6 = iom_varid( numror, 'dissl', ldstop = .FALSE. ) 808 id2 = iom_varid( numror, 'avt_k', ldstop = .FALSE. ) 809 id3 = iom_varid( numror, 'avm_k', ldstop = .FALSE. ) 810 id4 = iom_varid( numror, 'dissl', ldstop = .FALSE. ) 836 811 ! 837 812 IF( id1 > 0 ) THEN ! 'en' exists 838 813 CALL iom_get( numror, jpdom_autoglo, 'en', en ) 839 IF( MIN( id2, id3, id4, id5, id6 ) > 0 ) THEN ! all required arrays exist 840 CALL iom_get( numror, jpdom_autoglo, 'avt' , avt ) 841 CALL iom_get( numror, jpdom_autoglo, 'avm' , avm ) 842 CALL iom_get( numror, jpdom_autoglo, 'avmu' , avmu ) 843 CALL iom_get( numror, jpdom_autoglo, 'avmv' , avmv ) 814 IF( MIN( id2, id3, id4 ) > 0 ) THEN ! all required arrays exist 815 CALL iom_get( numror, jpdom_autoglo, 'avt_k', avt_k ) 816 CALL iom_get( numror, jpdom_autoglo, 'avm_k', avm_k ) 844 817 CALL iom_get( numror, jpdom_autoglo, 'dissl', dissl ) 845 818 ELSE ! one at least array is missing 846 CALL tke_avn ! compute avt, avm, a vmu, avmv and dissl (approximation)819 CALL tke_avn ! compute avt, avm, and dissl (approximation) 847 820 ENDIF 848 821 ELSE ! No TKE array found: initialisation 849 822 IF(lwp) WRITE(numout,*) ' ===>>>> : previous run without tke scheme, en computed by iterative loop' 850 en (:,:,:) = rn_emin * tmask(:,:,:) 851 CALL tke_avn ! recompute avt, avm, avmu, avmv and dissl (approximation) 852 ! 853 avt_k (:,:,:) = avt (:,:,:) 854 avm_k (:,:,:) = avm (:,:,:) 855 avmu_k(:,:,:) = avmu(:,:,:) 856 avmv_k(:,:,:) = avmv(:,:,:) 823 en (:,:,:) = rn_emin * wmask(:,:,:) 824 CALL tke_avn ! recompute avt, avm, and dissl (approximation) 825 avt_k(:,:,:) = avt(:,:,:) 826 avm_k(:,:,:) = avm(:,:,:) 857 827 ! 858 828 DO jit = nit000 + 1, nit000 + 10 ; CALL zdf_tke( jit ) ; END DO 829 avt_k(:,:,:) = avt(:,:,:) 830 avm_k(:,:,:) = avm(:,:,:) 859 831 ENDIF 860 832 ELSE !* Start from rest 861 833 en(:,:,:) = rn_emin * tmask(:,:,:) 862 834 DO jk = 1, jpk ! set the Kz to the background value 863 avt (:,:,jk) = avtb(jk) * wmask (:,:,jk) 864 avm (:,:,jk) = avmb(jk) * wmask (:,:,jk) 865 avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk) 866 avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk) 835 avt_k(:,:,jk) = avtb(jk) * wmask (:,:,jk) 836 avm_k(:,:,jk) = avmb(jk) * wmask (:,:,jk) 867 837 END DO 868 838 ENDIF … … 871 841 ! ! ------------------- 872 842 IF(lwp) WRITE(numout,*) '---- tke-rst ----' 873 CALL iom_rstput( kt, nitrst, numrow, 'en' , en ) 874 CALL iom_rstput( kt, nitrst, numrow, 'avt' , avt_k ) 875 CALL iom_rstput( kt, nitrst, numrow, 'avm' , avm_k ) 876 CALL iom_rstput( kt, nitrst, numrow, 'avmu' , avmu_k ) 877 CALL iom_rstput( kt, nitrst, numrow, 'avmv' , avmv_k ) 878 CALL iom_rstput( kt, nitrst, numrow, 'dissl', dissl ) 843 CALL iom_rstput( kt, nitrst, numrow, 'en' , en ) 844 CALL iom_rstput( kt, nitrst, numrow, 'avt_k', avt_k ) 845 CALL iom_rstput( kt, nitrst, numrow, 'avm_k', avm_k ) 846 CALL iom_rstput( kt, nitrst, numrow, 'dissl', dissl ) 879 847 ! 880 848 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
