feat: add prepared stationary frequency response

[oameye]

Supersedes the architectural prototype in #21 with a stationary response layer rebuilt on current main.

What this adds

• reusable prepared FrequencyResponse objects;
• dense shifted-Hessenberg response backend behind an abstract solver boundary;
• generic Kubo susceptibility;
• complete multiport normal/anomalous (Bogoliubov/Nambu) scattering response;
• single S-parameter convenience access;
• connected normally ordered output-emission spectra;
• vacuum-normalized quadrature/squeezing spectra;
• standard SLH output convention b_out = S b_in + L.

The expensive symbolic translation, Liouvillian construction, stationary-state validation and Hessenberg preparation are performed once and reused. Port response RHS vectors and trace-contraction weights are cached directly in Liouville space, so prepared frequency sweeps avoid reconstructing temporary QuantumOptics operators and avoid per-frequency matrix products for scalar observables.

For nontrivial static scattering matrices, S contributes both to the direct output and to the internal weak coherent drive through (S†L)_j. Symbolic scattering entries are lowered through the existing SecondQuantizedAlgebra/Symbolics interfaces, including canonical SecondQuantizedAlgebra.expim(ϕ) phases and ordinary symbolic sin/cos beam-splitter matrices.

Physical conventions

The output field follows the package SLH convention

b_out = S b_in + L.

Emission spectra are connected and normally ordered, with vacuum floor zero. The coherent stationary mean is removed from the smooth spectrum because its elastic contribution is a delta peak. Quadrature spectra are connected homodyne spectra with shot-noise floor one.

Parametrically driven systems expose both normal and anomalous scattering blocks rather than forcing a scalar-S description onto a Bogoliubov response.

Solver assumptions

DenseHessenberg targets small-to-medium Hilbert spaces. It explicitly assumes a unique stationary state and checks normalization and stationarity by default. Rank-one deflation removes the stationary null direction; it is not intended to regularize additional zero or undamped modes. The public response API is solver-agnostic so sparse/Krylov backends can be added later without changing the physical interface.

Prepared states and numeric susceptibility observables are also checked for compatibility with the translated QuantumOptics basis before the implementation drops to raw Liouville-space matrices.

Validation

The test suite includes:

• analytic passive two-port cavity S11/S21 and lossless power conservation;
• analytic Kubo susceptibility;
• one-port resonant reflection;
• symbolic phase and symbolic real rotation scattering matrices;
• observable internal loss;
• DPA normal and anomalous scattering against analytic formulas;
• Bogoliubov conservation;
• analytic squeezed and anti-squeezed quadrature spectra;
• analytic DPA emission;
• nonlinear Kerr-parametric-oscillator response, susceptibility and emission;
• prepared/one-shot API equivalence and public API contract/error paths;
• rejection of same-dimension but incompatible QuantumOptics bases for prepared states and numeric susceptibility observables.

Compiler acceptance includes package-wide JET plus explicit JET.@test_opt gates for scalar and vector-frequency prepared scattering, susceptibility, emission and quadrature APIs. These gates pass on the implementation. Package tests are validated on Julia 1.10 LTS and Julia 1.12.

Performance

The benchmark suite reproduces #21's prepared operator-reconstruction hot loop and checks numerical equality to < 1e-10 before timing it. On the final exact head 6fc7ba639314f1d6b297b1ff717638819bf8d3bf, a 201-point prepared S11 sweep measured:

implementation	time	memory	allocations
#21 reproduced prepared hot loop	20.347 ms	3.385 MB	4,425
#66 prepared Liouville-vector path	20.181 ms	2.221 MB	2,817

The final exact-head prepared path is therefore about 0.8% faster, uses about 34% less memory, and performs about 36% fewer allocations than the reproduced #21 hot loop on the same runner and resolved environment. A previous same-run measurement on a faster hosted runner gave 13.482 ms versus 13.663 ms with the same memory/allocation counts, so the relative result is reproducible across runner speed differences.

The final exact-head benchmark measurement step completed successfully. The overall Benchmark Tracking workflow is red because the repository-wide historical comparison flags numerous unrelated existing benchmark families (pulse coupling, translation, SLH algebra, and interaction-picture workloads) as slower than the stored main baseline on this hosted runner. #66 does not modify those code paths. The response-specific paired old/new benchmark is therefore the relevant no-regression acceptance check; benchmark alert/failure thresholds have not been weakened.

The benchmark workflow reports the BenchmarkTools time, memory and allocation counters explicitly.

[codecov]

Codecov Report

❌ Patch coverage is 95.91837% with 10 lines in your changes missing coverage. Please review.
✅ Project coverage is 97.00%. Comparing base (6f71914) to head (6fc7ba6).

Files with missing lines	Patch %	Lines
src/response.jl	95.91%	10 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main      #66      +/-   ##
==========================================
- Coverage   97.58%   97.00%   -0.59%     
==========================================
  Files           4        5       +1     
  Lines         455      700     +245     
==========================================
+ Hits          444      679     +235     
- Misses         11       21      +10     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[github-actions]

Benchmark Results

Details

Benchmark suite	Current: 6fc7ba6	Previous: 6f71914	Ratio
Correlations/two-time/single photon cavity	162215353 ns	150701969 ns	1.08
Interaction Picture/coefficient matrix M/analytical (2 equal modes)	11415 ns	10333 ns	1.10
Interaction Picture/coefficient matrix M/numerical (ODE)	53462 ns	46293 ns	1.15
Interaction Picture/coupling matrix evaluation/2 modes	83.2655601659751 ns	77.48098663926002 ns	1.07
Interaction Picture/coupling matrix evaluation/4 modes	162.96 ns	135.16782407407408 ns	1.21
Interaction Picture/operator substitution/TLS cascade	18345 ns	15344 ns	1.20
Pulse Couplings/multi-pulse/input 2 modes	310431 ns	265843 ns	1.17
Pulse Couplings/multi-pulse/output 2 modes	370430 ns	322080 ns	1.15
Pulse Couplings/single-pulse/input	241506 ns	210935 ns	1.14
Pulse Couplings/single-pulse/input Gaussian	63.710204081632654 ns	54.118781725888326 ns	1.18
Pulse Couplings/single-pulse/output	271193 ns	232318 ns	1.17
Pulse Couplings/single-pulse/output Gaussian	64.62346938775511 ns	55.17157360406091 ns	1.17
Response/S11 single prepared	77403 ns
Response/S11 sweep legacy operator	16084400 ns
Response/S11 sweep one-shot	18340409 ns
Response/S11 sweep prepared	15844151 ns
Response/emission sweep prepared	15856519 ns
Response/full Bogoliubov sweep prepared	22726502 ns
Response/prepare	1884248 ns
Response/quadrature sweep prepared	31732270 ns
SLH Algebra/closure evaluation/2-QD waveguide H(t)	28059 ns	26974 ns	1.04
SLH Algebra/closure evaluation/2-QD waveguide L(t)	2819.6666666666665 ns	3000.1111111111113 ns	0.94
SLH Algebra/numeric/2-QD waveguide composition	25163 ns	22400 ns	1.12
SLH Algebra/symbolic/3-cavity cascade	6269.4 ns	5143.714285714285 ns	1.22
SLH Algebra/symbolic/concatenate + cascade	2126.1 ns	1863.5 ns	1.14
SLH Algebra/symbolic/feedback OPO loop	5468.857142857143 ns	4480.875 ns	1.22
Translation/closure evaluation/3-cavity H(t)	8330.2 ns	7275.666666666667 ns	1.14
Translation/closure evaluation/3-cavity L(t)	5540.2 ns	4986.666666666667 ns	1.11
Translation/static/atom-cavity	17115 ns	14047 ns	1.22
Translation/time-dependent/3-cavity H+L	383703 ns	314809 ns	1.22
Translation/time-dependent/atom-cavity	52968 ns	43828 ns	1.21

This comment was automatically generated by workflow using github-action-benchmark.

[github-actions]

⚠️ Performance Alert ⚠️

Possible performance regression was detected for benchmark 'Benchmark Results'.
Benchmark result of this commit is worse than the previous benchmark result exceeding threshold 1.30.

Benchmark suite	Current: 6fc7ba6	Previous: 6f71914	Ratio
Interaction Picture/coefficient matrix M/analytical (2 equal modes)	14342 ns	10333 ns	1.39
Interaction Picture/coefficient matrix M/numerical (ODE)	67111 ns	46293 ns	1.45
Interaction Picture/coupling matrix evaluation/2 modes	102.44957081545064 ns	77.48098663926002 ns	1.32
Interaction Picture/coupling matrix evaluation/4 modes	217.9388888888889 ns	135.16782407407408 ns	1.61
Interaction Picture/operator substitution/TLS cascade	25418 ns	15344 ns	1.66
Pulse Couplings/multi-pulse/input 2 modes	442305 ns	265843 ns	1.66
Pulse Couplings/multi-pulse/output 2 modes	499371 ns	322080 ns	1.55
Pulse Couplings/single-pulse/input	358930 ns	210935 ns	1.70
Pulse Couplings/single-pulse/input Gaussian	94.52410901467505 ns	54.118781725888326 ns	1.75
Pulse Couplings/single-pulse/output	396316 ns	232318 ns	1.71
Pulse Couplings/single-pulse/output Gaussian	94.27148846960168 ns	55.17157360406091 ns	1.71
SLH Algebra/closure evaluation/2-QD waveguide H(t)	36846 ns	26974 ns	1.37
SLH Algebra/numeric/2-QD waveguide composition	30746 ns	22400 ns	1.37
SLH Algebra/symbolic/3-cavity cascade	8419.333333333334 ns	5143.714285714285 ns	1.64
SLH Algebra/symbolic/concatenate + cascade	2560.5555555555557 ns	1863.5 ns	1.37
SLH Algebra/symbolic/feedback OPO loop	6261.4 ns	4480.875 ns	1.40
Translation/closure evaluation/3-cavity L(t)	9684 ns	4986.666666666667 ns	1.94
Translation/static/atom-cavity	22934 ns	14047 ns	1.63
Translation/time-dependent/3-cavity H+L	514885 ns	314809 ns	1.64
Translation/time-dependent/atom-cavity	81293 ns	43828 ns	1.85

This comment was automatically generated by workflow using github-action-benchmark.