Coarse-grained kinetic energy equation
The CoarseGrainedKineticEnergyEquation module provides diagnostics for the kinetic energy budget of the coarse-grained (filtered) flow, in which a low-pass spatial filter $\widetilde{(\,\cdot\,)}$ separates a filtered from a subfilter scale. The section below derives that budget.
Deriving the coarse-grained kinetic energy budget
Oceananigans' NonhydrostaticModel evolves the velocity $v_i$ with the momentum equation (here without the background-flow, surface-wave, and forcing terms)
\[\partial_t v_i = - v_j \, \partial_j v_i - \epsilon_{ijk} \, f_j \, v_k - \partial_i p + b \, \hat g_i - \partial_j \tau_{ij} ,\]
with Coriolis parameter $f_i$, kinematic pressure $p$, buoyancy $b$ along the vertical $\hat g_i$, the permutation symbol $\epsilon_{ijk}$, and the kinematic stress tensor $\tau_{ij}$ supplied by the closure. The velocity components are $(v_1, v_2, v_3) = (u, v, w)$, and the resolved viscous dissipation $\varepsilon = -\partial_j v_i \, \tau_{ij}$ is what KineticEnergyDissipationRate computes. Incompressibility $\partial_i v_i = 0$ lets us write advection in flux form, $v_j \, \partial_j v_i = \partial_j (v_i v_j)$.
We can define the residual (subfilter) stress tensor
\[\tau^r_{ij} = \widetilde{v_i v_j} - \tilde v_i\,\tilde v_j ,\]
and re-write the filtered momentum equation as:
\[\partial_t \tilde v_i = - \tilde v_j \, \partial_j \tilde v_i - \partial_j \tau^r_{ij} - \epsilon_{ijk} \, f_j \, \tilde v_k - \partial_i \tilde p + \tilde b \, \hat g_i - \partial_j \tilde\tau_{ij} .\]
Multiplying by the filtered velocity $\tilde v_i$ gives the budget for the kinetic energy of the filtered flow $K^l = \tfrac{1}{2}\,\tilde v_i\,\tilde v_i$. Advection, pressure, and the two stress terms each split into a transport divergence plus a local term; Coriolis does no work ($\tilde v_i\,\epsilon_{ijk}\,f_j\,\tilde v_k = 0$), leaving
\[\partial_t K^l + \partial_j J_j = \underbrace{\tilde w\,\tilde b}_{\text{buoyancy production}} - \underbrace{\Pi_K}_{\text{cross-scale flux}} - \underbrace{\varepsilon^l}_{\text{dissipation}} ,\]
where $J_i$ collects the (advective, pressure, and stress) transport fluxes, which vanish when integrated over a closed or periodic domain. The two local exchange terms are
\[\Pi_K = -\tau^r_{ij}\,\widetilde{S}_{ij} , \qquad \varepsilon^l = -\partial_j \tilde v_i \, \tilde\tau_{ij} , \qquad \widetilde{S}_{ij} = \tfrac{1}{2}\left(\partial_j \tilde v_i + \partial_i \tilde v_j\right) .\]
The buoyancy production $\tilde w\,\tilde b$ converts between filtered kinetic and potential energy. $\Pi_K$ (KineticEnergyCrossScaleFlux) is the cross-scale kinetic energy flux: the rate at which the filter transfers kinetic energy from the filtered to the subfilter scales, following Aluie et al. (2018). A positive $\Pi_K$ denotes a forward (downscale) transfer, and $\widetilde{S}_{ij}$ is the strain rate tensor of the filtered velocity. The residual stress $\tau^r_{ij}$ itself is available as subfilter_stress_tensor.
$\varepsilon^l$ (CoarseGrainedKineticEnergyDissipationRate) is the viscous dissipation of the filtered flow: the filtered velocity gradient contracted with the filtered stress $\tilde\tau_{ij}$. The stress is filtered, not recomputed from $\tilde v_i$: $\widetilde{\tau_{ij}(v_i)}$ and $\tau_{ij}(\tilde v_i)$ agree only for a constant, uniform viscosity. When $\tau_{ij}$ is symmetric (as for an isotropic closure) the antisymmetric part of the gradient drops out and the dissipation can be written with the strain rate, $\varepsilon^l = -\tilde\tau_{ij}\,\widetilde{S}_{ij}$, reducing further to $2\nu\,\widetilde{S}_{ij}\widetilde{S}_{ij}$ for a constant-viscosity closure.
These diagnostics take a filter argument: any callable mapping a field to its low-pass-filtered counterpart, typically a GaussianFilter or BoxFilter. The directions the filter acts in (set inside filter) are independent of how each diagnostic contracts: the stress tensor and cross-scale flux take a dims argument selecting the directions they contract over — so you can filter horizontally yet contract the full 3D tensor — while the coarse-grained dissipation always forms the full viscous contraction.
Example
using Oceananigans, Oceanostics
grid = RectilinearGrid(size=(16, 16, 16), extent=(1, 1, 1), topology=(Periodic, Periodic, Bounded))
model = NonhydrostaticModel(grid; closure=ScalarDiffusivity(ν=1e-4))
ℓ = 0.2 # Gaussian filter scale (full width at half maximum) in all three directions
filter = GaussianFilter(; dims=(1, 2, 3), σ=ℓ / (2√(2log(2))), boundary=(left=0, right=0))
τ = subfilter_stress_tensor(model, filter) # the subfilter stress tensor components
Πₖ = KineticEnergyCrossScaleFlux(model, filter) # the cross-scale KE flux, at (Center, Center, Center)
εˡ = CoarseGrainedKineticEnergyDissipationRate(model, filter) # dissipation of the filtered flow
# equivalently, the convenience methods build the Gaussian filter from σ for you:
εˡ = CoarseGrainedKineticEnergyDissipationRate(model; σ=ℓ / (2√(2log(2))), boundary=(left=0, right=0))
# output
CoarseGrainedKineticEnergyDissipationRate KernelFunctionOperation at (Center, Center, Center)
├── grid: 16×16×16 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── kernel_function: coarse_grained_dissipation_rate_ccc (generic function with 1 method)
└── arguments: ("NamedTuple", "NamedTuple")
└── computes: coarse-grained kinetic energy dissipation rate εˡ = ∂ⱼūᵢ·F̄ᵢⱼSubfilter-scale stress tensor
Oceanostics.CoarseGrainedKineticEnergyEquation.subfilter_stress_tensor — Function
subfilter_stress_tensor(
model,
filter;
dims,
collocate_diagonals
)
Return the components of the subfilter-scale (SFS) stress tensor τ, the residual momentum flux that a low-pass filter removes from the filtered scales:
τⁱʲ = filter(uⁱuʲ) - ūⁱ ūʲ , ūⁱ = filter(uⁱ)(also called the sub-grid stress in LES, or the generalized central moment in the coarse-graining framework of Aluie et al., 2018, J. Phys. Oceanogr., doi:10.1175/JPO-D-17-0100.1). It is the quantity contracted with the filtered strain rate to form the cross-scale kinetic-energy flux — see KineticEnergyCrossScaleFlux.
filter is any callable that maps a field to its low-pass-filtered counterpart, e.g. a reusable GaussianFilter or BoxFilter:
using Oceananigans, Oceanostics
grid = RectilinearGrid(size=(4, 4, 4), extent=(1, 1, 1), topology=(Periodic, Periodic, Bounded))
model = NonhydrostaticModel(grid)
filter = GaussianFilter(; dims=(1, 2, 3), σ=0.1)
τ = subfilter_stress_tensor(model, filter)
keys(τ)
# output
(:τ₁₁, :τ₂₂, :τ₃₃, :τ₁₂, :τ₁₃, :τ₂₃)The result is a NamedTuple with the independent components, each living at the same staggered location as the corresponding StressTensor component; collocate_diagonals has the same meaning as there and is forwarded to it (use collocate_diagonals = true to put the diagonals at ccc, e.g. to form the subfilter kinetic energy ½(τ₁₁ + τ₂₂ + τ₃₃)). The filtered velocities ūⁱ and the filtered momentum fluxes filter(uⁱuʲ) are materialized as Fields internally (the filter's fast staged path only fires when wrapped directly in a Field), so each returned component is a lazy operation over those computed fields and recomputes correctly when written by an OutputWriter.
dims selects which spatial directions enter the tensor, exactly as in StressTensor: component τⁱʲ is kept only when both i and j are in dims, and only the velocities used by the kept components are filtered. The default dims = (1, 2, 3) returns the full tensor; dims = (1, 3) returns the x–z subset (τ₁₁, τ₃₃, τ₁₃).
A convenience method subfilter_stress_tensor(model; σ, dims, boundary, N, collocate_diagonals) builds the Gaussian filter for you from a standard deviation σ (a Gaussian of full width at half maximum ℓ has σ = ℓ / (2√(2 ln 2))).
Cross-scale kinetic energy flux
Oceanostics.CoarseGrainedKineticEnergyEquation.KineticEnergyCrossScaleFlux — Type
KineticEnergyCrossScaleFlux(model, filter; dims)
Return the cross-scale (scale-to-scale) kinetic-energy flux Πₖ, the rate at which a low-pass filter transfers kinetic energy from the filtered to the subfilter scales (Aluie et al., 2018, J. Phys. Oceanogr., doi:10.1175/JPO-D-17-0100.1):
Πₖ = -τⁱʲ S̄ⁱʲwhere τⁱʲ = filter(uⁱuʲ) - ūⁱūʲ is the subfilter-scale stress tensor (subfilter_stress_tensor) and S̄ⁱʲ = ½(∂ūⁱ/∂xʲ + ∂ūʲ/∂xⁱ) is the strain rate tensor of the filtered velocity (StrainRateTensor applied to ūⁱ). The contraction is evaluated at (Center, Center, Center), with off-diagonal components counted twice by symmetry. Πₖ > 0 is forward (downscale, filtered → subfilter) transfer. The result is per unit mass (units m² s⁻³); multiply by a reference density ρ₀ for a volumetric power.
filter is any callable mapping a field to its filtered counterpart, e.g. a reusable GaussianFilter:
using Oceananigans, Oceanostics
grid = RectilinearGrid(size=(4, 4, 4), extent=(1, 1, 1), topology=(Periodic, Periodic, Bounded))
model = NonhydrostaticModel(grid)
ℓ = 0.2 # filter scale (full width at half maximum)
filter = GaussianFilter(; dims=(1, 2, 3), σ=ℓ / (2√(2log(2))))
KineticEnergyCrossScaleFlux(model, filter)
# output
KineticEnergyCrossScaleFlux KernelFunctionOperation at (Center, Center, Center)
├── grid: 4×4×4 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── kernel_function: cross_scale_ke_flux_ccc (generic function with 1 method)
└── arguments: ("Oceananigans.AbstractOperations.UnaryOperation",)
└── computes: cross-scale kinetic energy flux Πₖ = -τⁱʲS̄ⁱʲThe returned object is a lazy operation over internally materialized filtered Fields, so it is ready for Field, Integral, and OutputWriters and recomputes as the simulation evolves.
dims selects which directions enter the tensors (only their i,j components are summed, and only the velocities they use are filtered): the default dims = (1, 2, 3) gives the full 3D flux, while dims = (1, 3) gives the 2D x–z flux Πₖ = -(τ₁₁S̄₁₁ + τ₃₃S̄₃₃ + 2τ₁₃S̄₁₃). The filter's own directions are set independently inside filter, so you can filter horizontally yet contract the full tensor.
A convenience method KineticEnergyCrossScaleFlux(model; σ, dims, boundary, N) builds the Gaussian filter for you from a standard deviation σ (with σ = ℓ / (2√(2 ln 2)) for a FWHM ℓ).
Coarse-grained kinetic energy dissipation
Oceanostics.CoarseGrainedKineticEnergyEquation.CoarseGrainedKineticEnergyDissipationRate — Type
CoarseGrainedKineticEnergyDissipationRate(model, filter)
Return the coarse-grained (filtered-flow) kinetic-energy dissipation rate εˡ, the rate at which viscosity removes kinetic energy from the filtered velocity field ūᵢ = filter(uᵢ):
εˡ = ∂ⱼūᵢ · F̄ᵢⱼ , F̄ᵢⱼ = filter(Fᵢⱼ(u))Here Fᵢⱼ(u) is the model's viscous momentum-flux tensor built from the full velocities and closure (the same fluxes KineticEnergyDissipationRate contracts), and F̄ᵢⱼ = filter(Fᵢⱼ(u)) is that flux low-pass filtered. Contracting the filtered flux with the filtered velocity gradient gives the viscous sink in the budget of the filtered kinetic energy Kˡ = ½ūᵢūᵢ (coarse-graining framework of Aluie et al., 2018, J. Phys. Oceanogr., doi:10.1175/JPO-D-17-0100.1).
Note the flux is filtered, filter(Fᵢⱼ(u)), not recomputed from the filtered velocity, Fᵢⱼ(ū). The two agree only for a constant, uniform viscosity, where the filter commutes with the flux; they differ once the viscosity varies in space (e.g. an eddy viscosity), and only the filtered-flux form is the dissipation that appears in the filtered KE budget. For a constant-viscosity ScalarDiffusivity it reduces to 2ν S̄ᵢⱼS̄ᵢⱼ, the dissipation of the resolved strain. It is evaluated at (Center, Center, Center), per unit mass (units m² s⁻³); multiply by a reference density ρ₀ for a volumetric power.
filter is any callable mapping a field to its filtered counterpart, e.g. a reusable GaussianFilter:
using Oceananigans, Oceanostics
grid = RectilinearGrid(size=(4, 4, 4), extent=(1, 1, 1), topology=(Periodic, Periodic, Bounded))
model = NonhydrostaticModel(grid; closure=ScalarDiffusivity(ν=1e-4))
ℓ = 0.2 # filter scale (full width at half maximum)
filter = GaussianFilter(; dims=(1, 2, 3), σ=ℓ / (2√(2log(2))))
CoarseGrainedKineticEnergyDissipationRate(model, filter)
# output
CoarseGrainedKineticEnergyDissipationRate KernelFunctionOperation at (Center, Center, Center)
├── grid: 4×4×4 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── kernel_function: coarse_grained_dissipation_rate_ccc (generic function with 1 method)
└── arguments: ("NamedTuple", "NamedTuple")
└── computes: coarse-grained kinetic energy dissipation rate εˡ = ∂ⱼūᵢ·F̄ᵢⱼThe viscosity and fluxes come from model.closure/model.closure_fields, exactly as in KineticEnergyDissipationRate, so the model needs a closure whose viscous fluxes are defined. The filtered velocities and the filtered fluxes are materialized as Fields internally (and refreshed on recompute), so the returned object is a lazy operation ready for Field, Integral, and OutputWriters and recomputes as the simulation evolves.
Unlike the cross-scale flux and the stress tensor, this diagnostic takes no dims argument: it always forms the full viscous contraction (matching KineticEnergyDissipationRate). The directions the filter acts in are set inside filter.
A convenience method CoarseGrainedKineticEnergyDissipationRate(model; σ, dims, boundary, N) builds the Gaussian filter for you from a standard deviation σ (with σ = ℓ / (2√(2 ln 2)) for a FWHM ℓ); dims here selects the directions the Gaussian filter acts in.