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

Generated by: LCOV version 1.13