Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 1 | namespace Eigen { |
| 2 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 3 | /** \page TutorialMatrixArithmetic Tutorial page 2 - %Matrix and vector arithmetic |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 4 | \ingroup Tutorial |
| 5 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 6 | \li \b Previous: \ref TutorialMatrixClass |
| 7 | \li \b Next: \ref TutorialArrayClass |
| 8 | |
| 9 | This tutorial aims to provide an overview and some details on how to perform arithmetic |
| 10 | between matrices, vectors and scalars with Eigen. |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 11 | |
| 12 | \b Table \b of \b contents |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 13 | - \ref TutorialArithmeticIntroduction |
| 14 | - \ref TutorialArithmeticMentionCommaInitializer |
| 15 | - \ref TutorialArithmeticAddSub |
| 16 | - \ref TutorialArithmeticScalarMulDiv |
| 17 | - \ref TutorialArithmeticMentionXprTemplates |
| 18 | - \ref TutorialArithmeticMatrixMul |
| 19 | - \ref TutorialArithmeticDotAndCross |
| 20 | - \ref TutorialArithmeticRedux |
| 21 | - \ref TutorialArithmeticValidity |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 22 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 23 | \section TutorialArithmeticIntroduction Introduction |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 24 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 25 | Eigen offers matrix/vector arithmetic operations either through overloads of common C++ arithmetic operators such as +, -, *, |
| 26 | or through special methods such as dot(), cross(), etc. |
| 27 | For the Matrix class (matrices and vectors), operators are only overloaded to support |
| 28 | linear-algebraic operations. For example, \c matrix1 \c * \c matrix2 means matrix-matrix product, |
| 29 | and \c vector \c + \c scalar is just not allowed. If you want to perform all kinds of array operations, |
| 30 | not linear algebra, see \ref TutorialArrayClass "next page". |
| 31 | |
| 32 | \section TutorialArithmeticMentionCommaInitializer A note about comma-initialization |
| 33 | |
| 34 | In examples below, we'll be initializing matrices with a very convenient device known as the \em comma-initializer. |
| 35 | For now, it is enough to know this example: |
| 36 | \include Tutorial_commainit_01.cpp |
| 37 | Output: \verbinclude Tutorial_commainit_01.out |
| 38 | The comma initializer is discussed in \ref TutorialAdvancedInitialization "this page". |
| 39 | |
| 40 | \section TutorialArithmeticAddSub Addition and substraction |
| 41 | |
| 42 | The left hand side and right hand side must, of course, have the same numbers of rows and of columns. They must |
| 43 | also have the same \c Scalar type, as Eigen doesn't do automatic type promotion. The operators at hand here are: |
| 44 | \li binary operator + as in \c a+b |
| 45 | \li binary operator - as in \c a-b |
| 46 | \li unary operator - as in \c -a |
| 47 | \li compound operator += as in \c a+=b |
| 48 | \li compound operator -= as in \c a-=b |
| 49 | |
| 50 | Example: \include tut_arithmetic_add_sub.cpp |
| 51 | Output: \verbinclude tut_arithmetic_add_sub.out |
| 52 | |
| 53 | \section TutorialArithmeticScalarMulDiv Scalar multiplication and division |
| 54 | |
| 55 | Multiplication and division by a scalar is very simple too. The operators at hand here are: |
| 56 | \li binary operator * as in \c matrix*scalar |
| 57 | \li binary operator * as in \c scalar*matrix |
| 58 | \li binary operator / as in \c matrix/scalar |
| 59 | \li compound operator *= as in \c matrix*=scalar |
| 60 | \li compound operator /= as in \c matrix/=scalar |
| 61 | |
| 62 | Example: \include tut_arithmetic_scalar_mul_div.cpp |
| 63 | Output: \verbinclude tut_arithmetic_scalar_mul_div.out |
| 64 | |
| 65 | \section TutorialArithmeticMentionXprTemplates A note about expression templates |
| 66 | |
| 67 | This is an advanced topic that we explain in \ref TopicEigenExpressionTemplates "this page", |
| 68 | but it is useful to just mention it now. In Eigen, arithmetic operators such as \c operator+ don't |
| 69 | perform any computation by themselves, they just return an "expression object" describing the computation to be |
| 70 | performed. The actual computation happens later, when the whole expression is evaluated, typically in \c operator=. |
| 71 | While this might sound heavy, any modern optimizing compiler is able to optimize away that abstraction and |
| 72 | the result is perfectly optimized code. For example, when you do: |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 73 | \code |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 74 | VectorXf a(50), b(50), c(50), d(50); |
| 75 | ... |
| 76 | a = 3*b + 4*c + 5*d; |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 77 | \endcode |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 78 | Eigen compiles it to just one for loop, so that the arrays are traversed only once. Simplyfying (e.g. ignoring |
| 79 | SIMD optimizations), this loop looks like this: |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 80 | \code |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 81 | for(int i = 0; i < 50; ++i) |
| 82 | a[i] = 3*b[i] + 4*c[i] + 5*d[i]; |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 83 | \endcode |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 84 | Thus, you should not be afraid of using relatively large arithmetic expressions with Eigen: it only gives Eigen |
| 85 | more opportunities for optimization. |
| 86 | |
| 87 | \section TutorialArithmeticMatrixMul Matrix-matrix and matrix-vector multiplication |
| 88 | |
| 89 | Matrix-matrix multiplication is again done with \c operator*. Since vectors are a special |
| 90 | case of matrices, they are implicitly handled there too, so matrix-vector product is really just a special |
| 91 | case of matrix-matrix product, and so is vector-vector outer product. Thus, all these cases are handled by just |
| 92 | two operators: |
| 93 | \li binary operator * as in \c a*b |
| 94 | \li compound operator *= as in \c a*=b |
| 95 | |
| 96 | Example: \include tut_arithmetic_matrix_mul.cpp |
| 97 | Output: \verbinclude tut_arithmetic_matrix_mul.out |
| 98 | |
| 99 | Note: if you read the above paragraph on expression templates and are worried that doing \c m=m*m might cause |
| 100 | aliasing issues, be reassured for now: Eigen treats matrix multiplication as a special case and takes care of |
| 101 | introducing a temporary here, so it will compile \c m=m*m as: |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 102 | \code |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 103 | tmp = m*m; |
| 104 | m = tmp; |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 105 | \endcode |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 106 | For more details on this topic, see \ref TopicEigenExpressionTemplates "this page". |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 107 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 108 | \section TutorialArithmeticDotAndCross Dot product and cross product |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 109 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 110 | The above-discussed \c operator* does not allow to compute dot and cross products. For that, you need the dot() and cross() methods. |
| 111 | Example: \include tut_arithmetic_dot_cross.cpp |
| 112 | Output: \verbinclude tut_arithmetic_dot_cross.out |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 113 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 114 | Remember that cross product is only for vectors of size 3. Dot product is for vectors of any sizes. |
| 115 | When using complex numbers, Eigen's dot product is conjugate-linear in the first variable and linear in the |
| 116 | second variable. |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 117 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 118 | \section TutorialArithmeticRedux Basic arithmetic reduction operations |
| 119 | Eigen also provides some reduction operations to obtain values such as the sum or the maximum |
| 120 | or minimum of all the coefficients in a given matrix or vector. |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 121 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 122 | TODO: convert this from table format to tutorial/examples format. |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 123 | |
| 124 | <table class="tutorial_code" align="center"> |
| 125 | <tr><td align="center">\b Reduction \b operation</td><td align="center">\b Usage \b example</td></tr> |
| 126 | <tr><td> |
| 127 | Sum of all the coefficients in a matrix</td><td>\code |
| 128 | MatrixXf m; |
| 129 | float totalSum = m.sum();\endcode</td></tr> |
| 130 | <tr><td> |
| 131 | Maximum coefficient in a matrix</td><td>\code |
| 132 | MatrixXf m; |
| 133 | int row, col; |
| 134 | |
| 135 | // minimum value will be stored in minValue |
| 136 | // and the row and column where it was found in row and col, |
| 137 | // (these two parameters are optional) |
| 138 | float minValue = m.minCoeff(&row,&col);\endcode</td></tr> |
| 139 | <tr><td> |
| 140 | Maximum coefficient in a matrix</td><td>\code |
| 141 | MatrixXf m; |
| 142 | int row, col; |
| 143 | |
| 144 | // maximum value will be stored in maxValue |
| 145 | // and the row and column where it was found in row and col, |
| 146 | // (these two parameters are optional) |
| 147 | float maxValue = m.maxCoeff(&row,&col);\endcode</td></tr> |
| 148 | <tr><td> |
| 149 | Product between all coefficients in a matrix</td><td>\code |
| 150 | MatrixXf m; |
| 151 | |
| 152 | float product = m.prod();\endcode</td></tr> |
| 153 | <tr><td> |
| 154 | Mean of coefficients in a matrix</td><td>\code |
| 155 | MatrixXf m; |
| 156 | |
| 157 | float mean = m.mean();\endcode</td></tr> |
| 158 | <tr><td> |
| 159 | Matrix's trace</td><td>\code |
| 160 | MatrixXf m; |
| 161 | |
| 162 | float trace = m.trace();\endcode</td></tr> |
| 163 | </table> |
| 164 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 165 | \subsection TutorialArithmeticValidity Validity of operations |
| 166 | Eigen checks the validity of the operations that you perform. When possible, |
| 167 | it checks them at compile-time, producing compilation errors. These error messages can be long and ugly, |
| 168 | but Eigen writes the important message in UPPERCASE_LETTERS_SO_IT_STANDS_OUT. For example: |
| 169 | \code |
| 170 | Matrix3f m; |
| 171 | Vector4f v; |
| 172 | v = m*v; // Compile-time error: YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES |
| 173 | \endcode |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 174 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 175 | Of course, in many cases, for example when checking dynamic sizes, the check cannot be performed at compile time. |
| 176 | Eigen then uses runtime assertions. This means that executing an illegal operation will result in a crash at runtime, |
| 177 | with an error message. |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 178 | |
Benoit Jacob | e078bb2 | 2010-06-26 14:00:00 -0400 | [diff] [blame] | 179 | \code |
| 180 | MatrixXf m(3,3); |
| 181 | VectorXf v(4); |
| 182 | v = m * v; // Run-time assertion failure here: "invalid matrix product" |
| 183 | \endcode |
| 184 | |
| 185 | For more details on this topic, see \ref TopicAssertions "this page". |
| 186 | |
| 187 | \li \b Next: \ref TutorialArrayClass |
Carlos Becker | 9d44005 | 2010-06-25 20:16:12 -0400 | [diff] [blame] | 188 | |
| 189 | */ |
| 190 | |
| 191 | } |