PostgreSQL Source Code  git master
float.c
Go to the documentation of this file.
1 /*-------------------------------------------------------------------------
2  *
3  * float.c
4  * Functions for the built-in floating-point types.
5  *
6  * Portions Copyright (c) 1996-2019, PostgreSQL Global Development Group
7  * Portions Copyright (c) 1994, Regents of the University of California
8  *
9  *
10  * IDENTIFICATION
11  * src/backend/utils/adt/float.c
12  *
13  *-------------------------------------------------------------------------
14  */
15 #include "postgres.h"
16 
17 #include <ctype.h>
18 #include <float.h>
19 #include <math.h>
20 #include <limits.h>
21 
22 #include "catalog/pg_type.h"
23 #include "common/int.h"
24 #include "common/shortest_dec.h"
25 #include "libpq/pqformat.h"
26 #include "miscadmin.h"
27 #include "utils/array.h"
28 #include "utils/float.h"
29 #include "utils/fmgrprotos.h"
30 #include "utils/sortsupport.h"
31 #include "utils/timestamp.h"
32 
33 
34 /*
35  * Configurable GUC parameter
36  *
37  * If >0, use shortest-decimal format for output; this is both the default and
38  * allows for compatibility with clients that explicitly set a value here to
39  * get round-trip-accurate results. If 0 or less, then use the old, slow,
40  * decimal rounding method.
41  */
43 
44 /* Cached constants for degree-based trig functions */
45 static bool degree_consts_set = false;
46 static float8 sin_30 = 0;
48 static float8 asin_0_5 = 0;
49 static float8 acos_0_5 = 0;
50 static float8 atan_1_0 = 0;
51 static float8 tan_45 = 0;
52 static float8 cot_45 = 0;
53 
54 /*
55  * These are intentionally not static; don't "fix" them. They will never
56  * be referenced by other files, much less changed; but we don't want the
57  * compiler to know that, else it might try to precompute expressions
58  * involving them. See comments for init_degree_constants().
59  */
65 
66 /* State for drandom() and setseed() */
67 static bool drandom_seed_set = false;
68 static unsigned short drandom_seed[3] = {0, 0, 0};
69 
70 /* Local function prototypes */
71 static double sind_q1(double x);
72 static double cosd_q1(double x);
73 static void init_degree_constants(void);
74 
75 #ifndef HAVE_CBRT
76 /*
77  * Some machines (in particular, some versions of AIX) have an extern
78  * declaration for cbrt() in <math.h> but fail to provide the actual
79  * function, which causes configure to not set HAVE_CBRT. Furthermore,
80  * their compilers spit up at the mismatch between extern declaration
81  * and static definition. We work around that here by the expedient
82  * of a #define to make the actual name of the static function different.
83  */
84 #define cbrt my_cbrt
85 static double cbrt(double x);
86 #endif /* HAVE_CBRT */
87 
88 
89 /*
90  * Returns -1 if 'val' represents negative infinity, 1 if 'val'
91  * represents (positive) infinity, and 0 otherwise. On some platforms,
92  * this is equivalent to the isinf() macro, but not everywhere: C99
93  * does not specify that isinf() needs to distinguish between positive
94  * and negative infinity.
95  */
96 int
98 {
99  int inf = isinf(val);
100 
101  if (inf == 0)
102  return 0;
103  else if (val > 0)
104  return 1;
105  else
106  return -1;
107 }
108 
109 
110 /* ========== USER I/O ROUTINES ========== */
111 
112 
113 /*
114  * float4in - converts "num" to float4
115  *
116  * Note that this code now uses strtof(), where it used to use strtod().
117  *
118  * The motivation for using strtof() is to avoid a double-rounding problem:
119  * for certain decimal inputs, if you round the input correctly to a double,
120  * and then round the double to a float, the result is incorrect in that it
121  * does not match the result of rounding the decimal value to float directly.
122  *
123  * One of the best examples is 7.038531e-26:
124  *
125  * 0xAE43FDp-107 = 7.03853069185120912085...e-26
126  * midpoint 7.03853100000000022281...e-26
127  * 0xAE43FEp-107 = 7.03853130814879132477...e-26
128  *
129  * making 0xAE43FDp-107 the correct float result, but if you do the conversion
130  * via a double, you get
131  *
132  * 0xAE43FD.7FFFFFF8p-107 = 7.03853099999999907487...e-26
133  * midpoint 7.03853099999999964884...e-26
134  * 0xAE43FD.80000000p-107 = 7.03853100000000022281...e-26
135  * 0xAE43FD.80000008p-107 = 7.03853100000000137076...e-26
136  *
137  * so the value rounds to the double exactly on the midpoint between the two
138  * nearest floats, and then rounding again to a float gives the incorrect
139  * result of 0xAE43FEp-107.
140  *
141  */
142 Datum
144 {
145  char *num = PG_GETARG_CSTRING(0);
146  char *orig_num;
147  float val;
148  char *endptr;
149 
150  /*
151  * endptr points to the first character _after_ the sequence we recognized
152  * as a valid floating point number. orig_num points to the original input
153  * string.
154  */
155  orig_num = num;
156 
157  /* skip leading whitespace */
158  while (*num != '\0' && isspace((unsigned char) *num))
159  num++;
160 
161  /*
162  * Check for an empty-string input to begin with, to avoid the vagaries of
163  * strtod() on different platforms.
164  */
165  if (*num == '\0')
166  ereport(ERROR,
167  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
168  errmsg("invalid input syntax for type %s: \"%s\"",
169  "real", orig_num)));
170 
171  errno = 0;
172  val = strtof(num, &endptr);
173 
174  /* did we not see anything that looks like a double? */
175  if (endptr == num || errno != 0)
176  {
177  int save_errno = errno;
178 
179  /*
180  * C99 requires that strtof() accept NaN, [+-]Infinity, and [+-]Inf,
181  * but not all platforms support all of these (and some accept them
182  * but set ERANGE anyway...) Therefore, we check for these inputs
183  * ourselves if strtof() fails.
184  *
185  * Note: C99 also requires hexadecimal input as well as some extended
186  * forms of NaN, but we consider these forms unportable and don't try
187  * to support them. You can use 'em if your strtof() takes 'em.
188  */
189  if (pg_strncasecmp(num, "NaN", 3) == 0)
190  {
191  val = get_float4_nan();
192  endptr = num + 3;
193  }
194  else if (pg_strncasecmp(num, "Infinity", 8) == 0)
195  {
196  val = get_float4_infinity();
197  endptr = num + 8;
198  }
199  else if (pg_strncasecmp(num, "+Infinity", 9) == 0)
200  {
201  val = get_float4_infinity();
202  endptr = num + 9;
203  }
204  else if (pg_strncasecmp(num, "-Infinity", 9) == 0)
205  {
206  val = -get_float4_infinity();
207  endptr = num + 9;
208  }
209  else if (pg_strncasecmp(num, "inf", 3) == 0)
210  {
211  val = get_float4_infinity();
212  endptr = num + 3;
213  }
214  else if (pg_strncasecmp(num, "+inf", 4) == 0)
215  {
216  val = get_float4_infinity();
217  endptr = num + 4;
218  }
219  else if (pg_strncasecmp(num, "-inf", 4) == 0)
220  {
221  val = -get_float4_infinity();
222  endptr = num + 4;
223  }
224  else if (save_errno == ERANGE)
225  {
226  /*
227  * Some platforms return ERANGE for denormalized numbers (those
228  * that are not zero, but are too close to zero to have full
229  * precision). We'd prefer not to throw error for that, so try to
230  * detect whether it's a "real" out-of-range condition by checking
231  * to see if the result is zero or huge.
232  *
233  * Use isinf() rather than HUGE_VALF on VS2013 because it
234  * generates a spurious overflow warning for -HUGE_VALF. Also use
235  * isinf() if HUGE_VALF is missing.
236  */
237  if (val == 0.0 ||
238 #if !defined(HUGE_VALF) || (defined(_MSC_VER) && (_MSC_VER < 1900))
239  isinf(val)
240 #else
241  (val >= HUGE_VALF || val <= -HUGE_VALF)
242 #endif
243  )
244  ereport(ERROR,
245  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
246  errmsg("\"%s\" is out of range for type real",
247  orig_num)));
248  }
249  else
250  ereport(ERROR,
251  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
252  errmsg("invalid input syntax for type %s: \"%s\"",
253  "real", orig_num)));
254  }
255 #ifdef HAVE_BUGGY_SOLARIS_STRTOD
256  else
257  {
258  /*
259  * Many versions of Solaris have a bug wherein strtod sets endptr to
260  * point one byte beyond the end of the string when given "inf" or
261  * "infinity".
262  */
263  if (endptr != num && endptr[-1] == '\0')
264  endptr--;
265  }
266 #endif /* HAVE_BUGGY_SOLARIS_STRTOD */
267 
268  /* skip trailing whitespace */
269  while (*endptr != '\0' && isspace((unsigned char) *endptr))
270  endptr++;
271 
272  /* if there is any junk left at the end of the string, bail out */
273  if (*endptr != '\0')
274  ereport(ERROR,
275  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
276  errmsg("invalid input syntax for type %s: \"%s\"",
277  "real", orig_num)));
278 
279  PG_RETURN_FLOAT4(val);
280 }
281 
282 /*
283  * float4out - converts a float4 number to a string
284  * using a standard output format
285  */
286 Datum
288 {
289  float4 num = PG_GETARG_FLOAT4(0);
290  char *ascii = (char *) palloc(32);
291  int ndig = FLT_DIG + extra_float_digits;
292 
293  if (extra_float_digits > 0)
294  {
295  float_to_shortest_decimal_buf(num, ascii);
296  PG_RETURN_CSTRING(ascii);
297  }
298 
299  (void) pg_strfromd(ascii, 32, ndig, num);
300  PG_RETURN_CSTRING(ascii);
301 }
302 
303 /*
304  * float4recv - converts external binary format to float4
305  */
306 Datum
308 {
310 
312 }
313 
314 /*
315  * float4send - converts float4 to binary format
316  */
317 Datum
319 {
320  float4 num = PG_GETARG_FLOAT4(0);
322 
323  pq_begintypsend(&buf);
324  pq_sendfloat4(&buf, num);
326 }
327 
328 /*
329  * float8in - converts "num" to float8
330  */
331 Datum
333 {
334  char *num = PG_GETARG_CSTRING(0);
335 
336  PG_RETURN_FLOAT8(float8in_internal(num, NULL, "double precision", num));
337 }
338 
339 /* Convenience macro: set *have_error flag (if provided) or throw error */
340 #define RETURN_ERROR(throw_error, have_error) \
341 do { \
342  if (have_error) { \
343  *have_error = true; \
344  return 0.0; \
345  } else { \
346  throw_error; \
347  } \
348 } while (0)
349 
350 /*
351  * float8in_internal_opt_error - guts of float8in()
352  *
353  * This is exposed for use by functions that want a reasonably
354  * platform-independent way of inputting doubles. The behavior is
355  * essentially like strtod + ereport on error, but note the following
356  * differences:
357  * 1. Both leading and trailing whitespace are skipped.
358  * 2. If endptr_p is NULL, we throw error if there's trailing junk.
359  * Otherwise, it's up to the caller to complain about trailing junk.
360  * 3. In event of a syntax error, the report mentions the given type_name
361  * and prints orig_string as the input; this is meant to support use of
362  * this function with types such as "box" and "point", where what we are
363  * parsing here is just a substring of orig_string.
364  *
365  * "num" could validly be declared "const char *", but that results in an
366  * unreasonable amount of extra casting both here and in callers, so we don't.
367  *
368  * When "*have_error" flag is provided, it's set instead of throwing an
369  * error. This is helpful when caller need to handle errors by itself.
370  */
371 double
372 float8in_internal_opt_error(char *num, char **endptr_p,
373  const char *type_name, const char *orig_string,
374  bool *have_error)
375 {
376  double val;
377  char *endptr;
378 
379  if (have_error)
380  *have_error = false;
381 
382  /* skip leading whitespace */
383  while (*num != '\0' && isspace((unsigned char) *num))
384  num++;
385 
386  /*
387  * Check for an empty-string input to begin with, to avoid the vagaries of
388  * strtod() on different platforms.
389  */
390  if (*num == '\0')
392  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
393  errmsg("invalid input syntax for type %s: \"%s\"",
394  type_name, orig_string))),
395  have_error);
396 
397  errno = 0;
398  val = strtod(num, &endptr);
399 
400  /* did we not see anything that looks like a double? */
401  if (endptr == num || errno != 0)
402  {
403  int save_errno = errno;
404 
405  /*
406  * C99 requires that strtod() accept NaN, [+-]Infinity, and [+-]Inf,
407  * but not all platforms support all of these (and some accept them
408  * but set ERANGE anyway...) Therefore, we check for these inputs
409  * ourselves if strtod() fails.
410  *
411  * Note: C99 also requires hexadecimal input as well as some extended
412  * forms of NaN, but we consider these forms unportable and don't try
413  * to support them. You can use 'em if your strtod() takes 'em.
414  */
415  if (pg_strncasecmp(num, "NaN", 3) == 0)
416  {
417  val = get_float8_nan();
418  endptr = num + 3;
419  }
420  else if (pg_strncasecmp(num, "Infinity", 8) == 0)
421  {
422  val = get_float8_infinity();
423  endptr = num + 8;
424  }
425  else if (pg_strncasecmp(num, "+Infinity", 9) == 0)
426  {
427  val = get_float8_infinity();
428  endptr = num + 9;
429  }
430  else if (pg_strncasecmp(num, "-Infinity", 9) == 0)
431  {
432  val = -get_float8_infinity();
433  endptr = num + 9;
434  }
435  else if (pg_strncasecmp(num, "inf", 3) == 0)
436  {
437  val = get_float8_infinity();
438  endptr = num + 3;
439  }
440  else if (pg_strncasecmp(num, "+inf", 4) == 0)
441  {
442  val = get_float8_infinity();
443  endptr = num + 4;
444  }
445  else if (pg_strncasecmp(num, "-inf", 4) == 0)
446  {
447  val = -get_float8_infinity();
448  endptr = num + 4;
449  }
450  else if (save_errno == ERANGE)
451  {
452  /*
453  * Some platforms return ERANGE for denormalized numbers (those
454  * that are not zero, but are too close to zero to have full
455  * precision). We'd prefer not to throw error for that, so try to
456  * detect whether it's a "real" out-of-range condition by checking
457  * to see if the result is zero or huge.
458  *
459  * On error, we intentionally complain about double precision not
460  * the given type name, and we print only the part of the string
461  * that is the current number.
462  */
463  if (val == 0.0 || val >= HUGE_VAL || val <= -HUGE_VAL)
464  {
465  char *errnumber = pstrdup(num);
466 
467  errnumber[endptr - num] = '\0';
469  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
470  errmsg("\"%s\" is out of range for type double precision",
471  errnumber))),
472  have_error);
473  }
474  }
475  else
477  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
478  errmsg("invalid input syntax for type "
479  "%s: \"%s\"",
480  type_name, orig_string))),
481  have_error);
482  }
483 #ifdef HAVE_BUGGY_SOLARIS_STRTOD
484  else
485  {
486  /*
487  * Many versions of Solaris have a bug wherein strtod sets endptr to
488  * point one byte beyond the end of the string when given "inf" or
489  * "infinity".
490  */
491  if (endptr != num && endptr[-1] == '\0')
492  endptr--;
493  }
494 #endif /* HAVE_BUGGY_SOLARIS_STRTOD */
495 
496  /* skip trailing whitespace */
497  while (*endptr != '\0' && isspace((unsigned char) *endptr))
498  endptr++;
499 
500  /* report stopping point if wanted, else complain if not end of string */
501  if (endptr_p)
502  *endptr_p = endptr;
503  else if (*endptr != '\0')
505  (errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
506  errmsg("invalid input syntax for type "
507  "%s: \"%s\"",
508  type_name, orig_string))),
509  have_error);
510 
511  return val;
512 }
513 
514 /*
515  * Interface to float8in_internal_opt_error() without "have_error" argument.
516  */
517 double
518 float8in_internal(char *num, char **endptr_p,
519  const char *type_name, const char *orig_string)
520 {
521  return float8in_internal_opt_error(num, endptr_p, type_name,
522  orig_string, NULL);
523 }
524 
525 
526 /*
527  * float8out - converts float8 number to a string
528  * using a standard output format
529  */
530 Datum
532 {
533  float8 num = PG_GETARG_FLOAT8(0);
534 
536 }
537 
538 /*
539  * float8out_internal - guts of float8out()
540  *
541  * This is exposed for use by functions that want a reasonably
542  * platform-independent way of outputting doubles.
543  * The result is always palloc'd.
544  */
545 char *
547 {
548  char *ascii = (char *) palloc(32);
549  int ndig = DBL_DIG + extra_float_digits;
550 
551  if (extra_float_digits > 0)
552  {
553  double_to_shortest_decimal_buf(num, ascii);
554  return ascii;
555  }
556 
557  (void) pg_strfromd(ascii, 32, ndig, num);
558  return ascii;
559 }
560 
561 /*
562  * float8recv - converts external binary format to float8
563  */
564 Datum
566 {
568 
570 }
571 
572 /*
573  * float8send - converts float8 to binary format
574  */
575 Datum
577 {
578  float8 num = PG_GETARG_FLOAT8(0);
580 
581  pq_begintypsend(&buf);
582  pq_sendfloat8(&buf, num);
584 }
585 
586 
587 /* ========== PUBLIC ROUTINES ========== */
588 
589 
590 /*
591  * ======================
592  * FLOAT4 BASE OPERATIONS
593  * ======================
594  */
595 
596 /*
597  * float4abs - returns |arg1| (absolute value)
598  */
599 Datum
601 {
602  float4 arg1 = PG_GETARG_FLOAT4(0);
603 
604  PG_RETURN_FLOAT4((float4) fabs(arg1));
605 }
606 
607 /*
608  * float4um - returns -arg1 (unary minus)
609  */
610 Datum
612 {
613  float4 arg1 = PG_GETARG_FLOAT4(0);
614  float4 result;
615 
616  result = -arg1;
617  PG_RETURN_FLOAT4(result);
618 }
619 
620 Datum
622 {
624 
625  PG_RETURN_FLOAT4(arg);
626 }
627 
628 Datum
630 {
631  float4 arg1 = PG_GETARG_FLOAT4(0);
632  float4 arg2 = PG_GETARG_FLOAT4(1);
633  float4 result;
634 
635  if (float4_gt(arg1, arg2))
636  result = arg1;
637  else
638  result = arg2;
639  PG_RETURN_FLOAT4(result);
640 }
641 
642 Datum
644 {
645  float4 arg1 = PG_GETARG_FLOAT4(0);
646  float4 arg2 = PG_GETARG_FLOAT4(1);
647  float4 result;
648 
649  if (float4_lt(arg1, arg2))
650  result = arg1;
651  else
652  result = arg2;
653  PG_RETURN_FLOAT4(result);
654 }
655 
656 /*
657  * ======================
658  * FLOAT8 BASE OPERATIONS
659  * ======================
660  */
661 
662 /*
663  * float8abs - returns |arg1| (absolute value)
664  */
665 Datum
667 {
668  float8 arg1 = PG_GETARG_FLOAT8(0);
669 
670  PG_RETURN_FLOAT8(fabs(arg1));
671 }
672 
673 
674 /*
675  * float8um - returns -arg1 (unary minus)
676  */
677 Datum
679 {
680  float8 arg1 = PG_GETARG_FLOAT8(0);
681  float8 result;
682 
683  result = -arg1;
684  PG_RETURN_FLOAT8(result);
685 }
686 
687 Datum
689 {
691 
692  PG_RETURN_FLOAT8(arg);
693 }
694 
695 Datum
697 {
698  float8 arg1 = PG_GETARG_FLOAT8(0);
699  float8 arg2 = PG_GETARG_FLOAT8(1);
700  float8 result;
701 
702  if (float8_gt(arg1, arg2))
703  result = arg1;
704  else
705  result = arg2;
706  PG_RETURN_FLOAT8(result);
707 }
708 
709 Datum
711 {
712  float8 arg1 = PG_GETARG_FLOAT8(0);
713  float8 arg2 = PG_GETARG_FLOAT8(1);
714  float8 result;
715 
716  if (float8_lt(arg1, arg2))
717  result = arg1;
718  else
719  result = arg2;
720  PG_RETURN_FLOAT8(result);
721 }
722 
723 
724 /*
725  * ====================
726  * ARITHMETIC OPERATORS
727  * ====================
728  */
729 
730 /*
731  * float4pl - returns arg1 + arg2
732  * float4mi - returns arg1 - arg2
733  * float4mul - returns arg1 * arg2
734  * float4div - returns arg1 / arg2
735  */
736 Datum
738 {
739  float4 arg1 = PG_GETARG_FLOAT4(0);
740  float4 arg2 = PG_GETARG_FLOAT4(1);
741 
742  PG_RETURN_FLOAT4(float4_pl(arg1, arg2));
743 }
744 
745 Datum
747 {
748  float4 arg1 = PG_GETARG_FLOAT4(0);
749  float4 arg2 = PG_GETARG_FLOAT4(1);
750 
751  PG_RETURN_FLOAT4(float4_mi(arg1, arg2));
752 }
753 
754 Datum
756 {
757  float4 arg1 = PG_GETARG_FLOAT4(0);
758  float4 arg2 = PG_GETARG_FLOAT4(1);
759 
760  PG_RETURN_FLOAT4(float4_mul(arg1, arg2));
761 }
762 
763 Datum
765 {
766  float4 arg1 = PG_GETARG_FLOAT4(0);
767  float4 arg2 = PG_GETARG_FLOAT4(1);
768 
769  PG_RETURN_FLOAT4(float4_div(arg1, arg2));
770 }
771 
772 /*
773  * float8pl - returns arg1 + arg2
774  * float8mi - returns arg1 - arg2
775  * float8mul - returns arg1 * arg2
776  * float8div - returns arg1 / arg2
777  */
778 Datum
780 {
781  float8 arg1 = PG_GETARG_FLOAT8(0);
782  float8 arg2 = PG_GETARG_FLOAT8(1);
783 
784  PG_RETURN_FLOAT8(float8_pl(arg1, arg2));
785 }
786 
787 Datum
789 {
790  float8 arg1 = PG_GETARG_FLOAT8(0);
791  float8 arg2 = PG_GETARG_FLOAT8(1);
792 
793  PG_RETURN_FLOAT8(float8_mi(arg1, arg2));
794 }
795 
796 Datum
798 {
799  float8 arg1 = PG_GETARG_FLOAT8(0);
800  float8 arg2 = PG_GETARG_FLOAT8(1);
801 
802  PG_RETURN_FLOAT8(float8_mul(arg1, arg2));
803 }
804 
805 Datum
807 {
808  float8 arg1 = PG_GETARG_FLOAT8(0);
809  float8 arg2 = PG_GETARG_FLOAT8(1);
810 
811  PG_RETURN_FLOAT8(float8_div(arg1, arg2));
812 }
813 
814 
815 /*
816  * ====================
817  * COMPARISON OPERATORS
818  * ====================
819  */
820 
821 /*
822  * float4{eq,ne,lt,le,gt,ge} - float4/float4 comparison operations
823  */
824 int
826 {
827  if (float4_gt(a, b))
828  return 1;
829  if (float4_lt(a, b))
830  return -1;
831  return 0;
832 }
833 
834 Datum
836 {
837  float4 arg1 = PG_GETARG_FLOAT4(0);
838  float4 arg2 = PG_GETARG_FLOAT4(1);
839 
840  PG_RETURN_BOOL(float4_eq(arg1, arg2));
841 }
842 
843 Datum
845 {
846  float4 arg1 = PG_GETARG_FLOAT4(0);
847  float4 arg2 = PG_GETARG_FLOAT4(1);
848 
849  PG_RETURN_BOOL(float4_ne(arg1, arg2));
850 }
851 
852 Datum
854 {
855  float4 arg1 = PG_GETARG_FLOAT4(0);
856  float4 arg2 = PG_GETARG_FLOAT4(1);
857 
858  PG_RETURN_BOOL(float4_lt(arg1, arg2));
859 }
860 
861 Datum
863 {
864  float4 arg1 = PG_GETARG_FLOAT4(0);
865  float4 arg2 = PG_GETARG_FLOAT4(1);
866 
867  PG_RETURN_BOOL(float4_le(arg1, arg2));
868 }
869 
870 Datum
872 {
873  float4 arg1 = PG_GETARG_FLOAT4(0);
874  float4 arg2 = PG_GETARG_FLOAT4(1);
875 
876  PG_RETURN_BOOL(float4_gt(arg1, arg2));
877 }
878 
879 Datum
881 {
882  float4 arg1 = PG_GETARG_FLOAT4(0);
883  float4 arg2 = PG_GETARG_FLOAT4(1);
884 
885  PG_RETURN_BOOL(float4_ge(arg1, arg2));
886 }
887 
888 Datum
890 {
891  float4 arg1 = PG_GETARG_FLOAT4(0);
892  float4 arg2 = PG_GETARG_FLOAT4(1);
893 
895 }
896 
897 static int
899 {
900  float4 arg1 = DatumGetFloat4(x);
901  float4 arg2 = DatumGetFloat4(y);
902 
903  return float4_cmp_internal(arg1, arg2);
904 }
905 
906 Datum
908 {
910 
911  ssup->comparator = btfloat4fastcmp;
912  PG_RETURN_VOID();
913 }
914 
915 /*
916  * float8{eq,ne,lt,le,gt,ge} - float8/float8 comparison operations
917  */
918 int
920 {
921  if (float8_gt(a, b))
922  return 1;
923  if (float8_lt(a, b))
924  return -1;
925  return 0;
926 }
927 
928 Datum
930 {
931  float8 arg1 = PG_GETARG_FLOAT8(0);
932  float8 arg2 = PG_GETARG_FLOAT8(1);
933 
934  PG_RETURN_BOOL(float8_eq(arg1, arg2));
935 }
936 
937 Datum
939 {
940  float8 arg1 = PG_GETARG_FLOAT8(0);
941  float8 arg2 = PG_GETARG_FLOAT8(1);
942 
943  PG_RETURN_BOOL(float8_ne(arg1, arg2));
944 }
945 
946 Datum
948 {
949  float8 arg1 = PG_GETARG_FLOAT8(0);
950  float8 arg2 = PG_GETARG_FLOAT8(1);
951 
952  PG_RETURN_BOOL(float8_lt(arg1, arg2));
953 }
954 
955 Datum
957 {
958  float8 arg1 = PG_GETARG_FLOAT8(0);
959  float8 arg2 = PG_GETARG_FLOAT8(1);
960 
961  PG_RETURN_BOOL(float8_le(arg1, arg2));
962 }
963 
964 Datum
966 {
967  float8 arg1 = PG_GETARG_FLOAT8(0);
968  float8 arg2 = PG_GETARG_FLOAT8(1);
969 
970  PG_RETURN_BOOL(float8_gt(arg1, arg2));
971 }
972 
973 Datum
975 {
976  float8 arg1 = PG_GETARG_FLOAT8(0);
977  float8 arg2 = PG_GETARG_FLOAT8(1);
978 
979  PG_RETURN_BOOL(float8_ge(arg1, arg2));
980 }
981 
982 Datum
984 {
985  float8 arg1 = PG_GETARG_FLOAT8(0);
986  float8 arg2 = PG_GETARG_FLOAT8(1);
987 
989 }
990 
991 static int
993 {
994  float8 arg1 = DatumGetFloat8(x);
995  float8 arg2 = DatumGetFloat8(y);
996 
997  return float8_cmp_internal(arg1, arg2);
998 }
999 
1000 Datum
1002 {
1004 
1005  ssup->comparator = btfloat8fastcmp;
1006  PG_RETURN_VOID();
1007 }
1008 
1009 Datum
1011 {
1012  float4 arg1 = PG_GETARG_FLOAT4(0);
1013  float8 arg2 = PG_GETARG_FLOAT8(1);
1014 
1015  /* widen float4 to float8 and then compare */
1016  PG_RETURN_INT32(float8_cmp_internal(arg1, arg2));
1017 }
1018 
1019 Datum
1021 {
1022  float8 arg1 = PG_GETARG_FLOAT8(0);
1023  float4 arg2 = PG_GETARG_FLOAT4(1);
1024 
1025  /* widen float4 to float8 and then compare */
1026  PG_RETURN_INT32(float8_cmp_internal(arg1, arg2));
1027 }
1028 
1029 /*
1030  * in_range support function for float8.
1031  *
1032  * Note: we needn't supply a float8_float4 variant, as implicit coercion
1033  * of the offset value takes care of that scenario just as well.
1034  */
1035 Datum
1037 {
1039  float8 base = PG_GETARG_FLOAT8(1);
1040  float8 offset = PG_GETARG_FLOAT8(2);
1041  bool sub = PG_GETARG_BOOL(3);
1042  bool less = PG_GETARG_BOOL(4);
1043  float8 sum;
1044 
1045  /*
1046  * Reject negative or NaN offset. Negative is per spec, and NaN is
1047  * because appropriate semantics for that seem non-obvious.
1048  */
1049  if (isnan(offset) || offset < 0)
1050  ereport(ERROR,
1051  (errcode(ERRCODE_INVALID_PRECEDING_OR_FOLLOWING_SIZE),
1052  errmsg("invalid preceding or following size in window function")));
1053 
1054  /*
1055  * Deal with cases where val and/or base is NaN, following the rule that
1056  * NaN sorts after non-NaN (cf float8_cmp_internal). The offset cannot
1057  * affect the conclusion.
1058  */
1059  if (isnan(val))
1060  {
1061  if (isnan(base))
1062  PG_RETURN_BOOL(true); /* NAN = NAN */
1063  else
1064  PG_RETURN_BOOL(!less); /* NAN > non-NAN */
1065  }
1066  else if (isnan(base))
1067  {
1068  PG_RETURN_BOOL(less); /* non-NAN < NAN */
1069  }
1070 
1071  /*
1072  * Deal with infinite offset (necessarily +inf, at this point). We must
1073  * special-case this because if base happens to be -inf, their sum would
1074  * be NaN, which is an overflow-ish condition we should avoid.
1075  */
1076  if (isinf(offset))
1077  {
1078  PG_RETURN_BOOL(sub ? !less : less);
1079  }
1080 
1081  /*
1082  * Otherwise it should be safe to compute base +/- offset. We trust the
1083  * FPU to cope if base is +/-inf or the true sum would overflow, and
1084  * produce a suitably signed infinity, which will compare properly against
1085  * val whether or not that's infinity.
1086  */
1087  if (sub)
1088  sum = base - offset;
1089  else
1090  sum = base + offset;
1091 
1092  if (less)
1093  PG_RETURN_BOOL(val <= sum);
1094  else
1095  PG_RETURN_BOOL(val >= sum);
1096 }
1097 
1098 /*
1099  * in_range support function for float4.
1100  *
1101  * We would need a float4_float8 variant in any case, so we supply that and
1102  * let implicit coercion take care of the float4_float4 case.
1103  */
1104 Datum
1106 {
1108  float4 base = PG_GETARG_FLOAT4(1);
1109  float8 offset = PG_GETARG_FLOAT8(2);
1110  bool sub = PG_GETARG_BOOL(3);
1111  bool less = PG_GETARG_BOOL(4);
1112  float8 sum;
1113 
1114  /*
1115  * Reject negative or NaN offset. Negative is per spec, and NaN is
1116  * because appropriate semantics for that seem non-obvious.
1117  */
1118  if (isnan(offset) || offset < 0)
1119  ereport(ERROR,
1120  (errcode(ERRCODE_INVALID_PRECEDING_OR_FOLLOWING_SIZE),
1121  errmsg("invalid preceding or following size in window function")));
1122 
1123  /*
1124  * Deal with cases where val and/or base is NaN, following the rule that
1125  * NaN sorts after non-NaN (cf float8_cmp_internal). The offset cannot
1126  * affect the conclusion.
1127  */
1128  if (isnan(val))
1129  {
1130  if (isnan(base))
1131  PG_RETURN_BOOL(true); /* NAN = NAN */
1132  else
1133  PG_RETURN_BOOL(!less); /* NAN > non-NAN */
1134  }
1135  else if (isnan(base))
1136  {
1137  PG_RETURN_BOOL(less); /* non-NAN < NAN */
1138  }
1139 
1140  /*
1141  * Deal with infinite offset (necessarily +inf, at this point). We must
1142  * special-case this because if base happens to be -inf, their sum would
1143  * be NaN, which is an overflow-ish condition we should avoid.
1144  */
1145  if (isinf(offset))
1146  {
1147  PG_RETURN_BOOL(sub ? !less : less);
1148  }
1149 
1150  /*
1151  * Otherwise it should be safe to compute base +/- offset. We trust the
1152  * FPU to cope if base is +/-inf or the true sum would overflow, and
1153  * produce a suitably signed infinity, which will compare properly against
1154  * val whether or not that's infinity.
1155  */
1156  if (sub)
1157  sum = base - offset;
1158  else
1159  sum = base + offset;
1160 
1161  if (less)
1162  PG_RETURN_BOOL(val <= sum);
1163  else
1164  PG_RETURN_BOOL(val >= sum);
1165 }
1166 
1167 
1168 /*
1169  * ===================
1170  * CONVERSION ROUTINES
1171  * ===================
1172  */
1173 
1174 /*
1175  * ftod - converts a float4 number to a float8 number
1176  */
1177 Datum
1179 {
1180  float4 num = PG_GETARG_FLOAT4(0);
1181 
1182  PG_RETURN_FLOAT8((float8) num);
1183 }
1184 
1185 
1186 /*
1187  * dtof - converts a float8 number to a float4 number
1188  */
1189 Datum
1191 {
1192  float8 num = PG_GETARG_FLOAT8(0);
1193 
1194  check_float4_val((float4) num, isinf(num), num == 0);
1195 
1196  PG_RETURN_FLOAT4((float4) num);
1197 }
1198 
1199 
1200 /*
1201  * dtoi4 - converts a float8 number to an int4 number
1202  */
1203 Datum
1205 {
1206  float8 num = PG_GETARG_FLOAT8(0);
1207 
1208  /*
1209  * Get rid of any fractional part in the input. This is so we don't fail
1210  * on just-out-of-range values that would round into range. Note
1211  * assumption that rint() will pass through a NaN or Inf unchanged.
1212  */
1213  num = rint(num);
1214 
1215  /*
1216  * Range check. We must be careful here that the boundary values are
1217  * expressed exactly in the float domain. We expect PG_INT32_MIN to be an
1218  * exact power of 2, so it will be represented exactly; but PG_INT32_MAX
1219  * isn't, and might get rounded off, so avoid using it.
1220  */
1221  if (unlikely(num < (float8) PG_INT32_MIN ||
1222  num >= -((float8) PG_INT32_MIN) ||
1223  isnan(num)))
1224  ereport(ERROR,
1225  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1226  errmsg("integer out of range")));
1227 
1228  PG_RETURN_INT32((int32) num);
1229 }
1230 
1231 
1232 /*
1233  * dtoi2 - converts a float8 number to an int2 number
1234  */
1235 Datum
1237 {
1238  float8 num = PG_GETARG_FLOAT8(0);
1239 
1240  /*
1241  * Get rid of any fractional part in the input. This is so we don't fail
1242  * on just-out-of-range values that would round into range. Note
1243  * assumption that rint() will pass through a NaN or Inf unchanged.
1244  */
1245  num = rint(num);
1246 
1247  /*
1248  * Range check. We must be careful here that the boundary values are
1249  * expressed exactly in the float domain. We expect PG_INT16_MIN to be an
1250  * exact power of 2, so it will be represented exactly; but PG_INT16_MAX
1251  * isn't, and might get rounded off, so avoid using it.
1252  */
1253  if (unlikely(num < (float8) PG_INT16_MIN ||
1254  num >= -((float8) PG_INT16_MIN) ||
1255  isnan(num)))
1256  ereport(ERROR,
1257  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1258  errmsg("smallint out of range")));
1259 
1260  PG_RETURN_INT16((int16) num);
1261 }
1262 
1263 
1264 /*
1265  * i4tod - converts an int4 number to a float8 number
1266  */
1267 Datum
1269 {
1270  int32 num = PG_GETARG_INT32(0);
1271 
1272  PG_RETURN_FLOAT8((float8) num);
1273 }
1274 
1275 
1276 /*
1277  * i2tod - converts an int2 number to a float8 number
1278  */
1279 Datum
1281 {
1282  int16 num = PG_GETARG_INT16(0);
1283 
1284  PG_RETURN_FLOAT8((float8) num);
1285 }
1286 
1287 
1288 /*
1289  * ftoi4 - converts a float4 number to an int4 number
1290  */
1291 Datum
1293 {
1294  float4 num = PG_GETARG_FLOAT4(0);
1295 
1296  /*
1297  * Get rid of any fractional part in the input. This is so we don't fail
1298  * on just-out-of-range values that would round into range. Note
1299  * assumption that rint() will pass through a NaN or Inf unchanged.
1300  */
1301  num = rint(num);
1302 
1303  /*
1304  * Range check. We must be careful here that the boundary values are
1305  * expressed exactly in the float domain. We expect PG_INT32_MIN to be an
1306  * exact power of 2, so it will be represented exactly; but PG_INT32_MAX
1307  * isn't, and might get rounded off, so avoid using it.
1308  */
1309  if (unlikely(num < (float4) PG_INT32_MIN ||
1310  num >= -((float4) PG_INT32_MIN) ||
1311  isnan(num)))
1312  ereport(ERROR,
1313  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1314  errmsg("integer out of range")));
1315 
1316  PG_RETURN_INT32((int32) num);
1317 }
1318 
1319 
1320 /*
1321  * ftoi2 - converts a float4 number to an int2 number
1322  */
1323 Datum
1325 {
1326  float4 num = PG_GETARG_FLOAT4(0);
1327 
1328  /*
1329  * Get rid of any fractional part in the input. This is so we don't fail
1330  * on just-out-of-range values that would round into range. Note
1331  * assumption that rint() will pass through a NaN or Inf unchanged.
1332  */
1333  num = rint(num);
1334 
1335  /*
1336  * Range check. We must be careful here that the boundary values are
1337  * expressed exactly in the float domain. We expect PG_INT16_MIN to be an
1338  * exact power of 2, so it will be represented exactly; but PG_INT16_MAX
1339  * isn't, and might get rounded off, so avoid using it.
1340  */
1341  if (unlikely(num < (float4) PG_INT16_MIN ||
1342  num >= -((float4) PG_INT16_MIN) ||
1343  isnan(num)))
1344  ereport(ERROR,
1345  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1346  errmsg("smallint out of range")));
1347 
1348  PG_RETURN_INT16((int16) num);
1349 }
1350 
1351 
1352 /*
1353  * i4tof - converts an int4 number to a float4 number
1354  */
1355 Datum
1357 {
1358  int32 num = PG_GETARG_INT32(0);
1359 
1360  PG_RETURN_FLOAT4((float4) num);
1361 }
1362 
1363 
1364 /*
1365  * i2tof - converts an int2 number to a float4 number
1366  */
1367 Datum
1369 {
1370  int16 num = PG_GETARG_INT16(0);
1371 
1372  PG_RETURN_FLOAT4((float4) num);
1373 }
1374 
1375 
1376 /*
1377  * =======================
1378  * RANDOM FLOAT8 OPERATORS
1379  * =======================
1380  */
1381 
1382 /*
1383  * dround - returns ROUND(arg1)
1384  */
1385 Datum
1387 {
1388  float8 arg1 = PG_GETARG_FLOAT8(0);
1389 
1390  PG_RETURN_FLOAT8(rint(arg1));
1391 }
1392 
1393 /*
1394  * dceil - returns the smallest integer greater than or
1395  * equal to the specified float
1396  */
1397 Datum
1399 {
1400  float8 arg1 = PG_GETARG_FLOAT8(0);
1401 
1402  PG_RETURN_FLOAT8(ceil(arg1));
1403 }
1404 
1405 /*
1406  * dfloor - returns the largest integer lesser than or
1407  * equal to the specified float
1408  */
1409 Datum
1411 {
1412  float8 arg1 = PG_GETARG_FLOAT8(0);
1413 
1414  PG_RETURN_FLOAT8(floor(arg1));
1415 }
1416 
1417 /*
1418  * dsign - returns -1 if the argument is less than 0, 0
1419  * if the argument is equal to 0, and 1 if the
1420  * argument is greater than zero.
1421  */
1422 Datum
1424 {
1425  float8 arg1 = PG_GETARG_FLOAT8(0);
1426  float8 result;
1427 
1428  if (arg1 > 0)
1429  result = 1.0;
1430  else if (arg1 < 0)
1431  result = -1.0;
1432  else
1433  result = 0.0;
1434 
1435  PG_RETURN_FLOAT8(result);
1436 }
1437 
1438 /*
1439  * dtrunc - returns truncation-towards-zero of arg1,
1440  * arg1 >= 0 ... the greatest integer less
1441  * than or equal to arg1
1442  * arg1 < 0 ... the least integer greater
1443  * than or equal to arg1
1444  */
1445 Datum
1447 {
1448  float8 arg1 = PG_GETARG_FLOAT8(0);
1449  float8 result;
1450 
1451  if (arg1 >= 0)
1452  result = floor(arg1);
1453  else
1454  result = -floor(-arg1);
1455 
1456  PG_RETURN_FLOAT8(result);
1457 }
1458 
1459 
1460 /*
1461  * dsqrt - returns square root of arg1
1462  */
1463 Datum
1465 {
1466  float8 arg1 = PG_GETARG_FLOAT8(0);
1467  float8 result;
1468 
1469  if (arg1 < 0)
1470  ereport(ERROR,
1471  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_POWER_FUNCTION),
1472  errmsg("cannot take square root of a negative number")));
1473 
1474  result = sqrt(arg1);
1475 
1476  check_float8_val(result, isinf(arg1), arg1 == 0);
1477  PG_RETURN_FLOAT8(result);
1478 }
1479 
1480 
1481 /*
1482  * dcbrt - returns cube root of arg1
1483  */
1484 Datum
1486 {
1487  float8 arg1 = PG_GETARG_FLOAT8(0);
1488  float8 result;
1489 
1490  result = cbrt(arg1);
1491  check_float8_val(result, isinf(arg1), arg1 == 0);
1492  PG_RETURN_FLOAT8(result);
1493 }
1494 
1495 
1496 /*
1497  * dpow - returns pow(arg1,arg2)
1498  */
1499 Datum
1501 {
1502  float8 arg1 = PG_GETARG_FLOAT8(0);
1503  float8 arg2 = PG_GETARG_FLOAT8(1);
1504  float8 result;
1505 
1506  /*
1507  * The POSIX spec says that NaN ^ 0 = 1, and 1 ^ NaN = 1, while all other
1508  * cases with NaN inputs yield NaN (with no error). Many older platforms
1509  * get one or more of these cases wrong, so deal with them via explicit
1510  * logic rather than trusting pow(3).
1511  */
1512  if (isnan(arg1))
1513  {
1514  if (isnan(arg2) || arg2 != 0.0)
1516  PG_RETURN_FLOAT8(1.0);
1517  }
1518  if (isnan(arg2))
1519  {
1520  if (arg1 != 1.0)
1522  PG_RETURN_FLOAT8(1.0);
1523  }
1524 
1525  /*
1526  * The SQL spec requires that we emit a particular SQLSTATE error code for
1527  * certain error conditions. Specifically, we don't return a
1528  * divide-by-zero error code for 0 ^ -1.
1529  */
1530  if (arg1 == 0 && arg2 < 0)
1531  ereport(ERROR,
1532  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_POWER_FUNCTION),
1533  errmsg("zero raised to a negative power is undefined")));
1534  if (arg1 < 0 && floor(arg2) != arg2)
1535  ereport(ERROR,
1536  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_POWER_FUNCTION),
1537  errmsg("a negative number raised to a non-integer power yields a complex result")));
1538 
1539  /*
1540  * pow() sets errno only on some platforms, depending on whether it
1541  * follows _IEEE_, _POSIX_, _XOPEN_, or _SVID_, so we try to avoid using
1542  * errno. However, some platform/CPU combinations return errno == EDOM
1543  * and result == NaN for negative arg1 and very large arg2 (they must be
1544  * using something different from our floor() test to decide it's
1545  * invalid). Other platforms (HPPA) return errno == ERANGE and a large
1546  * (HUGE_VAL) but finite result to signal overflow.
1547  */
1548  errno = 0;
1549  result = pow(arg1, arg2);
1550  if (errno == EDOM && isnan(result))
1551  {
1552  if ((fabs(arg1) > 1 && arg2 >= 0) || (fabs(arg1) < 1 && arg2 < 0))
1553  /* The sign of Inf is not significant in this case. */
1554  result = get_float8_infinity();
1555  else if (fabs(arg1) != 1)
1556  result = 0;
1557  else
1558  result = 1;
1559  }
1560  else if (errno == ERANGE && result != 0 && !isinf(result))
1561  result = get_float8_infinity();
1562 
1563  check_float8_val(result, isinf(arg1) || isinf(arg2), arg1 == 0);
1564  PG_RETURN_FLOAT8(result);
1565 }
1566 
1567 
1568 /*
1569  * dexp - returns the exponential function of arg1
1570  */
1571 Datum
1573 {
1574  float8 arg1 = PG_GETARG_FLOAT8(0);
1575  float8 result;
1576 
1577  errno = 0;
1578  result = exp(arg1);
1579  if (errno == ERANGE && result != 0 && !isinf(result))
1580  result = get_float8_infinity();
1581 
1582  check_float8_val(result, isinf(arg1), false);
1583  PG_RETURN_FLOAT8(result);
1584 }
1585 
1586 
1587 /*
1588  * dlog1 - returns the natural logarithm of arg1
1589  */
1590 Datum
1592 {
1593  float8 arg1 = PG_GETARG_FLOAT8(0);
1594  float8 result;
1595 
1596  /*
1597  * Emit particular SQLSTATE error codes for ln(). This is required by the
1598  * SQL standard.
1599  */
1600  if (arg1 == 0.0)
1601  ereport(ERROR,
1602  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_LOG),
1603  errmsg("cannot take logarithm of zero")));
1604  if (arg1 < 0)
1605  ereport(ERROR,
1606  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_LOG),
1607  errmsg("cannot take logarithm of a negative number")));
1608 
1609  result = log(arg1);
1610 
1611  check_float8_val(result, isinf(arg1), arg1 == 1);
1612  PG_RETURN_FLOAT8(result);
1613 }
1614 
1615 
1616 /*
1617  * dlog10 - returns the base 10 logarithm of arg1
1618  */
1619 Datum
1621 {
1622  float8 arg1 = PG_GETARG_FLOAT8(0);
1623  float8 result;
1624 
1625  /*
1626  * Emit particular SQLSTATE error codes for log(). The SQL spec doesn't
1627  * define log(), but it does define ln(), so it makes sense to emit the
1628  * same error code for an analogous error condition.
1629  */
1630  if (arg1 == 0.0)
1631  ereport(ERROR,
1632  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_LOG),
1633  errmsg("cannot take logarithm of zero")));
1634  if (arg1 < 0)
1635  ereport(ERROR,
1636  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_LOG),
1637  errmsg("cannot take logarithm of a negative number")));
1638 
1639  result = log10(arg1);
1640 
1641  check_float8_val(result, isinf(arg1), arg1 == 1);
1642  PG_RETURN_FLOAT8(result);
1643 }
1644 
1645 
1646 /*
1647  * dacos - returns the arccos of arg1 (radians)
1648  */
1649 Datum
1651 {
1652  float8 arg1 = PG_GETARG_FLOAT8(0);
1653  float8 result;
1654 
1655  /* Per the POSIX spec, return NaN if the input is NaN */
1656  if (isnan(arg1))
1658 
1659  /*
1660  * The principal branch of the inverse cosine function maps values in the
1661  * range [-1, 1] to values in the range [0, Pi], so we should reject any
1662  * inputs outside that range and the result will always be finite.
1663  */
1664  if (arg1 < -1.0 || arg1 > 1.0)
1665  ereport(ERROR,
1666  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1667  errmsg("input is out of range")));
1668 
1669  result = acos(arg1);
1670 
1671  check_float8_val(result, false, true);
1672  PG_RETURN_FLOAT8(result);
1673 }
1674 
1675 
1676 /*
1677  * dasin - returns the arcsin of arg1 (radians)
1678  */
1679 Datum
1681 {
1682  float8 arg1 = PG_GETARG_FLOAT8(0);
1683  float8 result;
1684 
1685  /* Per the POSIX spec, return NaN if the input is NaN */
1686  if (isnan(arg1))
1688 
1689  /*
1690  * The principal branch of the inverse sine function maps values in the
1691  * range [-1, 1] to values in the range [-Pi/2, Pi/2], so we should reject
1692  * any inputs outside that range and the result will always be finite.
1693  */
1694  if (arg1 < -1.0 || arg1 > 1.0)
1695  ereport(ERROR,
1696  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1697  errmsg("input is out of range")));
1698 
1699  result = asin(arg1);
1700 
1701  check_float8_val(result, false, true);
1702  PG_RETURN_FLOAT8(result);
1703 }
1704 
1705 
1706 /*
1707  * datan - returns the arctan of arg1 (radians)
1708  */
1709 Datum
1711 {
1712  float8 arg1 = PG_GETARG_FLOAT8(0);
1713  float8 result;
1714 
1715  /* Per the POSIX spec, return NaN if the input is NaN */
1716  if (isnan(arg1))
1718 
1719  /*
1720  * The principal branch of the inverse tangent function maps all inputs to
1721  * values in the range [-Pi/2, Pi/2], so the result should always be
1722  * finite, even if the input is infinite.
1723  */
1724  result = atan(arg1);
1725 
1726  check_float8_val(result, false, true);
1727  PG_RETURN_FLOAT8(result);
1728 }
1729 
1730 
1731 /*
1732  * atan2 - returns the arctan of arg1/arg2 (radians)
1733  */
1734 Datum
1736 {
1737  float8 arg1 = PG_GETARG_FLOAT8(0);
1738  float8 arg2 = PG_GETARG_FLOAT8(1);
1739  float8 result;
1740 
1741  /* Per the POSIX spec, return NaN if either input is NaN */
1742  if (isnan(arg1) || isnan(arg2))
1744 
1745  /*
1746  * atan2 maps all inputs to values in the range [-Pi, Pi], so the result
1747  * should always be finite, even if the inputs are infinite.
1748  */
1749  result = atan2(arg1, arg2);
1750 
1751  check_float8_val(result, false, true);
1752  PG_RETURN_FLOAT8(result);
1753 }
1754 
1755 
1756 /*
1757  * dcos - returns the cosine of arg1 (radians)
1758  */
1759 Datum
1761 {
1762  float8 arg1 = PG_GETARG_FLOAT8(0);
1763  float8 result;
1764 
1765  /* Per the POSIX spec, return NaN if the input is NaN */
1766  if (isnan(arg1))
1768 
1769  /*
1770  * cos() is periodic and so theoretically can work for all finite inputs,
1771  * but some implementations may choose to throw error if the input is so
1772  * large that there are no significant digits in the result. So we should
1773  * check for errors. POSIX allows an error to be reported either via
1774  * errno or via fetestexcept(), but currently we only support checking
1775  * errno. (fetestexcept() is rumored to report underflow unreasonably
1776  * early on some platforms, so it's not clear that believing it would be a
1777  * net improvement anyway.)
1778  *
1779  * For infinite inputs, POSIX specifies that the trigonometric functions
1780  * should return a domain error; but we won't notice that unless the
1781  * platform reports via errno, so also explicitly test for infinite
1782  * inputs.
1783  */
1784  errno = 0;
1785  result = cos(arg1);
1786  if (errno != 0 || isinf(arg1))
1787  ereport(ERROR,
1788  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1789  errmsg("input is out of range")));
1790 
1791  check_float8_val(result, false, true);
1792  PG_RETURN_FLOAT8(result);
1793 }
1794 
1795 
1796 /*
1797  * dcot - returns the cotangent of arg1 (radians)
1798  */
1799 Datum
1801 {
1802  float8 arg1 = PG_GETARG_FLOAT8(0);
1803  float8 result;
1804 
1805  /* Per the POSIX spec, return NaN if the input is NaN */
1806  if (isnan(arg1))
1808 
1809  /* Be sure to throw an error if the input is infinite --- see dcos() */
1810  errno = 0;
1811  result = tan(arg1);
1812  if (errno != 0 || isinf(arg1))
1813  ereport(ERROR,
1814  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1815  errmsg("input is out of range")));
1816 
1817  result = 1.0 / result;
1818  check_float8_val(result, true /* cot(0) == Inf */ , true);
1819  PG_RETURN_FLOAT8(result);
1820 }
1821 
1822 
1823 /*
1824  * dsin - returns the sine of arg1 (radians)
1825  */
1826 Datum
1828 {
1829  float8 arg1 = PG_GETARG_FLOAT8(0);
1830  float8 result;
1831 
1832  /* Per the POSIX spec, return NaN if the input is NaN */
1833  if (isnan(arg1))
1835 
1836  /* Be sure to throw an error if the input is infinite --- see dcos() */
1837  errno = 0;
1838  result = sin(arg1);
1839  if (errno != 0 || isinf(arg1))
1840  ereport(ERROR,
1841  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1842  errmsg("input is out of range")));
1843 
1844  check_float8_val(result, false, true);
1845  PG_RETURN_FLOAT8(result);
1846 }
1847 
1848 
1849 /*
1850  * dtan - returns the tangent of arg1 (radians)
1851  */
1852 Datum
1854 {
1855  float8 arg1 = PG_GETARG_FLOAT8(0);
1856  float8 result;
1857 
1858  /* Per the POSIX spec, return NaN if the input is NaN */
1859  if (isnan(arg1))
1861 
1862  /* Be sure to throw an error if the input is infinite --- see dcos() */
1863  errno = 0;
1864  result = tan(arg1);
1865  if (errno != 0 || isinf(arg1))
1866  ereport(ERROR,
1867  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
1868  errmsg("input is out of range")));
1869 
1870  check_float8_val(result, true /* tan(pi/2) == Inf */ , true);
1871  PG_RETURN_FLOAT8(result);
1872 }
1873 
1874 
1875 /* ========== DEGREE-BASED TRIGONOMETRIC FUNCTIONS ========== */
1876 
1877 
1878 /*
1879  * Initialize the cached constants declared at the head of this file
1880  * (sin_30 etc). The fact that we need those at all, let alone need this
1881  * Rube-Goldberg-worthy method of initializing them, is because there are
1882  * compilers out there that will precompute expressions such as sin(constant)
1883  * using a sin() function different from what will be used at runtime. If we
1884  * want exact results, we must ensure that none of the scaling constants used
1885  * in the degree-based trig functions are computed that way. To do so, we
1886  * compute them from the variables degree_c_thirty etc, which are also really
1887  * constants, but the compiler cannot assume that.
1888  *
1889  * Other hazards we are trying to forestall with this kluge include the
1890  * possibility that compilers will rearrange the expressions, or compute
1891  * some intermediate results in registers wider than a standard double.
1892  *
1893  * In the places where we use these constants, the typical pattern is like
1894  * volatile float8 sin_x = sin(x * RADIANS_PER_DEGREE);
1895  * return (sin_x / sin_30);
1896  * where we hope to get a value of exactly 1.0 from the division when x = 30.
1897  * The volatile temporary variable is needed on machines with wide float
1898  * registers, to ensure that the result of sin(x) is rounded to double width
1899  * the same as the value of sin_30 has been. Experimentation with gcc shows
1900  * that marking the temp variable volatile is necessary to make the store and
1901  * reload actually happen; hopefully the same trick works for other compilers.
1902  * (gcc's documentation suggests using the -ffloat-store compiler switch to
1903  * ensure this, but that is compiler-specific and it also pessimizes code in
1904  * many places where we don't care about this.)
1905  */
1906 static void
1908 {
1910  one_minus_cos_60 = 1.0 - cos(degree_c_sixty * RADIANS_PER_DEGREE);
1911  asin_0_5 = asin(degree_c_one_half);
1912  acos_0_5 = acos(degree_c_one_half);
1913  atan_1_0 = atan(degree_c_one);
1916  degree_consts_set = true;
1917 }
1918 
1919 #define INIT_DEGREE_CONSTANTS() \
1920 do { \
1921  if (!degree_consts_set) \
1922  init_degree_constants(); \
1923 } while(0)
1924 
1925 
1926 /*
1927  * asind_q1 - returns the inverse sine of x in degrees, for x in
1928  * the range [0, 1]. The result is an angle in the
1929  * first quadrant --- [0, 90] degrees.
1930  *
1931  * For the 3 special case inputs (0, 0.5 and 1), this
1932  * function will return exact values (0, 30 and 90
1933  * degrees respectively).
1934  */
1935 static double
1936 asind_q1(double x)
1937 {
1938  /*
1939  * Stitch together inverse sine and cosine functions for the ranges [0,
1940  * 0.5] and (0.5, 1]. Each expression below is guaranteed to return
1941  * exactly 30 for x=0.5, so the result is a continuous monotonic function
1942  * over the full range.
1943  */
1944  if (x <= 0.5)
1945  {
1946  volatile float8 asin_x = asin(x);
1947 
1948  return (asin_x / asin_0_5) * 30.0;
1949  }
1950  else
1951  {
1952  volatile float8 acos_x = acos(x);
1953 
1954  return 90.0 - (acos_x / acos_0_5) * 60.0;
1955  }
1956 }
1957 
1958 
1959 /*
1960  * acosd_q1 - returns the inverse cosine of x in degrees, for x in
1961  * the range [0, 1]. The result is an angle in the
1962  * first quadrant --- [0, 90] degrees.
1963  *
1964  * For the 3 special case inputs (0, 0.5 and 1), this
1965  * function will return exact values (0, 60 and 90
1966  * degrees respectively).
1967  */
1968 static double
1969 acosd_q1(double x)
1970 {
1971  /*
1972  * Stitch together inverse sine and cosine functions for the ranges [0,
1973  * 0.5] and (0.5, 1]. Each expression below is guaranteed to return
1974  * exactly 60 for x=0.5, so the result is a continuous monotonic function
1975  * over the full range.
1976  */
1977  if (x <= 0.5)
1978  {
1979  volatile float8 asin_x = asin(x);
1980 
1981  return 90.0 - (asin_x / asin_0_5) * 30.0;
1982  }
1983  else
1984  {
1985  volatile float8 acos_x = acos(x);
1986 
1987  return (acos_x / acos_0_5) * 60.0;
1988  }
1989 }
1990 
1991 
1992 /*
1993  * dacosd - returns the arccos of arg1 (degrees)
1994  */
1995 Datum
1997 {
1998  float8 arg1 = PG_GETARG_FLOAT8(0);
1999  float8 result;
2000 
2001  /* Per the POSIX spec, return NaN if the input is NaN */
2002  if (isnan(arg1))
2004 
2006 
2007  /*
2008  * The principal branch of the inverse cosine function maps values in the
2009  * range [-1, 1] to values in the range [0, 180], so we should reject any
2010  * inputs outside that range and the result will always be finite.
2011  */
2012  if (arg1 < -1.0 || arg1 > 1.0)
2013  ereport(ERROR,
2014  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2015  errmsg("input is out of range")));
2016 
2017  if (arg1 >= 0.0)
2018  result = acosd_q1(arg1);
2019  else
2020  result = 90.0 + asind_q1(-arg1);
2021 
2022  check_float8_val(result, false, true);
2023  PG_RETURN_FLOAT8(result);
2024 }
2025 
2026 
2027 /*
2028  * dasind - returns the arcsin of arg1 (degrees)
2029  */
2030 Datum
2032 {
2033  float8 arg1 = PG_GETARG_FLOAT8(0);
2034  float8 result;
2035 
2036  /* Per the POSIX spec, return NaN if the input is NaN */
2037  if (isnan(arg1))
2039 
2041 
2042  /*
2043  * The principal branch of the inverse sine function maps values in the
2044  * range [-1, 1] to values in the range [-90, 90], so we should reject any
2045  * inputs outside that range and the result will always be finite.
2046  */
2047  if (arg1 < -1.0 || arg1 > 1.0)
2048  ereport(ERROR,
2049  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2050  errmsg("input is out of range")));
2051 
2052  if (arg1 >= 0.0)
2053  result = asind_q1(arg1);
2054  else
2055  result = -asind_q1(-arg1);
2056 
2057  check_float8_val(result, false, true);
2058  PG_RETURN_FLOAT8(result);
2059 }
2060 
2061 
2062 /*
2063  * datand - returns the arctan of arg1 (degrees)
2064  */
2065 Datum
2067 {
2068  float8 arg1 = PG_GETARG_FLOAT8(0);
2069  float8 result;
2070  volatile float8 atan_arg1;
2071 
2072  /* Per the POSIX spec, return NaN if the input is NaN */
2073  if (isnan(arg1))
2075 
2077 
2078  /*
2079  * The principal branch of the inverse tangent function maps all inputs to
2080  * values in the range [-90, 90], so the result should always be finite,
2081  * even if the input is infinite. Additionally, we take care to ensure
2082  * than when arg1 is 1, the result is exactly 45.
2083  */
2084  atan_arg1 = atan(arg1);
2085  result = (atan_arg1 / atan_1_0) * 45.0;
2086 
2087  check_float8_val(result, false, true);
2088  PG_RETURN_FLOAT8(result);
2089 }
2090 
2091 
2092 /*
2093  * atan2d - returns the arctan of arg1/arg2 (degrees)
2094  */
2095 Datum
2097 {
2098  float8 arg1 = PG_GETARG_FLOAT8(0);
2099  float8 arg2 = PG_GETARG_FLOAT8(1);
2100  float8 result;
2101  volatile float8 atan2_arg1_arg2;
2102 
2103  /* Per the POSIX spec, return NaN if either input is NaN */
2104  if (isnan(arg1) || isnan(arg2))
2106 
2108 
2109  /*
2110  * atan2d maps all inputs to values in the range [-180, 180], so the
2111  * result should always be finite, even if the inputs are infinite.
2112  *
2113  * Note: this coding assumes that atan(1.0) is a suitable scaling constant
2114  * to get an exact result from atan2(). This might well fail on us at
2115  * some point, requiring us to decide exactly what inputs we think we're
2116  * going to guarantee an exact result for.
2117  */
2118  atan2_arg1_arg2 = atan2(arg1, arg2);
2119  result = (atan2_arg1_arg2 / atan_1_0) * 45.0;
2120 
2121  check_float8_val(result, false, true);
2122  PG_RETURN_FLOAT8(result);
2123 }
2124 
2125 
2126 /*
2127  * sind_0_to_30 - returns the sine of an angle that lies between 0 and
2128  * 30 degrees. This will return exactly 0 when x is 0,
2129  * and exactly 0.5 when x is 30 degrees.
2130  */
2131 static double
2132 sind_0_to_30(double x)
2133 {
2134  volatile float8 sin_x = sin(x * RADIANS_PER_DEGREE);
2135 
2136  return (sin_x / sin_30) / 2.0;
2137 }
2138 
2139 
2140 /*
2141  * cosd_0_to_60 - returns the cosine of an angle that lies between 0
2142  * and 60 degrees. This will return exactly 1 when x
2143  * is 0, and exactly 0.5 when x is 60 degrees.
2144  */
2145 static double
2146 cosd_0_to_60(double x)
2147 {
2148  volatile float8 one_minus_cos_x = 1.0 - cos(x * RADIANS_PER_DEGREE);
2149 
2150  return 1.0 - (one_minus_cos_x / one_minus_cos_60) / 2.0;
2151 }
2152 
2153 
2154 /*
2155  * sind_q1 - returns the sine of an angle in the first quadrant
2156  * (0 to 90 degrees).
2157  */
2158 static double
2159 sind_q1(double x)
2160 {
2161  /*
2162  * Stitch together the sine and cosine functions for the ranges [0, 30]
2163  * and (30, 90]. These guarantee to return exact answers at their
2164  * endpoints, so the overall result is a continuous monotonic function
2165  * that gives exact results when x = 0, 30 and 90 degrees.
2166  */
2167  if (x <= 30.0)
2168  return sind_0_to_30(x);
2169  else
2170  return cosd_0_to_60(90.0 - x);
2171 }
2172 
2173 
2174 /*
2175  * cosd_q1 - returns the cosine of an angle in the first quadrant
2176  * (0 to 90 degrees).
2177  */
2178 static double
2179 cosd_q1(double x)
2180 {
2181  /*
2182  * Stitch together the sine and cosine functions for the ranges [0, 60]
2183  * and (60, 90]. These guarantee to return exact answers at their
2184  * endpoints, so the overall result is a continuous monotonic function
2185  * that gives exact results when x = 0, 60 and 90 degrees.
2186  */
2187  if (x <= 60.0)
2188  return cosd_0_to_60(x);
2189  else
2190  return sind_0_to_30(90.0 - x);
2191 }
2192 
2193 
2194 /*
2195  * dcosd - returns the cosine of arg1 (degrees)
2196  */
2197 Datum
2199 {
2200  float8 arg1 = PG_GETARG_FLOAT8(0);
2201  float8 result;
2202  int sign = 1;
2203 
2204  /*
2205  * Per the POSIX spec, return NaN if the input is NaN and throw an error
2206  * if the input is infinite.
2207  */
2208  if (isnan(arg1))
2210 
2211  if (isinf(arg1))
2212  ereport(ERROR,
2213  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2214  errmsg("input is out of range")));
2215 
2217 
2218  /* Reduce the range of the input to [0,90] degrees */
2219  arg1 = fmod(arg1, 360.0);
2220 
2221  if (arg1 < 0.0)
2222  {
2223  /* cosd(-x) = cosd(x) */
2224  arg1 = -arg1;
2225  }
2226 
2227  if (arg1 > 180.0)
2228  {
2229  /* cosd(360-x) = cosd(x) */
2230  arg1 = 360.0 - arg1;
2231  }
2232 
2233  if (arg1 > 90.0)
2234  {
2235  /* cosd(180-x) = -cosd(x) */
2236  arg1 = 180.0 - arg1;
2237  sign = -sign;
2238  }
2239 
2240  result = sign * cosd_q1(arg1);
2241 
2242  check_float8_val(result, false, true);
2243  PG_RETURN_FLOAT8(result);
2244 }
2245 
2246 
2247 /*
2248  * dcotd - returns the cotangent of arg1 (degrees)
2249  */
2250 Datum
2252 {
2253  float8 arg1 = PG_GETARG_FLOAT8(0);
2254  float8 result;
2255  volatile float8 cot_arg1;
2256  int sign = 1;
2257 
2258  /*
2259  * Per the POSIX spec, return NaN if the input is NaN and throw an error
2260  * if the input is infinite.
2261  */
2262  if (isnan(arg1))
2264 
2265  if (isinf(arg1))
2266  ereport(ERROR,
2267  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2268  errmsg("input is out of range")));
2269 
2271 
2272  /* Reduce the range of the input to [0,90] degrees */
2273  arg1 = fmod(arg1, 360.0);
2274 
2275  if (arg1 < 0.0)
2276  {
2277  /* cotd(-x) = -cotd(x) */
2278  arg1 = -arg1;
2279  sign = -sign;
2280  }
2281 
2282  if (arg1 > 180.0)
2283  {
2284  /* cotd(360-x) = -cotd(x) */
2285  arg1 = 360.0 - arg1;
2286  sign = -sign;
2287  }
2288 
2289  if (arg1 > 90.0)
2290  {
2291  /* cotd(180-x) = -cotd(x) */
2292  arg1 = 180.0 - arg1;
2293  sign = -sign;
2294  }
2295 
2296  cot_arg1 = cosd_q1(arg1) / sind_q1(arg1);
2297  result = sign * (cot_arg1 / cot_45);
2298 
2299  /*
2300  * On some machines we get cotd(270) = minus zero, but this isn't always
2301  * true. For portability, and because the user constituency for this
2302  * function probably doesn't want minus zero, force it to plain zero.
2303  */
2304  if (result == 0.0)
2305  result = 0.0;
2306 
2307  check_float8_val(result, true /* cotd(0) == Inf */ , true);
2308  PG_RETURN_FLOAT8(result);
2309 }
2310 
2311 
2312 /*
2313  * dsind - returns the sine of arg1 (degrees)
2314  */
2315 Datum
2317 {
2318  float8 arg1 = PG_GETARG_FLOAT8(0);
2319  float8 result;
2320  int sign = 1;
2321 
2322  /*
2323  * Per the POSIX spec, return NaN if the input is NaN and throw an error
2324  * if the input is infinite.
2325  */
2326  if (isnan(arg1))
2328 
2329  if (isinf(arg1))
2330  ereport(ERROR,
2331  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2332  errmsg("input is out of range")));
2333 
2335 
2336  /* Reduce the range of the input to [0,90] degrees */
2337  arg1 = fmod(arg1, 360.0);
2338 
2339  if (arg1 < 0.0)
2340  {
2341  /* sind(-x) = -sind(x) */
2342  arg1 = -arg1;
2343  sign = -sign;
2344  }
2345 
2346  if (arg1 > 180.0)
2347  {
2348  /* sind(360-x) = -sind(x) */
2349  arg1 = 360.0 - arg1;
2350  sign = -sign;
2351  }
2352 
2353  if (arg1 > 90.0)
2354  {
2355  /* sind(180-x) = sind(x) */
2356  arg1 = 180.0 - arg1;
2357  }
2358 
2359  result = sign * sind_q1(arg1);
2360 
2361  check_float8_val(result, false, true);
2362  PG_RETURN_FLOAT8(result);
2363 }
2364 
2365 
2366 /*
2367  * dtand - returns the tangent of arg1 (degrees)
2368  */
2369 Datum
2371 {
2372  float8 arg1 = PG_GETARG_FLOAT8(0);
2373  float8 result;
2374  volatile float8 tan_arg1;
2375  int sign = 1;
2376 
2377  /*
2378  * Per the POSIX spec, return NaN if the input is NaN and throw an error
2379  * if the input is infinite.
2380  */
2381  if (isnan(arg1))
2383 
2384  if (isinf(arg1))
2385  ereport(ERROR,
2386  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2387  errmsg("input is out of range")));
2388 
2390 
2391  /* Reduce the range of the input to [0,90] degrees */
2392  arg1 = fmod(arg1, 360.0);
2393 
2394  if (arg1 < 0.0)
2395  {
2396  /* tand(-x) = -tand(x) */
2397  arg1 = -arg1;
2398  sign = -sign;
2399  }
2400 
2401  if (arg1 > 180.0)
2402  {
2403  /* tand(360-x) = -tand(x) */
2404  arg1 = 360.0 - arg1;
2405  sign = -sign;
2406  }
2407 
2408  if (arg1 > 90.0)
2409  {
2410  /* tand(180-x) = -tand(x) */
2411  arg1 = 180.0 - arg1;
2412  sign = -sign;
2413  }
2414 
2415  tan_arg1 = sind_q1(arg1) / cosd_q1(arg1);
2416  result = sign * (tan_arg1 / tan_45);
2417 
2418  /*
2419  * On some machines we get tand(180) = minus zero, but this isn't always
2420  * true. For portability, and because the user constituency for this
2421  * function probably doesn't want minus zero, force it to plain zero.
2422  */
2423  if (result == 0.0)
2424  result = 0.0;
2425 
2426  check_float8_val(result, true /* tand(90) == Inf */ , true);
2427  PG_RETURN_FLOAT8(result);
2428 }
2429 
2430 
2431 /*
2432  * degrees - returns degrees converted from radians
2433  */
2434 Datum
2436 {
2437  float8 arg1 = PG_GETARG_FLOAT8(0);
2438 
2440 }
2441 
2442 
2443 /*
2444  * dpi - returns the constant PI
2445  */
2446 Datum
2448 {
2450 }
2451 
2452 
2453 /*
2454  * radians - returns radians converted from degrees
2455  */
2456 Datum
2458 {
2459  float8 arg1 = PG_GETARG_FLOAT8(0);
2460 
2462 }
2463 
2464 
2465 /* ========== HYPERBOLIC FUNCTIONS ========== */
2466 
2467 
2468 /*
2469  * dsinh - returns the hyperbolic sine of arg1
2470  */
2471 Datum
2473 {
2474  float8 arg1 = PG_GETARG_FLOAT8(0);
2475  float8 result;
2476 
2477  errno = 0;
2478  result = sinh(arg1);
2479 
2480  /*
2481  * if an ERANGE error occurs, it means there is an overflow. For sinh,
2482  * the result should be either -infinity or infinity, depending on the
2483  * sign of arg1.
2484  */
2485  if (errno == ERANGE)
2486  {
2487  if (arg1 < 0)
2488  result = -get_float8_infinity();
2489  else
2490  result = get_float8_infinity();
2491  }
2492 
2493  check_float8_val(result, true, true);
2494  PG_RETURN_FLOAT8(result);
2495 }
2496 
2497 
2498 /*
2499  * dcosh - returns the hyperbolic cosine of arg1
2500  */
2501 Datum
2503 {
2504  float8 arg1 = PG_GETARG_FLOAT8(0);
2505  float8 result;
2506 
2507  errno = 0;
2508  result = cosh(arg1);
2509 
2510  /*
2511  * if an ERANGE error occurs, it means there is an overflow. As cosh is
2512  * always positive, it always means the result is positive infinity.
2513  */
2514  if (errno == ERANGE)
2515  result = get_float8_infinity();
2516 
2517  check_float8_val(result, true, false);
2518  PG_RETURN_FLOAT8(result);
2519 }
2520 
2521 /*
2522  * dtanh - returns the hyperbolic tangent of arg1
2523  */
2524 Datum
2526 {
2527  float8 arg1 = PG_GETARG_FLOAT8(0);
2528  float8 result;
2529 
2530  /*
2531  * For tanh, we don't need an errno check because it never overflows.
2532  */
2533  result = tanh(arg1);
2534 
2535  check_float8_val(result, false, true);
2536  PG_RETURN_FLOAT8(result);
2537 }
2538 
2539 /*
2540  * dasinh - returns the inverse hyperbolic sine of arg1
2541  */
2542 Datum
2544 {
2545  float8 arg1 = PG_GETARG_FLOAT8(0);
2546  float8 result;
2547 
2548  /*
2549  * For asinh, we don't need an errno check because it never overflows.
2550  */
2551  result = asinh(arg1);
2552 
2553  check_float8_val(result, true, true);
2554  PG_RETURN_FLOAT8(result);
2555 }
2556 
2557 /*
2558  * dacosh - returns the inverse hyperbolic cosine of arg1
2559  */
2560 Datum
2562 {
2563  float8 arg1 = PG_GETARG_FLOAT8(0);
2564  float8 result;
2565 
2566  /*
2567  * acosh is only defined for inputs >= 1.0. By checking this ourselves,
2568  * we need not worry about checking for an EDOM error, which is a good
2569  * thing because some implementations will report that for NaN. Otherwise,
2570  * no error is possible.
2571  */
2572  if (arg1 < 1.0)
2573  ereport(ERROR,
2574  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2575  errmsg("input is out of range")));
2576 
2577  result = acosh(arg1);
2578 
2579  check_float8_val(result, true, true);
2580  PG_RETURN_FLOAT8(result);
2581 }
2582 
2583 /*
2584  * datanh - returns the inverse hyperbolic tangent of arg1
2585  */
2586 Datum
2588 {
2589  float8 arg1 = PG_GETARG_FLOAT8(0);
2590  float8 result;
2591 
2592  /*
2593  * atanh is only defined for inputs between -1 and 1. By checking this
2594  * ourselves, we need not worry about checking for an EDOM error, which is
2595  * a good thing because some implementations will report that for NaN.
2596  */
2597  if (arg1 < -1.0 || arg1 > 1.0)
2598  ereport(ERROR,
2599  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2600  errmsg("input is out of range")));
2601 
2602  /*
2603  * Also handle the infinity cases ourselves; this is helpful because old
2604  * glibc versions may produce the wrong errno for this. All other inputs
2605  * cannot produce an error.
2606  */
2607  if (arg1 == -1.0)
2608  result = -get_float8_infinity();
2609  else if (arg1 == 1.0)
2610  result = get_float8_infinity();
2611  else
2612  result = atanh(arg1);
2613 
2614  check_float8_val(result, true, true);
2615  PG_RETURN_FLOAT8(result);
2616 }
2617 
2618 
2619 /*
2620  * drandom - returns a random number
2621  */
2622 Datum
2624 {
2625  float8 result;
2626 
2627  /* Initialize random seed, if not done yet in this process */
2628  if (unlikely(!drandom_seed_set))
2629  {
2630  /*
2631  * If possible, initialize the seed using high-quality random bits.
2632  * Should that fail for some reason, we fall back on a lower-quality
2633  * seed based on current time and PID.
2634  */
2636  {
2638  uint64 iseed;
2639 
2640  /* Mix the PID with the most predictable bits of the timestamp */
2641  iseed = (uint64) now ^ ((uint64) MyProcPid << 32);
2642  drandom_seed[0] = (unsigned short) iseed;
2643  drandom_seed[1] = (unsigned short) (iseed >> 16);
2644  drandom_seed[2] = (unsigned short) (iseed >> 32);
2645  }
2646  drandom_seed_set = true;
2647  }
2648 
2649  /* pg_erand48 produces desired result range [0.0 - 1.0) */
2650  result = pg_erand48(drandom_seed);
2651 
2652  PG_RETURN_FLOAT8(result);
2653 }
2654 
2655 
2656 /*
2657  * setseed - set seed for the random number generator
2658  */
2659 Datum
2661 {
2662  float8 seed = PG_GETARG_FLOAT8(0);
2663  uint64 iseed;
2664 
2665  if (seed < -1 || seed > 1 || isnan(seed))
2666  ereport(ERROR,
2667  (errcode(ERRCODE_INVALID_PARAMETER_VALUE),
2668  errmsg("setseed parameter %g is out of allowed range [-1,1]",
2669  seed)));
2670 
2671  /* Use sign bit + 47 fractional bits to fill drandom_seed[] */
2672  iseed = (int64) (seed * (float8) UINT64CONST(0x7FFFFFFFFFFF));
2673  drandom_seed[0] = (unsigned short) iseed;
2674  drandom_seed[1] = (unsigned short) (iseed >> 16);
2675  drandom_seed[2] = (unsigned short) (iseed >> 32);
2676  drandom_seed_set = true;
2677 
2678  PG_RETURN_VOID();
2679 }
2680 
2681 
2682 
2683 /*
2684  * =========================
2685  * FLOAT AGGREGATE OPERATORS
2686  * =========================
2687  *
2688  * float8_accum - accumulate for AVG(), variance aggregates, etc.
2689  * float4_accum - same, but input data is float4
2690  * float8_avg - produce final result for float AVG()
2691  * float8_var_samp - produce final result for float VAR_SAMP()
2692  * float8_var_pop - produce final result for float VAR_POP()
2693  * float8_stddev_samp - produce final result for float STDDEV_SAMP()
2694  * float8_stddev_pop - produce final result for float STDDEV_POP()
2695  *
2696  * The naive schoolbook implementation of these aggregates works by
2697  * accumulating sum(X) and sum(X^2). However, this approach suffers from
2698  * large rounding errors in the final computation of quantities like the
2699  * population variance (N*sum(X^2) - sum(X)^2) / N^2, since each of the
2700  * intermediate terms is potentially very large, while the difference is often
2701  * quite small.
2702  *
2703  * Instead we use the Youngs-Cramer algorithm [1] which works by accumulating
2704  * Sx=sum(X) and Sxx=sum((X-Sx/N)^2), using a numerically stable algorithm to
2705  * incrementally update those quantities. The final computations of each of
2706  * the aggregate values is then trivial and gives more accurate results (for
2707  * example, the population variance is just Sxx/N). This algorithm is also
2708  * fairly easy to generalize to allow parallel execution without loss of
2709  * precision (see, for example, [2]). For more details, and a comparison of
2710  * this with other algorithms, see [3].
2711  *
2712  * The transition datatype for all these aggregates is a 3-element array
2713  * of float8, holding the values N, Sx, Sxx in that order.
2714  *
2715  * Note that we represent N as a float to avoid having to build a special
2716  * datatype. Given a reasonable floating-point implementation, there should
2717  * be no accuracy loss unless N exceeds 2 ^ 52 or so (by which time the
2718  * user will have doubtless lost interest anyway...)
2719  *
2720  * [1] Some Results Relevant to Choice of Sum and Sum-of-Product Algorithms,
2721  * E. A. Youngs and E. M. Cramer, Technometrics Vol 13, No 3, August 1971.
2722  *
2723  * [2] Updating Formulae and a Pairwise Algorithm for Computing Sample
2724  * Variances, T. F. Chan, G. H. Golub & R. J. LeVeque, COMPSTAT 1982.
2725  *
2726  * [3] Numerically Stable Parallel Computation of (Co-)Variance, Erich
2727  * Schubert and Michael Gertz, Proceedings of the 30th International
2728  * Conference on Scientific and Statistical Database Management, 2018.
2729  */
2730 
2731 static float8 *
2732 check_float8_array(ArrayType *transarray, const char *caller, int n)
2733 {
2734  /*
2735  * We expect the input to be an N-element float array; verify that. We
2736  * don't need to use deconstruct_array() since the array data is just
2737  * going to look like a C array of N float8 values.
2738  */
2739  if (ARR_NDIM(transarray) != 1 ||
2740  ARR_DIMS(transarray)[0] != n ||
2741  ARR_HASNULL(transarray) ||
2742  ARR_ELEMTYPE(transarray) != FLOAT8OID)
2743  elog(ERROR, "%s: expected %d-element float8 array", caller, n);
2744  return (float8 *) ARR_DATA_PTR(transarray);
2745 }
2746 
2747 /*
2748  * float8_combine
2749  *
2750  * An aggregate combine function used to combine two 3 fields
2751  * aggregate transition data into a single transition data.
2752  * This function is used only in two stage aggregation and
2753  * shouldn't be called outside aggregate context.
2754  */
2755 Datum
2757 {
2758  ArrayType *transarray1 = PG_GETARG_ARRAYTYPE_P(0);
2759  ArrayType *transarray2 = PG_GETARG_ARRAYTYPE_P(1);
2760  float8 *transvalues1;
2761  float8 *transvalues2;
2762  float8 N1,
2763  Sx1,
2764  Sxx1,
2765  N2,
2766  Sx2,
2767  Sxx2,
2768  tmp,
2769  N,
2770  Sx,
2771  Sxx;
2772 
2773  transvalues1 = check_float8_array(transarray1, "float8_combine", 3);
2774  transvalues2 = check_float8_array(transarray2, "float8_combine", 3);
2775 
2776  N1 = transvalues1[0];
2777  Sx1 = transvalues1[1];
2778  Sxx1 = transvalues1[2];
2779 
2780  N2 = transvalues2[0];
2781  Sx2 = transvalues2[1];
2782  Sxx2 = transvalues2[2];
2783 
2784  /*--------------------
2785  * The transition values combine using a generalization of the
2786  * Youngs-Cramer algorithm as follows:
2787  *
2788  * N = N1 + N2
2789  * Sx = Sx1 + Sx2
2790  * Sxx = Sxx1 + Sxx2 + N1 * N2 * (Sx1/N1 - Sx2/N2)^2 / N;
2791  *
2792  * It's worth handling the special cases N1 = 0 and N2 = 0 separately
2793  * since those cases are trivial, and we then don't need to worry about
2794  * division-by-zero errors in the general case.
2795  *--------------------
2796  */
2797  if (N1 == 0.0)
2798  {
2799  N = N2;
2800  Sx = Sx2;
2801  Sxx = Sxx2;
2802  }
2803  else if (N2 == 0.0)
2804  {
2805  N = N1;
2806  Sx = Sx1;
2807  Sxx = Sxx1;
2808  }
2809  else
2810  {
2811  N = N1 + N2;
2812  Sx = float8_pl(Sx1, Sx2);
2813  tmp = Sx1 / N1 - Sx2 / N2;
2814  Sxx = Sxx1 + Sxx2 + N1 * N2 * tmp * tmp / N;
2815  check_float8_val(Sxx, isinf(Sxx1) || isinf(Sxx2), true);
2816  }
2817 
2818  /*
2819  * If we're invoked as an aggregate, we can cheat and modify our first
2820  * parameter in-place to reduce palloc overhead. Otherwise we construct a
2821  * new array with the updated transition data and return it.
2822  */
2823  if (AggCheckCallContext(fcinfo, NULL))
2824  {
2825  transvalues1[0] = N;
2826  transvalues1[1] = Sx;
2827  transvalues1[2] = Sxx;
2828 
2829  PG_RETURN_ARRAYTYPE_P(transarray1);
2830  }
2831  else
2832  {
2833  Datum transdatums[3];
2834  ArrayType *result;
2835 
2836  transdatums[0] = Float8GetDatumFast(N);
2837  transdatums[1] = Float8GetDatumFast(Sx);
2838  transdatums[2] = Float8GetDatumFast(Sxx);
2839 
2840  result = construct_array(transdatums, 3,
2841  FLOAT8OID,
2842  sizeof(float8), FLOAT8PASSBYVAL, 'd');
2843 
2844  PG_RETURN_ARRAYTYPE_P(result);
2845  }
2846 }
2847 
2848 Datum
2850 {
2851  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
2853  float8 *transvalues;
2854  float8 N,
2855  Sx,
2856  Sxx,
2857  tmp;
2858 
2859  transvalues = check_float8_array(transarray, "float8_accum", 3);
2860  N = transvalues[0];
2861  Sx = transvalues[1];
2862  Sxx = transvalues[2];
2863 
2864  /*
2865  * Use the Youngs-Cramer algorithm to incorporate the new value into the
2866  * transition values.
2867  */
2868  N += 1.0;
2869  Sx += newval;
2870  if (transvalues[0] > 0.0)
2871  {
2872  tmp = newval * N - Sx;
2873  Sxx += tmp * tmp / (N * transvalues[0]);
2874 
2875  /*
2876  * Overflow check. We only report an overflow error when finite
2877  * inputs lead to infinite results. Note also that Sxx should be NaN
2878  * if any of the inputs are infinite, so we intentionally prevent Sxx
2879  * from becoming infinite.
2880  */
2881  if (isinf(Sx) || isinf(Sxx))
2882  {
2883  if (!isinf(transvalues[1]) && !isinf(newval))
2884  ereport(ERROR,
2885  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2886  errmsg("value out of range: overflow")));
2887 
2888  Sxx = get_float8_nan();
2889  }
2890  }
2891 
2892  /*
2893  * If we're invoked as an aggregate, we can cheat and modify our first
2894  * parameter in-place to reduce palloc overhead. Otherwise we construct a
2895  * new array with the updated transition data and return it.
2896  */
2897  if (AggCheckCallContext(fcinfo, NULL))
2898  {
2899  transvalues[0] = N;
2900  transvalues[1] = Sx;
2901  transvalues[2] = Sxx;
2902 
2903  PG_RETURN_ARRAYTYPE_P(transarray);
2904  }
2905  else
2906  {
2907  Datum transdatums[3];
2908  ArrayType *result;
2909 
2910  transdatums[0] = Float8GetDatumFast(N);
2911  transdatums[1] = Float8GetDatumFast(Sx);
2912  transdatums[2] = Float8GetDatumFast(Sxx);
2913 
2914  result = construct_array(transdatums, 3,
2915  FLOAT8OID,
2916  sizeof(float8), FLOAT8PASSBYVAL, 'd');
2917 
2918  PG_RETURN_ARRAYTYPE_P(result);
2919  }
2920 }
2921 
2922 Datum
2924 {
2925  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
2926 
2927  /* do computations as float8 */
2929  float8 *transvalues;
2930  float8 N,
2931  Sx,
2932  Sxx,
2933  tmp;
2934 
2935  transvalues = check_float8_array(transarray, "float4_accum", 3);
2936  N = transvalues[0];
2937  Sx = transvalues[1];
2938  Sxx = transvalues[2];
2939 
2940  /*
2941  * Use the Youngs-Cramer algorithm to incorporate the new value into the
2942  * transition values.
2943  */
2944  N += 1.0;
2945  Sx += newval;
2946  if (transvalues[0] > 0.0)
2947  {
2948  tmp = newval * N - Sx;
2949  Sxx += tmp * tmp / (N * transvalues[0]);
2950 
2951  /*
2952  * Overflow check. We only report an overflow error when finite
2953  * inputs lead to infinite results. Note also that Sxx should be NaN
2954  * if any of the inputs are infinite, so we intentionally prevent Sxx
2955  * from becoming infinite.
2956  */
2957  if (isinf(Sx) || isinf(Sxx))
2958  {
2959  if (!isinf(transvalues[1]) && !isinf(newval))
2960  ereport(ERROR,
2961  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
2962  errmsg("value out of range: overflow")));
2963 
2964  Sxx = get_float8_nan();
2965  }
2966  }
2967 
2968  /*
2969  * If we're invoked as an aggregate, we can cheat and modify our first
2970  * parameter in-place to reduce palloc overhead. Otherwise we construct a
2971  * new array with the updated transition data and return it.
2972  */
2973  if (AggCheckCallContext(fcinfo, NULL))
2974  {
2975  transvalues[0] = N;
2976  transvalues[1] = Sx;
2977  transvalues[2] = Sxx;
2978 
2979  PG_RETURN_ARRAYTYPE_P(transarray);
2980  }
2981  else
2982  {
2983  Datum transdatums[3];
2984  ArrayType *result;
2985 
2986  transdatums[0] = Float8GetDatumFast(N);
2987  transdatums[1] = Float8GetDatumFast(Sx);
2988  transdatums[2] = Float8GetDatumFast(Sxx);
2989 
2990  result = construct_array(transdatums, 3,
2991  FLOAT8OID,
2992  sizeof(float8), FLOAT8PASSBYVAL, 'd');
2993 
2994  PG_RETURN_ARRAYTYPE_P(result);
2995  }
2996 }
2997 
2998 Datum
3000 {
3001  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3002  float8 *transvalues;
3003  float8 N,
3004  Sx;
3005 
3006  transvalues = check_float8_array(transarray, "float8_avg", 3);
3007  N = transvalues[0];
3008  Sx = transvalues[1];
3009  /* ignore Sxx */
3010 
3011  /* SQL defines AVG of no values to be NULL */
3012  if (N == 0.0)
3013  PG_RETURN_NULL();
3014 
3015  PG_RETURN_FLOAT8(Sx / N);
3016 }
3017 
3018 Datum
3020 {
3021  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3022  float8 *transvalues;
3023  float8 N,
3024  Sxx;
3025 
3026  transvalues = check_float8_array(transarray, "float8_var_pop", 3);
3027  N = transvalues[0];
3028  /* ignore Sx */
3029  Sxx = transvalues[2];
3030 
3031  /* Population variance is undefined when N is 0, so return NULL */
3032  if (N == 0.0)
3033  PG_RETURN_NULL();
3034 
3035  /* Note that Sxx is guaranteed to be non-negative */
3036 
3037  PG_RETURN_FLOAT8(Sxx / N);
3038 }
3039 
3040 Datum
3042 {
3043  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3044  float8 *transvalues;
3045  float8 N,
3046  Sxx;
3047 
3048  transvalues = check_float8_array(transarray, "float8_var_samp", 3);
3049  N = transvalues[0];
3050  /* ignore Sx */
3051  Sxx = transvalues[2];
3052 
3053  /* Sample variance is undefined when N is 0 or 1, so return NULL */
3054  if (N <= 1.0)
3055  PG_RETURN_NULL();
3056 
3057  /* Note that Sxx is guaranteed to be non-negative */
3058 
3059  PG_RETURN_FLOAT8(Sxx / (N - 1.0));
3060 }
3061 
3062 Datum
3064 {
3065  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3066  float8 *transvalues;
3067  float8 N,
3068  Sxx;
3069 
3070  transvalues = check_float8_array(transarray, "float8_stddev_pop", 3);
3071  N = transvalues[0];
3072  /* ignore Sx */
3073  Sxx = transvalues[2];
3074 
3075  /* Population stddev is undefined when N is 0, so return NULL */
3076  if (N == 0.0)
3077  PG_RETURN_NULL();
3078 
3079  /* Note that Sxx is guaranteed to be non-negative */
3080 
3081  PG_RETURN_FLOAT8(sqrt(Sxx / N));
3082 }
3083 
3084 Datum
3086 {
3087  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3088  float8 *transvalues;
3089  float8 N,
3090  Sxx;
3091 
3092  transvalues = check_float8_array(transarray, "float8_stddev_samp", 3);
3093  N = transvalues[0];
3094  /* ignore Sx */
3095  Sxx = transvalues[2];
3096 
3097  /* Sample stddev is undefined when N is 0 or 1, so return NULL */
3098  if (N <= 1.0)
3099  PG_RETURN_NULL();
3100 
3101  /* Note that Sxx is guaranteed to be non-negative */
3102 
3103  PG_RETURN_FLOAT8(sqrt(Sxx / (N - 1.0)));
3104 }
3105 
3106 /*
3107  * =========================
3108  * SQL2003 BINARY AGGREGATES
3109  * =========================
3110  *
3111  * As with the preceding aggregates, we use the Youngs-Cramer algorithm to
3112  * reduce rounding errors in the aggregate final functions.
3113  *
3114  * The transition datatype for all these aggregates is a 6-element array of
3115  * float8, holding the values N, Sx=sum(X), Sxx=sum((X-Sx/N)^2), Sy=sum(Y),
3116  * Syy=sum((Y-Sy/N)^2), Sxy=sum((X-Sx/N)*(Y-Sy/N)) in that order.
3117  *
3118  * Note that Y is the first argument to all these aggregates!
3119  *
3120  * It might seem attractive to optimize this by having multiple accumulator
3121  * functions that only calculate the sums actually needed. But on most
3122  * modern machines, a couple of extra floating-point multiplies will be
3123  * insignificant compared to the other per-tuple overhead, so I've chosen
3124  * to minimize code space instead.
3125  */
3126 
3127 Datum
3129 {
3130  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3131  float8 newvalY = PG_GETARG_FLOAT8(1);
3132  float8 newvalX = PG_GETARG_FLOAT8(2);
3133  float8 *transvalues;
3134  float8 N,
3135  Sx,
3136  Sxx,
3137  Sy,
3138  Syy,
3139  Sxy,
3140  tmpX,
3141  tmpY,
3142  scale;
3143 
3144  transvalues = check_float8_array(transarray, "float8_regr_accum", 6);
3145  N = transvalues[0];
3146  Sx = transvalues[1];
3147  Sxx = transvalues[2];
3148  Sy = transvalues[3];
3149  Syy = transvalues[4];
3150  Sxy = transvalues[5];
3151 
3152  /*
3153  * Use the Youngs-Cramer algorithm to incorporate the new values into the
3154  * transition values.
3155  */
3156  N += 1.0;
3157  Sx += newvalX;
3158  Sy += newvalY;
3159  if (transvalues[0] > 0.0)
3160  {
3161  tmpX = newvalX * N - Sx;
3162  tmpY = newvalY * N - Sy;
3163  scale = 1.0 / (N * transvalues[0]);
3164  Sxx += tmpX * tmpX * scale;
3165  Syy += tmpY * tmpY * scale;
3166  Sxy += tmpX * tmpY * scale;
3167 
3168  /*
3169  * Overflow check. We only report an overflow error when finite
3170  * inputs lead to infinite results. Note also that Sxx, Syy and Sxy
3171  * should be NaN if any of the relevant inputs are infinite, so we
3172  * intentionally prevent them from becoming infinite.
3173  */
3174  if (isinf(Sx) || isinf(Sxx) || isinf(Sy) || isinf(Syy) || isinf(Sxy))
3175  {
3176  if (((isinf(Sx) || isinf(Sxx)) &&
3177  !isinf(transvalues[1]) && !isinf(newvalX)) ||
3178  ((isinf(Sy) || isinf(Syy)) &&
3179  !isinf(transvalues[3]) && !isinf(newvalY)) ||
3180  (isinf(Sxy) &&
3181  !isinf(transvalues[1]) && !isinf(newvalX) &&
3182  !isinf(transvalues[3]) && !isinf(newvalY)))
3183  ereport(ERROR,
3184  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
3185  errmsg("value out of range: overflow")));
3186 
3187  if (isinf(Sxx))
3188  Sxx = get_float8_nan();
3189  if (isinf(Syy))
3190  Syy = get_float8_nan();
3191  if (isinf(Sxy))
3192  Sxy = get_float8_nan();
3193  }
3194  }
3195 
3196  /*
3197  * If we're invoked as an aggregate, we can cheat and modify our first
3198  * parameter in-place to reduce palloc overhead. Otherwise we construct a
3199  * new array with the updated transition data and return it.
3200  */
3201  if (AggCheckCallContext(fcinfo, NULL))
3202  {
3203  transvalues[0] = N;
3204  transvalues[1] = Sx;
3205  transvalues[2] = Sxx;
3206  transvalues[3] = Sy;
3207  transvalues[4] = Syy;
3208  transvalues[5] = Sxy;
3209 
3210  PG_RETURN_ARRAYTYPE_P(transarray);
3211  }
3212  else
3213  {
3214  Datum transdatums[6];
3215  ArrayType *result;
3216 
3217  transdatums[0] = Float8GetDatumFast(N);
3218  transdatums[1] = Float8GetDatumFast(Sx);
3219  transdatums[2] = Float8GetDatumFast(Sxx);
3220  transdatums[3] = Float8GetDatumFast(Sy);
3221  transdatums[4] = Float8GetDatumFast(Syy);
3222  transdatums[5] = Float8GetDatumFast(Sxy);
3223 
3224  result = construct_array(transdatums, 6,
3225  FLOAT8OID,
3226  sizeof(float8), FLOAT8PASSBYVAL, 'd');
3227 
3228  PG_RETURN_ARRAYTYPE_P(result);
3229  }
3230 }
3231 
3232 /*
3233  * float8_regr_combine
3234  *
3235  * An aggregate combine function used to combine two 6 fields
3236  * aggregate transition data into a single transition data.
3237  * This function is used only in two stage aggregation and
3238  * shouldn't be called outside aggregate context.
3239  */
3240 Datum
3242 {
3243  ArrayType *transarray1 = PG_GETARG_ARRAYTYPE_P(0);
3244  ArrayType *transarray2 = PG_GETARG_ARRAYTYPE_P(1);
3245  float8 *transvalues1;
3246  float8 *transvalues2;
3247  float8 N1,
3248  Sx1,
3249  Sxx1,
3250  Sy1,
3251  Syy1,
3252  Sxy1,
3253  N2,
3254  Sx2,
3255  Sxx2,
3256  Sy2,
3257  Syy2,
3258  Sxy2,
3259  tmp1,
3260  tmp2,
3261  N,
3262  Sx,
3263  Sxx,
3264  Sy,
3265  Syy,
3266  Sxy;
3267 
3268  transvalues1 = check_float8_array(transarray1, "float8_regr_combine", 6);
3269  transvalues2 = check_float8_array(transarray2, "float8_regr_combine", 6);
3270 
3271  N1 = transvalues1[0];
3272  Sx1 = transvalues1[1];
3273  Sxx1 = transvalues1[2];
3274  Sy1 = transvalues1[3];
3275  Syy1 = transvalues1[4];
3276  Sxy1 = transvalues1[5];
3277 
3278  N2 = transvalues2[0];
3279  Sx2 = transvalues2[1];
3280  Sxx2 = transvalues2[2];
3281  Sy2 = transvalues2[3];
3282  Syy2 = transvalues2[4];
3283  Sxy2 = transvalues2[5];
3284 
3285  /*--------------------
3286  * The transition values combine using a generalization of the
3287  * Youngs-Cramer algorithm as follows:
3288  *
3289  * N = N1 + N2
3290  * Sx = Sx1 + Sx2
3291  * Sxx = Sxx1 + Sxx2 + N1 * N2 * (Sx1/N1 - Sx2/N2)^2 / N
3292  * Sy = Sy1 + Sy2
3293  * Syy = Syy1 + Syy2 + N1 * N2 * (Sy1/N1 - Sy2/N2)^2 / N
3294  * Sxy = Sxy1 + Sxy2 + N1 * N2 * (Sx1/N1 - Sx2/N2) * (Sy1/N1 - Sy2/N2) / N
3295  *
3296  * It's worth handling the special cases N1 = 0 and N2 = 0 separately
3297  * since those cases are trivial, and we then don't need to worry about
3298  * division-by-zero errors in the general case.
3299  *--------------------
3300  */
3301  if (N1 == 0.0)
3302  {
3303  N = N2;
3304  Sx = Sx2;
3305  Sxx = Sxx2;
3306  Sy = Sy2;
3307  Syy = Syy2;
3308  Sxy = Sxy2;
3309  }
3310  else if (N2 == 0.0)
3311  {
3312  N = N1;
3313  Sx = Sx1;
3314  Sxx = Sxx1;
3315  Sy = Sy1;
3316  Syy = Syy1;
3317  Sxy = Sxy1;
3318  }
3319  else
3320  {
3321  N = N1 + N2;
3322  Sx = float8_pl(Sx1, Sx2);
3323  tmp1 = Sx1 / N1 - Sx2 / N2;
3324  Sxx = Sxx1 + Sxx2 + N1 * N2 * tmp1 * tmp1 / N;
3325  check_float8_val(Sxx, isinf(Sxx1) || isinf(Sxx2), true);
3326  Sy = float8_pl(Sy1, Sy2);
3327  tmp2 = Sy1 / N1 - Sy2 / N2;
3328  Syy = Syy1 + Syy2 + N1 * N2 * tmp2 * tmp2 / N;
3329  check_float8_val(Syy, isinf(Syy1) || isinf(Syy2), true);
3330  Sxy = Sxy1 + Sxy2 + N1 * N2 * tmp1 * tmp2 / N;
3331  check_float8_val(Sxy, isinf(Sxy1) || isinf(Sxy2), true);
3332  }
3333 
3334  /*
3335  * If we're invoked as an aggregate, we can cheat and modify our first
3336  * parameter in-place to reduce palloc overhead. Otherwise we construct a
3337  * new array with the updated transition data and return it.
3338  */
3339  if (AggCheckCallContext(fcinfo, NULL))
3340  {
3341  transvalues1[0] = N;
3342  transvalues1[1] = Sx;
3343  transvalues1[2] = Sxx;
3344  transvalues1[3] = Sy;
3345  transvalues1[4] = Syy;
3346  transvalues1[5] = Sxy;
3347 
3348  PG_RETURN_ARRAYTYPE_P(transarray1);
3349  }
3350  else
3351  {
3352  Datum transdatums[6];
3353  ArrayType *result;
3354 
3355  transdatums[0] = Float8GetDatumFast(N);
3356  transdatums[1] = Float8GetDatumFast(Sx);
3357  transdatums[2] = Float8GetDatumFast(Sxx);
3358  transdatums[3] = Float8GetDatumFast(Sy);
3359  transdatums[4] = Float8GetDatumFast(Syy);
3360  transdatums[5] = Float8GetDatumFast(Sxy);
3361 
3362  result = construct_array(transdatums, 6,
3363  FLOAT8OID,
3364  sizeof(float8), FLOAT8PASSBYVAL, 'd');
3365 
3366  PG_RETURN_ARRAYTYPE_P(result);
3367  }
3368 }
3369 
3370 
3371 Datum
3373 {
3374  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3375  float8 *transvalues;
3376  float8 N,
3377  Sxx;
3378 
3379  transvalues = check_float8_array(transarray, "float8_regr_sxx", 6);
3380  N = transvalues[0];
3381  Sxx = transvalues[2];
3382 
3383  /* if N is 0 we should return NULL */
3384  if (N < 1.0)
3385  PG_RETURN_NULL();
3386 
3387  /* Note that Sxx is guaranteed to be non-negative */
3388 
3389  PG_RETURN_FLOAT8(Sxx);
3390 }
3391 
3392 Datum
3394 {
3395  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3396  float8 *transvalues;
3397  float8 N,
3398  Syy;
3399 
3400  transvalues = check_float8_array(transarray, "float8_regr_syy", 6);
3401  N = transvalues[0];
3402  Syy = transvalues[4];
3403 
3404  /* if N is 0 we should return NULL */
3405  if (N < 1.0)
3406  PG_RETURN_NULL();
3407 
3408  /* Note that Syy is guaranteed to be non-negative */
3409 
3410  PG_RETURN_FLOAT8(Syy);
3411 }
3412 
3413 Datum
3415 {
3416  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3417  float8 *transvalues;
3418  float8 N,
3419  Sxy;
3420 
3421  transvalues = check_float8_array(transarray, "float8_regr_sxy", 6);
3422  N = transvalues[0];
3423  Sxy = transvalues[5];
3424 
3425  /* if N is 0 we should return NULL */
3426  if (N < 1.0)
3427  PG_RETURN_NULL();
3428 
3429  /* A negative result is valid here */
3430 
3431  PG_RETURN_FLOAT8(Sxy);
3432 }
3433 
3434 Datum
3436 {
3437  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3438  float8 *transvalues;
3439  float8 N,
3440  Sx;
3441 
3442  transvalues = check_float8_array(transarray, "float8_regr_avgx", 6);
3443  N = transvalues[0];
3444  Sx = transvalues[1];
3445 
3446  /* if N is 0 we should return NULL */
3447  if (N < 1.0)
3448  PG_RETURN_NULL();
3449 
3450  PG_RETURN_FLOAT8(Sx / N);
3451 }
3452 
3453 Datum
3455 {
3456  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3457  float8 *transvalues;
3458  float8 N,
3459  Sy;
3460 
3461  transvalues = check_float8_array(transarray, "float8_regr_avgy", 6);
3462  N = transvalues[0];
3463  Sy = transvalues[3];
3464 
3465  /* if N is 0 we should return NULL */
3466  if (N < 1.0)
3467  PG_RETURN_NULL();
3468 
3469  PG_RETURN_FLOAT8(Sy / N);
3470 }
3471 
3472 Datum
3474 {
3475  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3476  float8 *transvalues;
3477  float8 N,
3478  Sxy;
3479 
3480  transvalues = check_float8_array(transarray, "float8_covar_pop", 6);
3481  N = transvalues[0];
3482  Sxy = transvalues[5];
3483 
3484  /* if N is 0 we should return NULL */
3485  if (N < 1.0)
3486  PG_RETURN_NULL();
3487 
3488  PG_RETURN_FLOAT8(Sxy / N);
3489 }
3490 
3491 Datum
3493 {
3494  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3495  float8 *transvalues;
3496  float8 N,
3497  Sxy;
3498 
3499  transvalues = check_float8_array(transarray, "float8_covar_samp", 6);
3500  N = transvalues[0];
3501  Sxy = transvalues[5];
3502 
3503  /* if N is <= 1 we should return NULL */
3504  if (N < 2.0)
3505  PG_RETURN_NULL();
3506 
3507  PG_RETURN_FLOAT8(Sxy / (N - 1.0));
3508 }
3509 
3510 Datum
3512 {
3513  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3514  float8 *transvalues;
3515  float8 N,
3516  Sxx,
3517  Syy,
3518  Sxy;
3519 
3520  transvalues = check_float8_array(transarray, "float8_corr", 6);
3521  N = transvalues[0];
3522  Sxx = transvalues[2];
3523  Syy = transvalues[4];
3524  Sxy = transvalues[5];
3525 
3526  /* if N is 0 we should return NULL */
3527  if (N < 1.0)
3528  PG_RETURN_NULL();
3529 
3530  /* Note that Sxx and Syy are guaranteed to be non-negative */
3531 
3532  /* per spec, return NULL for horizontal and vertical lines */
3533  if (Sxx == 0 || Syy == 0)
3534  PG_RETURN_NULL();
3535 
3536  PG_RETURN_FLOAT8(Sxy / sqrt(Sxx * Syy));
3537 }
3538 
3539 Datum
3541 {
3542  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3543  float8 *transvalues;
3544  float8 N,
3545  Sxx,
3546  Syy,
3547  Sxy;
3548 
3549  transvalues = check_float8_array(transarray, "float8_regr_r2", 6);
3550  N = transvalues[0];
3551  Sxx = transvalues[2];
3552  Syy = transvalues[4];
3553  Sxy = transvalues[5];
3554 
3555  /* if N is 0 we should return NULL */
3556  if (N < 1.0)
3557  PG_RETURN_NULL();
3558 
3559  /* Note that Sxx and Syy are guaranteed to be non-negative */
3560 
3561  /* per spec, return NULL for a vertical line */
3562  if (Sxx == 0)
3563  PG_RETURN_NULL();
3564 
3565  /* per spec, return 1.0 for a horizontal line */
3566  if (Syy == 0)
3567  PG_RETURN_FLOAT8(1.0);
3568 
3569  PG_RETURN_FLOAT8((Sxy * Sxy) / (Sxx * Syy));
3570 }
3571 
3572 Datum
3574 {
3575  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3576  float8 *transvalues;
3577  float8 N,
3578  Sxx,
3579  Sxy;
3580 
3581  transvalues = check_float8_array(transarray, "float8_regr_slope", 6);
3582  N = transvalues[0];
3583  Sxx = transvalues[2];
3584  Sxy = transvalues[5];
3585 
3586  /* if N is 0 we should return NULL */
3587  if (N < 1.0)
3588  PG_RETURN_NULL();
3589 
3590  /* Note that Sxx is guaranteed to be non-negative */
3591 
3592  /* per spec, return NULL for a vertical line */
3593  if (Sxx == 0)
3594  PG_RETURN_NULL();
3595 
3596  PG_RETURN_FLOAT8(Sxy / Sxx);
3597 }
3598 
3599 Datum
3601 {
3602  ArrayType *transarray = PG_GETARG_ARRAYTYPE_P(0);
3603  float8 *transvalues;
3604  float8 N,
3605  Sx,
3606  Sxx,
3607  Sy,
3608  Sxy;
3609 
3610  transvalues = check_float8_array(transarray, "float8_regr_intercept", 6);
3611  N = transvalues[0];
3612  Sx = transvalues[1];
3613  Sxx = transvalues[2];
3614  Sy = transvalues[3];
3615  Sxy = transvalues[5];
3616 
3617  /* if N is 0 we should return NULL */
3618  if (N < 1.0)
3619  PG_RETURN_NULL();
3620 
3621  /* Note that Sxx is guaranteed to be non-negative */
3622 
3623  /* per spec, return NULL for a vertical line */
3624  if (Sxx == 0)
3625  PG_RETURN_NULL();
3626 
3627  PG_RETURN_FLOAT8((Sy - Sx * Sxy / Sxx) / N);
3628 }
3629 
3630 
3631 /*
3632  * ====================================
3633  * MIXED-PRECISION ARITHMETIC OPERATORS
3634  * ====================================
3635  */
3636 
3637 /*
3638  * float48pl - returns arg1 + arg2
3639  * float48mi - returns arg1 - arg2
3640  * float48mul - returns arg1 * arg2
3641  * float48div - returns arg1 / arg2
3642  */
3643 Datum
3645 {
3646  float4 arg1 = PG_GETARG_FLOAT4(0);
3647  float8 arg2 = PG_GETARG_FLOAT8(1);
3648 
3649  PG_RETURN_FLOAT8(float8_pl((float8) arg1, arg2));
3650 }
3651 
3652 Datum
3654 {
3655  float4 arg1 = PG_GETARG_FLOAT4(0);
3656  float8 arg2 = PG_GETARG_FLOAT8(1);
3657 
3658  PG_RETURN_FLOAT8(float8_mi((float8) arg1, arg2));
3659 }
3660 
3661 Datum
3663 {
3664  float4 arg1 = PG_GETARG_FLOAT4(0);
3665  float8 arg2 = PG_GETARG_FLOAT8(1);
3666 
3667  PG_RETURN_FLOAT8(float8_mul((float8) arg1, arg2));
3668 }
3669 
3670 Datum
3672 {
3673  float4 arg1 = PG_GETARG_FLOAT4(0);
3674  float8 arg2 = PG_GETARG_FLOAT8(1);
3675 
3676  PG_RETURN_FLOAT8(float8_div((float8) arg1, arg2));
3677 }
3678 
3679 /*
3680  * float84pl - returns arg1 + arg2
3681  * float84mi - returns arg1 - arg2
3682  * float84mul - returns arg1 * arg2
3683  * float84div - returns arg1 / arg2
3684  */
3685 Datum
3687 {
3688  float8 arg1 = PG_GETARG_FLOAT8(0);
3689  float4 arg2 = PG_GETARG_FLOAT4(1);
3690 
3691  PG_RETURN_FLOAT8(float8_pl(arg1, (float8) arg2));
3692 }
3693 
3694 Datum
3696 {
3697  float8 arg1 = PG_GETARG_FLOAT8(0);
3698  float4 arg2 = PG_GETARG_FLOAT4(1);
3699 
3700  PG_RETURN_FLOAT8(float8_mi(arg1, (float8) arg2));
3701 }
3702 
3703 Datum
3705 {
3706  float8 arg1 = PG_GETARG_FLOAT8(0);
3707  float4 arg2 = PG_GETARG_FLOAT4(1);
3708 
3709  PG_RETURN_FLOAT8(float8_mul(arg1, (float8) arg2));
3710 }
3711 
3712 Datum
3714 {
3715  float8 arg1 = PG_GETARG_FLOAT8(0);
3716  float4 arg2 = PG_GETARG_FLOAT4(1);
3717 
3718  PG_RETURN_FLOAT8(float8_div(arg1, (float8) arg2));
3719 }
3720 
3721 /*
3722  * ====================
3723  * COMPARISON OPERATORS
3724  * ====================
3725  */
3726 
3727 /*
3728  * float48{eq,ne,lt,le,gt,ge} - float4/float8 comparison operations
3729  */
3730 Datum
3732 {
3733  float4 arg1 = PG_GETARG_FLOAT4(0);
3734  float8 arg2 = PG_GETARG_FLOAT8(1);
3735 
3736  PG_RETURN_BOOL(float8_eq((float8) arg1, arg2));
3737 }
3738 
3739 Datum
3741 {
3742  float4 arg1 = PG_GETARG_FLOAT4(0);
3743  float8 arg2 = PG_GETARG_FLOAT8(1);
3744 
3745  PG_RETURN_BOOL(float8_ne((float8) arg1, arg2));
3746 }
3747 
3748 Datum
3750 {
3751  float4 arg1 = PG_GETARG_FLOAT4(0);
3752  float8 arg2 = PG_GETARG_FLOAT8(1);
3753 
3754  PG_RETURN_BOOL(float8_lt((float8) arg1, arg2));
3755 }
3756 
3757 Datum
3759 {
3760  float4 arg1 = PG_GETARG_FLOAT4(0);
3761  float8 arg2 = PG_GETARG_FLOAT8(1);
3762 
3763  PG_RETURN_BOOL(float8_le((float8) arg1, arg2));
3764 }
3765 
3766 Datum
3768 {
3769  float4 arg1 = PG_GETARG_FLOAT4(0);
3770  float8 arg2 = PG_GETARG_FLOAT8(1);
3771 
3772  PG_RETURN_BOOL(float8_gt((float8) arg1, arg2));
3773 }
3774 
3775 Datum
3777 {
3778  float4 arg1 = PG_GETARG_FLOAT4(0);
3779  float8 arg2 = PG_GETARG_FLOAT8(1);
3780 
3781  PG_RETURN_BOOL(float8_ge((float8) arg1, arg2));
3782 }
3783 
3784 /*
3785  * float84{eq,ne,lt,le,gt,ge} - float8/float4 comparison operations
3786  */
3787 Datum
3789 {
3790  float8 arg1 = PG_GETARG_FLOAT8(0);
3791  float4 arg2 = PG_GETARG_FLOAT4(1);
3792 
3793  PG_RETURN_BOOL(float8_eq(arg1, (float8) arg2));
3794 }
3795 
3796 Datum
3798 {
3799  float8 arg1 = PG_GETARG_FLOAT8(0);
3800  float4 arg2 = PG_GETARG_FLOAT4(1);
3801 
3802  PG_RETURN_BOOL(float8_ne(arg1, (float8) arg2));
3803 }
3804 
3805 Datum
3807 {
3808  float8 arg1 = PG_GETARG_FLOAT8(0);
3809  float4 arg2 = PG_GETARG_FLOAT4(1);
3810 
3811  PG_RETURN_BOOL(float8_lt(arg1, (float8) arg2));
3812 }
3813 
3814 Datum
3816 {
3817  float8 arg1 = PG_GETARG_FLOAT8(0);
3818  float4 arg2 = PG_GETARG_FLOAT4(1);
3819 
3820  PG_RETURN_BOOL(float8_le(arg1, (float8) arg2));
3821 }
3822 
3823 Datum
3825 {
3826  float8 arg1 = PG_GETARG_FLOAT8(0);
3827  float4 arg2 = PG_GETARG_FLOAT4(1);
3828 
3829  PG_RETURN_BOOL(float8_gt(arg1, (float8) arg2));
3830 }
3831 
3832 Datum
3834 {
3835  float8 arg1 = PG_GETARG_FLOAT8(0);
3836  float4 arg2 = PG_GETARG_FLOAT4(1);
3837 
3838  PG_RETURN_BOOL(float8_ge(arg1, (float8) arg2));
3839 }
3840 
3841 /*
3842  * Implements the float8 version of the width_bucket() function
3843  * defined by SQL2003. See also width_bucket_numeric().
3844  *
3845  * 'bound1' and 'bound2' are the lower and upper bounds of the
3846  * histogram's range, respectively. 'count' is the number of buckets
3847  * in the histogram. width_bucket() returns an integer indicating the
3848  * bucket number that 'operand' belongs to in an equiwidth histogram
3849  * with the specified characteristics. An operand smaller than the
3850  * lower bound is assigned to bucket 0. An operand greater than the
3851  * upper bound is assigned to an additional bucket (with number
3852  * count+1). We don't allow "NaN" for any of the float8 inputs, and we
3853  * don't allow either of the histogram bounds to be +/- infinity.
3854  */
3855 Datum
3857 {
3858  float8 operand = PG_GETARG_FLOAT8(0);
3859  float8 bound1 = PG_GETARG_FLOAT8(1);
3860  float8 bound2 = PG_GETARG_FLOAT8(2);
3861  int32 count = PG_GETARG_INT32(3);
3862  int32 result;
3863 
3864  if (count <= 0.0)
3865  ereport(ERROR,
3866  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_WIDTH_BUCKET_FUNCTION),
3867  errmsg("count must be greater than zero")));
3868 
3869  if (isnan(operand) || isnan(bound1) || isnan(bound2))
3870  ereport(ERROR,
3871  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_WIDTH_BUCKET_FUNCTION),
3872  errmsg("operand, lower bound, and upper bound cannot be NaN")));
3873 
3874  /* Note that we allow "operand" to be infinite */
3875  if (isinf(bound1) || isinf(bound2))
3876  ereport(ERROR,
3877  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_WIDTH_BUCKET_FUNCTION),
3878  errmsg("lower and upper bounds must be finite")));
3879 
3880  if (bound1 < bound2)
3881  {
3882  if (operand < bound1)
3883  result = 0;
3884  else if (operand >= bound2)
3885  {
3886  if (pg_add_s32_overflow(count, 1, &result))
3887  ereport(ERROR,
3888  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
3889  errmsg("integer out of range")));
3890  }
3891  else
3892  result = ((float8) count * (operand - bound1) / (bound2 - bound1)) + 1;
3893  }
3894  else if (bound1 > bound2)
3895  {
3896  if (operand > bound1)
3897  result = 0;
3898  else if (operand <= bound2)
3899  {
3900  if (pg_add_s32_overflow(count, 1, &result))
3901  ereport(ERROR,
3902  (errcode(ERRCODE_NUMERIC_VALUE_OUT_OF_RANGE),
3903  errmsg("integer out of range")));
3904  }
3905  else
3906  result = ((float8) count * (bound1 - operand) / (bound1 - bound2)) + 1;
3907  }
3908  else
3909  {
3910  ereport(ERROR,
3911  (errcode(ERRCODE_INVALID_ARGUMENT_FOR_WIDTH_BUCKET_FUNCTION),
3912  errmsg("lower bound cannot equal upper bound")));
3913  result = 0; /* keep the compiler quiet */
3914  }
3915 
3916  PG_RETURN_INT32(result);
3917 }
3918 
3919 /* ========== PRIVATE ROUTINES ========== */
3920 
3921 #ifndef HAVE_CBRT
3922 
3923 static double
3924 cbrt(double x)
3925 {
3926  int isneg = (x < 0.0);
3927  double absx = fabs(x);
3928  double tmpres = pow(absx, (double) 1.0 / (double) 3.0);
3929 
3930  /*
3931  * The result is somewhat inaccurate --- not really pow()'s fault, as the
3932  * exponent it's handed contains roundoff error. We can improve the
3933  * accuracy by doing one iteration of Newton's formula. Beware of zero
3934  * input however.
3935  */
3936  if (tmpres > 0.0)
3937  tmpres -= (tmpres - absx / (tmpres * tmpres)) / (double) 3.0;
3938 
3939  return isneg ? -tmpres : tmpres;
3940 }
3941 
3942 #endif /* !HAVE_CBRT */
Datum float48lt(PG_FUNCTION_ARGS)
Definition: float.c:3749
struct SortSupportData * SortSupport
Definition: sortsupport.h:58
Datum float8ne(PG_FUNCTION_ARGS)
Definition: float.c:938
Datum dtoi2(PG_FUNCTION_ARGS)
Definition: float.c:1236
#define PG_GETARG_FLOAT8(n)
Definition: fmgr.h:276
signed short int16
Definition: c.h:345
Datum dlog1(PG_FUNCTION_ARGS)
Definition: float.c:1591
static bool float8_ne(const float8 val1, const float8 val2)
Definition: float.h:302
Datum float4pl(PG_FUNCTION_ARGS)
Definition: float.c:737
Datum float4mi(PG_FUNCTION_ARGS)
Definition: float.c:746
Datum float84ge(PG_FUNCTION_ARGS)
Definition: float.c:3833
Datum dcos(PG_FUNCTION_ARGS)
Definition: float.c:1760
float8 degree_c_sixty
Definition: float.c:62
Datum dacos(PG_FUNCTION_ARGS)
Definition: float.c:1650
#define PG_GETARG_INT32(n)
Definition: fmgr.h:264
#define cbrt
Definition: float.c:84
Datum radians(PG_FUNCTION_ARGS)
Definition: float.c:2457
Datum float8_stddev_pop(PG_FUNCTION_ARGS)
Definition: float.c:3063
static float8 cot_45
Definition: float.c:52
Datum float4eq(PG_FUNCTION_ARGS)
Definition: float.c:835
#define RADIANS_PER_DEGREE
Definition: float.h:26
Datum float8larger(PG_FUNCTION_ARGS)
Definition: float.c:696
float8 degree_c_one
Definition: float.c:64
#define PG_RETURN_FLOAT4(x)
Definition: fmgr.h:355
Datum dcosd(PG_FUNCTION_ARGS)
Definition: float.c:2198
int MyProcPid
Definition: globals.c:40
static float4 float4_div(const float4 val1, const float4 val2)
Definition: float.h:244
Datum float8_regr_syy(PG_FUNCTION_ARGS)
Definition: float.c:3393
static float8 get_float8_infinity(void)
Definition: float.h:90
Datum float4le(PG_FUNCTION_ARGS)
Definition: float.c:862
Datum datand(PG_FUNCTION_ARGS)
Definition: float.c:2066
Datum float48mi(PG_FUNCTION_ARGS)
Definition: float.c:3653
#define INIT_DEGREE_CONSTANTS()
Definition: float.c:1919
float8 degree_c_one_half
Definition: float.c:63
Datum float4larger(PG_FUNCTION_ARGS)
Definition: float.c:629
Datum float8mi(PG_FUNCTION_ARGS)
Definition: float.c:788
Datum float8_regr_sxy(PG_FUNCTION_ARGS)
Definition: float.c:3414
Datum dsin(PG_FUNCTION_ARGS)
Definition: float.c:1827
Datum datanh(PG_FUNCTION_ARGS)
Definition: float.c:2587
Datum float8gt(PG_FUNCTION_ARGS)
Definition: float.c:965
Datum dpi(PG_FUNCTION_ARGS)
Definition: float.c:2447
Datum float48eq(PG_FUNCTION_ARGS)
Definition: float.c:3731
Datum btfloat8cmp(PG_FUNCTION_ARGS)
Definition: float.c:983
char * float8out_internal(double num)
Definition: float.c:546
Datum float8in(PG_FUNCTION_ARGS)
Definition: float.c:332
static double asind_q1(double x)
Definition: float.c:1936
TimestampTz GetCurrentTimestamp(void)
Definition: timestamp.c:1583
Datum float8_var_samp(PG_FUNCTION_ARGS)
Definition: float.c:3041
int64 TimestampTz
Definition: timestamp.h:39
static bool degree_consts_set
Definition: float.c:45
Datum i4tof(PG_FUNCTION_ARGS)
Definition: float.c:1356
Datum float4ge(PG_FUNCTION_ARGS)
Definition: float.c:880
void pq_begintypsend(StringInfo buf)
Definition: pqformat.c:328
Datum float8um(PG_FUNCTION_ARGS)
Definition: float.c:678
char * pstrdup(const char *in)
Definition: mcxt.c:1186
Datum ftod(PG_FUNCTION_ARGS)
Definition: float.c:1178
static unsigned short drandom_seed[3]
Definition: float.c:68
void pq_sendfloat8(StringInfo buf, float8 f)
Definition: pqformat.c:278
StringInfoData * StringInfo
Definition: stringinfo.h:43
#define PG_RETURN_FLOAT8(x)
Definition: fmgr.h:356
Datum dpow(PG_FUNCTION_ARGS)
Definition: float.c:1500
Datum float84div(PG_FUNCTION_ARGS)
Definition: float.c:3713
static double sind_q1(double x)
Definition: float.c:2159
Datum float8up(PG_FUNCTION_ARGS)
Definition: float.c:688
Datum dtrunc(PG_FUNCTION_ARGS)
Definition: float.c:1446
#define PG_RETURN_INT32(x)
Definition: fmgr.h:344
Datum float8_covar_pop(PG_FUNCTION_ARGS)
Definition: float.c:3473
Datum float48le(PG_FUNCTION_ARGS)
Definition: float.c:3758
ArrayType * construct_array(Datum *elems, int nelems, Oid elmtype, int elmlen, bool elmbyval, char elmalign)
Definition: arrayfuncs.c:3291
static bool float8_le(const float8 val1, const float8 val2)
Definition: float.h:326
Datum ascii(PG_FUNCTION_ARGS)
Datum dceil(PG_FUNCTION_ARGS)
Definition: float.c:1398
Datum float8eq(PG_FUNCTION_ARGS)
Definition: float.c:929
#define PG_RETURN_INT16(x)
Definition: fmgr.h:346
void pq_sendfloat4(StringInfo buf, float4 f)
Definition: pqformat.c:254
int errcode(int sqlerrcode)
Definition: elog.c:570
Datum btfloat4cmp(PG_FUNCTION_ARGS)
Definition: float.c:889
int scale
Definition: pgbench.c:151
Datum dround(PG_FUNCTION_ARGS)
Definition: float.c:1386
static float8 acos_0_5
Definition: float.c:49
Datum float8_regr_r2(PG_FUNCTION_ARGS)
Definition: float.c:3540
#define PG_GETARG_POINTER(n)
Definition: fmgr.h:271
Datum float84ne(PG_FUNCTION_ARGS)
Definition: float.c:3797
float8 degree_c_thirty
Definition: float.c:60
static void init_degree_constants(void)
Definition: float.c:1907
Datum float8ge(PG_FUNCTION_ARGS)
Definition: float.c:974
#define PG_GETARG_BOOL(n)
Definition: fmgr.h:269
#define PG_RETURN_BYTEA_P(x)
Definition: fmgr.h:360
static bool drandom_seed_set
Definition: float.c:67
static float8 float8_mul(const float8 val1, const float8 val2)
Definition: float.h:232
Datum float8_regr_sxx(PG_FUNCTION_ARGS)
Definition: float.c:3372
Datum drandom(PG_FUNCTION_ARGS)
Definition: float.c:2623
Datum dlog10(PG_FUNCTION_ARGS)
Definition: float.c:1620
bytea * pq_endtypsend(StringInfo buf)
Definition: pqformat.c:348
Datum float8_avg(PG_FUNCTION_ARGS)
Definition: float.c:2999
int float8_cmp_internal(float8 a, float8 b)
Definition: float.c:919
int float4_cmp_internal(float4 a, float4 b)
Definition: float.c:825
float strtof(const char *nptr, char **endptr)
Definition: strtof.c:28
signed int int32
Definition: c.h:346
Datum float48pl(PG_FUNCTION_ARGS)
Definition: float.c:3644
Datum float4_accum(PG_FUNCTION_ARGS)
Definition: float.c:2923
Datum float8send(PG_FUNCTION_ARGS)
Definition: float.c:576
Datum float4up(PG_FUNCTION_ARGS)
Definition: float.c:621
Datum float48ge(PG_FUNCTION_ARGS)
Definition: float.c:3776
int pg_strncasecmp(const char *s1, const char *s2, size_t n)
Definition: pgstrcasecmp.c:69
Datum dtanh(PG_FUNCTION_ARGS)
Definition: float.c:2525
#define PG_GETARG_ARRAYTYPE_P(n)
Definition: array.h:251
static void check_float8_val(const float8 val, const bool inf_is_valid, const bool zero_is_valid)
Definition: float.h:151
Datum float8_regr_slope(PG_FUNCTION_ARGS)
Definition: float.c:3573
Datum dcosh(PG_FUNCTION_ARGS)
Definition: float.c:2502
static float4 float4_mi(const float4 val1, const float4 val2)
Definition: float.h:198
Datum datan2(PG_FUNCTION_ARGS)
Definition: float.c:1735
Datum ftoi2(PG_FUNCTION_ARGS)
Definition: float.c:1324
Datum dtof(PG_FUNCTION_ARGS)
Definition: float.c:1190
Datum float8pl(PG_FUNCTION_ARGS)
Definition: float.c:779
double float8in_internal(char *num, char **endptr_p, const char *type_name, const char *orig_string)
Definition: float.c:518
#define ERROR
Definition: elog.h:43
double float8
Definition: c.h:491
Datum in_range_float4_float8(PG_FUNCTION_ARGS)
Definition: float.c:1105
Datum float84pl(PG_FUNCTION_ARGS)
Definition: float.c:3686
Datum float4div(PG_FUNCTION_ARGS)
Definition: float.c:764
static void check_float4_val(const float4 val, const bool inf_is_valid, const bool zero_is_valid)
Definition: float.h:136
#define ARR_DIMS(a)
Definition: array.h:282
static double acosd_q1(double x)
Definition: float.c:1969
Datum dcbrt(PG_FUNCTION_ARGS)
Definition: float.c:1485
Datum setseed(PG_FUNCTION_ARGS)
Definition: float.c:2660
static int btfloat8fastcmp(Datum x, Datum y, SortSupport ssup)
Definition: float.c:992
static float8 get_float8_nan(void)
Definition: float.h:119
int isinf(double x)
Datum float8abs(PG_FUNCTION_ARGS)
Definition: float.c:666
char sign
Definition: informix.c:688
int(* comparator)(Datum x, Datum y, SortSupport ssup)
Definition: sortsupport.h:106
static float8 atan_1_0
Definition: float.c:50
Datum dasin(PG_FUNCTION_ARGS)
Definition: float.c:1680
#define ARR_DATA_PTR(a)
Definition: array.h:310
Datum float8div(PG_FUNCTION_ARGS)
Definition: float.c:806
float8 degree_c_forty_five
Definition: float.c:61
static char * buf
Definition: pg_test_fsync.c:68
Datum dtoi4(PG_FUNCTION_ARGS)
Definition: float.c:1204
Datum float8_regr_combine(PG_FUNCTION_ARGS)
Definition: float.c:3241
Datum float84mul(PG_FUNCTION_ARGS)
Definition: float.c:3704
Datum datan(PG_FUNCTION_ARGS)
Definition: float.c:1710
static double cosd_0_to_60(double x)
Definition: float.c:2146
Datum float4abs(PG_FUNCTION_ARGS)
Definition: float.c:600
static bool float8_ge(const float8 val1, const float8 val2)
Definition: float.h:350
#define ARR_HASNULL(a)
Definition: array.h:279
static float8 float8_pl(const float8 val1, const float8 val2)
Definition: float.h:187
int float_to_shortest_decimal_buf(float f, char *result)
Definition: f2s.c:781
static float8 float8_mi(const float8 val1, const float8 val2)
Definition: float.h:209
static float8 sin_30
Definition: float.c:46
double rint(double x)
Definition: rint.c:21
#define PG_INT16_MIN
Definition: c.h:437
Datum ftoi4(PG_FUNCTION_ARGS)
Definition: float.c:1292
#define PG_INT32_MIN
Definition: c.h:440
static bool float4_gt(const float4 val1, const float4 val2)
Definition: float.h:332
#define PG_RETURN_ARRAYTYPE_P(x)
Definition: array.h:253
#define PG_GETARG_FLOAT4(n)
Definition: fmgr.h:275
Datum dsind(PG_FUNCTION_ARGS)
Definition: float.c:2316
#define ereport(elevel, rest)
Definition: elog.h:141
Datum float8mul(PG_FUNCTION_ARGS)
Definition: float.c:797
Datum dfloor(PG_FUNCTION_ARGS)
Definition: float.c:1410
Datum btfloat48cmp(PG_FUNCTION_ARGS)
Definition: float.c:1010
Datum float4in(PG_FUNCTION_ARGS)
Definition: float.c:143
Datum dsign(PG_FUNCTION_ARGS)
Definition: float.c:1423
Datum float4mul(PG_FUNCTION_ARGS)
Definition: float.c:755
static bool float8_lt(const float8 val1, const float8 val2)
Definition: float.h:314
Datum float8_covar_samp(PG_FUNCTION_ARGS)
Definition: float.c:3492
Datum float8lt(PG_FUNCTION_ARGS)
Definition: float.c:947
Datum float8_regr_intercept(PG_FUNCTION_ARGS)
Definition: float.c:3600
int pg_strfromd(char *str, size_t count, int precision, double value)
Definition: snprintf.c:1265
static float4 get_float4_nan(void)
Definition: float.h:107
static float8 one_minus_cos_60
Definition: float.c:47
Datum float4recv(PG_FUNCTION_ARGS)
Definition: float.c:307
Datum degrees(PG_FUNCTION_ARGS)
Definition: float.c:2435
Datum float4out(PG_FUNCTION_ARGS)
Definition: float.c:287
Datum dcot(PG_FUNCTION_ARGS)
Definition: float.c:1800
float float4
Definition: c.h:490
Datum float48ne(PG_FUNCTION_ARGS)
Definition: float.c:3740
Datum float8recv(PG_FUNCTION_ARGS)
Definition: float.c:565
Datum float8_regr_avgx(PG_FUNCTION_ARGS)
Definition: float.c:3435
static double cosd_q1(double x)
Definition: float.c:2179
static double sind_0_to_30(double x)
Definition: float.c:2132
static float8 asin_0_5
Definition: float.c:48
Datum dsinh(PG_FUNCTION_ARGS)
Definition: float.c:2472
Datum float8_accum(PG_FUNCTION_ARGS)
Definition: float.c:2849
#define DatumGetFloat8(X)
Definition: postgres.h:728
#define PG_RETURN_BOOL(x)
Definition: fmgr.h:349
uintptr_t Datum
Definition: postgres.h:367
Datum float4ne(PG_FUNCTION_ARGS)
Definition: float.c:844
Datum float4smaller(PG_FUNCTION_ARGS)
Definition: float.c:643
Datum dexp(PG_FUNCTION_ARGS)
Definition: float.c:1572
Datum float84lt(PG_FUNCTION_ARGS)
Definition: float.c:3806
#define PG_GETARG_INT16(n)
Definition: fmgr.h:266
Datum dtand(PG_FUNCTION_ARGS)
Definition: float.c:2370
double pg_erand48(unsigned short xseed[3])
Definition: erand48.c:88
Datum float84mi(PG_FUNCTION_ARGS)
Definition: float.c:3695
Datum width_bucket_float8(PG_FUNCTION_ARGS)
Definition: float.c:3856
int extra_float_digits
Definition: float.c:42
Datum float8_var_pop(PG_FUNCTION_ARGS)
Definition: float.c:3019
static bool float4_lt(const float4 val1, const float4 val2)
Definition: float.h:308
Datum float8_regr_accum(PG_FUNCTION_ARGS)
Definition: float.c:3128
static bool float8_gt(const float8 val1, const float8 val2)
Definition: float.h:338
#define PG_RETURN_VOID()
Definition: fmgr.h:339
#define Float8GetDatumFast(X)
Definition: postgres.h:776
bool pg_strong_random(void *buf, size_t len)
Datum btfloat8sortsupport(PG_FUNCTION_ARGS)
Definition: float.c:1001
static bool float4_ge(const float4 val1, const float4 val2)
Definition: float.h:344
float8 pq_getmsgfloat8(StringInfo msg)
Definition: pqformat.c:490
static bool float4_ne(const float4 val1, const float4 val2)
Definition: float.h:296
#define DatumGetFloat4(X)
Definition: postgres.h:680
static bool pg_add_s32_overflow(int32 a, int32 b, int32 *result)
Definition: int.h:104
#define PG_RETURN_CSTRING(x)
Definition: fmgr.h:352
Datum float8smaller(PG_FUNCTION_ARGS)
Definition: float.c:710
static float8 float8_div(const float8 val1, const float8 val2)
Definition: float.h:260
Datum float4lt(PG_FUNCTION_ARGS)
Definition: float.c:853
#define newval
Datum float84gt(PG_FUNCTION_ARGS)
Definition: float.c:3824
int double_to_shortest_decimal_buf(double f, char *result)
Definition: d2s.c:1053
int AggCheckCallContext(FunctionCallInfo fcinfo, MemoryContext *aggcontext)
Definition: nodeAgg.c:3572
float4 pq_getmsgfloat4(StringInfo msg)
Definition: pqformat.c:471
static bool float4_eq(const float4 val1, const float4 val2)
Definition: float.h:284
Datum i2tof(PG_FUNCTION_ARGS)
Definition: float.c:1368
Datum float48mul(PG_FUNCTION_ARGS)
Definition: float.c:3662
#define ARR_NDIM(a)
Definition: array.h:278
Datum dasind(PG_FUNCTION_ARGS)
Definition: float.c:2031
#define RETURN_ERROR(throw_error, have_error)
Definition: float.c:340
int is_infinite(double val)
Definition: float.c:97
Datum float8_corr(PG_FUNCTION_ARGS)
Definition: float.c:3511
static bool float8_eq(const float8 val1, const float8 val2)
Definition: float.h:290
Datum dacosh(PG_FUNCTION_ARGS)
Definition: float.c:2561
Datum dasinh(PG_FUNCTION_ARGS)
Definition: float.c:2543
Datum dsqrt(PG_FUNCTION_ARGS)
Definition: float.c:1464
Datum datan2d(PG_FUNCTION_ARGS)
Definition: float.c:2096
Datum i4tod(PG_FUNCTION_ARGS)
Definition: float.c:1268
Datum float8_combine(PG_FUNCTION_ARGS)
Definition: float.c:2756
static float8 * check_float8_array(ArrayType *transarray, const char *caller, int n)
Definition: float.c:2732
Datum float4gt(PG_FUNCTION_ARGS)
Definition: float.c:871
Datum btfloat84cmp(PG_FUNCTION_ARGS)
Definition: float.c:1020
void * palloc(Size size)
Definition: mcxt.c:949
int errmsg(const char *fmt,...)
Definition: elog.c:784
static bool float4_le(const float4 val1, const float4 val2)
Definition: float.h:320
static float4 float4_pl(const float4 val1, const float4 val2)
Definition: float.h:176
Datum float8_stddev_samp(PG_FUNCTION_ARGS)
Definition: float.c:3085
static int btfloat4fastcmp(Datum x, Datum y, SortSupport ssup)
Definition: float.c:898
Datum float4send(PG_FUNCTION_ARGS)
Definition: float.c:318
#define elog(elevel,...)
Definition: elog.h:226
double float8in_internal_opt_error(char *num, char **endptr_p, const char *type_name, const char *orig_string, bool *have_error)
Definition: float.c:372
Datum float84le(PG_FUNCTION_ARGS)
Definition: float.c:3815
static float4 float4_mul(const float4 val1, const float4 val2)
Definition: float.h:220
void * arg
#define unlikely(x)
Definition: c.h:208
#define PG_GETARG_CSTRING(n)
Definition: fmgr.h:272
Datum float8out(PG_FUNCTION_ARGS)
Definition: float.c:531
Datum float48div(PG_FUNCTION_ARGS)
Definition: float.c:3671
#define PG_FUNCTION_ARGS
Definition: fmgr.h:188
Datum dcotd(PG_FUNCTION_ARGS)
Definition: float.c:2251
static float8 tan_45
Definition: float.c:51
Datum float8le(PG_FUNCTION_ARGS)
Definition: float.c:956
Datum float8_regr_avgy(PG_FUNCTION_ARGS)
Definition: float.c:3454
Datum float84eq(PG_FUNCTION_ARGS)
Definition: float.c:3788
Datum dacosd(PG_FUNCTION_ARGS)
Definition: float.c:1996
#define ARR_ELEMTYPE(a)
Definition: array.h:280
Datum float48gt(PG_FUNCTION_ARGS)
Definition: float.c:3767
long val
Definition: informix.c:684
Datum now(PG_FUNCTION_ARGS)
Definition: timestamp.c:1547
Datum float4um(PG_FUNCTION_ARGS)
Definition: float.c:611
#define PG_RETURN_NULL()
Definition: fmgr.h:335
Datum i2tod(PG_FUNCTION_ARGS)
Definition: float.c:1280
#define M_PI
Definition: earthdistance.c:10
Datum dtan(PG_FUNCTION_ARGS)
Definition: float.c:1853
Datum in_range_float8_float8(PG_FUNCTION_ARGS)
Definition: float.c:1036
Datum btfloat4sortsupport(PG_FUNCTION_ARGS)
Definition: float.c:907
static float4 get_float4_infinity(void)
Definition: float.h:70