stochastic-rs: High-Performance Stochastic Process Simulation in Rust
by Daniel Boros
Jul 26
14 min read
1041 views
The crate stochastic-rs is an open-source Rust library for simulating and analyzing stochastic processes, designed with a focus on quantitative finance and high performance. The project was born out of a need to generate massive amounts of synthetic data for research. As a PhD candidate in financial mathematics (with a focus on fractional stochastic processes), I needed to simulate "tens to hundreds of millions of ... data points" to train deep neural networks. Using Python proved inadequate for such performance-intensive tasks, so Rust was chosen to build a faster solution. The result is a library that began as a university research project and has since expanded to include a broad range of stochastic process models and tools.
Why simulation speed matters: In quantitative finance and AI, one often needs to generate thousands of stochastic trajectories (random paths) repeatedly – for example, to price derivatives via Monte Carlo or to create training data for machine learning models. Every bit of performance improvement counts when each experiment may require millions of path simulations. Rust's promise of "zero-cost abstractions" (high-level code with C-like performance) and safety guarantees make it an attractive choice for this domain. The stochastic-rs crate leverages Rust’s strengths to achieve high throughput in simulations, enabling more complex models to be explored without prohibitive run times.
Performance Optimizations in Rust
From the outset, the library's development has been guided by a performance-oriented mindset. Here are some of the key optimizations and design choices that make stochastic-rs efficient:
-
Parallel Path Generation: The library uses multi-threading (via the Rayon crate) to generate multiple sample paths in parallel. In typical use, one might need to simulate, say, 10,000 independent paths of an asset price process; by distributing these across threads, the crate achieves near-linear speedups on multi-core CPUs. Thread-safe structures (e.g. using
Arc<Mutex<...>>
where needed) ensure data race freedom even when generating random variates concurrently. -
SIMD and Low-Level Optimizations: Many calculations (such as vectorized arithmetic for generating Gaussian noise or applying transformations) are accelerated using SIMD instructions. The crate taps into packages like
wide
(for explicit SIMD) and numeric libraries likendarray
andndrustfft
which internally leverage vectorization. Rust’sunsafe
capabilities are used judiciously to enable low-level optimizations without sacrificing safety in higher-level code. For example, generating fractional Gaussian noise (the increments of fractional Brownian motion) involves heavy linear algebra that benefits from FFTs and vectorized operations. -
Custom Memory Allocator: When dealing with large simulations, memory allocation patterns can impact performance. I experimented with using jemalloc (via the
tikv-jemallocator
crate) as the global allocator instead of Rust’s default. This, combined with compiler tuning (e.g. settingcodegen-units=1
to avoid code fragmentation), yielded a significant boost – cutting runtime for a large fBM simulation from ~1.12 seconds down to ~0.71 seconds. This illustrates how low-level tweaks (memory management and compile settings) can enhance throughput for computation-heavy tasks. -
GPU Offloading (CUDA): Certain processes can even be generated on the GPU. The library includes experimental support for NVIDIA GPUs using CUDA (via libraries like
cudarc
). For example, the fractional Gaussian noise generator can leverage GPU parallelism for the FFT-based algorithm (discussed below), which is beneficial when simulating extremely long paths or a very high number of paths simultaneously. In one of the project's blog posts, I detailed implementing a CUDA kernel to accelerate fractional noise generation for AI training data. (Support for other platforms, such as Apple’s Metal, is on the roadmap – but as of now only CUDA GPUs are supported.) -
Algorithmic Improvements: Beyond hardware considerations, the core algorithms themselves were refined for efficiency. Initially, a straightforward function was used to generate a fractional Brownian motion path, but this meant re-computing a large covariance matrix for each trajectory. By refactoring into a reusable struct (holding precomputed FFT components of the covariance), the library avoids redundant computations on each simulation. Now, generating subsequent paths (with the same parameters) is much faster – one only needs to sample new random numbers and do a few vector-matrix operations, rather than recompute the entire kernel from scratch each time.
All these optimizations underscore Rust's strength: one can fine-tune performance at every level, from high-level parallel iterators down to memory allocation and CPU/GPU instructions. In the words of the author, "the biggest breakthroughs came from careful parallelization, avoiding redundant computations, and tapping into Rust's compiler optimizations" – resulting in a library that pushes the limits of stochastic simulation speed.
Fractional Brownian Motion and Rough Paths in Finance
A major motivation for stochastic-rs was modeling fractional stochastic processes, especially fractional Brownian motion (fBM). Fractional Brownian motion is a generalization of the standard Brownian motion (the classical Wiener process). Whereas ordinary Brownian motion has independent increments and no "memory", fractional Brownian motion introduces a Hurst parameter that governs the roughness or persistence of the process. For , fBM reduces to standard Brownian motion. But if , the increments are no longer independent:
-
If , the process exhibits long-range dependence (also called persistence or "memory"). This means increments are positively correlated – if the process has been trending up, it tends to continue trending up, and vice versa. In practical terms, a high- fBM has a smoother path that wanders more steadily (displaying momentum).
-
If , the process exhibits anti-persistence, meaning increments are negatively correlated. An upward movement is likely to be followed by a downward movement (and vice versa), resulting in a more jagged, mean-reverting path. Low- fBM paths are rougher than standard Brownian motion.
-
yields the familiar memoryless Brownian motion (no correlation between increments).
For a more intuitive picture: one can imagine an fBM as a time series that "remembers" its past to a degree determined by . This makes fractional processes very appealing for financial modeling – empirical studies of asset prices and volatility have found evidence of long memory (e.g., in volatility clustering, where volatility levels are correlated over long periods). By using fractional Brownian motion or related rough processes, quants can build models that capture these long-memory effects better than classical Markovian models.
Mathematically, fBM has stationary increments and a covariance function:
which for yields correlated increments. Because fractional Brownian motion is not a semimartingale (except in the trivial case), one cannot apply the standard Itô calculus to it. This has spurred the development of alternative integration theories (e.g. fractional calculus, Wiener–Lie integration, or rough path theory in the context of rough volatility models). The rough path viewpoint, advocated in some recent finance research, posits that asset price dynamics might be better captured by processes with (very rough paths), as seen in the Rough Bergomi or Rough Heston models for volatility.
stochastic-rs provides tools to experiment with such models. In fact, fractional Brownian motion is something of a "hello world" for the library’s performance optimizations. I aimed to create "the fastest fractional Brownian motion generator available in the open-source world." Achieving this required choosing the right algorithm: the two classical exact simulation methods for fBM are (1) Cholesky decomposition of the covariance matrix (which has complexity and was deemed too slow), and (2) the Davies–Harte method, an FFT-based algorithm that runs in . The Davies–Harte algorithm was chosen for stochastic-rs. In essence, this method embeds the target covariance matrix into a larger circulant matrix, whose Fourier transform yields the eigenvalues (spectrum) of the covariance. By taking the square root of those eigenvalues, one can generate fractional Gaussian noise via frequency-domain sampling. In simpler terms: generate Gaussian noise in the frequency domain such that when transformed back to time domain, it has the desired covariance structure. The resulting fractional Gaussian noise (fGn) is then cumulated to produce a fractional Brownian motion path.
To put it concretely, if we want to simulate an fBM with Hurst and 1024 time steps, the library will:
- Compute the theoretical autocovariance vector
r
for fractional Gaussian noise increments of length 1024 (using the formula below for , with ). - Form a symmetric circulant vector by concatenating
r
with a mirrored copy of itself. - Perform an FFT on this vector to get eigenvalues, take their square roots (this gives the frequency-domain std-dev for each frequency).
- Generate a vector of independent standard Gaussians of the same length (frequency domain noise).
- Multiply (element-wise) the sqrt-eigenvalue vector with the Gaussian noise vector, then inverse FFT to get time-domain fractional Gaussian noise.
- Scale appropriately and take the cumulative sum to obtain the fractional Brownian motion path.
This approach is much faster than naive covariance matrix simulation for large n
. And because FFT operations are used, it also opens the door to GPU acceleration (as large FFTs can be offloaded to cuFFT on NVIDIA GPUs, for example).
Features and Modules of the Library
While fractional Brownian motion was a driving use case, stochastic-rs is a general toolkit that covers a wide variety of stochastic processes and models. It aims to be a one-stop library for quantitative finance simulations. Key features include:
-
Standard Processes (Diffusions): The library provides common continuous-time processes such as Gaussian white noise, correlated Gaussian noise, standard Brownian motion, correlated multi-dimensional Brownian motion, Geometric Brownian motion (the model underlying Black–Scholes), mean-reverting Ornstein–Uhlenbeck and Jacobi processes, and Cox–Ingersoll–Ross (CIR) processes. These can be used to model stock prices, interest rates, and other financial variables under various assumptions.
-
Jump and Lévy Processes: It includes processes with jumps (discontinuous paths) as well. For example, a Poisson process and Compound Poisson process (summing random jump sizes) are available, as building blocks for jump-diffusion models. There are also more advanced Lévy processes like Variance Gamma, Normal-Inverse-Gaussian, and others often used to capture heavy tails in asset returns. (Some of these are marked unstable as development is ongoing.)
-
Stochastic Volatility Models: Several well-known financial models that combine stochastic processes are implemented. These include the Heston model (stochastic volatility with mean reversion), the Merton and Bates models (which add jumps to the Black–Scholes framework), the Vasicek model (for interest rates), and even the SABR model for volatility smile dynamics. These models are useful for pricing derivatives and risk management, and having them in one library allows users to simulate and compare different scenarios easily.
-
Fractional (Rough) Processes: As expected, the library has a rich selection of fractional processes. Beyond fractional Brownian motion and fractional Gaussian noise, it includes fractional Ornstein–Uhlenbeck, fractional CIR, fractional Jacobi, and even fractional versions of the Heston model and other advanced models. The so-called rough Heston and rough Bergomi models are examples that incorporate fBM to model the persistence in volatility. With stochastic-rs, a researcher can simulate such rough volatility models out-of-the-box, which is quite cutting-edge in quantitative finance.
-
Discrete-time Models: For completeness, some discrete-time stochastic processes common in econometrics are included. For example, AR(1) (autoregressive) processes and GARCH-family models for volatility are available (these were mentioned in the talk but not explicitly listed in the docs snippet). Such models are useful for time series analysis and were added to extend the library’s utility beyond continuous SDEs.
-
Numerical Tools (SDE Solvers & Calibration): The library provides generic solvers for stochastic differential equations (SDEs), such as Euler–Maruyama method for simulating SDEs step by step. This means users can define a custom SDE (by specifying drift and diffusion functions) and simulate it without writing the solver logic themselves. Additionally, there are tools for optimization and calibration tasks. For instance, one can calibrate model parameters to market data by minimizing error metrics – the crate integrates with optimization libraries (the presence of
argmin
andlevenberg-marquardt
in the dependencies confirms this). There is also a module for implied volatility and option pricing, enabling the calculation of implied vol surfaces and pricing of derivative contracts given simulated underlying dynamics. This makes stochastic-rs not just about path generation, but also about connecting those simulations to practical financial outcomes (like pricing an option or measuring risk). -
Sensitivity Analysis: An interesting advanced feature mentioned is the use of Malliavin calculus techniques (the Malliavin derivative) to compute sensitivities (Greeks) or analyze how randomness affects outcomes. Malliavin calculus allows one to differentiate through random processes, providing a way to estimate sensitivities (like option price / noise). This is still a niche feature (and likely in early development), but it showcases the depth of the library’s offerings – going so far as to include research-level tools that are not commonly found in standard finance libraries.
All the above are built with performance in mind and composed in a modular way. The code is organized into modules (e.g. stochastic::processes
, stochastic::quant
, stochastic::stats
, etc.), with many processes implementing common traits for simulation. For example, there is a trait (e.g. Sampling
or similar) that processes implement to generate sample paths or random samples. This allows writing generic code over different processes.
Example: Simulating a Fractional Brownian Motion Path
To illustrate how one might use stochastic-rs, let's walk through a simple example in Rust. Suppose we want to simulate a single trajectory of fractional Brownian motion and examine its properties. We choose a Hurst parameter (which implies some long-range dependence) and simulate, say, 1000 time steps over the interval with (so ).
First, add stochastic-rs
to your Cargo.toml
dependencies (the crate is published on crates.io):
[dependencies]
stochastic-rs = "*"
Now in your Rust code:
use stochastic_rs::stochastic::fbm::FractionalBrownianMotion;
use stochastic_rs::traits::Sampler;
let fbm = FractionalBrownianMotion::new(0.7, 1000, 1.0);
let path = fbm.sample();
println!("Sample fBM path (H=0.7) first 5 points: {:?}", &path[..5]);
In this snippet, FractionalBrownianMotion::new(hurst, n, T)
constructs the process with the given Hurst exponent, number of steps n
, and total time horizon T
. Calling sample()
(provided by a common trait implemented for many process structs) generates a full path (for fBM, an array of length n+1
including the initial value). Internally, this will execute the FFT-based Davies–Harte algorithm described earlier to produce the path efficiently. We could then analyze the path
– for instance, check its mean, variance, plot it, etc., or use it as part of a larger simulation (like pricing a derivative that depends on an underlying following fBM).
For comparison, if we wanted to simulate a standard Brownian motion instead, the library offers a similar interface:
use stochastic_rs::stochastic::brownian::BrownianMotion;
let bm = BrownianMotion::new(0.0, 1.0);
let bm_path = bm.sample(n_steps, T);
This would generate a random walk with independent Gaussian increments. Many other processes (Geometric BM, Heston, etc.) can be created and simulated in an analogous way, with new(...)
constructors for parameters and typically a sample()
or simulate()
method to produce paths. The library's design tries to make switching out one model for another as easy as changing the type of process you instantiate.
Outlook and Community
Stochastic-RS is a young but actively developing project (open-sourced under MIT license). The latest versions have introduced more models and improvements, but also occasionally break the API as it matures – so expect some changes if you update the crate (the author notes that "breaking changes may occur" while the library is still under rapid development). Future plans include expanding GPU support (possibly adding Apple's Metal alongside CUDA), completing a copula module for modeling dependencies between multiple processes, and enhancing the Malliavin calculus tools for even more sophisticated sensitivity analysis. There is also an intention to provide Python bindings, which would open the library to a wider audience who may want to harness Rust's speed from Python environments.
The reception in industry circles so far has been encouraging. As Rust gains traction in performance-critical domains, there is growing interest in Rust-based quant finance tools. While many quantitative developers still use Python (for its rich ecosystem and quick prototyping) or C++ (for established libraries like QuantLib), Rust is emerging as an attractive middle-ground – offering the speed of C++ with safer concurrency and modern language ergonomics. stochastic-rs demonstrates that complex financial models can be implemented in Rust without sacrificing performance, and serves as a proof-of-concept that Rust can handle the heavy lifting in simulations, AI data generation, and risk calculations. Firms are beginning to explore Rust for these purposes, and having open-source libraries like this lowers the barrier to adoption.
Get involved: The project is open-source (available on GitHub at rust-dd/stochastic-rs and welcomes contributions, feedback, and feature requests. The developer is active in the Rust quantitative finance community (you can even find a Discord channel in the Rust–DD community for discussions). Users can refer to the documentation on Docs.rs and the examples provided there to get started. Given the continuous improvements (e.g., adding more tests and examples is on the to-do list), community involvement can help shape the library's direction.
In summary, stochastic-rs is an exciting development for computational finance. It brings together modern Rust engineering and advanced stochastic calculus, enabling high-performance simulation of everything from plain vanilla Brownian motions to exotic fractional volatility models. As the field moves towards more data-intensive and simulation-heavy approaches (like calibrating deep learning models on synthetic financial data or exploring millions of Monte Carlo scenarios), such tools will be invaluable. This project illustrates how choosing the right language and algorithms can push the frontier of what's possible in quantitative simulations – achieving in milliseconds what might otherwise take seconds or minutes – and it invites quant developers to ride the new wave of Rust in finance.