Index: /trunk/DOC/TexFiles/Chapters/Annex_A.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Annex_A.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Annex_A.tex	(revision 1223)
@@ -8,10 +8,11 @@
 
 
-In order to establish the set of Primitive Equation in curvilinear $s$-coordinates ($i.e.$ 
-orthogonal curvilinear coordinate in the horizontal and $s$-coordinate in the vertical), we 
-start from the set of equation established in \S\ref{PE_zco_Eq} for the special case 
-$k = z$ and thus $e_3 = 1$, and we introduce an arbitrary vertical coordinate 
-$s = s(i,j,z,t)$. Let us define a new vertical scale factor by $e_3 = \partial z / \partial s$ 
-(which now depends on $(i,j,z,t)$) and the horizontal slope of $s$-surfaces by :
+In order to establish the set of Primitive Equation in curvilinear $s$-coordinates 
+($i.e.$ an orthogonal curvilinear coordinate in the horizontal and $s$-coordinate 
+in the vertical), we start from the set of equations established in \S\ref{PE_zco_Eq} 
+for the special case $k = z$ and thus $e_3 = 1$, and we introduce an arbitrary 
+vertical coordinate $s = s(i,j,z,t)$. Let us define a new vertical scale factor by 
+$e_3 = \partial z / \partial s$ (which now depends on $(i,j,z,t)$) and the horizontal 
+slope of $s$-surfaces by :
 \begin{equation} \label{Apdx_A_s_slope}
 \sigma _1 =\frac{1}{e_1 }\;\left. {\frac{\partial z}{\partial i}} \right|_s 
@@ -20,6 +21,6 @@
 \end{equation}
 
-The chain rule to establish the model equations in the curvilinear $s$-coordinate system 
-is:
+The chain rule to establish the model equations in the curvilinear $s$-coordinate 
+system is:
 \begin{equation} \label{Apdx_A_s_chain_rule}
 \begin{aligned}
@@ -41,5 +42,6 @@
 \end{equation}
 
-In particular applying the time derivative chain rule to $z$ provide the expression of $w_s$,  the vertical velocity of the $s-$surfaces:
+In particular applying the time derivative chain rule to $z$ provides the 
+expression for $w_s$,  the vertical velocity of the $s-$surfaces:
 \begin{equation} \label{Apdx_A_w_in_s}
 w_s   =  \left.   \frac{\partial z }{\partial t}   \right|_s 
@@ -54,5 +56,7 @@
 \label{Apdx_B_continuity}
 
-Using (\ref{Apdx_A_s_chain_rule}) and the fact that the horizontal scale factors $e_1$ and $e_2$ do not depend on the vertical coordinate, the divergence of the velocity relative to the ($i$,$j$,$z$) coordinate system is transformed as follows:
+Using (\ref{Apdx_A_s_chain_rule}) and the fact that the horizontal scale factors 
+$e_1$ and $e_2$ do not depend on the vertical coordinate, the divergence of 
+the velocity relative to the ($i$,$j$,$z$) coordinate system is transformed as follows:
 
 \begin{align*}
@@ -109,9 +113,12 @@
  \end{align*} 
  
-Here, $w$ is the vertical velocity relative to the $z-$coordinate system. Introducing the dia-surface velocity component, $\omega $, defined as the velocity relative to the moving $s$-surfaces and normal to them:
+Here, $w$ is the vertical velocity relative to the $z-$coordinate system. 
+Introducing the dia-surface velocity component, $\omega $, defined as 
+the velocity relative to the moving $s$-surfaces and normal to them:
 \begin{equation} \label{Apdx_A_w_s}
 \omega  = w - w_s - \sigma _1 \,u - \sigma _2 \,v    \\
 \end{equation}
-with $w_s$ given by \eqref{Apdx_A_w_in_s}, we obtain the expression of the divergence of the velocity in the curvilinear $s$-coordinate system:
+with $w_s$ given by \eqref{Apdx_A_w_in_s}, we obtain the expression for 
+the divergence of the velocity in the curvilinear $s$-coordinate system:
 \begin{align*} \label{Apdx_A_A4}
 \nabla \cdot {\rm {\bf U}} 
@@ -142,5 +149,6 @@
 \end{align*}
 
-As a result, the continuity equation \eqref{Eq_PE_continuity} in $s$-coordinates becomes:
+As a result, the continuity equation \eqref{Eq_PE_continuity} in the 
+$s$-coordinates becomes:
 \begin{equation} \label{Apdx_A_A5}
 \frac{1}{e_3 } \frac{\partial e_3}{\partial t} 
@@ -157,7 +165,7 @@
 \label{Apdx_B_momentum}
 
-Let us consider \eqref{Eq_PE_dyn_vect}, the first component of the 
-momentum equation in the vector invariant form (similar manipulations can be performed on the second one). Its non linear term can be transformed 
-as follows:
+Let us consider \eqref{Eq_PE_dyn_vect}, the first component of the momentum 
+equation in the vector invariant form (similar manipulations can be performed 
+on the second component). Its non-linear term can be transformed as follows:
 
 \begin{align*}
@@ -205,6 +213,7 @@
 \end{align*}
 
-Therefore, the non-linear terms of the momentum equation have the same form 
-in $z-$ and $s-$coordinates but with the addition of the time derivative of the velocity: 
+Therefore, the non-linear terms of the momentum equation have the same 
+form in $z-$ and $s-$coordinates but with the addition of the time derivative 
+of the velocity: 
 \begin{multline}  \label{Apdx_A_momentum_NL}
 +\left. \zeta \right|_z v-\frac{1}{2e_1 }\left. {\frac{\partial (u^2+v^2)}{\partial i}} \right|_z 
@@ -224,8 +233,11 @@
 \end{equation}
 
-An additional term appears in (\ref{Apdx_A_grad_p}) which accounts for the tilt of model 
-levels.
-
-Introducing \eqref{Apdx_A_momentum_NL} and \eqref{Apdx_A_grad_p} in \eqref{Eq_PE_dyn_vect} and regrouping the time derivative terms in the left hand side, and performing the same manipulation on the second component, we obtain the vector invariant form of momentum equation in $s-$coordinate :
+An additional term appears in (\ref{Apdx_A_grad_p}) which accounts for the 
+tilt of model levels.
+
+Introducing \eqref{Apdx_A_momentum_NL} and \eqref{Apdx_A_grad_p} in \eqref{Eq_PE_dyn_vect} and regrouping the time derivative terms in the left 
+hand side, and performing the same manipulation on the second component, 
+we obtain the vector invariant form of the momentum equations in the 
+$s-$coordinate :
 \begin{subequations} \label{Apdx_A_dyn_vect}
 \begin{multline} \label{Apdx_A_PE_dyn_vect_u}
@@ -249,5 +261,8 @@
 \end{subequations}
 
-It has the same form as in $z-$coordinate but the vertical scale factor that has appeared inside the time derivative. The form of the vertical physics and forcing terms remain unchanged. The form of the lateral physics is discussed in appendix~\ref{Apdx_B}.  
+It has the same form as in the $z-$coordinate but for the vertical scale factor 
+that has appeared inside the time derivative. The form of the vertical physics 
+and forcing terms remains unchanged. The form of the lateral physics is 
+discussed in appendix~\ref{Apdx_B}.  
 
 % ================================================================
@@ -257,6 +272,6 @@
 \label{Apdx_B_tracer}
 
-The tracer equation is obtained using the same calculation as for the 
-continuity equation and then regrouping the time derivative terms in the left hand side :
+The tracer equation is obtained using the same calculation as for the continuity 
+equation and then regrouping the time derivative terms in the left hand side :
 
 \begin{multline} \label{Apdx_A_tracer}
@@ -269,5 +284,5 @@
 
 
-The expression of the advection term is a straight consequence of (A.4), the 
-expression of the 3D divergence in $s$-coordinates established above. 
-
+The expression for the advection term is a straight consequence of (A.4), the 
+expression of the 3D divergence in the $s$-coordinates established above. 
+
Index: /trunk/DOC/TexFiles/Chapters/Annex_B.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Annex_B.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Annex_B.tex	(revision 1223)
@@ -13,5 +13,6 @@
 
 
-In $z$-coordinate, the horizontal/vertical second order tracer diffusive operator is given by:
+In the $z$-coordinate, the horizontal/vertical second order tracer diffusion operator 
+is given by:
 \begin{multline} \label{Apdx_B1}
  D^T = \frac{1}{e_1 \, e_2}      \left[ 
@@ -22,5 +23,7 @@
 \end{multline}
 
-In $s$-coordinate, we defined the slopes of $s$-surfaces, $\sigma_1$ and $\sigma_2$ by (A.1) and the vertical/horizontal ratio of diffusive coefficient by $\epsilon = A^{vT} / A^{lT}$. The diffusive operator is given by:
+In the $s$-coordinate, we defined the slopes of $s$-surfaces, $\sigma_1$ and 
+$\sigma_2$ by (!!!A.1!!!) and the vertical/horizontal ratio of diffusion coefficient 
+by $\epsilon = A^{vT} / A^{lT}$. The diffusion operator is given by:
 
 \begin{equation} \label{Apdx_B2}
@@ -34,5 +37,5 @@
 \end{array} }} \right)
 \end{equation}
-or in expended form:
+or in expanded form:
 \begin{multline} \label{Apdx_B3}
 D^T=\frac{1}{e_1\,e_2\,e_3 }\;\left[ {\quad \; \; e_2\,e_3\,A^{lT} \;\left. 
@@ -45,6 +48,12 @@
 \end{multline}
 
-
-Equation \eqref{Apdx_B2} (or equivalently \eqref{Apdx_B3}) is obtained from \eqref{Apdx_B1} without any additional assumption. Indeed, for the special case $k=z$ and thus $e_3 =1$, we introduce an arbitrary vertical coordinate $s = s (i,j,z)$ as in Appendix~\ref{Apdx_A} and use \eqref{Apdx_A_s_slope} and \eqref{Apdx_A_s_chain_rule}. Since no cross horizontal derivate $\partial _i \partial _j $ appears in \eqref{Apdx_B1}, the ($i$,$z$) and ($j$,$z$) planes are independent. The demonstration can then be done for the ($i$,$z$)~$\to$~($j$,$s$) transformation without any loss of generality:
+Equation \eqref{Apdx_B2} (or equivalently \eqref{Apdx_B3}) is obtained from 
+\eqref{Apdx_B1} without any additional assumption. Indeed, for the special 
+case $k=z$ and thus $e_3 =1$, we introduce an arbitrary vertical coordinate 
+$s = s (i,j,z)$ as in Appendix~\ref{Apdx_A} and use \eqref{Apdx_A_s_slope} 
+and \eqref{Apdx_A_s_chain_rule}. Since no cross horizontal derivative 
+$\partial _i \partial _j $ appears in \eqref{Apdx_B1}, the ($i$,$z$) and ($j$,$z$) 
+planes are independent. The derivation can then be demonstrated for the 
+($i$,$z$)~$\to$~($j$,$s$) transformation without any loss of generality:
 
 \begin{equation*}
@@ -87,5 +96,7 @@
 \end{multline*}
 
-Since the horizontal scale factor do not depend on the vertical coordinate, the last term of the first line and the first term of the last line cancel, while the second line reduces to a single vertical derivative, so it becomes:
+Since the horizontal scale factors do not depend on the vertical coordinate, 
+the last term of the first line and the first term of the last line cancel, while 
+the second line reduces to a single vertical derivative, so it becomes:
 
 \begin{multline*}
@@ -94,5 +105,5 @@
 \end{multline*}
 
-in other words, the horizontal Laplacian operator in the ($i$,$s$)plane takes the following expression :
+in other words, the horizontal Laplacian operator in the ($i$,$s$) plane takes the following form :
 
 \begin{equation*}
@@ -123,5 +134,7 @@
 
 
-The iso/diapycnal diffusive tensor $\textbf {A}_{\textbf I}$ expressed in the  ($i$,$j$,$k$) curvilinear coordinate system in which the equations of the ocean circulation model are formulated, takes the following expression \citep{Redi_JPO82}:
+The iso/diapycnal diffusive tensor $\textbf {A}_{\textbf I}$ expressed in the ($i$,$j$,$k$) 
+curvilinear coordinate system in which the equations of the ocean circulation model are 
+formulated, takes the following form \citep{Redi_JPO82}:
 
 \begin{equation*}
@@ -141,5 +154,6 @@
 \end{equation*}
 
-In practice, the isopycnal slopes are generally less than $10^{-2}$ in the ocean, so $\textbf {A}_{\textbf I}$ can be simplified appreciably \citep{Cox1987}:
+In practice, the isopycnal slopes are generally less than $10^{-2}$ in the ocean, so 
+$\textbf {A}_{\textbf I}$ can be simplified appreciably \citep{Cox1987}:
 \begin{equation*}
 {\textbf{A}_{\textbf{I}}} \approx A^{lT}
@@ -151,5 +165,7 @@
 \end{equation*}
 
-The resulting isopycnal operator conserves the quantity and dissipates its square. The demonstration of the first property is trivial as \eqref{Apdx_B2} is the divergence of fluxes. Let us demonstrate the second one:
+The resulting isopycnal operator conserves the quantity and dissipates its square. 
+The demonstration of the first property is trivial as \eqref{Apdx_B2} is the divergence 
+of fluxes. Let us demonstrate the second one:
 \begin{equation*}
 \iiint\limits_D T\;\nabla .\left( {\textbf{A}}_{\textbf{I}} \nabla T \right)\,dv = -\iiint\limits_D \nabla T\;.\left( {\textbf{A}}_{\textbf{I}} \nabla T \right)\,dv
@@ -169,6 +185,5 @@
 the property becomes obvious. 
 
-The resulting diffusive operator in $z$-coordinates has the following 
-expression :
+The resulting diffusion operator in $z$-coordinate has the following form :
 \begin{multline*} \label{Apdx_B_ldfiso}
  D^T=\frac{1}{e_1 e_2 }\left\{ {\;\frac{\partial }{\partial i}\left[ {A_h \left( {\frac{e_2 }{e_1 }\frac{\partial T}{\partial i}-a_1 \frac{e_2}{e_3}\frac{\partial T}{\partial k}} \right)} \right]\;} \right.\;\; \\ 
@@ -177,5 +192,8 @@
 \end{multline*}
 
-It has to be emphasised that the simplification introduced leads to a decoupling between ($i$,$z$) and ($j$,$z$) planes. The operator has therefore the same expression as \eqref{Apdx_B3}, the diffusive operator obtained for geopotential diffusion in $s$-coordinate. 
+It has to be emphasised that the simplification introduced, leads to a decoupling 
+between ($i$,$z$) and ($j$,$z$) planes. The operator has therefore the same 
+expression as \eqref{Apdx_B3}, the diffusion operator obtained for geopotential 
+diffusion in the $s$-coordinate. 
 
 % ================================================================
@@ -185,5 +203,7 @@
 \label{Apdx_B_3}
 
-The second order momentum diffusive operator (Laplacian) in $z$-coordinate is found by applying \eqref{Eq_PE_lap_vector}, the expression of the Laplacian of a vector,  to the horizontal velocity vector : 
+The second order momentum diffusion operator (Laplacian) in the $z$-coordinate 
+is found by applying \eqref{Eq_PE_lap_vector}, the expression for the Laplacian 
+of a vector,  to the horizontal velocity vector : 
 \begin{align*}
 \Delta {\textbf{U}}_h 
@@ -221,13 +241,16 @@
 \end{array} }} \right)
 \end{align*}
-Using \eqref{Eq_PE_div}, the definition of the horizontal divergence, the third componant of the second vector is obviously zero and thus :
+Using \eqref{Eq_PE_div}, the definition of the horizontal divergence, the third 
+componant of the second vector is obviously zero and thus :
 \begin{equation*}
 \Delta {\textbf{U}}_h = \nabla _h \left( \chi \right) - \nabla _h \times \left( \zeta \right) + \frac {1}{e_3 } \frac {\partial }{\partial k} \left( {\frac {1}{e_3 } \frac{\partial {\textbf{ U}}_h }{\partial k}} \right)
 \end{equation*}
 
-Note that this operator ensures a full separation between the vorticity and horizontal divergence fields (see Appendix~\ref{Apdx_C}). It is only equal to a Laplacian applied on each component in Cartesian coordinate, not on the sphere.
+Note that this operator ensures a full separation between the vorticity and horizontal 
+divergence fields (see Appendix~\ref{Apdx_C}). It is only equal to a Laplacian 
+applied to each component in Cartesian coordinates, not on the sphere.
 
 The horizontal/vertical second order (Laplacian type) operator used to diffuse 
-horizontal momentum in $z$-coordinate takes therefore the following expression :
+horizontal momentum in the $z$-coordinate therefore takes the following form :
 \begin{equation} \label{Apdx_B_Lap_U}
  {\textbf{D}}^{\textbf{U}} = 
@@ -237,5 +260,5 @@
 				\frac{\partial {\rm {\bf U}}_h }{\partial k}} \right) \\ 
 \end{equation}
-that is in expanded form:
+that is, in expanded form:
 \begin{align*}
 D^{\textbf{U}}_u 
@@ -249,3 +272,8 @@
 \end{align*}
 
-Note Bene: introducing a rotation in \eqref{Apdx_B_Lap_U} does not lead to any useful expression for the iso/diapycnal Laplacian operator in $z$-coordinate. Similarly, we did not found an expression of practical use for the geopotential horizontal/vertical Laplacian operator in $s$-coordinate. Generally, \eqref{Apdx_B_Lap_U} is used in both $z$- and $s$-coordinate system, that is a Laplacian diffusion is applied on momentum along the coordinate directions.
+Note Bene: introducing a rotation in \eqref{Apdx_B_Lap_U} does not lead to a 
+useful expression for the iso/diapycnal Laplacian operator in the $z$-coordinate. 
+Similarly, we did not found an expression of practical use for the geopotential 
+horizontal/vertical Laplacian operator in the $s$-coordinate. Generally, 
+\eqref{Apdx_B_Lap_U} is used in both $z$- and $s$-coordinate systems, that is 
+a Laplacian diffusion is applied on momentum along the coordinate directions.
Index: /trunk/DOC/TexFiles/Chapters/Annex_C.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Annex_C.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Annex_C.tex	(revision 1223)
@@ -17,10 +17,16 @@
 \label{Apdx_C.1}
 
-
-First, the boundary condition on the vertical velocity (no flux through the surface and the bottom) is established for the discrete set of momentum equations. Then, it is shown that the non linear terms of the momentum equation are written such that the potential enstrophy of a horizontally non divergent flow is preserved while all the other non-diffusive terms preserve the kinetic energy: the energy is also preserved in practice. In addition, an option is also offer for the vorticity term discretization which provides 
-a total kinetic energy conserving discretization for that term. 
-
-Nota Bene: these properties are established here in the rigid-lid case and for the 2nd order centered scheme. A forthcoming update will be their generalisation to the free surface case
-and higher order scheme.
+First, the boundary condition on the vertical velocity (no flux through the surface 
+and the bottom) is established for the discrete set of momentum equations. 
+Then, it is shown that the non-linear terms of the momentum equation are written 
+such that the potential enstrophy of a horizontally non-divergent flow is preserved 
+while all the other non-diffusive terms preserve the kinetic energy; in practice the 
+energy is also preserved. In addition, an option is also offered for the vorticity term 
+discretization which provides a total kinetic energy conserving discretization for 
+that term. 
+
+Nota Bene: these properties are established here in the rigid-lid case and for the 
+2nd order centered scheme. A forthcoming update will be their generalisation to 
+the free surface case and higher order scheme.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -31,9 +37,9 @@
 
 
-The discrete set of momentum equations used in rigid lid approximation 
+The discrete set of momentum equations used in the rigid-lid approximation 
 automatically satisfies the surface and bottom boundary conditions 
 (no flux through the surface and the bottom: $w_{surface} =w_{bottom} =~0$). 
 Indeed, taking the discrete horizontal divergence of the vertical sum of the 
-horizontal momentum equations (Eqs. (II.2.1) and (II.2.2)~) wheighted by the 
+horizontal momentum equations (!!!Eqs. (II.2.1) and (II.2.2)!!!) weighted by the 
 vertical scale factors, it becomes:
 \begin{flalign*}
@@ -83,9 +89,15 @@
 
 
-The surface boundary condition associated with the rigid lid approximation ($w_{surface} =0)$ is imposed in the computation of the vertical velocity (II.2.5). Therefore, it turns out to be:
+The surface boundary condition associated with the rigid lid approximation 
+($w_{surface} =0)$ is imposed in the computation of the vertical velocity (!!! II.2.5!!!!). 
+Therefore, it turns out to be:
 \begin{equation*}
 \frac{\partial } {\partial t}w_{bottom} \equiv 0
 \end{equation*}
-As the bottom velocity is initially set to zero, it remains zero all the time. Symmetrically, if $w_{bottom} =0$ is used in the computation of the vertical velocity (upward integral of the horizontal divergence), the same computation leads to $w_{surface} =0$ as soon as the surface vertical velocity is initially set to zero.
+As the bottom velocity is initially set to zero, it remains zero all the time. 
+Symmetrically, if $w_{bottom} =0$ is used in the computation of the vertical 
+velocity (upward integral of the horizontal divergence), the same computation 
+leads to $w_{surface} =0$ as soon as the surface vertical velocity is initially 
+set to zero.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -101,5 +113,11 @@
 \label{Apdx_C_vor} 
 
-Potential vorticity is located at $f$-points and defined as: $\zeta / e_{3f}$. The standard discrete formulation of the relative vorticity term obviously conserves potential vorticity (ENS scheme). It also conserves the potential enstrophy for a horizontally non-divergent flow (i.e. $\chi $=0) but not the total kinetic energy. Indeed, using the symmetry or skew symmetry properties of the operators (Eqs \eqref{DOM_mi_adj} and \eqref{DOM_di_adj}), it can be shown that:
+Potential vorticity is located at $f$-points and defined as: $\zeta / e_{3f}$. 
+The standard discrete formulation of the relative vorticity term obviously 
+conserves potential vorticity (ENS scheme). It also conserves the potential 
+enstrophy for a horizontally non-divergent flow (i.e. $\chi $=0) but not the 
+total kinetic energy. Indeed, using the symmetry or skew symmetry properties 
+of the operators (Eqs \eqref{DOM_mi_adj} and \eqref{DOM_di_adj}), it can 
+be shown that:
 \begin{equation} \label{Apdx_C_1.1}
 \int_D {\zeta / e_3\,\;{\textbf{k}}\cdot \frac{1} {e_3} \nabla \times \left( {\zeta \;{\textbf{k}}\times {\textbf{U}}_h } \right)\;dv} \equiv 0
@@ -164,8 +182,8 @@
 \end{align*}
 
-Note that the demonstration is done here for the relative potential 
-vorticity but it still hold for the planetary ($f/e_3$) and the total 
-potential vorticity $((\zeta +f) /e_3 )$. Another formulation of 
-the two components of the vorticity term is optionally offered (ENE scheme) :
+Note that the derivation is demonstrated here for the relative potential 
+vorticity but it applies also to the planetary ($f/e_3$) and the total 
+potential vorticity $((\zeta +f) /e_3 )$. Another formulation of the two 
+components of the vorticity term is optionally offered (ENE scheme) :
 \begin{equation*}
 	- \zeta \;{\textbf{k}}\times {\textbf {U}}_h 
@@ -185,7 +203,8 @@
 \end{equation*}
 
-This formulation does not conserve the enstrophy but the total kinetic 
-energy. It is also possible to mix the two formulations in order to conserve 
-enstrophy on the relative vorticity term and energy on the Coriolis term.
+This formulation does not conserve the enstrophy but it does conserve the 
+total kinetic energy. It is also possible to mix the two formulations in order 
+to conserve enstrophy on the relative vorticity term and energy on the 
+Coriolis term.
 \begin{flalign*}
 &\int\limits_D - \textbf{U}_h \cdot   \left(  \zeta \;\textbf{k} \times \textbf{U}_h  \right)  \;  dv &&  \\
@@ -219,6 +238,10 @@
  = - \int_D \textbf{U}_h \cdot w \frac{\partial \textbf{U}_h} {\partial k}\;dv
 \end{equation*}
-Indeed, using successively \eqref{DOM_di_adj} ($i.e.$ the skew symmetry property of the $\delta$ operator) and the incompressibility, then again \eqref{DOM_di_adj}, then 
-the commutativity of operators $\overline {\,\cdot \,}$ and $\delta$, and finally \eqref{DOM_mi_adj} ($i.e.$ the  symmetry property of the $\overline {\,\cdot \,}$ operator) applied in the horizontal and vertical direction, it becomes:
+Indeed, using successively \eqref{DOM_di_adj} ($i.e.$ the skew symmetry 
+property of the $\delta$ operator) and the incompressibility, then 
+\eqref{DOM_di_adj} again, then the commutativity of operators 
+$\overline {\,\cdot \,}$ and $\delta$, and finally \eqref{DOM_mi_adj} 
+($i.e.$ the symmetry property of the $\overline {\,\cdot \,}$ operator) 
+applied in the horizontal and vertical directions, it becomes:
 \begin{flalign*}
 &\int_D \textbf{U}_h \cdot \nabla_h \left( \frac{1}{2}\;{\textbf{U}_h}^2 \right)\;dv	&&&\\
@@ -252,5 +275,5 @@
 	\biggl\{  \overline {e_{1T}\,e_{2T} \,w}^{\,i+1/2}\;2 	
 			\overline {u}^{\,k+1/2}\; \delta_{k+1/2}         \left[ u \right]     %&&&  \\
-	 + \overline {e_{1T}\,e_{2T} \,w}^{\,j+1/2}\;2 	\overline {v}^{\,k+1/2}\; \delta_{k+1/2} 			\left[ v \right]  \;
+	 + \overline {e_{1T}\,e_{2T} \,w}^{\,j+1/2}\;2 	\overline {v}^{\,k+1/2}\; \delta_{k+1/2} 	\left[ v \right]  \;
 	\biggr\} 
 	&&\displaybreak[0] \\ 
@@ -272,5 +295,9 @@
 \end{flalign*}
 
-The main point here is that the satisfaction of this property links the choice of the discrete formulation of vertical advection and of horizontal gradient of KE. Choosing one imposes the other. For example KE can also be discretized as $1/2\,({\overline u^{\,i}}^2 + {\overline v^{\,j}}^2)$. This leads to the following expression for the vertical advection:
+The main point here is that the satisfaction of this property links the choice of 
+the discrete formulation of the vertical advection and of the horizontal gradient 
+of KE. Choosing one imposes the other. For example KE can also be discretized 
+as $1/2\,({\overline u^{\,i}}^2 + {\overline v^{\,j}}^2)$. This leads to the following 
+expression for the vertical advection:
 \begin{equation*}
 \frac{1} {e_3 }\; w\; \frac{\partial \textbf{U}_h } {\partial k}
@@ -284,5 +311,7 @@
 \end{array}} } \right)
 \end{equation*}
-a formulation that requires a additional horizontal mean compare to the one used in NEMO. Nine velocity points have to be used instead of 3. This is the reason why it has not been choosen.
+a formulation that requires an additional horizontal mean in contrast with 
+the one used in NEMO. Nine velocity points have to be used instead of 3. 
+This is the reason why it has not been chosen.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -298,5 +327,7 @@
 \label{Apdx_C.1.3.1} 
 
-In flux from the vorticity term reduces to a Coriolis term in which the Coriolis parameter has been modified to account for the ``metric'' term. This altered Coriolis parameter is discretised at F-point. It is given by: 
+In flux from the vorticity term reduces to a Coriolis term in which the Coriolis 
+parameter has been modified to account for the ``metric'' term. This altered 
+Coriolis parameter is discretised at an f-point. It is given by: 
 \begin{equation*}
 f+\frac{1} {e_1 e_2 } 
@@ -308,5 +339,7 @@
 \end{equation*}
 
-The ENE scheme is then applied to obtain the vorticity term in flux form. It therefore conserves the total KE. The demonstration is same as for the vorticity term in vector invariant form (\S\ref{Apdx_C_vor}).
+The ENE scheme is then applied to obtain the vorticity term in flux form. 
+It therefore conserves the total KE. The derivation is the same as for the 
+vorticity term in the vector invariant form (\S\ref{Apdx_C_vor}).
 
 % -------------------------------------------------------------------------------------------------------------
@@ -316,5 +349,9 @@
 \label{Apdx_C.1.3.2} 
 
-The flux form operator of the momentum advection is evaluated using a centered second order finite difference scheme. Because of the flux form, the discrete operator does not contribute to the global budget of linear momentum. Because of the centered second order scheme, it conserves the horizontal kinetic energy, that is :
+The flux form operator of the momentum advection is evaluated using a 
+centered second order finite difference scheme. Because of the flux form, 
+the discrete operator does not contribute to the global budget of linear 
+momentum. Because of the centered second order scheme, it conserves 
+the horizontal kinetic energy, that is :
 
 \begin{equation} \label{Apdx_C_I.3.10}
@@ -326,5 +363,7 @@
 \end{equation}
 
-Let us demonstrate this property for the first term of the scalar product (i.e. considering just the the terms associated with the i-component of the advection):
+Let us demonstrate this property for the first term of the scalar product 
+($i.e.$ considering just the the terms associated with the i-component of 
+the advection):
 \begin{flalign*}
 &\int_D u \cdot \nabla \cdot \left(   \textbf{U}\,u   \right) \; dv    &&&\\
@@ -382,5 +421,8 @@
 \end{flalign*}
 
-When the UBS scheme is used to evaluate the flux form momentum advection, the discrete operator does not contribute to the global budget of linear momentum (flux form). The horizontal kinetic energy is not conserved, but forced to decrease (the scheme is diffusive). 
+When the UBS scheme is used to evaluate the flux form momentum advection, 
+the discrete operator does not contribute to the global budget of linear momentum 
+(flux form). The horizontal kinetic energy is not conserved, but forced to decay 
+($i.e.$ the scheme is diffusive). 
 
 % -------------------------------------------------------------------------------------------------------------
@@ -391,5 +433,11 @@
 
 
-A pressure gradient has no contribution to the evolution of the vorticity as the curl of a gradient is zero. In $z$-coordinate, this properties is satisfied locally on a C-grid with 2nd order finite differences (property \eqref{Eq_DOM_curl_grad}). When the equation of state is linear ($i.e.$ when an advective-diffusive equation for density can be derived from those of temperature and salinity) the change of KE due to the work of pressure forces is balanced by the change of potential energy due to buoyancy forces: 
+A pressure gradient has no contribution to the evolution of the vorticity as the 
+curl of a gradient is zero. In the $z$-coordinate, this property is satisfied locally 
+on a C-grid with 2nd order finite differences (property \eqref{Eq_DOM_curl_grad}). 
+When the equation of state is linear ($i.e.$ when an advection-diffusion equation 
+for density can be derived from those of temperature and salinity) the change of 
+KE due to the work of pressure forces is balanced by the change of potential 
+energy due to buoyancy forces: 
 \begin{equation*}
 \int_D - \frac{1} {\rho_o} \left. \nabla p^h \right|_z \cdot \textbf{U}_h \;dv 
@@ -397,5 +445,8 @@
 \end{equation*}
 
-This property can be satisfied in discrete sense for both $z$- and $s$-coordinates. Indeed, defining the depth of a $T$-point, $z_T$ defined as the sum of the vertical scale factors at $w$-points starting from the surface, the workof pressure forces can be written as:
+This property can be satisfied in a discrete sense for both $z$- and $s$-coordinates. 
+Indeed, defining the depth of a $T$-point, $z_T$, as the sum of the vertical scale 
+factors at $w$-points starting from the surface, the work of pressure forces can be 
+written as:
 \begin{flalign*}
 &\int_D - \frac{1} {\rho_o} \left. \nabla p^h \right|_z \cdot \textbf{U}_h \;dv   &&& \\
@@ -410,6 +461,7 @@
 \end{flalign*}
 
-Using  \eqref{DOM_di_adj}, $i.e.$ the skew symmetry property of the $\delta$ operator, \eqref{Eq_wzv}, the continuity equation), and \eqref{Eq_dynhpg_sco}, the hydrostatic 
-equation in $s$-coordinate, it turns out to be: 
+Using  \eqref{DOM_di_adj}, $i.e.$ the skew symmetry property of the $\delta$ 
+operator, \eqref{Eq_wzv}, the continuity equation), and \eqref{Eq_dynhpg_sco}, 
+the hydrostatic equation in the $s$-coordinate, it becomes: 
 \begin{flalign*} 
 \equiv& \frac{1} {\rho_o} \sum\limits_{i,j,k}    \biggl\{ 
@@ -467,8 +519,8 @@
 \end{flalign*}
 
-Note that this property strongly constraints the discrete expression of both 
-the depth of $T-$points and of the term added to the pressure gradient in 
-$s$-coordinate. Nevertheless, it is almost never satisfied as a linear equation of state 
-is rarely used.
+Note that this property strongly constrains the discrete expression of both 
+the depth of $T-$points and of the term added to the pressure gradient in the
+$s$-coordinate. Nevertheless, it is almost never satisfied since a linear equation 
+of state is rarely used.
 
 % -------------------------------------------------------------------------------------------------------------
@@ -479,8 +531,16 @@
 
 
-The surface pressure gradient has no contribution to the evolution of the vorticity. This property is trivially satisfied locally as the equation verified by $\psi $ has been derived from the discrete formulation of the momentum equation and of the curl. But it has to be noticed that since the elliptic equation verified by $\psi $ is solved numerically by an iterative solver (preconditioned conjugate gradient or successive over relaxation), the 
-property is only satisfied at the precision required on the solver used.
-
-With the rigid-lid approximation, the change of KE due to the work of surface pressure forces is exactly zero. This is satisfied in discrete form, at the precision required on the elliptic solver used to solve this equation. This can be demonstrated as follows:
+The surface pressure gradient has no contribution to the evolution of the vorticity. 
+This property is trivially satisfied locally since the equation verified by $\psi$ has 
+been derived from the discrete formulation of the momentum equation and of the curl. 
+But it has to be noted that since the elliptic equation satisfied by $\psi$ is solved 
+numerically by an iterative solver (preconditioned conjugate gradient or successive 
+over relaxation), the property is only satisfied at the precision requested for the 
+solver used.
+
+With the rigid-lid approximation, the change of KE due to the work of surface 
+pressure forces is exactly zero. This is satisfied in discrete form, at the precision 
+requested for the elliptic solver used to solve this equation. This can be 
+demonstrated as follows:
 \begin{flalign*}
 \int\limits_D  - \frac{1} {\rho_o} \nabla_h \left( p_s \right) \cdot \textbf{U}_h \;dv%	&&& \\
@@ -503,5 +563,5 @@
 	&& \\ 
 %
-\intertext{using the relation between \textit{$\psi $} and the vertically sum of the velocity, it becomes:}
+\intertext{using the relation between \textit{$\psi $} and the vertical sum of the velocity, it becomes:}
 %
 &\equiv \sum\limits_{i,j} 
@@ -537,5 +597,10 @@
 \end{flalign*}
 
-The last equality is obtained using \eqref{Eq_dynspg_rl}, the discrete barotropic streamfunction time evolution equation. By the way, this shows that \eqref{Eq_dynspg_rl} is the only way do compute the streamfunction, otherwise the surface pressure forces will work. Nevertheless, since the elliptic equation verified by $\psi $ is solved numerically by an iterative solver, the property is only satisfied at the precision required on the solver.
+The last equality is obtained using \eqref{Eq_dynspg_rl}, the discrete barotropic 
+streamfunction time evolution equation. By the way, this shows that 
+\eqref{Eq_dynspg_rl} is the only way to compute the streamfunction, 
+otherwise the surface pressure forces will do work. Nevertheless, since 
+the elliptic equation satisfied by $\psi $ is solved numerically by an iterative 
+solver, the property is only satisfied at the precision requested for the solver.
 
 % ================================================================
@@ -546,5 +611,13 @@
 
 
-All the numerical schemes used in NEMO are written such that the tracer content is conserved by the internal dynamics and physics (equations in flux form). For advection, only the CEN2 scheme ($i.e.$ 2nd order finite different scheme) conserves the global variance of tracer. Nevertheless the other schemes ensure that the global variance decreases ($i.e.$ they are at least slightly diffusive). For diffusion, all the schemes ensure the decrease of the total tracer variance, but the iso-neutral operator. There is generally no strict conservation of mass, as the equation of state is non linear with respect to $T$ and $S$. In practice, the mass is conserved with a very good accuracy. 
+All the numerical schemes used in NEMO are written such that the tracer content 
+is conserved by the internal dynamics and physics (equations in flux form). 
+For advection, only the CEN2 scheme ($i.e.$ $2^{nd}$ order finite different scheme) 
+conserves the global variance of tracer. Nevertheless the other schemes ensure 
+that the global variance decreases ($i.e.$ they are at least slightly diffusive). 
+For diffusion, all the schemes ensure the decrease of the total tracer variance, 
+except the iso-neutral operator. There is generally no strict conservation of mass, 
+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
@@ -553,5 +626,7 @@
 \label{Apdx_C.2.1}
 
-Whatever the advection scheme considered it conserves of the tracer content as all the scheme are written in flux form. Let $\tau$ be the tracer interpolated at velocity point (whatever the interpolation is). The conservation of the tracer content is obtained as follows: 
+Whatever the advection scheme considered it conserves of the tracer content as all 
+the scheme are written in flux form. Let $\tau$ be the tracer interpolated at velocity point 
+(whatever the interpolation is). The conservation of the tracer content is obtained as follows: 
 \begin{flalign*}
 &\int_D \nabla \cdot \left( T \textbf{U} \right)\;dv    &&&\\
@@ -571,5 +646,6 @@
 \end{flalign*}
 
-The conservation of the variance of tracer can be achieved only with the CEN2 scheme. It can be demonstarted as follows:
+The conservation of the variance of tracer can be achieved only with the CEN2 scheme. 
+It can be demonstarted as follows:
 \begin{flalign*}
 &\int\limits_D T\;\nabla \cdot \left( T\; \textbf{U} \right)\;dv &&&\\
@@ -614,5 +690,5 @@
 
 The discrete formulation of the horizontal diffusion of momentum ensures the 
-conservation of potential vorticity and horizontal divergence and the 
+conservation of potential vorticity and the horizontal divergence, and the 
 dissipation of the square of these quantities (i.e. enstrophy and the 
 variance of the horizontal divergence) as well as the dissipation of the 
@@ -623,6 +699,9 @@
 horizontal divergence) and \textit{vice versa}. 
 
-These properties of the horizontal diffusive operator are a direct 
-consequence of properties \eqref{Eq_DOM_curl_grad} and \eqref{Eq_DOM_div_curl}. When the vertical curl of the horizontal diffusion of momentum (discrete sense) is taken, the term associated to the horizontal gradient of the divergence is zero locally. 
+These properties of the horizontal diffusion operator are a direct consequence 
+of properties \eqref{Eq_DOM_curl_grad} and \eqref{Eq_DOM_div_curl}. 
+When the vertical curl of the horizontal diffusion of momentum (discrete sense) 
+is taken, the term associated with the horizontal gradient of the divergence is 
+locally zero. 
 
 % -------------------------------------------------------------------------------------------------------------
@@ -789,6 +868,6 @@
 
 When the horizontal divergence of the horizontal diffusion of momentum 
-(discrete sense) is taken, the term associated to the vertical curl of the 
-vorticity is zero locally, due to (II.1.8). The resulting term conserves the 
+(discrete sense) is taken, the term associated with the vertical curl of the 
+vorticity is zero locally, due to (!!! II.1.8  !!!!!). The resulting term conserves the 
 $\chi$ and dissipates $\chi^2$ when the eddy coefficients are 
 horizontally uniform.
@@ -873,5 +952,7 @@
 
 
-As for the lateral momentum physics, the continuous form of the vertical diffusion of momentum satisfies the several integral constraints. The first two are associated to the conservation of momentum and the dissipation of horizontal kinetic energy:
+As for the lateral momentum physics, the continuous form of the vertical diffusion 
+of momentum satisfies several integral constraints. The first two are associated 
+with the conservation of momentum and the dissipation of horizontal kinetic energy:
 \begin{align*}
  \int\limits_D  
@@ -923,5 +1004,5 @@
 	&&&\\ 
 %
-\intertext{as the horizontal scale factor do not depend on $k$, it follows:}
+\intertext{since the horizontal scale factor does not depend on $k$, it follows:}
 %
 &\equiv - \sum\limits_{i,j,k} 
@@ -982,5 +1063,5 @@
 \end{flalign*}
 If the vertical diffusion coefficient is uniform over the whole domain, the 
-enstrophy is dissipated, i.e.
+enstrophy is dissipated, $i.e.$
 \begin{flalign*}
 \int\limits_D \zeta \, \textbf{k} \cdot \nabla \times 
@@ -1053,5 +1134,8 @@
 	&&\\ 
 \end{flalign*}
-Using the fact that the vertical diffusive coefficients are uniform and that in $z$-coordinates, the vertical scale factors do not depends on $i$ and $j$ so that: $e_{3f} =e_{3u} =e_{3v} =e_{3T} $ and $e_{3w} =e_{3uw} =e_{3vw} $, it follows:
+Using the fact that the vertical diffusion coefficients are uniform, and that in 
+$z$-coordinate, the vertical scale factors do not depend on $i$ and $j$ so 
+that: $e_{3f} =e_{3u} =e_{3v} =e_{3T} $ and $e_{3w} =e_{3uw} =e_{3vw} $, 
+it follows:
 \begin{flalign*}
 \equiv A^{\,vm} \sum\limits_{i,j,k} \zeta \;\delta_k 
@@ -1090,5 +1174,7 @@
 	&&&\\
 \end{flalign*}
-and the square of the horizontal divergence decreases (i.e. the horizontal divergence is dissipated) if vertical diffusion coefficient is uniform over the whole domain:
+and the square of the horizontal divergence decreases ($i.e.$ the horizontal 
+divergence is dissipated) if the vertical diffusion coefficient is uniform over the 
+whole domain:
 
 \begin{flalign*}
@@ -1103,6 +1189,5 @@
 	&&&\\
 \end{flalign*}
-This property is only satisfied in $z$-coordinates:
-
+This property is only satisfied in the $z$-coordinate:
 \begin{flalign*}
 \int\limits_D \chi \;\nabla \cdot 
@@ -1209,7 +1294,11 @@
 \label{Apdx_C.5}
 
-
-
-The numerical schemes used for tracer subgridscale physics are written such that the heat and salt contents are conserved (equations in flux form, second order centered finite differences). As a form flux is used to compute the temperature and salinity, the quadratic form of these quantities (i.e. their variance) globally tends to diminish. As for the advection term, there is generally no strict conservation of mass even if, in practice, the mass is conserved with a very good accuracy. 
+The numerical schemes used for tracer subgridscale physics are written such 
+that the heat and salt contents are conserved (equations in flux form, second 
+order centered finite differences). Since a flux form is used to compute the 
+temperature and salinity, the quadratic form of these quantities (i.e. their variance) 
+globally tends to diminish. As for the advection term, there is generally no strict 
+conservation of mass, even if in practice the mass is conserved to a very high 
+accuracy. 
 
 % -------------------------------------------------------------------------------------------------------------
@@ -1245,5 +1334,5 @@
 \end{flalign*}
 
-In fact, this property is simply resulting from the flux form of the operator. 
+In fact, this property simply results from the flux form of the operator. 
 
 % -------------------------------------------------------------------------------------------------------------
@@ -1253,5 +1342,5 @@
 \label{Apdx_C.5.2}
 
-constraint of dissipation of tracer variance:
+constraint on the dissipation of tracer variance:
 \begin{flalign*}
 \int\limits_D T\;\nabla & \cdot \left( A\;\nabla T \right)\;dv	&&&\\ 
Index: /trunk/DOC/TexFiles/Chapters/Annex_D.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Annex_D.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Annex_D.tex	(revision 1223)
@@ -7,10 +7,10 @@
 
 
-A "model life" is more than ten years. Its software, composed by a few 
-hundred modules, is used by many people who are scientists or students and 
-do not necessary know very well all computer aspects. Moreover, a well 
-thought-out programme is easy to read and understand, less difficult to 
-modify, produces fewer bugs and is easier to maintain. Therefore, it is 
-essential that the model development follows some rules :
+A "model life" is more than ten years. Its software, composed of a few 
+hundred modules, is used by many people who are scientists or students 
+and do not necessarily know every aspect of computing very well. 
+Moreover, a well thought-out program is easier to read and understand, 
+less difficult to modify, produces fewer bugs and is easier to maintain. 
+Therefore, it is essential that the model development follows some rules :
 
 - well planned and designed
@@ -27,5 +27,6 @@
 
 To satisfy part of these aims, \NEMO is written with a coding standard which 
-is close to the ECMWF rules, named DOCTOR \citep{Gibson_TR86}. These rules present some advantages like :
+is close to the ECMWF rules, named DOCTOR \citep{Gibson_TR86}. 
+These rules present some advantages like :
 
 - to provide a well presented program
@@ -51,12 +52,12 @@
 - the detail of input and output interfaces
 
-- the external routines and functions used (if exist)
+- the external routines and functions used (if they exist)
 
-- references (if exist)
+- references (if they exist)
 
-- the author name (s), the date of creation and of updates.
+- the author name(s), the date of creation and any updates.
 
-- Each program is splitted into several well separated sections and 
-sub-sections with a underlined title and specific labelled statements.
+- Each program is split into several well separated sections and 
+sub-sections with an underlined title and specific labelled statements.
 
 - A program has not more than 200 to 300 lines.
@@ -68,46 +69,14 @@
 \label{Apdx_D_coding}
 
-- Use of the universal language \textsc{Fortran} 5 ANSI 77, with non 
-standard extensions, NAMELIST and \textsc{Fortran} 90 (matrix resolution 
-algorithm), and some well-identified particular statements or functions for 
-weak and massive parallelism, and vectorization 
+- Use of the universal language \textsc{Fortran} 90, and try to avoid obsolescent
+features like statement functions, do not use GO TO and EQUIVALENCE statements.
 
-- A comment line begins with a uppercase character C at the first column. A 
-space line must have a C. For the on-line documentation, comments are 
-classified into three levels :
+- A continuation line begins with the character {\$} indented by three spaces 
+compared to the previous line.
 
-Overview, triggered by CCC in columns 1 to 3. Only the title and the purpose 
-of the program are identified like that. This overview documentation can be 
-extracted by the UNIX function : grep -e '\^{}CCC' *
+- Always use a three spaces indentation in DO loop, CASE, or IF-ELSEIF-ELSE-ENDIF 
+statements.
 
-External, triggered by CC in columns 1 to 2, and which correspond to 
-headlines of each programme, extracted by : grep -e '\^{}CC' *
-
-Internal which are all the comments, extracted by : grep -e '\^{}C' *
-
-- Statements GO TO, EQUIVALENCE are forbidden.
-
-- A section is numbered with labels which are in agreement with the 
-paragraph label and increase from the begin to the end of routine. Labels of 
-a hundred ( 200, 201.. 220..) are reserved to a unique section. The 
-\textsc{Fortran} 90 extension syntax DO/ENDDO is used except for multitasked 
-do-loop. In this case labels 1000, 2000, ... are used. The FORMAT statement 
-are labelled with numbers in the range 9000 to 9999.
-
-- A continuation line begins with the character {\$} in column 6.
-
-- All statements begin in column 7 with the following gaps :
-
-2 spaces toward the right in a DO loop.
-
-4 spaces toward the right in IF, ELSEIF, ELSE and ENDIF statements, with 
-only 2 spaces for ELSE and ELSEIF lines. All IF statement must be followed 
-by a ELSE statement.
-
-Some spaces in the continuation line for alignment.
-
-- Use of different labels for each DO loop statement.
-
-- STOP must be well documented with the name of the subroutine or a number.
+- use call to ctl\_stop routine instead of just a STOP.
 
 % ================================================================
@@ -118,7 +87,6 @@
 
 The purpose of the naming conventions is to use prefix letters to classify 
-model variables. These conventions allow to know easily the variable type 
-and to identify them rapidly:
-
+model variables. These conventions allow the variable type to be easily known
+and rapidly identified:
 
 %--------------------------------------------------TABLE--------------------------------------------------
@@ -128,10 +96,10 @@
 \hline  Type \par / Status &   integer&   real&   logical &   character&   double \par precision&   complex \\  
 \hline
-public & 
-\textbf{m n} \par \textit{but not } \par \textbf{nam}& 
-\textbf{a b e f g h o} \textbf{q} \textit{to} \textbf{x} \par but not \par \textbf{sf}& 
-\textbf{l} \par \textit{but not} \par \textbf{lp ld ll}& 
-\textbf{c} \par \textit{but not} \par \textbf{cp cd cl } \par \textbf{com cim}& 
-\textbf{d} \par \textit{but not} \par \textbf{dp dd dl}& 
+public  \par or  \par module variable& 
+\textbf{m n} \par \textit{but not } \par \textbf{nn\_}& 
+\textbf{a b e f g h o} \textbf{q} \textit{to} \textbf{x} \par but not \par \textbf{fs rn\_}& 
+\textbf{l} \par \textit{but not} \par \textbf{lp ld ll ln\_}& 
+\textbf{c} \par \textit{but not} \par \textbf{cp cd cl cn\_}& 
+\textbf{d} \par \textit{but not} \par \textbf{dp dd dl dn\_}& 
 \textbf{y} \par \textit{but not} \par \textbf{yp yd yl} \\
 \hline
@@ -168,5 +136,14 @@
 \textbf{yp} \\
 \hline
-statement \par function& 
+
+namelist&
+\textbf{nn\_}& 
+\textbf{rn\_}& 
+\textbf{ln\_}& 
+\textbf{cn\_}& 
+\textbf{dn\_}& 
+\\
+\hline
+CPP \par macro& 
 \textbf{kf}& 
 \textbf{sf} \par & 
Index: /trunk/DOC/TexFiles/Chapters/Annex_E.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Annex_E.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Annex_E.tex	(revision 1223)
@@ -29,7 +29,11 @@
 - \frac{1}{2}\, |U|_{i+1/2} \;\frac{1}{6} \;\delta_{i+1/2}[\tau"_i]
 \end{equation}
-where $U_{i+1/2} = e_{1u}\,e_{3u}\,u_{i+1/2}$ and $\tau "_i =\delta _i \left[ {\delta _{i+1/2} \left[ \tau \right]} \right]$. By choosing this expression for $\tau "$ we consider a fourth order approximation of $\partial_i^2$ with a constant i-grid spacing ($\Delta i=1$). 
-
-Alternative choice: introduce the scale factors:  $\tau "_i =\frac{e_{1T}}{e_{2T}\,e_{3T}}\delta _i \left[ \frac{e_{2u} e_{3u} }{e_{1u} }\delta _{i+1/2}[\tau] \right]$.
+where $U_{i+1/2} = e_{1u}\,e_{3u}\,u_{i+1/2}$ and 
+$\tau "_i =\delta _i \left[ {\delta _{i+1/2} \left[ \tau \right]} \right]$. 
+By choosing this expression for $\tau "$ we consider a fourth order approximation 
+of $\partial_i^2$ with a constant i-grid spacing ($\Delta i=1$). 
+
+Alternative choice: introduce the scale factors:  
+$\tau "_i =\frac{e_{1T}}{e_{2T}\,e_{3T}}\delta _i \left[ \frac{e_{2u} e_{3u} }{e_{1u} }\delta _{i+1/2}[\tau] \right]$.
 
 
@@ -48,11 +52,12 @@
 scheme when \np{ln\_traadv\_ubs}=T.
 
-For stability reasons, in \eqref{Eq_tra_adv_ubs}, the first term which corresponds to a 
-second order centred scheme is evaluated using the \textit{now} velocity (centred in 
-time) while the second term which is the diffusive part of the scheme, is 
-evaluated using the \textit{before} velocity (forward in time. This is discussed by \citet{Webb1998} in the context of the Quick advection scheme. UBS and QUICK 
+For stability reasons, in \eqref{Eq_tra_adv_ubs}, the first term which corresponds 
+to a second order centred scheme is evaluated using the \textit{now} velocity 
+(centred in time) while the second term which is the diffusive part of the scheme, 
+is evaluated using the \textit{before} velocity (forward in time. This is discussed 
+by \citet{Webb1998} in the context of the Quick advection scheme. UBS and QUICK 
 schemes only differ by one coefficient. Substituting 1/6 with 1/8 in 
-(\ref{Eq_tra_adv_ubs}) leads to the QUICK advection scheme \citep{Webb1998}. This 
-option is not available through a namelist parameter, since the 1/6 
+(\ref{Eq_tra_adv_ubs}) leads to the QUICK advection scheme \citep{Webb1998}. 
+This option is not available through a namelist parameter, since the 1/6 
 coefficient is hard coded. Nevertheless it is quite easy to make the 
 substitution in \mdl{traadv\_ubs} module and obtain a QUICK scheme
@@ -67,7 +72,7 @@
 NB 2 : In a forthcoming release four options will be proposed for the 
 vertical component used in the UBS scheme. $\tau _w^{ubs}$ will be 
-evaluated using either \textit{(a)} a centred $2^{nd}$ order scheme , or  \textit{(b)} a TVD 
-scheme, or  \textit{(c)} an interpolation based on conservative parabolic splines 
-following \citet{Sacha2005} implementation of UBS in ROMS, 
+evaluated using either \textit{(a)} a centred $2^{nd}$ order scheme , 
+or  \textit{(b)} a TVD scheme, or  \textit{(c)} an interpolation based on conservative 
+parabolic splines following \citet{Sacha2005} implementation of UBS in ROMS, 
 or  \textit{(d)} an UBS. The $3^{rd}$ case has dispersion properties similar to an 
 eight-order accurate conventional scheme.
@@ -88,6 +93,11 @@
 \end{split}
 \end{equation}
-\eqref{Eq_tra_adv_ubs2} has several advantages. First it clearly evidence that the UBS scheme is based on the fourth order scheme to which is added an upstream biased diffusive term. Second, this emphasises that the $4^{th}$ order part have to be evaluated at \emph{now} time step, not only the $2^{th}$ order part as stated above using \eqref{Eq_tra_adv_ubs}. Third, the diffusive term is in fact a biharmonic operator with a eddy coefficient with is simply proportional to the velocity.
-
+\eqref{Eq_tra_adv_ubs2} has several advantages. First it clearly evidence that 
+the UBS scheme is based on the fourth order scheme to which is added an 
+upstream biased diffusive term. Second, this emphasises that the $4^{th}$ order 
+part have to be evaluated at \emph{now} time step, not only the $2^{th}$ order 
+part as stated above using \eqref{Eq_tra_adv_ubs}. Third, the diffusive term is 
+in fact a biharmonic operator with a eddy coefficient with is simply proportional 
+to the velocity.
 
 laplacian diffusion:
@@ -179,5 +189,7 @@
 \end{align}
 \end{subequations}
-As for space operator, the adjoint of the derivation and averaging time operators are $\delta_t^*=\delta_{t+\Delta t/2}$ and $\overline{\cdot}^{\,t\,*}= \overline{\cdot}^{\,t+\Delta/2}$, respectively. 
+As for space operator, the adjoint of the derivation and averaging time operators are 
+$\delta_t^*=\delta_{t+\Delta t/2}$ and $\overline{\cdot}^{\,t\,*}= \overline{\cdot}^{\,t+\Delta/2}$
+, respectively. 
 
 The Leap-frog time stepping given by \eqref{Eq_DOM_nxt} can be defined as:
@@ -187,6 +199,8 @@
 		=         \frac{q^{t+\Delta t}-q^{t-\Delta t}}{2\Delta t}
 \end{equation} 
-Note that \eqref{LF} shows that the leapfrog time step is $\Delta t$, not $2\Delta t$ as it can be found sometime in literature. 
-The leap-Frog time stepping is a second order centered scheme. As such it respects the quadratic invariant in integral forms, $i.e.$ the following continuous property,
+Note that \eqref{LF} shows that the leapfrog time step is $\Delta t$, not $2\Delta t$ 
+as it can be found sometime in literature. 
+The leap-Frog time stepping is a second order centered scheme. As such it respects 
+the quadratic invariant in integral forms, $i.e.$ the following continuous property,
 \begin{equation} \label{Energy}
 \int_{t_0}^{t_1} {q\, \frac{\partial q}{\partial t} \;dt} 
@@ -205,6 +219,9 @@
 	   \equiv \frac{1}{2} \left( {q_{t_1}}^2 - {q_{t_0}}^2 \right)
 \end{split} \end{equation}
-NB here pb of boundary condition when applying the adjoin! In space, setting to 0 the quantity in land area is sufficient to get rid of the boundary condition (equivalently of the boundary value of the integration by part). In time this boundary condition is not physical and \textbf{add something here!!!}
-
-
-
+NB here pb of boundary condition when applying the adjoin! In space, setting to 0 
+the quantity in land area is sufficient to get rid of the boundary condition 
+(equivalently of the boundary value of the integration by part). In time this boundary 
+condition is not physical and \textbf{add something here!!!}
+
+
+
Index: /trunk/DOC/TexFiles/Chapters/Chap_Conservation.tex
===================================================================
--- /trunk/DOC/TexFiles/Chapters/Chap_Conservation.tex	(revision 1222)
+++ /trunk/DOC/TexFiles/Chapters/Chap_Conservation.tex	(revision 1223)
@@ -3,5 +3,5 @@
 % Invariant of the Equations
 % ================================================================
-\chapter{Annex --- Invariants of the Primitive Equations}
+\chapter{Invariants of the Primitive Equations}
 \label{Invariant}
 \minitoc
