num_valid/kernels/

native64.rs

#![deny(rustdoc::broken_intra_doc_links)]

//! # Native 64-bit Floating-Point Kernel
//!
//! This module provides a numerical kernel implementation based on Rust's native
//! 64-bit floating-point numbers ([`f64`]) and complex numbers ([`num::Complex<f64>`]).
//! It serves as the standard, high-performance double-precision backend for the
//! [`num-valid`](crate) library.
//!
//! ## Purpose and Role
//!
//! The primary role of this module is to act as a **bridge** between the abstract traits
//! defined by [`num-valid`](crate) (like [`RealScalar`] and [`ComplexScalar`]) and Rust's
//! concrete, hardware-accelerated numeric types. It adapts [`f64`] and [`Complex<f64>`]
//! so they can be used seamlessly in generic code written against the library's traits.
//!
//! ### `f64` as a `RealScalar`
//!
//! By implementing [`RealScalar`] for [`f64`], this module makes Rust's native float
//! a "first-class citizen" in the [`num-valid`](crate)ecosystem. This implementation
//! provides concrete logic for all the operations required by [`RealScalar`] and its
//! super-traits (e.g., [`FpScalar`], [`Rounding`], [`Sign`], [`Constants`]).
//! Most of these implementations simply delegate to the highly optimized methods
//! available in Rust's standard library (e.g., `f64::sin`, `f64::sqrt`).
//!
//! This allows generic functions bounded by `T: RealScalar` to be instantiated with `f64`,
//! leveraging direct CPU floating-point support for maximum performance.
//!
//! ### `Complex<f64>` as a `ComplexScalar`
//!
//! Similarly, this module implements [`ComplexScalar`] for [`num::Complex<f64>`]. This
//! makes the standard complex number type from the `num` crate compatible with the
//! [`num-valid`](crate)abstraction for complex numbers. It implements all required traits,
//! including [`FpScalar`], [`ComplexScalarGetParts`], [`ComplexScalarSetParts`], and
//! [`ComplexScalarMutateParts`], by delegating to the methods of `num::Complex`.
//!
//! This enables generic code written against `T: ComplexScalar` to operate on
//! `Complex<f64>` values, providing a performant backend for complex arithmetic.
//!
//! ## Core Components
//!
//! - **[`Native64`]**: A marker struct that implements the [`NumKernel`](crate::NumKernel) trait. It is used
//!   in generic contexts to select this native kernel, specifying [`f64`] as its
//!   [`Real`](crate::NumKernel::Real) type and [`num::Complex<f64>`] as its
//!   [`Complex`](crate::NumKernel::Complex) type.
//! - **Trait Implementations**: This module is primarily composed of `impl` blocks that
//!   satisfy the contracts of core traits (`FpScalar`, `RealScalar`, `ComplexScalar`,
//!   `Arithmetic`, `MulAddRef`, etc.) for `f64` and `Complex<f64>`.
//!
//! ## Validation
//!
//! It is crucial to understand that the base types `f64` and `Complex<f64>` do not
//! inherently enforce properties like finiteness (they can be `NaN` or `Infinity`).
//!
//! The [`num-valid`](crate) library introduces validation at the **trait and wrapper level**,
//! not at the base type level. For instance, [`RealScalar::try_from_f64`] for `f64`
//! uses a [`StrictFinitePolicy`] to ensure the input is finite before creating an instance.
//! This design separates the raw, high-performance types from the validated, safer
//! abstractions built on top of them.
//!
//! For scenarios requiring types that are guaranteed to be valid by construction,
//! consider using the validated wrappers from the [`native64_validated`](crate::kernels::native64_validated) module.

use crate::{
    Arithmetic, ComplexScalar, Constants, FpScalar, MulAddRef, RealScalar,
    functions::{
        Clamp, Classify, ExpM1, Hypot, Ln1p, Rounding, Sign, TotalCmp,
        complex::{
            ComplexScalarConstructors, ComplexScalarGetParts, ComplexScalarMutateParts,
            ComplexScalarSetParts,
        },
    },
    kernels::RawKernel,
    scalar_kind,
    validation::{
        ErrorsTryFromf64, ErrorsValidationRawComplex, ErrorsValidationRawReal,
        Native64RawRealStrictFinitePolicy, StrictFinitePolicy, validate_complex,
    },
};
use duplicate::duplicate_item;
use num::Complex;
use num_traits::MulAddAssign;
use std::{cmp::Ordering, num::FpCategory};
use try_create::ValidationPolicy;

//----------------------------------------------------------------------------------------------
/// Numerical kernel specifier for Rust's native 64-bit floating-point types.
///
/// This struct is used as a generic argument or associated type to indicate
/// that the floating-point values for computations should be [`f64`] for real numbers
/// and [`num::Complex<f64>`] for complex numbers.
///
/// It implements the [`RawKernel`] trait, bridging the native raw types to the
/// [`num-valid`](crate)'s generic numerical abstraction.
#[derive(Debug, Clone, PartialEq, PartialOrd)]
pub struct Native64;

impl RawKernel for Native64 {
    /// The type for a real scalar value in this numerical kernel.
    type RawReal = f64;

    /// The type for a complex scalar value in this numerical kernel.
    type RawComplex = Complex<f64>;
}

//----------------------------------------------------------------------------------------------

//----------------------------------------------------------------------------------------------
/*
/// Implements the [`NumKernel`] trait for [`Native64`].
impl NumKernel for Native64 {
    /// The type for a real scalar value in this numerical kernel.
    type Real = f64;

    /// The type for a complex scalar value in this numerical kernel.
    type Complex = Complex<f64>;
}
*/
//----------------------------------------------------------------------------------------------

//----------------------------------------------------------------------------------------------
/// Marker implementation of [`Arithmetic`] for [`f64`].
///
/// This signifies that [`f64`] supports the standard arithmetic operations
/// (addition, subtraction, multiplication, division, remainder, negation)
/// as defined by the traits aggregated within [`Arithmetic`].
impl Arithmetic for f64 {}

/// Marker implementation of [`Arithmetic`] for [`num::Complex<f64>`].
///
/// This signifies that [`Complex<f64>`] supports the standard arithmetic operations
/// (addition, subtraction, multiplication, division, negation)
/// as defined by the traits aggregated within [`Arithmetic`].
impl Arithmetic for Complex<f64> {}
//----------------------------------------------------------------------------------------------

//----------------------------------------------------------------------------------------------
/// Implements the [`FpScalar`] trait for [`f64`].
///
/// This trait provides a common interface for floating-point scalar types within
/// the [`num-valid`](crate) ecosystem.
impl FpScalar for f64 {
    type Kind = scalar_kind::Real;

    /// Defines the associated real type, which is [`f64`] itself.
    type RealType = f64;

    /// Returns a reference to the underlying raw real value.
    fn as_raw_ref(&self) -> &Self::InnerType {
        self
    }
}

/// Implements the [`FpScalar`] trait for [`num::Complex<f64>`].
///
/// This trait provides a common interface for floating-point scalar types within
/// the [`num-valid`](crate) ecosystem.
impl FpScalar for Complex<f64> {
    type Kind = scalar_kind::Complex;

    /// Defines the associated real type for complex numbers, which is [`f64`].
    type RealType = f64;

    /// Returns a reference to the underlying raw complex value.
    fn as_raw_ref(&self) -> &Self::InnerType {
        self
    }
}
//----------------------------------------------------------------------------------------------

//----------------------------------------------------------------------------------------------
impl Sign for f64 {
    /// Returns a number with the magnitude of `self` and the sign of `sign`.
    #[inline(always)]
    fn kernel_copysign(self, sign: &Self) -> Self {
        self.copysign(*sign)
    }

    /// Returns `true` if the value is negative, −0 or NaN with a negative sign.
    #[inline(always)]
    fn kernel_is_sign_negative(&self) -> bool {
        self.is_sign_negative()
    }

    /// Returns `true` if the value is positive, +0 or NaN with a positive sign.
    #[inline(always)]
    fn kernel_is_sign_positive(&self) -> bool {
        self.is_sign_positive()
    }

    /// Computes the signum.
    ///
    /// The returned value is:
    /// - `1.0` if the value is positive, +0.0 or +∞
    /// - `−1.0` if the value is negative, −0.0 or −∞
    /// - `NaN` if the value is NaN
    #[inline(always)]
    fn kernel_signum(self) -> Self {
        self.signum()
    }
}

impl Rounding for f64 {
    /// Returns the smallest integer greater than or equal to `self`.
    #[inline(always)]
    fn kernel_ceil(self) -> Self {
        self.ceil()
    }

    /// Returns the largest integer smaller than or equal to `self`.
    #[inline(always)]
    fn kernel_floor(self) -> Self {
        self.floor()
    }

    /// Returns the fractional part of `self`.
    #[inline(always)]
    fn kernel_fract(self) -> Self {
        self.fract()
    }

    /// Rounds `self` to the nearest integer, rounding half-way cases away from zero.
    #[inline(always)]
    fn kernel_round(self) -> Self {
        self.round()
    }

    /// Returns the nearest integer to a number. Rounds half-way cases to the number with an even least significant digit.
    ///
    /// This function always returns the precise result.
    ///
    /// # Examples
    /// ```
    /// use num_valid::{RealScalar, functions::Rounding};
    ///
    /// let f = 3.3_f64;
    /// let g = -3.3_f64;
    /// let h = 3.5_f64;
    /// let i = 4.5_f64;
    ///
    /// assert_eq!(f.kernel_round_ties_even(), 3.0);
    /// assert_eq!(g.kernel_round_ties_even(), -3.0);
    /// assert_eq!(h.kernel_round_ties_even(), 4.0);
    /// assert_eq!(i.kernel_round_ties_even(), 4.0);
    /// ```
    #[inline(always)]
    fn kernel_round_ties_even(self) -> Self {
        self.round_ties_even()
    }

    /// Returns the integer part of `self`. This means that non-integer numbers are always truncated towards zero.
    ///    
    /// # Examples
    /// ```
    /// use num_valid::{RealScalar, functions::Rounding};
    ///
    /// let f = 3.7_f64;
    /// let g = 3.0_f64;
    /// let h = -3.7_f64;
    ///
    /// assert_eq!(f.kernel_trunc(), 3.0);
    /// assert_eq!(g.kernel_trunc(), 3.0);
    /// assert_eq!(h.kernel_trunc(), -3.0);
    /// ```
    #[inline(always)]
    fn kernel_trunc(self) -> Self {
        self.trunc()
    }
}

impl Constants for f64 {
    /// [Machine epsilon] value for `f64`.
    ///
    /// This is the difference between `1.0` and the next larger representable number.
    ///
    /// [Machine epsilon]: https://en.wikipedia.org/wiki/Machine_epsilon
    #[inline(always)]
    fn epsilon() -> Self {
        f64::EPSILON
    }

    /// Build and return the (floating point) value -1. represented by a `f64`.
    #[inline(always)]
    fn negative_one() -> Self {
        -1.
    }

    /// Build and return the (floating point) value 0.5 represented by the proper type.
    #[inline(always)]
    fn one_div_2() -> Self {
        0.5
    }

    /// Build and return the (floating point) value `2.0`.
    #[inline(always)]
    fn two() -> Self {
        2.
    }

    /// Build and return the maximum finite value allowed by the current floating point representation.
    #[inline(always)]
    fn max_finite() -> Self {
        f64::MAX
    }

    /// Build and return the minimum finite (i.e., the most negative) value allowed by the current floating point representation.
    #[inline(always)]
    fn min_finite() -> Self {
        f64::MIN
    }

    /// Build and return the (floating point) value `π`.
    #[inline(always)]
    fn pi() -> Self {
        std::f64::consts::PI
    }

    /// Build and return the (floating point) value `2 π`.
    #[inline(always)]
    fn two_pi() -> Self {
        std::f64::consts::PI * 2.
    }

    /// Build and return the (floating point) value `π/2`.
    #[inline(always)]
    fn pi_div_2() -> Self {
        std::f64::consts::FRAC_PI_2
    }

    /// Build and return the natural logarithm of 2, i.e. the (floating point) value `ln(2)`.
    #[inline(always)]
    fn ln_2() -> Self {
        std::f64::consts::LN_2
    }

    /// Build and return the natural logarithm of 10, i.e. the (floating point) value `ln(10)`.
    #[inline(always)]
    fn ln_10() -> Self {
        std::f64::consts::LN_10
    }

    /// Build and return the base-10 logarithm of 2, i.e. the (floating point) value `Log_10(2)`.
    #[inline(always)]
    fn log10_2() -> Self {
        std::f64::consts::LOG10_2
    }

    /// Build and return the base-2 logarithm of 10, i.e. the (floating point) value `Log_2(10)`.
    #[inline(always)]
    fn log2_10() -> Self {
        std::f64::consts::LOG2_10
    }

    /// Build and return the base-2 logarithm of `e`, i.e. the (floating point) value `Log_2(e)`.
    #[inline(always)]
    fn log2_e() -> Self {
        std::f64::consts::LOG2_E
    }

    /// Build and return the base-10 logarithm of `e`, i.e. the (floating point) value `Log_10(e)`.
    #[inline(always)]
    fn log10_e() -> Self {
        std::f64::consts::LOG10_E
    }

    /// Build and return the (floating point) value `e` represented by the proper type.
    #[inline(always)]
    fn e() -> Self {
        std::f64::consts::E
    }
}

impl Clamp for f64 {
    /// Clamp the value within the specified bounds.
    ///
    /// Returns `max` if `self` is greater than `max`, and `min` if `self` is less than `min`.
    /// Otherwise this returns `self`.
    ///
    /// Note that this function returns `NaN` if the initial value was `NaN` as well.
    ///
    /// # Panics
    /// Panics if `min` > `max`, `min` is `NaN`, or `max` is `NaN`.
    /// ```
    /// use num_valid::functions::Clamp;
    ///
    /// assert!(Clamp::clamp(-3.0f64, &-2.0, &1.0) == -2.0);
    /// assert!(Clamp::clamp(0.0f64, &-2.0, &1.0) == 0.0);
    /// assert!(Clamp::clamp(2.0f64, &-2.0, &1.0) == 1.0);
    /// assert!(Clamp::clamp(f64::NAN, &-2.0, &1.0).is_nan());
    /// ```
    #[inline(always)]
    fn clamp(self, min: &Self, max: &Self) -> Self {
        f64::clamp(self, *min, *max)
    }
}

impl Classify for f64 {
    /// Returns the floating point category of the number. If only one property is going to be tested,
    /// it is generally faster to use the specific predicate instead.
    /// ```
    /// use num_valid::functions::Classify;
    /// use std::num::FpCategory;
    ///
    /// let num = 12.4_f64;
    /// let inf = f64::INFINITY;
    ///
    /// assert_eq!(Classify::classify(&num), FpCategory::Normal);
    /// assert_eq!(Classify::classify(&inf), FpCategory::Infinite);
    /// ```
    #[inline(always)]
    fn classify(&self) -> FpCategory {
        f64::classify(*self)
    }
}

impl ExpM1 for f64 {
    /// Returns `e^(self) - 1`` in a way that is accurate even if the number is close to zero.
    #[inline(always)]
    fn exp_m1(self) -> Self {
        f64::exp_m1(self)
    }
}

impl Hypot for f64 {
    /// Compute the distance between the origin and a point (`self`, `other`) on the Euclidean plane.
    /// Equivalently, compute the length of the hypotenuse of a right-angle triangle with other sides having length `self.abs()` and `other.abs()`.
    #[inline(always)]
    fn hypot(self, other: &Self) -> Self {
        f64::hypot(self, *other)
    }
}

impl Ln1p for f64 {
    /// Returns `ln(1. + self)` (natural logarithm) more accurately than if the operations were performed separately.
    #[inline(always)]
    fn ln_1p(self) -> Self {
        f64::ln_1p(self)
    }
}

impl TotalCmp for f64 {
    #[inline(always)]
    fn total_cmp(&self, other: &Self) -> Ordering {
        self.total_cmp(other)
    }
}

impl RealScalar for f64 {
    type RawReal = f64;

    /// Multiplies two products and adds them in one fused operation, rounding to the nearest with only one rounding error.
    /// `a.kernel_mul_add_mul_mut(&b, &c, &d)` produces a result like `&a * &b + &c * &d`, but stores the result in `a` using its precision.
    #[inline(always)]
    fn kernel_mul_add_mul_mut(&mut self, mul: &Self, add_mul1: &Self, add_mul2: &Self) {
        self.mul_add_assign(*mul, add_mul1 * add_mul2);
    }

    /// Multiplies two products and subtracts them in one fused operation, rounding to the nearest with only one rounding error.
    /// `a.kernel_mul_sub_mul_mut(&b, &c, &d)` produces a result like `&a * &b - &c * &d`, but stores the result in `a` using its precision.
    #[inline(always)]
    fn kernel_mul_sub_mul_mut(&mut self, mul: &Self, sub_mul1: &Self, sub_mul2: &Self) {
        self.mul_add_assign(*mul, -sub_mul1 * sub_mul2);
    }

    /// Try to build a [`f64`] instance from a [`f64`].
    /// The returned value is `Ok` if the input `value` is finite,
    /// otherwise the returned value is `ErrorsTryFromf64`.
    #[inline(always)]
    fn try_from_f64(value: f64) -> Result<Self, ErrorsTryFromf64<f64>> {
        StrictFinitePolicy::<f64, 53>::validate(value)
            .map_err(|e| ErrorsTryFromf64::Output { source: e })
    }
}

impl ComplexScalarConstructors for Complex<f64> {
    type RawComplex = Complex<f64>;

    fn try_new_complex(
        real_part: f64,
        imag_part: f64,
    ) -> Result<Self, ErrorsValidationRawComplex<ErrorsValidationRawReal<f64>>> {
        validate_complex::<Native64RawRealStrictFinitePolicy>(&real_part, &imag_part)
            .map(|_| Complex::new(real_part, imag_part))
    }
}

/// Implement the [`ComplexScalarGetParts`] trait for the `Complex<f64>` type.
impl ComplexScalarGetParts for Complex<f64> {
    /// Get the real part of the complex number.
    #[inline(always)]
    fn real_part(&self) -> f64 {
        self.re
    }

    /// Get the imaginary part of the complex number.
    #[inline(always)]
    fn imag_part(&self) -> f64 {
        self.im
    }

    /// Returns a reference to the raw real part of the complex number.
    #[inline(always)]
    fn raw_real_part(&self) -> &f64 {
        &self.re
    }

    /// Returns a reference to the raw imaginary part of the complex number.
    #[inline(always)]
    fn raw_imag_part(&self) -> &f64 {
        &self.im
    }
}

/// Implement the [`ComplexScalarSetParts`] trait for the `Complex<f64>` type.
impl ComplexScalarSetParts for Complex<f64> {
    /// Set the value of the real part.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if `real_part` is not finite.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn set_real_part(&mut self, real_part: f64) {
        debug_assert!(
            real_part.is_finite(),
            "The real part is not finite (i.e. is infinite or NaN)."
        );
        self.re = real_part;
    }

    /// Set the value of the imaginary part.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if `imag_part` is not finite.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn set_imaginary_part(&mut self, imag_part: f64) {
        debug_assert!(
            imag_part.is_finite(),
            "The imaginary part is not finite (i.e. is infinite or NaN)."
        );
        self.im = imag_part;
    }
}

/// Implement the [`ComplexScalarMutateParts`] trait for the `Complex<f64>` type.
impl ComplexScalarMutateParts for Complex<f64> {
    /// Add the value of the the real coefficient `c` to real part of `self`.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if the real part of `self` is not finite after the addition.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn add_to_real_part(&mut self, c: &f64) {
        self.re += c;

        debug_assert!(
            self.re.is_finite(),
            "The real part is not finite (i.e. is infinite or NaN)."
        );
    }

    /// Add the value of the the real coefficient `c` to imaginary part of `self`.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if the imaginary part of `self` is not finite after the addition.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn add_to_imaginary_part(&mut self, c: &f64) {
        self.im += c;

        debug_assert!(
            self.im.is_finite(),
            "The imaginary part is not finite (i.e. is infinite or NaN)."
        );
    }

    /// Multiply the value of the real part by the real coefficient `c`.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if the real part of `self` is not finite after the multiplication.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn multiply_real_part(&mut self, c: &f64) {
        self.re *= c;

        debug_assert!(
            self.re.is_finite(),
            "The real part is not finite (i.e. is infinite or NaN)."
        );
    }

    /// Multiply the value of the imaginary part by the real coefficient `c`.
    ///
    /// # Panics
    /// In **debug builds**, this will panic if the imaginary part of `self` is not finite after the multiplication.
    /// This check is disabled in release builds for performance.
    #[inline(always)]
    fn multiply_imaginary_part(&mut self, c: &f64) {
        self.im *= c;

        debug_assert!(
            self.im.is_finite(),
            "The imaginary part is not finite (i.e. is infinite or NaN)."
        );
    }
}

/// Implement the [`ComplexScalar`] trait for the `Complex<f64>` type.
impl ComplexScalar for Complex<f64> {
    fn into_parts(self) -> (Self::RealType, Self::RealType) {
        (self.re, self.im)
    }
}

//----------------------------------------------------------------------------------------------

//----------------------------------------------------------------------------------------------
#[duplicate_item(
    T trait_comment;
    [f64] ["Implementation of the [`MulAddRef`] trait for `f64`."];
    [Complex<f64>] ["Implementation of the [`MulAddRef`] trait for `Complex<f64>`."];
)]
#[doc = trait_comment]
impl MulAddRef for T {
    /// Multiplies and adds in one fused operation, rounding to the nearest with only one rounding error.
    ///
    /// `a.mul_add(b, c)` produces a result like `a * &b + &c`.
    #[inline(always)]
    fn mul_add_ref(self, b: &Self, c: &Self) -> Self {
        <Self as num::traits::MulAdd>::mul_add(self, *b, *c)
    }
}
//----------------------------------------------------------------------------------------------

#[cfg(test)]
mod tests {
    use crate::{
        functions::TotalCmp,
        validation::{ErrorsValidationRawComplex, ErrorsValidationRawReal},
    };

    mod real {
        use super::*;
        use crate::Constants;

        #[test]
        fn test_constants() {
            assert_eq!(<f64 as Constants>::epsilon(), f64::EPSILON);
            assert_eq!(<f64 as Constants>::negative_one(), -1.0);
            assert_eq!(<f64 as Constants>::one_div_2(), 0.5);
            assert_eq!(<f64 as Constants>::two(), 2.0);
            assert_eq!(<f64 as Constants>::max_finite(), f64::MAX);
            assert_eq!(<f64 as Constants>::min_finite(), f64::MIN);
            assert_eq!(<f64 as Constants>::pi(), std::f64::consts::PI);
            assert_eq!(<f64 as Constants>::two_pi(), std::f64::consts::PI * 2.0);
            assert_eq!(<f64 as Constants>::pi_div_2(), std::f64::consts::FRAC_PI_2);
            assert_eq!(<f64 as Constants>::ln_2(), std::f64::consts::LN_2);
            assert_eq!(<f64 as Constants>::ln_10(), std::f64::consts::LN_10);
            assert_eq!(<f64 as Constants>::log10_2(), std::f64::consts::LOG10_2);
            assert_eq!(<f64 as Constants>::log2_10(), std::f64::consts::LOG2_10);
            assert_eq!(<f64 as Constants>::log2_e(), std::f64::consts::LOG2_E);
            assert_eq!(<f64 as Constants>::log10_e(), std::f64::consts::LOG10_E);
            assert_eq!(<f64 as Constants>::e(), std::f64::consts::E);
        }

        #[test]
        #[allow(clippy::op_ref)]
        fn multiply_ref() {
            let a = 2.0f64;
            let b = 3.0f64;
            let result = a * &b;
            assert_eq!(result, 6.0);
        }

        #[test]
        fn total_cmp() {
            let a = 2.0f64;
            let b = 3.0f64;
            assert_eq!(
                <f64 as TotalCmp>::total_cmp(&a, &b),
                std::cmp::Ordering::Less
            );
            assert_eq!(
                <f64 as TotalCmp>::total_cmp(&b, &a),
                std::cmp::Ordering::Greater
            );
            assert_eq!(
                <f64 as TotalCmp>::total_cmp(&a, &a),
                std::cmp::Ordering::Equal
            );
        }

        mod from_f64 {
            use crate::{RealNative64StrictFinite, RealScalar};

            #[test]
            fn test_from_f64_valid_constants() {
                // Test with mathematical constants (known valid values)
                let pi = RealNative64StrictFinite::from_f64(std::f64::consts::PI);
                assert_eq!(pi.as_ref(), &std::f64::consts::PI);

                let e = RealNative64StrictFinite::from_f64(std::f64::consts::E);
                assert_eq!(e.as_ref(), &std::f64::consts::E);

                let sqrt2 = RealNative64StrictFinite::from_f64(std::f64::consts::SQRT_2);
                assert_eq!(sqrt2.as_ref(), &std::f64::consts::SQRT_2);
            }

            #[test]
            fn test_from_f64_valid_values() {
                // Test with regular valid values
                let x = RealNative64StrictFinite::from_f64(42.0);
                assert_eq!(x.as_ref(), &42.0);

                let y = RealNative64StrictFinite::from_f64(-3.0);
                assert_eq!(y.as_ref(), &-3.0);

                let z = RealNative64StrictFinite::from_f64(0.0);
                assert_eq!(z.as_ref(), &0.0);
            }

            #[test]
            #[should_panic(expected = "RealScalar::from_f64() failed")]
            fn test_from_f64_nan_panics() {
                let _ = RealNative64StrictFinite::from_f64(f64::NAN);
            }

            #[test]
            #[should_panic(expected = "RealScalar::from_f64() failed")]
            fn test_from_f64_infinity_panics() {
                let _ = RealNative64StrictFinite::from_f64(f64::INFINITY);
            }

            #[test]
            #[should_panic(expected = "RealScalar::from_f64() failed")]
            fn test_from_f64_neg_infinity_panics() {
                let _ = RealNative64StrictFinite::from_f64(f64::NEG_INFINITY);
            }

            #[test]
            #[should_panic(expected = "RealScalar::from_f64() failed")]
            fn test_from_f64_subnormal_panics() {
                // Test with a subnormal value (very small but non-zero)
                let _ = RealNative64StrictFinite::from_f64(f64::MIN_POSITIVE / 2.0);
            }

            #[test]
            fn test_try_from_f64_error_handling() {
                // Verify try_from_f64 still works for error handling
                assert!(RealNative64StrictFinite::try_from_f64(f64::NAN).is_err());
                assert!(RealNative64StrictFinite::try_from_f64(f64::INFINITY).is_err());
                assert!(RealNative64StrictFinite::try_from_f64(f64::NEG_INFINITY).is_err());

                // Valid values work with try_from_f64
                assert!(RealNative64StrictFinite::try_from_f64(3.0).is_ok());
                assert!(RealNative64StrictFinite::try_from_f64(0.0).is_ok());
            }
        }

        mod truncate_to_usize {
            use crate::{kernels::RawRealTrait, validation::ErrorsRawRealToInteger};

            #[test]
            fn test_f64_truncate_to_usize_valid() {
                assert_eq!(42.0_f64.truncate_to_usize().unwrap(), 42);
                assert_eq!(42.9_f64.truncate_to_usize().unwrap(), 42);
                assert_eq!(0.0_f64.truncate_to_usize().unwrap(), 0);
            }

            #[test]
            fn test_f64_truncate_to_usize_not_finite() {
                assert!(matches!(
                    f64::NAN.truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::NotFinite { .. })
                ));
                assert!(matches!(
                    f64::INFINITY.truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::NotFinite { .. })
                ));
                assert!(matches!(
                    f64::NEG_INFINITY.truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::NotFinite { .. })
                ));
            }

            #[test]
            fn test_f64_truncate_to_usize_out_of_range() {
                // Negative value
                assert!(matches!(
                    (-1.0_f64).truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::OutOfRange { .. })
                ));
                // Value too large
                assert!(matches!(
                    ((usize::MAX as f64) + 1.0).truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::OutOfRange { .. })
                ));

                // this is exactly usize::MAX, which is out of range because f64 cannot represent all integers above 2^53 exactly
                assert!(matches!(
                    (usize::MAX as f64).truncate_to_usize(),
                    Err(ErrorsRawRealToInteger::OutOfRange { .. })
                ));
            }
        }
    }

    mod complex {
        use super::*;
        use crate::{
            ComplexScalarConstructors, ComplexScalarGetParts, ComplexScalarMutateParts,
            ComplexScalarSetParts,
        };
        use num::{Complex, Zero};

        #[test]
        fn real_part() {
            let c1 = Complex::new(1.23, 4.56);
            assert_eq!(c1.real_part(), 1.23);

            let c2 = Complex::new(-7.89, 0.12);
            assert_eq!(c2.real_part(), -7.89);

            let c3 = Complex::new(0.0, 10.0);
            assert_eq!(c3.real_part(), 0.0);

            let c_nan_re = Complex::new(f64::NAN, 5.0);
            assert!(c_nan_re.real_part().is_nan());

            let c_inf_re = Complex::new(f64::INFINITY, 5.0);
            assert!(c_inf_re.real_part().is_infinite());
            assert!(c_inf_re.real_part().is_sign_positive());

            let c_neg_inf_re = Complex::new(f64::NEG_INFINITY, 5.0);
            assert!(c_neg_inf_re.real_part().is_infinite());
            assert!(c_neg_inf_re.real_part().is_sign_negative());
        }

        #[test]
        fn imag_part() {
            let c1 = Complex::new(1.23, 4.56);
            assert_eq!(c1.imag_part(), 4.56);

            let c2 = Complex::new(7.89, -0.12);
            assert_eq!(c2.imag_part(), -0.12);

            let c3 = Complex::new(10.0, 0.0);
            assert_eq!(c3.imag_part(), 0.0);

            let c_nan_im = Complex::new(5.0, f64::NAN);
            assert!(c_nan_im.imag_part().is_nan());

            let c_inf_im = Complex::new(5.0, f64::INFINITY);
            assert!(c_inf_im.imag_part().is_infinite());
            assert!(c_inf_im.imag_part().is_sign_positive());

            let c_neg_inf_im = Complex::new(5.0, f64::NEG_INFINITY);
            assert!(c_neg_inf_im.imag_part().is_infinite());
            assert!(c_neg_inf_im.imag_part().is_sign_negative());
        }

        #[test]
        fn try_new_complex() {
            let r1 = 1.23;
            let i1 = 4.56;
            let c1 = Complex::<f64>::try_new_complex(r1, i1).unwrap();
            assert_eq!(c1.re, r1);
            assert_eq!(c1.im, i1);
            assert_eq!(c1.real_part(), r1);
            assert_eq!(c1.imag_part(), i1);

            let r2 = -7.89;
            let i2 = -0.12;
            let c2 = Complex::<f64>::try_new_complex(r2, i2).unwrap();
            assert_eq!(c2.re, r2);
            assert_eq!(c2.im, i2);
            assert_eq!(c2.real_part(), r2);
            assert_eq!(c2.imag_part(), i2);

            let r3 = 0.0;
            let i3 = 0.0;
            let c3 = Complex::<f64>::try_new_complex(r3, i3).unwrap();
            assert_eq!(c3.re, r3);
            assert_eq!(c3.im, i3);
            assert!(c3.is_zero()); // Assuming Zero trait is available and implemented

            let c_nan_re = Complex::<f64>::try_new_complex(f64::NAN, 5.0).unwrap_err();
            assert!(matches!(
                c_nan_re,
                ErrorsValidationRawComplex::InvalidRealPart { .. }
            ));

            let c_inf_im = Complex::<f64>::try_new_complex(10.0, f64::INFINITY).unwrap_err();
            assert!(matches!(
                c_inf_im,
                ErrorsValidationRawComplex::InvalidImaginaryPart { .. }
            ));

            let c_nan_re_inf_im =
                Complex::<f64>::try_new_complex(f64::NAN, f64::INFINITY).unwrap_err();
            assert!(matches!(
                c_nan_re_inf_im,
                ErrorsValidationRawComplex::InvalidBothParts { .. }
            ));
        }

        #[test]
        fn try_new_pure_real() {
            let r1 = 1.23;
            let c1 = Complex::<f64>::try_new_pure_real(r1).unwrap();
            assert_eq!(c1.re, r1);
            assert_eq!(c1.im, 0.0);

            let c_nan = Complex::<f64>::try_new_pure_real(f64::NAN).unwrap_err();
            assert!(matches!(
                c_nan,
                ErrorsValidationRawComplex::InvalidRealPart {
                    source: ErrorsValidationRawReal::IsNaN { .. }
                }
            ));
        }

        #[test]
        fn try_new_pure_imaginary() {
            let i1 = 1.23;
            let c1 = Complex::<f64>::try_new_pure_imaginary(i1).unwrap();
            assert_eq!(c1.re, 0.0);
            assert_eq!(c1.im, i1);

            let c_nan = Complex::<f64>::try_new_pure_imaginary(f64::NAN).unwrap_err();
            assert!(matches!(
                c_nan,
                ErrorsValidationRawComplex::InvalidImaginaryPart {
                    source: ErrorsValidationRawReal::IsNaN { .. }
                }
            ));
        }

        #[test]
        fn add_to_real_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.add_to_real_part(&3.0);
            assert_eq!(c.re, 4.0);
            assert_eq!(c.im, 2.0);

            c.add_to_real_part(&-5.0);
            assert_eq!(c.re, -1.0);
            assert_eq!(c.im, 2.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The real part is not finite (i.e. is infinite or NaN).")]
        fn add_to_real_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.add_to_real_part(&f64::NAN);
        }

        #[test]
        fn add_to_imaginary_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.add_to_imaginary_part(&3.0);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, 5.0);

            c.add_to_imaginary_part(&-4.0);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, 1.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The imaginary part is not finite (i.e. is infinite or NaN).")]
        fn add_to_imaginary_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.add_to_imaginary_part(&f64::NAN);
        }

        #[test]
        fn multiply_real_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.multiply_real_part(&3.0);
            assert_eq!(c.re, 3.0);
            assert_eq!(c.im, 2.0);

            c.multiply_real_part(&-2.0);
            assert_eq!(c.re, -6.0);
            assert_eq!(c.im, 2.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The real part is not finite (i.e. is infinite or NaN).")]
        fn multiply_real_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.multiply_real_part(&f64::NAN);
        }

        #[test]
        fn multiply_imaginary_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.multiply_imaginary_part(&3.0);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, 6.0);

            c.multiply_imaginary_part(&-0.5);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, -3.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The imaginary part is not finite (i.e. is infinite or NaN).")]
        fn multiply_imaginary_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.multiply_imaginary_part(&f64::NAN);
        }

        #[test]
        fn set_real_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.set_real_part(3.0);
            assert_eq!(c.re, 3.0);
            assert_eq!(c.im, 2.0);

            c.set_real_part(-4.0);
            assert_eq!(c.re, -4.0);
            assert_eq!(c.im, 2.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The real part is not finite (i.e. is infinite or NaN).")]
        fn set_real_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.set_real_part(f64::NAN);
        }

        #[test]
        fn set_imaginary_part() {
            let mut c = Complex::new(1.0, 2.0);
            c.set_imaginary_part(3.0);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, 3.0);

            c.set_imaginary_part(-4.0);
            assert_eq!(c.re, 1.0);
            assert_eq!(c.im, -4.0);
        }

        #[cfg(debug_assertions)]
        #[test]
        #[should_panic(expected = "The imaginary part is not finite (i.e. is infinite or NaN).")]
        fn set_imaginary_part_nan() {
            let mut c = Complex::new(1.0, 2.0);
            c.set_imaginary_part(f64::NAN);
        }

        #[test]
        #[allow(clippy::op_ref)]
        fn multiply_ref() {
            let c1 = Complex::new(1.0, 2.0);
            let c2 = Complex::new(3.0, 4.0);
            let result = c1 * &c2;
            assert_eq!(result, Complex::new(-5.0, 10.0)); // (1*3 - 2*4) + (1*4 + 2*3)i
        }
    }
}