Computing Andrew's slog solution

[jaydfox]

Here's my SAGE code preparing the coefficients, for rational and for integer matrices:

For the rational solution, you'll need to divide each coefficient by k factorial (for the kth coefficient). I did that step elsewhere, which is why it isn't shown below. But for comparison purposes, it's probably best to put it in this function. I'd do it myself, but SAGE is still crunching numbers, so I can't test a code change right now.

Code:

# For optimal speed and accuracy, try to use a rational LogB
def Prepare_Slog_Q(LogB, size):
    print "Preparing slog for base e^(%s) with %s terms"%(LogB, size)
    m = matrix(QQ, size, size, [[[((c+1)^r)*factorial(c+1) -
        (c==(r-1))/(LogB^r)] for c in xrange(size)]
        for r in xrange(size)])
    v = vector(QQ, size, [(k==0) for k in xrange(size)])
    s = m.solve_right(v)
    ret = vector(QQ, len(s)+1, [-1] + list(s))
    return ret

Code:

# For optimal speed and accuracy, try to use a rational LogB
def Prepare_Slog(LogB, size):
    print "Preparing slog for base e^(%s) with %s terms"%(LogB, size)
    m = matrix(QQ, size, size, [[[((c+1)^r) -
        factorial(c+1)*(c==(r-1))/(LogB^r)]
        for c in xrange(size)] for r in xrange(size)])
    v = vector(QQ, size, [(k==0) for k in xrange(size)])
    s = m.solve_right(v)
    ret = vector(QQ, len(s)+1, [-1] + list(s))
    return ret