1/* Test file for mpfr_frac.
2
3Copyright 2002-2023 Free Software Foundation, Inc.
4Contributed by the AriC and Caramba projects, INRIA.
5
6This file is part of the GNU MPFR Library.
7
8The GNU MPFR Library is free software; you can redistribute it and/or modify
9it under the terms of the GNU Lesser General Public License as published by
10the Free Software Foundation; either version 3 of the License, or (at your
11option) any later version.
12
13The GNU MPFR Library is distributed in the hope that it will be useful, but
14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
16License for more details.
17
18You should have received a copy of the GNU Lesser General Public License
19along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
20https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
2151 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
29static void
30check0 (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
111static void
112check1 (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
129static void
130special (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
172static void
173bug20090918 (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
240int
241main (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