Boundary Conditions

This page documents the boundary condition options for all field types in GeoDynamo.jl. Proper boundary conditions are essential for physically meaningful simulations of planetary core dynamics.

Overview

GeoDynamo.jl supports boundary conditions at two radial surfaces:

• Inner boundary (ICB): Inner core boundary or center (for ball geometry)
• Outer boundary (CMB): Core-mantle boundary or outer surface

Each field type (velocity, magnetic, temperature, composition) has specific BC options appropriate for its physics.

Specifying Boundary Conditions

Boundary conditions are passed to GeodynamoModel per field, wrapped in FieldBoundaryConditions(inner = ..., outer = ...):

t_bcs = FieldBoundaryConditions(
    inner = ValueBoundaryCondition(1.0),   # Dirichlet: fixed value
    outer = FluxBoundaryCondition(0.0),    # Neumann: fixed flux
)

BC type names

The canonical names follow Oceananigans.jl; the original GeoDynamo names remain as aliases and keep working:

Two equivalent ways to pass them

Per-field keywords:

model = GeodynamoModel(grid;
    velocity_bcs    = FieldBoundaryConditions(inner = NoSlip(), outer = NoSlip()),
    temperature_bcs = FieldBoundaryConditions(inner = ValueBoundaryCondition(1.0),
                                              outer = ValueBoundaryCondition(0.0)),
)

or one boundary_conditions NamedTuple (Oceananigans style):

model = GeodynamoModel(grid;
    boundary_conditions = (
        velocity    = FieldBoundaryConditions(inner = NoSlip(), outer = NoSlip()),
        temperature = FieldBoundaryConditions(inner = ValueBoundaryCondition(1.0),
                                              outer = ValueBoundaryCondition(0.0)),
    ))

Specifying the same field through both paths throws an ArgumentError.

Defaults

Velocity Boundary Conditions

Velocity fields use toroidal-poloidal decomposition: u = ∇×(Tr) + ∇×∇×(Pr)

Available BC Types

Physical Interpretation

No-Slip (Rigid Boundary)

All velocity components = 0 at boundary
u_r = u_θ = u_φ = 0

• Appropriate for solid boundaries (e.g., rigid inner core, solid mantle)
• Most common choice for geodynamo simulations
• Toroidal: T = 0 directly sets tangential velocity to zero
• Poloidal: ∂P/∂r = 0 ensures radial velocity gradient vanishes

Stress-Free (Free-Slip)

Zero tangential stress at boundary
σ_rθ = σ_rφ = 0, with v_r = 0

• Appropriate for fluid-fluid interfaces or idealized boundaries
• Toroidal condition derived from: σrθ = η·r·∂(vθ/r)/∂r = 0 → ∂T/∂r = T/r
• Poloidal condition: ∂²P/∂r² = 0 ensures stress-free radial component

Usage

grid = SphericalShellGrid(nr = 64, lmax = 31)
model = GeodynamoModel(
    grid;
    velocity_bcs = BoundaryConditions(inner = NoSlip(), outer = NoSlip()),
)
simulation = Simulation(model; Δt = 1e-5, stop_time = 0.02)

Magnetic Field Boundary Conditions

Magnetic fields also use toroidal-poloidal decomposition: B = ∇×(Tr) + ∇×∇×(Pr)

Available BC Types

Physical Interpretation

Insulating Boundary

Exterior region has σ = 0 (no electrical conductivity)
Current cannot flow across boundary: J_n = 0
Field matches potential field solution outside

The l-dependent conditions ensure proper matching to potential field solutions. Because Br = l(l+1)P/r², the interior potential field (Φ ∝ r^l, so Br ∝ r^{l-1}) corresponds to P ∝ r^{l+1}, and the exterior potential field (Φ ∝ r^{-(l+1)}, so B_r ∝ r^{-(l+2)}) corresponds to P ∝ r^{-l}; matching P′/P at each boundary gives:

• Inner: (∂/∂r - (l+1)/r)P = 0 → P behaves as r^{l+1} (regular interior vacuum)
• Outer: (∂/∂r + l/r)P = 0 → P decays as r^{-l} (vanishes at infinity)

This is the standard choice for:

• CMB (core-mantle boundary): Mantle is nearly insulating
• ICB with no inner core: Insulating condition at inner boundary

Conducting Inner Core

Inner core has finite electrical conductivity
Magnetic field diffuses through inner core
B and ∂B/∂r continuous at ICB

For simulations with a solid, electrically conducting inner core:

• Field and its derivative must match across the ICB interface
• Outer boundary remains insulating (CMB)
• Inner core field satisfies regularity at r = 0

Perfect Conductor

Boundary material has σ → ∞
Tangential magnetic field = 0
B_θ = B_φ = 0 at boundary

Rarely used in geodynamo context, and not yet implemented.

Usage

Insulating is the default — embedded in the implicit solver matrices and applied automatically whenever the magnetic field is enabled. A conducting inner core is opt-in via the magnetic_inner_bc parameter.

# Default: insulating magnetic BCs (applied automatically)
params = SolverParameters(; include_magnetic = true, geometry = :shell, ...)

# Opt in to a conducting inner core (shell geometry, CNAB2 timestepper):
params = SolverParameters(; include_magnetic = true, geometry = :shell,
                            magnetic_inner_bc = :conducting_inner_core, ...)

There is no enforce_magnetic_boundary_constraints! function. Perfect conductor is not yet selectable.

Temperature Boundary Conditions

Temperature is a scalar field with simpler BC structure.

Available BC Types

BC Combinations

Special Case: Neumann-Neumann

Solution: For Neumann-Neumann BCs, GeoDynamo.jl automatically uses Dirichlet for the l=0 (mean) mode at the inner boundary. This pins the average temperature while allowing flux conditions for all other modes.

l=0 mode (spherically symmetric mean): Dirichlet at inner boundary
l>0 modes (variations): Neumann at both boundaries

Usage

# Fixed temperature boundaries
temperature_bcs = BoundaryConditions(
    inner = FixedTemperature(1.0),
    outer = FixedTemperature(0.0),
)

# Heat flux boundaries
temperature_bcs = BoundaryConditions(
    inner = FixedFlux(0.1),
    outer = FixedFlux(0.0),
)

model = GeodynamoModel(grid; temperature_bcs)

Composition Boundary Conditions

Composition follows the same structure as temperature.

Available BC Types

Physical Interpretation

• Dirichlet: Fixed light element concentration (e.g., at ICB during inner core growth)
• Neumann: Prescribed mass flux (e.g., compositional release during freezing)
• Zero flux (∂C/∂r = 0): No compositional exchange across boundary

Special Case: Neumann-Neumann

Same as temperature: when both boundaries have flux conditions, the l=0 mode uses Dirichlet at the inner boundary to pin the mean composition level.

Usage

# Fixed composition at ICB, zero flux at CMB
composition_bcs = BoundaryConditions(
    inner = FixedTemperature(1.0),
    outer = FixedFlux(0.0),
)

model = GeodynamoModel(grid; include_composition = true, composition_bcs)

Summary Table

Model Reference

grid = SphericalShellGrid(nr = 64, lmax = 31)
model = GeodynamoModel(
    grid;
    velocity_bcs = BoundaryConditions(inner = NoSlip(), outer = StressFree()),
    temperature_bcs = BoundaryConditions(inner = FixedTemperature(1.0), outer = FixedFlux(0.0)),
    composition_bcs = BoundaryConditions(inner = FixedTemperature(0.0), outer = FixedFlux(0.0)),
)
simulation = Simulation(model; Δt = 1e-5, stop_time = 0.02)

Programmatic API

# Per-field keyword form
model = GeodynamoModel(grid;
    velocity_bcs = BoundaryConditions(inner = NoSlip(), outer = StressFree()),
    temperature_bcs = BoundaryConditions(inner = FixedTemperature(1.0),
                                         outer = FixedFlux(0.0)),
    composition_bcs = BoundaryConditions(inner = FixedFlux(0.01),
                                         outer = FixedTemperature(0.0)),
    magnetic_inner_bc = :insulating,
)

# Oceananigans-style NamedTuple form
model = GeodynamoModel(grid;
    boundary_conditions = (
        velocity = FieldBoundaryConditions(inner = NoSlip(), outer = NoSlip()),
        temperature = FieldBoundaryConditions(inner = ValueBoundaryCondition(1.0),
                                              outer = FluxBoundaryCondition(0.0)),
    ),
    magnetic_inner_bc = :conducting_inner_core,
    include_magnetic = true,
)

Next Steps