tgmpop.c revision 1.1.1.1.8.1 1 1.1.1.1.8.1 tls /* Test file for mpfr_add_[q,z], mpfr_sub_[q,z], mpfr_div_[q,z],
2 1.1.1.1.8.1 tls mpfr_mul_[q,z], mpfr_cmp_[f,q,z]
3 1.1 mrg
4 1.1.1.1.8.1 tls Copyright 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013 Free Software Foundation, Inc.
5 1.1.1.1.8.1 tls Contributed by the AriC and Caramel projects, INRIA.
6 1.1 mrg
7 1.1 mrg This file is part of the GNU MPFR Library.
8 1.1 mrg
9 1.1 mrg The GNU MPFR Library is free software; you can redistribute it and/or modify
10 1.1 mrg it under the terms of the GNU Lesser General Public License as published by
11 1.1 mrg the Free Software Foundation; either version 3 of the License, or (at your
12 1.1 mrg option) any later version.
13 1.1 mrg
14 1.1 mrg The GNU MPFR Library is distributed in the hope that it will be useful, but
15 1.1 mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16 1.1 mrg or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17 1.1 mrg License for more details.
18 1.1 mrg
19 1.1 mrg You should have received a copy of the GNU Lesser General Public License
20 1.1 mrg along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
21 1.1 mrg http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
22 1.1 mrg 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
23 1.1 mrg
24 1.1 mrg #include <stdio.h>
25 1.1 mrg #include <stdlib.h>
26 1.1 mrg #include "mpfr-test.h"
27 1.1 mrg
28 1.1.1.1.8.1 tls #define CHECK_FOR(str, cond) \
29 1.1.1.1.8.1 tls if ((cond) == 0) { \
30 1.1.1.1.8.1 tls printf ("Special case error %s. Ternary value = %d, flags = %u\n", \
31 1.1.1.1.8.1 tls str, res, __gmpfr_flags); \
32 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y); \
33 1.1.1.1.8.1 tls printf ("X = "); mpfr_dump (x); \
34 1.1.1.1.8.1 tls printf ("Q = "); mpz_dump (mpq_numref(q)); \
35 1.1.1.1.8.1 tls printf (" /"); mpz_dump (mpq_denref(q)); \
36 1.1.1.1.8.1 tls exit (1); \
37 1.1.1.1.8.1 tls }
38 1.1.1.1.8.1 tls
39 1.1.1.1.8.1 tls #define CHECK_FORZ(str, cond) \
40 1.1.1.1.8.1 tls if ((cond) == 0) { \
41 1.1.1.1.8.1 tls printf ("Special case error %s. Ternary value = %d, flags = %u\n", \
42 1.1.1.1.8.1 tls str, res, __gmpfr_flags); \
43 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y); \
44 1.1.1.1.8.1 tls printf ("X = "); mpfr_dump (x); \
45 1.1.1.1.8.1 tls printf ("Z = "); mpz_dump (z); \
46 1.1.1.1.8.1 tls exit (1); \
47 1.1.1.1.8.1 tls }
48 1.1 mrg
49 1.1 mrg static void
50 1.1 mrg special (void)
51 1.1 mrg {
52 1.1 mrg mpfr_t x, y;
53 1.1.1.1.8.1 tls mpq_t q;
54 1.1.1.1.8.1 tls mpz_t z;
55 1.1 mrg int res = 0;
56 1.1 mrg
57 1.1 mrg mpfr_init (x);
58 1.1 mrg mpfr_init (y);
59 1.1.1.1.8.1 tls mpq_init (q);
60 1.1.1.1.8.1 tls mpz_init (z);
61 1.1 mrg
62 1.1 mrg /* cancellation in mpfr_add_q */
63 1.1 mrg mpfr_set_prec (x, 60);
64 1.1 mrg mpfr_set_prec (y, 20);
65 1.1.1.1.8.1 tls mpz_set_str (mpq_numref (q), "-187207494", 10);
66 1.1.1.1.8.1 tls mpz_set_str (mpq_denref (q), "5721", 10);
67 1.1 mrg mpfr_set_str_binary (x, "11111111101001011011100101100011011110010011100010000100001E-44");
68 1.1.1.1.8.1 tls mpfr_add_q (y, x, q, MPFR_RNDN);
69 1.1 mrg CHECK_FOR ("cancelation in add_q", mpfr_cmp_ui_2exp (y, 256783, -64) == 0);
70 1.1 mrg
71 1.1 mrg mpfr_set_prec (x, 19);
72 1.1 mrg mpfr_set_str_binary (x, "0.1011110101110011100E0");
73 1.1.1.1.8.1 tls mpz_set_str (mpq_numref (q), "187207494", 10);
74 1.1.1.1.8.1 tls mpz_set_str (mpq_denref (q), "5721", 10);
75 1.1 mrg mpfr_set_prec (y, 29);
76 1.1.1.1.8.1 tls mpfr_add_q (y, x, q, MPFR_RNDD);
77 1.1 mrg mpfr_set_prec (x, 29);
78 1.1 mrg mpfr_set_str_binary (x, "11111111101001110011010001001E-14");
79 1.1 mrg CHECK_FOR ("cancelation in add_q", mpfr_cmp (x,y) == 0);
80 1.1 mrg
81 1.1 mrg /* Inf */
82 1.1 mrg mpfr_set_inf (x, 1);
83 1.1.1.1.8.1 tls mpz_set_str (mpq_numref (q), "395877315", 10);
84 1.1.1.1.8.1 tls mpz_set_str (mpq_denref (q), "3508975966", 10);
85 1.1 mrg mpfr_set_prec (y, 118);
86 1.1.1.1.8.1 tls mpfr_add_q (y, x, q, MPFR_RNDU);
87 1.1 mrg CHECK_FOR ("inf", mpfr_inf_p (y) && mpfr_sgn (y) > 0);
88 1.1.1.1.8.1 tls mpfr_sub_q (y, x, q, MPFR_RNDU);
89 1.1 mrg CHECK_FOR ("inf", mpfr_inf_p (y) && mpfr_sgn (y) > 0);
90 1.1 mrg
91 1.1 mrg /* Nan */
92 1.1 mrg MPFR_SET_NAN (x);
93 1.1.1.1.8.1 tls mpfr_add_q (y, x, q, MPFR_RNDU);
94 1.1 mrg CHECK_FOR ("nan", mpfr_nan_p (y));
95 1.1.1.1.8.1 tls mpfr_sub_q (y, x, q, MPFR_RNDU);
96 1.1 mrg CHECK_FOR ("nan", mpfr_nan_p (y));
97 1.1 mrg
98 1.1 mrg /* Exact value */
99 1.1 mrg mpfr_set_prec (x, 60);
100 1.1 mrg mpfr_set_prec (y, 60);
101 1.1 mrg mpfr_set_str1 (x, "0.5");
102 1.1.1.1.8.1 tls mpz_set_str (mpq_numref (q), "3", 10);
103 1.1.1.1.8.1 tls mpz_set_str (mpq_denref (q), "2", 10);
104 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDU);
105 1.1 mrg CHECK_FOR ("0.5+3/2", mpfr_cmp_ui(y, 2)==0 && res==0);
106 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDU);
107 1.1 mrg CHECK_FOR ("0.5-3/2", mpfr_cmp_si(y, -1)==0 && res==0);
108 1.1 mrg
109 1.1 mrg /* Inf Rationnal */
110 1.1.1.1.8.1 tls mpq_set_ui (q, 1, 0);
111 1.1 mrg mpfr_set_str1 (x, "0.5");
112 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
113 1.1 mrg CHECK_FOR ("0.5+1/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
114 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
115 1.1 mrg CHECK_FOR ("0.5-1/0", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
116 1.1.1.1.8.1 tls mpq_set_si (q, -1, 0);
117 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
118 1.1 mrg CHECK_FOR ("0.5+ -1/0", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
119 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
120 1.1 mrg CHECK_FOR ("0.5- -1/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
121 1.1.1.1.8.1 tls res = mpfr_div_q (y, x, q, MPFR_RNDN);
122 1.1 mrg CHECK_FOR ("0.5 / (-1/0)", mpfr_zero_p (y) && MPFR_SIGN (y) < 0 && res == 0);
123 1.1.1.1.8.1 tls mpq_set_ui (q, 1, 0);
124 1.1.1.1.8.1 tls mpfr_set_inf (x, 1);
125 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
126 1.1.1.1.8.1 tls CHECK_FOR ("+Inf + +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
127 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
128 1.1.1.1.8.1 tls CHECK_FOR ("+Inf - +Inf", MPFR_IS_NAN (y) && res == 0);
129 1.1.1.1.8.1 tls mpfr_set_inf (x, -1);
130 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
131 1.1.1.1.8.1 tls CHECK_FOR ("-Inf + +Inf", MPFR_IS_NAN (y) && res == 0);
132 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
133 1.1.1.1.8.1 tls CHECK_FOR ("-Inf - +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
134 1.1.1.1.8.1 tls mpq_set_si (q, -1, 0);
135 1.1.1.1.8.1 tls mpfr_set_inf (x, 1);
136 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
137 1.1.1.1.8.1 tls CHECK_FOR ("+Inf + -Inf", MPFR_IS_NAN (y) && res == 0);
138 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
139 1.1.1.1.8.1 tls CHECK_FOR ("+Inf - -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
140 1.1.1.1.8.1 tls mpfr_set_inf (x, -1);
141 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
142 1.1.1.1.8.1 tls CHECK_FOR ("-Inf + -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
143 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
144 1.1.1.1.8.1 tls CHECK_FOR ("-Inf - -Inf", MPFR_IS_NAN (y) && res == 0);
145 1.1 mrg
146 1.1 mrg /* 0 */
147 1.1.1.1.8.1 tls mpq_set_ui (q, 0, 1);
148 1.1 mrg mpfr_set_ui (x, 42, MPFR_RNDN);
149 1.1.1.1.8.1 tls res = mpfr_add_q (y, x, q, MPFR_RNDN);
150 1.1 mrg CHECK_FOR ("42+0/1", mpfr_cmp_ui (y, 42) == 0 && res == 0);
151 1.1.1.1.8.1 tls res = mpfr_sub_q (y, x, q, MPFR_RNDN);
152 1.1 mrg CHECK_FOR ("42-0/1", mpfr_cmp_ui (y, 42) == 0 && res == 0);
153 1.1.1.1.8.1 tls res = mpfr_mul_q (y, x, q, MPFR_RNDN);
154 1.1 mrg CHECK_FOR ("42*0/1", mpfr_zero_p (y) && MPFR_SIGN (y) > 0 && res == 0);
155 1.1.1.1.8.1 tls mpfr_clear_flags ();
156 1.1.1.1.8.1 tls res = mpfr_div_q (y, x, q, MPFR_RNDN);
157 1.1.1.1.8.1 tls CHECK_FOR ("42/(0/1)", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
158 1.1.1.1.8.1 tls && mpfr_divby0_p ());
159 1.1.1.1.8.1 tls mpz_set_ui (z, 0);
160 1.1.1.1.8.1 tls mpfr_clear_flags ();
161 1.1.1.1.8.1 tls res = mpfr_div_z (y, x, z, MPFR_RNDN);
162 1.1.1.1.8.1 tls CHECK_FORZ ("42/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
163 1.1.1.1.8.1 tls && mpfr_divby0_p ());
164 1.1 mrg
165 1.1.1.1.8.1 tls mpz_clear (z);
166 1.1.1.1.8.1 tls mpq_clear (q);
167 1.1 mrg mpfr_clear (x);
168 1.1 mrg mpfr_clear (y);
169 1.1 mrg }
170 1.1 mrg
171 1.1 mrg static void
172 1.1 mrg check_for_zero (void)
173 1.1 mrg {
174 1.1 mrg /* Check that 0 is unsigned! */
175 1.1 mrg mpq_t q;
176 1.1 mrg mpz_t z;
177 1.1 mrg mpfr_t x;
178 1.1 mrg int r;
179 1.1 mrg mpfr_sign_t i;
180 1.1 mrg
181 1.1 mrg mpfr_init (x);
182 1.1 mrg mpz_init (z);
183 1.1 mrg mpq_init (q);
184 1.1 mrg
185 1.1 mrg mpz_set_ui (z, 0);
186 1.1 mrg mpq_set_ui (q, 0, 1);
187 1.1 mrg
188 1.1 mrg MPFR_SET_ZERO (x);
189 1.1 mrg RND_LOOP (r)
190 1.1 mrg {
191 1.1 mrg for (i = MPFR_SIGN_NEG ; i <= MPFR_SIGN_POS ;
192 1.1 mrg i+=MPFR_SIGN_POS-MPFR_SIGN_NEG)
193 1.1 mrg {
194 1.1 mrg MPFR_SET_SIGN(x, i);
195 1.1 mrg mpfr_add_z (x, x, z, (mpfr_rnd_t) r);
196 1.1 mrg if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
197 1.1 mrg {
198 1.1 mrg printf("GMP Zero errors for add_z & rnd=%s & s=%d\n",
199 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
200 1.1 mrg mpfr_dump (x);
201 1.1 mrg exit (1);
202 1.1 mrg }
203 1.1 mrg mpfr_sub_z (x, x, z, (mpfr_rnd_t) r);
204 1.1 mrg if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
205 1.1 mrg {
206 1.1 mrg printf("GMP Zero errors for sub_z & rnd=%s & s=%d\n",
207 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
208 1.1 mrg mpfr_dump (x);
209 1.1 mrg exit (1);
210 1.1 mrg }
211 1.1 mrg mpfr_mul_z (x, x, z, (mpfr_rnd_t) r);
212 1.1 mrg if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
213 1.1 mrg {
214 1.1 mrg printf("GMP Zero errors for mul_z & rnd=%s & s=%d\n",
215 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
216 1.1 mrg mpfr_dump (x);
217 1.1 mrg exit (1);
218 1.1 mrg }
219 1.1 mrg mpfr_add_q (x, x, q, (mpfr_rnd_t) r);
220 1.1 mrg if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
221 1.1 mrg {
222 1.1 mrg printf("GMP Zero errors for add_q & rnd=%s & s=%d\n",
223 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
224 1.1 mrg mpfr_dump (x);
225 1.1 mrg exit (1);
226 1.1 mrg }
227 1.1 mrg mpfr_sub_q (x, x, q, (mpfr_rnd_t) r);
228 1.1 mrg if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
229 1.1 mrg {
230 1.1 mrg printf("GMP Zero errors for sub_q & rnd=%s & s=%d\n",
231 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
232 1.1 mrg mpfr_dump (x);
233 1.1 mrg exit (1);
234 1.1 mrg }
235 1.1 mrg }
236 1.1 mrg }
237 1.1 mrg
238 1.1 mrg mpq_clear (q);
239 1.1 mrg mpz_clear (z);
240 1.1 mrg mpfr_clear (x);
241 1.1 mrg }
242 1.1 mrg
243 1.1 mrg static void
244 1.1 mrg test_cmp_z (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
245 1.1 mrg {
246 1.1 mrg mpfr_t x, z;
247 1.1 mrg mpz_t y;
248 1.1 mrg mpfr_prec_t p;
249 1.1 mrg int res1, res2;
250 1.1 mrg int n;
251 1.1 mrg
252 1.1 mrg mpfr_init (x);
253 1.1 mrg mpfr_init2 (z, MPFR_PREC_MIN);
254 1.1 mrg mpz_init (y);
255 1.1.1.1.8.1 tls
256 1.1.1.1.8.1 tls /* check the erange flag when x is NaN */
257 1.1.1.1.8.1 tls mpfr_set_nan (x);
258 1.1.1.1.8.1 tls mpz_set_ui (y, 17);
259 1.1.1.1.8.1 tls mpfr_clear_erangeflag ();
260 1.1.1.1.8.1 tls res1 = mpfr_cmp_z (x, y);
261 1.1.1.1.8.1 tls if (res1 != 0 || mpfr_erangeflag_p () == 0)
262 1.1.1.1.8.1 tls {
263 1.1.1.1.8.1 tls printf ("Error for mpfr_cmp_z (NaN, 17)\n");
264 1.1.1.1.8.1 tls printf ("Return value: expected 0, got %d\n", res1);
265 1.1.1.1.8.1 tls printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
266 1.1.1.1.8.1 tls exit (1);
267 1.1.1.1.8.1 tls }
268 1.1.1.1.8.1 tls
269 1.1 mrg for(p=pmin ; p < pmax ; p++)
270 1.1 mrg {
271 1.1 mrg mpfr_set_prec (x, p);
272 1.1 mrg for ( n = 0; n < nmax ; n++)
273 1.1 mrg {
274 1.1 mrg mpfr_urandomb (x, RANDS);
275 1.1 mrg mpz_urandomb (y, RANDS, 1024);
276 1.1 mrg if (!MPFR_IS_SINGULAR (x))
277 1.1 mrg {
278 1.1 mrg mpfr_sub_z (z, x, y, MPFR_RNDN);
279 1.1 mrg res1 = mpfr_sgn (z);
280 1.1 mrg res2 = mpfr_cmp_z (x, y);
281 1.1 mrg if (res1 != res2)
282 1.1 mrg {
283 1.1 mrg printf("Error for mpfr_cmp_z: res=%d sub_z gives %d\n",
284 1.1 mrg res2, res1);
285 1.1 mrg exit (1);
286 1.1 mrg }
287 1.1 mrg }
288 1.1 mrg }
289 1.1 mrg }
290 1.1 mrg mpz_clear (y);
291 1.1 mrg mpfr_clear (x);
292 1.1 mrg mpfr_clear (z);
293 1.1 mrg }
294 1.1 mrg
295 1.1 mrg static void
296 1.1 mrg test_cmp_q (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
297 1.1 mrg {
298 1.1 mrg mpfr_t x, z;
299 1.1 mrg mpq_t y;
300 1.1 mrg mpfr_prec_t p;
301 1.1 mrg int res1, res2;
302 1.1 mrg int n;
303 1.1 mrg
304 1.1 mrg mpfr_init (x);
305 1.1 mrg mpfr_init2 (z, MPFR_PREC_MIN);
306 1.1 mrg mpq_init (y);
307 1.1.1.1.8.1 tls
308 1.1.1.1.8.1 tls /* check the erange flag when x is NaN */
309 1.1.1.1.8.1 tls mpfr_set_nan (x);
310 1.1.1.1.8.1 tls mpq_set_ui (y, 17, 1);
311 1.1.1.1.8.1 tls mpfr_clear_erangeflag ();
312 1.1.1.1.8.1 tls res1 = mpfr_cmp_q (x, y);
313 1.1.1.1.8.1 tls if (res1 != 0 || mpfr_erangeflag_p () == 0)
314 1.1.1.1.8.1 tls {
315 1.1.1.1.8.1 tls printf ("Error for mpfr_cmp_q (NaN, 17)\n");
316 1.1.1.1.8.1 tls printf ("Return value: expected 0, got %d\n", res1);
317 1.1.1.1.8.1 tls printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
318 1.1.1.1.8.1 tls exit (1);
319 1.1.1.1.8.1 tls }
320 1.1.1.1.8.1 tls
321 1.1 mrg for(p=pmin ; p < pmax ; p++)
322 1.1 mrg {
323 1.1 mrg mpfr_set_prec (x, p);
324 1.1 mrg for (n = 0 ; n < nmax ; n++)
325 1.1 mrg {
326 1.1 mrg mpfr_urandomb (x, RANDS);
327 1.1 mrg mpq_set_ui (y, randlimb (), randlimb() );
328 1.1 mrg if (!MPFR_IS_SINGULAR (x))
329 1.1 mrg {
330 1.1 mrg mpfr_sub_q (z, x, y, MPFR_RNDN);
331 1.1 mrg res1 = mpfr_sgn (z);
332 1.1 mrg res2 = mpfr_cmp_q (x, y);
333 1.1 mrg if (res1 != res2)
334 1.1 mrg {
335 1.1 mrg printf("Error for mpfr_cmp_q: res=%d sub_z gives %d\n",
336 1.1 mrg res2, res1);
337 1.1 mrg exit (1);
338 1.1 mrg }
339 1.1 mrg }
340 1.1 mrg }
341 1.1 mrg }
342 1.1 mrg mpq_clear (y);
343 1.1 mrg mpfr_clear (x);
344 1.1 mrg mpfr_clear (z);
345 1.1 mrg }
346 1.1 mrg
347 1.1 mrg static void
348 1.1 mrg test_cmp_f (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
349 1.1 mrg {
350 1.1 mrg mpfr_t x, z;
351 1.1 mrg mpf_t y;
352 1.1 mrg mpfr_prec_t p;
353 1.1 mrg int res1, res2;
354 1.1 mrg int n;
355 1.1 mrg
356 1.1 mrg mpfr_init (x);
357 1.1 mrg mpfr_init2 (z, pmax+GMP_NUMB_BITS);
358 1.1 mrg mpf_init2 (y, MPFR_PREC_MIN);
359 1.1.1.1.8.1 tls
360 1.1.1.1.8.1 tls /* check the erange flag when x is NaN */
361 1.1.1.1.8.1 tls mpfr_set_nan (x);
362 1.1.1.1.8.1 tls mpf_set_ui (y, 17);
363 1.1.1.1.8.1 tls mpfr_clear_erangeflag ();
364 1.1.1.1.8.1 tls res1 = mpfr_cmp_f (x, y);
365 1.1.1.1.8.1 tls if (res1 != 0 || mpfr_erangeflag_p () == 0)
366 1.1.1.1.8.1 tls {
367 1.1.1.1.8.1 tls printf ("Error for mpfr_cmp_f (NaN, 17)\n");
368 1.1.1.1.8.1 tls printf ("Return value: expected 0, got %d\n", res1);
369 1.1.1.1.8.1 tls printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
370 1.1.1.1.8.1 tls exit (1);
371 1.1.1.1.8.1 tls }
372 1.1.1.1.8.1 tls
373 1.1 mrg for(p=pmin ; p < pmax ; p+=3)
374 1.1 mrg {
375 1.1 mrg mpfr_set_prec (x, p);
376 1.1 mrg mpf_set_prec (y, p);
377 1.1 mrg for ( n = 0; n < nmax ; n++)
378 1.1 mrg {
379 1.1 mrg mpfr_urandomb (x, RANDS);
380 1.1 mrg mpf_urandomb (y, RANDS, p);
381 1.1 mrg if (!MPFR_IS_SINGULAR (x))
382 1.1 mrg {
383 1.1 mrg mpfr_set_f (z, y, MPFR_RNDN);
384 1.1 mrg mpfr_sub (z, x, z, MPFR_RNDN);
385 1.1 mrg res1 = mpfr_sgn (z);
386 1.1 mrg res2 = mpfr_cmp_f (x, y);
387 1.1 mrg if (res1 != res2)
388 1.1 mrg {
389 1.1 mrg printf("Error for mpfr_cmp_f: res=%d sub gives %d\n",
390 1.1 mrg res2, res1);
391 1.1 mrg exit (1);
392 1.1 mrg }
393 1.1 mrg }
394 1.1 mrg }
395 1.1 mrg }
396 1.1 mrg mpf_clear (y);
397 1.1 mrg mpfr_clear (x);
398 1.1 mrg mpfr_clear (z);
399 1.1 mrg }
400 1.1 mrg
401 1.1 mrg static void
402 1.1 mrg test_specialz (int (*mpfr_func)(mpfr_ptr, mpfr_srcptr, mpz_srcptr, mpfr_rnd_t),
403 1.1 mrg void (*mpz_func)(mpz_ptr, mpz_srcptr, mpz_srcptr),
404 1.1 mrg const char *op)
405 1.1 mrg {
406 1.1 mrg mpfr_t x1, x2;
407 1.1 mrg mpz_t z1, z2;
408 1.1 mrg int res;
409 1.1 mrg
410 1.1 mrg mpfr_inits2 (128, x1, x2, (mpfr_ptr) 0);
411 1.1 mrg mpz_init (z1); mpz_init(z2);
412 1.1 mrg mpz_fac_ui (z1, 19); /* 19!+1 fits perfectly in a 128 bits mantissa */
413 1.1 mrg mpz_add_ui (z1, z1, 1);
414 1.1 mrg mpz_fac_ui (z2, 20); /* 20!+1 fits perfectly in a 128 bits mantissa */
415 1.1 mrg mpz_add_ui (z2, z2, 1);
416 1.1 mrg
417 1.1 mrg res = mpfr_set_z(x1, z1, MPFR_RNDN);
418 1.1 mrg if (res)
419 1.1 mrg {
420 1.1 mrg printf("Specialz %s: set_z1 error\n", op);
421 1.1 mrg exit(1);
422 1.1 mrg }
423 1.1 mrg mpfr_set_z (x2, z2, MPFR_RNDN);
424 1.1 mrg if (res)
425 1.1 mrg {
426 1.1 mrg printf("Specialz %s: set_z2 error\n", op);
427 1.1 mrg exit(1);
428 1.1 mrg }
429 1.1 mrg
430 1.1 mrg /* (19!+1) * (20!+1) fits in a 128 bits number */
431 1.1 mrg res = mpfr_func(x1, x1, z2, MPFR_RNDN);
432 1.1 mrg if (res)
433 1.1 mrg {
434 1.1 mrg printf("Specialz %s: wrong inexact flag.\n", op);
435 1.1 mrg exit(1);
436 1.1 mrg }
437 1.1 mrg mpz_func(z1, z1, z2);
438 1.1 mrg res = mpfr_set_z (x2, z1, MPFR_RNDN);
439 1.1 mrg if (res)
440 1.1 mrg {
441 1.1 mrg printf("Specialz %s: set_z2 error\n", op);
442 1.1 mrg exit(1);
443 1.1 mrg }
444 1.1 mrg if (mpfr_cmp(x1, x2))
445 1.1 mrg {
446 1.1 mrg printf("Specialz %s: results differ.\nx1=", op);
447 1.1 mrg mpfr_print_binary(x1);
448 1.1 mrg printf("\nx2=");
449 1.1 mrg mpfr_print_binary(x2);
450 1.1 mrg printf ("\nZ2=");
451 1.1 mrg mpz_out_str (stdout, 2, z1);
452 1.1 mrg putchar('\n');
453 1.1 mrg exit(1);
454 1.1 mrg }
455 1.1 mrg
456 1.1 mrg mpz_set_ui (z1, 1);
457 1.1 mrg mpz_set_ui (z2, 0);
458 1.1 mrg mpfr_set_ui (x1, 1, MPFR_RNDN);
459 1.1 mrg mpz_func (z1, z1, z2);
460 1.1 mrg res = mpfr_func(x1, x1, z2, MPFR_RNDN);
461 1.1 mrg mpfr_set_z (x2, z1, MPFR_RNDN);
462 1.1 mrg if (mpfr_cmp(x1, x2))
463 1.1 mrg {
464 1.1 mrg printf("Specialz %s: results differ(2).\nx1=", op);
465 1.1 mrg mpfr_print_binary(x1);
466 1.1 mrg printf("\nx2=");
467 1.1 mrg mpfr_print_binary(x2);
468 1.1 mrg putchar('\n');
469 1.1 mrg exit(1);
470 1.1 mrg }
471 1.1 mrg
472 1.1 mrg mpz_clear (z1); mpz_clear(z2);
473 1.1 mrg mpfr_clears (x1, x2, (mpfr_ptr) 0);
474 1.1 mrg }
475 1.1 mrg
476 1.1 mrg static void
477 1.1.1.1.8.1 tls test_special2z (int (*mpfr_func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
478 1.1.1.1.8.1 tls void (*mpz_func)(mpz_ptr, mpz_srcptr, mpz_srcptr),
479 1.1.1.1.8.1 tls const char *op)
480 1.1.1.1.8.1 tls {
481 1.1.1.1.8.1 tls mpfr_t x1, x2;
482 1.1.1.1.8.1 tls mpz_t z1, z2;
483 1.1.1.1.8.1 tls int res;
484 1.1.1.1.8.1 tls
485 1.1.1.1.8.1 tls mpfr_inits2 (128, x1, x2, (mpfr_ptr) 0);
486 1.1.1.1.8.1 tls mpz_init (z1); mpz_init(z2);
487 1.1.1.1.8.1 tls mpz_fac_ui (z1, 19); /* 19!+1 fits perfectly in a 128 bits mantissa */
488 1.1.1.1.8.1 tls mpz_add_ui (z1, z1, 1);
489 1.1.1.1.8.1 tls mpz_fac_ui (z2, 20); /* 20!+1 fits perfectly in a 128 bits mantissa */
490 1.1.1.1.8.1 tls mpz_add_ui (z2, z2, 1);
491 1.1.1.1.8.1 tls
492 1.1.1.1.8.1 tls res = mpfr_set_z(x1, z1, MPFR_RNDN);
493 1.1.1.1.8.1 tls if (res)
494 1.1.1.1.8.1 tls {
495 1.1.1.1.8.1 tls printf("Special2z %s: set_z1 error\n", op);
496 1.1.1.1.8.1 tls exit(1);
497 1.1.1.1.8.1 tls }
498 1.1.1.1.8.1 tls mpfr_set_z (x2, z2, MPFR_RNDN);
499 1.1.1.1.8.1 tls if (res)
500 1.1.1.1.8.1 tls {
501 1.1.1.1.8.1 tls printf("Special2z %s: set_z2 error\n", op);
502 1.1.1.1.8.1 tls exit(1);
503 1.1.1.1.8.1 tls }
504 1.1.1.1.8.1 tls
505 1.1.1.1.8.1 tls /* (19!+1) * (20!+1) fits in a 128 bits number */
506 1.1.1.1.8.1 tls res = mpfr_func(x1, z1, x2, MPFR_RNDN);
507 1.1.1.1.8.1 tls if (res)
508 1.1.1.1.8.1 tls {
509 1.1.1.1.8.1 tls printf("Special2z %s: wrong inexact flag.\n", op);
510 1.1.1.1.8.1 tls exit(1);
511 1.1.1.1.8.1 tls }
512 1.1.1.1.8.1 tls mpz_func(z1, z1, z2);
513 1.1.1.1.8.1 tls res = mpfr_set_z (x2, z1, MPFR_RNDN);
514 1.1.1.1.8.1 tls if (res)
515 1.1.1.1.8.1 tls {
516 1.1.1.1.8.1 tls printf("Special2z %s: set_z2 error\n", op);
517 1.1.1.1.8.1 tls exit(1);
518 1.1.1.1.8.1 tls }
519 1.1.1.1.8.1 tls if (mpfr_cmp(x1, x2))
520 1.1.1.1.8.1 tls {
521 1.1.1.1.8.1 tls printf("Special2z %s: results differ.\nx1=", op);
522 1.1.1.1.8.1 tls mpfr_print_binary(x1);
523 1.1.1.1.8.1 tls printf("\nx2=");
524 1.1.1.1.8.1 tls mpfr_print_binary(x2);
525 1.1.1.1.8.1 tls printf ("\nZ2=");
526 1.1.1.1.8.1 tls mpz_out_str (stdout, 2, z1);
527 1.1.1.1.8.1 tls putchar('\n');
528 1.1.1.1.8.1 tls exit(1);
529 1.1.1.1.8.1 tls }
530 1.1.1.1.8.1 tls
531 1.1.1.1.8.1 tls mpz_set_ui (z1, 0);
532 1.1.1.1.8.1 tls mpz_set_ui (z2, 1);
533 1.1.1.1.8.1 tls mpfr_set_ui (x2, 1, MPFR_RNDN);
534 1.1.1.1.8.1 tls res = mpfr_func(x1, z1, x2, MPFR_RNDN);
535 1.1.1.1.8.1 tls mpz_func (z1, z1, z2);
536 1.1.1.1.8.1 tls mpfr_set_z (x2, z1, MPFR_RNDN);
537 1.1.1.1.8.1 tls if (mpfr_cmp(x1, x2))
538 1.1.1.1.8.1 tls {
539 1.1.1.1.8.1 tls printf("Special2z %s: results differ(2).\nx1=", op);
540 1.1.1.1.8.1 tls mpfr_print_binary(x1);
541 1.1.1.1.8.1 tls printf("\nx2=");
542 1.1.1.1.8.1 tls mpfr_print_binary(x2);
543 1.1.1.1.8.1 tls putchar('\n');
544 1.1.1.1.8.1 tls exit(1);
545 1.1.1.1.8.1 tls }
546 1.1.1.1.8.1 tls
547 1.1.1.1.8.1 tls mpz_clear (z1); mpz_clear(z2);
548 1.1.1.1.8.1 tls mpfr_clears (x1, x2, (mpfr_ptr) 0);
549 1.1.1.1.8.1 tls }
550 1.1.1.1.8.1 tls
551 1.1.1.1.8.1 tls static void
552 1.1 mrg test_genericz (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
553 1.1 mrg int (*func)(mpfr_ptr, mpfr_srcptr, mpz_srcptr, mpfr_rnd_t),
554 1.1 mrg const char *op)
555 1.1 mrg {
556 1.1 mrg mpfr_prec_t prec;
557 1.1 mrg mpfr_t arg1, dst_big, dst_small, tmp;
558 1.1 mrg mpz_t arg2;
559 1.1 mrg mpfr_rnd_t rnd;
560 1.1 mrg int inexact, compare, compare2;
561 1.1 mrg unsigned int n;
562 1.1 mrg
563 1.1 mrg mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
564 1.1 mrg mpz_init (arg2);
565 1.1 mrg
566 1.1 mrg for (prec = p0; prec <= p1; prec++)
567 1.1 mrg {
568 1.1 mrg mpfr_set_prec (arg1, prec);
569 1.1 mrg mpfr_set_prec (tmp, prec);
570 1.1 mrg mpfr_set_prec (dst_small, prec);
571 1.1 mrg
572 1.1 mrg for (n=0; n<N; n++)
573 1.1 mrg {
574 1.1 mrg mpfr_urandomb (arg1, RANDS);
575 1.1 mrg mpz_urandomb (arg2, RANDS, 1024);
576 1.1 mrg rnd = RND_RAND ();
577 1.1 mrg mpfr_set_prec (dst_big, 2*prec);
578 1.1 mrg compare = func(dst_big, arg1, arg2, rnd);
579 1.1 mrg if (mpfr_can_round (dst_big, 2*prec, rnd, rnd, prec))
580 1.1 mrg {
581 1.1 mrg mpfr_set (tmp, dst_big, rnd);
582 1.1 mrg inexact = func(dst_small, arg1, arg2, rnd);
583 1.1 mrg if (mpfr_cmp (tmp, dst_small))
584 1.1 mrg {
585 1.1 mrg printf ("Results differ for prec=%u rnd_mode=%s and %s_z:\n"
586 1.1 mrg "arg1=",
587 1.1 mrg (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
588 1.1 mrg mpfr_print_binary (arg1);
589 1.1 mrg printf("\narg2=");
590 1.1.1.1.8.1 tls mpz_out_str (stdout, 10, arg2);
591 1.1.1.1.8.1 tls printf ("\ngot ");
592 1.1.1.1.8.1 tls mpfr_dump (dst_small);
593 1.1.1.1.8.1 tls printf ("expected ");
594 1.1.1.1.8.1 tls mpfr_dump (tmp);
595 1.1.1.1.8.1 tls printf ("approx ");
596 1.1.1.1.8.1 tls mpfr_dump (dst_big);
597 1.1.1.1.8.1 tls exit (1);
598 1.1.1.1.8.1 tls }
599 1.1.1.1.8.1 tls compare2 = mpfr_cmp (tmp, dst_big);
600 1.1.1.1.8.1 tls /* if rounding to nearest, cannot know the sign of t - f(x)
601 1.1.1.1.8.1 tls because of composed rounding: y = o(f(x)) and t = o(y) */
602 1.1.1.1.8.1 tls if (compare * compare2 >= 0)
603 1.1.1.1.8.1 tls compare = compare + compare2;
604 1.1.1.1.8.1 tls else
605 1.1.1.1.8.1 tls compare = inexact; /* cannot determine sign(t-f(x)) */
606 1.1.1.1.8.1 tls if (((inexact == 0) && (compare != 0)) ||
607 1.1.1.1.8.1 tls ((inexact > 0) && (compare <= 0)) ||
608 1.1.1.1.8.1 tls ((inexact < 0) && (compare >= 0)))
609 1.1.1.1.8.1 tls {
610 1.1.1.1.8.1 tls printf ("Wrong inexact flag for rnd=%s and %s_z:\n"
611 1.1.1.1.8.1 tls "expected %d, got %d\n",
612 1.1.1.1.8.1 tls mpfr_print_rnd_mode (rnd), op, compare, inexact);
613 1.1.1.1.8.1 tls printf ("\narg1="); mpfr_print_binary (arg1);
614 1.1.1.1.8.1 tls printf ("\narg2="); mpz_out_str(stdout, 2, arg2);
615 1.1.1.1.8.1 tls printf ("\ndstl="); mpfr_print_binary (dst_big);
616 1.1.1.1.8.1 tls printf ("\ndsts="); mpfr_print_binary (dst_small);
617 1.1.1.1.8.1 tls printf ("\ntmp ="); mpfr_dump (tmp);
618 1.1.1.1.8.1 tls exit (1);
619 1.1.1.1.8.1 tls }
620 1.1.1.1.8.1 tls }
621 1.1.1.1.8.1 tls }
622 1.1.1.1.8.1 tls }
623 1.1.1.1.8.1 tls
624 1.1.1.1.8.1 tls mpz_clear (arg2);
625 1.1.1.1.8.1 tls mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
626 1.1.1.1.8.1 tls }
627 1.1.1.1.8.1 tls
628 1.1.1.1.8.1 tls static void
629 1.1.1.1.8.1 tls test_generic2z (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
630 1.1.1.1.8.1 tls int (*func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
631 1.1.1.1.8.1 tls const char *op)
632 1.1.1.1.8.1 tls {
633 1.1.1.1.8.1 tls mpfr_prec_t prec;
634 1.1.1.1.8.1 tls mpfr_t arg1, dst_big, dst_small, tmp;
635 1.1.1.1.8.1 tls mpz_t arg2;
636 1.1.1.1.8.1 tls mpfr_rnd_t rnd;
637 1.1.1.1.8.1 tls int inexact, compare, compare2;
638 1.1.1.1.8.1 tls unsigned int n;
639 1.1.1.1.8.1 tls
640 1.1.1.1.8.1 tls mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
641 1.1.1.1.8.1 tls mpz_init (arg2);
642 1.1.1.1.8.1 tls
643 1.1.1.1.8.1 tls for (prec = p0; prec <= p1; prec++)
644 1.1.1.1.8.1 tls {
645 1.1.1.1.8.1 tls mpfr_set_prec (arg1, prec);
646 1.1.1.1.8.1 tls mpfr_set_prec (tmp, prec);
647 1.1.1.1.8.1 tls mpfr_set_prec (dst_small, prec);
648 1.1.1.1.8.1 tls
649 1.1.1.1.8.1 tls for (n=0; n<N; n++)
650 1.1.1.1.8.1 tls {
651 1.1.1.1.8.1 tls mpfr_urandomb (arg1, RANDS);
652 1.1.1.1.8.1 tls mpz_urandomb (arg2, RANDS, 1024);
653 1.1.1.1.8.1 tls rnd = RND_RAND ();
654 1.1.1.1.8.1 tls mpfr_set_prec (dst_big, 2*prec);
655 1.1.1.1.8.1 tls compare = func(dst_big, arg2, arg1, rnd);
656 1.1.1.1.8.1 tls if (mpfr_can_round (dst_big, 2*prec, rnd, rnd, prec))
657 1.1.1.1.8.1 tls {
658 1.1.1.1.8.1 tls mpfr_set (tmp, dst_big, rnd);
659 1.1.1.1.8.1 tls inexact = func(dst_small, arg2, arg1, rnd);
660 1.1.1.1.8.1 tls if (mpfr_cmp (tmp, dst_small))
661 1.1.1.1.8.1 tls {
662 1.1.1.1.8.1 tls printf ("Results differ for prec=%u rnd_mode=%s and %s_z:\n"
663 1.1.1.1.8.1 tls "arg1=",
664 1.1.1.1.8.1 tls (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
665 1.1.1.1.8.1 tls mpfr_print_binary (arg1);
666 1.1.1.1.8.1 tls printf("\narg2=");
667 1.1.1.1.8.1 tls mpz_out_str (stdout, 10, arg2);
668 1.1 mrg printf ("\ngot ");
669 1.1 mrg mpfr_dump (dst_small);
670 1.1 mrg printf ("expected ");
671 1.1 mrg mpfr_dump (tmp);
672 1.1 mrg printf ("approx ");
673 1.1 mrg mpfr_dump (dst_big);
674 1.1 mrg exit (1);
675 1.1 mrg }
676 1.1 mrg compare2 = mpfr_cmp (tmp, dst_big);
677 1.1 mrg /* if rounding to nearest, cannot know the sign of t - f(x)
678 1.1 mrg because of composed rounding: y = o(f(x)) and t = o(y) */
679 1.1 mrg if (compare * compare2 >= 0)
680 1.1 mrg compare = compare + compare2;
681 1.1 mrg else
682 1.1 mrg compare = inexact; /* cannot determine sign(t-f(x)) */
683 1.1 mrg if (((inexact == 0) && (compare != 0)) ||
684 1.1 mrg ((inexact > 0) && (compare <= 0)) ||
685 1.1 mrg ((inexact < 0) && (compare >= 0)))
686 1.1 mrg {
687 1.1 mrg printf ("Wrong inexact flag for rnd=%s and %s_z:\n"
688 1.1 mrg "expected %d, got %d\n",
689 1.1 mrg mpfr_print_rnd_mode (rnd), op, compare, inexact);
690 1.1 mrg printf ("\narg1="); mpfr_print_binary (arg1);
691 1.1 mrg printf ("\narg2="); mpz_out_str(stdout, 2, arg2);
692 1.1 mrg printf ("\ndstl="); mpfr_print_binary (dst_big);
693 1.1 mrg printf ("\ndsts="); mpfr_print_binary (dst_small);
694 1.1 mrg printf ("\ntmp ="); mpfr_dump (tmp);
695 1.1 mrg exit (1);
696 1.1 mrg }
697 1.1 mrg }
698 1.1 mrg }
699 1.1 mrg }
700 1.1 mrg
701 1.1 mrg mpz_clear (arg2);
702 1.1 mrg mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
703 1.1 mrg }
704 1.1 mrg
705 1.1 mrg static void
706 1.1 mrg test_genericq (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
707 1.1 mrg int (*func)(mpfr_ptr, mpfr_srcptr, mpq_srcptr, mpfr_rnd_t),
708 1.1 mrg const char *op)
709 1.1 mrg {
710 1.1 mrg mpfr_prec_t prec;
711 1.1 mrg mpfr_t arg1, dst_big, dst_small, tmp;
712 1.1 mrg mpq_t arg2;
713 1.1 mrg mpfr_rnd_t rnd;
714 1.1 mrg int inexact, compare, compare2;
715 1.1 mrg unsigned int n;
716 1.1 mrg
717 1.1 mrg mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
718 1.1 mrg mpq_init (arg2);
719 1.1 mrg
720 1.1 mrg for (prec = p0; prec <= p1; prec++)
721 1.1 mrg {
722 1.1 mrg mpfr_set_prec (arg1, prec);
723 1.1 mrg mpfr_set_prec (tmp, prec);
724 1.1 mrg mpfr_set_prec (dst_small, prec);
725 1.1 mrg
726 1.1 mrg for (n=0; n<N; n++)
727 1.1 mrg {
728 1.1 mrg mpfr_urandomb (arg1, RANDS);
729 1.1 mrg mpq_set_ui (arg2, randlimb (), randlimb() );
730 1.1 mrg mpq_canonicalize (arg2);
731 1.1 mrg rnd = RND_RAND ();
732 1.1 mrg mpfr_set_prec (dst_big, prec+10);
733 1.1 mrg compare = func(dst_big, arg1, arg2, rnd);
734 1.1 mrg if (mpfr_can_round (dst_big, prec+10, rnd, rnd, prec))
735 1.1 mrg {
736 1.1 mrg mpfr_set (tmp, dst_big, rnd);
737 1.1 mrg inexact = func(dst_small, arg1, arg2, rnd);
738 1.1 mrg if (mpfr_cmp (tmp, dst_small))
739 1.1 mrg {
740 1.1 mrg printf ("Results differ for prec=%u rnd_mode=%s and %s_q:\n"
741 1.1 mrg "arg1=",
742 1.1 mrg (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
743 1.1 mrg mpfr_print_binary (arg1);
744 1.1 mrg printf("\narg2=");
745 1.1 mrg mpq_out_str(stdout, 2, arg2);
746 1.1 mrg printf ("\ngot ");
747 1.1 mrg mpfr_print_binary (dst_small);
748 1.1 mrg printf ("\nexpected ");
749 1.1 mrg mpfr_print_binary (tmp);
750 1.1 mrg printf ("\napprox ");
751 1.1 mrg mpfr_print_binary (dst_big);
752 1.1 mrg putchar('\n');
753 1.1 mrg exit (1);
754 1.1 mrg }
755 1.1 mrg compare2 = mpfr_cmp (tmp, dst_big);
756 1.1 mrg /* if rounding to nearest, cannot know the sign of t - f(x)
757 1.1 mrg because of composed rounding: y = o(f(x)) and t = o(y) */
758 1.1 mrg if (compare * compare2 >= 0)
759 1.1 mrg compare = compare + compare2;
760 1.1 mrg else
761 1.1 mrg compare = inexact; /* cannot determine sign(t-f(x)) */
762 1.1 mrg if (((inexact == 0) && (compare != 0)) ||
763 1.1 mrg ((inexact > 0) && (compare <= 0)) ||
764 1.1 mrg ((inexact < 0) && (compare >= 0)))
765 1.1 mrg {
766 1.1 mrg printf ("Wrong inexact flag for rnd=%s and %s_q:\n"
767 1.1 mrg "expected %d, got %d",
768 1.1 mrg mpfr_print_rnd_mode (rnd), op, compare, inexact);
769 1.1 mrg printf ("\narg1="); mpfr_print_binary (arg1);
770 1.1 mrg printf ("\narg2="); mpq_out_str(stdout, 2, arg2);
771 1.1 mrg printf ("\ndstl="); mpfr_print_binary (dst_big);
772 1.1 mrg printf ("\ndsts="); mpfr_print_binary (dst_small);
773 1.1 mrg printf ("\ntmp ="); mpfr_print_binary (tmp);
774 1.1 mrg putchar('\n');
775 1.1 mrg exit (1);
776 1.1 mrg }
777 1.1 mrg }
778 1.1 mrg }
779 1.1 mrg }
780 1.1 mrg
781 1.1 mrg mpq_clear (arg2);
782 1.1 mrg mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
783 1.1 mrg }
784 1.1 mrg
785 1.1 mrg static void
786 1.1 mrg test_specialq (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
787 1.1 mrg int (*mpfr_func)(mpfr_ptr, mpfr_srcptr, mpq_srcptr, mpfr_rnd_t),
788 1.1 mrg void (*mpq_func)(mpq_ptr, mpq_srcptr, mpq_srcptr),
789 1.1 mrg const char *op)
790 1.1 mrg {
791 1.1 mrg mpfr_t fra, frb, frq;
792 1.1 mrg mpq_t q1, q2, qr;
793 1.1 mrg unsigned int n;
794 1.1 mrg mpfr_prec_t prec;
795 1.1 mrg
796 1.1 mrg for (prec = p0 ; prec < p1 ; prec++)
797 1.1 mrg {
798 1.1 mrg mpfr_inits2 (prec, fra, frb, frq, (mpfr_ptr) 0);
799 1.1 mrg mpq_init (q1); mpq_init(q2); mpq_init (qr);
800 1.1 mrg
801 1.1 mrg for( n = 0 ; n < N ; n++)
802 1.1 mrg {
803 1.1 mrg mpq_set_ui(q1, randlimb(), randlimb() );
804 1.1 mrg mpq_set_ui(q2, randlimb(), randlimb() );
805 1.1 mrg mpq_canonicalize (q1);
806 1.1 mrg mpq_canonicalize (q2);
807 1.1 mrg mpq_func (qr, q1, q2);
808 1.1 mrg mpfr_set_q (fra, q1, MPFR_RNDD);
809 1.1 mrg mpfr_func (fra, fra, q2, MPFR_RNDD);
810 1.1 mrg mpfr_set_q (frb, q1, MPFR_RNDU);
811 1.1 mrg mpfr_func (frb, frb, q2, MPFR_RNDU);
812 1.1 mrg mpfr_set_q (frq, qr, MPFR_RNDN);
813 1.1 mrg /* We should have fra <= qr <= frb */
814 1.1 mrg if ( (mpfr_cmp(fra, frq) > 0) || (mpfr_cmp (frq, frb) > 0))
815 1.1 mrg {
816 1.1.1.1.8.1 tls printf("Range error for prec=%lu and %s",
817 1.1.1.1.8.1 tls (unsigned long) prec, op);
818 1.1 mrg printf ("\nq1="); mpq_out_str(stdout, 2, q1);
819 1.1 mrg printf ("\nq2="); mpq_out_str(stdout, 2, q2);
820 1.1 mrg printf ("\nfr_dn="); mpfr_print_binary (fra);
821 1.1 mrg printf ("\nfr_q ="); mpfr_print_binary (frq);
822 1.1 mrg printf ("\nfr_up="); mpfr_print_binary (frb);
823 1.1 mrg putchar('\n');
824 1.1 mrg exit (1);
825 1.1 mrg }
826 1.1 mrg }
827 1.1 mrg
828 1.1 mrg mpq_clear (q1); mpq_clear (q2); mpq_clear (qr);
829 1.1 mrg mpfr_clears (fra, frb, frq, (mpfr_ptr) 0);
830 1.1 mrg }
831 1.1 mrg }
832 1.1 mrg
833 1.1.1.1.8.1 tls static void
834 1.1.1.1.8.1 tls bug_mul_q_20100810 (void)
835 1.1.1.1.8.1 tls {
836 1.1.1.1.8.1 tls mpfr_t x;
837 1.1.1.1.8.1 tls mpfr_t y;
838 1.1.1.1.8.1 tls mpq_t q;
839 1.1.1.1.8.1 tls int inexact;
840 1.1.1.1.8.1 tls
841 1.1.1.1.8.1 tls mpfr_init (x);
842 1.1.1.1.8.1 tls mpfr_init (y);
843 1.1.1.1.8.1 tls mpq_init (q);
844 1.1.1.1.8.1 tls
845 1.1.1.1.8.1 tls /* mpfr_mul_q: the inexact value must be set in case of overflow */
846 1.1.1.1.8.1 tls mpq_set_ui (q, 4096, 3);
847 1.1.1.1.8.1 tls mpfr_set_inf (x, +1);
848 1.1.1.1.8.1 tls mpfr_nextbelow (x);
849 1.1.1.1.8.1 tls inexact = mpfr_mul_q (y, x, q, MPFR_RNDU);
850 1.1.1.1.8.1 tls
851 1.1.1.1.8.1 tls if (inexact <= 0)
852 1.1.1.1.8.1 tls {
853 1.1.1.1.8.1 tls printf ("Overflow error in mpfr_mul_q. ");
854 1.1.1.1.8.1 tls printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
855 1.1.1.1.8.1 tls
856 1.1.1.1.8.1 tls exit (1);
857 1.1.1.1.8.1 tls }
858 1.1.1.1.8.1 tls if (!mpfr_inf_p (y))
859 1.1.1.1.8.1 tls {
860 1.1.1.1.8.1 tls printf ("Overflow error in mpfr_mul_q (y, x, q, MPFR_RNDD). ");
861 1.1.1.1.8.1 tls printf ("\nx = ");
862 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
863 1.1.1.1.8.1 tls printf ("\nq = ");
864 1.1.1.1.8.1 tls mpq_out_str (stdout, 10, q);
865 1.1.1.1.8.1 tls printf ("\ny = ");
866 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
867 1.1.1.1.8.1 tls printf (" (should be +infinity)\n");
868 1.1.1.1.8.1 tls
869 1.1.1.1.8.1 tls exit (1);
870 1.1.1.1.8.1 tls }
871 1.1.1.1.8.1 tls
872 1.1.1.1.8.1 tls mpq_clear (q);
873 1.1.1.1.8.1 tls mpfr_clear (y);
874 1.1.1.1.8.1 tls mpfr_clear (x);
875 1.1.1.1.8.1 tls }
876 1.1.1.1.8.1 tls
877 1.1.1.1.8.1 tls static void
878 1.1.1.1.8.1 tls bug_div_q_20100810 (void)
879 1.1.1.1.8.1 tls {
880 1.1.1.1.8.1 tls mpfr_t x;
881 1.1.1.1.8.1 tls mpfr_t y;
882 1.1.1.1.8.1 tls mpq_t q;
883 1.1.1.1.8.1 tls int inexact;
884 1.1.1.1.8.1 tls
885 1.1.1.1.8.1 tls mpfr_init (x);
886 1.1.1.1.8.1 tls mpfr_init (y);
887 1.1.1.1.8.1 tls mpq_init (q);
888 1.1.1.1.8.1 tls
889 1.1.1.1.8.1 tls /* mpfr_div_q: the inexact value must be set in case of overflow */
890 1.1.1.1.8.1 tls mpq_set_ui (q, 3, 4096);
891 1.1.1.1.8.1 tls mpfr_set_inf (x, +1);
892 1.1.1.1.8.1 tls mpfr_nextbelow (x);
893 1.1.1.1.8.1 tls inexact = mpfr_div_q (y, x, q, MPFR_RNDU);
894 1.1.1.1.8.1 tls
895 1.1.1.1.8.1 tls if (inexact <= 0)
896 1.1.1.1.8.1 tls {
897 1.1.1.1.8.1 tls printf ("Overflow error in mpfr_div_q. ");
898 1.1.1.1.8.1 tls printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
899 1.1.1.1.8.1 tls
900 1.1.1.1.8.1 tls exit (1);
901 1.1.1.1.8.1 tls }
902 1.1.1.1.8.1 tls if (!mpfr_inf_p (y))
903 1.1.1.1.8.1 tls {
904 1.1.1.1.8.1 tls printf ("Overflow error in mpfr_div_q (y, x, q, MPFR_RNDD). ");
905 1.1.1.1.8.1 tls printf ("\nx = ");
906 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
907 1.1.1.1.8.1 tls printf ("\nq = ");
908 1.1.1.1.8.1 tls mpq_out_str (stdout, 10, q);
909 1.1.1.1.8.1 tls printf ("\ny = ");
910 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
911 1.1.1.1.8.1 tls printf (" (should be +infinity)\n");
912 1.1.1.1.8.1 tls
913 1.1.1.1.8.1 tls exit (1);
914 1.1.1.1.8.1 tls }
915 1.1.1.1.8.1 tls
916 1.1.1.1.8.1 tls mpq_clear (q);
917 1.1.1.1.8.1 tls mpfr_clear (y);
918 1.1.1.1.8.1 tls mpfr_clear (x);
919 1.1.1.1.8.1 tls }
920 1.1.1.1.8.1 tls
921 1.1.1.1.8.1 tls static void
922 1.1.1.1.8.1 tls bug_mul_div_q_20100818 (void)
923 1.1.1.1.8.1 tls {
924 1.1.1.1.8.1 tls mpq_t qa, qb;
925 1.1.1.1.8.1 tls mpfr_t x1, x2, y1, y2, y3;
926 1.1.1.1.8.1 tls mpfr_exp_t emin, emax, e;
927 1.1.1.1.8.1 tls int inex;
928 1.1.1.1.8.1 tls int rnd;
929 1.1.1.1.8.1 tls
930 1.1.1.1.8.1 tls emin = mpfr_get_emin ();
931 1.1.1.1.8.1 tls emax = mpfr_get_emax ();
932 1.1.1.1.8.1 tls set_emin (MPFR_EMIN_MIN);
933 1.1.1.1.8.1 tls set_emax (MPFR_EMAX_MAX);
934 1.1.1.1.8.1 tls
935 1.1.1.1.8.1 tls mpq_init (qa);
936 1.1.1.1.8.1 tls mpq_init (qb);
937 1.1.1.1.8.1 tls mpfr_inits2 (32, x1, x2, y1, y2, y3, (mpfr_ptr) 0);
938 1.1.1.1.8.1 tls
939 1.1.1.1.8.1 tls mpq_set_ui (qa, 3, 17);
940 1.1.1.1.8.1 tls mpq_set_ui (qb, 17, 3);
941 1.1.1.1.8.1 tls inex = mpfr_set_ui (x1, 7, MPFR_RNDN);
942 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
943 1.1.1.1.8.1 tls
944 1.1.1.1.8.1 tls e = MPFR_EMAX_MAX - 3;
945 1.1.1.1.8.1 tls inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN); /* x2 = x1 * 2^e */
946 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
947 1.1.1.1.8.1 tls
948 1.1.1.1.8.1 tls RND_LOOP(rnd)
949 1.1.1.1.8.1 tls {
950 1.1.1.1.8.1 tls mpfr_mul_q (y1, x1, qa, (mpfr_rnd_t) rnd);
951 1.1.1.1.8.1 tls mpfr_div_q (y3, x1, qb, (mpfr_rnd_t) rnd);
952 1.1.1.1.8.1 tls MPFR_ASSERTN (mpfr_equal_p (y1, y3));
953 1.1.1.1.8.1 tls inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
954 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
955 1.1.1.1.8.1 tls inex = mpfr_mul (y3, y3, y1, MPFR_RNDN); /* y3 = y1 * 2^e */
956 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
957 1.1.1.1.8.1 tls mpfr_mul_q (y2, x2, qa, (mpfr_rnd_t) rnd);
958 1.1.1.1.8.1 tls if (! mpfr_equal_p (y2, y3))
959 1.1.1.1.8.1 tls {
960 1.1.1.1.8.1 tls printf ("Error 1 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
961 1.1.1.1.8.1 tls printf ("Expected "); mpfr_dump (y3);
962 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y2);
963 1.1.1.1.8.1 tls exit (1);
964 1.1.1.1.8.1 tls }
965 1.1.1.1.8.1 tls mpfr_div_q (y2, x2, qb, (mpfr_rnd_t) rnd);
966 1.1.1.1.8.1 tls if (! mpfr_equal_p (y2, y3))
967 1.1.1.1.8.1 tls {
968 1.1.1.1.8.1 tls printf ("Error 2 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
969 1.1.1.1.8.1 tls printf ("Expected "); mpfr_dump (y3);
970 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y2);
971 1.1.1.1.8.1 tls exit (1);
972 1.1.1.1.8.1 tls }
973 1.1.1.1.8.1 tls }
974 1.1.1.1.8.1 tls
975 1.1.1.1.8.1 tls e = MPFR_EMIN_MIN;
976 1.1.1.1.8.1 tls inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN); /* x2 = x1 * 2^e */
977 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
978 1.1.1.1.8.1 tls
979 1.1.1.1.8.1 tls RND_LOOP(rnd)
980 1.1.1.1.8.1 tls {
981 1.1.1.1.8.1 tls mpfr_div_q (y1, x1, qa, (mpfr_rnd_t) rnd);
982 1.1.1.1.8.1 tls mpfr_mul_q (y3, x1, qb, (mpfr_rnd_t) rnd);
983 1.1.1.1.8.1 tls MPFR_ASSERTN (mpfr_equal_p (y1, y3));
984 1.1.1.1.8.1 tls inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
985 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
986 1.1.1.1.8.1 tls inex = mpfr_mul (y3, y3, y1, MPFR_RNDN); /* y3 = y1 * 2^e */
987 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
988 1.1.1.1.8.1 tls mpfr_div_q (y2, x2, qa, (mpfr_rnd_t) rnd);
989 1.1.1.1.8.1 tls if (! mpfr_equal_p (y2, y3))
990 1.1.1.1.8.1 tls {
991 1.1.1.1.8.1 tls printf ("Error 3 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
992 1.1.1.1.8.1 tls printf ("Expected "); mpfr_dump (y3);
993 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y2);
994 1.1.1.1.8.1 tls exit (1);
995 1.1.1.1.8.1 tls }
996 1.1.1.1.8.1 tls mpfr_mul_q (y2, x2, qb, (mpfr_rnd_t) rnd);
997 1.1.1.1.8.1 tls if (! mpfr_equal_p (y2, y3))
998 1.1.1.1.8.1 tls {
999 1.1.1.1.8.1 tls printf ("Error 4 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
1000 1.1.1.1.8.1 tls printf ("Expected "); mpfr_dump (y3);
1001 1.1.1.1.8.1 tls printf ("Got "); mpfr_dump (y2);
1002 1.1.1.1.8.1 tls exit (1);
1003 1.1.1.1.8.1 tls }
1004 1.1.1.1.8.1 tls }
1005 1.1.1.1.8.1 tls
1006 1.1.1.1.8.1 tls mpq_clear (qa);
1007 1.1.1.1.8.1 tls mpq_clear (qb);
1008 1.1.1.1.8.1 tls mpfr_clears (x1, x2, y1, y2, y3, (mpfr_ptr) 0);
1009 1.1.1.1.8.1 tls
1010 1.1.1.1.8.1 tls set_emin (emin);
1011 1.1.1.1.8.1 tls set_emax (emax);
1012 1.1.1.1.8.1 tls }
1013 1.1.1.1.8.1 tls
1014 1.1.1.1.8.1 tls static void
1015 1.1.1.1.8.1 tls reduced_expo_range (void)
1016 1.1.1.1.8.1 tls {
1017 1.1.1.1.8.1 tls mpfr_t x;
1018 1.1.1.1.8.1 tls mpz_t z;
1019 1.1.1.1.8.1 tls mpq_t q;
1020 1.1.1.1.8.1 tls mpfr_exp_t emin;
1021 1.1.1.1.8.1 tls int inex;
1022 1.1.1.1.8.1 tls
1023 1.1.1.1.8.1 tls emin = mpfr_get_emin ();
1024 1.1.1.1.8.1 tls set_emin (4);
1025 1.1.1.1.8.1 tls
1026 1.1.1.1.8.1 tls mpfr_init2 (x, 32);
1027 1.1.1.1.8.1 tls
1028 1.1.1.1.8.1 tls mpz_init (z);
1029 1.1.1.1.8.1 tls mpfr_clear_flags ();
1030 1.1.1.1.8.1 tls inex = mpfr_set_ui (x, 17, MPFR_RNDN);
1031 1.1.1.1.8.1 tls MPFR_ASSERTN (inex == 0);
1032 1.1.1.1.8.1 tls mpz_set_ui (z, 3);
1033 1.1.1.1.8.1 tls inex = mpfr_mul_z (x, x, z, MPFR_RNDN);
1034 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 51) != 0)
1035 1.1.1.1.8.1 tls {
1036 1.1.1.1.8.1 tls printf ("Error 1 in reduce_expo_range: expected 51 with inex = 0,"
1037 1.1.1.1.8.1 tls " got\n");
1038 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
1039 1.1.1.1.8.1 tls printf ("with inex = %d\n", inex);
1040 1.1.1.1.8.1 tls exit (1);
1041 1.1.1.1.8.1 tls }
1042 1.1.1.1.8.1 tls inex = mpfr_div_z (x, x, z, MPFR_RNDN);
1043 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
1044 1.1.1.1.8.1 tls {
1045 1.1.1.1.8.1 tls printf ("Error 2 in reduce_expo_range: expected 17 with inex = 0,"
1046 1.1.1.1.8.1 tls " got\n");
1047 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
1048 1.1.1.1.8.1 tls printf ("with inex = %d\n", inex);
1049 1.1.1.1.8.1 tls exit (1);
1050 1.1.1.1.8.1 tls }
1051 1.1.1.1.8.1 tls inex = mpfr_add_z (x, x, z, MPFR_RNDN);
1052 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 20) != 0)
1053 1.1.1.1.8.1 tls {
1054 1.1.1.1.8.1 tls printf ("Error 3 in reduce_expo_range: expected 20 with inex = 0,"
1055 1.1.1.1.8.1 tls " got\n");
1056 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
1057 1.1.1.1.8.1 tls printf ("with inex = %d\n", inex);
1058 1.1.1.1.8.1 tls exit (1);
1059 1.1.1.1.8.1 tls }
1060 1.1.1.1.8.1 tls inex = mpfr_sub_z (x, x, z, MPFR_RNDN);
1061 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
1062 1.1.1.1.8.1 tls {
1063 1.1.1.1.8.1 tls printf ("Error 4 in reduce_expo_range: expected 17 with inex = 0,"
1064 1.1.1.1.8.1 tls " got\n");
1065 1.1.1.1.8.1 tls mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
1066 1.1.1.1.8.1 tls printf ("with inex = %d\n", inex);
1067 1.1.1.1.8.1 tls exit (1);
1068 1.1.1.1.8.1 tls }
1069 1.1.1.1.8.1 tls MPFR_ASSERTN (__gmpfr_flags == 0);
1070 1.1.1.1.8.1 tls if (mpfr_cmp_z (x, z) <= 0)
1071 1.1.1.1.8.1 tls {
1072 1.1.1.1.8.1 tls printf ("Error 5 in reduce_expo_range: expected a positive value.\n");
1073 1.1.1.1.8.1 tls exit (1);
1074 1.1.1.1.8.1 tls }
1075 1.1.1.1.8.1 tls mpz_clear (z);
1076 1.1.1.1.8.1 tls
1077 1.1.1.1.8.1 tls mpq_init (q);
1078 1.1.1.1.8.1 tls mpq_set_ui (q, 1, 1);
1079 1.1.1.1.8.1 tls mpfr_set_ui (x, 16, MPFR_RNDN);
1080 1.1.1.1.8.1 tls inex = mpfr_add_q (x, x, q, MPFR_RNDN);
1081 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
1082 1.1.1.1.8.1 tls {
1083 1.1.1.1.8.1 tls printf ("Error in reduce_expo_range for 16 + 1/1,"
1084 1.1.1.1.8.1 tls " got inex = %d and\nx = ", inex);
1085 1.1.1.1.8.1 tls mpfr_dump (x);
1086 1.1.1.1.8.1 tls exit (1);
1087 1.1.1.1.8.1 tls }
1088 1.1.1.1.8.1 tls inex = mpfr_sub_q (x, x, q, MPFR_RNDN);
1089 1.1.1.1.8.1 tls if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 16) != 0)
1090 1.1.1.1.8.1 tls {
1091 1.1.1.1.8.1 tls printf ("Error in reduce_expo_range for 17 - 1/1,"
1092 1.1.1.1.8.1 tls " got inex = %d and\nx = ", inex);
1093 1.1.1.1.8.1 tls mpfr_dump (x);
1094 1.1.1.1.8.1 tls exit (1);
1095 1.1.1.1.8.1 tls }
1096 1.1.1.1.8.1 tls mpq_clear (q);
1097 1.1.1.1.8.1 tls
1098 1.1.1.1.8.1 tls mpfr_clear (x);
1099 1.1.1.1.8.1 tls
1100 1.1.1.1.8.1 tls set_emin (emin);
1101 1.1.1.1.8.1 tls }
1102 1.1.1.1.8.1 tls
1103 1.1.1.1.8.1 tls static void
1104 1.1.1.1.8.1 tls addsubq_overflow_aux (mpfr_exp_t e)
1105 1.1.1.1.8.1 tls {
1106 1.1.1.1.8.1 tls mpfr_t x, y;
1107 1.1.1.1.8.1 tls mpq_t q;
1108 1.1.1.1.8.1 tls mpfr_exp_t emax;
1109 1.1.1.1.8.1 tls int inex;
1110 1.1.1.1.8.1 tls int rnd;
1111 1.1.1.1.8.1 tls int sign, sub;
1112 1.1.1.1.8.1 tls
1113 1.1.1.1.8.1 tls MPFR_ASSERTN (e <= LONG_MAX);
1114 1.1.1.1.8.1 tls emax = mpfr_get_emax ();
1115 1.1.1.1.8.1 tls set_emax (e);
1116 1.1.1.1.8.1 tls mpfr_inits2 (16, x, y, (mpfr_ptr) 0);
1117 1.1.1.1.8.1 tls mpq_init (q);
1118 1.1.1.1.8.1 tls
1119 1.1.1.1.8.1 tls mpfr_set_inf (x, 1);
1120 1.1.1.1.8.1 tls mpfr_nextbelow (x);
1121 1.1.1.1.8.1 tls mpq_set_ui (q, 1, 1);
1122 1.1.1.1.8.1 tls
1123 1.1.1.1.8.1 tls for (sign = 0; sign <= 1; sign++)
1124 1.1.1.1.8.1 tls {
1125 1.1.1.1.8.1 tls for (sub = 0; sub <= 1; sub++)
1126 1.1.1.1.8.1 tls {
1127 1.1.1.1.8.1 tls RND_LOOP(rnd)
1128 1.1.1.1.8.1 tls {
1129 1.1.1.1.8.1 tls unsigned int flags, ex_flags;
1130 1.1.1.1.8.1 tls int inf;
1131 1.1.1.1.8.1 tls
1132 1.1.1.1.8.1 tls inf = rnd == MPFR_RNDA ||
1133 1.1.1.1.8.1 tls rnd == (sign ? MPFR_RNDD : MPFR_RNDU);
1134 1.1.1.1.8.1 tls ex_flags = MPFR_FLAGS_INEXACT | (inf ? MPFR_FLAGS_OVERFLOW : 0);
1135 1.1.1.1.8.1 tls mpfr_clear_flags ();
1136 1.1.1.1.8.1 tls inex = sub ?
1137 1.1.1.1.8.1 tls mpfr_sub_q (y, x, q, (mpfr_rnd_t) rnd) :
1138 1.1.1.1.8.1 tls mpfr_add_q (y, x, q, (mpfr_rnd_t) rnd);
1139 1.1.1.1.8.1 tls flags = __gmpfr_flags;
1140 1.1.1.1.8.1 tls if (inex == 0 || flags != ex_flags ||
1141 1.1.1.1.8.1 tls (inf ? ! mpfr_inf_p (y) : ! mpfr_equal_p (x, y)))
1142 1.1.1.1.8.1 tls {
1143 1.1.1.1.8.1 tls printf ("Error in addsubq_overflow_aux(%ld),"
1144 1.1.1.1.8.1 tls " sign = %d, %s\n", (long) e, sign,
1145 1.1.1.1.8.1 tls mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
1146 1.1.1.1.8.1 tls printf ("Got inex = %d, y = ", inex);
1147 1.1.1.1.8.1 tls mpfr_dump (y);
1148 1.1.1.1.8.1 tls printf ("Expected flags:");
1149 1.1.1.1.8.1 tls flags_out (ex_flags);
1150 1.1.1.1.8.1 tls printf ("Got flags: ");
1151 1.1.1.1.8.1 tls flags_out (flags);
1152 1.1.1.1.8.1 tls exit (1);
1153 1.1.1.1.8.1 tls }
1154 1.1.1.1.8.1 tls }
1155 1.1.1.1.8.1 tls mpq_neg (q, q);
1156 1.1.1.1.8.1 tls }
1157 1.1.1.1.8.1 tls mpfr_neg (x, x, MPFR_RNDN);
1158 1.1.1.1.8.1 tls mpq_neg (q, q);
1159 1.1.1.1.8.1 tls }
1160 1.1.1.1.8.1 tls
1161 1.1.1.1.8.1 tls mpq_clear (q);
1162 1.1.1.1.8.1 tls mpfr_clears (x, y, (mpfr_ptr) 0);
1163 1.1.1.1.8.1 tls set_emax (emax);
1164 1.1.1.1.8.1 tls }
1165 1.1.1.1.8.1 tls
1166 1.1.1.1.8.1 tls static void
1167 1.1.1.1.8.1 tls addsubq_overflow (void)
1168 1.1.1.1.8.1 tls {
1169 1.1.1.1.8.1 tls addsubq_overflow_aux (4913);
1170 1.1.1.1.8.1 tls addsubq_overflow_aux (MPFR_EMAX_MAX);
1171 1.1.1.1.8.1 tls }
1172 1.1.1.1.8.1 tls
1173 1.1.1.1.8.1 tls static void
1174 1.1.1.1.8.1 tls coverage_mpfr_mul_q_20110218 (void)
1175 1.1.1.1.8.1 tls {
1176 1.1.1.1.8.1 tls mpfr_t cmp, res, op1;
1177 1.1.1.1.8.1 tls mpq_t op2;
1178 1.1.1.1.8.1 tls int status;
1179 1.1.1.1.8.1 tls
1180 1.1.1.1.8.1 tls mpfr_init2 (cmp, MPFR_PREC_MIN);
1181 1.1.1.1.8.1 tls mpfr_init2 (res, MPFR_PREC_MIN);
1182 1.1.1.1.8.1 tls mpfr_init_set_si (op1, 1, MPFR_RNDN);
1183 1.1.1.1.8.1 tls
1184 1.1.1.1.8.1 tls mpq_init (op2);
1185 1.1.1.1.8.1 tls mpq_set_si (op2, 0, 0);
1186 1.1.1.1.8.1 tls mpz_set_si (mpq_denref (op2), 0);
1187 1.1.1.1.8.1 tls
1188 1.1.1.1.8.1 tls status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
1189 1.1.1.1.8.1 tls
1190 1.1.1.1.8.1 tls if ((status != 0) || (mpfr_cmp (cmp, res) != 0))
1191 1.1.1.1.8.1 tls {
1192 1.1.1.1.8.1 tls printf ("Results differ %d.\nres=", status);
1193 1.1.1.1.8.1 tls mpfr_print_binary (res);
1194 1.1.1.1.8.1 tls printf ("\ncmp=");
1195 1.1.1.1.8.1 tls mpfr_print_binary (cmp);
1196 1.1.1.1.8.1 tls putchar ('\n');
1197 1.1.1.1.8.1 tls exit (1);
1198 1.1.1.1.8.1 tls }
1199 1.1.1.1.8.1 tls
1200 1.1.1.1.8.1 tls mpfr_set_si (op1, 1, MPFR_RNDN);
1201 1.1.1.1.8.1 tls mpq_set_si (op2, -1, 0);
1202 1.1.1.1.8.1 tls
1203 1.1.1.1.8.1 tls status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
1204 1.1.1.1.8.1 tls
1205 1.1.1.1.8.1 tls mpfr_set_inf (cmp, -1);
1206 1.1.1.1.8.1 tls if ((status != 0) || (mpfr_cmp(res, cmp) != 0))
1207 1.1.1.1.8.1 tls {
1208 1.1.1.1.8.1 tls printf ("mpfr_mul_q 1 * (-1/0) returned a wrong value :\n waiting for ");
1209 1.1.1.1.8.1 tls mpfr_print_binary (cmp);
1210 1.1.1.1.8.1 tls printf (" got ");
1211 1.1.1.1.8.1 tls mpfr_print_binary (res);
1212 1.1.1.1.8.1 tls printf ("\n trinary value is %d\n", status);
1213 1.1.1.1.8.1 tls exit (1);
1214 1.1.1.1.8.1 tls }
1215 1.1.1.1.8.1 tls
1216 1.1.1.1.8.1 tls mpq_clear (op2);
1217 1.1.1.1.8.1 tls mpfr_clear (op1);
1218 1.1.1.1.8.1 tls mpfr_clear (res);
1219 1.1.1.1.8.1 tls mpfr_clear (cmp);
1220 1.1.1.1.8.1 tls }
1221 1.1.1.1.8.1 tls
1222 1.1 mrg int
1223 1.1 mrg main (int argc, char *argv[])
1224 1.1 mrg {
1225 1.1 mrg tests_start_mpfr ();
1226 1.1 mrg
1227 1.1 mrg special ();
1228 1.1 mrg
1229 1.1 mrg test_specialz (mpfr_add_z, mpz_add, "add");
1230 1.1 mrg test_specialz (mpfr_sub_z, mpz_sub, "sub");
1231 1.1 mrg test_specialz (mpfr_mul_z, mpz_mul, "mul");
1232 1.1 mrg test_genericz (2, 100, 100, mpfr_add_z, "add");
1233 1.1 mrg test_genericz (2, 100, 100, mpfr_sub_z, "sub");
1234 1.1 mrg test_genericz (2, 100, 100, mpfr_mul_z, "mul");
1235 1.1 mrg test_genericz (2, 100, 100, mpfr_div_z, "div");
1236 1.1.1.1.8.1 tls test_special2z (mpfr_z_sub, mpz_sub, "sub");
1237 1.1.1.1.8.1 tls test_generic2z (2, 100, 100, mpfr_z_sub, "sub");
1238 1.1 mrg
1239 1.1 mrg test_genericq (2, 100, 100, mpfr_add_q, "add");
1240 1.1 mrg test_genericq (2, 100, 100, mpfr_sub_q, "sub");
1241 1.1 mrg test_genericq (2, 100, 100, mpfr_mul_q, "mul");
1242 1.1 mrg test_genericq (2, 100, 100, mpfr_div_q, "div");
1243 1.1 mrg test_specialq (2, 100, 100, mpfr_mul_q, mpq_mul, "mul");
1244 1.1 mrg test_specialq (2, 100, 100, mpfr_div_q, mpq_div, "div");
1245 1.1 mrg test_specialq (2, 100, 100, mpfr_add_q, mpq_add, "add");
1246 1.1 mrg test_specialq (2, 100, 100, mpfr_sub_q, mpq_sub, "sub");
1247 1.1 mrg
1248 1.1 mrg test_cmp_z (2, 100, 100);
1249 1.1 mrg test_cmp_q (2, 100, 100);
1250 1.1 mrg test_cmp_f (2, 100, 100);
1251 1.1 mrg
1252 1.1 mrg check_for_zero ();
1253 1.1 mrg
1254 1.1.1.1.8.1 tls bug_mul_q_20100810 ();
1255 1.1.1.1.8.1 tls bug_div_q_20100810 ();
1256 1.1.1.1.8.1 tls bug_mul_div_q_20100818 ();
1257 1.1.1.1.8.1 tls reduced_expo_range ();
1258 1.1.1.1.8.1 tls addsubq_overflow ();
1259 1.1.1.1.8.1 tls
1260 1.1.1.1.8.1 tls coverage_mpfr_mul_q_20110218 ();
1261 1.1.1.1.8.1 tls
1262 1.1 mrg tests_end_mpfr ();
1263 1.1 mrg return 0;
1264 1.1 mrg }
1265 1.1 mrg
1266