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