| OLD | NEW |
| (Empty) |
| 1 /* | |
| 2 Copyright 2000 Free Software Foundation, Inc. | |
| 3 | |
| 4 This file is part of the GNU MP Library. | |
| 5 | |
| 6 The GNU MP Library is free software; you can redistribute it and/or modify | |
| 7 it under the terms of the GNU Lesser General Public License as published by | |
| 8 the Free Software Foundation; either version 3 of the License, or (at your | |
| 9 option) any later version. | |
| 10 | |
| 11 The GNU MP Library is distributed in the hope that it will be useful, but | |
| 12 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY | |
| 13 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public | |
| 14 License for more details. | |
| 15 | |
| 16 You should have received a copy of the GNU Lesser General Public License | |
| 17 along with the GNU MP Library. If not, see http://www.gnu.org/licenses/. */ | |
| 18 | |
| 19 #include <stdio.h> | |
| 20 #include <stdlib.h> | |
| 21 #include <unistd.h> | |
| 22 #include <signal.h> | |
| 23 #include <math.h> | |
| 24 #include "gmp.h" | |
| 25 #include "gmpstat.h" | |
| 26 | |
| 27 #define RCSID(msg) \ | |
| 28 static /**/const char *const rcsid[] = { (char *)rcsid, "\100(#)" msg } | |
| 29 | |
| 30 RCSID("$Id$"); | |
| 31 | |
| 32 int g_debug = 0; | |
| 33 | |
| 34 static mpz_t a; | |
| 35 | |
| 36 static void | |
| 37 sh_status (int sig) | |
| 38 { | |
| 39 printf ("sh_status: signal %d caught. dumping status.\n", sig); | |
| 40 | |
| 41 printf (" a = "); | |
| 42 mpz_out_str (stdout, 10, a); | |
| 43 printf ("\n"); | |
| 44 fflush (stdout); | |
| 45 | |
| 46 if (SIGSEGV == sig) /* remove SEGV handler */ | |
| 47 signal (SIGSEGV, SIG_DFL); | |
| 48 } | |
| 49 | |
| 50 /* Input is a modulus (m). We shall find multiplyer (a) and adder (c) | |
| 51 conforming to the rules found in the first comment block in file | |
| 52 mpz/urandom.c. | |
| 53 | |
| 54 Then run a spectral test on the generator and discard any | |
| 55 multipliers not passing. */ | |
| 56 | |
| 57 /* TODO: | |
| 58 | |
| 59 . find a better algorithm than a+=8; bigger jumps perhaps? | |
| 60 | |
| 61 */ | |
| 62 | |
| 63 void | |
| 64 mpz_true_random (mpz_t s, unsigned long int nbits) | |
| 65 { | |
| 66 #if __FreeBSD__ | |
| 67 FILE *fs; | |
| 68 char c[1]; | |
| 69 int i; | |
| 70 | |
| 71 mpz_set_ui (s, 0); | |
| 72 for (i = 0; i < nbits; i += 8) | |
| 73 { | |
| 74 for (;;) | |
| 75 { | |
| 76 int nread; | |
| 77 fs = fopen ("/dev/random", "r"); | |
| 78 nread = fread (c, 1, 1, fs); | |
| 79 fclose (fs); | |
| 80 if (nread != 0) | |
| 81 break; | |
| 82 sleep (1); | |
| 83 } | |
| 84 mpz_mul_2exp (s, s, 8); | |
| 85 mpz_add_ui (s, s, ((unsigned long int) c[0]) & 0xff); | |
| 86 printf ("%d random bits\n", i + 8); | |
| 87 } | |
| 88 if (nbits % 8 != 0) | |
| 89 mpz_mod_2exp (s, s, nbits); | |
| 90 #endif | |
| 91 } | |
| 92 | |
| 93 int | |
| 94 main (int argc, char *argv[]) | |
| 95 { | |
| 96 const char usage[] = "usage: findlc [-dv] m2exp [low_merit [high_merit]]\n"; | |
| 97 int f; | |
| 98 int v_lose, m_lose, v_best, m_best; | |
| 99 int c; | |
| 100 int debug = 1; | |
| 101 int cnt_high_merit; | |
| 102 mpz_t m; | |
| 103 unsigned long int m2exp; | |
| 104 #define DIMS 6 /* dimensions run in spectral test */ | |
| 105 mpf_t v[DIMS-1]; /* spectral test result (there's no v | |
| 106 for 1st dimension */ | |
| 107 mpf_t f_merit, low_merit, high_merit; | |
| 108 mpz_t acc, minus8; | |
| 109 mpz_t min, max; | |
| 110 mpz_t s; | |
| 111 | |
| 112 | |
| 113 mpz_init (m); | |
| 114 mpz_init (a); | |
| 115 for (f = 0; f < DIMS-1; f++) | |
| 116 mpf_init (v[f]); | |
| 117 mpf_init (f_merit); | |
| 118 mpf_init_set_d (low_merit, .1); | |
| 119 mpf_init_set_d (high_merit, .1); | |
| 120 | |
| 121 while ((c = getopt (argc, argv, "a:di:hv")) != -1) | |
| 122 switch (c) | |
| 123 { | |
| 124 case 'd': /* debug */ | |
| 125 g_debug++; | |
| 126 break; | |
| 127 | |
| 128 case 'v': /* print version */ | |
| 129 puts (rcsid[1]); | |
| 130 exit (0); | |
| 131 | |
| 132 case 'h': | |
| 133 case '?': | |
| 134 default: | |
| 135 fputs (usage, stderr); | |
| 136 exit (1); | |
| 137 } | |
| 138 | |
| 139 argc -= optind; | |
| 140 argv += optind; | |
| 141 | |
| 142 if (argc < 1) | |
| 143 { | |
| 144 fputs (usage, stderr); | |
| 145 exit (1); | |
| 146 } | |
| 147 | |
| 148 /* Install signal handler. */ | |
| 149 if (SIG_ERR == signal (SIGSEGV, sh_status)) | |
| 150 { | |
| 151 perror ("signal (SIGSEGV)"); | |
| 152 exit (1); | |
| 153 } | |
| 154 if (SIG_ERR == signal (SIGHUP, sh_status)) | |
| 155 { | |
| 156 perror ("signal (SIGHUP)"); | |
| 157 exit (1); | |
| 158 } | |
| 159 | |
| 160 printf ("findlc: version: %s\n", rcsid[1]); | |
| 161 m2exp = atol (argv[0]); | |
| 162 mpz_init_set_ui (m, 1); | |
| 163 mpz_mul_2exp (m, m, m2exp); | |
| 164 printf ("m = 0x"); | |
| 165 mpz_out_str (stdout, 16, m); | |
| 166 puts (""); | |
| 167 | |
| 168 if (argc > 1) /* have low_merit */ | |
| 169 mpf_set_str (low_merit, argv[1], 0); | |
| 170 if (argc > 2) /* have high_merit */ | |
| 171 mpf_set_str (high_merit, argv[2], 0); | |
| 172 | |
| 173 if (debug) | |
| 174 { | |
| 175 fprintf (stderr, "low_merit = "); | |
| 176 mpf_out_str (stderr, 10, 2, low_merit); | |
| 177 fprintf (stderr, "; high_merit = "); | |
| 178 mpf_out_str (stderr, 10, 2, high_merit); | |
| 179 fputs ("\n", stderr); | |
| 180 } | |
| 181 | |
| 182 mpz_init (minus8); | |
| 183 mpz_set_si (minus8, -8L); | |
| 184 mpz_init_set_ui (acc, 0); | |
| 185 mpz_init (s); | |
| 186 mpz_init_set_d (min, 0.01 * pow (2.0, (double) m2exp)); | |
| 187 mpz_init_set_d (max, 0.99 * pow (2.0, (double) m2exp)); | |
| 188 | |
| 189 mpz_true_random (s, m2exp); /* Start. */ | |
| 190 mpz_setbit (s, 0); /* Make it odd. */ | |
| 191 | |
| 192 v_best = m_best = 2*(DIMS-1); | |
| 193 for (;;) | |
| 194 { | |
| 195 mpz_add (acc, acc, s); | |
| 196 mpz_mod_2exp (acc, acc, m2exp); | |
| 197 #if later | |
| 198 mpz_and_si (a, acc, -8L); | |
| 199 #else | |
| 200 mpz_and (a, acc, minus8); | |
| 201 #endif | |
| 202 mpz_add_ui (a, a, 5); | |
| 203 if (mpz_cmp (a, min) <= 0 || mpz_cmp (a, max) >= 0) | |
| 204 continue; | |
| 205 | |
| 206 spectral_test (v, DIMS, a, m); | |
| 207 for (f = 0, v_lose = m_lose = 0, cnt_high_merit = DIMS-1; | |
| 208 f < DIMS-1; f++) | |
| 209 { | |
| 210 merit (f_merit, f + 2, v[f], m); | |
| 211 | |
| 212 if (mpf_cmp_ui (v[f], 1 << (30 / (f + 2) + (f == 2))) < 0) | |
| 213 v_lose++; | |
| 214 | |
| 215 if (mpf_cmp (f_merit, low_merit) < 0) | |
| 216 m_lose++; | |
| 217 | |
| 218 if (mpf_cmp (f_merit, high_merit) >= 0) | |
| 219 cnt_high_merit--; | |
| 220 } | |
| 221 | |
| 222 if (0 == v_lose && 0 == m_lose) | |
| 223 { | |
| 224 mpz_out_str (stdout, 10, a); puts (""); fflush (stdout); | |
| 225 if (0 == cnt_high_merit) | |
| 226 break; /* leave loop */ | |
| 227 } | |
| 228 if (v_lose < v_best) | |
| 229 { | |
| 230 v_best = v_lose; | |
| 231 printf ("best (v_lose=%d; m_lose=%d): ", v_lose, m_lose); | |
| 232 mpz_out_str (stdout, 10, a); puts (""); fflush (stdout); | |
| 233 } | |
| 234 if (m_lose < m_best) | |
| 235 { | |
| 236 m_best = m_lose; | |
| 237 printf ("best (v_lose=%d; m_lose=%d): ", v_lose, m_lose); | |
| 238 mpz_out_str (stdout, 10, a); puts (""); fflush (stdout); | |
| 239 } | |
| 240 } | |
| 241 | |
| 242 mpz_clear (m); | |
| 243 mpz_clear (a); | |
| 244 for (f = 0; f < DIMS-1; f++) | |
| 245 mpf_clear (v[f]); | |
| 246 mpf_clear (f_merit); | |
| 247 mpf_clear (low_merit); | |
| 248 mpf_clear (high_merit); | |
| 249 | |
| 250 printf ("done.\n"); | |
| 251 return 0; | |
| 252 } | |
| OLD | NEW |