Re: Patch for erfinv.m

[John W. Eaton]

On 23-Jan-2008, Alois Schloegl wrote:

| John W. Eaton wrote:
| > On 21-Jan-2008, Alois Schloegl wrote:
| >
| > | >> Consider the following patch for erfinv.
| > | >>
| > | >> It (1) replaces z_old and z_new by a single variable z, and
| > | >> (2) makes the initial checks simpler.
| > | >>
| > | >> The code becomes leaner, and a bit faster.
| >
| > I applied the patch but used
| >
| >   ## x < 1 or x > 1 ==> NaN
| >   y(abs (x) >= 1) = NaN;
| >
| > instead of
| >
| > | +  y(~(abs(x) < 1)) = NaN;  %% x<1, x>1, x=NaN
| >
| > Any reason not to make this change?
| >   
| 
| In the latter case, x=NaN is not properly handled.
| 
| The comment line should be
| 
| - ## x < 1 or x > 1 ==> NaN
| + ## x < 1 or x > 1 or isnan(x) ==> NaN
| 
| Either use the version of the patch, or add an extra check using isnan.
| The solution in the patch is simpler, and perhaps faster, too.
| 
| Checking it again, I see that an even better solution is
| 
| - y(~(abs(x) < 1)) = NaN;  
| + y(~(abs(x) <= 1)) = NaN;  
| 
| because no NaN should be assigned in case of abs(x)==1. 

I don't think I'm following.

As I understand it, given a vector like this:

  y = [NaN, -2, -1, 0, 1, 2];

we want to convert all elements strictly outside the range [-1, 1] to
NaN?  If so, then I think

  y(abs (y) > 1) = NaN;

should work.  I see this result:

  octave:1> y = [NaN, -2, -1, 0, 1, 2];
  octave:2> y(abs (y) > 1) = NaN
  y =

     NaN   NaN    -1     0     1   NaN

So I don't see why you are using the condition <= and then negating
the result.  Why not simply use > as the test?

jwe