diffsol_la/scalar/
mod.rs

1use std::{
2    fmt::Display,
3    ops::{Add, AddAssign, Mul, MulAssign, Sub, SubAssign},
4};
5
6#[cfg(feature = "cuda")]
7pub mod cuda;
8
9use crate::vector::VectorView;
10
11/// A scalar type suitable for numerical computations in ODE solvers.
12///
13/// This trait aggregates the crate-local numeric requirements shared across diffsol.
14///
15/// # Implementations
16/// DiffSol provides implementations for `f64` and `f32`.
17///
18/// # Examples
19/// ```
20/// use diffsol_la::Scalar;
21///
22/// fn compute<T: Scalar>(x: T, y: T) -> T {
23///     x * x + y
24/// }
25/// ```
26pub trait Scalar:
27    num_traits::Signed
28    + num_traits::Pow<Self, Output = Self>
29    + num_traits::Pow<i32, Output = Self>
30    + num_traits::FromPrimitive
31    + num_traits::ToPrimitive
32    + Display
33    + std::fmt::Debug
34    + Copy
35    + PartialEq
36    + PartialOrd
37    + AddAssign
38    + SubAssign
39    + MulAssign
40    + std::ops::DivAssign
41    + Send
42    + Sync
43    + 'static
44{
45    /// Machine epsilon for this scalar type (smallest representable positive value such that 1.0 + EPSILON != 1.0).
46    const EPSILON: Self;
47    /// Positive infinity value for this scalar type.
48    const INFINITY: Self;
49    /// Not-a-Number (NaN) value for this scalar type.
50    const NAN: Self;
51    /// Check if this value is NaN.
52    fn is_nan(self) -> bool;
53    /// Square root.
54    fn sqrt(self) -> Self;
55    /// Exponential.
56    fn exp(self) -> Self;
57    /// Sine.
58    fn sin(self) -> Self;
59    /// Cosine.
60    fn cos(self) -> Self;
61    /// Maximum of two values.
62    fn max(self, other: Self) -> Self;
63}
64
65/// A [`Scalar`] that also satisfies nalgebra's numeric field requirements.
66pub trait NalgebraScalar:
67    Scalar + nalgebra::Scalar + nalgebra::SimdRealField + nalgebra::ComplexField<RealField = Self>
68{
69}
70
71impl<T> NalgebraScalar for T where
72    T: Scalar
73        + nalgebra::Scalar
74        + nalgebra::SimdRealField
75        + nalgebra::ComplexField<RealField = Self>
76{
77}
78
79/// A [`Scalar`] that also satisfies faer's numeric field requirements.
80pub trait FaerScalar: Scalar + faer::traits::ComplexField + faer::traits::RealField {}
81
82impl<T> FaerScalar for T where T: Scalar + faer::traits::ComplexField + faer::traits::RealField {}
83
84/// The index type used throughout DiffSol for indexing vectors and matrices.
85pub type IndexType = usize;
86
87impl Scalar for f64 {
88    const EPSILON: Self = f64::EPSILON;
89    const INFINITY: Self = f64::INFINITY;
90    const NAN: Self = f64::NAN;
91    fn is_nan(self) -> bool {
92        self.is_nan()
93    }
94    fn sqrt(self) -> Self {
95        self.sqrt()
96    }
97    fn exp(self) -> Self {
98        self.exp()
99    }
100    fn sin(self) -> Self {
101        self.sin()
102    }
103    fn cos(self) -> Self {
104        self.cos()
105    }
106    fn max(self, other: Self) -> Self {
107        self.max(other)
108    }
109}
110
111impl Scalar for f32 {
112    const EPSILON: Self = f32::EPSILON;
113    const INFINITY: Self = f32::INFINITY;
114    const NAN: Self = f32::NAN;
115    fn is_nan(self) -> bool {
116        self.is_nan()
117    }
118    fn sqrt(self) -> Self {
119        self.sqrt()
120    }
121    fn exp(self) -> Self {
122        self.exp()
123    }
124    fn sin(self) -> Self {
125        self.sin()
126    }
127    fn cos(self) -> Self {
128        self.cos()
129    }
130    fn max(self, other: Self) -> Self {
131        self.max(other)
132    }
133}
134
135impl<T: FaerScalar> From<faer::Scale<T>> for Scale<T> {
136    fn from(s: faer::Scale<T>) -> Self {
137        Scale(s.0)
138    }
139}
140impl<T: FaerScalar> From<Scale<T>> for faer::Scale<T> {
141    fn from(s: Scale<T>) -> Self {
142        faer::Scale(s.value())
143    }
144}
145impl<T: Scalar> From<T> for Scale<T> {
146    fn from(s: T) -> Self {
147        Scale(s)
148    }
149}
150
151/// A wrapper for scalar values used when scaling vectors and matrices.
152#[derive(Copy, Clone, Debug)]
153pub struct Scale<E: Scalar>(pub E);
154
155impl<E: Scalar> Scale<E> {
156    /// Get the underlying scalar value.
157    #[inline]
158    pub fn value(self) -> E {
159        self.0
160    }
161}
162
163/// Create a `Scale` wrapper from a scalar value.
164///
165/// This is a convenience function equivalent to `Scale(value)`.
166#[inline]
167pub fn scale<E: Scalar>(value: E) -> Scale<E> {
168    Scale(value)
169}
170
171macro_rules! impl_bin_op {
172    ($trait:ident, $method:ident, $operator:tt) => {
173        impl<E: Scalar> $trait<Scale<E>> for Scale<E> {
174            type Output = Scale<E>;
175
176            #[inline]
177            fn $method(self, rhs: Scale<E>) -> Self::Output {
178                Scale(self.0 $operator rhs.0)
179            }
180        }
181    };
182}
183
184macro_rules! impl_assign_bin_op {
185    ($trait:ident, $method:ident, $operator:tt) => {
186        impl<E: Scalar> $trait<Scale<E>> for Scale<E> {
187            #[inline]
188            fn $method(&mut self, rhs: Scale<E>) {
189                self.0 = self.0 $operator rhs.0
190            }
191        }
192    };
193}
194
195impl_bin_op!(Mul, mul, *);
196impl_bin_op!(Add, add, +);
197impl_bin_op!(Sub, sub, -);
198
199impl_assign_bin_op!(MulAssign, mul_assign, *);
200impl_assign_bin_op!(AddAssign, add_assign, +);
201impl_assign_bin_op!(SubAssign, sub_assign, -);
202
203impl<E: Scalar> PartialEq for Scale<E> {
204    #[inline]
205    fn eq(&self, rhs: &Self) -> bool {
206        self.0 == rhs.0
207    }
208}
209
210impl<V: VectorView<'static, T = E>, E: Scalar> Mul<V> for Scale<E> {
211    type Output = V::Owned;
212    #[inline]
213    fn mul(self, rhs: V) -> Self::Output {
214        rhs * scale(self.0)
215    }
216}
217
218#[test]
219fn test_scale() {
220    assert_eq!(scale(2.0) * scale(3.0), scale(6.0));
221}
diffsol_la/scalar/mod.rs

diffsol_la/scalar/
mod.rs
