Edit on GitHub

Two dimensional turbulence example

In this example, we initialize a random velocity field and observe its viscous, 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

For this example, we need PyPlot for plotting and Statistics for setting up a random initial condition with zero mean velocity.

using Oceananigans, PyPlot, Statistics

In addition to importing plotting and statistics packages, we import some types and functions from Oceananigans that will aid in the calculation and visualization of voriticty.

using Oceananigans: Face, Cell
using Oceananigans.TurbulenceClosures: ∂x_faa, ∂y_afa

We instantiate the model with a simple isotropic diffusivity

model = Model(
        grid = RegularCartesianGrid(N=(128, 128, 1), L=(2π, 2π, 2π)),
    buoyancy = nothing,
     tracers = nothing,
     closure = ConstantIsotropicDiffusivity(ν=1e-3, κ=1e-3)
)
Model{Oceananigans.AdamsBashforthTimestepper{Float64,NamedTuple{(:u, :v, :w),Tuple{Field{Oceananigans.Face,Oceananigans.Cell,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Face,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Cell,Oceananigans.Face,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}}},ConstantIsotropicDiffusivity{Float64,NamedTuple{(),Tuple{}}},CPU,RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}},Float64,Nothing,Nothing,NamedTuple{(:u, :v, :w),Tuple{Field{Oceananigans.Face,Oceananigans.Cell,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Face,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Cell,Oceananigans.Face,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}},NamedTuple{(),Tuple{}},NamedTuple{(:pHY′, :pNHS),Tuple{Field{Oceananigans.Cell,Oceananigans.Cell,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Cell,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}},NamedTuple{(:u, :v, :w),Tuple{typeof(Oceananigans.zeroforcing),typeof(Oceananigans.zeroforcing),typeof(Oceananigans.zeroforcing)}},NamedTuple{(:solution, :tendency, :pressure),Tuple{NamedTuple{(:u, :v, :w),Tuple{NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}}},NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}}},NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Oceananigans.NoPenetration,Nothing},BoundaryCondition{Oceananigans.NoPenetration,Nothing}}}}}},NamedTuple{(:u, :v, :w),Tuple{NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}}},NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}}},NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Oceananigans.NoPenetration,Nothing},BoundaryCondition{Oceananigans.NoPenetration,Nothing}}}}}},NamedTuple{(:x, :y, :z),Tuple{CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}},CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}}}}},Oceananigans.PoissonSolverCPU{Oceananigans.PPN,Array{Float64,3},Array{Complex{Float64},3},FFTW.cFFTWPlan{Complex{Float64},-1,true,3},FFTW.r2rFFTWPlan{Complex{Float64},(5,),true,3},AbstractFFTs.ScaledPlan{Complex{Float64},FFTW.cFFTWPlan{Complex{Float64},1,true,3},Float64},FFTW.r2rFFTWPlan{Complex{Float64},(4,),true,3}},Nothing,OrderedCollections.OrderedDict{Symbol,Oceananigans.AbstractOutputWriter},OrderedCollections.OrderedDict{Symbol,Oceananigans.AbstractDiagnostic},Nothing}(CPU(), RegularCartesianGrid{Float64}
  resolution (Nx, Ny, Nz) = (128, 128, 1)
   halo size (Hx, Hy, Hz) = (1, 1, 1)
      domain (Lx, Ly, Lz) = (6.283185307179586, 6.283185307179586, 6.283185307179586)
grid spacing (Δx, Δy, Δz) = (0.04908738521234052, 0.04908738521234052, 6.283185307179586), Clock{Float64}(0.0, 0), nothing, nothing, (u = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], v = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], w = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]), NamedTuple(), (pHY′ = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], pNHS = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]), (u = Oceananigans.zeroforcing, v = Oceananigans.zeroforcing, w = Oceananigans.zeroforcing), ConstantIsotropicDiffusivity{Float64,NamedTuple{(),Tuple{}}}(0.001, NamedTuple()), (solution = (u = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}(BoundaryCondition{Flux,Nothing}(nothing), BoundaryCondition{Flux,Nothing}(nothing))), v = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}(BoundaryCondition{Flux,Nothing}(nothing), BoundaryCondition{Flux,Nothing}(nothing))), w = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Oceananigans.NoPenetration,Nothing},BoundaryCondition{Oceananigans.NoPenetration,Nothing}}(BoundaryCondition{Oceananigans.NoPenetration,Nothing}(nothing), BoundaryCondition{Oceananigans.NoPenetration,Nothing}(nothing)))), tendency = (u = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}(BoundaryCondition{Flux,Nothing}(nothing), BoundaryCondition{Flux,Nothing}(nothing))), v = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}(BoundaryCondition{Flux,Nothing}(nothing), BoundaryCondition{Flux,Nothing}(nothing))), w = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Oceananigans.NoPenetration,Nothing},BoundaryCondition{Oceananigans.NoPenetration,Nothing}}(BoundaryCondition{Oceananigans.NoPenetration,Nothing}(nothing), BoundaryCondition{Oceananigans.NoPenetration,Nothing}(nothing)))), pressure = (x = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), y = CoordinateBoundaryConditions{BoundaryCondition{Periodic,Nothing},BoundaryCondition{Periodic,Nothing}}(BoundaryCondition{Periodic,Nothing}(nothing), BoundaryCondition{Periodic,Nothing}(nothing)), z = CoordinateBoundaryConditions{BoundaryCondition{Flux,Nothing},BoundaryCondition{Flux,Nothing}}(BoundaryCondition{Flux,Nothing}(nothing), BoundaryCondition{Flux,Nothing}(nothing)))), Oceananigans.AdamsBashforthTimestepper{Float64,NamedTuple{(:u, :v, :w),Tuple{Field{Oceananigans.Face,Oceananigans.Cell,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Face,Oceananigans.Cell,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}},Field{Oceananigans.Cell,Oceananigans.Cell,Oceananigans.Face,OffsetArrays.OffsetArray{Float64,3,Array{Float64,3}},RegularCartesianGrid{Float64,StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}}}}}}}((u = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], v = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], w = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]), (u = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], v = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ], w = [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]

[0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0… ; ]), 0.125), Oceananigans.PoissonSolverCPU{Oceananigans.PPN,Array{Float64,3},Array{Complex{Float64},3},FFTW.cFFTWPlan{Complex{Float64},-1,true,3},FFTW.r2rFFTWPlan{Complex{Float64},(5,),true,3},AbstractFFTs.ScaledPlan{Complex{Float64},FFTW.cFFTWPlan{Complex{Float64},1,true,3},Float64},FFTW.r2rFFTWPlan{Complex{Float64},(4,),true,3}}(Oceananigans.PPN(), [0.0; 0.9997992185115969; … ; 3.996788270156917; 0.9997992185116003], [0.0 0.9997992185115969 … 3.996788270156917 0.9997992185116003], [0.0], Complex{Float64}[0.0 + 0.0im 0.0 + 0.0im … 0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.0 + 0.0im … 0.0 + 0.0im 0.0 + 0.0im; … ; 0.0 + 0.0im 0.0 + 0.0im … 0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.0 + 0.0im … 0.0 + 0.0im 0.0 + 0.0im], FFTW in-place forward plan for 128×128×1 array of Complex{Float64}
(dft-rank>=2/1
  (dft-direct-128-x128 "n1fv_128_avx2")
  (dft-direct-128-x128 "n1fv_128_avx2")), FFTW in-place r2r REDFT10 plan for 128×128×1 array of Complex{Float64}
(redft10e-r2hc-1-x32768
  (rdft-nop)), 6.103515625e-5 * FFTW in-place backward plan for 128×128×1 array of Complex{Float64}
(dft-rank>=2/1
  (dft-direct-128-x128 "n1bv_128_avx2")
  (dft-direct-128-x128 "n1bv_128_avx2")), FFTW in-place r2r REDFT01 plan for 128×128×1 array of Complex{Float64}
(rdft-nop)), nothing, OrderedCollections.OrderedDict{Symbol,Oceananigans.AbstractOutputWriter}(), OrderedCollections.OrderedDict{Symbol,Oceananigans.AbstractDiagnostic}(), nothing)

Our initial condition randomizes u and v. We also ensure that both have zero mean for purely aesthetic reasons.

u₀ = rand(size(model.grid)...)
u₀ .-= mean(u₀)

set!(model, u=u₀, v=u₀)

Next we define a function for calculating the vertical vorticity associated with the velocity fields u and v.

function vorticity!(ω, u, v)
    for j = 1:u.grid.Ny, i = 1:u.grid.Nx
        @inbounds ω.data[i, j, 1] = ∂x_faa(i, j, 1, u.grid, v.data) - ∂y_afa(i, j, 1, u.grid, u.data)
    end
    return nothing
end
vorticity! (generic function with 1 method)

Finally, we create the vorticity field for storing u and v, initialize a figure, and run the model forward

ω = Field(Face, Face, Cell, model.architecture, model.grid)

close("all")
fig, ax = subplots()

for i = 1:10
    time_step!(model, Nt=100, Δt=1e-1)

    vorticity!(ω, model.velocities.u, model.velocities.v)

    cla()
    imshow(data(ω)[:, :, 1])
    ax.axis("off")
    pause(0.1)
end

We can plot out the final vorticity field.

gcf()

This page was generated using Literate.jl.