Home | History | Annotate | Line # | Download | only in tests
tui_sub.c revision 1.1
      1  1.1  mrg /* Test file for mpfr_ui_sub.
      2  1.1  mrg 
      3  1.1  mrg Copyright 2000, 2001, 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 #include <float.h>
     26  1.1  mrg 
     27  1.1  mrg #include "mpfr-test.h"
     28  1.1  mrg 
     29  1.1  mrg static void
     30  1.1  mrg special (void)
     31  1.1  mrg {
     32  1.1  mrg   mpfr_t x, y, res;
     33  1.1  mrg   int inexact;
     34  1.1  mrg 
     35  1.1  mrg   mpfr_init (x);
     36  1.1  mrg   mpfr_init (y);
     37  1.1  mrg   mpfr_init (res);
     38  1.1  mrg 
     39  1.1  mrg   mpfr_set_prec (x, 24);
     40  1.1  mrg   mpfr_set_prec (y, 24);
     41  1.1  mrg   mpfr_set_str_binary (y, "0.111100110001011010111");
     42  1.1  mrg   inexact = mpfr_ui_sub (x, 1, y, MPFR_RNDN);
     43  1.1  mrg   if (inexact)
     44  1.1  mrg     {
     45  1.1  mrg       printf ("Wrong inexact flag: got %d, expected 0\n", inexact);
     46  1.1  mrg       exit (1);
     47  1.1  mrg     }
     48  1.1  mrg 
     49  1.1  mrg   mpfr_set_prec (x, 24);
     50  1.1  mrg   mpfr_set_prec (y, 24);
     51  1.1  mrg   mpfr_set_str_binary (y, "0.111100110001011010111");
     52  1.1  mrg   if ((inexact = mpfr_ui_sub (x, 38181761, y, MPFR_RNDN)) >= 0)
     53  1.1  mrg     {
     54  1.1  mrg       printf ("Wrong inexact flag: got %d, expected -1\n", inexact);
     55  1.1  mrg       exit (1);
     56  1.1  mrg     }
     57  1.1  mrg 
     58  1.1  mrg   mpfr_set_prec (x, 63);
     59  1.1  mrg   mpfr_set_prec (y, 63);
     60  1.1  mrg   mpfr_set_str_binary (y, "0.111110010010100100110101101010001001100101110001000101110111111E-1");
     61  1.1  mrg   if ((inexact = mpfr_ui_sub (x, 1541116494, y, MPFR_RNDN)) <= 0)
     62  1.1  mrg     {
     63  1.1  mrg       printf ("Wrong inexact flag: got %d, expected +1\n", inexact);
     64  1.1  mrg       exit (1);
     65  1.1  mrg     }
     66  1.1  mrg 
     67  1.1  mrg   mpfr_set_prec (x, 32);
     68  1.1  mrg   mpfr_set_prec (y, 32);
     69  1.1  mrg   mpfr_set_str_binary (y, "0.11011000110111010001011100011100E-1");
     70  1.1  mrg   if ((inexact = mpfr_ui_sub (x, 2000375416, y, MPFR_RNDN)) >= 0)
     71  1.1  mrg     {
     72  1.1  mrg       printf ("Wrong inexact flag: got %d, expected -1\n", inexact);
     73  1.1  mrg       exit (1);
     74  1.1  mrg     }
     75  1.1  mrg 
     76  1.1  mrg   mpfr_set_prec (x, 24);
     77  1.1  mrg   mpfr_set_prec (y, 24);
     78  1.1  mrg   mpfr_set_str_binary (y, "0.110011011001010011110111E-2");
     79  1.1  mrg   if ((inexact = mpfr_ui_sub (x, 927694848, y, MPFR_RNDN)) <= 0)
     80  1.1  mrg     {
     81  1.1  mrg       printf ("Wrong inexact flag: got %d, expected +1\n", inexact);
     82  1.1  mrg       exit (1);
     83  1.1  mrg     }
     84  1.1  mrg 
     85  1.1  mrg   /* bug found by Mathieu Dutour, 12 Apr 2001 */
     86  1.1  mrg   mpfr_set_prec (x, 5);
     87  1.1  mrg   mpfr_set_prec (y, 5);
     88  1.1  mrg   mpfr_set_prec (res, 5);
     89  1.1  mrg   mpfr_set_str_binary (x, "1e-12");
     90  1.1  mrg 
     91  1.1  mrg   mpfr_ui_sub (y, 1, x, MPFR_RNDD);
     92  1.1  mrg   mpfr_set_str_binary (res, "0.11111");
     93  1.1  mrg   if (mpfr_cmp (y, res))
     94  1.1  mrg     {
     95  1.1  mrg       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDD) for x=2^(-12)\nexpected 1.1111e-1, got ");
     96  1.1  mrg       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
     97  1.1  mrg       printf ("\n");
     98  1.1  mrg       exit (1);
     99  1.1  mrg     }
    100  1.1  mrg 
    101  1.1  mrg   mpfr_ui_sub (y, 1, x, MPFR_RNDU);
    102  1.1  mrg   mpfr_set_str_binary (res, "1.0");
    103  1.1  mrg   if (mpfr_cmp (y, res))
    104  1.1  mrg     {
    105  1.1  mrg       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDU) for x=2^(-12)\n"
    106  1.1  mrg               "expected 1.0, got ");
    107  1.1  mrg       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
    108  1.1  mrg       printf ("\n");
    109  1.1  mrg       exit (1);
    110  1.1  mrg     }
    111  1.1  mrg 
    112  1.1  mrg   mpfr_ui_sub (y, 1, x, MPFR_RNDN);
    113  1.1  mrg   mpfr_set_str_binary (res, "1.0");
    114  1.1  mrg   if (mpfr_cmp (y, res))
    115  1.1  mrg     {
    116  1.1  mrg       printf ("Error in mpfr_ui_sub (y, 1, x, MPFR_RNDN) for x=2^(-12)\n"
    117  1.1  mrg               "expected 1.0, got ");
    118  1.1  mrg       mpfr_out_str (stdout, 2, 0, y, MPFR_RNDN);
    119  1.1  mrg       printf ("\n");
    120  1.1  mrg       exit (1);
    121  1.1  mrg     }
    122  1.1  mrg 
    123  1.1  mrg   mpfr_set_prec (x, 10);
    124  1.1  mrg   mpfr_set_prec (y, 10);
    125  1.1  mrg   mpfr_urandomb (x, RANDS);
    126  1.1  mrg   mpfr_ui_sub (y, 0, x, MPFR_RNDN);
    127  1.1  mrg   if (MPFR_IS_ZERO(x))
    128  1.1  mrg     MPFR_ASSERTN(MPFR_IS_ZERO(y));
    129  1.1  mrg   else
    130  1.1  mrg     MPFR_ASSERTN(mpfr_cmpabs (x, y) == 0 && mpfr_sgn (x) != mpfr_sgn (y));
    131  1.1  mrg 
    132  1.1  mrg   mpfr_set_prec (x, 73);
    133  1.1  mrg   mpfr_set_str_binary (x, "0.1101111010101011011011100011010000000101110001011111001011011000101111101E-99");
    134  1.1  mrg   mpfr_ui_sub (x, 1, x, MPFR_RNDZ);
    135  1.1  mrg   mpfr_nextabove (x);
    136  1.1  mrg   MPFR_ASSERTN(mpfr_cmp_ui (x, 1) == 0);
    137  1.1  mrg 
    138  1.1  mrg   mpfr_clear (x);
    139  1.1  mrg   mpfr_clear (y);
    140  1.1  mrg   mpfr_clear (res);
    141  1.1  mrg }
    142  1.1  mrg 
    143  1.1  mrg /* checks that (y-x) gives the right results with 53 bits of precision */
    144  1.1  mrg static void
    145  1.1  mrg check (unsigned long y, const char *xs, mpfr_rnd_t rnd_mode, const char *zs)
    146  1.1  mrg {
    147  1.1  mrg   mpfr_t xx, zz;
    148  1.1  mrg 
    149  1.1  mrg   mpfr_inits2 (53, xx, zz, (mpfr_ptr) 0);
    150  1.1  mrg   mpfr_set_str1 (xx, xs);
    151  1.1  mrg   mpfr_ui_sub (zz, y, xx, rnd_mode);
    152  1.1  mrg   if (mpfr_cmp_str1 (zz, zs) )
    153  1.1  mrg     {
    154  1.1  mrg       printf ("expected difference is %s, got\n",zs);
    155  1.1  mrg       mpfr_out_str(stdout, 10, 0, zz, MPFR_RNDN);
    156  1.1  mrg       printf ("mpfr_ui_sub failed for y=%lu x=%s with rnd_mode=%s\n",
    157  1.1  mrg               y, xs, mpfr_print_rnd_mode (rnd_mode));
    158  1.1  mrg       exit (1);
    159  1.1  mrg     }
    160  1.1  mrg   mpfr_clears (xx, zz, (mpfr_ptr) 0);
    161  1.1  mrg }
    162  1.1  mrg 
    163  1.1  mrg /* if u = o(x-y), v = o(u-x), w = o(v+y), then x-y = u-w */
    164  1.1  mrg static void
    165  1.1  mrg check_two_sum (mpfr_prec_t p)
    166  1.1  mrg {
    167  1.1  mrg   unsigned int x;
    168  1.1  mrg   mpfr_t y, u, v, w;
    169  1.1  mrg   mpfr_rnd_t rnd;
    170  1.1  mrg   int inexact;
    171  1.1  mrg 
    172  1.1  mrg   mpfr_inits2 (p, y, u, v, w, (mpfr_ptr) 0);
    173  1.1  mrg   do
    174  1.1  mrg     {
    175  1.1  mrg       x = randlimb ();
    176  1.1  mrg     }
    177  1.1  mrg   while (x < 1);
    178  1.1  mrg   mpfr_urandomb (y, RANDS);
    179  1.1  mrg   rnd = MPFR_RNDN;
    180  1.1  mrg   inexact = mpfr_ui_sub (u, x, y, rnd);
    181  1.1  mrg   mpfr_sub_ui (v, u, x, rnd);
    182  1.1  mrg   mpfr_add (w, v, y, rnd);
    183  1.1  mrg   /* as u = (x-y) + w, we should have inexact and w of same sign */
    184  1.1  mrg   if (((inexact == 0) && mpfr_cmp_ui (w, 0)) ||
    185  1.1  mrg       ((inexact > 0) && (mpfr_cmp_ui (w, 0) <= 0)) ||
    186  1.1  mrg       ((inexact < 0) && (mpfr_cmp_ui (w, 0) >= 0)))
    187  1.1  mrg     {
    188  1.1  mrg       printf ("Wrong inexact flag for prec=%u, rnd=%s\n",
    189  1.1  mrg               (unsigned int) p, mpfr_print_rnd_mode (rnd));
    190  1.1  mrg       printf ("x=%u\n", x);
    191  1.1  mrg       printf ("y="); mpfr_print_binary(y); puts ("");
    192  1.1  mrg       printf ("u="); mpfr_print_binary(u); puts ("");
    193  1.1  mrg       printf ("v="); mpfr_print_binary(v); puts ("");
    194  1.1  mrg       printf ("w="); mpfr_print_binary(w); puts ("");
    195  1.1  mrg       printf ("inexact = %d\n", inexact);
    196  1.1  mrg       exit (1);
    197  1.1  mrg     }
    198  1.1  mrg   mpfr_clears (y, u, v, w, (mpfr_ptr) 0);
    199  1.1  mrg }
    200  1.1  mrg 
    201  1.1  mrg static void
    202  1.1  mrg check_nans (void)
    203  1.1  mrg {
    204  1.1  mrg   mpfr_t  x, y;
    205  1.1  mrg 
    206  1.1  mrg   mpfr_init2 (x, 123L);
    207  1.1  mrg   mpfr_init2 (y, 123L);
    208  1.1  mrg 
    209  1.1  mrg   /* 1 - nan == nan */
    210  1.1  mrg   mpfr_set_nan (x);
    211  1.1  mrg   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
    212  1.1  mrg   MPFR_ASSERTN (mpfr_nan_p (y));
    213  1.1  mrg 
    214  1.1  mrg   /* 1 - +inf == -inf */
    215  1.1  mrg   mpfr_set_inf (x, 1);
    216  1.1  mrg   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
    217  1.1  mrg   MPFR_ASSERTN (mpfr_inf_p (y));
    218  1.1  mrg   MPFR_ASSERTN (mpfr_sgn (y) < 0);
    219  1.1  mrg 
    220  1.1  mrg   /* 1 - -inf == +inf */
    221  1.1  mrg   mpfr_set_inf (x, -1);
    222  1.1  mrg   mpfr_ui_sub (y, 1L, x, MPFR_RNDN);
    223  1.1  mrg   MPFR_ASSERTN (mpfr_inf_p (y));
    224  1.1  mrg   MPFR_ASSERTN (mpfr_sgn (y) > 0);
    225  1.1  mrg 
    226  1.1  mrg   mpfr_clear (x);
    227  1.1  mrg   mpfr_clear (y);
    228  1.1  mrg }
    229  1.1  mrg 
    230  1.1  mrg /* Check mpfr_ui_sub with u = 0 (unsigned). */
    231  1.1  mrg static void check_neg (void)
    232  1.1  mrg {
    233  1.1  mrg   mpfr_t x, yneg, ysub;
    234  1.1  mrg   int i, s;
    235  1.1  mrg   int r;
    236  1.1  mrg 
    237  1.1  mrg   mpfr_init2 (x, 64);
    238  1.1  mrg   mpfr_init2 (yneg, 32);
    239  1.1  mrg   mpfr_init2 (ysub, 32);
    240  1.1  mrg 
    241  1.1  mrg   for (i = 0; i <= 25; i++)
    242  1.1  mrg     {
    243  1.1  mrg       mpfr_sqrt_ui (x, i, MPFR_RNDN);
    244  1.1  mrg       for (s = 0; s <= 1; s++)
    245  1.1  mrg         {
    246  1.1  mrg           RND_LOOP (r)
    247  1.1  mrg             {
    248  1.1  mrg               int tneg, tsub;
    249  1.1  mrg 
    250  1.1  mrg               tneg = mpfr_neg (yneg, x, (mpfr_rnd_t) r);
    251  1.1  mrg               tsub = mpfr_ui_sub (ysub, 0, x, (mpfr_rnd_t) r);
    252  1.1  mrg               MPFR_ASSERTN (mpfr_equal_p (yneg, ysub));
    253  1.1  mrg               MPFR_ASSERTN (!(MPFR_IS_POS (yneg) ^ MPFR_IS_POS (ysub)));
    254  1.1  mrg               MPFR_ASSERTN (tneg == tsub);
    255  1.1  mrg             }
    256  1.1  mrg           mpfr_neg (x, x, MPFR_RNDN);
    257  1.1  mrg         }
    258  1.1  mrg     }
    259  1.1  mrg 
    260  1.1  mrg   mpfr_clear (x);
    261  1.1  mrg   mpfr_clear (yneg);
    262  1.1  mrg   mpfr_clear (ysub);
    263  1.1  mrg }
    264  1.1  mrg 
    265  1.1  mrg int
    266  1.1  mrg main (int argc, char *argv[])
    267  1.1  mrg {
    268  1.1  mrg   mpfr_prec_t p;
    269  1.1  mrg   unsigned k;
    270  1.1  mrg 
    271  1.1  mrg   tests_start_mpfr ();
    272  1.1  mrg 
    273  1.1  mrg   check_nans ();
    274  1.1  mrg 
    275  1.1  mrg   special ();
    276  1.1  mrg   for (p=2; p<100; p++)
    277  1.1  mrg     for (k=0; k<100; k++)
    278  1.1  mrg       check_two_sum (p);
    279  1.1  mrg 
    280  1.1  mrg   check(1196426492, "1.4218093058435347e-3", MPFR_RNDN,
    281  1.1  mrg         "1.1964264919985781e9");
    282  1.1  mrg   check(1092583421, "-1.0880649218158844e9", MPFR_RNDN,
    283  1.1  mrg         "2.1806483428158845901e9");
    284  1.1  mrg   check(948002822, "1.22191250737771397120e+20", MPFR_RNDN,
    285  1.1  mrg         "-1.2219125073682338611e20");
    286  1.1  mrg   check(832100416, "4.68311314939691330000e-215", MPFR_RNDD,
    287  1.1  mrg         "8.3210041599999988079e8");
    288  1.1  mrg   check(1976245324, "1.25296395864546893357e+232", MPFR_RNDZ,
    289  1.1  mrg         "-1.2529639586454686577e232");
    290  1.1  mrg   check(2128997392, "-1.08496826129284207724e+187", MPFR_RNDU,
    291  1.1  mrg         "1.0849682612928422704e187");
    292  1.1  mrg   check(293607738, "-1.9967571564050541e-5", MPFR_RNDU,
    293  1.1  mrg         "2.9360773800002003e8");
    294  1.1  mrg   check(354270183, "2.9469161763489528e3", MPFR_RNDN,
    295  1.1  mrg         "3.5426723608382362e8");
    296  1.1  mrg 
    297  1.1  mrg   check_neg ();
    298  1.1  mrg 
    299  1.1  mrg   tests_end_mpfr ();
    300  1.1  mrg   return 0;
    301  1.1  mrg }
    302