New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13284 for NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF – NEMO

Ignore:
Timestamp:
2020-07-09T17:12:23+02:00 (4 years ago)
Author:
smasson
Message:

4.0-HEAD: merge 4.0-HEAD_r12713_clem_dan_fixcpl into 4.0-HEAD

Location:
NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF/zdfdrg.F90

    r13268 r13284  
    3232   USE lib_mpp        ! distributed memory computing 
    3333   USE prtctl         ! Print control 
     34   USE sbc_oce , ONLY : nn_ice  
    3435 
    3536   IMPLICIT NONE 
     
    4647   LOGICAL          ::   ln_loglayer  ! logarithmic drag: Cd = vkarmn/log(z/z0) 
    4748   LOGICAL , PUBLIC ::   ln_drgimp    ! implicit top/bottom friction flag 
    48  
     49   LOGICAL , PUBLIC ::   ln_drgice_imp ! implicit ice-ocean drag  
    4950   !                                 !!* Namelist namdrg_top & _bot: TOP or BOTTOM coefficient namelist * 
    5051   REAL(wp)         ::   rn_Cd0       !: drag coefficient                                           [ - ] 
     
    231232      INTEGER   ::   ios, ioptio   ! local integers 
    232233      !! 
    233       NAMELIST/namdrg/ ln_drg_OFF, ln_lin, ln_non_lin, ln_loglayer, ln_drgimp 
     234      NAMELIST/namdrg/ ln_drg_OFF, ln_lin, ln_non_lin, ln_loglayer, ln_drgimp, ln_drgice_imp 
    234235      !!---------------------------------------------------------------------- 
    235236      ! 
     
    254255         WRITE(numout,*) '      logarithmic drag: Cd = vkarmn/log(z/z0)   ln_loglayer = ', ln_loglayer 
    255256         WRITE(numout,*) '      implicit friction                         ln_drgimp   = ', ln_drgimp 
     257         WRITE(numout,*) '      implicit ice-ocean drag                   ln_drgice_imp  =', ln_drgice_imp 
    256258      ENDIF 
    257259      ! 
     
    264266      IF( ioptio /= 1 )   CALL ctl_stop( 'zdf_drg_init: Choose ONE type of drag coef in namdrg' ) 
    265267      ! 
     268      IF ( ln_drgice_imp.AND.(.NOT.ln_drgimp) ) &  
     269         &                CALL ctl_stop( 'zdf_drg_init: ln_drgice_imp=T requires ln_drgimp=T' ) 
     270      ! 
     271      IF ( ln_drgice_imp.AND.( nn_ice /=2 ) ) & 
     272         &  CALL ctl_stop( 'zdf_drg_init: ln_drgice_imp=T requires si3' ) 
    266273      ! 
    267274      !                     !==  BOTTOM drag setting  ==!   (applied at seafloor) 
     
    274281      !                     !==  TOP drag setting  ==!   (applied at the top of ocean cavities) 
    275282      ! 
    276       IF( ln_isfcav ) THEN              ! Ocean cavities: top friction setting 
    277          ALLOCATE( rCd0_top(jpi,jpj), rCdU_top(jpi,jpj) ) 
     283      IF( ln_isfcav.OR.ln_drgice_imp ) THEN              ! Ocean cavities: top friction setting 
     284         ALLOCATE( rCdU_top(jpi,jpj) ) 
     285      ENDIF 
     286      ! 
     287      IF( ln_isfcav ) THEN 
     288         ALLOCATE( rCd0_top(jpi,jpj)) 
    278289         CALL drg_init( 'TOP   '   , mikt       ,                                         &   ! <== in 
    279290            &           r_Cdmin_top, r_Cdmax_top, r_z0_top, r_ke0_top, rCd0_top, rCdU_top )   ! ==> out 
  • NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF/zdfgls.F90

    r13268 r13284  
    5454   INTEGER  ::   nn_bc_bot         ! bottom boundary condition (=0/1) 
    5555   INTEGER  ::   nn_z0_met         ! Method for surface roughness computation 
     56   INTEGER  ::   nn_z0_ice         ! Roughness accounting for sea ice 
    5657   INTEGER  ::   nn_stab_func      ! stability functions G88, KC or Canuto (=0/1/2) 
    5758   INTEGER  ::   nn_clos           ! closure 0/1/2/3 MY82/k-eps/k-w/gen 
     
    6263   REAL(wp) ::   rn_crban          ! Craig and Banner constant for surface breaking waves mixing 
    6364   REAL(wp) ::   rn_hsro           ! Minimum surface roughness 
     65   REAL(wp) ::   rn_hsri           ! Ice ocean roughness 
    6466   REAL(wp) ::   rn_frac_hs        ! Fraction of wave height as surface roughness (if nn_z0_met > 1)  
    6567 
     
    151153      REAL(wp), DIMENSION(jpi,jpj)     ::   zflxs       ! Turbulence fluxed induced by internal waves  
    152154      REAL(wp), DIMENSION(jpi,jpj)     ::   zhsro       ! Surface roughness (surface waves) 
     155      REAL(wp), DIMENSION(jpi,jpj)     ::   zice_fra    ! Tapering of wave breaking under sea ice 
    153156      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   eb          ! tke at time before 
    154157      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   hmxl_b      ! mixing length at time before 
     
    166169      ustar2_bot (:,:) = 0._wp 
    167170 
     171      SELECT CASE ( nn_z0_ice ) 
     172      CASE( 0 )   ;   zice_fra(:,:) = 0._wp 
     173      CASE( 1 )   ;   zice_fra(:,:) =        TANH( fr_i(:,:) * 10._wp ) 
     174      CASE( 2 )   ;   zice_fra(:,:) =              fr_i(:,:) 
     175      CASE( 3 )   ;   zice_fra(:,:) = MIN( 4._wp * fr_i(:,:) , 1._wp ) 
     176      END SELECT 
     177       
    168178      ! Compute surface, top and bottom friction at T-points 
    169179      DO jj = 2, jpjm1              !==  surface ocean friction 
     
    211221      END SELECT 
    212222      ! 
     223      ! adapt roughness where there is sea ice 
     224      zhsro(:,:) = ( (1._wp-zice_fra(:,:)) * zhsro(:,:) + zice_fra(:,:) * rn_hsri )*tmask(:,:,1)  + (1._wp - tmask(:,:,1))*rn_hsro 
     225      ! 
    213226      DO jk = 2, jpkm1              !==  Compute dissipation rate  ==! 
    214227         DO jj = 1, jpjm1 
     
    305318      CASE ( 0 )             ! Dirichlet boundary condition (set e at k=1 & 2)  
    306319      ! First level 
    307       en   (:,:,1) = MAX(  rn_emin , rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3  ) 
     320      en   (:,:,1) = MAX(  rn_emin , rc02r * ustar2_surf(:,:) * (1._wp + (1._wp-zice_fra(:,:))*rsbc_tke1)**r2_3  ) 
    308321      zd_lw(:,:,1) = en(:,:,1) 
    309322      zd_up(:,:,1) = 0._wp 
     
    311324      !  
    312325      ! One level below 
    313       en   (:,:,2) =  MAX(  rc02r * ustar2_surf(:,:) * (  1._wp + rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2))  & 
    314          &                 / zhsro(:,:) )**(1.5_wp*ra_sf)  )**(2._wp/3._wp)                      , rn_emin   ) 
     326      en   (:,:,2) =  MAX(  rc02r * ustar2_surf(:,:) * (  1._wp + (1._wp-zice_fra(:,:))*rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2)) & 
     327         &                 / zhsro(:,:) )**(1.5_wp*ra_sf)  )**(2._wp/3._wp) , rn_emin   ) 
    315328      zd_lw(:,:,2) = 0._wp  
    316329      zd_up(:,:,2) = 0._wp 
     
    321334      ! 
    322335      ! Dirichlet conditions at k=1 
    323       en   (:,:,1) = MAX(  rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3 , rn_emin  ) 
     336      en   (:,:,1) = MAX(  rc02r * ustar2_surf(:,:) * (1._wp + (1._wp-zice_fra(:,:))*rsbc_tke1)**r2_3 , rn_emin  ) 
    324337      zd_lw(:,:,1) = en(:,:,1) 
    325338      zd_up(:,:,1) = 0._wp 
     
    331344      zd_lw(:,:,2) = 0._wp 
    332345      zkar (:,:)   = (rl_sf + (vkarmn-rl_sf)*(1.-EXP(-rtrans*gdept_n(:,:,1)/zhsro(:,:)) )) 
    333       zflxs(:,:)   = rsbc_tke2 * ustar2_surf(:,:)**1.5_wp * zkar(:,:) & 
     346      zflxs(:,:)   = rsbc_tke2 * (1._wp-zice_fra(:,:)) * ustar2_surf(:,:)**1.5_wp * zkar(:,:) & 
    334347          &                    * (  ( zhsro(:,:)+gdept_n(:,:,1) ) / zhsro(:,:)  )**(1.5_wp*ra_sf) 
    335348!!gm why not   :                        * ( 1._wp + gdept_n(:,:,1) / zhsro(:,:) )**(1.5_wp*ra_sf) 
     
    582595         zkar (:,:)   = rl_sf + (vkarmn-rl_sf)*(1._wp-EXP(-rtrans*gdept_n(:,:,1)/zhsro(:,:) )) ! Lengh scale slope 
    583596         zdep (:,:)   = ((zhsro(:,:) + gdept_n(:,:,1)) / zhsro(:,:))**(rmm*ra_sf) 
    584          zflxs(:,:)   = (rnn + rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 
     597         zflxs(:,:)   = (rnn + (1._wp-zice_fra(:,:))*rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:)) & 
     598            &           *(1._wp + (1._wp-zice_fra(:,:))*rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp) 
    585599         zdep (:,:)   = rsbc_psi1 * (zwall_psi(:,:,1)*p_avm(:,:,1)+zwall_psi(:,:,2)*p_avm(:,:,2)) * & 
    586600            &           ustar2_surf(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + gdept_n(:,:,1))**(rnn-1.) 
     
    855869      REAL(wp)::   zcr   ! local scalar 
    856870      !! 
    857       NAMELIST/namzdf_gls/rn_emin, rn_epsmin, ln_length_lim, & 
    858          &            rn_clim_galp, ln_sigpsi, rn_hsro,      & 
    859          &            rn_crban, rn_charn, rn_frac_hs,        & 
    860          &            nn_bc_surf, nn_bc_bot, nn_z0_met,     & 
     871      NAMELIST/namzdf_gls/rn_emin, rn_epsmin, ln_length_lim,       & 
     872         &            rn_clim_galp, ln_sigpsi, rn_hsro, rn_hsri,   & 
     873         &            rn_crban, rn_charn, rn_frac_hs,              & 
     874         &            nn_bc_surf, nn_bc_bot, nn_z0_met, nn_z0_ice, & 
    861875         &            nn_stab_func, nn_clos 
    862876      !!---------------------------------------------------------- 
     
    886900         WRITE(numout,*) '      Charnock coefficient                          rn_charn       = ', rn_charn 
    887901         WRITE(numout,*) '      Surface roughness formula                     nn_z0_met      = ', nn_z0_met 
     902         WRITE(numout,*) '      surface wave breaking under ice               nn_z0_ice      = ', nn_z0_ice 
     903         SELECT CASE( nn_z0_ice ) 
     904         CASE( 0 )   ;   WRITE(numout,*) '   ==>>>   no impact of ice cover on surface wave breaking' 
     905         CASE( 1 )   ;   WRITE(numout,*) '   ==>>>   roughness uses rn_hsri and is weigthed by 1-TANH( fr_i(:,:) * 10 )' 
     906         CASE( 2 )   ;   WRITE(numout,*) '   ==>>>   roughness uses rn_hsri and is weighted by 1-fr_i(:,:)' 
     907         CASE( 3 )   ;   WRITE(numout,*) '   ==>>>   roughness uses rn_hsri and is weighted by 1-MIN( 1, 4 * fr_i(:,:) )' 
     908         CASE DEFAULT 
     909            CALL ctl_stop( 'zdf_gls_init: wrong value for nn_z0_ice, should be 0,1,2, or 3') 
     910         END SELECT 
    888911         WRITE(numout,*) '      Wave height frac. (used if nn_z0_met=2)       rn_frac_hs     = ', rn_frac_hs 
    889912         WRITE(numout,*) '      Stability functions                           nn_stab_func   = ', nn_stab_func 
    890913         WRITE(numout,*) '      Type of closure                               nn_clos        = ', nn_clos 
    891914         WRITE(numout,*) '      Surface roughness (m)                         rn_hsro        = ', rn_hsro 
     915         WRITE(numout,*) '      Ice-ocean roughness (used if nn_z0_ice/=0)    rn_hsri        = ', rn_hsri 
    892916         WRITE(numout,*) 
    893917         WRITE(numout,*) '   Namelist namdrg_top/_bot:   used values:' 
  • NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF/zdfphy.F90

    r11536 r13284  
    2828   USE sbc_oce        ! surface module (only for nn_isf in the option compatibility test) 
    2929   USE sbcrnf         ! surface boundary condition: runoff variables 
     30   USE sbc_ice        ! sea ice drag 
    3031#if defined key_agrif 
    3132   USE agrif_oce_interp   ! interpavm 
     
    252253      ENDIF 
    253254      ! 
     255#if defined key_si3 
     256      IF ( ln_drgice_imp) THEN 
     257         IF ( ln_isfcav ) THEN 
     258            rCdU_top(:,:) = rCdU_top(:,:) + ssmask(:,:) * tmask(:,:,1) * rCdU_ice(:,:) 
     259         ELSE 
     260            rCdU_top(:,:) = rCdU_ice(:,:) 
     261         ENDIF 
     262      ENDIF 
     263#endif 
     264      !  
    254265      !                       !==  Kz from chosen turbulent closure  ==!   (avm_k, avt_k) 
    255266      ! 
  • NEMO/releases/r4.0/r4.0-HEAD/src/OCE/ZDF/zdftke.F90

    r13268 r13284  
    4646   USE zdfmxl         ! vertical physics: mixed layer 
    4747   ! 
     48#if defined key_si3 
     49   USE ice, ONLY: hm_i, h_i 
     50#endif 
     51#if defined key_cice 
     52   USE sbc_ice, ONLY: h_i 
     53#endif 
    4854   USE in_out_manager ! I/O manager 
    4955   USE iom            ! I/O manager library 
     
    6268   !                      !!** Namelist  namzdf_tke  ** 
    6369   LOGICAL  ::   ln_mxl0   ! mixing length scale surface value as function of wind stress or not 
     70   INTEGER  ::   nn_mxlice ! type of scaling under sea-ice (=0/1/2/3) 
     71   REAL(wp) ::   rn_mxlice ! ice thickness value when scaling under sea-ice 
    6472   INTEGER  ::   nn_mxl    ! type of mixing length (=0/1/2/3) 
    6573   REAL(wp) ::   rn_mxl0   ! surface  min value of mixing length (kappa*z_o=0.4*0.1 m)  [m] 
     
    7482   INTEGER  ::      nn_htau   ! type of tke profile of penetration (=0/1) 
    7583   REAL(wp) ::      rn_efr    ! fraction of TKE surface value which penetrates in the ocean 
    76    REAL(wp) ::      rn_eice   ! =0 ON below sea-ice, =4 OFF when ice fraction > 1/4    
    7784   LOGICAL  ::   ln_lc     ! Langmuir cells (LC) as a source term of TKE or not 
    7885   REAL(wp) ::      rn_lc     ! coef to compute vertical velocity of Langmuir cells 
     86   INTEGER  ::   nn_eice   ! attenutaion of langmuir & surface wave breaking under ice (=0/1/2/3)    
    7987 
    8088   REAL(wp) ::   ri_cri    ! critic Richardson number (deduced from rn_ediff and rn_ediss values) 
     
    190198      REAL(wp), DIMENSION(:,:,:) , INTENT(in   ) ::   p_avm, p_avt   ! vertical eddy viscosity & diffusivity (w-points) 
    191199      ! 
    192       INTEGER ::   ji, jj, jk              ! dummy loop arguments 
     200      INTEGER ::   ji, jj, jk                  ! dummy loop arguments 
    193201      REAL(wp) ::   zetop, zebot, zmsku, zmskv ! local scalars 
    194202      REAL(wp) ::   zrhoa  = 1.22              ! Air density kg/m3 
    195203      REAL(wp) ::   zcdrag = 1.5e-3            ! drag coefficient 
    196       REAL(wp) ::   zbbrau, zri                ! local scalars 
    197       REAL(wp) ::   zfact1, zfact2, zfact3     !   -         - 
    198       REAL(wp) ::   ztx2  , zty2  , zcof       !   -         - 
    199       REAL(wp) ::   ztau  , zdif               !   -         - 
    200       REAL(wp) ::   zus   , zwlc  , zind       !   -         - 
    201       REAL(wp) ::   zzd_up, zzd_lw             !   -         - 
     204      REAL(wp) ::   zbbrau, zbbirau, zri       ! local scalars 
     205      REAL(wp) ::   zfact1, zfact2, zfact3     !   -      - 
     206      REAL(wp) ::   ztx2  , zty2  , zcof       !   -      - 
     207      REAL(wp) ::   ztau  , zdif               !   -      - 
     208      REAL(wp) ::   zus   , zwlc  , zind       !   -      - 
     209      REAL(wp) ::   zzd_up, zzd_lw             !   -      - 
    202210      INTEGER , DIMENSION(jpi,jpj)     ::   imlc 
    203       REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc, zfr_i 
     211      REAL(wp), DIMENSION(jpi,jpj)     ::   zice_fra, zhlc, zus3 
    204212      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zpelc, zdiag, zd_up, zd_lw 
    205213      !!-------------------------------------------------------------------- 
    206214      ! 
    207       zbbrau = rn_ebb / rau0       ! Local constant initialisation 
    208       zfact1 = -.5_wp * rdt  
    209       zfact2 = 1.5_wp * rdt * rn_ediss 
    210       zfact3 = 0.5_wp       * rn_ediss 
     215      zbbrau  =  rn_ebb / rau0       ! Local constant initialisation 
     216      zbbirau =  3.75_wp / rau0 
     217      zfact1  = -0.5_wp * rdt  
     218      zfact2  =  1.5_wp * rdt * rn_ediss 
     219      zfact3  =  0.5_wp       * rn_ediss 
     220      ! 
     221      ! ice fraction considered for attenuation of langmuir & wave breaking 
     222      SELECT CASE ( nn_eice ) 
     223      CASE( 0 )   ;   zice_fra(:,:) = 0._wp 
     224      CASE( 1 )   ;   zice_fra(:,:) =        TANH( fr_i(:,:) * 10._wp ) 
     225      CASE( 2 )   ;   zice_fra(:,:) =              fr_i(:,:) 
     226      CASE( 3 )   ;   zice_fra(:,:) = MIN( 4._wp * fr_i(:,:) , 1._wp ) 
     227      END SELECT 
    211228      ! 
    212229      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    213230      !                     !  Surface/top/bottom boundary condition on tke 
    214231      !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    215        
     232      ! 
    216233      DO jj = 2, jpjm1            ! en(1)   = rn_ebb taum / rau0  (min value rn_emin0) 
    217234         DO ji = fs_2, fs_jpim1   ! vector opt. 
    218             en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1) 
     235            en(ji,jj,1) = MAX( rn_emin0, ( ( 1._wp - fr_i(ji,jj) ) * zbbrau + & 
     236               &                                     fr_i(ji,jj)   * zbbirau ) * taum(ji,jj) ) * tmask(ji,jj,1) 
    219237         END DO 
    220238      END DO 
     
    248266                  zetop = - 0.001875_wp * rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2  & 
    249267                     &                                           + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  ) 
    250                   en(ji,jj,mikt(ji,jj)) = en(ji,jj,1)           * tmask(ji,jj,1) & 
     268                  en(ji,jj,mikt(ji,jj)) = en(ji,jj,1)           * tmask(ji,jj,1) &      
    251269                     &                  + MAX( zetop, rn_emin ) * (1._wp - tmask(ji,jj,1)) * ssmask(ji,jj) 
    252270               END DO 
     
    286304            DO ji = fs_2, fs_jpim1   ! vector opt. 
    287305               zus  = zcof * SQRT( taum(ji,jj) )           ! Stokes drift 
    288                zfr_i(ji,jj) = ( 1._wp - 4._wp * fr_i(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok 
    289                IF (zfr_i(ji,jj) < 0. ) zfr_i(ji,jj) = 0. 
     306               zus3(ji,jj) = ( 1._wp - zice_fra(ji,jj) ) * zus * zus * zus * tmask(ji,jj,1) ! zus > 0. ok 
    290307            END DO 
    291308         END DO          
     
    293310            DO jj = 2, jpjm1 
    294311               DO ji = fs_2, fs_jpim1   ! vector opt. 
    295                   IF ( zfr_i(ji,jj) /= 0. ) THEN                
     312                  IF ( zus3(ji,jj) /= 0. ) THEN                
    296313                     ! vertical velocity due to LC    
    297314                     IF ( pdepw(ji,jj,jk) - zhlc(ji,jj) < 0 .AND. wmask(ji,jj,jk) /= 0. ) THEN 
    298315                        !                                           ! vertical velocity due to LC 
    299                         zwlc = rn_lc * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) )   ! warning: optimization: zus^3 is in zfr_i 
     316                        zwlc = rn_lc * SIN( rpi * pdepw(ji,jj,jk) / zhlc(ji,jj) ) 
    300317                        !                                           ! TKE Langmuir circulation source term 
    301                         en(ji,jj,jk) = en(ji,jj,jk) + rdt * zfr_i(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj) 
     318                        en(ji,jj,jk) = en(ji,jj,jk) + rdt * zus3(ji,jj) * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj) 
    302319                     ENDIF 
    303320                  ENDIF 
     
    399416       
    400417      IF( nn_etau == 1 ) THEN           !* penetration below the mixed layer (rn_efr fraction) 
    401          DO jk = 2, jpkm1                       ! rn_eice =0 ON below sea-ice, =4 OFF when ice fraction > 0.25 
     418         DO jk = 2, jpkm1                       ! nn_eice=0 : ON below sea-ice ; nn_eice>0 : partly OFF 
    402419            DO jj = 2, jpjm1 
    403420               DO ji = fs_2, fs_jpim1   ! vector opt. 
    404421                  en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -pdepw(ji,jj,jk) / htau(ji,jj) )   & 
    405                      &                                 * MAX(0.,1._wp - rn_eice *fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     422                     &                                 * ( 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    406423               END DO 
    407424            END DO 
     
    412429               jk = nmln(ji,jj) 
    413430               en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -pdepw(ji,jj,jk) / htau(ji,jj) )   & 
    414                   &                                 * MAX(0.,1._wp - rn_eice *fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     431                  &                                 * ( 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    415432            END DO 
    416433         END DO 
     
    425442                  zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add )  ! apply some modifications... 
    426443                  en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -pdepw(ji,jj,jk) / htau(ji,jj) )   & 
    427                      &                        * MAX(0.,1._wp - rn_eice *fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     444                     &                                 * ( 1._wp - zice_fra(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    428445               END DO 
    429446            END DO 
     
    477494      REAL(wp) ::   zrn2, zraug, zcoef, zav   ! local scalars 
    478495      REAL(wp) ::   zdku,   zdkv, zsqen       !   -      - 
    479       REAL(wp) ::   zemxl, zemlm, zemlp       !   -      - 
     496      REAL(wp) ::   zemxl, zemlm, zemlp, zmaxice       !   -      - 
    480497      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmxlm, zmxld   ! 3D workspace 
    481498      !!-------------------------------------------------------------------- 
     
    490507      zmxlm(:,:,:)  = rmxl_min     
    491508      zmxld(:,:,:)  = rmxl_min 
    492       ! 
     509      !  
    493510      IF( ln_mxl0 ) THEN            ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rau0*g) 
     511         ! 
    494512         zraug = vkarmn * 2.e5_wp / ( rau0 * grav ) 
    495          DO jj = 2, jpjm1 
     513#if ! defined key_si3 && ! defined key_cice 
     514         DO jj = 2, jpjm1                     ! No sea-ice 
    496515            DO ji = fs_2, fs_jpim1 
    497                zmxlm(ji,jj,1) = MAX( rn_mxl0, zraug * taum(ji,jj) * tmask(ji,jj,1) ) 
    498             END DO 
    499          END DO 
    500       ELSE  
     516               zmxlm(ji,jj,1) =  zraug * taum(ji,jj) * tmask(ji,jj,1) 
     517            END DO 
     518         END DO 
     519#else 
     520 
     521         SELECT CASE( nn_mxlice )             ! Type of scaling under sea-ice 
     522         ! 
     523         CASE( 0 )                      ! No scaling under sea-ice 
     524            DO jj = 2, jpjm1 
     525               DO ji = fs_2, fs_jpim1 
     526                  zmxlm(ji,jj,1) = zraug * taum(ji,jj) * tmask(ji,jj,1) 
     527               END DO 
     528            END DO 
     529            ! 
     530         CASE( 1 )                      ! scaling with constant sea-ice thickness 
     531            DO jj = 2, jpjm1 
     532               DO ji = fs_2, fs_jpim1 
     533                  zmxlm(ji,jj,1) =  ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 
     534                     &                          fr_i(ji,jj)   * rn_mxlice           ) * tmask(ji,jj,1) 
     535               END DO 
     536            END DO 
     537            ! 
     538         CASE( 2 )                      ! scaling with mean sea-ice thickness 
     539            DO jj = 2, jpjm1 
     540               DO ji = fs_2, fs_jpim1 
     541#if defined key_si3 
     542                  zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 
     543                     &                         fr_i(ji,jj)   * hm_i(ji,jj) * 2._wp ) * tmask(ji,jj,1) 
     544#elif defined key_cice 
     545                  zmaxice = MAXVAL( h_i(ji,jj,:) ) 
     546                  zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 
     547                     &                         fr_i(ji,jj)   * zmaxice             ) * tmask(ji,jj,1) 
     548#endif 
     549               END DO 
     550            END DO 
     551            ! 
     552         CASE( 3 )                      ! scaling with max sea-ice thickness 
     553            DO jj = 2, jpjm1 
     554               DO ji = fs_2, fs_jpim1 
     555                  zmaxice = MAXVAL( h_i(ji,jj,:) ) 
     556                  zmxlm(ji,jj,1) = ( ( 1._wp - fr_i(ji,jj) ) * zraug * taum(ji,jj) + & 
     557                     &                         fr_i(ji,jj)   * zmaxice             ) * tmask(ji,jj,1) 
     558               END DO 
     559            END DO 
     560            ! 
     561         END SELECT 
     562#endif 
     563         ! 
     564         DO jj = 2, jpjm1 
     565            DO ji = fs_2, fs_jpim1 
     566               zmxlm(ji,jj,1) = MAX( rn_mxl0, zmxlm(ji,jj,1) ) 
     567            END DO 
     568         END DO 
     569         ! 
     570      ELSE 
    501571         zmxlm(:,:,1) = rn_mxl0 
    502572      ENDIF 
     
    643713      INTEGER ::   ios 
    644714      !! 
    645       NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin  ,  & 
    646          &                 rn_emin0, rn_bshear, nn_mxl , ln_mxl0  ,  & 
    647          &                 rn_mxl0 , nn_pdl   , ln_lc  , rn_lc,      & 
    648          &                 nn_etau , nn_htau  , rn_efr , rn_eice   
     715      NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb   , rn_emin  ,  & 
     716         &                 rn_emin0, rn_bshear, nn_mxl   , ln_mxl0  ,  & 
     717         &                 rn_mxl0 , nn_mxlice, rn_mxlice,             & 
     718         &                 nn_pdl  , ln_lc    , rn_lc,                 & 
     719         &                 nn_etau , nn_htau  , rn_efr   , nn_eice   
    649720      !!---------------------------------------------------------------------- 
    650721      ! 
     
    675746         WRITE(numout,*) '         surface mixing length = F(stress) or not    ln_mxl0   = ', ln_mxl0 
    676747         WRITE(numout,*) '         surface  mixing length minimum value        rn_mxl0   = ', rn_mxl0 
     748         IF( ln_mxl0 ) THEN 
     749            WRITE(numout,*) '      type of scaling under sea-ice               nn_mxlice = ', nn_mxlice 
     750            IF( nn_mxlice == 1 ) & 
     751            WRITE(numout,*) '      ice thickness when scaling under sea-ice    rn_mxlice = ', rn_mxlice 
     752            SELECT CASE( nn_mxlice )             ! Type of scaling under sea-ice 
     753            CASE( 0 )   ;   WRITE(numout,*) '   ==>>>   No scaling under sea-ice' 
     754            CASE( 1 )   ;   WRITE(numout,*) '   ==>>>   scaling with constant sea-ice thickness' 
     755            CASE( 2 )   ;   WRITE(numout,*) '   ==>>>   scaling with mean sea-ice thickness' 
     756            CASE( 3 )   ;   WRITE(numout,*) '   ==>>>   scaling with max sea-ice thickness' 
     757            CASE DEFAULT 
     758               CALL ctl_stop( 'zdf_tke_init: wrong value for nn_mxlice, should be 0,1,2,3 or 4') 
     759            END SELECT 
     760         ENDIF 
    677761         WRITE(numout,*) '      Langmuir cells parametrization              ln_lc     = ', ln_lc 
    678762         WRITE(numout,*) '         coef to compute vertical velocity of LC     rn_lc  = ', rn_lc 
     
    680764         WRITE(numout,*) '          type of tke penetration profile            nn_htau   = ', nn_htau 
    681765         WRITE(numout,*) '          fraction of TKE that penetrates            rn_efr    = ', rn_efr 
    682          WRITE(numout,*) '          below sea-ice:  =0 ON                      rn_eice   = ', rn_eice 
    683          WRITE(numout,*) '          =4 OFF when ice fraction > 1/4   ' 
     766         WRITE(numout,*) '      langmuir & surface wave breaking under ice  nn_eice = ', nn_eice 
     767         SELECT CASE( nn_eice )  
     768         CASE( 0 )   ;   WRITE(numout,*) '   ==>>>   no impact of ice cover on langmuir & surface wave breaking' 
     769         CASE( 1 )   ;   WRITE(numout,*) '   ==>>>   weigthed by 1-TANH( fr_i(:,:) * 10 )' 
     770         CASE( 2 )   ;   WRITE(numout,*) '   ==>>>   weighted by 1-fr_i(:,:)' 
     771         CASE( 3 )   ;   WRITE(numout,*) '   ==>>>   weighted by 1-MIN( 1, 4 * fr_i(:,:) )' 
     772         CASE DEFAULT 
     773            CALL ctl_stop( 'zdf_tke_init: wrong value for nn_eice, should be 0,1,2, or 3') 
     774         END SELECT       
    684775         IF( .NOT.ln_drg_OFF ) THEN 
    685776            WRITE(numout,*) 
Note: See TracChangeset for help on using the changeset viewer.