/* -*- gp-script -*- */ \\% Basic manipulations of vectors. \\ Author: Ralf Stephan, (reformatted by Joerg Arndt) \\ License: GPL version 3 or later \\ online at http://www.jjj.de/pari/ \\ version: 2014-October-16 (18:31) ggf(v)= { /* Guess generating function for supplied vector of integers */ /* Example: ggf(vector(10,j,fibonacci(j)^2)) ==> (-x + 1)/(x^3 - 2*x^2 - 2*x + 1) */ my(l,m,p,q,B); l = length(v); B = floor(l/2); if ( B<4, return(0) ); m = matrix(B, B, x, y, v[x-y+B+1]); q = qflll(m, 4)[1]; if ( length(q)==0, return(0) ); p = sum(k=1, B, 'x^(k-1)*q[k,1]); q = Pol( Pol( vector(l, n, v[l-n+1]) )*p + O('x^(B+1)) ); if ( polcoeff(p,0)<0, q=-q; p=-p); q = q/p; p = Ser(q+O('x^(l+1))); for (m=1, l, if( polcoeff(p,m-1)!=v[m], return(0)) ); return(q); } /* ----- */ \\ aux: strsgn(a)=if(a>=0, Str("+" a), Str(a)); \\prtpol(c)= \\{ /* print coefficients with ascending degrees */ \\ print1("("); \\ for (d=(poldegree(c)-...) \\ print1(")"); \\} /* ----- */ \\ \\prtgf(v)= \\{ /* Print Generating function */ \\ my(gf, t); \\ gf = ggf(v); \\ if ( gf==0, return(0) ); \\ t = numerator(gf); \\ print1 \\ return(1); \\} /* ----- */ prtrec(v)= { /* Print recurrence relation */ my(c,d,c0,t); c=ggf(v); if ( c==0, return(0) ); c = denominator(c); \\ print(c); \\ c = polrecip(c); \\ print(c); d = poldegree(c); c0 = polcoeff(c,0); print1("a(n) ="); for (k=1, d, t=polcoeff(c,k)/c0; if (t, print1(" "); print1( strsgn( -t ), "*a(n-",k,")" ) ); ); print(); return(1); } /* ----- */ \\ ==== end of file ====