Home | History | Annotate | Line # | Download | only in tests
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