aboutsummaryrefslogtreecommitdiff
path: root/coqprime/gencertif/factorize.c
diff options
context:
space:
mode:
Diffstat (limited to 'coqprime/gencertif/factorize.c')
-rw-r--r--coqprime/gencertif/factorize.c675
1 files changed, 0 insertions, 675 deletions
diff --git a/coqprime/gencertif/factorize.c b/coqprime/gencertif/factorize.c
deleted file mode 100644
index d42f72e68..000000000
--- a/coqprime/gencertif/factorize.c
+++ /dev/null
@@ -1,675 +0,0 @@
-/*
-(*************************************************************)
-(* This file is distributed under the terms of the *)
-(* GNU Lesser General Public License Version 2.1 *)
-(*************************************************************)
-(* Benjamin.Gregoire@inria.fr Laurent.Thery@inria.fr *)
-(*************************************************************)
-*/
-
-#include <stdlib.h>
-#include <stdio.h>
-#include <string.h>
-#include "gmp.h"
-#include "ecm.h"
-#include "certif.h"
-
-#if defined (__STDC__) \
- || defined (__cplusplus) \
- || defined (_AIX) \
- || defined (__DECC) \
- || (defined (__mips) && defined (_SYSTYPE_SVR4)) \
- || defined (_MSC_VER) \
- || defined (_WIN32)
-#define __ECM_HAVE_TOKEN_PASTE 1
-#else
-#define __ECM_HAVE_TOKEN_PASTE 0
-#endif
-
-#ifndef __ECM
-#if __ECM_HAVE_TOKEN_PASTE
-#define __ECM(x) __ecm_##x
-#else
-#define __ECM(x) __ecm_/**/x
-#endif
-#endif
-
-#define pp1_random_seed __ECM(pp1_random_seed)
-void pp1_random_seed (mpz_t, mpz_t, gmp_randstate_t);
-#define pm1_random_seed __ECM(pm1_random_seed)
-void pm1_random_seed (mpz_t, mpz_t, gmp_randstate_t);
-#define get_random_ul __ECM(get_random_ul)
-unsigned long get_random_ul (void);
-
-
-static unsigned add[] = {4, 2, 4, 2, 4, 6, 2, 6};
-
-void
-factor_using_division (mpz_t t, pock_certif_t c)
-{
- mpz_t q, r;
- unsigned long int f;
- int ai;
- unsigned *addv = add;
- unsigned int failures;
- unsigned int limit;
-
- /* Set the trial division limit according the size of n. */
- limit = mpz_sizeinbase (t, 2);
- if (limit > 1000)
- limit = 1000 * 1000;
- else
- limit = limit * limit;
-
-
- if (flag_verbose)
- {
- printf ("[using trivial division (%u)] ", limit);
- fflush (stdout);
- }
-
- mpz_init (q);
- mpz_init (r);
-
- f = mpz_scan1 (t, 0);
- mpz_div_2exp (t, t, f);
- while (f)
- {
- if (flag_verbose) { printf ("2 "); fflush (stdout);}
- dec_add_ui(c, 2);
- f--;
- }
-
- for (;;)
- {
- mpz_tdiv_qr_ui (q, r, t, 3);
- if (mpz_cmp_ui (r, 0) != 0) break;
- mpz_set (t, q);
- if (flag_verbose) { printf ("3 "); fflush (stdout); }
- dec_add_ui(c,3);
- }
-
- for (;;)
- {
- mpz_tdiv_qr_ui (q, r, t, 5);
- if (mpz_cmp_ui (r, 0) != 0) break;
- mpz_set (t, q);
- if (flag_verbose) { printf ("5 "); fflush (stdout); }
- dec_add_ui(c,5);
- }
-
- failures = 0;
- f = 7;
- ai = 0;
- while (mpz_cmp_ui (t, 1) != 0)
- {
- mpz_tdiv_qr_ui (q, r, t, f);
-
- if (mpz_cmp_ui (r, 0) != 0)
- {
- f += addv[ai];
- if (mpz_cmp_ui (q, f) < 0) break;
- ai = (ai + 1) & 7;
- failures++;
- if (failures > limit) break;
- }
- else
- {
- mpz_swap (t, q);
- if (flag_verbose) { printf ("%lu ", f); fflush (stdout); }
- dec_add_ui(c,f);
- failures = 0;
- }
- }
-
- if (flag_verbose) fprintf(stdout,"\n");
-
- mpz_clear (q);
- mpz_clear (r);
- return;
-}
-
-void out_factor(mpz_t f,pock_certif_t c)
-{
- mpz_out_str (stdout, 10, f);
- fprintf(stdout," (%lu digits, F1 %lu digits) \n", mpz_sizeinbase(f,10),
- mpz_sizeinbase(c->_F1,10));
- fflush(stdout);
- return;
-}
-
-
-int
-factor_using_pollard_rho (mpz_t n, int a_int, unsigned long p,
- pock_certif_t pc)
-{
- mpz_t x, x1, y, P;
- mpz_t a;
- mpz_t g;
- mpz_t t1, t2;
- int k, l, c, i, res;
-
- if (flag_verbose)
- {
- printf ("[pollard-rho (%d)] ", a_int);
- fflush (stdout);
- }
-
- mpz_init (g);
- mpz_init (t1);
- mpz_init (t2);
-
- mpz_init_set_si (a, a_int);
- mpz_init_set_si (y, 2);
- mpz_init_set_si (x, 2);
- mpz_init_set_si (x1, 2);
- k = 1;
- l = 1;
- mpz_init_set_ui (P, 1);
- c = 0;
-
- res = 0;
-
- while (!res)
- {
-S2:
- if (p != 0)
- {
- mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
- }
- else
- {
- mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
- }
- mpz_sub (t1, x1, x); mpz_mul (t2, P, t1); mpz_mod (P, t2, n);
- c++;
- if (c == 20)
- {
- c = 0;
- mpz_gcd (g, P, n);
- if (mpz_cmp_ui (g, 1) != 0)
- goto S4;
- mpz_set (y, x);
- }
- /*S3: */
- k--;
- if (k > 0)
- goto S2;
-
- mpz_gcd (g, P, n);
- if (mpz_cmp_ui (g, 1) != 0)
- goto S4;
-
- mpz_set (x1, x);
- k = l;
- l = 2 * l;
- for (i = 0; i < k; i++)
- {
- if (p != 0)
- {
- mpz_powm_ui (x, x, p, n); mpz_add (x, x, a);
- }
- else
- {
- mpz_mul (x, x, x); mpz_add (x, x, a); mpz_mod (x, x, n);
- }
- }
- mpz_set (y, x);
- c = 0;
- goto S2;
-S4:
- do
- {
- if (p != 0)
- {
- mpz_powm_ui (y, y, p, n); mpz_add (y, y, a);
- }
- else
- {
- mpz_mul (y, y, y); mpz_add (y, y, a); mpz_mod (y, y, n);
- }
- mpz_sub (t1, x1, y); mpz_gcd (g, t1, n);
- }
- while (mpz_cmp_ui (g, 1) == 0);
-
- if (!mpz_probab_prime_p (g, 3))
- {
- do
- {
- mp_limb_t a_limb;
- mpn_random (&a_limb, (mp_size_t) 1);
- a_int = (int) a_limb;
- }
- while (a_int == -2 || a_int == 0);
-
- if (flag_verbose)
- {
- printf ("[composite factor--restarting pollard-rho] ");
- fflush (stdout);
- }
- res = factor_using_pollard_rho (g, a_int, p, pc);
- break;
- }
- else
- {
- dec_add_mpz(pc, g);
- if (flag_verbose) out_factor(g,pc);
- res = check_pock (pc);
- if (res) break;
- }
- mpz_div (n, n, g);
- mpz_mod (x, x, n);
- mpz_mod (x1, x1, n);
- mpz_mod (y, y, n);
- if (mpz_probab_prime_p (n, 3))
- {
- dec_add_mpz(pc, n);
- if (flag_verbose) out_factor(n,pc);
- res = check_pock (pc);
- break;
- }
- }
-
- mpz_clear (g);
- mpz_clear (P);
- mpz_clear (t2);
- mpz_clear (t1);
- mpz_clear (a);
- mpz_clear (x1);
- mpz_clear (x);
- mpz_clear (y);
-
- return res;
-}
-
-static double B1_table[] =
- { 11000, 50000, 250000, 1000000, 3000000,
- 11000000, 43000000, 110000000, 260000000, 850000000 };
-
-static int it_table[] =
- { 200, 214, 422, 30, 30,
- 20, 20, 20, 20, 20};
-
-static int size_table[] =
- { 20, 25, 30, 35, 40,
- 45, 50, 55, 60, 65 };
-
-int ecm_factorize(mpz_t n, pock_certif_t c);
-
-int my_ecm_factor(mpz_t n, pock_certif_t c, double B1, int iterate)
-{
- int i, res, found;
- mpz_t f;
- gmp_randstate_t randstate;
- ecm_params params;
-
- mpz_init(f);
- ecm_init(params);
- /* if (flag_verbose) params->verbose = 1; */
- gmp_randinit_default (randstate);
- gmp_randseed_ui (randstate, get_random_ul ());
- if (B1 > 11000) params->B1done = 11000;
-
- res = 0;
- i = 0;
- iterate += 5;
- while (i < iterate && !res && mpz_cmp_ui (n, 1) != 0) {
- if (i == 0) { /* start with pm1 */
- if (flag_verbose) {
- printf("using pm1 with B1 = %1.0f ", B1);
- fflush(stdout);
- }
-
- params->method = ECM_PM1;
- pm1_random_seed (params->x, n, randstate);
- found = ecm_factor(f, n, B1, params);
- if (found) {
- mpz_tdiv_q(n, n, f);
- if (mpz_probab_prime_p (f, 3))
- {
- dec_add_mpz(c,f);
- if (flag_verbose) out_factor(f,c);
- res = check_pock(c);
- }
- else
- {
- if (flag_verbose)
- {
- fprintf(stdout,"composite factor ");
- mpz_out_str (stdout, 10, f);
- fprintf(stdout,"(%lu digits)\n",mpz_sizeinbase(f,10));
- fflush(stdout);
- }
- if (B1 == 11000) res = factor_using_pollard_rho(f, 1, 0,c);
- else res = ecm_factorize(f, c);
- }
-
- if (!res && mpz_cmp_ui (n, 1) != 0 && mpz_probab_prime_p (n, 3)) {
- dec_add_mpz(c,n);
- if (flag_verbose) out_factor(n,c);
- mpz_tdiv_q(n, n, n);
- res = check_pock(c);
- }
- } else i++;
- if (flag_verbose) printf("\n");
- } else if (0 < i && i <= 3) {
- /* do 3 time pp1 */
- params->method = ECM_PP1;
- mpz_set_ui (params->x, 0);
- if (flag_verbose && i == 1)
- { printf("using pp1 with B1 = %1.0f ", B1);
- fflush(stdout);
- }
- pp1_random_seed (params->x, n, randstate);
- found = ecm_factor(f, n, B1, params);
- if (found) {
- mpz_tdiv_q(n, n, f);
- if (mpz_probab_prime_p (f, 3))
- {
- dec_add_mpz(c,f);
- if (flag_verbose) out_factor(f,c);
- res = check_pock(c);
- }
- else
- {
- if (flag_verbose)
- {
- fprintf(stdout,"composite factor ");
- mpz_out_str (stdout, 10, f);
- fprintf(stdout,"(%lu digits)\n",mpz_sizeinbase(f,10));
- fflush(stdout);
- }
- if (B1 == 11000) res = factor_using_pollard_rho(f, 1, 0,c);
- else res = ecm_factorize(f, c);
- }
-
- if (!res && mpz_cmp_ui (n, 1) != 0 && mpz_probab_prime_p (n, 3)) {
- dec_add_mpz(c,n);
- if (flag_verbose) out_factor(n,c);
- mpz_tdiv_q(n, n, n);
- res = check_pock(c);
- }
- i = 0; /* restarting to factorize */
- } else i++;
- if (flag_verbose && i == 3) printf("\n");
- } else { /* continue with ecm */
- params->method = ECM_ECM;
- mpz_set_ui (params->x, 0);
- if (flag_verbose && i == 4) {
- printf("using ecm with B1 = %1.0f ", B1);
- fflush(stdout);
- } else {printf("#%i ", i-4); fflush(stdout);}
-
- mpz_urandomb (params->sigma, randstate, 32);
- mpz_add_ui (params->sigma, params->sigma, 6);
-
- found = ecm_factor (f, n, B1, params);
-
- if (found > 0) { /* found a factor */
- mpz_tdiv_q(n, n, f);
- if (mpz_probab_prime_p (f, 3))
- {
- dec_add_mpz(c,f);
- if (flag_verbose) out_factor(f,c);
- res = check_pock(c);
- }
- else
- {
- if (flag_verbose)
- {
- fprintf(stdout,"composite factor ");
- mpz_out_str (stdout, 10, f);
- fprintf(stdout,"(%lu digits)\n",mpz_sizeinbase(f,10));
- fflush(stdout);
- }
- if (B1 == 11000) res = factor_using_pollard_rho(f, 1, 0,c);
- else res = ecm_factorize(f, c);
- }
- if (!res && mpz_cmp_ui (n, 1) != 0 && mpz_probab_prime_p (n, 3)) {
- dec_add_mpz(c,n);
- if (flag_verbose) out_factor(n,c);
- mpz_tdiv_q(n, n, n);
- res = check_pock(c);
- }
- i = 0; /* restarting to factorize */
- } else i++;
- }
- }
- if (flag_verbose) printf("\n");
-
- ecm_clear(params);
- mpz_clear(f);
- return res;
-}
-
-
-int ecm_factorize(mpz_t n, pock_certif_t c)
-{
-
- int iB1, res = 0;
-
- /*res = my_ecm_factor(n, c, B1_table[0], 4);
- if (!res) res = my_ecm_factor(n, c, B1_table[1], 4);
- if (!res) res = my_ecm_factor(n, c, B1_table[2], 4);
- if (!res) res = my_ecm_factor(n, c, B1_table[3], 3);
- if (!res) res = my_ecm_factor(n, c, B1_table[4], 2); */
-
- for (iB1 = 0; !res && iB1 < 10 && mpz_cmp_ui (n, 1) != 0; iB1++)
- {
- if (flag_verbose)
- printf("Searching factor of %i digits\n", size_table[iB1]);
- res = my_ecm_factor(n, c, B1_table[iB1], it_table[iB1]);
- }
-
- return res;
-}
-
-int factorize_no_small(mpz_t n, pock_certif_t c)
-{
- int res;
-
- if (mpz_probab_prime_p (n, 3))
- {
- if (flag_verbose) mpz_out_str (stdout, 10, n);
- dec_add_mpz(c,n);
- res = check_pock(c);
- }
- else
- res = ecm_factorize (n, c);
-
- if (flag_verbose) { fprintf(stdout,"\n");fflush(stdout); }
-
- return res;
-}
-
-int factorize(mpz_t n, pock_certif_t c)
-{
- int res;
-
- /* compute the factorization */
- if (flag_verbose) {
- fprintf(stdout," factorize ");
- mpz_out_str (stdout, 10, n);fflush(stdout);
- fprintf(stdout,"\n ");
- fprintf(stdout," of %lu digits\n", mpz_sizeinbase(n,10));
- fflush(stdout);
- }
-
- factor_using_division (n, c);
-
- res = check_pock(c);
-
- if (!res) res = factorize_no_small(n, c);
-
- if (flag_verbose) { fprintf(stdout,"\n");fflush(stdout); }
-
- return res;
-}
-
-int factorize_mersenne (unsigned long int p, pock_certif_t c)
-{
- unsigned long int q;
- int i,iB1,res,used;
- mpz_t n;
- __mpz_struct dec[100];
-
- if (flag_verbose) {
- fprintf(stdout, "\nfactorize mersenne %lu\n", p);
- fflush(stdout);
- }
-
- used = 0;
- q = p;
- mpz_init (n);
-
- while (q > 3) {
-
- if (q % 2 == 0)
- {
- q = q / 2;
-
- mpz_set_ui(n, 1); /* n = 1 */
- mpz_mul_2exp(n, n, q); /* n = 2^q */
- mpz_add_ui(n, n, 1); /* n = 2^q + 1 */
-
- factor_using_division(n,c);
-
- mpz_init_set(&(dec[used]), n);
- used++;
- }
- else if (q % 3 == 0 )
- {
- q = q /3;
-
- mpz_set_ui (n,1); /* n = 1 */
- mpz_mul_2exp (n, n, q); /* n = 2^q */
- mpz_add_ui (n, n, 1); /* n = 2^q + 1 */
- mpz_mul_2exp (n, n, q); /* n = 2^(2q) + 2^q */
- mpz_add_ui (n, n, 1); /* n = 2^(2q) + 2^q + 1 */
-
- factor_using_division (n,c);
- mpz_init_set(&(dec[used]), n);
- used++;
- }
- else break;
-
- }
-
- switch (q) {
- case 1:
- break;
- case 2:
- dec_add_ui(c,3);
- break;
- case 3:
- dec_add_ui(c,7);
- break;
- default:
- mpz_set_ui(n, 1);
- mpz_mul_2exp(n, n, q);
- mpz_sub_ui(n, n, 1);
- factor_using_division (n,c);
- mpz_init_set(&(dec[used]), n);
- used++;
- break;
- }
-
- res = check_pock(c);
- iB1 = 0;
- while (!res && iB1 < 10) {
- for (i = 0; i < used && !res; i++) {
- if (mpz_cmp_ui (&(dec[i]), 1) != 0)
- res = my_ecm_factor(&(dec[i]), c, B1_table[iB1], it_table[iB1]);
- }
- iB1++;
- }
-
- mpz_clear(n);
- return res;
-}
-
-
-pock_certif_t mersenne_certif (mpz_t t, unsigned long int p)
-{
- pock_certif_t c;
- c = pock_init(t);
- dec_add_ui(c, 2);
- factorize_mersenne (p-1, c);
- finalize_pock(c);
- return c;
-}
-
-
-
-pock_certif_t pock_certif (mpz_t t)
-{
-
- mpz_t tm1;
-
- pock_certif_t c;
-
- if (flag_verbose) {
- fprintf(stdout,"pocklington ");
- mpz_out_str (stdout, 10, t);fflush(stdout);
- fprintf(stdout,"\n");fflush(stdout);
- }
-
- /* initialize the decompostion */
- c = pock_init(t);
-
- /* compute t - 1 */
- mpz_init_set (tm1, c->_R1);
-
- /* compute the factorisation */
- factorize(tm1, c);
-
- mpz_clear(tm1);
-
- finalize_pock(c);
-
- return c;
-}
-
-
-
-int MAXPROOFPRIMES = 48611; /* 5000 first ones */
-
-pre_certif_t certif_2;
-
-void extend_lc (certif_t lc, pock_certif_t c, unsigned long int min,
- unsigned long int max )
-{
- int i, size;
- mpz_ptr *ptr;
- mpz_t t;
- mpz_init (t);
-
- ptr = c->_dec;
- size = c->_used;
-
- if (c->_pow2 > 0 && !_2_is_in(lc)) {
- mpz_t t2;
- pre_certif_t ct;
- mpz_init_set_ui (t2, 2);
- ct = mk_proof_certif(t2);
- add_pre(ct, lc);
- mpz_clear(t2);
- }
-
- for(i = size - 1; i >= 0; i--)
- {
- mpz_set(t, ptr[i]);
- if (!is_in(t, lc)) {
- pre_certif_t ct;
- if (mpz_cmp_ui(t, MAXPROOFPRIMES) <= 0 ||
- (mpz_cmp_ui (t, min) >= 0 && (mpz_cmp_ui (t, max) <= 0)))
- ct = mk_proof_certif(t);
- else {
- ct = mk_pock_certif (pock_certif(t));
- extend_lc(lc, ct->_certif._pock, min, max);
- }
- add_pre(ct, lc);
- }
- }
- return;
-}
generated by cgit v1.2.3 (git 2.39.1) at 2024-04-30 12:32:58 +0000