Mercurial > hg > octave-lyh
diff scripts/specfun/gammai.m @ 1586:54350c98055c
[project @ 1995-10-31 05:22:53 by jwe]
author | jwe |
---|---|
date | Tue, 31 Oct 1995 05:24:40 +0000 |
parents | 100413a7e8a2 |
children | 5d29638dd524 |
line wrap: on
line diff
--- a/scripts/specfun/gammai.m +++ b/scripts/specfun/gammai.m @@ -80,7 +80,7 @@ error ("gammai: all entries of x must be nonnegative"); endif - y = zeros(1, n); + y = zeros (1, n); # For x < a + 1, use summation. The below choice of k should ensure # that the overall error is less than eps ... @@ -90,8 +90,8 @@ if (s > 0) k = ceil (- max ([a(S), x(S)]) * log (eps)); K = (1:k)'; - M = ones(k, 1); - A = cumprod((M * x(S)) ./ (M * a(S) + K * ones(1, s))); + M = ones (k, 1); + A = cumprod ((M * x(S)) ./ (M * a(S) + K * ones(1, s))); y(S) = exp (-x(S) + a(S) .* log (x(S))) .* (1 + sum (A)) ./ gamma (a(S)+1); endif @@ -99,11 +99,13 @@ # Note, however, that this converges MUCH slower than the series # expansion for small a and x not too large! - S = find((x >= a + 1) & (x < Inf)); - s = length(S); + S = find ((x >= a + 1) & (x < Inf)); + s = length (S); if (s > 0) - u = [zeros(1, s); ones(1, s)]; - v = [ones(1, s); x(S)]; + t1 = zeros (1, s); + t2 = ones (1, s); + u = [t1; t2]; + v = [t2; x(S)]; c_old = 0; c_new = v(1,:) ./ v(2,:); n = 1;