Skip to content

Scattering API

scattering

Numerically-stable scattering-matrix backend.

The default engine multiplies per-layer transfer matrices Mᵢ = Vᵢ·Pᵢ·Vᵢ⁻¹ with Pᵢ = diag(exp(-i·kz·k₀·d)) (waves.get_matrix). For thick / lossy / strongly-evanescent layers the propagation diagonal mixes exponentially growing and decaying terms, so the products overflow and lose precision (energy conservation drifts, eventually NaN).

This module reformulates the same 4×4 anisotropic problem as a scattering matrix cascaded with the Redheffer star product, where every propagation factor is a decaying exponential (exp(-|Im(kz)|·k₀·d) ≤ 1) — the growing matrices are never materialised — so it is stable in all cases. It is the opt-in backend behind Structure.execute(payload, backend="scattering") and produces the same reflection/transmission coefficients as the transfer method where both are well-conditioned, and correct coefficients where the transfer method fails.

Algorithm and conventions follow PyLlama (Bay, Lafait & Lequime / Bay et al., J. Opt. Soc. Am. A 39, 1431 (2022); arXiv:2012.05945), eqs 18-28, adapted to this codebase's field-vector ordering [Ex, Ey, Hx, Hy] and partial-wave sort (WaveProfile.tangential_modes: columns [t0, t1, r0, r1] = forward, forward, backward, backward). Per-medium eigenvectors come from each layer's profile (interior layers and a semi-infinite crystal exit) or from the ambient dynamical matrices A_inc / A_exit (prism, isotropic exit).

Functions:

scattering_coefficients 

scattering_coefficients(layers, k_0)

Reflection/transmission coefficients via the stable scattering-matrix method.

Parameters:

  • layers (list) –

    The executed structure's layers (structure.layers), each exposing eigenvectors via profile or an ambient dynamical matrix.

  • k_0 (ndarray) –

    Canonical free-space wavenumber [1, 1, F, 1].

Returns:

  • dict[str, ndarray] –

    Dict of the eight complex coefficients ``r_pp, r_ss, r_ps, r_sp, t_pp,

  • dict[str, ndarray] –

    t_ss, t_ps, t_spin canonical[A, B, F, T]`` layout (un-presented).

Note

For an isotropic exit the transmission coefficients are in the clean s/p basis; for a semi-infinite anisotropic exit they are in that crystal's eigenmode basis (as in the transfer method), while the reflection coefficients are always in the prism's s/p basis.

Source code in hyperbolic_optics/scattering.py 
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
def scattering_coefficients(layers: list, k_0: np.ndarray) -> dict[str, np.ndarray]:
    """Reflection/transmission coefficients via the stable scattering-matrix method.

    Args:
        layers: The executed structure's layers (``structure.layers``), each
            exposing eigenvectors via ``profile`` or an ambient dynamical matrix.
        k_0: Canonical free-space wavenumber ``[1, 1, F, 1]``.

    Returns:
        Dict of the eight complex coefficients ``r_pp, r_ss, r_ps, r_sp, t_pp,
        t_ss, t_ps, t_sp`` in canonical ``[A, B, F, T]`` layout (un-presented).

    Note:
        For an isotropic exit the transmission coefficients are in the clean s/p
        basis; for a semi-infinite anisotropic exit they are in that crystal's
        eigenmode basis (as in the transfer method), while the reflection
        coefficients are always in the prism's s/p basis.
    """
    media = [_medium(layer) for layer in layers]
    scattering = None
    for (p_left, q_left, d_left), (p_right, _, _) in zip(media[:-1], media[1:], strict=True):
        interface = _interface_scattering(p_left, q_left, d_left, p_right, k_0)
        scattering = interface if scattering is None else _star(scattering, interface)

    assert_canonical(scattering, matrix_ndim=2, name="scattering_matrix")
    s = scattering
    # PyLlama r_kj is (out k, in j); this codebase uses r_{in->out}, so the cross
    # terms transpose (verified against the transfer backend).
    return {
        "r_pp": s[..., 2, 0],
        "r_ss": s[..., 3, 1],
        "r_ps": s[..., 3, 0],
        "r_sp": s[..., 2, 1],
        "t_pp": s[..., 0, 0],
        "t_ss": s[..., 1, 1],
        "t_ps": s[..., 1, 0],
        "t_sp": s[..., 0, 1],
    }