Home | History | Annotate | Line # | Download | only in tests
tgmpop.c revision 1.1.1.3
      1  1.1.1.2  mrg /* Test file for mpfr_add_[q,z], mpfr_sub_[q,z], mpfr_div_[q,z],
      2  1.1.1.2  mrg    mpfr_mul_[q,z], mpfr_cmp_[f,q,z]
      3      1.1  mrg 
      4  1.1.1.3  mrg Copyright 2004-2016 Free Software Foundation, Inc.
      5  1.1.1.3  mrg Contributed by the AriC and Caramba projects, INRIA.
      6      1.1  mrg 
      7      1.1  mrg This file is part of the GNU MPFR Library.
      8      1.1  mrg 
      9      1.1  mrg The GNU MPFR Library is free software; you can redistribute it and/or modify
     10      1.1  mrg it under the terms of the GNU Lesser General Public License as published by
     11      1.1  mrg the Free Software Foundation; either version 3 of the License, or (at your
     12      1.1  mrg option) any later version.
     13      1.1  mrg 
     14      1.1  mrg The GNU MPFR Library is distributed in the hope that it will be useful, but
     15      1.1  mrg WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
     16      1.1  mrg or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
     17      1.1  mrg License for more details.
     18      1.1  mrg 
     19      1.1  mrg You should have received a copy of the GNU Lesser General Public License
     20      1.1  mrg along with the GNU MPFR Library; see the file COPYING.LESSER.  If not, see
     21      1.1  mrg http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
     22      1.1  mrg 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
     23      1.1  mrg 
     24      1.1  mrg #include <stdio.h>
     25      1.1  mrg #include <stdlib.h>
     26      1.1  mrg #include "mpfr-test.h"
     27      1.1  mrg 
     28  1.1.1.2  mrg #define CHECK_FOR(str, cond)                                            \
     29  1.1.1.2  mrg   if ((cond) == 0) {                                                    \
     30  1.1.1.2  mrg     printf ("Special case error %s. Ternary value = %d, flags = %u\n",  \
     31  1.1.1.2  mrg             str, res, __gmpfr_flags);                                   \
     32  1.1.1.2  mrg     printf ("Got "); mpfr_dump (y);                                     \
     33  1.1.1.2  mrg     printf ("X = "); mpfr_dump (x);                                     \
     34  1.1.1.2  mrg     printf ("Q = "); mpz_dump (mpq_numref(q));                          \
     35  1.1.1.2  mrg     printf ("   /"); mpz_dump (mpq_denref(q));                          \
     36  1.1.1.2  mrg     exit (1);                                                           \
     37  1.1.1.2  mrg   }
     38  1.1.1.2  mrg 
     39  1.1.1.2  mrg #define CHECK_FORZ(str, cond)                                           \
     40  1.1.1.2  mrg   if ((cond) == 0) {                                                    \
     41  1.1.1.2  mrg     printf ("Special case error %s. Ternary value = %d, flags = %u\n",  \
     42  1.1.1.2  mrg             str, res, __gmpfr_flags);                                   \
     43  1.1.1.2  mrg     printf ("Got "); mpfr_dump (y);                                     \
     44  1.1.1.2  mrg     printf ("X = "); mpfr_dump (x);                                     \
     45  1.1.1.2  mrg     printf ("Z = "); mpz_dump (z);                                      \
     46  1.1.1.2  mrg     exit (1);                                                           \
     47  1.1.1.2  mrg   }
     48      1.1  mrg 
     49      1.1  mrg static void
     50      1.1  mrg special (void)
     51      1.1  mrg {
     52      1.1  mrg   mpfr_t x, y;
     53  1.1.1.2  mrg   mpq_t q;
     54  1.1.1.2  mrg   mpz_t z;
     55      1.1  mrg   int res = 0;
     56      1.1  mrg 
     57      1.1  mrg   mpfr_init (x);
     58      1.1  mrg   mpfr_init (y);
     59  1.1.1.2  mrg   mpq_init (q);
     60  1.1.1.2  mrg   mpz_init (z);
     61      1.1  mrg 
     62      1.1  mrg   /* cancellation in mpfr_add_q */
     63      1.1  mrg   mpfr_set_prec (x, 60);
     64      1.1  mrg   mpfr_set_prec (y, 20);
     65  1.1.1.2  mrg   mpz_set_str (mpq_numref (q), "-187207494", 10);
     66  1.1.1.2  mrg   mpz_set_str (mpq_denref (q), "5721", 10);
     67      1.1  mrg   mpfr_set_str_binary (x, "11111111101001011011100101100011011110010011100010000100001E-44");
     68  1.1.1.2  mrg   mpfr_add_q (y, x, q, MPFR_RNDN);
     69      1.1  mrg   CHECK_FOR ("cancelation in add_q", mpfr_cmp_ui_2exp (y, 256783, -64) == 0);
     70      1.1  mrg 
     71      1.1  mrg   mpfr_set_prec (x, 19);
     72      1.1  mrg   mpfr_set_str_binary (x, "0.1011110101110011100E0");
     73  1.1.1.2  mrg   mpz_set_str (mpq_numref (q), "187207494", 10);
     74  1.1.1.2  mrg   mpz_set_str (mpq_denref (q), "5721", 10);
     75      1.1  mrg   mpfr_set_prec (y, 29);
     76  1.1.1.2  mrg   mpfr_add_q (y, x, q, MPFR_RNDD);
     77      1.1  mrg   mpfr_set_prec (x, 29);
     78      1.1  mrg   mpfr_set_str_binary (x, "11111111101001110011010001001E-14");
     79      1.1  mrg   CHECK_FOR ("cancelation in add_q", mpfr_cmp (x,y) == 0);
     80      1.1  mrg 
     81      1.1  mrg   /* Inf */
     82      1.1  mrg   mpfr_set_inf (x, 1);
     83  1.1.1.2  mrg   mpz_set_str (mpq_numref (q), "395877315", 10);
     84  1.1.1.2  mrg   mpz_set_str (mpq_denref (q), "3508975966", 10);
     85      1.1  mrg   mpfr_set_prec (y, 118);
     86  1.1.1.2  mrg   mpfr_add_q (y, x, q, MPFR_RNDU);
     87      1.1  mrg   CHECK_FOR ("inf", mpfr_inf_p (y) && mpfr_sgn (y) > 0);
     88  1.1.1.2  mrg   mpfr_sub_q (y, x, q, MPFR_RNDU);
     89      1.1  mrg   CHECK_FOR ("inf", mpfr_inf_p (y) && mpfr_sgn (y) > 0);
     90      1.1  mrg 
     91      1.1  mrg   /* Nan */
     92      1.1  mrg   MPFR_SET_NAN (x);
     93  1.1.1.2  mrg   mpfr_add_q (y, x, q, MPFR_RNDU);
     94      1.1  mrg   CHECK_FOR ("nan", mpfr_nan_p (y));
     95  1.1.1.2  mrg   mpfr_sub_q (y, x, q, MPFR_RNDU);
     96      1.1  mrg   CHECK_FOR ("nan", mpfr_nan_p (y));
     97      1.1  mrg 
     98      1.1  mrg   /* Exact value */
     99      1.1  mrg   mpfr_set_prec (x, 60);
    100      1.1  mrg   mpfr_set_prec (y, 60);
    101      1.1  mrg   mpfr_set_str1 (x, "0.5");
    102  1.1.1.2  mrg   mpz_set_str (mpq_numref (q), "3", 10);
    103  1.1.1.2  mrg   mpz_set_str (mpq_denref (q), "2", 10);
    104  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDU);
    105      1.1  mrg   CHECK_FOR ("0.5+3/2", mpfr_cmp_ui(y, 2)==0 && res==0);
    106  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDU);
    107      1.1  mrg   CHECK_FOR ("0.5-3/2", mpfr_cmp_si(y, -1)==0 && res==0);
    108      1.1  mrg 
    109      1.1  mrg   /* Inf Rationnal */
    110  1.1.1.2  mrg   mpq_set_ui (q, 1, 0);
    111      1.1  mrg   mpfr_set_str1 (x, "0.5");
    112  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    113      1.1  mrg   CHECK_FOR ("0.5+1/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    114  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    115      1.1  mrg   CHECK_FOR ("0.5-1/0", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    116  1.1.1.2  mrg   mpq_set_si (q, -1, 0);
    117  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    118      1.1  mrg   CHECK_FOR ("0.5+ -1/0", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    119  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    120      1.1  mrg   CHECK_FOR ("0.5- -1/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    121  1.1.1.2  mrg   res = mpfr_div_q (y, x, q, MPFR_RNDN);
    122      1.1  mrg   CHECK_FOR ("0.5 / (-1/0)", mpfr_zero_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    123  1.1.1.2  mrg   mpq_set_ui (q, 1, 0);
    124  1.1.1.2  mrg   mpfr_set_inf (x, 1);
    125  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    126  1.1.1.2  mrg   CHECK_FOR ("+Inf + +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    127  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    128  1.1.1.2  mrg   CHECK_FOR ("+Inf - +Inf", MPFR_IS_NAN (y) && res == 0);
    129  1.1.1.2  mrg   mpfr_set_inf (x, -1);
    130  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    131  1.1.1.2  mrg   CHECK_FOR ("-Inf + +Inf", MPFR_IS_NAN (y) && res == 0);
    132  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    133  1.1.1.2  mrg   CHECK_FOR ("-Inf - +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    134  1.1.1.2  mrg   mpq_set_si (q, -1, 0);
    135  1.1.1.2  mrg   mpfr_set_inf (x, 1);
    136  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    137  1.1.1.2  mrg   CHECK_FOR ("+Inf + -Inf", MPFR_IS_NAN (y) && res == 0);
    138  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    139  1.1.1.2  mrg   CHECK_FOR ("+Inf - -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    140  1.1.1.2  mrg   mpfr_set_inf (x, -1);
    141  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    142  1.1.1.2  mrg   CHECK_FOR ("-Inf + -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    143  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    144  1.1.1.2  mrg   CHECK_FOR ("-Inf - -Inf", MPFR_IS_NAN (y) && res == 0);
    145      1.1  mrg 
    146      1.1  mrg   /* 0 */
    147  1.1.1.2  mrg   mpq_set_ui (q, 0, 1);
    148      1.1  mrg   mpfr_set_ui (x, 42, MPFR_RNDN);
    149  1.1.1.2  mrg   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    150      1.1  mrg   CHECK_FOR ("42+0/1", mpfr_cmp_ui (y, 42) == 0 && res == 0);
    151  1.1.1.2  mrg   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    152      1.1  mrg   CHECK_FOR ("42-0/1", mpfr_cmp_ui (y, 42) == 0 && res == 0);
    153  1.1.1.2  mrg   res = mpfr_mul_q (y, x, q, MPFR_RNDN);
    154      1.1  mrg   CHECK_FOR ("42*0/1", mpfr_zero_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    155  1.1.1.2  mrg   mpfr_clear_flags ();
    156  1.1.1.2  mrg   res = mpfr_div_q (y, x, q, MPFR_RNDN);
    157  1.1.1.2  mrg   CHECK_FOR ("42/(0/1)", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
    158  1.1.1.2  mrg              && mpfr_divby0_p ());
    159  1.1.1.2  mrg   mpz_set_ui (z, 0);
    160  1.1.1.2  mrg   mpfr_clear_flags ();
    161  1.1.1.2  mrg   res = mpfr_div_z (y, x, z, MPFR_RNDN);
    162  1.1.1.2  mrg   CHECK_FORZ ("42/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
    163  1.1.1.2  mrg               && mpfr_divby0_p ());
    164      1.1  mrg 
    165  1.1.1.2  mrg   mpz_clear (z);
    166  1.1.1.2  mrg   mpq_clear (q);
    167      1.1  mrg   mpfr_clear (x);
    168      1.1  mrg   mpfr_clear (y);
    169      1.1  mrg }
    170      1.1  mrg 
    171      1.1  mrg static void
    172      1.1  mrg check_for_zero (void)
    173      1.1  mrg {
    174      1.1  mrg   /* Check that 0 is unsigned! */
    175      1.1  mrg   mpq_t q;
    176      1.1  mrg   mpz_t z;
    177      1.1  mrg   mpfr_t x;
    178      1.1  mrg   int r;
    179      1.1  mrg   mpfr_sign_t i;
    180      1.1  mrg 
    181      1.1  mrg   mpfr_init (x);
    182      1.1  mrg   mpz_init (z);
    183      1.1  mrg   mpq_init (q);
    184      1.1  mrg 
    185      1.1  mrg   mpz_set_ui (z, 0);
    186      1.1  mrg   mpq_set_ui (q, 0, 1);
    187      1.1  mrg 
    188      1.1  mrg   MPFR_SET_ZERO (x);
    189      1.1  mrg   RND_LOOP (r)
    190      1.1  mrg     {
    191      1.1  mrg       for (i = MPFR_SIGN_NEG ; i <= MPFR_SIGN_POS ;
    192      1.1  mrg            i+=MPFR_SIGN_POS-MPFR_SIGN_NEG)
    193      1.1  mrg         {
    194      1.1  mrg           MPFR_SET_SIGN(x, i);
    195      1.1  mrg           mpfr_add_z (x, x, z, (mpfr_rnd_t) r);
    196      1.1  mrg           if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
    197      1.1  mrg             {
    198      1.1  mrg               printf("GMP Zero errors for add_z & rnd=%s & s=%d\n",
    199      1.1  mrg                      mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
    200      1.1  mrg               mpfr_dump (x);
    201      1.1  mrg               exit (1);
    202      1.1  mrg             }
    203      1.1  mrg           mpfr_sub_z (x, x, z, (mpfr_rnd_t) r);
    204      1.1  mrg           if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
    205      1.1  mrg             {
    206      1.1  mrg               printf("GMP Zero errors for sub_z & rnd=%s & s=%d\n",
    207      1.1  mrg                      mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
    208      1.1  mrg               mpfr_dump (x);
    209      1.1  mrg               exit (1);
    210      1.1  mrg             }
    211      1.1  mrg           mpfr_mul_z (x, x, z, (mpfr_rnd_t) r);
    212      1.1  mrg           if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
    213      1.1  mrg             {
    214      1.1  mrg               printf("GMP Zero errors for mul_z & rnd=%s & s=%d\n",
    215      1.1  mrg                      mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
    216      1.1  mrg               mpfr_dump (x);
    217      1.1  mrg               exit (1);
    218      1.1  mrg             }
    219      1.1  mrg           mpfr_add_q (x, x, q, (mpfr_rnd_t) r);
    220      1.1  mrg           if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
    221      1.1  mrg             {
    222      1.1  mrg               printf("GMP Zero errors for add_q & rnd=%s & s=%d\n",
    223      1.1  mrg                      mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
    224      1.1  mrg               mpfr_dump (x);
    225      1.1  mrg               exit (1);
    226      1.1  mrg             }
    227      1.1  mrg           mpfr_sub_q (x, x, q, (mpfr_rnd_t) r);
    228      1.1  mrg           if (!MPFR_IS_ZERO(x) || MPFR_SIGN(x)!=i)
    229      1.1  mrg             {
    230      1.1  mrg               printf("GMP Zero errors for sub_q & rnd=%s & s=%d\n",
    231      1.1  mrg                      mpfr_print_rnd_mode ((mpfr_rnd_t) r), i);
    232      1.1  mrg               mpfr_dump (x);
    233      1.1  mrg               exit (1);
    234      1.1  mrg              }
    235      1.1  mrg         }
    236      1.1  mrg     }
    237      1.1  mrg 
    238      1.1  mrg   mpq_clear (q);
    239      1.1  mrg   mpz_clear (z);
    240      1.1  mrg   mpfr_clear (x);
    241      1.1  mrg }
    242      1.1  mrg 
    243      1.1  mrg static void
    244      1.1  mrg test_cmp_z (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
    245      1.1  mrg {
    246      1.1  mrg   mpfr_t x, z;
    247      1.1  mrg   mpz_t  y;
    248      1.1  mrg   mpfr_prec_t p;
    249      1.1  mrg   int res1, res2;
    250      1.1  mrg   int n;
    251      1.1  mrg 
    252      1.1  mrg   mpfr_init (x);
    253      1.1  mrg   mpfr_init2 (z, MPFR_PREC_MIN);
    254      1.1  mrg   mpz_init (y);
    255  1.1.1.2  mrg 
    256  1.1.1.2  mrg   /* check the erange flag when x is NaN */
    257  1.1.1.2  mrg   mpfr_set_nan (x);
    258  1.1.1.2  mrg   mpz_set_ui (y, 17);
    259  1.1.1.2  mrg   mpfr_clear_erangeflag ();
    260  1.1.1.2  mrg   res1 = mpfr_cmp_z (x, y);
    261  1.1.1.2  mrg   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    262  1.1.1.2  mrg     {
    263  1.1.1.2  mrg       printf ("Error for mpfr_cmp_z (NaN, 17)\n");
    264  1.1.1.2  mrg       printf ("Return value: expected 0, got %d\n", res1);
    265  1.1.1.2  mrg       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    266  1.1.1.2  mrg       exit (1);
    267  1.1.1.2  mrg     }
    268  1.1.1.2  mrg 
    269      1.1  mrg   for(p=pmin ; p < pmax ; p++)
    270      1.1  mrg     {
    271      1.1  mrg       mpfr_set_prec (x, p);
    272      1.1  mrg       for ( n = 0; n < nmax ; n++)
    273      1.1  mrg         {
    274      1.1  mrg           mpfr_urandomb (x, RANDS);
    275      1.1  mrg           mpz_urandomb  (y, RANDS, 1024);
    276      1.1  mrg           if (!MPFR_IS_SINGULAR (x))
    277      1.1  mrg             {
    278      1.1  mrg               mpfr_sub_z (z, x, y, MPFR_RNDN);
    279      1.1  mrg               res1 = mpfr_sgn (z);
    280      1.1  mrg               res2 = mpfr_cmp_z (x, y);
    281      1.1  mrg               if (res1 != res2)
    282      1.1  mrg                 {
    283      1.1  mrg                   printf("Error for mpfr_cmp_z: res=%d sub_z gives %d\n",
    284      1.1  mrg                          res2, res1);
    285      1.1  mrg                   exit (1);
    286      1.1  mrg                 }
    287      1.1  mrg             }
    288      1.1  mrg         }
    289      1.1  mrg     }
    290      1.1  mrg   mpz_clear (y);
    291      1.1  mrg   mpfr_clear (x);
    292      1.1  mrg   mpfr_clear (z);
    293      1.1  mrg }
    294      1.1  mrg 
    295      1.1  mrg static void
    296      1.1  mrg test_cmp_q (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
    297      1.1  mrg {
    298      1.1  mrg   mpfr_t x, z;
    299      1.1  mrg   mpq_t  y;
    300      1.1  mrg   mpfr_prec_t p;
    301      1.1  mrg   int res1, res2;
    302      1.1  mrg   int n;
    303      1.1  mrg 
    304      1.1  mrg   mpfr_init (x);
    305      1.1  mrg   mpfr_init2 (z, MPFR_PREC_MIN);
    306      1.1  mrg   mpq_init (y);
    307  1.1.1.2  mrg 
    308  1.1.1.2  mrg   /* check the erange flag when x is NaN */
    309  1.1.1.2  mrg   mpfr_set_nan (x);
    310  1.1.1.2  mrg   mpq_set_ui (y, 17, 1);
    311  1.1.1.2  mrg   mpfr_clear_erangeflag ();
    312  1.1.1.2  mrg   res1 = mpfr_cmp_q (x, y);
    313  1.1.1.2  mrg   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    314  1.1.1.2  mrg     {
    315  1.1.1.2  mrg       printf ("Error for mpfr_cmp_q (NaN, 17)\n");
    316  1.1.1.2  mrg       printf ("Return value: expected 0, got %d\n", res1);
    317  1.1.1.2  mrg       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    318  1.1.1.2  mrg       exit (1);
    319  1.1.1.2  mrg     }
    320  1.1.1.2  mrg 
    321      1.1  mrg   for(p=pmin ; p < pmax ; p++)
    322      1.1  mrg     {
    323      1.1  mrg       mpfr_set_prec (x, p);
    324      1.1  mrg       for (n = 0 ; n < nmax ; n++)
    325      1.1  mrg         {
    326      1.1  mrg           mpfr_urandomb (x, RANDS);
    327      1.1  mrg           mpq_set_ui (y, randlimb (), randlimb() );
    328      1.1  mrg           if (!MPFR_IS_SINGULAR (x))
    329      1.1  mrg             {
    330      1.1  mrg               mpfr_sub_q (z, x, y, MPFR_RNDN);
    331      1.1  mrg               res1 = mpfr_sgn (z);
    332      1.1  mrg               res2 = mpfr_cmp_q (x, y);
    333      1.1  mrg               if (res1 != res2)
    334      1.1  mrg                 {
    335      1.1  mrg                   printf("Error for mpfr_cmp_q: res=%d sub_z gives %d\n",
    336      1.1  mrg                          res2, res1);
    337      1.1  mrg                   exit (1);
    338      1.1  mrg                 }
    339      1.1  mrg             }
    340      1.1  mrg         }
    341      1.1  mrg     }
    342      1.1  mrg   mpq_clear (y);
    343      1.1  mrg   mpfr_clear (x);
    344      1.1  mrg   mpfr_clear (z);
    345      1.1  mrg }
    346      1.1  mrg 
    347      1.1  mrg static void
    348      1.1  mrg test_cmp_f (mpfr_prec_t pmin, mpfr_prec_t pmax, int nmax)
    349      1.1  mrg {
    350      1.1  mrg   mpfr_t x, z;
    351      1.1  mrg   mpf_t  y;
    352      1.1  mrg   mpfr_prec_t p;
    353      1.1  mrg   int res1, res2;
    354      1.1  mrg   int n;
    355      1.1  mrg 
    356      1.1  mrg   mpfr_init (x);
    357      1.1  mrg   mpfr_init2 (z, pmax+GMP_NUMB_BITS);
    358      1.1  mrg   mpf_init2 (y, MPFR_PREC_MIN);
    359  1.1.1.2  mrg 
    360  1.1.1.2  mrg   /* check the erange flag when x is NaN */
    361  1.1.1.2  mrg   mpfr_set_nan (x);
    362  1.1.1.2  mrg   mpf_set_ui (y, 17);
    363  1.1.1.2  mrg   mpfr_clear_erangeflag ();
    364  1.1.1.2  mrg   res1 = mpfr_cmp_f (x, y);
    365  1.1.1.2  mrg   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    366  1.1.1.2  mrg     {
    367  1.1.1.2  mrg       printf ("Error for mpfr_cmp_f (NaN, 17)\n");
    368  1.1.1.2  mrg       printf ("Return value: expected 0, got %d\n", res1);
    369  1.1.1.2  mrg       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    370  1.1.1.2  mrg       exit (1);
    371  1.1.1.2  mrg     }
    372  1.1.1.2  mrg 
    373      1.1  mrg   for(p=pmin ; p < pmax ; p+=3)
    374      1.1  mrg     {
    375      1.1  mrg       mpfr_set_prec (x, p);
    376      1.1  mrg       mpf_set_prec (y, p);
    377      1.1  mrg       for ( n = 0; n < nmax ; n++)
    378      1.1  mrg         {
    379      1.1  mrg           mpfr_urandomb (x, RANDS);
    380      1.1  mrg           mpf_urandomb  (y, RANDS, p);
    381      1.1  mrg           if (!MPFR_IS_SINGULAR (x))
    382      1.1  mrg             {
    383      1.1  mrg               mpfr_set_f (z, y, MPFR_RNDN);
    384      1.1  mrg               mpfr_sub   (z, x, z, MPFR_RNDN);
    385      1.1  mrg               res1 = mpfr_sgn (z);
    386      1.1  mrg               res2 = mpfr_cmp_f (x, y);
    387      1.1  mrg               if (res1 != res2)
    388      1.1  mrg                 {
    389      1.1  mrg                   printf("Error for mpfr_cmp_f: res=%d sub gives %d\n",
    390      1.1  mrg                          res2, res1);
    391      1.1  mrg                   exit (1);
    392      1.1  mrg                 }
    393      1.1  mrg             }
    394      1.1  mrg         }
    395      1.1  mrg     }
    396      1.1  mrg   mpf_clear (y);
    397      1.1  mrg   mpfr_clear (x);
    398      1.1  mrg   mpfr_clear (z);
    399      1.1  mrg }
    400      1.1  mrg 
    401      1.1  mrg static void
    402      1.1  mrg test_specialz (int (*mpfr_func)(mpfr_ptr, mpfr_srcptr, mpz_srcptr, mpfr_rnd_t),
    403      1.1  mrg                void (*mpz_func)(mpz_ptr, mpz_srcptr, mpz_srcptr),
    404      1.1  mrg                const char *op)
    405      1.1  mrg {
    406      1.1  mrg   mpfr_t x1, x2;
    407      1.1  mrg   mpz_t  z1, z2;
    408      1.1  mrg   int res;
    409      1.1  mrg 
    410      1.1  mrg   mpfr_inits2 (128, x1, x2, (mpfr_ptr) 0);
    411      1.1  mrg   mpz_init (z1); mpz_init(z2);
    412      1.1  mrg   mpz_fac_ui (z1, 19); /* 19!+1 fits perfectly in a 128 bits mantissa */
    413      1.1  mrg   mpz_add_ui (z1, z1, 1);
    414      1.1  mrg   mpz_fac_ui (z2, 20); /* 20!+1 fits perfectly in a 128 bits mantissa */
    415      1.1  mrg   mpz_add_ui (z2, z2, 1);
    416      1.1  mrg 
    417      1.1  mrg   res = mpfr_set_z(x1, z1, MPFR_RNDN);
    418      1.1  mrg   if (res)
    419      1.1  mrg     {
    420      1.1  mrg       printf("Specialz %s: set_z1 error\n", op);
    421      1.1  mrg       exit(1);
    422      1.1  mrg     }
    423      1.1  mrg   mpfr_set_z (x2, z2, MPFR_RNDN);
    424      1.1  mrg   if (res)
    425      1.1  mrg     {
    426      1.1  mrg       printf("Specialz %s: set_z2 error\n", op);
    427      1.1  mrg       exit(1);
    428      1.1  mrg     }
    429      1.1  mrg 
    430      1.1  mrg   /* (19!+1) * (20!+1) fits in a 128 bits number */
    431      1.1  mrg   res = mpfr_func(x1, x1, z2, MPFR_RNDN);
    432      1.1  mrg   if (res)
    433      1.1  mrg     {
    434      1.1  mrg       printf("Specialz %s: wrong inexact flag.\n", op);
    435      1.1  mrg       exit(1);
    436      1.1  mrg     }
    437      1.1  mrg   mpz_func(z1, z1, z2);
    438      1.1  mrg   res = mpfr_set_z (x2, z1, MPFR_RNDN);
    439      1.1  mrg   if (res)
    440      1.1  mrg     {
    441      1.1  mrg       printf("Specialz %s: set_z2 error\n", op);
    442      1.1  mrg       exit(1);
    443      1.1  mrg     }
    444      1.1  mrg   if (mpfr_cmp(x1, x2))
    445      1.1  mrg     {
    446      1.1  mrg       printf("Specialz %s: results differ.\nx1=", op);
    447      1.1  mrg       mpfr_print_binary(x1);
    448      1.1  mrg       printf("\nx2=");
    449      1.1  mrg       mpfr_print_binary(x2);
    450      1.1  mrg       printf ("\nZ2=");
    451      1.1  mrg       mpz_out_str (stdout, 2, z1);
    452      1.1  mrg       putchar('\n');
    453      1.1  mrg       exit(1);
    454      1.1  mrg     }
    455      1.1  mrg 
    456      1.1  mrg   mpz_set_ui (z1, 1);
    457      1.1  mrg   mpz_set_ui (z2, 0);
    458      1.1  mrg   mpfr_set_ui (x1, 1, MPFR_RNDN);
    459      1.1  mrg   mpz_func (z1, z1, z2);
    460      1.1  mrg   res = mpfr_func(x1, x1, z2, MPFR_RNDN);
    461      1.1  mrg   mpfr_set_z (x2, z1, MPFR_RNDN);
    462      1.1  mrg   if (mpfr_cmp(x1, x2))
    463      1.1  mrg     {
    464      1.1  mrg       printf("Specialz %s: results differ(2).\nx1=", op);
    465      1.1  mrg       mpfr_print_binary(x1);
    466      1.1  mrg       printf("\nx2=");
    467      1.1  mrg       mpfr_print_binary(x2);
    468      1.1  mrg       putchar('\n');
    469      1.1  mrg       exit(1);
    470      1.1  mrg     }
    471      1.1  mrg 
    472      1.1  mrg   mpz_clear (z1); mpz_clear(z2);
    473      1.1  mrg   mpfr_clears (x1, x2, (mpfr_ptr) 0);
    474      1.1  mrg }
    475      1.1  mrg 
    476      1.1  mrg static void
    477  1.1.1.2  mrg test_special2z (int (*mpfr_func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
    478  1.1.1.2  mrg                void (*mpz_func)(mpz_ptr, mpz_srcptr, mpz_srcptr),
    479  1.1.1.2  mrg                const char *op)
    480  1.1.1.2  mrg {
    481  1.1.1.2  mrg   mpfr_t x1, x2;
    482  1.1.1.2  mrg   mpz_t  z1, z2;
    483  1.1.1.2  mrg   int res;
    484  1.1.1.2  mrg 
    485  1.1.1.2  mrg   mpfr_inits2 (128, x1, x2, (mpfr_ptr) 0);
    486  1.1.1.2  mrg   mpz_init (z1); mpz_init(z2);
    487  1.1.1.2  mrg   mpz_fac_ui (z1, 19); /* 19!+1 fits perfectly in a 128 bits mantissa */
    488  1.1.1.2  mrg   mpz_add_ui (z1, z1, 1);
    489  1.1.1.2  mrg   mpz_fac_ui (z2, 20); /* 20!+1 fits perfectly in a 128 bits mantissa */
    490  1.1.1.2  mrg   mpz_add_ui (z2, z2, 1);
    491  1.1.1.2  mrg 
    492  1.1.1.2  mrg   res = mpfr_set_z(x1, z1, MPFR_RNDN);
    493  1.1.1.2  mrg   if (res)
    494  1.1.1.2  mrg     {
    495  1.1.1.2  mrg       printf("Special2z %s: set_z1 error\n", op);
    496  1.1.1.2  mrg       exit(1);
    497  1.1.1.2  mrg     }
    498  1.1.1.2  mrg   mpfr_set_z (x2, z2, MPFR_RNDN);
    499  1.1.1.2  mrg   if (res)
    500  1.1.1.2  mrg     {
    501  1.1.1.2  mrg       printf("Special2z %s: set_z2 error\n", op);
    502  1.1.1.2  mrg       exit(1);
    503  1.1.1.2  mrg     }
    504  1.1.1.2  mrg 
    505  1.1.1.2  mrg   /* (19!+1) * (20!+1) fits in a 128 bits number */
    506  1.1.1.2  mrg   res = mpfr_func(x1, z1, x2, MPFR_RNDN);
    507  1.1.1.2  mrg   if (res)
    508  1.1.1.2  mrg     {
    509  1.1.1.2  mrg       printf("Special2z %s: wrong inexact flag.\n", op);
    510  1.1.1.2  mrg       exit(1);
    511  1.1.1.2  mrg     }
    512  1.1.1.2  mrg   mpz_func(z1, z1, z2);
    513  1.1.1.2  mrg   res = mpfr_set_z (x2, z1, MPFR_RNDN);
    514  1.1.1.2  mrg   if (res)
    515  1.1.1.2  mrg     {
    516  1.1.1.2  mrg       printf("Special2z %s: set_z2 error\n", op);
    517  1.1.1.2  mrg       exit(1);
    518  1.1.1.2  mrg     }
    519  1.1.1.2  mrg   if (mpfr_cmp(x1, x2))
    520  1.1.1.2  mrg     {
    521  1.1.1.2  mrg       printf("Special2z %s: results differ.\nx1=", op);
    522  1.1.1.2  mrg       mpfr_print_binary(x1);
    523  1.1.1.2  mrg       printf("\nx2=");
    524  1.1.1.2  mrg       mpfr_print_binary(x2);
    525  1.1.1.2  mrg       printf ("\nZ2=");
    526  1.1.1.2  mrg       mpz_out_str (stdout, 2, z1);
    527  1.1.1.2  mrg       putchar('\n');
    528  1.1.1.2  mrg       exit(1);
    529  1.1.1.2  mrg     }
    530  1.1.1.2  mrg 
    531  1.1.1.2  mrg   mpz_set_ui (z1, 0);
    532  1.1.1.2  mrg   mpz_set_ui (z2, 1);
    533  1.1.1.2  mrg   mpfr_set_ui (x2, 1, MPFR_RNDN);
    534  1.1.1.2  mrg   res = mpfr_func(x1, z1, x2, MPFR_RNDN);
    535  1.1.1.2  mrg   mpz_func (z1, z1, z2);
    536  1.1.1.2  mrg   mpfr_set_z (x2, z1, MPFR_RNDN);
    537  1.1.1.2  mrg   if (mpfr_cmp(x1, x2))
    538  1.1.1.2  mrg     {
    539  1.1.1.2  mrg       printf("Special2z %s: results differ(2).\nx1=", op);
    540  1.1.1.2  mrg       mpfr_print_binary(x1);
    541  1.1.1.2  mrg       printf("\nx2=");
    542  1.1.1.2  mrg       mpfr_print_binary(x2);
    543  1.1.1.2  mrg       putchar('\n');
    544  1.1.1.2  mrg       exit(1);
    545  1.1.1.2  mrg     }
    546  1.1.1.2  mrg 
    547  1.1.1.2  mrg   mpz_clear (z1); mpz_clear(z2);
    548  1.1.1.2  mrg   mpfr_clears (x1, x2, (mpfr_ptr) 0);
    549  1.1.1.2  mrg }
    550  1.1.1.2  mrg 
    551  1.1.1.2  mrg static void
    552      1.1  mrg test_genericz (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
    553      1.1  mrg                int (*func)(mpfr_ptr, mpfr_srcptr, mpz_srcptr, mpfr_rnd_t),
    554      1.1  mrg                const char *op)
    555      1.1  mrg {
    556      1.1  mrg   mpfr_prec_t prec;
    557      1.1  mrg   mpfr_t arg1, dst_big, dst_small, tmp;
    558      1.1  mrg   mpz_t  arg2;
    559      1.1  mrg   mpfr_rnd_t rnd;
    560      1.1  mrg   int inexact, compare, compare2;
    561      1.1  mrg   unsigned int n;
    562      1.1  mrg 
    563      1.1  mrg   mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    564      1.1  mrg   mpz_init (arg2);
    565      1.1  mrg 
    566      1.1  mrg   for (prec = p0; prec <= p1; prec++)
    567      1.1  mrg     {
    568      1.1  mrg       mpfr_set_prec (arg1, prec);
    569      1.1  mrg       mpfr_set_prec (tmp, prec);
    570      1.1  mrg       mpfr_set_prec (dst_small, prec);
    571      1.1  mrg 
    572      1.1  mrg       for (n=0; n<N; n++)
    573      1.1  mrg         {
    574      1.1  mrg           mpfr_urandomb (arg1, RANDS);
    575      1.1  mrg           mpz_urandomb (arg2, RANDS, 1024);
    576      1.1  mrg           rnd = RND_RAND ();
    577      1.1  mrg           mpfr_set_prec (dst_big, 2*prec);
    578      1.1  mrg           compare = func(dst_big, arg1, arg2, rnd);
    579      1.1  mrg           if (mpfr_can_round (dst_big, 2*prec, rnd, rnd, prec))
    580      1.1  mrg             {
    581      1.1  mrg               mpfr_set (tmp, dst_big, rnd);
    582      1.1  mrg               inexact = func(dst_small, arg1, arg2, rnd);
    583      1.1  mrg               if (mpfr_cmp (tmp, dst_small))
    584      1.1  mrg                 {
    585      1.1  mrg                   printf ("Results differ for prec=%u rnd_mode=%s and %s_z:\n"
    586      1.1  mrg                           "arg1=",
    587      1.1  mrg                           (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
    588      1.1  mrg                   mpfr_print_binary (arg1);
    589      1.1  mrg                   printf("\narg2=");
    590  1.1.1.2  mrg                   mpz_out_str (stdout, 10, arg2);
    591  1.1.1.2  mrg                   printf ("\ngot      ");
    592  1.1.1.2  mrg                   mpfr_dump (dst_small);
    593  1.1.1.2  mrg                   printf ("expected ");
    594  1.1.1.2  mrg                   mpfr_dump (tmp);
    595  1.1.1.2  mrg                   printf ("approx   ");
    596  1.1.1.2  mrg                   mpfr_dump (dst_big);
    597  1.1.1.2  mrg                   exit (1);
    598  1.1.1.2  mrg                 }
    599  1.1.1.2  mrg               compare2 = mpfr_cmp (tmp, dst_big);
    600  1.1.1.2  mrg               /* if rounding to nearest, cannot know the sign of t - f(x)
    601  1.1.1.2  mrg                  because of composed rounding: y = o(f(x)) and t = o(y) */
    602  1.1.1.2  mrg               if (compare * compare2 >= 0)
    603  1.1.1.2  mrg                 compare = compare + compare2;
    604  1.1.1.2  mrg               else
    605  1.1.1.2  mrg                 compare = inexact; /* cannot determine sign(t-f(x)) */
    606  1.1.1.2  mrg               if (((inexact == 0) && (compare != 0)) ||
    607  1.1.1.2  mrg                   ((inexact > 0) && (compare <= 0)) ||
    608  1.1.1.2  mrg                   ((inexact < 0) && (compare >= 0)))
    609  1.1.1.2  mrg                 {
    610  1.1.1.2  mrg                   printf ("Wrong inexact flag for rnd=%s and %s_z:\n"
    611  1.1.1.2  mrg                           "expected %d, got %d\n",
    612  1.1.1.2  mrg                           mpfr_print_rnd_mode (rnd), op, compare, inexact);
    613  1.1.1.2  mrg                   printf ("\narg1="); mpfr_print_binary (arg1);
    614  1.1.1.2  mrg                   printf ("\narg2="); mpz_out_str(stdout, 2, arg2);
    615  1.1.1.2  mrg                   printf ("\ndstl="); mpfr_print_binary (dst_big);
    616  1.1.1.2  mrg                   printf ("\ndsts="); mpfr_print_binary (dst_small);
    617  1.1.1.2  mrg                   printf ("\ntmp ="); mpfr_dump (tmp);
    618  1.1.1.2  mrg                   exit (1);
    619  1.1.1.2  mrg                 }
    620  1.1.1.2  mrg             }
    621  1.1.1.2  mrg         }
    622  1.1.1.2  mrg     }
    623  1.1.1.2  mrg 
    624  1.1.1.2  mrg   mpz_clear (arg2);
    625  1.1.1.2  mrg   mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    626  1.1.1.2  mrg }
    627  1.1.1.2  mrg 
    628  1.1.1.2  mrg static void
    629  1.1.1.2  mrg test_generic2z (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
    630  1.1.1.2  mrg                int (*func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
    631  1.1.1.2  mrg                const char *op)
    632  1.1.1.2  mrg {
    633  1.1.1.2  mrg   mpfr_prec_t prec;
    634  1.1.1.2  mrg   mpfr_t arg1, dst_big, dst_small, tmp;
    635  1.1.1.2  mrg   mpz_t  arg2;
    636  1.1.1.2  mrg   mpfr_rnd_t rnd;
    637  1.1.1.2  mrg   int inexact, compare, compare2;
    638  1.1.1.2  mrg   unsigned int n;
    639  1.1.1.2  mrg 
    640  1.1.1.2  mrg   mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    641  1.1.1.2  mrg   mpz_init (arg2);
    642  1.1.1.2  mrg 
    643  1.1.1.2  mrg   for (prec = p0; prec <= p1; prec++)
    644  1.1.1.2  mrg     {
    645  1.1.1.2  mrg       mpfr_set_prec (arg1, prec);
    646  1.1.1.2  mrg       mpfr_set_prec (tmp, prec);
    647  1.1.1.2  mrg       mpfr_set_prec (dst_small, prec);
    648  1.1.1.2  mrg 
    649  1.1.1.2  mrg       for (n=0; n<N; n++)
    650  1.1.1.2  mrg         {
    651  1.1.1.2  mrg           mpfr_urandomb (arg1, RANDS);
    652  1.1.1.2  mrg           mpz_urandomb (arg2, RANDS, 1024);
    653  1.1.1.2  mrg           rnd = RND_RAND ();
    654  1.1.1.2  mrg           mpfr_set_prec (dst_big, 2*prec);
    655  1.1.1.2  mrg           compare = func(dst_big, arg2, arg1, rnd);
    656  1.1.1.2  mrg           if (mpfr_can_round (dst_big, 2*prec, rnd, rnd, prec))
    657  1.1.1.2  mrg             {
    658  1.1.1.2  mrg               mpfr_set (tmp, dst_big, rnd);
    659  1.1.1.2  mrg               inexact = func(dst_small, arg2, arg1, rnd);
    660  1.1.1.2  mrg               if (mpfr_cmp (tmp, dst_small))
    661  1.1.1.2  mrg                 {
    662  1.1.1.2  mrg                   printf ("Results differ for prec=%u rnd_mode=%s and %s_z:\n"
    663  1.1.1.2  mrg                           "arg1=",
    664  1.1.1.2  mrg                           (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
    665  1.1.1.2  mrg                   mpfr_print_binary (arg1);
    666  1.1.1.2  mrg                   printf("\narg2=");
    667  1.1.1.2  mrg                   mpz_out_str (stdout, 10, arg2);
    668      1.1  mrg                   printf ("\ngot      ");
    669      1.1  mrg                   mpfr_dump (dst_small);
    670      1.1  mrg                   printf ("expected ");
    671      1.1  mrg                   mpfr_dump (tmp);
    672      1.1  mrg                   printf ("approx   ");
    673      1.1  mrg                   mpfr_dump (dst_big);
    674      1.1  mrg                   exit (1);
    675      1.1  mrg                 }
    676      1.1  mrg               compare2 = mpfr_cmp (tmp, dst_big);
    677      1.1  mrg               /* if rounding to nearest, cannot know the sign of t - f(x)
    678      1.1  mrg                  because of composed rounding: y = o(f(x)) and t = o(y) */
    679      1.1  mrg               if (compare * compare2 >= 0)
    680      1.1  mrg                 compare = compare + compare2;
    681      1.1  mrg               else
    682      1.1  mrg                 compare = inexact; /* cannot determine sign(t-f(x)) */
    683      1.1  mrg               if (((inexact == 0) && (compare != 0)) ||
    684      1.1  mrg                   ((inexact > 0) && (compare <= 0)) ||
    685      1.1  mrg                   ((inexact < 0) && (compare >= 0)))
    686      1.1  mrg                 {
    687      1.1  mrg                   printf ("Wrong inexact flag for rnd=%s and %s_z:\n"
    688      1.1  mrg                           "expected %d, got %d\n",
    689      1.1  mrg                           mpfr_print_rnd_mode (rnd), op, compare, inexact);
    690      1.1  mrg                   printf ("\narg1="); mpfr_print_binary (arg1);
    691      1.1  mrg                   printf ("\narg2="); mpz_out_str(stdout, 2, arg2);
    692      1.1  mrg                   printf ("\ndstl="); mpfr_print_binary (dst_big);
    693      1.1  mrg                   printf ("\ndsts="); mpfr_print_binary (dst_small);
    694      1.1  mrg                   printf ("\ntmp ="); mpfr_dump (tmp);
    695      1.1  mrg                   exit (1);
    696      1.1  mrg                 }
    697      1.1  mrg             }
    698      1.1  mrg         }
    699      1.1  mrg     }
    700      1.1  mrg 
    701      1.1  mrg   mpz_clear (arg2);
    702      1.1  mrg   mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    703      1.1  mrg }
    704      1.1  mrg 
    705      1.1  mrg static void
    706      1.1  mrg test_genericq (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
    707      1.1  mrg                int (*func)(mpfr_ptr, mpfr_srcptr, mpq_srcptr, mpfr_rnd_t),
    708      1.1  mrg                const char *op)
    709      1.1  mrg {
    710      1.1  mrg   mpfr_prec_t prec;
    711      1.1  mrg   mpfr_t arg1, dst_big, dst_small, tmp;
    712      1.1  mrg   mpq_t  arg2;
    713      1.1  mrg   mpfr_rnd_t rnd;
    714      1.1  mrg   int inexact, compare, compare2;
    715      1.1  mrg   unsigned int n;
    716      1.1  mrg 
    717      1.1  mrg   mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    718      1.1  mrg   mpq_init (arg2);
    719      1.1  mrg 
    720      1.1  mrg   for (prec = p0; prec <= p1; prec++)
    721      1.1  mrg     {
    722      1.1  mrg       mpfr_set_prec (arg1, prec);
    723      1.1  mrg       mpfr_set_prec (tmp, prec);
    724      1.1  mrg       mpfr_set_prec (dst_small, prec);
    725      1.1  mrg 
    726      1.1  mrg       for (n=0; n<N; n++)
    727      1.1  mrg         {
    728      1.1  mrg           mpfr_urandomb (arg1, RANDS);
    729      1.1  mrg           mpq_set_ui (arg2, randlimb (), randlimb() );
    730      1.1  mrg           mpq_canonicalize (arg2);
    731      1.1  mrg           rnd = RND_RAND ();
    732      1.1  mrg           mpfr_set_prec (dst_big, prec+10);
    733      1.1  mrg           compare = func(dst_big, arg1, arg2, rnd);
    734      1.1  mrg           if (mpfr_can_round (dst_big, prec+10, rnd, rnd, prec))
    735      1.1  mrg             {
    736      1.1  mrg               mpfr_set (tmp, dst_big, rnd);
    737      1.1  mrg               inexact = func(dst_small, arg1, arg2, rnd);
    738      1.1  mrg               if (mpfr_cmp (tmp, dst_small))
    739      1.1  mrg                 {
    740      1.1  mrg                   printf ("Results differ for prec=%u rnd_mode=%s and %s_q:\n"
    741      1.1  mrg                           "arg1=",
    742      1.1  mrg                           (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
    743      1.1  mrg                   mpfr_print_binary (arg1);
    744      1.1  mrg                   printf("\narg2=");
    745      1.1  mrg                   mpq_out_str(stdout, 2, arg2);
    746      1.1  mrg                   printf ("\ngot      ");
    747      1.1  mrg                   mpfr_print_binary (dst_small);
    748      1.1  mrg                   printf ("\nexpected ");
    749      1.1  mrg                   mpfr_print_binary (tmp);
    750      1.1  mrg                   printf ("\napprox  ");
    751      1.1  mrg                   mpfr_print_binary (dst_big);
    752      1.1  mrg                   putchar('\n');
    753      1.1  mrg                   exit (1);
    754      1.1  mrg                 }
    755      1.1  mrg               compare2 = mpfr_cmp (tmp, dst_big);
    756      1.1  mrg               /* if rounding to nearest, cannot know the sign of t - f(x)
    757      1.1  mrg                  because of composed rounding: y = o(f(x)) and t = o(y) */
    758      1.1  mrg               if (compare * compare2 >= 0)
    759      1.1  mrg                 compare = compare + compare2;
    760      1.1  mrg               else
    761      1.1  mrg                 compare = inexact; /* cannot determine sign(t-f(x)) */
    762      1.1  mrg               if (((inexact == 0) && (compare != 0)) ||
    763      1.1  mrg                   ((inexact > 0) && (compare <= 0)) ||
    764      1.1  mrg                   ((inexact < 0) && (compare >= 0)))
    765      1.1  mrg                 {
    766      1.1  mrg                   printf ("Wrong inexact flag for rnd=%s and %s_q:\n"
    767      1.1  mrg                           "expected %d, got %d",
    768      1.1  mrg                           mpfr_print_rnd_mode (rnd), op, compare, inexact);
    769      1.1  mrg                   printf ("\narg1="); mpfr_print_binary (arg1);
    770      1.1  mrg                   printf ("\narg2="); mpq_out_str(stdout, 2, arg2);
    771      1.1  mrg                   printf ("\ndstl="); mpfr_print_binary (dst_big);
    772      1.1  mrg                   printf ("\ndsts="); mpfr_print_binary (dst_small);
    773      1.1  mrg                   printf ("\ntmp ="); mpfr_print_binary (tmp);
    774      1.1  mrg                   putchar('\n');
    775      1.1  mrg                   exit (1);
    776      1.1  mrg                 }
    777      1.1  mrg             }
    778      1.1  mrg         }
    779      1.1  mrg     }
    780      1.1  mrg 
    781      1.1  mrg   mpq_clear (arg2);
    782      1.1  mrg   mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    783      1.1  mrg }
    784      1.1  mrg 
    785      1.1  mrg static void
    786      1.1  mrg test_specialq (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
    787      1.1  mrg                int (*mpfr_func)(mpfr_ptr, mpfr_srcptr, mpq_srcptr, mpfr_rnd_t),
    788      1.1  mrg                void (*mpq_func)(mpq_ptr, mpq_srcptr, mpq_srcptr),
    789      1.1  mrg                const char *op)
    790      1.1  mrg {
    791      1.1  mrg   mpfr_t fra, frb, frq;
    792      1.1  mrg   mpq_t  q1, q2, qr;
    793      1.1  mrg   unsigned int n;
    794      1.1  mrg   mpfr_prec_t prec;
    795      1.1  mrg 
    796      1.1  mrg   for (prec = p0 ; prec < p1 ; prec++)
    797      1.1  mrg     {
    798      1.1  mrg       mpfr_inits2 (prec, fra, frb, frq, (mpfr_ptr) 0);
    799      1.1  mrg       mpq_init (q1); mpq_init(q2); mpq_init (qr);
    800      1.1  mrg 
    801      1.1  mrg       for( n = 0 ; n < N ; n++)
    802      1.1  mrg         {
    803      1.1  mrg           mpq_set_ui(q1, randlimb(), randlimb() );
    804      1.1  mrg           mpq_set_ui(q2, randlimb(), randlimb() );
    805      1.1  mrg           mpq_canonicalize (q1);
    806      1.1  mrg           mpq_canonicalize (q2);
    807      1.1  mrg           mpq_func (qr, q1, q2);
    808      1.1  mrg           mpfr_set_q (fra, q1, MPFR_RNDD);
    809      1.1  mrg           mpfr_func (fra, fra, q2, MPFR_RNDD);
    810      1.1  mrg           mpfr_set_q (frb, q1, MPFR_RNDU);
    811      1.1  mrg           mpfr_func (frb, frb, q2, MPFR_RNDU);
    812      1.1  mrg           mpfr_set_q (frq, qr, MPFR_RNDN);
    813      1.1  mrg           /* We should have fra <= qr <= frb */
    814      1.1  mrg           if ( (mpfr_cmp(fra, frq) > 0) || (mpfr_cmp (frq, frb) > 0))
    815      1.1  mrg             {
    816  1.1.1.2  mrg               printf("Range error for prec=%lu and %s",
    817  1.1.1.2  mrg                      (unsigned long) prec, op);
    818      1.1  mrg               printf ("\nq1="); mpq_out_str(stdout, 2, q1);
    819      1.1  mrg               printf ("\nq2="); mpq_out_str(stdout, 2, q2);
    820      1.1  mrg               printf ("\nfr_dn="); mpfr_print_binary (fra);
    821      1.1  mrg               printf ("\nfr_q ="); mpfr_print_binary (frq);
    822      1.1  mrg               printf ("\nfr_up="); mpfr_print_binary (frb);
    823      1.1  mrg               putchar('\n');
    824      1.1  mrg               exit (1);
    825      1.1  mrg             }
    826      1.1  mrg         }
    827      1.1  mrg 
    828      1.1  mrg       mpq_clear (q1); mpq_clear (q2); mpq_clear (qr);
    829      1.1  mrg       mpfr_clears (fra, frb, frq, (mpfr_ptr) 0);
    830      1.1  mrg     }
    831      1.1  mrg }
    832      1.1  mrg 
    833  1.1.1.2  mrg static void
    834  1.1.1.2  mrg bug_mul_q_20100810 (void)
    835  1.1.1.2  mrg {
    836  1.1.1.2  mrg   mpfr_t x;
    837  1.1.1.2  mrg   mpfr_t y;
    838  1.1.1.2  mrg   mpq_t q;
    839  1.1.1.2  mrg   int inexact;
    840  1.1.1.2  mrg 
    841  1.1.1.2  mrg   mpfr_init (x);
    842  1.1.1.2  mrg   mpfr_init (y);
    843  1.1.1.2  mrg   mpq_init (q);
    844  1.1.1.2  mrg 
    845  1.1.1.2  mrg   /* mpfr_mul_q: the inexact value must be set in case of overflow */
    846  1.1.1.2  mrg   mpq_set_ui (q, 4096, 3);
    847  1.1.1.2  mrg   mpfr_set_inf (x, +1);
    848  1.1.1.2  mrg   mpfr_nextbelow (x);
    849  1.1.1.2  mrg   inexact = mpfr_mul_q (y, x, q, MPFR_RNDU);
    850  1.1.1.2  mrg 
    851  1.1.1.2  mrg   if (inexact <= 0)
    852  1.1.1.2  mrg     {
    853  1.1.1.2  mrg       printf ("Overflow error in mpfr_mul_q. ");
    854  1.1.1.2  mrg       printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
    855  1.1.1.2  mrg 
    856  1.1.1.2  mrg       exit (1);
    857  1.1.1.2  mrg     }
    858  1.1.1.2  mrg   if (!mpfr_inf_p (y))
    859  1.1.1.2  mrg     {
    860  1.1.1.2  mrg       printf ("Overflow error in mpfr_mul_q (y, x, q, MPFR_RNDD). ");
    861  1.1.1.2  mrg       printf ("\nx = ");
    862  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
    863  1.1.1.2  mrg       printf ("\nq = ");
    864  1.1.1.2  mrg       mpq_out_str (stdout, 10, q);
    865  1.1.1.2  mrg       printf ("\ny = ");
    866  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
    867  1.1.1.2  mrg       printf (" (should be +infinity)\n");
    868  1.1.1.2  mrg 
    869  1.1.1.2  mrg       exit (1);
    870  1.1.1.2  mrg     }
    871  1.1.1.2  mrg 
    872  1.1.1.2  mrg   mpq_clear (q);
    873  1.1.1.2  mrg   mpfr_clear (y);
    874  1.1.1.2  mrg   mpfr_clear (x);
    875  1.1.1.2  mrg }
    876  1.1.1.2  mrg 
    877  1.1.1.2  mrg static void
    878  1.1.1.2  mrg bug_div_q_20100810 (void)
    879  1.1.1.2  mrg {
    880  1.1.1.2  mrg   mpfr_t x;
    881  1.1.1.2  mrg   mpfr_t y;
    882  1.1.1.2  mrg   mpq_t q;
    883  1.1.1.2  mrg   int inexact;
    884  1.1.1.2  mrg 
    885  1.1.1.2  mrg   mpfr_init (x);
    886  1.1.1.2  mrg   mpfr_init (y);
    887  1.1.1.2  mrg   mpq_init (q);
    888  1.1.1.2  mrg 
    889  1.1.1.2  mrg   /* mpfr_div_q: the inexact value must be set in case of overflow */
    890  1.1.1.2  mrg   mpq_set_ui (q, 3, 4096);
    891  1.1.1.2  mrg   mpfr_set_inf (x, +1);
    892  1.1.1.2  mrg   mpfr_nextbelow (x);
    893  1.1.1.2  mrg   inexact = mpfr_div_q (y, x, q, MPFR_RNDU);
    894  1.1.1.2  mrg 
    895  1.1.1.2  mrg   if (inexact <= 0)
    896  1.1.1.2  mrg     {
    897  1.1.1.2  mrg       printf ("Overflow error in mpfr_div_q. ");
    898  1.1.1.2  mrg       printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
    899  1.1.1.2  mrg 
    900  1.1.1.2  mrg       exit (1);
    901  1.1.1.2  mrg     }
    902  1.1.1.2  mrg   if (!mpfr_inf_p (y))
    903  1.1.1.2  mrg     {
    904  1.1.1.2  mrg       printf ("Overflow error in mpfr_div_q (y, x, q, MPFR_RNDD). ");
    905  1.1.1.2  mrg       printf ("\nx = ");
    906  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
    907  1.1.1.2  mrg       printf ("\nq = ");
    908  1.1.1.2  mrg       mpq_out_str (stdout, 10, q);
    909  1.1.1.2  mrg       printf ("\ny = ");
    910  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
    911  1.1.1.2  mrg       printf (" (should be +infinity)\n");
    912  1.1.1.2  mrg 
    913  1.1.1.2  mrg       exit (1);
    914  1.1.1.2  mrg     }
    915  1.1.1.2  mrg 
    916  1.1.1.2  mrg   mpq_clear (q);
    917  1.1.1.2  mrg   mpfr_clear (y);
    918  1.1.1.2  mrg   mpfr_clear (x);
    919  1.1.1.2  mrg }
    920  1.1.1.2  mrg 
    921  1.1.1.2  mrg static void
    922  1.1.1.2  mrg bug_mul_div_q_20100818 (void)
    923  1.1.1.2  mrg {
    924  1.1.1.2  mrg   mpq_t qa, qb;
    925  1.1.1.2  mrg   mpfr_t x1, x2, y1, y2, y3;
    926  1.1.1.2  mrg   mpfr_exp_t emin, emax, e;
    927  1.1.1.2  mrg   int inex;
    928  1.1.1.2  mrg   int rnd;
    929  1.1.1.2  mrg 
    930  1.1.1.2  mrg   emin = mpfr_get_emin ();
    931  1.1.1.2  mrg   emax = mpfr_get_emax ();
    932  1.1.1.2  mrg   set_emin (MPFR_EMIN_MIN);
    933  1.1.1.2  mrg   set_emax (MPFR_EMAX_MAX);
    934  1.1.1.2  mrg 
    935  1.1.1.2  mrg   mpq_init (qa);
    936  1.1.1.2  mrg   mpq_init (qb);
    937  1.1.1.2  mrg   mpfr_inits2 (32, x1, x2, y1, y2, y3, (mpfr_ptr) 0);
    938  1.1.1.2  mrg 
    939  1.1.1.2  mrg   mpq_set_ui (qa, 3, 17);
    940  1.1.1.2  mrg   mpq_set_ui (qb, 17, 3);
    941  1.1.1.2  mrg   inex = mpfr_set_ui (x1, 7, MPFR_RNDN);
    942  1.1.1.2  mrg   MPFR_ASSERTN (inex == 0);
    943  1.1.1.2  mrg 
    944  1.1.1.2  mrg   e = MPFR_EMAX_MAX - 3;
    945  1.1.1.2  mrg   inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN);  /* x2 = x1 * 2^e */
    946  1.1.1.2  mrg   MPFR_ASSERTN (inex == 0);
    947  1.1.1.2  mrg 
    948  1.1.1.2  mrg   RND_LOOP(rnd)
    949  1.1.1.2  mrg     {
    950  1.1.1.2  mrg       mpfr_mul_q (y1, x1, qa, (mpfr_rnd_t) rnd);
    951  1.1.1.2  mrg       mpfr_div_q (y3, x1, qb, (mpfr_rnd_t) rnd);
    952  1.1.1.2  mrg       MPFR_ASSERTN (mpfr_equal_p (y1, y3));
    953  1.1.1.2  mrg       inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
    954  1.1.1.2  mrg       MPFR_ASSERTN (inex == 0);
    955  1.1.1.2  mrg       inex = mpfr_mul (y3, y3, y1, MPFR_RNDN);  /* y3 = y1 * 2^e */
    956  1.1.1.2  mrg       MPFR_ASSERTN (inex == 0);
    957  1.1.1.2  mrg       mpfr_mul_q (y2, x2, qa, (mpfr_rnd_t) rnd);
    958  1.1.1.2  mrg       if (! mpfr_equal_p (y2, y3))
    959  1.1.1.2  mrg         {
    960  1.1.1.2  mrg           printf ("Error 1 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    961  1.1.1.2  mrg           printf ("Expected "); mpfr_dump (y3);
    962  1.1.1.2  mrg           printf ("Got      "); mpfr_dump (y2);
    963  1.1.1.2  mrg           exit (1);
    964  1.1.1.2  mrg         }
    965  1.1.1.2  mrg       mpfr_div_q (y2, x2, qb, (mpfr_rnd_t) rnd);
    966  1.1.1.2  mrg       if (! mpfr_equal_p (y2, y3))
    967  1.1.1.2  mrg         {
    968  1.1.1.2  mrg           printf ("Error 2 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    969  1.1.1.2  mrg           printf ("Expected "); mpfr_dump (y3);
    970  1.1.1.2  mrg           printf ("Got      "); mpfr_dump (y2);
    971  1.1.1.2  mrg           exit (1);
    972  1.1.1.2  mrg         }
    973  1.1.1.2  mrg     }
    974  1.1.1.2  mrg 
    975  1.1.1.2  mrg   e = MPFR_EMIN_MIN;
    976  1.1.1.2  mrg   inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN);  /* x2 = x1 * 2^e */
    977  1.1.1.2  mrg   MPFR_ASSERTN (inex == 0);
    978  1.1.1.2  mrg 
    979  1.1.1.2  mrg   RND_LOOP(rnd)
    980  1.1.1.2  mrg     {
    981  1.1.1.2  mrg       mpfr_div_q (y1, x1, qa, (mpfr_rnd_t) rnd);
    982  1.1.1.2  mrg       mpfr_mul_q (y3, x1, qb, (mpfr_rnd_t) rnd);
    983  1.1.1.2  mrg       MPFR_ASSERTN (mpfr_equal_p (y1, y3));
    984  1.1.1.2  mrg       inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
    985  1.1.1.2  mrg       MPFR_ASSERTN (inex == 0);
    986  1.1.1.2  mrg       inex = mpfr_mul (y3, y3, y1, MPFR_RNDN);  /* y3 = y1 * 2^e */
    987  1.1.1.2  mrg       MPFR_ASSERTN (inex == 0);
    988  1.1.1.2  mrg       mpfr_div_q (y2, x2, qa, (mpfr_rnd_t) rnd);
    989  1.1.1.2  mrg       if (! mpfr_equal_p (y2, y3))
    990  1.1.1.2  mrg         {
    991  1.1.1.2  mrg           printf ("Error 3 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    992  1.1.1.2  mrg           printf ("Expected "); mpfr_dump (y3);
    993  1.1.1.2  mrg           printf ("Got      "); mpfr_dump (y2);
    994  1.1.1.2  mrg           exit (1);
    995  1.1.1.2  mrg         }
    996  1.1.1.2  mrg       mpfr_mul_q (y2, x2, qb, (mpfr_rnd_t) rnd);
    997  1.1.1.2  mrg       if (! mpfr_equal_p (y2, y3))
    998  1.1.1.2  mrg         {
    999  1.1.1.2  mrg           printf ("Error 4 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
   1000  1.1.1.2  mrg           printf ("Expected "); mpfr_dump (y3);
   1001  1.1.1.2  mrg           printf ("Got      "); mpfr_dump (y2);
   1002  1.1.1.2  mrg           exit (1);
   1003  1.1.1.2  mrg         }
   1004  1.1.1.2  mrg     }
   1005  1.1.1.2  mrg 
   1006  1.1.1.2  mrg   mpq_clear (qa);
   1007  1.1.1.2  mrg   mpq_clear (qb);
   1008  1.1.1.2  mrg   mpfr_clears (x1, x2, y1, y2, y3, (mpfr_ptr) 0);
   1009  1.1.1.2  mrg 
   1010  1.1.1.2  mrg   set_emin (emin);
   1011  1.1.1.2  mrg   set_emax (emax);
   1012  1.1.1.2  mrg }
   1013  1.1.1.2  mrg 
   1014  1.1.1.2  mrg static void
   1015  1.1.1.2  mrg reduced_expo_range (void)
   1016  1.1.1.2  mrg {
   1017  1.1.1.2  mrg   mpfr_t x;
   1018  1.1.1.2  mrg   mpz_t z;
   1019  1.1.1.2  mrg   mpq_t q;
   1020  1.1.1.2  mrg   mpfr_exp_t emin;
   1021  1.1.1.2  mrg   int inex;
   1022  1.1.1.2  mrg 
   1023  1.1.1.2  mrg   emin = mpfr_get_emin ();
   1024  1.1.1.2  mrg   set_emin (4);
   1025  1.1.1.2  mrg 
   1026  1.1.1.2  mrg   mpfr_init2 (x, 32);
   1027  1.1.1.2  mrg 
   1028  1.1.1.2  mrg   mpz_init (z);
   1029  1.1.1.2  mrg   mpfr_clear_flags ();
   1030  1.1.1.2  mrg   inex = mpfr_set_ui (x, 17, MPFR_RNDN);
   1031  1.1.1.2  mrg   MPFR_ASSERTN (inex == 0);
   1032  1.1.1.2  mrg   mpz_set_ui (z, 3);
   1033  1.1.1.2  mrg   inex = mpfr_mul_z (x, x, z, MPFR_RNDN);
   1034  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 51) != 0)
   1035  1.1.1.2  mrg     {
   1036  1.1.1.2  mrg       printf ("Error 1 in reduce_expo_range: expected 51 with inex = 0,"
   1037  1.1.1.2  mrg               " got\n");
   1038  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1039  1.1.1.2  mrg       printf ("with inex = %d\n", inex);
   1040  1.1.1.2  mrg       exit (1);
   1041  1.1.1.2  mrg     }
   1042  1.1.1.2  mrg   inex = mpfr_div_z (x, x, z, MPFR_RNDN);
   1043  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1044  1.1.1.2  mrg     {
   1045  1.1.1.2  mrg       printf ("Error 2 in reduce_expo_range: expected 17 with inex = 0,"
   1046  1.1.1.2  mrg               " got\n");
   1047  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1048  1.1.1.2  mrg       printf ("with inex = %d\n", inex);
   1049  1.1.1.2  mrg       exit (1);
   1050  1.1.1.2  mrg     }
   1051  1.1.1.2  mrg   inex = mpfr_add_z (x, x, z, MPFR_RNDN);
   1052  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 20) != 0)
   1053  1.1.1.2  mrg     {
   1054  1.1.1.2  mrg       printf ("Error 3 in reduce_expo_range: expected 20 with inex = 0,"
   1055  1.1.1.2  mrg               " got\n");
   1056  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1057  1.1.1.2  mrg       printf ("with inex = %d\n", inex);
   1058  1.1.1.2  mrg       exit (1);
   1059  1.1.1.2  mrg     }
   1060  1.1.1.2  mrg   inex = mpfr_sub_z (x, x, z, MPFR_RNDN);
   1061  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1062  1.1.1.2  mrg     {
   1063  1.1.1.2  mrg       printf ("Error 4 in reduce_expo_range: expected 17 with inex = 0,"
   1064  1.1.1.2  mrg               " got\n");
   1065  1.1.1.2  mrg       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1066  1.1.1.2  mrg       printf ("with inex = %d\n", inex);
   1067  1.1.1.2  mrg       exit (1);
   1068  1.1.1.2  mrg     }
   1069  1.1.1.2  mrg   MPFR_ASSERTN (__gmpfr_flags == 0);
   1070  1.1.1.2  mrg   if (mpfr_cmp_z (x, z) <= 0)
   1071  1.1.1.2  mrg     {
   1072  1.1.1.2  mrg       printf ("Error 5 in reduce_expo_range: expected a positive value.\n");
   1073  1.1.1.2  mrg       exit (1);
   1074  1.1.1.2  mrg     }
   1075  1.1.1.2  mrg   mpz_clear (z);
   1076  1.1.1.2  mrg 
   1077  1.1.1.2  mrg   mpq_init (q);
   1078  1.1.1.2  mrg   mpq_set_ui (q, 1, 1);
   1079  1.1.1.2  mrg   mpfr_set_ui (x, 16, MPFR_RNDN);
   1080  1.1.1.2  mrg   inex = mpfr_add_q (x, x, q, MPFR_RNDN);
   1081  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1082  1.1.1.2  mrg     {
   1083  1.1.1.2  mrg       printf ("Error in reduce_expo_range for 16 + 1/1,"
   1084  1.1.1.2  mrg               " got inex = %d and\nx = ", inex);
   1085  1.1.1.2  mrg       mpfr_dump (x);
   1086  1.1.1.2  mrg       exit (1);
   1087  1.1.1.2  mrg     }
   1088  1.1.1.2  mrg   inex = mpfr_sub_q (x, x, q, MPFR_RNDN);
   1089  1.1.1.2  mrg   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 16) != 0)
   1090  1.1.1.2  mrg     {
   1091  1.1.1.2  mrg       printf ("Error in reduce_expo_range for 17 - 1/1,"
   1092  1.1.1.2  mrg               " got inex = %d and\nx = ", inex);
   1093  1.1.1.2  mrg       mpfr_dump (x);
   1094  1.1.1.2  mrg       exit (1);
   1095  1.1.1.2  mrg     }
   1096  1.1.1.2  mrg   mpq_clear (q);
   1097  1.1.1.2  mrg 
   1098  1.1.1.2  mrg   mpfr_clear (x);
   1099  1.1.1.2  mrg 
   1100  1.1.1.2  mrg   set_emin (emin);
   1101  1.1.1.2  mrg }
   1102  1.1.1.2  mrg 
   1103  1.1.1.2  mrg static void
   1104  1.1.1.2  mrg addsubq_overflow_aux (mpfr_exp_t e)
   1105  1.1.1.2  mrg {
   1106  1.1.1.2  mrg   mpfr_t x, y;
   1107  1.1.1.2  mrg   mpq_t q;
   1108  1.1.1.2  mrg   mpfr_exp_t emax;
   1109  1.1.1.2  mrg   int inex;
   1110  1.1.1.2  mrg   int rnd;
   1111  1.1.1.2  mrg   int sign, sub;
   1112  1.1.1.2  mrg 
   1113  1.1.1.2  mrg   MPFR_ASSERTN (e <= LONG_MAX);
   1114  1.1.1.2  mrg   emax = mpfr_get_emax ();
   1115  1.1.1.2  mrg   set_emax (e);
   1116  1.1.1.2  mrg   mpfr_inits2 (16, x, y, (mpfr_ptr) 0);
   1117  1.1.1.2  mrg   mpq_init (q);
   1118  1.1.1.2  mrg 
   1119  1.1.1.2  mrg   mpfr_set_inf (x, 1);
   1120  1.1.1.2  mrg   mpfr_nextbelow (x);
   1121  1.1.1.2  mrg   mpq_set_ui (q, 1, 1);
   1122  1.1.1.2  mrg 
   1123  1.1.1.2  mrg   for (sign = 0; sign <= 1; sign++)
   1124  1.1.1.2  mrg     {
   1125  1.1.1.2  mrg       for (sub = 0; sub <= 1; sub++)
   1126  1.1.1.2  mrg         {
   1127  1.1.1.2  mrg           RND_LOOP(rnd)
   1128  1.1.1.2  mrg             {
   1129  1.1.1.2  mrg               unsigned int flags, ex_flags;
   1130  1.1.1.2  mrg               int inf;
   1131  1.1.1.2  mrg 
   1132  1.1.1.2  mrg               inf = rnd == MPFR_RNDA ||
   1133  1.1.1.2  mrg                     rnd == (sign ? MPFR_RNDD : MPFR_RNDU);
   1134  1.1.1.2  mrg               ex_flags = MPFR_FLAGS_INEXACT | (inf ? MPFR_FLAGS_OVERFLOW : 0);
   1135  1.1.1.2  mrg               mpfr_clear_flags ();
   1136  1.1.1.2  mrg               inex = sub ?
   1137  1.1.1.2  mrg                 mpfr_sub_q (y, x, q, (mpfr_rnd_t) rnd) :
   1138  1.1.1.2  mrg                 mpfr_add_q (y, x, q, (mpfr_rnd_t) rnd);
   1139  1.1.1.2  mrg               flags = __gmpfr_flags;
   1140  1.1.1.2  mrg               if (inex == 0 || flags != ex_flags ||
   1141  1.1.1.2  mrg                   (inf ? ! mpfr_inf_p (y) : ! mpfr_equal_p (x, y)))
   1142  1.1.1.2  mrg                 {
   1143  1.1.1.2  mrg                   printf ("Error in addsubq_overflow_aux(%ld),"
   1144  1.1.1.2  mrg                           " sign = %d, %s\n", (long) e, sign,
   1145  1.1.1.2  mrg                           mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
   1146  1.1.1.2  mrg                   printf ("Got inex = %d, y = ", inex);
   1147  1.1.1.2  mrg                   mpfr_dump (y);
   1148  1.1.1.2  mrg                   printf ("Expected flags:");
   1149  1.1.1.2  mrg                   flags_out (ex_flags);
   1150  1.1.1.2  mrg                   printf ("Got flags:     ");
   1151  1.1.1.2  mrg                   flags_out (flags);
   1152  1.1.1.2  mrg                   exit (1);
   1153  1.1.1.2  mrg                 }
   1154  1.1.1.2  mrg             }
   1155  1.1.1.2  mrg           mpq_neg (q, q);
   1156  1.1.1.2  mrg         }
   1157  1.1.1.2  mrg       mpfr_neg (x, x, MPFR_RNDN);
   1158  1.1.1.2  mrg       mpq_neg (q, q);
   1159  1.1.1.2  mrg     }
   1160  1.1.1.2  mrg 
   1161  1.1.1.2  mrg   mpq_clear (q);
   1162  1.1.1.2  mrg   mpfr_clears (x, y, (mpfr_ptr) 0);
   1163  1.1.1.2  mrg   set_emax (emax);
   1164  1.1.1.2  mrg }
   1165  1.1.1.2  mrg 
   1166  1.1.1.2  mrg static void
   1167  1.1.1.2  mrg addsubq_overflow (void)
   1168  1.1.1.2  mrg {
   1169  1.1.1.2  mrg   addsubq_overflow_aux (4913);
   1170  1.1.1.2  mrg   addsubq_overflow_aux (MPFR_EMAX_MAX);
   1171  1.1.1.2  mrg }
   1172  1.1.1.2  mrg 
   1173  1.1.1.2  mrg static void
   1174  1.1.1.2  mrg coverage_mpfr_mul_q_20110218 (void)
   1175  1.1.1.2  mrg {
   1176  1.1.1.2  mrg   mpfr_t cmp, res, op1;
   1177  1.1.1.2  mrg   mpq_t op2;
   1178  1.1.1.2  mrg   int status;
   1179  1.1.1.2  mrg 
   1180  1.1.1.2  mrg   mpfr_init2 (cmp, MPFR_PREC_MIN);
   1181  1.1.1.2  mrg   mpfr_init2 (res, MPFR_PREC_MIN);
   1182  1.1.1.2  mrg   mpfr_init_set_si (op1, 1, MPFR_RNDN);
   1183  1.1.1.2  mrg 
   1184  1.1.1.2  mrg   mpq_init (op2);
   1185  1.1.1.2  mrg   mpq_set_si (op2, 0, 0);
   1186  1.1.1.2  mrg   mpz_set_si (mpq_denref (op2), 0);
   1187  1.1.1.2  mrg 
   1188  1.1.1.2  mrg   status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
   1189  1.1.1.2  mrg 
   1190  1.1.1.2  mrg   if ((status != 0) || (mpfr_cmp (cmp, res) != 0))
   1191  1.1.1.2  mrg     {
   1192  1.1.1.2  mrg       printf ("Results differ %d.\nres=", status);
   1193  1.1.1.2  mrg       mpfr_print_binary (res);
   1194  1.1.1.2  mrg       printf ("\ncmp=");
   1195  1.1.1.2  mrg       mpfr_print_binary (cmp);
   1196  1.1.1.2  mrg       putchar ('\n');
   1197  1.1.1.2  mrg       exit (1);
   1198  1.1.1.2  mrg     }
   1199  1.1.1.2  mrg 
   1200  1.1.1.2  mrg   mpfr_set_si (op1, 1, MPFR_RNDN);
   1201  1.1.1.2  mrg   mpq_set_si (op2, -1, 0);
   1202  1.1.1.2  mrg 
   1203  1.1.1.2  mrg   status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
   1204  1.1.1.2  mrg 
   1205  1.1.1.2  mrg   mpfr_set_inf (cmp, -1);
   1206  1.1.1.2  mrg   if ((status != 0) || (mpfr_cmp(res, cmp) != 0))
   1207  1.1.1.2  mrg     {
   1208  1.1.1.2  mrg       printf ("mpfr_mul_q 1 * (-1/0) returned a wrong value :\n waiting for ");
   1209  1.1.1.2  mrg       mpfr_print_binary (cmp);
   1210  1.1.1.2  mrg       printf (" got ");
   1211  1.1.1.2  mrg       mpfr_print_binary (res);
   1212  1.1.1.2  mrg       printf ("\n trinary value is %d\n", status);
   1213  1.1.1.2  mrg       exit (1);
   1214  1.1.1.2  mrg     }
   1215  1.1.1.2  mrg 
   1216  1.1.1.2  mrg   mpq_clear (op2);
   1217  1.1.1.2  mrg   mpfr_clear (op1);
   1218  1.1.1.2  mrg   mpfr_clear (res);
   1219  1.1.1.2  mrg   mpfr_clear (cmp);
   1220  1.1.1.2  mrg }
   1221  1.1.1.2  mrg 
   1222      1.1  mrg int
   1223      1.1  mrg main (int argc, char *argv[])
   1224      1.1  mrg {
   1225      1.1  mrg   tests_start_mpfr ();
   1226      1.1  mrg 
   1227      1.1  mrg   special ();
   1228      1.1  mrg 
   1229      1.1  mrg   test_specialz (mpfr_add_z, mpz_add, "add");
   1230      1.1  mrg   test_specialz (mpfr_sub_z, mpz_sub, "sub");
   1231      1.1  mrg   test_specialz (mpfr_mul_z, mpz_mul, "mul");
   1232      1.1  mrg   test_genericz (2, 100, 100, mpfr_add_z, "add");
   1233      1.1  mrg   test_genericz (2, 100, 100, mpfr_sub_z, "sub");
   1234      1.1  mrg   test_genericz (2, 100, 100, mpfr_mul_z, "mul");
   1235      1.1  mrg   test_genericz (2, 100, 100, mpfr_div_z, "div");
   1236  1.1.1.2  mrg   test_special2z (mpfr_z_sub, mpz_sub, "sub");
   1237  1.1.1.2  mrg   test_generic2z (2, 100, 100, mpfr_z_sub, "sub");
   1238      1.1  mrg 
   1239      1.1  mrg   test_genericq (2, 100, 100, mpfr_add_q, "add");
   1240      1.1  mrg   test_genericq (2, 100, 100, mpfr_sub_q, "sub");
   1241      1.1  mrg   test_genericq (2, 100, 100, mpfr_mul_q, "mul");
   1242      1.1  mrg   test_genericq (2, 100, 100, mpfr_div_q, "div");
   1243      1.1  mrg   test_specialq (2, 100, 100, mpfr_mul_q, mpq_mul, "mul");
   1244      1.1  mrg   test_specialq (2, 100, 100, mpfr_div_q, mpq_div, "div");
   1245      1.1  mrg   test_specialq (2, 100, 100, mpfr_add_q, mpq_add, "add");
   1246      1.1  mrg   test_specialq (2, 100, 100, mpfr_sub_q, mpq_sub, "sub");
   1247      1.1  mrg 
   1248      1.1  mrg   test_cmp_z (2, 100, 100);
   1249      1.1  mrg   test_cmp_q (2, 100, 100);
   1250      1.1  mrg   test_cmp_f (2, 100, 100);
   1251      1.1  mrg 
   1252      1.1  mrg   check_for_zero ();
   1253      1.1  mrg 
   1254  1.1.1.2  mrg   bug_mul_q_20100810 ();
   1255  1.1.1.2  mrg   bug_div_q_20100810 ();
   1256  1.1.1.2  mrg   bug_mul_div_q_20100818 ();
   1257  1.1.1.2  mrg   reduced_expo_range ();
   1258  1.1.1.2  mrg   addsubq_overflow ();
   1259  1.1.1.2  mrg 
   1260  1.1.1.2  mrg   coverage_mpfr_mul_q_20110218 ();
   1261  1.1.1.2  mrg 
   1262      1.1  mrg   tests_end_mpfr ();
   1263      1.1  mrg   return 0;
   1264      1.1  mrg }
   1265      1.1  mrg 
   1266