1 /*
2
3 Copyright 2012, 2013, 2018 Free Software Foundation, Inc.
4
5 This file is part of the GNU MP Library test suite.
6
7 The GNU MP Library test suite is free software; you can redistribute it
8 and/or modify it under the terms of the GNU General Public License as
9 published by the Free Software Foundation; either version 3 of the License,
10 or (at your option) any later version.
11
12 The GNU MP Library test suite is distributed in the hope that it will be
13 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
15 Public License for more details.
16
17 You should have received a copy of the GNU General Public License along with
18 the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */
19
20 #include <limits.h>
21 #include <math.h>
22 #include <stdlib.h>
23 #include <stdio.h>
24 #include <string.h>
25 #include <float.h>
26
27 #include "testutils.h"
28
29 #define GMP_LIMB_BITS (sizeof(mp_limb_t) * CHAR_BIT)
30
31 mp_bitcnt_t
32 mpz_mantissasizeinbits (const mpz_t z)
33 {
34 return ! mpz_cmp_ui (z, 0) ? 0 :
35 mpz_sizeinbase (z, 2) - mpz_scan1 (z, 0);
36 }
37
38 #if defined(DBL_MANT_DIG) && FLT_RADIX == 2
39 int
40 mpz_get_d_exact_p (const mpz_t z)
41 {
42 return mpz_mantissasizeinbits (z) <= DBL_MANT_DIG;
43 }
44 #define HAVE_EXACT_P 1
45 #endif
46
47 #define COUNT 10000
48
49 void
50 test_matissa (void)
51 {
52 mpz_t x, y;
53 int i, c;
54
55 mpz_init (x);
56 mpz_init (y);
57
58 mini_urandomb (y, 4);
59 c = i = mpz_get_ui (y);
60
61 do {
62 double d;
63 int cmp;
64
65 mpz_setbit (x, c);
66 d = mpz_get_d (x);
67 mpz_set_d (y, d);
68 if (mpz_cmp_d (y, d) != 0)
69 {
70 fprintf (stderr, "mpz_cmp_d (y, d) failed:\n"
71 "d = %.20g\n"
72 "i = %i\n"
73 "c = %i\n",
74 d, i, c);
75 abort ();
76 }
77
78 cmp = mpz_cmp (x, y);
79
80 #if defined(HAVE_EXACT_P)
81 if ((mpz_get_d_exact_p (x) != 0) != (cmp == 0))
82 {
83 fprintf (stderr, "Not all bits converted:\n"
84 "d = %.20g\n"
85 "i = %i\n"
86 "c = %i\n",
87 d, i, c);
88 abort ();
89 }
90 #endif
91
92 if (cmp < 0)
93 {
94 fprintf (stderr, "mpz_get_d failed:\n"
95 "d = %.20g\n"
96 "i = %i\n"
97 "c = %i\n",
98 d, i, c);
99 abort ();
100 }
101 else if (cmp > 0)
102 {
103 if (mpz_cmp_d (x, d) <= 0)
104 {
105 fprintf (stderr, "mpz_cmp_d (x, d) failed:\n"
106 "d = %.20g\n"
107 "i = %i\n"
108 "c = %i\n",
109 d, i, c);
110 abort ();
111 }
112 break;
113 }
114 ++c;
115 } while (1);
116
117 mpz_clear (x);
118 mpz_clear (y);
119 }
120
121 #ifndef M_PI
122 #define M_PI 3.141592653589793238462643383279502884
123 #endif
124
125 static const struct
126 {
127 double d;
128 const char *s;
129 } values[] = {
130 { 0.0, "0" },
131 { 0.3, "0" },
132 { -0.3, "0" },
133 { M_PI, "3" },
134 { M_PI*1e15, "b29430a256d21" },
135 { -M_PI*1e15, "-b29430a256d21" },
136 /* 17 * 2^{200} =
137 27317946752402834684213355569799764242877450894307478200123392 */
138 {0.2731794675240283468421335556979976424288e62,
139 "1100000000000000000000000000000000000000000000000000" },
140 { 0.0, NULL }
141 };
142
143 void
144 testmain (int argc, char **argv)
145 {
146 unsigned i;
147 mpz_t x;
148
149 for (i = 0; values[i].s; i++)
150 {
151 char *s;
152 mpz_init_set_d (x, values[i].d);
153 s = mpz_get_str (NULL, 16, x);
154 if (strcmp (s, values[i].s) != 0)
155 {
156 fprintf (stderr, "mpz_set_d failed:\n"
157 "d = %.20g\n"
158 "s = %s\n"
159 "r = %s\n",
160 values[i].d, s, values[i].s);
161 abort ();
162 }
163 testfree (s, strlen(s) + 1);
164 mpz_clear (x);
165 }
166
167 mpz_init (x);
168
169 for (i = 0; i < COUNT; i++)
170 {
171 /* Use volatile, to avoid extended precision in floating point
172 registers, e.g., on m68k and 80387. */
173 volatile double d, f;
174 unsigned long m;
175 int e;
176
177 mini_rrandomb (x, GMP_LIMB_BITS);
178 m = mpz_get_ui (x);
179 mini_urandomb (x, 8);
180 e = mpz_get_ui (x) - 100;
181
182 d = ldexp ((double) m, e);
183 mpz_set_d (x, d);
184 f = mpz_get_d (x);
185 if (f != floor (d))
186 {
187 fprintf (stderr, "mpz_set_d/mpz_get_d failed:\n");
188 goto dumperror;
189 }
190 if ((f == d) ? (mpz_cmp_d (x, d) != 0) : (mpz_cmp_d (x, d) >= 0))
191 {
192 fprintf (stderr, "mpz_cmp_d (x, d) failed:\n");
193 goto dumperror;
194 }
195 f = d + 1.0;
196 if (f > d && ! (mpz_cmp_d (x, f) < 0))
197 {
198 fprintf (stderr, "mpz_cmp_d (x, f) failed:\n");
199 goto dumperror;
200 }
201
202 d = - d;
203
204 mpz_set_d (x, d);
205 f = mpz_get_d (x);
206 if (f != ceil (d))
207 {
208 fprintf (stderr, "mpz_set_d/mpz_get_d failed:\n");
209 dumperror:
210 dump ("x", x);
211 fprintf (stderr, "m = %lx, e = %i\n", m, e);
212 fprintf (stderr, "d = %.15g\n", d);
213 fprintf (stderr, "f = %.15g\n", f);
214 fprintf (stderr, "f - d = %.5g\n", f - d);
215 abort ();
216 }
217 if ((f == d) ? (mpz_cmp_d (x, d) != 0) : (mpz_cmp_d (x, d) <= 0))
218 {
219 fprintf (stderr, "mpz_cmp_d (x, d) failed:\n");
220 goto dumperror;
221 }
222 f = d - 1.0;
223 if (f < d && ! (mpz_cmp_d (x, f) > 0))
224 {
225 fprintf (stderr, "mpz_cmp_d (x, f) failed:\n");
226 goto dumperror;
227 }
228 }
229
230 mpz_clear (x);
231 test_matissa();
232 }