Pre-release: benchmarks, new solvers, and docs update#6
Merged
Conversation
Introduces @fdomain_component decorator to define components via an admittance matrix Y(f) rather than time-domain physics F(y), Q(y). - `@fdomain_component(ports=...)` wraps a function `(f, **params) -> Y_matrix` into a CircuitComponent subclass with `_is_fdomain=True` and no internal states. Parameters (excluding `f`) must have defaults. - `ComponentGroup` gains an `is_fdomain` flag set during `compile_netlist`. - Assembly: f-domain groups are skipped in `assemble_residual_only_*` (used by HB time-domain vmap) and evaluated at f=0 in `assemble_system_*` (DC). - HB solver: adds `Y(k·f₀) @ V_k` contributions to R_k for each f-domain group before the irfft, covering all harmonics simultaneously via vmap. - Transient: raises RuntimeError at setup time when f-domain groups are present. - Tests: 15 new tests covering decorator validation, DC, HB accuracy, JIT, grad, and the transient guard. - Notebook: Part 3 demonstrates that CapacitorFD/InductorFD defined via @fdomain_component give numerically identical HB results to the standard time-domain Capacitor/Inductor (max error ~1e-10 V).
- Add quick-reference table comparing @component, @source, @fdomain_component, and sax_component with solver support matrix - Add @fdomain_component section explaining: when to use it, solver behaviour at DC/HB/transient, SkinEffectResistor example, and the equivalence with time-domain reactive components for HB - Expand @component section with an Inductor example showing internal state variables - Reorganise Photonic Components into two subsections: manual @component with s_to_y() vs sax_component for importing SAX library models; clarify why photonic components are time-domain (optical wavelength is a parameter, not a solver variable) - Explain that sax_component and @fdomain_component serve orthogonal physical regimes and cannot be unified
…ify pixi tasks - Switch mkdocs.yml palette default to light mode (scheme: default first) - Fix logo SVG to use fill="currentColor" for automatic dark/light adaptation - Update README and index.md to use logo.svg instead of logo_white.svg - Eliminate docs/nbrun.py and docs/nbconvert_all.py Python helpers; replace with inline pixi tasks (nbrun, nbdocs) using /bin/sh and nbconvert_config.py - Add docs/nbconvert_config.py (Jupyter traitlets config) for flat notebook conversion to docs/examples/ without shell variable quoting issues - Regenerate docs/examples/ as flat structure (no electrical/photonic subdirs) - Update notebooks: replace transparent backgrounds with white, remove dark_background style cells - Simplify GitHub Actions docs.yaml to single pixi run docs-deploy step - Remove unused extra.css logo filter rules and dead mkdocs-jupyter CSS
Adds two JAX-traceable rescue strategies for difficult DC operating points, alongside API improvements: **Homotopy methods on CircuitLinearSolver:** - `_run_newton`: refactored inner Newton loop returning `(y, converged)` - `solve_dc_checked`: like `solve_dc` but exposes the boolean converged flag - `solve_dc_gmin`: GMIN stepping — steps g_leak log-uniformly from g_start down to nominal, using jax.lax.scan - `solve_dc_source`: source stepping — ramps source amplitudes 10%→100% via jax.lax.scan; requires amplitude_param tagging on sources - `solve_dc_auto`: automatic fallback via jax.lax.cond — tries direct Newton, falls back to GMIN + source stepping if it diverges **Source amplitude tagging:** - `amplitude_param` field added to `@component`/`@source` decorators and stored on generated CircuitComponent subclasses - `VoltageSource`, `SmoothPulse`, `VoltageSourceAC` tagged with `"V"`; `CurrentSource` tagged with `"I"` - `ComponentGroup` stores `amplitude_param`; assembly scales the tagged param leaf via `eqx.tree_at` when `source_scale != 1.0` **API improvements:** - `g_leak` now settable via `from_component_groups` and `analyze_circuit` - All tests migrated from `solver.from_component_groups` to `analyze_circuit` as the single canonical entry point
KLUSplitFactorSolver — "factor once per time step, solve N times per Newton iteration" (Modified Newton / frozen-Jacobian): - Remove duplicate field declarations inherited from KLUSplitSolver - Fix factor_jacobian: self.symbolic_handle → self._handle_wrapper.handle - Remove stale commented-out __del__/cleanup block - Add from_component_groups override with correct return type annotation - Add to backends dict as "klu_split_factor" (falls back to KLUSolver when klujax split interface is unavailable) FactorizedTransientSolver — uncommented and activated: - Assembles and factors the Jacobian once per time step at y_pred - Each Newton iteration uses assemble_residual_only + solve_with_frozen_jacobian - free_numeric called after Newton loop; import guarded for older klujax - test_factorized_transient_matches_vectorized added (skips when split interface unavailable) Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Rename KLUSplitFactorSolver → KLUSplitLinear (backward-compat alias kept). Add KLUSplitQuadratic(KLUSplitLinear) with refactor_jacobian(), which updates the numeric LU factorization in-place via klujax.refactor — reusing the existing memory and fill-reducing permutation for cheaper per-iteration Jacobian updates. Add RefactoringTransientSolver(FactorizedTransientSolver): factors once at the predicted state, then calls refactor_jacobian() + solve each Newton iteration. The returned handle is threaded through the computation graph to prevent XLA dead code elimination. Gives full Newton (quadratic) convergence at lower cost than re-factoring from scratch each iteration. backends["klu_split"] now resolves to KLUSplitQuadratic when klujax.refactor is available, falling back to KLUSplitLinear otherwise. split_refactor_available flag added for conditional availability checks. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Introduce three BDF2 solver classes — BDF2VectorizedTransientSolver, BDF2FactorizedTransientSolver, and BDF2RefactoringTransientSolver — each pairing with its Backward Euler counterpart and adding variable-step 2nd-order accuracy. Key changes: - assembly.py: add `alpha` parameter to assemble_system_real/complex so callers can pass α/h as the charge-derivative scaling factor; default α=1.0 preserves all existing callers (DC, harmonic balance, BE) - transient.py: add _bdf2_preamble() helper that computes the variable- step BDF2 coefficients (ω, α₀, α₁, α₂) and history charges, branches between BE (step 0) and BDF2 (step ≥ 1) via jnp.where, and returns a make_residual() closure; add three BDF2 solver classes with solver_state = (y_nm1, h_nm1, step_count: int32) - __init__.py: export the three new BDF2 solver classes - tests/test_transient.py: add BDF2 smoke tests (float + complex) and BDF2Factorized vs BDF2Vectorized agreement test
Introduce three SDIRK3 solver classes — SDIRK3VectorizedTransientSolver, SDIRK3FactorizedTransientSolver, and SDIRK3RefactoringTransientSolver — using Alexander's 3-stage L-stable tableau (γ ≈ 0.4359). The key efficiency property of SDIRK is that the diagonal element γ is identical for all three stages, so J_eff = dF/dy + (1/(γh))·dQ/dy is the same matrix across all stage solves — factor once, reuse across all three stages and their Newton iterations. Companion residual per stage: R_i(Y_i) = F(Y_i) + (Q(Y_i) - Q_hist_i) / (γh) = 0 with accumulated history charges Q_hist_2 and Q_hist_3 built from the previous stage solutions. solver_state = (y_prev, dt_prev) — identical to Backward Euler, since SDIRK3 is a one-step method requiring no cross-step history.
BDF2 is now the out-of-the-box solver for setup_transient(), replacing Backward Euler. Users who need explicit 1st-order behaviour can still pass transient_solver=VectorizedTransientSolver explicitly.
Two bugs in `_junction_charge` caused a negative capacitance (dQ/dV < 0) that destabilised the Newton solver when junction caps were non-zero: 1. Stray `- 1.0` negated the depletion charge formula, giving dQ/dV = -Cj 2. Linear extrapolation above fc*Vj used q_normal(v) evaluated at runtime v instead of the precomputed q at the threshold, so the slope diverged Fixed formula: Q = Cj0*Vj/(1-m) * [1 - (1 - v/Vj)^(1-m)] for v < fc*Vj Q = q_thresh + C_linear * (v - v_thresh) for v >= fc*Vj where q_thresh and C_linear are constants. dQ/dv = +Cj(v) > 0 everywhere. Also add `benchmark` pixi feature/environment with ngspice for circuit accuracy benchmarking against NGSpice reference simulations.
Three transient benchmarks comparing Circulax vs NGSpice reference: - rc_pulse_vs_ngspice.py: RC + PULSE source, τ=1ms, 1s simulation - diode_cascade_vs_ngspice.py: 20-stage rectifier chain, ~138V peak - fullwave_rect_vs_ngspice.py: full-wave bridge rectifier, D1N4007, 100µF load Shared infrastructure in bench_utils/: - ngspice.py: subprocess runner + output parser - metrics.py: waveform comparison (max/RMS error) - plotting.py: multi-panel comparison plots NGSpice netlists in circuits/. Run via: pixi run -e benchmark python benchmarking/<name>_vs_ngspice.py
- Rename *_vs_ngspice.py → *_testbench.py - Add bench_utils/runner.py: SolverResult dataclass + run_benchmark() callable dict interface — adding a new solver is a one-liner in SOLVERS - Warmup reduced to 2 steps (JAX traces on first call regardless of count) - RC pulse: backend dense; diode cascade: backend klu_split - Add pixi tasks bench_cascade and bench_rectifier
Sweeps N sections (default: 100, 1000, 10000, 50000) of an LC transmission line (L=10nH, C=4pF, Z₀=50Ω) with adaptive PID stepping. - klu_split backend throughout - SaveAt uses fn= projection to save only 2 I/O node voltages at [t0, T_MAX], avoiding O(N) memory per save point at large sizes - Reports compile, warmup, sim time and µs/step per size in a table - Pixi task: bench_lc_ladder
Introduces circulax/solvers/circuit_diffeq.py — a focused replacement for diffrax.diffeqsolve that removes all code paths not used by circuit simulation: - Events, dense output, progress meter removed from State and loop body - SDE/CDE support and backward-compat shims removed - made_jump treated as static False (never in carry state) - Assumes t0 < t1 (forward-only time direction) - saveat.steps saving mode not supported CircuitState shrinks from 16 → 9 fields, reducing pytree tracing overhead. Compatible with PIDController, SaveAt(ts=...), and produces diffrax.Solution output unchanged. setup_transient() updated to call circuit_diffeqsolve instead of diffrax.diffeqsolve; exposes checkpoints= kwarg for AD tuning. All 110 tests pass. Accuracy vs NGSpice unchanged across all benchmarks. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
- benchmarking/circuits/diode_clipper.cir: NGSpice netlist for D1N4007
half-wave clipper (Rs=1kΩ, RL=10kΩ, Vs=2V 1kHz, 20-cycle transient)
- benchmarking/diode_clipper_hb_testbench.py: three-way comparison of
NGSpice transient, Circulax HB, and Circulax transient solvers
- single-frequency mode: timing + accuracy + harmonic spectrum table
- --sweep mode: 100-frequency vmap sweep vs NGSpice serial (11× timed
speedup, ~1 mV RMS accuracy); uses eqx.tree_at to override the
VoltageSourceAC freq param per frequency so HB converges correctly
Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
…ncy sweep guide ReadMe.md / docs/index.md: - Fix import paths (circulax.base_component → circulax.components.base_component, circulax.components → circulax.components.electronic) - Remove removed diffrax.TqdmProgressMeter kwarg from example - Add Harmonic Balance and AC Sweep bullets to Analysis section docs/SUMMARY.md: - Expand Examples to explicit list including new ac_sweep page docs/examples/ac_sweep.md + ac_sweep_files/: - Execute ac_sweep.ipynb via papermill and convert to Markdown docs/harmonic_balance.md: - Add "Frequency Sweep" section documenting the jax.vmap + eqx.tree_at pattern for sweeping HB over 100 frequencies in one XLA call
Covers Dense, Sparse, KLU, and KLU Split (marked experimental) with a decision table and guidance on when to switch backends. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Remove incorrect claim that KLU backends don't support vmap/grad. All four backends are tested with jax.vmap and jax.grad in test_harmonic_balance.py. Drop the vmap/grad column from the summary table and remove related incorrect statements from KLU and KLU Split sections. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
…ontrol Documents BDF2/SDIRK3 solver order, A-stability, and adaptive step sizing via diffrax.PIDController with tuning guidance. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Adds a linux-64 activation script that sources .env from the project root (if present) using set -a / set +a so all variables are exported. No-op on macOS or when .env doesn't exist. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Move activation script from linux-64 target to workspace-level so it also runs on osx-arm64. Both platforms are POSIX so the same sh script works. Windows is not a supported platform.
Restructure the guide to lead with KLU as the default and industry-standard choice (used by SPICE, Spectre, HSPICE), demote Dense to a niche option for very small circuits or GPU/HB use, and update the summary table accordingly. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
Re-execute and convert all example notebooks; update pixi.lock to match current dependency resolution. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
circuit_diffeqsolve: stripped-down diffrax integration loop removing unused features (events, dense output, progress meter) for cleaner circuit simulation_junction_chargesign error and linear extrapolation anchor in BJT/diode models.envloading on pixi environment activation; pixi task cleanup