Annotation of OpenXM_contrib/gmp/tests/rand/findlc.c, Revision 1.1
1.1 ! maekawa 1: /*
! 2: Copyright (C) 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 2.1 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; see the file COPYING.LIB. If not, write to
! 18: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
! 19: MA 02111-1307, USA.
! 20: */
! 21:
! 22: #include <stdio.h>
! 23: #include <stdlib.h>
! 24: #include <unistd.h>
! 25: #include <signal.h>
! 26: #include <math.h>
! 27: #include "gmp.h"
! 28: #include "gmpstat.h"
! 29:
! 30: #define RCSID(msg) \
! 31: static /**/const char *const rcsid[] = { (char *)rcsid, "\100(#)" msg }
! 32:
! 33: RCSID("$Id: findlc.c,v 1.3 2000/07/24 17:04:39 tege Exp $");
! 34:
! 35: int g_debug = 0;
! 36:
! 37: static mpz_t a;
! 38:
! 39: static void
! 40: sh_status (int sig)
! 41: {
! 42: printf ("sh_status: signal %d caught. dumping status.\n", sig);
! 43:
! 44: printf (" a = ");
! 45: mpz_out_str (stdout, 10, a);
! 46: printf ("\n");
! 47: fflush (stdout);
! 48:
! 49: if (SIGSEGV == sig) /* remove SEGV handler */
! 50: signal (SIGSEGV, SIG_DFL);
! 51: }
! 52:
! 53: /* Input is a modulus (m). We shall find multiplyer (a) and adder (c)
! 54: conforming to the rules found in the first comment block in file
! 55: mpz/urandom.c.
! 56:
! 57: Then run a spectral test on the generator and discard any
! 58: multipliers not passing. */
! 59:
! 60: /* TODO:
! 61:
! 62: . find a better algorithm than a+=8; bigger jumps perhaps?
! 63:
! 64: */
! 65:
! 66: void
! 67: mpz_true_random (mpz_t s, unsigned long int nbits)
! 68: {
! 69: #if __FreeBSD__
! 70: FILE *fs;
! 71: char c[1];
! 72: int i;
! 73:
! 74: mpz_set_ui (s, 0);
! 75: for (i = 0; i < nbits; i += 8)
! 76: {
! 77: for (;;)
! 78: {
! 79: int nread;
! 80: fs = fopen ("/dev/random", "r");
! 81: nread = fread (c, 1, 1, fs);
! 82: fclose (fs);
! 83: if (nread != 0)
! 84: break;
! 85: sleep (1);
! 86: }
! 87: mpz_mul_2exp (s, s, 8);
! 88: mpz_add_ui (s, s, ((unsigned long int) c[0]) & 0xff);
! 89: printf ("%d random bits\n", i + 8);
! 90: }
! 91: if (nbits % 8 != 0)
! 92: mpz_mod_2exp (s, s, nbits);
! 93: #endif
! 94: }
! 95:
! 96: int
! 97: main (int argc, char *argv[])
! 98: {
! 99: const char usage[] = "usage: findlc [-dv] m2exp [low_merit [high_merit]]\n";
! 100: int f;
! 101: int v_lose, m_lose, v_best, m_best;
! 102: int c;
! 103: int debug = 1;
! 104: int cnt_high_merit;
! 105: mpz_t m;
! 106: unsigned long int m2exp;
! 107: #define DIMS 6 /* dimensions run in spectral test */
! 108: mpf_t v[DIMS-1]; /* spectral test result (there's no v
! 109: for 1st dimension */
! 110: mpf_t f_merit, low_merit, high_merit;
! 111: mpz_t acc, minus8;
! 112: mpz_t min, max;
! 113: mpz_t s;
! 114:
! 115:
! 116: mpz_init (m);
! 117: mpz_init (a);
! 118: for (f = 0; f < DIMS-1; f++)
! 119: mpf_init (v[f]);
! 120: mpf_init (f_merit);
! 121: mpf_init_set_d (low_merit, .1);
! 122: mpf_init_set_d (high_merit, .1);
! 123:
! 124: while ((c = getopt (argc, argv, "a:di:hv")) != -1)
! 125: switch (c)
! 126: {
! 127: case 'd': /* debug */
! 128: g_debug++;
! 129: break;
! 130:
! 131: case 'v': /* print version */
! 132: puts (rcsid[1]);
! 133: exit (0);
! 134:
! 135: case 'h':
! 136: case '?':
! 137: default:
! 138: fputs (usage, stderr);
! 139: exit (1);
! 140: }
! 141:
! 142: argc -= optind;
! 143: argv += optind;
! 144:
! 145: if (argc < 1)
! 146: {
! 147: fputs (usage, stderr);
! 148: exit (1);
! 149: }
! 150:
! 151: /* Install signal handler. */
! 152: if (SIG_ERR == signal (SIGSEGV, sh_status))
! 153: {
! 154: perror ("signal (SIGSEGV)");
! 155: exit (1);
! 156: }
! 157: if (SIG_ERR == signal (SIGHUP, sh_status))
! 158: {
! 159: perror ("signal (SIGHUP)");
! 160: exit (1);
! 161: }
! 162:
! 163: printf ("findlc: version: %s\n", rcsid[1]);
! 164: m2exp = atol (argv[0]);
! 165: mpz_init_set_ui (m, 1);
! 166: mpz_mul_2exp (m, m, m2exp);
! 167: printf ("m = 0x");
! 168: mpz_out_str (stdout, 16, m);
! 169: puts ("");
! 170:
! 171: if (argc > 1) /* have low_merit */
! 172: mpf_set_str (low_merit, argv[1], 0);
! 173: if (argc > 2) /* have high_merit */
! 174: mpf_set_str (high_merit, argv[2], 0);
! 175:
! 176: if (debug)
! 177: {
! 178: fprintf (stderr, "low_merit = ");
! 179: mpf_out_str (stderr, 10, 2, low_merit);
! 180: fprintf (stderr, "; high_merit = ");
! 181: mpf_out_str (stderr, 10, 2, high_merit);
! 182: fputs ("\n", stderr);
! 183: }
! 184:
! 185: mpz_init (minus8);
! 186: mpz_set_si (minus8, -8L);
! 187: mpz_init_set_ui (acc, 0);
! 188: mpz_init (s);
! 189: mpz_init_set_d (min, 0.01 * pow (2.0, (double) m2exp));
! 190: mpz_init_set_d (max, 0.99 * pow (2.0, (double) m2exp));
! 191:
! 192: mpz_true_random (s, m2exp); /* Start. */
! 193: mpz_setbit (s, 0); /* Make it odd. */
! 194:
! 195: v_best = m_best = 2*(DIMS-1);
! 196: for (;;)
! 197: {
! 198: mpz_add (acc, acc, s);
! 199: mpz_mod_2exp (acc, acc, m2exp);
! 200: #if later
! 201: mpz_and_si (a, acc, -8L);
! 202: #else
! 203: mpz_and (a, acc, minus8);
! 204: #endif
! 205: mpz_add_ui (a, a, 5);
! 206: if (mpz_cmp (a, min) <= 0 || mpz_cmp (a, max) >= 0)
! 207: continue;
! 208:
! 209: spectral_test (v, DIMS, a, m);
! 210: for (f = 0, v_lose = m_lose = 0, cnt_high_merit = DIMS-1;
! 211: f < DIMS-1; f++)
! 212: {
! 213: merit (f_merit, f + 2, v[f], m);
! 214:
! 215: if (mpf_cmp_ui (v[f], 1 << (30 / (f + 2) + (f == 2))) < 0)
! 216: v_lose++;
! 217:
! 218: if (mpf_cmp (f_merit, low_merit) < 0)
! 219: m_lose++;
! 220:
! 221: if (mpf_cmp (f_merit, high_merit) >= 0)
! 222: cnt_high_merit--;
! 223: }
! 224:
! 225: if (0 == v_lose && 0 == m_lose)
! 226: {
! 227: mpz_out_str (stdout, 10, a); puts (""); fflush (stdout);
! 228: if (0 == cnt_high_merit)
! 229: break; /* leave loop */
! 230: }
! 231: if (v_lose < v_best)
! 232: {
! 233: v_best = v_lose;
! 234: printf ("best (v_lose=%d; m_lose=%d): ", v_lose, m_lose);
! 235: mpz_out_str (stdout, 10, a); puts (""); fflush (stdout);
! 236: }
! 237: if (m_lose < m_best)
! 238: {
! 239: m_best = m_lose;
! 240: printf ("best (v_lose=%d; m_lose=%d): ", v_lose, m_lose);
! 241: mpz_out_str (stdout, 10, a); puts (""); fflush (stdout);
! 242: }
! 243: }
! 244:
! 245: mpz_clear (m);
! 246: mpz_clear (a);
! 247: for (f = 0; f < DIMS-1; f++)
! 248: mpf_clear (v[f]);
! 249: mpf_clear (f_merit);
! 250: mpf_clear (low_merit);
! 251: mpf_clear (high_merit);
! 252:
! 253: printf ("done.\n");
! 254: return 0;
! 255: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to findlc.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / gmp / tests / rand