2019-06-22 21:20:28 +03: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.
|
|
|
|
|
|
|
|
module math
|
|
|
|
|
2019-07-21 14:29:32 +03:00
|
|
|
#include <math.h>
|
2019-07-12 08:01:12 +03:00
|
|
|
|
|
|
|
// NOTE
|
|
|
|
// When adding a new function, please make sure it's in the right place.
|
2019-07-17 01:03:51 +03:00
|
|
|
// All functions are sorted alphabetically.
|
2019-07-12 08:01:12 +03:00
|
|
|
|
2019-06-22 21:20:28 +03:00
|
|
|
const (
|
2019-06-25 09:31:35 +03:00
|
|
|
E = 2.71828182845904523536028747135266249775724709369995957496696763
|
|
|
|
Pi = 3.14159265358979323846264338327950288419716939937510582097494459
|
|
|
|
Phi = 1.61803398874989484820458683436563811772030917980576286213544862
|
2019-06-27 20:16:02 +03:00
|
|
|
Tau = 6.28318530717958647692528676655900576839433879875021164194988918
|
2019-06-25 09:31:35 +03:00
|
|
|
|
|
|
|
Sqrt2 = 1.41421356237309504880168872420969807856967187537694807317667974
|
|
|
|
SqrtE = 1.64872127070012814684865078781416357165377610071014801157507931
|
|
|
|
SqrtPi = 1.77245385090551602729816748334114518279754945612238712821380779
|
2019-07-02 13:50:33 +03:00
|
|
|
SqrtTau = 2.50662827463100050241576528481104525300698674060993831662992357
|
2019-06-25 09:31:35 +03:00
|
|
|
SqrtPhi = 1.27201964951406896425242246173749149171560804184009624861664038
|
|
|
|
|
|
|
|
Ln2 = 0.693147180559945309417232121458176568075500134360255254120680009
|
|
|
|
Log2E = 1.0 / Ln2
|
|
|
|
Ln10 = 2.30258509299404568401799145468436420760110148862877297603332790
|
|
|
|
Log10E = 1.0 / Ln10
|
2019-06-22 21:20:28 +03:00
|
|
|
)
|
|
|
|
|
2019-07-10 17:05:39 +03:00
|
|
|
const (
|
|
|
|
MaxI8 = (1<<7) - 1
|
|
|
|
MinI8 = -1 << 7
|
|
|
|
MaxI16 = (1<<15) - 1
|
|
|
|
MinI16 = -1 << 15
|
|
|
|
MaxI32 = (1<<31) - 1
|
|
|
|
MinI32 = -1 << 31
|
2019-07-17 01:03:51 +03:00
|
|
|
// MaxI64 = ((1<<63) - 1)
|
|
|
|
// MinI64 = (-(1 << 63) )
|
|
|
|
MaxU8 = (1<<8) - 1
|
2019-07-10 17:05:39 +03:00
|
|
|
MaxU16 = (1<<16) - 1
|
|
|
|
MaxU32 = (1<<32) - 1
|
|
|
|
MaxU64 = (1<<64) - 1
|
2019-07-17 01:03:51 +03:00
|
|
|
)
|
2019-07-10 17:05:39 +03:00
|
|
|
|
2019-06-30 16:33:37 +03:00
|
|
|
// Returns the absolute value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn abs(a f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
if a < 0 {
|
|
|
|
return -a
|
|
|
|
}
|
|
|
|
return a
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// acos calculates inversed cosine (arccosine).
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn acos(a f64) f64 {
|
2019-06-24 17:05:30 +03:00
|
|
|
return C.acos(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// asin calculates inversed sine (arcsine).
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn asin(a f64) f64 {
|
2019-06-24 17:05:30 +03:00
|
|
|
return C.asin(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// atan calculates inversed tangent (arctangent).
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn atan(a f64) f64 {
|
2019-06-24 17:05:30 +03:00
|
|
|
return C.atan(a)
|
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// atan2 calculates inversed tangent with two arguments, returns angle between the X axis and the point.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn atan2(a, b f64) f64 {
|
2019-06-24 17:05:30 +03:00
|
|
|
return C.atan2(a, b)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// cbrt calculates cubic root.
|
2019-06-30 16:33:37 +03:00
|
|
|
pub fn cbrt(a f64) f64 {
|
2019-07-02 13:50:33 +03:00
|
|
|
return C.cbrt(a)
|
2019-06-30 16:33:37 +03:00
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// ceil returns the nearest integer greater or equal to the provided value.
|
2019-07-17 01:03:51 +03:00
|
|
|
pub fn ceil(a f64) int {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.ceil(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// cos calculates cosine.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn cos(a f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
return C.cos(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// cosh calculates hyperbolic cosine.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn cosh(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.cosh(a)
|
|
|
|
}
|
|
|
|
|
2019-07-12 08:46:40 +03:00
|
|
|
// degrees convert from degrees to radians.
|
|
|
|
pub fn degrees(radians f64) f64 {
|
|
|
|
return radians * (180.0 / Pi)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// exp calculates exponement of the number (math.pow(math.E, a)).
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn exp(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.exp(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// digits returns an array of the digits of n in the given base.
|
|
|
|
pub fn digits(n, base int) []int {
|
|
|
|
mut sign := 1
|
|
|
|
if n < 0 {
|
|
|
|
sign = -1
|
|
|
|
n = -n
|
|
|
|
}
|
|
|
|
mut res := []int
|
|
|
|
for n != 0 {
|
|
|
|
res << (n % base) * sign
|
|
|
|
n /= base
|
|
|
|
}
|
|
|
|
return res
|
|
|
|
}
|
|
|
|
|
2019-07-12 08:46:40 +03:00
|
|
|
// erf computes the error funtion value
|
|
|
|
pub fn erf(a f64) f64 {
|
|
|
|
return C.erf(a)
|
|
|
|
}
|
|
|
|
|
|
|
|
// erfc computes the complimentary error function value
|
|
|
|
pub fn erfc(a f64) f64 {
|
|
|
|
return C.erfc(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// exp2 returns the base-2 exponential function of a (math.pow(2, a)).
|
2019-06-30 16:33:37 +03:00
|
|
|
pub fn exp2(a f64) f64 {
|
2019-07-02 13:50:33 +03:00
|
|
|
return C.exp2(a)
|
2019-06-30 16:33:37 +03:00
|
|
|
}
|
|
|
|
|
2019-07-12 08:46:40 +03:00
|
|
|
// factorial calculates the factorial of the provided value.
|
2019-07-17 01:03:51 +03:00
|
|
|
fn recursive_product( n int, current_number_ptr &int) int{
|
|
|
|
mut m := n / 2
|
|
|
|
if (m == 0){
|
|
|
|
return *current_number_ptr += 2
|
|
|
|
}
|
|
|
|
if (n == 2){
|
|
|
|
return (*current_number_ptr += 2) * (*current_number_ptr += 2)
|
|
|
|
}
|
|
|
|
return recursive_product((n - m), *current_number_ptr) * recursive_product(m, *current_number_ptr)
|
|
|
|
}
|
|
|
|
|
|
|
|
pub fn factorial(n int) i64 {
|
|
|
|
if n < 0 {
|
|
|
|
panic('factorial: Cannot find factorial of negative number')
|
|
|
|
}
|
|
|
|
if n < 2 {
|
|
|
|
return i64(1)
|
|
|
|
}
|
|
|
|
mut r := 1
|
|
|
|
mut p := 1
|
|
|
|
mut current_number := 1
|
|
|
|
mut h := 0
|
|
|
|
mut shift := 0
|
|
|
|
mut high := 1
|
|
|
|
mut len := high
|
|
|
|
mut log2n := int(floor(log2(n)))
|
|
|
|
for ;h != n; {
|
|
|
|
shift += h
|
|
|
|
h = n >> log2n
|
|
|
|
log2n -= 1
|
|
|
|
len = high
|
|
|
|
high = (h - 1) | 1
|
|
|
|
len = (high - len)/2
|
|
|
|
if (len > 0){
|
|
|
|
p *= recursive_product(len, ¤t_number)
|
|
|
|
r *= p
|
|
|
|
}
|
|
|
|
}
|
|
|
|
return i64((r << shift))
|
2019-07-12 08:46:40 +03:00
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// floor returns the nearest integer lower or equal of the provided value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn floor(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.floor(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// fmod returns the floating-point remainder of number / denom (rounded towards zero):
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn fmod(a, b f64) f64 {
|
2019-06-24 14:08:38 +03:00
|
|
|
return C.fmod(a, b)
|
2019-06-24 11:50:45 +03:00
|
|
|
}
|
|
|
|
|
2019-07-12 08:46:40 +03:00
|
|
|
// gamma computes the gamma function value
|
|
|
|
pub fn gamma(a f64) f64 {
|
|
|
|
return C.tgamma(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// gcd calculates greatest common (positive) divisor (or zero if a and b are both zero).
|
2019-07-03 19:51:03 +03:00
|
|
|
pub fn gcd(a, b i64) i64 {
|
2019-06-29 18:24:55 +03:00
|
|
|
if a < 0 {
|
|
|
|
a = -a
|
|
|
|
}
|
|
|
|
if b < 0 {
|
|
|
|
b = -b
|
|
|
|
}
|
|
|
|
for b != 0 {
|
|
|
|
a %= b
|
|
|
|
if a == 0 {
|
|
|
|
return b
|
|
|
|
}
|
|
|
|
b %= a
|
|
|
|
}
|
|
|
|
return a
|
|
|
|
}
|
|
|
|
|
|
|
|
// lcm calculates least common (non-negative) multiple.
|
2019-07-03 19:51:03 +03:00
|
|
|
pub fn lcm(a, b i64) i64 {
|
2019-06-29 18:24:55 +03:00
|
|
|
if a == 0 {
|
|
|
|
return a
|
|
|
|
}
|
|
|
|
res := a * (b / gcd(b, a))
|
|
|
|
if res < 0 {
|
|
|
|
return -res
|
|
|
|
}
|
|
|
|
return res
|
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// log calculates natural (base-e) logarithm of the provided value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn log(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.log(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// log2 calculates base-2 logarithm of the provided value.
|
2019-06-30 16:33:37 +03:00
|
|
|
pub fn log2(a f64) f64 {
|
2019-07-11 09:26:31 +03:00
|
|
|
return C.log2(a)
|
2019-06-30 16:33:37 +03:00
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// log10 calculates the common (base-10) logarithm of the provided value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn log10(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.log10(a)
|
|
|
|
}
|
|
|
|
|
2019-07-12 08:46:40 +03:00
|
|
|
// log_gamma computes the log-gamma function value
|
|
|
|
pub fn log_gamma(a f64) f64 {
|
|
|
|
return C.lgamma(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// log_n calculates base-N logarithm of the provided value.
|
2019-06-30 16:33:37 +03:00
|
|
|
pub fn log_n(a, b f64) f64 {
|
|
|
|
return C.log(a) / C.log(b)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// max returns the maximum value of the two provided.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn max(a, b f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
if a > b {
|
|
|
|
return a
|
|
|
|
}
|
|
|
|
return b
|
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// min returns the minimum value of the two provided.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn min(a, b f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
if a < b {
|
|
|
|
return a
|
|
|
|
}
|
|
|
|
return b
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// pow returns base raised to the provided power.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn pow(a, b f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
return C.pow(a, b)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// radians convert from radians to degrees.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn radians(degrees f64) f64 {
|
2019-06-25 09:37:23 +03:00
|
|
|
return degrees * (Pi / 180.0)
|
2019-06-22 21:20:28 +03:00
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// round returns the integer nearest to the provided value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn round(f f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
return C.round(f)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// sin calculates sine.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn sin(a f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
return C.sin(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// sinh calculates hyperbolic sine.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn sinh(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.sinh(a)
|
|
|
|
}
|
|
|
|
|
2019-07-12 21:45:56 +03:00
|
|
|
// sqrt calculates square-root of the provided value.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn sqrt(a f64) f64 {
|
2019-06-22 21:20:28 +03:00
|
|
|
return C.sqrt(a)
|
|
|
|
}
|
2019-07-02 13:50:33 +03:00
|
|
|
// tan calculates tangent.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn tan(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.tan(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// tanh calculates hyperbolic tangent.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn tanh(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.tanh(a)
|
|
|
|
}
|
|
|
|
|
2019-07-02 13:50:33 +03:00
|
|
|
// trunc rounds a toward zero, returning the nearest integral value that is not
|
2019-06-30 16:33:37 +03:00
|
|
|
// larger in magnitude than a.
|
2019-06-26 18:49:50 +03:00
|
|
|
pub fn trunc(a f64) f64 {
|
2019-06-23 09:19:37 +03:00
|
|
|
return C.trunc(a)
|
|
|
|
}
|