xref: /netbsd-src/external/lgpl3/gmp/dist/tests/rand/spect.c (revision 72c7faa4dbb41dbb0238d6b4a109da0d4b236dd4)
1 /* spect.c -- the spectral test */
2 
3 /*
4 Copyright 1999 Free Software Foundation, Inc.
5 
6 This file is part of the GNU MP Library test suite.
7 
8 The GNU MP Library test suite is free software; you can redistribute it
9 and/or modify it under the terms of the GNU General Public License as
10 published by the Free Software Foundation; either version 3 of the License,
11 or (at your option) any later version.
12 
13 The GNU MP Library test suite is distributed in the hope that it will be
14 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General
16 Public License for more details.
17 
18 You should have received a copy of the GNU General Public License along with
19 the GNU MP Library test suite.  If not, see https://www.gnu.org/licenses/.  */
20 
21 /* T is upper dimension.  Z_A is the LC multiplier, which is
22    relatively prime to Z_M, the LC modulus.  The result is put in
23    rop[] with v[t] in rop[t-2]. */
24 
25 /* BUGS: Due to lazy allocation scheme, maximum T is hard coded to MAXT. */
26 
27 #include <stdio.h>
28 #include <stdlib.h>
29 #include <unistd.h>
30 #include <math.h>
31 
32 #include "gmpstat.h"
33 
34 int g_debug = 0;
35 
36 int
main(int argc,char * argv[])37 main (int argc, char *argv[])
38 {
39   const char usage[] = "usage: spect [-d] a m n\n";
40   int c;
41   unsigned int n;
42   mpz_t a, m;
43   mpf_t res[GMP_SPECT_MAXT], res_min[GMP_SPECT_MAXT], f_tmp;
44   int f;
45 
46 
47   mpz_init (a);
48   mpz_init (m);
49   for (f = 0; f < GMP_SPECT_MAXT; f++)
50     {
51       mpf_init (res[f]);
52       mpf_init (res_min[f]);
53     }
54   mpf_init (f_tmp);
55   mpf_set_ui (res_min[0], 32768); /* 2^15 */
56   mpf_set_ui (res_min[1], 1024); /* 2^10 */
57   mpf_set_ui (res_min[2], 256); /* 2^8 */
58   mpf_set_ui (res_min[3], 64); /* 2^6 */
59   mpf_set_ui (res_min[4], 32); /* 2^5 */
60 
61   while ((c = getopt (argc, argv, "dh")) != -1)
62     switch (c)
63       {
64       case 'd':			/* debug */
65 	g_debug++;
66 	break;
67       case 'h':
68       default:
69 	fputs (usage, stderr);
70 	exit (1);
71       }
72   argc -= optind;
73   argv += optind;
74 
75   if (argc < 3)
76     {
77       fputs (usage, stderr);
78       exit (1);
79     }
80 
81   mpz_set_str (a, argv[0], 0);
82   mpz_set_str (m, argv[1], 0);
83   n = (unsigned int) atoi (argv[2]);
84   if (n + 1 > GMP_SPECT_MAXT)
85     n = GMP_SPECT_MAXT + 1;
86 
87   spectral_test (res, n, a, m);
88 
89   for (f = 0; f < n - 1; f++)
90     {
91       /* print v */
92       printf ("%d: v = ", f + 2);
93       mpf_out_str (stdout, 10, 4, res[f]);
94 
95 #ifdef PRINT_RAISED_BY_TWO_AS_WELL
96       printf (" (^2 = ");
97       mpf_mul (f_tmp, res[f], res[f]);
98       mpf_out_str (stdout, 10, 4, f_tmp);
99       printf (")");
100 #endif /* PRINT_RAISED_BY_TWO_AS_WELL */
101 
102       /* print merit */
103       printf (" m = ");
104       merit (f_tmp, f + 2, res[f], m);
105       mpf_out_str (stdout, 10, 4, f_tmp);
106 
107       if (mpf_cmp (res[f], res_min[f]) < 0)
108 	printf ("\t*** v too low ***");
109       if (mpf_get_d (f_tmp) < .1)
110 	printf ("\t*** merit too low ***");
111 
112       puts ("");
113     }
114 
115   mpz_clear (a);
116   mpz_clear (m);
117   for (f = 0; f < GMP_SPECT_MAXT; f++)
118     {
119       mpf_clear (res[f]);
120       mpf_clear (res_min[f]);
121     }
122   mpf_clear (f_tmp);
123 
124   return 0;
125 }
126 
127 
128 void
debug_foo()129 debug_foo()
130 {
131   if (0)
132     {
133       mpz_dump (0);
134       mpf_dump (0);
135     }
136 }
137