Merge pull request #8 from dancixx/feat/rework-api-and-performance

feat: add rbergomi

Cargo.toml

@@ -1,3 +1,8 @@
 [workspace]
-members = ["stochastic-rs-core", "stochastic-rs-quant", "stochastic-rs-stats"]
+members = [
+    "stochastic-rs-core",
+    "stochastic-rs-ml",
+    "stochastic-rs-quant",
+    "stochastic-rs-stats",
+]
 resolver = "2"

stochastic-rs-core/src/lib.rs

@@ -49,7 +49,7 @@ pub trait Sampling<T: Clone + Send + Sync + Zero>: Send + Sync {
       panic!("m must be specified for parallel sampling");
     }
 
-    let mut xs = Array2::zeros((self.m().unwrap(), self.n()));
+    let mut xs = Array2::zeros((self.m().unwrap(), self.n() + 1));
 
     xs.axis_iter_mut(Axis(0)).into_par_iter().for_each(|mut x| {
       x.assign(&self.sample());

stochastic-rs-core/src/main.rs

@@ -1,24 +1,31 @@
-use stochastic_rs::process::fbm::Fbm;
+use stochastic_rs::{process::fbm::Fbm, Sampling};
 
 fn main() {
-  // let fbm = Fbm::new(0.75, 10000, Some(1.0), None);
-  // let mut runs = Vec::new();
+  let fbm = Fbm::new(&Fbm {
+    hurst: 0.9,
+    n: 10000,
+    t: None,
+    m: Some(100),
+    ..Default::default()
+  });
+  println!("{:?}", fbm.fgn.hurst);
+  let mut runs = Vec::new();
 
-  // for _ in 0..20 {
-  //   let start = std::time::Instant::now();
-  //   for _ in 0..1000 {
-  //     let _ = fbm.sample();
-  //   }
+  for _ in 0..20 {
+    let start = std::time::Instant::now();
+    for _ in 0..1000 {
+      let _ = fbm.sample_par();
+    }
 
-  //   let duration = start.elapsed();
-  //   println!(
-  //     "Time elapsed in expensive_function() is: {:?}",
-  //     duration.as_secs_f32()
-  //   );
-  //   runs.push(duration.as_secs_f32());
-  // }
+    let duration = start.elapsed();
+    println!(
+      "Time elapsed in expensive_function() is: {:?}",
+      duration.as_secs_f32()
+    );
+    runs.push(duration.as_secs_f32());
+  }
 
-  // let sum: f32 = runs.iter().sum();
-  // let average = sum / runs.len() as f32;
-  // println!("Average time: {}", average);
+  let sum: f32 = runs.iter().sum();
+  let average = sum / runs.len() as f32;
+  println!("Average time: {}", average);
 }

stochastic-rs-core/src/noise/fgn.rs

@@ -21,7 +21,7 @@ pub struct Fgn {
 
 impl Default for Fgn {
   fn default() -> Self {
-    Self::new(0.5, 1024, None, None)
+    Self::new(0.7, 1024, None, None)
   }
 }
 
@@ -112,7 +112,7 @@ mod tests {
 
   #[test]
   fn plot() {
-    let fgn = Fgn::new(0.7, 1024, Some(1.0), None);
+    let fgn = Fgn::new(0.7, 1024, Some(1.0), Some(1));
     let mut plot = Plot::new();
     let d = fgn.sample_par();
     for data in d.axis_iter(Axis(0)) {

stochastic-rs-core/src/process/fbm.rs

@@ -17,7 +17,7 @@ impl Fbm {
       panic!("Hurst parameter must be in (0, 1)")
     }
 
-    let fgn = Fgn::new(params.hurst, params.n, params.t, params.m);
+    let fgn = Fgn::new(params.hurst, params.n, params.t, None);
 
     Self {
       hurst: params.hurst,
@@ -33,10 +33,10 @@ impl Sampling<f64> for Fbm {
   fn sample(&self) -> Array1<f64> {
     let fgn = self.fgn.sample();
 
-    let mut fbm = Array1::<f64>::zeros(self.n);
+    let mut fbm = Array1::<f64>::zeros(self.n + 1);
     fbm.slice_mut(s![1..]).assign(&fgn);
 
-    for i in 1..self.n {
+    for i in 1..(self.n + 1) {
       fbm[i] += fbm[i - 1];
     }
 
@@ -62,10 +62,10 @@ mod tests {
   #[test]
   fn plot() {
     let fbm = Fbm::new(&Fbm {
-      hurst: 0.5,
+      hurst: 0.9,
       n: 1000,
       t: Some(1.0),
-      m: None,
+      m: Some(1),
       ..Default::default()
     });
     let mut plot = Plot::new();

stochastic-rs-core/src/volatility.rs

@@ -1,33 +1,4 @@
-//! This module contains the implementations of various diffusion processes.
-//!
-//! The following diffusion processes are implemented:
-//!
-//! - Duffie-Kan Model
-//!   - The Duffie-Kan model is a multifactor interest rate model incorporating correlated Brownian motions.
-//!   - SDE: `dr(t) = (a1 * r(t) + b1 * x(t) + c1) * dt + sigma1 * (alpha * r(t) + beta * x(t) + gamma) * dW_r(t)`
-//!   - SDE: `dx(t) = (a2 * r(t) + b2 * x(t) + c2) * dt + sigma2 * (alpha * r(t) + beta * x(t) + gamma) * dW_x(t)`
-//!   - where Corr(W_r(t), W_x(t)) = rho
-//!
-//! - **Heston Model**
-//!   - A stochastic volatility model used to describe the evolution of the volatility of an underlying asset.
-//!   - SDE: `dS(t) = mu * S(t) * dt + S(t) * sqrt(V(t)) * dW_1(t)`
-//!   - SDE: `dV(t) = kappa * (theta - V(t)) * dt + eta * sqrt(V(t)) * dW_2(t)`
-//!
-//! - SABR (Stochastic Alpha, Beta, Rho) Model
-//! - Widely used in financial mathematics for modeling stochastic volatility.
-//!   - SDE: `dF(t) = V(t) * F(t)^beta * dW_F(t)`
-//!   - SDE: `dV(t) = alpha * V(t) * dW_V(t)`
-//!   - where Corr(W_F(t), W_V(t)) = rho
-//!
-//! - **Vasicek Model**
-//!   - An Ornstein-Uhlenbeck process used to model interest rates.
-//!   - SDE: `dX(t) = theta * (mu - X(t)) * dt + sigma * dW(t)`
-//!
-//! - **Fractional Vasicek (fVasicek) Model**
-//!   - Incorporates fractional Brownian motion into the Vasicek model.
-//!   - SDE: `dX(t) = theta * (mu - X(t)) * dt + sigma * dW^H(t)`
-//!
-//! Each process has its own module and functions to generate sample paths.
-
+pub mod fheston;
 pub mod heston;
+pub mod rbergomi;
 pub mod sabr;

stochastic-rs-core/src/volatility/rbergomi.rs

@@ -0,0 +1,71 @@
+use ndarray::{s, Array1};
+
+use crate::{noise::cgns::Cgns, Sampling2D};
+
+#[derive(Default)]
+pub struct RoughBergomi {
+  pub hurst: f64,
+  pub nu: f64,
+  pub sigma_0: f64,
+  pub r: f64,
+  pub rho: f64,
+  pub n: usize,
+  pub t: Option<f64>,
+  pub m: Option<usize>,
+  pub cgns: Cgns,
+}
+
+impl RoughBergomi {
+  #[must_use]
+  pub fn new(params: &Self) -> Self {
+    let cgns = Cgns::new(&Cgns {
+      rho: params.rho,
+      n: params.n,
+      t: params.t,
+      m: params.m,
+    });
+
+    Self {
+      hurst: params.hurst,
+      nu: params.nu,
+      sigma_0: params.sigma_0,
+      r: params.r,
+      rho: params.rho,
+      n: params.n,
+      t: params.t,
+      m: params.m,
+      cgns,
+    }
+  }
+}
+
+impl Sampling2D<f64> for RoughBergomi {
+  fn sample(&self) -> [Array1<f64>; 2] {
+    let dt = self.t.unwrap_or(1.0) / self.n as f64;
+    let [cgn1, z] = self.cgns.sample();
+
+    let mut s = Array1::<f64>::zeros(self.n + 1);
+    let mut sigma2 = Array1::<f64>::zeros(self.n + 1);
+
+    for i in 1..(self.n + 1) {
+      s[i] = s[i - 1] + self.r * s[i - 1] + sigma2[i - 1] * cgn1[i - 1];
+
+      let sum_z = z.slice(s![..i]).sum();
+      let t = i as f64 * dt;
+      sigma2[i] = self.sigma_0.powi(2)
+        * (self.nu * (2.0 * self.hurst).sqrt() * t.powf(self.hurst - 0.5) * sum_z
+          - 0.5 * self.nu.powi(2) * t.powf(2.0 * self.hurst))
+        .exp();
+    }
+
+    [s, sigma2]
+  }
+
+  fn n(&self) -> usize {
+    self.n
+  }
+
+  fn m(&self) -> Option<usize> {
+    self.m
+  }
+}

stochastic-rs-ml/Cargo.toml

@@ -0,0 +1,6 @@
+[package]
+name = "stochastic-rs-ml"
+version = "0.1.0"
+edition = "2021"
+
+[dependencies]

stochastic-rs-ml/src/lib.rs

@@ -0,0 +1,14 @@
+pub fn add(left: u64, right: u64) -> u64 {
+  left + right
+}
+
+#[cfg(test)]
+mod tests {
+  use super::*;
+
+  #[test]
+  fn it_works() {
+    let result = add(2, 2);
+    assert_eq!(result, 4);
+  }
+}