Thread Rating:
  • 0 Vote(s) - 0 Average
  • 1
  • 2
  • 3
  • 4
  • 5
SAGE code for computing flow matrix for exp(z)-1
#1
(Update: Look two posts down, where I have updated code that is much faster, probably on par with the Mathematica code that Andrew posted.)


I don't recall what the exact term is for the matrix which this code computes. I call it the FMA, which I believe stands roughly for "flow matrix". You see, "fma" was the name of a variable in some code which Andrew provided two years ago:

http://math.eretrandre.org/tetrationforu...147#pid147

I don't remember what the matrix is called, but I know what it represents. The columns (or rows, depending on transposition) are the coefficients of Lagrange interpolation polynomials in n, and each polynomial interpolates a coefficient of the power series of the nth iterate of the desired function, provided it's centered at a parabolic fixed point. For example, I have been using it to calculate for exp(z)-1.

Well, the code apparently runs much slower in SAGE than it does in Mathematica. Andrew converted the code directly, it seems, so efficiencies in Mathematica are apparently not conserved in going to SAGE. Fear not, this merely means we must find an equivalent form that is efficient in SAGE, or at least relatively so.

If you look at the code, there are three subfunctions that Andrew wrote. The first calculates the coefficients of the integer iterates of a function centered at its parabolic fixed point (can't recall if the fixed point needs to be 0, but the change is trivial if so). The second two functions create the interpolation polynomials.

Well, the first change I made was in the parabolic_idm function. Rather than computing the taylor series of ser(expr), I changed it to compute expr(ser). It turns out that internal optimizations make this latter form much faster, about four times faster for size=100, five times faster for size=150, and generally with the speed increase getting slightly better with increasing size.

Using the matrix approach that Gottfried helped me figure out, I was able to reduce the second and third function calls to a matrix multiplication, which takes only a matter of seconds for matrices this small. These latter two function calls were taking almost as much time as the first, so reducing this to a few seconds helps dramatically.

Without further ago, here is the code, which runs in SAGE 4.0.1 and 4.1, but is probably backwards compatible as well. (Depending on how you convert Andrew's fma list to a matrix, it might be a transposition of the one I calculate below. So don't freak out if you run into this when trying to compare.)

Code:
# Define the custom IDM:
def parabolic_idm(expr, x, size):
    ret = [];
    ser = x;
    for i in xrange(size):
        ser = taylor(expr(ser), x, 0, size);
        cof = ser.coeffs(x);
        ret.append([0] + [cof[int(k)][0] for k in xrange(len(cof))]);
    top = [0, 1] + [0 for i in xrange(size-1)];
    return [top] + ret;

# Create a Lagrange Interpolation polynomial for each
# coefficient of the iteration function. Note that
# multiplying on the left by a row-vector of powers
# of n (where n is the desired iteration) will give
# a row vector which represents the coefficients of
# a power series in z. Alternatively, multiplying
# on the right by a column vector of powers of z
# will give the coefficients of a power series in n,
# which allows direct calculation of the nth iterate,
# centered at z.
def LagrangePolynomialMatrix(pidm, size):
    XI = matrix(ZZ, size+1);
    XI[0,0] = factorial(size);
    for ccol in xrange(1, size+1):
        XI[0, ccol] = (1-ccol)/ccol * XI[0, ccol-1]
    for rrow in xrange(1,size+1):
        for ccol in xrange(rrow, size+1):
            XI[rrow, ccol] =  (XI[rrow-1, ccol-1] + (1-ccol) * XI[rrow, ccol-1])/ccol
    PI = matrix(ZZ, size+1);
    PI[0,0] = 1;
    for rrow in xrange(1, size+1):
        for ccol in xrange(rrow+1):
            PI[rrow, ccol] = PI[rrow-1, ccol-1] - PI[rrow-1, ccol];
    return (XI*PI)*matrix(QQ, pidm)/factorial(size);

# test data: change or delete as needed
size = 10
time mtx = parabolic_idm(exp(x)-1, x, size)
time fma = LagrangePolynomialMatrix(mtx, size)

Note: the code gives a deprecation warning when calculating the Taylor series, but I'm not sure about the best way to avoid it. (I'll post an update if I can work around it without hurting performance.)

PS: If anyone tries this code and it doesn't work, let me know. If anyone can suggest further ways to optimize, I'd like to know, as learning to optimize SAGE code can be helpful for me in general.
~ Jay Daniel Fox
Reply


Messages In This Thread
SAGE code for computing flow matrix for exp(z)-1 - by jaydfox - 08/19/2009, 08:39 PM

Possibly Related Threads...
Thread Author Replies Views Last Post
  Revisting my accelerated slog solution using Abel matrix inversion jaydfox 21 6,123 02/09/2019, 02:25 PM
Last Post: sheldonison
  C++ code for tet, ate and hexp MorgothV8 0 2,831 07/10/2014, 04:24 PM
Last Post: MorgothV8
  "Kneser"/Riemann mapping method code for *complex* bases mike3 2 6,345 08/15/2011, 03:14 PM
Last Post: Gottfried
  An incremental method to compute (Abel) matrix inverses bo198214 3 8,780 07/20/2010, 12:13 PM
Last Post: Gottfried
  Sage Question? rsgerard 1 3,780 05/09/2010, 11:40 AM
Last Post: bo198214
  Single-exp series computation code mike3 0 2,820 04/20/2010, 08:59 PM
Last Post: mike3
  computing teh last digits without computing the number deepinlife 3 5,477 02/24/2009, 09:09 AM
Last Post: deepinlife
  How to free memory in SAGE? jaydfox 2 5,213 12/21/2007, 06:38 PM
Last Post: andydude
  Convergence of matrix solution for base e jaydfox 6 8,373 12/18/2007, 12:14 AM
Last Post: jaydfox
  Computing Abel function at a given center jaydfox 10 11,463 11/30/2007, 06:44 PM
Last Post: andydude



Users browsing this thread: 1 Guest(s)