moved function to own file + added sample variation

Cargo.toml

@@ -1,5 +1,5 @@
 [package]
-name = "MonteCarlo-Rust"
+name = "montecarlo-rust"
 version = "0.1.0"
 edition = "2021"
 

src/bounds.rs

@@ -12,19 +12,15 @@ impl Bounds {
         return Bounds {x, y};
     }
 
-    pub fn get_area(self) -> f64 {
+    pub fn get_area(&self) -> f64 {
         return self.x.get_length() * self.y.get_length();
     }
 
-    pub fn get_random_point(self) -> Point {
+    pub fn get_random_point(&self) -> Point {
         return Point::new(self.x.get_random_value(), self.y.get_random_value());
     }
 
-    pub fn get_x(self) -> LinearBounds {
+    pub fn get_x(&self) -> LinearBounds {
         return self.x;
     }
-
-    pub fn get_y(self) -> LinearBounds {
-        return self.y;
-    }
 }

src/function.rs

@@ -0,0 +1,61 @@
+use std::sync::Arc;
+use std::thread::JoinHandle;
+use crate::bounds::Bounds;
+use crate::linear_bounds::LinearBounds;
+use crate::point::Point;
+
+#[derive(Clone, Copy)]
+pub struct Function{}
+
+fn f(x: &f64) -> f64 {
+    return  x.powi(4) - 4_f64 * x.powi(3) + 18_f64 * x.powi(2) - 12_f64 * x - 69_f64
+}
+
+fn big_f(x: &f64) -> f64 {
+    return 0.2 * x.powi(5) - x.powi(4) + 6_f64 * x.powi(3) - 6_f64 * x.powi(2) - 69_f64 * x;
+}
+
+impl Function {
+    pub fn new() -> Function {
+        return Function {};
+    }
+
+    pub fn integrate(&self, bounds: &LinearBounds) -> f64 {
+        return big_f(&bounds.get_higher()) - big_f(&bounds.get_lower());
+    }
+
+    fn includes(&self, point: &Point) -> bool {
+        let y_0: f64 = f(&point.get_x());
+        return y_0.abs() > point.get_y().abs() && y_0 * point.get_y() >= 0_f64;
+    }
+
+    pub fn approximate(&self, bounds: &Bounds, samples: usize, thread_cnt: usize) -> f64 {
+        let mut hits: u64 = 0;
+        let mut threads: Vec<JoinHandle<u64>> = vec![];
+
+        let function: Arc<Function> = Arc::new(*self);
+        let bounds: Arc<Bounds> = Arc::new(*bounds);
+
+        for _ in 0..thread_cnt {
+            let function = function.clone();
+            let bounds = bounds.clone();
+
+            threads.push(std::thread::spawn(move || {
+                let mut hits: u64 = 0;
+                for _ in 0..samples {
+                    let point: Point = bounds.get_random_point();
+                    if function.includes(&point) {
+                        hits += 1;
+                    }
+                }
+                return hits;
+            }))
+        }
+
+        for thread in threads {
+            hits += thread.join().unwrap();
+        }
+
+        return (hits as f64 / (thread_cnt * samples) as f64) * bounds.get_area();
+    }
+}

src/linear_bounds.rs

@@ -12,20 +12,20 @@ impl LinearBounds {
         return LinearBounds {lower, higher};
     }
 
-    pub fn get_length(self) -> f64 {
+    pub fn get_length(&self) -> f64 {
         return (self.higher - self.lower).abs();
     }
 
-    pub fn get_random_value(self) -> f64 {
+    pub fn get_random_value(&self) -> f64 {
         let mut rand: ThreadRng = thread_rng();
         return self.lower + self.get_length() * rand.gen::<f64>();
     }
 
-    pub fn get_lower(self) -> f64 {
+    pub fn get_lower(&self) -> f64 {
         return self.lower;
     }
 
-    pub fn get_higher(self) -> f64 {
+    pub fn get_higher(&self) -> f64 {
         return self.higher;
     }
 }

src/main.rs

@@ -1,67 +1,35 @@
-use std::thread::JoinHandle;
-use std::time::Instant;
 use crate::bounds::Bounds;
+use crate::function::Function;
 use crate::linear_bounds::LinearBounds;
-use crate::point::Point;
-
-const SAMPLES: u64 = 1_000_000_000;
-const THREAD_CNT: u64 = 8;
+use crate::plot::Plot;
 
 mod point;
 mod linear_bounds;
 mod bounds;
+mod function;
+mod plot;
 
-fn f(x: f64) -> f64 {
-    return  x.powi(4) - 4_f64 * x.powi(3) + 18_f64 * x.powi(2) - 12_f64 * x - 69_f64
-}
-
-fn is_inside(point: Point) -> bool {
-    let y_0: f64 = f(point.get_x());
-    return y_0.abs() > point.get_y().abs() && y_0 * point.get_y() >= 0_f64;
-}
-
-fn get_real_integral(bounds: LinearBounds) -> f64 {
-    return big_f(bounds.get_higher()) - big_f(bounds.get_lower());
-}
-
-fn big_f(x: f64) -> f64 {
-    return 0.2 * x.powi(5) - x.powi(4) + 6_f64 * x.powi(3) - 6_f64 * x.powi(2) - 69_f64 * x;
-}
+const THREAD_CNT: usize = 8;
 
 fn main() {
-    let start: Instant = Instant::now();
-    let mut points_inside: u64 = 0;
     let bounds: Bounds = Bounds::new(LinearBounds::new(0_f64, 20_f64), LinearBounds::new(-100_f64, 150000_f64));
+    let func: Function = Function::new();
 
-    let mut threads: Vec<JoinHandle<u64>> = vec![];
+    let sample_limit: usize = 1_000;
+    let sub_samples: Vec<usize> = vec![1, 2, 5];
+    let samples_per_iteration: usize = 10;
 
-    for _ in 0..THREAD_CNT {
-        threads.push(std::thread::spawn(move || {
-            let mut points_inside: u64 = 0;
-            for _ in 0..SAMPLES {
-                let point: Point = bounds.get_random_point();
-                if is_inside(point) {
-                    points_inside += 1;
-                }
+    let mut result: Plot = Plot::new();
+
+    let mut sample_cnt: usize = 1;
+    while sample_cnt <= sample_limit {
+        for sub_sample in sub_samples.clone() {
+            let samples: usize = sub_sample * sample_cnt;
+            for _ in 0..samples_per_iteration {
+                result.add_sample(samples, func.approximate(&bounds, samples, THREAD_CNT));
             }
-            return points_inside;
-        }))
+        }
+        sample_cnt *= 10;
     }
-
-    for thread in threads {
-        points_inside += thread.join().unwrap();
-    }
-
-    let total_samples: u64 = SAMPLES * THREAD_CNT;
-
-    let approximated_integral: f64 = (points_inside as f64 / total_samples as f64) * bounds.get_area();
-    println!("The approximated Integral of the function is: {:.2}", approximated_integral);
-
-    let real_integral: f64 = get_real_integral(bounds.get_x());
-    println!("The real Integral of the function is: {:.2}", real_integral);
-
-    let error: f64 = (real_integral - approximated_integral).abs();
-    println!("That's an error of {:.2} or {:.5}%", error, (error/real_integral) * 100_f64);
-
-    println!("And the whole thing took {:.5} Seconds for {} samples", start.elapsed().as_micros() as f64 / 1000000_f64, total_samples);
-}
+    result.print();
+}

src/plot.rs

@@ -0,0 +1,29 @@
+use std::iter::zip;
+
+#[derive(Debug)]
+pub struct Plot {
+    x: Vec<usize>,
+    y: Vec<Vec<f64>>,
+}
+
+impl Plot {
+    pub fn new() -> Plot {
+        return Plot { x: vec![], y: vec![] }
+    }
+
+    pub fn add_sample(&mut self, sample: usize, value: f64) -> () {
+        if let Some(index) = self.x.iter().position(|&x| x == sample) {
+            self.y[index].push(value);
+        }
+        else {
+            self.x.push(sample);
+            self.y.push(vec![value]);
+        }
+    }
+
+    pub fn print(self) -> () {
+        for (x, y) in zip(self.x, self.y) {
+            println!("{}: {:?}", x, y);
+        }
+    }
+}

src/point.rs

@@ -9,11 +9,11 @@ impl Point {
         return Point {x, y};
     }
 
-    pub fn get_x(self) -> f64 {
+    pub fn get_x(&self) -> f64 {
         return self.x;
     }
 
-    pub fn get_y(self) -> f64 {
+    pub fn get_y(&self) -> f64 {
         return self.y;
     }
 }