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