MA522 Final Project (due Tue, Dec 15, 5pm).

----------------------------------------------------------------------------------------------
Your final project is to implement in Maple a complete factorization algorithm
for polynomials modulo p.

myfactor := proc(f   # f is a univariate polynomial with integer coefficients
                ,x   # the variable in f
                ,p   # p is prime number (you need not test that this input is actually prime)
                )
# returns a list
# 
# [c, [f1, e1],..., [fr,er]]
# 
# where c is an integer, e1,...,er are positive integers
# and f1,...,fr are monic non-constant polynomials in x
# such that
# 
# f congruent c * f1^e1 * ... * fr^er module p
# 
# and f1,...,fr are irreducible and pairwise relatively prime modulo p
# (the latter is to ensure that no two polynomials are congruent modulo p,
# as would be if f1 = x^2+2 and f2 = x^2-1 and p = 3.)

... # your code

end; #myfactor

This is essentially the Maple function

   Factor(f) mod p; 

Your procedure first computes c and a monic modular image of f,
then a squarefree factorization of f mod p via Algorithm 14.21 (GG)
with help of Exercise 14.27 (GG) and finally factors the squarefree
factors via Algorithm 14.3 and 14.8 (GG---Cantor/Zassenhaus).

Your implementation should be moderately efficient, e.g., should
use repeated squaring.

You may use polynomial and matrix operations provided for basic
problems in Z_p, for example, Gcd, Powmod, etc.
but of course not Factor or DistDeg.

You can test your program by multiplying polynomials together (Expand(f) mod p).
Also try to factor  x^4 + 1 mod ithprime(i)  for i = 1,2,...
and  numtheory[cyclotomic](7^i,x) mod 2  for i=1,2,3.

----------------------------------------------------------------------------------------------

Problem 2:
As an application of your algorithm, you are to factor primes, namely

myGIfactor = proc(p # p is a prime integer congruent 1 modulo 4
                    # (you need not verify primality here)
                 )
# returns the Gaussian integer a+I*b such that (a+I*b)*(a-I*b) = p
# E.g., myGIfactor(5) returns 2+I.
end;

This is the essential part of the Maple function GaussInt[GIfactor].

The algorithm is based on factoring x^2 + 1 mod p to computing
a residue c mod p such that c^2 mod p = -1.
Then the GCD(p,c+I) over the Gaussian integers (cf. Exercise 3.19 [GG p. 59])
is a+I*b.  Note the Maple function GaussInt[GIgcd], which you may use.

Please proof that your algorithm will always work.

----------------------------------------------------------------------------------------------

Good luck!
Yours, Erich