diff --git a/docs/conf.py b/docs/conf.py index d705e19878..5d84b3c37a 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -159,7 +159,7 @@ def latexPassthru(name, rawtext, text, lineno, inliner, options={}, content=[]): # General information about the project. project = u'MOM6' -copyright = u'2017-2020, MOM6 developers' +copyright = u'2017-2021, MOM6 developers' # The version info for the project you're documenting, acts as replacement for # |version| and |release|, also used in various other places throughout the diff --git a/docs/parameterizations_vertical.rst b/docs/parameterizations_vertical.rst index 4705cf6c48..ff0784b698 100644 --- a/docs/parameterizations_vertical.rst +++ b/docs/parameterizations_vertical.rst @@ -14,9 +14,13 @@ K-profile parameterization (KPP) Energetic Planetary Boundary Layer (ePBL) A energetically constrained boundary layer scheme following :cite:`reichl2018`. Implemented in MOM_energetic_PBL. + :ref:`EPBL` + Bulk mixed layer (BML) A 2-layer bulk mixed layer used in pure-isopycnal model. Implemented in MOM_bulk_mixed_layer. + :ref:`BML` + Interior and bottom-driven mixing --------------------------------- diff --git a/docs/zotero.bib b/docs/zotero.bib index bb400542b8..c0c7ee3bd9 100644 --- a/docs/zotero.bib +++ b/docs/zotero.bib @@ -2684,3 +2684,57 @@ @techreport{griffies2015a pages = {98 pp}, institution = {NOAA GFDL} } + +@inbook{niiler1977, + author = {P. P. Niiler and E. B. Kraus}, + chapter = {One-dimesional models of the upper ocean}, + title = {Modelling and Prediction of the Upper Layers of the Ocean}, + year = {1977}, + editor = {E. B. Kraus}, + publisher = {Pergamon Press} +} + +@article{oberhuber1993, + doi = {10.1175/1520-0485(1993)023<0830:sotacw>2.0.co;2}, + year = 1993, + publisher = {American Meteorological Society}, + volume = {23}, + number = {5}, + pages = {830--845}, + author = {J. M. Oberhuber}, + title = {Simulation of the Atlantic Circulation with a Coupled Sea Ice-Mixed Layer-Isopycnal General Circulation Model. Part {II}: Model Experiment}, + journal = {J. Phys. Oceanography} +} + +@techreport{muller2003, + doi = {10.21236/ada618366}, + year = 2003, + publisher = {Defense Technical Information Center}, + author = {P. Muller}, + institution = {School of Ocean and Earth Science and Technology}, + title = {A{\textquotesingle}ha Huliko{\textquotesingle}a Workshop Series} +} + +@article{wang2003, + doi = {10.1029/2003gl017869}, + year = 2003, + publisher = {American Geophysical Union ({AGU})}, + volume = {30}, + number = {18}, + author = {D. Wang}, + title = {Entrainment laws and a bulk mixed layer model of rotating convection derived from large-eddy simulations}, + journal = {Geophys. Res. Lett.} +} + +@article{kraus1967, + doi = {10.3402/tellusa.v19i1.9753}, + year = 1967, + publisher = {Informa {UK} Limited}, + volume = {19}, + number = {1}, + pages = {98--106}, + author = {E. B. Kraus and J. S. Turner}, + title = {A one-dimensional model of the seasonal thermocline {II}. The general theory and its consequences}, + journal = {Tellus} +} + diff --git a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 index 4a9d428807..137294eda1 100644 --- a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 +++ b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 @@ -155,36 +155,8 @@ module MOM_bulk_mixed_layer contains -!> This subroutine partially steps the bulk mixed layer model. -!! The following processes are executed, in the order listed. -!! 1. Undergo convective adjustment into mixed layer. -!! 2. Apply surface heating and cooling. -!! 3. Starting from the top, entrain whatever fluid the TKE budget -!! permits. Penetrating shortwave radiation is also applied at -!! this point. -!! 4. If there is any unentrained fluid that was formerly in the -!! mixed layer, detrain this fluid into the buffer layer. This -!! is equivalent to the mixed layer detraining to the Monin- -!! Obukhov depth. -!! 5. Divide the fluid in the mixed layer evenly into CS%nkml pieces. -!! 6. Split the buffer layer if appropriate. -!! Layers 1 to nkml are the mixed layer, nkml+1 to nkml+nkbl are the -!! buffer layers. The results of this subroutine are mathematically -!! identical if there are multiple pieces of the mixed layer with -!! the same density or if there is just a single layer. There is no -!! stability limit on the time step. -!! -!! The key parameters for the mixed layer are found in the control structure. -!! These include mstar, nstar, nstar2, pen_SW_frac, pen_SW_scale, and TKE_decay. -!! For the Oberhuber (1993) mixed layer, the values of these are: -!! pen_SW_frac = 0.42, pen_SW_scale = 15.0 m, mstar = 1.25, -!! nstar = 1, TKE_decay = 2.5, conv_decay = 0.5 -!! TKE_decay is 1/kappa in eq. 28 of Oberhuber (1993), while conv_decay is 1/mu. -!! Conv_decay has been eliminated in favor of the well-calibrated form for the -!! efficiency of penetrating convection from Wang (2003). -!! For a traditional Kraus-Turner mixed layer, the values are: -!! pen_SW_frac = 0.0, pen_SW_scale = 0.0 m, mstar = 1.25, -!! nstar = 0.4, TKE_decay = 0.0, conv_decay = 0.0 +!> This subroutine partially steps the bulk mixed layer model. +!! See \ref BML for more details. subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, CS, & optics, Hml, aggregate_FW_forcing, dt_diag, last_call) type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure. @@ -3708,16 +3680,15 @@ end function EF4 !! !! This file contains the subroutine (bulkmixedlayer) that !! implements a Kraus-Turner-like bulk mixed layer, based on the work -!! of various people, as described in the review paper by Niiler and -!! Kraus (1979), with particular attention to the form proposed by -!! Oberhuber (JPO, 1993, 808-829), with an extension to a refied bulk -!! mixed layer as described in Hallberg (Aha Huliko'a, 2003). The -!! physical processes portrayed in this subroutine include convective -!! adjustment and mixed layer entrainment and detrainment. -!! Penetrating shortwave radiation and an exponential decay of TKE -!! fluxes are also supported by this subroutine. Several constants +!! of various people, as described in the review paper by \cite Niiler1977, +!! with particular attention to the form proposed by \cite Oberhuber1993, +!! with an extension to a refined bulk mixed layer as described in +!! Hallberg (\cite muller2003). The physical processes portrayed in +!! this subroutine include convective adjustment and mixed layer entrainment +!! and detrainment. Penetrating shortwave radiation and an exponential decay +!! of TKE fluxes are also supported by this subroutine. Several constants !! can alternately be set to give a traditional Kraus-Turner mixed !! layer scheme, although that is not the preferred option. The -!! physical processes and arguments are described in detail below. +!! physical processes and arguments are described in detail in \ref BML. end module MOM_bulk_mixed_layer diff --git a/src/parameterizations/vertical/MOM_energetic_PBL.F90 b/src/parameterizations/vertical/MOM_energetic_PBL.F90 index 5a9e67bfd9..6920b8dd22 100644 --- a/src/parameterizations/vertical/MOM_energetic_PBL.F90 +++ b/src/parameterizations/vertical/MOM_energetic_PBL.F90 @@ -1030,7 +1030,7 @@ subroutine ePBL_column(h, u, v, T0, S0, dSV_dT, dSV_dS, TKE_forcing, B_flux, abs dt_h = (GV%Z_to_H**2*dt) / max(0.5*(h(k-1)+h(k)), 1e-15*h_sum) ! This tests whether the layers above and below this interface are in - ! a convetively stable configuration, without considering any effects of + ! a convectively stable configuration, without considering any effects of ! mixing at higher interfaces. It is an approximation to the more ! complete test dPEc_dKd_Kd0 >= 0.0, that would include the effects of ! mixing across interface K-1. The dT_to_dColHt here are effectively @@ -2079,7 +2079,7 @@ subroutine energetic_PBL_init(Time, G, GV, US, param_file, diag, CS) call log_param(param_file, mdl, "EPBL_MSTAR_SCHEME", tmpstr, & "EPBL_MSTAR_SCHEME selects the method for setting mstar. Valid values are: \n"//& "\t CONSTANT - Use a fixed mstar given by MSTAR \n"//& - "\t OM4 - Use L_Ekman/L_Obukhov in the sabilizing limit, as in OM4 \n"//& + "\t OM4 - Use L_Ekman/L_Obukhov in the stabilizing limit, as in OM4 \n"//& "\t REICHL_H18 - Use the scheme documented in Reichl & Hallberg, 2018.", & default=CONSTANT_STRING) tmpstr = uppercase(tmpstr) @@ -2468,10 +2468,10 @@ end subroutine energetic_PBL_end !! simple enough that it requires only a single vertical pass to !! determine the diffusivity. The development of bulk mixed layer !! models stems from the work of various people, as described in the -!! review paper by Niiler and Kraus (1979). The work here draws in -!! with particular on the form for TKE decay proposed by Oberhuber -!! (JPO, 1993, 808-829), with an extension to a refined bulk mixed -!! layer as described in Hallberg (Aha Huliko'a, 2003). The physical +!! review paper by \cite niiler1977. The work here draws in +!! with particular on the form for TKE decay proposed by +!! \cite oberhuber1993, with an extension to a refined bulk mixed +!! layer as described in Hallberg (\cite muller2003). The physical !! processes portrayed in this subroutine include convectively driven !! mixing and mechanically driven mixing. Unlike boundary-layer !! mixing, stratified shear mixing is not a one-directional turbulent diff --git a/src/parameterizations/vertical/_BML.dox b/src/parameterizations/vertical/_BML.dox new file mode 100644 index 0000000000..2786a26851 --- /dev/null +++ b/src/parameterizations/vertical/_BML.dox @@ -0,0 +1,49 @@ +/*! \page BML Bulk Surface Mixed Layer + +This bulk surface mixed layer scheme was designed to be used with a +purely isopycnal model. Following \cite niiler1977, \cite oberhuber1993, +and Hallberg (\cite muller2003) the TKE budget is used to construct a +time-evolving homogeneous mixed layer. A buffer layer sits between +the mixed layer and the interior ocean to mediate between the two. + + The following processes are executed, in the order listed. + +\li 1. Undergo convective adjustment into mixed layer. +\li 2. Apply surface heating and cooling. +\li 3. Starting from the top, entrain whatever fluid the TKE budget + permits. Penetrating shortwave radiation is also applied at + this point. +\li 4. If there is any unentrained fluid that was formerly in the + mixed layer, detrain this fluid into the buffer layer. This + is equivalent to the mixed layer detraining to the Monin- + Obukhov depth. +\li 5. Divide the fluid in the mixed layer evenly into CS\%nkml pieces. +\li 6. Split the buffer layer if appropriate. + +Layers 1 to nkml are the mixed layer, nkml+1 to nkml+nkbl are the +buffer layers. The results of this subroutine are mathematically +identical if there are multiple pieces of the mixed layer with +the same density or if there is just a single layer. There is no +stability limit on the time step. + +The key parameters for the mixed layer are found in the control structure. +These include mstar, nstar, nstar2, pen\_SW\_frac, pen\_SW\_scale, and TKE\_decay. +For the \cite oberhuber1993 and \cite kraus1967 mixed layers, the values of these are: + + + +
Model variables used in the bulk mixed layer
Symbol Value in Oberhuber (1993) Value in Kraus-Turner (1967) +
pen\_SW\_frac 0.42 0.0 +
pen\_SW\_scale 15.0 m 0.0 m +
mstar 1.25 1.25 +
nstar 1 0.4 +
TKE\_decay 2.5 0.0 +
conv\_decay 0.5 0.0 +
+ +TKE\_decay is \f$1/\kappa\f$ in eq. 28 of \cite oberhuber1993, while +conv\_decay is \f$1/\mu\f$. Conv\_decay has been eliminated in favor of +the well-calibrated form for the efficiency of penetrating convection +from \cite wang2003. + +*/ diff --git a/src/parameterizations/vertical/_EPBL.dox b/src/parameterizations/vertical/_EPBL.dox new file mode 100644 index 0000000000..d531c9ad9a --- /dev/null +++ b/src/parameterizations/vertical/_EPBL.dox @@ -0,0 +1,254 @@ +/*! \page EPBL Energetically-constrained Planetary Boundary Layer + +We here describe a scheme for modeling the ocean surface boundary layer +(OSBL) suitable for use in global climate models. It builds on the ideas in +\ref BML, bringing in some of the ideas from \ref subsection_kappa_shear, to +make an energetically consistent boundary layer suitable for use with +a generalized vertical coordinate. Unlike in \ref BML, variables are +allowed to have vertical structure within the boundary layer. The downward +turbulent flux of buoyant water by OSBL turbulence converts mechanical +energy into potential energy as it mixes with less buoyant water at the +base of the OSBL. As described in \cite reichl2018, we focus on OSBL +parameterizations that constrain this integrated potential energy +conversion due to turbulent mixing. + +The leading-order mean OSBL equation for arbitrary scalar \f$\phi\f$ is: + +\f[ + \frac{\partial \overline{\phi}}{\partial t} = - \frac{\partial}{\partial z} + \overline{w^\prime \phi^\prime} + \nu_\phi \frac{\partial^2 \overline{\phi}}{\partial z^2} +\f] + +where the symbols are as follows: + + + +
Symbols used in TKE equation
Symbol Meaning +
\f$u_i\f$ horizontal components of the velocity +
\f$\phi\f$ arbitrary scalar (tracer) quantity +
\f$w\f$ vertical component of the velocity +
\f$\overline{w}\f$ ensemble average \f$w\f$ +
\f$w^\prime\f$ fluctuations from \f$\overline{w}\f$ +
\f$k\f$ turbulent kinetic energy (TKE) +
\f$K_M\f$ turbulent mixing coefficient for momentum +
\f$K_\phi\f$ turbulent mixing coefficient for \f$\phi\f$ +
\f$\sigma_k\f$ turbulent Schmidt number +
\f$b\f$ buoyancy +
\f$\epsilon\f$ buoyancy turbulent dissipation rate +
+ +This equation describes the evolution of mean quantity \f$\overline{\phi}\f$ +due to vertical processes, including the often negligible molecular +mixing. We would like to parameterize the vertical mixing since we won't be +resolving all the relevant time and space scales. + +We use the Boussinesq hypothesis for turbulence closure. This approximates +the Reynolds stress terms using an eddy viscosity (eddy diffusivity for +turbulent scalar fluxes): + +\f[ + \overline{u_i^\prime w^\prime} = - K_M \frac{\partial \overline{u_i}}{\partial z} , +\f] + +Similarly, the eddy diffusivity is used to parameterize turbulent scalar fluxes as: + +\f[ + \overline{\phi^\prime w^\prime} = - K_\phi \frac{\partial \overline{\phi}}{\partial z} , +\f] + +The parameters needed to close the system of equations are then reduced to the turbulent +mixing coefficients, \f$K_\phi\f$ and \f$K_M\f$. + +We start with an equation for the turbulent kinetic energy (TKE): + +\f[ + \frac{\partial k}{\partial t} = \frac{\partial}{\partial z} \left( \frac{K_M}{\sigma_k} + \frac{\partial k}{\partial z} \right) - \overline{u_i^\prime w^\prime} \frac{\partial \overline{u_i}} + {\partial z} + \overline{w^\prime b^\prime} - \epsilon +\f] + +Terms in this equation represent TKE storage (LHS), TKE flux convergence, +shear production, buoyancy production, and dissipation. + +\section section_WMBL Well-mixed Boundary Layers (WMBL) + +Assuming steady state and other parameterizations, integrating vertically +over the surface boundary layer, \cite reichl2018 obtains the form: + +\f[ + \frac{1}{2} H_{bl} w_e \Delta b = m_\ast u_\ast^3 - n_\ast \frac{H_{bl}}{2} + B(H_{bl}) , +\f] + +with the following variables: + + + +
Symbols used in integrated TKE equation
Symbol Meaning +
\f$H_{bl}\f$ boundary layer thickness +
\f$w_e\f$ entrainment velocity +
\f$\Delta b\f$ change in buoyancy at base of mixed layer +
\f$m_\ast\f$ sum of mechanical coefficients +
\f$u_\ast\f$ friction velocity (\f$u_\ast = (|\tau| / \rho_0)^{1/2}\f$) +
\f$\tau\f$ wind stress +
\f$n_\ast\f$ convective proportionality coefficient +
1 for stabilizing surface buoyancy flux, less otherwise +
\f$B(H_{bl})\f$ surface buoyancy flux +
+ +\section section_ePBL Energetics-based Planetary Boundary Layer + +Once again, the goal is to formulate a surface mixing scheme to find the +turbulent eddy diffusivity (and viscosity) in a way that is suitable for use +in a global climate model, using long timesteps and large grid spacing. +After evaluating a well-mixed boundary layer (WMBL), the shear mixing of +\cite jackson2008 (JHL, \ref subsection_kappa_shear), as well as a more complete +boundary layer scheme, it was decided to combine a number of these ideas +into a new scheme: + +\f[ + K(z) = F_x(K_{ePBL}(z), K_{JHL}(z), K_n(z)) +\f] + +where \f$F_x\f$ is some unknown function of a new \f$K_{ePBL}\f$, +\f$K_{JHL}\f$, the diffusivity due to shear as determined by +\cite jackson2008, and \f$K_n\f$, the diffusivity from other ideas. +We start by specifying the form of \f$K_{ePBL}\f$ as being: + +\f[ + K_{ePBL}(z) = C_K w_t l , +\f] + +where \f$w_t\f$ is a turbulent velocity scale, \f$C_K\f$ is a coefficient, and +\f$l\f$ is a length scale. + +\subsection subsection_lengthscale Turbulent length scale + +We propose a form for the length scale as follows: + +\f[ + l = (z_0 + |z|) \times \max \left[ \frac{l_b}{H_{bl}} , \left( + \frac{H_{bl} - |z|}{H_{bl}} \right)^\gamma \, \right] , +\f] + +where we have the following variables: + + + +
Symbols used in ePBL length scale
Symbol Meaning +
\f$H_{bl}\f$ boundary layer thickness +
\f$z_0\f$ roughness length +
\f$\gamma\f$ coefficient, 2 is as in KPP, \cite large1994 +
\f$l_b\f$ bottom length scale +
+ +\subsection subsection_velocityscale Turbulent velocity scale + +We do not predict the TKE prognostically and therefore approximate the vertical TKE +profile to estimate \f$w_t\f$. An estimate for the mechanical contribution to the velocity +scale follows the standard two-equation approach. In one and two-equation second-order +\f$K\f$ parameterizations the boundary condition for the TKE is typically employed as a +flux boundary condition. + +\f[ + K \left. \frac{\partial k}{\partial z} \right|_{z=0} = c_\mu^0 u_\ast^3 . +\f] + +The profile of \f$k\f$ decays in the vertical from \f$k \propto (c_\mu^0)^{2/3} +u_\ast^2\f$ toward the base of the OSBL. Here we assume a similar relationship to estimate +the mechanical contribution to the TKE profile. The value of \f$w_t\f$ due to mechanical +sources, \f$v_\ast\f$, is estimate as \f$v_\ast (z=0) \propto (c_\mu^0)^{1/3} u_\ast\f$ at +the surface. Since we only parameterize OSBL turbulent mixing due to surface forcing, the +value of the velocity scale is assumed to decay moving away from the surface. For +simplicity we employ a linear decay in depth: + +\f[ + v_\ast (z) = (c_\mu^0)^{1/3} u_\ast \left( 1 - a \cdot \min \left[ 1, + \frac{|z|}{H_{bl}} \right] \right) , +\f] + +where \f$1 > a > 0\f$ has the effect of making \f$v_\ast(z=H_{bl}) > 0\f$. +Making the constant coefficient \f$a\f$ close to one has the effect of reducing the mixing +rate near the base of the boundary layer, thus producing a more diffuse entrainment +region. Making \f$a\f$ close to zero has the effect of increasing the mixing at the base +of the boundary layer, producing a more 'step-like' entrainment region. + +An estimate for the buoyancy contribution is found utilizing the convective velocity +scale: + +\f[ + w_\ast (z) = C_{w_\ast} \left( \int_z^0 \overline{w^\prime b^\prime} dz \right)^{1/3} , +\f] + +where \f$C_{w_\ast}\f$ is a non-dimensional empirical coefficient. Convection in one and +two-equation closure causes a TKE profile that peaks below the surface. The quantity +\f$\overline{w^\prime b^\prime}\f$ is solved for in ePBL as \f$KN^2\f$. + +These choices for the convective and mechanical components of the velocity scale in the +OSBL are then added together to get an estimate for the total turbulent velocity scale: + +\f[ + w_t (z) = w_\ast (z) + v_\ast (z) . +\f] + +The value of \f$a\f$ is arbitrarily chosen to be 0.95 here. + +\subsection subsection_ePBL_summary Summarizing the ePBL implementation + +The ePBL mixing coefficient is found by multiplying a velocity scale +(\ref subsection_velocityscale) by a length scale (\ref subsection_lengthscale). The +precise value of the coefficient \f$C_K\f$ used does not significantly alter the +prescribed potential energy change constraint. A reasonable value is \f$C_K \approx 0.55\f$ to +be consistent with other approaches (e.g. \cite umlauf2005). + +The boundary layer thickness (\f$H_{bl}\f$) within ePBL is based on +the depth where the energy requirement for turbulent mixing of density +exceeds the available energy (\ref section_WMBL). \f$H_{bl}\f$ is +determined by the energetic constraint imposed using the value of +\f$m_\ast\f$ and \f$n_\ast\f$. An iterative solver is required because +\f$m_\ast\f$ and the mixing length are dependent on \f$H_{bl}\f$. + +We use a constant value for convectively driven TKE of \f$n_\ast = 0.066\f$. The +parameterizations for \f$m_\ast\f$ are formulated specifically for the regimes where +\f$K_{JHL}\f$ is sensitive to model numerics \f$(|f| \Delta t \approx +1)\f$ (\cite reichl2018). + +\subsection subsection_ePBL_JHL Combining ePBL and JHL mixing coefficients + +We now address the combination of the ePBL mixing coefficient and the JHL mixing +coefficient. The function \f$F_x\f$ above cannot be the linear sum of \f$K_{ePBL}\f$ and +\f$K_{JHL}\f$. One reason this sum is not valid is because the JHL mixing coefficient is +determined by resolved current shear, including that driven by the surface wind. The +wind-driven current is also included in the ePBL mixing coefficient formulation. An +alternative approach is therefore needed to avoid double counting. + +\f$K_{ePBL}\f$ is not used at the equator as scalings are only investigated when \f$|f| > +0\f$. The solution we employ is to use the maximum mixing coefficient of the two +contributions, + +\f[ + K (z) = \max (K_{ePBL} (z), K_{JHL} (z)), +\f] + +where \f$m_\ast\f$ (and hence \f$K_{ePBL}\f$) is constrained to be small as \f$|f| +\rightarrow 0\f$. This form uses the JHL mixing coefficient when the ePBL coefficient is +small. + +This approach is reasonable when the wind-driven mixing dominates, since both JHL and ePBL +give a similar solution when deployed optimally. One weakness of this approach is the +tropical region, where the shear-driven ePBL \f$m_\ast\f$ coefficient is not formulated. +The JHL parameterization is skillful to simulate this mixing, but does not include the +contribution of convection. The convective portion of \f$K_{ePBL}\f$ should be combined +with \f$K_{JHL}\f$ in the equatorial region when shear and convection occur together. +Future research is warranted. + +Finally, one should note that the mixing coefficient here (\f$K\f$) is used for both +diffusivity and viscosity, implying a turbulent Prandtl number of 1.0. + +\subsection subsection_Langmuir Langmuir circulation + +While only briefly alluded to in \cite reichl2018, the MOM6 code implementing ePBL does +support the option to add a Langmuir parameterization. There are in fact two options, both +adjusting \f$m_\ast\f$. + +*/