Initial commit.

1 files changed, 194 insertions, 0 deletions

diff --git a/gmp-6.3.0/gen-sieve.c b/gmp-6.3.0/gen-sieve.c
new file mode 100644
index 0000000..7133918
--- /dev/null
+++ b/gmp-6.3.0/gen-sieve.c
@@ -0,0 +1,194 @@
+/* Generate primesieve data.
+
+   Contributed to the GNU project by Marco Bodrato.
+
+Copyright 2021, 2022 Free Software Foundation, Inc.
+
+This file is part of the GNU MP Library.
+
+The GNU MP Library is free software; you can redistribute it and/or modify
+it under the terms of either:
+
+  * the GNU Lesser General Public License as published by the Free
+    Software Foundation; either version 3 of the License, or (at your
+    option) any later version.
+
+or
+
+  * the GNU General Public License as published by the Free Software
+    Foundation; either version 2 of the License, or (at your option) any
+    later version.
+
+or both in parallel, as here.
+
+The GNU MP Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+for more details.
+
+You should have received copies of the GNU General Public License and the
+GNU Lesser General Public License along with the GNU MP Library.  If not,
+see https://www.gnu.org/licenses/.  */
+
+#include <stdio.h>
+#include "bootstrap.c"
+
+static int
+bit_to_n (int bit) { return (bit*3+4)|1; }
+
+int
+generate (int limb_bits, int limit)
+{
+  mpz_t  limb;
+  int i, lb, pc, c, totpc, maxprime;
+
+  mpz_init (limb);
+
+  printf ("/* This file generated by gen-sieve.c - DO NOT EDIT. */\n");
+  printf ("\n");
+  printf ("#if GMP_LIMB_BITS != %d\n", limb_bits);
+  printf ("Error, error, this data is for %d bits\n", limb_bits);
+  printf ("#endif\n");
+  printf ("\n");
+  printf ("#define PRIMESIEVE_INIT_TABLE ");
+
+  maxprime = 3;
+  lb = pc = c = totpc = 0;
+  for (i = 0; i < limit; i++)
+    {
+      if (! isprime (bit_to_n (i)))
+	mpz_setbit (limb, lb);
+      else
+	maxprime = bit_to_n (i), ++pc;
+      ++lb;
+      if (lb == limb_bits)
+	{
+	  ++c;
+	  printf ("\\\n\tCNST_LIMB (0x");
+	  mpz_out_str (stdout, -16, limb);
+	  printf ("),\t/* %d - %d (%d primes) */\t", bit_to_n (i + 1 - limb_bits),
+		  bit_to_n (i + 1) - 1, pc);
+	  totpc += pc;
+	  lb = pc = 0;
+	  mpz_set_ui (limb, 0);
+	}
+    }
+
+  if ((mpz_sgn (limb) | lb | pc) != 0)
+    {
+      printf ("\ngen-sieve: Internal error, during generate (%d, %d).\n", limb_bits,  limit);
+      abort();
+    }
+
+  mpz_clear (limb);
+
+  printf ("\n");
+  printf ("#define PRIMESIEVE_NUMBEROF_TABLE %d\n", c);
+
+  printf ("/* #define PRIMESIEVE_PRIMES_IN_TABLE %d */\n", totpc);
+  printf ("#define PRIMESIEVE_HIGHEST_PRIME %d\n", maxprime);
+  printf ("/* #define PRIMESIEVE_FIRST_UNCHECKED %d */\n", bit_to_n (limit));
+
+  return c;
+}
+
+void
+setmask (mpz_t mask, int a, int b)
+{
+  mpz_set_ui (mask, 0);
+  for (unsigned i = 0; i < 2 * a * b; ++i)
+    if ((bit_to_n (i) % a == 0) || (bit_to_n (i) % b == 0))
+      mpz_setbit (mask, i);
+}
+
+void
+gen_sieve_masks (int limb_bits) {
+  mpz_t  mask, limb;
+
+  mpz_init (mask);
+  mpz_init (limb);
+
+  printf ("\n");
+  if (limb_bits > 60 && limb_bits < 91)
+    {
+      setmask (mask, 5, 11);
+
+      mpz_tdiv_r_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_MASK1 CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+      mpz_tdiv_q_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_MASKT CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+
+      setmask (mask, 7, 13);
+
+      mpz_tdiv_r_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_2MSK1 CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+      mpz_tdiv_q_2exp (mask, mask, limb_bits);
+      mpz_tdiv_r_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_2MSK2 CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+      mpz_tdiv_q_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_2MSKT CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+    }
+  else if (limb_bits > 23 && limb_bits < 36)
+    {
+      setmask (mask, 5, 7);
+
+      mpz_tdiv_r_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_MASK1 CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+      mpz_tdiv_q_2exp (mask, mask, limb_bits);
+      mpz_tdiv_r_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_MASK2 CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+      mpz_tdiv_q_2exp (limb, mask, limb_bits);
+      printf ("#define SIEVE_MASKT CNST_LIMB(0x");
+      mpz_out_str (stdout, -16, limb);
+      printf (")\n");
+    }
+  printf ("\n");
+
+  mpz_clear (limb);
+  mpz_clear (mask);
+}
+
+/* 5*2 = 10
+   7*2 = 14
+   5*7*2 = 70 (2*35, 3*24, 4*18, 5*14...)
+   5*11*2 = 110 (2*55, 3*37, 4*28, 5*22...)
+   5*13*2 = 130 (2*65, 3*44, 4*33, 5*26...)
+   7*11*2 = 154 (2*77, 3*52, 4*39, 5*31...)
+   7*13*2 = 182 (2*91, 3*61, 4*46, 5*37...)
+*/
+
+int
+main (int argc, char *argv[])
+{
+  int  limb_bits, limit;
+
+  if (argc != 2)
+    {
+      fprintf (stderr, "Usage: gen-sieve <limbbits>\n");
+      exit (1);
+    }
+
+  limb_bits = atoi (argv[1]);
+
+  limit = 64 * 28; /* bits in the presieved sieve */
+  if (limit % limb_bits != 0)
+    limit += limb_bits - limit % limb_bits;
+  generate (limb_bits, limit);
+  gen_sieve_masks (limb_bits);
+
+  return 0;
+}