Mercurial > octave-nkf
annotate src/DLD-FUNCTIONS/gammainc.cc @ 7789:82be108cc558
First attempt at single precision tyeps
* * *
corrections to qrupdate single precision routines
* * *
prefer demotion to single over promotion to double
* * *
Add single precision support to log2 function
* * *
Trivial PROJECT file update
* * *
Cache optimized hermitian/transpose methods
* * *
Add tests for tranpose/hermitian and ChangeLog entry for new transpose code
| author | David Bateman <dbateman@free.fr> |
|---|---|
| date | Sun, 27 Apr 2008 22:34:17 +0200 |
| parents | a1dbe9d80eee |
| children | a41df65f3f00 |
| rev | line source |
|---|---|
| 3254 | 1 /* |
| 2 | |
| 7017 | 3 Copyright (C) 1997, 1999, 2000, 2004, 2005, 2006, 2007 John W. Eaton |
| 3254 | 4 |
| 5 This file is part of Octave. | |
| 6 | |
| 7 Octave is free software; you can redistribute it and/or modify it | |
| 8 under the terms of the GNU General Public License as published by the | |
| 7016 | 9 Free Software Foundation; either version 3 of the License, or (at your |
| 10 option) any later version. | |
| 3254 | 11 |
| 12 Octave is distributed in the hope that it will be useful, but WITHOUT | |
| 13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
| 14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
| 15 for more details. | |
| 16 | |
| 17 You should have received a copy of the GNU General Public License | |
| 7016 | 18 along with Octave; see the file COPYING. If not, see |
| 19 <http://www.gnu.org/licenses/>. | |
| 3254 | 20 |
| 21 */ | |
| 22 | |
| 23 #ifdef HAVE_CONFIG_H | |
| 24 #include <config.h> | |
| 25 #endif | |
| 26 | |
| 27 #include "lo-specfun.h" | |
| 28 | |
| 29 #include "defun-dld.h" | |
| 30 #include "error.h" | |
| 31 #include "gripes.h" | |
| 32 #include "oct-obj.h" | |
| 33 #include "utils.h" | |
| 34 | |
| 35 DEFUN_DLD (gammainc, args, , | |
| 3548 | 36 "-*- texinfo -*-\n\ |
| 3444 | 37 @deftypefn {Mapping Function} {} gammainc (@var{x}, @var{a})\n\ |
| 5517 | 38 Compute the normalized incomplete gamma function,\n\ |
| 3444 | 39 @iftex\n\ |
| 40 @tex\n\ | |
| 41 $$\n\ | |
| 3683 | 42 \\gamma (x, a) = {\\displaystyle\\int_0^x e^{-t} t^{a-1} dt \\over \\Gamma (a)}\n\ |
| 3444 | 43 $$\n\ |
| 44 @end tex\n\ | |
| 45 @end iftex\n\ | |
| 46 @ifinfo\n\ | |
| 3254 | 47 \n\ |
| 3444 | 48 @smallexample\n\ |
| 3683 | 49 x\n\ |
| 50 1 /\n\ | |
| 51 gammainc (x, a) = --------- | exp (-t) t^(a-1) dt\n\ | |
| 52 gamma (a) /\n\ | |
| 53 t=0\n\ | |
| 3444 | 54 @end smallexample\n\ |
| 5517 | 55 \n\ |
| 3444 | 56 @end ifinfo\n\ |
| 5517 | 57 with the limiting value of 1 as @var{x} approaches infinity.\n\ |
| 58 The standard notation is @math{P(a,x)}, e.g. Abramowitz and Stegun (6.5.1).\n\ | |
| 3254 | 59 \n\ |
| 3683 | 60 If @var{a} is scalar, then @code{gammainc (@var{x}, @var{a})} is returned\n\ |
| 3444 | 61 for each element of @var{x} and vice versa.\n\ |
| 62 \n\ | |
| 3683 | 63 If neither @var{x} nor @var{a} is scalar, the sizes of @var{x} and\n\ |
| 64 @var{a} must agree, and @var{gammainc} is applied element-by-element.\n\ | |
| 5642 | 65 @seealso{gamma, lgamma}\n\ |
| 66 @end deftypefn") | |
| 3254 | 67 { |
| 68 octave_value retval; | |
| 69 | |
| 70 int nargin = args.length (); | |
| 71 | |
| 72 if (nargin == 2) | |
| 73 { | |
| 74 octave_value x_arg = args(0); | |
| 75 octave_value a_arg = args(1); | |
| 76 | |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
77 // FIXME Can we make a template version of the duplicated code below |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
78 if (x_arg.is_single_type () || a_arg.is_single_type ()) |
| 3254 | 79 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
80 if (x_arg.is_scalar_type ()) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
81 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
82 float x = x_arg.float_value (); |
| 3254 | 83 |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
84 if (! error_state) |
| 3254 | 85 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
86 if (a_arg.is_scalar_type ()) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
87 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
88 float a = a_arg.float_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
89 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
90 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
91 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
92 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
93 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
94 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
95 FloatNDArray a = a_arg.float_array_value (); |
| 3254 | 96 |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
97 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
98 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
99 } |
| 3254 | 100 } |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
101 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
102 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
103 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
104 FloatNDArray x = x_arg.float_array_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
105 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
106 if (! error_state) |
| 3254 | 107 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
108 if (a_arg.is_scalar_type ()) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
109 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
110 float a = a_arg.float_value (); |
| 3254 | 111 |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
112 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
113 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
114 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
115 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
116 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
117 FloatNDArray a = a_arg.float_array_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
118 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
119 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
120 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
121 } |
| 3254 | 122 } |
| 123 } | |
| 124 } | |
| 125 else | |
| 126 { | |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
127 if (x_arg.is_scalar_type ()) |
| 3254 | 128 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
129 double x = x_arg.double_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
130 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
131 if (! error_state) |
| 3254 | 132 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
133 if (a_arg.is_scalar_type ()) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
134 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
135 double a = a_arg.double_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
136 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
137 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
138 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
139 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
140 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
141 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
142 NDArray a = a_arg.array_value (); |
| 3254 | 143 |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
144 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
145 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
146 } |
| 3254 | 147 } |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
148 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
149 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
150 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
151 NDArray x = x_arg.array_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
152 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
153 if (! error_state) |
| 3254 | 154 { |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
155 if (a_arg.is_scalar_type ()) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
156 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
157 double a = a_arg.double_value (); |
| 3254 | 158 |
|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
159 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
160 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
161 } |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
162 else |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
163 { |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
164 NDArray a = a_arg.array_value (); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
165 |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
166 if (! error_state) |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
167 retval = gammainc (x, a); |
|
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
7017
diff
changeset
|
168 } |
| 3254 | 169 } |
| 170 } | |
| 171 } | |
| 172 } | |
| 173 else | |
| 5823 | 174 print_usage (); |
| 3254 | 175 |
| 176 return retval; | |
| 177 } | |
| 178 | |
| 179 /* | |
| 180 ;;; Local Variables: *** | |
| 181 ;;; mode: C++ *** | |
| 182 ;;; End: *** | |
| 183 */ | |
| 184 |