Il y a quelques années, lorsque nous penchions sur le Green Code, nous avions trouvé une méthode de factorisation intéressante appelée communément appelé ECM (pour Elliptic Curve Method). Etant curieux, nous avons voulu percer les zones ombres de cette méthode. Notre article précédent nous a permis de définir les bases sur les courbes elliptiques. Cela nous sera utile pour cet article.
Hendrik Lenstra a, en 1985, découvert/inventé un algorithme de factorisation qui utilise des courbes elliptiques. Cet algorithme a la particularité de permettre de trouver les petits facteurs premiers d’un nombre aussi grand qu’il soit. Et là où la plus part des algorithmes de factorisation dépendent de la taille du nombre à factoriser celui-ce dépend principalement de la taille du facteur à trouver.
L’idée de la factorisation par courbe elliptique consiste à choisir une courbe elliptique modulo 
Comme nous l’avons vu, pour le calcul de la somme de deux points sur une courbe elliptique il faut calculer la pente m. Cette pente contient une division qui en arithmétique modulaire n’est rien d’autre que le calcul d’un inverse. Or cet inverse n’existe pas toujours dans le cas où
Petit exemple avec 91
Plus concrêtement, essayons de factoriser le nombre 
On choisit un point de la courbe : par exemple P (37,37) est bien sur la courbe car :
Puis on l’additionne à lui-même plusieurs fois : Au début tout se passe bien : P = (37,37) ; 2P = (56,15) ; 3P = (49,20) ; 4P = (86,40)
Mais pour 5P = 4P + P = (86,40)+(37,37), on doit calculer 
Comme le calcul échoue, on sait que 
Notez que bien que le point 5P n’existe pas, le point 6P existe bel et bien sur la courbe :(23,58).
Les calculs échouent pour k = 5 ( et tous ces multiples) car 7 est un diviseur de 91 et que 5 est le nombre d’éléments de E(1,1) modulo 7.
Projetons-nous
Si l’on veut comprendre d’où vient ce lien, il faut comprendre les liens des courbes E(a,b) modulo 91 et modulo 7.
Chaque point de la courbe modulo 91 peut être « projeter » sur la courbe modulo 7 simplement en réduisant ces coordonnées modulo 7.
Tout d’abord remarquons que sur la courbe E(1,1) modulo 7, il n’y a que 5 éléments.
sage: Z = Zmod(7) sage: E = EllipticCurve(Z, (1,1)) sage: E.points() [(0 : 1 : 0), (0 : 1 : 1), (0 : 6 : 1), (2 : 2 : 1), (2 : 5 : 1)] sage: E.cardinality() 5 sage:
Si l’on reprend le point 
Mais c’est Absurde
On peut continuer de vérifier les additions pour les multiples suivants. Vient alors le cas 
Si nous supposons qu’un tel point existe, alors on pourrait calculer 
Or on sait que l’ordre de cette courbe est 5 et donc d’après Lagrange, 5 fois n’importe quel point est égal à l’élément neutre. Dans le cas des courbes elliptiques, l’élément neutre a des coordonnées infinies.
D’où l’absurdité.
Quels paramètres choisir ?
Revenons au cas général. Si nous souhaitons utiliser cette technique pour factoriser un entier quelconque, plusieurs questions se posent :
Y-a-t-il un point P à choisir en particulier ?
La réponse à cette question est assez simple : pour l’instant nous n’avons pas de moyen de déterminer à l’avance si un point sera meilleurs qu’un autre. Le point est donc souvent choisi au hazard.
Quelle est la valeur optimale pour k (le nombre de fois qu’on multiplie le point P) ?
Pour répondre à cette question, posons-nous en une autre quelle sont les propriétés que nous souhaitons avoir pour k ? De plus, nous avons vu plus haut dans notre exemple que 5*P n’existait pas (ni ses multiples) mais que 6*P lui existait. Autrement dit, nous cherchons un k qui soit au moins égal à l’ordre de la courbe modulo le facteur cherché. Bien entendu, nous n’avons aucune idée de l’ordre en question.
Soyons naïfs :
Du coup une approche naïve pourrait être de prendre 
Si nous cherchons les facteurs plus petit que 1000000 on sait que l’ordre maximum sera 
Mais ça fait quand même beaucoup…
Soyons attentifs !
Nous pourrions donc utiliser non plus la factorielle mais la primorielle. Cela réduit significitivement notre nombre. Cependant cela reste incomplet et nous pouvons « probablement » mieux faire.
Cela reste incomplet car si l’ordre pour le facteur chercher est par exemple 9. La multiplication par k = 11# (primorielle 11 = 2*3*5*7*11), ne donnera pas le résultat escompté car nous n’avons pas de puissance des nombres premiers.
Idéalement, il faut choisir une borne B et calculer le produit de toutes les puissances de nombres premiers inférieurs à B. Pour 27, nous aurions 
Mais le nombre k est toujours très (trop) grand pour une utilisation en pratique
Par ailleurs, les lecteurs assidus de ce blog auront peut-être reconnu un problème se rapproche d’une notion déjà abordée ici. En effet, les plus attentifs auront probablement fait le lien avec la friabilité introduite lorsque nous parlions de l’algorithme p-1 de Pollard et l’optimisation de k.
Les valeurs de a et b varient chacune de 0 à N, Il y a un peu moins de 
Cependant on doit aussi s’intéresser aux ordres super-friables. Grâce à un petit script, nous obtenons les résultats suivants qui donnent une tendance.
| P | Ordres min et max | Nombre d’ordres différents | Valeurs B-smooth avec B = Sqrt(P) | Ratio |
| 11 | (5, 18) | 13 | 0 | 0 |
| 101 | (81, 122) | 41 | 4 | 0.097 |
| 1009 | (946, 1073) | 127 | 22 | 0.173 |
| 10009 | (9809, 10210) | 401 | 83 | 0.206 |
| 100003 | (99371, 100636) | 1265 | 314 | 0.248 |
| 1000003 | (998003, 1002004) | 4001 | 1036 | 0.258 |
| 10000019 | (9993695, 10006344) | 12649 | 3371 | 0.266 |
| 100000007 | (99980007, 100020008) | 40001 | 10961 | 0.274 |
| 1000000007 | (999936762, 1000063253) | 126491 | 35247 | 0.278 |
En effet, il semble que plus P augmente, plus le nombre de nombre B-superfriable (avec B = sqrt(B)) augmente. Nos calculs montrent que si nous avons un P proche de 
Pour la suite, soyons probabilistes !
En utilisant la notion de friabilité et le fait que nous pouvons tester plusieurs courbes, il est peut-être possible de drastiquement réduire le nombre k. En effet, avec un peu de chance, nous pouvons tomber sur une courbe d’ordre B-super-friable, avec un B de taille raisonnable.
Puisque le nombre total de courbes est très grand, il n’est pas possible d’effectuer des recherches sur toutes les courbes, nous devons nous limiter à effectuer des calculs sur un sous-ensemble des courbes mais combien précisement ? Il y a deux aspects à prendre en compte : il nous faut d’un coté avoir suffisament de courbes pour avoir une chance (probabilité) raisonnable de trouver un facteur de la taille choisie et d’un autre coté, il nous faut avoir un minimum de courbe pour que les temps de calculs soient compétitifs par rapport aux autres méthodes de recherche de facteur.
Nous avons donc réalisé une implémentation de l’algorithme ECM en utilisant l’outils SageMath. Le script donc le code est ci-dessous nous a permis d’étudier les valeurs optimales du nombre de courbe et de B en fonction de la taille du facteur à trouver. Le résultat est donné dans le tableau ci-aprés.
Taille du B Nb Moyen Nb de courbe pour les 10 tests facteur Nb Médian réalisés pour chaque B 1.00E+05 20 1322.8 1086 [849, 3276, 1932, 73, 1323, 2471, 599, 414, 187, 2104] 1.00E+05 50 41.5 24 [72, 25, 3, 23, 13, 18, 23, 140, 65, 33] 1.00E+05 100 18.3 15.5 [15, 11, 35, 5, 16, 35, 9, 14, 20, 23] 1.00E+05 200 5.8 3.5 [13, 3, 1, 12, 14, 3, 1, 5, 4, 2] 1.00E+05 500 3.1 2.5 [3, 2, 9, 1, 2, 3, 5, 4, 1, 1] 1.00E+05 1000 2 1 [5, 1, 1, 1, 2, 1, 2, 5, 1, 1] 1.00E+05 1500 2.2 2 [2, 1, 1, 2, 4, 5, 2, 2, 2, 1] 1.00E+06 20 5263.1 5603.5 [416, 26, 10000, 10000, 2149, 10000, 7829, 5929, 1004, 5278] 1.00E+06 50 156.4 128.5 [376, 31, 97, 160, 69, 168, 97, 358, 28, 180] 1.00E+06 100 21.6 12.5 [14, 10, 6, 72, 11, 1, 1, 25, 36, 40] 1.00E+06 200 12.8 7 [29, 21, 1, 4, 15, 1, 8, 4, 6, 39] 1.00E+06 500 4 3.5 [4, 3, 4, 3, 2, 7, 2, 4, 2, 9] 1.00E+06 1000 3.5 3.5 [4, 3, 2, 5, 4, 6, 6, 2, 1, 2] 1.00E+06 1500 3.1 2.5 [5, 3, 1, 1, 2, 5, 7, 1, 2, 4] 1.00E+07 20 10000 10000 [10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000] 1.00E+07 50 1204.3 802 [183, 794, 930, 2831, 372, 810, 55, 1747, 3974, 347] 1.00E+07 100 210.1 204 [67, 315, 71, 373, 10, 8, 555, 275, 133, 294] 1.00E+07 200 32.7 14.5 [10, 18, 3, 123, 1, 35, 83, 32, 11, 11] 1.00E+07 500 15.8 12 [14, 14, 6, 2, 7, 40, 5, 43, 17, 10] 1.00E+07 1000 3.6 3 [1, 1, 3, 6, 6, 4, 3, 3, 3, 6] 1.00E+07 1500 6.7 5 [5, 1, 8, 16, 5, 13, 3, 5, 4, 7] 1.00E+08 20 10000 10000 [10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000] 1.00E+08 50 2408.1 1555 [2199, 1213, 2309, 1724, 3939, 1386, 10000, 506, 385, 420] 1.00E+08 100 460.3 521 [599, 603, 127, 203, 868, 549, 1020, 90, 493, 51] 1.00E+08 200 53.2 50 [9, 51, 14, 49, 147, 67, 91, 3, 84, 17] 1.00E+08 500 23.1 16.5 [8, 25, 6, 18, 83, 7, 16, 13, 38, 17] 1.00E+08 1000 15.6 15 [4, 38, 17, 15, 12, 3, 22, 27, 3, 15] 1.00E+08 1500 9.4 3.5 [6, 1, 1, 2, 1, 4, 54, 13, 9, 3] 1.00E+09 20 10000 10000 [10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000] 1.00E+09 50 10000 10000 [10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000] 1.00E+09 100 3262.9 2229 [2496, 1305, 3913, 272, 282, 9391, 2991, 1962, 10000, 17] 1.00E+09 200 284.2 186.5 [92, 631, 123, 799, 169, 399, 168, 204, 233, 24] 1.00E+09 500 43 43.5 [41, 46, 49, 19, 4, 56, 88, 5, 119, 3] 1.00E+09 1000 25.7 18 [85, 15, 3, 21, 9, 36, 35, 13, 1, 39] 1.00E+09 1500 18.3 11.5 [62, 8, 15, 5, 19, 5, 5, 30, 26, 8] 1.00E+11 200 4203.5 2883.5 [479, 6265, 10000, 10000, 491, 6260, 3211, 2556, 1825, 948] 1.00E+11 500 188.9 203 [263, 238, 93, 362, 22, 229, 342, 57, 177, 106] 1.00E+11 1000 168.2 137.5 [134, 144, 223, 230, 567, 84, 71, 141, 80, 8] 1.00E+11 1500 109 58.5 [8, 40, 319, 35, 252, 57, 160, 60, 37, 122] 1.00E+12 200 7652.1 10000 [10000, 10000, 1126, 1425, 10000, 10000, 6572, 10000, 10000, 7398] 1.00E+12 500 1156.9 872 [403, 1414, 3639, 954, 1297, 303, 487, 663, 1619, 790] 1.00E+12 1000 288.3 186 [165, 207, 72, 728, 140, 324, 398, 628, 125, 96] 1.00E+12 1500 82.4 81 [74, 12, 12, 253, 91, 93, 1, 35, 88, 165] 1.00E+13 200 9562.1 10000 [10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 10000, 5621] 1.00E+13 500 3742.1 3837.5 [1038, 8005, 5419, 4126, 4256, 1254, 3449, 5933, 3549, 392] 1.00E+13 1000 746.6 645 [577, 1545, 389, 713, 1412, 1110, 85, 522, 921, 192]
Ici notons que nous avons une approche très naïve. Nous utilisons beaucoup de fonction aléatoire pour générer nos points et nos courbes. Cela donne que des résultats bien loin de standard moderne dont le raisonnement est hors de notre portée. Mais les mathématiciens sont des gens plein de ressources et ont trouvé des moyens pour encore rafiner et abaisser la borne B et le nombre de courbes en remarquant que certains courbes notamment ont des propriétés intéressantes.
Y-a-t-il une courbe (des coefficients a et b) à choisir en particulier ?
Cette dernière question est en effet importante. Les coefficients a et b permettent d’influer sur l’ordre de la courbe et comme le rappelle le théorème de Haas, l’ordre d’une courbe modulo p varie entre 
17: [4, 8, 24, 8, 16, 24, 28, 8, 32, 8, 28, 24, 16, 8, 24, 8, 4] 17 17 53: [39, 26, 104, 52, 104, 130, 52, 26, 260, 52, 117, 104, 156, 78, 156, 78, 156, 104, 117, 52, 260, 26, 52, 130, 104, 52, 104, 26, 39] 29 29 79: [52, 156, 78, 156, 91, 156, 156, 312, 78, 390, 78, 156, 156, 364, 117, 156, 234, 390, 234, 156, 117, 364, 156, 156, 78, 390, 78, 312, 156, 156, 91, 156, 78, 156, 52] 35 35 101: [25, 50, 300, 100, 100, 250, 300, 100, 400, 150, 500, 200, 200, 200, 600, 150, 200, 400, 375, 100, 700, 100, 375, 400, 200, 150, 600, 200, 200, 200, 500, 150, 400, 100, 300, 250, 100, 100, 300, 50, 25] 41 41 127: [126, 63, 336, 147, 252, 252, 630, 189, 252, 378, 504, 315, 252, 126, 1008, 504, 378, 252, 504, 189, 756, 273, 630, 273, 756, 189, 504, 252, 378, 504, 1008, 126, 252, 315, 504, 378, 252, 189, 630, 252, 252, 147, 336, 63, 126] 45 45 173 [129, 86, 516, 86, 516, 602, 344, 258, 1032, 172, 516, 602, 1032, 430, 688, 430, 516, 860, 516, 258, 2064, 344, 559, 430, 860, 430, 1204, 430, 860, 430, 559, 344, 2064, 258, 516, 860, 516, 430, 688, 430, 1032, 602, 516, 172, 1032, 258, 344, 602, 516, 86, 516, 86, 129] 53 53 257: [64, 128, 896, 256, 768, 1024, 768, 256, 1024, 384, 1536, 896, 768, 512, 2816, 640, 512, 1280, 1792, 896, 2048, 384, 768, 640, 1536, 1024, 3072, 512, 512, 1664, 1984, 512, 2048, 512, 1984, 1664, 512, 512, 3072, 1024, 1536, 640, 768, 384, 2048, 896, 1792, 1280, 512, 640, 2816, 512, 768, 896, 1536, 384, 1024, 256, 768, 1024, 768, 256, 896, 128, 64] 65 65
Même si aucune formule actuellement ne permet de déterminer l’ordre d’une courbe en fonction de a,b et p, En autre des mathématiciens en utilisant des notions que nous détaillerons pas ici (comprendre qu’on n’a pas vraiment compris) ont trouvé des sous-ensembles de courbes dont l’ordre est au moins toujours divisibles par 12 (Paramétrage de Sumaya). Cela donne un certain avantage et permet de réduire le nombre de courbes à tester.
Mais l’algorithme a encore et toujours été raffiné : tant sur l’implémentation (pour gagner en vitesse de calcul) que sur les choix des courbes et que celui de k.
Ainsi les dernieres valeurs recommandées lors de l’utilisation de GMP-ECM (la référence pour la factorisation par courbes elliptiques) sont les suivantes :
| Digits | optimal B1 | expected curves (default parameters for GMP-ECM 6) |
expected curves (default parameters for GMP-ECM 7) |
| 20 | 11,000 | 86 | 107 |
| 25 | 50,000 | 214 | 261 |
| 30 | 250,000 | 430 | 513 |
| 35 | 1,000,000 | 910 | 1,071 |
| 40 | 3,000,000 | 2,351 | 2,753 |
| 45 | 11,000,000 | 4,482 | 5,208 |
| 50 | 43,000,000 | 7,557 | 8,704 |
| 55 | 110,000,000 | 17,884 | 20,479 |
| 60 | 260,000,000 | 42,057 | 47,888 |
| 65 | 850,000,000 | 69,471 | 78,923 |
| 70 | 2,900,000,000 | 102,212 | 115,153 |
| 75 | 7,600,000,000 | 188,056 | 211,681 |
| 80 | 25,000,000,000 | 265,557 | 296,479 |
Voici quelques lectures intéressantes pour ce qui voudraient aller plus loin :
- FINDING ECM-FRIENDLY CURVES THROUGH A STUDY OF GALOIS PROPERTIES
- Implementing the Elliptic Curve Method of Factoring in Reconfigurable Hardware
- A PRACTICAL ANALYSIS OF THE ELLIPTIC CURVE FACTORING ALGORITHM
- 20 years of ECM
En espérant que vous avez réussi à suivre, sinon n’hésitez pas à laisser un commentaire .
Implémentation naïve d’ECM avec Sage:
#!/usr/bin/env sage
from numpy import median
from itertools import groupby
from sage.all import *
'''
Naive implemenation of ECM method
for test and education purposes
'''
def factorECM(k,m,verbose):
condition1 = True
condition2 = True
Z = Zmod(m)
curvecount = 0
x = 0
y = 1
while curvecount 0):
curvecount += 1
#The curve is defined modulo M
E = EllipticCurve(Z, (a,b))
if(verbose):
print "Curve", a ,b
'''For some reasons, Sage does not allow us to
use E.random_point() when M is composite
so we have to compute a point on the curve
We simply try various value of x and hope to find one
'''
while condition2 :
x = ZZ.random_element(0,m)
y2 = x**3 + a*x + b % m
try:
y = mod(y2,m).sqrt(extend=False)
except:
continue
break
P = E(x,y)
#We compute k*P on catch the error if any, if not we'll try another curve
try:
Q = k*P
except ZeroDivisionError as err:
d = Integer(err.args[0].split()[2])
if(verbose):
print "After testing ", curvecount, " curve(s), we have a winner : "
print "\t Point P :" , P
print "\t on E(",a,",",b,")"
print "\t In K*P computation, d=", d ,"has no inverse modulo", m
print "\t So we found one factor gcd(m,d) =", gcd(m,d)
print "m =", m , " = ", gcd(m,d),"x", m/gcd(m,d)
return curvecount
if(verbose):
print "No factor found with the given parameter : try increase B"
return curvecount
def main(argv):
if len(sys.argv) != 4:
print("Usage: %s " % sys.argv[0])
print("Naive implementation of ECM for factorization")
print("v : Try to factor with ")
print("q : Do some benchmark with random composite and predefine B")
sys.exit(1)
#The number we want to factor
m = sage_eval(sys.argv[2])
#super-smooth border B
B = sage_eval(sys.argv[3])
if sys.argv[1] == 'v':
#We compute the number k
k = LCM(2..previous_prime(B))
factorECM(k,m,True)
return 0
if sys.argv[1] == 'q':
Blist = (200,500,1000,1500)
for i in range (11,16):
rg = 10**i
for B in Blist:
nbcurve = list()
k = LCM(2..previous_prime(B))
for l in range(10):
p = next_prime(ZZ.random_element(rg,rg*10))
q = next_prime(ZZ.random_element(rg*321,rg*14311))
nbcurve.append(factorECM(k,p*q,False))
print rg, B,
print float(sum(nbcurve) / len(nbcurve)), float(median(nbcurve)),
print nbcurve
if __name__ == "__main__":
main(sys.argv[1:])
<span id="mce_SELREST_start" style="overflow:hidden;line-height:0;"></span>
Cryptobourrin 