Skip to content

Implementation status

What is operational

The repository is an executable 3D reference solver and inverse-design loop, not a code skeleton.

  • Three-dimensional first-order Nedelec edge elements on tetrahedra, with PEC truncation and anisotropic transformation-optics PML.
  • A regularized electric-current dipole, emitted-power LDOS, vacuum normalization, material interpolation, compact density filtering, and smooth projection.
  • Assembled CPU and CUDA truth solves, assembled GMRES, and element-local CPU/CUDA matrix-free GMRES. The CUDA action is a fused six-edge-per-tetrahedron kernel.
  • A complex-shifted ILU preconditioner. The shift affects only the approximate inverse; convergence is checked against the unshifted physical Maxwell equation.
  • General conjugate-transpose LDOS adjoints, analytic frequency partials, reciprocal complex-symmetric QNM sensitivities, and exact filter/projection transpose propagation.
  • Explicit non-Hermitian generalized shift-invert QNM acquisition, source-coupled deterministic starts, warm-started Arnoldi, overlap/coupling scoring, PML rejection, and physical pencil residuals.
  • Fixed-frequency acquisition followed by automatic shifted pole tracking, symmetric frequency-band penalty, move limits, backtracking, mode-identity guards, and rejected-solve recovery.
  • Restartable checkpoints that retain density, stage, pole, QNM field, and pole diagnostics.
  • PML, independent-mesh, and driven-linewidth validation commands; equal-discretization backend benchmarks; NPZ, JSON, YAML, and ParaView VTU output.

Verified on the included cases

The automated suite currently contains 16 tests. It checks CPU/CUDA agreement, assembled/matrix-free agreement, complex transpose symmetry, forward and adjoint residuals, QNM backward error, restart continuity, VTU output, and central directional derivatives.

On configs/tiny.yaml:

  • driven and adjoint residuals are approximately 1e-15;
  • the selected QNM is approximately 6.3886946 - 0.1658318i, with Q = 19.2626 and a generalized-pencil residual of 3.5e-15;
  • the best LDOS directional-gradient error is below 1e-7;
  • the analytic frequency-partial error is below 1e-6;
  • QNM and total shifted derivatives pass the 3e-3 simple-mode threshold;
  • a shifted checkpoint resumes with the same tracked pole and continues improving the objective.

On configs/benchmark.yaml (1,296 tetrahedra and 1,206 free edge unknowns), all four selected truth/scale paths agree in emitted power to better than 2e-9 at the configured iterative tolerance. The lighter complex-shifted factors require about 133 CPU or 180 CUDA Krylov steps while maintaining physical residuals below 1e-9. This is a correctness and integration benchmark, not a claim that a consumer RTX GPU is faster at this small size.

Validation results that intentionally fail

The tiny mesh is too coarse to support a scientific device claim. Its PML sweep changes normalized LDOS by about 12%, above the 5% test tolerance. A two-times driven mesh refinement changes normalized LDOS by about 422%. Its driven-response Lorentzian fit gives a Q about 17.7% below the QNM value and R^2 about 0.937. These commands return failure as designed.

A two-times refined tiny solve with a near-exact factor reaches the sparse-factorization wall at only a few thousand edge unknowns. Switching to the lighter complex-shifted factor in configs/validation.yaml finishes the same forward system in roughly 1.2 seconds on the development CPU, at the cost of hundreds of Krylov steps. The GPU matrix-free action itself is not the blocker; controlling iteration count without a heavy assembled coarse inverse is. This is the boundary between the implemented reference backend and the production backend below.

Production-scale boundary

The following are not implemented and must not be inferred from the current backend:

  • MFEM/libCEED partial assembly and an auxiliary-space H(curl) multigrid hierarchy;
  • GPU-resident multilevel coarse solves, QNM correction solves, and optimization state;
  • mode deflation/recycled FGMRES beyond a one-pole initial guess;
  • scalable Jacobi-Davidson or SLEPc QNM continuation;
  • second-order elements, adaptive/unstructured CAD meshing, MPI, or multi-GPU execution;
  • robust eroded/nominal/dilated fabrication optimization.

The current matrix-free operator still builds assembled matrices for the truth residual, QNM solve, and ILU preconditioner. It therefore validates the operator and GPU execution model, but it does not yet satisfy the final near-linear-memory scaling gate. A generic scalar algebraic multigrid experiment was rejected because it did not converge on the indefinite curl-curl system; the next backend needs a real Maxwell auxiliary-space method, not a renamed scalar AMG.

Commands and gates

uv sync --extra test --extra gpu
uv run pytest -q
uv run fem-qnm gradient-check configs/tiny.yaml --include-qnm
uv run fem-qnm benchmark configs/benchmark.yaml --backends cpu_direct cpu_matrix_free_gmres cuda_direct cuda_matrix_free_gmres
uv run fem-qnm validate-pml configs/tiny.yaml
uv run fem-qnm validate-mesh configs/validation.yaml --refinement-factors 1 2 --skip-qnm
uv run fem-qnm validate-linewidth configs/tiny.yaml

Do not report an optimized design until PML, mesh, and linewidth commands all accept on a verification mesh not used for optimization.