- Timestamp:
- 2019-07-31T18:05:50+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_STOPACK/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
r5407 r11384 53 53 USE timing ! Timing 54 54 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 55 USE stopack 55 56 56 57 IMPLICIT NONE … … 85 86 REAL(wp) :: rhftau_scl = 1.0_wp ! scale factor applied to HF part of taum (nn_etau=3) 86 87 87 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: en !: now turbulent kinetic energy [m2/s2] 88 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: rn_lc0, rn_ediff0, rn_ediss0, rn_ebb0, rn_efr0 89 88 90 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:) :: htau ! depth of tke penetration (nn_htau) 89 91 REAL(wp) , ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dissl ! now mixing lenght of dissipation 90 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: avt_k , avm_k ! not enhanced Kz 91 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: avmu_k, avmv_k ! not enhanced Kz 92 92 93 #if defined key_c1d 93 94 ! !!** 1D cfg only ** ('key_c1d') … … 115 116 & e_pdl(jpi,jpj,jpk) , e_ric(jpi,jpj,jpk) , & 116 117 #endif 117 & en (jpi,jpj,jpk) , htau (jpi,jpj) , dissl(jpi,jpj,jpk) , & 118 & avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk), & 119 & avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk), STAT= zdf_tke_alloc ) 118 & htau (jpi,jpj) , dissl(jpi,jpj,jpk) , STAT= zdf_tke_alloc ) 120 119 ! 121 120 IF( lk_mpp ) CALL mpp_sum ( zdf_tke_alloc ) … … 178 177 avmu(:,:,:) = avmu_k(:,:,:) 179 178 avmv(:,:,:) = avmv_k(:,:,:) 180 ENDIF 179 ENDIF 180 ! 181 IF( nn_spp_tkelc > 0 ) THEN 182 rn_lc0 = rn_lc 183 CALL spp_gen(kt,rn_lc0,nn_spp_tkelc,rn_tkelc_sd, jk_spp_tkelc ) 184 ENDIF 185 IF( nn_spp_tkedf > 0 ) THEN 186 rn_ediff0 = rn_ediff 187 CALL spp_gen(kt,rn_ediff0,nn_spp_tkedf,rn_tkedf_sd, jk_spp_tkedf ) 188 ENDIF 189 IF( nn_spp_tkeds > 0 ) THEN 190 rn_ediss0 = rn_ediss 191 CALL spp_gen(kt,rn_ediss0,nn_spp_tkeds,rn_tkeds_sd, jk_spp_tkeds ) 192 ENDIF 193 IF( nn_spp_tkebb > 0 ) THEN 194 rn_ebb0 = rn_ebb 195 CALL spp_gen(kt,rn_ebb0,nn_spp_tkebb,rn_tkebb_sd, jk_spp_tkebb ) 196 ENDIF 197 IF( nn_spp_tkefr > 0 ) THEN 198 rn_efr0 = rn_efr 199 CALL spp_gen(kt,rn_efr0,nn_spp_tkefr,rn_tkefr_sd, jk_spp_tkefr ) 200 ENDIF 181 201 ! 182 202 CALL tke_tke ! now tke (en) … … 188 208 avmu_k(:,:,:) = avmu(:,:,:) 189 209 avmv_k(:,:,:) = avmv(:,:,:) 210 ! 211 IF ( kt .eq. nitend ) THEN 212 DEALLOCATE ( rn_lc0, rn_ediff0, rn_ediss0, rn_ebb0, rn_efr0 ) 213 ENDIF 190 214 ! 191 215 END SUBROUTINE zdf_tke … … 214 238 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 215 239 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 216 REAL(wp) :: z bbrau, zesh2! temporary scalars217 REAL(wp) :: zfact1 , zfact2, zfact3! - -240 REAL(wp) :: zesh2 ! temporary scalars 241 REAL(wp) :: zfact1 ! - - 218 242 REAL(wp) :: ztx2 , zty2 , zcof ! - - 219 243 REAL(wp) :: ztau , zdif ! - - … … 222 246 !!bfr REAL(wp) :: zebot ! - - 223 247 INTEGER , POINTER, DIMENSION(:,: ) :: imlc 224 REAL(wp), POINTER, DIMENSION(:,: ) :: zhlc 248 REAL(wp), POINTER, DIMENSION(:,: ) :: zhlc, zbbrau,zfact2,zfact3 225 249 REAL(wp), POINTER, DIMENSION(:,:,:) :: zpelc, zdiag, zd_up, zd_lw 226 250 !!-------------------------------------------------------------------- … … 229 253 ! 230 254 CALL wrk_alloc( jpi,jpj, imlc ) ! integer 231 CALL wrk_alloc( jpi,jpj, zhlc ) 255 CALL wrk_alloc( jpi,jpj, zhlc ) 256 CALL wrk_alloc( jpi,jpj, zbbrau ) 257 CALL wrk_alloc( jpi,jpj, zfact2 ) 258 CALL wrk_alloc( jpi,jpj, zfact3 ) 232 259 CALL wrk_alloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw ) 233 260 ! 234 zbbrau = rn_ebb / rau0 ! Local constant initialisation261 zbbrau = rn_ebb0 / rau0 ! Local constant initialisation 235 262 zfact1 = -.5_wp * rdt 236 zfact2 = 1.5_wp * rdt * rn_ediss 237 zfact3 = 0.5_wp * rn_ediss 263 zfact2 = 1.5_wp * rdt * rn_ediss0 264 zfact3 = 0.5_wp * rn_ediss0 238 265 ! 239 266 ! … … 250 277 DO jj = 2, jpjm1 ! en(1) = rn_ebb taum / rau0 (min value rn_emin0) 251 278 DO ji = fs_2, fs_jpim1 ! vector opt. 252 en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1)279 en(ji,jj,1) = MAX( rn_emin0, zbbrau(ji,jj) * taum(ji,jj) ) * tmask(ji,jj,1) 253 280 END DO 254 281 END DO … … 315 342 ! ! vertical velocity due to LC 316 343 zind = 0.5 - SIGN( 0.5, fsdepw(ji,jj,jk) - zhlc(ji,jj) ) 317 zwlc = zind * rn_lc * zus * SIN( rpi * fsdepw(ji,jj,jk) / zhlc(ji,jj) )344 zwlc = zind * rn_lc0(ji,jj) * zus * SIN( rpi * fsdepw(ji,jj,jk) / zhlc(ji,jj) ) 318 345 ! ! TKE Langmuir circulation source term 319 en(ji,jj,jk) = en(ji,jj,jk) + rdt * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1) 346 en(ji,jj,jk) = en(ji,jj,jk) + rdt * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc ) / & 347 & zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1) 348 320 349 END DO 321 350 END DO … … 360 389 zd_up(ji,jj,jk) = zzd_up ! Matrix (zdiag, zd_up, zd_lw) 361 390 zd_lw(ji,jj,jk) = zzd_lw 362 zdiag(ji,jj,jk) = 1._wp - zzd_lw - zzd_up + zfact2 * dissl(ji,jj,jk) * tmask(ji,jj,jk)391 zdiag(ji,jj,jk) = 1._wp - zzd_lw - zzd_up + zfact2(ji,jj) * dissl(ji,jj,jk) * tmask(ji,jj,jk) 363 392 ! 364 393 ! ! right hand side in en 365 394 en(ji,jj,jk) = en(ji,jj,jk) + rdt * ( zesh2 - avt(ji,jj,jk) * rn2(ji,jj,jk) & 366 & + zfact3 * dissl(ji,jj,jk) * en (ji,jj,jk) ) &395 & + zfact3(ji,jj) * dissl(ji,jj,jk) * en (ji,jj,jk) ) & 367 396 & * wmask(ji,jj,jk) 368 397 END DO … … 420 449 DO jj = 2, jpjm1 421 450 DO ji = fs_2, fs_jpim1 ! vector opt. 422 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) &423 & * ( 1._wp -fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1)451 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr0(ji,jj) * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 452 & * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 424 453 END DO 425 454 END DO … … 429 458 DO ji = fs_2, fs_jpim1 ! vector opt. 430 459 jk = nmln(ji,jj) 431 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) &432 & * ( 1._wp -fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1)460 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr0(ji,jj) * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 461 & * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 433 462 END DO 434 463 END DO … … 445 474 zdif = taum(ji,jj) - ztau ! mean of modulus - modulus of the mean 446 475 zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add ) ! apply some modifications... 447 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) &448 & * ( 1._wp -fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1)476 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau(ji,jj) * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 477 & * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 449 478 END DO 450 479 END DO … … 455 484 CALL wrk_dealloc( jpi,jpj, imlc ) ! integer 456 485 CALL wrk_dealloc( jpi,jpj, zhlc ) 486 CALL wrk_dealloc( jpi,jpj, zbbrau ) 487 CALL wrk_dealloc( jpi,jpj, zfact2 ) 488 CALL wrk_dealloc( jpi,jpj, zfact3 ) 457 489 CALL wrk_dealloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw ) 458 490 ! … … 637 669 DO ji = fs_2, fs_jpim1 ! vector opt. 638 670 zsqen = SQRT( en(ji,jj,jk) ) 639 zav = rn_ediff * zmxlm(ji,jj,jk) * zsqen671 zav = rn_ediff0(ji,jj) * zmxlm(ji,jj,jk) * zsqen 640 672 avm (ji,jj,jk) = MAX( zav, avmb(jk) ) * wmask(ji,jj,jk) 641 673 avt (ji,jj,jk) = MAX( zav, avtb_2d(ji,jj) * avtb(jk) ) * wmask(ji,jj,jk) … … 643 675 END DO 644 676 END DO 677 IF(nn_spp_avt > 0 ) CALL spp_gen(1 ,avt(:,:,jk),nn_spp_avt,rn_avt_sd, jk_spp_avt, jk) 678 IF(nn_spp_avm > 0 ) CALL spp_gen(1 ,avm(:,:,jk),nn_spp_avm,rn_avm_sd, jk_spp_avm, jk) 645 679 END DO 646 680 CALL lbc_lnk( avm, 'W', 1. ) ! Lateral boundary conditions (sign unchanged) … … 710 744 !!---------------------------------------------------------------------- 711 745 INTEGER :: ji, jj, jk ! dummy loop indices 712 INTEGER :: ios 746 INTEGER :: ios, ierr 713 747 !! 714 748 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & … … 767 801 rn_mxl0 = rmxl_min 768 802 ENDIF 803 804 ALLOCATE( rn_lc0 (jpi,jpj) ) ; rn_lc0 = rn_lc 805 ALLOCATE( rn_ediff0(jpi,jpj) ) ; rn_ediff0 = rn_ediff 806 ALLOCATE( rn_ediss0(jpi,jpj) ) ; rn_ediss0 = rn_ediss 807 ALLOCATE( rn_ebb0 (jpi,jpj) ) ; rn_ebb0 = rn_ebb 808 ALLOCATE( rn_efr0 (jpi,jpj) ) ; rn_efr0 = rn_efr 769 809 770 IF( nn_etau == 2 ) CALL zdf_mxl( nit000 ) ! Initialization of nmln 810 IF( nn_etau == 2 ) THEN 811 ierr = zdf_mxl_alloc() 812 nmln(:,:) = nlb10 ! Initialization of nmln 813 ENDIF 771 814 772 815 ! !* depth of penetration of surface tke … … 836 879 CALL tke_avn ! recompute avt, avm, avmu, avmv and dissl (approximation) 837 880 ! 838 avt_k (:,:,:) = avt (:,:,:)839 avm_k (:,:,:) = avm (:,:,:)840 avmu_k(:,:,:) = avmu(:,:,:)841 avmv_k(:,:,:) = avmv(:,:,:)842 !843 881 DO jit = nit000 + 1, nit000 + 10 ; CALL zdf_tke( jit ) ; END DO 844 882 ENDIF 845 883 ELSE !* Start from rest 846 884 en(:,:,:) = rn_emin * tmask(:,:,:) 847 DO jk = 1, jpk ! set the Kz to the background value848 avt (:,:,jk) = avtb(jk) * wmask (:,:,jk)849 avm (:,:,jk) = avmb(jk) * wmask (:,:,jk)850 avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk)851 avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk)852 END DO853 885 ENDIF 854 886 ! 887 avt_k (:,:,:) = avt (:,:,:) 888 avm_k (:,:,:) = avm (:,:,:) 889 avmu_k(:,:,:) = avmu(:,:,:) 890 avmv_k(:,:,:) = avmv(:,:,:) 891 855 892 ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN ! Create restart file 856 893 ! ! -------------------
Note: See TracChangeset
for help on using the changeset viewer.