Changeset 11123 for NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex

Timestamp:

2019-06-17T14:22:27+02:00 (5 years ago)

Author:

nicolasmartin

Message:

Modification of LaTeX subfiles accordingly to new citations keys

Location:

NEMO/trunk/doc/latex/NEMO/subfiles

Files:

• . (modified) (1 prop)
• chap_ZDF.tex (modified) (46 diffs)

NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex

@@ -87 +87 @@
 a dependency between the vertical eddy coefficients and the local Richardson number
 (\ie the ratio of stratification to vertical shear).
-Following \citet{Pacanowski_Philander_JPO81}, the following formulation has been implemented:
+Following \citet{pacanowski.philander_JPO81}, the following formulation has been implemented:
 \[
   % \label{eq:zdfric}
@@ -124 +124 @@
 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}.
+the empirical values \np{rn\_wtmix} and \np{rn\_wvmix} \citep{lermusiaux_JMS01}.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -140 +140 @@
 a prognostic equation for $\bar{e}$, the turbulent kinetic energy,
 and a closure assumption for the turbulent length scales.
-This turbulent closure model has been developed by \citet{Bougeault1989} in the atmospheric case,
-adapted by \citet{Gaspar1990} for the oceanic case, and embedded in OPA, the ancestor of NEMO,
-by \citet{Blanke1993} for equatorial Atlantic simulations.
-Since then, significant modifications have been introduced by \citet{Madec1998} in both the implementation and
+This turbulent closure model has been developed by \citet{bougeault.lacarrere_MWR89} in the atmospheric case,
+adapted by \citet{gaspar.gregoris.ea_JGR90} for the oceanic case, and embedded in OPA, the ancestor of NEMO,
+by \citet{blanke.delecluse_JPO93} for equatorial Atlantic simulations.
+Since then, significant modifications have been introduced by \citet{madec.delecluse.ea_NPM98} in both the implementation and
 the formulation of the mixing length scale.
 The time evolution of $\bar{e}$ is the result of the production of $\bar{e}$ through vertical shear,
-its destruction through stratification, its vertical diffusion, and its dissipation of \citet{Kolmogorov1942} type:
+its destruction through stratification, its vertical diffusion, and its dissipation of \citet{kolmogorov_IANS42} type:
 \begin{equation}
   \label{eq:zdftke_e}
@@ -168 +168 @@
 $P_{rt}$ is the Prandtl number, $K_m$ and $K_\rho$ are the vertical eddy viscosity and diffusivity coefficients.
 The constants $C_k =  0.1$ and $C_\epsilon = \sqrt {2} /2$ $\approx 0.7$ are designed to deal with
-vertical mixing at any depth \citep{Gaspar1990}. 
+vertical mixing at any depth \citep{gaspar.gregoris.ea_JGR90}. 
 They are set through namelist parameters \np{nn\_ediff} and \np{nn\_ediss}.
-$P_{rt}$ can be set to unity or, following \citet{Blanke1993}, be a function of the local Richardson number, $R_i$:
+$P_{rt}$ can be set to unity or, following \citet{blanke.delecluse_JPO93}, be a function of the local Richardson number, $R_i$:
 \begin{align*}
   % \label{eq:prt}
@@ -185 +185 @@
 At the sea surface, the value of $\bar{e}$ is prescribed from the wind stress field as
 $\bar{e}_o = e_{bb} |\tau| / \rho_o$, with $e_{bb}$ the \np{rn\_ebb} namelist parameter.
-The default value of $e_{bb}$ is 3.75. \citep{Gaspar1990}), however a much larger value can be used when
+The default value of $e_{bb}$ is 3.75. \citep{gaspar.gregoris.ea_JGR90}), however a much larger value can be used when
 taking into account the surface wave breaking (see below Eq. \autoref{eq:ZDF_Esbc}).
 The bottom value of TKE is assumed to be equal to the value of the level just above.
@@ -191 +191 @@
 the numerical scheme does not ensure its positivity.
 To overcome this problem, a cut-off in the minimum value of $\bar{e}$ is used (\np{rn\_emin} namelist parameter).
-Following \citet{Gaspar1990}, the cut-off value is set to $\sqrt{2}/2~10^{-6}~m^2.s^{-2}$.
-This allows the subsequent formulations to match that of \citet{Gargett1984} for the diffusion in
+Following \citet{gaspar.gregoris.ea_JGR90}, the cut-off value is set to $\sqrt{2}/2~10^{-6}~m^2.s^{-2}$.
+This allows the subsequent formulations to match that of \citet{gargett_JMR84} for the diffusion in
 the thermocline and deep ocean :  $K_\rho = 10^{-3} / N$.
 In addition, a cut-off is applied on $K_m$ and $K_\rho$ to avoid numerical instabilities associated with
@@ -202 +202 @@
 
 For computational efficiency, the original formulation of the turbulent length scales proposed by
-\citet{Gaspar1990} has been simplified.
+\citet{gaspar.gregoris.ea_JGR90} has been simplified.
 Four formulations are proposed, the choice of which is controlled by the \np{nn\_mxl} namelist parameter.
-The first two are based on the following first order approximation \citep{Blanke1993}:
+The first two are based on the following first order approximation \citep{blanke.delecluse_JPO93}:
 \begin{equation}
   \label{eq:tke_mxl0_1}
@@ -212 +212 @@
 The resulting length scale is bounded by the distance to the surface or to the bottom
 (\np{nn\_mxl}\forcode{ = 0}) or by the local vertical scale factor (\np{nn\_mxl}\forcode{ = 1}).
-\citet{Blanke1993} notice that this simplification has two major drawbacks:
+\citet{blanke.delecluse_JPO93} notice that this simplification has two major drawbacks:
 it makes no sense for locally unstable stratification and the computation no longer uses all
 the information contained in the vertical density profile.
-To overcome these drawbacks, \citet{Madec1998} introduces the \np{nn\_mxl}\forcode{ = 2..3} cases,
+To overcome these drawbacks, \citet{madec.delecluse.ea_NPM98} introduces the \np{nn\_mxl}\forcode{ = 2..3} cases,
 which add an extra assumption concerning the vertical gradient of the computed length scale.
 So, the length scales are first evaluated as in \autoref{eq:tke_mxl0_1} and then bounded such that:
@@ -225 +225 @@
 \autoref{eq:tke_mxl_constraint} means that the vertical variations of the length scale cannot be larger than
 the variations of depth.
-It provides a better approximation of the \citet{Gaspar1990} formulation while being much less 
+It provides a better approximation of the \citet{gaspar.gregoris.ea_JGR90} formulation while being much less 
 time consuming.
 In particular, it allows the length scale to be limited not only by the distance to the surface or
@@ -258 +258 @@
 In the \np{nn\_mxl}\forcode{ = 2} case, the dissipation and mixing length scales take the same value:
 $ l_k=  l_\epsilon = \min \left(\ l_{up} \;,\;  l_{dwn}\ \right)$, while in the \np{nn\_mxl}\forcode{ = 3} case,
-the dissipation and mixing turbulent length scales are give as in \citet{Gaspar1990}:
+the dissipation and mixing turbulent length scales are give as in \citet{gaspar.gregoris.ea_JGR90}:
 \[
   % \label{eq:tke_mxl_gaspar}
@@ -270 +270 @@
 Usually the surface scale is given by $l_o = \kappa \,z_o$ where $\kappa = 0.4$ is von Karman's constant and
 $z_o$ the roughness parameter of the surface.
-Assuming $z_o=0.1$~m \citep{Craig_Banner_JPO94} leads to a 0.04~m, the default value of \np{rn\_mxl0}.
+Assuming $z_o=0.1$~m \citep{craig.banner_JPO94} leads to a 0.04~m, the default value of \np{rn\_mxl0}.
 In the ocean interior a minimum length scale is set to recover the molecular viscosity when
 $\bar{e}$ reach its minimum value ($1.10^{-6}= C_k\, l_{min} \,\sqrt{\bar{e}_{min}}$ ).
@@ -277 +277 @@
 %-----------------------------------------------------------------------%
 
-Following \citet{Mellor_Blumberg_JPO04}, the TKE turbulence closure model has been modified to
+Following \citet{mellor.blumberg_JPO04}, the TKE turbulence closure model has been modified to
 include the effect of surface wave breaking energetics.
 This results in a reduction of summertime surface temperature when the mixed layer is relatively shallow.
-The \citet{Mellor_Blumberg_JPO04} modifications acts on surface length scale and TKE values and
+The \citet{mellor.blumberg_JPO04} modifications acts on surface length scale and TKE values and
 air-sea drag coefficient. 
 The latter concerns the bulk formulea and is not discussed here. 
 
-Following \citet{Craig_Banner_JPO94}, the boundary condition on surface TKE value is :
+Following \citet{craig.banner_JPO94}, the boundary condition on surface TKE value is :
 \begin{equation}
   \label{eq:ZDF_Esbc}
   \bar{e}_o = \frac{1}{2}\,\left(  15.8\,\alpha_{CB} \right)^{2/3} \,\frac{|\tau|}{\rho_o}
 \end{equation}
-where $\alpha_{CB}$ is the \citet{Craig_Banner_JPO94} constant of proportionality which depends on the ''wave age'',
-ranging from 57 for mature waves to 146 for younger waves \citep{Mellor_Blumberg_JPO04}. 
+where $\alpha_{CB}$ is the \citet{craig.banner_JPO94} constant of proportionality which depends on the ''wave age'',
+ranging from 57 for mature waves to 146 for younger waves \citep{mellor.blumberg_JPO04}. 
 The boundary condition on the turbulent length scale follows the Charnock's relation:
 \begin{equation}
@@ -297 +297 @@
 \end{equation}
 where $\kappa=0.40$ is the von Karman constant, and $\beta$ is the Charnock's constant.
-\citet{Mellor_Blumberg_JPO04} suggest $\beta = 2.10^{5}$ the value chosen by
-\citet{Stacey_JPO99} citing observation evidence, and
+\citet{mellor.blumberg_JPO04} suggest $\beta = 2.10^{5}$ the value chosen by
+\citet{stacey_JPO99} citing observation evidence, and
 $\alpha_{CB} = 100$ the Craig and Banner's value.
 As the surface boundary condition on TKE is prescribed through $\bar{e}_o = e_{bb} |\tau| / \rho_o$,
@@ -315 +315 @@
 Although LC have nothing to do with convection, the circulation pattern is rather similar to
 so-called convective rolls in the atmospheric boundary layer.
-The detailed physics behind LC is described in, for example, \citet{Craik_Leibovich_JFM76}.
+The detailed physics behind LC is described in, for example, \citet{craik.leibovich_JFM76}.
 The prevailing explanation is that LC arise from a nonlinear interaction between the Stokes drift and
 wind drift currents. 
 
 Here we introduced in the TKE turbulent closure the simple parameterization of Langmuir circulations proposed by
-\citep{Axell_JGR02} for a $k-\epsilon$ turbulent closure.
+\citep{axell_JGR02} for a $k-\epsilon$ turbulent closure.
 The parameterization, tuned against large-eddy simulation, includes the whole effect of LC in
 an extra source terms of TKE, $P_{LC}$.
@@ -326 +326 @@
 \forcode{.true.} in the namtke namelist.
  
-By making an analogy with the characteristic convective velocity scale (\eg, \citet{D'Alessio_al_JPO98}),
+By making an analogy with the characteristic convective velocity scale (\eg, \citet{dalessio.abdella.ea_JPO98}),
 $P_{LC}$ is assumed to be : 
 \[
@@ -334 +334 @@
 With no information about the wave field, $w_{LC}$ is assumed to be proportional to 
 the Stokes drift $u_s = 0.377\,\,|\tau|^{1/2}$, where $|\tau|$ is the surface wind stress module 
-\footnote{Following \citet{Li_Garrett_JMR93}, the surface Stoke drift velocity may be expressed as
+\footnote{Following \citet{li.garrett_JMR93}, the surface Stoke drift velocity may be expressed as
   $u_s =  0.016 \,|U_{10m}|$.
   Assuming an air density of $\rho_a=1.22 \,Kg/m^3$ and a drag coefficient of
@@ -350 +350 @@
   \end{cases}
 \]
-where $c_{LC} = 0.15$ has been chosen by \citep{Axell_JGR02} as a good compromise to fit LES data.
+where $c_{LC} = 0.15$ has been chosen by \citep{axell_JGR02} as a good compromise to fit LES data.
 The chosen value yields maximum vertical velocities $w_{LC}$ of the order of a few centimeters per second.
 The value of $c_{LC}$ is set through the \np{rn\_lc} namelist parameter,
-having in mind that it should stay between 0.15 and 0.54 \citep{Axell_JGR02}. 
+having in mind that it should stay between 0.15 and 0.54 \citep{axell_JGR02}. 
 
 The $H_{LC}$ is estimated in a similar way as the turbulent length scale of TKE equations:
@@ -368 +368 @@
 produce mixed-layer depths that are too shallow during summer months and windy conditions.
 This bias is particularly acute over the Southern Ocean.
-To overcome this systematic bias, an ad hoc parameterization is introduced into the TKE scheme \cite{Rodgers_2014}. 
+To overcome this systematic bias, an ad hoc parameterization is introduced into the TKE scheme \cite{rodgers.aumont.ea_B14}. 
 The parameterization is an empirical one, \ie not derived from theoretical considerations,
 but rather is meant to account for observed processes that affect the density structure of 
@@ -427 +427 @@
 (first line in \autoref{eq:PE_zdf}).
 To do so a special care have to be taken for both the time and space discretization of
-the TKE equation \citep{Burchard_OM02,Marsaleix_al_OM08}.
+the TKE equation \citep{burchard_OM02,marsaleix.auclair.ea_OM08}.
 
 Let us first address the time stepping issue. \autoref{fig:TKE_time_scheme} shows how
@@ -524 +524 @@
 The Generic Length Scale (GLS) scheme is a turbulent closure scheme based on two prognostic equations:
 one for the turbulent kinetic energy $\bar {e}$, and another for the generic length scale,
-$\psi$ \citep{Umlauf_Burchard_JMS03, Umlauf_Burchard_CSR05}.
+$\psi$ \citep{umlauf.burchard_JMR03, umlauf.burchard_CSR05}.
 This later variable is defined as: $\psi = {C_{0\mu}}^{p} \ {\bar{e}}^{m} \ l^{n}$, 
 where the triplet $(p, m, n)$ value given in Tab.\autoref{tab:GLS} allows to recover a number of
-well-known turbulent closures ($k$-$kl$ \citep{Mellor_Yamada_1982}, $k$-$\epsilon$ \citep{Rodi_1987},
-$k$-$\omega$ \citep{Wilcox_1988} among others \citep{Umlauf_Burchard_JMS03,Kantha_Carniel_CSR05}). 
+well-known turbulent closures ($k$-$kl$ \citep{mellor.yamada_RG82}, $k$-$\epsilon$ \citep{rodi_JGR87},
+$k$-$\omega$ \citep{wilcox_AJ88} among others \citep{umlauf.burchard_JMR03,kantha.carniel_JMR03}). 
 The GLS scheme is given by the following set of equations:
 \begin{equation}
@@ -577 +577 @@
     \begin{tabular}{ccccc}
       &   $k-kl$   & $k-\epsilon$ & $k-\omega$ &   generic   \\
-      % & \citep{Mellor_Yamada_1982} &  \citep{Rodi_1987}       & \citep{Wilcox_1988} &                 \\
+      % & \citep{mellor.yamada_RG82} &  \citep{rodi_JGR87}       & \citep{wilcox_AJ88} &                 \\
       \hline
       \hline
@@ -604 +604 @@
 the mixing length towards $K z_b$ ($K$: Kappa and $z_b$: rugosity length) value near physical boundaries
 (logarithmic boundary layer law).
-$C_{\mu}$ and $C_{\mu'}$ are calculated from stability function proposed by \citet{Galperin_al_JAS88},
-or by \citet{Kantha_Clayson_1994} or one of the two functions suggested by \citet{Canuto_2001}
+$C_{\mu}$ and $C_{\mu'}$ are calculated from stability function proposed by \citet{galperin.kantha.ea_JAS88},
+or by \citet{kantha.clayson_JGR94} or one of the two functions suggested by \citet{canuto.howard.ea_JPO01}
 (\np{nn\_stab\_func}\forcode{ = 0..3}, resp.). 
 The value of $C_{0\mu}$ depends of the choice of the stability function.
@@ -612 +612 @@
 Neumann condition through \np{nn\_tkebc\_surf} and \np{nn\_tkebc\_bot}, resp.
 As for TKE closure, the wave effect on the mixing is considered when
-\np{ln\_crban}\forcode{ = .true.} \citep{Craig_Banner_JPO94, Mellor_Blumberg_JPO04}.
+\np{ln\_crban}\forcode{ = .true.} \citep{craig.banner_JPO94, mellor.blumberg_JPO04}.
 The \np{rn\_crban} namelist parameter is $\alpha_{CB}$ in \autoref{eq:ZDF_Esbc} and
 \np{rn\_charn} provides the value of $\beta$ in \autoref{eq:ZDF_Lsbc}. 
@@ -619 +619 @@
 almost all authors apply a clipping of the length scale as an \textit{ad hoc} remedy.
 With this clipping, the maximum permissible length scale is determined by $l_{max} = c_{lim} \sqrt{2\bar{e}}/ N$.
-A value of $c_{lim} = 0.53$ is often used \citep{Galperin_al_JAS88}.
-\cite{Umlauf_Burchard_CSR05} show that the value of the clipping factor is of crucial importance for
+A value of $c_{lim} = 0.53$ is often used \citep{galperin.kantha.ea_JAS88}.
+\cite{umlauf.burchard_CSR05} show that the value of the clipping factor is of crucial importance for
 the entrainment depth predicted in stably stratified situations,
 and that its value has to be chosen in accordance with the algebraic model for the turbulent fluxes.
@@ -627 +627 @@
 
 The time and space discretization of the GLS equations follows the same energetic consideration as for
-the TKE case described in \autoref{subsec:ZDF_tke_ene} \citep{Burchard_OM02}.
-Examples of performance of the 4 turbulent closure scheme can be found in \citet{Warner_al_OM05}.
+the TKE case described in \autoref{subsec:ZDF_tke_ene} \citep{burchard_OM02}.
+Examples of performance of the 4 turbulent closure scheme can be found in \citet{warner.sherwood.ea_OM05}.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -700 +700 @@
 the water column, but only until the density structure becomes neutrally stable
 (\ie until the mixed portion of the water column has \textit{exactly} the density of the water just below)
-\citep{Madec_al_JPO91}.
+\citep{madec.delecluse.ea_JPO91}.
 The associated algorithm is an iterative process used in the following way (\autoref{fig:npc}):
 starting from the top of the ocean, the first instability is found.
@@ -718 +718 @@
 the algorithm used in \NEMO converges for any profile in a number of iterations which is less than
 the number of vertical levels.
-This property is of paramount importance as pointed out by \citet{Killworth1989}:
+This property is of paramount importance as pointed out by \citet{killworth_iprc89}:
 it avoids the existence of permanent and unrealistic static instabilities at the sea surface.
 This non-penetrative convective algorithm has been proved successful in studies of the deep water formation in
-the north-western Mediterranean Sea \citep{Madec_al_JPO91, Madec_al_DAO91, Madec_Crepon_Bk91}.
+the north-western Mediterranean Sea \citep{madec.delecluse.ea_JPO91, madec.chartier.ea_DAO91, madec.crepon_iprc91}.
 
 The current implementation has been modified in order to deal with any non linear equation of seawater
@@ -748 +748 @@
 In this case, the vertical eddy mixing coefficients are assigned very large values
 (a typical value is $10\;m^2s^{-1})$ in regions where the stratification is unstable
-(\ie when $N^2$ the Brunt-Vais\"{a}l\"{a} frequency is negative) \citep{Lazar_PhD97, Lazar_al_JPO99}.
+(\ie when $N^2$ the Brunt-Vais\"{a}l\"{a} frequency is negative) \citep{lazar_phd97, lazar.madec.ea_JPO99}.
 This is done either on tracers only (\np{nn\_evdm}\forcode{ = 0}) or
 on both momentum and tracers (\np{nn\_evdm}\forcode{ = 1}).
@@ -764 +764 @@
 Note that the stability test is performed on both \textit{before} and \textit{now} values of $N^2$.
 This removes a potential source of divergence of odd and even time step in
-a leapfrog environment \citep{Leclair_PhD2010} (see \autoref{sec:STP_mLF}).
+a leapfrog environment \citep{leclair_phd10} (see \autoref{sec:STP_mLF}).
 
 % -------------------------------------------------------------------------------------------------------------
@@ -807 +807 @@
 The former condition leads to salt fingering and the latter to diffusive convection.
 Double-diffusive phenomena contribute to diapycnal mixing in extensive regions of the ocean.
-\citet{Merryfield1999} include a parameterisation of such phenomena in a global ocean model and show that 
+\citet{merryfield.holloway.ea_JPO99} include a parameterisation of such phenomena in a global ocean model and show that 
 it leads to relatively minor changes in circulation but exerts significant regional influences on
 temperature and salinity.
@@ -842 +842 @@
     \caption{
       \protect\label{fig:zdfddm}
-      From \citet{Merryfield1999} :
+      From \citet{merryfield.holloway.ea_JPO99} :
       (a) Diapycnal diffusivities $A_f^{vT}$ and $A_f^{vS}$ for temperature and salt in regions of salt fingering.
       Heavy curves denote $A^{\ast v} = 10^{-3}~m^2.s^{-1}$ and thin curves $A^{\ast v} = 10^{-4}~m^2.s^{-1}$;
@@ -855 +855 @@
 
 The factor 0.7 in \autoref{eq:zdfddm_f_T} reflects the measured ratio $\alpha F_T /\beta F_S \approx  0.7$ of
-buoyancy flux of heat to buoyancy flux of salt (\eg, \citet{McDougall_Taylor_JMR84}).
-Following  \citet{Merryfield1999}, we adopt $R_c = 1.6$, $n = 6$, and $A^{\ast v} = 10^{-4}~m^2.s^{-1}$.
+buoyancy flux of heat to buoyancy flux of salt (\eg, \citet{mcdougall.taylor_JMR84}).
+Following  \citet{merryfield.holloway.ea_JPO99}, we adopt $R_c = 1.6$, $n = 6$, and $A^{\ast v} = 10^{-4}~m^2.s^{-1}$.
 
 To represent mixing of S and T by diffusive layering,  the diapycnal diffusivities suggested by
@@ -963 +963 @@
 This coefficient is generally estimated by setting a typical decay time $\tau$ in the deep ocean, 
 and setting $r = H / \tau$, where $H$ is the ocean depth.
-Commonly accepted values of $\tau$ are of the order of 100 to 200 days \citep{Weatherly_JMR84}.
+Commonly accepted values of $\tau$ are of the order of 100 to 200 days \citep{weatherly_JMR84}.
 A value $\tau^{-1} = 10^{-7}$~s$^{-1}$ equivalent to 115 days, is usually used in quasi-geostrophic models.
 One may consider the linear friction as an approximation of quadratic friction, $r \approx 2\;C_D\;U_{av}$
-(\citet{Gill1982}, Eq. 9.6.6).
+(\citet{gill_bk82}, Eq. 9.6.6).
 For example, with a drag coefficient $C_D = 0.002$, a typical speed of tidal currents of $U_{av} =0.1$~m\;s$^{-1}$,
 and assuming an ocean depth $H = 4000$~m, the resulting friction coefficient is $r = 4\;10^{-4}$~m\;s$^{-1}$.
@@ -1005 +1005 @@
 internal waves breaking and other short time scale currents.
 A typical value of the drag coefficient is $C_D = 10^{-3} $.
-As an example, the CME experiment \citep{Treguier_JGR92} uses $C_D = 10^{-3}$ and
-$e_b = 2.5\;10^{-3}$m$^2$\;s$^{-2}$, while the FRAM experiment \citep{Killworth1992} uses $C_D = 1.4\;10^{-3}$ and
+As an example, the CME experiment \citep{treguier_JGR92} uses $C_D = 10^{-3}$ and
+$e_b = 2.5\;10^{-3}$m$^2$\;s$^{-2}$, while the FRAM experiment \citep{killworth_JPO92} uses $C_D = 1.4\;10^{-3}$ and
 $e_b =2.5\;\;10^{-3}$m$^2$\;s$^{-2}$.
 The CME choices have been set as default values (\np{rn\_bfri2} and \np{rn\_bfeb2} namelist parameters).
@@ -1235 +1235 @@
 Options are defined through the  \ngn{namzdf\_tmx} namelist variables.
 The parameterization of tidal mixing follows the general formulation for the vertical eddy diffusivity proposed by
-\citet{St_Laurent_al_GRL02} and first introduced in an OGCM by \citep{Simmons_al_OM04}. 
+\citet{st-laurent.simmons.ea_GRL02} and first introduced in an OGCM by \citep{simmons.jayne.ea_OM04}. 
 In this formulation an additional vertical diffusivity resulting from internal tide breaking,
 $A^{vT}_{tides}$ is expressed as a function of $E(x,y)$,
@@ -1252 +1252 @@
 with the remaining $1-q$ radiating away as low mode internal waves and
 contributing to the background internal wave field.
-A value of $q=1/3$ is typically used \citet{St_Laurent_al_GRL02}.
+A value of $q=1/3$ is typically used \citet{st-laurent.simmons.ea_GRL02}.
 The vertical structure function $F(z)$ models the distribution of the turbulent mixing in the vertical.
 It is implemented as a simple exponential decaying upward away from the bottom,
 with a vertical scale of $h_o$ (\np{rn\_htmx} namelist parameter,
-with a typical value of $500\,m$) \citep{St_Laurent_Nash_DSR04}, 
+with a typical value of $500\,m$) \citep{st-laurent.nash_DSR04}, 
 \[
   % \label{eq:Fz}
@@ -1274 +1274 @@
 the unrepresented internal waves induced by the tidal flow over rough topography in a stratified ocean.
 In the current version of \NEMO, the map is built from the output of
-the barotropic global ocean tide model MOG2D-G \citep{Carrere_Lyard_GRL03}.
+the barotropic global ocean tide model MOG2D-G \citep{carrere.lyard_GRL03}.
 This model provides the dissipation associated with internal wave energy for the M2 and K1 tides component
 (\autoref{fig:ZDF_M2_K1_tmx}).
@@ -1280 +1280 @@
 The internal wave energy is thus : $E(x, y) = 1.25 E_{M2} + E_{K1}$.
 Its global mean value is $1.1$ TW,
-in agreement with independent estimates \citep{Egbert_Ray_Nat00, Egbert_Ray_JGR01}. 
+in agreement with independent estimates \citep{egbert.ray_N00, egbert.ray_JGR01}. 
 
 %>>>>>>>>>>>>>>>>>>>>>>>>>>>>
@@ -1288 +1288 @@
     \caption{
       \protect\label{fig:ZDF_M2_K1_tmx}
-      (a) M2 and (b) K1 internal wave drag energy from \citet{Carrere_Lyard_GRL03} ($W/m^2$).
+      (a) M2 and (b) K1 internal wave drag energy from \citet{carrere.lyard_GRL03} ($W/m^2$).
     }
   \end{center}
@@ -1306 +1306 @@
 
 When \np{ln\_tmx\_itf}\forcode{ = .true.}, the two key parameters $q$ and $F(z)$ are adjusted following
-the parameterisation developed by \citet{Koch-Larrouy_al_GRL07}:
+the parameterisation developed by \citet{koch-larrouy.madec.ea_GRL07}:
 
 First, the Indonesian archipelago is a complex geographic region with a series of
@@ -1318 +1318 @@
 Second, the vertical structure function, $F(z)$, is no more associated with a bottom intensification of the mixing,
 but with a maximum of energy available within the thermocline.
-\citet{Koch-Larrouy_al_GRL07} have suggested that the vertical distribution of
+\citet{koch-larrouy.madec.ea_GRL07} have suggested that the vertical distribution of
 the energy dissipation proportional to $N^2$ below the core of the thermocline and to $N$ above. 
 The resulting $F(z)$ is:
@@ -1335 +1335 @@
 Introduced in a regional OGCM, the parameterization improves the water mass characteristics in
 the different Indonesian seas, suggesting that the horizontal and vertical distributions of
-the mixing are adequately prescribed \citep{Koch-Larrouy_al_GRL07, Koch-Larrouy_al_OD08a, Koch-Larrouy_al_OD08b}.
+the mixing are adequately prescribed \citep{koch-larrouy.madec.ea_GRL07, koch-larrouy.madec.ea_OD08*a, koch-larrouy.madec.ea_OD08*b}.
 Note also that such a parameterisation has a significant impact on the behaviour of
-global coupled GCMs \citep{Koch-Larrouy_al_CD10}.
+global coupled GCMs \citep{koch-larrouy.lengaigne.ea_CD10}.
 
 % ================================================================
@@ -1351 +1351 @@
 
 The parameterization of mixing induced by breaking internal waves is a generalization of
-the approach originally proposed by \citet{St_Laurent_al_GRL02}.
+the approach originally proposed by \citet{st-laurent.simmons.ea_GRL02}.
 A three-dimensional field of internal wave energy dissipation $\epsilon(x,y,z)$ is first constructed,
 and the resulting diffusivity is obtained as 
@@ -1361 +1361 @@
 the energy available for mixing.
 If the \np{ln\_mevar} namelist parameter is set to false, the mixing efficiency is taken as constant and
-equal to 1/6 \citep{Osborn_JPO80}.
+equal to 1/6 \citep{osborn_JPO80}.
 In the opposite (recommended) case, $R_f$ is instead a function of
 the turbulence intensity parameter $Re_b = \frac{ \epsilon}{\nu \, N^2}$,
-with $\nu$ the molecular viscosity of seawater, following the model of \cite{Bouffard_Boegman_DAO2013} and
-the implementation of \cite{de_lavergne_JPO2016_efficiency}.
+with $\nu$ the molecular viscosity of seawater, following the model of \cite{bouffard.boegman_DAO13} and
+the implementation of \cite{de-lavergne.madec.ea_JPO16}.
 Note that $A^{vT}_{wave}$ is bounded by $10^{-2}\,m^2/s$, a limit that is often reached when
 the mixing efficiency is constant.
@@ -1371 +1371 @@
 In addition to the mixing efficiency, the ratio of salt to heat diffusivities can chosen to vary 
 as a function of $Re_b$ by setting the \np{ln\_tsdiff} parameter to true, a recommended choice. 
-This parameterization of differential mixing, due to \cite{Jackson_Rehmann_JPO2014},
-is implemented as in \cite{de_lavergne_JPO2016_efficiency}.
+This parameterization of differential mixing, due to \cite{jackson.rehmann_JPO14},
+is implemented as in \cite{de-lavergne.madec.ea_JPO16}.
 
 The three-dimensional distribution of the energy available for mixing, $\epsilon(i,j,k)$,
@@ -1395 +1395 @@
 $h_{cri}$ is related to the large-scale topography of the ocean (etopo2) and
 $h_{bot}$ is a function of the energy flux $E_{bot}$, the characteristic horizontal scale of
-the abyssal hill topography \citep{Goff_JGR2010} and the latitude.
+the abyssal hill topography \citep{goff_JGR10} and the latitude.
 
 % ================================================================