Solving unsteady problems

The main entry point of the package is solve_unsteady:

state, outputs = solve_unsteady(;
    setup,
    start = (; u),
    tlims = (0.0, 1.0),
    params = (; viscosity = 1e-3),
)

It steps the initial start state (a named tuple of fields, e.g. (; u) or (; u, temp)) through the time interval tlims and returns the final state and the processor outputs.

Forces and parameters

The right-hand side is a function force!(force, state, t; setup, cache, params...) passed as the force! keyword. The built-in forces are navierstokes! (convection and diffusion, the default) and boussinesq! (adds a coupled temperature equation). Physical parameters are passed in the params named tuple and forwarded to force! as keyword arguments, e.g. params = (; viscosity = 1e-3) for navierstokes!. Custom forces follow the same signature; see examples/Kolmogorov2D.jl for a custom body force and examples/ChannelFlow.jl for adding an eddy-viscosity closure model.

Time step and method

By default, the time step is chosen adaptively from a CFL condition on the convective and diffusive limits (tune with the cfl keyword); pass a fixed Δt to disable this. The time integration method defaults to the low-storage third-order Runge-Kutta method LMWray3; any explicit Runge-Kutta tableau from RKMethods can be passed via the method keyword. The theory is discussed in Spatial and temporal discretization.

For finer control than solve_unsteady, the stepping primitives create_stepper and timestep/timestep! can be called directly, e.g. to write a custom (differentiable) simulation loop.

API

IncompressibleNavierStokes.boussinesq! Method

boussinesq!(
    force,
    state,
    t;
    setup,
    cache,
    viscosity,
    conductivity,
    gdir,
    gravity,
    dodissipation
)

Boussinesq forcing (Navier-Stokes + gravity for u, convection-diffusion for temp).

source

IncompressibleNavierStokes.boussinesq Method

boussinesq(
    state,
    t;
    setup,
    viscosity,
    conductivity,
    gdir,
    gravity,
    dodissipation
)

Boussinesq forcing (Navier-Stokes + gravity for u, convection-diffusion for temp).

source

IncompressibleNavierStokes.get_state Method

get_state(
    stepper
) -> Union{NamedTuple, DataStructures.SortedMultiDict}

Get state (; u, temp, t, n) from stepper.

source

IncompressibleNavierStokes.navierstokes! Method

navierstokes!(force, state, t; setup, cache, viscosity)

Navier-Stokes momentum forcing (convection + diffusion).

source

IncompressibleNavierStokes.navierstokes Method

navierstokes(state, t; setup, viscosity)

Navier-Stokes momentum forcing (convection + diffusion).

source

IncompressibleNavierStokes.solve_unsteady Method

solve_unsteady(
;
    setup,
    tlims,
    start,
    force!,
    docopy,
    method,
    psolver,
    Δt,
    Δt_min,
    cfl,
    n_adapt_Δt,
    processors,
    params,
    ode_cache,
    force_cache
)

Solve unsteady problem using method.

The initial start state is a named tuple of fields, e.g. (; u) or (; u, temp). The right-hand side force! is called as force!(force, state, t; setup, cache, params...), where params is a named tuple of parameters passed as keyword arguments (e.g. params = (; viscosity) for the default navierstokes!, or (; viscosity, conductivity, gdir, gravity, dodissipation) for boussinesq!).

If Δt is a real number, it is rounded such that (t_end - t_start) / Δt is an integer. If Δt = nothing, the time step is chosen every n_adapt_Δt iteration with CFL-number cfl. If Δt_min is given, the adaptive time step never goes below it.

The processors are called after every time step.

Note that the state observable passed to the processor.initialize function contains fields living on the device, and you may have to move them back to the host using Array in the processor.

Return (; state..., t), outputs, where outputs is a named tuple with the outputs of processors with the same field names.

source

Time stepping

IncompressibleNavierStokes.AbstractODEMethod Type

abstract type AbstractODEMethod{T}

Abstract ODE method.

Fields

source

IncompressibleNavierStokes.AbstractRungeKuttaMethod Type

abstract type AbstractRungeKuttaMethod{T} <: IncompressibleNavierStokes.AbstractODEMethod{T}

Abstract Runge Kutta method.

Fields

source

IncompressibleNavierStokes.ExplicitRungeKuttaMethod Type

struct ExplicitRungeKuttaMethod{T} <: IncompressibleNavierStokes.AbstractRungeKuttaMethod{T}

Explicit Runge Kutta method. See Sanderse [6].

Consider the velocity field $u_0$ at a certain time $t_0$. We will now perform one time step to $t=t_0+\mathrm{\Delta }t$. For explicit Runge-Kutta methods, this time step is divided into $s$ sub-steps $t_i=t_0+\mathrm{\Delta }{t}_i$ with increment $\mathrm{\Delta }{t}_i=c_i\mathrm{\Delta }t$. The final substep performs the full time step $\mathrm{\Delta }{t}_s=\mathrm{\Delta }t$ such that $t_s=t$.

For $i=1,\dots ,s$, the intermediate velocity $u_i$ and pressure $p_i$ are computed as follows:

$$\begin{array}{rl}{k}_i& =F(u_{i-1},t_{i-1})-y_G(t_{i-1})\\ v_i& =u_0+\mathrm{\Delta }t\sum _{j=1}^i{a}_{ij}{k}_j\\ L{p}_i& =WM\frac{1}{c_i}\sum _{j=1}^i{a}_{ij}{k}_j+W\frac{y_M(t_i)-y_M(t_0)}{\mathrm{\Delta }{t}_i}\\ & =W\frac{(M{v}_i+y_M(t_i))-(M{u}_0+y_M(t_0))}{\mathrm{\Delta }{t}_i^n}\\ & =W\frac{M{v}_i+y_M(t_i)}{\mathrm{\Delta }{t}_i^n}\\ u_i& =v_i-\mathrm{\Delta }{t}_{i}G{p}_i,\end{array}$$

where $(a_{ij}{)}_{ij}$ are the Butcher tableau coefficients of the RK-method, with the convention $c_i=\sum _{j=1}^i{a}_{ij}$.

Finally, we return $u_s$. If $u_0=u(t_0)$, we get the accuracy $u_s=u(t)+\mathcal{O}(\mathrm{\Delta }{t}^{r+1})$, where $r$ is the order of the RK-method. If we perform $n$ RK time steps instead of one, starting at exact initial conditions $u^0=u(0)$, then $u^n=u(t^n)+\mathcal{O}(\mathrm{\Delta }{t}^r)$ for all $n\in \{1,\dots ,N\}$. Note that for a given $u$, the corresponding pressure $p$ can be calculated to the same accuracy as $u$ by doing an additional pressure projection after each outer time step $\mathrm{\Delta }t$ (if we know $\frac{\mathrm{d}{y}_M}{\mathrm{d}t}(t)$), or to first order accuracy by simply returning $p_s$.

Note that each of the sub-step velocities $u_i$ is divergence free, after projecting the tentative velocities $v_i$. This is ensured due to the judiciously chosen replacement of $\frac{\mathrm{d}{y}_M}{\mathrm{d}t}(t_i)$ with $(y_M(t_i)-y_M(t_0))/\mathrm{\Delta }{t}_i$. The space-discrete divergence-freeness is thus perfectly preserved, even though the time discretization introduces other errors.

Fields

• A

• b

• c

• r

• p_add_solve

source

IncompressibleNavierStokes.LMWray3 Type

struct LMWray3{T} <: IncompressibleNavierStokes.AbstractRungeKuttaMethod{T}

Low memory Wray 3rd order scheme [8]. Uses 3 vector fields and one scalar field.

Fields

source

IncompressibleNavierStokes.runge_kutta_method Method

runge_kutta_method(
    A,
    b,
    c,
    r;
    T,
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod

Get Runge Kutta method. The function checks whether the method is explicit.

p_add_solve: whether to add a pressure solve step to the method.

For implicit RK methods: newton_type, maxiter, abstol, reltol.

source

IncompressibleNavierStokes.create_stepper Function

create_stepper(method; setup, psolver, state, t, n = 0)

Create time stepper.

source

IncompressibleNavierStokes.timestep Function

timestep(method, force, stepper, Δt; params = nothing)

Perform one time step.

Non-mutating/allocating/out-of-place version.

See also timestep!.

source

IncompressibleNavierStokes.timestep! Function

timestep!(method, force!, stepper, Δt; params = nothing, ode_cache, force_cache)

Perform one time step.

Mutating/non-allocating/in-place version.

See also timestep.

source

IncompressibleNavierStokes.get_cache Function

get_cache(method, state, setup)

Get time stepper cache for the given ODE method.

The method get_cache(force!, setup) returns the cache for a right-hand-side function force! instead.

source

Runge-Kutta methods

IncompressibleNavierStokes.RKMethods Module

Set up Butcher arrays A, b, and c, as well as and SSP coefficient r. For families of methods, optional input s is the number of stages.

Original (MATLAB) by David Ketcheson, extended by Benjamin Sanderse.

Exports

• DOPRI6

• FE11

• HEM3

• HEM3BS

• HEM5

• Heun33

• MTE22

• Mid22

• NSSP21

• NSSP32

• NSSP33

• NSSP53

• RK33C2

• RK33P2

• RK44

• RK44C2

• RK44C23

• RK44P2

• RK56

• SSP104

• SSP22

• SSP33

• SSP42

• SSP43

• Wray3

• rSSPs2

• rSSPs3

source

IncompressibleNavierStokes.RKMethods.DOPRI6 Method

DOPRI6(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Dormand-Prince pair.

source

IncompressibleNavierStokes.RKMethods.FE11 Method

FE11(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

FE11 (Forward Euler).

source

IncompressibleNavierStokes.RKMethods.HEM3 Method

HEM3(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Brasey and Hairer.

source

IncompressibleNavierStokes.RKMethods.HEM3BS Method

HEM3BS(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

HEM3BS.

source

IncompressibleNavierStokes.RKMethods.HEM5 Method

HEM5(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Brasey and Hairer, 5 stage, 4th order.

source

IncompressibleNavierStokes.RKMethods.Heun33 Method

Heun33(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Heun33.

source

IncompressibleNavierStokes.RKMethods.MTE22 Method

MTE22(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Minimal truncation error 22 method (Heun).

source

IncompressibleNavierStokes.RKMethods.Mid22 Method

Mid22(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Midpoint 22 method.

source

IncompressibleNavierStokes.RKMethods.NSSP21 Method

NSSP21(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

NSSP21.

source

IncompressibleNavierStokes.RKMethods.NSSP32 Method

NSSP32(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

NSSP32.

source

IncompressibleNavierStokes.RKMethods.NSSP33 Method

NSSP33(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

NSSP33.

source

IncompressibleNavierStokes.RKMethods.NSSP53 Method

NSSP53(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

NSSP53.

source

IncompressibleNavierStokes.RKMethods.RK33C2 Method

RK33C2(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK3 satisfying C(2) for i=3.

source

IncompressibleNavierStokes.RKMethods.RK33P2 Method

RK33P2(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK3 satisfying the second order condition for the pressure.

source

IncompressibleNavierStokes.RKMethods.RK44 Method

RK44(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Classical fourth order.

source

IncompressibleNavierStokes.RKMethods.RK44C2 Method

RK44C2(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK4 satisfying C(2) for i=3.

source

IncompressibleNavierStokes.RKMethods.RK44C23 Method

RK44C23(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK4 satisfying C(2) for i=3 and c2=c3.

source

IncompressibleNavierStokes.RKMethods.RK44P2 Method

RK44P2(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK4 satisfying the second order condition for the pressure (but not third order).

source

IncompressibleNavierStokes.RKMethods.RK56 Method

RK56(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

RK56.

source

IncompressibleNavierStokes.RKMethods.SSP104 Method

SSP104(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

SSP104.

source

IncompressibleNavierStokes.RKMethods.SSP22 Method

SSP22(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

SSP22.

source

IncompressibleNavierStokes.RKMethods.SSP33 Method

SSP33(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

SSP33.

source

IncompressibleNavierStokes.RKMethods.SSP42 Method

SSP42(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

SSP42.

source

IncompressibleNavierStokes.RKMethods.SSP43 Method

SSP43(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

SSP43.

source

IncompressibleNavierStokes.RKMethods.Wray3 Method

Wray3(
;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod{Float64}

Wray's RK3.

source

IncompressibleNavierStokes.RKMethods.rSSPs2 Function

rSSPs2(
;
    ...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod
rSSPs2(
    s;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod

Rational (optimal, low-storage) s-stage 2nd order SSP.

source

IncompressibleNavierStokes.RKMethods.rSSPs3 Function

rSSPs3(
;
    ...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod
rSSPs3(
    s;
    kwargs...
) -> IncompressibleNavierStokes.ExplicitRungeKuttaMethod

Rational (optimal, low-storage) s^2-stage 3rd order SSP.

source