Baby/giant step strategy

I wrote in my previous post that Cantor and Zassenhaus' algorithm of polynomial factorisation over finite fields can be divided into three stages:

• squarefree factorisation
• distinct-degree factorisation
• equal-degree factorisation

(in the standard Cantor-Zassenhaus algorithm the squarefree factorisation step is merged with the algorithm itself)

According to Modern Computer Algebra, "the second stage consumes the bulk of the computing time". In 1998 Erich Kaltofen and Victor Shoup designed a new approach to distinct-degree factorisation which allows to decrease it's cost. This strategy is called "baby/giant step" and exploits the fact from number theory:

Lemma. For nonnegative integers $i$ and $j$, the polynomial $x^{q^i} - x^{q^j} \in F_q [x]$ is divisible by precisely those irreducible polynomials in $F_q [x]$ whose degree divides $i−j$.

The proof is easy: if $i \geq j$ then  $x^{q^i} - x^{q^j} = (x^{q^{i-j}} - x)^{q^j}$. From the Theorem from the previous post we know that the factorisation of $x^{q^k} - x$ consists of all irreducible factors whose degree is a divisor of $k$. 

Now it's possible to formulate the algorithm:

Fast distinct-degree factorisation (DDF)

Input and output

This algorithm takes as input a square-free polynomial $f \in F_q [x]$ of degree $n$. 

The output is $f ^{[1]} , . . . , f^{[n]} \in F_q [x]$ such that for $1 \leq d \leq n$, $f ^{[d]}$ is the product of the monic irreducible factors of f of degree $d$. The algorithm is parameterized by a constant $\beta$, with $0 \leq \beta \leq 1$.

Step 1 (compute baby steps)

Let $l = \lceil n\beta \rceil$ . For $0\leq i\leq l$, compute $h_i = x^{q^i} \mod f$.
Step 2 (compute giant steps) 

Let $m = \lceil \frac{n}{2l} \rceil$ . For $1\leq j\leq m$, compute $H_j = x^{q^{lj}} \mod f$. 

Step 3 (compute interval polynomials) 

For $1\leq j\leq m$, compute $I_j = \prod\limits_{0\leq i\leq l}(H_j - h_i) \mod f$ .

(by Lemma the polynomial $I_j$ is divisible by those irreducible factors of $f$ whose degree divides an integer $k$ with  $(j-1)l < k \leq jl$)

Step 4 (compute coarse DDF) 

In this step, we compute polynomials $F_1, ... , F_m$ where $F_j = f^{[(j-1)l+1]} f^{[(j-1)l+2]} ... f^{[jl]}$.  This is done as follows.
$s\leftarrow f;$
for $j \leftarrow 1$ to $m$ do
$\left\{F_j \leftarrow \gcd(s, I_j); s\leftarrow \frac{s}{F_j}\right\}$
Step 5 (compute fine DDF) 

In this step, we compute the output polynomials $f^{[1]}, ... , f ^{[n]}$. 

First, initialize $f^{[1]}, ... , f^{[n]}$ to 1. Then do the following.

for $j \leftarrow1$ to $m$ do

{

$g\leftarrow F_j ;$

for $i\leftarrow l - 1$ down to 0 do
    $\left\{ f^{[lj-i]} \leftarrow \gcd(g, H_j - h_i); g\leftarrow \frac{g}{f^{[lj-i]}} \right\}$

}

if $s \neq 1$ then $f^{[\deg(s)]}\leftarrow s$

Some technical remarks

1.  Lemma. For any positive integer $r$ if we are given $h = x^{q^r} \mod f \in F_q [x]$, then for any $g \in F_q [x]$, we can compute  $g^{q^r} \mod f$ as $g(h) \mod f$.

Now to perform Step 2 we set $H_1 = h_l$ and then compute each $H_j$ as $H_{j−1} (H_1) \mod f$.

2. Computing $f(g) \mod h$ is a so-called "modular composition" problem. For solving this problem one can use e.g. Horner scheme or Brent and Kung's algorithm (this also requires fast matrix multiplication).

By now I implemented the simplest version of the fast CZ algorithm: it uses Horner scheme in Step 2.

The next step might be to implement matrix multiplication and Brent and Kung's algorithm.

Cantor-Zassenhaus algorithm

Last Saturday I finished the first big part of my project: I ported Cantor-Zassenhaus algorithm with all helper functions from nmod_poly module to fmpz_mod_poly module.

I think it was an important part of my job because it helped me to familiarise myself with the code: I found some important differences between nmod and fmpz modules (I spent a lot of time debugging my stupid errors relates to these differences).

During the last month I got a good impression of the code of FLINT. I like it's clearness: every function has it's own file and detailed description, the function names are easy-to-understand (I mean that one can guess what the function should do according to it's name), pieces of code doing one distinct thing are organised in separate functions (so the complex algorithms don't look ugly).


Now I'm going to describe Cantor-Zassenhaus algorithm (the complete description can be found in Modern Computer Algebra book).

Input and output
The algorithm gets on input a nonconstant polynomial $P$ over a field $Z_q$ ($q$ is an odd prime number). On output it gives all monic irreducible factors of P with their multipliсities.

Short description
In brief the algorithm can be described as follows: in a cycle for $i=1,2,...$ it computes distinct-degree polynomial $g_i$ (a product of all squarefree irreducible factors of $P$ of degree $i$), then it decomposes $g_i$ into equal degree factors $f_{ij}$ (each of degree $i$) and for all $j$ removes the highest possible degree of $f_{ij}$ (say $e_{ij}$) from $P$ simply dividing $P$ by $f_{ij}$. The cycle continues for $P=P'$ ($P'$ is $P$ without $f_{ij}^{e_{ij}}$) while $\deg P' > 1$.

Distinct-degree factorisation
This algorithm is based on Fermat's little theorem, more precisely on the consequence of the following


Theorem. For any $d\geq 1 (x^{q^d}-x) \in F_q[x]$ is the product of all monic irreducible polynomials in $F_q[x]$ whose degree divides $d$.


Now we can easily construct an algorithm: for $d = 1, 2, ...$ compute $h = x^{q^d} - x$ and find $g=\gcd(h,f)$, continue this procedure for $f=\frac{f}{g}$ until $f=1$.


Equal-degree factorisation
It is a recursive algorithm: it calls a probabilistic procedure of finding a factor $g$ of a given degree of a polynomial $P$, then it removes the found factor $g$ from $P$ and calls itself for $P=\frac{P}{f}$.

Equal-degree splitting
It is a probabilistic procedure of finding a factor $g$ of a given degree of a polynomial $P$. It is based on the following 


Lemma. Let $q$ be an odd prime number. Then $a^{(q-1)/2}\in\{1,-1\} \,\,\, \forall a\in F_q^{*}$ ($F_q^{*}$  is a multiplicative group of  $F_q$). 


Let the polynomial $P$ have the degree $n=rd$ and every it's single factor have the degree $d$. Then the probabilistic algorithm is as follows: choose $a\in F_q[x]$ with $\deg a < n$ at random. If $\deg a = 0$ return "failure". Else compute $g_1=\gcd(a,f)$. If $\deg g_1 > 0$ return $g_1$. Else (if $\deg g_1 = 0$) compute $b=a^{(q^d-1)/2}\mathrm{rem} f$ (see Lemma) and $g_2=\gcd(b-1, f)$. If $g_2\neq 1, g_2\neq f$ return $g_2$, else return "failure".


Some technical remarks
To compute powers of q the binary exponentiation method can be applied.
To compute the gcd  euclidean method can be used.


According to my project proposal I'm out of schedule, so I would like to change the order of the next two algorithms and to continue now with "baby/giant step" algorithm (and to implement Berlekamp after it if there is time left).

Adding helper functions, part 1

Now I'm working on the second step of porting code for Cantor-Zassenhaus factorisation. It is necessary to move functions for polynomial arithmetic such as mulmod, remove, pow_trunc, pow_trunc_binexp, powmod_ui_binexp. For now I moved the first two functions.