1/* tfmod -- test file for mpfr_fmod 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 28#if MPFR_VERSION >= MPFR_VERSION_NUM(2,4,0) 29 30#define TEST_FUNCTION mpfr_fmod 31#define TWO_ARGS 32#include "tgeneric.c" 33 34/* compute remainder as in definition: 35 r = x - n * y, where n = trunc(x/y). 36 warning: may change flags. */ 37static int 38slow_fmod (mpfr_ptr r, mpfr_srcptr x, mpfr_srcptr y, mpfr_rnd_t rnd) 39{ 40 mpfr_t q; 41 int inexact; 42 if (MPFR_UNLIKELY (MPFR_IS_SINGULAR (x) || MPFR_IS_SINGULAR (y))) 43 { 44 if (MPFR_IS_NAN (x) || MPFR_IS_NAN (y) || MPFR_IS_INF (x) 45 || MPFR_IS_ZERO (y)) 46 { 47 MPFR_SET_NAN (r); 48 MPFR_RET_NAN; 49 } 50 else /* either y is Inf and x is 0 or non-special, 51 or x is 0 and y is non-special, 52 in both cases the quotient is zero. */ 53 return mpfr_set (r, x, rnd); 54 } 55 /* regular cases */ 56 /* if 2^(ex-1) <= |x| < 2^ex, and 2^(ey-1) <= |y| < 2^ey, 57 then |x/y| < 2^(ex-ey+1) */ 58 mpfr_init2 (q, 59 MAX (MPFR_PREC_MIN, mpfr_get_exp (x) - mpfr_get_exp (y) + 1)); 60 mpfr_div (q, x, y, MPFR_RNDZ); 61 mpfr_trunc (q, q); /* may change inexact flag */ 62 mpfr_prec_round (q, mpfr_get_prec (q) + mpfr_get_prec (y), MPFR_RNDZ); 63 inexact = mpfr_mul (q, q, y, MPFR_RNDZ); /* exact */ 64 inexact = mpfr_sub (r, x, q, rnd); 65 mpfr_clear (q); 66 return inexact; 67} 68 69static void 70test_failed (mpfr_t erem, mpfr_t grem, int eret, int gret, mpfr_t x, mpfr_t y, 71 mpfr_rnd_t rnd) 72{ 73 printf ("error: mpfr_fmod (r, x, y, rnd)\n x = "); 74 mpfr_out_str (stdout, 10, 0, x, MPFR_RNDD); 75 printf ("\n y = "); 76 mpfr_out_str (stdout, 10, 0, y, MPFR_RNDD); 77 printf ("\nrnd = %s", mpfr_print_rnd_mode (rnd)); 78 if (eret != gret) 79 printf ("\nexpected %s return value, got %d", 80 (eret < 0 ? "negative" : eret > 0 ? "positive" : "zero"), gret); 81 printf ("\n expected r = "); 82 mpfr_out_str (stdout, 10, 0, erem, MPFR_RNDD); 83 printf ("\n got r = "); 84 mpfr_out_str (stdout, 10, 0, grem, MPFR_RNDD); 85 putchar ('\n'); 86 87 exit (1); 88} 89 90static void 91check (mpfr_t r0, mpfr_t x, mpfr_t y, mpfr_rnd_t rnd) 92{ 93 int inex0, inex1; 94 mpfr_t r1; 95 mpfr_init2 (r1, mpfr_get_prec (r0)); 96 97 inex0 = mpfr_fmod (r0, x, y, rnd); 98 inex1 = slow_fmod (r1, x, y, rnd); 99 if (!mpfr_equal_p (r0, r1) || inex0 != inex1) 100 test_failed (r1, r0, inex1, inex0, x, y, rnd); 101 mpfr_clear (r1); 102} 103 104static void 105regular (void) 106{ 107 mpfr_t x, y, r; 108 mpfr_inits (x, y, r, (mpfr_ptr) 0); 109 110 /* remainder = 0 */ 111 mpfr_set_str (y, "FEDCBA987654321p-64", 16, MPFR_RNDN); 112 mpfr_pow_ui (x, y, 42, MPFR_RNDN); 113 check (r, x, y, MPFR_RNDN); 114 115 /* x < y */ 116 mpfr_set_ui_2exp (x, 64723, -19, MPFR_RNDN); 117 mpfr_mul (x, x, y, MPFR_RNDN); 118 check (r, x, y, MPFR_RNDN); 119 120 /* sign(x) = sign (r) */ 121 mpfr_set_ui (x, 123798, MPFR_RNDN); 122 mpfr_set_ui (y, 10, MPFR_RNDN); 123 check (r, x, y, MPFR_RNDN); 124 125 /* huge difference between precisions */ 126 mpfr_set_prec (x, 314); 127 mpfr_set_prec (y, 8); 128 mpfr_set_prec (r, 123); 129 mpfr_const_pi (x, MPFR_RNDD); /* x = pi */ 130 mpfr_set_ui_2exp (y, 1, 3, MPFR_RNDD); /* y = 1/8 */ 131 check (r, x, y, MPFR_RNDD); 132 133 mpfr_clears (x, y, r, (mpfr_ptr) 0); 134} 135 136static void 137special (void) 138{ 139 int inexact; 140 mpfr_t x, y, r, nan; 141 mpfr_inits (x, y, r, nan, (mpfr_ptr) 0); 142 143 mpfr_set_nan (nan); 144 145 /* fmod (NaN, NaN) is NaN */ 146 mpfr_set_nan (x); 147 mpfr_set_nan (y); 148 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 149 if (!mpfr_nan_p (r) || inexact != 0) 150 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 151 152 /* fmod (NaN, +0) is NaN */ 153 mpfr_set_ui (y, 0, MPFR_RNDN); 154 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 155 if (!mpfr_nan_p (r) || inexact != 0) 156 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 157 158 /* fmod (+1, 0) is NaN */ 159 mpfr_set_ui (x, 1, MPFR_RNDN); 160 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 161 if (!mpfr_nan_p (r) || inexact != 0) 162 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 163 164 /* fmod (0, 0) is NaN */ 165 mpfr_set_ui (x, 0, MPFR_RNDN); 166 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 167 if (!mpfr_nan_p (r) || inexact != 0) 168 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 169 170 /* fmod (+inf, +0) is NaN */ 171 mpfr_set_inf (x, +1); 172 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 173 if (!mpfr_nan_p (r) || inexact != 0) 174 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 175 176 /* fmod (-inf, +0) is NaN */ 177 mpfr_set_inf (x, -1); 178 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 179 if (!mpfr_nan_p (r) || inexact != 0) 180 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 181 182 /* fmod (-inf, -0) is NaN */ 183 mpfr_neg (x, x, MPFR_RNDN); 184 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 185 if (!mpfr_nan_p (r) || inexact != 0) 186 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 187 188 /* fmod (-inf, +1) is NaN */ 189 mpfr_set_ui (y, +1, MPFR_RNDN); 190 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 191 if (!mpfr_nan_p (r) || inexact != 0) 192 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 193 194 /* fmod (+inf, +1) is NaN */ 195 mpfr_neg (x, x, MPFR_RNDN); 196 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 197 if (!mpfr_nan_p (r) || inexact != 0) 198 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 199 200 /* fmod (+inf, -inf) is NaN */ 201 mpfr_set_inf (y, -1); 202 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 203 if (!mpfr_nan_p (r) || inexact != 0) 204 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 205 206 /* fmod (-inf, -inf) is NaN */ 207 mpfr_neg (x, x, MPFR_RNDN); 208 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 209 if (!mpfr_nan_p (r) || inexact != 0) 210 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 211 212 /* fmod (-inf, +inf) is NaN */ 213 mpfr_neg (y, y, MPFR_RNDN); 214 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 215 if (!mpfr_nan_p (r) || inexact != 0) 216 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 217 218 /* fmod (+inf, +inf) is NaN */ 219 mpfr_neg (x, x, MPFR_RNDN); 220 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 221 if (!mpfr_nan_p (r) || inexact != 0) 222 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 223 224 /* fmod (x, +inf) = x, if x is finite */ 225 mpfr_set_ui (x, 1, MPFR_RNDN); 226 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 227 if (!mpfr_equal_p (r, x) || inexact != 0) 228 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 229 mpfr_neg (x, x, MPFR_RNDN); 230 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 231 if (!mpfr_equal_p (r, x) || inexact != 0) 232 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 233 234 /* fmod (+0, +inf) = +0 */ 235 mpfr_set_ui (x, 0, MPFR_RNDN); 236 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 237 if (!mpfr_equal_p (r, x) || inexact != 0) 238 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 239 240 /* fmod (-0, +inf) = -0 */ 241 mpfr_neg (x, x, MPFR_RNDN); 242 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 243 if (!mpfr_equal_p (r, x) || inexact != 0) 244 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 245 246 /* fmod (x, -inf) = x, if x is finite */ 247 mpfr_set_inf (y, -1); 248 mpfr_set_ui (x, 1, MPFR_RNDN); 249 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 250 if (!mpfr_equal_p (r, x) || inexact != 0) 251 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 252 mpfr_neg (x, x, MPFR_RNDN); 253 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 254 if (!mpfr_equal_p (r, x) || inexact != 0) 255 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 256 257 /* fmod (+0, -inf) = +0 */ 258 mpfr_set_ui (x, 0, MPFR_RNDN); 259 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 260 if (!mpfr_equal_p (r, x) || inexact != 0) 261 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 262 263 /* fmod (-0, -inf) = -0 */ 264 mpfr_neg (x, x, MPFR_RNDN); 265 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 266 if (!mpfr_equal_p (r, x) || inexact != 0) 267 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 268 269 /* fmod (+0, +0) is NaN */ 270 mpfr_set_ui (x, 0, MPFR_RNDN); 271 mpfr_set_ui (y, 0, MPFR_RNDN); 272 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 273 if (!mpfr_nan_p (r) || inexact != 0) 274 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 275 276 /* fmod (+0, -0) is NaN */ 277 mpfr_neg (y, y, MPFR_RNDN); 278 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 279 if (!mpfr_nan_p (r) || inexact != 0) 280 test_failed (r, nan, 0, inexact, x, y, MPFR_RNDN); 281 282 /* fmod (+0, +1) = +0 */ 283 mpfr_set_ui (y, 1, MPFR_RNDN); 284 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 285 if (!mpfr_equal_p (r, x) || inexact != 0) 286 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 287 288 /* fmod (+0, -1) = +0 */ 289 mpfr_neg (y, y, MPFR_RNDN); 290 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 291 if (!mpfr_equal_p (r, x) || inexact != 0) 292 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 293 294 /* fmod (-0, -1) = -0 */ 295 mpfr_neg (x, x, MPFR_RNDN); 296 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 297 if (!mpfr_equal_p (r, x) || inexact != 0) 298 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 299 300 /* fmod (-0, +1) = -0 */ 301 mpfr_neg (y, y, MPFR_RNDN); 302 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 303 if (!mpfr_equal_p (r, x) || inexact != 0) 304 test_failed (r, x, 0, inexact, x, y, MPFR_RNDN); 305 306 mpfr_clears (x, y, r, nan, (mpfr_ptr) 0); 307 return; 308} 309 310/* bug reported by Eric Veach */ 311static void 312bug20090519 (void) 313{ 314 mpfr_t x, y, r; 315 int inexact; 316 mpfr_inits2 (100, x, y, r, (mpfr_ptr) 0); 317 318 mpfr_set_prec (x, 3); 319 mpfr_set_prec (y, 3); 320 mpfr_set_prec (r, 3); 321 mpfr_set_si (x, 8, MPFR_RNDN); 322 mpfr_set_si (y, 7, MPFR_RNDN); 323 check (r, x, y, MPFR_RNDN); 324 325 mpfr_set_prec (x, 10); 326 mpfr_set_prec (y, 10); 327 mpfr_set_prec (r, 10); 328 mpfr_set_ui_2exp (x, 3, 26, MPFR_RNDN); 329 mpfr_set_si (y, (1 << 9) - 1, MPFR_RNDN); 330 check (r, x, y, MPFR_RNDN); 331 332 mpfr_set_prec (x, 100); 333 mpfr_set_prec (y, 100); 334 mpfr_set_prec (r, 100); 335 mpfr_set_str (x, "3.5", 10, MPFR_RNDN); 336 mpfr_set_str (y, "1.1", 10, MPFR_RNDN); 337 check (r, x, y, MPFR_RNDN); 338 /* double check, with a pre-computed value */ 339 { 340 mpfr_t er; 341 mpfr_init2 (er, 100); 342 mpfr_set_str (er, "CCCCCCCCCCCCCCCCCCCCCCCC8p-102", 16, MPFR_RNDN); 343 344 inexact = mpfr_fmod (r, x, y, MPFR_RNDN); 345 if (!mpfr_equal_p (r, er) || inexact != 0) 346 test_failed (er, r, 0, inexact, x, y, MPFR_RNDN); 347 348 mpfr_clear (er); 349 } 350 351 mpfr_set_si (x, 20, MPFR_RNDN); 352 mpfr_set_ui_2exp (y, 1, 1, MPFR_RNDN); /* exact */ 353 mpfr_sin (y, y, MPFR_RNDN); 354 check (r, x, y, MPFR_RNDN); 355 356 mpfr_clears(r, x, y, (mpfr_ptr) 0); 357} 358 359int 360main (int argc, char *argv[]) 361{ 362 tests_start_mpfr (); 363 364 bug20090519 (); 365 366 test_generic (2, 100, 100); 367 368 special (); 369 regular (); 370 371 tests_end_mpfr (); 372 return 0; 373} 374 375#else 376 377int 378main (void) 379{ 380 printf ("Warning! Test disabled for this MPFR version.\n"); 381 return 0; 382} 383 384#endif 385