Changeset 2992 for branches

Timestamp:

2011-10-25T16:15:52+02:00 (13 years ago)

Author:

poddo

Message:

ticket #879. Step 2: Add in changes from the 2011/dev_r2802_INGV4_zdfRic branch

Location:

branches/2011/dev_INGV_2011

Files:

• DOC/TexFiles/Biblio/Biblio.bib (modified) (9 diffs)
• DOC/TexFiles/Chapters/Chap_SBC.tex (modified) (5 diffs)
• DOC/TexFiles/Chapters/Chap_ZDF.tex (modified) (1 diff)
• DOC/TexFiles/Namelist/namsbc (modified) (1 diff)
• DOC/TexFiles/Namelist/namsbc_ecmwf (copied) (copied from branches/2011/dev_r2802_INGV4_zdfRic/DOC/TexFiles/Namelist/namsbc_ecmwf)
• DOC/TexFiles/Namelist/namzdf_ric (modified) (1 diff)
• NEMOGCM/CONFIG/GYRE/EXP00/namelist (modified) (1 diff)
• NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist (modified) (1 diff)
• NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist (modified) (1 diff)
• NEMOGCM/CONFIG/POMME/EXP00/namelist (modified) (1 diff)
• NEMOGCM/NEMO/OPA_SRC/ZDF/zdfric.F90 (modified) (8 diffs)

branches/2011/dev_INGV_2011/DOC/TexFiles/Biblio/Biblio.bib

@@ -355 +355 @@
 }
 
+@ARTICLE{Bignami_al_JGR95,
+  author = {F. Bignami and S. Marullo and R. Santoleri and M. E. Schiano},
+  title = {Longwave radiation budget in the Mediterranean Sea}, 
+  journal = JGR,
+  year = {1995},
+  volume = {100},  number = {C2},
+  pages = {2501--2514},
+  doi = {10.1029/94JC02496},
+}
+
 @ARTICLE{Blanke_al_JPO99,
   author = {B. Blanke and M. Arhan and G. Madec and S. Roche},
@@ -590 +600 @@
   doi = {10.1029/2002GL016473},
   url = {http://dx.doi.org/10.1029/2002GL016473}
+}
+
+@ARTICLE{Castellari_al_JMS1998,
+  author = {S> Castellari and  N. Pinardi and K. Leaman },
+  title = {A model study of air-sea interactions in the Mediterranean Sea.},
+  journal = JMS,
+  year = {1998},
+  volume = {18},
+  pages = {89--114}
 }
 
@@ -1247 +1266 @@
 }
 
+@ARTICLE{Hellerman_Rosenstein_JPO83,
+  author = {S. Hellerman and  M. Rosenstein },
+  title = {Normal monthly wind stress over the world ocean with error estimates},
+  journal = JPO,
+  year = {1983},
+  volume = {13},
+  pages = {1093--1104},
+}
+
 @ARTICLE{He_Ding_JSC01,
   author = {Y. He and C. H. Q. Ding},
@@ -1511 +1539 @@
   volume = {6},
   pages = {56--58}
+}
+@ARTICLE{Kondo1975,
+  author = {J. Kondo},
+  title = {Air-sea bulk transfer coefficients in diabatic conditions},
+  journal = {Boundary-Layer Meteorol},
+  year = {1975},
+  volume = {9},
+  pages = {91--112}
+}
+
+@ARTICLE{Lermusiaux2001,
+  author = {P. F. J. Lermusiaux},
+  title = {Evolving the subspace of three-dimensional miltiscale ocean variability: Massachusetts Bay},
+  journal = JMS,
+  year = {2001},
+  volume = {29},
+  pages = {385--422}
 }
 
@@ -1906 +1951 @@
   volume = {125},  number = {5},
   pages = {958--971}
+}
+
+@ARTICLE{Maggiore_al_PCE98,
+  author = {A. Maggiore and  M. Zavatarelli and M. G. Angelucci and N. Pinardi},
+  title = {Surface heat and water fluxes in the Adriatic Sea: seasonal and interannual variability},
+  journal = {Phys Chem Earth},
+  year = {1998},
+  volume = {23},
+  pages = {561--567}
 }
 
@@ -2115 +2169 @@
 }
 
+@ARTICLE{Oddo_al_OS09,
+  author = {P. Oddo and M. Adani and N. Pinardi and C. Fratianni and M. Tonani and D. Pettenuzzo},
+  title = {A nested Atlantic-Mediterranean Sea general circulation model for operational forecasting},
+  journal = OS,
+  year = {2009},
+  volume = {5},
+  pages = {1--13},
+}
+
 @PHDTHESIS{Olivier_PhD01,
   author = {F. Olivier},
@@ -2170 +2233 @@
 }
 
+@ARTICLE{Payne_JAS72,
+  author = {R. E. Payne},
+  title = {Albedo of the Sea Surface},
+  journal = JAS,
+  year = {1972},
+  volume = {29}
+  pages = {959--970}
+}
+
 @ARTICLE{Penduff_al_OM06,
   author = {T. Penduff and B. Barnier and J.-M. Molines and G. Madec},
@@ -2234 +2306 @@
   volume = {13},
   pages = {1154--1158}
+}
+
+@ARTICLE{Reed_JPO77,
+  author = {R. K. Reed},
+  title = {On estimating insolation over the ocean},
+  journal = JPO,
+  year = {1977},
+  volume = {1},
+  pages = {874--971}
 }
 
@@ -2502 +2583 @@
   volume = {8},
   pages = {175--201}
+}
+
+@ARTICLE{Tonani_al_OS08,
+  author = {M. Tonani and N. Pinardi and S. Dobricic and I. Pujol and C. Fratianni},
+  title = {A high-resolution free-surface model of the Mediterranean Sea},
+  journal = OS,
+  year = {2008},
+  volume = {4},
+  pages = {1--14}
 }
 

branches/2011/dev_INGV_2011/DOC/TexFiles/Chapters/Chap_SBC.tex

@@ -24 +24 @@
 \end{itemize}
 
-Four different ways to provide the first six fields to the ocean are available which 
+Five different ways to provide the first six fields to the ocean are available which 
 are controlled by namelist variables: an analytical formulation (\np{ln\_ana}~=~true), 
 a flux formulation (\np{ln\_flx}~=~true), a bulk formulae formulation (CORE 
-(\np{ln\_core}~=~true) or CLIO (\np{ln\_clio}~=~true) bulk formulae) and a coupled 
+(\np{ln\_core}~=~true), CLIO (\np{ln\_clio}~=~true) or MFS
+\footnote { Note that MFS bulk formulae compute fluxes only for the ocean component}
+(\np{ln\_ecmwf}~=~true) bulk formulae) and a coupled 
 formulation (exchanges with a atmospheric model via the OASIS coupler) 
 (\np{ln\_cpl}~=~true). When used, the atmospheric pressure forces both 
-ocean and ice dynamics (\np{ln\_apr\_dyn}~=~true)
-\footnote{The surface pressure field could be use in bulk formulae, nevertheless 
-none of the current bulk formulea (CLIO and CORE) uses the it.}. 
+ocean and ice dynamics (\np{ln\_apr\_dyn}~=~true).
 The frequency at which the six or seven fields have to be updated is the \np{nn\_fsbc} 
 namelist parameter. 
@@ -51 +51 @@
 
 In this chapter, we first discuss where the surface boundary condition appears in the
-model equations. Then we present the four ways of providing the surface boundary condition, 
+model equations. Then we present the five ways of providing the surface boundary condition, 
 followed by the description of the atmospheric pressure and the river runoff. 
 Next the scheme for interpolation on the fly is described.
@@ -480 +480 @@
 % Bulk formulation
 % ================================================================
-\section  [Bulk formulation (\textit{sbcblk\_core} or \textit{sbcblk\_clio}) ]
-      {Bulk formulation \small{(\mdl{sbcblk\_core} or \mdl{sbcblk\_clio} module)} }
+\section  [Bulk formulation (\textit{sbcblk\_core}, \textit{sbcblk\_clio} or \textit{sbcblk\_ecmwf}) ]
+      {Bulk formulation \small{(\mdl{sbcblk\_core} \mdl{sbcblk\_clio} \mdl{sbcblk\_ecmwf} modules)} }
 \label{SBC_blk}
 
@@ -487 +487 @@
 using bulk formulae and atmospheric fields and ocean (and ice) variables. 
 
-The atmospheric fields used depend on the bulk formulae used. Two bulk formulations 
-are available : the CORE and CLIO bulk formulea. The choice is made by setting to true
-one of the following namelist variable : \np{ln\_core} and \np{ln\_clio}.
-
-Note : in forced mode, when a sea-ice model is used, a bulk formulation have to be used. 
-Therefore the two bulk formulea provided include the computation of the fluxes over both 
+The atmospheric fields used depend on the bulk formulae used. Three bulk formulations 
+are available : the CORE, the CLIO and the MFS bulk formulea. The choice is made by setting to true
+one of the following namelist variable : \np{ln\_core} ; \np{ln\_clio} or  \np{ln\_ecmwf}.
+
+Note : in forced mode, when a sea-ice model is used, a bulk formulation (CLIO or CORE) have to be used. 
+Therefore the two bulk (CLIO and CORE) formulea include the computation of the fluxes over both 
 an ocean and an ice surface. 
 
@@ -583 +583 @@
 namelist (see \S\ref{SBC_fldread}). 
 
+% -------------------------------------------------------------------------------------------------------------
+%        ECMWF Bulk formulea
+% -------------------------------------------------------------------------------------------------------------
+\subsection    [MFS Bulk formulea (\np{ln\_ecmwf}=true)]
+            {MFS Bulk formulea (\np{ln\_ecmwf}=true, \mdl{sbcblk\_ecmwf})}
+\label{SBC_blk_ecmwf}
+%------------------------------------------namsbc_ecmwf----------------------------------------------------
+\namdisplay{namsbc_ecmwf} 
+%----------------------------------------------------------------------------------------------------------
+
+The MFS (Mediterranean Forecasting System) bulk formulae have been developed by
+ \citet{Castellari_al_JMS1998}. 
+They have been designed to handle the ECMWF operational data and are currently 
+in use in the MFS operational system \citep{Tonani_al_OS08}, \citep{Oddo_al_OS09}.
+The wind stress computation uses a drag coefficient computed according to \citet{Hellerman_Rosenstein_JPO83}.
+The surface boundary condition for temperature involves the balance between surface solar radiation,
+net long-wave radiation, the latent and sensible heat fluxes.
+Solar radiation is dependent on cloud cover and is computed by means of
+an astronomical formula \citep{Reed_JPO77}. Albedo monthly values are from \citet{Payne_JAS72} 
+as means of the values at $40^{o}N$ and $30^{o}N$ for the Atlantic Ocean (hence the same latitudinal
+band of the Mediterranean Sea). The net long-wave radiation flux
+\citep{Bignami_al_JGR95} is a function of
+air temperature, sea-surface temperature, cloud cover and relative humidity.
+Sensible heat and latent heat fluxes are computed by classical
+bulk formulae parameterized according to \citet{Kondo1975}.
+Details on the bulk formulae used can be found in \citet{Maggiore_al_PCE98} and \citet{Castellari_al_JMS1998}.
+
+The required 7 input fields must be provided on the model Grid-T and  are:
+\begin{itemize}
+\item          Zonal Component of the 10m wind ($ms^{-1}$)  (\np{sn\_windi})
+\item          Meridional Component of the 10m wind ($ms^{-1}$)  (\np{sn\_windj})
+\item          Total Claud Cover (\%)  (\np{sn\_clc})
+\item          2m Air Temperature ($K$) (\np{sn\_tair})
+\item          2m Dew Point Temperature ($K$)  (\np{sn\_rhm})
+\item          Total Precipitation ${Kg} m^{-2} s^{-1}$ (\np{sn\_prec})
+\item          Mean Sea Level Pressure (${Pa}) (\np{sn\_msl})
+\end{itemize}
+% -------------------------------------------------------------------------------------------------------------
 % ================================================================
 % Coupled formulation

branches/2011/dev_INGV_2011/DOC/TexFiles/Chapters/Chap_ZDF.tex

@@ -100 +100 @@
 $a=5$ and $n=2$. The last three values can be modified by setting the 
 \np{rn\_avmri}, \np{rn\_alp} and \np{nn\_ric} namelist parameters, respectively.
+
+A simple mixing-layer model to transfer and dissipate the atmospheric
+ forcings (wind-stress and buoyancy fluxes) can be activated setting 
+the \np{ln\_mldw} =.true. in the namelist.
+
+In this case, the local depth of turbulent wind-mixing or "Ekman depth"
+ $h_{e}(x,y,t)$ is evaluated and the vertical eddy coefficients prescribed within this layer.
+
+This depth is assumed proportional to the "depth of frictional influence" that is limited by rotation:
+\begin{equation}
+         h_{e} = Ek \frac {u^{*}} {f_{0}}    \\
+\end{equation}
+where, $Ek$ is an empirical parameter, $u^{*}$ is the friction velocity and $f_{0}$ is the Coriolis 
+parameter.
+
+In this similarity height relationship, the turbulent friction velocity:
+\begin{equation}
+         u^{*} = \sqrt \frac {|\tau|} {\rho_o}     \\
+\end{equation}
+
+is computed from the wind stress vector $|\tau|$ and the reference dendity $ \rho_o$.
+The final $h_{e}$ is further constrained by the adjustable bounds \np{rn\_mldmin} and \np{rn\_mldmax}.
+Once $h_{e}$ is computed, the vertical eddy coefficients within $h_{e}$ are set to 
+the empirical values \np{rn\_wtmix} and \np{rn\_wvmix} \citep{Lermusiaux2001}.
 
 % -------------------------------------------------------------------------------------------------------------

branches/2011/dev_INGV_2011/DOC/TexFiles/Namelist/namsbc

@@ -8 +8 @@
    ln_blk_clio = .false.   !  CLIO bulk formulation                     (T => fill namsbc_clio) 
    ln_blk_core = .true.    !  CORE bulk formulation                     (T => fill namsbc_core) 
+   ln_blk_ecmwf= .false.   !  MFS bulk formulation                      (T => fill namsbc_ecmwf)
    ln_cpl      = .false.   !  Coupled formulation                       (T => fill namsbc_cpl )
    ln_apr_dyn  = .false.   !  Patm gradient added in ocean & ice Eqs.   (T => fill namsbc_apr )

branches/2011/dev_INGV_2011/DOC/TexFiles/Namelist/namzdf_ric

@@ -5 +5 @@
    rn_alp      =   5.      !  coefficient of the parameterization
    nn_ric      =   2       !  coefficient of the parameterization
+   rn_ekmfc    =   0.7     !  Factor in the Ekman depth Equation
+   rn_mldmin   =   1.0     !  minimum allowable mixed-layer depth estimate (m)
+   rn_mldmax   =1000.0     !  maximum allowable mixed-layer depth estimate (m)
+   rn_wtmix    =  10.0     !  vertical eddy viscosity coeff [m2/s] in the mixed-layer
+   rn_wvmix    =  10.0     !  vertical eddy diffusion coeff [m2/s] in the mixed-layer
+   ln_mldw     = .true.    !  Flag to use or not the mized layer depth param.
 /

branches/2011/dev_INGV_2011/NEMOGCM/CONFIG/GYRE/EXP00/namelist

@@ -592 +592 @@
    rn_alp      =   5.      !  coefficient of the parameterization
    nn_ric      =   2       !  coefficient of the parameterization
+   rn_ekmfc    =   0.7     !  Factor in the Ekman depth Equation
+   rn_mldmin   =   1.0     !  minimum allowable mixed-layer depth estimate (m)
+   rn_mldmax   =1000.0     !  maximum allowable mixed-layer depth estimate (m)
+   rn_wtmix    =  10.0     !  vertical eddy viscosity coeff [m2/s] in the mixed-layer
+   rn_wvmix    =  10.0     !  vertical eddy diffusion coeff [m2/s] in the mixed-layer
+   ln_mldw     = .true.    !  Flag to use or not the mized layer depth param.
 /
 !-----------------------------------------------------------------------

branches/2011/dev_INGV_2011/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist

@@ -592 +592 @@
    rn_alp      =   5.      !  coefficient of the parameterization
    nn_ric      =   2       !  coefficient of the parameterization
+   rn_ekmfc    =   0.7     !  Factor in the Ekman depth Equation
+   rn_mldmin   =   1.0     !  minimum allowable mixed-layer depth estimate (m)
+   rn_mldmax   =1000.0     !  maximum allowable mixed-layer depth estimate (m)
+   rn_wtmix    =  10.0     !  vertical eddy viscosity coeff [m2/s] in the mixed-layer
+   rn_wvmix    =  10.0     !  vertical eddy diffusion coeff [m2/s] in the mixed-layer
+   ln_mldw     = .true.    !  Flag to use or not the mized layer depth param.
 /
 !-----------------------------------------------------------------------

branches/2011/dev_INGV_2011/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist

@@ -606 +606 @@
    rn_alp      =   5.      !  coefficient of the parameterization
    nn_ric      =   2       !  coefficient of the parameterization
+   rn_ekmfc    =   0.7     !  Factor in the Ekman depth Equation
+   rn_mldmin   =   1.0     !  minimum allowable mixed-layer depth estimate (m)
+   rn_mldmax   =1000.0     !  maximum allowable mixed-layer depth estimate (m)
+   rn_wtmix    =  10.0     !  vertical eddy viscosity coeff [m2/s] in the mixed-layer
+   rn_wvmix    =  10.0     !  vertical eddy diffusion coeff [m2/s] in the mixed-layer
+   ln_mldw     = .true.    !  Flag to use or not the mized layer depth param.
 /
 !-----------------------------------------------------------------------

branches/2011/dev_INGV_2011/NEMOGCM/CONFIG/POMME/EXP00/namelist

@@ -592 +592 @@
    rn_alp      =   5.      !  coefficient of the parameterization
    nn_ric      =   2       !  coefficient of the parameterization
+   rn_ekmfc    =   0.7     !  Factor in the Ekman depth Equation
+   rn_mldmin   =   1.0     !  minimum allowable mixed-layer depth estimate (m)
+   rn_mldmax   =1000.0     !  maximum allowable mixed-layer depth estimate (m)
+   rn_wtmix    =  10.0     !  vertical eddy viscosity coeff [m2/s] in the mixed-layer
+   rn_wvmix    =  10.0     !  vertical eddy diffusion coeff [m2/s] in the mixed-layer
+   ln_mldw     = .true.    !  Flag to use or not the mized layer depth param.
 /
 !-----------------------------------------------------------------------

branches/2011/dev_INGV_2011/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfric.F90

@@ -7 +7 @@
    !! History :  OPA  ! 1987-09  (P. Andrich)  Original code
    !!            4.0  ! 1991-11  (G. Madec)
-   !!            7.0  ! 1996-01  (G. Madec)  complet rewriting of multitasking suppression of common work arrays
-   !!            8.0  ! 1997-06 (G. Madec)  complete rewriting of zdfmix
+   !!            7.0  ! 1996-01  (G. Madec)  complete rewriting of multitasking suppression of common work arrays
+   !!            8.0  ! 1997-06  (G. Madec)  complete rewriting of zdfmix
    !!   NEMO     1.0  ! 2002-06  (G. Madec)  F90: Free form and module
    !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
+   !!            3.3.1! 2011-09  (P. Oddo) Mixed layer depth parameterization
    !!----------------------------------------------------------------------
 #if defined key_zdfric   ||   defined key_esopa
@@ -20 +21 @@
    !!   zdf_ric_init : initialization, namelist read, & parameters control
    !!----------------------------------------------------------------------
-   USE oce             ! ocean dynamics and tracers variables
-   USE dom_oce         ! ocean space and time domain variables
-   USE zdf_oce         ! ocean vertical physics
-   USE in_out_manager  ! I/O manager
-   USE lbclnk          ! ocean lateral boundary condition (or mpp link)
-   USE lib_mpp         ! MPP library
+   USE oce                   ! ocean dynamics and tracers variables
+   USE dom_oce               ! ocean space and time domain variables
+   USE zdf_oce               ! ocean vertical physics
+   USE in_out_manager        ! I/O manager
+   USE lbclnk                ! ocean lateral boundary condition (or mpp link)
+   USE lib_mpp               ! MPP library
+   USE eosbn2, ONLY : nn_eos
 
    IMPLICIT NONE
@@ -39 +41 @@
    REAL(wp) ::   rn_avmri = 100.e-4_wp   ! maximum value of the vertical eddy viscosity
    REAL(wp) ::   rn_alp   =   5._wp      ! coefficient of the parameterization
+   REAL(wp) ::   rn_ekmfc =   0.7_wp     ! Ekman Factor Coeff
+   REAL(wp) ::   rn_mldmin=   1.0_wp     ! minimum mixed layer (ML) depth    
+   REAL(wp) ::   rn_mldmax=1000.0_wp     ! maximum mixed layer depth
+   REAL(wp) ::   rn_wtmix =  10.0_wp     ! Vertical eddy Diff. in the ML
+   REAL(wp) ::   rn_wvmix =  10.0_wp     ! Vertical eddy Visc. in the ML
+   LOGICAL  ::   ln_mldw  = .TRUE.       ! Use or not the MLD parameters
 
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tmric   !: coef. for the horizontal mean at t-point
@@ -67 +75 @@
       !!                    
       !! ** Purpose :   Compute the before eddy viscosity and diffusivity as
-      !!              a function of the local richardson number.
+      !!                a function of the local richardson number.
       !!
       !! ** Method  :   Local richardson number dependent formulation of the 
-      !!              vertical eddy viscosity and diffusivity coefficients. 
+      !!                vertical eddy viscosity and diffusivity coefficients. 
       !!                The eddy coefficients are given by:
       !!                    avm = avm0 + avmb
       !!                    avt = avm0 / (1 + rn_alp*ri)
-      !!              with ri  = N^2 / dz(u)**2
-      !!                       = e3w**2 * rn2/[ mi( dk(ub) )+mj( dk(vb) ) ]
-      !!                   avm0= rn_avmri / (1 + rn_alp*ri)**nn_ric
+      !!                with ri  = N^2 / dz(u)**2
+      !!                         = e3w**2 * rn2/[ mi( dk(ub) )+mj( dk(vb) ) ]
+      !!                    avm0= rn_avmri / (1 + rn_alp*ri)**nn_ric
       !!      Where ri is the before local Richardson number,
       !!            rn_avmri is the maximum value reaches by avm and avt 
@@ -84 +92 @@
       !!      avtb=1.e-7 m2/s, rn_alp=5. and nn_ric=2.
       !!      a numerical threshold is impose on the vertical shear (1.e-20)
+      !!      As second step compute Ekman depth from wind stress forcing
+      !!      and apply namelist provided vertical coeff within this depth.
+      !!      The Ekman depth is:
+      !!              Ustar = SQRT(Taum/rho0)
+      !!              ekd= rn_ekmfc * Ustar / f0
+      !!      Large et al. (1994, eq.29) suggest rn_ekmfc=0.7; however, the derivation
+      !!      of the above equation indicates the value is somewhat arbitrary; therefore
+      !!      we allow the freedom to increase or decrease this value, if the
+      !!      Ekman depth estimate appears too shallow or too deep, respectively.
+      !!      Ekd is then limited by rn_mldmin and rn_mldmax provided in the
+      !!      namelist
       !!        N.B. the mask are required for implicit scheme, and surface
       !!      and bottom value already set in zdfini.F90
       !!
       !! References : Pacanowski & Philander 1981, JPO, 1441-1451.
+      !!              PFJ Lermusiaux 2001.
       !!----------------------------------------------------------------------
       USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
-      USE wrk_nemo, ONLY:   zwx => wrk_2d_1     ! 2D workspace
-      !!
-      INTEGER, INTENT( in ) ::   kt         ! ocean time-step indexocean time step
-      !!
-      INTEGER  ::   ji, jj, jk               ! dummy loop indices
-      REAL(wp) ::   zcoef, zdku, zdkv, zri, z05alp     ! temporary scalars
-      !!----------------------------------------------------------------------
-
-      IF( wrk_in_use(2, 1) ) THEN
+      USE wrk_nemo, ONLY:   zwx => wrk_2d_1         ! 2D workspace
+      USE wrk_nemo, ONLY:   ekm_dep => wrk_2d_2     ! 2D workspace
+      USE phycst,   ONLY:   rsmall,rau0
+      USE sbc_oce,  ONLY:   taum
+      !!
+      INTEGER, INTENT( in ) ::   kt                           ! ocean time-step
+      !!
+      INTEGER  ::   ji, jj, jk                                ! dummy loop indices
+      REAL(wp) ::   zcoef, zdku, zdkv, zri, z05alp, zflageos  ! temporary scalars
+      REAL(wp) ::   zrhos, zustar
+      !!----------------------------------------------------------------------
+
+      IF( wrk_in_use(2, 1,2) ) THEN
          CALL ctl_stop('zdf_ric : requested workspace array unavailable')   ;   RETURN
       ENDIF
@@ -145 +169 @@
       !                                                ! ===============
       !
+      IF( ln_mldw ) THEN
+
+      !  Compute Ekman depth from wind stress forcing.
+      ! -------------------------------------------------------
+      zflageos = ( 0.5 + SIGN( 0.5, nn_eos - 1. ) ) * rau0
+      DO jj = 1, jpj
+         DO ji = 1, jpi
+            zrhos          = rhop(ji,jj,1) + zflageos * ( 1. - tmask(ji,jj,1) )
+            zustar         = SQRT( taum(ji,jj) / ( zrhos +  rsmall ) )
+            ekm_dep(ji,jj) = rn_ekmfc * zustar / ( ABS( ff(ji,jj) ) + rsmall )
+            ekm_dep(ji,jj) = MAX(ekm_dep(ji,jj),rn_mldmin) ! Minimun allowed
+            ekm_dep(ji,jj) = MIN(ekm_dep(ji,jj),rn_mldmax) ! Maximum allowed
+         END DO
+      END DO
+
+      ! In the first model level vertical diff/visc coeff.s 
+      ! are always equal to the namelist values rn_wtmix/rn_wvmix
+      ! -------------------------------------------------------
+      DO jj = 1, jpj
+         DO ji = 1, jpi
+            avmv(ji,jj,1) = MAX( avmv(ji,jj,1), rn_wvmix )
+            avmu(ji,jj,1) = MAX( avmu(ji,jj,1), rn_wvmix )
+            avt( ji,jj,1) = MAX(  avt(ji,jj,1), rn_wtmix )
+         END DO
+      END DO
+
+      !  Force the vertical mixing coef within the Ekman depth
+      ! -------------------------------------------------------
+      DO jk = 2, jpkm1
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               IF( fsdept(ji,jj,jk) < ekm_dep(ji,jj) ) THEN
+                  avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk), rn_wvmix )
+                  avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk), rn_wvmix )
+                  avt( ji,jj,jk) = MAX(  avt(ji,jj,jk), rn_wtmix )
+               ENDIF
+            END DO
+         END DO
+      END DO
+
+      DO jk = 1, jpkm1                
+         DO jj = 1, jpj
+            DO ji = 1, jpi
+               avmv(ji,jj,jk) = avmv(ji,jj,jk) * vmask(ji,jj,jk)
+               avmu(ji,jj,jk) = avmu(ji,jj,jk) * umask(ji,jj,jk)
+               avt( ji,jj,jk) = avt( ji,jj,jk) * tmask(ji,jj,jk)
+            END DO
+         END DO
+      END DO
+
+     ENDIF
+
       CALL lbc_lnk( avt , 'W', 1. )                         ! Boundary conditions   (unchanged sign)
       CALL lbc_lnk( avmu, 'U', 1. )   ;   CALL lbc_lnk( avmv, 'V', 1. )
       !
-      IF( wrk_not_released(2, 1) )   CALL ctl_stop('zdf_ric: failed to release workspace array')
+      IF( wrk_not_released(2, 1,2) )   CALL ctl_stop('zdf_ric: failed to release workspace array')
       !
    END SUBROUTINE zdf_ric
@@ -168 +244 @@
       INTEGER :: ji, jj, jk   ! dummy loop indices
       !!
-      NAMELIST/namzdf_ric/ rn_avmri, rn_alp, nn_ric
+      NAMELIST/namzdf_ric/ rn_avmri, rn_alp   , nn_ric  , rn_ekmfc,  &
+         &                rn_mldmin, rn_mldmax, rn_wtmix, rn_wvmix, ln_mldw
       !!----------------------------------------------------------------------
       !
@@ -179 +256 @@
          WRITE(numout,*) '~~~~~~~'
          WRITE(numout,*) '   Namelist namzdf_ric : set Kz(Ri) parameters'
-         WRITE(numout,*) '      maximum vertical viscosity     rn_avmri = ', rn_avmri
-         WRITE(numout,*) '      coefficient                    rn_alp   = ', rn_alp
-         WRITE(numout,*) '      coefficient                    nn_ric   = ', nn_ric
+         WRITE(numout,*) '      maximum vertical viscosity     rn_avmri  = ', rn_avmri
+         WRITE(numout,*) '      coefficient                    rn_alp    = ', rn_alp
+         WRITE(numout,*) '      coefficient                    nn_ric    = ', nn_ric
+         WRITE(numout,*) '      Ekman Factor Coeff             rn_ekmfc  = ', rn_ekmfc
+         WRITE(numout,*) '      minimum mixed layer depth      rn_mldmin = ', rn_mldmin
+         WRITE(numout,*) '      maximum mixed layer depth      rn_mldmax = ', rn_mldmax
+         WRITE(numout,*) '      Vertical eddy Diff. in the ML  rn_wtmix  = ', rn_wtmix
+         WRITE(numout,*) '      Vertical eddy Visc. in the ML  rn_wvmix  = ', rn_wvmix
+         WRITE(numout,*) '      Use the MLD parameterization   ln_mldw   = ', ln_mldw
       ENDIF
       !