[Octave-bug-tracker] [bug #36607] norminv doesn't work for very small values

[Rik]

Update of bug #36607 (project octave):

                Category:                    None => Libraries              
              Item Group:                    None => Inaccurate Result      
                  Status:                    None => Confirmed              

    _______________________________________________________

Follow-up Comment #1:

Confirmed.  The problem can be traced to the erfinv() function.

Octave attempts to compute the inverse of the Gaussian using


inv = sqrt (2) * erfinv (2 * x - 1);


On machines using IEEE-854 8-byte double precision floating point numbers the
machine precision (eps()) is 2.2e-16.  When you have a calculation, such as
subtraction, where one nomber is of order 1 and the other is below eps, such
as 1e-17, then the smaller number is effectively rounded to zero.  This leads
to


erfinv (-1) = -Inf


Could you run the following code in Matlab to help elucidate things?


erfinv (-1 + [1e-16 1e-17 1e-18 1e-19 1e-20])
norminv ([1e-16 1e-17 1e-18 1e-19 1e-20])




    _______________________________________________________

Reply to this item at:

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

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