Two dimensional turbulence example
In this example, we initialize a random velocity field and observe its turbulent decay in a two-dimensional domain. This example demonstrates:
- How to run a model with no buoyancy equation or tracers;
- How to create user-defined fields
- How to use differentiation functions
- How to save a computed field
Model setup
We instantiate the model with a simple isotropic diffusivity. We also use a 4-th order advection scheme and Runge-Kutta 3rd order time-stepping scheme.
using Oceananigans, Oceananigans.Advection
model = IncompressibleModel(
grid = RegularCartesianGrid(size=(128, 128, 1), halo=(2, 2, 2), extent=(2π, 2π, 2π)),
timestepper = :RungeKutta3,
advection = CenteredFourthOrder(),
buoyancy = nothing,
tracers = nothing,
closure = IsotropicDiffusivity(ν=1e-4)
)
Setting initial conditions
Our initial condition randomizes u
and v
. We also ensure that both have zero mean for purely aesthetic reasons.
using Statistics
u₀ = rand(size(model.grid)...)
u₀ .-= mean(u₀)
set!(model, u=u₀, v=u₀)
Calculating vorticity
using Oceananigans.Fields, Oceananigans.AbstractOperations
Next we create an object called an ComputedField
that calculates vorticity. We'll use this object to calculate vorticity on-line and output it as the simulation progresses.
u, v, w = model.velocities
ω = ComputedField(∂x(v) - ∂y(u))
Now we construct a simulation.
simulation = Simulation(model, Δt=0.1, stop_iteration=2000)
Simulation{IncompressibleModel{CPU, Float64}} ├── Model clock: time = 0 seconds, iteration = 0 ├── Next time step (Float64): 100 ms ├── Iteration interval: 1 ├── Stop criteria: Any[Oceananigans.Simulations.iteration_limit_exceeded, Oceananigans.Simulations.stop_time_exceeded, Oceananigans.Simulations.wall_time_limit_exceeded] ├── Run time: 0 seconds, wall time limit: Inf ├── Stop time: Inf years, stop iteration: 2000 ├── Diagnostics: OrderedCollections.OrderedDict with no entries └── Output writers: OrderedCollections.OrderedDict with no entries
Output
We set up an output writer for the simulation that saves the vorticity every 20 iterations.
using Oceananigans.OutputWriters
simulation.output_writers[:fields] =
JLD2OutputWriter(model, (ω = ω,), schedule=IterationInterval(20),
prefix = "2d_turbulence_vorticity", force = true)
JLD2OutputWriter scheduled on IterationInterval(20): ├── filepath: ./2d_turbulence_vorticity.jld2 ├── 1 outputs: (:ω,) ├── field slicer: FieldSlicer(:, :, :, with_halos=false) ├── array type: Array{Float32} ├── including: [:grid, :coriolis, :buoyancy, :closure] └── max filesize: Inf YiB
Running the simulation
Finally, we run the simulation.
run!(simulation)
[ Info: Simulation is stopping. Model iteration 2000 has hit or exceeded simulation stop iteration 2000.
Visualizing the results
We load the output and make a movie.
using JLD2
file = jldopen(simulation.output_writers[:fields].filepath)
iterations = parse.(Int, keys(file["timeseries/t"]))
Construct the $x$, $y$ grid for plotting purposes,
using Oceananigans.Grids
x, y = xnodes(ω)[:], ynodes(ω)[:]
and animate the vorticity.
using Plots
anim = @animate for iteration in iterations
ω_snapshot = file["timeseries/ω/$iteration"][:, :, 1]
heatmap(x, y, ω_snapshot',
xlabel="x", ylabel="y",
aspectratio=1,
color=:balance,
clims=(-0.2, 0.2),
xlims=(0, model.grid.Lx),
ylims=(0, model.grid.Ly))
end
This page was generated using Literate.jl.