Re: gamma() and lgamma() functions

On Fri, 23 Aug 2024 at 13:40, Peter Eisentraut <peter@eisentraut.org> wrote:
>
> What are examples of where this would be useful in a database context?

gamma() and lgamma() are the kinds of functions that are generally
useful for a variety of tasks like statistical analysis and
combinatorial computations, and having them in the database allows
those sort of computations to be performed without the need to export
the data to an external tool. We discussed adding them in a thread
last year [1], and there has been at least a little prior interest
[2].

[1] https://www.postgresql.org/message-id/CA%2BhUKGKJAcB8Q5qziKTTSnkA4Mnv_6f%2B7-_XUgbh9jFjSdEFQg%40mail.gmail.com
[2] https://stackoverflow.com/questions/58884066/how-can-i-run-the-equivalent-of-excels-gammaln-function-in-postgres

Of course, there's a somewhat fuzzy line between what is generally
useful enough, and what is too specialised for core Postgres, but I
would argue that these qualify, since they are part of C99, and
commonly appear in other general-purpose math libraries like the
Python math module.

> > The error-handling for these functions seems to be a little trickier
> > than most, so that might need further discussion.
>
> Yeah, this is quite something.
>
> I'm not sure why you are doing the testing for special values (NaN etc.)
> yourself when the C library function already does it.  For example, if I
> remove the isnan(arg1) check in your dlgamma(), then it still behaves
> the same in all tests.

It's useful to do that so that we don't need to assume that every
platform conforms to the POSIX standard, and because it simplifies the
later overflow checks. This is consistent with the approach taken in
other functions, such as dexp(), dsin(), dcos(), etc.

> The overflow checks after the C library call are written
> differently for the two functions.  dgamma() does not check errno for
> ERANGE for example.  It might also be good if dgamma() checked errno for
> EDOM, because other the result of gamma(-1) being "overflow" seems a bit
> wrong.

They're intentionally different because the functions themselves are
different. In this case:

select gamma(-1);
ERROR:  value out of range: overflow

it is correct to throw an error, because gamma(-1) is undefined (it
goes to -Inf as x goes to -1 from above, and +Inf as x goes to -1 from
below, so there is no well-defined limit).

I've updated the patch to give a more specific error message for
negative integer inputs, as opposed to other overflow cases.

Relying on errno being ERANGE or EDOM doesn't seem possible though,
because the docs indicate that, while its behaviour is one thing
today, per POSIX, that will change in the future.

By contrast, lgamma() does not raise an error for such inputs:

select lgamma(-1);
  lgamma
----------
 Infinity

This is correct because lgamma() is the log of the absolute value of
the gamma function, so the limit is +Inf as x goes to -1 from both
sides.

> I'm also not clear why you are turning a genuine result overflow into
> infinity in lgamma().

OK, I've changed that to only return infinity if the input is
infinite, zero, or a negative integer. Otherwise, it now throws an
overflow error.

> Btw., I'm reading that lgamma() in the C library is not necessarily
> thread-safe.  Is that a possible problem?

It's not quite clear what to do about that. Some platforms may define
the lgamma_r() function, but not all. Do we have a standard pattern
for dealing with functions for which there is no thread-safe
alternative?

Regards,
Dean