Cornacchia algorithm

Soyez le premier à donner votre avis sur cette source.

Snippet vu 8 220 fois - Téléchargée 37 fois

Contenu du snippet

Ce code permet de résoudre l'équation x^2 + y^2 = p, p premier. Très pratique en cryptographie ;o)

Source / Exemple :


/* CORNACCHIA Algorithm
 *

  • GOAL: given d and p prime, find (x,y) such that x^2 + y^2 = d * p
*
  • in this implementation, d = 1. This code does not accept d <> 1 !!!
  • see "A Course in Computational Algebraic Number Theory" by Henri Cohen
  • coded by malik@hammoutene.com, 2004,the 5th of August
*
  • /
#pragma comment(lib, "gmp.lib") #include <stdlib.h> #include <stdio.h> #include <string.h> #include "gmp.h" #define _WIN32_WINNT 0x0400 #include <windows.h> #include <wincrypt.h> int genPrime(){ // GENERATION DE NOMBRE PREMIER mpz_t t,t2, prime_n; int d2; d2=getrand(); mpz_init(t); mpz_set_ui(t,d2); mpz_mul_ui(t,t,4); mpz_add_ui(t,t,1); if (mpz_probab_prime_p(t,20)==0){ while (mpz_probab_prime_p(t,20)==0){ d2=getrand(); mpz_set_ui(t,d2); mpz_mul_ui(t,t,4); mpz_add_ui(t,t,1); } } mpz_init(t2); mpz_set(t2,t); mpz_set(t,t2); // petite finte pour le return mpz_init(prime_n); mpz_set(prime_n,t); mpz_clear(t); mpz_clear(t2); return (prime_n); } int getrand() // GENERATION DE NOMBRE ALEATOIRE { HCRYPTPROV hProv = 0; int rnd; CryptAcquireContext(&hProv, 0, 0, PROV_RSA_FULL, CRYPT_VERIFYCONTEXT); CryptGenRandom(hProv, sizeof(rnd), (BYTE*)&rnd); CryptReleaseContext(hProv, 0); return abs(rnd/1000); // pour être dans un ordre honnête (~ 2^20 et plus) } int power(int val, int pow) { int ret_val = 1; int i; for(i = 0; i < pow; i++) ret_val *= val; return(ret_val); } //************* SHANKS-TONELLI Algorithm int shanks(mpz_t a2, mpz_t prime_number){ int random_n; mpz_t a,a_prime,b,x,y,z,q,t,temp,temp2; mpz_t rand_n; int e = 0; int compteur=1; int r,e2,m,mp,r2; mpz_init(a); mpz_sub_ui(a,prime_number,1); mpz_init(q); mpz_init(a_prime); mpz_mod_ui(q,a,2); while(mpz_cmp_ui(q,0) == 0) { mpz_div_ui(a_prime, a, 2); e = e+1; mpz_set(a,a_prime); mpz_mod_ui(q,a,2); } e2 = power(2,e); mpz_div_ui(q, prime_number, e2); // 1 Find generator random_n = getrand(); mpz_init(rand_n); mpz_set_ui(rand_n,random_n); if (mpz_legendre(rand_n, prime_number) != -1){ while (mpz_legendre(rand_n, prime_number) != -1){ random_n = getrand(); mpz_set_ui(rand_n,random_n); } } mpz_init(z); mpz_powm(z,rand_n,q,prime_number); // 2 Initialize mpz_init(y); mpz_set(y,z); r = e; mpz_init(temp); mpz_sub_ui(temp, q, 1); mpz_div_ui(temp, temp, 2); mpz_init(x); mpz_powm(x,a2,temp,prime_number); mpz_pow_ui(temp,x,2); mpz_mul(temp, temp, a2); mpz_init(b); mpz_mod(b,temp,prime_number); mpz_mul(temp, a2, x); mpz_mod(x,temp,prime_number); // 3 Find exponent mpz_mod(temp, b, prime_number); mpz_init(t); if (mpz_cmp_ui(temp,1) != 0){ mpz_init(temp2); while(mpz_cmp_ui(temp,1) != 0){ m=1; mp=power(2,m); mpz_powm_ui(temp2,b,mp,prime_number); while(mpz_cmp_ui(temp2,1) !=0) { m++; mp=power(2,m); mpz_powm_ui(temp2,b,mp,prime_number); } if (m==r){ printf("a is not a quadratic residue mod p...\n"); exit(0);} // 4 Reduce exponent r2 = r-m-1; r2 = power(2,r2); mpz_powm_ui(t,y,r2,prime_number); mpz_powm_ui(y,t,2,prime_number); r=m; mpz_mul(x,x,t); mpz_mod(x,x,prime_number); mpz_mul(b,b,y); mpz_mod(b,b,prime_number); mpz_mod(temp, b, prime_number); } } mpz_clear(a); mpz_clear(a2); mpz_clear(a_prime); mpz_clear(b); mpz_clear(rand_n); mpz_clear(y); mpz_clear(z); mpz_clear(q); mpz_clear(t); mpz_clear(temp); mpz_clear(temp2); return (x); } struct cornAlgo{ mpz_t x_sol; mpz_t y_sol; }; struct cornAlgo cornacchia(mpz_t prime){ struct cornAlgo resultats; int d = 1; int k2; mpz_t temp,temp3,a2,a,x01,b,l,c,r; mpz_t prime_number; mpz_init(prime_number); mpz_set(prime_number,prime); //************* CORNACCHIA ALGORITHM // 1 Test if residue mpz_init(temp); mpz_set_ui(temp,d); mpz_neg(temp,temp); k2= mpz_legendre(temp,prime_number); mpz_clear(temp); if (k2 == -1) { printf("the equation has no solution\n"); exit (0);} // 2 Compute square root mpz_init(a2); mpz_set_ui(a2,d); mpz_neg(a2,a2); mpz_init(x01); mpz_set(x01, shanks(a2,prime_number)); if (mpz_cmp(x01,prime_number)>0){ mpz_neg(x01,x01); } mpz_init(temp3); mpz_div_ui(temp3, prime_number,2); mpz_add_ui(temp3, temp3, 1); if (mpz_cmp(x01,temp3)>=0){ while (!(mpz_cmp(x01,temp3)>=0)){ mpz_add(x01, x01, prime_number); } } mpz_init(a); mpz_set(a,prime_number); mpz_init(b); mpz_set(b,x01); mpz_init(l); mpz_sqrt(l,prime_number); // 3 Euclidean algorithm mpz_init(r); if (mpz_cmp(b,l)>0){ while ((mpz_cmp(b,l)>0)){ mpz_mod(r,a,b); mpz_set(a,b); mpz_set(b,r); } } // 4 Test solution: ici, comme d = 1, les tests ne sont pas utiles mpz_pow_ui(temp3,b,2); mpz_neg(temp3,temp3); mpz_add(temp3,temp3,prime_number); mpz_init(c); mpz_set(c, temp3); mpz_sqrt(c,c); mpz_init(resultats.x_sol); mpz_set(resultats.x_sol,b); mpz_init(resultats.y_sol); mpz_set(resultats.y_sol,c); mpz_clear(temp3); mpz_clear(a); mpz_clear(x01); mpz_clear(l); mpz_clear(r); mpz_clear(prime_number); return resultats; } //************* MAIN main(int argc, char *argv[]) { struct cornAlgo resultatsXY; mpz_t prime_number; // génération d'un nombre premier tel que p = 1 + 4k mpz_init(prime_number); mpz_set(prime_number,genPrime()); // exécution de l'algorithme de Cornacchia resultatsXY = cornacchia(prime_number); // affichage du résultat mpz_out_str (stdout, 10, prime_number); printf(" = "); mpz_out_str (stdout, 10, resultatsXY.x_sol); printf("^2 + "); mpz_out_str (stdout, 10, resultatsXY.y_sol); printf("^2\n"); }

Conclusion :


N'oubliez pas d'installer GMP pour faire fonctionner ce code (j'ai mis une marche à suivre dans le forum).

Merci de me tenir au courant en cas d'amélioration (par exemple, si quelqu'un traite le cas d>1)

A voir également

Ajouter un commentaire

Commentaires

Messages postés
147
Date d'inscription
samedi 1 août 2009
Statut
Membre
Dernière intervention
5 novembre 2019

Malik je te pose la meme question qu'au autre user de gmp. As tu utilisé gmp avec le module gmpxx.h(cpp) et vc6++ (ou plus) si oui pourrais tu me faire parvenir le source de gmpxx.h sil te plai. Bonne route.Christophe
Messages postés
1154
Date d'inscription
mardi 9 septembre 2003
Statut
Membre
Dernière intervention
15 août 2009
14
J'aurais tendance à dire oui: les unsigned on un bit de moins (le bit de poids fort) donc ils sont plus courts, donc traîtés plus vite!
Messages postés
117
Date d'inscription
mercredi 3 septembre 2003
Statut
Membre
Dernière intervention
17 février 2007

J'ai cherché dans le manuel de GMP si les fonctions unsigned integer (ui) étaient plus rapides que les fonctions signed integer (si) ou que les fonctions "normales" (indéterminées). Sans résultat. Quelqu'un saurait m'informer?

(je ne travaille qu'avec des positifs. vaut-il mieux que j'utilise les fonctions "ui", les fonctions "si" ou les fonctions normales (indéterminées) ?)
Messages postés
117
Date d'inscription
mercredi 3 septembre 2003
Statut
Membre
Dernière intervention
17 février 2007

OK!!! unsigned integer = entier relatif! il suffisait d'avoir un bon traducteur.

Donc, logiquement, on emploie les fonctions qui se terminent par "ui" lorsqu'on veut manipuler des entiers relatifs (négatifs et positifs) et les fonctions de base (sans "ui") lorsqu'on ne manipule que des entiers positifs (ce qui est mon cas)...

Pourquoi le français n'est-il pas la première langue parlée au monde??

ScelW
Messages postés
1154
Date d'inscription
mardi 9 septembre 2003
Statut
Membre
Dernière intervention
15 août 2009
14
Il n'y est pas écrit "unsigned integer is...", mais au travers des fonctions on comprends ce que c'est!

Mais si ça suffit pas: http://www.rwc.uc.edu/koehler/comath/13.html

google is your friend ;o)
Afficher les 14 commentaires

Vous n'êtes pas encore membre ?

inscrivez-vous, c'est gratuit et ça prend moins d'une minute !

Les membres obtiennent plus de réponses que les utilisateurs anonymes.

Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.

Le fait d'être membre vous permet d'avoir des options supplémentaires.