feat: add Sobol quasi-Monte Carlo sampler - #244
Conversation
Adds a scrambled Sobol low-discrepancy sampler, shaped like RNG, that plugs into MonteCarlo.integrate via the existing rng= slot. For smooth integrands the error shrinks close to O(1/N) instead of Monte Carlo's O(1/sqrt(N)). - torch backend uses torch.quasirandom.SobolEngine (native, moved onto the default device so GPU domains line up); numpy/jax/tensorflow use scipy.stats.qmc.Sobol converted to the backend. - Sample points are constants, so autodiff through the integrand and domain is preserved exactly as with plain Monte Carlo. - Exposed as torchquad.Sobol. Tests (all backends): analytic accuracy, Sobol beating plain MC at equal N, seed reproducibility, and a torch gradient-flow check. Addresses #140 (MonteCarlo path; VEGAS and JIT support to follow). Motivated by the QMC request in #217.
Coverage Report
|
||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
- Tighten the accuracy-test tolerance to 1e-5 so it fails if Sobol degrades to plain-MC quality (measured error ~2e-8 SciPy / ~5e-7 torch). - Warn uniformly on non-power-of-two point counts: torch's SobolEngine would otherwise silently drop the balance property that SciPy warns about. - Document the rng.py/qmc.py per-backend exception in CLAUDE.md rule 3. - Add the Sobol sampler to the CHANGELOG.
|
Claude: Thanks — addressed in 4bb604d:
|
|
Claude: ### Accuracy comparison: Monte Carlo vs Sobol QMC vs VEGAS Relative error against the analytic value at an equal sample budget of N = 2¹⁴ = 16384, numpy backend, float64, fixed
Takeaways
Caveat: VEGAS's |
|
Claude: Added a whole-collection accuracy test (aa5a5a6) — thanks for the nudge, the previous single cos-integrand check was too thin.
Also updated |
- Replace the single-integrand accuracy check with a collection test that runs every function in integration_test_functions (real + complex, incl. multi-dim) in 1-D/3-D/10-D via compute_integration_test_errors, with per-dimension bounds far tighter than plain Monte Carlo's (torch's native 10-D Sobol is weaker, so that bound is looser). - REVIEW.md: require validating integrators/samplers/error-estimates against the whole analytic collection, not one or two ad-hoc integrands.
# Conflicts: # CHANGELOG.md
What
First-class quasi-Monte Carlo via a scrambled Sobol sampler that plugs into
MonteCarlo.integratethrough the existingrng=slot — no changes toMonteCarloitself:For smooth integrands the error shrinks close to
O(1/N)instead of plain MonteCarlo's
O(1/sqrt(N)). In the tests, at N=2¹³ the Sobol error is ~3e-10 vs~1e-4…1e-3 for pseudo-random MC across all four backends.
torch.quasirandom.SobolEngine(native; points moved onto thecurrent default device so GPU domains line up, matching
RNG'storch.rand).scipy.stats.qmc.Sobol(SciPy is already a harddependency) converted to the backend.
is preserved exactly as for plain MC.
RNG's per-backend shape (the established exception to thebackend-agnostic rule for RNG-like samplers) and is exposed as
torchquad.Sobol.Notes / scope
integratepath; the JIT-compiled path builds its own RNGand is out of scope here.
reproducible for a fixed seed within a backend but not bit-for-bit across
backends (documented on the class).
Test plan
sobol_test.py, all backends: analytic accuracy, Sobol < plain-MC error atequal N (deterministic, fixed seed), seed reproducibility, torch gradient flow — 13 passed
monte_carlo_test.pystill green (19 passed together)ruff/pydoclint/vultureclean;sphinx-build -Wbuilds (Sobol auto-documented via__all__)Addresses #140 (MonteCarlo path; VEGAS + JIT support to follow). Motivated by #217.
Roadmap F1.