Large Eddy Simulation (LES) Models

This page provides a comprehensive introduction to Large Eddy Simulation (LES) and the subgrid-scale (SGS) closure models available in Tarang.jl.

What is Large Eddy Simulation?

The Turbulence Challenge

Turbulent flows contain a vast range of length scales, from the largest energy-containing eddies down to the smallest dissipative scales (Kolmogorov scales). In a Direct Numerical Simulation (DNS), we resolve all these scales, which requires:

\[N \sim Re^{9/4}\]

grid points in 3D, where $Re$ is the Reynolds number. For atmospheric or oceanic flows with $Re \sim 10^9$, this is computationally impossible.

The LES Approach

Large Eddy Simulation offers a practical alternative:

1. Resolve the large, energy-containing eddies directly
2. Model the effect of small, unresolved eddies

┌─────────────────────────────────────────────────────────────┐
│                    Energy Spectrum E(k)                      │
│                                                              │
│    E(k)                                                      │
│     │                                                        │
│     │   ╭──╮                                                 │
│     │  ╱    ╲          k^(-5/3)                              │
│     │ ╱      ╲___                                            │
│     │╱            ╲___                                       │
│     │                  ╲___                                  │
│     │                       ╲___                             │
│     └──────────────────────────────────────────────► k       │
│         │                    │                               │
│         │◄── Resolved ──────►│◄── Modeled (SGS) ──►│        │
│         │   (large eddies)   │   (small eddies)    │        │
│                              │                               │
│                          filter cutoff (Δ)                   │
└─────────────────────────────────────────────────────────────┘

The filter width $\Delta$ (typically the grid spacing) separates resolved from unresolved scales.

Filtering the Navier-Stokes Equations

We apply a spatial filter to decompose velocity into resolved ($\bar{u}$) and subgrid ($u'$) parts:

\[u = \bar{u} + u'\]

Filtering the incompressible Navier-Stokes equations gives:

\[\frac{\partial \bar{u}_i}{\partial t} + \bar{u}_j \frac{\partial \bar{u}_i}{\partial x_j} = -\frac{1}{\rho}\frac{\partial \bar{p}}{\partial x_i} + \nu \nabla^2 \bar{u}_i - \frac{\partial \tau_{ij}}{\partial x_j}\]

where the subgrid-scale (SGS) stress tensor appears:

\[\tau_{ij} = \overline{u_i u_j} - \bar{u}_i \bar{u}_j\]

This tensor represents the effect of unresolved turbulent motions on the resolved flow. Since we cannot compute $\tau_{ij}$ directly (it involves unresolved velocities), we must model it.

The Closure Problem

Why We Need Models

The SGS stress $\tau_{ij}$ contains information about scales we don't resolve. This creates the closure problem: our filtered equations have more unknowns than equations.

The Eddy Viscosity Hypothesis

Most SGS models use the Boussinesq hypothesis, which assumes the SGS stress is proportional to the resolved strain rate:

\[\tau_{ij} - \frac{1}{3}\tau_{kk}\delta_{ij} = -2 \nu_e \bar{S}_{ij}\]

where:

• \[\nu_e\]

  is the eddy viscosity (to be modeled)

• \[\bar{S}_{ij}\]

  is the resolved strain rate tensor:

\[\bar{S}_{ij} = \frac{1}{2}\left(\frac{\partial \bar{u}_i}{\partial x_j} + \frac{\partial \bar{u}_j}{\partial x_i}\right)\]

This transforms the filtered momentum equation into:

\[\frac{\partial \bar{u}_i}{\partial t} + \bar{u}_j \frac{\partial \bar{u}_i}{\partial x_j} = -\frac{1}{\rho}\frac{\partial \bar{p}^*}{\partial x_i} + (\nu + \nu_e) \nabla^2 \bar{u}_i\]

where we've absorbed the isotropic part into a modified pressure $\bar{p}^*$.

Key insight: The SGS model effectively adds a spatially-varying viscosity $\nu_e(x,t)$ to the molecular viscosity $\nu$.

Energy Cascade and Dissipation

Forward Energy Cascade

In 3D turbulence, energy flows from large scales to small scales (forward cascade):

Large eddies ──► Medium eddies ──► Small eddies ──► Dissipation
  (production)                                        (ε = 2νSᵢⱼSᵢⱼ)

The SGS model must drain energy from resolved scales at the correct rate to maintain physical behavior.

SGS Dissipation

The rate at which energy is transferred from resolved to unresolved scales is:

\[\varepsilon_{sgs} = -\tau_{ij} \bar{S}_{ij} = 2\nu_e \bar{S}_{ij}\bar{S}_{ij} = \nu_e |\bar{S}|^2\]

where $|\bar{S}| = \sqrt{2\bar{S}_{ij}\bar{S}_{ij}}$ is the strain rate magnitude.

There is no extra factor of 2 in the final expression: the $\sqrt{2}$ in the definition of $|\bar{S}|$ already carries it. This is what sgs_dissipation returns — writing $2\nu_e|\bar{S}|^2$ would double-count the dissipation.

A good SGS model ensures $\varepsilon_{sgs}$ matches the actual energy transfer rate.

Available Models in Tarang.jl

Tarang.jl provides two eddy viscosity models:

Smagorinsky Model

Physical Basis

Joseph Smagorinsky (1963) proposed the first and most widely-used SGS model. It's based on dimensional analysis and mixing length theory.

Key assumptions:

1. The subgrid scales are in equilibrium (production = dissipation)
2. The characteristic length scale is proportional to the filter width $\Delta$
3. The characteristic velocity scale is $\Delta |\bar{S}|$

Mathematical Formulation

The eddy viscosity is computed as:

\[\nu_e = (C_s \Delta)^2 |\bar{S}|\]

where:

• \[C_s \approx 0.17\]

  is the Smagorinsky constant

• \[\Delta = (\Delta_x \Delta_y \Delta_z)^{1/3}\]

  is the effective filter width

• \[|\bar{S}| = \sqrt{2 \bar{S}_{ij} \bar{S}_{ij}}\]

  is the strain rate magnitude

Expanded form for 3D:

\[|\bar{S}| = \sqrt{2\left[\left(\frac{\partial \bar{u}}{\partial x}\right)^2 + \left(\frac{\partial \bar{v}}{\partial y}\right)^2 + \left(\frac{\partial \bar{w}}{\partial z}\right)^2\right] + \left(\frac{\partial \bar{u}}{\partial y} + \frac{\partial \bar{v}}{\partial x}\right)^2 + \left(\frac{\partial \bar{u}}{\partial z} + \frac{\partial \bar{w}}{\partial x}\right)^2 + \left(\frac{\partial \bar{v}}{\partial z} + \frac{\partial \bar{w}}{\partial y}\right)^2}\]

Usage in Tarang.jl

using Tarang

# Grid parameters
N = 128                     # Grid points per direction
L = 2π                      # Domain size
Δ = L / N                   # Grid spacing

# Create the Smagorinsky model
sgs_model = SmagorinskyModel(
    C_s = 0.17,                    # Smagorinsky constant
    filter_width = (Δ, Δ, Δ),      # Grid spacing in each direction
    field_size = (N, N, N)         # Number of grid points
)

# Compute eddy viscosity from velocity gradients
# You need all 9 components of the velocity gradient tensor
compute_eddy_viscosity!(sgs_model,
    ∂u∂x, ∂u∂y, ∂u∂z,    # Gradients of u
    ∂v∂x, ∂v∂y, ∂v∂z,    # Gradients of v
    ∂w∂x, ∂w∂y, ∂w∂z     # Gradients of w
)

# Retrieve the computed eddy viscosity field
νₑ = get_eddy_viscosity(sgs_model)

Choosing the Smagorinsky Constant

The constant $C_s$ is not universal and depends on the flow:

Rule of thumb: Start with $C_s = 0.17$ and reduce if you observe over-dissipation.

Limitations

1. Over-dissipation: The model is always "on", even in laminar regions
2. No backscatter: Cannot represent energy transfer from small to large scales
3. Wall behavior: Needs damping functions near solid walls
4. Isotropic assumption: Assumes same behavior in all directions

Anisotropic Minimum Dissipation (AMD) Model

Motivation

The AMD model (Rozema et al., 2015) addresses key limitations of Smagorinsky:

1. Automatic switch-off: $\nu_e = 0$ in laminar/transitional regions
2. Anisotropic grids: Properly handles $\Delta_x \neq \Delta_y \neq \Delta_z$
3. Minimum dissipation: Adds only the dissipation needed for stability

Mathematical Formulation

The AMD model computes eddy viscosity as:

\[\nu_e = \max\left(0, -C \frac{\Delta_k^2 \frac{\partial \bar{u}_i}{\partial x_k} \frac{\partial \bar{u}_j}{\partial x_k} \bar{S}_{ij}}{\frac{\partial \bar{u}_m}{\partial x_n} \frac{\partial \bar{u}_m}{\partial x_n}}\right)\]

with summation over $i$, $j$ and $k$. The filter width $\Delta_k$ is indexed by the derivative direction $k$, not by the velocity-component index — that is what makes the model anisotropy-aware, and it is the term most easily got wrong. Tying $\Delta$ to the component index instead (e.g. $\hat{\delta}_{ij} = \Delta_i \delta_{ij}$) leaves only the diagonal strain contributing, and gives a different sign as well as a different magnitude.

Simplified form:

\[\nu_e = \max\left(0, -C \frac{\text{numerator}}{\text{denominator}}\right)\]

• Numerator: Measures alignment between velocity gradients and strain rate
• Denominator: Total velocity gradient magnitude squared
• max(0, ...): Ensures non-negative eddy viscosity

Physical Interpretation

The numerator can be negative (giving $\nu_e > 0$) when:

• Velocity gradients are aligned with the strain rate
• Energy is being transferred to smaller scales

The numerator is positive (giving $\nu_e = 0$) when:

• Flow is laminar or transitional
• No SGS dissipation is needed

Usage in Tarang.jl

# Create the AMD model
sgs_model = AMDModel(
    C = 1/12,                      # Poincaré constant (for spectral methods)
    filter_width = (Δx, Δy, Δz),   # Can be different in each direction
    field_size = (Nx, Ny, Nz),
    clip_negative = true           # Ensure νₑ ≥ 0 (recommended)
)

# Compute eddy viscosity (same interface as Smagorinsky)
compute_eddy_viscosity!(sgs_model,
    ∂u∂x, ∂u∂y, ∂u∂z,
    ∂v∂x, ∂v∂y, ∂v∂z,
    ∂w∂x, ∂w∂y, ∂w∂z
)

νₑ = get_eddy_viscosity(sgs_model)

Choosing the AMD Constant

The constant $C$ depends on the numerical discretization:

Scalar Transport (Buoyancy)

For buoyancy-driven flows, AMD also provides eddy diffusivity for scalar transport:

\[\kappa_e = \max\left(0, -C \frac{\Delta_k^2 \frac{\partial \bar{u}_i}{\partial x_k} \frac{\partial \bar{\theta}}{\partial x_k} \frac{\partial \bar{\theta}}{\partial x_i}}{\frac{\partial \bar{\theta}}{\partial x_n} \frac{\partial \bar{\theta}}{\partial x_n}}\right)\]

The inner sum runs over every velocity component $u_i$, so the method needs the complete velocity-gradient tensor (9 components in 3D, 4 in 2D) before the scalar gradients — not the gradient of a single velocity component:

# Compute eddy diffusivity for buoyancy/temperature
compute_eddy_diffusivity!(sgs_model,
    ∂u∂x, ∂u∂y, ∂u∂z,    # Gradients of u
    ∂v∂x, ∂v∂y, ∂v∂z,    # Gradients of v
    ∂w∂x, ∂w∂y, ∂w∂z,    # Gradients of w
    ∂b∂x, ∂b∂y, ∂b∂z     # Buoyancy/scalar gradients
)

κₑ = get_eddy_diffusivity(sgs_model)

Choosing Between Models

Decision Flowchart

Start
  │
  ▼
Is your grid anisotropic (Δx ≠ Δy ≠ Δz)?
  │
  ├── Yes ──► Use AMD (handles anisotropy naturally)
  │
  └── No
       │
       ▼
Does your flow have laminar or transitional regions?
       │
       ├── Yes ──► Use AMD (automatically switches off)
       │
       └── No
            │
            ▼
Is computational cost a primary concern?
            │
            ├── Yes ──► Use Smagorinsky (slightly cheaper)
            │
            └── No ──► Either works; AMD is generally more accurate

Summary Comparison

Complete Example: LES of Decaying Turbulence

This example demonstrates the SGS side of an LES of decaying homogeneous isotropic turbulence: building the model, evaluating the velocity gradients, and refreshing νₑ each step. Advancing the momentum equation is not shown, because Tarang does not do it for you — see the scope note at the top of this page.

using Tarang
using Statistics

# ============================================================
# 1. Physical and Numerical Parameters
# ============================================================

N = 128                     # Grid points per direction
L = 2π                      # Domain size [m]
Δ = L / N                   # Grid spacing [m]
ν = 1e-4                    # Molecular (kinematic) viscosity [m²/s]
dt = 0.001                  # Time step [s]
nsteps = 1000               # Number of time steps

# ============================================================
# 2. Create the SGS Model
# ============================================================

# Option A: Smagorinsky (simple, robust)
sgs = SmagorinskyModel(
    C_s = 0.17,
    filter_width = (Δ, Δ, Δ),
    field_size = (N, N, N),
    architecture = CPU()       # GPU() keeps νₑ on the device
)

# Option B: AMD (recommended for most applications)
# sgs = AMDModel(
#     C = 1/12,
#     filter_width = (Δ, Δ, Δ),
#     field_size = (N, N, N)
# )

# ============================================================
# 3. Setup Computational Domain
# ============================================================

coords = CartesianCoordinates("x", "y", "z")
dist = Distributor(coords; mesh=(1, 1, 1))  # Single processor

# Fourier bases for periodic domain
xbasis = RealFourier(coords["x"]; size=N, bounds=(0.0, L))
ybasis = RealFourier(coords["y"]; size=N, bounds=(0.0, L))
zbasis = RealFourier(coords["z"]; size=N, bounds=(0.0, L))

domain = Domain(dist, (xbasis, ybasis, zbasis))

# Create velocity field
u = VectorField(dist, coords, "u", (xbasis, ybasis, zbasis))

# ============================================================
# 4. Initialize with a Turbulent Velocity Field
# ============================================================

# A Taylor-Green vortex keeps this example self-contained; substitute your
# own IC (e.g. random phases with a prescribed energy spectrum).
x, y, z = local_grids(dist, xbasis, ybasis, zbasis)
ensure_layout!(u, :g)
get_grid_data(u.components[1]) .=  sin.(x) .* cos.(y') .* cos.(reshape(z, 1, 1, :))
get_grid_data(u.components[2]) .= .-cos.(x) .* sin.(y') .* cos.(reshape(z, 1, 1, :))
get_grid_data(u.components[3]) .= 0.0
ensure_layout!(u, :c)

# ============================================================
# 5. Helper Function: Compute All Velocity Gradients
# ============================================================

function compute_velocity_gradients(u)
    # grad(u) is a lazy Gradient — it is not itself iterable, so it cannot be
    # destructured directly. evaluate() turns it into a TensorField whose
    # components are already in the component-major order the models expect:
    # (∂u∂x, ∂u∂y, ∂u∂z, ∂v∂x, ∂v∂y, ∂v∂z, ∂w∂x, ∂w∂y, ∂w∂z).
    G = evaluate(grad(u))
    return Tuple(get_grid_data(g) for g in G.components)
end

# ============================================================
# 6. Per-Step SGS Update
# ============================================================
#
# Nothing in the solver stack calls the SGS model for you. What follows is the
# SGS half of a timestep; advancing ū is code you supply.

for step in 1:nsteps
    # --- Step 1: Compute velocity gradients ---
    grads = compute_velocity_gradients(u)

    # --- Step 2: Update SGS eddy viscosity ---
    compute_eddy_viscosity!(sgs, grads...)

    # --- Step 3: Get effective viscosity ---
    νₑ = get_eddy_viscosity(sgs)
    ν_eff = ν .+ νₑ  # Total viscosity = molecular + SGS, one value per grid point

    # --- Step 4: Advance the momentum equation (your code) ---
    # The filtered Navier-Stokes with SGS model:
    #   ∂ū/∂t + (ū·∇)ū = -∇p̄/ρ + ∇·((ν + νₑ)∇ū)
    #
    # ν_eff is a plain array and no timestepper keyword accepts it, so the
    # variable-viscosity term must be assembled and applied explicitly.

    # --- Step 5: Diagnostics ---
    if step % 100 == 0
        mean_νₑ = mean_eddy_viscosity(sgs)
        max_νₑ = max_eddy_viscosity(sgs)

        println("Step $step:")
        println("  Mean eddy viscosity: $(mean_νₑ)")
        println("  Max eddy viscosity:  $(max_νₑ)")
        println("  Ratio νₑ/ν (mean):   $(mean_νₑ/ν)")
    end
end

Complete Example: LES of Rayleigh-Bénard Convection

For buoyancy-driven turbulence we need both momentum and scalar closures. As above, only the SGS update is shown — the momentum and temperature updates are yours to write.

using Tarang

# ============================================================
# Physical Parameters
# ============================================================

Nx, Nz = 256, 128           # Grid resolution
Lx, Lz = 4.0, 1.0           # Domain size

Ra = 1e8                    # Rayleigh number
Pr = 1.0                    # Prandtl number
ν = sqrt(Pr / Ra)           # Kinematic viscosity
κ = ν / Pr                  # Thermal diffusivity
nsteps = 1000               # Number of time steps

# ============================================================
# Domain and Fields
# ============================================================

coords = CartesianCoordinates("x", "z")
dist   = Distributor(coords; dtype=Float64, device=CPU())
xbasis = RealFourier(coords["x"]; size=Nx, bounds=(0.0, Lx))  # periodic
zbasis = ChebyshevT(coords["z"];  size=Nz, bounds=(0.0, Lz))  # wall-bounded
domain = Domain(dist, (xbasis, zbasis))

u = Field(dist; name="u", bases=(xbasis, zbasis))   # horizontal velocity
w = Field(dist; name="w", bases=(xbasis, zbasis))   # vertical velocity
T = Field(dist; name="T", bases=(xbasis, zbasis))   # temperature

# (Initialise u, w and T with your own IC and boundary conditions.)

# ============================================================
# Create AMD Model (recommended for RBC)
# ============================================================

# One filter width per axis. `grid_spacing` reports the *smallest* spacing on
# each axis, which is the conservative choice on the Chebyshev-clustered z axis.
Δx, Δz = grid_spacing(domain)

sgs = AMDModel(
    C = 1/12,
    filter_width = (Δx, Δz),
    field_size = (Nx, Nz)
)

# ============================================================
# Per-Step SGS Update
# ============================================================

for step in 1:nsteps
    # Velocity gradients. grad() is lazy and not iterable; evaluate() realises
    # it, and for a scalar field the result has one component per axis.
    ∂u∂x, ∂u∂z = (get_grid_data(g) for g in evaluate(grad(u)).components)
    ∂w∂x, ∂w∂z = (get_grid_data(g) for g in evaluate(grad(w)).components)

    # Temperature gradients
    ∂T∂x, ∂T∂z = (get_grid_data(g) for g in evaluate(grad(T)).components)

    # --- Momentum closure ---
    # 2D needs the complete 2×2 velocity-gradient tensor, component-major.
    compute_eddy_viscosity!(sgs, ∂u∂x, ∂u∂z, ∂w∂x, ∂w∂z)
    νₑ = get_eddy_viscosity(sgs)
    ν_eff = ν .+ νₑ

    # --- Scalar (temperature) closure ---
    # The same velocity-gradient tensor first, then the scalar gradients.
    compute_eddy_diffusivity!(sgs, ∂u∂x, ∂u∂z, ∂w∂x, ∂w∂z, ∂T∂x, ∂T∂z)
    κₑ = get_eddy_diffusivity(sgs)
    κ_eff = κ .+ κₑ

    # ν_eff and κ_eff are plain arrays. Applying them to the momentum and
    # temperature equations is your code — see the scope note at the top.
end

Diagnostics and Analysis

SGS Energy Dissipation

The SGS dissipation rate tells you how much energy is being drained by the model:

# Get strain magnitude (computed during eddy viscosity calculation)
S_mag = sgs.strain_magnitude  # Only for Smagorinsky

# Or compute it yourself — note the broadcasting dots (sqrt is not vectorised)
S_mag = sqrt.(2 .* (S11.^2 .+ S22.^2 .+ S33.^2 .+ 2 .* S12.^2 .+ 2 .* S13.^2 .+ 2 .* S23.^2))

# SGS dissipation field: ε_sgs = νₑ |S̄|²  (no extra factor of 2 — |S̄| carries it)
ε_sgs = sgs_dissipation(sgs, S_mag)

# Domain-averaged SGS dissipation
ε_sgs_mean = mean_sgs_dissipation(sgs, S_mag)

Monitoring Model Behavior

# Check ratio of SGS to molecular viscosity
νₑ_mean = mean_eddy_viscosity(sgs)
println("νₑ/ν = $(νₑ_mean/ν)")

# If νₑ/ν >> 1: SGS model is dominant (typical for high Re)
# If νₑ/ν << 1: Flow is nearly DNS-resolved
# If νₑ/ν ~ 1: Well-resolved LES

# For AMD: check how often νₑ = 0
νₑ = get_eddy_viscosity(sgs)
fraction_zero = sum(νₑ .== 0) / length(νₑ)
println("Fraction with νₑ = 0: $(fraction_zero)")
# High fraction → flow is largely laminar/transitional

Tips for Successful LES

Resolution Requirements

LES still requires adequate resolution:

┌────────────────────────────────────────────────────────────┐
│  Resolution Quality for LES                                 │
├────────────────────────────────────────────────────────────┤
│  80% of TKE resolved      →  Minimum acceptable LES        │
│  90% of TKE resolved      →  Good quality LES              │
│  95%+ of TKE resolved     →  Nearly DNS quality            │
└────────────────────────────────────────────────────────────┘

Rule of thumb: Grid spacing should resolve the inertial range, typically $\Delta \lesssim L_{integral}/10$.

Common Pitfalls

When LES May Not Be Appropriate

• Very low Reynolds numbers: DNS may be feasible
• Strongly anisotropic turbulence: May need special treatment
• Flows with strong backscatter: Standard models don't capture this
• Near-wall regions: May need wall models or finer grids

2D Flows

Both models support 2D simulations:

# 2D Smagorinsky
sgs_2d = SmagorinskyModel(
    C_s = 0.17,
    filter_width = (Δx, Δy),
    field_size = (Nx, Ny)
)

# 2D velocity gradients (4 components)
compute_eddy_viscosity!(sgs_2d, ∂u∂x, ∂u∂y, ∂v∂x, ∂v∂y)

Note on 2D turbulence: 2D turbulence has an inverse energy cascade (energy flows to large scales), which is fundamentally different from 3D. Standard SGS models may not be appropriate for 2D flows.

Choosing an Architecture (CPU / GPU)

Both models take an architecture keyword, defaulting to CPU(). It decides where the model's internal buffers (eddy_viscosity, strain_magnitude, and for AMD eddy_diffusivity) are allocated: plain Arrays on CPU(), CuArrays on GPU().

# CPU (the default)
sgs_cpu = SmagorinskyModel(
    C_s = 0.17,
    filter_width = (Δ, Δ, Δ),
    field_size = (N, N, N),
    architecture = CPU()
)

# AMD takes the same keyword
amd_cpu = AMDModel(
    C = 1/12,
    filter_width = (Δx, Δy, Δz),
    field_size = (Nx, Ny, Nz),
    architecture = CPU()
)

GPU() requires CUDA.jl to be loaded first, and a CUDA device to be present — it throws otherwise:

using CUDA        # must be loaded before GPU() can be constructed
using Tarang

sgs_gpu = SmagorinskyModel(
    C_s = 0.17,
    filter_width = (Δ, Δ, Δ),
    field_size = (N, N, N),
    architecture = GPU()
)

compute_eddy_viscosity! and compute_eddy_diffusivity! dispatch on the model's architecture: broadcast kernels on the GPU, SIMD loops on the CPU. Gradient arrays that are not already on the model's architecture are copied there on every call, so a CPU array handed to a GPU() model works but pays a host-to-device transfer each step — keep the gradients on the same device as the model. get_eddy_viscosity returns an array of the model's own kind (Array on CPU(), CuArray on GPU()).

API Reference

Constructors

SmagorinskyModel(;
    C_s = 0.17,                    # Smagorinsky constant
    filter_width::NTuple{N, Real}, # (Δx, Δy) or (Δx, Δy, Δz)
    field_size::NTuple{N, Int},    # (Nx, Ny) or (Nx, Ny, Nz)
    dtype = Float64,               # Precision
    architecture = CPU()           # CPU() or GPU()
)

AMDModel(;
    C = 1/12,                      # Poincaré constant
    filter_width::NTuple{N, Real}, # Can be anisotropic
    field_size::NTuple{N, Int},
    clip_negative = true,          # Ensure νₑ ≥ 0
    dtype = Float64,               # Precision
    architecture = CPU()           # CPU() or GPU()
)

Core Functions

Analysis Functions

Utility Functions

Troubleshooting

References

Original Papers

1. Smagorinsky, J. (1963). "General circulation experiments with the primitive equations: I. The basic experiment." Monthly Weather Review, 91(3), 99-164.

2. Rozema, W., Bae, H. J., Moin, P., & Verstappen, R. (2015). "Minimum-dissipation models for large-eddy simulation." Physics of Fluids, 27(8), 085107. PDF

3. Abkar, M., Bae, H. J., & Moin, P. (2016). "Minimum-dissipation scalar transport model for large-eddy simulation of turbulent flows." Physical Review Fluids, 1(4), 041701. Link

Textbooks and Reviews

1. Pope, S. B. (2000). Turbulent Flows. Cambridge University Press. (Chapter 13: Large-Eddy Simulation)

2. Sagaut, P. (2006). Large Eddy Simulation for Incompressible Flows. Springer.

Related Documentation

• dedaLES AMD Documentation

See Also

• Solvers - Time integration methods
• Operators - Gradient and differential operators
• Stochastic Forcing - Forcing for turbulence
• API: LES Models - Complete API reference