• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1namespace Eigen {
2
3/** \page QuickRefPage Quick reference guide
4
5\b Table \b of \b contents
6  - \ref QuickRef_Headers
7  - \ref QuickRef_Types
8  - \ref QuickRef_Map
9  - \ref QuickRef_ArithmeticOperators
10  - \ref QuickRef_Coeffwise
11  - \ref QuickRef_Reductions
12  - \ref QuickRef_Blocks
13  - \ref QuickRef_Misc
14  - \ref QuickRef_DiagTriSymm
15\n
16
17<hr>
18
19<a href="#" class="top">top</a>
20\section QuickRef_Headers Modules and Header files
21
22The Eigen library is divided in a Core module and several additional modules. Each module has a corresponding header file which has to be included in order to use the module. The \c %Dense and \c Eigen header files are provided to conveniently gain access to several modules at once.
23
24<table class="manual">
25<tr><th>Module</th><th>Header file</th><th>Contents</th></tr>
26<tr><td>\link Core_Module Core \endlink</td><td>\code#include <Eigen/Core>\endcode</td><td>Matrix and Array classes, basic linear algebra (including triangular and selfadjoint products), array manipulation</td></tr>
27<tr class="alt"><td>\link Geometry_Module Geometry \endlink</td><td>\code#include <Eigen/Geometry>\endcode</td><td>Transform, Translation, Scaling, Rotation2D and 3D rotations (Quaternion, AngleAxis)</td></tr>
28<tr><td>\link LU_Module LU \endlink</td><td>\code#include <Eigen/LU>\endcode</td><td>Inverse, determinant, LU decompositions with solver (FullPivLU, PartialPivLU)</td></tr>
29<tr><td>\link Cholesky_Module Cholesky \endlink</td><td>\code#include <Eigen/Cholesky>\endcode</td><td>LLT and LDLT Cholesky factorization with solver</td></tr>
30<tr class="alt"><td>\link Householder_Module Householder \endlink</td><td>\code#include <Eigen/Householder>\endcode</td><td>Householder transformations; this module is used by several linear algebra modules</td></tr>
31<tr><td>\link SVD_Module SVD \endlink</td><td>\code#include <Eigen/SVD>\endcode</td><td>SVD decomposition with least-squares solver (JacobiSVD)</td></tr>
32<tr class="alt"><td>\link QR_Module QR \endlink</td><td>\code#include <Eigen/QR>\endcode</td><td>QR decomposition with solver (HouseholderQR, ColPivHouseholderQR, FullPivHouseholderQR)</td></tr>
33<tr><td>\link Eigenvalues_Module Eigenvalues \endlink</td><td>\code#include <Eigen/Eigenvalues>\endcode</td><td>Eigenvalue, eigenvector decompositions (EigenSolver, SelfAdjointEigenSolver, ComplexEigenSolver)</td></tr>
34<tr class="alt"><td>\link Sparse_Module Sparse \endlink</td><td>\code#include <Eigen/Sparse>\endcode</td><td>%Sparse matrix storage and related basic linear algebra (SparseMatrix, DynamicSparseMatrix, SparseVector)</td></tr>
35<tr><td></td><td>\code#include <Eigen/Dense>\endcode</td><td>Includes Core, Geometry, LU, Cholesky, SVD, QR, and Eigenvalues header files</td></tr>
36<tr class="alt"><td></td><td>\code#include <Eigen/Eigen>\endcode</td><td>Includes %Dense and %Sparse header files (the whole Eigen library)</td></tr>
37</table>
38
39<a href="#" class="top">top</a>
40\section QuickRef_Types Array, matrix and vector types
41
42
43\b Recall: Eigen provides two kinds of dense objects: mathematical matrices and vectors which are both represented by the template class Matrix, and general 1D and 2D arrays represented by the template class Array:
44\code
45typedef Matrix<Scalar, RowsAtCompileTime, ColsAtCompileTime, Options> MyMatrixType;
46typedef Array<Scalar, RowsAtCompileTime, ColsAtCompileTime, Options> MyArrayType;
47\endcode
48
49\li \c Scalar is the scalar type of the coefficients (e.g., \c float, \c double, \c bool, \c int, etc.).
50\li \c RowsAtCompileTime and \c ColsAtCompileTime are the number of rows and columns of the matrix as known at compile-time or \c Dynamic.
51\li \c Options can be \c ColMajor or \c RowMajor, default is \c ColMajor. (see class Matrix for more options)
52
53All combinations are allowed: you can have a matrix with a fixed number of rows and a dynamic number of columns, etc. The following are all valid:
54\code
55Matrix<double, 6, Dynamic>                  // Dynamic number of columns (heap allocation)
56Matrix<double, Dynamic, 2>                  // Dynamic number of rows (heap allocation)
57Matrix<double, Dynamic, Dynamic, RowMajor>  // Fully dynamic, row major (heap allocation)
58Matrix<double, 13, 3>                       // Fully fixed (static allocation)
59\endcode
60
61In most cases, you can simply use one of the convenience typedefs for \ref matrixtypedefs "matrices" and \ref arraytypedefs "arrays". Some examples:
62<table class="example">
63<tr><th>Matrices</th><th>Arrays</th></tr>
64<tr><td>\code
65Matrix<float,Dynamic,Dynamic>   <=>   MatrixXf
66Matrix<double,Dynamic,1>        <=>   VectorXd
67Matrix<int,1,Dynamic>           <=>   RowVectorXi
68Matrix<float,3,3>               <=>   Matrix3f
69Matrix<float,4,1>               <=>   Vector4f
70\endcode</td><td>\code
71Array<float,Dynamic,Dynamic>    <=>   ArrayXXf
72Array<double,Dynamic,1>         <=>   ArrayXd
73Array<int,1,Dynamic>            <=>   RowArrayXi
74Array<float,3,3>                <=>   Array33f
75Array<float,4,1>                <=>   Array4f
76\endcode</td></tr>
77</table>
78
79Conversion between the matrix and array worlds:
80\code
81Array44f a1, a1;
82Matrix4f m1, m2;
83m1 = a1 * a2;                     // coeffwise product, implicit conversion from array to matrix.
84a1 = m1 * m2;                     // matrix product, implicit conversion from matrix to array.
85a2 = a1 + m1.array();             // mixing array and matrix is forbidden
86m2 = a1.matrix() + m1;            // and explicit conversion is required.
87ArrayWrapper<Matrix4f> m1a(m1);   // m1a is an alias for m1.array(), they share the same coefficients
88MatrixWrapper<Array44f> a1m(a1);
89\endcode
90
91In the rest of this document we will use the following symbols to emphasize the features which are specifics to a given kind of object:
92\li <a name="matrixonly"><a/>\matrixworld linear algebra matrix and vector only
93\li <a name="arrayonly"><a/>\arrayworld array objects only
94
95\subsection QuickRef_Basics Basic matrix manipulation
96
97<table class="manual">
98<tr><th></th><th>1D objects</th><th>2D objects</th><th>Notes</th></tr>
99<tr><td>Constructors</td>
100<td>\code
101Vector4d  v4;
102Vector2f  v1(x, y);
103Array3i   v2(x, y, z);
104Vector4d  v3(x, y, z, w);
105
106VectorXf  v5; // empty object
107ArrayXf   v6(size);
108\endcode</td><td>\code
109Matrix4f  m1;
110
111
112
113
114MatrixXf  m5; // empty object
115MatrixXf  m6(nb_rows, nb_columns);
116\endcode</td><td class="note">
117By default, the coefficients \n are left uninitialized</td></tr>
118<tr class="alt"><td>Comma initializer</td>
119<td>\code
120Vector3f  v1;     v1 << x, y, z;
121ArrayXf   v2(4);  v2 << 1, 2, 3, 4;
122
123\endcode</td><td>\code
124Matrix3f  m1;   m1 << 1, 2, 3,
125                      4, 5, 6,
126                      7, 8, 9;
127\endcode</td><td></td></tr>
128
129<tr><td>Comma initializer (bis)</td>
130<td colspan="2">
131\include Tutorial_commainit_02.cpp
132</td>
133<td>
134output:
135\verbinclude Tutorial_commainit_02.out
136</td>
137</tr>
138
139<tr class="alt"><td>Runtime info</td>
140<td>\code
141vector.size();
142
143vector.innerStride();
144vector.data();
145\endcode</td><td>\code
146matrix.rows();          matrix.cols();
147matrix.innerSize();     matrix.outerSize();
148matrix.innerStride();   matrix.outerStride();
149matrix.data();
150\endcode</td><td class="note">Inner/Outer* are storage order dependent</td></tr>
151<tr><td>Compile-time info</td>
152<td colspan="2">\code
153ObjectType::Scalar              ObjectType::RowsAtCompileTime
154ObjectType::RealScalar          ObjectType::ColsAtCompileTime
155ObjectType::Index               ObjectType::SizeAtCompileTime
156\endcode</td><td></td></tr>
157<tr class="alt"><td>Resizing</td>
158<td>\code
159vector.resize(size);
160
161
162vector.resizeLike(other_vector);
163vector.conservativeResize(size);
164\endcode</td><td>\code
165matrix.resize(nb_rows, nb_cols);
166matrix.resize(Eigen::NoChange, nb_cols);
167matrix.resize(nb_rows, Eigen::NoChange);
168matrix.resizeLike(other_matrix);
169matrix.conservativeResize(nb_rows, nb_cols);
170\endcode</td><td class="note">no-op if the new sizes match,<br/>otherwise data are lost<br/><br/>resizing with data preservation</td></tr>
171
172<tr><td>Coeff access with \n range checking</td>
173<td>\code
174vector(i)     vector.x()
175vector[i]     vector.y()
176              vector.z()
177              vector.w()
178\endcode</td><td>\code
179matrix(i,j)
180\endcode</td><td class="note">Range checking is disabled if \n NDEBUG or EIGEN_NO_DEBUG is defined</td></tr>
181
182<tr class="alt"><td>Coeff access without \n range checking</td>
183<td>\code
184vector.coeff(i)
185vector.coeffRef(i)
186\endcode</td><td>\code
187matrix.coeff(i,j)
188matrix.coeffRef(i,j)
189\endcode</td><td></td></tr>
190
191<tr><td>Assignment/copy</td>
192<td colspan="2">\code
193object = expression;
194object_of_float = expression_of_double.cast<float>();
195\endcode</td><td class="note">the destination is automatically resized (if possible)</td></tr>
196
197</table>
198
199\subsection QuickRef_PredefMat Predefined Matrices
200
201<table class="manual">
202<tr>
203  <th>Fixed-size matrix or vector</th>
204  <th>Dynamic-size matrix</th>
205  <th>Dynamic-size vector</th>
206</tr>
207<tr style="border-bottom-style: none;">
208  <td>
209\code
210typedef {Matrix3f|Array33f} FixedXD;
211FixedXD x;
212
213x = FixedXD::Zero();
214x = FixedXD::Ones();
215x = FixedXD::Constant(value);
216x = FixedXD::Random();
217x = FixedXD::LinSpaced(size, low, high);
218
219x.setZero();
220x.setOnes();
221x.setConstant(value);
222x.setRandom();
223x.setLinSpaced(size, low, high);
224\endcode
225  </td>
226  <td>
227\code
228typedef {MatrixXf|ArrayXXf} Dynamic2D;
229Dynamic2D x;
230
231x = Dynamic2D::Zero(rows, cols);
232x = Dynamic2D::Ones(rows, cols);
233x = Dynamic2D::Constant(rows, cols, value);
234x = Dynamic2D::Random(rows, cols);
235N/A
236
237x.setZero(rows, cols);
238x.setOnes(rows, cols);
239x.setConstant(rows, cols, value);
240x.setRandom(rows, cols);
241N/A
242\endcode
243  </td>
244  <td>
245\code
246typedef {VectorXf|ArrayXf} Dynamic1D;
247Dynamic1D x;
248
249x = Dynamic1D::Zero(size);
250x = Dynamic1D::Ones(size);
251x = Dynamic1D::Constant(size, value);
252x = Dynamic1D::Random(size);
253x = Dynamic1D::LinSpaced(size, low, high);
254
255x.setZero(size);
256x.setOnes(size);
257x.setConstant(size, value);
258x.setRandom(size);
259x.setLinSpaced(size, low, high);
260\endcode
261  </td>
262</tr>
263
264<tr><td colspan="3">Identity and \link MatrixBase::Unit basis vectors \endlink \matrixworld</td></tr>
265<tr style="border-bottom-style: none;">
266  <td>
267\code
268x = FixedXD::Identity();
269x.setIdentity();
270
271Vector3f::UnitX() // 1 0 0
272Vector3f::UnitY() // 0 1 0
273Vector3f::UnitZ() // 0 0 1
274\endcode
275  </td>
276  <td>
277\code
278x = Dynamic2D::Identity(rows, cols);
279x.setIdentity(rows, cols);
280
281
282
283N/A
284\endcode
285  </td>
286  <td>\code
287N/A
288
289
290VectorXf::Unit(size,i)
291VectorXf::Unit(4,1) == Vector4f(0,1,0,0)
292                    == Vector4f::UnitY()
293\endcode
294  </td>
295</tr>
296</table>
297
298
299
300\subsection QuickRef_Map Mapping external arrays
301
302<table class="manual">
303<tr>
304<td>Contiguous \n memory</td>
305<td>\code
306float data[] = {1,2,3,4};
307Map<Vector3f> v1(data);       // uses v1 as a Vector3f object
308Map<ArrayXf>  v2(data,3);     // uses v2 as a ArrayXf object
309Map<Array22f> m1(data);       // uses m1 as a Array22f object
310Map<MatrixXf> m2(data,2,2);   // uses m2 as a MatrixXf object
311\endcode</td>
312</tr>
313<tr>
314<td>Typical usage \n of strides</td>
315<td>\code
316float data[] = {1,2,3,4,5,6,7,8,9};
317Map<VectorXf,0,InnerStride<2> >  v1(data,3);                      // = [1,3,5]
318Map<VectorXf,0,InnerStride<> >   v2(data,3,InnerStride<>(3));     // = [1,4,7]
319Map<MatrixXf,0,OuterStride<3> >  m2(data,2,3);                    // both lines     |1,4,7|
320Map<MatrixXf,0,OuterStride<> >   m1(data,2,3,OuterStride<>(3));   // are equal to:  |2,5,8|
321\endcode</td>
322</tr>
323</table>
324
325
326<a href="#" class="top">top</a>
327\section QuickRef_ArithmeticOperators Arithmetic Operators
328
329<table class="manual">
330<tr><td>
331add \n subtract</td><td>\code
332mat3 = mat1 + mat2;           mat3 += mat1;
333mat3 = mat1 - mat2;           mat3 -= mat1;\endcode
334</td></tr>
335<tr class="alt"><td>
336scalar product</td><td>\code
337mat3 = mat1 * s1;             mat3 *= s1;           mat3 = s1 * mat1;
338mat3 = mat1 / s1;             mat3 /= s1;\endcode
339</td></tr>
340<tr><td>
341matrix/vector \n products \matrixworld</td><td>\code
342col2 = mat1 * col1;
343row2 = row1 * mat1;           row1 *= mat1;
344mat3 = mat1 * mat2;           mat3 *= mat1; \endcode
345</td></tr>
346<tr class="alt"><td>
347transposition \n adjoint \matrixworld</td><td>\code
348mat1 = mat2.transpose();      mat1.transposeInPlace();
349mat1 = mat2.adjoint();        mat1.adjointInPlace();
350\endcode
351</td></tr>
352<tr><td>
353\link MatrixBase::dot() dot \endlink product \n inner product \matrixworld</td><td>\code
354scalar = vec1.dot(vec2);
355scalar = col1.adjoint() * col2;
356scalar = (col1.adjoint() * col2).value();\endcode
357</td></tr>
358<tr class="alt"><td>
359outer product \matrixworld</td><td>\code
360mat = col1 * col2.transpose();\endcode
361</td></tr>
362
363<tr><td>
364\link MatrixBase::norm() norm \endlink \n \link MatrixBase::normalized() normalization \endlink \matrixworld</td><td>\code
365scalar = vec1.norm();         scalar = vec1.squaredNorm()
366vec2 = vec1.normalized();     vec1.normalize(); // inplace \endcode
367</td></tr>
368
369<tr class="alt"><td>
370\link MatrixBase::cross() cross product \endlink \matrixworld</td><td>\code
371#include <Eigen/Geometry>
372vec3 = vec1.cross(vec2);\endcode</td></tr>
373</table>
374
375<a href="#" class="top">top</a>
376\section QuickRef_Coeffwise Coefficient-wise \& Array operators
377Coefficient-wise operators for matrices and vectors:
378<table class="manual">
379<tr><th>Matrix API \matrixworld</th><th>Via Array conversions</th></tr>
380<tr><td>\code
381mat1.cwiseMin(mat2)
382mat1.cwiseMax(mat2)
383mat1.cwiseAbs2()
384mat1.cwiseAbs()
385mat1.cwiseSqrt()
386mat1.cwiseProduct(mat2)
387mat1.cwiseQuotient(mat2)\endcode
388</td><td>\code
389mat1.array().min(mat2.array())
390mat1.array().max(mat2.array())
391mat1.array().abs2()
392mat1.array().abs()
393mat1.array().sqrt()
394mat1.array() * mat2.array()
395mat1.array() / mat2.array()
396\endcode</td></tr>
397</table>
398
399It is also very simple to apply any user defined function \c foo using DenseBase::unaryExpr together with std::ptr_fun:
400\code mat1.unaryExpr(std::ptr_fun(foo))\endcode
401
402Array operators:\arrayworld
403
404<table class="manual">
405<tr><td>Arithmetic operators</td><td>\code
406array1 * array2     array1 / array2     array1 *= array2    array1 /= array2
407array1 + scalar     array1 - scalar     array1 += scalar    array1 -= scalar
408\endcode</td></tr>
409<tr><td>Comparisons</td><td>\code
410array1 < array2     array1 > array2     array1 < scalar     array1 > scalar
411array1 <= array2    array1 >= array2    array1 <= scalar    array1 >= scalar
412array1 == array2    array1 != array2    array1 == scalar    array1 != scalar
413\endcode</td></tr>
414<tr><td>Trigo, power, and \n misc functions \n and the STL variants</td><td>\code
415array1.min(array2)
416array1.max(array2)
417array1.abs2()
418array1.abs()                  std::abs(array1)
419array1.sqrt()                 std::sqrt(array1)
420array1.log()                  std::log(array1)
421array1.exp()                  std::exp(array1)
422array1.pow(exponent)          std::pow(array1,exponent)
423array1.square()
424array1.cube()
425array1.inverse()
426array1.sin()                  std::sin(array1)
427array1.cos()                  std::cos(array1)
428array1.tan()                  std::tan(array1)
429array1.asin()                 std::asin(array1)
430array1.acos()                 std::acos(array1)
431\endcode
432</td></tr>
433</table>
434
435<a href="#" class="top">top</a>
436\section QuickRef_Reductions Reductions
437
438Eigen provides several reduction methods such as:
439\link DenseBase::minCoeff() minCoeff() \endlink, \link DenseBase::maxCoeff() maxCoeff() \endlink,
440\link DenseBase::sum() sum() \endlink, \link DenseBase::prod() prod() \endlink,
441\link MatrixBase::trace() trace() \endlink \matrixworld,
442\link MatrixBase::norm() norm() \endlink \matrixworld, \link MatrixBase::squaredNorm() squaredNorm() \endlink \matrixworld,
443\link DenseBase::all() all() \endlink, and \link DenseBase::any() any() \endlink.
444All reduction operations can be done matrix-wise,
445\link DenseBase::colwise() column-wise \endlink or
446\link DenseBase::rowwise() row-wise \endlink. Usage example:
447<table class="manual">
448<tr><td rowspan="3" style="border-right-style:dashed;vertical-align:middle">\code
449      5 3 1
450mat = 2 7 8
451      9 4 6 \endcode
452</td> <td>\code mat.minCoeff(); \endcode</td><td>\code 1 \endcode</td></tr>
453<tr class="alt"><td>\code mat.colwise().minCoeff(); \endcode</td><td>\code 2 3 1 \endcode</td></tr>
454<tr style="vertical-align:middle"><td>\code mat.rowwise().minCoeff(); \endcode</td><td>\code
4551
4562
4574
458\endcode</td></tr>
459</table>
460
461Special versions of \link DenseBase::minCoeff(Index*,Index*) minCoeff \endlink and \link DenseBase::maxCoeff(Index*,Index*) maxCoeff \endlink:
462\code
463int i, j;
464s = vector.minCoeff(&i);        // s == vector[i]
465s = matrix.maxCoeff(&i, &j);    // s == matrix(i,j)
466\endcode
467Typical use cases of all() and any():
468\code
469if((array1 > 0).all()) ...      // if all coefficients of array1 are greater than 0 ...
470if((array1 < array2).any()) ... // if there exist a pair i,j such that array1(i,j) < array2(i,j) ...
471\endcode
472
473
474<a href="#" class="top">top</a>\section QuickRef_Blocks Sub-matrices
475
476Read-write access to a \link DenseBase::col(Index) column \endlink
477or a \link DenseBase::row(Index) row \endlink of a matrix (or array):
478\code
479mat1.row(i) = mat2.col(j);
480mat1.col(j1).swap(mat1.col(j2));
481\endcode
482
483Read-write access to sub-vectors:
484<table class="manual">
485<tr>
486<th>Default versions</th>
487<th>Optimized versions when the size \n is known at compile time</th></tr>
488<th></th>
489
490<tr><td>\code vec1.head(n)\endcode</td><td>\code vec1.head<n>()\endcode</td><td>the first \c n coeffs </td></tr>
491<tr><td>\code vec1.tail(n)\endcode</td><td>\code vec1.tail<n>()\endcode</td><td>the last \c n coeffs </td></tr>
492<tr><td>\code vec1.segment(pos,n)\endcode</td><td>\code vec1.segment<n>(pos)\endcode</td>
493    <td>the \c n coeffs in \n the range [\c pos : \c pos + \c n [</td></tr>
494<tr class="alt"><td colspan="3">
495
496Read-write access to sub-matrices:</td></tr>
497<tr>
498  <td>\code mat1.block(i,j,rows,cols)\endcode
499      \link DenseBase::block(Index,Index,Index,Index) (more) \endlink</td>
500  <td>\code mat1.block<rows,cols>(i,j)\endcode
501      \link DenseBase::block(Index,Index) (more) \endlink</td>
502  <td>the \c rows x \c cols sub-matrix \n starting from position (\c i,\c j)</td></tr>
503<tr><td>\code
504 mat1.topLeftCorner(rows,cols)
505 mat1.topRightCorner(rows,cols)
506 mat1.bottomLeftCorner(rows,cols)
507 mat1.bottomRightCorner(rows,cols)\endcode
508 <td>\code
509 mat1.topLeftCorner<rows,cols>()
510 mat1.topRightCorner<rows,cols>()
511 mat1.bottomLeftCorner<rows,cols>()
512 mat1.bottomRightCorner<rows,cols>()\endcode
513 <td>the \c rows x \c cols sub-matrix \n taken in one of the four corners</td></tr>
514 <tr><td>\code
515 mat1.topRows(rows)
516 mat1.bottomRows(rows)
517 mat1.leftCols(cols)
518 mat1.rightCols(cols)\endcode
519 <td>\code
520 mat1.topRows<rows>()
521 mat1.bottomRows<rows>()
522 mat1.leftCols<cols>()
523 mat1.rightCols<cols>()\endcode
524 <td>specialized versions of block() \n when the block fit two corners</td></tr>
525</table>
526
527
528
529<a href="#" class="top">top</a>\section QuickRef_Misc Miscellaneous operations
530
531\subsection QuickRef_Reverse Reverse
532Vectors, rows, and/or columns of a matrix can be reversed (see DenseBase::reverse(), DenseBase::reverseInPlace(), VectorwiseOp::reverse()).
533\code
534vec.reverse()           mat.colwise().reverse()   mat.rowwise().reverse()
535vec.reverseInPlace()
536\endcode
537
538\subsection QuickRef_Replicate Replicate
539Vectors, matrices, rows, and/or columns can be replicated in any direction (see DenseBase::replicate(), VectorwiseOp::replicate())
540\code
541vec.replicate(times)                                          vec.replicate<Times>
542mat.replicate(vertical_times, horizontal_times)               mat.replicate<VerticalTimes, HorizontalTimes>()
543mat.colwise().replicate(vertical_times, horizontal_times)     mat.colwise().replicate<VerticalTimes, HorizontalTimes>()
544mat.rowwise().replicate(vertical_times, horizontal_times)     mat.rowwise().replicate<VerticalTimes, HorizontalTimes>()
545\endcode
546
547
548<a href="#" class="top">top</a>\section QuickRef_DiagTriSymm Diagonal, Triangular, and Self-adjoint matrices
549(matrix world \matrixworld)
550
551\subsection QuickRef_Diagonal Diagonal matrices
552
553<table class="example">
554<tr><th>Operation</th><th>Code</th></tr>
555<tr><td>
556view a vector \link MatrixBase::asDiagonal() as a diagonal matrix \endlink \n </td><td>\code
557mat1 = vec1.asDiagonal();\endcode
558</td></tr>
559<tr><td>
560Declare a diagonal matrix</td><td>\code
561DiagonalMatrix<Scalar,SizeAtCompileTime> diag1(size);
562diag1.diagonal() = vector;\endcode
563</td></tr>
564<tr><td>Access the \link MatrixBase::diagonal() diagonal \endlink and \link MatrixBase::diagonal(Index) super/sub diagonals \endlink of a matrix as a vector (read/write)</td>
565 <td>\code
566vec1 = mat1.diagonal();        mat1.diagonal() = vec1;      // main diagonal
567vec1 = mat1.diagonal(+n);      mat1.diagonal(+n) = vec1;    // n-th super diagonal
568vec1 = mat1.diagonal(-n);      mat1.diagonal(-n) = vec1;    // n-th sub diagonal
569vec1 = mat1.diagonal<1>();     mat1.diagonal<1>() = vec1;   // first super diagonal
570vec1 = mat1.diagonal<-2>();    mat1.diagonal<-2>() = vec1;  // second sub diagonal
571\endcode</td>
572</tr>
573
574<tr><td>Optimized products and inverse</td>
575 <td>\code
576mat3  = scalar * diag1 * mat1;
577mat3 += scalar * mat1 * vec1.asDiagonal();
578mat3 = vec1.asDiagonal().inverse() * mat1
579mat3 = mat1 * diag1.inverse()
580\endcode</td>
581</tr>
582
583</table>
584
585\subsection QuickRef_TriangularView Triangular views
586
587TriangularView gives a view on a triangular part of a dense matrix and allows to perform optimized operations on it. The opposite triangular part is never referenced and can be used to store other information.
588
589\note The .triangularView() template member function requires the \c template keyword if it is used on an
590object of a type that depends on a template parameter; see \ref TopicTemplateKeyword for details.
591
592<table class="example">
593<tr><th>Operation</th><th>Code</th></tr>
594<tr><td>
595Reference to a triangular with optional \n
596unit or null diagonal (read/write):
597</td><td>\code
598m.triangularView<Xxx>()
599\endcode \n
600\c Xxx = ::Upper, ::Lower, ::StrictlyUpper, ::StrictlyLower, ::UnitUpper, ::UnitLower
601</td></tr>
602<tr><td>
603Writing to a specific triangular part:\n (only the referenced triangular part is evaluated)
604</td><td>\code
605m1.triangularView<Eigen::Lower>() = m2 + m3 \endcode
606</td></tr>
607<tr><td>
608Conversion to a dense matrix setting the opposite triangular part to zero:
609</td><td>\code
610m2 = m1.triangularView<Eigen::UnitUpper>()\endcode
611</td></tr>
612<tr><td>
613Products:
614</td><td>\code
615m3 += s1 * m1.adjoint().triangularView<Eigen::UnitUpper>() * m2
616m3 -= s1 * m2.conjugate() * m1.adjoint().triangularView<Eigen::Lower>() \endcode
617</td></tr>
618<tr><td>
619Solving linear equations:\n
620\f$ M_2 := L_1^{-1} M_2 \f$ \n
621\f$ M_3 := {L_1^*}^{-1} M_3 \f$ \n
622\f$ M_4 := M_4 U_1^{-1} \f$
623</td><td>\n \code
624L1.triangularView<Eigen::UnitLower>().solveInPlace(M2)
625L1.triangularView<Eigen::Lower>().adjoint().solveInPlace(M3)
626U1.triangularView<Eigen::Upper>().solveInPlace<OnTheRight>(M4)\endcode
627</td></tr>
628</table>
629
630\subsection QuickRef_SelfadjointMatrix Symmetric/selfadjoint views
631
632Just as for triangular matrix, you can reference any triangular part of a square matrix to see it as a selfadjoint
633matrix and perform special and optimized operations. Again the opposite triangular part is never referenced and can be
634used to store other information.
635
636\note The .selfadjointView() template member function requires the \c template keyword if it is used on an
637object of a type that depends on a template parameter; see \ref TopicTemplateKeyword for details.
638
639<table class="example">
640<tr><th>Operation</th><th>Code</th></tr>
641<tr><td>
642Conversion to a dense matrix:
643</td><td>\code
644m2 = m.selfadjointView<Eigen::Lower>();\endcode
645</td></tr>
646<tr><td>
647Product with another general matrix or vector:
648</td><td>\code
649m3  = s1 * m1.conjugate().selfadjointView<Eigen::Upper>() * m3;
650m3 -= s1 * m3.adjoint() * m1.selfadjointView<Eigen::Lower>();\endcode
651</td></tr>
652<tr><td>
653Rank 1 and rank K update: \n
654\f$ upper(M_1) \mathrel{{+}{=}} s_1 M_2 M_2^* \f$ \n
655\f$ lower(M_1) \mathbin{{-}{=}} M_2^* M_2 \f$
656</td><td>\n \code
657M1.selfadjointView<Eigen::Upper>().rankUpdate(M2,s1);
658M1.selfadjointView<Eigen::Lower>().rankUpdate(M2.adjoint(),-1); \endcode
659</td></tr>
660<tr><td>
661Rank 2 update: (\f$ M \mathrel{{+}{=}} s u v^* + s v u^* \f$)
662</td><td>\code
663M.selfadjointView<Eigen::Upper>().rankUpdate(u,v,s);
664\endcode
665</td></tr>
666<tr><td>
667Solving linear equations:\n(\f$ M_2 := M_1^{-1} M_2 \f$)
668</td><td>\code
669// via a standard Cholesky factorization
670m2 = m1.selfadjointView<Eigen::Upper>().llt().solve(m2);
671// via a Cholesky factorization with pivoting
672m2 = m1.selfadjointView<Eigen::Lower>().ldlt().solve(m2);
673\endcode
674</td></tr>
675</table>
676
677*/
678
679/*
680<table class="tutorial_code">
681<tr><td>
682\link MatrixBase::asDiagonal() make a diagonal matrix \endlink \n from a vector </td><td>\code
683mat1 = vec1.asDiagonal();\endcode
684</td></tr>
685<tr><td>
686Declare a diagonal matrix</td><td>\code
687DiagonalMatrix<Scalar,SizeAtCompileTime> diag1(size);
688diag1.diagonal() = vector;\endcode
689</td></tr>
690<tr><td>Access \link MatrixBase::diagonal() the diagonal and super/sub diagonals of a matrix \endlink as a vector (read/write)</td>
691 <td>\code
692vec1 = mat1.diagonal();            mat1.diagonal() = vec1;      // main diagonal
693vec1 = mat1.diagonal(+n);          mat1.diagonal(+n) = vec1;    // n-th super diagonal
694vec1 = mat1.diagonal(-n);          mat1.diagonal(-n) = vec1;    // n-th sub diagonal
695vec1 = mat1.diagonal<1>();         mat1.diagonal<1>() = vec1;   // first super diagonal
696vec1 = mat1.diagonal<-2>();        mat1.diagonal<-2>() = vec1;  // second sub diagonal
697\endcode</td>
698</tr>
699
700<tr><td>View on a triangular part of a matrix (read/write)</td>
701 <td>\code
702mat2 = mat1.triangularView<Xxx>();
703// Xxx = Upper, Lower, StrictlyUpper, StrictlyLower, UnitUpper, UnitLower
704mat1.triangularView<Upper>() = mat2 + mat3; // only the upper part is evaluated and referenced
705\endcode</td></tr>
706
707<tr><td>View a triangular part as a symmetric/self-adjoint matrix (read/write)</td>
708 <td>\code
709mat2 = mat1.selfadjointView<Xxx>();     // Xxx = Upper or Lower
710mat1.selfadjointView<Upper>() = mat2 + mat2.adjoint();  // evaluated and write to the upper triangular part only
711\endcode</td></tr>
712
713</table>
714
715Optimized products:
716\code
717mat3 += scalar * vec1.asDiagonal() * mat1
718mat3 += scalar * mat1 * vec1.asDiagonal()
719mat3.noalias() += scalar * mat1.triangularView<Xxx>() * mat2
720mat3.noalias() += scalar * mat2 * mat1.triangularView<Xxx>()
721mat3.noalias() += scalar * mat1.selfadjointView<Upper or Lower>() * mat2
722mat3.noalias() += scalar * mat2 * mat1.selfadjointView<Upper or Lower>()
723mat1.selfadjointView<Upper or Lower>().rankUpdate(mat2);
724mat1.selfadjointView<Upper or Lower>().rankUpdate(mat2.adjoint(), scalar);
725\endcode
726
727Inverse products: (all are optimized)
728\code
729mat3 = vec1.asDiagonal().inverse() * mat1
730mat3 = mat1 * diag1.inverse()
731mat1.triangularView<Xxx>().solveInPlace(mat2)
732mat1.triangularView<Xxx>().solveInPlace<OnTheRight>(mat2)
733mat2 = mat1.selfadjointView<Upper or Lower>().llt().solve(mat2)
734\endcode
735
736*/
737}
738