M1 Cryptographie - Arithmétique modulaire

Arithmétique élémentaire

L'arithmétique est l'étude des nombres entiers $${\mathbb N}=\{0,1,2,3,4,\ldots\}$$ $${\mathbb Z}=\{\ldots,-3,-2,-1,0,1,2,3,\ldots\}$$

Division euclidienne

La propriété fondamentale est la division euclidienne : Étant donnés deux entiers $a\in{\mathbb Z}$ et $b>0$, il existe un unique couple d'entiers $(q,r)$ (quotient et reste) tel que $$a =bq+r\ {\rm avec}\ 0\le r<b.$$

Lorsque le reste $r$ est nul, on dit que $b$ divise $a$ ou que $a$ est multiple de $b$, et on écrit $b|a$.

Un entier $>1$ qui n'a pas de diviseurs non-triviaux, c'est à dire autres que 1 et lui-même est dit premier. Par exemple, $2,3,5,7,11,13,17,19,23$ sont premiers.

Le p.g.c.d. de deux entiers $a$ et $b$ est le plus grand entier positif qui divise à la fois $a$ et $b$. On le note $d=a\wedge b$. Par exemple $(-12)\wedge 18=6$.

Le plus petit entier positif divisible à la fois par $a$ et $b$ est appelé p.p.c.m. de $a$ et $b$. On le note $m=a\vee b$. On a $$a\vee b=\frac{ab}{a\wedge b}$$

Algorithme d'Euclide

Pour calculer $a\wedge b$, (avec $a\ge b$) on écrit la division euclidienne $a=bq+r$. Alors $r=a-bq$ est divisible par tout entier divisant à la fois $a$ et $b$ donc par le p.g.c.d. de $a$ et $b$, et tout entier divisant $r$ et $b$ divise aussi $a$. Donc, $$a\wedge b = b\wedge r$$ Comme $r<b$, $(b,r)<(a,b)$ et on itère le procédé jusqu'à arriver à un reste nul. Alors, $d\wedge 0=d$ et on afini.

In [1]:
def gcd(a,b):
    if b==0: return a
    return gcd(b,a%b)

gcd(192,36)    
Out[1]:
12
In [2]:
gcd(36,192)
Out[2]:
12
In [3]:
gcd(-36,192)
Out[3]:
12
In [4]:
gcd(192,-36)
Out[4]:
-12

On voit qu'il faut faire attention aux signes. On évitera aussi la récursivité, car on va travailler sur des entiers très longs.

In [5]:
def gcd(a,b):
    a,b = abs(a), abs(b)
    r = a%b
    while r:
        a,b=b,r
        r = a%b
    return b

Algorithme d'Euclide étendu et théorème de Bézout

Dans la version itérative, on calcule donc une suite de quotients et de restes qu'il est commode de numéroter ainsi. On pose $a_{-1}=a$, $a_0=b$, $q_0=q$, $a_1=r$, et les itérations calculent $$a_{-1}=a_0q_0+a_1$$ $$a_0=a_1q_1+a_2$$ $$\cdots=\cdots$$ $$a_{k-1}=a_kq_k+a_{k+1}$$ et si le premier reste nul est $a_{k+1}$, alors le p.g.c.d. est $d=a_k$.

Ceci implique qu'à chaque étape, on peut calculer des entiers $u_i,v_i$ tels que $a_i=u_ia+v_ib$ En effet, au départ, $a_1=1a+0b$, $a_0=0a+1b$, et par récurrence, comme $a_{i-2}=a_{i-1}q_{i-1}+a_i$, on peut écrire $$a_i = (u_{i-2}a+v_{i-2}b) - q_{i-1}(u_{i-1}a+v_{i-1}b) = (u_{i-2}-q_{i-1}u_{i-1})a+(v_{i-2}-q_{i-1}v_{i-1})b$$ de sorte que $$u_i = u_{i-2}-q_{i-1}u_{i-1}$$ $$v_i = v_{i-2}-q_{i-1}v_{i-1}$$ Les deux suites $u_i$ et $v_i$ vérifient la même récurrence d'ordre 2, avec des conditions initiales différentes $${u_{-1}\choose v_{-1}}={1\choose 0}$$ $${u_0\choose v_0}={0\choose 1}$$

Si le p.g.c.d. est $a_k$, on aura donc $d=u_ka+v_kb$.

Cette propriété, connue sous le nom de théorème de Bézout, sera très importante dans la suite :

Si $a\wedge b=d$, il existe deux entiers $u$ et $v$ tels que $d=ua+vb$, et l'agorithme ci-dessus, appelé algorithme d'Euclide étendu, permet de les calculer.

Voir ici pour des explications complémentaires.

On programmera cet algorithme en TD. Voilà à quoi devrait ressembler le résultat (le fichier ent.py sera fourni au TD suivant).

In [6]:
from ent import *
In [7]:
xgcd(36,196)
Out[7]:
(4, 11, -2)
In [8]:
11*36-2*196
Out[8]:
4
In [9]:
xgcd(12345678910111213141516171819202122232425,3232323232323232323232)
Out[9]:
(1L, -1328660027339171737607L, 5074743117952096996736521329039888644043L)

Congruences

On dit que $a$ est congru à $b$ modulo $n$ si $a$ et $b$ ont même reste dans la division par $n$, ou encore si $a-b$ est divisible par $n$.

On écrit $a\equiv b \mod n$ ou $a\equiv b [n]$.

En python, le test est a%n == b%n.

La congruence modulo $n$ est une relation d'équivalence. On peut manipuler les congruences comme des égalités (d'où la notation), car $a\equiv b [n]$ et $a'\equiv b' [n]$ entraîne $a+a'\equiv b+b' [n]$, $aa'\equiv bb' [n]$, etc.

D'où l'idée d'introduire des symboles $\bar a$ ($a\in{\mathbb Z}$) qui vérifient $$\bar a =  \bar b\ \Leftrightarrow \ a\equiv b [n]$$

L'addition et la multiplication modulo $n$ possèdent clairement les mêmes propriétés d'associativité et de ditributivité que leurs versions usuelles. On dit que l'ensemble $${\mathbb Z}/n{\mathbb Z}=\{\bar 0,\bar 1, \ldots,\overline{n-1}\},$$ muni de ces opérations, est un anneau commutatif unitaire.

Voici la table d'addition modulo 12, familière grâce aux horloges à cadran.

In [10]:
def tabadd(n):
    h = ' + ' + ''.join(["%3d"%i for i in range(n)])+'\n'
    for i in range(n):     
        l = ('%2d '% i ) + ''.join(["%3d"% ((i+j)%n) for j in range(n)])+'\n'
        h+=l
    return h

print tabadd(12)
 +   0  1  2  3  4  5  6  7  8  9 10 11
 0   0  1  2  3  4  5  6  7  8  9 10 11
 1   1  2  3  4  5  6  7  8  9 10 11  0
 2   2  3  4  5  6  7  8  9 10 11  0  1
 3   3  4  5  6  7  8  9 10 11  0  1  2
 4   4  5  6  7  8  9 10 11  0  1  2  3
 5   5  6  7  8  9 10 11  0  1  2  3  4
 6   6  7  8  9 10 11  0  1  2  3  4  5
 7   7  8  9 10 11  0  1  2  3  4  5  6
 8   8  9 10 11  0  1  2  3  4  5  6  7
 9   9 10 11  0  1  2  3  4  5  6  7  8
10  10 11  0  1  2  3  4  5  6  7  8  9
11  11  0  1  2  3  4  5  6  7  8  9 10


La table de multiplication est peut-être moins familière. Elle possède quelques propriétés étranges.

In [11]:
def tabmul(n):
    h = ' * ' + ''.join(["%3d"%i for i in range(1,n)])+'\n'
    for i in range(1,n):     
        l = ('%2d '% i ) + ''.join(["%3d"% ((i*j)%n) for j in range(1,n)])+'\n'
        h+=l
    return h

print tabmul(12) 
 *   1  2  3  4  5  6  7  8  9 10 11
 1   1  2  3  4  5  6  7  8  9 10 11
 2   2  4  6  8 10  0  2  4  6  8 10
 3   3  6  9  0  3  6  9  0  3  6  9
 4   4  8  0  4  8  0  4  8  0  4  8
 5   5 10  3  8  1  6 11  4  9  2  7
 6   6  0  6  0  6  0  6  0  6  0  6
 7   7  2  9  4 11  6  1  8  3 10  5
 8   8  4  0  8  4  0  8  4  0  8  4
 9   9  6  3  0  9  6  3  0  9  6  3
10  10  8  6  4  2  0 10  8  6  4  2
11  11 10  9  8  7  6  5  4  3  2  1


On constate qu'à la différence de la table d'addition, toutes les lignes ou colonnes ne sont pas des permutations les unes des autres. Dans la table d'addition, chaque ligne contient l'élément neutre 0, $({\mathbb Z}/n{\mathbb Z},+)$ est un groupe.

Par contre, modulo 12, certains éléments n'ont pas d'inverse, et sont des diviseurs de zéro. Par exemple, 2 n'a pas d'inverse, et $\bar 2\times \bar 6 = \bar 0$.

En fait, un diviseur de 0 ne peut jamais être inversible : si on suppose $ab=0$ et $a'a=1$, alors $a'ab=(a'a)b=(1)b=b$ par associativité, et d'un autre côté, $a'ab = a'(ab)=a'0=0$. Donc $b=0$, et $a$ n'est pas un diviseur de 0.

Un corps est un anneau dans lequel tout élément non nul est inversible. Par exemple, ${\mathbb R}$ est un corps, mais ${\mathbb Z}$ n'en est pas un, car par exemple, 2 n'est pas inversible ($\frac12$ n'est pas un entier).

Un anneau unitaire fini sans diviseurs de 0 est toujours un corps : si $a_1=1,\ldots, a_n$ sont ses éléments non nuls, multiplions les tous par l'un d'entre eux $a_i$. Alors, les $n$ éléments $a_1a_i,\ldots,a_na_i$ sont tous différents, car si on avait $a_ja_i=a_ka_i$, on aurait aussi $(a_j-a_k)a_i=0$ et donc des divseurs de zéro. Donc, l'un de ces éléments est égal à 1, et $a_i$ possède donc un inverse.

Dans ${\mathbb Z}/n{\mathbb Z}$, un élément $\bar a$ est inversible si et seulement si $a\wedge n=1$. En effet, si $a$ et $n$ ont un diviseur commun non trivial $d$, $a=da'$, $n=dn'$, alors $\bar a\bar n'=\bar 0$, bien que $\bar a$ et $\bar n'$ soient non nuls.

Inversement, si $a$ est premier avec $n$, $d=a\wedge n=1$, et le théorème de Bézout nous fournit $u$ et $v$ tels que $1=au+vn$. Donc, dans ${\mathbb Z}/n{\mathbb Z}$, $\bar 1 = \bar a\bar u$, et $\bar a$ est inversible.

Le théorème de Bézout fournit donc aussi un algorithme pour calculer les inverses.

In [12]:
def inversemod(x,n):
    d,u,v = xgcd(x,n)
    if d!=1: raise Exception('Element non inversible')
    return u if u>0 else n+u 
In [13]:
inversemod(5,12)
Out[13]:
5
In [14]:
inversemod(6,12)
---------------------------------------------------------------------------
Exception                                 Traceback (most recent call last)
<ipython-input-14-ec6c8b5f9548> in <module>()
----> 1 inversemod(6,12)

<ipython-input-12-bc211e74865f> in inversemod(x, n)
      1 def inversemod(x,n):
      2     d,u,v = xgcd(x,n)
----> 3     if d!=1: raise Exception('Element non inversible')
      4     return u if u>0 else n+u

Exception: Element non inversible
In [15]:
p = 2**607-1
p
Out[15]:
531137992816767098689588206552468627329593117727031923199444138200403559860852242739162502265229285668889329486246501015346579337652707239409519978766587351943831270835393219031728127L
In [16]:
inversemod(123456,p)
Out[16]:
513244954713111266948320083892964698633832864870737164187435907164443554632582377543852619821920830842093461149889424828504089515648935779061669784432652696708489171177179763242192120L
In [17]:
(_*123456)%p
Out[17]:
1L

On vient donc de démontrer :

"Lorsque $p$ est un nombre premier, ${\mathbb Z}/p{\mathbb Z}$ est un corps."

Cette propriété est fondamentale. Tout ce qu'on apprend en premier cycle sur les polynômes ou les matrices reste vrai pourvu que les coefficients soient pris dans un corps. L'existence de corps finis tels que les ${\mathbb Z}/p{\mathbb Z}$ permet d'utiliser ces notions en cryptographie ou pour la construction de codes correcteurs d'erreurs.

Le (petit) théorème de Fermat

Il est donc important d'étudier les nombres premiers, et en particulier d'être capable de les reconnaître. L'observation suivante sera fondamentale pour la suite.

Rappelons la formule du binôme de Newton $$(x+y)^n = \sum_{k=0}^n{n\choose k}x^{n-k}y^k,\quad {n\choose k}=\frac{n!}{k!(n-k)!}=\frac{n(n-1)\cdots(n-k+1)}{1\cdot 2\cdots k}$$

Regardons le triangle de Pascal modulo divers entiers.

In [18]:
def binomial(n,k):
    if n<k: return 0
    if k==0: return 1
    return binomial(n-1,k-1)+binomial(n-1,k)


def triangle(n,r):
    for m in range(n+1):
        print "%2d" % m,
        for k in range(m+1): 
            print "%3d" % (binomial(m,k) % r),
        print

# Modulo 2        
triangle(15,2)
 0   1
 1   1   1
 2   1   0   1
 3   1   1   1   1
 4   1   0   0   0   1
 5   1   1   0   0   1   1
 6   1   0   1   0   1   0   1
 7   1   1   1   1   1   1   1   1
 8   1   0   0   0   0   0   0   0   1
 9   1   1   0   0   0   0   0   0   1   1
10   1   0   1   0   0   0   0   0   1   0   1
11   1   1   1   1   0   0   0   0   1   1   1   1
12   1   0   0   0   1   0   0   0   1   0   0   0   1
13   1   1   0   0   1   1   0   0   1   1   0   0   1   1
14   1   0   1   0   1   0   1   0   1   0   1   0   1   0   1
15   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1   1

In [19]:
# Modulo 3
triangle(15,3)
 0   1
 1   1   1
 2   1   2   1
 3   1   0   0   1
 4   1   1   0   1   1
 5   1   2   1   1   2   1
 6   1   0   0   2   0   0   1
 7   1   1   0   2   2   0   1   1
 8   1   2   1   2   1   2   1   2   1
 9   1   0   0   0   0   0   0   0   0   1
10   1   1   0   0   0   0   0   0   0   1   1
11   1   2   1   0   0   0   0   0   0   1   2   1
12   1   0   0   1   0   0   0   0   0   1   0   0   1
13   1   1   0   1   1   0   0   0   0   1   1   0   1   1
14   1   2   1   1   2   1   0   0   0   1   2   1   1   2   1
15   1   0   0   2   0   0   1   0   0   1   0   0   2   0   0   1

In [20]:
# Modulo 4
triangle(15,4)
 0   1
 1   1   1
 2   1   2   1
 3   1   3   3   1
 4   1   0   2   0   1
 5   1   1   2   2   1   1
 6   1   2   3   0   3   2   1
 7   1   3   1   3   3   1   3   1
 8   1   0   0   0   2   0   0   0   1
 9   1   1   0   0   2   2   0   0   1   1
10   1   2   1   0   2   0   2   0   1   2   1
11   1   3   3   1   2   2   2   2   1   3   3   1
12   1   0   2   0   3   0   0   0   3   0   2   0   1
13   1   1   2   2   3   3   0   0   3   3   2   2   1   1
14   1   2   3   0   1   2   3   0   3   2   1   0   3   2   1
15   1   3   1   3   1   3   1   3   3   1   3   1   3   1   3   1

In [21]:
# Modulo 5
triangle(15,5)
 0   1
 1   1   1
 2   1   2   1
 3   1   3   3   1
 4   1   4   1   4   1
 5   1   0   0   0   0   1
 6   1   1   0   0   0   1   1
 7   1   2   1   0   0   1   2   1
 8   1   3   3   1   0   1   3   3   1
 9   1   4   1   4   1   1   4   1   4   1
10   1   0   0   0   0   2   0   0   0   0   1
11   1   1   0   0   0   2   2   0   0   0   1   1
12   1   2   1   0   0   2   4   2   0   0   1   2   1
13   1   3   3   1   0   2   1   1   2   0   1   3   3   1
14   1   4   1   4   1   2   3   2   3   2   1   4   1   4   1
15   1   0   0   0   0   3   0   0   0   0   3   0   0   0   0   1

On observe que lorsque $n$ est premier, apparemment, tous les ${n\choose k}$ sont divisbles par $n$, sauf le premier et le dernier qui valent toujours 1. De fait, la formule avec les factorielles montre que ce nombre, qui est un entier, est le résultat de la simplification d'une fraction dont le numérateur est divisible par $n$, mais pas le dénominateur si on suppose que $n$ est premier. Donc le résultat est divisible par $n$.

Pour $p$ premier, la formule du binôme dans ${\mathbb Z}/p{\mathbb Z}$ devient ainsi $$(\bar x + \bar y)^p = \bar x^p + \bar y^p.$$ On peut itérer : $(\bar x + \bar y+ \bar z)^p = \bar x^p + (\bar y+\bar z)^p$, etc., de sorte qu'en écrivant un entier quelconque $a=1+1+\cdots+1$ ($a$ fois), on a $$\bar a^p = (\bar 1 +\bar 1 + \cdots + \bar 1)^p = \bar 1^p + \bar 1^p + \cdots+\bar 1^p = \bar 1 +\bar 1 + \cdots + \bar 1=\bar a$$

Ainsi, tout élément de ${\mathbb Z}/p{\mathbb Z}$ vérifie $$\bar a^p=\bar a$$ et s'il est inversible, $$\bar a^{p-1}=\bar 1.$$

Cet énoncé constitue le petit théorème de Fermat, qu'on énonce habituellement avec des congruences :

"Si $p$ est premier, tout entier $a$ vérifie $a^p\equiv a\mod p$, et si de plus $a$ est premier avec $p$, alors $a^{p-1}\equiv 1\mod p$."

In [22]:
# Par exemple, 2^607-1 est premier
p
Out[22]:
531137992816767098689588206552468627329593117727031923199444138200403559860852242739162502265229285668889329486246501015346579337652707239409519978766587351943831270835393219031728127L
In [23]:
# On écrira une fonction powermod(a,m,p) qui calcule a**m % p
powermod(42,p-1,p)
Out[23]:
1L
In [24]:
# Celui ci n'est pas premier
n = 2**647-1
n
Out[24]:
583992399055640987986069965529637289586333248927815671114136642291107221402710705472756839848623539171666215625420084135768154204336056063776340648924443416096255318318113913610607896607565283327L
In [25]:
# Ce calcul le prouve
powermod(3,n-1,n)
Out[25]:
21428223903520356654322755377445578958208675410619789103596222568646209820633205931060590298644292314996507952112174748503760913459299420318597663375422240564022629209691880886344754719021729800L

Exponentiation modulaire

La fonction powermod calcule efficacement $a^m\mod n$ par le procédé suivant (voir ici). On écrit (en pensée) l'exposant $m$ en binaire $m=b_k\cdots b_1b_0$. Les bits $b_i$ sont calculés au fur est à mesure. Au départ, le résultat intermédiare res est 1. On calcule $b_0=m\mod 2$. Si $b_0=1$ on multiplie res par $a$ modulo $n$. On divise $m$ par 2, on remplace $a$ par $a^2\mod n$, et on recommence. Ceci calcule $$\bar a^m = (\bar a)^{b_0}(\bar a^2)^{b_1})(\bar a^{(2^2)})^{b_2}\cdots (\bar a^{(2^k)})^{b_k}$$ en remarquant que $(x^{(2^i)})^2 = x^{(2^{i+1})}$ (programme de maths de 4ème ...).

Tests de primalité

Le théorème de Fermat permet de prouver qu'un nombre $n$ n'est pas premier en exhibant un $a$ tel que $a^{n-1}\not\equiv 1\mod n$. C'est le test de Fermat. On peut montrer que chaque fois qu'on obtient 1, la probabilité que $n$ ne soit pas premier est presque divisée par 4. Malheureusement, il existe des nombres non premiers (les nombres de Carmichael) qui passent tous les tests de Fermat $a^n\equiv a \mod n$. Ils sont très rares, mais suffisent à faire capoter l'entreprise.

In [26]:
c = 1590231231043178376951698401
for a in range(2,12): print powermod(3,c-1,c),
1 1 1 1 1 1 1 1 1 1

In [27]:
print factor(c)
[(17, 1), (19, 1), (23, 1), (29, 1), (31, 1), (37, 1), (41, 1), (43, 1), (61, 1), (67, 1), (71, 1), (73, 1), (79, 1), (97, 1), (113, 1), (199L, 1)]

Le test de Miller-Rabin permet de remédier à ce problème. C'est une petite modification du calcul de $a^{n-1}\mod n$ dans laquelle on rajoute un test.

L'idée est que dans un corps, l'équation $x^2-1=0$ ne peut pas avoir plus de deux solutions. Dans les entiers modulo $n$, on a toujours $\bar 1$ et $-\bar 1$, et $x^2-\bar 1=(x-\bar 1)(x+\bar 1)$. Si on en trouve une autre $\bar b$, alors $(\bar b-\bar 1)(\bar b+\bar 1)=\bar b^2-\bar 1=\bar 0$, donc on a des diviseurs de 0 et $n$ n'est pas premier.

Pour cela, on écrit $n-1=2^km$ avec $m$ impair, et on commence par calculer $b=a^m \mod n$. Si $b=1$, on renvoie True (probablement premier). Sinon, on calcule les $b^{(2^i)}\mod n$ et la première fois que le résultat vaut 1, on teste si la valeur précédente était $\equiv -1\mod n$. Si ce n'est pas le cas, on renvoie False (composé). Si on arrive à $b^{(2^k)}$ sans avoir rencontré de racine carrée non triviale de 1, on renvoie True.

In [28]:
miller_rabin(c)
Out[28]:
False

Le théorème des restes chinois

Il permet de résoudre les système de conguences

$$ \left\{ \begin{matrix} x &\equiv& x_1& [m_1] \\ x &\equiv& x_2& [m_2]\\ \vdots & \cdots &\vdots \\ x&\equiv &x_n&[m_n]\end{matrix}\right.$$

lorsque les modules $m_i$ sont deux à deux premiers entre eux.

Pour chaque $i$, il est facile de trouver $a_i$ qui vérifie $a_i\equiv 0\ [m_j]$ pour tous les $j\not=i$. Il suffit de prendre $a_i=\prod_{j\not=i}m_j$. Ce nombre est premier avec $m_i$, donc inversible modulo $m_i$. Calculons les coefficients de Bézout tels que $1=u_ia_i+v_im_i$. Alors, $y_i=u_ia_i$ vérifie $y_i\equiv 1 \ [m_i]$ et $y_i\equiv 0\ [m_j]$ pour $j\not=i$, ce qui entraîne que $$X := y_1x_1+y_2x_2+\cdots y_nx_n$$ est solution du système. Il y en a d'autres : en ajoutant le produit $m=m_1m_2\cdots m_n$ à une solution, on obtient une autre solution, et deux solutions quelconques $x',x''$ diffèrent forcément d'un multiple de de $m$, puisque leur différence est divisible par tous les $m_i$.

La plus petite solution positive du système est donc $$x = y_1x_1+y_2x_2+\cdots y_nx_n \mod m.$$

Voir ici.

Par exemple, avec $(x_1,x_2,x_3,x_4)=(1,2,3,4)$ et (m_1,m2,m_3,m4)=(1,13,17,19)$, on trouve $x= 41823$ :

solve_chinese([1,2,3,4], [11,13,17,19])
41823
In []: