/** Program to calculate e to arbitrary precision. This uses binary splitting for speed. You usually use this by calling E.getE[digits] which will return e to the number of digits you specify. This version has been updated to be fastest with Frink: The Next Generation. Due to limitations in Java's BigDecimal class, this is limited to a little over 101 million digits. This differs from the basic algorithm in eBinarySplitting.frink because when you ask for more digits of e than you previously had, this *resumes* the calculation rather than starting from scratch. See: http://numbers.computation.free.fr/Constants/Algorithms/splitting.html http://www.ginac.de/CLN/binsplit.pdf */ class E { class var largestDigits = -1 class var cacheE = undef /** These variables help us resume the state of the binary-splitting algorithm. */ class var largestP = 1 class var largestQ = 1 class var largestK = 1 class var logFactorial = 0. /** This is the main public method to get the value of e to a certain number of digits, calculating it if need be. If e has already been calculated to a sufficient number of digits, this returns it from the cache. */ class getE[digits = getPrecision[]] := { origPrec = getPrecision[] try { setPrecision[digits] if (largestDigits >= digits) return 1. * cacheE else return 1. * calcE[digits] } finally setPrecision[origPrec] } /** This is an internal method that calculates the digits of e if necessary. */ class calcE[digits] := { origPrec = getPrecision[] try { setPrecision[15] // Find number of terms to calculate. // ln[x!] = ln[1] + ln[2] + ... + ln[x] logMax = digits * ln[10] k = largestK while (logFactorial < logMax) { logFactorial = logFactorial + ln[k] k = k + 1 } setPrecision[digits+3] if largestK < k // We may have already calculated enough! { if (largestK == 1) // This is the first iteration { p = P[0,k] q = Q[0,k] } else { // Continuing iterations pmb = P[largestK, k] qmb = Q[largestK, k] p = largestP qmb + pmb q = largestQ qmb } // Store the biggest values back in the cache largestP = p largestQ = q largestK = k cacheE = 1 + (1. * p)/q } setPrecision[digits] return 1. * cacheE } finally setPrecision[origPrec] } /** Private method for binary splitting. */ class P[a,b] := { if (b-a) == 1 return 1 m = (a+b) div 2 r = P[a,m] Q[m,b] + P[m,b] return r } /** Private method for binary splitting. */ class Q[a,b] := { if (b-a) == 1 return b m = (a+b) div 2 return Q[a,m] Q[m,b] } }