solprop.modes.connection#

Derive the connection of a basis whose material interfaces move along z.

A guide whose transverse profile moves is written in a moving basis, and a moving basis costs a term. With the field expanded on the local modes at every position, \(E(r, z) = \sum_n c_n(z)\,\psi_n(r; z)\), the coefficients obey

\[\frac{dc}{dz} = \bigl(\text{the local generator}\bigr) c - \Gamma c , \qquad \Gamma = W^{-1} S , \qquad S_{mn} = \langle \psi_m \mid \partial_z \psi_n \rangle ,\]

where \(\langle\cdot\mid\cdot\rangle\) is this package’s own flux metric – the Hermitian Poynting pairing of solprop.modes.vector.flux_product() – and \(W\) is that pairing’s Gram matrix, the same flux_metric the source projection map Q is built from. Using one pairing for the dual, the projection and the connection is what makes the three consistent; the bilinear reciprocity pairing never enters, exactly as for the source path.

The flux self-terms are kept, and this module exists because they must be. Differentiating the Gram matrix gives the one identity the whole certificate rests on:

\[W' = \frac{dW}{dz} = \langle \partial_z\psi_m \mid \psi_n\rangle + \langle \psi_m \mid \partial_z\psi_n\rangle = S^{H} + S , \qquad\text{equivalently}\qquad W\Gamma + \Gamma^{H} W = W' .\]

So \(\Gamma\) is anti-Hermitian with respect to the flux metric, up to that metric’s own drift, and the consequence is exact: with \(dc/dz = -\Gamma c\) the flux norm obeys \(d(c^{H} W c)/dz = c^{H}(W' - W\Gamma - \Gamma^{H}W)c = 0\). The connection alone moves no power. Zeroing the diagonal – which the Bessel formula for the capillary does, correctly, for its basis – breaks that identity. For diagonal \(W\), the resulting norm defect is \(\sum_n |c_n|^2 W'_{nn}\), which is why solprop.modes.capillary.wall_connection()’s formula is not transplanted here.

The diagonal splits into the two pieces the identity names:

for diagonal \(W\), its real part is \(W'_{nn}/(2W_{nn})\), reducing to \(\tfrac12 W'_{nn}\) where \(W_{nn}=1\). This flux self-term describes the rate of change of the mode’s power normalisation. For a dense metric the full matrix identity must be used;
its imaginary part is the gauge phase rate, a geometric (Berry) connection. It vanishes only in a gauge that parallel-transports the phase, which is what solprop.translate.modal_tables.align_modal_sample() does to first order by making the overlap with the previous knot real and positive.

Neither is zero in general for a moving interface, and a declared z-varying change of coordinates \(c = T(z) d\) adds \(T^{H}T'\) to the connection (and its negative to the generator) in the declared basis, whose diagonal is pure imaginary and is exactly the second piece.

How it is computed. \(\partial_z\psi_n\) is a three-point Lagrange derivative of the provider’s own solves, tracked into one gauge, on one shared node set whose panel edges land on every stencil position’s interfaces (connection_nodes()). That last point is not a detail: a rule built at one position and reused at another leaves its panel edge off the moving interface, and the quadrature of a Poynting integrand with a jump at the interface is then wrong by a term first order in the displacement – which is first order in the step and therefore survives the difference undiminished. Measured on a step-index taper it puts a spurious 0.98 1/m on the diagonal and reads |W'| = 2.8 1/m where the true value is below 1e-8. Splitting the panels at every stencil radius restores piecewise-smooth quadrature. Finite quadrature, mode-solve and derivative errors remain and require independent refinement.

The step is chosen so that the difference’s truncation, which falls as h^2, and its round-off, which grows as 1/h, are balanced (STEP_FRACTION). The measured ModalConnection.antisymmetry is the certificate, and it is the residual of the identity above with W' differenced independently of the pairing, in 1/m.

What this is not. The connection describes the basis the provider retains. What lies outside the retained set does not vanish at a material interface while every retained basis function’s expansion of it converges only in \(L^2\), so a propagated coefficient departs from the direct projection of the same frozen field at second order in the total contraction – as the capillary’s wall connection does in the stated comparisons. This observed scaling is not a general error bound for arbitrary taper shapes or contractions.

Module Attributes

`STEP_FRACTION`	Difference step as a fraction of the z scale, dimensionless: 1.53e-5.
`MAX_STENCIL_SPAN`	Largest share of the z scale a difference stencil may span, dimensionless.
`CERTIFICATE_WARNING_RTOL`	Relative antisymmetry above which a measured connection warns, dimensionless.
`NEGLIGIBLE_CONVERSION`	Below this much conversion over the z domain a connection is rounding, dimensionless.

Functions

`band_connection`(provider, omega, z, *[, ...])	Measure the moving-basis connection at one position over a band of frequencies.
`basis_connection`(transform, derivative)	Return the connection a declared z-varying change of coordinates costs.
`connection_nodes`(provider, omega, positions, *)	Build a shared node set split at every stencil position's panel edges.
`modal_connection`(provider, omega, z, *[, ...])	Measure the moving-basis connection of a provider at one point.

Classes

ModalConnection(matrix, metric, metric_rate, ...)

Carry one position's moving-basis connection with the certificate it met.

solprop.modes.connection.CERTIFICATE_WARNING_RTOL: float = 0.0001#

Relative antisymmetry above which a measured connection warns, dimensionless.

antisymmetry / max|Gamma|. The default step meets 7.4e-11 and 5.7e-9 on the two qualified providers and the tests admit 1e-8 there; a stencil shifted onto one side of a declared kink reads 1e-5, round-off limited and harmless; a caller declaring a step of a fifth of the scale reads 2.6e-2 with the connection itself 1.2e-2 wrong while still inside MAX_STENCIL_SPAN. The warning sits at 1e-4, where the connection is itself wrong at a part in ten thousand, so the certificate is what says a declared step was too coarse. modal_connection() warns by name above it and does not refuse: the certificate is the caller’s to read, and a noisy provider may legitimately sit above it.

solprop.modes.connection.MAX_STENCIL_SPAN: float = 0.25#

Largest share of the z scale a difference stencil may span, dimensionless.

The derivative is local, so a stencil that reaches a quarter of the way across the scale it was sized from is not measuring a local slope. A caller that asks for a step this large is refused rather than given a secant.

It is also what makes the three positions fit: at a quarter of the scale a stencil clamped against either end of the interval still lies inside it, so modal_connection() needs no second guard for that. A value above a half would, which is why this one is not a tuning knob.

class solprop.modes.connection.ModalConnection(matrix, metric, metric_rate, antisymmetry, step, z, omega)[source]#

Bases: Displayed

Carry one position’s moving-basis connection with the certificate it met.

Attributes:

matrixnumpy.ndarray: The connection Gamma = W^-1 S in 1/m, complex, shape (n_modes_connection_out, n_modes_connection_in). The coefficients obey dc/dz = ... - Gamma c.
metricnumpy.ndarray: The flux Gram W at the same position, dimensionless and Hermitian, shape (n_modes_metric_out, n_modes_metric_in).
metric_ratenumpy.ndarray: Its derivative W' in 1/m, same shape, differenced from the stencil’s own Gram matrices. The derivation makes it S + S^H, twice the Hermitian part of the pairing, and antisymmetry is what that costs to meet.
antisymmetryfloat: The certificate, in 1/m: norm(W Gamma + Gamma^H W - W', 2) with W' differenced independently of the pairing, which the derivation makes exactly zero. What is measured is the stencil’s own truncation and round-off – and a quadrature whose panel edge has left a moving interface, which is what it exists to catch.
stepfloat: The difference step actually used, m, positive.
zfloat: The position the connection was taken at, m.
omegafloat: The frequency it was taken at, rad/s, positive.

Raises:

solprop.exceptions.ModeError: If the matrices are not square and equally wide, or a scalar is out of range.

Examples

>>> import numpy as np
>>> from solprop.modes.connection import ModalConnection
>>> link = ModalConnection(
...     np.array([[0.0, 1.0], [-1.0, 0.0]], dtype=complex),
...     np.eye(2, dtype=complex), np.zeros((2, 2), dtype=complex),
...     0.0, 1e-9, 0.5, 2.3e15,
... )
>>> link.n_modes, float(np.max(abs(link.self_terms)))
(2, 0.0)

antisymmetry: float#

matrix: Complex[ndarray, 'n_modes_connection_out n_modes_connection_in']#

metric: Complex[ndarray, 'n_modes_metric_out n_modes_metric_in']#

metric_rate: Complex[ndarray, 'n_modes_rate_out n_modes_rate_in']#

property n_modes: int[source]#

Return how many channels the connection couples.

Returns:

int: The channel count.

omega: float#

property self_terms: Complex[ndarray, 'n_modes_self'][source]#

Return the diagonal: the flux self-term and the gauge phase rate.

The real part is W'_nn/2, the rate the channel’s own power normalisation moves; the imaginary part is the gauge’s phase rate. Neither is zero in general for a moving material interface.

Returns:

numpy.ndarray: Complex, shape (n_modes,), in 1/m.

step: float#

z: float#

solprop.modes.connection.NEGLIGIBLE_CONVERSION: float = 1e-06#

Below this much conversion over the z domain a connection is rounding, dimensionless.

max|Gamma| * (z_high - z_low): a degenerate pair’s connection is zero by symmetry and reads ~1e-8 1/m of solve noise, so its relative certificate is meaningless and modal_connection() does not warn on it.

solprop.modes.connection.STEP_FRACTION: float = 1.525878906250001e-05#

Difference step as a fraction of the z scale, dimensionless: 1.53e-5.

The three-point derivative’s truncation falls as (h/L)^2 and its round-off grows as delta L/h, where L is the scale the sampled fields vary over and delta is the relative noise of a solve; the two are equal at h/L = (delta/2)^(1/3). delta = 16 eps – the rounding of a short reduction, not of one operation – is what a provider that solves an eigenproblem and integrates a quadrature actually carries, and the measurement agrees. On a step-index taper of L = 1e-3 m and a connection whose own scale is 3.05e2 1/m, the certificate ModalConnection.antisymmetry reads 5.3e-10, 2.8e-10, 7.4e-11, 5.6e-10, 2.4e-9 and 9.5e-9 relative at h/L of 6.1e-6, 1.0e-5, 1.53e-5, 3.0e-5, 6.0e-5 and 1.2e-4, and the connection itself agrees to ten digits across the whole scan. A caller whose provider is noisier declares its own step and reads the certificate.

solprop.modes.connection.band_connection(provider, omega, z, *, nodes=None, step=None, degenerate_blocks=(), basis=None, n_r=None, n_theta=None)[source]#

Measure the moving-basis connection at one position over a band of frequencies.

modal_connection() at every frequency of omega, with one difference that is a cost and not a physics choice: the stencil and its union node rule (connection_nodes()) are built once, at the lowest frequency of the band, where a bound field reaches furthest and the rule’s certified tail is a bound over the whole band – the same reason solprop.modes.guide.ModalGuide.from_provider() asks for its rule there. A rule with a certified tail costs a provider solve of its own, so building three of them per frequency would cost more than the connection itself. Each frequency is then measured and certified exactly as modal_connection() measures it: the three solves tracked into the centre’s gauge, the declared basis applied after the tracking, Gamma = W^-1 S on the flux pairing and the antisymmetry certificate against W' differenced independently.

This is what a moving guide’s path needs: the connection at every sampled position and every frequency knot of the guide, interpolated in frequency onto the run’s bins as the fields are.

Parameters:

providersolprop.modes.vector.ModalProvider: The guide, solved at three positions per frequency.
omegasequence of float or numpy.ndarray: Positive angular frequencies in rad/s inside the provider’s band, shape (n_omega,), at least one.
zfloat: Position in m inside the provider’s z domain.
nodessolprop.modes.vector.ModalNodeRule or None, optional: A further rule joining the union, as connection_nodes() takes it.
stepfloat or None, optional: The difference step in m, as modal_connection() takes it.
degenerate_blockssequence of sequence of int, optional: Channel blocks the tracking may rotate inside.
basissolprop.modes.guide.ChannelBasis or None, optional: A declared constant unitary change of the channel coordinates.
n_r, n_thetaint or None, optional: The resolution each stencil position’s rule is asked for (connection_nodes()).

Returns:

tuple of ModalConnection: One per frequency, in the order given.

Raises:

solprop.exceptions.ModeError: If omega is empty, or on anything modal_connection() refuses.

Examples

>>> import numpy as np
>>> from solprop.modes.connection import band_connection
>>> from solprop.modes.graded_index import (
...     GradedIndexProvider, PowerLawPermittivity,
... )
>>> from solprop.modes.material import IndexMaterial
>>> from solprop.profiles import Linear
>>> guide = GradedIndexProvider(
...     permittivity=PowerLawPermittivity(IndexMaterial(1.45), 0.005),
...     exterior=IndexMaterial(1.45 * np.sqrt(1 - 2 * 0.005)),
...     radius=Linear(25e-6, 26e-6, 1e-3), family="LP", count=2,
...     elements=16, extent=2.0,
... )
>>> links = band_connection(guide, [1.80e15, 1.84e15], 5e-4)
>>> len(links), links[0].n_modes
(2, 2)

solprop.modes.connection.basis_connection(transform, derivative)[source]#

Return the connection a declared z-varying change of coordinates costs.

In general, for \(c = T(z) d\) the declared coordinates carry \(-T^{-1}T'\) beside the transformed generator, so the connection in the declared basis gains \(T^{-1}T'\). For a unitary T – the only kind ChannelBasis declares – that is \(T^{H}T'\), which needs no inverse and is exactly anti-Hermitian because differentiating \(T^{H}T = I\) gives \(T^{H}T' + (T^{H}T')^{H} = 0\). Its diagonal is therefore pure imaginary and is the gauge phase rate of this module’s opening derivation: a relabelling that produces no physical conversion, which is the pure change-of-basis test.

Parameters:

transformnumpy.ndarray: The unitary T(z), shape (n_modes, n_modes), dimensionless.
derivativenumpy.ndarray: T'(z) in 1/m, the same shape.

Returns:

numpy.ndarray: T^H T' in 1/m, anti-Hermitian, shape (n_modes, n_modes).

Raises:

solprop.exceptions.ModeError: If the two are not square matrices of the same width.

Examples

>>> import numpy as np
>>> from solprop.modes.connection import basis_connection
>>> angle, rate = 0.3, 2.0
>>> turn = np.array([[np.cos(angle), -np.sin(angle)],
...                  [np.sin(angle), np.cos(angle)]], dtype=complex)
>>> slope = rate * np.array([[-np.sin(angle), -np.cos(angle)],
...                          [np.cos(angle), -np.sin(angle)]], dtype=complex)
>>> link = basis_connection(turn, slope)
>>> bool(np.allclose(link, -link.conj().T)), round(float(np.real(link[0, 1])), 12)
(True, -2.0)

solprop.modes.connection.connection_nodes(provider, omega, positions, *, nodes=None, n_r=None, n_theta=None)[source]#

Build a shared node set split at every stencil position’s panel edges.

The rules returned share their coordinates and their area weights exactly – that is what makes a difference of the sampled fields meaningful – and differ only in region_index, which each position computes from its own material radii. The radial panels are split at the union of every position’s own declared panel edges, so each solve’s integrand is smooth on every panel and quadrature can converge on its smooth pieces, which a single frozen rule cannot ensure (see this module’s own docstring for what a frozen rule costs).

Taking the union of the declared edges rather than of the material radii alone is what makes this work for a discretised provider: a graded-index sector is a P1 finite-element field whose slope breaks at every element edge, and its whole mesh scales with the core radius, so a rule split at the interface alone integrates a kinked function with a Gauss rule and reads a flux metric 0.39 from the identity (measured).

Parameters:

providersolprop.modes.vector.ModalProvider: The guide, whose regions declare the radii that move and whose own node_rule declares the panel edges at each position.
omegafloat: Positive angular frequency in rad/s, at which each position’s rule is asked for.
positionssequence of float or numpy.ndarray: The stencil positions in m, inside the provider’s z domain, shape (n_z,).
nodessolprop.modes.vector.ModalNodeRule or None, optional: A further rule whose panel edges and resolution join the union, so that a caller can measure the connection on a superset of the quadrature its run uses. None takes the provider’s own rules alone.
n_r, n_thetaint or None, optional: The resolution each position’s own rule is asked for, forwarded to the provider’s node_rule; None takes the provider’s own default, which is what modal_connection() measures at. A moving guide’s path passes the resolution its own sampling rule declares, so that the connection and the maps it is propagated beside are resolved alike.

Returns:

tuple of ModalNodeRule: One rule per position, in the same order.

Raises:

solprop.exceptions.ModeError: If a rule does not declare the panel edges, the radial resolution and the truncation radius this needs, or none of its edges lies inside the domain.

solprop.modes.connection.modal_connection(provider, omega, z, *, nodes=None, step=None, degenerate_blocks=(), basis=None)[source]#

Measure the moving-basis connection of a provider at one point.

The derivation is this module’s own docstring: Gamma = W^-1 S with S_mn = <psi_m | d_z psi_n> on the package’s flux metric, with the flux self-terms kept. The three samples of the difference stencil are put in one gauge by solprop.translate.modal_tables.align_modal_sample() against the centre, so the phase the derivative measures is the parallel-transported one and not the solver’s arbitrary one, and a declared basis is applied after it in the same fixed order the tables compose.

The stencil is centred where it fits and shifted where it does not, at the ends of the provider’s z domain and against a declared kink; the derivative weights are the exact Lagrange ones of the three positions actually used, so the order is the same either way.

Parameters:

providersolprop.modes.vector.ModalProvider: The guide, solved at three positions.
omegafloat: Positive angular frequency in rad/s inside the provider’s band.
zfloat: Position in m inside the provider’s z domain.
nodessolprop.modes.vector.ModalNodeRule or None, optional: The rule to refine onto the moving interfaces (connection_nodes()); None asks the provider for one at z.
stepfloat or None, optional: The difference step in m; None takes STEP_FRACTION of the z scale, which is the provider’s domain narrowed to the interval between the kinks z lies between.
degenerate_blockssequence of sequence of int, optional: Channel blocks the tracking may rotate inside, as solprop.translate.modal_tables.sample_modal_table() takes them.
basissolprop.modes.guide.ChannelBasis or None, optional: A declared unitary change of the channel coordinates, applied after the tracking gauge. It is a constant declaration: a transform that varies with z adds basis_connection() on top, which the caller composes.

Returns:

ModalConnection: The connection, the metric, its measured rate and the certificate.

Raises:

solprop.exceptions.ModeError: If step is not positive and finite, if the stencil spans more than MAX_STENCIL_SPAN of the z scale or does not fit inside the z domain at all, or on anything the provider, the tracking or the declared basis refuses.

Examples

>>> import numpy as np
>>> from solprop.modes.connection import modal_connection
>>> from solprop.modes.graded_index import (
...     GradedIndexProvider, PowerLawPermittivity,
... )
>>> from solprop.modes.material import IndexMaterial
>>> from solprop.profiles import Linear
>>> guide = GradedIndexProvider(
...     permittivity=PowerLawPermittivity(IndexMaterial(1.45), 0.005),
...     exterior=IndexMaterial(1.45 * np.sqrt(1 - 2 * 0.005)),
...     radius=Linear(25e-6, 26e-6, 1e-3), family="LP", count=2,
...     elements=16, extent=2.0,
... )
>>> link = modal_connection(guide, 2 * np.pi * 3e8 / 1.03e-6, 5e-4)
>>> link.n_modes, bool(link.antisymmetry < 1e-6 * np.max(abs(link.matrix)))
(2, True)