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