Home | History | Annotate | Line # | Download | only in tests
tgmpop.c revision 1.1.1.1.8.1
      1  1.1.1.1.8.1  tls /* Test file for mpfr_add_[q,z], mpfr_sub_[q,z], mpfr_div_[q,z],
      2  1.1.1.1.8.1  tls    mpfr_mul_[q,z], mpfr_cmp_[f,q,z]
      3          1.1  mrg 
      4  1.1.1.1.8.1  tls Copyright 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012, 2013 Free Software Foundation, Inc.
      5  1.1.1.1.8.1  tls Contributed by the AriC and Caramel 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.1.8.1  tls #define CHECK_FOR(str, cond)                                            \
     29  1.1.1.1.8.1  tls   if ((cond) == 0) {                                                    \
     30  1.1.1.1.8.1  tls     printf ("Special case error %s. Ternary value = %d, flags = %u\n",  \
     31  1.1.1.1.8.1  tls             str, res, __gmpfr_flags);                                   \
     32  1.1.1.1.8.1  tls     printf ("Got "); mpfr_dump (y);                                     \
     33  1.1.1.1.8.1  tls     printf ("X = "); mpfr_dump (x);                                     \
     34  1.1.1.1.8.1  tls     printf ("Q = "); mpz_dump (mpq_numref(q));                          \
     35  1.1.1.1.8.1  tls     printf ("   /"); mpz_dump (mpq_denref(q));                          \
     36  1.1.1.1.8.1  tls     exit (1);                                                           \
     37  1.1.1.1.8.1  tls   }
     38  1.1.1.1.8.1  tls 
     39  1.1.1.1.8.1  tls #define CHECK_FORZ(str, cond)                                           \
     40  1.1.1.1.8.1  tls   if ((cond) == 0) {                                                    \
     41  1.1.1.1.8.1  tls     printf ("Special case error %s. Ternary value = %d, flags = %u\n",  \
     42  1.1.1.1.8.1  tls             str, res, __gmpfr_flags);                                   \
     43  1.1.1.1.8.1  tls     printf ("Got "); mpfr_dump (y);                                     \
     44  1.1.1.1.8.1  tls     printf ("X = "); mpfr_dump (x);                                     \
     45  1.1.1.1.8.1  tls     printf ("Z = "); mpz_dump (z);                                      \
     46  1.1.1.1.8.1  tls     exit (1);                                                           \
     47  1.1.1.1.8.1  tls   }
     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.1.8.1  tls   mpq_t q;
     54  1.1.1.1.8.1  tls   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.1.8.1  tls   mpq_init (q);
     60  1.1.1.1.8.1  tls   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.1.8.1  tls   mpz_set_str (mpq_numref (q), "-187207494", 10);
     66  1.1.1.1.8.1  tls   mpz_set_str (mpq_denref (q), "5721", 10);
     67          1.1  mrg   mpfr_set_str_binary (x, "11111111101001011011100101100011011110010011100010000100001E-44");
     68  1.1.1.1.8.1  tls   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.1.8.1  tls   mpz_set_str (mpq_numref (q), "187207494", 10);
     74  1.1.1.1.8.1  tls   mpz_set_str (mpq_denref (q), "5721", 10);
     75          1.1  mrg   mpfr_set_prec (y, 29);
     76  1.1.1.1.8.1  tls   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.1.8.1  tls   mpz_set_str (mpq_numref (q), "395877315", 10);
     84  1.1.1.1.8.1  tls   mpz_set_str (mpq_denref (q), "3508975966", 10);
     85          1.1  mrg   mpfr_set_prec (y, 118);
     86  1.1.1.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   mpfr_add_q (y, x, q, MPFR_RNDU);
     94          1.1  mrg   CHECK_FOR ("nan", mpfr_nan_p (y));
     95  1.1.1.1.8.1  tls   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.1.8.1  tls   mpz_set_str (mpq_numref (q), "3", 10);
    103  1.1.1.1.8.1  tls   mpz_set_str (mpq_denref (q), "2", 10);
    104  1.1.1.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   mpq_set_ui (q, 1, 0);
    111          1.1  mrg   mpfr_set_str1 (x, "0.5");
    112  1.1.1.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   mpq_set_si (q, -1, 0);
    117  1.1.1.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   mpq_set_ui (q, 1, 0);
    124  1.1.1.1.8.1  tls   mpfr_set_inf (x, 1);
    125  1.1.1.1.8.1  tls   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    126  1.1.1.1.8.1  tls   CHECK_FOR ("+Inf + +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    127  1.1.1.1.8.1  tls   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    128  1.1.1.1.8.1  tls   CHECK_FOR ("+Inf - +Inf", MPFR_IS_NAN (y) && res == 0);
    129  1.1.1.1.8.1  tls   mpfr_set_inf (x, -1);
    130  1.1.1.1.8.1  tls   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    131  1.1.1.1.8.1  tls   CHECK_FOR ("-Inf + +Inf", MPFR_IS_NAN (y) && res == 0);
    132  1.1.1.1.8.1  tls   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    133  1.1.1.1.8.1  tls   CHECK_FOR ("-Inf - +Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    134  1.1.1.1.8.1  tls   mpq_set_si (q, -1, 0);
    135  1.1.1.1.8.1  tls   mpfr_set_inf (x, 1);
    136  1.1.1.1.8.1  tls   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    137  1.1.1.1.8.1  tls   CHECK_FOR ("+Inf + -Inf", MPFR_IS_NAN (y) && res == 0);
    138  1.1.1.1.8.1  tls   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    139  1.1.1.1.8.1  tls   CHECK_FOR ("+Inf - -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0);
    140  1.1.1.1.8.1  tls   mpfr_set_inf (x, -1);
    141  1.1.1.1.8.1  tls   res = mpfr_add_q (y, x, q, MPFR_RNDN);
    142  1.1.1.1.8.1  tls   CHECK_FOR ("-Inf + -Inf", mpfr_inf_p (y) && MPFR_SIGN (y) < 0 && res == 0);
    143  1.1.1.1.8.1  tls   res = mpfr_sub_q (y, x, q, MPFR_RNDN);
    144  1.1.1.1.8.1  tls   CHECK_FOR ("-Inf - -Inf", MPFR_IS_NAN (y) && res == 0);
    145          1.1  mrg 
    146          1.1  mrg   /* 0 */
    147  1.1.1.1.8.1  tls   mpq_set_ui (q, 0, 1);
    148          1.1  mrg   mpfr_set_ui (x, 42, MPFR_RNDN);
    149  1.1.1.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   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.1.8.1  tls   mpfr_clear_flags ();
    156  1.1.1.1.8.1  tls   res = mpfr_div_q (y, x, q, MPFR_RNDN);
    157  1.1.1.1.8.1  tls   CHECK_FOR ("42/(0/1)", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
    158  1.1.1.1.8.1  tls              && mpfr_divby0_p ());
    159  1.1.1.1.8.1  tls   mpz_set_ui (z, 0);
    160  1.1.1.1.8.1  tls   mpfr_clear_flags ();
    161  1.1.1.1.8.1  tls   res = mpfr_div_z (y, x, z, MPFR_RNDN);
    162  1.1.1.1.8.1  tls   CHECK_FORZ ("42/0", mpfr_inf_p (y) && MPFR_SIGN (y) > 0 && res == 0
    163  1.1.1.1.8.1  tls               && mpfr_divby0_p ());
    164          1.1  mrg 
    165  1.1.1.1.8.1  tls   mpz_clear (z);
    166  1.1.1.1.8.1  tls   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.1.8.1  tls 
    256  1.1.1.1.8.1  tls   /* check the erange flag when x is NaN */
    257  1.1.1.1.8.1  tls   mpfr_set_nan (x);
    258  1.1.1.1.8.1  tls   mpz_set_ui (y, 17);
    259  1.1.1.1.8.1  tls   mpfr_clear_erangeflag ();
    260  1.1.1.1.8.1  tls   res1 = mpfr_cmp_z (x, y);
    261  1.1.1.1.8.1  tls   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    262  1.1.1.1.8.1  tls     {
    263  1.1.1.1.8.1  tls       printf ("Error for mpfr_cmp_z (NaN, 17)\n");
    264  1.1.1.1.8.1  tls       printf ("Return value: expected 0, got %d\n", res1);
    265  1.1.1.1.8.1  tls       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    266  1.1.1.1.8.1  tls       exit (1);
    267  1.1.1.1.8.1  tls     }
    268  1.1.1.1.8.1  tls 
    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.1.8.1  tls 
    308  1.1.1.1.8.1  tls   /* check the erange flag when x is NaN */
    309  1.1.1.1.8.1  tls   mpfr_set_nan (x);
    310  1.1.1.1.8.1  tls   mpq_set_ui (y, 17, 1);
    311  1.1.1.1.8.1  tls   mpfr_clear_erangeflag ();
    312  1.1.1.1.8.1  tls   res1 = mpfr_cmp_q (x, y);
    313  1.1.1.1.8.1  tls   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    314  1.1.1.1.8.1  tls     {
    315  1.1.1.1.8.1  tls       printf ("Error for mpfr_cmp_q (NaN, 17)\n");
    316  1.1.1.1.8.1  tls       printf ("Return value: expected 0, got %d\n", res1);
    317  1.1.1.1.8.1  tls       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    318  1.1.1.1.8.1  tls       exit (1);
    319  1.1.1.1.8.1  tls     }
    320  1.1.1.1.8.1  tls 
    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.1.8.1  tls 
    360  1.1.1.1.8.1  tls   /* check the erange flag when x is NaN */
    361  1.1.1.1.8.1  tls   mpfr_set_nan (x);
    362  1.1.1.1.8.1  tls   mpf_set_ui (y, 17);
    363  1.1.1.1.8.1  tls   mpfr_clear_erangeflag ();
    364  1.1.1.1.8.1  tls   res1 = mpfr_cmp_f (x, y);
    365  1.1.1.1.8.1  tls   if (res1 != 0 || mpfr_erangeflag_p () == 0)
    366  1.1.1.1.8.1  tls     {
    367  1.1.1.1.8.1  tls       printf ("Error for mpfr_cmp_f (NaN, 17)\n");
    368  1.1.1.1.8.1  tls       printf ("Return value: expected 0, got %d\n", res1);
    369  1.1.1.1.8.1  tls       printf ("Erange flag: expected set, got %d\n", mpfr_erangeflag_p ());
    370  1.1.1.1.8.1  tls       exit (1);
    371  1.1.1.1.8.1  tls     }
    372  1.1.1.1.8.1  tls 
    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.1.8.1  tls test_special2z (int (*mpfr_func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
    478  1.1.1.1.8.1  tls                void (*mpz_func)(mpz_ptr, mpz_srcptr, mpz_srcptr),
    479  1.1.1.1.8.1  tls                const char *op)
    480  1.1.1.1.8.1  tls {
    481  1.1.1.1.8.1  tls   mpfr_t x1, x2;
    482  1.1.1.1.8.1  tls   mpz_t  z1, z2;
    483  1.1.1.1.8.1  tls   int res;
    484  1.1.1.1.8.1  tls 
    485  1.1.1.1.8.1  tls   mpfr_inits2 (128, x1, x2, (mpfr_ptr) 0);
    486  1.1.1.1.8.1  tls   mpz_init (z1); mpz_init(z2);
    487  1.1.1.1.8.1  tls   mpz_fac_ui (z1, 19); /* 19!+1 fits perfectly in a 128 bits mantissa */
    488  1.1.1.1.8.1  tls   mpz_add_ui (z1, z1, 1);
    489  1.1.1.1.8.1  tls   mpz_fac_ui (z2, 20); /* 20!+1 fits perfectly in a 128 bits mantissa */
    490  1.1.1.1.8.1  tls   mpz_add_ui (z2, z2, 1);
    491  1.1.1.1.8.1  tls 
    492  1.1.1.1.8.1  tls   res = mpfr_set_z(x1, z1, MPFR_RNDN);
    493  1.1.1.1.8.1  tls   if (res)
    494  1.1.1.1.8.1  tls     {
    495  1.1.1.1.8.1  tls       printf("Special2z %s: set_z1 error\n", op);
    496  1.1.1.1.8.1  tls       exit(1);
    497  1.1.1.1.8.1  tls     }
    498  1.1.1.1.8.1  tls   mpfr_set_z (x2, z2, MPFR_RNDN);
    499  1.1.1.1.8.1  tls   if (res)
    500  1.1.1.1.8.1  tls     {
    501  1.1.1.1.8.1  tls       printf("Special2z %s: set_z2 error\n", op);
    502  1.1.1.1.8.1  tls       exit(1);
    503  1.1.1.1.8.1  tls     }
    504  1.1.1.1.8.1  tls 
    505  1.1.1.1.8.1  tls   /* (19!+1) * (20!+1) fits in a 128 bits number */
    506  1.1.1.1.8.1  tls   res = mpfr_func(x1, z1, x2, MPFR_RNDN);
    507  1.1.1.1.8.1  tls   if (res)
    508  1.1.1.1.8.1  tls     {
    509  1.1.1.1.8.1  tls       printf("Special2z %s: wrong inexact flag.\n", op);
    510  1.1.1.1.8.1  tls       exit(1);
    511  1.1.1.1.8.1  tls     }
    512  1.1.1.1.8.1  tls   mpz_func(z1, z1, z2);
    513  1.1.1.1.8.1  tls   res = mpfr_set_z (x2, z1, MPFR_RNDN);
    514  1.1.1.1.8.1  tls   if (res)
    515  1.1.1.1.8.1  tls     {
    516  1.1.1.1.8.1  tls       printf("Special2z %s: set_z2 error\n", op);
    517  1.1.1.1.8.1  tls       exit(1);
    518  1.1.1.1.8.1  tls     }
    519  1.1.1.1.8.1  tls   if (mpfr_cmp(x1, x2))
    520  1.1.1.1.8.1  tls     {
    521  1.1.1.1.8.1  tls       printf("Special2z %s: results differ.\nx1=", op);
    522  1.1.1.1.8.1  tls       mpfr_print_binary(x1);
    523  1.1.1.1.8.1  tls       printf("\nx2=");
    524  1.1.1.1.8.1  tls       mpfr_print_binary(x2);
    525  1.1.1.1.8.1  tls       printf ("\nZ2=");
    526  1.1.1.1.8.1  tls       mpz_out_str (stdout, 2, z1);
    527  1.1.1.1.8.1  tls       putchar('\n');
    528  1.1.1.1.8.1  tls       exit(1);
    529  1.1.1.1.8.1  tls     }
    530  1.1.1.1.8.1  tls 
    531  1.1.1.1.8.1  tls   mpz_set_ui (z1, 0);
    532  1.1.1.1.8.1  tls   mpz_set_ui (z2, 1);
    533  1.1.1.1.8.1  tls   mpfr_set_ui (x2, 1, MPFR_RNDN);
    534  1.1.1.1.8.1  tls   res = mpfr_func(x1, z1, x2, MPFR_RNDN);
    535  1.1.1.1.8.1  tls   mpz_func (z1, z1, z2);
    536  1.1.1.1.8.1  tls   mpfr_set_z (x2, z1, MPFR_RNDN);
    537  1.1.1.1.8.1  tls   if (mpfr_cmp(x1, x2))
    538  1.1.1.1.8.1  tls     {
    539  1.1.1.1.8.1  tls       printf("Special2z %s: results differ(2).\nx1=", op);
    540  1.1.1.1.8.1  tls       mpfr_print_binary(x1);
    541  1.1.1.1.8.1  tls       printf("\nx2=");
    542  1.1.1.1.8.1  tls       mpfr_print_binary(x2);
    543  1.1.1.1.8.1  tls       putchar('\n');
    544  1.1.1.1.8.1  tls       exit(1);
    545  1.1.1.1.8.1  tls     }
    546  1.1.1.1.8.1  tls 
    547  1.1.1.1.8.1  tls   mpz_clear (z1); mpz_clear(z2);
    548  1.1.1.1.8.1  tls   mpfr_clears (x1, x2, (mpfr_ptr) 0);
    549  1.1.1.1.8.1  tls }
    550  1.1.1.1.8.1  tls 
    551  1.1.1.1.8.1  tls 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.1.8.1  tls                   mpz_out_str (stdout, 10, arg2);
    591  1.1.1.1.8.1  tls                   printf ("\ngot      ");
    592  1.1.1.1.8.1  tls                   mpfr_dump (dst_small);
    593  1.1.1.1.8.1  tls                   printf ("expected ");
    594  1.1.1.1.8.1  tls                   mpfr_dump (tmp);
    595  1.1.1.1.8.1  tls                   printf ("approx   ");
    596  1.1.1.1.8.1  tls                   mpfr_dump (dst_big);
    597  1.1.1.1.8.1  tls                   exit (1);
    598  1.1.1.1.8.1  tls                 }
    599  1.1.1.1.8.1  tls               compare2 = mpfr_cmp (tmp, dst_big);
    600  1.1.1.1.8.1  tls               /* if rounding to nearest, cannot know the sign of t - f(x)
    601  1.1.1.1.8.1  tls                  because of composed rounding: y = o(f(x)) and t = o(y) */
    602  1.1.1.1.8.1  tls               if (compare * compare2 >= 0)
    603  1.1.1.1.8.1  tls                 compare = compare + compare2;
    604  1.1.1.1.8.1  tls               else
    605  1.1.1.1.8.1  tls                 compare = inexact; /* cannot determine sign(t-f(x)) */
    606  1.1.1.1.8.1  tls               if (((inexact == 0) && (compare != 0)) ||
    607  1.1.1.1.8.1  tls                   ((inexact > 0) && (compare <= 0)) ||
    608  1.1.1.1.8.1  tls                   ((inexact < 0) && (compare >= 0)))
    609  1.1.1.1.8.1  tls                 {
    610  1.1.1.1.8.1  tls                   printf ("Wrong inexact flag for rnd=%s and %s_z:\n"
    611  1.1.1.1.8.1  tls                           "expected %d, got %d\n",
    612  1.1.1.1.8.1  tls                           mpfr_print_rnd_mode (rnd), op, compare, inexact);
    613  1.1.1.1.8.1  tls                   printf ("\narg1="); mpfr_print_binary (arg1);
    614  1.1.1.1.8.1  tls                   printf ("\narg2="); mpz_out_str(stdout, 2, arg2);
    615  1.1.1.1.8.1  tls                   printf ("\ndstl="); mpfr_print_binary (dst_big);
    616  1.1.1.1.8.1  tls                   printf ("\ndsts="); mpfr_print_binary (dst_small);
    617  1.1.1.1.8.1  tls                   printf ("\ntmp ="); mpfr_dump (tmp);
    618  1.1.1.1.8.1  tls                   exit (1);
    619  1.1.1.1.8.1  tls                 }
    620  1.1.1.1.8.1  tls             }
    621  1.1.1.1.8.1  tls         }
    622  1.1.1.1.8.1  tls     }
    623  1.1.1.1.8.1  tls 
    624  1.1.1.1.8.1  tls   mpz_clear (arg2);
    625  1.1.1.1.8.1  tls   mpfr_clears (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    626  1.1.1.1.8.1  tls }
    627  1.1.1.1.8.1  tls 
    628  1.1.1.1.8.1  tls static void
    629  1.1.1.1.8.1  tls test_generic2z (mpfr_prec_t p0, mpfr_prec_t p1, unsigned int N,
    630  1.1.1.1.8.1  tls                int (*func)(mpfr_ptr, mpz_srcptr, mpfr_srcptr, mpfr_rnd_t),
    631  1.1.1.1.8.1  tls                const char *op)
    632  1.1.1.1.8.1  tls {
    633  1.1.1.1.8.1  tls   mpfr_prec_t prec;
    634  1.1.1.1.8.1  tls   mpfr_t arg1, dst_big, dst_small, tmp;
    635  1.1.1.1.8.1  tls   mpz_t  arg2;
    636  1.1.1.1.8.1  tls   mpfr_rnd_t rnd;
    637  1.1.1.1.8.1  tls   int inexact, compare, compare2;
    638  1.1.1.1.8.1  tls   unsigned int n;
    639  1.1.1.1.8.1  tls 
    640  1.1.1.1.8.1  tls   mpfr_inits (arg1, dst_big, dst_small, tmp, (mpfr_ptr) 0);
    641  1.1.1.1.8.1  tls   mpz_init (arg2);
    642  1.1.1.1.8.1  tls 
    643  1.1.1.1.8.1  tls   for (prec = p0; prec <= p1; prec++)
    644  1.1.1.1.8.1  tls     {
    645  1.1.1.1.8.1  tls       mpfr_set_prec (arg1, prec);
    646  1.1.1.1.8.1  tls       mpfr_set_prec (tmp, prec);
    647  1.1.1.1.8.1  tls       mpfr_set_prec (dst_small, prec);
    648  1.1.1.1.8.1  tls 
    649  1.1.1.1.8.1  tls       for (n=0; n<N; n++)
    650  1.1.1.1.8.1  tls         {
    651  1.1.1.1.8.1  tls           mpfr_urandomb (arg1, RANDS);
    652  1.1.1.1.8.1  tls           mpz_urandomb (arg2, RANDS, 1024);
    653  1.1.1.1.8.1  tls           rnd = RND_RAND ();
    654  1.1.1.1.8.1  tls           mpfr_set_prec (dst_big, 2*prec);
    655  1.1.1.1.8.1  tls           compare = func(dst_big, arg2, arg1, rnd);
    656  1.1.1.1.8.1  tls           if (mpfr_can_round (dst_big, 2*prec, rnd, rnd, prec))
    657  1.1.1.1.8.1  tls             {
    658  1.1.1.1.8.1  tls               mpfr_set (tmp, dst_big, rnd);
    659  1.1.1.1.8.1  tls               inexact = func(dst_small, arg2, arg1, rnd);
    660  1.1.1.1.8.1  tls               if (mpfr_cmp (tmp, dst_small))
    661  1.1.1.1.8.1  tls                 {
    662  1.1.1.1.8.1  tls                   printf ("Results differ for prec=%u rnd_mode=%s and %s_z:\n"
    663  1.1.1.1.8.1  tls                           "arg1=",
    664  1.1.1.1.8.1  tls                           (unsigned) prec, mpfr_print_rnd_mode (rnd), op);
    665  1.1.1.1.8.1  tls                   mpfr_print_binary (arg1);
    666  1.1.1.1.8.1  tls                   printf("\narg2=");
    667  1.1.1.1.8.1  tls                   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.1.8.1  tls               printf("Range error for prec=%lu and %s",
    817  1.1.1.1.8.1  tls                      (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.1.8.1  tls static void
    834  1.1.1.1.8.1  tls bug_mul_q_20100810 (void)
    835  1.1.1.1.8.1  tls {
    836  1.1.1.1.8.1  tls   mpfr_t x;
    837  1.1.1.1.8.1  tls   mpfr_t y;
    838  1.1.1.1.8.1  tls   mpq_t q;
    839  1.1.1.1.8.1  tls   int inexact;
    840  1.1.1.1.8.1  tls 
    841  1.1.1.1.8.1  tls   mpfr_init (x);
    842  1.1.1.1.8.1  tls   mpfr_init (y);
    843  1.1.1.1.8.1  tls   mpq_init (q);
    844  1.1.1.1.8.1  tls 
    845  1.1.1.1.8.1  tls   /* mpfr_mul_q: the inexact value must be set in case of overflow */
    846  1.1.1.1.8.1  tls   mpq_set_ui (q, 4096, 3);
    847  1.1.1.1.8.1  tls   mpfr_set_inf (x, +1);
    848  1.1.1.1.8.1  tls   mpfr_nextbelow (x);
    849  1.1.1.1.8.1  tls   inexact = mpfr_mul_q (y, x, q, MPFR_RNDU);
    850  1.1.1.1.8.1  tls 
    851  1.1.1.1.8.1  tls   if (inexact <= 0)
    852  1.1.1.1.8.1  tls     {
    853  1.1.1.1.8.1  tls       printf ("Overflow error in mpfr_mul_q. ");
    854  1.1.1.1.8.1  tls       printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
    855  1.1.1.1.8.1  tls 
    856  1.1.1.1.8.1  tls       exit (1);
    857  1.1.1.1.8.1  tls     }
    858  1.1.1.1.8.1  tls   if (!mpfr_inf_p (y))
    859  1.1.1.1.8.1  tls     {
    860  1.1.1.1.8.1  tls       printf ("Overflow error in mpfr_mul_q (y, x, q, MPFR_RNDD). ");
    861  1.1.1.1.8.1  tls       printf ("\nx = ");
    862  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
    863  1.1.1.1.8.1  tls       printf ("\nq = ");
    864  1.1.1.1.8.1  tls       mpq_out_str (stdout, 10, q);
    865  1.1.1.1.8.1  tls       printf ("\ny = ");
    866  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
    867  1.1.1.1.8.1  tls       printf (" (should be +infinity)\n");
    868  1.1.1.1.8.1  tls 
    869  1.1.1.1.8.1  tls       exit (1);
    870  1.1.1.1.8.1  tls     }
    871  1.1.1.1.8.1  tls 
    872  1.1.1.1.8.1  tls   mpq_clear (q);
    873  1.1.1.1.8.1  tls   mpfr_clear (y);
    874  1.1.1.1.8.1  tls   mpfr_clear (x);
    875  1.1.1.1.8.1  tls }
    876  1.1.1.1.8.1  tls 
    877  1.1.1.1.8.1  tls static void
    878  1.1.1.1.8.1  tls bug_div_q_20100810 (void)
    879  1.1.1.1.8.1  tls {
    880  1.1.1.1.8.1  tls   mpfr_t x;
    881  1.1.1.1.8.1  tls   mpfr_t y;
    882  1.1.1.1.8.1  tls   mpq_t q;
    883  1.1.1.1.8.1  tls   int inexact;
    884  1.1.1.1.8.1  tls 
    885  1.1.1.1.8.1  tls   mpfr_init (x);
    886  1.1.1.1.8.1  tls   mpfr_init (y);
    887  1.1.1.1.8.1  tls   mpq_init (q);
    888  1.1.1.1.8.1  tls 
    889  1.1.1.1.8.1  tls   /* mpfr_div_q: the inexact value must be set in case of overflow */
    890  1.1.1.1.8.1  tls   mpq_set_ui (q, 3, 4096);
    891  1.1.1.1.8.1  tls   mpfr_set_inf (x, +1);
    892  1.1.1.1.8.1  tls   mpfr_nextbelow (x);
    893  1.1.1.1.8.1  tls   inexact = mpfr_div_q (y, x, q, MPFR_RNDU);
    894  1.1.1.1.8.1  tls 
    895  1.1.1.1.8.1  tls   if (inexact <= 0)
    896  1.1.1.1.8.1  tls     {
    897  1.1.1.1.8.1  tls       printf ("Overflow error in mpfr_div_q. ");
    898  1.1.1.1.8.1  tls       printf ("Wrong inexact flag: got %d, should be positive.\n", inexact);
    899  1.1.1.1.8.1  tls 
    900  1.1.1.1.8.1  tls       exit (1);
    901  1.1.1.1.8.1  tls     }
    902  1.1.1.1.8.1  tls   if (!mpfr_inf_p (y))
    903  1.1.1.1.8.1  tls     {
    904  1.1.1.1.8.1  tls       printf ("Overflow error in mpfr_div_q (y, x, q, MPFR_RNDD). ");
    905  1.1.1.1.8.1  tls       printf ("\nx = ");
    906  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD);
    907  1.1.1.1.8.1  tls       printf ("\nq = ");
    908  1.1.1.1.8.1  tls       mpq_out_str (stdout, 10, q);
    909  1.1.1.1.8.1  tls       printf ("\ny = ");
    910  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD);
    911  1.1.1.1.8.1  tls       printf (" (should be +infinity)\n");
    912  1.1.1.1.8.1  tls 
    913  1.1.1.1.8.1  tls       exit (1);
    914  1.1.1.1.8.1  tls     }
    915  1.1.1.1.8.1  tls 
    916  1.1.1.1.8.1  tls   mpq_clear (q);
    917  1.1.1.1.8.1  tls   mpfr_clear (y);
    918  1.1.1.1.8.1  tls   mpfr_clear (x);
    919  1.1.1.1.8.1  tls }
    920  1.1.1.1.8.1  tls 
    921  1.1.1.1.8.1  tls static void
    922  1.1.1.1.8.1  tls bug_mul_div_q_20100818 (void)
    923  1.1.1.1.8.1  tls {
    924  1.1.1.1.8.1  tls   mpq_t qa, qb;
    925  1.1.1.1.8.1  tls   mpfr_t x1, x2, y1, y2, y3;
    926  1.1.1.1.8.1  tls   mpfr_exp_t emin, emax, e;
    927  1.1.1.1.8.1  tls   int inex;
    928  1.1.1.1.8.1  tls   int rnd;
    929  1.1.1.1.8.1  tls 
    930  1.1.1.1.8.1  tls   emin = mpfr_get_emin ();
    931  1.1.1.1.8.1  tls   emax = mpfr_get_emax ();
    932  1.1.1.1.8.1  tls   set_emin (MPFR_EMIN_MIN);
    933  1.1.1.1.8.1  tls   set_emax (MPFR_EMAX_MAX);
    934  1.1.1.1.8.1  tls 
    935  1.1.1.1.8.1  tls   mpq_init (qa);
    936  1.1.1.1.8.1  tls   mpq_init (qb);
    937  1.1.1.1.8.1  tls   mpfr_inits2 (32, x1, x2, y1, y2, y3, (mpfr_ptr) 0);
    938  1.1.1.1.8.1  tls 
    939  1.1.1.1.8.1  tls   mpq_set_ui (qa, 3, 17);
    940  1.1.1.1.8.1  tls   mpq_set_ui (qb, 17, 3);
    941  1.1.1.1.8.1  tls   inex = mpfr_set_ui (x1, 7, MPFR_RNDN);
    942  1.1.1.1.8.1  tls   MPFR_ASSERTN (inex == 0);
    943  1.1.1.1.8.1  tls 
    944  1.1.1.1.8.1  tls   e = MPFR_EMAX_MAX - 3;
    945  1.1.1.1.8.1  tls   inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN);  /* x2 = x1 * 2^e */
    946  1.1.1.1.8.1  tls   MPFR_ASSERTN (inex == 0);
    947  1.1.1.1.8.1  tls 
    948  1.1.1.1.8.1  tls   RND_LOOP(rnd)
    949  1.1.1.1.8.1  tls     {
    950  1.1.1.1.8.1  tls       mpfr_mul_q (y1, x1, qa, (mpfr_rnd_t) rnd);
    951  1.1.1.1.8.1  tls       mpfr_div_q (y3, x1, qb, (mpfr_rnd_t) rnd);
    952  1.1.1.1.8.1  tls       MPFR_ASSERTN (mpfr_equal_p (y1, y3));
    953  1.1.1.1.8.1  tls       inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
    954  1.1.1.1.8.1  tls       MPFR_ASSERTN (inex == 0);
    955  1.1.1.1.8.1  tls       inex = mpfr_mul (y3, y3, y1, MPFR_RNDN);  /* y3 = y1 * 2^e */
    956  1.1.1.1.8.1  tls       MPFR_ASSERTN (inex == 0);
    957  1.1.1.1.8.1  tls       mpfr_mul_q (y2, x2, qa, (mpfr_rnd_t) rnd);
    958  1.1.1.1.8.1  tls       if (! mpfr_equal_p (y2, y3))
    959  1.1.1.1.8.1  tls         {
    960  1.1.1.1.8.1  tls           printf ("Error 1 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    961  1.1.1.1.8.1  tls           printf ("Expected "); mpfr_dump (y3);
    962  1.1.1.1.8.1  tls           printf ("Got      "); mpfr_dump (y2);
    963  1.1.1.1.8.1  tls           exit (1);
    964  1.1.1.1.8.1  tls         }
    965  1.1.1.1.8.1  tls       mpfr_div_q (y2, x2, qb, (mpfr_rnd_t) rnd);
    966  1.1.1.1.8.1  tls       if (! mpfr_equal_p (y2, y3))
    967  1.1.1.1.8.1  tls         {
    968  1.1.1.1.8.1  tls           printf ("Error 2 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    969  1.1.1.1.8.1  tls           printf ("Expected "); mpfr_dump (y3);
    970  1.1.1.1.8.1  tls           printf ("Got      "); mpfr_dump (y2);
    971  1.1.1.1.8.1  tls           exit (1);
    972  1.1.1.1.8.1  tls         }
    973  1.1.1.1.8.1  tls     }
    974  1.1.1.1.8.1  tls 
    975  1.1.1.1.8.1  tls   e = MPFR_EMIN_MIN;
    976  1.1.1.1.8.1  tls   inex = mpfr_set_ui_2exp (x2, 7, e, MPFR_RNDN);  /* x2 = x1 * 2^e */
    977  1.1.1.1.8.1  tls   MPFR_ASSERTN (inex == 0);
    978  1.1.1.1.8.1  tls 
    979  1.1.1.1.8.1  tls   RND_LOOP(rnd)
    980  1.1.1.1.8.1  tls     {
    981  1.1.1.1.8.1  tls       mpfr_div_q (y1, x1, qa, (mpfr_rnd_t) rnd);
    982  1.1.1.1.8.1  tls       mpfr_mul_q (y3, x1, qb, (mpfr_rnd_t) rnd);
    983  1.1.1.1.8.1  tls       MPFR_ASSERTN (mpfr_equal_p (y1, y3));
    984  1.1.1.1.8.1  tls       inex = mpfr_set_ui_2exp (y3, 1, e, MPFR_RNDN);
    985  1.1.1.1.8.1  tls       MPFR_ASSERTN (inex == 0);
    986  1.1.1.1.8.1  tls       inex = mpfr_mul (y3, y3, y1, MPFR_RNDN);  /* y3 = y1 * 2^e */
    987  1.1.1.1.8.1  tls       MPFR_ASSERTN (inex == 0);
    988  1.1.1.1.8.1  tls       mpfr_div_q (y2, x2, qa, (mpfr_rnd_t) rnd);
    989  1.1.1.1.8.1  tls       if (! mpfr_equal_p (y2, y3))
    990  1.1.1.1.8.1  tls         {
    991  1.1.1.1.8.1  tls           printf ("Error 3 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
    992  1.1.1.1.8.1  tls           printf ("Expected "); mpfr_dump (y3);
    993  1.1.1.1.8.1  tls           printf ("Got      "); mpfr_dump (y2);
    994  1.1.1.1.8.1  tls           exit (1);
    995  1.1.1.1.8.1  tls         }
    996  1.1.1.1.8.1  tls       mpfr_mul_q (y2, x2, qb, (mpfr_rnd_t) rnd);
    997  1.1.1.1.8.1  tls       if (! mpfr_equal_p (y2, y3))
    998  1.1.1.1.8.1  tls         {
    999  1.1.1.1.8.1  tls           printf ("Error 4 in bug_mul_div_q_20100818 (rnd = %d)\n", rnd);
   1000  1.1.1.1.8.1  tls           printf ("Expected "); mpfr_dump (y3);
   1001  1.1.1.1.8.1  tls           printf ("Got      "); mpfr_dump (y2);
   1002  1.1.1.1.8.1  tls           exit (1);
   1003  1.1.1.1.8.1  tls         }
   1004  1.1.1.1.8.1  tls     }
   1005  1.1.1.1.8.1  tls 
   1006  1.1.1.1.8.1  tls   mpq_clear (qa);
   1007  1.1.1.1.8.1  tls   mpq_clear (qb);
   1008  1.1.1.1.8.1  tls   mpfr_clears (x1, x2, y1, y2, y3, (mpfr_ptr) 0);
   1009  1.1.1.1.8.1  tls 
   1010  1.1.1.1.8.1  tls   set_emin (emin);
   1011  1.1.1.1.8.1  tls   set_emax (emax);
   1012  1.1.1.1.8.1  tls }
   1013  1.1.1.1.8.1  tls 
   1014  1.1.1.1.8.1  tls static void
   1015  1.1.1.1.8.1  tls reduced_expo_range (void)
   1016  1.1.1.1.8.1  tls {
   1017  1.1.1.1.8.1  tls   mpfr_t x;
   1018  1.1.1.1.8.1  tls   mpz_t z;
   1019  1.1.1.1.8.1  tls   mpq_t q;
   1020  1.1.1.1.8.1  tls   mpfr_exp_t emin;
   1021  1.1.1.1.8.1  tls   int inex;
   1022  1.1.1.1.8.1  tls 
   1023  1.1.1.1.8.1  tls   emin = mpfr_get_emin ();
   1024  1.1.1.1.8.1  tls   set_emin (4);
   1025  1.1.1.1.8.1  tls 
   1026  1.1.1.1.8.1  tls   mpfr_init2 (x, 32);
   1027  1.1.1.1.8.1  tls 
   1028  1.1.1.1.8.1  tls   mpz_init (z);
   1029  1.1.1.1.8.1  tls   mpfr_clear_flags ();
   1030  1.1.1.1.8.1  tls   inex = mpfr_set_ui (x, 17, MPFR_RNDN);
   1031  1.1.1.1.8.1  tls   MPFR_ASSERTN (inex == 0);
   1032  1.1.1.1.8.1  tls   mpz_set_ui (z, 3);
   1033  1.1.1.1.8.1  tls   inex = mpfr_mul_z (x, x, z, MPFR_RNDN);
   1034  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 51) != 0)
   1035  1.1.1.1.8.1  tls     {
   1036  1.1.1.1.8.1  tls       printf ("Error 1 in reduce_expo_range: expected 51 with inex = 0,"
   1037  1.1.1.1.8.1  tls               " got\n");
   1038  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1039  1.1.1.1.8.1  tls       printf ("with inex = %d\n", inex);
   1040  1.1.1.1.8.1  tls       exit (1);
   1041  1.1.1.1.8.1  tls     }
   1042  1.1.1.1.8.1  tls   inex = mpfr_div_z (x, x, z, MPFR_RNDN);
   1043  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1044  1.1.1.1.8.1  tls     {
   1045  1.1.1.1.8.1  tls       printf ("Error 2 in reduce_expo_range: expected 17 with inex = 0,"
   1046  1.1.1.1.8.1  tls               " got\n");
   1047  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1048  1.1.1.1.8.1  tls       printf ("with inex = %d\n", inex);
   1049  1.1.1.1.8.1  tls       exit (1);
   1050  1.1.1.1.8.1  tls     }
   1051  1.1.1.1.8.1  tls   inex = mpfr_add_z (x, x, z, MPFR_RNDN);
   1052  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 20) != 0)
   1053  1.1.1.1.8.1  tls     {
   1054  1.1.1.1.8.1  tls       printf ("Error 3 in reduce_expo_range: expected 20 with inex = 0,"
   1055  1.1.1.1.8.1  tls               " got\n");
   1056  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1057  1.1.1.1.8.1  tls       printf ("with inex = %d\n", inex);
   1058  1.1.1.1.8.1  tls       exit (1);
   1059  1.1.1.1.8.1  tls     }
   1060  1.1.1.1.8.1  tls   inex = mpfr_sub_z (x, x, z, MPFR_RNDN);
   1061  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1062  1.1.1.1.8.1  tls     {
   1063  1.1.1.1.8.1  tls       printf ("Error 4 in reduce_expo_range: expected 17 with inex = 0,"
   1064  1.1.1.1.8.1  tls               " got\n");
   1065  1.1.1.1.8.1  tls       mpfr_out_str (stdout, 10, 0, x, MPFR_RNDN);
   1066  1.1.1.1.8.1  tls       printf ("with inex = %d\n", inex);
   1067  1.1.1.1.8.1  tls       exit (1);
   1068  1.1.1.1.8.1  tls     }
   1069  1.1.1.1.8.1  tls   MPFR_ASSERTN (__gmpfr_flags == 0);
   1070  1.1.1.1.8.1  tls   if (mpfr_cmp_z (x, z) <= 0)
   1071  1.1.1.1.8.1  tls     {
   1072  1.1.1.1.8.1  tls       printf ("Error 5 in reduce_expo_range: expected a positive value.\n");
   1073  1.1.1.1.8.1  tls       exit (1);
   1074  1.1.1.1.8.1  tls     }
   1075  1.1.1.1.8.1  tls   mpz_clear (z);
   1076  1.1.1.1.8.1  tls 
   1077  1.1.1.1.8.1  tls   mpq_init (q);
   1078  1.1.1.1.8.1  tls   mpq_set_ui (q, 1, 1);
   1079  1.1.1.1.8.1  tls   mpfr_set_ui (x, 16, MPFR_RNDN);
   1080  1.1.1.1.8.1  tls   inex = mpfr_add_q (x, x, q, MPFR_RNDN);
   1081  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 17) != 0)
   1082  1.1.1.1.8.1  tls     {
   1083  1.1.1.1.8.1  tls       printf ("Error in reduce_expo_range for 16 + 1/1,"
   1084  1.1.1.1.8.1  tls               " got inex = %d and\nx = ", inex);
   1085  1.1.1.1.8.1  tls       mpfr_dump (x);
   1086  1.1.1.1.8.1  tls       exit (1);
   1087  1.1.1.1.8.1  tls     }
   1088  1.1.1.1.8.1  tls   inex = mpfr_sub_q (x, x, q, MPFR_RNDN);
   1089  1.1.1.1.8.1  tls   if (inex != 0 || MPFR_IS_NAN (x) || mpfr_cmp_ui (x, 16) != 0)
   1090  1.1.1.1.8.1  tls     {
   1091  1.1.1.1.8.1  tls       printf ("Error in reduce_expo_range for 17 - 1/1,"
   1092  1.1.1.1.8.1  tls               " got inex = %d and\nx = ", inex);
   1093  1.1.1.1.8.1  tls       mpfr_dump (x);
   1094  1.1.1.1.8.1  tls       exit (1);
   1095  1.1.1.1.8.1  tls     }
   1096  1.1.1.1.8.1  tls   mpq_clear (q);
   1097  1.1.1.1.8.1  tls 
   1098  1.1.1.1.8.1  tls   mpfr_clear (x);
   1099  1.1.1.1.8.1  tls 
   1100  1.1.1.1.8.1  tls   set_emin (emin);
   1101  1.1.1.1.8.1  tls }
   1102  1.1.1.1.8.1  tls 
   1103  1.1.1.1.8.1  tls static void
   1104  1.1.1.1.8.1  tls addsubq_overflow_aux (mpfr_exp_t e)
   1105  1.1.1.1.8.1  tls {
   1106  1.1.1.1.8.1  tls   mpfr_t x, y;
   1107  1.1.1.1.8.1  tls   mpq_t q;
   1108  1.1.1.1.8.1  tls   mpfr_exp_t emax;
   1109  1.1.1.1.8.1  tls   int inex;
   1110  1.1.1.1.8.1  tls   int rnd;
   1111  1.1.1.1.8.1  tls   int sign, sub;
   1112  1.1.1.1.8.1  tls 
   1113  1.1.1.1.8.1  tls   MPFR_ASSERTN (e <= LONG_MAX);
   1114  1.1.1.1.8.1  tls   emax = mpfr_get_emax ();
   1115  1.1.1.1.8.1  tls   set_emax (e);
   1116  1.1.1.1.8.1  tls   mpfr_inits2 (16, x, y, (mpfr_ptr) 0);
   1117  1.1.1.1.8.1  tls   mpq_init (q);
   1118  1.1.1.1.8.1  tls 
   1119  1.1.1.1.8.1  tls   mpfr_set_inf (x, 1);
   1120  1.1.1.1.8.1  tls   mpfr_nextbelow (x);
   1121  1.1.1.1.8.1  tls   mpq_set_ui (q, 1, 1);
   1122  1.1.1.1.8.1  tls 
   1123  1.1.1.1.8.1  tls   for (sign = 0; sign <= 1; sign++)
   1124  1.1.1.1.8.1  tls     {
   1125  1.1.1.1.8.1  tls       for (sub = 0; sub <= 1; sub++)
   1126  1.1.1.1.8.1  tls         {
   1127  1.1.1.1.8.1  tls           RND_LOOP(rnd)
   1128  1.1.1.1.8.1  tls             {
   1129  1.1.1.1.8.1  tls               unsigned int flags, ex_flags;
   1130  1.1.1.1.8.1  tls               int inf;
   1131  1.1.1.1.8.1  tls 
   1132  1.1.1.1.8.1  tls               inf = rnd == MPFR_RNDA ||
   1133  1.1.1.1.8.1  tls                     rnd == (sign ? MPFR_RNDD : MPFR_RNDU);
   1134  1.1.1.1.8.1  tls               ex_flags = MPFR_FLAGS_INEXACT | (inf ? MPFR_FLAGS_OVERFLOW : 0);
   1135  1.1.1.1.8.1  tls               mpfr_clear_flags ();
   1136  1.1.1.1.8.1  tls               inex = sub ?
   1137  1.1.1.1.8.1  tls                 mpfr_sub_q (y, x, q, (mpfr_rnd_t) rnd) :
   1138  1.1.1.1.8.1  tls                 mpfr_add_q (y, x, q, (mpfr_rnd_t) rnd);
   1139  1.1.1.1.8.1  tls               flags = __gmpfr_flags;
   1140  1.1.1.1.8.1  tls               if (inex == 0 || flags != ex_flags ||
   1141  1.1.1.1.8.1  tls                   (inf ? ! mpfr_inf_p (y) : ! mpfr_equal_p (x, y)))
   1142  1.1.1.1.8.1  tls                 {
   1143  1.1.1.1.8.1  tls                   printf ("Error in addsubq_overflow_aux(%ld),"
   1144  1.1.1.1.8.1  tls                           " sign = %d, %s\n", (long) e, sign,
   1145  1.1.1.1.8.1  tls                           mpfr_print_rnd_mode ((mpfr_rnd_t) rnd));
   1146  1.1.1.1.8.1  tls                   printf ("Got inex = %d, y = ", inex);
   1147  1.1.1.1.8.1  tls                   mpfr_dump (y);
   1148  1.1.1.1.8.1  tls                   printf ("Expected flags:");
   1149  1.1.1.1.8.1  tls                   flags_out (ex_flags);
   1150  1.1.1.1.8.1  tls                   printf ("Got flags:     ");
   1151  1.1.1.1.8.1  tls                   flags_out (flags);
   1152  1.1.1.1.8.1  tls                   exit (1);
   1153  1.1.1.1.8.1  tls                 }
   1154  1.1.1.1.8.1  tls             }
   1155  1.1.1.1.8.1  tls           mpq_neg (q, q);
   1156  1.1.1.1.8.1  tls         }
   1157  1.1.1.1.8.1  tls       mpfr_neg (x, x, MPFR_RNDN);
   1158  1.1.1.1.8.1  tls       mpq_neg (q, q);
   1159  1.1.1.1.8.1  tls     }
   1160  1.1.1.1.8.1  tls 
   1161  1.1.1.1.8.1  tls   mpq_clear (q);
   1162  1.1.1.1.8.1  tls   mpfr_clears (x, y, (mpfr_ptr) 0);
   1163  1.1.1.1.8.1  tls   set_emax (emax);
   1164  1.1.1.1.8.1  tls }
   1165  1.1.1.1.8.1  tls 
   1166  1.1.1.1.8.1  tls static void
   1167  1.1.1.1.8.1  tls addsubq_overflow (void)
   1168  1.1.1.1.8.1  tls {
   1169  1.1.1.1.8.1  tls   addsubq_overflow_aux (4913);
   1170  1.1.1.1.8.1  tls   addsubq_overflow_aux (MPFR_EMAX_MAX);
   1171  1.1.1.1.8.1  tls }
   1172  1.1.1.1.8.1  tls 
   1173  1.1.1.1.8.1  tls static void
   1174  1.1.1.1.8.1  tls coverage_mpfr_mul_q_20110218 (void)
   1175  1.1.1.1.8.1  tls {
   1176  1.1.1.1.8.1  tls   mpfr_t cmp, res, op1;
   1177  1.1.1.1.8.1  tls   mpq_t op2;
   1178  1.1.1.1.8.1  tls   int status;
   1179  1.1.1.1.8.1  tls 
   1180  1.1.1.1.8.1  tls   mpfr_init2 (cmp, MPFR_PREC_MIN);
   1181  1.1.1.1.8.1  tls   mpfr_init2 (res, MPFR_PREC_MIN);
   1182  1.1.1.1.8.1  tls   mpfr_init_set_si (op1, 1, MPFR_RNDN);
   1183  1.1.1.1.8.1  tls 
   1184  1.1.1.1.8.1  tls   mpq_init (op2);
   1185  1.1.1.1.8.1  tls   mpq_set_si (op2, 0, 0);
   1186  1.1.1.1.8.1  tls   mpz_set_si (mpq_denref (op2), 0);
   1187  1.1.1.1.8.1  tls 
   1188  1.1.1.1.8.1  tls   status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
   1189  1.1.1.1.8.1  tls 
   1190  1.1.1.1.8.1  tls   if ((status != 0) || (mpfr_cmp (cmp, res) != 0))
   1191  1.1.1.1.8.1  tls     {
   1192  1.1.1.1.8.1  tls       printf ("Results differ %d.\nres=", status);
   1193  1.1.1.1.8.1  tls       mpfr_print_binary (res);
   1194  1.1.1.1.8.1  tls       printf ("\ncmp=");
   1195  1.1.1.1.8.1  tls       mpfr_print_binary (cmp);
   1196  1.1.1.1.8.1  tls       putchar ('\n');
   1197  1.1.1.1.8.1  tls       exit (1);
   1198  1.1.1.1.8.1  tls     }
   1199  1.1.1.1.8.1  tls 
   1200  1.1.1.1.8.1  tls   mpfr_set_si (op1, 1, MPFR_RNDN);
   1201  1.1.1.1.8.1  tls   mpq_set_si (op2, -1, 0);
   1202  1.1.1.1.8.1  tls 
   1203  1.1.1.1.8.1  tls   status = mpfr_mul_q (res, op1, op2, MPFR_RNDN);
   1204  1.1.1.1.8.1  tls 
   1205  1.1.1.1.8.1  tls   mpfr_set_inf (cmp, -1);
   1206  1.1.1.1.8.1  tls   if ((status != 0) || (mpfr_cmp(res, cmp) != 0))
   1207  1.1.1.1.8.1  tls     {
   1208  1.1.1.1.8.1  tls       printf ("mpfr_mul_q 1 * (-1/0) returned a wrong value :\n waiting for ");
   1209  1.1.1.1.8.1  tls       mpfr_print_binary (cmp);
   1210  1.1.1.1.8.1  tls       printf (" got ");
   1211  1.1.1.1.8.1  tls       mpfr_print_binary (res);
   1212  1.1.1.1.8.1  tls       printf ("\n trinary value is %d\n", status);
   1213  1.1.1.1.8.1  tls       exit (1);
   1214  1.1.1.1.8.1  tls     }
   1215  1.1.1.1.8.1  tls 
   1216  1.1.1.1.8.1  tls   mpq_clear (op2);
   1217  1.1.1.1.8.1  tls   mpfr_clear (op1);
   1218  1.1.1.1.8.1  tls   mpfr_clear (res);
   1219  1.1.1.1.8.1  tls   mpfr_clear (cmp);
   1220  1.1.1.1.8.1  tls }
   1221  1.1.1.1.8.1  tls 
   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.1.8.1  tls   test_special2z (mpfr_z_sub, mpz_sub, "sub");
   1237  1.1.1.1.8.1  tls   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.1.8.1  tls   bug_mul_q_20100810 ();
   1255  1.1.1.1.8.1  tls   bug_div_q_20100810 ();
   1256  1.1.1.1.8.1  tls   bug_mul_div_q_20100818 ();
   1257  1.1.1.1.8.1  tls   reduced_expo_range ();
   1258  1.1.1.1.8.1  tls   addsubq_overflow ();
   1259  1.1.1.1.8.1  tls 
   1260  1.1.1.1.8.1  tls   coverage_mpfr_mul_q_20110218 ();
   1261  1.1.1.1.8.1  tls 
   1262          1.1  mrg   tests_end_mpfr ();
   1263          1.1  mrg   return 0;
   1264          1.1  mrg }
   1265          1.1  mrg 
   1266