Skip to content

Note: Output is not generated for this example (to save resources on GitHub).

Decaying Homogeneous Isotropic Turbulence - 3D

In this example we consider decaying homogeneous isotropic turbulence, similar to the cases considered in [1] and [2]. The initial velocity field is created randomly, but with a specific energy spectrum. Due to viscous dissipation, the turbulent features eventually group to form larger visible eddies.

julia
using CairoMakie
using IncompressibleNavierStokes

Output directory

julia
outdir = joinpath(@__DIR__, "output", "DecayingTurbulence3D")
ispath(outdir) || mkpath(outdir)

Floating point precision

julia
T = Float64

Array type

julia
ArrayType = Array
# using CUDA; ArrayType = CuArray
# using AMDGPU; ArrayType = ROCArray
# using oneAPI; ArrayType = oneArray
# using Metal; ArrayType = MtlArray

Reynolds number

julia
Re = T(6_000)

A 3D grid is a Cartesian product of three vectors

julia
n = 32
lims = T(0), T(1)
x = LinRange(lims..., n + 1), LinRange(lims..., n + 1), LinRange(lims..., n + 1)

Build setup and assemble operators

julia
setup = Setup(; x, Re, ArrayType);
nothing #hide

Since the grid is uniform and identical for x, y, and z, we may use a specialized spectral pressure solver

julia
psolver = psolver_spectral(setup);
nothing #hide

Initial conditions

julia
ustart = random_field(setup; psolver);
nothing #hide

Solve unsteady problem

julia
(; u, t), outputs = solve_unsteady(;
    setup,
    ustart,
    tlims = (T(0), T(1)),
    Δt = T(1e-3),
    psolver,
    processors = (
        # rtp = realtimeplotter(; setup, plot = fieldplot, nupdate = 10),
        ehist = realtimeplotter(;
            setup,
            plot = energy_history_plot,
            nupdate = 10,
            displayfig = false,
        ),
        espec = realtimeplotter(; setup, plot = energy_spectrum_plot, nupdate = 10),
        # anim = animator(; setup, path = "$outdir/solution.mkv", nupdate = 20),
        # vtk = vtk_writer(; setup, nupdate = 10, dir = outdir, filename = "solution"),
        # field = fieldsaver(; setup, nupdate = 10),
        log = timelogger(; nupdate = 100),
    ),
);
nothing #hide

Post-process

We may visualize or export the computed fields (u, p)

Energy history

julia
outputs.ehist

Energy spectrum

julia
outputs.espec

Export to VTK

julia
save_vtk(state; setup, filename = joinpath(outdir, "solution"), psolver)

This page was generated using Literate.jl.