Changeset 10478 for branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF
- Timestamp:
- 2019-01-09T12:11:55+01:00 (5 years ago)
- Location:
- branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90
r10473 r10478 10 10 USE in_out_manager ! I/O manager 11 11 USE lib_mpp ! MPP library 12 USE sbc_oce, ONLY : ln_zdfqiao 12 13 13 14 IMPLICIT NONE … … 35 36 INTEGER , PUBLIC :: nn_npc !: non penetrative convective scheme call frequency 36 37 INTEGER , PUBLIC :: nn_npcp !: non penetrative convective scheme print frequency 37 LOGICAL , PUBLIC :: ln_zdfqiao !: Enhanced wave vertical mixing Qiao(2010) formulation flag38 38 39 39 -
branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90
r10473 r10478 24 24 USE phycst ! physical constants 25 25 USE zdfmxl ! mixed layer 26 USE sbcwave , ONLY: hsw ! significant wave height26 USE sbcwave 27 27 ! 28 28 USE lbclnk ! ocean lateral boundary conditions (or mpp link) … … 48 48 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ustars2 !: Squared surface velocity scale at T-points 49 49 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ustarb2 !: Squared bottom velocity scale at T-points 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 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: rtrans 50 55 51 56 ! !! ** Namelist namzdf_gls ** … … 61 66 REAL(wp) :: rn_emin ! minimum value of TKE (m2/s2) 62 67 REAL(wp) :: rn_charn ! Charnock constant for surface breaking waves mixing : 1400. (standard) or 2.e5 (Stacey value) 63 REAL(wp) :: rn_crban 68 REAL(wp) :: rn_crban_default ! Craig and Banner constant for surface breaking waves mixing 64 69 REAL(wp) :: rn_hsro ! Minimum surface roughness 65 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, 3) 66 71 67 72 REAL(wp) :: rcm_sf = 0.73_wp ! Shear free turbulence parameters … … 92 97 REAL(wp) :: rm7 = 0.0_wp 93 98 REAL(wp) :: rm8 = 0.318_wp 94 REAL(wp) :: rtrans 99 REAL(wp) :: rtrans_default = 0.1_wp 95 100 REAL(wp) :: rc02, rc02r, rc03, rc04 ! coefficients deduced from above parameters 96 REAL(wp) :: rsbc_tke1 , rsbc_tke2, rfact_tke ! - - - -97 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 ! - - - - 98 103 REAL(wp) :: rsbc_zs1, rsbc_zs2 ! - - - - 99 104 REAL(wp) :: rc0, rc2, rc3, rf6, rcff, rc_diff ! - - - - … … 117 122 !! *** FUNCTION zdf_gls_alloc *** 118 123 !!---------------------------------------------------------------------- 119 ALLOCATE( mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) , & 120 & 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) , & 126 & rsbc_tke1(jpi,jpj), rsbc_tke3(jpi,jpj) , & 127 & rsbc_psi1(jpi,jpj), rtrans(jpi,jpj), STAT= zdf_gls_alloc ) 121 128 ! 122 129 IF( lk_mpp ) CALL mpp_sum ( zdf_gls_alloc ) … … 163 170 ustars2 = 0._wp ; ustarb2 = 0._wp ; psi = 0._wp ; zwall_psi = 0._wp 164 171 172 ! variable initialization 173 IF( ln_phioc ) THEN 174 IF( nn_wmix==jp_janssen ) THEN 175 rsbc_tke1(:,:) = (-rsc_tke*rn_crban(:,:)/(rcm_sf*ra_sf*rl_sf))**(2._wp/3._wp) ! k_eps = 53.Dirichlet + Wave breaking 176 rsbc_tke3(:,:) = rdt * rn_crban(:,:) ! Neumann + Wave breaking 177 rsbc_psi1(:,:) = rc0**rpp * rsbc_tke1(:,:)**rmm * rl_sf**rnn ! Dirichlet + Wave breaking 178 ELSE IF( nn_wmix==jp_craigbanner ) THEN 179 rsbc_tke1(:,:) = -3._wp/2._wp*rn_crban(:,:)*ra_sf*rl_sf 180 rsbc_tke3(:,:) = rdt * rn_crban(:,:) / rl_sf 181 rtrans(:,:) = 0.2_wp / MAX( rsmall, rsbc_tke1(:,:)**(1./(-ra_sf*3._wp/2._wp))-1._wp ) 182 ENDIF 183 ENDIF 184 165 185 IF( kt /= nit000 ) THEN ! restore before value to compute tke 166 186 avt (:,:,:) = avt_k (:,:,:) … … 200 220 zhsro(:,:) = MAX(rsbc_zs2 * ustars2(:,:) * zdep(:,:)**1.5, rn_hsro) ! zhsro = rn_frac_hs * Hsw (eq. 11) 201 221 CASE ( 3 ) ! Roughness given by the wave model (coupled or read in file) 202 zhsro (:,:) = hsw(:,:)222 zhsro = MAX(rn_frac_hs * hsw, rn_hsro) 203 223 END SELECT 204 224 … … 314 334 ! --------------------------------- 315 335 ! 336 IF( ln_phioc ) THEN 337 SELECT CASE ( nn_bc_surf ) 338 ! 339 CASE ( 0 ) ! Dirichlet case 340 IF( nn_wmix==jp_janssen ) THEN 341 ! First level 342 en(:,:,1) = MAX( rsbc_tke1(:,:) * ustars2(:,:), rn_emin ) 343 z_elem_a(:,:,1) = en(:,:,1) 344 z_elem_c(:,:,1) = 0._wp 345 z_elem_b(:,:,1) = 1._wp 346 347 ! One level below 348 en(:,:,2) = MAX( rsbc_tke1(:,:) * ustars2(:,:) * ((zhsro(:,:)+fsdepw(:,:,2))/zhsro(:,:) )**ra_sf, rn_emin ) 349 z_elem_a(:,:,2) = 0._wp 350 z_elem_c(:,:,2) = 0._wp 351 z_elem_b(:,:,2) = 1._wp 352 ELSE IF( nn_wmix==jp_craigbanner ) THEN 353 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1(:,:))**(2._wp/3._wp) 354 en(:,:,1) = MAX(en(:,:,1), rn_emin) 355 z_elem_a(:,:,1) = en(:,:,1) 356 z_elem_c(:,:,1) = 0._wp 357 z_elem_b(:,:,1) = 1._wp 358 ! 359 ! One level below 360 en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1(:,:) * ((zhsro(:,:)+fsdepw(:,:,2)) & 361 & / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 362 en(:,:,2) = MAX(en(:,:,2), rn_emin ) 363 z_elem_a(:,:,2) = 0._wp 364 z_elem_c(:,:,2) = 0._wp 365 z_elem_b(:,:,2) = 1._wp 366 ! 367 ENDIF 368 CASE ( 1 ) ! Neumann boundary condition on d(e)/dz 369 IF( nn_wmix==jp_janssen ) THEN 370 ! Dirichlet conditions at k=1 371 en(:,:,1) = MAX( rsbc_tke1(:,:) * ustars2(:,:), rn_emin ) 372 z_elem_a(:,:,1) = en(:,:,1) 373 z_elem_c(:,:,1) = 0._wp 374 z_elem_b(:,:,1) = 1._wp 375 ! at k=2, set de/dz=Fw 376 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 377 z_elem_a(:,:,2) = 0._wp 378 zflxs(:,:) = rsbc_tke3(:,:) * ustars2(:,:)**1.5_wp * ((zhsro(:,:)+fsdept(:,:,1) ) / zhsro(:,:) )**(1.5*ra_sf) 379 en(:,:,2) = en(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 380 ELSE IF( nn_wmix==jp_craigbanner ) THEN 381 ! Dirichlet conditions at k=1 382 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1(:,:))**(2._wp/3._wp) 383 en(:,:,1) = MAX(en(:,:,1), rn_emin) 384 z_elem_a(:,:,1) = en(:,:,1) 385 z_elem_c(:,:,1) = 0._wp 386 z_elem_b(:,:,1) = 1._wp 387 ! 388 ! at k=2, set de/dz=Fw 389 !cbr 390 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 391 z_elem_a(:,:,2) = 0._wp 392 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans(:,:)*fsdept(:,:,1)/zhsro(:,:)) )) 393 zflxs(:,:) = rsbc_tke3(:,:) * ustars2(:,:)**1.5_wp * zkar(:,:) & 394 & * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 395 396 en(:,:,2) = en(:,:,2) + zflxs(:,:)/fse3w(:,:,2) 397 ! 398 ENDIF 399 END SELECT 400 ELSE 316 401 SELECT CASE ( nn_bc_surf ) 317 402 ! 318 403 CASE ( 0 ) ! Dirichlet case 319 404 ! First level 320 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 )**(2._wp/3._wp)405 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default)**(2._wp/3._wp) 321 406 en(:,:,1) = MAX(en(:,:,1), rn_emin) 322 407 z_elem_a(:,:,1) = en(:,:,1) … … 325 410 ! 326 411 ! One level below 327 en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2)) &412 en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default * ((zhsro(:,:)+fsdepw(:,:,2)) & 328 413 & / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 329 414 en(:,:,2) = MAX(en(:,:,2), rn_emin ) … … 334 419 ! 335 420 CASE ( 1 ) ! Neumann boundary condition on d(e)/dz 336 !337 421 ! Dirichlet conditions at k=1 338 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 )**(2._wp/3._wp)422 en(:,:,1) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1_default)**(2._wp/3._wp) 339 423 en(:,:,1) = MAX(en(:,:,1), rn_emin) 340 424 z_elem_a(:,:,1) = en(:,:,1) … … 346 430 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 347 431 z_elem_a(:,:,2) = 0._wp 348 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans *fsdept(:,:,1)/zhsro(:,:)) ))432 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans_default*fsdept(:,:,1)/zhsro(:,:)) )) 349 433 zflxs(:,:) = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) & 350 434 & * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) … … 354 438 ! 355 439 END SELECT 440 ENDIF ! ln_phioc 356 441 357 442 ! Bottom boundary condition on tke … … 539 624 ! 540 625 CASE ( 0 ) ! Dirichlet boundary conditions 541 ! 542 ! Surface value 543 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 544 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 545 z_elem_a(:,:,1) = psi(:,:,1) 546 z_elem_c(:,:,1) = 0._wp 547 z_elem_b(:,:,1) = 1._wp 548 ! 549 ! One level below 550 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdepw(:,:,2)/zhsro(:,:) ))) 551 zdep(:,:) = (zhsro(:,:) + fsdepw(:,:,2)) * zkar(:,:) 552 psi (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 553 z_elem_a(:,:,2) = 0._wp 554 z_elem_c(:,:,2) = 0._wp 555 z_elem_b(:,:,2) = 1._wp 556 ! 557 ! 626 IF( ln_phioc ) THEN ! Wave induced mixing case 627 ! en(1) = q2(1) = 0.5 * (15.8 * Ccb)^(2/3) * u*^2 628 ! balance between the production and the 629 ! dissipation terms including the wave effect 630 IF( nn_wmix==jp_janssen ) THEN 631 ! Surface value 632 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 633 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 634 z_elem_a(:,:,1) = psi(:,:,1) 635 z_elem_c(:,:,1) = 0._wp 636 z_elem_b(:,:,1) = 1._wp 637 ! 638 ! One level below 639 zex1 = (rmm*ra_sf+rnn) 640 zex2 = (rmm*ra_sf) 641 zdep(:,:) = ( (zhsro(:,:) + fsdepw(:,:,2))**zex1 ) / zhsro(:,:)**zex2 642 psi (:,:,2) = rsbc_psi1(:,:) * ustars2(:,:)**rmm * zdep(:,:) * tmask(:,:,1) 643 z_elem_a(:,:,2) = 0._wp 644 z_elem_c(:,:,2) = 0._wp 645 z_elem_b(:,:,2) = 1._wp 646 ! 647 ! 648 ELSE IF( nn_wmix==jp_craigbanner ) THEN 649 ! Surface value 650 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 651 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 652 z_elem_a(:,:,1) = psi(:,:,1) 653 z_elem_c(:,:,1) = 0._wp 654 z_elem_b(:,:,1) = 1._wp 655 ! 656 ! One level below 657 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans(:,:)*fsdepw(:,:,2)/zhsro(:,:) ))) 658 zdep(:,:) = (zhsro(:,:) + fsdepw(:,:,2)) * zkar(:,:) 659 psi (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 660 z_elem_a(:,:,2) = 0._wp 661 z_elem_c(:,:,2) = 0._wp 662 z_elem_b(:,:,2) = 1._wp 663 ! 664 ! 665 ENDIF 666 ELSE ! Wave induced mixing case with default values 667 ! Surface value 668 zdep(:,:) = zhsro(:,:) * rl_sf ! Cosmetic 669 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 670 z_elem_a(:,:,1) = psi(:,:,1) 671 z_elem_c(:,:,1) = 0._wp 672 z_elem_b(:,:,1) = 1._wp 673 ! 674 ! One level below 675 zkar(:,:) = (rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans_default*fsdepw(:,:,2)/zhsro(:,:) ))) 676 zdep(:,:) = (zhsro(:,:) + fsdepw(:,:,2)) * zkar(:,:) 677 psi (:,:,2) = rc0**rpp * en(:,:,2)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 678 z_elem_a(:,:,2) = 0._wp 679 z_elem_c(:,:,2) = 0._wp 680 z_elem_b(:,:,2) = 1._wp 681 ! 682 ! 683 ENDIF 558 684 CASE ( 1 ) ! Neumann boundary condition on d(psi)/dz 559 ! 560 ! Surface value: Dirichlet 561 zdep(:,:) = zhsro(:,:) * rl_sf 562 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 563 z_elem_a(:,:,1) = psi(:,:,1) 564 z_elem_c(:,:,1) = 0._wp 565 z_elem_b(:,:,1) = 1._wp 566 ! 567 ! Neumann condition at k=2 568 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 569 z_elem_a(:,:,2) = 0._wp 570 ! 571 ! Set psi vertical flux at the surface: 572 zkar(:,:) = rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 573 zdep(:,:) = ((zhsro(:,:) + fsdept(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 574 zflxs(:,:) = (rnn + rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 575 zdep(:,:) = rsbc_psi1 * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * & 576 & ustars2(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + fsdept(:,:,1))**(rnn-1.) 577 zflxs(:,:) = zdep(:,:) * zflxs(:,:) 578 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 579 580 ! 581 ! 685 IF( ln_phioc ) THEN ! Wave induced mixing case with forced/coupled fields 686 IF( nn_wmix==jp_janssen ) THEN 687 ! 688 zdep(:,:) = rl_sf * zhsro(:,:) 689 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 690 z_elem_a(:,:,1) = psi(:,:,1) 691 z_elem_c(:,:,1) = 0._wp 692 z_elem_b(:,:,1) = 1._wp 693 ! 694 ! Neumann condition at k=2 695 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 696 z_elem_a(:,:,2) = 0._wp 697 ! 698 ! Set psi vertical flux at the surface: 699 zdep(:,:) = (zhsro(:,:) + fsdept(:,:,1))**(rmm*ra_sf+rnn-1._wp) / zhsro(:,:)**(rmm*ra_sf) 700 zflxs(:,:) = rsbc_psi3 * ( zwall_psi(:,:,1)*avm(:,:,1) + zwall_psi(:,:,2)*avm(:,:,2) ) & 701 & * en(:,:,1)**rmm * zdep 702 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 703 ! 704 ELSE IF( nn_wmix==jp_craigbanner ) THEN 705 ! Surface value: Dirichlet 706 zdep(:,:) = zhsro(:,:) * rl_sf 707 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 708 z_elem_a(:,:,1) = psi(:,:,1) 709 z_elem_c(:,:,1) = 0._wp 710 z_elem_b(:,:,1) = 1._wp 711 ! 712 ! Neumann condition at k=2 713 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 714 z_elem_a(:,:,2) = 0._wp 715 ! 716 ! Set psi vertical flux at the surface: 717 zkar(:,:) = rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans(:,:)*fsdept(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 718 zdep(:,:) = ((zhsro(:,:) + fsdept(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 719 zflxs(:,:) = (rnn + rsbc_tke1(:,:) * (rnn + rmm*ra_sf) * zdep(:,:)) * & 720 (1._wp + rsbc_tke1(:,:) * zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 721 zdep(:,:) = rsbc_psi1(:,:) * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * & 722 & ustars2(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + fsdept(:,:,1))**(rnn-1.) 723 zflxs(:,:) = zdep(:,:) * zflxs(:,:) 724 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 725 ! 726 ENDIF 727 ELSE ! Wave induced mixing case with default values 728 ! Surface value: Dirichlet 729 zdep(:,:) = zhsro(:,:) * rl_sf 730 psi (:,:,1) = rc0**rpp * en(:,:,1)**rmm * zdep(:,:)**rnn * tmask(:,:,1) 731 z_elem_a(:,:,1) = psi(:,:,1) 732 z_elem_c(:,:,1) = 0._wp 733 z_elem_b(:,:,1) = 1._wp 734 ! 735 ! Neumann condition at k=2 736 z_elem_b(:,:,2) = z_elem_b(:,:,2) + z_elem_a(:,:,2) ! Remove z_elem_a from z_elem_b 737 z_elem_a(:,:,2) = 0._wp 738 ! 739 ! Set psi vertical flux at the surface: 740 zkar(:,:) = rl_sf + (vkarmn-rl_sf)*(1._wp-exp(-rtrans_default*fsdept(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 741 zdep(:,:) = ((zhsro(:,:) + fsdept(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 742 zflxs(:,:) = (rnn + rsbc_tke1_default * (rnn + rmm*ra_sf) * zdep(:,:)) * & 743 (1._wp + rsbc_tke1_default * zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 744 zdep(:,:) = rsbc_psi1(:,:) * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * & 745 & ustars2(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + fsdept(:,:,1))**(rnn-1.) 746 zflxs(:,:) = zdep(:,:) * zflxs(:,:) 747 psi(:,:,2) = psi(:,:,2) + zflxs(:,:) / fse3w(:,:,2) 748 ! 749 ! 750 ENDIF 582 751 END SELECT 583 752 … … 870 1039 NAMELIST/namzdf_gls/rn_emin, rn_epsmin, ln_length_lim, & 871 1040 & rn_clim_galp, ln_sigpsi, rn_hsro, & 872 & rn_crban , rn_charn, rn_frac_hs,&1041 & rn_crban_default, rn_charn, rn_frac_hs,& 873 1042 & nn_bc_surf, nn_bc_bot, nn_z0_met, & 874 1043 & nn_stab_func, nn_clos … … 898 1067 WRITE(numout,*) ' TKE Bottom boundary condition nn_bc_bot = ', nn_bc_bot 899 1068 WRITE(numout,*) ' Modify psi Schmidt number (wb case) ln_sigpsi = ', ln_sigpsi 900 WRITE(numout,*) ' Craig and Banner coefficient rn_crban = ', rn_crban1069 WRITE(numout,*) ' Craig and Banner coefficient (default) rn_crban = ', rn_crban_default 901 1070 WRITE(numout,*) ' Charnock coefficient rn_charn = ', rn_charn 902 1071 WRITE(numout,*) ' Surface roughness formula nn_z0_met = ', nn_z0_met … … 915 1084 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' ) 916 1085 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' ) 917 IF( nn_z0_met == 3 .AND. .NOT.ln_wave ) CALL ctl_stop( 'zdf_gls_init: nn_z0_met=3 requires ln_wave=T' ) 1086 IF( nn_z0_met == 3 .AND. .NOT.ln_rough ) CALL ctl_stop( 'zdf_gls_init: nn_z0_met=3 requires ln_rough=T' ) 1087 IF( nn_z0_met .NE. 3 .AND. ln_rough ) THEN 1088 CALL ctl_warn('W A R N I N G: ln_rough=.TRUE. but nn_z0_met is not 3 - resetting nn_z0_met to 3') 1089 nn_z0_met = 3 1090 ENDIF 918 1091 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' ) 919 1092 IF( nn_clos < 0 .OR. nn_clos > 3 ) CALL ctl_stop( 'zdf_gls_init: bad flag: nn_clos is 0, 1, 2 or 3' ) … … 1089 1262 & - SQRT(rsc_tke + 24._wp*rsc_psi0*rpsi2 ) ) 1090 1263 1091 IF ( rn_crban==0._wp ) THEN1264 IF( .NOT. ln_phioc .AND. rn_crban_default==0._wp ) THEN 1092 1265 rl_sf = vkarmn 1093 1266 ELSE … … 1128 1301 rc03 = rc02 * rc0 1129 1302 rc04 = rc03 * rc0 1130 rsbc_tke1 = -3._wp/2._wp*rn_crban*ra_sf*rl_sf ! Dirichlet + Wave breaking1131 rsbc_tke2 = rdt * rn_crban / rl_sf ! Neumann + Wave breaking1132 zcr = MAX(rsmall, rsbc_tke1**(1./(-ra_sf*3._wp/2._wp))-1._wp )1133 rtrans = 0.2_wp / zcr ! Ad. inverse transition length between log and wave layer1303 rsbc_tke1_default = -3._wp/2._wp*rn_crban_default*ra_sf*rl_sf 1304 rsbc_tke2 = rdt * rn_crban_default / rl_sf 1305 zcr = MAX( rsmall, rsbc_tke1_default**(1./(-ra_sf*3._wp/2._wp))-1._wp ) 1306 rtrans_default = 0.2_wp / zcr 1134 1307 rsbc_zs1 = rn_charn/grav ! Charnock formula for surface roughness 1135 1308 rsbc_zs2 = rn_frac_hs / 0.85_wp / grav * 665._wp ! Rascle formula for surface roughness -
branches/UKMO/AMM15_v3_6_STABLE_package_collate_coupling/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90
r10473 r10478 54 54 !! 55 55 NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp, & 56 & ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp, & 57 & ln_zdfqiao 56 & ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp 58 57 !!---------------------------------------------------------------------- 59 58 … … 84 83 WRITE(numout,*) ' npc call frequency nn_npc = ', nn_npc 85 84 WRITE(numout,*) ' npc print frequency nn_npcp = ', nn_npcp 86 WRITE(numout,*) ' Qiao formulation flag ln_zdfqiao=', ln_zdfqiao87 85 ENDIF 88 86
Note: See TracChangeset
for help on using the changeset viewer.