Metadata-Version: 2.4
Name: photonfdtd
Version: 0.12.0
Summary: Local-only FDTD + waveguide mode solver for integrated photonics (NumPy core; optional JAX and compiled Rust/CUDA backends).
Project-URL: Homepage, https://github.com/ngpaladi/photonfdtd
Project-URL: Issues, https://github.com/ngpaladi/photonfdtd/issues
Author: ngpaladi
License: MIT
License-File: LICENSE
Keywords: electromagnetics,fdtd,mode-solver,photonics,simulation
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: License :: OSI Approved :: MIT License
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3 :: Only
Classifier: Topic :: Scientific/Engineering :: Physics
Requires-Python: >=3.10
Requires-Dist: matplotlib>=3.7
Requires-Dist: numpy>=1.24
Requires-Dist: scipy>=1.10
Provides-Extra: compress
Requires-Dist: zstandard>=0.21; extra == 'compress'
Provides-Extra: docs
Requires-Dist: sphinx-rtd-theme>=3; extra == 'docs'
Requires-Dist: sphinx>=7; extra == 'docs'
Provides-Extra: jax
Requires-Dist: jax>=0.4; extra == 'jax'
Provides-Extra: test
Requires-Dist: numba>=0.59; extra == 'test'
Requires-Dist: pytest>=7.4; extra == 'test'
Requires-Dist: zstandard>=0.21; extra == 'test'
Description-Content-Type: text/markdown

# photonfdtd

A small, fully-local open-source FDTD + waveguide mode solver written in
Python/NumPy. It implements a Yee-grid FDTD time-stepper with CPML
absorbing boundaries, dispersive materials, anisotropic subpixel smoothing,
and a full-vectorial waveguide mode solver, with an API intentionally similar
in spirit to Tidy3D.

This is alpha-stage software (v0.1). The pieces that exist are tested and
correct; the pieces that don't, don't. See *Status* below.

## Install

```bash
pip install photonfdtd
```

or from source, either directly from GitHub:

```bash
pip install "git+https://github.com/ngpaladi/photonfdtd"
```

or from a checkout:

```bash
git clone https://github.com/ngpaladi/photonfdtd
cd photonfdtd
pip install -e .
```

The PyPI wheel is pure Python (NumPy core; `photonfdtd[jax]` for the
differentiable XLA backend). The optional compiled Rust / CUDA backends
(`backend="rust"`, `backend="rust-cuda"`) are built separately from the Rust
sources shipped in the sdist — see [`benchmarks/README.md`](benchmarks/README.md).

### Releasing (maintainers)

Releases publish to PyPI automatically via
[`.github/workflows/publish.yml`](.github/workflows/publish.yml) using PyPI
**Trusted Publishing** (OIDC — no API token stored). To cut a release:

1. Bump the version in `pyproject.toml` **and** `src/photonfdtd/__init__.py`
   (kept in sync by `tests/test_version.py`).
2. Publish a GitHub Release tagged `vX.Y.Z` matching that version. The workflow
   runs the tests, guards tag-vs-version, builds the sdist + wheel, and
   publishes.

One-time setup before the first release: on PyPI, add a *trusted publisher* for
the project (Publishing → add a GitHub Actions publisher) with owner `ngpaladi`,
repo `photonfdtd`, workflow `publish.yml`, environment `pypi`. Until the project
exists, use PyPI's *pending publisher* form with the same values.

## Quick start

### A point dipole radiating in 2D vacuum

```python
import photonfdtd as pf

lam0 = 1.0e-6
freq0 = pf.C_0 / lam0
dx = lam0 / 20

grid = pf.Grid(size=(4e-6, 4e-6), cell_size=dx, pml_layers=(12, 12, 0))
src = pf.PointDipole(
    position=(0.0, 0.0),
    component="Ez",
    waveform=pf.GaussianPulse(freq0=freq0, fwhm=10e-15),
)
mon = pf.FieldMonitor(name="snap", components=("Ez",), interval=20)

sim = pf.Simulation(grid, sources=[src], monitors=[mon], run_time=200e-15)
result = sim.run()
```

`result.fields["snap"]["Ez"]` is a `(n_frames, ny, nz)` array of field snapshots.

### Solving a slab waveguide mode

```python
import photonfdtd as pf

clad = pf.Medium.from_index(1.0)
core = pf.Medium.from_index(2.0)
slab = pf.Box(center=(0.0, 0.0), size=(0.5e-6, 0.3e-6), medium=core)

ms = pf.ModeSolver(
    size=(4e-6, 3e-6),
    cell_size=20e-9,
    structures=[slab],
    wavelength=1.55e-6,
    num_modes=2,
)
result = ms.solve()
print(result.n_eff)            # array of effective indices
```

## What v0.2 actually does

- **Yee-grid FDTD** in 1D / 2D / 3D, with the Courant-stable time step
  selected automatically.
- **CPML** absorbing boundaries (Roden & Gedney 2000, kappa = 1) on any
  number of axes.
- **Isotropic dielectric media, non-dispersive or dispersive.** Non-dispersive
  media stamp a single `eps_r` per cell; dispersive media (`DispersiveMedium`,
  with `.lorentz` / `.drude` / `.sellmeier` constructors) carry Lorentz, Drude
  and Sellmeier poles advanced by the auxiliary-differential-equation (ADE)
  method during stepping. A cited material library (`pf.silica()`,
  `pf.silicon()`, `pf.silicon_nitride()`, `pf.lithium_niobate()`, `pf.gold()`,
  `pf.silver()`) reproduces each source's published index — see
  [`docs/material_data.md`](docs/material_data.md).
- **Anisotropic subpixel smoothing** (`Simulation(subpixel=True)`): partially
  filled interface cells get an effective diagonal permittivity tensor
  (harmonic mean normal to the interface, arithmetic mean tangential), removing
  the staircase error and restoring second-order boundary accuracy
  (Farjadpour et al. 2006). See [`docs/accuracy.rst`](docs/accuracy.rst).
- **Geometry primitives**: axis-aligned `Box` and arbitrary-polygon
  `PolySlab` (a polygon in xy extruded between two z bounds).
- **Sources**: soft point-dipole (`PointDipole`), distributed line/area
  mode injection (`ModeSource`), an energy-normalised `SinglePhotonSource`
  whose amplitude is set so the launched wavepacket carries exactly
  $h\,\nu$ of total electromagnetic energy, and a moving-charge
  `ChargedParticle` current source that emits **Cherenkov radiation** when
  it outruns the local phase velocity (see `examples/04_cherenkov.py`).
- **Monitors**: time-domain field snapshots (`FieldMonitor`), Poynting flux
  (`FluxMonitor`), and a frequency-domain `DFTMonitor` that accumulates a
  running Fourier transform at chosen frequencies so storage scales with the
  number of frequencies rather than the number of timesteps (routinely
  50-1000x smaller than an equivalent time-domain monitor, exact at those
  frequencies).
- **Memory efficiency for large volumes**: per-array `float32`/`float64`
  precision, the permittivity grid released during stepping, a shared curl
  scratch buffer that avoids per-step full-domain temporaries, and
  `FieldMonitor(compression=...)`, which streams snapshots to disk as
  per-frame-scaled, quantised, compressed blocks - keeping RAM flat regardless
  of recording length and shrinking stored data ~10-30x (8-bit) or ~5-7x
  (16-bit) versus an uncompressed float64 monitor.
- **Out-of-core stepping** (`sim.run(out_of_core=True, tile_cells=...)`): for a
  volume whose fields do not fit in RAM, the six field arrays, `ce_field` and
  the CPML state are memory-mapped to disk and each timestep sweeps the domain
  in slabs along one axis, so **peak RAM is bounded by the tile size, not the
  grid**. Reproduces the in-core result to machine precision. NumPy backend,
  point sources and `FieldMonitor` (incl. `compression=`) supported; see
  `photonfdtd.outofcore`.
- **Differentiable PIC inverse design.** A substrate/etched-core/cladding stack
  driven by a 2-D etch density `ρ(x,y)` (`EtchedCore`) maps — conic filter →
  `tanh` projection → anisotropic subpixel sidewalls → 3-D ε — entirely in JAX,
  so `value_and_grad_density` returns `d(loss)/d(ρ)` through the full time
  evolution (topology optimization; gradient matches finite differences to
  ~1e-8). Ports are launched one-way with `UniModeSource` (equivalence-principle
  current sheets). See [`docs/accuracy.rst`](docs/accuracy.rst).
- **Mode-decomposition S-parameters** (`s_parameters`, `mode_amplitudes`):
  project the frequency-domain fields on a memory-light port plane (any axis via
  `DFTMonitor(plane_axis=..., plane_position=...)`) onto a solved mode to get
  complex forward/backward modal amplitudes and hence S-parameters. A mode down
  a straight lossless guide gives `|S21|≈1`; the overlap is backend-agnostic
  (NumPy/CuPy/JAX) so `|S_ij|²` is differentiable through the adjoint — the PIC
  inverse-design objective. See [`docs/accuracy.rst`](docs/accuracy.rst).
- **Full-vectorial waveguide mode solver** (`ModeSolver`): a finite-difference
  eigenproblem for the transverse fields on the Yee grid that distinguishes TE
  from TM and captures high-index-contrast boundary effects. Validated against
  the exact symmetric-slab dispersion to better than 0.01% on both TE0 and TM0
  (Zhu & Brown 2002; Fallahkhair et al. 2008).
- **gdsfactory adapter** (`from_gdsfactory`) that reads a layout
  `Component`, maps its layers onto user-supplied materials, and returns a
  pre-built `Simulation`.
- **Optional acceleration backends** for time-stepping: a Numba JIT path and
  a CuPy GPU path (`use_gpu=True`). The GPU path uses only generic CuPy array
  ops, so it runs on either an **NVIDIA** GPU (CuPy's CUDA build,
  `cupy-cuda12x`) or an **AMD** GPU (CuPy's ROCm build, `cupy-rocm-5-0`).
  Both are optional; without them the plain NumPy core runs.
- **Differentiable JAX backend** (`use_jax=True`): a pure-functional
  reimplementation of the Yee + CPML step under `jax.lax.scan`, JIT-compiled
  through XLA (CPU/GPU/TPU from one path) and matching the in-core solver to
  machine precision. Because the stepper is a pure function of the
  permittivity, `pf.jax_value_and_grad_eps(sim, loss)` returns the gradient of
  any monitor-based scalar w.r.t. `eps_r` — a **time-domain adjoint for
  gradient-based inverse design / topology optimization**. See
  `photonfdtd.jaxbackend`. (`pip install "photonfdtd[jax]"`.)

### Running on GPU

**JAX is the GPU backend.** With `use_jax=True` the forward run and the
differentiable / reversible adjoints run on the GPU automatically once a CUDA
`jaxlib` is present — no code change: `pip install "jax[cuda12]"`. Validated on
an NVIDIA RTX 4080 (CUDA 12): forward matches NumPy to ~1e-16, and both adjoints
run on-device. If JAX falls back to CPU with a "cuSPARSE not found"-style error,
the nvidia CUDA libraries just aren't on the loader path — add them, e.g.
`export LD_LIBRARY_PATH=$(python -c "import nvidia,glob,os;print(':'.join(glob.glob(os.path.dirname(nvidia.__file__)+'/*/lib')))")`.
Large differentiable runs fit in GPU memory thanks to gradient checkpointing and
the O(1)-memory reversible adjoint (`photonfdtd.reversible`).

The **CuPy backend** (`use_gpu=True`) is a deprecated legacy path — superseded by
JAX for GPU work, but retained (optional; `pip install cupy-cuda12x`) for
AMD/ROCm GPUs and for **GPU/host/disk out-of-core** (`use_gpu=True` +
`run(out_of_core=True)`), where fields live on disk and each tile is processed on
the GPU so peak device memory is one tile — a volume larger than GPU memory still
runs (JAX cannot stream to disk).

## What v0.2 does *not* do (yet)

Items below are intentional out-of-scope and tracked in the issue
tracker:

- Dispersive media (Lorentz / Drude / Debye) - planned next, via the
  auxiliary-differential-equation method.
- Full-vectorial 2D mode solving with E and H couplings.
- Anisotropic media.
- Sub-cell averaging for polygon edges.
- Total-field / scattered-field plane wave injection.

## Roadmap

- v0.3: dispersive media (Lorentz pole-residue), full-vectorial mode
  solver, anisotropic materials.
- v0.4: Numba/JAX backend for ~10-100x speed-up.

## License

MIT.
