solprop.backends.low_storage

The low-storage interaction-picture RK4 march: replay-only, accepted in place.

SEP 24 (docs/dev/seps/archive/, sep-0024-low-storage-interaction-picture-runge-kutta.md). Classical fourth-order Runge-Kutta in the interaction picture anchored at the step’s midpoint, Lawson’s transformation [1] in Hult’s arrangement [2]. For a step [z_n, z_n + h] with z_m = z_n + h/2, P = U(z_m, z_n), Q = U(z_n + h, z_m) and U(z_2, z_1) = exp(Phi(z_2) - Phi(z_1)), the step is

g1 = P F(z_n, A_n), F2 = F(z_m, P A_n + h/2 g1), F3 = F(z_m, P A_n + h/2 F2), F4 = F(z_n + h, Q (P A_n + h F3)), A_{n+1} = Q (P A_n + h/6 g1 + h/3 F2 + h/3 F3) + h/6 F4,

with F the complete source (the response, the in-stage absorbers and the modified noise draw at the evaluation’s position). Only forward flows are applied, four per element, and on a constant operator P = Q = exp(h D / 2).

The register form (SEP 24 section 4.3) holds three field-sized registers beside the source pipeline’s scratch S: the state A, overwritten in place by P A and then by the accepted step, the accumulator B and the stage register C. A step is four passes, each forming an input, evaluating F into S and folding it into the registers:

pass	evaluation	finish
1	F1 at z_n	A <- P A; g <- P S; B <- A + h/6 g; C <- A + h/2 g
2	F2 at z_m	B <- B + h/3 S; C <- A + h/2 S
3	F3 at z_m	B <- B + h/3 S; C <- Q (A + h S)
4	F4 at end	A <- Q B + h/6 S

The march (LowStorageMarch) replays a prescribed schedule and never rejects, so it accepts each step in place: from the step’s start to its commit the backend is in flight (StepDriver._in_flight) and every entry that reads or advances the state refuses it, until StepDriver.set_initial loads a state again (SEP 24 section 6.2, ADR-050). It records no error estimate: the step history’s norm is 0.0, “not estimated” (section 6.3). Divergence is detected on the accepted state’s norm instead.

The engine holds the registers and makes the passes (LowStorageEngine); a backend may keep its own where its state lives (StepDriver._low_storage_engine). HostLowStorage is the host engine and the oracle, over the backend’s _adams_lab_source and _adams_propagator.

References

[1]

J. D. Lawson, “Generalized Runge-Kutta processes for stable systems with large Lipschitz constants”, SIAM Journal on Numerical Analysis 4, 372 (1967). doi:10.1137/0704033

[2]

J. Hult, “A fourth-order Runge-Kutta in the interaction picture method for simulating supercontinuum generation in optical fibers”, J. Lightwave Technol. 25, 3770 (2007). doi:10.1109/JLT.2007.909373

Functions

schedule_key(keys, edges, z, dz)	Return the table key of the schedule step from z, from table_steps().
table_steps(schedule)	Return the step each of a schedule's steps builds and keys the flow table on.

Classes

HostLowStorage(backend, method)	The host engine: the four register passes over host arrays (SEP 24 6.5).
LowStorageEngine(*args, **kwargs)	What LowStorageMarch asks of the engine that holds its registers.
LowStorageMarch(backend, method)	Replay a schedule with the low-storage RK4IP march, accepting in place.

class solprop.backends.low_storage.HostLowStorage(backend, method)

Bases: object

The host engine: the four register passes over host arrays (SEP 24 6.5).

The backend’s state is the register A and is updated in place; B and C are this engine’s host arrays, allocated at the first step. Each pass makes one source evaluation through StepDriver._adams_lab_source (the complete source, the noise draw advanced to the pass’s position and the in-stage absorbers included) and takes P and Q from StepDriver._adams_propagator or the table, allocating transients as the reference does. The step’s three positions are known before its first pass, so their host scalars are read together (StepDriver._stage_bundles) and the hooks take the kept bundles.

The flows. "tabulate" on a constant operator holds one volume, exp(h_t D / 2), built at the first step after a load and reused while the steps keep the same key h_t, their run’s mean step (table_steps(), SEP 24 Q10): a uniform schedule’s steps are differences of rounded positions and differ in their last bits. The table is built from the key alone, from the operator’s origin to h_t / 2 – a constant operator’s increment is the distance itself – so its exponent carries one rounding of D h_t / 2 and none of the step’s position (fl(z + h/2) - z would be off by up to half an ulp of z, an error that accumulates). Each reused flow’s exponent is then off by |D| |h - h_t| / 2, a few ulps of the position times |D|, which does not accumulate over the run. A varying operator’s P and Q change every step, and "tabulate" forms the pair once per step; "recompute" forms each flow in the finish that applies it, from the positions, whose differences telescope. The treatments agree to rounding, not bit for bit.

The stiffness indicator (SEP 24 section 6.3, Q1): h rho = 2 ||F3 - F2|| / ||F2 - g1|| over the norm mask estimates |h lambda| of the linearised source on the content the step integrates (for u' = lambda u it is exactly |h lambda|), the quantity the stability interval of section 4.6 bounds. It is kept in stiffness and never enters the norm column. It reads the content’s margin, not the worst case: measured on an ionising capillary it passed 2.785 on the step before the replay diverged, and beside an in-stage collar it read 0.17 to 0.65 of h alpha_max / 2 as the collar held less or more of the field.

Parameters:

backendsolprop.backends.driver.StepDriver

The backend whose state, hooks and boundary treatments are used; its state is a host array in the canonical layout.

methodsolprop.stepper.low_storage.LowStorageRungeKutta

The march.

Attributes:

stiffnessfloat

h rho of the last completed step, dimensionless; NaN before one, or when F2 - g1 vanished on the norm mask.

atomicbool

Always False: the passes are Python calls a signal may interrupt.

Notes

What the engine holds between steps is the registers B and C and, on a constant operator tabulating, the table: held() names them for the backend’s memory ledger, and close() releases them.

close()

Release the registers and the flow table.

held()

Return the bytes of every array the engine holds between steps, by name.

Returns:

dict

register_accumulator, register_stage and propagators (the table), each zero while it is not allocated.

reset()

Drop the flow table, which a resume builds again from its own first step.

step(z, dz)

Take one RK4IP step in the four passes of the register form.

Parameters:

zfloat

Where the step starts, m.

dzfloat

The step, m, positive.

Returns:

float

The accepted state’s squared norm over the norm mask, before the boundary treatments.

class solprop.backends.low_storage.LowStorageEngine(*args, **kwargs)

Bases: Protocol

What LowStorageMarch asks of the engine that holds its registers.

Part of the kernel seam (schema 29): the march keeps the schedule, the commit and the records; an engine keeps the registers and the flow table and does every state-sized operation of a step. HostLowStorage is the reference.

atomic: bool

Whether step() is one native call that no signal interrupts (schema 32).

A signal’s Python handler runs at the first bytecode after such a call returns, which is before the commit, so the march then takes the step, the backend’s _finish and the position move in one committing() block (SEP 24 section 6.7): only a native exception can abandon the step in flight. An engine whose step runs Python between its passes is not atomic, and a signal there abandons the step at once rather than waiting up to a step for the commit (section 6.2, Q6).

close()

Release the registers and the flow table; the engine is not used again.

reset()

Forget what a step left behind: the state was loaded outside the march.

step(z, dz)

Take one step in place, leaving the raw accepted state in the backend.

Parameters:

zfloat

Where the step starts, m.

dzfloat

The step, m, positive.

Returns:

float

The accepted state’s squared norm over the norm mask, before the boundary treatments, in the state’s units squared; not finite when the field has diverged.

stiffness: float

h rho of the last completed step (SEP 24 section 6.3), or NaN where the engine does not form it: the march warns when it passes RK4’s stability interval.

class solprop.backends.low_storage.LowStorageMarch(backend, method)

Bases: object

Replay a schedule with the low-storage RK4IP march, accepting in place.

Follows solprop.backends.multistep.AdamsMarch: it owns its registers through an engine and replays a prescribed schedule, which is the only thing it does (SEP 24 section 6.2). Each step calls the engine, applies the backend’s _finish (the hard band mask, and the split-step collar where chosen) and moves the position in one commit, then records the step through the backend.

Parameters:

backendsolprop.backends.driver.StepDriver

The backend whose state is marched; it owns the step history and the statistics rows.

methodsolprop.stepper.low_storage.LowStorageRungeKutta

The march.

Attributes:

source_evaluationsint

Source evaluations of the committed steps since the state was last loaded (StepDriver.set_initial, which sets it to zero), four a step: a physical count, never a logical response call (SEP 24 section 5). A step that did not commit – abandoned in flight, or diverged – adds none.

Notes

The march has no error estimate to warn on. It logs a warning, once per call of replay_to(), when the engine’s stiffness indicator passes RK4’s negative-real stability interval (LOW_STORAGE_REAL_INTERVAL): an explicitly integrated damping term the translation cannot bound in advance – plasma loss – is then near the step’s stability limit.

close()

Release the engine’s registers and table, with the backend’s buffers.

property engine: LowStorageEngine

Return the engine that holds the registers.

Returns:

LowStorageEngine

The host engine, or the backend’s own.

replay_step(dz, end, slack)

Take one prescribed step from the backend’s position and accept it.

The unit a lockstep batch advances each member by (solprop.backends.batch.LockstepRunner), with solprop.backends.multistep.AdamsMarch.replay_step()’s signature.

Parameters:

dzfloat

The step, m.

endfloat

The schedule’s edge it ends on, m.

slackfloat

How far an edge may sit from a target it lands on, m; unused, since the march restarts nothing at a kink.

Returns:

float

The norm the history records: 0.0, “not estimated”.

Raises:

BackendError

If an earlier step was abandoned in flight.

ConvergenceError

If the accepted state is not finite: the field has diverged, and a replay cannot reject the step. Raised in flight, so the state is lost with it.

replay_to(schedule, target)

Replay the schedule’s steps to target, mirroring the Adams replay.

Parameters:

schedulesolprop.stepper.schedule.FixedSchedule

The schedule to replay.

targetfloat

Position to reach, m.

Raises:

BackendError

If the schedule runs out before the target is reached.

reset()

Forget what the engine holds: the state was loaded outside the march.

snapshot()

Return nothing for a resume to continue: the march keeps no history.

A resume at a schedule edge needs no history and no recycled stage (SEP 24 section 5), as under RungeKutta().

solprop.backends.low_storage.schedule_key(keys, edges, z, dz)

Return the table key of the schedule step from z, from table_steps().

Found by position rather than by a step counter, which a lockstep member does not advance (LockstepRunner keeps the batch’s).

Parameters:

keysnumpy.ndarray

table_steps() of the schedule.

edgesnumpy.ndarray

The schedule’s edges, m.

zfloat

Where the step starts, m.

dzfloat

The step, m.

Returns:

float

The key, m.

solprop.backends.low_storage.table_steps(schedule)

Return the step each of a schedule’s steps builds and keys the flow table on.

SEP 24 Q10, decided in phase 2: the mean step of the step’s run, the longest stretch of consecutive steps within LANDING_ULPS ulps of the run’s first, the ulps being of the run’s farthest position (or of the step, near the origin). uniform() takes a segment’s steps from differences of rounded positions, which differ in their last bits by up to an ulp of the position: far from the origin that is many ulps of the step, and a table keyed on one step’s own value within ulps of the step was rebuilt nearly every step. Keyed on the run’s mean h_bar instead, a uniform segment builds one table and reuses it; each flow’s exponent is then off by |D| |h - h_bar| / 2, a few ulps of the position times |D|, and the error does not drift: the steps of a run sum to n h_bar to rounding, so the phase the flows apply over a run is the run’s own. A schedule whose steps all differ (a discovered one) has runs of one step and keys each step on itself.

Parameters:

schedulesolprop.stepper.schedule.FixedSchedule

The schedule replayed.

Returns:

numpy.ndarray

One key per step, m, shape (n_steps,): constant over a run.

Examples

>>> from solprop.backends.low_storage import table_steps
>>> from solprop.stepper import FixedSchedule
>>> schedule = FixedSchedule.uniform(0.0, 1.0, 0.1, save_planes=[0.25, 1.0])
>>> keys = table_steps(schedule)
>>> len(set(keys.tolist())), round(float(keys[0]), 12), round(float(keys[-1]), 12)
(2, 0.083333333333, 0.09375)