Biogeochemistry in submesoscale eddies in the Eady model
In this example we setup a 3D model with a constant background buoyancy gradient with associated thermal wind (the Eady model) with the LOBSTER biogeochemical model. This example demonstrates how to use biogeochemistry in a more complicated physical model. The parameters in this example correspond roughly to those used by Taylor (2016) and result to the generation of a single submesoscale eddy.
Install dependencies
First we ensure we have the required dependencies installed
using Pkg
pkg "add OceanBioME, Oceananigans, CairoMakie"Model setup
We load the required packages. Although not required, we also set the random seed to ensure reproducibility of the results.
using OceanBioME, Oceananigans, Printf, CairoMakie
using Oceananigans.Units
using Random
Random.seed!(11)Random.TaskLocalRNG()Construct a grid with uniform grid spacing on CPU
grid = RectilinearGrid(size = (32, 32, 8),
extent = (1kilometer, 1kilometer, 100meters))32×32×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── Periodic x ∈ [0.0, 1000.0) regularly spaced with Δx=31.25
├── Periodic y ∈ [0.0, 1000.0) regularly spaced with Δy=31.25
└── Bounded z ∈ [-100.0, 0.0] regularly spaced with Δz=12.5Set the Coriolis and buoyancy models.
coriolis = FPlane(f = 1e-4) # [s⁻¹]
buoyancy = SeawaterBuoyancy()SeawaterBuoyancy{Float64}:
├── gravitational_acceleration: 9.80665
└── equation_of_state: LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078)Specify parameters that are used to construct the background state.
background_state_parameters = (M = 1e-4, # s⁻¹, geostrophic shear
f = coriolis.f, # s⁻¹, Coriolis parameter
N = 1e-4, # s⁻¹, buoyancy frequency
H = grid.Lz,
g = buoyancy.gravitational_acceleration,
α = buoyancy.equation_of_state.thermal_expansion)(M = 0.0001, f = 0.0001, N = 0.0001, H = 100.0, g = 9.80665, α = 0.000167)We assume a background buoyancy $B$ with a constant stratification and also a constant lateral gradient (in the zonal direction). The background velocity components $U$ and $V$ are prescribed so that the thermal wind relationship is satisfied, that is, $f \partial_z U = - \partial_y B$ and $f \partial_z V = \partial_x B$.
T(x, y, z, t, p) = (p.M^2 * x + p.N^2 * (z + p.H)) / (p.g * p.α)
V(x, y, z, t, p) = p.M^2 / p.f * (z + p.H)
V_field = BackgroundField(V, parameters = background_state_parameters)
T_field = BackgroundField(T, parameters = background_state_parameters)BackgroundField{typeof(Main.var"##142".T), @NamedTuple{M::Float64, f::Float64, N::Float64, H::Float64, g::Float64, α::Float64}}
├── func: T (generic function with 1 method)
└── parameters: (M = 0.0001, f = 0.0001, N = 0.0001, H = 100.0, g = 9.80665, α = 0.000167)Specify some horizontal and vertical viscosity and diffusivity.
νᵥ = κᵥ = 1e-4 # [m² s⁻¹]
vertical_diffusivity = VerticalScalarDiffusivity(ν = νᵥ, κ = κᵥ)VerticalScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0001, κ=0.0001)Setup the biogeochemical model with optional carbonate chemistry turned on.
biogeochemistry = LOBSTER(grid;
inorganic_carbon = CarbonateSystem(),
scale_negatives = true)
DIC_bcs = FieldBoundaryConditions(top = CarbonDioxideGasExchangeBoundaryCondition())Oceananigans.FieldBoundaryConditions, with boundary conditions
├── west: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)
├── east: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)
├── bottom: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)
├── top: FluxBoundaryCondition: DiscreteBoundaryFunction with GasExchange
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: nothing)Model instantiation
model = NonhydrostaticModel(grid;
biogeochemistry,
boundary_conditions = (DIC = DIC_bcs, ),
advection = WENO(),
coriolis,
tracers = (:T, :S),
buoyancy,
background_fields = (T = T_field, v = V_field),
closure = vertical_diffusivity)NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 32×32×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme:
│ ├── momentum: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── T: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── S: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── NO₃: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── NH₄: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── P: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── Z: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── DOM: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── sPOM: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── bPOM: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ ├── DIC: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
│ └── Alk: WENO{3, Float64, Oceananigans.Utils.BackendOptimizedDivision}(order=5)
├── tracers: (T, S, NO₃, NH₄, P, Z, DOM, sPOM, bPOM, DIC, Alk)
├── closure: VerticalScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0001, κ=(T=0.0001, S=0.0001, NO₃=0.0001, NH₄=0.0001, P=0.0001, Z=0.0001, DOM=0.0001, sPOM=0.0001, bPOM=0.0001, DIC=0.0001, Alk=0.0001))
├── buoyancy: SeawaterBuoyancy with g=9.80665 and LinearEquationOfState(thermal_expansion=0.000167, haline_contraction=0.00078) with ĝ = NegativeZDirection()
└── coriolis: FPlane{Float64}(f=0.0001)Initial conditions
Start with a bit of random noise added to the background thermal wind and an arbitary biogeochemical state.
Ξ(z) = randn() * z / grid.Lz * (z / grid.Lz + 1)
Ũ = 1e-3
uᵢ(x, y, z) = Ũ * Ξ(z)
vᵢ(x, y, z) = Ũ * Ξ(z)
set!(model, u=uᵢ, v=vᵢ, P = 0.03, Z = 0.03, NO₃ = 4.0, NH₄ = 0.05, DIC = 2200.0, Alk = 2409.0, S = 35, T = 20)Setup the simulation
Choose an appropriate initial timestep for this resolution and set up the simulation
Δx = minimum_xspacing(grid, Center(), Center(), Center())
Δy = minimum_yspacing(grid, Center(), Center(), Center())
Δz = minimum_zspacing(grid, Center(), Center(), Center())
Δt₀ = 0.5 * min(Δx, Δy, Δz) / V(0, 0, 0, 0, background_state_parameters)
simulation = Simulation(model, Δt = Δt₀, stop_time = 10days)Simulation of NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── Next time step: 10.417 minutes
├── run_wall_time: 0 seconds
├── run_wall_time / iteration: NaN days
├── stop_time: 10 days
├── stop_iteration: Inf
├── wall_time_limit: Inf
├── minimum_relative_step: 0.0
├── callbacks: OrderedDict with 4 entries:
│ ├── stop_time_exceeded => Callback of stop_time_exceeded on IterationInterval(1)
│ ├── stop_iteration_exceeded => Callback of stop_iteration_exceeded on IterationInterval(1)
│ ├── wall_time_limit_exceeded => Callback of wall_time_limit_exceeded on IterationInterval(1)
│ └── nan_checker => Callback of NaNChecker for u on IterationInterval(100)
└── output_writers: OrderedDict with no entriesAdapt the time step while keeping the CFL number fixed.
wizard = TimeStepWizard(cfl = 0.5, diffusive_cfl = 0.5, max_Δt = 30minutes)
simulation.callbacks[:wizard] = Callback(wizard, IterationInterval(5))Create a progress message.
progress(sim) = @printf("i: % 6d, sim time: % 10s, wall time: % 10s, Δt: % 10s, CFL: %.2e\n",
sim.model.clock.iteration,
prettytime(sim.model.clock.time),
prettytime(sim.run_wall_time),
prettytime(sim.Δt),
AdvectiveCFL(sim.Δt)(sim.model))
simulation.callbacks[:progress] = Callback(progress, IterationInterval(20))Callback of progress on IterationInterval(20)Here, we add some diagnostics to calculate and output.
u, v, w = model.velocities # unpack velocity `Field`sNamedTuple with 3 Fields on 32×32×8 RectilinearGrid{Float64, Periodic, Periodic, Bounded} on CPU with 3×3×3 halo:
├── u: 32×32×8 Field{Face, Center, Center} on RectilinearGrid on CPU
├── v: 32×32×8 Field{Center, Face, Center} on RectilinearGrid on CPU
└── w: 32×32×9 Field{Center, Center, Face} on RectilinearGrid on CPUand also calculate the vertical vorticity.
ζ = Field(∂x(v) - ∂y(u))32×32×8 Field{Face, Face, Center} on RectilinearGrid on CPU
├── grid: 32×32×8 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: Nothing
├── operand: BinaryOperation at (Face, Face, Center)
├── status: time=0.0
└── data: 38×38×14 OffsetArray(::Array{Float64, 3}, -2:35, -2:35, -2:11) with eltype Float64 with indices -2:35×-2:35×-2:11
└── max=4.88217e-5, min=-4.88463e-5, mean=7.65142e-24Periodically save the velocities and vorticity to a file.
simulation.output_writers[:fields] = JLD2Writer(model, merge(model.tracers, (; u, v, w, ζ));
schedule = TimeInterval(2hours),
filename = "eady_turbulence_bgc",
overwrite_files = true)
Run the simulation
run!(simulation)[ Info: Initializing simulation...
i: 0, sim time: 0 seconds, wall time: 0 seconds, Δt: 11.458 minutes, CFL: 2.17e-01
[ Info: ... simulation initialization complete (15.761 seconds)
[ Info: Executing initial time step...
[ Info: ... initial time step complete (1.655 seconds).
i: 20, sim time: 4.254 hours, wall time: 18.050 seconds, Δt: 16.776 minutes, CFL: 3.13e-01
i: 40, sim time: 10 hours, wall time: 18.700 seconds, Δt: 24.562 minutes, CFL: 4.66e-01
i: 60, sim time: 18 hours, wall time: 19.383 seconds, Δt: 25.213 minutes, CFL: 5.00e-01
i: 80, sim time: 1.082 days, wall time: 20.086 seconds, Δt: 23.606 minutes, CFL: 5.00e-01
i: 100, sim time: 1.350 days, wall time: 20.763 seconds, Δt: 23.733 minutes, CFL: 5.00e-01
i: 120, sim time: 1.632 days, wall time: 21.424 seconds, Δt: 23.224 minutes, CFL: 5.00e-01
i: 140, sim time: 1.911 days, wall time: 22.066 seconds, Δt: 21.906 minutes, CFL: 5.00e-01
i: 160, sim time: 2.182 days, wall time: 22.731 seconds, Δt: 21.011 minutes, CFL: 5.00e-01
i: 180, sim time: 2.443 days, wall time: 23.379 seconds, Δt: 18.832 minutes, CFL: 5.00e-01
i: 200, sim time: 2.679 days, wall time: 24.018 seconds, Δt: 16.676 minutes, CFL: 5.00e-01
i: 220, sim time: 2.887 days, wall time: 24.630 seconds, Δt: 14.988 minutes, CFL: 5.00e-01
i: 240, sim time: 3.078 days, wall time: 25.255 seconds, Δt: 13.396 minutes, CFL: 5.00e-01
i: 260, sim time: 3.244 days, wall time: 25.886 seconds, Δt: 11.870 minutes, CFL: 5.00e-01
i: 280, sim time: 3.393 days, wall time: 26.520 seconds, Δt: 10.587 minutes, CFL: 5.00e-01
i: 300, sim time: 3.527 days, wall time: 27.150 seconds, Δt: 9.913 minutes, CFL: 5.00e-01
i: 320, sim time: 3.655 days, wall time: 27.784 seconds, Δt: 9.371 minutes, CFL: 5.00e-01
i: 340, sim time: 3.775 days, wall time: 28.577 seconds, Δt: 8.650 minutes, CFL: 5.00e-01
i: 360, sim time: 3.892 days, wall time: 29.183 seconds, Δt: 8.412 minutes, CFL: 5.00e-01
i: 380, sim time: 4 days, wall time: 29.785 seconds, Δt: 8.320 minutes, CFL: 5.00e-01
i: 400, sim time: 4.111 days, wall time: 30.406 seconds, Δt: 7.907 minutes, CFL: 5.00e-01
i: 420, sim time: 4.214 days, wall time: 31.007 seconds, Δt: 7.480 minutes, CFL: 5.00e-01
i: 440, sim time: 4.318 days, wall time: 31.608 seconds, Δt: 7.526 minutes, CFL: 5.00e-01
i: 460, sim time: 4.422 days, wall time: 32.228 seconds, Δt: 7.871 minutes, CFL: 5.00e-01
i: 480, sim time: 4.527 days, wall time: 32.831 seconds, Δt: 8.097 minutes, CFL: 5.00e-01
i: 500, sim time: 4.639 days, wall time: 33.424 seconds, Δt: 7.984 minutes, CFL: 5.00e-01
i: 520, sim time: 4.750 days, wall time: 34.020 seconds, Δt: 8.011 minutes, CFL: 5.00e-01
i: 540, sim time: 4.861 days, wall time: 34.645 seconds, Δt: 8.154 minutes, CFL: 5.00e-01
i: 560, sim time: 4.974 days, wall time: 35.240 seconds, Δt: 8.200 minutes, CFL: 5.00e-01
i: 580, sim time: 5.083 days, wall time: 35.835 seconds, Δt: 8.269 minutes, CFL: 5.00e-01
i: 600, sim time: 5.196 days, wall time: 36.460 seconds, Δt: 8.354 minutes, CFL: 5.00e-01
i: 620, sim time: 5.310 days, wall time: 37.059 seconds, Δt: 8.920 minutes, CFL: 5.00e-01
i: 640, sim time: 5.422 days, wall time: 37.679 seconds, Δt: 8.347 minutes, CFL: 5.00e-01
i: 660, sim time: 5.534 days, wall time: 38.272 seconds, Δt: 8.411 minutes, CFL: 5.00e-01
i: 680, sim time: 5.649 days, wall time: 38.863 seconds, Δt: 8.843 minutes, CFL: 5.00e-01
i: 700, sim time: 5.768 days, wall time: 39.481 seconds, Δt: 8.812 minutes, CFL: 5.00e-01
i: 720, sim time: 5.891 days, wall time: 40.080 seconds, Δt: 9.249 minutes, CFL: 5.00e-01
i: 740, sim time: 6.006 days, wall time: 40.701 seconds, Δt: 9.183 minutes, CFL: 5.00e-01
i: 760, sim time: 6.128 days, wall time: 41.302 seconds, Δt: 8.891 minutes, CFL: 5.00e-01
i: 780, sim time: 6.244 days, wall time: 41.890 seconds, Δt: 8.231 minutes, CFL: 5.00e-01
i: 800, sim time: 6.350 days, wall time: 42.514 seconds, Δt: 7.611 minutes, CFL: 5.00e-01
i: 820, sim time: 6.452 days, wall time: 43.114 seconds, Δt: 6.961 minutes, CFL: 5.00e-01
i: 840, sim time: 6.541 days, wall time: 43.706 seconds, Δt: 6.371 minutes, CFL: 5.00e-01
i: 860, sim time: 6.625 days, wall time: 44.299 seconds, Δt: 5.858 minutes, CFL: 5.00e-01
i: 880, sim time: 6.702 days, wall time: 44.889 seconds, Δt: 5.618 minutes, CFL: 5.00e-01
i: 900, sim time: 6.777 days, wall time: 45.488 seconds, Δt: 5.491 minutes, CFL: 5.00e-01
i: 920, sim time: 6.848 days, wall time: 46.084 seconds, Δt: 5.380 minutes, CFL: 5.00e-01
i: 940, sim time: 6.920 days, wall time: 46.681 seconds, Δt: 5.594 minutes, CFL: 5.00e-01
i: 960, sim time: 6.999 days, wall time: 47.255 seconds, Δt: 5.838 minutes, CFL: 5.00e-01
i: 980, sim time: 7.079 days, wall time: 47.853 seconds, Δt: 6.216 minutes, CFL: 5.00e-01
i: 1000, sim time: 7.166 days, wall time: 48.446 seconds, Δt: 6.438 minutes, CFL: 5.00e-01
i: 1020, sim time: 7.250 days, wall time: 49.044 seconds, Δt: 6.484 minutes, CFL: 5.00e-01
i: 1040, sim time: 7.333 days, wall time: 49.639 seconds, Δt: 6.003 minutes, CFL: 5.00e-01
i: 1060, sim time: 7.414 days, wall time: 50.238 seconds, Δt: 5.493 minutes, CFL: 5.00e-01
i: 1080, sim time: 7.488 days, wall time: 50.832 seconds, Δt: 5.218 minutes, CFL: 5.00e-01
i: 1100, sim time: 7.557 days, wall time: 51.433 seconds, Δt: 5.007 minutes, CFL: 5.00e-01
i: 1120, sim time: 7.624 days, wall time: 52.031 seconds, Δt: 4.821 minutes, CFL: 5.00e-01
i: 1140, sim time: 7.690 days, wall time: 52.630 seconds, Δt: 4.893 minutes, CFL: 5.00e-01
i: 1160, sim time: 7.757 days, wall time: 53.226 seconds, Δt: 5.068 minutes, CFL: 5.00e-01
i: 1180, sim time: 7.825 days, wall time: 53.803 seconds, Δt: 4.714 minutes, CFL: 5.00e-01
i: 1200, sim time: 7.888 days, wall time: 54.396 seconds, Δt: 4.634 minutes, CFL: 5.00e-01
i: 1220, sim time: 7.951 days, wall time: 54.993 seconds, Δt: 4.362 minutes, CFL: 5.00e-01
i: 1240, sim time: 8.009 days, wall time: 55.584 seconds, Δt: 4.225 minutes, CFL: 5.00e-01
i: 1260, sim time: 8.068 days, wall time: 56.159 seconds, Δt: 4.295 minutes, CFL: 5.00e-01
i: 1280, sim time: 8.126 days, wall time: 56.757 seconds, Δt: 4.490 minutes, CFL: 5.00e-01
i: 1300, sim time: 8.186 days, wall time: 57.351 seconds, Δt: 4.547 minutes, CFL: 5.00e-01
i: 1320, sim time: 8.249 days, wall time: 57.926 seconds, Δt: 4.673 minutes, CFL: 5.00e-01
i: 1340, sim time: 8.312 days, wall time: 58.519 seconds, Δt: 4.785 minutes, CFL: 5.00e-01
i: 1360, sim time: 8.377 days, wall time: 59.118 seconds, Δt: 4.861 minutes, CFL: 5.00e-01
i: 1380, sim time: 8.444 days, wall time: 59.711 seconds, Δt: 4.913 minutes, CFL: 5.00e-01
i: 1400, sim time: 8.510 days, wall time: 1.005 minutes, Δt: 5.043 minutes, CFL: 5.00e-01
i: 1420, sim time: 8.581 days, wall time: 1.015 minutes, Δt: 5.255 minutes, CFL: 5.00e-01
i: 1440, sim time: 8.655 days, wall time: 1.025 minutes, Δt: 5.641 minutes, CFL: 5.00e-01
i: 1460, sim time: 8.731 days, wall time: 1.034 minutes, Δt: 5.806 minutes, CFL: 5.00e-01
i: 1480, sim time: 8.812 days, wall time: 1.044 minutes, Δt: 5.957 minutes, CFL: 5.00e-01
i: 1500, sim time: 8.890 days, wall time: 1.054 minutes, Δt: 5.752 minutes, CFL: 5.00e-01
i: 1520, sim time: 8.967 days, wall time: 1.064 minutes, Δt: 5.507 minutes, CFL: 5.00e-01
i: 1540, sim time: 9.042 days, wall time: 1.074 minutes, Δt: 5.519 minutes, CFL: 5.00e-01
i: 1560, sim time: 9.117 days, wall time: 1.084 minutes, Δt: 5.435 minutes, CFL: 5.00e-01
i: 1580, sim time: 9.189 days, wall time: 1.094 minutes, Δt: 5.387 minutes, CFL: 5.00e-01
i: 1600, sim time: 9.261 days, wall time: 1.104 minutes, Δt: 5.404 minutes, CFL: 5.00e-01
i: 1620, sim time: 9.333 days, wall time: 1.113 minutes, Δt: 5.404 minutes, CFL: 5.00e-01
i: 1640, sim time: 9.408 days, wall time: 1.123 minutes, Δt: 5.412 minutes, CFL: 5.00e-01
i: 1660, sim time: 9.481 days, wall time: 1.133 minutes, Δt: 5.496 minutes, CFL: 5.00e-01
i: 1680, sim time: 9.554 days, wall time: 1.143 minutes, Δt: 5.575 minutes, CFL: 5.00e-01
i: 1700, sim time: 9.631 days, wall time: 1.153 minutes, Δt: 5.960 minutes, CFL: 5.00e-01
i: 1720, sim time: 9.712 days, wall time: 1.163 minutes, Δt: 5.975 minutes, CFL: 5.00e-01
i: 1740, sim time: 9.791 days, wall time: 1.173 minutes, Δt: 5.856 minutes, CFL: 5.00e-01
i: 1760, sim time: 9.869 days, wall time: 1.183 minutes, Δt: 5.720 minutes, CFL: 5.00e-01
i: 1780, sim time: 9.944 days, wall time: 1.193 minutes, Δt: 5.459 minutes, CFL: 5.00e-01
[ Info: Simulation is stopping after running for 1.200 minutes.
[ Info: Simulation time 10 days equals or exceeds stop time 10 days.
Now load the saved output,
ζ = FieldTimeSeries("eady_turbulence_bgc.jld2", "ζ")
P = FieldTimeSeries("eady_turbulence_bgc.jld2", "P")
NO₃ = FieldTimeSeries("eady_turbulence_bgc.jld2", "NO₃")
NH₄ = FieldTimeSeries("eady_turbulence_bgc.jld2", "NH₄")
DIC = FieldTimeSeries("eady_turbulence_bgc.jld2", "DIC")
times = ζ.times
xζ, yζ, zζ = nodes(ζ)
xc, yc, zc = nodes(P)and build the frames,
n = Observable(1)
Nz = grid.Nz
ζₙ = @lift interior( ζ[$n], :, :, Nz)
Nₙ = @lift interior(NO₃[$n], :, :, Nz) .+ interior(NH₄[$n], :, :, Nz)
Pₙ = @lift interior( P[$n], :, :, Nz)
DICₙ = @lift interior(DIC[$n], :, :, Nz)Observable([2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0; 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0 2200.0])
now plot
fig = Figure(size = (1600, 1600), fontsize = 20)
lims = [(minimum(T), maximum(T)) for T in ( ζ[:, :, Nz, :],
NO₃[:, :, Nz, :] .+ NH₄[:, :, Nz, :],
P[:, :, Nz, :],
DIC[:, :, Nz, :])]
axis_kwargs = (xlabel = "x (m)", ylabel = "y (m)", aspect = DataAspect())
ax1 = Axis(fig[1, 1]; title = "Vertical vorticity (1 / s)", axis_kwargs...)
hm1 = heatmap!(ax1, xζ, yζ, ζₙ, colormap = :balance, colorrange = lims[1])
Colorbar(fig[1, 2], hm1)
ax2 = Axis(fig[1, 3]; title = "Nutrient (NO₃ + NH₄) concentration (mmol N / m³)", axis_kwargs...)
hm2 = heatmap!(ax2, xc, yc, Nₙ, colormap = Reverse(:bamako), colorrange = lims[2])
Colorbar(fig[1, 4], hm2)
ax3 = Axis(fig[2, 1]; title = "Phytoplankton concentration (mmol N / m³)", axis_kwargs...)
hm3 = heatmap!(ax3, xc, yc, Pₙ, colormap = Reverse(:batlow), colorrange = lims[3])
Colorbar(fig[2, 2], hm3)
ax4 = Axis(fig[2, 3]; title = "Dissolved inorganic carbon (mmol C / m³)", axis_kwargs...)
hm4 = heatmap!(ax4, xc, yc, DICₙ, colormap = Reverse(:devon), colorrange = lims[4])
Colorbar(fig[2, 4], hm4)
title = @lift "t = $(prettytime(times[$n]))"
Label(fig[0, :], title, fontsize = 30)Label()and record the movie
record(fig, "eady.mp4", 1:length(times), framerate = 12) do i
n[] = i
end"eady.mp4"This page was generated using Literate.jl.