1use core::fmt;
2use core::ops::*;
3use core::write;
4
5#[allow(dead_code)]
6type Polar = crate::complex::polar::ComplexPolar<FT>;
7type Rectangular = crate::complex::rectangular::Complex<FT>;
8type FT = f64;
9
10impl Add<Rectangular> for FT {
11 type Output = Rectangular;
12 fn add(self, z: Self::Output) -> Self::Output {
13 z + self
14 }
15}
16
17impl Sub<Rectangular> for FT {
18 type Output = Rectangular;
19 fn sub(self, z: Self::Output) -> Self::Output {
20 Rectangular::new(self - z.re, -z.im)
21 }
22}
23
24impl Mul<Rectangular> for FT {
25 type Output = Rectangular;
26 fn mul(self, z: Self::Output) -> Self::Output {
27 z * self
28 }
29}
30
31impl Div<Rectangular> for FT {
32 type Output = Rectangular;
33 fn div(self, z: Self::Output) -> Self::Output {
34 self * z.recip()
35 }
36}
37
38impl Rectangular {
39 pub const NEG_ONE: Self = Self::new(-1.0, 0.0);
40 pub const NEG_I: Self = Self::new(0.0, -1.0);
41
42 #[cfg(feature = "glam")]
44 pub fn as_vec2(self) -> glam::Vec2 {
45 glam::vec2(self.re as f32, self.im as f32)
46 }
47
48 #[cfg(feature = "glam")]
50 pub fn as_dvec2(self) -> glam::DVec2 {
51 glam::dvec2(self.re, self.im)
52 }
53}
54
55impl fmt::Display for Rectangular {
56 fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
57 fn fmt_x(f: &mut fmt::Formatter, x: FT, sign_plus: bool) -> fmt::Result {
58 match (f.precision(), sign_plus) {
59 (None, false) => write!(f, "{}", x),
60 (None, true) => write!(f, "{:+}", x),
61 (Some(p), false) => write!(f, "{:.*}", p, x),
62 (Some(p), true) => write!(f, "{:+.*}", p, x),
63 }
64 }
65 match (self.re, self.im, f.sign_plus()) {
66 (re, 0.0, sp) => fmt_x(f, re, sp),
67 (0.0, 1.0, false) => write!(f, "i"),
68 (0.0, 1.0, true) => write!(f, "+i"),
69 (0.0, -1.0, _) => write!(f, "-i"),
70 (0.0, im, sp) => {
71 fmt_x(f, im, sp)?;
72 write!(f, "i")
73 }
74 (re, 1.0, sp) => {
75 fmt_x(f, re, sp)?;
76 write!(f, "+i")
77 }
78 (re, -1.0, sp) => {
79 fmt_x(f, re, sp)?;
80 write!(f, "-i")
81 }
82 (re, im, sp) => {
83 fmt_x(f, re, sp)?;
84 fmt_x(f, im, true)?;
85 write!(f, "i")
86 }
87 }
88 }
89}
90
91#[cfg(feature = "rand")]
92impl rand::distr::Distribution<Rectangular> for rand::distr::StandardUniform {
93 fn sample<R: rand::Rng + ?Sized>(&self, rng: &mut R) -> Rectangular {
94 rng.sample::<Polar, _>(self).to_rectangular()
95 }
96}
97
98#[cfg(feature = "glam")]
99impl From<glam::Vec2> for Rectangular {
100 fn from(v: glam::Vec2) -> Self {
101 Self::new(v.x as FT, v.y as FT)
102 }
103}
104
105#[cfg(feature = "glam")]
106impl From<glam::DVec2> for Rectangular {
107 fn from(v: glam::DVec2) -> Self {
108 Self::new(v.x, v.y)
109 }
110}
111
112#[cfg(feature = "glam")]
113impl From<Rectangular> for glam::Vec2 {
114 fn from(z: Rectangular) -> Self {
115 z.as_vec2()
116 }
117}
118
119#[cfg(feature = "glam")]
120impl From<Rectangular> for glam::DVec2 {
121 fn from(z: Rectangular) -> Self {
122 z.as_dvec2()
123 }
124}
125
126#[cfg(test)]
127mod tests {
128 use super::*;
129 use approx::*;
130 use core::f64::consts::{E, FRAC_PI_2, PI, SQRT_2};
131 use core::iter::Iterator;
132 use rand::{
133 Rng, SeedableRng,
134 distr::{Distribution, StandardUniform},
135 rngs::StdRng,
136 };
137
138 const NUM_SAMPLES: usize = 100;
139 const SEED: u64 = 21;
140
141 fn random_samples<T>() -> impl Iterator<Item = T>
142 where
143 StandardUniform: Distribution<T>,
144 {
145 StdRng::seed_from_u64(SEED)
146 .sample_iter(StandardUniform)
147 .take(NUM_SAMPLES)
148 }
149
150 #[test]
151 fn addition() {
152 for z0 in random_samples::<Rectangular>() {
153 for z1 in random_samples::<Rectangular>() {
154 let z = z0 + z1;
155 assert_eq!(z.re, z0.re + z1.re);
156 assert_eq!(z.im, z0.im + z1.im);
157
158 let z = z0 + z1.re;
159 assert_eq!(z.re, z0.re + z1.re);
160 assert_eq!(z.im, z0.im);
161
162 let z = z0.re + z1;
163 assert_eq!(z.re, z0.re + z1.re);
164 assert_eq!(z.im, z1.im);
165
166 let mut z = z0;
167 z += z1;
168 assert_eq!(z, z0 + z1);
169
170 let mut z = z0;
171 z += z1.re;
172 assert_eq!(z, z0 + z1.re);
173 }
174 assert_eq!(z0 + Rectangular::ZERO, z0);
175 }
176 }
177
178 #[test]
179 fn subtraction() {
180 for z0 in random_samples::<Rectangular>() {
181 for z1 in random_samples::<Rectangular>() {
182 let z = z0 - z1;
183 assert_eq!(z.re, z0.re - z1.re);
184 assert_eq!(z.im, z0.im - z1.im);
185
186 let z = z0 - z1.re;
187 assert_eq!(z.re, z0.re - z1.re);
188 assert_eq!(z.im, z0.im);
189
190 let z = z0.re - z1;
191 assert_eq!(z.re, z0.re - z1.re);
192 assert_eq!(z.im, -z1.im);
193
194 let mut z = z0;
195 z -= z1;
196 assert_eq!(z, z0 - z1);
197
198 let mut z = z0;
199 z -= z1.re;
200 assert_eq!(z, z0 - z1.re);
201 }
202 assert_eq!(z0 - z0, Rectangular::ZERO);
203 assert_eq!(z0 - Rectangular::ZERO, z0);
204 }
205 }
206
207 #[test]
208 fn multiplication() {
209 for z0 in random_samples::<Rectangular>() {
210 for z1 in random_samples::<Rectangular>() {
211 let z = z0 * z1;
212 assert_ulps_eq!(z.abs(), z0.abs() * z1.abs());
213 assert_ulps_eq!(
214 z.arg().sin(),
215 (z0.arg() + z1.arg()).sin(),
216 epsilon = 4.0 * FT::EPSILON
217 );
218
219 let z = z0 * z1.re;
220 assert_eq!(z.re, z0.re * z1.re);
221 assert_eq!(z.im, z0.im * z1.re);
222
223 let z = z0.re * z1;
224 assert_eq!(z.re, z0.re * z1.re);
225 assert_eq!(z.im, z0.re * z1.im);
226
227 let mut z = z0;
228 z *= z1;
229 assert_eq!(z, z0 * z1);
230
231 let mut z = z0;
232 z *= z1.re;
233 assert_eq!(z, z0 * z1.re);
234 }
235 assert_eq!(z0 * Rectangular::ONE, z0);
236 assert_eq!(z0 * Rectangular::ZERO, Rectangular::ZERO);
237 assert_eq!(z0 * 0.0, Rectangular::ZERO);
238 }
239 }
240
241 #[test]
242 fn division() {
243 for z0 in random_samples::<Rectangular>() {
244 for z1 in random_samples::<Rectangular>() {
245 let z = z0 / z1;
246 assert_relative_eq!(
247 z.abs(),
248 z0.abs() / z1.abs(),
249 max_relative = 3.0 * FT::EPSILON
250 );
251 assert_ulps_eq!(
252 z.arg().sin(),
253 (z0.arg() - z1.arg()).sin(),
254 epsilon = 4.0 * FT::EPSILON
255 );
256
257 let z = z0 / z1.re;
258 assert_eq!(z.re, z0.re / z1.re);
259 assert_eq!(z.im, z0.im / z1.re);
260
261 let z = z0.re / z1;
262 assert_ulps_eq!(z.abs(), z0.re.abs() / z1.abs());
263 assert_ulps_eq!(
264 z.arg().sin(),
265 (-z0.re.signum() * z1.arg()).sin(),
266 epsilon = 2.0 * FT::EPSILON
267 );
268
269 let mut z = z0;
270 z /= z1;
271 assert_eq!(z, z0 / z1);
272
273 let mut z = z0;
274 z /= z1.re;
275 assert_eq!(z, z0 / z1.re);
276 }
277 assert_ulps_eq!(z0 / z0, Rectangular::ONE);
278 assert_eq!(Rectangular::ZERO / z0, Rectangular::ZERO);
279 }
280 }
281
282 #[test]
283 fn negation() {
284 for z in random_samples::<Rectangular>() {
285 assert_eq!(-z, Rectangular::new(-z.re, -z.im));
286 }
287 assert_eq!(-Rectangular::ONE, Rectangular::NEG_ONE);
288 assert_eq!(-Rectangular::I, Rectangular::NEG_I);
289 assert_eq!(-Rectangular::NEG_ONE, Rectangular::ONE);
290 assert_eq!(-Rectangular::NEG_I, Rectangular::I);
291 }
292
293 #[test]
294 fn reciprocal() {
295 for z in random_samples::<Rectangular>() {
296 assert_eq!(z.recip(), 1.0 / z);
297 assert_ulps_eq!(z * z.recip(), Rectangular::ONE);
298 }
299 assert_eq!(Rectangular::ONE.recip(), Rectangular::ONE);
300 assert_eq!(Rectangular::I.recip(), Rectangular::NEG_I);
301 assert_eq!(Rectangular::NEG_ONE.recip(), Rectangular::NEG_ONE);
302 assert_eq!(Rectangular::NEG_I.recip(), Rectangular::I);
303 }
304
305 #[test]
306 fn sqrt() {
307 for z in random_samples::<Rectangular>() {
308 assert_ulps_eq!(z.sqrt().abs(), z.abs().sqrt());
309 assert_ulps_eq!(
310 z.sqrt().arg(),
311 z.arg() / 2.0,
312 epsilon = 1400.0 * FT::EPSILON
313 );
314 }
315 assert_eq!(Rectangular::ONE.sqrt(), Rectangular::ONE);
316 assert_eq!(Rectangular::NEG_ONE.sqrt(), Rectangular::I);
317 assert_eq!(
318 Rectangular::new(0.0, 2.0).sqrt(),
319 Rectangular::new(1.0, 1.0)
320 );
321 assert_eq!(
322 Rectangular::new(0.0, -2.0).sqrt(),
323 Rectangular::new(1.0, -1.0)
324 );
325 }
326
327 #[test]
328 fn abs() {
329 for z in random_samples::<Rectangular>() {
330 assert_ulps_eq!(z.abs_sq(), z.abs() * z.abs());
331 assert_eq!(z.abs_sq(), z.re * z.re + z.im * z.im);
332 }
333 assert_eq!(Rectangular::ONE.abs(), 1.0);
334 assert_eq!(Rectangular::I.abs(), 1.0);
335 assert_eq!(Rectangular::NEG_ONE.abs(), 1.0);
336 assert_eq!(Rectangular::NEG_I.abs(), 1.0);
337 assert_eq!(Rectangular::new(1.0, 1.0).abs(), SQRT_2);
338 assert_eq!(Rectangular::new(-1.0, 1.0).abs(), SQRT_2);
339 assert_eq!(Rectangular::new(-1.0, -1.0).abs(), SQRT_2);
340 assert_eq!(Rectangular::new(1.0, -1.0).abs(), SQRT_2);
341 }
342
343 #[test]
344 fn conjugate() {
345 for z in random_samples::<Rectangular>() {
346 assert_eq!(z.conjugate().re, z.re);
347 assert_eq!(z.conjugate().im, -z.im);
348 assert_eq!(z.conjugate().conjugate(), z);
349 }
350 assert_eq!(Rectangular::ONE.conjugate(), Rectangular::ONE);
351 assert_eq!(Rectangular::I.conjugate(), Rectangular::NEG_I);
352 assert_eq!(Rectangular::NEG_ONE.conjugate(), Rectangular::NEG_ONE);
353 assert_eq!(Rectangular::NEG_I.conjugate(), Rectangular::I);
354 }
355
356 #[test]
357 fn arg() {
358 assert_eq!(Rectangular::ONE.arg(), 0.0);
359 assert_eq!(Rectangular::I.arg(), FRAC_PI_2);
360 assert_eq!(Rectangular::NEG_ONE.arg(), PI);
361 assert_eq!(Rectangular::NEG_I.arg(), -FRAC_PI_2);
362 }
363
364 #[test]
365 fn exp() {
366 for z in random_samples::<Rectangular>() {
367 assert_eq!(z.exp().abs, z.re.exp());
368 assert_eq!(z.exp().arg, z.im);
369 assert_ulps_eq!(z.exp().ln(), z);
370 }
371 assert_eq!(Rectangular::ONE.exp(), Polar::new(E, 0.0));
372 assert_eq!(Rectangular::I.exp(), Polar::new(1.0, 1.0));
373 assert_eq!(Rectangular::NEG_ONE.exp(), Polar::new(E.recip(), 0.0));
374 assert_eq!(Rectangular::NEG_I.exp(), Polar::new(1.0, -1.0));
375 }
376
377 #[test]
378 fn log() {
379 for z in random_samples::<Rectangular>() {
380 assert_eq!(z.ln().re, z.abs().ln());
381 assert_eq!(z.ln().im, z.arg());
382 assert_ulps_eq!(z.ln().exp(), z.to_polar());
383 }
384 assert_eq!(Rectangular::ONE.ln(), Rectangular::ZERO);
385 assert_eq!(Rectangular::I.ln(), Rectangular::I * FRAC_PI_2);
386 assert_eq!(Rectangular::NEG_ONE.ln(), Rectangular::I * PI);
387 assert_eq!(Rectangular::NEG_I.ln(), Rectangular::I * -FRAC_PI_2);
388
389 assert_ulps_eq!(Rectangular::new(E, 0.0).ln(), Rectangular::ONE);
390 assert_ulps_eq!(Rectangular::new(2.0, 0.0).log2(), Rectangular::ONE);
391 assert_ulps_eq!(Rectangular::new(10.0, 0.0).log10(), Rectangular::ONE);
392 }
393
394 #[test]
395 fn powi() {
396 for z in random_samples::<Rectangular>() {
397 assert_eq!(z.powi(0), Polar::ONE);
398 assert_eq!(z.powi(1), z.to_polar());
399 for n in random_samples::<i32>() {
400 assert_eq!(z.powi(n).abs, z.abs().powi(n));
401 assert_eq!(z.powi(n).arg, z.arg() * n as FT);
402 }
403 }
404 for n in random_samples::<i32>() {
405 assert_eq!(Rectangular::ZERO.powi(n.abs()), Polar::ZERO);
406 assert_eq!(Rectangular::ONE.powi(n), Polar::ONE);
407 }
408 }
409
410 #[test]
411 fn powf() {
412 for z in random_samples::<Rectangular>() {
413 assert_eq!(z.powf(0.0), Polar::ONE);
414 assert_eq!(z.powf(1.0), z.to_polar());
415 for n in random_samples::<i32>() {
416 let x = n as FT * 0.01;
417 assert_eq!(z.powf(x).abs, z.abs().powf(x));
418 assert_eq!(z.powf(x).arg, z.arg() * x);
419 }
420 }
421 for n in random_samples::<i32>() {
422 let x = n as FT * 0.01;
423 assert_eq!(Rectangular::ZERO.powf(x.abs()), Polar::ZERO);
424 assert_eq!(Rectangular::ONE.powf(x), Polar::ONE);
425 }
426 }
427}