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 StaticStridesOffsetArray(::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 (11.440 seconds)
[ Info: Executing initial time step...
[ Info: ... initial time step complete (1.491 seconds).
i: 20, sim time: 4.254 hours, wall time: 13.520 seconds, Δt: 16.776 minutes, CFL: 3.13e-01
i: 40, sim time: 10 hours, wall time: 14.121 seconds, Δt: 24.562 minutes, CFL: 4.66e-01
i: 60, sim time: 18 hours, wall time: 14.734 seconds, Δt: 25.280 minutes, CFL: 5.00e-01
i: 80, sim time: 1.082 days, wall time: 15.397 seconds, Δt: 23.707 minutes, CFL: 5.00e-01
i: 100, sim time: 1.350 days, wall time: 15.998 seconds, Δt: 23.859 minutes, CFL: 5.00e-01
i: 120, sim time: 1.632 days, wall time: 16.585 seconds, Δt: 23.274 minutes, CFL: 5.00e-01
i: 140, sim time: 1.911 days, wall time: 17.336 seconds, Δt: 21.946 minutes, CFL: 5.00e-01
i: 160, sim time: 2.182 days, wall time: 17.945 seconds, Δt: 21.129 minutes, CFL: 5.00e-01
i: 180, sim time: 2.444 days, wall time: 18.531 seconds, Δt: 19.091 minutes, CFL: 5.00e-01
i: 200, sim time: 2.679 days, wall time: 19.107 seconds, Δt: 16.970 minutes, CFL: 5.00e-01
i: 220, sim time: 2.888 days, wall time: 19.677 seconds, Δt: 15.357 minutes, CFL: 5.00e-01
i: 240, sim time: 3.083 days, wall time: 20.246 seconds, Δt: 13.579 minutes, CFL: 5.00e-01
i: 260, sim time: 3.258 days, wall time: 20.820 seconds, Δt: 11.634 minutes, CFL: 5.00e-01
i: 280, sim time: 3.408 days, wall time: 21.389 seconds, Δt: 10.441 minutes, CFL: 5.00e-01
i: 300, sim time: 3.542 days, wall time: 21.953 seconds, Δt: 9.585 minutes, CFL: 5.00e-01
i: 320, sim time: 3.667 days, wall time: 22.508 seconds, Δt: 9.238 minutes, CFL: 5.00e-01
i: 340, sim time: 3.786 days, wall time: 23.017 seconds, Δt: 8.458 minutes, CFL: 5.00e-01
i: 360, sim time: 3.896 days, wall time: 23.524 seconds, Δt: 8.215 minutes, CFL: 5.00e-01
i: 380, sim time: 4.006 days, wall time: 24.031 seconds, Δt: 8.101 minutes, CFL: 5.00e-01
i: 400, sim time: 4.110 days, wall time: 24.568 seconds, Δt: 7.649 minutes, CFL: 5.00e-01
i: 420, sim time: 4.213 days, wall time: 25.094 seconds, Δt: 7.328 minutes, CFL: 5.00e-01
i: 440, sim time: 4.310 days, wall time: 25.645 seconds, Δt: 7.187 minutes, CFL: 5.00e-01
i: 460, sim time: 4.410 days, wall time: 26.173 seconds, Δt: 7.514 minutes, CFL: 5.00e-01
i: 480, sim time: 4.511 days, wall time: 26.685 seconds, Δt: 7.669 minutes, CFL: 5.00e-01
i: 500, sim time: 4.616 days, wall time: 27.208 seconds, Δt: 7.961 minutes, CFL: 5.00e-01
i: 520, sim time: 4.721 days, wall time: 27.773 seconds, Δt: 7.883 minutes, CFL: 5.00e-01
i: 540, sim time: 4.827 days, wall time: 28.340 seconds, Δt: 7.854 minutes, CFL: 5.00e-01
i: 560, sim time: 4.928 days, wall time: 28.967 seconds, Δt: 7.969 minutes, CFL: 5.00e-01
i: 580, sim time: 5.040 days, wall time: 29.522 seconds, Δt: 8.241 minutes, CFL: 5.00e-01
i: 600, sim time: 5.153 days, wall time: 30.087 seconds, Δt: 8.516 minutes, CFL: 5.00e-01
i: 620, sim time: 5.262 days, wall time: 30.644 seconds, Δt: 8.634 minutes, CFL: 5.00e-01
i: 640, sim time: 5.382 days, wall time: 31.190 seconds, Δt: 8.580 minutes, CFL: 5.00e-01
i: 660, sim time: 5.498 days, wall time: 31.748 seconds, Δt: 8.308 minutes, CFL: 5.00e-01
i: 680, sim time: 5.608 days, wall time: 32.300 seconds, Δt: 8.775 minutes, CFL: 5.00e-01
i: 700, sim time: 5.727 days, wall time: 32.849 seconds, Δt: 8.687 minutes, CFL: 5.00e-01
i: 720, sim time: 5.846 days, wall time: 33.406 seconds, Δt: 9.247 minutes, CFL: 5.00e-01
i: 740, sim time: 5.967 days, wall time: 33.960 seconds, Δt: 9.005 minutes, CFL: 5.00e-01
i: 760, sim time: 6.083 days, wall time: 34.513 seconds, Δt: 9.042 minutes, CFL: 5.00e-01
i: 780, sim time: 6.203 days, wall time: 35.070 seconds, Δt: 8.595 minutes, CFL: 5.00e-01
i: 800, sim time: 6.316 days, wall time: 35.626 seconds, Δt: 7.535 minutes, CFL: 5.00e-01
i: 820, sim time: 6.413 days, wall time: 36.178 seconds, Δt: 6.694 minutes, CFL: 5.00e-01
i: 840, sim time: 6.500 days, wall time: 36.727 seconds, Δt: 6.237 minutes, CFL: 5.00e-01
i: 860, sim time: 6.583 days, wall time: 37.282 seconds, Δt: 5.910 minutes, CFL: 5.00e-01
i: 880, sim time: 6.664 days, wall time: 37.830 seconds, Δt: 5.584 minutes, CFL: 5.00e-01
i: 900, sim time: 6.739 days, wall time: 38.378 seconds, Δt: 5.412 minutes, CFL: 5.00e-01
i: 920, sim time: 6.813 days, wall time: 38.929 seconds, Δt: 5.358 minutes, CFL: 5.00e-01
i: 940, sim time: 6.885 days, wall time: 39.477 seconds, Δt: 5.413 minutes, CFL: 5.00e-01
i: 960, sim time: 6.958 days, wall time: 40.033 seconds, Δt: 5.527 minutes, CFL: 5.00e-01
i: 980, sim time: 7.035 days, wall time: 40.585 seconds, Δt: 5.660 minutes, CFL: 5.00e-01
i: 1000, sim time: 7.111 days, wall time: 41.136 seconds, Δt: 5.716 minutes, CFL: 5.00e-01
i: 1020, sim time: 7.191 days, wall time: 41.684 seconds, Δt: 5.679 minutes, CFL: 5.00e-01
i: 1040, sim time: 7.265 days, wall time: 42.239 seconds, Δt: 5.426 minutes, CFL: 5.00e-01
i: 1060, sim time: 7.337 days, wall time: 42.788 seconds, Δt: 5.186 minutes, CFL: 5.00e-01
i: 1080, sim time: 7.409 days, wall time: 43.332 seconds, Δt: 5.159 minutes, CFL: 5.00e-01
i: 1100, sim time: 7.477 days, wall time: 43.887 seconds, Δt: 5.055 minutes, CFL: 5.00e-01
i: 1120, sim time: 7.545 days, wall time: 44.442 seconds, Δt: 4.879 minutes, CFL: 5.00e-01
i: 1140, sim time: 7.611 days, wall time: 45.003 seconds, Δt: 4.869 minutes, CFL: 5.00e-01
i: 1160, sim time: 7.677 days, wall time: 45.560 seconds, Δt: 4.938 minutes, CFL: 5.00e-01
i: 1180, sim time: 7.745 days, wall time: 46.109 seconds, Δt: 4.811 minutes, CFL: 5.00e-01
i: 1200, sim time: 7.810 days, wall time: 46.661 seconds, Δt: 4.736 minutes, CFL: 5.00e-01
i: 1220, sim time: 7.873 days, wall time: 47.217 seconds, Δt: 4.840 minutes, CFL: 5.00e-01
i: 1240, sim time: 7.941 days, wall time: 47.773 seconds, Δt: 5.074 minutes, CFL: 5.00e-01
i: 1260, sim time: 8.011 days, wall time: 48.327 seconds, Δt: 5.320 minutes, CFL: 5.00e-01
i: 1280, sim time: 8.083 days, wall time: 48.899 seconds, Δt: 5.755 minutes, CFL: 5.00e-01
i: 1300, sim time: 8.166 days, wall time: 49.458 seconds, Δt: 6.170 minutes, CFL: 5.00e-01
i: 1320, sim time: 8.249 days, wall time: 50.008 seconds, Δt: 6.415 minutes, CFL: 5.00e-01
i: 1340, sim time: 8.333 days, wall time: 50.564 seconds, Δt: 6.305 minutes, CFL: 5.00e-01
i: 1360, sim time: 8.421 days, wall time: 51.127 seconds, Δt: 6.887 minutes, CFL: 5.00e-01
i: 1380, sim time: 8.514 days, wall time: 51.681 seconds, Δt: 6.719 minutes, CFL: 5.00e-01
i: 1400, sim time: 8.601 days, wall time: 52.243 seconds, Δt: 6.141 minutes, CFL: 5.00e-01
i: 1420, sim time: 8.683 days, wall time: 52.800 seconds, Δt: 6.043 minutes, CFL: 5.00e-01
i: 1440, sim time: 8.768 days, wall time: 53.354 seconds, Δt: 6.361 minutes, CFL: 5.00e-01
i: 1460, sim time: 8.851 days, wall time: 53.906 seconds, Δt: 6.192 minutes, CFL: 5.00e-01
i: 1480, sim time: 8.934 days, wall time: 54.453 seconds, Δt: 6.146 minutes, CFL: 5.00e-01
i: 1500, sim time: 9.017 days, wall time: 54.999 seconds, Δt: 6.105 minutes, CFL: 5.00e-01
i: 1520, sim time: 9.100 days, wall time: 55.557 seconds, Δt: 6.203 minutes, CFL: 5.00e-01
i: 1540, sim time: 9.184 days, wall time: 56.114 seconds, Δt: 6.372 minutes, CFL: 5.00e-01
i: 1560, sim time: 9.272 days, wall time: 56.668 seconds, Δt: 6.359 minutes, CFL: 5.00e-01
i: 1580, sim time: 9.360 days, wall time: 57.218 seconds, Δt: 6.218 minutes, CFL: 5.00e-01
i: 1600, sim time: 9.442 days, wall time: 57.768 seconds, Δt: 6.117 minutes, CFL: 5.00e-01
i: 1620, sim time: 9.526 days, wall time: 58.318 seconds, Δt: 6.200 minutes, CFL: 5.00e-01
i: 1640, sim time: 9.609 days, wall time: 58.873 seconds, Δt: 6.193 minutes, CFL: 5.00e-01
i: 1660, sim time: 9.692 days, wall time: 59.423 seconds, Δt: 6.071 minutes, CFL: 5.00e-01
i: 1680, sim time: 9.775 days, wall time: 59.977 seconds, Δt: 5.839 minutes, CFL: 5.00e-01
i: 1700, sim time: 9.853 days, wall time: 1.009 minutes, Δt: 5.456 minutes, CFL: 5.00e-01
i: 1720, sim time: 9.924 days, wall time: 1.018 minutes, Δt: 5.230 minutes, CFL: 5.00e-01
i: 1740, sim time: 9.997 days, wall time: 1.027 minutes, Δt: 5.334 minutes, CFL: 5.00e-01
[ Info: Simulation is stopping after running for 1.028 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.