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