diff options
Diffstat (limited to 'doc')
113 files changed, 2923 insertions, 720 deletions
diff --git a/doc/A05_PortingFrom2To3.dox b/doc/A05_PortingFrom2To3.dox index d885b4f6d..51555f996 100644 --- a/doc/A05_PortingFrom2To3.dox +++ b/doc/A05_PortingFrom2To3.dox @@ -9,11 +9,8 @@ and gives tips to help porting your application from Eigen2 to Eigen3. \section CompatibilitySupport Eigen2 compatibility support -In order to ease the switch from Eigen2 to Eigen3, Eigen3 features \subpage Eigen2SupportModes "Eigen2 support modes". - -The quick way to enable this is to define the \c EIGEN2_SUPPORT preprocessor token \b before including any Eigen header (typically it should be set in your project options). - -A more powerful, \em staged migration path is also provided, which may be useful to migrate larger projects from Eigen2 to Eigen3. This is explained in the \ref Eigen2SupportModes "Eigen 2 support modes" page. +Up to version 3.2 %Eigen provides <a href="http://eigen.tuxfamily.org/dox-3.2/Eigen2SupportModes.html">Eigen2 support modes</a>. These are removed now, because they were barely used anymore and became hard to maintain after internal re-designs. +You can still use them by first <a href="http://eigen.tuxfamily.org/dox-3.2/Eigen2ToEigen3.html">porting your code to Eigen 3.2</a>. \section Using The USING_PART_OF_NAMESPACE_EIGEN macro @@ -226,7 +223,7 @@ triangular part to work on</td></tr> \section GeometryModule Changes in the Geometry module -The Geometry module is the one that changed the most. If you rely heavily on it, it's probably a good idea to use the \ref Eigen2SupportModes "Eigen 2 support modes" to perform your migration. +The Geometry module is the one that changed the most. If you rely heavily on it, it's probably a good idea to use the <a href="http://eigen.tuxfamily.org/dox-3.2/Eigen2SupportModes.html">"Eigen 2 support modes"</a> to perform your migration. \section Transform The Transform class @@ -264,7 +261,7 @@ use it unless you are sure of what you are doing, i.e., you have rigourosly meas The EIGEN_ALIGN_128 macro has been renamed to EIGEN_ALIGN16. Don't be surprised, it's just that we switched to counting in bytes ;-) -The EIGEN_DONT_ALIGN option still exists in Eigen 3, but it has a new cousin: EIGEN_DONT_ALIGN_STATICALLY. It allows to get rid of all static alignment issues while keeping alignment of dynamic-size heap-allocated arrays, thus keeping vectorization for dynamic-size objects. +The \link TopicPreprocessorDirectivesPerformance EIGEN_DONT_ALIGN \endlink option still exists in Eigen 3, but it has a new cousin: \link TopicPreprocessorDirectivesPerformance EIGEN_DONT_ALIGN_STATICALLY.\endlink It allows to get rid of all static alignment issues while keeping alignment of dynamic-size heap-allocated arrays. Vectorization of statically allocated arrays is still preserved (unless you define \link TopicPreprocessorDirectivesPerformance EIGEN_UNALIGNED_VECTORIZE \endlink =0), at the cost of unaligned memory stores. \section AlignedMap Aligned Map objects @@ -281,7 +278,7 @@ result = Vector4f::MapAligned(some_aligned_array); \section StdContainers STL Containers -In Eigen2, <tt>#include<Eigen/StdVector></tt> tweaked std::vector to automatically align elements. The problem was that that was quite invasive. In Eigen3, we only override standard behavior if you use Eigen::aligned_allocator<T> as your allocator type. So for example, if you use std::vector<Matrix4f>, you need to do the following change (note that aligned_allocator is under namespace Eigen): +In Eigen2, <tt>\#include\<Eigen/StdVector\></tt> tweaked std::vector to automatically align elements. The problem was that that was quite invasive. In Eigen3, we only override standard behavior if you use Eigen::aligned_allocator<T> as your allocator type. So for example, if you use std::vector<Matrix4f>, you need to do the following change (note that aligned_allocator is under namespace Eigen): <table class="manual"> <tr><th>Eigen 2</th><th>Eigen 3</th></tr> diff --git a/doc/A10_Eigen2SupportModes.dox b/doc/A10_Eigen2SupportModes.dox deleted file mode 100644 index abfdb4683..000000000 --- a/doc/A10_Eigen2SupportModes.dox +++ /dev/null @@ -1,93 +0,0 @@ -namespace Eigen { - -/** \page Eigen2SupportModes Eigen 2 support modes - -This page documents the Eigen2 support modes, a powerful tool to help migrating your project from Eigen 2 to Eigen 3. -Don't miss our page on \ref Eigen2ToEigen3 "API changes" between Eigen 2 and Eigen 3. - -\eigenAutoToc - -\section EIGEN2_SUPPORT_Macro The quick way: define EIGEN2_SUPPORT - -By defining EIGEN2_SUPPORT before including any Eigen 3 header, you get back a large part of the Eigen 2 API, while keeping the Eigen 3 API and ABI unchanged. - -This defaults to the \ref Stage30 "stage 30" described below. - -The rest of this page describes an optional, more powerful \em staged migration path. - -\section StagedMigrationPathOverview Overview of the staged migration path - -The primary reason why EIGEN2_SUPPORT alone may not be enough to migrate a large project from Eigen 2 to Eigen 3 is that some of the Eigen 2 API is inherently incompatible with the Eigen 3 API. This happens when the same identifier is used in Eigen 2 and in Eigen 3 with different meanings. To help migrate projects that rely on such API, we provide a staged migration path allowing to perform the migration \em incrementally. - -It goes as follows: -\li Step 0: start with a project using Eigen 2. -\li Step 1: build your project against Eigen 3 with \ref Stage10 "Eigen 2 support stage 10". This mode enables maximum compatibility with the Eigen 2 API, with just a few exceptions. -\li Step 2: build your project against Eigen 3 with \ref Stage20 "Eigen 2 support stage 20". This mode forces you to add eigen2_ prefixes to the Eigen2 identifiers that conflict with Eigen 3 API. -\li Step 3: build your project against Eigen 3 with \ref Stage30 "Eigen 2 support stage 30". This mode enables the full Eigen 3 API. -\li Step 4: build your project against Eigen 3 with \ref Stage40 "Eigen 2 support stage 40". This mode enables the full Eigen 3 strictness on matters, such as const-correctness, where Eigen 2 was looser. -\li Step 5: build your project against Eigen 3 without any Eigen 2 support mode. - -\section Stage10 Stage 10: define EIGEN2_SUPPORT_STAGE10_FULL_EIGEN2_API - -Enable this mode by defining the EIGEN2_SUPPORT_STAGE10_FULL_EIGEN2_API preprocessor macro before including any Eigen 3 header. - -This mode maximizes support for the Eigen 2 API. As a result, it does not offer the full Eigen 3 API. Also, it doesn't offer quite 100% of the Eigen 2 API. - -The part of the Eigen 3 API that is not present in this mode, is Eigen 3's Geometry module. Indeed, this mode completely replaces it by a copy of Eigen 2's Geometry module. - -The parts of the API that are still not 100% Eigen 2 compatible in this mode are: -\li Dot products over complex numbers. Eigen 2's dot product was linear in the first variable. Eigen 3's dot product is linear in the second variable. In other words, the Eigen 2 code \code x.dot(y) \endcode is equivalent to the Eigen 3 code \code y.dot(x) \endcode In yet other words, dot products are complex-conjugated in Eigen 3 compared to Eigen 2. The switch to the new convention was commanded by common usage, especially with the notation \f$ x^Ty \f$ for dot products of column-vectors. -\li The Sparse module. -\li Certain fine details of linear algebraic decompositions. For example, LDLT decomposition is now pivoting in Eigen 3 whereas it wasn't in Eigen 2, so code that was relying on its underlying matrix structure will break. -\li Usage of Eigen types in STL containers, \ref Eigen2ToEigen3 "as explained on this page". - -\section Stage20 Stage 20: define EIGEN2_SUPPORT_STAGE20_RESOLVE_API_CONFLICTS - -Enable this mode by defining the EIGEN2_SUPPORT_STAGE10_FULL_EIGEN2_API preprocessor macro before including any Eigen 3 header. - -This mode removes the Eigen 2 API that is directly conflicting with Eigen 3 API. Instead, these bits of Eigen 2 API remain available with eigen2_ prefixes. The main examples of such API are: -\li the whole Geometry module. For example, replace \c Quaternion by \c eigen2_Quaternion, replace \c Transform3f by \c eigen2_Transform3f, etc. -\li the lu() method to obtain a LU decomposition. Replace by eigen2_lu(). - -There is also one more eigen2_-prefixed identifier that you should know about, even though its use is not checked at compile time by this mode: the dot() method. As was discussed above, over complex numbers, its meaning is different between Eigen 2 and Eigen 3. You can use eigen2_dot() to get the Eigen 2 behavior. - -\section Stage30 Stage 30: define EIGEN2_SUPPORT_STAGE30_FULL_EIGEN3_API - -Enable this mode by defining the EIGEN2_SUPPORT_STAGE30_FULL_EIGEN3_API preprocessor macro before including any Eigen 3 header. Also, this mode is what you get by default when you just define EIGEN2_SUPPORT. - -This mode gives you the full unaltered Eigen 3 API, while still keeping as much support as possible for the Eigen 2 API. - -The eigen2_-prefixed identifiers are still available, but at this stage you should now replace them by Eigen 3 identifiers. Have a look at our page on \ref Eigen2ToEigen3 "API changes" between Eigen 2 and Eigen 3. - -\section Stage40 Stage 40: define EIGEN2_SUPPORT_STAGE40_FULL_EIGEN3_STRICTNESS - -Enable this mode by defining the EIGEN2_SUPPORT_STAGE40_FULL_EIGEN3_STRICTNESS preprocessor macro before including any Eigen 3 header. - -This mode tightens the last bits of strictness, especially const-correctness, that had to be loosened to support what Eigen 2 allowed. For example, this code compiled in Eigen 2: -\code -const float array[4]; -x = Map<Vector4f>(array); -\endcode -That allowed to circumvent constness. This is no longer allowed in Eigen 3. If you have to map const data in Eigen 3, map it as a const-qualified type. However, rather than explictly constructing Map objects, we strongly encourage you to use the static Map methods instead, as they take care of all of this for you: -\code -const float array[4]; -x = Vector4f::Map(array); -\endcode -This lets Eigen do the right thing for you and works equally well in Eigen 2 and in Eigen 3. - -\section FinallyDropAllEigen2Support Finally drop all Eigen 2 support - -Stage 40 is the first where it's "comfortable" to stay for a little longer period, since it preserves 100% Eigen 3 compatibility. However, we still encourage you to complete your migration as quickly as possible. While we do run the Eigen 2 test suite against Eigen 3's stage 10 support mode, we can't guarantee the same level of support and quality assurance for Eigen 2 support as we do for Eigen 3 itself, especially not in the long term. \ref Eigen2ToEigen3 "This page" describes a large part of the changes that you may need to perform. - -\section ABICompatibility What about ABI compatibility? - -It goes as follows: -\li Stage 10 already is ABI compatible with Eigen 3 for the basic (Matrix, Array, SparseMatrix...) types. However, since this stage uses a copy of Eigen 2's Geometry module instead of Eigen 3's own Geometry module, the ABI in the Geometry module is not Eigen 3 compatible. -\li Stage 20 removes the Eigen 3-incompatible Eigen 2 Geometry module (it remains available with eigen2_ prefix). So at this stage, all the identifiers that exist in Eigen 3 have the Eigen 3 ABI (and API). -\li Stage 30 introduces the remaining Eigen 3 identifiers. So at this stage, you have the full Eigen 3 ABI. -\li Stage 40 is no different than Stage 30 in these matters. - - -*/ - -} diff --git a/doc/AsciiQuickReference.txt b/doc/AsciiQuickReference.txt index d3bfd439c..8409f8850 100644 --- a/doc/AsciiQuickReference.txt +++ b/doc/AsciiQuickReference.txt @@ -32,17 +32,19 @@ A << 1, 2, 3, // Initialize A. The elements can also be B << A, A, A; // B is three horizontally stacked A's. A.fill(10); // Fill A with all 10's. -// Eigen // Matlab -MatrixXd::Identity(rows,cols) // eye(rows,cols) -C.setIdentity(rows,cols) // C = eye(rows,cols) -MatrixXd::Zero(rows,cols) // zeros(rows,cols) -C.setZero(rows,cols) // C = ones(rows,cols) -MatrixXd::Ones(rows,cols) // ones(rows,cols) -C.setOnes(rows,cols) // C = ones(rows,cols) -MatrixXd::Random(rows,cols) // rand(rows,cols)*2-1 // MatrixXd::Random returns uniform random numbers in (-1, 1). -C.setRandom(rows,cols) // C = rand(rows,cols)*2-1 -VectorXd::LinSpace(size,low,high) // linspace(low,high,size)' -v.setLinSpace(size,low,high) // v = linspace(low,high,size)' +// Eigen // Matlab +MatrixXd::Identity(rows,cols) // eye(rows,cols) +C.setIdentity(rows,cols) // C = eye(rows,cols) +MatrixXd::Zero(rows,cols) // zeros(rows,cols) +C.setZero(rows,cols) // C = zeros(rows,cols) +MatrixXd::Ones(rows,cols) // ones(rows,cols) +C.setOnes(rows,cols) // C = ones(rows,cols) +MatrixXd::Random(rows,cols) // rand(rows,cols)*2-1 // MatrixXd::Random returns uniform random numbers in (-1, 1). +C.setRandom(rows,cols) // C = rand(rows,cols)*2-1 +VectorXd::LinSpaced(size,low,high) // linspace(low,high,size)' +v.setLinSpaced(size,low,high) // v = linspace(low,high,size)' +VectorXi::LinSpaced(((hi-low)/step)+1, // low:step:hi + low,low+step*(size-1)) // // Matrix slicing and blocks. All expressions listed here are read/write. @@ -67,10 +69,10 @@ P.rightCols<cols>() // P(:, end-cols+1:end) P.rightCols(cols) // P(:, end-cols+1:end) P.topRows<rows>() // P(1:rows, :) P.topRows(rows) // P(1:rows, :) -P.middleRows<rows>(i) // P(:, i+1:i+rows) -P.middleRows(i, rows) // P(:, i+1:i+rows) -P.bottomRows<rows>() // P(:, end-rows+1:end) -P.bottomRows(rows) // P(:, end-rows+1:end) +P.middleRows<rows>(i) // P(i+1:i+rows, :) +P.middleRows(i, rows) // P(i+1:i+rows, :) +P.bottomRows<rows>() // P(end-rows+1:end, :) +P.bottomRows(rows) // P(end-rows+1:end, :) P.topLeftCorner(rows, cols) // P(1:rows, 1:cols) P.topRightCorner(rows, cols) // P(1:rows, end-cols+1:end) P.bottomLeftCorner(rows, cols) // P(end-rows+1:end, 1:cols) @@ -82,16 +84,20 @@ P.bottomRightCorner<rows,cols>() // P(end-rows+1:end, end-cols+1:end) // Of particular note is Eigen's swap function which is highly optimized. // Eigen // Matlab -R.row(i) = P.col(j); // R(i, :) = P(:, i) +R.row(i) = P.col(j); // R(i, :) = P(:, j) R.col(j1).swap(mat1.col(j2)); // R(:, [j1 j2]) = R(:, [j2, j1]) -// Views, transpose, etc; all read-write except for .adjoint(). +// Views, transpose, etc; // Eigen // Matlab R.adjoint() // R' -R.transpose() // R.' or conj(R') -R.diagonal() // diag(R) +R.transpose() // R.' or conj(R') // Read-write +R.diagonal() // diag(R) // Read-write x.asDiagonal() // diag(x) -R.transpose().colwise().reverse(); // rot90(R) +R.transpose().colwise().reverse() // rot90(R) // Read-write +R.rowwise().reverse() // fliplr(R) +R.colwise().reverse() // flipud(R) +R.replicate(i,j) // repmat(P,i,j) + // All the same as Matlab, but matlab doesn't have *= style operators. // Matrix-vector. Matrix-matrix. Matrix-scalar. @@ -103,37 +109,40 @@ a *= M; R = P + Q; R = P/s; R -= Q; R /= s; // Vectorized operations on each element independently -// Eigen // Matlab -R = P.cwiseProduct(Q); // R = P .* Q -R = P.array() * s.array();// R = P .* s -R = P.cwiseQuotient(Q); // R = P ./ Q -R = P.array() / Q.array();// R = P ./ Q -R = P.array() + s.array();// R = P + s -R = P.array() - s.array();// R = P - s -R.array() += s; // R = R + s -R.array() -= s; // R = R - s -R.array() < Q.array(); // R < Q -R.array() <= Q.array(); // R <= Q -R.cwiseInverse(); // 1 ./ P -R.array().inverse(); // 1 ./ P -R.array().sin() // sin(P) -R.array().cos() // cos(P) -R.array().pow(s) // P .^ s -R.array().square() // P .^ 2 -R.array().cube() // P .^ 3 -R.cwiseSqrt() // sqrt(P) -R.array().sqrt() // sqrt(P) -R.array().exp() // exp(P) -R.array().log() // log(P) -R.cwiseMax(P) // max(R, P) -R.array().max(P.array()) // max(R, P) -R.cwiseMin(P) // min(R, P) -R.array().min(P.array()) // min(R, P) -R.cwiseAbs() // abs(P) -R.array().abs() // abs(P) -R.cwiseAbs2() // abs(P.^2) -R.array().abs2() // abs(P.^2) -(R.array() < s).select(P,Q); // (R < s ? P : Q) +// Eigen // Matlab +R = P.cwiseProduct(Q); // R = P .* Q +R = P.array() * s.array(); // R = P .* s +R = P.cwiseQuotient(Q); // R = P ./ Q +R = P.array() / Q.array(); // R = P ./ Q +R = P.array() + s.array(); // R = P + s +R = P.array() - s.array(); // R = P - s +R.array() += s; // R = R + s +R.array() -= s; // R = R - s +R.array() < Q.array(); // R < Q +R.array() <= Q.array(); // R <= Q +R.cwiseInverse(); // 1 ./ P +R.array().inverse(); // 1 ./ P +R.array().sin() // sin(P) +R.array().cos() // cos(P) +R.array().pow(s) // P .^ s +R.array().square() // P .^ 2 +R.array().cube() // P .^ 3 +R.cwiseSqrt() // sqrt(P) +R.array().sqrt() // sqrt(P) +R.array().exp() // exp(P) +R.array().log() // log(P) +R.cwiseMax(P) // max(R, P) +R.array().max(P.array()) // max(R, P) +R.cwiseMin(P) // min(R, P) +R.array().min(P.array()) // min(R, P) +R.cwiseAbs() // abs(P) +R.array().abs() // abs(P) +R.cwiseAbs2() // abs(P.^2) +R.array().abs2() // abs(P.^2) +(R.array() < s).select(P,Q ); // (R < s ? P : Q) +R = (Q.array()==0).select(P,A) // R(Q==0) = P(Q==0) +R = P.unaryExpr(ptr_fun(func)) // R = arrayfun(func, P) // with: scalar func(const scalar &x); + // Reductions. int r, c; @@ -163,13 +172,27 @@ x.squaredNorm() // dot(x, x) Note the equivalence is not true for co x.dot(y) // dot(x, y) x.cross(y) // cross(x, y) Requires #include <Eigen/Geometry> +//// Type conversion +// Eigen // Matlab +A.cast<double>(); // double(A) +A.cast<float>(); // single(A) +A.cast<int>(); // int32(A) +A.real(); // real(A) +A.imag(); // imag(A) +// if the original type equals destination type, no work is done + +// Note that for most operations Eigen requires all operands to have the same type: +MatrixXf F = MatrixXf::Zero(3,3); +A += F; // illegal in Eigen. In Matlab A = A+F is allowed +A += F.cast<double>(); // F converted to double and then added (generally, conversion happens on-the-fly) + // Eigen can map existing memory into Eigen matrices. float array[3]; -Map<Vector3f>(array, 3).fill(10); -int data[4] = 1, 2, 3, 4; -Matrix2i mat2x2(data); -MatrixXi mat2x2 = Map<Matrix2i>(data); -MatrixXi mat2x2 = Map<MatrixXi>(data, 2, 2); +Vector3f::Map(array).fill(10); // create a temporary Map over array and sets entries to 10 +int data[4] = {1, 2, 3, 4}; +Matrix2i mat2x2(data); // copies data into mat2x2 +Matrix2i::Map(data) = 2*mat2x2; // overwrite elements of data with 2*mat2x2 +MatrixXi::Map(data, 2, 2) += mat2x2; // adds mat2x2 to elements of data (alternative syntax if size is not know at compile time) // Solve Ax = b. Result stored in x. Matlab: x = A \ b. x = A.ldlt().solve(b)); // A sym. p.s.d. #include <Eigen/Cholesky> diff --git a/doc/B01_Experimental.dox b/doc/B01_Experimental.dox index 5fc0ccd60..e1f031db8 100644 --- a/doc/B01_Experimental.dox +++ b/doc/B01_Experimental.dox @@ -4,7 +4,7 @@ namespace Eigen { \eigenAutoToc -\section summary Summary +\section Experimental_summary Summary With the 2.0 release, Eigen's API is, to a large extent, stable. However, we wish to retain the freedom to make API incompatible changes. To that effect, we call many parts of Eigen "experimental" which means that they are not subject to API stability guarantee. @@ -17,7 +17,7 @@ Experimental features may at any time: \li be subject to an API incompatible change; \li introduce API or ABI incompatible changes in your own code if you let them affect your API or ABI. -\section modules Experimental modules +\section Experimental_modules Experimental modules The following modules are considered entirely experimental, and we make no firm API stability guarantee about them for the time being: \li SVD @@ -26,7 +26,7 @@ The following modules are considered entirely experimental, and we make no firm \li Sparse \li Geometry (this one should be mostly stable, but it's a little too early to make a formal guarantee) -\section core Experimental parts of the Core module +\section Experimental_core Experimental parts of the Core module In the Core module, the only classes subject to ABI stability guarantee (meaning that you can use it for data members in your public ABI) is: \li Matrix diff --git a/doc/CMakeLists.txt b/doc/CMakeLists.txt index 8a493031c..db413bc65 100644 --- a/doc/CMakeLists.txt +++ b/doc/CMakeLists.txt @@ -10,12 +10,20 @@ if(CMAKE_COMPILER_IS_GNUCXX) endif(CMAKE_SYSTEM_NAME MATCHES Linux) endif(CMAKE_COMPILER_IS_GNUCXX) +option(EIGEN_INTERNAL_DOCUMENTATION "Build internal documentation" OFF) + + # Set some Doxygen flags set(EIGEN_DOXY_PROJECT_NAME "Eigen") set(EIGEN_DOXY_OUTPUT_DIRECTORY_SUFFIX "") set(EIGEN_DOXY_INPUT "\"${Eigen_SOURCE_DIR}/Eigen\" \"${Eigen_SOURCE_DIR}/doc\"") set(EIGEN_DOXY_HTML_COLORSTYLE_HUE "220") set(EIGEN_DOXY_TAGFILES "") +if(EIGEN_INTERNAL_DOCUMENTATION) + set(EIGEN_DOXY_INTERNAL "YES") +else(EIGEN_INTERNAL_DOCUMENTATION) + set(EIGEN_DOXY_INTERNAL "NO") +endif(EIGEN_INTERNAL_DOCUMENTATION) configure_file( ${CMAKE_CURRENT_SOURCE_DIR}/Doxyfile.in @@ -70,6 +78,8 @@ add_custom_target( COMMAND ${CMAKE_COMMAND} -E make_directory ${CMAKE_CURRENT_BINARY_DIR}/html/ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/eigen_navtree_hacks.js ${CMAKE_CURRENT_BINARY_DIR}/html/ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/Eigen_Silly_Professor_64x64.png ${CMAKE_CURRENT_BINARY_DIR}/html/ + COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/ftv2pnode.png ${CMAKE_CURRENT_BINARY_DIR}/html/ + COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/ftv2node.png ${CMAKE_CURRENT_BINARY_DIR}/html/ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/AsciiQuickReference.txt ${CMAKE_CURRENT_BINARY_DIR}/html/ WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} ) @@ -80,6 +90,8 @@ add_custom_target( COMMAND ${CMAKE_COMMAND} -E make_directory ${Eigen_BINARY_DIR}/doc/html/unsupported COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/eigen_navtree_hacks.js ${CMAKE_CURRENT_BINARY_DIR}/html/unsupported/ COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/Eigen_Silly_Professor_64x64.png ${CMAKE_CURRENT_BINARY_DIR}/html/unsupported/ + COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/ftv2pnode.png ${CMAKE_CURRENT_BINARY_DIR}/html/unsupported/ + COMMAND ${CMAKE_COMMAND} -E copy ${CMAKE_CURRENT_SOURCE_DIR}/ftv2node.png ${CMAKE_CURRENT_BINARY_DIR}/html/unsupported/ WORKING_DIRECTORY ${Eigen_BINARY_DIR}/doc ) @@ -89,9 +101,11 @@ add_dependencies(doc-unsupported-prerequisites unsupported_snippets unsupported_ add_custom_target(doc ALL COMMAND doxygen COMMAND doxygen Doxyfile-unsupported + COMMAND ${CMAKE_COMMAND} -E copy ${Eigen_BINARY_DIR}/doc/html/group__TopicUnalignedArrayAssert.html ${Eigen_BINARY_DIR}/doc/html/TopicUnalignedArrayAssert.html COMMAND ${CMAKE_COMMAND} -E rename html eigen-doc COMMAND ${CMAKE_COMMAND} -E remove eigen-doc/eigen-doc.tgz - COMMAND ${CMAKE_COMMAND} -E tar cfz eigen-doc/eigen-doc.tgz eigen-doc + COMMAND ${CMAKE_COMMAND} -E tar cfz eigen-doc.tgz eigen-doc + COMMAND ${CMAKE_COMMAND} -E rename eigen-doc.tgz eigen-doc/eigen-doc.tgz COMMAND ${CMAKE_COMMAND} -E rename eigen-doc html WORKING_DIRECTORY ${Eigen_BINARY_DIR}/doc) diff --git a/doc/CoeffwiseMathFunctionsTable.dox b/doc/CoeffwiseMathFunctionsTable.dox new file mode 100644 index 000000000..3ae9420dc --- /dev/null +++ b/doc/CoeffwiseMathFunctionsTable.dox @@ -0,0 +1,525 @@ +namespace Eigen { + +/** \eigenManualPage CoeffwiseMathFunctions Catalog of coefficient-wise math functions + + +<!-- <span style="font-size:300%; color:red; font-weight: 900;">!WORK IN PROGRESS!</span> --> + +This table presents a catalog of the coefficient-wise math functions supported by %Eigen. +In this table, \c a, \c b, refer to Array objects or expressions, and \c m refers to a linear algebra Matrix/Vector object. Standard scalar types are abbreviated as follows: + - \c int: \c i32 + - \c float: \c f + - \c double: \c d + - \c std::complex<float>: \c cf + - \c std::complex<double>: \c cd + +For each row, the first column list the equivalent calls for arrays, and matrices when supported. Of course, all functions are available for matrices by first casting it as an array: \c m.array(). + +The third column gives some hints in the underlying scalar implementation. In most cases, %Eigen does not implement itself the math function but relies on the STL for standard scalar types, or user-provided functions for custom scalar types. +For instance, some simply calls the respective function of the STL while preserving <a href="http://en.cppreference.com/w/cpp/language/adl">argument-dependent lookup</a> for custom types. +The following: +\code +using std::foo; +foo(a[i]); +\endcode +means that the STL's function \c std::foo will be potentially called if it is compatible with the underlying scalar type. If not, then the user must ensure that an overload of the function foo is available for the given scalar type (usually defined in the same namespace as the given scalar type). +This also means that, unless specified, if the function \c std::foo is available only in some recent c++ versions (e.g., c++11), then the respective %Eigen's function/method will be usable on standard types only if the compiler support the required c++ version. + +<table class="manual-hl"> +<tr> +<th>API</th><th>Description</th><th>Default scalar implementation</th><th>SIMD</th> +</tr> +<tr><td colspan="4"></td></tr> +<tr><th colspan="4">Basic operations</th></tr> +<tr> + <td class="code"> + \anchor cwisetable_abs + a.\link ArrayBase::abs abs\endlink(); \n + \link Eigen::abs abs\endlink(a); \n + m.\link MatrixBase::cwiseAbs cwiseAbs\endlink(); + </td> + <td>absolute value (\f$ |a_i| \f$) </td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/fabs">std::abs</a>; \n + abs(a[i]); + </td> + <td>SSE2, AVX (i32,f,d)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_inverse + a.\link ArrayBase::inverse inverse\endlink(); \n + \link Eigen::inverse inverse\endlink(a); \n + m.\link MatrixBase::cwiseInverse cwiseInverse\endlink(); + </td> + <td>inverse value (\f$ 1/a_i \f$) </td> + <td class="code"> + 1/a[i]; + </td> + <td>All engines (f,d,fc,fd)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_conj + a.\link ArrayBase::conjugate conjugate\endlink(); \n + \link Eigen::conj conj\endlink(a); \n + m.\link MatrixBase::conjugate conjugate(); + </td> + <td><a href="https://en.wikipedia.org/wiki/Complex_conjugate">complex conjugate</a> (\f$ \bar{a_i} \f$),\n + no-op for real </td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/complex/conj">std::conj</a>; \n + conj(a[i]); + </td> + <td>All engines (fc,fd)</td> +</tr> +<tr> +<th colspan="4">Exponential functions</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_exp + a.\link ArrayBase::exp exp\endlink(); \n + \link Eigen::exp exp\endlink(a); + </td> + <td>\f$ e \f$ raised to the given power (\f$ e^{a_i} \f$) </td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/exp">std::exp</a>; \n + exp(a[i]); + </td> + <td>SSE2, AVX (f,d)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_log + a.\link ArrayBase::log log\endlink(); \n + \link Eigen::log log\endlink(a); + </td> + <td>natural (base \f$ e \f$) logarithm (\f$ \ln({a_i}) \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/log">std::log</a>; \n + log(a[i]); + </td> + <td>SSE2, AVX (f)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_log1p + a.\link ArrayBase::log1p log1p\endlink(); \n + \link Eigen::log1p log1p\endlink(a); + </td> + <td>natural (base \f$ e \f$) logarithm of 1 plus \n the given number (\f$ \ln({1+a_i}) \f$)</td> + <td>built-in generic implementation based on \c log,\n + plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/log1p">\c std::log1p </a>; \cpp11</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_log10 + a.\link ArrayBase::log10 log10\endlink(); \n + \link Eigen::log10 log10\endlink(a); + </td> + <td>base 10 logarithm (\f$ \log_{10}({a_i}) \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/log10">std::log10</a>; \n + log10(a[i]); + </td> + <td></td> +</tr> +<tr> +<th colspan="4">Power functions</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_pow + a.\link ArrayBase::pow pow\endlink(b); \n + \link Eigen::pow pow\endlink(a,b); + </td> + <td>raises a number to the given power (\f$ a_i ^ {b_i} \f$) \n \c a and \c b can be either an array or scalar.</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/pow">std::pow</a>; \n + pow(a[i],b[i]);\n + (plus builtin for integer types)</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_sqrt + a.\link ArrayBase::sqrt sqrt\endlink(); \n + \link Eigen::sqrt sqrt\endlink(a);\n + m.\link MatrixBase::cwiseSqrt cwiseSqrt\endlink(); + </td> + <td>computes square root (\f$ \sqrt a_i \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/sqrt">std::sqrt</a>; \n + sqrt(a[i]);</td> + <td>SSE2, AVX (f,d)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_rsqrt + a.\link ArrayBase::rsqrt rsqrt\endlink(); \n + \link Eigen::rsqrt rsqrt\endlink(a); + </td> + <td><a href="https://en.wikipedia.org/wiki/Fast_inverse_square_root">reciprocal square root</a> (\f$ 1/{\sqrt a_i} \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/sqrt">std::sqrt</a>; \n + 1/sqrt(a[i]); \n + </td> + <td>SSE2, AVX, AltiVec, ZVector (f,d)\n + (approx + 1 Newton iteration)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_square + a.\link ArrayBase::square square\endlink(); \n + \link Eigen::square square\endlink(a); + </td> + <td>computes square power (\f$ a_i^2 \f$)</td> + <td class="code"> + a[i]*a[i]</td> + <td>All (i32,f,d,cf,cd)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_cube + a.\link ArrayBase::cube cube\endlink(); \n + \link Eigen::cube cube\endlink(a); + </td> + <td>computes cubic power (\f$ a_i^3 \f$)</td> + <td class="code"> + a[i]*a[i]*a[i]</td> + <td>All (i32,f,d,cf,cd)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_abs2 + a.\link ArrayBase::abs2 abs2\endlink(); \n + \link Eigen::abs2 abs2\endlink(a);\n + m.\link MatrixBase::cwiseAbs2 cwiseAbs2\endlink(); + </td> + <td>computes the squared absolute value (\f$ |a_i|^2 \f$)</td> + <td class="code"> + real: a[i]*a[i] \n + complex: real(a[i])*real(a[i]) \n + + imag(a[i])*imag(a[i])</td> + <td>All (i32,f,d)</td> +</tr> +<tr> +<th colspan="4">Trigonometric functions</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_sin + a.\link ArrayBase::sin sin\endlink(); \n + \link Eigen::sin sin\endlink(a); + </td> + <td>computes sine</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/sin">std::sin</a>; \n + sin(a[i]);</td> + <td>SSE2, AVX (f)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_cos + a.\link ArrayBase::cos cos\endlink(); \n + \link Eigen::cos cos\endlink(a); + </td> + <td>computes cosine</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/cos">std::cos</a>; \n + cos(a[i]);</td> + <td>SSE2, AVX (f)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_tan + a.\link ArrayBase::tan tan\endlink(); \n + \link Eigen::tan tan\endlink(a); + </td> + <td>computes tangent</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/tan">std::tan</a>; \n + tan(a[i]);</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_asin + a.\link ArrayBase::asin asin\endlink(); \n + \link Eigen::asin asin\endlink(a); + </td> + <td>computes arc sine (\f$ \sin^{-1} a_i \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/asin">std::asin</a>; \n + asin(a[i]);</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_acos + a.\link ArrayBase::acos acos\endlink(); \n + \link Eigen::acos acos\endlink(a); + </td> + <td>computes arc cosine (\f$ \cos^{-1} a_i \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/acos">std::acos</a>; \n + acos(a[i]);</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_atan + a.\link ArrayBase::atan tan\endlink(); \n + \link Eigen::atan atan\endlink(a); + </td> + <td>computes arc tangent (\f$ \tan^{-1} a_i \f$)</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/atan">std::atan</a>; \n + atan(a[i]);</td> + <td></td> +</tr> +<tr> +<th colspan="4">Hyperbolic functions</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_sinh + a.\link ArrayBase::sinh sinh\endlink(); \n + \link Eigen::sinh sinh\endlink(a); + </td> + <td>computes hyperbolic sine</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/sinh">std::sinh</a>; \n + sinh(a[i]);</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_cosh + a.\link ArrayBase::cosh cohs\endlink(); \n + \link Eigen::cosh cosh\endlink(a); + </td> + <td>computes hyperbolic cosine</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/cosh">std::cosh</a>; \n + cosh(a[i]);</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_tanh + a.\link ArrayBase::tanh tanh\endlink(); \n + \link Eigen::tanh tanh\endlink(a); + </td> + <td>computes hyperbolic tangent</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/tanh">std::tanh</a>; \n + tanh(a[i]);</td> + <td></td> +</tr> +<tr> +<th colspan="4">Nearest integer floating point operations</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_ceil + a.\link ArrayBase::ceil ceil\endlink(); \n + \link Eigen::ceil ceil\endlink(a); + </td> + <td>nearest integer not less than the given value</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/ceil">std::ceil</a>; \n + ceil(a[i]);</td> + <td>SSE4,AVX,ZVector (f,d)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_floor + a.\link ArrayBase::floor floor\endlink(); \n + \link Eigen::floor floor\endlink(a); + </td> + <td>nearest integer not greater than the given value</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/floor">std::floor</a>; \n + floor(a[i]);</td> + <td>SSE4,AVX,ZVector (f,d)</td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_round + a.\link ArrayBase::round round\endlink(); \n + \link Eigen::round round\endlink(a); + </td> + <td>nearest integer, \n rounding away from zero in halfway cases</td> + <td>built-in generic implementation \n based on \c floor and \c ceil,\n + plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/round">\c std::round </a>; \cpp11</td> + <td>SSE4,AVX,ZVector (f,d)</td> +</tr> +<tr> +<th colspan="4">Floating point manipulation functions</th> +</tr> +<tr> +<th colspan="4">Classification and comparison</th> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_isfinite + a.\link ArrayBase::isFinite isFinite\endlink(); \n + \link Eigen::isfinite isfinite\endlink(a); + </td> + <td>checks if the given number has finite value</td> + <td>built-in generic implementation,\n + plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/isfinite">\c std::isfinite </a>; \cpp11</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_isinf + a.\link ArrayBase::isInf isInf\endlink(); \n + \link Eigen::isinf isinf\endlink(a); + </td> + <td>checks if the given number is infinite</td> + <td>built-in generic implementation,\n + plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/isinf">\c std::isinf </a>; \cpp11</td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_isnan + a.\link ArrayBase::isNaN isNaN\endlink(); \n + \link Eigen::isnan isnan\endlink(a); + </td> + <td>checks if the given number is not a number</td> + <td>built-in generic implementation,\n + plus \c using <a href="http://en.cppreference.com/w/cpp/numeric/math/isnan">\c std::isnan </a>; \cpp11</td> + <td></td> +</tr> +<tr> +<th colspan="4">Error and gamma functions</th> +</tr> +<tr> <td colspan="4"> Require \c \#include \c <unsupported/Eigen/SpecialFunctions> </td></tr> +<tr> + <td class="code"> + \anchor cwisetable_erf + a.\link ArrayBase::erf erf\endlink(); \n + \link Eigen::erf erf\endlink(a); + </td> + <td>error function</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/erf">std::erf</a>; \cpp11 \n + erf(a[i]); + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_erfc + a.\link ArrayBase::erfc erfc\endlink(); \n + \link Eigen::erfc erfc\endlink(a); + </td> + <td>complementary error function</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/erfc">std::erfc</a>; \cpp11 \n + erfc(a[i]); + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_lgamma + a.\link ArrayBase::lgamma lgamma\endlink(); \n + \link Eigen::lgamma lgamma\endlink(a); + </td> + <td>natural logarithm of the gamma function</td> + <td class="code"> + using <a href="http://en.cppreference.com/w/cpp/numeric/math/lgamma">std::lgamma</a>; \cpp11 \n + lgamma(a[i]); + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_digamma + a.\link ArrayBase::digamma digamma\endlink(); \n + \link Eigen::digamma digamma\endlink(a); + </td> + <td><a href="https://en.wikipedia.org/wiki/Digamma_function">logarithmic derivative of the gamma function</a></td> + <td> + built-in for float and double + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_igamma + \link Eigen::igamma igamma\endlink(a,x); + </td> + <td><a href="https://en.wikipedia.org/wiki/Incomplete_gamma_function">lower incomplete gamma integral</a> + \n \f$ \gamma(a_i,x_i)= \frac{1}{|a_i|} \int_{0}^{x_i}e^{\text{-}t} t^{a_i-1} \mathrm{d} t \f$</td> + <td> + built-in for float and double,\n but requires \cpp11 + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_igammac + \link Eigen::igammac igammac\endlink(a,x); + </td> + <td><a href="https://en.wikipedia.org/wiki/Incomplete_gamma_function">upper incomplete gamma integral</a> + \n \f$ \Gamma(a_i,x_i) = \frac{1}{|a_i|} \int_{x_i}^{\infty}e^{\text{-}t} t^{a_i-1} \mathrm{d} t \f$</td> + <td> + built-in for float and double,\n but requires \cpp11 + </td> + <td></td> +</tr> +<tr> +<th colspan="4">Special functions</th> +</tr> +<tr> <td colspan="4"> Require \c \#include \c <unsupported/Eigen/SpecialFunctions> </td></tr> +<tr> + <td class="code"> + \anchor cwisetable_polygamma + \link Eigen::polygamma polygamma\endlink(n,x); + </td> + <td><a href="https://en.wikipedia.org/wiki/Polygamma_function">n-th derivative of digamma at x</a></td> + <td> + built-in generic based on\n <a href="#cwisetable_lgamma">\c lgamma </a>, + <a href="#cwisetable_digamma"> \c digamma </a> + and <a href="#cwisetable_zeta">\c zeta </a>. + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_betainc + \link Eigen::betainc betainc\endlink(a,b,x); + </td> + <td><a href="https://en.wikipedia.org/wiki/Beta_function#Incomplete_beta_function">Incomplete beta function</a></td> + <td> + built-in for float and double,\n but requires \cpp11 + </td> + <td></td> +</tr> +<tr> + <td class="code"> + \anchor cwisetable_zeta + \link Eigen::zeta zeta\endlink(a,b); + </td> + <td><a href="https://en.wikipedia.org/wiki/Hurwitz_zeta_function">Hurwitz zeta function</a> + \n \f$ \zeta(a_i,b_i)=\sum_{k=0}^{\infty}(b_i+k)^{\text{-}a_i} \f$</td> + <td> + built-in for float and double + </td> + <td></td> +</tr> +<tr><td colspan="4"></td></tr> +</table> + +\n + +*/ + +} diff --git a/doc/CustomizingEigen.dox b/doc/CustomizingEigen.dox deleted file mode 100644 index 0850863aa..000000000 --- a/doc/CustomizingEigen.dox +++ /dev/null @@ -1,165 +0,0 @@ -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. - -\eigenAutoToc - -\section ExtendingMatrixBase Extending MatrixBase (and other classes) - -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 your 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. - -You can extend many of the other classes used in Eigen by defining similarly named preprocessor symbols. For instance, define \c EIGEN_ARRAYBASE_PLUGIN if you want to extend the ArrayBase class. A full list of classes that can be extended in this way and the corresponding preprocessor symbols can be found on our page \ref TopicPreprocessorDirectives. - -Here is an example of an extension file for adding methods to MatrixBase: \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 internal::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<internal::scalar_add_op<Scalar>, Derived> -operator+(const Scalar& scalar) const -{ return CwiseUnaryOp<internal::scalar_add_op<Scalar>, Derived>(derived(), internal::scalar_add_op<Scalar>(scalar)); } - -friend const CwiseUnaryOp<internal::scalar_add_op<Scalar>, Derived> -operator+(const Scalar& scalar, const MatrixBase<Derived>& mat) -{ return CwiseUnaryOp<internal::scalar_add_op<Scalar>, Derived>(mat.derived(), internal::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, MyVerySpecificVector2 -> 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: - -\include CustomizingEigen_Inheritance.cpp - -Output: \verbinclude CustomizingEigen_Inheritance.out - -This is the kind of error you can get if you don't provide those methods -\verbatim -error: no match for ‘operator=’ in ‘v = Eigen::operator*( -const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1, 0, -0x000000001, 1> >::Scalar&, -const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType&) -(((const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType&) -((const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType*)(& v))))’ -\endverbatim - -\anchor user_defined_scalars \section CustomScalarType Using custom scalar types - -By default, Eigen currently supports standard floating-point types (\c float, \c double, \c std::complex<float>, \c std::complex<double>, \c long \c double), as well as all native integer types (e.g., \c int, \c unsigned \c int, \c short, etc.), and \c bool. -On x86-64 systems, \c long \c double permits to locally enforces the use of x87 registers with extended accuracy (in comparison to SSE). - -In order to add support for a custom type \c T you need: --# make sure the common operator (+,-,*,/,etc.) are supported by the type \c T --# add a specialization of struct Eigen::NumTraits<T> (see \ref NumTraits) --# define the math functions that makes sense for your type. This includes standard ones like sqrt, pow, sin, tan, conj, real, imag, etc, as well as abs2 which is Eigen specific. - (see the file Eigen/src/Core/MathFunctions.h) - -The math function should be defined in the same namespace than \c T, or in the \c std namespace though that second approach is not recommended. - -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 ADOLCSUPPORT_H -#define ADOLCSUPPORT_H - -#define ADOLC_TAPELESS -#include <adolc/adouble.h> -#include <Eigen/Core> - -namespace Eigen { - -template<> struct NumTraits<adtl::adouble> - : NumTraits<double> // permits to get the epsilon, dummy_precision, lowest, highest functions -{ - typedef adtl::adouble Real; - typedef adtl::adouble NonInteger; - typedef adtl::adouble Nested; - - enum { - IsComplex = 0, - IsInteger = 0, - IsSigned = 1, - RequireInitialization = 1, - ReadCost = 1, - AddCost = 3, - MulCost = 3 - }; -}; - -} - -namespace adtl { - -inline const adouble& conj(const adouble& x) { return x; } -inline const adouble& real(const adouble& x) { return x; } -inline adouble imag(const adouble&) { return 0.; } -inline adouble abs(const adouble& x) { return fabs(x); } -inline adouble abs2(const adouble& x) { return x*x; } - -} - -#endif // ADOLCSUPPORT_H -\endcode - - -\sa \ref TopicPreprocessorDirectives - -*/ - -} diff --git a/doc/CustomizingEigen_CustomScalar.dox b/doc/CustomizingEigen_CustomScalar.dox new file mode 100644 index 000000000..1ee78cbe5 --- /dev/null +++ b/doc/CustomizingEigen_CustomScalar.dox @@ -0,0 +1,120 @@ +namespace Eigen { + +/** \page TopicCustomizing_CustomScalar Using custom scalar types +\anchor user_defined_scalars + +By default, Eigen currently supports standard floating-point types (\c float, \c double, \c std::complex<float>, \c std::complex<double>, \c long \c double), as well as all native integer types (e.g., \c int, \c unsigned \c int, \c short, etc.), and \c bool. +On x86-64 systems, \c long \c double permits to locally enforces the use of x87 registers with extended accuracy (in comparison to SSE). + +In order to add support for a custom type \c T you need: +-# make sure the common operator (+,-,*,/,etc.) are supported by the type \c T +-# add a specialization of struct Eigen::NumTraits<T> (see \ref NumTraits) +-# define the math functions that makes sense for your type. This includes standard ones like sqrt, pow, sin, tan, conj, real, imag, etc, as well as abs2 which is Eigen specific. + (see the file Eigen/src/Core/MathFunctions.h) + +The math function should be defined in the same namespace than \c T, or in the \c std namespace though that second approach is not recommended. + +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 ADOLCSUPPORT_H +#define ADOLCSUPPORT_H + +#define ADOLC_TAPELESS +#include <adolc/adouble.h> +#include <Eigen/Core> + +namespace Eigen { + +template<> struct NumTraits<adtl::adouble> + : NumTraits<double> // permits to get the epsilon, dummy_precision, lowest, highest functions +{ + typedef adtl::adouble Real; + typedef adtl::adouble NonInteger; + typedef adtl::adouble Nested; + + enum { + IsComplex = 0, + IsInteger = 0, + IsSigned = 1, + RequireInitialization = 1, + ReadCost = 1, + AddCost = 3, + MulCost = 3 + }; +}; + +} + +namespace adtl { + +inline const adouble& conj(const adouble& x) { return x; } +inline const adouble& real(const adouble& x) { return x; } +inline adouble imag(const adouble&) { return 0.; } +inline adouble abs(const adouble& x) { return fabs(x); } +inline adouble abs2(const adouble& x) { return x*x; } + +} + +#endif // ADOLCSUPPORT_H +\endcode + +This other example adds support for the \c mpq_class type from <a href="https://gmplib.org/">GMP</a>. It shows in particular how to change the way Eigen picks the best pivot during LU factorization. It selects the coefficient with the highest score, where the score is by default the absolute value of a number, but we can define a different score, for instance to prefer pivots with a more compact representation (this is an example, not a recommendation). Note that the scores should always be non-negative and only zero is allowed to have a score of zero. Also, this can interact badly with thresholds for inexact scalar types. + +\code +#include <gmpxx.h> +#include <Eigen/Core> +#include <boost/operators.hpp> + +namespace Eigen { + template<> struct NumTraits<mpq_class> : GenericNumTraits<mpq_class> + { + typedef mpq_class Real; + typedef mpq_class NonInteger; + typedef mpq_class Nested; + + static inline Real epsilon() { return 0; } + static inline Real dummy_precision() { return 0; } + static inline Real digits10() { return 0; } + + enum { + IsInteger = 0, + IsSigned = 1, + IsComplex = 0, + RequireInitialization = 1, + ReadCost = 6, + AddCost = 150, + MulCost = 100 + }; + }; + + namespace internal { + + template<> struct scalar_score_coeff_op<mpq_class> { + struct result_type : boost::totally_ordered1<result_type> { + std::size_t len; + result_type(int i = 0) : len(i) {} // Eigen uses Score(0) and Score() + result_type(mpq_class const& q) : + len(mpz_size(q.get_num_mpz_t())+ + mpz_size(q.get_den_mpz_t())-1) {} + friend bool operator<(result_type x, result_type y) { + // 0 is the worst possible pivot + if (x.len == 0) return y.len > 0; + if (y.len == 0) return false; + // Prefer a pivot with a small representation + return x.len > y.len; + } + friend bool operator==(result_type x, result_type y) { + // Only used to test if the score is 0 + return x.len == y.len; + } + }; + result_type operator()(mpq_class const& x) const { return x; } + }; + } +} +\endcode + +*/ + +} diff --git a/doc/CustomizingEigen_InheritingMatrix.dox b/doc/CustomizingEigen_InheritingMatrix.dox new file mode 100644 index 000000000..b21e55433 --- /dev/null +++ b/doc/CustomizingEigen_InheritingMatrix.dox @@ -0,0 +1,34 @@ +namespace Eigen { + +/** \page TopicCustomizing_InheritingMatrix 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, MyVerySpecificVector2 -> 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: + +\include CustomizingEigen_Inheritance.cpp + +Output: \verbinclude CustomizingEigen_Inheritance.out + +This is the kind of error you can get if you don't provide those methods +\verbatim +error: no match for ‘operator=’ in ‘v = Eigen::operator*( +const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1, 0, -0x000000001, 1> >::Scalar&, +const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType&) +(((const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType&) +((const Eigen::MatrixBase<Eigen::Matrix<double, -0x000000001, 1> >::StorageBaseType*)(& v))))’ +\endverbatim + +*/ + +} diff --git a/doc/CustomizingEigen_NullaryExpr.dox b/doc/CustomizingEigen_NullaryExpr.dox new file mode 100644 index 000000000..37c8dcd89 --- /dev/null +++ b/doc/CustomizingEigen_NullaryExpr.dox @@ -0,0 +1,86 @@ +namespace Eigen { + +/** \page TopicCustomizing_NullaryExpr Matrix manipulation via nullary-expressions + + +The main purpose of the class CwiseNullaryOp is to define \em procedural matrices such as constant or random matrices as returned by the Ones(), Zero(), Constant(), Identity() and Random() methods. +Nevertheless, with some imagination it is possible to accomplish very sophisticated matrix manipulation with minimal efforts such that \ref TopicNewExpressionType "implementing new expression" is rarely needed. + +\section NullaryExpr_Circulant Example 1: circulant matrix + +To explore these possibilities let us start with the \em circulant example of the \ref TopicNewExpressionType "implementing new expression" topic. +Let us recall that a circulant matrix is a matrix where each column is the same as the +column to the left, except that it is cyclically shifted downwards. +For example, here is a 4-by-4 circulant matrix: +\f[ \begin{bmatrix} + 1 & 8 & 4 & 2 \\ + 2 & 1 & 8 & 4 \\ + 4 & 2 & 1 & 8 \\ + 8 & 4 & 2 & 1 +\end{bmatrix} \f] +A circulant matrix is uniquely determined by its first column. We wish +to write a function \c makeCirculant which, given the first column, +returns an expression representing the circulant matrix. + +For this exercise, the return type of \c makeCirculant will be a CwiseNullaryOp that we need to instantiate with: +1 - a proper \c circulant_functor storing the input vector and implementing the adequate coefficient accessor \c operator(i,j) +2 - a template instantiation of class Matrix conveying compile-time information such as the scalar type, sizes, and preferred storage layout. + +Calling \c ArgType the type of the input vector, we can construct the equivalent squared Matrix type as follows: + +\snippet make_circulant2.cpp square + +This little helper structure will help us to implement our \c makeCirculant function as follows: + +\snippet make_circulant2.cpp makeCirculant + +As usual, our function takes as argument a \c MatrixBase (see this \ref TopicFunctionTakingEigenTypes "page" for more details). +Then, the CwiseNullaryOp object is constructed through the DenseBase::NullaryExpr static method with the adequate runtime sizes. + +Then, we need to implement our \c circulant_functor, which is a straightforward exercise: + +\snippet make_circulant2.cpp circulant_func + +We are now all set to try our new feature: + +\snippet make_circulant2.cpp main + + +If all the fragments are combined, the following output is produced, +showing that the program works as expected: + +\include make_circulant2.out + +This implementation of \c makeCirculant is much simpler than \ref TopicNewExpressionType "defining a new expression" from scratch. + + +\section NullaryExpr_Indexing Example 2: indexing rows and columns + +The goal here is to mimic MatLab's ability to index a matrix through two vectors of indices referencing the rows and columns to be picked respectively, like this: + +\snippet nullary_indexing.out main1 + +To this end, let us first write a nullary-functor storing references to the input matrix and to the two arrays of indices, and implementing the required \c operator()(i,j): + +\snippet nullary_indexing.cpp functor + +Then, let's create an \c indexing(A,rows,cols) function creating the nullary expression: + +\snippet nullary_indexing.cpp function + +Finally, here is an example of how this function can be used: + +\snippet nullary_indexing.cpp main1 + +This straightforward implementation is already quite powerful as the row or column index arrays can also be expressions to perform offsetting, modulo, striding, reverse, etc. + +\snippet nullary_indexing.cpp main2 + +and the output is: + +\snippet nullary_indexing.out main2 + +*/ + +} + diff --git a/doc/CustomizingEigen_Plugins.dox b/doc/CustomizingEigen_Plugins.dox new file mode 100644 index 000000000..d88f2409b --- /dev/null +++ b/doc/CustomizingEigen_Plugins.dox @@ -0,0 +1,69 @@ +namespace Eigen { + +/** \page TopicCustomizing_Plugins Extending MatrixBase (and other classes) + +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 your 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. + +You can extend many of the other classes used in Eigen by defining similarly named preprocessor symbols. For instance, define \c EIGEN_ARRAYBASE_PLUGIN if you want to extend the ArrayBase class. A full list of classes that can be extended in this way and the corresponding preprocessor symbols can be found on our page \ref TopicPreprocessorDirectives. + +Here is an example of an extension file for adding methods to MatrixBase: \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 internal::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 CwiseBinaryOp<internal::scalar_sum_op<Scalar>, const Derived, const ConstantReturnType> +operator+(const Scalar& scalar) const +{ return CwiseBinaryOp<internal::scalar_sum_op<Scalar>, const Derived, const ConstantReturnType>(derived(), Constant(rows(),cols(),scalar)); } + +friend const CwiseBinaryOp<internal::scalar_sum_op<Scalar>, const ConstantReturnType, Derived> +operator+(const Scalar& scalar, const MatrixBase<Derived>& mat) +{ return CwiseBinaryOp<internal::scalar_sum_op<Scalar>, const ConstantReturnType, Derived>(Constant(rows(),cols(),scalar), mat.derived()); } +\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 + +*/ + +} diff --git a/doc/DenseDecompositionBenchmark.dox b/doc/DenseDecompositionBenchmark.dox new file mode 100644 index 000000000..7be9c70cd --- /dev/null +++ b/doc/DenseDecompositionBenchmark.dox @@ -0,0 +1,42 @@ +namespace Eigen { + +/** \eigenManualPage DenseDecompositionBenchmark Benchmark of dense decompositions + +This page presents a speed comparison of the dense matrix decompositions offered by %Eigen for a wide range of square matrices and overconstrained problems. + +For a more general overview on the features and numerical robustness of linear solvers and decompositions, check this \link TopicLinearAlgebraDecompositions table \endlink. + +This benchmark has been run on a laptop equipped with an Intel core i7 \@ 2,6 GHz, and compiled with clang with \b AVX and \b FMA instruction sets enabled but without multi-threading. +It uses \b single \b precision \b float numbers. For double, you can get a good estimate by multiplying the timings by a factor 2. + +The square matrices are symmetric, and for the overconstrained matrices, the reported timmings include the cost to compute the symmetric covariance matrix \f$ A^T A \f$ for the first four solvers based on Cholesky and LU, as denoted by the \b * symbol (top-right corner part of the table). +Timings are in \b milliseconds, and factors are relative to the LLT decomposition which is the fastest but also the least general and robust. + +<table class="manual"> +<tr><th>solver/size</th> + <th>8x8</th> <th>100x100</th> <th>1000x1000</th> <th>4000x4000</th> <th>10000x8</th> <th>10000x100</th> <th>10000x1000</th> <th>10000x4000</th></tr> +<tr><td>LLT</td><td>0.05</td><td>0.42</td><td>5.83</td><td>374.55</td><td>6.79 <sup><a href="#note_ls">*</a></sup></td><td>30.15 <sup><a href="#note_ls">*</a></sup></td><td>236.34 <sup><a href="#note_ls">*</a></sup></td><td>3847.17 <sup><a href="#note_ls">*</a></sup></td></tr> +<tr class="alt"><td>LDLT</td><td>0.07 (x1.3)</td><td>0.65 (x1.5)</td><td>26.86 (x4.6)</td><td>2361.18 (x6.3)</td><td>6.81 (x1) <sup><a href="#note_ls">*</a></sup></td><td>31.91 (x1.1) <sup><a href="#note_ls">*</a></sup></td><td>252.61 (x1.1) <sup><a href="#note_ls">*</a></sup></td><td>5807.66 (x1.5) <sup><a href="#note_ls">*</a></sup></td></tr> +<tr><td>PartialPivLU</td><td>0.08 (x1.5)</td><td>0.69 (x1.6)</td><td>15.63 (x2.7)</td><td>709.32 (x1.9)</td><td>6.81 (x1) <sup><a href="#note_ls">*</a></sup></td><td>31.32 (x1) <sup><a href="#note_ls">*</a></sup></td><td>241.68 (x1) <sup><a href="#note_ls">*</a></sup></td><td>4270.48 (x1.1) <sup><a href="#note_ls">*</a></sup></td></tr> +<tr class="alt"><td>FullPivLU</td><td>0.1 (x1.9)</td><td>4.48 (x10.6)</td><td>281.33 (x48.2)</td><td>-</td><td>6.83 (x1) <sup><a href="#note_ls">*</a></sup></td><td>32.67 (x1.1) <sup><a href="#note_ls">*</a></sup></td><td>498.25 (x2.1) <sup><a href="#note_ls">*</a></sup></td><td>-</td></tr> +<tr><td>HouseholderQR</td><td>0.19 (x3.5)</td><td>2.18 (x5.2)</td><td>23.42 (x4)</td><td>1337.52 (x3.6)</td><td>34.26 (x5)</td><td>129.01 (x4.3)</td><td>377.37 (x1.6)</td><td>4839.1 (x1.3)</td></tr> +<tr class="alt"><td>ColPivHouseholderQR</td><td>0.23 (x4.3)</td><td>2.23 (x5.3)</td><td>103.34 (x17.7)</td><td>9987.16 (x26.7)</td><td>36.05 (x5.3)</td><td>163.18 (x5.4)</td><td>2354.08 (x10)</td><td>37860.5 (x9.8)</td></tr> +<tr><td>CompleteOrthogonalDecomposition</td><td>0.23 (x4.3)</td><td>2.22 (x5.2)</td><td>99.44 (x17.1)</td><td>10555.3 (x28.2)</td><td>35.75 (x5.3)</td><td>169.39 (x5.6)</td><td>2150.56 (x9.1)</td><td>36981.8 (x9.6)</td></tr> +<tr class="alt"><td>FullPivHouseholderQR</td><td>0.23 (x4.3)</td><td>4.64 (x11)</td><td>289.1 (x49.6)</td><td>-</td><td>69.38 (x10.2)</td><td>446.73 (x14.8)</td><td>4852.12 (x20.5)</td><td>-</td></tr> +<tr><td>JacobiSVD</td><td>1.01 (x18.6)</td><td>71.43 (x168.4)</td><td>-</td><td>-</td><td>113.81 (x16.7)</td><td>1179.66 (x39.1)</td><td>-</td><td>-</td></tr> +<tr class="alt"><td>BDCSVD</td><td>1.07 (x19.7)</td><td>21.83 (x51.5)</td><td>331.77 (x56.9)</td><td>18587.9 (x49.6)</td><td>110.53 (x16.3)</td><td>397.67 (x13.2)</td><td>2975 (x12.6)</td><td>48593.2 (x12.6)</td></tr> +</table> + +<a name="note_ls">\b *: </a> This decomposition do not support direct least-square solving for over-constrained problems, and the reported timing include the cost to form the symmetric covariance matrix \f$ A^T A \f$. + +\b Observations: + + LLT is always the fastest solvers. + + For largely over-constrained problems, the cost of Cholesky/LU decompositions is dominated by the computation of the symmetric covariance matrix. + + For large problem sizes, only the decomposition implementing a cache-friendly blocking strategy scale well. Those include LLT, PartialPivLU, HouseholderQR, and BDCSVD. This explain why for a 4k x 4k matrix, HouseholderQR is faster than LDLT. In the future, LDLT and ColPivHouseholderQR will also implement blocking strategies. + + CompleteOrthogonalDecomposition is based on ColPivHouseholderQR and they thus achieve the same level of performance. + +The above table has been generated by the <a href="https://bitbucket.org/eigen/eigen/raw/default/bench/dense_solvers.cpp">bench/dense_solvers.cpp</a> file, feel-free to hack it to generate a table matching your hardware, compiler, and favorite problem sizes. + +*/ + +} diff --git a/doc/Doxyfile.in b/doc/Doxyfile.in index 85af9f1d4..2109978fe 100644 --- a/doc/Doxyfile.in +++ b/doc/Doxyfile.in @@ -125,7 +125,7 @@ ALWAYS_DETAILED_SEC = NO # members were ordinary class members. Constructors, destructors and assignment # operators of the base classes will not be shown. -INLINE_INHERITED_MEMB = YES +INLINE_INHERITED_MEMB = NO # If the FULL_PATH_NAMES tag is set to YES then Doxygen will prepend the full # path before files name in the file list and in the header files. If set @@ -206,7 +206,7 @@ TAB_SIZE = 8 # You can put \n's in the value part of an alias to insert newlines. ALIASES = "only_for_vectors=This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column." \ - "not_reentrant=\warning This function is not re-entrant." \ + "not_reentrant=\warning This function is not re-entrant." \ "array_module=This is defined in the %Array module. \code #include <Eigen/Array> \endcode" \ "cholesky_module=This is defined in the %Cholesky module. \code #include <Eigen/Cholesky> \endcode" \ "eigenvalues_module=This is defined in the %Eigenvalues module. \code #include <Eigen/Eigenvalues> \endcode" \ @@ -216,6 +216,7 @@ ALIASES = "only_for_vectors=This is only for vectors (either row- "lu_module=This is defined in the %LU module. \code #include <Eigen/LU> \endcode" \ "qr_module=This is defined in the %QR module. \code #include <Eigen/QR> \endcode" \ "svd_module=This is defined in the %SVD module. \code #include <Eigen/SVD> \endcode" \ + "specialfunctions_module=This is defined in the \b unsupported SpecialFunctions module. \code #include <Eigen/SpecialFunctions> \endcode" \ "label=\bug" \ "matrixworld=<a href='#matrixonly' style='color:green;text-decoration: none;'>*</a>" \ "arrayworld=<a href='#arrayonly' style='color:blue;text-decoration: none;'>*</a>" \ @@ -223,7 +224,13 @@ ALIASES = "only_for_vectors=This is only for vectors (either row- "note_about_using_kernel_to_study_multiple_solutions=If you need a complete analysis of the space of solutions, take the one solution obtained by this method and add to it elements of the kernel, as determined by kernel()." \ "note_about_checking_solutions=This method just tries to find as good a solution as possible. If you want to check whether a solution exists or if it is accurate, just call this function to get a result and then compute the error of this result, or use MatrixBase::isApprox() directly, for instance like this: \code bool a_solution_exists = (A*result).isApprox(b, precision); \endcode This method avoids dividing by zero, so that the non-existence of a solution doesn't by itself mean that you'll get \c inf or \c nan values." \ "note_try_to_help_rvo=This function returns the result by value. In order to make that efficient, it is implemented as just a return statement using a special constructor, hopefully allowing the compiler to perform a RVO (return value optimization)." \ - "nonstableyet=\warning This is not considered to be part of the stable public API yet. Changes may happen in future releases. See \ref Experimental \"Experimental parts of Eigen\" + "nonstableyet=\warning This is not considered to be part of the stable public API yet. Changes may happen in future releases. See \ref Experimental \"Experimental parts of Eigen\"" \ + "implsparsesolverconcept=This class follows the \link TutorialSparseSolverConcept sparse solver concept \endlink." \ + "blank= " \ + "cpp11=<span class='cpp11'>[c++11]</span>" \ + "cpp14=<span class='cpp14'>[c++14]</span>" \ + "cpp17=<span class='cpp17'>[c++17]</span>" \ + "newin{1}=<span class='newin3x'>New in %Eigen \1.</span>" ALIASES += "eigenAutoToc= " @@ -272,7 +279,7 @@ OPTIMIZE_OUTPUT_VHDL = NO # (default is Fortran), use: inc=Fortran f=C. Note that for custom extensions # you also need to set FILE_PATTERNS otherwise the files are not read by doxygen. -EXTENSION_MAPPING = +EXTENSION_MAPPING = .h=C++ no_extension=C++ # If MARKDOWN_SUPPORT is enabled (the default) then doxygen pre-processes all # comments according to the Markdown format, which allows for more readable @@ -317,7 +324,7 @@ IDL_PROPERTY_SUPPORT = YES # member in the group (if any) for the other members of the group. By default # all members of a group must be documented explicitly. -DISTRIBUTE_GROUP_DOC = NO +DISTRIBUTE_GROUP_DOC = YES # Set the SUBGROUPING tag to YES (the default) to allow class member groups of # the same type (for instance a group of public functions) to be put as a @@ -367,7 +374,7 @@ TYPEDEF_HIDES_STRUCT = NO # 2^(16+SYMBOL_CACHE_SIZE). The valid range is 0..9, the default is 0, # corresponding to a cache size of 2^16 = 65536 symbols. -SYMBOL_CACHE_SIZE = 0 +# SYMBOL_CACHE_SIZE = 0 # Similar to the SYMBOL_CACHE_SIZE the size of the symbol lookup cache can be # set using LOOKUP_CACHE_SIZE. This cache is used to resolve symbols given @@ -403,7 +410,7 @@ EXTRACT_PACKAGE = NO # If the EXTRACT_STATIC tag is set to YES all static members of a file # will be included in the documentation. -EXTRACT_STATIC = NO +EXTRACT_STATIC = YES # If the EXTRACT_LOCAL_CLASSES tag is set to YES classes (and structs) # defined locally in source files will be included in the documentation. @@ -460,7 +467,7 @@ HIDE_IN_BODY_DOCS = NO # to NO (the default) then the documentation will be excluded. # Set it to YES to include the internal documentation. -INTERNAL_DOCS = NO +INTERNAL_DOCS = ${EIGEN_DOXY_INTERNAL} # If the CASE_SENSE_NAMES tag is set to NO then Doxygen will only generate # file names in lower-case letters. If set to YES upper-case letters are also @@ -474,13 +481,13 @@ CASE_SENSE_NAMES = YES # will show members with their full class and namespace scopes in the # documentation. If set to YES the scope will be hidden. -HIDE_SCOPE_NAMES = YES +HIDE_SCOPE_NAMES = NO # If the SHOW_INCLUDE_FILES tag is set to YES (the default) then Doxygen # will put a list of the files that are included by a file in the documentation # of that file. -SHOW_INCLUDE_FILES = NO +SHOW_INCLUDE_FILES = ${EIGEN_DOXY_INTERNAL} # If the FORCE_LOCAL_INCLUDES tag is set to YES then Doxygen # will list include files with double quotes in the documentation @@ -546,7 +553,7 @@ STRICT_PROTO_MATCHING = NO # disable (NO) the todo list. This list is created by putting \todo # commands in the documentation. -GENERATE_TODOLIST = NO +GENERATE_TODOLIST = ${EIGEN_DOXY_INTERNAL} # The GENERATE_TESTLIST tag can be used to enable (YES) or # disable (NO) the test list. This list is created by putting \test @@ -558,13 +565,13 @@ GENERATE_TESTLIST = NO # disable (NO) the bug list. This list is created by putting \bug # commands in the documentation. -GENERATE_BUGLIST = NO +GENERATE_BUGLIST = ${EIGEN_DOXY_INTERNAL} # The GENERATE_DEPRECATEDLIST tag can be used to enable (YES) or # disable (NO) the deprecated list. This list is created by putting # \deprecated commands in the documentation. -GENERATE_DEPRECATEDLIST= NO +GENERATE_DEPRECATEDLIST= YES # The ENABLED_SECTIONS tag can be used to enable conditional # documentation sections, marked by \if sectionname ... \endif. @@ -721,7 +728,8 @@ RECURSIVE = YES # Note that relative paths are relative to the directory from which doxygen is # run. -EXCLUDE = "${Eigen_SOURCE_DIR}/Eigen/Eigen2Support" \ +EXCLUDE = "${Eigen_SOURCE_DIR}/Eigen/src/Core/products" \ + "${Eigen_SOURCE_DIR}/Eigen/Eigen2Support" \ "${Eigen_SOURCE_DIR}/Eigen/src/Eigen2Support" \ "${Eigen_SOURCE_DIR}/doc/examples" \ "${Eigen_SOURCE_DIR}/doc/special_examples" \ @@ -802,7 +810,7 @@ EXAMPLE_RECURSIVE = NO # directories that contain image that are included in the documentation (see # the \image command). -IMAGE_PATH = +IMAGE_PATH = ${Eigen_BINARY_DIR}/doc/html # The INPUT_FILTER tag can be used to specify a program that doxygen should # invoke to filter for each input file. Doxygen will invoke the filter program @@ -866,13 +874,13 @@ STRIP_CODE_COMMENTS = YES # then for each documented function all documented # functions referencing it will be listed. -REFERENCED_BY_RELATION = YES +REFERENCED_BY_RELATION = NO # If the REFERENCES_RELATION tag is set to YES # then for each documented function all documented entities # called/used by that function will be listed. -REFERENCES_RELATION = YES +REFERENCES_RELATION = NO # If the REFERENCES_LINK_SOURCE tag is set to YES (the default) # and SOURCE_BROWSER tag is set to YES, then the hyperlinks from @@ -1467,13 +1475,13 @@ XML_OUTPUT = xml # which can be used by a validating XML parser to check the # syntax of the XML files. -XML_SCHEMA = +# XML_SCHEMA = # The XML_DTD tag can be used to specify an XML DTD, # which can be used by a validating XML parser to check the # syntax of the XML files. -XML_DTD = +# XML_DTD = # If the XML_PROGRAMLISTING tag is set to YES Doxygen will # dump the program listings (including syntax highlighting @@ -1583,9 +1591,14 @@ PREDEFINED = EIGEN_EMPTY_STRUCT \ EIGEN_VECTORIZE \ EIGEN_QT_SUPPORT \ EIGEN_STRONG_INLINE=inline \ - "EIGEN2_SUPPORT_STAGE=99" \ + EIGEN_DEVICE_FUNC= \ "EIGEN_MAKE_CWISE_BINARY_OP(METHOD,FUNCTOR)=template<typename OtherDerived> const CwiseBinaryOp<FUNCTOR<Scalar>, const Derived, const OtherDerived> METHOD(const EIGEN_CURRENT_STORAGE_BASE_CLASS<OtherDerived> &other) const;" \ - "EIGEN_CWISE_PRODUCT_RETURN_TYPE(LHS,RHS)=CwiseBinaryOp<internal::scalar_product_op<typename LHS::Scalar, typename RHS::Scalar >, const LHS, const RHS>" + "EIGEN_CWISE_PRODUCT_RETURN_TYPE(LHS,RHS)=CwiseBinaryOp<internal::scalar_product_op<LHS::Scalar,RHS::Scalar>, const LHS, const RHS>"\ + "EIGEN_CAT2(a,b)= a ## b"\ + "EIGEN_CAT(a,b)=EIGEN_CAT2(a,b)"\ + "EIGEN_CWISE_BINARY_RETURN_TYPE(LHS,RHS,OPNAME)=CwiseBinaryOp<EIGEN_CAT(EIGEN_CAT(internal::scalar_,OPNAME),_op)<LHS::Scalar, RHS::Scalar>, const LHS, const RHS>"\ + DOXCOMMA=, + # If the MACRO_EXPANSION and EXPAND_ONLY_PREDEF tags are set to YES then # this tag can be used to specify a list of macro names that should be expanded. @@ -1601,7 +1614,15 @@ EXPAND_AS_DEFINED = EIGEN_MAKE_TYPEDEFS \ EIGEN_CURRENT_STORAGE_BASE_CLASS \ EIGEN_MATHFUNC_IMPL \ _EIGEN_GENERIC_PUBLIC_INTERFACE \ - EIGEN2_SUPPORT + EIGEN_ARRAY_DECLARE_GLOBAL_UNARY \ + EIGEN_EMPTY \ + EIGEN_EULER_ANGLES_TYPEDEFS \ + EIGEN_EULER_ANGLES_SINGLE_TYPEDEF \ + EIGEN_EULER_SYSTEM_TYPEDEF \ + EIGEN_DOC_UNARY_ADDONS \ + EIGEN_DOC_BLOCK_ADDONS_NOT_INNER_PANEL \ + EIGEN_DOC_BLOCK_ADDONS_INNER_PANEL_IF + # If the SKIP_FUNCTION_MACROS tag is set to YES (the default) then # doxygen's preprocessor will remove all references to function-like macros @@ -1701,7 +1722,7 @@ DOT_NUM_THREADS = 0 # the DOTFONTPATH environment variable or by setting DOT_FONTPATH to the # directory containing the font. -DOT_FONTNAME = FreeSans +DOT_FONTNAME = # The DOT_FONTSIZE tag can be used to set the size of the font of dot graphs. # The default size is 10pt. diff --git a/doc/FixedSizeVectorizable.dox b/doc/FixedSizeVectorizable.dox index 8ae135173..0012465ca 100644 --- a/doc/FixedSizeVectorizable.dox +++ b/doc/FixedSizeVectorizable.dox @@ -1,10 +1,10 @@ namespace Eigen { -/** \eigenManualPage TopicFixedSizeVectorizable Fixed-size vectorizable Eigen objects +/** \eigenManualPage TopicFixedSizeVectorizable Fixed-size vectorizable %Eigen objects The goal of this page is to explain what we mean by "fixed-size vectorizable". -\section summary Executive Summary +\section FixedSizeVectorizable_summary Executive Summary An Eigen object is called "fixed-size vectorizable" if it has fixed size and that size is a multiple of 16 bytes. @@ -21,17 +21,17 @@ Examples include: \li Eigen::Quaterniond \li Eigen::Quaternionf -\section explanation Explanation +\section FixedSizeVectorizable_explanation Explanation -First, "fixed-size" should be clear: an Eigen object has fixed size if its number of rows and its number of columns are fixed at compile-time. So for example Matrix3f has fixed size, but MatrixXf doesn't (the opposite of fixed-size is dynamic-size). +First, "fixed-size" should be clear: an %Eigen object has fixed size if its number of rows and its number of columns are fixed at compile-time. So for example \ref Matrix3f has fixed size, but \ref MatrixXf doesn't (the opposite of fixed-size is dynamic-size). -The array of coefficients of a fixed-size Eigen object is a plain "static array", it is not dynamically allocated. For example, the data behind a Matrix4f is just a "float array[16]". +The array of coefficients of a fixed-size %Eigen object is a plain "static array", it is not dynamically allocated. For example, the data behind a \ref Matrix4f is just a "float array[16]". Fixed-size objects are typically very small, which means that we want to handle them with zero runtime overhead -- both in terms of memory usage and of speed. -Now, vectorization (both SSE and AltiVec) works with 128-bit packets. Moreover, for performance reasons, these packets need to be have 128-bit alignment. +Now, vectorization works with 128-bit packets (e.g., SSE, AltiVec, NEON), 256-bit packets (e.g., AVX), or 512-bit packets (e.g., AVX512). Moreover, for performance reasons, these packets are most efficiently read and written if they have the same alignment as the packet size, that is 16 bytes, 32 bytes, and 64 bytes respectively. -So it turns out that the only way that fixed-size Eigen objects can be vectorized, is if their size is a multiple of 128 bits, or 16 bytes. Eigen will then request 16-byte alignment for these objects, and henceforth rely on these objects being aligned so no runtime check for alignment is performed. +So it turns out that the best way that fixed-size %Eigen objects can be vectorized, is if their size is a multiple of 16 bytes (or more). %Eigen will then request 16-byte alignment (or more) for these objects, and henceforth rely on these objects being aligned to achieve maximal efficiency. */ diff --git a/doc/InplaceDecomposition.dox b/doc/InplaceDecomposition.dox new file mode 100644 index 000000000..cb1c6d413 --- /dev/null +++ b/doc/InplaceDecomposition.dox @@ -0,0 +1,115 @@ +namespace Eigen { + +/** \eigenManualPage InplaceDecomposition Inplace matrix decompositions + +Starting from %Eigen 3.3, the LU, Cholesky, and QR decompositions can operate \em inplace, that is, directly within the given input matrix. +This feature is especially useful when dealing with huge matrices, and or when the available memory is very limited (embedded systems). + +To this end, the respective decomposition class must be instantiated with a Ref<> matrix type, and the decomposition object must be constructed with the input matrix as argument. As an example, let us consider an inplace LU decomposition with partial pivoting. + +Let's start with the basic inclusions, and declaration of a 2x2 matrix \c A: + +<table class="example"> +<tr><th>code</th><th>output</th></tr> +<tr> + <td>\snippet TutorialInplaceLU.cpp init + </td> + <td>\snippet TutorialInplaceLU.out init + </td> +</tr> +</table> + +No surprise here! Then, let's declare our inplace LU object \c lu, and check the content of the matrix \c A: + +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp declaration + </td> + <td>\snippet TutorialInplaceLU.out declaration + </td> +</tr> +</table> + +Here, the \c lu object computes and stores the \c L and \c U factors within the memory held by the matrix \c A. +The coefficients of \c A have thus been destroyed during the factorization, and replaced by the L and U factors as one can verify: + +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp matrixLU + </td> + <td>\snippet TutorialInplaceLU.out matrixLU + </td> +</tr> +</table> + +Then, one can use the \c lu object as usual, for instance to solve the Ax=b problem: +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp solve + </td> + <td>\snippet TutorialInplaceLU.out solve + </td> +</tr> +</table> + +Here, since the content of the original matrix \c A has been lost, we had to declared a new matrix \c A0 to verify the result. + +Since the memory is shared between \c A and \c lu, modifying the matrix \c A will make \c lu invalid. +This can easily be verified by modifying the content of \c A and trying to solve the initial problem again: + +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp modifyA + </td> + <td>\snippet TutorialInplaceLU.out modifyA + </td> +</tr> +</table> + +Note that there is no shared pointer under the hood, it is the \b responsibility \b of \b the \b user to keep the input matrix \c A in life as long as \c lu is living. + +If one wants to update the factorization with the modified A, one has to call the compute method as usual: +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp recompute + </td> + <td>\snippet TutorialInplaceLU.out recompute + </td> +</tr> +</table> + +Note that calling compute does not change the memory which is referenced by the \c lu object. Therefore, if the compute method is called with another matrix \c A1 different than \c A, then the content of \c A1 won't be modified. This is still the content of \c A that will be used to store the L and U factors of the matrix \c A1. +This can easily be verified as follows: +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp recompute_bis0 + </td> + <td>\snippet TutorialInplaceLU.out recompute_bis0 + </td> +</tr> +</table> +The matrix \c A1 is unchanged, and one can thus solve A1*x=b, and directly check the residual without any copy of \c A1: +<table class="example"> +<tr> + <td>\snippet TutorialInplaceLU.cpp recompute_bis1 + </td> + <td>\snippet TutorialInplaceLU.out recompute_bis1 + </td> +</tr> +</table> + + +Here is the list of matrix decompositions supporting this inplace mechanism: + +- class LLT +- class LDLT +- class PartialPivLU +- class FullPivLU +- class HouseholderQR +- class ColPivHouseholderQR +- class FullPivHouseholderQR +- class CompleteOrthogonalDecomposition + +*/ + +}
\ No newline at end of file diff --git a/doc/LeastSquares.dox b/doc/LeastSquares.dox new file mode 100644 index 000000000..e2191a22f --- /dev/null +++ b/doc/LeastSquares.dox @@ -0,0 +1,70 @@ +namespace Eigen { + +/** \eigenManualPage LeastSquares Solving linear least squares systems + +This page describes how to solve linear least squares systems using %Eigen. An overdetermined system +of equations, say \a Ax = \a b, has no solutions. In this case, it makes sense to search for the +vector \a x which is closest to being a solution, in the sense that the difference \a Ax - \a b is +as small as possible. This \a x is called the least square solution (if the Euclidean norm is used). + +The three methods discussed on this page are the SVD decomposition, the QR decomposition and normal +equations. Of these, the SVD decomposition is generally the most accurate but the slowest, normal +equations is the fastest but least accurate, and the QR decomposition is in between. + +\eigenAutoToc + + +\section LeastSquaresSVD Using the SVD decomposition + +The \link JacobiSVD::solve() solve() \endlink method in the JacobiSVD class can be directly used to +solve linear squares systems. It is not enough to compute only the singular values (the default for +this class); you also need the singular vectors but the thin SVD decomposition suffices for +computing least squares solutions: + +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr> + <td>\include TutorialLinAlgSVDSolve.cpp </td> + <td>\verbinclude TutorialLinAlgSVDSolve.out </td> +</tr> +</table> + +This is example from the page \link TutorialLinearAlgebra Linear algebra and decompositions \endlink. + + +\section LeastSquaresQR Using the QR decomposition + +The solve() method in QR decomposition classes also computes the least squares solution. There are +three QR decomposition classes: HouseholderQR (no pivoting, so fast but unstable), +ColPivHouseholderQR (column pivoting, thus a bit slower but more accurate) and FullPivHouseholderQR +(full pivoting, so slowest and most stable). Here is an example with column pivoting: + +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr> + <td>\include LeastSquaresQR.cpp </td> + <td>\verbinclude LeastSquaresQR.out </td> +</tr> +</table> + + +\section LeastSquaresNormalEquations Using normal equations + +Finding the least squares solution of \a Ax = \a b is equivalent to solving the normal equation +<i>A</i><sup>T</sup><i>Ax</i> = <i>A</i><sup>T</sup><i>b</i>. This leads to the following code + +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr> + <td>\include LeastSquaresNormalEquations.cpp </td> + <td>\verbinclude LeastSquaresNormalEquations.out </td> +</tr> +</table> + +If the matrix \a A is ill-conditioned, then this is not a good method, because the condition number +of <i>A</i><sup>T</sup><i>A</i> is the square of the condition number of \a A. This means that you +lose twice as many digits using normal equation than if you use the other methods. + +*/ + +}
\ No newline at end of file diff --git a/doc/Manual.dox b/doc/Manual.dox index 3367982ca..342b145fd 100644 --- a/doc/Manual.dox +++ b/doc/Manual.dox @@ -3,18 +3,32 @@ namespace Eigen { +/** \page UserManual_CustomizingEigen Extending/Customizing Eigen + %Eigen can be extended in several ways, for instance, by defining global methods, by inserting custom methods within main %Eigen's classes through the \ref TopicCustomizing_Plugins "plugin" mechanism, by adding support to \ref TopicCustomizing_CustomScalar "custom scalar types" etc. See below for the respective sub-topics. + - \subpage TopicCustomizing_Plugins + - \subpage TopicCustomizing_InheritingMatrix + - \subpage TopicCustomizing_CustomScalar + - \subpage TopicCustomizing_NullaryExpr + - \subpage TopicNewExpressionType + \sa \ref TopicPreprocessorDirectives +*/ + + /** \page UserManual_Generalities General topics - \subpage Eigen2ToEigen3 - \subpage TopicFunctionTakingEigenTypes - \subpage TopicPreprocessorDirectives - \subpage TopicAssertions - - \subpage TopicCustomizingEigen - \subpage TopicMultiThreading + - \subpage TopicUsingBlasLapack - \subpage TopicUsingIntelMKL + - \subpage TopicCUDA + - \subpage TopicPitfalls - \subpage TopicTemplateKeyword - \subpage UserManual_UnderstandingEigen + - \subpage TopicCMakeGuide */ - + /** \page UserManual_UnderstandingEigen Understanding Eigen - \subpage TopicInsideEigenExample - \subpage TopicClassHierarchy @@ -56,6 +70,8 @@ namespace Eigen { \ingroup DenseMatrixManipulation_chapter */ /** \addtogroup TutorialMapClass \ingroup DenseMatrixManipulation_chapter */ +/** \addtogroup TutorialReshapeSlicing + \ingroup DenseMatrixManipulation_chapter */ /** \addtogroup TopicAliasing \ingroup DenseMatrixManipulation_chapter */ /** \addtogroup TopicStorageOrders @@ -85,6 +101,9 @@ namespace Eigen { /** \addtogroup Householder_Module \ingroup DenseMatrixManipulation_Reference */ +/** \addtogroup CoeffwiseMathFunctions + \ingroup DenseMatrixManipulation_chapter */ + /** \addtogroup QuickRefPage \ingroup DenseMatrixManipulation_chapter */ @@ -96,6 +115,12 @@ namespace Eigen { \ingroup DenseLinearSolvers_chapter */ /** \addtogroup TopicLinearAlgebraDecompositions \ingroup DenseLinearSolvers_chapter */ +/** \addtogroup LeastSquares + \ingroup DenseLinearSolvers_chapter */ +/** \addtogroup InplaceDecomposition + \ingroup DenseLinearSolvers_chapter */ +/** \addtogroup DenseDecompositionBenchmark + \ingroup DenseLinearSolvers_chapter */ /** \addtogroup DenseLinearSolvers_Reference \ingroup DenseLinearSolvers_chapter */ @@ -120,6 +145,8 @@ namespace Eigen { \ingroup Sparse_chapter */ /** \addtogroup TopicSparseSystems \ingroup Sparse_chapter */ +/** \addtogroup MatrixfreeSolverExample + \ingroup Sparse_chapter */ /** \addtogroup Sparse_Reference \ingroup Sparse_chapter */ @@ -156,4 +183,7 @@ namespace Eigen { \ingroup Geometry_Reference */ /** \addtogroup Splines_Module \ingroup Geometry_Reference */ + +/** \internal \brief Namespace containing low-level routines from the %Eigen library. */ +namespace internal {} } diff --git a/doc/MatrixfreeSolverExample.dox b/doc/MatrixfreeSolverExample.dox new file mode 100644 index 000000000..3efa292b5 --- /dev/null +++ b/doc/MatrixfreeSolverExample.dox @@ -0,0 +1,20 @@ + +namespace Eigen { + +/** + +\eigenManualPage MatrixfreeSolverExample Matrix-free solvers + +Iterative solvers such as ConjugateGradient and BiCGSTAB can be used in a matrix free context. To this end, user must provide a wrapper class inheriting EigenBase<> and implementing the following methods: + - \c Index \c rows() and \c Index \c cols(): returns number of rows and columns respectively + - \c operator* with your type and an %Eigen dense column vector (its actual implementation goes in a specialization of the internal::generic_product_impl class) + +\c Eigen::internal::traits<> must also be specialized for the wrapper type. + +Here is a complete example wrapping an Eigen::SparseMatrix: +\include matrixfree_cg.cpp +Output: \verbinclude matrixfree_cg.out + +*/ + +}
\ No newline at end of file diff --git a/doc/NewExpressionType.dox b/doc/NewExpressionType.dox new file mode 100644 index 000000000..c2f243312 --- /dev/null +++ b/doc/NewExpressionType.dox @@ -0,0 +1,143 @@ +namespace Eigen { + +/** \page TopicNewExpressionType Adding a new expression type + +<!--<span style="font-size:130%; color:red; font-weight: 900;"></span>--> +\warning +Disclaimer: this page is tailored to very advanced users who are not afraid of dealing with some %Eigen's internal aspects. +In most cases, a custom expression can be avoided by either using custom \ref MatrixBase::unaryExpr "unary" or \ref MatrixBase::binaryExpr "binary" functors, +while extremely complex matrix manipulations can be achieved by a nullary functors as described in the \ref TopicCustomizing_NullaryExpr "previous page". + +This page describes with the help of an example how to implement a new +light-weight expression type in %Eigen. This consists of three parts: +the expression type itself, a traits class containing compile-time +information about the expression, and the evaluator class which is +used to evaluate the expression to a matrix. + +\b TO \b DO: Write a page explaining the design, with details on +vectorization etc., and refer to that page here. + + +\eigenAutoToc + +\section TopicSetting The setting + +A circulant matrix is a matrix where each column is the same as the +column to the left, except that it is cyclically shifted downwards. +For example, here is a 4-by-4 circulant matrix: +\f[ \begin{bmatrix} + 1 & 8 & 4 & 2 \\ + 2 & 1 & 8 & 4 \\ + 4 & 2 & 1 & 8 \\ + 8 & 4 & 2 & 1 +\end{bmatrix} \f] +A circulant matrix is uniquely determined by its first column. We wish +to write a function \c makeCirculant which, given the first column, +returns an expression representing the circulant matrix. + +For simplicity, we restrict the \c makeCirculant function to dense +matrices. It may make sense to also allow arrays, or sparse matrices, +but we will not do so here. We also do not want to support +vectorization. + + +\section TopicPreamble Getting started + +We will present the file implementing the \c makeCirculant function +part by part. We start by including the appropriate header files and +forward declaring the expression class, which we will call +\c Circulant. The \c makeCirculant function will return an object of +this type. The class \c Circulant is in fact a class template; the +template argument \c ArgType refers to the type of the vector passed +to the \c makeCirculant function. + +\include make_circulant.cpp.preamble + + +\section TopicTraits The traits class + +For every expression class \c X, there should be a traits class +\c Traits<X> in the \c Eigen::internal namespace containing +information about \c X known as compile time. + +As explained in \ref TopicSetting, we designed the \c Circulant +expression class to refer to dense matrices. The entries of the +circulant matrix have the same type as the entries of the vector +passed to the \c makeCirculant function. The type used to index the +entries is also the same. Again for simplicity, we will only return +column-major matrices. Finally, the circulant matrix is a square +matrix (number of rows equals number of columns), and the number of +rows equals the number of rows of the column vector passed to the +\c makeCirculant function. If this is a dynamic-size vector, then the +size of the circulant matrix is not known at compile-time. + +This leads to the following code: + +\include make_circulant.cpp.traits + + +\section TopicExpression The expression class + +The next step is to define the expression class itself. In our case, +we want to inherit from \c MatrixBase in order to expose the interface +for dense matrices. In the constructor, we check that we are passed a +column vector (see \ref TopicAssertions) and we store the vector from +which we are going to build the circulant matrix in the member +variable \c m_arg. Finally, the expression class should compute the +size of the corresponding circulant matrix. As explained above, this +is a square matrix with as many columns as the vector used to +construct the matrix. + +\b TO \b DO: What about the \c Nested typedef? It seems to be +necessary; is this only temporary? + +\include make_circulant.cpp.expression + + +\section TopicEvaluator The evaluator + +The last big fragment implements the evaluator for the \c Circulant +expression. The evaluator computes the entries of the circulant +matrix; this is done in the \c .coeff() member function. The entries +are computed by finding the corresponding entry of the vector from +which the circulant matrix is constructed. Getting this entry may +actually be non-trivial when the circulant matrix is constructed from +a vector which is given by a complicated expression, so we use the +evaluator which corresponds to the vector. + +The \c CoeffReadCost constant records the cost of computing an entry +of the circulant matrix; we ignore the index computation and say that +this is the same as the cost of computing an entry of the vector from +which the circulant matrix is constructed. + +In the constructor, we save the evaluator for the column vector which +defined the circulant matrix. We also save the size of that vector; +remember that we can query an expression object to find the size but +not the evaluator. + +\include make_circulant.cpp.evaluator + + +\section TopicEntry The entry point + +After all this, the \c makeCirculant function is very simple. It +simply creates an expression object and returns it. + +\include make_circulant.cpp.entry + + +\section TopicMain A simple main function for testing + +Finally, a short \c main function that shows how the \c makeCirculant +function can be called. + +\include make_circulant.cpp.main + +If all the fragments are combined, the following output is produced, +showing that the program works as expected: + +\include make_circulant.out + +*/ +} + diff --git a/doc/Overview.dox b/doc/Overview.dox index 9ab96233a..dbb49bd21 100644 --- a/doc/Overview.dox +++ b/doc/Overview.dox @@ -17,7 +17,9 @@ You're a MatLab user? There is also a <a href="AsciiQuickReference.txt">short AS The \b main \b documentation is organized into \em chapters covering different domains of features. They are themselves composed of \em user \em manual pages describing the different features in a comprehensive way, and \em reference pages that gives you access to the API documentation through the related Eigen's \em modules and \em classes. -Under the \subpage UserManual_Generalities section, you will find documentation on more general topics such as preprocessor directives, controlling assertions, multi-threading, MKL support, some Eigen's internal insights, and much more... +Under the \subpage UserManual_CustomizingEigen section, you will find discussions and examples on extending %Eigen's features and supporting custom scalar types. + +Under the \subpage UserManual_Generalities section, you will find documentation on more general topics such as preprocessor directives, controlling assertions, multi-threading, MKL support, some Eigen's internal insights, and much more... Finally, do not miss the search engine, useful to quickly get to the documentation of a given class or function. diff --git a/doc/PassingByValue.dox b/doc/PassingByValue.dox index bf4d0ef4b..9254fe6d8 100644 --- a/doc/PassingByValue.dox +++ b/doc/PassingByValue.dox @@ -4,21 +4,21 @@ namespace Eigen { Passing objects by value is almost always a very bad idea in C++, as this means useless copies, and one should pass them by reference instead. -With Eigen, this is even more important: passing \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen objects" by value is not only inefficient, it can be illegal or make your program crash! And the reason is that these Eigen objects have alignment modifiers that aren't respected when they are passed by value. +With %Eigen, this is even more important: passing \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen objects" by value is not only inefficient, it can be illegal or make your program crash! And the reason is that these %Eigen objects have alignment modifiers that aren't respected when they are passed by value. -So for example, a function like this, where v is passed by value: +For example, a function like this, where \c v is passed by value: \code void my_function(Eigen::Vector2d v); \endcode -needs to be rewritten as follows, passing v by reference: +needs to be rewritten as follows, passing \c v by const reference: \code void my_function(const Eigen::Vector2d& v); \endcode -Likewise if you have a class having a Eigen object as member: +Likewise if you have a class having an %Eigen object as member: \code struct Foo diff --git a/doc/Pitfalls.dox b/doc/Pitfalls.dox new file mode 100644 index 000000000..cf42effef --- /dev/null +++ b/doc/Pitfalls.dox @@ -0,0 +1,38 @@ +namespace Eigen { + +/** \page TopicPitfalls Common pitfalls + +\section TopicPitfalls_template_keyword Compilation error with template methods + +See this \link TopicTemplateKeyword page \endlink. + +\section TopicPitfalls_auto_keyword C++11 and the auto keyword + +In short: do not use the auto keywords with Eigen's expressions, unless you are 100% sure about what you are doing. In particular, do not use the auto keyword as a replacement for a Matrix<> type. Here is an example: + +\code +MatrixXd A, B; +auto C = A*B; +for(...) { ... w = C * v; ...} +\endcode + +In this example, the type of C is not a MatrixXd but an abstract expression representing a matrix product and storing references to A and B. Therefore, the product of A*B will be carried out multiple times, once per iteration of the for loop. Moreover, if the coefficients of A or B change during the iteration, then C will evaluate to different values. + +Here is another example leading to a segfault: +\code +auto C = ((A+B).eval()).transpose(); +// do something with C +\endcode +The problem is that eval() returns a temporary object (in this case a MatrixXd) which is then referenced by the Transpose<> expression. However, this temporary is deleted right after the first line, and there the C expression reference a dead object. The same issue might occur when sub expressions are automatically evaluated by Eigen as in the following example: +\code +VectorXd u, v; +auto C = u + (A*v).normalized(); +// do something with C +\endcode +where the normalized() method has to evaluate the expensive product A*v to avoid evaluating it twice. On the other hand, the following example is perfectly fine: +\code +auto C = (u + (A*v).normalized()).eval(); +\endcode +In this case, C will be a regular VectorXd object. +*/ +} diff --git a/doc/PreprocessorDirectives.dox b/doc/PreprocessorDirectives.dox index 6e40f919f..f01b39aec 100644 --- a/doc/PreprocessorDirectives.dox +++ b/doc/PreprocessorDirectives.dox @@ -5,7 +5,7 @@ namespace Eigen { You can control some aspects of %Eigen by defining the preprocessor tokens using \c \#define. These macros should be defined before any %Eigen headers are included. Often they are best set in the project options. -This page lists the preprocesor tokens recognised by %Eigen. +This page lists the preprocessor tokens recognized by %Eigen. \eigenAutoToc @@ -18,25 +18,67 @@ one option, and other parts (or libraries that you use) are compiled with anothe fail to link or exhibit subtle bugs. Nevertheless, these options can be useful for people who know what they are doing. - - \b EIGEN2_SUPPORT - if defined, enables the Eigen2 compatibility mode. This is meant to ease the transition - of Eigen2 to Eigen3 (see \ref Eigen2ToEigen3). Not defined by default. - - \b EIGEN2_SUPPORT_STAGEnn_xxx (for various values of nn and xxx) - staged migration path from Eigen2 to - Eigen3; see \ref Eigen2SupportModes. + - \b EIGEN2_SUPPORT and \b EIGEN2_SUPPORT_STAGEnn_xxx are disabled starting from the 3.3 release. + Defining one of these will raise a compile-error. If you need to compile Eigen2 code, + <a href="http://eigen.tuxfamily.org/index.php?title=Eigen2">check this site</a>. - \b EIGEN_DEFAULT_DENSE_INDEX_TYPE - the type for column and row indices in matrices, vectors and array (DenseBase::Index). Set to \c std::ptrdiff_t by default. - \b EIGEN_DEFAULT_IO_FORMAT - the IOFormat to use when printing a matrix if no %IOFormat is specified. Defaults to the %IOFormat constructed by the default constructor IOFormat::IOFormat(). - \b EIGEN_INITIALIZE_MATRICES_BY_ZERO - if defined, all entries of newly constructed matrices and arrays are initialized to zero, as are new entries in matrices and arrays after resizing. Not defined by default. + \warning The unary (resp. binary) constructor of \c 1x1 (resp. \c 2x1 or \c 1x2) fixed size matrices is + always interpreted as an initialization constructor where the argument(s) are the coefficient values + and not the sizes. For instance, \code Vector2d v(2,1); \endcode will create a vector with coeficients [2,1], + and \b not a \c 2x1 vector initialized with zeros (i.e., [0,0]). If such cases might occur, then it is + recommended to use the default constructor with a explicit call to resize: + \code + Matrix<?,SizeAtCompileTime,1> v; + v.resize(size); + Matrix<?,RowsAtCompileTime,ColsAtCompileTime> m; + m.resize(rows,cols); + \endcode - \b EIGEN_INITIALIZE_MATRICES_BY_NAN - if defined, all entries of newly constructed matrices and arrays are initialized to NaN, as are new entries in matrices and arrays after resizing. This option is especially useful for debugging purpose, though a memory tool like <a href="http://valgrind.org/">valgrind</a> is preferable. Not defined by default. + \warning See the documentation of \c EIGEN_INITIALIZE_MATRICES_BY_ZERO for a discussion on a limitations + of these macros when applied to \c 1x1, \c 1x2, and \c 2x1 fixed-size matrices. - \b EIGEN_NO_AUTOMATIC_RESIZING - if defined, the matrices (or arrays) on both sides of an assignment <tt>a = b</tt> have to be of the same size; otherwise, %Eigen automatically resizes \c a so that it is of the correct size. Not defined by default. +\section TopicPreprocessorDirectivesCppVersion C++ standard features + +By default, %Eigen strive to automatically detect and enable langage features at compile-time based on +the information provided by the compiler. + + - \b EIGEN_MAX_CPP_VER - disables usage of C++ features requiring a version greater than EIGEN_MAX_CPP_VER. + Possible values are: 03, 11, 14, 17, etc. If not defined (the default), %Eigen enables all features supported + by the compiler. + +Individual features can be explicitly enabled or disabled by defining the following token to 0 or 1 respectively. +For instance, one might limit the C++ version to C++03 by defining EIGEN_MAX_CPP_VER=03, but still enable C99 math +functions by defining EIGEN_HAS_C99_MATH=1. + + - \b EIGEN_HAS_C99_MATH - controls the usage of C99 math functions such as erf, erfc, lgamma, etc. + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_CXX11_MATH - controls the implementation of some functions such as round, logp1, isinf, isnan, etc. + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_RVALUE_REFERENCES - defines whetehr rvalue references are supported + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_STD_RESULT_OF - defines whether std::result_of is supported + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_VARIADIC_TEMPLATES - defines whether variadic templates are supported + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_CONSTEXPR - defines whether relaxed const expression are supported + Automatic detection disabled if EIGEN_MAX_CPP_VER<14. + - \b EIGEN_HAS_CXX11_CONTAINERS - defines whether STL's containers follows C++11 specifications + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + - \b EIGEN_HAS_CXX11_NOEXCEPT - defines whether noexcept is supported + Automatic detection disabled if EIGEN_MAX_CPP_VER<11. + \section TopicPreprocessorDirectivesAssertions Assertions The %Eigen library contains many assertions to guard against programming errors, both at compile time and at @@ -55,30 +97,39 @@ run time. However, these assertions do cost time and can thus be turned off. \section TopicPreprocessorDirectivesPerformance Alignment, vectorization and performance tweaking - - \b EIGEN_MALLOC_ALREADY_ALIGNED - Can be set to 0 or 1 to tell whether default system malloc already + - \b \c EIGEN_MALLOC_ALREADY_ALIGNED - Can be set to 0 or 1 to tell whether default system \c malloc already returns aligned buffers. In not defined, then this information is automatically deduced from the compiler and system preprocessor tokens. - - \b EIGEN_DONT_ALIGN - disables alignment completely. %Eigen will not try to align its objects and does not - expect that any objects passed to it are aligned. This will turn off vectorization. Not defined by default. - - \b EIGEN_DONT_ALIGN_STATICALLY - disables alignment of arrays on the stack. Not defined by default, unless - \c EIGEN_DONT_ALIGN is defined. + - \b \c EIGEN_MAX_ALIGN_BYTES - Must be a power of two, or 0. Defines an upper bound on the memory boundary in bytes on which dynamically and statically allocated data may be aligned by %Eigen. If not defined, a default value is automatically computed based on architecture, compiler, and OS. + This option is typically used to enforce binary compatibility between code/libraries compiled with different SIMD options. For instance, one may compile AVX code and enforce ABI compatibility with existing SSE code by defining \c EIGEN_MAX_ALIGN_BYTES=16. In the other way round, since by default AVX implies 32 bytes alignment for best performance, one can compile SSE code to be ABI compatible with AVX code by defining \c EIGEN_MAX_ALIGN_BYTES=32. + - \b \c EIGEN_MAX_STATIC_ALIGN_BYTES - Same as \c EIGEN_MAX_ALIGN_BYTES but for statically allocated data only. By default, if only \c EIGEN_MAX_ALIGN_BYTES is defined, then \c EIGEN_MAX_STATIC_ALIGN_BYTES == \c EIGEN_MAX_ALIGN_BYTES, otherwise a default value is automatically computed based on architecture, compiler, and OS (can be smaller than the default value of EIGEN_MAX_ALIGN_BYTES on architectures that do not support stack alignment). + Let us emphasize that \c EIGEN_MAX_*_ALIGN_BYTES define only a diserable upper bound. In practice data is aligned to largest power-of-two common divisor of \c EIGEN_MAX_STATIC_ALIGN_BYTES and the size of the data, such that memory is not wasted. + - \b \c EIGEN_DONT_PARALLELIZE - if defined, this disables multi-threading. This is only relevant if you enabled OpenMP. + See \ref TopicMultiThreading for details. - \b EIGEN_DONT_VECTORIZE - disables explicit vectorization when defined. Not defined by default, unless alignment is disabled by %Eigen's platform test or the user defining \c EIGEN_DONT_ALIGN. - - \b EIGEN_FAST_MATH - enables some optimizations which might affect the accuracy of the result. This currently + - \b \c EIGEN_UNALIGNED_VECTORIZE - disables/enables vectorization with unaligned stores. Default is 1 (enabled). + If set to 0 (disabled), then expression for which the destination cannot be aligned are not vectorized (e.g., unaligned + small fixed size vectors or matrices) + - \b \c EIGEN_FAST_MATH - enables some optimizations which might affect the accuracy of the result. This currently enables the SSE vectorization of sin() and cos(), and speedups sqrt() for single precision. Defined to 1 by default. Define it to 0 to disable. - - \b EIGEN_UNROLLING_LIMIT - defines the size of a loop to enable meta unrolling. Set it to zero to disable + - \b \c EIGEN_UNROLLING_LIMIT - defines the size of a loop to enable meta unrolling. Set it to zero to disable unrolling. The size of a loop here is expressed in %Eigen's own notion of "number of FLOPS", it does not correspond to the number of iterations or the number of instructions. The default is value 100. - - \b EIGEN_STACK_ALLOCATION_LIMIT - defines the maximum bytes for a buffer to be allocated on the stack. For internal + - \b \c EIGEN_STACK_ALLOCATION_LIMIT - defines the maximum bytes for a buffer to be allocated on the stack. For internal temporary buffers, dynamic memory allocation is employed as a fall back. For fixed-size matrices or arrays, exceeding this threshold raises a compile time assertion. Use 0 to set no limit. Default is 128 KB. + - \c EIGEN_DONT_ALIGN - Deprecated, it is a synonym for \c EIGEN_MAX_ALIGN_BYTES=0. It disables alignment completely. %Eigen will not try to align its objects and does not expect that any objects passed to it are aligned. This will turn off vectorization if \b EIGEN_UNALIGNED_VECTORIZE=1. Not defined by default. + - \c EIGEN_DONT_ALIGN_STATICALLY - Deprecated, it is a synonym for \c EIGEN_MAX_STATIC_ALIGN_BYTES=0. It disables alignment of arrays on the stack. Not defined by default, unless \c EIGEN_DONT_ALIGN is defined. + + \section TopicPreprocessorDirectivesPlugins Plugins It is possible to add new methods to many fundamental classes in %Eigen by writing a plugin. As explained in -the section \ref ExtendingMatrixBase, the plugin is specified by defining a \c EIGEN_xxx_PLUGIN macro. The +the section \ref TopicCustomizing_Plugins, the plugin is specified by defining a \c EIGEN_xxx_PLUGIN macro. The following macros are supported; none of them are defined by default. - \b EIGEN_ARRAY_PLUGIN - filename of plugin for extending the Array class. @@ -89,6 +140,8 @@ following macros are supported; none of them are defined by default. - \b EIGEN_MATRIX_PLUGIN - filename of plugin for extending the Matrix class. - \b EIGEN_MATRIXBASE_PLUGIN - filename of plugin for extending the MatrixBase class. - \b EIGEN_PLAINOBJECTBASE_PLUGIN - filename of plugin for extending the PlainObjectBase class. + - \b EIGEN_MAPBASE_PLUGIN - filename of plugin for extending the MapBase class. + - \b EIGEN_QUATERNION_PLUGIN - filename of plugin for extending the Quaternion class. - \b EIGEN_QUATERNIONBASE_PLUGIN - filename of plugin for extending the QuaternionBase class. - \b EIGEN_SPARSEMATRIX_PLUGIN - filename of plugin for extending the SparseMatrix class. - \b EIGEN_SPARSEMATRIXBASE_PLUGIN - filename of plugin for extending the SparseMatrixBase class. diff --git a/doc/QuickReference.dox b/doc/QuickReference.dox index a4be0f68a..44f5410db 100644 --- a/doc/QuickReference.dox +++ b/doc/QuickReference.dox @@ -13,17 +13,17 @@ The Eigen library is divided in a Core module and several additional modules. Ea <table class="manual"> <tr><th>Module</th><th>Header file</th><th>Contents</th></tr> -<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> +<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> <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> -<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> -<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> -<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> -<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> -<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> -<tr><td>\link Eigenvalues_Module Eigenvalues \endlink</td><td>\code#include <Eigen/Eigenvalues>\endcode</td><td>Eigenvalue, eigenvector decompositions (EigenSolver, SelfAdjointEigenSolver, ComplexEigenSolver)</td></tr> -<tr class="alt"><td>\link Sparse_modules Sparse \endlink</td><td>\code#include <Eigen/Sparse>\endcode</td><td>%Sparse matrix storage and related basic linear algebra (SparseMatrix, DynamicSparseMatrix, SparseVector)</td></tr> -<tr><td></td><td>\code#include <Eigen/Dense>\endcode</td><td>Includes Core, Geometry, LU, Cholesky, SVD, QR, and Eigenvalues header files</td></tr> -<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> +<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> +<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> +<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> +<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> +<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> +<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> +<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> +<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> +<tr ><td></td><td>\code#include <Eigen/Eigen>\endcode</td><td>Includes %Dense and %Sparse header files (the whole Eigen library)</td></tr> </table> <a href="#" class="top">top</a> @@ -340,7 +340,7 @@ mat1 = mat2.adjoint(); mat1.adjointInPlace(); \endcode </td></tr> <tr><td> -\link MatrixBase::dot() dot \endlink product \n inner product \matrixworld</td><td>\code +\link MatrixBase::dot dot \endlink product \n inner product \matrixworld</td><td>\code scalar = vec1.dot(vec2); scalar = col1.adjoint() * col2; scalar = (col1.adjoint() * col2).value();\endcode @@ -364,32 +364,10 @@ vec3 = vec1.cross(vec2);\endcode</td></tr> <a href="#" class="top">top</a> \section QuickRef_Coeffwise Coefficient-wise \& Array operators -Coefficient-wise operators for matrices and vectors: -<table class="manual"> -<tr><th>Matrix API \matrixworld</th><th>Via Array conversions</th></tr> -<tr><td>\code -mat1.cwiseMin(mat2) -mat1.cwiseMax(mat2) -mat1.cwiseAbs2() -mat1.cwiseAbs() -mat1.cwiseSqrt() -mat1.cwiseProduct(mat2) -mat1.cwiseQuotient(mat2)\endcode -</td><td>\code -mat1.array().min(mat2.array()) -mat1.array().max(mat2.array()) -mat1.array().abs2() -mat1.array().abs() -mat1.array().sqrt() -mat1.array() * mat2.array() -mat1.array() / mat2.array() -\endcode</td></tr> -</table> - -It is also very simple to apply any user defined function \c foo using DenseBase::unaryExpr together with std::ptr_fun: -\code mat1.unaryExpr(std::ptr_fun(foo))\endcode -Array operators:\arrayworld +In addition to the aforementioned operators, Eigen supports numerous coefficient-wise operator and functions. +Most of them unambiguously makes sense in array-world\arrayworld. The following operators are readily available for arrays, +or available through .array() for vectors and matrices: <table class="manual"> <tr><td>Arithmetic operators</td><td>\code @@ -400,28 +378,108 @@ array1 + scalar array1 - scalar array1 += scalar array1 -= scalar array1 < array2 array1 > array2 array1 < scalar array1 > scalar array1 <= array2 array1 >= array2 array1 <= scalar array1 >= scalar array1 == array2 array1 != array2 array1 == scalar array1 != scalar +array1.min(array2) array1.max(array2) array1.min(scalar) array1.max(scalar) \endcode</td></tr> -<tr><td>Trigo, power, and \n misc functions \n and the STL variants</td><td>\code -array1.min(array2) -array1.max(array2) +<tr><td>Trigo, power, and \n misc functions \n and the STL-like variants</td><td>\code array1.abs2() array1.abs() abs(array1) array1.sqrt() sqrt(array1) array1.log() log(array1) +array1.log10() log10(array1) array1.exp() exp(array1) -array1.pow(exponent) pow(array1,exponent) +array1.pow(array2) pow(array1,array2) +array1.pow(scalar) pow(array1,scalar) + pow(scalar,array2) array1.square() array1.cube() array1.inverse() + array1.sin() sin(array1) array1.cos() cos(array1) array1.tan() tan(array1) array1.asin() asin(array1) array1.acos() acos(array1) +array1.atan() atan(array1) +array1.sinh() sinh(array1) +array1.cosh() cosh(array1) +array1.tanh() tanh(array1) +array1.arg() arg(array1) + +array1.floor() floor(array1) +array1.ceil() ceil(array1) +array1.round() round(aray1) + +array1.isFinite() isfinite(array1) +array1.isInf() isinf(array1) +array1.isNaN() isnan(array1) +\endcode +</td></tr> +</table> + + +The following coefficient-wise operators are available for all kind of expressions (matrices, vectors, and arrays), and for both real or complex scalar types: + +<table class="manual"> +<tr><th>Eigen's API</th><th>STL-like APIs\arrayworld </th><th>Comments</th></tr> +<tr><td>\code +mat1.real() +mat1.imag() +mat1.conjugate() +\endcode +</td><td>\code +real(array1) +imag(array1) +conj(array1) +\endcode +</td><td> +\code + // read-write, no-op for real expressions + // read-only for real, read-write for complexes + // no-op for real expressions \endcode </td></tr> </table> +Some coefficient-wise operators are readily available for for matrices and vectors through the following cwise* methods: +<table class="manual"> +<tr><th>Matrix API \matrixworld</th><th>Via Array conversions</th></tr> +<tr><td>\code +mat1.cwiseMin(mat2) mat1.cwiseMin(scalar) +mat1.cwiseMax(mat2) mat1.cwiseMax(scalar) +mat1.cwiseAbs2() +mat1.cwiseAbs() +mat1.cwiseSqrt() +mat1.cwiseInverse() +mat1.cwiseProduct(mat2) +mat1.cwiseQuotient(mat2) +mat1.cwiseEqual(mat2) mat1.cwiseEqual(scalar) +mat1.cwiseNotEqual(mat2) +\endcode +</td><td>\code +mat1.array().min(mat2.array()) mat1.array().min(scalar) +mat1.array().max(mat2.array()) mat1.array().max(scalar) +mat1.array().abs2() +mat1.array().abs() +mat1.array().sqrt() +mat1.array().inverse() +mat1.array() * mat2.array() +mat1.array() / mat2.array() +mat1.array() == mat2.array() mat1.array() == scalar +mat1.array() != mat2.array() +\endcode</td></tr> +</table> +The main difference between the two API is that the one based on cwise* methods returns an expression in the matrix world, +while the second one (based on .array()) returns an array expression. +Recall that .array() has no cost, it only changes the available API and interpretation of the data. + +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), <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): +\code +mat1.unaryExpr(std::ptr_fun(foo)); +mat1.unaryExpr(std::ref(foo)); +mat1.unaryExpr([](double x) { return foo(x); }); +\endcode + + <a href="#" class="top">top</a> \section QuickRef_Reductions Reductions diff --git a/doc/SparseLinearSystems.dox b/doc/SparseLinearSystems.dox index c00be10d3..fc33b93e7 100644 --- a/doc/SparseLinearSystems.dox +++ b/doc/SparseLinearSystems.dox @@ -4,52 +4,87 @@ In Eigen, there are several methods available to solve linear systems when the c \eigenAutoToc -\section TutorialSparseDirectSolvers Sparse solvers +\section TutorialSparseSolverList List of sparse solvers -%Eigen currently provides a limited set of built-in solvers, as well as wrappers to external solver libraries. -They are summarized in the following table: +%Eigen currently provides a wide set of built-in solvers, as well as wrappers to external solver libraries. +They are summarized in the following tables: + +\subsection TutorialSparseSolverList_Direct Built-in direct solvers <table class="manual"> -<tr><th>Class</th><th>Module</th><th>Solver kind</th><th>Matrix kind</th><th>Features related to performance</th> - <th>Dependencies,License</th><th class="width20em"><p>Notes</p></th></tr> -<tr><td>SimplicialLLT </td><td>\link SparseCholesky_Module SparseCholesky \endlink</td><td>Direct LLt factorization</td><td>SPD</td><td>Fill-in reducing</td> - <td>built-in, LGPL</td> +<tr><th>Class</th><th>Solver kind</th><th>Matrix kind</th><th>Features related to performance</th> + <th>License</th><th class="width20em"><p>Notes</p></th></tr> + +<tr><td>SimplicialLLT \n <tt>\#include<Eigen/\link SparseCholesky_Module SparseCholesky\endlink></tt></td><td>Direct LLt factorization</td><td>SPD</td><td>Fill-in reducing</td> + <td>LGPL</td> <td>SimplicialLDLT is often preferable</td></tr> -<tr><td>SimplicialLDLT </td><td>\link SparseCholesky_Module SparseCholesky \endlink</td><td>Direct LDLt factorization</td><td>SPD</td><td>Fill-in reducing</td> - <td>built-in, LGPL</td> + +<tr><td>SimplicialLDLT \n <tt>\#include<Eigen/\link SparseCholesky_Module SparseCholesky\endlink></tt></td><td>Direct LDLt factorization</td><td>SPD</td><td>Fill-in reducing</td> + <td>LGPL</td> <td>Recommended for very sparse and not too large problems (e.g., 2D Poisson eq.)</td></tr> -<tr><td>ConjugateGradient</td><td>\link IterativeLinearSolvers_Module IterativeLinearSolvers \endlink</td><td>Classic iterative CG</td><td>SPD</td><td>Preconditionning</td> - <td>built-in, MPL2</td> - <td>Recommended for large symmetric problems (e.g., 3D Poisson eq.)</td></tr> -<tr><td>BiCGSTAB</td><td>\link IterativeLinearSolvers_Module IterativeLinearSolvers \endlink</td><td>Iterative stabilized bi-conjugate gradient</td><td>Square</td><td>Preconditionning</td> - <td>built-in, MPL2</td> - <td>To speedup the convergence, try it with the \ref IncompleteLUT preconditioner.</td></tr> -<tr><td>SparseLU</td> <td>\link SparseLU_Module SparseLU \endlink </td> <td>LU factorization </td> + +<tr><td>SparseLU \n <tt>\#include<Eigen/\link SparseLU_Module SparseLU\endlink></tt></td> <td>LU factorization </td> <td>Square </td><td>Fill-in reducing, Leverage fast dense algebra</td> - <td> built-in, MPL2</td> <td>optimized for small and large problems with irregular patterns </td></tr> -<tr><td>SparseQR</td> <td>\link SparseQR_Module SparseQR \endlink</td> <td> QR factorization</td> + <td>MPL2</td> + <td>optimized for small and large problems with irregular patterns </td></tr> + +<tr><td>SparseQR \n <tt>\#include<Eigen/\link SparseQR_Module SparseQR\endlink></tt></td> <td> QR factorization</td> <td>Any, rectangular</td><td> Fill-in reducing</td> - <td>built-in, MPL2</td><td>recommended for least-square problems, has a basic rank-revealing feature</td></tr> -<tr> <th colspan="7"> Wrappers to external solvers </th></tr> + <td>MPL2</td> + <td>recommended for least-square problems, has a basic rank-revealing feature</td></tr> + </table> + +\subsection TutorialSparseSolverList_Iterative Built-in iterative solvers + +<table class="manual"> +<tr><th>Class</th><th>Solver kind</th><th>Matrix kind</th><th>Supported preconditioners, [default]</th> + <th>License</th><th class="width20em"><p>Notes</p></th></tr> + +<tr><td>ConjugateGradient \n <tt>\#include<Eigen/\link IterativeLinearSolvers_Module IterativeLinearSolvers\endlink></tt></td> <td>Classic iterative CG</td><td>SPD</td> + <td>IdentityPreconditioner, [DiagonalPreconditioner], IncompleteCholesky</td> + <td>MPL2</td> + <td>Recommended for large symmetric problems (e.g., 3D Poisson eq.)</td></tr> + +<tr><td>LeastSquaresConjugateGradient \n <tt>\#include<Eigen/\link IterativeLinearSolvers_Module IterativeLinearSolvers\endlink></tt></td><td>CG for rectangular least-square problem</td><td>Rectangular</td> + <td>IdentityPreconditioner, [LeastSquareDiagonalPreconditioner]</td> + <td>MPL2</td> + <td>Solve for min |A'Ax-b|^2 without forming A'A</td></tr> + +<tr><td>BiCGSTAB \n <tt>\#include<Eigen/\link IterativeLinearSolvers_Module IterativeLinearSolvers\endlink></tt></td><td>Iterative stabilized bi-conjugate gradient</td><td>Square</td> + <td>IdentityPreconditioner, [DiagonalPreconditioner], IncompleteLUT</td> + <td>MPL2</td> + <td>To speedup the convergence, try it with the \ref IncompleteLUT preconditioner.</td></tr> +</table> + +\subsection TutorialSparseSolverList_Wrapper Wrappers to external solvers + +<table class="manual"> +<tr><th>Class</th><th>Module</th><th>Solver kind</th><th>Matrix kind</th><th>Features related to performance</th> + <th>Dependencies,License</th><th class="width20em"><p>Notes</p></th></tr> <tr><td>PastixLLT \n PastixLDLT \n PastixLU</td><td>\link PaStiXSupport_Module PaStiXSupport \endlink</td><td>Direct LLt, LDLt, LU factorizations</td><td>SPD \n SPD \n Square</td><td>Fill-in reducing, Leverage fast dense algebra, Multithreading</td> <td>Requires the <a href="http://pastix.gforge.inria.fr">PaStiX</a> package, \b CeCILL-C </td> <td>optimized for tough problems and symmetric patterns</td></tr> <tr><td>CholmodSupernodalLLT</td><td>\link CholmodSupport_Module CholmodSupport \endlink</td><td>Direct LLt factorization</td><td>SPD</td><td>Fill-in reducing, Leverage fast dense algebra</td> - <td>Requires the <a href="http://www.cise.ufl.edu/research/sparse/SuiteSparse/">SuiteSparse</a> package, \b GPL </td> + <td>Requires the <a href="http://www.suitesparse.com">SuiteSparse</a> package, \b GPL </td> <td></td></tr> <tr><td>UmfPackLU</td><td>\link UmfPackSupport_Module UmfPackSupport \endlink</td><td>Direct LU factorization</td><td>Square</td><td>Fill-in reducing, Leverage fast dense algebra</td> - <td>Requires the <a href="http://www.cise.ufl.edu/research/sparse/SuiteSparse/">SuiteSparse</a> package, \b GPL </td> + <td>Requires the <a href="http://www.suitesparse.com">SuiteSparse</a> package, \b GPL </td> <td></td></tr> <tr><td>SuperLU</td><td>\link SuperLUSupport_Module SuperLUSupport \endlink</td><td>Direct LU factorization</td><td>Square</td><td>Fill-in reducing, Leverage fast dense algebra</td> <td>Requires the <a href="http://crd-legacy.lbl.gov/~xiaoye/SuperLU/">SuperLU</a> library, (BSD-like)</td> <td></td></tr> <tr><td>SPQR</td><td>\link SPQRSupport_Module SPQRSupport \endlink </td> <td> QR factorization </td> <td> Any, rectangular</td><td>fill-in reducing, multithreaded, fast dense algebra</td> - <td> requires the <a href="http://www.cise.ufl.edu/research/sparse/SuiteSparse/">SuiteSparse</a> package, \b GPL </td><td>recommended for linear least-squares problems, has a rank-revealing feature</tr> + <td> requires the <a href="http://www.suitesparse.com">SuiteSparse</a> package, \b GPL </td><td>recommended for linear least-squares problems, has a rank-revealing feature</tr> +<tr><td>PardisoLLT \n PardisoLDLT \n PardisoLU</td><td>\link PardisoSupport_Module PardisoSupport \endlink</td><td>Direct LLt, LDLt, LU factorizations</td><td>SPD \n SPD \n Square</td><td>Fill-in reducing, Leverage fast dense algebra, Multithreading</td> + <td>Requires the <a href="http://eigen.tuxfamily.org/Counter/redirect_to_mkl.php">Intel MKL</a> package, \b Proprietary </td> + <td>optimized for tough problems patterns, see also \link TopicUsingIntelMKL using MKL with Eigen \endlink</td></tr> </table> Here \c SPD means symmetric positive definite. +\section TutorialSparseSolverConcept Sparse solver concept + All these solvers follow the same general concept. Here is a typical and general example: \code @@ -101,8 +136,10 @@ x2 = solver.solve(b2); \endcode The compute() method is equivalent to calling both analyzePattern() and factorize(). -Finally, each solver provides some specific features, such as determinant, access to the factors, controls of the iterations, and so on. -More details are availble in the documentations of the respective classes. +Each solver provides some specific features, such as determinant, access to the factors, controls of the iterations, and so on. +More details are available in the documentations of the respective classes. + +Finally, most of the iterative solvers, can also be used in a \b matrix-free context, see the following \link MatrixfreeSolverExample example \endlink. \section TheSparseCompute The Compute Step In the compute() function, the matrix is generally factorized: LLT for self-adjoint matrices, LDLT for general hermitian matrices, LU for non hermitian matrices and QR for rectangular matrices. These are the results of using direct solvers. For this class of solvers precisely, the compute step is further subdivided into analyzePattern() and factorize(). @@ -140,7 +177,16 @@ x2 = solver.solve(b2); For direct methods, the solution are computed at the machine precision. Sometimes, the solution need not be too accurate. In this case, the iterative methods are more suitable and the desired accuracy can be set before the solve step using \b setTolerance(). For all the available functions, please, refer to the documentation of the \link IterativeLinearSolvers_Module Iterative solvers module \endlink. \section BenchmarkRoutine -Most of the time, all you need is to know how much time it will take to qolve your system, and hopefully, what is the most suitable solver. In Eigen, we provide a benchmark routine that can be used for this purpose. It is very easy to use. In the build directory, navigate to bench/spbench and compile the routine by typing \b make \e spbenchsolver. Run it with --help option to get the list of all available options. Basically, the matrices to test should be in <a href="http://math.nist.gov/MatrixMarket/formats.html">MatrixMarket Coordinate format</a>, and the routine returns the statistics from all available solvers in Eigen. +Most of the time, all you need is to know how much time it will take to solve your system, and hopefully, what is the most suitable solver. In Eigen, we provide a benchmark routine that can be used for this purpose. It is very easy to use. In the build directory, navigate to bench/spbench and compile the routine by typing \b make \e spbenchsolver. Run it with --help option to get the list of all available options. Basically, the matrices to test should be in <a href="http://math.nist.gov/MatrixMarket/formats.html">MatrixMarket Coordinate format</a>, and the routine returns the statistics from all available solvers in Eigen. + +To export your matrices and right-hand-side vectors in the matrix-market format, you can the the unsupported SparseExtra module: +\code +#include <unsupported/Eigen/SparseExtra> +... +Eigen::saveMarket(A, "filename.mtx"); +Eigen::saveMarket(A, "filename_SPD.mtx", Eigen::Symmetric); // if A is symmetric-positive-definite +Eigen::saveMarketVector(B, "filename_b.mtx"); +\endcode The following table gives an example of XML statistics from several Eigen built-in and external solvers. <TABLE border="1"> diff --git a/doc/SparseQuickReference.dox b/doc/SparseQuickReference.dox index 4a33d0cc9..a25622e80 100644 --- a/doc/SparseQuickReference.dox +++ b/doc/SparseQuickReference.dox @@ -21,7 +21,7 @@ i.e either row major or column major. The default is column major. Most arithmet <td> Resize/Reserve</td> <td> \code - sm1.resize(m,n); //Change sm1 to a m x n matrix. + sm1.resize(m,n); // Change sm1 to a m x n matrix. sm1.reserve(nnz); // Allocate room for nnz nonzeros elements. \endcode </td> @@ -71,11 +71,10 @@ i.e either row major or column major. The default is column major. Most arithmet <td> Constant or Random Insertion</td> <td> \code -sm1.setZero(); // Set the matrix with zero elements -sm1.setConstant(val); //Replace all the nonzero values with val +sm1.setZero(); \endcode </td> -<td> The matrix sm1 should have been created before ???</td> +<td>Remove all non-zero coefficients</td> </tr> </table> @@ -152,10 +151,10 @@ It is easy to perform arithmetic operations on sparse matrices provided that the <td> Permutation </td> <td> \code -perm.indices(); // Reference to the vector of indices +perm.indices(); // Reference to the vector of indices sm1.twistedBy(perm); // Permute rows and columns -sm2 = sm1 * perm; //Permute the columns -sm2 = perm * sm1; // Permute the columns +sm2 = sm1 * perm; // Permute the columns +sm2 = perm * sm1; // Permute the columns \endcode </td> <td> @@ -182,9 +181,9 @@ sm2 = perm * sm1; // Permute the columns \section sparseotherops Other supported operations <table class="manual"> -<tr><th>Operations</th> <th> Code </th> <th> Notes</th> </tr> +<tr><th style="min-width:initial"> Code </th> <th> Notes</th> </tr> +<tr><td colspan="2">Sub-matrices</td></tr> <tr> -<td>Sub-matrices</td> <td> \code sm1.block(startRow, startCol, rows, cols); @@ -194,25 +193,31 @@ sm2 = perm * sm1; // Permute the columns sm1.bottomLeftCorner( rows, cols); sm1.bottomRightCorner( rows, cols); \endcode -</td> <td> </td> +</td><td> +Contrary to dense matrices, here <strong>all these methods are read-only</strong>.\n +See \ref TutorialSparse_SubMatrices and below for read-write sub-matrices. +</td> </tr> -<tr> -<td> Range </td> +<tr class="alt"><td colspan="2"> Range </td></tr> +<tr class="alt"> <td> \code - sm1.innerVector(outer); - sm1.innerVectors(start, size); - sm1.leftCols(size); - sm2.rightCols(size); - sm1.middleRows(start, numRows); - sm1.middleCols(start, numCols); - sm1.col(j); + sm1.innerVector(outer); // RW + sm1.innerVectors(start, size); // RW + sm1.leftCols(size); // RW + sm2.rightCols(size); // RO because sm2 is row-major + sm1.middleRows(start, numRows); // RO because sm1 is column-major + sm1.middleCols(start, numCols); // RW + sm1.col(j); // RW \endcode </td> -<td>A inner vector is either a row (for row-major) or a column (for column-major). As stated earlier, the evaluation can be done in a matrix with different storage order </td> +<td> +A inner vector is either a row (for row-major) or a column (for column-major).\n +As stated earlier, for a read-write sub-matrix (RW), the evaluation can be done in a matrix with different storage order. +</td> </tr> +<tr><td colspan="2"> Triangular and selfadjoint views</td></tr> <tr> -<td> Triangular and selfadjoint views</td> <td> \code sm2 = sm1.triangularview<Lower>(); @@ -223,26 +228,44 @@ sm2 = perm * sm1; // Permute the columns \code \endcode </td> </tr> -<tr> -<td>Triangular solve </td> +<tr class="alt"><td colspan="2">Triangular solve </td></tr> +<tr class="alt"> <td> \code dv2 = sm1.triangularView<Upper>().solve(dv1); - dv2 = sm1.topLeftCorner(size, size).triangularView<Lower>().solve(dv1); + dv2 = sm1.topLeftCorner(size, size) + .triangularView<Lower>().solve(dv1); \endcode </td> <td> For general sparse solve, Use any suitable module described at \ref TopicSparseSystems </td> </tr> +<tr><td colspan="2"> Low-level API</td></tr> <tr> -<td> Low-level API</td> <td> \code -sm1.valuePtr(); // Pointer to the values -sm1.innerIndextr(); // Pointer to the indices. -sm1.outerIndexPtr(); //Pointer to the beginning of each inner vector +sm1.valuePtr(); // Pointer to the values +sm1.innerIndextr(); // Pointer to the indices. +sm1.outerIndexPtr(); // Pointer to the beginning of each inner vector +\endcode +</td> +<td> +If the matrix is not in compressed form, makeCompressed() should be called before.\n +Note that these functions are mostly provided for interoperability purposes with external libraries.\n +A better access to the values of the matrix is done by using the InnerIterator class as described in \link TutorialSparse the Tutorial Sparse \endlink section</td> +</tr> +<tr class="alt"><td colspan="2">Mapping external buffers</td></tr> +<tr class="alt"> +<td> +\code +int outerIndexPtr[cols+1]; +int innerIndices[nnz]; +double values[nnz]; +Map<SparseMatrix<double> > sm1(rows,cols,nnz,outerIndexPtr, // read-write + innerIndices,values); +Map<const SparseMatrix<double> > sm2(...); // read-only \endcode </td> -<td> If the matrix is not in compressed form, makeCompressed() should be called before. Note that these functions are mostly provided for interoperability purposes with external libraries. A better access to the values of the matrix is done by using the InnerIterator class as described in \link TutorialSparse the Tutorial Sparse \endlink section</td> +<td>As for dense matrices, class Map<SparseMatrixType> can be used to see external buffers as an %Eigen's SparseMatrix object. </td> </tr> </table> */ diff --git a/doc/StlContainers.dox b/doc/StlContainers.dox index d8d0d529c..e0f8714a9 100644 --- a/doc/StlContainers.dox +++ b/doc/StlContainers.dox @@ -4,7 +4,7 @@ namespace Eigen { \eigenAutoToc -\section summary Executive summary +\section StlContainers_summary Executive summary Using STL containers on \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen types", or classes having members of such types, requires taking the following two steps: @@ -28,7 +28,7 @@ std::map<int, Eigen::Vector4f, std::less<int>, \endcode Note that the third parameter "std::less<int>" is just the default value, but we have to include it because we want to specify the fourth parameter, which is the allocator type. -\section vector The case of std::vector +\section StlContainers_vector The case of std::vector The situation with std::vector was even worse (explanation below) so we had to specialize it for the Eigen::aligned_allocator type. In practice you \b must use the Eigen::aligned_allocator (not another aligned allocator), \b and \#include <Eigen/StdVector>. diff --git a/doc/StructHavingEigenMembers.dox b/doc/StructHavingEigenMembers.dox index 74a8d5217..7fbed0eb0 100644 --- a/doc/StructHavingEigenMembers.dox +++ b/doc/StructHavingEigenMembers.dox @@ -4,11 +4,11 @@ namespace Eigen { \eigenAutoToc -\section summary Executive Summary +\section StructHavingEigenMembers_summary Executive Summary -If you define a structure having members of \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen types", you must overload its "operator new" so that it generates 16-bytes-aligned pointers. Fortunately, Eigen provides you with a macro EIGEN_MAKE_ALIGNED_OPERATOR_NEW that does that for you. +If you define a structure having members of \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen types", you must overload its "operator new" so that it generates 16-bytes-aligned pointers. Fortunately, %Eigen provides you with a macro EIGEN_MAKE_ALIGNED_OPERATOR_NEW that does that for you. -\section what What kind of code needs to be changed? +\section StructHavingEigenMembers_what What kind of code needs to be changed? The kind of code that needs to be changed is this: @@ -27,7 +27,7 @@ Foo *foo = new Foo; In other words: you have a class that has as a member a \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen object", and then you dynamically create an object of that class. -\section how How should such code be modified? +\section StructHavingEigenMembers_how How should such code be modified? Very easy, you just need to put a EIGEN_MAKE_ALIGNED_OPERATOR_NEW macro in a public part of your class, like this: @@ -48,9 +48,9 @@ Foo *foo = new Foo; This macro makes "new Foo" always return an aligned pointer. -If this approach is too intrusive, see also the \ref othersolutions. +If this approach is too intrusive, see also the \ref StructHavingEigenMembers_othersolutions "other solutions". -\section why Why is this needed? +\section StructHavingEigenMembers_why Why is this needed? OK let's say that your code looks like this: @@ -67,7 +67,7 @@ class Foo Foo *foo = new Foo; \endcode -A Eigen::Vector2d consists of 2 doubles, which is 128 bits. Which is exactly the size of a SSE packet, which makes it possible to use SSE for all sorts of operations on this vector. But SSE instructions (at least the ones that Eigen uses, which are the fast ones) require 128-bit alignment. Otherwise you get a segmentation fault. +A Eigen::Vector2d consists of 2 doubles, which is 128 bits. Which is exactly the size of a SSE packet, which makes it possible to use SSE for all sorts of operations on this vector. But SSE instructions (at least the ones that %Eigen uses, which are the fast ones) require 128-bit alignment. Otherwise you get a segmentation fault. For this reason, Eigen takes care by itself to require 128-bit alignment for Eigen::Vector2d, by doing two things: \li Eigen requires 128-bit alignment for the Eigen::Vector2d's array (of 2 doubles). With GCC, this is done with a __attribute__ ((aligned(16))). @@ -81,7 +81,7 @@ The alignment attribute of the member v is then relative to the start of the cla The solution is to let class Foo have an aligned "operator new", as we showed in the previous section. -\section movetotop Should I then put all the members of Eigen types at the beginning of my class? +\section StructHavingEigenMembers_movetotop Should I then put all the members of Eigen types at the beginning of my class? That's not required. Since Eigen takes care of declaring 128-bit alignment, all members that need it are automatically 128-bit aligned relatively to the class. So code like this works fine: @@ -95,15 +95,15 @@ public: }; \endcode -\section dynamicsize What about dynamic-size matrices and vectors? +\section StructHavingEigenMembers_dynamicsize What about dynamic-size matrices and vectors? Dynamic-size matrices and vectors, such as Eigen::VectorXd, allocate dynamically their own array of coefficients, so they take care of requiring absolute alignment automatically. So they don't cause this issue. The issue discussed here is only with \ref TopicFixedSizeVectorizable "fixed-size vectorizable matrices and vectors". -\section bugineigen So is this a bug in Eigen? +\section StructHavingEigenMembers_bugineigen So is this a bug in Eigen? No, it's not our bug. It's more like an inherent problem of the C++98 language specification, and seems to be taken care of in the upcoming language revision: <a href="http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2007/n2341.pdf">see this document</a>. -\section conditional What if I want to do this conditionnally (depending on template parameters) ? +\section StructHavingEigenMembers_conditional What if I want to do this conditionnally (depending on template parameters) ? For this situation, we offer the macro EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF(NeedsToAlign). It will generate aligned operators like EIGEN_MAKE_ALIGNED_OPERATOR_NEW if NeedsToAlign is true. It will generate operators with the default alignment if NeedsToAlign is false. @@ -128,7 +128,7 @@ Foo<3> *foo3 = new Foo<3>; // foo3 has only the system default alignment guarant \endcode -\section othersolutions Other solutions +\section StructHavingEigenMembers_othersolutions Other solutions In case putting the EIGEN_MAKE_ALIGNED_OPERATOR_NEW macro everywhere is too intrusive, there exists at least two other solutions. diff --git a/doc/TemplateKeyword.dox b/doc/TemplateKeyword.dox index c9944ae05..b84cfdae9 100644 --- a/doc/TemplateKeyword.dox +++ b/doc/TemplateKeyword.dox @@ -73,23 +73,23 @@ for operator<". The reason that the \c template keyword is necessary in the last example has to do with the rules for how templates are supposed to be compiled in C++. The compiler has to check the code for correct syntax at the point where the template is defined, without knowing the actual value of the template arguments (\c Derived1 -and \c Derived2 in the example). That means that the compiler cannot know that <tt>dst.triangularPart</tt> is +and \c Derived2 in the example). That means that the compiler cannot know that <tt>dst.triangularView</tt> is a member template and that the following < symbol is part of the delimiter for the template -parameter. Another possibility would be that <tt>dst.triangularPart</tt> is a member variable with the < +parameter. Another possibility would be that <tt>dst.triangularView</tt> is a member variable with the < symbol refering to the <tt>operator<()</tt> function. In fact, the compiler should choose the second -possibility, according to the standard. If <tt>dst.triangularPart</tt> is a member template (as in our case), +possibility, according to the standard. If <tt>dst.triangularView</tt> is a member template (as in our case), the programmer should specify this explicitly with the \c template keyword and write <tt>dst.template -triangularPart</tt>. +triangularView</tt>. The precise rules are rather complicated, but ignoring some subtleties we can summarize them as follows: - A <em>dependent name</em> is name that depends (directly or indirectly) on a template parameter. In the example, \c dst is a dependent name because it is of type <tt>MatrixBase<Derived1></tt> which depends on the template parameter \c Derived1. -- If the code contains either one of the contructions <tt>xxx.yyy</tt> or <tt>xxx->yyy</tt> and \c xxx is a +- If the code contains either one of the constructs <tt>xxx.yyy</tt> or <tt>xxx->yyy</tt> and \c xxx is a dependent name and \c yyy refers to a member template, then the \c template keyword must be used before \c yyy, leading to <tt>xxx.template yyy</tt> or <tt>xxx->template yyy</tt>. -- If the code contains the contruction <tt>xxx::yyy</tt> and \c xxx is a dependent name and \c yyy refers to a - member typedef, then the \c typename keyword must be used before the whole construction, leading to +- If the code contains the construct <tt>xxx::yyy</tt> and \c xxx is a dependent name and \c yyy refers to a + member typedef, then the \c typename keyword must be used before the whole construct, leading to <tt>typename xxx::yyy</tt>. As an example where the \c typename keyword is required, consider the following code in \ref TutorialSparse diff --git a/doc/TopicAliasing.dox b/doc/TopicAliasing.dox index c2654aed2..a8f164428 100644 --- a/doc/TopicAliasing.dox +++ b/doc/TopicAliasing.dox @@ -153,10 +153,11 @@ not necessary to evaluate the right-hand side explicitly. \section TopicAliasingMatrixMult Aliasing and matrix multiplication -Matrix multiplication is the only operation in %Eigen that assumes aliasing by default. Thus, if \c matA is a -matrix, then the statement <tt>matA = matA * matA;</tt> is safe. All other operations in %Eigen assume that -there are no aliasing problems, either because the result is assigned to a different matrix or because it is a -component-wise operation. +Matrix multiplication is the only operation in %Eigen that assumes aliasing by default, <strong>under the +condition that the destination matrix is not resized</strong>. +Thus, if \c matA is a \b squared matrix, then the statement <tt>matA = matA * matA;</tt> is safe. +All other operations in %Eigen assume that there are no aliasing problems, +either because the result is assigned to a different matrix or because it is a component-wise operation. <table class="example"> <tr><th>Example</th><th>Output</th></tr> @@ -198,6 +199,27 @@ may get wrong results: \verbinclude TopicAliasing_mult3.out </td></tr></table> +Moreover, starting in Eigen 3.3, aliasing is \b not assumed if the destination matrix is resized and the product is not directly assigned to the destination. +Therefore, the following example is also wrong: + +<table class="example"> +<tr><th>Example</th><th>Output</th></tr> +<tr><td> +\include TopicAliasing_mult4.cpp +</td> +<td> +\verbinclude TopicAliasing_mult4.out +</td></tr></table> + +As for any aliasing issue, you can resolve it by explicitly evaluating the expression prior to assignment: +<table class="example"> +<tr><th>Example</th><th>Output</th></tr> +<tr><td> +\include TopicAliasing_mult5.cpp +</td> +<td> +\verbinclude TopicAliasing_mult5.out +</td></tr></table> \section TopicAliasingSummary Summary diff --git a/doc/TopicAssertions.dox b/doc/TopicAssertions.dox index 4ead40174..c8b4d84f2 100644 --- a/doc/TopicAssertions.dox +++ b/doc/TopicAssertions.dox @@ -16,7 +16,7 @@ Both eigen_assert and eigen_plain_assert are defined in Macros.h. Defining eigen #include <stdexcept> #undef eigen_assert #define eigen_assert(x) \ - if (!x) { throw (std::runtime_error("Put your message here")); } + if (!(x)) { throw (std::runtime_error("Put your message here")); } \endcode \subsection DisableAssert Disabling assertions diff --git a/doc/TopicCMakeGuide.dox b/doc/TopicCMakeGuide.dox new file mode 100644 index 000000000..cf767d0dd --- /dev/null +++ b/doc/TopicCMakeGuide.dox @@ -0,0 +1,56 @@ +namespace Eigen { + +/** + +\page TopicCMakeGuide Using %Eigen in CMake Projects + +%Eigen provides native CMake support which allows the library to be easily +used in CMake projects. + +\note %CMake 3.0 (or later) is required to enable this functionality. + +%Eigen exports a CMake target called `Eigen3::Eigen` which can be imported +using the `find_package` CMake command and used by calling +`target_link_libraries` as in the following example: +\code{.cmake} +cmake_minimum_required (VERSION 3.0) +project (myproject) + +find_package (Eigen3 3.3 REQUIRED NO_MODULE) + +add_executable (example example.cpp) +target_link_libraries (example Eigen3::Eigen) +\endcode + +The above code snippet must be placed in a file called `CMakeLists.txt` alongside +`example.cpp`. After running +\code{.sh} +$ cmake path-to-example-directory +\endcode +CMake will produce project files that generate an executable called `example` +which requires at least version 3.3 of %Eigen. Here, `path-to-example-directory` +is the path to the directory that contains both `CMakeLists.txt` and +`example.cpp`. + +Do not forget to set the <a href="https://cmake.org/cmake/help/v3.7/variable/CMAKE_PREFIX_PATH.html">\c CMAKE_PREFIX_PATH </a> variable if Eigen is not installed in a default location or if you want to pick a specific version. For instance: +\code{.sh} +$ cmake path-to-example-directory -DCMAKE_PREFIX_PATH=$HOME/mypackages +\endcode +An alternative is to set the \c Eigen3_DIR cmake's variable to the respective path containing the \c Eigen3*.cmake files. For instance: +\code{.sh} +$ cmake path-to-example-directory -DEigen3_DIR=$HOME/mypackages/share/eigen3/cmake/ +\endcode + +If the `REQUIRED` option is omitted when locating %Eigen using +`find_package`, one can check whether the package was found as follows: +\code{.cmake} +find_package (Eigen3 3.3 NO_MODULE) + +if (TARGET Eigen3::Eigen) + # Use the imported target +endif (TARGET Eigen3::Eigen) +\endcode + +*/ + +} diff --git a/doc/TopicLazyEvaluation.dox b/doc/TopicLazyEvaluation.dox index 393bc41d8..101ef8c72 100644 --- a/doc/TopicLazyEvaluation.dox +++ b/doc/TopicLazyEvaluation.dox @@ -36,7 +36,7 @@ Here is now a more involved example: Eigen chooses lazy evaluation at every stage in that example, which is clearly the correct choice. In fact, lazy evaluation is the "default choice" and Eigen will choose it except in a few circumstances. -<b>The first circumstance</b> in which Eigen chooses immediate evaluation, is when it sees an assignment <tt>a = b;</tt> and the expression \c b has the evaluate-before-assigning \link flags flag\endlink. The most important example of such an expression is the \link GeneralProduct matrix product expression\endlink. For example, when you do +<b>The first circumstance</b> in which Eigen chooses immediate evaluation, is when it sees an assignment <tt>a = b;</tt> and the expression \c b has the evaluate-before-assigning \link flags flag\endlink. The most important example of such an expression is the \link Product matrix product expression\endlink. For example, when you do \code matrix = matrix * matrix; \endcode @@ -48,7 +48,7 @@ What if you know that the result does no alias the operand of the product and wa Here, since we know that matrix2 is not the same matrix as matrix1, we know that lazy evaluation is not dangerous, so we may force lazy evaluation. Concretely, the effect of noalias() here is to bypass the evaluate-before-assigning \link flags flag\endlink. -<b>The second circumstance</b> in which Eigen chooses immediate evaluation, is when it sees a nested expression such as <tt>a + b</tt> where \c b is already an expression having the evaluate-before-nesting \link flags flag\endlink. Again, the most important example of such an expression is the \link GeneralProduct matrix product expression\endlink. For example, when you do +<b>The second circumstance</b> in which Eigen chooses immediate evaluation, is when it sees a nested expression such as <tt>a + b</tt> where \c b is already an expression having the evaluate-before-nesting \link flags flag\endlink. Again, the most important example of such an expression is the \link Product matrix product expression\endlink. For example, when you do \code matrix1 = matrix2 + matrix3 * matrix4; \endcode diff --git a/doc/TopicLinearAlgebraDecompositions.dox b/doc/TopicLinearAlgebraDecompositions.dox index 8649cc27b..491470627 100644 --- a/doc/TopicLinearAlgebraDecompositions.dox +++ b/doc/TopicLinearAlgebraDecompositions.dox @@ -4,6 +4,7 @@ namespace Eigen { This page presents a catalogue of the dense matrix decompositions offered by Eigen. For an introduction on linear solvers and decompositions, check this \link TutorialLinearAlgebra page \endlink. +To get an overview of the true relative speed of the different decomposition, check this \link DenseDecompositionBenchmark benchmark \endlink. \section TopicLinAlgBigTable Catalogue of decompositions offered by Eigen @@ -116,7 +117,7 @@ For an introduction on linear solvers and decompositions, check this \link Tutor <td>JacobiSVD (two-sided)</td> <td>-</td> <td>Slow (but fast for small matrices)</td> - <td>Excellent-Proven<sup><a href="#note3">3</a></sup></td> + <td>Proven<sup><a href="#note3">3</a></sup></td> <td>Yes</td> <td>Singular values/vectors, least squares</td> <td>Yes (and does least squares)</td> @@ -132,7 +133,7 @@ For an introduction on linear solvers and decompositions, check this \link Tutor <td>Yes</td> <td>Eigenvalues/vectors</td> <td>-</td> - <td>Good</td> + <td>Excellent</td> <td><em>Closed forms for 2x2 and 3x3</em></td> </tr> @@ -249,13 +250,14 @@ For an introduction on linear solvers and decompositions, check this \link Tutor <dt><b>Implicit Multi Threading (MT)</b></dt> <dd>Means the algorithm can take advantage of multicore processors via OpenMP. "Implicit" means the algortihm itself is not parallelized, but that it relies on parallelized matrix-matrix product rountines.</dd> <dt><b>Explicit Multi Threading (MT)</b></dt> - <dd>Means the algorithm is explicitely parallelized to take advantage of multicore processors via OpenMP.</dd> + <dd>Means the algorithm is explicitly parallelized to take advantage of multicore processors via OpenMP.</dd> <dt><b>Meta-unroller</b></dt> <dd>Means the algorithm is automatically and explicitly unrolled for very small fixed size matrices.</dd> <dt><b></b></dt> <dd></dd> </dl> + */ } diff --git a/doc/TopicMultithreading.dox b/doc/TopicMultithreading.dox index fb944af29..47c9b261f 100644 --- a/doc/TopicMultithreading.dox +++ b/doc/TopicMultithreading.dox @@ -8,22 +8,26 @@ Some Eigen's algorithms can exploit the multiple cores present in your hardware. * GCC: \c -fopenmp * ICC: \c -openmp * MSVC: check the respective option in the build properties. -You can control the number of thread that will be used using either the OpenMP API or Eiegn's API using the following priority: +You can control the number of thread that will be used using either the OpenMP API or Eigen's API using the following priority: \code OMP_NUM_THREADS=n ./my_program omp_set_num_threads(n); Eigen::setNbThreads(n); \endcode -Unless setNbThreads has been called, Eigen uses the number of threads specified by OpenMP. You can restore this bahavior by calling \code setNbThreads(0); \endcode +Unless setNbThreads has been called, Eigen uses the number of threads specified by OpenMP. You can restore this behavior by calling \code setNbThreads(0); \endcode You can query the number of threads that will be used with: \code -n = Eigen::nbThreads(n); +n = Eigen::nbThreads( ); \endcode You can disable Eigen's multi threading at compile time by defining the EIGEN_DONT_PARALLELIZE preprocessor token. Currently, the following algorithms can make use of multi-threading: - * general matrix - matrix products - * PartialPivLU + - general dense matrix - matrix products + - PartialPivLU + - row-major-sparse * dense vector/matrix products + - ConjugateGradient with \c Lower|Upper as the \c UpLo template parameter. + - BiCGSTAB with a row-major sparse matrix format. + - LeastSquaresConjugateGradient \section TopicMultiThreading_UsingEigenWithMT Using Eigen in a multi-threaded application @@ -39,7 +43,9 @@ int main(int argc, char** argv) } \endcode -\warning note that all functions generating random matrices are \b not re-entrant nor thread-safe. Those include DenseBase::Random(), and DenseBase::setRandom() despite a call to Eigen::initParallel(). This is because these functions are based on std::rand which is not re-entrant. For thread-safe random generator, we recommend the use of boost::random of c++11 random feature. +\note With Eigen 3.3, and a fully C++11 compliant compiler (i.e., <a href="http://en.cppreference.com/w/cpp/language/storage_duration#Static_local_variables">thread-safe static local variable initialization</a>), then calling \c initParallel() is optional. + +\warning note that all functions generating random matrices are \b not re-entrant nor thread-safe. Those include DenseBase::Random(), and DenseBase::setRandom() despite a call to Eigen::initParallel(). This is because these functions are based on std::rand which is not re-entrant. For thread-safe random generator, we recommend the use of boost::random or c++11 random feature. In the case your application is parallelized with OpenMP, you might want to disable Eigen's own parallization as detailed in the previous section. diff --git a/doc/TutorialArrayClass.dox b/doc/TutorialArrayClass.dox index 6432684aa..f6f351091 100644 --- a/doc/TutorialArrayClass.dox +++ b/doc/TutorialArrayClass.dox @@ -157,7 +157,7 @@ The following example shows how to use array operations on a Matrix object by em * to multiply them coefficient-wise and assigns the result to the matrix variable \c result (this is legal because Eigen allows assigning array expressions to matrix variables). -As a matter of fact, this usage case is so common that Eigen provides a \link MatrixBase::cwiseProduct() const +As a matter of fact, this usage case is so common that Eigen provides a \link MatrixBase::cwiseProduct const .cwiseProduct(.) \endlink method for matrices to compute the coefficient-wise product. This is also shown in the example program. diff --git a/doc/TutorialGeometry.dox b/doc/TutorialGeometry.dox index 372a275de..2e1420f98 100644 --- a/doc/TutorialGeometry.dox +++ b/doc/TutorialGeometry.dox @@ -126,11 +126,12 @@ Apply the transformation to a \b vector </td><td>\code VectorNf vec1, vec2; vec2 = t.linear() * vec1;\endcode</td></tr> <tr><td> -Apply a \em general transformation \n to a \b normal \b vector -(<a href="http://femto.cs.uiuc.edu/faqs/cga-faq.html#S5.27">explanations</a>)</td><td>\code +Apply a \em general transformation \n to a \b normal \b vector \n +</td><td>\code VectorNf n1, n2; MatrixNf normalMatrix = t.linear().inverse().transpose(); n2 = (normalMatrix * n1).normalized();\endcode</td></tr> +<tr><td colspan="2">(See subject 5.27 of this <a href="http://www.faqs.org/faqs/graphics/algorithms-faq">faq</a> for the explanations)</td></tr> <tr class="alt"><td> Apply a transformation with \em pure \em rotation \n to a \b normal \b vector (no scaling, no shear)</td><td>\code @@ -231,8 +232,8 @@ On the other hand, since there exist 24 different conventions, they are pretty c to create a rotation matrix according to the 2-1-2 convention.</td><td>\code Matrix3f m; m = AngleAxisf(angle1, Vector3f::UnitZ()) -* * AngleAxisf(angle2, Vector3f::UnitY()) -* * AngleAxisf(angle3, Vector3f::UnitZ()); + * AngleAxisf(angle2, Vector3f::UnitY()) + * AngleAxisf(angle3, Vector3f::UnitZ()); \endcode</td></tr> </table> diff --git a/doc/TutorialLinearAlgebra.dox b/doc/TutorialLinearAlgebra.dox index b09f3543e..cb92ceeae 100644 --- a/doc/TutorialLinearAlgebra.dox +++ b/doc/TutorialLinearAlgebra.dox @@ -40,8 +40,9 @@ depending on your matrix and the trade-off you want to make: <tr> <th>Decomposition</th> <th>Method</th> - <th>Requirements on the matrix</th> - <th>Speed</th> + <th>Requirements<br/>on the matrix</th> + <th>Speed<br/> (small-to-medium)</th> + <th>Speed<br/> (large)</th> <th>Accuracy</th> </tr> <tr> @@ -49,6 +50,7 @@ depending on your matrix and the trade-off you want to make: <td>partialPivLu()</td> <td>Invertible</td> <td>++</td> + <td>++</td> <td>+</td> </tr> <tr class="alt"> @@ -56,6 +58,7 @@ depending on your matrix and the trade-off you want to make: <td>fullPivLu()</td> <td>None</td> <td>-</td> + <td>- -</td> <td>+++</td> </tr> <tr> @@ -63,20 +66,23 @@ depending on your matrix and the trade-off you want to make: <td>householderQr()</td> <td>None</td> <td>++</td> + <td>++</td> <td>+</td> </tr> <tr class="alt"> <td>ColPivHouseholderQR</td> <td>colPivHouseholderQr()</td> <td>None</td> - <td>+</td> <td>++</td> + <td>-</td> + <td>+++</td> </tr> <tr> <td>FullPivHouseholderQR</td> <td>fullPivHouseholderQr()</td> <td>None</td> <td>-</td> + <td>- -</td> <td>+++</td> </tr> <tr class="alt"> @@ -84,21 +90,31 @@ depending on your matrix and the trade-off you want to make: <td>llt()</td> <td>Positive definite</td> <td>+++</td> + <td>+++</td> <td>+</td> </tr> <tr> <td>LDLT</td> <td>ldlt()</td> - <td>Positive or negative semidefinite</td> + <td>Positive or negative<br/> semidefinite</td> <td>+++</td> + <td>+</td> <td>++</td> </tr> + <tr class="alt"> + <td>JacobiSVD</td> + <td>jacobiSvd()</td> + <td>None</td> + <td>- -</td> + <td>- - -</td> + <td>+++</td> + </tr> </table> All of these decompositions offer a solve() method that works as in the above example. For example, if your matrix is positive definite, the above table says that a very good -choice is then the LDLT decomposition. Here's an example, also demonstrating that using a general +choice is then the LLT or LDLT decomposition. Here's an example, also demonstrating that using a general matrix (not a vector) as right hand side is possible. <table class="example"> @@ -167,8 +183,8 @@ Here is an example: \section TutorialLinAlgLeastsquares Least squares solving -The best way to do least squares solving is with a SVD decomposition. Eigen provides one as the JacobiSVD class, and its solve() -is doing least-squares solving. +The most accurate method to do least squares solving is with a SVD decomposition. Eigen provides one +as the JacobiSVD class, and its solve() is doing least-squares solving. Here is an example: <table class="example"> @@ -179,9 +195,10 @@ Here is an example: </tr> </table> -Another way, potentially faster but less reliable, is to use a LDLT decomposition -of the normal matrix. In any case, just read any reference text on least squares, and it will be very easy for you -to implement any linear least squares computation on top of Eigen. +Another methods, potentially faster but less reliable, are to use a Cholesky decomposition of the +normal matrix or a QR decomposition. Our page on \link LeastSquares least squares solving \endlink +has more details. + \section TutorialLinAlgSeparateComputation Separating the computation from the construction diff --git a/doc/TutorialReductionsVisitorsBroadcasting.dox b/doc/TutorialReductionsVisitorsBroadcasting.dox index 992cf6f34..f5322b4a6 100644 --- a/doc/TutorialReductionsVisitorsBroadcasting.dox +++ b/doc/TutorialReductionsVisitorsBroadcasting.dox @@ -32,7 +32,7 @@ Eigen also provides the \link MatrixBase::norm() norm() \endlink method, which r These operations can also operate on matrices; in that case, a n-by-p matrix is seen as a vector of size (n*p), so for example the \link MatrixBase::norm() norm() \endlink method returns the "Frobenius" or "Hilbert-Schmidt" norm. We refrain from speaking of the \f$\ell^2\f$ norm of a matrix because that can mean different things. -If you want other \f$\ell^p\f$ norms, use the \link MatrixBase::lpNorm() lpNnorm<p>() \endlink method. The template parameter \a p can take the special value \a Infinity if you want the \f$\ell^\infty\f$ norm, which is the maximum of the absolute values of the coefficients. +If you want other coefficient-wise \f$\ell^p\f$ norms, use the \link MatrixBase::lpNorm lpNorm<p>() \endlink method. The template parameter \a p can take the special value \a Infinity if you want the \f$\ell^\infty\f$ norm, which is the maximum of the absolute values of the coefficients. The following example demonstrates these methods. @@ -45,6 +45,17 @@ The following example demonstrates these methods. \verbinclude Tutorial_ReductionsVisitorsBroadcasting_reductions_norm.out </td></tr></table> +\b Operator \b norm: The 1-norm and \f$\infty\f$-norm <a href="https://en.wikipedia.org/wiki/Operator_norm">matrix operator norms</a> can easily be computed as follows: +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr><td> +\include Tutorial_ReductionsVisitorsBroadcasting_reductions_operatornorm.cpp +</td> +<td> +\verbinclude Tutorial_ReductionsVisitorsBroadcasting_reductions_operatornorm.out +</td></tr></table> +See below for more explanations on the syntax of these expressions. + \subsection TutorialReductionsVisitorsBroadcastingReductionsBool Boolean reductions The following reductions operate on boolean values: @@ -79,7 +90,7 @@ Array. The arguments passed to a visitor are pointers to the variables where the row and column position are to be stored. These variables should be of type -\link DenseBase::Index Index \endlink, as shown below: +\link Eigen::Index Index \endlink, as shown below: <table class="example"> <tr><th>Example:</th><th>Output:</th></tr> @@ -90,17 +101,16 @@ row and column position are to be stored. These variables should be of type \verbinclude Tutorial_ReductionsVisitorsBroadcasting_visitors.out </td></tr></table> -Note that both functions also return the value of the minimum or maximum coefficient if needed, -as if it was a typical reduction operation. +Both functions also return the value of the minimum or maximum coefficient. \section TutorialReductionsVisitorsBroadcastingPartialReductions Partial reductions Partial reductions are reductions that can operate column- or row-wise on a Matrix or Array, applying the reduction operation on each column or row and -returning a column or row-vector with the corresponding values. Partial reductions are applied +returning a column or row vector with the corresponding values. Partial reductions are applied with \link DenseBase::colwise() colwise() \endlink or \link DenseBase::rowwise() rowwise() \endlink. A simple example is obtaining the maximum of the elements -in each column in a given matrix, storing the result in a row-vector: +in each column in a given matrix, storing the result in a row vector: <table class="example"> <tr><th>Example:</th><th>Output:</th></tr> @@ -122,8 +132,7 @@ The same operation can be performed row-wise: \verbinclude Tutorial_ReductionsVisitorsBroadcasting_rowwise.out </td></tr></table> -<b>Note that column-wise operations return a 'row-vector' while row-wise operations -return a 'column-vector'</b> +<b>Note that column-wise operations return a row vector, while row-wise operations return a column vector.</b> \subsection TutorialReductionsVisitorsBroadcastingPartialReductionsCombined Combining partial reductions with other operations It is also possible to use the result of a partial reduction to do further processing. @@ -165,7 +174,7 @@ The concept behind broadcasting is similar to partial reductions, with the diffe constructs an expression where a vector (column or row) is interpreted as a matrix by replicating it in one direction. -A simple example is to add a certain column-vector to each column in a matrix. +A simple example is to add a certain column vector to each column in a matrix. This can be accomplished with: <table class="example"> @@ -242,7 +251,7 @@ is a new matrix whose size is the same as matrix <tt>m</tt>: \f[ \f] - <tt>(m.colwise() - v).colwise().squaredNorm()</tt> is a partial reduction, computing the squared norm column-wise. The result of -this operation is a row-vector where each coefficient is the squared Euclidean distance between each column in <tt>m</tt> and <tt>v</tt>: \f[ +this operation is a row vector where each coefficient is the squared Euclidean distance between each column in <tt>m</tt> and <tt>v</tt>: \f[ \mbox{(m.colwise() - v).colwise().squaredNorm()} = \begin{bmatrix} 1 & 505 & 32 & 50 diff --git a/doc/TutorialReshapeSlicing.dox b/doc/TutorialReshapeSlicing.dox new file mode 100644 index 000000000..3730a5de6 --- /dev/null +++ b/doc/TutorialReshapeSlicing.dox @@ -0,0 +1,65 @@ +namespace Eigen { + +/** \eigenManualPage TutorialReshapeSlicing Reshape and Slicing + +%Eigen does not expose convenient methods to take slices or to reshape a matrix yet. +Nonetheless, such features can easily be emulated using the Map class. + +\eigenAutoToc + +\section TutorialReshape Reshape + +A reshape operation consists in modifying the sizes of a matrix while keeping the same coefficients. +Instead of modifying the input matrix itself, which is not possible for compile-time sizes, the approach consist in creating a different \em view on the storage using class Map. +Here is a typical example creating a 1D linear view of a matrix: + +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr><td> +\include Tutorial_ReshapeMat2Vec.cpp +</td> +<td> +\verbinclude Tutorial_ReshapeMat2Vec.out +</td></tr></table> + +Remark how the storage order of the input matrix modifies the order of the coefficients in the linear view. +Here is another example reshaping a 2x6 matrix to a 6x2 one: +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr><td> +\include Tutorial_ReshapeMat2Mat.cpp +</td> +<td> +\verbinclude Tutorial_ReshapeMat2Mat.out +</td></tr></table> + + + +\section TutorialSlicing Slicing + +Slicing consists in taking a set of rows, columns, or elements, uniformly spaced within a matrix. +Again, the class Map allows to easily mimic this feature. + +For instance, one can skip every P elements in a vector: +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr><td> +\include Tutorial_SlicingVec.cpp +</td> +<td> +\verbinclude Tutorial_SlicingVec.out +</td></tr></table> + +One can olso take one column over three using an adequate outer-stride or inner-stride depending on the actual storage order: +<table class="example"> +<tr><th>Example:</th><th>Output:</th></tr> +<tr><td> +\include Tutorial_SlicingCol.cpp +</td> +<td> +\verbinclude Tutorial_SlicingCol.out +</td></tr></table> + +*/ + +} diff --git a/doc/TutorialSparse.dox b/doc/TutorialSparse.dox index 835c59354..352907408 100644 --- a/doc/TutorialSparse.dox +++ b/doc/TutorialSparse.dox @@ -83,7 +83,7 @@ There is no notion of compressed/uncompressed mode for a SparseVector. \section TutorialSparseExample First example -Before describing each individual class, let's start with the following typical example: solving the Laplace equation \f$ \nabla u = 0 \f$ on a regular 2D grid using a finite difference scheme and Dirichlet boundary conditions. +Before describing each individual class, let's start with the following typical example: solving the Laplace equation \f$ \Delta u = 0 \f$ on a regular 2D grid using a finite difference scheme and Dirichlet boundary conditions. Such problem can be mathematically expressed as a linear problem of the form \f$ Ax=b \f$ where \f$ x \f$ is the vector of \c m unknowns (in our case, the values of the pixels), \f$ b \f$ is the right hand side vector resulting from the boundary conditions, and \f$ A \f$ is an \f$ m \times m \f$ matrix containing only a few non-zero elements resulting from the discretization of the Laplacian operator. <table class="manual"> @@ -241,11 +241,11 @@ In the following \em sm denotes a sparse matrix, \em sv a sparse vector, \em dm sm1.real() sm1.imag() -sm1 0.5*sm1 sm1+sm2 sm1-sm2 sm1.cwiseProduct(sm2) \endcode -However, a strong restriction is that the storage orders must match. For instance, in the following example: +However, <strong>a strong restriction is that the storage orders must match</strong>. For instance, in the following example: \code sm4 = sm1 + sm2 + sm3; \endcode -sm1, sm2, and sm3 must all be row-major or all column major. +sm1, sm2, and sm3 must all be row-major or all column-major. On the other hand, there is no restriction on the target matrix sm4. For instance, this means that for computing \f$ A^T + A \f$, the matrix \f$ A^T \f$ must be evaluated into a temporary matrix of compatible storage order: \code @@ -257,7 +257,14 @@ Binary coefficient wise operators can also mix sparse and dense expressions: \code sm2 = sm1.cwiseProduct(dm1); dm2 = sm1 + dm1; +dm2 = dm1 - sm1; \endcode +Performance-wise, the adding/subtracting sparse and dense matrices is better performed in two steps. For instance, instead of doing <tt>dm2 = sm1 + dm1</tt>, better write: +\code +dm2 = dm1; +dm2 += sm1; +\endcode +This version has the advantage to fully exploit the higher performance of dense storage (no indirection, SIMD, etc.), and to pay the cost of slow sparse evaluation on the few non-zeros of the sparse matrix only. %Sparse expressions also support transposition: @@ -304,6 +311,26 @@ sm2 = sm1.transpose() * P; \endcode +\subsection TutorialSparse_SubMatrices Block operations + +Regarding read-access, sparse matrices expose the same API than for dense matrices to access to sub-matrices such as blocks, columns, and rows. See \ref TutorialBlockOperations for a detailed introduction. +However, for performance reasons, writing to a sub-sparse-matrix is much more limited, and currently only contiguous sets of columns (resp. rows) of a column-major (resp. row-major) SparseMatrix are writable. Moreover, this information has to be known at compile-time, leaving out methods such as <tt>block(...)</tt> and <tt>corner*(...)</tt>. The available API for write-access to a SparseMatrix are summarized below: +\code +SparseMatrix<double,ColMajor> sm1; +sm1.col(j) = ...; +sm1.leftCols(ncols) = ...; +sm1.middleCols(j,ncols) = ...; +sm1.rightCols(ncols) = ...; + +SparseMatrix<double,RowMajor> sm2; +sm2.row(i) = ...; +sm2.topRows(nrows) = ...; +sm2.middleRows(i,nrows) = ...; +sm2.bottomRows(nrows) = ...; +\endcode + +In addition, sparse matrices expose the SparseMatrixBase::innerVector() and SparseMatrixBase::innerVectors() methods, which are aliases to the col/middleCols methods for a column-major storage, and to the row/middleRows methods for a row-major storage. + \subsection TutorialSparse_TriangularSelfadjoint Triangular and selfadjoint views Just as with dense matrices, the triangularView() function can be used to address a triangular part of the matrix, and perform triangular solves with a dense right hand side: diff --git a/doc/UnalignedArrayAssert.dox b/doc/UnalignedArrayAssert.dox index 8c97d7874..0f7022973 100644 --- a/doc/UnalignedArrayAssert.dox +++ b/doc/UnalignedArrayAssert.dox @@ -7,8 +7,8 @@ Hello! You are seeing this webpage because your program terminated on an asserti my_program: path/to/eigen/Eigen/src/Core/DenseStorage.h:44: Eigen::internal::matrix_array<T, Size, MatrixOptions, Align>::internal::matrix_array() [with T = double, int Size = 2, int MatrixOptions = 2, bool Align = true]: -Assertion `(reinterpret_cast<size_t>(array) & 0xf) == 0 && "this assertion -is explained here: http://eigen.tuxfamily.org/dox/UnalignedArrayAssert.html +Assertion `(reinterpret_cast<size_t>(array) & (sizemask)) == 0 && "this assertion +is explained here: http://eigen.tuxfamily.org/dox-devel/group__TopicUnalignedArrayAssert.html **** READ THIS WEB PAGE !!! ****"' failed. </pre> @@ -46,9 +46,9 @@ then you need to read this separate page: \ref TopicStructHavingEigenMembers "St Note that here, Eigen::Vector2d is only used as an example, more generally the issue arises for all \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen types". -\section c2 Cause 2: STL Containers +\section c2 Cause 2: STL Containers or manual memory allocation -If you use STL Containers such as std::vector, std::map, ..., with Eigen objects, or with classes containing Eigen objects, like this, +If you use STL Containers such as std::vector, std::map, ..., with %Eigen objects, or with classes containing %Eigen objects, like this, \code std::vector<Eigen::Matrix2f> my_vector; @@ -60,6 +60,8 @@ then you need to read this separate page: \ref TopicStlContainers "Using STL Con Note that here, Eigen::Matrix2f is only used as an example, more generally the issue arises for all \ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen types" and \ref TopicStructHavingEigenMembers "structures having such Eigen objects as member". +The same issue will be exhibited by any classes/functions by-passing operator new to allocate memory, that is, by performing custom memory allocation followed by calls to the placement new operator. This is for instance typically the case of \c std::make_shared or \c std::allocate_shared for which is the solution is to use an \ref aligned_allocator "aligned allocator" as detailed in the \ref TopicStlContainers "solution for STL containers". + \section c3 Cause 3: Passing Eigen objects by value If some function in your code is getting an Eigen object passed by value, like this, @@ -90,24 +92,37 @@ Note that here, Eigen::Quaternionf is only used as an example, more generally th \section explanation General explanation of this assertion -\ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen objects" must absolutely be created at 16-byte-aligned locations, otherwise SIMD instructions adressing them will crash. +\ref TopicFixedSizeVectorizable "fixed-size vectorizable Eigen objects" must absolutely be created at 16-byte-aligned locations, otherwise SIMD instructions addressing them will crash. Eigen normally takes care of these alignment issues for you, by setting an alignment attribute on them and by overloading their "operator new". However there are a few corner cases where these alignment settings get overridden: they are the possible causes for this assertion. -\section getrid I don't care about vectorization, how do I get rid of that stuff? +\section getrid I don't care about optimal vectorization, how do I get rid of that stuff? -Two possibilities: +Three possibilities: <ul> - <li>Define EIGEN_DONT_ALIGN_STATICALLY. That disables all 128-bit static alignment code, while keeping 128-bit heap alignment. This has the effect of - disabling vectorization for fixed-size objects (like Matrix4d) while keeping vectorization of dynamic-size objects - (like MatrixXd). But do note that this breaks ABI compatibility with the default behavior of 128-bit static alignment.</li> - <li>Or define both EIGEN_DONT_VECTORIZE and EIGEN_DISABLE_UNALIGNED_ARRAY_ASSERT. This keeps the - 128-bit alignment code and thus preserves ABI compatibility, but completely disables vectorization.</li> + <li>Use the \c DontAlign option to Matrix, Array, Quaternion, etc. objects that gives you trouble. This way Eigen won't try to align them, and thus won"t assume any special alignment. On the down side, you will pay the cost of unaligned loads/stores for them, but on modern CPUs, the overhead is either null or marginal. See \link StructHavingEigenMembers_othersolutions here \endlink for an example.</li> + <li>Define \link TopicPreprocessorDirectivesPerformance EIGEN_DONT_ALIGN_STATICALLY \endlink. That disables all 16-byte (and above) static alignment code, while keeping 16-byte (or above) heap alignment. This has the effect of + vectorizing fixed-size objects (like Matrix4d) through unaligned stores (as controlled by \link TopicPreprocessorDirectivesPerformance EIGEN_UNALIGNED_VECTORIZE \endlink), while keeping unchanged the vectorization of dynamic-size objects + (like MatrixXd). But do note that this breaks ABI compatibility with the default behavior of static alignment.</li> + <li>Or define both \link TopicPreprocessorDirectivesPerformance EIGEN_DONT_VECTORIZE \endlink and EIGEN_DISABLE_UNALIGNED_ARRAY_ASSERT. This keeps the + 16-byte alignment code and thus preserves ABI compatibility, but completely disables vectorization.</li> </ul> -For more information, see <a href="http://eigen.tuxfamily.org/index.php?title=FAQ#I_disabled_vectorization.2C_but_I.27m_still_getting_annoyed_about_alignment_issues.21">this FAQ</a>. +If you want to know why defining EIGEN_DONT_VECTORIZE does not by itself disable 16-byte alignment and the assertion, here's the explanation: + +It doesn't disable the assertion, because otherwise code that runs fine without vectorization would suddenly crash when enabling vectorization. +It doesn't disable 16-byte alignment, because that would mean that vectorized and non-vectorized code are not mutually ABI-compatible. This ABI compatibility is very important, even for people who develop only an in-house application, as for instance one may want to have in the same application a vectorized path and a non-vectorized path. + +\section checkmycode How can I check my code is safe regarding alignment issues? + +Unfortunately, there is no possibility in C++ to detect any of the aformentioned shortcoming at compile time (though static analysers are becoming more and more powerful and could detect some of them). +Even at runtime, all we can do is to catch invalid unaligned allocation and trigger the explicit assertion mentioned at the begining of this page. +Therefore, if your program runs fine on a given system with some given compilation flags, then this does not guarantee that your code is safe. For instance, on most 64 bits systems buffer are aligned on 16 bytes boundary and so, if you do not enable AVX instruction set, then your code will run fine. On the other hand, the same code may assert if moving to a more exotic platform, or enabling AVX instructions that required 32 bytes alignment by default. + +The situation is not hopeless though. Assuming your code is well covered by unit test, then you can check its alignment safety by linking it to a custom malloc library returning 8 bytes aligned buffers only. This way all alignment shortcomings should pop-up. To this end, you must also compile your program with \link TopicPreprocessorDirectivesPerformance EIGEN_MALLOC_ALREADY_ALIGNED=0 \endlink. + */ diff --git a/doc/UsingBlasLapackBackends.dox b/doc/UsingBlasLapackBackends.dox new file mode 100644 index 000000000..caa597122 --- /dev/null +++ b/doc/UsingBlasLapackBackends.dox @@ -0,0 +1,133 @@ +/* + Copyright (c) 2011, Intel Corporation. All rights reserved. + Copyright (C) 2011-2016 Gael Guennebaud <gael.guennebaud@inria.fr> + + Redistribution and use in source and binary forms, with or without modification, + are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright notice, this + list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright notice, + this list of conditions and the following disclaimer in the documentation + and/or other materials provided with the distribution. + * Neither the name of Intel Corporation nor the names of its contributors may + be used to endorse or promote products derived from this software without + specific prior written permission. + + THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND + ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED + WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE + DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR + ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES + (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; + LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON + ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT + (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS + SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + + ******************************************************************************** + * Content : Documentation on the use of BLAS/LAPACK libraries through Eigen + ******************************************************************************** +*/ + +namespace Eigen { + +/** \page TopicUsingBlasLapack Using BLAS/LAPACK from %Eigen + + +Since %Eigen version 3.3 and later, any F77 compatible BLAS or LAPACK libraries can be used as backends for dense matrix products and dense matrix decompositions. +For instance, one can use <a href="http://eigen.tuxfamily.org/Counter/redirect_to_mkl.php">Intel® MKL</a>, Apple's Accelerate framework on OSX, <a href="http://www.openblas.net/">OpenBLAS</a>, <a href="http://www.netlib.org/lapack">Netlib LAPACK</a>, etc. + +Do not miss this \link TopicUsingIntelMKL page \endlink for further discussions on the specific use of Intel® MKL (also includes VML, PARDISO, etc.) + +In order to use an external BLAS and/or LAPACK library, you must link you own application to the respective libraries and their dependencies. +For LAPACK, you must also link to the standard <a href="http://www.netlib.org/lapack/lapacke.html">Lapacke</a> library, which is used as a convenient think layer between %Eigen's C++ code and LAPACK F77 interface. Then you must activate their usage by defining one or multiple of the following macros (\b before including any %Eigen's header): + +\note For Mac users, in order to use the lapack version shipped with the Accelerate framework, you also need the lapacke library. +Using <a href="https://www.macports.org/">MacPorts</a>, this is as easy as: +\code +sudo port install lapack +\endcode +and then use the following link flags: \c -framework \c Accelerate \c /opt/local/lib/lapack/liblapacke.dylib + +<table class="manual"> +<tr><td>\c EIGEN_USE_BLAS </td><td>Enables the use of external BLAS level 2 and 3 routines (compatible with any F77 BLAS interface)</td></tr> +<tr class="alt"><td>\c EIGEN_USE_LAPACKE </td><td>Enables the use of external Lapack routines via the <a href="http://www.netlib.org/lapack/lapacke.html">Lapacke</a> C interface to Lapack (compatible with any F77 LAPACK interface)</td></tr> +<tr><td>\c EIGEN_USE_LAPACKE_STRICT </td><td>Same as \c EIGEN_USE_LAPACKE but algorithms of lower numerical robustness are disabled. \n This currently concerns only JacobiSVD which otherwise would be replaced by \c gesvd that is less robust than Jacobi rotations.</td></tr> +</table> + +When doing so, a number of %Eigen's algorithms are silently substituted with calls to BLAS or LAPACK routines. +These substitutions apply only for \b Dynamic \b or \b large enough objects with one of the following four standard scalar types: \c float, \c double, \c complex<float>, and \c complex<double>. +Operations on other scalar types or mixing reals and complexes will continue to use the built-in algorithms. + +The breadth of %Eigen functionality that can be substituted is listed in the table below. +<table class="manual"> +<tr><th>Functional domain</th><th>Code example</th><th>BLAS/LAPACK routines</th></tr> +<tr><td>Matrix-matrix operations \n \c EIGEN_USE_BLAS </td><td>\code +m1*m2.transpose(); +m1.selfadjointView<Lower>()*m2; +m1*m2.triangularView<Upper>(); +m1.selfadjointView<Lower>().rankUpdate(m2,1.0); +\endcode</td><td>\code +?gemm +?symm/?hemm +?trmm +dsyrk/ssyrk +\endcode</td></tr> +<tr class="alt"><td>Matrix-vector operations \n \c EIGEN_USE_BLAS </td><td>\code +m1.adjoint()*b; +m1.selfadjointView<Lower>()*b; +m1.triangularView<Upper>()*b; +\endcode</td><td>\code +?gemv +?symv/?hemv +?trmv +\endcode</td></tr> +<tr><td>LU decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code +v1 = m1.lu().solve(v2); +\endcode</td><td>\code +?getrf +\endcode</td></tr> +<tr class="alt"><td>Cholesky decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code +v1 = m2.selfadjointView<Upper>().llt().solve(v2); +\endcode</td><td>\code +?potrf +\endcode</td></tr> +<tr><td>QR decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code +m1.householderQr(); +m1.colPivHouseholderQr(); +\endcode</td><td>\code +?geqrf +?geqp3 +\endcode</td></tr> +<tr class="alt"><td>Singular value decomposition \n \c EIGEN_USE_LAPACKE </td><td>\code +JacobiSVD<MatrixXd> svd; +svd.compute(m1, ComputeThinV); +\endcode</td><td>\code +?gesvd +\endcode</td></tr> +<tr><td>Eigen-value decompositions \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code +EigenSolver<MatrixXd> es(m1); +ComplexEigenSolver<MatrixXcd> ces(m1); +SelfAdjointEigenSolver<MatrixXd> saes(m1+m1.transpose()); +GeneralizedSelfAdjointEigenSolver<MatrixXd> + gsaes(m1+m1.transpose(),m2+m2.transpose()); +\endcode</td><td>\code +?gees +?gees +?syev/?heev +?syev/?heev, +?potrf +\endcode</td></tr> +<tr class="alt"><td>Schur decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code +RealSchur<MatrixXd> schurR(m1); +ComplexSchur<MatrixXcd> schurC(m1); +\endcode</td><td>\code +?gees +\endcode</td></tr> +</table> +In the examples, m1 and m2 are dense matrices and v1 and v2 are dense vectors. + +*/ + +} diff --git a/doc/UsingIntelMKL.dox b/doc/UsingIntelMKL.dox index 4b624a156..a1a3a18f2 100644 --- a/doc/UsingIntelMKL.dox +++ b/doc/UsingIntelMKL.dox @@ -32,106 +32,45 @@ namespace Eigen { -/** \page TopicUsingIntelMKL Using Intel® Math Kernel Library from Eigen +/** \page TopicUsingIntelMKL Using Intel® MKL from %Eigen -\section TopicUsingIntelMKL_Intro Eigen and Intel® Math Kernel Library (Intel® MKL) +<!-- \section TopicUsingIntelMKL_Intro Eigen and Intel® Math Kernel Library (Intel® MKL) --> + +Since %Eigen version 3.1 and later, users can benefit from built-in Intel® Math Kernel Library (MKL) optimizations with an installed copy of Intel MKL 10.3 (or later). -Since Eigen version 3.1 and later, users can benefit from built-in Intel MKL optimizations with an installed copy of Intel MKL 10.3 (or later). <a href="http://eigen.tuxfamily.org/Counter/redirect_to_mkl.php"> Intel MKL </a> provides highly optimized multi-threaded mathematical routines for x86-compatible architectures. Intel MKL is available on Linux, Mac and Windows for both Intel64 and IA32 architectures. -\warning Be aware that Intel® MKL is a proprietary software. It is the responsibility of the users to buy MKL licenses for their products. Moreover, the license of the user product has to allow linking to proprietary software that excludes any unmodified versions of the GPL. +\note +Intel® MKL is a proprietary software and it is the responsibility of users to buy or register for community (free) Intel MKL licenses for their products. Moreover, the license of the user product has to allow linking to proprietary software that excludes any unmodified versions of the GPL. -Using Intel MKL through Eigen is easy: --# define the \c EIGEN_USE_MKL_ALL macro before including any Eigen's header +Using Intel MKL through %Eigen is easy: +-# define the \c EIGEN_USE_MKL_ALL macro before including any %Eigen's header -# link your program to MKL libraries (see the <a href="http://software.intel.com/en-us/articles/intel-mkl-link-line-advisor/">MKL linking advisor</a>) -# on a 64bits system, you must use the LP64 interface (not the ILP64 one) -When doing so, a number of Eigen's algorithms are silently substituted with calls to Intel MKL routines. +When doing so, a number of %Eigen's algorithms are silently substituted with calls to Intel MKL routines. These substitutions apply only for \b Dynamic \b or \b large enough objects with one of the following four standard scalar types: \c float, \c double, \c complex<float>, and \c complex<double>. Operations on other scalar types or mixing reals and complexes will continue to use the built-in algorithms. -In addition you can coarsely select choose which parts will be substituted by defining one or multiple of the following macros: +In addition you can choose which parts will be substituted by defining one or multiple of the following macros: <table class="manual"> -<tr><td>\c EIGEN_USE_BLAS </td><td>Enables the use of external BLAS level 2 and 3 routines (currently works with Intel MKL only)</td></tr> -<tr class="alt"><td>\c EIGEN_USE_LAPACKE </td><td>Enables the use of external Lapack routines via the <a href="http://www.netlib.org/lapack/lapacke.html">Intel Lapacke</a> C interface to Lapack (currently works with Intel MKL only)</td></tr> -<tr><td>\c EIGEN_USE_LAPACKE_STRICT </td><td>Same as \c EIGEN_USE_LAPACKE but algorithm of lower robustness are disabled. This currently concerns only JacobiSVD which otherwise would be replaced by \c gesvd that is less robust than Jacobi rotations.</td></tr> +<tr><td>\c EIGEN_USE_BLAS </td><td>Enables the use of external BLAS level 2 and 3 routines</td></tr> +<tr class="alt"><td>\c EIGEN_USE_LAPACKE </td><td>Enables the use of external Lapack routines via the <a href="http://www.netlib.org/lapack/lapacke.html">Lapacke</a> C interface to Lapack</td></tr> +<tr><td>\c EIGEN_USE_LAPACKE_STRICT </td><td>Same as \c EIGEN_USE_LAPACKE but algorithm of lower robustness are disabled. \n This currently concerns only JacobiSVD which otherwise would be replaced by \c gesvd that is less robust than Jacobi rotations.</td></tr> <tr class="alt"><td>\c EIGEN_USE_MKL_VML </td><td>Enables the use of Intel VML (vector operations)</td></tr> <tr><td>\c EIGEN_USE_MKL_ALL </td><td>Defines \c EIGEN_USE_BLAS, \c EIGEN_USE_LAPACKE, and \c EIGEN_USE_MKL_VML </td></tr> </table> -Finally, the PARDISO sparse solver shipped with Intel MKL can be used through the \ref PardisoLU, \ref PardisoLLT and \ref PardisoLDLT classes of the \ref PardisoSupport_Module. - +Note that the BLAS and LAPACKE backends can be enabled for any F77 compatible BLAS and LAPACK libraries. See this \link TopicUsingBlasLapack page \endlink for the details. -\section TopicUsingIntelMKL_SupportedFeatures List of supported features +Finally, the PARDISO sparse solver shipped with Intel MKL can be used through the \ref PardisoLU, \ref PardisoLLT and \ref PardisoLDLT classes of the \ref PardisoSupport_Module. -The breadth of Eigen functionality covered by Intel MKL is listed in the table below. +The following table summarizes the list of functions covered by \c EIGEN_USE_MKL_VML: <table class="manual"> -<tr><th>Functional domain</th><th>Code example</th><th>MKL routines</th></tr> -<tr><td>Matrix-matrix operations \n \c EIGEN_USE_BLAS </td><td>\code -m1*m2.transpose(); -m1.selfadjointView<Lower>()*m2; -m1*m2.triangularView<Upper>(); -m1.selfadjointView<Lower>().rankUpdate(m2,1.0); -\endcode</td><td>\code -?gemm -?symm/?hemm -?trmm -dsyrk/ssyrk -\endcode</td></tr> -<tr class="alt"><td>Matrix-vector operations \n \c EIGEN_USE_BLAS </td><td>\code -m1.adjoint()*b; -m1.selfadjointView<Lower>()*b; -m1.triangularView<Upper>()*b; -\endcode</td><td>\code -?gemv -?symv/?hemv -?trmv -\endcode</td></tr> -<tr><td>LU decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code -v1 = m1.lu().solve(v2); -\endcode</td><td>\code -?getrf -\endcode</td></tr> -<tr class="alt"><td>Cholesky decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code -v1 = m2.selfadjointView<Upper>().llt().solve(v2); -\endcode</td><td>\code -?potrf -\endcode</td></tr> -<tr><td>QR decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code -m1.householderQr(); -m1.colPivHouseholderQr(); -\endcode</td><td>\code -?geqrf -?geqp3 -\endcode</td></tr> -<tr class="alt"><td>Singular value decomposition \n \c EIGEN_USE_LAPACKE </td><td>\code -JacobiSVD<MatrixXd> svd; -svd.compute(m1, ComputeThinV); -\endcode</td><td>\code -?gesvd -\endcode</td></tr> -<tr><td>Eigen-value decompositions \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code -EigenSolver<MatrixXd> es(m1); -ComplexEigenSolver<MatrixXcd> ces(m1); -SelfAdjointEigenSolver<MatrixXd> saes(m1+m1.transpose()); -GeneralizedSelfAdjointEigenSolver<MatrixXd> - gsaes(m1+m1.transpose(),m2+m2.transpose()); -\endcode</td><td>\code -?gees -?gees -?syev/?heev -?syev/?heev, -?potrf -\endcode</td></tr> -<tr class="alt"><td>Schur decomposition \n \c EIGEN_USE_LAPACKE \n \c EIGEN_USE_LAPACKE_STRICT </td><td>\code -RealSchur<MatrixXd> schurR(m1); -ComplexSchur<MatrixXcd> schurC(m1); -\endcode</td><td>\code -?gees -\endcode</td></tr> -<tr><td>Vector Math \n \c EIGEN_USE_MKL_VML </td><td>\code +<tr><th>Code example</th><th>MKL routines</th></tr> +<tr><td>\code v2=v1.array().sin(); v2=v1.array().asin(); v2=v1.array().cos(); @@ -155,7 +94,7 @@ v?Sqr v?Powx \endcode</td></tr> </table> -In the examples, m1 and m2 are dense matrices and v1 and v2 are dense vectors. +In the examples, v1 and v2 are dense vectors. \section TopicUsingIntelMKL_Links Links diff --git a/doc/UsingNVCC.dox b/doc/UsingNVCC.dox new file mode 100644 index 000000000..f8e755b79 --- /dev/null +++ b/doc/UsingNVCC.dox @@ -0,0 +1,32 @@ + +namespace Eigen { + +/** \page TopicCUDA Using Eigen in CUDA kernels + +\b Disclaimer: this page is about an \b experimental feature in %Eigen. + +Staring from CUDA 5.0, the CUDA compiler, \c nvcc, is able to properly parse %Eigen's code (almost). +A few adaptations of the %Eigen's code already allows to use some parts of %Eigen in your own CUDA kernels. +To this end you need the devel branch of %Eigen, CUDA 5.0 or greater with GCC. + +Known issues: + + - \c nvcc with MS Visual Studio does not work (patch welcome) + + - \c nvcc with \c clang does not work (patch welcome) + + - \c nvcc 5.5 with gcc-4.7 (or greater) has issues with the standard \c \<limits\> header file. To workaround this, you can add the following before including any other files: + \code + // workaround issue between gcc >= 4.7 and cuda 5.5 + #if (defined __GNUC__) && (__GNUC__>4 || __GNUC_MINOR__>=7) + #undef _GLIBCXX_ATOMIC_BUILTINS + #undef _GLIBCXX_USE_INT128 + #endif + \endcode + + - On 64bits system Eigen uses \c long \c int as the default type for indexes and sizes. On CUDA device, it would make sense to default to 32 bits \c int. + However, to keep host and CUDA code compatible, this cannot be done automatically by %Eigen, and the user is thus required to define \c EIGEN_DEFAULT_DENSE_INDEX_TYPE to \c int throughout his code (or only for CUDA code if there is no interaction between host and CUDA code through %Eigen's object). + +*/ + +} diff --git a/doc/eigendoxy.css b/doc/eigendoxy.css index 60243d870..6ce2b839b 100644 --- a/doc/eigendoxy.css +++ b/doc/eigendoxy.css @@ -45,7 +45,7 @@ pre.fragment { /* Common style for all Eigen's tables */ -table.example, table.manual, table.manual-vl { +table.example, table.manual, table.manual-vl, table.manual-hl { max-width:100%; border-collapse: collapse; border-style: solid; @@ -58,7 +58,7 @@ table.example, table.manual, table.manual-vl { -webkit-box-shadow: 5px 5px 5px rgba(0, 0, 0, 0.15); } -table.example th, table.manual th, table.manual-vl th { +table.example th, table.manual th, table.manual-vl th, table.manual-hl th { padding: 0.5em 0.5em 0.5em 0.5em; text-align: left; padding-right: 1em; @@ -70,7 +70,7 @@ table.example th, table.manual th, table.manual-vl th { filter: progid:DXImageTransform.Microsoft.gradient(startColorstr='#FFFFFF', endColorstr='#F4F4E5'); } -table.example td, table.manual td, table.manual-vl td { +table.example td, table.manual td, table.manual-vl td, table.manual-hl td { vertical-align:top; border-width: 1px; border-color: #cccccc; @@ -108,15 +108,15 @@ table.example td { /* standard class for the manual */ -table.manual, table.manual-vl { +table.manual, table.manual-vl, table.manual-hl { padding: 0.2em 0em 0.5em 0em; } -table.manual th, table.manual-vl th { +table.manual th, table.manual-vl th, table.manual-hl th { margin: 0em 0em 0.3em 0em; } -table.manual td, table.manual-vl td { +table.manual td, table.manual-vl td, table.manual-hl td { padding: 0.3em 0.5em 0.3em 0.5em; vertical-align:top; border-width: 1px; @@ -136,6 +136,16 @@ table.manual-vl th.inter { border-style: solid solid solid solid; } +table.manual-hl td { + border-color: #cccccc; + border-width: 1px; + border-style: solid none solid none; +} + +table td.code { + font-family: monospace; +} + h2 { margin-top:2em; border-style: none none solid none; @@ -166,6 +176,16 @@ div.toc ul { margin: 0.2em 0 0.4em 0.5em; } +span.cpp11,span.cpp14,span.cpp17 { + color: #119911; + font-weight: bold; +} + +.newin3x { + color: #a37c1a; + font-weight: bold; +} + /**** old Eigen's styles ****/ @@ -177,8 +197,8 @@ table.tutorial_code td { /* Whenever doxygen meets a '\n' or a '<BR/>', it will put - * the text containing the characted into a <p class="starttd">. - * This little hack togehter with table.tutorial_code td.note + * the text containing the character into a <p class="starttd">. + * This little hack together with table.tutorial_code td.note * aims at fixing this issue. */ table.tutorial_code td.note p.starttd { margin: 0px; diff --git a/doc/examples/CMakeLists.txt b/doc/examples/CMakeLists.txt index 71b272a15..f7a19055f 100644 --- a/doc/examples/CMakeLists.txt +++ b/doc/examples/CMakeLists.txt @@ -6,13 +6,16 @@ foreach(example_src ${examples_SRCS}) if(EIGEN_STANDARD_LIBRARIES_TO_LINK_TO) target_link_libraries(${example} ${EIGEN_STANDARD_LIBRARIES_TO_LINK_TO}) endif() - get_target_property(example_executable - ${example} LOCATION) add_custom_command( TARGET ${example} POST_BUILD - COMMAND ${example_executable} + COMMAND ${example} ARGS >${CMAKE_CURRENT_BINARY_DIR}/${example}.out ) add_dependencies(all_examples ${example}) endforeach(example_src) + +check_cxx_compiler_flag("-std=c++11" EIGEN_COMPILER_SUPPORT_CPP11) +if(EIGEN_COMPILER_SUPPORT_CPP11) +ei_add_target_property(nullary_indexing COMPILE_FLAGS "-std=c++11") +endif()
\ No newline at end of file diff --git a/doc/examples/Cwise_erf.cpp b/doc/examples/Cwise_erf.cpp new file mode 100644 index 000000000..e7cd2c1c0 --- /dev/null +++ b/doc/examples/Cwise_erf.cpp @@ -0,0 +1,9 @@ +#include <Eigen/Core> +#include <unsupported/Eigen/SpecialFunctions> +#include <iostream> +using namespace Eigen; +int main() +{ + Array4d v(-0.5,2,0,-7); + std::cout << v.erf() << std::endl; +} diff --git a/doc/examples/Cwise_erfc.cpp b/doc/examples/Cwise_erfc.cpp new file mode 100644 index 000000000..d8bb04c30 --- /dev/null +++ b/doc/examples/Cwise_erfc.cpp @@ -0,0 +1,9 @@ +#include <Eigen/Core> +#include <unsupported/Eigen/SpecialFunctions> +#include <iostream> +using namespace Eigen; +int main() +{ + Array4d v(-0.5,2,0,-7); + std::cout << v.erfc() << std::endl; +} diff --git a/doc/examples/Cwise_lgamma.cpp b/doc/examples/Cwise_lgamma.cpp new file mode 100644 index 000000000..f1c4f503e --- /dev/null +++ b/doc/examples/Cwise_lgamma.cpp @@ -0,0 +1,9 @@ +#include <Eigen/Core> +#include <unsupported/Eigen/SpecialFunctions> +#include <iostream> +using namespace Eigen; +int main() +{ + Array4d v(0.5,10,0,-1); + std::cout << v.lgamma() << std::endl; +}
\ No newline at end of file diff --git a/doc/examples/MatrixBase_cwise_const.cpp b/doc/examples/MatrixBase_cwise_const.cpp deleted file mode 100644 index 23700e0b6..000000000 --- a/doc/examples/MatrixBase_cwise_const.cpp +++ /dev/null @@ -1,18 +0,0 @@ -#define EIGEN2_SUPPORT -#include <Eigen/Core> -#include <iostream> - -using namespace Eigen; -using namespace std; - -int main() -{ - Matrix3i m = Matrix3i::Random(); - cout << "Here is the matrix m:" << endl << m << endl; - Matrix3i n = Matrix3i::Random(); - cout << "And here is the matrix n:" << endl << n << endl; - cout << "The coefficient-wise product of m and n is:" << endl; - cout << m.cwise() * n << endl; - cout << "Taking the cube of the coefficients of m yields:" << endl; - cout << m.cwise().pow(3) << endl; -} diff --git a/doc/examples/TutorialInplaceLU.cpp b/doc/examples/TutorialInplaceLU.cpp new file mode 100644 index 000000000..cb9c59b60 --- /dev/null +++ b/doc/examples/TutorialInplaceLU.cpp @@ -0,0 +1,61 @@ +#include <iostream> +struct init { + init() { std::cout << "[" << "init" << "]" << std::endl; } +}; +init init_obj; +// [init] +#include <iostream> +#include <Eigen/Dense> + +using namespace std; +using namespace Eigen; + +int main() +{ + MatrixXd A(2,2); + A << 2, -1, 1, 3; + cout << "Here is the input matrix A before decomposition:\n" << A << endl; +cout << "[init]" << endl; + +cout << "[declaration]" << endl; + PartialPivLU<Ref<MatrixXd> > lu(A); + cout << "Here is the input matrix A after decomposition:\n" << A << endl; +cout << "[declaration]" << endl; + +cout << "[matrixLU]" << endl; + cout << "Here is the matrix storing the L and U factors:\n" << lu.matrixLU() << endl; +cout << "[matrixLU]" << endl; + +cout << "[solve]" << endl; + MatrixXd A0(2,2); A0 << 2, -1, 1, 3; + VectorXd b(2); b << 1, 2; + VectorXd x = lu.solve(b); + cout << "Residual: " << (A0 * x - b).norm() << endl; +cout << "[solve]" << endl; + +cout << "[modifyA]" << endl; + A << 3, 4, -2, 1; + x = lu.solve(b); + cout << "Residual: " << (A0 * x - b).norm() << endl; +cout << "[modifyA]" << endl; + +cout << "[recompute]" << endl; + A0 = A; // save A + lu.compute(A); + x = lu.solve(b); + cout << "Residual: " << (A0 * x - b).norm() << endl; +cout << "[recompute]" << endl; + +cout << "[recompute_bis0]" << endl; + MatrixXd A1(2,2); + A1 << 5,-2,3,4; + lu.compute(A1); + cout << "Here is the input matrix A1 after decomposition:\n" << A1 << endl; +cout << "[recompute_bis0]" << endl; + +cout << "[recompute_bis1]" << endl; + x = lu.solve(b); + cout << "Residual: " << (A1 * x - b).norm() << endl; +cout << "[recompute_bis1]" << endl; + +} diff --git a/doc/examples/TutorialLinAlgInverseDeterminant.cpp b/doc/examples/TutorialLinAlgInverseDeterminant.cpp index 43970ff05..14dde5b35 100644 --- a/doc/examples/TutorialLinAlgInverseDeterminant.cpp +++ b/doc/examples/TutorialLinAlgInverseDeterminant.cpp @@ -13,4 +13,4 @@ int main() cout << "Here is the matrix A:\n" << A << endl; cout << "The determinant of A is " << A.determinant() << endl; cout << "The inverse of A is:\n" << A.inverse() << endl; -}
\ No newline at end of file +} diff --git a/doc/examples/Tutorial_ReductionsVisitorsBroadcasting_reductions_operatornorm.cpp b/doc/examples/Tutorial_ReductionsVisitorsBroadcasting_reductions_operatornorm.cpp new file mode 100644 index 000000000..62e28fc31 --- /dev/null +++ b/doc/examples/Tutorial_ReductionsVisitorsBroadcasting_reductions_operatornorm.cpp @@ -0,0 +1,18 @@ +#include <Eigen/Dense> +#include <iostream> + +using namespace Eigen; +using namespace std; + +int main() +{ + MatrixXf m(2,2); + m << 1,-2, + -3,4; + + cout << "1-norm(m) = " << m.cwiseAbs().colwise().sum().maxCoeff() + << " == " << m.colwise().lpNorm<1>().maxCoeff() << endl; + + cout << "infty-norm(m) = " << m.cwiseAbs().rowwise().sum().maxCoeff() + << " == " << m.rowwise().lpNorm<1>().maxCoeff() << endl; +} diff --git a/doc/examples/make_circulant.cpp b/doc/examples/make_circulant.cpp new file mode 100644 index 000000000..92e6aaa2b --- /dev/null +++ b/doc/examples/make_circulant.cpp @@ -0,0 +1,11 @@ +/* +This program is presented in several fragments in the doc page. +Every fragment is in its own file; this file simply combines them. +*/ + +#include "make_circulant.cpp.preamble" +#include "make_circulant.cpp.traits" +#include "make_circulant.cpp.expression" +#include "make_circulant.cpp.evaluator" +#include "make_circulant.cpp.entry" +#include "make_circulant.cpp.main" diff --git a/doc/examples/make_circulant.cpp.entry b/doc/examples/make_circulant.cpp.entry new file mode 100644 index 000000000..f9d2eb8a9 --- /dev/null +++ b/doc/examples/make_circulant.cpp.entry @@ -0,0 +1,5 @@ +template <class ArgType> +Circulant<ArgType> makeCirculant(const Eigen::MatrixBase<ArgType>& arg) +{ + return Circulant<ArgType>(arg.derived()); +} diff --git a/doc/examples/make_circulant.cpp.evaluator b/doc/examples/make_circulant.cpp.evaluator new file mode 100644 index 000000000..2ba79e783 --- /dev/null +++ b/doc/examples/make_circulant.cpp.evaluator @@ -0,0 +1,32 @@ +namespace Eigen { + namespace internal { + template<typename ArgType> + struct evaluator<Circulant<ArgType> > + : evaluator_base<Circulant<ArgType> > + { + typedef Circulant<ArgType> XprType; + typedef typename nested_eval<ArgType, XprType::ColsAtCompileTime>::type ArgTypeNested; + typedef typename remove_all<ArgTypeNested>::type ArgTypeNestedCleaned; + typedef typename XprType::CoeffReturnType CoeffReturnType; + + enum { + CoeffReadCost = evaluator<ArgTypeNestedCleaned>::CoeffReadCost, + Flags = Eigen::ColMajor + }; + + evaluator(const XprType& xpr) + : m_argImpl(xpr.m_arg), m_rows(xpr.rows()) + { } + + CoeffReturnType coeff(Index row, Index col) const + { + Index index = row - col; + if (index < 0) index += m_rows; + return m_argImpl.coeff(index); + } + + evaluator<ArgTypeNestedCleaned> m_argImpl; + const Index m_rows; + }; + } +} diff --git a/doc/examples/make_circulant.cpp.expression b/doc/examples/make_circulant.cpp.expression new file mode 100644 index 000000000..380cd4450 --- /dev/null +++ b/doc/examples/make_circulant.cpp.expression @@ -0,0 +1,20 @@ +template <class ArgType> +class Circulant : public Eigen::MatrixBase<Circulant<ArgType> > +{ +public: + Circulant(const ArgType& arg) + : m_arg(arg) + { + EIGEN_STATIC_ASSERT(ArgType::ColsAtCompileTime == 1, + YOU_TRIED_CALLING_A_VECTOR_METHOD_ON_A_MATRIX); + } + + typedef typename Eigen::internal::ref_selector<Circulant>::type Nested; + + typedef Eigen::Index Index; + Index rows() const { return m_arg.rows(); } + Index cols() const { return m_arg.rows(); } + + typedef typename Eigen::internal::ref_selector<ArgType>::type ArgTypeNested; + ArgTypeNested m_arg; +}; diff --git a/doc/examples/make_circulant.cpp.main b/doc/examples/make_circulant.cpp.main new file mode 100644 index 000000000..877f97f62 --- /dev/null +++ b/doc/examples/make_circulant.cpp.main @@ -0,0 +1,8 @@ +int main() +{ + Eigen::VectorXd vec(4); + vec << 1, 2, 4, 8; + Eigen::MatrixXd mat; + mat = makeCirculant(vec); + std::cout << mat << std::endl; +} diff --git a/doc/examples/make_circulant.cpp.preamble b/doc/examples/make_circulant.cpp.preamble new file mode 100644 index 000000000..e575cce14 --- /dev/null +++ b/doc/examples/make_circulant.cpp.preamble @@ -0,0 +1,4 @@ +#include <Eigen/Core> +#include <iostream> + +template <class ArgType> class Circulant; diff --git a/doc/examples/make_circulant.cpp.traits b/doc/examples/make_circulant.cpp.traits new file mode 100644 index 000000000..4e04535d3 --- /dev/null +++ b/doc/examples/make_circulant.cpp.traits @@ -0,0 +1,19 @@ +namespace Eigen { + namespace internal { + template <class ArgType> + struct traits<Circulant<ArgType> > + { + typedef Eigen::Dense StorageKind; + typedef Eigen::MatrixXpr XprKind; + typedef typename ArgType::StorageIndex StorageIndex; + typedef typename ArgType::Scalar Scalar; + enum { + Flags = Eigen::ColMajor, + RowsAtCompileTime = ArgType::RowsAtCompileTime, + ColsAtCompileTime = ArgType::RowsAtCompileTime, + MaxRowsAtCompileTime = ArgType::MaxRowsAtCompileTime, + MaxColsAtCompileTime = ArgType::MaxRowsAtCompileTime + }; + }; + } +} diff --git a/doc/examples/make_circulant2.cpp b/doc/examples/make_circulant2.cpp new file mode 100644 index 000000000..95d3dd31a --- /dev/null +++ b/doc/examples/make_circulant2.cpp @@ -0,0 +1,52 @@ +#include <Eigen/Core> +#include <iostream> + +using namespace Eigen; + +// [circulant_func] +template<class ArgType> +class circulant_functor { + const ArgType &m_vec; +public: + circulant_functor(const ArgType& arg) : m_vec(arg) {} + + const typename ArgType::Scalar& operator() (Index row, Index col) const { + Index index = row - col; + if (index < 0) index += m_vec.size(); + return m_vec(index); + } +}; +// [circulant_func] + +// [square] +template<class ArgType> +struct circulant_helper { + typedef Matrix<typename ArgType::Scalar, + ArgType::SizeAtCompileTime, + ArgType::SizeAtCompileTime, + ColMajor, + ArgType::MaxSizeAtCompileTime, + ArgType::MaxSizeAtCompileTime> MatrixType; +}; +// [square] + +// [makeCirculant] +template <class ArgType> +CwiseNullaryOp<circulant_functor<ArgType>, typename circulant_helper<ArgType>::MatrixType> +makeCirculant(const Eigen::MatrixBase<ArgType>& arg) +{ + typedef typename circulant_helper<ArgType>::MatrixType MatrixType; + return MatrixType::NullaryExpr(arg.size(), arg.size(), circulant_functor<ArgType>(arg.derived())); +} +// [makeCirculant] + +// [main] +int main() +{ + Eigen::VectorXd vec(4); + vec << 1, 2, 4, 8; + Eigen::MatrixXd mat; + mat = makeCirculant(vec); + std::cout << mat << std::endl; +} +// [main] diff --git a/doc/examples/matrixfree_cg.cpp b/doc/examples/matrixfree_cg.cpp new file mode 100644 index 000000000..6a205aea3 --- /dev/null +++ b/doc/examples/matrixfree_cg.cpp @@ -0,0 +1,128 @@ +#include <iostream> +#include <Eigen/Core> +#include <Eigen/Dense> +#include <Eigen/IterativeLinearSolvers> +#include <unsupported/Eigen/IterativeSolvers> + +class MatrixReplacement; +using Eigen::SparseMatrix; + +namespace Eigen { +namespace internal { + // MatrixReplacement looks-like a SparseMatrix, so let's inherits its traits: + template<> + struct traits<MatrixReplacement> : public Eigen::internal::traits<Eigen::SparseMatrix<double> > + {}; +} +} + +// Example of a matrix-free wrapper from a user type to Eigen's compatible type +// For the sake of simplicity, this example simply wrap a Eigen::SparseMatrix. +class MatrixReplacement : public Eigen::EigenBase<MatrixReplacement> { +public: + // Required typedefs, constants, and method: + typedef double Scalar; + typedef double RealScalar; + typedef int StorageIndex; + enum { + ColsAtCompileTime = Eigen::Dynamic, + MaxColsAtCompileTime = Eigen::Dynamic, + IsRowMajor = false + }; + + Index rows() const { return mp_mat->rows(); } + Index cols() const { return mp_mat->cols(); } + + template<typename Rhs> + Eigen::Product<MatrixReplacement,Rhs,Eigen::AliasFreeProduct> operator*(const Eigen::MatrixBase<Rhs>& x) const { + return Eigen::Product<MatrixReplacement,Rhs,Eigen::AliasFreeProduct>(*this, x.derived()); + } + + // Custom API: + MatrixReplacement() : mp_mat(0) {} + + void attachMyMatrix(const SparseMatrix<double> &mat) { + mp_mat = &mat; + } + const SparseMatrix<double> my_matrix() const { return *mp_mat; } + +private: + const SparseMatrix<double> *mp_mat; +}; + + +// Implementation of MatrixReplacement * Eigen::DenseVector though a specialization of internal::generic_product_impl: +namespace Eigen { +namespace internal { + + template<typename Rhs> + struct generic_product_impl<MatrixReplacement, Rhs, SparseShape, DenseShape, GemvProduct> // GEMV stands for matrix-vector + : generic_product_impl_base<MatrixReplacement,Rhs,generic_product_impl<MatrixReplacement,Rhs> > + { + typedef typename Product<MatrixReplacement,Rhs>::Scalar Scalar; + + template<typename Dest> + static void scaleAndAddTo(Dest& dst, const MatrixReplacement& lhs, const Rhs& rhs, const Scalar& alpha) + { + // This method should implement "dst += alpha * lhs * rhs" inplace, + // however, for iterative solvers, alpha is always equal to 1, so let's not bother about it. + assert(alpha==Scalar(1) && "scaling is not implemented"); + + // Here we could simply call dst.noalias() += lhs.my_matrix() * rhs, + // but let's do something fancier (and less efficient): + for(Index i=0; i<lhs.cols(); ++i) + dst += rhs(i) * lhs.my_matrix().col(i); + } + }; + +} +} + +int main() +{ + int n = 10; + Eigen::SparseMatrix<double> S = Eigen::MatrixXd::Random(n,n).sparseView(0.5,1); + S = S.transpose()*S; + + MatrixReplacement A; + A.attachMyMatrix(S); + + Eigen::VectorXd b(n), x; + b.setRandom(); + + // Solve Ax = b using various iterative solver with matrix-free version: + { + Eigen::ConjugateGradient<MatrixReplacement, Eigen::Lower|Eigen::Upper, Eigen::IdentityPreconditioner> cg; + cg.compute(A); + x = cg.solve(b); + std::cout << "CG: #iterations: " << cg.iterations() << ", estimated error: " << cg.error() << std::endl; + } + + { + Eigen::BiCGSTAB<MatrixReplacement, Eigen::IdentityPreconditioner> bicg; + bicg.compute(A); + x = bicg.solve(b); + std::cout << "BiCGSTAB: #iterations: " << bicg.iterations() << ", estimated error: " << bicg.error() << std::endl; + } + + { + Eigen::GMRES<MatrixReplacement, Eigen::IdentityPreconditioner> gmres; + gmres.compute(A); + x = gmres.solve(b); + std::cout << "GMRES: #iterations: " << gmres.iterations() << ", estimated error: " << gmres.error() << std::endl; + } + + { + Eigen::DGMRES<MatrixReplacement, Eigen::IdentityPreconditioner> gmres; + gmres.compute(A); + x = gmres.solve(b); + std::cout << "DGMRES: #iterations: " << gmres.iterations() << ", estimated error: " << gmres.error() << std::endl; + } + + { + Eigen::MINRES<MatrixReplacement, Eigen::Lower|Eigen::Upper, Eigen::IdentityPreconditioner> minres; + minres.compute(A); + x = minres.solve(b); + std::cout << "MINRES: #iterations: " << minres.iterations() << ", estimated error: " << minres.error() << std::endl; + } +} diff --git a/doc/examples/nullary_indexing.cpp b/doc/examples/nullary_indexing.cpp new file mode 100644 index 000000000..e27c3585a --- /dev/null +++ b/doc/examples/nullary_indexing.cpp @@ -0,0 +1,66 @@ +#include <Eigen/Core> +#include <iostream> + +using namespace Eigen; + +// [functor] +template<class ArgType, class RowIndexType, class ColIndexType> +class indexing_functor { + const ArgType &m_arg; + const RowIndexType &m_rowIndices; + const ColIndexType &m_colIndices; +public: + typedef Matrix<typename ArgType::Scalar, + RowIndexType::SizeAtCompileTime, + ColIndexType::SizeAtCompileTime, + ArgType::Flags&RowMajorBit?RowMajor:ColMajor, + RowIndexType::MaxSizeAtCompileTime, + ColIndexType::MaxSizeAtCompileTime> MatrixType; + + indexing_functor(const ArgType& arg, const RowIndexType& row_indices, const ColIndexType& col_indices) + : m_arg(arg), m_rowIndices(row_indices), m_colIndices(col_indices) + {} + + const typename ArgType::Scalar& operator() (Index row, Index col) const { + return m_arg(m_rowIndices[row], m_colIndices[col]); + } +}; +// [functor] + +// [function] +template <class ArgType, class RowIndexType, class ColIndexType> +CwiseNullaryOp<indexing_functor<ArgType,RowIndexType,ColIndexType>, typename indexing_functor<ArgType,RowIndexType,ColIndexType>::MatrixType> +indexing(const Eigen::MatrixBase<ArgType>& arg, const RowIndexType& row_indices, const ColIndexType& col_indices) +{ + typedef indexing_functor<ArgType,RowIndexType,ColIndexType> Func; + typedef typename Func::MatrixType MatrixType; + return MatrixType::NullaryExpr(row_indices.size(), col_indices.size(), Func(arg.derived(), row_indices, col_indices)); +} +// [function] + + +int main() +{ + std::cout << "[main1]\n"; + Eigen::MatrixXi A = Eigen::MatrixXi::Random(4,4); + Array3i ri(1,2,1); + ArrayXi ci(6); ci << 3,2,1,0,0,2; + Eigen::MatrixXi B = indexing(A, ri, ci); + std::cout << "A =" << std::endl; + std::cout << A << std::endl << std::endl; + std::cout << "A([" << ri.transpose() << "], [" << ci.transpose() << "]) =" << std::endl; + std::cout << B << std::endl; + std::cout << "[main1]\n"; + + std::cout << "[main2]\n"; + B = indexing(A, ri+1, ci); + std::cout << "A(ri+1,ci) =" << std::endl; + std::cout << B << std::endl << std::endl; +#if __cplusplus >= 201103L + B = indexing(A, ArrayXi::LinSpaced(13,0,12).unaryExpr([](int x){return x%4;}), ArrayXi::LinSpaced(4,0,3)); + std::cout << "A(ArrayXi::LinSpaced(13,0,12).unaryExpr([](int x){return x%4;}), ArrayXi::LinSpaced(4,0,3)) =" << std::endl; + std::cout << B << std::endl << std::endl; +#endif + std::cout << "[main2]\n"; +} + diff --git a/doc/ftv2node.png b/doc/ftv2node.png Binary files differnew file mode 100644 index 000000000..63c605bb4 --- /dev/null +++ b/doc/ftv2node.png diff --git a/doc/ftv2pnode.png b/doc/ftv2pnode.png Binary files differnew file mode 100644 index 000000000..c6ee22f93 --- /dev/null +++ b/doc/ftv2pnode.png diff --git a/doc/snippets/CMakeLists.txt b/doc/snippets/CMakeLists.txt index 92a22ea61..1baf32fba 100644 --- a/doc/snippets/CMakeLists.txt +++ b/doc/snippets/CMakeLists.txt @@ -14,17 +14,13 @@ foreach(snippet_src ${snippets_SRCS}) if(EIGEN_STANDARD_LIBRARIES_TO_LINK_TO) target_link_libraries(${compile_snippet_target} ${EIGEN_STANDARD_LIBRARIES_TO_LINK_TO}) endif() - get_target_property(compile_snippet_executable - ${compile_snippet_target} LOCATION) add_custom_command( TARGET ${compile_snippet_target} POST_BUILD - COMMAND ${compile_snippet_executable} + COMMAND ${compile_snippet_target} ARGS >${CMAKE_CURRENT_BINARY_DIR}/${snippet}.out ) add_dependencies(all_snippets ${compile_snippet_target}) set_source_files_properties(${CMAKE_CURRENT_BINARY_DIR}/${compile_snippet_src} PROPERTIES OBJECT_DEPENDS ${snippet_src}) endforeach(snippet_src) - -ei_add_target_property(compile_tut_arithmetic_transpose_aliasing COMPILE_FLAGS -DEIGEN_NO_DEBUG)
\ No newline at end of file diff --git a/doc/snippets/Cwise_arg.cpp b/doc/snippets/Cwise_arg.cpp new file mode 100644 index 000000000..3f45133b6 --- /dev/null +++ b/doc/snippets/Cwise_arg.cpp @@ -0,0 +1,3 @@ +ArrayXcf v = ArrayXcf::Random(3); +cout << v << endl << endl; +cout << arg(v) << endl; diff --git a/doc/snippets/Cwise_array_power_array.cpp b/doc/snippets/Cwise_array_power_array.cpp new file mode 100644 index 000000000..432a76ee5 --- /dev/null +++ b/doc/snippets/Cwise_array_power_array.cpp @@ -0,0 +1,4 @@ +Array<double,1,3> x(8,25,3), + e(1./3.,0.5,2.); +cout << "[" << x << "]^[" << e << "] = " << x.pow(e) << endl; // using ArrayBase::pow +cout << "[" << x << "]^[" << e << "] = " << pow(x,e) << endl; // using Eigen::pow diff --git a/doc/snippets/Cwise_atan.cpp b/doc/snippets/Cwise_atan.cpp new file mode 100644 index 000000000..446844726 --- /dev/null +++ b/doc/snippets/Cwise_atan.cpp @@ -0,0 +1,2 @@ +ArrayXd v = ArrayXd::LinSpaced(5,0,1); +cout << v.atan() << endl; diff --git a/doc/snippets/Cwise_boolean_not.cpp b/doc/snippets/Cwise_boolean_not.cpp new file mode 100644 index 000000000..40009f15a --- /dev/null +++ b/doc/snippets/Cwise_boolean_not.cpp @@ -0,0 +1,5 @@ +Array3d v(1,2,3); +v(1) *= 0.0/0.0; +v(2) /= 0.0; +cout << v << endl << endl; +cout << !isfinite(v) << endl; diff --git a/doc/snippets/Cwise_boolean_xor.cpp b/doc/snippets/Cwise_boolean_xor.cpp new file mode 100644 index 000000000..fafbec806 --- /dev/null +++ b/doc/snippets/Cwise_boolean_xor.cpp @@ -0,0 +1,2 @@ +Array3d v(-1,2,1), w(-3,2,3); +cout << ((v<w) ^ (v<0)) << endl; diff --git a/doc/snippets/Cwise_ceil.cpp b/doc/snippets/Cwise_ceil.cpp new file mode 100644 index 000000000..76cf661f4 --- /dev/null +++ b/doc/snippets/Cwise_ceil.cpp @@ -0,0 +1,3 @@ +ArrayXd v = ArrayXd::LinSpaced(7,-2,2); +cout << v << endl << endl; +cout << ceil(v) << endl; diff --git a/doc/snippets/Cwise_cosh.cpp b/doc/snippets/Cwise_cosh.cpp new file mode 100644 index 000000000..80ee75da5 --- /dev/null +++ b/doc/snippets/Cwise_cosh.cpp @@ -0,0 +1,2 @@ +ArrayXd v = ArrayXd::LinSpaced(5,0,1); +cout << cosh(v) << endl; diff --git a/doc/snippets/Cwise_floor.cpp b/doc/snippets/Cwise_floor.cpp new file mode 100644 index 000000000..73756b41c --- /dev/null +++ b/doc/snippets/Cwise_floor.cpp @@ -0,0 +1,3 @@ +ArrayXd v = ArrayXd::LinSpaced(7,-2,2); +cout << v << endl << endl; +cout << floor(v) << endl; diff --git a/doc/snippets/Cwise_isFinite.cpp b/doc/snippets/Cwise_isFinite.cpp new file mode 100644 index 000000000..1da55fd16 --- /dev/null +++ b/doc/snippets/Cwise_isFinite.cpp @@ -0,0 +1,5 @@ +Array3d v(1,2,3); +v(1) *= 0.0/0.0; +v(2) /= 0.0; +cout << v << endl << endl; +cout << isfinite(v) << endl; diff --git a/doc/snippets/Cwise_isInf.cpp b/doc/snippets/Cwise_isInf.cpp new file mode 100644 index 000000000..be793081c --- /dev/null +++ b/doc/snippets/Cwise_isInf.cpp @@ -0,0 +1,5 @@ +Array3d v(1,2,3); +v(1) *= 0.0/0.0; +v(2) /= 0.0; +cout << v << endl << endl; +cout << isinf(v) << endl; diff --git a/doc/snippets/Cwise_isNaN.cpp b/doc/snippets/Cwise_isNaN.cpp new file mode 100644 index 000000000..7b2a93082 --- /dev/null +++ b/doc/snippets/Cwise_isNaN.cpp @@ -0,0 +1,5 @@ +Array3d v(1,2,3); +v(1) *= 0.0/0.0; +v(2) /= 0.0; +cout << v << endl << endl; +cout << isnan(v) << endl; diff --git a/doc/snippets/Cwise_log10.cpp b/doc/snippets/Cwise_log10.cpp new file mode 100644 index 000000000..b7ae4a834 --- /dev/null +++ b/doc/snippets/Cwise_log10.cpp @@ -0,0 +1,2 @@ +Array4d v(-1,0,1,2); +cout << log10(v) << endl; diff --git a/doc/snippets/Cwise_round.cpp b/doc/snippets/Cwise_round.cpp new file mode 100644 index 000000000..e5c88230b --- /dev/null +++ b/doc/snippets/Cwise_round.cpp @@ -0,0 +1,3 @@ +ArrayXd v = ArrayXd::LinSpaced(7,-2,2); +cout << v << endl << endl; +cout << round(v) << endl; diff --git a/doc/snippets/Cwise_scalar_power_array.cpp b/doc/snippets/Cwise_scalar_power_array.cpp new file mode 100644 index 000000000..c968b2c84 --- /dev/null +++ b/doc/snippets/Cwise_scalar_power_array.cpp @@ -0,0 +1,2 @@ +Array<double,1,3> e(2,-3,1./3.); +cout << "10^[" << e << "] = " << pow(10,e) << endl; diff --git a/doc/snippets/Cwise_sign.cpp b/doc/snippets/Cwise_sign.cpp new file mode 100644 index 000000000..49920e4f1 --- /dev/null +++ b/doc/snippets/Cwise_sign.cpp @@ -0,0 +1,2 @@ +Array3d v(-3,5,0); +cout << v.sign() << endl; diff --git a/doc/snippets/Cwise_sinh.cpp b/doc/snippets/Cwise_sinh.cpp new file mode 100644 index 000000000..fac9b19a8 --- /dev/null +++ b/doc/snippets/Cwise_sinh.cpp @@ -0,0 +1,2 @@ +ArrayXd v = ArrayXd::LinSpaced(5,0,1); +cout << sinh(v) << endl; diff --git a/doc/snippets/Cwise_tanh.cpp b/doc/snippets/Cwise_tanh.cpp new file mode 100644 index 000000000..30cd0450d --- /dev/null +++ b/doc/snippets/Cwise_tanh.cpp @@ -0,0 +1,2 @@ +ArrayXd v = ArrayXd::LinSpaced(5,0,1); +cout << tanh(v) << endl; diff --git a/doc/snippets/DenseBase_LinSpacedInt.cpp b/doc/snippets/DenseBase_LinSpacedInt.cpp new file mode 100644 index 000000000..0d7ae068e --- /dev/null +++ b/doc/snippets/DenseBase_LinSpacedInt.cpp @@ -0,0 +1,8 @@ +cout << "Even spacing inputs:" << endl; +cout << VectorXi::LinSpaced(8,1,4).transpose() << endl; +cout << VectorXi::LinSpaced(8,1,8).transpose() << endl; +cout << VectorXi::LinSpaced(8,1,15).transpose() << endl; +cout << "Uneven spacing inputs:" << endl; +cout << VectorXi::LinSpaced(8,1,7).transpose() << endl; +cout << VectorXi::LinSpaced(8,1,9).transpose() << endl; +cout << VectorXi::LinSpaced(8,1,16).transpose() << endl; diff --git a/doc/snippets/DirectionWise_hnormalized.cpp b/doc/snippets/DirectionWise_hnormalized.cpp new file mode 100644 index 000000000..3410790a8 --- /dev/null +++ b/doc/snippets/DirectionWise_hnormalized.cpp @@ -0,0 +1,7 @@ +typedef Matrix<double,4,Dynamic> Matrix4Xd; +Matrix4Xd M = Matrix4Xd::Random(4,5); +Projective3d P(Matrix4d::Random()); +cout << "The matrix M is:" << endl << M << endl << endl; +cout << "M.colwise().hnormalized():" << endl << M.colwise().hnormalized() << endl << endl; +cout << "P*M:" << endl << P*M << endl << endl; +cout << "(P*M).colwise().hnormalized():" << endl << (P*M).colwise().hnormalized() << endl << endl;
\ No newline at end of file diff --git a/doc/snippets/EigenSolver_eigenvectors.cpp b/doc/snippets/EigenSolver_eigenvectors.cpp index 0fad4dadb..8355f76c9 100644 --- a/doc/snippets/EigenSolver_eigenvectors.cpp +++ b/doc/snippets/EigenSolver_eigenvectors.cpp @@ -1,4 +1,4 @@ MatrixXd ones = MatrixXd::Ones(3,3); EigenSolver<MatrixXd> es(ones); -cout << "The first eigenvector of the 3x3 matrix of ones is:" - << endl << es.eigenvectors().col(1) << endl; +cout << "The first eigenvector of the 3x3 matrix of ones is:" + << endl << es.eigenvectors().col(0) << endl; diff --git a/doc/snippets/LeastSquaresNormalEquations.cpp b/doc/snippets/LeastSquaresNormalEquations.cpp new file mode 100644 index 000000000..997cf1715 --- /dev/null +++ b/doc/snippets/LeastSquaresNormalEquations.cpp @@ -0,0 +1,4 @@ +MatrixXf A = MatrixXf::Random(3, 2); +VectorXf b = VectorXf::Random(3); +cout << "The solution using normal equations is:\n" + << (A.transpose() * A).ldlt().solve(A.transpose() * b) << endl; diff --git a/doc/snippets/LeastSquaresQR.cpp b/doc/snippets/LeastSquaresQR.cpp new file mode 100644 index 000000000..6c9704547 --- /dev/null +++ b/doc/snippets/LeastSquaresQR.cpp @@ -0,0 +1,4 @@ +MatrixXf A = MatrixXf::Random(3, 2); +VectorXf b = VectorXf::Random(3); +cout << "The solution using the QR decomposition is:\n" + << A.colPivHouseholderQr().solve(b) << endl; diff --git a/doc/snippets/MatrixBase_applyOnTheLeft.cpp b/doc/snippets/MatrixBase_applyOnTheLeft.cpp new file mode 100644 index 000000000..6398c873a --- /dev/null +++ b/doc/snippets/MatrixBase_applyOnTheLeft.cpp @@ -0,0 +1,7 @@ +Matrix3f A = Matrix3f::Random(3,3), B; +B << 0,1,0, + 0,0,1, + 1,0,0; +cout << "At start, A = " << endl << A << endl; +A.applyOnTheLeft(B); +cout << "After applyOnTheLeft, A = " << endl << A << endl; diff --git a/doc/snippets/MatrixBase_applyOnTheRight.cpp b/doc/snippets/MatrixBase_applyOnTheRight.cpp new file mode 100644 index 000000000..e4b71b2d8 --- /dev/null +++ b/doc/snippets/MatrixBase_applyOnTheRight.cpp @@ -0,0 +1,9 @@ +Matrix3f A = Matrix3f::Random(3,3), B; +B << 0,1,0, + 0,0,1, + 1,0,0; +cout << "At start, A = " << endl << A << endl; +A *= B; +cout << "After A *= B, A = " << endl << A << endl; +A.applyOnTheRight(B); // equivalent to A *= B +cout << "After applyOnTheRight, A = " << endl << A << endl; diff --git a/doc/snippets/MatrixBase_cwiseSign.cpp b/doc/snippets/MatrixBase_cwiseSign.cpp new file mode 100644 index 000000000..efd717955 --- /dev/null +++ b/doc/snippets/MatrixBase_cwiseSign.cpp @@ -0,0 +1,4 @@ +MatrixXd m(2,3); +m << 2, -4, 6, + -5, 1, 0; +cout << m.cwiseSign() << endl; diff --git a/doc/snippets/MatrixBase_hnormalized.cpp b/doc/snippets/MatrixBase_hnormalized.cpp new file mode 100644 index 000000000..652cd77c0 --- /dev/null +++ b/doc/snippets/MatrixBase_hnormalized.cpp @@ -0,0 +1,6 @@ +Vector4d v = Vector4d::Random(); +Projective3d P(Matrix4d::Random()); +cout << "v = " << v.transpose() << "]^T" << endl; +cout << "v.hnormalized() = " << v.hnormalized().transpose() << "]^T" << endl; +cout << "P*v = " << (P*v).transpose() << "]^T" << endl; +cout << "(P*v).hnormalized() = " << (P*v).hnormalized().transpose() << "]^T" << endl;
\ No newline at end of file diff --git a/doc/snippets/MatrixBase_homogeneous.cpp b/doc/snippets/MatrixBase_homogeneous.cpp new file mode 100644 index 000000000..457c28f91 --- /dev/null +++ b/doc/snippets/MatrixBase_homogeneous.cpp @@ -0,0 +1,6 @@ +Vector3d v = Vector3d::Random(), w; +Projective3d P(Matrix4d::Random()); +cout << "v = [" << v.transpose() << "]^T" << endl; +cout << "h.homogeneous() = [" << v.homogeneous().transpose() << "]^T" << endl; +cout << "(P * v.homogeneous()) = [" << (P * v.homogeneous()).transpose() << "]^T" << endl; +cout << "(P * v.homogeneous()).hnormalized() = [" << (P * v.homogeneous()).eval().hnormalized().transpose() << "]^T" << endl;
\ No newline at end of file diff --git a/doc/snippets/MatrixBase_marked.cpp b/doc/snippets/MatrixBase_marked.cpp deleted file mode 100644 index f60712178..000000000 --- a/doc/snippets/MatrixBase_marked.cpp +++ /dev/null @@ -1,14 +0,0 @@ -#ifndef _MSC_VER - #warning deprecated -#endif -/* -Matrix3d m = Matrix3d::Zero(); -m.part<Eigen::UpperTriangular>().setOnes(); -cout << "Here is the matrix m:" << endl << m << endl; -Matrix3d n = Matrix3d::Ones(); -n.part<Eigen::LowerTriangular>() *= 2; -cout << "Here is the matrix n:" << endl << n << endl; -cout << "And now here is m.inverse()*n, taking advantage of the fact that" - " m is upper-triangular:" << endl - << m.marked<Eigen::UpperTriangular>().solveTriangular(n); -*/
\ No newline at end of file diff --git a/doc/snippets/MatrixBase_part.cpp b/doc/snippets/MatrixBase_part.cpp deleted file mode 100644 index d3e7f482e..000000000 --- a/doc/snippets/MatrixBase_part.cpp +++ /dev/null @@ -1,13 +0,0 @@ -#ifndef _MSC_VER - #warning deprecated -#endif -/* -Matrix3d m = Matrix3d::Zero(); -m.part<Eigen::StrictlyUpperTriangular>().setOnes(); -cout << "Here is the matrix m:" << endl << m << endl; -cout << "And let us now compute m*m.adjoint() in a very optimized way" << endl - << "taking advantage of the symmetry." << endl; -Matrix3d n; -n.part<Eigen::SelfAdjoint>() = (m*m.adjoint()).lazy(); -cout << "The result is:" << endl << n << endl; -*/
\ No newline at end of file diff --git a/doc/snippets/MatrixBase_selfadjointView.cpp b/doc/snippets/MatrixBase_selfadjointView.cpp new file mode 100644 index 000000000..4bd3c7eeb --- /dev/null +++ b/doc/snippets/MatrixBase_selfadjointView.cpp @@ -0,0 +1,6 @@ +Matrix3i m = Matrix3i::Random(); +cout << "Here is the matrix m:" << endl << m << endl; +cout << "Here is the symmetric matrix extracted from the upper part of m:" << endl + << Matrix3i(m.selfadjointView<Upper>()) << endl; +cout << "Here is the symmetric matrix extracted from the lower part of m:" << endl + << Matrix3i(m.selfadjointView<Lower>()) << endl; diff --git a/doc/snippets/MatrixBase_extract.cpp b/doc/snippets/MatrixBase_triangularView.cpp index c96220f72..03aa303f0 100644 --- a/doc/snippets/MatrixBase_extract.cpp +++ b/doc/snippets/MatrixBase_triangularView.cpp @@ -1,13 +1,9 @@ -#ifndef _MSC_VER - #warning deprecated -#endif -/* deprecated Matrix3i m = Matrix3i::Random(); cout << "Here is the matrix m:" << endl << m << endl; cout << "Here is the upper-triangular matrix extracted from m:" << endl - << m.part<Eigen::UpperTriangular>() << endl; + << Matrix3i(m.triangularView<Eigen::Upper>()) << endl; cout << "Here is the strictly-upper-triangular matrix extracted from m:" << endl - << m.part<Eigen::StrictlyUpperTriangular>() << endl; + << Matrix3i(m.triangularView<Eigen::StrictlyUpper>()) << endl; cout << "Here is the unit-lower-triangular matrix extracted from m:" << endl - << m.part<Eigen::UnitLowerTriangular>() << endl; -*/
\ No newline at end of file + << Matrix3i(m.triangularView<Eigen::UnitLower>()) << endl; +// FIXME need to implement output for triangularViews (Bug 885) diff --git a/doc/snippets/SparseMatrix_coeffs.cpp b/doc/snippets/SparseMatrix_coeffs.cpp new file mode 100644 index 000000000..f71a69b07 --- /dev/null +++ b/doc/snippets/SparseMatrix_coeffs.cpp @@ -0,0 +1,9 @@ +SparseMatrix<double> A(3,3); +A.insert(1,2) = 0; +A.insert(0,1) = 1; +A.insert(2,0) = 2; +A.makeCompressed(); +cout << "The matrix A is:" << endl << MatrixXd(A) << endl; +cout << "it has " << A.nonZeros() << " stored non zero coefficients that are: " << A.coeffs().transpose() << endl; +A.coeffs() += 10; +cout << "After adding 10 to every stored non zero coefficient, the matrix A is:" << endl << MatrixXd(A) << endl; diff --git a/doc/snippets/TopicAliasing_mult4.cpp b/doc/snippets/TopicAliasing_mult4.cpp new file mode 100644 index 000000000..8a8992f6c --- /dev/null +++ b/doc/snippets/TopicAliasing_mult4.cpp @@ -0,0 +1,5 @@ +MatrixXf A(2,2), B(3,2); +B << 2, 0, 0, 3, 1, 1; +A << 2, 0, 0, -2; +A = (B * A).cwiseAbs(); +cout << A;
\ No newline at end of file diff --git a/doc/snippets/TopicAliasing_mult5.cpp b/doc/snippets/TopicAliasing_mult5.cpp new file mode 100644 index 000000000..1a36defde --- /dev/null +++ b/doc/snippets/TopicAliasing_mult5.cpp @@ -0,0 +1,5 @@ +MatrixXf A(2,2), B(3,2); +B << 2, 0, 0, 3, 1, 1; +A << 2, 0, 0, -2; +A = (B * A).eval().cwiseAbs(); +cout << A; diff --git a/doc/snippets/Triangular_solve.cpp b/doc/snippets/Triangular_solve.cpp new file mode 100644 index 000000000..548442467 --- /dev/null +++ b/doc/snippets/Triangular_solve.cpp @@ -0,0 +1,11 @@ +Matrix3d m = Matrix3d::Zero(); +m.triangularView<Eigen::Upper>().setOnes(); +cout << "Here is the matrix m:\n" << m << endl; +Matrix3d n = Matrix3d::Ones(); +n.triangularView<Eigen::Lower>() *= 2; +cout << "Here is the matrix n:\n" << n << endl; +cout << "And now here is m.inverse()*n, taking advantage of the fact that" + " m is upper-triangular:\n" + << m.triangularView<Eigen::Upper>().solve(n) << endl; +cout << "And this is n*m.inverse():\n" + << m.triangularView<Eigen::Upper>().solve<Eigen::OnTheRight>(n); diff --git a/doc/snippets/Tutorial_AdvancedInitialization_Join.cpp b/doc/snippets/Tutorial_AdvancedInitialization_Join.cpp index 84e8715cb..55a21539d 100644 --- a/doc/snippets/Tutorial_AdvancedInitialization_Join.cpp +++ b/doc/snippets/Tutorial_AdvancedInitialization_Join.cpp @@ -3,7 +3,7 @@ vec1 << 1, 2, 3; std::cout << "vec1 = " << vec1 << std::endl; RowVectorXd vec2(4); -vec2 << 1, 4, 9, 16;; +vec2 << 1, 4, 9, 16; std::cout << "vec2 = " << vec2 << std::endl; RowVectorXd joined(7); diff --git a/doc/snippets/Tutorial_ReshapeMat2Mat.cpp b/doc/snippets/Tutorial_ReshapeMat2Mat.cpp new file mode 100644 index 000000000..f84d6e76d --- /dev/null +++ b/doc/snippets/Tutorial_ReshapeMat2Mat.cpp @@ -0,0 +1,6 @@ +MatrixXf M1(2,6); // Column-major storage +M1 << 1, 2, 3, 4, 5, 6, + 7, 8, 9, 10, 11, 12; + +Map<MatrixXf> M2(M1.data(), 6,2); +cout << "M2:" << endl << M2 << endl;
\ No newline at end of file diff --git a/doc/snippets/Tutorial_ReshapeMat2Vec.cpp b/doc/snippets/Tutorial_ReshapeMat2Vec.cpp new file mode 100644 index 000000000..95bd4e0e6 --- /dev/null +++ b/doc/snippets/Tutorial_ReshapeMat2Vec.cpp @@ -0,0 +1,11 @@ +MatrixXf M1(3,3); // Column-major storage +M1 << 1, 2, 3, + 4, 5, 6, + 7, 8, 9; + +Map<RowVectorXf> v1(M1.data(), M1.size()); +cout << "v1:" << endl << v1 << endl; + +Matrix<float,Dynamic,Dynamic,RowMajor> M2(M1); +Map<RowVectorXf> v2(M2.data(), M2.size()); +cout << "v2:" << endl << v2 << endl;
\ No newline at end of file diff --git a/doc/snippets/Tutorial_SlicingCol.cpp b/doc/snippets/Tutorial_SlicingCol.cpp new file mode 100644 index 000000000..f667ff689 --- /dev/null +++ b/doc/snippets/Tutorial_SlicingCol.cpp @@ -0,0 +1,11 @@ +MatrixXf M1 = MatrixXf::Random(3,8); +cout << "Column major input:" << endl << M1 << "\n"; +Map<MatrixXf,0,OuterStride<> > M2(M1.data(), M1.rows(), (M1.cols()+2)/3, OuterStride<>(M1.outerStride()*3)); +cout << "1 column over 3:" << endl << M2 << "\n"; + +typedef Matrix<float,Dynamic,Dynamic,RowMajor> RowMajorMatrixXf; +RowMajorMatrixXf M3(M1); +cout << "Row major input:" << endl << M3 << "\n"; +Map<RowMajorMatrixXf,0,Stride<Dynamic,3> > M4(M3.data(), M3.rows(), (M3.cols()+2)/3, + Stride<Dynamic,3>(M3.outerStride(),3)); +cout << "1 column over 3:" << endl << M4 << "\n";
\ No newline at end of file diff --git a/doc/snippets/Tutorial_SlicingVec.cpp b/doc/snippets/Tutorial_SlicingVec.cpp new file mode 100644 index 000000000..07e10bf69 --- /dev/null +++ b/doc/snippets/Tutorial_SlicingVec.cpp @@ -0,0 +1,4 @@ +RowVectorXf v = RowVectorXf::LinSpaced(20,0,19); +cout << "Input:" << endl << v << endl; +Map<RowVectorXf,0,InnerStride<2> > v2(v.data(), v.size()/2); +cout << "Even:" << v2 << endl;
\ No newline at end of file diff --git a/doc/snippets/VectorwiseOp_homogeneous.cpp b/doc/snippets/VectorwiseOp_homogeneous.cpp new file mode 100644 index 000000000..aba4fed0e --- /dev/null +++ b/doc/snippets/VectorwiseOp_homogeneous.cpp @@ -0,0 +1,7 @@ +typedef Matrix<double,3,Dynamic> Matrix3Xd; +Matrix3Xd M = Matrix3Xd::Random(3,5); +Projective3d P(Matrix4d::Random()); +cout << "The matrix M is:" << endl << M << endl << endl; +cout << "M.colwise().homogeneous():" << endl << M.colwise().homogeneous() << endl << endl; +cout << "P * M.colwise().homogeneous():" << endl << P * M.colwise().homogeneous() << endl << endl; +cout << "P * M.colwise().homogeneous().hnormalized(): " << endl << (P * M.colwise().homogeneous()).colwise().hnormalized() << endl << endl;
\ No newline at end of file diff --git a/doc/snippets/compile_snippet.cpp.in b/doc/snippets/compile_snippet.cpp.in index 82ae89162..d63f371a3 100644 --- a/doc/snippets/compile_snippet.cpp.in +++ b/doc/snippets/compile_snippet.cpp.in @@ -1,5 +1,13 @@ -#include <Eigen/Eigen> +static bool eigen_did_assert = false; +#define eigen_assert(X) if(!eigen_did_assert && !(X)){ std::cout << "### Assertion raised in " << __FILE__ << ":" << __LINE__ << ":\n" #X << "\n### The following would happen without assertions:\n"; eigen_did_assert = true;} + #include <iostream> +#include <Eigen/Eigen> + +#ifndef M_PI +#define M_PI 3.1415926535897932384626433832795 +#endif + using namespace Eigen; using namespace std; diff --git a/doc/special_examples/CMakeLists.txt b/doc/special_examples/CMakeLists.txt index 45e2339e3..101fbc5f9 100644 --- a/doc/special_examples/CMakeLists.txt +++ b/doc/special_examples/CMakeLists.txt @@ -1,4 +1,3 @@ - if(NOT EIGEN_TEST_NOQT) find_package(Qt4) if(QT4_FOUND) @@ -6,17 +5,17 @@ if(NOT EIGEN_TEST_NOQT) endif() endif(NOT EIGEN_TEST_NOQT) - if(QT4_FOUND) add_executable(Tutorial_sparse_example Tutorial_sparse_example.cpp Tutorial_sparse_example_details.cpp) target_link_libraries(Tutorial_sparse_example ${EIGEN_STANDARD_LIBRARIES_TO_LINK_TO} ${QT_QTCORE_LIBRARY} ${QT_QTGUI_LIBRARY}) - + add_custom_command( - TARGET Tutorial_sparse_example - POST_BUILD - COMMAND Tutorial_sparse_example ARGS ${CMAKE_CURRENT_BINARY_DIR}/../html/Tutorial_sparse_example.jpeg + TARGET Tutorial_sparse_example + POST_BUILD + COMMAND ${CMAKE_COMMAND} -E make_directory ${CMAKE_CURRENT_BINARY_DIR}/../html/ + COMMAND Tutorial_sparse_example ARGS ${CMAKE_CURRENT_BINARY_DIR}/../html/Tutorial_sparse_example.jpeg ) - + add_dependencies(all_examples Tutorial_sparse_example) endif(QT4_FOUND) @@ -26,15 +25,11 @@ if(EIGEN_COMPILER_SUPPORT_CPP11) target_link_libraries(random_cpp11 ${EIGEN_STANDARD_LIBRARIES_TO_LINK_TO}) add_dependencies(all_examples random_cpp11) ei_add_target_property(random_cpp11 COMPILE_FLAGS "-std=c++11") - - get_target_property(random_cpp11_exec - random_cpp11 LOCATION) + add_custom_command( TARGET random_cpp11 POST_BUILD - COMMAND ${random_cpp11_exec} + COMMAND random_cpp11 ARGS >${CMAKE_CURRENT_BINARY_DIR}/random_cpp11.out ) - - -endif()
\ No newline at end of file +endif() diff --git a/doc/special_examples/Tutorial_sparse_example.cpp b/doc/special_examples/Tutorial_sparse_example.cpp index 002f19f01..830e196ea 100644 --- a/doc/special_examples/Tutorial_sparse_example.cpp +++ b/doc/special_examples/Tutorial_sparse_example.cpp @@ -9,6 +9,8 @@ void saveAsBitmap(const Eigen::VectorXd& x, int n, const char* filename); int main(int argc, char** argv) { + assert(argc==2); + int n = 300; // size of the image int m = n*n; // number of unknows (=number of pixels) diff --git a/doc/special_examples/Tutorial_sparse_example_details.cpp b/doc/special_examples/Tutorial_sparse_example_details.cpp index 7d820b44a..bc18b0188 100644 --- a/doc/special_examples/Tutorial_sparse_example_details.cpp +++ b/doc/special_examples/Tutorial_sparse_example_details.cpp @@ -8,7 +8,7 @@ typedef Eigen::Triplet<double> T; void insertCoefficient(int id, int i, int j, double w, std::vector<T>& coeffs, Eigen::VectorXd& b, const Eigen::VectorXd& boundary) { - int n = boundary.size(); + int n = int(boundary.size()); int id1 = i+j*n; if(i==-1 || i==n) b(id) -= w * boundary(j); // constrained coefficient diff --git a/doc/special_examples/random_cpp11.cpp b/doc/special_examples/random_cpp11.cpp index ccd7c77d0..33744c051 100644 --- a/doc/special_examples/random_cpp11.cpp +++ b/doc/special_examples/random_cpp11.cpp @@ -7,7 +7,7 @@ using namespace Eigen; int main() { std::default_random_engine generator; std::poisson_distribution<int> distribution(4.1); - auto poisson = [&] (int) {return distribution(generator);}; + auto poisson = [&] () {return distribution(generator);}; RowVectorXi v = RowVectorXi::NullaryExpr(10, poisson ); std::cout << v << "\n"; |
