embedded-complex-f32 0.1.0

// Copyright (C) 2026 Jorge Andre Castro
//
// Ce programme est un logiciel libre : vous pouvez le redistribuer et/ou le modifier
// selon les termes de la Licence Publique Générale GNU telle que publiée par la
// Free Software Foundation, soit la version 2 de la licence, soit (à votre gré)
// n'importe quelle version ultérieure.

//! # embedded-complex-f32
//!
//! Nombres complexes `f32` pour systèmes embarqués `no_std`.
//!
//! Zéro dépendance C · Sans `unsafe` · Racine carrée via `[`embedded-f32-sqrt`]`
//!
//! ## Fonctionnalités
//!
//! - Arithmétique complète : `+`, `-`, `*`, `/`
//! - Module (norme), carré de la norme, conjugué
//! - Puissance entière (`powi`)
//! - Racine carrée complexe (`csqrt`)
//! - Inverse et division vérifiée
//! - Gestion robuste de NaN / Infinity (IEEE 754)
//! - `Display` et `Debug` via `core::fmt`
//!
//! ## Hors périmètre
//!
//! `arg()`, `to_polar()`, `from_polar()`, `exp()`, `ln()` requièrent des
//! fonctions trigonométriques précises. Utilisez `libm` ou `micromath` selon
//! votre target et construisez par-dessus ce type.
//!
//! ## Exemple rapide
//!
//! ```rust
//! use embedded_complex_f32::Complex;
//!
//! let a = Complex::new(3.0, 4.0);
//! assert!((a.norm() - 5.0).abs() < 1e-4);  // |3 + 4i| = 5
//!
//! let b = Complex::new(1.0, 0.0);
//! let c = a + b;
//! assert_eq!(c.re(), 4.0);
//!
//! let conj = a.conj();
//! assert_eq!(conj.im(), -4.0);
//! ```

#![no_std]
#![forbid(unsafe_code)]
#![warn(missing_docs)]

use core::fmt;
use core::ops::{Add, AddAssign, Div, DivAssign, Mul, MulAssign, Neg, Sub, SubAssign};

use embedded_f32_sqrt::sqrt;

// Erreurs

/// Erreurs possibles dans les opérations sur les complexes.
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum ComplexError {
    /// Division par zéro (module nul du dénominateur).
    DivisionByZero,
    /// Argument impossible : entrée négative pour une racine réelle.
    NegativeInput,
    /// Valeur non définie (NaN propagé).
    Undefined,
}

impl fmt::Display for ComplexError {
    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
        match self {
            ComplexError::DivisionByZero => write!(f, "ComplexError: division by zero"),
            ComplexError::NegativeInput  => write!(f, "ComplexError: negative input"),
            ComplexError::Undefined      => write!(f, "ComplexError: undefined (NaN)"),
        }
    }
}

// Type principal

/// Nombre complexe en virgule flottante simple précision.
///
/// Représenté sous forme cartésienne `re + im·i`.
///
/// # Invariants
///
/// Aucun invariant fort n'est imposé : NaN et Infinity sont autorisés
/// (ils se propagent naturellement comme en IEEE 754).
#[derive(Clone, Copy, PartialEq)]
pub struct Complex {
    re: f32,
    im: f32,
}

// Constructeurs & accesseurs
impl Complex {
    /// Crée un nouveau nombre complexe `re + im·i`.
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let z = Complex::new(1.0, -2.0);
    /// assert_eq!(z.re(), 1.0);
    /// assert_eq!(z.im(), -2.0);
    /// ```
    #[inline]
    pub const fn new(re: f32, im: f32) -> Self {
        Self { re, im }
    }

    /// Le zéro complexe `0 + 0i`.
    pub const ZERO: Self = Self::new(0.0, 0.0);

    /// L'unité réelle `1 + 0i`.
    pub const ONE: Self = Self::new(1.0, 0.0);

    /// L'unité imaginaire `0 + 1i`.
    pub const I: Self = Self::new(0.0, 1.0);

    /// Partie réelle.
    #[inline]
    pub const fn re(self) -> f32 { self.re }

    /// Partie imaginaire.
    #[inline]
    pub const fn im(self) -> f32 { self.im }

    /// Retourne `true` si l'un des composants est NaN.
    #[inline]
    pub fn is_nan(self) -> bool {
        self.re.is_nan() || self.im.is_nan()
    }

    /// Retourne `true` si l'un des composants est infini.
    #[inline]
    pub fn is_infinite(self) -> bool {
        self.re.is_infinite() || self.im.is_infinite()
    }

    /// Retourne `true` si les deux composants sont finis.
    #[inline]
    pub fn is_finite(self) -> bool {
        self.re.is_finite() && self.im.is_finite()
    }
}


// Opérations algébriques

impl Complex {
    /// Conjugué : `re - im·i`.
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let z = Complex::new(3.0, -4.0);
    /// assert_eq!(z.conj(), Complex::new(3.0, 4.0));
    /// ```
    #[inline]
    pub fn conj(self) -> Self {
        Self::new(self.re, -self.im)
    }

    /// Norme (module) : `√(re² + im²)`.
    ///
    /// Utilise [`embedded_f32_sqrt::sqrt`] aucune dépendance C requise.
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let z = Complex::new(3.0, 4.0);
    /// assert!((z.norm() - 5.0).abs() < 1e-4);
    /// ```
    pub fn norm(self) -> f32 {
        sqrt(self.re * self.re + self.im * self.im).unwrap_or(f32::NAN)
    }

    /// Carré de la norme : `re² + im²` (évite la racine carrée).
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let z = Complex::new(3.0, 4.0);
    /// assert!((z.norm_sq() - 25.0).abs() < 1e-4);
    /// ```
    #[inline]
    pub fn norm_sq(self) -> f32 {
        self.re * self.re + self.im * self.im
    }

    /// Division vérifiée — retourne `Err(ComplexError::DivisionByZero)` si `|rhs| == 0`.
    ///
    /// ```rust
    /// use embedded_complex_f32::{Complex, ComplexError};
    /// let z = Complex::new(1.0, 0.0);
    /// assert_eq!(z.checked_div(Complex::ZERO), Err(ComplexError::DivisionByZero));
    /// ```
    pub fn checked_div(self, rhs: Self) -> Result<Self, ComplexError> {
        let denom = rhs.norm_sq();
        if denom == 0.0 {
            return Err(ComplexError::DivisionByZero);
        }
        Ok(Self::new(
            (self.re * rhs.re + self.im * rhs.im) / denom,
            (self.im * rhs.re - self.re * rhs.im) / denom,
        ))
    }

    /// Inverse : `1 / self`.
    ///
    /// Retourne `Err(ComplexError::DivisionByZero)` si `|self| == 0`.
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let r = Complex::new(4.0, 0.0).inv().unwrap();
    /// assert!((r.re() - 0.25).abs() < 1e-5);
    /// ```
    pub fn inv(self) -> Result<Self, ComplexError> {
        Self::ONE.checked_div(self)
    }

    /// Racine carrée complexe.
    ///
    /// Pour `z = r·e^{iθ}`, retourne `√r · e^{iθ/2}` via les identités
    /// trigonométriques demi-angle  sans libm, sans appel trig approché.
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// let z = Complex::new(3.0, 4.0);
    /// let s = z.csqrt().unwrap();
    /// let back = s * s;
    /// assert!((back.re() - 3.0).abs() < 1e-4);
    /// assert!((back.im() - 4.0).abs() < 1e-4);
    /// ```
    pub fn csqrt(self) -> Result<Self, ComplexError> {
        if self.is_nan() {
            return Err(ComplexError::Undefined);
        }
        let r = self.norm();
        let sqrt_r = sqrt(r).map_err(|_| ComplexError::NegativeInput)?;

        if sqrt_r == 0.0 {
            return Ok(Self::ZERO);
        }

        // cos(θ/2) = √((1 + cosθ)/2),  sin(θ/2) = sign(sinθ)·√((1 − cosθ)/2)
        let cos_theta = self.re / r;
        let half_cos = sqrt(((1.0 + cos_theta) / 2.0).max(0.0))
            .map_err(|_| ComplexError::NegativeInput)?;
        let half_sin_abs = sqrt(((1.0 - cos_theta) / 2.0).max(0.0))
            .map_err(|_| ComplexError::NegativeInput)?;
        let half_sin = if self.im < 0.0 { -half_sin_abs } else { half_sin_abs };

        Ok(Self::new(sqrt_r * half_cos, sqrt_r * half_sin))
    }

    /// Puissance entière `self^n` par exponentiation rapide.
    ///
    /// `n` peut être négatif (utilise [`inv`](Self::inv) si `n < 0`).
    ///
    /// ```rust
    /// use embedded_complex_f32::Complex;
    /// // i⁴ = 1
    /// let r = Complex::I.powi(4).unwrap();
    /// assert!((r.re() - 1.0).abs() < 1e-5);
    /// ```
    pub fn powi(self, n: i32) -> Result<Self, ComplexError> {
        let base = if n < 0 { self.inv()? } else { self };
        let mut exp = n.unsigned_abs();
        let mut result = Self::ONE;
        let mut b = base;
        while exp > 0 {
            if exp & 1 == 1 { result = result * b; }
            b = b * b;
            exp >>= 1;
        }
        Ok(result)
    }
}

// Opérateurs Rust (core::ops)

impl Add for Complex {
    type Output = Self;
    #[inline]
    fn add(self, rhs: Self) -> Self { Self::new(self.re + rhs.re, self.im + rhs.im) }
}
impl AddAssign for Complex {
    #[inline]
    fn add_assign(&mut self, rhs: Self) { *self = *self + rhs; }
}

impl Sub for Complex {
    type Output = Self;
    #[inline]
    fn sub(self, rhs: Self) -> Self { Self::new(self.re - rhs.re, self.im - rhs.im) }
}
impl SubAssign for Complex {
    #[inline]
    fn sub_assign(&mut self, rhs: Self) { *self = *self - rhs; }
}

/// `(a+bi)(c+di) = (ac−bd) + (ad+bc)i`
impl Mul for Complex {
    type Output = Self;
    #[inline]
    fn mul(self, rhs: Self) -> Self {
        Self::new(
            self.re * rhs.re - self.im * rhs.im,
            self.re * rhs.im + self.im * rhs.re,
        )
    }
}
impl MulAssign for Complex {
    #[inline]
    fn mul_assign(&mut self, rhs: Self) { *self = *self * rhs; }
}

/// Division via conjugué : `(a+bi)/(c+di) = [(a+bi)(c−di)] / (c²+d²)`.
///
/// Retourne `NaN + NaN·i` si `|rhs| == 0`. Préférez [`Complex::checked_div`]
/// pour une gestion explicite de l'erreur.
impl Div for Complex {
    type Output = Self;
    fn div(self, rhs: Self) -> Self {
        self.checked_div(rhs).unwrap_or(Self::new(f32::NAN, f32::NAN))
    }
}
impl DivAssign for Complex {
    #[inline]
    fn div_assign(&mut self, rhs: Self) { *self = *self / rhs; }
}

impl Neg for Complex {
    type Output = Self;
    #[inline]
    fn neg(self) -> Self { Self::new(-self.re, -self.im) }
}

// Opérateurs scalaires f32
impl Add<f32> for Complex {
    type Output = Self;
    #[inline]
    fn add(self, rhs: f32) -> Self { Self::new(self.re + rhs, self.im) }
}
impl Sub<f32> for Complex {
    type Output = Self;
    #[inline]
    fn sub(self, rhs: f32) -> Self { Self::new(self.re - rhs, self.im) }
}
impl Mul<f32> for Complex {
    type Output = Self;
    #[inline]
    fn mul(self, rhs: f32) -> Self { Self::new(self.re * rhs, self.im * rhs) }
}
impl Div<f32> for Complex {
    type Output = Self;
    #[inline]
    fn div(self, rhs: f32) -> Self { Self::new(self.re / rhs, self.im / rhs) }
}

// fmt

impl fmt::Display for Complex {
    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
        if self.im >= 0.0 || self.im.is_nan() {
            write!(f, "{} + {}i", self.re, self.im)
        } else {
            write!(f, "{} - {}i", self.re, -self.im)
        }
    }
}

impl fmt::Debug for Complex {
    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
        write!(f, "Complex {{ re: {}, im: {} }}", self.re, self.im)
    }
}

// Conversions

impl From<f32> for Complex {
    /// Embed un scalaire réel comme `x + 0i`.
    #[inline]
    fn from(x: f32) -> Self { Self::new(x, 0.0) }
}

impl From<(f32, f32)> for Complex {
    #[inline]
    fn from((re, im): (f32, f32)) -> Self { Self::new(re, im) }
}

impl From<Complex> for (f32, f32) {
    #[inline]
    fn from(z: Complex) -> Self { (z.re, z.im) }
}

// Tests
#[cfg(test)]
mod tests {
    use super::*;

    const EPS: f32 = 1e-4;

    fn approx_eq(a: f32, b: f32) -> bool { (a - b).abs() < EPS }
    fn complex_approx_eq(a: Complex, b: Complex) -> bool {
        approx_eq(a.re, b.re) && approx_eq(a.im, b.im)
    }

    //  Constantes 

    #[test]
    fn constants() {
        assert_eq!(Complex::ZERO, Complex::new(0.0, 0.0));
        assert_eq!(Complex::ONE,  Complex::new(1.0, 0.0));
        assert_eq!(Complex::I,    Complex::new(0.0, 1.0));
    }

    //  Arithmétique 

    #[test]
    fn addition() {
        assert_eq!(Complex::new(1.0, 2.0) + Complex::new(3.0, -1.0), Complex::new(4.0, 1.0));
    }

    #[test]
    fn subtraction() {
        assert_eq!(Complex::new(5.0, 3.0) - Complex::new(2.0, 1.0), Complex::new(3.0, 2.0));
    }

    #[test]
    fn multiplication() {
        // (1+i)(1−i) = 2
        let r = Complex::new(1.0, 1.0) * Complex::new(1.0, -1.0);
        assert!(approx_eq(r.re, 2.0));
        assert!(approx_eq(r.im, 0.0));
    }

    #[test]
    fn i_squared_is_minus_one() {
        let r = Complex::I * Complex::I;
        assert!(approx_eq(r.re, -1.0));
        assert!(approx_eq(r.im,  0.0));
    }

    #[test]
    fn division() {
        // (4+2i) / (1+i) = 3−i
        let r = Complex::new(4.0, 2.0) / Complex::new(1.0, 1.0);
        assert!(approx_eq(r.re,  3.0));
        assert!(approx_eq(r.im, -1.0));
    }

    #[test]
    fn division_by_zero_returns_nan() {
        let r = Complex::ONE / Complex::ZERO;
        assert!(r.is_nan());
    }

    #[test]
    fn checked_div_by_zero_returns_err() {
        assert_eq!(Complex::ONE.checked_div(Complex::ZERO), Err(ComplexError::DivisionByZero));
    }

    #[test]
    fn negation() {
        assert_eq!(-Complex::new(1.0, -2.0), Complex::new(-1.0, 2.0));
    }

    #[test]
    fn scalar_ops() {
        let z = Complex::new(3.0, 4.0);
        assert_eq!(z * 2.0, Complex::new(6.0, 8.0));
        assert_eq!(z / 2.0, Complex::new(1.5, 2.0));
        assert_eq!(z + 1.0, Complex::new(4.0, 4.0));
        assert_eq!(z - 1.0, Complex::new(2.0, 4.0));
    }

    #[test]
    fn assign_ops() {
        let mut z = Complex::new(1.0, 2.0);
        z += Complex::new(0.5, 0.5);
        assert!(approx_eq(z.re, 1.5));
        z *= Complex::new(2.0, 0.0);
        assert!(approx_eq(z.re, 3.0));
    }

    // Norme & conjugué 

    #[test]
    fn norm_pythagorean() {
        assert!(approx_eq(Complex::new(3.0,  4.0).norm(),  5.0));
        assert!(approx_eq(Complex::new(5.0, 12.0).norm(), 13.0));
    }

    #[test]
    fn norm_sq() {
        assert!(approx_eq(Complex::new(3.0, 4.0).norm_sq(), 25.0));
    }

    #[test]
    fn conjugate() {
        let z = Complex::new(3.0, -4.0);
        let c = z.conj();
        assert_eq!(c.im(), 4.0);
        let prod = z * c;
        assert!(approx_eq(prod.im, 0.0));
        assert!(prod.re > 0.0);
    }

    //  Racine carrée 

    #[test]
    fn csqrt_real_positive() {
        let r = Complex::new(9.0, 0.0).csqrt().unwrap();
        assert!(approx_eq(r.re, 3.0));
        assert!(approx_eq(r.im, 0.0));
    }

    #[test]
    fn csqrt_minus_one_gives_i() {
        let r = Complex::new(-1.0, 0.0).csqrt().unwrap();
        assert!(approx_eq(r.norm(), 1.0));
    }

    #[test]
    fn csqrt_general() {
        // s² == z est la vraie vérification, indépendante de toute trig
        let z = Complex::new(3.0, 4.0);
        let back = z.csqrt().unwrap() * z.csqrt().unwrap();
        assert!(approx_eq(back.re, 3.0));
        assert!(approx_eq(back.im, 4.0));
    }

    //Puissance entière 

    #[test]
    fn powi_zero_exp() {
        assert!(complex_approx_eq(Complex::new(5.0, 3.0).powi(0).unwrap(), Complex::ONE));
    }

    #[test]
    fn powi_i4_is_one() {
        let r = Complex::I.powi(4).unwrap();
        assert!(approx_eq(r.re, 1.0));
        assert!(approx_eq(r.im, 0.0));
    }

    #[test]
    fn powi_negative_exp() {
        let r = Complex::new(2.0, 0.0).powi(-1).unwrap();
        assert!(approx_eq(r.re, 0.5));
    }

    //  Inverse 

    #[test]
    fn inv_real() {
        let r = Complex::new(4.0, 0.0).inv().unwrap();
        assert!(approx_eq(r.re, 0.25));
    }

    #[test]
    fn inv_zero_returns_err() {
        assert_eq!(Complex::ZERO.inv(), Err(ComplexError::DivisionByZero));
    }

    //  Conversions 

    #[test]
    fn from_f32() {
        assert_eq!(Complex::from(3.0f32), Complex::new(3.0, 0.0));
    }

    #[test]
    fn from_tuple() {
        let z: Complex = (2.0f32, -1.0f32).into();
        assert_eq!(z, Complex::new(2.0, -1.0));
    }

    #[test]
    fn into_tuple() {
        let (re, im): (f32, f32) = Complex::new(7.0, -3.0).into();
        assert_eq!(re, 7.0);
        assert_eq!(im, -3.0);
    }

    // IEEE 754 

    #[test]
    fn nan_propagation() {
        let z = Complex::new(f32::NAN, 0.0);
        assert!(z.is_nan());
        assert!((z + Complex::ONE).is_nan());
    }

    #[test]
    fn finite_check() {
        assert!( Complex::new(1.0, 2.0).is_finite());
        assert!(!Complex::new(f32::INFINITY, 0.0).is_finite());
    }
}