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