Skip to main content

abels_complex/complex/
polar.rs

1pub type ComplexPolar32 = ComplexPolar<f32>;
2pub type ComplexPolar64 = ComplexPolar<f64>;
3use crate::traits::Number;
4use core::fmt;
5
6use super::Complex as Rectangular;
7use core::ops::*;
8
9/// Creates a complex number in polar form.
10#[inline(always)]
11#[must_use]
12pub const fn complex_polar<FT>(abs: FT, arg: FT) -> ComplexPolar<FT> {
13    ComplexPolar::new(abs, arg)
14}
15
16/// A complex number in polar form.
17#[derive(Clone, Copy, PartialEq, Debug, Default)]
18#[repr(C)]
19pub struct ComplexPolar<FT> {
20    pub abs: FT,
21    pub arg: FT,
22}
23
24impl<FT> ComplexPolar<FT> {
25    /// Creates a complex number.
26    pub const fn new(abs: FT, arg: FT) -> Self {
27        Self { abs, arg }
28    }
29}
30
31impl<FT: Number> ComplexPolar<FT> {
32    pub const ZERO: Self = Self::new(FT::ZERO, FT::ZERO);
33    pub const ONE: Self = Self::new(FT::ONE, FT::ZERO);
34
35    /// Computes the conjugate.
36    pub fn conjugate(self) -> Self {
37        Self::new(self.abs, -self.arg)
38    }
39
40    /// Computes the real component.
41    pub fn re(self) -> FT {
42        self.abs * self.arg.cos()
43    }
44
45    /// Computes the imaginary component.
46    pub fn im(self) -> FT {
47        self.abs * self.arg.sin()
48    }
49
50    /// Computes the squared absolute value.
51    pub fn abs_sq(self) -> FT {
52        self.abs * self.abs
53    }
54
55    /// Computes the reciprocal.
56    pub fn recip(self) -> Self {
57        Self::new(self.abs.recip(), -self.arg)
58    }
59
60    /// Computes the principal square root.
61    pub fn sqrt(self) -> Self {
62        let two = FT::ONE + FT::ONE;
63        Self::new(self.abs.sqrt(), self.arg / two)
64    }
65
66    /// Convert to rectangular form.
67    pub fn to_rectangular(self) -> Rectangular<FT> {
68        let (sin, cos) = self.arg.sin_cos();
69        Rectangular::new(cos, sin) * self.abs
70    }
71
72    /// Computes `e^self` where `e` is the base of the natural logarithm.
73    pub fn exp(self) -> Self {
74        self.to_rectangular().exp()
75    }
76
77    /// Computes `2^self`.
78    pub fn exp2(self) -> Self {
79        self.to_rectangular().exp2()
80    }
81
82    /// Computes the principal natural logarithm.
83    pub fn ln(self) -> Rectangular<FT> {
84        Rectangular::new(self.abs.ln(), self.arg)
85    }
86
87    /// Computes the principal natural logarithm of `1 + self`.
88    ///
89    /// More numerically stable than `(self + 1).ln()` when `self ≈ 0`.
90    pub fn ln_1p(self) -> Rectangular<FT> {
91        self.to_rectangular().ln_1p()
92    }
93
94    /// Computes the principal logarithm in base 2.
95    pub fn log2(self) -> Rectangular<FT> {
96        self.ln() / FT::LN_2()
97    }
98
99    /// Computes the principal logarithm in base 10.
100    pub fn log10(self) -> Rectangular<FT> {
101        self.ln() / FT::LN_10()
102    }
103
104    /// Raises `self` to a floating point power.
105    pub fn powf(self, x: FT) -> Self {
106        if x < FT::ZERO && self.abs == FT::ZERO {
107            return Self::ZERO;
108        }
109        Self::new(self.abs.powf(x), self.arg * x)
110    }
111
112    /// Raises `self` to an integer power.
113    pub fn powi(self, n: i32) -> Self {
114        if n < 0 && self.abs == FT::ZERO {
115            return Self::ZERO;
116        }
117        Self::new(self.abs.powi(n), self.arg * FT::from_i32(n))
118    }
119
120    /// Normalizes the absolute value and the argument into the range `[0, ∞)` and `(-π, +π]` respectively.
121    pub fn normalize(mut self) -> Self {
122        self.arg = self.arg.rem_euclid(&FT::TAU());
123        if self.abs < FT::ZERO {
124            self.abs = -self.abs;
125            if self.arg <= FT::ZERO {
126                self.arg += FT::PI();
127            } else {
128                self.arg -= FT::PI();
129            }
130        } else if self.arg > FT::PI() {
131            self.arg -= FT::TAU();
132        } else if self.arg <= -FT::PI() {
133            self.arg += FT::TAU();
134        }
135        self
136    }
137}
138
139impl<FT: Number> Mul for ComplexPolar<FT> {
140    type Output = Self;
141    fn mul(mut self, other: Self) -> Self {
142        self *= other;
143        self
144    }
145}
146
147impl<FT: Number> Mul<FT> for ComplexPolar<FT> {
148    type Output = Self;
149    fn mul(mut self, re: FT) -> Self::Output {
150        self *= re;
151        self
152    }
153}
154
155impl<FT: Number> MulAssign for ComplexPolar<FT> {
156    fn mul_assign(&mut self, other: Self) {
157        self.abs *= other.abs;
158        self.arg += other.arg;
159    }
160}
161
162impl<FT: Number> MulAssign<FT> for ComplexPolar<FT> {
163    fn mul_assign(&mut self, re: FT) {
164        self.abs *= re;
165    }
166}
167
168impl<FT: Number> Div for ComplexPolar<FT> {
169    type Output = Self;
170    fn div(mut self, other: Self) -> Self {
171        self /= other;
172        self
173    }
174}
175
176impl<FT: Number> Div<FT> for ComplexPolar<FT> {
177    type Output = Self;
178    fn div(mut self, re: FT) -> Self {
179        self /= re;
180        self
181    }
182}
183
184impl<FT: Number> DivAssign for ComplexPolar<FT> {
185    fn div_assign(&mut self, other: Self) {
186        *self *= other.recip();
187    }
188}
189
190impl<FT: Number> DivAssign<FT> for ComplexPolar<FT> {
191    fn div_assign(&mut self, re: FT) {
192        self.abs /= re;
193    }
194}
195
196impl<FT: Number> Neg for ComplexPolar<FT> {
197    type Output = Self;
198    fn neg(mut self) -> Self {
199        self.abs = -self.abs;
200        self
201    }
202}
203
204impl<FT: Number + fmt::Display> fmt::Display for ComplexPolar<FT> {
205    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
206        fn fmt_x<FT: fmt::Display>(f: &mut fmt::Formatter, x: FT) -> fmt::Result {
207            if let Some(p) = f.precision() {
208                write!(f, "{x:.*}", p)
209            } else {
210                write!(f, "{x}")
211            }
212        }
213        let pi_radians = self.arg / FT::PI();
214        fmt_x(f, self.abs)?;
215        if pi_radians == FT::ZERO || self.abs == FT::ZERO {
216            Ok(())
217        } else if pi_radians == FT::ONE {
218            write!(f, "e^iπ")
219        } else {
220            write!(f, "e^")?;
221            fmt_x(f, pi_radians)?;
222            write!(f, "iπ")
223        }
224    }
225}
226
227impl<FT: Number> From<FT> for ComplexPolar<FT> {
228    fn from(value: FT) -> Self {
229        Self::new(value, FT::ZERO)
230    }
231}
232
233#[cfg(feature = "approx")]
234use approx::{AbsDiffEq, RelativeEq, UlpsEq};
235
236#[cfg(feature = "approx")]
237impl<FT: AbsDiffEq + Copy> AbsDiffEq for ComplexPolar<FT>
238where
239    <FT as AbsDiffEq>::Epsilon: Copy,
240{
241    type Epsilon = <FT as AbsDiffEq>::Epsilon;
242    fn default_epsilon() -> Self::Epsilon {
243        FT::default_epsilon()
244    }
245    fn abs_diff_eq(&self, other: &Self, epsilon: Self::Epsilon) -> bool {
246        FT::abs_diff_eq(&self.abs, &other.abs, epsilon)
247            && FT::abs_diff_eq(&self.arg, &other.arg, epsilon)
248    }
249}
250
251#[cfg(feature = "approx")]
252impl<FT: RelativeEq + Copy> RelativeEq for ComplexPolar<FT>
253where
254    <FT as AbsDiffEq>::Epsilon: Copy,
255{
256    fn default_max_relative() -> Self::Epsilon {
257        FT::default_max_relative()
258    }
259    fn relative_eq(
260        &self,
261        other: &Self,
262        epsilon: Self::Epsilon,
263        max_relative: Self::Epsilon,
264    ) -> bool {
265        FT::relative_eq(&self.abs, &other.abs, epsilon, max_relative)
266            && FT::relative_eq(&self.arg, &other.arg, epsilon, max_relative)
267    }
268}
269
270#[cfg(feature = "approx")]
271impl<FT: UlpsEq + Copy> UlpsEq for ComplexPolar<FT>
272where
273    <FT as AbsDiffEq>::Epsilon: Copy,
274{
275    fn default_max_ulps() -> u32 {
276        FT::default_max_ulps()
277    }
278    fn ulps_eq(&self, other: &Self, epsilon: Self::Epsilon, max_ulps: u32) -> bool {
279        FT::ulps_eq(&self.abs, &other.abs, epsilon, max_ulps)
280            && FT::ulps_eq(&self.arg, &other.arg, epsilon, max_ulps)
281    }
282}