use num::{FromPrimitive, Num, NumAssign, One, Signed, Zero};
use std::any::Any;
use std::fmt::{Debug, Display};
use std::ops::Neg;
use std::{f32, f64};
use crate::general::{Field, JoinSemilattice, MeetSemilattice, RealField, SubsetOf, SupersetOf};
#[cfg(not(feature = "std"))]
use num::Float;
#[allow(missing_docs)]
pub trait ComplexField:
SubsetOf<Self>
+ SupersetOf<f64>
+ Field
+ Copy
+ Num
+ NumAssign
+ FromPrimitive
+ Neg<Output = Self>
+ MeetSemilattice
+ JoinSemilattice
+ Send
+ Sync
+ Any
+ 'static
+ Debug
+ Display
{
type RealField: RealField;
fn from_real(re: Self::RealField) -> Self;
fn real(self) -> Self::RealField;
fn imaginary(self) -> Self::RealField;
fn modulus(self) -> Self::RealField;
fn modulus_squared(self) -> Self::RealField;
fn argument(self) -> Self::RealField;
fn norm1(self) -> Self::RealField;
fn scale(self, factor: Self::RealField) -> Self;
fn unscale(self, factor: Self::RealField) -> Self;
fn to_polar(self) -> (Self::RealField, Self::RealField) {
(self.modulus(), self.argument())
}
fn to_exp(self) -> (Self::RealField, Self) {
let m = self.modulus();
if !m.is_zero() {
(m, self.unscale(m))
} else {
(Self::RealField::zero(), Self::one())
}
}
fn signum(self) -> Self {
self.to_exp().1
}
fn floor(self) -> Self;
fn ceil(self) -> Self;
fn round(self) -> Self;
fn trunc(self) -> Self;
fn fract(self) -> Self;
fn mul_add(self, a: Self, b: Self) -> Self;
fn abs(self) -> Self::RealField;
fn hypot(self, other: Self) -> Self::RealField;
fn recip(self) -> Self;
fn conjugate(self) -> Self;
fn sin(self) -> Self;
fn cos(self) -> Self;
fn sin_cos(self) -> (Self, Self);
#[inline]
fn sinh_cosh(self) -> (Self, Self) {
(self.sinh(), self.cosh())
}
fn tan(self) -> Self;
fn asin(self) -> Self;
fn acos(self) -> Self;
fn atan(self) -> Self;
fn sinh(self) -> Self;
fn cosh(self) -> Self;
fn tanh(self) -> Self;
fn asinh(self) -> Self;
fn acosh(self) -> Self;
fn atanh(self) -> Self;
fn is_finite(&self) -> bool;
#[inline]
fn sinc(self) -> Self {
if self.is_zero() {
Self::one()
} else {
self.sin() / self
}
}
#[inline]
fn sinhc(self) -> Self {
if self.is_zero() {
Self::one()
} else {
self.sinh() / self
}
}
#[inline]
fn cosc(self) -> Self {
if self.is_zero() {
Self::one()
} else {
self.cos() / self
}
}
#[inline]
fn coshc(self) -> Self {
if self.is_zero() {
Self::one()
} else {
self.cosh() / self
}
}
fn log(self, base: Self::RealField) -> Self;
fn log2(self) -> Self;
fn log10(self) -> Self;
fn ln(self) -> Self;
fn ln_1p(self) -> Self;
fn sqrt(self) -> Self;
fn try_sqrt(self) -> Option<Self>;
fn exp(self) -> Self;
fn exp2(self) -> Self;
fn exp_m1(self) -> Self;
fn powi(self, n: i32) -> Self;
fn powf(self, n: Self::RealField) -> Self;
fn powc(self, n: Self) -> Self;
fn cbrt(self) -> Self;
}
macro_rules! impl_complex(
($($T:ty, $M:ident, $libm: ident);*) => ($(
impl ComplexField for $T {
type RealField = $T;
#[inline]
fn from_real(re: Self::RealField) -> Self {
re
}
#[inline]
fn real(self) -> Self::RealField {
self
}
#[inline]
fn imaginary(self) -> Self::RealField {
Self::zero()
}
#[inline]
fn norm1(self) -> Self::RealField {
$libm::abs(self)
}
#[inline]
fn modulus(self) -> Self::RealField {
$libm::abs(self)
}
#[inline]
fn modulus_squared(self) -> Self::RealField {
self * self
}
#[inline]
fn argument(self) -> Self::RealField {
if self >= Self::zero() {
Self::zero()
} else {
Self::pi()
}
}
#[inline]
fn to_exp(self) -> (Self, Self) {
if self >= Self::zero() {
(self, Self::one())
} else {
(-self, -Self::one())
}
}
#[inline]
fn recip(self) -> Self {
$M::recip(self)
}
#[inline]
fn conjugate(self) -> Self {
self
}
#[inline]
fn scale(self, factor: Self::RealField) -> Self {
self * factor
}
#[inline]
fn unscale(self, factor: Self::RealField) -> Self {
self / factor
}
#[inline]
fn floor(self) -> Self {
$libm::floor(self)
}
#[inline]
fn ceil(self) -> Self {
$libm::ceil(self)
}
#[inline]
fn round(self) -> Self {
$libm::round(self)
}
#[inline]
fn trunc(self) -> Self {
$libm::trunc(self)
}
#[inline]
fn fract(self) -> Self {
$libm::fract(self)
}
#[inline]
fn abs(self) -> Self {
$libm::abs(self)
}
#[inline]
fn signum(self) -> Self {
Signed::signum(&self)
}
#[inline]
fn mul_add(self, a: Self, b: Self) -> Self {
$libm::mul_add(self, a, b)
}
#[cfg(feature = "std")]
#[inline]
fn powi(self, n: i32) -> Self {
self.powi(n)
}
#[cfg(not(feature = "std"))]
#[inline]
fn powi(self, n: i32) -> Self {
$libm::powf(self, n as $T)
}
#[inline]
fn powf(self, n: Self) -> Self {
$libm::powf(self, n)
}
#[inline]
fn powc(self, n: Self) -> Self {
$libm::powf(self, n)
}
#[inline]
fn sqrt(self) -> Self {
$libm::sqrt(self)
}
#[inline]
fn try_sqrt(self) -> Option<Self> {
if self >= Self::zero() {
Some($libm::sqrt(self))
} else {
None
}
}
#[inline]
fn exp(self) -> Self {
$libm::exp(self)
}
#[inline]
fn exp2(self) -> Self {
$libm::exp2(self)
}
#[inline]
fn exp_m1(self) -> Self {
$libm::exp_m1(self)
}
#[inline]
fn ln_1p(self) -> Self {
$libm::ln_1p(self)
}
#[inline]
fn ln(self) -> Self {
$libm::ln(self)
}
#[inline]
fn log(self, base: Self) -> Self {
$libm::log(self, base)
}
#[inline]
fn log2(self) -> Self {
$libm::log2(self)
}
#[inline]
fn log10(self) -> Self {
$libm::log10(self)
}
#[inline]
fn cbrt(self) -> Self {
$libm::cbrt(self)
}
#[inline]
fn hypot(self, other: Self) -> Self::RealField {
$libm::hypot(self, other)
}
#[inline]
fn sin(self) -> Self {
$libm::sin(self)
}
#[inline]
fn cos(self) -> Self {
$libm::cos(self)
}
#[inline]
fn tan(self) -> Self {
$libm::tan(self)
}
#[inline]
fn asin(self) -> Self {
$libm::asin(self)
}
#[inline]
fn acos(self) -> Self {
$libm::acos(self)
}
#[inline]
fn atan(self) -> Self {
$libm::atan(self)
}
#[inline]
fn sin_cos(self) -> (Self, Self) {
$libm::sin_cos(self)
}
#[inline]
fn sinh(self) -> Self {
$libm::sinh(self)
}
#[inline]
fn cosh(self) -> Self {
$libm::cosh(self)
}
#[inline]
fn tanh(self) -> Self {
$libm::tanh(self)
}
#[inline]
fn asinh(self) -> Self {
$libm::asinh(self)
}
#[inline]
fn acosh(self) -> Self {
$libm::acosh(self)
}
#[inline]
fn atanh(self) -> Self {
$libm::atanh(self)
}
#[inline]
fn is_finite(&self) -> bool {
$M::is_finite(*self)
}
}
)*)
);
#[cfg(not(feature = "std"))]
impl_complex!(
f32, f32, Float;
f64, f64, Float
);
#[cfg(feature = "std")]
impl_complex!(
f32,f32,f32;
f64,f64,f64
);
impl<N: RealField> ComplexField for num_complex::Complex<N> {
type RealField = N;
#[inline]
fn from_real(re: Self::RealField) -> Self {
Self::new(re, Self::RealField::zero())
}
#[inline]
fn real(self) -> Self::RealField {
self.re
}
#[inline]
fn imaginary(self) -> Self::RealField {
self.im
}
#[inline]
fn argument(self) -> Self::RealField {
self.im.atan2(self.re)
}
#[inline]
fn modulus(self) -> Self::RealField {
self.re.hypot(self.im)
}
#[inline]
fn modulus_squared(self) -> Self::RealField {
self.re * self.re + self.im * self.im
}
#[inline]
fn norm1(self) -> Self::RealField {
self.re.abs() + self.im.abs()
}
#[inline]
fn recip(self) -> Self {
Self::one() / self
}
#[inline]
fn conjugate(self) -> Self {
self.conj()
}
#[inline]
fn scale(self, factor: Self::RealField) -> Self {
self * factor
}
#[inline]
fn unscale(self, factor: Self::RealField) -> Self {
self / factor
}
#[inline]
fn floor(self) -> Self {
Self::new(self.re.floor(), self.im.floor())
}
#[inline]
fn ceil(self) -> Self {
Self::new(self.re.ceil(), self.im.ceil())
}
#[inline]
fn round(self) -> Self {
Self::new(self.re.round(), self.im.round())
}
#[inline]
fn trunc(self) -> Self {
Self::new(self.re.trunc(), self.im.trunc())
}
#[inline]
fn fract(self) -> Self {
Self::new(self.re.fract(), self.im.fract())
}
#[inline]
fn mul_add(self, a: Self, b: Self) -> Self {
self * a + b
}
#[inline]
fn abs(self) -> Self::RealField {
self.modulus()
}
#[inline]
fn exp2(self) -> Self {
let _2 = N::one() + N::one();
num_complex::Complex::new(_2, N::zero()).powc(self)
}
#[inline]
fn exp_m1(self) -> Self {
self.exp() - Self::one()
}
#[inline]
fn ln_1p(self) -> Self {
(Self::one() + self).ln()
}
#[inline]
fn log2(self) -> Self {
let _2 = N::one() + N::one();
self.log(_2)
}
#[inline]
fn log10(self) -> Self {
let _10 = N::from_subset(&10.0f64);
self.log(_10)
}
#[inline]
fn cbrt(self) -> Self {
let one_third = N::from_subset(&(1.0 / 3.0));
self.powf(one_third)
}
#[inline]
fn powi(self, n: i32) -> Self {
let n = N::from_subset(&(n as f64));
self.powf(n)
}
#[inline]
fn is_finite(&self) -> bool {
self.re.is_finite() && self.im.is_finite()
}
#[inline]
fn exp(self) -> Self {
complex_from_polar(self.re.exp(), self.im)
}
#[inline]
fn ln(self) -> Self {
let (r, theta) = self.to_polar();
Self::new(r.ln(), theta)
}
#[inline]
fn sqrt(self) -> Self {
let two = N::one() + N::one();
let (r, theta) = self.to_polar();
complex_from_polar(r.sqrt(), theta / two)
}
#[inline]
fn try_sqrt(self) -> Option<Self> {
Some(self.sqrt())
}
#[inline]
fn hypot(self, b: Self) -> Self::RealField {
(self.modulus_squared() + b.modulus_squared()).sqrt()
}
#[inline]
fn powf(self, exp: Self::RealField) -> Self {
let (r, theta) = self.to_polar();
complex_from_polar(r.powf(exp), theta * exp)
}
#[inline]
fn log(self, base: N) -> Self {
let (r, theta) = self.to_polar();
Self::new(r.log(base), theta / base.ln())
}
#[inline]
fn powc(self, exp: Self) -> Self {
let (r, theta) = self.to_polar();
complex_from_polar(
r.powf(exp.re) * (-exp.im * theta).exp(),
exp.re * theta + exp.im * r.ln(),
)
}
#[inline]
fn sin(self) -> Self {
Self::new(
self.re.sin() * self.im.cosh(),
self.re.cos() * self.im.sinh(),
)
}
#[inline]
fn cos(self) -> Self {
Self::new(
self.re.cos() * self.im.cosh(),
-self.re.sin() * self.im.sinh(),
)
}
#[inline]
fn sin_cos(self) -> (Self, Self) {
let (rsin, rcos) = self.re.sin_cos();
let (isinh, icosh) = self.im.sinh_cosh();
let sin = Self::new(rsin * icosh, rcos * isinh);
let cos = Self::new(rcos * icosh, -rsin * isinh);
(sin, cos)
}
#[inline]
fn tan(self) -> Self {
let (two_re, two_im) = (self.re + self.re, self.im + self.im);
Self::new(two_re.sin(), two_im.sinh()).unscale(two_re.cos() + two_im.cosh())
}
#[inline]
fn asin(self) -> Self {
let i = Self::i();
-i * ((Self::one() - self * self).sqrt() + i * self).ln()
}
#[inline]
fn acos(self) -> Self {
let i = Self::i();
-i * (i * (Self::one() - self * self).sqrt() + self).ln()
}
#[inline]
fn atan(self) -> Self {
let i = Self::i();
let one = Self::one();
let two = one + one;
if self == i {
return Self::new(N::zero(), N::one() / N::zero());
} else if self == -i {
return Self::new(N::zero(), -N::one() / N::zero());
}
((one + i * self).ln() - (one - i * self).ln()) / (two * i)
}
#[inline]
fn sinh(self) -> Self {
Self::new(
self.re.sinh() * self.im.cos(),
self.re.cosh() * self.im.sin(),
)
}
#[inline]
fn cosh(self) -> Self {
Self::new(
self.re.cosh() * self.im.cos(),
self.re.sinh() * self.im.sin(),
)
}
#[inline]
fn sinh_cosh(self) -> (Self, Self) {
let (rsinh, rcosh) = self.re.sinh_cosh();
let (isin, icos) = self.im.sin_cos();
let sin = Self::new(rsinh * icos, rcosh * isin);
let cos = Self::new(rcosh * icos, rsinh * isin);
(sin, cos)
}
#[inline]
fn tanh(self) -> Self {
let (two_re, two_im) = (self.re + self.re, self.im + self.im);
Self::new(two_re.sinh(), two_im.sin()).unscale(two_re.cosh() + two_im.cos())
}
#[inline]
fn asinh(self) -> Self {
let one = Self::one();
(self + (one + self * self).sqrt()).ln()
}
#[inline]
fn acosh(self) -> Self {
let one = Self::one();
let two = one + one;
two * (((self + one) / two).sqrt() + ((self - one) / two).sqrt()).ln()
}
#[inline]
fn atanh(self) -> Self {
let one = Self::one();
let two = one + one;
if self == one {
return Self::new(N::one() / N::zero(), N::zero());
} else if self == -one {
return Self::new(-N::one() / N::zero(), N::zero());
}
((one + self).ln() - (one - self).ln()) / two
}
}
#[inline]
fn complex_from_polar<N: RealField>(r: N, theta: N) -> num_complex::Complex<N> {
num_complex::Complex::new(r * theta.cos(), r * theta.sin())
}