- Timestamp:
- 2017-03-17T10:44:05+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/r6232_HZG_WAVE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90
r6204 r7807 24 24 USE phycst ! physical constants 25 25 USE zdfmxl ! mixed layer 26 USE sbcwave 27 ! 26 28 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 27 29 USE lib_mpp ! MPP manager … … 47 49 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ustarb2 !: Squared bottom velocity scale at T-points 48 50 51 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rsbc_tke1 52 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rsbc_tke3 53 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rsbc_psi1 54 49 55 ! !! ** Namelist namzdf_gls ** 56 LOGICAL :: ln_crban = .FALSE. !! WAVE2NEMO added ln_crban (=.TRUE. use Craig and Banner scheme) 50 57 LOGICAL :: ln_length_lim ! use limit on the dissipation rate under stable stratification (Galperin et al. 1988) 51 58 LOGICAL :: ln_sigpsi ! Activate Burchard (2003) modification for k-eps closure & wave breaking mixing … … 59 66 REAL(wp) :: rn_emin ! minimum value of TKE (m2/s2) 60 67 REAL(wp) :: rn_charn ! Charnock constant for surface breaking waves mixing : 1400. (standard) or 2.e5 (Stacey value) 61 REAL(wp) :: rn_crban 68 REAL(wp) :: rn_crban_default ! Craig and Banner constant for surface breaking waves mixing 62 69 REAL(wp) :: rn_hsro ! Minimum surface roughness 63 REAL(wp) :: rn_frac_hs ! Fraction of wave height as surface roughness (if nn_z0_met >1)70 REAL(wp) :: rn_frac_hs ! Fraction of wave height as surface roughness (if nn_z0_met = 1) 64 71 65 72 REAL(wp) :: rcm_sf = 0.73_wp ! Shear free turbulence parameters … … 92 99 REAL(wp) :: rtrans = 0.1_wp 93 100 REAL(wp) :: rc02, rc02r, rc03, rc04 ! coefficients deduced from above parameters 94 REAL(wp) :: rsbc_tke1 , rsbc_tke2, rfact_tke! - - - -95 REAL(wp) :: rsbc_psi 1, rsbc_psi2, rfact_psi ! - - - -101 REAL(wp) :: rsbc_tke1_default, rsbc_tke2, rfact_tke ! - - - - 102 REAL(wp) :: rsbc_psi2, rsbc_psi3, rfact_psi ! - - - - 96 103 REAL(wp) :: rsbc_zs1, rsbc_zs2 ! - - - - 97 104 REAL(wp) :: rc0, rc2, rc3, rf6, rcff, rc_diff ! - - - - … … 115 122 !! *** FUNCTION zdf_gls_alloc *** 116 123 !!---------------------------------------------------------------------- 117 ALLOCATE( mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) , & 118 & ustars2(jpi,jpj) , ustarb2(jpi,jpj) , STAT= zdf_gls_alloc ) 124 ALLOCATE( mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) , & 125 & ustars2(jpi,jpj), ustarb2(jpi,jpj),rsbc_tke1(jpi,jpj), & 126 & rsbc_tke3(jpi,jpj), rsbc_psi1(jpi,jpj) , STAT=zdf_gls_alloc ) 119 127 ! 120 128 IF( lk_mpp ) CALL mpp_sum ( zdf_gls_alloc ) … … 161 169 ustars2 = 0._wp ; ustarb2 = 0._wp ; psi = 0._wp ; zwall_psi = 0._wp 162 170 171 172 ! variable initialization 173 IF( ln_crban ) THEN 174 rsbc_tke1(:,:) = (-rsc_tke*rn_crban(:,:)/(rcm_sf*ra_sf*rl_sf))**(2._wp/3._wp) ! k_eps = 53.Dirichlet + Wave breaking 175 rsbc_tke3(:,:) = rdt * rn_crban(:,:) ! Neumann + Wave breaking 176 rsbc_psi1(:,:) = rc0**rpp * rsbc_tke1(:,:)**rmm * rl_sf**rnn ! Dirichlet + Wave breaking 177 ENDIF 178 163 179 IF( kt /= nit000 ) THEN ! restore before value to compute tke 164 180 avt (:,:,:) = avt_k (:,:,:) … … 197 213 zdep(:,:) = 30.*TANH(2.*0.3/(28.*SQRT(MAX(ustars2(:,:),rsmall)))) ! Wave age (eq. 10) 198 214 zhsro(:,:) = MAX(rsbc_zs2 * ustars2(:,:) * zdep(:,:)**1.5, rn_hsro) ! zhsro = rn_frac_hs * Hsw (eq. 11) 199 ! 215 CASE ( 3 ) ! Roughness given by the wave model (coupled or read in file) 216 WHERE( swh_wavepar == 0._wp ) ! surface roughness length according to Charnock formula when sign. wave height 0 217 zhsro = MAX(rn_charn / grav * ustars2, rn_hsro) 218 ELSEWHERE 219 zhsro = MAX(swh_wavepar, rn_hsro) 220 END WHERE 200 221 END SELECT 201 222 … … 314 335 ! 315 336 CASE ( 0 ) ! Dirichlet case 316 ! First level 317 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1)**(2._wp/3._wp) 318 en(:,:,1) = MAX(en(:,:,1), rn_emin) 319 z_elem_a(:,:,1) = en(:,:,1) 320 z_elem_c(:,:,1) = 0._wp 321 z_elem_b(:,:,1) = 1._wp 322 ! 323 ! One level below 324 en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2)) & 325 & / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 326 en(:,:,2) = MAX(en(:,:,2), rn_emin ) 327 z_elem_a(:,:,2) = 0._wp 328 z_elem_c(:,:,2) = 0._wp 329 z_elem_b(:,:,2) = 1._wp 330 ! 331 ! 337 IF( ln_crban ) THEN ! wave induced mixing case with forced/coupled fields 338 ! First level 339 en(:,:,1) = MAX( rsbc_tke1(:,:) * ustars2(:,:), rn_emin ) 340 z_elem_a(:,:,1) = en(:,:,1) 341 z_elem_c(:,:,1) = 0._wp 342 z_elem_b(:,:,1) = 1._wp 343 344 ! One level below 345 en(:,:,2) = MAX( rsbc_tke1(:,:) * ustars2(:,:) * ((zhsro(:,:)+fsdepw(:,:,2))/zhsro(:,:) )**ra_sf, rn_emin ) 346 z_elem_a(:,:,2) = 0._wp 347 z_elem_c(:,:,2) = 0._wp 348 z_elem_b(:,:,2) = 1._wp 349 ELSE ! wave induced mixing case with default values 350 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default)**(2._wp/3._wp) 351 en(:,:,1) = MAX(en(:,:,1), rn_emin) 352 z_elem_a(:,:,1) = en(:,:,1) 353 z_elem_c(:,:,1) = 0._wp 354 z_elem_b(:,:,1) = 1._wp 355 ! 356 ! One level below 357 en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default * ((zhsro(:,:)+fsdepw(:,:,2)) & 358 & / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 359 en(:,:,2) = MAX(en(:,:,2), rn_emin ) 360 z_elem_a(:,:,2) = 0._wp 361 z_elem_c(:,:,2) = 0._wp 362 z_elem_b(:,:,2) = 1._wp 363 ! 364 ! 365 ENDIF 332 366 CASE ( 1 ) ! Neumann boundary condition on d(e)/dz 333 ! 334 ! Dirichlet conditions at k=1 335 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1)**(2._wp/3._wp) 336 en(:,:,1) = MAX(en(:,:,1), rn_emin) 337 z_elem_a(:,:,1) = en(:,:,1) 338 z_elem_c(:,:,1) = 0._wp 339 z_elem_b(:,:,1) = 1._wp 340 ! 341 ! at k=2, set de/dz=Fw 342 !cbr 343 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 344 z_elem_a(:,:,2) = 0._wp 345 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:)) )) 346 zflxs(:,:) = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) & 347 & * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 348 349 en(:,:,2) = en(:,:,2) + zflxs(:,:)/fse3w(:,:,2) 350 ! 351 ! 367 IF( ln_crban ) THEN ! Shear free case: d(e)/dz=Fw with forced/coupled fields 368 ! Dirichlet conditions at k=1 369 en(:,:,1) = MAX( rsbc_tke1(:,:) * ustars2(:,:), rn_emin ) 370 z_elem_a(:,:,1) = en(:,:,1) 371 z_elem_c(:,:,1) = 0._wp 372 z_elem_b(:,:,1) = 1._wp 373 ! at k=2, set de/dz=Fw 374 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 375 z_elem_a(:,:,2) = 0._wp 376 zflxs(:,:) = rsbc_tke3(:,:) * ustars2(:,:)**1.5_wp * ((zhsro(:,:)+fsdept(:,:,1) ) / zhsro(:,:) )**(1.5*ra_sf) 377 en(:,:,2) = en(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 378 ELSE ! Shear free case: d(e)/dz=Fw with default values 379 ! Dirichlet conditions at k=1 380 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default)**(2._wp/3._wp) 381 en(:,:,1) = MAX(en(:,:,1), rn_emin) 382 z_elem_a(:,:,1) = en(:,:,1) 383 z_elem_c(:,:,1) = 0._wp 384 z_elem_b(:,:,1) = 1._wp 385 ! 386 ! at k=2, set de/dz=Fw 387 !cbr 388 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 389 z_elem_a(:,:,2) = 0._wp 390 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:)) )) 391 zflxs(:,:) = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) & 392 & * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 393 394 en(:,:,2) = en(:,:,2) + zflxs(:,:)/fse3w(:,:,2) 395 ! 396 ! 397 ENDIF 352 398 END SELECT 353 399 … … 536 582 ! 537 583 CASE ( 0 ) ! Dirichlet boundary conditions 538 ! 539 ! Surface value 540 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 541 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 542 z_elem_a(:,:,1) = psi(:,:,1) 543 z_elem_c(:,:,1) = 0._wp 544 z_elem_b(:,:,1) = 1._wp 545 ! 546 ! One level below 547 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdepw(:,:,2)/zhsro(:,:) ))) 548 zdep(:,:) = (zhsro(:,:) + fsdepw(:,:,2)) * zkar(:,:) 549 psi (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 550 z_elem_a(:,:,2) = 0._wp 551 z_elem_c(:,:,2) = 0._wp 552 z_elem_b(:,:,2) = 1._wp 553 ! 554 ! 584 IF( ln_crban ) THEN ! Wave induced mixing case 585 ! en(1) = q2(1) = 0.5 * (15.8 * Ccb)^(2/3) * u*^2 586 ! balance between the production and the 587 ! dissipation terms including the wave effect 588 ! Surface value 589 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 590 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 591 z_elem_a(:,:,1) = psi(:,:,1) 592 z_elem_c(:,:,1) = 0._wp 593 z_elem_b(:,:,1) = 1._wp 594 ! 595 ! One level below 596 zex1 = (rmm*ra_sf+rnn) 597 zex2 = (rmm*ra_sf) 598 zdep(:,:) = ( (zhsro(:,:) + fsdepw(:,:,2))**zex1 ) / zhsro(:,:)**zex2 599 psi (:,:,2) = rsbc_psi1(:,:) * ustars2(:,:)**rmm * zdep(:,:) * tmask(:,:,1) 600 z_elem_a(:,:,2) = 0._wp 601 z_elem_c(:,:,2) = 0._wp 602 z_elem_b(:,:,2) = 1._wp 603 ! 604 ! 605 ELSE ! Wave induced mixing case with default values 606 ! Surface value 607 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 608 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 609 z_elem_a(:,:,1) = psi(:,:,1) 610 z_elem_c(:,:,1) = 0._wp 611 z_elem_b(:,:,1) = 1._wp 612 ! 613 ! One level below 614 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdepw(:,:,2)/zhsro(:,:) ))) 615 zdep(:,:) = (zhsro(:,:) + fsdepw(:,:,2)) * zkar(:,:) 616 psi (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 617 z_elem_a(:,:,2) = 0._wp 618 z_elem_c(:,:,2) = 0._wp 619 z_elem_b(:,:,2) = 1._wp 620 ! 621 ! 622 ENDIF 555 623 CASE ( 1 ) ! Neumann boundary condition on d(psi)/dz 556 ! 557 ! Surface value: Dirichlet 558 zdep(:,:) = zhsro(:,:) * rl_sf 559 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 560 z_elem_a(:,:,1) = psi(:,:,1) 561 z_elem_c(:,:,1) = 0._wp 562 z_elem_b(:,:,1) = 1._wp 563 ! 564 ! Neumann condition at k=2 565 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 566 z_elem_a(:,:,2) = 0._wp 567 ! 568 ! Set psi vertical flux at the surface: 569 zkar(:,:) = rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 570 zdep(:,:) = ((zhsro(:,:) + fsdept(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 571 zflxs(:,:) = (rnn + rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 572 zdep(:,:) = rsbc_psi1 * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * & 573 & ustars2(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + fsdept(:,:,1))**(rnn-1.) 574 zflxs(:,:) = zdep(:,:) * zflxs(:,:) 575 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 576 577 ! 578 ! 624 IF( ln_crban ) THEN ! Wave induced mixing case with forced/coupled fields 625 ! 626 zdep(:,:) = rl_sf * zhsro(:,:) 627 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 628 z_elem_a(:,:,1) = psi(:,:,1) 629 z_elem_c(:,:,1) = 0._wp 630 z_elem_b(:,:,1) = 1._wp 631 ! 632 ! Neumann condition at k=2 633 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 634 z_elem_a(:,:,2) = 0._wp 635 ! 636 ! Set psi vertical flux at the surface: 637 zdep(:,:) = (zhsro(:,:) + fsdept(:,:,1))**(rmm*ra_sf+rnn-1._wp) / zhsro(:,:)**(rmm*ra_sf) 638 zflxs(:,:) = rsbc_psi3 * ( zwall_psi(:,:,1)*avm(:,:,1) + zwall_psi(:,:,2)*avm(:,:,2) ) & 639 & * en(:,:,1)**rmm * zdep 640 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 641 ! 642 ELSE ! Wave induced mixing case with default values 643 ! Surface value: Dirichlet 644 zdep(:,:) = zhsro(:,:) * rl_sf 645 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 646 z_elem_a(:,:,1) = psi(:,:,1) 647 z_elem_c(:,:,1) = 0._wp 648 z_elem_b(:,:,1) = 1._wp 649 ! 650 ! Neumann condition at k=2 651 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 652 z_elem_a(:,:,2) = 0._wp 653 ! 654 ! Set psi vertical flux at the surface: 655 zkar(:,:) = rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 656 zdep(:,:) = ((zhsro(:,:) + fsdept(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 657 zflxs(:,:) = (rnn + rsbc_tke1_default * (rnn + rmm*ra_sf) * zdep(:,:)) * & 658 (1._wp + rsbc_tke1_default * zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 659 zdep(:,:) = rsbc_psi1(:,:) * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * & 660 & ustars2(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + fsdept(:,:,1))**(rnn-1.) 661 zflxs(:,:) = zdep(:,:) * zflxs(:,:) 662 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 663 ! 664 ! 665 ENDIF 579 666 END SELECT 580 667 … … 866 953 !! 867 954 NAMELIST/namzdf_gls/rn_emin, rn_epsmin, ln_length_lim, & 868 & rn_clim_galp, ln_ sigpsi, rn_hsro, &869 & rn_crban , rn_charn, rn_frac_hs,&955 & rn_clim_galp, ln_crban, ln_sigpsi, rn_hsro, & 956 & rn_crban_default, rn_charn, rn_frac_hs,& 870 957 & nn_bc_surf, nn_bc_bot, nn_z0_met, & 871 958 & nn_stab_func, nn_clos … … 895 982 WRITE(numout,*) ' TKE Bottom boundary condition nn_bc_bot = ', nn_bc_bot 896 983 WRITE(numout,*) ' Modify psi Schmidt number (wb case) ln_sigpsi = ', ln_sigpsi 897 WRITE(numout,*) ' Craig and Banner coefficient rn_crban = ', rn_crban 984 WRITE(numout,*) ' Craig and Banner coefficient (default) rn_crban = ', rn_crban_default 985 WRITE(numout,*) ' Craig and Banner scheme ln_crban = ', ln_crban 898 986 WRITE(numout,*) ' Charnock coefficient rn_charn = ', rn_charn 899 987 WRITE(numout,*) ' Surface roughness formula nn_z0_met = ', nn_z0_met … … 909 997 910 998 ! !* Check of some namelist values 911 IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'bad flag: nn_bc_surf is 0 or 1' ) 912 IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'bad flag: nn_bc_surf is 0 or 1' ) 913 IF( nn_z0_met < 0 .OR. nn_z0_met > 2 ) CALL ctl_stop( 'bad flag: nn_z0_met is 0, 1 or 2' ) 914 IF( nn_stab_func < 0 .OR. nn_stab_func > 3 ) CALL ctl_stop( 'bad flag: nn_stab_func is 0, 1, 2 and 3' ) 915 IF( nn_clos < 0 .OR. nn_clos > 3 ) CALL ctl_stop( 'bad flag: nn_clos is 0, 1, 2 or 3' ) 999 IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_bc_surf is 0 or 1' ) 1000 IF( nn_bc_surf < 0 .OR. nn_bc_surf > 1 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_bc_surf is 0 or 1' ) 1001 IF( nn_z0_met < 0 .OR. nn_z0_met > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_z0_met is 0, 1, 2 or 3' ) 1002 ! IF( nn_z0_met == 3 .AND. .NOT.ln_sdw ) CALL ctl_stop( 'zdf_gls_init: nn_z0_met=3 requires ln_sdw=T' ) 1003 IF( nn_stab_func < 0 .OR. nn_stab_func > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_stab_func is 0, 1, 2 and 3' ) 1004 IF( nn_clos < 0 .OR. nn_clos > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_clos is 0, 1, 2 or 3' ) 916 1005 917 1006 SELECT CASE ( nn_clos ) !* set the parameters for the chosen closure … … 1085 1174 & - SQRT(rsc_tke + 24._wp*rsc_psi0*rpsi2 ) ) 1086 1175 1087 IF ( rn_crban==0._wp ) THEN1176 IF( .NOT. ln_crban .AND. rn_crban_default==0._wp ) THEN 1088 1177 rl_sf = vkarmn 1089 1178 ELSE … … 1124 1213 rc03 = rc02 * rc0 1125 1214 rc04 = rc03 * rc0 1126 rsbc_tke1 = -3._wp/2._wp*rn_crban*ra_sf*rl_sf ! Dirichlet + Wave breaking1127 rsbc_tke2 = rdt * rn_crban / rl_sf ! Neumann + Wave breaking1128 zcr = MAX(rsmall, rsbc_tke1**(1./(-ra_sf*3._wp/2._wp))-1._wp )1129 rtrans = 0.2_wp / zcr ! Ad. inverse transition length between log and wave layer1215 rsbc_tke1_default = -3._wp/2._wp*rn_crban_default*ra_sf*rl_sf 1216 rsbc_tke2 = rdt * rn_crban_default / rl_sf 1217 zcr = MAX( rsmall, rsbc_tke1_default**(1./(-ra_sf*3._wp/2._wp))-1._wp ) 1218 rtrans = 0.2_wp / zcr 1130 1219 rsbc_zs1 = rn_charn/grav ! Charnock formula for surface roughness 1131 1220 rsbc_zs2 = rn_frac_hs / 0.85_wp / grav * 665._wp ! Rascle formula for surface roughness
Note: See TracChangeset
for help on using the changeset viewer.