1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
|
namespace Eigen {
/** \page TopicCustomizingEigen Customizing/Extending Eigen
Eigen can be extended in several ways, for instance, by defining global methods, \ref ExtendingMatrixBase "by adding custom methods to MatrixBase", adding support to \ref CustomScalarType "custom types" etc.
\b Table \b of \b contents
- \ref ExtendingMatrixBase
- \ref InheritingFromMatrix
- \ref CustomScalarType
- \ref PreprocessorDirectives
\section ExtendingMatrixBase Extending MatrixBase
In this section we will see how to add custom methods to MatrixBase. Since all expressions and matrix types inherit MatrixBase, adding a method to MatrixBase make it immediately available to all expressions ! A typical use case is, for instance, to make Eigen compatible with another API.
You certainly know that in C++ it is not possible to add methods to an existing class. So how that's possible ? Here the trick is to include in the declaration of MatrixBase a file defined by the preprocessor token \c EIGEN_MATRIXBASE_PLUGIN:
\code
class MatrixBase {
// ...
#ifdef EIGEN_MATRIXBASE_PLUGIN
#include EIGEN_MATRIXBASE_PLUGIN
#endif
};
\endcode
Therefore to extend MatrixBase with you own methods you just have to create a file with your method declaration and define EIGEN_MATRIXBASE_PLUGIN before you include any Eigen's header file.
Here is an example of such an extension file: \n
\b MatrixBaseAddons.h
\code
inline Scalar at(uint i, uint j) const { return this->operator()(i,j); }
inline Scalar& at(uint i, uint j) { return this->operator()(i,j); }
inline Scalar at(uint i) const { return this->operator[](i); }
inline Scalar& at(uint i) { return this->operator[](i); }
inline RealScalar squaredLength() const { return squaredNorm(); }
inline RealScalar length() const { return norm(); }
inline RealScalar invLength(void) const { return fast_inv_sqrt(squaredNorm()); }
template<typename OtherDerived>
inline Scalar squaredDistanceTo(const MatrixBase<OtherDerived>& other) const
{ return (derived() - other.derived()).squaredNorm(); }
template<typename OtherDerived>
inline RealScalar distanceTo(const MatrixBase<OtherDerived>& other) const
{ return ei_sqrt(derived().squaredDistanceTo(other)); }
inline void scaleTo(RealScalar l) { RealScalar vl = norm(); if (vl>1e-9) derived() *= (l/vl); }
inline Transpose<Derived> transposed() {return this->transpose();}
inline const Transpose<Derived> transposed() const {return this->transpose();}
inline uint minComponentId(void) const { int i; this->minCoeff(&i); return i; }
inline uint maxComponentId(void) const { int i; this->maxCoeff(&i); return i; }
template<typename OtherDerived>
void makeFloor(const MatrixBase<OtherDerived>& other) { derived() = derived().cwiseMin(other.derived()); }
template<typename OtherDerived>
void makeCeil(const MatrixBase<OtherDerived>& other) { derived() = derived().cwiseMax(other.derived()); }
const CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>
operator+(const Scalar& scalar) const
{ return CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>(derived(), ei_scalar_add_op<Scalar>(scalar)); }
friend const CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>
operator+(const Scalar& scalar, const MatrixBase<Derived>& mat)
{ return CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>(mat.derived(), ei_scalar_add_op<Scalar>(scalar)); }
\endcode
Then one can the following declaration in the config.h or whatever prerequisites header file of his project:
\code
#define EIGEN_MATRIXBASE_PLUGIN "MatrixBaseAddons.h"
\endcode
\section InheritingFromMatrix Inheriting from Matrix
Before inheriting from Matrix, be really, i mean REALLY sure that using
EIGEN_MATRIX_PLUGIN is not what you really want (see previous section).
If you just need to add few members to Matrix, this is the way to go.
An example of when you actually need to inherit Matrix, is when you have
several layers of heritage such as MyVerySpecificVector1,MyVerySpecificVector1 -> MyVector1 -> Matrix and.
MyVerySpecificVector3,MyVerySpecificVector4 -> MyVector2 -> Matrix.
In order for your object to work within the Eigen framework, you need to
define a few members in your inherited class.
Here is a minimalistic example:\n
\code
class MyVectorType : public Eigen::VectorXd
{
public:
MyVectorType(void):Eigen::VectorXd() {}
// You need to define this for your object to work
typedef Eigen::VectorXd Base;
template<typename OtherDerived>
MyVectorType & operator= (const Eigen::MatrixBase <OtherDerived>& other)
{
this->Base::operator=(other);
return *this;
}
};
\endcode
This is the kind of error you can get if you don't provide those methods
\code
error: no match for ‘operator=’ in ‘delta =
(((Eigen::MatrixBase<Eigen::Matrix<std::complex<float>, 10000, 1, 2, 10000,
1> >*)(& delta)) + 8u)->Eigen::MatrixBase<Derived>::cwise [with Derived =
Eigen::Matrix<std::complex<float>, 10000, 1, 2, 10000,
1>]().Eigen::Cwise<ExpressionType>::operator* [with OtherDerived =
Eigen::Matrix<std::complex<float>, 10000, 1, 2, 10000, 1>, ExpressionType =
Eigen::Matrix<std::complex<float>, 10000, 1, 2, 10000, 1>](((const
Eigen::MatrixBase<Eigen::Matrix<std::complex<float>, 10000, 1, 2, 10000, 1>
>&)(((const Eigen::MatrixBase<Eigen::Matrix<std::complex<float>, 10000, 1,
>2, 10000, 1> >*)((const spectral1d*)where)) + 8u)))’
\endcode
\anchor user_defined_scalars \section CustomScalarType Using custom scalar types
By default, Eigen currently supports the following scalar types: \c int, \c float, \c double, \c std::complex<float>, \c std::complex<double>, \c long \c double, \c long \c long \c int (64 bits integers), and \c bool. The \c long \c double is especially useful on x86-64 systems or when the SSE2 instruction set is enabled because it enforces the use of x87 registers with extended accuracy.
In order to add support for a custom type \c T you need:
1 - make sure the common operator (+,-,*,/,etc.) are supported by the type \c T
2 - add a specialization of struct Eigen::NumTraits<T> (see \ref NumTraits)
3 - define a couple of math functions for your type such as: ei_sqrt, ei_abs, etc...
(see the file Eigen/src/Core/MathFunctions.h)
Here is a concrete example adding support for the Adolc's \c adouble type. <a href="https://projects.coin-or.org/ADOL-C">Adolc</a> is an automatic differentiation library. The type \c adouble is basically a real value tracking the values of any number of partial derivatives.
\code
#ifndef ADLOCSUPPORT_H
#define ADLOCSUPPORT_H
#define ADOLC_TAPELESS
#include <adolc/adouble.h>
#include <Eigen/Core>
namespace Eigen {
template<> struct NumTraits<adtl::adouble>
{
typedef adtl::adouble Real;
typedef adtl::adouble NonInteger;
typedef adtl::adouble Nested;
enum {
IsComplex = 0,
IsInteger = 0,
IsSigned,
ReadCost = 1,
AddCost = 1,
MulCost = 1
};
};
}
// the Adolc's type adouble is defined in the adtl namespace
// therefore, the following ei_* functions *must* be defined
// in the same namespace
namespace adtl {
inline const adouble& ei_conj(const adouble& x) { return x; }
inline const adouble& ei_real(const adouble& x) { return x; }
inline adouble ei_imag(const adouble&) { return 0.; }
inline adouble ei_abs(const adouble& x) { return fabs(x); }
inline adouble ei_abs2(const adouble& x) { return x*x; }
inline adouble ei_sqrt(const adouble& x) { return sqrt(x); }
inline adouble ei_exp(const adouble& x) { return exp(x); }
inline adouble ei_log(const adouble& x) { return log(x); }
inline adouble ei_sin(const adouble& x) { return sin(x); }
inline adouble ei_cos(const adouble& x) { return cos(x); }
inline adouble ei_pow(const adouble& x, adouble y) { return pow(x, y); }
}
#endif // ADLOCSUPPORT_H
\endcode
\section PreprocessorDirectives Preprocessor directives
You can control some aspects of Eigen by defining the following preprocessor tokens them before including any of Eigen's headers.
- \b #EIGEN_NO_DEBUG disables Eigen assertions. Like NDEBUG but only affects Eigen's assertions.
- \b EIGEN_DONT_VECTORIZE disables explicit vectorization when defined.
- \b EIGEN_UNROLLING_LIMIT defines the maximal instruction counts to enable meta unrolling of loops. Set it to zero to disable unrolling. The default is 100.
- \b EIGEN_DEFAULT_TO_ROW_MAJOR the default storage order for matrices becomes row-major instead of column-major.
- \b EIGEN_TUNE_FOR_CPU_CACHE_SIZE represents the maximal size in Bytes of L2 blocks. Since several blocks have to stay concurently in L2 cache, this value should correspond to at most 1/4 of the size of L2 cache.
- \b EIGEN_NO_STATIC_ASSERT replaces compile time static assertions by runtime assertions
- \b EIGEN_MATRIXBASE_PLUGIN see \ref ExtendingMatrixBase
*/
}
|
