From 6dc3fdf024bc598274ced97843df6ab326d1afc4 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Tue, 28 Jul 2026 10:59:26 +0900 Subject: [PATCH 1/7] Describe the formulation of a slab ocean model --- doc/description/ocean.tex | 57 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 57 insertions(+) create mode 100644 doc/description/ocean.tex diff --git a/doc/description/ocean.tex b/doc/description/ocean.tex new file mode 100644 index 00000000..a30d6c09 --- /dev/null +++ b/doc/description/ocean.tex @@ -0,0 +1,57 @@ +\section{Ocean model} +{\bf \Large +\begin{tabular}{ccc} +\hline + Corresponding author & : & Yuta Kawai\\ +\hline +\end{tabular} +} + +\subsection{A single-layered slab ocean model} +%% +This model only treats a thermodynamics in a single-layer ocean, +which calculate the evolution of sea temperature. + +The governing equation is described using the internal energy of ocean as +%% +\begin{align} + \dfrac{\partial (C_o e_{{\rm in},o})}{\partial t} =& + R_{\rm SW,dn} - R_{\rm SW,up} + R_{\rm LW,dn} - R_{\rm LW,up} - F_{\rm SH} - F_{\rm LH} \nonumber \\ + +& e_{\rm in,precip} + Q_{\rm ext}, \label{eq:goevern_eq_slab_ocean} +\end{align} +%% +where $C_o$ is the heat capacity of the slab ocean +and +$e_{{\rm in},o}$ is the internal energy. +For the radiative fluxes at the surface $R_{*,*}$, +the lower subscripts SW, LW, dn, and up denotes shortwave, longwave, downward, and upward, respectively. +$F_{\rm SH}$ and $F_{\rm LH}$ are the sensible heat and latent heat fluxes, respectively. +$e_{\rm precip}$ is the internal energy flux due to the precipitation. +Finally, $Q_{\rm ext}$ is the heat flux with sources of external heating. + +$C_o$ is calculated as +\begin{align} + C_o = \rho_o c_{p,o} h_o, +\end{align} +where $\rho_o$ and $c_{p,o}$ are the density and specific heat of ocean, respectively, +and $h_o$ is the depth of the slab ocean. + + +$R_{\rm SW,dn}$ and $R_{\rm LW,dn}$ are provided by the radiation schemes described in Sec.\,\ref{sec:radiation}. +On the other hand, using the surface albedo $\alpha_o$ and the sea surface temperature $T_s$ diagnosed by the ocean model, +$R_{\rm SW,up}$ is calculated as +%% +\begin{align} + R_{\rm SW,up} = \alpha_o R_{\rm SW,dn}, +\end{align} +%% +and $R_{\rm LW,up}$ is calculated based on the Stefan-Boltzmann law as +\begin{align} + R_{\rm LW,up} = \sigma T_s^4, +\end{align} +where $\sigma$ is the Stefan-Boltzmann constant. +$F_{\rm SH}$ and $F_{\rm LH}$ are provided by the surface flux schemes, +and the simple bulk-flux scheme described in Sec.\,\ref{ssec:simple_sfcflux_bulk} is used in the current implementation. + +In the temporal discretization in Eq.\,\eqref{eq:goevern_eq_slab_ocean}, the forward Euler method is used +because the temporal variation of $e_{{\rm in},o}$ is relatively slow due to the large heat capacity of the slab ocean. \ No newline at end of file From 27677431f469639ff9ed5695a4b2756966c139c2 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Wed, 29 Jul 2026 15:49:39 +0900 Subject: [PATCH 2/7] Add radiation section with detailed formulations --- doc/description/ocean.tex | 58 ++++--- doc/description/radiation.tex | 300 ++++++++++++++++++++++++++++++++++ doc/reference.bib | 82 ++++++++++ 3 files changed, 415 insertions(+), 25 deletions(-) create mode 100644 doc/description/radiation.tex diff --git a/doc/description/ocean.tex b/doc/description/ocean.tex index a30d6c09..b3c59c51 100644 --- a/doc/description/ocean.tex +++ b/doc/description/ocean.tex @@ -7,10 +7,13 @@ \section{Ocean model} \end{tabular} } -\subsection{A single-layered slab ocean model} +\subsection{A single-layer slab ocean model} %% -This model only treats a thermodynamics in a single-layer ocean, -which calculate the evolution of sea temperature. +The current version of SCALE-DG provides a simple single-layer slab ocean model. +The model represents only the thermodynamics of the upper ocean and predicts the evolution of the sea surface temperature (SST). +The slab ocean is treated independently at each DG node, without horizontal ocean circulation or vertical mixing below the mixed layer. +This model is intended for idealized experiments and for studies in which ocean dynamics play only a secondary role. +This model only treats a thermodynamics in a single-layer ocean, which calculate the evolution of sea temperature. The governing equation is described using the internal energy of ocean as %% @@ -20,38 +23,43 @@ \subsection{A single-layered slab ocean model} +& e_{\rm in,precip} + Q_{\rm ext}, \label{eq:goevern_eq_slab_ocean} \end{align} %% -where $C_o$ is the heat capacity of the slab ocean -and -$e_{{\rm in},o}$ is the internal energy. -For the radiative fluxes at the surface $R_{*,*}$, -the lower subscripts SW, LW, dn, and up denotes shortwave, longwave, downward, and upward, respectively. -$F_{\rm SH}$ and $F_{\rm LH}$ are the sensible heat and latent heat fluxes, respectively. -$e_{\rm precip}$ is the internal energy flux due to the precipitation. -Finally, $Q_{\rm ext}$ is the heat flux with sources of external heating. +where $C_o$ is the heat capacity of the slab ocean and $e_{{\rm in},o}$ is its internal energy. +For the surface radiation fluxes $R_{*,*}$, the subscripts SW, LW, dn, and up denote shortwave radiation, +longwave radiation, downward, and upward fluxes, respectively. +$F_{\rm SH}$ and $F_{\rm LH}$ are the sensible and latent heat fluxes exchanged between the ocean and the atmosphere. +$e_{\rm in,precip}$ represents the internal energy supplied by precipitation. +Finally, $Q_{\rm ext}$ denotes an externally prescribed heat source, +such as an idealized ocean heat transport or a prescribed Q-flux. -$C_o$ is calculated as +The heat capacity of the slab ocean is assumed to be constant and is calculated as +%% \begin{align} C_o = \rho_o c_{p,o} h_o, \end{align} -where $\rho_o$ and $c_{p,o}$ are the density and specific heat of ocean, respectively, -and $h_o$ is the depth of the slab ocean. - +%% +where $\rho_o$ and $c_{p,o}$ are the density and specific heat capacity of seawater, +respectively, and $h_o$ is the prescribed depth of the slab ocean. +Increasing $h_o$ increases the thermal inertia of the ocean +and therefore lengthens the response time of the SST to atmospheric forcing. -$R_{\rm SW,dn}$ and $R_{\rm LW,dn}$ are provided by the radiation schemes described in Sec.\,\ref{sec:radiation}. -On the other hand, using the surface albedo $\alpha_o$ and the sea surface temperature $T_s$ diagnosed by the ocean model, -$R_{\rm SW,up}$ is calculated as +The downward radiative fluxes, $R_{\rm SW,dn}$ and $R_{\rm LW,dn}$, are provided by the radiation schemes described in Sec.\,\ref{sec:radiation}. +The upward shortwave radiation is calculated using the prescribed ocean surface albedo $\alpha_o$ as %% \begin{align} R_{\rm SW,up} = \alpha_o R_{\rm SW,dn}, \end{align} %% -and $R_{\rm LW,up}$ is calculated based on the Stefan-Boltzmann law as +whereas the upward longwave radiation is diagnosed from the SST according to the Stefan--Boltzmann law, +%% \begin{align} - R_{\rm LW,up} = \sigma T_s^4, + R_{\rm LW,up} = \sigma_{\rm SB} T_s^4, \end{align} -where $\sigma$ is the Stefan-Boltzmann constant. -$F_{\rm SH}$ and $F_{\rm LH}$ are provided by the surface flux schemes, -and the simple bulk-flux scheme described in Sec.\,\ref{ssec:simple_sfcflux_bulk} is used in the current implementation. +%% +where $\sigma_{\rm SB}$ is the Stefan--Boltzmann constant and +$T_s$ is the sea surface temperature diagnosed from the internal energy of the slab ocean. +The sensible and latent heat fluxes are supplied by the surface-flux scheme. +In the current implementation, the simple bulk-flux scheme described in Sec.\,\ref{ssec:simple_sfcflux_bulk} is used. -In the temporal discretization in Eq.\,\eqref{eq:goevern_eq_slab_ocean}, the forward Euler method is used -because the temporal variation of $e_{{\rm in},o}$ is relatively slow due to the large heat capacity of the slab ocean. \ No newline at end of file +Because the slab ocean possesses a much larger heat capacity than the atmosphere, +its characteristic evolution time is much longer than the atmospheric time step. +Therefore, the governing equation (Eq.\,\eqref{eq:goevern_eq_slab_ocean}) is integrated using the forward Euler method. diff --git a/doc/description/radiation.tex b/doc/description/radiation.tex new file mode 100644 index 00000000..a47e4f1d --- /dev/null +++ b/doc/description/radiation.tex @@ -0,0 +1,300 @@ +\section{Radiation}\label{sec:radiation} + +{\bf \Large +\begin{tabular}{ccc} +\hline + Corresponding author & : & Yuta Kawai\\ +\hline +\end{tabular} +} +\\ + +Radiative transfer describes the propagation of electromagnetic radiation through the atmosphere and determines the radiative heating and cooling resulting from the absorption, emission, and, in general, scattering of radiation. + +\subsection{Basic formulation of radiative transfer} +%% +The following formulation follows the standard theory of radiative transfer~\citep[e.g.,][]{Liou2002}. +In general, the specific intensity depends on the position, propagation direction, wavelength, and time. +For a given wavelength and propagation direction, the radiative transfer equation is written as +%% +\begin{equation} + \frac{dI_\lambda(\bm{\sigma})}{ds} = + -\chi_\lambda I_\lambda(\bm{\sigma}) + +\kappa_\lambda B_\lambda + +\sigma_{s,\lambda}\mathcal{S}_\lambda(\bm{\sigma}), +\label{eq:radiation_transfer_equation_general_form} +\end{equation} +%% +where $I_\lambda(\bm{\sigma})$ is the specific intensity at wavelength $\lambda$ in the propagation direction $\bm{\sigma}$, +$s$ is the distance along the radiation path, +$\kappa_\lambda$ and $\sigma_{s,\lambda}$ are the absorption and scattering coefficients, respectively, +$\chi_\lambda=\kappa_\lambda+\sigma_{s,\lambda}$ is the extinction coefficient, +and $B_\lambda$ is the Planck function. $\mathcal{S}_\lambda$ is the scattering source function defined as +%% +\begin{equation} + \mathcal{S}_\lambda(\bm{\Omega}) = \frac{1}{4\pi}\int_{4\pi} I_\lambda(\bm{\Omega}') P_\lambda(\bm{\Omega},\bm{\Omega}') d\bm{\Omega}', +\end{equation} +%% +where $P_\lambda(\bm{\Omega},\bm{\Omega}')$ is the scattering phase function that describes the angular distribution of scattered radiation. +Introducing the optical thickness $\tau_\lambda$ defined as +%% +\begin{equation} + d\tau_\lambda = \chi_\lambda ds, +\end{equation} +%% +Eq.\,\eqref{eq:radiation_transfer_equation_general_form} becomes +%% +\begin{equation} + \frac{dI_\lambda}{d\tau_\lambda} = - I_\lambda + S_\lambda. +\label{eq:radiation_transfer_equation_general_form_tau} +\end{equation} +Here, $S_\lambda$ is the source function defined as +%% +\begin{align} + S_\lambda &= \frac{\kappa_\lambda B_\lambda + \sigma_{s,\lambda} J_\lambda}{\chi_\lambda} \\ + &= (1-\omega_\lambda) B_\lambda + \omega_\lambda J_\lambda, +\end{align} +%% +where $\omega_\lambda = \sigma_{s,\lambda}/\chi_\lambda$ is the single scattering albedo. + + +In atmospheric simulations, the radiative transfer is typically solved independently in each atmospheric column under the plane-parallel approximation, +whereby horizontal radiative transport is neglected. +We introduce the optical thickness along the vertical direction $z$ as +%% +\begin{equation} + d\tau_{\lambda,z}= -\chi_\lambda dz= - \chi_\lambda \mu\, ds = -\mu\, d\tau_\lambda, +\end{equation} +%% +where $\mu$ is the cosine of the zenith angle $\theta$. +Thus, we can rewrite Eq.\,\eqref{eq:radiation_transfer_equation_general_form_tau} as +%% +\begin{equation} + \mu \frac{dI_\lambda}{d\tau_{\lambda,z}} = I_\lambda - S_\lambda. +\label{eq:radiation_transfer_equation_general_form_plane_parallel} +\end{equation} +The equation for the radiative fluxes is obtained by integrating Eq.\,\eqref{eq:radiation_transfer_equation_general_form_plane_parallel} over the semi-sphere as +%% +\begin{align} + \dfrac{dF_\lambda^\uparrow}{d\tau_{\lambda,z}} =& + 4\pi J_\lambda^\uparrow + - 2\pi (1-\omega_\lambda) B_\lambda + - \dfrac{\omega_\lambda}{2}\int_0^1 d\mu \int_{4\pi} P_\lambda(\bm{\Omega},\bm{\Omega}') I_\lambda(\bm{\Omega}') d\Omega', \\ + \dfrac{dF_\lambda^\downarrow}{d\tau_{\lambda,z}} =& + - 4\pi J_\lambda^\uparrow + + 2\pi (1-\omega_\lambda) B_\lambda + + \dfrac{\omega_\lambda}{2}\int_{-1}^0 d\mu \int_{4\pi} P_\lambda(\bm{\Omega},\bm{\Omega}') I_\lambda(\bm{\Omega}') d\Omega'. \label{eq:rte_general_form_plane_parallel_flux} +\end{align} +%% +Here, superscripts $\uparrow$ and $\downarrow$ indicate the quantities integrated over the upward and downward hemispheres, respectively, +$J_\lambda$ is the mean intensity defined as +\begin{align} + J_\lambda = \frac{1}{4\pi}\int_{4\pi} I_\lambda(\bm{\Omega}) d\Omega, +\end{align} +%% +and the upward and downward fluxes are defined as +%% +\begin{align} +F_\lambda^\uparrow =2\pi \int_0^1 I_\lambda(\mu)\mu\,d\mu,\;\; +F_\lambda^\downarrow =2\pi \int_{-1}^0 I_\lambda(\mu)|\mu|\,d\mu, +\end{align} +%% +where both fluxes are defined as positive quantities. + + + +The hemispherically integrated equations are not closed because they still contain the angular moments of the radiation field and the scattering integral. +In the two-stream approximation~\citep{Toon1989}, the angular dependence within each hemisphere is represented by a single effective propagation direction. +The resulting angular integration is approximated using a diffusivity factor $D_\lambda$, +which accounts for the increased optical path length of obliquely propagating radiation. +Then, the two-stream equations can therefore be written schematically as +\begin{align} +\frac{dF_\lambda^\uparrow}{d\tau_{\lambda,z}} +&= D_\lambda^\uparrow \left( F_\lambda^\uparrow -\pi S_\lambda^\uparrow \right), \\ +\frac{dF_\lambda^\downarrow}{d\tau_{\lambda,z}} +&= +-D_\lambda^\downarrow \left( F_\lambda^\downarrow -\pi S_\lambda^\downarrow \right), +\label{eq:rte_twostream_approx} +\end{align} +%% +where $D_\lambda^\uparrow$ and $D_\lambda^\downarrow$ the diffusivity factors determined by the adopted angular quadrature, +$S_\lambda^\uparrow$ and $S_\lambda^\downarrow$ denote the corresponding hemispherically averaged source functions. +Note that the source term $\pi S_\lambda$ is not obtained by directly evaluating $2\pi\int_0^1 S_\lambda\,d\mu$. +Rather, it appears as part of the angular closure introduced in the two-stream approximation. + + +In the current version of SCALE-DG, +the spectral dependence is treated separately in each radiation scheme by dividing the spectrum into one or more bands. +Furthermore, the scattering is neglected in all radiation schemes. +Consequently, +$\sigma_{s,\lambda} = 0$ and $\chi_\lambda = \kappa_\lambda$, and the source function is simplified as $S_\lambda = B_\lambda$. +The radiation transfer equation in Eq.\,\eqref{eq:rte_twostream_approx} using the two-stream approximation therefore becomes +%% +\begin{align} + \frac{dF_\lambda^\uparrow}{d\tau_{\lambda,z}} = D^\uparrow_\lambda \left( F_\lambda^\uparrow - \pi B_\lambda \right), \\ + \frac{dF_\lambda^\downarrow}{d\tau_{\lambda,z}} = -D^\downarrow_\lambda \left( F_\lambda^\downarrow - \pi B_\lambda \right). +\label{eq:rte_general_form_plane_parallel_flux_no_scattering} +\end{align} + + +The following sections describe the radiation schemes currently available in SCALE-DG. +Although the parameterizations of optical properties differ among the schemes, +all are based on Eq.\,\eqref{eq:rte_general_form_plane_parallel_flux_no_scattering}. + + +\subsection{Gray radiation scheme}\label{sec:gray_radiation} +%% +As the simplest radiation model incorporating moisture feedback, a gray radiation model is provided. +This scheme is gray in infrared so that a single optical thickness is defined for the entire longwave spectrum, +which includes a parameterization of long-wave absorption by CO$_2$. + + +Following Eq.\,(5) in \cite{Vallis2018}, +the vertical gradient of longwave optical thickness is parameterized as +%% +\begin{equation} +\frac{\partial \tau_{\mathrm{LW}}}{\partial \sigma} += +a_{\mathrm{LW}}\mu + b_{\mathrm{LW}} q_v + c_{\mathrm{LW}} \ln\!\left(\frac{C_{\mathrm{CO_2}}}{360}\right), +\label{eq:gray_radiation_tau} +\end{equation} +%% +where $\sigma=p/P_0$, $q_v$ is the specific humidity, +$C_{\mathrm{CO_2}}$ is the atmospheric CO$_2$ concentration [ppm], +$\mu$ is a prescribed scaling factor, +and $a_{\mathrm{LW}}$, $b_{\mathrm{LW}}$, and $c_{\mathrm{LW}}$ are empirical coefficients. + +\cite{Byrne2013} used $(a_{\mathrm{LW}}, b_{\mathrm{LW}}, c_{\mathrm{LW}})= (0.8678,\,1997.9,\,0)$ +with the coefficients obtained by fitting this parameterization to the longwave optical depths employed by ~\cite{Frierson2006}. +Their idealized experiments used a planetary albedo of approximately $0.38$. +For experiments employing an albedo closer to that of Earth, approximately $0.30$, +\cite{Vallis2018} suggested $(a_{\mathrm{LW}}, b_{\mathrm{LW}}, c_{\mathrm{LW}})= (0.1627,\,1997.9,\,0.17)$. +The latter parameter set also includes the explicit parameterization of longwave absorption by CO$_2$, +which was derived using output from the Santa Barbara DISORT Atmospheric Radiative Transfer model (SBDART)~\citep{Ricchiazzi1998}. + + +SCALE-DG uses the coefficient set of \cite{Vallis2018} by default, +while the \cite{Byrne2013} set is also available as an alternative. + + +\subsection{Two-band infrared and one-band solar radiation scheme}\label{sec:two_band_radiation} +%% +This scheme has two infrared bands and one shortwave band, as described in \cite{Geen2016}, +which provides an intermediate complexity between gray radiation and more sophisticated radiation schemes. +All bands were originally parameterized by fitting to data from SBDART for a range of atmospheric profiles. +\cite{Vallis2018} include the CO2 absorption in each band and changes the functional form of the non-window optical depth. +The upward flux of shortwave radiation is assumed to be transparent. + +Following Eq.\,(6) in \cite{Vallis2018}, +the optical thickness in the shortwave band ($< 4~\mu\mathrm{m}$) is parameterized as +%% +\begin{equation} +\frac{\partial \tau_{\mathrm{SW}}}{\partial \sigma} += a_{\mathrm{SW}} + b_{\mathrm{SW}}(\tau_{\mathrm{SW}}) q_v ++ c_{\mathrm{SW}} \ln\!\left(\frac{C_{\mathrm{CO_2}}}{360}\right), +\label{eq:2bandLW_1bandSW_rad_tau_sw} +\end{equation} +%% +where the humidity coefficient depends on the accumulated optical thickness as, +%% +\begin{equation} +\ln b_{\mathrm{SW}} (\tau_{\mathrm{SW}}) = \frac{0.01887}{\tau_{\mathrm{SW}}+0.009522} ++ \frac{1.603}{(\tau_{\mathrm{SW}}+0.5194)^2}. +\label{eq:twoband_bsw} +\end{equation} +%% +This dependence approximately represents the saturation of water-vapor absorption in the shortwave spectrum. + +The longwave optical thickness is divided into a non-window band in 8-14$~\mu\mathrm{m}$ and an atmospheric window band in other long-wave wavelengths ($> 4~\mu\mathrm{m}$). +Following Eq.\,(7) in \cite{Vallis2018}, +the non-window component is parameterized as +%% +\begin{equation} +\frac{\partial \tau_{\mathrm{LW}}^{\mathrm{nw}}}{\partial \sigma} += +a_{\mathrm{LW}} + b_{\mathrm{LW}} \ln\!\left(1+c_{\mathrm{LW}} q_v\right) ++ d_{\mathrm{LW}} \ln\!\left(\frac{C_{\mathrm{CO_2}}}{360~\mathrm{ppm}}\right), +\label{eq:twoband_tau_lw} +\end{equation} +%% +while the atmospheric window band is expressed as +%% +\begin{equation} +\frac{\partial \tau_{\mathrm{LW}}^{\mathrm{win}}}{\partial \sigma} += a_{\mathrm{win}} + b_{\mathrm{win}} q_v + c_{\mathrm{win}} q_v^2 ++ d_{\mathrm{win}} \ln\!\left(\frac{C_{\mathrm{CO_2}}}{360~\mathrm{ppm}}\right). +\label{eq:twoband_tau_window} +\end{equation} +%% +The total longwave flux is obtained as the weighted sum of the two bands. +Fixed weighting factors $R_{\mathrm{nw}}$ and $R_{\mathrm{win}}$ are prescribed for the non-window and window components, respectively. + +The default parameter values adopted in the Isca implementation of Vallis et al.~\cite{Vallis2018} are +%% +\begin{align} +(a_{\mathrm{SW}}, c_{\mathrm{SW}}) &= (0.0596,\,0.0029), \\ +(a_{\mathrm{LW}}, b_{\mathrm{LW}}, c_{\mathrm{LW}}, d_{\mathrm{LW}}) &= (0.1000,\,23.8,\,254.0,\,0.2023), \\ +(a_{\mathrm{win}}, b_{\mathrm{win}}, c_{\mathrm{win}}, d_{\mathrm{win}}) &= (0.2150,\,147.11,\,10814.0,\,0.0954). +\end{align} +%% +The partitioning factors for the longwave radiation are +$R_{\mathrm{nw}}=0.6268$ and $R_{\mathrm{win}}=0.3732$. + + +\subsection{Common formulation} +%% +The gray, one-band, and two-band radiation schemes share the same numerical procedure for evaluating radiative fluxes. +The parameterizations described in the preceding sections determine the optical-thickness gradient $\partial\tau/\partial\sigma$, +whereas the radiative transfer equations are solved using a common discretization and integration procedure. + +Radiative transfer is solved independently for each vertical column. +The vertical grid used for the radiation calculation coincides with the DG nodal grid. +For each interval between adjacent DG nodes, the optical thickness for each band is +approximated by +%% +\begin{equation} +\Delta\tau = \left(\frac{\partial\tau}{\partial\sigma}\right) \Delta\sigma, +\label{eq:dtau} +\end{equation} +%% +where +$\Delta\sigma=\Delta p/P_0$, +and the optical-thickness gradient is evaluated using the arithmetic mean of the thermodynamic variables at the two adjacent nodes. +For the one-band solar radiation scheme, +the coefficient $b_{\mathrm{SW}}$ depends on the accumulated shortwave optical thickness and is therefore evaluated sequentially from the model top downward. + +Radiative fluxes are calculated using the two-stream approximation without scattering. +Assuming that the optical thickness and the Planck source function are constant within each vertical interval, +the analytical solution of the radiative-transfer equation gives the upward and downward longwave fluxes as +%% +\begin{align} +F^{\uparrow}_{k+1} &= F^{\uparrow}_{k} e^{-D\Delta\tau} ++ \mathcal{B}\left(1-e^{-D\Delta\tau}\right), \\ +F^{\downarrow}_{k} &= F^{\downarrow}_{k+1} e^{-D\Delta\tau} +\mathcal{B} \left(1-e^{-D\Delta\tau}\right), +\end{align} +%% +where $k$ is the vertical index, +$D$ is the diffusivity factor, +and $\mathcal{B}=\sigma_{\mathrm{SB}}T^4$ is the blackbody emission. +When the radiation schemes described in Secs.\,\ref{sec:gray_radiation} and \ref{sec:two_band_radiation} are used, we set $D=1$. +For the two-band scheme, the above calculation is applied independently to the non-window and window bands, +and the total longwave flux is obtained as their weighted sum. + +The downward shortwave flux is evaluated only by atmospheric attenuation, +%% +\begin{equation} +F^{\downarrow}_{k} = F^{\downarrow}_{k+1} e^{-\Delta\tau}, +\end{equation} +%% +whereas the upward shortwave flux at the surface is determined by the prescribed surface albedo. +In the present implementation, the atmosphere is assumed to be +transparent to upward shortwave radiation. + + +Finally, the net radiative flux, +%% +\begin{equation} +F_{\mathrm{net}} = F^{\uparrow} - F^{\downarrow}, +\end{equation} +%% +is differentiated in the vertical direction using the DG derivative operator to obtain the radiative heating rate. \ No newline at end of file diff --git a/doc/reference.bib b/doc/reference.bib index 07796e38..68fdf365 100644 --- a/doc/reference.bib +++ b/doc/reference.bib @@ -639,4 +639,86 @@ @article{tiedtke1987parametrization author={Tiedtke, M}, journal={Lecture Notes NWP Course Parametrization of Diabatic Processes, European Centre for Medium-Range Weather Forecasts, Reading, UK}, year={1987} +} + +% Radiation +@book{Liou2002, + author = {Kuo-Nan Liou}, + title = {An Introduction to Atmospheric Radiation}, + edition = {2}, + publisher = {Academic Press}, + address = {San Diego}, + year = {2002}, + isbn = {9780124514515} +} +@article{Toon1989, + author = {Owen B. Toon and Carole P. McKay and Terence P. Ackerman and K. Santhanam}, + title = {Rapid Calculation of Radiative Heating Rates and Photodissociation Rates in Inhomogeneous Multiple Scattering Atmospheres}, + journal = {Journal of Geophysical Research}, + volume = {94}, + number = {D13}, + pages = {16287--16301}, + year = {1989}, +} + +@article{Frierson2006, + author = {Frierson, Dargan M. W. and Held, Isaac M. and Zurita-Gotor, Pablo}, + title = {A Gray-Radiation Aquaplanet Moist {GCM}. + Part I: Static Stability and Eddy Scale}, + journal = {Journal of the Atmospheric Sciences}, + year = {2006}, + volume = {63}, + number = {10}, + pages = {2548--2566}, +} +@article{Byrne2013, + author = {Byrne, Michael P. and O'Gorman, Paul A.}, + title = {Land--Ocean Warming Contrast over a Wide Range of Climates: + Convective Quasi-Equilibrium Theory and Idealized Simulations}, + journal = {Journal of Climate}, + year = {2013}, + volume = {26}, + number = {12}, + pages = {4000--4016}, +} +@article{Geen2016, + author = {Geen, Ruth and Czaja, Arnaud and Haigh, Joanna D.}, + title = {The Effects of Increasing Humidity on Heat Transport by Extratropical Waves}, + journal = {Geophysical Research Letters}, + year = {2016}, + volume = {43}, + number = {15}, + pages = {8314--8321}, + doi = {10.1002/2016GL070214} +} +@article{Vallis2018, + author = { + Vallis, Geoffrey K. + and Colyer, Greg + and Geen, Ruth + and Gerber, Edwin + and Jucker, Martin + and Maher, Penelope + and Paterson, Alexander + and Pietschnig, Marianne + and Penn, James + and Thomson, Stephen I. + }, + title = { {Isca}, v1.0: a framework for the global modelling of the atmospheres of Earth and other planets at varying levels of complexity}, + journal = {Geoscientific Model Development}, + year = {2018}, + volume = {11}, + number = {3}, + pages = {843--859}, + doi = {10.5194/gmd-11-843-2018} +} +@article{Ricchiazzi1998, + author = {Ricchiazzi, Paul and Yang, Shiren and Gautier, Catherine and Sowle, Leslie}, + title = {{SBDART}: A Research and Teaching Software Tool for + Plane-Parallel Radiative Transfer in the Earth's Atmosphere}, + journal = {Bulletin of the American Meteorological Society}, + year = {1998}, + volume = {79}, + number = {10}, + pages = {2101--2114}, } \ No newline at end of file From a79a7837b22cae2c7796116ed8aab80769cb6359 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Fri, 31 Jul 2026 11:28:44 +0900 Subject: [PATCH 3/7] Add boundary layer turbulence model section with PBL parameterizations --- doc/description/pbl_turbulence.tex | 212 +++++++++++++++++++++++++++++ 1 file changed, 212 insertions(+) create mode 100644 doc/description/pbl_turbulence.tex diff --git a/doc/description/pbl_turbulence.tex b/doc/description/pbl_turbulence.tex new file mode 100644 index 00000000..f07f4fd7 --- /dev/null +++ b/doc/description/pbl_turbulence.tex @@ -0,0 +1,212 @@ + +\section{Boundary layer turbulence model} +{\bf \Large +\begin{tabular}{ccc} +\hline + Corresponding author & : & Yuta Kawai\\ +\hline +\end{tabular} +} + +\subsection{Basic formulation of PBL parameterizations} +%% +Unlike the SGS turbulence models described in the previous section, +PBL parameterizations represent the turbulent transport associated with the entire planetary boundary layer (PBL), +including energetic turbulent eddies that are not explicitly resolved by mesoscale atmospheric simulations. +Under the assumption of local horizontal homogeneity, the resolved mean fields are considered to vary more weakly in the horizontal directions than in the vertical direction over the characteristic scales of boundary-layer turbulence. +Conventional PBL parameterizations are therefore commonly formulated as one-dimensional vertical-column models. +The principal effects of unresolved turbulence are represented through vertical fluxes of momentum, heat, moisture, and other scalars. +Horizontal subgrid-scale mixing, when required, is generally treated separately. + + +To describe the evolution of the resolved mean state, +the governing equations are Reynolds-averaged using the Reynolds decomposition as +$\phi = \overline{\phi} + \phi'$ +where the overbar and prime denote the ensemble-mean and turbulent components, respectively. +The mean horizontal momentum equations contain the vertical turbulent momentum fluxes as +% +\begin{align} + \frac{\partial \overline{u}}{\partial t} &= -\frac{\partial \overline{u'w'}}{\partial z} + \mathcal{F}_u, \\ + \frac{\partial \overline{v}}{\partial t} &= -\frac{\partial \overline{v'w'}}{\partial z} + \mathcal{F}_v, +\end{align} +% +where $\mathcal{F}_u$ and $\mathcal{F}_v$ denote the resolved dynamical forcing terms. +Similarly, the mean equation for a scalar quantity $\phi$ is written as +% +\begin{equation} + \frac{\partial \overline{\phi}}{\partial t} = -\frac{\partial \overline{\phi'w'}}{\partial z} + \mathcal{F}_{\phi}. +\end{equation} +%% +These equations are not closed because the turbulent fluxes, such as $\overline{u'w'}$ and $\overline{\phi'w'}$, are additional unknowns, +which constitutes the central problem of turbulence closure. + +In first-order closure models, the turbulent fluxes are directly assumed to be proportional to the corresponding gradients of the resolved mean fields. +In higher-order closure models, including the Mellor--Yamada (MY) family, the turbulent fluxes are instead determined from prognostic or diagnostic equations for second-order turbulent moments. +These equations describe the balances among shear production, buoyancy production or destruction, pressure correlations, turbulent transport, and viscous dissipation. +After suitable closure approximations are introduced, the second-moment equations form a closed system. + + +In the MY level-2 closure, the local time tendencies and turbulent transport terms in the second-moment equations are neglected. +The second moments are therefore assumed to be in local equilibrium, and their governing equations reduce to a set of algebraic relations. +By solving these relations, the vertical turbulent fluxes can be represented by the local gradients of the resolved mean fields as +% +\begin{align} + \overline{u'w'} + &= -K_M \frac{\partial \overline{u}}{\partial z}, \label{eq:eddy_visc_hypo}\\ + \overline{\phi'w'} + &= -K_H \frac{\partial \overline{\phi}}{\partial z} \label{eq:eddy_diff_hypo}, +\end{align} +% +where $K_M$ and $K_H$ denote the eddy viscosity and eddy diffusivity, respectively. +Thus, in the MY framework, the gradient-diffusion form is not merely imposed as an independent empirical assumption, +but arises as a constitutive relation obtained from the locally equilibrated second-moment equations. +The eddy viscosity and diffusivity are written as +% +\begin{align} + K_M &= l q S_M, \\ + K_H &= l q S_H, +\label{eq:eddy_viscdiff_1.5ord_tb_closure} +\end{align} +% +where $l$ is the turbulent mixing length, +$q=\sqrt{\overline{u'^2}+\overline{v'^2}+\overline{w'^2}}=\sqrt{2e}$ is the turbulent velocity scale, +and $S_M$ and $S_H$ are stability functions for momentum and scalar transport, respectively. +The factors $lq$ provide the characteristic turbulent diffusivity scale, +whereas the stability functions arise from the algebraic solution of the second-moment equations and represent the effects of shear and buoyancy on the turbulent fluxes. +Under unstable stratification, buoyancy production enhances turbulent motions and generally increases the mixing coefficients. +Under stable stratification, buoyancy suppresses vertical motions and reduces turbulent transport. +Different PBL parameterization schemes therefore differ mainly in the closure assumptions used to determine $l$, $q$, $S_M$, and $S_H$, +and in whether these quantities are diagnosed or predicted. + + +\subsection{Mellor--Yamada Nakanishi--Niino (MYNN) level 2.0 model} +%% +This model employs the local-equilibrium form of the MYNN second-moment closure. +In this model, the turbulent velocity scale $q$ is diagnosed from the local balance of turbulent kinetic energy (TKE), +rather than predicted from a time-dependent TKE equation. +The local TKE balance is written as +%% +\begin{equation} +P + B - \varepsilon = 0, +\end{equation} +%% +where $P$ is the shear production, $B$ is the buoyancy production, and $\varepsilon$ is the dissipation rate of TKE. +For a horizontally homogeneous flow, these terms are expressed as +%% +\begin{align} + P = \dfrac{\partial \bar{u}}{\partial z} \overline{u'w'} + + \dfrac{\partial \bar{v}}{\partial z} \overline{v'w'}, \;\;\; + B = - \dfrac{g}{\theta_0} \overline{\theta'w'}, \;\;\; + \varepsilon = \dfrac{q^3}{B_1 l}, +\label{eq:mynn_l2_tke_terms} +\end{align} +%% +where $\theta_0$ is a reference potential temperature and $B_1$ is a closure constant. +Using the gradient-diffusion relations in Eqs.\,\eqref{eq:eddy_visc_hypo} and \eqref{eq:eddy_diff_hypo}, +the production terms become +% +\begin{align} + P \sim K_M S^2, \;\; B \sim -K_H N^2, +\label{eq:mynn_l2_production} +\end{align} +% +where $S^2$ and $N^2$ are the squared magnitude of the vertical wind shear and the squared Brunt--V\"ais\"al\"a frequency, respectively. + + +Substituting $K_M=lqS_M$ and $K_H=lqS_H$ into the local TKE balance provides +%% +\begin{equation} +lq\left(S_M S^2-S_H N^2\right) - \frac{q^3}{B_1l} = 0. +\end{equation} +%% +The diagnosed turbulent velocity scale therefore satisfies +%% +\begin{equation} +q^2 = B_1 l^2 S_M S^2(1-R_f) +\label{eq:TKE_local_equilibrium} +\end{equation} +%% +where $R_f$ is the flux Richardson number as +% +\begin{equation} +R_f = \frac{S_HN^2}{S_MS^2} = \frac{K_HN^2}{K_MS^2}. +\label{eq:mynn_flux_richardson} +\end{equation} +% + +The stability functions are obtained by solving the locally equilibrated algebraic equations for the Reynolds stresses, turbulent heat fluxes, and scalar variance. +Although the detailed elimination of the individual second moments is omitted here, the resulting functions can be expressed in the form +% +\begin{align} +S_H = S_H(R_f), \;\; +S_M = S_M(R_f). +\label{eq:mynn_stability_functions_rf} +\end{align} +% +Thus, the effects of shear and stratification on momentum and heat transport are represented through a single nondimensional parameter, $R_f$. + +The flux Richardson number is related to the gradient Richardson number $Ri=N^2/S^2$ through the stability functions as +% +\begin{equation} +Ri = R_f\frac{S_M(R_f)}{S_H(R_f)}. +\label{eq:mynn_ri_rf_relation} +\end{equation} +% +Because the ratio $S_M/S_H$ is itself a known function of $R_f$, this algebraic relation can be inverted to obtain +% +\begin{equation} +R_f=R_f(Ri), +\label{eq:mynn_ri_rf_relation_quadratic} +\end{equation} +% +which corresponds to solving a quadratic equation for $R_f$. + + +The diagnostic procedure of the MYNN level 2.0 closure may therefore be summarized as +% +\begin{equation} +Ri +;\longrightarrow; +R_f +;\longrightarrow; +\left(S_M,S_H\right) +;\longrightarrow; +q +;\longrightarrow; +\left(K_M,K_H\right). +\label{eq:mynn_l2_diagnostic_sequence} +\end{equation} +% +$Ri$ is first calculated from the resolved vertical gradients. +The flux Richardson number is then obtained from the algebraic MYNN relation in Eq.\,\eqref{eq:mynn_ri_rf_relation_quadratic}, +followed by the evaluation of the stability functions in Eq.\,\eqref{eq:mynn_stability_functions_rf}. +The turbulent velocity scale is diagnosed from Eq.\,\eqref{eq:TKE_local_equilibrium}, +and the eddy viscosity and diffusivity are finally evaluated using Eq.\,\eqref{eq:eddy_viscdiff_1.5ord_tb_closure}. +The mixing length $l$ is specified through a separate MYNN mixing-length formulation. + +The complete algebraic formulation, including the closure constants and analytical expressions for the stability functions, +is given in Appendix A of \cite{Nakanishi2006}. + + +\subsection{Common formulation} +%% +Although different PBL parameterization schemes employ different formulations for the turbulent mixing coefficients, +their effects are incorporated into the governing equations in SCALE-DG through a common vertical diffusion operator. +For momentum, potential temperature, moisture, and other prognostic scalars, the turbulent fluxes are represented by the gradient-diffusion relations, +and the corresponding tendency is written in the generic form +% +\begin{equation} +\frac{\partial \rho\phi}{\partial t} = \frac{\partial}{\partial z}\left(\rho K_\phi \frac{\partial \phi}{\partial z}\right), +\end{equation} +% +where $K_\phi$ denotes the turbulent diffusion coefficient determined by the selected PBL parameterization. + +The vertical diffusion operators are spatially discretized using the symmetric interior penalty (SIP) formulation\citep{Arnold2002}. +This formulation weakly couples neighboring DG elements through consistent numerical fluxes and penalty terms, +providing a stable and high-order discretization of the diffusion operator. +We evaluate the turbulent diffusion coefficients at the LGL nodes and enter the SIP discretization as spatially varying diffusion coefficients. + + +The characteristic time scale associated with vertical turbulent diffusion is often much shorter than that of the resolved atmospheric dynamics, especially under strongly mixed boundary-layer conditions. +To avoid the severe timestep restriction imposed by explicit integration, +the vertical diffusion terms are integrated implicitly using the backward Euler method. From 78fd51bc5d9073069baab3874aac8070f68c388b Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Fri, 31 Jul 2026 11:29:08 +0900 Subject: [PATCH 4/7] Add labels to sections for improved referencing --- doc/description/dynamics_spatial_discretization.tex | 2 +- doc/description/surface_flux.tex | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/doc/description/dynamics_spatial_discretization.tex b/doc/description/dynamics_spatial_discretization.tex index 0802c1b6..eb61488b 100644 --- a/doc/description/dynamics_spatial_discretization.tex +++ b/doc/description/dynamics_spatial_discretization.tex @@ -6,7 +6,7 @@ \end{tabular} } -\section{Spatial discretization} +\section{Spatial discretization}\label{sec:dynamics_spatial_discretization} %% We perform the spatial discretization for Eq.\,\eqref{eq:GovernEq_NonhydrostaticDynCore} based on a nodal DGM \citep[e.g.,][]{hesthaven2007nodal}. diff --git a/doc/description/surface_flux.tex b/doc/description/surface_flux.tex index 1fb23154..72941c32 100644 --- a/doc/description/surface_flux.tex +++ b/doc/description/surface_flux.tex @@ -35,7 +35,7 @@ \subsection{Constant flux} The minimum values of the wind speed and $C_m$ can be modified from the default values (i.e., zero) to avoid numerical instability. -\subsection{Simple bulk flux scheme} +\subsection{Simple bulk flux scheme}\label{ssec:simple_sfcflux_bulk} %% This option calculates the surface fluxes based on constant values of the bulk drag coefficient $C_m$, the heat-transfer coefficient $C_h$, and the moisture-transfer coefficient $C_e$. From a668c80d59fdee0d94442900c80f14e6e598be68 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Fri, 31 Jul 2026 11:29:27 +0900 Subject: [PATCH 5/7] Update turbulence section --- doc/description/turbulence.tex | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/doc/description/turbulence.tex b/doc/description/turbulence.tex index c08c8528..067e119a 100644 --- a/doc/description/turbulence.tex +++ b/doc/description/turbulence.tex @@ -7,7 +7,11 @@ \section{Turbulence}\label{sec:turbulence} \end{tabular} } -\subsection{Smagorinsky-type model} +This section describes the subgrid-scale (SGS) turbulence model used for large-eddy simulation (LES) in SCALE-DG. +In LES, the large energy-containing turbulent eddies are explicitly resolved by the numerical model, +while the effects of unresolved SGS motions are represented through eddy-viscosity models. + +\subsection{Smagorinsky-type model}\label{sec:turbulence_model_SMG} %% As a turbulent model, this subsection describes @@ -126,3 +130,13 @@ \subsection{Smagorinsky-type model} For further details of the turbulent model, refer to Sect.\,2.2 of \cite{nishizawa2015influence}. +\subsection{Common formulation} +%% +Although the formulations of the SGS turbulence models may differ, +they all provide the parameterized SGS momentum and scalar fluxes required by the governing equations. +Once evaluated, these fluxes are incorporated into the governing equations as viscous and diffusive flux terms. +Their spatial discretization follows the common DG formulation for the diffusion terms described in Sect.~\ref{sec:dynamics_spatial_discretization}. + +The SGS flux terms are treated explicitly in time. +In the present implementation, their tendencies are integrated using the forward Euler method. +Therefore, the timestep must satisfy the stability constraint associated with the SGS eddy viscosity and diffusivity. \ No newline at end of file From 41dc903791787ffc423db29dc3cfe343caf00df9 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Fri, 31 Jul 2026 11:30:13 +0900 Subject: [PATCH 6/7] Add new sections for PBL turbulence, radiation, and ocean and update reference.bib --- doc/description/scale_dg_description.tex | 11 ++++++++++ doc/reference.bib | 26 ++++++++++++++++++++++-- 2 files changed, 35 insertions(+), 2 deletions(-) diff --git a/doc/description/scale_dg_description.tex b/doc/description/scale_dg_description.tex index 12fe3fed..eefc2f08 100644 --- a/doc/description/scale_dg_description.tex +++ b/doc/description/scale_dg_description.tex @@ -55,8 +55,19 @@ \chapter{Discretization of dynamics} \newpage \chapter{Physical parameterization} \input{turbulence} +\input{pbl_turbulence} +% +\newpage \input{microphysics} +% +\newpage +\input{radiation} +% +\newpage \input{surface_flux} +% +\newpage +\input{ocean} %%% diff --git a/doc/reference.bib b/doc/reference.bib index 68fdf365..c942c717 100644 --- a/doc/reference.bib +++ b/doc/reference.bib @@ -146,6 +146,7 @@ @book{hesthaven2007nodal publisher={Springer Science \& Business Media}, doi = {10.1007/978-0-387-72067-8} } +% Diffusion operator in DGM @article{BASSI1997267, author = {F. Bassi and S. Rebay}, journal = {Journal of Computational Physics}, @@ -155,6 +156,16 @@ @article{BASSI1997267 volume = {131}, year = {1997}, } +@article{Arnold2002, + author = {Douglas N. Arnold and Franco Brezzi and Bernardo Cockburn and L. Donatella Marini}, + title = {Unified Analysis of Discontinuous Galerkin Methods for Elliptic Problems}, + journal = {SIAM Journal on Numerical Analysis}, + volume = {39}, + number = {5}, + pages = {1749--1779}, + year = {2002}, + doi = {10.1137/S0036142901384162} +} % Stabilization in high-order FEM @article{WINTERS20181, author = {Andrew R. Winters and Rodrigo C. Moura and Gianmarco Mengaldo and Gregor J. Gassner and Stefanie Walch and Joaquim Peiro and Spencer J. Sherwin}, @@ -509,8 +520,6 @@ @article{Scotti_1993 volume = {5}, pages = {2306--2308}, } - - % Formulation of eddy viscous terms in cubed sphere coordinate @article{Rancic2017Nonhydro, author = {Miodrag Ran{\v c}i{\'c} and R. James Purser and Du{\v s}an Jovi{\'c} and Ratko Vasic and Thomas Black}, @@ -524,6 +533,19 @@ @article{Rancic2017Nonhydro doi = {10.1175/MWR-D-16-0178.1} } +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% RANS, turbulent models +@article{Nakanishi2006, + author = {Nakanishi, Mikio and Niino, Hiroshi}, + title = {An Improved Mellor--Yamada Level-3 Model: Its Numerical Stability and Application to a Regional Prediction of Advection Fog}, + journal = {Boundary-Layer Meteorology}, + year = {2006}, + volume = {119}, + number = {2}, + pages = {397--407}, + doi = {10.1007/s10546-005-9030-8} +} + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % SCALE @article{nishizawa2015influence, From 30b0448bbacd0633b0cd3ab31829862c95c399c6 Mon Sep 17 00:00:00 2001 From: Yuta Kawai Date: Fri, 31 Jul 2026 11:47:40 +0900 Subject: [PATCH 7/7] Fix typos --- doc/description/ocean.tex | 6 +++--- doc/description/pbl_turbulence.tex | 2 +- doc/description/radiation.tex | 8 ++++---- 3 files changed, 8 insertions(+), 8 deletions(-) diff --git a/doc/description/ocean.tex b/doc/description/ocean.tex index b3c59c51..b44693dc 100644 --- a/doc/description/ocean.tex +++ b/doc/description/ocean.tex @@ -13,14 +13,14 @@ \subsection{A single-layer slab ocean model} The model represents only the thermodynamics of the upper ocean and predicts the evolution of the sea surface temperature (SST). The slab ocean is treated independently at each DG node, without horizontal ocean circulation or vertical mixing below the mixed layer. This model is intended for idealized experiments and for studies in which ocean dynamics play only a secondary role. -This model only treats a thermodynamics in a single-layer ocean, which calculate the evolution of sea temperature. +This model treats only the thermodynamics of a single-layer ocean and predicts the evolution of the sea surface temperature. The governing equation is described using the internal energy of ocean as %% \begin{align} \dfrac{\partial (C_o e_{{\rm in},o})}{\partial t} =& R_{\rm SW,dn} - R_{\rm SW,up} + R_{\rm LW,dn} - R_{\rm LW,up} - F_{\rm SH} - F_{\rm LH} \nonumber \\ - +& e_{\rm in,precip} + Q_{\rm ext}, \label{eq:goevern_eq_slab_ocean} + +& e_{\rm in,precip} + Q_{\rm ext}, \label{eq:govern_eq_slab_ocean} \end{align} %% where $C_o$ is the heat capacity of the slab ocean and $e_{{\rm in},o}$ is its internal energy. @@ -62,4 +62,4 @@ \subsection{A single-layer slab ocean model} Because the slab ocean possesses a much larger heat capacity than the atmosphere, its characteristic evolution time is much longer than the atmospheric time step. -Therefore, the governing equation (Eq.\,\eqref{eq:goevern_eq_slab_ocean}) is integrated using the forward Euler method. +Therefore, the governing equation Eq.\,\eqref{eq:govern_eq_slab_ocean} is integrated using the forward Euler method. diff --git a/doc/description/pbl_turbulence.tex b/doc/description/pbl_turbulence.tex index f07f4fd7..4bedf9f9 100644 --- a/doc/description/pbl_turbulence.tex +++ b/doc/description/pbl_turbulence.tex @@ -201,7 +201,7 @@ \subsection{Common formulation} % where $K_\phi$ denotes the turbulent diffusion coefficient determined by the selected PBL parameterization. -The vertical diffusion operators are spatially discretized using the symmetric interior penalty (SIP) formulation\citep{Arnold2002}. +The vertical diffusion operators are spatially discretized using the symmetric interior penalty (SIP) formulation~\citep{Arnold2002}. This formulation weakly couples neighboring DG elements through consistent numerical fluxes and penalty terms, providing a stable and high-order discretization of the diffusion operator. We evaluate the turbulent diffusion coefficients at the LGL nodes and enter the SIP discretization as spatially varying diffusion coefficients. diff --git a/doc/description/radiation.tex b/doc/description/radiation.tex index a47e4f1d..53ae6da8 100644 --- a/doc/description/radiation.tex +++ b/doc/description/radiation.tex @@ -81,7 +81,7 @@ \subsection{Basic formulation of radiative transfer} - 2\pi (1-\omega_\lambda) B_\lambda - \dfrac{\omega_\lambda}{2}\int_0^1 d\mu \int_{4\pi} P_\lambda(\bm{\Omega},\bm{\Omega}') I_\lambda(\bm{\Omega}') d\Omega', \\ \dfrac{dF_\lambda^\downarrow}{d\tau_{\lambda,z}} =& - - 4\pi J_\lambda^\uparrow + - 4\pi J_\lambda^\downarrow + 2\pi (1-\omega_\lambda) B_\lambda + \dfrac{\omega_\lambda}{2}\int_{-1}^0 d\mu \int_{4\pi} P_\lambda(\bm{\Omega},\bm{\Omega}') I_\lambda(\bm{\Omega}') d\Omega'. \label{eq:rte_general_form_plane_parallel_flux} \end{align} @@ -117,7 +117,7 @@ \subsection{Basic formulation of radiative transfer} \label{eq:rte_twostream_approx} \end{align} %% -where $D_\lambda^\uparrow$ and $D_\lambda^\downarrow$ the diffusivity factors determined by the adopted angular quadrature, +where $D_\lambda^\uparrow$ and $D_\lambda^\downarrow$ are the diffusivity factors determined by the adopted angular quadrature, $S_\lambda^\uparrow$ and $S_\lambda^\downarrow$ denote the corresponding hemispherically averaged source functions. Note that the source term $\pi S_\lambda$ is not obtained by directly evaluating $2\pi\int_0^1 S_\lambda\,d\mu$. Rather, it appears as part of the angular closure introduced in the two-stream approximation. @@ -182,8 +182,8 @@ \subsection{Two-band infrared and one-band solar radiation scheme}\label{sec:two This scheme has two infrared bands and one shortwave band, as described in \cite{Geen2016}, which provides an intermediate complexity between gray radiation and more sophisticated radiation schemes. All bands were originally parameterized by fitting to data from SBDART for a range of atmospheric profiles. -\cite{Vallis2018} include the CO2 absorption in each band and changes the functional form of the non-window optical depth. -The upward flux of shortwave radiation is assumed to be transparent. + \cite{Vallis2018} includes the CO$_2$ absorption in each band and changes the functional form of the non-window optical depth. + The upward flux of shortwave radiation is assumed to be transparent. Following Eq.\,(6) in \cite{Vallis2018}, the optical thickness in the shortwave band ($< 4~\mu\mathrm{m}$) is parameterized as