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