From 4613ddc7267343e74703991267a1adb45eff3403 Mon Sep 17 00:00:00 2001 From: Dean Rasheed Date: Sun, 2 Mar 2025 09:23:20 +0000 Subject: [PATCH v5] Add support for gamma() and lgamma() functions. These are useful general-purpose math functions which are included in POSIX and C99, and are commonly included in other math libraries, so expose them as SQL-callable functions. Author: Dean Rasheed Reviewed-by: Stepan Neretin Reviewed-by: Tom Lane Reviewed-by: Peter Eisentraut Reviewed-by: Dmitry Koval Discussion: https://postgr.es/m/CAEZATCXpGyfjXCirFk9au+FvM0y2Ah+2-0WSJx7MO368ysNUPA@mail.gmail.com --- doc/src/sgml/func.sgml | 38 ++++++++++++ src/backend/utils/adt/float.c | 87 ++++++++++++++++++++++++++++ src/include/catalog/pg_proc.dat | 7 +++ src/test/regress/expected/float8.out | 57 ++++++++++++++++++ src/test/regress/sql/float8.sql | 23 ++++++++ 5 files changed, 212 insertions(+) diff --git a/doc/src/sgml/func.sgml b/doc/src/sgml/func.sgml index 0e6c534965..143031e03e 100644 --- a/doc/src/sgml/func.sgml +++ b/doc/src/sgml/func.sgml @@ -1396,6 +1396,27 @@ SELECT NOT(ROW(table.*) IS NOT NULL) FROM TABLE; -- detect at least one null in + + + + gamma + + gamma ( double precision ) + double precision + + + Gamma function + + + gamma(0.5) + 1.772453850905516 + + + gamma(6) + 120 + + + @@ -1436,6 +1457,23 @@ SELECT NOT(ROW(table.*) IS NOT NULL) FROM TABLE; -- detect at least one null in + + + + lgamma + + lgamma ( double precision ) + double precision + + + Natural logarithm of the absolute value of the gamma function + + + lgamma(1000) + 5905.220423209181 + + + diff --git a/src/backend/utils/adt/float.c b/src/backend/utils/adt/float.c index 2bc31eabf2..a46ebffae4 100644 --- a/src/backend/utils/adt/float.c +++ b/src/backend/utils/adt/float.c @@ -2786,6 +2786,93 @@ derfc(PG_FUNCTION_ARGS) } +/* ========== GAMMA FUNCTIONS ========== */ + + +/* + * dgamma - returns the gamma function of arg1 + */ +Datum +dgamma(PG_FUNCTION_ARGS) +{ + float8 arg1 = PG_GETARG_FLOAT8(0); + float8 result; + + /* + * Handle NaN and Inf cases explicitly. This simplifies the overflow + * checks on platforms that do not set errno. + */ + if (isnan(arg1)) + result = arg1; + else if (isinf(arg1)) + { + /* Per POSIX, an input of -Inf causes a domain error */ + if (arg1 < 0) + { + float_overflow_error(); + result = get_float8_nan(); /* keep compiler quiet */ + } + else + result = arg1; + } + else + { + /* + * Note: the POSIX/C99 gamma function is called "tgamma", not "gamma". + * + * On some platforms, tgamma() will not set errno but just return Inf, + * NaN, or zero to report overflow/underflow; therefore, test those + * cases explicitly (note that, like the exponential function, the + * gamma function has no zeros). + */ + errno = 0; + result = tgamma(arg1); + + if (errno != 0 || isinf(result) || isnan(result)) + { + if (result != 0.0) + float_overflow_error(); + else + float_underflow_error(); + } + else if (result == 0.0) + float_underflow_error(); + } + + PG_RETURN_FLOAT8(result); +} + + +/* + * dlgamma - natural logarithm of absolute value of gamma of arg1 + */ +Datum +dlgamma(PG_FUNCTION_ARGS) +{ + float8 arg1 = PG_GETARG_FLOAT8(0); + float8 result; + + /* + * Note: lgamma may not be thread-safe because it may write to a global + * variable signgam, which may not be thread-local. However, this doesn't + * matter to us, since we don't use signgam. + */ + errno = 0; + result = lgamma(arg1); + + /* + * If an ERANGE error occurs, it means there is an overflow. + * + * On some platforms, lgamma() will not set errno but just return infinity + * to report overflow, but it should never underflow. + */ + if (errno == ERANGE || (isinf(result) && !isinf(arg1))) + float_overflow_error(); + + PG_RETURN_FLOAT8(result); +} + + /* * ========================= diff --git a/src/include/catalog/pg_proc.dat b/src/include/catalog/pg_proc.dat index cd9422d0ba..2149fa4281 100644 --- a/src/include/catalog/pg_proc.dat +++ b/src/include/catalog/pg_proc.dat @@ -3542,6 +3542,13 @@ proname => 'erfc', prorettype => 'float8', proargtypes => 'float8', prosrc => 'derfc' }, +{ oid => '8702', descr => 'gamma function', + proname => 'gamma', prorettype => 'float8', proargtypes => 'float8', + prosrc => 'dgamma' }, +{ oid => '8703', descr => 'natural logarithm of absolute value of gamma function', + proname => 'lgamma', prorettype => 'float8', proargtypes => 'float8', + prosrc => 'dlgamma' }, + { oid => '1618', proname => 'interval_mul', prorettype => 'interval', proargtypes => 'interval float8', prosrc => 'interval_mul' }, diff --git a/src/test/regress/expected/float8.out b/src/test/regress/expected/float8.out index 9ef9793fe9..10a5a6e1b6 100644 --- a/src/test/regress/expected/float8.out +++ b/src/test/regress/expected/float8.out @@ -829,6 +829,63 @@ FROM (VALUES (float8 '-infinity'), NaN | NaN | NaN (22 rows) +RESET extra_float_digits; +-- gamma functions +-- we run these with extra_float_digits = -1, to get consistently rounded +-- results on all platforms. +SET extra_float_digits = -1; +SELECT x, + gamma(x), + lgamma(x) +FROM (VALUES (0.5), (1), (2), (3), (4), (5), + (float8 'infinity'), (float8 'nan')) AS t(x); + x | gamma | lgamma +----------+-----------------+------------------ + 0.5 | 1.7724538509055 | 0.5723649429247 + 1 | 1 | 0 + 2 | 1 | 0 + 3 | 2 | 0.69314718055995 + 4 | 6 | 1.7917594692281 + 5 | 24 | 3.1780538303479 + Infinity | Infinity | Infinity + NaN | NaN | NaN +(8 rows) + +-- test overflow/underflow handling +SELECT gamma(float8 '-infinity'); +ERROR: value out of range: overflow +SELECT lgamma(float8 '-infinity'); + lgamma +---------- + Infinity +(1 row) + +SELECT gamma(float8 '-1000.5'); +ERROR: value out of range: underflow +SELECT lgamma(float8 '-1000.5'); + lgamma +------------------ + -5914.4377011169 +(1 row) + +SELECT gamma(float8 '-1'); +ERROR: value out of range: overflow +SELECT lgamma(float8 '-1'); +ERROR: value out of range: overflow +SELECT gamma(float8 '0'); +ERROR: value out of range: overflow +SELECT lgamma(float8 '0'); +ERROR: value out of range: overflow +SELECT gamma(float8 '1000'); +ERROR: value out of range: overflow +SELECT lgamma(float8 '1000'); + lgamma +----------------- + 5905.2204232092 +(1 row) + +SELECT lgamma(float8 '1e308'); +ERROR: value out of range: overflow RESET extra_float_digits; -- test for over- and underflow INSERT INTO FLOAT8_TBL(f1) VALUES ('10e400'); diff --git a/src/test/regress/sql/float8.sql b/src/test/regress/sql/float8.sql index 81a35e0bf1..db8d5724c2 100644 --- a/src/test/regress/sql/float8.sql +++ b/src/test/regress/sql/float8.sql @@ -245,6 +245,29 @@ FROM (VALUES (float8 '-infinity'), RESET extra_float_digits; +-- gamma functions +-- we run these with extra_float_digits = -1, to get consistently rounded +-- results on all platforms. +SET extra_float_digits = -1; +SELECT x, + gamma(x), + lgamma(x) +FROM (VALUES (0.5), (1), (2), (3), (4), (5), + (float8 'infinity'), (float8 'nan')) AS t(x); +-- test overflow/underflow handling +SELECT gamma(float8 '-infinity'); +SELECT lgamma(float8 '-infinity'); +SELECT gamma(float8 '-1000.5'); +SELECT lgamma(float8 '-1000.5'); +SELECT gamma(float8 '-1'); +SELECT lgamma(float8 '-1'); +SELECT gamma(float8 '0'); +SELECT lgamma(float8 '0'); +SELECT gamma(float8 '1000'); +SELECT lgamma(float8 '1000'); +SELECT lgamma(float8 '1e308'); +RESET extra_float_digits; + -- test for over- and underflow INSERT INTO FLOAT8_TBL(f1) VALUES ('10e400'); -- 2.43.0