defmodule Math do
@vsn "0.2.0"
@moduledoc """
Mathematical functions and constants.
"""
# For practical uses floats can be considered equal if their difference is less than this value. See <~>.
@epsilon 1.0e-15
# Theoretical limit is 1.80e308, but Erlang errors at that value, so the practical limit is slightly below that one.
@max_value 1.79_769_313_486_231_580_793e308
@type x :: number
@type y :: number
@doc """
The mathematical constant *π* (pi).
The ratio of a circle's circumference to its diameter.
The returned number is a floating-point approximation (as π is irrational)
"""
@spec pi :: float
defdelegate pi, to: :math
@rad_in_deg 180 / :math.pi()
@doc """
The mathematical constant *τ* (tau).
The ratio of a circle's circumference to its radius.
Defined as 2 * π, but preferred by some mathematicians.
The returned number is a floating-point approximation (as τ is irrational)
"""
@spec tau :: float
def tau, do: pi() * 2
@doc """
The mathematical constant *ℯ* (e).
"""
@spec e :: float
def e, do: 2.718281828459045
@doc """
Equality-test for whether *x* and *y* are _nearly_ equal.
This is useful when working with floating-point numbers, as these introduce small rounding errors.
## Examples
iex> 2.3 - 0.3 == 2.0
false
iex> 2.3 - 0.3 <~> 2.0
true
"""
@spec number <~> number :: boolean
def x <~> y do
abs_x = abs(x)
abs_y = abs(y)
diff = abs(x - y)
# Hacky comparison for floats that are nearly equal.
cond do
x == y ->
true
x == 0 or y == 0 ->
diff < @epsilon
true ->
diff / min(abs_x + abs_y, @max_value) < @epsilon
end
end
# General
@doc """
Arithmetic exponentiation. Calculates *x* to the *n* -th power.
When both *x* and *n* are integers and *n* is positive, returns an `integer`.
When *n* is a negative integer, returns a `float`.
When working with integers, the Exponentiation by Squaring algorithm is used, to allow for a fast and precise result.
When one of the numbers is a float, returns a `float` by using erlang's `:math.pow/2` function.
It is possible to calculate roots by choosing *n* between 0.0 and 1.0 (To calculate the *p* -th-root, pass 1/*p* to the function)
## Examples
iex> Math.pow(2, 4)
16
iex> Math.pow(2.0, 4)
16.0
iex> Math.pow(2, 4.0)
16.0
iex> Math.pow(5, 100)
7888609052210118054117285652827862296732064351090230047702789306640625
iex> Math.pow(5.0, 100)
7.888609052210118e69
iex> Math.pow(2, (1 / 2))
1.4142135623730951
"""
@spec pow(number, number) :: number
def pow(x, n)
def pow(x, n) when is_integer(x) and is_integer(n), do: _pow(x, n)
# Float implementation. Uses erlang's math library.
def pow(x, n) do
:math.pow(x, n)
end
# Integer implementation. Uses Exponentiation by Squaring.
defp _pow(x, n, y \\ 1)
defp _pow(_x, 0, y), do: y
defp _pow(x, 1, y), do: x * y
defp _pow(x, n, y) when n < 0, do: _pow(1 / x, -n, y)
defp _pow(x, n, y) when rem(n, 2) == 0, do: _pow(x * x, div(n, 2), y)
defp _pow(x, n, y), do: _pow(x * x, div(n - 1, 2), x * y)
@doc """
Calculates the non-negative square root of *x*.
"""
@spec sqrt(x) :: float
defdelegate sqrt(x), to: :math
@doc """
Calculates the non-negative nth-root of *x*.
## Examples
iex> Math.nth_root(27, 3)
3.0
iex> Math.nth_root(65536, 8)
4.0
"""
@spec nth_root(x, number) :: float
def nth_root(x, n)
def nth_root(x, n), do: pow(x, 1 / n)
@doc """
Calculates the non-negative integer square root of *x* (rounded towards zero)
Does not accept negative numbers as input.
## Examples
iex> Math.isqrt(100)
10
iex> Math.isqrt(16)
4
iex> Math.isqrt(65536)
256
iex> Math.isqrt(10)
3
"""
@spec isqrt(integer) :: integer
def isqrt(x)
def isqrt(x) when x < 0, do: raise(ArithmeticError)
def isqrt(x), do: _isqrt(x, 1, div(1 + x, 2))
defp _isqrt(x, m, n) when abs(m - n) <= 1 and n * n <= x, do: n
defp _isqrt(_x, m, n) when abs(m - n) <= 1, do: n - 1
defp _isqrt(x, _, n) do
_isqrt(x, n, div(n + div(x, n), 2))
end
#
@doc """
Calculates the Greatest Common divisor of two numbers.
This is the largest positive integer that divides both *a* and *b* without leaving a remainder.
Also see `Math.lcm/2`
## Examples
iex> Math.gcd(2, 4)
2
iex> Math.gcd(2, 3)
1
iex> Math.gcd(12, 8)
4
iex> Math.gcd(54, 24)
6
iex> Math.gcd(-54, 24)
6
"""
@spec gcd(integer, integer) :: non_neg_integer
def gcd(a, b) when is_integer(a) and is_integer(b), do: egcd(a, b) |> elem(0)
@doc """
Calculates integers `gcd`, `s`, and `t` for `as + bt = gcd(a, b)`
Also see `Math.gcd/2`.
Returns a tuple: `{gcd, s, t}`
## Examples
iex> Math.egcd(2, 4)
{2, 1, 0}
iex> Math.egcd(2, 3)
{1, -1, 1}
iex> Math.egcd(12, 8)
{4, 1, -1}
iex> Math.egcd(54, 24)
{6, 1, -2}
iex> Math.egcd(-54, 24)
{6, 1, -2}
"""
@spec egcd(integer, integer) :: non_neg_integer
def egcd(a, b) when is_integer(a) and is_integer(b), do: _egcd(abs(a), abs(b), 0, 1, 1, 0)
defp _egcd(0, b, s, t, _u, _v), do: {b, s, t}
defp _egcd(a, b, s, t, u, v) do
q = div(b, a)
r = rem(b, a)
m = s - u * q
n = t - v * q
_egcd(r, a, u, v, m, n)
end
@doc """
Calculates the Least Common Multiple of two numbers.
This is the smallest positive integer that can be divided by both *a* by *b* without leaving a remainder.
Also see `Math.gcd/2`
## Examples
iex> Math.lcm(4, 6)
12
iex> Math.lcm(3, 7)
21
iex> Math.lcm(21, 6)
42
"""
@spec lcm(integer, integer) :: non_neg_integer
def lcm(a, b)
def lcm(0, 0), do: 0
def lcm(a, b) do
abs(Kernel.div(a * b, gcd(a, b)))
end
@precompute_factorials_up_to 1000
@doc """
Calculates the factorial of *n*: 1 * 2 * 3 * ... * *n*
To make this function faster, values of *n* up to `#{@precompute_factorials_up_to}` are precomputed at compile time.
## Examples
iex> Math.factorial(1)
1
iex> Math.factorial(5)
120
iex> Math.factorial(20)
2432902008176640000
"""
@spec factorial(non_neg_integer) :: pos_integer
def factorial(n)
def factorial(0), do: 1
for {n, fact} <-
1..@precompute_factorials_up_to
|> Enum.scan({0, 1}, fn n, {_prev_n, prev_fact} -> {n, n * prev_fact} end) do
def factorial(unquote(n)), do: unquote(fact)
end
def factorial(n) when n >= 0 do
n * factorial(n - 1)
end
@doc """
Calculates the k-permutations of *n*.
This is the number of distinct ways to create groups of size *k* from *n* distinct elements.
Notice that *n* is the first parameter, for easier piping.
## Examples
iex> Math.k_permutations(10, 2)
90
iex> Math.k_permutations(5, 5)
120
iex> Math.k_permutations(3, 4)
0
"""
@spec k_permutations(non_neg_integer, non_neg_integer) :: non_neg_integer
def k_permutations(n, k)
def k_permutations(n, k) when k > n, do: 0
def k_permutations(n, k) do
div(factorial(n), factorial(n - k))
end
@doc """
Calculates the k-combinations of *n*.
## Examples
iex> Math.k_combinations(10, 2)
45
iex> Math.k_combinations(5, 5)
1
iex> Math.k_combinations(3, 4)
0
"""
@spec k_combinations(non_neg_integer, non_neg_integer) :: non_neg_integer
def k_combinations(n, k)
def k_combinations(n, k) when k > n, do: 0
def k_combinations(n, k) do
div(factorial(n), factorial(k) * factorial(n - k))
end
# Logarithms and exponentiation
@doc """
Calculates ℯ to the xth power.
"""
@spec exp(x) :: float
defdelegate exp(x), to: :math
@doc """
Calculates the natural logarithm (base `ℯ`) of *x*.
See also `Math.e/0`.
"""
@spec log(x) :: float
defdelegate log(x), to: :math
@doc """
Calculates the base-*b* logarithm of *x*
Note that variants for the most common logarithms exist that are faster and more precise.
See also `Math.log/1`, `Math.log2/1` and `Math.log10/1`.
## Examples
iex> Math.log(5, 5)
1.0
iex> Math.log(20, 2) <~> Math.log2(20)
true
iex> Math.log(20, 10) <~> Math.log10(20)
true
iex> Math.log(2, 4)
0.5
iex> Math.log(10, 4)
1.6609640474436813
"""
@spec log(x, number) :: float
def log(x, x), do: 1.0
def log(x, b) do
:math.log(x) / :math.log(b)
end
@doc """
Calculates the binary logarithm (base `2`) of *x*.
See also `Math.log/2`.
"""
@spec log2(x) :: float
defdelegate log2(x), to: :math
@doc """
Calculates the common logarithm (base `10`) of *x*.
See also `Math.log/2`.
"""
@spec log10(x) :: float
defdelegate log10(x), to: :math
# Trigonometry
@doc """
Converts degrees to radians
## Examples
iex>Math.deg2rad(180)
3.141592653589793
"""
@spec deg2rad(x) :: float
def deg2rad(x) do
x / @rad_in_deg
end
@doc """
Converts radians to degrees
## Examples
iex>Math.rad2deg(Math.pi)
180.0
"""
@spec rad2deg(x) :: float
def rad2deg(x) do
x * @rad_in_deg
end
@doc """
Computes the sine of *x*.
*x* is interpreted as a value in radians.
"""
@spec sin(x) :: float
defdelegate sin(x), to: :math
@doc """
Computes the cosine of *x*.
*x* is interpreted as a value in radians.
"""
@spec cos(x) :: float
defdelegate cos(x), to: :math
@doc """
Computes the tangent of *x* (expressed in radians).
"""
@spec tan(x) :: float
defdelegate tan(x), to: :math
@doc """
Computes the arc sine of *x*. (expressed in radians)
"""
@spec asin(x) :: float
defdelegate asin(x), to: :math
@doc """
Computes the arc cosine of *x*. (expressed in radians)
"""
@spec acos(x) :: float
defdelegate acos(x), to: :math
@doc """
Computes the arc tangent of *x*. (expressed in radians)
"""
@spec atan(x) :: float
defdelegate atan(x), to: :math
@doc """
Computes the arc tangent given *y* and *x*. (expressed in radians)
This variant returns the inverse tangent in the correct quadrant, as the signs of both *x* and *y* are known.
"""
@spec atan2(y, x) :: float
defdelegate atan2(y, x), to: :math
# Advanced Trigonometry
@doc """
Computes the hyperbolic sine of *x* (expressed in radians).
"""
@spec sinh(x) :: float
defdelegate sinh(x), to: :math
@doc """
Computes the hyperbolic cosine of *x* (expressed in radians).
"""
@spec cosh(x) :: float
defdelegate cosh(x), to: :math
@doc """
Computes the hyperbolic tangent of *x* (expressed in radians).
"""
@spec tanh(x) :: float
defdelegate tanh(x), to: :math
@doc """
Computes the inverse hyperbolic sine of *x*.
"""
@spec asinh(x) :: float
defdelegate asinh(x), to: :math
@doc """
Computes the inverse hyperbolic cosine of *x*.
"""
@spec acosh(x) :: float
defdelegate acosh(x), to: :math
@doc """
Computes the inverse hyperbolic tangent of *x*.
"""
@spec atanh(x) :: float
defdelegate atanh(x), to: :math
@doc """
Computes the modular multiplicatibe inverse of `a` under modulo `m`
In other words, given integers `a` and `m` calculate a value `b` for `ab = 1 (mod m)`
Returns an `{:ok, b}` tuple or `{:error, :not_coprime}` when `b` cannot be calculated for the inputs.
## Examples
iex> Math.mod_inv 3, 11
{:ok, 4}
iex> Math.mod_inv 10, 17
{:ok, 12}
iex> Math.mod_inv 123, 455
{:ok, 37}
iex> Math.mod_inv 123, 456
{:error, :not_coprime}
"""
@spec mod_inv(a :: integer, m :: integer) :: {:ok, integer} | {:error, :not_coprime}
def mod_inv(a, m)
def mod_inv(a, m) when is_float(a) or is_float(m),
do: raise(ArgumentError, "Inputs cannot be of type float")
def mod_inv(a, m) when is_integer(a) and is_integer(m) do
case egcd(a, m) do
{1, s, _t} -> {:ok, rem(s + m, m)}
_ -> {:error, :not_coprime}
end
end
@doc """
Computes the modular multiplicatibe inverse of `a` under modulo `m`
Similar to `mod_in/2`, but returns only the value or raises an error.
## Examples
iex> Math.mod_inv! 3, 11
4
iex> Math.mod_inv! 10, 17
12
iex> Math.mod_inv! 123, 455
37
iex> Math.mod_inv! 123, 456
** (ArithmeticError) Inputs are not coprime!
"""
@spec mod_inv!(a :: integer, m :: integer) :: integer
def mod_inv!(a, m)
def mod_inv!(a, m) do
case mod_inv(a, m) do
{:ok, b} -> b
_ -> raise ArithmeticError, "Inputs are not coprime!"
end
end
# Interpolation
@doc """
Computes the y value of a given t point on a bezier_curve
If t is not in range [0,1] returns the interpolated point outside control points bounding box
## Examples
iex> Math.bezier_curve(0.5, [{0,0}, {1,1}])
{0.5,0.5}
"""
@spec bezier_curve(t :: number, control_points :: maybe_improper_list) :: tuple
def bezier_curve(t, control_points) when is_list(control_points) and is_number(t) do
new_points =
control_points
|> Enum.with_index()
|> Enum.slice(0..-2)
|> Enum.map(
fn {p0, index} ->
next_index = index + 1
p1 = Enum.at(control_points, next_index)
linear_interpolation(t, p0, p1)
end
)
if Enum.count(new_points) == 1 do
Enum.at(new_points, 0)
else
bezier_curve(t, new_points)
end
end
@doc """
Computes the y value of a given t point on a bezier_curve
Similar to `bezier_curve/2`, but raises an error if not on range [0,1].
## Examples
iex> Math.bezier_curve!(0.5, [{0,0}, {1,1}])
{0.5,0.5}
iex> Math.bezier_curve!(1.5, [{0,0}, {1,1}])
** (ArgumentError) t is not beetween 0 and 1
"""
@spec bezier_curve!(t :: number, control_points :: maybe_improper_list) :: tuple
def bezier_curve!(t, control_points)
when is_number(t)
and 0.0 <= t
and t <= 1.0
and is_list(control_points) do
bezier_curve(t, control_points)
end
def bezier_curve!(t, control_points)
when is_number(t)
and is_list(control_points) do
raise ArgumentError, "t is not beetween 0 and 1"
end
@doc """
Computes the y value of a given t point on a linear_interpolation between 2 points
If t is not in range [0,1] returns the interpolated point outside control points bounding box
## Examples
iex> Math.linear_interpolation(0.5, {0,0}, {1,1})
{0.5, 0.5}
"""
@spec linear_interpolation(t :: number, p0 :: tuple, p1 :: tuple) :: tuple
def linear_interpolation(t, p0, p1)
when is_number(t)
and is_tuple(p0)
and is_tuple(p1) do
last = tuple_size(p0) - 1
0..last
|> Enum.reduce({},
fn index, tuple ->
first_elem = p0 |> elem(index)
second_elem = p1 |> elem(index)
tuple
|> Tuple.append(
(1 - t) * first_elem + t * second_elem
)
end
)
end
@doc """
Computes the point of a given t on a linear_interpolation between 2 points
Similar to `linear_interpolation/3`, but raises an error if not on range [0,1].
## Examples
iex> Math.linear_interpolation!(0.5, {0,0}, {1,1})
{0.5,0.5}
iex> Math.linear_interpolation!(1.5, {0,0}, {1,1})
** (ArgumentError) t is not beetween 0 and 1
"""
@spec linear_interpolation!(t :: number, p0 :: tuple, p1 :: tuple) :: tuple
def linear_interpolation!(t, p0, p1)
when is_number(t)
and 0.0 <= t
and t <= 1.0
and is_tuple(p0)
and is_tuple(p1) do
linear_interpolation(t, p0, p1)
end
def linear_interpolation!(t, p0, p1)
when is_number(t)
and is_tuple(p0)
and is_tuple(p1) do
raise ArgumentError, "t is not beetween 0 and 1"
end
end
