granicus.if.org Git - postgresql/blob - src/backend/utils/adt/float.c

   1 /*-------------------------------------------------------------------------
   2  *
   3  * float.c
   4  *        Functions for the built-in floating-point types.
   5  *
   6  * Portions Copyright (c) 1996-2000, PostgreSQL, Inc
   7  * Portions Copyright (c) 1994, Regents of the University of California
   8  *
   9  *
  10  * IDENTIFICATION
  11  *        $Header: /cvsroot/pgsql/src/backend/utils/adt/float.c,v 1.61 2000/06/14 18:17:42 petere Exp $
  12  *
  13  *-------------------------------------------------------------------------
  14  */
  15 /*
  16  * OLD COMMENTS
  17  *              Basic float4 ops:
  18  *               float4in, float4out, float4abs, float4um
  19  *              Basic float8 ops:
  20  *               float8in, float8inAd, float8out, float8outAd, float8abs, float8um
  21  *              Arithmetic operators:
  22  *               float4pl, float4mi, float4mul, float4div
  23  *               float8pl, float8mi, float8mul, float8div
  24  *              Comparison operators:
  25  *               float4eq, float4ne, float4lt, float4le, float4gt, float4ge
  26  *               float8eq, float8ne, float8lt, float8le, float8gt, float8ge
  27  *              Conversion routines:
  28  *               ftod, dtof, i4tod, dtoi4, i2tod, dtoi2, itof, ftoi, i2tof, ftoi2
  29  *
  30  *              Random float8 ops:
  31  *               dround, dtrunc, dsqrt, dcbrt, dpow, dexp, dlog1
  32  *              Arithmetic operators:
  33  *               float48pl, float48mi, float48mul, float48div
  34  *               float84pl, float84mi, float84mul, float84div
  35  *              Comparison operators:
  36  *               float48eq, float48ne, float48lt, float48le, float48gt, float48ge
  37  *               float84eq, float84ne, float84lt, float84le, float84gt, float84ge
  38  *
  39  *              (You can do the arithmetic and comparison stuff using conversion
  40  *               routines, but then you pay the overhead of converting...)
  41  *
  42  * XXX GLUESOME STUFF. FIX IT! -AY '94
  43  *
  44  *              Added some additional conversion routines and cleaned up
  45  *               a bit of the existing code. Need to change the error checking
  46  *               for calls to pow(), exp() since on some machines (my Linux box
  47  *               included) these routines do not set errno. - tgl 97/05/10
  48  */
  49 #include <ctype.h>
  50 #include <errno.h>
  51 #include <float.h>                              /* faked on sunos4 */
  52 #include <math.h>
  53
  54 #include "postgres.h"
  55
  56 #ifdef HAVE_LIMITS_H
  57 #include <limits.h>
  58 #ifndef MAXINT
  59 #define MAXINT            INT_MAX
  60 #endif
  61 #else
  62 #ifdef HAVE_VALUES_H
  63 #include <values.h>
  64 #endif
  65 #endif
  66
  67 /* for finite() on Solaris */
  68 #ifdef HAVE_IEEEFP_H
  69 # include <ieeefp.h>
  70 #endif
  71
  72 #include "fmgr.h"
  73 #include "utils/builtins.h"
  74
  75 static void CheckFloat8Val(double val);
  76
  77 #ifndef NAN
  78 #define NAN             (0.0/0.0)
  79 #endif
  80
  81 #ifndef SHRT_MAX
  82 #define SHRT_MAX 32767
  83 #endif
  84 #ifndef SHRT_MIN
  85 #define SHRT_MIN (-32768)
  86 #endif
  87
  88 #define FORMAT                  'g'             /* use "g" output format as standard
  89                                                                  * format */
  90 /* not sure what the following should be, but better to make it over-sufficient */
  91 #define MAXFLOATWIDTH   64
  92 #define MAXDOUBLEWIDTH  128
  93
  94 #if !(NeXT && NX_CURRENT_COMPILER_RELEASE > NX_COMPILER_RELEASE_3_2)
  95  /* NS3.3 has conflicting declarations of these in <math.h> */
  96
  97 #ifndef atof
  98 extern double atof(const char *p);
  99
 100 #endif
 101
 102 #ifndef HAVE_CBRT
 103 #define cbrt my_cbrt
 104 static double cbrt(double x);
 105
 106 #else
 107 #if !defined(nextstep)
 108 extern double cbrt(double x);
 109
 110 #endif
 111 #endif
 112
 113 #ifndef HAVE_RINT
 114 #define rint my_rint
 115 static double rint(double x);
 116
 117 #else
 118 extern double rint(double x);
 119
 120 #endif
 121
 122 #endif
 123
 124 /* ========== USER I/O ROUTINES ========== */
 125
 126
 127 #define FLOAT4_MAX               FLT_MAX
 128 #define FLOAT4_MIN               FLT_MIN
 129 #define FLOAT8_MAX               DBL_MAX
 130 #define FLOAT8_MIN               DBL_MIN
 131
 132
 133 /*
 134    check to see if a float4 val is outside of
 135    the FLOAT4_MIN, FLOAT4_MAX bounds.
 136
 137    raise an elog warning if it is
 138 */
 139 static void
 140 CheckFloat4Val(double val)
 141 {
 142
 143         /*
 144          * defining unsafe floats's will make float4 and float8 ops faster at
 145          * the cost of safety, of course!
 146          */
 147 #ifdef UNSAFE_FLOATS
 148         return;
 149 #else
 150         if (fabs(val) > FLOAT4_MAX)
 151                 elog(ERROR, "Bad float4 input format -- overflow");
 152         if (val != 0.0 && fabs(val) < FLOAT4_MIN)
 153                 elog(ERROR, "Bad float4 input format -- underflow");
 154         return;
 155 #endif   /* UNSAFE_FLOATS */
 156 }
 157
 158 /*
 159    check to see if a float8 val is outside of
 160    the FLOAT8_MIN, FLOAT8_MAX bounds.
 161
 162    raise an elog warning if it is
 163 */
 164 static void
 165 CheckFloat8Val(double val)
 166 {
 167
 168         /*
 169          * defining unsafe floats's will make float4 and float8 ops faster at
 170          * the cost of safety, of course!
 171          */
 172 #ifdef UNSAFE_FLOATS
 173         return;
 174 #else
 175         if (fabs(val) > FLOAT8_MAX)
 176                 elog(ERROR, "Bad float8 input format -- overflow");
 177         if (val != 0.0 && fabs(val) < FLOAT8_MIN)
 178                 elog(ERROR, "Bad float8 input format -- underflow");
 179         return;
 180 #endif   /* UNSAFE_FLOATS */
 181 }
 182
 183 /*
 184  *              float4in                - converts "num" to float
 185  *                                                restricted syntax:
 186  *                                                {<sp>} [+|-] {digit} [.{digit}] [<exp>]
 187  *                                                where <sp> is a space, digit is 0-9,
 188  *                                                <exp> is "e" or "E" followed by an integer.
 189  */
 190 float32
 191 float4in(char *num)
 192 {
 193         float32         result = (float32) palloc(sizeof(float32data));
 194         double          val;
 195         char       *endptr;
 196
 197         errno = 0;
 198         val = strtod(num, &endptr);
 199         if (*endptr != '\0')
 200         {
 201                 /* Should we accept "NaN" or "Infinity" for float4? */
 202                 elog(ERROR, "Bad float4 input format '%s'", num);
 203         }
 204         else
 205         {
 206                 if (errno == ERANGE)
 207                         elog(ERROR, "Input '%s' is out of range for float4", num);
 208         }
 209
 210         /*
 211          * if we get here, we have a legal double, still need to check to see
 212          * if it's a legal float
 213          */
 214
 215         CheckFloat4Val(val);
 216
 217         *result = val;
 218         return result;
 219 }
 220
 221 /*
 222  *              float4out               - converts a float4 number to a string
 223  *                                                using a standard output format
 224  */
 225 char *
 226 float4out(float32 num)
 227 {
 228         char       *ascii = (char *) palloc(MAXFLOATWIDTH + 1);
 229
 230         if (!num)
 231                 return strcpy(ascii, "(null)");
 232
 233         sprintf(ascii, "%.*g", FLT_DIG, *num);
 234         return ascii;
 235 }
 236
 237
 238 /*
 239  *              float8in                - converts "num" to float8
 240  *                                                restricted syntax:
 241  *                                                {<sp>} [+|-] {digit} [.{digit}] [<exp>]
 242  *                                                where <sp> is a space, digit is 0-9,
 243  *                                                <exp> is "e" or "E" followed by an integer.
 244  */
 245 float64
 246 float8in(char *num)
 247 {
 248         float64         result = (float64) palloc(sizeof(float64data));
 249         double          val;
 250         char       *endptr;
 251
 252         errno = 0;
 253         val = strtod(num, &endptr);
 254         if (*endptr != '\0')
 255         {
 256                 if (strcasecmp(num, "NaN") == 0)
 257                         val = NAN;
 258                 else if (strcasecmp(num, "Infinity") == 0)
 259                         val = HUGE_VAL;
 260                 else
 261                         elog(ERROR, "Bad float8 input format '%s'", num);
 262         }
 263         else
 264         {
 265                 if (errno == ERANGE)
 266                         elog(ERROR, "Input '%s' is out of range for float8", num);
 267         }
 268
 269         CheckFloat8Val(val);
 270
 271         *result = val;
 272         return result;
 273 }
 274
 275
 276 /*
 277  *              float8out               - converts float8 number to a string
 278  *                                                using a standard output format
 279  */
 280 char *
 281 float8out(float64 num)
 282 {
 283         char       *ascii = (char *) palloc(MAXDOUBLEWIDTH + 1);
 284
 285         if (!num)
 286                 return strcpy(ascii, "(null)");
 287
 288         if (isnan(*num))
 289                 return strcpy(ascii, "NaN");
 290         if (isinf(*num))
 291                 return strcpy(ascii, "Infinity");
 292
 293         sprintf(ascii, "%.*g", DBL_DIG, *num);
 294         return ascii;
 295 }
 296
 297 /* ========== PUBLIC ROUTINES ========== */
 298
 299
 300 /*
 301  *              ======================
 302  *              FLOAT4 BASE OPERATIONS
 303  *              ======================
 304  */
 305
 306 /*
 307  *              float4abs               - returns a pointer to |arg1| (absolute value)
 308  */
 309 float32
 310 float4abs(float32 arg1)
 311 {
 312         float32         result;
 313         double          val;
 314
 315         if (!arg1)
 316                 return (float32) NULL;
 317
 318         val = fabs(*arg1);
 319
 320         CheckFloat4Val(val);
 321
 322         result = (float32) palloc(sizeof(float32data));
 323         *result = val;
 324         return result;
 325 }
 326
 327 /*
 328  *              float4um                - returns a pointer to -arg1 (unary minus)
 329  */
 330 float32
 331 float4um(float32 arg1)
 332 {
 333         float32         result;
 334         double          val;
 335
 336         if (!arg1)
 337                 return (float32) NULL;
 338
 339         val = ((*arg1 != 0) ? -(*arg1) : *arg1);
 340         CheckFloat4Val(val);
 341
 342         result = (float32) palloc(sizeof(float32data));
 343         *result = val;
 344         return result;
 345 }
 346
 347 float32
 348 float4larger(float32 arg1, float32 arg2)
 349 {
 350         float32         result;
 351
 352         if (!arg1 || !arg2)
 353                 return (float32) NULL;
 354
 355         result = (float32) palloc(sizeof(float32data));
 356
 357         *result = ((*arg1 > *arg2) ? *arg1 : *arg2);
 358         return result;
 359 }
 360
 361 float32
 362 float4smaller(float32 arg1, float32 arg2)
 363 {
 364         float32         result;
 365
 366         if (!arg1 || !arg2)
 367                 return (float32) NULL;
 368
 369         result = (float32) palloc(sizeof(float32data));
 370
 371         *result = ((*arg1 > *arg2) ? *arg2 : *arg1);
 372         return result;
 373 }
 374
 375 /*
 376  *              ======================
 377  *              FLOAT8 BASE OPERATIONS
 378  *              ======================
 379  */
 380
 381 /*
 382  *              float8abs               - returns a pointer to |arg1| (absolute value)
 383  */
 384 float64
 385 float8abs(float64 arg1)
 386 {
 387         float64         result;
 388         double          val;
 389
 390         if (!arg1)
 391                 return (float64) NULL;
 392
 393         result = (float64) palloc(sizeof(float64data));
 394
 395         val = fabs(*arg1);
 396         CheckFloat8Val(val);
 397         *result = val;
 398         return result;
 399 }
 400
 401
 402 /*
 403  *              float8um                - returns a pointer to -arg1 (unary minus)
 404  */
 405 float64
 406 float8um(float64 arg1)
 407 {
 408         float64         result;
 409         double          val;
 410
 411         if (!arg1)
 412                 return (float64) NULL;
 413
 414         val = ((*arg1 != 0) ? -(*arg1) : *arg1);
 415
 416         CheckFloat8Val(val);
 417         result = (float64) palloc(sizeof(float64data));
 418         *result = val;
 419         return result;
 420 }
 421
 422 float64
 423 float8larger(float64 arg1, float64 arg2)
 424 {
 425         float64         result;
 426
 427         if (!arg1 || !arg2)
 428                 return (float64) NULL;
 429
 430         result = (float64) palloc(sizeof(float64data));
 431
 432         *result = ((*arg1 > *arg2) ? *arg1 : *arg2);
 433         return result;
 434 }
 435
 436 float64
 437 float8smaller(float64 arg1, float64 arg2)
 438 {
 439         float64         result;
 440
 441         if (!arg1 || !arg2)
 442                 return (float64) NULL;
 443
 444         result = (float64) palloc(sizeof(float64data));
 445
 446         *result = ((*arg1 > *arg2) ? *arg2 : *arg1);
 447         return result;
 448 }
 449
 450
 451 /*
 452  *              ====================
 453  *              ARITHMETIC OPERATORS
 454  *              ====================
 455  */
 456
 457 /*
 458  *              float4pl                - returns a pointer to arg1 + arg2
 459  *              float4mi                - returns a pointer to arg1 - arg2
 460  *              float4mul               - returns a pointer to arg1 * arg2
 461  *              float4div               - returns a pointer to arg1 / arg2
 462  *              float4inc               - returns a poniter to arg1 + 1.0
 463  */
 464 float32
 465 float4pl(float32 arg1, float32 arg2)
 466 {
 467         float32         result;
 468         double          val;
 469
 470         if (!arg1 || !arg2)
 471                 return (float32) NULL;
 472
 473         val = *arg1 + *arg2;
 474         CheckFloat4Val(val);
 475
 476         result = (float32) palloc(sizeof(float32data));
 477         *result = val;
 478
 479         return result;
 480 }
 481
 482 float32
 483 float4mi(float32 arg1, float32 arg2)
 484 {
 485         float32         result;
 486         double          val;
 487
 488         if (!arg1 || !arg2)
 489                 return (float32) NULL;
 490
 491         val = *arg1 - *arg2;
 492
 493         CheckFloat4Val(val);
 494         result = (float32) palloc(sizeof(float32data));
 495         *result = val;
 496         return result;
 497 }
 498
 499 float32
 500 float4mul(float32 arg1, float32 arg2)
 501 {
 502         float32         result;
 503         double          val;
 504
 505         if (!arg1 || !arg2)
 506                 return (float32) NULL;
 507
 508         val = *arg1 * *arg2;
 509
 510         CheckFloat4Val(val);
 511         result = (float32) palloc(sizeof(float32data));
 512         *result = val;
 513         return result;
 514 }
 515
 516 float32
 517 float4div(float32 arg1, float32 arg2)
 518 {
 519         float32         result;
 520         double          val;
 521
 522         if (!arg1 || !arg2)
 523                 return (float32) NULL;
 524
 525         if (*arg2 == 0.0)
 526                 elog(ERROR, "float4div: divide by zero error");
 527
 528         val = *arg1 / *arg2;
 529
 530         CheckFloat4Val(val);
 531         result = (float32) palloc(sizeof(float32data));
 532         *result = val;
 533         return result;
 534 }
 535
 536 float32
 537 float4inc(float32 arg1)
 538 {
 539         float32         result;
 540         double          val;
 541
 542         if (!arg1)
 543                 return (float32) NULL;
 544
 545         val = *arg1 + (float32data) 1.0;
 546
 547         CheckFloat4Val(val);
 548         result = (float32) palloc(sizeof(float32data));
 549         *result = val;
 550         return result;
 551 }
 552
 553 /*
 554  *              float8pl                - returns a pointer to arg1 + arg2
 555  *              float8mi                - returns a pointer to arg1 - arg2
 556  *              float8mul               - returns a pointer to arg1 * arg2
 557  *              float8div               - returns a pointer to arg1 / arg2
 558  *              float8inc               - returns a pointer to arg1 + 1.0
 559  */
 560 float64
 561 float8pl(float64 arg1, float64 arg2)
 562 {
 563         float64         result;
 564         double          val;
 565
 566         if (!arg1 || !arg2)
 567                 return (float64) NULL;
 568
 569         result = (float64) palloc(sizeof(float64data));
 570
 571         val = *arg1 + *arg2;
 572         CheckFloat8Val(val);
 573         *result = val;
 574         return result;
 575 }
 576
 577 float64
 578 float8mi(float64 arg1, float64 arg2)
 579 {
 580         float64         result;
 581         double          val;
 582
 583         if (!arg1 || !arg2)
 584                 return (float64) NULL;
 585
 586         result = (float64) palloc(sizeof(float64data));
 587
 588         val = *arg1 - *arg2;
 589         CheckFloat8Val(val);
 590         *result = val;
 591         return result;
 592 }
 593
 594 float64
 595 float8mul(float64 arg1, float64 arg2)
 596 {
 597         float64         result;
 598         double          val;
 599
 600         if (!arg1 || !arg2)
 601                 return (float64) NULL;
 602
 603         result = (float64) palloc(sizeof(float64data));
 604
 605         val = *arg1 * *arg2;
 606         CheckFloat8Val(val);
 607         *result = val;
 608         return result;
 609 }
 610
 611 float64
 612 float8div(float64 arg1, float64 arg2)
 613 {
 614         float64         result;
 615         double          val;
 616
 617         if (!arg1 || !arg2)
 618                 return (float64) NULL;
 619
 620         result = (float64) palloc(sizeof(float64data));
 621
 622         if (*arg2 == 0.0)
 623                 elog(ERROR, "float8div: divide by zero error");
 624
 625         val = *arg1 / *arg2;
 626         CheckFloat8Val(val);
 627         *result = val;
 628         return result;
 629 }
 630
 631 float64
 632 float8inc(float64 arg1)
 633 {
 634         float64         result;
 635         double          val;
 636
 637         if (!arg1)
 638                 return (float64) NULL;
 639
 640         val = *arg1 + (float64data) 1.0;
 641         CheckFloat8Val(val);
 642         result = (float64) palloc(sizeof(float64data));
 643         *result = val;
 644         return result;
 645 }
 646
 647
 648 /*
 649  *              ====================
 650  *              COMPARISON OPERATORS
 651  *              ====================
 652  */
 653
 654 /*
 655  *              float4{eq,ne,lt,le,gt,ge}               - float4/float4 comparison operations
 656  */
 657 bool
 658 float4eq(float32 arg1, float32 arg2)
 659 {
 660         if (!arg1 || !arg2)
 661                 return 0;
 662
 663         return *arg1 == *arg2;
 664 }
 665
 666 bool
 667 float4ne(float32 arg1, float32 arg2)
 668 {
 669         if (!arg1 || !arg2)
 670                 return 0;
 671
 672         return *arg1 != *arg2;
 673 }
 674
 675 bool
 676 float4lt(float32 arg1, float32 arg2)
 677 {
 678         if (!arg1 || !arg2)
 679                 return 0;
 680
 681         return *arg1 < *arg2;
 682 }
 683
 684 bool
 685 float4le(float32 arg1, float32 arg2)
 686 {
 687         if (!arg1 || !arg2)
 688                 return 0;
 689
 690         return *arg1 <= *arg2;
 691 }
 692
 693 bool
 694 float4gt(float32 arg1, float32 arg2)
 695 {
 696         if (!arg1 || !arg2)
 697                 return 0;
 698
 699         return *arg1 > *arg2;
 700 }
 701
 702 bool
 703 float4ge(float32 arg1, float32 arg2)
 704 {
 705         if (!arg1 || !arg2)
 706                 return 0;
 707
 708         return *arg1 >= *arg2;
 709 }
 710
 711 /*
 712  *              float8{eq,ne,lt,le,gt,ge}               - float8/float8 comparison operations
 713  */
 714 bool
 715 float8eq(float64 arg1, float64 arg2)
 716 {
 717         if (!arg1 || !arg2)
 718                 return 0;
 719
 720         return *arg1 == *arg2;
 721 }
 722
 723 bool
 724 float8ne(float64 arg1, float64 arg2)
 725 {
 726         if (!arg1 || !arg2)
 727                 return 0;
 728
 729         return *arg1 != *arg2;
 730 }
 731
 732 bool
 733 float8lt(float64 arg1, float64 arg2)
 734 {
 735         if (!arg1 || !arg2)
 736                 return 0;
 737
 738         return *arg1 < *arg2;
 739 }
 740
 741 bool
 742 float8le(float64 arg1, float64 arg2)
 743 {
 744         if (!arg1 || !arg2)
 745                 return 0;
 746
 747         return *arg1 <= *arg2;
 748 }
 749
 750 bool
 751 float8gt(float64 arg1, float64 arg2)
 752 {
 753         if (!arg1 || !arg2)
 754                 return 0;
 755
 756         return *arg1 > *arg2;
 757 }
 758
 759 bool
 760 float8ge(float64 arg1, float64 arg2)
 761 {
 762         if (!arg1 || !arg2)
 763                 return 0;
 764
 765         return *arg1 >= *arg2;
 766 }
 767
 768
 769 /*
 770  *              ===================
 771  *              CONVERSION ROUTINES
 772  *              ===================
 773  */
 774
 775 /*
 776  *              ftod                    - converts a float4 number to a float8 number
 777  */
 778 float64
 779 ftod(float32 num)
 780 {
 781         float64         result;
 782
 783         if (!num)
 784                 return (float64) NULL;
 785
 786         result = (float64) palloc(sizeof(float64data));
 787
 788         *result = *num;
 789         return result;
 790 }
 791
 792
 793 /*
 794  *              dtof                    - converts a float8 number to a float4 number
 795  */
 796 float32
 797 dtof(float64 num)
 798 {
 799         float32         result;
 800
 801         if (!num)
 802                 return (float32) NULL;
 803
 804         CheckFloat4Val(*num);
 805
 806         result = (float32) palloc(sizeof(float32data));
 807
 808         *result = *num;
 809         return result;
 810 }
 811
 812
 813 /*
 814  *              dtoi4                   - converts a float8 number to an int4 number
 815  */
 816 int32
 817 dtoi4(float64 num)
 818 {
 819         int32           result;
 820
 821         if (!num)
 822                 return 0;                               /* fmgr will return NULL anyway */
 823
 824         if ((*num < INT_MIN) || (*num > INT_MAX))
 825                 elog(ERROR, "dtoi4: integer out of range");
 826
 827         result = rint(*num);
 828         return result;
 829 }
 830
 831
 832 /*
 833  *              dtoi2                   - converts a float8 number to an int2 number
 834  */
 835 Datum
 836 dtoi2(PG_FUNCTION_ARGS)
 837 {
 838         float8          num = PG_GETARG_FLOAT8(0);
 839         int16           result;
 840
 841         if ((num < SHRT_MIN) || (num > SHRT_MAX))
 842                 elog(ERROR, "dtoi2: integer out of range");
 843
 844         result = (int16) rint(num);
 845         PG_RETURN_INT16(result);
 846 }
 847
 848
 849 /*
 850  *              i4tod                   - converts an int4 number to a float8 number
 851  */
 852 Datum
 853 i4tod(PG_FUNCTION_ARGS)
 854 {
 855         int32           num = PG_GETARG_INT32(0);
 856         float8          result;
 857
 858         result = num;
 859         PG_RETURN_FLOAT8(result);
 860 }
 861
 862
 863 /*
 864  *              i2tod                   - converts an int2 number to a float8 number
 865  */
 866 Datum
 867 i2tod(PG_FUNCTION_ARGS)
 868 {
 869         int16           num = PG_GETARG_INT16(0);
 870         float8          result;
 871
 872         result = num;
 873         PG_RETURN_FLOAT8(result);
 874 }
 875
 876
 877 /*
 878  *              ftoi4                   - converts a float8 number to an int4 number
 879  */
 880 int32
 881 ftoi4(float32 num)
 882 {
 883         int32           result;
 884
 885         if (!num)
 886                 return 0;                               /* fmgr will return NULL anyway */
 887
 888         if ((*num < INT_MIN) || (*num > INT_MAX))
 889                 elog(ERROR, "ftoi4: integer out of range");
 890
 891         result = rint(*num);
 892         return result;
 893 }
 894
 895
 896 /*
 897  *              ftoi2                   - converts a float4 number to an int2 number
 898  */
 899 Datum
 900 ftoi2(PG_FUNCTION_ARGS)
 901 {
 902         float4          num = PG_GETARG_FLOAT4(0);
 903         int16           result;
 904
 905         if ((num < SHRT_MIN) || (num > SHRT_MAX))
 906                 elog(ERROR, "ftoi2: integer out of range");
 907
 908         result = (int16) rint(num);
 909         PG_RETURN_INT16(result);
 910 }
 911
 912
 913 /*
 914  *              i4tof                   - converts an int4 number to a float8 number
 915  */
 916 Datum
 917 i4tof(PG_FUNCTION_ARGS)
 918 {
 919         int32           num = PG_GETARG_INT32(0);
 920         float4          result;
 921
 922         result = num;
 923         PG_RETURN_FLOAT4(result);
 924 }
 925
 926
 927 /*
 928  *              i2tof                   - converts an int2 number to a float4 number
 929  */
 930 Datum
 931 i2tof(PG_FUNCTION_ARGS)
 932 {
 933         int16           num = PG_GETARG_INT16(0);
 934         float4          result;
 935
 936         result = num;
 937         PG_RETURN_FLOAT4(result);
 938 }
 939
 940
 941 /*
 942  *              float8_text             - converts a float8 number to a text string
 943  */
 944 text *
 945 float8_text(float64 num)
 946 {
 947         text       *result;
 948         int                     len;
 949         char       *str;
 950
 951         str = float8out(num);
 952         len = (strlen(str) + VARHDRSZ);
 953
 954         result = palloc(len);
 955
 956         VARSIZE(result) = len;
 957         memmove(VARDATA(result), str, (len - VARHDRSZ));
 958
 959         pfree(str);
 960         return result;
 961 }       /* float8_text() */
 962
 963
 964 /*
 965  *              text_float8             - converts a text string to a float8 number
 966  */
 967 float64
 968 text_float8(text *string)
 969 {
 970         float64         result;
 971         int                     len;
 972         char       *str;
 973
 974         len = (VARSIZE(string) - VARHDRSZ);
 975         str = palloc(len + 1);
 976         memmove(str, VARDATA(string), len);
 977         *(str + len) = '\0';
 978
 979         result = float8in(str);
 980         pfree(str);
 981
 982         return result;
 983 }       /* text_float8() */
 984
 985
 986 /*
 987  *              float4_text             - converts a float4 number to a text string
 988  */
 989 text *
 990 float4_text(float32 num)
 991 {
 992         text       *result;
 993         int                     len;
 994         char       *str;
 995
 996         str = float4out(num);
 997         len = (strlen(str) + VARHDRSZ);
 998
 999         result = palloc(len);
1000
1001         VARSIZE(result) = len;
1002         memmove(VARDATA(result), str, (len - VARHDRSZ));
1003
1004         pfree(str);
1005         return result;
1006 }       /* float4_text() */
1007
1008
1009 /*
1010  *              text_float4             - converts a text string to a float4 number
1011  */
1012 float32
1013 text_float4(text *string)
1014 {
1015         float32         result;
1016         int                     len;
1017         char       *str;
1018
1019         len = (VARSIZE(string) - VARHDRSZ);
1020         str = palloc(len + 1);
1021         memmove(str, VARDATA(string), len);
1022         *(str + len) = '\0';
1023
1024         result = float4in(str);
1025         pfree(str);
1026
1027         return result;
1028 }       /* text_float4() */
1029
1030
1031 /*
1032  *              =======================
1033  *              RANDOM FLOAT8 OPERATORS
1034  *              =======================
1035  */
1036
1037 /*
1038  *              dround                  - returns a pointer to  ROUND(arg1)
1039  */
1040 float64
1041 dround(float64 arg1)
1042 {
1043         float64         result;
1044         double          tmp;
1045
1046         if (!arg1)
1047                 return (float64) NULL;
1048
1049         result = (float64) palloc(sizeof(float64data));
1050
1051         tmp = *arg1;
1052         *result = (float64data) rint(tmp);
1053         return result;
1054 }
1055
1056
1057 /*
1058  *              dtrunc                  - returns a pointer to  truncation of arg1,
1059  *                                                arg1 >= 0 ... the greatest integer as float8 less
1060  *                                                                              than or equal to arg1
1061  *                                                arg1 < 0      ... the greatest integer as float8 greater
1062  *                                                                              than or equal to arg1
1063  */
1064 float64
1065 dtrunc(float64 arg1)
1066 {
1067         float64         result;
1068         double          tmp;
1069
1070         if (!arg1)
1071                 return (float64) NULL;
1072
1073         result = (float64) palloc(sizeof(float64data));
1074
1075         tmp = *arg1;
1076         if (*arg1 >= 0)
1077                 *result = (float64data) floor(tmp);
1078         else
1079                 *result = (float64data) -(floor(-tmp));
1080         return result;
1081 }
1082
1083
1084 /*
1085  *              dsqrt                   - returns a pointer to square root of arg1
1086  */
1087 float64
1088 dsqrt(float64 arg1)
1089 {
1090         float64         result;
1091         double          tmp;
1092
1093         if (!arg1)
1094                 return (float64) NULL;
1095
1096         result = (float64) palloc(sizeof(float64data));
1097
1098         tmp = *arg1;
1099         *result = (float64data) sqrt(tmp);
1100         return result;
1101 }
1102
1103
1104 /*
1105  *              dcbrt                   - returns a pointer to cube root of arg1
1106  */
1107 float64
1108 dcbrt(float64 arg1)
1109 {
1110         float64         result;
1111         double          tmp;
1112
1113         if (!arg1)
1114                 return (float64) NULL;
1115
1116         result = (float64) palloc(sizeof(float64data));
1117
1118         tmp = *arg1;
1119         *result = (float64data) cbrt(tmp);
1120         return result;
1121 }
1122
1123
1124 /*
1125  *              dpow                    - returns a pointer to pow(arg1,arg2)
1126  */
1127 float64
1128 dpow(float64 arg1, float64 arg2)
1129 {
1130         float64         result;
1131         double          tmp1,
1132                                 tmp2;
1133
1134         if (!arg1 || !arg2)
1135                 return (float64) NULL;
1136
1137         result = (float64) palloc(sizeof(float64data));
1138
1139         tmp1 = *arg1;
1140         tmp2 = *arg2;
1141
1142         /*
1143          * We must check both for errno getting set and for a NaN result, in
1144          * order to deal with the vagaries of different platforms...
1145          */
1146         errno = 0;
1147         *result = (float64data) pow(tmp1, tmp2);
1148         if (errno != 0
1149 #ifdef HAVE_FINITE
1150                 || !finite(*result)
1151 #endif
1152                 )
1153                 elog(ERROR, "pow() result is out of range");
1154
1155         CheckFloat8Val(*result);
1156         return result;
1157 }
1158
1159
1160 /*
1161  *              dexp                    - returns a pointer to the exponential function of arg1
1162  */
1163 float64
1164 dexp(float64 arg1)
1165 {
1166         float64         result;
1167         double          tmp;
1168
1169         if (!arg1)
1170                 return (float64) NULL;
1171
1172         result = (float64) palloc(sizeof(float64data));
1173
1174         tmp = *arg1;
1175
1176         /*
1177          * We must check both for errno getting set and for a NaN result, in
1178          * order to deal with the vagaries of different platforms. Also, a
1179          * zero result implies unreported underflow.
1180          */
1181         errno = 0;
1182         *result = (float64data) exp(tmp);
1183         if (errno != 0 || *result == 0.0
1184 #ifdef HAVE_FINITE
1185                 || !finite(*result)
1186 #endif
1187                 )
1188                 elog(ERROR, "exp() result is out of range");
1189
1190         CheckFloat8Val(*result);
1191         return result;
1192 }
1193
1194
1195 /*
1196  *              dlog1                   - returns a pointer to the natural logarithm of arg1
1197  *                                                ("dlog" is already a logging routine...)
1198  */
1199 float64
1200 dlog1(float64 arg1)
1201 {
1202         float64         result;
1203         double          tmp;
1204
1205         if (!PointerIsValid(arg1))
1206                 return (float64) NULL;
1207
1208         result = (float64) palloc(sizeof(float64data));
1209
1210         tmp = *arg1;
1211         if (tmp == 0.0)
1212                 elog(ERROR, "can't take log of zero");
1213         if (tmp < 0)
1214                 elog(ERROR, "can't take log of a negative number");
1215         *result = (float64data) log(tmp);
1216
1217         CheckFloat8Val(*result);
1218         return result;
1219 }       /* dlog1() */
1220
1221
1222 /*
1223  *              dlog10                  - returns a pointer to the base 10 logarithm of arg1
1224  */
1225 float64
1226 dlog10(float64 arg1)
1227 {
1228         float64         result;
1229         double          tmp;
1230
1231         if (!PointerIsValid(arg1))
1232                 return (float64) NULL;
1233
1234         result = (float64) palloc(sizeof(float64data));
1235
1236         tmp = *arg1;
1237         if (tmp == 0.0)
1238                 elog(ERROR, "can't take log of zero");
1239         if (tmp < 0)
1240                 elog(ERROR, "can't take log of a negative number");
1241         *result = (float64data) log10(tmp);
1242
1243         CheckFloat8Val(*result);
1244         return result;
1245 }       /* dlog10() */
1246
1247
1248 /*
1249  *              dacos                   - returns a pointer to the arccos of arg1 (radians)
1250  */
1251 float64
1252 dacos(float64 arg1)
1253 {
1254         float64         result;
1255         double          tmp;
1256
1257         if (!PointerIsValid(arg1))
1258                 return (float64) NULL;
1259
1260         result = (float64) palloc(sizeof(float64data));
1261
1262         tmp = *arg1;
1263         errno = 0;
1264         *result = (float64data) acos(tmp);
1265         if (errno != 0
1266 #ifdef HAVE_FINITE
1267                 || !finite(*result)
1268 #endif
1269                 )
1270                 elog(ERROR, "dacos(%f) input is out of range", *arg1);
1271
1272         CheckFloat8Val(*result);
1273         return result;
1274 }       /* dacos() */
1275
1276
1277 /*
1278  *              dasin                   - returns a pointer to the arcsin of arg1 (radians)
1279  */
1280 float64
1281 dasin(float64 arg1)
1282 {
1283         float64         result;
1284         double          tmp;
1285
1286         if (!PointerIsValid(arg1))
1287                 return (float64) NULL;
1288
1289         result = (float64) palloc(sizeof(float64data));
1290
1291         tmp = *arg1;
1292         errno = 0;
1293         *result = (float64data) asin(tmp);
1294         if (errno != 0
1295 #ifdef HAVE_FINITE
1296                 || !finite(*result)
1297 #endif
1298                 )
1299                 elog(ERROR, "dasin(%f) input is out of range", *arg1);
1300
1301         CheckFloat8Val(*result);
1302         return result;
1303 }       /* dasin() */
1304
1305
1306 /*
1307  *              datan                   - returns a pointer to the arctan of arg1 (radians)
1308  */
1309 float64
1310 datan(float64 arg1)
1311 {
1312         float64         result;
1313         double          tmp;
1314
1315         if (!PointerIsValid(arg1))
1316                 return (float64) NULL;
1317
1318         result = (float64) palloc(sizeof(float64data));
1319
1320         tmp = *arg1;
1321         errno = 0;
1322         *result = (float64data) atan(tmp);
1323         if (errno != 0
1324 #ifdef HAVE_FINITE
1325                 || !finite(*result)
1326 #endif
1327                 )
1328                 elog(ERROR, "atan(%f) input is out of range", *arg1);
1329
1330         CheckFloat8Val(*result);
1331         return result;
1332 }       /* datan() */
1333
1334
1335 /*
1336  *              atan2                   - returns a pointer to the arctan2 of arg1 (radians)
1337  */
1338 float64
1339 datan2(float64 arg1, float64 arg2)
1340 {
1341         float64         result;
1342
1343         if (!PointerIsValid(arg1) || !PointerIsValid(arg1))
1344                 return (float64) NULL;
1345
1346         result = (float64) palloc(sizeof(float64data));
1347
1348         errno = 0;
1349         *result = (float64data) atan2(*arg1, *arg2);
1350         if (errno != 0
1351 #ifdef HAVE_FINITE
1352                 || !finite(*result)
1353 #endif
1354                 )
1355                 elog(ERROR, "atan2(%f,%f) input is out of range", *arg1, *arg2);
1356
1357         CheckFloat8Val(*result);
1358         return result;
1359 }       /* datan2() */
1360
1361
1362 /*
1363  *              dcos                    - returns a pointer to the cosine of arg1 (radians)
1364  */
1365 float64
1366 dcos(float64 arg1)
1367 {
1368         float64         result;
1369         double          tmp;
1370
1371         if (!PointerIsValid(arg1))
1372                 return (float64) NULL;
1373
1374         result = (float64) palloc(sizeof(float64data));
1375
1376         tmp = *arg1;
1377         errno = 0;
1378         *result = (float64data) cos(tmp);
1379         if (errno != 0
1380 #ifdef HAVE_FINITE
1381                 || !finite(*result)
1382 #endif
1383                 )
1384                 elog(ERROR, "dcos(%f) input is out of range", *arg1);
1385
1386         CheckFloat8Val(*result);
1387         return result;
1388 }       /* dcos() */
1389
1390
1391 /*
1392  *              dcot                    - returns a pointer to the cotangent of arg1 (radians)
1393  */
1394 float64
1395 dcot(float64 arg1)
1396 {
1397         float64         result;
1398         double          tmp;
1399
1400         if (!PointerIsValid(arg1))
1401                 return (float64) NULL;
1402
1403         result = (float64) palloc(sizeof(float64data));
1404
1405         tmp = *arg1;
1406         errno = 0;
1407         *result = (float64data) tan(tmp);
1408         if ((errno != 0) || (*result == 0.0)
1409 #ifdef HAVE_FINITE
1410                 || !finite(*result)
1411 #endif
1412                 )
1413                 elog(ERROR, "dcot(%f) input is out of range", *arg1);
1414
1415         *result = 1.0 / (*result);
1416         CheckFloat8Val(*result);
1417         return result;
1418 }       /* dcot() */
1419
1420
1421 /*
1422  *              dsin                    - returns a pointer to the sine of arg1 (radians)
1423  */
1424 float64
1425 dsin(float64 arg1)
1426 {
1427         float64         result;
1428         double          tmp;
1429
1430         if (!PointerIsValid(arg1))
1431                 return (float64) NULL;
1432
1433         result = (float64) palloc(sizeof(float64data));
1434
1435         tmp = *arg1;
1436         errno = 0;
1437         *result = (float64data) sin(tmp);
1438         if (errno != 0
1439 #ifdef HAVE_FINITE
1440                 || !finite(*result)
1441 #endif
1442                 )
1443                 elog(ERROR, "dsin(%f) input is out of range", *arg1);
1444
1445         CheckFloat8Val(*result);
1446         return result;
1447 }       /* dsin() */
1448
1449
1450 /*
1451  *              dtan                    - returns a pointer to the tangent of arg1 (radians)
1452  */
1453 float64
1454 dtan(float64 arg1)
1455 {
1456         float64         result;
1457         double          tmp;
1458
1459         if (!PointerIsValid(arg1))
1460                 return (float64) NULL;
1461
1462         result = (float64) palloc(sizeof(float64data));
1463
1464         tmp = *arg1;
1465         errno = 0;
1466         *result = (float64data) tan(tmp);
1467         if (errno != 0
1468 #ifdef HAVE_FINITE
1469                 || !finite(*result)
1470 #endif
1471                 )
1472                 elog(ERROR, "dtan(%f) input is out of range", *arg1);
1473
1474         CheckFloat8Val(*result);
1475         return result;
1476 }       /* dtan() */
1477
1478
1479 #ifndef M_PI
1480 /* from my RH5.2 gcc math.h file - thomas 2000-04-03 */
1481 #define M_PI 3.14159265358979323846
1482 #endif
1483
1484
1485 /*
1486  *              degrees         - returns a pointer to degrees converted from radians
1487  */
1488 float64
1489 degrees(float64 arg1)
1490 {
1491         float64         result;
1492
1493         if (!arg1)
1494                 return (float64) NULL;
1495
1496         result = (float64) palloc(sizeof(float64data));
1497
1498         *result = ((*arg1) * (180.0 / M_PI));
1499
1500         CheckFloat8Val(*result);
1501         return result;
1502 }       /* degrees() */
1503
1504
1505 /*
1506  *              dpi                             - returns a pointer to degrees converted to radians
1507  */
1508 float64
1509 dpi(void)
1510 {
1511         float64         result;
1512
1513         result = (float64) palloc(sizeof(float64data));
1514
1515         *result = (M_PI);
1516
1517         return result;
1518 }       /* dpi() */
1519
1520
1521 /*
1522  *              radians         - returns a pointer to radians converted from degrees
1523  */
1524 float64
1525 radians(float64 arg1)
1526 {
1527         float64         result;
1528
1529         if (!arg1)
1530                 return (float64) NULL;
1531
1532         result = (float64) palloc(sizeof(float64data));
1533
1534         *result = ((*arg1) * (M_PI / 180.0));
1535
1536         CheckFloat8Val(*result);
1537         return result;
1538 }       /* radians() */
1539
1540
1541 /*
1542  *              drandom         - returns a random number
1543  */
1544 float64
1545 drandom(void)
1546 {
1547         float64         result;
1548
1549         result = (float64) palloc(sizeof(float64data));
1550
1551         /* result 0.0-1.0 */
1552         *result = (((double) random()) / RAND_MAX);
1553
1554         CheckFloat8Val(*result);
1555         return result;
1556 }       /* drandom() */
1557
1558
1559 /*
1560  *              setseed         - set seed for the random number generator
1561  */
1562 int32
1563 setseed(float64 seed)
1564 {
1565         int                     iseed = ((*seed) * RAND_MAX);
1566
1567         srandom((unsigned int) ((*seed) * RAND_MAX));
1568
1569         return iseed;
1570 }       /* setseed() */
1571
1572
1573 /*
1574  *              ====================
1575  *              ARITHMETIC OPERATORS
1576  *              ====================
1577  */
1578
1579 /*
1580  *              float48pl               - returns a pointer to arg1 + arg2
1581  *              float48mi               - returns a pointer to arg1 - arg2
1582  *              float48mul              - returns a pointer to arg1 * arg2
1583  *              float48div              - returns a pointer to arg1 / arg2
1584  */
1585 float64
1586 float48pl(float32 arg1, float64 arg2)
1587 {
1588         float64         result;
1589
1590         if (!arg1 || !arg2)
1591                 return (float64) NULL;
1592
1593         result = (float64) palloc(sizeof(float64data));
1594
1595         *result = *arg1 + *arg2;
1596         CheckFloat8Val(*result);
1597         return result;
1598 }
1599
1600 float64
1601 float48mi(float32 arg1, float64 arg2)
1602 {
1603         float64         result;
1604
1605         if (!arg1 || !arg2)
1606                 return (float64) NULL;
1607
1608         result = (float64) palloc(sizeof(float64data));
1609
1610         *result = *arg1 - *arg2;
1611         CheckFloat8Val(*result);
1612         return result;
1613 }
1614
1615 float64
1616 float48mul(float32 arg1, float64 arg2)
1617 {
1618         float64         result;
1619
1620         if (!arg1 || !arg2)
1621                 return (float64) NULL;
1622
1623         result = (float64) palloc(sizeof(float64data));
1624
1625         *result = *arg1 * *arg2;
1626         CheckFloat8Val(*result);
1627         return result;
1628 }
1629
1630 float64
1631 float48div(float32 arg1, float64 arg2)
1632 {
1633         float64         result;
1634
1635         if (!arg1 || !arg2)
1636                 return (float64) NULL;
1637
1638         result = (float64) palloc(sizeof(float64data));
1639
1640         if (*arg2 == 0.0)
1641                 elog(ERROR, "float48div: divide by zero");
1642
1643         *result = *arg1 / *arg2;
1644         CheckFloat8Val(*result);
1645         return result;
1646 }
1647
1648 /*
1649  *              float84pl               - returns a pointer to arg1 + arg2
1650  *              float84mi               - returns a pointer to arg1 - arg2
1651  *              float84mul              - returns a pointer to arg1 * arg2
1652  *              float84div              - returns a pointer to arg1 / arg2
1653  */
1654 float64
1655 float84pl(float64 arg1, float32 arg2)
1656 {
1657         float64         result;
1658
1659         if (!arg1 || !arg2)
1660                 return (float64) NULL;
1661
1662         result = (float64) palloc(sizeof(float64data));
1663
1664         *result = *arg1 + *arg2;
1665         CheckFloat8Val(*result);
1666         return result;
1667 }
1668
1669 float64
1670 float84mi(float64 arg1, float32 arg2)
1671 {
1672         float64         result;
1673
1674         if (!arg1 || !arg2)
1675                 return (float64) NULL;
1676
1677         result = (float64) palloc(sizeof(float64data));
1678
1679         *result = *arg1 - *arg2;
1680         CheckFloat8Val(*result);
1681         return result;
1682 }
1683
1684 float64
1685 float84mul(float64 arg1, float32 arg2)
1686 {
1687
1688         float64         result;
1689
1690         if (!arg1 || !arg2)
1691                 return (float64) NULL;
1692
1693         result = (float64) palloc(sizeof(float64data));
1694
1695         *result = *arg1 * *arg2;
1696         CheckFloat8Val(*result);
1697         return result;
1698 }
1699
1700 float64
1701 float84div(float64 arg1, float32 arg2)
1702 {
1703         float64         result;
1704
1705         if (!arg1 || !arg2)
1706                 return (float64) NULL;
1707
1708         result = (float64) palloc(sizeof(float64data));
1709
1710         if (*arg2 == 0.0)
1711                 elog(ERROR, "float48div: divide by zero");
1712
1713         *result = *arg1 / *arg2;
1714         CheckFloat8Val(*result);
1715         return result;
1716 }
1717
1718 /*
1719  *              ====================
1720  *              COMPARISON OPERATORS
1721  *              ====================
1722  */
1723
1724 /*
1725  *              float48{eq,ne,lt,le,gt,ge}              - float4/float8 comparison operations
1726  */
1727 bool
1728 float48eq(float32 arg1, float64 arg2)
1729 {
1730         if (!arg1 || !arg2)
1731                 return 0;
1732
1733         return *arg1 == *arg2;
1734 }
1735
1736 bool
1737 float48ne(float32 arg1, float64 arg2)
1738 {
1739         if (!arg1 || !arg2)
1740                 return 0;
1741
1742         return *arg1 != *arg2;
1743 }
1744
1745 bool
1746 float48lt(float32 arg1, float64 arg2)
1747 {
1748         if (!arg1 || !arg2)
1749                 return 0;
1750
1751         return *arg1 < *arg2;
1752 }
1753
1754 bool
1755 float48le(float32 arg1, float64 arg2)
1756 {
1757         if (!arg1 || !arg2)
1758                 return 0;
1759
1760         return *arg1 <= *arg2;
1761 }
1762
1763 bool
1764 float48gt(float32 arg1, float64 arg2)
1765 {
1766         if (!arg1 || !arg2)
1767                 return 0;
1768
1769         return *arg1 > *arg2;
1770 }
1771
1772 bool
1773 float48ge(float32 arg1, float64 arg2)
1774 {
1775         if (!arg1 || !arg2)
1776                 return 0;
1777
1778         return *arg1 >= *arg2;
1779 }
1780
1781 /*
1782  *              float84{eq,ne,lt,le,gt,ge}              - float4/float8 comparison operations
1783  */
1784 bool
1785 float84eq(float64 arg1, float32 arg2)
1786 {
1787         if (!arg1 || !arg2)
1788                 return 0;
1789
1790         return *arg1 == *arg2;
1791 }
1792
1793 bool
1794 float84ne(float64 arg1, float32 arg2)
1795 {
1796         if (!arg1 || !arg2)
1797                 return 0;
1798
1799         return *arg1 != *arg2;
1800 }
1801
1802 bool
1803 float84lt(float64 arg1, float32 arg2)
1804 {
1805         if (!arg1 || !arg2)
1806                 return 0;
1807
1808         return *arg1 < *arg2;
1809 }
1810
1811 bool
1812 float84le(float64 arg1, float32 arg2)
1813 {
1814         if (!arg1 || !arg2)
1815                 return 0;
1816
1817         return *arg1 <= *arg2;
1818 }
1819
1820 bool
1821 float84gt(float64 arg1, float32 arg2)
1822 {
1823         if (!arg1 || !arg2)
1824                 return 0;
1825
1826         return *arg1 > *arg2;
1827 }
1828
1829 bool
1830 float84ge(float64 arg1, float32 arg2)
1831 {
1832         if (!arg1 || !arg2)
1833                 return 0;
1834
1835         return *arg1 >= *arg2;
1836 }
1837
1838 /* ========== PRIVATE ROUTINES ========== */
1839
1840 /* From "fdlibm" @ netlib.att.com */
1841
1842 #ifndef HAVE_RINT
1843
1844 /* @(#)s_rint.c 5.1 93/09/24 */
1845 /*
1846  * ====================================================
1847  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
1848  *
1849  * Developed at SunPro, a Sun Microsystems, Inc. business.
1850  * Permission to use, copy, modify, and distribute this
1851  * software is freely granted, provided that this notice
1852  * is preserved.
1853  * ====================================================
1854  */
1855
1856 /*
1857  * rint(x)
1858  * Return x rounded to integral value according to the prevailing
1859  * rounding mode.
1860  * Method:
1861  *              Using floating addition.
1862  * Exception:
1863  *              Inexact flag raised if x not equal to rint(x).
1864  */
1865
1866 #ifdef __STDC__
1867 static const double
1868 #else
1869 static double
1870 #endif
1871                         one = 1.0,
1872                         TWO52[2] = {
1873         4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
1874         -4.50359962737049600000e+15,/* 0xC3300000, 0x00000000 */
1875 };
1876
1877 #ifdef __STDC__
1878 static double
1879 rint(double x)
1880 #else
1881 static double
1882 rint(x)
1883 double          x;
1884
1885 #endif
1886 {
1887         int                     i0,
1888                                 n0,
1889                                 j0,
1890                                 sx;
1891         unsigned        i,
1892                                 i1;
1893         double          w,
1894                                 t;
1895
1896         n0 = (*((int *) &one) >> 29) ^ 1;
1897         i0 = *(n0 + (int *) &x);
1898         sx = (i0 >> 31) & 1;
1899         i1 = *(1 - n0 + (int *) &x);
1900         j0 = ((i0 >> 20) & 0x7ff) - 0x3ff;
1901         if (j0 < 20)
1902         {
1903                 if (j0 < 0)
1904                 {
1905                         if (((i0 & 0x7fffffff) | i1) == 0)
1906                                 return x;
1907                         i1 |= (i0 & 0x0fffff);
1908                         i0 &= 0xfffe0000;
1909                         i0 |= ((i1 | -i1) >> 12) & 0x80000;
1910                         *(n0 + (int *) &x) = i0;
1911                         w = TWO52[sx] + x;
1912                         t = w - TWO52[sx];
1913                         i0 = *(n0 + (int *) &t);
1914                         *(n0 + (int *) &t) = (i0 & 0x7fffffff) | (sx << 31);
1915                         return t;
1916                 }
1917                 else
1918                 {
1919                         i = (0x000fffff) >> j0;
1920                         if (((i0 & i) | i1) == 0)
1921                                 return x;               /* x is integral */
1922                         i >>= 1;
1923                         if (((i0 & i) | i1) != 0)
1924                         {
1925                                 if (j0 == 19)
1926                                         i1 = 0x40000000;
1927                                 else
1928                                         i0 = (i0 & (~i)) | ((0x20000) >> j0);
1929                         }
1930                 }
1931         }
1932         else if (j0 > 51)
1933         {
1934                 if (j0 == 0x400)
1935                         return x + x;           /* inf or NaN */
1936                 else
1937                         return x;                       /* x is integral */
1938         }
1939         else
1940         {
1941                 i = ((unsigned) (0xffffffff)) >> (j0 - 20);
1942                 if ((i1 & i) == 0)
1943                         return x;                       /* x is integral */
1944                 i >>= 1;
1945                 if ((i1 & i) != 0)
1946                         i1 = (i1 & (~i)) | ((0x40000000) >> (j0 - 20));
1947         }
1948         *(n0 + (int *) &x) = i0;
1949         *(1 - n0 + (int *) &x) = i1;
1950         w = TWO52[sx] + x;
1951         return w - TWO52[sx];
1952 }
1953
1954 #endif   /* !HAVE_RINT */
1955
1956 #ifndef HAVE_CBRT
1957
1958 static
1959 double
1960 cbrt(x)
1961 double          x;
1962 {
1963         int                     isneg = (x < 0.0);
1964         double          tmpres = pow(fabs(x), (double) 1.0 / (double) 3.0);
1965
1966         return isneg ? -tmpres : tmpres;
1967 }
1968
1969 #endif   /* !HAVE_CBRT */