Skip to content

Buoyancy models and equations of state ​

The buoyancy option selects how buoyancy is treated in NonhydrostaticModels and HydrostaticFreeSurfaceModels (ShallowWaterModels do not have that option given the physics of the model). There are currently three alternatives:

  1. No buoyancy (and no gravity).

  2. Evolve buoyancy as a tracer.

  3. Seawater buoyancy: evolve temperature and salinity as tracers with a value for the gravitational acceleration and an equation of state of your choosing.

No buoyancy ​

To turn off buoyancy (and gravity) you can simply pass buoyancy = nothing to the model constructor. For example to create a NonhydrostaticModel:

julia
julia> grid = RectilinearGrid(size=(8, 8, 8), extent=(1, 1, 1));

julia> model = NonhydrostaticModel(grid)
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   └── momentum: Centered(order=2)
├── tracers: ()
├── closure: Nothing
├── buoyancy: Nothing
└── coriolis: Nothing

The option buoyancy = nothing is the default for NonhydrostaticModel, so omitting the buoyancy keyword argument from the NonhydrostaticModel constructor yields the same:

julia
julia> model = NonhydrostaticModel(grid)
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   └── momentum: Centered(order=2)
├── tracers: ()
├── closure: Nothing
├── buoyancy: Nothing
└── coriolis: Nothing

The same is true for HydrostaticFreeSurfaceModel,

julia
julia> model = HydrostaticFreeSurfaceModel(grid)
HydrostaticFreeSurfaceModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: QuasiAdamsBashforth2TimeStepper
├── tracers: ()
├── closure: Nothing
├── buoyancy: Nothing
├── free surface: ImplicitFreeSurface with gravitational acceleration 9.80665 m s⁻²
│   └── solver: FFTImplicitFreeSurfaceSolver
├── advection scheme:
│   └── momentum: VectorInvariant
├── vertical_coordinate: ZCoordinate
└── coriolis: Nothing

Buoyancy as a tracer ​

Both NonhydrostaticModel and HydrostaticFreeSurfaceModel support evolving a buoyancy tracer by including :b in tracers and specifying buoyancy = BuoyancyTracer():

julia
julia> model = NonhydrostaticModel(grid; buoyancy=BuoyancyTracer(), tracers=:b)
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   ├── momentum: Centered(order=2)
│   └── b: Centered(order=2)
├── tracers: b
├── closure: Nothing
├── buoyancy: BuoyancyTracer with ĝ = NegativeZDirection()
└── coriolis: Nothing

Similarly for a HydrostaticFreeSurfaceModel with buoyancy as a tracer:

julia
julia> model = HydrostaticFreeSurfaceModel(grid; buoyancy=BuoyancyTracer(), tracers=:b)
HydrostaticFreeSurfaceModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: QuasiAdamsBashforth2TimeStepper
├── tracers: b
├── closure: Nothing
├── buoyancy: BuoyancyTracer with ĝ = NegativeZDirection()
├── free surface: ImplicitFreeSurface with gravitational acceleration 9.80665 m s⁻²
│   └── solver: FFTImplicitFreeSurfaceSolver
├── advection scheme:
│   ├── momentum: VectorInvariant
│   └── b: Centered(order=2)
├── vertical_coordinate: ZCoordinate
└── coriolis: Nothing

Seawater buoyancy ​

NonhydrostaticModel and HydrostaticFreeSurfaceModel support modeling the buoyancy of seawater as a function of the gravitational acceleration, the conservative temperature , and the absolute salinity . The relationship between , , the geopotential height, and the density perturbation from a reference value is called the equation_of_state.

Specifying buoyancy = SeawaterBuoyancy() returns a buoyancy model with a linear equation of state, Earth standard gravitational_acceleration = 9.80665 (in S.I. units  ) and requires to add :T and :S as tracers:

julia
julia> model = NonhydrostaticModel(grid; buoyancy=SeawaterBuoyancy(), tracers=(:T, :S))
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   ├── momentum: Centered(order=2)
│   ├── T: Centered(order=2)
│   └── S: Centered(order=2)
├── tracers: (T, S)
├── closure: Nothing
├── buoyancy: SeawaterBuoyancy with g=9.80665 and LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078) with ĝ = NegativeZDirection()
└── coriolis: Nothing

and the same is true for HydrostaticFreeSurfaceModel,

julia
julia> model = HydrostaticFreeSurfaceModel(grid; buoyancy=SeawaterBuoyancy(), tracers=(:T, :S))
HydrostaticFreeSurfaceModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: QuasiAdamsBashforth2TimeStepper
├── tracers: (T, S)
├── closure: Nothing
├── buoyancy: SeawaterBuoyancy with g=9.80665 and LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078) with ĝ = NegativeZDirection()
├── free surface: ImplicitFreeSurface with gravitational acceleration 9.80665 m s⁻²
│   └── solver: FFTImplicitFreeSurfaceSolver
├── advection scheme:
│   ├── momentum: VectorInvariant
│   ├── T: Centered(order=2)
│   └── S: Centered(order=2)
├── vertical_coordinate: ZCoordinate
└── coriolis: Nothing

To model flows near the surface of Europa where gravitational_acceleration = 1.3  , we might alternatively specify

julia
julia> buoyancy = SeawaterBuoyancy(gravitational_acceleration=1.3)
SeawaterBuoyancy{Float64}:
├── gravitational_acceleration: 1.3
└── equation_of_state: LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078)

julia> model = NonhydrostaticModel(grid; buoyancy, tracers=(:T, :S))
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   ├── momentum: Centered(order=2)
│   ├── T: Centered(order=2)
│   └── S: Centered(order=2)
├── tracers: (T, S)
├── closure: Nothing
├── buoyancy: SeawaterBuoyancy with g=1.3 and LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078) with ĝ = NegativeZDirection()
└── coriolis: Nothing

for example.

Linear equation of state ​

To specify the thermal expansion and haline contraction coefficients     and   ,

julia
julia> buoyancy = SeawaterBuoyancy(equation_of_state=LinearEquationOfState(thermal_expansion=2e-3, haline_contraction=5e-4))
SeawaterBuoyancy{Float64}:
├── gravitational_acceleration: 9.80665
└── equation_of_state: LinearEquationOfState(thermal_expansion=0.002, haline_contraction=0.0005)

Idealized nonlinear equations of state ​

Instead of a linear equation of state, six idealized (second-order) nonlinear equations of state as described by DocumenterCitations.CitationSiteNode("Roquet15Idealized-cite-1")

may be used. These equations of state are provided via the SeawaterPolynomials.jl package.

julia
julia> using SeawaterPolynomials.SecondOrderSeawaterPolynomials

julia> eos = RoquetEquationOfState(:Freezing)
BoussinesqEquationOfState{Float64}:
├── seawater_polynomial: SecondOrderSeawaterPolynomial{Float64}
└── reference_density: 1024.6

julia> eos.seawater_polynomial # the density anomaly
ρ' = 0.7718 Sᴬ - 0.0491 Θ - 0.005027 Θ² - 2.5681e-5 Θ Z + 0.0 Sᴬ² + 0.0 Sᴬ Z + 0.0 Sᴬ Θ

julia> buoyancy = SeawaterBuoyancy(equation_of_state=eos)
SeawaterBuoyancy{Float64}:
├── gravitational_acceleration: 9.80665
└── equation_of_state: BoussinesqEquationOfState{Float64}

TEOS-10 equation of state ​

A high-accuracy 55-term polynomial approximation to the TEOS-10 equation of state suitable for use in Boussinesq models as described by DocumenterCitations.CitationSiteNode("Roquet15TEOS-cite-1")

is implemented in the SeawaterPolynomials.jl package and may be used.

julia
julia> using SeawaterPolynomials.TEOS10

julia> eos = TEOS10EquationOfState()
BoussinesqEquationOfState{Float64}:
├── seawater_polynomial: TEOS10SeawaterPolynomial{Float64}
└── reference_density: 1020.0

The direction of gravitational acceleration ​

To simulate gravitational accelerations that don't align with the vertical (z) coordinate, we use BuoyancyForce(formulation; gravity_unit_vector), wherein the buoyancy formulation can be BuoyancyTracer, SeawaterBuoyancy, etc, in addition to the gravity_unit_vector. For example,

julia
julia> grid = RectilinearGrid(size=(8, 8, 8), extent=(1, 1, 1));

julia> θ = 45; # degrees

julia> g̃ = (0, sind(θ), cosd(θ));

julia> buoyancy = BuoyancyForce(BuoyancyTracer(), gravity_unit_vector=g̃)
BuoyancyForce:
├── formulation: BuoyancyTracer
└── gravity_unit_vector: (0.0, 0.707107, 0.707107)

julia> model = NonhydrostaticModel(grid; buoyancy, tracers=:b)
NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 8×8×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│   ├── momentum: Centered(order=2)
│   └── b: Centered(order=2)
├── tracers: b
├── closure: Nothing
├── buoyancy: BuoyancyTracer with ĝ = (0.0, 0.707107, 0.707107)
└── coriolis: Nothing