SuiteSparse cleanup.
1. CreateSparseMatrixTransposeView now returns a struct instead
of a pointer.
2. Add AnalyzeCholeskyWithNaturalOrdering.
Change-Id: If27a5502949c3994edd95be0d25ec7a0d1fa1ae1
diff --git a/internal/ceres/sparse_normal_cholesky_solver.cc b/internal/ceres/sparse_normal_cholesky_solver.cc
index aacba86..7ece965 100644
--- a/internal/ceres/sparse_normal_cholesky_solver.cc
+++ b/internal/ceres/sparse_normal_cholesky_solver.cc
@@ -199,9 +199,7 @@
VectorRef(x, num_cols).setZero();
- scoped_ptr<cholmod_sparse> lhs(ss_.CreateSparseMatrixTransposeView(A));
- CHECK_NOTNULL(lhs.get());
-
+ cholmod_sparse lhs = ss_.CreateSparseMatrixTransposeView(A);
cholmod_dense* rhs = ss_.CreateDenseVector(Atb.data(), num_cols, num_cols);
event_logger.AddEvent("Setup");
@@ -216,7 +214,7 @@
}
event_logger.AddEvent("Analysis");
- cholmod_dense* sol = ss_.SolveCholesky(lhs.get(), factor_, rhs);
+ cholmod_dense* sol = ss_.SolveCholesky(&lhs, factor_, rhs);
event_logger.AddEvent("Solve");
ss_.Free(rhs);
diff --git a/internal/ceres/suitesparse.cc b/internal/ceres/suitesparse.cc
index 87fae75..9063fa0 100644
--- a/internal/ceres/suitesparse.cc
+++ b/internal/ceres/suitesparse.cc
@@ -87,23 +87,23 @@
return cholmod_triplet_to_sparse(&triplet, triplet.nnz, &cc_);
}
-cholmod_sparse* SuiteSparse::CreateSparseMatrixTransposeView(
+cholmod_sparse SuiteSparse::CreateSparseMatrixTransposeView(
CompressedRowSparseMatrix* A) {
- cholmod_sparse* m = new cholmod_sparse_struct;
- m->nrow = A->num_cols();
- m->ncol = A->num_rows();
- m->nzmax = A->num_nonzeros();
-
- m->p = reinterpret_cast<void*>(A->mutable_rows());
- m->i = reinterpret_cast<void*>(A->mutable_cols());
- m->x = reinterpret_cast<void*>(A->mutable_values());
-
- m->stype = 0; // Matrix is not symmetric.
- m->itype = CHOLMOD_INT;
- m->xtype = CHOLMOD_REAL;
- m->dtype = CHOLMOD_DOUBLE;
- m->sorted = 1;
- m->packed = 1;
+ cholmod_sparse m;
+ m.nrow = A->num_cols();
+ m.ncol = A->num_rows();
+ m.nzmax = A->num_nonzeros();
+ m.nz = NULL;
+ m.p = reinterpret_cast<void*>(A->mutable_rows());
+ m.i = reinterpret_cast<void*>(A->mutable_cols());
+ m.x = reinterpret_cast<void*>(A->mutable_values());
+ m.z = NULL;
+ m.stype = 0; // Matrix is not symmetric.
+ m.itype = CHOLMOD_INT;
+ m.xtype = CHOLMOD_REAL;
+ m.dtype = CHOLMOD_DOUBLE;
+ m.sorted = 1;
+ m.packed = 1;
return m;
}
@@ -172,6 +172,23 @@
return factor;
}
+cholmod_factor* SuiteSparse::AnalyzeCholeskyWithNaturalOrdering(cholmod_sparse* A) {
+ cc_.nmethods = 1;
+ cc_.method[0].ordering = CHOLMOD_NATURAL;
+ cc_.postorder = 0;
+
+ cholmod_factor* factor = cholmod_analyze(A, &cc_);
+ CHECK_EQ(cc_.status, CHOLMOD_OK)
+ << "Cholmod symbolic analysis failed " << cc_.status;
+ CHECK_NOTNULL(factor);
+
+ if (VLOG_IS_ON(2)) {
+ cholmod_print_common(const_cast<char*>("Symbolic Analysis"), &cc_);
+ }
+
+ return factor;
+}
+
bool SuiteSparse::BlockAMDOrdering(const cholmod_sparse* A,
const vector<int>& row_blocks,
const vector<int>& col_blocks,
diff --git a/internal/ceres/suitesparse.h b/internal/ceres/suitesparse.h
index 5a8ef63..131716d 100644
--- a/internal/ceres/suitesparse.h
+++ b/internal/ceres/suitesparse.h
@@ -70,10 +70,8 @@
// Create a cholmod_sparse wrapper around the contents of A. This is
// a shallow object, which refers to the contents of A and does not
- // use the SuiteSparse machinery to allocate memory, this object
- // should be disposed off with a delete and not a call to Free as is
- // the case for objects returned by CreateSparseMatrixTranspose.
- cholmod_sparse* CreateSparseMatrixTransposeView(CompressedRowSparseMatrix* A);
+ // use the SuiteSparse machinery to allocate memory.
+ cholmod_sparse CreateSparseMatrixTransposeView(CompressedRowSparseMatrix* A);
// Given a vector x, build a cholmod_dense vector of size out_size
// with the first in_size entries copied from x. If x is NULL, then
@@ -135,6 +133,12 @@
cholmod_factor* AnalyzeCholeskyWithUserOrdering(cholmod_sparse* A,
const vector<int>& ordering);
+ // Perform a symbolic factorization of A without re-ordering A. No
+ // postordering of the elimination tree is performed. This ensures
+ // that the symbolic factor does not introduce an extra permutation
+ // on the matrix. See the documentation for CHOLMOD for more details.
+ cholmod_factor* AnalyzeCholeskyWithNaturalOrdering(cholmod_sparse* A);
+
// Use the symbolic factorization in L, to find the numerical
// factorization for the matrix A or AA^T. Return true if
// successful, false otherwise. L contains the numeric factorization
