Eigenmode expansion (EME)¶
Mode-matching propagation through a sequence of uniform sections. Use the convergence helpers to assess the retained modes and interfaces.
Run¶
run_eme
¶
run_eme(
sections: Sequence[Section],
n_modes: Optional[int] = None,
*,
interface_rcond: float = 1e-10,
enforce_passivity: bool = False,
interface: str = "galerkin",
) -> EMEResult
Cascade a sequence of Sections into one device S-matrix.
Parameters¶
sections:
Cross-sections in propagation order. The first and last define the
input/output ports (typically length-0 leads). All must share the
same transverse grid (ny, nx) and spacing.
n_modes:
Optional cap applied independently to every section. Adjacent sections
may retain different counts; their interface blocks are rectangular.
A cap may not split a degenerate set of modes (a validated Yee beta
group, or the TE/TM pair of a square or round core).
enforce_passivity:
Clip the gain of each ordinary real guided interface (see
interface_smatrix). Off by default: the clip biases the
fundamental transmission low by 0.1 to 0.4 %, and more modes or
slabs do not remove the bias. Without it, a cascade of ordinary
guided ports whose energy_balance(0) exceeds 1 by more than 1e-3
warns: that gain is a basis that has not converged.
interface
"galerkin" (default) or "projection", the interface model of
interface_smatrix. Use "projection" for an open device
modelled with a truncated (for example guided-only) basis, where power
that a section's modes cannot represent radiates away: gentle bends,
lateral offsets and other junctions whose scattered field does not
return to the guide. Keep "galerkin" with a complete basis; the two
agree once the basis is complete, which is a useful check.
Returns¶
EMEResult The device scattering matrix and convenience power ratios.
run_eme_band
¶
run_eme_band(
base_sections: Sequence[Tuple[VectorModeSolver, float]],
wlens_um: Sequence[float],
*,
num_modes: int,
neff_margin: float = 0.0,
n_clad: Optional[float] = None,
n_modes: Optional[int] = None,
solve_kwargs: Optional[dict] = None,
basis_kwargs: Optional[dict] = None,
workers: Optional[int] = None,
backend: str = "auto",
enforce_passivity: bool = False,
interface: str = "galerkin",
) -> Dict[float, EMEResult]
Cascade a staircased device across a band of wavelengths, solving the independent per-(wavelength × section) FDE eigenmodes in parallel.
Parameters¶
base_sections:
Cross-sections in propagation order as (base_solver, length_um)
pairs, e.g. from rectangular_base_section. Each cross-section is
sampled on the mesh once and reused across the whole band.
length_um == 0 marks a port
lead (as in Section). All sections must share one transverse grid.
wlens_um:
The free-space wavelengths (microns) to sweep.
num_modes:
Guided modes solved per cross-section. Passed to
VectorModeSolver.solve on the ordinary path, or used as
num_guided when basis_kwargs selects a complete EME basis.
neff_margin:
Drop modes whose n_eff is within this margin of n_clad (near-cutoff, see waveguide_section). Requires n_clad.
n_clad:
Cladding index, needed only when neff_margin > 0.
n_modes:
Cap on the cascade basis size, forwarded to run_eme.
solve_kwargs:
Extra keyword arguments forwarded to VectorModeSolver.solve
(e.g. n_guess); must be picklable for the process backend.
basis_kwargs:
When provided, solve a complete basis with
VectorModeSolver.solve_eme_basis instead. num_modes supplies
num_guided and this dictionary supplies options such as
num_radiation, num_evanescent and pml_cells_xy. It is
mutually exclusive with solve_kwargs and must be picklable for the
process backend.
workers:
Worker count for the eigensolves. None → min(#solves, os.cpu_count).
1 (or a single solve) runs serially with no pool.
backend:
"auto" (default), "process" on fork platforms (Linux/CI/cloud,
~3.5x), "thread" on spawn platforms (macOS/Windows, ~1.7x without the
process-startup tax or __main__-guard requirement). Force "process"
for a large spawn-platform sweep, or "serial" to disable parallelism.
enforce_passivity:
Forwarded to run_eme (off by default).
interface:
Forwarded to run_eme ("galerkin" by default).
Returns¶
dict[float, EMEResult]
{wavelength_um: EMEResult}, one cascaded device S-matrix per
wavelength, identical to a per-wavelength run_eme cascade.
EMEResult
dataclass
¶
EMEResult(
s11: ndarray,
s12: ndarray,
s21: ndarray,
s22: ndarray,
n_modes: int,
n_left_modes: int = 0,
n_right_modes: int = 0,
section_mode_counts: Tuple[int, ...] = (),
left_modes: Tuple[VectorMode, ...] = (),
right_modes: Tuple[VectorMode, ...] = (),
left_power_metric: Optional[ndarray] = None,
right_power_metric: Optional[ndarray] = None,
interface_diagnostics: Tuple[
InterfaceDiagnostics, ...
] = (),
interface: str = "galerkin",
)
Device scattering matrix from run_eme.
Blocks may be rectangular: S21 is N_right x N_left. Physical power
is evaluated with Hermitian Poynting metrics, so evanescent/PML amplitudes
are never mistaken for |a|^2 power channels. n_modes retains its
legacy meaning as the minimum section-basis count; use n_left_modes,
n_right_modes, and section_mode_counts for the rectangular model.
right_power_metric
class-attribute
instance-attribute
¶
interface_diagnostics
class-attribute
instance-attribute
¶
interface_diagnostics: Tuple[InterfaceDiagnostics, ...] = ()
transmission
property
¶
Fundamental-to-fundamental physical power transmission.
Index-based, so it follows n_eff ordering and NOT polarization, see _fundamental_index. Warns when the two ports disagree on
the polarization at index 0. Prefer transmission_of
(result.transmission_of("TE")) whenever the ports can reorder,
which is the normal case for a taper.
transmission_of
¶
Power transmission between the fundamental polarization modes
of the two ports, resolved by mode identity rather than by index.
Metric-aware like transmission (physical power ratio);
polarization=None reproduces it exactly.
transmitted_power
¶
Total physical output-port power from input_mode.
reflected_power
¶
Total physical reflected power from input_mode.
energy_balance
¶
T_total + R_total for input_mode.
Unlike the former forced-unitary interface, this is a physical flux
diagnostic: radiation/PML attenuation or an incomplete trace basis may
produce a deficit. A value above one by more than about 1e-3 is a
passivity/convergence failure (run_eme warns); below that it
is the residual of guided traces that are not exactly orthogonal.
Sections¶
Section
dataclass
¶
Section(
modes: Sequence[VectorMode], length_um: float = 0.0
)
One z-invariant cross-section of an EME device.
Attributes¶
modes:
The eigenmodes of this cross-section. A guided-only tuple from
VectorModeSolver.solve remains valid; a complete EME basis from
VectorModeSolver.solve_eme_basis may additionally contain
experimental radiation and evanescent modes.
length_um:
Physical length of the section along propagation (microns). 0.0 marks
a port / semi-infinite lead, it contributes its interface with its
neighbour but no propagation phase.
waveguide_section
¶
waveguide_section(
*,
wlen_um: float,
dl_um: float,
core_w_um: float,
core_h_um: float,
n_core: float,
n_clad: float,
window_w_um: float,
window_h_um: float,
num_modes: int,
num_radiation_modes: int = 0,
num_evanescent_modes: int = 0,
pml_cells_xy: Optional[Tuple[int, int]] = None,
pml_imaginary_thickness_um: Tuple[float, float] = (
0.1,
0.1,
),
length_um: float = 0.0,
neff_margin: float = 0.0,
subpixel: bool = True,
subpixel_method: str = "tensor",
) -> Section
Solve a centered rectangular-core cross-section and wrap it as a
Section.
Pass the same wlen_um, dl_um, window_w_um and
window_h_um for every section of a device so they share one transverse
grid (required by run_eme). Only core_w_um / core_h_um should
vary between sections.
neff_margin drops modes whose n_eff is within this margin of the
cladding index, i.e. near-cutoff modes. Those are poorly resolved (large
evanescent tails into the window walls) and break the within-section
orthonormality the interface relies on, so a small margin (e.g. 0.05)
keeps the modal basis clean. Default 0.0 keeps every guided mode the
solver returns.
rectangular_base_section
¶
rectangular_base_section(
*,
core_w_um: float,
length_um: float = 0.0,
dl_um: float,
core_h_um: float,
n_core: float,
n_clad: float,
window_w_um: float,
window_h_um: float,
ref_wlen_um: float = 1.31,
subpixel: bool = True,
subpixel_method: str = "tensor",
) -> Tuple[VectorModeSolver, float]
Sample a centered rectangular-core cross-section once and pair it with
its length, ready for a band/length sweep via run_eme_band.
Mirrors the geometry arguments of waveguide_section, but returns the
unsolved base solver instead of solved modes: the permittivity sampled grid is
wavelength-independent (non-dispersive prototype), so ref_wlen_um is
only a placeholder, run_eme_band rebases every solve to the swept
wavelength with VectorModeSolver.at_wavelength. Pass the same
dl_um / window_* for every section of a device so they share one
transverse grid (required by run_eme).
Scattering algebra¶
interface_smatrix
¶
interface_smatrix(
left_modes: Sequence[VectorMode],
right_modes: Sequence[VectorMode],
dl_x_um: float,
dl_y_um: float,
*,
rcond: float = 1e-10,
enforce_passivity: bool = False,
return_diagnostics: bool = False,
method: str = "galerkin",
) -> Union[SMatrix, Tuple[SMatrix, InterfaceDiagnostics]]
Generalized two-sided Petrov-Galerkin match at a waveguide step.
Independent E- and H-continuity projections are tested from both sides with
reaction-adjoint modal traces. The resulting tall system supports unequal
mode counts and complex/PML bases. Uniform node-collocated and validated
full-box uniform Yee-staggered bases are both supported, provided every mode
on both sides uses the same registration. (On the Yee lattice the two
reaction-product pairs Ex*Hy and Ey*Hx are natively co-located.)
A rank-revealing SVD exposes truncation through continuity residuals; it
never fabricates exact unitarity.
method chooses how the scattering matrix is built from those
projections:
"galerkin"(default) solves the tall continuity system in the least-squares sense. With a complete basis this is the exact interface. With a truncated basis it cannot place power in channels the basis lacks, so it returns that power to the represented modes: a guided-only basis reports almost no loss at a junction whose modes are offset or reshaped."projection"matches each incidence only on the far side's modes: from the left, the E- and H-continuity equations tested with the right modes, solved for the transmitted and reflected amplitudes together (for one mode per side, t = 2cb / (c + b) and r = (b - c) / (b + c), with c and b the E and H projections); from the right, the same with the left modes. With a complete basis this is the exact interface, like the default. With a truncated basis, power that the far basis cannot represent is lost from the cascade, as radiation leaves an open device: a guided-only basis then keeps the mode-overlap loss of a lateral offset or a straight-to-bend junction. It drops radiation that would return to the guide further on, so it reads low wherever that matters, such as a sharp curvature reversal, a closed box or a very short abrupt step; there a complete basis is needed. Identical bases give the identity in both methods.
enforce_passivity=True clips any gain of an ordinary real guided
interface (singular values of the power-normalized S-matrix above 1).
It is off by default: the guided traces are not exactly orthogonal, so
raw interfaces are slightly active (about 1e-3 on an adiabatic taper,
1e-2 at an abrupt step), and the clip does not remove only that gain.
It also pulls the fundamental channel down, by 0.1 % (num_modes=2)
to 0.4 % (num_modes=4) on an adiabatic taper at dl 0.04 um, a bias
that does not shrink with more modes or slabs (it does shrink with a
finer cross-section mesh), while the raw cascade stays lossless to 3e-4. The
raw violation is always reported in the diagnostics.
propagation_smatrix
¶
propagation_smatrix(
modes: Sequence[VectorMode], length_um: float
) -> SMatrix
Diagonal propagation S-matrix for a uniform section of length_um.
Forward and backward both pick up Φ = diag(exp(-i k0 n_eff L - k0 k_eff L))
with no inter-mode coupling; reflection blocks are zero. Equivalently this
is exp(-i*beta_eme*L) for
beta_eme = k0*mode.n_eff_complex.conjugate().
The generalized interface transforms only its dual test basis, so these
raw primal eigenmodes retain their diagonal propagation constants. For an
evanescent mode n_eff≈0, k_eff>0 the same expression is the causal decay.
cascade
¶
Cascade a left-to-right sequence of S-matrices with the star product.
star_product
¶
Redheffer star product Sa ⋆ Sb (Sa on the left, Sb on the
right), connecting Sa's right port to Sb's left port.
Sa may connect n0 -> n1 and Sb n1 -> n2; only their internal
port counts must match. Linear solves are used rather than explicit matrix
inverses.
Convergence¶
eme_convergence_report
¶
eme_convergence_report(
results: Sequence[EMEResult],
*,
labels: Optional[Sequence[str]] = None,
port_left: int = 1,
port_right: int = 1,
tolerance: float = 0.001,
stability_steps: int = 2,
) -> EMEConvergenceReport
Compare a nested sweep on one fixed propagating external-port subspace.
Radiation eigenvectors reorder when the window or PML changes, so comparing their individual columns is meaningless. Propagating eigenvectors can likewise phase-flip or rotate inside a degenerate subspace. On a common external-port grid this report least-squares aligns each current port basis to the previous one before comparing the four complex S blocks. A requested guided subset may not split a degenerate beta cluster.
energy_balances is the worst-case generalized output/input power ratio
over the selected left port subspace, not the coordinate-dependent result
of exciting whichever degenerate vector happened to be returned first.
Field-continuity (both all-channel and guided-input maxima), reaction
orthogonality, and generalized reciprocity are exposed separately; they are
not folded into smatrix_stable because useful thresholds depend on the
problem and retained channel types.
stability_steps prevents one accidental flat pair from being declared
stable. The default requires two consecutive sub-tolerance deltas, hence
at least three results. Use stability_steps=1 only for an exact
algebraic identity check, not a numerical convergence claim.
EMEConvergenceReport
dataclass
¶
EMEConvergenceReport(
labels: Tuple[str, ...],
section_mode_counts: Tuple[Tuple[int, ...], ...],
port_smatrix_deltas: Tuple[float, ...],
port_subspace_errors: Tuple[float, ...],
energy_balances: Tuple[float, ...],
max_interface_residuals: Tuple[float, ...],
max_reciprocity_errors: Tuple[float, ...],
tolerance: float,
stability_steps: int,
smatrix_stable: bool,
max_guided_interface_residuals: Tuple[float, ...] = (),
max_reaction_orthogonality_errors: Tuple[
float, ...
] = (),
)
Successive fixed-port stability evidence for EME basis/PML sweeps.
smatrix_stable means that the requested number of consecutive,
basis-aligned propagating-port S-matrix deltas and port-subspace errors are
within tolerance. It deliberately does not mean that the underlying
field expansion is physically converged. Interface continuity,
reciprocity, passivity/energy balance, mesh, window, and PML sweeps remain
independent gates.
InterfaceDiagnostics
dataclass
¶
InterfaceDiagnostics(
n_left: int,
n_right: int,
rank: int,
equations: int,
condition_number: float,
raw_residual_by_input: ndarray,
residual_by_input: ndarray,
raw_passivity_violation: float,
passivity_violation: float,
passivity_correction_applied: bool,
reciprocity_error: float,
left_gram_condition: float,
right_gram_condition: float,
left_gram_asymmetry: float,
right_gram_asymmetry: float,
left_reaction_orthogonality_error: float = 0.0,
right_reaction_orthogonality_error: float = 0.0,
max_guided_input_residual: float = 0.0,
method: str = "galerkin",
)
Numerical evidence for one generalized mode-matched interface.
*_reaction_orthogonality_error is the largest self-normalized reaction
overlap between modes with distinct propagation constants. It complements
Gram condition/asymmetry: a large symmetric off-diagonal overlap can expose
an inaccurate high-order E/H reconstruction while leaving both of those
older diagnostics deceptively modest. max_guided_input_residual filters
the continuity residual to incident guided channels on either side, so a
poor continuum test channel cannot hide the accuracy of the physical ports.
method records the interface model that built the S-matrix; rank,
condition_number and equations describe that model's own system.