MODULE zdfphy !!====================================================================== !! *** MODULE zdfphy *** !! Ocean physics : manager of vertical mixing parametrizations !!====================================================================== !! History : 4.0 ! 2017-04 (G. Madec) original code !!---------------------------------------------------------------------- !!---------------------------------------------------------------------- !! zdf_phy_init : initialization of all vertical physics pakages !! zdf_phy : upadate at each time-step the vertical mixing coeff. !!---------------------------------------------------------------------- USE par_oce ! mesh and scale factors USE zdf_oce ! TKE vertical mixing USE sbc_oce ! surface module (only for nn_isf in the option compatibility test) USE zdfbfr ! bottom friction USE zdftke ! TKE vertical mixing USE zdfgls ! GLS vertical mixing USE zdfric ! Richardson vertical mixing USE zdfddm ! double diffusion mixing USE zdfevd ! enhanced vertical diffusion USE zdftmx ! internal tide-induced mixing USE zdfqiao !Qiao module wave induced mixing (zdf_qiao routine) USE zdfmxl ! Mixed-layer depth (zdf_mxl routine) USE tranpc ! convection: non penetrative adjustment USE sbcrnf ! surface boundary condition: runoff variables ! USE in_out_manager ! I/O manager USE iom ! IOM library USE lib_mpp ! distribued memory computing IMPLICIT NONE PRIVATE PUBLIC zdf_phy_init ! routine called by nemogcm.F90 PUBLIC zdf_phy ! routine called by step.F90 !!---------------------------------------------------------------------- !! NEMO/OPA 4.0 , NEMO Consortium (2011) !! $Id$ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE zdf_phy_init !!---------------------------------------------------------------------- !! *** ROUTINE zdf_phy_init *** !! !! ** Purpose : initializations of the vertical ocean physics !! !! ** Method : Read namelist namzdf, control logicals !!---------------------------------------------------------------------- INTEGER :: ioptio, ios ! local integers !! NAMELIST/namzdf/ ln_zdfcst, ln_zdfric, ln_zdftke, ln_zdfgls, & ! type of closure scheme & ln_zdfevd, nn_evdm, rn_evd , & ! convection : evd & ln_zdfnpc, nn_npc , nn_npcp, & ! convection : npc & ln_zdfddm, rn_avts, rn_hsbfr, & ! double diffusion & ln_zdftmx, & ! tidal mixing & ln_zdfqiao, & ! surface wave-induced mixing & ln_zdfexp, nn_zdfexp, & ! time-stepping & rn_avm0, rn_avt0, nn_avb, nn_havtb ! coefficients !!org NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp, & !!org & ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp, & !!org & ln_zdfqiao !!---------------------------------------------------------------------- REWIND( numnam_ref ) ! Namelist namzdf in reference namelist : Vertical mixing parameters READ ( numnam_ref, namzdf, IOSTAT = ios, ERR = 901) 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf in reference namelist', lwp ) REWIND( numnam_cfg ) ! Namelist namzdf in reference namelist : Vertical mixing parameters READ ( numnam_cfg, namzdf, IOSTAT = ios, ERR = 902 ) 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf in configuration namelist', lwp ) IF(lwm) WRITE ( numond, namzdf ) IF(lwp) THEN !* Parameter print WRITE(numout,*) WRITE(numout,*) 'zdf_phy_init : vertical physics' WRITE(numout,*) '~~~~~~~~' WRITE(numout,*) ' Namelist namzdf : set vertical mixing mixing parameters' WRITE(numout,*) ' vertical closure scheme' WRITE(numout,*) ' constant vertical mixing coefficient ln_zdfcst = ', ln_zdfcst WRITE(numout,*) ' constant vertical mixing coefficient ln_zdfric = ', ln_zdfric WRITE(numout,*) ' constant vertical mixing coefficient ln_zdftke = ', ln_zdftke WRITE(numout,*) ' constant vertical mixing coefficient ln_zdfgls = ', ln_zdfgls WRITE(numout,*) ' convection: ' WRITE(numout,*) ' enhanced vertical diffusion ln_zdfevd = ', ln_zdfevd WRITE(numout,*) ' applied on momentum (=1/0) nn_evdm = ', nn_evdm WRITE(numout,*) ' vertical coefficient for evd rn_evd = ', rn_evd WRITE(numout,*) ' non-penetrative convection (npc) ln_zdfnpc = ', ln_zdfnpc WRITE(numout,*) ' npc call frequency nn_npc = ', nn_npc WRITE(numout,*) ' npc print frequency nn_npcp = ', nn_npcp WRITE(numout,*) ' double diffusive mixing ln_zdfddm = ', ln_zdfddm WRITE(numout,*) ' maximum avs for dd mixing rn_avts = ', rn_avts WRITE(numout,*) ' heat/salt buoyancy flux ratio rn_hsbfr= ', rn_hsbfr WRITE(numout,*) ' surface wave-induced mixing ln_zdfqiao= ', ln_zdfqiao ! surface wave induced mixing WRITE(numout,*) ' tidal mixing ln_zdftmx = ', ln_zdftmx WRITE(numout,*) ' time splitting / backward scheme ln_zdfexp = ', ln_zdfexp WRITE(numout,*) ' number of sub-time step (ln_zdfexp=T) nn_zdfexp = ', nn_zdfexp WRITE(numout,*) ' coefficients : ' WRITE(numout,*) ' vertical eddy viscosity rn_avm0 = ', rn_avm0 WRITE(numout,*) ' vertical eddy diffusivity rn_avt0 = ', rn_avt0 WRITE(numout,*) ' constant background or profile nn_avb = ', nn_avb WRITE(numout,*) ' horizontal variation for avtb nn_havtb = ', nn_havtb ENDIF !!gm IF(ln_zdfddm) THEN ! double diffusive mixing' ! avs(:,:,:) = rn_avt0 * wmask(:,:,:) !!gm ENDIF ! !* Parameter & logical controls ! ! ---------------------------- ! ! ! ... check of vertical mixing scheme on tracers ! ==> will be done in trazdf module ! ! ! ... check of mixing coefficient IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' vertical mixing option :' ioptio = 0 IF( ln_zdfcst ) THEN IF(lwp) WRITE(numout,*) ' constant eddy diffusion coefficients' ioptio = ioptio+1 ENDIF IF( ln_zdfric ) THEN IF(lwp) WRITE(numout,*) ' Richardson dependent eddy coefficients' ioptio = ioptio+1 ENDIF IF( ln_zdftke ) THEN IF(lwp) WRITE(numout,*) ' TKE dependent eddy coefficients' ioptio = ioptio+1 ENDIF IF( ln_zdfgls ) THEN IF(lwp) WRITE(numout,*) ' GLS dependent eddy coefficients' ioptio = ioptio+1 ENDIF IF( ioptio == 0 .OR. ioptio > 1 ) & & CALL ctl_stop( ' one and only one vertical diffusion option has to be defined ' ) IF( ( ln_zdfric .OR. ln_zdfgls ) .AND. ln_isfcav ) & & CALL ctl_stop( ' only zdfcst and zdftke were tested with ice shelves cavities ' ) ! ! ! ... Convection IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' convection :' ! #if defined key_top IF( ln_zdfnpc ) CALL ctl_stop( ' zdf_phy_init: npc scheme is not working with key_top' ) #endif ! ioptio = 0 IF( ln_zdfnpc ) THEN IF(lwp) WRITE(numout,*) ' use non penetrative convective scheme' ioptio = ioptio+1 ENDIF IF( ln_zdfevd ) THEN IF(lwp) WRITE(numout,*) ' use enhanced vertical dif. scheme' ioptio = ioptio+1 ENDIF IF( ln_zdftke ) THEN IF(lwp) WRITE(numout,*) ' use the 1.5 turbulent closure' ENDIF IF( ln_zdfgls ) THEN IF(lwp) WRITE(numout,*) ' use the GLS closure scheme' ENDIF IF ( ioptio > 1 ) CALL ctl_stop( ' chose between ln_zdfnpc and ln_zdfevd' ) IF( ioptio == 0 .AND. .NOT.( ln_zdftke .OR. ln_zdfgls ) ) & CALL ctl_stop( ' except for TKE or GLS physics, a convection scheme is', & & ' required: ln_zdfevd or ln_zdfnpc logicals' ) ! !* Background eddy viscosity and diffusivity profil IF( nn_avb == 0 ) THEN ! Define avmb, avtb from namelist parameter avmb(:) = rn_avm0 avtb(:) = rn_avt0 ELSE ! Background profile of avt (fit a theoretical/observational profile (Krauss 1990) avmb(:) = rn_avm0 avtb(:) = rn_avt0 + ( 3.e-4_wp - 2._wp * rn_avt0 ) * 1.e-4_wp * gdepw_1d(:) ! m2/s IF(ln_sco .AND. lwp) CALL ctl_warn( 'avtb profile not valid in sco' ) ENDIF ! IF( ln_rstart ) THEN ! Read avmb, avtb in restart (if exist) ! if ln_traadv_cen, avmb, avtb have been modified in traadv_cen2 module. ! To ensure the restartability, avmb & avtb are written in the restart ! file in traadv_cen2 end read here. IF( iom_varid( numror, 'avmb', ldstop = .FALSE. ) > 0 ) THEN CALL iom_get( numror, jpdom_unknown, 'avmb', avmb ) CALL iom_get( numror, jpdom_unknown, 'avtb', avtb ) ENDIF ENDIF ! ! 2D shape of the avtb avtb_2d(:,:) = 1.e0 ! uniform ! IF( nn_havtb == 1 ) THEN ! decrease avtb in the equatorial band ! -15S -5S : linear decrease from avt0 to avt0/10. ! -5S +5N : cst value avt0/10. ! 5N 15N : linear increase from avt0/10, to avt0 WHERE(-15. <= gphit .AND. gphit < -5 ) avtb_2d = (1. - 0.09 * (gphit + 15.)) WHERE( -5. <= gphit .AND. gphit < 5 ) avtb_2d = 0.1 WHERE( 5. <= gphit .AND. gphit < 15 ) avtb_2d = (0.1 + 0.09 * (gphit - 5.)) ENDIF ! !!gm moved into zdf_phy_init ! CALL zdf_bfr_init ! bottom friction ioptio = 0 !== type of vertical turbulent closure ==! (set nzdfphy) ! ! IF( ln_zdfcst ) THEN ; ioptio = ioptio + 1 ; nzdf_phy = np_CST ; ENDIF ! IF( ln_zdfric ) THEN ; ioptio = ioptio + 1 ; nzdf_phy = np_RIC ; CALL zdf_ric_init ; ENDIF ! IF( ln_zdftke ) THEN ; ioptio = ioptio + 1 ; nzdf_phy = np_TKE ; CALL zdf_tke_init ; ENDIF ! IF( ln_zdfgls ) THEN ; ioptio = ioptio + 1 ; nzdf_phy = np_GLS ; CALL zdf_gls_init ; ENDIF ! IF( ln_zdfric ) CALL zdf_ric_init ! Richardson number dependent Kz IF( ln_zdftke ) CALL zdf_tke_init ! TKE closure scheme IF( ln_zdfgls ) CALL zdf_gls_init ! GLS closure scheme IF( ln_zdftmx ) CALL zdf_tmx_init ! tidal vertical mixing !!gm ! END SUBROUTINE zdf_phy_init SUBROUTINE zdf_phy( kstp ) !!---------------------------------------------------------------------- !! *** ROUTINE zdf_phy *** !! !! ** Purpose : Update ocean physics at each time-step !! !! ** Method : !! !! ** Action : avm, avt vertical eddy viscosity and diffusivity at w-points !! nmld ??? mixed layer depth in level and meters <<<<====verifier ! !! bottom stress..... <<<<====verifier ! !!---------------------------------------------------------------------- INTEGER, INTENT(in) :: kstp ! ocean time-step index ! INTEGER :: ji, jj, jk ! dummy loop indice !!---------------------------------------------------------------------- ! CALL zdf_bfr( kstp ) ! bottom friction (if quadratic) ! ! Vertical eddy viscosity and diffusivity coefficients IF( ln_zdfric ) CALL zdf_ric ( kstp ) ! Richardson number dependent Kz IF( ln_zdftke ) CALL zdf_tke ( kstp ) ! TKE closure scheme for Kz IF( ln_zdfgls ) CALL zdf_gls ( kstp ) ! GLS closure scheme for Kz IF( ln_zdfqiao ) CALL zdf_qiao( kstp ) ! Qiao vertical mixing ! IF( ln_zdfcst ) THEN ! Constant Kz (reset avt, avm[uv] to the background value) avt (:,:,:) = rn_avt0 * wmask (:,:,:) avm (:,:,:) = rn_avm0 * wmask (:,:,:) avmu(:,:,:) = rn_avm0 * wumask(:,:,:) avmv(:,:,:) = rn_avm0 * wvmask(:,:,:) ENDIF ! IF( ln_rnf_mouth ) THEN ! increase diffusivity at rivers mouths DO jk = 2, nkrnf ; avt(:,:,jk) = avt(:,:,jk) + 2._wp * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk) ; END DO ENDIF ! IF( ln_zdfevd ) CALL zdf_evd( kstp ) ! enhanced vertical eddy diffusivity ! IF( ln_zdfddm ) THEN ! double diffusive mixing CALL zdf_ddm( kstp ) ELSE ! avs=avt DO jk = 2, jpkm1 ; avs(:,:,jk) = avt(:,:,jk) ; END DO ENDIF ! IF( ln_zdftmx ) CALL zdf_tmx( kstp ) ! tidal vertical mixing CALL zdf_mxl( kstp ) ! mixed layer depth ! write TKE or GLS information in the restart file IF( lrst_oce .AND. ln_zdftke ) CALL tke_rst( kstp, 'WRITE' ) IF( lrst_oce .AND. ln_zdfgls ) CALL gls_rst( kstp, 'WRITE' ) ! END SUBROUTINE zdf_phy !!====================================================================== END MODULE zdfphy