Simulations

struct Simulation{M, E, S, R, RS, X, EX, DX, TS, FT, TR, C, H, O, U}

Simulation(mesh; eos,
           dust = nothing,
           boundary_conditions = BoundaryConditions(),
           reconstruction = PLM(),
           riemann_solver = HLLE(),
           stepper = RK2(FT),
           cfl = FT(3//10),
           sources = (),
           components = (),
           hooks = (;),
           outputs = (),
           units = nothing,
           gpu_aware_mpi = false)

A simulation's mesh, physics configuration, state arrays, and Clock. step! advances one cycle and run! advances to a stopping criterion.

Keyword arguments (FT is the mesh's float type; every scheme argument must carry it too):

• eos: the required equation of state, which selects hydrodynamics or constrained-transport MHD and the conserved-variable count.
• boundary_conditions: the six-face BoundaryConditions.
• reconstruction: DonorCell, PLM, PPM, PPM5, WENOZ, WENOZPP, WENOAO, or WENOAOPP; its stencil must fit the mesh's ghost width.
• riemann_solver: Rusanov, HLLE, HLLC / LHLLC (ideal hydrodynamics), or HLLD / LHLLD (MHD). Multidimensional LHLLD runs need two ghost cells. The default is the one solver defined for every equation of state.
• dust: a DustFluids system, or nothing. Coupled drag requires the IMEX2P stepper.
• stepper: RK2– RK4, or IMEX2P for stiff components.
• cfl: the Courant number as the fraction of the stepper's stability limit to run at, in (0, 1]. The crossing time underneath is a per-direction minimum, so multidimensional runs want roughly 1 / ndimensions or less; the default is safe in every dimensionality.
• sources: one or more cell-centered rate functions or SourceTerms.
• components: tuple of physics components extending the TimeSteppers protocol (TurbulenceDriving, EMFSource, user structs); within-slot execution order is the tuple order.
• hooks: NamedTuple with any subset of before_cycle / after_update / after_stage / after_cycle (see default_hooks); missing names stay empty.
• outputs: output writers (HistoryOutput, FieldOutput, RestartOutput), each with its own cadence.
• units: a CodeUnits naming what one code unit is in cgs, or nothing. Physics with dimensional inputs and the output metadata read from here.
• gpu_aware_mpi: hand device buffers straight to MPI (BoundaryExchange).

After construction the clock stands at t = 0 with dt = 0: set the state with set_initial_condition! — or fill state.w0 and call initialize! — before stepping.

• mesh::Any: domain geometry and block decomposition

• eos::Any: equation of state

• state::Any: state arrays: HydroState or MHDState by the EOS, inside a composite State when a subsystem is present

• reconstruction::Any: reconstruction scheme

• riemann_solver::Any: Riemann solver

• exchange::Any: the ghost exchange, holding the physical boundary conditions — for MHD the fused round carrying u0 and b0 together

• edge_exchange::Any: the mid-stage EMF edge round of sync_edges!; nothing for hydrodynamics

• dust_exchange::Any: ghost-cell exchange for the stacked dust species; nothing without dust

• stepper::Any: time integrator

• clock::Clock: simulation time, timestep, and cycle count

• cfl::Any: fraction of the stepper's stability limit to run at

• timestep_reduction::Any: in-flight slot of the final stage's overlapped timestep reduction (nothing when serial)

• cycle_seam::Bool: whether once-per-cycle injections defer the final-stage ghost exchange

• components::Any: physics components, executed in tuple order within each pipeline slot

• hooks::Any: user callables at the named cycle points

• outputs::Any: output writers, armed at construction and scheduled between cycles by run!

• units::Any: what one code unit is in cgs (CodeUnits), or nothing

using Aether

mesh = Mesh(CPU(); size = (64, 64, 1), extent = (1, 1, 1), cells_per_block = (32, 32, 1))
Simulation(mesh; eos = IdealHydro(1.4, 1e-12, 1e-10), cfl = 0.4)

# output
Simulation on CPU()
├── mesh: Mesh{Float64}(64 × 64 × 1 cells, 4 blocks)
├── eos: IdealHydro{Float64}(1.4, 1.0e-12, 1.0e-10, Inf)
├── state: HydroState{Float64}(36 × 36 × 1 cells, 5 variables, 4 blocks)
├── reconstruction: PLM()
├── riemann solver: HLLE()
├── boundary conditions: BoundaryConditions(PeriodicBC()/PeriodicBC(), PeriodicBC()/PeriodicBC(), PeriodicBC()/PeriodicBC())
├── stepper: SSPRK{Float64}(2 stages, cfl_limit = 1.0)
├── cfl: 0.4
├── components: (none)
├── hooks: (none)
├── outputs: (none)
└── clock: Clock{Float64}(t = 0.0, dt = 0.0, cycle = 0)

set_initial_condition!(
    problem,
    simulation;
    magnetic_field,
    vector_potential,
    dust
)

Set the initial condition from a point function of position and leave the simulation ready to step. problem is evaluated at every interior cell center as problem(x1, x2, x3) and returns the primitive variables — for the ideal equations of state (ρ, v1, v2, v3, e) with e the internal energy density (an ideal-gas pressure enters as e = p / (γ - 1)), for the isothermal ones (ρ, v1, v2, v3). Evaluation happens on the host, block by block, followed by one transfer; for kernel-filled or file-based setups fill simulation.state.w0 (and for MHD the faces of simulation.state.b0) directly and finish with initialize!, which this function ends with.

An MHD simulation takes its face-centered field from one of two keywords:

• vector_potential = A: A(x1, x2, x3) returns (A1, A2, A3); the field is the discrete curl of A sampled at edge centers, so $\nabla \cdot B$ vanishes to round-off identically, whatever A is. Preferred.
• magnetic_field = B: B(x1, x2, x3) returns (B1, B2, B3), sampled directly at face centers. The discrete divergence is then whatever the samples imply — exact for fields with each component constant along its own direction (uniform fields, Brio & Wu, Orszag–Tang), O(Δx^2) otherwise.

Positions along singleton directions evaluate at the domain center, and both faces of a singleton layer receive the same value. Omitting both keywords starts from B = 0; giving either with a hydrodynamic EOS throws.

A dusty simulation takes one point function per species through the dust keyword: dust = (f₁, …, f_N), each evaluated at the interior cell centers as f(x1, x2, x3) -> (ρ, v1, v2, v3). Omitting it leaves the dust at the zero-initialized state (floored to the dust density floor on the first refresh).

using Aether

mesh = Mesh(CPU(); size = (64, 1, 1), extent = (1, 1, 1))
simulation = Simulation(mesh; eos = IdealHydro(1.4, 1e-12, 1e-10),
                        boundary_conditions = BoundaryConditions(ix1 = OutflowBC(),
                                                                 ox1 = OutflowBC()),
                        cfl = 0.8)
set_initial_condition!(simulation) do x1, x2, x3
    ρ, p = x1 < 1/2 ? (1.0, 1.0) : (0.125, 0.1)
    return (ρ, 0, 0, 0, p / (1.4 - 1))
end
simulation.clock

# output
Clock{Float64}(t = 0.0, dt = 0.010564428184106458, cycle = 0)

initialize!(simulation)

Take a simulation whose interior primitive state w0 (and, for MHD, interior face field b0) is set and make it ready to step: conserved variables over the interior, one refresh_primitives! — ghost exchanges, bcc for MHD, primitives over the padded box with the floors applied — and the first clock.dt. set_initial_condition! ends here; call it directly after filling the state yourself.

run!(simulation; stop_time, stop_cycle, progress)

March the simulation with step! until clock.t reaches stop_time or clock.cycle reaches stop_cycle — at least one must be given. The last step is shortened to land on stop_time exactly, so calling run! again with a later stop_time continues the run in exact segments.

The simulation's outputs fire between cycles (write_outputs!): anything due at entry — a fresh run dumps its initial condition — then after every step, and a forced final write so the last state is never lost when stop_time is not a multiple of an output's dt. Outputs never touch the timestep: the trajectory is bitwise independent of the output configuration.

progress selects the run's progress display, handled entirely on rank 0 with no communication — the trajectory is also independent of it:

• :auto (default): a LiveProgress in an interactive session (REPL or Jupyter), a default LogProgress when a script runs on a terminal or under MPI on several ranks (mpiexec pipes the streams, but a distributed run still earns its log), and silent when the streams are redirected on a single rank — doctests and captured logs stay clean.
• a LiveProgress: an in-place animated status line with progress bar, cycle, t, dt, zones per second, and the estimated time remaining.
• a LogProgress: an appended log line every every_cycles cycles with walltime, cycle, t, dt, and zones per second — for batch runs.
• nothing: silent.

using Aether

mesh = Mesh(CPU(); size = (64, 1, 1), extent = (1, 1, 1))
simulation = Simulation(mesh; eos = IdealHydro(1.4, 1e-12, 1e-10),
                        boundary_conditions = BoundaryConditions(ix1 = OutflowBC(),
                                                                 ox1 = OutflowBC()),
                        cfl = 0.8)
set_initial_condition!(simulation) do x1, x2, x3
    ρ, p = x1 < 1/2 ? (1.0, 1.0) : (0.125, 0.1)
    return (ρ, 0, 0, 0, p / (1.4 - 1))
end
run!(simulation; stop_time = 0.2)
simulation.clock

# output
Clock{Float64}(t = 0.2, dt = 0.005655450061433209, cycle = 34)

struct LiveProgress

LiveProgress(; interval = 0.1, bar_width = 16, io = stderr)

In-place progress display for interactive sessions (REPL, Jupyter): one line, rewritten in place at most every interval seconds of walltime, showing an animated progress bar with the current cycle, simulation time t, timestep dt, the smoothed update rate in zones per second, and the estimated time to completion. The bar fills a starfield: the unsimulated region twinkles, a shooting star streaks across it every few seconds, and rare sights — a long-tailed comet, a supernova, a passing satellite — reward the long waits. The sky quickens for a few seconds when dt drops sharply, so the display reacts to the physics it is marching; a shimmer drifts through the simulated region, and output writes pin a brief note to the line. A one-time opening header names the mesh, architecture, rank count, and the run's span. A wordmark wave whose speed follows the update rate rolls in the header on a terminal, and leads the line itself in Jupyter.

Pass as run!'s progress keyword; it is also what progress = :auto (the default) selects when an interactive display is available. Only rank 0 displays, and no communication is involved: progress never perturbs the run.

using Aether

LiveProgress(interval = 0.5)

# output
LiveProgress(interval = 0.5 s)

struct LogProgress

LogProgress(; every_cycles = 100, io = stdout)

Appending progress log for batch runs (julia script.jl, job schedulers): one line every every_cycles cycles showing the walltime since the application opened (since Aether was loaded), the current cycle, simulation time t, timestep dt, and the update rate in zones per second averaged since the previous line. An opening header line names the mesh, architecture, rank count, and the run's span; a final line marked run complete reports the whole run's average rate.

Pass as run!'s progress keyword. Only rank 0 logs, and no communication is involved: progress never perturbs the run.

using Aether

LogProgress(every_cycles = 200)

# output
LogProgress(every_cycles = 200)