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

   1 /* ----------
   2  * numeric.c
   3  *
   4  *      An exact numeric data type for the Postgres database system
   5  *
   6  *      1998 Jan Wieck
   7  *
   8  * $Header: /cvsroot/pgsql/src/backend/utils/adt/numeric.c,v 1.30 2000/06/14 18:17:44 petere Exp $
   9  *
  10  * ----------
  11  */
  12
  13 #include "postgres.h"
  14
  15 #include <ctype.h>
  16 #include <float.h>
  17 #include <math.h>
  18 #include <errno.h>
  19 #include <sys/types.h>
  20
  21 #include "utils/builtins.h"
  22 #include "utils/numeric.h"
  23
  24 /* ----------
  25  * Uncomment the following to enable compilation of dump_numeric()
  26  * and dump_var() and to get a dump of any result produced by make_result().
  27  * ----------
  28 #define NUMERIC_DEBUG
  29  */
  30
  31
  32 /* ----------
  33  * Local definitions
  34  * ----------
  35  */
  36 #ifndef MIN
  37 #define MIN(a,b) (((a)<(b)) ? (a) : (b))
  38 #endif
  39 #ifndef MAX
  40 #define MAX(a,b) (((a)>(b)) ? (a) : (b))
  41 #endif
  42
  43 #ifndef NAN
  44 #define NAN             (0.0/0.0)
  45 #endif
  46
  47
  48 /* ----------
  49  * Local data types
  50  *
  51  * Note: the first digit of a NumericVar's value is assumed to be multiplied
  52  * by 10 ** weight.  Another way to say it is that there are weight+1 digits
  53  * before the decimal point.  It is possible to have weight < 0.
  54  *
  55  * The value represented by a NumericVar is determined by the sign, weight,
  56  * ndigits, and digits[] array.  The rscale and dscale are carried along,
  57  * but they are just auxiliary information until rounding is done before
  58  * final storage or display.  (Scales are the number of digits wanted
  59  * *after* the decimal point.  Scales are always >= 0.)
  60  *
  61  * buf points at the physical start of the palloc'd digit buffer for the
  62  * NumericVar.  digits points at the first digit in actual use (the one
  63  * with the specified weight).  We normally leave an unused byte or two
  64  * (preset to zeroes) between buf and digits, so that there is room to store
  65  * a carry out of the top digit without special pushups.  We just need to
  66  * decrement digits (and increment weight) to make room for the carry digit.
  67  *
  68  * If buf is NULL then the digit buffer isn't actually palloc'd and should
  69  * not be freed --- see the constants below for an example.
  70  *
  71  * NB: All the variable-level functions are written in a style that makes it
  72  * possible to give one and the same variable as argument and destination.
  73  * This is feasible because the digit buffer is separate from the variable.
  74  * ----------
  75  */
  76 typedef unsigned char NumericDigit;
  77
  78 typedef struct NumericVar
  79 {
  80         int                     ndigits;                /* number of digits in digits[] - can be
  81                                                                  * 0! */
  82         int                     weight;                 /* weight of first digit */
  83         int                     rscale;                 /* result scale */
  84         int                     dscale;                 /* display scale */
  85         int                     sign;                   /* NUMERIC_POS, NUMERIC_NEG, or
  86                                                                  * NUMERIC_NAN */
  87         NumericDigit *buf;                      /* start of palloc'd space for digits[] */
  88         NumericDigit *digits;           /* decimal digits */
  89 } NumericVar;
  90
  91
  92 /* ----------
  93  * Local data
  94  * ----------
  95  */
  96 static int      global_rscale = NUMERIC_MIN_RESULT_SCALE;
  97
  98 /* ----------
  99  * Some preinitialized variables we need often
 100  * ----------
 101  */
 102 static NumericDigit const_zero_data[1] = {0};
 103 static NumericVar const_zero =
 104 {0, 0, 0, 0, NUMERIC_POS, NULL, const_zero_data};
 105
 106 static NumericDigit const_one_data[1] = {1};
 107 static NumericVar const_one =
 108 {1, 0, 0, 0, NUMERIC_POS, NULL, const_one_data};
 109
 110 static NumericDigit const_two_data[1] = {2};
 111 static NumericVar const_two =
 112 {1, 0, 0, 0, NUMERIC_POS, NULL, const_two_data};
 113
 114 static NumericVar const_nan =
 115 {0, 0, 0, 0, NUMERIC_NAN, NULL, NULL};
 116
 117
 118
 119 /* ----------
 120  * Local functions
 121  * ----------
 122  */
 123
 124 #ifdef NUMERIC_DEBUG
 125 static void dump_numeric(char *str, Numeric num);
 126 static void dump_var(char *str, NumericVar *var);
 127
 128 #else
 129 #define dump_numeric(s,n)
 130 #define dump_var(s,v)
 131 #endif
 132
 133 #define digitbuf_alloc(size)  ((NumericDigit *) palloc(size))
 134 #define digitbuf_free(buf)      \
 135         do { \
 136                  if ((buf) != NULL) \
 137                          pfree(buf); \
 138         } while (0)
 139
 140 #define init_var(v)             memset(v,0,sizeof(NumericVar))
 141 static void alloc_var(NumericVar *var, int ndigits);
 142 static void free_var(NumericVar *var);
 143 static void zero_var(NumericVar *var);
 144
 145 static void set_var_from_str(char *str, NumericVar *dest);
 146 static void set_var_from_num(Numeric value, NumericVar *dest);
 147 static void set_var_from_var(NumericVar *value, NumericVar *dest);
 148 static char *get_str_from_var(NumericVar *var, int dscale);
 149
 150 static Numeric make_result(NumericVar *var);
 151
 152 static void apply_typmod(NumericVar *var, int32 typmod);
 153
 154 static int      cmp_var(NumericVar *var1, NumericVar *var2);
 155 static void add_var(NumericVar *var1, NumericVar *var2, NumericVar *result);
 156 static void sub_var(NumericVar *var1, NumericVar *var2, NumericVar *result);
 157 static void mul_var(NumericVar *var1, NumericVar *var2, NumericVar *result);
 158 static void div_var(NumericVar *var1, NumericVar *var2, NumericVar *result);
 159 static void mod_var(NumericVar *var1, NumericVar *var2, NumericVar *result);
 160 static void ceil_var(NumericVar *var, NumericVar *result);
 161 static void floor_var(NumericVar *var, NumericVar *result);
 162
 163 static void sqrt_var(NumericVar *arg, NumericVar *result);
 164 static void exp_var(NumericVar *arg, NumericVar *result);
 165 static void ln_var(NumericVar *arg, NumericVar *result);
 166 static void log_var(NumericVar *base, NumericVar *num, NumericVar *result);
 167 static void power_var(NumericVar *base, NumericVar *exp, NumericVar *result);
 168
 169 static int      cmp_abs(NumericVar *var1, NumericVar *var2);
 170 static void add_abs(NumericVar *var1, NumericVar *var2, NumericVar *result);
 171 static void sub_abs(NumericVar *var1, NumericVar *var2, NumericVar *result);
 172
 173
 174
 175 /* ----------------------------------------------------------------------
 176  *
 177  * Input-, output- and rounding-functions
 178  *
 179  * ----------------------------------------------------------------------
 180  */
 181
 182
 183 /* ----------
 184  * numeric_in() -
 185  *
 186  *      Input function for numeric data type
 187  * ----------
 188  */
 189 Datum
 190 numeric_in(PG_FUNCTION_ARGS)
 191 {
 192         char       *str = PG_GETARG_CSTRING(0);
 193 #ifdef NOT_USED
 194         Oid                     typelem = PG_GETARG_OID(1);
 195 #endif
 196         int32           typmod = PG_GETARG_INT32(2);
 197         NumericVar      value;
 198         Numeric         res;
 199
 200         /* ----------
 201          * Check for NaN
 202          * ----------
 203          */
 204         if (strcmp(str, "NaN") == 0)
 205                 PG_RETURN_NUMERIC(make_result(&const_nan));
 206
 207         /* ----------
 208          * Use set_var_from_str() to parse the input string
 209          * and return it in the packed DB storage format
 210          * ----------
 211          */
 212         init_var(&value);
 213         set_var_from_str(str, &value);
 214
 215         apply_typmod(&value, typmod);
 216
 217         res = make_result(&value);
 218         free_var(&value);
 219
 220         PG_RETURN_NUMERIC(res);
 221 }
 222
 223
 224 /* ----------
 225  * numeric_out() -
 226  *
 227  *      Output function for numeric data type
 228  * ----------
 229  */
 230 Datum
 231 numeric_out(PG_FUNCTION_ARGS)
 232 {
 233         Numeric         num = PG_GETARG_NUMERIC(0);
 234         NumericVar      x;
 235         char       *str;
 236
 237         /* ----------
 238          * Handle NaN
 239          * ----------
 240          */
 241         if (NUMERIC_IS_NAN(num))
 242                 PG_RETURN_CSTRING(pstrdup("NaN"));
 243
 244         /* ----------
 245          * Get the number in the variable format.
 246          *
 247          * Even if we didn't need to change format, we'd still need to copy
 248          * the value to have a modifiable copy for rounding.  set_var_from_num()
 249          * also guarantees there is extra digit space in case we produce a
 250          * carry out from rounding.
 251          * ----------
 252          */
 253         init_var(&x);
 254         set_var_from_num(num, &x);
 255
 256         str = get_str_from_var(&x, x.dscale);
 257
 258         free_var(&x);
 259
 260         PG_RETURN_CSTRING(str);
 261 }
 262
 263
 264 /* ----------
 265  * numeric() -
 266  *
 267  *      This is a special function called by the Postgres database system
 268  *      before a value is stored in a tuples attribute. The precision and
 269  *      scale of the attribute have to be applied on the value.
 270  * ----------
 271  */
 272 Datum
 273 numeric(PG_FUNCTION_ARGS)
 274 {
 275         Numeric         num = PG_GETARG_NUMERIC(0);
 276         int32           typmod = PG_GETARG_INT32(1);
 277         Numeric         new;
 278         int32           tmp_typmod;
 279         int                     precision;
 280         int                     scale;
 281         int                     maxweight;
 282         NumericVar      var;
 283
 284         /* ----------
 285          * Handle NaN
 286          * ----------
 287          */
 288         if (NUMERIC_IS_NAN(num))
 289                 PG_RETURN_NUMERIC(make_result(&const_nan));
 290
 291         /* ----------
 292          * If the value isn't a valid type modifier, simply return a
 293          * copy of the input value
 294          * ----------
 295          */
 296         if (typmod < (int32) (VARHDRSZ))
 297         {
 298                 new = (Numeric) palloc(num->varlen);
 299                 memcpy(new, num, num->varlen);
 300                 PG_RETURN_NUMERIC(new);
 301         }
 302
 303         /* ----------
 304          * Get the precision and scale out of the typmod value
 305          * ----------
 306          */
 307         tmp_typmod = typmod - VARHDRSZ;
 308         precision = (tmp_typmod >> 16) & 0xffff;
 309         scale = tmp_typmod & 0xffff;
 310         maxweight = precision - scale;
 311
 312         /* ----------
 313          * If the number is in bounds and due to the present result scale
 314          * no rounding could be necessary, just make a copy of the input
 315          * and modify its scale fields.
 316          * ----------
 317          */
 318         if (num->n_weight < maxweight && scale >= num->n_rscale)
 319         {
 320                 new = (Numeric) palloc(num->varlen);
 321                 memcpy(new, num, num->varlen);
 322                 new->n_rscale = scale;
 323                 new->n_sign_dscale = NUMERIC_SIGN(new) |
 324                         ((uint16) scale & NUMERIC_DSCALE_MASK);
 325                 PG_RETURN_NUMERIC(new);
 326         }
 327
 328         /* ----------
 329          * We really need to fiddle with things - unpack the number into
 330          * a variable and let apply_typmod() do it.
 331          * ----------
 332          */
 333         init_var(&var);
 334
 335         set_var_from_num(num, &var);
 336         apply_typmod(&var, typmod);
 337         new = make_result(&var);
 338
 339         free_var(&var);
 340
 341         PG_RETURN_NUMERIC(new);
 342 }
 343
 344
 345 /* ----------------------------------------------------------------------
 346  *
 347  * Sign manipulation, rounding and the like
 348  *
 349  * ----------------------------------------------------------------------
 350  */
 351
 352
 353 Numeric
 354 numeric_abs(Numeric num)
 355 {
 356         Numeric         res;
 357
 358         /* ----------
 359          * Handle NULL
 360          * ----------
 361          */
 362         if (num == NULL)
 363                 return NULL;
 364
 365         /* ----------
 366          * Handle NaN
 367          * ----------
 368          */
 369         if (NUMERIC_IS_NAN(num))
 370                 return make_result(&const_nan);
 371
 372         /* ----------
 373          * Do it the easy way directly on the packed format
 374          * ----------
 375          */
 376         res = (Numeric) palloc(num->varlen);
 377         memcpy(res, num, num->varlen);
 378
 379         res->n_sign_dscale = NUMERIC_POS | NUMERIC_DSCALE(num);
 380
 381         return res;
 382 }
 383
 384
 385 Numeric
 386 numeric_uminus(Numeric num)
 387 {
 388         Numeric         res;
 389
 390         /* ----------
 391          * Handle NULL
 392          * ----------
 393          */
 394         if (num == NULL)
 395                 return NULL;
 396
 397         /* ----------
 398          * Handle NaN
 399          * ----------
 400          */
 401         if (NUMERIC_IS_NAN(num))
 402                 return make_result(&const_nan);
 403
 404         /* ----------
 405          * Do it the easy way directly on the packed format
 406          * ----------
 407          */
 408         res = (Numeric) palloc(num->varlen);
 409         memcpy(res, num, num->varlen);
 410
 411         /* ----------
 412          * The packed format is known to be totally zero digit trimmed
 413          * always. So we can identify a ZERO by the fact that there
 414          * are no digits at all.  Do nothing to a zero.
 415          * ----------
 416          */
 417         if (num->varlen != NUMERIC_HDRSZ)
 418         {
 419                 /* Else, flip the sign */
 420                 if (NUMERIC_SIGN(num) == NUMERIC_POS)
 421                         res->n_sign_dscale = NUMERIC_NEG | NUMERIC_DSCALE(num);
 422                 else
 423                         res->n_sign_dscale = NUMERIC_POS | NUMERIC_DSCALE(num);
 424         }
 425
 426         return res;
 427 }
 428
 429
 430 Numeric
 431 numeric_sign(Numeric num)
 432 {
 433         Numeric         res;
 434         NumericVar      result;
 435
 436         /* ----------
 437          * Handle NULL
 438          * ----------
 439          */
 440         if (num == NULL)
 441                 return NULL;
 442
 443         /* ----------
 444          * Handle NaN
 445          * ----------
 446          */
 447         if (NUMERIC_IS_NAN(num))
 448                 return make_result(&const_nan);
 449
 450         init_var(&result);
 451
 452         /* ----------
 453          * The packed format is known to be totally zero digit trimmed
 454          * always. So we can identify a ZERO by the fact that there
 455          * are no digits at all.
 456          * ----------
 457          */
 458         if (num->varlen == NUMERIC_HDRSZ)
 459                 set_var_from_var(&const_zero, &result);
 460         else
 461         {
 462                 /* ----------
 463                  * And if there are some, we return a copy of ONE
 464                  * with the sign of our argument
 465                  * ----------
 466                  */
 467                 set_var_from_var(&const_one, &result);
 468                 result.sign = NUMERIC_SIGN(num);
 469         }
 470
 471         res = make_result(&result);
 472         free_var(&result);
 473
 474         return res;
 475 }
 476
 477
 478 /* ----------
 479  * numeric_round() -
 480  *
 481  *      Round a value to have 'scale' digits after the decimal point.
 482  *      We allow negative 'scale', implying rounding before the decimal
 483  *      point --- Oracle interprets rounding that way.
 484  * ----------
 485  */
 486 Datum
 487 numeric_round(PG_FUNCTION_ARGS)
 488 {
 489         Numeric         num = PG_GETARG_NUMERIC(0);
 490         int32           scale = PG_GETARG_INT32(1);
 491         Numeric         res;
 492         NumericVar      arg;
 493         int                     i;
 494
 495         /* ----------
 496          * Handle NaN
 497          * ----------
 498          */
 499         if (NUMERIC_IS_NAN(num))
 500                 PG_RETURN_NUMERIC(make_result(&const_nan));
 501
 502         /* ----------
 503          * Limit the scale value to avoid possible overflow in calculations below.
 504          * ----------
 505          */
 506         scale = MIN(NUMERIC_MAX_RESULT_SCALE,
 507                                 MAX(-NUMERIC_MAX_RESULT_SCALE, scale));
 508
 509         /* ----------
 510          * Unpack the argument and round it at the proper digit position
 511          * ----------
 512          */
 513         init_var(&arg);
 514         set_var_from_num(num, &arg);
 515
 516         i = arg.weight + scale + 1;
 517
 518         if (i < arg.ndigits)
 519         {
 520
 521                 /*
 522                  * If i = 0, the value loses all digits, but could round up if its
 523                  * first digit is more than 4.  If i < 0 the result must be 0.
 524                  */
 525                 if (i < 0)
 526                         arg.ndigits = 0;
 527                 else
 528                 {
 529                         int                     carry = (arg.digits[i] > 4) ? 1 : 0;
 530
 531                         arg.ndigits = i;
 532
 533                         while (carry)
 534                         {
 535                                 carry += arg.digits[--i];
 536                                 arg.digits[i] = carry % 10;
 537                                 carry /= 10;
 538                         }
 539
 540                         if (i < 0)
 541                         {
 542                                 Assert(i == -1);/* better not have added more than 1 digit */
 543                                 Assert(arg.digits > arg.buf);
 544                                 arg.digits--;
 545                                 arg.ndigits++;
 546                                 arg.weight++;
 547                         }
 548                 }
 549         }
 550
 551         /* ----------
 552          * Set result's scale to something reasonable.
 553          * ----------
 554          */
 555         scale = MIN(NUMERIC_MAX_DISPLAY_SCALE, MAX(0, scale));
 556         arg.rscale = scale;
 557         arg.dscale = scale;
 558
 559         /* ----------
 560          * Return the rounded result
 561          * ----------
 562          */
 563         res = make_result(&arg);
 564
 565         free_var(&arg);
 566         PG_RETURN_NUMERIC(res);
 567 }
 568
 569
 570 /* ----------
 571  * numeric_trunc() -
 572  *
 573  *      Truncate a value to have 'scale' digits after the decimal point.
 574  *      We allow negative 'scale', implying a truncation before the decimal
 575  *      point --- Oracle interprets truncation that way.
 576  * ----------
 577  */
 578 Datum
 579 numeric_trunc(PG_FUNCTION_ARGS)
 580 {
 581         Numeric         num = PG_GETARG_NUMERIC(0);
 582         int32           scale = PG_GETARG_INT32(1);
 583         Numeric         res;
 584         NumericVar      arg;
 585
 586         /* ----------
 587          * Handle NaN
 588          * ----------
 589          */
 590         if (NUMERIC_IS_NAN(num))
 591                 PG_RETURN_NUMERIC(make_result(&const_nan));
 592
 593         /* ----------
 594          * Limit the scale value to avoid possible overflow in calculations below.
 595          * ----------
 596          */
 597         scale = MIN(NUMERIC_MAX_RESULT_SCALE,
 598                                 MAX(-NUMERIC_MAX_RESULT_SCALE, scale));
 599
 600         /* ----------
 601          * Unpack the argument and truncate it at the proper digit position
 602          * ----------
 603          */
 604         init_var(&arg);
 605         set_var_from_num(num, &arg);
 606
 607         arg.ndigits = MIN(arg.ndigits, MAX(0, arg.weight + scale + 1));
 608
 609         /* ----------
 610          * Set result's scale to something reasonable.
 611          * ----------
 612          */
 613         scale = MIN(NUMERIC_MAX_DISPLAY_SCALE, MAX(0, scale));
 614         arg.rscale = scale;
 615         arg.dscale = scale;
 616
 617         /* ----------
 618          * Return the truncated result
 619          * ----------
 620          */
 621         res = make_result(&arg);
 622
 623         free_var(&arg);
 624         PG_RETURN_NUMERIC(res);
 625 }
 626
 627
 628 /* ----------
 629  * numeric_ceil() -
 630  *
 631  *      Return the smallest integer greater than or equal to the argument
 632  * ----------
 633  */
 634 Numeric
 635 numeric_ceil(Numeric num)
 636 {
 637         Numeric         res;
 638         NumericVar      result;
 639
 640         if (num == NULL)
 641                 return NULL;
 642
 643         if (NUMERIC_IS_NAN(num))
 644                 return make_result(&const_nan);
 645
 646         init_var(&result);
 647
 648         set_var_from_num(num, &result);
 649         ceil_var(&result, &result);
 650
 651         result.dscale = 0;
 652
 653         res = make_result(&result);
 654         free_var(&result);
 655
 656         return res;
 657 }
 658
 659
 660 /* ----------
 661  * numeric_floor() -
 662  *
 663  *      Return the largest integer equal to or less than the argument
 664  * ----------
 665  */
 666 Numeric
 667 numeric_floor(Numeric num)
 668 {
 669         Numeric         res;
 670         NumericVar      result;
 671
 672         if (num == NULL)
 673                 return NULL;
 674
 675         if (NUMERIC_IS_NAN(num))
 676                 return make_result(&const_nan);
 677
 678         init_var(&result);
 679
 680         set_var_from_num(num, &result);
 681         floor_var(&result, &result);
 682
 683         result.dscale = 0;
 684
 685         res = make_result(&result);
 686         free_var(&result);
 687
 688         return res;
 689 }
 690
 691
 692 /* ----------------------------------------------------------------------
 693  *
 694  * Comparison functions
 695  *
 696  * ----------------------------------------------------------------------
 697  */
 698
 699
 700 int32
 701 numeric_cmp(Numeric num1, Numeric num2)
 702 {
 703         int                     result;
 704         NumericVar      arg1;
 705         NumericVar      arg2;
 706
 707         if (num1 == NULL || num2 == NULL)
 708                 return (int32) 0;
 709
 710         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 711                 return (int32) 0;
 712
 713         init_var(&arg1);
 714         init_var(&arg2);
 715
 716         set_var_from_num(num1, &arg1);
 717         set_var_from_num(num2, &arg2);
 718
 719         result = cmp_var(&arg1, &arg2);
 720
 721         free_var(&arg1);
 722         free_var(&arg2);
 723
 724         return (int32) ((result == 0) ? 0 : ((result < 0) ? -1 : 1));
 725 }
 726
 727
 728 bool
 729 numeric_eq(Numeric num1, Numeric num2)
 730 {
 731         int                     result;
 732         NumericVar      arg1;
 733         NumericVar      arg2;
 734
 735         if (num1 == NULL || num2 == NULL)
 736                 return FALSE;
 737
 738         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 739                 return FALSE;
 740
 741         init_var(&arg1);
 742         init_var(&arg2);
 743
 744         set_var_from_num(num1, &arg1);
 745         set_var_from_num(num2, &arg2);
 746
 747         result = cmp_var(&arg1, &arg2);
 748
 749         free_var(&arg1);
 750         free_var(&arg2);
 751
 752         return (result == 0);
 753 }
 754
 755
 756 bool
 757 numeric_ne(Numeric num1, Numeric num2)
 758 {
 759         int                     result;
 760         NumericVar      arg1;
 761         NumericVar      arg2;
 762
 763         if (num1 == NULL || num2 == NULL)
 764                 return FALSE;
 765
 766         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 767                 return FALSE;
 768
 769         init_var(&arg1);
 770         init_var(&arg2);
 771
 772         set_var_from_num(num1, &arg1);
 773         set_var_from_num(num2, &arg2);
 774
 775         result = cmp_var(&arg1, &arg2);
 776
 777         free_var(&arg1);
 778         free_var(&arg2);
 779
 780         return (result != 0);
 781 }
 782
 783
 784 bool
 785 numeric_gt(Numeric num1, Numeric num2)
 786 {
 787         int                     result;
 788         NumericVar      arg1;
 789         NumericVar      arg2;
 790
 791         if (num1 == NULL || num2 == NULL)
 792                 return FALSE;
 793
 794         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 795                 return FALSE;
 796
 797         init_var(&arg1);
 798         init_var(&arg2);
 799
 800         set_var_from_num(num1, &arg1);
 801         set_var_from_num(num2, &arg2);
 802
 803         result = cmp_var(&arg1, &arg2);
 804
 805         free_var(&arg1);
 806         free_var(&arg2);
 807
 808         return (result > 0);
 809 }
 810
 811
 812 bool
 813 numeric_ge(Numeric num1, Numeric num2)
 814 {
 815         int                     result;
 816         NumericVar      arg1;
 817         NumericVar      arg2;
 818
 819         if (num1 == NULL || num2 == NULL)
 820                 return FALSE;
 821
 822         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 823                 return FALSE;
 824
 825         init_var(&arg1);
 826         init_var(&arg2);
 827
 828         set_var_from_num(num1, &arg1);
 829         set_var_from_num(num2, &arg2);
 830
 831         result = cmp_var(&arg1, &arg2);
 832
 833         free_var(&arg1);
 834         free_var(&arg2);
 835
 836         return (result >= 0);
 837 }
 838
 839
 840 bool
 841 numeric_lt(Numeric num1, Numeric num2)
 842 {
 843         int                     result;
 844         NumericVar      arg1;
 845         NumericVar      arg2;
 846
 847         if (num1 == NULL || num2 == NULL)
 848                 return FALSE;
 849
 850         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 851                 return FALSE;
 852
 853         init_var(&arg1);
 854         init_var(&arg2);
 855
 856         set_var_from_num(num1, &arg1);
 857         set_var_from_num(num2, &arg2);
 858
 859         result = cmp_var(&arg1, &arg2);
 860
 861         free_var(&arg1);
 862         free_var(&arg2);
 863
 864         return (result < 0);
 865 }
 866
 867
 868 bool
 869 numeric_le(Numeric num1, Numeric num2)
 870 {
 871         int                     result;
 872         NumericVar      arg1;
 873         NumericVar      arg2;
 874
 875         if (num1 == NULL || num2 == NULL)
 876                 return FALSE;
 877
 878         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 879                 return FALSE;
 880
 881         init_var(&arg1);
 882         init_var(&arg2);
 883
 884         set_var_from_num(num1, &arg1);
 885         set_var_from_num(num2, &arg2);
 886
 887         result = cmp_var(&arg1, &arg2);
 888
 889         free_var(&arg1);
 890         free_var(&arg2);
 891
 892         return (result <= 0);
 893 }
 894
 895
 896 /* ----------------------------------------------------------------------
 897  *
 898  * Arithmetic base functions
 899  *
 900  * ----------------------------------------------------------------------
 901  */
 902
 903
 904 /* ----------
 905  * numeric_add() -
 906  *
 907  *      Add two numerics
 908  * ----------
 909  */
 910 Numeric
 911 numeric_add(Numeric num1, Numeric num2)
 912 {
 913         NumericVar      arg1;
 914         NumericVar      arg2;
 915         NumericVar      result;
 916         Numeric         res;
 917
 918         /* ----------
 919          * Handle NULL
 920          * ----------
 921          */
 922         if (num1 == NULL || num2 == NULL)
 923                 return NULL;
 924
 925         /* ----------
 926          * Handle NaN
 927          * ----------
 928          */
 929         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 930                 return make_result(&const_nan);
 931
 932         /* ----------
 933          * Unpack the values, let add_var() compute the result
 934          * and return it. The internals of add_var() will automatically
 935          * set the correct result and display scales in the result.
 936          * ----------
 937          */
 938         init_var(&arg1);
 939         init_var(&arg2);
 940         init_var(&result);
 941
 942         set_var_from_num(num1, &arg1);
 943         set_var_from_num(num2, &arg2);
 944
 945         add_var(&arg1, &arg2, &result);
 946         res = make_result(&result);
 947
 948         free_var(&arg1);
 949         free_var(&arg2);
 950         free_var(&result);
 951
 952         return res;
 953 }
 954
 955
 956 /* ----------
 957  * numeric_sub() -
 958  *
 959  *      Subtract one numeric from another
 960  * ----------
 961  */
 962 Numeric
 963 numeric_sub(Numeric num1, Numeric num2)
 964 {
 965         NumericVar      arg1;
 966         NumericVar      arg2;
 967         NumericVar      result;
 968         Numeric         res;
 969
 970         /* ----------
 971          * Handle NULL
 972          * ----------
 973          */
 974         if (num1 == NULL || num2 == NULL)
 975                 return NULL;
 976
 977         /* ----------
 978          * Handle NaN
 979          * ----------
 980          */
 981         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
 982                 return make_result(&const_nan);
 983
 984         /* ----------
 985          * Unpack the two arguments, let sub_var() compute the
 986          * result and return it.
 987          * ----------
 988          */
 989         init_var(&arg1);
 990         init_var(&arg2);
 991         init_var(&result);
 992
 993         set_var_from_num(num1, &arg1);
 994         set_var_from_num(num2, &arg2);
 995
 996         sub_var(&arg1, &arg2, &result);
 997         res = make_result(&result);
 998
 999         free_var(&arg1);
1000         free_var(&arg2);
1001         free_var(&result);
1002
1003         return res;
1004 }
1005
1006
1007 /* ----------
1008  * numeric_mul() -
1009  *
1010  *      Calculate the product of two numerics
1011  * ----------
1012  */
1013 Numeric
1014 numeric_mul(Numeric num1, Numeric num2)
1015 {
1016         NumericVar      arg1;
1017         NumericVar      arg2;
1018         NumericVar      result;
1019         Numeric         res;
1020
1021         /* ----------
1022          * Handle NULL
1023          * ----------
1024          */
1025         if (num1 == NULL || num2 == NULL)
1026                 return NULL;
1027
1028         /* ----------
1029          * Handle NaN
1030          * ----------
1031          */
1032         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1033                 return make_result(&const_nan);
1034
1035         /* ----------
1036          * Unpack the arguments, let mul_var() compute the result
1037          * and return it. Unlike add_var() and sub_var(), mul_var()
1038          * will round the result to the scale stored in global_rscale.
1039          * In the case of numeric_mul(), which is invoked for the *
1040          * operator on numerics, we set it to the exact representation
1041          * for the product (rscale = sum(rscale of arg1, rscale of arg2)
1042          * and the same for the dscale).
1043          * ----------
1044          */
1045         init_var(&arg1);
1046         init_var(&arg2);
1047         init_var(&result);
1048
1049         set_var_from_num(num1, &arg1);
1050         set_var_from_num(num2, &arg2);
1051
1052         global_rscale = arg1.rscale + arg2.rscale;
1053
1054         mul_var(&arg1, &arg2, &result);
1055
1056         result.dscale = arg1.dscale + arg2.dscale;
1057
1058         res = make_result(&result);
1059
1060         free_var(&arg1);
1061         free_var(&arg2);
1062         free_var(&result);
1063
1064         return res;
1065 }
1066
1067
1068 /* ----------
1069  * numeric_div() -
1070  *
1071  *      Divide one numeric into another
1072  * ----------
1073  */
1074 Numeric
1075 numeric_div(Numeric num1, Numeric num2)
1076 {
1077         NumericVar      arg1;
1078         NumericVar      arg2;
1079         NumericVar      result;
1080         Numeric         res;
1081         int                     res_dscale;
1082
1083         /* ----------
1084          * Handle NULL
1085          * ----------
1086          */
1087         if (num1 == NULL || num2 == NULL)
1088                 return NULL;
1089
1090         /* ----------
1091          * Handle NaN
1092          * ----------
1093          */
1094         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1095                 return make_result(&const_nan);
1096
1097         /* ----------
1098          * Unpack the arguments
1099          * ----------
1100          */
1101         init_var(&arg1);
1102         init_var(&arg2);
1103         init_var(&result);
1104
1105         set_var_from_num(num1, &arg1);
1106         set_var_from_num(num2, &arg2);
1107
1108         /* ----------
1109          * The result scale of a division isn't specified in any
1110          * SQL standard. For Postgres it is the following (where
1111          * SR, DR are the result- and display-scales of the returned
1112          * value, S1, D1, S2 and D2 are the scales of the two arguments,
1113          * The minimum and maximum scales are compile time options from
1114          * numeric.h):
1115          *
1116          *      DR = MIN(MAX(D1 + D2, MIN_DISPLAY_SCALE), MAX_DISPLAY_SCALE)
1117          *      SR = MIN(MAX(MAX(S1 + S2, MIN_RESULT_SCALE), DR + 4), MAX_RESULT_SCALE)
1118          *
1119          * By default, any result is computed with a minimum of 34 digits
1120          * after the decimal point or at least with 4 digits more than
1121          * displayed.
1122          * ----------
1123          */
1124         res_dscale = MAX(arg1.dscale + arg2.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1125         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1126         global_rscale = MAX(arg1.rscale + arg2.rscale,
1127                                                 NUMERIC_MIN_RESULT_SCALE);
1128         global_rscale = MAX(global_rscale, res_dscale + 4);
1129         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1130
1131         /* ----------
1132          * Do the divide, set the display scale and return the result
1133          * ----------
1134          */
1135         div_var(&arg1, &arg2, &result);
1136
1137         result.dscale = res_dscale;
1138
1139         res = make_result(&result);
1140
1141         free_var(&arg1);
1142         free_var(&arg2);
1143         free_var(&result);
1144
1145         return res;
1146 }
1147
1148
1149 /* ----------
1150  * numeric_mod() -
1151  *
1152  *      Calculate the modulo of two numerics
1153  * ----------
1154  */
1155 Numeric
1156 numeric_mod(Numeric num1, Numeric num2)
1157 {
1158         Numeric         res;
1159         NumericVar      arg1;
1160         NumericVar      arg2;
1161         NumericVar      result;
1162
1163         if (num1 == NULL || num2 == NULL)
1164                 return NULL;
1165
1166         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1167                 return make_result(&const_nan);
1168
1169         init_var(&arg1);
1170         init_var(&arg2);
1171         init_var(&result);
1172
1173         set_var_from_num(num1, &arg1);
1174         set_var_from_num(num2, &arg2);
1175
1176         mod_var(&arg1, &arg2, &result);
1177
1178         result.dscale = result.rscale;
1179         res = make_result(&result);
1180
1181         free_var(&result);
1182         free_var(&arg2);
1183         free_var(&arg1);
1184
1185         return res;
1186 }
1187
1188
1189 /* ----------
1190  * numeric_inc() -
1191  *
1192  *      Increment a number by one
1193  * ----------
1194  */
1195 Numeric
1196 numeric_inc(Numeric num)
1197 {
1198         NumericVar      arg;
1199         Numeric         res;
1200
1201         /* ----------
1202          * Handle NULL
1203          * ----------
1204          */
1205         if (num == NULL)
1206                 return NULL;
1207
1208         /* ----------
1209          * Handle NaN
1210          * ----------
1211          */
1212         if (NUMERIC_IS_NAN(num))
1213                 return make_result(&const_nan);
1214
1215         /* ----------
1216          * Compute the result and return it
1217          * ----------
1218          */
1219         init_var(&arg);
1220
1221         set_var_from_num(num, &arg);
1222
1223         add_var(&arg, &const_one, &arg);
1224         res = make_result(&arg);
1225
1226         free_var(&arg);
1227
1228         return res;
1229 }
1230
1231
1232 /* ----------
1233  * numeric_dec() -
1234  *
1235  *      Decrement a number by one
1236  * ----------
1237  */
1238 Numeric
1239 numeric_dec(Numeric num)
1240 {
1241         NumericVar      arg;
1242         Numeric         res;
1243
1244         /* ----------
1245          * Handle NULL
1246          * ----------
1247          */
1248         if (num == NULL)
1249                 return NULL;
1250
1251         /* ----------
1252          * Handle NaN
1253          * ----------
1254          */
1255         if (NUMERIC_IS_NAN(num))
1256                 return make_result(&const_nan);
1257
1258         /* ----------
1259          * Compute the result and return it
1260          * ----------
1261          */
1262         init_var(&arg);
1263
1264         set_var_from_num(num, &arg);
1265
1266         sub_var(&arg, &const_one, &arg);
1267         res = make_result(&arg);
1268
1269         free_var(&arg);
1270
1271         return res;
1272 }
1273
1274
1275 /* ----------
1276  * numeric_smaller() -
1277  *
1278  *      Return the smaller of two numbers
1279  * ----------
1280  */
1281 Numeric
1282 numeric_smaller(Numeric num1, Numeric num2)
1283 {
1284         NumericVar      arg1;
1285         NumericVar      arg2;
1286         Numeric         res;
1287
1288         /* ----------
1289          * Handle NULL
1290          * ----------
1291          */
1292         if (num1 == NULL || num2 == NULL)
1293                 return NULL;
1294
1295         /* ----------
1296          * Handle NaN
1297          * ----------
1298          */
1299         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1300                 return make_result(&const_nan);
1301
1302         /* ----------
1303          * Unpack the values, and decide which is the smaller one
1304          * ----------
1305          */
1306         init_var(&arg1);
1307         init_var(&arg2);
1308
1309         set_var_from_num(num1, &arg1);
1310         set_var_from_num(num2, &arg2);
1311
1312         if (cmp_var(&arg1, &arg2) <= 0)
1313                 res = make_result(&arg1);
1314         else
1315                 res = make_result(&arg2);
1316
1317         free_var(&arg1);
1318         free_var(&arg2);
1319
1320         return res;
1321 }
1322
1323
1324 /* ----------
1325  * numeric_larger() -
1326  *
1327  *      Return the larger of two numbers
1328  * ----------
1329  */
1330 Numeric
1331 numeric_larger(Numeric num1, Numeric num2)
1332 {
1333         NumericVar      arg1;
1334         NumericVar      arg2;
1335         Numeric         res;
1336
1337         /* ----------
1338          * Handle NULL
1339          * ----------
1340          */
1341         if (num1 == NULL || num2 == NULL)
1342                 return NULL;
1343
1344         /* ----------
1345          * Handle NaN
1346          * ----------
1347          */
1348         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1349                 return make_result(&const_nan);
1350
1351         /* ----------
1352          * Unpack the values, and decide which is the larger one
1353          * ----------
1354          */
1355         init_var(&arg1);
1356         init_var(&arg2);
1357
1358         set_var_from_num(num1, &arg1);
1359         set_var_from_num(num2, &arg2);
1360
1361         if (cmp_var(&arg1, &arg2) >= 0)
1362                 res = make_result(&arg1);
1363         else
1364                 res = make_result(&arg2);
1365
1366         free_var(&arg1);
1367         free_var(&arg2);
1368
1369         return res;
1370 }
1371
1372
1373 /* ----------------------------------------------------------------------
1374  *
1375  * Complex math functions
1376  *
1377  * ----------------------------------------------------------------------
1378  */
1379
1380
1381 /* ----------
1382  * numeric_sqrt() -
1383  *
1384  *      Compute the square root of a numeric.
1385  * ----------
1386  */
1387 Numeric
1388 numeric_sqrt(Numeric num)
1389 {
1390         Numeric         res;
1391         NumericVar      arg;
1392         NumericVar      result;
1393         int                     res_dscale;
1394
1395         /* ----------
1396          * Handle NULL
1397          * ----------
1398          */
1399         if (num == NULL)
1400                 return NULL;
1401
1402         /* ----------
1403          * Handle NaN
1404          * ----------
1405          */
1406         if (NUMERIC_IS_NAN(num))
1407                 return make_result(&const_nan);
1408
1409         /* ----------
1410          * Unpack the argument, determine the scales like for divide,
1411          * let sqrt_var() do the calculation and return the result.
1412          * ----------
1413          */
1414         init_var(&arg);
1415         init_var(&result);
1416
1417         set_var_from_num(num, &arg);
1418
1419         res_dscale = MAX(arg.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1420         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1421         global_rscale = MAX(arg.rscale, NUMERIC_MIN_RESULT_SCALE);
1422         global_rscale = MAX(global_rscale, res_dscale + 4);
1423         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1424
1425         sqrt_var(&arg, &result);
1426
1427         result.dscale = res_dscale;
1428
1429         res = make_result(&result);
1430
1431         free_var(&result);
1432         free_var(&arg);
1433
1434         return res;
1435 }
1436
1437
1438 /* ----------
1439  * numeric_exp() -
1440  *
1441  *      Raise e to the power of x
1442  * ----------
1443  */
1444 Numeric
1445 numeric_exp(Numeric num)
1446 {
1447         Numeric         res;
1448         NumericVar      arg;
1449         NumericVar      result;
1450         int                     res_dscale;
1451
1452         /* ----------
1453          * Handle NULL
1454          * ----------
1455          */
1456         if (num == NULL)
1457                 return NULL;
1458
1459         /* ----------
1460          * Handle NaN
1461          * ----------
1462          */
1463         if (NUMERIC_IS_NAN(num))
1464                 return make_result(&const_nan);
1465
1466         /* ----------
1467          * Same procedure like for sqrt().
1468          * ----------
1469          */
1470         init_var(&arg);
1471         init_var(&result);
1472         set_var_from_num(num, &arg);
1473
1474         res_dscale = MAX(arg.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1475         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1476         global_rscale = MAX(arg.rscale, NUMERIC_MIN_RESULT_SCALE);
1477         global_rscale = MAX(global_rscale, res_dscale + 4);
1478         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1479
1480         exp_var(&arg, &result);
1481
1482         result.dscale = res_dscale;
1483
1484         res = make_result(&result);
1485
1486         free_var(&result);
1487         free_var(&arg);
1488
1489         return res;
1490 }
1491
1492
1493 /* ----------
1494  * numeric_ln() -
1495  *
1496  *      Compute the natural logarithm of x
1497  * ----------
1498  */
1499 Numeric
1500 numeric_ln(Numeric num)
1501 {
1502         Numeric         res;
1503         NumericVar      arg;
1504         NumericVar      result;
1505         int                     res_dscale;
1506
1507         /* ----------
1508          * Handle NULL
1509          * ----------
1510          */
1511         if (num == NULL)
1512                 return NULL;
1513
1514         /* ----------
1515          * Handle NaN
1516          * ----------
1517          */
1518         if (NUMERIC_IS_NAN(num))
1519                 return make_result(&const_nan);
1520
1521         /* ----------
1522          * Same procedure like for sqrt()
1523          * ----------
1524          */
1525         init_var(&arg);
1526         init_var(&result);
1527         set_var_from_num(num, &arg);
1528
1529         res_dscale = MAX(arg.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1530         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1531         global_rscale = MAX(arg.rscale, NUMERIC_MIN_RESULT_SCALE);
1532         global_rscale = MAX(global_rscale, res_dscale + 4);
1533         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1534
1535         ln_var(&arg, &result);
1536
1537         result.dscale = res_dscale;
1538
1539         res = make_result(&result);
1540
1541         free_var(&result);
1542         free_var(&arg);
1543
1544         return res;
1545 }
1546
1547
1548 /* ----------
1549  * numeric_log() -
1550  *
1551  *      Compute the logarithm of x in a given base
1552  * ----------
1553  */
1554 Numeric
1555 numeric_log(Numeric num1, Numeric num2)
1556 {
1557         Numeric         res;
1558         NumericVar      arg1;
1559         NumericVar      arg2;
1560         NumericVar      result;
1561         int                     res_dscale;
1562
1563         /* ----------
1564          * Handle NULL
1565          * ----------
1566          */
1567         if (num1 == NULL || num2 == NULL)
1568                 return NULL;
1569
1570         /* ----------
1571          * Handle NaN
1572          * ----------
1573          */
1574         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1575                 return make_result(&const_nan);
1576
1577         /* ----------
1578          * Initialize things and calculate scales
1579          * ----------
1580          */
1581         init_var(&arg1);
1582         init_var(&arg2);
1583         init_var(&result);
1584         set_var_from_num(num1, &arg1);
1585         set_var_from_num(num2, &arg2);
1586
1587         res_dscale = MAX(arg1.dscale + arg2.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1588         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1589         global_rscale = MAX(arg1.rscale + arg2.rscale, NUMERIC_MIN_RESULT_SCALE);
1590         global_rscale = MAX(global_rscale, res_dscale + 4);
1591         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1592
1593         /* ----------
1594          * Call log_var() to compute and return the result
1595          * ----------
1596          */
1597         log_var(&arg1, &arg2, &result);
1598
1599         result.dscale = res_dscale;
1600
1601         res = make_result(&result);
1602
1603         free_var(&result);
1604         free_var(&arg2);
1605         free_var(&arg1);
1606
1607         return res;
1608 }
1609
1610
1611 /* ----------
1612  * numeric_power() -
1613  *
1614  *      Raise m to the power of x
1615  * ----------
1616  */
1617 Numeric
1618 numeric_power(Numeric num1, Numeric num2)
1619 {
1620         Numeric         res;
1621         NumericVar      arg1;
1622         NumericVar      arg2;
1623         NumericVar      result;
1624         int                     res_dscale;
1625
1626         /* ----------
1627          * Handle NULL
1628          * ----------
1629          */
1630         if (num1 == NULL || num2 == NULL)
1631                 return NULL;
1632
1633         /* ----------
1634          * Handle NaN
1635          * ----------
1636          */
1637         if (NUMERIC_IS_NAN(num1) || NUMERIC_IS_NAN(num2))
1638                 return make_result(&const_nan);
1639
1640         /* ----------
1641          * Initialize things and calculate scales
1642          * ----------
1643          */
1644         init_var(&arg1);
1645         init_var(&arg2);
1646         init_var(&result);
1647         set_var_from_num(num1, &arg1);
1648         set_var_from_num(num2, &arg2);
1649
1650         res_dscale = MAX(arg1.dscale + arg2.dscale, NUMERIC_MIN_DISPLAY_SCALE);
1651         res_dscale = MIN(res_dscale, NUMERIC_MAX_DISPLAY_SCALE);
1652         global_rscale = MAX(arg1.rscale + arg2.rscale, NUMERIC_MIN_RESULT_SCALE);
1653         global_rscale = MAX(global_rscale, res_dscale + 4);
1654         global_rscale = MIN(global_rscale, NUMERIC_MAX_RESULT_SCALE);
1655
1656         /* ----------
1657          * Call log_var() to compute and return the result
1658          * ----------
1659          */
1660         power_var(&arg1, &arg2, &result);
1661
1662         result.dscale = res_dscale;
1663
1664         res = make_result(&result);
1665
1666         free_var(&result);
1667         free_var(&arg2);
1668         free_var(&arg1);
1669
1670         return res;
1671 }
1672
1673
1674 /* ----------------------------------------------------------------------
1675  *
1676  * Type conversion functions
1677  *
1678  * ----------------------------------------------------------------------
1679  */
1680
1681
1682 Datum
1683 int4_numeric(PG_FUNCTION_ARGS)
1684 {
1685         int32           val = PG_GETARG_INT32(0);
1686         Numeric         res;
1687         NumericVar      result;
1688         char       *tmp;
1689
1690         init_var(&result);
1691
1692         tmp = DatumGetCString(DirectFunctionCall1(int4out,
1693                                                                                           Int32GetDatum(val)));
1694         set_var_from_str(tmp, &result);
1695         res = make_result(&result);
1696
1697         free_var(&result);
1698         pfree(tmp);
1699
1700         PG_RETURN_NUMERIC(res);
1701 }
1702
1703
1704 int32
1705 numeric_int4(Numeric num)
1706 {
1707         NumericVar      x;
1708         char       *str;
1709         Datum           result;
1710
1711         if (num == NULL)
1712                 return 0;
1713
1714         if (NUMERIC_IS_NAN(num))
1715                 elog(ERROR, "Cannot convert NaN to int4");
1716
1717         /* ----------
1718          * Get the number in the variable format so we can round to integer.
1719          * ----------
1720          */
1721         init_var(&x);
1722         set_var_from_num(num, &x);
1723
1724         str = get_str_from_var(&x, 0);          /* dscale = 0 produces rounding */
1725
1726         free_var(&x);
1727
1728         result = DirectFunctionCall1(int4in, CStringGetDatum(str));
1729         pfree(str);
1730
1731         return result;
1732 }
1733
1734
1735 Numeric
1736 int8_numeric(int64 *val)
1737 {
1738         Numeric         res;
1739         NumericVar      result;
1740         char       *tmp;
1741
1742         init_var(&result);
1743
1744         tmp = DatumGetCString(DirectFunctionCall1(int8out,
1745                                                                                           PointerGetDatum(val)));
1746         set_var_from_str(tmp, &result);
1747         res = make_result(&result);
1748
1749         free_var(&result);
1750         pfree(tmp);
1751
1752         return res;
1753 }
1754
1755
1756 int64 *
1757 numeric_int8(Numeric num)
1758 {
1759         NumericVar      x;
1760         char       *str;
1761         Datum           result;
1762
1763         if (num == NULL)
1764                 return NULL;
1765
1766         if (NUMERIC_IS_NAN(num))
1767                 elog(ERROR, "Cannot convert NaN to int8");
1768
1769         /* ----------
1770          * Get the number in the variable format so we can round to integer.
1771          * ----------
1772          */
1773         init_var(&x);
1774         set_var_from_num(num, &x);
1775
1776         str = get_str_from_var(&x, 0);          /* dscale = 0 produces rounding */
1777
1778         free_var(&x);
1779
1780         result = DirectFunctionCall1(int8in, CStringGetDatum(str));
1781
1782         pfree(str);
1783
1784         return (int64 *) (result);
1785 }
1786
1787
1788 Datum
1789 int2_numeric(PG_FUNCTION_ARGS)
1790 {
1791         int16           val = PG_GETARG_INT16(0);
1792         Numeric         res;
1793         NumericVar      result;
1794         char       *tmp;
1795
1796         init_var(&result);
1797
1798         tmp = DatumGetCString(DirectFunctionCall1(int2out,
1799                                                                                           Int16GetDatum(val)));
1800         set_var_from_str(tmp, &result);
1801         res = make_result(&result);
1802
1803         free_var(&result);
1804         pfree(tmp);
1805
1806         PG_RETURN_NUMERIC(res);
1807 }
1808
1809
1810 Datum
1811 numeric_int2(PG_FUNCTION_ARGS)
1812 {
1813         Numeric         num = PG_GETARG_NUMERIC(0);
1814         NumericVar      x;
1815         char       *str;
1816         Datum           result;
1817
1818         if (NUMERIC_IS_NAN(num))
1819                 elog(ERROR, "Cannot convert NaN to int2");
1820
1821         /* ----------
1822          * Get the number in the variable format so we can round to integer.
1823          * ----------
1824          */
1825         init_var(&x);
1826         set_var_from_num(num, &x);
1827
1828         str = get_str_from_var(&x, 0);          /* dscale = 0 produces rounding */
1829
1830         free_var(&x);
1831
1832         result = DirectFunctionCall1(int2in, CStringGetDatum(str));
1833         pfree(str);
1834
1835         return result;
1836 }
1837
1838
1839 Numeric
1840 float8_numeric(float64 val)
1841 {
1842         Numeric         res;
1843         NumericVar      result;
1844         char            buf[DBL_DIG + 100];
1845
1846         if (val == NULL)
1847                 return NULL;
1848
1849         if (isnan(*val))
1850                 return make_result(&const_nan);
1851
1852         sprintf(buf, "%.*g", DBL_DIG, *val);
1853
1854         init_var(&result);
1855
1856         set_var_from_str(buf, &result);
1857         res = make_result(&result);
1858
1859         free_var(&result);
1860
1861         return res;
1862 }
1863
1864
1865 float64
1866 numeric_float8(Numeric num)
1867 {
1868         char       *tmp;
1869         float64         result;
1870
1871         if (num == NULL)
1872                 return NULL;
1873
1874         if (NUMERIC_IS_NAN(num))
1875         {
1876                 result = (float64) palloc(sizeof(float64data));
1877                 *result = NAN;
1878                 return result;
1879         }
1880
1881         tmp = DatumGetCString(DirectFunctionCall1(numeric_out,
1882                                                                                           NumericGetDatum(num)));
1883         result = float8in(tmp);
1884         pfree(tmp);
1885
1886         return result;
1887 }
1888
1889
1890 Numeric
1891 float4_numeric(float32 val)
1892 {
1893         Numeric         res;
1894         NumericVar      result;
1895         char            buf[FLT_DIG + 100];
1896
1897         if (val == NULL)
1898                 return NULL;
1899
1900         if (isnan(*val))
1901                 return make_result(&const_nan);
1902
1903         sprintf(buf, "%.*g", FLT_DIG, *val);
1904
1905         init_var(&result);
1906
1907         set_var_from_str(buf, &result);
1908         res = make_result(&result);
1909
1910         free_var(&result);
1911
1912         return res;
1913 }
1914
1915
1916 float32
1917 numeric_float4(Numeric num)
1918 {
1919         char       *tmp;
1920         float32         result;
1921
1922         if (num == NULL)
1923                 return NULL;
1924
1925         if (NUMERIC_IS_NAN(num))
1926         {
1927                 result = (float32) palloc(sizeof(float32data));
1928                 *result = NAN;
1929                 return result;
1930         }
1931
1932         tmp = DatumGetCString(DirectFunctionCall1(numeric_out,
1933                                                                                           NumericGetDatum(num)));
1934         result = float4in(tmp);
1935         pfree(tmp);
1936
1937         return result;
1938 }
1939
1940
1941 /* ----------------------------------------------------------------------
1942  *
1943  * Local functions follow
1944  *
1945  * ----------------------------------------------------------------------
1946  */
1947
1948
1949 #ifdef NUMERIC_DEBUG
1950
1951 /* ----------
1952  * dump_numeric() - Dump a value in the db storage format for debugging
1953  * ----------
1954  */
1955 static void
1956 dump_numeric(char *str, Numeric num)
1957 {
1958         int                     i;
1959
1960         printf("%s: NUMERIC w=%d r=%d d=%d ", str, num->n_weight, num->n_rscale,
1961                    NUMERIC_DSCALE(num));
1962         switch (NUMERIC_SIGN(num))
1963         {
1964                 case NUMERIC_POS:
1965                         printf("POS");
1966                         break;
1967                 case NUMERIC_NEG:
1968                         printf("NEG");
1969                         break;
1970                 case NUMERIC_NAN:
1971                         printf("NaN");
1972                         break;
1973                 default:
1974                         printf("SIGN=0x%x", NUMERIC_SIGN(num));
1975                         break;
1976         }
1977
1978         for (i = 0; i < num->varlen - NUMERIC_HDRSZ; i++)
1979                 printf(" %d %d", (num->n_data[i] >> 4) & 0x0f, num->n_data[i] & 0x0f);
1980         printf("\n");
1981 }
1982
1983
1984 /* ----------
1985  * dump_var() - Dump a value in the variable format for debugging
1986  * ----------
1987  */
1988 static void
1989 dump_var(char *str, NumericVar *var)
1990 {
1991         int                     i;
1992
1993         printf("%s: VAR w=%d r=%d d=%d ", str, var->weight, var->rscale,
1994                    var->dscale);
1995         switch (var->sign)
1996         {
1997                 case NUMERIC_POS:
1998                         printf("POS");
1999                         break;
2000                 case NUMERIC_NEG:
2001                         printf("NEG");
2002                         break;
2003                 case NUMERIC_NAN:
2004                         printf("NaN");
2005                         break;
2006                 default:
2007                         printf("SIGN=0x%x", var->sign);
2008                         break;
2009         }
2010
2011         for (i = 0; i < var->ndigits; i++)
2012                 printf(" %d", var->digits[i]);
2013
2014         printf("\n");
2015 }
2016
2017 #endif   /* NUMERIC_DEBUG */
2018
2019
2020 /* ----------
2021  * alloc_var() -
2022  *
2023  *      Allocate a digit buffer of ndigits digits (plus a spare digit for rounding)
2024  * ----------
2025  */
2026 static void
2027 alloc_var(NumericVar *var, int ndigits)
2028 {
2029         digitbuf_free(var->buf);
2030         var->buf = digitbuf_alloc(ndigits + 1);
2031         var->buf[0] = 0;
2032         var->digits = var->buf + 1;
2033         var->ndigits = ndigits;
2034 }
2035
2036
2037 /* ----------
2038  * free_var() -
2039  *
2040  *      Return the digit buffer of a variable to the free pool
2041  * ----------
2042  */
2043 static void
2044 free_var(NumericVar *var)
2045 {
2046         digitbuf_free(var->buf);
2047         var->buf = NULL;
2048         var->digits = NULL;
2049         var->sign = NUMERIC_NAN;
2050 }
2051
2052
2053 /* ----------
2054  * zero_var() -
2055  *
2056  *      Set a variable to ZERO.
2057  *      Note: rscale and dscale are not touched.
2058  * ----------
2059  */
2060 static void
2061 zero_var(NumericVar *var)
2062 {
2063         digitbuf_free(var->buf);
2064         var->buf = NULL;
2065         var->digits = NULL;
2066         var->ndigits = 0;
2067         var->weight = 0;                        /* by convention; doesn't really matter */
2068         var->sign = NUMERIC_POS;        /* anything but NAN... */
2069 }
2070
2071
2072 /* ----------
2073  * set_var_from_str()
2074  *
2075  *      Parse a string and put the number into a variable
2076  * ----------
2077  */
2078 static void
2079 set_var_from_str(char *str, NumericVar *dest)
2080 {
2081         char       *cp = str;
2082         bool            have_dp = FALSE;
2083         int                     i = 0;
2084
2085         while (*cp)
2086         {
2087                 if (!isspace((int) *cp))
2088                         break;
2089                 cp++;
2090         }
2091
2092         alloc_var(dest, strlen(cp));
2093         dest->weight = -1;
2094         dest->dscale = 0;
2095         dest->sign = NUMERIC_POS;
2096
2097         switch (*cp)
2098         {
2099                 case '+':
2100                         dest->sign = NUMERIC_POS;
2101                         cp++;
2102                         break;
2103
2104                 case '-':
2105                         dest->sign = NUMERIC_NEG;
2106                         cp++;
2107                         break;
2108         }
2109
2110         if (*cp == '.')
2111         {
2112                 have_dp = TRUE;
2113                 cp++;
2114         }
2115
2116         if (!isdigit((int) *cp))
2117                 elog(ERROR, "Bad numeric input format '%s'", str);
2118
2119         while (*cp)
2120         {
2121                 if (isdigit((int) *cp))
2122                 {
2123                         dest->digits[i++] = *cp++ - '0';
2124                         if (!have_dp)
2125                                 dest->weight++;
2126                         else
2127                                 dest->dscale++;
2128                 }
2129                 else if (*cp == '.')
2130                 {
2131                         if (have_dp)
2132                                 elog(ERROR, "Bad numeric input format '%s'", str);
2133                         have_dp = TRUE;
2134                         cp++;
2135                 }
2136                 else
2137                         break;
2138         }
2139         dest->ndigits = i;
2140
2141         /* Handle exponent, if any */
2142         if (*cp == 'e' || *cp == 'E')
2143         {
2144                 long            exponent;
2145                 char       *endptr;
2146
2147                 cp++;
2148                 exponent = strtol(cp, &endptr, 10);
2149                 if (endptr == cp)
2150                         elog(ERROR, "Bad numeric input format '%s'", str);
2151                 cp = endptr;
2152                 if (exponent > NUMERIC_MAX_PRECISION ||
2153                         exponent < -NUMERIC_MAX_PRECISION)
2154                         elog(ERROR, "Bad numeric input format '%s'", str);
2155                 dest->weight += (int) exponent;
2156                 dest->dscale -= (int) exponent;
2157                 if (dest->dscale < 0)
2158                         dest->dscale = 0;
2159         }
2160
2161         /* Should be nothing left but spaces */
2162         while (*cp)
2163         {
2164                 if (!isspace((int) *cp))
2165                         elog(ERROR, "Bad numeric input format '%s'", str);
2166                 cp++;
2167         }
2168
2169         /* Strip any leading zeroes */
2170         while (dest->ndigits > 0 && *(dest->digits) == 0)
2171         {
2172                 (dest->digits)++;
2173                 (dest->weight)--;
2174                 (dest->ndigits)--;
2175         }
2176         if (dest->ndigits == 0)
2177                 dest->weight = 0;
2178
2179         dest->rscale = dest->dscale;
2180 }
2181
2182
2183 /*
2184  * set_var_from_num() -
2185  *
2186  *      Parse back the packed db format into a variable
2187  *
2188  */
2189 static void
2190 set_var_from_num(Numeric num, NumericVar *dest)
2191 {
2192         NumericDigit *digit;
2193         int                     i;
2194         int                     n;
2195
2196         n = num->varlen - NUMERIC_HDRSZ;        /* number of digit-pairs in packed
2197                                                                                  * fmt */
2198
2199         alloc_var(dest, n * 2);
2200
2201         dest->weight = num->n_weight;
2202         dest->rscale = num->n_rscale;
2203         dest->dscale = NUMERIC_DSCALE(num);
2204         dest->sign = NUMERIC_SIGN(num);
2205
2206         digit = dest->digits;
2207
2208         for (i = 0; i < n; i++)
2209         {
2210                 unsigned char digitpair = num->n_data[i];
2211
2212                 *digit++ = (digitpair >> 4) & 0x0f;
2213                 *digit++ = digitpair & 0x0f;
2214         }
2215 }
2216
2217
2218 /* ----------
2219  * set_var_from_var() -
2220  *
2221  *      Copy one variable into another
2222  * ----------
2223  */
2224 static void
2225 set_var_from_var(NumericVar *value, NumericVar *dest)
2226 {
2227         NumericDigit *newbuf;
2228
2229         newbuf = digitbuf_alloc(value->ndigits + 1);
2230         newbuf[0] = 0;                          /* spare digit for rounding */
2231         memcpy(newbuf + 1, value->digits, value->ndigits);
2232
2233         digitbuf_free(dest->buf);
2234
2235         memcpy(dest, value, sizeof(NumericVar));
2236         dest->buf = newbuf;
2237         dest->digits = newbuf + 1;
2238 }
2239
2240
2241 /* ----------
2242  * get_str_from_var() -
2243  *
2244  *      Convert a var to text representation (guts of numeric_out).
2245  *      CAUTION: var's contents may be modified by rounding!
2246  *      Caller must have checked for NaN case.
2247  *      Returns a palloc'd string.
2248  * ----------
2249  */
2250 static char *
2251 get_str_from_var(NumericVar *var, int dscale)
2252 {
2253         char       *str;
2254         char       *cp;
2255         int                     i;
2256         int                     d;
2257
2258         /* ----------
2259          * Check if we must round up before printing the value and
2260          * do so.
2261          * ----------
2262          */
2263         i = dscale + var->weight + 1;
2264         if (i >= 0 && var->ndigits > i)
2265         {
2266                 int                     carry = (var->digits[i] > 4) ? 1 : 0;
2267
2268                 var->ndigits = i;
2269
2270                 while (carry)
2271                 {
2272                         carry += var->digits[--i];
2273                         var->digits[i] = carry % 10;
2274                         carry /= 10;
2275                 }
2276
2277                 if (i < 0)
2278                 {
2279                         Assert(i == -1);        /* better not have added more than 1 digit */
2280                         Assert(var->digits > var->buf);
2281                         var->digits--;
2282                         var->ndigits++;
2283                         var->weight++;
2284                 }
2285         }
2286         else
2287                 var->ndigits = MAX(0, MIN(i, var->ndigits));
2288
2289         /* ----------
2290          * Allocate space for the result
2291          * ----------
2292          */
2293         str = palloc(MAX(0, dscale) + MAX(0, var->weight) + 4);
2294         cp = str;
2295
2296         /* ----------
2297          * Output a dash for negative values
2298          * ----------
2299          */
2300         if (var->sign == NUMERIC_NEG)
2301                 *cp++ = '-';
2302
2303         /* ----------
2304          * Output all digits before the decimal point
2305          * ----------
2306          */
2307         i = MAX(var->weight, 0);
2308         d = 0;
2309
2310         while (i >= 0)
2311         {
2312                 if (i <= var->weight && d < var->ndigits)
2313                         *cp++ = var->digits[d++] + '0';
2314                 else
2315                         *cp++ = '0';
2316                 i--;
2317         }
2318
2319         /* ----------
2320          * If requested, output a decimal point and all the digits
2321          * that follow it.
2322          * ----------
2323          */
2324         if (dscale > 0)
2325         {
2326                 *cp++ = '.';
2327                 while (i >= -dscale)
2328                 {
2329                         if (i <= var->weight && d < var->ndigits)
2330                                 *cp++ = var->digits[d++] + '0';
2331                         else
2332                                 *cp++ = '0';
2333                         i--;
2334                 }
2335         }
2336
2337         /* ----------
2338          * terminate the string and return it
2339          * ----------
2340          */
2341         *cp = '\0';
2342         return str;
2343 }
2344
2345
2346 /* ----------
2347  * make_result() -
2348  *
2349  *      Create the packed db numeric format in palloc()'d memory from
2350  *      a variable.  The var's rscale determines the number of digits kept.
2351  * ----------
2352  */
2353 static Numeric
2354 make_result(NumericVar *var)
2355 {
2356         Numeric         result;
2357         NumericDigit *digit = var->digits;
2358         int                     weight = var->weight;
2359         int                     sign = var->sign;
2360         int                     n;
2361         int                     i,
2362                                 j;
2363
2364         if (sign == NUMERIC_NAN)
2365         {
2366                 result = (Numeric) palloc(NUMERIC_HDRSZ);
2367
2368                 result->varlen = NUMERIC_HDRSZ;
2369                 result->n_weight = 0;
2370                 result->n_rscale = 0;
2371                 result->n_sign_dscale = NUMERIC_NAN;
2372
2373                 dump_numeric("make_result()", result);
2374                 return result;
2375         }
2376
2377         n = MAX(0, MIN(var->ndigits, var->weight + var->rscale + 1));
2378
2379         /* truncate leading zeroes */
2380         while (n > 0 && *digit == 0)
2381         {
2382                 digit++;
2383                 weight--;
2384                 n--;
2385         }
2386         /* truncate trailing zeroes */
2387         while (n > 0 && digit[n - 1] == 0)
2388                 n--;
2389
2390         /* If zero result, force to weight=0 and positive sign */
2391         if (n == 0)
2392         {
2393                 weight = 0;
2394                 sign = NUMERIC_POS;
2395         }
2396
2397         result = (Numeric) palloc(NUMERIC_HDRSZ + (n + 1) / 2);
2398         result->varlen = NUMERIC_HDRSZ + (n + 1) / 2;
2399         result->n_weight = weight;
2400         result->n_rscale = var->rscale;
2401         result->n_sign_dscale = sign |
2402                 ((uint16) var->dscale & NUMERIC_DSCALE_MASK);
2403
2404         i = 0;
2405         j = 0;
2406         while (j < n)
2407         {
2408                 unsigned char digitpair = digit[j++] << 4;
2409
2410                 if (j < n)
2411                         digitpair |= digit[j++];
2412                 result->n_data[i++] = digitpair;
2413         }
2414
2415         dump_numeric("make_result()", result);
2416         return result;
2417 }
2418
2419
2420 /* ----------
2421  * apply_typmod() -
2422  *
2423  *      Do bounds checking and rounding according to the attributes
2424  *      typmod field.
2425  * ----------
2426  */
2427 static void
2428 apply_typmod(NumericVar *var, int32 typmod)
2429 {
2430         int                     precision;
2431         int                     scale;
2432         int                     maxweight;
2433         int                     i;
2434
2435         /* Do nothing if we have a default typmod (-1) */
2436         if (typmod < (int32) (VARHDRSZ))
2437                 return;
2438
2439         typmod -= VARHDRSZ;
2440         precision = (typmod >> 16) & 0xffff;
2441         scale = typmod & 0xffff;
2442         maxweight = precision - scale;
2443
2444         /* Round to target scale */
2445         i = scale + var->weight + 1;
2446         if (i >= 0 && var->ndigits > i)
2447         {
2448                 int                     carry = (var->digits[i] > 4) ? 1 : 0;
2449
2450                 var->ndigits = i;
2451
2452                 while (carry)
2453                 {
2454                         carry += var->digits[--i];
2455                         var->digits[i] = carry % 10;
2456                         carry /= 10;
2457                 }
2458
2459                 if (i < 0)
2460                 {
2461                         Assert(i == -1);        /* better not have added more than 1 digit */
2462                         Assert(var->digits > var->buf);
2463                         var->digits--;
2464                         var->ndigits++;
2465                         var->weight++;
2466                 }
2467         }
2468         else
2469                 var->ndigits = MAX(0, MIN(i, var->ndigits));
2470
2471         /* ----------
2472          * Check for overflow - note we can't do this before rounding,
2473          * because rounding could raise the weight.  Also note that the
2474          * var's weight could be inflated by leading zeroes, which will
2475          * be stripped before storage but perhaps might not have been yet.
2476          * In any case, we must recognize a true zero, whose weight doesn't
2477          * mean anything.
2478          * ----------
2479          */
2480         if (var->weight >= maxweight)
2481         {
2482                 /* Determine true weight; and check for all-zero result */
2483                 int                     tweight = var->weight;
2484
2485                 for (i = 0; i < var->ndigits; i++)
2486                 {
2487                         if (var->digits[i])
2488                                 break;
2489                         tweight--;
2490                 }
2491
2492                 if (tweight >= maxweight && i < var->ndigits)
2493                         elog(ERROR, "overflow on numeric "
2494                           "ABS(value) >= 10^%d for field with precision %d scale %d",
2495                                  tweight, precision, scale);
2496         }
2497
2498         var->rscale = scale;
2499         var->dscale = scale;
2500 }
2501
2502
2503 /* ----------
2504  * cmp_var() -
2505  *
2506  *      Compare two values on variable level
2507  * ----------
2508  */
2509 static int
2510 cmp_var(NumericVar *var1, NumericVar *var2)
2511 {
2512         if (var1->ndigits == 0)
2513         {
2514                 if (var2->ndigits == 0)
2515                         return 0;
2516                 if (var2->sign == NUMERIC_NEG)
2517                         return 1;
2518                 return -1;
2519         }
2520         if (var2->ndigits == 0)
2521         {
2522                 if (var1->sign == NUMERIC_POS)
2523                         return 1;
2524                 return -1;
2525         }
2526
2527         if (var1->sign == NUMERIC_POS)
2528         {
2529                 if (var2->sign == NUMERIC_NEG)
2530                         return 1;
2531                 return cmp_abs(var1, var2);
2532         }
2533
2534         if (var2->sign == NUMERIC_POS)
2535                 return -1;
2536
2537         return cmp_abs(var2, var1);
2538 }
2539
2540
2541 /* ----------
2542  * add_var() -
2543  *
2544  *      Full version of add functionality on variable level (handling signs).
2545  *      result might point to one of the operands too without danger.
2546  * ----------
2547  */
2548 static void
2549 add_var(NumericVar *var1, NumericVar *var2, NumericVar *result)
2550 {
2551         /* ----------
2552          * Decide on the signs of the two variables what to do
2553          * ----------
2554          */
2555         if (var1->sign == NUMERIC_POS)
2556         {
2557                 if (var2->sign == NUMERIC_POS)
2558                 {
2559                         /* ----------
2560                          * Both are positive
2561                          * result = +(ABS(var1) + ABS(var2))
2562                          * ----------
2563                          */
2564                         add_abs(var1, var2, result);
2565                         result->sign = NUMERIC_POS;
2566                 }
2567                 else
2568                 {
2569                         /* ----------
2570                          * var1 is positive, var2 is negative
2571                          * Must compare absolute values
2572                          * ----------
2573                          */
2574                         switch (cmp_abs(var1, var2))
2575                         {
2576                                 case 0: /* ----------
2577                                                                                                                                  * ABS(var1) == ABS(var2)
2578                                                                                                                                  * result = ZERO
2579                                                                                                                                  * ----------
2580                                                                                                                                  */
2581                                         zero_var(result);
2582                                         result->rscale = MAX(var1->rscale, var2->rscale);
2583                                         result->dscale = MAX(var1->dscale, var2->dscale);
2584                                         break;
2585
2586                                 case 1: /* ----------
2587                                                                                                                                  * ABS(var1) > ABS(var2)
2588                                                                                                                                  * result = +(ABS(var1) - ABS(var2))
2589                                                                                                                                  * ----------
2590                                                                                                                                  */
2591                                         sub_abs(var1, var2, result);
2592                                         result->sign = NUMERIC_POS;
2593                                         break;
2594
2595                                 case -1:                /* ----------
2596                                                                  * ABS(var1) < ABS(var2)
2597                                                                  * result = -(ABS(var2) - ABS(var1))
2598                                                                  * ----------
2599                                                                  */
2600                                         sub_abs(var2, var1, result);
2601                                         result->sign = NUMERIC_NEG;
2602                                         break;
2603                         }
2604                 }
2605         }
2606         else
2607         {
2608                 if (var2->sign == NUMERIC_POS)
2609                 {
2610                         /* ----------
2611                          * var1 is negative, var2 is positive
2612                          * Must compare absolute values
2613                          * ----------
2614                          */
2615                         switch (cmp_abs(var1, var2))
2616                         {
2617                                 case 0: /* ----------
2618                                                                                                                                  * ABS(var1) == ABS(var2)
2619                                                                                                                                  * result = ZERO
2620                                                                                                                                  * ----------
2621                                                                                                                                  */
2622                                         zero_var(result);
2623                                         result->rscale = MAX(var1->rscale, var2->rscale);
2624                                         result->dscale = MAX(var1->dscale, var2->dscale);
2625                                         break;
2626
2627                                 case 1: /* ----------
2628                                                                                                                                  * ABS(var1) > ABS(var2)
2629                                                                                                                                  * result = -(ABS(var1) - ABS(var2))
2630                                                                                                                                  * ----------
2631                                                                                                                                  */
2632                                         sub_abs(var1, var2, result);
2633                                         result->sign = NUMERIC_NEG;
2634                                         break;
2635
2636                                 case -1:                /* ----------
2637                                                                  * ABS(var1) < ABS(var2)
2638                                                                  * result = +(ABS(var2) - ABS(var1))
2639                                                                  * ----------
2640                                                                  */
2641                                         sub_abs(var2, var1, result);
2642                                         result->sign = NUMERIC_POS;
2643                                         break;
2644                         }
2645                 }
2646                 else
2647                 {
2648                         /* ----------
2649                          * Both are negative
2650                          * result = -(ABS(var1) + ABS(var2))
2651                          * ----------
2652                          */
2653                         add_abs(var1, var2, result);
2654                         result->sign = NUMERIC_NEG;
2655                 }
2656         }
2657 }
2658
2659
2660 /* ----------
2661  * sub_var() -
2662  *
2663  *      Full version of sub functionality on variable level (handling signs).
2664  *      result might point to one of the operands too without danger.
2665  * ----------
2666  */
2667 static void
2668 sub_var(NumericVar *var1, NumericVar *var2, NumericVar *result)
2669 {
2670         /* ----------
2671          * Decide on the signs of the two variables what to do
2672          * ----------
2673          */
2674         if (var1->sign == NUMERIC_POS)
2675         {
2676                 if (var2->sign == NUMERIC_NEG)
2677                 {
2678                         /* ----------
2679                          * var1 is positive, var2 is negative
2680                          * result = +(ABS(var1) + ABS(var2))
2681                          * ----------
2682                          */
2683                         add_abs(var1, var2, result);
2684                         result->sign = NUMERIC_POS;
2685                 }
2686                 else
2687                 {
2688                         /* ----------
2689                          * Both are positive
2690                          * Must compare absolute values
2691                          * ----------
2692                          */
2693                         switch (cmp_abs(var1, var2))
2694                         {
2695                                 case 0: /* ----------
2696                                                                                                                                  * ABS(var1) == ABS(var2)
2697                                                                                                                                  * result = ZERO
2698                                                                                                                                  * ----------
2699                                                                                                                                  */
2700                                         zero_var(result);
2701                                         result->rscale = MAX(var1->rscale, var2->rscale);
2702                                         result->dscale = MAX(var1->dscale, var2->dscale);
2703                                         break;
2704
2705                                 case 1: /* ----------
2706                                                                                                                                  * ABS(var1) > ABS(var2)
2707                                                                                                                                  * result = +(ABS(var1) - ABS(var2))
2708                                                                                                                                  * ----------
2709                                                                                                                                  */
2710                                         sub_abs(var1, var2, result);
2711                                         result->sign = NUMERIC_POS;
2712                                         break;
2713
2714                                 case -1:                /* ----------
2715                                                                  * ABS(var1) < ABS(var2)
2716                                                                  * result = -(ABS(var2) - ABS(var1))
2717                                                                  * ----------
2718                                                                  */
2719                                         sub_abs(var2, var1, result);
2720                                         result->sign = NUMERIC_NEG;
2721                                         break;
2722                         }
2723                 }
2724         }
2725         else
2726         {
2727                 if (var2->sign == NUMERIC_NEG)
2728                 {
2729                         /* ----------
2730                          * Both are negative
2731                          * Must compare absolute values
2732                          * ----------
2733                          */
2734                         switch (cmp_abs(var1, var2))
2735                         {
2736                                 case 0: /* ----------
2737                                                                                                                                  * ABS(var1) == ABS(var2)
2738                                                                                                                                  * result = ZERO
2739                                                                                                                                  * ----------
2740                                                                                                                                  */
2741                                         zero_var(result);
2742                                         result->rscale = MAX(var1->rscale, var2->rscale);
2743                                         result->dscale = MAX(var1->dscale, var2->dscale);
2744                                         break;
2745
2746                                 case 1: /* ----------
2747                                                                                                                                  * ABS(var1) > ABS(var2)
2748                                                                                                                                  * result = -(ABS(var1) - ABS(var2))
2749                                                                                                                                  * ----------
2750                                                                                                                                  */
2751                                         sub_abs(var1, var2, result);
2752                                         result->sign = NUMERIC_NEG;
2753                                         break;
2754
2755                                 case -1:                /* ----------
2756                                                                  * ABS(var1) < ABS(var2)
2757                                                                  * result = +(ABS(var2) - ABS(var1))
2758                                                                  * ----------
2759                                                                  */
2760                                         sub_abs(var2, var1, result);
2761                                         result->sign = NUMERIC_POS;
2762                                         break;
2763                         }
2764                 }
2765                 else
2766                 {
2767                         /* ----------
2768                          * var1 is negative, var2 is positive
2769                          * result = -(ABS(var1) + ABS(var2))
2770                          * ----------
2771                          */
2772                         add_abs(var1, var2, result);
2773                         result->sign = NUMERIC_NEG;
2774                 }
2775         }
2776 }
2777
2778
2779 /* ----------
2780  * mul_var() -
2781  *
2782  *      Multiplication on variable level. Product of var1 * var2 is stored
2783  *      in result.
2784  * ----------
2785  */
2786 static void
2787 mul_var(NumericVar *var1, NumericVar *var2, NumericVar *result)
2788 {
2789         NumericDigit *res_buf;
2790         NumericDigit *res_digits;
2791         int                     res_ndigits;
2792         int                     res_weight;
2793         int                     res_sign;
2794         int                     i,
2795                                 ri,
2796                                 i1,
2797                                 i2;
2798         long            sum = 0;
2799
2800         res_weight = var1->weight + var2->weight + 2;
2801         res_ndigits = var1->ndigits + var2->ndigits + 1;
2802         if (var1->sign == var2->sign)
2803                 res_sign = NUMERIC_POS;
2804         else
2805                 res_sign = NUMERIC_NEG;
2806
2807         res_buf = digitbuf_alloc(res_ndigits);
2808         res_digits = res_buf;
2809         memset(res_digits, 0, res_ndigits);
2810
2811         ri = res_ndigits;
2812         for (i1 = var1->ndigits - 1; i1 >= 0; i1--)
2813         {
2814                 sum = 0;
2815                 i = --ri;
2816
2817                 for (i2 = var2->ndigits - 1; i2 >= 0; i2--)
2818                 {
2819                         sum = sum + res_digits[i] + var1->digits[i1] * var2->digits[i2];
2820                         res_digits[i--] = sum % 10;
2821                         sum /= 10;
2822                 }
2823                 res_digits[i] = sum;
2824         }
2825
2826         i = res_weight + global_rscale + 2;
2827         if (i >= 0 && i < res_ndigits)
2828         {
2829                 sum = (res_digits[i] > 4) ? 1 : 0;
2830                 res_ndigits = i;
2831                 i--;
2832                 while (sum)
2833                 {
2834                         sum += res_digits[i];
2835                         res_digits[i--] = sum % 10;
2836                         sum /= 10;
2837                 }
2838         }
2839
2840         while (res_ndigits > 0 && *res_digits == 0)
2841         {
2842                 res_digits++;
2843                 res_weight--;
2844                 res_ndigits--;
2845         }
2846         while (res_ndigits > 0 && res_digits[res_ndigits - 1] == 0)
2847                 res_ndigits--;
2848
2849         if (res_ndigits == 0)
2850         {
2851                 res_sign = NUMERIC_POS;
2852                 res_weight = 0;
2853         }
2854
2855         digitbuf_free(result->buf);
2856         result->buf = res_buf;
2857         result->digits = res_digits;
2858         result->ndigits = res_ndigits;
2859         result->weight = res_weight;
2860         result->rscale = global_rscale;
2861         result->sign = res_sign;
2862 }
2863
2864
2865 /* ----------
2866  * div_var() -
2867  *
2868  *      Division on variable level.
2869  * ----------
2870  */
2871 static void
2872 div_var(NumericVar *var1, NumericVar *var2, NumericVar *result)
2873 {
2874         NumericDigit *res_digits;
2875         int                     res_ndigits;
2876         int                     res_sign;
2877         int                     res_weight;
2878         NumericVar      dividend;
2879         NumericVar      divisor[10];
2880         int                     ndigits_tmp;
2881         int                     weight_tmp;
2882         int                     rscale_tmp;
2883         int                     ri;
2884         int                     i;
2885         long            guess;
2886         long            first_have;
2887         long            first_div;
2888         int                     first_nextdigit;
2889         int                     stat = 0;
2890
2891         /* ----------
2892          * First of all division by zero check
2893          * ----------
2894          */
2895         ndigits_tmp = var2->ndigits + 1;
2896         if (ndigits_tmp == 1)
2897                 elog(ERROR, "division by zero on numeric");
2898
2899         /* ----------
2900          * Determine the result sign, weight and number of digits to calculate
2901          * ----------
2902          */
2903         if (var1->sign == var2->sign)
2904                 res_sign = NUMERIC_POS;
2905         else
2906                 res_sign = NUMERIC_NEG;
2907         res_weight = var1->weight - var2->weight + 1;
2908         res_ndigits = global_rscale + res_weight;
2909         if (res_ndigits <= 0)
2910                 res_ndigits = 1;
2911
2912         /* ----------
2913          * Now result zero check
2914          * ----------
2915          */
2916         if (var1->ndigits == 0)
2917         {
2918                 zero_var(result);
2919                 result->rscale = global_rscale;
2920                 return;
2921         }
2922
2923         /* ----------
2924          * Initialize local variables
2925          * ----------
2926          */
2927         init_var(&dividend);
2928         for (i = 1; i < 10; i++)
2929                 init_var(&divisor[i]);
2930
2931         /* ----------
2932          * Make a copy of the divisor which has one leading zero digit
2933          * ----------
2934          */
2935         divisor[1].ndigits = ndigits_tmp;
2936         divisor[1].rscale = var2->ndigits;
2937         divisor[1].sign = NUMERIC_POS;
2938         divisor[1].buf = digitbuf_alloc(ndigits_tmp);
2939         divisor[1].digits = divisor[1].buf;
2940         divisor[1].digits[0] = 0;
2941         memcpy(&(divisor[1].digits[1]), var2->digits, ndigits_tmp - 1);
2942
2943         /* ----------
2944          * Make a copy of the dividend
2945          * ----------
2946          */
2947         dividend.ndigits = var1->ndigits;
2948         dividend.weight = 0;
2949         dividend.rscale = var1->ndigits;
2950         dividend.sign = NUMERIC_POS;
2951         dividend.buf = digitbuf_alloc(var1->ndigits);
2952         dividend.digits = dividend.buf;
2953         memcpy(dividend.digits, var1->digits, var1->ndigits);
2954
2955         /* ----------
2956          * Setup the result
2957          * ----------
2958          */
2959         digitbuf_free(result->buf);
2960         result->buf = digitbuf_alloc(res_ndigits + 2);
2961         res_digits = result->buf;
2962         result->digits = res_digits;
2963         result->ndigits = res_ndigits;
2964         result->weight = res_weight;
2965         result->rscale = global_rscale;
2966         result->sign = res_sign;
2967         res_digits[0] = 0;
2968
2969         first_div = divisor[1].digits[1] * 10;
2970         if (ndigits_tmp > 2)
2971                 first_div += divisor[1].digits[2];
2972
2973         first_have = 0;
2974         first_nextdigit = 0;
2975
2976         weight_tmp = 1;
2977         rscale_tmp = divisor[1].rscale;
2978
2979         for (ri = 0; ri <= res_ndigits; ri++)
2980         {
2981                 first_have = first_have * 10;
2982                 if (first_nextdigit >= 0 && first_nextdigit < dividend.ndigits)
2983                         first_have += dividend.digits[first_nextdigit];
2984                 first_nextdigit++;
2985
2986                 guess = (first_have * 10) / first_div + 1;
2987                 if (guess > 9)
2988                         guess = 9;
2989
2990                 while (guess > 0)
2991                 {
2992                         if (divisor[guess].buf == NULL)
2993                         {
2994                                 int                     i;
2995                                 long            sum = 0;
2996
2997                                 memcpy(&divisor[guess], &divisor[1], sizeof(NumericVar));
2998                                 divisor[guess].buf = digitbuf_alloc(divisor[guess].ndigits);
2999                                 divisor[guess].digits = divisor[guess].buf;
3000                                 for (i = divisor[1].ndigits - 1; i >= 0; i--)
3001                                 {
3002                                         sum += divisor[1].digits[i] * guess;
3003                                         divisor[guess].digits[i] = sum % 10;
3004                                         sum /= 10;
3005                                 }
3006                         }
3007
3008                         divisor[guess].weight = weight_tmp;
3009                         divisor[guess].rscale = rscale_tmp;
3010
3011                         stat = cmp_abs(&dividend, &divisor[guess]);
3012                         if (stat >= 0)
3013                                 break;
3014
3015                         guess--;
3016                 }
3017
3018                 res_digits[ri + 1] = guess;
3019                 if (stat == 0)
3020                 {
3021                         ri++;
3022                         break;
3023                 }
3024
3025                 weight_tmp--;
3026                 rscale_tmp++;
3027
3028                 if (guess == 0)
3029                         continue;
3030
3031                 sub_abs(&dividend, &divisor[guess], &dividend);
3032
3033                 first_nextdigit = dividend.weight - weight_tmp;
3034                 first_have = 0;
3035                 if (first_nextdigit >= 0 && first_nextdigit < dividend.ndigits)
3036                         first_have = dividend.digits[first_nextdigit];
3037                 first_nextdigit++;
3038         }
3039
3040         result->ndigits = ri + 1;
3041         if (ri == res_ndigits + 1)
3042         {
3043                 int                     carry = (res_digits[ri] > 4) ? 1 : 0;
3044
3045                 result->ndigits = ri;
3046                 res_digits[ri] = 0;
3047
3048                 while (carry && ri > 0)
3049                 {
3050                         carry += res_digits[--ri];
3051                         res_digits[ri] = carry % 10;
3052                         carry /= 10;
3053                 }
3054         }
3055
3056         while (result->ndigits > 0 && *(result->digits) == 0)
3057         {
3058                 (result->digits)++;
3059                 (result->weight)--;
3060                 (result->ndigits)--;
3061         }
3062         while (result->ndigits > 0 && result->digits[result->ndigits - 1] == 0)
3063                 (result->ndigits)--;
3064         if (result->ndigits == 0)
3065                 result->sign = NUMERIC_POS;
3066
3067         /*
3068          * Tidy up
3069          *
3070          */
3071         digitbuf_free(dividend.buf);
3072         for (i = 1; i < 10; i++)
3073                 digitbuf_free(divisor[i].buf);
3074 }
3075
3076
3077 /* ----------
3078  * mod_var() -
3079  *
3080  *      Calculate the modulo of two numerics at variable level
3081  * ----------
3082  */
3083 static void
3084 mod_var(NumericVar *var1, NumericVar *var2, NumericVar *result)
3085 {
3086         NumericVar      tmp;
3087         int                     save_global_rscale;
3088
3089         init_var(&tmp);
3090
3091         /* ----------
3092          * We do it by fiddling around with global_rscale and truncating
3093          * the result of the division.
3094          * ----------
3095          */
3096         save_global_rscale = global_rscale;
3097         global_rscale = var2->rscale + 2;
3098
3099         div_var(var1, var2, &tmp);
3100         tmp.rscale = var2->rscale;
3101         tmp.ndigits = MAX(0, MIN(tmp.ndigits, tmp.weight + tmp.rscale + 1));
3102
3103         global_rscale = var2->rscale;
3104         mul_var(var2, &tmp, &tmp);
3105
3106         sub_var(var1, &tmp, result);
3107
3108         global_rscale = save_global_rscale;
3109         free_var(&tmp);
3110 }
3111
3112
3113 /* ----------
3114  * ceil_var() -
3115  *
3116  *      Return the smallest integer greater than or equal to the argument
3117  *      on variable level
3118  * ----------
3119  */
3120 static void
3121 ceil_var(NumericVar *var, NumericVar *result)
3122 {
3123         NumericVar      tmp;
3124
3125         init_var(&tmp);
3126         set_var_from_var(var, &tmp);
3127
3128         tmp.rscale = 0;
3129         tmp.ndigits = MIN(tmp.ndigits, MAX(0, tmp.weight + 1));
3130         if (tmp.sign == NUMERIC_POS && cmp_var(var, &tmp) != 0)
3131                 add_var(&tmp, &const_one, &tmp);
3132
3133         set_var_from_var(&tmp, result);
3134         free_var(&tmp);
3135 }
3136
3137
3138 /* ----------
3139  * floor_var() -
3140  *
3141  *      Return the largest integer equal to or less than the argument
3142  *      on variable level
3143  * ----------
3144  */
3145 static void
3146 floor_var(NumericVar *var, NumericVar *result)
3147 {
3148         NumericVar      tmp;
3149
3150         init_var(&tmp);
3151         set_var_from_var(var, &tmp);
3152
3153         tmp.rscale = 0;
3154         tmp.ndigits = MIN(tmp.ndigits, MAX(0, tmp.weight + 1));
3155         if (tmp.sign == NUMERIC_NEG && cmp_var(var, &tmp) != 0)
3156                 sub_var(&tmp, &const_one, &tmp);
3157
3158         set_var_from_var(&tmp, result);
3159         free_var(&tmp);
3160 }
3161
3162
3163 /* ----------
3164  * sqrt_var() -
3165  *
3166  *      Compute the square root of x using Newtons algorithm
3167  * ----------
3168  */
3169 static void
3170 sqrt_var(NumericVar *arg, NumericVar *result)
3171 {
3172         NumericVar      tmp_arg;
3173         NumericVar      tmp_val;
3174         NumericVar      last_val;
3175         int                     res_rscale;
3176         int                     save_global_rscale;
3177         int                     stat;
3178
3179         save_global_rscale = global_rscale;
3180         global_rscale += 8;
3181         res_rscale = global_rscale;
3182
3183         stat = cmp_var(arg, &const_zero);
3184         if (stat == 0)
3185         {
3186                 set_var_from_var(&const_zero, result);
3187                 result->rscale = res_rscale;
3188                 result->sign = NUMERIC_POS;
3189                 return;
3190         }
3191
3192         if (stat < 0)
3193                 elog(ERROR, "math error on numeric - cannot compute SQRT of negative value");
3194
3195         init_var(&tmp_arg);
3196         init_var(&tmp_val);
3197         init_var(&last_val);
3198
3199         set_var_from_var(arg, &tmp_arg);
3200         set_var_from_var(result, &last_val);
3201
3202         /* ----------
3203          * Initialize the result to the first guess
3204          * ----------
3205          */
3206         digitbuf_free(result->buf);
3207         result->buf = digitbuf_alloc(1);
3208         result->digits = result->buf;
3209         result->digits[0] = tmp_arg.digits[0] / 2;
3210         if (result->digits[0] == 0)
3211                 result->digits[0] = 1;
3212         result->ndigits = 1;
3213         result->weight = tmp_arg.weight / 2;
3214         result->rscale = res_rscale;
3215         result->sign = NUMERIC_POS;
3216
3217         for (;;)
3218         {
3219                 div_var(&tmp_arg, result, &tmp_val);
3220
3221                 add_var(result, &tmp_val, result);
3222                 div_var(result, &const_two, result);
3223
3224                 if (cmp_var(&last_val, result) == 0)
3225                         break;
3226                 set_var_from_var(result, &last_val);
3227         }
3228
3229         free_var(&last_val);
3230         free_var(&tmp_val);
3231         free_var(&tmp_arg);
3232
3233         global_rscale = save_global_rscale;
3234         div_var(result, &const_one, result);
3235 }
3236
3237
3238 /* ----------
3239  * exp_var() -
3240  *
3241  *      Raise e to the power of x
3242  * ----------
3243  */
3244 static void
3245 exp_var(NumericVar *arg, NumericVar *result)
3246 {
3247         NumericVar      x;
3248         NumericVar      xpow;
3249         NumericVar      ifac;
3250         NumericVar      elem;
3251         NumericVar      ni;
3252         int                     d;
3253         int                     i;
3254         int                     ndiv2 = 0;
3255         bool            xneg = FALSE;
3256         int                     save_global_rscale;
3257
3258         init_var(&x);
3259         init_var(&xpow);
3260         init_var(&ifac);
3261         init_var(&elem);
3262         init_var(&ni);
3263
3264         set_var_from_var(arg, &x);
3265
3266         if (x.sign == NUMERIC_NEG)
3267         {
3268                 xneg = TRUE;
3269                 x.sign = NUMERIC_POS;
3270         }
3271
3272         save_global_rscale = global_rscale;
3273         global_rscale = 0;
3274         for (i = x.weight, d = 0; i >= 0; i--, d++)
3275         {
3276                 global_rscale *= 10;
3277                 if (d < x.ndigits)
3278                         global_rscale += x.digits[d];
3279                 if (global_rscale >= 1000)
3280                         elog(ERROR, "argument for EXP() too big");
3281         }
3282
3283         global_rscale = global_rscale / 2 + save_global_rscale + 8;
3284
3285         while (cmp_var(&x, &const_one) > 0)
3286         {
3287                 ndiv2++;
3288                 global_rscale++;
3289                 div_var(&x, &const_two, &x);
3290         }
3291
3292         add_var(&const_one, &x, result);
3293         set_var_from_var(&x, &xpow);
3294         set_var_from_var(&const_one, &ifac);
3295         set_var_from_var(&const_one, &ni);
3296
3297         for (i = 2; TRUE; i++)
3298         {
3299                 add_var(&ni, &const_one, &ni);
3300                 mul_var(&xpow, &x, &xpow);
3301                 mul_var(&ifac, &ni, &ifac);
3302                 div_var(&xpow, &ifac, &elem);
3303
3304                 if (elem.ndigits == 0)
3305                         break;
3306
3307                 add_var(result, &elem, result);
3308         }
3309
3310         while (ndiv2-- > 0)
3311                 mul_var(result, result, result);
3312
3313         global_rscale = save_global_rscale;
3314         if (xneg)
3315                 div_var(&const_one, result, result);
3316         else
3317                 div_var(result, &const_one, result);
3318
3319         result->sign = NUMERIC_POS;
3320
3321         free_var(&x);
3322         free_var(&xpow);
3323         free_var(&ifac);
3324         free_var(&elem);
3325         free_var(&ni);
3326 }
3327
3328
3329 /* ----------
3330  * ln_var() -
3331  *
3332  *      Compute the natural log of x
3333  * ----------
3334  */
3335 static void
3336 ln_var(NumericVar *arg, NumericVar *result)
3337 {
3338         NumericVar      x;
3339         NumericVar      xx;
3340         NumericVar      ni;
3341         NumericVar      elem;
3342         NumericVar      fact;
3343         int                     i;
3344         int                     save_global_rscale;
3345
3346         if (cmp_var(arg, &const_zero) <= 0)
3347                 elog(ERROR, "math error on numeric - cannot compute LN of value <= zero");
3348
3349         save_global_rscale = global_rscale;
3350         global_rscale += 8;
3351
3352         init_var(&x);
3353         init_var(&xx);
3354         init_var(&ni);
3355         init_var(&elem);
3356         init_var(&fact);
3357
3358         set_var_from_var(&const_two, &fact);
3359         set_var_from_var(arg, &x);
3360
3361         while (cmp_var(&x, &const_two) >= 0)
3362         {
3363                 sqrt_var(&x, &x);
3364                 mul_var(&fact, &const_two, &fact);
3365         }
3366         set_var_from_str("0.5", &elem);
3367         while (cmp_var(&x, &elem) <= 0)
3368         {
3369                 sqrt_var(&x, &x);
3370                 mul_var(&fact, &const_two, &fact);
3371         }
3372
3373         sub_var(&x, &const_one, result);
3374         add_var(&x, &const_one, &elem);
3375         div_var(result, &elem, result);
3376         set_var_from_var(result, &xx);
3377         mul_var(result, result, &x);
3378
3379         set_var_from_var(&const_one, &ni);
3380
3381         for (i = 2; TRUE; i++)
3382         {
3383                 add_var(&ni, &const_two, &ni);
3384                 mul_var(&xx, &x, &xx);
3385                 div_var(&xx, &ni, &elem);
3386
3387                 if (cmp_var(&elem, &const_zero) == 0)
3388                         break;
3389
3390                 add_var(result, &elem, result);
3391         }
3392
3393         global_rscale = save_global_rscale;
3394         mul_var(result, &fact, result);
3395
3396         free_var(&x);
3397         free_var(&xx);
3398         free_var(&ni);
3399         free_var(&elem);
3400         free_var(&fact);
3401 }
3402
3403
3404 /* ----------
3405  * log_var() -
3406  *
3407  *      Compute the logarithm of x in a given base
3408  * ----------
3409  */
3410 static void
3411 log_var(NumericVar *base, NumericVar *num, NumericVar *result)
3412 {
3413         NumericVar      ln_base;
3414         NumericVar      ln_num;
3415
3416         global_rscale += 8;
3417
3418         init_var(&ln_base);
3419         init_var(&ln_num);
3420
3421         ln_var(base, &ln_base);
3422         ln_var(num, &ln_num);
3423
3424         global_rscale -= 8;
3425
3426         div_var(&ln_num, &ln_base, result);
3427
3428         free_var(&ln_num);
3429         free_var(&ln_base);
3430 }
3431
3432
3433 /* ----------
3434  * power_var() -
3435  *
3436  *      Raise base to the power of exp
3437  * ----------
3438  */
3439 static void
3440 power_var(NumericVar *base, NumericVar *exp, NumericVar *result)
3441 {
3442         NumericVar      ln_base;
3443         NumericVar      ln_num;
3444         int                     save_global_rscale;
3445
3446         save_global_rscale = global_rscale;
3447         global_rscale += global_rscale / 3 + 8;
3448
3449         init_var(&ln_base);
3450         init_var(&ln_num);
3451
3452         ln_var(base, &ln_base);
3453         mul_var(&ln_base, exp, &ln_num);
3454
3455         global_rscale = save_global_rscale;
3456
3457         exp_var(&ln_num, result);
3458
3459         free_var(&ln_num);
3460         free_var(&ln_base);
3461
3462 }
3463
3464
3465 /* ----------------------------------------------------------------------
3466  *
3467  * Following are the lowest level functions that operate unsigned
3468  * on the variable level
3469  *
3470  * ----------------------------------------------------------------------
3471  */
3472
3473
3474 /* ----------
3475  * cmp_abs() -
3476  *
3477  *      Compare the absolute values of var1 and var2
3478  *      Returns:        -1 for ABS(var1) < ABS(var2)
3479  *                              0  for ABS(var1) == ABS(var2)
3480  *                              1  for ABS(var1) > ABS(var2)
3481  * ----------
3482  */
3483 static int
3484 cmp_abs(NumericVar *var1, NumericVar *var2)
3485 {
3486         int                     i1 = 0;
3487         int                     i2 = 0;
3488         int                     w1 = var1->weight;
3489         int                     w2 = var2->weight;
3490         int                     stat;
3491
3492         while (w1 > w2 && i1 < var1->ndigits)
3493         {
3494                 if (var1->digits[i1++] != 0)
3495                         return 1;
3496                 w1--;
3497         }
3498         while (w2 > w1 && i2 < var2->ndigits)
3499         {
3500                 if (var2->digits[i2++] != 0)
3501                         return -1;
3502                 w2--;
3503         }
3504
3505         if (w1 == w2)
3506         {
3507                 while (i1 < var1->ndigits && i2 < var2->ndigits)
3508                 {
3509                         stat = var1->digits[i1++] - var2->digits[i2++];
3510                         if (stat)
3511                         {
3512                                 if (stat > 0)
3513                                         return 1;
3514                                 return -1;
3515                         }
3516                 }
3517         }
3518
3519         while (i1 < var1->ndigits)
3520         {
3521                 if (var1->digits[i1++] != 0)
3522                         return 1;
3523         }
3524         while (i2 < var2->ndigits)
3525         {
3526                 if (var2->digits[i2++] != 0)
3527                         return -1;
3528         }
3529
3530         return 0;
3531 }
3532
3533
3534 /* ----------
3535  * add_abs() -
3536  *
3537  *      Add the absolute values of two variables into result.
3538  *      result might point to one of the operands without danger.
3539  * ----------
3540  */
3541 static void
3542 add_abs(NumericVar *var1, NumericVar *var2, NumericVar *result)
3543 {
3544         NumericDigit *res_buf;
3545         NumericDigit *res_digits;
3546         int                     res_ndigits;
3547         int                     res_weight;
3548         int                     res_rscale;
3549         int                     res_dscale;
3550         int                     i,
3551                                 i1,
3552                                 i2;
3553         int                     carry = 0;
3554
3555         res_weight = MAX(var1->weight, var2->weight) + 1;
3556         res_rscale = MAX(var1->rscale, var2->rscale);
3557         res_dscale = MAX(var1->dscale, var2->dscale);
3558         res_ndigits = res_rscale + res_weight + 1;
3559         if (res_ndigits <= 0)
3560                 res_ndigits = 1;
3561
3562         res_buf = digitbuf_alloc(res_ndigits);
3563         res_digits = res_buf;
3564
3565         i1 = res_rscale + var1->weight + 1;
3566         i2 = res_rscale + var2->weight + 1;
3567         for (i = res_ndigits - 1; i >= 0; i--)
3568         {
3569                 i1--;
3570                 i2--;
3571                 if (i1 >= 0 && i1 < var1->ndigits)
3572                         carry += var1->digits[i1];
3573                 if (i2 >= 0 && i2 < var2->ndigits)
3574                         carry += var2->digits[i2];
3575
3576                 res_digits[i] = carry % 10;
3577                 carry /= 10;
3578         }
3579
3580         while (res_ndigits > 0 && *res_digits == 0)
3581         {
3582                 res_digits++;
3583                 res_weight--;
3584                 res_ndigits--;
3585         }
3586         while (res_ndigits > 0 && res_digits[res_ndigits - 1] == 0)
3587                 res_ndigits--;
3588
3589         if (res_ndigits == 0)
3590                 res_weight = 0;
3591
3592         digitbuf_free(result->buf);
3593         result->ndigits = res_ndigits;
3594         result->buf = res_buf;
3595         result->digits = res_digits;
3596         result->weight = res_weight;
3597         result->rscale = res_rscale;
3598         result->dscale = res_dscale;
3599 }
3600
3601
3602 /* ----------
3603  * sub_abs() -
3604  *
3605  *      Subtract the absolute value of var2 from the absolute value of var1
3606  *      and store in result. result might point to one of the operands
3607  *      without danger.
3608  *
3609  *      ABS(var1) MUST BE GREATER OR EQUAL ABS(var2) !!!
3610  * ----------
3611  */
3612 static void
3613 sub_abs(NumericVar *var1, NumericVar *var2, NumericVar *result)
3614 {
3615         NumericDigit *res_buf;
3616         NumericDigit *res_digits;
3617         int                     res_ndigits;
3618         int                     res_weight;
3619         int                     res_rscale;
3620         int                     res_dscale;
3621         int                     i,
3622                                 i1,
3623                                 i2;
3624         int                     borrow = 0;
3625
3626         res_weight = var1->weight;
3627         res_rscale = MAX(var1->rscale, var2->rscale);
3628         res_dscale = MAX(var1->dscale, var2->dscale);
3629         res_ndigits = res_rscale + res_weight + 1;
3630         if (res_ndigits <= 0)
3631                 res_ndigits = 1;
3632
3633         res_buf = digitbuf_alloc(res_ndigits);
3634         res_digits = res_buf;
3635
3636         i1 = res_rscale + var1->weight + 1;
3637         i2 = res_rscale + var2->weight + 1;
3638         for (i = res_ndigits - 1; i >= 0; i--)
3639         {
3640                 i1--;
3641                 i2--;
3642                 if (i1 >= 0 && i1 < var1->ndigits)
3643                         borrow += var1->digits[i1];
3644                 if (i2 >= 0 && i2 < var2->ndigits)
3645                         borrow -= var2->digits[i2];
3646
3647                 if (borrow < 0)
3648                 {
3649                         res_digits[i] = borrow + 10;
3650                         borrow = -1;
3651                 }
3652                 else
3653                 {
3654                         res_digits[i] = borrow;
3655                         borrow = 0;
3656                 }
3657         }
3658
3659         while (res_ndigits > 0 && *res_digits == 0)
3660         {
3661                 res_digits++;
3662                 res_weight--;
3663                 res_ndigits--;
3664         }
3665         while (res_ndigits > 0 && res_digits[res_ndigits - 1] == 0)
3666                 res_ndigits--;
3667
3668         if (res_ndigits == 0)
3669                 res_weight = 0;
3670
3671         digitbuf_free(result->buf);
3672         result->ndigits = res_ndigits;
3673         result->buf = res_buf;
3674         result->digits = res_digits;
3675         result->weight = res_weight;
3676         result->rscale = res_rscale;
3677         result->dscale = res_dscale;
3678 }