1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2016 Gael Guennebaud <[email protected]>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9
10
11 #ifndef EIGEN_SPECIALFUNCTIONS_ARRAYAPI_H
12 #define EIGEN_SPECIALFUNCTIONS_ARRAYAPI_H
13
14 namespace Eigen {
15
16 /** \cpp11 \returns an expression of the coefficient-wise igamma(\a a, \a x) to the given arrays.
17 *
18 * This function computes the coefficient-wise incomplete gamma function.
19 *
20 * \note This function supports only float and double scalar types in c++11 mode. To support other scalar types,
21 * or float/double in non c++11 mode, the user has to provide implementations of igammac(T,T) for any scalar
22 * type T to be supported.
23 *
24 * \sa Eigen::igammac(), Eigen::lgamma()
25 */
26 template<typename Derived,typename ExponentDerived>
27 EIGEN_STRONG_INLINE 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)28 igamma(const Eigen::ArrayBase<Derived>& a, const Eigen::ArrayBase<ExponentDerived>& x)
29 {
30 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_igamma_op<typename Derived::Scalar>, const Derived, const ExponentDerived>(
31 a.derived(),
32 x.derived()
33 );
34 }
35
36 /** \cpp11 \returns an expression of the coefficient-wise igamma_der_a(\a a, \a x) to the given arrays.
37 *
38 * This function computes the coefficient-wise derivative of the incomplete
39 * gamma function with respect to the parameter a.
40 *
41 * \note This function supports only float and double scalar types in c++11
42 * mode. To support other scalar types,
43 * or float/double in non c++11 mode, the user has to provide implementations
44 * of igamma_der_a(T,T) for any scalar
45 * type T to be supported.
46 *
47 * \sa Eigen::igamma(), Eigen::lgamma()
48 */
49 template <typename Derived, typename ExponentDerived>
50 EIGEN_STRONG_INLINE const Eigen::CwiseBinaryOp<Eigen::internal::scalar_igamma_der_a_op<typename Derived::Scalar>, const Derived, const ExponentDerived>
igamma_der_a(const Eigen::ArrayBase<Derived> & a,const Eigen::ArrayBase<ExponentDerived> & x)51 igamma_der_a(const Eigen::ArrayBase<Derived>& a, const Eigen::ArrayBase<ExponentDerived>& x) {
52 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_igamma_der_a_op<typename Derived::Scalar>, const Derived, const ExponentDerived>(
53 a.derived(),
54 x.derived());
55 }
56
57 /** \cpp11 \returns an expression of the coefficient-wise gamma_sample_der_alpha(\a alpha, \a sample) to the given arrays.
58 *
59 * This function computes the coefficient-wise derivative of the sample
60 * of a Gamma(alpha, 1) random variable with respect to the parameter alpha.
61 *
62 * \note This function supports only float and double scalar types in c++11
63 * mode. To support other scalar types,
64 * or float/double in non c++11 mode, the user has to provide implementations
65 * of gamma_sample_der_alpha(T,T) for any scalar
66 * type T to be supported.
67 *
68 * \sa Eigen::igamma(), Eigen::lgamma()
69 */
70 template <typename AlphaDerived, typename SampleDerived>
71 EIGEN_STRONG_INLINE const Eigen::CwiseBinaryOp<Eigen::internal::scalar_gamma_sample_der_alpha_op<typename AlphaDerived::Scalar>, const AlphaDerived, const SampleDerived>
gamma_sample_der_alpha(const Eigen::ArrayBase<AlphaDerived> & alpha,const Eigen::ArrayBase<SampleDerived> & sample)72 gamma_sample_der_alpha(const Eigen::ArrayBase<AlphaDerived>& alpha, const Eigen::ArrayBase<SampleDerived>& sample) {
73 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_gamma_sample_der_alpha_op<typename AlphaDerived::Scalar>, const AlphaDerived, const SampleDerived>(
74 alpha.derived(),
75 sample.derived());
76 }
77
78 /** \cpp11 \returns an expression of the coefficient-wise igammac(\a a, \a x) to the given arrays.
79 *
80 * This function computes the coefficient-wise complementary incomplete gamma function.
81 *
82 * \note This function supports only float and double scalar types in c++11 mode. To support other scalar types,
83 * or float/double in non c++11 mode, the user has to provide implementations of igammac(T,T) for any scalar
84 * type T to be supported.
85 *
86 * \sa Eigen::igamma(), Eigen::lgamma()
87 */
88 template<typename Derived,typename ExponentDerived>
89 EIGEN_STRONG_INLINE 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)90 igammac(const Eigen::ArrayBase<Derived>& a, const Eigen::ArrayBase<ExponentDerived>& x)
91 {
92 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_igammac_op<typename Derived::Scalar>, const Derived, const ExponentDerived>(
93 a.derived(),
94 x.derived()
95 );
96 }
97
98 /** \cpp11 \returns an expression of the coefficient-wise polygamma(\a n, \a x) to the given arrays.
99 *
100 * It returns the \a n -th derivative of the digamma(psi) evaluated at \c x.
101 *
102 * \note This function supports only float and double scalar types in c++11 mode. To support other scalar types,
103 * or float/double in non c++11 mode, the user has to provide implementations of polygamma(T,T) for any scalar
104 * type T to be supported.
105 *
106 * \sa Eigen::digamma()
107 */
108 // * \warning Be careful with the order of the parameters: x.polygamma(n) is equivalent to polygamma(n,x)
109 // * \sa ArrayBase::polygamma()
110 template<typename DerivedN,typename DerivedX>
111 EIGEN_STRONG_INLINE 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)112 polygamma(const Eigen::ArrayBase<DerivedN>& n, const Eigen::ArrayBase<DerivedX>& x)
113 {
114 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_polygamma_op<typename DerivedX::Scalar>, const DerivedN, const DerivedX>(
115 n.derived(),
116 x.derived()
117 );
118 }
119
120 /** \cpp11 \returns an expression of the coefficient-wise betainc(\a x, \a a, \a b) to the given arrays.
121 *
122 * This function computes the regularized incomplete beta function (integral).
123 *
124 * \note This function supports only float and double scalar types in c++11 mode. To support other scalar types,
125 * or float/double in non c++11 mode, the user has to provide implementations of betainc(T,T,T) for any scalar
126 * type T to be supported.
127 *
128 * \sa Eigen::betainc(), Eigen::lgamma()
129 */
130 template<typename ArgADerived, typename ArgBDerived, typename ArgXDerived>
131 EIGEN_STRONG_INLINE const Eigen::CwiseTernaryOp<Eigen::internal::scalar_betainc_op<typename ArgXDerived::Scalar>, const ArgADerived, const ArgBDerived, const ArgXDerived>
betainc(const Eigen::ArrayBase<ArgADerived> & a,const Eigen::ArrayBase<ArgBDerived> & b,const Eigen::ArrayBase<ArgXDerived> & x)132 betainc(const Eigen::ArrayBase<ArgADerived>& a, const Eigen::ArrayBase<ArgBDerived>& b, const Eigen::ArrayBase<ArgXDerived>& x)
133 {
134 return Eigen::CwiseTernaryOp<Eigen::internal::scalar_betainc_op<typename ArgXDerived::Scalar>, const ArgADerived, const ArgBDerived, const ArgXDerived>(
135 a.derived(),
136 b.derived(),
137 x.derived()
138 );
139 }
140
141
142 /** \returns an expression of the coefficient-wise zeta(\a x, \a q) to the given arrays.
143 *
144 * It returns the Riemann zeta function of two arguments \a x and \a q:
145 *
146 * \param x is the exponent, it must be > 1
147 * \param q is the shift, it must be > 0
148 *
149 * \note This function supports only float and double scalar types. To support other scalar types, the user has
150 * to provide implementations of zeta(T,T) for any scalar type T to be supported.
151 *
152 * \sa ArrayBase::zeta()
153 */
154 template<typename DerivedX,typename DerivedQ>
155 EIGEN_STRONG_INLINE const Eigen::CwiseBinaryOp<Eigen::internal::scalar_zeta_op<typename DerivedX::Scalar>, const DerivedX, const DerivedQ>
zeta(const Eigen::ArrayBase<DerivedX> & x,const Eigen::ArrayBase<DerivedQ> & q)156 zeta(const Eigen::ArrayBase<DerivedX>& x, const Eigen::ArrayBase<DerivedQ>& q)
157 {
158 return Eigen::CwiseBinaryOp<Eigen::internal::scalar_zeta_op<typename DerivedX::Scalar>, const DerivedX, const DerivedQ>(
159 x.derived(),
160 q.derived()
161 );
162 }
163
164
165 } // end namespace Eigen
166
167 #endif // EIGEN_SPECIALFUNCTIONS_ARRAYAPI_H
168