Mercurial > hg > octave-lyh
changeset 7410:8a3b2ccc4e11
[project @ 2008-01-22 21:34:24 by jwe]
author | jwe |
---|---|
date | Tue, 22 Jan 2008 21:34:24 +0000 |
parents | 73036cdd855d |
children | 83a8781b529d |
files | scripts/ChangeLog scripts/specfun/erfinv.m |
diffstat | 2 files changed, 15 insertions(+), 17 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ChangeLog +++ b/scripts/ChangeLog @@ -1,3 +1,8 @@ +2008-01-22 Schloegl Alois <alois.schloegl@tugraz.at> + + * specfun/erfinv.m: Replace z_old and z_new by a single variable z. + Simplify initial checks on argument values. + 2008-01-22 Michael Goffioul <michael.goffioul@gmail.com> * plot/gnuplot_drawnow.m: New function corresponding to the
--- a/scripts/specfun/erfinv.m +++ b/scripts/specfun/erfinv.m @@ -1,5 +1,6 @@ ## Copyright (C) 1995, 1996, 1997, 1999, 2000, 2002, 2004, 2005, 2006, ## 2007 Kurt Hornik +## Copyright (C) 2008 Alois Schloegl ## ## This file is part of Octave. ## @@ -44,31 +45,23 @@ x = reshape (x, nel, 1); y = zeros (nel, 1); - i = find ((x < -1) | (x > 1) | isnan(x)); - if any (i) - y(i) = NaN * ones (length (i), 1); - endif - - t = find (x == -1); - y (t) = (-Inf) * ones (size (t)); - - t = find (x == 1); - y (t) = Inf * ones (size (t)); + ## x < 1 or x > 1 ==> NaN + y(abs (x) >= 1) = NaN; + y(x == -1) = -Inf; + y(x == +1) = +Inf; i = find ((x > -1) & (x < 1)); - if any (i) + if (any (i)) s = sqrt (pi) / 2; - z_old = ones (length (i), 1); - z_new = sqrt (-log (1 - abs (x(i)))) .* sign (x(i)); - while (any (abs (erf (z_new) - x(i)) > tol * abs (x(i)))) - z_old = z_new; - z_new = z_old - (erf (z_old) - x(i)) .* exp (z_old.^2) * s; + z = sqrt (-log (1 - abs (x(i)))) .* sign (x(i)); + while (any (abs (erf (z) - x(i)) > tol * abs (x(i)))) + z = z - (erf (z) - x(i)) .* exp (z.^2) * s; if (++iterations > maxit) warning ("erfinv: iteration limit exceeded"); break; endif endwhile - y(i) = z_new; + y(i) = z; endif y = reshape (y, sz);