Changeset 1705 for trunk/NEMO/OPA_SRC/ZDF/zdftke.F90
- Timestamp:
- 2009-11-03T17:37:00+01:00 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/ZDF/zdftke.F90
r1695 r1705 65 65 66 66 ! !!! ** Namelist namzdf_tke ** 67 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 68 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 69 REAL(wp) :: rn_lmin0 = 0.4_wp ! surface min value of mixing length [m] 70 REAL(wp) :: rn_lmin = 0.1_wp ! interior min value of mixing length [m] 71 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 72 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 73 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 74 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 75 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 76 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 77 REAL(wp) :: rn_bshear= 1.e-20 ! background shear (>0) 78 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2) 79 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 80 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 81 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 82 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 67 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 68 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 69 REAL(wp) :: rn_lmin0 = 0.4_wp ! surface min value of mixing length [m] 70 REAL(wp) :: rn_lmin = 0.1_wp ! interior min value of mixing length [m] 71 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 72 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 73 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 74 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 75 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 76 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 77 REAL(wp) :: rn_bshear = 1.e-20 ! background shear (>0) 78 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2) 79 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 80 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 81 REAL(wp) :: rn_addhft = 0.0_wp ! add offset applied to HF tau 82 REAL(wp) :: rn_sclhft = 1.0_wp ! scale factor applied to HF tau 83 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 84 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 83 85 84 86 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) … … 177 179 USE oce, zd_lw => ta ! use ta as workspace 178 180 !! 179 INTEGER :: ji, jj, jk ! dummy loop arguments 180 REAL(wp) :: zbbrau, zesh2 ! temporary scalars 181 REAL(wp) :: zfact1, zfact2, zfact3 ! - - 182 REAL(wp) :: ztx2 , zty2 , zcof ! - - 183 REAL(wp) :: zus , zwlc , zind ! - - 184 REAL(wp) :: zzd_up, zzd_lw ! - - 181 INTEGER :: ji, jj, jk ! dummy loop arguments 185 182 INTEGER :: ikbu, ikbv, ikbum1, ikbvm1 ! temporary scalar 186 183 INTEGER :: ikbt, ikbumm1, ikbvmm1 ! temporary scalar 187 REAL(wp) :: zebot ! temporary scalars 184 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 185 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 186 REAL(wp) :: zbbrau, zesh2 ! temporary scalars 187 REAL(wp) :: zfact1, zfact2, zfact3 ! - - 188 REAL(wp) :: ztx2 , zty2 , zcof ! - - 189 REAL(wp) :: ztau , zdif ! - - 190 REAL(wp) :: zus , zwlc , zind ! - - 191 REAL(wp) :: zzd_up, zzd_lw ! - - 192 REAL(wp) :: zebot ! - - 188 193 INTEGER , DIMENSION(jpi,jpj) :: imlc ! 2D workspace 189 194 REAL(wp), DIMENSION(jpi,jpj) :: zhlc ! - - … … 242 247 END DO 243 248 ! !* finite Langmuir Circulation depth 249 zcof = 0.5 * 0.016 * 0.016 / ( zrhoa * zcdrag ) 244 250 imlc(:,:) = mbathy(:,:) ! Initialization to the number of w ocean point mbathy 245 251 DO jk = jpkm1, 2, -1 246 252 DO jj = 1, jpj ! Last w-level at which zpelc>=0.5*us*us 247 253 DO ji = 1, jpi ! with us=0.016*wind(starting from jpk-1) 248 zus = 0.000128 * wndm(ji,jj) * wndm(ji,jj)254 zus = zcof * taum(ji,jj) 249 255 IF( zpelc(ji,jj,jk) > zus ) imlc(ji,jj) = jk 250 256 END DO … … 265 271 hlc(:,:) = zhlc(:,:) * tmask(:,:,1) ! c1d configuration: save finite Langmuir Circulation depth 266 272 # endif 273 zcof = 0.016 / SQRT( zrhoa * zcdrag ) 267 274 !CDIR NOVERRCHK 268 275 DO jk = 2, jpkm1 !* TKE Langmuir circulation source term added to en … … 271 278 !CDIR NOVERRCHK 272 279 DO ji = fs_2, fs_jpim1 ! vector opt. 273 !!gm replace here wndn by a formulation with the stress module 274 zus = 0.016 * wndm(ji,jj) ! Stokes drift 280 zus = zcof * SQRT( taum(ji,jj) ) ! Stokes drift 275 281 ! ! vertical velocity due to LC 276 282 zind = 0.5 - SIGN( 0.5, fsdepw(ji,jj,jk) - zhlc(ji,jj) ) … … 371 377 ! ! TKE due to surface and internal wave breaking 372 378 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 373 IF( nn_etau == 1 ) THEN !* 379 IF( nn_etau == 1 ) THEN !* penetration throughout the water column 374 380 DO jk = 2, jpkm1 375 381 DO jj = 2, jpjm1 … … 386 392 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 387 393 & * ( 1.e0 - fr_i(ji,jj) ) * tmask(ji,jj,jk) 394 END DO 395 END DO 396 ELSEIF( nn_etau == 3 ) THEN !* penetration throughout the water column 397 !CDIR NOVERRCHK 398 DO jk = 2, jpkm1 399 !CDIR NOVERRCHK 400 DO jj = 2, jpjm1 401 !CDIR NOVERRCHK 402 DO ji = fs_2, fs_jpim1 ! vector opt. 403 ztx2 = utau(ji-1,jj ) + utau(ji,jj) 404 zty2 = vtau(ji ,jj-1) + vtau(ji,jj) 405 ztau = 0.5 * SQRT( ztx2 * ztx2 + zty2 * zty2 ) ! module of the mean stress 406 zdif = taum(ji,jj) - ztau ! mean of the module - module of the mean 407 zdif = rn_sclhft * MAX( 0.e0, zdif + rn_addhft ) ! apply some modifications... 408 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 409 & * ( 1.e0 - fr_i(ji,jj) ) * tmask(ji,jj,jk) 410 END DO 388 411 END DO 389 412 END DO … … 642 665 INTEGER :: ji, jj, jk ! dummy loop indices 643 666 !! 644 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb, rn_emin, & 645 & rn_emin0, rn_bshear, nn_mxl, ln_mxl0, & 646 & rn_lmin , rn_lmin0 , nn_pdl, nn_etau, & 647 & nn_htau , rn_efr , ln_lc , rn_lc 667 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 668 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 669 & rn_lmin , rn_lmin0 , nn_pdl , nn_etau , & 670 & nn_htau , rn_efr , rn_addhft, rn_sclhft, & 671 & ln_lc , rn_lc 648 672 !!---------------------------------------------------------------------- 649 673 … … 658 682 WRITE(numout,*) '~~~~~~~~' 659 683 WRITE(numout,*) ' Namelist namzdf_tke : set tke mixing parameters' 660 WRITE(numout,*) ' coef. to compute avt rn_ediff = ', rn_ediff 661 WRITE(numout,*) ' Kolmogoroff dissipation coef. rn_ediss = ', rn_ediss 662 WRITE(numout,*) ' tke surface input coef. rn_ebb = ', rn_ebb 663 WRITE(numout,*) ' minimum value of tke rn_emin = ', rn_emin 664 WRITE(numout,*) ' surface minimum value of tke rn_emin0 = ', rn_emin0 665 WRITE(numout,*) ' background shear (>0) rn_bshear= ', rn_bshear 666 WRITE(numout,*) ' mixing length type nn_mxl = ', nn_mxl 667 WRITE(numout,*) ' prandl number flag nn_pdl = ', nn_pdl 668 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 669 WRITE(numout,*) ' surface mixing length minimum value rn_lmin0 = ', rn_lmin0 670 WRITE(numout,*) ' interior mixing length minimum value rn_lmin0 = ', rn_lmin0 671 WRITE(numout,*) ' test param. to add tke induced by wind nn_etau = ', nn_etau 672 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 673 WRITE(numout,*) ' % of rn_emin0 which pene. the thermocline rn_efr = ', rn_efr 674 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc 675 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc 684 WRITE(numout,*) ' coef. to compute avt rn_ediff = ', rn_ediff 685 WRITE(numout,*) ' Kolmogoroff dissipation coef. rn_ediss = ', rn_ediss 686 WRITE(numout,*) ' tke surface input coef. rn_ebb = ', rn_ebb 687 WRITE(numout,*) ' minimum value of tke rn_emin = ', rn_emin 688 WRITE(numout,*) ' surface minimum value of tke rn_emin0 = ', rn_emin0 689 WRITE(numout,*) ' background shear (>0) rn_bshear = ', rn_bshear 690 WRITE(numout,*) ' mixing length type nn_mxl = ', nn_mxl 691 WRITE(numout,*) ' prandl number flag nn_pdl = ', nn_pdl 692 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 693 WRITE(numout,*) ' surface mixing length minimum value rn_lmin0 = ', rn_lmin0 694 WRITE(numout,*) ' interior mixing length minimum value rn_lmin0 = ', rn_lmin0 695 WRITE(numout,*) ' test param. to add tke induced by wind nn_etau = ', nn_etau 696 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 697 WRITE(numout,*) ' fraction of en which pene. the thermocline rn_efr = ', rn_efr 698 WRITE(numout,*) ' add offset applied to HF tau rn_addhft = ', rn_addhft 699 WRITE(numout,*) ' scale factor applied to HF tau rn_sclhft = ', rn_sclhft 700 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc 701 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc 676 702 WRITE(numout,*) 677 703 WRITE(numout,*) ' critical Richardson nb with your parameters ri_cri = ', ri_cri … … 679 705 680 706 ! !* Check of some namelist values 681 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 682 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 683 IF( nn_htau < 0 .OR. nn_htau > 1 ) CALL ctl_stop( 'bad flag: nn_htau is 0 or 1 ' ) 684 IF( rn_lc < 0.15 .OR. rn_lc > 0.2 ) CALL ctl_stop( 'bad value: rn_lc must be between 0.15 and 0.2 ' ) 707 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 708 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 709 IF( nn_htau < 0 .OR. nn_htau > 2 ) CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' ) 710 IF( rn_lc < 0.15 .OR. rn_lc > 0.2 ) CALL ctl_stop( 'bad value: rn_lc must be between 0.15 and 0.2 ' ) 711 IF( rn_sclhft == 0. .AND. nn_etau == 3 ) CALL ctl_stop( 'force null HF tau to penetrate the thermocline...' ) 712 IF( .NOT. lhftau .AND. nn_etau == 3 ) CALL ctl_stop( 'bad flag: nn_etau == 3 must be used with HF tau' ) 685 713 686 714 IF( nn_etau == 2 ) CALL zdf_mxl( nit000 ) ! Initialization of nmln … … 697 725 END DO 698 726 END DO 727 CASE( 2 ) ! constant depth penetration (here 30 meters) 728 htau(:,:) = 30.e0 699 729 END SELECT 700 730 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.