# BinarySplitting.frink

``` /** This library implements a general algorithm for performing binary splitting     calculations of arbitrary-precision calculations.  It follows the algorithm     described in:     "Fast multiprecision evaluation of series of rational numbers" by Bruno     Hible and Thomas Papanikolaou:     https://ginac.de/CLN/binsplit.pdf     To make this work for a particular binary splitting algorithm, you will     need to implement the interface ParamProvider below in one of your own     classes.  See BinarySplittingPi.frink for an example of calculating pi     using binary splitting or BinarySplittingExp.frink for an (incomplete)     implementation of the exp[x] function. */ class BinarySplitting {    /** This is usually the method you want to call.  It performs the binary        splitting and the final "fixup" calculation S = T / (B Q) to obtain the        final result.  If you need to manage significant digits better, you        may need to call the split[] method yourself and then do that fixup        calculation on the result.    */    class calc[arg is ParamProvider, n1, n2, digits = getPrecision[]] :=    {       result = split[arg, n1, n2]       origPrec = getPrecision[]       try       {          setPrecision[digits+3]          // This is a final floating-point result          S = result.T / (result.B * result.Q)          setPrecision[digits]          return 1. * S   // Round to the specified number of digits       }       finally          setPrecision[origPrec]    }    /** This performs recursive binary splitting on particular arguments.        You usually don't want to call this directly, but call the calc[]        function.        returns:  a SplitResult object containing the results of the calculation.    */    class split[arg is ParamProvider, n1, n2] :=    {       diff = n2 - n1       if diff == 0          return "split:  Error.  Difference n2 - n1 should be greater than 0."              result = new SplitResult              if diff == 1       {          // The result at the point n1          p = arg.p[n1]          result.P = p          result.Q = arg.q[n1]          result.B = arg.b[n1]          result.T = arg.a[n1] * p          return result       }       /* The paper above *very helpfully* says "cases 2, 3, 4 left out for          simplicity".  Augh.  But these are the hard-to-derive and tricky parts.          I've dug the implementations out of their source code for the CLN          libraries (written in C) but haven't implemented them here yet.          The libraries are available at:          https://ginac.de/CLN/          From their source, you probably want the very-obviously-named file          src/float/transcendental/cl_LF_ratseries_pqab.cc                 and in that, maybe the function eval_pqab_series_aux          to see how more of the recursive cases are defined.       */       // Default recursive case       // find the middle of the interval       nm = (n1 + n2) div 2       // sum left side       L = split[arg, n1, nm]       // sum right side       R = split[arg, nm, n2]       // Put together       result.P = L.P * R.P       result.Q = L.Q * R.Q       result.B = L.B * R.B       result.T = R.B * R.Q * L.T + L.B * L.P * R.T       return result    } } /** This is an interface that you will implement in a class that provides     the values of a[n], b[n], p[n], and q[n].  Since these are usually a     simple function of n (and maybe some other parameters,) it's implemented     as an interface with basic functions rather than by forcing the user to     create a (possibly infinite) array of all the possible values of each     variable for each n.  These functions should all return integers. */ interface ParamProvider {    // You must implement these functions to return values for integer values    // of n, n>=0    a[n]    b[n]    p[n]    q[n] } /** This class is a simple data structure that holds the current values of     P, Q, B, and T.  It may someday be used to *resume* (or somday checkpoint or     parallelize) calculations. */ class SplitResult {    var P    var Q    var B    var T } ```