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 13403 for NEMO/branches/UKMO/NEMO_4.0.1_FKOSM_m11715/src/OCE/ZDF/zdfosm.F90 – NEMO

Ignore:
Timestamp:
2020-08-14T17:08:06+02:00 (4 years ago)
Author:
dancopsey
Message:

Merge in George Nurser's improvements to Stokes Drift in OSMOSIS.
Custom merge from revision 13392 to revision 13400 of:
http://forge.ipsl.jussieu.fr/nemo/log/NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser_4.0

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0.1_FKOSM_m11715/src/OCE/ZDF/zdfosm.F90

    r13402 r13403  
    103103   REAL(wp) ::   rn_osm_la          ! Turbulent Langmuir number 
    104104   REAL(wp) ::   rn_osm_dstokes     ! Depth scale of Stokes drift 
     105   REAL(wp) ::   rn_zdfosm_adjust_sd = 1.0 ! factor to reduce Stokes drift by 
     106   REAL(wp) ::   rn_osm_hblfrac = 0.1! for nn_osm_wave = 3/4 specify fraction in top of hbl 
     107   LOGICAL  ::   ln_zdfosm_ice_shelter      ! flag to activate ice sheltering 
    105108   REAL(wp) ::   rn_osm_hbl0 = 10._wp       ! Initial value of hbl for 1D runs 
    106109   INTEGER  ::   nn_ave             ! = 0/1 flag for horizontal average on avt 
    107110   INTEGER  ::   nn_osm_wave = 0    ! = 0/1/2 flag for getting stokes drift from La# / PM wind-waves/Inputs into sbcwave 
     111   INTEGER  ::   nn_osm_SD_reduce   ! = 0/1/2 flag for getting effective stokes drift from surface value 
    108112   LOGICAL  ::   ln_dia_osm         ! Use namelist  rn_osm_la 
    109  
    110113 
    111114   LOGICAL  ::   ln_kpprimix  = .true.  ! Shear instability mixing 
     
    123126   !                                        ! parameters used in nn_osm_mle = 1 case 
    124127   REAL(wp) ::   rn_osm_mle_lat              ! reference latitude for a 5 km scale of ML front 
     128   LOGICAL  ::   ln_osm_hmle_limit           ! If true arbitrarily restrict hmle to rn_osm_hmle_limit*zmld 
     129   REAL(wp) ::   rn_osm_hmle_limit           ! If ln_osm_hmle_limit true arbitrarily restrict hmle to rn_osm_hmle_limit*zmld 
    125130   REAL(wp) ::   rn_osm_mle_rho_c        ! Density criterion for definition of MLD used by FK 
    126131   REAL(wp) ::   r5_21 = 5.e0 / 21.e0   ! factor used in mle streamfunction computation 
     
    321326      REAL(wp) :: zzeta_v = 0.46 
    322327      REAL(wp) :: zabsstke 
     328      REAL(wp) :: zsqrtpi, z_two_thirds, zproportion, ztransp, zthickness 
     329      REAL(wp) :: z2k_times_thickness, zsqrt_depth, zexp_depth, zdstokes0, zf, zexperfc 
    323330 
    324331      ! For debugging 
     
    407414           zwbav(ji,jj) = grav  * zthermal * zwthav(ji,jj) - grav  * zbeta * zwsav(ji,jj) 
    408415           ! Surface upward velocity fluxes 
    409            zuw0(ji,jj) = -utau(ji,jj) * r1_rau0 * tmask(ji,jj,1) 
    410            zvw0(ji,jj) = -vtau(ji,jj) * r1_rau0 * tmask(ji,jj,1) 
     416           zuw0(ji,jj) = - 0.5 * (utau(ji-1,jj) + utau(ji,jj)) * r1_rau0 * tmask(ji,jj,1) 
     417           zvw0(ji,jj) = - 0.5 * (vtau(ji,jj-1) + vtau(ji,jj)) * r1_rau0 * tmask(ji,jj,1) 
    411418           ! Friction velocity (zustar), at T-point : LMD94 eq. 2 
    412419           zustar(ji,jj) = MAX( SQRT( SQRT( zuw0(ji,jj) * zuw0(ji,jj) + zvw0(ji,jj) * zvw0(ji,jj) ) ), 1.0e-8 ) 
     
    423430              zus_x = zcos_wind(ji,jj) * zustar(ji,jj) / 0.3**2 
    424431              zus_y = zsin_wind(ji,jj) * zustar(ji,jj) / 0.3**2 
     432              ! Linearly 
    425433              zustke(ji,jj) = MAX ( SQRT( zus_x*zus_x + zus_y*zus_y), 1.0e-8 ) 
    426               ! dstokes(ji,jj) set to constant value rn_osm_dstokes from namelist in zdf_osm_init 
     434              dstokes(ji,jj) = rn_osm_dstokes 
    427435           END DO 
    428436        END DO 
     
    432440           DO ji = 2, jpim1 
    433441              ! Use wind speed wndm included in sbc_oce module 
    434               zustke(ji,jj) = MAX ( 0.016 * wndm(ji,jj), 1.0e-8 ) 
    435               dstokes(ji,jj) = MAX( 0.12 * wndm(ji,jj)**2 / grav, 5.e-1) 
     442              zustke(ji,jj) =  MAX ( 0.016 * wndm(ji,jj), 1.0e-8 ) 
     443              dstokes(ji,jj) = MAX ( 0.12 * wndm(ji,jj)**2 / grav, 5.e-1) 
    436444           END DO 
    437445        END DO 
     
    439447     CASE(2) 
    440448        zfac =  2.0_wp * rpi / 16.0_wp 
     449 
    441450        DO jj = 2, jpjm1 
    442451           DO ji = 2, jpim1 
    443               ! The Langmur number from the ECMWF model appears to give La<0.3 for wind-driven seas. 
    444               !    The coefficient 0.8 gives La=0.3  in this situation. 
    445               ! It could represent the effects of the spread of wave directions 
    446               ! around the mean wind. The effect of this adjustment needs to be tested. 
    447               zabsstke = SQRT(ut0sd(ji,jj)**2 + vt0sd(ji,jj)**2) 
    448               zustke(ji,jj) = MAX (0.8 * ( zcos_wind(ji,jj) * ut0sd(ji,jj) + zsin_wind(ji,jj)  * vt0sd(ji,jj) ), 1.0e-8) 
    449               dstokes(ji,jj) = MAX(zfac * hsw(ji,jj)*hsw(ji,jj) / ( MAX(zabsstke*wmp(ji,jj), 1.0e-7 ) ), 5.0e-1) !rn_osm_dstokes ! 
     452              IF (hsw(ji,jj) > 1.e-4) THEN 
     453                 ! Use  wave fields 
     454                 zabsstke = SQRT(ut0sd(ji,jj)**2 + vt0sd(ji,jj)**2) 
     455                 zustke(ji,jj) = MAX ( ( zcos_wind(ji,jj) * ut0sd(ji,jj) + zsin_wind(ji,jj)  * vt0sd(ji,jj) ), 1.0e-8) 
     456                 dstokes(ji,jj) = MAX (zfac * hsw(ji,jj)*hsw(ji,jj) / ( MAX(zabsstke * wmp(ji,jj), 1.0e-7 ) ), 5.0e-1) 
     457              ELSE 
     458                 ! Assume masking issue (e.g. ice in ECMWF reanalysis but not in model run) 
     459                 ! .. so default to Pierson-Moskowitz 
     460                 zustke(ji,jj) = MAX ( 0.016 * wndm(ji,jj), 1.0e-8 ) 
     461                 dstokes(ji,jj) = MAX ( 0.12 * wndm(ji,jj)**2 / grav, 5.e-1) 
     462              END IF 
     463            END DO 
     464        END DO 
     465     END SELECT 
     466 
     467     IF (ln_zdfosm_ice_shelter) THEN 
     468        ! Reduce both Stokes drift and its depth scale by ocean fraction to represent sheltering by ice 
     469        DO jj = 2, jpjm1 
     470           DO ji = 2, jpim1 
     471              zustke(ji,jj) = zustke(ji,jj) * (1.0_wp - fr_i(ji,jj)) 
     472              dstokes(ji,jj) = dstokes(ji,jj) * (1.0_wp - fr_i(ji,jj)) 
     473           END DO 
     474        END DO 
     475     END IF 
     476 
     477     SELECT CASE (nn_osm_SD_reduce) 
     478     ! Reduce surface Stokes drift by a constant factor or following Breivik (2016) + van  Roekel (2012) or Grant (2020). 
     479     CASE(0) 
     480        ! The Langmur number from the ECMWF model (or from PM)  appears to give La<0.3 for wind-driven seas. 
     481        !    The coefficient rn_zdfosm_adjust_sd = 0.8 gives La=0.3  in this situation. 
     482        ! It could represent the effects of the spread of wave directions 
     483        ! around the mean wind. The effect of this adjustment needs to be tested. 
     484        IF(nn_osm_wave > 0) THEN 
     485           zustke(2:jpim1,2:jpjm1) = rn_zdfosm_adjust_sd * zustke(2:jpim1,2:jpjm1) 
     486        END IF 
     487     CASE(1) 
     488        ! van  Roekel (2012): consider average SD over top 10% of boundary layer 
     489        ! assumes approximate depth profile of SD from Breivik (2016) 
     490        zsqrtpi = SQRT(rpi) 
     491        z_two_thirds = 2.0_wp / 3.0_wp 
     492 
     493        DO jj = 2, jpjm1 
     494           DO ji = 2, jpim1 
     495              zthickness = rn_osm_hblfrac*hbl(ji,jj) 
     496              z2k_times_thickness =  zthickness * 2.0_wp / MAX( ABS( 5.97_wp * dstokes(ji,jj) ), 0.0000001_wp ) 
     497              zsqrt_depth = SQRT(z2k_times_thickness) 
     498              zexp_depth  = EXP(-z2k_times_thickness) 
     499              zustke(ji,jj) = zustke(ji,jj) * (1.0_wp - zexp_depth  & 
     500                   &              - z_two_thirds * ( zsqrtpi*zsqrt_depth*z2k_times_thickness * ERFC(zsqrt_depth) & 
     501                   &              + 1.0_wp - (1.0_wp + z2k_times_thickness)*zexp_depth ) ) / z2k_times_thickness 
     502 
     503           END DO 
     504        END DO 
     505     CASE(2) 
     506        ! Grant (2020): Match to exponential with same SD and d/dz(Sd) at depth 10% of boundary layer 
     507        ! assumes approximate depth profile of SD from Breivik (2016) 
     508        zsqrtpi = SQRT(rpi) 
     509 
     510        DO jj = 2, jpjm1 
     511           DO ji = 2, jpim1 
     512              zthickness = rn_osm_hblfrac*hbl(ji,jj) 
     513              z2k_times_thickness =  zthickness * 2.0_wp / MAX( ABS( 5.97_wp * dstokes(ji,jj) ), 0.0000001_wp ) 
     514 
     515              IF(z2k_times_thickness < 50._wp) THEN 
     516                 zsqrt_depth = SQRT(z2k_times_thickness) 
     517                 zexperfc = zsqrtpi * zsqrt_depth * ERFC(zsqrt_depth) * EXP(z2k_times_thickness) 
     518              ELSE 
     519                 ! asymptotic expansion of sqrt(pi)*zsqrt_depth*EXP(z2k_times_thickness)*ERFC(zsqrt_depth) for large z2k_times_thickness 
     520                 ! See Abramowitz and Stegun, Eq. 7.1.23 
     521                 ! zexperfc = 1._wp - (1/2)/(z2k_times_thickness)  + (3/4)/(z2k_times_thickness**2) - (15/8)/(z2k_times_thickness**3) 
     522                 zexperfc = ((- 1.875_wp/z2k_times_thickness + 0.75_wp)/z2k_times_thickness - 0.5_wp)/z2k_times_thickness + 1.0_wp 
     523              END IF 
     524              zf = z2k_times_thickness*(1.0_wp/zexperfc - 1.0_wp) 
     525              dstokes(ji,jj) = 5.97 * zf * dstokes(ji,jj) 
     526              zustke(ji,jj) = zustke(ji,jj) * EXP(z2k_times_thickness * ( 1.0_wp / (2. * zf) - 1.0_wp )) * ( 1.0_wp - zexperfc) 
    450527           END DO 
    451528        END DO 
     
    486563     ! zdf_osm_zmld_horizontal_gradients need halo values for ibld, so must 
    487564     ! previously exist for hbl also. 
     565 
     566     ! agn 23/6/20: not clear all this is needed, as hbl checked after it is re-calculated anyway 
     567     ! ########################################################################## 
    488568      hbl(:,:) = MAX(hbl(:,:), gdepw_n(:,:,4) ) 
    489569      ibld(:,:) = 4 
     
    497577         END DO 
    498578      END DO 
     579     ! ########################################################################## 
    499580 
    500581      DO jj = 2, jpjm1 
     
    10701151 
    10711152 
    1072   
     1153 
    10731154       IF ( ln_osm_mle ) THEN  ! set up diffusivity and non-gradient mixing 
    10741155          DO jj = 2 , jpjm1 
     
    11631244            IF ( iom_use("wind_wave_abs_power") ) CALL iom_put( "wind_wave_abs_power", 1000.*rau0*tmask(:,:,1)*zustar**2*zustke ) 
    11641245         ! Stokes drift read in from sbcwave  (=2). 
    1165          CASE(2) 
     1246         CASE(2:3) 
    11661247            IF ( iom_use("us_x") ) CALL iom_put( "us_x", ut0sd*umask(:,:,1) )               ! x surface Stokes drift 
    11671248            IF ( iom_use("us_y") ) CALL iom_put( "us_y", vt0sd*vmask(:,:,1) )               ! y surface Stokes drift 
     
    15181599             zrf_shear = TANH( ( zustar(ji,jj) / zwcor )**0.69 ) 
    15191600             zrf_langmuir = TANH( ( zwstrl(ji,jj) / zwcor )**0.69 ) 
    1520              zr_stokes = 1.0 - EXP( -25.0 * dstokes(ji,jj) / hbl(ji,jj) & 
    1521                   &                  * ( 1.0 + 4.0 * dstokes(ji,jj) / hbl(ji,jj) ) ) 
     1601             IF (nn_osm_SD_reduce > 0 ) THEN 
     1602                ! Effective Stokes drift already reduced from surface value 
     1603                zr_stokes = 1.0_wp 
     1604             ELSE 
     1605                ! Effective Stokes drift only reduced by factor rn_zdfodm_adjust_sd, 
     1606                ! requires further reduction where BL is deep 
     1607                zr_stokes = 1.0 - EXP( -25.0 * dstokes(ji,jj) / hbl(ji,jj) & 
     1608                     &                  * ( 1.0 + 4.0 * dstokes(ji,jj) / hbl(ji,jj) ) ) 
     1609             END IF 
    15221610 
    15231611             zwb_ent(ji,jj) = - 2.0 * 0.2 * zrf_conv * zwbav(ji,jj) & 
     
    15561644 
    15571645                  zvel_max = - ( 1.0 + 1.0 * ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird * rn_rdt / hbl(ji,jj) ) * zwb_ent(ji,jj) / & 
    1558                   &        ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird 
     1646                  &        MAX((zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird, epsln) 
    15591647                  zdhdt(ji,jj) = -zwb_ent(ji,jj) / ( zvel_max + MAX(zdb_bl(ji,jj), 1.0e-15) ) 
    15601648                ! added ajgn 23 July as temporay fix 
     
    15891677                    zpert = MAX( 2.0 * ( 1.0 + 0.0 * 2.0 * zvstr(ji,jj) * rn_rdt / hbl(ji,jj) ) * zvstr(ji,jj)**2 / hbl(ji,jj), zdb_bl(ji,jj) ) 
    15901678                ENDIF 
    1591                 zdhdt(ji,jj) = 2.0 * zdhdt(ji,jj) / zpert 
     1679                zdhdt(ji,jj) = 2.0 * zdhdt(ji,jj) / MAX(zpert, epsln) 
    15921680            ENDIF 
    15931681         END DO 
     
    17221810                     ! OSBL is deepening. Note wb_fk_b is zero if ln_osm_mle=F 
    17231811                     IF ( zdb_bl(ji,jj) > 0._wp .and. zdbdz_ext(ji,jj) > 0._wp)THEN 
    1724                         IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN  ! near neutral stability 
     1812                        IF ( ( zwstrc(ji,jj) / MAX(zvstr(ji,jj), epsln) )**3 <= 0.5 ) THEN  ! near neutral stability 
    17251813                           zari = MIN( 1.5 * ( zdb_bl(ji,jj) / ( zdbdz_ext(ji,jj) * zhbl(ji,jj) ) ) / & 
    1726                                 (1.0 + zdb_bl(ji,jj)**2 / ( 4.5 * zvstr(ji,jj)**2 * zdbdz_ext(ji,jj) ) ), 0.2 ) 
     1814                                (1.0 + zdb_bl(ji,jj)**2 / MAX( 4.5 * zvstr(ji,jj)**2 * zdbdz_ext(ji,jj), epsln ) ), 0.2 ) 
    17271815                        ELSE                                                     ! unstable 
    17281816                           zari = MIN( 1.5 * ( zdb_bl(ji,jj) / ( zdbdz_ext(ji,jj) * zhbl(ji,jj) ) ) / & 
    17291817                                (1.0 + zdb_bl(ji,jj)**2 / ( 4.5 * zwstrc(ji,jj)**2 * zdbdz_ext(ji,jj) ) ), 0.2 ) 
    17301818                        ENDIF 
    1731                         ztau = 0.2 * hbl(ji,jj) / (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird 
     1819                        ztau = 0.2 * hbl(ji,jj) / MAX(epsln, (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird) 
    17321820                        zddhdt = zari * hbl(ji,jj) 
    17331821                     ELSE 
    1734                         ztau = 0.2 * hbl(ji,jj) / ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird 
     1822                        ztau = 0.2 * hbl(ji,jj) / MAX(epsln, (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird) 
    17351823                        zddhdt = 0.2 * hbl(ji,jj) 
    17361824                     ENDIF 
    17371825                  ELSE 
    1738                      ztau = 0.2 * hbl(ji,jj) / ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird 
     1826                     ztau = 0.2 * hbl(ji,jj) /  MAX(epsln, (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird) 
    17391827                     zddhdt = 0.2 * hbl(ji,jj) 
    17401828                  ENDIF 
    17411829               ELSE ! ln_osm_mle 
    17421830                  IF ( zdb_bl(ji,jj) > 0._wp .and. zdbdz_ext(ji,jj) > 0._wp)THEN 
    1743                      IF ( ( zwstrc(ji,jj) / zvstr(ji,jj) )**3 <= 0.5 ) THEN  ! near neutral stability 
    1744                         zari = MIN( 1.5 * ( zdb_bl(ji,jj) / ( zdbdz_ext(ji,jj) * zhbl(ji,jj) ) ) / & 
    1745                              (1.0 + zdb_bl(ji,jj)**2 / ( 4.5 * zvstr(ji,jj)**2 * zdbdz_ext(ji,jj) ) ), 0.2 ) 
     1831                     IF ( ( zwstrc(ji,jj) / MAX(zvstr(ji,jj), epsln) )**3 <= 0.5 ) THEN  ! near neutral stability 
     1832                        zari = MIN( 1.5 * ( zdb_bl(ji,jj) / MAX( zdbdz_ext(ji,jj) * zhbl(ji,jj), epsln ) ) / & 
     1833                             (1.0 + zdb_bl(ji,jj)**2 /  MAX(4.5 * zvstr(ji,jj)**2 * zdbdz_ext(ji,jj), epsln ) ), 0.2 ) 
    17461834                     ELSE                                                     ! unstable 
    1747                         zari = MIN( 1.5 * ( zdb_bl(ji,jj) / ( zdbdz_ext(ji,jj) * zhbl(ji,jj) ) ) / & 
    1748                              (1.0 + zdb_bl(ji,jj)**2 / ( 4.5 * zwstrc(ji,jj)**2 * zdbdz_ext(ji,jj) ) ), 0.2 ) 
     1835                        zari = MIN( 1.5 * ( zdb_bl(ji,jj) / MAX( zdbdz_ext(ji,jj) * zhbl(ji,jj), epsln ) ) / & 
     1836                             (1.0 + zdb_bl(ji,jj)**2 / MAX(4.5 * zvstr(ji,jj)**2 * zdbdz_ext(ji,jj), epsln ) ), 0.2 ) 
    17491837                     ENDIF 
    1750                      ztau = hbl(ji,jj) / (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird 
     1838                     ztau = hbl(ji,jj) / MAX(epsln, (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird) 
    17511839                     zddhdt = zari * hbl(ji,jj) 
    17521840                  ELSE 
    1753                      ztau = hbl(ji,jj) / ( zvstr(ji,jj)**3 + 0.5 * zwstrc(ji,jj)**3 )**pthird 
     1841                     ztau = hbl(ji,jj) / MAX(epsln, (zvstr(ji,jj)**3 + 0.5 *zwstrc(ji,jj)**3)**pthird) 
    17541842                     zddhdt = 0.2 * hbl(ji,jj) 
    17551843                  ENDIF 
     
    17601848               ! Alan: this hml is never defined or used 
    17611849            ELSE   ! IF (lconv) 
    1762                ztau = hbl(ji,jj) / zvstr(ji,jj) 
     1850               ztau = hbl(ji,jj) / MAX(zvstr(ji,jj), epsln) 
    17631851               IF ( zdhdt(ji,jj) >= 0.0 ) THEN    ! probably shouldn't include wm here 
    17641852                  ! boundary layer deepening 
     
    17661854                     ! pycnocline thickness set by stratification - use same relationship as for neutral conditions. 
    17671855                     zari = MIN( 4.5 * ( zvstr(ji,jj)**2 ) & 
    1768                           & / ( zdb_bl(ji,jj) * zhbl(ji,jj) ) + 0.01  , 0.2 ) 
     1856                          & /  MAX(zdb_bl(ji,jj) * zhbl(ji,jj), epsln ) + 0.01  , 0.2 ) 
    17691857                     zddhdt = MIN( zari, 0.2 ) * hbl(ji,jj) 
    17701858                  ELSE 
     
    18211909      zN2_c = grav * rn_osm_mle_rho_c * r1_rau0   ! convert density criteria into N^2 criteria 
    18221910      DO jk = nlb10, jpkm1 
    1823          DO jj = 1, jpj                ! Mixed layer level: w-level  
     1911         DO jj = 1, jpj                ! Mixed layer level: w-level 
    18241912            DO ji = 1, jpi 
    18251913               ikt = mbkt(ji,jj) 
     
    19402028          ENDIF 
    19412029          hmle(ji,jj) = MIN(hmle(ji,jj), ht_n(ji,jj)) 
    1942           hmle(ji,jj) = MIN(hmle(ji,jj), 1.2*zmld(ji,jj))  
     2030          IF(ln_osm_hmle_limit) hmle(ji,jj) = MIN(hmle(ji,jj), rn_osm_hmle_limit*zmld(ji,jj)) 
    19432031        END DO 
    19442032      END DO 
     
    19702058                  & + dbdx_mle(ji-1,jj) * dbdx_mle(ji-1,jj) + dbdy_mle(ji,jj-1) * dbdy_mle(ji,jj-1) ) ) 
    19712059             zwb_fk(ji,jj) = rn_osm_mle_ce * hmle(ji,jj) * hmle(ji,jj) *ztmp * zdbds_mle * zdbds_mle 
    1972       ! This vbelocity scale, defined in Fox-Kemper et al (2008), is needed for calculating dhdt. 
    1973              zvel_mle(ji,jj) = zdbds_mle * ztmp * hmle(ji,jj) * tmask(ji,jj,1)  
     2060      ! This velocity scale, defined in Fox-Kemper et al (2008), is needed for calculating dhdt. 
     2061             zvel_mle(ji,jj) = zdbds_mle * ztmp * hmle(ji,jj) * tmask(ji,jj,1) 
    19742062             zdiff_mle(ji,jj) = 1.e-4_wp * ztmp * zdbds_mle * zhmle(ji,jj)**3  / rn_osm_mle_lf 
    19752063          ENDIF 
     
    19982086     !! 
    19992087     NAMELIST/namzdf_osm/ ln_use_osm_la, rn_osm_la, rn_osm_dstokes, nn_ave & 
    2000           & ,nn_osm_wave, ln_dia_osm, rn_osm_hbl0 & 
     2088          & ,nn_osm_wave, ln_dia_osm, rn_osm_hbl0, rn_zdfosm_adjust_sd & 
    20012089          & ,ln_kpprimix, rn_riinfty, rn_difri, ln_convmix, rn_difconv, nn_osm_wave & 
    2002           & ,ln_osm_mle 
     2090          & ,nn_osm_SD_reduce, ln_osm_mle, rn_osm_hblfrac, ln_zdfosm_ice_shelter 
    20032091! Namelist for Fox-Kemper parametrization. 
    20042092      NAMELIST/namosm_mle/ nn_osm_mle, rn_osm_mle_ce, rn_osm_mle_lf, rn_osm_mle_time, rn_osm_mle_lat,& 
    2005            & rn_osm_mle_rho_c,rn_osm_mle_thresh,rn_osm_mle_tau 
     2093           & rn_osm_mle_rho_c, rn_osm_mle_thresh, rn_osm_mle_tau, ln_osm_hmle_limit, rn_osm_hmle_limit 
    20062094 
    20072095     !!---------------------------------------------------------------------- 
     
    20242112        WRITE(numout,*) '     Use  MLE in OBL, i.e. Fox-Kemper param        ln_osm_mle = ', ln_osm_mle 
    20252113        WRITE(numout,*) '     Turbulent Langmuir number                     rn_osm_la   = ', rn_osm_la 
     2114        WRITE(numout,*) '     Stokes drift reduction factor                 rn_zdfosm_adjust_sd   = ', rn_zdfosm_adjust_sd 
    20262115        WRITE(numout,*) '     Initial hbl for 1D runs                       rn_osm_hbl0   = ', rn_osm_hbl0 
    20272116        WRITE(numout,*) '     Depth scale of Stokes drift                   rn_osm_dstokes = ', rn_osm_dstokes 
     
    20352124        CASE(2) 
    20362125           WRITE(numout,*) '     calculated from ECMWF wave fields' 
     2126         END SELECT 
     2127        WRITE(numout,*) '     Stokes drift reduction                        nn_osm_SD_reduce', nn_osm_SD_reduce 
     2128        WRITE(numout,*) '     fraction of hbl to average SD over/fit' 
     2129        WRITE(numout,*) '     exponential with nn_osm_SD_reduce = 1 or 2    rn_osm_hblfrac =  ', rn_osm_hblfrac 
     2130        SELECT CASE (nn_osm_SD_reduce) 
     2131        CASE(0) 
     2132           WRITE(numout,*) '     No reduction' 
     2133        CASE(1) 
     2134           WRITE(numout,*) '     Average SD over upper rn_osm_hblfrac of BL' 
     2135        CASE(2) 
     2136           WRITE(numout,*) '     Fit exponential to slope rn_osm_hblfrac of BL' 
    20372137        END SELECT 
     2138        WRITE(numout,*) '     reduce surface SD and depth scale under ice   ln_zdfosm_ice_shelter=', ln_zdfosm_ice_shelter 
    20382139        WRITE(numout,*) '     Output osm diagnostics                       ln_dia_osm  = ',  ln_dia_osm 
    20392140        WRITE(numout,*) '     Use KPP-style shear instability mixing       ln_kpprimix = ', ln_kpprimix 
     
    20722173            WRITE(numout,*) '         Threshold used to define ML for FK                      rn_osm_mle_thresh  = ', rn_osm_mle_thresh, 'm^2/s' 
    20732174            WRITE(numout,*) '         Timescale for OSM-FK                                         rn_osm_mle_tau  = ', rn_osm_mle_tau, 's' 
     2175            WRITE(numout,*) '         switch to limit hmle                                      ln_osm_hmle_limit  = ', ln_osm_hmle_limit 
     2176            WRITE(numout,*) '         fraction of zmld to limit hmle to if ln_osm_hmle_limit =.T.  rn_osm_hmle_limit = ', rn_osm_hmle_limit 
    20742177         ENDIF         ! 
    20752178     ENDIF 
Note: See TracChangeset for help on using the changeset viewer.