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-9/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.9e-04 m s⁻¹, wall time: 8.808 seconds
[ Info: ... simulation initialization complete (3.233 seconds)
[ Info: Executing initial time step...
[ Info: ... initial time step complete (5.668 seconds).
Iteration: 0200, time: 26.889 minutes, Δt: 7.713 seconds, max|w|: 4.3e-04 m s⁻¹, wall time: 17.516 seconds
Iteration: 0400, time: 52.078 minutes, Δt: 7.512 seconds, max|w|: 3.3e-04 m s⁻¹, wall time: 18.956 seconds
Iteration: 0600, time: 1.282 hours, Δt: 7.443 seconds, max|w|: 2.9e-04 m s⁻¹, wall time: 20.385 seconds
Iteration: 0800, time: 1.691 hours, Δt: 7.389 seconds, max|w|: 2.7e-04 m s⁻¹, wall time: 21.822 seconds
Iteration: 1000, time: 2.100 hours, Δt: 7.365 seconds, max|w|: 2.2e-04 m s⁻¹, wall time: 23.273 seconds
Iteration: 1200, time: 2.509 hours, Δt: 7.343 seconds, max|w|: 3.0e-04 m s⁻¹, wall time: 24.742 seconds
Iteration: 1400, time: 2.915 hours, Δt: 7.333 seconds, max|w|: 7.4e-04 m s⁻¹, wall time: 26.215 seconds
Iteration: 1600, time: 3.322 hours, Δt: 7.322 seconds, max|w|: 2.6e-03 m s⁻¹, wall time: 27.690 seconds
Iteration: 1800, time: 3.728 hours, Δt: 7.380 seconds, max|w|: 5.5e-03 m s⁻¹, wall time: 29.535 seconds
Iteration: 2000, time: 4.142 hours, Δt: 7.562 seconds, max|w|: 6.2e-03 m s⁻¹, wall time: 31.021 seconds
Iteration: 2200, time: 4.557 hours, Δt: 7.443 seconds, max|w|: 5.1e-03 m s⁻¹, wall time: 32.472 seconds
Iteration: 2400, time: 4.970 hours, Δt: 7.414 seconds, max|w|: 4.6e-03 m s⁻¹, wall time: 33.917 seconds
Iteration: 2600, time: 5.378 hours, Δt: 7.297 seconds, max|w|: 4.8e-03 m s⁻¹, wall time: 35.388 seconds
Iteration: 2800, time: 5.782 hours, Δt: 7.270 seconds, max|w|: 4.5e-03 m s⁻¹, wall time: 36.847 seconds
Iteration: 3000, time: 6.183 hours, Δt: 7.241 seconds, max|w|: 5.4e-03 m s⁻¹, wall time: 38.341 seconds
Iteration: 3200, time: 6.585 hours, Δt: 7.218 seconds, max|w|: 4.9e-03 m s⁻¹, wall time: 39.827 seconds
Iteration: 3400, time: 6.986 hours, Δt: 7.251 seconds, max|w|: 4.4e-03 m s⁻¹, wall time: 41.309 seconds
Iteration: 3600, time: 7.388 hours, Δt: 7.256 seconds, max|w|: 4.0e-03 m s⁻¹, wall time: 42.802 seconds
Iteration: 3800, time: 7.790 hours, Δt: 7.267 seconds, max|w|: 6.8e-03 m s⁻¹, wall time: 44.262 seconds
Iteration: 4000, time: 8.192 hours, Δt: 7.287 seconds, max|w|: 5.3e-03 m s⁻¹, wall time: 45.752 seconds
Iteration: 4200, time: 8.597 hours, Δt: 7.328 seconds, max|w|: 5.9e-03 m s⁻¹, wall time: 47.234 seconds
Iteration: 4400, time: 9.000 hours, Δt: 7.279 seconds, max|w|: 5.2e-03 m s⁻¹, wall time: 48.711 seconds
Iteration: 4600, time: 9.397 hours, Δt: 7.234 seconds, max|w|: 5.5e-03 m s⁻¹, wall time: 50.192 seconds
Iteration: 4800, time: 9.804 hours, Δt: 7.387 seconds, max|w|: 6.9e-03 m s⁻¹, wall time: 51.666 seconds
Iteration: 5000, time: 10.213 hours, Δt: 7.366 seconds, max|w|: 5.8e-03 m s⁻¹, wall time: 53.138 seconds
Iteration: 5200, time: 10.619 hours, Δt: 7.361 seconds, max|w|: 7.0e-03 m s⁻¹, wall time: 54.609 seconds
Iteration: 5400, time: 11.028 hours, Δt: 7.271 seconds, max|w|: 5.5e-03 m s⁻¹, wall time: 56.084 seconds
Iteration: 5600, time: 11.430 hours, Δt: 7.273 seconds, max|w|: 6.6e-03 m s⁻¹, wall time: 57.550 seconds
Iteration: 5800, time: 11.836 hours, Δt: 7.324 seconds, max|w|: 6.7e-03 m s⁻¹, wall time: 59.024 seconds
Iteration: 6000, time: 12.240 hours, Δt: 7.295 seconds, max|w|: 8.4e-03 m s⁻¹, wall time: 1.008 minutes
Iteration: 6200, time: 12.641 hours, Δt: 7.190 seconds, max|w|: 5.4e-03 m s⁻¹, wall time: 1.033 minutes
Iteration: 6400, time: 13.038 hours, Δt: 7.245 seconds, max|w|: 5.0e-03 m s⁻¹, wall time: 1.058 minutes
Iteration: 6600, time: 13.435 hours, Δt: 7.151 seconds, max|w|: 5.1e-03 m s⁻¹, wall time: 1.083 minutes
Iteration: 6800, time: 13.833 hours, Δt: 7.247 seconds, max|w|: 7.6e-03 m s⁻¹, wall time: 1.107 minutes
Iteration: 7000, time: 14.235 hours, Δt: 7.218 seconds, max|w|: 8.6e-03 m s⁻¹, wall time: 1.132 minutes
Iteration: 7200, time: 14.634 hours, Δt: 7.252 seconds, max|w|: 6.4e-03 m s⁻¹, wall time: 1.157 minutes
Iteration: 7400, time: 15.032 hours, Δt: 7.190 seconds, max|w|: 7.3e-03 m s⁻¹, wall time: 1.182 minutes
Iteration: 7600, time: 15.435 hours, Δt: 7.297 seconds, max|w|: 8.1e-03 m s⁻¹, wall time: 1.206 minutes
Iteration: 7800, time: 15.838 hours, Δt: 7.285 seconds, max|w|: 7.3e-03 m s⁻¹, wall time: 1.230 minutes
Iteration: 8000, time: 16.239 hours, Δt: 7.210 seconds, max|w|: 7.2e-03 m s⁻¹, wall time: 1.254 minutes
Iteration: 8200, time: 16.637 hours, Δt: 7.164 seconds, max|w|: 6.6e-03 m s⁻¹, wall time: 1.278 minutes
Iteration: 8400, time: 17.034 hours, Δt: 7.192 seconds, max|w|: 9.4e-03 m s⁻¹, wall time: 1.303 minutes
Iteration: 8600, time: 17.433 hours, Δt: 7.234 seconds, max|w|: 8.4e-03 m s⁻¹, wall time: 1.327 minutes
Iteration: 8800, time: 17.830 hours, Δt: 7.212 seconds, max|w|: 6.2e-03 m s⁻¹, wall time: 1.351 minutes
Iteration: 9000, time: 18.227 hours, Δt: 7.215 seconds, max|w|: 7.5e-03 m s⁻¹, wall time: 1.375 minutes
Iteration: 9200, time: 18.627 hours, Δt: 7.197 seconds, max|w|: 7.0e-03 m s⁻¹, wall time: 1.400 minutes
Iteration: 9400, time: 19.024 hours, Δt: 7.270 seconds, max|w|: 7.0e-03 m s⁻¹, wall time: 1.425 minutes
Iteration: 9600, time: 19.428 hours, Δt: 7.251 seconds, max|w|: 8.1e-03 m s⁻¹, wall time: 1.450 minutes
Iteration: 9800, time: 19.827 hours, Δt: 7.164 seconds, max|w|: 7.6e-03 m s⁻¹, wall time: 1.475 minutes
Iteration: 10000, time: 20.218 hours, Δt: 7.217 seconds, max|w|: 7.2e-03 m s⁻¹, wall time: 1.499 minutes
Iteration: 10200, time: 20.617 hours, Δt: 7.172 seconds, max|w|: 7.3e-03 m s⁻¹, wall time: 1.524 minutes
Iteration: 10400, time: 21.012 hours, Δt: 6.965 seconds, max|w|: 8.0e-03 m s⁻¹, wall time: 1.549 minutes
Iteration: 10600, time: 21.407 hours, Δt: 7.129 seconds, max|w|: 8.3e-03 m s⁻¹, wall time: 1.574 minutes
Iteration: 10800, time: 21.801 hours, Δt: 6.991 seconds, max|w|: 8.0e-03 m s⁻¹, wall time: 1.598 minutes
Iteration: 11000, time: 22.179 hours, Δt: 6.628 seconds, max|w|: 7.8e-03 m s⁻¹, wall time: 1.623 minutes
Iteration: 11200, time: 22.565 hours, Δt: 7.085 seconds, max|w|: 5.6e-03 m s⁻¹, wall time: 1.648 minutes
Iteration: 11400, time: 22.961 hours, Δt: 7.121 seconds, max|w|: 5.4e-03 m s⁻¹, wall time: 1.673 minutes
Iteration: 11600, time: 23.359 hours, Δt: 7.190 seconds, max|w|: 6.2e-03 m s⁻¹, wall time: 1.698 minutes
Iteration: 11800, time: 23.757 hours, Δt: 7.207 seconds, max|w|: 7.0e-03 m s⁻¹, wall time: 1.723 minutes
[ Info: Simulation is stopping after running for 1.619 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.