Index: NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.sty
===================================================================
--- NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.sty	(revision 11041)
+++ NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.sty	(revision 11043)
@@ -6,5 +6,4 @@
 %% LaTeX packages
 %% ==============================================================================
-
 \usepackage{natbib}           %% bib
 \usepackage{caption}          %% caption
@@ -24,15 +23,8 @@
 
 %% Extensions in bundle package
-
 \usepackage{amssymb, graphicx, makeidx, tabularx}
 
-
 %% Configuration
-
 \captionsetup{margin=10pt, font={small}, labelsep=colon, labelfont={bf}}
-\hypersetup{
-	pdftitle={NEMO ocean engine}, pdfauthor={Gurvan Madec, and NEMO System Team},
-	colorlinks
-}
 \idxlayout{font=footnotesize, columns=3}
 \renewcommand{\bibfont}{\footnotesize}
@@ -42,15 +34,12 @@
 %% Styles
 %% ==============================================================================
-
 \pagestyle{fancy}
 \bibliographystyle{../../NEMO/main/ametsoc}
 
 %% Additionnal fonts
-
 \DeclareMathAlphabet{\mathpzc}{OT1}{pzc}{m}{it}
 
 
 %% Page layout
-
 \fancyhf{}
 \fancyhead[LE,RO]{\bfseries\thepage}
@@ -65,5 +54,4 @@
 
 %% Catcodes
-
 \makeatletter
 \def\LigneVerticale{\vrule height 5cm depth 2cm\hspace{0.1cm}\relax}
Index: NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.tex
===================================================================
--- NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/NEMO/main/NEMO_manual.tex	(revision 11043)
@@ -12,4 +12,8 @@
 %% Custom style (.sty)
 \usepackage{../main/NEMO_manual}
+\hypersetup{
+	pdftitle={NEMO ocean engine}, pdfauthor={Gurvan Madec, and NEMO System Team},
+	colorlinks
+}
 
 %% Include references and index for single subfile compilation
@@ -47,5 +51,4 @@
 %
 %  }                                                        \\
-%                                                           \\
   \textit{Issue 27, Notes du P\^{o}le de mod\'{e}lisation} \\
   \textit{Institut Pierre-Simon Laplace (IPSL)}            \\
@@ -56,5 +59,4 @@
 \maketitle
 \frontmatter
-
 
 %% ToC i.e. Table of Contents
@@ -123,4 +125,3 @@
 \printindex
 
-
 \end{document}
Index: NEMO/trunk/doc/latex/SI3/main/SI3_manual.bib
===================================================================
--- NEMO/trunk/doc/latex/SI3/main/SI3_manual.bib	(revision 11041)
+++ NEMO/trunk/doc/latex/SI3/main/SI3_manual.bib	(revision 11043)
@@ -1,11 +1,11 @@
 
-@Article{         assur_1958,
-  author        = {Assur, A},
-  year          = {1958},
-  month         = {01},
-  pages         = {106-138},
-  title         = {Composition of sea ice and its tensile strength},
-  volume        = {598},
-  journal       = {Arctic Sea Ice}
+@Article{	  assur_1958,
+  author	= {Assur, A},
+  year		= {1958},
+  month		= {01},
+  pages		= {106-138},
+  title		= {Composition of sea ice and its tensile strength},
+  volume	= {598},
+  journal	= {Arctic Sea Ice}
 }
 
@@ -234,5 +234,5 @@
 }
 
-@Article{	  h_yland_2002,
+@Article{	  hoyland_2002,
   author	= {Høyland, Knut V.},
   title		= {Consolidation of first-year sea ice ridges},
@@ -308,5 +308,5 @@
 }
 
-@Article{	  lepp_ranta_1995,
+@Article{	  lepparanta_1995,
   author	= {Leppäranta, Matti and Lensu, Mikko and Kosloff, Pekka and
 		  Veitch, Brian},
@@ -324,11 +324,11 @@
 }
 
-@Article{	  lepp_ranta_2011,
+@Article{	  lepparanta_2011,
   author	= {Leppäranta, Matti},
   title		= {Drift ice material},
   year		= 2011,
   pages		= {11–63},
-  doi		= {10.1007/978-3-642-04683-4_2},
-  url		= {http://dx.doi.org/10.1007/978-3-642-04683-4_2},
+  doi		= {10.1007/978-3-642-04683-4\_2},
+  url		= {http://dx.doi.org/10.1007/978-3-642-04683-4\_2},
   isbn		= 9783642046834,
   journal	= {The Drift of Sea Ice},
@@ -367,4 +367,13 @@
 }
 
+@Article{	  massonnet_2018,
+  author	= {Massonnet, F. and Barth\'el\'emy, A. and Worou, K. and
+		  Fichefet, T. and Vancoppenolle, M. and Rousset, C.},
+  title		= {Insights on the discretization of the ice thickness
+		  distribution in large-scale sea ice models},
+  journal	= {submitted},
+  year		= {2018}
+}
+
 @Article{	  maykut_1971,
   author	= {Maykut, Gary A. and Untersteiner, Norbert},
@@ -383,4 +392,14 @@
 }
 
+@Article{	  maykut_1973,
+  author	= {Maykut, G. A. and Thorndike, A. S.},
+  title		= {An approach to coupling the dynamics and thermodynamics of
+		  Arctic sea ice},
+  journal	= {AIDJEX Bulletin},
+  year		= {1973},
+  volume	= {21},
+  pages		= {23--29}
+}
+
 @Article{	  maykut_1986,
   author	= {Maykut, Gary A.},
@@ -388,6 +407,6 @@
   year		= 1986,
   pages		= {395–463},
-  doi		= {10.1007/978-1-4899-5352-0_6},
-  url		= {http://dx.doi.org/10.1007/978-1-4899-5352-0_6},
+  doi		= {10.1007/978-1-4899-5352-0\_6},
+  url		= {http://dx.doi.org/10.1007/978-1-4899-5352-0\_6},
   isbn		= 9781489953520,
   journal	= {The Geophysics of Sea Ice},
@@ -546,4 +565,14 @@
 }
 
+@Book{		  teos-10_2010,
+  title		= {{The international thermodynamic equation of seawater -
+		  2010: Calculation and use of thermodynamic properties}},
+  publisher	= {UNESCO (English)},
+  year		= {2010},
+  author	= {{IOC, SCOR and IAPSO}},
+  series	= {Intergovernmental Oceanographic Commission, Manuals and
+		  Guides No. 56}
+}
+
 @Article{	  thorndike_1975,
   author	= {Thorndike, A. S. and Rothrock, D. A. and Maykut, G. A. and
@@ -610,6 +639,6 @@
   year		= 1992,
   pages		= {113–138},
-  doi		= {10.1007/978-94-011-2809-4_20},
-  url		= {http://dx.doi.org/10.1007/978-94-011-2809-4_20},
+  doi		= {10.1007/978-94-011-2809-4\_20},
+  url		= {http://dx.doi.org/10.1007/978-94-011-2809-4\_20},
   isbn		= 9789401128094,
   journal	= {Interactive Dynamics of Convection and Solidification},
Index: NEMO/trunk/doc/latex/SI3/main/SI3_manual.tex
===================================================================
--- NEMO/trunk/doc/latex/SI3/main/SI3_manual.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/SI3/main/SI3_manual.tex	(revision 11043)
@@ -12,4 +12,9 @@
 %% Custom style (.sty)
 \usepackage{../../NEMO/main/NEMO_manual}
+\hypersetup{
+  pdftitle={SI³ – Sea Ice modelling Integrated Initiative – The NEMO Sea Ice engine},
+  pdfauthor={NEMO Sea Ice Working Group},
+  colorlinks
+}
 
 %% Include references and index for single subfile compilation
Index: NEMO/trunk/doc/latex/SI3/subfiles/chap_model_basics.tex
===================================================================
--- NEMO/trunk/doc/latex/SI3/subfiles/chap_model_basics.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/SI3/subfiles/chap_model_basics.tex	(revision 11043)
@@ -28,5 +28,5 @@
 
 \subsection{Scales, thermodynamics and dynamics}
-Because sea ice is much wider -- $\mathcal{O}$(100-1000 km) -- than thick -- $\mathcal{O}$(1 m) -- ice drift can be considered as purely horizontal: vertical motions around the hydrostatic equilibrium position are negligible. The same scaling argument justifies the assumption that heat exchanges are purely vertical\footnote{The latter assumption is probably less valid, because the horizontal scales of temperature variations are $\mathcal{O}$(10-100 m)}. It is on this basis that thermodynamics and dynamics are separated and rely upon different frameworks and sets of hypotheses: thermodynamics use the ice thickness distribution \citep{thorndike_1975} and the mushy-layer \citep{worster_1992} frameworks, whereas dynamics assume continuum mechanics \citep[e.g.,][]{lepp_ranta_2011}. Thermodynamics and dynamics interact by two means: first, advection impacts state variables; second, the horizontal momentum equation depends, among other things, on the ice state.
+Because sea ice is much wider -- $\mathcal{O}$(100-1000 km) -- than thick -- $\mathcal{O}$(1 m) -- ice drift can be considered as purely horizontal: vertical motions around the hydrostatic equilibrium position are negligible. The same scaling argument justifies the assumption that heat exchanges are purely vertical\footnote{The latter assumption is probably less valid, because the horizontal scales of temperature variations are $\mathcal{O}$(10-100 m)}. It is on this basis that thermodynamics and dynamics are separated and rely upon different frameworks and sets of hypotheses: thermodynamics use the ice thickness distribution \citep{thorndike_1975} and the mushy-layer \citep{worster_1992} frameworks, whereas dynamics assume continuum mechanics \citep[e.g.,][]{lepparanta_2011}. Thermodynamics and dynamics interact by two means: first, advection impacts state variables; second, the horizontal momentum equation depends, among other things, on the ice state.
 
 \subsection{Subgrid scale variations}
@@ -70,5 +70,5 @@
  & Description & Value & Units & Ref \\ \hline
 $c_i$ (cpic) & Pure ice specific heat & 2067 & J/kg/K & ? \\
-$c_w$ (rcp) & Seawater specific heat & 3991 & J/kg/K & \cite{TEOS_2010} \\
+$c_w$ (rcp) & Seawater specific heat & 3991 & J/kg/K & \cite{teos-10_2010} \\
 $L$ (lfus) & Latent heat of fusion (0$^\circ$C) & 334000 & J/kg/K & \cite{bitz_1999} \\
 $\rho_i$ (rhoic) & Sea ice density & 917 & kg/m$^3$ & \cite{bitz_1999} \\
@@ -154,5 +154,5 @@
 \subsection{Dynamic formulation}
 
-The formulation of ice dynamics is based on the continuum approach. The latter holds provided the drift ice particles are much larger than single ice floes, and much smaller than typical gradient scales. This compromise is rarely achieved in practice \citep{lepp_ranta_2011}. Yet the continuum approach generates a convenient momentum equation for the horizontal ice velocity vector $\mathbf{u}=(u,v)$, which can be solved with classical numerical methods (here, finite differences on the NEMO C-grid). The most important term in the momentum equation is internal stress. We follow the viscous-plastic (VP) rheological framework \citep{hibler_1979}, assuming that sea ice has no tensile strength but responds to compressive and shear deformations in a plastic way. In practice, the elastic-viscous-plastic (EVP) technique of  \citep{bouillon_2013} is used, more convient numerically than VP.  It is well accepted that the VP rheology and its relatives are the minimum complexity to get reasonable ice drift patterns \citep{kreyscher_2000}, but fail at generating the observed deformation patterns \citep{girard_2009}. This is a long-lasting problem: what is the ideal rheological model for sea ice and how it should be applied are still being debated \citep[see, e.g.][]{weiss_2013}. 
+The formulation of ice dynamics is based on the continuum approach. The latter holds provided the drift ice particles are much larger than single ice floes, and much smaller than typical gradient scales. This compromise is rarely achieved in practice \citep{lepparanta_2011}. Yet the continuum approach generates a convenient momentum equation for the horizontal ice velocity vector $\mathbf{u}=(u,v)$, which can be solved with classical numerical methods (here, finite differences on the NEMO C-grid). The most important term in the momentum equation is internal stress. We follow the viscous-plastic (VP) rheological framework \citep{hibler_1979}, assuming that sea ice has no tensile strength but responds to compressive and shear deformations in a plastic way. In practice, the elastic-viscous-plastic (EVP) technique of  \citep{bouillon_2013} is used, more convient numerically than VP.  It is well accepted that the VP rheology and its relatives are the minimum complexity to get reasonable ice drift patterns \citep{kreyscher_2000}, but fail at generating the observed deformation patterns \citep{girard_2009}. This is a long-lasting problem: what is the ideal rheological model for sea ice and how it should be applied are still being debated \citep[see, e.g.][]{weiss_2013}. 
 
 %-------------------------------------------------------------------------------------------------------------------------
@@ -296,5 +296,5 @@
 $C$ (rn\_crhg) & ice strength concentration param. & 20 & - & \citep{hibler_1979} \\
 $H^*$ (rn\_hstar) & maximum ridged ice thickness param. & 25 & m & \citep{lipscomb_2007} \\
-$p$ (rn\_por\_rdg) & porosity of new ridges & 0.3 & - & \citep{lepp_ranta_1995} \\
+$p$ (rn\_por\_rdg) & porosity of new ridges & 0.3 & - & \citep{lepparanta_1995} \\
 $amax$ (rn\_amax) & maximum ice concentration & 0.999 & - & -\\
 $h_0$ (rn\_hnewice) & thickness of newly formed ice & 0.1 & m & - \\
@@ -313,5 +313,5 @@
 Transport connects the horizontal velocity fields and the rest of the ice properties. LIM assumes that the ice properties in the different thickness categories are transported at the same velocity. The scheme of \cite{prather_1986}, based on the conservation of 0, 1$^{st}$ and 2$^{nd}$ order moments in $x-$ and $y-$directions,  is used, with some numerical diffusion if desired. Whereas this scheme is accurate, nearly conservative, it is also quite expensive since, for each advected field, five moments need to be advected, which proves CPU consuming, in particular when multiple categories are used. Other solutions are currently explored.
 
-The dissipation of energy associated with plastic failure under convergence and shear is accomplished by rafting (overriding of two ice plates) and ridging (breaking of an ice plate and subsequent piling of the broken ice blocks into pressure ridges). Thin ice preferentially rafts whereas thick ice preferentially ridges \citep{tuhkuri_2002}. Because observations of these processes are limited, their representation in LIM is rather heuristic. The amount of ice that rafts/ridges depends on the strain rate tensor invariants (shear and divergence) as in \citep{flato_1995}, while the ice categories involved are determined by a participation function favouring thin ice \citep{lipscomb_2007}. The thickness of ice being deformed ($h'$) determines whether ice rafts ($h'<$ 0.75 m) or ridges ($h'>$ 0.75 m), following \cite{haapala_2000}. The deformed ice thickness is $2h'$ after rafting, and is distributed between $2h'$ and $2 \sqrt{H^*h'}$ after ridging, where $H^* = 25$ m \citep{lipscomb_2007}. Newly ridged ice is highly porous, effectively trapping seawater. To represent this, a prescribed volume fraction (30\%) of newly ridged ice \citep{lepp_ranta_1995} incorporates mass, salt and heat are extracted from the ocean. Hence, in contrast with other models, the net thermodynamic ice production during convergence is not zero in LIM, since mass is added to sea ice during ridging. Consequently, simulated new ridges have high temperature and salinity as observed \citep{h_yland_2002}. A fraction of snow (50 \%) falls into the ocean during deformation.
+The dissipation of energy associated with plastic failure under convergence and shear is accomplished by rafting (overriding of two ice plates) and ridging (breaking of an ice plate and subsequent piling of the broken ice blocks into pressure ridges). Thin ice preferentially rafts whereas thick ice preferentially ridges \citep{tuhkuri_2002}. Because observations of these processes are limited, their representation in LIM is rather heuristic. The amount of ice that rafts/ridges depends on the strain rate tensor invariants (shear and divergence) as in \citep{flato_1995}, while the ice categories involved are determined by a participation function favouring thin ice \citep{lipscomb_2007}. The thickness of ice being deformed ($h'$) determines whether ice rafts ($h'<$ 0.75 m) or ridges ($h'>$ 0.75 m), following \cite{haapala_2000}. The deformed ice thickness is $2h'$ after rafting, and is distributed between $2h'$ and $2 \sqrt{H^*h'}$ after ridging, where $H^* = 25$ m \citep{lipscomb_2007}. Newly ridged ice is highly porous, effectively trapping seawater. To represent this, a prescribed volume fraction (30\%) of newly ridged ice \citep{lepparanta_1995} incorporates mass, salt and heat are extracted from the ocean. Hence, in contrast with other models, the net thermodynamic ice production during convergence is not zero in LIM, since mass is added to sea ice during ridging. Consequently, simulated new ridges have high temperature and salinity as observed \citep{hoyland_2002}. A fraction of snow (50 \%) falls into the ocean during deformation.
 
 \section{Ice thermodynamics}
Index: NEMO/trunk/doc/latex/SI3/subfiles/chap_ridging_rafting.tex
===================================================================
--- NEMO/trunk/doc/latex/SI3/subfiles/chap_ridging_rafting.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/SI3/subfiles/chap_ridging_rafting.tex	(revision 11043)
@@ -70,5 +70,5 @@
 \textbf{Rafting} is the piling of two ice sheets on top of each other. Rafting doubles the participating ice thickness and is a volume-conserving process. \cite{babko_2002} concluded that rafting plays a significant role during initial ice growth in fall, therefore we included it into the model. 
 
-\textbf{Ridging} is the piling of a series of broken ice blocks into pressure ridges. Ridging redistributes participating ice on a various range of thicknesses. Ridging does not conserve ice volume, as pressure ridges are porous. Therefore, the volume of ridged ice is larger than the volume of new ice being ridged. In the model, newly ridged is has a prescribed porosity $p=30\%$ (\textit{ridge\_por} in \textit{namelist\_ice}), following observations \citep{lepp_ranta_1995,h_yland_2002}. The importance of ridging is now since the early works of \citep{thorndike_1975}.
+\textbf{Ridging} is the piling of a series of broken ice blocks into pressure ridges. Ridging redistributes participating ice on a various range of thicknesses. Ridging does not conserve ice volume, as pressure ridges are porous. Therefore, the volume of ridged ice is larger than the volume of new ice being ridged. In the model, newly ridged is has a prescribed porosity $p=30\%$ (\textit{ridge\_por} in \textit{namelist\_ice}), following observations \citep{lepparanta_1995,hoyland_2002}. The importance of ridging is now since the early works of \citep{thorndike_1975}.
 
 The deformation modes are formulated using \textbf{participation} and \textbf{transfer} functions with specific contributions from ridging and rafting:
@@ -115,5 +115,5 @@
 \label{eq:nri}
 \end{equation}
-The redistributor $\gamma(h',h)$ specifies how area of thickness $h'$ is redistributed on area of thickness $h$. We follow \citep{hibler_1980} who constructed a rule, based on observations, that forces all ice participating in ridging with thickness $h'$ to be linearly distributed between ice that is between $2h'$ and $2\sqrt{H^*h'}$ thick, where $H^\star=100$ m (\textit{Hstar} in \textit{namelist\_ice}). This in turn determines how to construct the ice volume redistribution function $\Psi^v$. Volumes equal to participating area times thickness are removed from thin ice. They are redistributed following Hibler's rule. The factor $(1+p)$ accounts for initial ridge porosity $p$ (\textit{ridge\_por} in \textit{namelist\_ice}, defined as the fractional volume of seawater initially included into ridges. In many previous models, the initial ridge porosity has been assumed to be 0, which is not the case in reality since newly formed ridges are porous, as indicated by in-situ observations \citep{lepp_ranta_1995,h_yland_2002}. In other words, LIM3 creates a higher volume of ridged ice with the same participating ice.
+The redistributor $\gamma(h',h)$ specifies how area of thickness $h'$ is redistributed on area of thickness $h$. We follow \citep{hibler_1980} who constructed a rule, based on observations, that forces all ice participating in ridging with thickness $h'$ to be linearly distributed between ice that is between $2h'$ and $2\sqrt{H^*h'}$ thick, where $H^\star=100$ m (\textit{Hstar} in \textit{namelist\_ice}). This in turn determines how to construct the ice volume redistribution function $\Psi^v$. Volumes equal to participating area times thickness are removed from thin ice. They are redistributed following Hibler's rule. The factor $(1+p)$ accounts for initial ridge porosity $p$ (\textit{ridge\_por} in \textit{namelist\_ice}, defined as the fractional volume of seawater initially included into ridges. In many previous models, the initial ridge porosity has been assumed to be 0, which is not the case in reality since newly formed ridges are porous, as indicated by in-situ observations \citep{lepparanta_1995,hoyland_2002}. In other words, LIM3 creates a higher volume of ridged ice with the same participating ice.
 
 For the numerical computation of the integrals, we have to compute several temporary values:
@@ -152,5 +152,5 @@
 \section{Mechanical redistribution for other global ice variables}
 
-The other global ice state variables redistribution functions $\Psi^X$ are computed based on $\Psi^g$ for the ice age content and on $\Psi^{v^i}$ for the remainder (ice enthalpy and salt content, snow volume and enthalpy). The general principles behind this derivation are described in Appendix A of \cite{bitz_2001}. A fraction $f_s=0.5$ (\textit{fsnowrdg} and \textit{fsnowrft} in \textit{namelist\_ice}) of the snow volume and enthalpy is assumed to be lost during ridging and rafting and transferred to the ocean. The contribution of the seawater trapped into the porous ridges is included in the computation of the redistribution of ice enthalpy and salt content (i.e., $\Psi^{e^i}$ and $\Psi^{M^s}$). During this computation, seawater is supposed to be in thermal equilibrium with the surrounding ice blocks. Ridged ice desalination induces an implicit decrease in internal brine volume, and heat supply to the ocean, which accounts for ridge consolidation as described by \cite{h_yland_2002}. The inclusion of seawater in ridges does not imply any net change in ocean salinity. The energy used to cool down the seawater trapped in porous ridges until the seawater freezing point is rejected into the ocean.
+The other global ice state variables redistribution functions $\Psi^X$ are computed based on $\Psi^g$ for the ice age content and on $\Psi^{v^i}$ for the remainder (ice enthalpy and salt content, snow volume and enthalpy). The general principles behind this derivation are described in Appendix A of \cite{bitz_2001}. A fraction $f_s=0.5$ (\textit{fsnowrdg} and \textit{fsnowrft} in \textit{namelist\_ice}) of the snow volume and enthalpy is assumed to be lost during ridging and rafting and transferred to the ocean. The contribution of the seawater trapped into the porous ridges is included in the computation of the redistribution of ice enthalpy and salt content (i.e., $\Psi^{e^i}$ and $\Psi^{M^s}$). During this computation, seawater is supposed to be in thermal equilibrium with the surrounding ice blocks. Ridged ice desalination induces an implicit decrease in internal brine volume, and heat supply to the ocean, which accounts for ridge consolidation as described by \cite{hoyland_2002}. The inclusion of seawater in ridges does not imply any net change in ocean salinity. The energy used to cool down the seawater trapped in porous ridges until the seawater freezing point is rejected into the ocean.
 
 \end{document}
Index: NEMO/trunk/doc/latex/SI3/subfiles/introduction.tex
===================================================================
--- NEMO/trunk/doc/latex/SI3/subfiles/introduction.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/SI3/subfiles/introduction.tex	(revision 11043)
@@ -15,5 +15,5 @@
 
 % Limitations & scope
-%There are limitations to the applicability of models such as SI$^3$. The continuum approach is not invalid for grid cell size above at least 1 km, below which sea ice particles may include just a few floes, which is not sufficient \citep{lepp_ranta_2011}. Second, one must remember that our current knowledge of sea ice is not as complete as for the ocean: there are no fundamental equations such as Navier Stokes equations for sea ice. Besides, important features and processes span widely different scales, such as brine inclusions (1 $\mu$m-1 mm) \citep{perovich_1996}, horizontal thickness variations (1 m-100 km) \citep{percival_2008}, deformation and fracturing (10 m-1000 km) \citep{marsan_2004}. These impose complicated and often subjective subgrid-scale treatments. All in all, there is more empirism in sea ice models than in ocean models. 
+%There are limitations to the applicability of models such as SI$^3$. The continuum approach is not invalid for grid cell size above at least 1 km, below which sea ice particles may include just a few floes, which is not sufficient \citep{lepparanta_2011}. Second, one must remember that our current knowledge of sea ice is not as complete as for the ocean: there are no fundamental equations such as Navier Stokes equations for sea ice. Besides, important features and processes span widely different scales, such as brine inclusions (1 $\mu$m-1 mm) \citep{perovich_1996}, horizontal thickness variations (1 m-100 km) \citep{percival_2008}, deformation and fracturing (10 m-1000 km) \citep{marsan_2004}. These impose complicated and often subjective subgrid-scale treatments. All in all, there is more empirism in sea ice models than in ocean models. 
 
 In order to handle all the subsequent required subjective choices, we applied the following guidelines or principles:
Index: NEMO/trunk/doc/latex/TOP/main/TOP_manual.bib
===================================================================
--- NEMO/trunk/doc/latex/TOP/main/TOP_manual.bib	(revision 11041)
+++ NEMO/trunk/doc/latex/TOP/main/TOP_manual.bib	(revision 11043)
@@ -87,5 +87,5 @@
 		  CFC-113, CCl4, SF6 and N2O (NCEI Accession 0164584)},
   year		= 2017,
-  doi		= {10.3334/cdiac/otg.cfc_atm_hist_2015},
+  doi		= {10.3334/cdiac/otg.cfc\_atm\_hist\_2015},
   url		= {https://accession.nodc.noaa.gov/0164584},
   publisher	= {NOAA National Centers for Environmental Information}
@@ -349,6 +349,6 @@
   number	= {3–4},
   issn		= {0033-8222},
-  doi		= {10.2458/azu_js_rc.55.16402},
-  url		= {http://dx.doi.org/10.2458/azu_js_rc.55.16402},
+  doi		= {10.2458/azu\_js\_rc.55.16402},
+  url		= {http://dx.doi.org/10.2458/azu\_js\_rc.55.16402},
   journal	= {Radiocarbon},
   publisher	= {Cambridge University Press (CUP)}
@@ -422,5 +422,5 @@
 		  Béranger, K. and Schneider, A. and Beuvier, J. and Somot,
 		  S.},
-  title		= {Simulated anthropogenic CO&lt;sub&gt;2&lt;/sub&gt; storage
+  title		= {Simulated anthropogenic CO$_{2}$ storage
 		  and acidification of the Mediterranean Sea},
   year		= 2015,
@@ -448,6 +448,6 @@
   pages		= {1869–1887},
   issn		= {1945-5755},
-  doi		= {10.2458/azu_js_rc.55.16947},
-  url		= {http://dx.doi.org/10.2458/azu_js_rc.55.16947},
+  doi		= {10.2458/azu\_js\_rc.55.16947},
+  url		= {http://dx.doi.org/10.2458/azu\_js\_rc.55.16947},
   journal	= {Radiocarbon},
   publisher	= {Cambridge University Press (CUP)}
@@ -495,5 +495,5 @@
 }
 
-@Article{	  toggweiler_1989,
+@Article{	  toggweiler_1989a,
   author	= {Toggweiler, J. R. and Dixon, K. and Bryan, K.},
   title		= {Simulations of radiocarbon in a coarse-resolution world
@@ -510,5 +510,5 @@
 }
 
-@Article{	  toggweiler_1989,
+@Article{	  toggweiler_1989b,
   author	= {Toggweiler, J. R. and Dixon, K. and Bryan, K.},
   title		= {Simulations of radiocarbon in a coarse-resolution world
@@ -615,5 +615,5 @@
   doi		= {10.1016/j.tree.2012.10.021},
   url		= {http://dx.doi.org/10.1016/j.tree.2012.10.021},
-  journal	= {Trends in Ecology & Evolution},
-  publisher	= {Elsevier BV}
-}
+  journal	= {Trends in Ecology \& Evolution},
+  publisher	= {Elsevier BV}
+}
Index: NEMO/trunk/doc/latex/TOP/main/TOP_manual.tex
===================================================================
--- NEMO/trunk/doc/latex/TOP/main/TOP_manual.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/TOP/main/TOP_manual.tex	(revision 11043)
@@ -12,4 +12,9 @@
 %% Custom style (.sty)
 \usepackage{../../NEMO/main/NEMO_manual}
+\hypersetup{
+  pdftitle={TOP – Tracers in Ocean Paradigm – The NEMO Tracers engine},
+  pdfauthor={NEMO TOP Working Group},
+  colorlinks
+}
 
 %% Include references and index for single subfile compilation
Index: NEMO/trunk/doc/latex/TOP/subfiles/model_description.tex
===================================================================
--- NEMO/trunk/doc/latex/TOP/subfiles/model_description.tex	(revision 11041)
+++ NEMO/trunk/doc/latex/TOP/subfiles/model_description.tex	(revision 11043)
@@ -26,5 +26,5 @@
 \end{equation}
 
-where expressions of $D^{lC}$ and $D^{vC}$ depend on the choice for the lateral and vertical subgrid scale parameterizations, see equations 5.10 and 5.11 in \citep{Madec_Bk2008} 
+where expressions of $D^{lC}$ and $D^{vC}$ depend on the choice for the lateral and vertical subgrid scale parameterizations, see equations 5.10 and 5.11 in \citep{nemo_manual} 
 
 {S(C)} , the first term on the right hand side of \ref{Eq_tracer}; is the SMS - Source Minus Sink - inherent to the tracer.  In the case of biological tracer such as phytoplankton, {S(C)} is the balance between phytoplankton growth and its decay through mortality and grazing. In the case of a tracer comprising carbon,  {S(C)} accounts for gas exchange, river discharge, flux to the sediments, gravitational sinking and other biological processes. In the case of a radioactive tracer, {S(C)} is simply loss due to radioactive decay.
@@ -61,5 +61,5 @@
         \item \textbf{AGE}	    :    Water age tracking
         \item \textbf{MY\_TRC}  :	Template for creation of new modules and external BGC models coupling
-        \item \textbf{PISCES}    :	 Built in BGC model. See \citep{Aumont_al_2015} for a throughout description.
+        \item \textbf{PISCES}    :	 Built in BGC model. See \citep{aumont_2015} for a throughout description.
 \end{itemize}
 %  ----------------------------------------------------------
@@ -73,5 +73,5 @@
 \nlst{namtrc_adv}
 %-------------------------------------------------------------------------------------------------------------
-The advection schemes used for the passive tracers are the same than the ones for $T$ and $S$ and described in section 5.1 of \citep{Madec_Bk2008}. The choice of an advection scheme  can be selected independently and  can differ from the ones used for active tracers. This choice is made in the \textit{namtrc\_adv} namelist, by  setting to \textit{true} one and only one of the logicals \textit{ln\_trcadv\_xxx}, the same way of what is done for dynamics.
+The advection schemes used for the passive tracers are the same than the ones for $T$ and $S$ and described in section 5.1 of \citep{nemo_manual}. The choice of an advection scheme  can be selected independently and  can differ from the ones used for active tracers. This choice is made in the \textit{namtrc\_adv} namelist, by  setting to \textit{true} one and only one of the logicals \textit{ln\_trcadv\_xxx}, the same way of what is done for dynamics.
 cen2, MUSCL2, and UBS are not \textit{positive} schemes meaning that 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
 
@@ -80,5 +80,5 @@
 \nlst{namtrc_ldf}
 %-------------------------------------------------------------------------------------------------------------
-In NEMO v4.0, the passive tracer diffusion has necessarily the same form as the active tracer diffusion, meaning that the numerical scheme must be the same. However the passive tracer mixing coefficient can be chosen as a multiple of the active ones by changing the value of \textit{rn\_ldf\_multi} in namelist \textit{namtrc\_ldf}. The choice of numerical scheme is then set  in the \ngn{namtra\_ldf} namelist for the dynamic described in section 5.2 of \citep{Madec_Bk2008}. 
+In NEMO v4.0, the passive tracer diffusion has necessarily the same form as the active tracer diffusion, meaning that the numerical scheme must be the same. However the passive tracer mixing coefficient can be chosen as a multiple of the active ones by changing the value of \textit{rn\_ldf\_multi} in namelist \textit{namtrc\_ldf}. The choice of numerical scheme is then set  in the \ngn{namtra\_ldf} namelist for the dynamic described in section 5.2 of \citep{nemo_manual}. 
 
 
@@ -145,5 +145,5 @@
 
 
- This implementation was first used in the CORE-II intercomparison runs described e.g.\ in \citet{Danabasoglu_al_2014}.
+ This implementation was first used in the CORE-II intercomparison runs described e.g.\ in \citet{danabasoglu_2014}.
 
 \subsection{Inert carbons tracer}
@@ -178,5 +178,5 @@
 Measuring the dissolved concentrations of the gases -- as well as the mixing ratios between them -- shows circulation pathways within the ocean as well as water mass ages (i.e. the time since last contact with the
 atmosphere). This feature of the gases has made them valuable across a wide range of oceanographic problems. One use lies in ocean modelling, where they can be used to evaluate the realism of the circulation and 
-ventilation of models, key for understanding the behaviour of wider modelled marine biogeochemistry (e.g. \citep{Dutay_al_2002,Palmieri_2015}). \\
+ventilation of models, key for understanding the behaviour of wider modelled marine biogeochemistry (e.g. \citep{dutay_2002,palmieri_2015}). \\
 
 Modelling these gases (henceforth CFCs) in NEMO is done within the passive tracer transport module, TOP, using the conservation state equation \ref{Eq_tracer}
@@ -187,5 +187,5 @@
 stable within the ocean, we assume that there are no sinks (i.e. no loss processes) within the ocean interior. 
 Consequently, the sinks-minus-sources term for CFCs consists only of their air-sea fluxes, $F_{cfc}$, as
-described in the Ocean Model Inter-comparison Project (OMIP) protocol \citep{Orr_al_2017}:
+described in the Ocean Model Inter-comparison Project (OMIP) protocol \citep{orr_2017}:
 
 % Because CFCs being stable in the ocean, we consider that there is no CFCs sink.
@@ -213,5 +213,5 @@
 Where $Sol$ is the gas solubility in mol~m$^{-3}$~pptv$^{-1}$, as defined in Equation \ref{equ_Sol_CFC}; 
 and $P_{cfc}$ is the atmosphere concentration of the CFC (in parts per trillion by volume, pptv).
-This latter concentration is provided to the model by the historical time-series of \citet{Bullister_2015}.
+This latter concentration is provided to the model by the historical time-series of \citet{bullister_2017}.
 This includes bulk atmospheric concentrations of the CFCs for both hemispheres -- this is necessary because of 
 the geographical asymmetry in the production and release of CFCs to the atmosphere. 
@@ -220,5 +220,5 @@
 
 The piston velocity $K_{w}$ is a function of 10~m wind speed (in m~s$^{-1}$) and sea surface temperature, 
-$T$ (in $^{\circ}$C), and is calculated here following \citet{Wanninkhof_1992}:
+$T$ (in $^{\circ}$C), and is calculated here following \citet{wanninkhof_1992}:
 
 \begin{eqnarray}
@@ -229,7 +229,7 @@
 Where $X_{conv}$ = $\frac{0.01}{3600}$, a conversion factor that changes the piston velocity 
 from cm~h$^{-1}$ to m~s$^{-1}$; 
-$a$ is a constant re-estimated by \citet{Wanninkhof_2014} to 0.251 (in $\frac{cm~h^{-1}}{(m~s^{-1})^{2}}$);
+$a$ is a constant re-estimated by \citet{wanninkhof_2014} to 0.251 (in $\frac{cm~h^{-1}}{(m~s^{-1})^{2}}$);
 and $u$ is the 10~m wind speed in m~s$^{-1}$ from either an atmosphere model or reanalysis atmospheric forcing.
-$Sc$ is the Schmidt number, and is calculated as follow, using coefficients from \citet{Wanninkhof_2014} (see Table \ref{tab_Sc}).
+$Sc$ is the Schmidt number, and is calculated as follow, using coefficients from \citet{wanninkhof_2014} (see Table \ref{tab_Sc}).
 
 \begin{eqnarray}
@@ -240,5 +240,5 @@
 The solubility, $Sol$, used in Equation \ref{equ_C_sat} is calculated in mol~l$^{-1}$~atm$^{-1}$, 
 and is specific for each gas. 
-It has been experimentally estimated by \citet{Warner_Weiss_1985} as a function of temperature 
+It has been experimentally estimated by \citet{warner_1985} as a function of temperature 
 and salinity:
 
@@ -363,10 +363,10 @@
 where $\Rq_{\textrm{ref}}$ is a reference ratio. For the purpose of ocean ventilation studies $\Rq_{\textrm{ref}}$ is set to one.
 
-Here we adopt the approach of \cite{Fiadeiro_1982} and \cite{Toggweiler_al_1989a,Toggweiler_al_1989b} in which  the ratio $\Rq$ is transported rather than the individual concentrations C and $\cq$.
-This approach calls for a strong assumption, i.e., that of a homogeneous and constant dissolved inorganic carbon (DIC) field \citep{Toggweiler_al_1989a,Mouchet_2013}. While in terms of
-oceanic $\Dcq$, it yields similar results to approaches involving carbonate chemistry, it underestimates the bomb radiocarbon inventory because it assumes a constant air-sea $\cd$ disequilibrium (Mouchet, 2013). Yet, field reconstructions of the ocean bomb $\cq$ inventory are also biased low \citep{Naegler_2009} since they assume that the anthropogenic perturbation did not affect ocean DIC since the pre-bomb epoch. For these reasons, bomb $\cq$ inventories obtained with the present method are directly comparable to reconstructions based on field measurements.
-
-This simplified approach also neglects the effects of fractionation (e.g.,  air-sea exchange) and of biological processes. Previous studies by \cite{Bacastow_MaierReimer_1990} and \cite{Joos_al_1997} resulted in nearly identical $\Dcq$ distributions among experiments considering biology or not.
-Since observed $\Rq$ ratios are corrected for the isotopic fractionation when converted to the standard $\Dcq$ notation \citep{Stuiver_Polach_1977} the model results are directly comparable to observations. 
+Here we adopt the approach of \cite{fiadeiro_1982} and \cite{toggweiler_1989a,toggweiler_1989b} in which  the ratio $\Rq$ is transported rather than the individual concentrations C and $\cq$.
+This approach calls for a strong assumption, i.e., that of a homogeneous and constant dissolved inorganic carbon (DIC) field \citep{toggweiler_1989a,mouchet_2013}. While in terms of
+oceanic $\Dcq$, it yields similar results to approaches involving carbonate chemistry, it underestimates the bomb radiocarbon inventory because it assumes a constant air-sea $\cd$ disequilibrium (Mouchet, 2013). Yet, field reconstructions of the ocean bomb $\cq$ inventory are also biased low \citep{naegler_2009} since they assume that the anthropogenic perturbation did not affect ocean DIC since the pre-bomb epoch. For these reasons, bomb $\cq$ inventories obtained with the present method are directly comparable to reconstructions based on field measurements.
+
+This simplified approach also neglects the effects of fractionation (e.g.,  air-sea exchange) and of biological processes. Previous studies by \cite{bacastow_1990} and \cite{joos_1997} resulted in nearly identical $\Dcq$ distributions among experiments considering biology or not.
+Since observed $\Rq$ ratios are corrected for the isotopic fractionation when converted to the standard $\Dcq$ notation \citep{stuiver_1977} the model results are directly comparable to observations. 
 
 Therefore the simplified approach is justified for the purpose of assessing the circulation and ventilation of OGCMs.
@@ -378,5 +378,5 @@
 where $\lambda$ is the radiocarbon decay rate, ${\mathbf{u}}$ the 3-D velocity field, and $\mathbf{K}$ the diffusivity tensor.
 
-At the air-sea interface a Robin boundary condition \citep{Haine_2006} is applied to \eqref{eq:quick}, i.e., the flux
+At the air-sea interface a Robin boundary condition \citep{haine_2006} is applied to \eqref{eq:quick}, i.e., the flux
 through the interface is proportional to the difference in the ratios between
 the ocean and the atmosphere
@@ -391,5 +391,5 @@
 
 
-The $\cd$ transfer velocity is based on the empirical formulation of \cite{Wanninkhof_1992} with chemical enhancement \citep{Wanninkhof_Knox_1996,Wanninkhof_2014}. The original formulation is modified to account for the reduction of the  air-sea exchange rate in the presence of sea ice. Hence
+The $\cd$ transfer velocity is based on the empirical formulation of \cite{wanninkhof_1992} with chemical enhancement \citep{wanninkhof_1996,wanninkhof_2014}. The original formulation is modified to account for the reduction of the  air-sea exchange rate in the presence of sea ice. Hence
 \begin{equation}
 \kappa_\cd=\left( K_W\,\mathrm{w}^2 + b  \right)\, (1-f_\mathrm{ice})\,\sqrt{660/Sc}, \label{eq:wanc14}
@@ -397,5 +397,5 @@
 with $\mathrm{w}$ the wind magnitude, $f_\mathrm{ice}$ the fractional ice cover, and $Sc$ the Schmidt number.
 $K_W$ in \eqref{eq:wanc14} is an empirical coefficient with dimension of an inverse velocity.
-The chemical enhancement term $b$ is represented as a function of temperature $T$ \citep{Wanninkhof_1992}
+The chemical enhancement term $b$ is represented as a function of temperature $T$ \citep{wanninkhof_1992}
 \begin{equation}
 b=2.5 ( 0.5246 + 0.016256 T+ 0.00049946  * T^2 ). \label{eq:wanchem}
@@ -413,14 +413,14 @@
 \label{sec:param}
  %
-The radiocarbon decay rate (\CODE{rlam14}; in \texttt{trcnam\_c14} module) is set to $\lambda=(1/8267)$ yr$^{-1}$ \citep{Stuiver_Polach_1977}, which corresponds to a half-life of 5730 yr.\\[1pt]
-%
-The Schmidt number $Sc$, Eq. \eqref{eq:wanc14}, is calculated with the help of the formulation of \cite{Wanninkhof_2014}. The $\cd$ solubility $K_0$ in \eqref{eq:Rspeed} is taken from \cite{Weiss_1974}. $K_0$ and $Sc$ are computed with the OGCM temperature and salinity fields (\texttt{trcsms\_c14} module).\\[1pt]
+The radiocarbon decay rate (\CODE{rlam14}; in \texttt{trcnam\_c14} module) is set to $\lambda=(1/8267)$ yr$^{-1}$ \citep{stuiver_1977}, which corresponds to a half-life of 5730 yr.\\[1pt]
+%
+The Schmidt number $Sc$, Eq. \eqref{eq:wanc14}, is calculated with the help of the formulation of \cite{wanninkhof_2014}. The $\cd$ solubility $K_0$ in \eqref{eq:Rspeed} is taken from \cite{weiss_1974}. $K_0$ and $Sc$ are computed with the OGCM temperature and salinity fields (\texttt{trcsms\_c14} module).\\[1pt]
 %
 The following parameters intervening in the air-sea exchange rate are set in \texttt{namelist\_c14}:
 \begin{itemize}
-\item The reference DIC concentration $\overline{\Ct}$ (\CODE{xdicsur}) intervening in \eqref{eq:Rspeed} is classically set to 2 mol m$^{-3}$ \citep{Toggweiler_al_1989a,Orr_al_2001,Butzin_al_2005}.
-%
-\item The value of the empirical coefficient $K_W$ (\CODE{xkwind}) in \eqref{eq:wanc14} depends on the wind field and on the model upper ocean mixing rate \citep{Toggweiler_al_1989a,Wanninkhof_1992,Naegler_2009,Wanninkhof_2014}.
-It should be adjusted so that the globally averaged $\cd$ piston velocity is $\kappa_\cd = 16.5\pm 3.2$ cm/h \citep{Naegler_2009}.
+\item The reference DIC concentration $\overline{\Ct}$ (\CODE{xdicsur}) intervening in \eqref{eq:Rspeed} is classically set to 2 mol m$^{-3}$ \citep{toggweiler_1989a,orr_2001,butzin_2005}.
+%
+\item The value of the empirical coefficient $K_W$ (\CODE{xkwind}) in \eqref{eq:wanc14} depends on the wind field and on the model upper ocean mixing rate \citep{toggweiler_1989a,wanninkhof_1992,naegler_2009,wanninkhof_2014}.
+It should be adjusted so that the globally averaged $\cd$ piston velocity is $\kappa_\cd = 16.5\pm 3.2$ cm/h \citep{naegler_2009}.
 %The sensitivity to this parametrization is discussed in section \ref{sec:result}.
 %
@@ -440,7 +440,7 @@
 \CODE{kc14typ}=0
 
-Unless otherwise specified in \texttt{namelist\_c14}, the atmospheric $\Rq_a$ (\CODE{rc14at}) is set to one, the atmospheric $\cd$ (\CODE{pco2at}) to 280 ppm, and the ocean $\Rq$ is initialized with \CODE{rc14init=0.85}, i.e., $\Dcq=$-150\textperthousand  \cite[typical for deep-ocean, Fig 6 in][]{Key_al_2004}.
-
-Equilibrium experiment should last until 98\% of the ocean volume exhibit a drift of less than 0.001\textperthousand/year \citep{Orr_al_2000}; this is usually achieved after few kyr (Fig. \ref{fig:drift}).
+Unless otherwise specified in \texttt{namelist\_c14}, the atmospheric $\Rq_a$ (\CODE{rc14at}) is set to one, the atmospheric $\cd$ (\CODE{pco2at}) to 280 ppm, and the ocean $\Rq$ is initialized with \CODE{rc14init=0.85}, i.e., $\Dcq=$-150\textperthousand  \cite[typical for deep-ocean, Fig 6 in][]{key_2004}.
+
+Equilibrium experiment should last until 98\% of the ocean volume exhibit a drift of less than 0.001\textperthousand/year \citep{orr_2000}; this is usually achieved after few kyr (Fig. \ref{fig:drift}).
 %
 \begin{figure}[!h]
@@ -469,6 +469,6 @@
 
 The model  is integrated from a given initial date following the observed records provided from 1765 AD on ( Fig. \ref{fig:bomb}).
-The file \texttt{atmc14.dat}  \cite[][\& I. Levin, personal comm.]{Enting_al_1994} provides atmospheric $\Dcq$ for three latitudinal bands: 90S-20S,    20S-20N \&    20N-90N.
-Atmospheric $\cd$ in the file \texttt{splco2.dat} is obtained from a spline fit through ice core data and direct atmospheric measurements \cite[][\& J. Orr, personal comm.]{Orr_al_2000}.
+The file \texttt{atmc14.dat}  \cite[][\& I. Levin, personal comm.]{enting_1994} provides atmospheric $\Dcq$ for three latitudinal bands: 90S-20S,    20S-20N \&    20N-90N.
+Atmospheric $\cd$ in the file \texttt{splco2.dat} is obtained from a spline fit through ice core data and direct atmospheric measurements \cite[][\& J. Orr, personal comm.]{orr_2000}.
 Dates in these forcing files are expressed as yr AD.
 
@@ -496,6 +496,6 @@
 Atmospheric $\Rq_a$ and $\cd$ are prescribed from forcing files. The ocean $\Rq$ is initialized with the value attributed to \CODE{rc14init} in \texttt{namelist\_c14}.
 
-The file \texttt{intcal13.14c} \citep{Reimer_al_2013} contains atmospheric $\Dcq$ from 0 to 50 kyr cal BP\footnote{cal BP: number of years before 1950 AD}.
-The $\cd$ forcing is provided in file \texttt{ByrdEdcCO2.txt}. The content of this file is based on  the high resolution record from EPICA Dome C \citep{Monnin_al_2004} for the Holocene and the Transition, and on Byrd Ice Core CO2 Data for 20--90 kyr BP  \citep{Ahn_Brook_2008}. These atmospheric values are reproduced in Fig. \ref{fig:paleo}. Dates in these files are expressed as yr BP. 
+The file \texttt{intcal13.14c} \citep{reimer_2013} contains atmospheric $\Dcq$ from 0 to 50 kyr cal BP\footnote{cal BP: number of years before 1950 AD}.
+The $\cd$ forcing is provided in file \texttt{ByrdEdcCO2.txt}. The content of this file is based on  the high resolution record from EPICA Dome C \citep{monnin_2004} for the Holocene and the Transition, and on Byrd Ice Core CO2 Data for 20--90 kyr BP  \citep{ahn_2008}. These atmospheric values are reproduced in Fig. \ref{fig:paleo}. Dates in these files are expressed as yr BP. 
 
 To ensure that the atmospheric forcing is applied properly as well as that output files contain consistent dates and inventories the experiment should be set up carefully.
@@ -539,5 +539,5 @@
 The radiocarbon age is computed as  $(-1/\lambda) \ln{ \left( \Rq \right)}$, with zero age corresponding to $\Rq=1$. 
 
-The reservoir age is the age difference between the ocean uppermost layer and the atmosphere. It is usually reported as conventional radiocarbon age; i.e., computed by means of the Libby radiocarbon mean life \cite[8033 yr;][]{Stuiver_Polach_1977}
+The reservoir age is the age difference between the ocean uppermost layer and the atmosphere. It is usually reported as conventional radiocarbon age; i.e., computed by means of the Libby radiocarbon mean life \cite[8033 yr;][]{stuiver_1977}
 \begin{align}
 {^{14}\tau_\mathrm{c}}= -8033 \; \ln \left(1 + \frac{\Dcq}{10^3}\right), \label{eq:convage}
@@ -549,5 +549,5 @@
 N_A \Rq_\mathrm{oxa} \overline{\Ct} \left( \int_\Omega \Rq d\Omega \right) /10^{26}, \label{eq:inv}
 \end{equation}
-where $N_A$ is the Avogadro's number ($N_A=6.022\times10^{23}$ at/mol), $\Rq_\mathrm{oxa}$ is the oxalic acid radiocarbon standard \cite[$\Rq_\mathrm{oxa}=1.176\times10^{-12}$;][]{Stuiver_Polach_1977}, and $\Omega$ is the ocean volume.  Bomb $\cq$ inventories are traditionally reported in units of $10^{26}$ atoms, hence the denominator in \eqref{eq:inv}.
+where $N_A$ is the Avogadro's number ($N_A=6.022\times10^{23}$ at/mol), $\Rq_\mathrm{oxa}$ is the oxalic acid radiocarbon standard \cite[$\Rq_\mathrm{oxa}=1.176\times10^{-12}$;][]{stuiver_1977}, and $\Omega$ is the ocean volume.  Bomb $\cq$ inventories are traditionally reported in units of $10^{26}$ atoms, hence the denominator in \eqref{eq:inv}.
 
 All transformations from second to year, and inversely, are performed with the help of the physical constant \CODE{rsiyea} the sideral year length expressed in seconds\footnote{The variable (\CODE{nyear\_len}) which reports the length in days of the previous/current/future year (see \textrm{oce\_trc.F90}) is not a constant. }.
@@ -564,9 +564,9 @@
 Two versions of PISCES are available in NEMO v4.0 :
 
-PISCES-v2, by setting in namelist\_pisces\_ref  \np{ln\_p4z} to true,  can be seen as one of the many Monod models \citep{Monod_1942}. It assumes a constant Redfield ratio and phytoplankton growth depends on the external concentration in nutrients. There are twenty-four prognostic variables (tracers) including two phytoplankton compartments  (diatoms and nanophytoplankton), two zooplankton size-classes (microzooplankton and  mesozooplankton) and a description of the carbonate chemistry. Formulations in PISCES-v2 are based on a mixed Monod/Quota formalism: On one hand, stoichiometry of C/N/P is fixed and growth rate of phytoplankton is limited by the external availability in N, P and Si. On the other hand, the iron and silicium quotas are variable and growth rate of phytoplankton is limited by the internal availability in Fe. Various parameterizations can be activated in PISCES-v2, setting for instance the complexity of iron chemistry or the description of particulate organic materials. 
-
-PISCES-QUOTA has been built on the PISCES-v2 model described in \citet{Aumont_al_2015}. PISCES-QUOTA has thirty-nine prognostic compartments. Phytoplankton growth can be controlled by five modeled limiting nutrients: Nitrate and Ammonium, Phosphate, Silicate and Iron. Five living compartments are represented: Three phytoplankton size classes/groups corresponding to picophytoplankton, nanophytoplankton and diatoms, and two zooplankton size classes which are microzooplankton and mesozooplankton. For phytoplankton, the prognostic variables are the carbon, nitrogen, phosphorus,  iron, chlorophyll and silicon biomasses (the latter only for diatoms). This means that the N/C, P/C, Fe/C and Chl/C ratios of both phytoplankton groups as well as the Si/C ratio of diatoms are prognostically predicted  by the model. Zooplankton are assumed to be strictly homeostatic \citep[e.g.,][]{Sterner_2002,Woods_Wilson_2013,Meunier_al_2014}. As a consequence, the C/N/P/Fe ratios of these groups are maintained constant and are not allowed to vary. In PISCES, the Redfield ratios C/N/P are set to 122/16/1 \citep{Takahashi_al_1985} and the -O/C ratio is set to 1.34 \citep{Kortzinger_al_2001}. No silicified zooplankton is assumed. The bacterial pool is not yet explicitly modeled.
-
-There are three non-living compartments: Semi-labile dissolved organic matter, small sinking particles, and large sinking particles. As a consequence of the variable stoichiometric ratios of phytoplankton and of the stoichiometric regulation of zooplankton, elemental ratios in organic matter cannot be supposed constant anymore as that was the case in PISCES-v2. Indeed, the nitrogen, phosphorus, iron, silicon and calcite pools of the particles are now all explicitly modeled. The sinking speed of the particles is not altered by their content in calcite and biogenic silicate (''The ballast effect'', \citep{Honjo_1996,Armstrong_al_2002}). The latter particles are assumed to sink at the same speed as the large organic matter particles. All the non-living compartments experience aggregation due to turbulence and differential settling as well as Brownian coagulation for DOM.
+PISCES-v2, by setting in namelist\_pisces\_ref  \np{ln\_p4z} to true,  can be seen as one of the many Monod models \citep{monod_1958}. It assumes a constant Redfield ratio and phytoplankton growth depends on the external concentration in nutrients. There are twenty-four prognostic variables (tracers) including two phytoplankton compartments  (diatoms and nanophytoplankton), two zooplankton size-classes (microzooplankton and  mesozooplankton) and a description of the carbonate chemistry. Formulations in PISCES-v2 are based on a mixed Monod/Quota formalism: On one hand, stoichiometry of C/N/P is fixed and growth rate of phytoplankton is limited by the external availability in N, P and Si. On the other hand, the iron and silicium quotas are variable and growth rate of phytoplankton is limited by the internal availability in Fe. Various parameterizations can be activated in PISCES-v2, setting for instance the complexity of iron chemistry or the description of particulate organic materials. 
+
+PISCES-QUOTA has been built on the PISCES-v2 model described in \citet{aumont_2015}. PISCES-QUOTA has thirty-nine prognostic compartments. Phytoplankton growth can be controlled by five modeled limiting nutrients: Nitrate and Ammonium, Phosphate, Silicate and Iron. Five living compartments are represented: Three phytoplankton size classes/groups corresponding to picophytoplankton, nanophytoplankton and diatoms, and two zooplankton size classes which are microzooplankton and mesozooplankton. For phytoplankton, the prognostic variables are the carbon, nitrogen, phosphorus,  iron, chlorophyll and silicon biomasses (the latter only for diatoms). This means that the N/C, P/C, Fe/C and Chl/C ratios of both phytoplankton groups as well as the Si/C ratio of diatoms are prognostically predicted  by the model. Zooplankton are assumed to be strictly homeostatic \citep[e.g.,][]{sterner_2003,woods_2013,meunier_2014}. As a consequence, the C/N/P/Fe ratios of these groups are maintained constant and are not allowed to vary. In PISCES, the Redfield ratios C/N/P are set to 122/16/1 \citep{takahashi_1985} and the -O/C ratio is set to 1.34 \citep{kortzinger_2001}. No silicified zooplankton is assumed. The bacterial pool is not yet explicitly modeled.
+
+There are three non-living compartments: Semi-labile dissolved organic matter, small sinking particles, and large sinking particles. As a consequence of the variable stoichiometric ratios of phytoplankton and of the stoichiometric regulation of zooplankton, elemental ratios in organic matter cannot be supposed constant anymore as that was the case in PISCES-v2. Indeed, the nitrogen, phosphorus, iron, silicon and calcite pools of the particles are now all explicitly modeled. The sinking speed of the particles is not altered by their content in calcite and biogenic silicate (''The ballast effect'', \citep{honjo_1996,armstrong_2001}). The latter particles are assumed to sink at the same speed as the large organic matter particles. All the non-living compartments experience aggregation due to turbulence and differential settling as well as Brownian coagulation for DOM.
 
  
Index: NEMO/trunk/doc/manual_build.sh
===================================================================
--- NEMO/trunk/doc/manual_build.sh	(revision 11041)
+++ NEMO/trunk/doc/manual_build.sh	(revision 11043)
@@ -45,4 +45,5 @@
 	printf "\t  The export should be available at root\n"
 	printf "\t  If not check LaTeX log in ./latex/$model/main/${model}_manual.log\n"
+	echo
 done
 
