Proof of Generalized Little Theorem of Fermat

Every linear recurrence (with integer coefficients) is effectively 
a pseudo-primality test.  This is simply a slight generalization of 
Fermat's Little Theorem,  n^p = n (mod p).  Consider the general 
linear recurrence of order d

      s[n] = c_1 s[n-1] + c_2 s[n-2] + ... + c_d s[n-d]            (1)

The characteristic polynomial of this recurrence is

     f(x)  =  x^d - c_1 x^(d-1) - c_2 x^(d-2) - ... - c_d          (2)

If r is any root of this polynomial then clearly the sequence

          1   r   r^2   r^3   r^4  ...

satisfies recurrence (1).  The general solution of (1) (assuming
the roots of (2) are distinct) is just a linear combination of 
powers of the roots r_1, r_2,..,r_d of f(x).  In other words,

   s[n]  =  b_1 (r_1)^n  +  b_2 (r_2)^n  +  ...  +  b_d (r_d)^n

where the coefficients b_j are constants determined by the arbitrary 
initial values of the sequence.  For example, if we set b_j = 1 for 
all j, then the value of s[n] is just the nth "Newton Sum" for f(x), 
i.e., the sum of the nth powers of the roots.  In other words,

         s[n]  =  (r_1)^n  + (r_2)^n  +  ...  + (r_d)^n

In this case s[0]=d and s[1]=c_1.  By the multinomial theorem we 
know that if p is a prime then the coefficient of every cross-product
term in the expansion of (x+y+..+z)^p is divisible by p.  (The standard 
"additive" proof of Fermat's Little Theorem uses this same fact.)  
It follows immediately that in the expansion of

                 (r1 + r2 + .. + rd)^p

the coefficient of every term is divisible by p, with the exception
of the terms r1^p, r2^p,.., rd^p.  Thus we have

 (r1 + r2 + .. + rd)^p   =   r1^p + r2^p + ... + rd^p

                                (p-1)!
                    + p  SUM ------------- (r1^a1)(r2^a2)...(rd^ad)
                             a1! a2!...ad!

where the summation is evaluated over all sets of non-negative 
integers {aj} where each aj is less than p and SUM aj = p.  This
sum is a symmetrical sum & product function of the roots of f(x),
so it can be expressed in terms of the coefficients of f, which
implies that the summation is an integer.  As a result, we have

 (r1 + r2 + .. + rd)^p   =   r1^p + r2^p + ... + rd^p    (mod p)

which establishes that s[p]=s[1] (mod p) for every prime p.  In fact,
it proves that s[kp]=s[k] (mod p) for any integer k and prime p.  

Any composite integer c such that s[kc]=s[k] (mod c) for all positive 
integers k is called a symmetric pseudoprime relative to f(x).  (Note 
that if the congruence is satisfied for k=1,2,..,d, then it will be 
satisfied for all k, so these d congruences are a generalization of 
the "signature" referred to in Shanks and Adams 1982 paper.)  The 
larger the Galois group of a polynomial, the more rare are symmetric 
pseudoprimes relative to that polynomial.  The smallest such composite 
relative to the Fibonacci polynomial x^2 - x - 1 is 705.  The smallest 
such composite relative to Perrin's cubic x^3 - x - 1 is 27664033.  
The smallest relative to the quintic  x^5 - x^3 - 2x^2 + 1  is 
2258745004684033.