diff options
Diffstat (limited to 'v_windows/v/old/vlib/math')
26 files changed, 6107 insertions, 0 deletions
diff --git a/v_windows/v/old/vlib/math/big/big.v b/v_windows/v/old/vlib/math/big/big.v new file mode 100644 index 0000000..a93aa9e --- /dev/null +++ b/v_windows/v/old/vlib/math/big/big.v @@ -0,0 +1,322 @@ +module big + +// Wrapper for https://github.com/kokke/tiny-bignum-c +#flag -I @VEXEROOT/thirdparty/bignum +#flag @VEXEROOT/thirdparty/bignum/bn.o +#include "bn.h" + +struct C.bn { +mut: + array [32]u32 +} + +// Big unsigned integer number. +type Number = C.bn + +fn C.bignum_init(n &Number) + +fn C.bignum_from_int(n &Number, i u64) + +fn C.bignum_to_int(n &Number) int + +fn C.bignum_from_string(n &Number, s &char, nbytes int) + +fn C.bignum_to_string(n &Number, s &char, maxsize int) + +// c = a + b +fn C.bignum_add(a &Number, b &Number, c &Number) + +// c = a - b +fn C.bignum_sub(a &Number, b &Number, c &Number) + +// c = a * b +fn C.bignum_mul(a &Number, b &Number, c &Number) + +// c = a / b +fn C.bignum_div(a &Number, b &Number, c &Number) + +// c = a % b +fn C.bignum_mod(a &Number, b &Number, c &Number) + +// c = a/b d=a%b +fn C.bignum_divmod(a &Number, b &Number, c &Number, d &Number) + +// c = a & b +fn C.bignum_and(a &Number, b &Number, c &Number) + +// c = a | b +fn C.bignum_or(a &Number, b &Number, c &Number) + +// c = a xor b +fn C.bignum_xor(a &Number, b &Number, c &Number) + +// b = a << nbits +fn C.bignum_lshift(a &Number, b &Number, nbits int) + +// b = a >> nbits +fn C.bignum_rshift(a &Number, b &Number, nbits int) + +fn C.bignum_cmp(a &Number, b &Number) int + +fn C.bignum_is_zero(a &Number) int + +// n++ +fn C.bignum_inc(n &Number) + +// n-- +fn C.bignum_dec(n &Number) + +// c = a ^ b +fn C.bignum_pow(a &Number, b &Number, c &Number) + +// b = integer_square_root_of(a) +fn C.bignum_isqrt(a &Number, b &Number) + +// copy src number to dst number +fn C.bignum_assign(dst &Number, src &Number) + +// new returns a bignum, initialized to 0 +pub fn new() Number { + return Number{} +} + +// conversion actions to/from big numbers: +// from_int converts an ordinary int number `i` to big.Number +pub fn from_int(i int) Number { + n := Number{} + C.bignum_from_int(&n, i) + return n +} + +// from_u64 converts an ordinary u64 number `u` to big.Number +pub fn from_u64(u u64) Number { + n := Number{} + C.bignum_from_int(&n, u) + return n +} + +// from_hex_string converts a hex string to big.Number +pub fn from_hex_string(input string) Number { + mut s := input.trim_prefix('0x') + if s.len == 0 { + s = '0' + } + padding := '0'.repeat((8 - s.len % 8) % 8) + s = padding + s + n := Number{} + C.bignum_from_string(&n, &char(s.str), s.len) + return n +} + +// from_string converts a decimal string to big.Number +pub fn from_string(input string) Number { + mut n := from_int(0) + for _, c in input { + d := from_int(int(c - `0`)) + n = (n * big.ten) + d + } + return n +} + +// .int() converts (a small) big.Number `n` to an ordinary integer. +pub fn (n &Number) int() int { + r := C.bignum_to_int(n) + return r +} + +const ( + ten = from_int(10) +) + +// .str returns a decimal representation of the big unsigned integer number n. +pub fn (n &Number) str() string { + if n.is_zero() { + return '0' + } + mut digits := []byte{} + mut x := n.clone() + + for !x.is_zero() { + // changes to reflect new api + div, mod := divmod(&x, &big.ten) + digits << byte(mod.int()) + `0` + x = div + } + return digits.reverse().bytestr() +} + +// .hexstr returns a hexadecimal representation of the bignum `n` +pub fn (n &Number) hexstr() string { + mut buf := [8192]byte{} + mut s := '' + unsafe { + bp := &buf[0] + // NB: C.bignum_to_string(), returns the HEXADECIMAL representation of the bignum n + C.bignum_to_string(n, &char(bp), 8192) + s = tos_clone(bp) + } + if s.len == 0 { + return '0' + } + return s +} + +// ////////////////////////////////////////////////////////// +// overloaded ops for the numbers: +pub fn (a &Number) + (b &Number) Number { + c := Number{} + C.bignum_add(a, b, &c) + return c +} + +pub fn (a &Number) - (b &Number) Number { + c := Number{} + C.bignum_sub(a, b, &c) + return c +} + +pub fn (a &Number) * (b &Number) Number { + c := Number{} + C.bignum_mul(a, b, &c) + return c +} + +pub fn (a &Number) / (b &Number) Number { + c := Number{} + C.bignum_div(a, b, &c) + return c +} + +pub fn (a &Number) % (b &Number) Number { + c := Number{} + C.bignum_mod(a, b, &c) + return c +} + +// divmod returns a pair of quotient and remainder from div modulo operation +// between two bignums `a` and `b` +pub fn divmod(a &Number, b &Number) (Number, Number) { + c := Number{} + d := Number{} + C.bignum_divmod(a, b, &c, &d) + return c, d +} + +// ////////////////////////////////////////////////////////// +pub fn cmp(a &Number, b &Number) int { + return C.bignum_cmp(a, b) +} + +pub fn (a &Number) is_zero() bool { + return C.bignum_is_zero(a) != 0 +} + +pub fn (mut a Number) inc() { + C.bignum_inc(&a) +} + +pub fn (mut a Number) dec() { + C.bignum_dec(&a) +} + +pub fn pow(a &Number, b &Number) Number { + c := Number{} + C.bignum_pow(a, b, &c) + return c +} + +pub fn (a &Number) isqrt() Number { + b := Number{} + C.bignum_isqrt(a, &b) + return b +} + +// ////////////////////////////////////////////////////////// +pub fn b_and(a &Number, b &Number) Number { + c := Number{} + C.bignum_and(a, b, &c) + return c +} + +pub fn b_or(a &Number, b &Number) Number { + c := Number{} + C.bignum_or(a, b, &c) + return c +} + +pub fn b_xor(a &Number, b &Number) Number { + c := Number{} + C.bignum_xor(a, b, &c) + return c +} + +pub fn (a &Number) lshift(nbits int) Number { + b := Number{} + C.bignum_lshift(a, &b, nbits) + return b +} + +pub fn (a &Number) rshift(nbits int) Number { + b := Number{} + C.bignum_rshift(a, &b, nbits) + return b +} + +pub fn (a &Number) clone() Number { + b := Number{} + C.bignum_assign(&b, a) + return b +} + +// ////////////////////////////////////////////////////////// +pub fn factorial(nn &Number) Number { + mut n := nn.clone() + mut a := nn.clone() + n.dec() + mut i := 1 + for !n.is_zero() { + res := a * n + n.dec() + a = res + i++ + } + return a +} + +pub fn fact(n int) Number { + return factorial(from_int(n)) +} + +// bytes returns an array of the bytes for the number `n`, +// in little endian format, where .bytes()[0] is the least +// significant byte. The result is NOT trimmed, and will contain 0s, even +// after the significant bytes. +// This method is faster than .bytes_trimmed(), but may be less convenient. +// Example: assert big.from_int(1).bytes()[0] == byte(0x01) +// Example: assert big.from_int(1024).bytes()[1] == byte(0x04) +// Example: assert big.from_int(1048576).bytes()[2] == byte(0x10) +pub fn (n &Number) bytes() []byte { + mut res := []byte{len: 128, init: 0} + unsafe { C.memcpy(res.data, n, 128) } + return res +} + +// bytes_trimmed returns an array of the bytes for the number `n`, +// in little endian format, where .bytes_trimmed()[0] is the least +// significant byte. The result is trimmed, so that *the last* byte +// of the result is also the the last meaningfull byte, != 0 . +// Example: assert big.from_int(1).bytes_trimmed() == [byte(0x01)] +// Example: assert big.from_int(1024).bytes_trimmed() == [byte(0x00), 0x04] +// Example: assert big.from_int(1048576).bytes_trimmed() == [byte(0x00), 0x00, 0x10] +pub fn (n &Number) bytes_trimmed() []byte { + mut res := []byte{len: 128, init: 0} + unsafe { C.memcpy(res.data, n, 128) } + mut non_zero_idx := 127 + for ; non_zero_idx >= 0; non_zero_idx-- { + if res[non_zero_idx] != 0 { + break + } + } + res.trim(non_zero_idx + 1) + return res +} diff --git a/v_windows/v/old/vlib/math/big/big_test.v b/v_windows/v/old/vlib/math/big/big_test.v new file mode 100644 index 0000000..0742b81 --- /dev/null +++ b/v_windows/v/old/vlib/math/big/big_test.v @@ -0,0 +1,167 @@ +import math.big + +fn test_new_big() { + n := big.new() + assert sizeof(big.Number) == 128 + assert n.hexstr() == '0' +} + +fn test_from_int() { + assert big.from_int(255).hexstr() == 'ff' + assert big.from_int(127).hexstr() == '7f' + assert big.from_int(1024).hexstr() == '400' + assert big.from_int(2147483647).hexstr() == '7fffffff' + assert big.from_int(-1).hexstr() == 'ffffffffffffffff' +} + +fn test_from_u64() { + assert big.from_u64(255).hexstr() == 'ff' + assert big.from_u64(127).hexstr() == '7f' + assert big.from_u64(1024).hexstr() == '400' + assert big.from_u64(4294967295).hexstr() == 'ffffffff' + assert big.from_u64(4398046511104).hexstr() == '40000000000' + assert big.from_u64(-1).hexstr() == 'ffffffffffffffff' +} + +fn test_plus() { + mut a := big.from_u64(2) + b := big.from_u64(3) + c := a + b + assert c.hexstr() == '5' + assert (big.from_u64(1024) + big.from_u64(1024)).hexstr() == '800' + a += b + assert a.hexstr() == '5' + a.inc() + assert a.hexstr() == '6' + a.dec() + a.dec() + assert a.hexstr() == '4' +} + +fn test_minus() { + a := big.from_u64(2) + mut b := big.from_u64(3) + c := b - a + assert c.hexstr() == '1' + e := big.from_u64(1024) + ee := e - e + assert ee.hexstr() == '0' + b -= a + assert b.hexstr() == '1' +} + +fn test_divide() { + a := big.from_u64(2) + mut b := big.from_u64(3) + c := b / a + assert c.hexstr() == '1' + assert (b % a).hexstr() == '1' + e := big.from_u64(1024) // dec(1024) == hex(0x400) + ee := e / e + assert ee.hexstr() == '1' + assert (e / a).hexstr() == '200' + assert (e / (a * a)).hexstr() == '100' + b /= a + assert b.hexstr() == '1' +} + +fn test_multiply() { + mut a := big.from_u64(2) + b := big.from_u64(3) + c := b * a + assert c.hexstr() == '6' + e := big.from_u64(1024) + e2 := e * e + e4 := e2 * e2 + e8 := e2 * e2 * e2 * e2 + e9 := e8 + big.from_u64(1) + d := ((e9 * e9) + b) * c + assert e4.hexstr() == '10000000000' + assert e8.hexstr() == '100000000000000000000' + assert e9.hexstr() == '100000000000000000001' + assert d.hexstr() == '60000000000000000000c00000000000000000018' + a *= b + assert a.hexstr() == '6' +} + +fn test_mod() { + assert (big.from_u64(13) % big.from_u64(10)).int() == 3 + assert (big.from_u64(13) % big.from_u64(9)).int() == 4 + assert (big.from_u64(7) % big.from_u64(5)).int() == 2 +} + +fn test_divmod() { + x, y := big.divmod(big.from_u64(13), big.from_u64(10)) + assert x.int() == 1 + assert y.int() == 3 + p, q := big.divmod(big.from_u64(13), big.from_u64(9)) + assert p.int() == 1 + assert q.int() == 4 + c, d := big.divmod(big.from_u64(7), big.from_u64(5)) + assert c.int() == 1 + assert d.int() == 2 +} + +fn test_from_str() { + assert big.from_string('9870123').str() == '9870123' + assert big.from_string('').str() == '0' + assert big.from_string('0').str() == '0' + assert big.from_string('1').str() == '1' + for i := 1; i < 307; i += 61 { + input := '9'.repeat(i) + out := big.from_string(input).str() + // eprintln('>> i: $i input: $input.str()') + // eprintln('>> i: $i out: $out.str()') + assert input == out + } +} + +fn test_from_hex_str() { + assert big.from_hex_string('0x123').hexstr() == '123' + for i in 1 .. 33 { + input := 'e'.repeat(i) + out := big.from_hex_string(input).hexstr() + assert input == out + } + assert big.from_string('0').hexstr() == '0' +} + +fn test_str() { + assert big.from_u64(255).str() == '255' + assert big.from_u64(127).str() == '127' + assert big.from_u64(1024).str() == '1024' + assert big.from_u64(4294967295).str() == '4294967295' + assert big.from_u64(4398046511104).str() == '4398046511104' + assert big.from_int(int(4294967295)).str() == '18446744073709551615' + assert big.from_int(-1).str() == '18446744073709551615' + assert big.from_hex_string('e'.repeat(80)).str() == '1993587900192849410235353592424915306962524220866209251950572167300738410728597846688097947807470' +} + +fn test_factorial() { + f5 := big.factorial(big.from_u64(5)) + assert f5.hexstr() == '78' + f100 := big.factorial(big.from_u64(100)) + assert f100.hexstr() == '1b30964ec395dc24069528d54bbda40d16e966ef9a70eb21b5b2943a321cdf10391745570cca9420c6ecb3b72ed2ee8b02ea2735c61a000000000000000000000000' +} + +fn trimbytes(n int, x []byte) []byte { + mut res := x.clone() + res.trim(n) + return res +} + +fn test_bytes() { + assert big.from_int(0).bytes().len == 128 + assert big.from_hex_string('e'.repeat(100)).bytes().len == 128 + assert trimbytes(3, big.from_int(1).bytes()) == [byte(0x01), 0x00, 0x00] + assert trimbytes(3, big.from_int(1024).bytes()) == [byte(0x00), 0x04, 0x00] + assert trimbytes(3, big.from_int(1048576).bytes()) == [byte(0x00), 0x00, 0x10] +} + +fn test_bytes_trimmed() { + assert big.from_int(0).bytes_trimmed().len == 0 + assert big.from_hex_string('AB'.repeat(50)).bytes_trimmed().len == 50 + assert big.from_int(1).bytes_trimmed() == [byte(0x01)] + assert big.from_int(1024).bytes_trimmed() == [byte(0x00), 0x04] + assert big.from_int(1048576).bytes_trimmed() == [byte(0x00), 0x00, 0x10] +} diff --git a/v_windows/v/old/vlib/math/bits.v b/v_windows/v/old/vlib/math/bits.v new file mode 100644 index 0000000..0626504 --- /dev/null +++ b/v_windows/v/old/vlib/math/bits.v @@ -0,0 +1,59 @@ +// Copyright (c) 2019-2021 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 + +const ( + uvnan = u64(0x7FF8000000000001) + uvinf = u64(0x7FF0000000000000) + uvneginf = u64(0xFFF0000000000000) + uvone = u64(0x3FF0000000000000) + mask = 0x7FF + shift = 64 - 11 - 1 + bias = 1023 + sign_mask = (u64(1) << 63) + frac_mask = ((u64(1) << u64(shift)) - u64(1)) +) + +// inf returns positive infinity if sign >= 0, negative infinity if sign < 0. +pub fn inf(sign int) f64 { + v := if sign >= 0 { math.uvinf } else { math.uvneginf } + return f64_from_bits(v) +} + +// nan returns an IEEE 754 ``not-a-number'' value. +pub fn nan() f64 { + return f64_from_bits(math.uvnan) +} + +// is_nan reports whether f is an IEEE 754 ``not-a-number'' value. +pub fn is_nan(f f64) bool { + // IEEE 754 says that only NaNs satisfy f != f. + // To avoid the floating-point hardware, could use: + // x := f64_bits(f); + // return u32(x>>shift)&mask == mask && x != uvinf && x != uvneginf + return f != f +} + +// is_inf reports whether f is an infinity, according to sign. +// If sign > 0, is_inf reports whether f is positive infinity. +// If sign < 0, is_inf reports whether f is negative infinity. +// If sign == 0, is_inf reports whether f is either infinity. +pub fn is_inf(f f64, sign int) bool { + // Test for infinity by comparing against maximum float. + // To avoid the floating-point hardware, could use: + // x := f64_bits(f); + // return sign >= 0 && x == uvinf || sign <= 0 && x == uvneginf; + return (sign >= 0 && f > max_f64) || (sign <= 0 && f < -max_f64) +} + +// NOTE: (joe-c) exponent notation is borked +// normalize returns a normal number y and exponent exp +// satisfying x == y × 2**exp. It assumes x is finite and non-zero. +// pub fn normalize(x f64) (f64, int) { +// smallest_normal := 2.2250738585072014e-308 // 2**-1022 +// if abs(x) < smallest_normal { +// return x * (1 << 52), -52 +// } +// return x, 0 +// } diff --git a/v_windows/v/old/vlib/math/bits/bits.v b/v_windows/v/old/vlib/math/bits/bits.v new file mode 100644 index 0000000..a3d2220 --- /dev/null +++ b/v_windows/v/old/vlib/math/bits/bits.v @@ -0,0 +1,478 @@ +// Copyright (c) 2019-2021 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 bits + +const ( + // See http://supertech.csail.mit.edu/papers/debruijn.pdf + de_bruijn32 = u32(0x077CB531) + de_bruijn32tab = [byte(0), 1, 28, 2, 29, 14, 24, 3, 30, 22, 20, 15, 25, 17, 4, 8, 31, 27, 13, + 23, 21, 19, 16, 7, 26, 12, 18, 6, 11, 5, 10, 9] + de_bruijn64 = u64(0x03f79d71b4ca8b09) + de_bruijn64tab = [byte(0), 1, 56, 2, 57, 49, 28, 3, 61, 58, 42, 50, 38, 29, 17, 4, 62, 47, + 59, 36, 45, 43, 51, 22, 53, 39, 33, 30, 24, 18, 12, 5, 63, 55, 48, 27, 60, 41, 37, 16, + 46, 35, 44, 21, 52, 32, 23, 11, 54, 26, 40, 15, 34, 20, 31, 10, 25, 14, 19, 9, 13, 8, 7, + 6, + ] +) + +const ( + m0 = u64(0x5555555555555555) // 01010101 ... + m1 = u64(0x3333333333333333) // 00110011 ... + m2 = u64(0x0f0f0f0f0f0f0f0f) // 00001111 ... + m3 = u64(0x00ff00ff00ff00ff) // etc. + m4 = u64(0x0000ffff0000ffff) +) + +const ( + // save importing math mod just for these + max_u32 = u32(4294967295) + max_u64 = u64(18446744073709551615) +) + +// --- LeadingZeros --- +// leading_zeros_8 returns the number of leading zero bits in x; the result is 8 for x == 0. +pub fn leading_zeros_8(x byte) int { + return 8 - len_8(x) +} + +// leading_zeros_16 returns the number of leading zero bits in x; the result is 16 for x == 0. +pub fn leading_zeros_16(x u16) int { + return 16 - len_16(x) +} + +// leading_zeros_32 returns the number of leading zero bits in x; the result is 32 for x == 0. +pub fn leading_zeros_32(x u32) int { + return 32 - len_32(x) +} + +// leading_zeros_64 returns the number of leading zero bits in x; the result is 64 for x == 0. +pub fn leading_zeros_64(x u64) int { + return 64 - len_64(x) +} + +// --- TrailingZeros --- +// trailing_zeros_8 returns the number of trailing zero bits in x; the result is 8 for x == 0. +pub fn trailing_zeros_8(x byte) int { + return int(ntz_8_tab[x]) +} + +// trailing_zeros_16 returns the number of trailing zero bits in x; the result is 16 for x == 0. +pub fn trailing_zeros_16(x u16) int { + if x == 0 { + return 16 + } + // see comment in trailing_zeros_64 + return int(bits.de_bruijn32tab[u32(x & -x) * bits.de_bruijn32 >> (32 - 5)]) +} + +// trailing_zeros_32 returns the number of trailing zero bits in x; the result is 32 for x == 0. +pub fn trailing_zeros_32(x u32) int { + if x == 0 { + return 32 + } + // see comment in trailing_zeros_64 + return int(bits.de_bruijn32tab[(x & -x) * bits.de_bruijn32 >> (32 - 5)]) +} + +// trailing_zeros_64 returns the number of trailing zero bits in x; the result is 64 for x == 0. +pub fn trailing_zeros_64(x u64) int { + if x == 0 { + return 64 + } + // If popcount is fast, replace code below with return popcount(^x & (x - 1)). + // + // x & -x leaves only the right-most bit set in the word. Let k be the + // index of that bit. Since only a single bit is set, the value is two + // to the power of k. Multiplying by a power of two is equivalent to + // left shifting, in this case by k bits. The de Bruijn (64 bit) constant + // is such that all six bit, consecutive substrings are distinct. + // Therefore, if we have a left shifted version of this constant we can + // find by how many bits it was shifted by looking at which six bit + // substring ended up at the top of the word. + // (Knuth, volume 4, section 7.3.1) + return int(bits.de_bruijn64tab[(x & -x) * bits.de_bruijn64 >> (64 - 6)]) +} + +// --- OnesCount --- +// ones_count_8 returns the number of one bits ("population count") in x. +pub fn ones_count_8(x byte) int { + return int(pop_8_tab[x]) +} + +// ones_count_16 returns the number of one bits ("population count") in x. +pub fn ones_count_16(x u16) int { + return int(pop_8_tab[x >> 8] + pop_8_tab[x & u16(0xff)]) +} + +// ones_count_32 returns the number of one bits ("population count") in x. +pub fn ones_count_32(x u32) int { + return int(pop_8_tab[x >> 24] + pop_8_tab[x >> 16 & 0xff] + pop_8_tab[x >> 8 & 0xff] + + pop_8_tab[x & u32(0xff)]) +} + +// ones_count_64 returns the number of one bits ("population count") in x. +pub fn ones_count_64(x u64) int { + // Implementation: Parallel summing of adjacent bits. + // See "Hacker's Delight", Chap. 5: Counting Bits. + // The following pattern shows the general approach: + // + // x = x>>1&(m0&m) + x&(m0&m) + // x = x>>2&(m1&m) + x&(m1&m) + // x = x>>4&(m2&m) + x&(m2&m) + // x = x>>8&(m3&m) + x&(m3&m) + // x = x>>16&(m4&m) + x&(m4&m) + // x = x>>32&(m5&m) + x&(m5&m) + // return int(x) + // + // Masking (& operations) can be left away when there's no + // danger that a field's sum will carry over into the next + // field: Since the result cannot be > 64, 8 bits is enough + // and we can ignore the masks for the shifts by 8 and up. + // Per "Hacker's Delight", the first line can be simplified + // more, but it saves at best one instruction, so we leave + // it alone for clarity. + mut y := (x >> u64(1) & (bits.m0 & bits.max_u64)) + (x & (bits.m0 & bits.max_u64)) + y = (y >> u64(2) & (bits.m1 & bits.max_u64)) + (y & (bits.m1 & bits.max_u64)) + y = ((y >> 4) + y) & (bits.m2 & bits.max_u64) + y += y >> 8 + y += y >> 16 + y += y >> 32 + return int(y) & ((1 << 7) - 1) +} + +// --- RotateLeft --- +// rotate_left_8 returns the value of x rotated left by (k mod 8) bits. +// To rotate x right by k bits, call rotate_left_8(x, -k). +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn rotate_left_8(x byte, k int) byte { + n := byte(8) + s := byte(k) & (n - byte(1)) + return ((x << s) | (x >> (n - s))) +} + +// rotate_left_16 returns the value of x rotated left by (k mod 16) bits. +// To rotate x right by k bits, call rotate_left_16(x, -k). +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn rotate_left_16(x u16, k int) u16 { + n := u16(16) + s := u16(k) & (n - u16(1)) + return ((x << s) | (x >> (n - s))) +} + +// rotate_left_32 returns the value of x rotated left by (k mod 32) bits. +// To rotate x right by k bits, call rotate_left_32(x, -k). +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn rotate_left_32(x u32, k int) u32 { + n := u32(32) + s := u32(k) & (n - u32(1)) + return ((x << s) | (x >> (n - s))) +} + +// rotate_left_64 returns the value of x rotated left by (k mod 64) bits. +// To rotate x right by k bits, call rotate_left_64(x, -k). +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn rotate_left_64(x u64, k int) u64 { + n := u64(64) + s := u64(k) & (n - u64(1)) + return ((x << s) | (x >> (n - s))) +} + +// --- Reverse --- +// reverse_8 returns the value of x with its bits in reversed order. +[inline] +pub fn reverse_8(x byte) byte { + return rev_8_tab[x] +} + +// reverse_16 returns the value of x with its bits in reversed order. +[inline] +pub fn reverse_16(x u16) u16 { + return u16(rev_8_tab[x >> 8]) | (u16(rev_8_tab[x & u16(0xff)]) << 8) +} + +// reverse_32 returns the value of x with its bits in reversed order. +[inline] +pub fn reverse_32(x u32) u32 { + mut y := ((x >> u32(1) & (bits.m0 & bits.max_u32)) | ((x & (bits.m0 & bits.max_u32)) << 1)) + y = ((y >> u32(2) & (bits.m1 & bits.max_u32)) | ((y & (bits.m1 & bits.max_u32)) << u32(2))) + y = ((y >> u32(4) & (bits.m2 & bits.max_u32)) | ((y & (bits.m2 & bits.max_u32)) << u32(4))) + return reverse_bytes_32(u32(y)) +} + +// reverse_64 returns the value of x with its bits in reversed order. +[inline] +pub fn reverse_64(x u64) u64 { + mut y := ((x >> u64(1) & (bits.m0 & bits.max_u64)) | ((x & (bits.m0 & bits.max_u64)) << 1)) + y = ((y >> u64(2) & (bits.m1 & bits.max_u64)) | ((y & (bits.m1 & bits.max_u64)) << 2)) + y = ((y >> u64(4) & (bits.m2 & bits.max_u64)) | ((y & (bits.m2 & bits.max_u64)) << 4)) + return reverse_bytes_64(y) +} + +// --- ReverseBytes --- +// reverse_bytes_16 returns the value of x with its bytes in reversed order. +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn reverse_bytes_16(x u16) u16 { + return (x >> 8) | (x << 8) +} + +// reverse_bytes_32 returns the value of x with its bytes in reversed order. +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn reverse_bytes_32(x u32) u32 { + y := ((x >> u32(8) & (bits.m3 & bits.max_u32)) | ((x & (bits.m3 & bits.max_u32)) << u32(8))) + return u32((y >> 16) | (y << 16)) +} + +// reverse_bytes_64 returns the value of x with its bytes in reversed order. +// +// This function's execution time does not depend on the inputs. +[inline] +pub fn reverse_bytes_64(x u64) u64 { + mut y := ((x >> u64(8) & (bits.m3 & bits.max_u64)) | ((x & (bits.m3 & bits.max_u64)) << u64(8))) + y = ((y >> u64(16) & (bits.m4 & bits.max_u64)) | ((y & (bits.m4 & bits.max_u64)) << u64(16))) + return (y >> 32) | (y << 32) +} + +// --- Len --- +// len_8 returns the minimum number of bits required to represent x; the result is 0 for x == 0. +pub fn len_8(x byte) int { + return int(len_8_tab[x]) +} + +// len_16 returns the minimum number of bits required to represent x; the result is 0 for x == 0. +pub fn len_16(x u16) int { + mut y := x + mut n := 0 + if y >= 1 << 8 { + y >>= 8 + n = 8 + } + return n + int(len_8_tab[y]) +} + +// len_32 returns the minimum number of bits required to represent x; the result is 0 for x == 0. +pub fn len_32(x u32) int { + mut y := x + mut n := 0 + if y >= (1 << 16) { + y >>= 16 + n = 16 + } + if y >= (1 << 8) { + y >>= 8 + n += 8 + } + return n + int(len_8_tab[y]) +} + +// len_64 returns the minimum number of bits required to represent x; the result is 0 for x == 0. +pub fn len_64(x u64) int { + mut y := x + mut n := 0 + if y >= u64(1) << u64(32) { + y >>= 32 + n = 32 + } + if y >= u64(1) << u64(16) { + y >>= 16 + n += 16 + } + if y >= u64(1) << u64(8) { + y >>= 8 + n += 8 + } + return n + int(len_8_tab[y]) +} + +// --- Add with carry --- +// Add returns the sum with carry of x, y and carry: sum = x + y + carry. +// The carry input must be 0 or 1; otherwise the behavior is undefined. +// The carryOut output is guaranteed to be 0 or 1. +// +// add_32 returns the sum with carry of x, y and carry: sum = x + y + carry. +// The carry input must be 0 or 1; otherwise the behavior is undefined. +// The carryOut output is guaranteed to be 0 or 1. +// +// This function's execution time does not depend on the inputs. +pub fn add_32(x u32, y u32, carry u32) (u32, u32) { + sum64 := u64(x) + u64(y) + u64(carry) + sum := u32(sum64) + carry_out := u32(sum64 >> 32) + return sum, carry_out +} + +// add_64 returns the sum with carry of x, y and carry: sum = x + y + carry. +// The carry input must be 0 or 1; otherwise the behavior is undefined. +// The carryOut output is guaranteed to be 0 or 1. +// +// This function's execution time does not depend on the inputs. +pub fn add_64(x u64, y u64, carry u64) (u64, u64) { + sum := x + y + carry + // The sum will overflow if both top bits are set (x & y) or if one of them + // is (x | y), and a carry from the lower place happened. If such a carry + // happens, the top bit will be 1 + 0 + 1 = 0 (&^ sum). + carry_out := ((x & y) | ((x | y) & ~sum)) >> 63 + return sum, carry_out +} + +// --- Subtract with borrow --- +// Sub returns the difference of x, y and borrow: diff = x - y - borrow. +// The borrow input must be 0 or 1; otherwise the behavior is undefined. +// The borrowOut output is guaranteed to be 0 or 1. +// +// sub_32 returns the difference of x, y and borrow, diff = x - y - borrow. +// The borrow input must be 0 or 1; otherwise the behavior is undefined. +// The borrowOut output is guaranteed to be 0 or 1. +// +// This function's execution time does not depend on the inputs. +pub fn sub_32(x u32, y u32, borrow u32) (u32, u32) { + diff := x - y - borrow + // The difference will underflow if the top bit of x is not set and the top + // bit of y is set (^x & y) or if they are the same (^(x ^ y)) and a borrow + // from the lower place happens. If that borrow happens, the result will be + // 1 - 1 - 1 = 0 - 0 - 1 = 1 (& diff). + borrow_out := ((~x & y) | (~(x ^ y) & diff)) >> 31 + return diff, borrow_out +} + +// sub_64 returns the difference of x, y and borrow: diff = x - y - borrow. +// The borrow input must be 0 or 1; otherwise the behavior is undefined. +// The borrowOut output is guaranteed to be 0 or 1. +// +// This function's execution time does not depend on the inputs. +pub fn sub_64(x u64, y u64, borrow u64) (u64, u64) { + diff := x - y - borrow + // See Sub32 for the bit logic. + borrow_out := ((~x & y) | (~(x ^ y) & diff)) >> 63 + return diff, borrow_out +} + +// --- Full-width multiply --- +const ( + two32 = u64(0x100000000) + mask32 = two32 - 1 + overflow_error = 'Overflow Error' + divide_error = 'Divide Error' +) + +// mul_32 returns the 64-bit product of x and y: (hi, lo) = x * y +// with the product bits' upper half returned in hi and the lower +// half returned in lo. +// +// This function's execution time does not depend on the inputs. +pub fn mul_32(x u32, y u32) (u32, u32) { + tmp := u64(x) * u64(y) + hi := u32(tmp >> 32) + lo := u32(tmp) + return hi, lo +} + +// mul_64 returns the 128-bit product of x and y: (hi, lo) = x * y +// with the product bits' upper half returned in hi and the lower +// half returned in lo. +// +// This function's execution time does not depend on the inputs. +pub fn mul_64(x u64, y u64) (u64, u64) { + x0 := x & bits.mask32 + x1 := x >> 32 + y0 := y & bits.mask32 + y1 := y >> 32 + w0 := x0 * y0 + t := x1 * y0 + (w0 >> 32) + mut w1 := t & bits.mask32 + w2 := t >> 32 + w1 += x0 * y1 + hi := x1 * y1 + w2 + (w1 >> 32) + lo := x * y + return hi, lo +} + +// --- Full-width divide --- +// div_32 returns the quotient and remainder of (hi, lo) divided by y: +// quo = (hi, lo)/y, rem = (hi, lo)%y with the dividend bits' upper +// half in parameter hi and the lower half in parameter lo. +// div_32 panics for y == 0 (division by zero) or y <= hi (quotient overflow). +pub fn div_32(hi u32, lo u32, y u32) (u32, u32) { + if y != 0 && y <= hi { + panic(bits.overflow_error) + } + z := (u64(hi) << 32) | u64(lo) + quo := u32(z / u64(y)) + rem := u32(z % u64(y)) + return quo, rem +} + +// div_64 returns the quotient and remainder of (hi, lo) divided by y: +// quo = (hi, lo)/y, rem = (hi, lo)%y with the dividend bits' upper +// half in parameter hi and the lower half in parameter lo. +// div_64 panics for y == 0 (division by zero) or y <= hi (quotient overflow). +pub fn div_64(hi u64, lo u64, y1 u64) (u64, u64) { + mut y := y1 + if y == 0 { + panic(bits.overflow_error) + } + if y <= hi { + panic(bits.overflow_error) + } + s := u32(leading_zeros_64(y)) + y <<= s + yn1 := y >> 32 + yn0 := y & bits.mask32 + un32 := (hi << s) | (lo >> (64 - s)) + un10 := lo << s + un1 := un10 >> 32 + un0 := un10 & bits.mask32 + mut q1 := un32 / yn1 + mut rhat := un32 - q1 * yn1 + for q1 >= bits.two32 || q1 * yn0 > bits.two32 * rhat + un1 { + q1-- + rhat += yn1 + if rhat >= bits.two32 { + break + } + } + un21 := un32 * bits.two32 + un1 - q1 * y + mut q0 := un21 / yn1 + rhat = un21 - q0 * yn1 + for q0 >= bits.two32 || q0 * yn0 > bits.two32 * rhat + un0 { + q0-- + rhat += yn1 + if rhat >= bits.two32 { + break + } + } + return q1 * bits.two32 + q0, (un21 * bits.two32 + un0 - q0 * y) >> s +} + +// rem_32 returns the remainder of (hi, lo) divided by y. Rem32 panics +// for y == 0 (division by zero) but, unlike Div32, it doesn't panic +// on a quotient overflow. +pub fn rem_32(hi u32, lo u32, y u32) u32 { + return u32(((u64(hi) << 32) | u64(lo)) % u64(y)) +} + +// rem_64 returns the remainder of (hi, lo) divided by y. Rem64 panics +// for y == 0 (division by zero) but, unlike div_64, it doesn't panic +// on a quotient overflow. +pub fn rem_64(hi u64, lo u64, y u64) u64 { + // We scale down hi so that hi < y, then use div_64 to compute the + // rem with the guarantee that it won't panic on quotient overflow. + // Given that + // hi ≡ hi%y (mod y) + // we have + // hi<<64 + lo ≡ (hi%y)<<64 + lo (mod y) + _, rem := div_64(hi % y, lo, y) + return rem +} diff --git a/v_windows/v/old/vlib/math/bits/bits_tables.v b/v_windows/v/old/vlib/math/bits/bits_tables.v new file mode 100644 index 0000000..830e1e8 --- /dev/null +++ b/v_windows/v/old/vlib/math/bits/bits_tables.v @@ -0,0 +1,79 @@ +// Copyright (c) 2019-2021 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 bits + +const ( + ntz_8_tab = [byte(0x08), 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, + 0x02, 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, + 0x00, 0x02, 0x00, 0x01, 0x00, 0x05, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, + 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, + 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x06, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, + 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, + 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x05, 0x00, 0x01, 0x00, 0x02, 0x00, + 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, + 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x07, 0x00, 0x01, 0x00, + 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, + 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x05, 0x00, + 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x04, + 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, + 0x06, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, + 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, 0x00, + 0x01, 0x00, 0x05, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, 0x02, + 0x00, 0x01, 0x00, 0x04, 0x00, 0x01, 0x00, 0x02, 0x00, 0x01, 0x00, 0x03, 0x00, 0x01, 0x00, + 0x02, 0x00, 0x01, 0x00] + pop_8_tab = [byte(0x00), 0x01, 0x01, 0x02, 0x01, 0x02, 0x02, 0x03, 0x01, 0x02, 0x02, 0x03, + 0x02, 0x03, 0x03, 0x04, 0x01, 0x02, 0x02, 0x03, 0x02, 0x03, 0x03, 0x04, 0x02, 0x03, 0x03, + 0x04, 0x03, 0x04, 0x04, 0x05, 0x01, 0x02, 0x02, 0x03, 0x02, 0x03, 0x03, 0x04, 0x02, 0x03, + 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x03, + 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x01, 0x02, 0x02, 0x03, 0x02, 0x03, 0x03, 0x04, + 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, 0x04, + 0x05, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, + 0x04, 0x05, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x03, 0x04, 0x04, 0x05, 0x04, + 0x05, 0x05, 0x06, 0x04, 0x05, 0x05, 0x06, 0x05, 0x06, 0x06, 0x07, 0x01, 0x02, 0x02, 0x03, + 0x02, 0x03, 0x03, 0x04, 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x02, 0x03, 0x03, + 0x04, 0x03, 0x04, 0x04, 0x05, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x02, 0x03, + 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x03, + 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x04, 0x05, 0x05, 0x06, 0x05, 0x06, 0x06, 0x07, + 0x02, 0x03, 0x03, 0x04, 0x03, 0x04, 0x04, 0x05, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, + 0x06, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x04, 0x05, 0x05, 0x06, 0x05, 0x06, + 0x06, 0x07, 0x03, 0x04, 0x04, 0x05, 0x04, 0x05, 0x05, 0x06, 0x04, 0x05, 0x05, 0x06, 0x05, + 0x06, 0x06, 0x07, 0x04, 0x05, 0x05, 0x06, 0x05, 0x06, 0x06, 0x07, 0x05, 0x06, 0x06, 0x07, + 0x06, 0x07, 0x07, 0x08] + rev_8_tab = [byte(0x00), 0x80, 0x40, 0xc0, 0x20, 0xa0, 0x60, 0xe0, 0x10, 0x90, 0x50, 0xd0, + 0x30, 0xb0, 0x70, 0xf0, 0x08, 0x88, 0x48, 0xc8, 0x28, 0xa8, 0x68, 0xe8, 0x18, 0x98, 0x58, + 0xd8, 0x38, 0xb8, 0x78, 0xf8, 0x04, 0x84, 0x44, 0xc4, 0x24, 0xa4, 0x64, 0xe4, 0x14, 0x94, + 0x54, 0xd4, 0x34, 0xb4, 0x74, 0xf4, 0x0c, 0x8c, 0x4c, 0xcc, 0x2c, 0xac, 0x6c, 0xec, 0x1c, + 0x9c, 0x5c, 0xdc, 0x3c, 0xbc, 0x7c, 0xfc, 0x02, 0x82, 0x42, 0xc2, 0x22, 0xa2, 0x62, 0xe2, + 0x12, 0x92, 0x52, 0xd2, 0x32, 0xb2, 0x72, 0xf2, 0x0a, 0x8a, 0x4a, 0xca, 0x2a, 0xaa, 0x6a, + 0xea, 0x1a, 0x9a, 0x5a, 0xda, 0x3a, 0xba, 0x7a, 0xfa, 0x06, 0x86, 0x46, 0xc6, 0x26, 0xa6, + 0x66, 0xe6, 0x16, 0x96, 0x56, 0xd6, 0x36, 0xb6, 0x76, 0xf6, 0x0e, 0x8e, 0x4e, 0xce, 0x2e, + 0xae, 0x6e, 0xee, 0x1e, 0x9e, 0x5e, 0xde, 0x3e, 0xbe, 0x7e, 0xfe, 0x01, 0x81, 0x41, 0xc1, + 0x21, 0xa1, 0x61, 0xe1, 0x11, 0x91, 0x51, 0xd1, 0x31, 0xb1, 0x71, 0xf1, 0x09, 0x89, 0x49, + 0xc9, 0x29, 0xa9, 0x69, 0xe9, 0x19, 0x99, 0x59, 0xd9, 0x39, 0xb9, 0x79, 0xf9, 0x05, 0x85, + 0x45, 0xc5, 0x25, 0xa5, 0x65, 0xe5, 0x15, 0x95, 0x55, 0xd5, 0x35, 0xb5, 0x75, 0xf5, 0x0d, + 0x8d, 0x4d, 0xcd, 0x2d, 0xad, 0x6d, 0xed, 0x1d, 0x9d, 0x5d, 0xdd, 0x3d, 0xbd, 0x7d, 0xfd, + 0x03, 0x83, 0x43, 0xc3, 0x23, 0xa3, 0x63, 0xe3, 0x13, 0x93, 0x53, 0xd3, 0x33, 0xb3, 0x73, + 0xf3, 0x0b, 0x8b, 0x4b, 0xcb, 0x2b, 0xab, 0x6b, 0xeb, 0x1b, 0x9b, 0x5b, 0xdb, 0x3b, 0xbb, + 0x7b, 0xfb, 0x07, 0x87, 0x47, 0xc7, 0x27, 0xa7, 0x67, 0xe7, 0x17, 0x97, 0x57, 0xd7, 0x37, + 0xb7, 0x77, 0xf7, 0x0f, 0x8f, 0x4f, 0xcf, 0x2f, 0xaf, 0x6f, 0xef, 0x1f, 0x9f, 0x5f, 0xdf, + 0x3f, 0xbf, 0x7f, 0xff] + len_8_tab = [byte(0x00), 0x01, 0x02, 0x02, 0x03, 0x03, 0x03, 0x03, 0x04, 0x04, 0x04, 0x04, + 0x04, 0x04, 0x04, 0x04, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, 0x05, + 0x05, 0x05, 0x05, 0x05, 0x05, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, + 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, + 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x06, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, + 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, + 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, + 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, + 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x07, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, 0x08, + 0x08, 0x08, 0x08, 0x08] +) diff --git a/v_windows/v/old/vlib/math/bits/bits_test.v b/v_windows/v/old/vlib/math/bits/bits_test.v new file mode 100644 index 0000000..23a01b3 --- /dev/null +++ b/v_windows/v/old/vlib/math/bits/bits_test.v @@ -0,0 +1,288 @@ +// +// test suite for bits and bits math functions +// +module bits + +fn test_bits() { + mut i := 0 + mut i1 := u64(0) + + // + // --- LeadingZeros --- + // + + // 8 bit + i = 1 + for x in 0 .. 8 { + // C.printf("x:%02x lz: %d cmp: %d\n", i << x, leading_zeros_8(i << x), 7-x) + assert leading_zeros_8(byte(i << x)) == 7 - x + } + + // 16 bit + i = 1 + for x in 0 .. 16 { + // C.printf("x:%04x lz: %d cmp: %d\n", u16(i) << x, leading_zeros_16(u16(i) << x), 15-x) + assert leading_zeros_16(u16(i) << x) == 15 - x + } + + // 32 bit + i = 1 + for x in 0 .. 32 { + // C.printf("x:%08x lz: %d cmp: %d\n", u32(i) << x, leading_zeros_32(u32(i) << x), 31-x) + assert leading_zeros_32(u32(i) << x) == 31 - x + } + + // 64 bit + i = 1 + for x in 0 .. 64 { + // C.printf("x:%016llx lz: %llu cmp: %d\n", u64(i) << x, leading_zeros_64(u64(i) << x), 63-x) + assert leading_zeros_64(u64(i) << x) == 63 - x + } + + // + // --- ones_count --- + // + + // 8 bit + i = 0 + for x in 0 .. 9 { + // C.printf("x:%02x lz: %llu cmp: %d\n", byte(i), ones_count_8(byte(i)), x) + assert ones_count_8(byte(i)) == x + i = (i << 1) + 1 + } + + // 16 bit + i = 0 + for x in 0 .. 17 { + // C.printf("x:%04x lz: %llu cmp: %d\n", u16(i), ones_count_16(u16(i)), x) + assert ones_count_16(u16(i)) == x + i = (i << 1) + 1 + } + + // 32 bit + i = 0 + for x in 0 .. 33 { + // C.printf("x:%08x lz: %llu cmp: %d\n", u32(i), ones_count_32(u32(i)), x) + assert ones_count_32(u32(i)) == x + i = (i << 1) + 1 + } + + // 64 bit + i1 = 0 + for x in 0 .. 65 { + // C.printf("x:%016llx lz: %llu cmp: %d\n", u64(i1), ones_count_64(u64(i1)), x) + assert ones_count_64(i1) == x + i1 = (i1 << 1) + 1 + } + + // + // --- rotate_left/right --- + // + assert rotate_left_8(0x12, 4) == 0x21 + assert rotate_left_16(0x1234, 8) == 0x3412 + assert rotate_left_32(0x12345678, 16) == 0x56781234 + assert rotate_left_64(0x1234567887654321, 32) == 0x8765432112345678 + + // + // --- reverse --- + // + + // 8 bit + i = 0 + for _ in 0 .. 9 { + mut rv := byte(0) + mut bc := 0 + mut n := i + for bc < 8 { + rv = (rv << 1) | (byte(n) & 0x01) + bc++ + n = n >> 1 + } + // C.printf("x:%02x lz: %llu cmp: %d\n", byte(i), reverse_8(byte(i)), rv) + assert reverse_8(byte(i)) == rv + i = (i << 1) + 1 + } + + // 16 bit + i = 0 + for _ in 0 .. 17 { + mut rv := u16(0) + mut bc := 0 + mut n := i + for bc < 16 { + rv = (rv << 1) | (u16(n) & 0x01) + bc++ + n = n >> 1 + } + // C.printf("x:%04x lz: %llu cmp: %d\n", u16(i), reverse_16(u16(i)), rv) + assert reverse_16(u16(i)) == rv + i = (i << 1) + 1 + } + + // 32 bit + i = 0 + for _ in 0 .. 33 { + mut rv := u32(0) + mut bc := 0 + mut n := i + for bc < 32 { + rv = (rv << 1) | (u32(n) & 0x01) + bc++ + n = n >> 1 + } + // C.printf("x:%08x lz: %llu cmp: %d\n", u32(i), reverse_32(u32(i)), rv) + assert reverse_32(u32(i)) == rv + i = (i << 1) + 1 + } + + // 64 bit + i1 = 0 + for _ in 0 .. 64 { + mut rv := u64(0) + mut bc := 0 + mut n := i1 + for bc < 64 { + rv = (rv << 1) | (n & 0x01) + bc++ + n = n >> 1 + } + // C.printf("x:%016llx lz: %016llx cmp: %016llx\n", u64(i1), reverse_64(u64(i1)), rv) + assert reverse_64(i1) == rv + i1 = (i1 << 1) + 1 + } + + // + // --- add --- + // + + // 32 bit + i = 1 + for x in 0 .. 32 { + v := u32(i) << x + sum, carry := add_32(v, v, u32(0)) + // C.printf("x:%08x [%llu,%llu] %llu\n", u32(i) << x, sum, carry, u64(v) + u64(v)) + assert ((u64(carry) << 32) | u64(sum)) == u64(v) + u64(v) + } + mut sum_32t, mut carry_32t := add_32(0x8000_0000, 0x8000_0000, u32(0)) + assert sum_32t == u32(0) + assert carry_32t == u32(1) + + sum_32t, carry_32t = add_32(0xFFFF_FFFF, 0xFFFF_FFFF, u32(1)) + assert sum_32t == 0xFFFF_FFFF + assert carry_32t == u32(1) + + // 64 bit + i = 1 + for x in 0 .. 63 { + v := u64(i) << x + sum, carry := add_64(v, v, u64(0)) + // C.printf("x:%16x [%llu,%llu] %llu\n", u64(i) << x, sum, carry, u64(v >> 32) + u64(v >> 32)) + assert ((carry << 32) | sum) == v + v + } + mut sum_64t, mut carry_64t := add_64(0x8000_0000_0000_0000, 0x8000_0000_0000_0000, + u64(0)) + assert sum_64t == u64(0) + assert carry_64t == u64(1) + + sum_64t, carry_64t = add_64(0xFFFF_FFFF_FFFF_FFFF, 0xFFFF_FFFF_FFFF_FFFF, u64(1)) + assert sum_64t == 0xFFFF_FFFF_FFFF_FFFF + assert carry_64t == u64(1) + + // + // --- sub --- + // + + // 32 bit + i = 1 + for x in 1 .. 32 { + v0 := u32(i) << x + v1 := v0 >> 1 + mut diff, mut borrow_out := sub_32(v0, v1, u32(0)) + // C.printf("x:%08x [%llu,%llu] %08x\n", u32(i) << x, diff, borrow_out, v0 - v1) + assert diff == v1 + + diff, borrow_out = sub_32(v0, v1, u32(1)) + // C.printf("x:%08x [%llu,%llu] %08x\n", u32(i) << x, diff, borrow_out, v0 - v1) + assert diff == (v1 - 1) + assert borrow_out == u32(0) + + diff, borrow_out = sub_32(v1, v0, u32(1)) + // C.printf("x:%08x [%llu,%llu] %08x\n", u32(i) << x, diff, borrow_out, v1 - v0) + assert borrow_out == u32(1) + } + + // 64 bit + i = 1 + for x in 1 .. 64 { + v0 := u64(i) << x + v1 := v0 >> 1 + mut diff, mut borrow_out := sub_64(v0, v1, u64(0)) + // C.printf("x:%08x [%llu,%llu] %08x\n", u64(i) << x, diff, borrow_out, v0 - v1) + assert diff == v1 + + diff, borrow_out = sub_64(v0, v1, u64(1)) + // C.printf("x:%08x [%llu,%llu] %08x\n", u64(i) << x, diff, borrow_out, v0 - v1) + assert diff == (v1 - 1) + assert borrow_out == u64(0) + + diff, borrow_out = sub_64(v1, v0, u64(1)) + // C.printf("x:%08x [%llu,%llu] %08x\n",u64(i) << x, diff, borrow_out, v1 - v0) + assert borrow_out == u64(1) + } + + // + // --- mul --- + // + + // 32 bit + i = 1 + for x in 0 .. 32 { + v0 := u32(i) << x + v1 := v0 - 1 + hi, lo := mul_32(v0, v1) + assert (u64(hi) << 32) | (u64(lo)) == u64(v0) * u64(v1) + } + + // 64 bit + i = 1 + for x in 0 .. 64 { + v0 := u64(i) << x + v1 := v0 - 1 + hi, lo := mul_64(v0, v1) + // C.printf("v0: %llu v1: %llu [%llu,%llu] tt: %llu\n", v0, v1, hi, lo, (v0 >> 32) * (v1 >> 32)) + assert (hi & 0xFFFF_FFFF_0000_0000) == (((v0 >> 32) * (v1 >> 32)) & 0xFFFF_FFFF_0000_0000) + assert (lo & 0x0000_0000_FFFF_FFFF) == (((v0 & 0x0000_0000_FFFF_FFFF) * (v1 & 0x0000_0000_FFFF_FFFF)) & 0x0000_0000_FFFF_FFFF) + } + + // + // --- div --- + // + + // 32 bit + i = 1 + for x in 0 .. 31 { + hi := u32(i) << x + lo := hi - 1 + y := u32(3) << x + quo, rem := div_32(hi, lo, y) + // C.printf("[%08x_%08x] %08x (%08x,%08x)\n", hi, lo, y, quo, rem) + tst := ((u64(hi) << 32) | u64(lo)) + assert quo == (tst / u64(y)) + assert rem == (tst % u64(y)) + assert rem == rem_32(hi, lo, y) + } + + // 64 bit + i = 1 + for x in 0 .. 62 { + hi := u64(i) << x + lo := u64(2) // hi - 1 + y := u64(0x4000_0000_0000_0000) + quo, rem := div_64(hi, lo, y) + // C.printf("[%016llx_%016llx] %016llx (%016llx,%016llx)\n", hi, lo, y, quo, rem) + assert quo == u64(2) << (x + 1) + _, rem1 := div_64(hi % y, lo, y) + assert rem == rem1 + assert rem == rem_64(hi, lo, y) + } +} diff --git a/v_windows/v/old/vlib/math/complex/complex.v b/v_windows/v/old/vlib/math/complex/complex.v new file mode 100644 index 0000000..9fd7cc8 --- /dev/null +++ b/v_windows/v/old/vlib/math/complex/complex.v @@ -0,0 +1,374 @@ +// Copyright (c) 2019-2021 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 complex + +import math + +pub struct Complex { +pub: + re f64 + im f64 +} + +pub fn complex(re f64, im f64) Complex { + return Complex{re, im} +} + +// To String method +pub fn (c Complex) str() string { + mut out := '${c.re:.6f}' + out += if c.im >= 0 { '+${c.im:.6f}' } else { '${c.im:.6f}' } + out += 'i' + return out +} + +// Complex Modulus value +// mod() and abs() return the same +pub fn (c Complex) abs() f64 { + return C.hypot(c.re, c.im) +} + +pub fn (c Complex) mod() f64 { + return c.abs() +} + +// Complex Angle +pub fn (c Complex) angle() f64 { + return math.atan2(c.im, c.re) +} + +// Complex Addition c1 + c2 +pub fn (c1 Complex) + (c2 Complex) Complex { + return Complex{c1.re + c2.re, c1.im + c2.im} +} + +// Complex Substraction c1 - c2 +pub fn (c1 Complex) - (c2 Complex) Complex { + return Complex{c1.re - c2.re, c1.im - c2.im} +} + +// Complex Multiplication c1 * c2 +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 +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{((c1.re * c2.re) + ((c1.im * -c2.im) * -1)) / denom, ((c1.re * -c2.im) + + (c1.im * c2.re)) / denom} +} + +// Complex Conjugate +pub fn (c Complex) conjugate() Complex { + return Complex{c.re, -c.im} +} + +// Complex Additive Inverse +// Based on +// http://tutorial.math.lamar.edu/Extras/ComplexPrimer/Arithmetic.aspx +pub fn (c Complex) addinv() Complex { + return Complex{-c.re, -c.im} +} + +// Complex Multiplicative Inverse +// Based on +// http://tutorial.math.lamar.edu/Extras/ComplexPrimer/Arithmetic.aspx +pub fn (c Complex) mulinv() Complex { + return Complex{c.re / (c.re * c.re + c.im * c.im), -c.im / (c.re * c.re + c.im * c.im)} +} + +// 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 +pub fn (c Complex) pow(n f64) Complex { + r := math.pow(c.abs(), n) + angle := c.angle() + return Complex{r * math.cos(n * angle), r * math.sin(n * angle)} +} + +// Complex nth root +pub fn (c Complex) root(n f64) Complex { + return c.pow(1.0 / n) +} + +// Complex Exponential +// Using Euler's Identity +// Based on +// https://www.math.wisc.edu/~angenent/Free-Lecture-Notes/freecomplexnumbers.pdf +pub fn (c Complex) exp() Complex { + a := math.exp(c.re) + return Complex{a * math.cos(c.im), a * math.sin(c.im)} +} + +// Complex Natural Logarithm +// Based on +// http://www.chemistrylearning.com/logarithm-of-complex-number/ +pub fn (c Complex) ln() Complex { + return Complex{math.log(c.abs()), c.angle()} +} + +// Complex Log Base Complex +// Based on +// http://www.milefoot.com/math/complex/summaryops.htm +pub fn (c Complex) log(base Complex) Complex { + return base.ln().divide(c.ln()) +} + +// Complex Argument +// Based on +// http://mathworld.wolfram.com/ComplexArgument.html +pub fn (c Complex) arg() f64 { + return math.atan2(c.im, c.re) +} + +// Complex raised to Complex Power +// Based on +// http://mathworld.wolfram.com/ComplexExponentiation.html +pub fn (c Complex) cpow(p Complex) Complex { + a := c.arg() + b := math.pow(c.re, 2) + math.pow(c.im, 2) + d := p.re * a + (1.0 / 2) * p.im * math.log(b) + t1 := math.pow(b, p.re / 2) * math.exp(-p.im * a) + return Complex{t1 * math.cos(d), t1 * math.sin(d)} +} + +// Complex Sin +// Based on +// http://www.milefoot.com/math/complex/functionsofi.htm +pub fn (c Complex) sin() Complex { + return Complex{math.sin(c.re) * math.cosh(c.im), math.cos(c.re) * math.sinh(c.im)} +} + +// Complex Cosine +// Based on +// http://www.milefoot.com/math/complex/functionsofi.htm +pub fn (c Complex) cos() Complex { + return Complex{math.cos(c.re) * math.cosh(c.im), -(math.sin(c.re) * math.sinh(c.im))} +} + +// Complex Tangent +// Based on +// http://www.milefoot.com/math/complex/functionsofi.htm +pub fn (c Complex) tan() Complex { + return c.sin().divide(c.cos()) +} + +// Complex Cotangent +// Based on +// http://www.suitcaseofdreams.net/Trigonometric_Functions.htm +pub fn (c Complex) cot() Complex { + return c.cos().divide(c.sin()) +} + +// Complex Secant +// Based on +// http://www.suitcaseofdreams.net/Trigonometric_Functions.htm +pub fn (c Complex) sec() Complex { + return complex(1, 0).divide(c.cos()) +} + +// Complex Cosecant +// Based on +// http://www.suitcaseofdreams.net/Trigonometric_Functions.htm +pub fn (c Complex) csc() Complex { + return complex(1, 0).divide(c.sin()) +} + +// 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()) +} + +// Complex Arc Cotangent / Cotangent Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse_Functions.htm +pub fn (c Complex) acot() Complex { + return complex(1, 0).divide(c).atan() +} + +// Complex Arc Secant / Secant Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse_Functions.htm +pub fn (c Complex) asec() Complex { + return complex(1, 0).divide(c).acos() +} + +// Complex Arc Cosecant / Cosecant Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse_Functions.htm +pub fn (c Complex) acsc() Complex { + return complex(1, 0).divide(c).asin() +} + +// Complex Hyperbolic Sin +// Based on +// http://www.milefoot.com/math/complex/functionsofi.htm +pub fn (c Complex) sinh() Complex { + return Complex{math.cos(c.im) * math.sinh(c.re), math.sin(c.im) * math.cosh(c.re)} +} + +// Complex Hyperbolic Cosine +// Based on +// http://www.milefoot.com/math/complex/functionsofi.htm +pub fn (c Complex) cosh() Complex { + return Complex{math.cos(c.im) * math.cosh(c.re), math.sin(c.im) * math.sinh(c.re)} +} + +// 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()) +} + +// Complex Hyperbolic Cotangent +// Based on +// http://www.suitcaseofdreams.net/Hyperbolic_Functions.htm +pub fn (c Complex) coth() Complex { + return c.cosh().divide(c.sinh()) +} + +// Complex Hyperbolic Secant +// Based on +// http://www.suitcaseofdreams.net/Hyperbolic_Functions.htm +pub fn (c Complex) sech() Complex { + return complex(1, 0).divide(c.cosh()) +} + +// Complex Hyperbolic Cosecant +// Based on +// http://www.suitcaseofdreams.net/Hyperbolic_Functions.htm +pub fn (c Complex) csch() Complex { + return complex(1, 0).divide(c.sinh()) +} + +// 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 { + one := complex(1, 0) + if c.re < 1 { + return complex(1.0 / 2, 0).multiply(one.add(c).divide(one.subtract(c)).ln()) + } else { + return complex(1.0 / 2, 0).multiply(one.add(c).ln().subtract(one.subtract(c).ln())) + } +} + +// Complex Hyperbolic Arc Cotangent / Cotangent Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm +pub fn (c Complex) acoth() Complex { + one := complex(1, 0) + if c.re < 0 || c.re > 1 { + return complex(1.0 / 2, 0).multiply(c.add(one).divide(c.subtract(one)).ln()) + } else { + div := one.divide(c) + return complex(1.0 / 2, 0).multiply(one.add(div).ln().subtract(one.subtract(div).ln())) + } +} + +// Complex Hyperbolic Arc Secant / Secant Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm +// For certain scenarios, Result mismatch in crossverification with Wolfram Alpha - analysis pending +// pub fn (c Complex) asech() Complex { +// one := complex(1,0) +// if(c.re < -1.0) { +// return one.subtract( +// one.subtract( +// c.pow(2) +// ) +// .root(2) +// ) +// .divide(c) +// .ln() +// } +// else { +// return one.add( +// one.subtract( +// c.pow(2) +// ) +// .root(2) +// ) +// .divide(c) +// .ln() +// } +// } + +// Complex Hyperbolic Arc Cosecant / Cosecant Inverse +// Based on +// http://www.suitcaseofdreams.net/Inverse__Hyperbolic_Functions.htm +pub fn (c Complex) acsch() Complex { + one := complex(1, 0) + if c.re < 0 { + return one.subtract(one.add(c.pow(2)).root(2)).divide(c).ln() + } else { + return one.add(one.add(c.pow(2)).root(2)).divide(c).ln() + } +} + +// Complex Equals +pub fn (c1 Complex) equals(c2 Complex) bool { + return (c1.re == c2.re) && (c1.im == c2.im) +} diff --git a/v_windows/v/old/vlib/math/complex/complex_test.v b/v_windows/v/old/vlib/math/complex/complex_test.v new file mode 100644 index 0000000..ccd448e --- /dev/null +++ b/v_windows/v/old/vlib/math/complex/complex_test.v @@ -0,0 +1,797 @@ +import math +import math.complex as cmplx + +fn tst_res(str1 string, str2 string) bool { + if (math.abs(str1.f64() - str2.f64())) < 1e-5 { + return true + } + return false +} + +fn test_complex_addition() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c1 := cmplx.complex(0, -10) + mut c2 := cmplx.complex(-40, 8) + mut result := c1 + c2 + assert result.equals(cmplx.complex(-40, -2)) + c1 = cmplx.complex(-71, 2) + c2 = cmplx.complex(88, -12) + result = c1 + c2 + assert result.equals(cmplx.complex(17, -10)) + c1 = cmplx.complex(0, -30) + c2 = cmplx.complex(52, -30) + result = c1 + c2 + assert result.equals(cmplx.complex(52, -60)) + c1 = cmplx.complex(12, -9) + c2 = cmplx.complex(32, -6) + result = c1 + c2 + assert result.equals(cmplx.complex(44, -15)) +} + +fn test_complex_subtraction() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c1 := cmplx.complex(-8, 0) + mut c2 := cmplx.complex(6, 30) + mut result := c1 - c2 + assert result.equals(cmplx.complex(-14, -30)) + c1 = cmplx.complex(-19, 7) + c2 = cmplx.complex(29, 32) + result = c1 - c2 + assert result.equals(cmplx.complex(-48, -25)) + c1 = cmplx.complex(12, 0) + c2 = cmplx.complex(23, 13) + result = c1 - c2 + assert result.equals(cmplx.complex(-11, -13)) + c1 = cmplx.complex(-14, 3) + c2 = cmplx.complex(0, 14) + result = c1 - c2 + assert result.equals(cmplx.complex(-14, -11)) +} + +fn test_complex_multiplication() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c1 := cmplx.complex(1, 2) + mut c2 := cmplx.complex(1, -4) + mut result := c1 * c2 + assert result.equals(cmplx.complex(9, -2)) + c1 = cmplx.complex(-4, -4) + c2 = cmplx.complex(-5, -3) + result = c1 * c2 + assert result.equals(cmplx.complex(8, 32)) + c1 = cmplx.complex(4, 4) + c2 = cmplx.complex(-2, -5) + result = c1 * c2 + assert result.equals(cmplx.complex(12, -28)) + c1 = cmplx.complex(2, -2) + c2 = cmplx.complex(4, -4) + result = c1 * c2 + assert result.equals(cmplx.complex(0, -16)) +} + +fn test_complex_division() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c1 := cmplx.complex(-9, -6) + mut c2 := cmplx.complex(-3, -2) + mut result := c1 / c2 + assert result.equals(cmplx.complex(3, 0)) + c1 = cmplx.complex(-23, 11) + c2 = cmplx.complex(5, 1) + result = c1 / c2 + assert result.equals(cmplx.complex(-4, 3)) + c1 = cmplx.complex(8, -2) + c2 = cmplx.complex(-4, 1) + result = c1 / c2 + assert result.equals(cmplx.complex(-2, 0)) + c1 = cmplx.complex(11, 24) + c2 = cmplx.complex(-4, -1) + result = c1 / c2 + assert result.equals(cmplx.complex(-4, -5)) +} + +fn test_complex_conjugate() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c1 := cmplx.complex(0, 8) + mut result := c1.conjugate() + assert result.equals(cmplx.complex(0, -8)) + c1 = cmplx.complex(7, 3) + result = c1.conjugate() + assert result.equals(cmplx.complex(7, -3)) + c1 = cmplx.complex(2, 2) + result = c1.conjugate() + assert result.equals(cmplx.complex(2, -2)) + c1 = cmplx.complex(7, 0) + result = c1.conjugate() + assert result.equals(cmplx.complex(7, 0)) +} + +fn test_complex_equals() { + mut c1 := cmplx.complex(0, 8) + mut c2 := cmplx.complex(0, 8) + assert c1.equals(c2) + c1 = cmplx.complex(-3, 19) + c2 = cmplx.complex(-3, 19) + assert c1.equals(c2) +} + +fn test_complex_abs() { + mut c1 := cmplx.complex(3, 4) + assert c1.abs() == 5 + c1 = cmplx.complex(1, 2) + assert c1.abs() == math.sqrt(5) + assert c1.abs() == c1.conjugate().abs() + c1 = cmplx.complex(7, 0) + assert c1.abs() == 7 +} + +fn test_complex_angle() { + // Test is based on and verified from practice examples of Khan Academy + // https://www.khanacademy.org/math/precalculus/imaginary-and-complex-numbers + mut c := cmplx.complex(1, 0) + assert c.angle() * 180 / math.pi == 0 + c = cmplx.complex(1, 1) + assert c.angle() * 180 / math.pi == 45 + c = cmplx.complex(0, 1) + assert c.angle() * 180 / math.pi == 90 + c = cmplx.complex(-1, 1) + assert c.angle() * 180 / math.pi == 135 + c = cmplx.complex(-1, -1) + assert c.angle() * 180 / math.pi == -135 + cc := c.conjugate() + assert cc.angle() + c.angle() == 0 +} + +fn test_complex_addinv() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-5, -7) + mut result := c1.addinv() + assert result.equals(c2) + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(3, -4) + result = c1.addinv() + assert result.equals(c2) + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(1, 2) + result = c1.addinv() + assert result.equals(c2) +} + +fn test_complex_mulinv() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.067568, -0.094595) + mut result := c1.mulinv() + // Some issue with precision comparison in f64 using == operator hence serializing to string + println(c2.str()) + println(result.str()) + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.12, -0.16) + result = c1.mulinv() + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.2, 0.4) + result = c1.mulinv() + assert result.equals(c2) +} + +fn test_complex_mod() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut result := c1.mod() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(result.str(), '8.602325') + c1 = cmplx.complex(-3, 4) + result = c1.mod() + assert result == 5 + c1 = cmplx.complex(-1, -2) + result = c1.mod() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(result.str(), '2.236068') +} + +fn test_complex_pow() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-24.0, 70.0) + mut result := c1.pow(2) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(117, 44) + result = c1.pow(3) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-7, -24) + result = c1.pow(4) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_root() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(2.607904, 1.342074) + mut result := c1.root(2) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(1.264953, 1.150614) + result = c1.root(3) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(1.068059, -0.595482) + result = c1.root(4) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_exp() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(111.889015, 97.505457) + mut result := c1.exp() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.032543, -0.037679) + result = c1.exp() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.153092, -0.334512) + result = c1.exp() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_ln() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(2.152033, 0.950547) + mut result := c1.ln() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(1.609438, 2.214297) + result = c1.ln() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(0.804719, -2.034444) + result = c1.ln() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_arg() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(2.152033, 0.950547) + mut result := c1.arg() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(result.str(), '0.950547') + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(1.609438, 2.214297) + result = c1.arg() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(result.str(), '2.214297') + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(0.804719, -2.034444) + result = c1.arg() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(result.str(), '-2.034444') +} + +fn test_complex_log() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut b1 := cmplx.complex(-6, -2) + mut c2 := cmplx.complex(0.232873, -1.413175) + mut result := c1.log(b1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + b1 = cmplx.complex(3, -1) + c2 = cmplx.complex(0.152198, -0.409312) + result = c1.log(b1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + b1 = cmplx.complex(0, 9) + c2 = cmplx.complex(-0.298243, 1.197981) + result = c1.log(b1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_cpow() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut r1 := cmplx.complex(2, 2) + mut c2 := cmplx.complex(11.022341, -0.861785) + mut result := c1.cpow(r1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + r1 = cmplx.complex(-4, -2) + c2 = cmplx.complex(0.118303, 0.063148) + result = c1.cpow(r1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + r1 = cmplx.complex(8, -9) + c2 = cmplx.complex(-0.000000, 0.000007) + result = c1.cpow(r1) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_sin() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-525.794515, 155.536550) + mut result := c1.sin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-3.853738, -27.016813) + result = c1.sin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-3.165779, -1.959601) + result = c1.sin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_cos() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(155.536809, 525.793641) + mut result := c1.cos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-27.034946, 3.851153) + result = c1.cos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(2.032723, -3.051898) + result = c1.cos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_tan() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-0.000001, 1.000001) + mut result := c1.tan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(0.000187, 0.999356) + result = c1.tan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.033813, -1.014794) + result = c1.tan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_cot() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-0.000001, -0.999999) + mut result := c1.cot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(0.000188, -1.000644) + result = c1.cot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.032798, 0.984329) + result = c1.cot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_sec() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.000517, -0.001749) + mut result := c1.sec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.036253, -0.005164) + result = c1.sec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(0.151176, 0.226974) + result = c1.sec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_csc() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(-0.001749, -0.000517) + mut result := c1.csc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.005174, 0.036276) + result = c1.csc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.228375, 0.141363) + result = c1.csc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_asin() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.617064, 2.846289) + mut result := c1.asin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.633984, 2.305509) + result = c1.asin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.427079, -1.528571) + result = c1.asin() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_acos() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.953732, -2.846289) + mut result := c1.acos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(2.204780, -2.305509) + result = c1.acos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(1.997875, 1.528571) + result = c1.acos() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_atan() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(1.502727, 0.094441) + mut result := c1.atan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-1.448307, 0.158997) + result = c1.atan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-1.338973, -0.402359) + result = c1.atan() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_acot() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.068069, -0.094441) + mut result := c1.acot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.122489, -0.158997) + result = c1.acot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.231824, 0.402359) + result = c1.acot() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_asec() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(1.503480, 0.094668) + mut result := c1.asec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(1.689547, 0.160446) + result = c1.asec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(1.757114, -0.396568) + result = c1.asec() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_acsc() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.067317, -0.094668) + mut result := c1.acsc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.118751, -0.160446) + result = c1.acsc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.186318, 0.396568) + result = c1.acsc() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_sinh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(55.941968, 48.754942) + mut result := c1.sinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(6.548120, -7.619232) + result = c1.sinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(0.489056, -1.403119) + result = c1.sinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_cosh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(55.947047, 48.750515) + mut result := c1.cosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-6.580663, 7.581553) + result = c1.cosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.642148, 1.068607) + result = c1.cosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_tanh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.999988, 0.000090) + mut result := c1.tanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-1.000710, 0.004908) + result = c1.tanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-1.166736, 0.243458) + result = c1.tanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_coth() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(1.000012, -0.000090) + mut result := c1.coth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.999267, -0.004901) + result = c1.coth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.821330, -0.171384) + result = c1.coth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_sech() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.010160, -0.008853) + mut result := c1.sech() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.065294, -0.075225) + result = c1.sech() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.413149, -0.687527) + result = c1.sech() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_csch() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.010159, -0.008854) + mut result := c1.csch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(0.064877, 0.075490) + result = c1.csch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(0.221501, 0.635494) + result = c1.csch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_asinh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(2.844098, 0.947341) + mut result := c1.asinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-2.299914, 0.917617) + result = c1.asinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-1.469352, -1.063440) + result = c1.asinh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_acosh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(2.846289, 0.953732) + mut result := c1.acosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(2.305509, 2.204780) + result = c1.acosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(1.528571, -1.997875) + result = c1.acosh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_atanh() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.067066, 1.476056) + mut result := c1.atanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.117501, 1.409921) + result = c1.atanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.173287, -1.178097) + result = c1.atanh() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_acoth() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.067066, -0.094740) + mut result := c1.acoth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.117501, -0.160875) + result = c1.acoth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.173287, 0.392699) + result = c1.acoth() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +// fn test_complex_asech() { +// // Tests were also verified on Wolfram Alpha +// mut c1 := cmplx.complex(5,7) +// mut c2 := cmplx.complex(0.094668,-1.503480) +// mut result := c1.asech() +// // Some issue with precision comparison in f64 using == operator hence serializing to string +// assert result.str() == c2.str() +// c1 = cmplx.complex(-3,4) +// c2 = cmplx.complex(0.160446,-1.689547) +// result = c1.asech() +// // Some issue with precision comparison in f64 using == operator hence serializing to string +// assert result.str() c2.str() +// c1 = cmplx.complex(-1,-2) +// c2 = cmplx.complex(0.396568,1.757114) +// result = c1.asech() +// // Some issue with precision comparison in f64 using == operator hence serializing to string +// assert result.str() == c2.str() +// } + +fn test_complex_acsch() { + // Tests were also verified on Wolfram Alpha + mut c1 := cmplx.complex(5, 7) + mut c2 := cmplx.complex(0.067819, -0.094518) + mut result := c1.acsch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-3, 4) + c2 = cmplx.complex(-0.121246, -0.159507) + result = c1.acsch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() + c1 = cmplx.complex(-1, -2) + c2 = cmplx.complex(-0.215612, 0.401586) + result = c1.acsch() + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert result.str() == c2.str() +} + +fn test_complex_re_im() { + c := cmplx.complex(2.1, 9.05) + assert c.re == 2.1 + assert c.im == 9.05 +} diff --git a/v_windows/v/old/vlib/math/const.v b/v_windows/v/old/vlib/math/const.v new file mode 100644 index 0000000..5a831a9 --- /dev/null +++ b/v_windows/v/old/vlib/math/const.v @@ -0,0 +1,49 @@ +// Copyright (c) 2019-2021 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 + +pub const ( + e = 2.71828182845904523536028747135266249775724709369995957496696763 + pi = 3.14159265358979323846264338327950288419716939937510582097494459 + phi = 1.61803398874989484820458683436563811772030917980576286213544862 + tau = 6.28318530717958647692528676655900576839433879875021164194988918 + sqrt2 = 1.41421356237309504880168872420969807856967187537694807317667974 + sqrt_e = 1.64872127070012814684865078781416357165377610071014801157507931 + sqrt_pi = 1.77245385090551602729816748334114518279754945612238712821380779 + sqrt_tau = 2.50662827463100050241576528481104525300698674060993831662992357 + sqrt_phi = 1.27201964951406896425242246173749149171560804184009624861664038 + ln2 = 0.693147180559945309417232121458176568075500134360255254120680009 + log2_e = 1.0 / ln2 + ln10 = 2.30258509299404568401799145468436420760110148862877297603332790 + log10_e = 1.0 / ln10 +) + +// Floating-point limit values +// max is the largest finite value representable by the type. +// smallest_non_zero is the smallest positive, non-zero value representable by the type. +pub const ( + max_f32 = 3.40282346638528859811704183484516925440e+38 // 2**127 * (2**24 - 1) / 2**23 + smallest_non_zero_f32 = 1.401298464324817070923729583289916131280e-45 // 1 / 2**(127 - 1 + 23) + max_f64 = 1.797693134862315708145274237317043567981e+308 // 2**1023 * (2**53 - 1) / 2**52 + smallest_non_zero_f64 = 4.940656458412465441765687928682213723651e-324 // 1 / 2**(1023 - 1 + 52) +) + +// Integer limit values +pub const ( + max_i8 = 127 + min_i8 = -128 + max_i16 = 32767 + min_i16 = -32768 + max_i32 = 2147483647 + min_i32 = -2147483648 + // -9223372036854775808 is wrong because C compilers parse literal values + // without sign first, and 9223372036854775808 overflows i64, hence the + // consecutive subtraction by 1 + min_i64 = i64(-9223372036854775807 - 1) + max_i64 = i64(9223372036854775807) + max_u8 = 255 + max_u16 = 65535 + max_u32 = u32(4294967295) + max_u64 = u64(18446744073709551615) +) diff --git a/v_windows/v/old/vlib/math/factorial/factorial.v b/v_windows/v/old/vlib/math/factorial/factorial.v new file mode 100644 index 0000000..9668d5d --- /dev/null +++ b/v_windows/v/old/vlib/math/factorial/factorial.v @@ -0,0 +1,80 @@ +// Copyright (c) 2019-2021 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 created by Ulises Jeremias Cornejo Fandos based on +// the definitions provided in https://scientificc.github.io/cmathl/ + +module factorial + +import math + +// factorial calculates the factorial of the provided value. +pub fn factorial(n f64) f64 { + // For a large postive argument (n >= FACTORIALS.len) return max_f64 + + if n >= factorials_table.len { + return math.max_f64 + } + + // Otherwise return n!. + if n == f64(i64(n)) && n >= 0.0 { + return factorials_table[i64(n)] + } + + return math.gamma(n + 1.0) +} + +// log_factorial calculates the log-factorial of the provided value. +pub fn log_factorial(n f64) f64 { + // For a large postive argument (n < 0) return max_f64 + + if n < 0 { + return -math.max_f64 + } + + // If n < N then return ln(n!). + + if n != f64(i64(n)) { + return math.log_gamma(n + 1) + } else if n < log_factorials_table.len { + return log_factorials_table[i64(n)] + } + + // Otherwise return asymptotic expansion of ln(n!). + + return log_factorial_asymptotic_expansion(int(n)) +} + +fn log_factorial_asymptotic_expansion(n int) f64 { + m := 6 + mut term := []f64{} + xx := f64((n + 1) * (n + 1)) + mut xj := f64(n + 1) + + log_factorial := log_sqrt_2pi - xj + (xj - 0.5) * math.log(xj) + + mut i := 0 + + for i = 0; i < m; i++ { + term << b_numbers[i] / xj + xj *= xx + } + + mut sum := term[m - 1] + + for i = m - 2; i >= 0; i-- { + if math.abs(sum) <= math.abs(term[i]) { + break + } + + sum = term[i] + } + + for i >= 0 { + sum += term[i] + i-- + } + + return log_factorial + sum +} diff --git a/v_windows/v/old/vlib/math/factorial/factorial_tables.v b/v_windows/v/old/vlib/math/factorial/factorial_tables.v new file mode 100644 index 0000000..a669b09 --- /dev/null +++ b/v_windows/v/old/vlib/math/factorial/factorial_tables.v @@ -0,0 +1,375 @@ +// Copyright (c) 2019-2021 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 factorial + +const ( + log_sqrt_2pi = 9.18938533204672741780329736e-1 + + b_numbers = [ + /* + Bernoulli numbers B(2),B(4),B(6),...,B(20). Only B(2),...,B(10) currently + * used. + */ + f64(1.0 / (6.0 * 2.0 * 1.0)), + -1.0 / (30.0 * 4.0 * 3.0), + 1.0 / (42.0 * 6.0 * 5.0), + -1.0 / (30.0 * 8.0 * 7.0), + 5.0 / (66.0 * 10.0 * 9.0), + -691.0 / (2730.0 * 12.0 * 11.0), + 7.0 / (6.0 * 14.0 * 13.0), + -3617.0 / (510.0 * 16.0 * 15.0), + 43867.0 / (796.0 * 18.0 * 17.0), + -174611.0 / (330.0 * 20.0 * 19.0), + ] + + factorials_table = [ + f64(1.000000000000000000000e+0), /* 0! */ + 1.000000000000000000000e+0, /* 1! */ + 2.000000000000000000000e+0, /* 2! */ + 6.000000000000000000000e+0, /* 3! */ + 2.400000000000000000000e+1, /* 4! */ + 1.200000000000000000000e+2, /* 5! */ + 7.200000000000000000000e+2, /* 6! */ + 5.040000000000000000000e+3, /* 7! */ + 4.032000000000000000000e+4, /* 8! */ + 3.628800000000000000000e+5, /* 9! */ + 3.628800000000000000000e+6, /* 10! */ + 3.991680000000000000000e+7, /* 11! */ + 4.790016000000000000000e+8, /* 12! */ + 6.227020800000000000000e+9, /* 13! */ + 8.717829120000000000000e+10, /* 14! */ + 1.307674368000000000000e+12, /* 15! */ + 2.092278988800000000000e+13, /* 16! */ + 3.556874280960000000000e+14, /* 17! */ + 6.402373705728000000000e+15, /* 18! */ + 1.216451004088320000000e+17, /* 19! */ + 2.432902008176640000000e+18, /* 20! */ + 5.109094217170944000000e+19, /* 21! */ + 1.124000727777607680000e+21, /* 22! */ + 2.585201673888497664000e+22, /* 23! */ + 6.204484017332394393600e+23, /* 24! */ + 1.551121004333098598400e+25, /* 25! */ + 4.032914611266056355840e+26, /* 26! */ + 1.088886945041835216077e+28, /* 27! */ + 3.048883446117138605015e+29, /* 28! */ + 8.841761993739701954544e+30, /* 29! */ + 2.652528598121910586363e+32, /* 30! */ + 8.222838654177922817726e+33, /* 31! */ + 2.631308369336935301672e+35, /* 32! */ + 8.683317618811886495518e+36, /* 33! */ + 2.952327990396041408476e+38, /* 34! */ + 1.033314796638614492967e+40, /* 35! */ + 3.719933267899012174680e+41, /* 36! */ + 1.376375309122634504632e+43, /* 37! */ + 5.230226174666011117600e+44, /* 38! */ + 2.039788208119744335864e+46, /* 39! */ + 8.159152832478977343456e+47, /* 40! */ + 3.345252661316380710817e+49, /* 41! */ + 1.405006117752879898543e+51, /* 42! */ + 6.041526306337383563736e+52, /* 43! */ + 2.658271574788448768044e+54, /* 44! */ + 1.196222208654801945620e+56, /* 45! */ + 5.502622159812088949850e+57, /* 46! */ + 2.586232415111681806430e+59, /* 47! */ + 1.241391559253607267086e+61, /* 48! */ + 6.082818640342675608723e+62, /* 49! */ + 3.041409320171337804361e+64, /* 50! */ + 1.551118753287382280224e+66, /* 51! */ + 8.065817517094387857166e+67, /* 52! */ + 4.274883284060025564298e+69, /* 53! */ + 2.308436973392413804721e+71, /* 54! */ + 1.269640335365827592597e+73, /* 55! */ + 7.109985878048634518540e+74, /* 56! */ + 4.052691950487721675568e+76, /* 57! */ + 2.350561331282878571829e+78, /* 58! */ + 1.386831185456898357379e+80, /* 59! */ + 8.320987112741390144276e+81, /* 60! */ + 5.075802138772247988009e+83, /* 61! */ + 3.146997326038793752565e+85, /* 62! */ + 1.982608315404440064116e+87, /* 63! */ + 1.268869321858841641034e+89, /* 64! */ + 8.247650592082470666723e+90, /* 65! */ + 5.443449390774430640037e+92, /* 66! */ + 3.647111091818868528825e+94, /* 67! */ + 2.480035542436830599601e+96, /* 68! */ + 1.711224524281413113725e+98, /* 69! */ + 1.197857166996989179607e+100, /* 70! */ + 8.504785885678623175212e+101, /* 71! */ + 6.123445837688608686152e+103, /* 72! */ + 4.470115461512684340891e+105, /* 73! */ + 3.307885441519386412260e+107, /* 74! */ + 2.480914081139539809195e+109, /* 75! */ + 1.885494701666050254988e+111, /* 76! */ + 1.451830920282858696341e+113, /* 77! */ + 1.132428117820629783146e+115, /* 78! */ + 8.946182130782975286851e+116, /* 79! */ + 7.156945704626380229481e+118, /* 80! */ + 5.797126020747367985880e+120, /* 81! */ + 4.753643337012841748421e+122, /* 82! */ + 3.945523969720658651190e+124, /* 83! */ + 3.314240134565353266999e+126, /* 84! */ + 2.817104114380550276949e+128, /* 85! */ + 2.422709538367273238177e+130, /* 86! */ + 2.107757298379527717214e+132, /* 87! */ + 1.854826422573984391148e+134, /* 88! */ + 1.650795516090846108122e+136, /* 89! */ + 1.485715964481761497310e+138, /* 90! */ + 1.352001527678402962552e+140, /* 91! */ + 1.243841405464130725548e+142, /* 92! */ + 1.156772507081641574759e+144, /* 93! */ + 1.087366156656743080274e+146, /* 94! */ + 1.032997848823905926260e+148, /* 95! */ + 9.916779348709496892096e+149, /* 96! */ + 9.619275968248211985333e+151, /* 97! */ + 9.426890448883247745626e+153, /* 98! */ + 9.332621544394415268170e+155, /* 99! */ + 9.332621544394415268170e+157, /* 100! */ + 9.425947759838359420852e+159, /* 101! */ + 9.614466715035126609269e+161, /* 102! */ + 9.902900716486180407547e+163, /* 103! */ + 1.029901674514562762385e+166, /* 104! */ + 1.081396758240290900504e+168, /* 105! */ + 1.146280563734708354534e+170, /* 106! */ + 1.226520203196137939352e+172, /* 107! */ + 1.324641819451828974500e+174, /* 108! */ + 1.443859583202493582205e+176, /* 109! */ + 1.588245541522742940425e+178, /* 110! */ + 1.762952551090244663872e+180, /* 111! */ + 1.974506857221074023537e+182, /* 112! */ + 2.231192748659813646597e+184, /* 113! */ + 2.543559733472187557120e+186, /* 114! */ + 2.925093693493015690688e+188, /* 115! */ + 3.393108684451898201198e+190, /* 116! */ + 3.969937160808720895402e+192, /* 117! */ + 4.684525849754290656574e+194, /* 118! */ + 5.574585761207605881323e+196, /* 119! */ + 6.689502913449127057588e+198, /* 120! */ + 8.094298525273443739682e+200, /* 121! */ + 9.875044200833601362412e+202, /* 122! */ + 1.214630436702532967577e+205, /* 123! */ + 1.506141741511140879795e+207, /* 124! */ + 1.882677176888926099744e+209, /* 125! */ + 2.372173242880046885677e+211, /* 126! */ + 3.012660018457659544810e+213, /* 127! */ + 3.856204823625804217357e+215, /* 128! */ + 4.974504222477287440390e+217, /* 129! */ + 6.466855489220473672507e+219, /* 130! */ + 8.471580690878820510985e+221, /* 131! */ + 1.118248651196004307450e+224, /* 132! */ + 1.487270706090685728908e+226, /* 133! */ + 1.992942746161518876737e+228, /* 134! */ + 2.690472707318050483595e+230, /* 135! */ + 3.659042881952548657690e+232, /* 136! */ + 5.012888748274991661035e+234, /* 137! */ + 6.917786472619488492228e+236, /* 138! */ + 9.615723196941089004197e+238, /* 139! */ + 1.346201247571752460588e+241, /* 140! */ + 1.898143759076170969429e+243, /* 141! */ + 2.695364137888162776589e+245, /* 142! */ + 3.854370717180072770522e+247, /* 143! */ + 5.550293832739304789551e+249, /* 144! */ + 8.047926057471991944849e+251, /* 145! */ + 1.174997204390910823948e+254, /* 146! */ + 1.727245890454638911203e+256, /* 147! */ + 2.556323917872865588581e+258, /* 148! */ + 3.808922637630569726986e+260, /* 149! */ + 5.713383956445854590479e+262, /* 150! */ + 8.627209774233240431623e+264, /* 151! */ + 1.311335885683452545607e+267, /* 152! */ + 2.006343905095682394778e+269, /* 153! */ + 3.089769613847350887959e+271, /* 154! */ + 4.789142901463393876336e+273, /* 155! */ + 7.471062926282894447084e+275, /* 156! */ + 1.172956879426414428192e+278, /* 157! */ + 1.853271869493734796544e+280, /* 158! */ + 2.946702272495038326504e+282, /* 159! */ + 4.714723635992061322407e+284, /* 160! */ + 7.590705053947218729075e+286, /* 161! */ + 1.229694218739449434110e+289, /* 162! */ + 2.004401576545302577600e+291, /* 163! */ + 3.287218585534296227263e+293, /* 164! */ + 5.423910666131588774984e+295, /* 165! */ + 9.003691705778437366474e+297, /* 166! */ + 1.503616514864999040201e+300, /* 167! */ + 2.526075744973198387538e+302, /* 168! */ + 4.269068009004705274939e+304, /* 169! */ + 7.257415615307998967397e+306, /* 170! */ + ] + + log_factorials_table = [ + f64(0.000000000000000000000e+0), /* 0! */ + 0.000000000000000000000e+0, /* 1! */ + 6.931471805599453094172e-1, /* 2! */ + 1.791759469228055000812e+0, /* 3! */ + 3.178053830347945619647e+0, /* 4! */ + 4.787491742782045994248e+0, /* 5! */ + 6.579251212010100995060e+0, /* 6! */ + 8.525161361065414300166e+0, /* 7! */ + 1.060460290274525022842e+1, /* 8! */ + 1.280182748008146961121e+1, /* 9! */ + 1.510441257307551529523e+1, /* 10! */ + 1.750230784587388583929e+1, /* 11! */ + 1.998721449566188614952e+1, /* 12! */ + 2.255216385312342288557e+1, /* 13! */ + 2.519122118273868150009e+1, /* 14! */ + 2.789927138384089156609e+1, /* 15! */ + 3.067186010608067280376e+1, /* 16! */ + 3.350507345013688888401e+1, /* 17! */ + 3.639544520803305357622e+1, /* 18! */ + 3.933988418719949403622e+1, /* 19! */ + 4.233561646075348502966e+1, /* 20! */ + 4.538013889847690802616e+1, /* 21! */ + 4.847118135183522387964e+1, /* 22! */ + 5.160667556776437357045e+1, /* 23! */ + 5.478472939811231919009e+1, /* 24! */ + 5.800360522298051993929e+1, /* 25! */ + 6.126170176100200198477e+1, /* 26! */ + 6.455753862700633105895e+1, /* 27! */ + 6.788974313718153498289e+1, /* 28! */ + 7.125703896716800901007e+1, /* 29! */ + 7.465823634883016438549e+1, /* 30! */ + 7.809222355331531063142e+1, /* 31! */ + 8.155795945611503717850e+1, /* 32! */ + 8.505446701758151741396e+1, /* 33! */ + 8.858082754219767880363e+1, /* 34! */ + 9.213617560368709248333e+1, /* 35! */ + 9.571969454214320248496e+1, /* 36! */ + 9.933061245478742692933e+1, /* 37! */ + 1.029681986145138126988e+2, /* 38! */ + 1.066317602606434591262e+2, /* 39! */ + 1.103206397147573954291e+2, /* 40! */ + 1.140342117814617032329e+2, /* 41! */ + 1.177718813997450715388e+2, /* 42! */ + 1.215330815154386339623e+2, /* 43! */ + 1.253172711493568951252e+2, /* 44! */ + 1.291239336391272148826e+2, /* 45! */ + 1.329525750356163098828e+2, /* 46! */ + 1.368027226373263684696e+2, /* 47! */ + 1.406739236482342593987e+2, /* 48! */ + 1.445657439463448860089e+2, /* 49! */ + 1.484777669517730320675e+2, /* 50! */ + 1.524095925844973578392e+2, /* 51! */ + 1.563608363030787851941e+2, /* 52! */ + 1.603311282166309070282e+2, /* 53! */ + 1.643201122631951814118e+2, /* 54! */ + 1.683274454484276523305e+2, /* 55! */ + 1.723527971391628015638e+2, /* 56! */ + 1.763958484069973517152e+2, /* 57! */ + 1.804562914175437710518e+2, /* 58! */ + 1.845338288614494905025e+2, /* 59! */ + 1.886281734236715911873e+2, /* 60! */ + 1.927390472878449024360e+2, /* 61! */ + 1.968661816728899939914e+2, /* 62! */ + 2.010093163992815266793e+2, /* 63! */ + 2.051681994826411985358e+2, /* 64! */ + 2.093425867525368356464e+2, /* 65! */ + 2.135322414945632611913e+2, /* 66! */ + 2.177369341139542272510e+2, /* 67! */ + 2.219564418191303339501e+2, /* 68! */ + 2.261905483237275933323e+2, /* 69! */ + 2.304390435657769523214e+2, /* 70! */ + 2.347017234428182677427e+2, /* 71! */ + 2.389783895618343230538e+2, /* 72! */ + 2.432688490029827141829e+2, /* 73! */ + 2.475729140961868839366e+2, /* 74! */ + 2.518904022097231943772e+2, /* 75! */ + 2.562211355500095254561e+2, /* 76! */ + 2.605649409718632093053e+2, /* 77! */ + 2.649216497985528010421e+2, /* 78! */ + 2.692910976510198225363e+2, /* 79! */ + 2.736731242856937041486e+2, /* 80! */ + 2.780675734403661429141e+2, /* 81! */ + 2.824742926876303960274e+2, /* 82! */ + 2.868931332954269939509e+2, /* 83! */ + 2.913239500942703075662e+2, /* 84! */ + 2.957666013507606240211e+2, /* 85! */ + 3.002209486470141317540e+2, /* 86! */ + 3.046868567656687154726e+2, /* 87! */ + 3.091641935801469219449e+2, /* 88! */ + 3.136528299498790617832e+2, /* 89! */ + 3.181526396202093268500e+2, /* 90! */ + 3.226634991267261768912e+2, /* 91! */ + 3.271852877037752172008e+2, /* 92! */ + 3.317178871969284731381e+2, /* 93! */ + 3.362611819791984770344e+2, /* 94! */ + 3.408150588707990178690e+2, /* 95! */ + 3.453794070622668541074e+2, /* 96! */ + 3.499541180407702369296e+2, /* 97! */ + 3.545390855194408088492e+2, /* 98! */ + 3.591342053695753987760e+2, /* 99! */ + 3.637393755555634901441e+2, /* 100! */ + 3.683544960724047495950e+2, /* 101! */ + 3.729794688856890206760e+2, /* 102! */ + 3.776141978739186564468e+2, /* 103! */ + 3.822585887730600291111e+2, /* 104! */ + 3.869125491232175524822e+2, /* 105! */ + 3.915759882173296196258e+2, /* 106! */ + 3.962488170517915257991e+2, /* 107! */ + 4.009309482789157454921e+2, /* 108! */ + 4.056222961611448891925e+2, /* 109! */ + 4.103227765269373054205e+2, /* 110! */ + 4.150323067282496395563e+2, /* 111! */ + 4.197508055995447340991e+2, /* 112! */ + 4.244781934182570746677e+2, /* 113! */ + 4.292143918666515701285e+2, /* 114! */ + 4.339593239950148201939e+2, /* 115! */ + 4.387129141861211848399e+2, /* 116! */ + 4.434750881209189409588e+2, /* 117! */ + 4.482457727453846057188e+2, /* 118! */ + 4.530248962384961351041e+2, /* 119! */ + 4.578123879812781810984e+2, /* 120! */ + 4.626081785268749221865e+2, /* 121! */ + 4.674121995716081787447e+2, /* 122! */ + 4.722243839269805962399e+2, /* 123! */ + 4.770446654925856331047e+2, /* 124! */ + 4.818729792298879342285e+2, /* 125! */ + 4.867092611368394122258e+2, /* 126! */ + 4.915534482232980034989e+2, /* 127! */ + 4.964054784872176206648e+2, /* 128! */ + 5.012652908915792927797e+2, /* 129! */ + 5.061328253420348751997e+2, /* 130! */ + 5.110080226652360267439e+2, /* 131! */ + 5.158908245878223975982e+2, /* 132! */ + 5.207811737160441513633e+2, /* 133! */ + 5.256790135159950627324e+2, /* 134! */ + 5.305842882944334921812e+2, /* 135! */ + 5.354969431801695441897e+2, /* 136! */ + 5.404169241059976691050e+2, /* 137! */ + 5.453441777911548737966e+2, /* 138! */ + 5.502786517242855655538e+2, /* 139! */ + 5.552202941468948698523e+2, /* 140! */ + 5.601690540372730381305e+2, /* 141! */ + 5.651248810948742988613e+2, /* 142! */ + 5.700877257251342061414e+2, /* 143! */ + 5.750575390247102067619e+2, /* 144! */ + 5.800342727671307811636e+2, /* 145! */ + 5.850178793888391176022e+2, /* 146! */ + 5.900083119756178539038e+2, /* 147! */ + 5.950055242493819689670e+2, /* 148! */ + 6.000094705553274281080e+2, /* 149! */ + 6.050201058494236838580e+2, /* 150! */ + 6.100373856862386081868e+2, /* 151! */ + 6.150612662070848845750e+2, /* 152! */ + 6.200917041284773200381e+2, /* 153! */ + 6.251286567308909491967e+2, /* 154! */ + 6.301720818478101958172e+2, /* 155! */ + 6.352219378550597328635e+2, /* 156! */ + 6.402781836604080409209e+2, /* 157! */ + 6.453407786934350077245e+2, /* 158! */ + 6.504096828956552392500e+2, /* 159! */ + 6.554848567108890661717e+2, /* 160! */ + 6.605662610758735291676e+2, /* 161! */ + 6.656538574111059132426e+2, /* 162! */ + 6.707476076119126755767e+2, /* 163! */ + 6.758474740397368739994e+2, /* 164! */ + 6.809534195136374546094e+2, /* 165! */ + 6.860654073019939978423e+2, /* 166! */ + 6.911834011144107529496e+2, /* 167! */ + 6.963073650938140118743e+2, /* 168! */ + 7.014372638087370853465e+2, /* 169! */ + 7.065730622457873471107e+2, /* 170! */ + 7.117147258022900069535e+2, /* 171! */ + ] +) diff --git a/v_windows/v/old/vlib/math/factorial/factorial_test.v b/v_windows/v/old/vlib/math/factorial/factorial_test.v new file mode 100644 index 0000000..6c2b575 --- /dev/null +++ b/v_windows/v/old/vlib/math/factorial/factorial_test.v @@ -0,0 +1,14 @@ +import math +import math.factorial as fact + +fn test_factorial() { + assert fact.factorial(12) == 479001600 + assert fact.factorial(5) == 120 + assert fact.factorial(0) == 1 +} + +fn test_log_factorial() { + assert fact.log_factorial(12) == math.log(479001600) + assert fact.log_factorial(5) == math.log(120) + assert fact.log_factorial(0) == math.log(1) +} diff --git a/v_windows/v/old/vlib/math/fractions/approximations.v b/v_windows/v/old/vlib/math/fractions/approximations.v new file mode 100644 index 0000000..dd4f855 --- /dev/null +++ b/v_windows/v/old/vlib/math/fractions/approximations.v @@ -0,0 +1,119 @@ +// Copyright (c) 2019-2021 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 fractions + +import math + +const ( + default_eps = 1.0e-4 + max_iterations = 50 + zero = fraction(0, 1) +) + +// ------------------------------------------------------------------------ +// Unwrapped evaluation methods for fast evaluation of continued fractions. +// ------------------------------------------------------------------------ +// We need these functions because the evaluation of continued fractions +// always has to be done from the end. Also, the numerator-denominator pairs +// are generated from front to end. This means building a result from a +// previous one isn't possible. So we need unrolled versions to ensure that +// we don't take too much of a performance penalty by calling eval_cf +// several times. +// ------------------------------------------------------------------------ +// eval_1 returns the result of evaluating a continued fraction series of length 1 +fn eval_1(whole i64, d []i64) Fraction { + return fraction(whole * d[0] + 1, d[0]) +} + +// eval_2 returns the result of evaluating a continued fraction series of length 2 +fn eval_2(whole i64, d []i64) Fraction { + den := d[0] * d[1] + 1 + return fraction(whole * den + d[1], den) +} + +// eval_3 returns the result of evaluating a continued fraction series of length 3 +fn eval_3(whole i64, d []i64) Fraction { + d1d2_plus_n2 := d[1] * d[2] + 1 + den := d[0] * d1d2_plus_n2 + d[2] + return fraction(whole * den + d1d2_plus_n2, den) +} + +// eval_cf evaluates a continued fraction series and returns a Fraction. +fn eval_cf(whole i64, den []i64) Fraction { + count := den.len + // Offload some small-scale calculations + // to dedicated functions + match count { + 1 { + return eval_1(whole, den) + } + 2 { + return eval_2(whole, den) + } + 3 { + return eval_3(whole, den) + } + else { + last := count - 1 + mut n := i64(1) + mut d := den[last] + // The calculations are done from back to front + for index := count - 2; index >= 0; index-- { + t := d + d = den[index] * d + n + n = t + } + return fraction(d * whole + n, d) + } + } +} + +// approximate returns a Fraction that approcimates the given value to +// within the default epsilon value (1.0e-4). This means the result will +// be accurate to 3 places after the decimal. +pub fn approximate(val f64) Fraction { + return approximate_with_eps(val, fractions.default_eps) +} + +// approximate_with_eps returns a Fraction +pub fn approximate_with_eps(val f64, eps f64) Fraction { + if val == 0.0 { + return fractions.zero + } + if eps < 0.0 { + panic('Epsilon value cannot be negative.') + } + if math.fabs(val) > math.max_i64 { + panic('Value out of range.') + } + // The integer part is separated first. Then we process the fractional + // part to generate numerators and denominators in tandem. + whole := i64(val) + mut frac := val - f64(whole) + // Quick exit for integers + if frac == 0.0 { + return fraction(whole, 1) + } + mut d := []i64{} + mut partial := fractions.zero + // We must complete the approximation within the maximum number of + // itertations allowed. If we can't panic. + // Empirically tested: the hardest constant to approximate is the + // golden ratio (math.phi) and for f64s, it only needs 38 iterations. + for _ in 0 .. fractions.max_iterations { + // We calculate the reciprocal. That's why the numerator is + // always 1. + frac = 1.0 / frac + den := i64(frac) + d << den + // eval_cf is called often so it needs to be performant + partial = eval_cf(whole, d) + // Check if we're done + if math.fabs(val - partial.f64()) < eps { + return partial + } + frac -= f64(den) + } + panic("Couldn't converge. Please create an issue on https://github.com/vlang/v") +} diff --git a/v_windows/v/old/vlib/math/fractions/approximations_test.v b/v_windows/v/old/vlib/math/fractions/approximations_test.v new file mode 100644 index 0000000..5ee92bf --- /dev/null +++ b/v_windows/v/old/vlib/math/fractions/approximations_test.v @@ -0,0 +1,189 @@ +// Copyright (c) 2019-2021 Alexander Medvednikov. All rights reserved. +// Use of this source code is governed by an MIT license +// that can be found in the LICENSE file. +import math.fractions +import math + +fn test_half() { + float_val := 0.5 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(1, 2)) +} + +fn test_third() { + float_val := 1.0 / 3.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(1, 3)) +} + +fn test_minus_one_twelfth() { + float_val := -1.0 / 12.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(-1, 12)) +} + +fn test_zero() { + float_val := 0.0 + println('Pre') + fract_val := fractions.approximate(float_val) + println('Post') + assert fract_val.equals(fractions.fraction(0, 1)) +} + +fn test_minus_one() { + float_val := -1.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(-1, 1)) +} + +fn test_thirty_three() { + float_val := 33.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(33, 1)) +} + +fn test_millionth() { + float_val := 1.0 / 1000000.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(1, 1000000)) +} + +fn test_minus_27_by_57() { + float_val := -27.0 / 57.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(-27, 57)) +} + +fn test_29_by_104() { + float_val := 29.0 / 104.0 + fract_val := fractions.approximate(float_val) + assert fract_val.equals(fractions.fraction(29, 104)) +} + +fn test_140710_232() { + float_val := 140710.232 + fract_val := fractions.approximate(float_val) + // Approximation will match perfectly for upto 3 places after the decimal + // The result will be within default_eps of original value + assert fract_val.f64() == float_val +} + +fn test_pi_1_digit() { + assert fractions.approximate_with_eps(math.pi, 5.0e-2).equals(fractions.fraction(22, + 7)) +} + +fn test_pi_2_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-3).equals(fractions.fraction(22, + 7)) +} + +fn test_pi_3_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-4).equals(fractions.fraction(333, + 106)) +} + +fn test_pi_4_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-5).equals(fractions.fraction(355, + 113)) +} + +fn test_pi_5_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-6).equals(fractions.fraction(355, + 113)) +} + +fn test_pi_6_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-7).equals(fractions.fraction(355, + 113)) +} + +fn test_pi_7_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-8).equals(fractions.fraction(103993, + 33102)) +} + +fn test_pi_8_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-9).equals(fractions.fraction(103993, + 33102)) +} + +fn test_pi_9_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-10).equals(fractions.fraction(104348, + 33215)) +} + +fn test_pi_10_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-11).equals(fractions.fraction(312689, + 99532)) +} + +fn test_pi_11_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-12).equals(fractions.fraction(1146408, + 364913)) +} + +fn test_pi_12_digits() { + assert fractions.approximate_with_eps(math.pi, 5.0e-13).equals(fractions.fraction(4272943, + 1360120)) +} + +fn test_phi_1_digit() { + assert fractions.approximate_with_eps(math.phi, 5.0e-2).equals(fractions.fraction(5, + 3)) +} + +fn test_phi_2_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-3).equals(fractions.fraction(21, + 13)) +} + +fn test_phi_3_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-4).equals(fractions.fraction(55, + 34)) +} + +fn test_phi_4_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-5).equals(fractions.fraction(233, + 144)) +} + +fn test_phi_5_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-6).equals(fractions.fraction(610, + 377)) +} + +fn test_phi_6_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-7).equals(fractions.fraction(1597, + 987)) +} + +fn test_phi_7_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-8).equals(fractions.fraction(6765, + 4181)) +} + +fn test_phi_8_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-9).equals(fractions.fraction(17711, + 10946)) +} + +fn test_phi_9_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-10).equals(fractions.fraction(75025, + 46368)) +} + +fn test_phi_10_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-11).equals(fractions.fraction(196418, + 121393)) +} + +fn test_phi_11_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-12).equals(fractions.fraction(514229, + 317811)) +} + +fn test_phi_12_digits() { + assert fractions.approximate_with_eps(math.phi, 5.0e-13).equals(fractions.fraction(2178309, + 1346269)) +} diff --git a/v_windows/v/old/vlib/math/fractions/fraction.v b/v_windows/v/old/vlib/math/fractions/fraction.v new file mode 100644 index 0000000..2e0b7bd --- /dev/null +++ b/v_windows/v/old/vlib/math/fractions/fraction.v @@ -0,0 +1,259 @@ +// Copyright (c) 2019-2021 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 fractions + +import math +import math.bits + +// Fraction Struct +// --------------- +// A Fraction has a numerator (n) and a denominator (d). If the user uses +// the helper functions in this module, then the following are guaranteed: +// 1. If the user provides n and d with gcd(n, d) > 1, the fraction will +// not be reduced automatically. +// 2. d cannot be set to zero. The factory function will panic. +// 3. If provided d is negative, it will be made positive. n will change as well. +struct Fraction { + n i64 + d i64 +pub: + is_reduced bool +} + +// A factory function for creating a Fraction, adds a boundary condition +// to ensure that the denominator is non-zero. It automatically converts +// the negative denominator to positive and adjusts the numerator. +// NOTE: Fractions created are not reduced by default. +pub fn fraction(n i64, d i64) Fraction { + if d == 0 { + panic('Denominator cannot be zero') + } + // The denominator is always guaranteed to be positive (and non-zero). + if d < 0 { + return fraction(-n, -d) + } + return Fraction{ + n: n + d: d + is_reduced: math.gcd(n, d) == 1 + } +} + +// To String method +pub fn (f Fraction) str() string { + return '$f.n/$f.d' +} + +// +// + ---------------------+ +// | Arithmetic functions.| +// + ---------------------+ +// +// These are implemented from Knuth, TAOCP Vol 2. Section 4.5 +// +// Returns a correctly reduced result for both addition and subtraction +// NOTE: requires reduced inputs +fn general_addition_result(f1 Fraction, f2 Fraction, addition bool) Fraction { + d1 := math.gcd(f1.d, f2.d) + // d1 happens to be 1 around 600/(pi)^2 or 61 percent of the time (Theorem 4.5.2D) + if d1 == 1 { + num1n2d := f1.n * f2.d + num1d2n := f1.d * f2.n + n := if addition { num1n2d + num1d2n } else { num1n2d - num1d2n } + return Fraction{ + n: n + d: f1.d * f2.d + is_reduced: true + } + } + // Here d1 > 1. + f1den := f1.d / d1 + f2den := f2.d / d1 + term1 := f1.n * f2den + term2 := f2.n * f1den + t := if addition { term1 + term2 } else { term1 - term2 } + d2 := math.gcd(t, d1) + return Fraction{ + n: t / d2 + d: f1den * (f2.d / d2) + is_reduced: true + } +} + +// Fraction add using operator overloading +pub fn (f1 Fraction) + (f2 Fraction) Fraction { + return general_addition_result(f1.reduce(), f2.reduce(), true) +} + +// Fraction subtract using operator overloading +pub fn (f1 Fraction) - (f2 Fraction) Fraction { + return general_addition_result(f1.reduce(), f2.reduce(), false) +} + +// Returns a correctly reduced result for both multiplication and division +// NOTE: requires reduced inputs +fn general_multiplication_result(f1 Fraction, f2 Fraction, multiplication bool) Fraction { + // * Theorem: If f1 and f2 are reduced i.e. gcd(f1.n, f1.d) == 1 and gcd(f2.n, f2.d) == 1, + // then gcd(f1.n * f2.n, f1.d * f2.d) == gcd(f1.n, f2.d) * gcd(f1.d, f2.n) + // * Knuth poses this an exercise for 4.5.1. - Exercise 2 + // * Also, note that: + // The terms are flipped for multiplication and division, so the gcds must be calculated carefully + // We do multiple divisions in order to prevent any possible overflows. + // * One more thing: + // if d = gcd(a, b) for example, then d divides both a and b + if multiplication { + d1 := math.gcd(f1.n, f2.d) + d2 := math.gcd(f1.d, f2.n) + return Fraction{ + n: (f1.n / d1) * (f2.n / d2) + d: (f2.d / d1) * (f1.d / d2) + is_reduced: true + } + } else { + d1 := math.gcd(f1.n, f2.n) + d2 := math.gcd(f1.d, f2.d) + return Fraction{ + n: (f1.n / d1) * (f2.d / d2) + d: (f2.n / d1) * (f1.d / d2) + is_reduced: true + } + } +} + +// Fraction multiply using operator overloading +pub fn (f1 Fraction) * (f2 Fraction) Fraction { + return general_multiplication_result(f1.reduce(), f2.reduce(), true) +} + +// Fraction divide using operator overloading +pub fn (f1 Fraction) / (f2 Fraction) Fraction { + if f2.n == 0 { + panic('Cannot divide by zero') + } + // If the second fraction is negative, it will + // mess up the sign. We need positive denominator + if f2.n < 0 { + return f1.negate() / f2.negate() + } + return general_multiplication_result(f1.reduce(), f2.reduce(), false) +} + +// Fraction negate method +pub fn (f Fraction) negate() Fraction { + return Fraction{ + n: -f.n + d: f.d + is_reduced: f.is_reduced + } +} + +// Fraction reciprocal method +pub fn (f Fraction) reciprocal() Fraction { + if f.n == 0 { + panic('Denominator cannot be zero') + } + return Fraction{ + n: f.d + d: f.n + is_reduced: f.is_reduced + } +} + +// Fraction method which reduces the fraction +pub fn (f Fraction) reduce() Fraction { + if f.is_reduced { + return f + } + cf := math.gcd(f.n, f.d) + return Fraction{ + n: f.n / cf + d: f.d / cf + is_reduced: true + } +} + +// f64 converts the Fraction to 64-bit floating point +pub fn (f Fraction) f64() f64 { + return f64(f.n) / f64(f.d) +} + +// +// + ------------------+ +// | Utility functions.| +// + ------------------+ +// +// Returns the absolute value of an i64 +fn abs(num i64) i64 { + if num < 0 { + return -num + } else { + return num + } +} + +// cmp_i64s compares the two arguments, returns 0 when equal, 1 when the first is bigger, -1 otherwise +fn cmp_i64s(a i64, b i64) int { + if a == b { + return 0 + } else if a > b { + return 1 + } else { + return -1 + } +} + +// cmp_f64s compares the two arguments, returns 0 when equal, 1 when the first is bigger, -1 otherwise +fn cmp_f64s(a f64, b f64) int { + // V uses epsilon comparison internally + if a == b { + return 0 + } else if a > b { + return 1 + } else { + return -1 + } +} + +// Two integers are safe to multiply when their bit lengths +// sum up to less than 64 (conservative estimate). +fn safe_to_multiply(a i64, b i64) bool { + return (bits.len_64(u64(abs(a))) + bits.len_64(u64(abs(b)))) < 64 +} + +// cmp compares the two arguments, returns 0 when equal, 1 when the first is bigger, -1 otherwise +fn cmp(f1 Fraction, f2 Fraction) int { + if safe_to_multiply(f1.n, f2.d) && safe_to_multiply(f2.n, f1.d) { + return cmp_i64s(f1.n * f2.d, f2.n * f1.d) + } else { + return cmp_f64s(f1.f64(), f2.f64()) + } +} + +// +-----------------------------+ +// | Public comparison functions | +// +-----------------------------+ +// equals returns true if both the Fractions are equal +pub fn (f1 Fraction) equals(f2 Fraction) bool { + return cmp(f1, f2) == 0 +} + +// ge returns true if f1 >= f2 +pub fn (f1 Fraction) ge(f2 Fraction) bool { + return cmp(f1, f2) >= 0 +} + +// gt returns true if f1 > f2 +pub fn (f1 Fraction) gt(f2 Fraction) bool { + return cmp(f1, f2) > 0 +} + +// le returns true if f1 <= f2 +pub fn (f1 Fraction) le(f2 Fraction) bool { + return cmp(f1, f2) <= 0 +} + +// lt returns true if f1 < f2 +pub fn (f1 Fraction) lt(f2 Fraction) bool { + return cmp(f1, f2) < 0 +} diff --git a/v_windows/v/old/vlib/math/fractions/fraction_test.v b/v_windows/v/old/vlib/math/fractions/fraction_test.v new file mode 100644 index 0000000..4928b7c --- /dev/null +++ b/v_windows/v/old/vlib/math/fractions/fraction_test.v @@ -0,0 +1,269 @@ +// Copyright (c) 2019-2021 Alexander Medvednikov. All rights reserved. +// Use of this source code is governed by an MIT license +// that can be found in the LICENSE file. +import math.fractions + +// (Old) results are verified using https://www.calculatorsoup.com/calculators/math/fractions.php +// Newer ones are contrived for corner cases or prepared by hand. +fn test_4_by_8_f64_and_str() { + f := fractions.fraction(4, 8) + assert f.f64() == 0.5 + assert f.str() == '4/8' +} + +fn test_10_by_5_f64_and_str() { + f := fractions.fraction(10, 5) + assert f.f64() == 2.0 + assert f.str() == '10/5' +} + +fn test_9_by_3_f64_and_str() { + f := fractions.fraction(9, 3) + assert f.f64() == 3.0 + assert f.str() == '9/3' +} + +fn test_4_by_minus_5_f64_and_str() { + f := fractions.fraction(4, -5) + assert f.f64() == -0.8 + assert f.str() == '-4/5' +} + +fn test_minus_7_by_minus_92_str() { + f := fractions.fraction(-7, -5) + assert f.str() == '7/5' +} + +fn test_4_by_8_plus_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(5, 10) + sum := f1 + f2 + assert sum.f64() == 1.0 + assert sum.str() == '1/1' + assert sum.equals(fractions.fraction(1, 1)) +} + +fn test_5_by_5_plus_8_by_8() { + f1 := fractions.fraction(5, 5) + f2 := fractions.fraction(8, 8) + sum := f1 + f2 + assert sum.f64() == 2.0 + assert sum.str() == '2/1' + assert sum.equals(fractions.fraction(2, 1)) +} + +fn test_9_by_3_plus_1_by_3() { + f1 := fractions.fraction(9, 3) + f2 := fractions.fraction(1, 3) + sum := f1 + f2 + assert sum.str() == '10/3' + assert sum.equals(fractions.fraction(10, 3)) +} + +fn test_3_by_7_plus_1_by_4() { + f1 := fractions.fraction(3, 7) + f2 := fractions.fraction(1, 4) + sum := f1 + f2 + assert sum.str() == '19/28' + assert sum.equals(fractions.fraction(19, 28)) +} + +fn test_36529_by_12409100000_plus_418754901_by_9174901000() { + f1 := fractions.fraction(i64(36529), i64(12409100000)) + f2 := fractions.fraction(i64(418754901), i64(9174901000)) + sum := f1 + f2 + assert sum.str() == '5196706591957729/113852263999100000' +} + +fn test_4_by_8_plus_minus_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(-5, 10) + diff := f2 + f1 + assert diff.f64() == 0 + assert diff.str() == '0/1' +} + +fn test_4_by_8_minus_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(5, 10) + diff := f2 - f1 + assert diff.f64() == 0 + assert diff.str() == '0/1' +} + +fn test_5_by_5_minus_8_by_8() { + f1 := fractions.fraction(5, 5) + f2 := fractions.fraction(8, 8) + diff := f2 - f1 + assert diff.f64() == 0 + assert diff.str() == '0/1' +} + +fn test_9_by_3_minus_1_by_3() { + f1 := fractions.fraction(9, 3) + f2 := fractions.fraction(1, 3) + diff := f1 - f2 + assert diff.str() == '8/3' +} + +fn test_3_by_7_minus_1_by_4() { + f1 := fractions.fraction(3, 7) + f2 := fractions.fraction(1, 4) + diff := f1 - f2 + assert diff.str() == '5/28' +} + +fn test_36529_by_12409100000_minus_418754901_by_9174901000() { + f1 := fractions.fraction(i64(36529), i64(12409100000)) + f2 := fractions.fraction(i64(418754901), i64(9174901000)) + sum := f1 - f2 + assert sum.str() == '-5196036292040471/113852263999100000' +} + +fn test_4_by_8_times_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(5, 10) + product := f1 * f2 + assert product.f64() == 0.25 + assert product.str() == '1/4' +} + +fn test_5_by_5_times_8_by_8() { + f1 := fractions.fraction(5, 5) + f2 := fractions.fraction(8, 8) + product := f1 * f2 + assert product.f64() == 1.0 + assert product.str() == '1/1' +} + +fn test_9_by_3_times_1_by_3() { + f1 := fractions.fraction(9, 3) + f2 := fractions.fraction(1, 3) + product := f1 * f2 + assert product.f64() == 1.0 + assert product.str() == '1/1' +} + +fn test_3_by_7_times_1_by_4() { + f1 := fractions.fraction(3, 7) + f2 := fractions.fraction(1, 4) + product := f2 * f1 + assert product.f64() == (3.0 / 28.0) + assert product.str() == '3/28' +} + +fn test_4_by_8_over_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(5, 10) + q := f1 / f2 + assert q.f64() == 1.0 + assert q.str() == '1/1' +} + +fn test_5_by_5_over_8_by_8() { + f1 := fractions.fraction(5, 5) + f2 := fractions.fraction(8, 8) + q := f1 / f2 + assert q.f64() == 1.0 + assert q.str() == '1/1' +} + +fn test_9_by_3_over_1_by_3() { + f1 := fractions.fraction(9, 3) + f2 := fractions.fraction(1, 3) + q := f1 / f2 + assert q.f64() == 9.0 + assert q.str() == '9/1' +} + +fn test_3_by_7_over_1_by_4() { + f1 := fractions.fraction(3, 7) + f2 := fractions.fraction(1, 4) + q := f1 / f2 + assert q.str() == '12/7' +} + +fn test_reciprocal_4_by_8() { + f := fractions.fraction(4, 8) + assert f.reciprocal().str() == '8/4' +} + +fn test_reciprocal_5_by_10() { + f := fractions.fraction(5, 10) + assert f.reciprocal().str() == '10/5' +} + +fn test_reciprocal_5_by_5() { + f := fractions.fraction(5, 5) + assert f.reciprocal().str() == '5/5' +} + +fn test_reciprocal_8_by_8() { + f := fractions.fraction(8, 8) + assert f.reciprocal().str() == '8/8' +} + +fn test_reciprocal_9_by_3() { + f := fractions.fraction(9, 3) + assert f.reciprocal().str() == '3/9' +} + +fn test_reciprocal_1_by_3() { + f := fractions.fraction(1, 3) + assert f.reciprocal().str() == '3/1' +} + +fn test_reciprocal_7_by_3() { + f := fractions.fraction(7, 3) + assert f.reciprocal().str() == '3/7' +} + +fn test_reciprocal_1_by_4() { + f := fractions.fraction(1, 4) + assert f.reciprocal().str() == '4/1' +} + +fn test_4_by_8_equals_5_by_10() { + f1 := fractions.fraction(4, 8) + f2 := fractions.fraction(5, 10) + assert f1.equals(f2) +} + +fn test_1_by_2_does_not_equal_3_by_4() { + f1 := fractions.fraction(1, 2) + f2 := fractions.fraction(3, 4) + assert !f1.equals(f2) +} + +fn test_reduce_3_by_9() { + f := fractions.fraction(3, 9) + assert f.reduce().equals(fractions.fraction(1, 3)) +} + +fn test_1_by_3_less_than_2_by_4() { + f1 := fractions.fraction(1, 3) + f2 := fractions.fraction(2, 4) + assert f1.lt(f2) + assert f1.le(f2) +} + +fn test_2_by_3_greater_than_2_by_4() { + f1 := fractions.fraction(2, 3) + f2 := fractions.fraction(2, 4) + assert f1.gt(f2) + assert f1.ge(f2) +} + +fn test_5_by_7_not_less_than_2_by_4() { + f1 := fractions.fraction(5, 7) + f2 := fractions.fraction(2, 4) + assert !f1.lt(f2) + assert !f1.le(f2) +} + +fn test_49_by_75_not_greater_than_2_by_3() { + f1 := fractions.fraction(49, 75) + f2 := fractions.fraction(2, 3) + assert !f1.gt(f2) + assert !f1.ge(f2) +} diff --git a/v_windows/v/old/vlib/math/math.c.v b/v_windows/v/old/vlib/math/math.c.v new file mode 100644 index 0000000..689e648 --- /dev/null +++ b/v_windows/v/old/vlib/math/math.c.v @@ -0,0 +1,296 @@ +// Copyright (c) 2019-2021 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 + +#include <math.h> + +$if windows { + $if tinyc { + #flag @VEXEROOT/thirdparty/tcc/lib/openlibm.o + } +} $else { + #flag -lm +} + +fn C.acos(x f64) f64 + +fn C.asin(x f64) f64 + +fn C.atan(x f64) f64 + +fn C.atan2(y f64, x f64) f64 + +fn C.cbrt(x f64) f64 + +fn C.ceil(x f64) f64 + +fn C.cos(x f64) f64 + +fn C.cosf(x f32) f32 + +fn C.cosh(x f64) f64 + +fn C.erf(x f64) f64 + +fn C.erfc(x f64) f64 + +fn C.exp(x f64) f64 + +fn C.exp2(x f64) f64 + +fn C.fabs(x f64) f64 + +fn C.floor(x f64) f64 + +fn C.fmod(x f64, y f64) f64 + +fn C.hypot(x f64, y f64) f64 + +fn C.log(x f64) f64 + +fn C.log2(x f64) f64 + +fn C.log10(x f64) f64 + +fn C.lgamma(x f64) f64 + +fn C.pow(x f64, y f64) f64 + +fn C.powf(x f32, y f32) f32 + +fn C.round(x f64) f64 + +fn C.sin(x f64) f64 + +fn C.sinf(x f32) f32 + +fn C.sinh(x f64) f64 + +fn C.sqrt(x f64) f64 + +fn C.sqrtf(x f32) f32 + +fn C.tgamma(x f64) f64 + +fn C.tan(x f64) f64 + +fn C.tanf(x f32) f32 + +fn C.tanh(x f64) f64 + +fn C.trunc(x f64) f64 + +// NOTE +// When adding a new function, please make sure it's in the right place. +// All functions are sorted alphabetically. +// Returns the absolute value. +[inline] +pub fn abs(a f64) f64 { + return C.fabs(a) +} + +// acos calculates inverse cosine (arccosine). +[inline] +pub fn acos(a f64) f64 { + return C.acos(a) +} + +// asin calculates inverse sine (arcsine). +[inline] +pub fn asin(a f64) f64 { + return C.asin(a) +} + +// atan calculates inverse tangent (arctangent). +[inline] +pub fn atan(a f64) f64 { + return C.atan(a) +} + +// atan2 calculates inverse tangent with two arguments, returns the angle between the X axis and the point. +[inline] +pub fn atan2(a f64, b f64) f64 { + return C.atan2(a, b) +} + +// cbrt calculates cubic root. +[inline] +pub fn cbrt(a f64) f64 { + return C.cbrt(a) +} + +// ceil returns the nearest f64 greater or equal to the provided value. +[inline] +pub fn ceil(a f64) f64 { + return C.ceil(a) +} + +// cos calculates cosine. +[inline] +pub fn cos(a f64) f64 { + return C.cos(a) +} + +// cosf calculates cosine. (float32) +[inline] +pub fn cosf(a f32) f32 { + return C.cosf(a) +} + +// cosh calculates hyperbolic cosine. +[inline] +pub fn cosh(a f64) f64 { + return C.cosh(a) +} + +// exp calculates exponent of the number (math.pow(math.E, a)). +[inline] +pub fn exp(a f64) f64 { + return C.exp(a) +} + +// erf computes the error function value +[inline] +pub fn erf(a f64) f64 { + return C.erf(a) +} + +// erfc computes the complementary error function value +[inline] +pub fn erfc(a f64) f64 { + return C.erfc(a) +} + +// exp2 returns the base-2 exponential function of a (math.pow(2, a)). +[inline] +pub fn exp2(a f64) f64 { + return C.exp2(a) +} + +// floor returns the nearest f64 lower or equal of the provided value. +[inline] +pub fn floor(a f64) f64 { + return C.floor(a) +} + +// fmod returns the floating-point remainder of number / denom (rounded towards zero): +[inline] +pub fn fmod(a f64, b f64) f64 { + return C.fmod(a, b) +} + +// gamma computes the gamma function value +[inline] +pub fn gamma(a f64) f64 { + return C.tgamma(a) +} + +// Returns hypotenuse of a right triangle. +[inline] +pub fn hypot(a f64, b f64) f64 { + return C.hypot(a, b) +} + +// log calculates natural (base-e) logarithm of the provided value. +[inline] +pub fn log(a f64) f64 { + return C.log(a) +} + +// log2 calculates base-2 logarithm of the provided value. +[inline] +pub fn log2(a f64) f64 { + return C.log2(a) +} + +// log10 calculates the common (base-10) logarithm of the provided value. +[inline] +pub fn log10(a f64) f64 { + return C.log10(a) +} + +// log_gamma computes the log-gamma function value +[inline] +pub fn log_gamma(a f64) f64 { + return C.lgamma(a) +} + +// log_n calculates base-N logarithm of the provided value. +[inline] +pub fn log_n(a f64, b f64) f64 { + return C.log(a) / C.log(b) +} + +// pow returns base raised to the provided power. +[inline] +pub fn pow(a f64, b f64) f64 { + return C.pow(a, b) +} + +// powf returns base raised to the provided power. (float32) +[inline] +pub fn powf(a f32, b f32) f32 { + return C.powf(a, b) +} + +// round returns the integer nearest to the provided value. +[inline] +pub fn round(f f64) f64 { + return C.round(f) +} + +// sin calculates sine. +[inline] +pub fn sin(a f64) f64 { + return C.sin(a) +} + +// sinf calculates sine. (float32) +[inline] +pub fn sinf(a f32) f32 { + return C.sinf(a) +} + +// sinh calculates hyperbolic sine. +[inline] +pub fn sinh(a f64) f64 { + return C.sinh(a) +} + +// sqrt calculates square-root of the provided value. +[inline] +pub fn sqrt(a f64) f64 { + return C.sqrt(a) +} + +// sqrtf calculates square-root of the provided value. (float32) +[inline] +pub fn sqrtf(a f32) f32 { + return C.sqrtf(a) +} + +// tan calculates tangent. +[inline] +pub fn tan(a f64) f64 { + return C.tan(a) +} + +// tanf calculates tangent. (float32) +[inline] +pub fn tanf(a f32) f32 { + return C.tanf(a) +} + +// tanh calculates hyperbolic tangent. +[inline] +pub fn tanh(a f64) f64 { + return C.tanh(a) +} + +// trunc rounds a toward zero, returning the nearest integral value that is not +// larger in magnitude than a. +[inline] +pub fn trunc(a f64) f64 { + return C.trunc(a) +} diff --git a/v_windows/v/old/vlib/math/math.js.v b/v_windows/v/old/vlib/math/math.js.v new file mode 100644 index 0000000..437e8ab --- /dev/null +++ b/v_windows/v/old/vlib/math/math.js.v @@ -0,0 +1,260 @@ +// Copyright (c) 2019-2021 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 + +// TODO : The commented out functions need either a native V implementation, a +// JS specific implementation, or use some other JS math library, such as +// https://github.com/josdejong/mathjs + +// Replaces C.fabs +fn JS.Math.abs(x f64) f64 + +fn JS.Math.acos(x f64) f64 +fn JS.Math.asin(x f64) f64 +fn JS.Math.atan(x f64) f64 +fn JS.Math.atan2(y f64, x f64) f64 +fn JS.Math.cbrt(x f64) f64 +fn JS.Math.ceil(x f64) f64 +fn JS.Math.cos(x f64) f64 +fn JS.Math.cosh(x f64) f64 + +// fn JS.Math.erf(x f64) f64 // Not in standard JS Math object +// fn JS.Math.erfc(x f64) f64 // Not in standard JS Math object +fn JS.Math.exp(x f64) f64 + +// fn JS.Math.exp2(x f64) f64 // Not in standard JS Math object +fn JS.Math.floor(x f64) f64 + +// fn JS.Math.fmod(x f64, y f64) f64 // Not in standard JS Math object +// fn JS.Math.hypot(x f64, y f64) f64 // Not in standard JS Math object +fn JS.Math.log(x f64) f64 + +// fn JS.Math.log2(x f64) f64 // Not in standard JS Math object +// fn JS.Math.log10(x f64) f64 // Not in standard JS Math object +// fn JS.Math.lgamma(x f64) f64 // Not in standard JS Math object +fn JS.Math.pow(x f64, y f64) f64 +fn JS.Math.round(x f64) f64 +fn JS.Math.sin(x f64) f64 +fn JS.Math.sinh(x f64) f64 +fn JS.Math.sqrt(x f64) f64 + +// fn JS.Math.tgamma(x f64) f64 // Not in standard JS Math object +fn JS.Math.tan(x f64) f64 +fn JS.Math.tanh(x f64) f64 +fn JS.Math.trunc(x f64) f64 + +// NOTE +// When adding a new function, please make sure it's in the right place. +// All functions are sorted alphabetically. + +// Returns the absolute value. +[inline] +pub fn abs(a f64) f64 { + return JS.Math.abs(a) +} + +// acos calculates inverse cosine (arccosine). +[inline] +pub fn acos(a f64) f64 { + return JS.Math.acos(a) +} + +// asin calculates inverse sine (arcsine). +[inline] +pub fn asin(a f64) f64 { + return JS.Math.asin(a) +} + +// atan calculates inverse tangent (arctangent). +[inline] +pub fn atan(a f64) f64 { + return JS.Math.atan(a) +} + +// atan2 calculates inverse tangent with two arguments, returns the angle between the X axis and the point. +[inline] +pub fn atan2(a f64, b f64) f64 { + return JS.Math.atan2(a, b) +} + +// cbrt calculates cubic root. +[inline] +pub fn cbrt(a f64) f64 { + return JS.Math.cbrt(a) +} + +// ceil returns the nearest f64 greater or equal to the provided value. +[inline] +pub fn ceil(a f64) f64 { + return JS.Math.ceil(a) +} + +// cos calculates cosine. +[inline] +pub fn cos(a f64) f64 { + return JS.Math.cos(a) +} + +// cosf calculates cosine. (float32). This doesn't exist in JS +[inline] +pub fn cosf(a f32) f32 { + return f32(JS.Math.cos(a)) +} + +// cosh calculates hyperbolic cosine. +[inline] +pub fn cosh(a f64) f64 { + return JS.Math.cosh(a) +} + +// exp calculates exponent of the number (math.pow(math.E, a)). +[inline] +pub fn exp(a f64) f64 { + return JS.Math.exp(a) +} + +// erf computes the error function value +[inline] +pub fn erf(a f64) f64 { + return JS.Math.erf(a) +} + +// erfc computes the complementary error function value +[inline] +pub fn erfc(a f64) f64 { + return JS.Math.erfc(a) +} + +// exp2 returns the base-2 exponential function of a (math.pow(2, a)). +[inline] +pub fn exp2(a f64) f64 { + return JS.Math.exp2(a) +} + +// floor returns the nearest f64 lower or equal of the provided value. +[inline] +pub fn floor(a f64) f64 { + return JS.Math.floor(a) +} + +// fmod returns the floating-point remainder of number / denom (rounded towards zero): +[inline] +pub fn fmod(a f64, b f64) f64 { + return JS.Math.fmod(a, b) +} + +// gamma computes the gamma function value +[inline] +pub fn gamma(a f64) f64 { + return JS.Math.tgamma(a) +} + +// Returns hypotenuse of a right triangle. +[inline] +pub fn hypot(a f64, b f64) f64 { + return JS.Math.hypot(a, b) +} + +// log calculates natural (base-e) logarithm of the provided value. +[inline] +pub fn log(a f64) f64 { + return JS.Math.log(a) +} + +// log2 calculates base-2 logarithm of the provided value. +[inline] +pub fn log2(a f64) f64 { + return JS.Math.log2(a) +} + +// log10 calculates the common (base-10) logarithm of the provided value. +[inline] +pub fn log10(a f64) f64 { + return JS.Math.log10(a) +} + +// log_gamma computes the log-gamma function value +[inline] +pub fn log_gamma(a f64) f64 { + return JS.Math.lgamma(a) +} + +// log_n calculates base-N logarithm of the provided value. +[inline] +pub fn log_n(a f64, b f64) f64 { + return JS.Math.log(a) / JS.Math.log(b) +} + +// pow returns base raised to the provided power. +[inline] +pub fn pow(a f64, b f64) f64 { + return JS.Math.pow(a, b) +} + +// powf returns base raised to the provided power. (float32) +[inline] +pub fn powf(a f32, b f32) f32 { + return f32(JS.Math.pow(a, b)) +} + +// round returns the integer nearest to the provided value. +[inline] +pub fn round(f f64) f64 { + return JS.Math.round(f) +} + +// sin calculates sine. +[inline] +pub fn sin(a f64) f64 { + return JS.Math.sin(a) +} + +// sinf calculates sine. (float32) +[inline] +pub fn sinf(a f32) f32 { + return f32(JS.Math.sin(a)) +} + +// sinh calculates hyperbolic sine. +[inline] +pub fn sinh(a f64) f64 { + return JS.Math.sinh(a) +} + +// sqrt calculates square-root of the provided value. +[inline] +pub fn sqrt(a f64) f64 { + return JS.Math.sqrt(a) +} + +// sqrtf calculates square-root of the provided value. (float32) +[inline] +pub fn sqrtf(a f32) f32 { + return f32(JS.Math.sqrt(a)) +} + +// tan calculates tangent. +[inline] +pub fn tan(a f64) f64 { + return JS.Math.tan(a) +} + +// tanf calculates tangent. (float32) +[inline] +pub fn tanf(a f32) f32 { + return f32(JS.Math.tan(a)) +} + +// tanh calculates hyperbolic tangent. +[inline] +pub fn tanh(a f64) f64 { + return JS.Math.tanh(a) +} + +// trunc rounds a toward zero, returning the nearest integral value that is not +// larger in magnitude than a. +[inline] +pub fn trunc(a f64) f64 { + return JS.Math.trunc(a) +} diff --git a/v_windows/v/old/vlib/math/math.v b/v_windows/v/old/vlib/math/math.v new file mode 100644 index 0000000..f102f39 --- /dev/null +++ b/v_windows/v/old/vlib/math/math.v @@ -0,0 +1,148 @@ +// Copyright (c) 2019-2021 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 + +// aprox_sin returns an approximation of sin(a) made using lolremez +pub fn aprox_sin(a f64) f64 { + a0 := 1.91059300966915117e-31 + a1 := 1.00086760103908896 + a2 := -1.21276126894734565e-2 + a3 := -1.38078780785773762e-1 + a4 := -2.67353392911981221e-2 + a5 := 2.08026600266304389e-2 + a6 := -3.03996055049204407e-3 + a7 := 1.38235642404333740e-4 + return a0 + a * (a1 + a * (a2 + a * (a3 + a * (a4 + a * (a5 + a * (a6 + a * a7)))))) +} + +// aprox_cos returns an approximation of sin(a) made using lolremez +pub fn aprox_cos(a f64) f64 { + a0 := 9.9995999154986614e-1 + a1 := 1.2548995793001028e-3 + a2 := -5.0648546280678015e-1 + a3 := 1.2942246466519995e-2 + a4 := 2.8668384702547972e-2 + a5 := 7.3726485210586547e-3 + a6 := -3.8510875386947414e-3 + a7 := 4.7196604604366623e-4 + a8 := -1.8776444013090451e-5 + return a0 + a * (a1 + a * (a2 + a * (a3 + a * (a4 + a * (a5 + a * (a6 + a * (a7 + a * a8))))))) +} + +// copysign returns a value with the magnitude of x and the sign of y +[inline] +pub fn copysign(x f64, y f64) f64 { + return f64_from_bits((f64_bits(x) & ~sign_mask) | (f64_bits(y) & sign_mask)) +} + +// degrees convert from degrees to radians. +[inline] +pub fn degrees(radians f64) f64 { + return radians * (180.0 / pi) +} + +// digits returns an array of the digits of n in the given base. +pub fn digits(_n int, base int) []int { + if base < 2 { + panic('digits: Cannot find digits of n with base $base') + } + mut n := _n + mut sign := 1 + if n < 0 { + sign = -1 + n = -n + } + mut res := []int{} + for n != 0 { + res << (n % base) * sign + n /= base + } + return res +} + +[inline] +pub fn fabs(x f64) f64 { + if x < 0.0 { + return -x + } + return x +} + +// gcd calculates greatest common (positive) divisor (or zero if a and b are both zero). +pub fn gcd(a_ i64, b_ i64) i64 { + mut a := a_ + mut b := b_ + 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. +pub fn lcm(a i64, b i64) i64 { + if a == 0 { + return a + } + res := a * (b / gcd(b, a)) + if res < 0 { + return -res + } + return res +} + +// max returns the maximum value of the two provided. +[inline] +pub fn max(a f64, b f64) f64 { + if a > b { + return a + } + return b +} + +// min returns the minimum value of the two provided. +[inline] +pub fn min(a f64, b f64) f64 { + if a < b { + return a + } + return b +} + +// sign returns the corresponding sign -1.0, 1.0 of the provided number. +// if n is not a number, its sign is nan too. +[inline] +pub fn sign(n f64) f64 { + if is_nan(n) { + return nan() + } + return copysign(1.0, n) +} + +// signi returns the corresponding sign -1.0, 1.0 of the provided number. +[inline] +pub fn signi(n f64) int { + return int(copysign(1.0, n)) +} + +// radians convert from radians to degrees. +[inline] +pub fn radians(degrees f64) f64 { + return degrees * (pi / 180.0) +} + +// signbit returns a value with the boolean representation of the sign for x +[inline] +pub fn signbit(x f64) bool { + return f64_bits(x) & sign_mask != 0 +} diff --git a/v_windows/v/old/vlib/math/math_test.v b/v_windows/v/old/vlib/math/math_test.v new file mode 100644 index 0000000..c0cd7d8 --- /dev/null +++ b/v_windows/v/old/vlib/math/math_test.v @@ -0,0 +1,806 @@ +module math + +struct Fi { + f f64 + i int +} + +const ( + vf_ = [f64(4.9790119248836735e+00), 7.7388724745781045e+00, -2.7688005719200159e-01, + -5.0106036182710749e+00, 9.6362937071984173e+00, 2.9263772392439646e+00, + 5.2290834314593066e+00, 2.7279399104360102e+00, 1.8253080916808550e+00, + -8.6859247685756013e+00, + ] + // The expected results below were computed by the high precision calculators + // at https://keisan.casio.com/. More exact input values (array vf_[], above) + // were obtained by printing them with "%.26f". The answers were calculated + // to 26 digits (by using the "Digit number" drop-down control of each + // calculator). + acos_ = [f64(1.0496193546107222142571536e+00), 6.8584012813664425171660692e-01, + 1.5984878714577160325521819e+00, 2.0956199361475859327461799e+00, + 2.7053008467824138592616927e-01, 1.2738121680361776018155625e+00, + 1.0205369421140629186287407e+00, 1.2945003481781246062157835e+00, + 1.3872364345374451433846657e+00, 2.6231510803970463967294145e+00] + asin_ = [f64(5.2117697218417440497416805e-01), 8.8495619865825236751471477e-01, + -2.769154466281941332086016e-02, -5.2482360935268931351485822e-01, + 1.3002662421166552333051524e+00, 2.9698415875871901741575922e-01, + 5.5025938468083370060258102e-01, 2.7629597861677201301553823e-01, + 1.83559892257451475846656e-01, -1.0523547536021497774980928e+00] + atan_ = [f64(1.372590262129621651920085e+00), 1.442290609645298083020664e+00, + -2.7011324359471758245192595e-01, -1.3738077684543379452781531e+00, + 1.4673921193587666049154681e+00, 1.2415173565870168649117764e+00, + 1.3818396865615168979966498e+00, 1.2194305844639670701091426e+00, + 1.0696031952318783760193244e+00, -1.4561721938838084990898679e+00] + atan2_ = [f64(1.1088291730037004444527075e+00), 9.1218183188715804018797795e-01, + 1.5984772603216203736068915e+00, 2.0352918654092086637227327e+00, + 8.0391819139044720267356014e-01, 1.2861075249894661588866752e+00, + 1.0889904479131695712182587e+00, 1.3044821793397925293797357e+00, + 1.3902530903455392306872261e+00, 2.2859857424479142655411058e+00] + ceil_ = [f64(5.0000000000000000e+00), 8.0000000000000000e+00, copysign(0, -1), + -5.0000000000000000e+00, 1.0000000000000000e+01, 3.0000000000000000e+00, + 6.0000000000000000e+00, 3.0000000000000000e+00, 2.0000000000000000e+00, + -8.0000000000000000e+00, + ] + cos_ = [f64(2.634752140995199110787593e-01), 1.148551260848219865642039e-01, + 9.6191297325640768154550453e-01, 2.938141150061714816890637e-01, + -9.777138189897924126294461e-01, -9.7693041344303219127199518e-01, + 4.940088096948647263961162e-01, -9.1565869021018925545016502e-01, + -2.517729313893103197176091e-01, -7.39241351595676573201918e-01] + // Results for 100000 * pi + vf_[i] + cos_large_ = [f64(2.634752141185559426744e-01), 1.14855126055543100712e-01, + 9.61912973266488928113e-01, 2.9381411499556122552e-01, -9.777138189880161924641e-01, + -9.76930413445147608049e-01, 4.940088097314976789841e-01, -9.15658690217517835002e-01, + -2.51772931436786954751e-01, -7.3924135157173099849e-01] + cosh_ = [f64(7.2668796942212842775517446e+01), 1.1479413465659254502011135e+03, + 1.0385767908766418550935495e+00, 7.5000957789658051428857788e+01, + 7.655246669605357888468613e+03, 9.3567491758321272072888257e+00, + 9.331351599270605471131735e+01, 7.6833430994624643209296404e+00, + 3.1829371625150718153881164e+00, 2.9595059261916188501640911e+03] + exp_ = [f64(1.4533071302642137507696589e+02), 2.2958822575694449002537581e+03, + 7.5814542574851666582042306e-01, 6.6668778421791005061482264e-03, + 1.5310493273896033740861206e+04, 1.8659907517999328638667732e+01, + 1.8662167355098714543942057e+02, 1.5301332413189378961665788e+01, + 6.2047063430646876349125085e+00, 1.6894712385826521111610438e-04] + exp2_ = [f64(3.1537839463286288034313104e+01), 2.1361549283756232296144849e+02, + 8.2537402562185562902577219e-01, 3.1021158628740294833424229e-02, + 7.9581744110252191462569661e+02, 7.6019905892596359262696423e+00, + 3.7506882048388096973183084e+01, 6.6250893439173561733216375e+00, + 3.5438267900243941544605339e+00, 2.4281533133513300984289196e-03] + fabs_ = [f64(4.9790119248836735e+00), 7.7388724745781045e+00, 2.7688005719200159e-01, + 5.0106036182710749e+00, 9.6362937071984173e+00, 2.9263772392439646e+00, + 5.2290834314593066e+00, 2.7279399104360102e+00, 1.8253080916808550e+00, + 8.6859247685756013e+00, + ] + floor_ = [f64(4.0000000000000000e+00), 7.0000000000000000e+00, -1.0000000000000000e+00, + -6.0000000000000000e+00, 9.0000000000000000e+00, 2.0000000000000000e+00, + 5.0000000000000000e+00, 2.0000000000000000e+00, 1.0000000000000000e+00, + -9.0000000000000000e+00, + ] + fmod_ = [f64(4.197615023265299782906368e-02), 2.261127525421895434476482e+00, + 3.231794108794261433104108e-02, 4.989396381728925078391512e+00, + 3.637062928015826201999516e-01, 1.220868282268106064236690e+00, + 4.770916568540693347699744e+00, 1.816180268691969246219742e+00, + 8.734595415957246977711748e-01, 1.314075231424398637614104e+00] + gamma_ = [f64(2.3254348370739963835386613898e+01), 2.991153837155317076427529816e+03, + -4.561154336726758060575129109e+00, 7.719403468842639065959210984e-01, + 1.6111876618855418534325755566e+05, 1.8706575145216421164173224946e+00, + 3.4082787447257502836734201635e+01, 1.579733951448952054898583387e+00, + 9.3834586598354592860187267089e-01, -2.093995902923148389186189429e-05] + log_ = [f64(1.605231462693062999102599e+00), 2.0462560018708770653153909e+00, + -1.2841708730962657801275038e+00, 1.6115563905281545116286206e+00, + 2.2655365644872016636317461e+00, 1.0737652208918379856272735e+00, + 1.6542360106073546632707956e+00, 1.0035467127723465801264487e+00, + 6.0174879014578057187016475e-01, 2.161703872847352815363655e+00] + logb_ = [f64(2.0000000000000000e+00), 2.0000000000000000e+00, -2.0000000000000000e+00, + 2.0000000000000000e+00, 3.0000000000000000e+00, 1.0000000000000000e+00, + 2.0000000000000000e+00, 1.0000000000000000e+00, 0.0000000000000000e+00, + 3.0000000000000000e+00, + ] + log10_ = [f64(6.9714316642508290997617083e-01), 8.886776901739320576279124e-01, + -5.5770832400658929815908236e-01, 6.998900476822994346229723e-01, + 9.8391002850684232013281033e-01, 4.6633031029295153334285302e-01, + 7.1842557117242328821552533e-01, 4.3583479968917773161304553e-01, + 2.6133617905227038228626834e-01, 9.3881606348649405716214241e-01] + log1p_ = [f64(4.8590257759797794104158205e-02), 7.4540265965225865330849141e-02, + -2.7726407903942672823234024e-03, -5.1404917651627649094953380e-02, + 9.1998280672258624681335010e-02, 2.8843762576593352865894824e-02, + 5.0969534581863707268992645e-02, 2.6913947602193238458458594e-02, + 1.8088493239630770262045333e-02, -9.0865245631588989681559268e-02] + log2_ = [f64(2.3158594707062190618898251e+00), 2.9521233862883917703341018e+00, + -1.8526669502700329984917062e+00, 2.3249844127278861543568029e+00, + 3.268478366538305087466309e+00, 1.5491157592596970278166492e+00, + 2.3865580889631732407886495e+00, 1.447811865817085365540347e+00, + 8.6813999540425116282815557e-01, 3.118679457227342224364709e+00] + modf_ = [[f64(4.0000000000000000e+00), 9.7901192488367350108546816e-01], + [f64(7.0000000000000000e+00), 7.3887247457810456552351752e-01], + [f64(-0.0), -2.7688005719200159404635997e-01], + [f64(-5.0000000000000000e+00), + -1.060361827107492160848778e-02, + ], + [f64(9.0000000000000000e+00), 6.3629370719841737980004837e-01], + [f64(2.0000000000000000e+00), 9.2637723924396464525443662e-01], + [f64(5.0000000000000000e+00), 2.2908343145930665230025625e-01], + [f64(2.0000000000000000e+00), 7.2793991043601025126008608e-01], + [f64(1.0000000000000000e+00), 8.2530809168085506044576505e-01], + [f64(-8.0000000000000000e+00), -6.8592476857560136238589621e-01], + ] + nextafter32_ = [4.979012489318848e+00, 7.738873004913330e+00, -2.768800258636475e-01, + -5.010602951049805e+00, 9.636294364929199e+00, 2.926377534866333e+00, 5.229084014892578e+00, + 2.727940082550049e+00, 1.825308203697205e+00, -8.685923576354980e+00] + nextafter64_ = [f64(4.97901192488367438926388786e+00), 7.73887247457810545370193722e+00, + -2.7688005719200153853520874e-01, -5.01060361827107403343006808e+00, + 9.63629370719841915615688777e+00, 2.92637723924396508934364647e+00, + 5.22908343145930754047867595e+00, 2.72793991043601069534929593e+00, + 1.82530809168085528249036997e+00, -8.68592476857559958602905681e+00] + pow_ = [f64(9.5282232631648411840742957e+04), 5.4811599352999901232411871e+07, + 5.2859121715894396531132279e-01, 9.7587991957286474464259698e-06, + 4.328064329346044846740467e+09, 8.4406761805034547437659092e+02, + 1.6946633276191194947742146e+05, 5.3449040147551939075312879e+02, + 6.688182138451414936380374e+01, 2.0609869004248742886827439e-09] + remainder_ = [f64(4.197615023265299782906368e-02), 2.261127525421895434476482e+00, + 3.231794108794261433104108e-02, -2.120723654214984321697556e-02, + 3.637062928015826201999516e-01, 1.220868282268106064236690e+00, + -4.581668629186133046005125e-01, -9.117596417440410050403443e-01, + 8.734595415957246977711748e-01, 1.314075231424398637614104e+00] + round_ = [f64(5), 8, -0.0, -5, 10, 3, 5, 3, 2, -9] + signbit_ = [false, false, true, true, false, false, false, false, false, true] + sin_ = [f64(-9.6466616586009283766724726e-01), 9.9338225271646545763467022e-01, + -2.7335587039794393342449301e-01, 9.5586257685042792878173752e-01, + -2.099421066779969164496634e-01, 2.135578780799860532750616e-01, + -8.694568971167362743327708e-01, 4.019566681155577786649878e-01, + 9.6778633541687993721617774e-01, -6.734405869050344734943028e-01] + // Results for 100000 * pi + vf_[i] + sin_large_ = [f64(-9.646661658548936063912e-01), 9.933822527198506903752e-01, + -2.7335587036246899796e-01, 9.55862576853689321268e-01, -2.099421066862688873691e-01, + 2.13557878070308981163e-01, -8.694568970959221300497e-01, 4.01956668098863248917e-01, + 9.67786335404528727927e-01, -6.7344058693131973066e-01] + sinh_ = [f64(7.2661916084208532301448439e+01), 1.1479409110035194500526446e+03, + -2.8043136512812518927312641e-01, -7.499429091181587232835164e+01, + 7.6552466042906758523925934e+03, 9.3031583421672014313789064e+00, + 9.330815755828109072810322e+01, 7.6179893137269146407361477e+00, + 3.021769180549615819524392e+00, -2.95950575724449499189888e+03] + sqrt_ = [f64(2.2313699659365484748756904e+00), 2.7818829009464263511285458e+00, + 5.2619393496314796848143251e-01, 2.2384377628763938724244104e+00, + 3.1042380236055381099288487e+00, 1.7106657298385224403917771e+00, + 2.286718922705479046148059e+00, 1.6516476350711159636222979e+00, + 1.3510396336454586262419247e+00, 2.9471892997524949215723329e+00] + tan_ = [f64(-3.661316565040227801781974e+00), 8.64900232648597589369854e+00, + -2.8417941955033612725238097e-01, 3.253290185974728640827156e+00, + 2.147275640380293804770778e-01, -2.18600910711067004921551e-01, + -1.760002817872367935518928e+00, -4.389808914752818126249079e-01, + -3.843885560201130679995041e+00, 9.10988793377685105753416e-01] + // Results for 100000 * pi + vf_[i] + tan_large_ = [f64(-3.66131656475596512705e+00), 8.6490023287202547927e+00, + -2.841794195104782406e-01, 3.2532901861033120983e+00, 2.14727564046880001365e-01, + -2.18600910700688062874e-01, -1.760002817699722747043e+00, -4.38980891453536115952e-01, + -3.84388555942723509071e+00, 9.1098879344275101051e-01] + tanh_ = [f64(9.9990531206936338549262119e-01), 9.9999962057085294197613294e-01, + -2.7001505097318677233756845e-01, -9.9991110943061718603541401e-01, + 9.9999999146798465745022007e-01, 9.9427249436125236705001048e-01, + 9.9994257600983138572705076e-01, 9.9149409509772875982054701e-01, + 9.4936501296239685514466577e-01, -9.9999994291374030946055701e-01] + trunc_ = [f64(4.0000000000000000e+00), 7.0000000000000000e+00, copysign(0, -1), + -5.0000000000000000e+00, 9.0000000000000000e+00, 2.0000000000000000e+00, + 5.0000000000000000e+00, 2.0000000000000000e+00, 1.0000000000000000e+00, + -8.0000000000000000e+00, + ] +) + +fn tolerance(a f64, b f64, tol f64) bool { + mut ee := tol + // Multiplying by ee here can underflow denormal values to zero. + // Check a==b so that at least if a and b are small and identical + // we say they match. + if a == b { + return true + } + mut d := a - b + if d < 0 { + d = -d + } + // note: b is correct (expected) value, a is actual value. + // make error tolerance a fraction of b, not a. + if b != 0 { + ee = ee * b + if ee < 0 { + ee = -ee + } + } + return d < ee +} + +fn close(a f64, b f64) bool { + return tolerance(a, b, 1e-14) +} + +fn veryclose(a f64, b f64) bool { + return tolerance(a, b, 4e-16) +} + +fn soclose(a f64, b f64, e f64) bool { + return tolerance(a, b, e) +} + +fn alike(a f64, b f64) bool { + if is_nan(a) && is_nan(b) { + return true + } else if a == b { + return signbit(a) == signbit(b) + } + return false +} + +fn test_nan() { + nan_f64 := nan() + assert nan_f64 != nan_f64 + nan_f32 := f32(nan_f64) + assert nan_f32 != nan_f32 +} + +fn test_acos() { + for i := 0; i < math.vf_.len; i++ { + a := math.vf_[i] / 10 + f := acos(a) + assert soclose(math.acos_[i], f, 1e-7) + } + vfacos_sc_ := [-pi, 1, pi, nan()] + acos_sc_ := [nan(), 0, nan(), nan()] + for i := 0; i < vfacos_sc_.len; i++ { + f := acos(vfacos_sc_[i]) + assert alike(acos_sc_[i], f) + } +} + +fn test_asin() { + for i := 0; i < math.vf_.len; i++ { + a := math.vf_[i] / 10 + f := asin(a) + assert veryclose(math.asin_[i], f) + } + vfasin_sc_ := [-pi, copysign(0, -1), 0, pi, nan()] + asin_sc_ := [nan(), copysign(0, -1), 0, nan(), nan()] + for i := 0; i < vfasin_sc_.len; i++ { + f := asin(vfasin_sc_[i]) + assert alike(asin_sc_[i], f) + } +} + +fn test_atan() { + for i := 0; i < math.vf_.len; i++ { + f := atan(math.vf_[i]) + assert veryclose(math.atan_[i], f) + } + vfatan_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + atan_sc_ := [f64(-pi / 2), copysign(0, -1), 0, pi / 2, nan()] + for i := 0; i < vfatan_sc_.len; i++ { + f := atan(vfatan_sc_[i]) + assert alike(atan_sc_[i], f) + } +} + +fn test_atan2() { + for i := 0; i < math.vf_.len; i++ { + f := atan2(10, math.vf_[i]) + assert veryclose(math.atan2_[i], f) + } + vfatan2_sc_ := [[inf(-1), inf(-1)], [inf(-1), -pi], [inf(-1), 0], + [inf(-1), pi], [inf(-1), inf(1)], [inf(-1), nan()], [-pi, inf(-1)], + [-pi, 0], [-pi, inf(1)], [-pi, nan()], [f64(-0.0), inf(-1)], + [f64(-0.0), -pi], [f64(-0.0), -0.0], [f64(-0.0), 0], [f64(-0.0), pi], + [f64(-0.0), inf(1)], [f64(-0.0), nan()], [f64(0), inf(-1)], + [f64(0), -pi], [f64(0), -0.0], [f64(0), 0], [f64(0), pi], + [f64(0), inf(1)], [f64(0), nan()], [pi, inf(-1)], [pi, 0], + [pi, inf(1)], [pi, nan()], [inf(1), inf(-1)], [inf(1), -pi], + [inf(1), 0], [inf(1), pi], [inf(1), inf(1)], [inf(1), nan()], + [nan(), nan()], + ] + atan2_sc_ := [f64(-3.0) * pi / 4.0, /* atan2(-inf, -inf) */ -pi / 2, /* atan2(-inf, -pi) */ + -pi / 2, + /* atan2(-inf, +0) */ -pi / 2, /* atan2(-inf, pi) */ -pi / 4, /* atan2(-inf, +inf) */ + nan(), /* atan2(-inf, nan) */ -pi, /* atan2(-pi, -inf) */ -pi / 2, /* atan2(-pi, +0) */ + -0.0, + /* atan2(-pi, inf) */ nan(), /* atan2(-pi, nan) */ -pi, /* atan2(-0, -inf) */ -pi, + /* atan2(-0, -pi) */ -pi, /* atan2(-0, -0) */ -0.0, /* atan2(-0, +0) */ -0.0, /* atan2(-0, pi) */ + -0.0, + /* atan2(-0, +inf) */ nan(), /* atan2(-0, nan) */ pi, /* atan2(+0, -inf) */ pi, /* atan2(+0, -pi) */ + pi, /* atan2(+0, -0) */ 0, /* atan2(+0, +0) */ 0, /* atan2(+0, pi) */ 0, /* atan2(+0, +inf) */ + nan(), /* atan2(+0, nan) */ pi, /* atan2(pi, -inf) */ pi / 2, /* atan2(pi, +0) */ 0, + /* atan2(pi, +inf) */ nan(), /* atan2(pi, nan) */ 3.0 * pi / 4, /* atan2(+inf, -inf) */ + pi / 2, /* atan2(+inf, -pi) */ pi / 2, /* atan2(+inf, +0) */ pi / 2, /* atan2(+inf, pi) */ + pi / 4, /* atan2(+inf, +inf) */ nan(), /* atan2(+inf, nan) */ + nan(), /* atan2(nan, nan) */ + ] + for i := 0; i < vfatan2_sc_.len; i++ { + f := atan2(vfatan2_sc_[i][0], vfatan2_sc_[i][1]) + assert alike(atan2_sc_[i], f) + } +} + +fn test_ceil() { + // for i := 0; i < vf_.len; i++ { + // f := ceil(vf_[i]) + // assert alike(ceil_[i], f) + // } + vfceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + ceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + for i := 0; i < vfceil_sc_.len; i++ { + f := ceil(vfceil_sc_[i]) + assert alike(ceil_sc_[i], f) + } +} + +fn test_cos() { + for i := 0; i < math.vf_.len; i++ { + f := cos(math.vf_[i]) + assert veryclose(math.cos_[i], f) + } + vfcos_sc_ := [inf(-1), inf(1), nan()] + cos_sc_ := [nan(), nan(), nan()] + for i := 0; i < vfcos_sc_.len; i++ { + f := cos(vfcos_sc_[i]) + assert alike(cos_sc_[i], f) + } +} + +fn test_cosh() { + for i := 0; i < math.vf_.len; i++ { + f := cosh(math.vf_[i]) + assert close(math.cosh_[i], f) + } + vfcosh_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + cosh_sc_ := [inf(1), 1, 1, inf(1), nan()] + for i := 0; i < vfcosh_sc_.len; i++ { + f := cosh(vfcosh_sc_[i]) + assert alike(cosh_sc_[i], f) + } +} + +fn test_abs() { + for i := 0; i < math.vf_.len; i++ { + f := abs(math.vf_[i]) + assert math.fabs_[i] == f + } +} + +fn test_floor() { + for i := 0; i < math.vf_.len; i++ { + f := floor(math.vf_[i]) + assert alike(math.floor_[i], f) + } + vfceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + ceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + for i := 0; i < vfceil_sc_.len; i++ { + f := floor(vfceil_sc_[i]) + assert alike(ceil_sc_[i], f) + } +} + +fn test_max() { + for i := 0; i < math.vf_.len; i++ { + f := max(math.vf_[i], math.ceil_[i]) + assert math.ceil_[i] == f + } +} + +fn test_min() { + for i := 0; i < math.vf_.len; i++ { + f := min(math.vf_[i], math.floor_[i]) + assert math.floor_[i] == f + } +} + +fn test_signi() { + assert signi(inf(-1)) == -1 + assert signi(-72234878292.4586129) == -1 + assert signi(-10) == -1 + assert signi(-pi) == -1 + assert signi(-1) == -1 + assert signi(-0.000000000001) == -1 + assert signi(-0.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001) == -1 + assert signi(-0.0) == -1 + // + assert signi(inf(1)) == 1 + assert signi(72234878292.4586129) == 1 + assert signi(10) == 1 + assert signi(pi) == 1 + assert signi(1) == 1 + assert signi(0.000000000001) == 1 + assert signi(0.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001) == 1 + assert signi(0.0) == 1 + assert signi(nan()) == 1 +} + +fn test_sign() { + assert sign(inf(-1)) == -1.0 + assert sign(-72234878292.4586129) == -1.0 + assert sign(-10) == -1.0 + assert sign(-pi) == -1.0 + assert sign(-1) == -1.0 + assert sign(-0.000000000001) == -1.0 + assert sign(-0.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001) == -1.0 + assert sign(-0.0) == -1.0 + // + assert sign(inf(1)) == 1.0 + assert sign(72234878292.4586129) == 1 + assert sign(10) == 1.0 + assert sign(pi) == 1.0 + assert sign(1) == 1.0 + assert sign(0.000000000001) == 1.0 + assert sign(0.000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001) == 1.0 + assert sign(0.0) == 1.0 + assert is_nan(sign(nan())) + assert is_nan(sign(-nan())) +} + +fn test_exp() { + for i := 0; i < math.vf_.len; i++ { + f := exp(math.vf_[i]) + assert veryclose(math.exp_[i], f) + } + vfexp_sc_ := [inf(-1), -2000, 2000, inf(1), nan(), /* smallest f64 that overflows Exp(x) */ + 7.097827128933841e+02, 1.48852223e+09, 1.4885222e+09, 1, /* near zero */ + 3.725290298461915e-09, + /* denormal */ -740] + exp_sc_ := [f64(0), 0, inf(1), inf(1), nan(), inf(1), inf(1), + inf(1), 2.718281828459045, 1.0000000037252903, 4.2e-322] + for i := 0; i < vfexp_sc_.len; i++ { + f := exp(vfexp_sc_[i]) + assert alike(exp_sc_[i], f) + } +} + +fn test_exp2() { + for i := 0; i < math.vf_.len; i++ { + f := exp2(math.vf_[i]) + assert soclose(math.exp2_[i], f, 1e-9) + } + vfexp2_sc_ := [f64(-2000), 2000, inf(1), nan(), /* smallest f64 that overflows Exp2(x) */ + 1024, /* near underflow */ -1.07399999999999e+03, /* near zero */ 3.725290298461915e-09] + exp2_sc_ := [f64(0), inf(1), inf(1), nan(), inf(1), 5e-324, 1.0000000025821745] + for i := 0; i < vfexp2_sc_.len; i++ { + f := exp2(vfexp2_sc_[i]) + assert alike(exp2_sc_[i], f) + } +} + +fn test_gamma() { + vfgamma_ := [[inf(1), inf(1)], [inf(-1), nan()], [f64(0), inf(1)], + [f64(-0.0), inf(-1)], [nan(), nan()], [f64(-1), nan()], + [f64(-2), nan()], [f64(-3), nan()], [f64(-1e+16), nan()], + [f64(-1e+300), nan()], [f64(1.7e+308), inf(1)], /* Test inputs inspi_red by Python test suite. */ + // Outputs computed at high precision by PARI/GP. + // If recomputing table entries), be careful to use + // high-precision (%.1000g) formatting of the f64 inputs. + // For example), -2.0000000000000004 is the f64 with exact value + //-2.00000000000000044408920985626161695), and + // gamma(-2.0000000000000004) = -1249999999999999.5386078562728167651513), while + // gamma(-2.00000000000000044408920985626161695) = -1125899906826907.2044875028130093136826. + // Thus the table lists -1.1258999068426235e+15 as the answer. + [f64(0.5), 1.772453850905516], [f64(1.5), 0.886226925452758], + [f64(2.5), 1.329340388179137], [f64(3.5), 3.3233509704478426], + [f64(-0.5), -3.544907701811032], [f64(-1.5), 2.363271801207355], + [f64(-2.5), -0.9453087204829419], [f64(-3.5), 0.2700882058522691], + [f64(0.1), 9.51350769866873], [f64(0.01), 99.4325851191506], + [f64(1e-08), 9.999999942278434e+07], [f64(1e-16), 1e+16], + [f64(0.001), 999.4237724845955], [f64(1e-16), 1e+16], + [f64(1e-308), 1e+308], [f64(5.6e-309), 1.7857142857142864e+308], + [f64(5.5e-309), inf(1)], [f64(1e-309), inf(1)], [f64(1e-323), inf(1)], + [f64(5e-324), inf(1)], [f64(-0.1), -10.686287021193193], + [f64(-0.01), -100.58719796441078], [f64(-1e-08), -1.0000000057721567e+08], + [f64(-1e-16), -1e+16], [f64(-0.001), -1000.5782056293586], + [f64(-1e-16), -1e+16], [f64(-1e-308), -1e+308], [f64(-5.6e-309), -1.7857142857142864e+308], + [f64(-5.5e-309), inf(-1)], [f64(-1e-309), inf(-1)], [f64(-1e-323), inf(-1)], + [f64(-5e-324), inf(-1)], [f64(-0.9999999999999999), -9.007199254740992e+15], + [f64(-1.0000000000000002), 4.5035996273704955e+15], + [f64(-1.9999999999999998), + 2.2517998136852485e+15, + ], + [f64(-2.0000000000000004), -1.1258999068426235e+15], + [f64(-100.00000000000001), + -7.540083334883109e-145, + ], + [f64(-99.99999999999999), 7.540083334884096e-145], [f64(17), 2.0922789888e+13], + [f64(171), 7.257415615307999e+306], [f64(171.6), 1.5858969096672565e+308], + [f64(171.624), 1.7942117599248104e+308], [f64(171.625), inf(1)], + [f64(172), inf(1)], [f64(2000), inf(1)], [f64(-100.5), -3.3536908198076787e-159], + [f64(-160.5), -5.255546447007829e-286], [f64(-170.5), -3.3127395215386074e-308], + [f64(-171.5), 1.9316265431712e-310], [f64(-176.5), -1.196e-321], + [f64(-177.5), 5e-324], [f64(-178.5), -0.0], [f64(-179.5), 0], + [f64(-201.0001), 0], [f64(-202.9999), -0.0], [f64(-1000.5), -0.0], + [f64(-1.0000000003e+09), -0.0], [f64(-4.5035996273704955e+15), 0], + [f64(-63.349078729022985), 4.177797167776188e-88], + [f64(-127.45117632943295), + 1.183111089623681e-214, + ], + ] + _ := vfgamma_[0][0] + // for i := 0; i < math.vf_.len; i++ { + // f := gamma(math.vf_[i]) + // assert veryclose(math.gamma_[i], f) + // } + // for _, g in vfgamma_ { + // f := gamma(g[0]) + // if is_nan(g[1]) || is_inf(g[1], 0) || g[1] == 0 || f == 0 { + // assert alike(g[1], f) + // } else if g[0] > -50 && g[0] <= 171 { + // assert veryclose(g[1], f) + // } else { + // assert soclose(g[1], f, 1e-9) + // } + // } +} + +fn test_hypot() { + for i := 0; i < math.vf_.len; i++ { + a := abs(1e+200 * math.tanh_[i] * sqrt(2.0)) + f := hypot(1e+200 * math.tanh_[i], 1e+200 * math.tanh_[i]) + assert veryclose(a, f) + } + vfhypot_sc_ := [[inf(-1), inf(-1)], [inf(-1), 0], [inf(-1), + inf(1), + ], + [inf(-1), nan()], [f64(-0.0), -0.0], [f64(-0.0), 0], [f64(0), -0.0], + [f64(0), 0], /* +0,0 */ [f64(0), inf(-1)], [f64(0), inf(1)], + [f64(0), nan()], [inf(1), inf(-1)], [inf(1), 0], [inf(1), + inf(1), + ], + [inf(1), nan()], [nan(), inf(-1)], [nan(), 0], [nan(), + inf(1), + ], + [nan(), nan()], + ] + hypot_sc_ := [inf(1), inf(1), inf(1), inf(1), 0, 0, 0, 0, inf(1), + inf(1), nan(), inf(1), inf(1), inf(1), inf(1), inf(1), + nan(), inf(1), nan()] + for i := 0; i < vfhypot_sc_.len; i++ { + f := hypot(vfhypot_sc_[i][0], vfhypot_sc_[i][1]) + assert alike(hypot_sc_[i], f) + } +} + +fn test_log() { + for i := 0; i < math.vf_.len; i++ { + a := abs(math.vf_[i]) + f := log(a) + assert math.log_[i] == f + } + vflog_sc_ := [inf(-1), -pi, copysign(0, -1), 0, 1, inf(1), + nan(), + ] + log_sc_ := [nan(), nan(), inf(-1), inf(-1), 0, inf(1), nan()] + f := log(10) + assert f == ln10 + for i := 0; i < vflog_sc_.len; i++ { + g := log(vflog_sc_[i]) + assert alike(log_sc_[i], g) + } +} + +fn test_log10() { + for i := 0; i < math.vf_.len; i++ { + a := abs(math.vf_[i]) + f := log10(a) + assert veryclose(math.log10_[i], f) + } + vflog_sc_ := [inf(-1), -pi, copysign(0, -1), 0, 1, inf(1), + nan(), + ] + log_sc_ := [nan(), nan(), inf(-1), inf(-1), 0, inf(1), nan()] + for i := 0; i < vflog_sc_.len; i++ { + f := log10(vflog_sc_[i]) + assert alike(log_sc_[i], f) + } +} + +fn test_pow() { + for i := 0; i < math.vf_.len; i++ { + f := pow(10, math.vf_[i]) + assert close(math.pow_[i], f) + } + vfpow_sc_ := [[inf(-1), -pi], [inf(-1), -3], [inf(-1), -0.0], + [inf(-1), 0], [inf(-1), 1], [inf(-1), 3], [inf(-1), pi], + [inf(-1), 0.5], [inf(-1), nan()], [-pi, inf(-1)], [-pi, -pi], + [-pi, -0.0], [-pi, 0], [-pi, 1], [-pi, pi], [-pi, inf(1)], + [-pi, nan()], [f64(-1), inf(-1)], [f64(-1), inf(1)], [f64(-1), nan()], + [f64(-1 / 2), inf(-1)], [f64(-1 / 2), inf(1)], [f64(-0.0), inf(-1)], + [f64(-0.0), -pi], [f64(-0.0), -0.5], [f64(-0.0), -3], + [f64(-0.0), 3], [f64(-0.0), pi], [f64(-0.0), 0.5], [f64(-0.0), inf(1)], + [f64(0), inf(-1)], [f64(0), -pi], [f64(0), -3], [f64(0), -0.0], + [f64(0), 0], [f64(0), 3], [f64(0), pi], [f64(0), inf(1)], + [f64(0), nan()], [f64(1 / 2), inf(-1)], [f64(1 / 2), inf(1)], + [f64(1), inf(-1)], [f64(1), inf(1)], [f64(1), nan()], + [pi, inf(-1)], [pi, -0.0], [pi, 0], [pi, 1], [pi, inf(1)], + [pi, nan()], [inf(1), -pi], [inf(1), -0.0], [inf(1), 0], + [inf(1), 1], [inf(1), pi], [inf(1), nan()], [nan(), -pi], + [nan(), -0.0], [nan(), 0], [nan(), 1], [nan(), pi], [nan(), + nan(), + ]] + pow_sc_ := [f64(0), /* pow(-inf, -pi) */ -0.0, /* pow(-inf, -3) */ 1, /* pow(-inf, -0) */ 1, /* pow(-inf, +0) */ + inf(-1), /* pow(-inf, 1) */ inf(-1), /* pow(-inf, 3) */ + inf(1), /* pow(-inf, pi) */ inf(1), /* pow(-inf, 0.5) */ + nan(), /* pow(-inf, nan) */ 0, /* pow(-pi, -inf) */ nan(), /* pow(-pi, -pi) */ + 1, /* pow(-pi, -0) */ 1, /* pow(-pi, +0) */ -pi, /* pow(-pi, 1) */ nan(), /* pow(-pi, pi) */ + inf(1), /* pow(-pi, +inf) */ nan(), /* pow(-pi, nan) */ 1, /* pow(-1, -inf) IEEE 754-2008 */ + 1, /* pow(-1, +inf) IEEE 754-2008 */ nan(), /* pow(-1, nan) */ + inf(1), /* pow(-1/2, -inf) */ 0, /* pow(-1/2, +inf) */ inf(1), /* pow(-0, -inf) */ + inf(1), /* pow(-0, -pi) */ inf(1), /* pow(-0, -0.5) */ + inf(-1), /* pow(-0, -3) IEEE 754-2008 */ -0.0, /* pow(-0, 3) IEEE 754-2008 */ 0, /* pow(-0, pi) */ + 0, /* pow(-0, 0.5) */ 0, /* pow(-0, +inf) */ inf(1), /* pow(+0, -inf) */ + inf(1), /* pow(+0, -pi) */ inf(1), /* pow(+0, -3) */ 1, /* pow(+0, -0) */ 1, /* pow(+0, +0) */ + 0, /* pow(+0, 3) */ 0, + /* pow(+0, pi) */ 0, /* pow(+0, +inf) */ nan(), /* pow(+0, nan) */ + inf(1), /* pow(1/2, -inf) */ 0, /* pow(1/2, +inf) */ 1, /* pow(1, -inf) IEEE 754-2008 */ + 1, /* pow(1, +inf) IEEE 754-2008 */ 1, /* pow(1, nan) IEEE 754-2008 */ 0, /* pow(pi, -inf) */ + 1, /* pow(pi, -0) */ 1, /* pow(pi, +0) */ pi, /* pow(pi, 1) */ inf(1), /* pow(pi, +inf) */ + nan(), /* pow(pi, nan) */ 0, /* pow(+inf, -pi) */ 1, /* pow(+inf, -0) */ 1, /* pow(+inf, +0) */ + inf(1), /* pow(+inf, 1) */ inf(1), /* pow(+inf, pi) */ + nan(), /* pow(+inf, nan) */ nan(), /* pow(nan, -pi) */ 1, /* pow(nan, -0) */ 1, /* pow(nan, +0) */ + nan(), /* pow(nan, 1) */ nan(), /* pow(nan, pi) */ nan(), /* pow(nan, nan) */] + for i := 0; i < vfpow_sc_.len; i++ { + f := pow(vfpow_sc_[i][0], vfpow_sc_[i][1]) + assert alike(pow_sc_[i], f) + } +} + +fn test_round() { + for i := 0; i < math.vf_.len; i++ { + f := round(math.vf_[i]) + assert alike(math.round_[i], f) + } + vfround_sc_ := [[f64(0), 0], [nan(), nan()], [inf(1), inf(1)]] + vfround_even_sc_ := [[f64(0), 0], [f64(1.390671161567e-309), 0], /* denormal */ + [f64(0.49999999999999994), 0], /* 0.5-epsilon */ [f64(0.5), 0], + [f64(0.5000000000000001), 1], /* 0.5+epsilon */ [f64(-1.5), -2], + [f64(-2.5), -2], [nan(), nan()], [inf(1), inf(1)], + [f64(2251799813685249.5), 2251799813685250], + /* 1 bit fractian */ [f64(2251799813685250.5), 2251799813685250], + [f64(4503599627370495.5), 4503599627370496], /* 1 bit fraction, rounding to 0 bit fractian */ + [f64(4503599627370497), 4503599627370497], /* large integer */ + ] + _ := vfround_even_sc_[0][0] + for i := 0; i < vfround_sc_.len; i++ { + f := round(vfround_sc_[i][0]) + assert alike(vfround_sc_[i][1], f) + } +} + +fn test_sin() { + for i := 0; i < math.vf_.len; i++ { + f := sin(math.vf_[i]) + assert veryclose(math.sin_[i], f) + } + vfsin_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + sin_sc_ := [nan(), copysign(0, -1), 0, nan(), nan()] + for i := 0; i < vfsin_sc_.len; i++ { + f := sin(vfsin_sc_[i]) + assert alike(sin_sc_[i], f) + } +} + +fn test_sinh() { + for i := 0; i < math.vf_.len; i++ { + f := sinh(math.vf_[i]) + assert close(math.sinh_[i], f) + } + vfsinh_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + sinh_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + for i := 0; i < vfsinh_sc_.len; i++ { + f := sinh(vfsinh_sc_[i]) + assert alike(sinh_sc_[i], f) + } +} + +fn test_sqrt() { + for i := 0; i < math.vf_.len; i++ { + mut a := abs(math.vf_[i]) + mut f := sqrt(a) + assert veryclose(math.sqrt_[i], f) + a = abs(math.vf_[i]) + f = sqrt(a) + assert veryclose(math.sqrt_[i], f) + } + vfsqrt_sc_ := [inf(-1), -pi, copysign(0, -1), 0, inf(1), nan()] + sqrt_sc_ := [nan(), nan(), copysign(0, -1), 0, inf(1), nan()] + for i := 0; i < vfsqrt_sc_.len; i++ { + mut f := sqrt(vfsqrt_sc_[i]) + assert alike(sqrt_sc_[i], f) + f = sqrt(vfsqrt_sc_[i]) + assert alike(sqrt_sc_[i], f) + } +} + +fn test_tan() { + for i := 0; i < math.vf_.len; i++ { + f := tan(math.vf_[i]) + assert veryclose(math.tan_[i], f) + } + vfsin_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + sin_sc_ := [nan(), copysign(0, -1), 0, nan(), nan()] + // same special cases as sin + for i := 0; i < vfsin_sc_.len; i++ { + f := tan(vfsin_sc_[i]) + assert alike(sin_sc_[i], f) + } +} + +fn test_tanh() { + for i := 0; i < math.vf_.len; i++ { + f := tanh(math.vf_[i]) + assert veryclose(math.tanh_[i], f) + } + vftanh_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + tanh_sc_ := [f64(-1), copysign(0, -1), 0, 1, nan()] + for i := 0; i < vftanh_sc_.len; i++ { + f := tanh(vftanh_sc_[i]) + assert alike(tanh_sc_[i], f) + } +} + +fn test_trunc() { + // for i := 0; i < vf_.len; i++ { + // f := trunc(vf_[i]) + // assert alike(trunc_[i], f) + // } + vfceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + ceil_sc_ := [inf(-1), copysign(0, -1), 0, inf(1), nan()] + for i := 0; i < vfceil_sc_.len; i++ { + f := trunc(vfceil_sc_[i]) + assert alike(ceil_sc_[i], f) + } +} + +fn test_gcd() { + assert gcd(6, 9) == 3 + assert gcd(6, -9) == 3 + assert gcd(-6, -9) == 3 + assert gcd(0, 0) == 0 +} + +fn test_lcm() { + assert lcm(2, 3) == 6 + assert lcm(-2, 3) == 6 + assert lcm(-2, -3) == 6 + assert lcm(0, 0) == 0 +} + +fn test_digits() { + digits_in_10th_base := digits(125, 10) + assert digits_in_10th_base[0] == 5 + assert digits_in_10th_base[1] == 2 + assert digits_in_10th_base[2] == 1 + digits_in_16th_base := digits(15, 16) + assert digits_in_16th_base[0] == 15 + negative_digits := digits(-4, 2) + assert negative_digits[2] == -1 +} + +// Check that math functions of high angle values +// return accurate results. [since (vf_[i] + large) - large != vf_[i], +// testing for Trig(vf_[i] + large) == Trig(vf_[i]), where large is +// a multiple of 2 * pi, is misleading.] +fn test_large_cos() { + large := 100000.0 * pi + for i := 0; i < math.vf_.len; i++ { + f1 := math.cos_large_[i] + f2 := cos(math.vf_[i] + large) + assert soclose(f1, f2, 4e-9) + } +} + +fn test_large_sin() { + large := 100000.0 * pi + for i := 0; i < math.vf_.len; i++ { + f1 := math.sin_large_[i] + f2 := sin(math.vf_[i] + large) + assert soclose(f1, f2, 4e-9) + } +} + +fn test_large_tan() { + large := 100000.0 * pi + for i := 0; i < math.vf_.len; i++ { + f1 := math.tan_large_[i] + f2 := tan(math.vf_[i] + large) + assert soclose(f1, f2, 4e-9) + } +} diff --git a/v_windows/v/old/vlib/math/mathutil/mathutil.v b/v_windows/v/old/vlib/math/mathutil/mathutil.v new file mode 100644 index 0000000..0930e26 --- /dev/null +++ b/v_windows/v/old/vlib/math/mathutil/mathutil.v @@ -0,0 +1,19 @@ +// Copyright (c) 2019-2021 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 mathutil + +[inline] +pub fn min<T>(a T, b T) T { + return if a < b { a } else { b } +} + +[inline] +pub fn max<T>(a T, b T) T { + return if a > b { a } else { b } +} + +[inline] +pub fn abs<T>(a T) T { + return if a > 0 { a } else { -a } +} diff --git a/v_windows/v/old/vlib/math/mathutil/mathutil_test.v b/v_windows/v/old/vlib/math/mathutil/mathutil_test.v new file mode 100644 index 0000000..6b3b68a --- /dev/null +++ b/v_windows/v/old/vlib/math/mathutil/mathutil_test.v @@ -0,0 +1,22 @@ +import math.mathutil as mu + +fn test_min() { + assert mu.min(42, 13) == 13 + assert mu.min(5, -10) == -10 + assert mu.min(7.1, 7.3) == 7.1 + assert mu.min(u32(32), u32(17)) == 17 +} + +fn test_max() { + assert mu.max(42, 13) == 42 + assert mu.max(5, -10) == 5 + assert mu.max(7.1, 7.3) == 7.3 + assert mu.max(u32(60), u32(17)) == 60 +} + +fn test_abs() { + assert mu.abs(99) == 99 + assert mu.abs(-10) == 10 + assert mu.abs(1.2345) == 1.2345 + assert mu.abs(-5.5) == 5.5 +} diff --git a/v_windows/v/old/vlib/math/stats/stats.v b/v_windows/v/old/vlib/math/stats/stats.v new file mode 100644 index 0000000..d7317bf --- /dev/null +++ b/v_windows/v/old/vlib/math/stats/stats.v @@ -0,0 +1,249 @@ +module stats + +import math + +// TODO: Implement all of them with generics + +// This module defines the following statistical operations on f64 array +// --------------------------- +// | Summary of Functions | +// --------------------------- +// ----------------------------------------------------------------------- +// freq - Frequency +// mean - Mean +// geometric_mean - Geometric Mean +// harmonic_mean - Harmonic Mean +// median - Median +// mode - Mode +// rms - Root Mean Square +// population_variance - Population Variance +// sample_variance - Sample Variance +// population_stddev - Population Standard Deviation +// sample_stddev - Sample Standard Deviation +// mean_absdev - Mean Absolute Deviation +// min - Minimum of the Array +// max - Maximum of the Array +// range - Range of the Array ( max - min ) +// ----------------------------------------------------------------------- + +// Measure of Occurance +// Frequency of a given number +// Based on +// https://www.mathsisfun.com/data/frequency-distribution.html +pub fn freq(arr []f64, val f64) int { + if arr.len == 0 { + return 0 + } + mut count := 0 + for v in arr { + if v == val { + count++ + } + } + return count +} + +// Measure of Central Tendancy +// Mean of the given input array +// Based on +// https://www.mathsisfun.com/data/central-measures.html +pub fn mean(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut sum := f64(0) + for v in arr { + sum += v + } + return sum / f64(arr.len) +} + +// Measure of Central Tendancy +// Geometric Mean of the given input array +// Based on +// https://www.mathsisfun.com/numbers/geometric-mean.html +pub fn geometric_mean(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut sum := f64(1) + for v in arr { + sum *= v + } + return math.pow(sum, f64(1) / arr.len) +} + +// Measure of Central Tendancy +// Harmonic Mean of the given input array +// Based on +// https://www.mathsisfun.com/numbers/harmonic-mean.html +pub fn harmonic_mean(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut sum := f64(0) + for v in arr { + sum += f64(1) / v + } + return f64(arr.len) / sum +} + +// Measure of Central Tendancy +// Median of the given input array ( input array is assumed to be sorted ) +// Based on +// https://www.mathsisfun.com/data/central-measures.html +pub fn median(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + if arr.len % 2 == 0 { + mid := (arr.len / 2) - 1 + return (arr[mid] + arr[mid + 1]) / f64(2) + } else { + return arr[((arr.len - 1) / 2)] + } +} + +// Measure of Central Tendancy +// Mode of the given input array +// Based on +// https://www.mathsisfun.com/data/central-measures.html +pub fn mode(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut freqs := []int{} + for v in arr { + freqs << freq(arr, v) + } + mut max := 0 + for i in 0 .. freqs.len { + if freqs[i] > freqs[max] { + max = i + } + } + return arr[max] +} + +// Root Mean Square of the given input array +// Based on +// https://en.wikipedia.org/wiki/Root_mean_square +pub fn rms(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut sum := f64(0) + for v in arr { + sum += math.pow(v, 2) + } + return math.sqrt(sum / f64(arr.len)) +} + +// Measure of Dispersion / Spread +// Population Variance of the given input array +// Based on +// https://www.mathsisfun.com/data/standard-deviation.html +pub fn population_variance(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + m := mean(arr) + mut sum := f64(0) + for v in arr { + sum += math.pow(v - m, 2) + } + return sum / f64(arr.len) +} + +// Measure of Dispersion / Spread +// Sample Variance of the given input array +// Based on +// https://www.mathsisfun.com/data/standard-deviation.html +pub fn sample_variance(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + m := mean(arr) + mut sum := f64(0) + for v in arr { + sum += math.pow(v - m, 2) + } + return sum / f64(arr.len - 1) +} + +// Measure of Dispersion / Spread +// Population Standard Deviation of the given input array +// Based on +// https://www.mathsisfun.com/data/standard-deviation.html +pub fn population_stddev(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + return math.sqrt(population_variance(arr)) +} + +// Measure of Dispersion / Spread +// Sample Standard Deviation of the given input array +// Based on +// https://www.mathsisfun.com/data/standard-deviation.html +pub fn sample_stddev(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + return math.sqrt(sample_variance(arr)) +} + +// Measure of Dispersion / Spread +// Mean Absolute Deviation of the given input array +// Based on +// https://en.wikipedia.org/wiki/Average_absolute_deviation +pub fn mean_absdev(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + amean := mean(arr) + mut sum := f64(0) + for v in arr { + sum += math.abs(v - amean) + } + return sum / f64(arr.len) +} + +// Minimum of the given input array +pub fn min(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut min := arr[0] + for v in arr { + if v < min { + min = v + } + } + return min +} + +// Maximum of the given input array +pub fn max(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + mut max := arr[0] + for v in arr { + if v > max { + max = v + } + } + return max +} + +// Measure of Dispersion / Spread +// Range ( Maximum - Minimum ) of the given input array +// Based on +// https://www.mathsisfun.com/data/range.html +pub fn range(arr []f64) f64 { + if arr.len == 0 { + return f64(0) + } + return max(arr) - min(arr) +} diff --git a/v_windows/v/old/vlib/math/stats/stats_test.v b/v_windows/v/old/vlib/math/stats/stats_test.v new file mode 100644 index 0000000..c18daff --- /dev/null +++ b/v_windows/v/old/vlib/math/stats/stats_test.v @@ -0,0 +1,269 @@ +import math.stats +import math + +fn test_freq() { + // Tests were also verified on Wolfram Alpha + data := [f64(10.0), f64(10.0), f64(5.9), f64(2.7)] + mut o := stats.freq(data, 10.0) + assert o == 2 + o = stats.freq(data, 2.7) + assert o == 1 + o = stats.freq(data, 15) + assert o == 0 +} + +fn tst_res(str1 string, str2 string) bool { + if (math.abs(str1.f64() - str2.f64())) < 1e-5 { + return true + } + return false +} + +fn test_mean() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '5.762500') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '17.650000') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '37.708000') +} + +fn test_geometric_mean() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.geometric_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '5.15993') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.geometric_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert o.str() == 'nan' || o.str() == '-nan' || o.str() == '-1.#IND00' || o == f64(0) + || o.str() == '-nan(ind)' // Because in math it yields a complex number + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.geometric_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '25.064496') +} + +fn test_harmonic_mean() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.harmonic_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '4.626519') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.harmonic_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '9.134577') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.harmonic_mean(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '16.555477') +} + +fn test_median() { + // Tests were also verified on Wolfram Alpha + // Assumes sorted array + + // Even + mut data := [f64(2.7), f64(4.45), f64(5.9), f64(10.0)] + mut o := stats.median(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '5.175000') + data = [f64(-3.0), f64(1.89), f64(4.4), f64(67.31)] + o = stats.median(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '3.145000') + data = [f64(7.88), f64(12.0), f64(54.83), f64(76.122)] + o = stats.median(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '33.415000') + + // Odd + data = [f64(2.7), f64(4.45), f64(5.9), f64(10.0), f64(22)] + o = stats.median(data) + assert o == f64(5.9) + data = [f64(-3.0), f64(1.89), f64(4.4), f64(9), f64(67.31)] + o = stats.median(data) + assert o == f64(4.4) + data = [f64(7.88), f64(3.3), f64(12.0), f64(54.83), f64(76.122)] + o = stats.median(data) + assert o == f64(12.0) +} + +fn test_mode() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(2.7), f64(2.7), f64(4.45), f64(5.9), f64(10.0)] + mut o := stats.mode(data) + assert o == f64(2.7) + data = [f64(-3.0), f64(1.89), f64(1.89), f64(1.89), f64(9), f64(4.4), f64(4.4), f64(9), + f64(67.31), + ] + o = stats.mode(data) + assert o == f64(1.89) + // Testing greedy nature + data = [f64(2.0), f64(4.0), f64(2.0), f64(4.0)] + o = stats.mode(data) + assert o == f64(2.0) +} + +fn test_rms() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.rms(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '6.362046') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.rms(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '33.773393') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.rms(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '47.452561') +} + +fn test_population_variance() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.population_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '7.269219') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.population_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '829.119550') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.population_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '829.852282') +} + +fn test_sample_variance() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.sample_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '9.692292') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.sample_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '1105.492733') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.sample_variance(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '1106.469709') +} + +fn test_population_stddev() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.population_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '2.696149') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.population_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '28.794436') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.population_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '28.807157') +} + +fn test_sample_stddev() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.sample_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '3.113245') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.sample_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '33.248951') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.sample_stddev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '33.263639') +} + +fn test_mean_absdev() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.mean_absdev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '2.187500') + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.mean_absdev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '24.830000') + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.mean_absdev(data) + // Some issue with precision comparison in f64 using == operator hence serializing to string + assert tst_res(o.str(), '27.768000') +} + +fn test_min() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.min(data) + assert o == f64(2.7) + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.min(data) + assert o == f64(-3.0) + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.min(data) + assert o == f64(7.88) +} + +fn test_max() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.max(data) + assert o == f64(10.0) + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.max(data) + assert o == f64(67.31) + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.max(data) + assert o == f64(76.122) +} + +fn test_range() { + // Tests were also verified on Wolfram Alpha + mut data := [f64(10.0), f64(4.45), f64(5.9), f64(2.7)] + mut o := stats.range(data) + assert o == f64(7.3) + data = [f64(-3.0), f64(67.31), f64(4.4), f64(1.89)] + o = stats.range(data) + assert o == f64(70.31) + data = [f64(12.0), f64(7.88), f64(76.122), f64(54.83)] + o = stats.range(data) + assert o == f64(68.242) +} + +fn test_passing_empty() { + data := []f64{} + assert stats.freq(data, 0) == 0 + assert stats.mean(data) == f64(0) + assert stats.geometric_mean(data) == f64(0) + assert stats.harmonic_mean(data) == f64(0) + assert stats.median(data) == f64(0) + assert stats.mode(data) == f64(0) + assert stats.rms(data) == f64(0) + assert stats.population_variance(data) == f64(0) + assert stats.sample_variance(data) == f64(0) + assert stats.population_stddev(data) == f64(0) + assert stats.sample_stddev(data) == f64(0) + assert stats.mean_absdev(data) == f64(0) + assert stats.min(data) == f64(0) + assert stats.max(data) == f64(0) + assert stats.range(data) == f64(0) +} diff --git a/v_windows/v/old/vlib/math/unsafe.v b/v_windows/v/old/vlib/math/unsafe.v new file mode 100644 index 0000000..e6ebc6f --- /dev/null +++ b/v_windows/v/old/vlib/math/unsafe.v @@ -0,0 +1,38 @@ +// Copyright (c) 2019-2021 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 + +// f32_bits returns the IEEE 754 binary representation of f, +// with the sign bit of f and the result in the same bit position. +// f32_bits(f32_from_bits(x)) == x. +pub fn f32_bits(f f32) u32 { + p := *unsafe { &u32(&f) } + return p +} + +// f32_from_bits returns the floating-point number corresponding +// to the IEEE 754 binary representation b, with the sign bit of b +// and the result in the same bit position. +// f32_from_bits(f32_bits(x)) == x. +pub fn f32_from_bits(b u32) f32 { + p := *unsafe { &f32(&b) } + return p +} + +// f64_bits returns the IEEE 754 binary representation of f, +// with the sign bit of f and the result in the same bit position, +// and f64_bits(f64_from_bits(x)) == x. +pub fn f64_bits(f f64) u64 { + p := *unsafe { &u64(&f) } + return p +} + +// f64_from_bits returns the floating-point number corresponding +// to the IEEE 754 binary representation b, with the sign bit of b +// and the result in the same bit position. +// f64_from_bits(f64_bits(x)) == x. +pub fn f64_from_bits(b u64) f64 { + p := *unsafe { &f64(&b) } + return p +} diff --git a/v_windows/v/old/vlib/math/util/util.v b/v_windows/v/old/vlib/math/util/util.v new file mode 100644 index 0000000..0802d28 --- /dev/null +++ b/v_windows/v/old/vlib/math/util/util.v @@ -0,0 +1,82 @@ +// Copyright (c) 2019-2021 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 util + +// imin returns the smallest of two integer values +[inline] +pub fn imin(a int, b int) int { + return if a < b { a } else { b } +} + +// imin returns the biggest of two integer values +[inline] +pub fn imax(a int, b int) int { + return if a > b { a } else { b } +} + +// iabs returns an integer as absolute value +[inline] +pub fn iabs(v int) int { + return if v > 0 { v } else { -v } +} + +// umin returns the smallest of two u32 values +[inline] +pub fn umin(a u32, b u32) u32 { + return if a < b { a } else { b } +} + +// umax returns the biggest of two u32 values +[inline] +pub fn umax(a u32, b u32) u32 { + return if a > b { a } else { b } +} + +// uabs returns an u32 as absolute value +[inline] +pub fn uabs(v u32) u32 { + return if v > 0 { v } else { -v } +} + +// fmin_32 returns the smallest `f32` of input `a` and `b`. +// Example: assert fmin_32(2.0,3.0) == 2.0 +[inline] +pub fn fmin_32(a f32, b f32) f32 { + return if a < b { a } else { b } +} + +// fmax_32 returns the largest `f32` of input `a` and `b`. +// Example: assert fmax_32(2.0,3.0) == 3.0 +[inline] +pub fn fmax_32(a f32, b f32) f32 { + return if a > b { a } else { b } +} + +// fabs_32 returns the absolute value of `a` as a `f32` value. +// Example: assert fabs_32(-2.0) == 2.0 +[inline] +pub fn fabs_32(v f32) f32 { + return if v > 0 { v } else { -v } +} + +// fmin_64 returns the smallest `f64` of input `a` and `b`. +// Example: assert fmin_64(2.0,3.0) == 2.0 +[inline] +pub fn fmin_64(a f64, b f64) f64 { + return if a < b { a } else { b } +} + +// fmax_64 returns the largest `f64` of input `a` and `b`. +// Example: assert fmax_64(2.0,3.0) == 3.0 +[inline] +pub fn fmax_64(a f64, b f64) f64 { + return if a > b { a } else { b } +} + +// fabs_64 returns the absolute value of `a` as a `f64` value. +// Example: assert fabs_64(-2.0) == f64(2.0) +[inline] +pub fn fabs_64(v f64) f64 { + return if v > 0 { v } else { -v } +} |