Library
Documenting the public user interface.
Oceananigans.jl
Oceananigans.Oceananigans
— ModuleMain module for Oceananigans.jl
– a Julia software for fast, friendly, flexible, data-driven, ocean-flavored fluid dynamics on CPUs and GPUs.
Abstract operations
Oceananigans.AbstractOperations.volume
— Constantvolume = VolumeMetric()
Instance of VolumeMetric
that generates BinaryOperation
s between AbstractField
s and their cell volumes. Summing this BinaryOperation
yields an integral of AbstractField
over the domain.
Example
julia> using Oceananigans
julia> using Oceananigans.AbstractOperations: volume
julia> c = CenterField(RectilinearGrid(size=(2, 2, 2), extent=(1, 2, 3)));
julia> c .= 1;
julia> c_dV = c * volume
BinaryOperation at (Center, Center, Center)
├── grid: 2×2×2 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
* at (Center, Center, Center)
├── 2×2×2 Field{Center, Center, Center} on RectilinearGrid on CPU
└── Vᶜᶜᶜ at (Center, Center, Center)
julia> c_dV[1, 1, 1]
0.75
julia> sum(c_dV)
6.0
Oceananigans.AbstractOperations.Δz
— ConstantΔz = ZSpacingMetric()
Instance of ZSpacingMetric
that generates BinaryOperation
s between AbstractField
s and the vertical grid spacing evaluated at the same location as the AbstractField
.
Δx
and Δy
play a similar role for horizontal grid spacings.
Example
julia> using Oceananigans
julia> using Oceananigans.AbstractOperations: Δz
julia> c = CenterField(RectilinearGrid(size=(1, 1, 1), extent=(1, 2, 3)));
julia> c_dz = c * Δz # returns BinaryOperation between Field and GridMetricOperation
BinaryOperation at (Center, Center, Center)
├── grid: 1×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
* at (Center, Center, Center)
├── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
└── Δzᶜᶜᶜ at (Center, Center, Center)
julia> c .= 1;
julia> c_dz[1, 1, 1]
3.0
Oceananigans.AbstractOperations.Average
— MethodAverage(field::AbstractField; condition = nothing, mask = 0, dims=:)
Return Reduction
representing a spatial average of field
over dims
.
Over regularly-spaced dimensions this is equivalent to a numerical mean!
.
Over dimensions of variable spacing, field
is multiplied by the appropriate grid length, area or volume, and divided by the total spatial extent of the interval.
Oceananigans.AbstractOperations.BinaryOperation
— MethodBinaryOperation{LX, LY, LZ}(op, a, b, ▶a, ▶b, grid)
Return an abstract representation of the binary operation op(▶a(a), ▶b(b))
on grid
, where ▶a
and ▶b
interpolate a
and b
to locations (LX, LY, LZ)
.
Oceananigans.AbstractOperations.ConditionalOperation
— MethodConditionalOperation(operand::AbstractField;
func = identity,
condition = nothing,
mask = 0)
Return an abstract representation of a masking procedure applied when condition
is satisfied on a field described by func(operand)
.
Positional arguments
operand
: TheAbstractField
to be masked (it must have agrid
property!)
Keyword arguments
func
: A unary transformation applied element-wise to the fieldoperand
at locations wherecondition == true
. Default isidentity
.condition
: either a function of(i, j, k, grid, operand)
returning a Boolean, or a 3-dimensional BooleanAbstractArray
. At locations wherecondition == false
, operand will be masked bymask
mask
: the scalar mask
condition_operand
is a convenience function used to construct a ConditionalOperation
condition_operand(func::Function, operand::AbstractField, condition, mask) = ConditionalOperation(operand; func, condition, mask)
Example
julia> using Oceananigans
julia> using Oceananigans.Fields: condition_operand
julia> c = CenterField(RectilinearGrid(size=(2, 1, 1), extent=(1, 1, 1)));
julia> f(i, j, k, grid, c) = i < 2; d = condition_operand(cos, c, f, 10)
ConditionalOperation at (Center, Center, Center)
├── operand: 2×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
├── grid: 2×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── func: cos (generic function with 30 methods)
├── condition: f (generic function with 1 method)
└── mask: 10
julia> d[1, 1, 1]
1.0
julia> d[2, 1, 1]
10
Oceananigans.AbstractOperations.Derivative
— MethodDerivative{LX, LY, LZ}(∂, arg, ▶, grid)
Return an abstract representation of the derivative ∂
on arg
, and subsequent interpolation by ▶
on grid
.
Oceananigans.AbstractOperations.Integral
— MethodIntegral(field; dims=:)
Return a Reduction
representing a spatial integral of field
over dims
.
Oceananigans.AbstractOperations.KernelFunctionOperation
— MethodKernelFunctionOperation{LX, LY, LZ}(kernel_function, grid, arguments...)
Construct a KernelFunctionOperation
at location (LX, LY, LZ)
on grid
with arguments
.
kernel_function
is called with
kernel_function(i, j, k, grid, arguments...)
Note that compute!(kfo::KernelFunctionOperation)
calls compute!
on all kfo.arguments
.
Examples
Construct a KernelFunctionOperation
that returns random numbers:
using Oceananigans
grid = RectilinearGrid(size=(1, 8, 8), extent=(1, 1, 1));
random_kernel_function(i, j, k, grid) = rand(); # use CUDA.rand on the GPU
kernel_op = KernelFunctionOperation{Center, Center, Center}(random_kernel_function, grid)
# output
KernelFunctionOperation at (Center, Center, Center)
├── grid: 1×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── kernel_function: random_kernel_function (generic function with 1 method)
└── arguments: ()
Construct a KernelFunctionOperation
using the vertical vorticity operator used internally to compute vertical vorticity on all grids:
using Oceananigans.Operators: ζ₃ᶠᶠᶜ # called with signature ζ₃ᶠᶠᶜ(i, j, k, grid, u, v)
model = HydrostaticFreeSurfaceModel(; grid);
u, v, w = model.velocities;
ζ_op = KernelFunctionOperation{Face, Face, Center}(ζ₃ᶠᶠᶜ, grid, u, v)
# output
KernelFunctionOperation at (Face, Face, Center)
├── grid: 1×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── kernel_function: ζ₃ᶠᶠᶜ (generic function with 1 method)
└── arguments: ("1×8×8 Field{Face, Center, Center} on RectilinearGrid on CPU", "1×8×8 Field{Center, Face, Center} on RectilinearGrid on CPU")
Oceananigans.AbstractOperations.UnaryOperation
— MethodUnaryOperation{LX, LY, LZ}(op, arg, ▶, grid)
Returns an abstract UnaryOperation
representing the action of op
on arg
, and subsequent interpolation by ▶
on grid
.
Oceananigans.AbstractOperations.∂x
— MethodReturn the $x$-derivative function acting at (X
, Y
, Any
).
Oceananigans.AbstractOperations.∂x
— Method∂x(L::Tuple, arg::AbstractField)
Return an abstract representation of an $x$-derivative acting on field arg
followed by interpolation to L
, where L
is a 3-tuple of Face
s and Center
s.
Oceananigans.AbstractOperations.∂x
— Method∂x(arg::AbstractField)
Return an abstract representation of a $x$-derivative acting on field arg
.
Oceananigans.AbstractOperations.∂y
— MethodReturn the $y$-derivative function acting at (X
, Y
, Any
).
Oceananigans.AbstractOperations.∂y
— Method∂y(L::Tuple, arg::AbstractField)
Return an abstract representation of a $y$-derivative acting on field arg
followed by interpolation to L
, where L
is a 3-tuple of Face
s and Center
s.
Oceananigans.AbstractOperations.∂y
— Method∂y(arg::AbstractField)
Return an abstract representation of a $y$-derivative acting on field arg
.
Oceananigans.AbstractOperations.∂z
— MethodReturn the $z$-derivative function acting at (Any
, Any
, Z
).
Oceananigans.AbstractOperations.∂z
— Method∂z(L::Tuple, arg::AbstractField)
Return an abstract representation of a $z$-derivative acting on field arg
followed by interpolation to L
, where L
is a 3-tuple of Face
s and Center
s.
Oceananigans.AbstractOperations.∂z
— Method∂z(arg::AbstractField)
Return an abstract representation of a $z$-derivative acting on field arg
.
Oceananigans.AbstractOperations.@at
— Macro@at location abstract_operation
Modify the abstract_operation
so that it returns values at location
, where location
is a 3-tuple of Face
s and Center
s.
Oceananigans.AbstractOperations.@binary
— Macro@binary op1 op2 op3...
Turn each binary function in the list (op1, op2, op3...)
into a binary operator on Oceananigans.Fields
for use in AbstractOperations
.
Note: a binary function is a function with two arguments: for example, +(x, y)
is a binary function.
Also note: a binary function in Base
must be imported to be extended: use import Base: op; @binary op
.
Example
julia> using Oceananigans, Oceananigans.AbstractOperations
julia> using Oceananigans.AbstractOperations: BinaryOperation, AbstractGridMetric, choose_location
julia> plus_or_times(x, y) = x < 0 ? x + y : x * y
plus_or_times (generic function with 1 method)
julia> @binary plus_or_times
Set{Any} with 6 elements:
:+
:/
:^
:-
:*
:plus_or_times
julia> c, d = (CenterField(RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1))) for i = 1:2);
julia> plus_or_times(c, d)
BinaryOperation at (Center, Center, Center)
├── grid: 1×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
plus_or_times at (Center, Center, Center)
├── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
└── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
Oceananigans.AbstractOperations.@multiary
— Macro@multiary op1 op2 op3...
Turn each multiary operator in the list (op1, op2, op3...)
into a multiary operator on Oceananigans.Fields
for use in AbstractOperations
.
Note that a multiary operator:
- is a function with two or more arguments: for example,
+(x, y, z)
is a multiary function; - must be imported to be extended if part of
Base
: useimport Base: op; @multiary op
; - can only be called on
Oceananigans.Field
s if the "location" is noted explicitly; see example.
Example
julia> using Oceananigans, Oceananigans.AbstractOperations
julia> harmonic_plus(a, b, c) = 1/3 * (1/a + 1/b + 1/c)
harmonic_plus (generic function with 1 method)
julia> c, d, e = Tuple(CenterField(RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1))) for i = 1:3);
julia> harmonic_plus(c, d, e) # before magic @multiary transformation
BinaryOperation at (Center, Center, Center)
├── grid: 1×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
* at (Center, Center, Center)
├── 0.3333333333333333
└── + at (Center, Center, Center)
├── / at (Center, Center, Center)
│ ├── 1
│ └── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
├── / at (Center, Center, Center)
│ ├── 1
│ └── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
└── / at (Center, Center, Center)
├── 1
└── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
julia> @multiary harmonic_plus
Set{Any} with 3 elements:
:+
:harmonic_plus
:*
julia> harmonic_plus(c, d, e)
MultiaryOperation at (Center, Center, Center)
├── grid: 1×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
harmonic_plus at (Center, Center, Center)
├── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
├── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
└── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
Oceananigans.AbstractOperations.@unary
— Macro@unary op1 op2 op3...
Turn each unary function in the list (op1, op2, op3...)
into a unary operator on Oceananigans.Fields
for use in AbstractOperations
.
Note: a unary function is a function with one argument: for example, sin(x)
is a unary function.
Also note: a unary function in Base
must be imported to be extended: use import Base: op; @unary op
.
Example
julia> using Oceananigans, Oceananigans.Grids, Oceananigans.AbstractOperations
julia> square_it(x) = x^2
square_it (generic function with 1 method)
julia> @unary square_it
Set{Any} with 10 elements:
:+
:sqrt
:square_it
:cos
:exp
:interpolate_identity
:-
:tanh
:sin
:abs
julia> c = CenterField(RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1)));
julia> square_it(c)
UnaryOperation at (Center, Center, Center)
├── grid: 1×1×1 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
└── tree:
square_it at (Center, Center, Center) via identity
└── 1×1×1 Field{Center, Center, Center} on RectilinearGrid on CPU
Advection
Oceananigans.Advection.Centered
— Typestruct Centered{N, FT, XT, YT, ZT, CA} <: AbstractCenteredAdvectionScheme{N, FT}
Centered reconstruction scheme.
Oceananigans.Advection.UpwindBiased
— Typestruct UpwindBiasedFifthOrder <: AbstractUpwindBiasedAdvectionScheme{3}
Upwind-biased fifth-order advection scheme.
Oceananigans.Advection.VectorInvariant
— MethodVectorInvariant(; vorticity_scheme::AbstractAdvectionScheme{N, FT} = EnstrophyConservingScheme(),
divergence_scheme = nothing,
vorticity_stencil = VelocityStencil(),
divergence_stencil = DefaultStencil(),
vertical_scheme = EnergyConservingScheme()) where {N, FT}
Construct a vector invariant momentum advection scheme of order N * 2 - 1
.
Keyword arguments
vorticity_scheme
: Scheme used forCenter
reconstruction of vorticity, options are upwind advection schemes -UpwindBiased
andWENO
- in addition to anEnergyConservingScheme
and anEnstrophyConservingScheme
(defaults toEnstrophyConservingScheme
)divergence_scheme
: Scheme used forFace
reconstruction of divergence. Options are upwind advection schemes -UpwindBiased
andWENO
- ornothing
. In casenothing
is specified, divergence flux is absorbed into the vertical advection term (defaults tonothing
). Ifvertical_scheme
isaEnergyConservingScheme
, divergence flux is absorbed in vertical advection and this keyword argument has no effectvorticity_stencil
: Stencil used for smoothness indicators in case of aWENO
upwind reconstruction. Choices are betweenVelocityStencil
which uses the horizontal velocity field to diagnose smoothness andDefaultStencil
which uses the variable being transported (defaults toVelocityStencil
)divergence_stencil
: same asvorticity_stencil
but for divergence reconstruction (defaults toDefaultStencil
)vertical_scheme
: Scheme used for vertical advection of horizontal momentum. It has to be consistent with the choice ofdivergence_stencil
. If the latter is aNothing
, onlyEnergyConservingScheme
is available (this keyword argument has no effect). In casedivergence_scheme
is anAbstractUpwindBiasedAdvectionScheme
,vertical_scheme
describes a flux form reconstruction of vertical momentum advection, and any advection scheme can be used -Centered
,UpwindBiased
andWENO
(defaults toEnergyConservingScheme
)
Examples
julia> using Oceananigans
julia> VectorInvariant()
Vector Invariant reconstruction, maximum order 1
Vorticity flux scheme:
└── EnstrophyConservingScheme{Float64}
Divergence flux scheme:
└── Nothing
Vertical advection scheme:
└── EnergyConservingScheme{Float64}
julia> using Oceananigans
julia> VectorInvariant(vorticity_scheme = WENO(), divergence_scheme = WENO(), vertical_scheme = WENO(order = 3))
Vector Invariant reconstruction, maximum order 5
Vorticity flux scheme:
└── WENO reconstruction order 5 with smoothness stencil Oceananigans.Advection.VelocityStencil()
Divergence flux scheme:
└── WENO reconstruction order 5 with smoothness stencil Oceananigans.Advection.DefaultStencil()
Vertical advection scheme:
└── WENO reconstruction order 3
Oceananigans.Advection.WENO
— TypeWENO([FT=Float64;]
order = 5,
grid = nothing,
zweno = true,
bounds = nothing)
Construct a weigthed essentially non-oscillatory advection scheme of order order
.
Keyword arguments
order
: The order of the WENO advection scheme. Default: 5grid
: (defaults tonothing
)zweno
: Whentrue
implement a Z-WENO formulation for the WENO weights calculation. (defaults totrue
)
Examples
julia> using Oceananigans
julia> WENO()
WENO reconstruction order 5
Smoothness formulation:
└── Z-weno
Boundary scheme:
└── WENO reconstruction order 3
Symmetric scheme:
└── Centered reconstruction order 4
Directions:
├── X regular
├── Y regular
└── Z regular
julia> using Oceananigans
julia> Nx, Nz = 16, 10;
julia> Lx, Lz = 1e4, 1e3;
julia> chebychev_spaced_z_faces(k) = - Lz/2 - Lz/2 * cos(π * (k - 1) / Nz);
julia> grid = RectilinearGrid(size = (Nx, Nz), halo = (4, 4), topology=(Periodic, Flat, Bounded),
x = (0, Lx), z = chebychev_spaced_z_faces);
julia> WENO(grid; order=7)
WENO reconstruction order 7
Smoothness formulation:
└── Z-weno
Boundary scheme:
└── WENO reconstruction order 5
Symmetric scheme:
└── Centered reconstruction order 6
Directions:
├── X regular
├── Y regular
└── Z stretched
Oceananigans.Advection.div_Uc
— Methoddiv_uc(i, j, k, grid, advection, U, c)
Calculate the divergence of the flux of a tracer quantity $c$ being advected by a velocity field, $𝛁⋅(𝐯 c)$,
1/V * [δxᶜᵃᵃ(Ax * u * ℑxᶠᵃᵃ(c)) + δyᵃᶜᵃ(Ay * v * ℑyᵃᶠᵃ(c)) + δzᵃᵃᶜ(Az * w * ℑzᵃᵃᶠ(c))]
which ends up at the location ccc
.
Oceananigans.Advection.div_𝐯u
— Methoddiv_𝐯u(i, j, k, grid, advection, U, u)
Calculate the advection of momentum in the $x$-direction using the conservative form, $𝛁⋅(𝐯 u)$,
1/Vᵘ * [δxᶠᵃᵃ(ℑxᶜᵃᵃ(Ax * u) * ℑxᶜᵃᵃ(u)) + δy_fca(ℑxᶠᵃᵃ(Ay * v) * ℑyᵃᶠᵃ(u)) + δz_fac(ℑxᶠᵃᵃ(Az * w) * ℑzᵃᵃᶠ(u))]
which ends up at the location fcc
.
Oceananigans.Advection.div_𝐯v
— Methoddiv_𝐯v(i, j, k, grid, advection, U, v)
Calculate the advection of momentum in the $y$-direction using the conservative form, $𝛁⋅(𝐯 v)$,
1/Vʸ * [δx_cfa(ℑyᵃᶠᵃ(Ax * u) * ℑxᶠᵃᵃ(v)) + δyᵃᶠᵃ(ℑyᵃᶜᵃ(Ay * v) * ℑyᵃᶜᵃ(v)) + δz_afc(ℑxᶠᵃᵃ(Az * w) * ℑzᵃᵃᶠ(w))]
which ends up at the location cfc
.
Oceananigans.Advection.div_𝐯w
— Methoddiv_𝐯w(i, j, k, grid, advection, U, w)
Calculate the advection of momentum in the $z$-direction using the conservative form, $𝛁⋅(𝐯 w)$,
1/Vʷ * [δx_caf(ℑzᵃᵃᶠ(Ax * u) * ℑxᶠᵃᵃ(w)) + δy_acf(ℑzᵃᵃᶠ(Ay * v) * ℑyᵃᶠᵃ(w)) + δzᵃᵃᶠ(ℑzᵃᵃᶜ(Az * w) * ℑzᵃᵃᶜ(w))]
which ends up at the location ccf
.
Architectures
Oceananigans.Architectures.AbstractArchitecture
— TypeAbstractArchitecture
Abstract supertype for architectures supported by Oceananigans.
Oceananigans.Architectures.CPU
— TypeCPU <: AbstractArchitecture
Run Oceananigans on one CPU node. Uses multiple threads if the environment variable JULIA_NUM_THREADS
is set.
Oceananigans.Architectures.GPU
— TypeGPU <: AbstractArchitecture
Run Oceananigans on a single NVIDIA CUDA GPU.
Boundary conditions
Oceananigans.BoundaryConditions.BoundaryCondition
— Typestruct BoundaryCondition{C<:AbstractBoundaryConditionClassification, T}
Container for boundary conditions.
Oceananigans.BoundaryConditions.BoundaryCondition
— MethodBoundaryCondition(Classification::DataType, condition)
Construct a boundary condition of type BC
with a number or array as a condition
.
Boundary condition types include Periodic
, Flux
, Value
, Gradient
, and Open
.
Oceananigans.BoundaryConditions.BoundaryCondition
— MethodBoundaryCondition(Classification::DataType, condition::Function;
parameters = nothing,
discrete_form = false,
field_dependencies=())
Construct a boundary condition of type Classification
with a function boundary condition
.
By default, the function boudnary condition
is assumed to have the 'continuous form' condition(ξ, η, t)
, where t
is time and ξ
and η
vary along the boundary. In particular:
- On
x
-boundaries,condition(y, z, t)
. - On
y
-boundaries,condition(x, z, t)
. - On
z
-boundaries,condition(x, y, t)
.
If parameters
is not nothing
, then function boundary conditions have the form func(ξ, η, t, parameters)
, where ξ
and η
are spatial coordinates varying along the boundary as explained above.
If discrete_form = true
, the function condition
is assumed to have the "discrete form",
condition(i, j, grid, clock, model_fields)
where i
, and j
are indices that vary along the boundary. If discrete_form = true
and parameters
is not nothing
, the function condition
is called with
condition(i, j, grid, clock, model_fields, parameters)
Oceananigans.BoundaryConditions.FieldBoundaryConditions
— TypeFieldBoundaryConditions(grid, location, indices=(:, :, :);
west = default_auxiliary_bc(topology(grid, 1)(), location[1]()),
east = default_auxiliary_bc(topology(grid, 1)(), location[1]()),
south = default_auxiliary_bc(topology(grid, 2)(), location[2]()),
north = default_auxiliary_bc(topology(grid, 2)(), location[2]()),
bottom = default_auxiliary_bc(topology(grid, 3)(), location[3]()),
top = default_auxiliary_bc(topology(grid, 3)(), location[3]()),
immersed = NoFluxBoundaryCondition())
Return boundary conditions for auxiliary fields (fields whose values are derived from a model's prognostic fields) on grid
and at location
.
Keyword arguments
Keyword arguments specify boundary conditions on the 6 possible boundaries:
west
, left end point in thex
-direction wherei = 1
east
, right end point in thex
-direction wherei = grid.Nx
south
, left end point in they
-direction wherej = 1
north
, right end point in they
-direction wherej = grid.Ny
bottom
, right end point in thez
-direction wherek = 1
top
, right end point in thez
-direction wherek = grid.Nz
immersed
: boundary between solid and fluid for immersed boundaries
If a boundary condition is unspecified, the default for auxiliary fields and the topology in the boundary-normal direction is used:
PeriodicBoundaryCondition
forPeriodic
directionsGradientBoundaryCondition(0)
forBounded
directions andCentered
-located fieldsnothing
forBounded
directions andFace
-located fieldsnothing
forFlat
directions and/orNothing
-located fields
Oceananigans.BoundaryConditions.FieldBoundaryConditions
— TypeFieldBoundaryConditions(; kwargs...)
Return a template for boundary conditions on prognostic fields.
Keyword arguments
Keyword arguments specify boundary conditions on the 7 possible boundaries:
west
: left end point in thex
-direction wherei = 1
east
: right end point in thex
-direction wherei = grid.Nx
south
: left end point in they
-direction wherej = 1
north
: right end point in they
-direction wherej = grid.Ny
bottom
: right end point in thez
-direction wherek = 1
top
: right end point in thez
-direction wherek = grid.Nz
immersed
: boundary between solid and fluid for immersed boundaries
If a boundary condition is unspecified, the default for prognostic fields and the topology in the boundary-normal direction is used:
PeriodicBoundaryCondition
forPeriodic
directionsNoFluxBoundaryCondition
forBounded
directions andCentered
-located fieldsImpenetrableBoundaryCondition
forBounded
directions andFace
-located fieldsnothing
forFlat
directions and/orNothing
-located fields
Oceananigans.BoundaryConditions.Flux
— Typestruct Flux <: AbstractBoundaryConditionClassification
A classification specifying a boundary condition on the flux of a field.
The sign convention is such that a positive flux represents the flux of a quantity in the positive direction. For example, a positive vertical flux implies a quantity is fluxed upwards, in the $+z$ direction.
Due to this convention, a positive flux applied to the top boundary specifies that a quantity is fluxed upwards across the top boundary and thus out of the domain. As a result, a positive flux applied to a top boundary leads to a reduction of that quantity in the interior of the domain; for example, a positive, upwards flux of heat at the top of the domain acts to cool the interior of the domain. Conversely, a positive flux applied to the bottom boundary leads to an increase of the quantity in the interior of the domain. The same logic holds for east, west, north, and south boundaries.
Oceananigans.BoundaryConditions.Gradient
— Typestruct Gradient <: AbstractBoundaryConditionClassification
A classification specifying a boundary condition on the derivative or gradient of a field. Also called a Neumann boundary condition.
Oceananigans.BoundaryConditions.Open
— Typestruct Open <: AbstractBoundaryConditionClassification
A classification that specifies the halo regions of a field directly.
For fields located at Faces, Open also specifies field value on the boundary.
Open boundary conditions are used to specify the component of a velocity field normal to a boundary and can also be used to describe nested or linked simulation domains.
Oceananigans.BoundaryConditions.Value
— Typestruct Value <: AbstractBoundaryConditionClassification
A classification specifying a boundary condition on the value of a field. Also called a Dirchlet boundary condition.
Oceananigans.BoundaryConditions.apply_x_bcs!
— MethodApply flux boundary conditions to a field c
by adding the associated flux divergence to the source term Gc
at the left and right.
Oceananigans.BoundaryConditions.apply_y_bcs!
— MethodApply flux boundary conditions to a field c
by adding the associated flux divergence to the source term Gc
at the left and right.
Oceananigans.BoundaryConditions.apply_z_bcs!
— MethodApply flux boundary conditions to a field c
by adding the associated flux divergence to the source term Gc
at the top and bottom.
Oceananigans.BoundaryConditions.fill_halo_regions!
— MethodFill halo regions in $x$, $y$, and $z$ for a given field's data.
Buoyancy models
Oceananigans.BuoyancyModels.Buoyancy
— MethodBuoyancy(; model, gravity_unit_vector=NegativeZDirection())
Construct a buoyancy
given a buoyancy model
. Optional keyword argument gravity_unit_vector
can be used to specify the direction of gravity (default NegativeZDirection()
). The buoyancy acceleration acts in the direction opposite to gravity.
Example
using Oceananigans
grid = RectilinearGrid(size=(1, 8, 8), extent=(1, 1, 1))
θ = 45 # degrees
g̃ = (0, sind(θ), cosd(θ));
buoyancy = Buoyancy(model=BuoyancyTracer(), gravity_unit_vector=g̃)
model = NonhydrostaticModel(grid=grid, buoyancy=buoyancy, tracers=:b)
# output
┌ Warning: The meaning of `gravity_unit_vector` changed in version 0.80.0.
│ In versions 0.79 and earlier, `gravity_unit_vector` indicated the direction _opposite_ to gravity.
│ In versions 0.80.0 and later, `gravity_unit_vector` indicates the direction of gravitational acceleration.
└ @ Oceananigans.BuoyancyModels ~/builds/tartarus-16/clima/oceananigans/src/BuoyancyModels/buoyancy.jl:48
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 1×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: QuasiAdamsBashforth2TimeStepper
├── tracers: b
├── closure: Nothing
├── buoyancy: BuoyancyTracer with ĝ = Tuple{Float64, Float64, Float64}
└── coriolis: Nothing
Oceananigans.BuoyancyModels.BuoyancyTracer
— TypeBuoyancyTracer <: AbstractBuoyancyModel{Nothing}
Type indicating that the tracer b
represents buoyancy.
Oceananigans.BuoyancyModels.LinearEquationOfState
— TypeLinearEquationOfState([FT=Float64;] thermal_expansion=1.67e-4, haline_contraction=7.80e-4)
Return LinearEquationOfState
for SeawaterBuoyancy
with thermal_expansion
coefficient and haline_contraction
coefficient. The buoyancy perturbation $b$ for LinearEquationOfState
is
\[ b = g (α T - β S),\]
where $g$ is gravitational acceleration, $α$ is thermal_expansion
, $β$ is haline_contraction
, $T$ is temperature, and $S$ is practical salinity units.
Default constants in units inverse Kelvin and practical salinity units for thermal_expansion
and haline_contraction
, respectively, are taken from Table 1.2 (page 33) of Vallis, "Atmospheric and Oceanic Fluid Dynamics: Fundamentals and Large-Scale Circulation" (2nd ed, 2017).
Oceananigans.BuoyancyModels.LinearEquationOfState
— TypeLinearEquationOfState{FT} <: AbstractEquationOfState
Linear equation of state for seawater.
Oceananigans.BuoyancyModels.SeawaterBuoyancy
— TypeSeawaterBuoyancy([FT = Float64;]
gravitational_acceleration = g_Earth,
equation_of_state = LinearEquationOfState(FT),
constant_temperature = false,
constant_salinity = false)
Returns parameters for a temperature- and salt-stratified seawater buoyancy model with a gravitational_acceleration
constant (typically called $g$), and an equation_of_state
that related temperature and salinity (or conservative temperature and absolute salinity) to density anomalies and buoyancy.
constant_temperature
indicates that buoyancy depends only on salinity. For a nonlinear equation of state, constant_temperature
is used as the temperature of the system. The same logic, with the roles of salinity and temperature reversed, holds when constant_salinity
is provided.
For a linear equation of state, the values of constant_temperature
or constant_salinity
are irrelevant; in this case, constant_temperature=true
(and similar for constant_salinity
) is valid input.
Oceananigans.BuoyancyModels.SeawaterBuoyancy
— TypeSeawaterBuoyancy{FT, EOS, T, S} <: AbstractBuoyancyModel{EOS}
BuoyancyModels model for seawater. T
and S
are either nothing
if both temperature and salinity are active, or of type FT
if temperature or salinity are constant, respectively.
Oceananigans.BuoyancyModels.∂x_b
— Method∂x_b(i, j, k, grid, b::SeawaterBuoyancy, C)
Returns the $x$-derivative of buoyancy for temperature and salt-stratified water,
\[∂_x b = g ( α ∂_x T - β ∂_x S ) ,\]
where $g$ is gravitational acceleration, $α$ is the thermal expansion coefficient, $β$ is the haline contraction coefficient, $T$ is conservative temperature, and $S$ is absolute salinity.
Note: In Oceananigans, model.tracers.T
is conservative temperature and model.tracers.S
is absolute salinity.
Note that $∂_x T$ (∂x_T
), $∂_x S$ (∂x_S
), $α$, and $β$ are all evaluated at cell interfaces in x
and cell centers in y
and z
.
Oceananigans.BuoyancyModels.∂y_b
— Method∂y_b(i, j, k, grid, b::SeawaterBuoyancy, C)
Returns the $y$-derivative of buoyancy for temperature and salt-stratified water,
\[∂_y b = g ( α ∂_y T - β ∂_y S ) ,\]
where $g$ is gravitational acceleration, $α$ is the thermal expansion coefficient, $β$ is the haline contraction coefficient, $T$ is conservative temperature, and $S$ is absolute salinity.
Note: In Oceananigans, model.tracers.T
is conservative temperature and model.tracers.S
is absolute salinity.
Note that $∂_y T$ (∂y_T
), $∂_y S$ (∂y_S
), $α$, and $β$ are all evaluated at cell interfaces in y
and cell centers in x
and z
.
Oceananigans.BuoyancyModels.∂z_b
— Method∂z_b(i, j, k, grid, b::SeawaterBuoyancy, C)
Returns the vertical derivative of buoyancy for temperature and salt-stratified water,
\[∂_z b = N^2 = g ( α ∂_z T - β ∂_z S ) ,\]
where $g$ is gravitational acceleration, $α$ is the thermal expansion coefficient, $β$ is the haline contraction coefficient, $T$ is conservative temperature, and $S$ is absolute salinity.
Note: In Oceananigans, model.tracers.T
is conservative temperature and model.tracers.S
is absolute salinity.
Note that $∂_z T$ (∂z_T
), $∂_z S$ (∂z_S
), $α$, and $β$ are all evaluated at cell interfaces in z
and cell centers in x
and y
.
Coriolis
Oceananigans.Coriolis.ActiveCellEnstrophyConservingScheme
— Typestruct ActiveCellEnstrophyConservingScheme
A parameter object for an enstrophy-conserving Coriolis scheme that excludes inactive (dry/land) edges (indices for which peripheral_node == true
) from the velocity interpolation.
Oceananigans.Coriolis.BetaPlane
— TypeBetaPlane([T=Float64;] f₀=nothing, β=nothing,
rotation_rate=Ω_Earth, latitude=nothing, radius=R_Earth)
Return a $β$-plane Coriolis parameter, $f = f₀ + β y$.
The user may specify both f₀
and β
, or the three parameters rotation_rate
, latitude
(in degrees), and radius
that specify the rotation rate and radius of a planet, and the central latitude (where $y = 0$) at which the β
-plane approximation is to be made.
If f₀
and β
are not specified, they are calculated from rotation_rate
, latitude
, and radius
according to the relations f₀ = 2 * rotation_rate * sind(latitude)
and β = 2 * rotation_rate * cosd(latitude) / radius
.
By default, the rotation_rate
and planet radius
are assumed to be Earth's.
Oceananigans.Coriolis.BetaPlane
— Typestruct BetaPlane{T} <: AbstractRotation
A parameter object for meridionally increasing Coriolis parameter (f = f₀ + β y
) that accounts for the variation of the locally vertical component of the rotation vector with latitude.
Oceananigans.Coriolis.ConstantCartesianCoriolis
— TypeConstantCartesianCoriolis([FT=Float64;] fx=nothing, fy=nothing, fz=nothing,
f=nothing, rotation_axis=ZDirection(),
rotation_rate=Ω_Earth, latitude=nothing)
Return a parameter object for a constant rotation decomposed into the x
, y
, and z
directions. In oceanography the components x
, y
, z
correspond to the directions east, north, and up. This constant rotation can be specified in three different ways:
- Specifying all components
fx
,fy
andfz
directly. - Specifying the Coriolis parameter
f
and (optionally) arotation_axis
(which defaults to thez
direction if not specified). - Specifying
latitude
(in degrees) and (optionally) arotation_rate
in radians per second (which defaults to Earth's rotation rate).
Oceananigans.Coriolis.ConstantCartesianCoriolis
— Typestruct ConstantCartesianCoriolis{FT} <: AbstractRotation
A Coriolis implementation that accounts for the locally vertical and possibly both local horizontal components of a constant rotation vector. This is a more general implementation of FPlane
, which only accounts for the locally vertical component.
Oceananigans.Coriolis.FPlane
— TypeFPlane([FT=Float64;] f=nothing, rotation_rate=Ω_Earth, latitude=nothing)
Return a parameter object for constant rotation at the angular frequency f/2
, and therefore with background vorticity f
, around a vertical axis. If f
is not specified, it is calculated from rotation_rate
and latitude
(in degrees) according to the relation f = 2 * rotation_rate * sind(latitude)
.
By default, rotation_rate
is assumed to be Earth's.
Also called FPlane
, after the "f-plane" approximation for the local effect of a planet's rotation in a planar coordinate system tangent to the planet's surface.
Oceananigans.Coriolis.FPlane
— Typestruct FPlane{FT} <: AbstractRotation
A parameter object for constant rotation around a vertical axis.
Oceananigans.Coriolis.HydrostaticSphericalCoriolis
— Typestruct HydrostaticSphericalCoriolis{S, FT} <: AbstractRotation
A parameter object for constant rotation around a vertical axis on the sphere.
Oceananigans.Coriolis.HydrostaticSphericalCoriolis
— MethodHydrostaticSphericalCoriolis([FT=Float64;]
rotation_rate = Ω_Earth,
scheme = EnergyConservingScheme())
Return a parameter object for Coriolis forces on a sphere rotating at rotation_rate
. By default, rotation_rate
is assumed to be Earth's.
Keyword arguments
scheme
: EitherEnergyConservingScheme()
(default),EnstrophyConservingScheme()
, orActiveCellEnstrophyConservingScheme()
.
Oceananigans.Coriolis.NonTraditionalBetaPlane
— TypeNonTraditionalBetaPlane(FT=Float64;
fz=nothing, fy=nothing, β=nothing, γ=nothing,
rotation_rate=Ω_Earth, latitude=nothing, radius=R_Earth)
The user may directly specify fz
, fy
, β
, γ
, and radius
or the three parameters rotation_rate
, latitude
(in degrees), and radius
that specify the rotation rate and radius of a planet, and the central latitude (where $y = 0$) at which the non-traditional β
-plane approximation is to be made.
If fz
, fy
, β
, and γ
are not specified, they are calculated from rotation_rate
, latitude
, and radius
according to the relations fz = 2 * rotation_rate * sind(latitude)
, fy = 2 * rotation_rate * cosd(latitude)
, β = 2 * rotation_rate * cosd(latitude) / radius
, and γ = - 4 * rotation_rate * sind(latitude) / radius
.
By default, the rotation_rate
and planet radius
is assumed to be Earth's.
Oceananigans.Coriolis.NonTraditionalBetaPlane
— Typestruct NonTraditionalBetaPlane{FT} <: AbstractRotation
A Coriolis implementation that accounts for the latitudinal variation of both the locally vertical and the locally horizontal components of the rotation vector. The "traditional" approximation in ocean models accounts for only the locally vertical component of the rotation vector (see BetaPlane
).
This implementation is based off of section 5 of Dellar (2011). It conserve energy, angular momentum, and potential vorticity.
References
Dellar, P. (2011). Variations on a beta-plane: Derivation of non-traditional beta-plane equations from Hamilton's principle on a sphere. Journal of Fluid Mechanics, 674, 174-195. doi:10.1017/S0022112010006464
Diagnostics
Oceananigans.Diagnostics.CFL
— Typestruct CFL{D, S}
An object for computing the Courant-Freidrichs-Lewy (CFL) number.
Oceananigans.Diagnostics.CFL
— MethodCFL(Δt [, timescale = Oceananigans.Advection.cell_advection_timescale])
Return an object for computing the Courant-Freidrichs-Lewy (CFL) number associated with time step Δt
or TimeStepWizard
and timescale
.
See also AdvectiveCFL
and DiffusiveCFL
.
Oceananigans.Diagnostics.StateChecker
— MethodStateChecker(; schedule, fields)
Returns a StateChecker
that logs field information (minimum, maximum, mean) for each field in a named tuple of fields
when schedule
actuates.
Oceananigans.Diagnostics.AdvectiveCFL
— MethodAdvectiveCFL(Δt)
Return an object for computing the Courant-Freidrichs-Lewy (CFL) number associated with time step Δt
or TimeStepWizard
and the time scale for advection across a cell. The advective CFL is, e.g., $U Δt / Δx$.
Example
julia> using Oceananigans
julia> model = NonhydrostaticModel(grid = RectilinearGrid(size=(16, 16, 16), extent=(8, 8, 8)));
julia> Δt = 1.0;
julia> cfl = AdvectiveCFL(Δt);
julia> model.velocities.u .= π;
julia> cfl(model)
6.283185307179586
Oceananigans.Diagnostics.DiffusiveCFL
— MethodDiffusiveCFL(Δt)
Returns an object for computing the diffusive Courant-Freidrichs-Lewy (CFL) number associated with time step Δt
or TimeStepWizard
and the time scale for diffusion across a cell associated with model.closure
. The diffusive CFL, e.g., for viscosity is $ν Δt / Δx²$.
The maximum diffusive CFL number among viscosity and all tracer diffusivities is returned.
Example
julia> using Oceananigans
julia> model = NonhydrostaticModel(grid = RectilinearGrid(size=(16, 16, 16), extent=(1, 1, 1)),
closure = ScalarDiffusivity(; ν = 1e-2));
julia> Δt = 0.1;
julia> dcfl = DiffusiveCFL(Δt);
julia> dcfl(model)
0.256
Distributed
Oceananigans.Distributed.DistributedArch
— TypeDistributedArch(child_architecture = CPU();
topology = (Periodic, Periodic, Periodic),
ranks,
use_buffers = false,
devices = nothing,
communicator = MPI.COMM_WORLD)
Constructor for a distributed architecture that uses MPI for communications
Positional arguments
child_architecture
: Specifies whether the computation is performed on CPUs or GPUs. Default:child_architecture = CPU()
.
Keyword arguments
topology
: the topology we want the grid to have. It is used to establish connectivity. Default:topology = (Periodic, Periodic, Periodic)
.ranks
(required): A 3-tuple(Rx, Ry, Rz)
specifying the total processors in thex
,y
andz
direction. NOTE: support for distributed z direction is limited, soRz = 1
is strongly suggested.use_buffers
: iftrue
, buffered halo communication is implemented. Iffalse
, halos will be exchanged through views. Buffered communication is not necessary in case ofCPU
execution, but it is necessary forGPU
execution without CUDA-aware MPIdevices
:GPU
device linked to local rank. The GPU will be assigned based on the local node rank as suchdevices[node_rank]
. Make sure to run--ntasks-per-node
<=--gres=gpu
. Ifnothing
, the devices will be assigned automatically based on the available resourcescommunicator
: the MPI communicator,MPI.COMM_WORLD
. This keyword argument should not be tampered with if not for testing or developing. Change at your own risk!
Oceananigans.Distributed.DistributedFFTBasedPoissonSolver
— MethodDistributedFFTBasedPoissonSolver(global_grid, local_grid)
Return a FFT-based solver for the Poisson equation,
\[∇²φ = b\]
for DistributedArch
itectures.
Supported configurations
We support two "modes":
Vertical pencil decompositions: two-dimensional decompositions in $(x, y)$ for three dimensional problems that satisfy either
Nz > Rx
orNz > Ry
.One-dimensional decompositions in either $x$ or $y$.
Above, Nz = size(global_grid, 3)
and Rx, Ry, Rz = architecture(local_grid).ranks
.
Other configurations that are decomposed in $(x, y)$ but have too few Nz
, or any configuration decomposed in $z$, are not supported.
Algorithm for two-dimensional decompositions
For two-dimensional decompositions (useful for three-dimensional problems), there are three forward transforms, three backward transforms, and four transpositions requiring MPI communication. In the schematic below, the first dimension is always the local dimension. In our implementation of the PencilFFTs algorithm, we require either Nz >= Rx
, or Nz >= Ry
, where Nz
is the number of vertical cells, Rx
is the number of ranks in $x$, and Ry
is the number of ranks in $y$. Below, we outline the algorithm for the case Nz >= Rx
. If Nz < Rx
, but Nz > Ry
, a similar algorithm applies with $x$ and $y$ swapped:
first(storage)
is initialized with layout $(z, x, y)$.- Transform along $z$.
3 Transpose + communicate to storage[2]
in layout $(x, z, y)$, which is distributed into (Rx, Ry)
processes in $(z, y)$.
- Transform along $x$.
- Transpose + communicate to
last(storage)
in layout $(y, x, z)$, which is distributed into(Rx, Ry)
processes in $(x, z)$. - Transform in $y$.
At this point the three in-place forward transforms are complete, and we solve the Poisson equation by updating last(storage)
. Then the process is reversed to obtain first(storage)
in physical space and with the layout $(z, x, y)$.
Restrictions
The algorithm for two-dimensional decompositions requires that Nz = size(global_grid, 3)
is larger than either Rx = ranks[1]
or Ry = ranks[2]
, where ranks
are configured when building DistributedArch
. If Nz
does not satisfy this condition, we can only support a one-dimensional decomposition.
Algorithm for one-dimensional decompositions
This algorithm requires a one-dimensional decomposition with either Rx = 1
or Ry = 1
, and is important to support two-dimensional transforms.
For one-dimensional decompositions, we place the decomposed direction last. If the number of ranks is Rh = max(Rx, Ry)
, this algorithm requires that both Nx > Rh
and Ny > Rh
. The resulting flow of transposes and transforms is similar to the two-dimensional case. It remains somewhat of a mystery why this succeeds (i.e., why the last transform is correctly decomposed).
Oceananigans.Distributed.reconstruct_global_grid
— Methodreconstruct_global_grid(grid::DistributedGrid)
Return the global grid on child_architecture(grid)
Fields
Oceananigans.Fields.AbstractField
— TypeAbstractField{LX, LY, LZ, G, T, N}
Abstract supertype for fields located at (LX, LY, LZ)
and defined on a grid G
with eltype T
and N
dimensions.
Note: we need the parameter T
to subtype AbstractArray.
Oceananigans.Fields.BackgroundField
— TypeBackgroundField{F, P}
Temporary container for storing information about BackgroundFields
.
Oceananigans.Fields.BackgroundField
— MethodBackgroundField(func; parameters=nothing)
Returns a BackgroundField
to be passed to NonhydrostaticModel
for use as a background velocity or tracer field.
If parameters
is not provided, func
must be callable with the signature
func(x, y, z, t)
If parameters
is provided, func
must be callable with the signature
func(x, y, z, t, parameters)
Oceananigans.Fields.Field
— MethodField{LX, LY, LZ}(grid::AbstractGrid,
T::DataType=eltype(grid); kw...) where {LX, LY, LZ}
Construct a Field
on grid
with data type T
at the location (LX, LY, LZ)
. Each of (LX, LY, LZ)
is either Center
or Face
and determines the field's location in (x, y, z)
respectively.
Keyword arguments
data :: OffsetArray
: An offset array with the fields data. If nothing is provided the field is filled with zeros.boundary_conditions
: If nothing is provided, then field is created using the default
FieldBoundaryConditions
.indices
: Used to prescribe where a reduced field lives on. For example, at whichk
index does a two-dimensional $x$-$y$ field lives on. Default:(:, :, :)
.
Example
A field at location (Face, Face, Center)
.
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(2, 3, 4), extent=(1, 1, 1));
julia> ω = Field{Face, Face, Center}(grid)
2×3×4 Field{Face, Face, Center} on RectilinearGrid on CPU
├── grid: 2×3×4 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── boundary conditions: FieldBoundaryConditions
│ └── west: Periodic, east: Periodic, south: Periodic, north: Periodic, bottom: ZeroFlux, top: ZeroFlux, immersed: ZeroFlux
└── data: 8×9×10 OffsetArray(::Array{Float64, 3}, -2:5, -2:6, -2:7) with eltype Float64 with indices -2:5×-2:6×-2:7
└── max=0.0, min=0.0, mean=0.0
Now, using indices
we can create a two dimensional $x$-$y$ field at location (Face, Face, Center)
to compute, e.g., the vertical vorticity $∂v/∂x - ∂u/∂y$ at the fluid's surface $z = 0$, which for Center
corresponds to k = Nz
.
julia> u = XFaceField(grid); v = YFaceField(grid);
julia> ωₛ = Field(∂x(v) - ∂y(u), indices=(:, :, grid.Nz))
2×3×1 Field{Face, Face, Center} on RectilinearGrid on CPU
├── grid: 2×3×4 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── boundary conditions: FieldBoundaryConditions
│ └── west: Periodic, east: Periodic, south: Periodic, north: Periodic, bottom: Nothing, top: Nothing, immersed: ZeroFlux
├── operand: BinaryOperation at (Face, Face, Center)
├── status: time=0.0
└── data: 8×9×1 OffsetArray(::Array{Float64, 3}, -2:5, -2:6, 4:4) with eltype Float64 with indices -2:5×-2:6×4:4
└── max=0.0, min=0.0, mean=0.0
julia> compute!(ωₛ)
2×3×1 Field{Face, Face, Center} on RectilinearGrid on CPU
├── grid: 2×3×4 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── boundary conditions: FieldBoundaryConditions
│ └── west: Periodic, east: Periodic, south: Periodic, north: Periodic, bottom: Nothing, top: Nothing, immersed: ZeroFlux
├── operand: BinaryOperation at (Face, Face, Center)
├── status: time=0.0
└── data: 8×9×1 OffsetArray(::Array{Float64, 3}, -2:5, -2:6, 4:4) with eltype Float64 with indices -2:5×-2:6×4:4
└── max=0.0, min=0.0, mean=0.0
Oceananigans.Fields.Reduction
— MethodReduction(reduce!, operand; dims)
Return a Reduction
of operand
with reduce!
, along dims
. Note that Reduction
expects reduce!
to operate in-place.
Example
using Oceananigans
Nx, Ny, Nz = 3, 3, 3
grid = RectilinearGrid(size=(Nx, Ny, Nz), x=(0, 1), y=(0, 1), z=(0, 1),
topology=(Periodic, Periodic, Periodic))
c = CenterField(grid)
set!(c, (x, y, z) -> x + y + z)
max_c² = Field(Reduction(maximum!, c^2, dims=3))
compute!(max_c²)
max_c²[1:Nx, 1:Ny]
# output
3×3 Matrix{Float64}:
1.36111 2.25 3.36111
2.25 3.36111 4.69444
3.36111 4.69444 6.25
Oceananigans.Fields.CenterField
— FunctionCenterField(grid; kw...)
Return a Field{Center, Center, Center}
on grid
. Additional keyword arguments are passed to the Field
constructor.
Oceananigans.Fields.PressureFields
— FunctionPressureFields(grid, bcs::NamedTuple)
Return a NamedTuple
with pressure fields pHY′
and pNHS
initialized as CenterField
s on grid
. Boundary conditions bcs
may be specified via a named tuple of FieldBoundaryCondition
s.
Oceananigans.Fields.PressureFields
— MethodPressureFields(proposed_pressures::NamedTuple{(:pHY′, :pNHS)}, grid, bcs)
Return a NamedTuple
of pressure fields with, overwriting boundary conditions in proposed_tracer_fields
with corresponding fields in the NamedTuple
bcs
.
Oceananigans.Fields.TendencyFields
— MethodTendencyFields(grid, tracer_names;
u = XFaceField(grid),
v = YFaceField(grid),
w = ZFaceField(grid),
kwargs...)
Return a NamedTuple
with tendencies for all solution fields (velocity fields and tracer fields), initialized on grid
. Optional kwargs
can be specified to assign data arrays to each tendency field.
Oceananigans.Fields.TracerFields
— MethodTracerFields(tracer_names, grid, user_bcs)
Return a NamedTuple
with tracer fields specified by tracer_names
initialized as CenterField
s on grid
. Boundary conditions user_bcs
may be specified via a named tuple of FieldBoundaryCondition
s.
Oceananigans.Fields.TracerFields
— MethodTracerFields(tracer_names, grid; kwargs...)
Return a NamedTuple
with tracer fields specified by tracer_names
initialized as CenterField
s on grid
. Fields may be passed via optional keyword arguments kwargs
for each field.
This function is used by OutputWriters.Checkpointer
and TendencyFields
. ```
Oceananigans.Fields.TracerFields
— MethodTracerFields(proposed_tracers::NamedTuple, grid, bcs)
Return a NamedTuple
of tracers, overwriting boundary conditions in proposed_tracers
with corresponding fields in the NamedTuple
bcs
.
Oceananigans.Fields.TracerFields
— MethodShortcut constructor for empty tracer fields.
Oceananigans.Fields.VelocityFields
— FunctionVelocityFields(grid, user_bcs = NamedTuple())
Return a NamedTuple
with fields u
, v
, w
initialized on grid
. Boundary conditions bcs
may be specified via a named tuple of FieldBoundaryCondition
s.
Oceananigans.Fields.VelocityFields
— MethodVelocityFields(proposed_velocities::NamedTuple{(:u, :v, :w)}, grid, bcs)
Return a NamedTuple
of velocity fields, overwriting boundary conditions in proposed_velocities
with corresponding fields in the NamedTuple
bcs
.
Oceananigans.Fields.XFaceField
— FunctionXFaceField(grid; kw...)
Return a Field{Face, Center, Center}
on grid
. Additional keyword arguments are passed to the Field
constructor.
Oceananigans.Fields.YFaceField
— FunctionYFaceField(grid; kw...)
Return a Field{Center, Face, Center}
on grid
. Additional keyword arguments are passed to the Field
constructor.
Oceananigans.Fields.ZFaceField
— FunctionZFaceField(grid; kw...)
Return a Field{Center, Center, Face}
on grid
. Additional keyword arguments are passed to the Field
constructor.
Oceananigans.Fields.compute!
— Functioncompute!(field)
Computes field.data
from field.operand
.
Oceananigans.Fields.field
— Methodfield(loc, a, grid)
Build a field from a
at loc
and on grid
.
Oceananigans.Fields.interior
— Methodinterior(f::Field)
Returns a view of f
that excludes halo points."
Oceananigans.Fields.interpolate
— Methodinterpolate(field, x, y, z)
Interpolate field
to the physical point (x, y, z)
using trilinear interpolation.
Oceananigans.Fields.regrid!
— Methodregrid!(a, b)
Regrid field b
onto the grid of field a
.
Example
Generate a tracer field on a vertically stretched grid and regrid it on a regular grid.
using Oceananigans
Nz, Lz = 2, 1.0
topology = (Flat, Flat, Bounded)
input_grid = RectilinearGrid(size=Nz, z = [0, Lz/3, Lz], topology=topology, halo=1)
input_field = CenterField(input_grid)
input_field[1, 1, 1:Nz] = [2, 3]
output_grid = RectilinearGrid(size=Nz, z=(0, Lz), topology=topology, halo=1)
output_field = CenterField(output_grid)
regrid!(output_field, input_field)
output_field[1, 1, :]
# output
4-element OffsetArray(::Vector{Float64}, 0:3) with eltype Float64 with indices 0:3:
0.0
2.333333333333333
3.0
0.0
Oceananigans.location
— MethodReturns the location (LX, LY, LZ)
of an AbstractField{LX, LY, LZ}
.
Oceananigans.Fields.@compute
— Macro@compute(exprs...)
Call compute!
on fields after defining them.
Forcings
Oceananigans.Forcings.AdvectiveForcing
— TypeAdvectiveForcing(scheme=UpwindBiasedFifthOrder(), u=ZeroField(), v=ZeroField(), w=ZeroField())
Build a forcing term representing advection by the velocity field u, v, w
with an advection scheme
.
Example
Using a tracer field to model sinking particles
using Oceananigans
# Physical parameters
gravitational_acceleration = 9.81 # m s⁻²
ocean_density = 1026 # kg m⁻³
mean_particle_density = 2000 # kg m⁻³
mean_particle_radius = 1e-3 # m
ocean_molecular_kinematic_viscosity = 1.05e-6 # m² s⁻¹
# Terminal velocity of a sphere in viscous flow
Δb = gravitational_acceleration * (mean_particle_density - ocean_density) / ocean_density
ν = ocean_molecular_kinematic_viscosity
R = mean_particle_radius
w_Stokes = - 2/9 * Δb / ν * R^2 # m s⁻¹
settling = AdvectiveForcing(UpwindBiasedFifthOrder(), w=w_Stokes)
# output
AdvectiveForcing with the UpwindBiased scheme:
├── u: ZeroField{Int64}
├── v: ZeroField{Int64}
└── w: ConstantField(-1.97096)
Oceananigans.Forcings.ContinuousForcing
— TypeContinuousForcing{LX, LY, LZ, P, F, D, I, ℑ}
A callable object that implements a "continuous form" forcing function on a field at the location LX, LY, LZ
with optional parameters.
Oceananigans.Forcings.ContinuousForcing
— MethodContinuousForcing(func; parameters=nothing, field_dependencies=())
Construct a "continuous form" forcing with optional parameters
and optional field_dependencies
on other fields in a model.
If neither parameters
nor field_dependencies
are provided, then func
must be callable with the signature
func(x, y, z, t)
where x, y, z
are the east-west, north-south, and vertical spatial coordinates, and t
is time.
If field_dependencies
are provided, the signature of func
must include them. For example, if field_dependencies=(:u, :S)
(and parameters
are not provided), then func
must be callable with the signature
func(x, y, z, t, u, S)
where u
is assumed to be the u
-velocity component, and S
is a tracer. Note that any field which does not have the name u
, v
, or w
is assumed to be a tracer and must be present in model.tracers
.
If parameters
are provided, then the last argument to func
must be parameters
. For example, if func
has no field_dependencies
but does depend on parameters
, then it must be callable with the signature
func(x, y, z, t, parameters)
With field_dependencies=(:u, :v, :w, :c)
and parameters
, then func
must be callable with the signature
func(x, y, z, t, u, v, w, c, parameters)
Oceananigans.Forcings.DiscreteForcing
— Typestruct DiscreteForcing{P, F}
Wrapper for "discrete form" forcing functions with optional parameters
.
Oceananigans.Forcings.DiscreteForcing
— MethodDiscreteForcing(func; parameters=nothing)
Construct a "discrete form" forcing function with optional parameters. The forcing function is applied at grid point i, j, k
.
When parameters
are not specified, func
must be callable with the signature
func(i, j, k, grid, clock, model_fields)
where grid
is model.grid
, clock.time
is the current simulation time and clock.iteration
is the current model iteration, and model_fields
is a NamedTuple
with u, v, w
and the fields in model.tracers
.
Note that the index end
does not access the final physical grid point of a model field in any direction. The final grid point must be explicitly specified, as in model_fields.u[i, j, grid.Nz]
.
When parameters
is specified, func
must be callable with the signature.
func(i, j, k, grid, clock, model_fields, parameters)
Above, parameters
is, in principle, arbitrary. Note, however, that GPU compilation can place constraints on typeof(parameters)
.
Oceananigans.Forcings.GaussianMask
— TypeGaussianMask{D}(center, width)
Callable object that returns a Gaussian masking function centered on center
, with width
, and varying along direction D
, i.e.,
exp(-(D - center)^2 / (2 * width^2))
Examples
- Create a Gaussian mask centered on
z=0
with width1
meter.
julia> mask = GaussianMask{:z}(center=0, width=1)
Oceananigans.Forcings.LinearTarget
— TypeLinearTarget{D}(intercept, gradient)
Callable object that returns a Linear target function with intercept
and gradient
, and varying along direction D
, i.e.,
intercept + D * gradient
Examples
Create a linear target function varying in
z
, equal to0
atz=0
and with gradient 10⁻⁶:julia> target = LinearTarget{:z}(intercept=0, gradient=1e-6)
Oceananigans.Forcings.Relaxation
— Typestruct Relaxation{R, M, T}
Callable object for restoring fields to a target
at some rate
and within a mask
ed region in x, y, z
.
Oceananigans.Forcings.Relaxation
— MethodRelaxation(; rate, mask=onefunction, target=zerofunction)
Returns a Forcing
that restores a field to target(x, y, z, t)
at the specified rate
, in the region mask(x, y, z)
.
The functions onefunction
and zerofunction
always return 1 and 0, respectively. Thus the default mask
leaves the whole domain uncovered, and the default target
is zero.
Example
- Restore a field to zero on a timescale of "3600" (equal to one hour if the time units of the simulation are seconds).
using Oceananigans
damping = Relaxation(rate = 1/3600)
# output
Relaxation{Float64, typeof(Oceananigans.Forcings.onefunction), typeof(Oceananigans.Forcings.zerofunction)}
├── rate: 0.0002777777777777778
├── mask: 1
└── target: 0
- Restore a field to a linear z-gradient within the bottom 1/4 of a domain on a timescale of "60" (equal to one minute if the time units of the simulation are seconds).
dTdz = 0.001 # ⁰C m⁻¹, temperature gradient
T₀ = 20 # ⁰C, surface temperature at z=0
Lz = 100 # m, depth of domain
bottom_sponge_layer = Relaxation(; rate = 1/60,
target = LinearTarget{:z}(intercept=T₀, gradient=dTdz),
mask = GaussianMask{:z}(center=-Lz, width=Lz/4))
# output
Relaxation{Float64, GaussianMask{:z, Float64}, LinearTarget{:z, Float64}}
├── rate: 0.016666666666666666
├── mask: exp(-(z + 100.0)^2 / (2 * 25.0^2))
└── target: 20.0 + 0.001 * z
Oceananigans.Forcings.Forcing
— MethodForcing(func; parameters=nothing, field_dependencies=(), discrete_form=false)
Returns a forcing function added to the tendency of an Oceananigans model field.
If discrete_form=false
(the default), and neither parameters
nor field_dependencies
are provided, then func
must be callable with the signature
func(x, y, z, t)
where x, y, z
are the east-west, north-south, and vertical spatial coordinates, and t
is time. Note that this form is also default in the constructor for NonhydrostaticModel
, so that Forcing
is not needed.
If discrete_form=false
(the default), and field_dependencies
are provided, the signature of func
must include them. For example, if field_dependencies=(:u, :S)
(and parameters
are not provided), then func
must be callable with the signature
func(x, y, z, t, u, S)`
where u
is assumed to be the u
-velocity component, and S
is a tracer. Note that any field which does not have the name u
, v
, or w
is assumed to be a tracer and must be present in model.tracers
.
If discrete_form=false
(the default) and parameters
are provided, then the last argument to func
must be parameters
. For example, if func
has no field_dependencies
but does depend on parameters
, then it must be callable with the signature
func(x, y, z, t, parameters)
The object parameters
is arbitrary in principle, however GPU compilation can place constraints on typeof(parameters)
.
With field_dependencies=(:u, :v, :w, :c)
and parameters
, then func
must be callable with the signature
func(x, y, z, t, u, v, w, c, parameters)
If discrete_form=true
then func
must be callable with the "discrete form"
func(i, j, k, grid, clock, model_fields)
where i, j, k
is the grid point at which the forcing is applied, grid
is model.grid
, clock.time
is the current simulation time and clock.iteration
is the current model iteration, and model_fields
is a NamedTuple
with u, v, w
, the fields in model.tracers
, and the fields in model.diffusivity_fields
, each of which is an OffsetArray
s (or NamedTuple
s of OffsetArray
s depending on the turbulence closure) of field data.
When discrete_form=true
and parameters
is specified, func
must be callable with the signature
func(i, j, k, grid, clock, model_fields, parameters)
Examples
using Oceananigans
# Parameterized forcing
parameterized_func(x, y, z, t, p) = p.μ * exp(z / p.λ) * cos(p.ω * t)
v_forcing = Forcing(parameterized_func, parameters = (μ=42, λ=0.1, ω=π))
# output
ContinuousForcing{NamedTuple{(:μ, :λ, :ω), Tuple{Int64, Float64, Irrational{:π}}}}
├── func: parameterized_func (generic function with 1 method)
├── parameters: (μ = 42, λ = 0.1, ω = π)
└── field dependencies: ()
Note that because forcing locations are regularized within the NonhydrostaticModel
constructor:
grid = RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1))
model = NonhydrostaticModel(grid=grid, forcing=(v=v_forcing,))
model.forcing.v
# output
ContinuousForcing{NamedTuple{(:μ, :λ, :ω), Tuple{Int64, Float64, Irrational{:π}}}} at (Center, Face, Center)
├── func: parameterized_func (generic function with 1 method)
├── parameters: (μ = 42, λ = 0.1, ω = π)
└── field dependencies: ()
After passing through the constructor for NonhydrostaticModel
, the v
-forcing location information is available and set to Center, Face, Center
.
# Field-dependent forcing
growth_in_sunlight(x, y, z, t, P) = exp(z) * P
plankton_forcing = Forcing(growth_in_sunlight, field_dependencies=:P)
# output
ContinuousForcing{Nothing}
├── func: growth_in_sunlight (generic function with 1 method)
├── parameters: nothing
└── field dependencies: (:P,)
# Parameterized, field-dependent forcing
tracer_relaxation(x, y, z, t, c, p) = p.μ * exp((z + p.H) / p.λ) * (p.dCdz * z - c)
c_forcing = Forcing(tracer_relaxation,
field_dependencies = :c,
parameters = (μ=1/60, λ=10, H=1000, dCdz=1))
# output
ContinuousForcing{NamedTuple{(:μ, :λ, :H, :dCdz), Tuple{Float64, Int64, Int64, Int64}}}
├── func: tracer_relaxation (generic function with 1 method)
├── parameters: (μ = 0.016666666666666666, λ = 10, H = 1000, dCdz = 1)
└── field dependencies: (:c,)
# Unparameterized discrete-form forcing function
filtered_relaxation(i, j, k, grid, clock, model_fields) =
@inbounds - (model_fields.c[i-1, j, k] + model_fields.c[i, j, k] + model_fields.c[i+1, j, k]) / 3
filtered_forcing = Forcing(filtered_relaxation, discrete_form=true)
# output
DiscreteForcing{Nothing}
├── func: filtered_relaxation (generic function with 1 method)
└── parameters: nothing
# Discrete-form forcing function with parameters
masked_damping(i, j, k, grid, clock, model_fields, parameters) =
@inbounds - parameters.μ * exp(grid.zᵃᵃᶜ[k] / parameters.λ) * model_fields.u[i, j, k]
masked_damping_forcing = Forcing(masked_damping, parameters=(μ=42, λ=π), discrete_form=true)
# output
DiscreteForcing{NamedTuple{(:μ, :λ), Tuple{Int64, Irrational{:π}}}}
├── func: masked_damping (generic function with 1 method)
└── parameters: (μ = 42, λ = π)
Grids
Oceananigans.Grids.AbstractCurvilinearGrid
— TypeAbstractCurvilinearGrid{FT, TX, TY, TZ}
Abstract supertype for curvilinear grids with elements of type FT
and topology {TX, TY, TZ}
.
Oceananigans.Grids.AbstractGrid
— TypeAbstractGrid{FT, TX, TY, TZ}
Abstract supertype for grids with elements of type FT
and topology {TX, TY, TZ}
.
Oceananigans.Grids.AbstractHorizontallyCurvilinearGrid
— TypeAbstractHorizontallyCurvilinearGrid{FT, TX, TY, TZ}
Abstract supertype for horizontally-curvilinear grids with elements of type FT
and topology {TX, TY, TZ}
.
Oceananigans.Grids.AbstractRectilinearGrid
— TypeAbstractRectilinearGrid{FT, TX, TY, TZ}
Abstract supertype for rectilinear grids with elements of type FT
and topology {TX, TY, TZ}
.
Oceananigans.Grids.AbstractTopology
— TypeAbstractTopology
Abstract supertype for grid topologies.
Oceananigans.Grids.AbstractUnderlyingGrid
— TypeAbstractUnderlyingGrid{FT, TX, TY, TZ}
Abstract supertype for "primary" grids (as opposed to grids with immersed boundaries) with elements of type FT
and topology {TX, TY, TZ}
.
Oceananigans.Grids.Bounded
— TypeBounded
Grid topology for bounded dimensions, e.g., wall-bounded dimensions.
Oceananigans.Grids.Center
— TypeCenter
A type describing the location at the center of a grid cell.
Oceananigans.Grids.Face
— TypeFace
A type describing the location at the face of a grid cell.
Oceananigans.Grids.Flat
— TypeFlat
Grid topology for flat dimensions, generally with one grid point, along which the solution is uniform and does not vary.
Oceananigans.Grids.FullyConnected
— TypeFullyConnected
Grid topology for dimensions that are connected to other models or domains.
Oceananigans.Grids.LatitudeLongitudeGrid
— TypeLatitudeLongitudeGrid([architecture = CPU(), FT = Float64];
size,
longitude,
latitude,
z = nothing,
radius = R_Earth,
topology = nothing,
precompute_metrics = true,
halo = nothing)
Creates a LatitudeLongitudeGrid
with coordinates (λ, φ, z)
denoting longitude, latitude, and vertical coordinate respectively.
Positional arguments
architecture
: Specifies whether arrays of coordinates and spacings are stored on the CPU or GPU. Default:CPU()
.FT
: Floating point data type. Default:Float64
.
Keyword arguments
size
(required): A 3-tuple prescribing the number of grid points each direction.longitude
(required),latitude
(required),z
(default:nothing
): Each is either a:- 2-tuple that specify the end points of the domain,
- one-dimensional array specifying the cell interface locations, or
- a single-argument function that takes an index and returns cell interface location.
Note: the latitude and longitude coordinates extents are expected in degrees.
radius
: The radius of the sphere the grid lives on. By default is equal to the radius of Earth.topology
: Tuple of topologies (Flat
,Bounded
,Periodic
) for each direction. The verticaltopology[3]
must beBounded
, while the latitude-longitude topologies can beBounded
,Periodic
, orFlat
. If no topology is provided then, by default, the topology is (Periodic
,Bounded
,Bounded
) if the latitudinal extent is 360 degrees or (Bounded
,Bounded
,Bounded
) otherwise.precompute_metrics
: Boolean specifying whether to precompute horizontal spacings and areas. Default:true
. Whenfalse
, horizontal spacings and areas are computed on-the-fly during a simulation.halo
: A 3-tuple of integers specifying the size of the halo region of cells surrounding the physical interior. The default is 3 halo cells in every direction.
Examples
- A default grid with
Float64
type:
julia> using Oceananigans
julia> grid = LatitudeLongitudeGrid(size=(36, 34, 25),
longitude = (-180, 180),
latitude = (-85, 85),
z = (-1000, 0))
36×34×25 LatitudeLongitudeGrid{Float64, Periodic, Bounded, Bounded} on CPU with 3×3×3 halo and with precomputed metrics
├── longitude: Periodic λ ∈ [-180.0, 180.0) regularly spaced with Δλ=10.0
├── latitude: Bounded φ ∈ [-85.0, 85.0] regularly spaced with Δφ=5.0
└── z: Bounded z ∈ [-1000.0, 0.0] regularly spaced with Δz=40.0
- A bounded spherical sector with cell interfaces stretched hyperbolically near the top:
julia> using Oceananigans
julia> σ = 1.1; # stretching factor
julia> Nz = 24; # vertical resolution
julia> Lz = 1000; # depth (m)
julia> hyperbolically_spaced_faces(k) = - Lz * (1 - tanh(σ * (k - 1) / Nz) / tanh(σ));
julia> grid = LatitudeLongitudeGrid(size=(36, 34, Nz),
longitude = (-180, 180),
latitude = (-20, 20),
z = hyperbolically_spaced_faces,
topology = (Bounded, Bounded, Bounded))
36×34×24 LatitudeLongitudeGrid{Float64, Bounded, Bounded, Bounded} on CPU with 3×3×3 halo and with precomputed metrics
├── longitude: Bounded λ ∈ [-180.0, 180.0] regularly spaced with Δλ=10.0
├── latitude: Bounded φ ∈ [-20.0, 20.0] regularly spaced with Δφ=1.17647
└── z: Bounded z ∈ [-1000.0, -0.0] variably spaced with min(Δz)=21.3342, max(Δz)=57.2159
Oceananigans.Grids.LeftConnected
— TypeLeftConnected
Grid topology for dimensions that are connected to other models or domains only on the left (the other direction is bounded)
Oceananigans.Grids.OrthogonalSphericalShellGrid
— TypeOrthogonalSphericalShellGrid(architecture::AbstractArchitecture = CPU(),
FT::DataType = Float64;
size,
z,
topology = (Bounded, Bounded, Bounded),
ξ = (-1, 1),
η = (-1, 1),
radius = R_Earth,
halo = (1, 1, 1),
rotation = nothing)
Create a OrthogonalSphericalShellGrid
that represents a section of a sphere after it has been mapped from the face of a cube. The cube's coordinates are ξ
and η
(which, by default, take values in the range $[-1, 1]$.
The mapping from the face of the cube to the sphere is done via the CubedSphere.jl package.
Positional arguments
architecture
: Specifies whether arrays of coordinates and spacings are stored on the CPU or GPU. Default:CPU()
.FT
: Floating point data type. Default:Float64
.
Keyword arguments
size
(required): A 3-tuple prescribing the number of grid points each direction.z
(required): Either a- 2-tuple that specify the end points of the $z$-domain,
- one-dimensional array specifying the cell interface locations, or
- a single-argument function that takes an index and returns cell interface location.
radius
: The radius of the sphere the grid lives on. By default is equal to the radius of Earth.halo
: A 3-tuple of integers specifying the size of the halo region of cells surrounding the physical interior. The default is 1 halo cells in every direction.rotation
: Rotation of the spherical shell grid about some axis that passes through the center of the sphere. Ifnothing
is provided (default), then the spherical shell includes the North Pole of the sphere in its center.
Examples
- A default grid with
Float64
type:
julia> using Oceananigans
julia> grid = OrthogonalSphericalShellGrid(size=(36, 34, 25), z=(-1000, 0))
36×34×25 OrthogonalSphericalShellGrid{Float64, Bounded, Bounded, Bounded} on CPU with 1×1×1 halo and with precomputed metrics
├── longitude: Bounded λ ∈ [-176.397, 180.0] variably spaced with min(Δλ)=48351.7, max(Δλ)=2.87833e5
├── latitude: Bounded φ ∈ [35.2644, 90.0] variably spaced with min(Δφ)=50632.2, max(Δφ)=3.04768e5
└── z: Bounded z ∈ [-1000.0, 0.0] regularly spaced with Δz=40.0
Oceananigans.Grids.Periodic
— TypePeriodic
Grid topology for periodic dimensions.
Oceananigans.Grids.RectilinearGrid
— TypeRectilinearGrid([architecture = CPU(), FT = Float64];
size,
x = nothing,
y = nothing,
z = nothing,
halo = nothing,
extent = nothing,
topology = (Periodic, Periodic, Bounded))
Create a RectilinearGrid
with size = (Nx, Ny, Nz)
grid points.
Positional arguments
architecture
: Specifies whether arrays of coordinates and spacings are stored on the CPU or GPU. Default:CPU()
.FT
: Floating point data type. Default:Float64
.
Keyword arguments
size
(required): A tuple prescribing the number of grid points in non-Flat
directions.size
is a 3-tuple for 3D models, a 2-tuple for 2D models, and either a scalar or 1-tuple for 1D models.topology
: A 3-tuple(TX, TY, TZ)
specifying the topology of the domain.TX
,TY
, andTZ
specify whether thex
-,y
-, andz
directions arePeriodic
,Bounded
, orFlat
. The topologyFlat
indicates that a model does not vary in those directions so that derivatives and interpolation are zero. The default istopology = (Periodic, Periodic, Bounded)
.extent
: A tuple prescribing the physical extent of the grid in non-Flat
directions, e.g.,(Lx, Ly, Lz)
. All directions are contructed with regular grid spacing and the domain (in the case that no direction isFlat
) is $0 ≤ x ≤ L_x$, $0 ≤ y ≤ L_y$, and $-L_z ≤ z ≤ 0$, which is most appropriate for oceanic applications in which $z = 0$ usually is the ocean's surface.x
,y
, andz
: Each ofx, y, z
are either (i) 2-tuples that specify the end points of the domain in their respect directions (in which case scalar values may be used inFlat
directions), or (ii) arrays or functions of the corresponding indicesi
,j
, ork
that specify the locations of cell faces in thex
-,y
-, orz
-direction, respectively. For example, to prescribe the cell faces inz
we need to provide a function that takesk
as argument and retuns the location of the faces for indicesk = 1
throughk = Nz + 1
, whereNz
is thesize
of the stretchedz
dimension.
Note: Either extent
, or all of x
, y
, and z
must be specified.
halo
: A tuple of integers that specifies the size of the halo region of cells surrounding the physical interior for each non-Flat
direction. The default is 3 halo cells in every direction.
The physical extent of the domain can be specified either via x
, y
, and z
keyword arguments indicating the left and right endpoints of each dimensions, e.g., x = (-π, π)
or via the extent
argument, e.g., extent = (Lx, Ly, Lz)
, which specifies the extent of each dimension in which case $0 ≤ x ≤ L_x$, $0 ≤ y ≤ L_y$, and $-L_z ≤ z ≤ 0$.
A grid topology may be specified via a tuple assigning one of Periodic
, Bounded
, and, Flat
to each dimension. By default, a horizontally periodic grid topology (Periodic, Periodic, Bounded)
is assumed.
Constants are stored using floating point values of type FT
. By default this is Float64
. Make sure to specify the desired FT
if not using Float64
.
Grid properties
(Nx, Ny, Nz) :: Int
: Number of physical points in the $(x, y, z)$-direction.(Hx, Hy, Hz) :: Int
: Number of halo points in the $(x, y, z)$-direction.(Lx, Ly, Lz) :: FT
: Physical extent of the grid in the $(x, y, z)$-direction.(Δxᶜᵃᵃ, Δyᵃᶜᵃ, Δzᵃᵃᶜ)
: Grid spacing in the $(x, y, z)$-direction between cell centers. Defined at cell centers in $x$, $y$, and $z$.(Δxᶠᵃᵃ, Δyᵃᶠᵃ, Δzᵃᵃᶠ)
: Grid spacing in the $(x, y, z)$-direction between cell faces. Defined at cell faces in $x$, $y$, and $z$.(xᶜᵃᵃ, yᵃᶜᵃ, zᵃᵃᶜ)
: $(x, y, z)$ coordinates of cell centers.(xᶠᵃᵃ, yᵃᶠᵃ, zᵃᵃᶠ)
: $(x, y, z)$ coordinates of cell faces.
Examples
- A grid with the default
Float64
type:
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(32, 32, 32), extent=(1, 2, 3))
32×32×32 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── Periodic x ∈ [0.0, 1.0) regularly spaced with Δx=0.03125
├── Periodic y ∈ [0.0, 2.0) regularly spaced with Δy=0.0625
└── Bounded z ∈ [-3.0, 0.0] regularly spaced with Δz=0.09375
- A grid with
Float32
type:
julia> using Oceananigans
julia> grid = RectilinearGrid(Float32; size=(32, 32, 16), x=(0, 8), y=(-10, 10), z=(-π, π))
32×32×16 RectilinearGrid{Float32, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── Periodic x ∈ [0.0, 8.0) regularly spaced with Δx=0.25
├── Periodic y ∈ [-10.0, 10.0) regularly spaced with Δy=0.625
└── Bounded z ∈ [-3.14159, 3.14159] regularly spaced with Δz=0.392699
- A two-dimenisional, horizontally-periodic grid:
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(32, 32), extent=(2π, 4π), topology=(Periodic, Periodic, Flat))
32×32×1 RectilinearGrid{Float64, Periodic, Periodic, Flat} on CPU with 3×3×0 halo
├── Periodic x ∈ [3.60072e-17, 6.28319) regularly spaced with Δx=0.19635
├── Periodic y ∈ [7.20145e-17, 12.5664) regularly spaced with Δy=0.392699
└── Flat z
- A one-dimensional "column" grid:
julia> using Oceananigans
julia> grid = RectilinearGrid(size=256, z=(-128, 0), topology=(Flat, Flat, Bounded))
1×1×256 RectilinearGrid{Float64, Flat, Flat, Bounded} on CPU with 0×0×3 halo
├── Flat x
├── Flat y
└── Bounded z ∈ [-128.0, 0.0] regularly spaced with Δz=0.5
- A horizontally-periodic regular grid with cell interfaces stretched hyperbolically near the top:
julia> using Oceananigans
julia> σ = 1.1; # stretching factor
julia> Nz = 24; # vertical resolution
julia> Lz = 32; # depth (m)
julia> hyperbolically_spaced_faces(k) = - Lz * (1 - tanh(σ * (k - 1) / Nz) / tanh(σ));
julia> grid = RectilinearGrid(size = (32, 32, Nz),
x = (0, 64), y = (0, 64),
z = hyperbolically_spaced_faces)
32×32×24 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── Periodic x ∈ [0.0, 64.0) regularly spaced with Δx=2.0
├── Periodic y ∈ [0.0, 64.0) regularly spaced with Δy=2.0
└── Bounded z ∈ [-32.0, -0.0] variably spaced with min(Δz)=0.682695, max(Δz)=1.83091
- A three-dimensional grid with regular spacing in $x$, cell interfaces at Chebyshev nodes in $y$, and cell interfaces hyperbolically stretched in $z$ near the top:
julia> using Oceananigans
julia> Nx, Ny, Nz = 32, 30, 24;
julia> Lx, Ly, Lz = 200, 100, 32; # (m)
julia> chebychev_nodes(j) = - Ly/2 * cos(π * (j - 1) / Ny);
julia> σ = 1.1; # stretching factor
julia> hyperbolically_spaced_faces(k) = - Lz * (1 - tanh(σ * (k - 1) / Nz) / tanh(σ));
julia> grid = RectilinearGrid(size = (Nx, Ny, Nz),
topology = (Periodic, Bounded, Bounded),
x = (0, Lx),
y = chebychev_nodes,
z = hyperbolically_spaced_faces)
32×30×24 RectilinearGrid{Float64, Periodic, Bounded, Bounded} on CPU with 3×3×3 halo
├── Periodic x ∈ [0.0, 200.0) regularly spaced with Δx=6.25
├── Bounded y ∈ [-50.0, 50.0] variably spaced with min(Δy)=0.273905, max(Δy)=5.22642
└── Bounded z ∈ [-32.0, -0.0] variably spaced with min(Δz)=0.682695, max(Δz)=1.83091
Oceananigans.Grids.RightConnected
— TypeRightConnected
Grid topology for dimensions that are connected to other models or domains only on the right (the other direction is bounded)
Oceananigans.Grids.halo_size
— Methodhalo_size(grid)
Return a 3-tuple with the number of halo cells on either side of the domain in (x, y, z).
Oceananigans.Grids.minimum_xspacing
— Methodminimum_xspacing(grid, ℓx, ℓy, ℓz)
minimum_xspacing(grid) = minimum_xspacing(grid, Center(), Center(), Center())
Return the minimum spacing for grid
in $x$ direction at location ℓx, ℓy, ℓz
.
Examples
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(2, 4, 8), extent=(1, 1, 1));
julia> minimum_xspacing(grid, Center(), Center(), Center())
0.5
Oceananigans.Grids.minimum_yspacing
— Methodminimum_yspacing(grid, ℓx, ℓy, ℓz)
minimum_yspacing(grid) = minimum_yspacing(grid, Center(), Center(), Center())
Return the minimum spacing for grid
in $y$ direction at location ℓx, ℓy, ℓz
.
Examples
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(2, 4, 8), extent=(1, 1, 1));
julia> minimum_yspacing(grid, Center(), Center(), Center())
0.25
Oceananigans.Grids.minimum_zspacing
— Methodminimum_zspacing(grid, ℓx, ℓy, ℓz)
minimum_zspacing(grid) = minimum_zspacing(grid, Center(), Center(), Center())
Return the minimum spacing for grid
in $z$ direction at location ℓx, ℓy, ℓz
.
Examples
julia> using Oceananigans
julia> grid = RectilinearGrid(size=(2, 4, 8), extent=(1, 1, 1));
julia> minimum_zspacing(grid, Center(), Center(), Center())
0.125
Oceananigans.Grids.new_data
— Functionnew_data(FT, arch, loc, topo, sz, halo_sz, indices)
Returns an OffsetArray
of zeros of float type FT
on arch
itecture, with indices corresponding to a field on a grid
of size(grid)
and located at loc
.
Oceananigans.Grids.nodes
— Methodnodes(grid, (ℓx, ℓy, ℓz); reshape=false, with_halos=false)
nodes(grid, ℓx, ℓy, ℓz; reshape=false, with_halos=false)
Return a 3-tuple of views over the interior nodes of the grid
's native coordinates at the locations in loc=(ℓx, ℓy, ℓz)
in x, y, z
.
If reshape=true
, the views are reshaped to 3D arrays with non-singleton dimensions 1, 2, 3 for x, y, z
, respectively. These reshaped arrays can then be used in broadcast operations with 3D fields or arrays.
For RectilinearGrid
s the native coordinates are x, y, z
; for curvilinear grids, like LatitudeLongitudeGrid
or OrthogonalSphericalShellGrid
the native coordinates are λ, φ, z
.
Oceananigans.Grids.offset_data
— Functionoffset_data(underlying_data, grid::AbstractGrid, loc)
Returns an OffsetArray
that maps to underlying_data
in memory, with offset indices appropriate for the data
of a field on a grid
of size(grid)
and located at loc
.
Oceananigans.Grids.topology
— Methodtopology(grid, dim)
Return the topology of the grid
for the dim
-th dimension.
Oceananigans.Grids.topology
— Methodtopology(grid)
Return a tuple with the topology of the grid
for each dimension.
Oceananigans.Grids.total_size
— Functiontotal_size(grid, loc)
Return the "total" size of a grid
at loc
. This is a 3-tuple of integers corresponding to the number of grid points along x, y, z
.
Oceananigans.Grids.xnodes
— Methodxnodes(grid, ℓx, ℓy, ℓz, with_halos=false)
Return the positions over the interior nodes on grid
in the $x$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
See znodes
for examples.
Oceananigans.Grids.xspacings
— Methodxspacings(grid, ℓx, ℓy, ℓz; with_halos=true)
Return the spacings over the interior nodes on grid
in the $x$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
julia> using Oceananigans
julia> grid = LatitudeLongitudeGrid(size=(8, 15, 10), longitude=(-20, 60), latitude=(-10, 50), z=(-100, 0));
julia> xspacings(grid, Center(), Face(), Center())
16-element view(OffsetArray(::Vector{Float64}, -2:18), 1:16) with eltype Float64:
1.0950562585518518e6
1.1058578920188267e6
1.1112718969963323e6
1.1112718969963323e6
1.1058578920188267e6
1.0950562585518518e6
1.0789196210678827e6
1.0575265956426917e6
1.0309814069457315e6
999413.38046802
962976.3124613502
921847.720658409
876227.979424229
826339.3435524226
772424.8654621692
714747.2110712599
Oceananigans.Grids.ynodes
— Methodynodes(grid, ℓx, ℓy, ℓz, with_halos=false)
Return the positions over the interior nodes on grid
in the $y$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
See znodes
for examples.
Oceananigans.Grids.yspacings
— Methodyspacings(grid, ℓx, ℓy, ℓz; with_halos=true)
Return the spacings over the interior nodes on grid
in the $y$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
julia> using Oceananigans
julia> grid = LatitudeLongitudeGrid(size=(20, 15, 10), longitude=(0, 20), latitude=(-15, 15), z=(-100, 0));
julia> yspacings(grid, Center(), Center(), Center())
222389.85328911748
Oceananigans.Grids.znodes
— Methodznodes(grid, ℓx, ℓy, ℓz; with_halos=false)
Return the positions over the interior nodes on grid
in the $z$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
julia> using Oceananigans
julia> horz_periodic_grid = RectilinearGrid(size=(3, 3, 3), extent=(2π, 2π, 1), halo=(1, 1, 1),
topology=(Periodic, Periodic, Bounded));
julia> zC = znodes(horz_periodic_grid, Center())
3-element view(OffsetArray(::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, 0:4), 1:3) with eltype Float64:
-0.8333333333333334
-0.5
-0.16666666666666666
julia> zC = znodes(horz_periodic_grid, Center(), Center(), Center())
3-element view(OffsetArray(::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, 0:4), 1:3) with eltype Float64:
-0.8333333333333334
-0.5
-0.16666666666666666
julia> zC = znodes(horz_periodic_grid, Center(), Center(), Center(), with_halos=true)
-1.1666666666666667:0.3333333333333333:0.16666666666666666 with indices 0:4
Oceananigans.Grids.zspacings
— Methodzspacings(grid, ℓx, ℓy, ℓz; with_halos=true)
Return the spacings over the interior nodes on grid
in the $z$-direction for the location ℓx
, ℓy
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
julia> using Oceananigans
julia> grid = LatitudeLongitudeGrid(size=(20, 15, 10), longitude=(0, 20), latitude=(-15, 15), z=(-100, 0));
julia> zspacings(grid, Center(), Center(), Center())
10.0
Oceananigans.Grids.λnodes
— Methodλnodes(grid::AbstractCurvilinearGrid, ℓx, ℓy, ℓz, with_halos=false)
Return the positions over the interior nodes on a curvilinear grid
in the $λ$-direction for the location ℓλ
, ℓφ
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
See znodes
for examples.
Oceananigans.Grids.φnodes
— Methodφnodes(grid::AbstractCurvilinearGrid, ℓx, ℓy, ℓz, with_halos=false)
Return the positions over the interior nodes on a curvilinear grid
in the $φ$-direction for the location ℓλ
, ℓφ
, ℓz
. For Bounded
directions, Face
nodes include the boundary points.
See znodes
for examples.
Immersed boundaries
Oceananigans.ImmersedBoundaries.GridFittedBottom
— TypeGridFittedBottom(bottom_height, [immersed_condition=CenterImmersedCondition()])
Return an immersed boundary with an irregular bottom fit to the underlying grid.
Oceananigans.ImmersedBoundaries.ImmersedBoundaryCondition
— MethodImmersedBoundaryCondition(; interfaces...)
Return an ImmersedBoundaryCondition with conditions on individual cell interfaces ∈ (west, east, south, north, bottom, top)
between the fluid and immersed boundary.
Oceananigans.ImmersedBoundaries.ImmersedBoundaryGrid
— MethodImmersedBoundaryGrid(grid, ib::GridFittedBottom)
Return a grid with GridFittedBottom
immersed boundary.
Computes ib.bottom_height and wraps in an array.
Lagrangian particle tracking
Oceananigans.LagrangianParticleTracking.LagrangianParticles
— MethodLagrangianParticles(particles::StructArray; restitution=1.0, tracked_fields::NamedTuple=NamedTuple(), dynamics=no_dynamics)
Construct some LagrangianParticles
that can be passed to a model. The particles
should be a StructArray
and can contain custom fields. The coefficient of restitution for particle-wall collisions is specified by restitution
.
A number of tracked_fields
may be passed in as a NamedTuple
of fields. Each particle will track the value of each field. Each tracked field must have a corresponding particle property. So if T
is a tracked field, then T
must also be a custom particle property.
dynamics
is a function of (lagrangian_particles, model, Δt)
that is called prior to advecting particles. parameters
can be accessed inside the dynamics
function.
Oceananigans.LagrangianParticleTracking.LagrangianParticles
— MethodLagrangianParticles(; x, y, z, restitution=1.0, dynamics=no_dynamics, parameters=nothing)
Construct some LagrangianParticles
that can be passed to a model. The particles will have initial locations x
, y
, and z
. The coefficient of restitution for particle-wall collisions is specified by restitution
.
dynamics
is a function of (lagrangian_particles, model, Δt)
that is called prior to advecting particles. parameters
can be accessed inside the dynamics
function.
Logger
Oceananigans.Logger.OceananigansLogger
— TypeOceananigansLogger(stream::IO=stdout, level=Logging.Info; show_info_source=false)
Based on Logging.SimpleLogger, it tries to log all messages in the following format:
[yyyy/mm/dd HH:MM:SS.sss] log_level message [-@-> source_file:line_number]
where the source of the message between the square brackets is included only if show_info_source=true
or if the message is not an info level message.
Models
Non-hydrostatic models
Oceananigans.Models.NonhydrostaticModels.NonhydrostaticModel
— MethodNonhydrostaticModel(; grid,
clock = Clock{eltype(grid)}(0, 0, 1),
advection = CenteredSecondOrder(),
buoyancy = nothing,
coriolis = nothing,
stokes_drift = nothing,
forcing::NamedTuple = NamedTuple(),
closure = nothing,
boundary_conditions::NamedTuple = NamedTuple(),
tracers = (),
timestepper = :QuasiAdamsBashforth2,
background_fields::NamedTuple = NamedTuple(),
particles::ParticlesOrNothing = nothing,
biogeochemistry::AbstractBGCOrNothing = nothing, velocities = nothing, pressures = nothing, diffusivityfields = nothing, pressuresolver = nothing, immersedboundary = nothing, auxiliaryfields = NamedTuple())
Construct a model for a non-hydrostatic, incompressible fluid on grid
, using the Boussinesq approximation when buoyancy != nothing
. By default, all Bounded directions are rigid and impenetrable.
Keyword arguments
grid
: (required) The resolution and discrete geometry on which themodel
is solved. The architecture (CPU/GPU) that the model is solved on is inferred from the architecture of thegrid
. Note that the grid needs to be regularly spaced in the horizontal dimensions, $x$ and $y$.advection
: The scheme that advects velocities and tracers. SeeOceananigans.Advection
.buoyancy
: The buoyancy model. SeeOceananigans.BuoyancyModels
.coriolis
: Parameters for the background rotation rate of the model.stokes_drift
: Parameters for Stokes drift fields associated with surface waves. Default:nothing
.forcing
:NamedTuple
of user-defined forcing functions that contribute to solution tendencies.closure
: The turbulence closure formodel
. SeeOceananigans.TurbulenceClosures
.boundary_conditions
:NamedTuple
containing field boundary conditions.tracers
: A tuple of symbols defining the names of the modeled tracers, or aNamedTuple
of preallocatedCenterField
s.timestepper
: A symbol that specifies the time-stepping method. Either:QuasiAdamsBashforth2
or:RungeKutta3
.background_fields
:NamedTuple
with background fields (e.g., background flow). Default:nothing
.particles
: Lagrangian particles to be advected with the flow. Default:nothing
.biogeochemistry
: Biogeochemical model fortracers
.velocities
: The model velocities. Default:nothing
.pressures
: Hydrostatic and non-hydrostatic pressure fields. Default:nothing
.diffusivity_fields
: Diffusivity fields. Default:nothing
.pressure_solver
: Pressure solver to be used in the model. Ifnothing
(default), the model constructor chooses the default based on thegrid
provide.immersed_boundary
: The immersed boundary. Default:nothing
.auxiliary_fields
:NamedTuple
of auxiliary fields. Default:nothing
Hydrostatic free-surface models
Oceananigans.Models.HydrostaticFreeSurfaceModels.ExplicitFreeSurface
— Typestruct ExplicitFreeSurface{E, T}
The explicit free surface solver.
η::Any
: free surface elevationgravitational_acceleration::Any
: gravitational accelerations
Oceananigans.Models.HydrostaticFreeSurfaceModels.HydrostaticFreeSurfaceModel
— MethodHydrostaticFreeSurfaceModel(; grid,
clock = Clock{eltype(grid)}(0, 0, 1),
momentum_advection = CenteredSecondOrder(),
tracer_advection = CenteredSecondOrder(),
buoyancy = SeawaterBuoyancy(eltype(grid)),
coriolis = nothing,
free_surface = ImplicitFreeSurface(gravitational_acceleration=g_Earth),
forcing::NamedTuple = NamedTuple(),
closure = nothing,
boundary_conditions::NamedTuple = NamedTuple(),
tracers = (:T, :S),
particles::ParticlesOrNothing = nothing,
biogeochemistry::AbstractBGCOrNothing = nothing,
velocities = nothing,
pressure = nothing,
diffusivity_fields = nothing,
auxiliary_fields = NamedTuple(),
calculate_only_active_cells_tendencies = false
)
Construct a hydrostatic model with a free surface on grid
.
Keyword arguments
grid
: (required) The resolution and discrete geometry on whichmodel
is solved. The architecture (CPU/GPU) that the model is solve is inferred from the architecture of the grid.momentum_advection
: The scheme that advects velocities. SeeOceananigans.Advection
.tracer_advection
: The scheme that advects tracers. SeeOceananigans.Advection
.buoyancy
: The buoyancy model. SeeOceananigans.BuoyancyModels
.coriolis
: Parameters for the background rotation rate of the model.forcing
:NamedTuple
of user-defined forcing functions that contribute to solution tendencies.free_surface
: The free surface model.closure
: The turbulence closure formodel
. SeeOceananigans.TurbulenceClosures
.boundary_conditions
:NamedTuple
containing field boundary conditions.tracers
: A tuple of symbols defining the names of the modeled tracers, or aNamedTuple
of preallocatedCenterField
s.particles
: Lagrangian particles to be advected with the flow. Default:nothing
.biogeochemistry
: Biogeochemical model fortracers
.velocities
: The model velocities. Default:nothing
.pressure
: Hydrostatic pressure field. Default:nothing
.diffusivity_fields
: Diffusivity fields. Default:nothing
.auxiliary_fields
:NamedTuple
of auxiliary fields. Default:nothing
Oceananigans.Models.HydrostaticFreeSurfaceModels.ImplicitFreeSurface
— MethodImplicitFreeSurface(; solver_method=:Default, gravitational_acceleration=g_Earth, solver_settings...)
Return an implicit free-surface solver. The implicit free-surface equation is
\[\left [ 𝛁_h ⋅ (H 𝛁_h) - \frac{1}{g Δt^2} \right ] η^{n+1} = \frac{𝛁_h ⋅ 𝐐_⋆}{g Δt} - \frac{η^{n}}{g Δt^2} ,\]
where $η^n$ is the free-surface elevation at the $n$-th time step, $H$ is depth, $g$ is the gravitational acceleration, $Δt$ is the time step, $𝛁_h$ is the horizontal gradient operator, and $𝐐_⋆$ is the barotropic volume flux associated with the predictor velocity field $𝐮_⋆$, i.e.,
\[𝐐_⋆ = \int_{-H}^0 𝐮_⋆ \, 𝖽 z ,\]
where
\[𝐮_⋆ = 𝐮^n + \int_{t_n}^{t_{n+1}} 𝐆ᵤ \, 𝖽t .\]
This equation can be solved, in general, using the PreconditionedConjugateGradientSolver
but other solvers can be invoked in special cases.
If $H$ is constant, we divide through out to obtain
\[\left ( ∇^2_h - \frac{1}{g H Δt^2} \right ) η^{n+1} = \frac{1}{g H Δt} \left ( 𝛁_h ⋅ 𝐐_⋆ - \frac{η^{n}}{Δt} \right ) .\]
Thus, for constant $H$ and on grids with regular spacing in $x$ and $y$ directions, the free surface can be obtained using the FFTBasedPoissonSolver
.
solver_method
can be either of:
:FastFourierTransform
forFFTBasedPoissonSolver
:HeptadiagonalIterativeSolver
forHeptadiagonalIterativeSolver
:PreconditionedConjugateGradient
forPreconditionedConjugateGradientSolver
By default, if the grid has regular spacing in the horizontal directions then the :FastFourierTransform
is chosen, otherwise the :HeptadiagonalIterativeSolver
.
Oceananigans.Models.HydrostaticFreeSurfaceModels.PrescribedVelocityFields
— MethodPrescribedVelocityFields(; u=zerofunc, v=zerofunc, w=zerofunc, parameters=nothing)
Builds PrescribedVelocityFields
with prescribed functions u
, v
, and w
.
If isnothing(parameters)
, then u, v, w
are called with the signature
u(x, y, z, t) = # something interesting
If !isnothing(parameters)
, then u, v, w
are called with the signature
u(x, y, z, t, parameters) = # something parameterized and interesting
In the constructor for HydrostaticFreeSurfaceModel
, the functions u, v, w
are wrapped in FunctionField
and associated with the model's grid
and clock
.
Oceananigans.Models.HydrostaticFreeSurfaceModels.SplitExplicitFreeSurface
— Typestruct SplitExplicitFreeSurface
The split-explicit free surface solver.
η
: The instantaneous free surface (ReducedField
)state
: The entire state for the split-explicit solver (SplitExplicitState
)auxiliary
: Parameters for timestepping split-explicit solver (NamedTuple
)gravitational_acceleration
: Gravitational accelerationsettings
: Settings for the split-explicit scheme (NamedTuple
)
Oceananigans.Models.HydrostaticFreeSurfaceModels.SplitExplicitFreeSurface
— MethodSplitExplicitFreeSurface(; gravitational_acceleration = g_Earth, kwargs...)
Return a SplitExplicitFreeSurface
representing an explicit time discretization of oceanic free surface dynamics with gravitational_acceleration
.
Keyword Arguments
substeps
: The number of substeps that divide the range(t, t + 2Δt)
. Note that some averaging functions do not require substepping until2Δt
. The number of substeps is reduced automatically to the last index ofaveraging_weights
for whichaveraging_weights > 0
.barotropic_averaging_kernel
: function ofτ
used to average the barotropic transportU
and free surfaceη
within the barotropic advancement.τ
is the fractional substep going from 0 to 2 with the baroclinic time stept + Δt
located atτ = 1
. This function should be centered atτ = 1
, that is, $∑ (aₘ m /M) = 1$.timestepper
: Time stepping scheme used, eitherForwardBackwardScheme()
orAdamsBashforth3Scheme()
.
Shallow-water models
Oceananigans.Models.ShallowWaterModels.ShallowWaterModel
— MethodShallowWaterModel(; grid,
gravitational_acceleration,
clock = Clock{eltype(grid)}(0, 0, 1),
momentum_advection = UpwindBiasedFifthOrder(),
tracer_advection = WENO(),
mass_advection = WENO(),
coriolis = nothing,
forcing::NamedTuple = NamedTuple(),
closure = nothing,
bathymetry = nothing,
tracers = (),
diffusivity_fields = nothing,
boundary_conditions::NamedTuple = NamedTuple(),
timestepper::Symbol = :RungeKutta3,
formulation = ConservativeFormulation())
Construct a shallow water model on grid
with gravitational_acceleration
constant.
Keyword arguments
grid
: (required) The resolution and discrete geometry on whichmodel
is solved. The architecture (CPU/GPU) that the model is solve is inferred from the architecture of the grid.gravitational_acceleration
: (required) The gravitational acceleration constant.clock
: Theclock
for the model.momentum_advection
: The scheme that advects velocities. SeeOceananigans.Advection
. Default:UpwindBiasedFifthOrder()
.tracer_advection
: The scheme that advects tracers. SeeOceananigans.Advection
. Default:WENO()
.mass_advection
: The scheme that advects the mass equation. SeeOceananigans.Advection
. Default:WENO()
.coriolis
: Parameters for the background rotation rate of the model.forcing
:NamedTuple
of user-defined forcing functions that contribute to solution tendencies.closure
: The turbulence closure formodel
. SeeOceananigans.TurbulenceClosures
.bathymetry
: The bottom bathymetry.tracers
: A tuple of symbols defining the names of the modeled tracers, or aNamedTuple
of preallocatedCenterField
s.diffusivity_fields
: Stores diffusivity fields when the closures require a diffusivity to be calculated at each timestep.boundary_conditions
:NamedTuple
containing field boundary conditions.timestepper
: A symbol that specifies the time-stepping method. Either:QuasiAdamsBashforth2
or:RungeKutta3
(default).formulation
: Whether the dynamics are expressed in conservative form (ConservativeFormulation()
; default) or in non-conservative form with a vector-invariant formulation for the non-linear terms (VectorInvariantFormulation()
).
The ConservativeFormulation()
requires RectilinearGrid
. Use VectorInvariantFormulation()
with LatitudeLongitudeGrid
.
Oceananigans.Models.ShallowWaterModels.ShallowWaterScalarDiffusivity
— TypeShallowWaterScalarDiffusivity([FT::DataType=Float64;]
ν=0, ξ=0, discrete_form=false)
Return a scalar diffusivity for the shallow water model.
The diffusivity for the shallow water model is calculated as h * ν
so that we get a viscous term in the form $h^{-1} 𝛁 ⋅ (h ν t)$, where $t$ is the 2D stress tensor plus a trace, i.e., $t = 𝛁𝐮 + (𝛁𝐮)^T - ξ I ⋅ (𝛁 ⋅ 𝐮)$.
With the VectorInvariantFormulation()
(that evolves $u$ and $v$) we compute $h^{-1} 𝛁(ν h 𝛁 t)$, while with the ConservativeFormulation()
(that evolves $u h$ and $v h$) we compute $𝛁 (ν h 𝛁 t)$.
MultiRegion
Oceananigans.MultiRegion.MultiRegionGrid
— MethodMultiRegionGrid(global_grid; partition = XPartition(2), devices = nothing)
Split a global_grid
into different regions handled by devices
.
Positional Arguments
global_grid
: the grid to be divided into regions
Keyword Arguments
partition
: the partitioning required. The implemented partitioning areXPartition
(division along the x direction) andYPartition
(division along the y direction)devices
: the devices to allocate memory on.nothing
will allocate memory on theCPU
. ForGPU
computation it is possible to specify the total number ofGPU
s or the specificGPU
s to allocate memory on. The number of devices does not have to match the number of regions
Operators
Oceananigans.Operators.div_xyᶜᶜᶜ
— Methoddiv_xyᶜᶜᵃ(i, j, k, grid, u, v)
Return the discrete div_xy = ∂x u + ∂y v
of velocity field u, v
defined as
1 / Azᶜᶜᵃ * [δxᶜᵃᵃ(Δyᵃᶜᵃ * u) + δyᵃᶜᵃ(Δxᶜᵃᵃ * v)]
at i, j, k
, where Azᶜᶜᵃ
is the area of the cell centered on (Center, Center, Any) –- a tracer cell, Δy
is the length of the cell centered on (Face, Center, Any) in y
(a u
cell), and Δx
is the length of the cell centered on (Center, Face, Any) in x
(a v
cell). div_xyᶜᶜᵃ
ends up at the location cca
.
Oceananigans.Operators.divᶜᶜᶜ
— Methoddivᶜᶜᶜ(i, j, k, grid, u, v, w)
Calculate the divergence $𝛁·𝐕$ of a vector field $𝐕 = (u, v, w)$,
1/V * [δxᶜᵃᵃ(Ax * u) + δxᵃᶜᵃ(Ay * v) + δzᵃᵃᶜ(Az * w)]
which ends up at the cell centers ccc
.
Oceananigans.Operators.ζ₃ᶠᶠᶜ
— Methodζ₃ᶠᶠᶜ(i, j, k, grid, u, v)
The vertical vorticity associated with horizontal velocities $u$ and $v$.
Oceananigans.Operators.∇²ᶜᶜᶜ
— Method∇²ᶜᶜᶜ(i, j, k, grid, c)
Calculate the Laplacian of $c$ via
1/V * [δxᶜᵃᵃ(Ax * ∂xᶠᵃᵃ(c)) + δyᵃᶜᵃ(Ay * ∂yᵃᶠᵃ(c)) + δzᵃᵃᶜ(Az * ∂zᵃᵃᶠ(c))]
which ends up at the location ccc
.
Output readers
Oceananigans.OutputReaders.FieldDataset
— MethodFieldDataset(filepath;
architecture=CPU(), grid=nothing, backend=InMemory(), metadata_paths=["metadata"])
Returns a Dict
containing a FieldTimeSeries
for each field in the JLD2 file located at filepath
. Note that model output must have been saved with halos.
Keyword arguments
backend
: EitherInMemory()
(default) orOnDisk()
. TheInMemory
backend will
load the data fully in memory as a 4D multi-dimensional array while the OnDisk()
backend will lazily load field time snapshots when the FieldTimeSeries
is indexed linearly.
metadata_paths
: A list of JLD2 paths to look for metadata. By default it looks infile["metadata"]
.grid
: May be specified to override the grid used in the JLD2 file.
Oceananigans.OutputReaders.FieldTimeSeries
— MethodFieldTimeSeries(path, name;
backend = InMemory(),
grid = nothing,
iterations = nothing,
times = nothing)
Return a FieldTimeSeries
containing a time-series of the field name
load from JLD2 output located at path
.
Keyword arguments
backend
:InMemory()
to load data into a 4D array orOnDisk()
to lazily load data from disk when indexing intoFieldTimeSeries
.grid
: A grid to associated with data, in the case that the native grid was not serialized properly.iterations
: Iterations to load. Defaults to all iterations found in the file.times
: Save times to load, as determined through an approximate floating point comparison to recorded save times. Defaults to times associated withiterations
. Takes precedence overiterations
iftimes
is specified.
Oceananigans.OutputReaders.FieldTimeSeries
— MethodFieldTimeSeries{LX, LY, LZ}(grid, times, [FT=eltype(grid);]
indices = (:, :, :),
boundary_conditions = nothing)
Return a FieldTimeSeries
at location (LX, LY, LZ)
, on grid
, at times
.
Output writers
Oceananigans.OutputWriters.AveragedTimeInterval
— Typemutable struct AveragedTimeInterval <: AbstractSchedule
Container for parameters that configure and handle time-averaged output.
Oceananigans.OutputWriters.AveragedTimeInterval
— MethodAveragedTimeInterval(interval; window=interval, stride=1)
Returns a schedule
that specifies periodic time-averaging of output. The time window
specifies the extent of the time-average, which reoccurs every interval
.
output
is computed and accumulated into the average every stride
iterations during the averaging window. For example, stride=1
computs output every iteration, whereas stride=2
computes output every other iteration. Time-averages with longer stride
s are faster to compute, but less accurate.
The time-average of $a$ is a left Riemann sum corresponding to
\[⟨a⟩ = T⁻¹ \int_{tᵢ-T}^{tᵢ} a \mathrm{d} t \, ,\]
where $⟨a⟩$ is the time-average of $a$, $T$ is the time-window for averaging, and the $tᵢ$ are discrete times separated by the time interval
. The $tᵢ$ specify both the end of the averaging window and the time at which output is written.
Example
using Oceananigans.OutputWriters: AveragedTimeInterval
using Oceananigans.Utils: days
schedule = AveragedTimeInterval(4days, window=2days)
# output
AveragedTimeInterval(window=2 days, stride=1, interval=4 days)
An AveragedTimeInterval
schedule directs an output writer to time-average its outputs before writing them to disk:
using Oceananigans
using Oceananigans.OutputWriters: JLD2OutputWriter
using Oceananigans.Utils: minutes
model = NonhydrostaticModel(grid=RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1)))
simulation = Simulation(model, Δt=10minutes, stop_time=30days)
simulation.output_writers[:velocities] = JLD2OutputWriter(model, model.velocities,
filename= "averaged_velocity_data.jld2",
schedule = AveragedTimeInterval(4days, window=2days, stride=2))
# output
JLD2OutputWriter scheduled on TimeInterval(4 days):
├── filepath: ./averaged_velocity_data.jld2
├── 3 outputs: (u, v, w) averaged on AveragedTimeInterval(window=2 days, stride=2, interval=4 days)
├── array type: Array{Float64}
├── including: [:grid, :coriolis, :buoyancy, :closure]
└── max filesize: Inf YiB
Oceananigans.OutputWriters.Checkpointer
— MethodCheckpointer(model; schedule,
dir = ".",
prefix = "checkpoint",
overwrite_existing = false,
cleanup = false,
additional_kwargs...)
Construct a Checkpointer
that checkpoints the model to a JLD2 file on schedule.
The model.clock.iteration
is included in the filename to distinguish between multiple checkpoint files.
To restart or "pickup" a model from a checkpoint, specify pickup = true
when calling run!
, ensuring that the checkpoint file is in directory dir
. See run!
for more details.
Note that extra model properties
can be safely specified, but removing crucial properties such as :velocities
will make restoring from the checkpoint impossible.
The checkpointer attempts to serialize as much of the model to disk as possible, but functions or objects containing functions cannot be serialized at this time.
Keyword arguments
schedule
(required): Schedule that determines when to checkpoint.dir
: Directory to save output to. Default:"."
(current working directory).prefix
: Descriptive filename prefixed to all output files. Default: "checkpoint".overwrite_existing
: Remove existing files if their filenames conflict. Default:false
.verbose
: Log what the output writer is doing with statistics on compute/write times and file sizes. Default:false
.cleanup
: Previous checkpoint files will be deleted once a new checkpoint file is written. Default:false
.properties
: List of model properties to checkpoint. This list must contain[:grid, :architecture, :timestepper, :particles]
. Default: [:architecture, :grid, :clock, :coriolis, :buoyancy, :closure, :velocities, :tracers, :timestepper, :particles]
Oceananigans.OutputWriters.JLD2OutputWriter
— MethodJLD2OutputWriter(model, outputs; filename, schedule,
dir = ".",
indices = (:, :, :),
with_halos = false,
array_type = Array{Float64},
max_filesize = Inf,
overwrite_existing = false,
init = noinit,
including = [:grid, :coriolis, :buoyancy, :closure],
verbose = false,
part = 1,
jld2_kw = Dict{Symbol, Any}())
Construct a JLD2OutputWriter
for an Oceananigans model
that writes label, output
pairs in outputs
to a JLD2 file.
The argument outputs
may be a Dict
or NamedTuple
. The keys of outputs
are symbols or strings that "name" output data. The values of outputs
are either AbstractField
s, objects that are called with the signature output(model)
, or WindowedTimeAverage
s of AbstractFields
s, functions, or callable objects.
Keyword arguments
Filenaming
filename
(required): Descriptive filename.".jld2"
is appended tofilename
in the file path iffilename
does not end in".jld2"
.dir
: Directory to save output to. Default:"."
(current working directory).
Output frequency and time-averaging
schedule
(required):AbstractSchedule
that determines when output is saved.
Slicing and type conversion prior to output
indices
: Specifies the indices to write to disk with aTuple
ofColon
,UnitRange
, orInt
elements. Indices must beColon
,Int
, or contiguousUnitRange
. Defaults to(:, :, :)
or "all indices". If!with_halos
, halo regions are removed fromindices
. For example,indices = (:, :, 1)
will save xy-slices of the bottom-most index.with_halos
(Bool): Whether or not to slice halo regions from fields before writing output. Note, that to postprocess saved output (e.g., compute derivatives, etc) information about the boundary conditions is often crucial. In that case you might need to setwith_halos = true
.array_type
: The array type to which output arrays are converted to prior to saving. Default:Array{Float64}
.
File management
max_filesize
: The writer will stop writing to the output file once the file size exceedsmax_filesize
, and write to a new one with a consistent naming scheme ending inpart1
,part2
, etc. Defaults toInf
.overwrite_existing
: Remove existing files if their filenames conflict. Default:false
.
Output file metadata management
init
: A function of the forminit(file, model)
that runs when a JLD2 output file is initialized. Default:noinit(args...) = nothing
.including
: List of model properties to save with every file. Default:[:grid, :coriolis, :buoyancy, :closure]
Miscellaneous keywords
verbose
: Log what the output writer is doing with statistics on compute/write times and file sizes. Default:false
.part
: The starting part number used ifmax_filesize
is finite. Default: 1.jld2_kw
: Dict of kwargs to be passed tojldopen
when data is written.
Example
Write out 3D fields for $u$, $v$, $w$, and a tracer $c$, along with a horizontal average:
using Oceananigans
using Oceananigans.Utils: hour, minute
model = NonhydrostaticModel(grid=RectilinearGrid(size=(1, 1, 1), extent=(1, 1, 1)), tracers=(:c,))
simulation = Simulation(model, Δt=12, stop_time=1hour)
function init_save_some_metadata!(file, model)
file["author"] = "Chim Riggles"
file["parameters/coriolis_parameter"] = 1e-4
file["parameters/density"] = 1027
return nothing
end
c_avg = Field(Average(model.tracers.c, dims=(1, 2)))
# Note that model.velocities is NamedTuple
simulation.output_writers[:velocities] = JLD2OutputWriter(model, model.velocities,
filename = "some_data.jld2",
schedule = TimeInterval(20minute),
init = init_save_some_metadata!)
# output
JLD2OutputWriter scheduled on TimeInterval(20 minutes):
├── filepath: ./some_data.jld2
├── 3 outputs: (u, v, w)
├── array type: Array{Float64}
├── including: [:grid, :coriolis, :buoyancy, :closure]
└── max filesize: Inf YiB
and a time- and horizontal-average of tracer $c$ every 20 minutes of simulation time to a file called some_averaged_data.jld2
simulation.output_writers[:avg_c] = JLD2OutputWriter(model, (; c=c_avg),
filename = "some_averaged_data.jld2",
schedule = AveragedTimeInterval(20minute, window=5minute))
# output
JLD2OutputWriter scheduled on TimeInterval(20 minutes):
├── filepath: ./some_averaged_data.jld2
├── 1 outputs: c averaged on AveragedTimeInterval(window=5 minutes, stride=1, interval=20 minutes)
├── array type: Array{Float64}
├── including: [:grid, :coriolis, :buoyancy, :closure]
└── max filesize: Inf YiB
Oceananigans.OutputWriters.NetCDFOutputWriter
— MethodNetCDFOutputWriter(model, outputs; filename, schedule
dir = ".",
array_type = Array{Float64},
indices = nothing,
with_halos = false,
global_attributes = Dict(),
output_attributes = Dict(),
dimensions = Dict(),
overwrite_existing = false,
compression = 0,
verbose = false)
Construct a NetCDFOutputWriter
that writes (label, output)
pairs in outputs
(which should be a Dict
) to a NetCDF file, where label
is a string that labels the output and output
is either a Field
(e.g. model.velocities.u
) or a function f(model)
that returns something to be written to disk. Custom output requires the spatial dimensions
(a Dict
) to be manually specified (see examples).
Keyword arguments
Filenaming
filename
(required): Descriptive filename.".nc"
is appended tofilename
iffilename
does not end in".nc"
.dir
: Directory to save output to.
Output frequency and time-averaging
schedule
(required):AbstractSchedule
that determines when output is saved.
Slicing and type conversion prior to output
indices
: Tuple of indices of the output variables to include. Default is(:, :, :)
, which includes the full fields.with_halos
: Boolean defining whether or not to include halos in the outputs. Default:false
. Note, that to postprocess saved output (e.g., compute derivatives, etc) information about the boundary conditions is often crucial. In that case you might need to setwith_halos = true
.array_type
: The array type to which output arrays are converted to prior to saving. Default:Array{Float64}
.dimensions
: ADict
of dimension tuples to apply to outputs (required for function outputs).
File management
overwrite_existing
: Iffalse
,NetCDFOutputWriter
will be set to append tofilepath
. Iftrue
,NetCDFOutputWriter
will overwritefilepath
if it exists or create it if not. Default:false
. See NCDatasets.jl documentation for more information about itsmode
option.compression
: Determines the compression level of data (0-9; default: 0).
Miscellaneous keywords
global_attributes
: Dict of model properties to save with every file. Default:Dict()
.output_attributes
: Dict of attributes to be saved with each field variable (reasonable defaults are provided for velocities, buoyancy, temperature, and salinity; otherwiseoutput_attributes
must be user-provided).
Examples
Saving the $u$ velocity field and temperature fields, the full 3D fields and surface 2D slices to separate NetCDF files:
using Oceananigans
grid = RectilinearGrid(size=(16, 16, 16), extent=(1, 1, 1))
model = NonhydrostaticModel(grid=grid, tracers=:c)
simulation = Simulation(model, Δt=12, stop_time=3600)
fields = Dict("u" => model.velocities.u, "c" => model.tracers.c)
simulation.output_writers[:field_writer] =
NetCDFOutputWriter(model, fields, filename="fields.nc", schedule=TimeInterval(60))
# output
NetCDFOutputWriter scheduled on TimeInterval(1 minute):
├── filepath: ./fields.nc
├── dimensions: zC(16), zF(17), xC(16), yF(16), xF(16), yC(16), time(0)
├── 2 outputs: (c, u)
└── array type: Array{Float64}
simulation.output_writers[:surface_slice_writer] =
NetCDFOutputWriter(model, fields, filename="surface_xy_slice.nc",
schedule=TimeInterval(60), indices=(:, :, grid.Nz))
# output
NetCDFOutputWriter scheduled on TimeInterval(1 minute):
├── filepath: ./surface_xy_slice.nc
├── dimensions: zC(1), zF(1), xC(16), yF(16), xF(16), yC(16), time(0)
├── 2 outputs: (c, u)
└── array type: Array{Float64}
simulation.output_writers[:averaged_profile_writer] =
NetCDFOutputWriter(model, fields,
filename = "averaged_z_profile.nc",
schedule = AveragedTimeInterval(60, window=20),
indices = (1, 1, :))
# output
NetCDFOutputWriter scheduled on TimeInterval(1 minute):
├── filepath: ./averaged_z_profile.nc
├── dimensions: zC(16), zF(17), xC(1), yF(1), xF(1), yC(1), time(0)
├── 2 outputs: (c, u) averaged on AveragedTimeInterval(window=20 seconds, stride=1, interval=1 minute)
└── array type: Array{Float64}
NetCDFOutputWriter
also accepts output functions that write scalars and arrays to disk, provided that their dimensions
are provided:
using Oceananigans
Nx, Ny, Nz = 16, 16, 16
grid = RectilinearGrid(size=(Nx, Ny, Nz), extent=(1, 2, 3))
model = NonhydrostaticModel(grid=grid)
simulation = Simulation(model, Δt=1.25, stop_iteration=3)
f(model) = model.clock.time^2; # scalar output
g(model) = model.clock.time .* exp.(znodes(Center, grid)) # vector/profile output
xC, yF = xnodes(grid, Center()), ynodes(grid, Face())
XC = [xC[i] for i in 1:Nx, j in 1:Ny]
YF = [yF[j] for i in 1:Nx, j in 1:Ny]
h(model) = @. model.clock.time * sin(XC) * cos(YF) # xy slice output
outputs = Dict("scalar" => f, "profile" => g, "slice" => h)
dims = Dict("scalar" => (), "profile" => ("zC",), "slice" => ("xC", "yC"))
output_attributes = Dict(
"scalar" => Dict("longname" => "Some scalar", "units" => "bananas"),
"profile" => Dict("longname" => "Some vertical profile", "units" => "watermelons"),
"slice" => Dict("longname" => "Some slice", "units" => "mushrooms")
)
global_attributes = Dict("location" => "Bay of Fundy", "onions" => 7)
simulation.output_writers[:things] =
NetCDFOutputWriter(model, outputs,
schedule=IterationInterval(1), filename="things.nc", dimensions=dims, verbose=true,
global_attributes=global_attributes, output_attributes=output_attributes)
# output
NetCDFOutputWriter scheduled on IterationInterval(1):
├── filepath: ./things.nc
├── dimensions: zC(16), zF(17), xC(16), yF(16), xF(16), yC(16), time(0)
├── 3 outputs: (profile, slice, scalar)
└── array type: Array{Float64}
Oceananigans.OutputWriters.WindowedTimeAverage
— TypeWindowedTimeAverage(operand, model=nothing; schedule)
Returns an object for computing running averages of operand
over schedule.window
and recurring on schedule.interval
, where schedule
is an AveragedTimeInterval
. During the collection period, averages are computed every schedule.stride
iteration.
operand
may be a Oceananigans.Field
or a function that returns an array or scalar.
Calling wta(model)
for wta::WindowedTimeAverage
object returns wta.result
.
Simulations
Oceananigans.Simulations.Callback
— TypeCallback(func, schedule=IterationInterval(1); parameters=nothing)
Return Callback
that executes func
on schedule
with optional parameters
. schedule = IterationInterval(1)
by default.
If isnothing(parameters)
, func(sim::Simulation)
is called. Otherwise, func
is called via func(sim::Simulation, parameters)
.
Oceananigans.Simulations.Simulation
— MethodSimulation(model; Δt,
verbose = true,
stop_iteration = Inf,
stop_time = Inf,
wall_time_limit = Inf)
Construct a Simulation
for a model
with time step Δt
.
Keyword arguments
Δt
: Required keyword argument specifying the simulation time step. Can be aNumber
for constant time steps or aTimeStepWizard
for adaptive time-stepping.stop_iteration
: Stop the simulation after this many iterations.stop_time
: Stop the simulation once this much model clock time has passed.wall_time_limit
: Stop the simulation if it's been running for longer than this many seconds of wall clock time.
Oceananigans.Simulations.TimeStepWizard
— TypeTimeStepWizard([FT=Float64;]
cfl = 0.2,
diffusive_cfl = Inf,
max_change = 1.1,
min_change = 0.5,
max_Δt = Inf,
min_Δt = 0.0,
cell_advection_timescale = cell_advection_timescale,
cell_diffusion_timescale = infinite_diffusion_timescale)
Callback for adapting simulation to maintain the advective Courant-Freidrichs-Lewy (CFL) number to cfl
, the diffusive_cfl
, while also maintaining max_Δt
, min_Δt
, and satisfying max_change
and min_change
criteria so that the simulation's timestep simulation.Δt
is not adapted "too quickly".
For more information on the CFL number, see its wikipedia entry.
Example
To use TimeStepWizard
, adapt in a Callback
and add it to a Simulation
:
julia> simulation = Simulation(model, Δt=0.9, stop_iteration=100)
julia> wizard = TimeStepWizard(cfl=0.2)
julia> simulation.callbacks[:wizard] = Callback(wizard, IterationInterval(4))
Then when run!(simulation)
is invoked, the time-step simulation.Δt
will be updated every 4 iterations.
(Note that the name :wizard
is unimportant.)
Oceananigans.Simulations.erroring_NaNChecker!
— Methoderroring_NaNChecker!(simulation)
Toggles simulation
's NaNChecker
to throw an error when a NaN
is detected.
Oceananigans.Simulations.iteration
— Methoditeration(sim::Simulation)
Return the current simulation iteration.
Oceananigans.Simulations.run!
— Methodrun!(simulation; pickup=false)
Run a simulation
until one of simulation.stop_criteria
evaluates true
. The simulation will then stop.
Picking simulations up from a checkpoint
Simulations are "picked up" from a checkpoint if pickup
is either true
, a String
, or an Integer
greater than 0.
Picking up a simulation sets field and tendency data to the specified checkpoint, leaving all other model properties unchanged.
Possible values for pickup
are:
pickup=true
picks a simulation up from the latest checkpoint associated with theCheckpointer
insimulation.output_writers
.pickup=iteration::Int
picks a simulation up from the checkpointed file associated withiteration
and theCheckpointer
insimulation.output_writers
.pickup=filepath::String
picks a simulation up from checkpointer data infilepath
.
Note that pickup=true
and pickup=iteration
fails if simulation.output_writers
contains more than one checkpointer.
Solvers
Oceananigans.Solvers.BatchedTridiagonalSolver
— TypeBatchedTridiagonalSolver
A batched solver for large numbers of triadiagonal systems.
Oceananigans.Solvers.BatchedTridiagonalSolver
— MethodBatchedTridiagonalSolver(grid; lower_diagonal, diagonal, upper_diagonal, parameters=nothing, tridiagonal_direction=ZDirection())
Construct a solver for batched tridiagonal systems on grid
of the form
bⁱʲ¹ ϕⁱʲ¹ + cⁱʲ¹ ϕⁱʲ² = fⁱʲ¹, k = 1
aⁱʲᵏ⁻¹ ϕⁱʲᵏ⁻¹ + bⁱʲᵏ ϕⁱʲᵏ + cⁱʲᵏ ϕⁱʲᵏ⁺¹ = fⁱʲᵏ, k = 2, ..., N-1
aⁱʲᴺ⁻¹ ϕⁱʲᴺ⁻¹ + bⁱʲᴺ ϕⁱʲᴺ = fⁱʲᴺ, k = N
where a
is the lower_diagonal
, b
is the diagonal
, and c
is the upper_diagonal
. ϕ
is the solution and f
is the right hand side source term passed to solve!(ϕ, tridiagonal_solver, f)
a
, b
, c
, and f
can be specified in three ways:
A 1D array means that
aⁱʲᵏ = a[k]
.A 3D array means that
aⁱʲᵏ = a[i, j, k]
.Otherwise,
a
is assumed to be callable:- If
isnothing(parameters)
thenaⁱʲᵏ = a(i, j, k, grid, args...)
. - If
!isnothing(parameters)
thenaⁱʲᵏ = a(i, j, k, grid, parameters, args...)
.
where
args...
areVarargs
passed tosolve_batched_tridiagonal_system!(ϕ, solver, args...)
.- If
Oceananigans.Solvers.FFTBasedPoissonSolver
— TypeFFTBasedPoissonSolver(grid, planner_flag=FFTW.PATIENT)
Return an FFTBasedPoissonSolver
that solves the "generalized" Poisson equation,
\[(∇² + m) ϕ = b,\]
where $m$ is a number, using a eigenfunction expansion of the discrete Poisson operator on a staggered grid and for periodic or Neumann boundary conditions.
In-place transforms are applied to $b$, which means $b$ must have complex-valued elements (typically the same type as solver.storage
).
See solve!
for more information about the FFT-based Poisson solver algorithm.
Oceananigans.Solvers.HeptadiagonalIterativeSolver
— MethodHeptadiagonalIterativeSolver(coeffs;
grid,
iterative_solver = cg!,
maximum_iterations = prod(size(grid)),
tolerance = 1e-13,
reduced_dim = (false, false, false),
placeholder_timestep = -1.0,
preconditioner_method = :Default,
preconditioner_settings = nothing,
template = arch_array(architecture(grid), zeros(prod(size(grid)))),
verbose = false)
Return a HeptadiagonalIterativeSolver
to solve the problem A * x = b
, provided that A
is a symmetric matrix.
The solver relies on a sparse version of the matrix A
that is stored in matrix_constructors
.
In particular, given coefficients Ax
, Ay
, Az
, C
, D
, the solved problem is
Axᵢ₊₁ ηᵢ₊₁ + Axᵢ ηᵢ₋₁ + Ayⱼ₊₁ ηⱼ₊₁ + Ayⱼ ηⱼ₋₁ + Azₖ₊₁ ηₖ₊₁ + Azₖ ηₖ₋₁
- 2 ( Axᵢ₊₁ + Axᵢ + Ayⱼ₊₁ + Ayⱼ + Azₖ₊₁ + Azₖ ) ηᵢⱼₖ
+ ( Cᵢⱼₖ + Dᵢⱼₖ/Δt^2 ) ηᵢⱼₖ = b
To have the equation solved at location {Center, Center, Center}
, the coefficients must be specified at:
Ax
->{Face, Center, Center}
Ay
->{Center, Face, Center}
Az
->{Center, Center, Face}
C
->{Center, Center, Center}
D
->{Center, Center, Center}
solver.matrix
is precomputed with a placeholder timestep value of placeholder_timestep = -1.0
.
The sparse matrix A
can be constructed with:
SparseMatrixCSC(constructors...)
for CPUCuSparseMatrixCSC(constructors...)
for GPU
The matrix constructors are calculated based on the pentadiagonal coeffients passed as an input to matrix_from_coefficients
function.
To allow for variable time step, the diagonal term - Az / (g * Δt²)
is only added later on and it is updated only when the previous time step changes (previous_Δt != Δt
).
Preconditioning is done through the various methods implemented in Solvers/sparse_preconditioners.jl
.
The iterative_solver
used can is to be chosen from the IterativeSolvers.jl package. The default solver is a Conjugate Gradient (cg
):
solver = HeptadiagonalIterativeSolver((Ax, Ay, Az, C, D); grid)
Oceananigans.Solvers.PreconditionedConjugateGradientSolver
— MethodPreconditionedConjugateGradientSolver(linear_operation;
template_field,
maxiter = size(template_field.grid),
reltol = sqrt(eps(eltype(template_field.grid))),
abstol = 0,
preconditioner = nothing)
Returns a PreconditionedConjugateGradientSolver
that solves the linear equation $A x = b$ using a iterative conjugate gradient method with optional preconditioning.
The solver is used by calling
solve!(x, solver::PreconditionedConjugateGradientOperator, b, args...)
for solver
, right-hand side b
, solution x
, and optional arguments args...
.
Arguments
linear_operation
: Function with signaturelinear_operation!(p, y, args...)
that calculatesA * y
and stores the result inp
for a "candidate solution"y
.args...
are optional positional arguments passed fromsolve!(x, solver, b, args...)
.template_field
: Dummy field that is the same type and size asx
andb
, which is used to infer thearchitecture
,grid
, and to create work arrays that are used internally by the solver.maxiter
: Maximum number of iterations the solver may perform before exiting.reltol, abstol
: Relative and absolute tolerance for convergence of the algorithm. The iteration stops whennorm(A * x - b) < tolerance
.preconditioner
: Object for whichprecondition!(z, preconditioner, r, args...)
computesz = P * r
, wherer
is the residual. TypicallyP
is approximatelyA⁻¹
.
See solve!
for more information about the preconditioned conjugate-gradient algorithm.
Oceananigans.Solvers.solve!
— Functionsolve!(ϕ, solver::FFTBasedPoissonSolver, b, m=0)
Solves the "generalized" Poisson equation,
\[(∇² + m) ϕ = b,\]
where $m$ is a number, using a eigenfunction expansion of the discrete Poisson operator on a staggered grid and for periodic or Neumann boundary conditions.
In-place transforms are applied to $b$, which means $b$ must have complex-valued elements (typically the same type as solver.storage
).
Equation $(∇² + m) ϕ = b$ is sometimes referred to as the "screened Poisson" equation when $m < 0$, or the Helmholtz equation when $m > 0$.
Oceananigans.Solvers.solve!
— Methodsolve!(ϕ, solver::BatchedTridiagonalSolver, rhs, args...)
Solve the batched tridiagonal system of linear equations with right hand side rhs
and lower diagonal, diagonal, and upper diagonal coefficients described by the BatchedTridiagonalSolver
solver
. BatchedTridiagonalSolver
uses a modified TriDiagonal Matrix Algorithm (TDMA).
The result is stored in ϕ
which must have size (grid.Nx, grid.Ny, grid.Nz)
.
Reference implementation per Numerical Recipes, Press et. al 1992 (§ 2.4).
Oceananigans.Solvers.solve!
— Methodsolve!(x, solver::PreconditionedConjugateGradientSolver, b, args...)
Solve A * x = b
using an iterative conjugate-gradient method, where A * x
is determined by solver.linear_operation
See figure 2.5 in
The Preconditioned Conjugate Gradient Method in "Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods" Barrett et. al, 2nd Edition.
Given:
- Linear Preconditioner operator
M!(solution, x, other_args...)
that computesM * x = solution
- A matrix operator
A
as a functionA()
; - A dot product function
norm()
; - A right-hand side
b
; - An initial guess
x
; and - Local vectors:
z
,r
,p
,q
This function executes the psuedocode algorithm
β = 0
r = b - A(x)
iteration = 0
Loop:
if iteration > maxiter
break
end
ρ = r ⋅ z
z = M(r)
β = ρⁱ⁻¹ / ρ
p = z + β * p
q = A(p)
α = ρ / (p ⋅ q)
x = x + α * p
r = r - α * q
if |r| < tolerance
break
end
iteration += 1
ρⁱ⁻¹ = ρ
Stokes drift
Oceananigans.StokesDrift.UniformStokesDrift
— TypeUniformStokesDrift{P, UZ, VZ, UT, VT} <: AbstractStokesDrift
Parameter struct for Stokes drift fields associated with surface waves.
Oceananigans.StokesDrift.UniformStokesDrift
— MethodUniformStokesDrift(; ∂z_uˢ=addzero, ∂z_vˢ=addzero, ∂t_uˢ=addzero, ∂t_vˢ=addzero, parameters=nothing)
Construct a set of functions that describes the Stokes drift field beneath a horizontally-uniform surface gravity wave field.
Time steppers
Oceananigans.TimeSteppers.Clock
— Typemutable struct Clock{T<:Number}
Keeps track of the current time
, iteration
number, and time-stepping stage
. The stage
is updated only for multi-stage time-stepping methods. The time::T
is either a number or a DateTime
object.
Oceananigans.TimeSteppers.Clock
— MethodClock(; time, iteration=0, stage=1)
Returns a Clock
object. By default, Clock
is initialized to the zeroth iteration
and first time step stage
.
Oceananigans.TimeSteppers.QuasiAdamsBashforth2TimeStepper
— MethodQuasiAdamsBashforth2TimeStepper(grid, tracers,
χ = 0.1;
implicit_solver = nothing,
Gⁿ = TendencyFields(grid, tracers),
G⁻ = TendencyFields(grid, tracers))
Return a 2nd-order quasi Adams-Bashforth (AB2) time stepper (QuasiAdamsBashforth2TimeStepper
) on grid
, with tracers
, and AB2 parameter χ
. The tendency fields Gⁿ
and G⁻
can be specified via optional kwargs
.
The 2nd-order quasi Adams-Bashforth timestepper steps forward the state Uⁿ
by Δt
via
Uⁿ⁺¹ = Uⁿ + Δt * [(3/2 + χ) * Gⁿ - (1/2 + χ) * Gⁿ⁻¹]
where Uⁿ
is the state at the $n$-th timestep, Gⁿ
is the tendency at the $n$-th timestep, and Gⁿ⁻¹
is the tendency at the previous timestep (G⁻
).
For the first timestep, since there are no saved tendencies from the previous timestep, the QuasiAdamsBashforth2TimeStepper
performs an Euler timestep:
Uⁿ⁺¹ = Uⁿ + Δt * Gⁿ
Oceananigans.TimeSteppers.RungeKutta3TimeStepper
— TypeRungeKutta3TimeStepper{FT, TG} <: AbstractTimeStepper
Holds parameters and tendency fields for a low storage, third-order Runge-Kutta-Wray time-stepping scheme described by Le and Moin (1991).
Oceananigans.TimeSteppers.RungeKutta3TimeStepper
— MethodRungeKutta3TimeStepper(grid, tracers;
implicit_solver = nothing,
Gⁿ = TendencyFields(grid, tracers),
G⁻ = TendencyFields(grid, tracers))
Return a 3rd-order Runge0Kutta timestepper (RungeKutta3TimeStepper
) on grid
and with tracers
. The tendency fields Gⁿ
and G⁻
can be specified via optional kwargs
.
The scheme described by Le and Moin (1991) (see H. Le, P. Moin (1991)). In a nutshel, the 3rd-order Runge Kutta timestepper steps forward the state Uⁿ
by Δt
via 3 substeps. A pressure correction step is applied after at each substep.
The state U
after each substep m
is
Uᵐ⁺¹ = Uᵐ + Δt * (γᵐ * Gᵐ + ζᵐ * Gᵐ⁻¹)
where Uᵐ
is the state at the $m$-th substep, Gᵐ
is the tendency at the $m$-th substep, Gᵐ⁻¹
is the tendency at the previous substep, and constants $γ¹ = 8/15$, $γ² = 5/12$, $γ³ = 3/4$, $ζ¹ = 0$, $ζ² = -17/60$, $ζ³ = -5/12$.
The state at the first substep is taken to be the one that corresponds to the $n$-th timestep, U¹ = Uⁿ
, and the state after the third substep is then the state at the Uⁿ⁺¹ = U⁴
.
Oceananigans.TimeSteppers.time_step!
— Methodtime_step!(model::AbstractModel{<:QuasiAdamsBashforth2TimeStepper}, Δt; euler=false)
Step forward model
one time step Δt
with a 2nd-order Adams-Bashforth method and pressure-correction substep. Setting euler=true
will take a forward Euler time step.
Oceananigans.TimeSteppers.time_step!
— Methodtime_step!(model::AbstractModel{<:RungeKutta3TimeStepper}, Δt)
Step forward model
one time step Δt
with a 3rd-order Runge-Kutta method. The 3rd-order Runge-Kutta method takes three intermediate substep stages to achieve a single timestep. A pressure correction step is applied at each intermediate stage.
Turbulence closures
Oceananigans.TurbulenceClosures.AnisotropicMinimumDissipation
— TypeAnisotropicMinimumDissipation{FT} <: AbstractTurbulenceClosure
Parameters for the "anisotropic minimum dissipation" turbulence closure for large eddy simulation proposed originally by Wybe Rozema, Hyun J. Bae, Parviz Moin, Roel Verstappen (2015) and Mahdi Abkar, Hyun J. Bae, Parviz Moin (2016), then modified by Roel Verstappen (2018), and finally described and validated for by Catherine A. Vreugdenhil, John R. Taylor (2018).
Oceananigans.TurbulenceClosures.AnisotropicMinimumDissipation
— MethodAnisotropicMinimumDissipation([time_discretization = ExplicitTimeDiscretization, FT = Float64;]
C = 1/12, Cν = nothing, Cκ = nothing, Cb = nothing)
Return parameters of type FT
for the AnisotropicMinimumDissipation
turbulence closure.
Keyword arguments
C
: Poincaré constant for both eddy viscosity and eddy diffusivities.C
is overridden for eddy viscosity or eddy diffusivity ifCν
orCκ
are set, respecitvely.Cν
: Poincaré constant for momentum eddy viscosity.Cκ
: Poincaré constant for tracer eddy diffusivities. If one number or function, the same number or function is applied to all tracers. If aNamedTuple
, it must possess a field specifying the Poncaré constant for every tracer.Cb
: Buoyancy modification multiplier (Cb = nothing
turns it off,Cb = 1
was used by Mahdi Abkar, Hyun J. Bae, Parviz Moin (2016)). Note: that we do not subtract the horizontally-average component before computing this buoyancy modification term. This implementation differs from Mahdi Abkar, Hyun J. Bae, Parviz Moin (2016)'s proposal and the impact of this approximation has not been tested or validated.time_discretization
: EitherExplicitTimeDiscretization()
orVerticallyImplicitTimeDiscretization()
, which integrates the terms involving only z-derivatives in the viscous and diffusive fluxes with an implicit time discretization.
By default: C = Cν = Cκ
= 1/12, which is appropriate for a finite-volume method employing a second-order advection scheme, Cb = nothing
, which terms off the buoyancy modification term.
Cν
or Cκ
may be constant numbers, or functions of x, y, z
.
Examples
julia> using Oceananigans
julia> pretty_diffusive_closure = AnisotropicMinimumDissipation(C=1/2)
AnisotropicMinimumDissipation{ExplicitTimeDiscretization} turbulence closure with:
Poincaré constant for momentum eddy viscosity Cν: 0.5
Poincaré constant for tracer(s) eddy diffusivit(ies) Cκ: 0.5
Buoyancy modification multiplier Cb: nothing
julia> using Oceananigans
julia> const Δz = 0.5; # grid resolution at surface
julia> surface_enhanced_tracer_C(x, y, z) = 1/12 * (1 + exp((z + Δz/2) / 8Δz));
julia> fancy_closure = AnisotropicMinimumDissipation(Cκ=surface_enhanced_tracer_C)
AnisotropicMinimumDissipation{ExplicitTimeDiscretization} turbulence closure with:
Poincaré constant for momentum eddy viscosity Cν: 0.08333333333333333
Poincaré constant for tracer(s) eddy diffusivit(ies) Cκ: surface_enhanced_tracer_C
Buoyancy modification multiplier Cb: nothing
julia> using Oceananigans
julia> tracer_specific_closure = AnisotropicMinimumDissipation(Cκ=(c₁=1/12, c₂=1/6))
AnisotropicMinimumDissipation{ExplicitTimeDiscretization} turbulence closure with:
Poincaré constant for momentum eddy viscosity Cν: 0.08333333333333333
Poincaré constant for tracer(s) eddy diffusivit(ies) Cκ: (c₁ = 0.08333333333333333, c₂ = 0.16666666666666666)
Buoyancy modification multiplier Cb: nothing
References
Vreugdenhil C., and Taylor J. (2018), "Large-eddy simulations of stratified plane Couette flow using the anisotropic minimum-dissipation model", Physics of Fluids 30, 085104.
Verstappen, R. (2018), "How much eddy dissipation is needed to counterbalance the nonlinear production of small, unresolved scales in a large-eddy simulation of turbulence?", Computers & Fluids 176, pp. 276-284.
Oceananigans.TurbulenceClosures.ConvectiveAdjustmentVerticalDiffusivity
— TypeConvectiveAdjustmentVerticalDiffusivity([time_discretization = VerticallyImplicitTimeDiscretization(), FT=Float64;]
convective_κz = 0,
convective_νz = 0,
background_κz = 0,
background_νz = 0)
Return a convective adjustment vertical diffusivity closure that applies different values of diffusivity and/or viscosity depending whether the region is statically stable (positive or zero buoyancy gradient) or statically unstable (negative buoyancy gradient).
Arguments
time_discretization
: EitherExplicitTimeDiscretization()
orVerticallyImplicitTimeDiscretization()
; defaultVerticallyImplicitTimeDiscretization()
.FT
: Float type; defaultFloat64
.
Keyword arguments
convective_κz
: Vertical tracer diffusivity in regions with negative (unstable) buoyancy gradients. Either a single number, function, array, field, or tuple of diffusivities for each tracer.background_κz
: Vertical tracer diffusivity in regions with zero or positive (stable) buoyancy gradients.convective_νz
: Vertical viscosity in regions with negative (unstable) buoyancy gradients. Either a number, function, array, or field.background_κz
: Vertical viscosity in regions with zero or positive (stable) buoyancy gradients.
Example
julia> using Oceananigans
julia> cavd = ConvectiveAdjustmentVerticalDiffusivity(convective_κz = 1)
ConvectiveAdjustmentVerticalDiffusivity{VerticallyImplicitTimeDiscretization}(background_κz=0.0 convective_κz=1 background_νz=0.0 convective_νz=0.0)
Oceananigans.TurbulenceClosures.ExplicitTimeDiscretization
— Typestruct ExplicitTimeDiscretization <: AbstractTimeDiscretization
A fully-explicit time-discretization of a TurbulenceClosure
.
Oceananigans.TurbulenceClosures.HorizontalDivergenceScalarDiffusivity
— TypeHorizontalDivergenceScalarDiffusivity([time_discretization=ExplicitTimeDiscretization(),
FT::DataType=Float64;]
kwargs...)
Shorthand for a ScalarDiffusivity
with HorizontalDivergenceFormulation()
. See ScalarDiffusivity
.
Oceananigans.TurbulenceClosures.HorizontalScalarDiffusivity
— TypeHorizontalScalarDiffusivity([time_discretization=ExplicitTimeDiscretization(),
FT::DataType=Float64;]
kwargs...)
Shorthand for a ScalarDiffusivity
with HorizontalFormulation()
. See ScalarDiffusivity
.
Oceananigans.TurbulenceClosures.IsopycnalSkewSymmetricDiffusivity
— MethodIsopycnalSkewSymmetricDiffusivity([time_disc=VerticallyImplicitTimeDiscretization(), FT=Float64;]
κ_skew = 0,
κ_symmetric = 0,
isopycnal_tensor = SmallSlopeIsopycnalTensor(),
slope_limiter = FluxTapering(1e-2))
Return parameters for an isopycnal skew-symmetric tracer diffusivity with skew diffusivity κ_skew
and symmetric diffusivity κ_symmetric
that uses an isopycnal_tensor
model for for calculating the isopycnal slopes, and (optionally) applying a slope_limiter
to the calculated isopycnal slope values.
Both κ_skew
and κ_symmetric
may be constants, arrays, fields, or functions of (x, y, z, t)
.
Oceananigans.TurbulenceClosures.RiBasedVerticalDiffusivity
— TypeRiBasedVerticalDiffusivity([time_discretization = VerticallyImplicitTimeDiscretization(),
FT = Float64;]
Ri_dependent_tapering = HyperbolicTangentRiDependentTapering(),
ν₀ = 0.7,
κ₀ = 0.5,
κᶜᵃ = 1.7,
Cᵉⁿ = 0.1,
Cᵃᵛ = 0.6,
Ri₀ = 0.1,
Riᵟ = 0.4,
warning = true)
Return a closure that estimates the vertical viscosity and diffusivity from "convective adjustment" coefficients ν₀
and κ₀
multiplied by a decreasing function of the Richardson number, $Ri$.
Arguments
time_discretization
: EitherExplicitTimeDiscretization()
orVerticallyImplicitTimeDiscretization()
, which integrates the terms involving only $z$-derivatives in the viscous and diffusive fluxes with an implicit time discretization. DefaultVerticallyImplicitTimeDiscretization()
.FT
: Float type; defaultFloat64
.
Keyword arguments
Ri_dependent_tapering
: The $Ri$-dependent tapering. Options are:PiecewiseLinearRiDependentTapering()
,HyperbolicTangentRiDependentTapering()
(default), andExponentialRiDependentTapering()
.ν₀
: Non-convective viscosity.κ₀
: Non-convective diffusivity for tracers.κᶜᵃ
: Convective adjustment diffusivity for tracers.Cᵉⁿ
: Entrainment coefficient for tracers.Cᵃᵛ
: Time-averaging coefficient for viscosity and diffusivity.Ri₀
: $Ri$ threshold for decreasing viscosity and diffusivity.Riᵟ
: $Ri$-width over which viscosity and diffusivity decreases to 0.
Oceananigans.TurbulenceClosures.ScalarBiharmonicDiffusivity
— TypeScalarBiharmonicDiffusivity([formulation=ThreeDimensionalFormulation(), FT=Float64;]
ν=0, κ=0,
discrete_form = false)
Return a scalar biharmonic diffusivity turbulence closure with viscosity coefficient ν
and tracer diffusivities κ
for each tracer field in tracers
. If a single κ
is provided, it is applied to all tracers. Otherwise κ
must be a NamedTuple
with values for every tracer individually.
Arguments
formulation
:HorizontalFormulation()
for diffusivity applied in the horizontal direction(s)VerticalFormulation()
for diffusivity applied in the vertical direction,ThreeDimensionalFormulation()
(default) for diffusivity applied isotropically to all directions
FT
: the float datatype (default:Float64
)
Keyword arguments
ν
: Viscosity.Number
,AbstractArray
,Field
, orFunction
.κ
: Diffusivity.Number
, three-dimensionalAbstractArray
,Field
,Function
, orNamedTuple
of diffusivities with entries for each tracer.discrete_form
:Boolean
; default:False
.
When prescribing the viscosities or diffusivities as functions, depending on the value of keyword argument discrete_form
, the constructor expects:
discrete_form = false
(default): functions of the grid's native coordinates and time, e.g.,(x, y, z, t)
for aRectilinearGrid
or(λ, φ, z, t)
for aLatitudeLongitudeGrid
.discrete_form = true
: functions of(i, j, k, grid, ℓx, ℓy, ℓz)
withℓx
,ℓy
andℓz
eitherFace()
orCenter()
.
For examples see ScalarDiffusivity
.
Oceananigans.TurbulenceClosures.ScalarBiharmonicDiffusivity
— Typestruct ScalarBiharmonicDiffusivity{F, N, K} <: AbstractScalarBiharmonicDiffusivity{F}
Holds viscosity and diffusivities for models with prescribed isotropic diffusivities.
Oceananigans.TurbulenceClosures.ScalarDiffusivity
— TypeScalarDiffusivity([time_discretization=ExplicitTimeDiscretization(),
formulation=ThreeDimensionalFormulation(), FT=Float64];
ν=0, κ=0,
discrete_form = false)
Return ScalarDiffusivity
turbulence closure with viscosity ν
and tracer diffusivities κ
for each tracer field in tracers
. If a single κ
is provided, it is applied to all tracers. Otherwise κ
must be a NamedTuple
with values for every tracer individually.
Arguments
time_discretization
: eitherExplicitTimeDiscretization()
(default) orVerticallyImplicitTimeDiscretization()
.formulation
:HorizontalFormulation()
for diffusivity applied in the horizontal direction(s)VerticalFormulation()
for diffusivity applied in the vertical direction,ThreeDimensionalFormulation()
(default) for diffusivity applied isotropically to all directions
FT
: the float datatype (default:Float64
)
Keyword arguments
ν
: Viscosity.Number
, three-dimensionalAbstractArray
,Field
, orFunction
.κ
: Diffusivity.Number
,AbstractArray
,Field
,Function
, orNamedTuple
of diffusivities with entries for each tracer.discrete_form
:Boolean
; default:False
.
When prescribing the viscosities or diffusivities as functions, depending on the value of keyword argument discrete_form
, the constructor expects:
discrete_form = false
(default): functions of the grid's native coordinates and time, e.g.,(x, y, z, t)
for aRectilinearGrid
or(λ, φ, z, t)
for aLatitudeLongitudeGrid
.discrete_form = true
: functions of(i, j, k, grid, ℓx, ℓy, ℓz)
withℓx
,ℓy
andℓz
eitherFace()
orCenter()
.
Examples
julia> using Oceananigans
julia> ScalarDiffusivity(ν = 1000, κ=2000)
ScalarDiffusivity{ExplicitTimeDiscretization}(ν=1000.0, κ=2000.0)
julia> const depth_scale = 100;
julia> @inline ν(x, y, z) = 1000 * exp(z / depth_scale)
ν (generic function with 1 method)
julia> ScalarDiffusivity(ν = ν)
ScalarDiffusivity{ExplicitTimeDiscretization}(ν=ν (generic function with 1 method), κ=0.0)
julia> using Oceananigans.Grids: znode
julia> @inline function κ(i, j, k, grid, ℓx, ℓy, ℓz)
z = znode(i, j, k, grid, ℓx, ℓy, ℓz)
return 2000 * exp(z / depth_scale)
end
κ (generic function with 1 method)
julia> ScalarDiffusivity(κ = κ)
ScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0, κ=κ (generic function with 1 method))
julia> ScalarDiffusivity(κ = κ, discrete_form = true)
ScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0, κ=Oceananigans.TurbulenceClosures.DiscreteDiffusionFunction{Nothing, Nothing, Nothing, Nothing, typeof(κ)})
Oceananigans.TurbulenceClosures.SmagorinskyLilly
— MethodSmagorinskyLilly([time_discretization::TD = ExplicitTimeDiscretization(), FT=Float64;] C=0.16, Cb=1.0, Pr=1.0)
Return a SmagorinskyLilly
type associated with the turbulence closure proposed by Lilly (1962) and Smagorinsky (1958, 1963), which has an eddy viscosity of the form
νₑ = (C * Δᶠ)² * √(2Σ²) * √(1 - Cb * N² / Σ²)
and an eddy diffusivity of the form
κₑ = νₑ / Pr
where Δᶠ
is the filter width, Σ² = ΣᵢⱼΣᵢⱼ
is the double dot product of the strain tensor Σᵢⱼ
, Pr
is the turbulent Prandtl number, and N²
is the total buoyancy gradient, and Cb
is a constant the multiplies the Richardson number modification to the eddy viscosity.
Arguments
time_discretization
: EitherExplicitTimeDiscretization()
orVerticallyImplicitTimeDiscretization()
, which integrates the terms involving only $z$-derivatives in the viscous and diffusive fluxes with an implicit time discretization. DefaultExplicitTimeDiscretization()
.FT
: Float type; defaultFloat64
.
Keyword arguments
C
: Smagorinsky constant. Default value is 0.16 as obtained by Lilly (1966).Cb
: Buoyancy term multipler based on Lilly (1962) (Cb = 0
turns it off,Cb ≠ 0
turns it on. Typically, and according to the original work by Lilly (1962),Cb = 1 / Pr
.)Pr
: Turbulent Prandtl numbers for each tracer. Either a constant applied to every tracer, or aNamedTuple
with fields for each tracer individually.
References
Smagorinsky, J. "On the numerical integration of the primitive equations of motion for baroclinic flow in a closed region." Monthly Weather Review (1958)
Lilly, D. K. "On the numerical simulation of buoyant convection." Tellus (1962)
Smagorinsky, J. "General circulation experiments with the primitive equations: I. The basic experiment." Monthly weather review (1963)
Lilly, D. K. "The representation of small-scale turbulence in numerical simulation experiments." NCAR Manuscript No. 281, 0, 1966.
Oceananigans.TurbulenceClosures.TwoDimensionalLeith
— TypeTwoDimensionalLeith(FT=Float64;
C=0.3, C_Redi=1, C_GM=1,
isopycnal_model=SmallSlopeIsopycnalTensor())
Return a TwoDimensionalLeith
type associated with the turbulence closure proposed by Leith (1965) and Fox-Kemper & Menemenlis (2008) which has an eddy viscosity of the form
νₑ = (C * Δᶠ)³ * √(|∇ₕ ζ|² + |∇ₕ ∂w/∂z|²)
and an eddy diffusivity of the form...
where Δᶠ
is the filter width, ζ = ∂v/∂x - ∂u/∂y
is the vertical vorticity, and C
is a model constant.
Keyword arguments
C
: Model constantC_Redi
: Coefficient for down-gradient tracer diffusivity for each tracer. Either a constant applied to every tracer, or aNamedTuple
with fields for each tracer individually.C_GM
: Coefficient for down-gradient tracer diffusivity for each tracer. Either a constant applied to every tracer, or aNamedTuple
with fields for each tracer individually.
References
Leith, C. E. (1968). "Diffusion Approximation for Two‐Dimensional Turbulence", The Physics of Fluids 11, 671. doi: 10.1063/1.1691968
Fox‐Kemper, B., & D. Menemenlis (2008), "Can large eddy simulation techniques improve mesoscale rich ocean models?", in Ocean Modeling in an Eddying Regime, Geophys. Monogr. Ser., vol. 177, pp. 319–337. doi: 10.1029/177GM19
Pearson, B. et al. (2017) , "Evaluation of scale-aware subgrid mesoscale eddy models in a global eddy rich model", Ocean Modelling 115, 42-58. doi: 10.1016/j.ocemod.2017.05.007
Oceananigans.TurbulenceClosures.VerticalScalarDiffusivity
— TypeVerticalScalarDiffusivity([time_discretization=ExplicitTimeDiscretization(),
FT::DataType=Float64;]
kwargs...)
Shorthand for a ScalarDiffusivity
with VerticalFormulation()
. See ScalarDiffusivity
.
Oceananigans.TurbulenceClosures.VerticallyImplicitTimeDiscretization
— Typestruct VerticallyImplicitTimeDiscretization <: AbstractTimeDiscretization
A vertically-implicit time-discretization of a TurbulenceClosure
.
This implies that a flux divergence such as $𝛁 ⋅ 𝐪$ at the $n$-th timestep is time-discretized as
[∇ ⋅ q]ⁿ = [explicit_flux_divergence]ⁿ + [∂z (κ ∂z c)]ⁿ⁺¹
Oceananigans.TurbulenceClosures.diffusivity
— Functiondiffusivity(closure, tracer_index, diffusivity_fields)
Returns the scalar diffusivity associated with closure
and tracer_index
.
Oceananigans.TurbulenceClosures.viscosity
— Functionviscosity(closure, diffusivities)
Returns the scalar viscosity associated with closure
.
Utilities
Oceananigans.Utils.AndSchedule
— MethodAndSchedule(schedules...)
Return a schedule that actuates when all child_schedule
s actuate.
Oceananigans.Utils.IterationInterval
— MethodIterationInterval(interval; offset=0)
Return a callable IterationInterval
that "actuates" (schedules output or callback execution) whenever the model iteration (modified by offset
) is a multiple of interval
.
For example,
IterationInterval(100)
actuates at iterations[100, 200, 300, ...]
.IterationInterval(100, offset=-1)
actuates at iterations[99, 199, 299, ...]
.
Oceananigans.Utils.MultiRegionObject
— MethodMultiRegionObject(regional_objects::Tuple; devices)
Return a MultiRegionObject
Oceananigans.Utils.OrSchedule
— MethodOrSchedule(schedules...)
Return a schedule that actuates when any of the child_schedule
s actuates.
Oceananigans.Utils.SpecifiedTimes
— MethodSpecifiedTimes(times)
Return a callable TimeInterval
that "actuates" (schedules output or callback execution) whenever the model's clock equals the specified values in times
. For example,
SpecifiedTimes([1, 15.3])
actuates whenmodel.clock.time
is1
and15.3
.
The specified times
need not be ordered as the SpecifiedTimes
constructor will check and order them in ascending order if needed.
Oceananigans.Utils.TimeInterval
— Typestruct TimeInterval <: AbstractSchedule
Callable TimeInterval
schedule for periodic output or diagnostic evaluation according to model.clock.time
.
Oceananigans.Utils.TimeInterval
— MethodTimeInterval(interval)
Return a callable TimeInterval
that schedules periodic output or diagnostic evaluation on a interval
of simulation time, as kept by model.clock
.
Oceananigans.Utils.WallTimeInterval
— MethodWallTimeInterval(interval; start_time = time_ns() * 1e-9)
Return a callable WallTimeInterval
that schedules periodic output or diagnostic evaluation on a interval
of "wall time" while a simulation runs, in units of seconds.
The "wall time" is the actual real world time in seconds, as kept by an actual or hypothetical clock hanging on your wall.
The keyword argument start_time
can be used to specify a starting wall time other than the moment WallTimeInterval
is constructed.
Oceananigans.Utils.launch!
— Methodlaunch!(arch, grid, layout, kernel!, args...; kwargs...)
Launches kernel!
, with arguments args
and keyword arguments kwargs
, over the dims
of grid
on the architecture arch
. kernels run on the defaul stream
Oceananigans.Utils.pretty_filesize
— Functionpretty_filesize(s, suffix="B")
Convert a floating point value s
representing a file size to a more human-friendly formatted string with one decimal places with a suffix
defaulting to "B". Depending on the value of s
the string will be formatted to show s
using an SI prefix from bytes, kiB (1024 bytes), MiB (1024² bytes), and so on up to YiB (1024⁸ bytes).
Oceananigans.Utils.prettytime
— Functionprettytime(t, longform=true)
Convert a floating point value t
representing an amount of time in SI units of seconds to a human-friendly string with three decimal places. Depending on the value of t
the string will be formatted to show t
in nanoseconds (ns), microseconds (μs), milliseconds (ms), seconds, minutes, hours, or days.
With longform=false
, we use s, m, hrs, and d in place of seconds, minutes, and hours.
Oceananigans.Utils.with_tracers
— Methodwith_tracers(tracer_names, initial_tuple, tracer_default)
Create a tuple corresponding to the solution variables u
, v
, w
, and tracer_names
. initial_tuple
is a NamedTuple
that at least has fields u
, v
, and w
, and may have some fields corresponding to the names in tracer_names
. tracer_default
is a function that produces a default tuple value for each tracer if not included in initial_tuple
.
Oceananigans.Utils.work_layout
— Methodwork_layout(grid, dims; include_right_boundaries=false, location=nothing)
Returns the workgroup
and worksize
for launching a kernel over dims
on grid
. The workgroup
is a tuple specifying the threads per block in each dimension. The worksize
specifies the range of the loop in each dimension.
Specifying include_right_boundaries=true
will ensure the work layout includes the right face end points along bounded dimensions. This requires the field location
to be specified.
For more information, see: https://github.com/CliMA/Oceananigans.jl/pull/308
Oceananigans.Utils.@apply_regionally
— Macro@apply_regionally expr
Use @apply_regionally
to distribute locally the function calls. Call compute_regionally
in case of a returning value and apply_regionally!
in case of no return.
Units
Oceananigans.Units.GiB
— ConstantGiB
A Float64
constant equal to 1024MiB
. Useful for increasing the clarity of scripts, e.g. max_filesize = 50GiB
.
Oceananigans.Units.KiB
— ConstantKiB
A Float64
constant equal to 1024.0. Useful for increasing the clarity of scripts, e.g. max_filesize = 250KiB
.
Oceananigans.Units.MiB
— ConstantMiB
A Float64
constant equal to 1024KiB
. Useful for increasing the clarity of scripts, e.g. max_filesize = 100MiB
.
Oceananigans.Units.TiB
— ConstantTiB
A Float64
constant equal to 1024GiB
. Useful for increasing the clarity of scripts, e.g. max_filesize = 2TiB
.
Oceananigans.Units.day
— Constantday
A Float64
constant equal to 24hours
. Useful for increasing the clarity of scripts, e.g. stop_time = 1day
.
Oceananigans.Units.days
— Constantdays
A Float64
constant equal to 24hours
. Useful for increasing the clarity of scripts, e.g. stop_time = 7days
.
Oceananigans.Units.hour
— Constanthour
A Float64
constant equal to 60minutes
. Useful for increasing the clarity of scripts, e.g. Δt = 1hour
.
Oceananigans.Units.hours
— Constanthours
A Float64
constant equal to 60minutes
. Useful for increasing the clarity of scripts, e.g. Δt = 3hours
.
Oceananigans.Units.kilometer
— Constantkilometer
A Float64
constant equal to 1000meters
. Useful for increasing the clarity of scripts, e.g. Lx = 1kilometer
.
Oceananigans.Units.kilometers
— Constantkilometers
A Float64
constant equal to 1000meters
. Useful for increasing the clarity of scripts, e.g. Lx = 5000kilometers
.
Oceananigans.Units.meter
— Constantmeter
A Float64
constant equal to 1.0. Useful for increasing the clarity of scripts, e.g. Lx = 1meter
.
Oceananigans.Units.meters
— Constantmeters
A Float64
constant equal to 1.0. Useful for increasing the clarity of scripts, e.g. Lx = 50meters
.
Oceananigans.Units.minute
— Constantminute
A Float64
constant equal to 60seconds
. Useful for increasing the clarity of scripts, e.g. Δt = 1minute
.
Oceananigans.Units.minutes
— Constantminutes
A Float64
constant equal to 60seconds
. Useful for increasing the clarity of scripts, e.g. Δt = 15minutes
.
Oceananigans.Units.second
— Constantsecond
A Float64
constant equal to 1.0. Useful for increasing the clarity of scripts, e.g. Δt = 1second
.
Oceananigans.Units.seconds
— Constantseconds
A Float64
constant equal to 1.0. Useful for increasing the clarity of scripts, e.g. Δt = 7seconds
.