Time stepping

Overview

Terrarium explicitly separates process computations in compute_auxiliary! and compute_tendencies! from the choice of time stepping scheme. As a general rule, only models can be configured for timestepping. A model can be initialized for timestepping via

Terrarium.initializeMethod
initialize(
    model::Terrarium.AbstractModel{NF, Grid} where Grid<:(Terrarium.AbstractLandGrid{NF}),
    timestepper::Terrarium.AbstractTimeStepper,
    inputs::InputSource...;
    clock,
    boundary_conditions,
    initializers,
    fields
) -> ModelIntegrator{_A, Arch, _B, _C, _D, StateVars, @NamedTuple{}} where {_A, Arch<:Oceananigans.Architectures.AbstractArchitecture, _B<:Terrarium.AbstractLandGrid{_A, Arch}, _C<:Terrarium.AbstractTimeStepper{_A}, _D<:Terrarium.AbstractModel{_A, _B}, StateVars<:(StateVariables{NF, _A, _B, _C, _D, _E, _F, _G, _H, _I, _J, ClockType} where {NF, _A, _B, _C, _D, _E, _F, _G, _H, _I, _J, ClockType<:(Clock{_A, Float64, Int64, Int64} where _A)})}

Creates and initializes a ModelIntegrator for the given model and timestepper with input variables populated by the given inputs. This method allocates all necessary Fields for the state variables and subsequently calls initialize!(::ModelIntegrator).

Note that this method is not type stable and thus should not be called from Enzyme autodiff. To reinitialize the model for an existing state, use initialize!(state, model).

See the docstring for initialize(::AbstractModel) for further details.

source

This will return a ModelIntegrator:

Terrarium.ModelIntegratorType
struct ModelIntegrator{NF, Arch<:Oceananigans.Architectures.AbstractArchitecture, Grid<:Terrarium.AbstractLandGrid{NF, Arch<:Oceananigans.Architectures.AbstractArchitecture}, TimeStepper<:Terrarium.AbstractTimeStepper{NF}, Model<:Terrarium.AbstractModel{NF, Grid<:Terrarium.AbstractLandGrid{NF, Arch<:Oceananigans.Architectures.AbstractArchitecture}}, StateVars<:Terrarium.AbstractStateVariables, Inits<:NamedTuple, Inputs<:InputSources} <: Oceananigans.AbstractModel{TimeStepper<:Terrarium.AbstractTimeStepper{NF}, Arch<:Oceananigans.Architectures.AbstractArchitecture}

Represents a "integrator" for a simulation of a given model. ModelIntegrator consists of a clock, a model, and an initialized StateVariables data structure, as well as a stage for the timestepper and any relevant inputs provided by a corresponding InputProvider. The ModelIntegrator implements the Oceananigans.AbstractModel interface and can thus be treated as a "model" in Oceananigans Simulations and output reading/writing utilities.

source

A single iteration of the ModelIntegrator roughly involves four steps:

  1. Invoke update_inputs! to populate all input Fields from their respective InputSources based on the current clock time,
  2. Invoke compute_auxiliary! on all model components to derive auxiliary state variables from the current prognostic state,
  3. Invoke compute_tendencies! on all model components to calculate tendencies for all prognostic variables,
  4. Apply the tendencies computed in step (3) to update the prognostic state based on the selected timestepping scheme. Higher order explicit or implicit timesteppers may repeat steps 1-3 multiple times within a single time step.

As an example, let's consider the construction and initialization of a SoilModel:

arch = CPU()
grid = ColumnGrid(arch, Float32, ExponentialSpacing(N=10))
model = SoilModel(grid)
integrator = initialize(model, ForwardEuler(Float32))
Integrator of SoilModel{Float32, ColumnGrid{Float32, CPU, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}}, SoilEnergyWaterCarbon{Float32, HomogeneousStratigraphy{Float32, ConstantSoilPorosity{Float32}}, SoilEnergyBalance{Float32, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float32, FreeWater, InverseQuadratic}}, SoilHydrology{Float32, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float32, BrooksCorey{FreezeCurves.SoilWaterVolume{Float32, Float32, Float32}, Float32, Float32}, UnsatKLinear{Float32}}, Nothing}, ConstantSoilCarbonDensity{Float32}}, PhysicalConstants{Float32}, SoilInitializer{Float32, QuasiThermalSteadyState{Float32}, SaturationWaterTable{Float32}, DefaultInitializer{Float32}}} with ForwardEuler{Float32}
├── Current time: 0.0
├── StateVariables{Float32}(clock = Clock{Float32, Float64}(time=0 seconds, iteration=0, last_Δt=Inf days), prognostic = (:internal_energy,), auxiliary = (:temperature, :liquid_water_fraction, :ground_temperature, :saturation_water_ice, :water_table, :hydraulic_conductivity), inputs = (), namespaces = ())

Here integrator corresponds to a ModelIntegrator configured for a ForwardEuler time stepping scheme. State variable Fields can be accessed via integrator.state:

integrator.state.temperature   # current temperature Field
1×1×10 Field{Center, Center, Center} on Oceananigans.Grids.RectilinearGrid on CPU
├── grid: 1×1×10 RectilinearGrid{Float32, Periodic, Flat, Bounded} on CPU with 1×0×3 halo
├── boundary conditions: FieldBoundaryConditions
│   └── west: Periodic, east: Periodic, south: Nothing, north: Nothing, bottom: ZeroFlux, top: ZeroFlux, immersed: Nothing
└── data: 3×1×16 OffsetArray(::Array{Float32, 3}, 0:2, 1:1, -2:13) with eltype Float32 with indices 0:2×1:1×-2:13
    └── max=2.50774, min=0.0005, mean=0.439117

Time stepping schemes

Terrarium currently provides two choices of explicit time steppers: ForwardEuler and Heun.

Terrarium.HeunType
struct Heun{NF, Stage} <: Terrarium.AbstractTimeStepper{NF}

Simple forward 2nd order Heun / improved Euler time stepping scheme.

source

Both are constructed with a default timestep Δt in seconds. Δt can be overridden at run time by specifying it when calling run! or timestep!.

Transient simulations

ModelIntegrator defines dispatches for timestep! and run! that allow for transient, open-ended simulations starting from the initialized state. We can use timestep! to take a single time step,

timestep!(integrator)

and run! to advance over a fixed number steps or a specified period,

# Advance by a fixed number of steps
run!(integrator; steps = 100)

# Advance for a given wall-clock period
run!(integrator; period = Day(30), Δt = 3600.0)
Integrator of SoilModel{Float32, ColumnGrid{Float32, CPU, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}}, SoilEnergyWaterCarbon{Float32, HomogeneousStratigraphy{Float32, ConstantSoilPorosity{Float32}}, SoilEnergyBalance{Float32, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float32, FreeWater, InverseQuadratic}}, SoilHydrology{Float32, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float32, BrooksCorey{FreezeCurves.SoilWaterVolume{Float32, Float32, Float32}, Float32, Float32}, UnsatKLinear{Float32}}, Nothing}, ConstantSoilCarbonDensity{Float32}}, PhysicalConstants{Float32}, SoilInitializer{Float32, QuasiThermalSteadyState{Float32}, SaturationWaterTable{Float32}, DefaultInitializer{Float32}}} with ForwardEuler{Float32}
├── Current time: 2.6223e6
├── StateVariables{Float32}(clock = Clock{Float32, Float64}(time=30.351 days, iteration=821, last_Δt=1 hour), prognostic = (:internal_energy,), auxiliary = (:temperature, :liquid_water_fraction, :ground_temperature, :saturation_water_ice, :water_table, :hydraulic_conductivity), inputs = (), namespaces = ())

Setting up a Simulation

The timestep! and run! methods allow us to control the integrate the model forward in time, but we only have access to the transient model state at each time step. In order for the model to be useful, we also need to be able to define a finite time period over which to run a simulation and save outputs during that time period.

Fortunately, Terrarium's ModelIntegrator implements the Oceananigans model interface. As a result, we can take advantage of the built-in Oceananigans infrastructure for configuring and running Simulations:

sim = Simulation(integrator; stop_time = 24*3600.0, Δt = 300.0)
run!(sim)
[ Info: Initializing simulation...
[ Info:     ... simulation initialization complete (162.104 μs)
[ Info: Executing initial time step...
[ Info: Simulation is stopping after running for 0 seconds.
[ Info: Simulation time 30.354 days equals or exceeds stop time 1 day.
[ Info:     ... initial time step complete (182.045 ms).

Note that stop_time is by default expressed in the same units as the clock, which is here seconds. However, Oceananigans.Units provides convenient conversion factors for other time units (e.g. days, hours):

using Oceananigans.Units: days

sim = Simulation(integrator; stop_time = 1days, Δt = 300.0)
Simulation of ModelIntegrator{Float32, CPU, ColumnGrid{Float32, CPU, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}}, ForwardEuler{Float32}, SoilModel{Float32, ColumnGrid{Float32, CPU, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}}, SoilEnergyWaterCarbon{Float32, HomogeneousStratigraphy{Float32, ConstantSoilPorosity{Float32}}, SoilEnergyBalance{Float32, Terrarium.ExplicitTwoPhaseHeatConduction, SoilEnergyTemperatureClosure, SoilThermalProperties{Float32, FreeWater, InverseQuadratic}}, SoilHydrology{Float32, NoFlow, SoilSaturationPressureClosure, SoilHydraulicsSURFEX{Float32, BrooksCorey{FreezeCurves.SoilWaterVolume{Float32, Float32, Float32}, Float32, Float32}, UnsatKLinear{Float32}}, Nothing}, ConstantSoilCarbonDensity{Float32}}, PhysicalConstants{Float32}, SoilInitializer{Float32, QuasiThermalSteadyState{Float32}, SaturationWaterTable{Float32}, DefaultInitializer{Float32}}}, StateVariables{Float32, (:internal_energy,), (:temperature, :liquid_water_fraction), (:temperature, :liquid_water_fraction, :ground_temperature, :saturation_water_ice, :water_table, :hydraulic_conductivity), (), (), Tuple{Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, NoFluxBoundaryCondition, NoFluxBoundaryCondition, Nothing, @NamedTuple{bottom_and_top::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_bottom_and_top_halo!)}, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 16)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:16)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{NoFluxBoundaryCondition, NoFluxBoundaryCondition}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}}, Tuple{Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, NoFluxBoundaryCondition, NoFluxBoundaryCondition, Nothing, @NamedTuple{bottom_and_top::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_bottom_and_top_halo!)}, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 16)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:16)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{NoFluxBoundaryCondition, NoFluxBoundaryCondition}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}}, Tuple{Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, NoFluxBoundaryCondition, NoFluxBoundaryCondition, Nothing, @NamedTuple{bottom_and_top::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_bottom_and_top_halo!)}, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 16)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:16)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{NoFluxBoundaryCondition, NoFluxBoundaryCondition}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}, Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, NoFluxBoundaryCondition, NoFluxBoundaryCondition, Nothing, @NamedTuple{bottom_and_top::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_bottom_and_top_halo!)}, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 16)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:16)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{NoFluxBoundaryCondition, NoFluxBoundaryCondition}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}, Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, UnitRange{Int64}}, OffsetArrays.OffsetArray{Float32, 3, SubArray{Float32, 3, Array{Float32, 3}, Tuple{Base.Slice{Base.OneTo{Int64}}, Base.Slice{Base.OneTo{Int64}}, UnitRange{Int64}}, true}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, Nothing, Nothing, Nothing, @NamedTuple{bottom_and_top::Nothing, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{Nothing, Nothing}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}, Field{Center, Center, Center, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, NoFluxBoundaryCondition, NoFluxBoundaryCondition, Nothing, @NamedTuple{bottom_and_top::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_bottom_and_top_halo!)}, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 16)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:16)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{NoFluxBoundaryCondition, NoFluxBoundaryCondition}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}, Field{Center, Center, Nothing, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, Nothing, Nothing, Nothing, @NamedTuple{bottom_and_top::Nothing, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 1)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:1)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{Nothing, Nothing}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}, Field{Center, Center, Face, Nothing, Oceananigans.Grids.RectilinearGrid{Float32, Oceananigans.Grids.Periodic, Oceananigans.Grids.Flat, Oceananigans.Grids.Bounded, Oceananigans.Grids.StaticVerticalDiscretization{OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}, OffsetArrays.OffsetVector{Float32, Vector{Float32}}}, Float32, Float32, OffsetArrays.OffsetVector{Float32, StepRangeLen{Float32, Float64, Float64, Int64}}, Nothing, CPU}, Tuple{Colon, Colon, Colon}, OffsetArrays.OffsetArray{Float32, 3, Array{Float32, 3}}, Float32, Oceananigans.BoundaryConditions.FieldBoundaryConditions{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Nothing, Nothing, Nothing, Nothing, Nothing, @NamedTuple{bottom_and_top::Nothing, south_and_north::Nothing, west_and_east::KernelAbstractions.Kernel{KernelAbstractions.CPU, KernelAbstractions.NDIteration.StaticSize{(1, 17)}, Oceananigans.Utils.OffsetStaticSize{(1:1, 1:17)}, typeof(Oceananigans.BoundaryConditions.cpu__fill_periodic_west_and_east_halo!)}}, @NamedTuple{bottom_and_top::Tuple{Nothing, Nothing}, south_and_north::Tuple{Nothing, Nothing}, west_and_east::Tuple{Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}, Oceananigans.BoundaryConditions.BoundaryCondition{Oceananigans.BoundaryConditions.Periodic, Nothing}}}}, Nothing, Nothing}}, Tuple{}, Tuple{}, Clock{Float32, Float64, Int64, Int64}}, @NamedTuple{}, InputSources{Tuple{}}}
├── Next time step: 5 minutes
├── run_wall_time: 0 seconds
├── run_wall_time / iteration: 0 seconds
├── stop_time: 1 day
├── stop_iteration: Inf
├── wall_time_limit: Inf
├── minimum_relative_step: 0.0
├── callbacks: OrderedDict with 3 entries:
│   ├── stop_time_exceeded => Callback of stop_time_exceeded on IterationInterval(1)
│   ├── stop_iteration_exceeded => Callback of stop_iteration_exceeded on IterationInterval(1)
│   └── wall_time_limit_exceeded => Callback of wall_time_limit_exceeded on IterationInterval(1)
└── output_writers: OrderedDict with no entries

Output writers

We can make use of Oceananigans.OutputWriters to flexibly save output to disk in a variety of formats.

The simplest choice output writer is JLD2Writer which saves selected Fields to a .jld2 file at a specified schedule:

using Oceananigans: JLD2Writer, TimeInterval
using Oceananigans.Units: hours

output_file = "$(tempname()).jld2"
println("Writing output to $(output_file)")

sim.output_writers[:soil] = JLD2Writer(
    integrator,
    (temperature = integrator.state.temperature,
     saturation  = integrator.state.saturation_water_ice);
    filename = output_file,
    overwrite_existing = true,
    including = [:grid], # save the grid with the output
    schedule = TimeInterval(2hours),
)

# Re-initialize integrator
Terrarium.initialize!(integrator)

# Run the simulation again
run!(sim)
Writing output to /tmp/jl_RWnEB2N01Q.jld2
[ Info: Initializing simulation...
[ Info:     ... simulation initialization complete (3.485 seconds)
[ Info: Executing initial time step...
[ Info:     ... initial time step complete (78.436 ms).
[ Info: Simulation is stopping after running for 3.672 seconds.
[ Info: Simulation time 1 day equals or exceeds stop time 1 day.
Note

Always call Terrarium.initialize!(integrator) before calling run!(sim) when re-running a simulation; otherwise, the simulation will not run due to the stopping condition already being satisfied.

The saved file can be read back as a FieldTimeSeries for post-processing:

using Oceananigans: FieldTimeSeries

temperature_series = FieldTimeSeries(output_file, "temperature")
temperature_series[end]  # Extract temperature Field at the last saved time
1×1×10 Field{Center, Center, Center} on Oceananigans.Grids.RectilinearGrid on CPU
├── grid: 1×1×10 RectilinearGrid{Float32, Periodic, Flat, Bounded} on CPU with 1×0×3 halo
├── boundary conditions: FieldBoundaryConditions
│   └── west: Periodic, east: Periodic, south: Nothing, north: Nothing, bottom: ZeroFlux, top: ZeroFlux, immersed: Nothing
└── data: 3×1×16 OffsetArray(view(::Array{Float32, 4}, :, :, :, 13), 0:2, 1:1, -2:13) with eltype Float32 with indices 0:2×1:1×-2:13
    └── max=2.50774, min=0.00546988, mean=0.440098

Output writers can also accept subtypes AbstractSchedule via the schedule keyword. Schedules determine the frequency and aggregation of the output Fields:

Schedule typeDescription
TimeInterval(Δt)Write every Δt seconds of simulation time
IterationInterval(n)Write every n timesteps
AveragedTimeInterval(Δt)Write time-averaged outputs over windows of Δt seconds

Multiple output writers can be added to the same simulation, e.g. to save different variables at different frequencies:

sim.output_writers[:fast]   = JLD2Writer(integrator, (temperature = ...,); schedule = TimeInterval(1hours),  ...)
sim.output_writers[:daily]  = JLD2Writer(integrator, (pressure_head = ...,); schedule = AveragedTimeInterval(1days), ...)

Callbacks

Terrarium simulations inherit the full Oceananigans Callback machinery. A callback is a function that is called by Simulation at a given schedule during run!:

using Oceananigans: Callback, IterationInterval

function print_progress(sim)
    t = sim.model.clock.time
    iter = sim.model.clock.iteration
    @info "Iteration $iter, time $t s"
end

sim.callbacks[:progress] = Callback(print_progress, IterationInterval(100))
Callback of print_progress on IterationInterval(100)

The callback function receives the Simulation object, giving it access to the full integrator via sim.model, and to the current state via sim.model.state.