• 0 Vote(s) - 0 Average
• 1
• 2
• 3
• 4
• 5
 SAGE code for computing flow matrix for exp(z)-1 jaydfox Long Time Fellow    Posts: 440 Threads: 31 Joined: Aug 2007 08/19/2009, 08:39 PM (This post was last modified: 08/19/2009, 11:59 PM by jaydfox.) (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( + [cof[int(k)] 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 « Next Oldest | Next Newest »

 Messages In This Thread SAGE code for computing flow matrix for exp(z)-1 - by jaydfox - 08/19/2009, 08:39 PM RE: SAGE code for computing flow matrix for exp(z)-1 - by jaydfox - 08/19/2009, 09:17 PM RE: SAGE code for computing flow matrix for exp(z)-1 - by jaydfox - 08/19/2009, 11:47 PM RE: SAGE code for computing flow matrix for exp(z)-1 - by bo198214 - 08/21/2009, 12:38 PM RE: SAGE code for computing flow matrix for exp(z)-1 - by jaydfox - 08/21/2009, 05:32 PM

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

Users browsing this thread: 1 Guest(s) 