1/* tremquo -- test file for mpfr_remquo and mpfr_remainder 2 3Copyright 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc. 4Contributed by the Arenaire and Cacao projects, INRIA. 5 6This file is part of the GNU MPFR Library. 7 8The GNU MPFR Library is free software; you can redistribute it and/or modify 9it under the terms of the GNU Lesser General Public License as published by 10the Free Software Foundation; either version 3 of the License, or (at your 11option) any later version. 12 13The GNU MPFR Library is distributed in the hope that it will be useful, but 14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 15or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 16License for more details. 17 18You should have received a copy of the GNU Lesser General Public License 19along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 20http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 2151 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23#include <stdio.h> 24#include <stdlib.h> 25 26#include "mpfr-test.h" 27 28static void 29bug20090227 (void) 30{ 31 mpfr_t x, y, r1, r2; 32 int inex1, inex2; 33 34 mpfr_init2 (x, 118); 35 mpfr_init2 (y, 181); 36 mpfr_init2 (r1, 140); 37 mpfr_init2 (r2, 140); 38 mpfr_set_si (x, -1, MPFR_RNDN); 39 mpfr_set_str_binary (y, "1.100100100001111110110101010001000100001011010001100001000110100110001001100011001100010100010111000000011011100000111001101000100101001000000100100111000001000100010100110011111010"); 40 inex1 = mpfr_remainder (r1, x, y, MPFR_RNDU); 41 /* since the quotient is -1, r1 is the rounding of x+y */ 42 inex2 = mpfr_add (r2, x, y, MPFR_RNDU); 43 if (mpfr_cmp (r1, r2)) 44 { 45 printf ("Error in mpfr_remainder (bug20090227)\n"); 46 printf ("Expected "); 47 mpfr_dump (r2); 48 printf ("Got "); 49 mpfr_dump (r1); 50 exit (1); 51 } 52 mpfr_clear (x); 53 mpfr_clear (y); 54 mpfr_clear (r1); 55 mpfr_clear (r2); 56} 57 58int 59main (int argc, char *argv[]) 60{ 61 mpfr_t x, y, r; 62 long q[1]; 63 64 if (argc == 3) /* usage: tremquo x y (rnd=MPFR_RNDN implicit) */ 65 { 66 mpfr_init2 (x, GMP_NUMB_BITS); 67 mpfr_init2 (y, GMP_NUMB_BITS); 68 mpfr_init2 (r, GMP_NUMB_BITS); 69 mpfr_set_str (x, argv[1], 10, MPFR_RNDN); 70 mpfr_set_str (y, argv[2], 10, MPFR_RNDN); 71 mpfr_remquo (r, q, x, y, MPFR_RNDN); 72 printf ("r="); 73 mpfr_out_str (stdout, 10, 0, r, MPFR_RNDN); 74 printf (" q=%ld\n", q[0]); 75 mpfr_clear (x); 76 mpfr_clear (y); 77 mpfr_clear (r); 78 return 0; 79 } 80 81 tests_start_mpfr (); 82 83 bug20090227 (); 84 85 mpfr_init (x); 86 mpfr_init (y); 87 mpfr_init (r); 88 89 /* special values */ 90 mpfr_set_nan (x); 91 mpfr_set_ui (y, 1, MPFR_RNDN); 92 mpfr_remquo (r, q, x, y, MPFR_RNDN); 93 MPFR_ASSERTN(mpfr_nan_p (r)); 94 95 mpfr_set_ui (x, 1, MPFR_RNDN); 96 mpfr_set_nan (y); 97 mpfr_remquo (r, q, x, y, MPFR_RNDN); 98 MPFR_ASSERTN(mpfr_nan_p (r)); 99 100 mpfr_set_inf (x, 1); /* +Inf */ 101 mpfr_set_ui (y, 1, MPFR_RNDN); 102 mpfr_remquo (r, q, x, y, MPFR_RNDN); 103 MPFR_ASSERTN (mpfr_nan_p (r)); 104 105 mpfr_set_inf (x, 1); /* +Inf */ 106 mpfr_set_ui (y, 0, MPFR_RNDN); 107 mpfr_remquo (r, q, x, y, MPFR_RNDN); 108 MPFR_ASSERTN (mpfr_nan_p (r)); 109 110 mpfr_set_inf (x, 1); /* +Inf */ 111 mpfr_set_inf (y, 1); 112 mpfr_remquo (r, q, x, y, MPFR_RNDN); 113 MPFR_ASSERTN (mpfr_nan_p (r)); 114 115 mpfr_set_ui (x, 0, MPFR_RNDN); 116 mpfr_set_inf (y, 1); 117 mpfr_remquo (r, q, x, y, MPFR_RNDN); 118 MPFR_ASSERTN (mpfr_cmp_ui (r, 0) == 0 && MPFR_IS_POS (r)); 119 MPFR_ASSERTN (q[0] == (long) 0); 120 121 mpfr_set_ui (x, 0, MPFR_RNDN); 122 mpfr_neg (x, x, MPFR_RNDN); /* -0 */ 123 mpfr_set_inf (y, 1); 124 mpfr_remquo (r, q, x, y, MPFR_RNDN); 125 MPFR_ASSERTN (mpfr_cmp_ui (r, 0) == 0 && MPFR_IS_NEG (r)); 126 MPFR_ASSERTN (q[0] == (long) 0); 127 128 mpfr_set_ui (x, 17, MPFR_RNDN); 129 mpfr_set_inf (y, 1); 130 mpfr_remquo (r, q, x, y, MPFR_RNDN); 131 MPFR_ASSERTN (mpfr_cmp (r, x) == 0); 132 MPFR_ASSERTN (q[0] == (long) 0); 133 134 mpfr_set_ui (x, 17, MPFR_RNDN); 135 mpfr_set_ui (y, 0, MPFR_RNDN); 136 mpfr_remquo (r, q, x, y, MPFR_RNDN); 137 MPFR_ASSERTN (mpfr_nan_p (r)); 138 139 mpfr_set_ui (x, 0, MPFR_RNDN); 140 mpfr_set_ui (y, 17, MPFR_RNDN); 141 mpfr_remquo (r, q, x, y, MPFR_RNDN); 142 MPFR_ASSERTN (mpfr_cmp_ui (r, 0) == 0 && MPFR_IS_POS (r)); 143 MPFR_ASSERTN (q[0] == (long) 0); 144 145 mpfr_set_ui (x, 0, MPFR_RNDN); 146 mpfr_neg (x, x, MPFR_RNDN); 147 mpfr_set_ui (y, 17, MPFR_RNDN); 148 mpfr_remquo (r, q, x, y, MPFR_RNDN); 149 MPFR_ASSERTN (mpfr_cmp_ui (r, 0) == 0 && MPFR_IS_NEG (r)); 150 MPFR_ASSERTN (q[0] == (long) 0); 151 152 mpfr_set_prec (x, 53); 153 mpfr_set_prec (y, 53); 154 155 /* check four possible sign combinations */ 156 mpfr_set_ui (x, 42, MPFR_RNDN); 157 mpfr_set_ui (y, 17, MPFR_RNDN); 158 mpfr_remquo (r, q, x, y, MPFR_RNDN); 159 MPFR_ASSERTN (mpfr_cmp_ui (r, 8) == 0); 160 MPFR_ASSERTN (q[0] == (long) 2); 161 mpfr_set_si (x, -42, MPFR_RNDN); 162 mpfr_set_ui (y, 17, MPFR_RNDN); 163 mpfr_remquo (r, q, x, y, MPFR_RNDN); 164 MPFR_ASSERTN (mpfr_cmp_si (r, -8) == 0); 165 MPFR_ASSERTN (q[0] == (long) -2); 166 mpfr_set_si (x, -42, MPFR_RNDN); 167 mpfr_set_si (y, -17, MPFR_RNDN); 168 mpfr_remquo (r, q, x, y, MPFR_RNDN); 169 MPFR_ASSERTN (mpfr_cmp_si (r, -8) == 0); 170 MPFR_ASSERTN (q[0] == (long) 2); 171 mpfr_set_ui (x, 42, MPFR_RNDN); 172 mpfr_set_si (y, -17, MPFR_RNDN); 173 mpfr_remquo (r, q, x, y, MPFR_RNDN); 174 MPFR_ASSERTN (mpfr_cmp_ui (r, 8) == 0); 175 MPFR_ASSERTN (q[0] == (long) -2); 176 177 mpfr_set_prec (x, 100); 178 mpfr_set_prec (y, 50); 179 mpfr_set_ui (x, 42, MPFR_RNDN); 180 mpfr_nextabove (x); /* 42 + 2^(-94) */ 181 mpfr_set_ui (y, 21, MPFR_RNDN); 182 mpfr_remquo (r, q, x, y, MPFR_RNDN); 183 MPFR_ASSERTN (mpfr_cmp_ui_2exp (r, 1, -94) == 0); 184 MPFR_ASSERTN (q[0] == (long) 2); 185 186 mpfr_set_prec (x, 50); 187 mpfr_set_prec (y, 100); 188 mpfr_set_ui (x, 42, MPFR_RNDN); 189 mpfr_nextabove (x); /* 42 + 2^(-44) */ 190 mpfr_set_ui (y, 21, MPFR_RNDN); 191 mpfr_remquo (r, q, x, y, MPFR_RNDN); 192 MPFR_ASSERTN (mpfr_cmp_ui_2exp (r, 1, -44) == 0); 193 MPFR_ASSERTN (q[0] == (long) 2); 194 195 mpfr_set_prec (x, 100); 196 mpfr_set_prec (y, 50); 197 mpfr_set_ui (x, 42, MPFR_RNDN); 198 mpfr_set_ui (y, 21, MPFR_RNDN); 199 mpfr_nextabove (y); /* 21 + 2^(-45) */ 200 mpfr_remquo (r, q, x, y, MPFR_RNDN); 201 /* r should be 42 - 2*(21 + 2^(-45)) = -2^(-44) */ 202 MPFR_ASSERTN (mpfr_cmp_si_2exp (r, -1, -44) == 0); 203 MPFR_ASSERTN (q[0] == (long) 2); 204 205 mpfr_set_prec (x, 50); 206 mpfr_set_prec (y, 100); 207 mpfr_set_ui (x, 42, MPFR_RNDN); 208 mpfr_set_ui (y, 21, MPFR_RNDN); 209 mpfr_nextabove (y); /* 21 + 2^(-95) */ 210 mpfr_remquo (r, q, x, y, MPFR_RNDN); 211 /* r should be 42 - 2*(21 + 2^(-95)) = -2^(-94) */ 212 MPFR_ASSERTN (mpfr_cmp_si_2exp (r, -1, -94) == 0); 213 MPFR_ASSERTN (q[0] == (long) 2); 214 215 /* exercise large quotient */ 216 mpfr_set_ui_2exp (x, 1, 65, MPFR_RNDN); 217 mpfr_set_ui (y, 1, MPFR_RNDN); 218 /* quotient is 2^65 */ 219 mpfr_remquo (r, q, x, y, MPFR_RNDN); 220 MPFR_ASSERTN (mpfr_cmp_si (r, 0) == 0); 221 MPFR_ASSERTN (q[0] % 1073741824L == 0L); 222 223 /* another large quotient */ 224 mpfr_set_prec (x, 65); 225 mpfr_set_prec (y, 65); 226 mpfr_const_pi (x, MPFR_RNDN); 227 mpfr_mul_2exp (x, x, 63, MPFR_RNDN); 228 mpfr_const_log2 (y, MPFR_RNDN); 229 mpfr_set_prec (r, 10); 230 mpfr_remquo (r, q, x, y, MPFR_RNDN); 231 /* q should be 41803643793084085130, r should be 605/2048 */ 232 MPFR_ASSERTN (mpfr_cmp_ui_2exp (r, 605, -11) == 0); 233 MPFR_ASSERTN ((q[0] > 0) && ((q[0] % 1073741824L) == 733836170L)); 234 235 /* check cases where quotient is 1.5 +/- eps */ 236 mpfr_set_prec (x, 65); 237 mpfr_set_prec (y, 65); 238 mpfr_set_prec (r, 63); 239 mpfr_set_ui (x, 3, MPFR_RNDN); 240 mpfr_set_ui (y, 2, MPFR_RNDN); 241 mpfr_remquo (r, q, x, y, MPFR_RNDN); 242 /* x/y = 1.5, quotient should be 2 (even rule), remainder should be -1 */ 243 MPFR_ASSERTN (mpfr_cmp_si (r, -1) == 0); 244 MPFR_ASSERTN (q[0] == 2L); 245 mpfr_set_ui (x, 3, MPFR_RNDN); 246 mpfr_nextabove (x); /* 3 + 2^(-63) */ 247 mpfr_set_ui (y, 2, MPFR_RNDN); 248 mpfr_remquo (r, q, x, y, MPFR_RNDN); 249 /* x/y = 1.5 + 2^(-64), quo should be 2, r should be -1 + 2^(-63) */ 250 MPFR_ASSERTN (mpfr_add_ui (r, r, 1, MPFR_RNDN) == 0); 251 MPFR_ASSERTN (mpfr_cmp_ui_2exp (r, 1, -63) == 0); 252 MPFR_ASSERTN (q[0] == 2L); 253 mpfr_set_ui (x, 3, MPFR_RNDN); 254 mpfr_set_ui (y, 2, MPFR_RNDN); 255 mpfr_nextabove (y); /* 2 + 2^(-63) */ 256 mpfr_remquo (r, q, x, y, MPFR_RNDN); 257 /* x/y = 1.5 - eps, quo should be 1, r should be 1 - 2^(-63) */ 258 MPFR_ASSERTN (mpfr_sub_ui (r, r, 1, MPFR_RNDN) == 0); 259 MPFR_ASSERTN (mpfr_cmp_si_2exp (r, -1, -63) == 0); 260 MPFR_ASSERTN (q[0] == 1L); 261 262 /* bug founds by Kaveh Ghazi, 3 May 2007 */ 263 mpfr_set_ui (x, 2, MPFR_RNDN); 264 mpfr_set_ui (y, 3, MPFR_RNDN); 265 mpfr_remainder (r, x, y, MPFR_RNDN); 266 MPFR_ASSERTN (mpfr_cmp_si (r, -1) == 0); 267 268 mpfr_set_si (x, -1, MPFR_RNDN); 269 mpfr_set_ui (y, 1, MPFR_RNDN); 270 mpfr_remainder (r, x, y, MPFR_RNDN); 271 MPFR_ASSERTN (mpfr_cmp_si (r, 0) == 0 && MPFR_SIGN (r) < 0); 272 273 /* check argument reuse */ 274 mpfr_set_si (x, -1, MPFR_RNDN); 275 mpfr_set_ui (y, 1, MPFR_RNDN); 276 mpfr_remainder (x, x, y, MPFR_RNDN); 277 MPFR_ASSERTN (mpfr_cmp_si (x, 0) == 0 && MPFR_SIGN (x) < 0); 278 279 mpfr_set_ui_2exp (x, 1, mpfr_get_emax () - 1, MPFR_RNDN); 280 mpfr_set_ui_2exp (y, 1, mpfr_get_emin (), MPFR_RNDN); 281 mpfr_remquo (r, q, x, y, MPFR_RNDN); 282 MPFR_ASSERTN (mpfr_zero_p (r) && MPFR_SIGN (r) > 0); 283 MPFR_ASSERTN (q[0] == 0); 284 285 mpfr_clear (x); 286 mpfr_clear (y); 287 mpfr_clear (r); 288 289 tests_end_mpfr (); 290 291 return 0; 292} 293