diff --git a/doc/C06_TutorialLinearAlgebra.dox b/doc/C06_TutorialLinearAlgebra.dox
new file mode 100644
index 000000000..f31fcfc0f
--- /dev/null
+++ b/doc/C06_TutorialLinearAlgebra.dox
@@ -0,0 +1,118 @@
+namespace Eigen {
+
+/** \page TutorialLinearAlgebra Tutorial page 6 - Linear algebra and decompositions
+ \ingroup Tutorial
+
+\li \b Previous: TODO
+\li \b Next: TODO
+
+This tutorial explains how to solve linear systems, compute various decompositions such as LU,
+QR, SVD, eigendecompositions... for more advanced topics, don't miss our special page on
+\ref TopicLinearAlgebraDecompositions "this topic".
+
+\section TutorialLinAlgBasicSolve How do I solve a system of linear equations?
+
+\b The \b problem: You have a system of equations, that you have written as a single matrix equation
+ \f[ Ax \: = \: b \f]
+Where \a A and \a b are matrices (\a b could be a vector, as a special case). You want to find a solution \a x.
+
+\b The \b solution: You can choose between various decompositions, depending on what your matrix \a A looks like,
+and depending on whether you favor speed or accuracy. However, let's start with an example that works in all cases,
+and is a good compromise:
+
+
+ | \include TutorialLinAlgExSolveColPivHouseholderQR.cpp |
+ output: \verbinclude TutorialLinAlgExSolveColPivHouseholderQR.out |
+
+
+
+In this example, the colPivHouseholderQr() method returns an object of class ColPivHouseholderQR. This line could
+have been replaced by:
+\code
+ColPivHouseholderQR dec(A);
+Vector3f x = dec.solve(b);
+\endcode
+
+Here, ColPivHouseholderQR is a QR decomposition with column pivoting. It's a good compromise for this tutorial, as it
+works for all matrices while being quite fast. Here is a table of some other decompositions that you can choose from,
+depending on your matrix and the trade-off you want to make:
+
+
+
+
+ | Decomposition |
+ Method |
+ Requirements on the matrix |
+ Speed |
+ Accuracy |
+
+
+
+ | PartialPivLU |
+ partialPivLu() |
+ Invertible |
+ ++ |
+ + |
+
+
+
+ | FullPivLU |
+ fullPivLu() |
+ None |
+ - |
+ +++ |
+
+
+
+ | HouseholderQR |
+ householderQr() |
+ None |
+ ++ |
+ + |
+
+
+
+ | ColPivHouseholderQR |
+ colPivHouseholderQr() |
+ None |
+ + |
+ ++ |
+
+
+
+ | FullPivHouseholderQR |
+ fullPivHouseholderQr() |
+ None |
+ - |
+ +++ |
+
+
+
+ | LLT |
+ llt() |
+ Positive definite |
+ +++ |
+ + |
+
+
+
+ | LDLT |
+ ldlt() |
+ Positive or negative semidefinite |
+ +++ |
+ ++ |
+
+
+
+
+All of these decompositions offer a solve() method that works as in the above example.
+
+For a much more complete table comparing all decompositions supported by Eigen (notice that Eigen
+supports many other decompositions), see our special page on
+\ref TopicLinearAlgebraDecompositions "this topic".
+
+
+
+*/
+
+}
diff --git a/doc/TopicLinearAlgebraDecompositions.dox b/doc/TopicLinearAlgebraDecompositions.dox
index 5065be70a..1def1e776 100644
--- a/doc/TopicLinearAlgebraDecompositions.dox
+++ b/doc/TopicLinearAlgebraDecompositions.dox
@@ -50,8 +50,8 @@ namespace Eigen {
HouseholderQR |
- |
Fast |
- Average |
Depends on condition number |
+ - |
Orthogonalization |
Yes |
Excellent |
diff --git a/doc/examples/TutorialLinAlgExSolveColPivHouseholderQR.cpp b/doc/examples/TutorialLinAlgExSolveColPivHouseholderQR.cpp
new file mode 100644
index 000000000..29c22be41
--- /dev/null
+++ b/doc/examples/TutorialLinAlgExSolveColPivHouseholderQR.cpp
@@ -0,0 +1,17 @@
+#include
+#include
+
+using namespace std;
+using namespace Eigen;
+
+int main()
+{
+ Matrix3f A;
+ Vector3f b;
+ A << 1,2,3, 4,5,6, 7,8,10;
+ b << 3, 3, 4;
+ cout << "Here is the matrix A:" << endl << A << endl;
+ cout << "Here is the vector b:" << endl << b << endl;
+ Vector3f x = A.colPivHouseholderQr().solve(b);
+ cout << "The solution is:" << endl << x << endl;
+}