summaryrefslogtreecommitdiffstats
path: root/libtommath/etc/mersenne.c
diff options
context:
space:
mode:
authorjan.nijtmans <nijtmans@users.sourceforge.net>2017-04-12 12:00:11 (GMT)
committerjan.nijtmans <nijtmans@users.sourceforge.net>2017-04-12 12:00:11 (GMT)
commitece45e7fb6469e3ee3ad49f168f8711fb36f93ce (patch)
treedb4a77927de2a4d6c6cf2bc672ebda4098b9b1a0 /libtommath/etc/mersenne.c
parent6f3388528ef453d29fbddba3f5a054d2f5268207 (diff)
parent473bfc0f18451046035f638732a609fc86d5a0aa (diff)
downloadtcl-ece45e7fb6469e3ee3ad49f168f8711fb36f93ce.zip
tcl-ece45e7fb6469e3ee3ad49f168f8711fb36f93ce.tar.gz
tcl-ece45e7fb6469e3ee3ad49f168f8711fb36f93ce.tar.bz2
merge trunkinitsubsystems
Diffstat (limited to 'libtommath/etc/mersenne.c')
-rw-r--r--libtommath/etc/mersenne.c140
1 files changed, 0 insertions, 140 deletions
diff --git a/libtommath/etc/mersenne.c b/libtommath/etc/mersenne.c
deleted file mode 100644
index 28ac834..0000000
--- a/libtommath/etc/mersenne.c
+++ /dev/null
@@ -1,140 +0,0 @@
-/* Finds Mersenne primes using the Lucas-Lehmer test
- *
- * Tom St Denis, tomstdenis@gmail.com
- */
-#include <time.h>
-#include <tommath.h>
-
-int
-is_mersenne (long s, int *pp)
-{
- mp_int n, u;
- int res, k;
-
- *pp = 0;
-
- if ((res = mp_init (&n)) != MP_OKAY) {
- return res;
- }
-
- if ((res = mp_init (&u)) != MP_OKAY) {
- goto LBL_N;
- }
-
- /* n = 2^s - 1 */
- if ((res = mp_2expt(&n, s)) != MP_OKAY) {
- goto LBL_MU;
- }
- if ((res = mp_sub_d (&n, 1, &n)) != MP_OKAY) {
- goto LBL_MU;
- }
-
- /* set u=4 */
- mp_set (&u, 4);
-
- /* for k=1 to s-2 do */
- for (k = 1; k <= s - 2; k++) {
- /* u = u^2 - 2 mod n */
- if ((res = mp_sqr (&u, &u)) != MP_OKAY) {
- goto LBL_MU;
- }
- if ((res = mp_sub_d (&u, 2, &u)) != MP_OKAY) {
- goto LBL_MU;
- }
-
- /* make sure u is positive */
- while (u.sign == MP_NEG) {
- if ((res = mp_add (&u, &n, &u)) != MP_OKAY) {
- goto LBL_MU;
- }
- }
-
- /* reduce */
- if ((res = mp_reduce_2k (&u, &n, 1)) != MP_OKAY) {
- goto LBL_MU;
- }
- }
-
- /* if u == 0 then its prime */
- if (mp_iszero (&u) == 1) {
- mp_prime_is_prime(&n, 8, pp);
- if (*pp != 1) printf("FAILURE\n");
- }
-
- res = MP_OKAY;
-LBL_MU:mp_clear (&u);
-LBL_N:mp_clear (&n);
- return res;
-}
-
-/* square root of a long < 65536 */
-long
-i_sqrt (long x)
-{
- long x1, x2;
-
- x2 = 16;
- do {
- x1 = x2;
- x2 = x1 - ((x1 * x1) - x) / (2 * x1);
- } while (x1 != x2);
-
- if (x1 * x1 > x) {
- --x1;
- }
-
- return x1;
-}
-
-/* is the long prime by brute force */
-int
-isprime (long k)
-{
- long y, z;
-
- y = i_sqrt (k);
- for (z = 2; z <= y; z++) {
- if ((k % z) == 0)
- return 0;
- }
- return 1;
-}
-
-
-int
-main (void)
-{
- int pp;
- long k;
- clock_t tt;
-
- k = 3;
-
- for (;;) {
- /* start time */
- tt = clock ();
-
- /* test if 2^k - 1 is prime */
- if (is_mersenne (k, &pp) != MP_OKAY) {
- printf ("Whoa error\n");
- return -1;
- }
-
- if (pp == 1) {
- /* count time */
- tt = clock () - tt;
-
- /* display if prime */
- printf ("2^%-5ld - 1 is prime, test took %ld ticks\n", k, tt);
- }
-
- /* goto next odd exponent */
- k += 2;
-
- /* but make sure its prime */
- while (isprime (k) == 0) {
- k += 2;
- }
- }
- return 0;
-}