[Octave-bug-tracker] [bug #47800] gammainc(x, a, "upper") rounds down to zero if output is below eps

[Nir Krakauer]

Follow-up Comment #33, bug #47800 (project octave):

Yes, unfortunately the new implementation can be very slow:

m = 5;
n = 100;
a = rand (m, n);
x = rand (m, n);
tic
for i = 1:n
  y = gammainc (x(:, i), a(:, i));
endfor
t = toc #old: 0.4 s, new: 19 s

Profiling shows that the main bottleneck is the ngammainc_l subfunction, where
the while loop can run thousands of times before the tolerance is reached. I
streamlined this loop a bit (attached), but it's seemingly impossible to relax
the tolerance without losing accuracy (tests start to fail).

(file #37329)
    _______________________________________________________

Additional Item Attachment:

File name: gammainc.m                     Size:10 KB


    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?47800>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/