solprop.backends.numpy_ref.passes#

The passes of one interaction-picture stage, one function each.

docs/dev/architecture.md section 8 lists the anatomy of a stage as seven numbered passes over the data, and docs/dev/backends.md section 2 says the reference backend writes each of them as a named function with its formula in the docstring. That is this module. Nothing here is fused, nothing is written in place and nothing allocates a buffer it keeps: every function takes arrays and returns a new one, so a reader can follow one pass without holding the rest of the step in mind. The optimised backends fuse these passes and are validated against the numbers this module produces.

The sequence is: pad, inverse transform, conversion to volts per metre, the responses, the temporal window, the forward transform, the crop, the prefactor and the spectral window.

Scales. The oversampling pad scale is already inside transform.synthesis and the crop scale is already inside responses.prefactor (docs/dev/kernel-contract.md sections 4 and 5), so to_time() and to_frequency() scatter and gather through grids.crop_index and apply no scale of their own. Using solprop.fourier.oversample() and solprop.fourier.crop() here instead would count both factors twice.

Convention. Every complex array is in the FFT convention of docs/dev/conventions.md section 2; nothing here conjugates. The Raman kernel spectrum is an internal-convention array whose imaginary part is what distinguishes gain from loss, so it is never passed through solprop.fourier.to_optics() (docs/guide/media/raman-numerics.md, “The sign that means gain”).

The one physical constant. solprop.constants.EPS0 is read here, by the Raman pass alone. Every other constant reaches a backend through a compiled coefficient – eps0 chi3 for the Kerr term, e^2/m_e for the plasma – but the Raman response writes eps0 outside the convolution integral, so folding it into a row would change the stated units of raman_kernel_undamped and cost the bitwise agreement with solprop.nonlinear.raman.RamanResponse.evaluate() that pins this pass.

Functions

absorber_term(modal_t, temporal_rate_fine, ...)

Return the temporal absorber as a contribution to the stage derivative.

accept_step(state, propagator, dz, accumulated)

Advance the state over a completed step.

accumulate(total, stage, weight)

Add one weighted stage derivative to a running accumulator.

apply_hard_mask(state, active_mask)

Zero the state outside the active band.

apply_temporal_rate(state, rate, dz, *, kind)

Apply the temporal absorber as an exact rate over one accepted step.

dispersive_kerr_polarisation(field_t, terms, ...)

Evaluate SEP 6's dispersive electronic polarisation at every node of a chunk.

error_norm(error, state, state_new, ...)

Reduce the embedded error estimate to the scaled norm the controller consumes.

held_plasma_spectrum(gated, edge, omega, ...)

Transform a current held at the temporal collar (SEP 27 section 4.4).

ionisation_rate(field_t, table)

Look the ionisation rate up on the magnitude of the field.

kerr_polarisation(field_t, kerr, density[, ...])

Evaluate the instantaneous Kerr polarisation at every node of a chunk.

mixture_raman_spectrum(components, rates, ...)

Form one amplitude-weighted kernel spectrum for all Raman components.

nonlinear_polarisation(field_t, responses, ...)

Run the fused response pass over a chunk of nodes and window it.

nonlinear_sources(field_t, responses, ...[, ...])

Separate the collared local polarisation and the ungated plasma current.

plasma_channel_fractions(field_t, channels, ...)

Evaluate independent exposure scans in channel order.

plasma_channels_polarisation(field_t, ...[, ...])

Sum channel weights before the two shared current integrals.

plasma_current(field_t, rate, fraction, ...)

Add the ionisation loss to the free current and integrate to a polarisation.

plasma_integration_diagnostics(field_t, ...)

Return the parts of the column's integration diagnostics (SEP 18).

plasma_integration_statistics(parts)

Form the three diagnostics from plasma_integration_diagnostics()'s parts.

plasma_scan(field_t, rate, coefficients, *, ...)

Integrate the ionised fraction and the free-electron current along t.

project(polarisation_t, projection, taper, n_pol)

Project a chunk's polarisation back onto the modes.

raman_polarisation(field_t, spectrum, *, ...)

Evaluate the delayed Raman polarisation at every node of a chunk.

raman_spectrum(undamped, rates, fraction, dt)

Damp each undamped Raman row at its own rate, sum, scale and transform.

stage_input(state, stage_sum, dz, ...)

Combine the stages, return to the lab frame and add the noise.

stage_output(polarisation, prefactor, ...[, ...])

Normalise the polarisation and return it to the interaction picture.

synthesise(modal_t, synthesis, taper, n_pol)

Synthesise the physical field on a chunk of quadrature nodes.

tensor_raman_spectrum(components, rates, ...)

Form effective trace and traceless spectra from independently damped terms.

to_frequency(polarisation_t, crop_index, *, kind)

Transform the modal polarisation to frequency and crop to the coarse grid.

to_time(spectrum, crop_index, n_to, *, kind)

Pad onto the fine grid and transform every modal column to time.

Classes

NonlinearSources(polarisation, current)

Carry the separate local polarisation and field-plasma current.

class solprop.backends.numpy_ref.passes.NonlinearSources(polarisation, current)[source]#

Bases: object

Carry the separate local polarisation and field-plasma current.

Attributes:
polarisationnumpy.ndarray

Collared Kerr and Raman polarisation, including the separately collared dispersive Kerr, in C/m^2, shape (n_nodes_chunk, n_pol, n_to).

currentnumpy.ndarray or None

Total ungated plasma current density in A/m^2 with the same shape; None when no field plasma is compiled. Every channel and both currents enter once.

current: Float[ndarray, 'n_nodes_chunk n_pol n_to'] | None#
polarisation: Num[ndarray, 'n_nodes_chunk n_pol n_to']#
solprop.backends.numpy_ref.passes.absorber_term(modal_t, temporal_rate_fine, crop_index, active_mask, *, kind)[source]#

Return the temporal absorber as a contribution to the stage derivative.

The "nonlinear" absorber scheme: -alpha(t) E(t) / 2 on the modal time-domain field, transformed and cropped so that it can be added to the stage in the same place the normalised polarisation is. alpha is a power coefficient, so the field carries the half (solprop.boundaries.nonlinear_absorber_term()).

It cannot ride the polarisation’s own transform, and this is the whole reason the two schemes exist: the polarisation is multiplied by the UPPE prefactor after the forward transform, and the absorber is already a rate on the field rather than a polarisation, so it would have to be divided by a prefactor that is exactly zero outside the active band. One extra modal forward transform per stage is its price, and it is what buys FSAL, since nothing then changes the state after the last stage of a step.

The fine grid is used rather than the coarse one because alpha(t) E(t) has more bandwidth than E(t) and the coarse grid would alias the product back into the band; the contract carries both axes for exactly this reason (kernel-contract.md section 5).

The result is masked to the active band, as the polarisation beside it already is through a prefactor that is exactly zero there. alpha(t) E(t) has content at every frequency, and the active mask is the statement that the equation is not solved outside that band: leaving it in would put content into the state that the accept pass’s hard mask then removes, which is both meaningless and enough to break the first-same-as-last identity (measured at 7.4e-7 before the mask and 1.8e-16 after it).

Parameters:
modal_tnumpy.ndarray

The lab-frame modal field on the fine time grid, shape (n_modes, n_pol, n_to), as to_time() produced it.

temporal_rate_finenumpy.ndarray

The collar rate on the fine time axis, 1/m, non-negative, shape (n_to,).

crop_indexnumpy.ndarray

Where each coarse bin sits in the fine spectrum, numpy.intp, shape (n_omega,).

active_masknumpy.ndarray

The band the equation is solved on, Boolean, shape (n_omega,).

kind{‘field’, ‘envelope’}

Spectral layout.

Returns:
numpy.ndarray

The absorber’s contribution on the coarse grid, complex128, shape (n_modes, n_pol, n_omega), exactly zero outside the active band, in the same units as the stage derivative.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import absorber_term
>>> field = np.ones((1, 1, 8))
>>> rate = np.full(8, 2.0)
>>> index = np.arange(2, dtype=np.intp)
>>> absorber_term(field, rate, index, np.ones(2, bool), kind="field").real
array([[[-8.,  0.]]])
solprop.backends.numpy_ref.passes.accept_step(state, propagator, dz, accumulated)[source]#

Advance the state over a completed step.

y <- exp(dPhi) (y + h sum_i b_i k_i). The interaction picture is anchored at the start of the step, so the whole interaction-variable solution – the state and the accumulated stages together – is carried to the end of the step by one exponential. The same factor applied to h sum_i (b_i - bhat_i) k_i gives the local error of the propagated state, which is what error_norm() scales.

Parameters:
statenumpy.ndarray

The state at the anchor, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts.

propagatornumpy.ndarray

exp(dPhi) over the whole step, dimensionless, broadcastable against the state.

dzfloat

Length of the step, m, positive.

accumulatednumpy.ndarray

sum_i b_i k_i, same shape as the state, in the square root of watts per metre.

Returns:
numpy.ndarray

The state at the end of the step, same shape and units.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import accept_step
>>> state = np.ones((1, 1, 1, 2), dtype=complex)
>>> accept_step(state, 2.0 * np.ones_like(state), 0.5, state).ravel()
array([3.+0.j, 3.+0.j])
solprop.backends.numpy_ref.passes.accumulate(total, stage, weight)[source]#

Add one weighted stage derivative to a running accumulator.

total + w k_i. The solution and the embedded error are running accumulators rather than materialised stage combinations, so nothing beyond the stages the next one needs is kept.

Parameters:
totalnumpy.ndarray

The accumulator so far, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts per metre.

stagenumpy.ndarray

The stage derivative k_i, same shape and units.

weightfloat

The tableau weight: b_i for the solution, b_i - bhat_i for the error.

Returns:
numpy.ndarray

The updated accumulator, same shape and units.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import accumulate
>>> total = np.zeros((1, 1, 1, 2), dtype=complex)
>>> accumulate(total, np.ones_like(total), 0.5).ravel()
array([0.5+0.j, 0.5+0.j])
solprop.backends.numpy_ref.passes.apply_hard_mask(state, active_mask)[source]#

Zero the state outside the active band.

y <- y * active_mask, once per accepted step. The hard band limit is the one boundary that stays a 0/1 mask rather than becoming a rate, because it is idempotent: the operator and the prefactor are already exactly zero outside the active band, so nothing regrows there and the total absorption does not depend on the step count. The active mask, not the narrower in-band one, which would erase the spectral collar the rate acts on.

Parameters:
statenumpy.ndarray

The state, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts.

active_masknumpy.ndarray

The active band, omega_min/2 < omega < 1.1 omega_max, Boolean, shape (n_omega,).

Returns:
numpy.ndarray

The masked state, same shape and units.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import apply_hard_mask
>>> state = np.ones((1, 1, 1, 3), dtype=complex)
>>> apply_hard_mask(state, np.array([True, True, False])).ravel()
array([1.+0.j, 1.+0.j, 0.+0.j])
solprop.backends.numpy_ref.passes.apply_temporal_rate(state, rate, dz, *, kind)[source]#

Apply the temporal absorber as an exact rate over one accepted step.

y <- F[F^-1[y] exp(-alpha(t) h / 2)], the split-step absorber scheme. The absorber cannot be folded into Phi because it is not diagonal in omega, and a per-step window would give total absorption W^n and no converged answer; as a rate the factors telescope, so the total is exp(-alpha(t) L / 2) however the steps fall. alpha is a power coefficient, hence the half on the field.

Parameters:
statenumpy.ndarray

The state, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts.

ratenumpy.ndarray

alpha(t) on the coarse time grid, float64, shape (n_t,), in 1/m, non-negative.

dzfloat

Length of the accepted step, m, positive.

kind{‘field’, ‘envelope’}

Spectral layout.

Returns:
numpy.ndarray

The absorbed state, same shape and units.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import apply_temporal_rate
>>> state = np.array([[[[4.0 + 0.0j, 0.0j, 0.0j]]]])
>>> absorbed = apply_temporal_rate(state, np.zeros(4), 1.0, kind="field")
>>> bool(np.allclose(absorbed, state))
True
solprop.backends.numpy_ref.passes.dispersive_kerr_polarisation(field_t, terms, densities, *, kind, carrier=None, window=None)[source]#

Evaluate SEP 6’s dispersive electronic polarisation at every node of a chunk.

The dispersive_kerr stage of the fused pass: per channel, the transforms of the field the evaluator’s forms need, the products of the filtered fields in the leaf’s sector, the temporal collar on each product, the output multipliers and one inverse transform, by solprop.nonlinear.responses.dispersive_cubic(), which is the arithmetic of the reference evaluate() on the leaf’s sampled multipliers. The polynomial’s channels add, P = sum_s c_s rho_s cubic_s(E); the symmetric evaluator’s mixture takes one filter built from its species’ densities at this stage, P = F G_mix[w (G_mix E)^3] (evaluations()). It reads the bare field and touches no other response’s source.

Parameters:
field_tnumpy.ndarray

The field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m: real on a field grid, complex on an envelope grid.

termssolprop.translate.responses.DispersiveKerrTerms

The compiled leaf.

densitiessequence of float

Each channel’s number density at this position, m^-3, in the leaf’s order.

kind{‘field’, ‘envelope’}

The spectral layout.

carriernumpy.ndarray or None, optional

exp(2 i omega_0 t) on the fine time grid, shape (n_to,); required when the leaf’s carrier is set.

windownumpy.ndarray or None, optional

The temporal collar, shape (n_to,), applied to each product before the output multipliers, G[w (G E)^3]; None is a unit collar.

Returns:
numpy.ndarray

The polarisation in C/m^2, shape (n_nodes_chunk, n_pol, n_to), real on a field grid and complex on an envelope grid, the collar already inside it.

Raises:
BackendError

If the densities do not match the channels, or the carrier is missing.

solprop.backends.numpy_ref.passes.error_norm(error, state, state_new, norm_mask, controller)[source]#

Reduce the embedded error estimate to the scaled norm the controller consumes.

||e||_2 / (rtol max(||y||_2, ||y_new||_2) + atol) over the in-band mask, one global scale – the weak norm of docs/dev/stepper.md section 3. When the controller carries weights the scale is per bin instead and the division happens before the reduction; the magnitude stays global either way, because a per-element magnitude collapses the step to zero on a spectrum whose in-band dynamic range reaches 1e17 (section 10, trap 1).

Parameters:
errornumpy.ndarray

The embedded difference of the propagated state, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts.

state, state_newnumpy.ndarray

The state before and after the trial step, same shape and units.

norm_masknumpy.ndarray

The bins the norm scores (grids.norm_mask, the controller’s region), Boolean, shape (n_omega,).

controllersolprop.stepper.controller.PIController

The controller, which carries the tolerances and the choice of norm.

Returns:
float

The scaled norm, dimensionless. One means the step sat exactly at tolerance.

Examples

>>> import numpy as np
>>> from solprop.stepper import PIController
>>> from solprop.backends.numpy_ref.passes import error_norm
>>> state = np.ones((1, 1, 1, 2), dtype=complex)
>>> mask = np.ones(2, dtype=bool)
>>> control = PIController(rtol=1e-2, atol=0.0)
>>> round(error_norm(1e-2 * state, state, state, mask, control), 12)
1.0
solprop.backends.numpy_ref.passes.held_plasma_spectrum(gated, edge, omega, crop_index, active)[source]#

Transform a current held at the temporal collar (SEP 27 section 4.4).

Form DFT[g J]/(i omega) - DFT[g' J]/(i omega)^2 only on active bins. Projection can precede or follow this linear map; crop and amplitude factors belong to the route. No running integral or numerical derivative enters it.

Parameters:
gated, edgenumpy.ndarray

Sampled g J and g' J, in A/m^2 and A/(m^2 s), with time last.

omeganumpy.ndarray

Coarse angular frequencies in rad/s, shape (n_omega,).

crop_indexnumpy.ndarray

Coarse positions in the fine half-spectrum, shape (n_omega,).

activenumpy.ndarray

Active frequency mask, shape (n_omega,); all selected frequencies positive.

Returns:
numpy.ndarray

Polarisation-equivalent spectrum in raw DFT units, with inactive bins zero.

solprop.backends.numpy_ref.passes.ionisation_rate(field_t, table)[source]#

Look the ionisation rate up on the magnitude of the field.

The intensity and ionisation_rate stages of the fused pass, together: the magnitude is what the intensity stage forms and the only thing the rate needs. For one component that is |E|; for two it is hypot(E_x, E_y) rather than sqrt(E_x^2 + E_y^2), because the tunnelling exponent is tens to hundreds and turns one unit in the last place of the magnitude into 1e-13 of the rate (solprop.nonlinear.responses.PlasmaResponse.evaluate() takes the same route).

The table is exactly zero below its field_min and holds its last value above field_max rather than extrapolating, which is what keeps the loss current of plasma_current() finite and a diverged step from becoming an infinity.

Parameters:
field_tnumpy.ndarray

The real field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m.

tablesolprop.translate.ionisation.RateTable

The compiled rate, responses.plasma_coefficients.rate.

Returns:
numpy.ndarray

The rate per atom in s^-1, shape (n_nodes_chunk, n_to), finite and non-negative.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import ionisation_rate
>>> from solprop.nonlinear.optmat_rates import OptmatRate
>>> from solprop.translate.ionisation import build_rate_table
>>> table = build_rate_table(OptmatRate("He", 800e-9), (2e9, 1e11), rtol=1e-4)
>>> field = np.array([[[5e10, 1e6]]])
>>> np.round(np.log10(ionisation_rate(field, table)[0, 0]), 3).item()
11.207
solprop.backends.numpy_ref.passes.kerr_polarisation(field_t, kerr, density, carrier=None)[source]#

Evaluate the instantaneous Kerr polarisation at every node of a chunk.

The kerr stage of the fused pass. The formulas are the two cubic forms of solprop.translate.responses.KerrCoefficients, which were read off solprop.nonlinear.responses.KerrResponse.evaluate() rather than re-derived:

real

P_p = c rho (sum_q E_q^2) E_p – the full cubic of a real field, which already carries the third harmonic and the counter-rotating terms.

analytic

P_p = c rho [s (sum_q |a_q|^2) a_p + k (sum_q a_q^2) conj(a_p)], plus c rho m exp(2 i omega_0 t) (sum_q a_q^2) a_p when the envelope grid radiates the harmonic, with s, k and m the self, conjugate and harmonic factors. For one component the first collapses to (3/4)|a|^2 a.

The sums run over the polarisation axis of one node, which is what the vector-isotropic coupling means: two components of the field at one point of the cross-section, never two different points.

Parameters:
field_tnumpy.ndarray

The field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m: real on a field grid, complex on an envelope grid.

kerrsolprop.translate.responses.KerrCoefficients

The compiled coefficients, whose coefficient is eps0 chi3 per unit of the density path.

densityfloat

The density path at this position, responses.density(z): a number density in m^-3 for a gas, 1.0 for a medium whose chi3 is already absolute. Non-negative.

carriernumpy.ndarray or None, optional

exp(2 i omega_0 t) on the fine time grid, shape (n_to,). Required when kerr.harmonic_factor is not zero; the harmonic term is added whenever a carrier is given, and None selects the fundamental cubic alone.

Returns:
numpy.ndarray

The polarisation in C/m^2, shape (n_nodes_chunk, n_pol, n_to), real on a field grid and complex on an envelope grid.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import kerr_polarisation
>>> from solprop.translate.responses import KerrCoefficients
>>> coefficients = KerrCoefficients(1.0, "analytic", 0.5, 0.25, 0.0, False, False)
>>> field = np.array([[[2.0 + 0.0j]]])
>>> polarisation = kerr_polarisation(field, coefficients, 1.0)
>>> bool(np.allclose(polarisation, 0.75 * abs(2.0) ** 2 * 2.0))
True
solprop.backends.numpy_ref.passes.mixture_raman_spectrum(components, rates, densities, dt)[source]#

Form one amplitude-weighted kernel spectrum for all Raman components.

Parameters:
componentstuple of RamanCoefficients

Raman-active populations in mixture order.

ratestuple of float

Decay rates in s^-1, flattened by component then kernel row.

densitiestuple of float

Component amplitude densities in m^-3.

dtfloat

Fine-grid sample spacing in s.

Returns:
numpy.ndarray

One damped kernel spectrum including densities and dt.

Raises:
BackendError

If the amplitude or decay-rate counts disagree with the components.

solprop.backends.numpy_ref.passes.nonlinear_polarisation(field_t, responses, t_window_fine, *, kind, density, preionisation, dt, carrier=None, raman_spectrum=None, scalars=None)[source]#

Run the fused response pass over a chunk of nodes and window it.

Pass 4 of architecture.md section 8: the traversal responses.response_order names, term by term, accumulating P(t) and applying the temporal collar at the end – but for SEP 6’s dispersive response, whose collar acts on its products before its output multipliers (dispersive_kerr_polarisation()) and which is added after. Each stage is one of the named functions of this module, so the order the translation layer compiled is the order that runs, and a term this backend cannot execute is refused at construction rather than dropped here.

Parameters:
field_tnumpy.ndarray

The field at the chunk’s nodes from synthesise(), shape (n_nodes_chunk, n_pol, n_to), in V/m.

responsessolprop.translate.responses.CompiledResponses

Section 5 of the kernel contract: the order, the coefficients and the tables.

t_window_finenumpy.ndarray

The temporal collar profile in [0, 1], float64, shape (n_to,).

kind{‘field’, ‘envelope’}

The spectral layout, which the Raman drive depends on.

densityfloat

The density path at this position, responses.density(z), non-negative.

preionisationfloat

The pre-ionised fraction at this position, in [0, 1]. Ignored when no plasma is compiled.

dtfloat

Sample spacing of the fine time grid, s, positive.

carriernumpy.ndarray or None, optional

exp(2 i omega_0 t) on the fine time grid, for an envelope-grid harmonic term.

raman_spectrumnumpy.ndarray or None, optional

The damped Raman kernel spectrum at this position from raman_spectrum(), shape (n_to + 1,). Required when the order names raman_drive; the caller forms it once per stage rather than once per chunk, because it does not depend on the nodes.

scalarsStageScalars or None, optional

Labelled component densities and preionisations at the stage. Required for multiple plasma channels; None retains the standalone single-response call.

Returns:
numpy.ndarray

The total polarisation in C/m^2, shape (n_nodes_chunk, n_pol, n_to), real on a field grid and complex on an envelope grid.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import nonlinear_polarisation
>>> from solprop.grid import FieldGrid
>>> from solprop.nonlinear.gnlse import SpecifiedMode
>>> from solprop.nonlinear.responses import KerrResponse
>>> from solprop.profiles import Constant
>>> from solprop.translate.responses import compile_responses
>>> grid = FieldGrid.from_limits(1.0, 800e-9, (400e-9, 2e-6), 400e-15)
>>> mode = SpecifiedMode.from_taylor(
...     [0.0], omega_ref=grid.omega0, aeff=1e-10
... )
>>> compiled = compile_responses(
...     [KerrResponse(chi3=1.0)],
...     grid,
...     mode,
...     density=Constant(1.0),
...     z_max=1.0,
...     prefactor_index=1.0,
... )
>>> field = np.full((1, 1, 4), 2.0)
>>> out = nonlinear_polarisation(
...     field, compiled, np.ones(4), kind="field", density=1.0, preionisation=0.0,
...     dt=grid.dt_fine,
... )
>>> bool(np.allclose(out, compiled.kerr_coefficients.coefficient * 8.0))
True
solprop.backends.numpy_ref.passes.nonlinear_sources(field_t, responses, t_window_fine, *, kind, density, preionisation, dt, carrier=None, raman_spectrum=None, scalars=None)[source]#

Separate the collared local polarisation and the ungated plasma current.

Pass 4 of architecture.md section 8: the traversal responses.response_order names, term by term, accumulating P(t) and applying the temporal collar at the end – but for SEP 6’s dispersive response, whose collar acts on its products before its output multipliers (dispersive_kerr_polarisation()) and which is added after. Each stage is one of the named functions of this module, so the order the translation layer compiled is the order that runs, and a term this backend cannot execute is refused at construction rather than dropped here.

Parameters:
field_tnumpy.ndarray

The field at the chunk’s nodes from synthesise(), shape (n_nodes_chunk, n_pol, n_to), in V/m.

responsessolprop.translate.responses.CompiledResponses

Section 5 of the kernel contract: the order, the coefficients and the tables.

t_window_finenumpy.ndarray

The temporal collar profile in [0, 1], float64, shape (n_to,).

kind{‘field’, ‘envelope’}

The spectral layout, which the Raman drive depends on.

densityfloat

The density path at this position, responses.density(z), non-negative.

preionisationfloat

The pre-ionised fraction at this position, in [0, 1]. Ignored when no plasma is compiled.

dtfloat

Sample spacing of the fine time grid, s, positive.

carriernumpy.ndarray or None, optional

exp(2 i omega_0 t) on the fine time grid, for an envelope-grid harmonic term.

raman_spectrumnumpy.ndarray or None, optional

The damped Raman kernel spectrum at this position from raman_spectrum(), shape (n_to + 1,). Required when the order names raman_drive; the caller forms it once per stage rather than once per chunk, because it does not depend on the nodes.

scalarsStageScalars or None, optional

Labelled component densities and preionisations at the stage. Required for multiple plasma channels; None retains the standalone single-response call.

Returns:
NonlinearSources

The collared non-plasma polarisation and total ungated current density.

solprop.backends.numpy_ref.passes.plasma_channel_fractions(field_t, channels, preionisation, dt, integral=None)[source]#

Evaluate independent exposure scans in channel order.

Parameters:
field_tnumpy.ndarray

Real node fields in V/m, shape (n_nodes_chunk, n_pol, n_to).

channelstuple of PlasmaChannel

The independent rate tables and population metadata.

preionisationtuple of float

Initial fractions, in channel order, in [0, 1].

dtfloat

Fine-grid sample spacing in s.

integralsolprop.nonlinear.ionisation.RunningIntegral or None, optional

The running-integral rule; None is the trapezoid.

Returns:
ratesnumpy.ndarray

Instantaneous rates in s^-1, shape (n_channels, n_nodes_chunk, n_to).

fractionsnumpy.ndarray

Ionised fractions with the same shape, dimensionless.

Raises:
BackendError

If the channel scalars do not match their declared count.

solprop.backends.numpy_ref.passes.plasma_channels_polarisation(field_t, channels, densities, preionisation, dt, integral=None)[source]#

Sum channel weights before the two shared current integrals.

Each fraction has its own exposure scan. The free-current weight is sum rho_c f_c e_c²/m_c; the loss numerator is sum rho_c I_p,c W_c (1-f_c). Their time integrals are shared. The original one-channel pass retains its density outside the integrals.

Parameters:
field_tnumpy.ndarray

Real node fields in V/m, shape (n_nodes_chunk, n_pol, n_to).

channelstuple of PlasmaChannel

Ordered independent channel coefficients.

densitiestuple of float

Number density of each population in m^-3.

preionisationtuple of float

Initial fraction of each population, dimensionless in [0, 1].

dtfloat

Fine-grid sample spacing in s.

integralsolprop.nonlinear.ionisation.RunningIntegral or None, optional

The running-integral rule shared by every channel; None is the trapezoid.

Returns:
numpy.ndarray

Plasma polarisation in C/m², with the field shape.

Raises:
BackendError

If there is not one density per channel.

solprop.backends.numpy_ref.passes.plasma_current(field_t, rate, fraction, free_current, coefficients, *, density, dt, integral=None)[source]#

Add the ionisation loss to the free current and integrate to a polarisation.

The plasma_current stage, and the third cumulative integral. J_loss = I_p W (1 - f) E / |E|^2 is added pointwise; its form is fixed by E . J_loss = I_p W (1 - f) = I_p df/dt, so the field pays exactly one ionisation potential for every electron it frees. For one component it is I_p W (1-f)/E, sign and all; for two it is that magnitude directed along the field, which is what makes a linearly polarised pulse give the scalar answer rotated. The polarisation is P = rho int (J_free + J_loss) dt.

The division is written with a where mask rather than a branch, as the oracle writes it, so it never happens at a sample where the rate or the field vanishes: the term is finite only because the rate is exactly zero below threshold, which is what solprop.translate.ionisation.RateTable guarantees below its field_min; a clipped but non-zero rate would divide by a tiny |E|^2.

Parameters:
field_tnumpy.ndarray

The real field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m.

ratenumpy.ndarray

The ionisation rate from ionisation_rate(), s^-1, shape (n_nodes_chunk, n_to).

fractionnumpy.ndarray

The fraction ionised from plasma_scan(), shape (n_nodes_chunk, n_to).

free_currentnumpy.ndarray

J_free from plasma_scan(), shape (n_nodes_chunk, n_pol, n_to).

coefficientssolprop.translate.responses.PlasmaCoefficients

The compiled plasma group; only ionisation_potential is read here.

densityfloat

The density path at this position, responses.density(z), in m^-3, non-negative. The currents are per atom; this is what turns them into a polarisation.

dtfloat

Sample spacing of the fine time grid, s, positive.

integralsolprop.nonlinear.ionisation.RunningIntegral or None, optional

The running-integral rule, responses.plasma_integral; None is the trapezoid.

Returns:
numpy.ndarray

The plasma polarisation in C/m^2, shape (n_nodes_chunk, n_pol, n_to).

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import plasma_current
>>> from solprop.nonlinear.optmat_rates import OptmatRate
>>> from solprop.nonlinear.responses import PlasmaResponse
>>> from solprop.translate.responses import plasma_coefficients
>>> compiled = plasma_coefficients(
...     PlasmaResponse("Ar", OptmatRate("Ar", 800e-9)), (1e9, 1e11), z_max=1.0
... )
>>> zeros = np.zeros((1, 1, 4))
>>> flat = np.zeros((1, 4))
>>> plasma_current(
...     zeros, flat, flat, zeros, compiled, density=1e25, dt=1e-16
... ).tolist()
[[[0.0, 0.0, 0.0, 0.0]]]
solprop.backends.numpy_ref.passes.plasma_integration_diagnostics(field_t, channels, preionisation, dt, integral=None)[source]#

Return the parts of the column’s integration diagnostics (SEP 18).

For every channel and node: the exposure S as the rule forms it, before the floor, and its increments; the fraction; the work the loss current pays per ionisation potential, int W (1 - f) dt by the trapezoid of the raw rate; the electrons that freed, f_end - f_0; and the free-current integrand |f E|. Every part is a maximum over nodes and channels, so chunks combine by taking the largest of each, and each ratio is normalised by the most-ionised node rather than node by node: a node the field barely ionises, whose rate is one sampled spike, would otherwise report the filter’s full undershoot on a rate at the level of rounding (plasma_integration_statistics() forms the ratios).

Parameters:
field_tnumpy.ndarray

Real node fields in V/m, shape (n_nodes_chunk, n_pol, n_to).

channelssequence of PlasmaCoefficients

The channels’ compiled coefficients, each with its rate table.

preionisationtuple of float

Each channel’s initial fraction at this plane.

dtfloat

Fine-grid sample spacing in s.

integralsolprop.nonlinear.ionisation.RunningIntegral or None, optional

The running-integral rule; None is the trapezoid.

Returns:
numpy.ndarray

The largest exposure dip -min dS and rise max dS; the largest |paid - freed| and freed; the largest |f E| in the first 2K samples and anywhere. Dimensionless but for the last two, in V/m.

solprop.backends.numpy_ref.passes.plasma_integration_statistics(parts)[source]#

Form the three diagnostics from plasma_integration_diagnostics()’s parts.

Parameters:
partsnumpy.ndarray

The six maxima, combined over every chunk of nodes.

Returns:
tuple of float

The undershoot max(-dS) / max(dS), the work defect max |paid - freed| / max freed and the boundary level, each zero where its denominator is.

solprop.backends.numpy_ref.passes.plasma_scan(field_t, rate, coefficients, *, preionisation, dt, integral=None)[source]#

Integrate the ionised fraction and the free-electron current along t.

The plasma_scan stage: two of the three cumulative integrals the response declares, which is why its domain is memory – the value at one sample depends on every earlier sample and on none of the later ones.

f = 1 - (1 - f_0) exp(-int W dt), with the pre-ionised fraction inside the exponential. Writing f_0 + 1 - exp(-int W dt) instead would let f exceed one, make 1 - f negative and stop df/dt being W (1 - f), so the loss term would no longer pay exactly one ionisation potential per electron. Then J_free = (e^2/m_e) int f E dt, the current of electrons already freed and still being accelerated: what lowers the refractive index and blue-shifts the pulse.

Both integrals are the problem’s running-integral rule (SEP 18), applied in the argument order of solprop.nonlinear.responses.PlasmaResponse.evaluate(), each a vectorised cumulative sum (scipy.integrate.cumulative_trapezoid); this pass and the oracle are the same computation, and other backends’ scans are held to the scan’s rounding bound (SEP 25), not to its order. Under the band rule the exposure is floored at zero, so the fraction stays in [preionisation, 1) and is non-decreasing to within the filter’s undershoot. Here the multipass form is the reference, and fusing the three sums into one sequential sweep belongs to the compiled backends (docs/guide/media/nonlinear.md).

Parameters:
field_tnumpy.ndarray

The real field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m.

ratenumpy.ndarray

The ionisation rate from ionisation_rate(), s^-1, shape (n_nodes_chunk, n_to).

coefficientssolprop.translate.responses.PlasmaCoefficients

The compiled plasma group; only charge_mass_ratio is read here.

preionisationfloat

The pre-ionised fraction at this position, coefficients.preionisation(z), in [0, 1].

dtfloat

Sample spacing of the fine time grid, s, positive.

integralsolprop.nonlinear.ionisation.RunningIntegral or None, optional

The running-integral rule, responses.plasma_integral; None is the trapezoid.

Returns:
fractionnumpy.ndarray

The fraction ionised, shape (n_nodes_chunk, n_to), in [preionisation, 1] and non-decreasing in time (under the band rule, to within the filter’s undershoot).

free_currentnumpy.ndarray

J_free per unit density, shape (n_nodes_chunk, n_pol, n_to), in C/(m^2 s) per m^-3.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import plasma_scan
>>> from solprop.nonlinear.optmat_rates import OptmatRate
>>> from solprop.nonlinear.responses import PlasmaResponse
>>> from solprop.translate.responses import plasma_coefficients
>>> compiled = plasma_coefficients(
...     PlasmaResponse("Ar", OptmatRate("Ar", 800e-9)), (1e9, 1e11), z_max=1.0
... )
>>> field = np.zeros((1, 1, 4))
>>> fraction, current = plasma_scan(
...     field, np.zeros((1, 4)), compiled, preionisation=0.25, dt=1e-16
... )
>>> fraction.tolist(), bool(np.all(current == 0.0))
([[0.25, 0.25, 0.25, 0.25]], True)
solprop.backends.numpy_ref.passes.project(polarisation_t, projection, taper, n_pol)[source]#

Project a chunk’s polarisation back onto the modes.

Pass 5 of architecture.md section 8: P[m, t] = sum_{i, p} W[m, i n_pol + p] P[i, p, t] / tau(z), one GEMM with the quadrature weights already folded into projection. The caller accumulates over chunks, which is what makes the space-domain buffers bounded by the chunk rather than by the node count (architecture.md section 13).

Parameters:
polarisation_tnumpy.ndarray

The polarisation at the chunk’s nodes, shape (n_nodes_chunk, n_pol_node, n_to), in C/m^2.

projectionnumpy.ndarray

The chunk’s columns of the projection matrix, float64, shape (n_modes, n_nodes_chunk * n_pol), weights folded in.

taperfloat

The rescaling at this position, transform.taper_scale(z), or 1.0. The projection is divided by it, as the synthesis is multiplied by it.

n_polint

Transverse components the transform resolves at each node, 1 or 2.

Returns:
numpy.ndarray

The modal polarisation on the fine time grid, shape (n_modes, n_pol, n_to).

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import project
>>> nodes = np.ones((2, 1, 4))
>>> project(nodes, np.array([[1.0, 3.0]]), 1.0, 1).shape
(1, 1, 4)
solprop.backends.numpy_ref.passes.raman_polarisation(field_t, spectrum, *, kind, thg, density)[source]#

Evaluate the delayed Raman polarisation at every node of a chunk.

The scalar reduction preserves the formula of solprop.nonlinear.raman.RamanResponse.evaluate() in its original order of operations:

\[P(t) = \varepsilon_0\,\rho\,E(t) \int h(t - t')\,D(t')\,\mathrm{d}t' ,\]

with D the real drive solprop.nonlinear.raman.raman_drive() chooses for the grid – E^2 on a field grid that keeps the harmonic, |E_a|^2/2 on one that does not, |A|^2/2 on an envelope grid.

The convolution runs on a doubled, zero-padded grid, because a transform product is a circular convolution and a rotational coherence outlives the window by orders of magnitude: hydrogen’s lives about 2.9 ns at one bar against a window of tens of picoseconds, so the wrapped tail would be almost the whole response. The drive and the kernel occupy the first half of an array of twice the length, which makes the product the linear convolution exactly (docs/guide/media/raman-numerics.md, “The doubled grid”).

For two Cartesian components, convolve the trace, difference and twice the off-diagonal entry of the real quadratic tensor, then apply the resulting symmetric matrix to the field. Isotropic-only declarations need just the trace history. Each history is real on both grids and uses one real transform pair. numpy.fft preserves the scalar oracle’s transform arithmetic.

Parameters:
field_tnumpy.ndarray

Field at the chunk’s nodes, shape (n_nodes_chunk, n_pol, n_to), in V/m: real on a field grid, complex on an envelope grid. One or two physical Cartesian components.

spectrumnumpy.ndarray

Damped complex128 kernels with dt included. Scalar fields use shape (n_to + 1,) from raman_spectrum(); vector fields use shape (n_raman_kernels, n_to + 1) from tensor_raman_spectrum(), with one trace kernel or the trace/traceless pair.

kind{‘field’, ‘envelope’}

The spectral layout the samples belong to.

thgbool

Whether a field grid’s drive keeps the 2 omega beat, and so radiates a third harmonic through the later multiplication by the field. Ignored on an envelope grid, which cannot carry it in a quadratic drive.

densityfloat

The density path at this position, responses.density(z), in m^-3 for a gas and 1.0 for a medium whose kernel is already absolute. Non-negative.

Returns:
numpy.ndarray

Raman polarisation in C/m^2, shape (n_nodes_chunk, n_pol, n_to), real on a field grid and complex on an envelope grid.

Raises:
BackendError

If the component count or effective kernel shape is incompatible with the doubled temporal grid.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import raman_polarisation
>>> field = np.ones((1, 1, 4))
>>> flat = np.zeros(5, dtype=complex)
>>> polarisation = raman_polarisation(
...     field, flat, kind="field", thg=True, density=1.0
... )
>>> bool(np.all(polarisation == 0.0))
True
solprop.backends.numpy_ref.passes.raman_spectrum(undamped, rates, fraction, dt)[source]#

Damp each undamped Raman row at its own rate, sum, scale and transform.

The half of the kernel that moves along z. raman_kernel_undamped holds one causal, tapered row per independently damped term and raman_decay_rates one scalar path per row (kernel-contract.md section 5); this forms f_R sum_k u_k(t) exp(-gamma_k t) on the doubled grid and hands it to solprop.translate.raman.kernel_spectrum(), which folds in the sample spacing the discrete sum needs to become an integral. There is one statement of that normalisation and this is not it.

Damping each term separately is the whole reason the contract carries rows rather than one array: a molecular gas’s rotational and vibrational coherences decay 3.5 times (H2) to 17 times (N2) apart at one bar, so damping their sum once would be wrong for every such gas (docs/guide/media/raman-numerics.md).

The terms are summed in the order of the rows and the fraction multiplies the sum, which is the arithmetic of solprop.translate.raman.sample_kernel() operation for operation, so the two spectra agree to the last bit. The padding half of each row is zero, so damping it costs nothing and changes nothing.

Parameters:
undampednumpy.ndarray

responses.raman_kernel_undamped, float64, shape (n_terms, 2 n_to), in m^5 V^-2 s^-1.

ratesarray_like

The decay rate of each term at this position, s^-1, non-negative, one per row: [table(z) for table in responses.raman_decay_rates].

fractionfloat

responses.raman_fraction, in (0, 1]; 1 for an absolute molecular response.

dtfloat

Sample spacing of the fine time grid, s, positive.

Returns:
numpy.ndarray

The damped kernel spectrum with the dt scaling folded in, complex128, shape (n_to + 1,), in the internal convention.

Raises:
BackendError

If there is not one rate per row.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import raman_spectrum
>>> rows = np.zeros((1, 8))
>>> rows[0, 1] = 1.0
>>> spectrum = raman_spectrum(rows, [0.0], 1.0, 2.0)
>>> spectrum.shape, bool(np.isclose(spectrum[0], 2.0))
((5,), True)
solprop.backends.numpy_ref.passes.stage_input(state, stage_sum, dz, propagator, noise)[source]#

Combine the stages, return to the lab frame and add the noise.

u = exp(dPhi_i) (y + h sum_j a_ij k_j) + N(z_n + c_i h), the first pass of architecture.md section 8. The interaction picture is anchored at the start of the step, so exp(dPhi_i) carries the stage input forward to where the nonlinear term is evaluated, and the shot-noise field is added there, in the frequency domain, where it costs one multiply-add and no transform.

Parameters:
statenumpy.ndarray

The state at the anchor, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts.

stage_sumnumpy.ndarray

sum_j a_ij k_j in the interaction variable, same shape and dtype, in the square root of watts per metre.

dzfloat

Length of the step, m, positive.

propagatornumpy.ndarray

exp(dPhi_i), dimensionless, broadcastable against the state; the batch axis is 1 when the members share a fibre.

noisenumpy.ndarray or None

The shot-noise field already advanced to this position, N0 exp(i Im Phi(omega, z)), in the square root of watts, or None when no noise model is active.

Returns:
numpy.ndarray

The lab-frame field the nonlinear term is evaluated at, same shape and dtype.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import stage_input
>>> state = np.ones((1, 1, 1, 2), dtype=complex)
>>> zero = np.zeros_like(state)
>>> stage_input(state, zero, 0.1, 2.0 * np.ones_like(state), None).ravel()
array([2.+0.j, 2.+0.j])
solprop.backends.numpy_ref.passes.stage_output(polarisation, prefactor, omega_window, inverse_propagator, absorber=None)[source]#

Normalise the polarisation and return it to the interaction picture.

k_i = exp(-dPhi_i) (p(omega, z) W(omega) crop(Ptilde) + A_i), pass 7 of architecture.md section 8, with A_i the temporal absorber’s own term when the "nonlinear" absorber scheme is active and nothing otherwise. The absorber sits beside the normalised polarisation rather than inside it because it is a rate on the field and not a polarisation: the prefactor does not apply to it. The prefactor is the UPPE factor -i omega / (2 eps0 c Re n_eff) (or its envelope form) with the crop scale folded in, and it is exactly zero outside the active band, so no unnormalised polarisation survives there for the stages and the error estimate to see.

Parameters:
polarisationnumpy.ndarray

The modal polarisation on the coarse grid, complex128, shape (n_modes, n_pol, n_omega), in C/m^2.

prefactornumpy.ndarray

responses.prefactor_at(responses, z), complex128, in m/(V s), broadcastable against the polarisation.

omega_windownumpy.ndarray

The spectral collar profile in [0, 1], float64, shape (n_omega,).

inverse_propagatornumpy.ndarray

exp(-dPhi_i), dimensionless, broadcastable against the state; the negation of the exponent stage_input() used, never an independent evaluation.

absorbernumpy.ndarray or None, optional

The temporal absorber’s contribution from absorber_term(), or None under the "split_step" scheme, where the rate is applied to the state once per accepted step instead.

Returns:
numpy.ndarray

The stage derivative in the interaction variable, complex128, shape (n_batch, n_modes, n_pol, n_omega), in the square root of watts per metre.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import stage_output
>>> polarisation = np.ones((1, 1, 2), dtype=complex)
>>> prefactor = np.array([[[-2.0j, 0.0j]]])
>>> unit = np.ones((1, 1, 1, 2), dtype=complex)
>>> stage_output(polarisation, prefactor, np.ones(2), unit).ravel().tolist()
[-2j, 0j]
solprop.backends.numpy_ref.passes.synthesise(modal_t, synthesis, taper, n_pol)[source]#

Synthesise the physical field on a chunk of quadrature nodes.

Pass 3 of architecture.md section 8: E[i, p, t] = tau(z) sum_m S[i n_pol + p, m] a[m, t], one GEMM over the mode axis. synthesis is the chunk’s rows of transform.synthesis, which already carry the oversampling pad scale, and taper is transform.taper_scale(z), the exact rescaling a scale-invariant waveguide needs (solprop.modes.quadrature.ModeMatrices.rescaled()).

The transverse components are carried by the rows of synthesis or by the state’s polarisation axis and never by both (kernel-contract.md section 5), so the returned polarisation axis has n_pol times the state’s length and one of the two factors is one.

Parameters:
modal_tnumpy.ndarray

The modal field on the fine time grid from to_time(), shape (n_modes, n_pol, n_to), in the square root of watts.

synthesisnumpy.ndarray

The chunk’s rows of the synthesis matrix, float64, shape (n_nodes_chunk * n_pol, n_modes), in V/(m sqrt(W)), already carrying the pad scale.

taperfloat

The rescaling at this position, transform.taper_scale(z), or 1.0 for a waveguide whose cross-section does not move.

n_polint

Transverse components the transform resolves at each node, 1 or 2.

Returns:
numpy.ndarray

The field at the chunk’s nodes, shape (n_nodes_chunk, n_pol_node, n_to), in V/m: real on a field grid, complex on an envelope grid.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import synthesise
>>> modal = np.ones((1, 1, 4))
>>> synthesise(modal, np.array([[2.0], [3.0]]), 1.0, 1).shape
(2, 1, 4)
solprop.backends.numpy_ref.passes.tensor_raman_spectrum(components, rates, densities, dt)[source]#

Form effective trace and traceless spectra from independently damped terms.

Parameters:
componentstuple of RamanCoefficients

Tensor-capable components in mixture order, with per-row numeric weights.

ratestuple of float

Decay rates in s^-1, flattened by component then kernel row.

densitiestuple of float

Amplitude densities in m^-3. Supply (1.0,) when a single component’s density is multiplied outside the convolution.

dtfloat

Fine-grid sample spacing in s.

Returns:
numpy.ndarray

Shape (n_raman_kernels, n_to + 1) with dt included. One trace kernel for isotropic-only declarations; otherwise trace and traceless rows.

Raises:
BackendError

If counts disagree or any term lacks a vector contraction.

solprop.backends.numpy_ref.passes.to_frequency(polarisation_t, crop_index, *, kind)[source]#

Transform the modal polarisation to frequency and crop to the coarse grid.

Pass 6 of architecture.md section 8: one batched forward transform of all n_modes * n_pol columns. The crop takes the coarse bins out of the fine spectrum at crop_index with no scale factor: the crop scale n_t/n_to is inside responses.prefactor (kernel-contract.md section 5), which stage_output() applies next.

Parameters:
polarisation_tnumpy.ndarray

The modal polarisation on the fine time grid, accumulated over the chunks by project(), shape (n_modes, n_pol, n_to), in C/m^2.

crop_indexnumpy.ndarray

Where each coarse bin sits in the fine spectrum, numpy.intp, shape (n_omega,).

kind{‘field’, ‘envelope’}

Spectral layout.

Returns:
numpy.ndarray

The modal polarisation on the coarse grid, complex128, shape (n_modes, n_pol, n_omega).

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import to_frequency
>>> samples = np.ones((1, 1, 8))
>>> index = np.arange(3, dtype=np.intp)
>>> to_frequency(samples, index, kind="field").real
array([[[8., 0., 0.]]])
solprop.backends.numpy_ref.passes.to_time(spectrum, crop_index, n_to, *, kind)[source]#

Pad onto the fine grid and transform every modal column to time.

Pass 2 of architecture.md section 8: one batched inverse transform of all n_modes * n_pol columns, whatever the quadrature does next. The coarse bins are scattered into the fine spectrum at crop_index with no scale factor: the pad scale n_to/n_t is inside transform.synthesis (kernel-contract.md section 4), so that every response sees the field in volts per metre whatever its degree in the field.

The mode and polarisation axes are folded together for the transform, so the columns a scalar problem transforms are exactly the ones it transformed before the polarisation axis existed.

Parameters:
spectrumnumpy.ndarray

Coarse modal spectrum, complex128, shape (n_modes, n_pol, n_omega), in the square root of watts.

crop_indexnumpy.ndarray

Where each coarse bin sits in the fine spectrum, numpy.intp, shape (n_omega,).

n_toint

Number of fine time samples, positive.

kind{‘field’, ‘envelope’}

Spectral layout, which selects the complex-to-real or complex-to-complex transform.

Returns:
numpy.ndarray

The modal field on the fine time grid, shape (n_modes, n_pol, n_to), in the square root of watts: real on a field grid, complex on an envelope grid.

Examples

>>> import numpy as np
>>> from solprop.backends.numpy_ref.passes import to_time
>>> spectrum = np.array([[[1.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j]]])
>>> index = np.arange(3, dtype=np.intp)
>>> to_time(spectrum, index, 8, kind="field")
array([[[0.125, 0.125, 0.125, 0.125, 0.125, 0.125, 0.125, 0.125]]])