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

Generated by: LCOV version 1.13