1 #ifndef NUMBER_H 2 #define NUMBER_H 3 4 #ifdef HAVE_CONFIG_H 5 # include "config.h" 6 #else 7 # ifndef HAVE_LIBMPFR 8 # define HAVE_LIBMPFR 9 # endif 10 #endif 11 12 #ifdef HAVE_LIBMPFR 13 14 # if HAVE_LIMITS_H 15 # include <limits.h> /* MPFR uses ULONG_MAX on some systems */ 16 # endif 17 18 # ifndef S_SPLINT_S /* splint barfs on the gmp.h header */ 19 # include <gmp.h> 20 # include <mpfr.h> 21 # endif 22 23 # define Number mpfr_t 24 # define NUM_PREC_MIN MPFR_PREC_MIN 25 # define NUM_PREC_MAX MPFR_PREC_MAX 26 27 # define num_init(n) mpfr_init(n) 28 # define num_init_set(n1, n2) mpfr_init_set((n1), (n2), GMP_RNDN) 29 # define num_init_set_d(n, d) mpfr_init_set_d((n), (d), GMP_RNDN) 30 # define num_init_set_ui(n, ui) mpfr_init_set_ui((n), (ui), GMP_RNDN) 31 # define num_init_set_str(n, str, base) mpfr_init_set_str((n), (str), (base), GMP_RNDN) 32 33 # define num_free(n) mpfr_clear(n) 34 35 typedef mp_exp_t num_exp_t; 36 # define num_get_d(n) mpfr_get_d((n), GMP_RNDN) 37 # define num_get_ui(n) mpfr_get_ui((n), GMP_RNDN) 38 # define num_get_str(a, b, c, d, e) mpfr_get_str(a, b, c, d, e, GMP_RNDN) 39 # define num_free_str(a) mpfr_free_str(a) 40 41 # define num_set(n1, n2) mpfr_set((n1), (n2), GMP_RNDN) 42 # define num_set_d(n, d) mpfr_set_d((n), (d), GMP_RNDN) 43 # define num_set_ui(n, ui) mpfr_set_ui((n), (ui), GMP_RNDN) 44 # define num_set_str(n, str, base) mpfr_set_str((n), (str), (base), GMP_RNDN) 45 # define num_set_nan(n) mpfr_set_nan(n) 46 # define num_set_inf(n, s) mpfr_set_inf(n, s) 47 48 # define num_add(ret, n1, n2) mpfr_add((ret), (n1), (n2), GMP_RNDN) 49 # define num_sub(ret, n1, n2) mpfr_sub((ret), (n1), (n2), GMP_RNDN) 50 # define num_add_ui(ret, n, ui) mpfr_add_ui((ret), (n), (ui), GMP_RNDN) 51 # define num_sub_ui(ret, n, ui) mpfr_sub_ui((ret), (n), (ui), GMP_RNDN) 52 53 # define num_mul(ret, n1, n2) mpfr_mul((ret), (n1), (n2), GMP_RNDN) 54 # define num_mul_si(ret, n, si) mpfr_mul_si((ret), (n), (si), GMP_RNDN) 55 # define num_mul_ui(ret, n, ui) mpfr_mul_ui((ret), (n), (ui), GMP_RNDN) 56 # define num_sqr(ret, n) mpfr_sqr((ret), (n), GMP_RNDN) 57 # define num_sqrt(ret, n) mpfr_sqrt((ret), (n), GMP_RNDN) 58 # define num_cbrt(ret, n) mpfr_cbrt((ret), (n), GMP_RNDN) 59 60 # define num_div(ret, n1, n2) mpfr_div((ret), (n1), (n2), GMP_RNDN) 61 # define num_div_ui(ret, n, ui) mpfr_div_ui((ret), (n), (ui), GMP_RNDN) 62 63 # define num_rint(ret, n) mpfr_rint((ret), (n), GMP_RNDN) 64 # define num_rintz(ret, n) mpfr_rint((ret), (n), GMP_RNDZ) 65 # define num_floor(ret, n) mpfr_floor((ret), (n)) 66 # define num_ceil(ret, n) mpfr_ceil((ret), (n)) 67 68 # define num_pow(ret, n, exp) mpfr_pow((ret), (n), (exp), GMP_RNDN) 69 # define num_pow_si(ret, n, exp) mpfr_pow_si((ret), (n), (exp), GMP_RNDN) 70 # define num_exp(ret, n) mpfr_exp((ret), (n), GMP_RNDN) 71 # define num_factorial(ret, op) mpfr_fac_ui((ret), (op), GMP_RNDN) 72 # define num_log10(ret, n) mpfr_log10((ret), (n), GMP_RNDN) 73 # define num_log2(ret, n) mpfr_log2((ret), (n), GMP_RNDN) 74 # define num_log(ret, n) mpfr_log((ret), (n), GMP_RNDN) 75 76 # define num_sin(ret, n) mpfr_sin((ret), (n), GMP_RNDN) 77 # define num_cos(ret, n) mpfr_cos((ret), (n), GMP_RNDN) 78 # define num_tan(ret, n) mpfr_tan((ret), (n), GMP_RNDN) 79 # define num_asin(ret, n) mpfr_asin((ret), (n), GMP_RNDN) 80 # define num_acos(ret, n) mpfr_acos((ret), (n), GMP_RNDN) 81 # define num_atan(ret, n) mpfr_atan((ret), (n), GMP_RNDN) 82 # define num_sinh(ret, n) mpfr_sinh((ret), (n), GMP_RNDN) 83 # define num_cosh(ret, n) mpfr_cosh((ret), (n), GMP_RNDN) 84 # define num_tanh(ret, n) mpfr_tanh((ret), (n), GMP_RNDN) 85 # define num_asinh(ret, n) mpfr_asinh((ret), (n), GMP_RNDN) 86 # define num_acosh(ret, n) mpfr_acosh((ret), (n), GMP_RNDN) 87 # define num_atanh(ret, n) mpfr_atanh((ret), (n), GMP_RNDN) 88 # define num_acoth(ret, n) do { \ 89 mpfr_pow_si((n), (n), -1, GMP_RNDN); \ 90 mpfr_atanh((ret), (n), GMP_RNDN); \ 91 } while (0) 92 # define num_asech(ret, n) do { \ 93 mpfr_pow_si((n), (n), -1, GMP_RNDN); \ 94 mpfr_acosh((ret), (n), GMP_RNDN); \ 95 } while (0) 96 # define num_acsch(ret, n) do { \ 97 mpfr_pow_si((n), (n), -1, GMP_RNDN); \ 98 mpfr_asinh((ret), (n), GMP_RNDN); \ 99 } while (0) 100 101 extern gmp_randstate_t randstate; 102 # define num_random(n) mpfr_urandomb((n), randstate) 103 # define num_zeta(ret, n) mpfr_zeta((ret), (n), GMP_RNDN); 104 # define num_sinc(ret, n) do { \ 105 if(mpfr_zero_p((n))) { \ 106 mpfr_set_ui((ret), 1, GMP_RNDN); \ 107 } else { \ 108 mpfr_sin((ret), (n), GMP_RNDN); \ 109 mpfr_div((ret), (ret), (n), GMP_RNDN); \ 110 } \ 111 } while (0) 112 113 # define num_const_pi(n) mpfr_const_pi((n), GMP_RNDN) 114 # define num_const_euler(n) mpfr_const_euler((n), GMP_RNDN) 115 116 # define num_out_str(fd, base, n) mpfr_out_str((fd), (base), 0, (n), GMP_RNDN) 117 # define num_is_nan(n) mpfr_nan_p(n) 118 # define num_is_inf(n) mpfr_inf_p(n) 119 # define num_is_zero(n) mpfr_zero_p(n) 120 # define num_is_equal(n1, n2) mpfr_equal_p((n1), (n2)) 121 # define num_is_greater(n1, n2) mpfr_greater_p((n1), (n2)) 122 # define num_is_greaterequal(n1, n2) mpfr_greaterequal_p((n1), (n2)) 123 # define num_is_less(n1, n2) mpfr_less_p((n1), (n2)) 124 # define num_is_lessequal(n1, n2) mpfr_lessequal_p((n1), (n2)) 125 # define num_sign(n) mpfr_sgn(n) 126 # define num_neg(ret, n) mpfr_neg((ret), (n), GMP_RNDN) 127 # define num_abs(ret, n) mpfr_abs((ret), (n), GMP_RNDN) 128 129 # define num_cmp_si(n, si) mpfr_cmp_si((n), (si)) 130 131 typedef mp_prec_t num_prec_t; 132 # define num_get_default_prec() mpfr_get_default_prec() 133 # define num_set_default_prec(prec) mpfr_set_default_prec(prec) 134 135 # ifdef HAVE_MPFR_22 136 # define num_const_catalan(n) mpfr_const_catalan((n), GMP_RNDN) 137 # define num_cot(ret, n) mpfr_cot((ret), (n), GMP_RNDN) 138 # define num_sec(ret, n) mpfr_sec((ret), (n), GMP_RNDN) 139 # define num_csc(ret, n) mpfr_csc((ret), (n), GMP_RNDN) 140 # define num_coth(ret, n) mpfr_coth((ret), (n), GMP_RNDN) 141 # define num_sech(ret, n) mpfr_sech((ret), (n), GMP_RNDN) 142 # define num_csch(ret, n) mpfr_csch((ret), (n), GMP_RNDN) 143 # define num_eint(ret, n) mpfr_eint((ret), (n), GMP_RNDN) 144 # define num_gamma(ret, n) mpfr_gamma((ret), (n), GMP_RNDN) 145 # define num_lngamma(ret, n) mpfr_lngamma((ret), (n), GMP_RNDN) 146 # else // ifdef HAVE_MPFR_22 147 # define num_const_catalan(n) num_set_str((n), W_CATALAN, 0) 148 # define num_cot(ret, n) do { \ 149 mpfr_tan((ret), (n), GMP_RNDN); \ 150 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 151 } while (0) 152 # define num_sec(ret, n) do { \ 153 mpfr_cos((ret), (n), GMP_RNDN); \ 154 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 155 } while (0) 156 # define num_csc(ret, n) do { \ 157 mpfr_sin((ret), (n), GMP_RNDN); \ 158 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 159 } while (0) 160 # define num_coth(ret, n) do { \ 161 mpfr_tanh((ret), (n), GMP_RNDN); \ 162 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 163 } while (0) 164 # define num_sech(ret, n) do { \ 165 mpfr_cosh((ret), (n), GMP_RNDN); \ 166 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 167 } while (0) 168 # define num_csch(ret, n) do { \ 169 mpfr_sinh((ret), (n), GMP_RNDN); \ 170 mpfr_pow_si((ret), (n), -1, GMP_RNDN); \ 171 } while (0) 172 # endif // ifdef HAVE_MPFR_22 173 174 # define num_bor(ret, n1, n2) do { \ 175 mpz_t intfirst, intsecond, intoutput; \ 176 mpz_init(intfirst); \ 177 mpz_init(intsecond); \ 178 mpz_init(intoutput); \ 179 mpfr_get_z(intfirst, (n1), GMP_RNDZ); \ 180 mpfr_get_z(intsecond, (n2), GMP_RNDZ); \ 181 mpz_ior(intoutput, intfirst, intsecond); \ 182 mpfr_set_z((ret), intoutput, GMP_RNDN); \ 183 mpz_clear(intfirst); \ 184 mpz_clear(intsecond); \ 185 mpz_clear(intoutput); \ 186 } while (0) 187 # define num_band(ret, n1, n2) do { \ 188 mpz_t intfirst, intsecond, intoutput; \ 189 mpz_init(intfirst); \ 190 mpz_init(intsecond); \ 191 mpz_init(intoutput); \ 192 mpfr_get_z(intfirst, (n1), GMP_RNDZ); \ 193 mpfr_get_z(intsecond, (n2), GMP_RNDZ); \ 194 mpz_and(intoutput, intfirst, intsecond); \ 195 mpfr_set_z((ret), intoutput, GMP_RNDN); \ 196 mpz_clear(intfirst); \ 197 mpz_clear(intsecond); \ 198 mpz_clear(intoutput); \ 199 } while (0) 200 # define num_bxor(ret, n1, n2) do { \ 201 mpz_t intfirst, intsecond, intoutput; \ 202 mpz_init(intfirst); \ 203 mpz_init(intsecond); \ 204 mpz_init(intoutput); \ 205 mpfr_get_z(intfirst, (n1), GMP_RNDZ); \ 206 mpfr_get_z(intsecond, (n2), GMP_RNDZ); \ 207 mpz_xor(intoutput, intfirst, intsecond); \ 208 mpfr_set_z((ret), intoutput, GMP_RNDN); \ 209 mpz_clear(intfirst); \ 210 mpz_clear(intsecond); \ 211 mpz_clear(intoutput); \ 212 } while (0) 213 # define num_trunc(ret, n) mpfr_trunc((ret), (n)) 214 # define num_comp(ret, n) do { \ 215 mpz_t integer, intoutput; \ 216 mpz_init(integer); \ 217 mpz_init(intoutput); \ 218 mpfr_get_z(integer, (n), GMP_RNDZ); \ 219 mpz_com(intoutput, integer); \ 220 mpfr_set_z((ret), intoutput, GMP_RNDN); \ 221 mpz_clear(integer); \ 222 mpz_clear(intoutput); \ 223 } while (0) 224 /* this is because binary not doesn't make sense when your integer can be 225 * arbitrarily big */ 226 # ifdef _MPFR_H_HAVE_INTMAX_T 227 # define num_bnot(ret, n) mpfr_set_uj((ret), ~mpfr_get_uj((n), GMP_RNDN), GMP_RNDN) 228 # else 229 # define num_bnot(ret, n) mpfr_set_ui((ret), ~mpfr_get_ui((n), GMP_RNDN), GMP_RNDN) 230 # endif 231 232 #else /* this is to reimplement mpfr using doubles */ 233 # include <float.h> /* for DBL_MANT_DIG, according to C89 */ 234 # include <stddef.h> /* for size_t, according to C89 */ 235 # include <stdio.h> /* for FILE* */ 236 struct numberstruct { 237 double value; 238 unsigned int nan; 239 }; 240 typedef struct numberstruct Number[1]; 241 242 # define NUM_PREC_MIN DBL_MANT_DIG 243 # define NUM_PREC_MAX DBL_MANT_DIG 244 245 typedef int num_exp_t; 246 typedef long int num_prec_t; 247 248 void num_init(Number n); 249 void num_init_set(Number n1, 250 const Number n2); 251 void num_init_set_d(Number n, 252 const double d); 253 void num_init_set_ui(Number n, 254 const unsigned int ui); 255 int num_init_set_str(Number n, 256 const char *str, 257 const int base); 258 259 void num_set(Number n1, 260 const Number n2); 261 void num_set_d(Number n, 262 const double d); 263 void num_set_ui(Number n, 264 const unsigned int ui); 265 int num_set_str(Number n, 266 const char *str, 267 const int base); 268 void num_set_nan(Number n); 269 270 void num_mul(Number ret, 271 const Number n1, 272 const Number n2); 273 void num_mul_si(Number ret, 274 const Number n, 275 const int si); 276 void num_mul_ui(Number ret, 277 const Number n, 278 const unsigned int ui); 279 void num_sqr(Number ret, 280 const Number n); 281 void num_sqrt(Number ret, 282 const Number n); 283 void num_cbrt(Number ret, 284 const Number n); 285 286 void num_pow(Number ret, 287 const Number n, 288 const Number exp); 289 void num_pow_si(Number ret, 290 const Number n, 291 const int exp); 292 void num_exp(Number ret, 293 const Number n); 294 void num_factorial(Number ret, 295 unsigned int op); 296 void num_log10(Number ret, 297 const Number n); 298 void num_log2(Number ret, 299 const Number n); 300 void num_log(Number ret, 301 const Number n); 302 303 void num_sin(Number ret, 304 const Number n); 305 void num_cos(Number ret, 306 const Number n); 307 void num_tan(Number ret, 308 const Number n); 309 void num_asin(Number ret, 310 const Number n); 311 void num_acos(Number ret, 312 const Number n); 313 void num_atan(Number ret, 314 const Number n); 315 void num_sinh(Number ret, 316 const Number n); 317 void num_cosh(Number ret, 318 const Number n); 319 void num_tanh(Number ret, 320 const Number n); 321 void num_asinh(Number ret, 322 const Number n); 323 void num_acosh(Number ret, 324 const Number n); 325 void num_atanh(Number ret, 326 const Number n); 327 void num_acoth(Number ret, 328 const Number n); 329 void num_asech(Number ret, 330 const Number n); 331 void num_acsch(Number ret, 332 const Number n); 333 334 void num_add(Number ret, 335 const Number n1, 336 const Number n2); 337 void num_sub(Number ret, 338 const Number n1, 339 const Number n2); 340 void num_add_ui(Number ret, 341 const Number n, 342 const unsigned int ui); 343 void num_sub_ui(Number ret, 344 const Number n, 345 const unsigned int ui); 346 347 void num_div(Number ret, 348 const Number n1, 349 const Number n2); 350 void num_div_ui(Number ret, 351 const Number n, 352 const unsigned int ui); 353 354 void num_rint(Number ret, 355 const Number n); 356 void num_rintz(Number ret, 357 const Number n); 358 void num_floor(Number ret, 359 const Number n); 360 void num_ceil(Number ret, 361 const Number n); 362 363 int num_random(Number n); 364 # define num_zeta(ret, n) num_unimplemented() 365 void num_sinc(Number ret, 366 const Number n); 367 368 void num_const_pi(Number n); 369 void num_const_euler(Number n); 370 void num_const_catalan(Number n); 371 372 void num_out_str(FILE *fd, 373 const int base, 374 const Number n); 375 int num_is_nan(const Number n); 376 int num_is_inf(const Number n); 377 int num_is_zero(const Number n); 378 int num_is_equal(const Number n1, 379 const Number n2); 380 int num_is_greater(const Number n1, 381 const Number n2); 382 int num_is_greaterequal(const Number n1, 383 const Number n2); 384 int num_is_less(const Number n1, 385 const Number n2); 386 int num_is_lessequal(const Number n1, 387 const Number n2); 388 int num_sign(const Number n); 389 void num_neg(Number ret, 390 const Number n); 391 void num_abs(Number ret, 392 const Number n); 393 int num_cmp_si(const Number n, 394 const int si); 395 396 void num_free(Number); 397 double num_get_d(const Number); 398 unsigned int num_get_ui(const Number); 399 char *num_get_str(char *str, 400 num_exp_t *expptr, 401 const int b, 402 const size_t n, 403 const Number op); 404 # define num_free_str(a) free(a); 405 num_prec_t num_get_default_prec(); 406 # define num_set_default_prec(prec) /* num_set_default_prec doesn't exist */ 407 void num_cot(Number ret, 408 const Number n); 409 void num_sec(Number ret, 410 const Number n); 411 void num_csc(Number ret, 412 const Number n); 413 void num_coth(Number ret, 414 const Number n); 415 void num_sech(Number ret, 416 const Number n); 417 void num_csch(Number ret, 418 const Number n); 419 # define num_lngamma(ret, n) num_unimplemented() 420 # define num_gamma(ret, n) num_unimplemented() 421 # define num_eint(ret, n) num_unimplemented() 422 void num_comp(Number ret, 423 const Number n); 424 void num_bnot(Number ret, 425 const Number n); 426 void num_band(Number ret, 427 const Number n1, 428 const Number n2); 429 void num_bxor(Number ret, 430 const Number n1, 431 const Number n2); 432 void num_bor(Number ret, 433 const Number n1, 434 const Number n2); 435 #endif // ifdef HAVE_LIBMPFR 436 437 void init_numbers(void); 438 int is_int(const Number); 439 #endif /* NUMBER_H */ 440 /* vim:set expandtab: */ 441