integral_lib/
lib.rs

1#![cfg_attr(feature="simd", feature(portable_simd))]
2
3use std::mem::swap;
4
5#[cfg(feature = "simd")]
6mod simd;
7
8
9pub fn trapezoidal(mut a: f64, mut b: f64, steps: usize, f: fn(f64) -> f64) -> f64 {
10    if a > b {
11        swap(&mut a, &mut b);
12    }
13
14    let h = (b - a) / steps as f64;
15
16    let mut sum = (f(a) + f(b)) / 2f64;
17    let mut x = a + h;
18    while x < b {
19        sum += f(x);
20        x += h;
21    }
22
23    return h * sum;
24}
25
26pub fn simpsons_one_third(mut a: f64, mut b: f64, steps: usize, f: fn(f64) -> f64) -> f64 {
27    let stepsf64 = steps as f64;
28
29    if a > b {
30        swap(&mut a, &mut b);
31    }
32
33    let h = (b - a) / (stepsf64);
34
35    let mut sum = (f(a) + f(b)) / 2f64;
36    let mut s_2 = 0f64;
37    let mut s_4 = 0f64;
38    let mut x = a + h;
39    let mut i = 0;
40    while x < b {
41        if (i + 1) % 2 == 1 {
42            s_4 += f(x);
43        } else {
44            s_2 += f(x);
45        }
46        x += h;
47        i += 1;
48    }
49    sum += (s_2 * 2f64) + (s_4 * 4f64);
50
51    return (h * sum) / 3f64;
52}
53
54pub fn simpsons_three_eights(mut a: f64, mut b: f64, steps: usize, f: fn(f64) -> f64) -> f64 {
55    let stepsf64 = steps as f64;
56
57    if a > b {
58        swap(&mut a, &mut b);
59    }
60
61    let h = (b - a) / stepsf64;
62
63    let mut sum = f(a) + f(b);
64    let mut s_2 = 0f64;
65    let mut s_3 = 0f64;
66    let mut x = a + h;
67    let mut i = 0;
68    while x < b {
69        if (i + 1) % 3 == 0 {
70            s_2 += f(x);
71        } else {
72            s_3 += f(x);
73        }
74        x += h;
75        i += 1;
76    }
77    sum += (s_2 * 2f64) + (s_3 * 3f64);
78
79    return h * sum * 3f64 / 8f64;
80}
81
82#[cfg(test)]
83mod tests {
84
85    use super::*;
86    // tests are made for aproximating equation:
87    // pi = 2 * integral (-1, 1) sqrt(1-x^2) dx
88
89    #[test]
90    fn test_trapezoidal() {
91        println!(
92            "trapezoidal: {}",
93            trapezoidal(-1f64, 1f64, 600, |x| (1f64 - (x * x)).powf(0.5))
94        );
95    }
96
97    #[test]
98    fn test_simpsons_one_third() {
99        println!(
100            "one third: {}",
101            simpsons_one_third(-1f64, 1f64, 600, |x| (1f64 - (x * x)).powf(0.5))
102        );
103    }
104
105    #[test]
106    fn test_simpsons_three_eights() {
107        println!(
108            "three eights: {}",
109            simpsons_three_eights(-1f64, 1f64, 600, |x| (1f64 - (x * x)).powf(0.5))
110        );
111    }
112}
integral_lib/lib.rs

integral_lib/
lib.rs
