Category talk:Wren-big: Difference between revisions

m
→‎Source code: Now uses Wren S/H lexer.
(→‎Source code: Added iroot and icbrt methods to BigInt class.)
m (→‎Source code: Now uses Wren S/H lexer.)
 
(33 intermediate revisions by the same user not shown)
Line 5:
Wren is unable to deal natively and precisely with integers outside a 53-bit range and so either has to limit such tasks to integers which it can deal with or is unable to offer a realistic solution at all.
 
This module aims to remedy that situation by providing classes for arbitrary precision arithmetic (integers, rational and rationaldecimal numbers) in an easy to use form. It is based on or guided by Peter Olson's public domain [https://github.com/peterolson/BigInteger.js BigInteger.js] and [https://github.com/peterolson/BigRational.js BigRational.js] libraries for Javascript, a language which until the recent addition of native BigInt support was in the same boat as Wren regarding the 53-bit limitation.
 
;BigInt class
Line 23:
;BigRat class
 
This follows the lines of the Rat class in the [https[://rosettacode.org/wiki/Category:Wren-rat |Wren-rat]] module though I have been guided by [https://github.com/peterolson/BigRational.js BigRational.js] when coding some of the trickier parts.
 
The ''BigRat.new'' constructor accepts a BigInt numerator and denominator (or just a numerator on its own). However, you can pass in Nums or Strings if a BigInt can be built from them
Line 29:
There are also separate methods to enable you to to generate a BigRat from a Rat object, from a string in rational form (i.e. "n/d"), from a mixed string (e.g. "1_2/3") from a decimal string or value (e.g. "1.234") or from a Num.
 
;BigDec class
As rational numbers can easily be created from or displayed as decimal numbers, there seems little point in adding a BigDecimal class to the module as well.
 
As rational numbers can easily be created from or displayed as decimal numbers, at first glance there seems little point in adding a BigDec class to the module as well.
 
However, a problem with BigRat is that it always works to maximum precision which uses a lot of memory and can be very slow in certain circumstances.
 
BigDec aims to fix this problem by working to a user defined precision, where precision is measured in digits after the decimal point. One can therefore limit precision having regard to the accuracy needed for a particular task though, to dove-tail with the maximum accuracy of ordinary Nums, a minimum (and default) precision of 16 is imposed.
 
Otherwise, BigDec is simply a thin wrapper over BigRat and defers to it for the implementation of most of its methods. Where all you need is decimal numbers, it is also somewhat easier to use than BigRat as the user never has to deal with the latter directly.
 
;General observations
 
The original Javascript libraries and this module are easy to use as they try to mimic native arithmetic by always producing immutable objects and by implicitly converting suitable operands into BigInts, BigRats or BigRatsBigDecs. This is accentuated further by Wren's support for operator overloading which means that all the familiar operations on numbers can be used.
 
Of course, such a system inevitably keeps the garbage collector busy and is not going to be as fast and efficient as systems where BigInts, BigRats or BigRatsBigDecs are mutable and where memory usage can therefore be managed by the user. On the other hand the latter systems can be difficult or unpleasant to use and it is therefore felt that the present approach is a good fit for Wren where ease of use is paramount.
 
===Source code===
<langsyntaxhighlight ecmascriptlang="wren">/* Module "big.wren" */
 
import "./trait" for Comparable
import "random" for Random
 
Line 59 ⟶ 67:
static four { BigInt.small_( 4) }
static five { BigInt.small_( 5) }
static six { BigInt.small_( 6) }
static seven { BigInt.small_( 7) }
static eight { BigInt.small_( 8) }
static nine { BigInt.small_( 9) }
static ten { BigInt.small_(10) }
 
Line 639 ⟶ 651:
if (!(b is BigInt)) b = BigInt.new(b)
return (a < b) ? a : b
}
 
// Returns the positive difference of two BigInts.
static dim(a, b) {
if (!(a is BigInt)) a = BigInt.new(a)
if (!(b is BigInt)) b = BigInt.new(b)
if (a >= b) return a - b
return zero
}
 
Line 675 ⟶ 695:
if (!(a is BigInt)) a = BigInt.new(a)
if (!(b is BigInt)) b = BigInt.new(b)
if (a.isZero && b.isZero) return BigInt.zero
return a / gcd(a, b) * b
}
 
// Returns the factorial of 'n'.
static factorial(n) {
if (!(n is Num && n >= 0)) Fiber.abort("Argument must be a non-negative integer")
if (n < 2) return BigInt.one
var fact = BigInt.one
for (i in 2..n) fact = fact * i
return fact
}
 
// Returns the multinomial coefficient of n over a list f where sum(f) == n.
static multinomial(n, f) {
if (!(n is Num && n >= 0)) Fiber.abort("First argument must be a non-negative integer")
if (!(f is List)) Fiber.abort("Second argument must be a list.")
var sum = f.reduce { |acc, i| acc + i }
if (n != sum) {
Fiber.abort("The elements of the list must sum to 'n'.")
}
var prod = BigInt.one
for (e in f) {
if (e < 0) Fiber.abort("The elements of the list must be non-negative integers.")
if (e > 1) prod = prod * BigInt.factorial(e)
}
return BigInt.factorial(n)/prod
}
 
// Returns the binomial coefficent of n over k.
static binomial(n, k) { multinomial(n, [k, n-k]) }
 
// Returns whether or not 'n' is an instance of BigInt.
Line 915 ⟶ 964:
isPrime { isPrime(false) }
isProbablePrime { isProbablePrime(5) }
 
// Returns the next probable prime number greater than 'this'.
nextPrime(iterations) {
if (this < BigInt.two) return BigInt.two
var n = isEven ? this + 1 : this + 2
while (true) {
if (n.isProbablePrime(iterations)) return n
n = n + 2
}
}
 
// Returns the previous probable prime number less than 'this' or null if there isn't one.
prevPrime(iterations) {
if (this < BigInt.three) return null
if (this == BigInt.three) return BigInt.two
var n = isEven ? this - 1 : this - 2
while (true) {
if (n.isProbablePrime(iterations)) return n
n = n - 2
}
}
 
// Convenience versions of the above methods which use 5 iterations.
nextPrime { nextPrime(5) }
prevPrime { prevPrime(5) }
 
// Negates a BigInt.
Line 935 ⟶ 1,009:
+ (n) {
if (!(n is BigInt)) n = BigInt.new(n)
if (n.isZero) return this.copy()
if (this.isSmall) return smallAdd_(n)
if (_signed != n.signed_) return this - (-n)
Line 958 ⟶ 1,032:
- (n) {
if (!(n is BigInt)) n = BigInt.new(n)
if (n.isZero) return this.copy()
if (this.isSmall) return smallSubtract_(n)
if (_signed != n.signed_) return this + (-n)
Line 976 ⟶ 1,050:
}
if (a.value_ == 0) return BigInt.zero
if (a.value_ == 1) return this.copy()
if (a.value_ == -1) return -this
return BigInt.multiplySmallAndList_(a.value_.abs, _value, _signed != a.signed_)
Line 990 ⟶ 1,064:
if (n.isSmall) {
if (b == 0) return BigInt.zero
if (b == 1) return this.copy()
if (b == -1) return -this
var ab = b.abs
Line 1,012 ⟶ 1,086:
}
 
// Returns the integer n'th root ofCube the current instance .
cube { square * this }
 
// Returns true if the current instance is a perfect square, false otherwise.
isSquare {
var root = isqrt
return root.square == this
}
 
// Returns true if the current instance is a perfect cube, false otherwise.
isCube {
var root = icbrt
return root.square * root == this
}
 
// Returns the integer n'th root of the current instance
// i.e. the precise n'th root truncated towards zero if not an integer.
// Throws an error if the current instance is negative and n is even.
Line 1,019 ⟶ 1,108:
Fiber.abort("Argument must be a positive integer.")
}
if (n == 1) return this.copy()
var t = copy()
var neg = t < 0
Line 1,039 ⟶ 1,128:
}
 
// Returns the cube root of the current instance i.e. the largest integer 'r' such that
icbrt// {r.cube <= this.
icbrt {
if (isSmall) return BigInt.small_(toSmall.cbrt.floor)
return iroot(3)
}
 
Line 1,057 ⟶ 1,147:
a = b
}
}
// Returns a list containing the quotient and the remainder after dividing the current instance
Line 1,262 ⟶ 1,352:
}
 
// Returns true if the 'n'th bit of the current instance is set or false otherwise.
testBit(n) {
if (n.type != Num || !n.isInteger || n < 0) Fiber.abort("Argument must be a non-negative integer.")
Line 1,305 ⟶ 1,395:
var sign = _signed ? "-" : ""
return sign + str
}
 
// Creates and returns a range of BigInts from 'this' to 'n' with a step of 1.
// 'this' and 'n' must both be 'small' integers. 'n' can be a Num or a BigInt.
..(n) { BigInts.range(this, n, 1, true) } // inclusive of 'n'
...(n) { BigInts.range(this, n, 1, false) } // exclusive of 'n'
 
/* Prime factorization methods. */
 
// Private worker method for Pollard's Rho algorithm.
static pollardRho_(m, seed, c) {
var g = Fn.new { |x| (x*x + c) % m }
var x = BigInt.new(seed)
var y = BigInt.new(seed)
var z = BigInt.one
var d = BigInt.one
var count = 0
while (true) {
x = g.call(x)
y = g.call(g.call(y))
d = (x - y).abs % m
z = z * d
count = count + 1
if (count == 100) {
d = BigInt.gcd(z, m)
if (d != BigInt.one) break
z = BigInt.one
count = 0
}
}
if (d == m) return BigInt.zero
return d
}
 
// Returns a factor (BigInt) of 'm' (a BigInt or an integral Num) using the
// Pollard's Rho algorithm. Both the 'seed' and 'c' can be set to integers.
// Returns BigInt.zero in the event of failure.
static pollardRho(m, seed, c) {
if (m < 2) return BigInt.zero
if (m is Num) m = BigInt.new(m)
if (m.isProbablePrime(10)) return m.copy()
if (m.isSquare) return m.isqrt
for (p in [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31]) {
if (m.isDivisibleBy(p)) return BigInt.new(p)
}
return pollardRho_(m, seed, c)
}
 
// Convenience version of the above method which uses a seed of 2 and a value for c of 1.
static pollardRho(m) { pollardRho(m, 2, 1) }
 
// Private method for factorizing smaller numbers (BigInt) using a wheel with basis [2, 3, 5].
static primeFactorsWheel_(m) {
var n = m.copy()
var inc = [4, 2, 4, 2, 4, 6, 2, 6]
var factors = []
var k = BigInt.new(37)
var i = 0
while (k * k <= n) {
if (n.isDivisibleBy(k)) {
factors.add(k.copy())
n = n / k
} else {
k = k + inc[i]
i = (i + 1) % 8
}
}
if (n > BigInt.one) factors.add(n)
return factors
}
 
// Private worker method (recursive) to obtain the prime factors of a number (BigInt).
static primeFactors_(m, trialDivs) {
if (m.isProbablePrime(10)) return [m.copy()]
var n = m.copy()
var factors = []
var seed = 2
var c = 1
var checkPrime = true
var threshold = 1e11 // from which using PR may be advantageous
if (trialDivs) {
for (p in [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31]) {
while (n.isDivisibleBy(p)) {
factors.add(BigInt.new(p))
n = n / p
}
}
}
while (n > BigInt.one) {
if (checkPrime && n.isProbablePrime(10)) {
factors.add(n)
break
}
if (n >= threshold) {
var d = pollardRho_(n, seed, c)
if (d != BigInt.zero) {
factors.addAll(primeFactors_(d, false))
n = n / d
checkPrime = true
} else if (c == 1) {
if (n.isSquare) {
n = n.isqrt
var pf = primeFactors_(n, false)
factors.addAll(pf)
factors.addAll(pf)
break
} else {
c = 2
checkPrime = false
}
} else if (c < 101) {
c = c + 1
} else if (seed < 101) {
seed = seed + 1
} else {
factors.addAll(primeFactorsWheel_(n))
break
}
} else {
factors.addAll(primeFactorsWheel_(n))
break
}
}
factors.sort()
return factors
}
 
// Returns a list of the primes factors (BigInt) of 'm' (a Bigint or an integral Num)
// using the wheel based factorization and/or Pollard's Rho algorithm as appropriate.
static primeFactors(m) {
if (m < 2) return []
if (m is Num) m = BigInt.new(m)
return primeFactors_(m, true)
}
 
// Returns all the divisors of 'n' including 1 and 'n' itself.
static divisors(n) {
if (n < 1) return []
if (n is Num) n = BigInt.new(n)
var divs = []
var divs2 = []
var i = one
var k = (n.isEven) ? one : two
var sqrt = n.isqrt
while (i <= sqrt) {
if (n.isDivisibleBy(i)) {
divs.add(i)
var j = n / i
if (j != i) divs2.add(j)
}
i = i + k
}
if (!divs2.isEmpty) divs = divs + divs2[-1..0]
return divs
}
 
// Returns all the divisors of 'n' excluding 'n'.
static properDivisors(n) {
var d = divisors(n)
var c = d.count
return (c <= 1) ? [] : d[0..-2]
}
 
// As 'divisors' method but uses a different algorithm.
// Better for larger numbers.
static divisors2(n) {
if (n is Num) n = BigInt.new(n)
var pf = primeFactors(n)
if (pf.count == 0) return (n == 1) ? [one] : pf
var arr = []
if (pf.count == 1) {
arr.add([pf[0].copy(), 1])
} else {
var prevPrime = pf[0]
var count = 1
for (i in 1...pf.count) {
if (pf[i] == prevPrime) {
count = count + 1
} else {
arr.add([prevPrime.copy(), count])
prevPrime = pf[i]
count = 1
}
}
arr.add([prevPrime.copy(), count])
}
var divisors = []
var generateDivs
generateDivs = Fn.new { |currIndex, currDivisor|
if (currIndex == arr.count) {
divisors.add(currDivisor.copy())
return
}
for (i in 0..arr[currIndex][1]) {
generateDivs.call(currIndex+1, currDivisor)
currDivisor = currDivisor * arr[currIndex][0]
}
}
generateDivs.call(0, one)
return divisors.sort()
}
 
// As 'properDivisors' but uses 'divisors2' method.
static properDivisors2(n) {
var d = divisors2(n)
var c = d.count
return (c <= 1) ? [] : d[0..-2]
}
 
// Returns the sum of all the divisors of 'n' including 1 and 'n' itself.
static divisorSum(n) {
if (n < 1) return zero
if (n is Num) n = BigInt.new(n)
var total = one
var power = two
while (n % 2 == 0) {
total = total + power
power = power * 2
n = n / 2
}
var i = three
while (i * i <= n) {
var sum = one
power = i
while (n % i == 0) {
sum = sum + power
power = power * i
n = n / i
}
total = total * sum
i = i + 2
}
if (n > 1) total = total * (n + 1)
return total
}
 
// Returns the number of divisors of 'n' including 1 and 'n' itself.
static divisorCount(n) {
if (n < 1) return 0
if (n is Num) n = BigInt.new(n)
var count = 0
var prod = 1
while (n % 2 == 0) {
count = count + 1
n = n / 2
}
prod = prod * (1 + count)
var i = three
while (i * i <= n) {
count = 0
while (n % i == 0) {
count = count + 1
n = n / i
}
prod = prod * (1 + count)
i = i + 2
}
if (n > 2) prod = prod * 2
return prod
}
}
 
/* BigInts contains various routines applicable to lists of big integers. */
and for creating and iterating though ranges of such numbers. */
class BigInts {
class BigInts is Sequence {
static sum(a) { a.reduce(BigInt.zero) { |acc, x| acc + x } }
static prod(a) { a.reduce(BigInt.one) { |acc, x| acc * x } }
static max(a) { a.reduce { |acc, x| (x > acc) ? x : acc } }
static min(a) { a.reduce { |acc, x| (x < acc) ? x : acc } }
 
// Private helper method for creating ranges.
static checkValue_(v, name) {
if (v is BigInt && v.isSmall) {
return v.toSmall
} else if (v is Num && v.isInteger && v.abs <= Num.maxSafeInteger) {
return v
} else {
Fiber.abort("Invalid value for '%(name)'.")
}
}
 
// Creates a range of 'small' BigInts analogous to the Range class but allowing for steps > 1.
// Use the 'ascending' parameter to check that the range's direction is as intended.
construct range(from, to, step, inclusive, ascending) {
from = BigInts.checkValue_(from, "from")
to = BigInts.checkValue_(to, "to")
step = BigInts.checkValue_(step, "step")
if (step < 1) Fiber.abort("'step' must be a positive integer.")
if (ascending && from > to) Fiber.abort("'from' cannot exceed 'to'.")
if (!ascending && from < to) Fiber.abort("'to' cannot exceed 'from'.")
_rng = inclusive ? from..to : from...to
_step = step
}
 
// Convenience versions of 'range' which use default values for some parameters.
static range(from, to, step, inclusive) { range(from, to, step, inclusive, from <= to) }
static range(from, to, step) { range(from, to, step, true, from <= to) }
static range(from, to) { range(from, to, 1, true, from <= to) }
 
// Self-explanatory read only properties.
from { _rng.from }
to { _rng.to }
step { _step }
min { from.min(to) }
max { from.max(to) }
isInclusive { _rng.isInclusive }
isAscending { from <= to }
 
// Iterator protocol methods.
iterate(iterator) {
if (!iterator || _step == 1) {
return _rng.iterate(iterator)
} else {
var count = _step
while (count > 0 && iterator) {
iterator = _rng.iterate(iterator)
count = count - 1
}
return iterator
}
}
 
iteratorValue(iterate) { BigInt.small_(_rng.iteratorValue(iterate)) }
}
 
Line 1,333 ⟶ 1,737:
 
// Constants.
static minusOne { BigRat.new(BigInt.minusOne, BigInt.one) }
static zero { BigRat.new(BigInt.zero, BigInt.one) }
static one { BigRat.new(BigInt.one, BigInt.one) }
static two { BigRat.new(BigInt.two, BigInt.one) }
static ten three { BigRat.new(BigInt.tenthree, BigInt.one) }
static halffour { BigRat.new(BigInt.onefour, BigInt.twoone) }
static tenthfive { BigRat.new(BigInt.onefive, BigInt.tenone) }
static six { BigRat.new(BigInt.six, BigInt.one) }
static seven { BigRat.new(BigInt.seven, BigInt.one) }
static eight { BigRat.new(BigInt.eight, BigInt.one) }
static nine { BigRat.new(BigInt.nine, BigInt.one) }
static ten { BigRat.new(BigInt.ten, BigInt.one) }
static half { BigRat.new(BigInt.one, BigInt.two) }
static third { BigRat.new(BigInt.one, BigInt.three) }
static quarter { BigRat.new(BigInt.one, BigInt.four) }
static fifth { BigRat.new(BigInt.one, BigInt.five) }
static sixth { BigRat.new(BigInt.one, BigInt.six) }
static seventh { BigRat.new(BigInt.one, BigInt.seven) }
static eighth { BigRat.new(BigInt.one, BigInt.eight) }
static ninth { BigRat.new(BigInt.one, BigInt.nine) }
static tenth { BigRat.new(BigInt.one, BigInt.ten) }
 
// Constructs a new BigRat object by passing it a numerator and a denominator.
Line 1,472 ⟶ 1,890:
 
// Rounding methods (similar to those in Num class).
ceil { // higher integer, towards zero
if (isInteger) return this.copy()
var div = _n/_d
if (!this.isNegative) div = div.inc
Line 1,480 ⟶ 1,898:
 
floor { // lower integer
if (isInteger) return this.copy()
var div = _n/_d
if (this.isNegative) div = div.dec
Line 1,489 ⟶ 1,907:
 
round { // nearer integer
if (isInteger) return this.copy()
var div = _n / _d
if (_d == 2) {
Line 1,496 ⟶ 1,914:
}
return (this + BigRat.half).floor
}
 
roundUp { this >= 0 ? ceil : floor } // rounds to higher integer, away from zero
 
round(digits) { // rounds to 'digits' decimal places
if (digits == 0) return round
return BigRat.fromDecimal(toDecimal(digits))
}
 
Line 1,519 ⟶ 1,944:
pow(i) {
if (!((i is Num) && i.isInteger)) Fiber.abort("Argument must be an integer.")
if (i == 0) return thisBigRat.one
if (i == 1) return this.copy()
if (i == -1) return this.inverse
var np = _n.pow(i.abs)
var dp = _d.pow(i.abs)
Line 1,527 ⟶ 1,954:
// Returns the square of the current instance.
square { BigRat.new(_n * _n , _d *_d) }
 
// Returns the cube of the current instance.
cube { square * this }
 
// Returns the square root of the current instance to 'digits' decimal places.
Line 1,536 ⟶ 1,966:
digits = digits + 5
var powd = BigInt.ten.pow(digits)
var sqtd = (powd.square * _n / _d).isqrtiroot(2)
return BigRat.new(sqtd, powd)
}
Line 1,542 ⟶ 1,972:
// Convenience version of the above method which uses 14 decimal places.
sqrt { sqrt(14) }
 
// Returns the cube root of the current instance to 'digits' decimal places.
// Five more decimals is used to try to ensure accuracy though this is not guaranteed.
cbrt(digits) {
if (!((digits is Num) && digits.isInteger && digits >= 0)) {
Fiber.abort("Digits must be a non-negative integer.")
}
digits = digits + 5
var powd = BigInt.ten.pow(digits)
var cbtd = (powd.cube * _n / _d).iroot(3)
return BigRat.new(cbtd, powd)
}
 
// Convenience version of the above method which uses 14 decimal places.
cbrt { cbrt(14) }
 
// Other methods.
inc { this + BigRat.one } // increment
dec { this - BigRat.one } // decrement
abs { (_n >= BigInt.zero) ? thiscopy() : -this } // absolute value
sign { _n.sign } // sign
 
// The inherited 'clone' method just returns 'this' as BigRat objects are immutable.
Line 1,569 ⟶ 2,014:
// Converts the current instance to a Num where possible.
// Will probably lose accuracy if the numerator and/or denominator are not 'small'.
toFloat { Num_n.fromString(thistoNum / _d.toDecimal(14))toNum }
 
// Converts the current instance to an integer where possible with any fractional part truncated.
Line 1,627 ⟶ 2,072:
// Returns a string represenation of this instance in the form "i_n/d" where 'i' is an integer.
toMixedString {
var qsign = _n.isNegative /? _d"-" : ""
var rnn = _n % _d.abs
ifvar (r.isNegative) rq = -rnn / _d
return q.toString + "_" +var r.toString += "/"nn +% _d.toString
return sign + q.toString + "_" + r.toString + "/" + _d.toString
}
 
Line 1,646 ⟶ 2,092:
}
 
/* BigDec represents an arbitrary length decimal number. Internally, it is simply a BigRat object
// Type aliases for classes in case of any name clashes with other modules.
but with the number of decimal places in its expansion limited to a given 'precision'.
var Big_BigInt = BigInt
The minimum and default precision is 16 decimal places. BigDec objects are immutable.
var Big_BigInts = BigInts
*/
var Big_BigRat = BigRat
class BigDec is Comparable {
var Big_BigRats = BigRats
// Private helper function to check that 'o' is a suitable type and throw an error otherwise.
var Big_Comparable = Comparable // in case imported indirectly
// Numbers, rational numbers and numeric strings are returned as BigDecs.
static check_(o) {
if (o is BigDec) return o
if (o is Num || o is String) return new(o)
if (o.type.toString == "Rat") return new_(BigRat.fromRat(o))
if (o is BigRat) return fromBigRat(o)
if (o is BigInt) return fromBigInt(o)
Fiber.abort("Argument must either be a number, a rational number or a numeric string.")
}
 
// Numeric constants (all have default precision)
static minusOne { fromInt(-1) }
static zero { fromInt(0) }
static one { fromInt(1) }
static two { fromInt(2) }
static three { fromInt(3) }
static four { fromInt(4) }
static five { fromInt(5) }
static six { fromInt(6) }
static seven { fromInt(7) }
static eight { fromInt(8) }
static nine { fromInt(9) }
static ten { fromInt(10) }
static half { new(0.5) }
static quarter { new(0.25) }
static fifth { new(0.2) }
static eighth { new(0.125) }
static tenth { new(0.1) }
 
// Transcendental constants (all have precision 64)
static e { new("2.7182818284590452353602874713526624977572470936999595749669676277", 64) }
static ln2 { new("0.6931471805599453094172321214581765680755001343602552541206800095", 64) }
static ln10 { new("2.3025850929940456840179914546843642076011014886287729760333279010", 64) }
static pi { new("3.1415926535897932384626433832795028841971693993751058209749445923", 64) }
static tau { new("6.2831853071795864769252867665590057683943387987502116419498891846", 64) }
static phi { new("1.6180339887498948482045868343656381177203091798057628621354486227", 64) }
static sqrt2 { new("1.4142135623730950488016887242096980785696718753769480731766797380", 64) }
 
// Returns 'pi' to a given precision using the Almqvest Berndt AGM method.
// See https://rosettacode.org/wiki/Arithmetic-geometric_mean/Calculate_Pi for details.
static pi(prec) {
if (prec < 16) prec = 16
var an = BigRat.one
var bn = BigRat.half.sqrt(prec)
var tn = BigRat.half.square
var pn = BigRat.one
while (pn <= prec) {
var prevAn = an
an = (bn + an) * BigRat.half
bn = (bn * prevAn).sqrt(prec)
prevAn = prevAn - an
tn = tn - (prevAn.square * pn)
pn = pn + pn
}
var res = (an + bn).square / (tn * 4)
return BigDec.fromBigRat(res, prec)
}
 
// Returns the greater of two BigDec objects.
static max(d1, d2) { (d1 < d2) ? d2 : d1 }
 
// Returns the smaller of two BigDec objects.
static min(d1, d2) { (d1 < d2) ? d1 : d2 }
 
// Constructs a BigDec object from a number or decimal numeric string.
construct new(s, prec) {
_prec = prec.max(16)
_br = BigRat.fromDecimal(s)
if (s is String) _br = _br.round(_prec)
}
 
// Constructs a BigDec object from an integer.
construct fromInt(i, prec) {
if (!(i is Num && i.isInteger)) Fiber.abort ("Argument must be an integer.")
_br = BigRat.new(i, 1)
_prec = prec.max(16)
}
 
// Constructs a BigDec object from a big integer.
construct fromBigInt(bi, prec) {
if (!(bi is BigInt)) Fiber.abort ("Argument must be a big integer.")
_br = BigRat.new(bi, BigInt.one)
_prec = prec.max(16)
}
 
// Constructs a BigDec object from a rational number.
construct fromRat(r, prec) {
_br = BigRat.fromRat(r)
_prec = prec.max(16)
}
 
// Constructs a BigDec object from a big rational number.
construct fromBigRat(br, prec) {
if (!(br is BigRat)) Fiber.abort ("Argument must be a big rational number.")
_br = br.copy()
_prec = prec.max(16)
}
 
// Convenience versions of the above constructors which always use a precision of 16.
static new(s) { new(s, 16) }
static fromInt(i) { fromInt(i, 16) }
static fromBigInt(bi) { fromBigInt(bi, 16) }
static fromRat(r) { fromRat(r, 16) }
static fromBigRat(br) { fromBigRat(br, 16) }
 
// Private constructor.
construct new_(br, prec) {
_br = br.round(prec)
_prec = prec
}
// Private method to return the internal BigRat without copying.
br_ { _br }
 
// Get or sets the precision for this object.
prec { _prec }
prec=(p) { _prec = p.floor.max(16) }
 
// Other properties.
isInteger { _br.isInteger } // checks if integral or not
isPositive { _br.isPositive } // checks if positive
isNegative { _br.isNegative } // checks if negative
isUnit { _br.isUnit } // checks if plus or minus one
isZero { _br.isZero } // checks if zero
 
/* Note that the result of unary operations has the same precision as the operand
and the result of binary operations has the greater precision of the operands */
 
// Rounding methods (similar to those in Num class).
ceil { BigDec.new_(_br.ceil, _prec) } // higher integer, towards zero
floor { BigDec.new_(_br.floor, _prec) } // lower integer
truncate { BigDec.new_(_br.truncate, _prec) } // lower integer, towards zero
round { BigDec.new_(_br.round, _prec) } // nearer integer
roundUp { BigDec.new_(_br.roundUp, _prec) } // higher integer, away from zero
 
round(digits) { BigDec.new(_br.round(digits), _prec) } // rounds to 'digits' decimal places
fraction { this - truncate } // fractional part (same sign as this)
 
// Negation.
- { BigDec.new (-_br, _prec) }
 
// Arithmetic operators (work with numbers, rational numbers and numeric strings
// as well as other decimal numbers).
+(o) { (o = BigDec.check_(o)) && BigDec.new_(_br + o.br_, _prec.max(o.prec)) }
-(o) { (o = BigDec.check_(o)) && BigDec.new_(_br - o.br_, _prec.max(o.prec)) }
*(o) { (o = BigDec.check_(o)) && BigDec.new_(_br * o.br_, _prec.max(o.prec)) }
/(o) { (o = BigDec.check_(o)) && BigDec.new_(_br / o.br_, _prec.max(o.prec)) }
%(o) { (o = BigDec.check_(o)) && BigDec.new_(_br % o.br_, _prec.max(o.prec)) }
 
// Integral power methods.
pow(i) { BigDec.new_(_br.pow(i), _prec) }
square { BigDec.new_(_br.square, _prec) }
cube { BigDec.new_(_br.cube, _prec) }
 
// Roots to 'digits' places or 16 places by default.
sqrt(digits) { BigDec.new_(_br.sqrt(digits), _prec) }
cbrt(digits) { BigDec.new_(_br.cbrt(digits), _prec) }
sqrt { sqrt(16) }
cbrt { cbrt(16) }
 
// Other methods.
inverse { BigDec.new_(_br.inverse, _prec) } // inverse
idiv(o) { (this/o).truncate } // integer division
inc { this + BigDec.one } // increment
dec { this - BigDec.one } // decrement
abs { BigDec.new_(_br.abs, _prec)} // absolute value
sign { _br.num.sign } // sign
 
// The inherited 'clone' method just returns 'this' as BigDec objects are immutable.
// If you need an actual copy use this method instead.
copy() { BigDec.new_(_br.copy(), _prec) }
 
// Compares this BigDec with 'other' to enable comparison operators via Comparable trait.
compare(other) { (other is BigDec) ? _br.compare(other.br_) : _br.compare(other) }
 
// As above but compares absolute values.
compareAbs(other) { (other is BigDec) ? _br.compareAbs(other.br_) : _br.compareAbs(other) }
 
// Returns this BigDec expressed as a BigInt with any fractional part truncated.
toBigInt { _br.num/_br.den }
 
// Returns this BigDec expressed as a BigRat.
toBigRat { _br.copy() }
 
// Converts the current instance to a Num where possible.
// Will probably lose accuracy for larger numbers.
toFloat { _br.toFloat }
 
// Converts the current instance to an integer where possible with any fractional part truncated.
// Will probably lose accuracy for integers >= 2^53.
toInt { this.toFloat.truncate }
 
// Returns the string representation of this BigDec object to 'digits' decimal places.
// If 'rounded' is true, the value is rounded to that number of places with halves
// being rounded away from zero. Otherwise the value is truncated to that number of places.
// If 'zfill' is true, any unfilled decimal places are filled with zeros.
toString(digits, rounded, zfill) { _br.toDecimal(digits, rounded, zfill) }
 
// Convenience versions of the above which use default values for some or all parameters
// and never add trailing zeros.
toString(digits, rounded) { toString(digits, rounded, false) }
toString(digits) { toString(digits, true, false) } // always rounded
toString { toString(_prec, true, false) } // precision digits, always rounded
toDefaultString { toString(16, true, false) } // 16 digits, always rounded
 
// Private worker method for exponentiation using Taylor series.
exp_ {
var y = _br
var x = _br
var z = _br
var res = x + BigRat.one
var eps = BigRat.new(BigInt.one, BigInt.ten.pow(_prec))
var k = 2
var fact = BigRat.two
while (x.abs > eps) {
z = z * y
x = z / fact
res = res + x.round(_prec+2)
k = k + 1
fact = fact * k
}
return BigDec.fromBigRat(res, _prec)
}
 
// Returns the exponential of this instance (e ^ this) to the same precision.
exp {
if (_br.isInteger) return exp_
return (truncate/2).exp_.square * fraction.exp_
}
 
// Private worker method for natural logarithm using Taylor series.
log_ {
var y = _br
var x = (y - BigRat.one) / (y + BigRat.one)
var z = x * x
var res = BigRat.zero
var eps = BigRat.new(BigInt.one, BigInt.ten.pow(_prec))
var k = 1
while (x.abs > eps) {
res = res + (BigRat.two * x / k).round(_prec+2)
x = (x * z).round(_prec+2)
k = k + 2
}
return BigDec.new_(res, _prec)
}
 
// Returns the natural logarithm of this instance to the same precision.
log {
if (sign <= 0) Fiber.abort("Instance must be positive.")
if (this < BigDec.two) return log_
var d = this.copy()
var pow = 0
var two = BigDec.new(2, _prec)
while (d >= two) {
d = d / 2
pow = pow + 1
}
return d.log_ + two.log_ * pow
}
 
// Returns this ^ f where f may be fractional.
powf(f) {
if ((f is Num) && f.isInteger) return pow(f)
return (log * f).exp
}
}
 
/* BigDecs contains various routines applicable to lists of big decimal numbers */
class BigDecs {
static sum(a) { a.reduce(BigDec.zero) { |acc, x| acc + x } }
static mean(a) { sum(a)/a.count }
static prod(a) { a.reduce(BigDec.one) { |acc, x| acc * x } }
static max(a) { a.reduce { |acc, x| (x > acc) ? x : acc } }
static min(a) { a.reduce { |acc, x| (x < acc) ? x : acc } }
}
 
// Initialize static fields.
BigInt.init_()</langsyntaxhighlight>
9,476

edits