Integrator Examples

The snippets below are written to match current oakcli usage. Where a block creates its own ctx, you can run the same pattern under oakcli run --script your_script.lua.

Fixed-Step RK4

The workhorse for smooth, well-conditioned fields. Four drift evaluations per step give 4th-order accuracy without any adaptive overhead:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.1,
    wavenumber = 0.35,
    omega = 0.2
})

local integrator = ooc.create_context_integrator(ctx, "rk4", {
    initial_dt = 0.01
})
ooc.set_integrator(ctx, integrator)

for i = 1, 1000 do
    ooc.step(ctx)
end

---

Adaptive RKF45

The Dormand–Prince embedded pair estimates the local truncation error and automatically shrinks or grows the timestep to stay within tolerance. Ideal for problems with intermittent sharp features:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.02)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.15,
    wavenumber = 0.6,
    omega = 0.3
})

local integrator = ooc.create_context_integrator(ctx, "rkf45", {
    initial_dt = 0.02,
    min_dt     = 1.0e-5,
    max_dt     = 0.1,
    tolerance  = 1.0e-4,
    safety     = 0.9,
    adaptive   = true
})
ooc.set_integrator(ctx, integrator)

for i = 1, 500 do
    ooc.step(ctx)
    local m = ooc.step_metrics_latest(ctx)
    ooc.log("step=%d dt=%.5f err=%.2e", m.step_index, m.accepted_dt, m.rms_error)
end

---

ETDRK4 for Supported Semilinear Complex Flows

Use ETDRK4 when the target field is complex and the plan decomposes into a supported exact-linear part plus same-field nonlinear or dt-scaled source terms. This is a strong fit for dispersive or damped spectral flows with additive forcing:

local ctx = ooc.create()
local dt = 0.0125

ooc.set_timestep(ctx, dt)

local field = ooc.add_field(ctx, {1024}, {
    type = "complex_double",
    fill = {0.0, 0.0}
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_chirp",
    amplitude = 0.16,
    wavenumber = 0.018,
    omega = -0.18,
    phase = math.pi / 5,
    kdot = -0.055,
    wdot = 0.014,
    rotation = 0.45,
    scale_by_dt = true
})

ooc.add_dispersion_operator(ctx, field, {
    coefficient = 0.78,
    order = 1.08,
    spacing = 0.02
})

ooc.add_linear_dissipative_operator(ctx, field, {
    coefficient = 0.018,
    order = 2.0,
    spacing = 0.02
})

local integrator = ooc.create_context_integrator(ctx, "etdrk4", {
    initial_dt = dt,
    target_field = field
})

ooc.set_integrator(ctx, integrator)

for i = 1, 200 do
    ooc.step(ctx)
end

Keep the plan on that single target field. For ETDRK4-supported stimulus operators, set scale_by_dt = true; time-dependent stimulus families also need their normal evolving clock behavior rather than fixed_clock = true.

---

Implicit Crank–Nicolson

A-stable trapezoidal rule — safe on stiff diffusion-dominated equations without sacrificing 2nd-order accuracy:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.05)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 1.0
})

ooc.add_linear_dissipative_operator(ctx, field, {
    coefficient = 0.04,
    order = 2.0,
    spacing = 0.02
})

local integrator = ooc.create_context_integrator(ctx, "crank_nicolson", {
    initial_dt = 0.05
})
ooc.set_integrator(ctx, integrator)

for i = 1, 200 do
    ooc.step(ctx)
end

---

Backward Euler — Maximum Stiffness Tolerance

L-stable implicit method; first-order but unconditionally stable on arbitrarily stiff linear operators. Use when Crank–Nicolson introduces oscillations near sharp transients:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.1)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 1.0
})

ooc.add_linear_dissipative_operator(ctx, field, {
    coefficient = 0.08,
    order = 2.0,
    spacing = 0.02
})

local integrator = ooc.create_context_integrator(ctx, "backward_euler", {
    initial_dt = 0.1
})
ooc.set_integrator(ctx, integrator)

for i = 1, 200 do
    ooc.step(ctx)
end

---

Subordination with Stochastic Noise

Two-stage subordinated integrator with Gaussian noise injection — suitable for fractional diffusion or stochastic PDEs:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.08,
    wavenumber = 0.4,
    omega = 0.15
})

local integrator = ooc.create_context_integrator(ctx, "subordination", {
    initial_dt          = 0.01,
    enable_stochastic   = true,
    stochastic_strength = 0.2,
    noise               = "gaussian",
    random_seed         = 1337
})
ooc.set_integrator(ctx, integrator)

for i = 1, 100 do
    ooc.step(ctx)
end

Switch to Laplace noise for heavier-tailed fluctuations:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.08,
    wavenumber = 0.4,
    omega = 0.15
})

local integrator = ooc.create_context_integrator(ctx, "subordination", {
    initial_dt          = 0.01,
    enable_stochastic   = true,
    stochastic_strength = 0.1,
    noise               = "laplace",
})
ooc.set_integrator(ctx, integrator)

for i = 1, 100 do
    ooc.step(ctx)
end

---

Multi-Integrator Sequence

set_integrator_sequence lets you assign different integration methods to different fields and step them together each tick. Each integrator advances at the same requested dt; the simulation clock advances by the minimum accepted step and reports the maximum error across all integrators.

A typical use: explicit RK4 for a reaction field (fast but smooth) and Crank–Nicolson for a diffusion field (stiff, needs implicit treatment):

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

-- Field 0: reaction term (explicit OK)
local reaction_field  = ooc.add_field(ctx, {256}, { fill = 0.0 })
-- Field 1: diffusion concentration (stiff)
local diffusion_field = ooc.add_field(ctx, {256}, { fill = 0.5 })

ooc.add_stimulus_operator(ctx, reaction_field, {
    type = "stimulus_sine",
    amplitude = 0.05,
    wavenumber = 0.08,
    omega = 0.1
})

ooc.add_linear_dissipative_operator(ctx, diffusion_field, {
    coefficient = 0.03,
    order = 2.0,
    spacing = 0.02
})

local rk4 = ooc.create_context_integrator(ctx, "rk4", {
    initial_dt   = 0.01,
    target_field = reaction_field
})

local cn = ooc.create_context_integrator(ctx, "crank_nicolson", {
    initial_dt   = 0.01,
    target_field = diffusion_field
})

ooc.set_integrator_sequence(ctx, { rk4, cn })

for i = 1, 500 do
    ooc.step(ctx)
end

Clear the sequence (detach all integrators):

ooc.set_integrator_sequence(ctx, {})  -- empty table clears the sequence
-- or equivalently:
ooc.detach_integrator(ctx)

---

Manual Step Loop with Integrator Step

Drive an integrator manually for custom multi-rate schemes or when you want explicit control over the requested dt:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.1,
    wavenumber = 0.5,
    omega = 0.2
})

local integrator = ooc.create_context_integrator(ctx, "rkf45", {
    initial_dt = 0.01,
    adaptive   = true,
    min_dt     = 1e-5,
    max_dt     = 0.1,
    tolerance  = 1e-4,
})

local dt = 0.01
for i = 1, 200 do
    local ok, err = ooc.integrator_step(ctx, integrator, dt)
    if not ok then
        ooc.log("step %d failed: error code %d", i, err)
        break
    end
    local m = ooc.step_metrics_latest(ctx)
    dt = m.next_dt  -- follow the adaptive suggestion for the next iteration
end

Use auto_advance = false only for probe steps where you intentionally do not want public time or step metrics advanced:

local ok, err = ooc.integrator_step(ctx, integrator, 0.01, false)
if not ok then
    ooc.log("probe step failed: code %d", err)
end

There is no later zero-dt “commit” call for that probe result. If you want the simulation clock and metrics history updated, call integrator_step with auto_advance = true.

---

Reading Integrator Info

Inspect the live configuration of an attached integrator:

local ctx = ooc.create()
local field = ooc.add_field(ctx, {128}, { fill = 0.0 })
ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.1,
    wavenumber = 0.25,
    omega = 0.1
})

local integrator = ooc.create_context_integrator(ctx, "rk4", {
    initial_dt = 0.01
})
ooc.set_integrator(ctx, integrator)

local info = ooc.integrator_info(ctx)
if info then
    ooc.log("method=%s  adaptive=%s  current_dt=%.4f  tolerance=%.2e",
        info.name,
        tostring(info.adaptive_enabled),
        info.current_dt,
        info.tolerance)
end

---

Tuning Parameters at Runtime

Switch adaptive mode and tighten tolerance mid-run without recreating the integrator:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.02)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.12,
    wavenumber = 0.5,
    omega = 0.25
})

local integrator = ooc.create_context_integrator(ctx, "rkf45", {
    initial_dt = 0.02,
    adaptive   = true,
    tolerance  = 1.0e-3
})
ooc.set_integrator(ctx, integrator)

-- First 100 steps: loose tolerance for fast spin-up
ooc.set_integrator_tolerance(ctx, 1.0e-2)
for i = 1, 100 do ooc.step(ctx) end

-- Tighten for the production run
ooc.set_integrator_tolerance(ctx, 1.0e-5)
ooc.set_integrator_safety(ctx, 0.85)
for i = 1, 1000 do ooc.step(ctx) end

Enable noise partway through a run:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.01)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.08,
    wavenumber = 0.35,
    omega = 0.15
})

local integrator = ooc.create_context_integrator(ctx, "subordination", {
    initial_dt = 0.01
})
ooc.set_integrator(ctx, integrator)

for i = 1, 200 do ooc.step(ctx) end

ooc.set_integrator_stochastic(ctx, true)
ooc.set_integrator_stochastic_strength(ctx, 0.05)
ooc.set_integrator_noise(ctx, "gaussian")

for i = 1, 800 do ooc.step(ctx) end

---

Tracking Step Metrics

Log the most recent step’s timing and error:

local ctx = ooc.create()
ooc.set_timestep(ctx, 0.02)

local field = ooc.add_field(ctx, {256}, {
    type = "double",
    fill = 0.0
})

ooc.add_stimulus_operator(ctx, field, {
    type = "stimulus_sine",
    amplitude = 0.1,
    wavenumber = 0.45,
    omega = 0.2
})

local integrator = ooc.create_context_integrator(ctx, "rkf45", {
    initial_dt = 0.02,
    adaptive   = true,
    tolerance  = 1.0e-4
})
ooc.set_integrator(ctx, integrator)

for i = 1, 500 do
    ooc.step(ctx)
end

local m = ooc.step_metrics_latest(ctx)
ooc.log("last step: index=%d accepted_dt=%.4f rms_err=%.2e wall=%.2fms",
    m.step_index, m.accepted_dt, m.rms_error, m.step_wall_ms)

Compute average wall time from the step history:

local history = ooc.step_metrics_history(ctx)
local total_ms = 0.0
local min_dt   = math.huge
local max_dt   = 0.0

for _, m in ipairs(history) do
    total_ms = total_ms + m.step_wall_ms
    if m.accepted_dt < min_dt then min_dt = m.accepted_dt end
    if m.accepted_dt > max_dt then max_dt = m.accepted_dt end
end

ooc.log("samples=%d  avg_wall=%.3fms  dt_range=[%.5f, %.5f]",
    #history, total_ms / #history, min_dt, max_dt)