test: use the compute_error method to set tolerance for tests (#20)

* remove mutable borrow from compute_error method * correct compute_error with abs of the derivatives * add an `is_within_tolerance` function that uses error approximation * update `generate_test!` macro * fix mistake? in compute_error
imrn99 · Jun 2, 2024 · dd304f7 · dd304f7
1 parent e0c102c
commit dd304f7
Show file tree

Hide file tree

Showing 3 changed files with 44 additions and 34 deletions.
diff --git a/integraal/src/structure/implementations.rs b/integraal/src/structure/implementations.rs
@@ -212,7 +212,7 @@ impl<'a, X: Scalar> Integraal<'a, X> {
     /// This method returns a `Result` taking the following values:
     /// - `Ok(X: Scalar)` -- The computation was successfuly done
     /// - `Err(IntegraalError)` -- The computation failed for the reason specified by the enum.
-    pub fn compute_error(&mut self) -> Result<X, IntegraalError> {
+    pub fn compute_error(&self) -> Result<X, IntegraalError> {
         if self.domain.is_none() | self.function.is_none() | self.method.is_none() {
             return Err(IntegraalError::MissingParameters(
                 "one or more parameter is missing",
@@ -224,8 +224,8 @@ impl<'a, X: Scalar> Integraal<'a, X> {
             n_step,
         }) = self.domain
         {
-            // ref: https://en.wikipedia.org/wiki/Riemann_sum#Riemann_summation_methods
             let res: X = match self.method {
+                // ref: https://en.wikipedia.org/wiki/Riemann_sum#Riemann_summation_methods
                 Some(ComputeMethod::RectangleLeft | ComputeMethod::RectangleRight) => {
                     let m1: X = (1..n_step)
                         .map(|step_id| match &self.function {
@@ -239,11 +239,12 @@ impl<'a, X: Scalar> Integraal<'a, X> {
                             }
                             None => unreachable!(),
                         })
-                        .max_by(|t1, t2| t1.partial_cmp(t2).unwrap())
+                        .max_by(|t1, t2| t1.abs().partial_cmp(&t2.abs()).unwrap())
                         .unwrap();
                     let end = start + step * X::from_usize(n_step).unwrap();
                     m1 * (end - start).powi(2) / X::from_usize(2 * n_step).unwrap()
                 }
+                // ref: https://en.wikipedia.org/wiki/Trapezoidal_rule#Error_analysis
                 Some(ComputeMethod::Trapezoid) => {
                     let d1: Vec<X> = (1..n_step)
                         .map(|step_id| match &self.function {
@@ -260,10 +261,10 @@ impl<'a, X: Scalar> Integraal<'a, X> {
                         .collect();
                     let m2: X = (1..n_step - 2)
                         .map(|step_id| d1[step_id] - d1[step_id - 1] / step)
-                        .max_by(|t1, t2| t1.partial_cmp(t2).unwrap())
+                        .max_by(|t1, t2| t1.abs().partial_cmp(&t2.abs()).unwrap())
                         .unwrap();
                     let end = start + step * X::from_usize(n_step).unwrap();
-                    m2 * (end - start).powi(3) / X::from_usize(24 * n_step.pow(2)).unwrap()
+                    -m2 * (end - start).powi(3) / X::from_usize(24 * n_step.pow(2)).unwrap()
                 }
                 #[cfg(feature = "montecarlo")]
                 Some(ComputeMethod::MonteCarlo { .. }) => {

diff --git a/integraal/src/structure/tests.rs b/integraal/src/structure/tests.rs
@@ -2,7 +2,9 @@
 
 // ------ IMPORTS
 
-use crate::{ComputeMethod, DomainDescriptor, FunctionDescriptor, Integraal, IntegraalError};
+use crate::{
+    ComputeMethod, DomainDescriptor, FunctionDescriptor, Integraal, IntegraalError, Scalar,
+};
 
 // ------ CONTENT
 
@@ -12,12 +14,6 @@ const RECTANGLE_TOLERANCE: f64 = 1e-5;
 const TRAPEZOID_TOLERANCE: f64 = 1e-5;
 const STEP: f64 = 0.001;
 
-macro_rules! almost_equal {
-    ($v1: expr, $v2: expr, $tol: ident) => {
-        ($v1 - $v2).abs() < $tol
-    };
-}
-
 macro_rules! generate_sample_descriptors {
     ($f: ident, $d: ident, $c: ident) => {
         let $f: FunctionDescriptor<f64> = FunctionDescriptor::Closure(Box::new(|x| x));
@@ -103,6 +99,26 @@ fn inconsistent_parameters() {
 
 // correct usages
 
+fn is_within_tolerance<T: Scalar>(
+    mut integraal: Integraal<T>,
+    expected_result: T,
+) -> (bool, String) {
+    let tolerance = integraal.compute_error().unwrap();
+    let computed_result = integraal.compute().unwrap();
+    let sum = expected_result.abs() + computed_result.abs();
+    let delta = (computed_result - expected_result).abs();
+    (
+        delta < tolerance + T::epsilon() * sum.min(T::max_value()),
+        format!("computed value: {computed_result:?}\nexpected value {expected_result:?}\ncomputed tolerance: {tolerance:?}"),
+    )
+}
+
+macro_rules! almost_equal {
+    ($v1: expr, $v2: expr, $tol: ident) => {
+        ($v1 - $v2).abs() < $tol
+    };
+}
+
 // test are groups per module according to the integral & the computation method
 // test names follow this pattern:
 // <FunctionDescriptorEnum><DomainDescriptorEnum>
@@ -134,25 +150,20 @@ macro_rules! generate_test {
             );
         }
     };
-    ($name: ident, $fnd: expr, $dmd: expr, $met: expr, $tol: ident) => {
+    ($name: ident, $fnd: expr, $dmd: expr, $met: expr) => {
         #[allow(non_snake_case)]
         #[test]
         fn $name() {
             let functiond = $fnd;
             let domaind = $dmd;
             let computem = $met;
             let mut integraal = Integraal::default();
-            let res = integraal
+            integraal
                 .function(functiond)
                 .domain(domaind)
-                .method(computem)
-                .compute();
-            assert!(res.is_ok());
-            assert!(
-                almost_equal!(res.unwrap(), 2.0, $tol),
-                "left: {} \nright: 2.0",
-                res.unwrap()
-            );
+                .method(computem);
+            let (res, msg) = is_within_tolerance(integraal, 2.0);
+            assert!(res, "{msg}");
         }
     };
 }
@@ -180,8 +191,7 @@ mod a_rectangleleft {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::RectangleLeft,
-        RECTANGLE_TOLERANCE
+        ComputeMethod::RectangleLeft
     );
 
     generate_test!(
@@ -207,8 +217,7 @@ mod a_rectangleleft {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::RectangleLeft,
-        RECTANGLE_TOLERANCE
+        ComputeMethod::RectangleLeft
     );
 }
 
@@ -235,8 +244,7 @@ mod a_rectangleright {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::RectangleRight,
-        RECTANGLE_TOLERANCE
+        ComputeMethod::RectangleRight
     );
 
     generate_test!(
@@ -262,8 +270,7 @@ mod a_rectangleright {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::RectangleRight,
-        RECTANGLE_TOLERANCE
+        ComputeMethod::RectangleRight
     );
 }
 
@@ -290,8 +297,7 @@ mod a_trapezoid {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::Trapezoid,
-        TRAPEZOID_TOLERANCE
+        ComputeMethod::Trapezoid
     );
 
     generate_test!(
@@ -317,8 +323,7 @@ mod a_trapezoid {
             step: STEP,
             n_step: (1000. * std::f64::consts::PI) as usize,
         },
-        ComputeMethod::Trapezoid,
-        TRAPEZOID_TOLERANCE
+        ComputeMethod::Trapezoid
     );
 }
 

diff --git a/integraal/src/traits.rs b/integraal/src/traits.rs
@@ -1,11 +1,14 @@
 //! Common traits implementation
 
+use std::fmt::Debug;
+
 /// Scalar value trait.
 ///
 /// This trait is automatically implemented for all types implementing its requirements.
 pub trait Scalar:
     Clone
     + Copy
+    + Debug
     + num_traits::Float
     + num_traits::Signed
     + num_traits::FromPrimitive
@@ -18,6 +21,7 @@ pub trait Scalar:
 impl<
         X: Clone
             + Copy
+            + Debug
             + num_traits::Float
             + num_traits::Signed
             + num_traits::FromPrimitive