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