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

Generated by: LCOV version 1.13