Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/DOC/TexFiles/Chapters/Chap_Model_Basics.tex
- Timestamp:
- 2016-01-08T10:35:19+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/DOC/TexFiles/Chapters/Chap_Model_Basics.tex
r3294 r6225 247 247 sufficient to solve a linearized version of (\ref{Eq_PE_ssh}), which still allows 248 248 to take into account freshwater fluxes applied at the ocean surface \citep{Roullet_Madec_JGR00}. 249 Nevertheless, with the linearization, an exact conservation of heat and salt contents is lost. 249 250 250 251 The filtering of EGWs in models with a free surface is usually a matter of discretisation 251 of the temporal derivatives, using the time splitting method \citep{Killworth_al_JPO91, Zhang_Endoh_JGR92} 252 or the implicit scheme \citep{Dukowicz1994}. In \NEMO, we use a slightly different approach 253 developed by \citet{Roullet_Madec_JGR00}: the damping of EGWs is ensured by introducing an 254 additional force in the momentum equation. \eqref{Eq_PE_dyn} becomes: 255 \begin{equation} \label{Eq_PE_flt} 256 \frac{\partial {\rm {\bf U}}_h }{\partial t}= {\rm {\bf M}} 257 - g \nabla \left( \tilde{\rho} \ \eta \right) 258 - g \ T_c \nabla \left( \widetilde{\rho} \ \partial_t \eta \right) 259 \end{equation} 260 where $T_c$, is a parameter with dimensions of time which characterizes the force, 261 $\widetilde{\rho} = \rho / \rho_o$ is the dimensionless density, and $\rm {\bf M}$ 262 represents the collected contributions of the Coriolis, hydrostatic pressure gradient, 263 non-linear and viscous terms in \eqref{Eq_PE_dyn}. 264 265 The new force can be interpreted as a diffusion of vertically integrated volume flux divergence. 266 The time evolution of $D$ is thus governed by a balance of two terms, $-g$ \textbf{A} $\eta$ 267 and $g \, T_c \,$ \textbf{A} $D$, associated with a propagative regime and a diffusive regime 268 in the temporal spectrum, respectively. In the diffusive regime, the EGWs no longer propagate, 269 $i.e.$ they are stationary and damped. The diffusion regime applies to the modes shorter than 270 $T_c$. For longer ones, the diffusion term vanishes. Hence, the temporally unresolved EGWs 271 can be damped by choosing $T_c > \rdt$. \citet{Roullet_Madec_JGR00} demonstrate that 272 (\ref{Eq_PE_flt}) can be integrated with a leap frog scheme except the additional term which 273 has to be computed implicitly. This is not surprising since the use of a large time step has a 274 necessarily numerical cost. Two gains arise in comparison with the previous formulations. 275 Firstly, the damping of EGWs can be quantified through the magnitude of the additional term. 276 Secondly, the numerical scheme does not need any tuning. Numerical stability is ensured as 277 soon as $T_c > \rdt$. 278 279 When the variations of free surface elevation are small compared to the thickness of the first 280 model layer, the free surface equation (\ref{Eq_PE_ssh}) can be linearized. As emphasized 281 by \citet{Roullet_Madec_JGR00} the linearization of (\ref{Eq_PE_ssh}) has consequences on the 282 conservation of salt in the model. With the nonlinear free surface equation, the time evolution 283 of the total salt content is 284 \begin{equation} \label{Eq_PE_salt_content} 285 \frac{\partial }{\partial t}\int\limits_{D\eta } {S\;dv} 286 =\int\limits_S {S\;(-\frac{\partial \eta }{\partial t}-D+P-E)\;ds} 287 \end{equation} 288 where $S$ is the salinity, and the total salt is integrated over the whole ocean volume 289 $D_\eta$ bounded by the time-dependent free surface. The right hand side (which is an 290 integral over the free surface) vanishes when the nonlinear equation (\ref{Eq_PE_ssh}) 291 is satisfied, so that the salt is perfectly conserved. When the free surface equation is 292 linearized, \citet{Roullet_Madec_JGR00} show that the total salt content integrated in the fixed 293 volume $D$ (bounded by the surface $z=0$) is no longer conserved: 294 \begin{equation} \label{Eq_PE_salt_content_linear} 295 \frac{\partial }{\partial t}\int\limits_D {S\;dv} 296 = - \int\limits_S {S\;\frac{\partial \eta }{\partial t}ds} 297 \end{equation} 298 299 The right hand side of (\ref{Eq_PE_salt_content_linear}) is small in equilibrium solutions 300 \citep{Roullet_Madec_JGR00}. It can be significant when the freshwater forcing is not balanced and 301 the globally averaged free surface is drifting. An increase in sea surface height \textit{$\eta $} 302 results in a decrease of the salinity in the fixed volume $D$. Even in that case though, 303 the total salt integrated in the variable volume $D_{\eta}$ varies much less, since 304 (\ref{Eq_PE_salt_content_linear}) can be rewritten as 305 \begin{equation} \label{Eq_PE_salt_content_corrected} 306 \frac{\partial }{\partial t}\int\limits_{D\eta } {S\;dv} 307 =\frac{\partial}{\partial t} \left[ \;{\int\limits_D {S\;dv} +\int\limits_S {S\eta \;ds} } \right] 308 =\int\limits_S {\eta \;\frac{\partial S}{\partial t}ds} 309 \end{equation} 310 311 Although the total salt content is not exactly conserved with the linearized free surface, 312 its variations are driven by correlations of the time variation of surface salinity with the 313 sea surface height, which is a negligible term. This situation contrasts with the case of 314 the rigid lid approximation in which case freshwater forcing is represented by a virtual 315 salt flux, leading to a spurious source of salt at the ocean surface 316 \citep{Huang_JPO93, Roullet_Madec_JGR00}. 317 318 \newpage 319 $\ $\newline % force a new ligne 252 of the temporal derivatives, using a split-explicit method \citep{Killworth_al_JPO91, Zhang_Endoh_JGR92} 253 or the implicit scheme \citep{Dukowicz1994} or the addition of a filtering force in the momentum equation 254 \citep{Roullet_Madec_JGR00}. With the present release, \NEMO offers the choice between 255 an explicit free surface (see \S\ref{DYN_spg_exp}) or a split-explicit scheme strongly 256 inspired the one proposed by \citet{Shchepetkin_McWilliams_OM05} (see \S\ref{DYN_spg_ts}). 257 258 %\newpage 259 %$\ $\newline % force a new ligne 320 260 321 261 % ================================================================ … … 655 595 the surface pressure, is given by: 656 596 \begin{equation} \label{Eq_PE_spg} 657 p_s = \left\{ \begin{split} 658 \rho \,g \,\eta & \qquad \qquad \; \qquad \text{ standard free surface} \\ 659 \rho \,g \,\eta &+ \rho_o \,\mu \,\frac{\partial \eta }{\partial t} \qquad \text{ filtered free surface} 660 \end{split} 661 \right. 597 p_s = \rho \,g \,\eta 662 598 \end{equation} 663 599 with $\eta$ is solution of \eqref{Eq_PE_ssh} … … 773 709 \end{equation} 774 710 775 The equations solved by the ocean model \eqref{Eq_PE} in $s-$coordinate can be written as follows :711 The equations solved by the ocean model \eqref{Eq_PE} in $s-$coordinate can be written as follows (see Appendix~\ref{Apdx_A_momentum}): 776 712 777 713 \vspace{0.5cm} 778 * momentum equation:714 $\bullet$ Vector invariant form of the momentum equation : 779 715 \begin{multline} \label{Eq_PE_sco_u} 780 \frac{ 1}{e_3} \frac{\partial \left( e_3\,u \right)}{\partial t}=716 \frac{\partial u }{\partial t}= 781 717 + \left( {\zeta +f} \right)\,v 782 718 - \frac{1}{2\,e_1} \frac{\partial}{\partial i} \left( u^2+v^2 \right) … … 787 723 \end{multline} 788 724 \begin{multline} \label{Eq_PE_sco_v} 789 \frac{ 1}{e_3} \frac{\partial \left( e_3\,v \right)}{\partial t}=725 \frac{\partial v }{\partial t}= 790 726 - \left( {\zeta +f} \right)\,u 791 727 - \frac{1}{2\,e_2 }\frac{\partial }{\partial j}\left( u^2+v^2 \right) … … 795 731 + D_v^{\vect{U}} + F_v^{\vect{U}} \quad 796 732 \end{multline} 733 734 \vspace{0.5cm} 735 $\bullet$ Vector invariant form of the momentum equation : 736 \begin{multline} \label{Eq_PE_sco_u} 737 \frac{1}{e_3} \frac{\partial \left( e_3\,u \right) }{\partial t}= 738 + \left( { f + \frac{1}{e_1 \; e_2 } 739 \left( v \frac{\partial e_2}{\partial i} 740 -u \frac{\partial e_1}{\partial j} \right)} \right) \, v \\ 741 - \frac{1}{e_1 \; e_2 \; e_3 } \left( 742 \frac{\partial \left( {e_2 \, e_3 \, u\,u} \right)}{\partial i} 743 + \frac{\partial \left( {e_1 \, e_3 \, v\,u} \right)}{\partial j} \right) 744 - \frac{1}{e_3 }\frac{\partial \left( { \omega\,u} \right)}{\partial k} \\ 745 - \frac{1}{e_1} \frac{\partial}{\partial i} \left( \frac{p_s + p_h}{\rho _o} \right) 746 + g\frac{\rho }{\rho _o}\sigma _1 747 + D_u^{\vect{U}} + F_u^{\vect{U}} \quad 748 \end{multline} 749 \begin{multline} \label{Eq_PE_sco_v} 750 \frac{1}{e_3} \frac{\partial \left( e_3\,v \right) }{\partial t}= 751 - \left( { f + \frac{1}{e_1 \; e_2} 752 \left( v \frac{\partial e_2}{\partial i} 753 -u \frac{\partial e_1}{\partial j} \right)} \right) \, u \\ 754 - \frac{1}{e_1 \; e_2 \; e_3 } \left( 755 \frac{\partial \left( {e_2 \; e_3 \,u\,v} \right)}{\partial i} 756 + \frac{\partial \left( {e_1 \; e_3 \,v\,v} \right)}{\partial j} \right) 757 - \frac{1}{e_3 } \frac{\partial \left( { \omega\,v} \right)}{\partial k} \\ 758 - \frac{1}{e_2 }\frac{\partial }{\partial j}\left( \frac{p_s+p_h }{\rho _o} \right) 759 + g\frac{\rho }{\rho _o }\sigma _2 760 + D_v^{\vect{U}} + F_v^{\vect{U}} \quad 761 \end{multline} 762 797 763 where the relative vorticity, \textit{$\zeta $}, the surface pressure gradient, and the hydrostatic 798 764 pressure have the same expressions as in $z$-coordinates although they do not represent 799 765 exactly the same quantities. $\omega$ is provided by the continuity equation 800 766 (see Appendix~\ref{Apdx_A}): 801 802 767 \begin{equation} \label{Eq_PE_sco_continuity} 803 768 \frac{\partial e_3}{\partial t} + e_3 \; \chi + \frac{\partial \omega }{\partial s} = 0 … … 809 774 810 775 \vspace{0.5cm} 811 *tracer equations:776 $\bullet$ tracer equations: 812 777 \begin{multline} \label{Eq_PE_sco_t} 813 778 \frac{1}{e_3} \frac{\partial \left( e_3\,T \right) }{\partial t}= … … 1024 989 1025 990 The $\tilde{z}$-coordinate has been developed by \citet{Leclair_Madec_OM10s}. 1026 It is not available in the current version of \NEMO. 991 It is available in \NEMO since the version 3.4. Nevertheless, it is currently not robust enough 992 to be used in all possible configurations. Its use is therefore not recommended. 993 We 1027 994 1028 995 \newpage … … 1180 1147 ocean (see Appendix~\ref{Apdx_B}). 1181 1148 1149 For \textit{iso-level} diffusion, $r_1$ and $r_2 $ are zero. $\Re $ reduces to the identity 1150 in the horizontal direction, no rotation is applied. 1151 1182 1152 For \textit{geopotential} diffusion, $r_1$ and $r_2 $ are the slopes between the 1183 geopotential and computational surfaces: in $z$-coordinates they are zero 1184 ($r_1 = r_2 = 0$) while in $s$-coordinate (including $\textit{z*}$ case) they are 1185 equal to $\sigma _1$ and $\sigma _2$, respectively (see \eqref{Eq_PE_sco_slope} ). 1153 geopotential and computational surfaces: they are equal to $\sigma _1$ and $\sigma _2$, 1154 respectively (see \eqref{Eq_PE_sco_slope} ). 1186 1155 1187 1156 For \textit{isoneutral} diffusion $r_1$ and $r_2$ are the slopes between the isoneutral … … 1235 1204 The lateral fourth order tracer diffusive operator is defined by: 1236 1205 \begin{equation} \label{Eq_PE_bilapT} 1237 D^{lT}=\Delta \left( {A^{lT}\;\Delta T}\right)1238 \qquad \text{where} \ D^{lT}=\Delta \left( {A^{lT}\;\Delta T} \right)1206 D^{lT}=\Delta \left( \;\Delta T \right) 1207 \qquad \text{where} \;\; \Delta \bullet = \nabla \left( {\sqrt{B^{lT}\,}\;\Re \;\nabla \bullet} \right) 1239 1208 \end{equation} 1240 1241 1209 It is the second order operator given by \eqref{Eq_PE_iso_tensor} applied twice with 1242 the eddy diffusion coefficient correctly placed.1210 the harmonic eddy diffusion coefficient set to the square root of the biharmonic one. 1243 1211 1244 1212 … … 1262 1230 1263 1231 Such a formulation ensures a complete separation between the vorticity and 1264 horizontal divergence fields (see Appendix~\ref{Apdx_C}). Unfortunately, it is not1265 available for geopotential diffusion in $s-$coordinates and for isoneutral1266 diffusion in both $z$- and $s$-coordinates ($i.e.$ when a rotation is required).1267 In these two cases, the $u$ and $v-$fields are considered as independent scalar1268 fields, so that the diffusive operator is given by:1232 horizontal divergence fields (see Appendix~\ref{Apdx_C}). 1233 Unfortunately, it is only available in \textit{iso-level} direction. 1234 When a rotation is required ($i.e.$ geopotential diffusion in $s-$coordinates 1235 or isoneutral diffusion in both $z$- and $s$-coordinates), the $u$ and $v-$fields 1236 are considered as independent scalar fields, so that the diffusive operator is given by: 1269 1237 \begin{equation} \label{Eq_PE_lapU_iso} 1270 1238 \begin{split} 1271 D_u^{l{\rm {\bf U}}} &= \nabla .\left( { \Re \;\nabla u} \right) \\1272 D_v^{l{\rm {\bf U}}} &= \nabla .\left( { \Re \;\nabla v} \right)1239 D_u^{l{\rm {\bf U}}} &= \nabla .\left( {A^{lm} \;\Re \;\nabla u} \right) \\ 1240 D_v^{l{\rm {\bf U}}} &= \nabla .\left( {A^{lm} \;\Re \;\nabla v} \right) 1273 1241 \end{split} 1274 1242 \end{equation} … … 1282 1250 1283 1251 As for tracers, the fourth order momentum diffusive operator along $z$ or $s$-surfaces 1284 is a re-entering second order operator \eqref{Eq_PE_lapU} or \eqref{Eq_PE_lapU} 1285 with the eddy viscosity coefficient correctly placed: 1286 1287 geopotential diffusion in $z$-coordinate: 1288 \begin{equation} \label{Eq_PE_bilapU} 1289 \begin{split} 1290 {\rm {\bf D}}^{l{\rm {\bf U}}} &=\nabla _h \left\{ {\;\nabla _h {\rm {\bf 1291 .}}\left[ {A^{lm}\,\nabla _h \left( \chi \right)} \right]\;} 1292 \right\}\; \\ 1293 &+\nabla _h \times \left\{ {\;{\rm {\bf k}}\cdot \nabla \times 1294 \left[ {A^{lm}\,\nabla _h \times \left( {\zeta \;{\rm {\bf k}}} \right)} 1295 \right]\;} \right\} 1296 \end{split} 1297 \end{equation} 1298 1299 \gmcomment{ change the position of the coefficient, both here and in the code} 1300 1301 geopotential diffusion in $s$-coordinate: 1302 \begin{equation} \label{Eq_bilapU_iso} 1303 \left\{ \begin{aligned} 1304 D_u^{l{\rm {\bf U}}} =\Delta \left( {A^{lm}\;\Delta u} \right) \\ 1305 D_v^{l{\rm {\bf U}}} =\Delta \left( {A^{lm}\;\Delta v} \right) 1306 \end{aligned} \right. 1307 \quad \text{where} \quad 1308 \Delta \left( \bullet \right) = \nabla \cdot \left( \Re \nabla(\bullet) \right) 1309 \end{equation} 1310 1252 is a re-entering second order operator \eqref{Eq_PE_lapU} or \eqref{Eq_PE_lapU_iso} 1253 with the harmonic eddy diffusion coefficient set to the square root of the biharmonic one. 1254
Note: See TracChangeset
for help on using the changeset viewer.