Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_DOMAINcfg.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_DOMAINcfg.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_DOMAINcfg.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Appendix DOMAINcfg : A brief guide to the DOMAINcfg tool
-% ================================================================
 \chapter{A brief guide to the DOMAINcfg tool}
 \label{apdx:DOMCFG}
@@ -19,6 +16,4 @@
 \end{figure}
 
-\newpage
-
 This appendix briefly describes some of the options available in the
 \forcode{DOMAINcfg} tool mentioned in \autoref{chap:DOM}.
@@ -35,7 +30,4 @@
 of those described elsewhere in this manual.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Choice of horizontal grid
-% -------------------------------------------------------------------------------------------------------------
 \section{Choice of horizontal grid}
 \label{sec:DOMCFG_hor}
@@ -55,9 +47,9 @@
 
 \begin{description}
- \item[{\np{jphgr_mesh}{jphgr\_mesh}=0}]  The most general curvilinear orthogonal grids.
+ \item [{\np{jphgr_mesh}{jphgr\_mesh}=0}]  The most general curvilinear orthogonal grids.
   The coordinates and their first derivatives with respect to $i$ and $j$ are provided
   in a input file (\ifile{coordinates}), read in \rou{hgr\_read} subroutine of the domhgr module.
   This is now the only option available within \NEMO\ itself from v4.0 onwards.
-\item[{\np{jphgr_mesh}{jphgr\_mesh}=1 to 5}] A few simple analytical grids are provided (see below).
+\item [{\np{jphgr_mesh}{jphgr\_mesh}=1 to 5}] A few simple analytical grids are provided (see below).
   For other analytical grids, the \mdl{domhgr} module (\texttt{DOMAINcfg} variant) must be
   modified by the user. In most cases, modifying the \mdl{usrdef\_hgr} module of \NEMO\ is
@@ -102,7 +94,4 @@
 (and the number of grid points).
 
-% -------------------------------------------------------------------------------------------------------------
-%        vertical reference coordinate transformation
-% -------------------------------------------------------------------------------------------------------------
 \section{Vertical grid}
 \label{sec:DOMCFG_vert}
@@ -111,5 +100,4 @@
 \label{sec:DOMCFG_zref}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tb]
   \centering
@@ -121,5 +109,4 @@
   \label{fig:DOMCFG_zgr}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The reference coordinate transformation $z_0(k)$ defines the arrays $gdept_0$ and
@@ -210,15 +197,11 @@
 999999$., in \nam{cfg}{cfg} namelist, and specifies instead the four following parameters:
 \begin{itemize}
-\item
-  \np{ppacr}{ppacr}~$= h_{cr}$: stretching factor (nondimensional).
+\item \np{ppacr}{ppacr}~$= h_{cr}$: stretching factor (nondimensional).
   The larger \np{ppacr}{ppacr}, the smaller the stretching.
   Values from $3$ to $10$ are usual.
-\item
-  \np{ppkth}{ppkth}~$= h_{th}$: is approximately the model level at which maximum stretching occurs
+\item \np{ppkth}{ppkth}~$= h_{th}$: is approximately the model level at which maximum stretching occurs
   (nondimensional, usually of order 1/2 or 2/3 of \jp{jpk})
-\item
-  \np{ppdzmin}{ppdzmin}: minimum thickness for the top layer (in meters).
-\item
-  \np{pphmax}{pphmax}: total depth of the ocean (meters).
+\item \np{ppdzmin}{ppdzmin}: minimum thickness for the top layer (in meters).
+\item \np{pphmax}{pphmax}: total depth of the ocean (meters).
 \end{itemize}
 
@@ -227,5 +210,4 @@
 \np{pphmax}{pphmax}~$= 5750~m$.
 
-%% %>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{table}
   \centering
@@ -304,5 +286,4 @@
 \end{table}
 %%%YY
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 %% % -------------------------------------------------------------------------------------------------------------
 %% %        Meter Bathymetry
@@ -314,12 +295,12 @@
 \np{nn_bathy}{nn\_bathy} (found in \nam{dom}{dom} namelist (\texttt{DOMAINCFG} variant) ):
 \begin{description}
-\item[{\np[=0]{nn_bathy}{nn\_bathy}}]:
+\item [{\np[=0]{nn_bathy}{nn\_bathy}}]:
   a flat-bottom domain is defined.
   The total depth $z_w (jpk)$ is given by the coordinate transformation.
   The domain can either be a closed basin or a periodic channel depending on the parameter \np{jperio}{jperio}.
-\item[{\np[=-1]{nn_bathy}{nn\_bathy}}]:
+\item [{\np[=-1]{nn_bathy}{nn\_bathy}}]:
   a domain with a bump of topography one third of the domain width at the central latitude.
   This is meant for the "EEL-R5" configuration, a periodic or open boundary channel with a seamount.
-\item[{\np[=1]{nn_bathy}{nn\_bathy}}]:
+\item [{\np[=1]{nn_bathy}{nn\_bathy}}]:
   read a bathymetry and ice shelf draft (if needed).
   The \ifile{bathy\_meter} file (Netcdf format) provides the ocean depth (positive, in meters) at
@@ -341,13 +322,13 @@
 After reading the bathymetry, the algorithm for vertical grid definition differs between the different options:
 \begin{description}
-\item[\forcode{ln_zco = .true.}]
+\item [\forcode{ln_zco = .true.}]
   set a reference coordinate transformation $z_0(k)$, and set $z(i,j,k,t) = z_0(k)$ where $z_0(k)$ is the closest match to the depth at $(i,j)$.
-\item[\forcode{ln_zps = .true.}]
+\item [\forcode{ln_zps = .true.}]
   set a reference coordinate transformation $z_0(k)$, and calculate the thickness of the deepest level at
   each $(i,j)$ point using the bathymetry, to obtain the final three-dimensional depth and scale factor arrays.
-\item[\forcode{ln_sco = .true.}]
+\item [\forcode{ln_sco = .true.}]
   smooth the bathymetry to fulfill the hydrostatic consistency criteria and
   set the three-dimensional transformation.
-\item[\forcode{s-z and s-zps}]
+\item [\forcode{s-z and s-zps}]
   smooth the bathymetry to fulfill the hydrostatic consistency criteria and
   set the three-dimensional transformation $z(i,j,k)$,
@@ -356,7 +337,4 @@
 %%%
 
-% -------------------------------------------------------------------------------------------------------------
-%        z-coordinate with constant thickness
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection[$Z$-coordinate with uniform thickness levels (\forcode{ln_zco})]{$Z$-coordinate with uniform thickness levels (\protect\np{ln_zco}{ln\_zco})}
 \label{subsec:DOMCFG_zco}
@@ -368,7 +346,4 @@
 rarely used in modern simulations but it can be useful for testing purposes.
 
-% -------------------------------------------------------------------------------------------------------------
-%        z-coordinate with partial step
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection[$Z$-coordinate with partial step (\forcode{ln_zps})]{$Z$-coordinate with partial step (\protect\np{ln_zps}{ln\_zps})}
 \label{subsec:DOMCFG_zps}
@@ -398,7 +373,4 @@
 the default thickness $e_{3t}(jk)$).
 
-% -------------------------------------------------------------------------------------------------------------
-%        s-coordinate
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection[$S$-coordinate (\forcode{ln_sco})]{$S$-coordinate (\protect\np{ln_sco}{ln\_sco})}
 \label{sec:DOMCFG_sco}
@@ -465,5 +437,4 @@
 \]
 
-%% %>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -474,5 +445,4 @@
   \label{fig:DOMCFG_sco_function}
 \end{figure}
-%% %>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 where $H_c$ is the critical depth (\np{rn_hc}{rn\_hc}) at which the coordinate transitions from pure $\sigma$ to
@@ -516,5 +486,4 @@
 where the namelist parameters \np{rn_zb_a}{rn\_zb\_a} and \np{rn_zb_b}{rn\_zb\_b} are $a$ and $b$ respectively.
 
-%% %>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -550,7 +519,4 @@
 and is output as part of the model mesh file at the start of the run.
 
-% -------------------------------------------------------------------------------------------------------------
-%        z*- or s*-coordinate
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection[\zstar- or \sstar-coordinate (\forcode{ln_linssh})]{\zstar- or \sstar-coordinate (\protect\np{ln_linssh}{ln\_linssh})}
 \label{subsec:DOMCFG_zgr_star}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_algos.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_algos.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_algos.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Appendix E : Note on some algorithms
-% ================================================================
 \chapter{Note on some algorithms}
 \label{apdx:ALGOS}
@@ -10,11 +7,6 @@
 \chaptertoc
 
-\newpage
-
 This appendix some on going consideration on algorithms used or planned to be used in \NEMO.
 
-% -------------------------------------------------------------------------------------------------------------
-%        UBS scheme
-% -------------------------------------------------------------------------------------------------------------
 \section{Upstream Biased Scheme (UBS) (\protect\np[=.true.]{ln_traadv_ubs}{ln\_traadv\_ubs})}
 \label{sec:ALGOS_tra_adv_ubs}
@@ -188,7 +180,4 @@
 which leads to ${A_u^{lT}} = \frac{1}{12} {e_{1u}}^3\ |u|$
 
-% -------------------------------------------------------------------------------------------------------------
-%        Leap-Frog energetic
-% -------------------------------------------------------------------------------------------------------------
 \section{Leapfrog energetic}
 \label{sec:ALGOS_LF}
@@ -244,13 +233,6 @@
 In time this boundary condition is not physical and \textbf{add something here!!!}
 
-% ================================================================
-% Iso-neutral diffusion :
-% ================================================================
-
 \section{Lateral diffusion operator}
 
-% ================================================================
-% Griffies' iso-neutral diffusion operator :
-% ================================================================
 \subsection{Griffies iso-neutral diffusion operator}
 
@@ -305,5 +287,4 @@
 the presence of partial cell at the ocean bottom (see \autoref{subsec:ALGOS_Gf_operator}).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -315,5 +296,4 @@
   \label{fig:ALGOS_ISO_triad}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The four iso-neutral fluxes associated with the triads are defined at $T$-point.
@@ -372,5 +352,5 @@
 This expression of the iso-neutral diffusion has been chosen in order to satisfy the following six properties:
 \begin{description}
-\item[$\bullet$ horizontal diffusion]
+\item [$\bullet$ horizontal diffusion]
   The discretization of the diffusion operator recovers the traditional five-point Laplacian in
   the limit of flat iso-neutral direction:
@@ -383,5 +363,5 @@
   \]
 
-\item[$\bullet$ implicit treatment in the vertical]
+\item [$\bullet$ implicit treatment in the vertical]
   In the diagonal term associated with the vertical divergence of the iso-neutral fluxes
   \ie\ the term associated with a second order vertical derivative)
@@ -397,5 +377,5 @@
   can be quite large.
 
-\item[$\bullet$ pure iso-neutral operator]
+\item [$\bullet$ pure iso-neutral operator]
   The iso-neutral flux of locally referenced potential density is zero, \ie
   \begin{align*}
@@ -413,5 +393,5 @@
   the definition of the triads' slopes \autoref{eq:ALGOS_Gf_slopes}.
 
-\item[$\bullet$ conservation of tracer]
+\item [$\bullet$ conservation of tracer]
   The iso-neutral diffusion term conserve the total tracer content, \ie
   \[
@@ -421,5 +401,5 @@
 This property is trivially satisfied since the iso-neutral diffusive operator is written in flux form.
 
-\item[$\bullet$ decrease of tracer variance]
+\item [$\bullet$ decrease of tracer variance]
   The iso-neutral diffusion term does not increase the total tracer variance, \ie
   \[
@@ -434,5 +414,5 @@
 the field on which it is applied become free of grid-point noise.
 
-\item[$\bullet$ self-adjoint operator]
+\item [$\bullet$ self-adjoint operator]
   The iso-neutral diffusion operator is self-adjoint, \ie
   \[
@@ -446,7 +426,4 @@
 \end{description}
 
-% ================================================================
-% Skew flux formulation for Eddy Induced Velocity :
-% ================================================================
 \subsection{Eddy induced velocity and skew flux formulation}
 
@@ -602,7 +579,4 @@
 
 $\ $\newpage      %force an empty line
-% ================================================================
-% Discrete Invariants of the iso-neutral diffrusion
-% ================================================================
 \subsection{Discrete invariants of the iso-neutral diffrusion}
 \label{subsec:ALGOS_Gf_operator}
@@ -767,9 +741,4 @@
 There is no need to develop a specific to obtain it.
 
-\newpage
-
-% ================================================================
-% Discrete Invariants of the skew flux formulation
-% ================================================================
 \subsection{Discrete invariants of the skew flux formulation}
 \label{subsec:ALGOS_eiv_skew}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_diff_opers.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_diff_opers.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_diff_opers.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter Appendix B : Diffusive Operators
-% ================================================================
 \chapter{Diffusive Operators}
 \label{apdx:DIFFOPERS}
@@ -10,9 +7,4 @@
 \chaptertoc
 
-\newpage
-
-% ================================================================
-% Horizontal/Vertical 2nd Order Tracer Diffusive Operators
-% ================================================================
 \section{Horizontal/Vertical $2^{nd}$ order tracer diffusive operators}
 \label{sec:DIFFOPERS_1}
@@ -156,7 +148,4 @@
 %\addtocounter{equation}{-2}
 
-% ================================================================
-% Isopycnal/Vertical 2nd Order Tracer Diffusive Operators
-% ================================================================
 \section{Iso/Diapycnal $2^{nd}$ order tracer diffusive operators}
 \label{sec:DIFFOPERS_2}
@@ -320,8 +309,4 @@
 \autoref{sec:DIFFOPERS_1} onto $s$-coordinates is exact, however steep the $s$-surfaces.
 
-
-% ================================================================
-% Lateral/Vertical Momentum Diffusive Operators
-% ================================================================
 \section{Lateral/Vertical momentum diffusive operators}
 \label{sec:DIFFOPERS_3}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_invariants.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_invariants.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_invariants.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter Ñ Appendix C : Discrete Invariants of the Equations
-% ================================================================
 \chapter{Discrete Invariants of the Equations}
 \label{apdx:INVARIANTS}
@@ -15,9 +12,4 @@
 %\gmcomment{
 
-\newpage
-
-% ================================================================
-% Introduction / Notations
-% ================================================================
 \section{Introduction / Notations}
 \label{sec:INVARIANTS_0}
@@ -93,7 +85,4 @@
 \end{flalign}
 
-% ================================================================
-% Continuous Total energy Conservation
-% ================================================================
 \section{Continuous conservation}
 \label{sec:INVARIANTS_1}
@@ -322,13 +311,7 @@
 %
 
-% ================================================================
-% Discrete Total energy Conservation : vector invariant form
-% ================================================================
 \section{Discrete total energy conservation: vector invariant form}
 \label{sec:INVARIANTS_2}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Total energy conservation
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Total energy conservation}
 \label{subsec:INVARIANTS_KE+PE_vect}
@@ -354,7 +337,4 @@
 leads to the discrete equivalent of the four equations \autoref{eq:INVARIANTS_E_tot_flux}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Vorticity term (coriolis + vorticity part of the advection)
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Vorticity term (coriolis + vorticity part of the advection)}
 \label{subsec:INVARIANTS_vor}
@@ -363,7 +343,4 @@
 or the planetary ($q=f/e_{3f}$), or the total potential vorticity ($q=(\zeta +f) /e_{3f}$).
 Two discretisation of the vorticity term (ENE and EEN) allows the conservation of the kinetic energy.
-% -------------------------------------------------------------------------------------------------------------
-%       Vorticity Term with ENE scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Vorticity term with ENE scheme (\protect\np[=.true.]{ln_dynvor_ene}{ln\_dynvor\_ene})}
 \label{subsec:INVARIANTS_vorENE}
@@ -403,7 +380,4 @@
 In other words, the domain averaged kinetic energy does not change due to the vorticity term.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Vorticity Term with EEN scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Vorticity term with EEN scheme (\protect\np[=.true.]{ln_dynvor_een}{ln\_dynvor\_een})}
 \label{subsec:INVARIANTS_vorEEN_vect}
@@ -475,7 +449,4 @@
 \end{flalign*}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Gradient of Kinetic Energy / Vertical Advection
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Gradient of kinetic energy / Vertical advection}
 \label{subsec:INVARIANTS_zad}
@@ -585,8 +556,4 @@
 Blah blah required on the the step representation of bottom topography.....
 
-
-% -------------------------------------------------------------------------------------------------------------
-%       Pressure Gradient Term
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Pressure gradient term}
 \label{subsec:INVARIANTS_2.6}
@@ -731,13 +698,7 @@
 Nevertheless, it is almost never satisfied since a linear equation of state is rarely used.
 
-% ================================================================
-% Discrete Total energy Conservation : flux form
-% ================================================================
 \section{Discrete total energy conservation: flux form}
 \label{sec:INVARIANTS_3}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Total energy conservation
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Total energy conservation}
 \label{subsec:INVARIANTS_KE+PE_flux}
@@ -760,14 +721,7 @@
 vector invariant or in flux form, leads to the discrete equivalent of the ????
 
-
-% -------------------------------------------------------------------------------------------------------------
-%       Coriolis and advection terms: flux form
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Coriolis and advection terms: flux form}
 \label{subsec:INVARIANTS_3.2}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Coriolis plus ``metric'' Term
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Coriolis plus ``metric'' term}
 \label{subsec:INVARIANTS_3.3}
@@ -788,7 +742,4 @@
 The derivation is the same as for the vorticity term in the vector invariant form (\autoref{subsec:INVARIANTS_vor}).
 
-% -------------------------------------------------------------------------------------------------------------
-%       Flux form advection
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Flux form advection}
 \label{subsec:INVARIANTS_3.4}
@@ -869,13 +820,7 @@
 The horizontal kinetic energy is not conserved, but forced to decay (\ie\ the scheme is diffusive).
 
-% ================================================================
-% Discrete Enstrophy Conservation
-% ================================================================
 \section{Discrete enstrophy conservation}
 \label{sec:INVARIANTS_4}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Vorticity Term with ENS scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Vorticity term with ENS scheme  (\protect\np[=.true.]{ln_dynvor_ens}{ln\_dynvor\_ens})}
 \label{subsec:INVARIANTS_vorENS}
@@ -944,7 +889,4 @@
 The later equality is obtain only when the flow is horizontally non-divergent, \ie\ $\chi$=$0$.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Vorticity Term with EEN scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Vorticity Term with EEN scheme (\protect\np[=.true.]{ln_dynvor_een}{ln\_dynvor\_een})}
 \label{subsec:INVARIANTS_vorEEN}
@@ -1017,7 +959,4 @@
 \end{flalign*}
 
-% ================================================================
-% Conservation Properties on Tracers
-% ================================================================
 \section{Conservation properties on tracers}
 \label{sec:INVARIANTS_5}
@@ -1033,7 +972,4 @@
 as the equation of state is non linear with respect to $T$ and $S$.
 In practice, the mass is conserved to a very high accuracy.
-% -------------------------------------------------------------------------------------------------------------
-%       Advection Term
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Advection term}
 \label{subsec:INVARIANTS_5.1}
@@ -1099,7 +1035,4 @@
 which is the discrete form of $ \frac{1}{2} \int_D {  T^2 \frac{1}{e_3} \frac{\partial  e_3 }{\partial t} \;dv }$.
 
-% ================================================================
-% Conservation Properties on Lateral Momentum Physics
-% ================================================================
 \section{Conservation properties on lateral momentum physics}
 \label{sec:INVARIANTS_dynldf_properties}
@@ -1120,7 +1053,4 @@
 the term associated with the horizontal gradient of the divergence is locally zero.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation of Potential Vorticity
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Conservation of potential vorticity}
 \label{subsec:INVARIANTS_6.1}
@@ -1154,7 +1084,4 @@
 \end{flalign*}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Dissipation of Horizontal Kinetic Energy
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Dissipation of horizontal kinetic energy}
 \label{subsec:INVARIANTS_6.2}
@@ -1206,7 +1133,4 @@
 \]
 
-% -------------------------------------------------------------------------------------------------------------
-%       Dissipation of Enstrophy
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Dissipation of enstrophy}
 \label{subsec:INVARIANTS_6.3}
@@ -1230,7 +1154,4 @@
 \end{flalign*}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation of Horizontal Divergence
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Conservation of horizontal divergence}
 \label{subsec:INVARIANTS_6.4}
@@ -1257,7 +1178,4 @@
 \end{flalign*}
 
-% -------------------------------------------------------------------------------------------------------------
-%       Dissipation of Horizontal Divergence Variance
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Dissipation of horizontal divergence variance}
 \label{subsec:INVARIANTS_6.5}
@@ -1283,7 +1201,4 @@
 \end{flalign*}
 
-% ================================================================
-% Conservation Properties on Vertical Momentum Physics
-% ================================================================
 \section{Conservation properties on vertical momentum physics}
 \label{sec:INVARIANTS_7}
@@ -1454,7 +1369,4 @@
 \end{flalign*}
 
-% ================================================================
-% Conservation Properties on Tracer Physics
-% ================================================================
 \section{Conservation properties on tracer physics}
 \label{sec:INVARIANTS_8}
@@ -1466,7 +1378,4 @@
 As for the advection term, there is conservation of mass only if the Equation Of Seawater is linear.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation of Tracers
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Conservation of tracers}
 \label{subsec:INVARIANTS_8.1}
@@ -1499,7 +1408,4 @@
 In fact, this property simply results from the flux form of the operator.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Dissipation of Tracer Variance
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Dissipation of tracer variance}
 \label{subsec:INVARIANTS_8.2}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_s_coord.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_s_coord.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_s_coord.tex	(revision 11596)
@@ -3,7 +3,4 @@
 \begin{document}
 
-% ================================================================
-% Chapter Appendix A : Curvilinear s-Coordinate Equations
-% ================================================================
 \chapter{Curvilinear $s-$Coordinate Equations}
 \label{apdx:SCOORD}
@@ -21,10 +18,4 @@
 \end{figure}
 
-
-\newpage
-
-% ================================================================
-% Chain rule
-% ================================================================
 \section{Chain rule for $s-$coordinates}
 \label{sec:SCOORD_chain}
@@ -121,8 +112,4 @@
 \end{equation}
 
-
-% ================================================================
-% continuity equation
-% ================================================================
 \section{Continuity equation in $s-$coordinates}
 \label{sec:SCOORD_continuity}
@@ -240,8 +227,4 @@
 the contribution of the time variation of the vertical coordinate to the volume budget.
 
-
-% ================================================================
-% momentum equation
-% ================================================================
 \section{Momentum equation in $s-$coordinate}
 \label{sec:SCOORD_momentum}
@@ -571,8 +554,4 @@
 \ie\ the volume flux across the moving $s$-surfaces per unit horizontal area.
 
-
-% ================================================================
-% Tracer equation
-% ================================================================
 \section{Tracer equation}
 \label{sec:SCOORD_tracer}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_triads.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_triads.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/apdx_triads.tex	(revision 11596)
@@ -12,13 +12,8 @@
 
 \begin{document}
-% ================================================================
-% Iso-neutral diffusion :
-% ================================================================
 \chapter{Iso-Neutral Diffusion and Eddy Advection using Triads}
 \label{apdx:TRIADS}
 
 \chaptertoc
-
-\newpage
 
 \section[Choice of \forcode{namtra\_ldf} namelist parameters]{Choice of \protect\nam{tra_ldf}{tra\_ldf} namelist parameters}
@@ -42,5 +37,5 @@
 The options specific to the Griffies scheme include:
 \begin{description}
-\item[{\np{ln_triad_iso}{ln\_triad\_iso}}]
+\item [{\np{ln_triad_iso}{ln\_triad\_iso}}]
   See \autoref{sec:TRIADS_taper}.
   If this is set false (the default),
@@ -53,5 +48,5 @@
   giving an almost pure horizontal diffusive tracer flux within the mixed layer.
   This is similar to the tapering suggested by \citet{gerdes.koberle.ea_CD91}. See \autoref{subsec:TRIADS_Gerdes-taper}
-\item[{\np{ln_botmix_triad}{ln\_botmix\_triad}}]
+\item [{\np{ln_botmix_triad}{ln\_botmix\_triad}}]
   See \autoref{sec:TRIADS_iso_bdry}.
   If this is set false (the default) then the lateral diffusive fluxes
@@ -59,11 +54,11 @@
   If it is set true, however, then these lateral diffusive fluxes are applied,
   giving smoother bottom tracer fields at the cost of introducing diapycnal mixing.
-\item[{\np{rn_sw_triad}{rn\_sw\_triad}}]
+\item [{\np{rn_sw_triad}{rn\_sw\_triad}}]
   blah blah to be added....
 \end{description}
 The options shared with the Standard scheme include:
 \begin{description}
-\item[{\np{ln_traldf_msc}{ln\_traldf\_msc}}]   blah blah to be added
-\item[{\np{rn_slpmax}{rn\_slpmax}}]  blah blah to be added
+\item [{\np{ln_traldf_msc}{ln\_traldf\_msc}}]   blah blah to be added
+\item [{\np{rn_slpmax}{rn\_slpmax}}]  blah blah to be added
 \end{description}
 
@@ -548,5 +543,5 @@
 The diffusion scheme satisfies the following six properties:
 \begin{description}
-\item[$\bullet$ horizontal diffusion]
+\item [$\bullet$ horizontal diffusion]
   The discretization of the diffusion operator recovers the traditional five-point Laplacian
   \autoref{eq:TRIADS_lat-normal} in the limit of flat iso-neutral direction:
@@ -559,5 +554,5 @@
   \]
 
-\item[$\bullet$ implicit treatment in the vertical]
+\item [$\bullet$ implicit treatment in the vertical]
   Only tracer values associated with a single water column appear in the expression \autoref{eq:TRIADS_i33} for
   the $_{33}$ fluxes, vertical fluxes driven by vertical gradients.
@@ -575,9 +570,9 @@
   (where $b_w= e_{1w}\,e_{2w}\,e_{3w}$ is the volume of $w$-cells) can be quite large.
 
-\item[$\bullet$ pure iso-neutral operator]
+\item [$\bullet$ pure iso-neutral operator]
   The iso-neutral flux of locally referenced potential density is zero.
   See \autoref{eq:TRIADS_latflux-rho} and \autoref{eq:TRIADS_vertflux-triad2}.
 
-\item[$\bullet$ conservation of tracer]
+\item [$\bullet$ conservation of tracer]
   The iso-neutral diffusion conserves tracer content, \ie
   \[
@@ -587,5 +582,5 @@
   This property is trivially satisfied since the iso-neutral diffusive operator is written in flux form.
 
-\item[$\bullet$ no increase of tracer variance]
+\item [$\bullet$ no increase of tracer variance]
   The iso-neutral diffusion does not increase the tracer variance, \ie
   \[
@@ -600,5 +595,5 @@
   the field on which it is applied becomes free of grid-point noise.
 
-\item[$\bullet$ self-adjoint operator]
+\item [$\bullet$ self-adjoint operator]
   The iso-neutral diffusion operator is self-adjoint, \ie
   \begin{equation}
@@ -753,6 +748,5 @@
 described above by \autoref{eq:TRIADS_Rtilde}.
 \begin{enumerate}
-\item
-  Mixed-layer depth is defined so as to avoid including regions of weak vertical stratification in
+\item Mixed-layer depth is defined so as to avoid including regions of weak vertical stratification in
   the slope definition.
   At each $i,j$ (simplified to $i$ in \autoref{fig:TRIADS_MLB_triad}),
@@ -766,6 +760,5 @@
   output the diagnosed mixed-layer depth $h_{\mathrm{ML}}=|z_{W}|_{k_{\mathrm{ML}}+1/2}$,
   the depth of the $w$-point above the $i,k_{\mathrm{ML}}$ tracer point.
-\item
-  We define `basal' triad slopes ${\:}_i{\mathbb{R}_{\mathrm{base}}}_{\,i_p}^{k_p}$ as
+\item We define `basal' triad slopes ${\:}_i{\mathbb{R}_{\mathrm{base}}}_{\,i_p}^{k_p}$ as
   the slopes of those triads whose vertical `arms' go down from the $i,k_{\mathrm{ML}}$ tracer point to
   the $i,k_{\mathrm{ML}}-1$ tracer point below.
@@ -790,6 +783,5 @@
 one gridbox deeper than the diagnosed ML depth $z_{\mathrm{ML}})$ that sets the $h$ used to taper the slopes in
 \autoref{eq:TRIADS_rmtilde}.
-\item
-  Finally, we calculate the adjusted triads ${\:}_i^k{\mathbb{R}_{\mathrm{ML}}}_{\,i_p}^{k_p}$ within
+\item Finally, we calculate the adjusted triads ${\:}_i^k{\mathbb{R}_{\mathrm{ML}}}_{\,i_p}^{k_p}$ within
   the mixed layer, by multiplying the appropriate ${\:}_i{\mathbb{R}_{\mathrm{base}}}_{\,i_p}^{k_p}$ by
   the ratio of the depth of the $w$-point ${z_w}_{k+k_p}$ to ${z_{\mathrm{base}}}_{\,i}$.
@@ -872,7 +864,4 @@
 % This may give strange looking results,
 % particularly where the mixed-layer depth varies strongly laterally.
-% ================================================================
-% Skew flux formulation for Eddy Induced Velocity :
-% ================================================================
 \section{Eddy induced advection formulated as a skew flux}
 \label{sec:TRIADS_skew-flux}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ASM.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ASM.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ASM.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter Assimilation increments (ASM)
-% ================================================================
 \chapter{Apply Assimilation Increments (ASM)}
 \label{chap:ASM}
@@ -19,6 +16,4 @@
 \end{tabular}
 \end{figure}
-
-\newpage
 
 The ASM code adds the functionality to apply increments to the model variables: temperature, salinity,
@@ -138,5 +133,4 @@
 This specifies the number of iterations of the divergence damping. Setting a value of the order of 100 will result in a significant reduction in the vertical velocity induced by the increments.
 
-
 %==========================================================================
 
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIA.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIA.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIA.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter I/O & Diagnostics
-% ================================================================
 \chapter{Output and Diagnostics (IOM, DIA, TRD, FLO)}
 \label{chap:DIA}
@@ -23,9 +20,4 @@
 \end{figure}
 
-\newpage
-
-% ================================================================
-%       Old Model Output
-% ================================================================
 \section{Model output}
 \label{sec:DIA_io_old}
@@ -53,7 +45,4 @@
 %\gmcomment{                    % start of gmcomment
 
-% ================================================================
-% Diagnostics
-% ================================================================
 \section{Standard model output (IOM)}
 \label{sec:DIA_iom}
@@ -64,9 +53,7 @@
 
 \begin{enumerate}
-\item
-  The complete and flexible control of the output files through external XML files adapted by
+\item The complete and flexible control of the output files through external XML files adapted by
   the user from standard templates.
-\item
-  To achieve high performance and scalable output through the optional distribution of
+\item To achieve high performance and scalable output through the optional distribution of
   all diagnostic output related tasks to dedicated processes.
 \end{enumerate}
@@ -76,17 +63,11 @@
 
 \begin{itemize}
-\item
-  The choice of output frequencies that can be different for each file (including real months and years).
-\item
-  The choice of file contents; includes complete flexibility over which data are written in which files
+\item The choice of output frequencies that can be different for each file (including real months and years).
+\item The choice of file contents; includes complete flexibility over which data are written in which files
   (the same data can be written in different files).
-\item
-  The possibility to split output files at a chosen frequency.
-\item
-  The possibility to extract a vertical or an horizontal subdomain.
-\item
-  The choice of the temporal operation to perform, \eg: average, accumulate, instantaneous, min, max and once.
-\item
-  Control over metadata via a large XML "database" of possible output fields.
+\item The possibility to split output files at a chosen frequency.
+\item The possibility to extract a vertical or an horizontal subdomain.
+\item The choice of the temporal operation to perform, \eg: average, accumulate, instantaneous, min, max and once.
+\item Control over metadata via a large XML "database" of possible output fields.
 \end{itemize}
 
@@ -150,8 +131,8 @@
 If an additional variable must be written to a restart file, the following steps are needed:
 \begin{description}
-   \item[step 1:] add variable name to a list of restart variables (in subroutine \rou{iom\_set\_rst\_vars,} \mdl{iom}) and
+   \item [step 1:] add variable name to a list of restart variables (in subroutine \rou{iom\_set\_rst\_vars,} \mdl{iom}) and
 define correct grid for the variable (\forcode{grid_N_3D} - 3D variable, \forcode{grid_N} - 2D variable, \forcode{grid_vector} -
 1D variable, \forcode{grid_scalar} - scalar),
-   \item[step 2:] add variable to the list of fields written by restart.  This can be done either in subroutine
+   \item [step 2:] add variable to the list of fields written by restart.  This can be done either in subroutine
 \rou{iom\_set\_rstw\_core} (\mdl{iom}) or by calling  \rou{iom\_set\_rstw\_active} (\mdl{iom}) with the name of a variable
 as an argument. This convention follows approach for writing restart using iom, where variables are
@@ -159,7 +140,5 @@
 \end{description}
 
-
 An older versions of XIOS do not support reading functionality. It's recommended to use at least XIOS2@1451.
-
 
 \subsection{XIOS: XML Inputs-Outputs Server}
@@ -276,10 +255,10 @@
 
 \begin{enumerate}
-\item[1.]
+\item [1.]
   in \NEMO\ code, add a \forcode{CALL iom_put( 'identifier', array )} where you want to output a 2D or 3D array.
-\item[2.]
+\item [2.]
   If necessary, add \forcode{USE iom ! I/O manager library} to the list of used modules in
   the upper part of your module.
-\item[3.]
+\item [3.]
   in the field\_def.xml file, add the definition of your variable using the same identifier you used in the f90 code
   (see subsequent sections for a details of the XML syntax and rules).
@@ -310,5 +289,5 @@
 \xmlcode{<field_group id="SBC" ...>} which has been defined with the correct frequency of operations
 (iom\_set\_field\_attr in \mdl{iom})
-\item[4.]
+\item [4.]
   add your field in one of the output files defined in iodef.xml
   (again see subsequent sections for syntax and rules)
@@ -737,6 +716,5 @@
 
 \begin{enumerate}
-\item
-  Simple computation: directly define the computation when refering to the variable in the file definition.
+\item Simple computation: directly define the computation when refering to the variable in the file definition.
 
 \begin{xmllines}
@@ -746,6 +724,5 @@
 \end{xmllines}
 
-\item
-  Simple computation: define a new variable and use it in the file definition.
+\item Simple computation: define a new variable and use it in the file definition.
 
 in field\_definition:
@@ -764,6 +741,5 @@
 sst2 won't be evaluated.
 
-\item
-  Change of variable precision:
+\item Change of variable precision:
 
 \begin{xmllines}
@@ -778,6 +754,5 @@
 Forcing double precision outputs with prec="8" (for example in the field\_definition) will avoid this problem.
 
-\item
-  add user defined attributes:
+\item add user defined attributes:
 
 \begin{xmllines}
@@ -794,6 +769,5 @@
 \end{xmllines}
 
-\item
-  use of the ``@'' function: example 1, weighted temporal average
+\item use of the ``@'' function: example 1, weighted temporal average
 
  - define a new variable in field\_definition
@@ -823,6 +797,5 @@
 Note that in this case, freq\_op must be equal to the file output\_freq.
 
-\item
-  use of the ``@'' function: example 2, monthly SSH standard deviation
+\item use of the ``@'' function: example 2, monthly SSH standard deviation
 
  - define a new variable in field\_definition
@@ -854,6 +827,5 @@
 Note that in this case, freq\_op must be equal to the file output\_freq.
 
-\item
-  use of the ``@'' function: example 3, monthly average of SST diurnal cycle
+\item use of the ``@'' function: example 3, monthly average of SST diurnal cycle
 
  - define 2 new variables in field\_definition
@@ -1326,8 +1298,4 @@
 This must be set to true if these metadata are to be included in the output files.
 
-
-% ================================================================
-%       NetCDF4 support
-% ================================================================
 \section[NetCDF4 support (\texttt{\textbf{key\_netcdf4}})]{NetCDF4 support (\protect\key{netcdf4})}
 \label{sec:DIA_nc4}
@@ -1446,7 +1414,4 @@
 the invidual processing regions and different chunking choices may be desired.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Tracer/Dynamics Trends
-% -------------------------------------------------------------------------------------------------------------
 \section[Tracer/Dynamics trends (\forcode{&namtrd})]{Tracer/Dynamics trends (\protect\nam{trd}{trd})}
 \label{sec:DIA_trd}
@@ -1470,21 +1435,21 @@
 
 \begin{description}
-\item[{\np{ln_glo_trd}{ln\_glo\_trd}}]:
+\item [{\np{ln_glo_trd}{ln\_glo\_trd}}]:
   at each \np{nn_trd}{nn\_trd} time-step a check of the basin averaged properties of
   the momentum and tracer equations is performed.
   This also includes a check of $T^2$, $S^2$, $\tfrac{1}{2} (u^2+v2)$,
   and potential energy time evolution equations properties;
-\item[{\np{ln_dyn_trd}{ln\_dyn\_trd}}]:
+\item [{\np{ln_dyn_trd}{ln\_dyn\_trd}}]:
   each 3D trend of the evolution of the two momentum components is output;
-\item[{\np{ln_dyn_mxl}{ln\_dyn\_mxl}}]:
+\item [{\np{ln_dyn_mxl}{ln\_dyn\_mxl}}]:
   each 3D trend of the evolution of the two momentum components averaged over the mixed layer is output;
-\item[{\np{ln_vor_trd}{ln\_vor\_trd}}]:
+\item [{\np{ln_vor_trd}{ln\_vor\_trd}}]:
   a vertical summation of the moment tendencies is performed,
   then the curl is computed to obtain the barotropic vorticity tendencies which are output;
-\item[{\np{ln_KE_trd}{ln\_KE\_trd}}] :
+\item [{\np{ln_KE_trd}{ln\_KE\_trd}}] :
   each 3D trend of the Kinetic Energy equation is output;
-\item[{\np{ln_tra_trd}{ln\_tra\_trd}}]:
+\item [{\np{ln_tra_trd}{ln\_tra\_trd}}]:
   each 3D trend of the evolution of temperature and salinity is output;
-\item[{\np{ln_tra_mxl}{ln\_tra\_mxl}}]:
+\item [{\np{ln_tra_mxl}{ln\_tra\_mxl}}]:
   each 2D trend of the evolution of temperature and salinity averaged over the mixed layer is output;
 \end{description}
@@ -1497,7 +1462,4 @@
 and none of the options have been tested with variable volume (\ie\ \np[=.true.]{ln_linssh}{ln\_linssh}).
 
-% -------------------------------------------------------------------------------------------------------------
-%       On-line Floats trajectories
-% -------------------------------------------------------------------------------------------------------------
 \section[FLO: On-Line Floats trajectories (\texttt{\textbf{key\_floats}})]{FLO: On-Line Floats trajectories (\protect\key{floats})}
 \label{sec:DIA_FLO}
@@ -1602,8 +1564,4 @@
 \end{xmllines}
 
-
-% -------------------------------------------------------------------------------------------------------------
-%       Harmonic analysis of tidal constituents
-% -------------------------------------------------------------------------------------------------------------
 \section[Harmonic analysis of tidal constituents (\texttt{\textbf{key\_diaharm}})]{Harmonic analysis of tidal constituents (\protect\key{diaharm})}
 \label{sec:DIA_diag_harm}
@@ -1654,7 +1612,4 @@
 We obtain in output $C_{j}$ and $S_{j}$ for each tidal wave.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Sections transports
-% -------------------------------------------------------------------------------------------------------------
 \section[Transports across sections (\texttt{\textbf{key\_diadct}})]{Transports across sections (\protect\key{diadct})}
 \label{sec:DIA_diag_dct}
@@ -1800,10 +1755,6 @@
 \end{table}
 
-% ================================================================
-% Steric effect in sea surface height
-% ================================================================
 \section{Diagnosing the steric effect in sea surface height}
 \label{sec:DIA_steric}
-
 
 Changes in steric sea level are caused when changes in the density of the water column imply an expansion or
@@ -1980,7 +1931,4 @@
 Both steric and thermosteric sea level are computed in \mdl{diaar5}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Other Diagnostics
-% -------------------------------------------------------------------------------------------------------------
 \section{Other diagnostics}
 \label{sec:DIA_diag_others}
@@ -2001,6 +1949,4 @@
 - the depth of the thermocline (maximum of the vertical temperature gradient) (\mdl{diahth})
 
-
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -2019,5 +1965,4 @@
   \label{fig:DIA_mask_subasins}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 % -----------------------------------------------------------
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIU.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIU.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DIU.tex	(revision 11596)
@@ -2,8 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Diurnal SST models (DIU)
-% Edited by James While
-% ================================================================
 \chapter{Diurnal SST Models (DIU)}
 \label{chap:DIU}
@@ -11,6 +7,4 @@
 \chaptertoc
 
-
-\newpage
 $\ $\newline % force a new line
 
@@ -19,11 +13,8 @@
 The skin temperature can be split into three parts:
 \begin{itemize}
-\item
-  A foundation SST which is free from diurnal warming.
-\item
-  A warm layer, typically ~3\,m thick,
+\item A foundation SST which is free from diurnal warming.
+\item A warm layer, typically ~3\,m thick,
   where heating from solar radiation can cause a warm stably stratified layer during the daytime
-\item
-  A cool skin, a thin layer, approximately ~1\, mm thick,
+\item A cool skin, a thin layer, approximately ~1\, mm thick,
   where long wave cooling is dominant and cools the immediate ocean surface.
 \end{itemize}
@@ -46,7 +37,7 @@
 This namelist contains only two variables:
 \begin{description}
-\item[{\np{ln_diurnal}{ln\_diurnal}}]
+\item [{\np{ln_diurnal}{ln\_diurnal}}]
   A logical switch for turning on/off both the cool skin and warm layer.
-\item[{\np{ln_diurnal_only}{ln\_diurnal\_only}}]
+\item [{\np{ln_diurnal_only}{ln\_diurnal\_only}}]
   A logical switch which if \forcode{.true.} will run the diurnal model without the other dynamical parts of \NEMO.
   \np{ln_diurnal_only}{ln\_diurnal\_only} must be \forcode{.false.} if \np{ln_diurnal}{ln\_diurnal} is \forcode{.false.}.
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DOM.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DOM.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DOM.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter 2 ——— Space and Time Domain (DOM)
-% ================================================================
 \chapter{Space Domain (DOM)}
 \label{chap:DOM}
@@ -36,6 +33,4 @@
 \end{table}
 
-\newpage
-
 Having defined the continuous equations in \autoref{chap:MB} and chosen a time discretisation \autoref{chap:TD},
 we need to choose a grid for spatial discretisation and related numerical algorithms.
@@ -43,17 +38,10 @@
 and other relevant information about the DOM (DOMain) source code modules.
 
-% ================================================================
-% Fundamentals of the Discretisation
-% ================================================================
 \section{Fundamentals of the discretisation}
 \label{sec:DOM_basics}
 
-% -------------------------------------------------------------------------------------------------------------
-%        Arrangement of Variables
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Arrangement of variables}
 \label{subsec:DOM_cell}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tb]
   \centering
@@ -68,5 +56,4 @@
   \label{fig:DOM_cell}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The numerical techniques used to solve the Primitive Equations in this model are based on the traditional,
@@ -99,5 +86,4 @@
 (see \autoref{eq:DOM_bar} in the next section).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{table}[!tb]
   \centering
@@ -130,5 +116,4 @@
   \label{tab:DOM_cell}
 \end{table}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 Note that the definition of the scale factors
@@ -145,5 +130,4 @@
 (rather than allowing the user to set arbitrary jumps in thickness between adjacent layers) \citep{treguier.dukowicz.ea_JGR96}.
 An example of the effect of such a choice is shown in \autoref{fig:DOM_zgr_e3}.
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -160,9 +144,5 @@
   \label{fig:DOM_zgr_e3}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        Vector Invariant Formulation
-% -------------------------------------------------------------------------------------------------------------
+
 \subsection{Discrete operators}
 \label{subsec:DOM_operators}
@@ -255,11 +235,7 @@
 demonstrate integral conservative properties of the discrete formulation chosen.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Numerical Indexing
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Numerical indexing}
 \label{subsec:DOM_Num_Index}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tb]
   \centering
@@ -271,5 +247,4 @@
   \label{fig:DOM_index_hor}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The array representation used in the \fortran\ code requires an integer indexing.
@@ -315,5 +290,4 @@
 accommodate the opposing vertical index directions in implementation and documentation.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!pt]
   \centering
@@ -326,9 +300,5 @@
   \label{fig:DOM_index_vert}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        Domain configuration
-% -------------------------------------------------------------------------------------------------------------
+
 \section{Spatial domain configuration}
 \label{subsec:DOM_config}
@@ -385,13 +355,7 @@
 See \autoref{sec:LBC_jperio} for details on the available options and the corresponding values for \jp{jperio}.
 
-% ================================================================
-% Domain: Horizontal Grid (mesh)
-% ================================================================
 \subsection[Horizontal grid mesh (\textit{domhgr.F90}]{Horizontal grid mesh (\protect\mdl{domhgr})}
 \label{subsec:DOM_hgr}
 
-% ================================================================
-% Domain: List of hgr-related fields needed
-% ================================================================
 \subsubsection{Required fields}
 \label{sec:DOM_hgr_fields}
@@ -453,8 +417,4 @@
 thus no specific arrays are defined at $w$ points.
 
-
-% ================================================================
-% Domain: Vertical Grid (domzgr)
-% ================================================================
 \subsection[Vertical grid (\textit{domzgr.F90})]{Vertical grid (\protect\mdl{domzgr})}
 \label{subsec:DOM_zgr}
@@ -476,5 +436,4 @@
 \end{enumerate}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tb]
   \centering
@@ -491,5 +450,4 @@
   \label{fig:DOM_z_zps_s_sps}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The choice of a vertical coordinate is made when setting up the configuration;
@@ -583,11 +541,6 @@
 With ice cavities, \jp{top\_level} determines the first wet point below the overlying ice shelf.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        level bathymetry and mask
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Level bathymetry and mask}
 \label{subsec:DOM_msk}
-
 
 From \jp{top\_level} and \jp{bottom\_level} fields, the mask fields are defined as follows:
@@ -619,5 +572,4 @@
 %% (see \autoref{fig:LBC_jperio}).
 
-
 %-------------------------------------------------------------------------------------------------
 %        Closed seas
@@ -643,7 +595,4 @@
 \end{clines}
 
-% -------------------------------------------------------------------------------------------------------------
-%        Grid files
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Output grid files}
 \label{subsec:DOM_meshmask}
@@ -664,7 +613,4 @@
 This file contains additional fields that can be useful for post-processing applications.
 
-% ================================================================
-% Domain: Initial State (dtatsd & istate)
-% ================================================================
 \section[Initial state (\textit{istate.F90} and \textit{dtatsd.F90})]{Initial state (\protect\mdl{istate} and \protect\mdl{dtatsd})}
 \label{sec:DOM_DTA_tsd}
@@ -682,5 +628,5 @@
 
 \begin{description}
-\item[{\np[=.true.]{ln_tsd_init}{ln\_tsd\_init}}]
+\item [{\np[=.true.]{ln_tsd_init}{ln\_tsd\_init}}]
   Use T and S input files that can be given on the model grid itself or on their native input data grids.
   In the latter case, the data will be interpolated on-the-fly both in the horizontal and the vertical to the model grid
@@ -688,5 +634,5 @@
   The information relating to the input files are specified in the \np{sn_tem}{sn\_tem} and \np{sn_sal}{sn\_sal} structures.
   The computation is done in the \mdl{dtatsd} module.
-\item[{\np[=.false.]{ln_tsd_init}{ln\_tsd\_init}}]
+\item [{\np[=.false.]{ln_tsd_init}{ln\_tsd\_init}}]
   Initial values for T and S are set via a user supplied \rou{usr\_def\_istate} routine contained in \mdl{userdef\_istate}.
   The default version sets horizontally uniform T and profiles as used in the GYRE configuration
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DYN.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DYN.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_DYN.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter ——— Ocean Dynamics (DYN)
-% ================================================================
 \chapter{Ocean Dynamics (DYN)}
 \label{chap:DYN}
@@ -56,7 +53,4 @@
 MISC correspond to "extracting tendency terms" or "vorticity balance"?}
 
-% ================================================================
-% Sea Surface Height evolution & Diagnostics variables
-% ================================================================
 \section{Sea surface height and diagnostic variables ($\eta$, $\zeta$, $\chi$, $w$)}
 \label{sec:DYN_divcur_wzv}
@@ -158,8 +152,4 @@
 (see \autoref{subsec:DOM_Num_Index_vertical}).
 
-
-% ================================================================
-% Coriolis and Advection terms: vector invariant form
-% ================================================================
 \section{Coriolis and advection: vector invariant form}
 \label{sec:DYN_adv_cor_vect}
@@ -182,7 +172,4 @@
 \autoref{chap:LBC}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Vorticity term
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Vorticity term (\textit{dynvor.F90})]{Vorticity term (\protect\mdl{dynvor})}
 \label{subsec:DYN_vor}
@@ -314,5 +301,4 @@
 \end{equation}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -414,8 +400,4 @@
 an option which is only available with a TVD scheme (see \np{ln_traadv_tvd_zts}{ln\_traadv\_tvd\_zts} in \autoref{subsec:TRA_adv_tvd}).
 
-
-% ================================================================
-% Coriolis and Advection : flux form
-% ================================================================
 \section{Coriolis and advection: flux form}
 \label{sec:DYN_adv_cor_flux}
@@ -430,5 +412,4 @@
 At the lateral boundaries either free slip,
 no slip or partial slip boundary conditions are applied following \autoref{chap:LBC}.
-
 
 %--------------------------------------------------------------------------------------------------------------
@@ -562,7 +543,4 @@
 %%%
 
-% ================================================================
-%           Hydrostatic pressure gradient term
-% ================================================================
 \section[Hydrostatic pressure gradient (\textit{dynhpg.F90})]{Hydrostatic pressure gradient (\protect\mdl{dynhpg})}
 \label{sec:DYN_hpg}
@@ -780,7 +758,4 @@
 This option is controlled by  \np{nn_dynhpg_rst}{nn\_dynhpg\_rst}, a namelist parameter.
 
-% ================================================================
-% Surface Pressure Gradient
-% ================================================================
 \section[Surface pressure gradient (\textit{dynspg.F90})]{Surface pressure gradient (\protect\mdl{dynspg})}
 \label{sec:DYN_spg}
@@ -809,5 +784,4 @@
 so that the update of the next velocities is done in module \mdl{dynspg\_flt} and not in \mdl{dynnxt}.
 
-
 The form of the surface pressure gradient term depends on how the user wants to
 handle the fast external gravity waves that are a solution of the analytical equation (\autoref{sec:MB_hor_pg}).
@@ -820,5 +794,4 @@
 The extra term introduced in the filtered method is calculated implicitly, so that a solver is used to compute it.
 As a consequence the update of the $next$ velocities is done in module \mdl{dynspg\_flt} and not in \mdl{dynnxt}.
-
 
 %--------------------------------------------------------------------------------------------------------------
@@ -1092,5 +1065,4 @@
 %>>>>>===============
 
-
 %--------------------------------------------------------------------------------------------------------------
 % Filtered free surface formulation
@@ -1121,7 +1093,4 @@
 It is computed once and for all and applies to all ocean time steps.
 
-% ================================================================
-% Lateral diffusion term
-% ================================================================
 \section[Lateral diffusion term and operators (\textit{dynldf.F90})]{Lateral diffusion term and operators (\protect\mdl{dynldf})}
 \label{sec:DYN_ldf}
@@ -1161,536 +1130,6 @@
 }
 
-% ================================================================
-\subsection[Iso-level laplacian (\forcode{ln_dynldf_lap})]{Iso-level laplacian operator (\protect\np{ln_dynldf_lap}{ln\_dynldf\_lap})}
-\label{subsec:DYN_ldf_lap}
-
-For lateral iso-level diffusion, the discrete operator is:
-\begin{equation}
-  \label{eq:DYN_ldf_lap}
-  \left\{
-    \begin{aligned}
-      D_u^{l{\mathrm {\mathbf U}}} =\frac{1}{e_{1u} }\delta_{i+1/2} \left[ {A_T^{lm}
-          \;\chi } \right]-\frac{1}{e_{2u} {\kern 1pt}e_{3u} }\delta_j \left[
-        {A_f^{lm} \;e_{3f} \zeta } \right] \\ \\
-      D_v^{l{\mathrm {\mathbf U}}} =\frac{1}{e_{2v} }\delta_{j+1/2} \left[ {A_T^{lm}
-          \;\chi } \right]+\frac{1}{e_{1v} {\kern 1pt}e_{3v} }\delta_i \left[
-        {A_f^{lm} \;e_{3f} \zeta } \right]
-    \end{aligned}
-  \right.
-\end{equation}
-
-As explained in \autoref{subsec:MB_ldf},
-this formulation (as the gradient of a divergence and curl of the vorticity) preserves symmetry and
-ensures a complete separation between the vorticity and divergence parts of the momentum diffusion.
-
-%--------------------------------------------------------------------------------------------------------------
-%           Rotated laplacian operator
-%--------------------------------------------------------------------------------------------------------------
-\subsection[Rotated laplacian (\forcode{ln_dynldf_iso})]{Rotated laplacian operator (\protect\np{ln_dynldf_iso}{ln\_dynldf\_iso})}
-\label{subsec:DYN_ldf_iso}
-
-A rotation of the lateral momentum diffusion operator is needed in several cases:
-for iso-neutral diffusion in the $z$-coordinate (\np[=.true.]{ln_dynldf_iso}{ln\_dynldf\_iso}) and
-for either iso-neutral (\np[=.true.]{ln_dynldf_iso}{ln\_dynldf\_iso}) or
-geopotential (\np[=.true.]{ln_dynldf_hor}{ln\_dynldf\_hor}) diffusion in the $s$-coordinate.
-In the partial step case, coordinates are horizontal except at the deepest level and
-no rotation is performed when \np[=.true.]{ln_dynldf_hor}{ln\_dynldf\_hor}.
-The diffusion operator is defined simply as the divergence of down gradient momentum fluxes on
-each momentum component.
-It must be emphasized that this formulation ignores constraints on the stress tensor such as symmetry.
-The resulting discrete representation is:
-\begin{equation}
-  \label{eq:DYN_ldf_iso}
-  \begin{split}
-    D_u^{l\textbf{U}} &= \frac{1}{e_{1u} \, e_{2u} \, e_{3u} }	\\
-    &  \left\{\quad  {\delta_{i+1/2} \left[ {A_T^{lm}  \left(
-              {\frac{e_{2t} \; e_{3t} }{e_{1t} } \,\delta_{i}[u]
-                -e_{2t} \; r_{1t} \,\overline{\overline {\delta_{k+1/2}[u]}}^{\,i,\,k}}
-            \right)} \right]} 	\right. \\
-    & \qquad +\ \delta_j \left[ {A_f^{lm} \left( {\frac{e_{1f}\,e_{3f} }{e_{2f}
-            }\,\delta_{j+1/2} [u] - e_{1f}\, r_{2f}
-            \,\overline{\overline {\delta_{k+1/2} [u]}} ^{\,j+1/2,\,k}}
-        \right)} \right] \\
-    &\qquad +\ \delta_k \left[ {A_{uw}^{lm} \left( {-e_{2u} \, r_{1uw} \,\overline{\overline
-              {\delta_{i+1/2} [u]}}^{\,i+1/2,\,k+1/2} }
-        \right.} \right. \\
-    &  \ \qquad \qquad \qquad \quad\
-    - e_{1u} \, r_{2uw} \,\overline{\overline {\delta_{j+1/2} [u]}} ^{\,j,\,k+1/2} \\
-    & \left. {\left. { \ \qquad \qquad \qquad \ \ \ \left. {\
-                +\frac{e_{1u}\, e_{2u} }{e_{3uw} }\,\left( {r_{1uw}^2+r_{2uw}^2}
-                \right)\,\delta_{k+1/2} [u]} \right)} \right]\;\;\;} \right\} \\ \\
-    D_v^{l\textbf{V}} &= \frac{1}{e_{1v} \, e_{2v} \, e_{3v} } \\
-    &  \left\{\quad  {\delta_{i+1/2} \left[ {A_f^{lm}  \left(
-              {\frac{e_{2f} \; e_{3f} }{e_{1f} } \,\delta_{i+1/2}[v]
-                -e_{2f} \; r_{1f} \,\overline{\overline {\delta_{k+1/2}[v]}}^{\,i+1/2,\,k}}
-            \right)} \right]} 	\right. \\
-    & \qquad +\ \delta_j \left[ {A_T^{lm} \left( {\frac{e_{1t}\,e_{3t} }{e_{2t}
-            }\,\delta_{j} [v] - e_{1t}\, r_{2t}
-            \,\overline{\overline {\delta_{k+1/2} [v]}} ^{\,j,\,k}}
-        \right)} \right] \\
-    & \qquad +\ \delta_k \left[ {A_{vw}^{lm} \left( {-e_{2v} \, r_{1vw} \,\overline{\overline
-              {\delta_{i+1/2} [v]}}^{\,i+1/2,\,k+1/2} }\right.} \right. \\
-    &  \ \qquad \qquad \qquad \quad\
-    - e_{1v} \, r_{2vw} \,\overline{\overline {\delta_{j+1/2} [v]}} ^{\,j+1/2,\,k+1/2} \\
-    & \left. {\left. { \ \qquad \qquad \qquad \ \ \ \left. {\
-                +\frac{e_{1v}\, e_{2v} }{e_{3vw} }\,\left( {r_{1vw}^2+r_{2vw}^2}
-                \right)\,\delta_{k+1/2} [v]} \right)} \right]\;\;\;} \right\}
-  \end{split}
-\end{equation}
-where $r_1$ and $r_2$ are the slopes between the surface along which the diffusion operator acts and
-the surface of computation ($z$- or $s$-surfaces).
-The way these slopes are evaluated is given in the lateral physics chapter (\autoref{chap:LDF}).
-
-%--------------------------------------------------------------------------------------------------------------
-%           Iso-level bilaplacian operator
-%--------------------------------------------------------------------------------------------------------------
-\subsection[Iso-level bilaplacian (\forcode{ln_dynldf_bilap})]{Iso-level bilaplacian operator (\protect\np{ln_dynldf_bilap}{ln\_dynldf\_bilap})}
-\label{subsec:DYN_ldf_bilap}
-
-The lateral fourth order operator formulation on momentum is obtained by applying \autoref{eq:DYN_ldf_lap} twice.
-It requires an additional assumption on boundary conditions:
-the first derivative term normal to the coast depends on the free or no-slip lateral boundary conditions chosen,
-while the third derivative terms normal to the coast are set to zero (see \autoref{chap:LBC}).
-%%%
-\gmcomment{add a remark on the the change in the position of the coefficient}
-%%%
-
-% ================================================================
 %           Vertical diffusion term
-% ================================================================
-\section[Vertical diffusion term (\textit{dynzdf.F90})]{Vertical diffusion term (\protect\mdl{dynzdf})}
-\label{sec:DYN_zdf}
-%----------------------------------------------namzdf------------------------------------------------------
-
-%-------------------------------------------------------------------------------------------------------------
-
-Options are defined through the \nam{zdf}{zdf} namelist variables.
-The large vertical diffusion coefficient found in the surface mixed layer together with high vertical resolution implies that in the case of explicit time stepping there would be too restrictive a constraint on the time step.
-Two time stepping schemes can be used for the vertical diffusion term:
-$(a)$ a forward time differencing scheme
-(\np[=.true.]{ln_zdfexp}{ln\_zdfexp}) using a time splitting technique (\np{nn_zdfexp}{nn\_zdfexp} $>$ 1) or
-$(b)$ a backward (or implicit) time differencing scheme (\np[=.false.]{ln_zdfexp}{ln\_zdfexp})
-(see \autoref{chap:TD}).
-Note that namelist variables \np{ln_zdfexp}{ln\_zdfexp} and \np{nn_zdfexp}{nn\_zdfexp} apply to both tracers and dynamics.
-
-The formulation of the vertical subgrid scale physics is the same whatever the vertical coordinate is.
-The vertical diffusion operators given by \autoref{eq:MB_zdf} take the following semi-discrete space form:
-\[
-  % \label{eq:DYN_zdf}
-  \left\{
-    \begin{aligned}
-      D_u^{vm} &\equiv \frac{1}{e_{3u}} \ \delta_k \left[ \frac{A_{uw}^{vm} }{e_{3uw} }
-        \ \delta_{k+1/2} [\,u\,]         \right]     \\
-      \\
-      D_v^{vm} &\equiv \frac{1}{e_{3v}} \ \delta_k \left[ \frac{A_{vw}^{vm} }{e_{3vw} }
-        \ \delta_{k+1/2} [\,v\,]         \right]
-    \end{aligned}
-  \right.
-\]
-where $A_{uw}^{vm} $ and $A_{vw}^{vm} $ are the vertical eddy viscosity and diffusivity coefficients.
-The way these coefficients are evaluated depends on the vertical physics used (see \autoref{chap:ZDF}).
-
-The surface boundary condition on momentum is the stress exerted by the wind.
-At the surface, the momentum fluxes are prescribed as the boundary condition on
-the vertical turbulent momentum fluxes,
-\begin{equation}
-  \label{eq:DYN_zdf_sbc}
-  \left.{\left( {\frac{A^{vm} }{e_3 }\ \frac{\partial \textbf{U}_h}{\partial k}} \right)} \right|_{z=1}
-  = \frac{1}{\rho_o} \binom{\tau_u}{\tau_v }
-\end{equation}
-where $\left( \tau_u ,\tau_v \right)$ are the two components of the wind stress vector in
-the (\textbf{i},\textbf{j}) coordinate system.
-The high mixing coefficients in the surface mixed layer ensure that the surface wind stress is distributed in
-the vertical over the mixed layer depth.
-If the vertical mixing coefficient is small (when no mixed layer scheme is used)
-the surface stress enters only the top model level, as a body force.
-The surface wind stress is calculated in the surface module routines (SBC, see \autoref{chap:SBC}).
-
-The turbulent flux of momentum at the bottom of the ocean is specified through a bottom friction parameterisation
-(see \autoref{sec:ZDF_drg})
-
-% ================================================================
 % External Forcing
-% ================================================================
-\section{External forcings}
-\label{sec:DYN_forcing}
-
-Besides the surface and bottom stresses (see the above section)
-which are introduced as boundary conditions on the vertical mixing,
-three other forcings may enter the dynamical equations by affecting the surface pressure gradient.
-
-(1) When \np[=.true.]{ln_apr_dyn}{ln\_apr\_dyn} (see \autoref{sec:SBC_apr}),
-the atmospheric pressure is taken into account when computing the surface pressure gradient.
-
-(2) When \np[=.true.]{ln_tide_pot}{ln\_tide\_pot} and \np[=.true.]{ln_tide}{ln\_tide} (see \autoref{sec:SBC_tide}),
-the tidal potential is taken into account when computing the surface pressure gradient.
-
-(3) When \np[=2]{nn_ice_embd}{nn\_ice\_embd} and LIM or CICE is used
-(\ie\ when the sea-ice is embedded in the ocean),
-the snow-ice mass is taken into account when computing the surface pressure gradient.
-
-
-\gmcomment{ missing : the lateral boundary condition !!!   another external forcing
- }
-
-% ================================================================
 % Wetting and drying
-% ================================================================
-\section{Wetting and drying }
-\label{sec:DYN_wetdry}
-
-There are two main options for wetting and drying code (wd):
-(a) an iterative limiter (il) and (b) a directional limiter (dl).
-The directional limiter is based on the scheme developed by \cite{warner.defne.ea_CG13} for RO
-MS
-which was in turn based on ideas developed for POM by \cite{oey_OM06}. The iterative
-limiter is a new scheme.  The iterative limiter is activated by setting $\mathrm{ln\_wd\_il} = \mathrm{.true.}$
-and $\mathrm{ln\_wd\_dl} = \mathrm{.false.}$. The directional limiter is activated
-by setting $\mathrm{ln\_wd\_dl} = \mathrm{.true.}$ and $\mathrm{ln\_wd\_il} = \mathrm{.false.}$.
-
-\begin{listing}
-  \nlst{namwad}
-  \caption{\forcode{&namwad}}
-  \label{lst:namwad}
-\end{listing}
-
-The following terminology is used. The depth of the topography (positive downwards)
-at each $(i,j)$ point is the quantity stored in array $\mathrm{ht\_wd}$ in the \NEMO\ code.
-The height of the free surface (positive upwards) is denoted by $ \mathrm{ssh}$. Given the sign
-conventions used, the water depth, $h$, is the height of the free surface plus the depth of the
-topography (i.e. $\mathrm{ssh} + \mathrm{ht\_wd}$).
-
-Both wd schemes take all points in the domain below a land elevation of $\mathrm{rn\_wdld}$ to be
-covered by water. They require the topography specified with a model
-configuration to have negative depths at points where the land is higher than the
-topography's reference sea-level. The vertical grid in \NEMO\ is normally computed relative to an
-initial state with zero sea surface height elevation.
-The user can choose to compute the vertical grid and heights in the model relative to
-a non-zero reference height for the free surface. This choice affects the calculation of the metrics and depths
-(i.e. the $\mathrm{e3t\_0, ht\_0}$ etc. arrays).
-
-Points where the water depth is less than $\mathrm{rn\_wdmin1}$ are interpreted as ``dry''.
-$\mathrm{rn\_wdmin1}$ is usually chosen to be of order $0.05$m but extreme topographies
-with very steep slopes require larger values for normal choices of time-step. Surface fluxes
-are also switched off for dry cells to prevent freezing, boiling etc. of very thin water layers.
-The fluxes are tappered down using a $\mathrm{tanh}$ weighting function
-to no flux as the dry limit $\mathrm{rn\_wdmin1}$ is approached. Even wet cells can be very shallow.
-The depth at which to start tapering is controlled by the user by setting $\mathrm{rn\_wd\_sbcdep}$.
-The fraction $(<1)$ of sufrace fluxes to use at this depth is set by $\mathrm{rn\_wd\_sbcfra}$.
-
-Both versions of the code have been tested in six test cases provided in the WAD\_TEST\_CASES configuration
-and in ``realistic'' configurations covering parts of the north-west European shelf.
-All these configurations have used pure sigma coordinates. It is expected that
-the wetting and drying code will work in domains with more general s-coordinates provided
-the coordinates are pure sigma in the region where wetting and drying actually occurs.
-
-The next sub-section descrbies the directional limiter and the following sub-section the iterative limiter.
-The final sub-section covers some additional considerations that are relevant to both schemes.
-
-
-%-----------------------------------------------------------------------------------------
-%   Iterative limiters
-%-----------------------------------------------------------------------------------------
-\subsection[Directional limiter (\textit{wet\_dry.F90})]{Directional limiter (\mdl{wet\_dry})}
-\label{subsec:DYN_wd_directional_limiter}
-
-The principal idea of the directional limiter is that
-water should not be allowed to flow out of a dry tracer cell (i.e. one whose water depth is less than \np{rn_wdmin1}{rn\_wdmin1}).
-
-All the changes associated with this option are made to the barotropic solver for the non-linear
-free surface code within dynspg\_ts.
-On each barotropic sub-step the scheme determines the direction of the flow across each face of all the tracer cells
-and sets the flux across the face to zero when the flux is from a dry tracer cell. This prevents cells
-whose depth is rn\_wdmin1 or less from drying out further. The scheme does not force $h$ (the water depth) at tracer cells
-to be at least the minimum depth and hence is able to conserve mass / volume.
-
-The flux across each $u$-face of a tracer cell is multiplied by a factor zuwdmask (an array which depends on ji and jj).
-If the user sets \np[=.false.]{ln_wd_dl_ramp}{ln\_wd\_dl\_ramp} then zuwdmask is 1 when the
-flux is from a cell with water depth greater than \np{rn_wdmin1}{rn\_wdmin1} and 0 otherwise. If the user sets
-\np[=.true.]{ln_wd_dl_ramp}{ln\_wd\_dl\_ramp} the flux across the face is ramped down as the water depth decreases
-from 2 * \np{rn_wdmin1}{rn\_wdmin1} to \np{rn_wdmin1}{rn\_wdmin1}. The use of this ramp reduced grid-scale noise in idealised test cases.
-
-At the point where the flux across a $u$-face is multiplied by zuwdmask , we have chosen
-also to multiply the corresponding velocity on the ``now'' step at that face by zuwdmask. We could have
-chosen not to do that and to allow fairly large velocities to occur in these ``dry'' cells.
-The rationale for setting the velocity to zero is that it is the momentum equations that are being solved
-and the total momentum of the upstream cell (treating it as a finite volume) should be considered
-to be its depth times its velocity. This depth is considered to be zero at ``dry'' $u$-points consistent with its
-treatment in the calculation of the flux of mass across the cell face.
-
-
-\cite{warner.defne.ea_CG13} state that in their scheme the velocity masks at the cell faces for the baroclinic
-timesteps are set to 0 or 1 depending on whether the average of the masks over the barotropic sub-steps is respectively less than
-or greater than 0.5. That scheme does not conserve tracers in integrations started from constant tracer
-fields (tracers independent of $x$, $y$ and $z$). Our scheme conserves constant tracers because
-the velocities used at the tracer cell faces on the baroclinic timesteps are carefully calculated by dynspg\_ts
-to equal their mean value during the barotropic steps. If the user sets \np[=.true.]{ln_wd_dl_bc}{ln\_wd\_dl\_bc}, the
-baroclinic velocities are also multiplied by a suitably weighted average of zuwdmask.
-
-%-----------------------------------------------------------------------------------------
-%   Iterative limiters
-%-----------------------------------------------------------------------------------------
-
-\subsection[Iterative limiter (\textit{wet\_dry.F90})]{Iterative limiter (\mdl{wet\_dry})}
-\label{subsec:DYN_wd_iterative_limiter}
-
-\subsubsection[Iterative flux limiter (\textit{wet\_dry.F90})]{Iterative flux limiter (\mdl{wet\_dry})}
-\label{subsec:DYN_wd_il_spg_limiter}
-
-The iterative limiter modifies the fluxes across the faces of cells that are either already ``dry''
-or may become dry within the next time-step using an iterative method.
-
-The flux limiter for the barotropic flow (devised by Hedong Liu) can be understood as follows:
-
-The continuity equation for the total water depth in a column
-\begin{equation}
-  \label{eq:DYN_wd_continuity}
-  \frac{\partial h}{\partial t} + \mathbf{\nabla.}(h\mathbf{u}) = 0 .
-\end{equation}
-can be written in discrete form  as
-
-\begin{align}
-  \label{eq:DYN_wd_continuity_2}
-  \frac{e_1 e_2}{\Delta t} ( h_{i,j}(t_{n+1}) - h_{i,j}(t_e) )
-  &= - ( \mathrm{flxu}_{i+1,j} - \mathrm{flxu}_{i,j}  + \mathrm{flxv}_{i,j+1} - \mathrm{flxv}_{i,j} ) \\
-  &= \mathrm{zzflx}_{i,j} .
-\end{align}
-
-In the above $h$ is the depth of the water in the column at point $(i,j)$,
-$\mathrm{flxu}_{i+1,j}$ is the flux out of the ``eastern'' face of the cell and
-$\mathrm{flxv}_{i,j+1}$ the flux out of the ``northern'' face of the cell; $t_{n+1}$ is
-the new timestep, $t_e$ is the old timestep (either $t_b$ or $t_n$) and $ \Delta t =
-t_{n+1} - t_e$; $e_1 e_2$ is the area of the tracer cells centred at $(i,j)$ and
-$\mathrm{zzflx}$ is the sum of the fluxes through all the faces.
-
-The flux limiter splits the flux $\mathrm{zzflx}$ into fluxes that are out of the cell
-(zzflxp) and fluxes that are into the cell (zzflxn).  Clearly
-
-\begin{equation}
-  \label{eq:DYN_wd_zzflx_p_n_1}
-  \mathrm{zzflx}_{i,j} = \mathrm{zzflxp}_{i,j} + \mathrm{zzflxn}_{i,j} .
-\end{equation}
-
-The flux limiter iteratively adjusts the fluxes $\mathrm{flxu}$ and $\mathrm{flxv}$ until
-none of the cells will ``dry out''. To be precise the fluxes are limited until none of the
-cells has water depth less than $\mathrm{rn\_wdmin1}$ on step $n+1$.
-
-Let the fluxes on the $m$th iteration step be denoted by $\mathrm{flxu}^{(m)}$ and
-$\mathrm{flxv}^{(m)}$.  Then the adjustment is achieved by seeking a set of coefficients,
-$\mathrm{zcoef}_{i,j}^{(m)}$ such that:
-
-\begin{equation}
-  \label{eq:DYN_wd_continuity_coef}
-  \begin{split}
-    \mathrm{zzflxp}^{(m)}_{i,j} =& \mathrm{zcoef}_{i,j}^{(m)} \mathrm{zzflxp}^{(0)}_{i,j} \\
-    \mathrm{zzflxn}^{(m)}_{i,j} =& \mathrm{zcoef}_{i,j}^{(m)} \mathrm{zzflxn}^{(0)}_{i,j}
-  \end{split}
-\end{equation}
-
-where the coefficients are $1.0$ generally but can vary between $0.0$ and $1.0$ around
-cells that would otherwise dry.
-
-The iteration is initialised by setting
-
-\begin{equation}
-  \label{eq:DYN_wd_zzflx_initial}
-  \mathrm{zzflxp^{(0)}}_{i,j} = \mathrm{zzflxp}_{i,j} , \quad  \mathrm{zzflxn^{(0)}}_{i,j} = \mathrm{zzflxn}_{i,j} .
-\end{equation}
-
-The fluxes out of cell $(i,j)$ are updated at the $m+1$th iteration if the depth of the
-cell on timestep $t_e$, namely $h_{i,j}(t_e)$, is less than the total flux out of the cell
-times the timestep divided by the cell area. Using (\autoref{eq:DYN_wd_continuity_2}) this
-condition is
-
-\begin{equation}
-  \label{eq:DYN_wd_continuity_if}
-  h_{i,j}(t_e)  - \mathrm{rn\_wdmin1} <  \frac{\Delta t}{e_1 e_2} ( \mathrm{zzflxp}^{(m)}_{i,j} + \mathrm{zzflxn}^{(m)}_{i,j} ) .
-\end{equation}
-
-Rearranging (\autoref{eq:DYN_wd_continuity_if}) we can obtain an expression for the maximum
-outward flux that can be allowed and still maintain the minimum wet depth:
-
-\begin{equation}
-  \label{eq:DYN_wd_max_flux}
-  \begin{split}
-    \mathrm{zzflxp}^{(m+1)}_{i,j} = \Big[ (h_{i,j}(t_e) & - \mathrm{rn\_wdmin1} - \mathrm{rn\_wdmin2})  \frac{e_1 e_2}{\Delta t} \phantom{]} \\
-    \phantom{[} & -  \mathrm{zzflxn}^{(m)}_{i,j} \Big]
-  \end{split}
-\end{equation}
-
-Note a small tolerance ($\mathrm{rn\_wdmin2}$) has been introduced here {\itshape [Q: Why is
-this necessary/desirable?]}. Substituting from (\autoref{eq:DYN_wd_continuity_coef}) gives an
-expression for the coefficient needed to multiply the outward flux at this cell in order
-to avoid drying.
-
-\begin{equation}
-  \label{eq:DYN_wd_continuity_nxtcoef}
-  \begin{split}
-    \mathrm{zcoef}^{(m+1)}_{i,j} = \Big[ (h_{i,j}(t_e) & - \mathrm{rn\_wdmin1} - \mathrm{rn\_wdmin2})  \frac{e_1 e_2}{\Delta t} \phantom{]} \\
-    \phantom{[} & -  \mathrm{zzflxn}^{(m)}_{i,j} \Big] \frac{1}{ \mathrm{zzflxp}^{(0)}_{i,j} }
-  \end{split}
-\end{equation}
-
-Only the outward flux components are altered but, of course, outward fluxes from one cell
-are inward fluxes to adjacent cells and the balance in these cells may need subsequent
-adjustment; hence the iterative nature of this scheme.  Note, for example, that the flux
-across the ``eastern'' face of the $(i,j)$th cell is only updated at the $m+1$th iteration
-if that flux at the $m$th iteration is out of the $(i,j)$th cell. If that is the case then
-the flux across that face is into the $(i+1,j)$ cell and that flux will not be updated by
-the calculation for the $(i+1,j)$th cell. In this sense the updates to the fluxes across
-the faces of the cells do not ``compete'' (they do not over-write each other) and one
-would expect the scheme to converge relatively quickly. The scheme is flux based so
-conserves mass. It also conserves constant tracers for the same reason that the
-directional limiter does.
-
-
-%----------------------------------------------------------------------------------------
-%      Surface pressure gradients
-%----------------------------------------------------------------------------------------
-\subsubsection[Modification of surface pressure gradients (\textit{dynhpg.F90})]{Modification of surface pressure gradients (\mdl{dynhpg})}
-\label{subsec:DYN_wd_il_spg}
-
-At ``dry'' points the water depth is usually close to $\mathrm{rn\_wdmin1}$. If the
-topography is sloping at these points the sea-surface will have a similar slope and there
-will hence be very large horizontal pressure gradients at these points. The WAD modifies
-the magnitude but not the sign of the surface pressure gradients (zhpi and zhpj) at such
-points by mulitplying them by positive factors (zcpx and zcpy respectively) that lie
-between $0$ and $1$.
-
-We describe how the scheme works for the ``eastward'' pressure gradient, zhpi, calculated
-at the $(i,j)$th $u$-point. The scheme uses the ht\_wd depths and surface heights at the
-neighbouring $(i+1,j)$ and $(i,j)$ tracer points.  zcpx is calculated using two logicals
-variables, $\mathrm{ll\_tmp1}$ and $\mathrm{ll\_tmp2}$ which are evaluated for each grid
-column.  The three possible combinations are illustrated in \autoref{fig:DYN_WAD_dynhpg}.
-
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-\begin{figure}[!ht]
-  \centering
-  \includegraphics[width=0.66\textwidth]{Fig_WAD_dynhpg}
-  \caption[Combinations controlling the limiting of the horizontal pressure gradient in
-  wetting and drying regimes]{
-    Three possible combinations of the logical variables controlling the
-    limiting of the horizontal pressure gradient in wetting and drying regimes}
-  \label{fig:DYN_WAD_dynhpg}
-\end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-The first logical, $\mathrm{ll\_tmp1}$, is set to true if and only if the water depth at
-both neighbouring points is greater than $\mathrm{rn\_wdmin1} + \mathrm{rn\_wdmin2}$ and
-the minimum height of the sea surface at the two points is greater than the maximum height
-of the topography at the two points:
-
-\begin{equation}
-  \label{eq:DYN_ll_tmp1}
-  \begin{split}
-    \mathrm{ll\_tmp1}  = & \mathrm{MIN(sshn(ji,jj), sshn(ji+1,jj))} > \\
-                     & \quad \mathrm{MAX(-ht\_wd(ji,jj), -ht\_wd(ji+1,jj))\  .and.} \\
-                     & \mathrm{MAX(sshn(ji,jj) + ht\_wd(ji,jj),} \\
-                     & \mathrm{\phantom{MAX(}sshn(ji+1,jj) + ht\_wd(ji+1,jj))} >\\
-                     & \quad\quad\mathrm{rn\_wdmin1 + rn\_wdmin2 }
-  \end{split}
-\end{equation}
-
-The second logical, $\mathrm{ll\_tmp2}$, is set to true if and only if the maximum height
-of the sea surface at the two points is greater than the maximum height of the topography
-at the two points plus $\mathrm{rn\_wdmin1} + \mathrm{rn\_wdmin2}$
-
-\begin{equation}
-  \label{eq:DYN_ll_tmp2}
-  \begin{split}
-    \mathrm{ ll\_tmp2 } = & \mathrm{( ABS( sshn(ji,jj) - sshn(ji+1,jj) ) > 1.E-12 )\ .AND.}\\
-    & \mathrm{( MAX(sshn(ji,jj), sshn(ji+1,jj)) > } \\
-    & \mathrm{\phantom{(} MAX(-ht\_wd(ji,jj), -ht\_wd(ji+1,jj)) + rn\_wdmin1 + rn\_wdmin2}) .
-  \end{split}
-\end{equation}
-
-If $\mathrm{ll\_tmp1}$ is true then the surface pressure gradient, zhpi at the $(i,j)$
-point is unmodified. If both logicals are false zhpi is set to zero.
-
-If $\mathrm{ll\_tmp1}$ is true and $\mathrm{ll\_tmp2}$ is false then the surface pressure
-gradient is multiplied through by zcpx which is the absolute value of the difference in
-the water depths at the two points divided by the difference in the surface heights at the
-two points. Thus the sign of the sea surface height gradient is retained but the magnitude
-of the pressure force is determined by the difference in water depths rather than the
-difference in surface height between the two points. Note that dividing by the difference
-between the sea surface heights can be problematic if the heights approach parity. An
-additional condition is applied to $\mathrm{ ll\_tmp2 }$ to ensure it is .false. in such
-conditions.
-
-\subsubsection[Additional considerations (\textit{usrdef\_zgr.F90})]{Additional considerations (\mdl{usrdef\_zgr})}
-\label{subsec:DYN_WAD_additional}
-
-In the very shallow water where wetting and drying occurs the parametrisation of
-bottom drag is clearly very important. In order to promote stability
-it is sometimes useful to calculate the bottom drag using an implicit time-stepping approach.
-
-Suitable specifcation of the surface heat flux in wetting and drying domains in forced and
-coupled simulations needs further consideration. In order to prevent freezing or boiling
-in uncoupled integrations the net surface heat fluxes need to be appropriately limited.
-
-%----------------------------------------------------------------------------------------
-%      The WAD test cases
-%----------------------------------------------------------------------------------------
-\subsection[The WAD test cases (\textit{usrdef\_zgr.F90})]{The WAD test cases (\mdl{usrdef\_zgr})}
-\label{subsec:DYN_WAD_test_cases}
-
-See the WAD tests MY\_DOC documention for details of the WAD test cases.
-
-
-
-% ================================================================
 % Time evolution term
-% ================================================================
-\section[Time evolution term (\textit{dynnxt.F90})]{Time evolution term (\protect\mdl{dynnxt})}
-\label{sec:DYN_nxt}
-
-%----------------------------------------------namdom----------------------------------------------------
-
-%-------------------------------------------------------------------------------------------------------------
-
-Options are defined through the \nam{dom}{dom} namelist variables.
-The general framework for dynamics time stepping is a leap-frog scheme,
-\ie\ a three level centred time scheme associated with an Asselin time filter (cf. \autoref{chap:TD}).
-The scheme is applied to the velocity, except when
-using the flux form of momentum advection (cf. \autoref{sec:DYN_adv_cor_flux})
-in the variable volume case (\texttt{vvl?} defined),
-where it has to be applied to the thickness weighted velocity (see \autoref{sec:SCOORD_momentum})
-
-$\bullet$ vector invariant form or linear free surface
-(\np[=.true.]{ln_dynhpg_vec}{ln\_dynhpg\_vec} ; \texttt{vvl?} not defined):
-\[
-  % \label{eq:DYN_nxt_vec}
-  \left\{
-    \begin{aligned}
-      &u^{t+\rdt} = u_f^{t-\rdt} + 2\rdt  \ \text{RHS}_u^t  	\\
-      &u_f^t \;\quad = u^t+\gamma \,\left[ {u_f^{t-\rdt} -2u^t+u^{t+\rdt}} \right]
-    \end{aligned}
-  \right.
-\]
-
-$\bullet$ flux form and nonlinear free surface
-(\np[=.false.]{ln_dynhpg_vec}{ln\_dynhpg\_vec} ; \texttt{vvl?} defined):
-\[
-  % \label{eq:DYN_nxt_flux}
-  \left\{
-    \begin{aligned}
-      &\left(e_{3u}\,u\right)^{t+\rdt} = \left(e_{3u}\,u\right)_f^{t-\rdt} + 2\rdt \; e_{3u} \;\text{RHS}_u^t  	\\
-      &\left(e_{3u}\,u\right)_f^t \;\quad = \left(e_{3u}\,u\right)^t
-      +\gamma \,\left[ {\left(e_{3u}\,u\right)_f^{t-\rdt} -2\left(e_{3u}\,u\right)^t+\left(e_{3u}\,u\right)^{t+\rdt}} \right]
-    \end{aligned}
-  \right.
-\]
-where RHS is the right hand side of the momentum equation,
-the subscript $f$ denotes filtered values and $\gamma$ is the Asselin coefficient.
-$\gamma$ is initialized as \np{nn_atfp}{nn\_atfp} (namelist parameter).
-Its default value is \np[=10.e-3]{nn_atfp}{nn\_atfp}.
-In both cases, the modified Asselin filter is not applied since perfect conservation is not an issue for
-the momentum equations.
-
-Note that with the filtered free surface,
-the update of the \textit{after} velocities is done in the \mdl{dynsp\_flt} module,
-and only array swapping and Asselin filtering is done in \mdl{dynnxt}.
-
-\onlyinsubfile{\input{../../global/epilogue}}
-
-\end{document}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LBC.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LBC.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LBC.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter — Lateral Boundary Condition (LBC)
-% ================================================================
 \chapter{Lateral Boundary Condition (LBC)}
 \label{chap:LBC}
@@ -10,11 +7,6 @@
 \chaptertoc
 
-\newpage
-
 %gm% add here introduction to this chapter
 
-% ================================================================
-% Boundary Condition at the Coast
-% ================================================================
 \section[Boundary condition at the coast (\forcode{rn_shlat})]{Boundary condition at the coast (\protect\np{rn_shlat}{rn\_shlat})}
 \label{sec:LBC_coast}
@@ -68,5 +60,4 @@
 (normal velocity $u$ remains zero at the coast) (\autoref{fig:LBC_uv}).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -77,5 +68,4 @@
   \label{fig:LBC_uv}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 For momentum the situation is a bit more complex as two boundary conditions must be provided along the coast
@@ -95,5 +85,4 @@
 These are:
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!p]
   \centering
@@ -108,9 +97,8 @@
   \label{fig:LBC_shlat}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 \begin{description}
 
-\item[free-slip boundary condition ({\np[=0]{rn_shlat}{rn\_shlat}})] the tangential velocity at
+\item [free-slip boundary condition ({\np[=0]{rn_shlat}{rn\_shlat}})] the tangential velocity at
   the coastline is equal to the offshore velocity,
   \ie\ the normal derivative of the tangential velocity is zero at the coast,
@@ -118,5 +106,5 @@
   (\autoref{fig:LBC_shlat}-a).
 
-\item[no-slip boundary condition ({\np[=2]{rn_shlat}{rn\_shlat}})] the tangential velocity vanishes at the coastline.
+\item [no-slip boundary condition ({\np[=2]{rn_shlat}{rn\_shlat}})] the tangential velocity vanishes at the coastline.
   Assuming that the tangential velocity decreases linearly from
   the closest ocean velocity grid point to the coastline,
@@ -139,10 +127,10 @@
   \]
 
-\item["partial" free-slip boundary condition (0$<$\np{rn_shlat}{rn\_shlat}$<$2)] the tangential velocity at
+\item ["partial" free-slip boundary condition (0$<$\np{rn_shlat}{rn\_shlat}$<$2)] the tangential velocity at
   the coastline is smaller than the offshore velocity, \ie\ there is a lateral friction but
   not strong enough to make the tangential velocity at the coast vanish (\autoref{fig:LBC_shlat}-c).
   This can be selected by providing a value of mask$_{f}$ strictly inbetween $0$ and $2$.
 
-\item["strong" no-slip boundary condition (2$<$\np{rn_shlat}{rn\_shlat})] the viscous boundary layer is assumed to
+\item ["strong" no-slip boundary condition (2$<$\np{rn_shlat}{rn\_shlat})] the viscous boundary layer is assumed to
   be smaller than half the grid size (\autoref{fig:LBC_shlat}-d).
   The friction is thus larger than in the no-slip case.
@@ -154,8 +142,4 @@
 it is only applied next to the coast where the minimum water depth can be quite shallow.
 
-
-% ================================================================
-% Boundary Condition around the Model Domain
-% ================================================================
 \section[Model domain boundary condition (\forcode{jperio})]{Model domain boundary condition (\protect\jp{jperio})}
 \label{sec:LBC_jperio}
@@ -166,7 +150,4 @@
 The north-fold boundary condition is associated with the 3-pole ORCA mesh.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Closed, cyclic (\jp{jperio}\forcode{ = 0..2})
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Closed, cyclic (\forcode{=0,1,2,7})]{Closed, cyclic (\protect\jp{jperio}\forcode{=0,1,2,7})}
 \label{subsec:LBC_jperio012}
@@ -182,9 +163,9 @@
 \begin{description}
 
-\item[For closed boundary (\jp{jperio}\forcode{=0})],
+\item [For closed boundary (\jp{jperio}\forcode{=0})],
   solid walls are imposed at all model boundaries:
   first and last rows and columns are set to zero.
 
-\item[For cyclic east-west boundary (\jp{jperio}\forcode{=1})],
+\item [For cyclic east-west boundary (\jp{jperio}\forcode{=1})],
   first and last rows are set to zero (closed) whilst the first column is set to
   the value of the last-but-one column and the last column to the value of the second one
@@ -192,5 +173,5 @@
   Whatever flows out of the eastern (western) end of the basin enters the western (eastern) end.
 
-\item[For cyclic north-south boundary (\jp{jperio}\forcode{=2})],
+\item [For cyclic north-south boundary (\jp{jperio}\forcode{=2})],
   first and last columns are set to zero (closed) whilst the first row is set to
   the value of the last-but-one row and the last row to the value of the second one
@@ -198,9 +179,8 @@
   Whatever flows out of the northern (southern) end of the basin enters the southern (northern) end.
 
-\item[Bi-cyclic east-west and north-south boundary (\jp{jperio}\forcode{=7})] combines cases 1 and 2.
+\item [Bi-cyclic east-west and north-south boundary (\jp{jperio}\forcode{=7})] combines cases 1 and 2.
 
 \end{description}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -210,9 +190,5 @@
   \label{fig:LBC_jperio}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        North fold (\textit{jperio = 3 }to $6)$
-% -------------------------------------------------------------------------------------------------------------
+
 \subsection[North-fold (\forcode{=3,6})]{North-fold (\protect\jp{jperio}\forcode{=3,6})}
 \label{subsec:LBC_north_fold}
@@ -224,5 +200,4 @@
 Further information can be found in \mdl{lbcnfd} module which applies the north fold boundary condition.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -234,9 +209,5 @@
   \label{fig:LBC_North_Fold_T}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% ====================================================================
-% Exchange with neighbouring processors
-% ====================================================================
+
 \section[Exchange with neighbouring processors (\textit{lbclnk.F90}, \textit{lib\_mpp.F90})]{Exchange with neighbouring processors (\protect\mdl{lbclnk}, \protect\mdl{lib\_mpp})}
 \label{sec:LBC_mpp}
@@ -284,5 +255,4 @@
 many communications during 1 time step of the model.\\
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -291,5 +261,4 @@
   \label{fig:LBC_mpp}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 In \NEMO, the splitting is regular and arithmetic.
@@ -339,5 +308,4 @@
 When land processors are eliminated, the value corresponding to these locations in the model output files is undefined. \np{ln_mskland}{ln\_mskland} must be activated in order avoid Not a Number values in output files. Note that it is better to not eliminate land processors when creating a meshmask file (\ie\ when setting a non-zero value to \np{nn_msh}{nn\_msh}).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -352,10 +320,5 @@
   \label{fig:LBC_mppini2}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-
-% ====================================================================
-% Unstructured open boundaries BDY
-% ====================================================================
+
 \section{Unstructured open boundary conditions (BDY)}
 \label{sec:LBC_bdy}
@@ -411,12 +374,12 @@
 
 \begin{description}
-\item[\forcode{'none'}:] No boundary condition applied.
+\item [\forcode{'none'}:] No boundary condition applied.
   So the solution will ``see'' the land points around the edge of the edge of the domain.
-\item[\forcode{'specified'}:] Specified boundary condition applied (only available for baroclinic velocity and tracer variables).
-\item[\forcode{'neumann'}:] Value at the boundary are duplicated (No gradient). Only available for baroclinic velocity and tracer variables.
-\item[\forcode{'frs'}:] Flow Relaxation Scheme (FRS) available for all variables.
-\item[\forcode{'Orlanski'}:] Orlanski radiation scheme (fully oblique) for barotropic, baroclinic and tracer variables.
-\item[\forcode{'Orlanski_npo'}:] Orlanski radiation scheme for barotropic, baroclinic and tracer variables.
-\item[\forcode{'flather'}:] Flather radiation scheme for the barotropic variables only.
+\item [\forcode{'specified'}:] Specified boundary condition applied (only available for baroclinic velocity and tracer variables).
+\item [\forcode{'neumann'}:] Value at the boundary are duplicated (No gradient). Only available for baroclinic velocity and tracer variables.
+\item [\forcode{'frs'}:] Flow Relaxation Scheme (FRS) available for all variables.
+\item [\forcode{'Orlanski'}:] Orlanski radiation scheme (fully oblique) for barotropic, baroclinic and tracer variables.
+\item [\forcode{'Orlanski_npo'}:] Orlanski radiation scheme for barotropic, baroclinic and tracer variables.
+\item [\forcode{'flather'}:] Flather radiation scheme for the barotropic variables only.
 \end{description}
 
@@ -620,5 +583,4 @@
 Only one mask file is used even if multiple boundary sets are defined.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -627,5 +589,4 @@
   \label{fig:LBC_bdy_geom}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 %----------------------------------------------
@@ -657,5 +618,4 @@
 will re-order the data in old BDY data files.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -665,5 +625,4 @@
   \label{fig:LBC_nc_header}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 %----------------------------------------------
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LDF.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LDF.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_LDF.tex	(revision 11596)
@@ -3,13 +3,8 @@
 \begin{document}
 
-% ================================================================
-% Chapter Lateral Ocean Physics (LDF)
-% ================================================================
 \chapter{Lateral Ocean Physics (LDF)}
 \label{chap:LDF}
 
 \chaptertoc
-
-\newpage
 
 The lateral physics terms in the momentum and tracer equations have been described in \autoref{eq:MB_zdf} and
@@ -31,7 +26,4 @@
 %--------------------------------------------------------------------------------------------------------------
 
-% ================================================================
-% Lateral Mixing Operator
-% ================================================================
 \section[Lateral mixing operators]{Lateral mixing operators}
 \label{sec:LDF_op}
@@ -55,7 +47,4 @@
 We stress again that from \NEMO\ 4, the simultaneous use Laplacian and Bilaplacian operators is not allowed.
 
-% ================================================================
-% Direction of lateral Mixing
-% ================================================================
 \section[Direction of lateral mixing (\textit{ldfslp.F90})]{Direction of lateral mixing (\protect\mdl{ldfslp})}
 \label{sec:LDF_slp}
@@ -151,5 +140,5 @@
 \begin{description}
 
-\item[$z$-coordinate with full step: ]
+\item [$z$-coordinate with full step: ]
   in \autoref{eq:LDF_slp_iso} the densities appearing in the $i$ and $j$ derivatives  are taken at the same depth,
   thus the $in situ$ density can be used.
@@ -158,9 +147,9 @@
   (see \autoref{subsec:TRA_bn2}).
 
-\item[$z$-coordinate with partial step: ]
+\item [$z$-coordinate with partial step: ]
   this case is identical to the full step case except that at partial step level,
   the \emph{horizontal} density gradient is evaluated as described in \autoref{sec:TRA_zpshde}.
 
-\item[$s$- or hybrid $s$-$z$- coordinate: ]
+\item [$s$- or hybrid $s$-$z$- coordinate: ]
   in the current release of \NEMO, iso-neutral mixing is only employed for $s$-coordinates if
   the Griffies scheme is used (\np[=.true.]{ln_traldf_triad}{ln\_traldf\_triad};
@@ -234,5 +223,4 @@
 contrary to the \citet{griffies.gnanadesikan.ea_JPO98} operator which has that property.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -241,5 +229,4 @@
   \label{fig:LDF_ZDF1}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 %There are three additional questions about the slope calculation.
@@ -249,6 +236,4 @@
 
 %from griffies: chapter 13.1....
-
-
 
 % In addition and also for numerical stability reasons \citep{cox_OM87, griffies_bk04},
@@ -264,5 +249,4 @@
 \colorbox{yellow}{The way slopes are tapered has be checked. Not sure that this is still what is actually done.}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -286,5 +270,4 @@
   \label{fig:LDF_eiv_slp}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 \colorbox{yellow}{add here a discussion about the flattening of the slopes, vs tapering the coefficient.}
@@ -316,8 +299,4 @@
 (see \autoref{sec:LBC_coast}).
 
-
-% ================================================================
-% Lateral Mixing Coefficients
-% ================================================================
 \section[Lateral mixing coefficient (\forcode{nn_aht_ijk_t} \& \forcode{nn_ahm_ijk_t})]{Lateral mixing coefficient (\protect\np{nn_aht_ijk_t}{nn\_aht\_ijk\_t} \& \protect\np{nn_ahm_ijk_t}{nn\_ahm\_ijk\_t})}
 \label{sec:LDF_coef}
@@ -467,7 +446,4 @@
 (\autoref{sec:INVARIANTS_dynldf_properties}).
 
-% ================================================================
-% Eddy Induced Mixing
-% ================================================================
 \section[Eddy induced velocity (\forcode{ln_ldfeiv})]{Eddy induced velocity (\protect\np{ln_ldfeiv}{ln\_ldfeiv})}
 
@@ -483,5 +459,4 @@
 
 %--------------------------------------------------------------------------------------------------------------
-
 
 %%gm  from Triad appendix  : to be incorporated....
@@ -538,7 +513,4 @@
 In case of setting \np[=.true.]{ln_traldf_triad}{ln\_traldf\_triad}, a skew form of the eddy induced advective fluxes is used, which is described in \autoref{apdx:TRIADS}.
 
-% ================================================================
-% Mixed layer eddies
-% ================================================================
 \section[Mixed layer eddies (\forcode{ln_mle})]{Mixed layer eddies (\protect\np{ln_mle}{ln\_mle})}
 \label{sec:LDF_mle}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_OBS.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_OBS.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_OBS.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter observation operator (OBS)
-% ================================================================
 \chapter{Observation and Model Comparison (OBS)}
 \label{chap:OBS}
@@ -22,6 +19,4 @@
 \end{figure}
 
-\newpage
-
 The observation and model comparison code, the observation operator (OBS), reads in observation files
 (profile temperature and salinity, sea surface temperature, sea level anomaly, sea ice concentration, and velocity) and calculates an interpolated model equivalent value at the observation location and nearest model time step.
@@ -60,7 +55,4 @@
 In \autoref{sec:OBS_obsutils} we describe some utilities to help work with the files produced by the OBS code.
 
-% ================================================================
-% Example
-% ================================================================
 \section{Running the observation operator code example}
 \label{sec:OBS_example}
@@ -610,5 +602,5 @@
 \begin{enumerate}
 
-\item[1.] {\bfseries Great-Circle distance-weighted interpolation.}
+\item [1.] {\bfseries Great-Circle distance-weighted interpolation.}
   The weights are computed as a function of the great-circle distance $s(P, \cdot)$ between $P$ and
   the model grid points $A$, $B$ etc.
@@ -655,5 +647,5 @@
    \end{alignat*}
 
-\item[2.] {\bfseries Great-Circle distance-weighted interpolation with small angle approximation.}
+\item [2.] {\bfseries Great-Circle distance-weighted interpolation with small angle approximation.}
   Similar to the previous interpolation but with the distance $s$ computed as
   \begin{alignat*}{2}
@@ -665,8 +657,8 @@
   where $M$ corresponds to $A$, $B$, $C$ or $D$.
 
-\item[3.] {\bfseries Bilinear interpolation for a regular spaced grid.}
+\item [3.] {\bfseries Bilinear interpolation for a regular spaced grid.}
   The interpolation is split into two 1D interpolations in the longitude and latitude directions, respectively.
 
-\item[4.] {\bfseries Bilinear remapping interpolation for a general grid.}
+\item [4.] {\bfseries Bilinear remapping interpolation for a general grid.}
   An iterative scheme that involves first mapping a quadrilateral cell into
   a cell with coordinates (0,0), (1,0), (0,1) and (1,1).
@@ -697,5 +689,4 @@
 \autoref{fig:OBS_avgrec} and~\autoref{fig:OBS_avgrad}.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -708,5 +699,4 @@
 % >>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -717,6 +707,4 @@
   \label{fig:OBS_avgrad}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
 
 \subsection{Grid search}
@@ -786,5 +774,4 @@
 \subsubsection{Geographical distribution of observations among processors}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -795,5 +782,4 @@
   \label{fig:OBS_local}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 This is the simplest option in which the observations are distributed according to
@@ -814,5 +800,4 @@
 \subsubsection{Round-robin distribution of observations among processors}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -823,5 +808,4 @@
   \label{fig:OBS_global}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 An alternative approach is to distribute the observations equally among processors and
@@ -843,10 +827,4 @@
 
 For profile observation types we do both vertical and horizontal interpolation. \NEMO\ has a generalised vertical coordinate system this means the vertical level depths can vary with location. Therefore, it is necessary first to perform vertical interpolation of the model value to the observation depths for each of the four surrounding grid points. After this the model values, at these points, at the observation depth, are horizontally interpolated to the observation location.
-
-\newpage
-
-% ================================================================
-% Standalone observation operator documentation
-% ================================================================
 
 %\usepackage{framed}
@@ -957,6 +935,4 @@
 climatologies with the same set of observations.
 This approach is referred to as \emph{Class 4} since it is the fourth metric defined by the GODAE intercomparison project. This requires multiple runs of the SAO and running an additional utility (not currently in the \NEMO\ repository) to combine the feedback files into one class 4 file.
-
-\newpage
 
 \section{Observation utilities}
@@ -1146,5 +1122,4 @@
 The rightmost group of buttons will print the plot window as a postscript, save it as png, or exit from dataplot.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -1153,10 +1128,8 @@
   \label{fig:OBS_dataplotmain}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 If a profile point is clicked with the mouse button a plot of the observation and background values as
 a function of depth (\autoref{fig:OBS_dataplotprofile}).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}
   \centering
@@ -1166,5 +1139,4 @@
   \label{fig:OBS_dataplotprofile}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 \onlyinsubfile{\input{../../global/epilogue}}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_SBC.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_SBC.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_SBC.tex	(revision 11596)
@@ -3,13 +3,8 @@
 \begin{document}
 
-% ================================================================
-% Chapter —— Surface Boundary Condition (SBC, SAS, ISF, ICB)
-% ================================================================
 \chapter{Surface Boundary Condition (SBC, SAS, ISF, ICB)}
 \label{chap:SBC}
 
 \chaptertoc
-
-\newpage
 
 %---------------------------------------namsbc--------------------------------------------------
@@ -25,14 +20,9 @@
 
 \begin{itemize}
-\item
-  the two components of the surface ocean stress $\left( {\tau_u \;,\;\tau_v} \right)$
-\item
-  the incoming solar and non solar heat fluxes $\left( {Q_{ns} \;,\;Q_{sr} } \right)$
-\item
-  the surface freshwater budget $\left( {\textit{emp}} \right)$
-\item
-  the surface salt flux associated with freezing/melting of seawater $\left( {\textit{sfx}} \right)$
-\item
-  the atmospheric pressure at the ocean surface $\left( p_a \right)$
+\item the two components of the surface ocean stress $\left( {\tau_u \;,\;\tau_v} \right)$
+\item the incoming solar and non solar heat fluxes $\left( {Q_{ns} \;,\;Q_{sr} } \right)$
+\item the surface freshwater budget $\left( {\textit{emp}} \right)$
+\item the surface salt flux associated with freezing/melting of seawater $\left( {\textit{sfx}} \right)$
+\item the atmospheric pressure at the ocean surface $\left( p_a \right)$
 \end{itemize}
 
@@ -41,13 +31,9 @@
 
 \begin{itemize}
-\item
-  a bulk formulation (\np[=.true.]{ln_blk}{ln\_blk} with four possible bulk algorithms),
-\item
-  a flux formulation (\np[=.true.]{ln_flx}{ln\_flx}),
-\item
-  a coupled or mixed forced/coupled formulation (exchanges with a atmospheric model via the OASIS coupler),
+\item a bulk formulation (\np[=.true.]{ln_blk}{ln\_blk} with four possible bulk algorithms),
+\item a flux formulation (\np[=.true.]{ln_flx}{ln\_flx}),
+\item a coupled or mixed forced/coupled formulation (exchanges with a atmospheric model via the OASIS coupler),
 (\np{ln_cpl}{ln\_cpl} or \np[=.true.]{ln_mixcpl}{ln\_mixcpl}),
-\item
-  a user defined formulation (\np[=.true.]{ln_usr}{ln\_usr}).
+\item a user defined formulation (\np[=.true.]{ln_usr}{ln\_usr}).
 \end{itemize}
 
@@ -66,45 +52,27 @@
 
 \begin{itemize}
-\item
-  the rotation of vector components supplied relative to an east-north coordinate system onto
+\item the rotation of vector components supplied relative to an east-north coordinate system onto
   the local grid directions in the model,
-\item
-  the use of a land/sea mask for input fields (\np[=.true.]{nn_lsm}{nn\_lsm}),
-\item
-  the addition of a surface restoring term to observed SST and/or SSS (\np[=.true.]{ln_ssr}{ln\_ssr}),
-\item
-  the modification of fluxes below ice-covered areas (using climatological ice-cover or a sea-ice model)
+\item the use of a land/sea mask for input fields (\np[=.true.]{nn_lsm}{nn\_lsm}),
+\item the addition of a surface restoring term to observed SST and/or SSS (\np[=.true.]{ln_ssr}{ln\_ssr}),
+\item the modification of fluxes below ice-covered areas (using climatological ice-cover or a sea-ice model)
   (\np[=0..3]{nn_ice}{nn\_ice}),
-\item
-  the addition of river runoffs as surface freshwater fluxes or lateral inflow (\np[=.true.]{ln_rnf}{ln\_rnf}),
-\item
-  the addition of ice-shelf melting as lateral inflow (parameterisation) or
+\item the addition of river runoffs as surface freshwater fluxes or lateral inflow (\np[=.true.]{ln_rnf}{ln\_rnf}),
+\item the addition of ice-shelf melting as lateral inflow (parameterisation) or
   as fluxes applied at the land-ice ocean interface (\np[=.true.]{ln_isf}{ln\_isf}),
-\item
-  the addition of a freshwater flux adjustment in order to avoid a mean sea-level drift
+\item the addition of a freshwater flux adjustment in order to avoid a mean sea-level drift
   (\np[=0..2]{nn_fwb}{nn\_fwb}),
-\item
-  the transformation of the solar radiation (if provided as daily mean) into an analytical diurnal cycle
+\item the transformation of the solar radiation (if provided as daily mean) into an analytical diurnal cycle
   (\np[=.true.]{ln_dm2dc}{ln\_dm2dc}),
-\item
-  the activation of wave effects from an external wave model  (\np[=.true.]{ln_wave}{ln\_wave}),
-\item
-  a neutral drag coefficient is read from an external wave model (\np[=.true.]{ln_cdgw}{ln\_cdgw}),
-\item
-  the Stokes drift from an external wave model is accounted for (\np[=.true.]{ln_sdw}{ln\_sdw}),
-\item
-  the choice of the Stokes drift profile parameterization (\np[=0..2]{nn_sdrift}{nn\_sdrift}),
-\item
-  the surface stress given to the ocean is modified by surface waves (\np[=.true.]{ln_tauwoc}{ln\_tauwoc}),
-\item
-  the surface stress given to the ocean is read from an external wave model (\np[=.true.]{ln_tauw}{ln\_tauw}),
-\item
-  the Stokes-Coriolis term is included (\np[=.true.]{ln_stcor}{ln\_stcor}),
-\item
-  the light penetration in the ocean (\np[=.true.]{ln_traqsr}{ln\_traqsr} with namelist \nam{tra_qsr}{tra\_qsr}),
-\item
-  the atmospheric surface pressure gradient effect on ocean and ice dynamics (\np[=.true.]{ln_apr_dyn}{ln\_apr\_dyn} with namelist \nam{sbc_apr}{sbc\_apr}),
-\item
-  the effect of sea-ice pressure on the ocean (\np[=.true.]{ln_ice_embd}{ln\_ice\_embd}).
+\item the activation of wave effects from an external wave model  (\np[=.true.]{ln_wave}{ln\_wave}),
+\item a neutral drag coefficient is read from an external wave model (\np[=.true.]{ln_cdgw}{ln\_cdgw}),
+\item the Stokes drift from an external wave model is accounted for (\np[=.true.]{ln_sdw}{ln\_sdw}),
+\item the choice of the Stokes drift profile parameterization (\np[=0..2]{nn_sdrift}{nn\_sdrift}),
+\item the surface stress given to the ocean is modified by surface waves (\np[=.true.]{ln_tauwoc}{ln\_tauwoc}),
+\item the surface stress given to the ocean is read from an external wave model (\np[=.true.]{ln_tauw}{ln\_tauw}),
+\item the Stokes-Coriolis term is included (\np[=.true.]{ln_stcor}{ln\_stcor}),
+\item the light penetration in the ocean (\np[=.true.]{ln_traqsr}{ln\_traqsr} with namelist \nam{tra_qsr}{tra\_qsr}),
+\item the atmospheric surface pressure gradient effect on ocean and ice dynamics (\np[=.true.]{ln_apr_dyn}{ln\_apr\_dyn} with namelist \nam{sbc_apr}{sbc\_apr}),
+\item the effect of sea-ice pressure on the ocean (\np[=.true.]{ln_ice_embd}{ln\_ice\_embd}).
 \end{itemize}
 
@@ -119,9 +87,4 @@
 which provides additional sources of fresh water.
 
-
-
-% ================================================================
-% Surface boundary condition for the ocean
-% ================================================================
 \section{Surface boundary condition for the ocean}
 \label{sec:SBC_ocean}
@@ -155,5 +118,4 @@
 $(ii)$ it changes the surface temperature and salinity through the heat and salt contents of
 the mass exchanged with atmosphere, sea-ice and ice shelves.
-
 
 %\colorbox{yellow}{Miss: }
@@ -185,5 +147,4 @@
 these averaged fields are used to compute the surface fluxes at the frequency of \np{nn_fsbc}{nn\_fsbc} time-steps.
 
-
 %-------------------------------------------------TABLE---------------------------------------------------
 \begin{table}[tb]
@@ -210,9 +171,4 @@
 %\colorbox{yellow}{Penser a} mettre dans le restant l'info nn\_fsbc ET nn\_fsbc*rdt de sorte de reinitialiser la moyenne si on change la frequence ou le pdt
 
-
-
-% ================================================================
-%       Input Data
-% ================================================================
 \section{Input data generic interface}
 \label{sec:SBC_input}
@@ -223,13 +179,9 @@
 The module is designed with four main objectives in mind:
 \begin{enumerate}
-\item
-  optionally provide a time interpolation of the input data every specified model time-step, whatever their input frequency is,
+\item optionally provide a time interpolation of the input data every specified model time-step, whatever their input frequency is,
   and according to the different calendars available in the model.
-\item
-  optionally provide an on-the-fly space interpolation from the native input data grid to the model grid.
-\item
-  make the run duration independent from the period cover by the input files.
-\item
-  provide a simple user interface and a rather simple developer interface by
+\item optionally provide an on-the-fly space interpolation from the native input data grid to the model grid.
+\item make the run duration independent from the period cover by the input files.
+\item provide a simple user interface and a rather simple developer interface by
   limiting the number of prerequisite informations.
 \end{enumerate}
@@ -251,8 +203,4 @@
 By default, the data are assumed to be in the same directory as the executable, so that cn\_dir='./'.
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Input Data specification (\mdl{fldread})
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Input data specification (\textit{fldread.F90})]{Input data specification (\protect\mdl{fldread})}
 \label{subsec:SBC_fldread}
@@ -265,5 +213,5 @@
 where
 \begin{description}
-\item[File name]:
+\item [File name]:
   the stem name of the NetCDF file to be opened.
   This stem will be completed automatically by the model, with the addition of a '.nc' at its end and
@@ -301,6 +249,5 @@
 %--------------------------------------------------------------------------------------------------------------
 
-
-\item[Record frequency]:
+\item [Record frequency]:
   the frequency of the records contained in the input file.
   Its unit is in hours if it is positive (for example 24 for daily forcing) or in months if negative
@@ -309,8 +256,8 @@
   On some computers, setting it to '24.' can be interpreted as 240!
 
-\item[Variable name]:
+\item [Variable name]:
   the name of the variable to be read in the input NetCDF file.
 
-\item[Time interpolation]:
+\item [Time interpolation]:
   a logical to activate, or not, the time interpolation.
   If set to 'false', the forcing will have a steplike shape remaining constant during each forcing period.
@@ -322,5 +269,5 @@
   linear interpolation will be performed between mid-day of two consecutive days.
 
-\item[Climatological forcing]:
+\item [Climatological forcing]:
   a logical to specify if a input file contains climatological forcing which can be cycle in time,
   or an interannual forcing which will requires additional files if
@@ -328,5 +275,5 @@
   See the above file naming strategy which impacts the expected name of the file to be opened.
 
-\item[Open/close frequency]:
+\item [Open/close frequency]:
   the frequency at which forcing files must be opened/closed.
   Four cases are coded:
@@ -337,5 +284,5 @@
   the experiment is not starting at the beginning of the year.
 
-\item[Others]:
+\item [Others]:
   'weights filename', 'pairing rotation' and 'land/sea mask' are associated with
   on-the-fly interpolation which is described in \autoref{subsec:SBC_iof}.
@@ -378,8 +325,4 @@
 a useful feature for user considering that it is too heavy to manipulate the complete file for year Y-1.
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Interpolation on the Fly
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Interpolation on-the-fly}
 \label{subsec:SBC_iof}
@@ -404,8 +347,4 @@
 Note that nn\_lsm=0 forces the code to not apply the procedure, even if a land/sea mask file is supplied.
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Bilinear interpolation
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Bilinear interpolation}
 \label{subsec:SBC_iof_bilinear}
@@ -429,8 +368,4 @@
 and wgt(1) corresponds to variable "wgt01" for example.
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Bicubic interpolation
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Bicubic interpolation}
 \label{subsec:SBC_iof_bicubic}
@@ -451,8 +386,4 @@
 the spatial dependency has been included into the weights.
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Implementation
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Implementation}
 \label{subsec:SBC_iof_imp}
@@ -490,24 +421,15 @@
 or is a copy of one from the first few columns on the opposite side of the grid (cyclical case).
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Limitations
-% -------------------------------------------------------------------------------------------------------------
 \subsubsection{Limitations}
 \label{subsec:SBC_iof_lim}
 
 \begin{enumerate}
-\item
-  The case where input data grids are not logically rectangular (irregular grid case) has not been tested.
-\item
-  This code is not guaranteed to produce positive definite answers from positive definite inputs when
+\item The case where input data grids are not logically rectangular (irregular grid case) has not been tested.
+\item This code is not guaranteed to produce positive definite answers from positive definite inputs when
   a bicubic interpolation method is used.
-\item
-  The cyclic condition is only applied on left and right columns, and not to top and bottom rows.
-\item
-  The gradients across the ends of a cyclical grid assume that the grid spacing between
+\item The cyclic condition is only applied on left and right columns, and not to top and bottom rows.
+\item The gradients across the ends of a cyclical grid assume that the grid spacing between
   the two columns involved are consistent with the weights used.
-\item
-  Neither interpolation scheme is conservative. (There is a conservative scheme available in SCRIP,
+\item Neither interpolation scheme is conservative. (There is a conservative scheme available in SCRIP,
   but this has not been implemented.)
 \end{enumerate}
@@ -520,8 +442,4 @@
 (see the directory NEMOGCM/TOOLS/WEIGHTS).
 
-
-% -------------------------------------------------------------------------------------------------------------
-% Standalone Surface Boundary Condition Scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Standalone surface boundary condition scheme (SAS)}
 \label{subsec:SBC_SAS}
@@ -541,15 +459,10 @@
 
 \begin{itemize}
-\item
-  Multiple runs of the model are required in code development to
+\item Multiple runs of the model are required in code development to
   see the effect of different algorithms in the bulk formulae.
-\item
-  The effect of different parameter sets in the ice model is to be examined.
-\item
-  Development of sea-ice algorithms or parameterizations.
-\item
-  Spinup of the iceberg floats
-\item
-  Ocean/sea-ice simulation with both models running in parallel (\np[=.true.]{ln_mixcpl}{ln\_mixcpl})
+\item The effect of different parameter sets in the ice model is to be examined.
+\item Development of sea-ice algorithms or parameterizations.
+\item Spinup of the iceberg floats
+\item Ocean/sea-ice simulation with both models running in parallel (\np[=.true.]{ln_mixcpl}{ln\_mixcpl})
 \end{itemize}
 
@@ -562,27 +475,21 @@
 
 \begin{itemize}
-\item
-  \mdl{nemogcm}:
+\item \mdl{nemogcm}:
   This routine initialises the rest of the model and repeatedly calls the stp time stepping routine (\mdl{step}).
   Since the ocean state is not calculated all associated initialisations have been removed.
-\item
-  \mdl{step}:
+\item \mdl{step}:
   The main time stepping routine now only needs to call the sbc routine (and a few utility functions).
-\item
-  \mdl{sbcmod}:
+\item \mdl{sbcmod}:
   This has been cut down and now only calculates surface forcing and the ice model required.
   New surface modules that can function when only the surface level of the ocean state is defined can also be added
   (\eg\ icebergs).
-\item
-  \mdl{daymod}:
+\item \mdl{daymod}:
   No ocean restarts are read or written (though the ice model restarts are retained),
   so calls to restart functions have been removed.
   This also means that the calendar cannot be controlled by time in a restart file,
   so the user must check that nn\_date0 in the model namelist is correct for his or her purposes.
-\item
-  \mdl{stpctl}:
+\item \mdl{stpctl}:
   Since there is no free surface solver, references to it have been removed from \rou{stp\_ctl} module.
-\item
-  \mdl{diawri}:
+\item \mdl{diawri}:
   All 3D data have been removed from the output.
   The surface temperature, salinity and velocity components (which have been read in) are written along with
@@ -593,6 +500,5 @@
 
 \begin{itemize}
-\item
-  \mdl{sbcsas}:
+\item \mdl{sbcsas}:
   This module initialises the input files needed for reading temperature, salinity and
   velocity arrays at the surface.
@@ -604,13 +510,7 @@
 \end{itemize}
 
-
 The user can also choose in the \nam{sbc_sas}{sbc\_sas} namelist to read the mean (nn\_fsbc time-step) fraction of solar net radiation absorbed in the 1st T level using
  (\np[=.true.]{ln_flx}{ln\_flx}) and to provide 3D oceanic velocities instead of 2D ones (\np{ln_flx}{ln\_flx}\forcode{=.true.}). In that last case, only the 1st level will be read in.
 
-
-
-% ================================================================
-% Flux formulation
-% ================================================================
 \section[Flux formulation (\textit{sbcflx.F90})]{Flux formulation (\protect\mdl{sbcflx})}
 \label{sec:SBC_flx}
@@ -634,9 +534,4 @@
 See \autoref{subsec:SBC_ssr} for its specification.
 
-
-
-% ================================================================
-% Bulk formulation
-% ================================================================
 \section[Bulk formulation (\textit{sbcblk.F90})]{Bulk formulation (\protect\mdl{sbcblk})}
 \label{sec:SBC_blk}
@@ -720,8 +615,4 @@
 the namsbc\_blk namelist (see \autoref{subsec:SBC_fldread}).
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Ocean-Atmosphere Bulk formulae
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Ocean-Atmosphere Bulk formulae (\textit{sbcblk\_algo\_coare.F90, sbcblk\_algo\_coare3p5.F90, sbcblk\_algo\_ecmwf.F90, sbcblk\_algo\_ncar.F90})]{Ocean-Atmosphere Bulk formulae (\mdl{sbcblk\_algo\_coare}, \mdl{sbcblk\_algo\_coare3p5}, \mdl{sbcblk\_algo\_ecmwf}, \mdl{sbcblk\_algo\_ncar})}
 \label{subsec:SBC_blk_ocean}
@@ -731,6 +622,5 @@
 their neutral transfer coefficients relationships with neutral wind.
 \begin{itemize}
-\item
-  NCAR (\np[=.true.]{ln_NCAR}{ln\_NCAR}):
+\item NCAR (\np[=.true.]{ln_NCAR}{ln\_NCAR}):
   The NCAR bulk formulae have been developed by \citet{large.yeager_rpt04}.
   They have been designed to handle the NCAR forcing, a mixture of NCEP reanalysis and satellite data.
@@ -741,19 +631,13 @@
   Note that substituting ERA40 to NCEP reanalysis fields does not require changes in the bulk formulea themself.
   This is the so-called DRAKKAR Forcing Set (DFS) \citep{brodeau.barnier.ea_OM10}.
-\item
-  COARE 3.0 (\np[=.true.]{ln_COARE_3p0}{ln\_COARE\_3p0}):
+\item COARE 3.0 (\np[=.true.]{ln_COARE_3p0}{ln\_COARE\_3p0}):
   See \citet{fairall.bradley.ea_JC03} for more details
-\item
-  COARE 3.5 (\np[=.true.]{ln_COARE_3p5}{ln\_COARE\_3p5}):
+\item COARE 3.5 (\np[=.true.]{ln_COARE_3p5}{ln\_COARE\_3p5}):
   See \citet{edson.jampana.ea_JPO13} for more details
-\item
-  ECMWF (\np[=.true.]{ln_ECMWF}{ln\_ECMWF}):
+\item ECMWF (\np[=.true.]{ln_ECMWF}{ln\_ECMWF}):
   Based on \href{https://www.ecmwf.int/node/9221}{IFS (Cy31)} implementation and documentation.
   Surface roughness lengths needed for the Obukhov length are computed following \citet{beljaars_QJRMS95}.
 \end{itemize}
 
-% -------------------------------------------------------------------------------------------------------------
-%        Ice-Atmosphere Bulk formulae
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Ice-Atmosphere Bulk formulae}
 \label{subsec:SBC_blk_ice}
@@ -762,9 +646,7 @@
 
 \begin{itemize}
-\item
-  Constant value (\np[ Cd_ice=1.4e-3 ]{constant value}{constant\ value}):
+\item Constant value (\np[ Cd_ice=1.4e-3 ]{constant value}{constant\ value}):
   default constant value used for momentum and heat neutral transfer coefficients
-\item
-  \citet{lupkes.gryanik.ea_JGR12} (\np[=.true.]{ln_Cd_L12}{ln\_Cd\_L12}):
+\item \citet{lupkes.gryanik.ea_JGR12} (\np[=.true.]{ln_Cd_L12}{ln\_Cd\_L12}):
   This scheme adds a dependency on edges at leads, melt ponds and flows
   of the constant neutral air-ice drag. After some approximations,
@@ -773,6 +655,5 @@
   starting at 1.5e-3 for A=0, reaching 1.97e-3 for A=0.5 and going down 1.4e-3 for A=1.
   It is theoretically applicable to all ice conditions (not only MIZ).
-\item
-  \citet{lupkes.gryanik_JGR15} (\np[=.true.]{ln_Cd_L15}{ln\_Cd\_L15}):
+\item \citet{lupkes.gryanik_JGR15} (\np[=.true.]{ln_Cd_L15}{ln\_Cd\_L15}):
   Alternative turbulent transfer coefficients formulation between sea-ice
   and atmosphere with distinct momentum and heat coefficients depending
@@ -784,9 +665,4 @@
 \end{itemize}
 
-
-
-% ================================================================
-% Coupled formulation
-% ================================================================
 \section[Coupled formulation (\textit{sbccpl.F90})]{Coupled formulation (\protect\mdl{sbccpl})}
 \label{sec:SBC_cpl}
@@ -826,9 +702,4 @@
 In cases where this is definitely not possible, the model should abort with an error message.
 
-
-
-% ================================================================
-%        Atmospheric pressure
-% ================================================================
 \section[Atmospheric pressure (\textit{sbcapr.F90})]{Atmospheric pressure (\protect\mdl{sbcapr})}
 \label{sec:SBC_apr}
@@ -868,9 +739,4 @@
 \np{ln_apr_obc}{ln\_apr\_obc}  might be set to true.
 
-
-
-% ================================================================
-%        Surface Tides Forcing
-% ================================================================
 \section[Surface tides (\textit{sbctide.F90})]{Surface tides (\protect\mdl{sbctide})}
 \label{sec:SBC_tide}
@@ -925,9 +791,4 @@
 \forcode{.false.} removes the SAL contribution.
 
-
-
-% ================================================================
-%        River runoffs
-% ================================================================
 \section[River runoffs (\textit{sbcrnf.F90})]{River runoffs (\protect\mdl{sbcrnf})}
 \label{sec:SBC_rnf}
@@ -950,5 +811,4 @@
 %required to properly represent the diurnal cycle \citep{bernie.woolnough.ea_JC05}. see also \autoref{fig:SBC_dcy}.}.
 
-
 %To do this we need to treat evaporation/precipitation fluxes and river runoff differently in the
 %\mdl{tra\_sbc} module.  We decided to separate them throughout the code, so that the variable
@@ -957,5 +817,4 @@
 %emp).  This meant many uses of emp and emps needed to be changed, a list of all modules which use
 %emp or emps and the changes made are below:
-
 
 %Rachel:
@@ -1035,5 +894,4 @@
 as the ocean water leaving the domain removes heat and salt (at the same concentration) with it.
 
-
 %\colorbox{yellow}{Nevertheless, Pb of vertical resolution and 3D input : increase vertical mixing near river mouths to mimic a 3D river
 
@@ -1055,9 +913,4 @@
 %To do this we need to treat evaporation/precipitation fluxes and river runoff differently in the tra_sbc module.  We decided to separate them throughout the code, so that the variable emp represented solely evaporation minus precipitation fluxes, and a new 2d variable rnf was added which represents the volume flux of river runoff (in kg/m2s to remain consistent with emp).  This meant many uses of emp and emps needed to be changed, a list of all modules which use emp or emps and the changes made are below:
 
-
-
-% ================================================================
-%        Ice shelf melting
-% ================================================================
 \section[Ice shelf melting (\textit{sbcisf.F90})]{Ice shelf melting (\protect\mdl{sbcisf})}
 \label{sec:SBC_isf}
@@ -1077,5 +930,5 @@
 \begin{description}
 
-  \item[{\np[=1]{nn_isf}{nn\_isf}}]:
+  \item [{\np[=1]{nn_isf}{nn\_isf}}]:
   The ice shelf cavity is represented (\np[=.true.]{ln_isfcav}{ln\_isfcav} needed).
   The fwf and heat flux are depending of the local water properties.
@@ -1084,8 +937,8 @@
 
 	\begin{description}
-	\item[{\np[=1]{nn_isfblk}{nn\_isfblk}}]:
+	\item [{\np[=1]{nn_isfblk}{nn\_isfblk}}]:
      The melt rate is based on a balance between the upward ocean heat flux and
      the latent heat flux at the ice shelf base. A complete description is available in \citet{hunter_rpt06}.
-	\item[{\np[=2]{nn_isfblk}{nn\_isfblk}}]:
+	\item [{\np[=2]{nn_isfblk}{nn\_isfblk}}]:
      The melt rate and the heat flux are based on a 3 equations formulation
      (a heat flux budget at the ice base, a salt flux budget at the ice base and a linearised freezing point temperature equation).
@@ -1104,5 +957,5 @@
      There are 3 different ways to compute the exchange coeficient:
 	\begin{description}
-        \item[{\np[=0]{nn_gammablk}{nn\_gammablk}}]:
+        \item [{\np[=0]{nn_gammablk}{nn\_gammablk}}]:
      The salt and heat exchange coefficients are constant and defined by \np{rn_gammas0}{rn\_gammas0} and \np{rn_gammat0}{rn\_gammat0}.
      \begin{gather*}
@@ -1112,5 +965,5 @@
      \end{gather*}
      This is the recommended formulation for ISOMIP.
-	\item[{\np[=1]{nn_gammablk}{nn\_gammablk}}]:
+	\item [{\np[=1]{nn_gammablk}{nn\_gammablk}}]:
      The salt and heat exchange coefficients are velocity dependent and defined as
      \begin{gather*}
@@ -1120,5 +973,5 @@
      where $u_{*}$ is the friction velocity in the top boundary layer (ie first \np{rn_hisf_tbl}{rn\_hisf\_tbl} meters).
      See \citet{jenkins.nicholls.ea_JPO10} for all the details on this formulation. It is the recommended formulation for realistic application.
-	\item[{\np[=2]{nn_gammablk}{nn\_gammablk}}]:
+	\item [{\np[=2]{nn_gammablk}{nn\_gammablk}}]:
      The salt and heat exchange coefficients are velocity and stability dependent and defined as:
 \[
@@ -1131,5 +984,5 @@
      This formulation has not been extensively tested in \NEMO\ (not recommended).
    \end{description}
-  \item[{\np[=2]{nn_isf}{nn\_isf}}]:
+  \item [{\np[=2]{nn_isf}{nn\_isf}}]:
    The ice shelf cavity is not represented.
    The fwf and heat flux are computed using the \citet{beckmann.goosse_OM03} parameterisation of isf melting.
@@ -1138,5 +991,5 @@
    (\np{sn_depmin_isf}{sn\_depmin\_isf}) as in (\np[=3]{nn_isf}{nn\_isf}).
    The effective melting length (\np{sn_Leff_isf}{sn\_Leff\_isf}) is read from a file.
-  \item[{\np[=3]{nn_isf}{nn\_isf}}]:
+  \item [{\np[=3]{nn_isf}{nn\_isf}}]:
    The ice shelf cavity is not represented.
    The fwf (\np{sn_rnfisf}{sn\_rnfisf}) is prescribed and distributed along the ice shelf edge between
@@ -1144,5 +997,5 @@
    the base of the ice shelf along the calving front (\np{sn_depmin_isf}{sn\_depmin\_isf}).
    The heat flux ($Q_h$) is computed as $Q_h = fwf \times L_f$.
-  \item[{\np[=4]{nn_isf}{nn\_isf}}]:
+  \item [{\np[=4]{nn_isf}{nn\_isf}}]:
    The ice shelf cavity is opened (\np[=.true.]{ln_isfcav}{ln\_isfcav} needed).
    However, the fwf is not computed but specified from file \np{sn_fwfisf}{sn\_fwfisf}).
@@ -1167,5 +1020,4 @@
 See the runoff section \autoref{sec:SBC_rnf} for all the details about the divergence correction.\\
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -1176,11 +1028,5 @@
   \label{fig:SBC_isf}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-
-
-% ================================================================
-%        Ice sheet coupling
-% ================================================================
+
 \section{Ice sheet coupling}
 \label{sec:SBC_iscpl}
@@ -1198,9 +1044,9 @@
 
 \begin{description}
-\item[Step 1]: the ice sheet model send a new bathymetry and ice shelf draft netcdf file.
-\item[Step 2]: a new domcfg.nc file is built using the DOMAINcfg tools.
-\item[Step 3]: \NEMO\ run for a specific period and output the average melt rate over the period.
-\item[Step 4]: the ice sheet model run using the melt rate outputed in step 4.
-\item[Step 5]: go back to 1.
+\item [Step 1]: the ice sheet model send a new bathymetry and ice shelf draft netcdf file.
+\item [Step 2]: a new domcfg.nc file is built using the DOMAINcfg tools.
+\item [Step 3]: \NEMO\ run for a specific period and output the average melt rate over the period.
+\item [Step 4]: the ice sheet model run using the melt rate outputed in step 4.
+\item [Step 5]: go back to 1.
 \end{description}
 
@@ -1210,19 +1056,19 @@
 
 \begin{description}
-\item[Thin a cell down]:
+\item [Thin a cell down]:
   T/S/ssh are unchanged and U/V in the top cell are corrected to keep the barotropic transport (bt) constant
   ($bt_b=bt_n$).
-\item[Enlarge  a cell]:
+\item [Enlarge  a cell]:
   See case "Thin a cell down"
-\item[Dry a cell]:
+\item [Dry a cell]:
   mask, T/S, U/V and ssh are set to 0.
   Furthermore, U/V into the water column are modified to satisfy ($bt_b=bt_n$).
-\item[Wet a cell]:
+\item [Wet a cell]:
   mask is set to 1, T/S is extrapolated from neighbours, $ssh_n = ssh_b$ and U/V set to 0.
   If no neighbours, T/S is extrapolated from old top cell value.
   If no neighbours along i,j and k (both previous test failed), T/S/U/V/ssh and mask are set to 0.
-\item[Dry a column]:
+\item [Dry a column]:
 	mask, T/S, U/V are set to 0 everywhere in the column and ssh set to 0.
-\item[Wet a column]:
+\item [Wet a column]:
   set mask to 1, T/S is extrapolated from neighbours, ssh is extrapolated from neighbours and U/V set to 0.
   If no neighbour, T/S/U/V and mask set to 0.
@@ -1247,9 +1093,4 @@
 The corrective increment is apply into the cell itself (if it is a wet cell), the neigbouring cells or the closest wet cell (if the cell is now dry).
 
-
-
-% ================================================================
-%        Handling of icebergs
-% ================================================================
 \section{Handling of icebergs (ICB)}
 \label{sec:SBC_ICB_icebergs}
@@ -1275,5 +1116,5 @@
 Two initialisation schemes are possible.
 \begin{description}
-\item[{\np{nn_test_icebergs}{nn\_test\_icebergs}~$>$~0}]
+\item [{\np{nn_test_icebergs}{nn\_test\_icebergs}~$>$~0}]
   In this scheme, the value of \np{nn_test_icebergs}{nn\_test\_icebergs} represents the class of iceberg to generate
   (so between 1 and 10), and \np{nn_test_icebergs}{nn\_test\_icebergs} provides a lon/lat box in the domain at each grid point of
@@ -1282,5 +1123,5 @@
   \np{nn_test_icebergs}{nn\_test\_icebergs} is defined by four numbers in \np{nn_test_box}{nn\_test\_box} representing the corners of
   the geographical box: lonmin,lonmax,latmin,latmax
-\item[{\np[=-1]{nn_test_icebergs}{nn\_test\_icebergs}}]
+\item [{\np[=-1]{nn_test_icebergs}{nn\_test\_icebergs}}]
   In this scheme, the model reads a calving file supplied in the \np{sn_icb}{sn\_icb} parameter.
   This should be a file with a field on the configuration grid (typically ORCA)
@@ -1307,8 +1148,8 @@
 The amount of information is controlled by two integer parameters:
 \begin{description}
-\item[{\np{nn_verbose_level}{nn\_verbose\_level}}] takes a value between one and four and
+\item [{\np{nn_verbose_level}{nn\_verbose\_level}}] takes a value between one and four and
   represents an increasing number of points in the code at which variables are written,
   and an increasing level of obscurity.
-\item[{\np{nn_verbose_write}{nn\_verbose\_write}}] is the number of timesteps between writes
+\item [{\np{nn_verbose_write}{nn\_verbose\_write}}] is the number of timesteps between writes
 \end{description}
 
@@ -1321,9 +1162,4 @@
 since its trajectory data may be spread across multiple files.
 
-
-
-% =============================================================================================================
-%        Interactions with waves (sbcwave.F90, ln_wave)
-% =============================================================================================================
 \section[Interactions with waves (\textit{sbcwave.F90}, \forcode{ln_wave})]{Interactions with waves (\protect\mdl{sbcwave}, \protect\np{ln_wave}{ln\_wave})}
 \label{sec:SBC_wave}
@@ -1349,11 +1185,10 @@
 Wave fields can be provided either in forced or coupled mode:
 \begin{description}
-\item[forced mode]: wave fields should be defined through the \nam{sbc_wave}{sbc\_wave} namelist
+\item [forced mode]: wave fields should be defined through the \nam{sbc_wave}{sbc\_wave} namelist
 for external data names, locations, frequency, interpolation and all the miscellanous options allowed by
 Input Data generic Interface (see \autoref{sec:SBC_input}).
-\item[coupled mode]: \NEMO\ and an external wave model can be coupled by setting \np[=.true.]{ln_cpl}{ln\_cpl}
+\item [coupled mode]: \NEMO\ and an external wave model can be coupled by setting \np[=.true.]{ln_cpl}{ln\_cpl}
 in \nam{sbc}{sbc} namelist and filling the \nam{sbc_cpl}{sbc\_cpl} namelist.
 \end{description}
-
 
 % ----------------------------------------------------------------
@@ -1369,5 +1204,4 @@
 the drag coefficient is computed according to the stable/unstable conditions of the
 air-sea interface following \citet{large.yeager_rpt04}.
-
 
 % ----------------------------------------------------------------
@@ -1409,5 +1243,5 @@
 
 \begin{description}
-\item[{\np{nn_sdrift}{nn\_sdrift} = 0}]: exponential integral profile parameterization proposed by
+\item [{\np{nn_sdrift}{nn\_sdrift} = 0}]: exponential integral profile parameterization proposed by
 \citet{breivik.janssen.ea_JPO14}:
 
@@ -1428,5 +1262,5 @@
 where $H_s$ is the significant wave height and $\omega$ is the wave frequency.
 
-\item[{\np{nn_sdrift}{nn\_sdrift} = 1}]: velocity profile based on the Phillips spectrum which is considered to be a
+\item [{\np{nn_sdrift}{nn\_sdrift} = 1}]: velocity profile based on the Phillips spectrum which is considered to be a
 reasonable estimate of the part of the spectrum mostly contributing to the Stokes drift velocity near the surface
 \citep{breivik.bidlot.ea_OM16}:
@@ -1440,5 +1274,5 @@
 where $erf$ is the complementary error function and $k_p$ is the peak wavenumber.
 
-\item[{\np{nn_sdrift}{nn\_sdrift} = 2}]: velocity profile based on the Phillips spectrum as for \np{nn_sdrift}{nn\_sdrift} = 1
+\item [{\np{nn_sdrift}{nn\_sdrift} = 2}]: velocity profile based on the Phillips spectrum as for \np{nn_sdrift}{nn\_sdrift} = 1
 but using the wave frequency from a wave model.
 
@@ -1465,5 +1299,4 @@
   - (\mathbf{U} + \mathbf{U}_{st}) \cdot \nabla{c}
 \]
-
 
 % ----------------------------------------------------------------
@@ -1479,5 +1312,4 @@
 approximations described in \autoref{subsec:SBC_wave_sdw}),
 \np[=.true.]{ln_stcor}{ln\_stcor} has to be set.
-
 
 % ----------------------------------------------------------------
@@ -1521,16 +1353,7 @@
 meridional stress components by setting \np[=.true.]{ln_tauw}{ln\_tauw}.
 
-
-
-% ================================================================
-% Miscellanea options
-% ================================================================
 \section{Miscellaneous options}
 \label{sec:SBC_misc}
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Diurnal cycle
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Diurnal cycle (\textit{sbcdcy.F90})]{Diurnal cycle (\protect\mdl{sbcdcy})}
 \label{subsec:SBC_dcy}
@@ -1540,5 +1363,4 @@
 %-------------------------------------------------------------------------------------------------------------
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -1553,5 +1375,4 @@
   \label{fig:SBC_diurnal}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 \cite{bernie.woolnough.ea_JC05} have shown that to capture 90$\%$ of the diurnal variability of SST requires a vertical resolution in upper ocean of 1~m or better and a temporal resolution of the surface fluxes of 3~h or less.
@@ -1576,5 +1397,4 @@
 one every 2~hours (from 1am to 11pm).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -1586,5 +1406,4 @@
   \label{fig:SBC_dcy}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 Note also that the setting a diurnal cycle in SWF is highly recommended when
@@ -1592,8 +1411,4 @@
 an inconsistency between the scale of the vertical resolution and the forcing acting on that scale.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Rotation of vector pairs onto the model grid directions
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Rotation of vector pairs onto the model grid directions}
 \label{subsec:SBC_rotation}
@@ -1612,8 +1427,4 @@
 The rot\_rep routine from the \mdl{geo2ocean} module is used to perform the rotation.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Surface restoring to observed SST and/or SSS
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Surface restoring to observed SST and/or SSS (\textit{sbcssr.F90})]{Surface restoring to observed SST and/or SSS (\protect\mdl{sbcssr})}
 \label{subsec:SBC_ssr}
@@ -1660,8 +1471,4 @@
 reduce the uncertainties we have on the observed freshwater budget.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Handling of ice-covered area
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Handling of ice-covered area  (\textit{sbcice\_...})}
 \label{subsec:SBC_ice-cover}
@@ -1671,10 +1478,10 @@
 the value of the \np{nn_ice}{nn\_ice} namelist parameter found in \nam{sbc}{sbc} namelist.
 \begin{description}
-\item[nn\_ice = 0]
+\item [nn\_ice = 0]
   there will never be sea-ice in the computational domain.
   This is a typical namelist value used for tropical ocean domain.
   The surface fluxes are simply specified for an ice-free ocean.
   No specific things is done for sea-ice.
-\item[nn\_ice = 1]
+\item [nn\_ice = 1]
   sea-ice can exist in the computational domain, but no sea-ice model is used.
   An observed ice covered area is read in a file.
@@ -1688,5 +1495,5 @@
   is usually referred as the \textit{ice-if} model.
   It can be found in the \mdl{sbcice\_if} module.
-\item[nn\_ice = 2 or more]
+\item [nn\_ice = 2 or more]
   A full sea ice model is used.
   This model computes the ice-ocean fluxes,
@@ -1701,8 +1508,4 @@
 %GS: ocean-ice (SI3) interface is not located in SBC directory anymore, so it should be included in SI3 doc
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        CICE-ocean Interface
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Interface to CICE (\textit{sbcice\_cice.F90})]{Interface to CICE (\protect\mdl{sbcice\_cice})}
 \label{subsec:SBC_cice}
@@ -1735,8 +1538,4 @@
 there is no sea ice.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        Freshwater budget control
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Freshwater budget control (\textit{sbcfwb.F90})]{Freshwater budget control (\protect\mdl{sbcfwb})}
 \label{subsec:SBC_fwb}
@@ -1747,12 +1546,12 @@
 
 \begin{description}
-\item[{\np[=0]{nn_fwb}{nn\_fwb}}]
+\item [{\np[=0]{nn_fwb}{nn\_fwb}}]
   no control at all.
   The mean sea level is free to drift, and will certainly do so.
-\item[{\np[=1]{nn_fwb}{nn\_fwb}}]
+\item [{\np[=1]{nn_fwb}{nn\_fwb}}]
   global mean \textit{emp} set to zero at each model time step.
   %GS: comment below still relevant ?
   %Note that with a sea-ice model, this technique only controls the mean sea level with linear free surface and no mass flux between ocean and ice (as it is implemented in the current ice-ocean coupling).
-\item[{\np[=2]{nn_fwb}{nn\_fwb}}]
+\item [{\np[=2]{nn_fwb}{nn\_fwb}}]
   freshwater budget is adjusted from the previous year annual mean budget which
   is read in the \textit{EMPave\_old.dat} file.
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_STO.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_STO.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_STO.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter stochastic parametrization of EOS (STO)
-% ================================================================
 \chapter{Stochastic Parametrization of EOS (STO)}
 \label{chap:STO}
@@ -24,6 +21,4 @@
 P.-A. Bouttier release 3.6 inital version
 
-\newpage
-
 As a result of the nonlinearity of the seawater equation of state, unresolved scales represent a major source of uncertainties in the computation of the large-scale horizontal density gradient from the large-scale temperature and salinity fields. Following  \cite{brankart_OM13}, the impact of these uncertainties can be simulated by random processes representing unresolved T/S fluctuations. The Stochastic Parametrization of EOS (STO) module implements this parametrization.
 
@@ -43,5 +38,4 @@
 a parametrized decorrelation time scale, and horizontal and vertical standard deviations $\sigma_s$.
 $\mathbf{\xi}$ are uncorrelated over the horizontal and fully correlated along the vertical.
-
 
 \section{Stochastic processes}
@@ -69,6 +63,5 @@
 
 \begin{itemize}
-\item
-  for order~1 processes, $w^{(i)}$ is a Gaussian white noise, with zero mean and standard deviation equal to~1,
+\item for order~1 processes, $w^{(i)}$ is a Gaussian white noise, with zero mean and standard deviation equal to~1,
   and the parameters $a^{(i)}$, $b^{(i)}$, $c^{(i)}$ are given by:
 
@@ -84,6 +77,5 @@
   \]
 
-\item
-  for order~$n>1$ processes, $w^{(i)}$ is an order~$n-1$ autoregressive process, with zero mean,
+\item for order~$n>1$ processes, $w^{(i)}$ is an order~$n-1$ autoregressive process, with zero mean,
   standard deviation equal to~$\sigma^{(i)}$;
   correlation timescale equal to~$\tau^{(i)}$;
@@ -124,5 +116,4 @@
 \section{Implementation details}
 \label{sec:STO_thech_details}
-
 
 The code implementing stochastic parametrisation is located in the src/OCE/STO directory.
@@ -186,11 +177,11 @@
 
 \begin{description}
-\item[{\np{nn_sto_eos}{nn\_sto\_eos}:}]   number of independent random walks
-\item[{\np{rn_eos_stdxy}{rn\_eos\_stdxy}:}] random walk horizontal standard deviation (in grid points)
-\item[{\np{rn_eos_stdz}{rn\_eos\_stdz}:}]  random walk vertical standard deviation (in grid points)
-\item[{\np{rn_eos_tcor}{rn\_eos\_tcor}:}]  random walk time correlation (in timesteps)
-\item[{\np{nn_eos_ord}{nn\_eos\_ord}:}]   order of autoregressive processes
-\item[{\np{nn_eos_flt}{nn\_eos\_flt}:}]   passes of Laplacian filter
-\item[{\np{rn_eos_lim}{rn\_eos\_lim}:}]   limitation factor (default = 3.0)
+\item [{\np{nn_sto_eos}{nn\_sto\_eos}:}]   number of independent random walks
+\item [{\np{rn_eos_stdxy}{rn\_eos\_stdxy}:}] random walk horizontal standard deviation (in grid points)
+\item [{\np{rn_eos_stdz}{rn\_eos\_stdz}:}]  random walk vertical standard deviation (in grid points)
+\item [{\np{rn_eos_tcor}{rn\_eos\_tcor}:}]  random walk time correlation (in timesteps)
+\item [{\np{nn_eos_ord}{nn\_eos\_ord}:}]   order of autoregressive processes
+\item [{\np{nn_eos_flt}{nn\_eos\_flt}:}]   passes of Laplacian filter
+\item [{\np{rn_eos_lim}{rn\_eos\_lim}:}]   limitation factor (default = 3.0)
 \end{description}
 
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_TRA.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_TRA.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_TRA.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter 1 ——— Ocean Tracers (TRA)
-% ================================================================
 \chapter{Ocean Tracers (TRA)}
 \label{chap:TRA}
@@ -57,7 +54,4 @@
 (\np{ln_tra_trd}{ln\_tra\_trd} or \np[=.true.]{ln_tra_mxl}{ln\_tra\_mxl}), as described in \autoref{chap:DIA}.
 
-% ================================================================
-% Tracer Advection
-% ================================================================
 \section[Tracer advection (\textit{traadv.F90})]{Tracer advection (\protect\mdl{traadv})}
 \label{sec:TRA_adv}
@@ -91,5 +85,4 @@
 In other words, by setting $\tau = 1$ in (\autoref{eq:TRA_adv}) we recover the discrete form of
 the continuity equation which is used to calculate the vertical velocity.
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -109,5 +102,4 @@
   \label{fig:TRA_adv_scheme}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The key difference between the advection schemes available in \NEMO\ is the choice made in space and
@@ -120,5 +112,5 @@
 
 \begin{description}
-\item[linear free surface:]
+\item [linear free surface:]
   (\np[=.true.]{ln_linssh}{ln\_linssh})
   the first level thickness is constant in time:
@@ -128,5 +120,5 @@
   $\tau_w|_{k = 1/2} = T_{k = 1}$, \ie\ the product of surface velocity (at $z = 0$) by
   the first level tracer value.
-\item[non-linear free surface:]
+\item [non-linear free surface:]
   (\np[=.false.]{ln_linssh}{ln\_linssh})
   convergence/divergence in the first ocean level moves the free surface up/down.
@@ -163,14 +155,11 @@
 
 \begin{enumerate}
-\item
-  CEN and FCT schemes require an explicit diffusion operator while the other schemes are diffusive enough so that
+\item CEN and FCT schemes require an explicit diffusion operator while the other schemes are diffusive enough so that
   they do not necessarily need additional diffusion;
-\item
-  CEN and UBS are not \textit{positive} schemes
+\item CEN and UBS are not \textit{positive} schemes
   \footnote{negative values can appear in an initially strictly positive tracer field which is advected},
   implying that false extrema are permitted.
   Their use is not recommended on passive tracers;
-\item
-  It is recommended that the same advection-diffusion scheme is used on both active and passive tracers.
+\item It is recommended that the same advection-diffusion scheme is used on both active and passive tracers.
 \end{enumerate}
 
@@ -182,7 +171,4 @@
 their results.
 
-% -------------------------------------------------------------------------------------------------------------
-%        2nd and 4th order centred schemes
-% -------------------------------------------------------------------------------------------------------------
 \subsection[CEN: Centred scheme (\forcode{ln_traadv_cen})]{CEN: Centred scheme (\protect\np{ln_traadv_cen}{ln\_traadv\_cen})}
 \label{subsec:TRA_adv_cen}
@@ -249,7 +235,4 @@
 these near boundary grid points.
 
-% -------------------------------------------------------------------------------------------------------------
-%        FCT scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[FCT: Flux Corrected Transport scheme (\forcode{ln_traadv_fct})]{FCT: Flux Corrected Transport scheme (\protect\np{ln_traadv_fct}{ln\_traadv\_fct})}
 \label{subsec:TRA_adv_tvd}
@@ -285,5 +268,4 @@
 A comparison of FCT-2 with MUSCL and a MPDATA scheme can be found in \citet{levy.estublier.ea_GRL01}.
 
-
 For stability reasons (see \autoref{chap:TD}),
 $\tau_u^{cen}$ is evaluated in (\autoref{eq:TRA_adv_fct}) using the \textit{now} tracer while
@@ -292,7 +274,4 @@
 while a forward scheme is used for the diffusive part.
 
-% -------------------------------------------------------------------------------------------------------------
-%        MUSCL scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[MUSCL: Monotone Upstream Scheme for Conservative Laws (\forcode{ln_traadv_mus})]{MUSCL: Monotone Upstream Scheme for Conservative Laws (\protect\np{ln_traadv_mus}{ln\_traadv\_mus})}
 \label{subsec:TRA_adv_mus}
@@ -328,7 +307,4 @@
 (\np[=.true.]{ln_mus_ups}{ln\_mus\_ups}).
 
-% -------------------------------------------------------------------------------------------------------------
-%        UBS scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[UBS a.k.a. UP3: Upstream-Biased Scheme (\forcode{ln_traadv_ubs})]{UBS a.k.a. UP3: Upstream-Biased Scheme (\protect\np{ln_traadv_ubs}{ln\_traadv\_ubs})}
 \label{subsec:TRA_adv_ubs}
@@ -400,7 +376,4 @@
 Note the current version of \NEMO\ uses the computationally more efficient formulation \autoref{eq:TRA_adv_ubs}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        QCK scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[QCK: QuiCKest scheme (\forcode{ln_traadv_qck})]{QCK: QuiCKest scheme (\protect\np{ln_traadv_qck}{ln\_traadv\_qck})}
 \label{subsec:TRA_adv_qck}
@@ -423,7 +396,4 @@
 %%%gmcomment   :  Cross term are missing in the current implementation....
 
-% ================================================================
-% Tracer Lateral Diffusion
-% ================================================================
 \section[Tracer lateral diffusion (\textit{traldf.F90})]{Tracer lateral diffusion (\protect\mdl{traldf})}
 \label{sec:TRA_ldf}
@@ -455,7 +425,4 @@
 the pure vertical component is split into an explicit and an implicit part \citep{lemarie.debreu.ea_OM12}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Type of operator
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Type of operator (\forcode{ln_traldf_}\{\forcode{OFF,lap,blp}\})]{Type of operator (\protect\np{ln_traldf_OFF}{ln\_traldf\_OFF}, \protect\np{ln_traldf_lap}{ln\_traldf\_lap}, or \protect\np{ln_traldf_blp}{ln\_traldf\_blp})}
 \label{subsec:TRA_ldf_op}
@@ -464,13 +431,13 @@
 
 \begin{description}
-\item[{\np[=.true.]{ln_traldf_OFF}{ln\_traldf\_OFF}}]
+\item [{\np[=.true.]{ln_traldf_OFF}{ln\_traldf\_OFF}}]
   no operator selected, the lateral diffusive tendency will not be applied to the tracer equation.
   This option can be used when the selected advection scheme is diffusive enough (MUSCL scheme for example).
-\item[{\np[=.true.]{ln_traldf_lap}{ln\_traldf\_lap}}]
+\item [{\np[=.true.]{ln_traldf_lap}{ln\_traldf\_lap}}]
   a laplacian operator is selected.
   This harmonic operator takes the following expression:  $\mathcal{L}(T) = \nabla \cdot A_{ht} \; \nabla T $,
   where the gradient operates along the selected direction (see \autoref{subsec:TRA_ldf_dir}),
   and $A_{ht}$ is the eddy diffusivity coefficient expressed in $m^2/s$ (see \autoref{chap:LDF}).
-\item[{\np[=.true.]{ln_traldf_blp}{ln\_traldf\_blp}}]:
+\item [{\np[=.true.]{ln_traldf_blp}{ln\_traldf\_blp}}]:
   a bilaplacian operator is selected.
   This biharmonic operator takes the following expression:
@@ -489,7 +456,4 @@
 whereas the laplacian damping time scales only like $\lambda^{-2}$.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Direction of action
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Action direction (\forcode{ln_traldf_}\{\forcode{lev,hor,iso,triad}\})]{Direction of action (\protect\np{ln_traldf_lev}{ln\_traldf\_lev}, \protect\np{ln_traldf_hor}{ln\_traldf\_hor}, \protect\np{ln_traldf_iso}{ln\_traldf\_iso}, or \protect\np{ln_traldf_triad}{ln\_traldf\_triad})}
 \label{subsec:TRA_ldf_dir}
@@ -515,7 +479,4 @@
 the next two sub-sections.
 
-% -------------------------------------------------------------------------------------------------------------
-%       iso-level operator
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Iso-level (bi-)laplacian operator (\forcode{ln_traldf_iso})]{Iso-level (bi-)laplacian operator ( \protect\np{ln_traldf_iso}{ln\_traldf\_iso})}
 \label{subsec:TRA_ldf_lev}
@@ -546,7 +507,4 @@
 They are calculated in the \mdl{zpshde} module, described in \autoref{sec:TRA_zpshde}.
 
-% -------------------------------------------------------------------------------------------------------------
-%         Rotated laplacian operator
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Standard and triad (bi-)laplacian operator}
 \label{subsec:TRA_ldf_iso_triad}
@@ -626,7 +584,4 @@
 \end{itemize}
 
-% ================================================================
-% Tracer Vertical Diffusion
-% ================================================================
 \section[Tracer vertical diffusion (\textit{trazdf.F90})]{Tracer vertical diffusion (\protect\mdl{trazdf})}
 \label{sec:TRA_zdf}
@@ -663,13 +618,7 @@
 it overcomes the stability constraint.
 
-% ================================================================
-% External Forcing
-% ================================================================
 \section{External forcing}
 \label{sec:TRA_sbc_qsr_bbc}
 
-% -------------------------------------------------------------------------------------------------------------
-%        surface boundary condition
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Surface boundary condition (\textit{trasbc.F90})]{Surface boundary condition (\protect\mdl{trasbc})}
 \label{subsec:TRA_sbc}
@@ -692,19 +641,14 @@
 
 \begin{itemize}
-\item
-  $Q_{ns}$, the non-solar part of the net surface heat flux that crosses the sea surface
+\item $Q_{ns}$, the non-solar part of the net surface heat flux that crosses the sea surface
   (\ie\ the difference between the total surface heat flux and the fraction of the short wave flux that
   penetrates into the water column, see \autoref{subsec:TRA_qsr})
   plus the heat content associated with of the mass exchange with the atmosphere and lands.
-\item
-  $\textit{sfx}$, the salt flux resulting from ice-ocean mass exchange (freezing, melting, ridging...)
-\item
-  \textit{emp}, the mass flux exchanged with the atmosphere (evaporation minus precipitation) and
+\item $\textit{sfx}$, the salt flux resulting from ice-ocean mass exchange (freezing, melting, ridging...)
+\item \textit{emp}, the mass flux exchanged with the atmosphere (evaporation minus precipitation) and
   possibly with the sea-ice and ice-shelves.
-\item
-  \textit{rnf}, the mass flux associated with runoff
+\item \textit{rnf}, the mass flux associated with runoff
   (see \autoref{sec:SBC_rnf} for further detail of how it acts on temperature and salinity tendencies)
-\item
-  \textit{fwfisf}, the mass flux associated with ice shelf melt,
+\item \textit{fwfisf}, the mass flux associated with ice shelf melt,
   (see \autoref{sec:SBC_isf} for further details on how the ice shelf melt is computed and applied).
 \end{itemize}
@@ -742,7 +686,4 @@
 This is the reason why the modified filter is not applied in the linear free surface case (see \autoref{chap:TD}).
 
-% -------------------------------------------------------------------------------------------------------------
-%        Solar Radiation Penetration
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Solar radiation penetration (\textit{traqsr.F90})]{Solar radiation penetration (\protect\mdl{traqsr})}
 \label{subsec:TRA_qsr}
@@ -823,13 +764,13 @@
 
 \begin{description}
-\item[{\np[=0]{nn_chldta}{nn\_chldta}}]
+\item [{\np[=0]{nn_chldta}{nn\_chldta}}]
   a constant 0.05 g.Chl/L value everywhere ;
-\item[{\np[=1]{nn_chldta}{nn\_chldta}}]
+\item [{\np[=1]{nn_chldta}{nn\_chldta}}]
   an observed time varying chlorophyll deduced from satellite surface ocean color measurement spread uniformly in
   the vertical direction;
-\item[{\np[=2]{nn_chldta}{nn\_chldta}}]
+\item [{\np[=2]{nn_chldta}{nn\_chldta}}]
   same as previous case except that a vertical profile of chlorophyl is used.
   Following \cite{morel.berthon_LO89}, the profile is computed from the local surface chlorophyll value;
-\item[{\np[=.true.]{ln_qsr_bio}{ln\_qsr\_bio}}]
+\item [{\np[=.true.]{ln_qsr_bio}{ln\_qsr\_bio}}]
   simulated time varying chlorophyll by TOP biogeochemical model.
   In this case, the RGB formulation is used to calculate both the phytoplankton light limitation in
@@ -849,5 +790,4 @@
 (\ie\ $I$ is masked).
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -863,9 +803,5 @@
   \label{fig:TRA_qsr_irradiance}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        Bottom Boundary Condition
-% -------------------------------------------------------------------------------------------------------------
+
 \subsection[Bottom boundary condition (\textit{trabbc.F90}) - \forcode{ln_trabbc})]{Bottom boundary condition (\protect\mdl{trabbc} - \protect\np{ln_trabbc}{ln\_trabbc})}
 \label{subsec:TRA_bbc}
@@ -878,5 +814,4 @@
 \end{listing}
 %--------------------------------------------------------------------------------------------------------------
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -887,5 +822,4 @@
   \label{fig:TRA_geothermal}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 Usually it is assumed that there is no exchange of heat or salt through the ocean bottom,
@@ -905,7 +839,4 @@
 the \ifile{geothermal\_heating} NetCDF file (\autoref{fig:TRA_geothermal}) \citep{emile-geay.madec_OS09}.
 
-% ================================================================
-% Bottom Boundary Layer
-% ================================================================
 \section[Bottom boundary layer (\textit{trabbl.F90} - \forcode{ln_trabbl})]{Bottom boundary layer (\protect\mdl{trabbl} - \protect\np{ln_trabbl}{ln\_trabbl})}
 \label{sec:TRA_bbl}
@@ -941,7 +872,4 @@
 \citet{campin.goosse_T99}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Diffusive BBL
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Diffusive bottom boundary layer (\forcode{nn_bbl_ldf=1})]{Diffusive bottom boundary layer (\protect\np[=1]{nn_bbl_ldf}{nn\_bbl\_ldf})}
 \label{subsec:TRA_bbl_diff}
@@ -980,7 +908,4 @@
 $\overline H^\sigma$, the along bottom mean temperature, salinity and depth, respectively.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Advective BBL
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Advective bottom boundary layer (\forcode{nn_bbl_adv=1,2})]{Advective bottom boundary layer (\protect\np[=1,2]{nn_bbl_adv}{nn\_bbl\_adv})}
 \label{subsec:TRA_bbl_adv}
@@ -991,5 +916,4 @@
 %}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -1006,5 +930,4 @@
   \label{fig:TRA_bbl}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 %!!      nn_bbl_adv = 1   use of the ocean velocity as bbl velocity
@@ -1071,7 +994,4 @@
 It has to be used to compute the effective velocity as well as the effective overturning circulation.
 
-% ================================================================
-% Tracer damping
-% ================================================================
 \section[Tracer damping (\textit{tradmp.F90})]{Tracer damping (\protect\mdl{tradmp})}
 \label{sec:TRA_dmp}
@@ -1130,7 +1050,4 @@
 \path{./tools/DMP_TOOLS}.
 
-% ================================================================
-% Tracer time evolution
-% ================================================================
 \section[Tracer time evolution (\textit{tranxt.F90})]{Tracer time evolution (\protect\mdl{tranxt})}
 \label{sec:TRA_nxt}
@@ -1165,7 +1082,4 @@
 $T^{t - \rdt} = T^t$ and $T^t = T_f$.
 
-% ================================================================
-% Equation of State (eosbn2)
-% ================================================================
 \section[Equation of state (\textit{eosbn2.F90})]{Equation of state (\protect\mdl{eosbn2})}
 \label{sec:TRA_eosbn2}
@@ -1179,10 +1093,6 @@
 %--------------------------------------------------------------------------------------------------------------
 
-% -------------------------------------------------------------------------------------------------------------
-%        Equation of State
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Equation of seawater (\forcode{ln_}\{\forcode{teos10,eos80,seos}\})]{Equation of seawater (\protect\np{ln_teos10}{ln\_teos10}, \protect\np{ln_teos80}{ln\_teos80}, or \protect\np{ln_seos}{ln\_seos})}
 \label{subsec:TRA_eos}
-
 
 The Equation Of Seawater (EOS) is an empirical nonlinear thermodynamic relationship linking seawater density,
@@ -1216,5 +1126,5 @@
 
 \begin{description}
-\item[{\np[=.true.]{ln_teos10}{ln\_teos10}}]
+\item [{\np[=.true.]{ln_teos10}{ln\_teos10}}]
   the polyTEOS10-bsq equation of seawater \citep{roquet.madec.ea_OM15} is used.
   The accuracy of this approximation is comparable to the TEOS-10 rational function approximation,
@@ -1235,5 +1145,5 @@
   either computing the air-sea and ice-sea fluxes (forced mode) or
   sending the SST field to the atmosphere (coupled mode).
-\item[{\np[=.true.]{ln_eos80}{ln\_eos80}}]
+\item [{\np[=.true.]{ln_eos80}{ln\_eos80}}]
   the polyEOS80-bsq equation of seawater is used.
   It takes the same polynomial form as the polyTEOS10, but the coefficients have been optimized to
@@ -1247,5 +1157,5 @@
   Nevertheless, a severe assumption is made in order to have a heat content ($C_p T_p$) which
   is conserved by the model: $C_p$ is set to a constant value, the TEOS10 value.
-\item[{\np[=.true.]{ln_seos}{ln\_seos}}]
+\item [{\np[=.true.]{ln_seos}{ln\_seos}}]
   a simplified EOS (S-EOS) inspired by \citet{vallis_bk06} is chosen,
   the coefficients of which has been optimized to fit the behavior of TEOS10
@@ -1277,5 +1187,4 @@
 \end{description}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{table}[!tb]
   \centering
@@ -1302,9 +1211,5 @@
   \label{tab:TRA_SEOS}
 \end{table}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        Brunt-V\"{a}is\"{a}l\"{a} Frequency
-% -------------------------------------------------------------------------------------------------------------
+
 \subsection[Brunt-V\"{a}is\"{a}l\"{a} frequency]{Brunt-V\"{a}is\"{a}l\"{a} frequency}
 \label{subsec:TRA_bn2}
@@ -1327,7 +1232,4 @@
 They are computed through \textit{eos\_rab}, a \fortran\ function that can be found in \mdl{eosbn2}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Freezing Point of Seawater
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Freezing point of seawater}
 \label{subsec:TRA_fzp}
@@ -1349,7 +1251,4 @@
 a \fortran\ function that can be found in \mdl{eosbn2}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Potential Energy
-% -------------------------------------------------------------------------------------------------------------
 %\subsection{Potential Energy anomalies}
 %\label{subsec:TRA_bn2}
@@ -1358,7 +1257,4 @@
 %
 
-% ================================================================
-% Horizontal Derivative in zps-coordinate
-% ================================================================
 \section[Horizontal derivative in \textit{zps}-coordinate (\textit{zpshde.F90})]{Horizontal derivative in \textit{zps}-coordinate (\protect\mdl{zpshde})}
 \label{sec:TRA_zpshde}
@@ -1380,5 +1276,4 @@
 For example, for temperature in the $i$-direction the needed interpolated temperature, $\widetilde T$, is:
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!p]
   \centering
@@ -1398,5 +1293,4 @@
   \label{fig:TRA_Partial_step_scheme}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \[
   \widetilde T = \lt\{
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_ZDF.tex	(revision 11596)
@@ -5,7 +5,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter  Vertical Ocean Physics (ZDF)
-% ================================================================
 \chapter{Vertical Ocean Physics (ZDF)}
 \label{chap:ZDF}
@@ -15,9 +12,4 @@
 %gm% Add here a small introduction to ZDF and naming of the different physics (similar to what have been written for TRA and DYN.
 
-\newpage
-
-% ================================================================
-% Vertical Mixing
-% ================================================================
 \section{Vertical mixing}
 \label{sec:ZDF}
@@ -55,7 +47,4 @@
 %--------------------------------------------------------------------------------------------------------------
 
-% -------------------------------------------------------------------------------------------------------------
-%        Constant
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Constant (\forcode{ln_zdfcst})]{Constant (\protect\np{ln_zdfcst}{ln\_zdfcst})}
 \label{subsec:ZDF_cst}
@@ -77,7 +66,4 @@
 $\sim10^{-9}~m^2.s^{-1}$ for salinity.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Richardson Number Dependent
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Richardson number dependent (\forcode{ln_zdfric})]{Richardson number dependent (\protect\np{ln_zdfric}{ln\_zdfric})}
 \label{subsec:ZDF_ric}
@@ -138,7 +124,4 @@
 the empirical values \np{rn_wtmix}{rn\_wtmix} and \np{rn_wvmix}{rn\_wvmix} \citep{lermusiaux_JMS01}.
 
-% -------------------------------------------------------------------------------------------------------------
-%        TKE Turbulent Closure Scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[TKE turbulent closure scheme (\forcode{ln_zdftke})]{TKE turbulent closure scheme (\protect\np{ln_zdftke}{ln\_zdftke})}
 \label{subsec:ZDF_tke}
@@ -248,5 +231,4 @@
 evaluate the dissipation and mixing length scales as
 (and note that here we use numerical indexing):
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -255,5 +237,4 @@
   \label{fig:ZDF_mixing_length}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \[
   % \label{eq:ZDF_tke_mxl2}
@@ -421,7 +402,4 @@
 % (\eg\ Mellor, 1989; Large et al., 1994; Meier, 2001; Axell, 2002; St. Laurent and Garrett, 2002).
 
-% -------------------------------------------------------------------------------------------------------------
-%        GLS Generic Length Scale Scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[GLS: Generic Length Scale (\forcode{ln_zdfgls})]{GLS: Generic Length Scale (\protect\np{ln_zdfgls}{ln\_zdfgls})}
 \label{subsec:ZDF_gls}
@@ -544,8 +522,4 @@
  in \citet{reffray.guillaume.ea_GMD15} for the \NEMO\ model.
 
-
-% -------------------------------------------------------------------------------------------------------------
-%        OSM OSMOSIS BL Scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[OSM: OSMosis boundary layer scheme (\forcode{ln_zdfosm})]{OSM: OSMosis boundary layer scheme (\protect\np{ln_zdfosm}{ln\_zdfosm})}
 \label{subsec:ZDF_osm}
@@ -561,11 +535,7 @@
 The OSMOSIS turbulent closure scheme is based on......   TBC
 
-% -------------------------------------------------------------------------------------------------------------
-%        TKE and GLS discretization considerations
-% -------------------------------------------------------------------------------------------------------------
 \subsection[ Discrete energy conservation for TKE and GLS schemes]{Discrete energy conservation for TKE and GLS schemes}
 \label{subsec:ZDF_tke_ene}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -576,5 +546,4 @@
   \label{fig:ZDF_TKE_time_scheme}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The production of turbulence by vertical shear (the first term of the right hand side of
@@ -666,7 +635,4 @@
 %For the latter, it is in fact the ratio $\sqrt{\bar{e}}/l_\epsilon$ which is stored.
 
-% ================================================================
-% Convection
-% ================================================================
 \section{Convection}
 \label{sec:ZDF_conv}
@@ -679,11 +645,7 @@
 or/and the use of a turbulent closure scheme.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Non-Penetrative Convective Adjustment
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Non-penetrative convective adjustment (\forcode{ln_tranpc})]{Non-penetrative convective adjustment (\protect\np{ln_tranpc}{ln\_tranpc})}
 \label{subsec:ZDF_npc}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!htb]
   \centering
@@ -704,5 +666,4 @@
   \label{fig:ZDF_npc}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 Options are defined through the \nam{zdf}{zdf} namelist variables.
@@ -744,7 +705,4 @@
 having to recompute the expansion coefficients at each mixing iteration.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Enhanced Vertical Diffusion
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Enhanced vertical diffusion (\forcode{ln_zdfevd})]{Enhanced vertical diffusion (\protect\np{ln_zdfevd}{ln\_zdfevd})}
 \label{subsec:ZDF_evd}
@@ -770,10 +728,6 @@
 a leapfrog environment \citep{leclair_phd10} (see \autoref{sec:TD_mLF}).
 
-% -------------------------------------------------------------------------------------------------------------
-%       Turbulent Closure Scheme
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Handling convection with turbulent closure schemes (\forcode{ln_zdf_}\{\forcode{tke,gls,osm}\})]{Handling convection with turbulent closure schemes (\forcode{ln_zdf{tke,gls,osm}})}
 \label{subsec:ZDF_tcs}
-
 
 The turbulent closure schemes presented in \autoref{subsec:ZDF_tke}, \autoref{subsec:ZDF_gls} and
@@ -798,10 +752,6 @@
 % gm%  + one word on non local flux with KPP scheme trakpp.F90 module...
 
-% ================================================================
-% Double Diffusion Mixing
-% ================================================================
 \section[Double diffusion mixing (\forcode{ln_zdfddm})]{Double diffusion mixing (\protect\np{ln_zdfddm}{ln\_zdfddm})}
 \label{subsec:ZDF_ddm}
-
 
 %-------------------------------------------namzdf_ddm-------------------------------------------------
@@ -818,5 +768,4 @@
 it leads to relatively minor changes in circulation but exerts significant regional influences on
 temperature and salinity.
-
 
 Diapycnal mixing of S and T are described by diapycnal diffusion coefficients
@@ -844,5 +793,4 @@
 \end{align}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -861,5 +809,4 @@
   \label{fig:ZDF_ddm}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The factor 0.7 in \autoref{eq:ZDF_ddm_f_T} reflects the measured ratio $\alpha F_T /\beta F_S \approx  0.7$ of
@@ -893,7 +840,4 @@
 This avoids duplication in the computation of $\alpha$ and $\beta$ (which is usually quite expensive).
 
-% ================================================================
-% Bottom Friction
-% ================================================================
 \section[Bottom and top friction (\textit{zdfdrg.F90})]{Bottom and top friction (\protect\mdl{zdfdrg})}
 \label{sec:ZDF_drg}
@@ -924,5 +868,4 @@
 As the friction processes at the top and the bottom are treated in and identical way,
 the description below considers mostly the bottom friction case, if not stated otherwise.
-
 
 Both the surface momentum flux (wind stress) and the bottom momentum flux (bottom friction) enter the equations as
@@ -973,7 +916,4 @@
 Note than from \NEMO\ 4.0, drag coefficients are only computed at cell centers (\ie\ at T-points) and refer to as $c_b^T$ in the following. These are then linearly interpolated in space to get $c_b^\textbf{U}$ at velocity points.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Linear Bottom Friction
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Linear top/bottom friction (\forcode{ln_lin})]{Linear top/bottom friction (\protect\np{ln_lin}{ln\_lin})}
 \label{subsec:ZDF_drg_linear}
@@ -1012,7 +952,4 @@
 $mask\_value$ * \np{rn_boost}{rn\_boost} * \np{rn_Cd0}{rn\_Cd0}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Non-Linear Bottom Friction
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Non-linear top/bottom friction (\forcode{ln_non_lin})]{Non-linear top/bottom friction (\protect\np{ln_non_lin}{ln\_non\_lin})}
 \label{subsec:ZDF_drg_nonlinear}
@@ -1047,7 +984,4 @@
 $mask\_value$ * \np{rn_boost}{rn\_boost} * \np{rn_Cd0}{rn\_Cd0}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Bottom Friction Log-layer
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Log-layer top/bottom friction (\forcode{ln_loglayer})]{Log-layer top/bottom friction (\protect\np{ln_loglayer}{ln\_loglayer})}
 \label{subsec:ZDF_drg_loglayer}
@@ -1073,7 +1007,4 @@
 %In this case, the relevant namelist parameters are \np{rn_tfrz0}{rn\_tfrz0}, \np{rn_tfri2}{rn\_tfri2} and \np{rn_tfri2_max}{rn\_tfri2\_max}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Explicit bottom Friction
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Explicit top/bottom friction (\forcode{ln_drgimp=.false.})]{Explicit top/bottom friction (\protect\np[=.false.]{ln_drgimp}{ln\_drgimp})}
 \label{subsec:ZDF_drg_stability}
@@ -1134,7 +1065,4 @@
 The number of potential breaches of the explicit stability criterion are still reported for information purposes.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Implicit Bottom Friction
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Implicit top/bottom friction (\forcode{ln_drgimp=.true.})]{Implicit top/bottom friction (\protect\np[=.true.]{ln_drgimp}{ln\_drgimp})}
 \label{subsec:ZDF_drg_imp}
@@ -1164,7 +1092,4 @@
 Superscript $n+1$ means the velocity used in the friction formula is to be calculated, so it is implicit.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Bottom Friction with split-explicit free surface
-% -------------------------------------------------------------------------------------------------------------
 \subsection[Bottom friction with split-explicit free surface]{Bottom friction with split-explicit free surface}
 \label{subsec:ZDF_drg_ts}
@@ -1180,8 +1105,4 @@
 Note that other strategies are possible, like considering vertical diffusion step in advance, \ie\ prior barotropic integration.
 
-
-% ================================================================
-% Internal wave-driven mixing
-% ================================================================
 \section[Internal wave-driven mixing (\forcode{ln_zdfiwm})]{Internal wave-driven mixing (\protect\np{ln_zdfiwm}{ln\_zdfiwm})}
 \label{subsec:ZDF_tmx_new}
@@ -1245,7 +1166,4 @@
 % Jc: input files names ?
 
-% ================================================================
-% surface wave-induced mixing
-% ================================================================
 \section[Surface wave-induced mixing (\forcode{ln_zdfswm})]{Surface wave-induced mixing (\protect\np{ln_zdfswm}{ln\_zdfswm})}
 \label{subsec:ZDF_swm}
@@ -1278,7 +1196,4 @@
 (for more information on wave parameters and settings see \autoref{sec:SBC_wave})
 
-% ================================================================
-% Adaptive-implicit vertical advection
-% ================================================================
 \section[Adaptive-implicit vertical advection (\forcode{ln_zad_Aimp})]{Adaptive-implicit vertical advection(\protect\np{ln_zad_Aimp}{ln\_zad\_Aimp})}
 \label{subsec:ZDF_aimp}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_cfgs.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_cfgs.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_cfgs.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter Configurations
-% ================================================================
 \chapter{Configurations}
 \label{chap:CFGS}
@@ -10,9 +7,4 @@
 \chaptertoc
 
-\newpage
-
-% ================================================================
-% Introduction
-% ================================================================
 \section{Introduction}
 \label{sec:CFGS_intro}
@@ -35,7 +27,4 @@
 %-------------------------------------------------------------------------------------------------------------
 
-% ================================================================
-% 1D model configuration
-% ================================================================
 \section[C1D: 1D Water column model (\texttt{\textbf{key\_c1d}})]{C1D: 1D Water column model (\protect\key{c1d})}
 \label{sec:CFGS_c1d}
@@ -60,11 +49,11 @@
 Therefore, defining \key{c1d} changes some things in the code behaviour:
 \begin{description}
-\item[(1)]
+\item [(1)]
   a simplified \rou{stp} routine is used (\rou{stp\_c1d}, see \mdl{step\_c1d} module) in which
   both lateral tendancy terms and lateral physics are not called;
-\item[(2)]
+\item [(2)]
   the vertical velocity is zero
   (so far, no attempt at introducing a Ekman pumping velocity has been made);
-\item[(3)]
+\item [(3)]
   a simplified treatment of the Coriolis term is performed as $U$- and $V$-points are the same
   (see \mdl{dyncor\_c1d}).
@@ -75,7 +64,4 @@
 % to be added:  a test case on the yearlong Ocean Weather Station (OWS) Papa dataset of Martin (1985)
 
-% ================================================================
-% ORCA family configurations
-% ================================================================
 \section{ORCA family: global ocean with tripolar grid}
 \label{sec:CFGS_orca}
@@ -90,5 +76,4 @@
 In this namelist\_cfg the name of domain input file is set in \nam{cfg}{cfg} block of namelist.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -104,9 +89,5 @@
   \label{fig:CFGS_ORCA_msh}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%       ORCA tripolar grid
-% -------------------------------------------------------------------------------------------------------------
+
 \subsection{ORCA tripolar grid}
 \label{subsec:CFGS_orca_grid}
@@ -121,5 +102,4 @@
 The resulting mesh presents no loss of continuity in either the mesh lines or the scale factors,
 or even the scale factor derivatives over the whole ocean domain, as the mesh is not a composite mesh.
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tbp]
   \centering
@@ -136,5 +116,4 @@
   \label{fig:CFGS_ORCA_e1e2}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 The method is applied to Mercator grid (\ie\ same zonal and meridional grid spacing) poleward of 20\deg{N},
@@ -149,7 +128,4 @@
 while the ratio of anisotropy remains close to one except near the Victoria Island in the Canadian Archipelago.
 
-% -------------------------------------------------------------------------------------------------------------
-%       ORCA-ICE(-PISCES) configurations
-% -------------------------------------------------------------------------------------------------------------
 \subsection{ORCA pre-defined resolution}
 \label{subsec:CFGS_orca_resolution}
@@ -182,5 +158,4 @@
 %--------------------------------------------------------------------------------------------------------------
 
-
 The ORCA\_R2 configuration has the following specificity: starting from a 2\deg\ ORCA mesh,
 local mesh refinements were applied to the Mediterranean, Red, Black and Caspian Seas,
@@ -221,7 +196,4 @@
 sponge layers at open boundaries.
 
-% -------------------------------------------------------------------------------------------------------------
-%       GYRE family: double gyre basin
-% -------------------------------------------------------------------------------------------------------------
 \section{GYRE family: double gyre basin}
 \label{sec:CFGS_gyre}
@@ -274,5 +246,4 @@
 namelist \path{./cfgs/GYRE_PISCES/EXPREF/namelist_cfg}.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -283,9 +254,5 @@
   \label{fig:CFGS_GYRE}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%       AMM configuration
-% -------------------------------------------------------------------------------------------------------------
+
 \section{AMM: atlantic margin configuration}
 \label{sec:CFGS_config_AMM}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_conservation.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_conservation.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_conservation.tex	(revision 11596)
@@ -3,7 +3,4 @@
 \begin{document}
 
-% ================================================================
-% Invariant of the Equations
-% ================================================================
 \chapter{Invariants of the Primitive Equations}
 \label{chap:CONS}
@@ -42,7 +39,4 @@
 \citep{Marti1992?, Levy1996?, Levy1998?}.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation Properties on Ocean Dynamics
-% -------------------------------------------------------------------------------------------------------------
 \section{Conservation properties on ocean dynamics}
 \label{sec:CONS_Invariant_dyn}
@@ -158,7 +152,4 @@
 otherwise there is no guarantee that the surface pressure force work vanishes.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation Properties on Ocean Thermodynamics
-% -------------------------------------------------------------------------------------------------------------
 \section{Conservation properties on ocean thermodynamics}
 \label{sec:CONS_Invariant_tra}
@@ -179,7 +170,4 @@
 In practice, the mass is conserved with a very good accuracy.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation Properties on Momentum Physics
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Conservation properties on momentum physics}
 \label{subsec:CONS_Invariant_dyn_physics}
@@ -223,5 +211,4 @@
           {A^{lm}\;\zeta \;{\mathrm {\mathbf k}}} \right)} \right]\;dv} \leqslant 0
 \]
-
 
 (II.4.6a) and (II.4.6b) means that the horizontal diffusion of momentum conserve both the potential vorticity and
@@ -286,7 +273,4 @@
 \ie\ the vertical momentum physics conserve momentum, potential vorticity, and horizontal divergence.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Conservation Properties on Tracer Physics
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Conservation properties on tracer physics}
 \label{subsec:CONS_Invariant_tra_physics}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_misc.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_misc.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_misc.tex	(revision 11596)
@@ -2,7 +2,4 @@
 
 \begin{document}
-% ================================================================
-% Chapter --- Miscellaneous Topics
-% ================================================================
 \chapter{Miscellaneous Topics}
 \label{chap:MISC}
@@ -10,9 +7,4 @@
 \chaptertoc
 
-\newpage
-
-% ================================================================
-% Representation of Unresolved Straits
-% ================================================================
 \section{Representation of unresolved straits}
 \label{sec:MISC_strait}
@@ -36,7 +28,4 @@
 lateral friction.
 
-% -------------------------------------------------------------------------------------------------------------
-%       Hand made geometry changes
-% -------------------------------------------------------------------------------------------------------------
 \subsection{Hand made geometry changes}
 \label{subsec:MISC_strait_hand}
@@ -83,5 +72,4 @@
 \texttt{fmask} for any other configuration.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!tbp]
   \centering
@@ -100,7 +88,5 @@
   \label{fig:MISC_strait_hand}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
+
 \begin{figure}[!tbp]
   \centering
@@ -116,9 +102,5 @@
   \label{fig:MISC_closea_mask_example}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% ================================================================
-% Closed seas
-% ================================================================
+
 \section[Closed seas (\textit{closea.F90})]{Closed seas (\protect\mdl{closea})}
 \label{sec:MISC_closea}
@@ -141,8 +123,8 @@
 
 \begin{enumerate}
-\item{{\bfseries No ``closea\_mask'' field is included in domain configuration
+\item {{\bfseries No ``closea\_mask'' field is included in domain configuration
   file.} In this case the closea module does nothing.}
 
-\item{{\bfseries A field called closea\_mask is included in the domain
+\item {{\bfseries A field called closea\_mask is included in the domain
 configuration file and ln\_closea=.false. in namelist namcfg.} In this
 case the inland seas defined by the closea\_mask field are filled in
@@ -150,5 +132,5 @@
 closea\_mask that is nonzero is set to be a land point.}
 
-\item{{\bfseries A field called closea\_mask is included in the domain
+\item {{\bfseries A field called closea\_mask is included in the domain
 configuration file and ln\_closea=.true. in namelist namcfg.} Each
 inland sea or group of inland seas is set to a positive integer value
@@ -159,5 +141,5 @@
 closea\_mask is zero).}
 
-\item{{\bfseries Fields called closea\_mask and closea\_mask\_rnf are
+\item {{\bfseries Fields called closea\_mask and closea\_mask\_rnf are
 included in the domain configuration file and ln\_closea=.true. in
 namelist namcfg.} This option works as for option 3, except that if
@@ -173,5 +155,5 @@
 ocean.}
 
-\item{{\bfseries Fields called closea\_mask and closea\_mask\_emp are
+\item {{\bfseries Fields called closea\_mask and closea\_mask\_emp are
 included in the domain configuration file and ln\_closea=.true. in
 namelist namcfg.} This option works the same as option 4 except that
@@ -185,7 +167,4 @@
 them to the domain configuration file in the utils/tools/DOMAINcfg directory.
 
-% ================================================================
-% Sub-Domain Functionality
-% ================================================================
 \section{Sub-domain functionality}
 \label{sec:MISC_zoom}
@@ -213,5 +192,5 @@
 
 \begin{itemize}
-\item  Add the new attribute to any input files requiring a j-row offset, i.e:
+\item Add the new attribute to any input files requiring a j-row offset, i.e:
 \begin{cmds}
 ncatted  -a open_ocean_jstart,global,a,d,41 eORCA1_domcfg.nc
@@ -251,7 +230,4 @@
 conditions. Experimenting with this remains an exercise for the user.
 
-% ================================================================
-% Accuracy and Reproducibility
-% ================================================================
 \section[Accuracy and reproducibility (\textit{lib\_fortran.F90})]{Accuracy and reproducibility (\protect\mdl{lib\_fortran})}
 \label{sec:MISC_fortran}
@@ -281,5 +257,4 @@
 We use a CPP key as the overwritting of a intrinsic function can present performance issues with
 some computers/compilers.
-
 
 \subsection{MPP reproducibility}
@@ -336,7 +311,4 @@
 non-reference configuration.
 
-% ================================================================
-% Model optimisation, Control Print and Benchmark
-% ================================================================
 \section{Model optimisation, control print and benchmark}
 \label{sec:MISC_opt}
@@ -367,9 +339,9 @@
 
 \begin{enumerate}
-\item{\np{ln_ctl}{ln\_ctl}: compute and print the trends averaged over the interior domain in all TRA, DYN, LDF and
+\item {\np{ln_ctl}{ln\_ctl}: compute and print the trends averaged over the interior domain in all TRA, DYN, LDF and
 ZDF modules.
 This option is very helpful when diagnosing the origin of an undesired change in model results. }
 
-\item{also \np{ln_ctl}{ln\_ctl} but using the nictl and njctl namelist parameters to check the source of differences between
+\item {also \np{ln_ctl}{ln\_ctl} but using the nictl and njctl namelist parameters to check the source of differences between
 mono and multi processor runs.}
 \end{enumerate}
@@ -415,7 +387,2 @@
 increment also applies to the time.step file which is otherwise updated every timestep.
 
-% ================================================================
-
-\onlyinsubfile{\input{../../global/epilogue}}
-
-\end{document}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics.tex	(revision 11596)
@@ -4,7 +4,4 @@
 \begin{document}
 
-% ================================================================
-% Chapter 1  Model Basics
-% ================================================================
 \chapter{Model Basics}
 \label{chap:MB}
@@ -12,16 +9,9 @@
 \chaptertoc
 
-\newpage
-
-% ================================================================
-% Primitive Equations
-% ================================================================
+%% =================================================================================================
 \section{Primitive equations}
 \label{sec:MB_PE}
 
-% -------------------------------------------------------------------------------------------------------------
-%        Vector Invariant Formulation
-% -------------------------------------------------------------------------------------------------------------
-
+%% =================================================================================================
 \subsection{Vector invariant formulation}
 \label{subsec:MB_PE_vector}
@@ -33,17 +23,13 @@
 
 \begin{enumerate}
-\item
-  \textit{spherical Earth approximation}: the geopotential surfaces are assumed to be oblate spheriods
+\item \textit{spherical Earth approximation}: the geopotential surfaces are assumed to be oblate spheriods
   that follow the Earth's bulge; these spheroids are approximated by spheres with
   gravity locally vertical (parallel to the Earth's radius) and independent of latitude
   \citep[][section 2]{white.hoskins.ea_QJRMS05}.
-\item
-  \textit{thin-shell approximation}: the ocean depth is neglected compared to the earth's radius
-\item
-  \textit{turbulent closure hypothesis}: the turbulent fluxes
+\item \textit{thin-shell approximation}: the ocean depth is neglected compared to the earth's radius
+\item \textit{turbulent closure hypothesis}: the turbulent fluxes
   (which represent the effect of small scale processes on the large-scale)
   are expressed in terms of large-scale features
-\item
-  \textit{Boussinesq hypothesis}: density variations are neglected except in their contribution to
+\item \textit{Boussinesq hypothesis}: density variations are neglected except in their contribution to
   the buoyancy force
   \begin{equation}
@@ -51,6 +37,5 @@
     \rho = \rho \ (T,S,p)
   \end{equation}
-\item
-  \textit{Hydrostatic hypothesis}: the vertical momentum equation is reduced to a balance between
+\item \textit{Hydrostatic hypothesis}: the vertical momentum equation is reduced to a balance between
   the vertical pressure gradient and the buoyancy force
   (this removes convective processes from the initial Navier-Stokes equations and so
@@ -60,6 +45,5 @@
     \pd[p]{z} = - \rho \ g
   \end{equation}
-\item
-  \textit{Incompressibility hypothesis}: the three dimensional divergence of the velocity vector $\vect U$
+\item \textit{Incompressibility hypothesis}: the three dimensional divergence of the velocity vector $\vect U$
   is assumed to be zero.
   \begin{equation}
@@ -67,6 +51,5 @@
     \nabla \cdot \vect U = 0
   \end{equation}
- \item
-  \textit{Neglect of additional Coriolis terms}: the Coriolis terms that vary with the cosine of latitude are neglected.
+\item \textit{Neglect of additional Coriolis terms}: the Coriolis terms that vary with the cosine of latitude are neglected.
   These terms may be non-negligible where the Brunt-Vaisala frequency $N$ is small, either in the deep ocean or
   in the sub-mesoscale motions of the mixed layer, or near the equator \citep[][section 1]{white.hoskins.ea_QJRMS05}.
@@ -108,7 +91,5 @@
 Their nature and formulation are discussed in \autoref{sec:MB_zdf_ldf} and \autoref{subsec:MB_boundary_condition}.
 
-% -------------------------------------------------------------------------------------------------------------
-% Boundary condition
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Boundary conditions}
 \label{subsec:MB_boundary_condition}
@@ -128,5 +109,4 @@
 the other components of the earth system.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!ht]
   \centering
@@ -138,9 +118,7 @@
   \label{fig:MB_ocean_bc}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 \begin{description}
-\item[Land - ocean interface:]
-  the major flux between continental margins and the ocean is a mass exchange of fresh water through river runoff.
+\item [Land - ocean interface:]  the major flux between continental margins and the ocean is a mass exchange of fresh water through river runoff.
   Such an exchange modifies the sea surface salinity especially in the vicinity of major river mouths.
   It can be neglected for short range integrations but has to be taken into account for long term integrations as
@@ -148,6 +126,5 @@
   It is required in order to close the water cycle of the climate system.
   It is usually specified as a fresh water flux at the air-sea interface in the vicinity of river mouths.
-\item[Solid earth - ocean interface:]
-  heat and salt fluxes through the sea floor are small, except in special areas of little extent.
+\item [Solid earth - ocean interface:]  heat and salt fluxes through the sea floor are small, except in special areas of little extent.
   They are usually neglected in the model
   \footnote{
@@ -171,6 +148,5 @@
   $\vect D^{\vect U}$ in \autoref{eq:MB_PE_dyn}.
   It is discussed in \autoref{eq:MB_zdf}.% and Chap. III.6 to 9.
-\item[Atmosphere - ocean interface:]
-  the kinematic surface condition plus the mass flux of fresh water PE (the precipitation minus evaporation budget)
+\item [Atmosphere - ocean interface:]  the kinematic surface condition plus the mass flux of fresh water PE (the precipitation minus evaporation budget)
   leads to:
   \[
@@ -181,6 +157,5 @@
   leads to the continuity of pressure across the interface $z = \eta$.
   The atmosphere and ocean also exchange horizontal momentum (wind stress), and heat.
-\item[Sea ice - ocean interface:]
-  the ocean and sea ice exchange heat, salt, fresh water and momentum.
+\item [Sea ice - ocean interface:]  the ocean and sea ice exchange heat, salt, fresh water and momentum.
   The sea surface temperature is constrained to be at the freezing point at the interface.
   Sea ice salinity is very low ($\sim4-6 \, psu$) compared to those of the ocean ($\sim34 \, psu$).
@@ -188,13 +163,9 @@
 \end{description}
 
-% ================================================================
-% The Horizontal Pressure Gradient
-% ================================================================
+%% =================================================================================================
 \section{Horizontal pressure gradient}
 \label{sec:MB_hor_pg}
 
-% -------------------------------------------------------------------------------------------------------------
-% Pressure Formulation
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Pressure formulation}
 \label{subsec:MB_p_formulation}
@@ -228,7 +199,5 @@
 Only the free surface formulation is now described in this document (see the next sub-section).
 
-% -------------------------------------------------------------------------------------------------------------
-% Free Surface Formulation
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Free surface formulation}
 \label{subsec:MB_free_surface}
@@ -280,13 +249,9 @@
 (see \autoref{subsec:DYN_spg_ts}).
 
-% ================================================================
-% Curvilinear z-coordinate System
-% ================================================================
+%% =================================================================================================
 \section{Curvilinear \textit{z-}coordinate system}
 \label{sec:MB_zco}
 
-% -------------------------------------------------------------------------------------------------------------
-% Tensorial Formalism
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Tensorial formalism}
 \label{subsec:MB_tensorial}
@@ -338,5 +303,4 @@
   \label{fig:MB_referential}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 Since the ocean depth is far smaller than the earth's radius, $a + z$, can be replaced by $a$ in
@@ -373,7 +337,5 @@
 where $q$ is a scalar quantity and $\vect A = (a_1,a_2,a_3)$ a vector in the $(i,j,k)$ coordinates system.
 
-% -------------------------------------------------------------------------------------------------------------
-% Continuous Model Equations
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Continuous model equations}
 \label{subsec:MB_zco_Eq}
@@ -496,6 +458,5 @@
 
 \begin{itemize}
-\item
-  \textbf{Vector invariant form of the momentum equations}:
+\item \textbf{Vector invariant form of the momentum equations}:
   \begin{equation}
     \label{eq:MB_dyn_vect}
@@ -510,6 +471,5 @@
     \end{split}
   \end{equation}
-\item
-  \textbf{flux form of the momentum equations}:
+\item \textbf{flux form of the momentum equations}:
   % \label{eq:MB_dyn_flux}
   \begin{multline*}
@@ -544,6 +504,5 @@
   where the divergence of the horizontal velocity, $\chi$ is given by \autoref{eq:MB_div_Uh}.
 
-\item
-  \textbf{tracer equations}:
+\item \textbf{tracer equations}:
   \begin{equation}
   \begin{split}
@@ -562,9 +521,5 @@
 are discussed in \autoref{chap:SBC}.
 
-\newpage
-
-% ================================================================
-% Curvilinear generalised vertical coordinate System
-% ================================================================
+%% =================================================================================================
 \section{Curvilinear generalised vertical coordinate system}
 \label{sec:MB_gco}
@@ -647,7 +602,5 @@
 %}
 
-% -------------------------------------------------------------------------------------------------------------
-% The s-coordinate Formulation
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{\textit{S}-coordinate formulation}
 
@@ -737,11 +690,8 @@
 }
 
-% -------------------------------------------------------------------------------------------------------------
-% Curvilinear \zstar-coordinate System
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Curvilinear \zstar-coordinate system}
 \label{subsec:MB_zco_star}
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!b]
   \centering
@@ -754,5 +704,4 @@
   \label{fig:MB_z_zstar}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 
 In this case, the free surface equation is nonlinear, and the variations of volume are fully taken into account.
@@ -827,15 +776,9 @@
 %end MOM doc %%%
 
-\newpage
-
-% -------------------------------------------------------------------------------------------------------------
-% Terrain following  coordinate System
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Curvilinear terrain-following \textit{s}--coordinate}
 \label{subsec:MB_sco}
 
-% -------------------------------------------------------------------------------------------------------------
-% Introduction
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsubsection{Introduction}
 
@@ -918,7 +861,5 @@
 It also offers a completely general transformation, $s=s(i,j,z)$ for the vertical coordinate.
 
-% -------------------------------------------------------------------------------------------------------------
-% Curvilinear z-tilde coordinate System
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{\texorpdfstring{Curvilinear \ztilde-coordinate}{}}
 \label{subsec:MB_zco_tilde}
@@ -929,9 +870,5 @@
 Its use is therefore not recommended.
 
-\newpage
-
-% ================================================================
-% Subgrid Scale Physics
-% ================================================================
+%% =================================================================================================
 \section{Subgrid scale physics}
 \label{sec:MB_zdf_ldf}
@@ -957,7 +894,5 @@
 The formulation of these terms and their underlying physics are briefly discussed in the next two subsections.
 
-% -------------------------------------------------------------------------------------------------------------
-% Vertical Subgrid Scale Physics
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Vertical subgrid scale physics}
 \label{subsec:MB_zdf}
@@ -991,7 +926,5 @@
 The choices available in \NEMO\ are discussed in \autoref{chap:ZDF}).
 
-% -------------------------------------------------------------------------------------------------------------
-% Lateral Diffusive and Viscous Operators Formulation
-% -------------------------------------------------------------------------------------------------------------
+%% =================================================================================================
 \subsection{Formulation of the lateral diffusive and viscous operators}
 \label{subsec:MB_ldf}
@@ -1047,4 +980,5 @@
 and UBS advection schemes when flux form is chosen for the momentum advection.
 
+%% =================================================================================================
 \subsubsection{Lateral laplacian tracer diffusive operator}
 
@@ -1088,4 +1022,5 @@
 while in $s$-coordinates $\pd[]{k}$ is replaced by $\pd[]{s}$.
 
+%% =================================================================================================
 \subsubsection{Eddy induced velocity}
 
@@ -1124,4 +1059,5 @@
 The latter strategy is used in \NEMO\ (cf. \autoref{chap:LDF}).
 
+%% =================================================================================================
 \subsubsection{Lateral bilaplacian tracer diffusive operator}
 
@@ -1135,4 +1071,5 @@
 the harmonic eddy diffusion coefficient set to the square root of the biharmonic one.
 
+%% =================================================================================================
 \subsubsection{Lateral Laplacian momentum diffusive operator}
 
@@ -1167,4 +1104,5 @@
 a geographical coordinate system \citep{lengaigne.madec.ea_JGR03}.
 
+%% =================================================================================================
 \subsubsection{Lateral bilaplacian momentum diffusive operator}
 
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics_zstar.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics_zstar.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_model_basics_zstar.tex	(revision 11596)
@@ -2,16 +2,6 @@
 
 \begin{document}
-% ================================================================
-% Chapter 1 Model Basics
-% ================================================================
-% ================================================================
-% Curvilinear \zstar- \sstar-coordinate System
-% ================================================================
 \chapter{ essai \zstar \sstar}
 \section{Curvilinear \zstar- or \sstar coordinate system}
-
-% -------------------------------------------------------------------------------------------------------------
-% ????
-% -------------------------------------------------------------------------------------------------------------
 
 \colorbox{yellow}{ to be updated }
@@ -70,7 +60,4 @@
 %%%
 
-% ================================================================
-% Surface Pressure Gradient and Sea Surface Height
-% ================================================================
 \section[Surface pressure gradient and sea surface heigth (\textit{dynspg.F90})]{Surface pressure gradient and sea surface heigth (\protect\mdl{dynspg})}
 \label{sec:MBZ_dyn_hpg_spg}
Index: /NEMO/trunk/doc/latex/NEMO/subfiles/chap_time_domain.tex
===================================================================
--- /NEMO/trunk/doc/latex/NEMO/subfiles/chap_time_domain.tex	(revision 11595)
+++ /NEMO/trunk/doc/latex/NEMO/subfiles/chap_time_domain.tex	(revision 11596)
@@ -3,7 +3,4 @@
 \begin{document}
 
-% ================================================================
-% Chapter 2 ——— Time Domain (step.F90)
-% ================================================================
 \chapter{Time Domain}
 \label{chap:TD}
@@ -16,6 +13,4 @@
   would help  ==> to be added}
 %%%%
-
-\newpage
 
 Having defined the continuous equations in \autoref{chap:MB}, we need now to choose a time discretization,
@@ -25,7 +20,4 @@
 the consequences for the order in which the equations are solved.
 
-% ================================================================
-% Time Discretisation
-% ================================================================
 \section{Time stepping environment}
 \label{sec:TD_environment}
@@ -55,7 +47,4 @@
 The time stepping itself is performed once at each time step where implicit vertical diffusion is computed, \ie\ in the \mdl{trazdf} and \mdl{dynzdf} modules.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Non-Diffusive Part---Leapfrog Scheme
-% -------------------------------------------------------------------------------------------------------------
 \section{Non-diffusive part --- Leapfrog scheme}
 \label{sec:TD_leap_frog}
@@ -98,7 +87,4 @@
 filter parameter and the viscosity and diffusion coefficients.
 
-% -------------------------------------------------------------------------------------------------------------
-%        Diffusive Part---Forward or Backward Scheme
-% -------------------------------------------------------------------------------------------------------------
 \section{Diffusive part --- Forward or backward scheme}
 \label{sec:TD_forward_imp}
@@ -165,7 +151,4 @@
 (see for example \citet{richtmyer.morton_bk67}).
 
-% -------------------------------------------------------------------------------------------------------------
-%        Surface Pressure gradient
-% -------------------------------------------------------------------------------------------------------------
 \section{Surface pressure gradient}
 \label{sec:TD_spg_ts}
@@ -184,5 +167,4 @@
 
 %\gmcomment{
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -197,10 +179,6 @@
   \label{fig:TD_TimeStep_flowchart}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 %}
 
-% -------------------------------------------------------------------------------------------------------------
-%        The Modified Leapfrog -- Asselin Filter scheme
-% -------------------------------------------------------------------------------------------------------------
 \section{Modified Leapfrog -- Asselin filter scheme}
 \label{sec:TD_mLF}
@@ -245,5 +223,4 @@
 even if separated by only $\rdt$ since the time filter is no longer applied to the forcing term.
 
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
 \begin{figure}[!t]
   \centering
@@ -261,9 +238,5 @@
   \label{fig:TD_MLF_forcing}
 \end{figure}
-%>>>>>>>>>>>>>>>>>>>>>>>>>>>>
-
-% -------------------------------------------------------------------------------------------------------------
-%        Start/Restart strategy
-% -------------------------------------------------------------------------------------------------------------
+
 \section{Start/Restart strategy}
 \label{sec:TD_rst}
@@ -315,76 +288,2 @@
 %-------------------------------------------------------------------------------------------------------------
 %        Time Domain
-% -------------------------------------------------------------------------------------------------------------
-\subsection{Time domain}
-\label{subsec:TD_time}
-%--------------------------------------------namrun-------------------------------------------
-
-%--------------------------------------------------------------------------------------------------------------
-
-Options are defined through the  \nam{dom}{dom} namelist variables.
- \colorbox{yellow}{add here a few word on nit000 and nitend}
-
- \colorbox{yellow}{Write documentation on the calendar and the key variable adatrj}
-
-add a description of daymod, and the model calandar (leap-year and co)
-
-}     	%% end add
-
-
-
-%%
-\gmcomment{       % add implicit in vvl case  and Crant-Nicholson scheme
-
-Implicit time stepping in case of variable volume thickness.
-
-Tracer case (NB for momentum in vector invariant form take care!)
-
-\begin{flalign*}
-  &\frac{\lt( e_{3t}\,T \rt)_k^{t+1}-\lt( e_{3t}\,T \rt)_k^{t-1}}{2\rdt}
-  \equiv \text{RHS}+ \delta_k \lt[ {\frac{A_w^{vt} }{e_{3w}^{t+1} }\delta_{k + 1/2} \lt[ {T^{t+1}} \rt]}
-  \rt]      \\
-  &\lt( e_{3t}\,T \rt)_k^{t+1}-\lt( e_{3t}\,T \rt)_k^{t-1}
-  \equiv {2\rdt} \ \text{RHS}+ {2\rdt} \ \delta_k \lt[ {\frac{A_w^{vt} }{e_{3w}^{t+1} }\delta_{k + 1/2} \lt[ {T^{t+1}} \rt]}
-  \rt]      \\
-  &\lt( e_{3t}\,T \rt)_k^{t+1}-\lt( e_{3t}\,T \rt)_k^{t-1}
-  \equiv 2\rdt \ \text{RHS}
-  + 2\rdt \ \lt\{ \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k + 1/2} [ T_{k +1}^{t+1} - T_k      ^{t+1} ]
-    - \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k - 1/2} [ T_k       ^{t+1} - T_{k -1}^{t+1} ]  \rt\}     \\
-  &\\
-  &\lt( e_{3t}\,T \rt)_k^{t+1}
-  -  {2\rdt} \           \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k + 1/2}                  T_{k +1}^{t+1}
-  + {2\rdt} \ \lt\{  \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k + 1/2}
-    +  \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k - 1/2}     \rt\}   T_{k    }^{t+1}
-  -  {2\rdt} \           \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k - 1/2}                  T_{k -1}^{t+1}      \\
-  &\equiv \lt( e_{3t}\,T \rt)_k^{t-1} + {2\rdt} \ \text{RHS}    \\
-  %
-\end{flalign*}
-\begin{flalign*}
-  \allowdisplaybreaks
-  \intertext{ Tracer case }
-  %
-  &  \qquad \qquad  \quad   -  {2\rdt}                  \ \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k + 1/2}
-  \qquad \qquad \qquad  \qquad  T_{k +1}^{t+1}   \\
-  &+ {2\rdt} \ \biggl\{  (e_{3t})_{k   }^{t+1}  \bigg. +    \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k + 1/2}
-  +   \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k - 1/2} \bigg. \biggr\}  \ \ \ T_{k   }^{t+1}  &&\\
-  & \qquad \qquad  \qquad \qquad \qquad \quad \ \ -  {2\rdt} \                          \lt[ \frac{A_w^{vt}}{e_{3w}^{t+1}} \rt]_{k - 1/2}                          \quad \ \ T_{k -1}^{t+1}
-  \ \equiv \ \lt( e_{3t}\,T \rt)_k^{t-1} + {2\rdt} \ \text{RHS}  \\
-  %
-\end{flalign*}
-\begin{flalign*}
-  \allowdisplaybreaks
-  \intertext{ Tracer content case }
-  %
-  & -  {2\rdt} \              & \frac{(A_w^{vt})_{k + 1/2}} {(e_{3w})_{k + 1/2}^{t+1}\;(e_{3t})_{k +1}^{t+1}}  && \  \lt( e_{3t}\,T \rt)_{k +1}^{t+1}   &\\
-  & + {2\rdt} \ \lt[ 1  \rt.+ & \frac{(A_w^{vt})_{k + 1/2}} {(e_{3w})_{k + 1/2}^{t+1}\;(e_{3t})_k^{t+1}}
-  + & \frac{(A_w^{vt})_{k - 1/2}} {(e_{3w})_{k - 1/2}^{t+1}\;(e_{3t})_k^{t+1}}  \lt.  \rt]  & \lt( e_{3t}\,T \rt)_{k   }^{t+1}  &\\
-  & -  {2\rdt} \               & \frac{(A_w^{vt})_{k - 1/2}} {(e_{3w})_{k - 1/2}^{t+1}\;(e_{3t})_{k -1}^{t+1}}     &\  \lt( e_{3t}\,T \rt)_{k -1}^{t+1}
-  \equiv \lt( e_{3t}\,T \rt)_k^{t-1} + {2\rdt} \ \text{RHS}  &
-\end{flalign*}
-
-%%
-}
-
-\onlyinsubfile{\input{../../global/epilogue}}
-
-\end{document}
