refactor!: replace f64 by generics Scalar (#13)

* add generic t(raits for domain & image values * add DomainValue to parameters enums * update sigs with the new trait * add image & integrated values generics to the main structure * add num dependency * fix implem code * fix tests * add minimal doc entry for the traits * merge all generic traits into on * update implem * fix warnings is this lint even working correctly?
imrn99 · May 21, 2024 · dbd2488 · dbd2488
1 parent 602f0aa
commit dbd2488
Show file tree

Hide file tree

Showing 8 changed files with 86 additions and 40 deletions.
diff --git a/Cargo.toml b/Cargo.toml
@@ -26,3 +26,4 @@ integraal-examples = { version = "0.0.2", path = "./examples" }
 # external
 rand = "0.9.0-alpha.1"
 rustversion = "1.0.15"
+num = "0.4.3"
diff --git a/integraal/Cargo.toml b/integraal/Cargo.toml
@@ -20,6 +20,7 @@ montecarlo = ["dep:rand"]
 # DEPS
 
 [dependencies]
+num.workspace = true
 rand = { workspace = true, features = ["small_rng"], optional = true }
 
 [build-dependencies]

diff --git a/integraal/src/lib.rs b/integraal/src/lib.rs
@@ -22,8 +22,10 @@
 
 mod parameters;
 mod structure;
+mod traits;
 
 // --- RE-EXPORTS
 
 pub use parameters::{ComputeMethod, DomainDescriptor, FunctionDescriptor};
 pub use structure::{Integraal, IntegraalError};
+pub use traits::Scalar;
diff --git a/integraal/src/parameters.rs b/integraal/src/parameters.rs
@@ -1,5 +1,7 @@
 //! integral parameterization code
 
+use crate::Scalar;
+
 /// Domain description enum
 ///
 /// This is essentially a discretization of the integrated space.
@@ -8,15 +10,15 @@
 /// `f64` values (i.e. the type used for further computations). In the future, adding support
 /// for higher dimension & generic value type can be considered.
 #[derive(Debug, Clone)]
-pub enum DomainDescriptor<'a> {
+pub enum DomainDescriptor<'a, T: Scalar> {
     /// List of values taken by the variable on which we integrate.
-    Explicit(&'a [f64]),
+    Explicit(&'a [T]),
     /// Description of a uniform discretization over a certain range of values.
     Uniform {
         /// First value of the range
-        start: f64,
+        start: T,
         /// Step between each value of the range
-        step: f64,
+        step: T,
         /// Total number of values
         n_step: usize,
     },
@@ -26,13 +28,16 @@ pub enum DomainDescriptor<'a> {
 ///
 /// This enum is used to provide either the values taken by the function or describe of to compute
 /// those.
-pub enum FunctionDescriptor {
+pub enum FunctionDescriptor<X>
+where
+    X: Scalar,
+{
     /// Direct expression of the function, taking a value of the domain as input & returning the
     /// image of that value through the function.
-    Closure(Box<dyn Fn(f64) -> f64>),
+    Closure(Box<dyn Fn(X) -> X>),
     /// List of values taken by the function. The coherence with the domain description must
     /// be ensured by the user in this case.
-    Values(Vec<f64>),
+    Values(Vec<X>),
 }
 
 /// Numerical integration method enum

diff --git a/integraal/src/structure/definitions.rs b/integraal/src/structure/definitions.rs
@@ -2,7 +2,7 @@
 
 // ------ IMPORTS
 
-use crate::{ComputeMethod, DomainDescriptor, FunctionDescriptor};
+use crate::{ComputeMethod, DomainDescriptor, FunctionDescriptor, Scalar};
 
 // ------ CONTENT
 
@@ -56,8 +56,8 @@ pub enum IntegraalError {
 /// # }
 /// ```
 #[derive(Default)]
-pub struct Integraal<'a> {
-    pub(crate) domain: Option<DomainDescriptor<'a>>,
-    pub(crate) function: Option<FunctionDescriptor>,
+pub struct Integraal<'a, X: Scalar> {
+    pub(crate) domain: Option<DomainDescriptor<'a, X>>,
+    pub(crate) function: Option<FunctionDescriptor<X>>,
     pub(crate) method: Option<ComputeMethod>,
 }
diff --git a/integraal/src/structure/implementations.rs b/integraal/src/structure/implementations.rs
@@ -2,19 +2,23 @@
 
 // ------ IMPORTS
 
-use crate::{ComputeMethod, DomainDescriptor, FunctionDescriptor, Integraal, IntegraalError};
+use crate::{
+    ComputeMethod, DomainDescriptor, FunctionDescriptor, Integraal, IntegraalError, Scalar,
+};
+use num::abs;
+use std::ops::Deref;
 
 // ------ CONTENT
 
-impl<'a> Integraal<'a> {
+impl<'a, X: Scalar> Integraal<'a, X> {
     /// Set the domain descriptor.
-    pub fn domain(&mut self, domain_descriptor: DomainDescriptor<'a>) -> &mut Self {
+    pub fn domain(&mut self, domain_descriptor: DomainDescriptor<'a, X>) -> &mut Self {
         self.domain = Some(domain_descriptor);
         self
     }
 
     /// Set the function descriptor.
-    pub fn function(&mut self, function_descriptor: FunctionDescriptor) -> &mut Self {
+    pub fn function(&mut self, function_descriptor: FunctionDescriptor<X>) -> &mut Self {
         self.function = Some(function_descriptor);
         self
     }
@@ -25,7 +29,11 @@ impl<'a> Integraal<'a> {
         self
     }
 
-    #[allow(clippy::missing_errors_doc, clippy::too_many_lines)]
+    #[allow(
+        clippy::missing_errors_doc,
+        clippy::missing_panics_doc,
+        clippy::too_many_lines
+    )]
     /// This method attempts to compute the integral. If it is successful, it will clear the
     /// internal [`FunctionDescriptor`] object before returning the result.
     ///
@@ -34,7 +42,7 @@ impl<'a> Integraal<'a> {
     /// This method returns a `Result` taking the following values:
     /// - `Ok(f64)` -- The computation was successfuly done
     /// - `Err(IntegraalError)` -- The computation failed for the reason specified by the enum.
-    pub fn compute(&mut self) -> Result<f64, IntegraalError> {
+    pub fn compute(&mut 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",
@@ -54,21 +62,21 @@ impl<'a> Integraal<'a> {
                     Some(ComputeMethod::RectangleLeft) => (1..n_sample)
                         .map(|idx| {
                             let step = args[idx] - args[idx - 1];
-                            step * vals[idx - 1]
+                            vals[idx - 1] * step
                         })
                         .sum(),
                     Some(ComputeMethod::RectangleRight) => (1..n_sample)
                         .map(|idx| {
                             let step = args[idx] - args[idx - 1];
-                            step * vals[idx]
+                            vals[idx] * step
                         })
                         .sum(),
                     Some(ComputeMethod::Trapezoid) => (1..n_sample)
                         .map(|idx| {
                             let step = args[idx] - args[idx - 1];
                             let y1 = vals[idx - 1];
                             let y2 = vals[idx];
-                            step * (y1.min(y2) + (y1 - y2).abs() / 2.0)
+                            (y1.min(y2) + abs(y1 - y2) / X::from_f32(2.0).unwrap()) * step
                         })
                         .sum(),
                     #[cfg(feature = "montecarlo")]
@@ -96,17 +104,17 @@ impl<'a> Integraal<'a> {
                 match &self.method {
                     Some(ComputeMethod::RectangleLeft) => {
                         // ignore the last value since its a left rule
-                        (0..*n_step - 1).map(|step_id| vals[step_id] * step).sum()
+                        (0..*n_step - 1).map(|step_id| vals[step_id] * *step).sum()
                     }
                     Some(ComputeMethod::RectangleRight) => {
                         // ignore the last value since its a left rule
-                        (1..*n_step).map(|step_id| vals[step_id] * step).sum()
+                        (1..*n_step).map(|step_id| vals[step_id] * *step).sum()
                     }
                     Some(ComputeMethod::Trapezoid) => (1..*n_step)
                         .map(|step_id| {
                             let y1 = vals[step_id - 1];
                             let y2 = vals[step_id];
-                            step * (y1.min(y2) + (y1 - y2).abs() / 2.0)
+                            (y1.min(y2) + (y1 - y2).abs() / X::from_f32(2.0).unwrap()) * *step
                         })
                         .sum(),
                     #[cfg(feature = "montecarlo")]
@@ -123,21 +131,21 @@ impl<'a> Integraal<'a> {
                 Some(ComputeMethod::RectangleLeft) => (1..args.len())
                     .map(|idx| {
                         let step = args[idx] - args[idx - 1];
-                        step * closure(args[idx - 1])
+                        closure(args[idx - 1]) * step
                     })
                     .sum(),
                 Some(ComputeMethod::RectangleRight) => (1..args.len())
                     .map(|idx| {
                         let step = args[idx] - args[idx - 1];
-                        step * closure(args[idx])
+                        closure(args[idx]) * step
                     })
                     .sum(),
                 Some(ComputeMethod::Trapezoid) => (1..args.len())
                     .map(|idx| {
                         let step = args[idx] - args[idx - 1];
-                        let y1 = closure(args[idx - 1]);
+                        let y1 = closure.deref()(args[idx - 1]);
                         let y2 = closure(args[idx]);
-                        step * (y1.min(y2) + (y1 - y2).abs() / 2.0)
+                        (y1.min(y2) + (y1 - y2).abs() / X::from_f32(2.0).unwrap()) * step
                     })
                     .sum(),
                 #[cfg(feature = "montecarlo")]
@@ -158,23 +166,23 @@ impl<'a> Integraal<'a> {
                 match &self.method {
                     Some(ComputeMethod::RectangleLeft) => (0..*n_step - 1)
                         .map(|step_id| {
-                            let x = start + step * step_id as f64;
-                            closure(x) * step
+                            let x = *start + *step * X::from_usize(step_id).unwrap();
+                            closure(x) * *step
                         })
                         .sum(),
                     Some(ComputeMethod::RectangleRight) => (1..*n_step)
                         .map(|step_id| {
-                            let x = start + step * step_id as f64;
-                            closure(x) * step
+                            let x = *start + *step * X::from_usize(step_id).unwrap();
+                            closure(x) * *step
                         })
                         .sum(),
                     Some(ComputeMethod::Trapezoid) => (1..*n_step)
                         .map(|step_id| {
-                            let x1 = start + step * (step_id - 1) as f64;
-                            let x2 = start + step * step_id as f64;
-                            let y1 = closure(x1);
+                            let x1 = *start + *step * X::from_usize(step_id - 1).unwrap();
+                            let x2 = *start + *step * X::from_usize(step_id).unwrap();
+                            let y1 = closure.deref()(x1);
                             let y2 = closure(x2);
-                            step * (y1.min(y2) + (y1 - y2).abs() / 2.0)
+                            (y1.min(y2) + (y1 - y2).abs() / X::from_f32(2.0).unwrap()) * *step
                         })
                         .sum(),
                     #[cfg(feature = "montecarlo")]

diff --git a/integraal/src/structure/tests.rs b/integraal/src/structure/tests.rs
@@ -20,15 +20,15 @@ macro_rules! almost_equal {
 
 macro_rules! generate_sample_descriptors {
     ($f: ident, $d: ident, $c: ident) => {
-        let $f = FunctionDescriptor::Closure(Box::new(|x| x));
-        let $d = DomainDescriptor::Explicit(&[]);
-        let $c = ComputeMethod::RectangleLeft;
+        let $f: FunctionDescriptor<f64> = FunctionDescriptor::Closure(Box::new(|x| x));
+        let $d: DomainDescriptor<'_, f64> = DomainDescriptor::Explicit(&[]);
+        let $c: ComputeMethod = ComputeMethod::RectangleLeft;
     };
 }
 
 macro_rules! generate_missing {
     ($a: ident, $b: ident) => {
-        let mut integral = Integraal::default();
+        let mut integral: Integraal<'_, f64> = Integraal::default();
         integral.$a($a).$b($b);
         assert_eq!(
             integral.compute(),
@@ -57,7 +57,7 @@ fn missing_parameters() {
     generate_missing!(function, domain);
 
     // missing all but one
-    let mut integral = Integraal::default();
+    let mut integral: Integraal<'_, f64> = Integraal::default();
     integral.method(method);
     assert_eq!(
         integral.compute(),

diff --git a/integraal/src/traits.rs b/integraal/src/traits.rs
@@ -0,0 +1,29 @@
+//! module doc
+
+/// Scalar value trait.
+///
+/// This trait is automatically implemented for all types implementing its requirements.
+pub trait Scalar:
+    Clone
+    + Copy
+    + num::Float
+    + num::Signed
+    + num::FromPrimitive
+    + std::ops::Sub<Output = Self>
+    + std::ops::Mul<Output = Self>
+    + std::iter::Sum
+{
+}
+
+impl<
+        X: Clone
+            + Copy
+            + num::Float
+            + num::Signed
+            + num::FromPrimitive
+            + std::ops::Sub<Output = Self>
+            + std::ops::Mul<Output = Self>
+            + std::iter::Sum,
+    > Scalar for X
+{
+}