qudit_core/
scalar.rs

1use crate::c32;
2use crate::c64;
3use crate::memory::Memorable;
4
5use faer_traits::ComplexField;
6use faer_traits::RealField;
7use num::FromPrimitive;
8use num::Signed;
9use num::ToPrimitive;
10use num::bigint::BigInt;
11use num::complex::ComplexFloat;
12use num::rational::Ratio;
13use num_traits::AsPrimitive;
14use num_traits::ConstOne;
15use num_traits::ConstZero;
16use num_traits::Float;
17use num_traits::FloatConst;
18use num_traits::NumAssign;
19use num_traits::NumAssignOps;
20use num_traits::NumAssignRef;
21use num_traits::NumOps;
22use rand_distr::{Distribution, StandardNormal, Uniform};
23
24use crate::bitwidth::BitWidthConvertible;
25
26/// A generic real number within the OpenQudit library.
27pub trait RealScalar:
28    Clone
29    + Copy
30    + Sync
31    + Send
32    + Sized
33    + std::any::Any
34    + std::fmt::Debug
35    + std::fmt::Display
36    + std::iter::Sum
37    + std::iter::Product
38    + Float
39    + FloatConst
40    + FromPrimitive
41    + ToPrimitive
42    + AsPrimitive<f32>
43    + AsPrimitive<f64>
44    + NumAssign
45    + NumAssignRef
46    + Default
47    + Signed
48    + 'static
49    + RealField
50    + ConstOne
51    + ConstZero
52    + Memorable
53    + BitWidthConvertible<Width32 = f32, Width64 = f64>
54{
55    /// The complex number type associated with this real number.
56    type C: ComplexScalar<R = Self>;
57
58    /// Convert this scalar to a rational number representation
59    fn into_ratio(self) -> Option<Ratio<BigInt>>;
60
61    /// Create a scalar from a rational number representation  
62    fn from_ratio(ratio: Ratio<BigInt>) -> Option<Self>;
63
64    /// Generate a random value from the standard normal distribution (mean=0, std=1)
65    fn standard_random() -> Self;
66
67    /// Generate a random value from a uniform distribution in the given range [min, max)
68    fn uniform_random(min: Self, max: Self) -> Self;
69
70    /// Check if two values are close using default tolerances
71    fn is_close(self, other: impl RealScalar) -> bool;
72
73    /// Check if two values are close with custom tolerances
74    /// Uses the formula: abs(a - b) <= (atol + rtol * abs(b))
75    fn is_close_with_tolerance(self, other: impl RealScalar, rtol: Self, atol: Self) -> bool;
76}
77
78impl RealScalar for f32 {
79    type C = c32;
80
81    fn into_ratio(self) -> Option<Ratio<BigInt>> {
82        Ratio::from_float(self as f64)
83    }
84
85    fn from_ratio(ratio: Ratio<BigInt>) -> Option<Self> {
86        ratio.to_f32()
87    }
88
89    fn standard_random() -> Self {
90        let mut rng = rand::rng();
91        StandardNormal.sample(&mut rng)
92    }
93
94    fn uniform_random(min: Self, max: Self) -> Self {
95        let mut rng = rand::rng();
96        Uniform::new(min, max).unwrap().sample(&mut rng)
97    }
98
99    fn is_close(self, other: impl RealScalar) -> bool {
100        self.is_close_with_tolerance(other, 1e-5, 1e-8)
101    }
102
103    fn is_close_with_tolerance(self, other: impl RealScalar, rtol: Self, atol: Self) -> bool {
104        (self - other.to32()).abs() <= (atol + rtol * other.to32().abs())
105    }
106}
107
108impl RealScalar for f64 {
109    type C = c64;
110
111    fn into_ratio(self) -> Option<Ratio<BigInt>> {
112        Ratio::from_float(self)
113    }
114
115    fn from_ratio(ratio: Ratio<BigInt>) -> Option<Self> {
116        ratio.to_f64()
117    }
118
119    fn standard_random() -> Self {
120        let mut rng = rand::rng();
121        StandardNormal.sample(&mut rng)
122    }
123
124    fn uniform_random(min: Self, max: Self) -> Self {
125        let mut rng = rand::rng();
126        Uniform::new(min, max).unwrap().sample(&mut rng)
127    }
128
129    fn is_close(self, other: impl RealScalar) -> bool {
130        self.is_close_with_tolerance(other, 1e-9, 1e-12)
131    }
132
133    fn is_close_with_tolerance(self, other: impl RealScalar, rtol: Self, atol: Self) -> bool {
134        (self - other.to64()).abs() <= (atol + rtol * other.to64().abs())
135    }
136}
137
138/// A generic complex number within the OpenQudit library.
139pub trait ComplexScalar:
140    Clone
141    + Copy
142    + Sync
143    + Send
144    + Sized
145    + std::any::Any
146    + std::fmt::Debug
147    + std::fmt::Display
148    + std::iter::Sum
149    + std::iter::Product
150    + NumAssign
151    + NumAssignRef
152    + NumOps<Self::R>
153    + NumAssignOps<Self::R>
154    + Default
155    + 'static
156    + ConstOne
157    + ConstZero
158    + ComplexField<Real = Self::R>
159    + ComplexFloat<Real = Self::R>
160    + Memorable
161    + BitWidthConvertible<Width32 = c32, Width64 = c64>
162{
163    /// The real number type associated with this complex number.
164    type R: RealScalar<C = Self>;
165
166    /// Create a complex number from real and imaginary parts
167    fn new(real: impl RealScalar, imag: impl RealScalar) -> Self;
168
169    /// Create a complex number from just the real part (imaginary = 0)
170    fn from_real(real: impl RealScalar) -> Self;
171
172    /// The real component of the complex number.
173    fn real(&self) -> Self::R;
174
175    /// The imaginary component of the complex number.
176    fn imag(&self) -> Self::R;
177
178    /// Generate a random complex number with both real and imaginary parts from standard normal
179    fn standard_random() -> Self;
180
181    /// Generate a random complex number with both parts uniform in given ranges
182    fn uniform_random(
183        real_min: Self::R,
184        real_max: Self::R,
185        imag_min: Self::R,
186        imag_max: Self::R,
187    ) -> Self;
188
189    /// Calculate the squared norm (|z|²) of the complex number
190    #[inline(always)]
191    fn norm_squared(self) -> Self::R {
192        self.re() * self.re() + self.im() * self.im()
193    }
194}
195
196impl ComplexScalar for c32 {
197    type R = f32;
198
199    fn new(real: impl RealScalar, imag: impl RealScalar) -> Self {
200        c32::new(real.to32(), imag.to32())
201    }
202
203    fn from_real(real: impl RealScalar) -> Self {
204        c32::new(real.to32(), 0.0)
205    }
206
207    #[inline(always)]
208    fn real(&self) -> Self::R {
209        self.re
210    }
211
212    #[inline(always)]
213    fn imag(&self) -> Self::R {
214        self.im
215    }
216
217    fn standard_random() -> Self {
218        c32::new(f32::standard_random(), f32::standard_random())
219    }
220
221    fn uniform_random(
222        real_min: Self::R,
223        real_max: Self::R,
224        imag_min: Self::R,
225        imag_max: Self::R,
226    ) -> Self {
227        c32::new(
228            f32::uniform_random(real_min, real_max),
229            f32::uniform_random(imag_min, imag_max),
230        )
231    }
232}
233
234impl ComplexScalar for c64 {
235    type R = f64;
236
237    fn new(real: impl RealScalar, imag: impl RealScalar) -> Self {
238        c64::new(real.to64(), imag.to64())
239    }
240
241    fn from_real(real: impl RealScalar) -> Self {
242        c64::new(real.to64(), 0.0)
243    }
244
245    #[inline(always)]
246    fn real(&self) -> Self::R {
247        self.re
248    }
249
250    #[inline(always)]
251    fn imag(&self) -> Self::R {
252        self.im
253    }
254
255    fn standard_random() -> Self {
256        c64::new(f64::standard_random(), f64::standard_random())
257    }
258
259    fn uniform_random(
260        real_min: Self::R,
261        real_max: Self::R,
262        imag_min: Self::R,
263        imag_max: Self::R,
264    ) -> Self {
265        c64::new(
266            f64::uniform_random(real_min, real_max),
267            f64::uniform_random(imag_min, imag_max),
268        )
269    }
270}