LCOV - code coverage report
Current view: top level - src/backend/utils/adt - float.c (source / functions) Hit Total Coverage
Test: PostgreSQL 15beta1 Lines: 1066 1328 80.3 %
Date: 2022-05-18 03:10:05 Functions: 137 155 88.4 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.14