S-matrix¶
Drive each port in turn and assemble the complex scattering matrix, with reciprocity and passivity checks, Touchstone files and wavelength-interpolated circuit models.
Plan and run¶
plan_smatrix
¶
plan_smatrix(
simulation,
ports: Sequence[Union[SMatrixPort, Any]],
*,
freqs_hz: Sequence[float],
source_time: Optional[SourceTimeType] = None,
power_watts: float = 1.0,
launch: str = "auto",
drive: Optional[Sequence[str]] = None,
keep_monitors: bool = False,
) -> SMatrixPlan
Prepare the N driven-port simulations + port readers for a full S-matrix, without running anything.
Parameters¶
simulation:
The device: a complete Simulation
(structures, grid, boundaries, run controls). Its sources are
placeholders, each driven copy replaces them with that port's mode
launch, and its monitors are dropped unless keep_monitors.
ports:
The port declarations (unique names): photonhub.Port values,
in your frame and at the same place as in Simulation(ports=...),
whose outward side and mode window resolve as they do there; or
SMatrixPort values, which name their frame on a simulation
fitted with domain=. A port in the half a symmetry plane mirrors
away raises ValueError, where ports= reads it through its
image: use the image, or drop the symmetry. A port read from
sim.ports of this same simulation is in the wire frame and stands
as it is, also after a rename or an edit of its centre, which is in
that frame too. Such an edit keeps the port's out_direction: moved
across the device, a port needs out_direction=None. A new
photonhub.Port built from its fields is in your frame.
freqs_hz:
Monitor frequencies for every port, the S-matrix's frequency axis.
source_time:
Drive pulse; None auto-tunes a
for_band over
freqs_hz.
power_watts / launch:
Forwarded to mode_launch.
drive:
Port names to actually drive (default: all). Undriven ports still get
monitors in every run (their rows fill; their columns stay NaN).
keep_monitors:
Also carry the base simulation's own monitors into every driven copy.
Returns¶
SMatrixPlan
With one prepared simulation per driven port and one
SPort reader per port.
SMatrixPlan
dataclass
¶
SMatrixPlan(
base: Any,
ports: Tuple[SMatrixPort, ...],
freqs_hz: Tuple[float, ...],
source_time: SourceTimeType,
simulations: Mapping[str, Any],
sports: Tuple[SPort, ...],
)
A fully-prepared S-matrix extraction: the driven simulations, the port
readers, and the run method. Built by plan_smatrix; inspect
simulations / sports (or quote()) before committing
compute, then call run. ports are the ports as placed, in the
wire frame on a simulation fitted with domain=.
run
¶
run(
runner: Union[
str, Callable[[Dict[str, Any]], Mapping[str, Any]]
] = "local",
*,
path_dir=None,
max_workers: int = 1,
colocate: bool = True,
allow_partial: bool = False,
**runner_kwargs: Any,
) -> SMatrixResult
Run every driven-port simulation and assemble the S-matrix.
Parameters¶
runner:
"local" (default) runs through
Batch (path_dir /
max_workers / solver_path / timeout / progress
forwarded); "cloud" (or the deprecated "web") submits via ph.cloud.Batch
(max_workers and extra kwargs such as device= and solver= forwarded,
with the cloud defaults when omitted: the service's only GPU and
the installed solver's commit, resolved once for all drives;
path_dir is not a web concept and must be None). Any other
callable is invoked as runner(simulations_dict) and must
return a mapping driven-port name -> RunResult (or any
name -> DataArray mapping per run), the injection point for
custom backends and tests.
colocate:
Forwarded to smatrix, keep True
for real (Yee-staggered) engine output.
allow_partial:
When some drives fail: False (default) raises with the per-drive
errors; True assembles the surviving columns (missing entries NaN)
and records the failures on the result.
run_smatrix
¶
run_smatrix(
simulation,
ports: Sequence[Union[SMatrixPort, Any]],
*,
freqs_hz: Sequence[float],
source_time: Optional[SourceTimeType] = None,
power_watts: float = 1.0,
launch: str = "auto",
drive: Optional[Sequence[str]] = None,
keep_monitors: bool = False,
runner: Union[
str, Callable[[Dict[str, Any]], Mapping[str, Any]]
] = "local",
path_dir=None,
max_workers: int = 1,
colocate: bool = True,
allow_partial: bool = False,
**runner_kwargs: Any,
) -> SMatrixResult
One call: plan_smatrix then SMatrixPlan.run. See both
for the parameters. For anything beyond a straight shot (cost preview,
inspecting the generated sims, custom per-drive tweaks) build the plan
explicitly.
SMatrixResult
dataclass
¶
SMatrixResult(
S: Any,
ports: Tuple[SMatrixPort, ...],
data: Mapping[str, Any],
errors: Mapping[str, Exception],
simulations: Mapping[str, Any],
metadata: Mapping[str, Any] = dict(),
)
The assembled S-matrix plus everything needed to audit it.
S is the complex xarray.DataArray from
assemble_smatrix, dims
(port_out, port_in, f), |S_ij|^2 a power ratio. data maps each
DRIVEN port name to its run's RunResult; errors carries
per-drive failures when allow_partial=True let the assembly proceed
(those columns are NaN).
On a simulation fitted with domain=, ports are the ports as placed,
in the wire frame, from the fitted box's low corner, while the field arrays
in data read in your frame, as any RunResult does.
metadata
class-attribute
instance-attribute
¶
reciprocity_error
¶
Maximum finite abs(S_ij - S_ji) across ports and frequencies.
Ignore non-finite differences. Return zero when none remain.
passivity_violation
¶
Largest eigenvalue of S.conj().T @ S minus one, across frequencies.
Skip matrices with non-finite entries. Return zero when no finite frequency matrix is available. A non-positive value passes passivity.
is_reciprocal
¶
True when reciprocity_error is at most atol.
is_passive
¶
True when passivity_violation is at most atol.
to_touchstone
¶
Write the matrix as a Touchstone file, see
write_touchstone.
to_model
¶
to_model(**kwargs) -> SMatrixModel
Build a checked wavelength-interpolated S-matrix model.
Ports¶
SMatrixPort
dataclass
¶
SMatrixPort(
name: str,
axis: str,
position_um: float,
out_direction: str,
half_w_um: float,
half_v_um: float,
polarization: str = "TE",
mode_index: int = 0,
center_um: Optional[Tuple[float, float]] = None,
dl_um: Optional[float] = None,
source_offset_um: Optional[float] = None,
supersample: int = 8,
num_modes: Optional[int] = None,
frame: Optional[str] = None,
)
One S-matrix port declaration: a plane, an outgoing direction, and a single mode channel.
Parameters¶
name:
Port label, the S-matrix index, the monitor name, and the batch entry
key (so it must be filesystem-safe: letters, digits, -_.).
axis:
Propagation axis of the port waveguide ('x' | 'y' | 'z').
position_um:
The port plane along axis, where the port monitor sits and where
S is referenced, in the frame frame names.
out_direction:
'+' or '-': the direction along axis that points OUT of the device
through this port (toward the nearby domain wall). The drive launches
the opposite way.
half_w_um / half_v_um:
Half-extents of the mode-solve window in the plane's natural
(horizontal, vertical) in-plane axes, same meaning as
solve_yee_mode. Choose them wide
enough that the guided mode's evanescent tail dies inside the window.
polarization / mode_index:
The mode channel, counted WITHIN the TE/TM family ('TE', 0 = TE0).
center_um:
Transverse waveguide location as (horizontal, vertical) in-plane
coordinates ((y, z) for an x-normal port, (x, z) for y-normal, (x, y)
for z-normal), in the frame frame names; None = the domain
centre. The port's mode is solved, read and launched at this centre.
dl_um:
Mode-solve transverse step; None = the simulation grid's dl_um.
source_offset_um:
Distance from the port plane to the drive's launch plane, measured
toward the wall (out_direction). None = a third of a
wavelength in the background, the standoff the declarative ports take.
supersample / num_modes:
Forwarded to the Yee mode solve (eigensolver frame controls).
frame:
What position_um and center_um are measured from on a
simulation fitted with domain=, whose stored coordinates start at
origin_um. "user" is your own frame, the one a
photonhub.Port reads in; the plan takes origin_um off.
"wire" is the fitted box's low corner, the frame of
solve_yee_mode, mode_monitor and mode_launch; the values
stand. None, the default, is for a simulation at the origin
(0, 0, 0), where the two are one place; on any other,
plan_smatrix raises ValueError.
in_direction
property
¶
The launch direction for this port's drive, into the device.
SPort
dataclass
¶
SPort(
name: str,
monitor: ModeMonitor,
out_direction: Optional[str] = None,
)
One port of an S-matrix: a ModeMonitor + its outgoing
direction.
Parameters¶
name:
Port label (the S-matrix index). Must be unique among the ports.
monitor:
The ModeMonitor recording the
port plane (its mode, axis, center_um, thickness_axis and
modes_by_freq drive the overlap).
out_direction:
The sign along the monitor's propagation axis ('+' or '-') that points
out of the device through this port, the direction a transmitted /
scattered wave leaves. The incident (incoming) direction is the opposite.
Defaults to the monitor's own direction (treated as the outgoing one).
in_direction
property
¶
The incoming (toward-device) direction, opposite the outgoing one.
outgoing
¶
{freq: b}, the complex amplitude of the wave scattered out of
the device through this port.
incoming
¶
{freq: a}, the complex amplitude of the wave incident into the
device at this port (the source side at the driven port).
Assemble¶
smatrix
¶
Multiport scattering-matrix (S-matrix) assembler, pure Python, no engine.
This assembles the multi-port scattering matrix, built entirely on top of the
mode-resolved complex modal amplitude that photonhub.analysis.mode_overlap.mode_amplitude extracts
from a recorded DFT plane. Where mode_devices gives you single-mode power
T(f) = |c|^2 per monitor, this module keeps the phase of c and arranges
the per-port amplitudes into a proper scattering matrix.
S-parameter definition¶
For an N-port device, mode the field on each port plane as a sum of an inward
(toward the device) and an outward (away from the device) guided wave. The
S-matrix relates the outgoing modal amplitudes b to the incoming ones a:
b = S a , S_ij(f) = b_i(f) / a_j(f) with port j the only one driven.
a_j is the amplitude incident on the device at the driven port; b_i
is the amplitude scattered out of the device at port i. Each is the
power-normalized wave amplitude c * sqrt(|P_mode| * 1e-12) (µm²→m²,
the engine's SI flux units), built from the
normalized modal coefficient c = a_pm / P_mode of
mode_amplitude and the port mode's own
plane power P_mode, so
|a|^2/|b|^2is the modal power the wave carries through the port plane (the same numbermode_powerreports), hence|S_ij|^2 = |b_i|^2 / |a_j|^2is a true power transmission / reflection even when ports i and j carry different modes (unequal widths / polarizations, whose per-modeP_modediffer) and a lossless reciprocal device hasS† S = I. (The barecratio would be biased byP_mode,j / P_mode,i ≈ n_eff,j / n_eff,iacross unequal ports and would break reciprocity there; for same-mode ports the two coincide.)- the phase is
c's phase (sqrt(|P_mode|)is real positive), soScarries the physical propagation phase a mode accumulates along the guide.
Incident vs scattered: the directional projection¶
Each SPort carries the outgoing direction, the sign along the
monitor's propagation axis that points out of the device through that port
("+" or "-"). The two directional projections of the recorded plane are
then:
- outgoing / scattered
b_i= project onto the mode travelling in the port's outgoing direction; - incoming / incident
a_i= project onto the mode travelling the opposite way (into the device).
At the driven port the recorded plane (placed on the total-field side, just
inside the source) carries incident + reflected: its incoming projection is the
incident a_j and its outgoing projection is the reflected b_j (so
S_jj is the port reflection). At every other port only scattered light
leaves, so the outgoing projection is the transmitted b_i and S_ij is the
port-i transmission. This directional split is exactly what the forward/backward
overlap of mode_amplitude provides.
Normalization¶
a and b are built from the source-spectrum-normalized complex modal
amplitudes already living in the run data: the recorded phasors are divided by
A0*S(f), so that normalization cancels in the ratio
S_ij = b_i / a_j and leaves a dimensionless scattering parameter. Each port's
amplitude is the power-normalized c * sqrt(|P_mode| * 1e-12) (see above, the µm²→m² factor keeps |b|^2 equal to the flux-commensurate
ModeMonitor.mode_power), so
|S_ij|^2 is a power ratio for same-mode AND unequal-mode port pairs alike
(a clean straight through-guide reads |S21| ≈ 1). No de-embedding of the
source-to-monitor or monitor-to-port reference plane is applied, S is
referenced to the monitor planes as placed (see "What's not handled").
Output¶
smatrix returns one column of S (the column for the driven port) as
{(port_out, port_in): {freq_hz: S}}, i.e. S_ij for every port i with
j fixed to the driven port. Run it once per driven port and feed the columns
to assemble_smatrix to get the full matrix as an xarray.DataArray
indexed (port_out, port_in, f). reciprocity_error /
passivity_violation (and their assert_* / is_* wrappers) check the
assembled matrix.
What's not handled (deferred)¶
- Multimode ports. One mode per port (the monitor's
mode). A multimode S-matrix would index ports by (port, mode_index); the per-mode amplitudes are available by re-projecting onto each mode, but the bookkeeping here is single-mode. - De-embedding. S is referenced to the monitor planes; no phase de-embedding
back to a common reference plane is performed. Do it externally if needed: in
the engine's
e^{-i omega t}convention an outgoing wave accumulatese^{+i beta L}over the extra lead lengthL, so multiply bye^{-i beta L}per port to remove it. - Source-overlap correction at the driven port. The incident
a_jis read from the recorded total field's incoming projection, which assumes the launched mode matches the monitor'smode; a mode-mismatched launch would biasa_j.
SPort
dataclass
¶
SPort(
name: str,
monitor: ModeMonitor,
out_direction: Optional[str] = None,
)
One port of an S-matrix: a ModeMonitor + its outgoing
direction.
Parameters¶
name:
Port label (the S-matrix index). Must be unique among the ports.
monitor:
The ModeMonitor recording the
port plane (its mode, axis, center_um, thickness_axis and
modes_by_freq drive the overlap).
out_direction:
The sign along the monitor's propagation axis ('+' or '-') that points
out of the device through this port, the direction a transmitted /
scattered wave leaves. The incident (incoming) direction is the opposite.
Defaults to the monitor's own direction (treated as the outgoing one).
in_direction
property
¶
The incoming (toward-device) direction, opposite the outgoing one.
outgoing
¶
{freq: b}, the complex amplitude of the wave scattered out of
the device through this port.
incoming
¶
{freq: a}, the complex amplitude of the wave incident into the
device at this port (the source side at the driven port).
smatrix
¶
smatrix(
ports: Sequence[SPort],
driven,
data,
*,
colocate: bool = True,
) -> Dict[Tuple[str, str], Dict[float, complex]]
Assemble one column of the S-matrix from a single driven-port run.
Drives port driven (only) and returns S_ij(f) = b_i / a_j for every
port i in ports, with j = driven fixed, i.e. the column of the
full S-matrix belonging to the driven port. a_j is the amplitude incident
on the device at the driven port (its incoming projection); b_i is the
amplitude scattered out of the device at port i (its outgoing projection).
Parameters¶
ports:
The N SPorts. Names must be unique.
driven:
The driven port, an SPort or its name. Must be in ports.
data:
The run's RunResult (or any mapping name -> DataArray) for the
run that drove driven. Every port monitor's data must be present.
Returns¶
dict[(port_out, port_in), dict[float, complex]]
{(i, driven): {freq_hz: S_ij}} for each port i. S is the
complex scattering parameter (|S|^2 is a power ratio).
assemble_smatrix
¶
assemble_smatrix(
columns: Sequence[
Mapping[Tuple[str, str], Mapping[float, complex]]
],
*,
port_order: Optional[Sequence[str]] = None,
)
Combine per-driven-port smatrix columns into a full S-matrix as an
xarray.DataArray indexed (port_out, port_in, f).
Parameters¶
columns:
An iterable of the dicts returned by smatrix (one per driven
port). Together they should cover every (i, j) entry; missing entries are
filled with nan and a clear coord layout is still returned.
port_order:
Optional explicit ordering of the port labels for both axes. Defaults to
first-seen order across the columns (driven ports first, then any extra
out-ports).
Returns¶
xarray.DataArray
Complex S with dims ("port_out", "port_in", "f") and string
port coords. S.sel(port_out=i, port_in=j) is S_ij(f).
reciprocity_error
¶
Max |S_ij - S_ji| over all off-diagonal entries and frequencies.
A reciprocal (non-magnetic, non-gyrotropic) device has a symmetric S-matrix,
so this should be ~0. nan entries (uncomputed) are ignored.
is_reciprocal
¶
True when reciprocity_error is within atol.
assert_reciprocal
¶
Raise AssertionError if S is not symmetric to atol.
passivity_violation
¶
Worst passivity excess max(eig(S† S)) - 1 over all frequencies.
For a passive device the scattering operator does not create energy:
S† S <= I, i.e. every eigenvalue of the Hermitian Gram matrix S† S is
<= 1. This returns the largest eigenvalue minus one (<= 0 means
passive; a small positive value is numerical slack). Frequencies with any
nan S entry are skipped.
is_passive
¶
True when passivity_violation is within atol (so the
largest singular value of S is <= 1 + atol).
assert_passive
¶
Raise AssertionError if S amplifies energy beyond atol (i.e.
max eig(S† S) > 1 + atol).
assemble_smatrix
¶
assemble_smatrix(
columns: Sequence[
Mapping[Tuple[str, str], Mapping[float, complex]]
],
*,
port_order: Optional[Sequence[str]] = None,
)
Combine per-driven-port smatrix columns into a full S-matrix as an
xarray.DataArray indexed (port_out, port_in, f).
Parameters¶
columns:
An iterable of the dicts returned by smatrix (one per driven
port). Together they should cover every (i, j) entry; missing entries are
filled with nan and a clear coord layout is still returned.
port_order:
Optional explicit ordering of the port labels for both axes. Defaults to
first-seen order across the columns (driven ports first, then any extra
out-ports).
Returns¶
xarray.DataArray
Complex S with dims ("port_out", "port_in", "f") and string
port coords. S.sel(port_out=i, port_in=j) is S_ij(f).
Physical checks¶
is_reciprocal
¶
True when reciprocity_error is within atol.
assert_reciprocal
¶
Raise AssertionError if S is not symmetric to atol.
reciprocity_error
¶
Max |S_ij - S_ji| over all off-diagonal entries and frequencies.
A reciprocal (non-magnetic, non-gyrotropic) device has a symmetric S-matrix,
so this should be ~0. nan entries (uncomputed) are ignored.
is_passive
¶
True when passivity_violation is within atol (so the
largest singular value of S is <= 1 + atol).
assert_passive
¶
Raise AssertionError if S amplifies energy beyond atol (i.e.
max eig(S† S) > 1 + atol).
passivity_violation
¶
Worst passivity excess max(eig(S† S)) - 1 over all frequencies.
For a passive device the scattering operator does not create energy:
S† S <= I, i.e. every eigenvalue of the Hermitian Gram matrix S† S is
<= 1. This returns the largest eigenvalue minus one (<= 0 means
passive; a small positive value is numerical slack). Frequencies with any
nan S entry are skipped.
Export¶
write_touchstone
¶
write_touchstone(
S,
path,
*,
z0: float = 50.0,
comments: Sequence[str] = (),
version: int = 1,
metadata: Optional[Mapping[str, Any]] = None,
strict: bool = False,
passivity_atol: float = 0.001,
reciprocity_atol: float = 0.001,
) -> Path
Write an assembled S-matrix as a Touchstone file.
S is the xarray.DataArray from
assemble_smatrix /
SMatrixResult.S, dims (port_out, port_in, f) with matching
port labels on both axes. Frequencies are written ascending in Hz,
real/imaginary (RI) format, one reference resistance z0 (the
Touchstone header requires one; it is dimensional bookkeeping only for
these modal wave amplitudes, power waves are already normalized so
|S_ij|^2 is a power ratio).
Version 1 layout follows the v1 spec: 1-port and 3+-port matrices in row-major
(port_out-major) order with at most four S-pairs per line and each
matrix row starting a new line; the 2-port special column-major order
S11 S21 S12 S22 on one line. NaN entries (an unassembled partial
matrix) are rejected.
Version 1's file extension is forced to the port count: device or
device.s3p for a 3-port; a mismatched suffix raises. Version 2 uses
.ts by default and records port count in its header. Failed checks
warn and are recorded, or raise before writing when strict=True.
Port labels, comments and metadata must contain ASCII text.
read_touchstone
¶
Read a Touchstone v1 or v2 S-parameter file.
Parameters¶
path:
Touchstone .sNp or .ts file. RI, MA and DB data and Hz, kHz,
MHz and GHz frequencies are accepted. Version 2 noise data is ignored.
Returns¶
xarray.DataArray
Complex matrix with dimensions (port_out, port_in, f) and frequency
coordinates in Hz.
smatrix_checks
¶
smatrix_checks(
S,
*,
passivity_atol: float = 0.001,
reciprocity_atol: float = 0.001,
) -> SMatrixChecks
Measure passivity and reciprocity with nonnegative finite tolerances.
Parameters¶
S:
S-matrix with axes (port_out, port_in, f) or an array accepted by
the S-matrix diagnostics.
passivity_atol, reciprocity_atol:
Allowed excess for each consistency check.
SMatrixChecks
dataclass
¶
SMatrixChecks(
passivity_violation: float,
passivity_atol: float,
passive: bool,
reciprocity_error: float,
reciprocity_atol: float,
reciprocal: bool,
num_freqs: int,
)
Passivity and reciprocity measurements; num_freqs counts samples.
Circuit models¶
smatrix_model
¶
smatrix_model(
S,
*,
metadata: Mapping[str, Any] | None = None,
strict: bool = False,
passivity_atol: float = 0.001,
reciprocity_atol: float = 0.001,
) -> SMatrixModel
Build a checked wavelength-interpolated model from S-parameters.
Parameters¶
S: Complex S-matrix with matching output and input port labels. metadata: JSON-serializable descriptive values with string keys. strict: Raise if passivity or reciprocity exceeds its tolerance. passivity_atol, reciprocity_atol: Nonnegative consistency check tolerances.
Polar interpolation is the default. An interval switches to whole-matrix
complex interpolation only when its polar excess exceeds the largest
sample gain (at least zero) plus the tolerance and a roundoff margin.
Noisy or uniformly gained samples keep polar and are reported by the
checks, whose tolerance is unchanged. Construction
rechecks the selected interpolant at 15 interior points (k/16, k=1..15),
in addition to the samples. Finite probes can miss unsampled polar peaks;
they do not certify continuous passivity or accuracy. Strict mode refuses
failures; non-strict mode warns and records them. Every switch warns with
its worst significant midpoint loss in dB: passivity is not accuracy.
See SMatrixModel for warning thresholds and sampling limitations.
Plain NumPy calls work in SAX circuits, but not with jax.jit.
SMatrixModel
dataclass
¶
SMatrixModel(
ports: tuple[str, ...],
wlens_um: ndarray,
s: ndarray,
checks: SMatrixChecks,
metadata: Mapping[str, Any],
max_phase_step_rad: float,
complex_intervals: tuple[int, ...] = (),
interpolation: Mapping[str, Any] = dict(),
)
A sampled S-matrix interpolated in wavelength for circuit use.
s uses [port_out, port_in, wavelength] order. Callable dictionary
keys use (port_from, port_to) order.
Polar (linear magnitude and unwrapped phase) interpolation is the default. At build time, an interval switches to whole-matrix complex interpolation only if a polar interior probe exceeds the largest sample gain (at least zero) plus the passivity tolerance and a roundoff margin. Samples and 15 interior points (k/16, k=1..15) are checked. These finite probes are diagnostics, not a bound on unsampled polar peaks or an accuracy certificate. Saved choices are used on every call without rechecking.
Any switch warns once per model with the count, worst interval bounds
and worst significant-entry midpoint amplitude drop in dB versus polar.
Significant means at least 10% of the larger row/column maximum there.
Phase steps above pi/2 warn when both endpoint magnitudes exceed 1e-6.
Sampling more densely resolves rapid phase evolution,
but cannot remove a zero's phase jump or detect aliased whole turns.
Plain NumPy calls work in SAX circuits, but not with jax.jit.