Use solveInPlace correctly.
When solving a linear system using Eigen's dense Cholesky factorization
if the right hand side of the linear system is the same vector
that will store the solution, call solveInPlace instead of solve.
Change-Id: I3e6d2f21ff420c25217cd87ee5d269fdfabbf19a
diff --git a/internal/ceres/eigen_dense_cholesky.cc b/internal/ceres/eigen_dense_cholesky.cc
index 1b7f87f..993f92c 100644
--- a/internal/ceres/eigen_dense_cholesky.cc
+++ b/internal/ceres/eigen_dense_cholesky.cc
@@ -35,7 +35,7 @@
SolveUpperTriangularUsingCholesky(int size,
const double* lhs_values,
const double* rhs_values,
- double* solution) {
+ double* solution_values) {
ConstMatrixRef lhs(lhs_values, size, size);
// On ARM we have experienced significant numerical problems with
@@ -51,8 +51,13 @@
#endif
if (cholesky.info() == Eigen::Success) {
- ConstVectorRef rhs(rhs_values, size);
- VectorRef(solution, size) = cholesky.solve(rhs);
+ VectorRef solution(solution_values, size);
+ if (solution_values == rhs_values) {
+ cholesky.solveInPlace(solution);
+ } else {
+ ConstVectorRef rhs(rhs_values, size);
+ solution = cholesky.solve(rhs);
+ }
}
return cholesky.info();
diff --git a/internal/ceres/eigen_dense_cholesky.h b/internal/ceres/eigen_dense_cholesky.h
index 79a091e..0d6b18a 100644
--- a/internal/ceres/eigen_dense_cholesky.h
+++ b/internal/ceres/eigen_dense_cholesky.h
@@ -57,8 +57,8 @@
// factorization. rhs_values and solution can point to the same array.
Eigen::ComputationInfo
SolveUpperTriangularUsingCholesky(int size,
- const double* lhs_values,
- const double* rhs_values,
+ const double* lhs,
+ const double* rhs,
double* solution);
} // namespace internal
