$\newcommand\ZZ{\mathbb Z}$
L'objectif de ce TP est d'implanter de l'arithmétique rapide sur les polynômes. Je vous fournis ci-dessous deux classes, qui permettent de représenter le corps $\ZZ/p\ZZ$ (classe ZpZ
) et les éléments de corps (classe ZpZ_elt
). L'idée est d'implanter les algorithmes de manière générique (sans tenir compte du corps de base) et qu'ils marchent aussi bien avec des polynômes à coefficients entiers (classe int
de Python) qu'avec des polynômes à coefficients dans $\ZZ/p\ZZ$ (classe ZpZ_elt
).
Il est bien sûr très intéressant d'aller regarder le code¹, mais vous pouvez vous contenter des explications ci-dessous et de l'exemple qui suit la définition des classes :
K = Zpz(17)
crée le corps $\ZZ/17\ZZ$ et le stocke dans la variable K
;K(12)
crée l'élément $12$ du corps K
;K(12) * K(2)
calcule $12\times 2$ dans K
(donc vaut $7$) ;+
, -
, /
(si on divise par autre chose que $0$...) et **
(puissance) ;K(12).inverse()
calcule l'inverse de $12$ dans le corps ;x
est un élément de $\ZZ/p\ZZ$, x.field()
renvoie $\ZZ/p\ZZ$ lui-même ;K.generator()
renvoie un générateur de K
, c'est-à-dire une racine primitive d'ordre $p-1$ de K
$=\ZZ/p\ZZ$ ;K.primitive_root(n)
renvoie une racine primitive d'ordre $n$ de K
, s'il en existe une (erreur sinon).Note. Je n'ai pas mis de vérification que le $p$ fourni est réellement premier. C'est-à-vous d'y faire attention !
¹ On pourra s'aider de la documentation de Python sur les opérateurs pour comprendre ce qui est fait.
class ZpZ:
def __init__(self, p):
self._p = p
def __repr__(self):
return "Z/{}Z".format(self._p)
def __call__(self, x):
return ZpZ_elt(x, self._p)
def generator(self):
x = 2
e = (self._p-1)//2
while (self(x)**e)._x == 1:
x += 1
return self(x)
def primitive_root(self, n):
if (self._p-1) % n != 0:
raise ValueError("{} n'a pas de racine primitive d'ordre {}".format(self, n))
return self.generator()**((self._p-1)//n)
class ZpZ_elt:
def __init__(self, x, p):
self._p = p
self._x = x % p
def __repr__(self):
return "{} (mod {})".format(self._x, self._p)
def __neg__(self):
return self.__class__((-self._x)%self._p, self._p)
def __add__(self, other):
if self._p != other._p:
raise ValueError("Eléments de corps différents.")
return self.__class__((self._x + other._x) % self._p, self._p)
def __sub__(self, other):
return self + (-other)
def __mul__(self, other):
if self._p != other._p:
raise ValueError("Eléments de corps différents.")
return self.__class__((self._x * other._x) % self._p, self._p)
def inverse(self):
if self._x == 0:
raise ZeroDivisionError("0 n'est pas inversible !")
r0 = self._p
r1 = self._x
u0 = v1 = 1
u1 = v0 = 0
while r1 > 1:
q = r0 // r1
r0, r1 = r1, r0 - q * r1
u0, u1 = u1, u0 - q * u1
v0, v1 = v1, v0 - q * v1
return self.__class__(v1 % self._p, self._p)
def __truediv__(self, other):
return self * other.inverse()
def __pow__(self, n):
r = 1
x = self._x
while n > 0:
if n % 2 == 1:
r = (r*x) % self._p
x = (x*x) % self._p
n //= 2
return self.__class__(r, self._p)
def field(self):
return ZpZ(self._p)
K = ZpZ(17)
K
x = K(12)
y = K(8)
x*y, x+y, x/y, x**3, x.inverse(), x*x.inverse()
x.field()
g = K.generator()
[g**i for i in range(1, 17)]
w = K.primitive_root(4)
[w**i for i in [1,2,3,4]]
w = K.primitive_root(3) # BOOM!
inverse
et __pow__
?Réponses :
On représente un polynôme par la liste de ses coefficients, le premier élément étant le coefficient constant. On n'impose pas que les représentations des polynômes soient normalisées, c'est-à-dire qu'on autorise des coefficients nuls en tête.
Tous les algorithmes de ce TP auront comme dernier paramètre le zéro du corps considéré. La raison est qu'on doit parfois initialiser des variables à $0$ : si on initialise une variable x
avec le $0$ des int
, on ne pourra pas ensuite faire (par exemple) x + K(1)
. Python ne saura pas faire l'addition entre un int
et un ZpZ_elt
. Ainsi, dans les algorithmes, on utilisera le paramètre zero
pour cela.
Je fournis une implantation possible pour l'addition de deux polynômes, pour illustrer le propos. Notez que la ligne 6 risquerait de provoquer une erreur si R
était initialisé par R = [0] * (...)
et qu'on appliquait l'algorithme avec des polynômes à coefficients dans $\ZZ/p\ZZ$.
def addition(P, Q, zero):
R = [zero] * max(len(P), len(Q)) # Polynôme résultat, de degré le max des degrés
for i in range(len(P)): # Parcours des coeffs de P
R[i] = P[i] # Recopiage de P dans R
for i in range(len(Q)):
R[i] += Q[i] # Ajout des coefficients de Q
return R
addition([1,2,3],[4,5,6,7], 0) # Ici, polynômes à coefficients de type int, donc zero est le 0 des int
P = [K(i) for i in [1,2,3]]
Q = [K(i) for i in [4,5,6,7]]
addition(P, Q, K(0)) # Polynômes à coefficients dans K, donc zero est K(0)
# Erreur, si on prend le mauvais 0 :
addition(P, Q, 0)
Compléter l'implantation de l'algorithme de Karatsuba. Pour se faciliter la vie, on se ramène au cas de deux polynômes de même taille, cette taille étant une puissance de $2$.
def Karatsuba(P, Q, zero):
n = (max(len(P), len(Q))-1).bit_length() # n est minimal tel que len(P) et len(Q) <= 2**n
if n == 0: return [P[0] * Q[0]] # cas de base
while len(P) < 2**n: P.append(zero) # on complète P...
while len(Q) < 2**n: Q.append(zero) # ... et Q
## COMPLETER ##
return
Karatsuba([1,2,4],[3,2], 0)[:4] == [3,8,16,8]
K = ZpZ(7)
P = [K(x) for x in [1,2,4]]
Q = [K(x) for x in [3,2]]
Karatsuba(P, Q, K(0))[:4] == [K(3), K(1), K(2), K(1)]
Implanter l'algorithme de FFT, puis la multiplication de polynômes par FFT. L'algorithme de multiplication ne sera utilisé qu'avec des polynômes à coefficients dans $\ZZ/p\ZZ$. On récupère le corps lui-même pour aller chercher une racine primitive de l'unité.
def FFT(P, w, zero): # w : racine primitive de l'unité d'ordre 2**(len(P))
return
def mul_FFT(P, Q, zero):
K = zero.field() # K: corps des coefficients de P et Q
return
Pour les tests, on pourra utiliser des corps dont qui sont assurés d'avoir des racines primitives de l'unité du bon ordre. La liste suivante contient des nombres premiers de la bonne forme.
Primes = [2**4 + 1, 3*2**5 + 1, 7*2**26 + 1, 5*2**55 + 1]
# Faire vos tests ici !
L'objectif est d'implanter la division euclidienne rapide des polynômes. Pour cela, on définit des opérations sur les séries tronquées, à partir des opérations sur les polynômes.
On écrit d'abord une fonction de multiplication de séries tronquées, qui prend (entre autres) une fonction de multiplication de polynômes. La fonction devra donc pouvoir être appelée avec comme valeurs de mul
soit Karatsuba
soit mul_FFT
.
def multiplication(P, Q, N, mul, zero):
"""
P, Q : séries tronquées à l'ordre N (le même pour les deux)
mul : une fonction de multiplication, qui prend 2 polys et zero en paramètres
"""
R = mul(P, Q, zero)
return R[:N]
## EXEMPLE : si on définit un algo de multiplication naive de polynomes, on peut appeler multiplication avec mulnaive
def mulnaive(P, Q, zero):
N = max(len(P), len(Q))
R = [zero] * (2*N-1)
for i in range(len(P)):
for j in range(len(Q)):
R[i+j] += P[i] * Q[j]
return R
multiplication([1,2,3], [4, 5, 6], 3, mulnaive, 0)
Le réciproque d'un polynôme se calcule très simplement :
def reciproque(P):
return list(reversed(P))
Écrire l'algorithme d'inversion de séries :
def inverse(P, N, mul, zero):
"""
Entrées : P (série tronquée), N (ordre de troncation), mul (fct de multiplication polynomiale), zero (comme d'hab)
Sortie : Inverse de P, à précision N
"""
return
K
Écrire l'algorithme de division euclidienne rapide :
def div_eucl(F, G, mul, zero):
"""
Entrées : F et G deux polynômes, mul et zero comme dans inverse
Sorties : Q, R les quotient et reste dans la division euclidienne de F par G
"""
return