New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 2226 for branches/DEV_r1826_DOC/DOC – NEMO

Ignore:
Timestamp:
2010-10-12T16:23:34+02:00 (14 years ago)
Author:
gm
Message:

ticket:#658 update ZDF with TKE time-stepping and GLS scheme

Location:
branches/DEV_r1826_DOC/DOC/TexFiles
Files:
1 added
2 edited

Legend:

Unmodified
Added
Removed
  • branches/DEV_r1826_DOC/DOC/TexFiles/Biblio/Biblio.bib

    r2213 r2226  
    1212@STRING{AW = {Addison-Wesley}} 
    1313 
     14@STRING{CP = {Clarendon Press}} 
     15 
    1416@STRING{CD = {Clim. Dyn.}} 
    1517 
    16 @STRING{CP = {Clarendon Press}} 
    17  
    1818@STRING{CUP = {Cambridge University Press}} 
    1919 
     20@STRING{CSR = {Cont. Shelf Res.}} 
     21 
    2022@STRING{D = {Dover Publications}} 
    2123 
     
    8385 
    8486@STRING{Recherche = {La Recherche}} 
     87 
     88@STRING{RGSP = {Rev. Geophys. Space Phys.}} 
    8589 
    8690@STRING{Science = {Science}} 
     
    387391  author = {D. Bernie and E. Guilyardi and G. Madec and J. M. Slingo and S. J. 
    388392   Woolnough}, 
    389   title = {Impact of resolving the diurnal cycle in an oceanatmosphere GCM. 
     393  title = {Impact of resolving the diurnal cycle in an ocean--atmosphere GCM. 
    390394   Part 2: A diurnally coupled CGCM}, 
    391395  journal = CD, 
     
    402406  author = {D. Bernie and E. Guilyardi and G. Madec and J. M. Slingo and S. J. 
    403407   Woolnough}, 
    404   title = {Impact of resolving the diurnal cycle in an oceanatmosphere GCM. 
     408  title = {Impact of resolving the diurnal cycle in an ocean--atmosphere GCM. 
    405409   Part 1: a diurnally forced OGCM}, 
    406410  journal = CD, 
     
    789793 
    790794@TECHREPORT{Chanut2005, 
    791  author = {J. Chanut} 
    792  year = {2005} 
    793  institution = {European Union: Marine Environment and Security for the European Area (MERSEA) Integrated Project} 
    794  title = {Nesting code for NEMO} 
    795  note = {MERSEA-WP09-MERCA-TASK-9.1.1} 
     795 author = {J. Chanut}, 
     796 title = {Nesting code for NEMO}, 
     797 year = {2005}, 
     798 institution = {European Union: Marine Environment and Security for the European Area (MERSEA) Integrated Project}, 
     799note = {MERSEA-WP09-MERCA-TASK-9.1.1} 
    796800}  
    797801 
     
    805809  volume = {33}, 
    806810  pages = {2504-2526}, 
    807   file = {:Users/mlelod/Documents/Biblio/vertical_coordinates/Chassignet_et_al_JPO_2003.pdf:PDF}, 
    808811  timestamp = {2010.02.01} 
    809812} 
     
    826829} 
    827830 
     831@ARTICLE{Canuto_2001, 
     832  author = {V. M. Canuto and A. Howard and Y. Cheng and M. S. Dubovikov}, 
     833  title = {Ocean turbulence. PartI: One-point closure model-momentum and heat vertical diffusivities}, 
     834  journal = JPO, 
     835  year = {2001}, 
     836  volume = {24, 12}, 
     837  pages = {2546--2559}, 
     838  owner = {gr}, 
     839  timestamp = {2010.09.09} 
     840} 
     841 
    828842@ARTICLE{Cox1987, 
    829843  author = {M. Cox}, 
     
    835849  owner = {gm}, 
    836850  timestamp = {2007.08.03} 
     851} 
     852 
     853@ARTICLE{Craig_Banner_1994, 
     854  author = {P. D. Banner and M. L. Banner}, 
     855  title = {Modeling wave-enhanced turbulence in the ocean surface layer}, 
     856  journal = JPO, 
     857  year = {1994}, 
     858  volume = {24, 12}, 
     859  pages = {2546--2559}, 
     860  owner = {g5}, 
     861  timestamp = {2010.09.09} 
    837862} 
    838863 
     
    863888@ARTICLE{Dandonneau_al_S04, 
    864889  author = {Y. Dandonneau and C. Menkes and T. Gorgues and G. Madec}, 
    865   title = {Reply to Peter Killworth, 2004 : « Comment on the Oceanic Rossby 
    866    Waves acting as a “Hay Rake” for ecosystem by-products »}, 
     890  title = {Reply to Peter Killworth, 2004 : '' Comment on the Oceanic Rossby 
     891   Waves acting as a “Hay Rake” for ecosystem by-products ''}, 
    867892  journal = {Science}, 
    868893  year = {2004}, 
     
    873898} 
    874899 
    875 @ARTICLE{Davies_QJRMetSoc76, 
    876  author = {H.C. Davies}  
    877  title = {A lateral boundary formulation for multi-level prediction models} 
    878  year = {1976} 
    879  journal = {Quarterly Journal of the Royal Meteorological Society} 
    880  volume = {102} 
     900@ARTICLE{Davies_QJRMS76, 
     901 author = {H.C. Davies}, 
     902 title = {A lateral boundary formulation for multi-level prediction models}, 
     903 year = {1976}, 
     904 journal = {QJRMS}, 
     905 volume = {102}, 
    881906 pages = {405--418} 
    882907} 
     
    898923  journal = {In \textit{Modeling the Earth's Climate and its Variability}, Les 
    899924   Houches, Session, LXVII 1997, 
    900     
    901925   Eds. W. R. Holland, S. Joussaume and F. David, Elsevier Science}, 
    902926  year = {2000}, 
     
    11251149} 
    11261150 
    1127 @ARTICLE{Engedahl1995, 
    1128  author = {H. Engerdahl} 
    1129  year = {1995} 
    1130  title = {Use of the flow relaxation scheme in a three-dimensional baroclinic ocean model with realistic topography} 
    1131  journal = {Tellus} 
    1132  volume = {47A} 
    1133  pages = {365--382} 
     1151@ARTICLE{Engerdahl_Tel95, 
     1152 author = {H. Engerdahl}, 
     1153 year = {1995}, 
     1154 title = {Use of the flow relaxation scheme in a three-dimensional baroclinic ocean model with realistic topography}, 
     1155 journal = {Tellus}, 
     1156 volume = {47A}, 
     1157 pages = {365--382}, 
    11341158} 
    11351159 
     
    11871211 
    11881212@ARTICLE{Flather1976, 
    1189  author = {R.A. Flather} 
    1190  year = {1976} 
    1191  title = {A tidal model of the north-west European continental shelf} 
    1192  journal = {Memoires de la Societ\'{e} Royale des Sciences de Li\`{e}ge} 
    1193  volume = {6} 
     1213 author = {R.A. Flather}, 
     1214 year = {1976}, 
     1215 title = {A tidal model of the north-west European continental shelf}, 
     1216 journal = {Memoires de la Societ\'{e} Royale des Sciences de Li\`{e}ge}, 
     1217 volume = {6}, 
    11941218 pages = {141--164} 
     1219} 
     1220 
     1221@ARTICLE{Flather_JPO94, 
     1222 author = {R.A. Flather}, 
     1223 year = {1994}, 
     1224 title = {A storm surge prediction model for the northern Bay of Bengal with application to the cyclone disaster in April 1991}, 
     1225 journal = {JPO}, 
     1226 volume = {24}, 
     1227 pages = {172--190} 
    11951228} 
    11961229 
     
    12061239  owner = {gm}, 
    12071240  timestamp = {2007.08.04} 
     1241} 
     1242 
     1243@ARTICLE{Galperin_al_JAS88, 
     1244  author = {B. Galperin and L. H. Kantha and S. Hassid and A. Rosati}, 
     1245  title = {A quasi-equilibrium turbulent energy model for geophysical flows}, 
     1246  journal = JAS, 
     1247  year = {1988}, 
     1248  volume = {45}, 
     1249  pages = {55--62}, 
     1250  owner = {gr}, 
     1251  timestamp = {2010.09.09} 
    12081252} 
    12091253 
     
    17791823  owner = {gm}, 
    17801824  timestamp = {2008.08.31} 
     1825} 
     1826 
     1827@ARTICLE{Kantha_Clayson_1994, 
     1828  author = {L. H. Kantha and C. A. Clayson}, 
     1829  title = {An improved mixed layer model for geophysical applications}, 
     1830  journal = JGR, 
     1831  year = {1994}, 
     1832  volume = {99}, 
     1833  pages = {25,235--25,266}, 
     1834  owner = {gr}, 
     1835  timestamp = {2010.09.09} 
     1836} 
     1837 
     1838@ARTICLE{Kantha_Carniel_CSR05, 
     1839  author = {L. Kantha and S. Carniel}, 
     1840  title = {Comment on ''Generic length-scale equation for geophysical turbulence models'' by L. Umlauf and H. Burchard}, 
     1841  journal = {Journal of Marine Systems},  
     1842  year = {2005}, 
     1843  volume = {61}, 
     1844  pages = {693--702}, 
     1845  owner = {gm}, 
     1846  timestamp = {2010.09.09} 
    17811847} 
    17821848 
     
    26152681} 
    26162682 
     2683@ARTICLE{Mellor_Yamada_1982, 
     2684  author = {G. L. Mellor and T. Yamada}, 
     2685  title = {Development of a turbulence closure model for geophysical fluid problems}, 
     2686  journal = RGSP, 
     2687  year = {1982}, 
     2688  volume = {20}, 
     2689  pages = {851-875}, 
     2690  owner = {gr}, 
     2691  timestamp = {2010.09.09} 
     2692} 
     2693 
    26172694@ARTICLE{Menkes_al_JPO06, 
    26182695  author = {C. Menkes and J. Vialard and S C. Kennan and J.-P. Boulanger and 
     
    29563033  timestamp = {2009.08.20}, 
    29573034  url = {http://dx.doi.org/10.1029/2002GL016003} 
     3035} 
     3036 
     3037@ARTICLE{Rodi_1987, 
     3038  author = {W. Rodi}, 
     3039  title = {Examples of calculation methods for flow and mixing in stratified fluids}, 
     3040  journal = JGR, 
     3041  year = {1987}, 
     3042  volume = {92, C5}, 
     3043  pages = {5305--5328}, 
     3044  owner = {gr}, 
     3045  timestamp = {2010.09.09} 
    29583046} 
    29593047 
     
    34773565} 
    34783566 
     3567@ARTICLE{Umlauf_Burchard_JMS03, 
     3568  author = {L. Umlauf and H. Burchard}, 
     3569  title = {A generic length-scale equation for geophysical turbulence models}, 
     3570  journal = {JMS},  
     3571  year = {2003}, 
     3572  volume = {61}, 
     3573  pages = {235--265}, 
     3574  number = {2}, 
     3575  owner = {gr}, 
     3576  timestamp = {2010.09.09} 
     3577} 
     3578 
     3579@ARTICLE{Umlauf_Burchard_CSR05, 
     3580  author = {L. Umlauf and H. Burchard}, 
     3581  title = {Second-order turbulence closure models for geophysical boundary layers. A review of recent work}, 
     3582  journal = {Journal of Marine Systems},  
     3583  year = {2005}, 
     3584  volume = {25}, 
     3585  pages = {795--827}, 
     3586  owner = {gm}, 
     3587  timestamp = {2010.09.09} 
     3588} 
     3589 
    34793590@BOOK{UNESCO1983, 
    34803591  title = {Algorithms for computation of fundamental property of sea water}, 
     
    36123723  owner = {gm}, 
    36133724  timestamp = {2010.04.14} 
     3725} 
     3726 
     3727@ARTICLE{Wilcox_1988, 
     3728  author = {D. C. Wilcox}, 
     3729  title = {Reassessment of the scale-determining equation for advanced turbulence models}, 
     3730  journal = {AIAA journal}, 
     3731  year = {1988}, 
     3732  volume = {26, 11}, 
     3733  pages = {1299--1310}, 
     3734  owner = {gr}, 
     3735  timestamp = {2010.09.09} 
    36143736} 
    36153737 
  • branches/DEV_r1826_DOC/DOC/TexFiles/Chapters/Chap_ZDF.tex

    r2164 r2226  
    5959\begin{align*}  
    6060A_u^{vm} = A_v^{vm} &= 1.2\ 10^{-4}~m^2.s^{-1}  \\ 
    61 \\ 
    6261A^{vT} = A^{vS} &= 1.2\ 10^{-5}~m^2.s^{-1} 
    6362\end{align*} 
     
    9190   \left\{      \begin{aligned} 
    9291         A^{vT} &= \frac {A_{ric}^{vT}}{\left( 1+a \; Ri \right)^n} + A_b^{vT}       \\ 
    93          \\ 
    9492         A^{vm} &= \frac{A^{vT}        }{\left( 1+ a \;Ri  \right)   } + A_b^{vm} 
    9593   \end{aligned}    \right. 
     
    126124\begin{equation} \label{Eq_zdftke_e} 
    127125\frac{\partial \bar{e}}{\partial t} =  
    128 \frac{A^{vm}}{{e_3}^2 }\;\left[ {\left( {\frac{\partial u}{\partial k}} \right)^2 
    129                      +\left( {\frac{\partial v}{\partial k}} \right)^2} \right] 
    130 -A^{vT}\,N^2 
     126\frac{K_m}{{e_3}^2 }\;\left[ {\left( {\frac{\partial u}{\partial k}} \right)^2 
     127                    +\left( {\frac{\partial v}{\partial k}} \right)^2} \right] 
     128-K_\rho\,N^2 
    131129+\frac{1}{e_3} \;\frac{\partial }{\partial k}\left[ {\frac{A^{vm}}{e_3 } 
    132130            \;\frac{\partial \bar{e}}{\partial k}} \right] 
     
    135133\begin{equation} \label{Eq_zdftke_kz} 
    136134   \begin{split} 
    137          A^{vm} &= C_k\  l_k\  \sqrt {\bar{e}\; }     \\ 
    138          A^{vT} &= A^{vm} / P_{rt} 
     135         K_m &= C_k\  l_k\  \sqrt {\bar{e}\; }     \\ 
     136         K_\rho &= A^{vm} / P_{rt} 
    139137   \end{split} 
    140138\end{equation} 
    141139where $N$ is the local Brunt-Vais\"{a}l\"{a} frequency (see \S\ref{TRA_bn2}),  
    142140$l_{\epsilon }$ and $l_{\kappa }$ are the dissipation and mixing length scales,  
    143 $P_{rt}$ is the Prandtl number. The constants $C_k =  0.1$ and  
    144 $C_\epsilon = \sqrt {2} /2$  $\approx 0.7$ are designed to deal with vertical mixing  
    145 at any depth \citep{Gaspar1990}. They are set through namelist parameters  
    146 \np{nn\_ediff} and \np{nn\_ediss}. $P_{rt}$ can be set to unity or, following  
    147 \citet{Blanke1993}, be a function of the local Richardson number, $R_i$: 
     141$P_{rt}$ is the Prandtl number, $K_m$ and $K_\rho$ are the vertical eddy viscosity  
     142and diffusivity coefficients. The constants $C_k =  0.1$ and $C_\epsilon = \sqrt {2} /2$   
     143$\approx 0.7$ are designed to deal with vertical mixing at any depth \citep{Gaspar1990}.  
     144They are set through namelist parameters \np{nn\_ediff} and \np{nn\_ediss}.  
     145$P_{rt}$ can be set to unity or, following \citet{Blanke1993}, be a function  
     146of the local Richardson number, $R_i$: 
    148147\begin{align*} \label{Eq_prt} 
    149148P_{rt} = \begin{cases} 
     
    165164which is valid in a stable stratified region with constant values of the Brunt- 
    166165Vais\"{a}l\"{a} frequency. The resulting length scale is bounded by the distance  
    167 to the surface or to the bottom (\np{nn\_mxl}=0) or by the local vertical scale factor  
    168 (\np{nn\_mxl}=1). \citet{Blanke1993} notice that this simplification has two major  
     166to the surface or to the bottom (\np{nn\_mxl} = 0) or by the local vertical scale factor  
     167(\np{nn\_mxl} = 1). \citet{Blanke1993} notice that this simplification has two major  
    169168drawbacks: it makes no sense for locally unstable stratification and the  
    170169computation no longer uses all the information contained in the vertical density  
    171170profile. To overcome these drawbacks, \citet{Madec1998} introduces the  
    172 \np{nn\_mxl}=2 or 3 cases, which add an extra assumption concerning the vertical  
     171\np{nn\_mxl} = 2 or 3 cases, which add an extra assumption concerning the vertical  
    173172gradient of the computed length scale. So, the length scales are first evaluated  
    174173as in \eqref{Eq_tke_mxl0_1} and then bounded such that: 
     
    227226to $\sqrt{2}/2~10^{-6}~m^2.s^{-2}$. This allows the subsequent formulations  
    228227to match that of \citet{Gargett1984} for the diffusion in the thermocline and  
    229 deep ocean :  $(A^{vT} = 10^{-3} / N)$.  
    230 In addition, a cut-off is applied on $A^{vm}$ and $A^{vT}$ to avoid numerical  
     228deep ocean :  $K_\rho = 10^{-3} / N$.  
     229In addition, a cut-off is applied on $K_m$ and $K_\rho$ to avoid numerical  
    231230instabilities associated with too weak vertical diffusion. They must be  
    232231specified at least larger than the molecular values, and are set through  
     
    236235%        TKE Turbulent Closure Scheme : new organization to energetic considerations 
    237236% ------------------------------------------------------------------------------------------------------------- 
    238 \subsection{TKE Turbulent Closure Scheme - time integration (\key{zdftke})} 
     237\subsection{TKE discretization considerations (\key{zdftke})} 
    239238\label{ZDF_tke_ene} 
    240239 
     
    249248The production of turbulence by vertical shear (the first term of the right hand side  
    250249of \eqref{Eq_zdftke_e}) should balance the loss of kinetic energy associated with 
    251 the vertical momentum diffusion (first line in \eqref{Eq_PE_zdf}).  
    252 The total loss of kinetic energy (in 1D for the demonstration)  
    253 due to this term is obtained by multiply this quantity by $u^n$ and verticaly integrating:    
    254  
     250the vertical momentum diffusion (first line in \eqref{Eq_PE_zdf}). To do so a special care  
     251have to be taken for both the time and space discretization of the TKE equation \citep{Burchard_OM02}. 
     252 
     253Let us first address the time stepping issue. Fig.~\ref{Fig_TKE_time_scheme} shows  
     254how the two-level Leap-Frog time stepping of the momentum and tracer equations interplays  
     255with the one-level forward time stepping of TKE equation. With this framework, the total loss  
     256of kinetic energy (in 1D for the demonstration) due to the vertical momentum diffusion is  
     257obtained by multiplying this quantity by $u^t$ and summing the result vertically:    
    255258\begin{equation} \label{Eq_energ1} 
    256 \int_{k_b}^{k_s} {u^t \frac{1}{e_3} 
    257                                  \frac{\partial   }  
    258                                         {\partial k} \left( \frac{A^{vm}}{e_3}  
    259                                                                  \frac{\partial{u^{t+1}}} 
    260                                                                         {\partial k           }     \right) \; e_3 \; dk } 
    261 = \left[  \frac{u^t}{e_3}   A^{vm} \frac{\partial{u^{t+1}}}{\partial k} \right]_{k_b}^{k_s} 
    262  - \int_{k_b}^{k_s}{\frac{A^{vm}}{{e_3}}\frac{\partial{u^t}}{\partial k}\frac{\partial{u^{t+1}}}{\partial k}} \ dk 
    263 \end{equation} 
    264  
    265 The first term of the right hand side of \eqref{Eq_energ1} represents the kinetic  
    266 energy transfer at the surface (atmospheric forcing) and at the bottom (friction effect).  
    267 The second term is always negative and have to balance the term of \eqref{Eq_zdftke_e}  
    268 previously identified. 
    269  
    270 The sink term (possibly a source term in statically unstable situations) of turbulence  
    271 by buoyancy (second term of the right hand side of \eqref{Eq_zdftke_e}) must balance 
    272 the source of potential energy associated with the vertical diffusion  
    273 in the density equation (second line in \eqref{Eq_PE_zdf}). The source of potential  
    274 energy (in 1D for the demonstration) due to this term is obtained by multiply this quantity  
    275 by $gz{\rho_r}^{-1}$ and verticaly integrating: 
    276  
     259\begin{split} 
     260\int_{-H}^{\eta}  u^t \,\partial_z &\left( {K_m}^t \,(\partial_z u)^{t+\rdt}  \right) \,dz   \\ 
     261&= \Bigl[  u^t \,{K_m}^t \,(\partial_z u)^{t+\rdt} \Bigr]_{-H}^{\eta}           
     262 - \int_{-H}^{\eta}{ {K_m}^t \,\partial_z{u^t} \,\partial_z u^{t+\rdt} \,dz } 
     263\end{split} 
     264\end{equation} 
     265Here, the vertical diffusion of momentum is discretized backward in time  
     266with a coefficient, $K_m$, known at time $t$ (Fig.~\ref{Fig_TKE_time_scheme}),  
     267as it is required when using the TKE scheme (see \S\ref{STP_forward_imp}).  
     268The first term of the right hand side of \eqref{Eq_energ1} represents the kinetic energy  
     269transfer at the surface (atmospheric forcing) and at the bottom (friction effect).  
     270The second term is always negative. It is the dissipation rate of kinetic energy,  
     271and thus minus the shear production rate of $\bar{e}$. \eqref{Eq_energ1}  
     272implies that, to be energetically consistent, the production rate of $\bar{e}$  
     273used to compute $(\bar{e})^t$ (and thus ${K_m}^t$) should be expressed as  
     274${K_m}^{t-\rdt}\,(\partial_z u)^{t-\rdt} \,(\partial_z u)^t$ (and not by the more straightforward  
     275$K_m \left( \partial_z u \right)^2$ expression taken at time $t$ or $t-\rdt$). 
     276 
     277A similar consideration applies on the destruction rate of $\bar{e}$ due to stratification  
     278(second term of the right hand side of \eqref{Eq_zdftke_e}). This term  
     279must balance the input of potential energy resulting from vertical mixing.  
     280The rate of change of potential energy (in 1D for the demonstration) due vertical  
     281mixing is obtained by multiplying vertical density diffusion  
     282tendency by $g\,z$ and and summing the result vertically: 
    277283\begin{equation} \label{Eq_energ2} 
    278 \begin{aligned} 
    279 \int_{k_b}^{k_s}{\frac{g\;z}{e_3} \frac{\partial }{\partial k} 
    280                         \left( \frac{A^{vT}}{e_3}  
    281                         \frac{\partial{\rho^{t+1}}}{\partial k}   \right)} \; e_3 \; dk 
    282 =\left[  g\;z \frac{A^{vT}}{e_3} 
    283                  \frac{\partial{\rho^{t+1}}}{\partial k} \right]_{k_b}^{k_s}   
    284 - \int_{k_b}^{k_s}{  \frac{A^{vT}}{e_3} g \frac{\partial{\rho^{t+1}}}{\partial k}} \; dk\\ 
    285 \\ 
    286 = - \left[  z\,A^{vT} {N_{t+1}}^2 \right]_{k_b}^{k_s} 
    287 + \int_{k_b}^{k_s}{  \rho^{t+1} \; A^{vT}{N_{t+1}}^2 \; e_3 \; dk  }\\ 
    288 \end{aligned} 
    289 \end{equation} 
    290 where $N^2_{t+1}$ is  the Brunt-Vais\"{a}l\"{a} frequency at $t+1$  
    291 and noting that $\frac{\partial z}{\partial k} = e_3$. 
    292  
    293 The first term is always zero because theBrunt-Vais\"{a}l\"{a} frequency is set to zero at the  
    294 surface and the bottom. The second term is of opposite sign than the buoyancy term 
    295 identified previously. 
    296  
    297 Under these energetic considerations, \eqref{Eq_zdftke_e} have to be written like this  
    298 to be consistant:  
    299  
     284\begin{split} 
     285\int_{-H}^{\eta} g\,z\,\partial_z &\left( {K_\rho}^t \,(\partial_k \rho)^{t+\rdt}   \right) \,dz    \\ 
     286&= \Bigl[  g\,z \,{K_\rho}^t \,(\partial_z \rho)^{t+\rdt} \Bigr]_{-H}^{\eta}   
     287   - \int_{-H}^{\eta}{ g \,{K_\rho}^t \,(\partial_k \rho)^{t+\rdt} } \,dz   \\ 
     288&= - \Bigl[  z\,{K_\rho}^t \,(N^2)^{t+\rdt} \Bigr]_{-H}^{\eta} 
     289+ \int_{-H}^{\eta}{  \rho^{t+\rdt} \, {K_\rho}^t \,(N^2)^{t+\rdt} \,dz  } 
     290\end{split} 
     291\end{equation} 
     292where we use $N^2 = -g \,\partial_k \rho / (e_3 \rho)$.  
     293The first term of the right hand side of \eqref{Eq_energ2}  is always zero  
     294because there is no diffusive flux through the ocean surface and bottom).  
     295The second term is minus the destruction rate of  $\bar{e}$ due to stratification.  
     296Therefore \eqref{Eq_energ1} implies that, to be energetically consistent, the product  
     297${K_\rho}^{t-\rdt}\,(N^2)^t$ should be used in \eqref{Eq_zdftke_e}, the TKE equation. 
     298 
     299Let us now address the space discretization issue.  
     300The vertical eddy coefficients are defined at $w$-point whereas the horizontal velocity  
     301components are in the centre of the side faces of a $t$-box in staggered C-grid  
     302(Fig.\ref{Fig_cell}). A space averaging is thus required to obtain the shear TKE production term. 
     303By redoing the \eqref{Eq_energ1} in the 3D case, it can be shown that the product of  
     304eddy coefficient by the shear at $t$ and $t-\rdt$ must be performed prior to the averaging. 
     305Furthermore, the possible time variation of $e_3$ (\key{vvl} case) have to be taken into  
     306account. 
     307 
     308The above energetic considerations leads to  
     309the following final discrete form for the TKE equation: 
    300310\begin{equation} \label{Eq_zdftke_ene} 
    301 \frac{\partial \bar{e}}{\partial t} = 
    302 \frac{A^{vm}}{{e_3}^2 }\;\left[ {\left( {\frac{\partial u^a}{\partial k}} \right)\left( {\frac{\partial u^n}{\partial k}} \right) 
    303                                                         +\left( {\frac{\partial v^a}{\partial k}} \right)\left( {\frac{\partial v^n}{\partial k}} \right) 
    304 } \right]-A^{vT}\,{N_n}^2+... 
    305 \end{equation} 
    306  
    307 Note that during a time step, the equation \eqref{Eq_zdftke_e} is resolved before those  
    308 of momentum and density. So, the indice "a" (after) become "n" (now) and the indice "n"  
    309 (now) become "b" (before). 
    310  
    311 Moreover, the vertical shear have to be multiply by the appropriate viscosity for numerical  
    312 stability. Thus, the vertical shear at U-point have to be multiply by the viscosity avmu and  
    313 the vertical shear at V-point have to be multiply by the viscosity avmv. Next, these two  
    314 quantities are averaged to obtain a production term by vertical shear at W-point :  
    315  
    316 \begin{equation} \label{Eq_zdftke_ene2} 
    317 \frac{\partial \bar{e}}{\partial t} = 
    318 \frac{1}{{e_3}^2 }\;\left[ {   
    319     A^{vmu}\left({\frac{\partial u^a}{\partial k}} \right)\left( {\frac{\partial u^n}{\partial k}} \right) 
    320    +A^{vmv}\left({\frac{\partial v^a}{\partial k}} \right)\left( {\frac{\partial v^n}{\partial k}} \right) 
    321                            } \right]-A^{vT}\,{N_n}^2+... 
    322 \end{equation}  
    323  
    324 The TKE equation is resolved before the mixing length, the viscosity and diffusivity. Two tabs  
    325 are then declared : dissl (dissipation length) (Remark : it's not only the dissipation lenght,  
    326 it's the root of the TKE divided by the dissipation lenght) and avmt (viscosity at the points T)  
    327 used for the vertical diffusion of the TKE.  
    328  
    329 This new organization needs also a reorganization of the routine step.F90.  
    330 The bigger change is the estimation of the Brunt-Vais\"{a}l\"{a} 
    331 frequency at "n" instead of "b". Moreover for energetic considerations, the call of tranxt.F90  
    332 is done after the density update. On the contrary, the density is updated with scalars fields  
    333 filtered by the Asselin filter.  
    334  
    335 This new organisation of the routine zdftke force to save five three dimensionnal tabs in  
    336 the restart : avmu, avmv, avt, avmt and dissl are needed to calculate $e_n$. At the end  
    337 of the run (time step = nitend), the alternative is to save only $e_n$ estimated at the  
    338 following time step (nitend+1). The next run using this restart file, the mixing length  
    339 and turbulents coefficients are directly calculated using $e_n$. It is the same thing  
    340 for the intermediate restart. 
    341  
    342  
    343 %%GM for figure of the time scheme: 
    344 \begin{equation}  
    345   \rho^{t+1} \; A^{vT}{N_{t+1}}^2 \; dk  \\ 
    346 \end{equation} 
    347  
    348 %%end GM  
    349  
     311\begin{split} 
     312\frac { (\bar{e})^t - (\bar{e})^{t-\rdt} } {\rdt}  \equiv   
     313\Biggl\{ \Biggr. 
     314  &\overline{ \left( \left(\overline{K_m}^{\,i+1/2}\right)^{t-\rdt} \,\frac{\delta_{k+1/2}[u^{t+\rdt}]}{{e_3u}^{t+\rdt} }  
     315                                                                              \ \frac{\delta_{k+1/2}[u^ t         ]}{{e_3u}^ t          }  \right) }^{\,i} \\ 
     316+&\overline{  \left( \left(\overline{K_m}^{\,j+1/2}\right)^{t-\rdt} \,\frac{\delta_{k+1/2}[v^{t+\rdt}]}{{e_3v}^{t+\rdt} }  
     317                                                                               \ \frac{\delta_{k+1/2}[v^ t         ]}{{e_3v}^ t          }  \right) }^{\,j}  
     318\Biggr. \Biggr\}   \\ 
     319% 
     320- &{K_\rho}^{t-\rdt}\,{(N^2)^t}    \\ 
     321% 
     322+&\frac{1}{{e_3w}^{t+\rdt}}  \;\delta_{k+1/2} \left[   {K_m}^{t-\rdt} \,\frac{\delta_{k}[(\bar{e})^{t+\rdt}]} {{e_3w}^{t+\rdt}}   \right]   \\ 
     323% 
     324- &c_\epsilon \; \left( \frac{\sqrt{\bar {e}}}{l_\epsilon}\right)^{t-\rdt}\,(\bar {e})^{t+\rdt} 
     325\end{split} 
     326\end{equation} 
     327where the last two terms in \eqref{Eq_zdftke_ene} (vertical diffusion and Kolmogorov dissipation)  
     328are time stepped using a backward scheme (see\S\ref{STP_forward_imp}).  
     329Note that the Kolmogorov term has been linearized in time in order to render  
     330the implicit computation possible. The restart of the TKE scheme  
     331requires the storage of $\bar {e}$, $K_m$, $K_\rho$ and $l_\epsilon$ as they all appear in  
     332the right hand side of \eqref{Eq_zdftke_ene}. For the latter, it is in fact  
     333the ratio $\sqrt{\bar{e}}/l_\epsilon$ which is stored.  
     334 
     335% ------------------------------------------------------------------------------------------------------------- 
     336%        GLS Generic Length Scale Scheme  
     337% ------------------------------------------------------------------------------------------------------------- 
     338\subsection{GLS Generic Length Scale (\key{zdfgls})} 
     339\label{ZDF_gls} 
     340 
     341%--------------------------------------------namgls--------------------------------------------------------- 
     342\namdisplay{namgls} 
     343%-------------------------------------------------------------------------------------------------------------- 
     344 
     345The Generic Length Scale (GLS) scheme is a turbulent closure scheme based on  
     346two prognostic equations: one for the turbulent kinetic energy $\bar {e}$, and another  
     347for the generic length scale, $\psi$ \citep{Umlauf_Burchard_JMS03, Umlauf_Burchard_CSR05}.  
     348This later variable is defined as : $\psi = {C_{0\mu}}^{p} \ {\bar{e}}^{m} \ l^{n}$,  
     349where the triplet $(p, m, n)$ value given in Tab.\ref{Tab_GLS} allows to recover  
     350a number of well-known turbulent closures ($k$-$kl$ \citep{Mellor_Yamada_1982},  
     351$k$-$\epsilon$ \citep{Rodi_1987}, $k$-$\omega$ \citep{Wilcox_1988}  
     352among others \citep{Umlauf_Burchard_JMS03,Kantha_Carniel_CSR05}).  
     353The GLS scheme is given by the following set of equations: 
     354\begin{equation} \label{Eq_zdfgls_e} 
     355\frac{\partial \bar{e}}{\partial t} =  
     356\frac{K_m}{\sigma_e e_3 }\;\left[ {\left( \frac{\partial u}{\partial k} \right)^2 
     357                                                   +\left( \frac{\partial v}{\partial k} \right)^2} \right] 
     358-K_\rho \,N^2 
     359+\frac{1}{e_3}\,\frac{\partial}{\partial k} \left[ \frac{K_m}{e_3}\,\frac{\partial \bar{e}}{\partial k} \right] 
     360- \epsilon 
     361\end{equation} 
     362 
     363\begin{equation} \label{Eq_zdfgls_psi} 
     364   \begin{split} 
     365\frac{\partial \psi}{\partial t} =& \frac{\psi}{\bar{e}} \left\{ 
     366\frac{C_1\,K_m}{\sigma_{\psi} {e_3}}\;\left[ {\left( \frac{\partial u}{\partial k} \right)^2 
     367                                                                   +\left( \frac{\partial v}{\partial k} \right)^2} \right] 
     368- C_3 \,K_\rho\,N^2   - C_2 \,\epsilon \,Fw   \right\}             \\ 
     369&+\frac{1}{e_3}  \;\frac{\partial }{\partial k}\left[ {\frac{K_m}{e_3 } 
     370                                  \;\frac{\partial \psi}{\partial k}} \right]\; 
     371   \end{split} 
     372\end{equation} 
     373 
     374\begin{equation} \label{Eq_zdfgls_kz} 
     375   \begin{split} 
     376         K_m    &= C_{\mu} \ \sqrt {\bar{e}} \ l         \\ 
     377         K_\rho &= C_{\mu'}\ \sqrt {\bar{e}} \ l 
     378   \end{split} 
     379\end{equation} 
     380 
     381\begin{equation} \label{Eq_zdfgls_eps} 
     382{\epsilon} = C_{0\mu} \,\frac{\bar {e}^{3/2}}{l} \; 
     383\end{equation} 
     384where $N$ is the local Brunt-Vais\"{a}l\"{a} frequency (see \S\ref{TRA_bn2})  
     385and $\epsilon$ the dissipation rate.  
     386The constants $C_1$, $C_2$, $C_3$, ${\sigma_e}$, ${\sigma_{\psi}}$ and the wall function ($Fw$)  
     387depends of the choice of the turbulence model. Four different turbulent models are pre-defined  
     388(Tab.\ref{Tab_GLS}). They are made available through th \np{gls} namelist parameter.  
     389 
     390%--------------------------------------------------TABLE-------------------------------------------------- 
     391\begin{table}[htbp]  \label{Tab_GLS} 
     392\begin{center} 
     393%\begin{tabular}{cp{70pt}cp{70pt}cp{70pt}cp{70pt}cp{70pt}cp{70pt}c} 
     394\begin{tabular}{ccccc} 
     395                         &   $k-kl$   & $k-\epsilon$ & $k-\omega$ &   generic   \\   
     396%                        & \citep{Mellor_Yamada_1982} &  \citep{Rodi_1987}       & \citep{Wilcox_1988} &                 \\   
     397\hline  \hline  
     398\np{nn\_clo}     & \textbf{0} &   \textbf{1}  &   \textbf{2}   &    \textbf{3}   \\   
     399\hline  
     400$( p , n , m )$          &   ( 0 , 1 , 1 )   & ( 3 , 1.5 , -1 )   & ( -1 , 0.5 , -1 )    &  ( 2 , 1 , -0.67 )  \\ 
     401$\sigma_k$      &    2.44         &     1.              &      2.                &      0.8          \\ 
     402$\sigma_\psi$  &    2.44         &     1.3            &      2.                 &       1.07       \\ 
     403$C_1$              &      0.9         &     1.44          &      0.555          &       1.           \\ 
     404$C_2$              &      0.5         &     1.92          &      0.833          &       1.22       \\ 
     405$C_3$              &      1.           &     1.              &      1.                &       1.           \\ 
     406$F_{wall}$        &      Yes        &       --             &     --                  &      --          \\ 
     407\hline 
     408\hline 
     409\end{tabular} 
     410\caption {Set of predefined GLS parameters, or equivalently predefined turbulence models available with \key{gls} and controlled by the \np{nn\_clos} namelist parameter.} 
     411\end{center} 
     412\end{table} 
     413%-------------------------------------------------------------------------------------------------------------- 
     414 
     415In the Mellor-Yamada model, the negativity of $n$ allows to use a wall function to force 
     416the convergence of the mixing length towards $K\,z_b$ ($K$: Kappa and $z_b$: rugosity length)  
     417value near physical boundaries (logarithmic boundary layer law). $C_{\mu}$ and $C_{\mu'}$  
     418are calculated from stability function proposed by \citet{Galperin_al_JAS88}, or by \citet{Kantha_Clayson_1994}  
     419or one of the two functions suggested by \citet{Canuto_2001}  (\np{nn\_stab\_func} = 0, 1, 2 or 3, resp.}). The value of $C_{0\mu}$ depends of the choice of the stability function. 
     420 
     421The surface and bottom boundary condition on both $\bar{e}$ and $\psi$ can be calculated  
     422thanks to Dirichlet or Neumann condition through \np{nn\_tkebc\_surf} and \np{nn\_tkebc\_bot}, resp.  
     423The wave effect on the mixing could be also being considered \citep{Craig_Banner_1994}. 
     424 
     425The $\psi$ equation is known to fail in stably stratified flows, and for this reason  
     426almost all authors apply a clipping of the length scale as an \textit{ad hoc} remedy.  
     427With this clipping, the maximum permissible length scale is determined by  
     428$l_{max} = c_{lim} \sqrt{2\bar{e}}/ N$. A value of $c_{lim} = 0.53$ is often used  
     429\citep{Galperin_al_JAS88}. \cite{Umlauf_Burchard_CSR05} show that the value of  
     430the clipping factor is of crucial importance for the entrainment depth predicted in  
     431stably stratified situations, and that its value has to be chosen in accordance  
     432with the algebraic model for the turbulent ßuxes. The clipping is only activated  
     433if \np{ln\_length\_lim}=true, and the $c_{lim}$ is set to the \np{clim\_galp} value. 
    350434 
    351435% ------------------------------------------------------------------------------------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.