Dear GAP Forum,

David Joyner wrote:

One example I've been trying to illustrate in my undergraduate
algebra class is the following:

(a) the coefficients of x^i in the power series expansion
of 1/(1+x^3+x^4) (mod 2),
[...]

GAP does all these, but (a) is the hardest. My code
for (a) is enclosed below. It is so slow that it is almost
impractical for classroom use (try computing 45 terms of the

If you look at your code you'll see that key, a and p are never changed,
only n is a real parameter.

This means, that you call your function again and again with the same
parameters -- this takes time.

A common technique from computer science (``dynamic programming'') for
resolving such a problem is to store results once they were obtained and
never (re)compute them again (in effect you're trading time for space, but
the space you need is minimal). This way an (exponential time) recursion
becomes polynomial time.

> QUESTION: Does anyone have any suggestions for speeding up my code?
I append one modified version (very quickly and clumsily written, I hope I
didn't introduce errors) of the code, which uses caching.
It does 1000 terms in under a minute on my computer.

The problem of asking for the same result several times comes up not too
infrequent in calculations. Some other systems, (for example Maple)
therefore extensively cache results. When we designed GAP we debated about
this, but eventually decided that this could be too memory consuming. GAP
thus (essentially) caches only results that are declared as an attribute.
(In a few cases results with respect to a parent group are also stored in an
attribute.)

For illustration of this, I also append a modified version of the polynomial
version of the code that uses attributes for caching. It is substantially
slower than the list-based version, as much more work has to be done in
constructing the derivatives.

I hope this is of some help.

All the best,

Alexander Hulpke

-- Colorado State University, Department of Mathematics,
Weber Building, Fort Collins, CO 80523, USA
email: hulpke@math.colostate.edu, Phone: ++1-970-4914288
http://www.math.colostate.edu/~hulpke

---------------------------------------------------------------------
# modified code, using caching of intermediate results

Dolfsr:=function(key,a,n,p)
local k,cache,lfsr,lfsrcache,offset;
k:=Size(key);
cache:=[];
offset:=1;

   # check whether value known
   lfsrcache:=function(n)
     if not IsBound(cache[n+offset]) then
cache[n+offset]:=lfsr(n);
     fi;
     return cache[n+offset];
   end;

   lfsr:=function(n)
   local i;
   if (n<k and p>0) then 
     for i in [1..k] do
if n=i-1 then return (key[i] mod p); fi;
     od;
   fi;
   if n<k then 
     for i in [1..k] do
if n=i-1 then return key[i]; fi; 
     od;
   fi;
   if n>=k and p>0 then  
  return Sum([1..k],i->a[i]*lfsrcache(n-i)) mod p; 
   fi;
   if n>=k then 
  return Sum([1..k],i->a[i]*lfsrcache(n-i)); 
   fi;
   end;

return lfsrcache(n);
end;

L:=[];
key:=[1,0,0,1];a:=[0,0,1,1];p:=2;
for i in [0..1000] do Add(L,Dolfsr(key,a,i,p)); od;
return L;

---------------------------------------------------------------------
# modified polynomial code

# derivative is declared as operation in the library, not as attribute
# (this -- as well as the missing Value method -- will be added in the nexet
# release.)
MyDerivative:=NewAttribute("MyDerivative",IsUnivariateRationalFunction);

InstallMethod(MyDerivative,"for univariate rational functions",
true,[IsUnivariateRationalFunction],0,
function(r)
#extends Derivative to rational fcns
local f,g;
# this is now done automatically by the method installation
#if not IsRationalFunction(r) then Print("wrong type\n"); return 0; fi;
f:=NumeratorOfRationalFunction(r);
g:=DenominatorOfRationalFunction(r);
return (Derivative(f)*g-f*Derivative(g))/g^2;
end);

higher_derivative:=function(r,k)
## returns k-th der of rat fcn r
if k=0 then return r; fi;
if k=1 then return MyDerivative(r); fi;
if k>1 then return higher_derivative(MyDerivative(r),k-1); fi;
end;

# replacement for `evaluate'
InstallOtherMethod(Value,"univariate rational functions",true,
[IsUnivariateRationalFunction,IsObject],0,
function(r,a)
#extends Value to rational fcns
#returns value of rat fcn r at a
local f,g;
f:=NumeratorOfRationalFunction(r);
g:=DenominatorOfRationalFunction(r);
return Value(f,a)/Value(g,a);
end);

series:=function(r,a,n)
##computes Taylor series of rat fcn r about x=a
local s,f,g,i;
s:=Value(r,a);
for i in [1..n] do
s:=s+Value(higher_derivative(r,i),a)*x^i/Factorial(i);
od;
return s;
end;