API

A stopping rule evaluated on the fine trajectories from the current starts.

AbstractPararealBackend

How the fine sweep of Parareal executes. The algorithm lives ONCE in parareal/parareal.jl; a backend supplies only these primitives:

• is_root(backend) — whether this process runs the serial parts (coarse init, correction, convergence check). true everywhere except non-root MPI ranks.
• sync_starts!(backend, cache, k) — make the chunk starting values U[n] available wherever chunk n will be solved. No-op off MPI.
• fine_map!(backend, cache, solution, problem, parareal, k; …) — solve chunks k:N with the fine solver; deposit boundary values into cache.F[n+1] on the root and chunk solutions into solution.lastiterate (MPI ranks keep their chunk local until collect_chunks!).
• snapshot!(backend, solution, k; …) — record iterate k when saveiterates is on. In memory by default; on MPI the ranks have already written their chunk to the run directory inside fine_map!.
• sync_converged(backend, flag) — agree the convergence decision.
• collect_chunks!(backend, solution, cache; …) — bring final chunks (and, under MPI with saveiterates, the per-iterate history) to the root.

An abstract type for parameters of Parareals.

An abstract type for caching intermediate computations in AbstractTimeParallelSolvers.

An abstract type for iterates in AbstractTimeParallelSolutions.

An abstract type for parameters in AbstractTimeParallelSolvers.

An abstract type for time-parallel solutions of NSDEBase.AbstractInitialValueProblems.

An abstract type for time-parallel solvers of NSDEBase.AbstractInitialValueProblems.

Immutable weight-rule settings; bases are per unit time.

Fine chunks via pmap over the default worker pool (Distributed.jl).

Use a fixed positive per-unit-time discount base.

Pipeline-only relative frontier test. It is not the bulk all-interface update test.

SPMD over MPI ranks: rank r owns chunk r + 1. Requires MPI.Init().

Legacy mean relative interface-update criterion.

Legacy neighbouring-chunk running maximum; prefer an explicit fixed or secant rule.

Disable tolerance acceptance while retaining budgets, failure checks, and exact completion.

Immutable Parareal settings. Per-run propagation plans and weights live in the cache.

PararealCache <: AbstractTimeParallelCache

Pre-allocated state for one Parareal run: the chunk boundary values U, the fine and coarse boundary results F and G, the previous-iterate snapshot U_ for the error function, the chunk time grid T, the chunk problems (built once, their u0 slots referencing U so the correction updates them in place), and reusable solver caches. Method-dependent fields can be replaced when switching solver, RHS, or weight-rule types; the state shape and chunk count stay fixed. Each run refreshes the initial state from the supplied problem. Moving a cache to a different span requires shiftwindow!. RHS changes rebuild the chunk caches; do not mutate RHS internals in place while retaining their solver caches.

Run accounting. Bulk-synchronous counts are unavailable (nothing) for a pipeline.

PararealIterate <: AbstractTimeParallelIterate

A composite type for a single iterate in a PararealSolution: one fine chunk solution per time chunk. The chunk storage is CONCRETE — chunks are built eagerly from the fine solver, so eltype(iterate.chunks) is the fine solver's concrete solution type, not an abstract box. This matters: these objects cross thread and process boundaries in the hot path.

Constructors

PararealIterate(chunks::AbstractVector{𝕊}) where 𝕊<:AbstractInitialValueSolution
PararealIterate(problem::AbstractInitialValueProblem, parareal::Parareal)

Functions

• firstindex/getindex/lastindex/setindex! : chunk access
• length/numchunks : number of chunks

Methods

(iterate::PararealIterate)(t::Real)

returns the value of iterate at t via interpolation of the owning chunk.

PararealParameters <: AbstractPararealParameters

A composite type for the basic parameters of Parareal.

Constructors

PararealParameters(N, K)
PararealParameters(; N=10, K=N)

Arguments

• N :: Integer : number of time chunks/processors.
• K :: Integer : maximum number of fine batches. A bulk-synchronous return uses at most K-1 correction sweeps because its fine trajectories are evaluated from the returned starts.

PararealSolution <: AbstractTimeParallelSolution

A composite type for an AbstractTimeParallelSolution obtained using Parareal.

Constructors

PararealSolution(lastiterate, errors)
PararealSolution(problem::AbstractInitialValueProblem, parareal::Parareal)

Arguments

• lastiterate :: PararealIterate
• errors :: AbstractVector{ℝ} where ℝ<:Real : iteration errors.

Functions

• firstindex : first index.
• getindex : get iterate.
• lastindex : last index.
• numiterates : number of iterates.
• numchunks : number of chunks of last iterate.
• setindex! : set iterate.

Methods

(solution::PararealSolution)(t::Real)

returns the value of solution at t via interpolation.

PipelinedMPIBackend()

MPI frontier scheduling with one rank per chunk. Supported criteria are WeightedDefect(p=Inf, weight=FixedWeight(...)), FrontierRelativeUpdate, and NoStop. Other criteria are rejected; the bulk all-interface update rule and changing secant weights require synchronized iterates.

Maximum-defect acceptance waits for the upstream start to freeze before accepting its outgoing defect. Thus the assembled accepted starts satisfy the selected defect tolerance. Converting that tolerance to a trajectory error additionally requires bounds on the fine-map amplification.

Call MPI.Init() first and launch at least N ranks. Global histories are unsupported. Diagnostics report actual total fine/coarse calls and per-chunk local_sweeps; sweeps and fine_batches are nothing because this schedule has no synchronized batches. residuals contains the maximum diagnostic among ranks active at each local iteration index, not a bulk residual.

The final gather returns each fine trajectory from its accepted start. With NoStop() and K >= N, completion reproduces reference(problem, solver).

All-interface mixed update test, with explicit absolute and relative tolerances.

Latest-pair same-chunk secant rule with a scaled denominator guard and a floor of one.

Runs everything in the calling task. The reference implementation.

Fine chunks via Threads.@threads with dynamic scheduling in shared memory.

Normalized weighted defect p-mean; p=1 is the average and p=Inf the maximum.

firstindex(iterate::PararealIterate)

returns the first index of iterate.

firstindex(solution::PararealSolution)

returns the first index of solution.

getindex(iterate::PararealIterate, n::Integer)

returns the n-th chunk of iterate.

getindex(solution::PararealSolution, n::Integer)

returns the n-th chunk of the last iteration of a PararealSolution.

lastindex(iterate::PararealIterate)

returns the last index of iterate.

lastindex(solution::PararealSolution)

returns the last index of solution.

length(iterate::PararealIterate)

returns the number of chunks of iterate.

length(solution::PararealSolution)

returns the number of chunks of solution.

setindex!(iterate::PararealIterate, value::AbstractInitialValueSolution, n::Integer)

stores value into the n-th chunk of iterate.

setindex!(solution::PararealSolution, chunk::AbstractInitialValueSolution, n::Integer)

stores an AbstractInitialValueSolution as the n-th chunk of the last iteration of a PararealSolution.

NSDEBase.initialize_cache(problem, parareal::Parareal) :: PararealCache

returns a reusable PararealCache for problem, honouring the same contract NSDERungeKutta honours for its solvers. Build once, then call parareal(cache, solution, problem; ...) repeatedly — e.g. in timing loops, to report setup costs separately from repeated solve costs. Reuse is not a zero-allocation guarantee.

NSDEBase.initialize_solution(problem, parareal::Parareal; history=:none) :: PararealSolution

returns a PararealSolution for problem. Reusable across solves.

Build changed Parareal settings without changing the original.

solve!(solution::AbstractTimeParallelSolution, problem, solver::AbstractTimeParallelSolver; kwargs...) :: AbstractTimeParallelSolution

computes the solution of problem using solver.

solve(problem, solver::AbstractTimeParallelSolver; kwargs...) :: AbstractTimeParallelSolution

computes the solution of problem using solver.

TimeParallelSolution(problem::AbstractInitialValueProblem, parareal::Parareal)

returns a PararealSolution constructor for the solution of problem with parareal.

Wnorm(iterate, reference, weight; base=:time)

Weighted Euclidean distance at the chunk interfaces. reference may be a solution or the vector returned by reference(problem, parareal). The base is per unit time, matching all weight rules. base=:chunk explicitly selects the per-chunk convention. Supply a scalar or FixedWeight; adaptive rule objects contain only their initial settings and are refused. To use a measured weight, pass solution.diagnostics.weights[end], or call Wnorm(solution, reference).

Wnorm(solution::PararealSolution, reference)

Weighted interface error using the last measured per-unit-time weight. A run with no finite fine batch has no measured weight and raises ArgumentError.

returns whether a criterion accepts its evaluated residual.

boundarytimes(iterate::PararealIterate)

returns the N + 1 chunk-boundary times of iterate: each chunk's start, plus the final chunk's end. These are the points Parareal iterates on — the grid every convergence measure (Wnorm, per-iterate error studies) is taken over.

boundaryvalues(iterate::PararealIterate)

returns the N + 1 chunk-boundary values of iterate: each chunk's first state, plus the final chunk's last state — the U vector of the Parareal iteration, read off the chunk solutions. The returned states ALIAS the chunk storage; copy them before mutating.

chunk_endpoint(sol, t)

returns the value of sol at the chunk endpoint t: the stored final value when the solver stopped exactly at t, otherwise the interpolant (for generic solvers that overshoot).

coarseinit!(cache, parareal; mpi=false)

runs the serial coarse pass over the chunk grid, filling G and — where makeGs allows — seeding the chunk starting values U. Chunks whose makeGs[n] is false keep whatever U[n] already holds (injected guesses, e.g. from a moving-window driver).

collect_iterates!(solution::PararealSolution; directory)

reads the per-iterate chunk files written by the MPI backend under saveiterates back into solution.iterates. Chunks below the diagonal (n < k) are final from earlier iterations and are copied forward.

contractionrate(errors::AbstractVector{<:Real}) :: Float64
contractionrate(solution::AbstractTimeParallelSolution) :: Float64

the observed per-iteration contraction rate β of an error trace: the least-squares slope of log(errors[k]) against k, exponentiated, so that errors[k] ≈ C βᵏ. Non-finite entries and EXACT ZEROS are excluded before fitting — finite termination drives the last error to bitwise zero, which is a property of the algorithm's fixed point, not of its rate. Returns NaN when fewer than two usable points remain: a run that converges in one iteration exhibits no observable rate, and downstream bounds (e.g. the outer radius R = r/βᴷ) must treat that as unbounded rather than fabricate a number.

correct!(cache, parareal, k; mpi=false)

the serial Parareal correction sweep after fine batch k: re-run the coarse solver from the freshest starts and update U[n+1] = F[n+1] + (Gnew − Gold), writing in place against the U_ snapshot maintained by the iteration loop.

Two details make finite termination exact in floating point, not just in real arithmetic:

• Order of operations. F + (Gnew − Gold) returns F bitwise whenever the two coarse values cancel bitwise; (Gnew + F) − Gold does not — it rounds F into Gnew first and loses low bits (with Gnew = Gold = 1e16 and F = 1 it returns 0). Do not "simplify" this expression.
• The newly exact interface is assigned, not computed. After batch k the start U[k] is already exact, so F[k+1] is the exact value of U[k+1] and U[k] has not moved since G[k+1] was last computed: the coarse re-run at n = k would reproduce G[k+1] bitwise. Skip it and copy F[k+1] across — one coarse solve saved per sweep, and no arithmetic on the exact prefix at all.

costratio(problem, parareal; repeats=3) :: Real

measures ζ, the coarse/fine cost ratio over ONE chunk of problem under parareal's chunking — the input theoretical_speedup needs. Both solvers are timed on the first chunk with @elapsed (best of repeats, after a compile warm-up). Measure ζ rather than guessing it from step-size ratios: those ignore per-step cost differences between the two solvers.

flatten(iterate::PararealIterate)

concatenates the chunk solutions into one pair of vectors (u, t), dropping each chunk's FIRST point after chunk 1: boundaries are stored twice (chunk n's last point and chunk n+1's corrected start), and keeping both would double-weight every boundary in any downstream statistic. The kept value is the fine endpoint; the dropped corrected start agrees with it only to the solve tolerance. t is monotone.

flatten(solution::PararealSolution)

flattens the last iterate: one (u, t) pair for the whole span, chunk boundaries deduplicated.

iteration_budget(S::Real, N::Integer, ζ::Real; estimate=:wallclock) :: Int

returns the largest fine-batch budget K ∈ 1:N that reaches the target modelled speed-up S; returns zero if none does. A returned K permits K-1 corrections. This is a cost budget, not a prediction that the iteration converges within that budget.

Integer bisection inverts theoretical_speedup with the same cost model and exact-prefix coarse count. Both :wallclock (default) and :work are supported; :naive and :improved are aliases for them, respectively.

maxseam(iterate::PararealIterate) :: Real

returns the largest chunk-boundary seam of iterate (see seams), or 0.0 for a single-chunk iterate.

maxseam(solution::PararealSolution) :: Real

maxseam of the last iterate.

numchunks(iterate::PararealIterate)

returns the number of chunks of iterate.

numchunks(solution::PararealSolution)

returns the number of chunks of solution.

numiterates(solution::PararealSolution)

returns the number of iterates of solution.

parareal!(cache, solution, problem, parareal, backend; saveiterates, directory, nocollect)

runs Parareal with the given AbstractPararealBackend. Under MPI this is SPMD: every rank calls it; is_root guards the serial parts.

prepare_chunks!(cache::PararealCache, parareal::Parareal)

rebuilds the chunk schedules and solver caches after a change of window or method.

prepare_problem!(cache::PararealCache, problem, parareal)

refreshes the initial data and RHS caches from problem; moving the time grid needs shiftwindow!.

reference(problem, parareal::Parareal) :: Vector

returns the serial fine solution at the N + 1 chunk interfaces, computed on exactly the same absolute-time chunk grid and schedule as parareal.

evaluates a criterion from the current fine batch and its cached starts.

seams(iterate::PararealIterate) :: Vector{<:Real}

returns the N − 1 chunk-boundary seam sizes of iterate: at each interior boundary, ‖(start of chunk n+1) − (end of chunk n)‖ — the distance between the corrected starting value U[n+1] and the fine endpoint F(U[n]) it should agree with at convergence. This is exactly the discontinuity that flatten hides by dropping the corrected start, so any statistical claim made on flattened output should report maxseam alongside it: under the weighted criterion ψ₂, late-window seams are unconstrained BY DESIGN (the discount forgives them), and only the seam sizes say which regime a run actually exercised. Empty for a single-chunk iterate.

seams(solution::PararealSolution) :: Vector{<:Real}

seams of the last iterate.

shiftwindow!(cache::PararealCache, τ0, τN)

re-targets a cache at the window [τ0, τN]: recomputes the chunk grid T and updates every chunk problem's tspan in place. This is the official seam for window drivers (e.g. NSDEMovingWindow) — chunk problems are built once at cache construction, so their time spans must be moved through this function, never by poking cache.T alone.

theoretical_speedup(K::Integer, N::Integer, ζ::Real; estimate=:wallclock) :: Real
theoretical_speedup(diagnostics::PararealDiagnostics, ζ::Real; estimate=:wallclock) :: Real
theoretical_speedup(solution::PararealSolution, ζ::Real; estimate=:wallclock) :: Real

returns the modelled speed-up over N serial fine chunks. K is the number of fine batches, not corrections: a finite bulk return has s = K-1, with 1 ≤ K ≤ N. The cost ratio ζ is coarse/fine time per chunk; see costratio.

The exact-prefix implementation uses K_G = K*(N-1) - K*(K-1)/2 coarse calls and K*N - K*(K-1)/2 fine calls, including the initial coarse pass and final fine batch. The diagnostics overload uses the recorded call counts directly.

• :wallclock (default) charges one fine-chunk time per batch and sequential coarse calls: N/(K + ζ*K_G). It assumes at least N workers and equal chunk costs, and excludes communication, synchronization, and load imbalance.
• :work divides the total fine-call cost by N, retaining the same sequential coarse cost. This assumes perfect redistribution of the shrinking fine work; it is not the critical path of the bulk backend.

These are cost models, not measured speed-ups or unconditional performance bounds. The compatibility keyword wallclock=false selects :work, while an explicit estimate takes precedence. Pipeline diagnostics, nonfinite or unstarted runs, and runs with no fine batch are refused: a pipeline needs its task schedule.