Category talk:Wren-math: Difference between revisions
Content added Content deleted
(→Source code: Added Math.lerp method) |
(Added modMul, modPow, modInv and isPrime2 (using Miller-Rabin) methods to Int class.) |
||
Line 86: | Line 86: | ||
} |
} |
||
// Returns the value of a polynomial when the |
// Returns the value of a polynomial when the unknown is 'x' using Horner's rule. |
||
// Polynomials are represented by an ordered list of coefficients |
// Polynomials are represented by an ordered list of coefficients |
||
// from the term with the highest degree down to the constant term. |
// from the term with the highest degree down to the constant term. |
||
Line 227: | Line 227: | ||
if (a.count == 2) return l |
if (a.count == 2) return l |
||
return lcm(a[2..-1] + [l]) |
return lcm(a[2..-1] + [l]) |
||
} |
} |
||
// Returns 'x' multiplied by 'n' modulo 'm'. |
|||
static modMul(x, n, m) { |
|||
if (m == 0) Fiber.abort("Cannot take modMul with modulus 0.") |
|||
var sign = x.sign * n.sign |
|||
var z = 0 |
|||
x = x.abs |
|||
n = n.abs |
|||
m = m.abs |
|||
while (n > 0) { |
|||
if (n % 2 == 1) z = (z + x) % m |
|||
x = (x * 2) % m |
|||
n = div(n, 2) |
|||
} |
|||
return z * sign |
|||
} |
|||
// Returns the remainder when 'b' raised to the power 'e' is divided by 'm'. |
// Returns the remainder when 'b' raised to the power 'e' is divided by 'm'. |
||
static modPow(b, e, m) { |
static modPow(b, e, m) { |
||
if (m == |
if (m == 0) Fiber.abort("Cannot take modPow with modulus 0.") |
||
if (m == 1 || m == -1) return 0 |
|||
var r = 1 |
var r = 1 |
||
b = b % m |
b = b % m |
||
if (e < 0) { |
|||
e = -e |
|||
b = modInv(b, m) |
|||
} |
|||
while (e > 0) { |
while (e > 0) { |
||
if ( |
if (b == 0) return 0 |
||
e |
if (e % 2 == 1) r = modMul(r, b, m) |
||
e = div(e, 2) |
|||
b = modMul(b, b, m) |
|||
} |
} |
||
return r |
return r |
||
} |
|||
// Returns the multiplicative inverse of 'x' modulo 'n'. |
|||
static modInv(x, n) { |
|||
var r = n |
|||
var newR = x.abs |
|||
var t = 0 |
|||
var newT = 1 |
|||
while (newR != 0) { |
|||
var q = quo(r, newR) |
|||
var lastT = t |
|||
var lastR = r |
|||
t = newT |
|||
r = newR |
|||
newT = lastT - q * newT |
|||
newR = lastR - q * newR |
|||
} |
|||
if (r != 1) Fiber.abort("%(x) and %(n) are not co-prime.") |
|||
if (t < 0) t = t + n |
|||
return (x < 0) ? -t : t |
|||
} |
} |
||
Line 285: | Line 327: | ||
if (n%d == 0) return false |
if (n%d == 0) return false |
||
d = d + 4 |
d = d + 4 |
||
} |
|||
return true |
|||
} |
|||
// Returns true if the current instance is prime, false otherwise but uses the |
|||
// Miller-Rabin test which should be faster for checking isolated larger integers. |
|||
static isPrime2(n) { |
|||
if (!n.isInteger || n < 2) return false |
|||
if (n%2 == 0) return n == 2 |
|||
if (n%3 == 0) return n == 3 |
|||
if (n%5 == 0) return n == 5 |
|||
if (n%7 == 0) return n == 7 |
|||
if (n < 121) return true |
|||
var a |
|||
if (n < 2047) { |
|||
a = [2] |
|||
} else if (n < 1373653) { |
|||
a = [2, 3] |
|||
} else if (n < 9080191) { |
|||
a = [31, 73] |
|||
} else if (n < 25326001) { |
|||
a = [2, 3, 5] |
|||
} else if (n < 3215031751) { |
|||
a = [2, 3, 5, 7] |
|||
} else if (n < 4759123141) { |
|||
a = [2, 7, 61] |
|||
} else if (n < 1122004669633) { |
|||
a = [2, 13, 23, 1662803] |
|||
} else if (n < 2152302898747) { |
|||
a = [2, 3, 5, 7, 11] |
|||
} else if (n < 3474749660383) { |
|||
a = [2, 3, 5, 7, 11, 13] |
|||
} else if (n < 341550071728321) { |
|||
a = [2, 3, 5, 7, 11, 13, 17] |
|||
} else { |
|||
a = [2, 3, 5, 7, 11, 13, 17, 19, 23] |
|||
} |
|||
var nPrev = n - 1 |
|||
var b = nPrev |
|||
var r = 0 |
|||
while (b % 2 == 0) { |
|||
b = div(b, 2) |
|||
r = r + 1 |
|||
} |
|||
for (i in 0...a.count) { |
|||
if (n >= a[i]) { |
|||
var x = modPow(a[i], b, n) |
|||
if (x != 1 && x != nPrev) { |
|||
var d = r - 1 |
|||
var next = false |
|||
while (d != 0) { |
|||
x = modMul(x, x, n) |
|||
if (x == 1) return false |
|||
if (x == nPrev) { |
|||
next = true |
|||
break |
|||
} |
|||
d = d - 1 |
|||
} |
|||
if (!next) return false |
|||
} |
|||
} |
|||
} |
} |
||
return true |
return true |