Changeset 12966 for NEMO/branches

Timestamp:

2020-05-25T12:04:00+02:00 (4 years ago)

Author:

jchanut

Message:

Switch to log-law bcs below ice in tke and gls - hard coded 3cm ice-ocean roughness: #2468

Location:

NEMO/branches/MOI/NEMO_4.03_IODRAG/src/OCE/ZDF

Files:

• zdfgls.F90 (modified) (8 diffs)
• zdftke.F90 (modified) (7 diffs)

NEMO/branches/MOI/NEMO_4.03_IODRAG/src/OCE/ZDF/zdfgls.F90

@@ -146 +146 @@
       REAL(wp) ::   gh, gm, shr, dif, zsqen, zavt, zavm !   -      -
       REAL(wp) ::   zmsku, zmskv                        !   -      -
+      REAL(wp) ::   rn_hsri ! ice ocean roughness
       REAL(wp), DIMENSION(jpi,jpj)     ::   zdep
       REAL(wp), DIMENSION(jpi,jpj)     ::   zkar
       REAL(wp), DIMENSION(jpi,jpj)     ::   zflxs       ! Turbulence fluxed induced by internal waves 
       REAL(wp), DIMENSION(jpi,jpj)     ::   zhsro       ! Surface roughness (surface waves)
+      REAL(wp), DIMENSION(jpi,jpj)     ::   ziconc      ! Tapering of wave breaking under sea ice
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   eb          ! tke at time before
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   hmxl_b      ! mixing length at time before
@@ -165 +167 @@
       ustar2_bot (:,:) = 0._wp
 
+      ziconc(:,:) = 1._wp - TANH(fr_i(:,:) * 10._wp)
       ! Compute surface, top and bottom friction at T-points
       DO jj = 2, jpjm1          
@@ -206 +209 @@
       END SELECT
       !
+      rn_hsri = 0.03_wp 
+      zhsro(:,:) =  (ziconc(:,:) * zhsro(:,:) + (1._wp - ziconc(:,:)) * rn_hsri)*tmask(:,:,1)  + (1._wp - tmask(:,:,1)) * rn_hsro
       DO jk = 2, jpkm1              !==  Compute dissipation rate  ==!
          DO jj = 1, jpjm1
@@ -300 +305 @@
       CASE ( 0 )             ! Dirichlet boundary condition (set e at k=1 & 2) 
       ! First level
-      en   (:,:,1) = MAX(  rn_emin , rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3  )
+      en   (:,:,1) = MAX(  rn_emin , rc02r * ustar2_surf(:,:) * (1._wp + ziconc(:,:)*rsbc_tke1)**r2_3  )
       zd_lw(:,:,1) = en(:,:,1)
       zd_up(:,:,1) = 0._wp
@@ -306 +311 @@
       ! 
       ! One level below
-      en   (:,:,2) =  MAX(  rc02r * ustar2_surf(:,:) * (  1._wp + rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2))   &
+      en   (:,:,2) =  MAX(  rc02r * ustar2_surf(:,:) * (  1._wp + ziconc(:,:)*rsbc_tke1 * ((zhsro(:,:)+gdepw_n(:,:,2))   &
          &                 / zhsro(:,:) )**(1.5_wp*ra_sf)  )**(2._wp/3._wp)                      , rn_emin   )
       zd_lw(:,:,2) = 0._wp 
@@ -316 +321 @@
       !
       ! Dirichlet conditions at k=1
-      en   (:,:,1) = MAX(  rc02r * ustar2_surf(:,:) * (1._wp + rsbc_tke1)**r2_3 , rn_emin  )
+      en   (:,:,1) = MAX(  rc02r * ustar2_surf(:,:) * (1._wp + ziconc(:,:)*rsbc_tke1)**r2_3 , rn_emin  )
       zd_lw(:,:,1) = en(:,:,1)
       zd_up(:,:,1) = 0._wp
@@ -326 +331 @@
       zd_lw(:,:,2) = 0._wp
       zkar (:,:)   = (rl_sf + (vkarmn-rl_sf)*(1.-EXP(-rtrans*gdept_n(:,:,1)/zhsro(:,:)) ))
-      zflxs(:,:)   = rsbc_tke2 * ustar2_surf(:,:)**1.5_wp * zkar(:,:) &
+      zflxs(:,:)   = rsbc_tke2 * ziconc(:,:) * ustar2_surf(:,:)**1.5_wp * zkar(:,:) &
           &                    * (  ( zhsro(:,:)+gdept_n(:,:,1) ) / zhsro(:,:)  )**(1.5_wp*ra_sf)
 !!gm why not   :                        * ( 1._wp + gdept_n(:,:,1) / zhsro(:,:) )**(1.5_wp*ra_sf)
@@ -577 +582 @@
          zkar (:,:)   = rl_sf + (vkarmn-rl_sf)*(1._wp-EXP(-rtrans*gdept_n(:,:,1)/zhsro(:,:) )) ! Lengh scale slope
          zdep (:,:)   = ((zhsro(:,:) + gdept_n(:,:,1)) / zhsro(:,:))**(rmm*ra_sf)
-         zflxs(:,:)   = (rnn + rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp)
+         zflxs(:,:)   = (rnn + ziconc(:,:)*rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + ziconc(:,:)*rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp)
          zdep (:,:)   = rsbc_psi1 * (zwall_psi(:,:,1)*p_avm(:,:,1)+zwall_psi(:,:,2)*p_avm(:,:,2)) * &
             &           ustar2_surf(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + gdept_n(:,:,1))**(rnn-1.)

NEMO/branches/MOI/NEMO_4.03_IODRAG/src/OCE/ZDF/zdftke.F90

@@ -195 +195 @@
       REAL(wp) ::   zrhoa  = 1.22              ! Air density kg/m3
       REAL(wp) ::   zcdrag = 1.5e-3            ! drag coefficient
-      REAL(wp) ::   zbbrau, zri                ! local scalars
+      REAL(wp) ::   zbbrau, zbbrau0, zri       ! local scalars
       REAL(wp) ::   zfact1, zfact2, zfact3     !   -         -
       REAL(wp) ::   ztx2  , zty2  , zcof       !   -         -
@@ -202 +202 @@
       REAL(wp) ::   zzd_up, zzd_lw             !   -         -
       INTEGER , DIMENSION(jpi,jpj)     ::   imlc
+      REAL(wp), DIMENSION(jpi,jpj)     ::   ziconc
       REAL(wp), DIMENSION(jpi,jpj)     ::   zhlc, zfr_i
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zpelc, zdiag, zd_up, zd_lw
@@ -207 +208 @@
       !
       zbbrau = rn_ebb / rau0       ! Local constant initialisation
+      zbbrau0 = 3.75_wp / rau0
       zfact1 = -.5_wp * rdt 
       zfact2 = 1.5_wp * rdt * rn_ediss
@@ -214 +216 @@
       !                     !  Surface/top/bottom boundary condition on tke
       !                     !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
-      
+      ziconc(:,:) = 1._wp - TANH(fr_i(:,:) * 10._wp)      
+      !
       DO jj = 2, jpjm1            ! en(1)   = rn_ebb taum / rau0  (min value rn_emin0)
          DO ji = fs_2, fs_jpim1   ! vector opt.
-            en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1)
+            en(ji,jj,1) = MAX( rn_emin0, (ziconc(ji,jj)*zbbrau + (1._wp-ziconc(ji,jj))*zbbrau0) &
+                             & * taum(ji,jj) ) * tmask(ji,jj,1)
          END DO
       END DO
@@ -249 +253 @@
                   zetop = - 0.001875_wp * rCdU_top(ji,jj) * SQRT(  ( zmsku*( ub(ji,jj,mikt(ji,jj))+ub(ji-1,jj,mikt(ji,jj)) ) )**2  &
                      &                                           + ( zmskv*( vb(ji,jj,mikt(ji,jj))+vb(ji,jj-1,mikt(ji,jj)) ) )**2  )
-                  en(ji,jj,mikt(ji,jj)) = en(ji,jj,1)           * tmask(ji,jj,1) &
+                  en(ji,jj,mikt(ji,jj)) = en(ji,jj,1)           * tmask(ji,jj,1) &     
                      &                  + MAX( zetop, rn_emin ) * (1._wp - tmask(ji,jj,1)) * ssmask(ji,jj)
                END DO
@@ -479 +483 @@
       REAL(wp) ::   zdku,   zdkv, zsqen       !   -      -
       REAL(wp) ::   zemxl, zemlm, zemlp       !   -      -
+      REAL(wp) ::   rn_hsri
+      REAL(wp), DIMENSION(jpi,jpj) :: ziconc
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmxlm, zmxld   ! 3D workspace
       !!--------------------------------------------------------------------
@@ -502 +508 @@
          zmxlm(:,:,1) = rn_mxl0
       ENDIF
+      ziconc(:,:) = 1._wp - TANH(fr_i(:,:) * 10._wp)
+      rn_hsri = vkarmn * 0.03_wp ! Sea ice roughness
+      zmxlm(:,:,1) =  (ziconc(:,:) * zmxlm(:,:,1) + (1._wp - ziconc(:,:)) * rn_hsri)*tmask(:,:,1)  &
+                 & + (1._wp - tmask(:,:,1)) * rn_mxl0
       !
       DO jk = 2, jpkm1              ! interior value : l=sqrt(2*e/n^2)