Tilted bottom boundary layer example
This example simulates a two-dimensional oceanic bottom boundary layer in a domain that's tilted with respect to gravity. We simulate the perturbation away from a constant along-slope (y-direction) velocity constant density stratification. This perturbation develops into a turbulent bottom boundary layer due to momentum loss at the bottom boundary modeled with a quadratic drag law.
This example illustrates
- changing the direction of gravitational acceleration in the buoyancy model;
- changing the axis of rotation for Coriolis forces.
Install dependencies
First let's make sure we have all required packages installed.
using Pkg
pkg"add Oceananigans, NCDatasets, CairoMakie"The domain
We create a grid with finer resolution near the bottom,
using Oceananigans
using Oceananigans.Units
Lx = 200meters
Lz = 100meters
Nx = 64
Nz = 64
# Creates a grid with near-constant spacing `refinement * Lz / Nz`
# near the bottom:
refinement = 1.8 # controls spacing near surface (higher means finer spaced)
stretching = 10 # controls rate of stretching at bottom
# "Warped" height coordinate
h(k) = (Nz + 1 - k) / Nz
# Linear near-surface generator
ζ(k) = 1 + (h(k) - 1) / refinement
# Bottom-intensified stretching function
Σ(k) = (1 - exp(-stretching * h(k))) / (1 - exp(-stretching))
# Generating function
z_faces(k) = - Lz * (ζ(k) * Σ(k) - 1)
grid = RectilinearGrid(topology = (Periodic, Flat, Bounded),
size = (Nx, Nz),
x = (0, Lx),
z = z_faces)64×1×64 RectilinearGrid{Float64, Periodic, Flat, Bounded} on CPU with 3×0×3 halo
├── Periodic x ∈ [0.0, 200.0) regularly spaced with Δx=3.125
├── Flat y
└── Bounded z ∈ [-0.0, 100.0] variably spaced with min(Δz)=0.868817, max(Δz)=6.55496Let's make sure the grid spacing is both finer and near-uniform at the bottom,
using CairoMakie
lines(zspacings(grid, Center()), znodes(grid, Center()),
axis = (ylabel = "Depth (m)",
xlabel = "Vertical spacing (m)"))
scatter!(zspacings(grid, Center()), znodes(grid, Center()))
Tilting the domain
We use a domain that's tilted with respect to gravity by
θ = 3 # degrees3so that $x$ is the along-slope direction, $z$ is the across-slope direction that is perpendicular to the bottom, and the unit vector anti-aligned with gravity is
ĝ = [sind(θ), 0, cosd(θ)]3-element Vector{Float64}:
0.052335956242943835
0.0
0.9986295347545738Changing the vertical direction impacts both the gravity_unit_vector for Buoyancy as well as the rotation_axis for Coriolis forces,
buoyancy = Buoyancy(model = BuoyancyTracer(), gravity_unit_vector = -ĝ)
coriolis = ConstantCartesianCoriolis(f = 1e-4, rotation_axis = ĝ)ConstantCartesianCoriolis{Float64}: fx = 5.23e-06, fy = 0.00e+00, fz = 9.99e-05where above we used a constant Coriolis parameter $f = 10^{-4} \, \rm{s}^{-1}$. The tilting also affects the kind of density stratified flows we can model. In particular, a constant density stratification in the tilted coordinate system
@inline constant_stratification(x, z, t, p) = p.N² * (x * p.ĝ[1] + z * p.ĝ[3])constant_stratification (generic function with 1 method)is not periodic in $x$. Thus we cannot explicitly model a constant stratification on an $x$-periodic grid such as the one used here. Instead, we simulate periodic perturbations away from the constant density stratification by imposing a constant stratification as a BackgroundField,
N² = 1e-5 # s⁻² # background vertical buoyancy gradient
B∞_field = BackgroundField(constant_stratification, parameters=(; ĝ, N² = N²))BackgroundField{typeof(Main.var"##269".constant_stratification), @NamedTuple{ĝ::Vector{Float64}, N²::Float64}}
├── func: constant_stratification (generic function with 1 method)
└── parameters: (ĝ = [0.052335956242943835, 0.0, 0.9986295347545738], N² = 1.0e-5)We choose to impose a bottom boundary condition of zero total diffusive buoyancy flux across the seafloor,
\[∂_z B = ∂_z b + N^{2} \cos{\theta} = 0.\]
This shows that to impose a no-flux boundary condition on the total buoyancy field $B$, we must apply a boundary condition to the perturbation buoyancy $b$, ```math ∂_z b = - N^{2} \cos{\theta}.
#```
∂z_b_bottom = - N² * cosd(θ)
negative_background_diffusive_flux = GradientBoundaryCondition(∂z_b_bottom)
b_bcs = FieldBoundaryConditions(bottom = negative_background_diffusive_flux)Oceananigans.FieldBoundaryConditions, with boundary conditions
├── west: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── east: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── bottom: GradientBoundaryCondition: -9.9863e-6
├── top: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)Bottom drag
We impose bottom drag that follows Monin–Obukhov theory. We include the background flow in the drag calculation, which is the only effect the background flow enters the problem,
V∞ = 0.1 # m s⁻¹
z₀ = 0.1 # m (roughness length)
κ = 0.4 # von Karman constant
z₁ = first(znodes(grid, Center())) # Closest grid center to the bottom
cᴰ = (κ / log(z₁ / z₀))^2 # Drag coefficient
@inline drag_u(x, t, u, v, p) = - p.cᴰ * √(u^2 + (v + p.V∞)^2) * u
@inline drag_v(x, t, u, v, p) = - p.cᴰ * √(u^2 + (v + p.V∞)^2) * (v + p.V∞)
drag_bc_u = FluxBoundaryCondition(drag_u, field_dependencies=(:u, :v), parameters=(; cᴰ, V∞))
drag_bc_v = FluxBoundaryCondition(drag_v, field_dependencies=(:u, :v), parameters=(; cᴰ, V∞))
u_bcs = FieldBoundaryConditions(bottom = drag_bc_u)
v_bcs = FieldBoundaryConditions(bottom = drag_bc_v)Oceananigans.FieldBoundaryConditions, with boundary conditions
├── west: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── east: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── bottom: FluxBoundaryCondition: ContinuousBoundaryFunction drag_v at (Nothing, Nothing, Nothing)
├── top: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)Create the NonhydrostaticModel
We are now ready to create the model. We create a NonhydrostaticModel with an UpwindBiasedFifthOrder advection scheme, a RungeKutta3 timestepper, and a constant viscosity and diffusivity. Here we use a smallish value of $10^{-4} \, \rm{m}^2\, \rm{s}^{-1}$.
ν = 1e-4
κ = 1e-4
closure = ScalarDiffusivity(ν=ν, κ=κ)
model = NonhydrostaticModel(; grid, buoyancy, coriolis, closure,
timestepper = :RungeKutta3,
advection = UpwindBiasedFifthOrder(),
tracers = :b,
boundary_conditions = (u=u_bcs, v=v_bcs, b=b_bcs),
background_fields = (; b=B∞_field))NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── grid: 64×1×64 RectilinearGrid{Float64, Periodic, Flat, Bounded} on CPU with 3×0×3 halo
├── timestepper: RungeKutta3TimeStepper
├── advection scheme: Upwind Biased reconstruction order 5
├── tracers: b
├── closure: ScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0001, κ=(b=0.0001,))
├── buoyancy: BuoyancyTracer with ĝ = Tuple{Float64, Float64, Float64}
└── coriolis: ConstantCartesianCoriolis{Float64}Let's introduce a bit of random noise at the bottom of the domain to speed up the onset of turbulence:
noise(x, z) = 1e-3 * randn() * exp(-(10z)^2 / grid.Lz^2)
set!(model, u=noise, w=noise)Create and run a simulation
We are now ready to create the simulation. We begin by setting the initial time step conservatively, based on the smallest grid size of our domain and either an advective or diffusive time scaling, depending on which is shorter.
Δt₀ = 0.5 * minimum([minimum_zspacing(grid) / V∞, minimum_zspacing(grid)^2/κ])
simulation = Simulation(model, Δt = Δt₀, stop_time = 1day)Simulation of NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── Next time step: 4.344 seconds
├── Elapsed wall time: 0 seconds
├── Wall time per iteration: NaN days
├── Stop time: 1 day
├── Stop iteration : Inf
├── Wall time limit: Inf
├── 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 entries
└── Diagnostics: OrderedDict with no entriesWe use a TimeStepWizard to adapt our time-step,
wizard = TimeStepWizard(max_change=1.1, cfl=0.7)
simulation.callbacks[:wizard] = Callback(wizard, IterationInterval(4))Callback of TimeStepWizard(cfl=0.7, max_Δt=Inf, min_Δt=0.0) on IterationInterval(4)and also we add another callback to print a progress message,
using Printf
start_time = time_ns() # so we can print the total elapsed wall time
progress_message(sim) =
@printf("Iteration: %04d, time: %s, Δt: %s, max|w|: %.1e m s⁻¹, wall time: %s\n",
iteration(sim), prettytime(time(sim)),
prettytime(sim.Δt), maximum(abs, sim.model.velocities.w),
prettytime((time_ns() - start_time) * 1e-9))
simulation.callbacks[:progress] = Callback(progress_message, IterationInterval(200))Callback of progress_message on IterationInterval(200)Add outputs to the simulation
We add outputs to our model using the NetCDFOutputWriter,
u, v, w = model.velocities
b = model.tracers.b
B∞ = model.background_fields.tracers.b
B = b + B∞
V = v + V∞
ωy = ∂z(u) - ∂x(w)
outputs = (; u, V, w, B, ωy)
simulation.output_writers[:fields] = NetCDFOutputWriter(model, outputs;
filename = joinpath(@__DIR__, "tilted_bottom_boundary_layer.nc"),
schedule = TimeInterval(20minutes),
overwrite_existing = true)NetCDFOutputWriter scheduled on TimeInterval(20 minutes):
├── filepath: /var/lib/buildkite-agent/builds/tartarus-8/clima/oceananigans/docs/src/literated/tilted_bottom_boundary_layer.nc
├── dimensions: zC(64), zF(65), xC(64), yF(1), xF(64), yC(1), time(0)
├── 5 outputs: (B, w, ωy, V, u)
└── array type: Array{Float64}
├── file_splitting: NoFileSplitting
└── file size: 19.2 KiBNow we just run it!
run!(simulation)[ Info: Initializing simulation...
Iteration: 0000, time: 0 seconds, Δt: 4.778 seconds, max|w|: 9.6e-04 m s⁻¹, wall time: 10.326 seconds
[ Info: ... simulation initialization complete (3.743 seconds)
[ Info: Executing initial time step...
[ Info: ... initial time step complete (6.613 seconds).
Iteration: 0200, time: 26.895 minutes, Δt: 7.728 seconds, max|w|: 4.1e-04 m s⁻¹, wall time: 20.449 seconds
Iteration: 0400, time: 52.205 minutes, Δt: 7.513 seconds, max|w|: 3.2e-04 m s⁻¹, wall time: 22.093 seconds
Iteration: 0600, time: 1.284 hours, Δt: 7.449 seconds, max|w|: 2.3e-04 m s⁻¹, wall time: 23.743 seconds
Iteration: 0800, time: 1.695 hours, Δt: 7.390 seconds, max|w|: 2.7e-04 m s⁻¹, wall time: 25.360 seconds
Iteration: 1000, time: 2.104 hours, Δt: 7.366 seconds, max|w|: 2.3e-04 m s⁻¹, wall time: 27.053 seconds
Iteration: 1200, time: 2.513 hours, Δt: 7.339 seconds, max|w|: 3.6e-04 m s⁻¹, wall time: 28.714 seconds
Iteration: 1400, time: 2.919 hours, Δt: 7.330 seconds, max|w|: 7.7e-04 m s⁻¹, wall time: 30.545 seconds
Iteration: 1600, time: 3.325 hours, Δt: 7.323 seconds, max|w|: 2.7e-03 m s⁻¹, wall time: 32.254 seconds
Iteration: 1800, time: 3.732 hours, Δt: 7.433 seconds, max|w|: 5.9e-03 m s⁻¹, wall time: 33.851 seconds
Iteration: 2000, time: 4.150 hours, Δt: 7.567 seconds, max|w|: 5.1e-03 m s⁻¹, wall time: 35.471 seconds
Iteration: 2200, time: 4.567 hours, Δt: 7.472 seconds, max|w|: 4.2e-03 m s⁻¹, wall time: 37.089 seconds
Iteration: 2400, time: 4.978 hours, Δt: 7.348 seconds, max|w|: 4.5e-03 m s⁻¹, wall time: 38.725 seconds
Iteration: 2600, time: 5.384 hours, Δt: 7.328 seconds, max|w|: 3.9e-03 m s⁻¹, wall time: 40.360 seconds
Iteration: 2800, time: 5.788 hours, Δt: 7.237 seconds, max|w|: 4.6e-03 m s⁻¹, wall time: 41.968 seconds
Iteration: 3000, time: 6.190 hours, Δt: 7.263 seconds, max|w|: 4.8e-03 m s⁻¹, wall time: 43.593 seconds
Iteration: 3200, time: 6.591 hours, Δt: 7.259 seconds, max|w|: 4.7e-03 m s⁻¹, wall time: 45.239 seconds
Iteration: 3400, time: 6.994 hours, Δt: 7.273 seconds, max|w|: 3.6e-03 m s⁻¹, wall time: 46.896 seconds
Iteration: 3600, time: 7.394 hours, Δt: 7.262 seconds, max|w|: 3.5e-03 m s⁻¹, wall time: 48.531 seconds
Iteration: 3800, time: 7.796 hours, Δt: 7.266 seconds, max|w|: 3.7e-03 m s⁻¹, wall time: 50.132 seconds
Iteration: 4000, time: 8.199 hours, Δt: 7.315 seconds, max|w|: 4.5e-03 m s⁻¹, wall time: 51.776 seconds
Iteration: 4200, time: 8.604 hours, Δt: 7.317 seconds, max|w|: 7.2e-03 m s⁻¹, wall time: 53.387 seconds
Iteration: 4400, time: 9.010 hours, Δt: 7.339 seconds, max|w|: 9.8e-03 m s⁻¹, wall time: 55.030 seconds
Iteration: 4600, time: 9.417 hours, Δt: 7.373 seconds, max|w|: 7.2e-03 m s⁻¹, wall time: 56.633 seconds
Iteration: 4800, time: 9.831 hours, Δt: 7.492 seconds, max|w|: 4.9e-03 m s⁻¹, wall time: 58.245 seconds
Iteration: 5000, time: 10.246 hours, Δt: 7.507 seconds, max|w|: 3.2e-03 m s⁻¹, wall time: 59.860 seconds
Iteration: 5200, time: 10.657 hours, Δt: 7.379 seconds, max|w|: 3.2e-03 m s⁻¹, wall time: 1.024 minutes
Iteration: 5400, time: 11.065 hours, Δt: 7.358 seconds, max|w|: 3.5e-03 m s⁻¹, wall time: 1.051 minutes
Iteration: 5600, time: 11.472 hours, Δt: 7.314 seconds, max|w|: 3.9e-03 m s⁻¹, wall time: 1.078 minutes
Iteration: 5800, time: 11.874 hours, Δt: 7.250 seconds, max|w|: 4.2e-03 m s⁻¹, wall time: 1.105 minutes
Iteration: 6000, time: 12.279 hours, Δt: 7.234 seconds, max|w|: 5.1e-03 m s⁻¹, wall time: 1.131 minutes
Iteration: 6200, time: 12.685 hours, Δt: 7.272 seconds, max|w|: 4.1e-03 m s⁻¹, wall time: 1.158 minutes
Iteration: 6400, time: 13.088 hours, Δt: 7.351 seconds, max|w|: 3.8e-03 m s⁻¹, wall time: 1.184 minutes
Iteration: 6600, time: 13.491 hours, Δt: 7.276 seconds, max|w|: 4.9e-03 m s⁻¹, wall time: 1.211 minutes
Iteration: 6800, time: 13.894 hours, Δt: 7.228 seconds, max|w|: 5.3e-03 m s⁻¹, wall time: 1.238 minutes
Iteration: 7000, time: 14.292 hours, Δt: 7.223 seconds, max|w|: 7.6e-03 m s⁻¹, wall time: 1.265 minutes
Iteration: 7200, time: 14.691 hours, Δt: 7.269 seconds, max|w|: 9.1e-03 m s⁻¹, wall time: 1.291 minutes
Iteration: 7400, time: 15.092 hours, Δt: 7.228 seconds, max|w|: 7.1e-03 m s⁻¹, wall time: 1.317 minutes
Iteration: 7600, time: 15.492 hours, Δt: 7.188 seconds, max|w|: 6.8e-03 m s⁻¹, wall time: 1.344 minutes
Iteration: 7800, time: 15.889 hours, Δt: 7.221 seconds, max|w|: 5.0e-03 m s⁻¹, wall time: 1.371 minutes
Iteration: 8000, time: 16.290 hours, Δt: 7.264 seconds, max|w|: 5.6e-03 m s⁻¹, wall time: 1.398 minutes
Iteration: 8200, time: 16.691 hours, Δt: 7.327 seconds, max|w|: 5.7e-03 m s⁻¹, wall time: 1.426 minutes
Iteration: 8400, time: 17.094 hours, Δt: 7.181 seconds, max|w|: 5.6e-03 m s⁻¹, wall time: 1.452 minutes
Iteration: 8600, time: 17.492 hours, Δt: 7.137 seconds, max|w|: 4.9e-03 m s⁻¹, wall time: 1.479 minutes
Iteration: 8800, time: 17.889 hours, Δt: 7.132 seconds, max|w|: 4.8e-03 m s⁻¹, wall time: 1.505 minutes
Iteration: 9000, time: 18.284 hours, Δt: 7.175 seconds, max|w|: 6.9e-03 m s⁻¹, wall time: 1.531 minutes
Iteration: 9200, time: 18.681 hours, Δt: 7.184 seconds, max|w|: 6.9e-03 m s⁻¹, wall time: 1.558 minutes
Iteration: 9400, time: 19.068 hours, Δt: 7.133 seconds, max|w|: 7.1e-03 m s⁻¹, wall time: 1.585 minutes
Iteration: 9600, time: 19.468 hours, Δt: 7.253 seconds, max|w|: 8.9e-03 m s⁻¹, wall time: 1.611 minutes
Iteration: 9800, time: 19.865 hours, Δt: 7.134 seconds, max|w|: 6.3e-03 m s⁻¹, wall time: 1.637 minutes
Iteration: 10000, time: 20.259 hours, Δt: 7.138 seconds, max|w|: 7.1e-03 m s⁻¹, wall time: 1.664 minutes
Iteration: 10200, time: 20.657 hours, Δt: 7.164 seconds, max|w|: 6.6e-03 m s⁻¹, wall time: 1.690 minutes
Iteration: 10400, time: 21.052 hours, Δt: 7.157 seconds, max|w|: 6.9e-03 m s⁻¹, wall time: 1.717 minutes
Iteration: 10600, time: 21.449 hours, Δt: 7.184 seconds, max|w|: 7.4e-03 m s⁻¹, wall time: 1.743 minutes
Iteration: 10800, time: 21.848 hours, Δt: 7.160 seconds, max|w|: 6.6e-03 m s⁻¹, wall time: 1.770 minutes
Iteration: 11000, time: 22.246 hours, Δt: 7.212 seconds, max|w|: 9.1e-03 m s⁻¹, wall time: 1.796 minutes
Iteration: 11200, time: 22.645 hours, Δt: 7.178 seconds, max|w|: 9.9e-03 m s⁻¹, wall time: 1.823 minutes
Iteration: 11400, time: 23.042 hours, Δt: 7.176 seconds, max|w|: 6.6e-03 m s⁻¹, wall time: 1.850 minutes
Iteration: 11600, time: 23.442 hours, Δt: 7.237 seconds, max|w|: 5.5e-03 m s⁻¹, wall time: 1.876 minutes
Iteration: 11800, time: 23.842 hours, Δt: 7.121 seconds, max|w|: 5.7e-03 m s⁻¹, wall time: 1.902 minutes
[ Info: Simulation is stopping after running for 1.773 minutes.
[ Info: Simulation time 1 day equals or exceeds stop time 1 day.
Visualize the results
First we load the required package to load NetCDF output files and define the coordinates for plotting using existing objects:
using NCDatasets, CairoMakie
xb, yb, zb = nodes(B)
xω, yω, zω = nodes(ωy)
xv, yv, zv = nodes(V)([1.5625, 4.6875, 7.8125, 10.9375, 14.0625, 17.1875, 20.3125, 23.4375, 26.5625, 29.6875, 32.8125, 35.9375, 39.0625, 42.1875, 45.3125, 48.4375, 51.5625, 54.6875, 57.8125, 60.9375, 64.0625, 67.1875, 70.3125, 73.4375, 76.5625, 79.6875, 82.8125, 85.9375, 89.0625, 92.1875, 95.3125, 98.4375, 101.5625, 104.6875, 107.8125, 110.9375, 114.0625, 117.1875, 120.3125, 123.4375, 126.5625, 129.6875, 132.8125, 135.9375, 139.0625, 142.1875, 145.3125, 148.4375, 151.5625, 154.6875, 157.8125, 160.9375, 164.0625, 167.1875, 170.3125, 173.4375, 176.5625, 179.6875, 182.8125, 185.9375, 189.0625, 192.1875, 195.3125, 198.4375], StepRangeLen(1.0, 0.0, 1), [0.4344083608847693, 1.303282217470314, 2.1722793473188617, 3.041419168936904, 3.9107241431416084, 4.780220246692868, 5.649937519111903, 6.519910693894487, 7.390179927024476, 8.26079163764491, 9.131799477986252, 10.00326545222724, 10.875261206921866, 11.747869519021842, 12.621186011421164, 13.495321130420056, 14.370402424632783, 15.24657717074116, 16.124015398230675, 17.00291337296054, 17.88349760825156, 18.766029482284154, 19.650810552161538, 20.538188668213625, 21.428565007224726, 22.322402160523495, 23.22023343257072, 24.122673528152244, 25.030430831905008, 25.944321513097137, 26.865285721823064, 27.794406180595843, 28.732929518320184, 29.682290742481506, 30.64414130083751, 31.62038124678409, 32.61319609381579, 33.62509902514256, 34.65897921569778, 35.71815712673464, 36.80644774933731, 37.92823290397964, 39.08854385039015, 40.29315562720595, 41.548694726134066, 42.86276191262984, 44.24407223660885, 45.702614534695954, 47.24983301231516, 48.89883381190526, 50.66461982505285, 52.56435739360966, 54.61767896995204, 56.84702627115933, 59.27803896724765, 61.93999449018695, 64.86630513774502, 68.09507927238374, 71.66975407682737, 75.63980801853194, 80.06156188607441, 84.99907797852119, 90.52516773625328, 96.72251877439135])Read in the simulation's output_writer for the two-dimensional fields and then create an animation showing the $y$-component of vorticity.
ds = NCDataset(simulation.output_writers[:fields].filepath, "r")
fig = Figure(size = (800, 600))
axis_kwargs = (xlabel = "Across-slope distance (m)",
ylabel = "Slope-normal\ndistance (m)",
limits = ((0, Lx), (0, Lz)),
)
ax_ω = Axis(fig[2, 1]; title = "Along-slope vorticity", axis_kwargs...)
ax_v = Axis(fig[3, 1]; title = "Along-slope velocity (v)", axis_kwargs...)
n = Observable(1)
ωy = @lift ds["ωy"][:, 1, :, $n]
B = @lift ds["B"][:, 1, :, $n]
hm_ω = heatmap!(ax_ω, xω, zω, ωy, colorrange = (-0.015, +0.015), colormap = :balance)
Colorbar(fig[2, 2], hm_ω; label = "s⁻¹")
ct_b = contour!(ax_ω, xb, zb, B, levels=-1e-3:0.5e-4:1e-3, color=:black)
V = @lift ds["V"][:, 1, :, $n]
V_max = @lift maximum(abs, ds["V"][:, 1, :, $n])
hm_v = heatmap!(ax_v, xv, zv, V, colorrange = (-V∞, +V∞), colormap = :balance)
Colorbar(fig[3, 2], hm_v; label = "m s⁻¹")
ct_b = contour!(ax_v, xb, zb, B, levels=-1e-3:0.5e-4:1e-3, color=:black)
times = collect(ds["time"])
title = @lift "t = " * string(prettytime(times[$n]))
fig[1, :] = Label(fig, title, fontsize=20, tellwidth=false)
fig
Finally, we record a movie.
frames = 1:length(times)
record(fig, "tilted_bottom_boundary_layer.mp4", frames, framerate=12) do i
n[] = i
endDon't forget to close the NetCDF file!
close(ds)closed DatasetThis page was generated using Literate.jl.