tmul.c revision 1.1.1.4 1 1.1 mrg /* Test file for mpfr_mul.
2 1.1 mrg
3 1.1.1.4 mrg Copyright 1999, 2001-2018 Free Software Foundation, Inc.
4 1.1.1.3 mrg Contributed by the AriC and Caramba projects, INRIA.
5 1.1 mrg
6 1.1 mrg This file is part of the GNU MPFR Library.
7 1.1 mrg
8 1.1 mrg The GNU MPFR Library is free software; you can redistribute it and/or modify
9 1.1 mrg it under the terms of the GNU Lesser General Public License as published by
10 1.1 mrg the Free Software Foundation; either version 3 of the License, or (at your
11 1.1 mrg option) any later version.
12 1.1 mrg
13 1.1 mrg The GNU MPFR Library is distributed in the hope that it will be useful, but
14 1.1 mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 1.1 mrg or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 1.1 mrg License for more details.
17 1.1 mrg
18 1.1 mrg You should have received a copy of the GNU Lesser General Public License
19 1.1 mrg along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
20 1.1 mrg http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 1.1 mrg 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
22 1.1 mrg
23 1.1 mrg #include "mpfr-test.h"
24 1.1 mrg
25 1.1 mrg #ifdef CHECK_EXTERNAL
26 1.1 mrg static int
27 1.1 mrg test_mul (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode)
28 1.1 mrg {
29 1.1 mrg int res;
30 1.1 mrg int ok = rnd_mode == MPFR_RNDN && mpfr_number_p (b) && mpfr_number_p (c);
31 1.1 mrg if (ok)
32 1.1 mrg {
33 1.1 mrg mpfr_print_raw (b);
34 1.1 mrg printf (" ");
35 1.1 mrg mpfr_print_raw (c);
36 1.1 mrg }
37 1.1 mrg res = mpfr_mul (a, b, c, rnd_mode);
38 1.1 mrg if (ok)
39 1.1 mrg {
40 1.1 mrg printf (" ");
41 1.1 mrg mpfr_print_raw (a);
42 1.1 mrg printf ("\n");
43 1.1 mrg }
44 1.1 mrg return res;
45 1.1 mrg }
46 1.1 mrg #else
47 1.1 mrg #define test_mul mpfr_mul
48 1.1 mrg #endif
49 1.1 mrg
50 1.1 mrg /* checks that xs * ys gives the expected result res */
51 1.1 mrg static void
52 1.1 mrg check (const char *xs, const char *ys, mpfr_rnd_t rnd_mode,
53 1.1 mrg unsigned int px, unsigned int py, unsigned int pz, const char *res)
54 1.1 mrg {
55 1.1 mrg mpfr_t xx, yy, zz;
56 1.1 mrg
57 1.1 mrg mpfr_init2 (xx, px);
58 1.1 mrg mpfr_init2 (yy, py);
59 1.1 mrg mpfr_init2 (zz, pz);
60 1.1 mrg mpfr_set_str1 (xx, xs);
61 1.1 mrg mpfr_set_str1 (yy, ys);
62 1.1 mrg test_mul(zz, xx, yy, rnd_mode);
63 1.1 mrg if (mpfr_cmp_str1 (zz, res) )
64 1.1 mrg {
65 1.1.1.4 mrg printf ("(1) mpfr_mul failed for x=%s y=%s with rnd=%s\n",
66 1.1 mrg xs, ys, mpfr_print_rnd_mode (rnd_mode));
67 1.1 mrg printf ("correct is %s, mpfr_mul gives ", res);
68 1.1.1.4 mrg mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN);
69 1.1.1.4 mrg putchar ('\n');
70 1.1 mrg exit (1);
71 1.1 mrg }
72 1.1.1.4 mrg mpfr_clear (xx);
73 1.1.1.4 mrg mpfr_clear (yy);
74 1.1.1.4 mrg mpfr_clear (zz);
75 1.1 mrg }
76 1.1 mrg
77 1.1 mrg static void
78 1.1 mrg check53 (const char *xs, const char *ys, mpfr_rnd_t rnd_mode, const char *zs)
79 1.1 mrg {
80 1.1 mrg mpfr_t xx, yy, zz;
81 1.1 mrg
82 1.1 mrg mpfr_inits2 (53, xx, yy, zz, (mpfr_ptr) 0);
83 1.1 mrg mpfr_set_str1 (xx, xs);
84 1.1 mrg mpfr_set_str1 (yy, ys);
85 1.1 mrg test_mul (zz, xx, yy, rnd_mode);
86 1.1 mrg if (mpfr_cmp_str1 (zz, zs) )
87 1.1 mrg {
88 1.1 mrg printf ("(2) mpfr_mul failed for x=%s y=%s with rnd=%s\n",
89 1.1 mrg xs, ys, mpfr_print_rnd_mode(rnd_mode));
90 1.1 mrg printf ("correct result is %s,\n mpfr_mul gives ", zs);
91 1.1.1.4 mrg mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN);
92 1.1.1.4 mrg putchar ('\n');
93 1.1 mrg exit (1);
94 1.1 mrg }
95 1.1 mrg mpfr_clears (xx, yy, zz, (mpfr_ptr) 0);
96 1.1 mrg }
97 1.1 mrg
98 1.1 mrg /* checks that x*y gives the right result with 24 bits of precision */
99 1.1 mrg static void
100 1.1 mrg check24 (const char *xs, const char *ys, mpfr_rnd_t rnd_mode, const char *zs)
101 1.1 mrg {
102 1.1 mrg mpfr_t xx, yy, zz;
103 1.1 mrg
104 1.1 mrg mpfr_inits2 (24, xx, yy, zz, (mpfr_ptr) 0);
105 1.1 mrg mpfr_set_str1 (xx, xs);
106 1.1 mrg mpfr_set_str1 (yy, ys);
107 1.1 mrg test_mul (zz, xx, yy, rnd_mode);
108 1.1 mrg if (mpfr_cmp_str1 (zz, zs) )
109 1.1 mrg {
110 1.1 mrg printf ("(3) mpfr_mul failed for x=%s y=%s with "
111 1.1 mrg "rnd=%s\n", xs, ys, mpfr_print_rnd_mode(rnd_mode));
112 1.1 mrg printf ("correct result is gives %s, mpfr_mul gives ", zs);
113 1.1.1.4 mrg mpfr_out_str (stdout, 10, 0, zz, MPFR_RNDN);
114 1.1.1.4 mrg putchar ('\n');
115 1.1 mrg exit (1);
116 1.1 mrg }
117 1.1 mrg mpfr_clears (xx, yy, zz, (mpfr_ptr) 0);
118 1.1 mrg }
119 1.1 mrg
120 1.1 mrg /* the following examples come from the paper "Number-theoretic Test
121 1.1 mrg Generation for Directed Rounding" from Michael Parks, Table 1 */
122 1.1 mrg static void
123 1.1 mrg check_float (void)
124 1.1 mrg {
125 1.1 mrg check24("8388609.0", "8388609.0", MPFR_RNDN, "70368760954880.0");
126 1.1 mrg check24("16777213.0", "8388609.0", MPFR_RNDN, "140737479966720.0");
127 1.1 mrg check24("8388611.0", "8388609.0", MPFR_RNDN, "70368777732096.0");
128 1.1 mrg check24("12582911.0", "8388610.0", MPFR_RNDN, "105553133043712.0");
129 1.1 mrg check24("12582914.0", "8388610.0", MPFR_RNDN, "105553158209536.0");
130 1.1 mrg check24("13981013.0", "8388611.0", MPFR_RNDN, "117281279442944.0");
131 1.1 mrg check24("11184811.0", "8388611.0", MPFR_RNDN, "93825028587520.0");
132 1.1 mrg check24("11184810.0", "8388611.0", MPFR_RNDN, "93825020198912.0");
133 1.1 mrg check24("13981014.0", "8388611.0", MPFR_RNDN, "117281287831552.0");
134 1.1 mrg
135 1.1 mrg check24("8388609.0", "8388609.0", MPFR_RNDZ, "70368760954880.0");
136 1.1 mrg check24("16777213.0", "8388609.0", MPFR_RNDZ, "140737471578112.0");
137 1.1 mrg check24("8388611.0", "8388609.0", MPFR_RNDZ, "70368777732096.0");
138 1.1 mrg check24("12582911.0", "8388610.0", MPFR_RNDZ, "105553124655104.0");
139 1.1 mrg check24("12582914.0", "8388610.0", MPFR_RNDZ, "105553158209536.0");
140 1.1 mrg check24("13981013.0", "8388611.0", MPFR_RNDZ, "117281271054336.0");
141 1.1 mrg check24("11184811.0", "8388611.0", MPFR_RNDZ, "93825028587520.0");
142 1.1 mrg check24("11184810.0", "8388611.0", MPFR_RNDZ, "93825011810304.0");
143 1.1 mrg check24("13981014.0", "8388611.0", MPFR_RNDZ, "117281287831552.0");
144 1.1 mrg
145 1.1 mrg check24("8388609.0", "8388609.0", MPFR_RNDU, "70368769343488.0");
146 1.1 mrg check24("16777213.0", "8388609.0", MPFR_RNDU, "140737479966720.0");
147 1.1 mrg check24("8388611.0", "8388609.0", MPFR_RNDU, "70368786120704.0");
148 1.1 mrg check24("12582911.0", "8388610.0", MPFR_RNDU, "105553133043712.0");
149 1.1 mrg check24("12582914.0", "8388610.0", MPFR_RNDU, "105553166598144.0");
150 1.1 mrg check24("13981013.0", "8388611.0", MPFR_RNDU, "117281279442944.0");
151 1.1 mrg check24("11184811.0", "8388611.0", MPFR_RNDU, "93825036976128.0");
152 1.1 mrg check24("11184810.0", "8388611.0", MPFR_RNDU, "93825020198912.0");
153 1.1 mrg check24("13981014.0", "8388611.0", MPFR_RNDU, "117281296220160.0");
154 1.1 mrg
155 1.1 mrg check24("8388609.0", "8388609.0", MPFR_RNDD, "70368760954880.0");
156 1.1 mrg check24("16777213.0", "8388609.0", MPFR_RNDD, "140737471578112.0");
157 1.1 mrg check24("8388611.0", "8388609.0", MPFR_RNDD, "70368777732096.0");
158 1.1 mrg check24("12582911.0", "8388610.0", MPFR_RNDD, "105553124655104.0");
159 1.1 mrg check24("12582914.0", "8388610.0", MPFR_RNDD, "105553158209536.0");
160 1.1 mrg check24("13981013.0", "8388611.0", MPFR_RNDD, "117281271054336.0");
161 1.1 mrg check24("11184811.0", "8388611.0", MPFR_RNDD, "93825028587520.0");
162 1.1 mrg check24("11184810.0", "8388611.0", MPFR_RNDD, "93825011810304.0");
163 1.1 mrg check24("13981014.0", "8388611.0", MPFR_RNDD, "117281287831552.0");
164 1.1 mrg }
165 1.1 mrg
166 1.1 mrg /* check sign of result */
167 1.1 mrg static void
168 1.1 mrg check_sign (void)
169 1.1 mrg {
170 1.1 mrg mpfr_t a, b;
171 1.1 mrg
172 1.1 mrg mpfr_init2 (a, 53);
173 1.1 mrg mpfr_init2 (b, 53);
174 1.1 mrg mpfr_set_si (a, -1, MPFR_RNDN);
175 1.1 mrg mpfr_set_ui (b, 2, MPFR_RNDN);
176 1.1 mrg test_mul(a, b, b, MPFR_RNDN);
177 1.1 mrg if (mpfr_cmp_ui (a, 4) )
178 1.1 mrg {
179 1.1 mrg printf ("2.0*2.0 gives \n");
180 1.1.1.4 mrg mpfr_out_str (stdout, 10, 0, a, MPFR_RNDN);
181 1.1.1.4 mrg putchar ('\n');
182 1.1 mrg exit (1);
183 1.1 mrg }
184 1.1 mrg mpfr_clear(a); mpfr_clear(b);
185 1.1 mrg }
186 1.1 mrg
187 1.1 mrg /* checks that the inexact return value is correct */
188 1.1 mrg static void
189 1.1 mrg check_exact (void)
190 1.1 mrg {
191 1.1 mrg mpfr_t a, b, c, d;
192 1.1 mrg mpfr_prec_t prec;
193 1.1 mrg int i, inexact;
194 1.1 mrg mpfr_rnd_t rnd;
195 1.1 mrg
196 1.1 mrg mpfr_init (a);
197 1.1 mrg mpfr_init (b);
198 1.1 mrg mpfr_init (c);
199 1.1 mrg mpfr_init (d);
200 1.1 mrg
201 1.1 mrg mpfr_set_prec (a, 17);
202 1.1 mrg mpfr_set_prec (b, 17);
203 1.1 mrg mpfr_set_prec (c, 32);
204 1.1 mrg mpfr_set_str_binary (a, "1.1000111011000100e-1");
205 1.1 mrg mpfr_set_str_binary (b, "1.0010001111100111e-1");
206 1.1 mrg if (test_mul (c, a, b, MPFR_RNDZ))
207 1.1 mrg {
208 1.1 mrg printf ("wrong return value (1)\n");
209 1.1 mrg exit (1);
210 1.1 mrg }
211 1.1 mrg
212 1.1.1.4 mrg for (prec = MPFR_PREC_MIN; prec < 100; prec++)
213 1.1 mrg {
214 1.1 mrg mpfr_set_prec (a, prec);
215 1.1 mrg mpfr_set_prec (b, prec);
216 1.1.1.4 mrg /* for prec=1, ensure PREC(c) >= 1 */
217 1.1.1.4 mrg mpfr_set_prec (c, 2 * prec - 2 + (prec == 1));
218 1.1 mrg mpfr_set_prec (d, 2 * prec);
219 1.1 mrg for (i = 0; i < 1000; i++)
220 1.1 mrg {
221 1.1 mrg mpfr_urandomb (a, RANDS);
222 1.1 mrg mpfr_urandomb (b, RANDS);
223 1.1 mrg rnd = RND_RAND ();
224 1.1 mrg inexact = test_mul (c, a, b, rnd);
225 1.1 mrg if (test_mul (d, a, b, rnd)) /* should be always exact */
226 1.1 mrg {
227 1.1 mrg printf ("unexpected inexact return value\n");
228 1.1 mrg exit (1);
229 1.1 mrg }
230 1.1.1.4 mrg if ((inexact == 0) && mpfr_cmp (c, d) && rnd != MPFR_RNDF)
231 1.1 mrg {
232 1.1.1.4 mrg printf ("rnd=%s: inexact=0 but results differ\n",
233 1.1.1.4 mrg mpfr_print_rnd_mode (rnd));
234 1.1.1.4 mrg printf ("a=");
235 1.1.1.4 mrg mpfr_out_str (stdout, 2, 0, a, rnd);
236 1.1.1.4 mrg printf ("\nb=");
237 1.1.1.4 mrg mpfr_out_str (stdout, 2, 0, b, rnd);
238 1.1.1.4 mrg printf ("\nc=");
239 1.1.1.4 mrg mpfr_out_str (stdout, 2, 0, c, rnd);
240 1.1.1.4 mrg printf ("\nd=");
241 1.1.1.4 mrg mpfr_out_str (stdout, 2, 0, d, rnd);
242 1.1.1.4 mrg printf ("\n");
243 1.1 mrg exit (1);
244 1.1 mrg }
245 1.1.1.4 mrg else if (inexact && (mpfr_cmp (c, d) == 0) && rnd != MPFR_RNDF)
246 1.1 mrg {
247 1.1 mrg printf ("inexact!=0 but results agree\n");
248 1.1 mrg printf ("prec=%u rnd=%s a=", (unsigned int) prec,
249 1.1 mrg mpfr_print_rnd_mode (rnd));
250 1.1 mrg mpfr_out_str (stdout, 2, 0, a, rnd);
251 1.1 mrg printf ("\nb=");
252 1.1 mrg mpfr_out_str (stdout, 2, 0, b, rnd);
253 1.1 mrg printf ("\nc=");
254 1.1 mrg mpfr_out_str (stdout, 2, 0, c, rnd);
255 1.1 mrg printf ("\nd=");
256 1.1 mrg mpfr_out_str (stdout, 2, 0, d, rnd);
257 1.1 mrg printf ("\n");
258 1.1 mrg exit (1);
259 1.1 mrg }
260 1.1 mrg }
261 1.1 mrg }
262 1.1 mrg
263 1.1 mrg mpfr_clear (a);
264 1.1 mrg mpfr_clear (b);
265 1.1 mrg mpfr_clear (c);
266 1.1 mrg mpfr_clear (d);
267 1.1 mrg }
268 1.1 mrg
269 1.1 mrg static void
270 1.1 mrg check_max(void)
271 1.1 mrg {
272 1.1 mrg mpfr_t xx, yy, zz;
273 1.1 mrg mpfr_exp_t emin;
274 1.1 mrg
275 1.1 mrg mpfr_init2(xx, 4);
276 1.1 mrg mpfr_init2(yy, 4);
277 1.1 mrg mpfr_init2(zz, 4);
278 1.1 mrg mpfr_set_str1 (xx, "0.68750");
279 1.1 mrg mpfr_mul_2si(xx, xx, MPFR_EMAX_DEFAULT/2, MPFR_RNDN);
280 1.1 mrg mpfr_set_str1 (yy, "0.68750");
281 1.1 mrg mpfr_mul_2si(yy, yy, MPFR_EMAX_DEFAULT - MPFR_EMAX_DEFAULT/2 + 1, MPFR_RNDN);
282 1.1 mrg mpfr_clear_flags();
283 1.1 mrg test_mul(zz, xx, yy, MPFR_RNDU);
284 1.1 mrg if (!(mpfr_overflow_p() && MPFR_IS_INF(zz)))
285 1.1 mrg {
286 1.1 mrg printf("check_max failed (should be an overflow)\n");
287 1.1 mrg exit(1);
288 1.1 mrg }
289 1.1 mrg
290 1.1 mrg mpfr_clear_flags();
291 1.1 mrg test_mul(zz, xx, yy, MPFR_RNDD);
292 1.1 mrg if (mpfr_overflow_p() || MPFR_IS_INF(zz))
293 1.1 mrg {
294 1.1 mrg printf("check_max failed (should NOT be an overflow)\n");
295 1.1 mrg exit(1);
296 1.1 mrg }
297 1.1 mrg mpfr_set_str1 (xx, "0.93750");
298 1.1 mrg mpfr_mul_2si(xx, xx, MPFR_EMAX_DEFAULT, MPFR_RNDN);
299 1.1 mrg if (!(MPFR_IS_FP(xx) && MPFR_IS_FP(zz)))
300 1.1 mrg {
301 1.1 mrg printf("check_max failed (internal error)\n");
302 1.1 mrg exit(1);
303 1.1 mrg }
304 1.1 mrg if (mpfr_cmp(xx, zz) != 0)
305 1.1 mrg {
306 1.1 mrg printf("check_max failed: got ");
307 1.1 mrg mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ);
308 1.1 mrg printf(" instead of ");
309 1.1 mrg mpfr_out_str(stdout, 2, 0, xx, MPFR_RNDZ);
310 1.1 mrg printf("\n");
311 1.1 mrg exit(1);
312 1.1 mrg }
313 1.1 mrg
314 1.1 mrg /* check underflow */
315 1.1 mrg emin = mpfr_get_emin ();
316 1.1 mrg set_emin (0);
317 1.1 mrg mpfr_set_str_binary (xx, "0.1E0");
318 1.1 mrg mpfr_set_str_binary (yy, "0.1E0");
319 1.1 mrg test_mul (zz, xx, yy, MPFR_RNDN);
320 1.1 mrg /* exact result is 0.1E-1, which should round to 0 */
321 1.1 mrg MPFR_ASSERTN(mpfr_cmp_ui (zz, 0) == 0 && MPFR_IS_POS(zz));
322 1.1 mrg set_emin (emin);
323 1.1 mrg
324 1.1 mrg /* coverage test for mpfr_powerof2_raw */
325 1.1 mrg emin = mpfr_get_emin ();
326 1.1 mrg set_emin (0);
327 1.1 mrg mpfr_set_prec (xx, mp_bits_per_limb + 1);
328 1.1 mrg mpfr_set_str_binary (xx, "0.1E0");
329 1.1 mrg mpfr_nextabove (xx);
330 1.1 mrg mpfr_set_str_binary (yy, "0.1E0");
331 1.1 mrg test_mul (zz, xx, yy, MPFR_RNDN);
332 1.1 mrg /* exact result is just above 0.1E-1, which should round to minfloat */
333 1.1 mrg MPFR_ASSERTN(mpfr_cmp (zz, yy) == 0);
334 1.1 mrg set_emin (emin);
335 1.1 mrg
336 1.1 mrg mpfr_clear(xx);
337 1.1 mrg mpfr_clear(yy);
338 1.1 mrg mpfr_clear(zz);
339 1.1 mrg }
340 1.1 mrg
341 1.1 mrg static void
342 1.1 mrg check_min(void)
343 1.1 mrg {
344 1.1 mrg mpfr_t xx, yy, zz;
345 1.1 mrg
346 1.1 mrg mpfr_init2(xx, 4);
347 1.1 mrg mpfr_init2(yy, 4);
348 1.1 mrg mpfr_init2(zz, 3);
349 1.1 mrg mpfr_set_str1(xx, "0.9375");
350 1.1 mrg mpfr_mul_2si(xx, xx, MPFR_EMIN_DEFAULT/2, MPFR_RNDN);
351 1.1 mrg mpfr_set_str1(yy, "0.9375");
352 1.1 mrg mpfr_mul_2si(yy, yy, MPFR_EMIN_DEFAULT - MPFR_EMIN_DEFAULT/2 - 1, MPFR_RNDN);
353 1.1 mrg test_mul(zz, xx, yy, MPFR_RNDD);
354 1.1.1.4 mrg if (MPFR_NOTZERO (zz))
355 1.1 mrg {
356 1.1 mrg printf("check_min failed: got ");
357 1.1 mrg mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ);
358 1.1 mrg printf(" instead of 0\n");
359 1.1 mrg exit(1);
360 1.1 mrg }
361 1.1 mrg
362 1.1 mrg test_mul(zz, xx, yy, MPFR_RNDU);
363 1.1 mrg mpfr_set_str1 (xx, "0.5");
364 1.1 mrg mpfr_mul_2si(xx, xx, MPFR_EMIN_DEFAULT, MPFR_RNDN);
365 1.1 mrg if (mpfr_sgn(xx) <= 0)
366 1.1 mrg {
367 1.1 mrg printf("check_min failed (internal error)\n");
368 1.1 mrg exit(1);
369 1.1 mrg }
370 1.1 mrg if (mpfr_cmp(xx, zz) != 0)
371 1.1 mrg {
372 1.1 mrg printf("check_min failed: got ");
373 1.1 mrg mpfr_out_str(stdout, 2, 0, zz, MPFR_RNDZ);
374 1.1 mrg printf(" instead of ");
375 1.1 mrg mpfr_out_str(stdout, 2, 0, xx, MPFR_RNDZ);
376 1.1 mrg printf("\n");
377 1.1 mrg exit(1);
378 1.1 mrg }
379 1.1 mrg
380 1.1 mrg mpfr_clear(xx);
381 1.1 mrg mpfr_clear(yy);
382 1.1 mrg mpfr_clear(zz);
383 1.1 mrg }
384 1.1 mrg
385 1.1 mrg static void
386 1.1 mrg check_nans (void)
387 1.1 mrg {
388 1.1 mrg mpfr_t p, x, y;
389 1.1 mrg
390 1.1 mrg mpfr_init2 (x, 123L);
391 1.1 mrg mpfr_init2 (y, 123L);
392 1.1 mrg mpfr_init2 (p, 123L);
393 1.1 mrg
394 1.1 mrg /* nan * 0 == nan */
395 1.1 mrg mpfr_set_nan (x);
396 1.1 mrg mpfr_set_ui (y, 0L, MPFR_RNDN);
397 1.1 mrg test_mul (p, x, y, MPFR_RNDN);
398 1.1 mrg MPFR_ASSERTN (mpfr_nan_p (p));
399 1.1 mrg
400 1.1 mrg /* 1 * nan == nan */
401 1.1 mrg mpfr_set_ui (x, 1L, MPFR_RNDN);
402 1.1 mrg mpfr_set_nan (y);
403 1.1 mrg test_mul (p, x, y, MPFR_RNDN);
404 1.1 mrg MPFR_ASSERTN (mpfr_nan_p (p));
405 1.1 mrg
406 1.1 mrg /* 0 * +inf == nan */
407 1.1 mrg mpfr_set_ui (x, 0L, MPFR_RNDN);
408 1.1 mrg mpfr_set_nan (y);
409 1.1 mrg test_mul (p, x, y, MPFR_RNDN);
410 1.1 mrg MPFR_ASSERTN (mpfr_nan_p (p));
411 1.1 mrg
412 1.1 mrg /* +1 * +inf == +inf */
413 1.1 mrg mpfr_set_ui (x, 1L, MPFR_RNDN);
414 1.1 mrg mpfr_set_inf (y, 1);
415 1.1 mrg test_mul (p, x, y, MPFR_RNDN);
416 1.1 mrg MPFR_ASSERTN (mpfr_inf_p (p));
417 1.1 mrg MPFR_ASSERTN (mpfr_sgn (p) > 0);
418 1.1 mrg
419 1.1 mrg /* -1 * +inf == -inf */
420 1.1 mrg mpfr_set_si (x, -1L, MPFR_RNDN);
421 1.1 mrg mpfr_set_inf (y, 1);
422 1.1 mrg test_mul (p, x, y, MPFR_RNDN);
423 1.1 mrg MPFR_ASSERTN (mpfr_inf_p (p));
424 1.1 mrg MPFR_ASSERTN (mpfr_sgn (p) < 0);
425 1.1 mrg
426 1.1 mrg mpfr_clear (x);
427 1.1 mrg mpfr_clear (y);
428 1.1 mrg mpfr_clear (p);
429 1.1 mrg }
430 1.1 mrg
431 1.1 mrg #define BUFSIZE 1552
432 1.1 mrg
433 1.1 mrg static void
434 1.1 mrg get_string (char *s, FILE *fp)
435 1.1 mrg {
436 1.1 mrg int c, n = BUFSIZE;
437 1.1 mrg
438 1.1 mrg while ((c = getc (fp)) != '\n')
439 1.1 mrg {
440 1.1 mrg if (c == EOF)
441 1.1 mrg {
442 1.1 mrg printf ("Error in get_string: end of file\n");
443 1.1 mrg exit (1);
444 1.1 mrg }
445 1.1 mrg *(unsigned char *)s++ = c;
446 1.1 mrg if (--n == 0)
447 1.1 mrg {
448 1.1 mrg printf ("Error in get_string: buffer is too small\n");
449 1.1 mrg exit (1);
450 1.1 mrg }
451 1.1 mrg }
452 1.1 mrg *s = '\0';
453 1.1 mrg }
454 1.1 mrg
455 1.1 mrg static void
456 1.1 mrg check_regression (void)
457 1.1 mrg {
458 1.1 mrg mpfr_t x, y, z;
459 1.1 mrg int i;
460 1.1 mrg FILE *fp;
461 1.1 mrg char s[BUFSIZE];
462 1.1 mrg
463 1.1 mrg mpfr_inits2 (6177, x, y, z, (mpfr_ptr) 0);
464 1.1 mrg /* we read long strings from a file since ISO C90 does not support strings of
465 1.1 mrg length > 509 */
466 1.1 mrg fp = src_fopen ("tmul.dat", "r");
467 1.1 mrg if (fp == NULL)
468 1.1 mrg {
469 1.1 mrg fprintf (stderr, "Error, cannot open tmul.dat in srcdir\n");
470 1.1 mrg exit (1);
471 1.1 mrg }
472 1.1 mrg get_string (s, fp);
473 1.1 mrg mpfr_set_str (y, s, 16, MPFR_RNDN);
474 1.1 mrg get_string (s, fp);
475 1.1 mrg mpfr_set_str (z, s, 16, MPFR_RNDN);
476 1.1 mrg i = mpfr_mul (x, y, z, MPFR_RNDN);
477 1.1 mrg get_string (s, fp);
478 1.1 mrg if (mpfr_cmp_str (x, s, 16, MPFR_RNDN) != 0 || i != -1)
479 1.1 mrg {
480 1.1 mrg printf ("Regression test 1 failed (i=%d, expected -1)\nx=", i);
481 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n');
482 1.1 mrg exit (1);
483 1.1 mrg }
484 1.1 mrg fclose (fp);
485 1.1 mrg
486 1.1 mrg mpfr_set_prec (x, 606);
487 1.1 mrg mpfr_set_prec (y, 606);
488 1.1 mrg mpfr_set_prec (z, 606);
489 1.1 mrg
490 1.1 mrg mpfr_set_str (y, "-f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff92daefc3f8052ca9f58736564d9e93e62d324@-1", 16, MPFR_RNDN);
491 1.1 mrg mpfr_set_str (z, "-f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff92daefc3f8052ca9f58736564d9e93e62d324@-1", 16, MPFR_RNDN);
492 1.1 mrg i = mpfr_mul (x, y, z, MPFR_RNDU);
493 1.1 mrg mpfr_set_str (y, "f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff25b5df87f00a5953eb0e6cac9b3d27cc5a64c@-1", 16, MPFR_RNDN);
494 1.1 mrg if (mpfr_cmp (x, y) || i <= 0)
495 1.1 mrg {
496 1.1 mrg printf ("Regression test (2) failed! (i=%d - Expected 1)\n", i);
497 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN); putchar ('\n');
498 1.1 mrg exit (1);
499 1.1 mrg }
500 1.1 mrg
501 1.1 mrg mpfr_set_prec (x, 184);
502 1.1 mrg mpfr_set_prec (y, 92);
503 1.1 mrg mpfr_set_prec (z, 1023);
504 1.1 mrg
505 1.1 mrg mpfr_set_str (y, "6.9b8c8498882770d8038c3b0@-1", 16, MPFR_RNDN);
506 1.1 mrg mpfr_set_str (z, "7.44e24b986e7fb296f1e936ce749fec3504cbf0d5ba769466b1c9f1578115efd5d29b4c79271191a920a99280c714d3a657ad6e3afbab77ffce9d697e9bb9110e26d676069afcea8b69f1d1541f2365042d80a97c21dcccd8ace4f1bb58b49922003e738e6f37bb82ef653cb2e87f763974e6ae50ae54e7724c38b80653e3289@255", 16, MPFR_RNDN);
507 1.1 mrg i = mpfr_mul (x, y, z, MPFR_RNDU);
508 1.1 mrg mpfr_set_prec (y, 184);
509 1.1 mrg mpfr_set_str (y, "3.0080038f2ac5054e3e71ccbb95f76aaab2221715025a28@255",
510 1.1 mrg 16, MPFR_RNDN);
511 1.1 mrg if (mpfr_cmp (x, y) || i <= 0)
512 1.1 mrg {
513 1.1 mrg printf ("Regression test (4) failed! (i=%d - expected 1)\n", i);
514 1.1 mrg printf ("Ref: 3.0080038f2ac5054e3e71ccbb95f76aaab2221715025a28@255\n"
515 1.1 mrg "Got: ");
516 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN);
517 1.1 mrg printf ("\n");
518 1.1 mrg exit (1);
519 1.1 mrg }
520 1.1 mrg
521 1.1 mrg mpfr_set_prec (x, 908);
522 1.1 mrg mpfr_set_prec (y, 908);
523 1.1 mrg mpfr_set_prec (z, 908);
524 1.1 mrg mpfr_set_str (y, "-f.fffffffffffffffffffffffffffffffffffffffffffffffffffffff"
525 1.1 mrg "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff"
526 1.1 mrg "ffffffffffffffffffffffffffffffffffffffffffffffffffffff99be91f83ec6f0ed28a3d42"
527 1.1 mrg "e6e9a327230345ea6@-1", 16, MPFR_RNDN);
528 1.1 mrg mpfr_set_str (z, "-f.fffffffffffffffffffffffffffffffffffffffffffffffffffffff"
529 1.1 mrg "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff"
530 1.1 mrg "ffffffffffffffffffffffffffffffffffffffffffffffffffffff99be91f83ec6f0ed28a3d42"
531 1.1 mrg "e6e9a327230345ea6@-1", 16, MPFR_RNDN);
532 1.1 mrg i = mpfr_mul (x, y, z, MPFR_RNDU);
533 1.1 mrg mpfr_set_str (y, "f.ffffffffffffffffffffffffffffffffffffffffffffffffffffffff"
534 1.1 mrg "fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff"
535 1.1 mrg "fffffffffffffffffffffffffffffffffffffffffffffffffffff337d23f07d8de1da5147a85c"
536 1.1 mrg "dd3464e46068bd4d@-1", 16, MPFR_RNDN);
537 1.1 mrg if (mpfr_cmp (x, y) || i <= 0)
538 1.1 mrg {
539 1.1 mrg printf ("Regression test (5) failed! (i=%d - expected 1)\n", i);
540 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN);
541 1.1 mrg printf ("\n");
542 1.1 mrg exit (1);
543 1.1 mrg }
544 1.1 mrg
545 1.1 mrg
546 1.1 mrg mpfr_set_prec (x, 50);
547 1.1 mrg mpfr_set_prec (y, 40);
548 1.1 mrg mpfr_set_prec (z, 53);
549 1.1 mrg mpfr_set_str (y, "4.1ffffffff8", 16, MPFR_RNDN);
550 1.1 mrg mpfr_set_str (z, "4.2000000ffe0000@-4", 16, MPFR_RNDN);
551 1.1 mrg i = mpfr_mul (x, y, z, MPFR_RNDN);
552 1.1 mrg if (mpfr_cmp_str (x, "1.104000041d6c0@-3", 16, MPFR_RNDN) != 0
553 1.1 mrg || i <= 0)
554 1.1 mrg {
555 1.1 mrg printf ("Regression test (6) failed! (i=%d - expected 1)\nx=", i);
556 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN);
557 1.1 mrg printf ("\nMore prec=");
558 1.1 mrg mpfr_set_prec (x, 93);
559 1.1 mrg mpfr_mul (x, y, z, MPFR_RNDN);
560 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN);
561 1.1 mrg printf ("\n");
562 1.1 mrg exit (1);
563 1.1 mrg }
564 1.1 mrg
565 1.1 mrg mpfr_set_prec (x, 439);
566 1.1 mrg mpfr_set_prec (y, 393);
567 1.1 mrg mpfr_set_str (y, "-1.921fb54442d18469898cc51701b839a252049c1114cf98e804177d"
568 1.1 mrg "4c76273644a29410f31c6809bbdf2a33679a748636600",
569 1.1 mrg 16, MPFR_RNDN);
570 1.1 mrg i = mpfr_mul (x, y, y, MPFR_RNDU);
571 1.1 mrg if (mpfr_cmp_str (x, "2.77a79937c8bbcb495b89b36602306b1c2159a8ff834288a19a08"
572 1.1 mrg "84094f1cda3dc426da61174c4544a173de83c2500f8bfea2e0569e3698",
573 1.1 mrg 16, MPFR_RNDN) != 0
574 1.1 mrg || i <= 0)
575 1.1 mrg {
576 1.1 mrg printf ("Regression test (7) failed! (i=%d - expected 1)\nx=", i);
577 1.1 mrg mpfr_out_str (stdout, 16, 0, x, MPFR_RNDN);
578 1.1 mrg printf ("\n");
579 1.1 mrg exit (1);
580 1.1 mrg }
581 1.1 mrg
582 1.1 mrg mpfr_set_prec (x, 1023);
583 1.1 mrg mpfr_set_prec (y, 1023);
584 1.1 mrg mpfr_set_prec (z, 511);
585 1.1 mrg mpfr_set_ui (x, 17, MPFR_RNDN);
586 1.1 mrg mpfr_set_ui (y, 42, MPFR_RNDN);
587 1.1 mrg i = mpfr_mul (z, x, y, MPFR_RNDN);
588 1.1 mrg if (mpfr_cmp_ui (z, 17*42) != 0 || i != 0)
589 1.1 mrg {
590 1.1 mrg printf ("Regression test (8) failed! (i=%d - expected 0)\nz=", i);
591 1.1 mrg mpfr_out_str (stdout, 16, 0, z, MPFR_RNDN);
592 1.1 mrg printf ("\n");
593 1.1 mrg exit (1);
594 1.1 mrg }
595 1.1 mrg
596 1.1 mrg mpfr_clears (x, y, z, (mpfr_ptr) 0);
597 1.1 mrg }
598 1.1 mrg
599 1.1 mrg #define TEST_FUNCTION test_mul
600 1.1 mrg #define TWO_ARGS
601 1.1 mrg #define RAND_FUNCTION(x) mpfr_random2(x, MPFR_LIMB_SIZE (x), randlimb () % 100, RANDS)
602 1.1 mrg #include "tgeneric.c"
603 1.1 mrg
604 1.1 mrg /* multiplies x by 53-bit approximation of Pi */
605 1.1 mrg static int
606 1.1 mrg mpfr_mulpi (mpfr_ptr y, mpfr_srcptr x, mpfr_rnd_t r)
607 1.1 mrg {
608 1.1 mrg mpfr_t z;
609 1.1 mrg int inex;
610 1.1 mrg
611 1.1 mrg mpfr_init2 (z, 53);
612 1.1 mrg mpfr_set_str_binary (z, "11.001001000011111101101010100010001000010110100011");
613 1.1 mrg inex = mpfr_mul (y, x, z, r);
614 1.1 mrg mpfr_clear (z);
615 1.1 mrg return inex;
616 1.1 mrg }
617 1.1 mrg
618 1.1.1.2 mrg static void
619 1.1.1.2 mrg valgrind20110503 (void)
620 1.1.1.2 mrg {
621 1.1.1.2 mrg mpfr_t a, b, c;
622 1.1.1.2 mrg
623 1.1.1.2 mrg mpfr_init2 (a, 2);
624 1.1.1.2 mrg mpfr_init2 (b, 2005);
625 1.1.1.2 mrg mpfr_init2 (c, 2);
626 1.1.1.2 mrg
627 1.1.1.2 mrg mpfr_set_ui (b, 5, MPFR_RNDN);
628 1.1.1.2 mrg mpfr_nextabove (b);
629 1.1.1.2 mrg mpfr_set_ui (c, 1, MPFR_RNDN);
630 1.1.1.2 mrg mpfr_mul (a, b, c, MPFR_RNDZ);
631 1.1.1.2 mrg /* After the call to mpfr_mulhigh_n, valgrind complains:
632 1.1.1.2 mrg Conditional jump or move depends on uninitialised value(s) */
633 1.1.1.2 mrg
634 1.1.1.2 mrg mpfr_clears (a, b, c, (mpfr_ptr) 0);
635 1.1.1.2 mrg }
636 1.1.1.2 mrg
637 1.1.1.4 mrg static void
638 1.1.1.4 mrg testall_rndf (mpfr_prec_t pmax)
639 1.1.1.4 mrg {
640 1.1.1.4 mrg mpfr_t a, b, c, d;
641 1.1.1.4 mrg mpfr_prec_t pa, pb, pc;
642 1.1.1.4 mrg
643 1.1.1.4 mrg for (pa = MPFR_PREC_MIN; pa <= pmax; pa++)
644 1.1.1.4 mrg {
645 1.1.1.4 mrg mpfr_init2 (a, pa);
646 1.1.1.4 mrg mpfr_init2 (d, pa);
647 1.1.1.4 mrg for (pb = MPFR_PREC_MIN; pb <= pmax; pb++)
648 1.1.1.4 mrg {
649 1.1.1.4 mrg mpfr_init2 (b, pb);
650 1.1.1.4 mrg mpfr_set_ui (b, 1, MPFR_RNDN);
651 1.1.1.4 mrg while (mpfr_cmp_ui (b, 2) < 0)
652 1.1.1.4 mrg {
653 1.1.1.4 mrg for (pc = MPFR_PREC_MIN; pc <= pmax; pc++)
654 1.1.1.4 mrg {
655 1.1.1.4 mrg mpfr_init2 (c, pc);
656 1.1.1.4 mrg mpfr_set_ui (c, 1, MPFR_RNDN);
657 1.1.1.4 mrg while (mpfr_cmp_ui (c, 2) < 0)
658 1.1.1.4 mrg {
659 1.1.1.4 mrg mpfr_mul (a, b, c, MPFR_RNDF);
660 1.1.1.4 mrg mpfr_mul (d, b, c, MPFR_RNDD);
661 1.1.1.4 mrg if (!mpfr_equal_p (a, d))
662 1.1.1.4 mrg {
663 1.1.1.4 mrg mpfr_mul (d, b, c, MPFR_RNDU);
664 1.1.1.4 mrg if (!mpfr_equal_p (a, d))
665 1.1.1.4 mrg {
666 1.1.1.4 mrg printf ("Error: mpfr_mul(a,b,c,RNDF) does not "
667 1.1.1.4 mrg "match RNDD/RNDU\n");
668 1.1.1.4 mrg printf ("b="); mpfr_dump (b);
669 1.1.1.4 mrg printf ("c="); mpfr_dump (c);
670 1.1.1.4 mrg printf ("a="); mpfr_dump (a);
671 1.1.1.4 mrg exit (1);
672 1.1.1.4 mrg }
673 1.1.1.4 mrg }
674 1.1.1.4 mrg mpfr_nextabove (c);
675 1.1.1.4 mrg }
676 1.1.1.4 mrg mpfr_clear (c);
677 1.1.1.4 mrg }
678 1.1.1.4 mrg mpfr_nextabove (b);
679 1.1.1.4 mrg }
680 1.1.1.4 mrg mpfr_clear (b);
681 1.1.1.4 mrg }
682 1.1.1.4 mrg mpfr_clear (a);
683 1.1.1.4 mrg mpfr_clear (d);
684 1.1.1.4 mrg }
685 1.1.1.4 mrg }
686 1.1.1.4 mrg
687 1.1.1.3 mrg /* Check underflow flag corresponds to *after* rounding.
688 1.1.1.3 mrg *
689 1.1.1.3 mrg * More precisely, we want to test mpfr_mul on inputs b and c such that
690 1.1.1.3 mrg * EXP(b*c) < emin but EXP(round(b*c, p, rnd)) = emin. Thus an underflow
691 1.1.1.3 mrg * must not be generated.
692 1.1.1.3 mrg */
693 1.1.1.3 mrg static void
694 1.1.1.3 mrg test_underflow (mpfr_prec_t pmax)
695 1.1.1.3 mrg {
696 1.1.1.3 mrg mpfr_exp_t emin;
697 1.1.1.3 mrg mpfr_prec_t p;
698 1.1.1.3 mrg mpfr_t a0, a, b, c;
699 1.1.1.3 mrg int inex;
700 1.1.1.3 mrg
701 1.1.1.3 mrg mpfr_init2 (a0, MPFR_PREC_MIN);
702 1.1.1.3 mrg emin = mpfr_get_emin ();
703 1.1.1.3 mrg mpfr_setmin (a0, emin); /* 0.5 * 2^emin */
704 1.1.1.3 mrg
705 1.1.1.3 mrg /* for RNDN, we want b*c < 0.5 * 2^emin but RNDN(b*c, p) = 0.5 * 2^emin,
706 1.1.1.3 mrg thus b*c >= (0.5 - 1/4 * ulp_p(0.5)) * 2^emin */
707 1.1.1.3 mrg for (p = MPFR_PREC_MIN; p <= pmax; p++)
708 1.1.1.3 mrg {
709 1.1.1.3 mrg mpfr_init2 (a, p + 1);
710 1.1.1.3 mrg mpfr_init2 (b, p + 10);
711 1.1.1.3 mrg mpfr_init2 (c, p + 10);
712 1.1.1.3 mrg do mpfr_urandomb (b, RANDS); while (MPFR_IS_ZERO (b));
713 1.1.1.3 mrg inex = mpfr_set_ui_2exp (a, 1, -1, MPFR_RNDZ); /* a = 0.5 */
714 1.1.1.3 mrg MPFR_ASSERTN (inex == 0);
715 1.1.1.3 mrg mpfr_nextbelow (a); /* 0.5 - 1/2*ulp_{p+1}(0.5) = 0.5 - 1/4*ulp_p(0.5) */
716 1.1.1.3 mrg inex = mpfr_div (c, a, b, MPFR_RNDU);
717 1.1.1.3 mrg /* 0.5 - 1/4 * ulp_p(0.5) = a <= b*c < 0.5 */
718 1.1.1.3 mrg mpfr_mul_2si (b, b, emin / 2, MPFR_RNDZ);
719 1.1.1.3 mrg mpfr_mul_2si (c, c, (emin - 1) / 2, MPFR_RNDZ);
720 1.1.1.3 mrg /* now (0.5 - 1/4 * ulp_p(0.5)) * 2^emin <= b*c < 0.5 * 2^emin,
721 1.1.1.3 mrg thus b*c should be rounded to 0.5 * 2^emin */
722 1.1.1.3 mrg mpfr_set_prec (a, p);
723 1.1.1.3 mrg mpfr_clear_underflow ();
724 1.1.1.3 mrg mpfr_mul (a, b, c, MPFR_RNDN);
725 1.1.1.3 mrg if (mpfr_underflow_p () || ! mpfr_equal_p (a, a0))
726 1.1.1.3 mrg {
727 1.1.1.3 mrg printf ("Error, b*c incorrect or underflow flag incorrectly set"
728 1.1.1.3 mrg " for emin=%" MPFR_EXP_FSPEC "d, rnd=%s\n",
729 1.1.1.3 mrg (mpfr_eexp_t) emin, mpfr_print_rnd_mode (MPFR_RNDN));
730 1.1.1.3 mrg printf ("b="); mpfr_dump (b);
731 1.1.1.3 mrg printf ("c="); mpfr_dump (c);
732 1.1.1.3 mrg printf ("a="); mpfr_dump (a);
733 1.1.1.3 mrg mpfr_set_prec (a, mpfr_get_prec (b) + mpfr_get_prec (c));
734 1.1.1.3 mrg mpfr_mul_2exp (b, b, 1, MPFR_RNDN);
735 1.1.1.3 mrg inex = mpfr_mul (a, b, c, MPFR_RNDN);
736 1.1.1.3 mrg MPFR_ASSERTN (inex == 0);
737 1.1.1.3 mrg printf ("Exact 2*a="); mpfr_dump (a);
738 1.1.1.3 mrg exit (1);
739 1.1.1.3 mrg }
740 1.1.1.3 mrg mpfr_clear (a);
741 1.1.1.3 mrg mpfr_clear (b);
742 1.1.1.3 mrg mpfr_clear (c);
743 1.1.1.3 mrg }
744 1.1.1.3 mrg
745 1.1.1.3 mrg /* for RNDU, we want b*c < 0.5*2^emin but RNDU(b*c, p) = 0.5*2^emin thus
746 1.1.1.3 mrg b*c > (0.5 - 1/2 * ulp_p(0.5)) * 2^emin */
747 1.1.1.3 mrg for (p = MPFR_PREC_MIN; p <= pmax; p++)
748 1.1.1.3 mrg {
749 1.1.1.3 mrg mpfr_init2 (a, p);
750 1.1.1.3 mrg mpfr_init2 (b, p + 10);
751 1.1.1.3 mrg mpfr_init2 (c, p + 10);
752 1.1.1.3 mrg do mpfr_urandomb (b, RANDS); while (MPFR_IS_ZERO (b));
753 1.1.1.3 mrg inex = mpfr_set_ui_2exp (a, 1, -1, MPFR_RNDZ); /* a = 0.5 */
754 1.1.1.3 mrg MPFR_ASSERTN (inex == 0);
755 1.1.1.3 mrg mpfr_nextbelow (a); /* 0.5 - 1/2 * ulp_p(0.5) */
756 1.1.1.3 mrg inex = mpfr_div (c, a, b, MPFR_RNDU);
757 1.1.1.3 mrg /* 0.5 - 1/2 * ulp_p(0.5) <= b*c < 0.5 */
758 1.1.1.3 mrg mpfr_mul_2si (b, b, emin / 2, MPFR_RNDZ);
759 1.1.1.3 mrg mpfr_mul_2si (c, c, (emin - 1) / 2, MPFR_RNDZ);
760 1.1.1.3 mrg if (inex == 0)
761 1.1.1.3 mrg mpfr_nextabove (c); /* ensures b*c > (0.5 - 1/2 * ulp_p(0.5)) * 2^emin.
762 1.1.1.3 mrg Warning: for p=1, 0.5 - 1/2 * ulp_p(0.5)
763 1.1.1.3 mrg = 0.25, thus b*c > 2^(emin-2), which should
764 1.1.1.3 mrg also be rounded up with p=1 to 0.5 * 2^emin
765 1.1.1.3 mrg with an unbounded exponent range. */
766 1.1.1.3 mrg mpfr_clear_underflow ();
767 1.1.1.3 mrg mpfr_mul (a, b, c, MPFR_RNDU);
768 1.1.1.3 mrg if (mpfr_underflow_p () || ! mpfr_equal_p (a, a0))
769 1.1.1.3 mrg {
770 1.1.1.3 mrg printf ("Error, b*c incorrect or underflow flag incorrectly set"
771 1.1.1.3 mrg " for emin=%" MPFR_EXP_FSPEC "d, rnd=%s\n",
772 1.1.1.3 mrg (mpfr_eexp_t) emin, mpfr_print_rnd_mode (MPFR_RNDU));
773 1.1.1.3 mrg printf ("b="); mpfr_dump (b);
774 1.1.1.3 mrg printf ("c="); mpfr_dump (c);
775 1.1.1.3 mrg printf ("a="); mpfr_dump (a);
776 1.1.1.3 mrg mpfr_set_prec (a, mpfr_get_prec (b) + mpfr_get_prec (c));
777 1.1.1.3 mrg mpfr_mul_2exp (b, b, 1, MPFR_RNDN);
778 1.1.1.3 mrg inex = mpfr_mul (a, b, c, MPFR_RNDN);
779 1.1.1.3 mrg MPFR_ASSERTN (inex == 0);
780 1.1.1.3 mrg printf ("Exact 2*a="); mpfr_dump (a);
781 1.1.1.3 mrg exit (1);
782 1.1.1.3 mrg }
783 1.1.1.3 mrg mpfr_clear (a);
784 1.1.1.3 mrg mpfr_clear (b);
785 1.1.1.3 mrg mpfr_clear (c);
786 1.1.1.3 mrg }
787 1.1.1.3 mrg
788 1.1.1.3 mrg mpfr_clear (a0);
789 1.1.1.3 mrg }
790 1.1.1.3 mrg
791 1.1.1.4 mrg /* checks special case where no underflow should occur */
792 1.1.1.4 mrg static void
793 1.1.1.4 mrg bug20161209 (void)
794 1.1.1.4 mrg {
795 1.1.1.4 mrg mpfr_exp_t emin;
796 1.1.1.4 mrg mpfr_t x, y, z;
797 1.1.1.4 mrg
798 1.1.1.4 mrg emin = mpfr_get_emin ();
799 1.1.1.4 mrg set_emin (-1);
800 1.1.1.4 mrg
801 1.1.1.4 mrg /* test for mpfr_mul_1 for 64-bit limb, mpfr_mul_2 for 32-bit limb */
802 1.1.1.4 mrg mpfr_init2 (x, 53);
803 1.1.1.4 mrg mpfr_init2 (y, 53);
804 1.1.1.4 mrg mpfr_init2 (z, 53);
805 1.1.1.4 mrg mpfr_set_str_binary (x, "0.101000001E-1"); /* x = 321/2^10 */
806 1.1.1.4 mrg mpfr_set_str_binary (y, "0.110011000010100101111000011011000111011000001E-1");
807 1.1.1.4 mrg /* y = 28059810762433/2^46 */
808 1.1.1.4 mrg /* x * y = (2^53+1)/2^56 = 0.001...000[1]000..., and should round to 0.25 */
809 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDN);
810 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0);
811 1.1.1.4 mrg
812 1.1.1.4 mrg /* test for mpfr_mul_2 for 64-bit limb */
813 1.1.1.4 mrg mpfr_set_prec (x, 65);
814 1.1.1.4 mrg mpfr_set_prec (y, 65);
815 1.1.1.4 mrg mpfr_set_prec (z, 65);
816 1.1.1.4 mrg mpfr_set_str_binary (x, "0.101101000010010110100001E-1"); /* 11806113/2^25 */
817 1.1.1.4 mrg mpfr_set_str_binary (y, "0.101101011110010101011010001111111001100001E-1");
818 1.1.1.4 mrg /* y = 3124947910241/2^43 */
819 1.1.1.4 mrg /* x * y = (2^65+1)/2^68 = 0.001...000[1]000..., and should round to 0.25 */
820 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDN);
821 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0);
822 1.1.1.4 mrg
823 1.1.1.4 mrg /* test for the generic code */
824 1.1.1.4 mrg mpfr_set_prec (x, 54);
825 1.1.1.4 mrg mpfr_set_prec (y, 55);
826 1.1.1.4 mrg mpfr_set_prec (z, 53);
827 1.1.1.4 mrg mpfr_set_str_binary (x, "0.101000001E-1");
828 1.1.1.4 mrg mpfr_set_str_binary (y, "0.110011000010100101111000011011000111011000001E-1");
829 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDN);
830 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0);
831 1.1.1.4 mrg
832 1.1.1.4 mrg /* another test for the generic code */
833 1.1.1.4 mrg mpfr_set_prec (x, 66);
834 1.1.1.4 mrg mpfr_set_prec (y, 67);
835 1.1.1.4 mrg mpfr_set_prec (z, 65);
836 1.1.1.4 mrg mpfr_set_str_binary (x, "0.101101000010010110100001E-1");
837 1.1.1.4 mrg mpfr_set_str_binary (y, "0.101101011110010101011010001111111001100001E-1");
838 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDN);
839 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -2) == 0);
840 1.1.1.4 mrg
841 1.1.1.4 mrg mpfr_clear (x);
842 1.1.1.4 mrg mpfr_clear (y);
843 1.1.1.4 mrg mpfr_clear (z);
844 1.1.1.4 mrg set_emin (emin);
845 1.1.1.4 mrg }
846 1.1.1.4 mrg
847 1.1.1.4 mrg /* test for case in mpfr_mul_1() where:
848 1.1.1.4 mrg ax = __gmpfr_emin - 1
849 1.1.1.4 mrg ap[0] == ~mask
850 1.1.1.4 mrg rnd_mode = MPFR_RNDZ.
851 1.1.1.4 mrg Whatever the values of rb and sb, we should round to zero (underflow). */
852 1.1.1.4 mrg static void
853 1.1.1.4 mrg bug20161209a (void)
854 1.1.1.4 mrg {
855 1.1.1.4 mrg mpfr_exp_t emin;
856 1.1.1.4 mrg mpfr_t x, y, z;
857 1.1.1.4 mrg
858 1.1.1.4 mrg emin = mpfr_get_emin ();
859 1.1.1.4 mrg set_emin (-1);
860 1.1.1.4 mrg
861 1.1.1.4 mrg mpfr_init2 (x, 53);
862 1.1.1.4 mrg mpfr_init2 (y, 53);
863 1.1.1.4 mrg mpfr_init2 (z, 53);
864 1.1.1.4 mrg
865 1.1.1.4 mrg /* case rb = sb = 0 */
866 1.1.1.4 mrg mpfr_set_str_binary (x, "0.11010010100110000110110011111E-1");
867 1.1.1.4 mrg mpfr_set_str_binary (y, "0.1001101110011000110100001");
868 1.1.1.4 mrg /* x = 441650591/2^30, y = 20394401/2^25 */
869 1.1.1.4 mrg /* x * y = (2^53-1)/2^55 = 0.00111...111[0]000..., and should round to 0 */
870 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDZ);
871 1.1.1.4 mrg MPFR_ASSERTN(mpfr_zero_p (z));
872 1.1.1.4 mrg
873 1.1.1.4 mrg /* case rb = 1, sb = 0 */
874 1.1.1.4 mrg mpfr_set_str_binary (x, "0.111111111000000000000000000111111111E-1");
875 1.1.1.4 mrg mpfr_set_str_binary (y, "0.1000000001000000001");
876 1.1.1.4 mrg /* x = 68585259519/2^37, y = 262657/2^19 */
877 1.1.1.4 mrg /* x * y = (2^54-1)/2^56 = 0.00111...111[1]000..., and should round to 0 */
878 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDZ);
879 1.1.1.4 mrg MPFR_ASSERTN(mpfr_zero_p (z));
880 1.1.1.4 mrg
881 1.1.1.4 mrg /* case rb = 0, sb = 1 */
882 1.1.1.4 mrg mpfr_set_str_binary (x, "0.110010011001011110001100100001000001E-1");
883 1.1.1.4 mrg mpfr_set_str_binary (y, "0.10100010100010111101");
884 1.1.1.4 mrg /* x = 541144371852^37, y = 665789/2^20 */
885 1.1.1.4 mrg /* x * y = (2^55-3)/2^57 = 0.00111...111[0]100..., and should round to 0 */
886 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDZ);
887 1.1.1.4 mrg MPFR_ASSERTN(mpfr_zero_p (z));
888 1.1.1.4 mrg
889 1.1.1.4 mrg /* case rb = sb = 1 */
890 1.1.1.4 mrg mpfr_set_str_binary (x, "0.10100110001001001010001111110010100111E-1");
891 1.1.1.4 mrg mpfr_set_str_binary (y, "0.110001010011101001");
892 1.1.1.4 mrg /* x = 178394823847/2^39, y = 201961/2^18 */
893 1.1.1.4 mrg /* x * y = (2^55-1)/2^57 = 0.00111...111[1]100..., and should round to 0 */
894 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDZ);
895 1.1.1.4 mrg MPFR_ASSERTN(mpfr_zero_p (z));
896 1.1.1.4 mrg
897 1.1.1.4 mrg /* similar test for mpfr_mul_2 (we only check rb = sb = 1 here) */
898 1.1.1.4 mrg mpfr_set_prec (x, 65);
899 1.1.1.4 mrg mpfr_set_prec (y, 65);
900 1.1.1.4 mrg mpfr_set_prec (z, 65);
901 1.1.1.4 mrg /* 2^67-1 = 193707721 * 761838257287 */
902 1.1.1.4 mrg mpfr_set_str_binary (x, "0.1011100010111011111011001001E-1");
903 1.1.1.4 mrg mpfr_set_str_binary (y, "0.1011000101100001000110010100010010000111");
904 1.1.1.4 mrg mpfr_mul (z, x, y, MPFR_RNDZ);
905 1.1.1.4 mrg MPFR_ASSERTN(mpfr_zero_p (z));
906 1.1.1.4 mrg
907 1.1.1.4 mrg mpfr_clear (x);
908 1.1.1.4 mrg mpfr_clear (y);
909 1.1.1.4 mrg mpfr_clear (z);
910 1.1.1.4 mrg set_emin (emin);
911 1.1.1.4 mrg }
912 1.1.1.4 mrg
913 1.1.1.4 mrg /* bug for RNDF */
914 1.1.1.4 mrg static void
915 1.1.1.4 mrg bug20170602 (void)
916 1.1.1.4 mrg {
917 1.1.1.4 mrg mpfr_t x, u, y, yd, yu;
918 1.1.1.4 mrg
919 1.1.1.4 mrg mpfr_init2 (x, 493);
920 1.1.1.4 mrg mpfr_init2 (u, 493);
921 1.1.1.4 mrg mpfr_init2 (y, 503);
922 1.1.1.4 mrg mpfr_init2 (yd, 503);
923 1.1.1.4 mrg mpfr_init2 (yu, 503);
924 1.1.1.4 mrg mpfr_set_str_binary (x, "0.1111100000000000001111111110000000001111111111111000000000000000000011111111111111111111111000000000000000000001111111111111111111111111111111111111111000000000011111111111111111111000000000011111111111111000000000000001110000000000000000000000000000000000000000011111111111110011111111111100000000000000011111111111111111110000000011111111111111111110011111111111110000000000001111111111111111000000000000000000000000000000000000111111111111111111111111111111011111111111111111111111111111100E44");
925 1.1.1.4 mrg mpfr_set_str_binary (u, "0.1110000000000000001111111111111111111111111111111111111000000000111111111111111111111111111111000000000000000000001111111000000000000000011111111111111111111111111111111111111111111111111111111000000000000000011111111111111000000011111111111111111110000000000000001111111111111111111111111111111111111110000000000001111111111111111111111111111111111111000000000000000000000000000000000001111111111111000000000000000000001111111111100000000000000011111111111111111111111111111111111111111111111E35");
926 1.1.1.4 mrg mpfr_mul (y, x, u, MPFR_RNDF);
927 1.1.1.4 mrg mpfr_mul (yd, x, u, MPFR_RNDD);
928 1.1.1.4 mrg mpfr_mul (yu, x, u, MPFR_RNDU);
929 1.1.1.4 mrg if (mpfr_cmp (y, yd) != 0 && mpfr_cmp (y, yu) != 0)
930 1.1.1.4 mrg {
931 1.1.1.4 mrg printf ("RNDF is neither RNDD nor RNDU\n");
932 1.1.1.4 mrg printf ("x"); mpfr_dump (x);
933 1.1.1.4 mrg printf ("u"); mpfr_dump (u);
934 1.1.1.4 mrg printf ("y(RNDF)="); mpfr_dump (y);
935 1.1.1.4 mrg printf ("y(RNDD)="); mpfr_dump (yd);
936 1.1.1.4 mrg printf ("y(RNDU)="); mpfr_dump (yu);
937 1.1.1.4 mrg exit (1);
938 1.1.1.4 mrg }
939 1.1.1.4 mrg mpfr_clear (x);
940 1.1.1.4 mrg mpfr_clear (u);
941 1.1.1.4 mrg mpfr_clear (y);
942 1.1.1.4 mrg mpfr_clear (yd);
943 1.1.1.4 mrg mpfr_clear (yu);
944 1.1.1.4 mrg }
945 1.1.1.4 mrg
946 1.1.1.4 mrg /* Test for 1 to 3 limbs. */
947 1.1.1.4 mrg static void
948 1.1.1.4 mrg small_prec (void)
949 1.1.1.4 mrg {
950 1.1.1.4 mrg mpfr_exp_t emin, emax;
951 1.1.1.4 mrg mpfr_t x, y, z1, z2, zz;
952 1.1.1.4 mrg int xq, yq, zq;
953 1.1.1.4 mrg mpfr_rnd_t rnd;
954 1.1.1.4 mrg mpfr_flags_t flags1, flags2;
955 1.1.1.4 mrg int inex1, inex2;
956 1.1.1.4 mrg int i, j, r, s, ediff;
957 1.1.1.4 mrg
958 1.1.1.4 mrg emin = mpfr_get_emin ();
959 1.1.1.4 mrg emax = mpfr_get_emax ();
960 1.1.1.4 mrg
961 1.1.1.4 mrg /* The mpfr_mul implementation doesn't extend the exponent range,
962 1.1.1.4 mrg so that it is OK to reduce it here for the test to make sure that
963 1.1.1.4 mrg mpfr_mul_2si can be used. */
964 1.1.1.4 mrg set_emin (-1000);
965 1.1.1.4 mrg set_emax (1000);
966 1.1.1.4 mrg
967 1.1.1.4 mrg mpfr_inits2 (3 * GMP_NUMB_BITS, x, y, z1, z2, (mpfr_ptr) 0);
968 1.1.1.4 mrg mpfr_init2 (zz, 6 * GMP_NUMB_BITS);
969 1.1.1.4 mrg for (i = 0; i < 3; i++)
970 1.1.1.4 mrg for (j = 0; j < 10000; j++)
971 1.1.1.4 mrg {
972 1.1.1.4 mrg xq = i * GMP_NUMB_BITS + 1 + randlimb () % GMP_NUMB_BITS;
973 1.1.1.4 mrg mpfr_set_prec (x, xq);
974 1.1.1.4 mrg yq = i * GMP_NUMB_BITS + 1 + randlimb () % GMP_NUMB_BITS;
975 1.1.1.4 mrg mpfr_set_prec (y, yq);
976 1.1.1.4 mrg zq = i * GMP_NUMB_BITS + 1 + randlimb () % (GMP_NUMB_BITS-1);
977 1.1.1.4 mrg mpfr_set_prec (z1, zq);
978 1.1.1.4 mrg mpfr_set_prec (z2, zq);
979 1.1.1.4 mrg s = r = randlimb () & 0x7f;
980 1.1.1.4 mrg do mpfr_urandomb (x, RANDS); while (mpfr_zero_p (x));
981 1.1.1.4 mrg if (s & 1)
982 1.1.1.4 mrg mpfr_neg (x, x, MPFR_RNDN);
983 1.1.1.4 mrg s >>= 1;
984 1.1.1.4 mrg if (s & 1)
985 1.1.1.4 mrg {
986 1.1.1.4 mrg do mpfr_urandomb (y, RANDS); while (mpfr_zero_p (y));
987 1.1.1.4 mrg }
988 1.1.1.4 mrg else
989 1.1.1.4 mrg {
990 1.1.1.4 mrg mpfr_ui_div (y, 1, x, MPFR_RNDN);
991 1.1.1.4 mrg mpfr_set_exp (y, 0);
992 1.1.1.4 mrg }
993 1.1.1.4 mrg s >>= 1;
994 1.1.1.4 mrg if (s & 1)
995 1.1.1.4 mrg mpfr_neg (y, y, MPFR_RNDN);
996 1.1.1.4 mrg s >>= 1;
997 1.1.1.4 mrg rnd = RND_RAND_NO_RNDF ();
998 1.1.1.4 mrg inex1 = mpfr_mul (zz, x, y, MPFR_RNDN);
999 1.1.1.4 mrg MPFR_ASSERTN (inex1 == 0);
1000 1.1.1.4 mrg if (s == 0)
1001 1.1.1.4 mrg {
1002 1.1.1.4 mrg ediff = __gmpfr_emin - MPFR_EXP (x);
1003 1.1.1.4 mrg mpfr_set_exp (x, __gmpfr_emin);
1004 1.1.1.4 mrg }
1005 1.1.1.4 mrg else if (s == 1)
1006 1.1.1.4 mrg {
1007 1.1.1.4 mrg ediff = __gmpfr_emax - MPFR_EXP (x) + 1;
1008 1.1.1.4 mrg mpfr_set_exp (x, __gmpfr_emax);
1009 1.1.1.4 mrg mpfr_mul_2ui (y, y, 1, MPFR_RNDN);
1010 1.1.1.4 mrg }
1011 1.1.1.4 mrg else
1012 1.1.1.4 mrg ediff = 0;
1013 1.1.1.4 mrg mpfr_clear_flags ();
1014 1.1.1.4 mrg inex1 = mpfr_mul_2si (z1, zz, ediff, rnd);
1015 1.1.1.4 mrg flags1 = __gmpfr_flags;
1016 1.1.1.4 mrg mpfr_clear_flags ();
1017 1.1.1.4 mrg inex2 = mpfr_mul (z2, x, y, rnd);
1018 1.1.1.4 mrg flags2 = __gmpfr_flags;
1019 1.1.1.4 mrg if (!(mpfr_equal_p (z1, z2) &&
1020 1.1.1.4 mrg SAME_SIGN (inex1, inex2) &&
1021 1.1.1.4 mrg flags1 == flags2))
1022 1.1.1.4 mrg {
1023 1.1.1.4 mrg printf ("Error in small_prec on i = %d, j = %d\n", i, j);
1024 1.1.1.4 mrg printf ("r = 0x%x, xq = %d, yq = %d, zq = %d, rnd = %s\n",
1025 1.1.1.4 mrg r, xq, yq, zq, mpfr_print_rnd_mode (rnd));
1026 1.1.1.4 mrg printf ("x = ");
1027 1.1.1.4 mrg mpfr_dump (x);
1028 1.1.1.4 mrg printf ("y = ");
1029 1.1.1.4 mrg mpfr_dump (y);
1030 1.1.1.4 mrg printf ("ediff = %d\n", ediff);
1031 1.1.1.4 mrg printf ("zz = ");
1032 1.1.1.4 mrg mpfr_dump (zz);
1033 1.1.1.4 mrg printf ("Expected ");
1034 1.1.1.4 mrg mpfr_dump (z1);
1035 1.1.1.4 mrg printf ("with inex = %d and flags =", inex1);
1036 1.1.1.4 mrg flags_out (flags1);
1037 1.1.1.4 mrg printf ("Got ");
1038 1.1.1.4 mrg mpfr_dump (z2);
1039 1.1.1.4 mrg printf ("with inex = %d and flags =", inex2);
1040 1.1.1.4 mrg flags_out (flags2);
1041 1.1.1.4 mrg exit (1);
1042 1.1.1.4 mrg }
1043 1.1.1.4 mrg }
1044 1.1.1.4 mrg mpfr_clears (x, y, z1, z2, zz, (mpfr_ptr) 0);
1045 1.1.1.4 mrg
1046 1.1.1.4 mrg set_emin (emin);
1047 1.1.1.4 mrg set_emax (emax);
1048 1.1.1.4 mrg }
1049 1.1.1.4 mrg
1050 1.1.1.4 mrg /* check ax < __gmpfr_emin with rnd_mode == MPFR_RNDN, rb = 0 and sb <> 0 */
1051 1.1.1.4 mrg static void
1052 1.1.1.4 mrg test_underflow2 (void)
1053 1.1.1.4 mrg {
1054 1.1.1.4 mrg mpfr_t x, y, z;
1055 1.1.1.4 mrg mpfr_exp_t emin;
1056 1.1.1.4 mrg
1057 1.1.1.4 mrg emin = mpfr_get_emin ();
1058 1.1.1.4 mrg mpfr_set_emin (0);
1059 1.1.1.4 mrg
1060 1.1.1.4 mrg mpfr_init2 (x, 24);
1061 1.1.1.4 mrg mpfr_init2 (y, 24);
1062 1.1.1.4 mrg mpfr_init2 (z, 24);
1063 1.1.1.4 mrg
1064 1.1.1.4 mrg mpfr_set_ui_2exp (x, 12913, -14, MPFR_RNDN);
1065 1.1.1.4 mrg mpfr_set_ui_2exp (y, 5197, -13, MPFR_RNDN);
1066 1.1.1.4 mrg mpfr_clear_underflow ();
1067 1.1.1.4 mrg /* x*y = 0.0111111111111111111111111[01] thus underflow */
1068 1.1.1.4 mrg mpfr_mul (z, y, x, MPFR_RNDN);
1069 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -1) == 0);
1070 1.1.1.4 mrg MPFR_ASSERTN(mpfr_underflow_p ());
1071 1.1.1.4 mrg
1072 1.1.1.4 mrg mpfr_set_prec (x, 65);
1073 1.1.1.4 mrg mpfr_set_prec (y, 65);
1074 1.1.1.4 mrg mpfr_set_prec (z, 65);
1075 1.1.1.4 mrg mpfr_set_str_binary (x, "0.10011101000110111011111100101001111");
1076 1.1.1.4 mrg mpfr_set_str_binary (y, "0.110100001001000111000011011110011");
1077 1.1.1.4 mrg mpfr_clear_underflow ();
1078 1.1.1.4 mrg /* x*y = 0.011111111111111111111111111111111111111111111111111111111111111111[01] thus underflow */
1079 1.1.1.4 mrg mpfr_mul (z, y, x, MPFR_RNDN);
1080 1.1.1.4 mrg MPFR_ASSERTN(mpfr_cmp_ui_2exp (z, 1, -1) == 0);
1081 1.1.1.4 mrg MPFR_ASSERTN(mpfr_underflow_p ());
1082 1.1.1.4 mrg
1083 1.1.1.4 mrg mpfr_clear (y);
1084 1.1.1.4 mrg mpfr_clear (x);
1085 1.1.1.4 mrg mpfr_clear (z);
1086 1.1.1.4 mrg
1087 1.1.1.4 mrg mpfr_set_emin (emin);
1088 1.1.1.4 mrg }
1089 1.1.1.4 mrg
1090 1.1 mrg int
1091 1.1 mrg main (int argc, char *argv[])
1092 1.1 mrg {
1093 1.1 mrg tests_start_mpfr ();
1094 1.1 mrg
1095 1.1.1.4 mrg testall_rndf (9);
1096 1.1 mrg check_nans ();
1097 1.1 mrg check_exact ();
1098 1.1 mrg check_float ();
1099 1.1 mrg
1100 1.1 mrg check53("6.9314718055994530941514e-1", "0.0", MPFR_RNDZ, "0.0");
1101 1.1 mrg check53("0.0", "6.9314718055994530941514e-1", MPFR_RNDZ, "0.0");
1102 1.1 mrg check_sign();
1103 1.1 mrg check53("-4.165000000e4", "-0.00004801920768307322868063274915", MPFR_RNDN,
1104 1.1 mrg "2.0");
1105 1.1 mrg check53("2.71331408349172961467e-08", "-6.72658901114033715233e-165",
1106 1.1 mrg MPFR_RNDZ, "-1.8251348697787782844e-172");
1107 1.1 mrg check53("2.71331408349172961467e-08", "-6.72658901114033715233e-165",
1108 1.1 mrg MPFR_RNDA, "-1.8251348697787786e-172");
1109 1.1 mrg check53("0.31869277231188065", "0.88642843322303122", MPFR_RNDZ,
1110 1.1 mrg "2.8249833483992453642e-1");
1111 1.1 mrg check("8.47622108205396074254e-01", "3.24039313247872939883e-01", MPFR_RNDU,
1112 1.1 mrg 28, 45, 2, "0.375");
1113 1.1 mrg check("8.47622108205396074254e-01", "3.24039313247872939883e-01", MPFR_RNDA,
1114 1.1 mrg 28, 45, 2, "0.375");
1115 1.1 mrg check("2.63978122803639081440e-01", "6.8378615379333496093e-1", MPFR_RNDN,
1116 1.1 mrg 34, 23, 31, "0.180504585267044603");
1117 1.1 mrg check("1.0", "0.11835170935876249132", MPFR_RNDU, 6, 41, 36,
1118 1.1 mrg "0.1183517093595583");
1119 1.1 mrg check53("67108865.0", "134217729.0", MPFR_RNDN, "9.007199456067584e15");
1120 1.1 mrg check("1.37399642157394197284e-01", "2.28877275604219221350e-01", MPFR_RNDN,
1121 1.1 mrg 49, 15, 32, "0.0314472340833162888");
1122 1.1 mrg check("4.03160720978664954828e-01", "5.854828e-1"
1123 1.1 mrg /*"5.85483042917246621073e-01"*/, MPFR_RNDZ,
1124 1.1 mrg 51, 22, 32, "0.2360436821472831");
1125 1.1 mrg check("3.90798504668055102229e-14", "9.85394674650308388664e-04", MPFR_RNDN,
1126 1.1 mrg 46, 22, 12, "0.385027296503914762e-16");
1127 1.1 mrg check("4.58687081072827851358e-01", "2.20543551472118792844e-01", MPFR_RNDN,
1128 1.1 mrg 49, 3, 2, "0.09375");
1129 1.1 mrg check_max();
1130 1.1 mrg check_min();
1131 1.1.1.4 mrg small_prec ();
1132 1.1 mrg
1133 1.1 mrg check_regression ();
1134 1.1.1.4 mrg test_generic (MPFR_PREC_MIN, 500, 100);
1135 1.1 mrg
1136 1.1 mrg data_check ("data/mulpi", mpfr_mulpi, "mpfr_mulpi");
1137 1.1 mrg
1138 1.1.1.2 mrg valgrind20110503 ();
1139 1.1.1.3 mrg test_underflow (128);
1140 1.1.1.4 mrg bug20161209 ();
1141 1.1.1.4 mrg bug20161209a ();
1142 1.1.1.4 mrg bug20170602 ();
1143 1.1.1.4 mrg test_underflow2 ();
1144 1.1.1.2 mrg
1145 1.1 mrg tests_end_mpfr ();
1146 1.1 mrg return 0;
1147 1.1 mrg }
1148