Math-Prime-Util-GMP

 view release on metacpan or  search on metacpan

bls75.c  view on Meta::CPAN

 *     - gcd(a^((n-1)/f)-1,n) = 1  for ALL factors f of A
 * then:
 *     if s = 0 or r*r - 8*s is not a perfect square
 *         n is prime
 *     else
 *         n is composite
 *
 * The generalized Pocklington test is also sometimes known as the
 * Pocklington-Lehmer test.  It's definitely an improvement over Lucas
 * since we only have to find factors up to sqrt(n), _and_ we can choose
 * a different 'a' value for each factor.  This is corollary 1 from BLS75.
 *
 * BLS T5 is the Brillhart-Lehmer-Selfridge 1975 theorem 5 (see link below).
 * We can factor even less of n, and the test lets us kick out some
 * composites early, without having to test n-3 different 'a' values.
 *
 * Once we've found the factors of n-1 (or enough of them), verification
 * usually happens really fast.  a=2 works for most, and few seem to require
 * more than ~ log2(n).  However all but BLS75 require testing all integers
 * 1 < a < n-1 before answering in the negative, which is impractical.
 *

bls75.c  view on Meta::CPAN

  mpz_sub(y, y, t);
  mpz_mul(y, y, F2);
  mpz_add_ui(y, y, 1);   /* y = 2F2^2 + (B2 - r)F2 + 1 */
  mpz_mul_ui(t, F2, B2);
  mpz_sub_ui(t, t, 1);
  mpz_mul(y, y, t);      /* times (B2F2-1) */

  return (mpz_cmp(n, y) < 0) ? 1 : 0;
}

static int bls_corollary11_limit(const mpz_t n, mpz_t R1, mpz_t F1, mpz_t F2, UV B,
                                 mpz_t t, mpz_t g, mpz_t r, mpz_t s)
{
  if (mpz_cmp(F1,F2) >= 0) {
    mpz_tdiv_q_2exp(t, F2, 1);
    mpz_mul(t, t, F1);
    mpz_mul(t, t, F1);
  } else {
    mpz_tdiv_q_2exp(t, F1, 1);
    mpz_mul(t, t, F2);
    mpz_mul(t, t, F2);

bls75.c  view on Meta::CPAN

  mpz_mul(t, t, g);
  return (mpz_cmp(n, t) < 0);
}

static int bls_theorem20_limit(const mpz_t n, mpz_t R1, mpz_t F1, mpz_t F2,
                               UV B, UV m,
                               mpz_t t, mpz_t g, mpz_t r, mpz_t s)
{
  int m_used = 0;

  if (bls_corollary11_limit(n,R1,F1,F2,B,t,g,r,s)) {
    mpz_set_ui(s,0);   /* No test for m needed */
    return 1;
  }

  mpz_mul_ui(t, F1, B);
  mpz_add_ui(g, t, 1);

  mpz_mul_ui(t, F2, B);
  mpz_sub_ui(t, t, 1);

bls75.c  view on Meta::CPAN

    if (success > 0) {
      mpz_mul(t, r, r);
      mpz_submul_ui(t, s, 8);   /* t = r^2 - 8s */
      /* N is prime if and only if s=0 OR t not a perfect square */
      success = (mpz_sgn(s) == 0 || !mpz_perfect_square_p(t))  ?  1  :  -1;
    }
    goto end_hybrid;    /* Theorem 5 or 7 */
  }

  /* Rather than looking at theorem 19 for an N+1 proof, we'll just go to
   * theorem 20 (or corollary 11).  T20 is faster. */

  /* Check N < B^3 F1*F2*F2/2  or  N < B^3 F1*F1*F2/2 */
  success = bls_theorem20_limit(n, R1, F1, F2, B1, m, t, u, r, s);

  if (get_verbose_level() > 0) printf("BLS75 proof using N-1 / N+1 (T20)\n");

  /* Trim some factors from f2stack if possible */
  if (nstack(&f2stack) > 1) {
    int i;
    mpz_set_ui(F2, 1);