Help:=proc(): print(`OneStep(p,q,r,y),FINDpqr(p,a,y) `): print(`CG(p,a,y),CGwp(p,a,y,PARA)`): end: #OneStep(p,q,r,y): Given a polynomial p and hyperexpon. a(y) #and a variable y, such that (a/p)'/(a/p) is a rational #function, finds, if possible a bigger p and a smaller a OneStep:=proc(p,q,r,y) local g,j,b: if not type(q,polynom) or not type(r, polynom) then ERROR(`not hyper`): fi: #gcd(r(y), q(y)-jr0(y)) for j from 1 to 1000 do g:=gcd(r,q-j*diff(r,y)): if degree(g,y)>0 then RETURN(p*g^j,normal(j*diff(r/g,y)+(q-j*diff(r,y))/g), normal(r/g), y): fi: od: p,q,r,y: end: #FINDpqr(p,a,y): inputs a polynomial p of y #hyperexp. a(y) and a variable y #and outputs the FINAL P,Q,R FINDpqr:=proc(p,a,y) local b,q,r,Paul,Paul1,var: b:=normal(diff(log(a),y)): q:=numer(b): r:=denom(b): Paul:=[p,q,r,y]: Paul1:=[OneStep(op(Paul))]: while Paul1<>Paul do Paul:=Paul1: Paul1:=[OneStep(op(Paul1))]: od: Paul1: end: #CG(p,a,y): performs the CG:=proc(p,a,y) local Paul,P,Q,R,k,L,f,c,i,G8,var,eq,var1: Paul:=FINDpqr(p,a,y): P:=Paul[1]: Q:=Paul[2]: R:=Paul[3]: L:=degree(Q+diff(R,y),y): if degree(R,y)<=L then k:=degree(P,y)-L: else #ignoring k0 for now k:=degree(P,y)-degree(R,y)+1: fi: if k<0 then RETURN(FAIL): fi: f:=add(c[i]*y^i,i=0..k): var:={seq(c[i],i=0..k)}: G8:=expand(P-(Q+diff(R,y))*f-R*diff(f,y)): eq:={coeffs(G8,y)}: var1:=solve(eq,var): if var1=NULL then RETURN(FAIL): fi: f:=subs(var1,f): normal(f*R*a/P): end: #CGwp(p,a,y,PARA): performs the CGwp:=proc(p,a,y,PARA) local Paul,P,Q,R,k,L,f,c,i,G8,var,eq,var1: Paul:=FINDpqr(p,a,y): P:=Paul[1]: Q:=Paul[2]: R:=Paul[3]: L:=degree(Q+diff(R,y),y): if degree(R,y)<=L then k:=degree(P,y)-L: else #ignoring k0 for now k:=degree(P,y)-degree(R,y)+1: fi: if k<0 then RETURN(FAIL): fi: f:=add(c[i]*y^i,i=0..k): var:={seq(c[i],i=0..k)}: G8:=expand(P-(Q+diff(R,y))*f-R*diff(f,y)): eq:={coeffs(G8,y)}: var1:=solve(eq,var union PARA): f:=subs(var1,f): if f=0 then RETURN(FAIL): fi: normal(f*R*a/P): end: