Update

Update

Module for integrating the prognostic equations.

Provides functions for updating the prognostic variables at the various stages of the semi-implicit time scheme.

See also

• PinCFlow.Types

• PinCFlow.Boundaries

Cartesian

Singleton for transformations from the terrain-following system to the Cartesian one.

KEK

Singleton for the eddy diffusion coefficient for turbulent kinetic energy.

KH

Singleton for the eddy diffusion coefficient for heat.

KM

Singleton for the eddy diffusion coefficient for momentum.

LHS

Singleton for the integration of the left-hand side of an equation.

RHS

Singleton for the integration of the right-hand side of an equation.

Transformed

Singleton for transformations from the Cartesian system to the terrain-following one.

X

Singleton for $\hat{x}$-axis along which a calculation should be performed.

Y

Singleton for $\hat{y}$-axis along which a calculation should be performed.

Z

Singleton for $\hat{z}$-axis along which a calculation should be performed.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::AbstractPredictand,
)

Perform an implicit substep to integrate the Rayleigh-damping term that represents the LHS sponge in the prognostic equation for variable by dispatching to the appropriate model-specific method.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::Rho,
    model::Val{:Boussinesq},
)

Return in Boussinesq mode (constant density).

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::Rho,
    model::Union{Val{:PseudoIncompressible}, Val{:Compressible}},
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the continuity equation.

The update is given by

\[\rho \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} \left(\rho + \alpha_\mathrm{R} \Delta t \bar{\rho}\right),\]

where $\alpha_\mathrm{R}$ is the Rayleigh-damping coefficient computed by PinCFlow.Update.compute_sponges! and $\Delta t$ is the time step given as input to this method.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::RhoP,
    model::Val{:Compressible},
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the auxiliary equation in compressible mode.

The update is given by

\[\rho' \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} \left[\rho' + \alpha_\mathrm{R} \Delta t \bar{\rho} \left(1 - \frac{P}{\rho \bar{\theta}}\right)\right].\]

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::RhoP,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the auxiliary equation in non-compressible modes.

The update is given by

\[\rho' \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} \rho'.\]

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::U,
    model::Union{
        Val{:Boussinesq},
        Val{:PseudoIncompressible},
        Val{:Compressible},
    },
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the zonal-momentum equation.

The update is given by

\[u_{i + 1 / 2} \rightarrow \left(1 + \alpha_{\mathrm{R}, i + 1 / 2} \Delta t\right)^{- 1} \left(u_{i + 1 / 2} + \alpha_{\mathrm{R}, i + 1 / 2} \Delta t u_{\mathrm{R}, i + 1 / 2}\right).\]

If state.namelists.sponge.relax_to_mean is false, $u_{\mathrm{R}, i + 1 / 2}$ is computed with state.namelist.sponge.relaxed_u. Otherwise, it is replaced with the average of $u_{i + 1 / 2}$ across the terrain-following coordinate surface.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::V,
    model::Union{
        Val{:Boussinesq},
        Val{:PseudoIncompressible},
        Val{:Compressible},
    },
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the meridional-momentum equation.

The update is given by

\[v_{j + 1 / 2} \rightarrow \left(1 + \alpha_{\mathrm{R}, j + 1 / 2} \Delta t\right)^{- 1} \left(v_{j + 1 / 2} + \alpha_{\mathrm{R}, j + 1 / 2} \Delta t v_{\mathrm{R}, j + 1 / 2}\right).\]

If state.namelists.sponge.relax_to_mean is false, $v_{\mathrm{R}, j + 1 / 2}$ is computed with state.namelist.sponge.relaxed_v. Otherwise, it is replaced with the average of $v_{j + 1 / 2}$ across the terrain-following coordinate surface.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::W,
    model::Union{
        Val{:Boussinesq},
        Val{:PseudoIncompressible},
        Val{:Compressible},
    },
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the transformed-vertical-momentum equation.

The update is given by

\[\hat{w}_{k + 1 / 2} \rightarrow \left(1 + \alpha_{\mathrm{R}, k + 1 / 2} \Delta t\right)^{- 1} \left(\hat{w}_{k + 1 / 2} + \alpha_{\mathrm{R}, k + 1 / 2} \Delta t \hat{w}_{\mathrm{R}, k + 1 / 2}\right).\]

If state.namelists.sponge.relax_to_mean is false, $\hat{w}_{\mathrm{R}, k + 1 / 2}$ is computed with the functions relaxed_u, relaxed_v and relaxed_w in state.namelists.sponge. Otherwise, it is replaced with the average of $\hat{w}_{k + 1 / 2}$ across the terrain-following coordinate surface.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::PiP,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)

Return in non-compressible modes (Exner-pressure fluctuations are only updated in the corrector step).

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::PiP,
    model::Val{:Compressible},
)

Update the Exner-pressure fluctuations to account for the Rayleigh damping applied to the mass-weighted potential temperature.

The update is given by

\[\pi' \rightarrow \pi' - \alpha_\mathrm{R} \Delta t P \frac{\partial \pi'}{\partial P} \left(1 - \frac{\bar{\rho}}{\rho}\right).\]

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::P,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)

Return in non-compressible modes (mass-weighted potential temperature is constant in time).

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::P,
    model::Val{:Compressible},
)

Integrate the Rayleigh-damping term that represents the LHS sponge in the thermodynamic-energy equation.

The update is given by

\[P \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} P \left(1 + \alpha_\mathrm{R} \Delta t \frac{\bar{\rho}}{\rho}\right).\]

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::Chi,
)

Integrate the Rayleigh-damping terms that represent the LHS sponge in the tracer equations by dispatching to the appropriate method.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::Chi,
    tracer_setup::Val{:NoTracer},
)

Return for configurations without tracer transport.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::Chi,
    tracer_setup::Val{:TracerOn},
)

Integrate the Rayleigh-damping terms that represent the LHS sponge in the tracer equations if state.namelists.tracer.apply_lhs_sponge_to_tracer is true.

In each tracer equation, the update is given by

\[\left(\rho \chi\right) \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} \left(\rho \chi + \alpha_\mathrm{R} \Delta t \rho \chi_\mathrm{R}\right),\]

where $\chi_\mathrm{R}$ is computed with the function relaxed_chi in state.namelists.tracer.

apply_lhs_sponge!(
    state::State,
    dt::AbstractFloat,
    time::AbstractFloat,
    variable::TKE,
)

Integrate the Rayleigh-damping terms that represent the LHS sponge in the turbulent kinetic energy equation.

In the equation for the turbulent kinetic energy, the update is given by

\[\left(\rho e_\mathrm{k}\right) \rightarrow \left(1 + \alpha_\mathrm{R} \Delta t\right)^{- 1} \left(\rho e_\mathrm{k} + \alpha_\mathrm{R} \Delta t \rho e_\mathrm{k,\mathrm{min}}\right),\]

where $e_\mathrm{k,\mathrm{min}}$ represents the minimum allowed TKE value stored in state.turbulence.turbulenceconstants.tkemin.

Arguments

• state: Model state.

• dt: Time step.

• time: Simulation time.

• variable: Variable to apply Rayleigh damping to.

• model: Dynamic equations.

• tracer_setup: General tracer-transport configuration.

check_tke!(state::State)

Enforce a minimum mass-specific turbulent kinetic energy value set by state.turbulence.turbulenceconstants.tkemin.

Arguments

• state: Model state.

compute_buoyancy_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::Union{RhoP, W},
)::AbstractFloat

Compute the factor by which the buoyancy term should be multiplied at $\left(i, j, k\right)$ or $\left(i, j, k + 1 / 2\right)$, by dispatching to a method specific for the dynamic equations and variable, and return the result.

In pseudo-incompressible mode, $\rho'$ are deviations of the total density from $\bar{\rho}$, which describes the reference atmosphere. However, in compressible mode, $\rho' = \rho - P / \bar{\theta}$ does not reduce to this, i.e. the density background has a spatiotemporal dependence. As a consequence, the right-hand side of the prognostic equation for $\rho'$ is given by

\[\left(\frac{\partial \rho'}{\partial t}\right)_{N^2} = f_{\rho'} \frac{N^2 \rho w}{g},\]

with $f_{\rho'} = \bar{\rho} / \rho$ in pseudo-incompressible mode and $f_{\rho'} = P / \left(\rho \bar{\theta}\right)$ in compressible mode. This method returns either $f_{\rho'}$ at $\left(i, j, k\right)$ or $f_w$, which is the interpolation of $f_{\rho'}$ to $\left(i, j, k + 1 / 2\right)$, based on the type of variable.

compute_buoyancy_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::RhoP,
    model::Val{:Compressible},
)::AbstractFloat

Return $f_{\rho'} = P / \left(\rho \bar{\theta}\right)$ as the factor by which the buoyancy term should be multiplied at $\left(i, j, k\right)$ in compressible mode.

compute_buoyancy_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::RhoP,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)::AbstractFloat

Return $f_{\rho'} = \bar{\rho} / \rho$ as the factor by which the buoyancy term should be multiplied at $\left(i, j, k\right)$ in pseudo-incompressible mode (this method is also used in Boussinesq mode).

compute_buoyancy_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    model::Val{:Compressible},
)::AbstractFloat

Return $f_w = \left(P / \bar{\theta}\right)_{k + 1 / 2} / \rho_{k + 1 / 2}$ as the factor by which the buoyancy term should be multiplied in compressible mode.

compute_buoyancy_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)::AbstractFloat

Return $f_w = \bar{\rho}_{k + 1 / 2} / \rho_{k + 1 / 2}$ as the factor by which the buoyancy term should be multiplied at $\left(i, j, k + 1 / 2\right)$ in pseudo-incompressible mode (this method is also used in Boussinesq mode).

Arguments

• state: Model state.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• variable: Variable for which the factor is needed.

• model: Dynamic equations.

compute_compressible_wind_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::Union{U, V, W},
)::AbstractFloat

Compute the factor by which the wind should be multiplied at $\left(i + 1 / 2, j, k\right)$, $\left(i, j + 1 / 2, k\right)$ or $\left(i, j, k + 1 / 2\right)$, by dispatching to a method specific for the dynamic equations and the component indicated by variable, and return the result.

In compressible mode, the Euler steps that are used to integrate the right-hand side of the momentum equation update $\left(J P\right)_{i + 1 / 2} u_{i + 1 / 2}$, $\left(J P\right)_{j + 1 / 2} v_{j + 1 / 2}$ and $\left(J P\right)_{k + 1 / 2} \hat{w}_{k + 1 / 2}$ instead of $u_{i + 1 / 2}$, $v_{j + 1 / 2}$ and $\hat{w}_{k + 1 / 2}$.

compute_compressible_wind_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::Union{U, V, W},
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)::AbstractFloat

Return $1$ as the factor by which the wind should be multiplied in non-compressible mode.

compute_compressible_wind_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
    model::Val{:Compressible},
)::AbstractFloat

Return $\left(J P\right)_{i + 1 / 2}$ as the factor by which the zonal wind should be multiplied in compressible mode.

compute_compressible_wind_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
    model::Val{:Compressible},
)::AbstractFloat

Return $\left(J P\right)_{j + 1 / 2}$ as the factor by which the meridional wind should be multiplied in compressible mode.

compute_compressible_wind_factor(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    model::Val{:Compressible},
)::AbstractFloat

Return $\left(J P\right)_{k + 1 / 2}$ as the factor by which the transformed vertical wind should be multiplied in compressible mode.

Arguments

• state: Model state.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• variable: Variable for which the factor is needed.

• model: Dynamic equations.

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
    direction::X,
)::AbstractFloat

Compute and return the diffusive zonal momentum fluxes in $\hat{x}$-direction, i.e.

\[\Xi_u^{\hat{x}} = \frac{u_{i + 1 / 2} - u_{i - 1 / 2}}{\Delta \hat{x}} + G^{13} \frac{u_{k + 1} - u_{k - 1}}{2\Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
    direction::Y,
)::AbstractFloat

Compute and return the diffusive zonal momentum fluxes in $\hat{y}$-direction, i.e.

\[\Xi_u^{\hat{y}} = \frac{u_{j + 1} - u_{j - 1}}{2\Delta \hat{y}} + G^{23} \frac{u_{k + 1} - u_{k - 1}}{2\Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
    direction::Z,
)::AbstractFloat

Compute and return the diffusive zonal momentum fluxes in $\hat{z}$-direction, i.e.

\[\Xi_u^{\hat{z}} = G^{13}\frac{u_{i + 1 / 2} - u_{i - 1 / 2}}{\Delta \hat{x}} + G^{23} \frac{u_{j + 1} - u_{j - 1}}{2 \Delta \hat{y}} + G^{33} \frac{u_{k + 1} - u_{k - 1}}{2 \Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
    direction::X,
)::AbstractFloat

Compute and return the diffusive meridional momentum fluxes in $\hat{x}$-direction, i.e.

\[\Xi_v^{\hat{x}} = \frac{v_{i + 1} - v_{i - 1}}{2 \Delta \hat{x}} + G^{13} \frac{v_{k + 1} - v_{k - 1}}{2 \Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
    direction::Y,
)::AbstractFloat

Compute and return the diffusive meridional momentum fluxes in $\hat{y}$-direction, i.e.

\[\Xi_v^{\hat{y}} = \frac{v_{j + 1 / 2} - v_{j - 1 / 2}}{\Delta \hat{y}} + G^{23} \frac{v_{k + 1} - v_{k - 1}}{2 \Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
    direction::Z,
)::AbstractFloat

Compute and return the diffusive meridional momentum fluxes in $\hat{z}$-direction, i.e.

\[\Xi_v^{\hat{z}} = G^{13}\frac{v_{i + 1} - v_{i - 1}}{2 \Delta \hat{x}} + G^{23} \frac{v_{j + 1 / 2} - v_{j - 1 / 2}}{\Delta \hat{y}} + G^{33} \frac{v_{k + 1} - v_{k - 1}}{2 \Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    direction::X,
)::AbstractFloat

Compute and return the diffusive vertical momentum fluxes in $\hat{x}$-direction, i.e.

\[\Xi_w^{\hat{x}} = \frac{w_{i + 1} - w_{i - 1}}{2 \Delta \hat{x}} + G^{13} \frac{w_{k + 1 / 2} - w_{k - 1 / 2}}{\Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    direction::Y,
)::AbstractFloat

Compute and return the diffusive vertical momentum fluxes in $\hat{y}$-direction, i.e.

\[\Xi_w^{\hat{y}} = \frac{w_{j + 1} - w_{j - 1}}{2 \Delta \hat{y}} + G^{23} \frac{w_{k + 1 / 2} - w_{k - 1 / 2}}{\Delta \hat{z}}.\]

compute_momentum_diffusion_terms(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    direction::Z,
)::AbstractFloat

Compute and return the diffusive vertical momentum fluxes in $\hat{z}$-direction, i.e.

\[\Xi_w^{\hat{z}} = G^{13} \frac{w_{i + 1} - w_{i - 1}}{2 \Delta \hat{x}} + G^{23} \frac{w_{j + 1} - w_{j - 1}}{2 \Delta \hat{y}} + G^{33} \frac{w_{k + 1 / 2} - w_{k - 1 / 2}}{\Delta \hat{z}}.\]

Arguments

• state: Model state.

• i: Grid-cell index on the $\hat{x}$-axis.

• j: Grid-cell index on the $\hat{y}$-axis.

• k: Grid-cell index on the $\hat{z}$-axis.

• variable: Wind direction.

• direction: Direction of the flux.

compute_pressure_gradient(
    state::State,
    pip::AbstractArray{<:AbstractFloat, 3},
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
)::AbstractFloat

Compute and return the pressure(-difference)-gradient term in the zonal-wind equation at $\left(i + 1 / 2, j, k\right)$, using the pressure(-difference) field pip.

The pressure-gradient component is given by

\[\mathcal{P}^{\rho u}_{i + 1 / 2} = \frac{\pi'_{i + 1} - \pi'}{\Delta \hat{x}} + G^{13}_{i + 1 / 2} \frac{\pi'_{i + 1 / 2, k + 1} - \pi'_{i + 1 / 2, k - 1}}{2 \Delta \hat{z}}.\]

Since the Exner-pressure is not known in the vertical ghost cells, a different discretization is needed at the vertical boundaries. At $k = k_0$ (in the first process in $\hat{z}$), the alternative second-order-accurate approximation

\[\mathcal{P}^{\rho u}_{i + 1 / 2} = \frac{\pi'_{i + 1} - \pi'}{\Delta \hat{x}} + G^{13}_{i + 1 / 2} \frac{- \pi'_{i + 1 / 2, k + 2} + 4 \pi'_{i + 1 / 2, k + 1} - 3 \pi'_{i + 1 / 2}}{2 \Delta \hat{z}}\]

is used and, in a similar manner, one has

\[\mathcal{P}^{\rho u}_{i + 1 / 2} = \frac{\pi'_{i + 1} - \pi'}{\Delta \hat{x}} + G^{13}_{i + 1 / 2} \frac{\pi'_{i + 1 / 2, k - 2} - 4 \pi'_{i + 1 / 2, k - 1} + 3 \pi'_{i + 1 / 2}}{2 \Delta \hat{z}}\]

at $k = k_1$ (in the last process in $\hat{z}$). The corresponding pressure-difference-gradient component $\mathcal{D}^{\rho u}_{i + 1 / 2}$ is obtained by replacing $\pi'$ with $\Delta \pi'$. The returned quantity also includes the factor $c_p \left(P_{i + 1 / 2} / \rho_{i + 1 / 2}\right)$.

compute_pressure_gradient(
    state::State,
    pip::AbstractArray{<:AbstractFloat, 3},
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
)::AbstractFloat

Compute and return the pressure-gradient term in the meridional-wind equation at $\left(i, j + 1 / 2, k\right)$, using the pressure(-difference) field pip.

The pressure-gradient component is given by

\[\mathcal{P}^{\rho v}_{j + 1 / 2} = \frac{\pi'_{j + 1} - \pi'}{\Delta \hat{y}} + G^{23}_{j + 1 / 2} \frac{\pi'_{j + 1 / 2, k + 1} - \pi'_{j + 1 / 2, k - 1}}{2 \Delta \hat{z}}.\]

Analogous to the component in the zonal-wind equation, one has

\[\mathcal{P}^{\rho v}_{j + 1 / 2} = \frac{\pi'_{j + 1} - \pi'}{\Delta \hat{y}} + G^{23}_{j + 1 / 2} \frac{- \pi'_{j + 1 / 2, k + 2} + 4 \pi'_{j + 1 / 2, k + 1} - 3 \pi'_{j + 1 / 2}}{2 \Delta \hat{z}}\]

at $k = k_0$ (in the first process in $\hat{z}$) and

\[\mathcal{P}^{\rho v}_{j + 1 / 2} = \frac{\pi'_{j + 1} - \pi'}{\Delta \hat{y}} + G^{23}_{j + 1 / 2} \frac{\pi'_{j + 1 / 2, k - 2} - 4 \pi'_{j + 1 / 2, k - 1} + 3 \pi'_{j + 1 / 2}}{2 \Delta \hat{z}}\]

at $k = k_1$ (in the last process in $\hat{z}$). The corresponding pressure-difference-gradient component $\mathcal{D}^{\rho v}_{j + 1 / 2}$ is obtained by replacing $\pi'$ with $\Delta \pi'$. The returned quantity also includes the factor $c_p \left(P_{j + 1 / 2} / \rho_{j + 1 / 2}\right)$.

compute_pressure_gradient(
    state::State,
    pip::AbstractArray{<:AbstractFloat, 3},
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
)::AbstractFloat

Compute and return the pressure(-difference)-gradient term in the transformed-vertical-wind equation at $\left(i, j, k + 1 / 2\right)$, using the pressure(-difference) field pip.

The pressure-gradient component is given by

\[\begin{align*} \mathcal{P}^{\rho \hat{w}}_{k + 1 / 2} & = G^{13}_{k + 1 / 2} \frac{\pi'_{i + 1, k + 1 / 2} - \pi'_{i - 1, k + 1 / 2}}{2 \Delta \hat{x}} + G^{23}_{k + 1 / 2} \frac{\pi'_{j + 1, k + 1 / 2} - \pi'_{j - 1, k + 1 / 2}}{2 \Delta \hat{y}}\\ & \quad + G^{33}_{k + 1 / 2} \frac{\pi'_{k + 1} - \pi'}{\Delta \hat{z}}. \end{align*}\]

At $k = k_0 - 1$ (in the first process in $\hat{z}$) and $k = k_1$ (in the last process in $\hat{z}$), it is set to zero. The corresponding pressure-difference-gradient component $\mathcal{D}^{\rho \hat{w}}_{k + 1 / 2}$ is obtained by replacing $\pi'$ with $\Delta \pi'$. The returned quantity also includes the factor $c_p \left(P_{k + 1 / 2} / \rho_{k + 1 / 2}\right)$.

Arguments

• state: Model state.

• pip: Pressure field.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• variable: Equation in which the respective pressure-gradient component is needed.

compute_sponges!(state::State, dt::AbstractFloat, time::AbstractFloat)

Compute the Rayleigh-damping coefficients of the two sponges.

The coefficients are computed with the functions lhs_sponge and rhs_sponge in state.namelists.sponge.

Arguments

• state: Model state.

• dt: Time step.

• time: Simulation time.

See also

• PinCFlow.Boundaries.set_zonal_boundaries_of_field!

• PinCFlow.Boundaries.set_meridional_boundaries_of_field!

compute_stress_tensor(
    i::Integer,
    j::Integer,
    k::Integer,
    mu::Integer,
    nu::Integer,
    state::State,
)::AbstractFloat

Compute and return the element $\left(\mu, \nu\right)$ of the Cartesian stress tensor at the grid point $\left(i, j, k\right)$.

The discretized elements of the Cartesian stress tensor are given by

\[\begin{align*} \Pi^{1 1} & = \frac{2}{\Delta \hat{x}} \left(u_{i + 1 / 2} - u_{i - 1 / 2}\right) + \frac{G^{1 3}}{\Delta \hat{z}} \left(u_{k + 1} - u_{k - 1}\right) - \frac{2}{3} \delta,\\ \Pi^{1 2} & = \frac{1}{2 \Delta \hat{y}} \left(u_{j + 1} - u_{j - 1}\right) + \frac{G^{2 3}}{2 \Delta \hat{z}} \left(u_{k + 1} - u_{k - 1}\right) + \frac{1}{2 \Delta \hat{x}} \left(v_{i + 1} - v_{i - 1}\right) + \frac{G^{1 3}}{2 \Delta \hat{z}} \left(v_{k + 1} - v_{k - 1}\right),\\ \Pi^{1 3} & = \frac{1}{2 J \Delta \hat{z}} \left(u_{k + 1} - u_{k - 1}\right) + \frac{1}{2 \Delta \hat{x}} \left(w_{i + 1} - w_{i - 1}\right) + \frac{G^{1 3}}{\Delta \hat{z}} \left(w_{k + 1 / 2} - w_{k - 1 / 2}\right),\\ \Pi^{2 2} & = \frac{2}{\Delta \hat{y}} \left(v_{j + 1 / 2} - v_{j - 1 / 2}\right) + \frac{G^{2 3}}{\Delta \hat{z}} \left(v_{k + 1} - v_{k - 1}\right) - \frac{2}{3} \delta,\\ \Pi^{2 3} & = \frac{1}{2 J \Delta \hat{z}} \left(v_{k + 1} - v_{k - 1}\right) + \frac{1}{2 \Delta \hat{y}} \left(w_{j + 1} - w_{j - 1}\right) + \frac{G^{2 3}}{\Delta \hat{z}} \left(w_{k + 1 / 2} - w_{k - 1 / 2}\right),\\ \Pi^{3 3} & = \frac{2}{J \Delta \hat{z}} \left(w_{k + 1 / 2} - w_{k - 1 / 2}\right) - \frac{2}{3} \delta, \end{align*}\]

where

\[\begin{align*} \delta & = \frac{1}{J} \left[\frac{1}{\Delta \hat{x}} \left(J_{i + 1 / 2} u_{i + 1 / 2} - J_{i - 1 / 2} u_{i - 1 / 2}\right) + \frac{1}{\Delta \hat{y}} \left(J_{j + 1 / 2} v_{j + 1 / 2} - J_{j - 1 / 2} v_{j - 1 / 2}\right)\right.\\ & \qquad \quad + \left.\frac{1}{\Delta \hat{z}} \left(J_{k + 1 / 2} \hat{w}_{k + 1 / 2} - J_{k - 1 / 2} \hat{w}_{k - 1 / 2}\right)\right]. \end{align*}\]

Arguments

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• mu: First contravariant tensor index.

• nu: Second contravariant tensor index.

• state: Model state.

See also

• PinCFlow.Update.compute_vertical_wind

compute_vertical_wind(
    i::Integer,
    j::Integer,
    k::Integer,
    state::State,
)::AbstractFloat

Compute and return the Cartesian vertical wind at the grid point $\left(i, j, k + 1 / 2\right)$.

Arguments

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• state: Model state.

See also

• PinCFlow.Update.transform

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::Union{U, V, W, P, Chi},
)::AbstractFloat

Return the volume force in the equation specified by variable, by dispatching to an equation-and-WKB-mode specific method.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::Union{U, V, W, Chi},
    wkb_mode::Val{:NoWKB},
)::AbstractFloat

Return $0$ as the volume force in non-WKB configurations (for all variables except the mass-weighted potential temperature).

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::U,
    wkb_mode::Union{Val{:SteadyState}, Val{:SingleColumn}, Val{:MultiColumn}},
)::AbstractFloat

Return the gravity-wave drag on the zonal momentum, interpolated to $\left(i + 1 / 2, j, k\right)$.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::V,
    wkb_mode::Union{Val{:SteadyState}, Val{:SingleColumn}, Val{:MultiColumn}},
)::AbstractFloat

Return the gravity-wave drag on the meridional momentum, interpolated to $\left(i, j + 1 / 2, k\right)$.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::W,
    wkb_mode::Union{Val{:SteadyState}, Val{:SingleColumn}, Val{:MultiColumn}},
)::AbstractFloat

Return the gravity-wave drag on the transformed vertical momentum, interpolated to $\left(i, j, k + 1 / 2\right)$, as given by

\[\left(\frac{\partial \hat{w}}{\partial t}\right)_\mathrm{w} = \left[G^{1 3} \left(\frac{\partial u}{\partial t}\right)_\mathrm{w}\right]_{k + 1 / 2} + \left[G^{2 3} \left(\frac{\partial v}{\partial t}\right)_\mathrm{w}\right]_{k + 1 / 2}.\]

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::P,
    wkb_mode::Val{:NoWKB},
)::AbstractFloat

Return the conductive heating.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variable::P,
    wkb_mode::Union{Val{:SteadyState}, Val{:SingleColumn}, Val{:MultiColumn}},
)::AbstractFloat

Return the sum of gravity-wave impact on the mass-weighted potential temperature and conductive heating.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variables::Chi,
    wkb_mode::Union{Val{:SteadyState}, Val{:SingleColumn}, Val{:MultiColumn}},
)::AbstractFloat

Return the tracer flux convergence due to gravity waves.

compute_volume_force(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    variables::TKE,
)::AbstractFloat

Return the mass-weighted impact of shear $\mathcal{S}$ and buoyancy $\mathcal{B}$ on the TKE, given by

\[\left(\frac{\partial \rho e_\mathrm{k}}{\partial t}\right) = \rho\mathcal{S} + \rho\mathcal{B}\]

where

\[\begin{align*} \mathcal{S} &= K_\mathrm{M}\left[\left(\frac{\partial u}{\partial \hat{z}}\right)^2 + \left(\frac{\partial v}{\partial \hat{z}}\right)^2\right] \;, \\ \mathcal{B} &= -K_\mathrm{H}\left(N^2 + \frac{\partial b}{\partial \hat{z}}\right) \;, \end{align*}\]

and $K_\mathrm{M}$ and $K_\mathrm{H}$ represent the eddy diffusion coefficients for momentum and heat, respectively.

Arguments

• state: Model state.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• variable: Variable (equation) of choice.

• wkb_mode: Approximations used by MS-GWaM.

See also

• PinCFlow.Update.conductive_heating

• PinCFlow.Update.compute_momentum_diffusion_terms

conductive_heating(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
)::AbstractFloat

Compute and return the conductive heating by dispatching to specialized methods dependent on the model.

conductive_heating(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    model::Val{:Boussinesq},
)::AbstractFloat

Return $0$ as conductive heating in Boussinesq mode.

conductive_heating(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    model::Val{:PseudoIncompressible},
)::AbstractFloat

Return $0$ as conductive heating in PseudoIncompressible mode.

conductive_heating(
    state::State,
    i::Integer,
    j::Integer,
    k::Integer,
    model::Val{:Compressible},
)::AbstractFloat

Compute and return the conductive heating as the convergence of potential temperature fluxes (weighted by the density), i.e.

\[\left(\frac{\partial P}{\partial t}\right)_\lambda = - \frac{\rho}{J} \left(\frac{\mathcal{F}_{i + 1 / 2}^{\theta, \hat{x}} - \mathcal{F}_{i - 1 / 2}^{\theta, \hat{x}}}{\Delta \hat{x}} + \frac{\mathcal{F}_{j + 1 / 2}^{\theta, \hat{y}} - \mathcal{F}_{j - 1 / 2}^{\theta, \hat{y}}}{\Delta \hat{y}} + \frac{\mathcal{F}_{k + 1 / 2}^{\theta, \hat{z}} - \mathcal{F}_{k - 1 / 2}^{\theta, \hat{z}}}{\Delta \hat{z}}\right).\]

Arguments

• state: Model state.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• model: Dynamic equations.

reset_thomas!(state::State)

Reset arrays needed for the Thomas tridiagonal algorithm to $0$.

Arguments

• state: Model state.

thomas_algorithm!(state::State)

Solves a tridiagonal system in $\hat{z}$-direction using the Thomas tridiagonal matrix algorithm (see Durran, 2010) . Since the Thomas algorithm consists of an upward elimination sweep and a downward pass, this method performs sequential one-way MPI communication if the domain is parallelized in the vertical.

The system is defined as:

\[a_k \phi_{k-1} + b_k\phi_k + c_k\phi_{k+1} = f_k\;.\]

The result is stored in state.variables.auxiliaries.fth.

Arguments

• state: Model state.

transform(
    i::Integer,
    j::Integer,
    k::Integer,
    uedger::AbstractFloat,
    uuedger::AbstractFloat,
    uedgel::AbstractFloat,
    uuedgel::AbstractFloat,
    vedgef::AbstractFloat,
    vuedgef::AbstractFloat,
    vedgeb::AbstractFloat,
    vuedgeb::AbstractFloat,
    wedgeu::AbstractFloat,
    coordinate::Cartesian,
    state::State,
)::AbstractFloat

Perform the transformation of a vertical-wind-like variable from the transformed system to the Cartesian one, given the wind-like components at the grid points surrounding $\left(i, j, k + 1 / 2\right)$, and return the result.

The discretized transformation rule for the vertical wind is given by

\[w_{k + 1 / 2} = J_{k + 1 / 2} \left[- \left(G^{1 3} u\right)_{k + 1 / 2} - \left(G^{2 3} v\right)_{k + 1 / 2} + \hat{w}_{k + 1 / 2}\right].\]

transform(
    i::Integer,
    j::Integer,
    k::Integer,
    uedger::AbstractFloat,
    uuedger::AbstractFloat,
    uedgel::AbstractFloat,
    uuedgel::AbstractFloat,
    vedgef::AbstractFloat,
    vuedgef::AbstractFloat,
    vedgeb::AbstractFloat,
    vuedgeb::AbstractFloat,
    wedgeu::AbstractFloat,
    coordinate::Transformed,
    state::State,
)::AbstractFloat

Perform the transformation of a vertical-wind-like variable from the Cartesian system to the transformed one, given the wind-like components at the grid points surrounding $\left(i, j, k + 1 / 2\right)$, and return the result.

The discretized transformation rule for the vertical wind is given by

\[\hat{w}_{k + 1 / 2} = \left(G^{1 3} u\right)_{k + 1 / 2} + \left(G^{2 3} v\right)_{k + 1 / 2} + \frac{w_{k + 1 / 2}}{J_{k + 1 / 2}}.\]

Arguments

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• uedger: Zonal-wind equivalent at $\left(i + 1 / 2, j, k\right)$.

• uuedger: Zonal-wind equivalent at $\left(i + 1 / 2, j, k + 1\right)$.

• uedgel: Zonal-wind equivalent at $\left(i - 1 / 2, j, k\right)$.

• uuedgel: Zonal-wind equivalent at $\left(i - 1 / 2, j, k + 1\right)$.

• vedgef: Meridional-wind equivalent at $\left(i, j + 1 / 2, k\right)$.

• vuedgef: Meridional-wind equivalent at $\left(i, j + 1 / 2, k + 1\right)$.

• vedgeb: Meridional-wind equivalent at $\left(i, j - 1 / 2, k\right)$.

• vuedgeb: Meridional-wind equivalent at $\left(i, j - 1 / 2, k + 1\right)$.

• wedgeu: Transformed-vertical-wind equivalent at $\left(i, j, k + 1 / 2\right)$

• coordinate: Coordinate system to transform to.

• state: Model state.

turbulence_diffusion_coefficient(state::State, i::Integer, j::Integer, k::Integer, variable::KM)

Compute the eddy diffusion coefficient for momentum at $(i, j, k)$.

The eddy diffusion coefficient for momentum is given by

\[ K_M = l_v \sqrt{2 e_\mathrm{k}} \;, \]

with turbulence mixing length $l_v$ stored in state.turbulence.turbulenceconstants.lv.

turbulence_diffusion_coefficient(state::State, i::Integer, j::Integer, k::Integer, variable::KH)

Compute the eddy diffusion coefficient for heat at $(i, j, k)$.

The eddy diffusion coefficient for heat is given by

\[ K_H = l_b \sqrt{2 e_\mathrm{k}} \;, \]

with turbulence mixing length $l_b$ stored in state.turbulence.turbulenceconstants.lb.

turbulence_diffusion_coefficient(state::State, i::Integer, j::Integer, k::Integer, variable::KEK)

Compute the eddy diffusion coefficient for turbulent kinetic energy at $(i, j, k)$.

The eddy diffusion coefficient for turbulent kinetic energy is given by

\[ K_{e_\mathrm{k}} = l_t \sqrt{2 e_\mathrm{k}} \;, \]

with turbulence mixing length $l_t$ stored in state.turbulence.turbulenceconstants.lt.

Arguments

• state: Model state.

• i: Zonal grid-cell index.

• j: Meridional grid-cell index.

• k: Vertical grid-cell index.

• variable: Eddy diffusion coefficient to be computed.

turbulent_diffusion!(state::State, dt::AbstractFloat)

Apply diffusion to the momentum, mass-weighted potential temperature, and tracers by dispatching to turbulence parameterization-specific method.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    turbulence_scheme::Val{:NoTurbulence},
)

Return for configurations without turbulence parameterization.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    turbulence_scheme::Val{:TKEScheme},
)

Apply diffusion by dispatching to specialized methods for momentum, mass-weighted potential temperature, and tracers, according to configurations set by state.namelist.turbulence.momentum_coupling, state.namelist.turbulence.entropy_coupling, and state.namelist.turbulence.tracer_coupling, respectively.

turbulent_diffusion!(state::State, dt::AbstractFloat, variable::U)

Apply diffusion to the zonal momentum.

The prognostic equation

\[\frac{\partial u}{\partial t} = \frac{1}{J}\frac{\partial}{\partial \hat{z}}\left(\frac{K_\mathrm{M}}{J}\frac{\partial u}{\partial \hat{z}}\right)\]

is solved using the Crank-Nicolson scheme, where the system

\[a_{i+1/2,k} u_{i+1/2,k-1}^{n+1} + b_{i+1/2,k} u_{i+1/2,k}^{n+1} + c_{i+1/2,k} u_{i+1/2,k+1}^{n+1} = f_{i+1/2,k}\]

is solved using a Thomas tridiagonal solver, with $\mathcal{K}_\mathrm{M} = \frac{K_\mathrm{M}}{J}$ and

\[\begin{align*} a_{i+1/2,k} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k-1/2}}{J_{i+1/2}} ,\\ b_{i+1/2,k} &= 1 + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k+1/2}}{J_{i+1/2}} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k-1/2}}{J_{i+1/2}} ,\\ c_{i+1/2,k} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k+1/2}}{J_{i+1/2}} ,\\ f_{i+1/2,k} &= \left[1 - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k+1/2}}{J_{i+1/2}} - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k-1/2}}{J_{i+1/2}}\right]v_{i+1/2}^{n} \\ &+\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k+1/2}}{J_{i+1/2}}v_{i+1/2,k+1}^{n} +\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},i+1/2,k-1/2}}{J_{i+1/2}}v_{i+1/2,k-1}^{n} . \end{align*}\]

turbulent_diffusion!(state::State, dt::AbstractFloat, variable::V)

Apply diffusion to the meridional momentum.

The prognostic equation

\[\frac{\partial v}{\partial t} = \frac{1}{J}\frac{\partial}{\partial \hat{z}}\left(\frac{K_\mathrm{M}}{J}\frac{\partial v}{\partial \hat{z}}\right)\]

is solved using the Crank-Nicolson scheme, where the system

\[a_{j+1/2,k} v_{j+1/2,k-1}^{n+1} + b_{j+1/2,k} v_{j+1/2,k}^{n+1} + c_{j+1/2,k} v_{j+1/2,k+1}^{n+1} = f_{j+1/2,k}\]

is solved using a Thomas tridiagonal solver, with $\mathcal{K}_\mathrm{M} = \frac{K_\mathrm{M}}{J}$ and

\[\begin{align*} a_{j+1/2,k} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k-1/2}}{J_{j+1/2}} ,\\ b_{j+1/2,k} &= 1 + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k+1/2}}{J_{j+1/2}} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k-1/2}}{J_{j+1/2}} ,\\ c_{j+1/2,k} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k+1/2}}{J_{j+1/2}} ,\\ f_{j+1/2,k} &= \left[1 - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k+1/2}}{J_{j+1/2}} - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k-1/2}}{J_{j+1/2}}\right]v_{j+1/2}^{n}\\ &+\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k+1/2}}{J_{j+1/2}}v_{j+1/2,k+1}^{n} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},j+1/2,k-1/2}}{J_{j+1/2}}v_{j+1/2,k-1}^{n} . \end{align*}\]

turbulent_diffusion!(state::State, dt::AbstractFloat, variable::W)

Apply diffusion to the vertical momentum.

The prognostic equation

\[\frac{\partial w}{\partial t} = \frac{1}{J}\frac{\partial}{\partial \hat{z}}\left(\frac{K_\mathrm{M}}{J}\frac{\partial w}{\partial \hat{z}}\right)\]

is solved using the Crank-Nicolson scheme, where the system

\[a_{k+1/2} w_{k-1/2}^{n+1} + b_{k+1/2} w_{k+1/2}^{n+1} + c_{k+1/2} w_{k+3/2}^{n+1} = f_{k+1/2}\]

is solved using a Thomas tridiagonal solver, with $\mathcal{K}_\mathrm{M} = \frac{K_\mathrm{M}}{J}$ and

\[\begin{align*} a_{k+1/2} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k}}{J_{k+1/2}^2} ,\\ b_{k+1/2} &= 1 + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k+1}}{J_{k+1/2}^2} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k}}{J_{k+1/2}^2} ,\\ c_{k+1/2} &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k+1}}{J_{k+1/2}^2} ,\\ f_{k+1/2} &= \left[1 - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k+1}}{J_{k+1/2}} - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k}}{J_{k+1/2}}\right]w_{k+1/2}^{n}\\ &+\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k+1}}{J_{k+1/2}}w_{k+1}^{n} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{M},k}}{J_{k+1/2}}w_{k-1}^{n} . \end{align*}\]

Using the relation

\[\frac{\partial \hat{w}}{\partial t} = G^{13}\frac{\partial u}{\partial t} + G^{23}\frac{\partial v}{\partial t} + \frac{1}{J}\frac{\partial w}{\partial t} ,\]

the transformed wind is calculated:

\[\begin{align*} \hat{w}^{n+1}_{k+1/2} &= G_{k+1/2}^{13}u^{n+1}_{k+1/2} + G_{k+1/2}^{23}v^{n+1}_{k+1/2} + \frac{1}{J_{k+1/2}}w^{n+1}_{k+1/2}. \end{align*}\]

turbulent_diffusion!(state::State, dt::AbstractFloat, variable::Theta)

Apply diffusion to the mass-weighted potential temperature by dispatching to model-specific methods.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    variable::Theta,
    model::Union{PseudoIncompressible, Boussinesq},
)

Return for configurations in Boussinesq and pseudo-incompressible mode.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    variable::Theta,
    model::Compressible,
)

Apply diffusion to the potential temperature for configurations in Compressible mode.

The prognostic equation

\[\frac{\partial \theta}{\partial t} = \frac{1}{J}\frac{\partial}{\partial \hat{z}}\left(\frac{K_\mathrm{H}}{J}\frac{\partial \theta}{\partial \hat{z}}\right)\]

is solved using the Crank-Nicolson scheme, where the system

\[a_k \theta_{k-1}^{n+1} + b_k \theta_k^{n+1} + c_k \theta_{k+1}^{n+1} = f_k\]

is solved using a Thomas tridiagonal solver, with $\mathcal{K}_\mathrm{H} = \frac{K_\mathrm{H}}{J}$ and

\[\begin{align*} a_k &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k} , \\ b_k &= 1 + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}{J_k}}, \\ c_k &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} , \\ f_k &= \left[1 - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k}\right] \left(\theta\right)_k^{n} \\ & + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} \left(\theta\right)_{k+1}^{n} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k} \left(\theta\right)_{k-1}^{n}. \end{align*}\]

turbulent_diffusion!(state::State, dt::AbstractFloat, variable::Chi)

Apply diffusion to tracers by dispatching to tracer-setup-specific configurations.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    variable::Chi,
    tracer_setup::NoTracer,
)

Return for configurations without tracer transport.

turbulent_diffusion!(
    state::State,
    dt::AbstractFloat,
    variable::Chi,
    tracer_setup::TracerOn,
)

Apply diffusion to the tracers variables.

The prognostic equation

\[\frac{\partial \chi}{\partial t} = \frac{1}{J}\frac{\partial}{\partial \hat{z}}\left(\frac{K_\mathrm{H}}{J}\frac{\partial \chi}{\partial \hat{z}}\right)\]

is solved using the Crank-Nicolson scheme, where the system

\[a_k \chi_{k-1}^{n+1} + b_k \chi_k^{n+1} + c_k \chi_{k+1}^{n+1} = f_k\]

is solved using a Thomas tridiagonal solver, with $\mathcal{K}_\mathrm{H} = \frac{K_\mathrm{H}}{J}$ and

\[\begin{align*} a_k &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k} , \\ b_k &= 1 + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k}, \\ c_k &= -\frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} , \\ f_k &= \left[1 - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} - \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k}\right] \chi_k^{n} \\ &+ \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k+1/2}}{J_k} \chi_{k+1}^{n} + \frac{\Delta t}{2(\Delta \hat{z})^2}\frac{\mathcal{K}_{\mathrm{H},k-1/2}}{J_k} \chi_{k-1}^{n}. \end{align*}\]

Arguments

• state: Model state.

• dt: Time step.

• turbulence_scheme: General turbulence parameterization configuration.

• variable: Variable (equation) of choice.

• model: Dynamic equations.

• tracer_setup: General tracer-transport configuration.

See also

• PinCFlow.Update.check_tke!

• PinCFlow.Update.reset_thomas!

• PinCFlow.Update.thomas_algorithm!

• PinCFlow.Update.turbulence_diffusion_coefficient

update!(state::State, dt::AbstractFloat, m::Integer, variable::Rho)

Update the density if the atmosphere is not Boussinesq by dispatching to the appropriate method.

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::Rho,
    model::Val{:Boussinesq},
)

Return in Boussinesq mode (the density is constant).

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::Rho,
    model::Union{Val{:PseudoIncompressible}, Val{:Compressible}},
)

Update the density with a Runge-Kutta step on the left-hand side of the equation (the right-hand side is zero).

The update is given by

\[\begin{align*} q^\rho & \rightarrow - \frac{\Delta t}{J} \left(\frac{\mathcal{F}^{\rho, \hat{x}}_{i + 1 / 2} - \mathcal{F}^{\rho, \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho, \hat{y}}_{j + 1 / 2} - \mathcal{F}^{\rho, \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{\rho, \hat{z}}_{k + 1 / 2} - \mathcal{F}^{\rho, \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + \alpha_\mathrm{RK} q^\rho,\\ \rho & \rightarrow \rho + \beta_\mathrm{RK} q^\rho, \end{align*}\]

where $\Delta t$ is the time step given as input to this method.

update!(state::State, dt::AbstractFloat, m::Integer, variable::RhoP, side::LHS)

Update the density fluctuations with a Runge-Kutta step on the left-hand-side of the equation.

The update is given by

\[\begin{align*} q^{\rho'} & \rightarrow - \frac{\Delta t}{J} \left(\frac{\mathcal{F}^{\rho', \hat{x}}_{i + 1 / 2} - \mathcal{F}^{\rho', \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho', \hat{y}}_{j + 1 / 2} - \mathcal{F}^{\rho', \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{\rho', \hat{z}}_{k + 1 / 2} - \mathcal{F}^{\rho', \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + \alpha_\mathrm{RK} q^{\rho'},\\ \rho' & \rightarrow \rho' + \beta_\mathrm{RK} q^{\rho'} \end{align*}\]

in Boussinesq/pseudo-incompressible mode and

\[\begin{align*} q^{\rho'} & \rightarrow \Delta t \left[- \frac{1}{J} \left(\frac{\mathcal{F}^{\rho', \hat{x}}_{i + 1 / 2} - \mathcal{F}^{\rho', \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho', \hat{y}}_{j + 1 / 2} - \mathcal{F}^{\rho', \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{\rho', \hat{z}}_{k + 1 / 2} - \mathcal{F}^{\rho', \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + \frac{F^P}{\bar{\theta}}\right] + \alpha_\mathrm{RK} q^{\rho'},\\ \rho' & \rightarrow \rho' + \beta_\mathrm{RK} q^{\rho'} \end{align*}\]

in compressible mode.

update!(
    state::State,
    dt::AbstractFloat,
    variable::RhoP,
    side::RHS,
    integration::Explicit,
)

Update the density fluctuations with an explicit Euler step the on right-hand side of the equation, without the Rayleigh-damping term.

The update is given by

\[\rho' \rightarrow - \frac{\rho}{g} \left(b' - \Delta t N^2 \frac{\bar{\rho}}{\rho} w\right)\]

in Boussinesq/pseudo-incompressible mode and

\[\rho' \rightarrow - \frac{\rho}{g} \left[b' - \Delta t N^2 \frac{P / \bar{\theta}}{\rho} \left(\frac{W_{k + 1 / 2}}{\left(J P\right)_{k + 1 / 2}}\right)\right]\]

in compressible mode, where $b' = - g \rho' / \rho$.

update!(
    state::State,
    dt::AbstractFloat,
    variable::RhoP,
    side::RHS,
    integration::Implicit,
    rayleigh_factor::AbstractFloat,
)

Update the density fluctuations with an implicit Euler step on the right-hand side of the equation.

The update is given by

\[\begin{align*} \rho' & \rightarrow - \frac{\rho}{g} \left[1 + \beta_\mathrm{R} \Delta t + \frac{\bar{\rho}}{\rho} \left(N \Delta t\right)^2\right]^{- 1}\\ & \quad \times \left\{- \frac{\bar{\rho}}{\rho} N^2 \Delta t J \left[\hat{w}_\mathrm{old} + \Delta t \left(- \left(c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}}\right) + \left(\frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)\right)\right] + \left(1 + \beta_\mathrm{R} \Delta t\right) b'\right.\\ & \qquad \quad + \left.\frac{\bar{\rho}}{\rho} N^2 \Delta t J \left(1 + \beta_\mathrm{R} \Delta t\right) \left(G^{13} u + G^{23} v\right)\vphantom{\left[\left(\frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)\right]}\right\}, \end{align*}\]

in Boussinesq/pseudo-incompressible mode and

\[\begin{align*} \rho' & \rightarrow - \frac{\rho}{g} \left[1 + \beta_\mathrm{R} \Delta t + \frac{P / \bar{\theta}}{\rho} \left(N \Delta t\right)^2\right]^{- 1}\\ & \quad \times \left\{- \frac{P / \bar{\theta}}{\rho} N^2 \Delta t J \left[\left(\frac{\hat{W}_{\mathrm{old}, k + 1 / 2}}{\left(J P\right)_{k + 1 / 2}}\right) + \Delta t \left(- \left(c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}}\right) + \left(\frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)\right)\right]\right.\\ & \qquad \quad + \left(1 + \beta_\mathrm{R} \Delta t\right) b' + \frac{P / \bar{\theta}}{\rho} N^2 \Delta t J \left(1 + \beta_\mathrm{R} \Delta t\right)\\ & \qquad \quad \times \left.\left[G^{13} \left(\frac{U_{i + 1 / 2}}{\left(J P\right)_{i + 1 / 2}}\right) + G^{23} \left(\frac{V_{j + 1 / 2}}{\left(J P\right)_{j + 1 / 2}}\right)\right]\right\}, \end{align*}\]

in compressible mode, where $\hat{w}_\mathrm{old}$ is the transformed vertical wind stored in state.variables.backups.

update!(state::State, dt::AbstractFloat, m::Integer, variable::U, side::LHS)

Update the zonal momentum with a Runge-Kutta step on the left-hand side of the equation.

The update is given by

\[\begin{align*} q^{\rho u}_{i + 1 / 2} & \rightarrow \Delta t \left[- \frac{1}{J_{i + 1 / 2}} \left(\frac{\mathcal{F}^{\rho u, \hat{x}}_{i + 1} - \mathcal{F}^{\rho u, \hat{x}}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho u, \hat{y}}_{i + 1 / 2, j + 1 / 2} - \mathcal{F}^{\rho u, \hat{y}}_{i + 1 / 2, j - 1 / 2}}{\Delta \hat{y}}\right.\right.\\ & \qquad \qquad \qquad \qquad + \left.\left.\frac{\mathcal{F}^{\rho u, \hat{z}}_{i + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho u, \hat{z}}_{i + 1 / 2, k - 1 / 2}}{\Delta \hat{z}}\right) + f \left(\rho_\mathrm{old} v\right)_{i + 1 / 2}\right] + \alpha_\mathrm{RK} q^{\rho u}_{i + 1 / 2},\\ u_{i + 1 / 2} & \rightarrow \rho_{i + 1 / 2}^{- 1} \left(\rho_{\mathrm{old}, i + 1 / 2} u_{i + 1 / 2} + \beta_\mathrm{RK} q^{\rho u}_{i + 1 / 2}\right), \end{align*}\]

where $\rho_\mathrm{old}$ is the density stored in state.variables.backups.

update!(
    state::State,
    dt::AbstractFloat,
    variable::U,
    side::RHS,
    integration::Explicit,
)

Update the zonal wind with an explicit Euler step on the right-hand side of the equation, without the Rayleigh-damping term.

The update is given by

\[u_{i + 1 / 2} \rightarrow u_{i + 1 / 2} + \Delta t \left(- c_p \frac{P_{i + 1 / 2}}{\rho_{i + 1 / 2}} \mathcal{P}_{i + 1 / 2}^{\rho u} + \frac{F_{i + 1 / 2}^{\rho u}}{\rho_{i + 1 / 2}}\right)\]

in Boussinesq/pseudo-incompressible mode and

\[U_{i + 1 / 2} \rightarrow U_{i + 1 / 2} + \Delta t \left(J P\right)_{i + 1 / 2} \left(- c_p \frac{P_{i + 1 / 2}}{\rho_{i + 1 / 2}} \mathcal{P}_{i + 1 / 2}^{\rho u} + \frac{F_{i + 1 / 2}^{\rho u}}{\rho_{i + 1 / 2}}\right)\]

in compressible mode.

update!(
    state::State,
    dt::AbstractFloat,
    variable::U,
    side::RHS,
    integration::Implicit,
    rayleigh_factor::AbstractFloat,
)

Update the zonal wind with an implicit Euler step on the right-hand side of the equation.

The update is given by

\[u_{i + 1 / 2} \rightarrow \left(1 + \beta_{\mathrm{R}, i + 1 / 2} \Delta t\right)^{- 1} \left[u_{i + 1 / 2} + \Delta t \left(- c_p \frac{P_{i + 1 / 2}}{\rho_{i + 1 / 2}} \mathcal{P}_{i + 1 / 2}^{\rho u} + \frac{F_{i + 1 / 2}^{\rho u}}{\rho_{i + 1 / 2}}\right)\right]\]

in Boussinesq/pseudo-incompressible mode and

\[U_{i + 1 / 2} \rightarrow \left(1 + \beta_{\mathrm{R}, i + 1 / 2} \Delta t\right)^{- 1} \left[U_{i + 1 / 2} + \Delta t \left(J P\right)_{i + 1 / 2} \left(- c_p \frac{P_{i + 1 / 2}}{\rho_{i + 1 / 2}} \mathcal{P}_{i + 1 / 2}^{\rho u} + \frac{F_{i + 1 / 2}^{\rho u}}{\rho_{i + 1 / 2}}\right)\right]\]

in compressible mode.

update!(state::State, dt::AbstractFloat, m::Integer, variable::V, side::LHS)

Update the meridional momentum with a Runge-Kutta step on the left-hand side of the equation.

The update is given by

\[\begin{align*} q^{\rho v}_{j + 1 / 2} & \rightarrow \Delta t \left[- \frac{1}{J_{j + 1 / 2}} \left(\frac{\mathcal{F}^{\rho v, \hat{x}}_{i + 1 / 2, j + 1 / 2} - \mathcal{F}^{\rho v, \hat{x}}_{i - 1 / 2, j + 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho v, \hat{y}}_{j + 1} - \mathcal{F}^{\rho v, \hat{y}}}{\Delta \hat{y}}\right.\right.\\ & \qquad \qquad \qquad \qquad + \left.\left.\frac{\mathcal{F}^{\rho v, \hat{z}}_{j + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho v, \hat{z}}_{j + 1 / 2, k - 1 / 2}}{\Delta \hat{z}}\right) - f \left(\rho_\mathrm{old} u_\mathrm{old}\right)_{j + 1 / 2}\right] + \alpha_\mathrm{RK} q^{\rho v}_{j + 1 / 2},\\ v_{j + 1 / 2} & \rightarrow \rho_{j + 1 / 2}^{- 1} \left(\rho_{\mathrm{old}, j + 1 / 2} v_{j + 1 / 2} + \beta_\mathrm{RK} q^{\rho v}_{j + 1 / 2}\right), \end{align*}\]

where $\rho_\mathrm{old}$ and $u_{\mathrm{old}, i + 1 / 2}$ are the density and zonal wind stored in state.variables.backups.

update!(
    state::State,
    dt::AbstractFloat,
    variable::V,
    side::RHS,
    integration::Explicit,
)

Update the meridional wind with an explicit Euler step on the right-hand side of the equation, without the Rayleigh-damping term.

The update is given by

\[v_{i + 1 / 2} \rightarrow v_{j + 1 / 2} + \Delta t \left(- c_p \frac{P_{j + 1 / 2}}{\rho_{j + 1 / 2}} \mathcal{P}_{j + 1 / 2}^{\rho v} + \frac{F_{j + 1 / 2}^{\rho v}}{\rho_{j + 1 / 2}}\right)\]

in Boussinesq/pseudo-incompressible mode and

\[V_{j + 1 / 2} \rightarrow V_{j + 1 / 2} + \Delta t \left(J P\right)_{j + 1 / 2} \left(- c_p \frac{P_{j + 1 / 2}}{\rho_{j + 1 / 2}} \mathcal{P}_{j + 1 / 2}^{\rho v} + \frac{F_{j + 1 / 2}^{\rho v}}{\rho_{j + 1 / 2}}\right)\]

in compressible mode.

update!(
    state::State,
    dt::AbstractFloat,
    variable::V,
    side::RHS,
    integration::Implicit,
    rayleigh_factor::AbstractFloat,
)

Update the meridional wind with an implicit Euler step on the right-hand side of the equation.

The update is given by

\[v_{j + 1 / 2} \rightarrow \left(1 + \beta_{\mathrm{R}, j + 1 / 2} \Delta t\right)^{- 1} \left[v_{j + 1 / 2} + \Delta t \left(- c_p \frac{P_{j + 1 / 2}}{\rho_{j + 1 / 2}} \mathcal{P}_{j + 1 / 2}^{\rho v} + \frac{F_{j + 1 / 2}^{\rho v}}{\rho_{j + 1 / 2}}\right)\right]\]

in Boussinesq/pseudo-incompressible mode and

\[V_{j + 1 / 2} \rightarrow \left(1 + \beta_{\mathrm{R}, j + 1 / 2} \Delta t\right)^{- 1} \left[V_{j + 1 / 2} + \Delta t \left(J P\right)_{j + 1 / 2} \left(- c_p \frac{P_{j + 1 / 2}}{\rho_{j + 1 / 2}} \mathcal{P}_{j + 1 / 2}^{\rho v} + \frac{F_{j + 1 / 2}^{\rho v}}{\rho_{j + 1 / 2}}\right)\right]\]

in compressible mode.

update!(state::State, dt::AbstractFloat, m::Integer, variable::W, side::LHS)

Update the transformed vertical momentum with a Runge-Kutta step on the left-hand side of the equation.

The update is given by

\[\begin{align*} q^{\rho \hat{w}}_{k + 1 / 2} & \rightarrow \Delta t \left\{- \left[G^{13} \left(\frac{1}{J_{i + 1 / 2}} \left(\frac{\mathcal{F}^{\rho u, \hat{x}}_{i + 1} - \mathcal{F}^{\rho u, \hat{x}}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho u, \hat{y}}_{i + 1 / 2, j + 1 / 2} - \mathcal{F}^{\rho u, \hat{y}}_{i + 1 / 2, j - 1 / 2}}{\Delta \hat{y}}\right.\right.\right.\right.\\ & \qquad \qquad \qquad \qquad \qquad \qquad + \left.\left.\left.\frac{\mathcal{F}^{\rho u, \hat{z}}_{i + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho u, \hat{z}}_{i + 1 / 2, k - 1 / 2}}{\Delta \hat{z}}\right)\right)\right]_{k + 1 / 2}\\ & \qquad \qquad - \left[G^{23} \left(\frac{1}{J_{j + 1 / 2}} \left(\frac{\mathcal{F}^{\rho v, \hat{x}}_{i + 1 / 2, j + 1 / 2} - \mathcal{F}^{\rho v, \hat{x}}_{i - 1 / 2, j + 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho v, \hat{y}}_{j + 1} - \mathcal{F}^{\rho v, \hat{y}}}{\Delta \hat{y}}\right.\right.\right.\\ & \qquad \qquad \qquad \qquad \qquad \qquad + \left.\left.\left.\frac{\mathcal{F}^{\rho v, \hat{z}}_{j + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho v, \hat{z}}_{j + 1 / 2, k - 1 / 2}}{\Delta \hat{z}}\right)\right)\right]_{k + 1 / 2}\\ & \qquad \qquad - \frac{1}{J_{k + 1 / 2}^2} \left(\frac{\mathcal{F}^{\rho w, \hat{x}}_{i + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho w, \hat{x}}_{i - 1 / 2, k + 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho w, \hat{y}}_{j + 1 / 2, k + 1 / 2} - \mathcal{F}^{\rho w, \hat{y}}_{j - 1 / 2, k + 1 / 2}}{\Delta \hat{y}}\right.\\ & \qquad \qquad \qquad \qquad \quad + \left.\frac{\mathcal{F}^{\rho w, \hat{z}}_{k + 1} - \mathcal{F}^{\rho w, \hat{z}}}{\Delta \hat{z}}\right)\\ & \qquad \qquad + \left.G^{13} f \left(\rho_\mathrm{old} v_\mathrm{old}\right)_{k + 1 / 2} - G^{23} f \left(\rho_\mathrm{old} u_\mathrm{old}\right)_{k + 1 / 2}\vphantom{- \frac{1}{J^2} \left(\frac{\mathcal{F}^{\rho w, \hat{z}}_{k + 1} - \mathcal{F}^{\rho w, \hat{z}}}{\Delta \hat{z}}\right)}\right\} + \alpha_\mathrm{RK} q^{\rho \hat{w}}_{k + 1 / 2},\\ \hat{w}_{k + 1 / 2} & \rightarrow \rho_{k + 1 / 2}^{- 1} \left(\rho_{\mathrm{old}, k + 1 / 2} \hat{w}_{k + 1 / 2} + \beta_\mathrm{RK} q^{\rho \hat{w}}_{k + 1 / 2}\right), \end{align*}\]

where $\rho_\mathrm{old}$, $u_{\mathrm{old}, i + 1 / 2}$ and $v_{\mathrm{old}, j + 1 / 2}$ are the density, zonal wind and meridional wind stored in state.variables.backups.

update!(
    state::State,
    dt::AbstractFloat,
    variable::W,
    side::RHS,
    integration::Explicit,
)

Update the transformed vertical wind with an explicit Euler step on the right-hand side of the equation, without the Rayleigh-damping term.

The update is given by

\[\hat{w}_{k + 1 / 2} \rightarrow \hat{w}_{k + 1 / 2} + \Delta t \left[- c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}} + \left(\frac{b'_\mathrm{old}}{J}\right)_{k + 1 / 2} + \frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right]\]

in Boussinesq/pseudo-incompressible mode and

\[\hat{W}_{k + 1 / 2} \rightarrow \hat{W}_{k + 1 / 2} + \Delta t \left(J P\right)_{k + 1 / 2} \left[- c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}} + \left(\frac{b'_\mathrm{old}}{J}\right)_{k + 1 / 2} + \frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right]\]

in compressible mode, where $b'_\mathrm{old} = - g \rho'_\mathrm{old} / \rho$, with $\rho'_\mathrm{old}$ being the density fluctuations stored in state.variables.backups.

update!(
    state::State,
    dt::AbstractFloat,
    variable::W,
    side::RHS,
    integration::Implicit,
    rayleigh_factor::AbstractFloat,
)

Update the transformed vertical wind with an implicit Euler step on the right-hand side of the equation.

The update is given by

\[\begin{align*} \hat{w}_{k + 1 / 2} & \rightarrow \left[1 + \beta_{\mathrm{R}, k + 1 / 2} \Delta t + \frac{\bar{\rho}_{k + 1 / 2}}{\rho_{k + 1 / 2}} N^2_{k + 1 / 2} \left(\Delta t\right)^2\right]^{- 1}\\ & \quad \times \left\{\hat{w}_{k + 1 / 2} + \Delta t \left(- c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}} + \left(\frac{b'}{J}\right)_{k + 1 / 2} + \frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)\right.\\ & \qquad \quad + \left.\frac{\bar{\rho}_{k + 1 / 2}}{\rho_{k + 1 / 2}} N^2_{k + 1 / 2} \left(\Delta t\right)^2 \left[\left(G^{13} u\right)_{k + 1 / 2} + \left(G^{2 3} v\right)_{k + 1 / 2}\right]\vphantom{\left(\frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)}\right\} \end{align*}\]

in Boussinesq/pseudo-incompressible mode and

\[\begin{align*} \hat{W}_{k + 1 / 2} & \rightarrow \left[1 + \beta_{\mathrm{R}, k + 1 / 2} \Delta t + \frac{\left(P / \bar{\theta}\right)_{k + 1 / 2}}{\rho_{k + 1 / 2}} N^2_{k + 1 / 2} \left(\Delta t\right)^2\right]^{- 1}\\ & \quad \times \left\{\hat{W}_{k + 1 / 2} + \Delta t \left(J P\right)_{k + 1 / 2} \left(- c_p \frac{P_{k + 1 / 2}}{\rho_{k + 1 / 2}} \mathcal{P}_{k + 1 / 2}^{\rho \hat{w}} + \left(\frac{b'}{J}\right)_{k + 1 / 2} + \frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)\right.\\ & \qquad \quad + \left(J P\right)_{k + 1 / 2} \frac{\left(P / \bar{\theta}\right)_{k + 1 / 2}}{\rho_{k + 1 / 2}} N^2_{k + 1 / 2} \left(\Delta t\right)^2\\ & \qquad \quad \times \left.\left[\left(G^{13} \left(\frac{U_{i + 1 / 2}}{\left(J P\right)_{i + 1 / 2}}\right)\right)_{k + 1 / 2} + \left(G^{2 3} \left(\frac{V_{j + 1 / 2}}{\left(J P\right)_{j + 1 / 2}}\right)\right)_{k + 1 / 2}\right]\vphantom{\left(\frac{F_{k + 1 / 2}^{\rho \hat{w}}}{\rho_{k + 1 / 2}}\right)}\right\} \end{align*}\]

in compressible mode.

update!(state::State, dt::AbstractFloat, variable::PiP)

Update the Exner-pressure if the atmosphere is compressible by dispatching to the appropriate method.

update!(
    state::State,
    dt::AbstractFloat,
    variable::PiP,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)

Return in non-compressible modes.

update!(
    state::State,
    dt::AbstractFloat,
    variable::PiP,
    model::Val{:Compressible},
)

Update the Exner-pressure such that it is synchronized with the updated mass-weighted potential temperature.

The update is given by

\[\begin{align*} \pi' & \rightarrow \pi' + \Delta t \left(\frac{\partial \pi'}{\partial P}\right) \left[- \frac{1}{J} \left(\frac{U_{\mathrm{old}, i + 1 / 2} - U_{\mathrm{old}, i - 1 / 2}}{\Delta \hat{x}} + \frac{V_{\mathrm{old}, j + 1 / 2} - V_{\mathrm{old}, j - 1 / 2}}{\Delta \hat{y}}\right.\right.\\ & \qquad \qquad \qquad \qquad \qquad \qquad + \left.\left.\frac{\hat{W}_{\mathrm{old}, k + 1 / 2} - \hat{W}_{\mathrm{old}, k - 1 / 2}}{\Delta \hat{z}}\right) + F^P\right], \end{align*}\]

where $U_{\mathrm{old}, i + 1 / 2}$, $V_{\mathrm{old}, j + 1 / 2}$ and $\hat{W}_{\mathrm{old}, k + 1 / 2}$ are the transformed wind components (including the factor $J P$) stored in state.variables.backups.

update!(state::State, dt::AbstractFloat, m::Integer, variable::P)

Update the mass-weighted potential temperature if the atmosphere is compressible by dispatching to the appropriate method.

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::P,
    model::Union{Val{:Boussinesq}, Val{:PseudoIncompressible}},
)

Return in non-compressible modes.

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::P,
    model::Val{:Compressible},
)

Update the mass-weighted potential temperature with a Runge-Kutta step on the left-hand side of the equation (the right-hand side is zero).

The update is given by

\[\begin{align*} q^P & \rightarrow \Delta t \left[- \frac{1}{J} \left(\frac{\mathcal{F}^{P, \hat{x}}_{i + 1 / 2} - \mathcal{F}^{P, \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{P, \hat{y}}_{j + 1 / 2} - \mathcal{F}^{P, \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{P, \hat{z}}_{k + 1 / 2} - \mathcal{F}^{P, \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + F^P\right] + \alpha_\mathrm{RK} q^P,\\ P & \rightarrow P + \beta_\mathrm{RK} q^P. \end{align*}\]

update!(state::State, dt::AbstractFloat, m::Integer, variable::Chi)

Update the tracers by dispatching to the appropriate method.

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::Chi,
    tracer_setup::Val{:NoTracer},
)

Return for configurations without tracer transport.

update!(
    state::State,
    dt::AbstractFloat,
    m::Integer,
    variable::Chi,
    tracer_setup::Val{:TracerOn},
)

Update the tracers with a Runge-Kutta step on the left-hand sides of the equations with WKB right-hand side terms according to namelists configuration.

The update is given by

\[\begin{align*} q^{\rho \chi} & \rightarrow \Delta t \left[- \frac{1}{J} \left(\frac{\mathcal{F}^{\rho \chi, \hat{x}}_{i + 1 / 2} - \mathcal{F}^{\rho \chi, \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho \chi, \hat{y}}_{j + 1 / 2} - \mathcal{F}^{\rho \chi, \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{\rho \chi, \hat{z}}_{k + 1 / 2} - \mathcal{F}^{\rho \chi, \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + F^{\rho \chi}\right] + \alpha_\mathrm{RK} q^{\rho \chi},\\ \left(\rho \chi\right) & \rightarrow \left(\rho \chi\right) + \beta_\mathrm{RK} q^{\rho \chi}. \end{align*}\]

update!(state::State, dt::AbstractFloat, m::Integer, variable::TKE)

Update the turbulent kinetic energy with a Runge-Kutta step on the left-hand sides of the equations with shear and buoyancy production terms.

The update is given by

\[\begin{align*} q^{\rho e_\mathrm{k}} & \rightarrow \Delta t \left[- \frac{1}{J} \left(\frac{\mathcal{F}^{\rho e_\mathrm{k}, \hat{x}}_{i + 1 / 2} - \mathcal{F}^{\rho e_\mathrm{k}, \hat{x}}_{i - 1 / 2}}{\Delta \hat{x}} + \frac{\mathcal{F}^{\rho e_\mathrm{k}, \hat{y}}_{j + 1 / 2} - \mathcal{F}^{\rho e_\mathrm{k}, \hat{y}}_{j - 1 / 2}}{\Delta \hat{y}} + \frac{\mathcal{F}^{\rho e_\mathrm{k}, \hat{z}}_{k + 1 / 2} - \mathcal{F}^{\rho e_\mathrm{k}, \hat{z}}_{k - 1 / 2}}{\Delta \hat{z}}\right) + F^{\rho e_\mathrm{k}}\right] + \alpha_\mathrm{RK} q^{\rho e_\mathrm{k}},\\ \left(\rho e_\mathrm{k}\right) & \rightarrow \left(\rho e_\mathrm{k}\right) + \beta_\mathrm{RK} q^{\rho e_\mathrm{k}}. \end{align*}\]

Arguments

• state: Model state.

• dt: Time step.

• m: Runge-Kutta-stage index.

• variable: Variable to update.

• model: Dynamic equations.

• side: Side of the equation.

• integration: Type of the Euler step.

• rayleigh_factor: Factor by which the Rayleigh-damping coefficient is multiplied.

• tracer_setup: General tracer-transport configuration.

See also

• PinCFlow.Update.compute_volume_force

• PinCFlow.Update.compute_compressible_wind_factor

• PinCFlow.Update.compute_vertical_wind

• PinCFlow.Update.compute_buoyancy_factor

• PinCFlow.Update.compute_pressure_gradient

• PinCFlow.Update.transform

• PinCFlow.Update.conductive_heating