Pour la théorie, voir ici (entre autres): http://fr.wikipedia.org/wiki/Th%C3%A9or%C3%A8me_des_restes_chinois

Sur le plan historique, ça vient des temps anciens de la Chine, et on cite plusieurs problèmes résolus. J'en prends un au hasard:

Pour vérifier, que les 1000 soldats sont tous là, on les place par rangée de 7, puis par 11 puis par 13. S'il manque un seul soldat pour compléter la dernière rangée dans les 3 cas, alors, il y a bien 1000 soldats.

En arithmétique modulaire, on va représenter le problème comme cela:

x =  6 mod  7  => c'est à dire que 6 est le reste de la division entière par 7
x = 10 mod 11  => c'est à dire que 10 est le reste de la division entière par 11
x = 12 mod 13  => c'est à dire que 12 est le reste de la division entière par 13

Et en appliquant le code de cette page, vous pourrez voir que ça marche: on obtient bien x = 1000!

Sur un plan plus général, il s'agit ici de résoudre (=trouver x et les conditions pour que x existe) des équations modulaires du genre:

x = a1 mod m1
x = a2 mod m2
...
x = ai mod mi
...
x = an mod mn

Rappelons que x = a1 mod m1 (mod=modulo) signifie:

• que x % m1 = a1 % m1 (même reste dans la division entière)
• ou que (x-a1) % m1 = 0
• ou que x = k*m1 + a1 (avec k entier quelconque)
• et si a1<m1: a1 est le reste de la division entière (ou “euclidienne”) de x par m1

Le signe '=' n'en est pas vraiment un. Il faudrait dire ici: “x congru à a1 modulo m1”: ce sont des relations de congruence.

On doit donc ici trouver x (s'il existe) qui satisfait ces 2 équations modulaire:

x = a1 mod m1
x = a2 mod m2

Au lieu de se plonger tout de suite dans des considérations mathématiques, on va étudier un cas pratique: deux engrenages qui engrènent ensemble!

Données de départ:

• 2 roues dentées qui engrènent ensemble
• l'une possède m1 dents, et l'autre m2 dents.
• on numérote leurs dents: la dent en contact au départ est notée '0'
• les autres dents sont notées 1, 2, 3, … sur chaque engrenage, dans leur sens de rotation (ce n'est pas le même!), de sorte qu'en commençant à tourner, la dent n°1 du 1er engrenage engrène avec la dent n°1 du 2ème engrenage.

A partir de la position de départ, on va faire tourner les 2 engrenages dent à dent, et on va compter le nombre de dents qui passent devant le point zéro initial.

On va se poser les questions suivantes:

• combien de dents (x) vais-je passer pour atteindre la position voulue pour les 2 roues, par exemple: la dent n°3 (=a1) de la roue 1 engrène avec la dent n°5 (=a2) de la roue 2?

• quelles sont les conditions pour qu'une telle solution soit possible?

Et bien, plutôt que de raisonner dans le vide, on va simuler la rotation avec une fonction.

def simul(a1, m1, a2, m2):
    c = 0 # compteur des dents qui engrènent au contact des 2 roues
    d1 = 0 # compteur de dents sur la roue 1
    d2 = 0 # compteur de dents sur la roue 2
    while True:
        c += 1
        if c >100000: 
            # test de bouclage infini = pas de solution (renvoyer -1)
            c = -1
            break
        d1 = (d1+1)%m1  # on avance d'une dent sur la roue 1
        d2 = (d2+1)%m2  # idem sur la roue 2
        if d1==a1%m1 and d2==a2%m2:
            # ça y est, on est arrivé à la position voulue au bout de c dents  
            break
    return c

Dans cette fonction, on part d'une position 0 et on fait tourner dent par dent jusqu'à obtenir, si elle est possible,la position voulue, c'est à dire quand la dent n°a1 de la roue 1 engrène avec la dent n°a2 de la roue 2. la fonction renvoie c, le nombre de dents qui ont engrenées depuis le départ jusqu'à la position voulue.

On essaye d'abord avec 2 roues dont les nombres de dents sont premiers entre eux (ils sont même premiers ici):

m1 = 3
m2 = 5

Voilà toutes les positions possibles:

Avec a1= 0 et a2= 0 => on trouve c=  15 dents engrenés: Roue 1:  5 tours + 0 dents et roue 2: 3 tours + 0 dents
Avec a1= 1 et a2= 0 => on trouve c=  10 dents engrenés: Roue 1:  3 tours + 1 dents et roue 2: 2 tours + 0 dents
Avec a1= 2 et a2= 0 => on trouve c=  5 dents engrenés: Roue 1:  1 tours + 2 dents et roue 2: 1 tours + 0 dents
Avec a1= 0 et a2= 1 => on trouve c=  6 dents engrenés: Roue 1:  2 tours + 0 dents et roue 2: 1 tours + 1 dents
Avec a1= 1 et a2= 1 => on trouve c=  1 dents engrenés: Roue 1:  0 tours + 1 dents et roue 2: 0 tours + 1 dents
Avec a1= 2 et a2= 1 => on trouve c=  11 dents engrenés: Roue 1:  3 tours + 2 dents et roue 2: 2 tours + 1 dents
Avec a1= 0 et a2= 2 => on trouve c=  12 dents engrenés: Roue 1:  4 tours + 0 dents et roue 2: 2 tours + 2 dents
Avec a1= 1 et a2= 2 => on trouve c=  7 dents engrenés: Roue 1:  2 tours + 1 dents et roue 2: 1 tours + 2 dents
Avec a1= 2 et a2= 2 => on trouve c=  2 dents engrenés: Roue 1:  0 tours + 2 dents et roue 2: 0 tours + 2 dents
Avec a1= 0 et a2= 3 => on trouve c=  3 dents engrenés: Roue 1:  1 tours + 0 dents et roue 2: 0 tours + 3 dents
Avec a1= 1 et a2= 3 => on trouve c=  13 dents engrenés: Roue 1:  4 tours + 1 dents et roue 2: 2 tours + 3 dents
Avec a1= 2 et a2= 3 => on trouve c=  8 dents engrenés: Roue 1:  2 tours + 2 dents et roue 2: 1 tours + 3 dents
Avec a1= 0 et a2= 4 => on trouve c=  9 dents engrenés: Roue 1:  3 tours + 0 dents et roue 2: 1 tours + 4 dents
Avec a1= 1 et a2= 4 => on trouve c=  4 dents engrenés: Roue 1:  1 tours + 1 dents et roue 2: 0 tours + 4 dents
Avec a1= 2 et a2= 4 => on trouve c=  14 dents engrenés: Roue 1:  4 tours + 2 dents et roue 2: 2 tours + 4 dents

Pour le premier cas: on retrouvé la position initiale au bout de 15 dents engrenés, c'est à dire quand la roue 1 a fait exactement 5 tours et la roue 2 exactement 3 tours.

En fait, ces 15 dents correspondent au PPCM de m1 et m2. Cela veut dire:

• que quand la solution existe, sa 1ère valeur se trouve entre 0 et ppcm(m1, m2).

• ou que la solution sera définie à k*ppcm(m1,m2) près, avec k entier quelconque.

Si on cherche toutes les solutions correspondant à toutes les positions possibles de a1 et a2, on les trouve toutes, grâce au fait que m1 et m2 sont premiers entre eux.

Prenons maintenant 2 roues ayant des nombres de dents non premiers entre eux, c'est à dire ayant au moins un des facteurs en commun.

m1 = 4  # 2*2
m2 = 6  # 2*3

Voilà toutes les positions possibles:

Avec a1= 0 et a2= 0 => on trouve c=  12 dents engrenés: Roue 1:  3 tours + 0 dents et roue 2: 2 tours + 0 dents
Avec a1= 1 et a2= 0 pas de solution
Avec a1= 2 et a2= 0 => on trouve c=  6 dents engrenés: Roue 1:  1 tours + 2 dents et roue 2: 1 tours + 0 dents
Avec a1= 3 et a2= 0 pas de solution
Avec a1= 0 et a2= 1 pas de solution
Avec a1= 1 et a2= 1 => on trouve c=  1 dents engrenés: Roue 1:  0 tours + 1 dents et roue 2: 0 tours + 1 dents
Avec a1= 2 et a2= 1 pas de solution
Avec a1= 3 et a2= 1 => on trouve c=  7 dents engrenés: Roue 1:  1 tours + 3 dents et roue 2: 1 tours + 1 dents
Avec a1= 0 et a2= 2 => on trouve c=  8 dents engrenés: Roue 1:  2 tours + 0 dents et roue 2: 1 tours + 2 dents
Avec a1= 1 et a2= 2 pas de solution
Avec a1= 2 et a2= 2 => on trouve c=  2 dents engrenés: Roue 1:  0 tours + 2 dents et roue 2: 0 tours + 2 dents
Avec a1= 3 et a2= 2 pas de solution
Avec a1= 0 et a2= 3 pas de solution
Avec a1= 1 et a2= 3 => on trouve c=  9 dents engrenés: Roue 1:  2 tours + 1 dents et roue 2: 1 tours + 3 dents
Avec a1= 2 et a2= 3 pas de solution
Avec a1= 3 et a2= 3 => on trouve c=  3 dents engrenés: Roue 1:  0 tours + 3 dents et roue 2: 0 tours + 3 dents
Avec a1= 0 et a2= 4 => on trouve c=  4 dents engrenés: Roue 1:  1 tours + 0 dents et roue 2: 0 tours + 4 dents
Avec a1= 1 et a2= 4 pas de solution
Avec a1= 2 et a2= 4 => on trouve c=  10 dents engrenés: Roue 1:  2 tours + 2 dents et roue 2: 1 tours + 4 dents
Avec a1= 3 et a2= 4 pas de solution
Avec a1= 0 et a2= 5 pas de solution
Avec a1= 1 et a2= 5 => on trouve c=  5 dents engrenés: Roue 1:  1 tours + 1 dents et roue 2: 0 tours + 5 dents
Avec a1= 2 et a2= 5 pas de solution
Avec a1= 3 et a2= 5 => on trouve c=  11 dents engrenés: Roue 1:  2 tours + 3 dents et roue 2: 1 tours + 5 dents

On voit bien que, même si m1 et m2 ne sont pas premiers entre eux, il y a tout de même des solutions possibles.

En fait, il y a solution quand on a: (a1-a2) = 0 mod pgcd(m1, m2) et seulement dans ce cas.

On retrouve d'ailleurs que quand m1 et m2 sont premiers entre eux, pgcd(m1,m2)=1 et ça marche quelque soit a1 et a2 (puisque quelque soit (a1-a2): (a1-a2) % 1 == 0).

Dernier point que nous n'avons pas vu à cause de la simulation mécanique des 2 engrenages:

• est-ce que a1, m1, a2, m2 pourraient être négatifs? Oui! Mais m1 et m2 ne peuvent pas être nuls.

Résumons:

Nous cherchons x qui satisfait:

x = a1 mod m1
x = a2 mod m2

• a1, m1, a2, m2 sont des entiers quelconques positifs ou négatifs, mais il faut m1<>0 et m2<>0.

• il y a une solution seulement quand (a1-a2) % pgcd(m1,m2) == 0

• cas particulier: si m1 et m2 sont premier entre eux, il y a toujours une solution (puisque pgcd(m1,m2)==1 et donc (a1-a2)%1==0)

• la solution renvoyée x est modulo ppcm(m1,m2), ce modulo étant renvoyé aussi

Voilà le code proposé pour calculer x et son modulo, optimisé pour 2 équations.

Ce code utilise le PGCD étendu pour calculer les coefficients de Bézout, ainsi que le PPCM (voir autre page sur ce site)

def rchinois2(a1, m1, a2, m2):
    """Calcul des restes chinois pour 2 équations modulaires
       calcul de x qui satisfait: x%m1=a1 et x%m2=a2
       a1, m1, a2 et m2 entiers quelconques >0 ou <0, mais m1<>0 et m2<>0
       la solution existe si (a1-a2)%pgcd(m1,m2)==0. Sinon: renvoi de None
       le résultat est modulo ppcm(m1, m2) (le modulo est renvoyé aussi)
       """
    r, y1, y2 = pgcde(m2, m1) # r = pgcd; y1 et y2 = coef de Bézout de m2 et m1
    if (a1-a2) % r != 0:
        return None, None  # il n'y a pas de solution dans ce cas: on renvoie None
    return ((a1*m2*y1 + a2*m1*y2) % (m1*m2))//r, ppcm(m1, m2)

Par rapport aux codes qu'on trouve sur le net, vous voyez qu'on corrige le dernier résultat en le divisant par r=pgcd(m1, m2): c'est pour que le résultat soit bon même quand m1 et m2 ne sont pas premiers entre eux. Dans ce cas, il faut vérifier qu'il y a bien une solution (si (a1-a2) % r == 0), et si ce n'est pas le cas, de renvoyer un signal à l'appelant (ici None, mais ça pourrait être une exception avec raise).

On passe aux choses sérieuses: résolution d'un nombre quelconque d'équations modulaires:

x = a1 mod m1
x = a2 mod m2
...
x = ai mod mi
...
x = an mod mn

On a ici les contraintes suivantes:

• les ai et mi sont des entiers quelconques positifs ou négatifs, mais les mi sont <>0.

• il y a une solution seulement quand pour tous les i,j (i<>j), on a: (ai-aj)%pgcd(mi,mj)==0

• cas particulier: si les mi,mj sont tous premiers entre eux 2 à 2, il y a toujours une solution quelques soient les ai.

• la solution renvoyée x est modulo ppcm(m1,m2, …, mn), ce modulo étant renvoyé aussi

On a plus de mal à imaginer n engrenages qui engrènent ensemble, alors on va simplifier un peu le modèle: on va avoir ici n engrenages qui vont tous engrener sur une même crémaillère de longueur infinie. A chaque fois que la crémaillère glissera d'une dent, tous les engrenages tourneront d'une dent.

On peut reprendre en le généralisant le code de simulation précédent, ce qui donnera:

def simul(arg):
 
    n = len(arg)  # n = nombre d'équations modulaires à résoudre (>=2)
 
    # vérification qu'il y aura bien une solution
    for i in xrange(1,n):
        ai, mi = arg[i]
        for j in xrange(0,i):
            aj, mj = arg[j]
            if (ai-aj) % pgcd(mi, mj) != 0:
                return None, None
 
    # calcul du ppcm(m1, m2, ..., mn) qui sera le modulo du résultat final
    pp = arg[0][1]
    for a, m in arg[1:]:
        pp = ppcm(pp, m)
 
    c = 0 # compteur des dents qui engrènent avec la crémaillère infinie
    d = [0 for i in xrange(0,n)] # d[i] = compteur de dents sur la roue i
    while True:
 
        # glissement d'une dent de la crémaillère 
        c += 1
 
        # on incrémente le numéro de dent sur chacune des roues 1 à n
        for i in xrange(0,n):
            d[i] = (d[i]+1) % arg[i][1]
 
        # on teste pour savoir si on est arrivé à la position demandée a1 à an
        ok = True
        for i in xrange(0,n):
            ai, mi = arg[i]
            if d[i] != ai % mi:
                ok = False
                break
 
        # on teste pour savoir si on a trouvé
        if ok:
            # oui: on a fini!
            break
 
    # retour de x et de son modulo
    return c%pp, pp 

On aura comme paramètre une liste de couples arg = [[a1,m1], [a2,m2], ..., [an,mn]] représentant les paramètres des équations modulaires à résoudre.

L'avantage de ce code est qu'il donne toutes les solutions possibles sans hypothèse mathématique particulière. C'est donc un excellent moyen pour vérifier un code plus “sérieux” comme le code suivant. Son seul défaut, c'est qu'il est très lent (environ 10000 fois moins rapide…).

Comme nous avons déjà un très bon algorithme pour résoudre 2 équations modulaires, on va s'en servir. Voilà comment.

Dans la liste:

x = a1 mod m1
x = a2 mod m2
...
x = ai mod mi
...
x = an mod mn

On a remplacer les 2 premières équations par leur résultat: a, m = rchinois2(a1,m1,a2,m2) donc x = a mod m

Puis, on remplacera aussi cette nouvelle équation combinée avec la suivante par son résultat a, m = rchinois2(a,m,a3,m3) donc x = a mod m

Etc… jusqu'à la fin!

Voilà le code proposé. Pour gagner encore en rapidité, j'ai intégré le calcul de 2 équations dedans:

def rchinois(arg):
    """Calcul des restes chinois pour n (>=2) équations modulaires
    Argument = liste de n couples: arg = [[a1,m1],[a2,m2],...,[an,mn]]
    Calcul de x qui résout: x%m1==a1,  x%m2==a2 ... x%mn==an
    Avec m1, m2, ..., mn entiers quelconques (>0 ou <0) mais <>0
    Condition d'existence de solution: (ai,aj)%pgcd(mi,mj)==0 (i<>j)
    Retour de x accompagné de son modulo
    """
 
    a1, m1 = arg[0]
    for a2, m2 in arg[1:]:        
        r, y1, y2 = pgcde(m2, m1) # r=pgcd; y1 et y2 =coef de Bézout de m2 et m1
        if (a1-a2) % r != 0:
            return None, None  # il n'y a pas de solution => on renvoie None
        a1 = ((a1*m2*y1 + a2*m1*y2) % (m1*m2))//r
        m1 = ppcm(m1,m2)
    return a1 % m1, m1  # renvoi de x et de son modulo, 

Le paramètre à donner à la fonction est une liste des couples [ai, mi] telle que: arg = [[a1, m1],[a2, m2], ..., [ai, mi], ..., [an, mn]] qui représentent les 2 paramètres ai, mi des n équations modulaires à résoudre.

Vous pouvez vérifier que le résultat x renvoyé est bon en vérifiant que toutes les équations modulaires sont effectivement satisfaites:

def verifsolution(x, arg):
    """Vérifie si x est bien une solution des équations modulaires arg"""
    for ai, mi in arg:
        if (x-ai) % mi != 0:
            return False
    return True

Exemple d'utilisation: reprenons l'exemple chinois cité au début de cette page

x =  6 mod  7 
x = 10 mod 11
x = 12 mod 13

arg = [[6, 7], [10, 11], [12, 13]]
print rchinois(arg)
1000 

Simple, non?

Si on a souvent la résolution de 3 équations modulaires, on peut faciliter l'écriture des paramètres:

def rchinois3(a1,m1,a2,m2,a3,m3):
    return rchinois([[a1, m1], [a2, m2], [a3, m3]])

Voilà un autre exemple de résolution avec n=10 équations modulaires:

x = -175 mod  181   
x =  162 mod -194   
x =  146 mod  -96   
x =  195 mod -139    
x =   65 mod -171    
x =   75 mod   47  
x =  174 mod  -34  
x =  133 mod   59   
x =  -86 mod  137  
x =  106 mod  173  

Résultat:

x = 534140579673186482

Allez, encore un petit pour la route: un exemple de 20 équations modulaires:

x =   41 mod   88  
x = -154 mod  -79   
x =  139 mod -111  
x = -194 mod -173  
x = -136 mod -193   
x =  -43 mod  -38   
x =  109 mod   51   
x =  151 mod  125  
x = -172 mod  -41   
x =  107 mod -109  
x = -172 mod  199    
x =   38 mod  -67   
x =  -36 mod -149   
x =  -38 mod  169  
x = -163 mod   86   
x =  154 mod -157    
x =   20 mod  -61  
x = -196 mod  -97  
x = -169 mod  -71    
x =   83 mod   47  

Ce qui donne:

x = 122955764519783709894623890189898495401

résultat obtenu en 2/10000 de seconde: essayez donc de faire ça à la main

Amusez-vous bien!