2019-07-06 20:19:09 +02:00
|
|
|
// Copyright (c) 2019 Alexander Medvednikov. All rights reserved.
|
|
|
|
// Use of this source code is governed by an MIT license
|
|
|
|
// that can be found in the LICENSE file.
|
|
|
|
|
2019-07-10 21:40:29 +02:00
|
|
|
module cmath
|
|
|
|
|
|
|
|
import math
|
2019-07-06 20:19:09 +02:00
|
|
|
|
|
|
|
struct Complex {
|
|
|
|
re f64
|
|
|
|
im f64
|
|
|
|
}
|
|
|
|
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn complex(re f64, im f64) Complex {
|
|
|
|
return Complex{re, im}
|
2019-07-06 20:19:09 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
// To String method
|
|
|
|
pub fn (c Complex) str() string {
|
|
|
|
mut out := '$c.re'
|
|
|
|
out += if c.im >= 0 {
|
|
|
|
'+$c.im'
|
|
|
|
}
|
|
|
|
else {
|
|
|
|
'$c.im'
|
|
|
|
}
|
|
|
|
out += 'i'
|
|
|
|
return out
|
|
|
|
}
|
|
|
|
|
2019-07-10 13:51:48 +02:00
|
|
|
// Complex Modulus value
|
|
|
|
// mod() and abs() return the same
|
2019-07-09 21:12:52 +02:00
|
|
|
pub fn (c Complex) abs() f64 {
|
2019-07-10 13:51:48 +02:00
|
|
|
return C.hypot(c.re, c.im)
|
2019-07-09 21:12:52 +02:00
|
|
|
}
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) mod() f64 {
|
|
|
|
return c.abs()
|
|
|
|
}
|
|
|
|
|
2019-07-09 21:12:52 +02:00
|
|
|
|
2019-07-07 16:13:17 +02:00
|
|
|
// Complex Angle
|
|
|
|
pub fn (c Complex) angle() f64 {
|
2019-07-10 21:40:29 +02:00
|
|
|
return math.atan2(c.im, c.re)
|
2019-07-07 16:13:17 +02:00
|
|
|
}
|
|
|
|
|
2019-07-06 20:19:09 +02:00
|
|
|
// Complex Addition c1 + c2
|
|
|
|
pub fn (c1 Complex) + (c2 Complex) Complex {
|
2019-07-10 13:51:48 +02:00
|
|
|
return Complex{c1.re + c2.re, c1.im + c2.im}
|
2019-07-06 20:19:09 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Substraction c1 - c2
|
|
|
|
pub fn (c1 Complex) - (c2 Complex) Complex {
|
2019-07-10 13:51:48 +02:00
|
|
|
return Complex{c1.re - c2.re, c1.im - c2.im}
|
2019-07-06 20:19:09 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Multiplication c1 * c2
|
|
|
|
// Currently Not Supported
|
|
|
|
// pub fn (c1 Complex) * (c2 Complex) Complex {
|
|
|
|
// return Complex{
|
|
|
|
// (c1.re * c2.re) + ((c1.im * c2.im) * -1),
|
|
|
|
// (c1.re * c2.im) + (c1.im * c2.re)
|
|
|
|
// }
|
|
|
|
// }
|
|
|
|
|
|
|
|
// Complex Division c1 / c2
|
|
|
|
// Currently Not Supported
|
|
|
|
// pub fn (c1 Complex) / (c2 Complex) Complex {
|
|
|
|
// denom := (c2.re * c2.re) + (c2.im * c2.im)
|
|
|
|
// return Complex {
|
|
|
|
// ((c1.re * c2.re) + ((c1.im * -c2.im) * -1))/denom,
|
|
|
|
// ((c1.re * -c2.im) + (c1.im * c2.re))/denom
|
|
|
|
// }
|
|
|
|
// }
|
|
|
|
|
|
|
|
// Complex Addition c1.add(c2)
|
|
|
|
pub fn (c1 Complex) add(c2 Complex) Complex {
|
|
|
|
return c1 + c2
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Subtraction c1.subtract(c2)
|
|
|
|
pub fn (c1 Complex) subtract(c2 Complex) Complex {
|
|
|
|
return c1 - c2
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Multiplication c1.multiply(c2)
|
|
|
|
pub fn (c1 Complex) multiply(c2 Complex) Complex {
|
|
|
|
return Complex{
|
|
|
|
(c1.re * c2.re) + ((c1.im * c2.im) * -1),
|
|
|
|
(c1.re * c2.im) + (c1.im * c2.re)
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Division c1.divide(c2)
|
|
|
|
pub fn (c1 Complex) divide(c2 Complex) Complex {
|
|
|
|
denom := (c2.re * c2.re) + (c2.im * c2.im)
|
|
|
|
return Complex {
|
2019-07-10 13:51:48 +02:00
|
|
|
((c1.re * c2.re) + ((c1.im * -c2.im) * -1)) / denom,
|
|
|
|
((c1.re * -c2.im) + (c1.im * c2.re)) / denom
|
2019-07-06 20:19:09 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Conjugate
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) conjugate() Complex{
|
|
|
|
return Complex{c.re, -c.im}
|
2019-07-06 20:19:09 +02:00
|
|
|
}
|
|
|
|
|
2019-07-08 19:41:37 +02:00
|
|
|
// Complex Additive Inverse
|
|
|
|
// Based on
|
|
|
|
// http://tutorial.math.lamar.edu/Extras/ComplexPrimer/Arithmetic.aspx
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) addinv() Complex {
|
|
|
|
return Complex{-c.re, -c.im}
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Multiplicative Inverse
|
|
|
|
// Based on
|
|
|
|
// http://tutorial.math.lamar.edu/Extras/ComplexPrimer/Arithmetic.aspx
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) mulinv() Complex {
|
2019-07-08 19:41:37 +02:00
|
|
|
return Complex {
|
2019-07-10 13:51:48 +02:00
|
|
|
c.re / (c.re * c.re + c.im * c.im),
|
|
|
|
-c.im / (c.re * c.re + c.im * c.im)
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Power
|
|
|
|
// Based on
|
|
|
|
// https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers/multiplying-and-dividing-complex-numbers-in-polar-form/a/complex-number-polar-form-review
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) pow(n f64) Complex {
|
2019-07-10 21:40:29 +02:00
|
|
|
r := math.pow(c.abs(), n)
|
2019-07-10 13:51:48 +02:00
|
|
|
angle := c.angle()
|
2019-07-08 19:41:37 +02:00
|
|
|
return Complex {
|
2019-07-10 21:40:29 +02:00
|
|
|
r * math.cos(n * angle),
|
|
|
|
r * math.sin(n * angle)
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex nth root
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) root(n f64) Complex {
|
|
|
|
return c.pow(1.0 / n)
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Exponential
|
|
|
|
// Using Euler's Identity
|
|
|
|
// Based on
|
|
|
|
// https://www.math.wisc.edu/~angenent/Free-Lecture-Notes/freecomplexnumbers.pdf
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) exp() Complex {
|
2019-07-10 21:40:29 +02:00
|
|
|
a := math.exp(c.re)
|
2019-07-08 19:41:37 +02:00
|
|
|
return Complex {
|
2019-07-10 21:40:29 +02:00
|
|
|
a * math.cos(c.im),
|
|
|
|
a * math.sin(c.im)
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Natural Logarithm
|
|
|
|
// Based on
|
|
|
|
// http://www.chemistrylearning.com/logarithm-of-complex-number/
|
2019-07-10 13:51:48 +02:00
|
|
|
pub fn (c Complex) ln() Complex {
|
2019-07-08 19:41:37 +02:00
|
|
|
return Complex {
|
2019-07-10 21:40:29 +02:00
|
|
|
math.log(c.abs()),
|
2019-07-10 13:51:48 +02:00
|
|
|
c.angle()
|
2019-07-08 19:41:37 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
2019-07-10 19:28:19 +02:00
|
|
|
// Complex Sin
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) sin() Complex {
|
|
|
|
return Complex{
|
2019-07-10 21:40:29 +02:00
|
|
|
math.sin(c.re) * math.cosh(c.im),
|
|
|
|
math.cos(c.re) * math.sinh(c.im)
|
2019-07-10 19:28:19 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Cosine
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) cos() Complex {
|
|
|
|
return Complex{
|
2019-07-10 21:40:29 +02:00
|
|
|
math.cos(c.re) * math.cosh(c.im),
|
|
|
|
-(math.sin(c.re) * math.sinh(c.im))
|
2019-07-10 19:28:19 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Tangent
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) tan() Complex {
|
|
|
|
return c.sin().divide(c.cos())
|
|
|
|
}
|
|
|
|
|
2019-07-11 15:35:06 +02:00
|
|
|
// Complex Arc Sin / Sin Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/summaryops.htm
|
|
|
|
pub fn (c Complex) asin() Complex {
|
|
|
|
return complex(0,-1).multiply(
|
|
|
|
complex(0,1)
|
|
|
|
.multiply(c)
|
|
|
|
.add(
|
|
|
|
complex(1,0)
|
|
|
|
.subtract(c.pow(2))
|
|
|
|
.root(2)
|
|
|
|
)
|
|
|
|
.ln()
|
|
|
|
)
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Arc Consine / Consine Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/summaryops.htm
|
|
|
|
pub fn (c Complex) acos() Complex {
|
|
|
|
return complex(0,-1).multiply(
|
|
|
|
c.add(
|
|
|
|
complex(0,1)
|
|
|
|
.multiply(
|
|
|
|
complex(1,0)
|
|
|
|
.subtract(c.pow(2))
|
|
|
|
.root(2)
|
|
|
|
)
|
|
|
|
)
|
|
|
|
.ln()
|
|
|
|
)
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Arc Tangent / Tangent Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/summaryops.htm
|
|
|
|
pub fn (c Complex) atan() Complex {
|
|
|
|
i := complex(0,1)
|
|
|
|
return complex(0,1.0/2).multiply(
|
|
|
|
i.add(c)
|
|
|
|
.divide(
|
|
|
|
i.subtract(c)
|
|
|
|
)
|
|
|
|
.ln()
|
|
|
|
)
|
|
|
|
}
|
|
|
|
|
2019-07-10 19:28:19 +02:00
|
|
|
// Complex Hyperbolic Sin
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) sinh() Complex {
|
|
|
|
return Complex{
|
2019-07-10 21:40:29 +02:00
|
|
|
math.cos(c.im) * math.sinh(c.re),
|
|
|
|
math.sin(c.im) * math.cosh(c.re)
|
2019-07-10 19:28:19 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Hyperbolic Cosine
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) cosh() Complex {
|
|
|
|
return Complex{
|
2019-07-10 21:40:29 +02:00
|
|
|
math.cos(c.im) * math.cosh(c.re),
|
|
|
|
math.sin(c.im) * math.sinh(c.re)
|
2019-07-10 19:28:19 +02:00
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Hyperbolic Tangent
|
|
|
|
// Based on
|
|
|
|
// http://www.milefoot.com/math/complex/functionsofi.htm
|
|
|
|
pub fn (c Complex) tanh() Complex {
|
|
|
|
return c.sinh().divide(c.cosh())
|
|
|
|
}
|
|
|
|
|
2019-07-11 15:35:06 +02:00
|
|
|
// Complex Hyperbolic Arc Sin / Sin Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm
|
|
|
|
pub fn (c Complex) asinh() Complex {
|
|
|
|
return c.add(
|
|
|
|
c.pow(2)
|
|
|
|
.add(complex(1,0))
|
|
|
|
.root(2)
|
|
|
|
).ln()
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Hyperbolic Arc Consine / Consine Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm
|
|
|
|
pub fn (c Complex) acosh() Complex {
|
|
|
|
if(c.re > 1) {
|
|
|
|
return c.add(
|
|
|
|
c.pow(2)
|
|
|
|
.subtract(complex(1,0))
|
|
|
|
.root(2)
|
|
|
|
).ln()
|
|
|
|
}
|
|
|
|
else {
|
|
|
|
one := complex(1,0)
|
|
|
|
return c.add(
|
|
|
|
c.add(one)
|
|
|
|
.root(2)
|
|
|
|
.multiply(
|
|
|
|
c.subtract(one)
|
|
|
|
.root(2)
|
|
|
|
)
|
|
|
|
).ln()
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
// Complex Hyperbolic Arc Tangent / Tangent Inverse
|
|
|
|
// Based on
|
|
|
|
// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm
|
|
|
|
pub fn (c Complex) atanh() Complex {
|
|
|
|
if(c.re < 1) {
|
|
|
|
one := complex(1,0)
|
|
|
|
return complex(1.0/2,0).multiply(
|
|
|
|
one
|
|
|
|
.add(c)
|
|
|
|
.divide(
|
|
|
|
one
|
|
|
|
.subtract(c)
|
|
|
|
)
|
|
|
|
.ln()
|
|
|
|
)
|
|
|
|
}
|
|
|
|
else {
|
|
|
|
one := complex(1,0)
|
|
|
|
return complex(1.0/2,0).multiply(
|
|
|
|
one
|
|
|
|
.add(c)
|
|
|
|
.ln()
|
|
|
|
.subtract(
|
|
|
|
one
|
|
|
|
.subtract(c)
|
|
|
|
.ln()
|
|
|
|
)
|
|
|
|
)
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
2019-07-06 20:19:09 +02:00
|
|
|
// Complex Equals
|
|
|
|
pub fn (c1 Complex) equals(c2 Complex) bool {
|
|
|
|
return (c1.re == c2.re) && (c1.im == c2.im)
|
2019-07-07 16:13:17 +02:00
|
|
|
}
|