11 #include "../Eigen/SpecialFunctions" 13 template<
typename X,
typename Y>
34 Scalar plusinf = std::numeric_limits<Scalar>::infinity();
35 Scalar nan = std::numeric_limits<Scalar>::quiet_NaN();
37 Index rows = internal::random<Index>(1,30);
42 ArrayType
m1 = ArrayType::Random(rows,cols);
43 #if EIGEN_HAS_C99_MATH 48 #endif // EIGEN_HAS_C99_MATH 52 #if EIGEN_HAS_C99_MATH 58 ArrayType
m1 = ArrayType::Random(rows,cols);
59 ArrayType
m2 = ArrayType::Random(rows,cols);
67 ArrayType
a = m1.abs() + 2;
68 ArrayType
x = m2.abs() + 2;
69 ArrayType
zero = ArrayType::Zero(rows, cols);
70 ArrayType one = ArrayType::Constant(rows, cols,
Scalar(1.0));
71 ArrayType a_m1 = a - one;
73 ArrayType Gamma_a_m1_x =
Eigen::igammac(a_m1, x) * a_m1.lgamma().exp();
75 ArrayType gamma_a_m1_x =
Eigen::igamma(a_m1, x) * a_m1.lgamma().exp();
96 Scalar igamma_s[][6] = {{0.0, nan, nan, nan, nan, nan},
97 {0.0, 0.6321205588285578, 0.7768698398515702,
98 0.9816843611112658, 9.999500016666262e-05, 1.0},
99 {0.0, 0.4275932955291202, 0.608374823728911,
100 0.9539882943107686, 7.522076445089201e-07, 1.0},
101 {0.0, 0.01898815687615381, 0.06564245437845008,
102 0.5665298796332909, 4.166333347221828e-18, 1.0},
103 {0.0, 0.9999780593618628, 0.9999899967080838,
104 0.9999996219837988, 0.9991370418689945, 1.0},
105 {0.0, 0.0, 0.0, 0.0, 0.0, 0.5042041932513908}};
106 Scalar igammac_s[][6] = {{nan, nan, nan, nan, nan, nan},
107 {1.0, 0.36787944117144233, 0.22313016014842982,
108 0.018315638888734182, 0.9999000049998333, 0.0},
109 {1.0, 0.5724067044708798, 0.3916251762710878,
110 0.04601170568923136, 0.9999992477923555, 0.0},
111 {1.0, 0.9810118431238462, 0.9343575456215499,
112 0.4334701203667089, 1.0, 0.0},
113 {1.0, 2.1940638138146658e-05, 1.0003291916285e-05,
114 3.7801620118431334e-07, 0.0008629581310054535,
116 {1.0, 1.0, 1.0, 1.0, 1.0, 0.49579580674813944}};
117 for (
int i = 0;
i < 6; ++
i) {
118 for (
int j = 0;
j < 6; ++
j) {
134 #endif // EIGEN_HAS_C99_MATH 138 ArrayType
x(7),
q(7),
res(7),
ref(7);
139 x << 1.5, 4, 10.5, 10000.5, 3, 1, 0.9;
140 q << 2, 1.5, 3, 1.0001, -2.5, 1.2345, 1.2345;
141 ref << 1.61237534869, 0.234848505667, 1.03086757337e-5, 0.367879440865, 0.054102025820864097, plusinf, nan;
149 ArrayType
x(7),
res(7),
ref(7);
150 x << 1, 1.5, 4, -10.5, 10000.5, 0, -1;
151 ref << -0.5772156649015329, 0.03648997397857645, 1.2561176684318, 2.398239129535781, 9.210340372392849, plusinf, plusinf;
159 #if EIGEN_HAS_C99_MATH 161 ArrayType
n(11),
x(11),
res(11),
ref(11);
162 n << 1, 1, 1, 1.5, 17, 31, 28, 8, 42, 147, 170;
163 x << 2, 3, 25.5, 1.5, 4.7, 11.8, 17.7, 30.2, 15.8, 54.1, 64;
164 ref << 0.644934066848, 0.394934066848, 0.0399946696496, nan, 293.334565435, 0.445487887616, -2.47810300902e-07, -8.29668781082e-09, -0.434562276666, 0.567742190178, -0.0108615497927;
167 if(
sizeof(RealScalar)>=8) {
179 #if EIGEN_HAS_C99_MATH 196 a << 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
197 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
198 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
199 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
200 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
201 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
202 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
203 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
204 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
205 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
206 0.03062277660168379, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999,
207 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999,
208 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999, 0.999,
209 31.62177660168379, 31.62177660168379, 31.62177660168379,
210 31.62177660168379, 31.62177660168379, 31.62177660168379,
211 31.62177660168379, 31.62177660168379, 31.62177660168379,
212 31.62177660168379, 31.62177660168379, 31.62177660168379,
213 31.62177660168379, 31.62177660168379, 31.62177660168379,
214 31.62177660168379, 31.62177660168379, 31.62177660168379,
215 31.62177660168379, 31.62177660168379, 31.62177660168379,
216 31.62177660168379, 31.62177660168379, 31.62177660168379,
217 31.62177660168379, 999.999, 999.999, 999.999, 999.999, 999.999, 999.999,
218 999.999, 999.999, 999.999, 999.999, 999.999, 999.999, 999.999, 999.999,
219 999.999, 999.999, 999.999, 999.999, 999.999, 999.999, 999.999, 999.999,
220 999.999, 999.999, 999.999;
222 b << 0.0, 0.0, 0.0, 0.0, 0.0, 0.03062277660168379, 0.03062277660168379,
223 0.03062277660168379, 0.03062277660168379, 0.03062277660168379, 0.999,
224 0.999, 0.999, 0.999, 0.999, 31.62177660168379, 31.62177660168379,
225 31.62177660168379, 31.62177660168379, 31.62177660168379, 999.999,
226 999.999, 999.999, 999.999, 999.999, 0.0, 0.0, 0.0, 0.0, 0.0,
227 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
228 0.03062277660168379, 0.03062277660168379, 0.999, 0.999, 0.999, 0.999,
229 0.999, 31.62177660168379, 31.62177660168379, 31.62177660168379,
230 31.62177660168379, 31.62177660168379, 999.999, 999.999, 999.999,
231 999.999, 999.999, 0.0, 0.0, 0.0, 0.0, 0.0, 0.03062277660168379,
232 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
233 0.03062277660168379, 0.999, 0.999, 0.999, 0.999, 0.999,
234 31.62177660168379, 31.62177660168379, 31.62177660168379,
235 31.62177660168379, 31.62177660168379, 999.999, 999.999, 999.999,
236 999.999, 999.999, 0.0, 0.0, 0.0, 0.0, 0.0, 0.03062277660168379,
237 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
238 0.03062277660168379, 0.999, 0.999, 0.999, 0.999, 0.999,
239 31.62177660168379, 31.62177660168379, 31.62177660168379,
240 31.62177660168379, 31.62177660168379, 999.999, 999.999, 999.999,
241 999.999, 999.999, 0.0, 0.0, 0.0, 0.0, 0.0, 0.03062277660168379,
242 0.03062277660168379, 0.03062277660168379, 0.03062277660168379,
243 0.03062277660168379, 0.999, 0.999, 0.999, 0.999, 0.999,
244 31.62177660168379, 31.62177660168379, 31.62177660168379,
245 31.62177660168379, 31.62177660168379, 999.999, 999.999, 999.999,
248 x << -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5,
249 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2,
250 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1,
251 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1,
252 -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8,
253 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5,
254 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2,
255 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1,
256 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5, 0.8, 1.1, -0.1, 0.2, 0.5,
259 v << nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
260 nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan, nan,
261 nan, nan, nan, 0.47972119876364683, 0.5, 0.5202788012363533, nan, nan,
262 0.9518683957740043, 0.9789663010413743, 0.9931729188073435, nan, nan,
263 0.999995949033062, 0.9999999999993698, 0.9999999999999999, nan, nan,
264 0.9999999999999999, 0.9999999999999999, 0.9999999999999999, nan, nan,
265 nan, nan, nan, nan, nan, 0.006827081192655869, 0.0210336989586256,
266 0.04813160422599567, nan, nan, 0.20014344256217678, 0.5000000000000001,
267 0.7998565574378232, nan, nan, 0.9991401428435834, 0.999999999698403,
268 0.9999999999999999, nan, nan, 0.9999999999999999, 0.9999999999999999,
269 0.9999999999999999, nan, nan, nan, nan, nan, nan, nan,
270 1.0646600232370887e-25, 6.301722877826246e-13, 4.050966937974938e-06,
271 nan, nan, 7.864342668429763e-23, 3.015969667594166e-10,
272 0.0008598571564165444, nan, nan, 6.031987710123844e-08,
273 0.5000000000000007, 0.9999999396801229, nan, nan, 0.9999999999999999,
274 0.9999999999999999, 0.9999999999999999, nan, nan, nan, nan, nan, nan,
275 nan, 0.0, 7.029920380986636e-306, 2.2450728208591345e-101, nan, nan,
276 0.0, 9.275871147869727e-302, 1.2232913026152827e-97, nan, nan, 0.0,
277 3.0891393081932924e-252, 2.9303043666183996e-60, nan, nan,
278 2.248913486879199e-196, 0.5000000000004947, 0.9999999999999999, nan;
286 ArrayType
m1 = ArrayType::Random(32);
287 ArrayType
m2 = ArrayType::Random(32);
288 ArrayType m3 = ArrayType::Random(32);
289 ArrayType one = ArrayType::Constant(32,
Scalar(1.0));
291 ArrayType
a = (m1 * 4.0).exp();
292 ArrayType
b = (m2 * 4.0).exp();
293 ArrayType
x = m3.abs();
304 ArrayType expected = one - (one -
x).
pow(b);
310 ArrayType expected = one;
315 ArrayType num = x.pow(a) * (one -
x).
pow(b);
316 ArrayType denom = a * (a.lgamma() + b.lgamma() - (a +
b).
lgamma()).
exp();
319 ArrayType expected =
betainc(a, b, x) - num / denom + eps;
321 if (
sizeof(Scalar) >= 8) {
332 ArrayType num = x.pow(a) * (one -
x).
pow(b);
333 ArrayType denom = b * (a.lgamma() + b.lgamma() - (a +
b).
lgamma()).
exp();
334 ArrayType expected =
betainc(a, b, x) + num / denom + eps;
335 ArrayType test =
betainc(a, b + one, x) + eps;
343 CALL_SUBTEST_1(array_special_functions<ArrayXf>());
344 CALL_SUBTEST_2(array_special_functions<ArrayXd>());
const CwiseBinaryOp< internal::scalar_zeta_op< Scalar >, const Derived, const DerivedQ > zeta(const EIGEN_CURRENT_STORAGE_BASE_CLASS< DerivedQ > &q) const
EIGEN_DEVICE_FUNC const ExpReturnType exp() const
EIGEN_DEVICE_FUNC const ErfReturnType erf() const
EIGEN_DONT_INLINE Scalar zero()
EIGEN_DEVICE_FUNC const SqrtReturnType sqrt() const
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
const Eigen::CwiseBinaryOp< Eigen::internal::scalar_igammac_op< typename Derived::Scalar >, const Derived, const ExponentDerived > igammac(const Eigen::ArrayBase< Derived > &a, const Eigen::ArrayBase< ExponentDerived > &x)
const Eigen::CwiseBinaryOp< Eigen::internal::scalar_igamma_op< typename Derived::Scalar >, const Derived, const ExponentDerived > igamma(const Eigen::ArrayBase< Derived > &a, const Eigen::ArrayBase< ExponentDerived > &x)
void verify_component_wise(const X &x, const Y &y)
cout<< "Here is the matrix m:"<< endl<< m<< endl;Matrix< ptrdiff_t, 3, 1 > res
void array_special_functions()
internal::enable_if< !(internal::is_same< typename Derived::Scalar, ScalarExponent >::value)&&EIGEN_SCALAR_BINARY_SUPPORTED(pow, typename Derived::Scalar, ScalarExponent), const EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(Derived, ScalarExponent, pow) >::type pow(const Eigen::ArrayBase< Derived > &x, const ScalarExponent &exponent)
#define VERIFY_IS_APPROX(a, b)
#define VERIFY_IS_EQUAL(a, b)
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
EIGEN_DEVICE_FUNC const Scalar & q
void test_special_functions()
NumTraits< Scalar >::Real RealScalar
EIGEN_DEVICE_FUNC const ErfcReturnType erfc() const
#define CALL_SUBTEST(FUNC)
const Eigen::CwiseBinaryOp< Eigen::internal::scalar_polygamma_op< typename DerivedX::Scalar >, const DerivedN, const DerivedX > polygamma(const Eigen::ArrayBase< DerivedN > &n, const Eigen::ArrayBase< DerivedX > &x)
EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE const TensorCwiseTernaryOp< internal::scalar_betainc_op< typename XDerived::Scalar >, const ADerived, const BDerived, const XDerived > betainc(const ADerived &a, const BDerived &b, const XDerived &x)
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
EIGEN_DEVICE_FUNC const DigammaReturnType digamma() const
EIGEN_DEVICE_FUNC const LgammaReturnType lgamma() const