solprop.translate.coupled#

The metric and the coupled-flow representation of a general modal guide.

Two leaves of section 3 of docs/dev/kernel-contract.md, added at kernel schema 12: ModalMetricTables, the explicit physical flux metric the controller, the linear flow and the output layer all read, and CoupledFlowTables, the frozen or ordered representation of a generator whose channels couple.

CoupledFlowTables is the record, its validation and its printing, built by constant_schur_flow() and frozen_residual_flow() for the routes whose channels couple – a capillary’s coherent wall and taper coupling, and a general modal guide whose declared basis or motion couples its channels. Everywhere else it is None and a guide propagates through the diagonal cumulative flow. ModalMetricTables admits all three kinds: "identity", "diagonal" and the "dense" kind a declared unitary change of basis produces (dense_metric()). Which of them a backend executes is its own admission, and numpy_ref refuses the diagonal one, which nothing here builds.

Units are SI throughout (docs/dev/conventions.md): omega in rad/s, z_edges in m, the generator and every term of it in 1/m, metric entries dimensionless, errors and bounds dimensionless. Every complex array is in the internal FFT convention, like the rest of the contract; there is no fourth conversion site.

Module Attributes

`MetricKind`	Which representation of the physical flux metric a problem carries.
`SectionMethod`	How one interval of a coupled flow is integrated.
`SECTION_TOL`	Relative coefficient tail a section's z expansion is resolved to.
`MAX_SECTION_DEGREE`	Highest Chebyshev degree along z tried on one section before it is refused.

Functions

`advisory_step`(flow[, scored])	Return the advisory step limit of every section and frequency, in m.
`constant_schur_flow`(omega, z_edges, ...[, ...])	Build the coupled flow of one constant section, with its certificate.
`coupled_flow_factors`(flow)	Freeze and decompose each section of a dense coupled flow on its active bins.
`dense_metric`(omega, z_edges, sampled, *, ...)	Build the dense metric of a basis from its own sampled flux Gram.
`eigenbasis_limit`(n_modes)	Return the largest eigenbasis condition the compiled flow actions accept.
`flow_edge_slack`(flow)	Return the tolerance a position is matched to a section edge within.
`flow_section`(flow, z)	Return the index of the section a position lies in.
`flow_weights`(flow, z[, section])	Return the weights that combine a coupled flow's coefficients at one position.
`frozen_generator`(flow, index)	Return the generator one section of a coupled flow is frozen at.
`frozen_residual_flow`(omega, z_edges, ...)	Build the certified coupled flow of a varying guide, frozen section by section.
`interval_count`(edges)	Return how many z intervals a set of stored positions declares.
`metric_at`(metric, z)	Evaluate a basis metric at one statistics or observation plane.
`orthonormal_metric`(omega, z_edges, *, ...)	Build the identity metric of a basis certified orthonormal in true power.
`residual_truncation`(x)	Return the frozen-residual truncation defect of one step, dimensionless.
`schur_blocks`(generator)	Return the complex Schur factors of a generator, one per frequency.
`schur_eigenbasis`(vectors, triangular)	Return the eigenbasis of each Schur factor, and where it is well conditioned.
`section_expansion`(sample, z_edges, *[, tol, ...])	Expand an array-valued function of z in Chebyshev coefficients per section.
`taper_identity`(flow, *, channels, active[, ...])	Return the modal identity of a capillary taper's coupling.
`wall_identity`(flow, limits, *, channels[, ...])	Return the modal identity of a capillary carrying a coherent wall.

Classes

`CoupledFlowFactors`(active, vectors, ...)	A coupled flow's frozen blocks and Schur factors on its active bins.
`CoupledFlowTables`(omega, z_edges, generator, ...)	A generator whose channels couple, frozen or ordered along z.
`ModalMetricTables`(kind, omega, z_edges, ...)	The explicit physical flux metric of a modal basis.

class solprop.translate.coupled.CoupledFlowFactors(active, vectors, triangular, frozen, condition_bound, residual_terms)[source]#

Bases: object

A coupled flow’s frozen blocks and Schur factors on its active bins.

What both backends execute a flow from, derived once from the leaf (coupled_flow_factors()); a decomposition is not a contract row, so it is formed from the stored generator rather than carried across the seam.

Attributes:

activenumpy.ndarray: Bool, (n_omega,): the bins the generator is non-zero at; the rest of the band has an identity action.
vectors, triangularnumpy.ndarray: Complex128, (n_intervals, n_active, n_modes, n_modes): each section’s Schur factors Z and T per active bin (schur_blocks()).
frozennumpy.ndarray: Complex128, same shape: each section’s frozen block (frozen_generator()), in 1/m.
condition_boundnumpy.ndarray: Float64, (n_intervals, n_active): the declared conditioning row.
residual_termsnumpy.ndarray or None: Complex128, (n_intervals * n_deg, n_active, n_modes, n_modes): a frozen-residual flow’s generator less its connection, flattened over (n_intervals, n_deg); None for a constant-Schur flow.

active: Bool[ndarray, 'n_omega']#

condition_bound: Float[ndarray, 'n_intervals n_active']#

frozen: Complex[ndarray, 'n_intervals n_active n_modes n_modes']#

residual_terms: Complex[ndarray, 'n_terms n_active n_modes n_modes'] | None#

triangular: Complex[ndarray, 'n_intervals n_active n_modes n_modes']#

vectors: Complex[ndarray, 'n_intervals n_active n_modes n_modes']#

class solprop.translate.coupled.CoupledFlowTables(omega, z_edges, generator, fixed_terms, term_weights, wall_vectors, wall_weights, connection, generator_error, residual_norm_bound, inverse_condition_bound, section_methods, gauge_id, certificate_id)[source]#

Bases: Displayed

A generator whose channels couple, frozen or ordered along z.

The leaf of contract section 3 that carries the coherent wall operator and the taper coupling. numpy_ref and native execute two combinations. The first is the dense generator of a single constant section, section_methods == ("constant_schur",), no connection and an exactly zero residual, built by constant_schur_flow() from a constant coherent block. The second is the dense generator of any number of "frozen_residual" sections at any degree along z, with the taper connection and an explicit residual the Runge-Kutta source carries, built by frozen_residual_flow() from the same block evaluated along the declared profile. Every other method, representation and mixture is refused by name by the problem’s validation (solprop.translate.admission.check_general_route()).

The dense generator is optional. The factored alternative combines fixed matrices with scalar weights, rank-one wall terms and connection coefficients without a mandatory dense (z, omega, modes, modes) allocation; exactly one of the two complete representations is stored, never an ambiguous mixture.

Attributes:

omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, shape (n_omega,), increasing and non-negative. A problem’s leaf carries grids.omega itself, whose first bin is zero on a field grid (solprop.translate.problem.PropagationProblem.validate()).
z_edgesnumpy.ndarray: Interval edges in m, shape (n_edges,), increasing.
generatornumpy.ndarray or None: Dense generator in 1/m, shape (n_intervals, n_deg, n_modes, n_modes, n_omega); the dense representation, mutually exclusive with the factored one.
fixed_termsnumpy.ndarray or None: Fixed matrices in 1/m, shape (n_terms, n_modes, n_modes, n_omega).
term_weightsnumpy.ndarray or None: Dimensionless scalar weights of those matrices, shape (n_intervals, n_deg, n_terms); supplied with fixed_terms.
wall_vectorsnumpy.ndarray or None: Rank-one wall vectors, dimensionless, shape (n_intervals, n_modes, n_omega); supplied with wall_weights.
wall_weightsnumpy.ndarray or None: Their scalar weights in 1/m, shape (n_intervals, n_deg, n_omega), so that w (v x v) carries the generator’s own 1/m, and constant_schur_flow() checks the reconstruction.
connectionnumpy.ndarray or None: Taper connection coefficients in 1/m, shape (n_intervals, n_deg, n_modes, n_modes, n_omega).
generator_errornumpy.ndarray: Non-negative dimensionless representation error, shape (n_intervals, n_omega).
residual_norm_boundnumpy.ndarray: Non-negative bound in 1/m on the explicit residual, same shape.
inverse_condition_boundnumpy.ndarray: Dimensionless bound of at least one on the conditioning of the inverse action, same shape.
section_methodstuple of str: One of "diagonal", "constant_schur" or "frozen_residual" per interval.
gauge_id, certificate_idstr: Non-empty identities; scalars, never arrays, enter resume identity.

Raises:

TranslationError: Naming the field that failed.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import CoupledFlowTables
>>> flow = CoupledFlowTables(
...     np.linspace(1e15, 2e15, 2), np.array([0.0, 1.0]),
...     np.zeros((1, 1, 2, 2, 2), dtype=complex), None, None, None, None, None,
...     np.zeros((1, 2)), np.zeros((1, 2)), np.ones((1, 2)), ("diagonal",),
...     "lab", "none",
... )
>>> flow.representation, flow.n_modes
('dense', 2)

certificate_id: str#

connection: Complex[ndarray, 'n_intervals_connection n_deg_connection n_modes_connection_out n_modes_connection_in n_omega_connection'] | None#

fixed_terms: Complex[ndarray, 'n_terms_fixed n_modes_fixed_out n_modes_fixed_in n_omega_fixed'] | None#

gauge_id: str#

generator: Complex[ndarray, 'n_intervals_generator n_deg_generator n_modes_generator_out n_modes_generator_in n_omega_generator'] | None#

generator_error: Float[ndarray, 'n_intervals_generator_error n_omega_generator_error']#

inverse_condition_bound: Float[ndarray, 'n_intervals_inverse_bound n_omega_inverse_bound']#

property n_deg: int[source]#

How many Chebyshev coefficients along z each section carries.

The degree axis of the dense generator, or one for the factored representation, whose weights carry the z dependence instead.

property n_modes: int[source]#: How many channels the flow couples.

omega: Float[ndarray, 'n_omega_flow']#

property representation: Literal['dense', 'factored'][source]#

Which of the two complete representations is stored.

"dense" when the generator array is present, "factored" when the fixed matrices and their weights are.

residual_norm_bound: Float[ndarray, 'n_intervals_residual_bound n_omega_residual_bound']#

section_methods: tuple[Literal['diagonal', 'constant_schur', 'frozen_residual'], ...]#

term_weights: Complex[ndarray, 'n_intervals_term_weights n_deg_term_weights n_terms_term_weights'] | None#

wall_vectors: Complex[ndarray, 'n_intervals_wall_vectors n_modes_wall_vectors n_omega_wall_vectors'] | None#

wall_weights: Complex[ndarray, 'n_intervals_wall_weights n_deg_wall_weights n_omega_wall_weights'] | None#

z_edges: Float[ndarray, 'n_edges_flow']#

solprop.translate.coupled.MAX_SECTION_DEGREE: int = 32#

Highest Chebyshev degree along z tried on one section before it is refused.

A section of a declared taper or pressure stage carries a smooth rational function of a(z) and rho(z); a degree this high on one section means the section spans structure the caller should have declared or should divide, which is what CoherentWall.sections is for. The search starts at solprop.translate.chebyshev.INITIAL_DEGREE and doubles.

solprop.translate.coupled.MetricKind#

Which representation of the physical flux metric a problem carries.

"identity" is an orthonormal power basis, whose metric owns no coefficient array at all; "diagonal" is a per-channel positive weight; "dense" is the full Hermitian matrix a nonorthogonal basis needs, which a declared unitary change of basis produces.

alias of Literal[‘identity’, ‘diagonal’, ‘dense’]

class solprop.translate.coupled.ModalMetricTables(kind, omega, z_edges, diagonal, dense, error, minimum_eigenvalue, maximum_eigenvalue, gauge_id, certificate_id)[source]#

Bases: Displayed

The explicit physical flux metric of a modal basis.

The documented leaf of contract section 3 that the step controller’s norm, the linear flow and every physical-energy statistic read, so that one W is shared by identity and never rebuilt. A coordinate change c = T d takes it to T.H @ W @ T; interference terms cannot be assigned unique positive per-channel energies, so coordinate intensities are not additive in a nonorthogonal basis.

Attributes:

kind{‘identity’, ‘diagonal’, ‘dense’}: Which representation is stored. "identity" owns no coefficient array; "diagonal" owns diagonal; "dense" owns dense, is Hermitian at every interval, degree and frequency, and is what a declared unitary change of basis produces (dense_metric()). Which of them a backend executes is its own admission: numpy_ref takes the identity and a dense metric that is a unitary transform of one, and refuses a diagonal metric, which nothing builds.
omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, shape (n_omega,), increasing and non-negative. A problem’s leaf carries grids.omega itself, whose first bin is zero on a field grid (solprop.translate.problem.PropagationProblem.validate()).
z_edgesnumpy.ndarray: Stored positions in m, shape (n_edges,), increasing. The interval-indexed arrays are counted by interval_count().
diagonalnumpy.ndarray or None: Positive dimensionless per-channel weights, shape (n_intervals, n_deg, n_modes, n_omega); present exactly for kind="diagonal".
densenumpy.ndarray or None: Hermitian dimensionless metric, shape (n_intervals, n_deg, n_modes, n_modes, n_omega); present exactly for kind="dense".
errornumpy.ndarray: Non-negative dimensionless bound on the stored representation, shape (n_intervals, n_omega). Declared by the caller; nothing here infers it.
minimum_eigenvalue, maximum_eigenvaluenumpy.ndarray: Non-negative dimensionless extreme eigenvalues of W, shape (n_intervals, n_omega), with the minimum never above the maximum. A minimum of zero is a semidefinite metric: it may be stored as flux data but certifies no positive error or noise norm.
gauge_id, certificate_idstr: Non-empty phase/coordinate identity and the identity of the certificate the bounds came from. Scalar identities, never arrays, enter resume identity.

Raises:

TranslationError: Naming the field that failed.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import ModalMetricTables
>>> zeros = np.zeros((1, 3))
>>> metric = ModalMetricTables(
...     "identity", np.linspace(1e15, 2e15, 3), np.array([0.0]), None, None,
...     zeros, np.ones((1, 3)), np.ones((1, 3)), "lab", "orthonormal",
... )
>>> metric.kind, metric.n_omega
('identity', 3)

certificate_id: str#

check_rows()[source]#

Check the certificate rows against the record’s own axes, by name.

Construction checks them; a leaf replaced past it is checked again by solprop.translate.problem.PropagationProblem.validate() before anything reads a row.

Raises:

TranslationError: If error, minimum_eigenvalue or maximum_eigenvalue is not one row per interval of z_edges (interval_count()) by one column per coarse frequency.

dense: Complex[ndarray, 'n_intervals_dense n_deg_dense n_modes_dense_out n_modes_dense_in n_omega_dense'] | None#

diagonal: Float[ndarray, 'n_intervals_diagonal n_deg_diagonal n_modes_diagonal n_omega_diagonal'] | None#

error: Float[ndarray, 'n_intervals_metric_error n_omega_metric_error']#

gauge_id: str#

kind: Literal['identity', 'diagonal', 'dense']#

maximum_eigenvalue: Float[ndarray, 'n_intervals_metric_max n_omega_metric_max']#

minimum_eigenvalue: Float[ndarray, 'n_intervals_metric_min n_omega_metric_min']#

property n_modes: int[source]#

How many channels the metric describes.

The channel count of whichever coefficient array is stored; 0 for "identity", which owns none and takes its size from the state.

property n_omega: int[source]#: How many coarse frequencies the metric is tabulated on.

omega: Float[ndarray, 'n_omega_metric']#

z_edges: Float[ndarray, 'n_edges_metric']#

solprop.translate.coupled.SECTION_TOL: float = 1e-12#

Relative coefficient tail a section’s z expansion is resolved to.

The same target solprop.translate.chebyshev.DEFAULT_TOL sets for the cumulative operator, and for the same reason: below it the interpolant is resolving the round-off of its own samples. What is not resolved to it is refused rather than accepted with a larger error, because the residual bound this expansion feeds is the whole certificate of the frozen split.

The samples’ own floor is well below it. The coherent wall’s diagonal k - beta_11 is formed from the two refractivities ((n_core - 1) - (n_eff - 1)), not as the difference of two wavenumbers that agree to five or six digits, which carried the absolute rounding of k (about 5e-9 1/m, up to 5e-10 of k - beta_11 itself) and now carries 1e-14 of it at worst. Forced to degree 64 across the transition of the three- and six-channel qualification tapers and the TE/TM taper, the generator’s tail floors at 4e-14 relative (8e-15 to 4.0e-14), 25 times below this target, so a tail above it is unresolved structure and never round-off.

solprop.translate.coupled.SectionMethod#

How one interval of a coupled flow is integrated.

"diagonal" is the existing cumulative commuting flow; "constant_schur" is the Schur matrix-function action of a constant coherent block; "frozen_residual" freezes the block inside the section and leaves its variation to the explicit source.

alias of Literal[‘diagonal’, ‘constant_schur’, ‘frozen_residual’]

solprop.translate.coupled.advisory_step(flow, scored=None)[source]#

Return the advisory step limit of every section and frequency, in m.

h_advisory = 0.9 / (K ||U_s|| ||U_s^-1||), with K the interval upper bound on the norm of the whole explicit linear residual – which is CoupledFlowTables.residual_norm_bound – and the conditioning factor its interaction-picture similarity transform costs, which is CoupledFlowTables.inverse_condition_bound. Both rows are already across the seam, so this needs no row of its own.

The factor 0.9 is a warning margin below the production tableau’s imaginary stability radius 0.997189, not an approved solver default and not a stability proof: a wall difference can have positive real eigenvalues, representing local growth relative to the frozen attenuating frame, and no positive step makes that residual alone contractive. In a lossy or nonnormal case this is advisory only and a growth certificate is still required (solprop.translate.problem.validate_coupled_flow() says so in the warning).

A section with an exactly zero residual – every constant one – returns infinity, because there is no explicit term to resolve. So does a bin outside scored: the explicit residual is applied only on the bins the step controller’s error norm scores, and the spectral collar beyond them is carried by the frozen flow alone, which has no stability limit.

Parameters:

flowCoupledFlowTables: The leaf.
scorednumpy.ndarray of bool, optional: The bins the run applies the residual on, shape (n_omega,): the error norm’s grids.inband_mask. None takes every bin.

Returns:

numpy.ndarray: Float64, shape (n_intervals, n_omega), positive, inf where the residual is exactly zero.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import advisory_step, constant_schur_flow
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 0] = block[1, 1] = -1.0
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> bool(np.all(np.isinf(advisory_step(flow))))
True

solprop.translate.coupled.constant_schur_flow(omega, z_edges, generator, *, dz_max, condition_limit, wall_vectors=None, wall_weights=None, gauge_id, certificate_id)[source]#

Build the coupled flow of one constant section, with its certificate.

The coupled flow at constant geometry: one section whose complete generator is frozen – exactly, because it does not vary – and integrated by a Schur matrix-function action, with the explicit residual identically zero. The dense generator is stored, and the rank-one wall pair beside it when the caller has one, so that a reader can see which part of the block is the wall without factoring it back out.

The generator is in the internal FFT convention, like every other array of the contract, and carries only what the diagonal cumulative flow does not: the shared scalar carrier of the section stays in linear.phi, where it is integrated analytically and commutes with the whole block.

Parameters:

omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, increasing.
z_edgesnumpy.ndarray: The section’s two edges in m, covering the run.
generatornumpy.ndarray: The constant internal-convention generator in 1/m, complex, shape (n_modes, n_modes, n_omega), exactly zero on every bin the model was not evaluated at.
dz_maxfloat: Largest step the run may take, m, positive: what the conditioning bound is declared over.
condition_limitfloat: Largest admitted conditioning, dimensionless and at least one.
wall_vectorsnumpy.ndarray or None, optional: The rank-one wall vector h, dimensionless, shape (n_modes,).
wall_weightsnumpy.ndarray or None, optional: Its weight in 1/m, complex, shape (n_omega,), supplied with the vector. The pair multiplies out to the generator’s own 1/m: the vector carries no dimension and the weight carries the whole of it.
gauge_idstr: The basis’s phase and coordinate identity.
certificate_idstr: Identity of the certificate the bounds came from.

Returns:

CoupledFlowTables: The validated leaf, section_methods == ("constant_schur",).

Raises:

TranslationError: If the generator varies in shape from the axes, if the wall pair does not reconstruct the generator’s off-diagonal, or if the conditioning bound exceeds condition_limit at any frequency.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import constant_schur_flow
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 0] = block[1, 1] = -1.0
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> flow.section_methods, float(flow.residual_norm_bound.max())
(('constant_schur',), 0.0)

solprop.translate.coupled.coupled_flow_factors(flow)[source]#

Freeze and decompose each section of a dense coupled flow on its active bins.

Parameters:

flowCoupledFlowTables: A leaf of the dense representation, as the problem’s validation admits it.

Returns:

CoupledFlowFactors: The frozen blocks, their Schur factors, the conditioning row and the residual.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import constant_schur_flow, coupled_flow_factors
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 1] = block[1, 0] = 1j
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> coupled_flow_factors(flow).frozen.shape
(1, 3, 2, 2)

solprop.translate.coupled.dense_metric(omega, z_edges, sampled, *, defect, gauge_id, certificate_id)[source]#

Build the dense metric of a basis from its own sampled flux Gram.

The companion of orthonormal_metric() for a basis whose flux metric is no longer the identity in any gauge it could be written in. The caller supplies W(omega) measured on the basis itself – the Hermitian Poynting pairing solprop.modes.vector.flux_product() integrates – and this stores it as the documented leaf of contract section 3, with the eigenvalue certificate a controller norm and a noise draw need.

Why it is genuinely dense. Under a declared unitary change of basis U a reference whose sampled Gram is I + delta becomes I + U.H delta U. The exact identity is invariant, but the measured defect is not: the off-diagonal of U.H delta U is not the off-diagonal of delta, and it is the measured object a run carries. Storing a literal identity here instead would be a decorated identity rather than a declaration of what the basis is.

The metric is constant along z in this constructor: one interval carrying one Chebyshev coefficient, which is degree zero, because the general modal route propagates constant geometry and its reference is a snapshot.

Parameters:

omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, increasing, shape (n_omega,).
z_edgesnumpy.ndarray: Stored positions in m, increasing, covering the run.
samplednumpy.ndarray: The Hermitian W(omega), shape (n_modes, n_modes, n_omega), dimensionless and positive definite at every frequency.
defectfloat: The measured norm(W - I, 2) of the reference basis this one was built from, non-negative and dimensionless: the table’s gram_defect, which a unitary transform leaves unchanged in the spectral norm.
gauge_idstr: The basis’s composed phase and coordinate identity.
certificate_idstr: Identity of the certificate defect came from.

Returns:

ModalMetricTables: The validated dense metric, with the extreme eigenvalues of the stored block per frequency.

Raises:

TranslationError: If the sampled metric is not a square positive-definite Hermitian block on the given frequencies.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import dense_metric
>>> sampled = np.broadcast_to(np.eye(2, dtype=complex)[:, :, None], (2, 2, 3))
>>> metric = dense_metric(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1.0]), sampled,
...     defect=1e-12, gauge_id="lab+circular(0,1)", certificate_id="sampled-gram",
... )
>>> metric.kind, metric.n_modes
('dense', 2)

solprop.translate.coupled.eigenbasis_limit(n_modes)[source]#

Return the largest eigenbasis condition the compiled flow actions accept.

The native backend forms a flow action Z exp(l T) Z^H from the eigenbasis of the Schur factor, W exp(l Lambda) Y (schur_eigenbasis()), where its rounding is bounded by kappa (n + 11 + 2 a) u relative to the action, with kappa = || |W| ||_2 || |Y| ||_2 and a = max |l lambda|. The reference’s scaling-and-squaring is bounded by 16 u + 2 sqrt(2) n (n + 2) u at its smallest – no squaring, a vanishing exponent: the sixteen-term series and the two products of the back-transform – and grows with the step from there. The limit is the kappa at which the two routes’ vanishing-step rounding terms meet – a heuristic gate, not a bound comparison at every step length (the eigenbasis bound also carries the measured ||WY - I|| and a drift that grows with the step); measured, the eigenbasis route is at least as accurate as the reference’s on every case tried; a bin beyond it (a nearly degenerate pair with a nonnormal coupling, where kappa grows as the inverse of the separation) is formed by the compiled scaling-and-squaring instead. The bounds are derived in tests/flow_parity.py and docs/dev/backends.md section 3.9.

Parameters:

n_modesint: Channels of the flow’s block, positive.

Returns:

float: The limit, dimensionless: 2.04 for one channel, rising with n; a one-channel block has kappa = 1 exactly.

Examples

>>> from solprop.translate.coupled import eigenbasis_limit
>>> round(eigenbasis_limit(3), 3)
4.173

solprop.translate.coupled.flow_edge_slack(flow)[source]#

Return the tolerance a position is matched to a section edge within.

A step schedule lands on a declared kink by summing its steps, so it arrives within a few ulp of the total span rather than on it exactly; chebyshev._EDGE_SLACK_FRACTION of the span – the same fraction solprop.translate.chebyshev.PiecewiseChebyshev matches its own breakpoints to – makes the section a stage is evaluated on depend on the schedule’s intent rather than on its last bits.

It is also what PropagationProblem._check_schedule of solprop.translate.problem holds a fixed schedule to beside a coupled flow. That check’s own window is wider (1e-9 of the run, which is 1e-9 m for anything under a metre), and a schedule edge between the two lands a step anchored in the closing section but running past the corner, against a frozen block that is not the one its stages read: fifth order collapses to first with nothing said. The two windows must therefore be one quantity, and a schedule that falls between them is refused rather than silently integrated.

Parameters:

flowCoupledFlowTables: The leaf.

Returns:

float: The tolerance in m, positive.

solprop.translate.coupled.flow_section(flow, z)[source]#

Return the index of the section a position lies in.

Parameters:

flowCoupledFlowTables: The leaf.
zfloat: Position in m.

Returns:

int: The section index, clamped to the ends: a position a rounding outside the span belongs to the nearest section, as a step landing on the last plane is. A position within 1e-10 of the span of an interior edge belongs to the section starting there, because a schedule that lands on a kink arrives a few ulp short of it as often as a few ulp past, and the two must not name different sections for the same step.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import constant_schur_flow, flow_section
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 0] = block[1, 1] = -1.0
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> flow_section(flow, 2e-3)
0

solprop.translate.coupled.flow_weights(flow, z, section=None)[source]#

Return the weights that combine a coupled flow’s coefficients at one position.

The stage-bundle half of the seam’s evaluation form, for the coupled flow’s own (n_intervals, n_deg) coefficient block: a backend evaluates G(z) = sum_j w_j terms_j over the flattened block, exactly as it evaluates dPhi over solprop.translate.stage.StageTerms.phi (docs/dev/kernel-contract.md section 1). Every weight outside the section containing z is exactly zero, so no section’s expansion is ever extrapolated, and the connection is combined over the same weights because it shares the axes.

Which section a position on an edge belongs to is the step’s to say. A section edge is a kink: the declared taper corner the sections were placed on, where a' jumps and the connection with it. A stage that lands exactly on the end of its own step’s section must read the closing side of that corner, not the opening one, or one stage of that step evaluates a source discontinuously larger than the rest and the whole run drops to first order – measured, before the tie was broken this way, as a global error of h times 1.0e2 per metre on the qualification taper. section is therefore the step’s own, and it wins whenever the position lies inside its closed interval; a position genuinely outside it – which only a step spanning an edge produces, and the driver lands on every hard kink – falls back to its own section rather than extrapolating a Chebyshev expansion past its end, which diverges as exp(2 n sqrt(eps)). That fallback takes the section the stage closes in rather than the one that opens at it, for the same reason: a stage that ran past its step’s section arrived from below, and the residual it is measured against is the one-sided limit on the side it came from.

Parameters:

flowCoupledFlowTables: The leaf.
zfloat: Position in m.
sectionint or None, optional: The step’s own section, from solprop.translate.stage.StageScalars.flow_section. None takes the section containing z, which is what a caller with no step to speak of wants.

Returns:

numpy.ndarray: Float64, shape (n_intervals * n_deg,).

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import constant_schur_flow, flow_weights
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 0] = block[1, 1] = -1.0
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> flow_weights(flow, 5e-4).tolist()
[1.0]

solprop.translate.coupled.frozen_generator(flow, index)[source]#

Return the generator one section of a coupled flow is frozen at.

The frozen point is the section’s midpoint, and it is read out of the stored expansion rather than carried beside it: the translation layer’s certificate and the backend’s Schur factors are then the same matrix to the last bit, because they come from the same coefficients through the same evaluation. At the midpoint x = 0, so the weights are T_k(0) – one, zero, minus one, zero, … – and a section of one degree returns its only coefficient unchanged, which keeps the constant-Schur flow bitwise the stored block.

Parameters:

flowCoupledFlowTables: The leaf.
indexint: Which section, in [0, n_intervals).

Returns:

numpy.ndarray: The frozen generator in 1/m, shape (n_modes, n_modes, n_omega).

Raises:

TranslationError: If the flow carries the factored representation, which has no dense generator to freeze.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import constant_schur_flow, frozen_generator
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 0] = block[1, 1] = -1.0
>>> flow = constant_schur_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1e-3]), block,
...     dz_max=1e-4, condition_limit=10.0, gauge_id="lab",
...     certificate_id="schur",
... )
>>> bool(np.array_equal(frozen_generator(flow, 0), block))
True

solprop.translate.coupled.frozen_residual_flow(omega, z_edges, generator, connection, *, dz_max, condition_limit, representation_error, gauge_id, certificate_id)[source]#

Build the certified coupled flow of a varying guide, frozen section by section.

The coupled flow at varying geometry. Each section carries the complete generator as a Chebyshev expansion along z; the block frozen inside the exponential is that expansion at the section midpoint (frozen_generator()), and everything else –

R_s(z) = D(z) - D(z_s) + H(z) - H(z_s) - C(z)

– is the explicit residual the interaction-picture Runge-Kutta source carries. With ``R`` evaluated exactly, freezing is a rearrangement of the original equation, not a second-order approximation; section refinement reduces the residual’s stiffness and the inverse amplification and is not permission to discard it.

The residual row is a true enclosure over each whole section rather than a sampled maximum: with R_s = sum_k r_k T_k and |T_k| <= 1, sum_k ||r_k||_2 bounds it at every position, extrema no sample saw included. The interpolant’s own error is the representation_error the caller measured, which joins the Schur diagnostic in generator_error.

Parameters:

omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, increasing.
z_edgesnumpy.ndarray: The section edges in m, shape (n_intervals + 1,), covering the run.
generatornumpy.ndarray: Chebyshev coefficients along z of the complete internal-convention generator in 1/m, shape (n_intervals, n_deg, n_modes, n_modes, n_omega), exactly zero on every bin the model was not evaluated at.
connectionnumpy.ndarray or None: The taper connection’s coefficients in the same layout and units, or None for a guide whose radius does not move.
dz_maxfloat: Largest step the run may take, m, positive: what the conditioning bound is declared over.
condition_limitfloat: Largest admitted conditioning, dimensionless and at least one.
representation_errorfloat: The relative Chebyshev tail of the expansions, dimensionless and non-negative, as section_expansion() measured it.
gauge_idstr: The basis’s phase and coordinate identity.
certificate_idstr: Identity of the certificate the bounds came from.

Returns:

CoupledFlowTables: The validated leaf, section_methods == ("frozen_residual",) * n_intervals.

Raises:

TranslationError: If the arrays disagree with the axes, if representation_error is negative, or if the conditioning bound exceeds condition_limit at any frequency.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import frozen_residual_flow
>>> block = np.zeros((2, 1, 2, 2, 3), dtype=complex)
>>> block[:, 0, 0, 0] = block[:, 0, 1, 1] = -1.0
>>> flow = frozen_residual_flow(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 5e-4, 1e-3]), block, None,
...     dz_max=1e-4, condition_limit=10.0, representation_error=0.0,
...     gauge_id="lab", certificate_id="frozen-residual",
... )
>>> flow.section_methods, float(flow.residual_norm_bound.max())
(('frozen_residual', 'frozen_residual'), 0.0)

solprop.translate.coupled.interval_count(edges)[source]#

Return how many z intervals a set of stored positions declares.

Parameters:

edgesint: Number of stored positions, at least one.

Returns:

int: edges - 1 for a path, and 1 for a single snapshot, whose one interval is the position it was taken at. A record’s interval-indexed arrays – the error rows, the coefficient representations – use this count, while its dense arrays are indexed by the positions themselves.

Examples

>>> from solprop.translate.coupled import interval_count
>>> interval_count(1), interval_count(4)
(1, 3)

solprop.translate.coupled.metric_at(metric, z)[source]#

Evaluate a basis metric at one statistics or observation plane.

The one arithmetic path from the stored leaf to the W a statistic is handed, so that the run-time reduction and a later check of it read the same numbers. None is the identity and is returned as None rather than as an array of ones: the statistics take metric=None as their fast, bitwise-unchanged path, and building an identity here would make every unweighted run pay for a matrix product it does not need.

Parameters:

metricModalMetricTables or None: The leaf, or None for a problem that carries none.
zfloat: Position in m inside the metric’s own span.

Returns:

numpy.ndarray or None: W(omega), shape (n_modes, n_modes, n_omega), dimensionless; None for an identity metric or no metric at all.

Raises:

TranslationError: If the metric is diagonal, which no route builds and nothing executes.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import dense_metric, metric_at
>>> sampled = np.broadcast_to(np.eye(2, dtype=complex)[:, :, None], (2, 2, 3))
>>> metric = dense_metric(
...     np.linspace(1e15, 2e15, 3), np.array([0.0, 1.0]), sampled,
...     defect=0.0, gauge_id="lab+circular(0,1)", certificate_id="sampled-gram",
... )
>>> metric_at(metric, 0.5).shape
(2, 2, 3)
>>> metric_at(None, 0.5) is None
True

solprop.translate.coupled.orthonormal_metric(omega, z_edges, *, defect, gauge_id, certificate_id)[source]#

Build the identity metric of a basis certified orthonormal in true power.

A Poynting-normalised bound basis whose sampled Gram passed the construction gate of solprop.translate.modal_tables.DenseModalTable has W = I to within that gate’s budget, so the metric owns no coefficient array at all and the measured defect is what the error row carries. The extreme eigenvalues are one: a certified positive metric, which is what the noise draw and the controller norm need.

Parameters:

omeganumpy.ndarray: The canonical coarse angular frequencies in rad/s, increasing.
z_edgesnumpy.ndarray: Stored positions in m, increasing, covering the run.
defectfloat: The measured norm(W - I, 2) of the basis, non-negative and dimensionless: the table’s gram_defect.
gauge_idstr: The basis’s phase and coordinate identity.
certificate_idstr: Identity of the certificate the defect came from.

Returns:

ModalMetricTables: The validated identity metric.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import orthonormal_metric
>>> metric = orthonormal_metric(
...     np.linspace(0.0, 2e15, 3), np.array([0.0, 1.0]), defect=1e-12,
...     gauge_id="lab", certificate_id="sampled-gram",
... )
>>> metric.kind, metric.diagonal is None
('identity', True)

solprop.translate.coupled.residual_truncation(x)[source]#

Return the frozen-residual truncation defect of one step, dimensionless.

||R(h R_s) - exp(h R_s)|| <= x^6/3600 + exp(x) x^7/5040 with x = h K, valid in any matrix norm and without an eigenvector assumption. The first term is |1/600 - 1/720| x^6, the difference between the production DOPRI5 tableau’s own sixth-order coefficient 1/600 and the exponential’s 1/720; the second bounds the exponential tail from the seventh order.

It is a local frozen defect, not a bound for a time-varying stage sequence, and it is charged to the linear share of the error budget, with the Schur action and the explicit-linear integration.

Parameters:

xfloat: h K, dimensionless and non-negative.

Returns:

float: The bound, dimensionless; inf where the exponential tail overflows a double, which a step so far above the advisory one reaches (a declared dz_max of half a capillary, for instance) and which the advisory warning then reports rather than failing on.

Raises:

TranslationError: If x is negative or not finite.

Examples

>>> from solprop.translate.coupled import residual_truncation
>>> f"{residual_truncation(0.1):.4e}"
'2.9971e-10'
>>> residual_truncation(0.0)
0.0
>>> residual_truncation(1e3)
inf

solprop.translate.coupled.schur_blocks(generator)[source]#

Return the complex Schur factors of a generator, one per frequency.

G = Z T Z^H with Z unitary and T upper triangular, which is the decomposition a matrix-function action of a nonnormal generator is entitled to: an eigenvector basis of a wall block near a degeneracy is arbitrarily ill-conditioned, and the Schur factor is not. A frequency whose block is exactly zero – an out-of-band bin the wall formula was never evaluated at – is not decomposed at all: its factors are the identity and zero, which is what the decomposition would give and costs nothing.

Parameters:

generatornumpy.ndarray: The generator in 1/m, complex, shape (n_modes, n_modes, n_omega).

Returns:

tuple of numpy.ndarray: Z and T, each complex of shape (n_omega, n_modes, n_modes): the frequency axis leads, because that is the axis a batched action broadcasts over.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import schur_blocks
>>> block = np.zeros((2, 2, 3), dtype=complex)
>>> block[0, 1, 1] = 3.0
>>> vectors, triangular = schur_blocks(block)
>>> vectors.shape, float(abs(triangular[1]).max())
((3, 2, 2), 3.0)

solprop.translate.coupled.schur_eigenbasis(vectors, triangular)[source]#

Return the eigenbasis of each Schur factor, and where it is well conditioned.

With G = Z T Z^H and T upper triangular with distinct diagonal lambda, T = V Lambda V^-1 with V unit upper triangular, so every action of the flow over a length l is W exp(l Lambda) Y with W = Z V and Y = V^-1 Z^H: n exponentials and one n^3 product, shared by the forward and inverse actions and by every stage length of a step. V is formed column by column by back substitution, V_ij = sum_k T_ik V_kj / (lambda_j - lambda_i), and its inverse by the division-free unit triangular recurrence, each sum taken from zero in increasing k.

Where two diagonal entries nearly coincide with a coupling between them, V is ill-conditioned and the route loses accuracy as the inverse of the separation; where they coincide exactly it does not exist. A bin is therefore marked for this route only where V is finite and || |W| ||_2 || |Y| ||_2 is within eigenbasis_limit(); the others are formed by scaling and squaring. Their W and Y are returned as zeros.

Parameters:

vectorsnumpy.ndarray: The unitary Schur factors Z, complex, shape (n_omega, n_modes, n_modes), as schur_blocks() returns them.
triangularnumpy.ndarray: The triangular factors T, 1/m, of the same shape.

Returns:

tuple of numpy.ndarray: W and Y, complex of the input’s shape; the eigenvalues lambda, 1/m, shape (n_omega, n_modes); and the boolean route, shape (n_omega,), true where the eigenbasis is used.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import schur_eigenbasis
>>> vectors = np.eye(2, dtype=complex)[None]
>>> triangular = np.array([[[-1.0, 0.1], [0.0, -2.0]]], dtype=complex)
>>> left, right, eigenvalues, route = schur_eigenbasis(vectors, triangular)
>>> bool(route[0]), bool(np.allclose(left[0] @ right[0], np.eye(2)))
(True, True)

solprop.translate.coupled.section_expansion(sample, z_edges, *, tol=1e-12, max_degree=32, name='the section expansion')[source]#

Expand an array-valued function of z in Chebyshev coefficients per section.

The fixed-section counterpart of solprop.translate.chebyshev.build(): the edges are given – they are the flow’s sections, which already land on every declared kink – so nothing is split and nothing is guessed, and the only search is the degree. One degree is resolved for every section, because the leaf’s n_deg axis is rectangular, and it is the largest any section needed.

The expansion is what makes the explicit residual a bounded quantity rather than a sampled one: |T_k| <= 1, so sum_k ||c_k|| encloses the function over the whole section including the extrema no sample sees (frozen_residual_flow()).

Parameters:

samplecallable: Takes positions in m, shape (n_z,), and returns the values with the position axis first, shape (n_z, *value_shape).
z_edgesnumpy.ndarray: The section edges in m, shape (n_intervals + 1,), strictly increasing.
tolfloat, optional: Relative coefficient tail to resolve to. Default SECTION_TOL.
max_degreefloat, optional: Highest degree tried. Default MAX_SECTION_DEGREE.
namestr, optional: What the refusal calls the quantity.

Returns:

tuple: The coefficients, shape (n_intervals, n_deg, *value_shape), and the worst relative tail actually reached, dimensionless.

Raises:

TranslationError: If a section’s tail is still above tol at max_degree.

Examples

>>> import numpy as np
>>> from solprop.translate.coupled import section_expansion
>>> coefficients, tail = section_expansion(
...     lambda z: np.exp(-z)[:, None], np.array([0.0, 0.5, 1.0])
... )
>>> coefficients.shape[0], bool(tail < 1e-12)
(2, True)

solprop.translate.coupled.taper_identity(flow, *, channels, active, scored=None, n_components=1, has_noise=False)[source]#

Return the modal identity of a capillary taper’s coupling.

A multimode capillary taper without a wall carries the taper coupling in a coupled flow, and a file must say so: a resume against the adiabatic route of the same geometry is a different equation. The record is wall_identity()’s with the wall’s own entries – its validity ratios, its radius and its impedance – absent, because no wall block was built; flow.model is "taper coupling" and the certificate is solprop.translate.problem.TAPER_FLOW_CERTIFICATE. A new key entry and not a new schema, as for wall_identity(). The basis is still the orthogonal "capillary HE_1m" set, so the ordinary accessors stay admitted.

Parameters:

flowCoupledFlowTables: The leaf the run carries.
channelssequence of solprop.modes.protocol.ChannelLabel: The channel labels in state order.
activenumpy.ndarray of bool: The run’s active band, whose ends are recorded.
scorednumpy.ndarray of bool, optional: The error norm’s bins, over which the recorded advisory step is taken (advisory_step()); None takes every bin.
n_componentsint, optional: The transverse components the synthesis resolves.
has_noisebool, optional: Whether the run carries an input-noise draw, recorded as wall_identity() records it.

Returns:

collections.abc.Mapping: The modal record of solprop.output.metadata.RunMetadata.

Examples

>>> import numpy as np
>>> from solprop.modes.protocol import ChannelLabel
>>> from solprop.translate.coupled import frozen_residual_flow, taper_identity
>>> omega = np.array([2.3e15, 2.4e15])
>>> block = np.zeros((1, 1, 2, 2, 2), dtype=complex)
>>> block[0, 0, 1, 1] = -1j * 1e3
>>> flow = frozen_residual_flow(
...     omega, np.array([0.0, 1e-3]), block, None, dz_max=1e-5,
...     condition_limit=100.0, representation_error=0.0,
...     gauge_id="capillary HE_1m", certificate_id="taper-frozen-residual",
... )
>>> labels = [ChannelLabel("HE", 1, m, None, None, f"HE1{m}") for m in (1, 2)]
>>> record = taper_identity(flow, channels=labels, active=np.ones(2, bool))
>>> record["flow"]["model"], record["basis"]
('taper coupling', 'capillary HE_1m')

solprop.translate.coupled.wall_identity(flow, limits, *, channels, impedance='mean', scored=None, n_components=1, has_noise=False)[source]#

Return the modal identity of a capillary carrying a coherent wall.

modal.flow joins the identity the general modal route writes, as modal.guide.sampling.basis does. A route without a coupled flow writes no flow entry.

The identity is what a resume compares, so no array enters it and every measured float is quantised (solprop.modes.material.quantised_measurement()) against the last bits of the BLAS that produced it. What is recorded of the block is what makes it the block it is: the section methods and how many of them there are, the degree of their expansions along z, whether a taper connection is carried, the two identity strings, the three validity ratios, the band they were taken over, the three certificate rows reduced to their extremes and the advisory step they imply (advisory_step()). The radius and the ratios are the worst section’s (solprop.modes.capillary.worst_wall_limits()), which for a constant guide is the only one there is, taken over the error norm’s scored band as the run’s admission takes them.

modal.basis is "capillary HE_1m", which is how solprop.output.capability.capability_of() tells this identity from a general guide’s: the coherent wall does not change the basis, which is the same orthogonal HE_1m set it always was, only the flow through it.

Parameters:

flowCoupledFlowTables: The leaf the run carries.
limitssolprop.modes.capillary.WallLimits: The validity ratios the block was admitted under.
channelssequence of solprop.modes.protocol.ChannelLabel: The channel labels in state order.
impedancestr, optional: Which wall the run declared. A "te_tm" block writes the model "te_tm impedance coherent wall" and the basis "capillary azimuthal order" rather than "capillary HE_1m", because its channels are two radial families of one azimuthal order and not one HE_1m set; it is still an orthogonal capillary mode set, so solprop.output.capability.capability_of() still reads it as the ordinary modal basis and the ordinary accessors stay admitted. A resume across a changed wall impedance is refused by this field.
scorednumpy.ndarray of bool, optional: The error norm’s bins, over which the recorded advisory step is taken, as the run’s warning takes it (advisory_step()); None takes every bin.
n_componentsint, optional: The transverse components the synthesis resolves, 1 or 2: two for a supplied TE_0m, TM_0m or EH set, one for the capillary’s own HE_1m set.
has_noisebool, optional: Whether the run carries an input-noise draw, which is made per channel in the orthonormal true-power basis and recorded as such.

Returns:

collections.abc.Mapping: The modal record of solprop.output.metadata.RunMetadata.

Examples

>>> import numpy as np
>>> from solprop.modes import MarcatiliMode, wall_operator
>>> from solprop.translate.coupled import constant_schur_flow, wall_identity
>>> modes = [MarcatiliMode(125e-6, "He", 1.0, m=index) for index in (1, 2)]
>>> omega = np.array([2.3e15, 2.4e15])
>>> block = wall_operator(modes, omega)
>>> flow = constant_schur_flow(
...     omega, np.array([0.0, 1e-3]), np.conj(block.generator), dz_max=1e-5,
...     condition_limit=100.0, gauge_id="capillary HE_1m",
...     certificate_id="constant-schur",
... )
>>> record = wall_identity(
...     flow, block.limits, channels=[mode.channel_label for mode in modes]
... )
>>> record["basis"], record["flow"]["section_methods"]
('capillary HE_1m', ['constant_schur'])