In fact, I just implemented yet another algorithm, this time by Beckermann and Labahn. It seems that this one is even faster. In fact, a *lot* faster.
If you are interested, here comes a preliminary version.
fffg.spad
Description: Binary data
usage: x: List POLY INT := [i for i in 1..30]; _ y: List POLY INT := [random(10) for i in 1..#x]; _ C:SquareMatrix(#x, POLY INT) := diagonalMatrix x; _ f:Vector SUP POLY INT := [-1, fun]; _ c2:(NNI, Vector SUP POLY INT, Vector SUP POLY INT) -> POLY INT c2(k, f, M) == eval(dot(f,M)(x.(k+1)),fun=y.(k+1)) r := fffg(C,c2,f,append([1 for i in 1..14],[2 for i in 15..#x])); computes a certain list of matrix, the interpolant of degree 14/15 then is r.1.(1,1)/r.1.(2,1) Another advantage is that this algorithm also computes other recurrence relations as needed in my guessing package :-) The corresponding article is available online: I implemented only the Algorithm from Table 1, the more general algorithm from Table 2 is not yet done. Furthermore, I did not optimize in any way. In particular, the representation of the elements of the vector space, the representation of the special multiplication rule and the computation of c_k are handled in a very unsatisfactory way. I'm very surprised that the algorithm is even in this setting so fast. Well, maybe I'm missing something. (Or maybe there is an efficiency bug in the implementation of Egecioglu-Koc-Gemignanis algorithm...) "Fraction-free Computation of Matrix Rational Interpolants and Matrix GCDs" http://math.univ-lille1.fr/~bbecker//ano/pub/1997/ano375.ps Other publications: http://math.univ-lille1.fr/~bbecker//bb_pub.html The following might also be useful: http://math.univ-lille1.fr/~bbecker//ano/pub/2002/ano445.ps If you feel like programming, how about implementing a domain for holonomic (= d-finite = P-recursive) functions? You might want to read http://www.risc.uni-linz.ac.at/research/combinat/software/GeneratingFunctions/pub/mallinger96.pdf to this end. Martin
_______________________________________________ Axiom-math mailing list [email protected] http://lists.nongnu.org/mailman/listinfo/axiom-math
