lib.rs source code [crates/num_complex/src/lib.rs]

1	// Copyright 2013 The Rust Project Developers. See the COPYRIGHT
2	// file at the top-level directory of this distribution and at
3	// http://rust-lang.org/COPYRIGHT.
4	//
5	// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
6	// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license
7	// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your
8	// option. This file may not be copied, modified, or distributed
9	// except according to those terms.
10
11	//! Complex numbers.
12	//!
13	//! ## Compatibility
14	//!
15	//! The `num-complex` crate is tested for rustc 1.60 and greater.
16
17	#![doc(html_root_url = "https://docs.rs/num-complex/0.4")]
18	#![no_std]
19
20	#[cfg(any(test, feature = "std"))]
21	#[cfg_attr(test, macro_use)]
22	extern crate std;
23
24	use core::fmt;
25	#[cfg(test)]
26	use core::hash;
27	use core::iter::{Product, Sum};
28	use core::ops::{Add, Div, Mul, Neg, Rem, Sub};
29	use core::str::FromStr;
30	#[cfg(feature = "std")]
31	use std::error::Error;
32
33	use num_traits::{ConstOne, ConstZero, Inv, MulAdd, Num, One, Pow, Signed, Zero};
34
35	use num_traits::float::FloatCore;
36	#[cfg(any(feature = "std", feature = "libm"))]
37	use num_traits::float::{Float, FloatConst};
38
39	mod cast;
40	mod pow;
41
42	#[cfg(any(feature = "std", feature = "libm"))]
43	mod complex_float;
44	#[cfg(any(feature = "std", feature = "libm"))]
45	pub use crate::complex_float::ComplexFloat;
46
47	#[cfg(feature = "rand")]
48	mod crand;
49	#[cfg(feature = "rand")]
50	pub use crate::crand::ComplexDistribution;
51
52	// FIXME #1284: handle complex NaN & infinity etc. This
53	// probably doesn't map to C's _Complex correctly.
54
55	/// A complex number in Cartesian form.
56	///
57	/// ## Representation and Foreign Function Interface Compatibility
58	///
59	/// `Complex<T>` is memory layout compatible with an array `[T; 2]`.
60	///
61	/// Note that `Complex<F>` where F is a floating point type is only* memory*
62	/// layout compatible with C's complex types, not* necessarily calling*
63	/// convention compatible. This means that for FFI you can only pass
64	/// `Complex<F>` behind a pointer, not as a value.
65	///
66	/// ## Examples
67	///
68	/// Example of extern function declaration.
69	///
70	/// ```
71	/// use num_complex::Complex;
72	/// use std::os::raw::c_int;
73	///
74	/// extern "C" {
75	/// fn zaxpy_(n: *const c_int, alpha: *const Complex<f64>,
76	/// x: *const Complex<f64>, incx: *const c_int,
77	/// y: *mut Complex<f64>, incy: *const c_int);
78	/// }
79	/// ```
80	#[derive(PartialEq, Eq, Copy, Clone, Hash, Debug, Default)]
81	#[repr(C)]
82	#[cfg_attr(
83	feature = "rkyv",
84	derive(rkyv::Archive, rkyv::Serialize, rkyv::Deserialize)
85	)]
86	#[cfg_attr(feature = "rkyv", archive(as = "Complex<T::Archived>"))]
87	#[cfg_attr(feature = "bytecheck", derive(bytecheck::CheckBytes))]
88	pub struct Complex<T> {
89	/// Real portion of the complex number
90	pub re: T,
91	/// Imaginary portion of the complex number
92	pub im: T,
93	}
94
95	/// Alias for a [`Complex<f32>`]
96	pub type Complex32 = Complex<f32>;
97
98	/// Create a new [`Complex<f32>`] with arguments that can convert [`Into<f32>`].
99	///
100	/// ```
101	/// use num_complex::{c32, Complex32};
102	/// assert_eq!(c32(`1u8`, `2`), Complex32::new(`1.0`, `2.0`));
103	/// ```
104	///
105	/// Note: ambiguous integer literals in Rust will [default] to `i32`, which does not* implement*
106	/// `Into<f32>`, so a call like `c32(1, 2)` will result in a type error. The example above uses a
107	/// suffixed `1u8` to set its type, and then the `2` can be inferred as the same type.
108	///
109	/// [default]: https://doc.rust-lang.org/reference/expressions/literal-expr.html#integer-literal-expressions
110	#[inline]
111	pub fn c32<T: Into<f32>>(re: T, im: T) -> Complex32 {
112	Complex::new(re.into(), im.into())
113	}
114
115	/// Alias for a [`Complex<f64>`]
116	pub type Complex64 = Complex<f64>;
117
118	/// Create a new [`Complex<f64>`] with arguments that can convert [`Into<f64>`].
119	///
120	/// ```
121	/// use num_complex::{c64, Complex64};
122	/// assert_eq!(c64(`1`, `2`), Complex64::new(`1.0`, `2.0`));
123	/// ```
124	#[inline]
125	pub fn c64<T: Into<f64>>(re: T, im: T) -> Complex64 {
126	Complex::new(re.into(), im.into())
127	}
128
129	impl<T> Complex<T> {
130	/// Create a new `Complex`
131	#[inline]
132	pub const fn new(re: T, im: T) -> Self {
133	Complex { re, im }
134	}
135	}
136
137	impl<T: Clone + Num> Complex<T> {
138	/// Returns the imaginary unit.
139	///
140	/// See also [`Complex::I`].
141	#[inline]
142	pub fn i() -> Self {
143	Self::new(T::zero(), T::one())
144	}
145
146	/// Returns the square of the norm (since `T` doesn't necessarily
147	/// have a sqrt function), i.e. `re^2 + im^2`.
148	#[inline]
149	pub fn norm_sqr(&self) -> T {
150	self.re.clone() * self.re.clone() + self.im.clone() * self.im.clone()
151	}
152
153	/// Multiplies `self` by the scalar `t`.
154	#[inline]
155	pub fn scale(&self, t: T) -> Self {
156	Self::new(self.re.clone() * t.clone(), self.im.clone() * t)
157	}
158
159	/// Divides `self` by the scalar `t`.
160	#[inline]
161	pub fn unscale(&self, t: T) -> Self {
162	Self::new(self.re.clone() / t.clone(), self.im.clone() / t)
163	}
164
165	/// Raises `self` to an unsigned integer power.
166	#[inline]
167	pub fn powu(&self, exp: u32) -> Self {
168	Pow::pow(self, exp)
169	}
170	}
171
172	impl<T: Clone + Num + Neg<Output = T>> Complex<T> {
173	/// Returns the complex conjugate. i.e. `re - i im`
174	#[inline]
175	pub fn conj(&self) -> Self {
176	Self::new(self.re.clone(), -self.im.clone())
177	}
178
179	/// Returns `1/self`
180	#[inline]
181	pub fn inv(&self) -> Self {
182	let norm_sqr: T = self.norm_sqr();
183	Self::new(
184	self.re.clone() / norm_sqr.clone(),
185	-self.im.clone() / norm_sqr,
186	)
187	}
188
189	/// Raises `self` to a signed integer power.
190	#[inline]
191	pub fn powi(&self, exp: i32) -> Self {
192	Pow::pow(self, rhs:exp)
193	}
194	}
195
196	impl<T: Clone + Signed> Complex<T> {
197	/// Returns the L1 norm `\|re\| + \|im\|` -- the [Manhattan distance] from the origin.
198	///
199	/// [Manhattan distance]: https://en.wikipedia.org/wiki/Taxicab_geometry
200	#[inline]
201	pub fn l1_norm(&self) -> T {
202	self.re.abs() + self.im.abs()
203	}
204	}
205
206	#[cfg(any(feature = "std", feature = "libm"))]
207	impl<T: Float> Complex<T> {
208	/// Create a new Complex with a given phase: `exp(i phase)`.*
209	/// See [cis (mathematics)](https://en.wikipedia.org/wiki/Cis_(mathematics)).
210	#[inline]
211	pub fn cis(phase: T) -> Self {
212	Self::new(phase.cos(), phase.sin())
213	}
214
215	/// Calculate \|self\|
216	#[inline]
217	pub fn norm(self) -> T {
218	self.re.hypot(self.im)
219	}
220	/// Calculate the principal Arg of self.
221	#[inline]
222	pub fn arg(self) -> T {
223	self.im.atan2(self.re)
224	}
225	/// Convert to polar form (r, theta), such that
226	/// `self = r exp(i * theta)`*
227	#[inline]
228	pub fn to_polar(self) -> (T, T) {
229	(self.norm(), self.arg())
230	}
231	/// Convert a polar representation into a complex number.
232	#[inline]
233	pub fn from_polar(r: T, theta: T) -> Self {
234	Self::new(r * theta.cos(), r * theta.sin())
235	}
236
237	/// Computes `e^(self)`, where `e` is the base of the natural logarithm.
238	#[inline]
239	pub fn exp(self) -> Self {
240	// formula: e^(a + bi) = e^a (cos(b) + isin(b)) = from_polar(e^a, b)*
241
242	let Complex { re, mut im } = self;
243	// Treat the corner cases +∞, -∞, and NaN
244	if re.is_infinite() {
245	if re < T::zero() {
246	if !im.is_finite() {
247	return Self::new(T::zero(), T::zero());
248	}
249	} else if im == T::zero() \|\| !im.is_finite() {
250	if im.is_infinite() {
251	im = T::nan();
252	}
253	return Self::new(re, im);
254	}
255	} else if re.is_nan() && im == T::zero() {
256	return self;
257	}
258
259	Self::from_polar(re.exp(), im)
260	}
261
262	/// Computes the principal value of natural logarithm of `self`.
263	///
264	/// This function has one branch cut:
265	///
266	/// `(-∞, 0]`, continuous from above.*
267	///
268	/// The branch satisfies `-π ≤ arg(ln(z)) ≤ π`.
269	#[inline]
270	pub fn ln(self) -> Self {
271	// formula: ln(z) = ln\|z\| + iarg(z)*
272	let (r, theta) = self.to_polar();
273	Self::new(r.ln(), theta)
274	}
275
276	/// Computes the principal value of the square root of `self`.
277	///
278	/// This function has one branch cut:
279	///
280	/// `(-∞, 0)`, continuous from above.*
281	///
282	/// The branch satisfies `-π/2 ≤ arg(sqrt(z)) ≤ π/2`.
283	#[inline]
284	pub fn sqrt(self) -> Self {
285	if self.im.is_zero() {
286	if self.re.is_sign_positive() {
287	// simple positive real √r, and copy `im` for its sign
288	Self::new(self.re.sqrt(), self.im)
289	} else {
290	// √(r e^(iπ)) = √r e^(iπ/2) = i√r
291	// √(r e^(-iπ)) = √r e^(-iπ/2) = -i√r
292	let re = T::zero();
293	let im = (-self.re).sqrt();
294	if self.im.is_sign_positive() {
295	Self::new(re, im)
296	} else {
297	Self::new(re, -im)
298	}
299	}
300	} else if self.re.is_zero() {
301	// √(r e^(iπ/2)) = √r e^(iπ/4) = √(r/2) + i√(r/2)
302	// √(r e^(-iπ/2)) = √r e^(-iπ/4) = √(r/2) - i√(r/2)
303	let one = T::one();
304	let two = one + one;
305	let x = (self.im.abs() / two).sqrt();
306	if self.im.is_sign_positive() {
307	Self::new(x, x)
308	} else {
309	Self::new(x, -x)
310	}
311	} else {
312	// formula: sqrt(r e^(it)) = sqrt(r) e^(it/2)
313	let one = T::one();
314	let two = one + one;
315	let (r, theta) = self.to_polar();
316	Self::from_polar(r.sqrt(), theta / two)
317	}
318	}
319
320	/// Computes the principal value of the cube root of `self`.
321	///
322	/// This function has one branch cut:
323	///
324	/// `(-∞, 0)`, continuous from above.*
325	///
326	/// The branch satisfies `-π/3 ≤ arg(cbrt(z)) ≤ π/3`.
327	///
328	/// Note that this does not match the usual result for the cube root of
329	/// negative real numbers. For example, the real cube root of `-8` is `-2`,
330	/// but the principal complex cube root of `-8` is `1 + i√3`.
331	#[inline]
332	pub fn cbrt(self) -> Self {
333	if self.im.is_zero() {
334	if self.re.is_sign_positive() {
335	// simple positive real ∛r, and copy `im` for its sign
336	Self::new(self.re.cbrt(), self.im)
337	} else {
338	// ∛(r e^(iπ)) = ∛r e^(iπ/3) = ∛r/2 + i∛r√3/2
339	// ∛(r e^(-iπ)) = ∛r e^(-iπ/3) = ∛r/2 - i∛r√3/2
340	let one = T::one();
341	let two = one + one;
342	let three = two + one;
343	let re = (-self.re).cbrt() / two;
344	let im = three.sqrt() * re;
345	if self.im.is_sign_positive() {
346	Self::new(re, im)
347	} else {
348	Self::new(re, -im)
349	}
350	}
351	} else if self.re.is_zero() {
352	// ∛(r e^(iπ/2)) = ∛r e^(iπ/6) = ∛r√3/2 + i∛r/2
353	// ∛(r e^(-iπ/2)) = ∛r e^(-iπ/6) = ∛r√3/2 - i∛r/2
354	let one = T::one();
355	let two = one + one;
356	let three = two + one;
357	let im = self.im.abs().cbrt() / two;
358	let re = three.sqrt() * im;
359	if self.im.is_sign_positive() {
360	Self::new(re, im)
361	} else {
362	Self::new(re, -im)
363	}
364	} else {
365	// formula: cbrt(r e^(it)) = cbrt(r) e^(it/3)
366	let one = T::one();
367	let three = one + one + one;
368	let (r, theta) = self.to_polar();
369	Self::from_polar(r.cbrt(), theta / three)
370	}
371	}
372
373	/// Raises `self` to a floating point power.
374	#[inline]
375	pub fn powf(self, exp: T) -> Self {
376	if exp.is_zero() {
377	return Self::one();
378	}
379	// formula: x^y = (ρ e^(i θ))^y = ρ^y e^(i θ y)
380	// = from_polar(ρ^y, θ y)
381	let (r, theta) = self.to_polar();
382	Self::from_polar(r.powf(exp), theta * exp)
383	}
384
385	/// Returns the logarithm of `self` with respect to an arbitrary base.
386	#[inline]
387	pub fn log(self, base: T) -> Self {
388	// formula: log_y(x) = log_y(ρ e^(i θ))
389	// = log_y(ρ) + log_y(e^(i θ)) = log_y(ρ) + ln(e^(i θ)) / ln(y)
390	// = log_y(ρ) + i θ / ln(y)
391	let (r, theta) = self.to_polar();
392	Self::new(r.log(base), theta / base.ln())
393	}
394
395	/// Raises `self` to a complex power.
396	#[inline]
397	pub fn powc(self, exp: Self) -> Self {
398	if exp.is_zero() {
399	return Self::one();
400	}
401	// formula: x^y = exp(y ln(x))*
402	(exp * self.ln()).exp()
403	}
404
405	/// Raises a floating point number to the complex power `self`.
406	#[inline]
407	pub fn expf(self, base: T) -> Self {
408	// formula: x^(a+bi) = x^a x^bi = x^a e^(b ln(x) i)
409	// = from_polar(x^a, b ln(x))
410	Self::from_polar(base.powf(self.re), self.im * base.ln())
411	}
412
413	/// Computes the sine of `self`.
414	#[inline]
415	pub fn sin(self) -> Self {
416	// formula: sin(a + bi) = sin(a)cosh(b) + icos(a)sinh(b)*
417	Self::new(
418	self.re.sin() * self.im.cosh(),
419	self.re.cos() * self.im.sinh(),
420	)
421	}
422
423	/// Computes the cosine of `self`.
424	#[inline]
425	pub fn cos(self) -> Self {
426	// formula: cos(a + bi) = cos(a)cosh(b) - isin(a)sinh(b)*
427	Self::new(
428	self.re.cos() * self.im.cosh(),
429	-self.re.sin() * self.im.sinh(),
430	)
431	}
432
433	/// Computes the tangent of `self`.
434	#[inline]
435	pub fn tan(self) -> Self {
436	// formula: tan(a + bi) = (sin(2a) + isinh(2b))/(cos(2a) + cosh(2b))*
437	let (two_re, two_im) = (self.re + self.re, self.im + self.im);
438	Self::new(two_re.sin(), two_im.sinh()).unscale(two_re.cos() + two_im.cosh())
439	}
440
441	/// Computes the principal value of the inverse sine of `self`.
442	///
443	/// This function has two branch cuts:
444	///
445	/// `(-∞, -1)`, continuous from above.*
446	/// `(1, ∞)`, continuous from below.*
447	///
448	/// The branch satisfies `-π/2 ≤ Re(asin(z)) ≤ π/2`.
449	#[inline]
450	pub fn asin(self) -> Self {
451	// formula: arcsin(z) = -i ln(sqrt(1-z^2) + iz)
452	let i = Self::i();
453	-i * ((Self::one() - self * self).sqrt() + i * self).ln()
454	}
455
456	/// Computes the principal value of the inverse cosine of `self`.
457	///
458	/// This function has two branch cuts:
459	///
460	/// `(-∞, -1)`, continuous from above.*
461	/// `(1, ∞)`, continuous from below.*
462	///
463	/// The branch satisfies `0 ≤ Re(acos(z)) ≤ π`.
464	#[inline]
465	pub fn acos(self) -> Self {
466	// formula: arccos(z) = -i ln(i sqrt(1-z^2) + z)
467	let i = Self::i();
468	-i * (i * (Self::one() - self * self).sqrt() + self).ln()
469	}
470
471	/// Computes the principal value of the inverse tangent of `self`.
472	///
473	/// This function has two branch cuts:
474	///
475	/// `(-∞i, -i]`, continuous from the left.*
476	/// `[i, ∞i)`, continuous from the right.*
477	///
478	/// The branch satisfies `-π/2 ≤ Re(atan(z)) ≤ π/2`.
479	#[inline]
480	pub fn atan(self) -> Self {
481	// formula: arctan(z) = (ln(1+iz) - ln(1-iz))/(2i)
482	let i = Self::i();
483	let one = Self::one();
484	let two = one + one;
485	if self == i {
486	return Self::new(T::zero(), T::infinity());
487	} else if self == -i {
488	return Self::new(T::zero(), -T::infinity());
489	}
490	((one + i * self).ln() - (one - i * self).ln()) / (two * i)
491	}
492
493	/// Computes the hyperbolic sine of `self`.
494	#[inline]
495	pub fn sinh(self) -> Self {
496	// formula: sinh(a + bi) = sinh(a)cos(b) + icosh(a)sin(b)*
497	Self::new(
498	self.re.sinh() * self.im.cos(),
499	self.re.cosh() * self.im.sin(),
500	)
501	}
502
503	/// Computes the hyperbolic cosine of `self`.
504	#[inline]
505	pub fn cosh(self) -> Self {
506	// formula: cosh(a + bi) = cosh(a)cos(b) + isinh(a)sin(b)*
507	Self::new(
508	self.re.cosh() * self.im.cos(),
509	self.re.sinh() * self.im.sin(),
510	)
511	}
512
513	/// Computes the hyperbolic tangent of `self`.
514	#[inline]
515	pub fn tanh(self) -> Self {
516	// formula: tanh(a + bi) = (sinh(2a) + isin(2b))/(cosh(2a) + cos(2b))*
517	let (two_re, two_im) = (self.re + self.re, self.im + self.im);
518	Self::new(two_re.sinh(), two_im.sin()).unscale(two_re.cosh() + two_im.cos())
519	}
520
521	/// Computes the principal value of inverse hyperbolic sine of `self`.
522	///
523	/// This function has two branch cuts:
524	///
525	/// `(-∞i, -i)`, continuous from the left.*
526	/// `(i, ∞i)`, continuous from the right.*
527	///
528	/// The branch satisfies `-π/2 ≤ Im(asinh(z)) ≤ π/2`.
529	#[inline]
530	pub fn asinh(self) -> Self {
531	// formula: arcsinh(z) = ln(z + sqrt(1+z^2))
532	let one = Self::one();
533	(self + (one + self * self).sqrt()).ln()
534	}
535
536	/// Computes the principal value of inverse hyperbolic cosine of `self`.
537	///
538	/// This function has one branch cut:
539	///
540	/// `(-∞, 1)`, continuous from above.*
541	///
542	/// The branch satisfies `-π ≤ Im(acosh(z)) ≤ π` and `0 ≤ Re(acosh(z)) < ∞`.
543	#[inline]
544	pub fn acosh(self) -> Self {
545	// formula: arccosh(z) = 2 ln(sqrt((z+1)/2) + sqrt((z-1)/2))
546	let one = Self::one();
547	let two = one + one;
548	two * (((self + one) / two).sqrt() + ((self - one) / two).sqrt()).ln()
549	}
550
551	/// Computes the principal value of inverse hyperbolic tangent of `self`.
552	///
553	/// This function has two branch cuts:
554	///
555	/// `(-∞, -1]`, continuous from above.*
556	/// `[1, ∞)`, continuous from below.*
557	///
558	/// The branch satisfies `-π/2 ≤ Im(atanh(z)) ≤ π/2`.
559	#[inline]
560	pub fn atanh(self) -> Self {
561	// formula: arctanh(z) = (ln(1+z) - ln(1-z))/2
562	let one = Self::one();
563	let two = one + one;
564	if self == one {
565	return Self::new(T::infinity(), T::zero());
566	} else if self == -one {
567	return Self::new(-T::infinity(), T::zero());
568	}
569	((one + self).ln() - (one - self).ln()) / two
570	}
571
572	/// Returns `1/self` using floating-point operations.
573	///
574	/// This may be more accurate than the generic `self.inv()` in cases
575	/// where `self.norm_sqr()` would overflow to ∞ or underflow to 0.
576	///
577	/// # Examples
578	///
579	/// ```
580	/// use num_complex::Complex64;
581	/// let c = Complex64::new(`1e300`, `1e300`);
582	///
583	/// // The generic `inv()` will overflow.
584	/// assert!(!c.inv().is_normal());
585	///
586	/// // But we can do better for `Float` types.
587	/// let inv = c.finv();
588	/// assert!(inv.is_normal());
589	/// println!("{:e}", inv);
590	///
591	/// let expected = Complex64::new(`5e-301`, `-5e-301`);
592	/// assert!((inv - expected).norm() < `1e-315`);
593	/// ```
594	#[inline]
595	pub fn finv(self) -> Complex<T> {
596	let norm = self.norm();
597	self.conj() / norm / norm
598	}
599
600	/// Returns `self/other` using floating-point operations.
601	///
602	/// This may be more accurate than the generic `Div` implementation in cases
603	/// where `other.norm_sqr()` would overflow to ∞ or underflow to 0.
604	///
605	/// # Examples
606	///
607	/// ```
608	/// use num_complex::Complex64;
609	/// let a = Complex64::new(`2.0`, `3.0`);
610	/// let b = Complex64::new(`1e300`, `1e300`);
611	///
612	/// // Generic division will overflow.
613	/// assert!(!(a / b).is_normal());
614	///
615	/// // But we can do better for `Float` types.
616	/// let quotient = a.fdiv(b);
617	/// assert!(quotient.is_normal());
618	/// println!("{:e}", quotient);
619	///
620	/// let expected = Complex64::new(`2.5e-300`, `5e-301`);
621	/// assert!((quotient - expected).norm() < `1e-315`);
622	/// ```
623	#[inline]
624	pub fn fdiv(self, other: Complex<T>) -> Complex<T> {
625	self * other.finv()
626	}
627	}
628
629	#[cfg(any(feature = "std", feature = "libm"))]
630	impl<T: Float + FloatConst> Complex<T> {
631	/// Computes `2^(self)`.
632	#[inline]
633	pub fn exp2(self) -> Self {
634	// formula: 2^(a + bi) = 2^a (cos(blog2) + isin(blog2))*
635	// = from_polar(2^a, blog2)*
636	Self::from_polar(self.re.exp2(), self.im * T::LN_2())
637	}
638
639	/// Computes the principal value of log base 2 of `self`.
640	#[inline]
641	pub fn log2(self) -> Self {
642	Self::ln(self) / T::LN_2()
643	}
644
645	/// Computes the principal value of log base 10 of `self`.
646	#[inline]
647	pub fn log10(self) -> Self {
648	Self::ln(self) / T::LN_10()
649	}
650	}
651
652	impl<T: FloatCore> Complex<T> {
653	/// Checks if the given complex number is NaN
654	#[inline]
655	pub fn is_nan(self) -> bool {
656	self.re.is_nan() \|\| self.im.is_nan()
657	}
658
659	/// Checks if the given complex number is infinite
660	#[inline]
661	pub fn is_infinite(self) -> bool {
662	!self.is_nan() && (self.re.is_infinite() \|\| self.im.is_infinite())
663	}
664
665	/// Checks if the given complex number is finite
666	#[inline]
667	pub fn is_finite(self) -> bool {
668	self.re.is_finite() && self.im.is_finite()
669	}
670
671	/// Checks if the given complex number is normal
672	#[inline]
673	pub fn is_normal(self) -> bool {
674	self.re.is_normal() && self.im.is_normal()
675	}
676	}
677
678	// Safety: `Complex<T>` is `repr(C)` and contains only instances of `T`, so we
679	// can guarantee it contains no added* padding. Thus, if `T: Zeroable`,*
680	// `Complex<T>` is also `Zeroable`
681	#[cfg(feature = "bytemuck")]
682	unsafe impl<T: bytemuck::Zeroable> bytemuck::Zeroable for Complex<T> {}
683
684	// Safety: `Complex<T>` is `repr(C)` and contains only instances of `T`, so we
685	// can guarantee it contains no added* padding. Thus, if `T: Pod`,*
686	// `Complex<T>` is also `Pod`
687	#[cfg(feature = "bytemuck")]
688	unsafe impl<T: bytemuck::Pod> bytemuck::Pod for Complex<T> {}
689
690	impl<T: Clone + Num> From<T> for Complex<T> {
691	#[inline]
692	fn from(re: T) -> Self {
693	Self::new(re, T::zero())
694	}
695	}
696
697	impl<'a, T: Clone + Num> From<&'a T> for Complex<T> {
698	#[inline]
699	fn from(re: &T) -> Self {
700	From::from(re.clone())
701	}
702	}
703
704	macro_rules! forward_ref_ref_binop {
705	(impl $imp:ident, $method:ident) => {
706	impl<'a, 'b, T: Clone + Num> $imp<&'b Complex<T>> for &'a Complex<T> {
707	type Output = Complex<T>;
708
709	#[inline]
710	fn $method(self, other: &Complex<T>) -> Self::Output {
711	self.clone().$method(other.clone())
712	}
713	}
714	};
715	}
716
717	macro_rules! forward_ref_val_binop {
718	(impl $imp:ident, $method:ident) => {
719	impl<'a, T: Clone + Num> $imp<Complex<T>> for &'a Complex<T> {
720	type Output = Complex<T>;
721
722	#[inline]
723	fn $method(self, other: Complex<T>) -> Self::Output {
724	self.clone().$method(other)
725	}
726	}
727	};
728	}
729
730	macro_rules! forward_val_ref_binop {
731	(impl $imp:ident, $method:ident) => {
732	impl<'a, T: Clone + Num> $imp<&'a Complex<T>> for Complex<T> {
733	type Output = Complex<T>;
734
735	#[inline]
736	fn $method(self, other: &Complex<T>) -> Self::Output {
737	self.$method(other.clone())
738	}
739	}
740	};
741	}
742
743	macro_rules! forward_all_binop {
744	(impl $imp:ident, $method:ident) => {
745	forward_ref_ref_binop!(impl $imp, $method);
746	forward_ref_val_binop!(impl $imp, $method);
747	forward_val_ref_binop!(impl $imp, $method);
748	};
749	}
750
751	// arithmetic
752	forward_all_binop!(impl Add, add);
753
754	// (a + i b) + (c + i d) == (a + c) + i (b + d)
755	impl<T: Clone + Num> Add<Complex<T>> for Complex<T> {
756	type Output = Self;
757
758	#[inline]
759	fn add(self, other: Self) -> Self::Output {
760	Self::Output::new(self.re + other.re, self.im + other.im)
761	}
762	}
763
764	forward_all_binop!(impl Sub, sub);
765
766	// (a + i b) - (c + i d) == (a - c) + i (b - d)
767	impl<T: Clone + Num> Sub<Complex<T>> for Complex<T> {
768	type Output = Self;
769
770	#[inline]
771	fn sub(self, other: Self) -> Self::Output {
772	Self::Output::new(self.re - other.re, self.im - other.im)
773	}
774	}
775
776	forward_all_binop!(impl Mul, mul);
777
778	// (a + i b) (c + i d) == (ac - bd) + i (ad + bc)*
779	impl<T: Clone + Num> Mul<Complex<T>> for Complex<T> {
780	type Output = Self;
781
782	#[inline]
783	fn mul(self, other: Self) -> Self::Output {
784	let re: T = self.re.clone() * other.re.clone() - self.im.clone() * other.im.clone();
785	let im: T = self.re * other.im + self.im * other.re;
786	Self::Output::new(re, im)
787	}
788	}
789
790	// (a + i b) (c + i d) + (e + i f) == ((ac + e) - bd) + i (ad + (bc + f))*
791	impl<T: Clone + Num + MulAdd<Output = T>> MulAdd<Complex<T>> for Complex<T> {
792	type Output = Complex<T>;
793
794	#[inline]
795	fn mul_add(self, other: Complex<T>, add: Complex<T>) -> Complex<T> {
796	let re: T = self.re.clone().mul_add(a:other.re.clone(), b:add.re)
797	- (self.im.clone() * other.im.clone()); // FIXME: use mulsub when available in rust
798	let im: T = self.re.mul_add(a:other.im, self.im.mul_add(a:other.re, b:add.im));
799	Complex::new(re, im)
800	}
801	}
802	impl<'a, 'b, T: Clone + Num + MulAdd<Output = T>> MulAdd<&'b Complex<T>> for &'a Complex<T> {
803	type Output = Complex<T>;
804
805	#[inline]
806	fn mul_add(self, other: &Complex<T>, add: &Complex<T>) -> Complex<T> {
807	self.clone().mul_add(a:other.clone(), b:add.clone())
808	}
809	}
810
811	forward_all_binop!(impl Div, div);
812
813	// (a + i b) / (c + i d) == [(a + i b) (c - i d)] / (cc + dd)*
814	// == [(ac + bd) / (cc + dd)] + i [(bc - ad) / (cc + dd)]
815	impl<T: Clone + Num> Div<Complex<T>> for Complex<T> {
816	type Output = Self;
817
818	#[inline]
819	fn div(self, other: Self) -> Self::Output {
820	let norm_sqr: T = other.norm_sqr();
821	let re: T = self.re.clone() * other.re.clone() + self.im.clone() * other.im.clone();
822	let im: T = self.im * other.re - self.re * other.im;
823	Self::Output::new(re:re / norm_sqr.clone(), im:im / norm_sqr)
824	}
825	}
826
827	forward_all_binop!(impl Rem, rem);
828
829	impl<T: Clone + Num> Complex<T> {
830	/// Find the gaussian integer corresponding to the true ratio rounded towards zero.
831	fn div_trunc(&self, divisor: &Self) -> Self {
832	let Complex { re: T, im: T } = self / divisor;
833	Complex::new(re:re.clone() - re % T::one(), im:im.clone() - im % T::one())
834	}
835	}
836
837	impl<T: Clone + Num> Rem<Complex<T>> for Complex<T> {
838	type Output = Self;
839
840	#[inline]
841	fn rem(self, modulus: Self) -> Self::Output {
842	let gaussian: Complex = self.div_trunc(&modulus);
843	self - modulus * gaussian
844	}
845	}
846
847	// Op Assign
848
849	mod opassign {
850	use core::ops::{AddAssign, DivAssign, MulAssign, RemAssign, SubAssign};
851
852	use num_traits::{MulAddAssign, NumAssign};
853
854	use crate::Complex;
855
856	impl<T: Clone + NumAssign> AddAssign for Complex<T> {
857	fn add_assign(&mut self, other: Self) {
858	self.re += other.re;
859	self.im += other.im;
860	}
861	}
862
863	impl<T: Clone + NumAssign> SubAssign for Complex<T> {
864	fn sub_assign(&mut self, other: Self) {
865	self.re -= other.re;
866	self.im -= other.im;
867	}
868	}
869
870	// (a + i b) (c + i d) == (ac - bd) + i (ad + bc)*
871	impl<T: Clone + NumAssign> MulAssign for Complex<T> {
872	fn mul_assign(&mut self, other: Self) {
873	let a = self.re.clone();
874
875	self.re *= other.re.clone();
876	self.re -= self.im.clone() * other.im.clone();
877
878	self.im *= other.re;
879	self.im += a * other.im;
880	}
881	}
882
883	// (a + i b) (c + i d) + (e + i f) == ((ac + e) - bd) + i (bc + (ad + f))*
884	impl<T: Clone + NumAssign + MulAddAssign> MulAddAssign for Complex<T> {
885	fn mul_add_assign(&mut self, other: Complex<T>, add: Complex<T>) {
886	let a = self.re.clone();
887
888	self.re.mul_add_assign(other.re.clone(), add.re); // (ac + e)*
889	self.re -= self.im.clone() * other.im.clone(); // ((ac + e) - bd)
890
891	let mut adf = a;
892	adf.mul_add_assign(other.im, add.im); // (ad + f)*
893	self.im.mul_add_assign(other.re, adf); // (bc + (ad + f))
894	}
895	}
896
897	impl<'a, 'b, T: Clone + NumAssign + MulAddAssign> MulAddAssign<&'a Complex<T>, &'b Complex<T>>
898	for Complex<T>
899	{
900	fn mul_add_assign(&mut self, other: &Complex<T>, add: &Complex<T>) {
901	self.mul_add_assign(other.clone(), add.clone());
902	}
903	}
904
905	// (a + i b) / (c + i d) == [(a + i b) (c - i d)] / (cc + dd)*
906	// == [(ac + bd) / (cc + dd)] + i [(bc - ad) / (cc + dd)]
907	impl<T: Clone + NumAssign> DivAssign for Complex<T> {
908	fn div_assign(&mut self, other: Self) {
909	let a = self.re.clone();
910	let norm_sqr = other.norm_sqr();
911
912	self.re *= other.re.clone();
913	self.re += self.im.clone() * other.im.clone();
914	self.re /= norm_sqr.clone();
915
916	self.im *= other.re;
917	self.im -= a * other.im;
918	self.im /= norm_sqr;
919	}
920	}
921
922	impl<T: Clone + NumAssign> RemAssign for Complex<T> {
923	fn rem_assign(&mut self, modulus: Self) {
924	let gaussian = self.div_trunc(&modulus);
925	self -= modulus gaussian;
926	}
927	}
928
929	impl<T: Clone + NumAssign> AddAssign<T> for Complex<T> {
930	fn add_assign(&mut self, other: T) {
931	self.re += other;
932	}
933	}
934
935	impl<T: Clone + NumAssign> SubAssign<T> for Complex<T> {
936	fn sub_assign(&mut self, other: T) {
937	self.re -= other;
938	}
939	}
940
941	impl<T: Clone + NumAssign> MulAssign<T> for Complex<T> {
942	fn mul_assign(&mut self, other: T) {
943	self.re *= other.clone();
944	self.im *= other;
945	}
946	}
947
948	impl<T: Clone + NumAssign> DivAssign<T> for Complex<T> {
949	fn div_assign(&mut self, other: T) {
950	self.re /= other.clone();
951	self.im /= other;
952	}
953	}
954
955	impl<T: Clone + NumAssign> RemAssign<T> for Complex<T> {
956	fn rem_assign(&mut self, other: T) {
957	self.re %= other.clone();
958	self.im %= other;
959	}
960	}
961
962	macro_rules! forward_op_assign {
963	(impl $imp:ident, $method:ident) => {
964	impl<'a, T: Clone + NumAssign> $imp<&'a Complex<T>> for Complex<T> {
965	#[inline]
966	fn $method(&mut self, other: &Self) {
967	self.$method(other.clone())
968	}
969	}
970	impl<'a, T: Clone + NumAssign> $imp<&'a T> for Complex<T> {
971	#[inline]
972	fn $method(&mut self, other: &T) {
973	self.$method(other.clone())
974	}
975	}
976	};
977	}
978
979	forward_op_assign!(impl AddAssign, add_assign);
980	forward_op_assign!(impl SubAssign, sub_assign);
981	forward_op_assign!(impl MulAssign, mul_assign);
982	forward_op_assign!(impl DivAssign, div_assign);
983	forward_op_assign!(impl RemAssign, rem_assign);
984	}
985
986	impl<T: Clone + Num + Neg<Output = T>> Neg for Complex<T> {
987	type Output = Self;
988
989	#[inline]
990	fn neg(self) -> Self::Output {
991	Self::Output::new(-self.re, -self.im)
992	}
993	}
994
995	impl<'a, T: Clone + Num + Neg<Output = T>> Neg for &'a Complex<T> {
996	type Output = Complex<T>;
997
998	#[inline]
999	fn neg(self) -> Self::Output {
1000	-self.clone()
1001	}
1002	}
1003
1004	impl<T: Clone + Num + Neg<Output = T>> Inv for Complex<T> {
1005	type Output = Self;
1006
1007	#[inline]
1008	fn inv(self) -> Self::Output {
1009	Complex::inv(&self)
1010	}
1011	}
1012
1013	impl<'a, T: Clone + Num + Neg<Output = T>> Inv for &'a Complex<T> {
1014	type Output = Complex<T>;
1015
1016	#[inline]
1017	fn inv(self) -> Self::Output {
1018	Complex::inv(self)
1019	}
1020	}
1021
1022	macro_rules! real_arithmetic {
1023	(@forward $imp:ident::$method:ident for $($real:ident),*) => (
1024	impl<'a, T: Clone + Num> $imp<&'a T> for Complex<T> {
1025	type Output = Complex<T>;
1026
1027	#[inline]
1028	fn $method(self, other: &T) -> Self::Output {
1029	self.$method(other.clone())
1030	}
1031	}
1032	impl<'a, T: Clone + Num> $imp<T> for &'a Complex<T> {
1033	type Output = Complex<T>;
1034
1035	#[inline]
1036	fn $method(self, other: T) -> Self::Output {
1037	self.clone().$method(other)
1038	}
1039	}
1040	impl<'a, 'b, T: Clone + Num> $imp<&'a T> for &'b Complex<T> {
1041	type Output = Complex<T>;
1042
1043	#[inline]
1044	fn $method(self, other: &T) -> Self::Output {
1045	self.clone().$method(other.clone())
1046	}
1047	}
1048	$(
1049	impl<'a> $imp<&'a Complex<$real>> for $real {
1050	type Output = Complex<$real>;
1051
1052	#[inline]
1053	fn $method(self, other: &Complex<$real>) -> Complex<$real> {
1054	self.$method(other.clone())
1055	}
1056	}
1057	impl<'a> $imp<Complex<$real>> for &'a $real {
1058	type Output = Complex<$real>;
1059
1060	#[inline]
1061	fn $method(self, other: Complex<$real>) -> Complex<$real> {
1062	self.clone().$method(other)
1063	}
1064	}
1065	impl<'a, 'b> $imp<&'a Complex<$real>> for &'b $real {
1066	type Output = Complex<$real>;
1067
1068	#[inline]
1069	fn $method(self, other: &Complex<$real>) -> Complex<$real> {
1070	self.clone().$method(other.clone())
1071	}
1072	}
1073	)*
1074	);
1075	($($real:ident),*) => (
1076	real_arithmetic!(@forward Add::add for $($real),*);
1077	real_arithmetic!(@forward Sub::sub for $($real),*);
1078	real_arithmetic!(@forward Mul::mul for $($real),*);
1079	real_arithmetic!(@forward Div::div for $($real),*);
1080	real_arithmetic!(@forward Rem::rem for $($real),*);
1081
1082	$(
1083	impl Add<Complex<$real>> for $real {
1084	type Output = Complex<$real>;
1085
1086	#[inline]
1087	fn add(self, other: Complex<$real>) -> Self::Output {
1088	Self::Output::new(self + other.re, other.im)
1089	}
1090	}
1091
1092	impl Sub<Complex<$real>> for $real {
1093	type Output = Complex<$real>;
1094
1095	#[inline]
1096	fn sub(self, other: Complex<$real>) -> Self::Output {
1097	Self::Output::new(self - other.re, $real::zero() - other.im)
1098	}
1099	}
1100
1101	impl Mul<Complex<$real>> for $real {
1102	type Output = Complex<$real>;
1103
1104	#[inline]
1105	fn mul(self, other: Complex<$real>) -> Self::Output {
1106	Self::Output::new(self * other.re, self * other.im)
1107	}
1108	}
1109
1110	impl Div<Complex<$real>> for $real {
1111	type Output = Complex<$real>;
1112
1113	#[inline]
1114	fn div(self, other: Complex<$real>) -> Self::Output {
1115	// a / (c + i d) == [a (c - i d)] / (cc + dd)*
1116	let norm_sqr = other.norm_sqr();
1117	Self::Output::new(self * other.re / norm_sqr.clone(),
1118	$real::zero() - self * other.im / norm_sqr)
1119	}
1120	}
1121
1122	impl Rem<Complex<$real>> for $real {
1123	type Output = Complex<$real>;
1124
1125	#[inline]
1126	fn rem(self, other: Complex<$real>) -> Self::Output {
1127	Self::Output::new(self, Self::zero()) % other
1128	}
1129	}
1130	)*
1131	);
1132	}
1133
1134	impl<T: Clone + Num> Add<T> for Complex<T> {
1135	type Output = Complex<T>;
1136
1137	#[inline]
1138	fn add(self, other: T) -> Self::Output {
1139	Self::Output::new(self.re + other, self.im)
1140	}
1141	}
1142
1143	impl<T: Clone + Num> Sub<T> for Complex<T> {
1144	type Output = Complex<T>;
1145
1146	#[inline]
1147	fn sub(self, other: T) -> Self::Output {
1148	Self::Output::new(self.re - other, self.im)
1149	}
1150	}
1151
1152	impl<T: Clone + Num> Mul<T> for Complex<T> {
1153	type Output = Complex<T>;
1154
1155	#[inline]
1156	fn mul(self, other: T) -> Self::Output {
1157	Self::Output::new(self.re * other.clone(), self.im * other)
1158	}
1159	}
1160
1161	impl<T: Clone + Num> Div<T> for Complex<T> {
1162	type Output = Self;
1163
1164	#[inline]
1165	fn div(self, other: T) -> Self::Output {
1166	Self::Output::new(self.re / other.clone(), self.im / other)
1167	}
1168	}
1169
1170	impl<T: Clone + Num> Rem<T> for Complex<T> {
1171	type Output = Complex<T>;
1172
1173	#[inline]
1174	fn rem(self, other: T) -> Self::Output {
1175	Self::Output::new(self.re % other.clone(), self.im % other)
1176	}
1177	}
1178
1179	real_arithmetic!(usize, u8, u16, u32, u64, u128, isize, i8, i16, i32, i64, i128, f32, f64);
1180
1181	// constants
1182	impl<T: ConstZero> Complex<T> {
1183	/// A constant `Complex` 0.
1184	pub const ZERO: Self = Self::new(T::ZERO, T::ZERO);
1185	}
1186
1187	impl<T: Clone + Num + ConstZero> ConstZero for Complex<T> {
1188	const ZERO: Self = Self::ZERO;
1189	}
1190
1191	impl<T: Clone + Num> Zero for Complex<T> {
1192	#[inline]
1193	fn zero() -> Self {
1194	Self::new(re:Zero::zero(), im:Zero::zero())
1195	}
1196
1197	#[inline]
1198	fn is_zero(&self) -> bool {
1199	self.re.is_zero() && self.im.is_zero()
1200	}
1201
1202	#[inline]
1203	fn set_zero(&mut self) {
1204	self.re.set_zero();
1205	self.im.set_zero();
1206	}
1207	}
1208
1209	impl<T: ConstOne + ConstZero> Complex<T> {
1210	/// A constant `Complex` 1.
1211	pub const ONE: Self = Self::new(T::ONE, T::ZERO);
1212
1213	/// A constant `Complex` _i_, the imaginary unit.
1214	pub const I: Self = Self::new(T::ZERO, T::ONE);
1215	}
1216
1217	impl<T: Clone + Num + ConstOne + ConstZero> ConstOne for Complex<T> {
1218	const ONE: Self = Self::ONE;
1219	}
1220
1221	impl<T: Clone + Num> One for Complex<T> {
1222	#[inline]
1223	fn one() -> Self {
1224	Self::new(re:One::one(), im:Zero::zero())
1225	}
1226
1227	#[inline]
1228	fn is_one(&self) -> bool {
1229	self.re.is_one() && self.im.is_zero()
1230	}
1231
1232	#[inline]
1233	fn set_one(&mut self) {
1234	self.re.set_one();
1235	self.im.set_zero();
1236	}
1237	}
1238
1239	macro_rules! write_complex {
1240	($f:ident, $t:expr, $prefix:expr, $re:expr, $im:expr, $T:ident) => {{
1241	let abs_re = if $re < Zero::zero() {
1242	$T::zero() - $re.clone()
1243	} else {
1244	$re.clone()
1245	};
1246	let abs_im = if $im < Zero::zero() {
1247	$T::zero() - $im.clone()
1248	} else {
1249	$im.clone()
1250	};
1251
1252	return if let Some(prec) = $f.precision() {
1253	fmt_re_im(
1254	$f,
1255	$re < $T::zero(),
1256	$im < $T::zero(),
1257	format_args!(concat!("{:.1$", $t, "}"), abs_re, prec),
1258	format_args!(concat!("{:.1$", $t, "}"), abs_im, prec),
1259	)
1260	} else {
1261	fmt_re_im(
1262	$f,
1263	$re < $T::zero(),
1264	$im < $T::zero(),
1265	format_args!(concat!("{:", $t, "}"), abs_re),
1266	format_args!(concat!("{:", $t, "}"), abs_im),
1267	)
1268	};
1269
1270	fn fmt_re_im(
1271	f: &mut fmt::Formatter<'_>,
1272	re_neg: bool,
1273	im_neg: bool,
1274	real: fmt::Arguments<'_>,
1275	imag: fmt::Arguments<'_>,
1276	) -> fmt::Result {
1277	let prefix = if f.alternate() { $prefix } else { "" };
1278	let sign = if re_neg {
1279	"-"
1280	} else if f.sign_plus() {
1281	"+"
1282	} else {
1283	""
1284	};
1285
1286	if im_neg {
1287	fmt_complex(
1288	f,
1289	format_args!(
1290	"{}{pre}{re}-{pre}{im}i",
1291	sign,
1292	re = real,
1293	im = imag,
1294	pre = prefix
1295	),
1296	)
1297	} else {
1298	fmt_complex(
1299	f,
1300	format_args!(
1301	"{}{pre}{re}+{pre}{im}i",
1302	sign,
1303	re = real,
1304	im = imag,
1305	pre = prefix
1306	),
1307	)
1308	}
1309	}
1310
1311	#[cfg(feature = "std")]
1312	// Currently, we can only apply width using an intermediate `String` (and thus `std`)
1313	fn fmt_complex(f: &mut fmt::Formatter<'_>, complex: fmt::Arguments<'_>) -> fmt::Result {
1314	use std::string::ToString;
1315	if let Some(width) = f.width() {
1316	write!(f, "{0: >1$}", complex.to_string(), width)
1317	} else {
1318	write!(f, "{}", complex)
1319	}
1320	}
1321
1322	#[cfg(not(feature = "std"))]
1323	fn fmt_complex(f: &mut fmt::Formatter<'_>, complex: fmt::Arguments<'_>) -> fmt::Result {
1324	write!(f, "{}", complex)
1325	}
1326	}};
1327	}
1328
1329	// string conversions
1330	impl<T> fmt::Display for Complex<T>
1331	where
1332	T: fmt::Display + Num + PartialOrd + Clone,
1333	{
1334	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1335	write_complex!(f, "", "", self.re, self.im, T)
1336	}
1337	}
1338
1339	impl<T> fmt::LowerExp for Complex<T>
1340	where
1341	T: fmt::LowerExp + Num + PartialOrd + Clone,
1342	{
1343	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1344	write_complex!(f, "e", "", self.re, self.im, T)
1345	}
1346	}
1347
1348	impl<T> fmt::UpperExp for Complex<T>
1349	where
1350	T: fmt::UpperExp + Num + PartialOrd + Clone,
1351	{
1352	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1353	write_complex!(f, "E", "", self.re, self.im, T)
1354	}
1355	}
1356
1357	impl<T> fmt::LowerHex for Complex<T>
1358	where
1359	T: fmt::LowerHex + Num + PartialOrd + Clone,
1360	{
1361	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1362	write_complex!(f, "x", "0x", self.re, self.im, T)
1363	}
1364	}
1365
1366	impl<T> fmt::UpperHex for Complex<T>
1367	where
1368	T: fmt::UpperHex + Num + PartialOrd + Clone,
1369	{
1370	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1371	write_complex!(f, "X", "0x", self.re, self.im, T)
1372	}
1373	}
1374
1375	impl<T> fmt::Octal for Complex<T>
1376	where
1377	T: fmt::Octal + Num + PartialOrd + Clone,
1378	{
1379	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1380	write_complex!(f, "o", "0o", self.re, self.im, T)
1381	}
1382	}
1383
1384	impl<T> fmt::Binary for Complex<T>
1385	where
1386	T: fmt::Binary + Num + PartialOrd + Clone,
1387	{
1388	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1389	write_complex!(f, "b", "0b", self.re, self.im, T)
1390	}
1391	}
1392
1393	fn from_str_generic<T, E, F>(s: &str, from: F) -> Result<Complex<T>, ParseComplexError<E>>
1394	where
1395	F: Fn(&str) -> Result<T, E>,
1396	T: Clone + Num,
1397	{
1398	let imag = match s.rfind('j') {
1399	None => 'i',
1400	_ => 'j',
1401	};
1402
1403	let mut neg_b = `false`;
1404	let mut a = s;
1405	let mut b = "";
1406
1407	for (i, w) in s.as_bytes().windows(`2`).enumerate() {
1408	let p = w[`0`];
1409	let c = w[`1`];
1410
1411	// ignore '+'/'-' if part of an exponent
1412	if (c == b'+' \|\| c == b'-') && !(p == b'e' \|\| p == b'E') {
1413	// trim whitespace around the separator
1414	a = s[..=i].trim_end_matches(char::is_whitespace);
1415	b = s[i + `2`..].trim_start_matches(char::is_whitespace);
1416	neg_b = c == b'-';
1417
1418	if b.is_empty() \|\| (neg_b && b.starts_with('-')) {
1419	return Err(ParseComplexError::expr_error());
1420	}
1421	break;
1422	}
1423	}
1424
1425	// split off real and imaginary parts
1426	if b.is_empty() {
1427	// input was either pure real or pure imaginary
1428	b = if a.ends_with(imag) { "0" } else { "0i" };
1429	}
1430
1431	let re;
1432	let neg_re;
1433	let im;
1434	let neg_im;
1435	if a.ends_with(imag) {
1436	im = a;
1437	neg_im = `false`;
1438	re = b;
1439	neg_re = neg_b;
1440	} else if b.ends_with(imag) {
1441	re = a;
1442	neg_re = `false`;
1443	im = b;
1444	neg_im = neg_b;
1445	} else {
1446	return Err(ParseComplexError::expr_error());
1447	}
1448
1449	// parse re
1450	let re = from(re).map_err(ParseComplexError::from_error)?;
1451	let re = if neg_re { T::zero() - re } else { re };
1452
1453	// pop imaginary unit off
1454	let mut im = &im[..im.len() - `1`];
1455	// handle im == "i" or im == "-i"
1456	if im.is_empty() \|\| im == "+" {
1457	im = "1";
1458	} else if im == "-" {
1459	im = "-1";
1460	}
1461
1462	// parse im
1463	let im = from(im).map_err(ParseComplexError::from_error)?;
1464	let im = if neg_im { T::zero() - im } else { im };
1465
1466	Ok(Complex::new(re, im))
1467	}
1468
1469	impl<T> FromStr for Complex<T>
1470	where
1471	T: FromStr + Num + Clone,
1472	{
1473	type Err = ParseComplexError<T::Err>;
1474
1475	/// Parses `a +/- bi`; `ai +/- b`; `a`; or `bi` where `a` and `b` are of type `T`
1476	fn from_str(s: &str) -> Result<Self, Self::Err> {
1477	from_str_generic(s, T::from_str)
1478	}
1479	}
1480
1481	impl<T: Num + Clone> Num for Complex<T> {
1482	type FromStrRadixErr = ParseComplexError<T::FromStrRadixErr>;
1483
1484	/// Parses `a +/- bi`; `ai +/- b`; `a`; or `bi` where `a` and `b` are of type `T`
1485	///
1486	/// `radix` must be <= 18; larger radix would include i* and j as digits,*
1487	/// which cannot be supported.
1488	///
1489	/// The conversion returns an error if 18 <= radix <= 36; it panics if radix > 36.
1490	///
1491	/// The elements of `T` are parsed using `Num::from_str_radix` too, and errors
1492	/// (or panics) from that are reflected here as well.
1493	fn from_str_radix(s: &str, radix: u32) -> Result<Self, Self::FromStrRadixErr> {
1494	assert!(
1495	radix <= `36`,
1496	"from_str_radix: radix is too high (maximum 36)"
1497	);
1498
1499	// larger radix would include 'i' and 'j' as digits, which cannot be supported
1500	if radix > `18` {
1501	return Err(ParseComplexError::unsupported_radix());
1502	}
1503
1504	from_str_generic(s, \|x\| -> Result<T, T::FromStrRadixErr> {
1505	T::from_str_radix(x, radix)
1506	})
1507	}
1508	}
1509
1510	impl<T: Num + Clone> Sum for Complex<T> {
1511	fn sum<I>(iter: I) -> Self
1512	where
1513	I: Iterator<Item = Self>,
1514	{
1515	iter.fold(Self::zero(), \|acc: Complex, c: Complex\| acc + c)
1516	}
1517	}
1518
1519	impl<'a, T: 'a + Num + Clone> Sum<&'a Complex<T>> for Complex<T> {
1520	fn sum<I>(iter: I) -> Self
1521	where
1522	I: Iterator<Item = &'a Complex<T>>,
1523	{
1524	iter.fold(Self::zero(), \|acc: Complex, c: &'a Complex\| acc + c)
1525	}
1526	}
1527
1528	impl<T: Num + Clone> Product for Complex<T> {
1529	fn product<I>(iter: I) -> Self
1530	where
1531	I: Iterator<Item = Self>,
1532	{
1533	iter.fold(Self::one(), \|acc: Complex, c: Complex\| acc * c)
1534	}
1535	}
1536
1537	impl<'a, T: 'a + Num + Clone> Product<&'a Complex<T>> for Complex<T> {
1538	fn product<I>(iter: I) -> Self
1539	where
1540	I: Iterator<Item = &'a Complex<T>>,
1541	{
1542	iter.fold(Self::one(), \|acc: Complex, c: &'a Complex\| acc * c)
1543	}
1544	}
1545
1546	#[cfg(feature = "serde")]
1547	impl<T> serde::Serialize for Complex<T>
1548	where
1549	T: serde::Serialize,
1550	{
1551	fn serialize<S>(&self, serializer: S) -> Result<S::Ok, S::Error>
1552	where
1553	S: serde::Serializer,
1554	{
1555	(&self.re, &self.im).serialize(serializer)
1556	}
1557	}
1558
1559	#[cfg(feature = "serde")]
1560	impl<'de, T> serde::Deserialize<'de> for Complex<T>
1561	where
1562	T: serde::Deserialize<'de>,
1563	{
1564	fn deserialize<D>(deserializer: D) -> Result<Self, D::Error>
1565	where
1566	D: serde::Deserializer<'de>,
1567	{
1568	let (re, im) = serde::Deserialize::deserialize(deserializer)?;
1569	Ok(Self::new(re, im))
1570	}
1571	}
1572
1573	#[derive(Debug, PartialEq)]
1574	pub struct ParseComplexError<E> {
1575	kind: ComplexErrorKind<E>,
1576	}
1577
1578	#[derive(Debug, PartialEq)]
1579	enum ComplexErrorKind<E> {
1580	ParseError(E),
1581	ExprError,
1582	UnsupportedRadix,
1583	}
1584
1585	impl<E> ParseComplexError<E> {
1586	fn expr_error() -> Self {
1587	ParseComplexError {
1588	kind: ComplexErrorKind::ExprError,
1589	}
1590	}
1591
1592	fn unsupported_radix() -> Self {
1593	ParseComplexError {
1594	kind: ComplexErrorKind::UnsupportedRadix,
1595	}
1596	}
1597
1598	fn from_error(error: E) -> Self {
1599	ParseComplexError {
1600	kind: ComplexErrorKind::ParseError(error),
1601	}
1602	}
1603	}
1604
1605	#[cfg(feature = "std")]
1606	impl<E: Error> Error for ParseComplexError<E> {
1607	#[allow(deprecated)]
1608	fn description(&self) -> &str {
1609	match self.kind {
1610	ComplexErrorKind::ParseError(ref e: &E) => e.description(),
1611	ComplexErrorKind::ExprError => "invalid or unsupported complex expression",
1612	ComplexErrorKind::UnsupportedRadix => "unsupported radix for conversion",
1613	}
1614	}
1615	}
1616
1617	impl<E: fmt::Display> fmt::Display for ParseComplexError<E> {
1618	fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
1619	match self.kind {
1620	ComplexErrorKind::ParseError(ref e: &E) => e.fmt(f),
1621	ComplexErrorKind::ExprError => "invalid or unsupported complex expression".fmt(f),
1622	ComplexErrorKind::UnsupportedRadix => "unsupported radix for conversion".fmt(f),
1623	}
1624	}
1625	}
1626
1627	#[cfg(test)]
1628	fn hash<T: hash::Hash>(x: &T) -> u64 {
1629	use std::collections::hash_map::RandomState;
1630	use std::hash::{BuildHasher, Hasher};
1631	let mut hasher = <RandomState as BuildHasher>::Hasher::new();
1632	x.hash(&mut hasher);
1633	hasher.finish()
1634	}
1635
1636	#[cfg(test)]
1637	pub(crate) mod test {
1638	#![allow(non_upper_case_globals)]
1639
1640	use super::{Complex, Complex64};
1641	use super::{ComplexErrorKind, ParseComplexError};
1642	use core::f64;
1643	use core::str::FromStr;
1644
1645	use std::string::{String, ToString};
1646
1647	use num_traits::{Num, One, Zero};
1648
1649	pub const _0_0i: Complex64 = Complex::new(`0.0`, `0.0`);
1650	pub const _1_0i: Complex64 = Complex::new(`1.0`, `0.0`);
1651	pub const _1_1i: Complex64 = Complex::new(`1.0`, `1.0`);
1652	pub const _0_1i: Complex64 = Complex::new(`0.0`, `1.0`);
1653	pub const _neg1_1i: Complex64 = Complex::new(`-1.0`, `1.0`);
1654	pub const _05_05i: Complex64 = Complex::new(`0.5`, `0.5`);
1655	pub const all_consts: [Complex64; `5`] = [_0_0i, _1_0i, _1_1i, _neg1_1i, _05_05i];
1656	pub const _4_2i: Complex64 = Complex::new(`4.0`, `2.0`);
1657	pub const _1_infi: Complex64 = Complex::new(`1.0`, f64::INFINITY);
1658	pub const _neg1_infi: Complex64 = Complex::new(`-1.0`, f64::INFINITY);
1659	pub const _1_nani: Complex64 = Complex::new(`1.0`, f64::NAN);
1660	pub const _neg1_nani: Complex64 = Complex::new(`-1.0`, f64::NAN);
1661	pub const _inf_0i: Complex64 = Complex::new(f64::INFINITY, `0.0`);
1662	pub const _neginf_1i: Complex64 = Complex::new(f64::NEG_INFINITY, `1.0`);
1663	pub const _neginf_neg1i: Complex64 = Complex::new(f64::NEG_INFINITY, `-1.0`);
1664	pub const _inf_1i: Complex64 = Complex::new(f64::INFINITY, `1.0`);
1665	pub const _inf_neg1i: Complex64 = Complex::new(f64::INFINITY, `-1.0`);
1666	pub const _neginf_infi: Complex64 = Complex::new(f64::NEG_INFINITY, f64::INFINITY);
1667	pub const _inf_infi: Complex64 = Complex::new(f64::INFINITY, f64::INFINITY);
1668	pub const _neginf_nani: Complex64 = Complex::new(f64::NEG_INFINITY, f64::NAN);
1669	pub const _inf_nani: Complex64 = Complex::new(f64::INFINITY, f64::NAN);
1670	pub const _nan_0i: Complex64 = Complex::new(f64::NAN, `0.0`);
1671	pub const _nan_1i: Complex64 = Complex::new(f64::NAN, `1.0`);
1672	pub const _nan_neg1i: Complex64 = Complex::new(f64::NAN, `-1.0`);
1673	pub const _nan_nani: Complex64 = Complex::new(f64::NAN, f64::NAN);
1674
1675	#[test]
1676	fn test_consts() {
1677	// check our constants are what Complex::new creates
1678	fn test(c: Complex64, r: f64, i: f64) {
1679	assert_eq!(c, Complex::new(r, i));
1680	}
1681	test(_0_0i, `0.0`, `0.0`);
1682	test(_1_0i, `1.0`, `0.0`);
1683	test(_1_1i, `1.0`, `1.0`);
1684	test(_neg1_1i, `-1.0`, `1.0`);
1685	test(_05_05i, `0.5`, `0.5`);
1686
1687	assert_eq!(_0_0i, Zero::zero());
1688	assert_eq!(_1_0i, One::one());
1689	}
1690
1691	#[test]
1692	fn test_scale_unscale() {
1693	assert_eq!(_05_05i.scale(`2.0`), _1_1i);
1694	assert_eq!(_1_1i.unscale(`2.0`), _05_05i);
1695	for &c in all_consts.iter() {
1696	assert_eq!(c.scale(`2.0`).unscale(`2.0`), c);
1697	}
1698	}
1699
1700	#[test]
1701	fn test_conj() {
1702	for &c in all_consts.iter() {
1703	assert_eq!(c.conj(), Complex::new(c.re, -c.im));
1704	assert_eq!(c.conj().conj(), c);
1705	}
1706	}
1707
1708	#[test]
1709	fn test_inv() {
1710	assert_eq!(_1_1i.inv(), _05_05i.conj());
1711	assert_eq!(_1_0i.inv(), _1_0i.inv());
1712	}
1713
1714	#[test]
1715	#[should_panic]
1716	fn test_divide_by_zero_natural() {
1717	let n = Complex::new(`2`, `3`);
1718	let d = Complex::new(`0`, `0`);
1719	let _x = n / d;
1720	}
1721
1722	#[test]
1723	fn test_inv_zero() {
1724	// FIXME #20: should this really fail, or just NaN?
1725	assert!(_0_0i.inv().is_nan());
1726	}
1727
1728	#[test]
1729	#[allow(clippy::float_cmp)]
1730	fn test_l1_norm() {
1731	assert_eq!(_0_0i.l1_norm(), `0.0`);
1732	assert_eq!(_1_0i.l1_norm(), `1.0`);
1733	assert_eq!(_1_1i.l1_norm(), `2.0`);
1734	assert_eq!(_0_1i.l1_norm(), `1.0`);
1735	assert_eq!(_neg1_1i.l1_norm(), `2.0`);
1736	assert_eq!(_05_05i.l1_norm(), `1.0`);
1737	assert_eq!(_4_2i.l1_norm(), `6.0`);
1738	}
1739
1740	#[test]
1741	fn test_pow() {
1742	for c in all_consts.iter() {
1743	assert_eq!(c.powi(`0`), _1_0i);
1744	let mut pos = _1_0i;
1745	let mut neg = _1_0i;
1746	for i in `1i32`..`20` {
1747	pos *= c;
1748	assert_eq!(pos, c.powi(i));
1749	if c.is_zero() {
1750	assert!(c.powi(-i).is_nan());
1751	} else {
1752	neg /= c;
1753	assert_eq!(neg, c.powi(-i));
1754	}
1755	}
1756	}
1757	}
1758
1759	#[cfg(any(feature = "std", feature = "libm"))]
1760	pub(crate) mod float {
1761
1762	use core::f64::INFINITY;
1763
1764	use super::*;
1765	use num_traits::{Float, Pow};
1766
1767	#[test]
1768	fn test_cis() {
1769	assert!(close(Complex::cis(`0.0` * f64::consts::PI), _1_0i));
1770	assert!(close(Complex::cis(`0.5` * f64::consts::PI), _0_1i));
1771	assert!(close(Complex::cis(`1.0` * f64::consts::PI), -_1_0i));
1772	assert!(close(Complex::cis(`1.5` * f64::consts::PI), -_0_1i));
1773	assert!(close(Complex::cis(`2.0` * f64::consts::PI), _1_0i));
1774	}
1775
1776	#[test]
1777	#[cfg_attr(target_arch = "x86", ignore)]
1778	// FIXME #7158: (maybe?) currently failing on x86.
1779	#[allow(clippy::float_cmp)]
1780	fn test_norm() {
1781	fn test(c: Complex64, ns: f64) {
1782	assert_eq!(c.norm_sqr(), ns);
1783	assert_eq!(c.norm(), ns.sqrt())
1784	}
1785	test(_0_0i, `0.0`);
1786	test(_1_0i, `1.0`);
1787	test(_1_1i, `2.0`);
1788	test(_neg1_1i, `2.0`);
1789	test(_05_05i, `0.5`);
1790	}
1791
1792	#[test]
1793	fn test_arg() {
1794	fn test(c: Complex64, arg: f64) {
1795	assert!((c.arg() - arg).abs() < `1.0e-6`)
1796	}
1797	test(_1_0i, `0.0`);
1798	test(_1_1i, `0.25` * f64::consts::PI);
1799	test(_neg1_1i, `0.75` * f64::consts::PI);
1800	test(_05_05i, `0.25` * f64::consts::PI);
1801	}
1802
1803	#[test]
1804	fn test_polar_conv() {
1805	fn test(c: Complex64) {
1806	let (r, theta) = c.to_polar();
1807	assert!((c - Complex::from_polar(r, theta)).norm() < `1e-6`);
1808	}
1809	for &c in all_consts.iter() {
1810	test(c);
1811	}
1812	}
1813
1814	pub(crate) fn close(a: Complex64, b: Complex64) -> bool {
1815	close_to_tol(a, b, `1e-10`)
1816	}
1817
1818	fn close_to_tol(a: Complex64, b: Complex64, tol: f64) -> bool {
1819	// returns true if a and b are reasonably close
1820	let close = (a == b) \|\| (a - b).norm() < tol;
1821	if !close {
1822	println!("{:?} != {:?}", a, b);
1823	}
1824	close
1825	}
1826
1827	// Version that also works if re or im are +inf, -inf, or nan
1828	fn close_naninf(a: Complex64, b: Complex64) -> bool {
1829	close_naninf_to_tol(a, b, `1.0e-10`)
1830	}
1831
1832	fn close_naninf_to_tol(a: Complex64, b: Complex64, tol: f64) -> bool {
1833	let mut close = `true`;
1834
1835	// Compare the real parts
1836	if a.re.is_finite() {
1837	if b.re.is_finite() {
1838	close = (a.re == b.re) \|\| (a.re - b.re).abs() < tol;
1839	} else {
1840	close = `false`;
1841	}
1842	} else if (a.re.is_nan() && !b.re.is_nan())
1843	\|\| (a.re.is_infinite()
1844	&& a.re.is_sign_positive()
1845	&& !(b.re.is_infinite() && b.re.is_sign_positive()))
1846	\|\| (a.re.is_infinite()
1847	&& a.re.is_sign_negative()
1848	&& !(b.re.is_infinite() && b.re.is_sign_negative()))
1849	{
1850	close = `false`;
1851	}
1852
1853	// Compare the imaginary parts
1854	if a.im.is_finite() {
1855	if b.im.is_finite() {
1856	close &= (a.im == b.im) \|\| (a.im - b.im).abs() < tol;
1857	} else {
1858	close = `false`;
1859	}
1860	} else if (a.im.is_nan() && !b.im.is_nan())
1861	\|\| (a.im.is_infinite()
1862	&& a.im.is_sign_positive()
1863	&& !(b.im.is_infinite() && b.im.is_sign_positive()))
1864	\|\| (a.im.is_infinite()
1865	&& a.im.is_sign_negative()
1866	&& !(b.im.is_infinite() && b.im.is_sign_negative()))
1867	{
1868	close = `false`;
1869	}
1870
1871	if close == `false` {
1872	println!("{:?} != {:?}", a, b);
1873	}
1874	close
1875	}
1876
1877	#[test]
1878	fn test_exp2() {
1879	assert!(close(_0_0i.exp2(), _1_0i));
1880	}
1881
1882	#[test]
1883	fn test_exp() {
1884	assert!(close(_1_0i.exp(), _1_0i.scale(f64::consts::E)));
1885	assert!(close(_0_0i.exp(), _1_0i));
1886	assert!(close(_0_1i.exp(), Complex::new(`1.0`.cos(), `1.0`.sin())));
1887	assert!(close(_05_05i.exp() * _05_05i.exp(), _1_1i.exp()));
1888	assert!(close(
1889	_0_1i.scale(-f64::consts::PI).exp(),
1890	_1_0i.scale(-`1.0`)
1891	));
1892	for &c in all_consts.iter() {
1893	// e^conj(z) = conj(e^z)
1894	assert!(close(c.conj().exp(), c.exp().conj()));
1895	// e^(z + 2 pi i) = e^z
1896	assert!(close(
1897	c.exp(),
1898	(c + _0_1i.scale(f64::consts::PI * `2.0`)).exp()
1899	));
1900	}
1901
1902	// The test values below were taken from https://en.cppreference.com/w/cpp/numeric/complex/exp
1903	assert!(close_naninf(_1_infi.exp(), _nan_nani));
1904	assert!(close_naninf(_neg1_infi.exp(), _nan_nani));
1905	assert!(close_naninf(_1_nani.exp(), _nan_nani));
1906	assert!(close_naninf(_neg1_nani.exp(), _nan_nani));
1907	assert!(close_naninf(_inf_0i.exp(), _inf_0i));
1908	assert!(close_naninf(_neginf_1i.exp(), `0.0` * Complex::cis(`1.0`)));
1909	assert!(close_naninf(_neginf_neg1i.exp(), `0.0` * Complex::cis(-`1.0`)));
1910	assert!(close_naninf(
1911	_inf_1i.exp(),
1912	f64::INFINITY * Complex::cis(`1.0`)
1913	));
1914	assert!(close_naninf(
1915	_inf_neg1i.exp(),
1916	f64::INFINITY * Complex::cis(-`1.0`)
1917	));
1918	assert!(close_naninf(_neginf_infi.exp(), _0_0i)); // Note: ±0±0i: signs of zeros are unspecified
1919	assert!(close_naninf(_inf_infi.exp(), _inf_nani)); // Note: ±∞+NaNi: sign of the real part is unspecified*
1920	assert!(close_naninf(_neginf_nani.exp(), _0_0i)); // Note: ±0±0i: signs of zeros are unspecified
1921	assert!(close_naninf(_inf_nani.exp(), _inf_nani)); // Note: ±∞+NaNi: sign of the real part is unspecified*
1922	assert!(close_naninf(_nan_0i.exp(), _nan_0i));
1923	assert!(close_naninf(_nan_1i.exp(), _nan_nani));
1924	assert!(close_naninf(_nan_neg1i.exp(), _nan_nani));
1925	assert!(close_naninf(_nan_nani.exp(), _nan_nani));
1926	}
1927
1928	#[test]
1929	fn test_ln() {
1930	assert!(close(_1_0i.ln(), _0_0i));
1931	assert!(close(_0_1i.ln(), _0_1i.scale(f64::consts::PI / `2.0`)));
1932	assert!(close(_0_0i.ln(), Complex::new(f64::neg_infinity(), `0.0`)));
1933	assert!(close(
1934	(_neg1_1i * _05_05i).ln(),
1935	_neg1_1i.ln() + _05_05i.ln()
1936	));
1937	for &c in all_consts.iter() {
1938	// ln(conj(z() = conj(ln(z))
1939	assert!(close(c.conj().ln(), c.ln().conj()));
1940	// for this branch, -pi <= arg(ln(z)) <= pi
1941	assert!(-f64::consts::PI <= c.ln().arg() && c.ln().arg() <= f64::consts::PI);
1942	}
1943	}
1944
1945	#[test]
1946	fn test_powc() {
1947	let a = Complex::new(`2.0`, `-3.0`);
1948	let b = Complex::new(`3.0`, `0.0`);
1949	assert!(close(a.powc(b), a.powf(b.re)));
1950	assert!(close(b.powc(a), a.expf(b.re)));
1951	let c = Complex::new(`1.0` / `3.0`, `0.1`);
1952	assert!(close_to_tol(
1953	a.powc(c),
1954	Complex::new(`1.65826`, -`0.33502`),
1955	`1e-5`
1956	));
1957	let z = Complex::new(`0.0`, `0.0`);
1958	assert!(close(z.powc(b), z));
1959	assert!(z.powc(Complex64::new(`0.`, INFINITY)).is_nan());
1960	assert!(z.powc(Complex64::new(`10.`, INFINITY)).is_nan());
1961	assert!(z.powc(Complex64::new(INFINITY, INFINITY)).is_nan());
1962	assert!(close(z.powc(Complex64::new(INFINITY, `0.`)), z));
1963	assert!(z.powc(Complex64::new(-`1.`, `0.`)).re.is_infinite());
1964	assert!(z.powc(Complex64::new(-`1.`, `0.`)).im.is_nan());
1965
1966	for c in all_consts.iter() {
1967	assert_eq!(c.powc(_0_0i), _1_0i);
1968	}
1969	assert_eq!(_nan_nani.powc(_0_0i), _1_0i);
1970	}
1971
1972	#[test]
1973	fn test_powf() {
1974	let c = Complex64::new(`2.0`, `-1.0`);
1975	let expected = Complex64::new(`-0.8684746`, `-16.695934`);
1976	assert!(close_to_tol(c.powf(`3.5`), expected, `1e-5`));
1977	assert!(close_to_tol(Pow::pow(c, `3.5_f64`), expected, `1e-5`));
1978	assert!(close_to_tol(Pow::pow(c, `3.5_f32`), expected, `1e-5`));
1979
1980	for c in all_consts.iter() {
1981	assert_eq!(c.powf(`0.0`), _1_0i);
1982	}
1983	assert_eq!(_nan_nani.powf(`0.0`), _1_0i);
1984	}
1985
1986	#[test]
1987	fn test_log() {
1988	let c = Complex::new(`2.0`, `-1.0`);
1989	let r = c.log(`10.0`);
1990	assert!(close_to_tol(r, Complex::new(`0.349485`, -`0.20135958`), `1e-5`));
1991	}
1992
1993	#[test]
1994	fn test_log2() {
1995	assert!(close(_1_0i.log2(), _0_0i));
1996	}
1997
1998	#[test]
1999	fn test_log10() {
2000	assert!(close(_1_0i.log10(), _0_0i));
2001	}
2002
2003	#[test]
2004	fn test_some_expf_cases() {
2005	let c = Complex::new(`2.0`, `-1.0`);
2006	let r = c.expf(`10.0`);
2007	assert!(close_to_tol(r, Complex::new(-`66.82015`, -`74.39803`), `1e-5`));
2008
2009	let c = Complex::new(`5.0`, `-2.0`);
2010	let r = c.expf(`3.4`);
2011	assert!(close_to_tol(r, Complex::new(-`349.25`, -`290.63`), `1e-2`));
2012
2013	let c = Complex::new(`-1.5`, `2.0` / `3.0`);
2014	let r = c.expf(`1.0` / `3.0`);
2015	assert!(close_to_tol(r, Complex::new(`3.8637`, -`3.4745`), `1e-2`));
2016	}
2017
2018	#[test]
2019	fn test_sqrt() {
2020	assert!(close(_0_0i.sqrt(), _0_0i));
2021	assert!(close(_1_0i.sqrt(), _1_0i));
2022	assert!(close(Complex::new(-`1.0`, `0.0`).sqrt(), _0_1i));
2023	assert!(close(Complex::new(-`1.0`, -`0.0`).sqrt(), _0_1i.scale(-`1.0`)));
2024	assert!(close(_0_1i.sqrt(), _05_05i.scale(`2.0`.sqrt())));
2025	for &c in all_consts.iter() {
2026	// sqrt(conj(z() = conj(sqrt(z))
2027	assert!(close(c.conj().sqrt(), c.sqrt().conj()));
2028	// for this branch, -pi/2 <= arg(sqrt(z)) <= pi/2
2029	assert!(
2030	-f64::consts::FRAC_PI_2 <= c.sqrt().arg()
2031	&& c.sqrt().arg() <= f64::consts::FRAC_PI_2
2032	);
2033	// sqrt(z) sqrt(z) = z*
2034	assert!(close(c.sqrt() * c.sqrt(), c));
2035	}
2036	}
2037
2038	#[test]
2039	fn test_sqrt_real() {
2040	for n in (`0`..`100`).map(f64::from) {
2041	// √(n² + 0i) = n + 0i
2042	let n2 = n * n;
2043	assert_eq!(Complex64::new(n2, `0.0`).sqrt(), Complex64::new(n, `0.0`));
2044	// √(-n² + 0i) = 0 + ni
2045	assert_eq!(Complex64::new(-n2, `0.0`).sqrt(), Complex64::new(`0.0`, n));
2046	// √(-n² - 0i) = 0 - ni
2047	assert_eq!(Complex64::new(-n2, -`0.0`).sqrt(), Complex64::new(`0.0`, -n));
2048	}
2049	}
2050
2051	#[test]
2052	fn test_sqrt_imag() {
2053	for n in (`0`..`100`).map(f64::from) {
2054	// √(0 + n²i) = n e^(iπ/4)
2055	let n2 = n * n;
2056	assert!(close(
2057	Complex64::new(`0.0`, n2).sqrt(),
2058	Complex64::from_polar(n, f64::consts::FRAC_PI_4)
2059	));
2060	// √(0 - n²i) = n e^(-iπ/4)
2061	assert!(close(
2062	Complex64::new(`0.0`, -n2).sqrt(),
2063	Complex64::from_polar(n, -f64::consts::FRAC_PI_4)
2064	));
2065	}
2066	}
2067
2068	#[test]
2069	fn test_cbrt() {
2070	assert!(close(_0_0i.cbrt(), _0_0i));
2071	assert!(close(_1_0i.cbrt(), _1_0i));
2072	assert!(close(
2073	Complex::new(-`1.0`, `0.0`).cbrt(),
2074	Complex::new(`0.5`, `0.75`.sqrt())
2075	));
2076	assert!(close(
2077	Complex::new(-`1.0`, -`0.0`).cbrt(),
2078	Complex::new(`0.5`, -(`0.75`.sqrt()))
2079	));
2080	assert!(close(_0_1i.cbrt(), Complex::new(`0.75`.sqrt(), `0.5`)));
2081	assert!(close(_0_1i.conj().cbrt(), Complex::new(`0.75`.sqrt(), -`0.5`)));
2082	for &c in all_consts.iter() {
2083	// cbrt(conj(z() = conj(cbrt(z))
2084	assert!(close(c.conj().cbrt(), c.cbrt().conj()));
2085	// for this branch, -pi/3 <= arg(cbrt(z)) <= pi/3
2086	assert!(
2087	-f64::consts::FRAC_PI_3 <= c.cbrt().arg()
2088	&& c.cbrt().arg() <= f64::consts::FRAC_PI_3
2089	);
2090	// cbrt(z) cbrt(z) cbrt(z) = z*
2091	assert!(close(c.cbrt() * c.cbrt() * c.cbrt(), c));
2092	}
2093	}
2094
2095	#[test]
2096	fn test_cbrt_real() {
2097	for n in (`0`..`100`).map(f64::from) {
2098	// ∛(n³ + 0i) = n + 0i
2099	let n3 = n * n * n;
2100	assert!(close(
2101	Complex64::new(n3, `0.0`).cbrt(),
2102	Complex64::new(n, `0.0`)
2103	));
2104	// ∛(-n³ + 0i) = n e^(iπ/3)
2105	assert!(close(
2106	Complex64::new(-n3, `0.0`).cbrt(),
2107	Complex64::from_polar(n, f64::consts::FRAC_PI_3)
2108	));
2109	// ∛(-n³ - 0i) = n e^(-iπ/3)
2110	assert!(close(
2111	Complex64::new(-n3, -`0.0`).cbrt(),
2112	Complex64::from_polar(n, -f64::consts::FRAC_PI_3)
2113	));
2114	}
2115	}
2116
2117	#[test]
2118	fn test_cbrt_imag() {
2119	for n in (`0`..`100`).map(f64::from) {
2120	// ∛(0 + n³i) = n e^(iπ/6)
2121	let n3 = n * n * n;
2122	assert!(close(
2123	Complex64::new(`0.0`, n3).cbrt(),
2124	Complex64::from_polar(n, f64::consts::FRAC_PI_6)
2125	));
2126	// ∛(0 - n³i) = n e^(-iπ/6)
2127	assert!(close(
2128	Complex64::new(`0.0`, -n3).cbrt(),
2129	Complex64::from_polar(n, -f64::consts::FRAC_PI_6)
2130	));
2131	}
2132	}
2133
2134	#[test]
2135	fn test_sin() {
2136	assert!(close(_0_0i.sin(), _0_0i));
2137	assert!(close(_1_0i.scale(f64::consts::PI * `2.0`).sin(), _0_0i));
2138	assert!(close(_0_1i.sin(), _0_1i.scale(`1.0`.sinh())));
2139	for &c in all_consts.iter() {
2140	// sin(conj(z)) = conj(sin(z))
2141	assert!(close(c.conj().sin(), c.sin().conj()));
2142	// sin(-z) = -sin(z)
2143	assert!(close(c.scale(-`1.0`).sin(), c.sin().scale(-`1.0`)));
2144	}
2145	}
2146
2147	#[test]
2148	fn test_cos() {
2149	assert!(close(_0_0i.cos(), _1_0i));
2150	assert!(close(_1_0i.scale(f64::consts::PI * `2.0`).cos(), _1_0i));
2151	assert!(close(_0_1i.cos(), _1_0i.scale(`1.0`.cosh())));
2152	for &c in all_consts.iter() {
2153	// cos(conj(z)) = conj(cos(z))
2154	assert!(close(c.conj().cos(), c.cos().conj()));
2155	// cos(-z) = cos(z)
2156	assert!(close(c.scale(-`1.0`).cos(), c.cos()));
2157	}
2158	}
2159
2160	#[test]
2161	fn test_tan() {
2162	assert!(close(_0_0i.tan(), _0_0i));
2163	assert!(close(_1_0i.scale(f64::consts::PI / `4.0`).tan(), _1_0i));
2164	assert!(close(_1_0i.scale(f64::consts::PI).tan(), _0_0i));
2165	for &c in all_consts.iter() {
2166	// tan(conj(z)) = conj(tan(z))
2167	assert!(close(c.conj().tan(), c.tan().conj()));
2168	// tan(-z) = -tan(z)
2169	assert!(close(c.scale(-`1.0`).tan(), c.tan().scale(-`1.0`)));
2170	}
2171	}
2172
2173	#[test]
2174	fn test_asin() {
2175	assert!(close(_0_0i.asin(), _0_0i));
2176	assert!(close(_1_0i.asin(), _1_0i.scale(f64::consts::PI / `2.0`)));
2177	assert!(close(
2178	_1_0i.scale(-`1.0`).asin(),
2179	_1_0i.scale(-f64::consts::PI / `2.0`)
2180	));
2181	assert!(close(_0_1i.asin(), _0_1i.scale((`1.0` + `2.0`.sqrt()).ln())));
2182	for &c in all_consts.iter() {
2183	// asin(conj(z)) = conj(asin(z))
2184	assert!(close(c.conj().asin(), c.asin().conj()));
2185	// asin(-z) = -asin(z)
2186	assert!(close(c.scale(-`1.0`).asin(), c.asin().scale(-`1.0`)));
2187	// for this branch, -pi/2 <= asin(z).re <= pi/2
2188	assert!(
2189	-f64::consts::PI / `2.0` <= c.asin().re && c.asin().re <= f64::consts::PI / `2.0`
2190	);
2191	}
2192	}
2193
2194	#[test]
2195	fn test_acos() {
2196	assert!(close(_0_0i.acos(), _1_0i.scale(f64::consts::PI / `2.0`)));
2197	assert!(close(_1_0i.acos(), _0_0i));
2198	assert!(close(
2199	_1_0i.scale(-`1.0`).acos(),
2200	_1_0i.scale(f64::consts::PI)
2201	));
2202	assert!(close(
2203	_0_1i.acos(),
2204	Complex::new(f64::consts::PI / `2.0`, (`2.0`.sqrt() - `1.0`).ln())
2205	));
2206	for &c in all_consts.iter() {
2207	// acos(conj(z)) = conj(acos(z))
2208	assert!(close(c.conj().acos(), c.acos().conj()));
2209	// for this branch, 0 <= acos(z).re <= pi
2210	assert!(`0.0` <= c.acos().re && c.acos().re <= f64::consts::PI);
2211	}
2212	}
2213
2214	#[test]
2215	fn test_atan() {
2216	assert!(close(_0_0i.atan(), _0_0i));
2217	assert!(close(_1_0i.atan(), _1_0i.scale(f64::consts::PI / `4.0`)));
2218	assert!(close(
2219	_1_0i.scale(-`1.0`).atan(),
2220	_1_0i.scale(-f64::consts::PI / `4.0`)
2221	));
2222	assert!(close(_0_1i.atan(), Complex::new(`0.0`, f64::infinity())));
2223	for &c in all_consts.iter() {
2224	// atan(conj(z)) = conj(atan(z))
2225	assert!(close(c.conj().atan(), c.atan().conj()));
2226	// atan(-z) = -atan(z)
2227	assert!(close(c.scale(-`1.0`).atan(), c.atan().scale(-`1.0`)));
2228	// for this branch, -pi/2 <= atan(z).re <= pi/2
2229	assert!(
2230	-f64::consts::PI / `2.0` <= c.atan().re && c.atan().re <= f64::consts::PI / `2.0`
2231	);
2232	}
2233	}
2234
2235	#[test]
2236	fn test_sinh() {
2237	assert!(close(_0_0i.sinh(), _0_0i));
2238	assert!(close(
2239	_1_0i.sinh(),
2240	_1_0i.scale((f64::consts::E - `1.0` / f64::consts::E) / `2.0`)
2241	));
2242	assert!(close(_0_1i.sinh(), _0_1i.scale(`1.0`.sin())));
2243	for &c in all_consts.iter() {
2244	// sinh(conj(z)) = conj(sinh(z))
2245	assert!(close(c.conj().sinh(), c.sinh().conj()));
2246	// sinh(-z) = -sinh(z)
2247	assert!(close(c.scale(-`1.0`).sinh(), c.sinh().scale(-`1.0`)));
2248	}
2249	}
2250
2251	#[test]
2252	fn test_cosh() {
2253	assert!(close(_0_0i.cosh(), _1_0i));
2254	assert!(close(
2255	_1_0i.cosh(),
2256	_1_0i.scale((f64::consts::E + `1.0` / f64::consts::E) / `2.0`)
2257	));
2258	assert!(close(_0_1i.cosh(), _1_0i.scale(`1.0`.cos())));
2259	for &c in all_consts.iter() {
2260	// cosh(conj(z)) = conj(cosh(z))
2261	assert!(close(c.conj().cosh(), c.cosh().conj()));
2262	// cosh(-z) = cosh(z)
2263	assert!(close(c.scale(-`1.0`).cosh(), c.cosh()));
2264	}
2265	}
2266
2267	#[test]
2268	fn test_tanh() {
2269	assert!(close(_0_0i.tanh(), _0_0i));
2270	assert!(close(
2271	_1_0i.tanh(),
2272	_1_0i.scale((f64::consts::E.powi(`2`) - `1.0`) / (f64::consts::E.powi(`2`) + `1.0`))
2273	));
2274	assert!(close(_0_1i.tanh(), _0_1i.scale(`1.0`.tan())));
2275	for &c in all_consts.iter() {
2276	// tanh(conj(z)) = conj(tanh(z))
2277	assert!(close(c.conj().tanh(), c.conj().tanh()));
2278	// tanh(-z) = -tanh(z)
2279	assert!(close(c.scale(-`1.0`).tanh(), c.tanh().scale(-`1.0`)));
2280	}
2281	}
2282
2283	#[test]
2284	fn test_asinh() {
2285	assert!(close(_0_0i.asinh(), _0_0i));
2286	assert!(close(_1_0i.asinh(), _1_0i.scale(`1.0` + `2.0`.sqrt()).ln()));
2287	assert!(close(_0_1i.asinh(), _0_1i.scale(f64::consts::PI / `2.0`)));
2288	assert!(close(
2289	_0_1i.asinh().scale(-`1.0`),
2290	_0_1i.scale(-f64::consts::PI / `2.0`)
2291	));
2292	for &c in all_consts.iter() {
2293	// asinh(conj(z)) = conj(asinh(z))
2294	assert!(close(c.conj().asinh(), c.conj().asinh()));
2295	// asinh(-z) = -asinh(z)
2296	assert!(close(c.scale(-`1.0`).asinh(), c.asinh().scale(-`1.0`)));
2297	// for this branch, -pi/2 <= asinh(z).im <= pi/2
2298	assert!(
2299	-f64::consts::PI / `2.0` <= c.asinh().im && c.asinh().im <= f64::consts::PI / `2.0`
2300	);
2301	}
2302	}
2303
2304	#[test]
2305	fn test_acosh() {
2306	assert!(close(_0_0i.acosh(), _0_1i.scale(f64::consts::PI / `2.0`)));
2307	assert!(close(_1_0i.acosh(), _0_0i));
2308	assert!(close(
2309	_1_0i.scale(-`1.0`).acosh(),
2310	_0_1i.scale(f64::consts::PI)
2311	));
2312	for &c in all_consts.iter() {
2313	// acosh(conj(z)) = conj(acosh(z))
2314	assert!(close(c.conj().acosh(), c.conj().acosh()));
2315	// for this branch, -pi <= acosh(z).im <= pi and 0 <= acosh(z).re
2316	assert!(
2317	-f64::consts::PI <= c.acosh().im
2318	&& c.acosh().im <= f64::consts::PI
2319	&& `0.0` <= c.cosh().re
2320	);
2321	}
2322	}
2323
2324	#[test]
2325	fn test_atanh() {
2326	assert!(close(_0_0i.atanh(), _0_0i));
2327	assert!(close(_0_1i.atanh(), _0_1i.scale(f64::consts::PI / `4.0`)));
2328	assert!(close(_1_0i.atanh(), Complex::new(f64::infinity(), `0.0`)));
2329	for &c in all_consts.iter() {
2330	// atanh(conj(z)) = conj(atanh(z))
2331	assert!(close(c.conj().atanh(), c.conj().atanh()));
2332	// atanh(-z) = -atanh(z)
2333	assert!(close(c.scale(-`1.0`).atanh(), c.atanh().scale(-`1.0`)));
2334	// for this branch, -pi/2 <= atanh(z).im <= pi/2
2335	assert!(
2336	-f64::consts::PI / `2.0` <= c.atanh().im && c.atanh().im <= f64::consts::PI / `2.0`
2337	);
2338	}
2339	}
2340
2341	#[test]
2342	fn test_exp_ln() {
2343	for &c in all_consts.iter() {
2344	// e^ln(z) = z
2345	assert!(close(c.ln().exp(), c));
2346	}
2347	}
2348
2349	#[test]
2350	fn test_exp2_log() {
2351	for &c in all_consts.iter() {
2352	// 2^log2(z) = z
2353	assert!(close(c.log2().exp2(), c));
2354	}
2355	}
2356
2357	#[test]
2358	fn test_trig_to_hyperbolic() {
2359	for &c in all_consts.iter() {
2360	// sin(iz) = i sinh(z)
2361	assert!(close((_0_1i * c).sin(), _0_1i * c.sinh()));
2362	// cos(iz) = cosh(z)
2363	assert!(close((_0_1i * c).cos(), c.cosh()));
2364	// tan(iz) = i tanh(z)
2365	assert!(close((_0_1i * c).tan(), _0_1i * c.tanh()));
2366	}
2367	}
2368
2369	#[test]
2370	fn test_trig_identities() {
2371	for &c in all_consts.iter() {
2372	// tan(z) = sin(z)/cos(z)
2373	assert!(close(c.tan(), c.sin() / c.cos()));
2374	// sin(z)^2 + cos(z)^2 = 1
2375	assert!(close(c.sin() * c.sin() + c.cos() * c.cos(), _1_0i));
2376
2377	// sin(asin(z)) = z
2378	assert!(close(c.asin().sin(), c));
2379	// cos(acos(z)) = z
2380	assert!(close(c.acos().cos(), c));
2381	// tan(atan(z)) = z
2382	// i and -i are branch points
2383	if c != _0_1i && c != _0_1i.scale(`-1.0`) {
2384	assert!(close(c.atan().tan(), c));
2385	}
2386
2387	// sin(z) = (e^(iz) - e^(-iz))/(2i)
2388	assert!(close(
2389	((_0_1i * c).exp() - (_0_1i * c).exp().inv()) / _0_1i.scale(`2.0`),
2390	c.sin()
2391	));
2392	// cos(z) = (e^(iz) + e^(-iz))/2
2393	assert!(close(
2394	((_0_1i * c).exp() + (_0_1i * c).exp().inv()).unscale(`2.0`),
2395	c.cos()
2396	));
2397	// tan(z) = i (1 - e^(2iz))/(1 + e^(2iz))
2398	assert!(close(
2399	_0_1i * (_1_0i - (_0_1i * c).scale(`2.0`).exp())
2400	/ (_1_0i + (_0_1i * c).scale(`2.0`).exp()),
2401	c.tan()
2402	));
2403	}
2404	}
2405
2406	#[test]
2407	fn test_hyperbolic_identites() {
2408	for &c in all_consts.iter() {
2409	// tanh(z) = sinh(z)/cosh(z)
2410	assert!(close(c.tanh(), c.sinh() / c.cosh()));
2411	// cosh(z)^2 - sinh(z)^2 = 1
2412	assert!(close(c.cosh() * c.cosh() - c.sinh() * c.sinh(), _1_0i));
2413
2414	// sinh(asinh(z)) = z
2415	assert!(close(c.asinh().sinh(), c));
2416	// cosh(acosh(z)) = z
2417	assert!(close(c.acosh().cosh(), c));
2418	// tanh(atanh(z)) = z
2419	// 1 and -1 are branch points
2420	if c != _1_0i && c != _1_0i.scale(`-1.0`) {
2421	assert!(close(c.atanh().tanh(), c));
2422	}
2423
2424	// sinh(z) = (e^z - e^(-z))/2
2425	assert!(close((c.exp() - c.exp().inv()).unscale(`2.0`), c.sinh()));
2426	// cosh(z) = (e^z + e^(-z))/2
2427	assert!(close((c.exp() + c.exp().inv()).unscale(`2.0`), c.cosh()));
2428	// tanh(z) = ( e^(2z) - 1)/(e^(2z) + 1)
2429	assert!(close(
2430	(c.scale(`2.0`).exp() - _1_0i) / (c.scale(`2.0`).exp() + _1_0i),
2431	c.tanh()
2432	));
2433	}
2434	}
2435	}
2436
2437	// Test both a + b and a += b
2438	macro_rules! test_a_op_b {
2439	($a:ident + $b:expr, $answer:expr) => {
2440	assert_eq!($a + $b, $answer);
2441	assert_eq!(
2442	{
2443	let mut x = $a;
2444	x += $b;
2445	x
2446	},
2447	$answer
2448	);
2449	};
2450	($a:ident - $b:expr, $answer:expr) => {
2451	assert_eq!($a - $b, $answer);
2452	assert_eq!(
2453	{
2454	let mut x = $a;
2455	x -= $b;
2456	x
2457	},
2458	$answer
2459	);
2460	};
2461	($a:ident * $b:expr, $answer:expr) => {
2462	assert_eq!($a * $b, $answer);
2463	assert_eq!(
2464	{
2465	let mut x = $a;
2466	x *= $b;
2467	x
2468	},
2469	$answer
2470	);
2471	};
2472	($a:ident / $b:expr, $answer:expr) => {
2473	assert_eq!($a / $b, $answer);
2474	assert_eq!(
2475	{
2476	let mut x = $a;
2477	x /= $b;
2478	x
2479	},
2480	$answer
2481	);
2482	};
2483	($a:ident % $b:expr, $answer:expr) => {
2484	assert_eq!($a % $b, $answer);
2485	assert_eq!(
2486	{
2487	let mut x = $a;
2488	x %= $b;
2489	x
2490	},
2491	$answer
2492	);
2493	};
2494	}
2495
2496	// Test both a + b and a + &b
2497	macro_rules! test_op {
2498	($a:ident $op:tt $b:expr, $answer:expr) => {
2499	test_a_op_b!($a $op $b, $answer);
2500	test_a_op_b!($a $op &$b, $answer);
2501	};
2502	}
2503
2504	mod complex_arithmetic {
2505	use super::{_05_05i, _0_0i, _0_1i, _1_0i, _1_1i, _4_2i, _neg1_1i, all_consts};
2506	use num_traits::{MulAdd, MulAddAssign, Zero};
2507
2508	#[test]
2509	fn test_add() {
2510	test_op!(_05_05i + _05_05i, _1_1i);
2511	test_op!(_0_1i + _1_0i, _1_1i);
2512	test_op!(_1_0i + _neg1_1i, _0_1i);
2513
2514	for &c in all_consts.iter() {
2515	test_op!(_0_0i + c, c);
2516	test_op!(c + _0_0i, c);
2517	}
2518	}
2519
2520	#[test]
2521	fn test_sub() {
2522	test_op!(_05_05i - _05_05i, _0_0i);
2523	test_op!(_0_1i - _1_0i, _neg1_1i);
2524	test_op!(_0_1i - _neg1_1i, _1_0i);
2525
2526	for &c in all_consts.iter() {
2527	test_op!(c - _0_0i, c);
2528	test_op!(c - c, _0_0i);
2529	}
2530	}
2531
2532	#[test]
2533	fn test_mul() {
2534	test_op!(_05_05i * _05_05i, _0_1i.unscale(`2.0`));
2535	test_op!(_1_1i * _0_1i, _neg1_1i);
2536
2537	// i^2 & i^4
2538	test_op!(_0_1i * _0_1i, -_1_0i);
2539	assert_eq!(_0_1i * _0_1i * _0_1i * _0_1i, _1_0i);
2540
2541	for &c in all_consts.iter() {
2542	test_op!(c * _1_0i, c);
2543	test_op!(_1_0i * c, c);
2544	}
2545	}
2546
2547	#[test]
2548	#[cfg(any(feature = "std", feature = "libm"))]
2549	fn test_mul_add_float() {
2550	assert_eq!(_05_05i.mul_add(_05_05i, _0_0i), _05_05i * _05_05i + _0_0i);
2551	assert_eq!(_05_05i * _05_05i + _0_0i, _05_05i.mul_add(_05_05i, _0_0i));
2552	assert_eq!(_0_1i.mul_add(_0_1i, _0_1i), _neg1_1i);
2553	assert_eq!(_1_0i.mul_add(_1_0i, _1_0i), _1_0i * _1_0i + _1_0i);
2554	assert_eq!(_1_0i * _1_0i + _1_0i, _1_0i.mul_add(_1_0i, _1_0i));
2555
2556	let mut x = _1_0i;
2557	x.mul_add_assign(_1_0i, _1_0i);
2558	assert_eq!(x, _1_0i * _1_0i + _1_0i);
2559
2560	for &a in &all_consts {
2561	for &b in &all_consts {
2562	for &c in &all_consts {
2563	let abc = a * b + c;
2564	assert_eq!(a.mul_add(b, c), abc);
2565	let mut x = a;
2566	x.mul_add_assign(b, c);
2567	assert_eq!(x, abc);
2568	}
2569	}
2570	}
2571	}
2572
2573	#[test]
2574	fn test_mul_add() {
2575	use super::Complex;
2576	const _0_0i: Complex<i32> = Complex { re: `0`, im: `0` };
2577	const _1_0i: Complex<i32> = Complex { re: `1`, im: `0` };
2578	const _1_1i: Complex<i32> = Complex { re: `1`, im: `1` };
2579	const _0_1i: Complex<i32> = Complex { re: `0`, im: `1` };
2580	const _neg1_1i: Complex<i32> = Complex { re: `-1`, im: `1` };
2581	const all_consts: [Complex<i32>; `5`] = [_0_0i, _1_0i, _1_1i, _0_1i, _neg1_1i];
2582
2583	assert_eq!(_1_0i.mul_add(_1_0i, _0_0i), _1_0i * _1_0i + _0_0i);
2584	assert_eq!(_1_0i * _1_0i + _0_0i, _1_0i.mul_add(_1_0i, _0_0i));
2585	assert_eq!(_0_1i.mul_add(_0_1i, _0_1i), _neg1_1i);
2586	assert_eq!(_1_0i.mul_add(_1_0i, _1_0i), _1_0i * _1_0i + _1_0i);
2587	assert_eq!(_1_0i * _1_0i + _1_0i, _1_0i.mul_add(_1_0i, _1_0i));
2588
2589	let mut x = _1_0i;
2590	x.mul_add_assign(_1_0i, _1_0i);
2591	assert_eq!(x, _1_0i * _1_0i + _1_0i);
2592
2593	for &a in &all_consts {
2594	for &b in &all_consts {
2595	for &c in &all_consts {
2596	let abc = a * b + c;
2597	assert_eq!(a.mul_add(b, c), abc);
2598	let mut x = a;
2599	x.mul_add_assign(b, c);
2600	assert_eq!(x, abc);
2601	}
2602	}
2603	}
2604	}
2605
2606	#[test]
2607	fn test_div() {
2608	test_op!(_neg1_1i / _0_1i, _1_1i);
2609	for &c in all_consts.iter() {
2610	if c != Zero::zero() {
2611	test_op!(c / c, _1_0i);
2612	}
2613	}
2614	}
2615
2616	#[test]
2617	fn test_rem() {
2618	test_op!(_neg1_1i % _0_1i, _0_0i);
2619	test_op!(_4_2i % _0_1i, _0_0i);
2620	test_op!(_05_05i % _0_1i, _05_05i);
2621	test_op!(_05_05i % _1_1i, _05_05i);
2622	assert_eq!((_4_2i + _05_05i) % _0_1i, _05_05i);
2623	assert_eq!((_4_2i + _05_05i) % _1_1i, _05_05i);
2624	}
2625
2626	#[test]
2627	fn test_neg() {
2628	assert_eq!(-_1_0i + _0_1i, _neg1_1i);
2629	assert_eq!((-_0_1i) * _0_1i, _1_0i);
2630	for &c in all_consts.iter() {
2631	assert_eq!(-(-c), c);
2632	}
2633	}
2634	}
2635
2636	mod real_arithmetic {
2637	use super::super::Complex;
2638	use super::{_4_2i, _neg1_1i};
2639
2640	#[test]
2641	fn test_add() {
2642	test_op!(_4_2i + `0.5`, Complex::new(`4.5`, `2.0`));
2643	assert_eq!(`0.5` + _4_2i, Complex::new(`4.5`, `2.0`));
2644	}
2645
2646	#[test]
2647	fn test_sub() {
2648	test_op!(_4_2i - `0.5`, Complex::new(`3.5`, `2.0`));
2649	assert_eq!(`0.5` - _4_2i, Complex::new(-`3.5`, -`2.0`));
2650	}
2651
2652	#[test]
2653	fn test_mul() {
2654	assert_eq!(_4_2i * `0.5`, Complex::new(`2.0`, `1.0`));
2655	assert_eq!(`0.5` * _4_2i, Complex::new(`2.0`, `1.0`));
2656	}
2657
2658	#[test]
2659	fn test_div() {
2660	assert_eq!(_4_2i / `0.5`, Complex::new(`8.0`, `4.0`));
2661	assert_eq!(`0.5` / _4_2i, Complex::new(`0.1`, -`0.05`));
2662	}
2663
2664	#[test]
2665	fn test_rem() {
2666	assert_eq!(_4_2i % `2.0`, Complex::new(`0.0`, `0.0`));
2667	assert_eq!(_4_2i % `3.0`, Complex::new(`1.0`, `2.0`));
2668	assert_eq!(`3.0` % _4_2i, Complex::new(`3.0`, `0.0`));
2669	assert_eq!(_neg1_1i % `2.0`, _neg1_1i);
2670	assert_eq!(-_4_2i % `3.0`, Complex::new(-`1.0`, -`2.0`));
2671	}
2672
2673	#[test]
2674	fn test_div_rem_gaussian() {
2675	// These would overflow with `norm_sqr` division.
2676	let max = Complex::new(`255u8`, `255u8`);
2677	assert_eq!(max / `200`, Complex::new(`1`, `1`));
2678	assert_eq!(max % `200`, Complex::new(`55`, `55`));
2679	}
2680	}
2681
2682	#[test]
2683	fn test_to_string() {
2684	fn test(c: Complex64, s: String) {
2685	assert_eq!(c.to_string(), s);
2686	}
2687	test(_0_0i, "0+0i".to_string());
2688	test(_1_0i, "1+0i".to_string());
2689	test(_0_1i, "0+1i".to_string());
2690	test(_1_1i, "1+1i".to_string());
2691	test(_neg1_1i, "-1+1i".to_string());
2692	test(-_neg1_1i, "1-1i".to_string());
2693	test(_05_05i, "0.5+0.5i".to_string());
2694	}
2695
2696	#[test]
2697	fn test_string_formatting() {
2698	let a = Complex::new(`1.23456`, `123.456`);
2699	assert_eq!(format!("{}", a), "1.23456+123.456i");
2700	assert_eq!(format!("{:.2}", a), "1.23+123.46i");
2701	assert_eq!(format!("{:.2e}", a), "1.23e0+1.23e2i");
2702	assert_eq!(format!("{:+.2E}", a), "+1.23E0+1.23E2i");
2703	#[cfg(feature = "std")]
2704	assert_eq!(format!("{:+20.2E}", a), " +1.23E0+1.23E2i");
2705
2706	let b = Complex::new(`0x80`, `0xff`);
2707	assert_eq!(format!("{:X}", b), "80+FFi");
2708	assert_eq!(format!("{:#x}", b), "0x80+0xffi");
2709	assert_eq!(format!("{:+#b}", b), "+0b10000000+0b11111111i");
2710	assert_eq!(format!("{:+#o}", b), "+0o200+0o377i");
2711	#[cfg(feature = "std")]
2712	assert_eq!(format!("{:+#16o}", b), " +0o200+0o377i");
2713
2714	let c = Complex::new(`-10`, `-10000`);
2715	assert_eq!(format!("{}", c), "-10-10000i");
2716	#[cfg(feature = "std")]
2717	assert_eq!(format!("{:16}", c), " -10-10000i");
2718	}
2719
2720	#[test]
2721	fn test_hash() {
2722	let a = Complex::new(`0i32`, `0i32`);
2723	let b = Complex::new(`1i32`, `0i32`);
2724	let c = Complex::new(`0i32`, `1i32`);
2725	assert!(crate::hash(&a) != crate::hash(&b));
2726	assert!(crate::hash(&b) != crate::hash(&c));
2727	assert!(crate::hash(&c) != crate::hash(&a));
2728	}
2729
2730	#[test]
2731	fn test_hashset() {
2732	use std::collections::HashSet;
2733	let a = Complex::new(`0i32`, `0i32`);
2734	let b = Complex::new(`1i32`, `0i32`);
2735	let c = Complex::new(`0i32`, `1i32`);
2736
2737	let set: HashSet<_> = [a, b, c].iter().cloned().collect();
2738	assert!(set.contains(&a));
2739	assert!(set.contains(&b));
2740	assert!(set.contains(&c));
2741	assert!(!set.contains(&(a + b + c)));
2742	}
2743
2744	#[test]
2745	fn test_is_nan() {
2746	assert!(!_1_1i.is_nan());
2747	let a = Complex::new(f64::NAN, f64::NAN);
2748	assert!(a.is_nan());
2749	}
2750
2751	#[test]
2752	fn test_is_nan_special_cases() {
2753	let a = Complex::new(`0f64`, f64::NAN);
2754	let b = Complex::new(f64::NAN, `0f64`);
2755	assert!(a.is_nan());
2756	assert!(b.is_nan());
2757	}
2758
2759	#[test]
2760	fn test_is_infinite() {
2761	let a = Complex::new(`2f64`, f64::INFINITY);
2762	assert!(a.is_infinite());
2763	}
2764
2765	#[test]
2766	fn test_is_finite() {
2767	assert!(_1_1i.is_finite())
2768	}
2769
2770	#[test]
2771	fn test_is_normal() {
2772	let a = Complex::new(`0f64`, f64::NAN);
2773	let b = Complex::new(`2f64`, f64::INFINITY);
2774	assert!(!a.is_normal());
2775	assert!(!b.is_normal());
2776	assert!(_1_1i.is_normal());
2777	}
2778
2779	#[test]
2780	fn test_from_str() {
2781	fn test(z: Complex64, s: &str) {
2782	assert_eq!(FromStr::from_str(s), Ok(z));
2783	}
2784	test(_0_0i, "0 + 0i");
2785	test(_0_0i, "0+0j");
2786	test(_0_0i, "0 - 0j");
2787	test(_0_0i, "0-0i");
2788	test(_0_0i, "0i + 0");
2789	test(_0_0i, "0");
2790	test(_0_0i, "-0");
2791	test(_0_0i, "0i");
2792	test(_0_0i, "0j");
2793	test(_0_0i, "+0j");
2794	test(_0_0i, "-0i");
2795
2796	test(_1_0i, "1 + 0i");
2797	test(_1_0i, "1+0j");
2798	test(_1_0i, "1 - 0j");
2799	test(_1_0i, "+1-0i");
2800	test(_1_0i, "-0j+1");
2801	test(_1_0i, "1");
2802
2803	test(_1_1i, "1 + i");
2804	test(_1_1i, "1+j");
2805	test(_1_1i, "1 + 1j");
2806	test(_1_1i, "1+1i");
2807	test(_1_1i, "i + 1");
2808	test(_1_1i, "1i+1");
2809	test(_1_1i, "+j+1");
2810
2811	test(_0_1i, "0 + i");
2812	test(_0_1i, "0+j");
2813	test(_0_1i, "-0 + j");
2814	test(_0_1i, "-0+i");
2815	test(_0_1i, "0 + 1i");
2816	test(_0_1i, "0+1j");
2817	test(_0_1i, "-0 + 1j");
2818	test(_0_1i, "-0+1i");
2819	test(_0_1i, "j + 0");
2820	test(_0_1i, "i");
2821	test(_0_1i, "j");
2822	test(_0_1i, "1j");
2823
2824	test(_neg1_1i, "-1 + i");
2825	test(_neg1_1i, "-1+j");
2826	test(_neg1_1i, "-1 + 1j");
2827	test(_neg1_1i, "-1+1i");
2828	test(_neg1_1i, "1i-1");
2829	test(_neg1_1i, "j + -1");
2830
2831	test(_05_05i, "0.5 + 0.5i");
2832	test(_05_05i, "0.5+0.5j");
2833	test(_05_05i, "5e-1+0.5j");
2834	test(_05_05i, "5E-1 + 0.5j");
2835	test(_05_05i, "5E-1i + 0.5");
2836	test(_05_05i, "0.05e+1j + 50E-2");
2837	}
2838
2839	#[test]
2840	fn test_from_str_radix() {
2841	fn test(z: Complex64, s: &str, radix: u32) {
2842	let res: Result<Complex64, <Complex64 as Num>::FromStrRadixErr> =
2843	Num::from_str_radix(s, radix);
2844	assert_eq!(res.unwrap(), z)
2845	}
2846	test(_4_2i, "4+2i", `10`);
2847	test(Complex::new(`15.0`, `32.0`), "F+20i", `16`);
2848	test(Complex::new(`15.0`, `32.0`), "1111+100000i", `2`);
2849	test(Complex::new(`-15.0`, `-32.0`), "-F-20i", `16`);
2850	test(Complex::new(`-15.0`, `-32.0`), "-1111-100000i", `2`);
2851
2852	fn test_error(s: &str, radix: u32) -> ParseComplexError<<f64 as Num>::FromStrRadixErr> {
2853	let res = Complex64::from_str_radix(s, radix);
2854
2855	res.expect_err(&format!("Expected failure on input {:?}", s))
2856	}
2857
2858	let err = test_error("1ii", `19`);
2859	if let ComplexErrorKind::UnsupportedRadix = err.kind {
2860	/ pass /
2861	} else {
2862	panic!("Expected failure on invalid radix, got {:?}", err);
2863	}
2864
2865	let err = test_error("1 + 0", `16`);
2866	if let ComplexErrorKind::ExprError = err.kind {
2867	/ pass /
2868	} else {
2869	panic!("Expected failure on expr error, got {:?}", err);
2870	}
2871	}
2872
2873	#[test]
2874	#[should_panic(expected = "radix is too high")]
2875	fn test_from_str_radix_fail() {
2876	// ensure we preserve the underlying panic on radix > 36
2877	let _complex = Complex64::from_str_radix("1", `37`);
2878	}
2879
2880	#[test]
2881	fn test_from_str_fail() {
2882	fn test(s: &str) {
2883	let complex: Result<Complex64, _> = FromStr::from_str(s);
2884	assert!(
2885	complex.is_err(),
2886	"complex {:?} -> {:?} should be an error",
2887	s,
2888	complex
2889	);
2890	}
2891	test("foo");
2892	test("6E");
2893	test("0 + 2.718");
2894	test("1 - -2i");
2895	test("314e-2ij");
2896	test("4.3j - i");
2897	test("1i - 2i");
2898	test("+ 1 - 3.0i");
2899	}
2900
2901	#[test]
2902	fn test_sum() {
2903	let v = vec![_0_1i, _1_0i];
2904	assert_eq!(v.iter().sum::<Complex64>(), _1_1i);
2905	assert_eq!(v.into_iter().sum::<Complex64>(), _1_1i);
2906	}
2907
2908	#[test]
2909	fn test_prod() {
2910	let v = vec![_0_1i, _1_0i];
2911	assert_eq!(v.iter().product::<Complex64>(), _0_1i);
2912	assert_eq!(v.into_iter().product::<Complex64>(), _0_1i);
2913	}
2914
2915	#[test]
2916	fn test_zero() {
2917	let zero = Complex64::zero();
2918	assert!(zero.is_zero());
2919
2920	let mut c = Complex::new(`1.23`, `4.56`);
2921	assert!(!c.is_zero());
2922	assert_eq!(c + zero, c);
2923
2924	c.set_zero();
2925	assert!(c.is_zero());
2926	}
2927
2928	#[test]
2929	fn test_one() {
2930	let one = Complex64::one();
2931	assert!(one.is_one());
2932
2933	let mut c = Complex::new(`1.23`, `4.56`);
2934	assert!(!c.is_one());
2935	assert_eq!(c * one, c);
2936
2937	c.set_one();
2938	assert!(c.is_one());
2939	}
2940
2941	#[test]
2942	#[allow(clippy::float_cmp)]
2943	fn test_const() {
2944	const R: f64 = `12.3`;
2945	const I: f64 = `-4.5`;
2946	const C: Complex64 = Complex::new(R, I);
2947
2948	assert_eq!(C.re, `12.3`);
2949	assert_eq!(C.im, -`4.5`);
2950	}
2951	}
2952