From 1eb0ec335381a38b88b281dea9bc5b866b9076a1 Mon Sep 17 00:00:00 2001 From: Yongha Kim Date: Wed, 4 Nov 2026 00:08:53 +0900 Subject: [PATCH 1/2] feat(rcwa): a polarization word takes the pure 1D TE/TM path A public rcwa.run of a one-row raster grating lit by PlaneWave(polarization='s') or 'p' now solves on the engine's pure TE/TM path (grating_type 0) instead of the conical one. No public field is added. The routing principle stays -- the declaration selects the path, and nothing else does -- and what counts as a declaration widens. The words 's' / 'p' declare the polarization fixed. A Jones pair, (1, 0) included, is a value and stays conical. The word takes the pure path only where that path answers the same light. Every test is structural, so jit and vmap decide it the same way: - every layer a one-row scalar raster on an axis-aligned cell (the routing's own one_dimensional test), with no tensor; - phi and ky absent, or a zero that nothing traces or differentiates; - the formulation is not first_order. Everywhere else the word solves bit for bit as the pair it names, silently. Unlike tetm it refuses nothing. Engine: - polarization accepts the words. The pair is stored as before, and the word in _polarization_word beside it. - _assign_grating_type keeps its refusals and assigns _routed_grating_type, which refuses and warns nothing, so a report can read it before any solve. - _polarization_word_holds is the eligibility test, and tetm_index reads a routed word. - _solve_path(planned=) names the path. - The JAX pytree carries the word in aux. - phi, kx and ky are stored under _eager(), as the material setters are. Without that a stated phi=0. was staged into a tracer inside jax.jit, and the jitted run of an eligible problem took the conical path (measured). The eligibility comparison runs under _eager() for the same reason. Adapter: the word goes to the engine as written, except under Outputs(scattering=True). The public matrix is 2N x 2N with channel 'sp', and the pure path's is N x N 'te'/'tm', so that run hands over the pair as before. The path is recorded: CheckReport.resolved['solve_path'] (planned) and Result.diagnostics['solve_path'] (per wavelength, as taken). The values are '1d_te', '1d_tm', '1d_conical' or '2d'. Nothing recorded the path before. Coherent-source runs report '2d', the layout their side-resolved first-order cascade is written on (CoherentSolver._solve_path). Measured through the public rcwa.run (NumPy, a 10-core laptop; docs/guide/rcwa/light.md "Cost"): 2.7x, 4.1x and 5.3x faster at orders (50, 0), (200, 0) and (800, 0), with a peak of 0.39 GB against 5.9 GB at the last; sum(R) 5e-13, 8e-12 and 4e-10 apart. Checks: - tests/api/test_unified_rcwa_pure_1d.py (new): the regression matrix holds the word against the pair to 1e-11 of each observable family; the worst measured gap is 1.9e-13. It also covers jax/torch backends, single precision, jax.grad, torch autograd and prepare -> jit / vmap / jit(vmap(grad)), with a kernel call counter. Fallback pins (skewed cell, nonzero phi and ky, first_order, scattering=True, traced or differentiated phi) are bit for bit against the pair. - The fallback cases ran on origin/main 93cf7cb1 and on this commit, through the same script: 184 arrays, all bit-identical. - tests/rcwa/test_pytree_round_trip.py lists the new aux key and checks that the word survives the round trip. Not fixed here: jax.jit of a public run whose Jones pair is exactly (0, 1) aborts the process. It does so on main too, under jaxlib 0.11.0: XLA's algebraic simplifier reads the runtime guard's complex constant as an iota and rewrites its abs. The new test writes that pair as (0, 2). Co-Authored-By: Claude Opus 5.5 --- src/wavebit/manifest.py | 5 +- src/wavebit/rcwa/_adapt.py | 19 +- src/wavebit/rcwa/_backends/jax_backend.py | 7 + src/wavebit/rcwa/_engine/engine.py | 203 +++++++++-- src/wavebit/rcwa/_engine/illumination.py | 8 + src/wavebit/rcwa/_lifecycle.py | 6 + src/wavebit/rcwa/results.py | 6 + tests/README.md | 1 + tests/api/test_unified_rcwa_pure_1d.py | 416 ++++++++++++++++++++++ tests/rcwa/test_pytree_round_trip.py | 21 +- 10 files changed, 656 insertions(+), 36 deletions(-) create mode 100644 tests/api/test_unified_rcwa_pure_1d.py diff --git a/src/wavebit/manifest.py b/src/wavebit/manifest.py index 1a780747..e27b9dd1 100644 --- a/src/wavebit/manifest.py +++ b/src/wavebit/manifest.py @@ -585,7 +585,10 @@ def _raster_normal_methods(): 'relative phase psi cannot spell — (1, -1j)/sqrt(2) and ' '(1, +1j)/sqrt(2) are the circular states. The pair is ' 'normalized to unit incident flux (overall scale is free, ' - 'overall phase is kept)', + "overall phase is kept). The words 's' and 'p' are the pairs " + '(1, 0) and (0, 1) declared fixed: the pure 1D TE/TM path solves ' + 'them on a one-row raster in an axis-aligned cell with phi and ky ' + "absent or a zero nothing traces, and the pair's path everywhere else", }, 'incidence_side': { 'values': ('top', 'bottom'), diff --git a/src/wavebit/rcwa/_adapt.py b/src/wavebit/rcwa/_adapt.py index b00c56ef..f500c8f5 100644 --- a/src/wavebit/rcwa/_adapt.py +++ b/src/wavebit/rcwa/_adapt.py @@ -177,6 +177,22 @@ def polarization(source, xp, execution): return pair / norm +def plane_wave_polarization(problem, source, xp, execution): + """What the engine is handed as one plane wave's polarization: the word, or the pair. + + A word, ``'s'`` or ``'p'``, goes over as written. The engine reads it as a fixed polarization + and decides from the structure and the incidence whether its pure 1D TE/TM path answers it + (``CoreSolveMixin._polarization_word_holds``); where that path does not, the word solves as + the pair it names. The one fact the engine cannot see is what the run asks for: a requested + scattering matrix keeps both polarizations' ports (``2N x 2N``, channel ``'sp'``), while the + pure path's has one polarization's (``N x N``, ``'te'`` / ``'tm'``). So a run that requests + one hands over the pair, and solves on the pair's path. A Jones pair is always a pair. + """ + if isinstance(source.polarization, str) and not problem.outputs.scattering: + return source.polarization + return polarization(source, xp, execution) + + def single_plane(problem): return len(problem.excitation.sources) == 1 and type(problem.excitation.sources[0]) is PlaneWave @@ -292,7 +308,8 @@ def build_solver(problem, numerics, execution, wavelength): kwargs['thickness'] = _array(xp, [0.], execution) if single_plane(problem): source = problem.excitation.sources[0] - kwargs.update(incidence_side=source.side, polarization=polarization(source, xp, execution)) + kwargs.update(incidence_side=source.side, + polarization=plane_wave_polarization(problem, source, xp, execution)) kwargs.update({k: getattr(source, k) for k in ('theta', 'phi', 'kx', 'ky') if getattr(source, k) is not None}) else: adapted = tuple(_engine_source(source, xp, execution) for source in problem.excitation.sources) diff --git a/src/wavebit/rcwa/_backends/jax_backend.py b/src/wavebit/rcwa/_backends/jax_backend.py index 5a7bd942..ff8254bf 100644 --- a/src/wavebit/rcwa/_backends/jax_backend.py +++ b/src/wavebit/rcwa/_backends/jax_backend.py @@ -795,6 +795,10 @@ def _tree_flatten(self): # Not the pair itself (that is a child now) but the fact that `psi` was given # explicitly, which the always-present psi child would otherwise overwrite. 'psi_given': self._psi_given, + # The word a pair was written as, `'s'` / `'p'` or `None`. A declaration, so static: + # the routing branches on it (`_polarization_word_holds`), and the pair child alone + # would come back a value, whose solve takes the conical path. + 'polarization_word': self._polarization_word, # `normal_grid` was missing beside the cell basis, and both were silent when they # were: the compiled `solve_field` flattens `self` on every call, so the solve ran # on whatever the *rebuilt* solver carries, not on what the caller configured. @@ -825,6 +829,8 @@ def _tree_unflatten(cls, aux_data, children): # legal when it was flattened. Rebuild without it, then restore the slot directly. aux_data = dict(aux_data) psi_given = aux_data.pop('psi_given') + # Restored beside the pair it names, which does not ride through `__init__` either. + polarization_word = aux_data.pop('polarization_word') # `formulation` is settable on the instance and on `solve()`, but not on `rcwa()`, # so it cannot go through the constructor either. Same treatment as `polarization`. formulation = aux_data.pop('formulation') @@ -862,6 +868,7 @@ def _tree_unflatten(cls, aux_data, children): # guards that with a try/except, but there is nothing left to validate here anyway -- # this pair was validated when the solver it came from was built. solver._polarization = polarization + solver._polarization_word = polarization_word solver._psi = psi solver._psi_given = psi_given solver._structure_drawn_on = drawn_on diff --git a/src/wavebit/rcwa/_engine/engine.py b/src/wavebit/rcwa/_engine/engine.py index 9065fc49..6da8319e 100644 --- a/src/wavebit/rcwa/_engine/engine.py +++ b/src/wavebit/rcwa/_engine/engine.py @@ -336,6 +336,26 @@ def _holds_autodiff(value): return False +def _absent_or_fixed_zero(value, eager): + """Whether an incidence component is absent, or a zero that nothing traces or differentiates. + + What a polarization word reads ``phi`` and ``ky`` with before it takes the pure 1D path + (:meth:`CoreSolveMixin._polarization_word_holds`). The test is structural before it is a + read: a JAX tracer -- ``jax.jit``, ``vmap`` and an eager ``jax.grad`` alike -- and a PyTorch + tensor under autodiff are an azimuth being swept or differentiated, which only the conical + path carries, so neither is read at all. A value with nothing attached is compared with zero + under ``eager`` (the solver's :meth:`~CoreSolveMixin._eager`), since inside ``jax.jit`` the + comparison of a concrete value is otherwise staged and has no truth value; one that cannot + be read (a ``torch.func.vmap`` batch) is not a fixed zero either. + """ + if value is None: + return True + if _holds_autodiff(value): + return False + with eager(): + return _concrete(value == 0) is True + + #: The text of :class:`~wavebit._common.backend.GainCutoffWarning`, constant apart from the side @@ -1330,6 +1350,9 @@ def _solve_fn(self, grating_type): _frozen_geometry = None _lattice_basis = None _polarization = None + #: ``'s'`` / ``'p'`` when :attr:`polarization` was written as that word, else ``None``. The + #: pair it names is in ``_polarization`` either way; this is the declaration beside it. + _polarization_word = None _tetm = None _chi = None #: ``None`` is a real state now, not an uninitialised one: it is what a solver whose @@ -1492,8 +1515,16 @@ def tetm_index(self): The core takes the index rather than the word because it indexes into the two formulations rather than naming them; the word is the user's spelling. + + Declared by :attr:`tetm`, or by a polarization word whose solve the routing put on the + pure path: ``'s'`` is TE and ``'p'`` TM, the same light (:attr:`polarization`). A word + routed elsewhere answers ``None``, as a pair does, since its solve carries both. """ - return None if self._tetm is None else int(self._tetm == 'tm') + if self._tetm is not None: + return int(self._tetm == 'tm') + if self._grating_type == 0 and self._polarization_word is not None: + return int(self._polarization_word == 'p') + return None def _refuse_y_orders_on_a_y_invariant_cell(self, one_dimensional): """A one-row raster cannot diffract out of plane, so it cannot keep ``y`` orders. @@ -1668,6 +1699,22 @@ def polarization(self): error, not an override. A 1D structure given a pair is solved on the conical path, which carries both polarizations; the pure 1D TE/TM path is single-polarization by construction. + + **The words ``'s'`` and ``'p'`` are a declaration, and a pair is a value** (since + 2026-11-03). A word reads back as the pair it names, ``(1, 0)`` or ``(0, 1)``, and says + beside it that the polarization is fixed: not a parameter anyone sweeps or + differentiates. That is what the pure 1D TE/TM path needs, so a word takes it wherever it + answers the same light -- a one-row raster on an axis-aligned cell, ``phi`` and ``ky`` + absent or a zero nothing traces (:meth:`_polarization_word_holds`) -- and everywhere else + solves exactly as the pair it names. Unlike :attr:`tetm` it refuses nothing: it states + the light, not the path. The two paths agree to round-off (2e-13 relative over + ``tests/api/test_unified_rcwa_pure_1d.py``), and through the public ``rcwa.run`` the pure + one was 2.7, 4.1 and 5.3 times faster at ``truncation_order`` ``(50, 0)``, ``(200, 0)`` + and ``(800, 0)``, at a peak of 0.39 GB against 5.9 GB at the last (measured 2026-11-03, + ``docs/guide/rcwa/light.md``). The public ``rcwa.PlaneWave`` hands + its own ``'s'`` / ``'p'`` here, unless a scattering matrix is asked for: the pure path's + is one polarization's, ``N x N`` with channel ``'te'`` / ``'tm'``, and the public one + keeps both, ``2N x 2N`` with channel ``'sp'``. """ return self._polarization @@ -1675,6 +1722,7 @@ def polarization(self): def polarization(self, polarization): if polarization is None: self._polarization = None + self._polarization_word = None return if self._psi_given: raise ValueError(_PSI_POLARIZATION_CLASH) @@ -1682,6 +1730,16 @@ def polarization(self, polarization): raise ValueError(_CHI_POLARIZATION_CLASH) if self._tetm is not None: raise ValueError(_TETM_INCIDENCE_CLASH) + word = None + if isinstance(polarization, str): + if polarization not in ('s', 'p'): + raise ValueError( + f"polarization must be 's', 'p' or a pair (a_s, a_p) of complex amplitudes, " + f'got {polarization!r}') + word = polarization + # The pair the word names, written as a caller writes it, so the conical and 2D + # paths are handed what `polarization=(1, 0)` / `(0, 1)` would hand them. + polarization = (1., 0.) if word == 's' else (0., 1.) pair = tuple(polarization) if len(pair) != 2: raise ValueError( @@ -1701,6 +1759,7 @@ def polarization(self, polarization): # inputs` while `(1.0+0j, -1j)` worked. A dtype is not something the caller should # have to know to get a gradient. self._polarization = tuple(self._as_complex(v) for v in pair) + self._polarization_word = word @property def lattice(self): @@ -2010,7 +2069,12 @@ def phi(self, phi): self._phi = None else: check_input_precision('phi', phi, self.complex_dtype) - self._phi = self._as_complex(phi) + # Under `_eager`, as the material setters are: whether the azimuth is a fixed zero + # decides a branch (`_polarization_word_holds`), and inside `jax.jit` a stated + # `phi=0.` was otherwise staged into a tracer that branch cannot tell from a swept + # one. A traced azimuth still stages. + with self._eager(): + self._phi = self._as_complex(phi) _kx = None _ky = None @@ -2026,11 +2090,14 @@ def _set_transverse_component(self, name, value): raise ValueError(problem) check_input_precision(name, value, self.complex_dtype) xp = self.xp - candidate = self._as_complex(value) - valid = xp.isfinite(candidate) & (xp.imag(candidate) == 0) - real = xp.real(candidate) - # Invalid traced values must not silently become a different real wave. - setattr(self, '_' + name, xp.where(valid, real, xp.full_like(real, float('nan')))) + # Under `_eager` for the reason `phi`'s setter gives: a stated `ky=0.` decides the + # polarization word's branch, and staged inside `jax.jit` it read as a traced one. + with self._eager(): + candidate = self._as_complex(value) + valid = xp.isfinite(candidate) & (xp.imag(candidate) == 0) + real = xp.real(candidate) + # Invalid traced values must not silently become a different real wave. + setattr(self, '_' + name, xp.where(valid, real, xp.full_like(real, float('nan')))) @property def kx(self): @@ -2971,12 +3038,17 @@ def _assign_grating_type(self): `3c7d166` fixed existed three times and had to be fixed three times. Same argument as the one that put the solve methods here (see this class's docstring). - * ``0`` -- pure 1D TE/TM: one grating vector, no azimuth, and a ``tetm`` declaration - that says the other polarization is not a degree of freedom of this problem. + * ``0`` -- pure 1D TE/TM: one grating vector, no azimuth, and a declaration that the + other polarization is not a degree of freedom of this problem -- ``tetm``, or a + polarization word ``'s'`` / ``'p'`` with no azimuth or ``ky`` left free + (:meth:`_polarization_word_holds`). * ``1`` -- 1D conical: the same one grating vector with the azimuth and the polarization both live. ``phi=None`` is what distinguishes it from ``phi=0``; a stated zero azimuth is still a stated one. * ``2`` -- 2D: everything else, and every scene whose k layout is not separable. + + The choice itself is :meth:`_routed_grating_type`, which refuses and warns nothing; this + puts the refusals and the warning a solve has to meet first in front of it, and assigns. """ self._refuse_a_structure_with_no_cell() self._refuse_a_raster_drawn_on_another_cell() @@ -2993,8 +3065,9 @@ def _assign_grating_type(self): one_dimensional = self.lattice_is_axis_aligned self._refuse_y_orders_on_a_y_invariant_cell(True) self._refuse_illegal_tetm(one_dimensional, one_row=True) - if self._tetm is None: - # Withdrawn since: the pure kernel cannot run undeclared (it raised + if not (one_dimensional and self._pure_path_declared()): + # Withdrawn since -- `tetm` cleared, or a polarization word whose phi, ky + # or cell has moved: the pure kernel cannot run undeclared (it raised # "tetm_index must be 0 (TE) or 1 (TM), got None"), so route as the raster # branch routes that raster undeclared -- conical on an axis-aligned cell, # the solve a never-declared solver gives bit for bit, and 2D on a skewed @@ -3002,13 +3075,30 @@ def _assign_grating_type(self): self._grating_type = 1 if one_dimensional else 2 return - if self._has_tensor_materials(): - self._grating_type = 2 - return - - if self.modeling_path == ModelingPath.RASTER: + if self.modeling_path == ModelingPath.RASTER and not self._has_tensor_materials(): one_row = self.structure.shape[1] == 1 self._refuse_y_orders_on_a_y_invariant_cell(one_row) + warn_if_xy_on_a_skewed_cell( + getattr(self, '_truncation_spelling', None), self.lattice_is_axis_aligned) + self._refuse_illegal_tetm(one_row and self.lattice_is_axis_aligned, one_row=one_row) + routed = self._routed_grating_type() + if routed is not None: + self._grating_type = routed + + def _routed_grating_type(self): + """The layout the structure and the declaration select -- ``0``, ``1`` or ``2``. + + :meth:`_assign_grating_type` assigns it once its refusals have passed, and a report + reads it before any solve (:meth:`_solve_path`). It refuses and warns nothing, so the + report cannot refuse a scene one call before the run would. ``None`` with no structure, + whose layout is kept from the last solve rather than chosen, and for a modeling path + neither branch knows. + """ + if self.structure is None: + return None + if self._has_tensor_materials(): + return 2 + if self.modeling_path == ModelingPath.RASTER: # One raster row means invariant along `a2`, which is a 1D grating -- but only on # an axis-aligned cell is it one the 1D operators can express. A skewed cell's # `b1` is perpendicular to `a2` rather than parallel to `a1`, so order `m` carries @@ -3017,21 +3107,65 @@ def _assign_grating_type(self): # problem: measured on a 50% binary grating, R out by 0.132 and T by 0.090 at a # 60-degree skew, energy conserving to 1e-15 throughout. See # `tests/rcwa/test_oblique_lattice.py::test_a_one_row_skewed_cell_is_the_grating_it_draws`. - one_dimensional = one_row and self.lattice_is_axis_aligned - warn_if_xy_on_a_skewed_cell( - getattr(self, '_truncation_spelling', None), self.lattice_is_axis_aligned) - self._refuse_illegal_tetm(one_dimensional, one_row=one_row) - if one_dimensional: - # The declaration selects the path, and nothing else does. A named incidence - # -- psi, or a pair, which mixes s and p the pure path cannot carry -- leaves - # the polarization a live parameter, and the conical formulation is the 1D - # one that carries it. - self._grating_type = 0 if self._tetm is not None else 1 - else: - self._grating_type = 2 + if not (self.structure.shape[1] == 1 and self.lattice_is_axis_aligned): + return 2 + # The declaration selects the path, and nothing else does. `tetm` declares the + # polarization and the azimuth away, and a polarization word declares the + # polarization fixed. A named incidence -- psi, chi, or a pair, which may mix s + # and p the pure path cannot carry, and whose value may be swept or + # differentiated -- leaves the polarization a live parameter, and the conical + # formulation is the 1D one that carries it. + return 0 if self._pure_path_declared() else 1 + if self.modeling_path == ModelingPath.PARAMETRIC: # staircase-free + return 2 + return None - elif self.modeling_path == ModelingPath.PARAMETRIC: # staircase-free - self._grating_type = 2 + def _pure_path_declared(self): + """Whether the incidence is declared to be one fixed polarization in the grating's plane. + + What the pure 1D TE/TM path needs beside a 1D structure: :attr:`tetm`, which declares + the polarization and the azimuth away and is refused where it cannot hold + (:meth:`_refuse_illegal_tetm`), or a polarization word that holds here. + """ + return self._tetm is not None or self._polarization_word_holds() + + def _polarization_word_holds(self): + """Whether ``polarization='s'`` / ``'p'`` takes the pure path on a 1D structure. + + The word fixes the polarization, and the pure path also needs the plane of incidence to + be the grating's own: ``phi`` and ``ky`` absent, or a zero that nothing traces or + differentiates (:func:`_absent_or_fixed_zero`). A nonzero one is off that plane, and a + swept or differentiated one is a degree of freedom whose derivative at zero only the + conical path has, so either keeps the conical path, silently -- the word states the + light, not the path. Every test is structural, so the answer is the same under + ``jax.jit`` and ``vmap`` as eagerly. The word against a pair: :attr:`polarization`. + """ + return (self._polarization_word is not None + and _absent_or_fixed_zero(self.phi, self._eager) + and _absent_or_fixed_zero(self.ky, self._eager)) + + #: :attr:`grating_type` by name, the pure layout split by its polarization (``tetm_index``): + #: what ``rcwa.CheckReport.resolved['solve_path']`` and ``rcwa.Result.diagnostics`` report. + _SOLVE_PATHS = {(0, 0): '1d_te', (0, 1): '1d_tm', (1, None): '1d_conical', (2, None): '2d'} + + def _solve_path(self, planned=False): + """The layout by name: the last solve's, or with ``planned`` the one a solve would take. + + ``'1d_te'`` / ``'1d_tm'`` (the pure path), ``'1d_conical'`` or ``'2d'``; ``None`` before + the first solve. ``planned`` reads the scene as it stands without solving or refusing: + :meth:`_routed_grating_type`, with ``formulation='first_order'`` on the 2D layout, as + :meth:`_solve_conv_signed` puts it (a z row is a tensor material's, already routed there). + """ + if not planned: + grating, index = self._grating_type, self.tetm_index + else: + grating = (2 if self._formulation == Formulation.FIRST_ORDER + else self._routed_grating_type()) + index = None + if grating == 0: + index = (int(self._tetm == 'tm') if self._tetm is not None + else int(self._polarization_word == 'p')) + return self._SOLVE_PATHS.get((grating, index if grating == 0 else None)) @property def grating_type(self): @@ -3088,10 +3222,11 @@ def _mu_1d(mu): def _solve_1d(self, epx_conv_all, epy_conv_all, epz_conv_i_all, mu=None, mu_top=None, mu_bot=None, eps_top=None, eps_bot=None, dng=None): - if self._polarization is not None: + if self._polarization is not None and self.tetm_index is None: # Unreachable through solve — `_assign_grating_type` sends a structure with a - # pair to the conical path, which carries both polarizations — so this only - # answers a direct call. + # pair, or with a word it did not route here, to the conical path, which carries + # both polarizations — so this only answers a direct call. A word routed here is + # read through `tetm_index`, and its pair is never used. raise ValueError( 'the pure 1D TE/TM path is single-polarization by construction and cannot ' 'take a polarization pair; use _solve_1d_conical, which solve() selects on ' @@ -3304,6 +3439,8 @@ def __init__(self, # end, for the reason above. polarization=None, # complex incident amplitude pair (a_s, a_p); # supersedes psi (naming both raises). See docs/reference/common/conventions.md. + # Or the word 's' / 'p': that pair, declared fixed, which takes the pure 1D + # path where it answers (see CoreSolveMixin.polarization). raster_normal_method=None, # raster NVF normal-field construction: # 'distance' (default) | 'contour' | 'spectral' | 'smooth'. Read only # where a raster normal field is built (factorization='nvf' on a raster); diff --git a/src/wavebit/rcwa/_engine/illumination.py b/src/wavebit/rcwa/_engine/illumination.py index c6accedb..cd79346d 100644 --- a/src/wavebit/rcwa/_engine/illumination.py +++ b/src/wavebit/rcwa/_engine/illumination.py @@ -465,6 +465,14 @@ class still land on the facade: the first is the facade's own slots, the second def formulation_assigned(self): return 'first_order' + def _solve_path(self, planned=False): + """``'2d'``, planned or taken: the side-resolved cascade is written on the 2D layout. + + The scene's own :meth:`~wavebit.rcwa._engine.engine.CoreSolveMixin._solve_path` would + answer ``None`` after a solve, since this cascade never routes it. + """ + return '2d' + @property def backend(self): return self._backend diff --git a/src/wavebit/rcwa/_lifecycle.py b/src/wavebit/rcwa/_lifecycle.py index 21681b10..11b1461f 100644 --- a/src/wavebit/rcwa/_lifecycle.py +++ b/src/wavebit/rcwa/_lifecycle.py @@ -145,6 +145,12 @@ def _report(problem, numerics, execution, solvers): 'layer_count': len(problem.stack.layers), 'fourier_coefficients': scene.fourier_coefficients.value, 'factorization': scene.factorization.value, 'algorithm': numerics.algorithm, 'formulation': numerics.formulation, + # The layout a solve takes, read from the solver before any has run: + # '1d_te' / '1d_tm' (the pure path a polarization word takes on a one-row + # raster), '1d_conical' or '2d' (also every coherent-source run, whose cascade + # is 2D). Result.diagnostics['solve_path'] has the layout each wavelength's + # solve took. + 'solve_path': solvers[0]._solve_path(planned=True), 'half_space_merge_tolerance': numerics.half_space_merge_tolerance, 'shape_nvf_grid_size': scene.shape_nvf_grid_size if scene.factorization.value == 'nvf' else None, 'backend': execution.backend, 'device': execution.device, 'precision': execution.precision, diff --git a/src/wavebit/rcwa/results.py b/src/wavebit/rcwa/results.py index 7c992c3a..789f302b 100644 --- a/src/wavebit/rcwa/results.py +++ b/src/wavebit/rcwa/results.py @@ -217,6 +217,7 @@ def pack(problem, numerics, execution, solves, checks, report): records, amplitudes, scattering = [], [], [] actual_wavelengths = [] actual_directions = [] + solve_paths = [] power_defined, efficiency_defined = [], [] merges = [] if single_plane(problem): @@ -228,6 +229,7 @@ def pack(problem, numerics, execution, solves, checks, report): for solver, result in solves: actual_wavelengths.append(_array(xp, solver.wavelength, execution)) actual_directions.append(solver.theta) + solve_paths.append(solver._solve_path()) # Which layers this wavelength's solve merged into a half-space, and which of them # a tolerance snapped: a snapped scene is a changed problem, recorded as detuning is. plan = solver._half_space_merged @@ -323,6 +325,7 @@ def pack(problem, numerics, execution, solves, checks, report): scene = solver._scene actual_wavelengths.append(_array(xp, scene.wavelength, execution)) actual_directions.append(None) + solve_paths.append(solver._solve_path()) region = result.power.region area = None if region == 'cell' if isinstance(region, str) else False: @@ -374,6 +377,9 @@ def pack(problem, numerics, execution, solves, checks, report): diagnostics={'report': report, 'effective_wavelengths': tuple(actual_wavelengths), 'source_group_counts': tuple(getattr(r, 'group_count', 1) for _, r in solves), 'effective_theta': tuple(actual_directions), + # The layout each wavelength's solve took: '1d_te' / '1d_tm' (the + # pure path), '1d_conical' or '2d'. report.resolved has the plan. + 'solve_path': tuple(solve_paths), 'port_power_defined': xp.stack(power_defined) if power_defined else None, # Says why power_entering and absorption are NaN when no # diffraction output carries the metadata flag. diff --git a/tests/README.md b/tests/README.md index 4240200e..0d28b0f3 100644 --- a/tests/README.md +++ b/tests/README.md @@ -123,6 +123,7 @@ compares the code to itself proves continuity, not correctness. | **Independent pair closure variation and native public geometry mapping** | `test_nvf_cartesian_rectangle_pair_api.py` | Frozen input indices preserve the two physical rectangles and material association; shared height/y/contact tangents, exact Fourier sets, uniform optical controls and package-only optional-backend behavior remain explicit. | | **An outside solver** | `test_reticolo_golden.py`, `test_fmmax_laurent_golden.py`, `test_oblique_lattice.py` (published Li 1997 table) | The strongest evidence here: a different group's code, in a different language. RETICOLO and fmmax's `FFT` formulation (against `factorization='laurent'`, which assembles the same operator) validate individual diffraction *orders*, not only totals. | | **Independent uniform interfaces** | `test_uniform_interfaces.py`, `fixtures/rcwa/uniform_interfaces/` | Do independent equal/neighboring material inputs preserve R/T, pointwise fields and derivatives against a 120-digit Fresnel reference, plus multilayer port and boundary continuity? | +| **A second formulation of the same problem** | `api/test_unified_rcwa_pure_1d.py` | Does the pure 1D TE/TM solve a polarization word takes reproduce every observable, field and derivative of the conical solve of the Jones pair it names, to round-off, under every backend and transform -- and does a word the pure path cannot answer solve bit for bit as that pair? The two kernels are written separately, so their agreement is evidence about both. They share the convolution matrices, so a factorization error is outside what it can see. | | **Physics** | `test_lossy_media.py`, `test_stability.py`, `test_fp_warnings.py`, `test_polarization.py` (circular states against the gyrotropic-slab Airy answer), `test_smm_field_state.py` (the plane wave inside a reflectionless uniform layer, every path's field, at complex angles and past the ETM's cap; SMM's field against the ETM's where its old states failed, at a layer light line, and against a label-free continuity-and-ports oracle), `test_slot_order.py` (every private builder splatted into `_solve_conv`), energy checks throughout | R+T = 1 on a lossless stack, Airy recursion on uniform layers. Independent of how we implement it. | | **Exact interface normals and material continuity** | `test_nvf_seam_materials.py` | Distinct equal-modulus materials retain the exact seam normal and a Hermitian rank-one projector; coefficient and optical changes shrink with the material perturbation, listings agree, and JAX/Torch gradients match secants. The study's opening creates a positive geometric gap. This does not certify optical convergence. | | **Painted material limits and geometric invariants** | `test_painted_material_continuity.py` | Real, absorbing and metallic equality ramps have no fixed step. Exact straight-interface normals, rank-one projectors, physical rotation/reflection, passivity, chunking, and JAX/Torch material derivatives constrain the contrast envelope independently of its construction. | diff --git a/tests/api/test_unified_rcwa_pure_1d.py b/tests/api/test_unified_rcwa_pure_1d.py new file mode 100644 index 00000000..2bbd350c --- /dev/null +++ b/tests/api/test_unified_rcwa_pure_1d.py @@ -0,0 +1,416 @@ +"""A polarization word on a one-row raster takes the pure 1D TE/TM path, and nothing else moves. + +``PlaneWave(polarization='s')`` / ``'p'`` is a declaration that the polarization is fixed, and +the engine's routing (``CoreSolveMixin._polarization_word_holds``) takes the pure TE/TM path +for it wherever that path answers the same light. A Jones pair is a value and keeps the conical +path. Since 2026-11-03; before, every public 1D run was conical. + +Two families of checks. The **regression matrix** runs each feature twice, once with the word +(the pure path) and once with the pair that names the same light (the conical path, which +every public 1D run took until this change), and holds every observable to round-off between +the two: two formulations of one problem, written independently, agree. The **fallback pins** +hold the cases the word does not take pure -- a skewed cell, a traced or differentiated azimuth, +a nonzero ``phi`` or ``ky``, ``formulation='first_order'``, a requested scattering matrix -- to +the pair *bit for bit*, which is what the word was lowered to before. Both also check the path +the run reports it took (``Result.diagnostics['solve_path']``) against what ``check`` planned +(``CheckReport.resolved['solve_path']``) and, where a transform could hide it, against the +kernel that was actually called. +""" +import collections + +import numpy as np +import pytest + +from wavebit import Material, Units, rcwa +from wavebit.rcwa._engine import engine as _engine + +NX = 40 +PERIOD, THICKNESS, WAVELENGTH = 1000.0, 200.0, 1363.0 +#: The word and the pair that names the same light, and the pure path the word takes. +WORDS = {'s': ((1.0, 0.0), '1d_te'), 'p': ((0.0, 1.0), '1d_tm')} +ORDERS = rcwa.Numerics(orders=(8, 0)) +X = np.linspace(-480.0, 480.0, 9) +Y = np.array([0.0, 130.0]) + + +def raster(epsilon=12.25): + """One row, 40 px: a bar of 18 px, off the cell centre so no pair of orders is mirror-equal.""" + grid = np.ones((1, NX), dtype=complex) + grid[0, 9:27] = epsilon + return grid + + +def grating(polarization='s', *, theta=0.2, phi=None, kx=None, ky=None, side='top', top=1.0, + bottom=2.25, wavelength=WAVELENGTH, layers=None, lattice=None, mu=None, + epsilon=None, thickness=THICKNESS, fields=True, scattering=False): + layer = rcwa.Layer(thickness=thickness, epsilon=raster() if epsilon is None else epsilon, mu=mu) + stack = dict(top=top if isinstance(top, Material) else Material(top), + bottom=bottom if isinstance(bottom, Material) else Material(bottom), + layers=layers or (layer,)) + stack.update({'period': (PERIOD, PERIOD)} if lattice is None else {'lattice': lattice}) + source = dict(polarization=polarization, side=side) + source.update(dict(kx=kx, ky=ky) if kx is not None or ky is not None else dict(theta=theta, phi=phi)) + return rcwa.Problem(units=Units(length='nm'), stack=rcwa.Stack(**stack), + excitation=rcwa.Excitation(wavelength=wavelength, + sources=(rcwa.PlaneWave(**source),)), + outputs=rcwa.Outputs(fields=fields, scattering=scattering)) + + +def host(value): + if hasattr(value, 'detach'): + return value.detach().cpu().numpy() + return np.asarray(value) + + +def observables(result, problem): + """Every number a single-plane-wave result carries, by name.""" + out = {} + for port in ('reflected', 'transmitted'): + o = result.outputs[port] + for name in ('amplitude_s', 'amplitude_p', 'amplitude_xyz', 'efficiency', 'power', + 'k_parallel', 'k_normal'): + out[f'{port}.{name}'] = host(getattr(o, name)) + out['absorption'] = host(result.absorption) + out['power_entering'] = host(result.outputs['power_entering']) + if problem.outputs.fields: + total = sum(float(host(layer.thickness)) for layer in problem.stack.layers) + z = np.array([0.0, 0.31*total, 0.77*total, total]) + out['field'] = host(result.field(x=X, y=Y, z=z).values) + for side in ('top', 'bottom'): + out[f'boundary.{side}'] = host(result.boundary_fields(side=side, x=X, y=Y).values) + return out + + +def _family(name): + """The observables one scale is read from: a quantity that is round-off on one path has no + scale of its own. A cross-polarized amplitude is an exact zero on the pure path and noise on + the conical one, measured against the co-polarized amplitude of the same port; a lossless + stack's absorption is measured against the power entering it.""" + port, _, quantity = name.rpartition('.') + if quantity in ('amplitude_s', 'amplitude_p'): + return f'{port}.amplitude_sp' + return 'budget' if name in ('absorption', 'power_entering') else name + + +def assert_round_off(pure, conical, rtol): + """Each observable to ``rtol`` of its family's largest entry, NaN where the other is NaN.""" + assert pure.keys() == conical.keys() + scales = collections.defaultdict(float) + for name in pure: + for values in (pure[name], conical[name]): + scales[_family(name)] = max(scales[_family(name)], + float(np.nanmax(np.abs(values), initial=0.0))) + for name in pure: + a, b = pure[name], conical[name] + assert a.shape == b.shape, name + assert np.array_equal(np.isnan(a), np.isnan(b)), name + gap = float(np.nanmax(np.abs(a - b), initial=0.0)) + scale = scales[_family(name)] + assert gap <= rtol*scale, f'{name}: {gap:.2e} of {scale:.2e}' + + +def assert_bits(word, pair): + assert word.keys() == pair.keys() + for name in word: + assert word[name].dtype == pair[name].dtype and np.array_equal(word[name], pair[name], + equal_nan=True), name + + +def paths(result): + return result.diagnostics['solve_path'], result.diagnostics['report'].resolved['solve_path'] + + +def run_both(word, make, numerics=ORDERS, backend='numpy', precision='double'): + """The word's run and the pair's, with the paths each reported.""" + pair, _ = WORDS[word] + out = [] + for polarization in (word, pair): + problem = make(polarization) + result = rcwa.run(problem, numerics=numerics, backend=backend, precision=precision) + out.append((observables(result, problem), result)) + return out + + +# ---------------------------------------------------------------- the regression matrix + +k0 = 2*np.pi/WAVELENGTH +_mu = np.ones((1, NX), dtype=complex) +_mu[0, 9:27] = 1.44 + +#: Each feature the pure path has to answer as the conical one does. +CASES = { + 'normal incidence': (lambda pol: grating(pol, theta=0.0), ORDERS), + 'theta 0.3': (lambda pol: grating(pol, theta=0.3), ORDERS), + 'theta -0.3': (lambda pol: grating(pol, theta=-0.3), ORDERS), + 'lit from below, glass below': (lambda pol: grating(pol, theta=0.3, side='bottom'), ORDERS), + 'a stated zero azimuth': (lambda pol: grating(pol, theta=0.3, phi=0.0), ORDERS), + 'kx spelling, ky = 0': (lambda pol: grating(pol, kx=0.3*k0, ky=0.0), ORDERS), + 'kx spelling alone': (lambda pol: grating(pol, kx=-0.2*k0), ORDERS), + 'raster mu': (lambda pol: grating(pol, mu=_mu), ORDERS), + 'lossy bar': (lambda pol: grating(pol, epsilon=raster(12.25+1j)), ORDERS), + 'magnetic half-spaces': (lambda pol: grating(pol, top=Material(2.25, mu=1.4), + bottom=Material(3.24, mu=1.1)), ORDERS), + 'a negative-index substrate': (lambda pol: grating(pol, bottom=Material(-2.0+0.05j, mu=-1.2+0.05j)), + ORDERS), + 'layers merged into the half-spaces': (lambda pol: grating(pol, layers=( + rcwa.Layer(thickness=30.0, epsilon=np.full((1, NX), 2.25+0j)), + rcwa.Layer(thickness=THICKNESS, epsilon=raster()), + rcwa.Layer(thickness=50.0, epsilon=np.full((1, NX), 1.0+0j)))), ORDERS), + 'a spectrum': (lambda pol: grating(pol, wavelength=(1300.0, 1363.0, 1420.0)), ORDERS), + 'ETM': (lambda pol: grating(pol), rcwa.Numerics(orders=(8, 0), algorithm='etm')), + 'detuning': (lambda pol: grating(pol), rcwa.Numerics(orders=(8, 0), detuning=1e-12)), + 'nvf': (lambda pol: grating(pol), rcwa.Numerics(orders=(8, 0), factorization='nvf')), + 'laurent, integrated': (lambda pol: grating(pol), rcwa.Numerics( + orders=(8, 0), factorization='laurent', fourier_coefficients='integrated')), + 'second_order stated': (lambda pol: grating(pol), rcwa.Numerics(orders=(8, 0), + formulation='second_order')), + # The -1 order exactly grazing in air: a Rayleigh point. + 'a Rayleigh point': (lambda pol: grating(pol, theta=float(np.arcsin(WAVELENGTH/PERIOD - 1.0)), + bottom=1.0), ORDERS), + # theta = pi/2 + 0.2j from glass: an evanescent incident wave, which has no power to + # normalize by, so efficiencies are NaN on both paths and the amplitudes carry the answer. + 'evanescent incidence': (lambda pol: grating(pol, theta=np.pi/2 + 0.2j, top=2.25, fields=False), + ORDERS), +} + + +@pytest.mark.parametrize('word', sorted(WORDS)) +@pytest.mark.parametrize('case', sorted(CASES)) +def test_the_pure_path_answers_as_the_conical_one(case, word): + make, numerics = CASES[case] + (pure, word_run), (conical, pair_run) = run_both(word, make, numerics) + expected = WORDS[word][1] + count = len(np.atleast_1d(make(word).excitation.wavelength)) + assert paths(word_run) == ((expected,)*count, expected) + assert paths(pair_run) == (('1d_conical',)*count, '1d_conical') + assert_round_off(pure, conical, rtol=1e-11) + + +@pytest.mark.parametrize('backend', ['jax', 'torch']) +@pytest.mark.parametrize('word', sorted(WORDS)) +def test_the_pure_path_answers_on_every_backend(backend, word): + pytest.importorskip(backend) + if backend == 'jax': + import jax + jax.config.update('jax_enable_x64', True) + (pure, word_run), (conical, _) = run_both(word, lambda pol: grating(pol, theta=0.25), + backend=backend) + assert paths(word_run) == ((WORDS[word][1],), WORDS[word][1]) + assert_round_off(pure, conical, rtol=1e-11) + + +@pytest.mark.parametrize('word', sorted(WORDS)) +def test_single_precision_agrees_to_its_own_round_off(word): + (pure, word_run), (conical, _) = run_both(word, lambda pol: grating(pol, theta=0.25), + precision='single') + assert paths(word_run) == ((WORDS[word][1],), WORDS[word][1]) + assert_round_off(pure, conical, rtol=1e-4) + + +def test_check_and_prepare_report_the_path_before_any_solve(): + for word, (pair, expected) in WORDS.items(): + assert rcwa.check(grating(word), numerics=ORDERS).resolved['solve_path'] == expected + prepared = rcwa.prepare(grating(word), numerics=ORDERS) + assert prepared.report.resolved['solve_path'] == expected + assert prepared.validate(grating(word)).resolved['solve_path'] == expected + assert rcwa.check(grating(pair), numerics=ORDERS).resolved['solve_path'] == '1d_conical' + skew = np.array([[PERIOD, 0.3*PERIOD], [0.0, PERIOD]]) + nvf = rcwa.Numerics(orders=(8, 0), factorization='nvf') + assert rcwa.check(grating('s', lattice=skew), numerics=nvf).resolved['solve_path'] == '2d' + assert rcwa.check(grating('s', phi=0.3), numerics=ORDERS).resolved['solve_path'] == '1d_conical' + assert rcwa.check(grating('s'), numerics=rcwa.Numerics(orders=(8, 0), formulation='first_order') + ).resolved['solve_path'] == '2d' + assert rcwa.check(grating('s', scattering=True), numerics=ORDERS + ).resolved['solve_path'] == '1d_conical' + + +# ---------------------------------------------------------------- derivatives and transforms + +def _jax(): + jax = pytest.importorskip('jax') + jax.config.update('jax_enable_x64', True) + return jax + + +@pytest.mark.parametrize('word', sorted(WORDS)) +def test_jax_derivatives_through_the_pure_path_are_the_conical_ones(word): + jax = _jax() + import jax.numpy as jnp + + def total_R(polarization, variable): + def R(v): + kw = {'thickness': dict(thickness=v), 'theta': dict(theta=v), + 'epsilon': dict(epsilon=jnp.where(jnp.asarray(raster()) == 1.0, 1.0+0j, v+0j))}[variable] + return jnp.sum(rcwa.run(grating(polarization, fields=False, **kw), numerics=ORDERS, + backend='jax').R) + return R + + pair = WORDS[word][0] + for variable, at in (('thickness', THICKNESS), ('theta', 0.2), ('epsilon', 12.25)): + pure = jax.grad(total_R(word, variable))(jnp.float64(at)) + conical = jax.grad(total_R(pair, variable))(jnp.float64(at)) + assert abs(float(pure) - float(conical)) <= 1e-11*max(abs(float(conical)), 1e-3), variable + + +@pytest.mark.parametrize('word', sorted(WORDS)) +def test_torch_autograd_through_the_pure_path_is_the_conical_one(word): + torch = pytest.importorskip('torch') + + def gradients(polarization): + thickness = torch.tensor(THICKNESS, dtype=torch.float64, requires_grad=True) + bar = torch.tensor(12.25, dtype=torch.float64, requires_grad=True) + epsilon = torch.where(torch.as_tensor(raster().real) == 1.0, torch.ones(1, NX, dtype=torch.float64), + bar).to(torch.complex128) + result = rcwa.run(grating(polarization, thickness=thickness, epsilon=epsilon, fields=False), + numerics=ORDERS, backend='torch') + torch.sum(result.R).backward() + return result.diagnostics['solve_path'], float(thickness.grad), float(bar.grad) + + path, *pure = gradients(word) + assert path == (WORDS[word][1],) + _, *conical = gradients(WORDS[word][0]) + np.testing.assert_allclose(pure, conical, rtol=1e-11) + + +@pytest.fixture +def kernels(monkeypatch): + """Which kernel each solve called: a transform traces the call even when it hides a result.""" + calls = collections.Counter() + + def counted(kind, kernel): + def call(*args, **kwargs): + calls[kind] += 1 + return kernel(*args, **kwargs) + return call + monkeypatch.setattr(_engine, '_SOLVES', { + algorithm: {kind: counted(kind, kernel) for kind, kernel in table.items()} + for algorithm, table in _engine._SOLVES.items()}) + return calls + + +@pytest.mark.parametrize('word', sorted(WORDS)) +def test_a_prepared_run_takes_the_pure_path_under_jit_and_vmap(word, kernels): + """A stated `phi=0.` is a fixed zero inside the trace too, so the word stays pure there.""" + jax = _jax() + import jax.numpy as jnp + make = lambda t: grating(word, phi=0.0, thickness=t, fields=False) # noqa: E731 + prepared = rcwa.prepare(make(THICKNESS), numerics=ORDERS, backend='jax') + R = lambda t: jnp.sum(prepared.run(make(t)).R) # noqa: E731 + thicknesses = jnp.asarray([190.0, 200.0, 210.0]) + eager = [float(R(t)) for t in thicknesses] + eager_grad = [float(jax.grad(R)(t)) for t in thicknesses] + kernels.clear() + jitted = jax.jit(lambda t: prepared.run(make(t)))(jnp.float64(THICKNESS)) + assert jitted.diagnostics['solve_path'] == (WORDS[word][1],) + np.testing.assert_allclose(jax.vmap(R)(thicknesses), eager, rtol=1e-12) + np.testing.assert_allclose(jax.jit(jax.vmap(jax.grad(R)))(thicknesses), eager_grad, rtol=1e-10) + assert kernels['1d'] > 0 and kernels['conical'] == 0, kernels + + +# ---------------------------------------------------------------- the fallback pins + +def _skew(): + return np.array([[PERIOD, 0.3*PERIOD], [0.0, PERIOD]]) + + +#: Where the word does not take the pure path, and the layout it takes instead. +FALLBACKS = { + # A one-row raster on a skewed cell is a 1D grating whose grooves do not run along y; the + # 2D layout is the one that expresses it (the scratch prototype that passed `tetm` from the + # adapter turned this solve into a ValidationError). + 'a skewed cell': (lambda pol: grating(pol, lattice=_skew()), + rcwa.Numerics(orders=(8, 0), factorization='nvf'), '2d'), + 'a nonzero azimuth': (lambda pol: grating(pol, phi=0.3), ORDERS, '1d_conical'), + 'a nonzero ky': (lambda pol: grating(pol, kx=0.2*k0, ky=0.1*k0), ORDERS, '1d_conical'), + 'first_order': (lambda pol: grating(pol), rcwa.Numerics(orders=(8, 0), formulation='first_order'), + '2d'), + # The public scattering matrix keeps both polarizations' ports, 2N x 2N with channel 'sp'; + # the pure path's is N x N, one polarization's, so a run that asks for one stays conical. + 'a scattering request': (lambda pol: grating(pol, scattering=True, fields=False), ORDERS, + '1d_conical'), +} + + +def _scattering_blocks(result): + if result._scattering is None: + return {} + return {f'S{i}.{block}': host(getattr(matrix, block)) + for i, matrix in enumerate(result.scattering) for block in ('S11', 'S12', 'S21', 'S22')} + + +@pytest.mark.parametrize('word', sorted(WORDS)) +@pytest.mark.parametrize('case', sorted(FALLBACKS)) +def test_a_word_the_pure_path_cannot_answer_solves_bit_for_bit_as_its_pair(case, word): + make, numerics, expected = FALLBACKS[case] + (a, word_run), (b, pair_run) = run_both(word, make, numerics) + assert paths(word_run) == paths(pair_run) == ((expected,), expected) + a.update(_scattering_blocks(word_run)) + b.update(_scattering_blocks(pair_run)) + assert_bits(a, b) + if case == 'a scattering request': + modes = 2*8 + 1 + assert [np.shape(m.S11) for m in word_run.scattering] == [(2*modes, 2*modes)] + assert word_run.scattering[0].entry_port.channel == 'sp' + + +def test_a_coherent_source_run_reports_the_2d_cascade_it_takes(): + """Several sources solve on the side-resolved first-order cascade, a word among them.""" + problem = grating('s', wavelength=700.0, fields=False) + problem = problem.replace( + outputs=rcwa.Outputs(diffraction=False), + excitation=problem.excitation.replace(sources=( + rcwa.PlaneWave(polarization='s', theta=0.1), rcwa.PlaneWave(polarization='p', theta=0.1)))) + assert rcwa.check(problem, numerics=ORDERS).resolved['solve_path'] == '2d' + assert paths(rcwa.run(problem, numerics=ORDERS)) == (('2d',), '2d') + + +def test_a_jones_pair_is_a_value_and_stays_conical(): + for pair in ((1.0, 0.0), (0.0, 1.0), (0.6, 0.8j)): + result = rcwa.run(grating(pair, fields=False), numerics=ORDERS) + assert paths(result) == (('1d_conical',), '1d_conical') + + +def test_a_traced_or_differentiated_azimuth_stays_conical(kernels): + """`phi` swept or differentiated is a degree of freedom only the conical path carries. + + The pair for ``'p'`` is written ``(0, 2)``, which the public door normalizes to exactly + ``(0, 1)``. Written ``(0, 1)``, ``jax.jit`` of the run aborts the process -- before this + change as after it, and on a pair alone: XLA's algebraic simplifier (jaxlib 0.11.0) reads the + runtime guard's constant ``(0, 1)`` as an iota and rewrites its complex ``abs`` to a value of + the wrong type. + """ + jax = _jax() + import jax.numpy as jnp + torch = pytest.importorskip('torch') + transforms = {'jit': jax.jit, 'grad': jax.grad, 'jit(grad)': lambda f: jax.jit(jax.grad(f)), + 'vmap': jax.vmap} + + def total_R(polarization): + prepared = rcwa.prepare(grating(polarization, phi=0.0, fields=False), numerics=ORDERS, + backend='jax') + return lambda ph: jnp.sum(prepared.run(grating(polarization, phi=ph, fields=False)).R) + + for word, (pair, _) in WORDS.items(): + pair = (0.0, 2.0) if word == 'p' else pair + R = {polarization: total_R(polarization) for polarization in (word, pair)} + kernels.clear() + for name, transform in transforms.items(): + at = jnp.asarray([0.0, 0.1]) if name == 'vmap' else jnp.float64(0.0) + np.testing.assert_array_equal(np.asarray(transform(R[word])(at)), + np.asarray(transform(R[pair])(at)), err_msg=name) + assert kernels['1d'] == 0 and kernels['conical'] > 0, kernels + # The plan is read from the concrete problem the preparation was given; the record is + # the path the traced run took. + prepared = rcwa.prepare(grating(word, phi=0.0, fields=False), numerics=ORDERS, backend='jax') + traced = jax.jit(lambda ph: prepared.run(grating(word, phi=ph, fields=False)))(jnp.float64(0.0)) + assert paths(traced) == (('1d_conical',), WORDS[word][1]) + + phi = {pol: torch.tensor(0.0, dtype=torch.float64, requires_grad=True) for pol in (word, pair)} + values = {} + for pol in (word, pair): + result = rcwa.run(grating(pol, phi=phi[pol], fields=False), numerics=ORDERS, backend='torch') + torch.sum(result.R).backward() + values[pol] = (result.diagnostics['solve_path'], host(result.R), float(phi[pol].grad)) + assert values[word][0] == ('1d_conical',) + np.testing.assert_array_equal(values[word][1], values[pair][1]) + assert values[word][2] == values[pair][2] diff --git a/tests/rcwa/test_pytree_round_trip.py b/tests/rcwa/test_pytree_round_trip.py index fda87d43..ac0c1fee 100644 --- a/tests/rcwa/test_pytree_round_trip.py +++ b/tests/rcwa/test_pytree_round_trip.py @@ -55,7 +55,9 @@ #: fact about the *flatten* (it compared the `Lattice.drawn` declaration with a concrete cell), #: not about the solver; `test_the_declaration_survives_every_backend_and_the_jit_boundary` and #: `test_the_flatten_check_is_carried_and_keyed_on_the_objects_it_compared` pin it. -AUX_NOT_IN_CARRIED = {'psi_given', 'device', 'drawn_checked'} +#: `polarization_word` is the word `polarization` was written as, beside the pair that +#: attribute reads back; `test_a_polarization_word_survives_beside_its_pair` checks it. +AUX_NOT_IN_CARRIED = {'psi_given', 'device', 'drawn_checked', 'polarization_word'} def _scene(**over): @@ -184,6 +186,23 @@ def test_psi_given_survives_because_the_psi_child_would_overwrite_it(): assert _round_trip(solver)._psi_given is given, kw +def test_a_polarization_word_survives_beside_its_pair(): + """`polarization='p'` is a declaration beside the pair `(0, 1)` it reads back as. + + The pair is a child and comes back as a value; without the word in aux the rebuilt solver + would hold a pair, and its solve would take the conical path where the original took the + pure one. A pair written as a pair has no word, and must not gain one. + """ + for given, word in (('p', 'p'), ((0.0, 1.0), None)): + solver = _one_row(polarization=given) + assert solver._polarization_word == word + rebuilt = _round_trip(solver) + assert rebuilt._polarization_word == word + assert _same(rebuilt.polarization, solver.polarization) + solver.solve(), rebuilt.solve() + assert rebuilt.grating_type == solver.grating_type == (0 if word else 1) + + def test_a_declaration_assigned_after_construction_comes_back(): """`tetm` set on a live 1D solver round-trips, as it does when it was given to `rcwa`. From 5eb84b9f5835789fa90c293f991676ddec8949a0 Mon Sep 17 00:00:00 2001 From: Yongha Kim Date: Wed, 4 Nov 2026 00:09:02 +0900 Subject: [PATCH 2/2] docs(rcwa): the words 's' and 'p' are a declaration, and a run names its path Public 1D runs are no longer always conical, so the statements that said they were are rewritten, in English and in Korean where a translation exists. - reference/keywords.md: the tetm row now points to the words. The closing paragraph no longer sends readers to the earlier interface, which was retired on 2026-10-15. - guide/rcwa/light.md: new section "The words are a declaration, the pair is a value", with the conditions, the fallbacks, the round-off agreement, the cost, and resolved/diagnostics['solve_path'], including what a traced run reports. - guide/rcwa/scene.md, fields.md (the pure path's exact zeros) and numerics.md (the path is not a setting). - guide/rcwa/anisotropic-materials.md section 8: the routing account and the "which path ran" requirement, which the result now meets. - learn/fundamentals light, overview and index, with their Korean pages. - development/source-map (and its Korean page) and tests/README name the new test. Co-Authored-By: Claude Opus 5.5 --- docs/development/source-map.md | 1 + docs/guide/rcwa/anisotropic-materials.md | 16 +++-- docs/guide/rcwa/fields.md | 11 +++- docs/guide/rcwa/light.md | 65 ++++++++++++++++--- docs/guide/rcwa/numerics.md | 2 + docs/guide/rcwa/scene.md | 4 +- docs/learn/fundamentals/index.md | 3 +- docs/learn/fundamentals/light.md | 4 +- docs/learn/fundamentals/overview.md | 3 +- .../ko/pages/development/source-map.md | 1 + .../ko/pages/learn/fundamentals/index.md | 2 +- .../ko/pages/learn/fundamentals/light.md | 4 +- .../ko/pages/learn/fundamentals/overview.md | 2 +- docs/reference/common/conventions.md | 2 +- docs/reference/keywords.md | 12 ++-- 15 files changed, 101 insertions(+), 31 deletions(-) diff --git a/docs/development/source-map.md b/docs/development/source-map.md index 2ea32094..304bb669 100644 --- a/docs/development/source-map.md +++ b/docs/development/source-map.md @@ -83,6 +83,7 @@ research cases; inspect its scope before launching a long calculation. | --- | --- | --- | | Material response or SI conversion | [test_common_material_models.py](https://github.com/fieldalgo/wavebit/blob/main/tests/materials/test_common_material_models.py), [test_lifecycle.py](https://github.com/fieldalgo/wavebit/blob/main/tests/fdtd/test_lifecycle.py) | Equivalent geometry units, SI conductivity/poles, supported causal material admission and intact derivative graphs. | | RCWA preparation, live values or result ownership | [test_unified_rcwa.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa.py) | Reuse versus fresh preparation, independent Fresnel fields/powers, invalid transformed inputs and earlier-result independence. | +| RCWA solve path (pure 1D TE/TM, conical or 2D) | [test_unified_rcwa_pure_1d.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa_pure_1d.py), [test_tetm_off_the_pure_path.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_tetm_off_the_pure_path.py) | A polarization word against the pair it names: every observable and derivative to round-off where the word takes the pure path, bit for bit where it falls back, and the reported path against the kernel actually called. | | Coherent sources and absolute power | [test_unified_rcwa_sources.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa_sources.py) | SI source strength, finite observation regions and independently evaluated fields. | | FDTD domain, source/monitor or reuse | [test_lifecycle.py](https://github.com/fieldalgo/wavebit/blob/main/tests/fdtd/test_lifecycle.py) | PML outside the physical domain, unchanged laboratory coordinates, live source/material reuse and resolved duration guards. | | RCWA modal or reconstruction algebra | [test_smm_field_state.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_smm_field_state.py), [test_lossy_media.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_lossy_media.py) | Independent interior fields and interface responses rather than a regenerated solver golden. | diff --git a/docs/guide/rcwa/anisotropic-materials.md b/docs/guide/rcwa/anisotropic-materials.md index 1f1eebda..29fd6307 100644 --- a/docs/guide/rcwa/anisotropic-materials.md +++ b/docs/guide/rcwa/anisotropic-materials.md @@ -492,14 +492,15 @@ on 2026-08-13 noted against each part, because the reasoning is what a future ax need and the code alone does not carry it. **The principle is already in the codebase.** `engine.py::_assign_grating_type` reads the `structure` -shape, `tetm` and `phi` and picks the cheapest formulation that can express the problem — its +shape, the incidence's declaration (`tetm`, or a polarization word `'s'` / `'p'`) and `phi` and `ky`, +and picks the cheapest formulation that can express the problem — its docstring opens *"Pick the formulation this scene needs: pure 1D TE/TM, 1D conical, or 2D."* 1D TE/TM splits into scalars, 1D conical does not, 2D couples everything. Anisotropy adds axes to that same decision rather than a new kind of decision. | axis | read from | branches | |---|---|---| -| grating dimensionality *(exists)* | `structure` shape, `tetm`, `phi` | 1D TE/TM, 1D conical, 2D | +| grating dimensionality *(exists)* | `structure` shape, `tetm` or a polarization word, `phi`, `ky` | 1D TE/TM, 1D conical, 2D | | **z row of the tensor** *(exists)* | whether the four z-row slots were given | second-order, first-order | | **permeability** | `mu != 1` | existing operator, magnetic operator | | **patterned or uniform** | layer has in-plane structure | factorization, or scaled identities | @@ -521,8 +522,11 @@ answered for s light whichever of the two it named. The message names the pair t the same light, `polarization=(1, 0)` / `(0, 1)` at `phi = 0`, and `tests/rcwa/test_tetm_off_the_pure_path.py` checks that it does: beside an all-zero z row, and under a forced `'first_order'`, the pair reproduces the pure `tetm` solve. `tetm` belongs to the -earlier interface only: `rcwa.PlaneWave` always names the polarization (`'s'`, `'p'` or a Jones -pair), so the current lifecycle never meets this refusal. +earlier interface only, so the current lifecycle never meets this refusal. It reaches the pure +path through a polarization word instead (`PlaneWave(polarization='s')` / `'p'`, since +2026-11-03), which declares the polarization fixed rather than absent and so is never refused: +beside a z row, under `'first_order'`, on a skewed cell or with a free azimuth it solves as the +pair it names, on the path that carries it (`engine.py::_polarization_word_holds`). ### Why it has to be automatic @@ -576,7 +580,9 @@ exactly this reason. A silent 8x is indistinguishable from a slow machine, and a cannot say which formulation it measured is not a benchmark. The attribute is `solver.formulation_assigned`, `'second_order'` or `'first_order'`, set by every 2D solve and `None` before the first one. Both are attributes of the earlier interface's solver. The current -`rcwa.Result` does not report the formulation or the path that ran yet, and the check report's +`rcwa.Result` reports the path each wavelength's solve took, `diagnostics['solve_path']` +(`'1d_te'`, `'1d_tm'`, `'1d_conical'` or `'2d'`), and the check report names the one a run will +take, `resolved['solve_path']`. It does not report which formulation a 2D solve assigned yet: `resolved['formulation']` is the requested value (`'auto'` by default), not the one assigned. ## References diff --git a/docs/guide/rcwa/fields.md b/docs/guide/rcwa/fields.md index 228ff2ac..a5a32f3b 100644 --- a/docs/guide/rcwa/fields.md +++ b/docs/guide/rcwa/fields.md @@ -125,12 +125,19 @@ The time convention is `e^{-iωt}`: the instantaneous real field is `real(E * exp(-1j*omega*t))`. `abs(E)` and `angle(E)` give each component's magnitude and phase. See [conventions](../../reference/common/conventions.md) for the source phase reference. -Every input is an electric amplitude: coupled s/p, `PlaneWave`, and the earlier interface's -pure 1D TM solve included. (Until 2026-09-27 the pure 1D TM solve was lit by a unit magnetic +Every input is an electric amplitude: coupled s/p, `PlaneWave`, and the pure 1D TM solve +included. (Until 2026-09-27 the pure 1D TM solve was lit by a unit magnetic input, so its field was `mu_inc/n_inc` of the conical p field for the same light.) The result API preserves the solved input's scale; it does not rescale the field by incident power. +A `PlaneWave` written `'s'` or `'p'` on a one-row raster solves on that pure TE/TM path +([the light guide](light.md#the-words-are-a-declaration-the-pair-is-a-value) says when), and its +`field` and `boundary_fields` carry the same six components as the conical path's, to round-off. +The three the polarization does not excite -- `E_x`, `E_z` and `H_y` for `'s'`, `E_y`, `H_x` and +`H_z` for `'p'` -- are exact zeros there, where the conical path returns round-off (1e-16 of the +excited components). + ## Boundary waves and polarization variants `side='top'` or `'bottom'` selects the outside half-space limit at that interface. diff --git a/docs/guide/rcwa/light.md b/docs/guide/rcwa/light.md index ed5ccee5..9c9e42e1 100644 --- a/docs/guide/rcwa/light.md +++ b/docs/guide/rcwa/light.md @@ -15,9 +15,10 @@ Plane-wave amplitude is in V/m. Polarization is `"s"`, `"p"`, or a complex `(s,p The [reference](../../reference/rcwa/index.md) lists their constructors. The sections below give the details of direction, polarization and incidence side. :::{note} -Three features below are not in the current API yet, so their passages use the earlier -`wavebit.compat.legacy` interface and say so: the pure 1D `tetm` declaration, the `psi`/`chi` -polarization coordinates, and the one-solve s/p responses `s_inc`/`p_inc`. +Two features below are not in the current API yet, so their passages use the earlier +`wavebit.compat.legacy` interface and say so: the `psi`/`chi` polarization coordinates, and the +one-solve s/p responses `s_inc`/`p_inc`. A third, the pure 1D `tetm` declaration, is now the +words `'s'` and `'p'` (below). ::: Specify the incident direction of each source with either angles or transverse wavevector @@ -178,8 +179,10 @@ that reaches the far side. ## The polarization comes in three kinds -The current `PlaneWave` takes the third kind, `polarization='s'`, `'p'` or `(a_s, a_p)`. The -first two are earlier-interface declarations that the current API does not have yet. +The current `PlaneWave` takes `polarization='s'`, `'p'` or `(a_s, a_p)`. The pair is the third +kind below, and the words now do the first kind's work: they declare the polarization fixed (see +"The words are a declaration" below). The second kind, the coordinates, is an earlier-interface +form the current API does not have yet. :::{admonition} Earlier interface :class: earlier-interface @@ -203,7 +206,8 @@ rcwa.PlaneWave(polarization=(1, -1j)) # the pair, straight from anothe Retired with `wavebit.compat` on 2026-10-15; this runs only on a historical checkout. ```python -# Earlier interface: the current API has no tetm or psi/chi declaration yet. +# Earlier interface: the current API has no psi/chi declaration yet, and tetm's pure path +# is taken by polarization='s' / 'p'. import numpy as np from wavebit.compat import legacy as wavebit @@ -213,9 +217,52 @@ wavebit.rcwa(psi=0.4, chi=np.pi / 4, ...) # circular ``` ::: -The current `polarization='s'` / `'p'` at `phi=0` runs the conical formulation: on a one-row -raster grating its amplitudes and efficiencies match `tetm='te'` / `'tm'` to 5e-15 (measured -2026-10-15), without the pure path's speed. +### The words are a declaration, the pair is a value + +`polarization='s'` or `'p'` says the polarization is fixed: not a parameter anyone sweeps or +differentiates. That is what the pure 1D TE/TM formulation needs, so since 2026-11-03 a word +takes it wherever it answers the same light: every layer a one-row raster (`Layer(epsilon=...)` +of shape `(1, Nx)`), an axis-aligned cell, and `phi` and `ky` absent or a zero that nothing +traces or differentiates. `'s'` is TE and `'p'` is TM. Before that date every one-row raster +solved on the conical formulation. + +Everywhere else a word solves exactly as the pair it names, `(1, 0)` or `(0, 1)`, and nothing is +refused: + +- a skewed cell, whose grooves do not run along y, solves on the 2D path; +- a nonzero `phi` or `ky`, or one that is traced or differentiated (`jax.jit`, `jax.grad`, + `vmap`, a `requires_grad` tensor), solves on the conical path, the one that carries an azimuth + and its derivative; +- `Numerics(formulation='first_order')` solves on the 2D layout it is written on; +- `Outputs(scattering=True)` stays conical, because the requested matrix keeps both + polarizations' ports (`2N x 2N`, channel `'sp'`) and the pure path's has one polarization's. + +A Jones pair is a value, and always keeps the conical path, `(1, 0)` included: it may be swept +or differentiated. The two formulations agree to round-off, 2e-13 relative or better on every +amplitude, efficiency, power, field and boundary field across the cases in +`tests/api/test_unified_rcwa_pure_1d.py`, including their JAX and PyTorch gradients. The pure +one is faster and smaller; see "Cost" below. `check(...).resolved['solve_path']` names the +path a run will take, and `result.diagnostics['solve_path']` the one each wavelength took: +`'1d_te'`, `'1d_tm'`, `'1d_conical'` or `'2d'`. The first is read from the concrete problem +`check` or `prepare` was given, so a prepared run that traces `phi` -- inside `jax.jit`, say -- +reports the plan `'1d_te'` and records `'1d_conical'`, the path it took. + +#### Cost + +Measured on 2026-11-03 through the public `rcwa.run`: NumPy on a 10-core laptop, a +1000 nm cell holding a 560 nm bar of permittivity 12.25 in air, 200 nm thick, at 1363 nm and +normal incidence, a 50-pixel raster with `fourier_coefficients='integrated'`, `'p'` against +`(0, 1)`. Times are the best of 7, 3 and 1 runs. Memory is the process's peak footprint above its +baseline during the run. + +| `orders` | `'p'` (pure TM) | `(0, 1)` (conical) | time ratio | peak memory, pure / conical | `sum(R)` apart | +| --- | --- | --- | --- | --- | --- | +| `(50, 0)` | 0.017 s | 0.045 s | 2.7x | 1.8 / 26 MB | 5e-13 | +| `(200, 0)` | 0.29 s | 1.22 s | 4.1x | 31 / 251 MB | 8e-12 | +| `(800, 0)` | 10.8 s | 57.8 s | 5.3x | 0.39 / 5.9 GB | 4e-10 | + +The conical layout solves matrices twice the side of the pure one's, so the gap in time and +memory widens with the truncation, and so does the last-digit difference between the two answers. :::{admonition} Earlier interface :class: earlier-interface diff --git a/docs/guide/rcwa/numerics.md b/docs/guide/rcwa/numerics.md index a54a90f7..fbee20e3 100644 --- a/docs/guide/rcwa/numerics.md +++ b/docs/guide/rcwa/numerics.md @@ -11,6 +11,8 @@ numerics = rcwa.Numerics( `orders=(Mu, Mv)` retains orders −Mu…Mu along the first lattice vector and −Mv…Mv along the second; it is always an explicit pair of nonnegative integers (`numerics.orders` refuses a bare integer). Result arrays place the wavelength axis first, then v, then u. Algorithms are `"smm"` or `"etm"`; formulations are `"auto"`, `"first_order"` or `"second_order"`. Factorization choices are explicit, with support checked against the material and geometry. The example's small settings are for a quick calculation, not a convergence recommendation. +The solve path is not a setting: the structure and the source decide whether a run takes the pure 1D TE/TM path, the 1D conical path or the 2D one ([the light guide](light.md#the-words-are-a-declaration-the-pair-is-a-value) says how). `rcwa.check(...).resolved['solve_path']` names the path a run will take, and `result.diagnostics['solve_path']` records the one each wavelength took: `'1d_te'`, `'1d_tm'`, `'1d_conical'` or `'2d'`. + Rayleigh detuning is opt-in: `detuning=None` requests the exact physical input. A positive `detuning` declares a relative step and the result records the changed coordinates. Never interpret a detuned result as the unmodified requested wavelength without reading that metadata. `topology_check='device'` uses available array guards and refuses other traced geometry. diff --git a/docs/guide/rcwa/scene.md b/docs/guide/rcwa/scene.md index b37a05a0..68c3270e 100644 --- a/docs/guide/rcwa/scene.md +++ b/docs/guide/rcwa/scene.md @@ -217,4 +217,6 @@ result = rcwa.run(problem, numerics=rcwa.Numerics(orders=(10, 0))) Each value has one place, so there is no second spelling for it to collide with: the cell, the half-spaces and the layers in `Stack`, the wavelength in `Excitation`, direction and polarization in the source, truncation in `Numerics`. The earlier `tetm='te'` is -`polarization='s'` at `phi=0.` here. +`polarization='s'` here: on a one-row `raster` it takes the same pure TE/TM path +(`result.diagnostics['solve_path']` reads `('1d_te',)`), where a Jones pair would take the +conical one. See [the light guide](light.md#the-words-are-a-declaration-the-pair-is-a-value). diff --git a/docs/learn/fundamentals/index.md b/docs/learn/fundamentals/index.md index 11c50979..67c21d07 100644 --- a/docs/learn/fundamentals/index.md +++ b/docs/learn/fundamentals/index.md @@ -3,7 +3,8 @@ :::{note} The downloadable teaching scripts, the text and the snippets of every chapter use the current `rcwa.Problem` / `rcwa.run` API. A few passages describe features the current API does not have -yet, such as the `psi`/`chi` polarization angles and the pure 1D `tetm` path; they sit in boxes +yet, such as the `psi`/`chi` polarization angles and the earlier `tetm` declaration (whose pure 1D +path the words `"s"` and `"p"` now take); they sit in boxes titled "Earlier interface" (retired with `wavebit.compat` on 2026-10-15). Some figures keep the earlier labels, and the chapters that show them name the current equivalents. Compare the [current getting-started guide](../../getting-started/index.md) and diff --git a/docs/learn/fundamentals/light.md b/docs/learn/fundamentals/light.md index da28ad7c..5a94cc18 100644 --- a/docs/learn/fundamentals/light.md +++ b/docs/learn/fundamentals/light.md @@ -257,12 +257,12 @@ te_result = rcwa.run(grating("s"), numerics=numerics) # TE tm_result = rcwa.run(grating("p"), numerics=numerics) # TM ``` -These run the general (conical) formulation, which also handles 2D patterns and nonzero `ky`. On a one-row raster grating at `phi=0` its amplitudes and efficiencies match a dedicated TE/TM solve to 5e-15 (measured 2026-10-15); it does not have that solve's speed. +On a one-row raster grating like this one, with `phi` and `ky` absent or zero, these run the dedicated TE/TM solve: the words `"s"` and `"p"` declare the polarization fixed, and that is what the solve needs (since 2026-11-03; `result.diagnostics["solve_path"]` reads `("1d_te",)` or `("1d_tm",)`). Written as a Jones pair, `(1, 0)` or `(0, 1)`, the same light runs the general (conical) formulation, which also handles 2D patterns and nonzero `ky`; the two agree to round-off, and the dedicated solve is several times faster on a fine truncation. A skewed cell, a nonzero or differentiated azimuth, and a requested scattering matrix keep the general formulation for the words too. :::{admonition} Earlier interface :class: earlier-interface -The earlier interface (retired with `wavebit.compat` on 2026-10-15) had a dedicated declaration for this case: `tetm="te"` or `tetm="tm"` solved TE or TM alone, 2 to 4 times faster than the conical formulation. It required `ky` to be zero (omitted, or explicitly `ky=0.`); nonzero `ky` raised an error. Because it declared a problem with no polarization degree of freedom, it could not be combined with `psi`, `chi` or `polarization`. The current API has no `tetm` declaration yet. +The earlier interface (retired with `wavebit.compat` on 2026-10-15) had a dedicated declaration for this case: `tetm="te"` or `tetm="tm"` solved TE or TM alone, 2 to 4 times faster than the conical formulation. It required `ky` to be zero (omitted, or explicitly `ky=0.`); nonzero `ky` raised an error. Because it declared a problem with no polarization degree of freedom, it could not be combined with `psi`, `chi` or `polarization`. In the current API the words `"s"` and `"p"` take that solve where it applies, and refuse nothing where it does not. ```python # Earlier interface only; runs on a historical checkout. diff --git a/docs/learn/fundamentals/overview.md b/docs/learn/fundamentals/overview.md index 582f58f2..2bc254ea 100644 --- a/docs/learn/fundamentals/overview.md +++ b/docs/learn/fundamentals/overview.md @@ -3,7 +3,8 @@ :::{note} The downloadable teaching scripts, the text and the snippets of every chapter use the current `rcwa.Problem` / `rcwa.run` API. A few passages describe features the current API does not have -yet, such as the `psi`/`chi` polarization angles and the pure 1D `tetm` path; they sit in boxes +yet, such as the `psi`/`chi` polarization angles and the earlier `tetm` declaration (whose pure 1D +path the words `"s"` and `"p"` now take); they sit in boxes titled "Earlier interface" (retired with `wavebit.compat` on 2026-10-15). Some figures keep the earlier labels, and the chapters that show them name the current equivalents. Compare the [current getting-started guide](../../getting-started/index.md) and diff --git a/docs/locales/ko/pages/development/source-map.md b/docs/locales/ko/pages/development/source-map.md index 6d68e083..aff92125 100644 --- a/docs/locales/ko/pages/development/source-map.md +++ b/docs/locales/ko/pages/development/source-map.md @@ -82,6 +82,7 @@ | --- | --- | --- | | 재료 응답 또는 SI 변환 | [test_common_material_models.py](https://github.com/fieldalgo/wavebit/blob/main/tests/materials/test_common_material_models.py), [test_lifecycle.py](https://github.com/fieldalgo/wavebit/blob/main/tests/fdtd/test_lifecycle.py) | 다른 기하 단위로 표현한 동일 문제, SI 전도도·극점, 지원하는 인과적 재료와 미분 연결을 확인합니다. | | RCWA 준비, 실행 시점 값 또는 결과 소유권 | [test_unified_rcwa.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa.py) | 재사용과 새 준비의 일치, 독립 Fresnel 장·전력, 변환 중 잘못된 입력, 이전 결과의 독립성을 검사합니다. | +| RCWA 풀이 경로(순수 1D TE/TM, conical, 2D) | [test_unified_rcwa_pure_1d.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa_pure_1d.py), [test_tetm_off_the_pure_path.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_tetm_off_the_pure_path.py) | 편광 낱말과 그 낱말이 가리키는 쌍을 비교합니다. 낱말이 순수 경로를 타면 모든 관측량과 미분이 반올림 오차 안에서 같고, 물러나면 비트 단위로 같으며, 보고한 경로가 실제로 호출된 커널과 맞는지 확인합니다. | | 결맞는 소스와 절대 전력 | [test_unified_rcwa_sources.py](https://github.com/fieldalgo/wavebit/blob/main/tests/api/test_unified_rcwa_sources.py) | SI 소스 세기, 유한한 관측 영역과 독립적으로 평가한 장을 확인합니다. | | FDTD 영역, 소스·모니터 또는 재사용 | [test_lifecycle.py](https://github.com/fieldalgo/wavebit/blob/main/tests/fdtd/test_lifecycle.py) | 물리 영역 바깥 PML, 고정된 실험실 좌표, 소스·재료 값 재사용, 실제 실행 시간 제약을 확인합니다. | | RCWA 모드 또는 장 복원 연산 | [test_smm_field_state.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_smm_field_state.py), [test_lossy_media.py](https://github.com/fieldalgo/wavebit/blob/main/tests/rcwa/test_lossy_media.py) | 솔버 golden을 다시 만드는 대신 독립적인 내부 장과 계면 응답으로 확인합니다. | diff --git a/docs/locales/ko/pages/learn/fundamentals/index.md b/docs/locales/ko/pages/learn/fundamentals/index.md index 0b7b3166..9db945df 100644 --- a/docs/locales/ko/pages/learn/fundamentals/index.md +++ b/docs/locales/ko/pages/learn/fundamentals/index.md @@ -2,7 +2,7 @@ :::{note} 다운로드 가능한 교육 스크립트와 모든 장의 본문·코드는 현재 `rcwa.Problem` / `rcwa.run` API를 씁니다. -`psi`/`chi` 편광 각도나 순수 1D `tetm` 경로처럼 현재 API에 아직 없는 기능을 설명하는 몇 대목은 +`psi`/`chi` 편광 각도나 이전 `tetm` 선언(그 순수 1D 경로는 이제 낱말 `"s"`·`"p"`가 탑니다)처럼 현재 API에 아직 없는 기능을 설명하는 몇 대목은 ‘이전 인터페이스’ 상자 안에 있습니다(2026-10-15에 `wavebit.compat`와 함께 퇴역). 일부 그림에는 이전 이름이 남아 있으며, 그런 그림이 있는 장에서 현재 이름을 알려 줍니다. [현재 시작 안내](../../getting-started/index.md)와 [API 이행 안내](../../migration/index.md)를 함께 참고하세요. 과거 측정 파일과 원본 출처는 그대로 둡니다. diff --git a/docs/locales/ko/pages/learn/fundamentals/light.md b/docs/locales/ko/pages/learn/fundamentals/light.md index 954a6847..ac46bc4a 100644 --- a/docs/locales/ko/pages/learn/fundamentals/light.md +++ b/docs/locales/ko/pages/learn/fundamentals/light.md @@ -257,12 +257,12 @@ te_result = rcwa.run(grating("s"), numerics=numerics) # TE tm_result = rcwa.run(grating("p"), numerics=numerics) # TM ``` -이 계산은 2D 패턴과 0이 아닌 `ky`도 다루는 일반(conical) 정식화로 풉니다. `phi=0`인 한 행 래스터 격자에서는 진폭과 효율이 TE/TM 전용 풀이와 5e-15 이내로 같습니다(2026-10-15 측정). 다만 전용 풀이의 속도는 얻지 못합니다. +이 격자처럼 한 행 래스터이고 `phi`와 `ky`가 없거나 0이면, 이 계산은 TE/TM 전용 풀이로 풉니다. 낱말 `"s"`와 `"p"`는 편광이 고정됐다는 선언이고, 전용 풀이에 필요한 것이 바로 그것입니다(2026-11-03부터이며, `result.diagnostics["solve_path"]`는 `("1d_te",)` 또는 `("1d_tm",)`입니다). 같은 빛을 Jones 쌍 `(1, 0)`이나 `(0, 1)`로 쓰면 2D 패턴과 0이 아닌 `ky`도 다루는 일반(conical) 정식화로 풉니다. 두 결과는 반올림 오차 안에서 같고, 차수를 많이 쓸수록 전용 풀이가 몇 배 빠릅니다. 기울어진 셀, 0이 아니거나 미분하는 방위각, 산란 행렬 요청에서는 낱말도 일반 정식화로 풉니다. :::{admonition} 이전 인터페이스 :class: earlier-interface -이전 인터페이스(2026-10-15에 `wavebit.compat`와 함께 퇴역)에는 이 경우를 위한 전용 선언이 있었습니다. `tetm="te"` 또는 `tetm="tm"`은 TE나 TM 하나만 풀었고, conical 정식화보다 2–4배 빨랐습니다. `ky`는 0이어야 했으며(생략하거나 `ky=0.`을 직접 적음), 0이 아닌 `ky`는 오류였습니다. 편광 자유도가 없는 문제를 선언하므로 `psi`, `chi`, `polarization`과 함께 쓸 수 없었습니다. 현재 API에는 아직 `tetm` 선언이 없습니다. +이전 인터페이스(2026-10-15에 `wavebit.compat`와 함께 퇴역)에는 이 경우를 위한 전용 선언이 있었습니다. `tetm="te"` 또는 `tetm="tm"`은 TE나 TM 하나만 풀었고, conical 정식화보다 2–4배 빨랐습니다. `ky`는 0이어야 했으며(생략하거나 `ky=0.`을 직접 적음), 0이 아닌 `ky`는 오류였습니다. 편광 자유도가 없는 문제를 선언하므로 `psi`, `chi`, `polarization`과 함께 쓸 수 없었습니다. 현재 API에서는 낱말 `"s"`와 `"p"`가 그 풀이를 쓸 수 있는 곳에서 쓰고, 쓸 수 없는 곳에서도 아무것도 거절하지 않습니다. ```python # 이전 인터페이스 전용. 과거 체크아웃에서만 실행됩니다. diff --git a/docs/locales/ko/pages/learn/fundamentals/overview.md b/docs/locales/ko/pages/learn/fundamentals/overview.md index f51ace36..28a8cb87 100644 --- a/docs/locales/ko/pages/learn/fundamentals/overview.md +++ b/docs/locales/ko/pages/learn/fundamentals/overview.md @@ -2,7 +2,7 @@ :::{note} 다운로드 가능한 교육 스크립트와 모든 장의 본문·코드는 현재 `rcwa.Problem` / `rcwa.run` API를 씁니다. -`psi`/`chi` 편광 각도나 순수 1D `tetm` 경로처럼 현재 API에 아직 없는 기능을 설명하는 몇 대목은 +`psi`/`chi` 편광 각도나 이전 `tetm` 선언(그 순수 1D 경로는 이제 낱말 `"s"`·`"p"`가 탑니다)처럼 현재 API에 아직 없는 기능을 설명하는 몇 대목은 ‘이전 인터페이스’ 상자 안에 있습니다(2026-10-15에 `wavebit.compat`와 함께 퇴역). 일부 그림에는 이전 이름이 남아 있으며, 그런 그림이 있는 장에서 현재 이름을 알려 줍니다. [현재 시작 안내](../../getting-started/index.md)와 [API 이행 안내](../../migration/index.md)를 함께 참고하세요. 과거 측정 파일과 원본 출처는 그대로 둡니다. diff --git a/docs/reference/common/conventions.md b/docs/reference/common/conventions.md index 55e94a74..bd734a16 100644 --- a/docs/reference/common/conventions.md +++ b/docs/reference/common/conventions.md @@ -274,7 +274,7 @@ All angles are **radians**. | `phi` | azimuth of the incident light's direction, from the `+x` axis | | `psi` | polarization: the angle of **E** from the p direction in the plane transverse to **k**. `psi = pi/2` is pure s (TE-like), `psi = 0` pure p (TM-like). **Real** — a complex `psi` is refused, see below | | `chi` | ellipticity angle: `tan(chi)` is the polarization ellipse's minor/major axis ratio and its sign is the rotation direction (`chi > 0` turns like `(1, -1j)`; see "Which way the circular states rotate" below). `chi = 0` is linear (the `psi` spelling), `|chi| = pi/4` circular. Real, like `psi`. Naming it selects the `(psi, chi)` chart — see below | -| `tetm` (1D only) | `'te'` / `'tm'` — a **declaration**, not an angle: it says this problem has no polarization and no azimuth degree of freedom, which is what makes the pure 1D TE/TM formulation legal (2-4x faster than conical, measured). Refused on a 2D cell, with `truncation_order[1] != 0`, or beside a nonzero `phi` / a `psi` / a `polarization`; and, by name, on every solve that cannot take the pure 1D path — a shape structure, a tensor material, a z row handed to `solve_conv`, `formulation='first_order'`, or no structure after anything but a pure solve (a solver holding none keeps the layout its last solve took, 2D if it never solved, and a kept pure layout re-checks `phi`, `ky`, the cell's skew and `truncation_order[1]` as the raster that chose it was checked). Since 2026-09-27; before, all but the tensor raised a `TypeError` from deep in the solve, or — with `tetm` assigned after construction — answered for s light whichever it named. `psi` reads back `None` under it — a declared-away degree of freedom is absent, not defaulted. It lights and reports exactly what `polarization=(1, 0)` / `(0, 1)` at `phi = 0` does, to ~1e-14: the same amplitudes in the same per-order s/p frame, the same `excitation` weights, the same unit-E incident field; its scattering matrix is the conical one's s-s / p-p *block*, in the same port variables. What it does not carry is the other polarization: `s_inc` / `p_inc` are `None`, and the cross-polarized amplitudes are exact zeros where the conical path has ~1e-16. Since 2026-09-27; before, it quoted `E_y` / `H_y` on every order (a sign apart on every order across `kx = 0` from the incident one), was lit by `+y` at any angle, and its TM input was a unit `H_y` | +| `tetm` (1D only) | `'te'` / `'tm'` — a **declaration**, not an angle: it says this problem has no polarization and no azimuth degree of freedom, which is what makes the pure 1D TE/TM formulation legal (2-4x faster than conical, measured). Refused on a 2D cell, with `truncation_order[1] != 0`, or beside a nonzero `phi` / a `psi` / a `polarization`; and, by name, on every solve that cannot take the pure 1D path — a shape structure, a tensor material, a z row handed to `solve_conv`, `formulation='first_order'`, or no structure after anything but a pure solve (a solver holding none keeps the layout its last solve took, 2D if it never solved, and a kept pure layout re-checks `phi`, `ky`, the cell's skew and `truncation_order[1]` as the raster that chose it was checked). Since 2026-09-27; before, all but the tensor raised a `TypeError` from deep in the solve, or — with `tetm` assigned after construction — answered for s light whichever it named. `psi` reads back `None` under it — a declared-away degree of freedom is absent, not defaulted. It lights and reports exactly what `polarization=(1, 0)` / `(0, 1)` at `phi = 0` does, to ~1e-14: the same amplitudes in the same per-order s/p frame, the same `excitation` weights, the same unit-E incident field; its scattering matrix is the conical one's s-s / p-p *block*, in the same port variables. What it does not carry is the other polarization: `s_inc` / `p_inc` are `None`, and the cross-polarized amplitudes are exact zeros where the conical path has ~1e-16. Since 2026-09-27; before, it quoted `E_y` / `H_y` on every order (a sign apart on every order across `kx = 0` from the incident one), was lit by `+y` at any angle, and its TM input was a unit `H_y`. The current `PlaneWave(polarization='s')` / `'p'` takes the same pure path wherever it answers the same light, and is never refused (since 2026-11-03; see [the words are a declaration](../../guide/rcwa/light.md#the-words-are-a-declaration-the-pair-is-a-value)) | ### Real-angle incidence diff --git a/docs/reference/keywords.md b/docs/reference/keywords.md index aa036b92..ca689131 100644 --- a/docs/reference/keywords.md +++ b/docs/reference/keywords.md @@ -393,9 +393,9 @@ amplitudes and scattering matrices agree to round-off (3e-13 or better). The | `incidence_side` | `PlaneWave(side=...)` | | | `theta`, `phi` | `PlaneWave(theta=..., phi=...)` | | | `kx`, `ky` | `PlaneWave(kx=..., ky=...)` | | -| `polarization` | `PlaneWave(polarization=...)` | also `'s'` and `'p'` | +| `polarization` | `PlaneWave(polarization=...)` | also `'s'` and `'p'`, which declare the polarization fixed (the `tetm` row) | | `psi`, `chi` | **none** | the same light is the Jones pair `polarization=(sin ψ cos χ + i cos ψ sin χ, cos ψ cos χ − i sin ψ sin χ)` (`chi = 0` is `(sin ψ, cos ψ)`), which agrees to round-off; the two-angle chart itself, without the pair's redundant scale and phase, has no field | -| `tetm` | **none** | no field declares the polarization away, so a one-row raster always solves on the conical path. The numbers are the same (to 3e-13), but in one measurement at `orders=(10, 0)` to `(30, 0)` the run was 2.2 to 2.9 times slower than the earlier pure TE/TM path | +| `tetm` | `PlaneWave(polarization='s')` / `'p'` | the word is the declaration: `'s'` is `'te'` and `'p'` is `'tm'`. On a one-row raster in an axis-aligned cell, with `phi` and `ky` absent or a zero nothing traces or differentiates, the run takes the pure TE/TM path, and `report.resolved['solve_path']` and `result.diagnostics['solve_path']` say `'1d_te'` / `'1d_tm'`. Unlike `tetm` the word refuses nothing: a skewed cell, a nonzero or traced `phi` or `ky`, `formulation='first_order'` and `Outputs(scattering=True)` solve as the Jones pair it names, on the conical or 2D path. A Jones pair always solves as a pair. Since 2026-11-03; before, every one-row raster solved on the conical path. The two paths agree to round-off (2e-13 relative over `tests/api/test_unified_rcwa_pure_1d.py`) | | `excitation`, an order map | `IncidentOrders(weights=(((m, n), weight), ...))` | `{0: 1.0}` on 1D becomes `((0, 0), 1.0)`. A coherent-source run: its total boundary power agrees with the earlier per-order sums, but **per-order efficiencies have no current form** (`rcwa.coherent_diffraction`) | | `excitation`, a side-keyed map or a source list | `Excitation(sources=(...))` | | | `stack=`, `incidence=` | `rcwa.Stack`, `rcwa.PlaneWave` | | @@ -409,6 +409,8 @@ amplitudes and scattering matrices agree to round-off (3e-13 or better). The | `power_region` | `Outputs(power_region=...)` | | | `length_unit` | `Problem(units=Units(length=...))` | one of `'m'`, `'mm'`, `'um'`, `'nm'` rather than a free factor | -Three earlier keywords have no current home: `psi`, `chi` and `tetm`. The order-map -`excitation` has one, without its per-order efficiencies. A problem that needs one of those -keeps the [earlier interface](../migration/legacy-api.md) for now. +Two earlier keywords have no current home: `psi` and `chi`, whose light is the Jones pair in +their row. The order-map `excitation` has one, without its per-order efficiencies. The earlier +interface that had them, `wavebit.compat`, was retired on 2026-10-15; it survives only in the +historical checkout the [earlier API reference](../migration/legacy-api.md) names, so it is a +record to read rather than a place to run a problem.