- Timestamp:
- 2021-06-11T10:22:51+02:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14608_AGRIF_domcfg/src/OCE/ZDF/zdfgls.F90
r14958 r14971 26 26 USE zdfmxl ! mixed layer 27 27 USE sbcwave , ONLY : hsw ! significant wave height 28 #if defined key_si3 29 USE ice, ONLY: hm_i, h_i 30 #endif 31 #if defined key_cice 32 USE sbc_ice, ONLY: h_i 33 #endif 28 34 ! 29 35 USE lbclnk ! ocean lateral boundary conditions (or mpp link) … … 51 57 LOGICAL :: ln_length_lim ! use limit on the dissipation rate under stable stratification (Galperin et al. 1988) 52 58 LOGICAL :: ln_sigpsi ! Activate Burchard (2003) modification for k-eps closure & wave breaking mixing 59 INTEGER :: nn_mxlice ! type of scaling under sea-ice (=0/1/2/3) 53 60 INTEGER :: nn_bc_surf ! surface boundary condition (=0/1) 54 61 INTEGER :: nn_bc_bot ! bottom boundary condition (=0/1) … … 224 231 ! 225 232 ! adapt roughness where there is sea ice 226 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 227 zhsro(ji,jj) = ( (1._wp-zice_fra(ji,jj)) * zhsro(ji,jj) + zice_fra(ji,jj) * rn_hsri )*tmask(ji,jj,1) + & 228 & (1._wp - tmask(ji,jj,1))*rn_hsro 229 END_2D 233 SELECT CASE( nn_mxlice ) ! Type of scaling under sea-ice 234 ! 235 CASE( 1 ) ! scaling with constant sea-ice roughness 236 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 237 zhsro(ji,jj) = ( (1._wp-zice_fra(ji,jj)) * zhsro(ji,jj) + zice_fra(ji,jj) * rn_hsri )*tmask(ji,jj,1) + (1._wp - tmask(ji,jj,1))*rn_hsro 238 END_2D 239 ! 240 CASE( 2 ) ! scaling with mean sea-ice thickness 241 #if defined key_si3 242 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 243 zhsro(ji,jj) = ( (1._wp-zice_fra(ji,jj)) * zhsro(ji,jj) + zice_fra(ji,jj) * hm_i(ji,jj) )*tmask(ji,jj,1) + (1._wp - tmask(ji,jj,1))*rn_hsro 244 END_2D 245 #endif 246 ! 247 CASE( 3 ) ! scaling with max sea-ice thickness 248 #if defined key_si3 || defined key_cice 249 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 250 zhsro(ji,jj) = ( (1._wp-zice_fra(ji,jj)) * zhsro(ji,jj) + zice_fra(ji,jj) * MAXVAL(h_i(ji,jj,:)) )*tmask(ji,jj,1) + (1._wp - tmask(ji,jj,1))*rn_hsro 251 END_2D 252 #endif 253 ! 254 END SELECT 230 255 ! 231 256 DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) !== Compute dissipation rate ==! … … 329 354 END_2D 330 355 ! 356 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 357 DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 358 IF( mikt(ji,jj) > 1 )THEN 359 itop = mikt(ji,jj) ! k top w-point 360 itopp1 = mikt(ji,jj) + 1 ! k+1 1st w-point below the top one 361 ! ! mask at the 362 ! ocean surface 363 ! points 364 z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 365 ! 366 ! Dirichlet condition applied at: 367 ! top level (itop) & Just below it (itopp1) 368 zd_lw(ji,jj,itop) = 0._wp ; zd_lw(ji,jj,itopp1) = 0._wp 369 zd_up(ji,jj,itop) = 0._wp ; zd_up(ji,jj,itopp1) = 0._wp 370 zdiag(ji,jj,itop) = 1._wp ; zdiag(ji,jj,itopp1) = 1._wp 371 en (ji,jj,itop) = z_en ; en (ji,jj,itopp1) = z_en 372 ENDIF 373 END_2D 374 ENDIF 331 375 ! 332 376 CASE ( 1 ) ! Neumann boundary condition (set d(e)/dz) … … 352 396 END_2D 353 397 ! 398 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 399 DO_2D_OVR( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 400 IF( mikt(ji,jj) > 1 )THEN 401 itop = mikt(ji,jj) ! k top w-point 402 itopp1 = mikt(ji,jj) + 1 ! k+1 1st w-point below the top one 403 ! ! mask at the 404 ! ocean surface 405 ! points 406 z_en = MAX( rc02r * ustar2_top(ji,jj), rn_emin ) * ( 1._wp - tmask(ji,jj,1) ) 407 ! 408 ! Bottom level Dirichlet condition: 409 ! Bottom level (ibot) & Just above it (ibotm1) 410 ! Dirichlet ! Neumann 411 zd_lw(ji,jj,itop) = 0._wp ! ! Remove zd_up from zdiag 412 zdiag(ji,jj,itop) = 1._wp ; zdiag(ji,jj,itopp1) = zdiag(ji,jj,itopp1) + zd_up(ji,jj,itopp1) 413 zd_up(ji,jj,itop) = 0._wp ; zd_up(ji,jj,itopp1) = 0._wp 414 en (ji,jj,itop) = z_en 415 ENDIF 416 END_2D 417 ENDIF 354 418 ! 355 419 END SELECT … … 607 671 END_2D 608 672 ! 673 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 674 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 675 IF ( mikt(ji,jj) > 1 ) THEN 676 itop = mikt(ji,jj) ! k top w-point 677 itopp1 = mikt(ji,jj) + 1 ! k+1 1st w-point below the top one 678 ! 679 zdep(ji,jj) = vkarmn * r_z0_top 680 psi (ji,jj,itop) = rc0**rpp * en(ji,jj,itop)**rmm *zdep(ji,jj)**rnn 681 zd_lw(ji,jj,itop) = 0._wp 682 zd_up(ji,jj,itop) = 0._wp 683 zdiag(ji,jj,itop) = 1._wp 684 ! 685 ! Just above last level, Dirichlet condition again (GOTM like) 686 zdep(ji,jj) = vkarmn * ( r_z0_top + e3t(ji,jj,itopp1,Kmm) ) 687 psi (ji,jj,itopp1) = rc0**rpp * en(ji,jj,itop )**rmm *zdep(ji,jj)**rnn 688 zd_lw(ji,jj,itopp1) = 0._wp 689 zd_up(ji,jj,itopp1) = 0._wp 690 zdiag(ji,jj,itopp1) = 1._wp 691 END IF 692 END_2D 693 END IF 694 ! 609 695 CASE ( 1 ) ! Neumman boundary condition 610 696 ! … … 631 717 psi(ji,jj,ibotm1) = psi(ji,jj,ibotm1) + zflxb / e3w(ji,jj,ibotm1,Kmm) 632 718 END_2D 719 ! 720 IF( ln_isfcav) THEN ! top boundary (ocean cavity) 721 DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 722 IF ( mikt(ji,jj) > 1 ) THEN 723 itop = mikt(ji,jj) ! k top w-point 724 itopp1 = mikt(ji,jj) + 1 ! k+1 1st w-point below the top one 725 ! 726 ! Bottom level Dirichlet condition: 727 zdep(ji,jj) = vkarmn * r_z0_top 728 psi (ji,jj,itop) = rc0**rpp * en(ji,jj,itop)**rmm *zdep(ji,jj)**rnn 729 ! 730 zd_lw(ji,jj,itop) = 0._wp 731 zd_up(ji,jj,itop) = 0._wp 732 zdiag(ji,jj,itop) = 1._wp 733 ! 734 ! Just below cavity level: Neumann condition with flux 735 ! injection 736 zdiag(ji,jj,itopp1) = zdiag(ji,jj,itopp1) + zd_up(ji,jj,itopp1) ! Remove zd_up from zdiag 737 zd_up(ji,jj,itopp1) = 0._wp 738 ! 739 ! Set psi vertical flux below cavity: 740 zdep(ji,jj) = r_z0_top + 0.5_wp*e3t(ji,jj,itopp1,Kmm) 741 zflxb = rsbc_psi2 * ( p_avm(ji,jj,itop) + p_avm(ji,jj,itopp1)) & 742 & * (0.5_wp*(en(ji,jj,itop)+en(ji,jj,itopp1)))**rmm * zdep(ji,jj)**(rnn-1._wp) 743 psi(ji,jj,itopp1) = psi(ji,jj,itopp1) + zflxb / e3w(ji,jj,itopp1,Kmm) 744 END IF 745 END_2D 746 END IF 747 633 748 ! 634 749 END SELECT … … 799 914 NAMELIST/namzdf_gls/rn_emin, rn_epsmin, ln_length_lim, & 800 915 & rn_clim_galp, ln_sigpsi, rn_hsro, rn_hsri, & 801 & rn_crban, rn_charn, rn_frac_hs,&916 & nn_mxlice, rn_crban, rn_charn, rn_frac_hs, & 802 917 & nn_bc_surf, nn_bc_bot, nn_z0_met, nn_z0_ice, & 803 918 & nn_stab_func, nn_clos … … 839 954 WRITE(numout,*) ' Type of closure nn_clos = ', nn_clos 840 955 WRITE(numout,*) ' Surface roughness (m) rn_hsro = ', rn_hsro 841 WRITE(numout,*) ' Ice-ocean roughness (used if nn_z0_ice/=0) rn_hsri = ', rn_hsri 956 WRITE(numout,*) ' type of scaling under sea-ice nn_mxlice = ', nn_mxlice 957 IF( nn_mxlice == 1 ) & 958 WRITE(numout,*) ' Ice-ocean roughness (used if nn_z0_ice/=0) rn_hsri = ', rn_hsri 959 SELECT CASE( nn_mxlice ) ! Type of scaling under sea-ice 960 CASE( 0 ) ; WRITE(numout,*) ' ==>>> No scaling under sea-ice' 961 CASE( 1 ) ; WRITE(numout,*) ' ==>>> scaling with constant sea-ice thickness' 962 CASE( 2 ) ; WRITE(numout,*) ' ==>>> scaling with mean sea-ice thickness' 963 CASE( 3 ) ; WRITE(numout,*) ' ==>>> scaling with max sea-ice thickness' 964 CASE DEFAULT 965 CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 ') 966 END SELECT 842 967 WRITE(numout,*) 843 968 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
