tfrac.c revision 1.1.1.1 1 1.1 mrg /* Test file for mpfr_frac.
2 1.1 mrg
3 1.1 mrg Copyright 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
4 1.1 mrg Contributed by the Arenaire and Cacao projects, INRIA.
5 1.1 mrg
6 1.1 mrg This file is part of the GNU MPFR Library.
7 1.1 mrg
8 1.1 mrg The GNU MPFR Library is free software; you can redistribute it and/or modify
9 1.1 mrg it under the terms of the GNU Lesser General Public License as published by
10 1.1 mrg the Free Software Foundation; either version 3 of the License, or (at your
11 1.1 mrg option) any later version.
12 1.1 mrg
13 1.1 mrg The GNU MPFR Library is distributed in the hope that it will be useful, but
14 1.1 mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 1.1 mrg or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 1.1 mrg License for more details.
17 1.1 mrg
18 1.1 mrg You should have received a copy of the GNU Lesser General Public License
19 1.1 mrg along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
20 1.1 mrg http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 1.1 mrg 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
22 1.1 mrg
23 1.1 mrg #include <stdio.h>
24 1.1 mrg #include <stdlib.h>
25 1.1 mrg
26 1.1 mrg #include "mpfr-test.h"
27 1.1 mrg
28 1.1 mrg #define PIP 70
29 1.1 mrg #define PFP 70
30 1.1 mrg #define PMAX (PIP+2*PFP)
31 1.1 mrg
32 1.1 mrg static void
33 1.1 mrg check0 (mpfr_ptr ip, mpfr_ptr fp, mpfr_prec_t prec, mpfr_rnd_t rnd)
34 1.1 mrg {
35 1.1 mrg mpfr_t sum, tmp, dst, fp2;
36 1.1 mrg int inex1, inex2;
37 1.1 mrg
38 1.1 mrg mpfr_init2 (sum, PMAX);
39 1.1 mrg mpfr_init2 (tmp, PMAX);
40 1.1 mrg mpfr_init2 (dst, prec);
41 1.1 mrg mpfr_init2 (fp2, prec);
42 1.1 mrg
43 1.1 mrg if (MPFR_SIGN (ip) != MPFR_SIGN (fp))
44 1.1 mrg {
45 1.1 mrg printf ("Internal error (1)\n");
46 1.1 mrg exit (1);
47 1.1 mrg }
48 1.1 mrg if (mpfr_add (sum, ip, fp, MPFR_RNDZ))
49 1.1 mrg {
50 1.1 mrg printf ("Wrong inexact flag in mpfr_add\n");
51 1.1 mrg exit (1);
52 1.1 mrg }
53 1.1 mrg if (MPFR_SIGN (sum) != MPFR_SIGN (fp))
54 1.1 mrg {
55 1.1 mrg printf ("Internal error (2)\n");
56 1.1 mrg exit (1);
57 1.1 mrg }
58 1.1 mrg
59 1.1 mrg inex1 = mpfr_frac (dst, sum, rnd);
60 1.1 mrg inex2 = mpfr_set (fp2, fp, rnd);
61 1.1 mrg if (inex1 != inex2)
62 1.1 mrg {
63 1.1 mrg printf ("Wrong inexact flag in mpfr_frac for\n");
64 1.1 mrg mpfr_out_str (stdout, 2, 0, sum, MPFR_RNDN);
65 1.1 mrg printf ("\nGot %d instead of %d\n", inex1, inex2);
66 1.1 mrg exit (1);
67 1.1 mrg }
68 1.1 mrg if (!mpfr_number_p (dst) ||
69 1.1 mrg MPFR_SIGN (dst) != MPFR_SIGN (fp2) ||
70 1.1 mrg mpfr_cmp (dst, fp2))
71 1.1 mrg {
72 1.1 mrg printf ("Error in mpfr_frac (y, x, %s) with\nx = ",
73 1.1 mrg mpfr_print_rnd_mode (rnd));
74 1.1 mrg mpfr_out_str (stdout, 2, 0, sum, MPFR_RNDN);
75 1.1 mrg printf ("\nGot ");
76 1.1 mrg mpfr_out_str (stdout, 2, 0, dst, MPFR_RNDN);
77 1.1 mrg printf ("\ninstead of ");
78 1.1 mrg mpfr_out_str (stdout, 2, 0, fp2, MPFR_RNDN);
79 1.1 mrg printf ("\n");
80 1.1 mrg exit (1);
81 1.1 mrg }
82 1.1 mrg
83 1.1 mrg if (prec == PMAX)
84 1.1 mrg {
85 1.1 mrg inex1 = mpfr_frac (sum, sum, rnd);
86 1.1 mrg if (inex1)
87 1.1 mrg {
88 1.1 mrg printf ("Wrong inexact flag in mpfr_frac\n");
89 1.1 mrg exit (1);
90 1.1 mrg }
91 1.1 mrg if (!mpfr_number_p (sum) ||
92 1.1 mrg MPFR_SIGN (sum) != MPFR_SIGN (fp) ||
93 1.1 mrg mpfr_cmp (sum, fp))
94 1.1 mrg {
95 1.1 mrg printf ("Error in mpfr_frac (x, x, %s) with\nx = ",
96 1.1 mrg mpfr_print_rnd_mode (rnd));
97 1.1 mrg mpfr_add (tmp, ip, fp, MPFR_RNDZ);
98 1.1 mrg mpfr_out_str (stdout, 2, 0, tmp, MPFR_RNDN);
99 1.1 mrg printf ("\nGot ");
100 1.1 mrg mpfr_out_str (stdout, 2, 0, sum, MPFR_RNDN);
101 1.1 mrg printf ("\ninstead of ");
102 1.1 mrg mpfr_out_str (stdout, 2, 0, fp, MPFR_RNDN);
103 1.1 mrg printf ("\n");
104 1.1 mrg exit (1);
105 1.1 mrg }
106 1.1 mrg }
107 1.1 mrg
108 1.1 mrg mpfr_clear (fp2);
109 1.1 mrg mpfr_clear (dst);
110 1.1 mrg mpfr_clear (tmp);
111 1.1 mrg mpfr_clear (sum);
112 1.1 mrg }
113 1.1 mrg
114 1.1 mrg static void
115 1.1 mrg check1 (mpfr_ptr ip, mpfr_ptr fp)
116 1.1 mrg {
117 1.1 mrg int rnd;
118 1.1 mrg
119 1.1 mrg for (rnd = 0; rnd < MPFR_RND_MAX ; rnd++)
120 1.1 mrg {
121 1.1 mrg check0 (ip, fp, PMAX, (mpfr_rnd_t) rnd);
122 1.1 mrg check0 (ip, fp, 70, (mpfr_rnd_t) rnd);
123 1.1 mrg mpfr_neg (fp, fp, MPFR_RNDN);
124 1.1 mrg mpfr_neg (ip, ip, MPFR_RNDN);
125 1.1 mrg check0 (ip, fp, PMAX, (mpfr_rnd_t) rnd);
126 1.1 mrg check0 (ip, fp, 70, (mpfr_rnd_t) rnd);
127 1.1 mrg mpfr_neg (fp, fp, MPFR_RNDN);
128 1.1 mrg mpfr_neg (ip, ip, MPFR_RNDN);
129 1.1 mrg }
130 1.1 mrg }
131 1.1 mrg
132 1.1 mrg static void
133 1.1 mrg special (void)
134 1.1 mrg {
135 1.1 mrg mpfr_t z, t;
136 1.1 mrg
137 1.1 mrg mpfr_init (z);
138 1.1 mrg mpfr_init (t);
139 1.1 mrg
140 1.1 mrg mpfr_set_nan (z);
141 1.1 mrg mpfr_frac (t, z, MPFR_RNDN);
142 1.1 mrg if (!mpfr_nan_p (t))
143 1.1 mrg {
144 1.1 mrg printf ("Error for frac(NaN)\n");
145 1.1 mrg exit (1);
146 1.1 mrg }
147 1.1 mrg
148 1.1 mrg mpfr_set_prec (z, 6);
149 1.1 mrg mpfr_set_prec (t, 3);
150 1.1 mrg
151 1.1 mrg mpfr_set_str_binary (z, "0.101101E3");
152 1.1 mrg mpfr_frac (t, z, MPFR_RNDN);
153 1.1 mrg mpfr_set_str_binary (z, "0.101");
154 1.1 mrg if (mpfr_cmp (t, z))
155 1.1 mrg {
156 1.1 mrg printf ("Error in frac(0.101101E3)\n");
157 1.1 mrg exit (1);
158 1.1 mrg }
159 1.1 mrg
160 1.1 mrg mpfr_set_prec (z, 34);
161 1.1 mrg mpfr_set_prec (t, 26);
162 1.1 mrg mpfr_set_str_binary (z, "0.101101010000010011110011001101E9");
163 1.1 mrg mpfr_frac (t, z, MPFR_RNDN);
164 1.1 mrg mpfr_set_str_binary (z, "0.000010011110011001101");
165 1.1 mrg if (mpfr_cmp (t, z))
166 1.1 mrg {
167 1.1 mrg printf ("Error in frac(0.101101010000010011110011001101E9)\n");
168 1.1 mrg exit (1);
169 1.1 mrg }
170 1.1 mrg
171 1.1 mrg mpfr_clear (z);
172 1.1 mrg mpfr_clear (t);
173 1.1 mrg }
174 1.1 mrg
175 1.1 mrg static void
176 1.1 mrg bug20090918 (void)
177 1.1 mrg {
178 1.1 mrg mpfr_t x, y, z;
179 1.1 mrg mp_limb_t y0;
180 1.1 mrg int inexy, inexz;
181 1.1 mrg int r, i;
182 1.1 mrg char *s[] = { "61680.352935791015625", "61680.999999" };
183 1.1 mrg mpfr_exp_t emin;
184 1.1 mrg
185 1.1 mrg emin = mpfr_get_emin ();
186 1.1 mrg mpfr_init2 (x, 32);
187 1.1 mrg mpfr_init2 (y, 13);
188 1.1 mrg
189 1.1 mrg for (i = 0; i <= 9; i++)
190 1.1 mrg {
191 1.1 mrg mpfr_set_str (x, s[i & 1], 10, MPFR_RNDZ);
192 1.1 mrg
193 1.1 mrg RND_LOOP(r)
194 1.1 mrg {
195 1.1 mrg set_emin ((i >> 1) - 3);
196 1.1 mrg inexy = mpfr_frac (y, x, (mpfr_rnd_t) r);
197 1.1 mrg set_emin (emin);
198 1.1 mrg y0 = MPFR_MANT(y)[0];
199 1.1 mrg while (y0 != 0 && (y0 >> 1) << 1 == y0)
200 1.1 mrg y0 >>= 1;
201 1.1 mrg if (y0 > 0x2000)
202 1.1 mrg {
203 1.1 mrg printf ("Error in bug20090918 (significand has more than"
204 1.1 mrg " 13 bits), i = %d, %s.\n", i,
205 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r));
206 1.1 mrg exit (1);
207 1.1 mrg }
208 1.1 mrg mpfr_init2 (z, 32);
209 1.1 mrg inexz = mpfr_frac (z, x, MPFR_RNDN);
210 1.1 mrg MPFR_ASSERTN (inexz == 0); /* exact */
211 1.1 mrg inexz = mpfr_prec_round (z, 13, (mpfr_rnd_t) r);
212 1.1 mrg set_emin ((i >> 1) - 3);
213 1.1 mrg inexz = mpfr_check_range (z, inexz, (mpfr_rnd_t) r);
214 1.1 mrg set_emin (emin);
215 1.1 mrg if (mpfr_cmp0 (y, z) != 0)
216 1.1 mrg {
217 1.1 mrg printf ("Error in bug20090918, i = %d, %s.\n", i,
218 1.1 mrg mpfr_print_rnd_mode ((mpfr_rnd_t) r));
219 1.1 mrg printf ("Expected ");
220 1.1 mrg mpfr_dump (z);
221 1.1 mrg printf ("Got ");
222 1.1 mrg mpfr_dump (y);
223 1.1 mrg exit (1);
224 1.1 mrg }
225 1.1 mrg if (! SAME_SIGN (inexy, inexz))
226 1.1 mrg {
227 1.1 mrg printf ("Incorrect ternary value in bug20090918, i = %d, %s.\n",
228 1.1 mrg i, mpfr_print_rnd_mode ((mpfr_rnd_t) r));
229 1.1 mrg printf ("Expected %d, got %d.\n", inexz, inexy);
230 1.1 mrg exit (1);
231 1.1 mrg }
232 1.1 mrg mpfr_clear (z);
233 1.1 mrg }
234 1.1 mrg }
235 1.1 mrg
236 1.1 mrg mpfr_clear (x);
237 1.1 mrg mpfr_clear (y);
238 1.1 mrg }
239 1.1 mrg
240 1.1 mrg #define TEST_FUNCTION mpfr_frac
241 1.1 mrg #include "tgeneric.c"
242 1.1 mrg
243 1.1 mrg int
244 1.1 mrg main (void)
245 1.1 mrg {
246 1.1 mrg mpfr_t ip, fp;
247 1.1 mrg int ni, nf1, nf2;
248 1.1 mrg
249 1.1 mrg tests_start_mpfr ();
250 1.1 mrg
251 1.1 mrg special ();
252 1.1 mrg
253 1.1 mrg mpfr_init2 (ip, PIP);
254 1.1 mrg mpfr_init2 (fp, PFP);
255 1.1 mrg
256 1.1 mrg for (ni = -1; ni < PIP; ni++)
257 1.1 mrg {
258 1.1 mrg if (ni <= 0)
259 1.1 mrg { /* ni + 1 */
260 1.1 mrg mpfr_set_si (ip, ni, MPFR_RNDN);
261 1.1 mrg mpfr_add_ui (ip, ip, 1, MPFR_RNDN);
262 1.1 mrg }
263 1.1 mrg else
264 1.1 mrg { /* 2^ni + 1 */
265 1.1 mrg mpfr_set_ui (ip, 1, MPFR_RNDN);
266 1.1 mrg mpfr_mul_2ui (ip, ip, ni, MPFR_RNDN);
267 1.1 mrg mpfr_add_ui (ip, ip, 1, MPFR_RNDN);
268 1.1 mrg }
269 1.1 mrg
270 1.1 mrg mpfr_set_ui (fp, 0, MPFR_RNDN);
271 1.1 mrg check1 (ip, fp);
272 1.1 mrg
273 1.1 mrg for (nf1 = 1; nf1 < PFP; nf1++)
274 1.1 mrg {
275 1.1 mrg mpfr_set_ui (fp, 1, MPFR_RNDN);
276 1.1 mrg mpfr_div_2ui (fp, fp, nf1, MPFR_RNDN);
277 1.1 mrg check1 (ip, fp);
278 1.1 mrg nf2 = 1 + (randlimb () % (PFP - 1));
279 1.1 mrg mpfr_set_ui (fp, 1, MPFR_RNDN);
280 1.1 mrg mpfr_div_2ui (fp, fp, nf2, MPFR_RNDN);
281 1.1 mrg mpfr_add_ui (fp, fp, 1, MPFR_RNDN);
282 1.1 mrg mpfr_div_2ui (fp, fp, nf1, MPFR_RNDN);
283 1.1 mrg check1 (ip, fp);
284 1.1 mrg }
285 1.1 mrg }
286 1.1 mrg
287 1.1 mrg mpfr_set_ui (ip, 1, MPFR_RNDN);
288 1.1 mrg mpfr_div_ui (ip, ip, 0, MPFR_RNDN);
289 1.1 mrg mpfr_set_ui (fp, 0, MPFR_RNDN);
290 1.1 mrg check1 (ip, fp); /* test infinities */
291 1.1 mrg
292 1.1 mrg mpfr_clear (ip);
293 1.1 mrg mpfr_clear (fp);
294 1.1 mrg
295 1.1 mrg bug20090918 ();
296 1.1 mrg
297 1.1 mrg test_generic (2, 1000, 10);
298 1.1 mrg
299 1.1 mrg tests_end_mpfr ();
300 1.1 mrg return 0;
301 1.1 mrg }
302