[Octave-bug-tracker] [bug #65230] Poor accuracy of betaln for high ratios of a/b

[anonymous]

Follow-up Comment #5, bug#65230 (group octave):

OP here again.  Thank you for the suggestions involving the symbolic package.

> The accuracy is limited not by the gammaln() function, but by the addition
operation using IEEE 8-byte doubles.  eps() is 2e-16 and if the difference in
scale between two arguments is greater than 1e16 then the smaller argument is
lost.

The accuracy of betaln(a, b) can be poor even if the addition operation a + b
is exact.  For example betaln(2^62, 2^10) returns -32768 for me, but the
actual value is about -37935.2477864453572569.



I am interested in the possibility of contributing a version of betaln.m that
improves accuracy at some cost to speed and simplicity.  I have a few
questions:
1. How likely is this to be accepted?  Assume for the sake of argument that
accuracy is at least 40 bits for the vast majority of the input space, and
that speed and simplicity are about 10 times worse than the existing
implementation.
2. The existing implementation has input validation that requires the
arguments to have equal size or at least one of them to be scalar.  Is there
any reason to retain this and not fully support broadcasting?
3. Some definitions of the beta function seem to require (the real parts of)
both arguments to be positive, but others relax this requirement.  What should
betaln return with one or more non-positive arguments?

Any input would be appreciated.  Thanks.



    _______________________________________________________

Reply to this item at:

  <https://savannah.gnu.org/bugs/?65230>

_______________________________________________
Message sent via Savannah
https://savannah.gnu.org/