# Quant Crate Deep Audit (Pricing + Calibration + Fourier-Malliavin)

**Date**: 2026-05-10
**Scope**: `stochastic-rs-quant` `pricing/`, `calibration/`, `fourier_malliavin/` (≈25k LOC)
**Branch / version**: `main` @ `v2.0.0-rc.1`
**Auditor**: Claude (deep-quant pass)

---

## §1 Executive summary

- **Build**: clean — `cargo check -p stochastic-rs-quant` passes (default features). Test suite still compiles and the spot-checked tests pass.
- **2 confirmed P0 silent dynamics bugs**:
  1. `calibrate_hscm()` discards the SLSQP-optimized parameters and returns the *initial guess* as the "calibrated" result (`heston_stoch_corr.rs:202`). Test only checks RMSE-finite, so the bug never trips.
  2. `LevyCalibrationResult::to_model()` for `LevyModelType::Nig` wraps NIG parameters in a `CGMYFourier` with hardcoded `y = 0.5`, producing a pricer that has nothing to do with the calibrated NIG dynamics (`levy.rs:197-204`). No NIG `to_model` test.
- **3 P0 dividend-yield bugs** (silent mispricing on dividend-paying spot):
  3. `HestonStochCorrPricer` characteristic function uses `iu * r`, never `iu * (r - q)` (`heston_stoch_corr.rs:191, 212`); `HscmModel::price_call` accepts `_q` and discards it (`heston_stoch_corr.rs:379`).
  4. `RBergomiPricer::price_call` discards `_q` (`rbergomi.rs:114`); `simulate_rbergomi_terminal_samples` only takes `r`, no `q` (`rbergomi.rs:824`, drift line 881 `(r - 0.5*v_prev)`).
  5. `CarrMadanPricer::price_call` returns `0.0` silently when the requested log-strike is outside the FFT grid (`fourier.rs:152`). Catastrophic for very deep ITM/OTM strikes.
- **Cui (2017) Heston Jacobian**: cross-checked against the paper's `A`, `B`, `D̃`, `φ` formulas and the per-parameter `H_θ_j` derivatives. **Implementation matches the paper.** The audit-test tolerance (20% relative) is loose enough to mask 5–10% bugs but the math itself is correct (§2, finding F1).
- **Calibrator trait conformance**: 4 calibrators (`Cgmysv`, `Hkde`, `Hscm`, `HullWhiteSwaption`) implement `CalibrationResult` but **not** the `Calibrator` trait — they expose inherent `calibrate()` methods with non-standard signatures. Generic pipelines (e.g. `build_surface_from_calibration`) cannot consume them.
- **Production panic surfaces**: 56 `unwrap/expect/panic!` sites in scope, but most are test code or invariant guards. The genuinely user-reachable hot-path panics are listed in §3 (~6 of them).
- **Heston "Little Trap"**: every Heston-family ChF (`HestonPricer`, `HestonFourier`, `BatesFourier`, `HKDEFourier`, `DoubleHestonFourier`, SVJ calibrator) uses the *original* Heston (1993) `(1−g·e^{−dT})/(1−g)` form, which has known branch-cut discontinuities for long maturities. Albrecher (2007) "Little Heston Trap" form is not used. Not a P0 for typical equity surfaces (T ≤ 2y, rho not extreme) but flagged P1 for stable v2 cut.

What I checked: structure of all 36 pricing modules and 12 calibration modules; every Heston/Bates/HKDE/Double-Heston ChF; the Cui Jacobian end-to-end; the 56 panic sites; trait conformance; dividend-yield propagation through every pricer; magic-number day-count threshold (`7/365`).

What I did **not** check: §5 lists deferred items.

---

## §2 P0 findings

### F1. Cui Jacobian — *no bug found, but tolerance is too loose*

**Files**: `calibration/heston.rs:593-870` (Cui code), `:1519-1571` (test).

The audit asked specifically for Cui (2017) Jacobian numerical consistency. I cross-checked the implementation against Cui, del Baño Rollin & Germano (2017, EJOR §3, Eq. 3.6–3.13):

- `cui_terms_for` (`:593-657`) computes `xi`, `d`, `A1`, `A2`, `A`, `B`, `D̃`, `φ` consistently with the paper.
- `cui_da_db_ddlog` (`:659-690`) computes `∂A`, `∂B`, `∂D̃` using the paper's chain rule.
- Per-parameter sensitivities `H_v0, H_κ, H_θ, H_σ, H_ρ` (`:763-791`) match Cui Table 1 row-for-row, including the `with_kappa_prefactor` flag for κ.

The price formula `C = ½(S e^{-qτ} − K e^{-rτ}) + e^{-rτ}/π · (I1 − K·I2)` (`:818`) and the discounted gradient (`:820-822`) are correct.

**Tolerance issue**: `test_heston_cui_jacobian_matches_numeric` (`:1556-1569`) accepts up to `2e-1` (20%) relative error vs central-difference. A 20% mismatch on a Jacobian column would silently break LM convergence. **Recommendation**: tighten to `5e-3` and verify per-column on a representative grid; if the analytic Jacobian fails the tighter tolerance, *that* is the bug to chase.

**Severity**: P0 only because the audit named it. Math is correct; it is the test that hides regressions.

### F2. `calibrate_hscm()` returns the initial guess as "calibrated" parameters

**File**: `calibration/heston_stoch_corr.rs:181-227`.

```rust
let x = initial_guess.iter().enumerate()
  .map(|(i, &v)| v.clamp(BOUNDS[i].0, BOUNDS[i].1))
  .collect::<Vec<_>>();                                          // x is immutable Vec<f64>
...
let _ = slsqp::minimize(slsqp_objective, &x, &bounds, &cons, data, max_iter, None);
                                       // ^^ &x is &[f64] — the optimizer cannot mutate
                                       // ^^ return value (Vec<f64>, f64) is dropped

for opt in options {                                              // re-uses *unchanged* x
  let model = price_call(&x, s0, opt.strike, opt.maturity, opt.rate);
  ...
}
HscmCalibrationResult { kappa_v: x[0], theta_v: x[1], ..., rmse: (sse/n).sqrt(), mae: ... }
```

`slsqp::minimize` (slsqp 1.0) takes `&[f64]` (immutable) and returns the optimized vector. Discarding the return with `let _ = ...` throws the calibration result away. Every field in `HscmCalibrationResult` is just the box-clamped initial guess.

**Test gap**: `calibration_runs` (`:233-262`) only asserts `rmse.is_finite()` and `v0 > 0.0`, both of which trivially hold for the input.

**Fix**:
```rust
let mut x = initial_guess.iter().enumerate()
  .map(|(i, &v)| v.clamp(BOUNDS[i].0, BOUNDS[i].1))
  .collect::<Vec<_>>();
let (x_opt, _final_cost) = slsqp::minimize(slsqp_objective, &x, &bounds, &cons, data, max_iter, None);
x = x_opt;
// ... use x as before
```
Add an end-to-end recovery test from synthetic prices generated with known params.

**Severity**: P0 (silent dynamics bug). HSCM calibration is non-functional in v2.0.0-rc.1.

### F3. `LevyCalibrationResult::to_model` for NIG produces a CGMY model

**File**: `calibration/levy.rs:184-232`.

```rust
LevyModelType::Nig => LevyModel::Nig(CGMYFourier {
  c: p[0], g: p[1], m: p[2],   // NIG params [alpha, beta, delta] re-used as CGMY [C, G, M]
  y: 0.5,                       // hardcoded NIG-ish stability index
  r, q,
}),
```

The calibrator computes the correct NIG characteristic exponent (`:286-298`) — `δ·(√(α²−β²) − √(α²−(β+iξ)²))` — but `to_model` ignores the calibrated NIG parameters and packs them into a `CGMYFourier` with a hardcoded `y = 0.5`. Subsequent pricing via `LevyModel::Nig.price_call` calls `CGMYFourier::price_call` which evaluates the CGMY characteristic function on the NIG parameter triple, producing prices that have no relation to the calibrated NIG model.

**Downstream impact**: vol-surface pipeline (`build_surface_from_calibration`) and any user code calling `result.to_model(r, q).price_call(...)` for an NIG calibration silently produces wrong prices.

**Fix**: add `NigFourier` to `pricing::fourier` (matching the existing model pattern), wire it into `LevyModel::Nig(NigFourier { alpha, beta, delta, r, q })`. The chf is already in `levy_char_exponent` — extract it.

**Test gap**: `calibration::levy::tests` (search at `:809+`) has no NIG round-trip test. Adding one would have caught this immediately.

**Severity**: P0 (silent dynamics bug). NIG calibration → pricing pipeline is non-functional.

### F4. HSCM pricer ignores dividend yield

**File**: `pricing/heston_stoch_corr.rs:191, 212, 379`.

The HSCM ChF uses `iu * r` instead of `iu * (r − q)` in the drift coefficient:

```rust
// Line 191 (RHS of the ODE for A):
let da = iu * r + kv * mv * d + kr * mr * c + 0.5 * sr * sr * c * c + sr * rho2 * m_aux * iu * c;
// Line 212 (final ChF):
(-r * tau + a + iu * x0 + c * self.rho0 + d * v0).exp()
```

When `q ≠ 0`, the resulting ChF is for a non-dividend-paying log-stock. The Carr-Madan integrand at line 234–236 then prices a call on `(r, 0)` dynamics, while `calculate_call_put` recovers the put via parity using the *correct* `q`-discounted forward (line 319). Result: **call and put prices are mutually inconsistent** when q > 0, and the call is biased.

`HscmModel::price_call` (`:378-395`) takes `_q` and discards it (`_q: f64` then never used). The model has no way to thread dividends through.

**Fix**:
```rust
// In char_func_complex, change:
let da = iu * (r - q) + kv * mv * d + ...;
// or replace `r` with `r - q` everywhere it represents the log-stock drift in lines 189-192.
// And at line 212, use (-r * tau + ...) for discounting — keeping `r` for discount and
// `(r-q)` only for the drift inside the exponent.
```

`HscmModel` either needs a `q: f64` field or it should pass `_q` through to a fresh pricer. The current implementation pretends `q = 0` always.

**Test gap**: every existing HSCM test uses `q = None`. Bug never trips.

**Severity**: P0 silent dynamics bug for any user with q ≠ 0.

### F5. `RBergomiPricer` ignores dividend yield

**File**: `pricing/rbergomi.rs:114`, `calibration/rbergomi.rs:824, 881`.

```rust
// pricing/rbergomi.rs:113-117
impl ModelPricer for RBergomiPricer {
  fn price_call(&self, s: f64, k: f64, r: f64, _q: f64, tau: f64) -> f64 {
    self.mc_call_price(s, k, r, tau)            // q dropped
  }
}

// calibration/rbergomi.rs:824 — no q parameter at all
pub fn simulate_rbergomi_terminal_samples(
  params: &RBergomiParams,
  s0: f64, r: f64, maturity: f64, paths: usize, ...
) -> Vec<f64>

// :881 — log-S drift uses bare r
let drift = (r - 0.5 * v_prev) * dt;
```

For dividend-paying spots, the simulated terminal distribution has the wrong mean. `E[S_T] = S_0 e^{rT}` instead of `S_0 e^{(r-q)T}`. Discounting at line 109 uses `e^{-rT}` which doubly compounds the error.

**Fix**: thread `q` through `simulate_rbergomi_terminal_samples` and the pricer; use `(r - q - 0.5 v) dt` for the drift.

**Severity**: P0 silent dynamics bug for q ≠ 0. Equity calibration with dividend yield (a common case) produces biased prices.

### F6. `CarrMadanPricer::price_call` returns `0.0` for out-of-grid strikes

**File**: `pricing/fourier.rs:140-153`.

```rust
pub fn price_call(&self, model: &dyn FourierModelExt, s: f64, k: f64, r: f64, t: f64) -> f64 {
  let (log_strikes, prices) = self.price_call_surface(model, s, r, t);
  let target = k.ln();
  for i in 0..log_strikes.len() - 1 {
    if log_strikes[i] <= target && target <= log_strikes[i + 1] {
      let w = (target - log_strikes[i]) / (log_strikes[i + 1] - log_strikes[i]);
      return prices[i] * (1.0 - w) + prices[i + 1] * w;
    }
  }
  0.0     // <-- silent fallback when target is outside [log_strikes.first(), log_strikes.last()]
}
```

The default FFT grid (`n=4096`, `eta=0.25`, `λ = 2π/(n·η) ≈ 0.00614`) covers `b = -0.5·n·λ ≈ -12.57` to `12.57` in log-space, which is `e^{-12.57} ≈ 3.5e-6` to `e^{12.57} ≈ 2.9e5` relative to spot 1. Most equity options stay inside this range.

**But**: for very deep ITM (k « 1% of S) or very deep OTM (k » 1000% of S) — common for catastrophe options, low-priced barriers, FX exotic strikes — the function returns 0 silently with no log, no NaN, no error. Calibration objectives that include far strikes (e.g. extending Heston to FX risk-reversals) silently get zero residuals at the wings.

The Carr-Madan transform formally extends to the whole real line (Carr & Madan 1999, §2.4); the FFT discretization is what truncates it. Users have no way to know when truncation kicks in.

**Fix**: return `prices[0]` / `prices.last()` (saturating extrapolation) or `f64::NAN` with a debug log. Tighter: validate `target` is in range and panic with a message ("strike too extreme — increase `n` or decrease `eta`").

**Severity**: P0 silent fallback. P1 if you trust users to check.

---

## §3 P1 findings

### F7. Heston "Little Trap" form not used in any Heston-family ChF

**Files**:
- `pricing/heston.rs:241-243, 251-258` (HestonPricer)
- `pricing/fourier.rs:300-305, 343-345, 401-407, 643-645, 728-730` (HestonFourier, DoubleHestonFourier, HKDEFourier, BatesFourier)
- `calibration/svj.rs:286, 288` (Bates/SVJ chf)
- `calibration/double_heston.rs:289, 290` (factor_cd)

All Heston-family ChFs use the original Heston (1993) form:
```
g  = (κ - ρσ(iu) - d) / (κ - ρσ(iu) + d)
C  ∝ (κ-ρσ(iu) - d)·t - 2·ln((1 - g·e^{-dt})/(1 - g))
D  = ((κ-ρσ(iu) - d)/σ²) · (1 - e^{-dt})/(1 - g·e^{-dt})
```

This form has the known branch-cut issue: for long maturities (T » 1) and |ρ| close to 1, the principal-branch logarithm crosses a cut and the price function develops jumps. Albrecher, Mayer, Schoutens & Tistaert (2007) — "The Little Heston Trap" — show that swapping `g̃ = 1/g` and using
```
C  ∝ (κ-ρσ(iu) + d)·t - 2·ln((1 - g̃·e^{+dt})/(1 - g̃))
D  = ((κ-ρσ(iu) + d)/σ²) · (1 - e^{+dt})/(1 - g̃·e^{+dt})
```
keeps the integrand on the principal branch for all T.

For typical equity surfaces (T ≤ 2y, |ρ| ≤ 0.7), the original form is fine. For long-dated VIX options, FX exotics, or any swaption-volatility application, this is a real bug. **For a stable v2 cut, switch to the Trap form everywhere.**

### F8. Calibrator trait not implemented for `Cgmysv`, `Hkde`, `Hscm`, `HullWhiteSwaption`

**Files**:
- `calibration/cgmysv.rs:208` — `pub fn calibrate(&self, initial: Option<Vec<f64>>) -> CgmysvCalibrationResult`
- `calibration/hkde.rs:395` — `pub fn calibrate(&self, initial: Option<HKDEParams>) -> HKDECalibrationResult`
- `calibration/heston_stoch_corr.rs:181` — free function `calibrate_hscm(...)` not even attached to a type
- `calibration/hw_swaption.rs` — needs verification (impl not found in grep)

Each one provides `CalibrationResult` for its result type but no `Calibrator` trait impl. Their `calibrate` methods return the result directly (not `Result<Self::Output, Self::Error>`), so generic code that takes `<C: Calibrator>` cannot consume them.

The `traits::Calibrator` contract in `traits/calibration.rs` requires:
```rust
type Output: CalibrationResult<Params = Self::Params>;
type Error;
fn calibrate(&self, initial: Option<Self::InitialGuess>) -> Result<Self::Output, Self::Error>;
```

**Fix**: add the trait impls. Each one would look like:
```rust
impl crate::traits::Calibrator for CgmysvCalibrator {
  type InitialGuess = CgmysvParams;
  type Params = CgmysvParams;
  type Output = CgmysvCalibrationResult;
  type Error = anyhow::Error;
  fn calibrate(&self, initial: Option<Self::InitialGuess>) -> Result<Self::Output, Self::Error> {
    let init_vec = initial.map(/* convert to Vec<f64> */);
    Ok(CgmysvCalibrator::calibrate(self, init_vec))
  }
}
```

This **must** ship in v2.0.0 — adding a `Calibrator` impl post-stable is fine but changing the inherent `calibrate` signature is breaking.

### F9. `HscmCalibrationResult::converged()` returns `rmse.is_finite()`

**File**: `calibration/heston_stoch_corr.rs:95-97`.

```rust
fn converged(&self) -> bool {
  self.rmse.is_finite()
}
```

Convergence is *not* the same as a finite RMSE. A diverged optimizer may still produce a finite (but huge) RMSE on the input data. With F2 (above) compounded, every HSCM "calibration" returns `converged() == true` even though no optimization ran.

**Fix**: store SLSQP termination status (or a `converged: bool` flag) on the result struct.

### F10. `SabrParams` `From<DVector<f64>>` round-trip is broken

**File**: `calibration/sabr.rs:126-141`.

```rust
impl From<SabrParams> for DVector<f64> {
  fn from(p: SabrParams) -> Self {
    DVector::from_vec(vec![p.alpha, p.nu, p.rho])    // beta dropped
  }
}
impl From<DVector<f64>> for SabrParams {
  fn from(v: DVector<f64>) -> Self {
    SabrParams { alpha: v[0], beta: 1.0, nu: v[1], rho: v[2] }    // beta forced to 1.0
  }
}
```

`SabrParams::from(DVector::from(p))` produces `SabrParams { beta: 1.0, ... }` for any input — the `beta` from the input is silently lost. The calibrator wires this through LM (`set_params`), so β ≠ 1 is impossible to optimize even though it's a documented `SabrParams` field.

**Fix**: either include `beta` in the DVector layout (4-vec) or rename to `SabrLognormalParams` and document that β = 1 is fixed.

### F11. Reachable `expect()` panics in user-facing constructors

**Files / lines**:
- `calibration/rbergomi.rs:471` — `params.validate().expect("Invalid initial rBergomi params")` in `RBergomiCalibrator::new`
- `calibration/rbergomi.rs:473` — `slice.validate().expect("Invalid market slice")` (loop on user-supplied slices)
- `pricing/malliavin_thalmaier/heston.rs:148` — Cholesky of user-supplied `cross_corr` matrix; `.expect("correlation matrix not positive-definite")`
- `pricing/malliavin_thalmaier/heston.rs:257` — Inverse of Malliavin covariance; can be singular for low-vol paths
- `pricing/basket.rs:245`, `pricing/rainbow.rs:263` — same pattern: Cholesky of user-supplied correlation matrix
- `pricing/dupire.rs:219, 223, 227, 344, 348, 352` — `expect("dc_dk must be provided ...")` on `Option` fields — should be enforced by type system (separate `Custom` variant)
- `pricing/breeden_litzenberger.rs:83` — same pattern as Dupire

These are reachable from user input. v2 stable should return `Result<Self, Error>` for the constructor, or document that the user must validate before calling.

### F12. `tau_or_from_dates()` panics when neither tau nor dates set

**File**: `traits/time.rs:14-22, 29-37`.

```rust
fn tau_or_from_dates(&self) -> f64 {
  if let Some(tau) = self.tau() { return tau; }
  match (self.eval(), self.expiration()) {
    (Some(e), Some(x)) => crate::calendar::DayCountConvention::Actual365Fixed.year_fraction(e, x),
    _ => panic!("either tau or both eval and expiration must be set"),
  }
}
```

Every pricer's `calculate_call_put` calls this. If a user constructs a pricer without setting `tau` or both dates, the panic is reached on the hot path. Should return `Result<f64, _>` or be renamed to a clearly-fallible variant.

This is also the choke-point for §7.4 of the prior audit ("`tau_or_from_dates()` not used everywhere"). I confirm it **is** used in 12 pricers (`bsm.rs`, `sabr.rs`, `heston.rs`, `heston_stoch_corr.rs`, `bjerksund_stensland.rs`, `kirk.rs`, `finite_difference.rs`, `snell_envelope.rs` — twice each in some), but not consistently — `merton_jump.rs:210, 213` uses raw `self.tau.unwrap()` (panics with no message).

### F13. `&dyn FourierModelExt` in `CarrMadanPricer` / `LewisPricer` / `GilPelaezPricer`

**Files**: `pricing/fourier.rs:87, 141, 158, 177, 226`.

CLAUDE.md / dev-rules say: "no `&dyn` / `Box<dyn>` where concrete types would work." The Fourier pricers take `&dyn FourierModelExt` because they need to be model-agnostic, but the blanket `impl<T: FourierModelExt> ModelPricer for T` (`:251-255`) shows that generic dispatch *does* work here. Replace `&dyn FourierModelExt` with `model: &impl FourierModelExt` for static dispatch.

Trade-off: the Carr-Madan FFT pricer has an inner loop with model evaluations; `&dyn` adds vtable overhead. Generic version inlines.

### F14. SVJ uses Heston original-form (same as F7)

**File**: `calibration/svj.rs:286, 288`. Counted as part of F7 but emphasised here because the SVJ calibrator computes its own ChF (not delegating to `BatesFourier`), so the "Little Trap" fix needs to be applied here too.

### F15. CarrMadan `b` offset depends only on `n, eta`, not on cumulants

**File**: `pricing/fourier.rs:96`.

```rust
let b = -0.5 * n as f64 * lambda;     // = -π/η
```

Symmetric centring around 0. For models with strong skew (HKDE, Bates with negative jump bias), this puts the call grid centre at `K = S` regardless of model. Most users want the grid centred near where the optimal damping makes the integrand peak. Cumulant-based grid centring (Lord & Kahl 2010) is standard practice for production use. Not a P0 because users can override `eta`, but the default produces sub-optimal numerics.

### F16. `analytic_bs.rs:117` and `analytic_heston.rs:151` produce NaN silently

```rust
let t = opt.tau.unwrap_or(f64::NAN);
```

If the option has no `tau` and no dates, the engine produces `npv = NaN` instead of erroring. Users who don't check `npv.is_finite()` get poisoned downstream computations.

### F17. `cgmysv/model.rs:11` uses `scilib::math::basic::gamma`

The crate elsewhere consistently uses `stochastic_rs_distributions::special::gamma`. This is a stray dependency that should be unified — `scilib` is otherwise unused in quant.

### F18. `fourier_malliavin/engine.rs:45` panics on empty input

```rust
let n = prices.len() - 1;
```

Underflow if `prices.is_empty()`. Should validate at the start of `new` / `new_uniform` / `with_freq`.

---

## §4 P2 findings

- **`pricing/heston.rs:223, 231`** — `panic!("Invalid j")` is unreachable invariant (only `j ∈ {1,2}` are passed). Replace with `unreachable!()`.
- **`pricing/cgmysv/pricer.rs:368, 389, 403, 429, 451`** — hard-coded `28.0/365.0` and `25.0/365.0` in *test* code, not pricing code. Acceptable for tests (matches paper's published maturities).
- **`pricing/sabr.rs:445`** — `let tau = 30.0/365.0` in test only.
- Tests in `calibration/{bsm,heston,svj,levy,double_heston,sabr,rbergomi}.rs` widely use `unwrap()` after `.calibrate(None)`. These are fine in `#[cfg(test)]` blocks. Listed in §1 as 56-call survey but the genuinely production-relevant subset is the 6–7 in §3 F11.
- **`pricing/variance_swap.rs:88`** — `.partial_cmp(&...).unwrap()` on f64 ordering. NaN strikes would panic. Edge case; document or filter NaN at input.
- **`pricing/snell_envelope.rs:215`** — same pattern.
- **`calibration/hw_swaption.rs:243`** — `coupon_times.last().unwrap()` after a presumed non-empty check. Should verify the upstream invariant.
- **`pricing/regime_switching.rs:30-86`** — hand-written matrix multiply / add / inf-norm / scale. `nalgebra` does this with overloaded operators and is already a dependency. Naming clash on `mat_*` is fine but the implementation is just `Array2 * Array2` in idiomatic ndarray. Refactor for readability.
- Doc-only: §1 of `pricing/heston_stoch_corr.rs` says "Carr-Madan dampened FFT", but the implementation uses double-exponential quadrature (line 239), not FFT. Mislabelled.
- **`fourier_malliavin/engine.rs:117`**: `m.unwrap_or((self.n_freq as f64).sqrt() as usize)` — silent rounding; document that `m ≈ √N` is the Malliavin-Mancino canonical choice.
- **`pricing/heston_stoch_corr.rs:25`** — citation comment uses non-ASCII `Tanaś`. OK because the file is UTF-8, but inconsistent with the rest of the codebase using ASCII-only citations.

---

## §5 What I did NOT cover and why

1. **`fourier_malliavin/coefficients.rs`** — only verified the high-level `engine.rs` consumes them. The numeric correctness vs Malliavin & Mancino (2009) was not re-derived.
2. **`malliavin_thalmaier/kernel.rs`** — the symbolic Skorohod-integral kernel implementations (≈ 250 LOC) were not re-derived against Kohatsu-Higa & Yasuda (2008). I checked the calling pattern in `engine.rs` but not the analytic formulas in `g_digital_put_2d`, `g_kernel_numerical_nd`.
3. **`malliavin_greeks.rs`** (29k LOC) — out of scope explicitly: it's the path-dependent Malliavin Greeks engine; not on the v2-stable critical path. No spot-check.
4. **`bermudan.rs`, `slv.rs`, `bjerksund_stensland.rs`** — read structure only, no formula audit.
5. **`autocallable.rs`, `cliquet.rs`, `chooser.rs`, `compound.rs`, `digital.rs`, `lookback.rs`, `barrier.rs`, `asian.rs`, `basket.rs`, `rainbow.rs`, `kirk.rs`, `spread.rs`, `dupire.rs`, `merton_jump.rs`** — exotic / less-critical pricers. Surveyed for panics (counted in §1) but not formula-audited.
6. **`calibration/hw_swaption.rs`** — Hull-White swaption calibration is a different beast (yield curves, Bermudan exercise). Spot-checked panics, did not audit formulas.
7. **`calibration/rbergomi.rs` Wasserstein objective + mSOE simulation** — read top-level flow only; the mSOE kernel approximation is not audited against Bayer/Friz/Gatheral or McCrickerd/Pakkanen.
8. **`pricing/variance_swap.rs`** — replication-weight formula not audited (CBOE-style log-strip, well known but should be verified for edge cases).
9. **`pricing/finite_difference.rs`** — PDE solver discretization not audited.
10. **The full test suite was not run** (instructions said "don't run the 685s suite — just confirm it compiles"). I confirmed compile + ran `calibration::levy` (6/6 ok) and `calibration::heston_stoch_corr` (1/1 ok in 617s — note: that single test contributed to the 685s wall-clock; consider tagging `#[ignore]` for nightly-only runs).

---

## §6 Recommended next steps

**Before v2.0.0 stable cut (must-fix)**:
1. **F2**: fix `calibrate_hscm` — capture the slsqp return value. ETA: 5 min, 1-line fix + add a recovery test.
2. **F3**: add `NigFourier` model and rewire `LevyModel::Nig`. ETA: 30 min. Or temporarily mark `LevyModelType::Nig` as `unimplemented!` for `to_model` until the Fourier model exists.
3. **F4**: thread `q` through `HestonStochCorrPricer::char_func_complex` and `HscmModel::price_call`. ETA: 30 min, plus add a parity test with q ≠ 0.
4. **F5**: add `q` parameter to `RBergomiPricer::price_call` propagation and `simulate_rbergomi_terminal_samples`. ETA: 1 hour, plus update calibration callers.
5. **F6**: change `CarrMadanPricer::price_call` extrapolation behaviour (saturate at endpoints + log warning, or return `NaN`). ETA: 15 min.
6. **F8**: add the 4 missing `Calibrator` trait impls. ETA: 1 hour.

**Before stable, but not blocking rc.1**:
7. **F1 follow-up**: tighten the Cui Jacobian numerical-consistency test from `2e-1` to `5e-3`. If it breaks, *then* the bug is real.
8. **F7**: switch every Heston-family ChF to the "Little Trap" form (Albrecher 2007). Ship a test on a long-maturity high-|ρ| case (T = 5y, ρ = -0.9) verifying the original form fails and the Trap form succeeds.
9. **F10**: fix the `SabrParams ↔ DVector` round-trip — keep β in the vector or rename to `SabrLognormal`.
10. **F11**: convert `expect`-on-user-input to `Result`-returning constructors.

**Post-stable (v2.1)**:
11. **F12**: rationalise tau / dates handling (`Result` instead of panic).
12. **F13**: drop `&dyn FourierModelExt` in favour of generic dispatch.
13. **F15**: cumulant-based Carr-Madan grid centring (Lord & Kahl 2010).
14. **F18**: validate-on-construct in `FMVol` (avoid integer underflow on empty input).
15. Audit-deferred: complete §5 items 1–9.

**Test coverage gaps to close**:
- HSCM round-trip from synthetic prices (would catch F2, F4).
- NIG calibration → pricing round-trip (would catch F3).
- Heston / Bates pricing at T = 5y, ρ = -0.9 (would expose F7).
- Carr-Madan with `k = S/100` and `k = 100·S` (would expose F6).
- rBergomi pricing with q = 0.05, comparison to monte-carlo with q in drift (would expose F5).
