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 |