Source-linked AI summary

Single- and Multilevel Quadrature with Error Control for Fourier Pricing under the Rough Heston Model

Chiheb Ben Hammouda, Abderrahmene Ben Romdhane, Michael Samet, Raul F. Tempone

arXiv:2609.00438v1q-fin.CPmath.NA

TL;DR

Fourier pricing under rough Heston is costly because every quadrature point requires a fractional Riccati solve, while a uniform time grid may waste work. The paper develops scaled single-level and multilevel Gauss-Laguerre methods that jointly allocate time-discretization and quadrature effort. Under stated assumptions, multilevel work scales as O(ε^-β/p), and experiments report lower total CPU time than BL2 in tested configurations.

  • Problem

    Rough-Heston Fourier pricing requires solving a fractional Riccati equation at every quadrature point, making a uniform time discretization potentially inefficient.

  • Method

    The paper develops scaled single-level and multilevel Gauss-Laguerre methods that jointly allocate fractional-Riccati time-discretization and Fourier-quadrature work.

  • Results

    O(ε^-β/p) computational work is proved for the multilevel method under stated assumptions, and the tested configurations report lower total CPU time than the BL2 Markovian approximation.

  • Takeaways & Limitations

    The multilevel construction provides computational savings over the single-level method while retaining the exact rough Heston fractional-Riccati representation.

  • Takeaways & Limitations

    The practical implementation replaces level-dependent parameters with empirically shared parameters, which is not part of the certified bounds; the relation between smoothness s_L and roughness α remains to be understood.

Abstract

from arXiv · show

Unlike the classical Heston model, Fourier pricing under the rough Heston model requires solving a fractional Riccati equation at every quadrature point. Since the required resolution varies with model parameters and quadrature point, a single uniform time discretization can be inefficient. We develop single- and multilevel Gauss-Laguerre quadrature methods that balance the time discretization and Fourier quadrature errors. Both methods scale the laguerre weight to the estimated Fourier integrand decay. The single-level method allocates a prescribed tolerance between the two errors. The multilevel method splits the integrand into a level-zero term and level differences, selecting quadrature points separately at each level. Suppose that the Fourier integrand discretization error is $O(Δt^p)$, that evaluating the characteristic function once costs $O(Δt^{-β})$, and that the algebraic Gauss-Laguerre quadrature error is $O(N^{-s_{SL}/2})$, where $s_{SL}$ is the smoothness index. Under this estimate and assumptions on the regularity and decay of level differences, we prove that the proposed single-level method requires $O(ε^{-(β/p+2/s_{SL})})$ computational work to achieve accuracy $ε$, whereas the proposed multilevel method requires $O(ε^{-β/p})$ computational work. We also study root-exponential Gauss-Laguerre error models for practical multilevel quadrature allocation. Numerical experiments support the observed fractional Riccati and Fourier integrand convergence rates and root-exponential quadrature behavior, and show substantial reductions in quadrature cost from the proposed scaling. The multilevel method provides clear computational savings over the single-level method. We further benchmark the multilevel fractional Riccati method against the BL2 Markovian approximation and report lower total CPU time in the tested configurations.

1 Introduction

Fourier pricing remains attractive for rough Heston, but each Fourier node requires a fractional Riccati solve whose cost and resolution can vary. The paper jointly controls time-discretization and quadrature errors through scaled single-level and multilevel Gauss-Laguerre methods, with theory and experiments indicating lower multilevel cost.

  • Motivation: Rough volatility improves short-maturity skew modeling but removes Markovianity, complicating PDE pricing and increasing the cost of path-dependent Monte Carlo.The rough Heston characteristic function nevertheless has a semi-explicit exponential-affine representation, preserving the appeal of Fourier pricing.
  • Open challenges: The paper addresses unresolved rough-Heston issues including unavailable explicit analyticity strips and insufficiently characterized large-Fourier-variable characteristic-function asymptotics.These gaps complicate contour selection, damping-parameter tuning, and quadrature design.
  • Contribution: The paper retains the exact rough Heston model and jointly allocates work between fractional-Riccati time discretization and Fourier quadrature.The required quadrature points are selected using an error estimate that accounts for the time step and characteristic-function evaluation cost.
  • Multilevel method: The multilevel method telescopically decomposes the discretized integrand and allocates quadrature points separately across a level-zero term and level differences.Its allocation exploits decay of the multilevel corrections.
  • Single-level method: The single-level method uses scaled Gauss-Laguerre quadrature, splitting the prescribed tolerance between time-discretization and Fourier quadrature errors.The Laguerre weight reflects estimated Fourier-integrand decay, and both algebraic and root-exponential error models are considered.
  • Numerical evidence: Numerical experiments support the predicted fractional-Riccati and Fourier-integrand convergence rates, root-exponential quadrature behavior, and substantial quadrature-cost reductions.The experiments compare single-level and multilevel work and CPU time, and include a BL2 benchmark comparison.

3 Methodology and Numerical Analysis

The paper develops error-controlled single-level and multilevel Gauss-Laguerre methods for rough Heston Fourier pricing, jointly managing fractional Riccati time discretization and Fourier quadrature. The multilevel construction distributes quadrature work across a level-zero term and level differences, while practical parameter selection uses asymptotic indicators and fitted quadrature models.

  • At each quadrature point, the characteristic function is computed by numerically solving the fractional Riccati-Volterra equation.
  • The fully discrete characteristic function combines nodal Riccati approximations with numerical time integration, producing the implemented Fourier integrand at each level.
  • The single-level method couples the time-discretization level and quadrature-point count because every quadrature point requires one fractional Riccati solve.
  • Scaled Gauss-Laguerre quadrature: Scaled Gauss-Laguerre quadrature replaces e^-u with e^-σ̃u, choosing σ̃ heuristically to match the estimated Fourier-integrand decay.
  • Practical parameter selection: Practical level selection uses Richardson-type indicators, but the asymptotic indicator is not certified unless early levels are asymptotic and their quadrature errors are sufficiently small.
  • Quadrature error models: The algebraic quadrature estimate is conservative relative to the root-exponential estimate, selecting more points for sufficiently small ε_quad when the root-exponential model applies.
  • Multilevel method: The multilevel method splits the integrand into a level-zero term and level differences, applying scaled quadrature separately across the hierarchy.

4 Numerical Results

Numerical diagnostics support the proposed convergence assumptions, scaling choices, and multilevel allocation across the tested rough Heston configurations. The multilevel method consistently reduces computational work and generally lowers CPU time, while several diagnostics remain finite-range evidence rather than proofs.

  • Experimental setup: The experiments use European calls with the EuRos and SPY rough Heston parameter sets as representative configurations.The parameter sets originate from calibrations to market implied-volatility data.
  • Time-discretization diagnostics: The fitted fractional-Riccati and Fourier-integrand rates are broadly consistent with 1 + α in the tested configurations.Reported rates include approximately 1.617–1.622 and 1.714–1.717 for one diagnostic, and 1.595–1.596 and 1.668–1.756 for another.
  • Time-discretization diagnostics: 2.5215 × 10^-7%–5.5594 × 10^-7% EuRos quadrature-to-level-difference ratios show negligible quadrature error relative to the level difference.The corresponding SPY ratios range from 0.1029% to 2.1979%.
  • Quadrature scaling: The scaled Gauss–Laguerre rule converges faster than the standard rule at short and long maturities and performs comparably at T = 1.This behavior is observed for both the level-zero integrand and the first level difference.
  • Quadrature error models: Root-exponential error estimates describe observed quadrature errors more closely than algebraic estimates over the tested quadrature orders.Algebraic estimates remain useful because they support explicit optimization and multilevel work-complexity analysis.
  • Single-level versus multilevel methods: 1.458 for SL and 1.237 for ML are the fitted work exponents, close to predicted exponents 1.485 and 1.235, respectively.Both methods satisfy the prescribed relative tolerance, and multilevel work grows more slowly as tolerance is reduced.
  • Single-level versus multilevel methods: More than ten times lower CPU time is achieved by the multilevel method at the smallest tested tolerance in the challenging short-maturity SPY case.The comparison provides numerical evidence rather than an asymptotic complexity result for root-exponential allocation.
  • Single-level versus multilevel methods: Across four configurations, ML is substantially faster in total CPU time and generally faster in pricing CPU time than SL.The pricing-time exception occurs at the tightest tolerance for the one-week SPY case.

5 Conclusion and Future Work

The paper develops a hierarchical Fourier-pricing framework for rough Heston that jointly controls Riccati time-discretization and Fourier quadrature costs. Multilevel allocation reduces expensive fine-grid solves, with lower total CPU time than BL2 in the reported configurations.

  • The framework evaluates the fractional Riccati equation hierarchically across Fourier quadrature levels, rather than using one finest grid for every node.Both methods use scaled Gauss-Laguerre quadrature matched to estimated Fourier-integrand decay.
  • The single-level complexity is algebraic, while the multilevel complexity reduces to O(ε^-β/p) when s0 ≥ 2p/β.The simplified multilevel estimate assumes the stated regularity and correction-decay conditions.
  • The empirical fractional Adams convergence rate is p ≈ 1 + α for the tested numerical configurations.The theoretical complexity analysis keeps p generic, while experiments use the observed rate.
  • The multilevel implementation has lower total CPU time than BL2 in the four reported SPY and EuRos benchmark configurations.The comparison uses the reference-assisted BL2 factor-count procedure.
  • Future work includes adaptive time stepping, alternative quadrature rules, and studying the relationship between Fourier-integrand smoothness and roughness.The proposed framework is also described as adaptable to other semi-explicit characteristic-function models.

Use of artificial intelligence

The authors used AI assistance for writing clarity, presentation, and reference identification while retaining responsibility for the paper's mathematical and empirical content.

  • AI assistance supported writing clarity, presentation, and identification of relevant references.
  • The authors retain responsibility for the mathematical content, numerical results, and arguments.

A Proof of Lemma 1

The proof fixes a Fourier variable, introduces the time-integration exponent based on the exact Riccati solution, and establishes bounds through recursive level differences.

  • The proof fixes u ≥ 0 and sets ξ = u + iR before applying the time-integration rule at grid points.
  • Because the initial exact and discrete Riccati values coincide, the initial level difference is zero.
  • The proof controls each level difference using linear and quadratic terms involving the Riccati coefficient and the bounded exact solution.
  • The resulting estimate is combined with the preceding relation to complete the level-difference bound.

B Proof of Proposition 2

The proof bounds characteristic-function error by the Fourier integrand's decay-weighted domination function and integrates this pointwise estimate over the Fourier variable.

  • Taking the supremum over levels establishes finiteness of Ψ(u), after which the characteristic-function error is expressed through Eℓ = G − Gℓ.
  • The characteristic-function error is bounded by |Φ(u + iR)| Ψ(u)∆t^pG.
  • The Fourier-integrand estimate follows from |Re(z)| ≤ |z| applied to the characteristic-function bound.
  • Integrability of the domination function then yields a finite integrated error bound.

C Proof of Proposition 3

The proof establishes Proposition 3 by invoking Assumption 9 and comparing the relevant equations.

  • Assumption 9 supplies the premise used in the proof.
  • The argument concludes Proposition 3 through this comparison.

D Proof of Proposition 6

The proof bounds single-level work by separating evaluation cost from quadrature-point cost, selecting the smallest admissible level, and enforcing both error constraints.

  • Single-level work is decomposed as W_SL = N_alg W_L.W_L denotes the cost of one evaluation at level L.
  • The smallest level L satisfying the discretization condition determines the time-step tolerance δ_ε.The coarser level L − 1 fails the discretization constraint.
  • The quadrature-point bound uses the estimated algebraic error and the relations A_L ≤ A and s_L ≥ s_SL.For sufficiently small ε, the proof also assumes A/ε_quad ≥ 1.
  • Choosing L and N_alg to satisfy both constraints yields |V − V_Nalg,L| ≤ ε.The conclusion follows by multiplying the two preceding estimates.

E Proof of Proposition 7

The proof of Proposition 7 splits the total tolerance equally between discretization and quadrature errors and uses the smallest admissible level.

  • The proof sets ε_disc = ε_quad = ε/2.
  • The selected level L is determined by the discretization constraint.
  • The asymptotic work estimate requires additional uniform bounds on A_L(ε) and B_L(ε).

F Proof of Proposition 8

The proof of Proposition 8 decouples level zero from correction levels, formulates the correction allocation as a convex optimization problem, and solves it using active constraints and stationarity.

  • The split constraints decouple the level-zero variable from the correction variables.At level zero, the proof formulates a separate optimization problem.
  • Because the objective increases with N_0, the level-zero constraint is active at the minimizer.
  • The correction-level problem is convex because N^(-s/2) is convex on (0, ∞).Strict feasibility allows the Karush-Kuhn-Tucker conditions to characterize the minimizer.
  • The correction constraint is active, and stationarity determines the optimal allocation.Positivity of c_ℓ and A_ℓ implies a positive multiplier λ.
  • Substitution into the active constraint gives the optimized correction-level expression.The proof then replaces c_ℓ with W_ℓ + W_{ℓ−1} and combines the resulting expression with N⋆.

G Proof of Proposition 9

The proof establishes the multilevel computational-work bound by combining level-wise quadrature estimates, refinement relations, and integer quadrature-point rounding. Under s0 ≥ 2p/β, the resulting complexity is O(ε^-β/p).

  • The condition s > 2p/β ensures the relevant exponent is negative under level refinement Δtℓ = Δt0 2^-ℓ.
  • The multilevel estimate combines the refinement relation Δtℓ−1 = 2Δtℓ with earlier bounds to prove the target result.
  • O(ε^-β/p) computational work follows when s0 ≥ 2p/β.The proof uses ε^-2/s0 = O(ε^-β/p) under this condition.

H Numerical evidence for Conjecture 1

Numerical experiments compare computable fine-grid reference quantities across levels and support the decay and integrability behavior required by Conjecture 1. The section also specifies an oracle-based BL2 benchmark whose factor count is selected by realized pricing error.

  • Numerical evidence for Conjecture 1: Observed stabilization of Dℓ,R and plateauing of Iℓ,R support the decay and integrability required in Conjecture 1.The experiment uses a finite Fourier interval, finitely many levels, and a finer-grid approximation rather than the exact characteristic function.
  • BL2 benchmark: The BL2 benchmark fixes each kernel rule before adaptively selecting Fourier truncation and Riccati and Fourier resolutions.Candidate rules are precomputed for n ∈ {1, 2, 3, 4}, while the internal half-tolerance excludes finite-factor kernel error.
  • BL2 benchmark: At most four BL2 factors were sufficient to achieve all tested tolerances.The minimum factor count is selected using the reference value Vref, making the procedure an oracle benchmark.
Loading 2609.00438v1…