feat: Add first draft of complex number support

This commit is contained in:
sigmasternchen 2024-12-20 19:44:58 +01:00
parent 46d2309268
commit aa9eae6177
2 changed files with 790 additions and 0 deletions

View file

@ -0,0 +1,251 @@
import gleam/bool
import gleam/float
import gleam/int
import gleam/io
import gleam/list
import gleam/order.{type Order, Eq, Gt, Lt}
import gleam/result
import gleam_community/maths
pub type Complex {
Complex(real: Float, imaginary: Float)
}
pub fn add(a: Complex, b: Complex) -> Complex {
Complex(a.real +. b.real, a.imaginary +. b.imaginary)
}
pub fn subtract(a: Complex, b: Complex) -> Complex {
Complex(a.real -. b.real, a.imaginary -. b.imaginary)
}
pub fn multiply(a: Complex, b: Complex) -> Complex {
Complex(
a.real *. b.real -. a.imaginary *. b.imaginary,
a.real *. b.imaginary +. a.imaginary *. b.real,
)
}
pub fn divide(a: Complex, b: Complex) -> Complex {
Complex(
{ a.real *. b.real +. a.imaginary *. b.imaginary }
/. { b.real *. b.real +. b.imaginary *. b.imaginary },
{ a.imaginary *. b.real -. a.real *. b.imaginary }
/. { b.real *. b.real +. b.imaginary *. b.imaginary },
)
}
pub fn absolute_value(z: Complex) -> Float {
let assert Ok(result) =
float.square_root(z.real *. z.real +. z.imaginary *. z.imaginary)
result
}
pub fn argument(z: Complex) -> Result(Float, Nil) {
case float.compare(z.real, 0.0), float.compare(z.imaginary, 0.0) {
Eq, Eq -> Error(Nil)
Eq, Lt -> Ok(0.0 -. maths.pi() /. 2.0)
Eq, Gt -> Ok(maths.pi() /. 2.0)
Lt, Lt -> Ok(maths.atan(z.imaginary /. z.real) -. maths.pi())
Lt, _ -> Ok(maths.atan(z.imaginary /. z.real) +. maths.pi())
Gt, _ -> Ok(maths.atan(z.imaginary /. z.real))
}
}
pub fn from_float(real: Float) -> Complex {
Complex(real, 0.0)
}
fn from_unit_angle(phi: Float) -> Complex {
Complex(maths.cos(phi), maths.sin(phi))
}
pub fn from_polar(r: Float, phi: Float) -> Complex {
multiply(from_float(r), from_unit_angle(phi))
}
pub fn to_polar(z: Complex) -> Result(#(Float, Float), Nil) {
case argument(z) {
Error(_) -> Error(Nil)
Ok(arg) -> Ok(#(absolute_value(z), arg))
}
}
pub fn exponential(z: Complex) -> Complex {
from_polar(maths.exponential(z.real), z.imaginary)
}
pub fn reciprocal(z: Complex) -> Complex {
let divisor = z.real *. z.real +. z.imaginary *. z.imaginary
Complex(z.real /. divisor, 0.0 -. z.imaginary /. divisor)
}
// De Moivres Theorem
pub fn power(z: Complex, n: Int) -> Result(Complex, Nil) {
let r = absolute_value(z)
case argument(z), int.compare(n, 0) {
// 0 ^ 0 -> undefined
Error(_), Eq -> Error(Nil)
// 0 ^ n -> 0
Error(_), _ -> Ok(Complex(0.0, 0.0))
// z ^ 0 -> 1
Ok(_), Eq -> Ok(multiplicative_identify())
// De Moivre's Theorem only works for positive integers
Ok(_), Lt -> Error(Nil)
Ok(arg), _ -> {
// n can't be 0
let assert Ok(r_to_the_n) = float.power(r, int.to_float(n))
Ok(from_polar(r_to_the_n, int.to_float(n) *. arg))
}
}
}
// De Moivres Theorem
pub fn nth_root(z: Complex, n: Int) -> Result(List(Complex), Nil) {
case int.compare(n, 0) {
Eq -> Error(Nil)
Lt -> {
let assert Ok(result_before_reciprocal) = nth_root(z, -n)
Ok(
result_before_reciprocal
|> list.map(reciprocal),
)
}
Gt ->
case argument(z) {
// any root of 0 = 0
Error(_) -> Ok([Complex(0.0, 0.0)])
Ok(arg) -> {
let r = absolute_value(z)
// r and n are always positive -> root is defined
let assert Ok(new_r) = maths.nth_root(r, n)
list.range(0, n - 1)
|> list.map(fn(k) {
from_polar(
new_r,
{ arg +. 2.0 *. int.to_float(k) *. maths.pi() } /. int.to_float(n),
)
})
|> Ok
}
}
}
}
fn zero() -> Complex {
Complex(0.0, 0.0)
}
fn multiplicative_identify() -> Complex {
from_float(1.0)
}
pub fn sum(arr: List(Complex)) -> Complex {
list.fold(arr, zero(), add)
}
pub fn product(arr: List(Complex)) -> Complex {
list.fold(arr, multiplicative_identify(), multiply)
}
pub fn weighted_sum(arr: List(#(Complex, Float))) -> Result(Complex, Nil) {
let weight_is_negative = list.any(arr, fn(tuple) { tuple.1 <. 0.0 })
case weight_is_negative {
True -> Error(Nil)
False ->
Ok(
list.fold(arr, zero(), fn(acc, tuple) {
add(acc, multiply(tuple.0, from_float(tuple.1)))
}),
)
}
}
pub fn weighted_product(arr: List(#(Complex, Int))) -> Result(Complex, Nil) {
let weight_is_negative = list.any(arr, fn(tuple) { tuple.1 < 0 })
case weight_is_negative {
True -> Error(Nil)
False ->
list.fold(arr, Ok(multiplicative_identify()), fn(acc_result, tuple) {
acc_result
|> result.then(fn(acc) {
power(tuple.0, tuple.1)
|> result.map(multiply(acc, _))
})
})
}
}
pub fn cumulative_sum(arr: List(Complex)) -> List(Complex) {
list.scan(arr, zero(), add)
}
pub fn cumulative_product(arr: List(Complex)) -> List(Complex) {
list.scan(arr, multiplicative_identify(), multiply)
}
pub fn absolute_difference(a: Complex, b: Complex) {
absolute_value(subtract(a, b))
}
pub fn mean(arr: List(Complex)) -> Result(Complex, Nil) {
case arr {
[] -> Error(Nil)
_ ->
sum(arr)
|> divide(arr |> list.length |> int.to_float |> from_float)
|> Ok
}
}
pub fn median(arr: List(Complex)) -> Result(Complex, Nil) {
use <- bool.guard(list.is_empty(arr), Error(Nil))
let length = list.length(arr)
let mid = length / 2
case length % 2 == 0 {
False -> arr |> list.drop(mid) |> list.first
True -> {
arr
|> list.drop(mid - 1)
|> list.take(2)
|> mean
}
}
}
pub fn is_close(x: Complex, y: Complex, rtol: Float, atol: Float) -> Bool {
let x = absolute_difference(x, y)
let y = atol +. rtol *. absolute_value(y)
x <=. y
}
pub fn all_close(
arr: List(#(Complex, Complex)),
rtol: Float,
atol: Float,
) -> Result(List(Bool), Nil) {
Ok(list.map(arr, fn(tuple) { is_close(tuple.0, tuple.1, rtol, atol) }))
}
pub fn compare_real(a: Complex, b: Complex) -> Order {
float.compare(a.real, b.real)
}
pub fn compare_imaginary(a: Complex, b: Complex) -> Order {
float.compare(a.imaginary, b.imaginary)
}
pub fn conjugate(z: Complex) -> Complex {
Complex(z.real, 0.0 -. z.imaginary)
}
pub fn to_string(z: Complex) -> String {
float.to_string(z.real)
<> case float.compare(z.imaginary, 0.0) {
Eq -> ""
Gt -> " + " <> float.to_string(z.imaginary) <> "i"
Lt -> " - " <> float.to_string(0.0 -. z.imaginary) <> "i"
}
}

View file

@ -0,0 +1,539 @@
import gleam/float
import gleam/function
import gleam/io
import gleam/list
import gleam/order.{Eq, Gt, Lt}
import gleam_community/complex.{Complex}
import gleam_community/maths
import gleeunit/should
pub fn is_close_relative_positive_test() {
complex.is_close(Complex(0.0, 0.9), Complex(0.0, 1.0), 0.1, 0.0)
|> should.be_true
}
pub fn is_close_relative_negative_test() {
complex.is_close(Complex(0.0, 0.9), Complex(0.0, 1.0), 0.09, 0.0)
|> should.be_false
}
pub fn is_close_absolute_positive_test() {
complex.is_close(Complex(0.0, 0.9), Complex(0.0, 1.0), 0.0, 0.1)
|> should.be_true
}
pub fn is_close_absolute_negative_test() {
complex.is_close(Complex(0.0, 0.9), Complex(0.0, 1.0), 0.0, 0.09)
|> should.be_false
}
pub fn add_test() {
complex.add(Complex(1.0, 2.0), Complex(3.0, 4.0))
|> should.equal(Complex(4.0, 6.0))
}
pub fn multiply_test() {
complex.multiply(Complex(1.0, 2.0), Complex(3.0, 4.0))
|> should.equal(Complex(-5.0, 10.0))
}
pub fn subtract_test() {
complex.subtract(Complex(1.0, 2.0), Complex(3.0, 4.0))
|> should.equal(Complex(-2.0, -2.0))
}
pub fn divide_test() {
complex.divide(Complex(1.0, 2.0), Complex(3.0, 4.0))
|> complex.is_close(Complex(0.44, 0.08), 0.0, 0.01)
|> should.be_true
}
pub fn absolute_value_test() {
complex.absolute_value(Complex(-4.0, 6.0))
|> float.loosely_equals(7.21, 0.05)
|> should.be_true
}
pub fn argument_origin_test() {
complex.argument(Complex(0.0, 0.0))
|> should.be_error
}
pub fn argument_real_positive_test() {
complex.argument(Complex(1.0, 0.0))
|> should.be_ok
|> float.loosely_equals(0.0, 0.01)
|> should.be_true
}
pub fn argument_real_negative_test() {
complex.argument(Complex(-1.0, 0.0))
|> should.be_ok
|> float.loosely_equals(maths.pi(), 0.01)
|> should.be_true
}
pub fn argument_imaginary_positive_test() {
complex.argument(Complex(0.0, 1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. 0.5, 0.01)
|> should.be_true
}
pub fn argument_imaginary_negative_test() {
complex.argument(Complex(0.0, -1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. -0.5, 0.01)
|> should.be_true
}
pub fn argument_1st_quadrant_test() {
complex.argument(Complex(1.0, 1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. 0.25, 0.01)
|> should.be_true
}
pub fn argument_2nd_quadrant_test() {
complex.argument(Complex(-1.0, 1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. 0.75, 0.01)
|> should.be_true
}
pub fn argument_3rd_quadrant_test() {
complex.argument(Complex(-1.0, -1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. -0.75, 0.01)
|> should.be_true
}
pub fn argument_4th_quadrant_test() {
complex.argument(Complex(1.0, -1.0))
|> should.be_ok
|> float.loosely_equals(maths.pi() *. -0.25, 0.01)
|> should.be_true
}
pub fn from_float_test() {
complex.from_float(69.42)
|> complex.is_close(Complex(69.42, 0.0), 0.0, 0.01)
|> should.be_true
}
pub fn from_polar_test() {
complex.from_polar(50.0, 4.0)
|> complex.is_close(Complex(-32.68, -37.84), 0.0, 0.01)
|> should.be_true
}
pub fn to_polar_origin_test() {
complex.to_polar(Complex(0.0, 0.0))
|> should.be_error
}
pub fn to_polar_default_test() {
let #(r, phi) =
complex.to_polar(Complex(123.0, -321.0))
|> should.be_ok
r
|> float.loosely_equals(343.76, 0.01)
|> should.be_true
phi
|> float.loosely_equals(-1.2, 0.01)
|> should.be_true
}
pub fn exponential_test() {
complex.exponential(Complex(-0.1, 10.0))
|> complex.is_close(Complex(-0.759, -0.492), 0.0, 0.01)
|> should.be_true
}
pub fn reciprocal_test() {
complex.reciprocal(Complex(0.19, -0.7))
|> complex.is_close(Complex(0.361, 1.33), 0.0, 0.1)
|> should.be_true
}
pub fn power_undefined_test() {
complex.power(Complex(0.0, 0.0), 0)
|> should.be_error
}
pub fn power_of_zero_test() {
complex.power(Complex(0.0, 0.0), 4)
|> should.be_ok
|> should.equal(Complex(0.0, 0.0))
}
pub fn power_to_the_zeroth_test() {
complex.power(Complex(42.0, 42.0), 0)
|> should.be_ok
|> should.equal(Complex(1.0, 0.0))
}
pub fn power_to_a_negative_test() {
complex.power(Complex(0.13, 0.65), -2)
|> should.be_error
}
pub fn power_default_test() {
complex.power(Complex(-1.5, 1.2), 10)
|> should.be_ok
|> complex.is_close(Complex(611.72, -306.3), 0.0, 0.1)
|> should.be_true
}
pub fn nth_root_zeroth_test() {
complex.nth_root(Complex(13.4, -16.4), 0)
|> should.be_error
}
pub fn nth_root_default_positive_test() {
complex.nth_root(Complex(13.4, -16.4), 5)
|> should.be_ok
|> list.zip([
Complex(1.81, -0.32),
Complex(0.868, 1.62),
Complex(-1.275, 1.328),
Complex(-1.657, -0.8),
Complex(0.25, -1.82),
])
|> list.all(fn(tuple) {
let #(actual, expected) = tuple
complex.is_close(actual, expected, 0.0, 0.01)
})
|> should.be_true
}
pub fn nth_root_default_negative_test() {
complex.nth_root(Complex(13.4, -16.4), -3)
|> should.be_ok
|> list.zip([
Complex(0.346, 0.105),
Complex(-0.08, -0.35),
Complex(-0.26, 0.247),
])
|> list.all(fn(tuple) {
let #(actual, expected) = tuple
complex.is_close(actual, expected, 0.0, 0.01)
})
|> should.be_true
}
pub fn sum_default_test() {
complex.sum([
Complex(1.0, 2.0),
Complex(-1.0, -1.0),
Complex(5.0, 4.0),
Complex(-5.0, -6.0),
])
|> complex.is_close(Complex(0.0, -1.0), 0.0, 0.01)
|> should.be_true
}
pub fn sum_empty_test() {
complex.sum([])
|> should.equal(Complex(0.0, 0.0))
}
pub fn product_default_test() {
complex.product([
Complex(1.0, 2.0),
Complex(-1.0, -1.0),
Complex(5.0, 4.0),
Complex(-5.0, -6.0),
])
|> complex.is_close(Complex(-151.0, -47.0), 0.0, 0.1)
|> should.be_true
}
pub fn product_empty_test() {
complex.product([])
|> should.equal(Complex(1.0, 0.0))
}
pub fn weighted_sum_default_test() {
complex.weighted_sum([
#(Complex(1.0, 2.0), 1.0),
#(Complex(-1.0, -1.0), 2.0),
#(Complex(5.0, 4.0), 2.0),
#(Complex(-5.0, -6.0), 0.5),
])
|> should.be_ok
|> complex.is_close(Complex(6.5, 5.0), 0.0, 0.01)
|> should.be_true
}
pub fn weighted_sum_empty_test() {
complex.weighted_sum([])
|> should.be_ok
|> should.equal(Complex(0.0, 0.0))
}
pub fn weighted_sum_negative_weights_test() {
complex.weighted_sum([
#(Complex(1.0, 2.0), -1.0),
#(Complex(-1.0, -1.0), 2.0),
#(Complex(5.0, 4.0), 2.0),
#(Complex(-5.0, -6.0), 0.5),
])
|> should.be_error
}
pub fn weighted_product_default_test() {
complex.weighted_product([
#(Complex(1.0, 2.0), 1),
#(Complex(-1.0, -1.0), 2),
#(Complex(5.0, 4.0), 2),
#(Complex(-5.0, -6.0), 1),
])
|> should.be_ok
|> complex.is_close(Complex(-272.0, 1406.0), 0.0, 0.01)
|> should.be_true
}
pub fn weighted_product_empty_test() {
complex.weighted_product([])
|> should.be_ok
|> should.equal(Complex(1.0, 0.0))
}
pub fn weighted_product_negative_weights_test() {
complex.weighted_product([
#(Complex(1.0, 2.0), -1),
#(Complex(-1.0, -1.0), 2),
#(Complex(5.0, 4.0), 2),
#(Complex(-5.0, -6.0), 1),
])
|> should.be_error
}
pub fn cumulative_sum_default_test() {
complex.cumulative_sum([
Complex(1.0, 2.0),
Complex(-1.0, -1.0),
Complex(5.0, 4.0),
Complex(-5.0, -6.0),
])
|> should.equal([
Complex(1.0, 2.0),
Complex(0.0, 1.0),
Complex(5.0, 5.0),
Complex(0.0, -1.0),
])
}
pub fn cumulative_sum_empty_test() {
complex.cumulative_sum([])
|> should.equal([])
}
pub fn cumulative_product_default_test() {
complex.cumulative_product([
Complex(1.0, 2.0),
Complex(-1.0, -1.0),
Complex(5.0, 4.0),
Complex(-5.0, -6.0),
])
|> should.equal([
Complex(1.0, 2.0),
Complex(1.0, -3.0),
Complex(17.0, -11.0),
Complex(-151.0, -47.0),
])
}
pub fn cumulative_product_empty_test() {
complex.cumulative_product([])
|> should.equal([])
}
pub fn mean_empty_test() {
complex.mean([])
|> should.be_error
}
pub fn mean_with_one_element_test() {
complex.mean([Complex(42.0, 1337.0)])
|> should.be_ok
|> should.equal(Complex(42.0, 1337.0))
}
pub fn mean_default_test() {
complex.mean([
Complex(1.0, 1.0),
Complex(3.0, 5.0),
Complex(-2.0, -3.0),
Complex(-3.0, -2.0),
])
|> should.be_ok
|> should.equal(Complex(-0.25, 0.25))
}
pub fn median_empty_test() {
complex.median([])
|> should.be_error
}
pub fn median_with_one_element_test() {
complex.mean([Complex(42.0, 1337.0)])
|> should.be_ok
|> should.equal(Complex(42.0, 1337.0))
}
pub fn median_with_two_elements_test() {
complex.median([Complex(42.0, 10.0), Complex(-42.0, 20.0)])
|> should.be_ok
|> should.equal(Complex(0.0, 15.0))
}
pub fn median_with_three_elements_test() {
complex.median([Complex(42.0, 10.0), Complex(-42.0, 20.0), Complex(0.0, 0.0)])
|> should.be_ok
|> should.equal(Complex(-42.0, 20.0))
}
pub fn median_default_odd_test() {
complex.median([
Complex(1.0, 1.0),
Complex(3.0, 5.0),
Complex(-2.0, -3.0),
Complex(-3.0, -2.0),
Complex(6.0, 9.0),
Complex(1.0, 1.0),
Complex(3.0, 5.0),
Complex(-2.0, -3.0),
Complex(-3.0, -2.0),
])
|> should.be_ok
|> should.equal(Complex(6.0, 9.0))
}
pub fn median_default_even_test() {
complex.median([
Complex(1.0, 1.0),
Complex(3.0, 5.0),
Complex(-2.0, -3.0),
Complex(-3.0, -2.0),
Complex(1.0, 1.0),
Complex(3.0, 5.0),
Complex(-2.0, -3.0),
Complex(-3.0, -2.0),
])
|> should.be_ok
|> should.equal(Complex(-1.0, -0.5))
}
pub fn all_close_absolute_positive_test() {
complex.all_close(
[
#(Complex(1.0, 1.0), Complex(1.1, 1.1)),
#(Complex(-4.0, -5.0), Complex(-4.0, -4.9)),
#(Complex(-3.0, 2.0), Complex(-2.9, 2.1)),
],
0.0,
0.2,
)
|> should.be_ok
|> list.all(function.identity)
|> should.be_true
}
pub fn all_close_absolute_negative_test() {
complex.all_close(
[
#(Complex(1.0, 1.0), Complex(1.1, 1.1)),
#(Complex(-4.0, -5.0), Complex(-4.0, -4.7)),
#(Complex(-3.0, 2.0), Complex(-2.9, 2.1)),
],
0.0,
0.2,
)
|> should.be_ok
|> list.all(function.identity)
|> should.be_false
}
pub fn all_close_relative_positive_test() {
complex.all_close(
[
#(Complex(1.0, 1.0), Complex(1.1, 1.1)),
#(Complex(-4.0, -5.0), Complex(-4.0, -4.9)),
#(Complex(-3.0, 2.0), Complex(-2.9, 2.1)),
],
0.1,
0.0,
)
|> should.be_ok
|> list.all(function.identity)
|> should.be_true
}
pub fn all_close_relative_negative_test() {
complex.all_close(
[
#(Complex(1.0, 1.0), Complex(1.1, 1.1)),
#(Complex(-4.0, -5.0), Complex(-4.0, -4.7)),
#(Complex(-3.0, 2.0), Complex(-2.9, 2.1)),
],
0.05,
0.0,
)
|> should.be_ok
|> list.all(function.identity)
|> should.be_false
}
pub fn compare_real_eq_test() {
complex.compare_real(Complex(5.0, 0.0), Complex(5.0, 42.0))
|> should.equal(Eq)
}
pub fn compare_real_lt_test() {
complex.compare_real(Complex(4.9, 0.0), Complex(5.0, 42.0))
|> should.equal(Lt)
}
pub fn compare_real_gt_test() {
complex.compare_real(Complex(5.1, 0.0), Complex(5.0, 42.0))
|> should.equal(Gt)
}
pub fn compare_imaginary_eq_test() {
complex.compare_imaginary(Complex(1.0, 42.0), Complex(5.0, 42.0))
|> should.equal(Eq)
}
pub fn compare_imaginary_lt_test() {
complex.compare_imaginary(Complex(1.0, 41.9), Complex(5.0, 42.0))
|> should.equal(Lt)
}
pub fn compare_imaginary_gt_test() {
complex.compare_imaginary(Complex(1.0, 42.1), Complex(5.0, 42.0))
|> should.equal(Gt)
}
pub fn conjugate_test() {
complex.conjugate(Complex(42.0, 69.0))
|> should.equal(Complex(42.0, -69.0))
}
pub fn to_string_only_real_test() {
complex.to_string(Complex(-42.0, 0.0))
|> should.equal("-42.0")
}
pub fn to_string_positive_imaginary_test() {
complex.to_string(Complex(4.0, 5.6))
|> should.equal("4.0 + 5.6i")
}
pub fn to_string_negative_imaginary_test() {
complex.to_string(Complex(0.0, -4.2))
|> should.equal("0.0 - 4.2i")
}