OpenEvolve for Quasi–Monte Carlo Design

Quasi-Monte Carlo (QMC) methods generate "low-discrepancy" point sets — mathematical constructions that fill a multi-dimensional space more uniformly and consistently than the independent and identically distributed points used in standard Monte Carlo simulations. With these structured point sets, QMC allows researchers to solve complex integration problems (such as pricing financial options or simulating physical systems) with significantly lower variance, greater accuracy, and faster convergence. OpenEvolve helped advance this field by treating QMC design as a program synthesis problem: instead of humans manually tweaking formulas, an LLM-driven evolutionary search was used to write and refine code that generates these point sets. This approach successfully discovered new, record-breaking 2D and 3D point sets with lower "star discrepancy" (a measure of uniformity) than previously known configurations. Additionally, OpenEvolve optimized Sobol' direction numbers, which are the "instructions" for generating high-quality sequences; the resulting AI-evolved parameters significantly reduced error in 32-dimensional financial option-pricing simulations compared to long-standing industry standards, such as Joe–Kuo. Full details in the paper and code.

Low-discrepancy point sets and digital sequences underpin quasi–Monte Carlo (QMC) methods for high-dimensional integration. We cast two long-standing QMC design problems as program synthesis. With OpenEvolve, we construct finite 2D/3D point sets with low star discrepancy and optimize Sobol' direction numbers that minimize randomized QMC (rQMC) error on downstream integrands, improving on the previous state of the art. On finite sets, we rediscover known optima and set new best-known 2D benchmarks for N≥40 and 3D benchmarks for N≥9. On digital sequences, evolved Sobol' parameters yield consistent reductions in randomized QMC integration error for several option-pricing tasks relative to widely used Joe–Kuo parameters.

I. Finite point sets: new low–star-discrepancy constructions

We first used OpenEvolve to construct point sets in [0, 1]^d with minimal star discrepancy. This objective is central in QMC as the Koksma–Hlawka inequality bounds deterministic QMC integration error by the product of the integrand's variation and the point set's star discrepancy.

Two-phase search

We used a two-phase procedure:

Phase I (Direct construction): the LLM proposes explicit constructions that directly output points (seeded in 2D by a shifted Fibonacci lattice and in 3D by a scrambled Sobol' construction).

Phase II (Iterative optimization): the LLM writes programs that generate an initialization and then refine it using iterative optimization (e.g., SLSQP), including restart logic and perturbation heuristics.

Figure 1: Visualization of N=16 point set generation in two dimensions. (A) Initial shifted Fibonacci lattice (Discrepancy: 0.0962). (B) Best direct construction found in Phase I (Discrepancy: 0.0924). (C) Final optimized point set from Phase II (Discrepancy: 0.0744), which is within 0.68% of the known optimal value of 0.0739.

A representative example appears at N=16 in 2D (Fig. 1). Starting from a shifted Fibonacci lattice with discrepancy 0.0962, Phase I finds a slightly improved construction at 0.0924; Phase II refines further to 0.0744, within 0.68% of the known optimum 0.0739.

Point set results

Figure 2: The Star Discrepancy of Sobol', Halton, Hammersley, Fibonacci, Rank-1-Lattice, MPMC (message passing Monte Carlo), and OpenEvolve (LLM) sets for increasing number of points N=100 to 1020 in 2D.

We benchmark against Fibonacci lattices, MPMC (message-passing Monte Carlo), and provably optimal sets where available. For 2D point sets with N≥40, the evolutionary search finds configurations with lower star discrepancy than best-known values in the literature (Fig. 2). For instance, at N=100 we obtain a star discrepancy 0.0150, improving over the previous best 0.0188, achieved by MPMC (Table 1).

Table 1: 2D Star Discrepancy Comparison for N=2⋯100 between Fibonacci, MPMC (Message-Passing Monte Carlo), LLM evolutionary search, and Clément et al. (provably optimal for N≤20). Best values are bolded.

In 3D, provably optimal solutions are known only up to N≤8. We match most known optima through this frontier and report improved benchmarks beyond it with explicit constructions reported up to N=16 (Table 2).

Table 2: 3D Star Discrepancy for N=1⋯8. Lower is better; best per N is bolded.

Ablation: beyond local optimization alone

To test whether Phase II achieves its gains by just applying standard SciPy minimization schemes, we ran SLSQP-only baselines initialized from strong seeds: Sobol', Fibonacci, and best Phase I sets with the same objective and similar iteration budgets. We found the LLM-evolved Phase II programs still achieve lower star discrepancy in 22 of 24 tested scenarios, with reductions up to 15% relative to the best SLSQP-only baseline. This suggests the gains are driven by nontrivial program structure—initializations and restarts—rather than simply plugging in an optimizer.

II. Sobol' direction numbers: improving rQMC error in 32D option pricing

Next, we turned to optimizing Sobol' direction numbers. Sobol' sequences are widely used due to their distribution properties and efficient generation, but performance depends critically on direction numbers, and finding strong parameters is a difficult combinatorial search problem.

Program synthesis for direction numbers

Here, each program returns a list of Sobol' parameters, initialized from the widely used Joe–Kuo (2008) direction numbers. We minimize the MSE of a rQMC estimator of a 32D Asian option price for N=8192, averaged over 1000 consistent randomizations. After several hundred iterations, the evolved solution makes sparse modifications: it updates only dimensions 4, 5, and 6, leaving all other dimensions (up to 32) identical to Joe–Kuo.

Results and validation

Figure 3: The % reduction in MSE (rQMC integration over 10000 random scrambles and shifts) using Sobol' direction numbers found via LLM evolutionary search vs. those of Joe & Kuo.

We validate on five common Asian option scenarios: out-of-the-money, in-the-money, at-the-money, low volatility, and high volatility. The evolved direction numbers yielded significantly lower MSE for N≥512 across scenarios. We found similar results for various other option types: Lookback, Barrier, Basket, and Bermudan.

We also compared our parameters against Sobol' digital nets constructed using LatNetBuilder (random-CBC with a projection-dependent t-value criterion in (d=32)). For all five Asian-style payoffs and N≤1024, the evolved direction numbers outperform LatNetBuilder constructions with differences diminishing at larger N where designs converge. Future work could integrate LatNetBuilder as a callable tool.

Takeaways

Unlike fixed approaches where a new, computationally expensive optimization must be performed from scratch for each value of N, a single set of discovered direction numbers can be used to generate a high-quality point set for any desired number of points. Furthermore, the optimized Sobol' sequence can make use of standard randomization techniques, such as Owen scrambling, to obtain unbiased error estimates of rQMC. Digital sequences extend naturally to higher dimensions. By contrast, direct point set construction becomes difficult to evaluate: exact star discrepancy computation becomes computationally intractable as d grows.

In practice, OpenEvolve could serve as a promising tool for optimizing digital nets for application-specific integrands. We found OpenEvolve to be surprisingly cheap: we used roughly 4,000 input and 2,000 output tokens per call with 2,000 calls to Gemini 2.0 Flash, roughly $3.60 in total. In our case, much of the computation went into verifying the solution: computing star discrepancy or integration error.