Solver::Options::ordering* are dead. Remove the old ordering API, and modify solver_impl.cc to use the new API everywhere. In the process also clean up the linear solver instantion logic in solver_impl.cc a bit too. Change-Id: Ia66898abc7f622070b184b21fce8cc6140c4cebf
diff --git a/internal/ceres/minimizer.h b/internal/ceres/minimizer.h index cfc98a3..667d80a 100644 --- a/internal/ceres/minimizer.h +++ b/internal/ceres/minimizer.h
@@ -74,7 +74,6 @@ lsqp_dump_directory = options.lsqp_dump_directory; lsqp_iterations_to_dump = options.lsqp_iterations_to_dump; lsqp_dump_format_type = options.lsqp_dump_format_type; - num_eliminate_blocks = options.num_eliminate_blocks; max_num_consecutive_invalid_steps = options.max_num_consecutive_invalid_steps; min_trust_region_radius = options.min_trust_region_radius; @@ -104,7 +103,6 @@ vector<int> lsqp_iterations_to_dump; DumpFormatType lsqp_dump_format_type; string lsqp_dump_directory; - int num_eliminate_blocks; int max_num_consecutive_invalid_steps; int min_trust_region_radius;
diff --git a/internal/ceres/ordering.cc b/internal/ceres/ordering.cc index 8b53ea5..4b3ec25 100644 --- a/internal/ceres/ordering.cc +++ b/internal/ceres/ordering.cc
@@ -70,6 +70,12 @@ return it->second; } +bool Ordering::ContainsParameterBlock(double* parameter_block) const { + const map<double*, int>::const_iterator it = + parameter_block_to_group_id_.find(parameter_block); + return (it != parameter_block_to_group_id_.end()); +} + int Ordering::NumParameterBlocks() const { return parameter_block_to_group_id_.size(); }
diff --git a/internal/ceres/solver.cc b/internal/ceres/solver.cc index 84c2756..5a84cfc 100644 --- a/internal/ceres/solver.cc +++ b/internal/ceres/solver.cc
@@ -82,8 +82,6 @@ num_parameters_reduced(-1), num_residual_blocks_reduced(-1), num_residuals_reduced(-1), - num_eliminate_blocks_given(-1), - num_eliminate_blocks_used(-1), num_threads_given(-1), num_threads_used(-1), num_linear_solver_threads_given(-1), @@ -91,7 +89,6 @@ linear_solver_type_given(SPARSE_NORMAL_CHOLESKY), linear_solver_type_used(SPARSE_NORMAL_CHOLESKY), preconditioner_type(IDENTITY), - ordering_type(NATURAL), trust_region_strategy_type(LEVENBERG_MARQUARDT), sparse_linear_algebra_library(SUITE_SPARSE) { } @@ -165,22 +162,7 @@ "N/A", "N/A"); } - internal::StringAppendF(&report, "Ordering %25s%25s\n", - OrderingTypeToString(ordering_type), - OrderingTypeToString(ordering_type)); - - if (IsSchurType(linear_solver_type_given)) { - if (ordering_type == SCHUR) { - internal::StringAppendF(&report, "num_eliminate_blocks%25s% 25d\n", - "N/A", - num_eliminate_blocks_used); - } else { - internal::StringAppendF(&report, "num_eliminate_blocks% 25d% 25d\n", - num_eliminate_blocks_given, - num_eliminate_blocks_used); - } - } - + // TODO(sameeragarwal): Add support for logging the ordering object. internal::StringAppendF(&report, "Threads: % 25d% 25d\n", num_threads_given, num_threads_used); internal::StringAppendF(&report, "Linear solver threads % 23d% 25d\n",
diff --git a/internal/ceres/solver_impl.cc b/internal/ceres/solver_impl.cc index 4bedbe3..7e90c9b 100644 --- a/internal/ceres/solver_impl.cc +++ b/internal/ceres/solver_impl.cc
@@ -52,6 +52,13 @@ #include "ceres/wall_time.h" namespace ceres { + +Solver::Options::~Options() { + if (ordering != NULL) { + delete ordering; + } +} + namespace internal { namespace { @@ -202,57 +209,11 @@ Solver::Summary* summary) { double solver_start_time = WallTimeInSeconds(); Solver::Options options(original_options); - - // Code for supporting both the old and the new ordering API. This - // code will get refactored and re-written as the old API is - // steadily deprecated. - - // Clobber all the input parameters from the old api parameters. - if (options.use_new_ordering_api) { - // For linear solvers which are not of Schur type, do nothing. - options.ordering.clear(); - options.num_eliminate_blocks = 0; - options.ordering_type = NATURAL; - - if (IsSchurType(options.linear_solver_type)) { - if (options.ordering_new_api == NULL) { - // User says we are free to find the independent set, and order - // any which way. - options.ordering_type = SCHUR; - } else { - // User has given an ordering and asked for a Schur type solver. - options.ordering_type = USER; - - // The lowest numbered group corresponds to - // num_eliminate_blocks e_blocks. - const map<int, set<double*> >& group_id_to_parameter_block - = options.ordering_new_api->group_id_to_parameter_blocks(); - - map<int, set<double*> >::const_iterator it = - group_id_to_parameter_block.begin(); - CHECK(it != group_id_to_parameter_block.end()); - - options.num_eliminate_blocks = it->second.size(); - for (; it != group_id_to_parameter_block.end(); ++it) { - options.ordering.insert(options.ordering.end(), - it->second.begin(), - it->second.end()); - } - CHECK_EQ(options.ordering.size(), - original_problem_impl->NumParameterBlocks()); - } - } - } else { - CHECK(options.ordering_new_api == NULL); - } - - Program* original_program = original_problem_impl->mutable_program(); ProblemImpl* problem_impl = original_problem_impl; // Reset the summary object to its default values. *CHECK_NOTNULL(summary) = Solver::Summary(); - #ifndef CERES_USE_OPENMP if (options.num_threads > 1) { LOG(WARNING) @@ -270,25 +231,10 @@ } #endif - if (IsSchurType(options.linear_solver_type) && - options.ordering_type != SCHUR && - options.num_eliminate_blocks < 1) { - summary->error = - StringPrintf("Using a Schur type solver with %s" - " ordering requires that" - " Solver::Options::num_eliminate_blocks" - " be set to a positive integer.", - OrderingTypeToString(options.ordering_type)); - LOG(WARNING) << summary->error; - return; - } - summary->linear_solver_type_given = options.linear_solver_type; - summary->num_eliminate_blocks_given = original_options.num_eliminate_blocks; summary->num_threads_given = original_options.num_threads; summary->num_linear_solver_threads_given = original_options.num_linear_solver_threads; - summary->ordering_type = original_options.ordering_type; summary->num_parameter_blocks = problem_impl->NumParameterBlocks(); summary->num_parameters = problem_impl->NumParameters(); @@ -301,6 +247,23 @@ summary->trust_region_strategy_type = options.trust_region_strategy_type; summary->dogleg_type = options.dogleg_type; + if (original_options.ordering != NULL) { + if (!IsOrderingValid(original_options, problem_impl, &summary->error)) { + LOG(WARNING) << summary->error; + return; + } + options.ordering = new Ordering(*original_options.ordering); + } else { + options.ordering = new Ordering; + const ProblemImpl::ParameterMap& parameter_map = + problem_impl->parameter_map(); + for (ProblemImpl::ParameterMap::const_iterator it = parameter_map.begin(); + it != parameter_map.end(); + ++it) { + options.ordering->AddParameterBlockToGroup(it->first, 0); + } + } + // Evaluate the initial cost, residual vector and the jacobian // matrix if requested by the user. The initial cost needs to be // computed on the original unpreprocessed problem, as it is used to @@ -313,13 +276,14 @@ options.return_initial_residuals ? &summary->initial_residuals : NULL, options.return_initial_gradient ? &summary->initial_gradient : NULL, options.return_initial_jacobian ? &summary->initial_jacobian : NULL); - original_program->SetParameterBlockStatePtrsToUserStatePtrs(); + original_program->SetParameterBlockStatePtrsToUserStatePtrs(); // If the user requests gradient checking, construct a new // ProblemImpl by wrapping the CostFunctions of problem_impl inside // GradientCheckingCostFunction and replacing problem_impl with // gradient_checking_problem_impl. scoped_ptr<ProblemImpl> gradient_checking_problem_impl; + // Save the original problem impl so we don't use the gradient // checking one when computing the residuals. if (options.check_gradients) { @@ -353,7 +317,6 @@ linear_solver(CreateLinearSolver(&options, &summary->error)); summary->linear_solver_type_used = options.linear_solver_type; summary->preconditioner_type = options.preconditioner_type; - summary->num_eliminate_blocks_used = options.num_eliminate_blocks; summary->num_linear_solver_threads_used = options.num_linear_solver_threads; if (linear_solver == NULL) { @@ -421,13 +384,67 @@ WallTimeInSeconds() - post_process_start_time; } +bool SolverImpl::IsOrderingValid(const Solver::Options& options, + const ProblemImpl* problem_impl, + string* error) { + if (options.ordering->NumParameterBlocks() != + problem_impl->NumParameterBlocks()) { + *error = "Number of parameter blocks in user supplied ordering " + "does not match the number of parameter blocks in the problem"; + return false; + } + + const Program& program = problem_impl->program(); + const vector<ParameterBlock*>& parameter_blocks = program.parameter_blocks(); + const vector<ResidualBlock*>& residual_blocks = program.residual_blocks(); + for (vector<ParameterBlock*>::const_iterator it = parameter_blocks.begin(); + it != parameter_blocks.end(); + ++it) { + if (!options.ordering->ContainsParameterBlock( + const_cast<double*>((*it)->user_state()))) { + *error = "Problem contains a parameter block that is not in " + "the user specified ordering."; + return false; + } + } + + if (IsSchurType(options.linear_solver_type) && + options.ordering->NumGroups() > 1) { + const set<double*>& e_blocks = + options.ordering->group_id_to_parameter_blocks().begin()->second; + + // Loop over each residual block and ensure that no two parameter + // blocks in the same residual block are part of the first + // elimination group, as that would violate the assumption that it + // is an independent set in the Hessian matrix. + for (vector<ResidualBlock*>::const_iterator it = residual_blocks.begin(); + it != residual_blocks.end(); + ++it) { + ParameterBlock* const* parameter_blocks = (*it)->parameter_blocks(); + const int num_parameter_blocks = (*it)->NumParameterBlocks(); + int count = 0; + for (int i = 0; i < num_parameter_blocks; ++i) { + count += e_blocks.count(const_cast<double*>( + parameter_blocks[i]->user_state())); + } + + if (count > 1) { + *error = "The user requested the use of a Schur type solver. " + "But the first elimination group in the ordering is not an " + "independent set."; + return false; + } + } + } + return true; +} + // Strips varying parameters and residuals, maintaining order, and updating // num_eliminate_blocks. bool SolverImpl::RemoveFixedBlocksFromProgram(Program* program, - int* num_eliminate_blocks, + Ordering* ordering, double* fixed_cost, string* error) { - int original_num_eliminate_blocks = *num_eliminate_blocks; vector<ParameterBlock*>* parameter_blocks = program->mutable_parameter_blocks(); @@ -484,7 +501,7 @@ } // Filter out unused or fixed parameter blocks, and update - // num_eliminate_blocks as necessary. + // the ordering. { vector<ParameterBlock*>* parameter_blocks = program->mutable_parameter_blocks(); @@ -493,8 +510,8 @@ ParameterBlock* parameter_block = (*parameter_blocks)[i]; if (parameter_block->index() == 1) { (*parameter_blocks)[j++] = parameter_block; - } else if (i < original_num_eliminate_blocks) { - (*num_eliminate_blocks)--; + } else { + ordering->RemoveParameterBlock(parameter_block->mutable_user_state()); } } parameter_blocks->resize(j); @@ -512,35 +529,17 @@ ProblemImpl* problem_impl, double* fixed_cost, string* error) { + CHECK_NOTNULL(options->ordering); Program* original_program = problem_impl->mutable_program(); scoped_ptr<Program> transformed_program(new Program(*original_program)); + Ordering* ordering = options->ordering; - if (options->ordering_type == USER && - !ApplyUserOrdering(*problem_impl, - options->ordering, - transformed_program.get(), - error)) { - return NULL; - } - - if (options->ordering_type == SCHUR && options->num_eliminate_blocks != 0) { - *error = "Can't specify SCHUR ordering and num_eliminate_blocks " - "at the same time; SCHUR ordering determines " - "num_eliminate_blocks automatically."; - return NULL; - } - - if (options->ordering_type == SCHUR && options->ordering.size() != 0) { - *error = "Can't specify SCHUR ordering type and the ordering " - "vector at the same time; SCHUR ordering determines " - "a suitable parameter ordering automatically."; - return NULL; - } - - int num_eliminate_blocks = options->num_eliminate_blocks; + const int min_group_id = + ordering->group_id_to_parameter_blocks().begin()->first; + const int original_num_groups = ordering->NumGroups(); if (!RemoveFixedBlocksFromProgram(transformed_program.get(), - &num_eliminate_blocks, + ordering, fixed_cost, error)) { return NULL; @@ -552,21 +551,70 @@ return transformed_program.release(); } - if (options->ordering_type == SCHUR) { + // If the user supplied an ordering with just one group, it is + // equivalent to the user supplying NULL as ordering. Ceres is + // completely free to choose the parameter block ordering as it sees + // fit. For Schur type solvers, this means that the user wishes for + // Ceres to identify the e_blocks, which we do by computing a + // maximal independent set. + if (original_num_groups == 1 && + IsSchurType(options->linear_solver_type)) { vector<ParameterBlock*> schur_ordering; - num_eliminate_blocks = ComputeSchurOrdering(*transformed_program, - &schur_ordering); + const int num_eliminate_blocks = ComputeSchurOrdering(*original_program, + &schur_ordering); CHECK_EQ(schur_ordering.size(), transformed_program->NumParameterBlocks()) << "Congratulations, you found a Ceres bug! Please report this error " << "to the developers."; - // Replace the transformed program's ordering with the schur ordering. - swap(*transformed_program->mutable_parameter_blocks(), schur_ordering); + for (int i = 0; i < schur_ordering.size(); ++i) { + ordering->AddParameterBlockToGroup(schur_ordering[i]->mutable_user_state(), + (i < num_eliminate_blocks) ? 0 : 1); + } } - options->num_eliminate_blocks = num_eliminate_blocks; - CHECK_GE(options->num_eliminate_blocks, 0) - << "Congratulations, you found a Ceres bug! Please report this error " - << "to the developers."; + + if (!ApplyUserOrdering( + *problem_impl, ordering, transformed_program.get(), error)) { + return NULL; + } + + // If the user requested the use of a Schur type solver, and + // supplied a non-NULL ordering object with more than one + // elimimation group, then it can happen that after all the + // parameter blocks which are fixed or unused have been removed from + // the program and the ordering, there are no more parameter blocks + // in the first elimination group. + // + // In such a case, the use of a Schur type solver is not possible, + // as they assume there is at least one e_block. Thus, we + // automatically switch to one of the other solvers, depending on + // the user's indicated preferences. + if (IsSchurType(options->linear_solver_type) && + original_num_groups > 1 && + ordering->GroupSize(min_group_id) == 0) { + string msg = "No e_blocks remaining. Switching from "; + if (options->linear_solver_type == SPARSE_SCHUR) { + options->linear_solver_type = SPARSE_NORMAL_CHOLESKY; + msg += "SPARSE_SCHUR to SPARSE_NORMAL_CHOLESKY."; + } else if (options->linear_solver_type == DENSE_SCHUR) { + // TODO(sameeragarwal): This is probably not a great choice. + // Ideally, we should have a DENSE_NORMAL_CHOLESKY, that can + // take a BlockSparseMatrix as input. + options->linear_solver_type = DENSE_QR; + msg += "DENSE_SCHUR to DENSE_QR."; + } else if (options->linear_solver_type == ITERATIVE_SCHUR) { + msg += StringPrintf("ITERATIVE_SCHUR with %s preconditioner " + "to CGNR with JACOBI preconditioner.", + PreconditionerTypeToString( + options->preconditioner_type)); + options->linear_solver_type = CGNR; + if (options->preconditioner_type != IDENTITY) { + // CGNR currently only supports the JACOBI preconditioner. + options->preconditioner_type = JACOBI; + } + } + + LOG(WARNING) << msg; + } // Since the transformed program is the "active" program, and it is mutated, // update the parameter offsets and indices. @@ -576,6 +624,10 @@ LinearSolver* SolverImpl::CreateLinearSolver(Solver::Options* options, string* error) { + CHECK_NOTNULL(options); + CHECK_NOTNULL(options->ordering); + CHECK_NOTNULL(error); + if (options->trust_region_strategy_type == DOGLEG) { if (options->linear_solver_type == ITERATIVE_SCHUR || options->linear_solver_type == CGNR) { @@ -593,6 +645,24 @@ "SuiteSparse was not enabled when Ceres was built."; return NULL; } + + if (linear_solver_options.preconditioner_type == SCHUR_JACOBI) { + *error = "SCHUR_JACOBI preconditioner not suppored. Please build Ceres " + "with SuiteSparse support."; + return NULL; + } + + if (linear_solver_options.preconditioner_type == CLUSTER_JACOBI) { + *error = "CLUSTER_JACOBI preconditioner not suppored. Please build Ceres " + "with SuiteSparse support."; + return NULL; + } + + if (linear_solver_options.preconditioner_type == CLUSTER_TRIDIAGONAL) { + *error = "CLUSTER_TRIDIAGONAL preconditioner not suppored. Please build " + "Ceres with SuiteSparse support."; + return NULL; + } #endif #ifdef CERES_NO_CXSPARSE @@ -604,6 +674,13 @@ } #endif +#if defined(CERES_NO_SUITESPARSE) && defined(CERES_NO_CXSPARSE) + if (linear_solver_options.type == SPARSE_SCHUR) { + *error = "Can't use SPARSE_SCHUR because neither SuiteSparse nor" + "CXSparse was enabled when Ceres was compiled."; + return NULL; + } +#endif if (options->linear_solver_max_num_iterations <= 0) { *error = "Solver::Options::linear_solver_max_num_iterations is 0."; @@ -629,61 +706,8 @@ linear_solver_options.preconditioner_type = options->preconditioner_type; linear_solver_options.sparse_linear_algebra_library = options->sparse_linear_algebra_library; - linear_solver_options.use_block_amd = options->use_block_amd; - -#ifdef CERES_NO_SUITESPARSE - if (linear_solver_options.preconditioner_type == SCHUR_JACOBI) { - *error = "SCHUR_JACOBI preconditioner not suppored. Please build Ceres " - "with SuiteSparse support."; - return NULL; - } - - if (linear_solver_options.preconditioner_type == CLUSTER_JACOBI) { - *error = "CLUSTER_JACOBI preconditioner not suppored. Please build Ceres " - "with SuiteSparse support."; - return NULL; - } - - if (linear_solver_options.preconditioner_type == CLUSTER_TRIDIAGONAL) { - *error = "CLUSTER_TRIDIAGONAL preconditioner not suppored. Please build " - "Ceres with SuiteSparse support."; - return NULL; - } -#endif linear_solver_options.num_threads = options->num_linear_solver_threads; - if (options->num_eliminate_blocks > 0) { - linear_solver_options - .elimination_groups.push_back(options->num_eliminate_blocks); - } - - // TODO(sameeragarwal): Fix this. Right now these are dummy values - // and violate the contract that elimination_groups should sum to - // the number of parameter blocks, but fixing this requires a bit - // more refactoring in solver_impl.cc, which is better done as we - // start deprecating the old API. - linear_solver_options.elimination_groups.push_back(1); - - if (linear_solver_options.elimination_groups.size() == 1 && - IsSchurType(linear_solver_options.type)) { -#if defined(CERES_NO_SUITESPARSE) && defined(CERES_NO_CXSPARSE) - LOG(INFO) << "No elimination block remaining switching to DENSE_QR."; - linear_solver_options.type = DENSE_QR; -#else - LOG(INFO) << "No elimination block remaining " - << "switching to SPARSE_NORMAL_CHOLESKY."; - linear_solver_options.type = SPARSE_NORMAL_CHOLESKY; -#endif - } - -#if defined(CERES_NO_SUITESPARSE) && defined(CERES_NO_CXSPARSE) - if (linear_solver_options.type == SPARSE_SCHUR) { - *error = "Can't use SPARSE_SCHUR because neither SuiteSparse nor" - "CXSparse was enabled when Ceres was compiled."; - return NULL; - } -#endif - // The matrix used for storing the dense Schur complement has a // single lock guarding the whole matrix. Running the // SchurComplementSolver with multiple threads leads to maximum @@ -698,56 +722,68 @@ << "switching to single-threaded."; linear_solver_options.num_threads = 1; } - - options->linear_solver_type = linear_solver_options.type; options->num_linear_solver_threads = linear_solver_options.num_threads; + linear_solver_options.use_block_amd = options->use_block_amd; + const map<int, set<double*> >& groups = + options->ordering->group_id_to_parameter_blocks(); + for (map<int, set<double*> >::const_iterator it = groups.begin(); + it != groups.end(); + ++it) { + linear_solver_options.elimination_groups.push_back(it->second.size()); + } + // Schur type solvers, expect at least two elimination groups. If + // there is only one elimination group, then CreateReducedProblem + // guarantees that this group only contains e_blocks. Thus we add a + // dummy elimination group with zero blocks in it. + if (IsSchurType(linear_solver_options.type) && + linear_solver_options.elimination_groups.size() == 1) { + linear_solver_options.elimination_groups.push_back(0); + } + return LinearSolver::Create(linear_solver_options); } bool SolverImpl::ApplyUserOrdering(const ProblemImpl& problem_impl, - vector<double*>& ordering, + const Ordering* ordering, Program* program, string* error) { - if (ordering.size() != program->NumParameterBlocks()) { + if (ordering->NumParameterBlocks() != program->NumParameterBlocks()) { *error = StringPrintf("User specified ordering does not have the same " "number of parameters as the problem. The problem" - "has %d blocks while the ordering has %ld blocks.", + "has %d blocks while the ordering has %d blocks.", program->NumParameterBlocks(), - ordering.size()); + ordering->NumParameterBlocks()); return false; } - // Ensure that there are no duplicates in the user's ordering. - { - vector<double*> ordering_copy(ordering); - sort(ordering_copy.begin(), ordering_copy.end()); - if (unique(ordering_copy.begin(), ordering_copy.end()) - != ordering_copy.end()) { - *error = "User specified ordering contains duplicates."; - return false; - } - } - vector<ParameterBlock*>* parameter_blocks = program->mutable_parameter_blocks(); + parameter_blocks->clear(); - fill(parameter_blocks->begin(), - parameter_blocks->end(), - static_cast<ParameterBlock*>(NULL)); + const ProblemImpl::ProblemImpl::ParameterMap& parameter_map = + problem_impl.parameter_map(); + const map<int, set<double*> >& groups = + ordering->group_id_to_parameter_blocks(); - const ProblemImpl::ParameterMap& parameter_map = problem_impl.parameter_map(); - for (int i = 0; i < ordering.size(); ++i) { - ProblemImpl::ParameterMap::const_iterator it = - parameter_map.find(ordering[i]); - if (it == parameter_map.end()) { - *error = StringPrintf("User specified ordering contains a pointer " - "to a double that is not a parameter block in the " - "problem. The invalid double is at position %d " - " in options.ordering.", i); - return false; + for (map<int, set<double*> >::const_iterator group_it = groups.begin(); + group_it != groups.end(); + ++group_it) { + const set<double*>& group = group_it->second; + for (set<double*>::const_iterator parameter_block_ptr_it = group.begin(); + parameter_block_ptr_it != group.end(); + ++parameter_block_ptr_it) { + ProblemImpl::ParameterMap::const_iterator parameter_block_it = + parameter_map.find(*parameter_block_ptr_it); + if (parameter_block_it == parameter_map.end()) { + *error = StringPrintf("User specified ordering contains a pointer " + "to a double that is not a parameter block in the " + "problem. The invalid double is in group: %d", + group_it->first); + return false; + } + parameter_blocks->push_back(parameter_block_it->second); } - (*parameter_blocks)[i] = it->second; } return true; } @@ -782,28 +818,26 @@ return true; } - CHECK_NE(0, options.num_eliminate_blocks) - << "Congratulations, you found a Ceres bug! Please report this error " - << "to the developers."; + const int num_eliminate_blocks = + options.ordering->group_id_to_parameter_blocks().begin()->second.size(); // Create a histogram of the number of residuals for each E block. There is an // extra bucket at the end to catch all non-eliminated F blocks. - vector<int> residual_blocks_per_e_block(options.num_eliminate_blocks + 1); + vector<int> residual_blocks_per_e_block(num_eliminate_blocks + 1); vector<ResidualBlock*>* residual_blocks = program->mutable_residual_blocks(); vector<int> min_position_per_residual(residual_blocks->size()); for (int i = 0; i < residual_blocks->size(); ++i) { ResidualBlock* residual_block = (*residual_blocks)[i]; - int position = MinParameterBlock(residual_block, - options.num_eliminate_blocks); + int position = MinParameterBlock(residual_block, num_eliminate_blocks); min_position_per_residual[i] = position; - DCHECK_LE(position, options.num_eliminate_blocks); + DCHECK_LE(position, num_eliminate_blocks); residual_blocks_per_e_block[position]++; } // Run a cumulative sum on the histogram, to obtain offsets to the start of // each histogram bucket (where each bucket is for the residuals for that // E-block). - vector<int> offsets(options.num_eliminate_blocks + 1); + vector<int> offsets(num_eliminate_blocks + 1); std::partial_sum(residual_blocks_per_e_block.begin(), residual_blocks_per_e_block.end(), offsets.begin()); @@ -842,7 +876,7 @@ // Sanity check #1: The difference in bucket offsets should match the // histogram sizes. - for (int i = 0; i < options.num_eliminate_blocks; ++i) { + for (int i = 0; i < num_eliminate_blocks; ++i) { CHECK_EQ(residual_blocks_per_e_block[i], offsets[i + 1] - offsets[i]) << "Congratulations, you found a Ceres bug! Please report this error " << "to the developers."; @@ -864,7 +898,12 @@ string* error) { Evaluator::Options evaluator_options; evaluator_options.linear_solver_type = options.linear_solver_type; - evaluator_options.num_eliminate_blocks = options.num_eliminate_blocks; + + evaluator_options.num_eliminate_blocks = + (options.ordering->NumGroups() > 0 && + IsSchurType(options.linear_solver_type)) + ? options.ordering->group_id_to_parameter_blocks().begin()->second.size() + : 0; evaluator_options.num_threads = options.num_threads; return Evaluator::Create(evaluator_options, program, error); }
diff --git a/internal/ceres/solver_impl.h b/internal/ceres/solver_impl.h index 11b44de..076a371 100644 --- a/internal/ceres/solver_impl.h +++ b/internal/ceres/solver_impl.h
@@ -37,6 +37,8 @@ #include "ceres/solver.h" namespace ceres { +class Ordering; + namespace internal { class Evaluator; @@ -76,7 +78,7 @@ // indicates an error was encountered whose cause is logged to // LOG(ERROR). static bool ApplyUserOrdering(const ProblemImpl& problem_impl, - vector<double*>& ordering, + const Ordering* ordering, Program* program, string* error); @@ -108,9 +110,13 @@ // If fixed_cost is not NULL, the residual blocks that are removed // are evaluated and the sum of their cost is returned in fixed_cost. static bool RemoveFixedBlocksFromProgram(Program* program, - int* num_eliminate_blocks, + Ordering* ordering, double* fixed_cost, string* error); + + static bool IsOrderingValid(const Solver::Options& options, + const ProblemImpl* problem_impl, + string* error); }; } // namespace internal
diff --git a/internal/ceres/solver_impl_test.cc b/internal/ceres/solver_impl_test.cc index 7d7740a..a19fca0 100644 --- a/internal/ceres/solver_impl_test.cc +++ b/internal/ceres/solver_impl_test.cc
@@ -31,6 +31,7 @@ #include "gtest/gtest.h" #include "ceres/autodiff_cost_function.h" #include "ceres/linear_solver.h" +#include "ceres/ordering.h" #include "ceres/parameter_block.h" #include "ceres/problem_impl.h" #include "ceres/program.h" @@ -87,29 +88,19 @@ string error; { - int num_eliminate_blocks = 0; + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 0); + Program program(*problem.mutable_program()); EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, NULL, &error)); EXPECT_EQ(program.NumParameterBlocks(), 3); EXPECT_EQ(program.NumResidualBlocks(), 3); - EXPECT_EQ(num_eliminate_blocks, 0); - } - - // Check that num_eliminate_blocks is preserved, when it contains - // all blocks. - { - int num_eliminate_blocks = 3; - Program program(problem.program()); - EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, - NULL, - &error)); - EXPECT_EQ(program.NumParameterBlocks(), 3); - EXPECT_EQ(program.NumResidualBlocks(), 3); - EXPECT_EQ(num_eliminate_blocks, 3); + EXPECT_EQ(ordering.NumParameterBlocks(), 3); } } @@ -121,16 +112,18 @@ problem.AddResidualBlock(new UnaryCostFunction(), NULL, &x); problem.SetParameterBlockConstant(&x); - int num_eliminate_blocks = 0; + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + Program program(problem.program()); string error; EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, NULL, &error)); EXPECT_EQ(program.NumParameterBlocks(), 0); EXPECT_EQ(program.NumResidualBlocks(), 0); - EXPECT_EQ(num_eliminate_blocks, 0); + EXPECT_EQ(ordering.NumParameterBlocks(), 0); } TEST(SolverImpl, RemoveFixedBlocksNoResidualBlocks) { @@ -143,16 +136,21 @@ problem.AddParameterBlock(&y, 1); problem.AddParameterBlock(&z, 1); - int num_eliminate_blocks = 0; + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 0); + + Program program(problem.program()); string error; EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, NULL, &error)); EXPECT_EQ(program.NumParameterBlocks(), 0); EXPECT_EQ(program.NumResidualBlocks(), 0); - EXPECT_EQ(num_eliminate_blocks, 0); + EXPECT_EQ(ordering.NumParameterBlocks(), 0); } TEST(SolverImpl, RemoveFixedBlocksOneParameterBlockConstant) { @@ -164,20 +162,26 @@ problem.AddParameterBlock(&x, 1); problem.AddParameterBlock(&y, 1); problem.AddParameterBlock(&z, 1); + + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 0); + problem.AddResidualBlock(new UnaryCostFunction(), NULL, &x); problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); problem.SetParameterBlockConstant(&x); - int num_eliminate_blocks = 0; + Program program(problem.program()); string error; EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, NULL, &error)); EXPECT_EQ(program.NumParameterBlocks(), 1); EXPECT_EQ(program.NumResidualBlocks(), 1); - EXPECT_EQ(num_eliminate_blocks, 0); + EXPECT_EQ(ordering.NumParameterBlocks(), 1); } TEST(SolverImpl, RemoveFixedBlocksNumEliminateBlocks) { @@ -194,16 +198,22 @@ problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); problem.SetParameterBlockConstant(&x); - int num_eliminate_blocks = 2; + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 1); + Program program(problem.program()); string error; EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, NULL, &error)); EXPECT_EQ(program.NumParameterBlocks(), 2); EXPECT_EQ(program.NumResidualBlocks(), 2); - EXPECT_EQ(num_eliminate_blocks, 1); + EXPECT_EQ(ordering.NumParameterBlocks(), 2); + EXPECT_EQ(ordering.GroupIdForParameterBlock(&y), 0); + EXPECT_EQ(ordering.GroupIdForParameterBlock(&z), 1); } TEST(SolverImpl, RemoveFixedBlocksFixedCost) { @@ -220,7 +230,11 @@ problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); problem.SetParameterBlockConstant(&x); - int num_eliminate_blocks = 2; + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 1); + double fixed_cost = 0.0; Program program(problem.program()); @@ -231,12 +245,14 @@ string error; EXPECT_TRUE(SolverImpl::RemoveFixedBlocksFromProgram(&program, - &num_eliminate_blocks, + &ordering, &fixed_cost, &error)); EXPECT_EQ(program.NumParameterBlocks(), 2); EXPECT_EQ(program.NumResidualBlocks(), 2); - EXPECT_EQ(num_eliminate_blocks, 1); + EXPECT_EQ(ordering.NumParameterBlocks(), 2); + EXPECT_EQ(ordering.GroupIdForParameterBlock(&y), 0); + EXPECT_EQ(ordering.GroupIdForParameterBlock(&z), 1); EXPECT_DOUBLE_EQ(fixed_cost, expected_fixed_cost); } @@ -253,6 +269,11 @@ problem.AddResidualBlock(new TernaryCostFunction(), NULL, &x, &y, &z); problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 0); + ordering.AddParameterBlockToGroup(&z, 0); + const vector<ResidualBlock*>& residual_blocks = problem.program().residual_blocks(); vector<ResidualBlock*> current_residual_blocks(residual_blocks); @@ -270,31 +291,6 @@ } } -TEST(SolverImpl, ReorderResidualBlockNumEliminateBlockDeathTest) { - ProblemImpl problem; - double x; - double y; - double z; - - problem.AddParameterBlock(&x, 1); - problem.AddParameterBlock(&y, 1); - problem.AddParameterBlock(&z, 1); - problem.AddResidualBlock(new UnaryCostFunction(), NULL, &x); - problem.AddResidualBlock(new TernaryCostFunction(), NULL, &x, &y, &z); - problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); - - Solver::Options options; - options.linear_solver_type = DENSE_SCHUR; - options.num_eliminate_blocks = 0; - string error; -#ifndef _WIN32 - EXPECT_DEATH( - SolverImpl::MaybeReorderResidualBlocks( - options, problem.mutable_program(), &error), - "Congratulations"); -#endif // _WIN32 -} - TEST(SolverImpl, ReorderResidualBlockNormalFunction) { ProblemImpl problem; double x; @@ -312,9 +308,14 @@ problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &y); problem.AddResidualBlock(new UnaryCostFunction(), NULL, &y); + Ordering* ordering = new Ordering; + ordering->AddParameterBlockToGroup(&x, 0); + ordering->AddParameterBlockToGroup(&y, 0); + ordering->AddParameterBlockToGroup(&z, 1); + Solver::Options options; options.linear_solver_type = DENSE_SCHUR; - options.num_eliminate_blocks = 2; + options.ordering = ordering; const vector<ResidualBlock*>& residual_blocks = problem.program().residual_blocks(); @@ -368,9 +369,14 @@ problem.AddResidualBlock(new BinaryCostFunction(), NULL, &x, &z); // 6 x problem.AddResidualBlock(new UnaryCostFunction(), NULL, &y); // 7 + Ordering* ordering = new Ordering; + ordering->AddParameterBlockToGroup(&x, 0); + ordering->AddParameterBlockToGroup(&z, 0); + ordering->AddParameterBlockToGroup(&y, 1); + Solver::Options options; options.linear_solver_type = DENSE_SCHUR; - options.num_eliminate_blocks = 2; + options.ordering = ordering; // Create the reduced program. This should remove the fixed block "z", // marking the index to -1 at the same time. x and y also get indices. @@ -423,42 +429,18 @@ problem.AddParameterBlock(&y, 1); problem.AddParameterBlock(&z, 1); - vector<double*> ordering; - ordering.push_back(&x); - ordering.push_back(&z); + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 1); Program program(problem.program()); string error; EXPECT_FALSE(SolverImpl::ApplyUserOrdering(problem, - ordering, + &ordering, &program, &error)); } -TEST(SolverImpl, ApplyUserOrderingHasDuplicates) { - ProblemImpl problem; - double x; - double y; - double z; - - problem.AddParameterBlock(&x, 1); - problem.AddParameterBlock(&y, 1); - problem.AddParameterBlock(&z, 1); - - vector<double*> ordering; - ordering.push_back(&x); - ordering.push_back(&z); - ordering.push_back(&z); - - Program program(problem.program()); - string error; - EXPECT_FALSE(SolverImpl::ApplyUserOrdering(problem, - ordering, - &program, - &error)); -} - - TEST(SolverImpl, ApplyUserOrderingNormal) { ProblemImpl problem; double x; @@ -469,16 +451,16 @@ problem.AddParameterBlock(&y, 1); problem.AddParameterBlock(&z, 1); - vector<double*> ordering; - ordering.push_back(&x); - ordering.push_back(&z); - ordering.push_back(&y); + Ordering ordering; + ordering.AddParameterBlockToGroup(&x, 0); + ordering.AddParameterBlockToGroup(&y, 2); + ordering.AddParameterBlockToGroup(&z, 1); Program* program = problem.mutable_program(); string error; EXPECT_TRUE(SolverImpl::ApplyUserOrdering(problem, - ordering, + &ordering, program, &error)); const vector<ParameterBlock*>& parameter_blocks = program->parameter_blocks(); @@ -502,6 +484,8 @@ Solver::Options options; options.linear_solver_type = DENSE_QR; options.linear_solver_max_num_iterations = -1; + // CreateLinearSolver assumes a non-empty ordering. + options.ordering = new Ordering; string error; EXPECT_EQ(SolverImpl::CreateLinearSolver(&options, &error), static_cast<LinearSolver*>(NULL)); @@ -511,6 +495,8 @@ Solver::Options options; options.linear_solver_type = DENSE_QR; options.linear_solver_min_num_iterations = -1; + // CreateLinearSolver assumes a non-empty ordering. + options.ordering = new Ordering; string error; EXPECT_EQ(SolverImpl::CreateLinearSolver(&options, &error), static_cast<LinearSolver*>(NULL)); @@ -521,44 +507,24 @@ options.linear_solver_type = DENSE_QR; options.linear_solver_min_num_iterations = 10; options.linear_solver_max_num_iterations = 5; + options.ordering = new Ordering; string error; EXPECT_EQ(SolverImpl::CreateLinearSolver(&options, &error), static_cast<LinearSolver*>(NULL)); } -TEST(SolverImpl, CreateLinearSolverZeroNumEliminateBlocks) { - Solver::Options options; - options.num_eliminate_blocks = 0; - options.linear_solver_type = DENSE_SCHUR; - string error; - scoped_ptr<LinearSolver> solver( - SolverImpl::CreateLinearSolver(&options, &error)); - EXPECT_TRUE(solver != NULL); - -#if defined(CERES_NO_SUITESPARSE) && defined(CERES_NO_CXSPARSE) - EXPECT_EQ(options.linear_solver_type, DENSE_QR); -#else - EXPECT_EQ(options.linear_solver_type, SPARSE_NORMAL_CHOLESKY); -#endif -} - -TEST(SolverImpl, ZeroNumEliminateBlocks) { - Solver::Options options; - options.num_eliminate_blocks = 0; - options.linear_solver_type = DENSE_SCHUR; - options.ordering_type = NATURAL; - ProblemImpl problem; - Solver::Summary summary; - SolverImpl::Solve(options, &problem, &summary); - EXPECT_EQ(summary.termination_type, DID_NOT_RUN); - EXPECT_EQ(summary.error.substr(0,25), "Using a Schur type solver"); -} - TEST(SolverImpl, CreateLinearSolverDenseSchurMultipleThreads) { Solver::Options options; - options.num_eliminate_blocks = 1; options.linear_solver_type = DENSE_SCHUR; options.num_linear_solver_threads = 2; + // The Schur type solvers can only be created with the Ordering + // contains at least one elimination group. + options.ordering = new Ordering; + double x; + double y; + options.ordering->AddParameterBlockToGroup(&x, 0); + options.ordering->AddParameterBlockToGroup(&y, 0); + string error; scoped_ptr<LinearSolver> solver( SolverImpl::CreateLinearSolver(&options, &error)); @@ -570,6 +536,8 @@ TEST(SolverImpl, CreateIterativeLinearSolverForDogleg) { Solver::Options options; options.trust_region_strategy_type = DOGLEG; + // CreateLinearSolver assumes a non-empty ordering. + options.ordering = new Ordering; string error; options.linear_solver_type = ITERATIVE_SCHUR; EXPECT_EQ(SolverImpl::CreateLinearSolver(&options, &error), @@ -584,6 +552,8 @@ Solver::Options options; scoped_ptr<LinearSolver> solver; options.linear_solver_type = DENSE_QR; + // CreateLinearSolver assumes a non-empty ordering. + options.ordering = new Ordering; string error; solver.reset(SolverImpl::CreateLinearSolver(&options, &error)); EXPECT_EQ(options.linear_solver_type, DENSE_QR); @@ -610,14 +580,17 @@ EXPECT_TRUE(solver.get() != NULL); #endif + double x; + double y; + options.ordering->AddParameterBlockToGroup(&x, 0); + options.ordering->AddParameterBlockToGroup(&y, 0); + options.linear_solver_type = DENSE_SCHUR; - options.num_eliminate_blocks = 2; solver.reset(SolverImpl::CreateLinearSolver(&options, &error)); EXPECT_EQ(options.linear_solver_type, DENSE_SCHUR); EXPECT_TRUE(solver.get() != NULL); options.linear_solver_type = SPARSE_SCHUR; - options.num_eliminate_blocks = 2; solver.reset(SolverImpl::CreateLinearSolver(&options, &error)); #if defined(CERES_NO_SUITESPARSE) && defined(CERES_NO_CXSPARSE) @@ -628,7 +601,6 @@ #endif options.linear_solver_type = ITERATIVE_SCHUR; - options.num_eliminate_blocks = 2; solver.reset(SolverImpl::CreateLinearSolver(&options, &error)); EXPECT_EQ(options.linear_solver_type, ITERATIVE_SCHUR); EXPECT_TRUE(solver.get() != NULL);
diff --git a/internal/ceres/system_test.cc b/internal/ceres/system_test.cc index 3dfbc00..4d78de7 100644 --- a/internal/ceres/system_test.cc +++ b/internal/ceres/system_test.cc
@@ -44,6 +44,7 @@ #include <string> #include "ceres/autodiff_cost_function.h" +#include "ceres/ordering.h" #include "ceres/problem.h" #include "ceres/rotation.h" #include "ceres/solver.h" @@ -57,28 +58,30 @@ namespace ceres { namespace internal { +const bool kAutomaticOrdering = true; +const bool kUserOrdering = false; + // Struct used for configuring the solver. struct SolverConfig { SolverConfig(LinearSolverType linear_solver_type, SparseLinearAlgebraLibraryType sparse_linear_algebra_library, - OrderingType ordering_type) + bool use_automatic_ordering) : linear_solver_type(linear_solver_type), sparse_linear_algebra_library(sparse_linear_algebra_library), - ordering_type(ordering_type), + use_automatic_ordering(use_automatic_ordering), preconditioner_type(IDENTITY), num_threads(1) { } SolverConfig(LinearSolverType linear_solver_type, SparseLinearAlgebraLibraryType sparse_linear_algebra_library, - OrderingType ordering_type, - PreconditionerType preconditioner_type, - int num_threads) + bool use_automatic_ordering, + PreconditionerType preconditioner_type) : linear_solver_type(linear_solver_type), sparse_linear_algebra_library(sparse_linear_algebra_library), - ordering_type(ordering_type), + use_automatic_ordering(use_automatic_ordering), preconditioner_type(preconditioner_type), - num_threads(num_threads) { + num_threads(1) { } string ToString() const { @@ -86,14 +89,14 @@ "(%s, %s, %s, %s, %d)", LinearSolverTypeToString(linear_solver_type), SparseLinearAlgebraLibraryTypeToString(sparse_linear_algebra_library), - OrderingTypeToString(ordering_type), + use_automatic_ordering ? "AUTOMATIC" : "USER", PreconditionerTypeToString(preconditioner_type), num_threads); } LinearSolverType linear_solver_type; SparseLinearAlgebraLibraryType sparse_linear_algebra_library; - OrderingType ordering_type; + bool use_automatic_ordering; PreconditionerType preconditioner_type; int num_threads; }; @@ -129,18 +132,14 @@ options.linear_solver_type = config.linear_solver_type; options.sparse_linear_algebra_library = config.sparse_linear_algebra_library; - options.ordering_type = config.ordering_type; options.preconditioner_type = config.preconditioner_type; options.num_threads = config.num_threads; options.num_linear_solver_threads = config.num_threads; options.return_final_residuals = true; - if (options.ordering_type == SCHUR || options.ordering_type == NATURAL) { - options.ordering.clear(); - } - - if (options.ordering_type == SCHUR) { - options.num_eliminate_blocks = 0; + if (config.use_automatic_ordering) { + delete options.ordering; + options.ordering = NULL; } LOG(INFO) << "Running solver configuration: " @@ -279,21 +278,19 @@ sparse_linear_algebra_library, \ ordering)) - CONFIGURE(DENSE_QR, SUITE_SPARSE, NATURAL); - CONFIGURE(DENSE_NORMAL_CHOLESKY, SUITE_SPARSE, NATURAL); - CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, SCHUR); + CONFIGURE(DENSE_QR, SUITE_SPARSE, kAutomaticOrdering); + CONFIGURE(DENSE_NORMAL_CHOLESKY, SUITE_SPARSE, kAutomaticOrdering); + CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, kAutomaticOrdering); #ifndef CERES_NO_SUITESPARSE - CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, NATURAL); - CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, SCHUR); + CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, kAutomaticOrdering); #endif // CERES_NO_SUITESPARSE #ifndef CERES_NO_CXSPARSE - CONFIGURE(SPARSE_NORMAL_CHOLESKY, CX_SPARSE, NATURAL); - CONFIGURE(SPARSE_NORMAL_CHOLESKY, CX_SPARSE, SCHUR); + CONFIGURE(SPARSE_NORMAL_CHOLESKY, CX_SPARSE, kAutomaticOrdering); #endif // CERES_NO_CXSPARSE - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, SCHUR); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kAutomaticOrdering); #undef CONFIGURE @@ -390,16 +387,17 @@ problem_.AddResidualBlock(cost_function, NULL, camera, point); } + options_.ordering = new Ordering; + // The points come before the cameras. for (int i = 0; i < num_points_; ++i) { - options_.ordering.push_back(points + 3 * i); + options_.ordering->AddParameterBlockToGroup(points + 3 * i, 0); } for (int i = 0; i < num_cameras_; ++i) { - options_.ordering.push_back(cameras + 9 * i); + options_.ordering->AddParameterBlockToGroup(cameras + 9 * i, 1); } - options_.num_eliminate_blocks = num_points(); options_.max_num_iterations = 25; options_.function_tolerance = 1e-10; options_.gradient_tolerance = 1e-10; @@ -479,46 +477,43 @@ TEST(SystemTest, BundleAdjustmentProblem) { vector<SolverConfig> configs; -#define CONFIGURE(linear_solver, sparse_linear_algebra_library, ordering, preconditioner, threads) \ +#define CONFIGURE(linear_solver, sparse_linear_algebra_library, ordering, preconditioner) \ configs.push_back(SolverConfig(linear_solver, \ sparse_linear_algebra_library, \ ordering, \ - preconditioner, \ - threads)) + preconditioner)) #ifndef CERES_NO_SUITESPARSE - CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, NATURAL, IDENTITY, 1); - CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, USER, IDENTITY, 1); - CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, SCHUR, IDENTITY, 1); + CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, kAutomaticOrdering, IDENTITY); + CONFIGURE(SPARSE_NORMAL_CHOLESKY, SUITE_SPARSE, kUserOrdering, IDENTITY); - CONFIGURE(SPARSE_SCHUR, SUITE_SPARSE, USER, IDENTITY, 1); - CONFIGURE(SPARSE_SCHUR, SUITE_SPARSE, SCHUR, IDENTITY, 1); + CONFIGURE(SPARSE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, IDENTITY); + CONFIGURE(SPARSE_SCHUR, SUITE_SPARSE, kUserOrdering, IDENTITY); #endif // CERES_NO_SUITESPARSE #ifndef CERES_NO_CXSPARSE - CONFIGURE(SPARSE_SCHUR, CX_SPARSE, USER, IDENTITY, 1); - CONFIGURE(SPARSE_SCHUR, CX_SPARSE, SCHUR, IDENTITY, 1); + CONFIGURE(SPARSE_SCHUR, CX_SPARSE, kAutomaticOrdering, IDENTITY); + CONFIGURE(SPARSE_SCHUR, CX_SPARSE, kUserOrdering, IDENTITY); #endif // CERES_NO_CXSPARSE - CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, USER, IDENTITY, 1); - CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, SCHUR, IDENTITY, 1); + CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, IDENTITY); + CONFIGURE(DENSE_SCHUR, SUITE_SPARSE, kUserOrdering, IDENTITY); - CONFIGURE(CGNR, SUITE_SPARSE, USER, JACOBI, 1); - - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, USER, JACOBI, 1); + CONFIGURE(CGNR, SUITE_SPARSE, kAutomaticOrdering, JACOBI); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kUserOrdering, JACOBI); #ifndef CERES_NO_SUITESPARSE - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, USER, SCHUR_JACOBI, 1); - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, USER, CLUSTER_JACOBI, 1); - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, USER, CLUSTER_TRIDIAGONAL, 1); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kUserOrdering, SCHUR_JACOBI); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kUserOrdering, CLUSTER_JACOBI); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kUserOrdering, CLUSTER_TRIDIAGONAL); #endif // CERES_NO_SUITESPARSE - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, SCHUR, JACOBI, 1); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, JACOBI); #ifndef CERES_NO_SUITESPARSE - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, SCHUR, SCHUR_JACOBI, 1); - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, SCHUR, CLUSTER_JACOBI, 1); - CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, SCHUR, CLUSTER_TRIDIAGONAL, 1); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, SCHUR_JACOBI); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, CLUSTER_JACOBI); + CONFIGURE(ITERATIVE_SCHUR, SUITE_SPARSE, kAutomaticOrdering, CLUSTER_TRIDIAGONAL); #endif // CERES_NO_SUITESPARSE #undef CONFIGURE
diff --git a/internal/ceres/trust_region_minimizer.cc b/internal/ceres/trust_region_minimizer.cc index 5e9fc04..9f93deb 100644 --- a/internal/ceres/trust_region_minimizer.cc +++ b/internal/ceres/trust_region_minimizer.cc
@@ -97,22 +97,14 @@ // moved inside TrustRegionStrategy, its not clear how we dump the // regularization vector/matrix anymore. // - // Doing this right requires either an API change to the - // TrustRegionStrategy and/or how LinearLeastSquares problems are - // stored on disk. + // Also num_eliminate_blocks is not visible to the trust region + // minimizer either. // - // For now, we will just not dump the regularizer. - return (!binary_search(options_.lsqp_iterations_to_dump.begin(), - options_.lsqp_iterations_to_dump.end(), - iteration) || - DumpLinearLeastSquaresProblem(options_.lsqp_dump_directory, - iteration, - options_.lsqp_dump_format_type, - jacobian, - NULL, - residuals, - step, - options_.num_eliminate_blocks)); + // Both of these indicate that this is the wrong place for this + // code, and going forward this should needs fixing/refactoring. + LOG(WARNING) << "Dumping linear least squares problem to disk is " + << "currently broken."; + return true; } void TrustRegionMinimizer::Minimize(const Minimizer::Options& options,
diff --git a/internal/ceres/types.cc b/internal/ceres/types.cc index 74c6af4..225384d 100644 --- a/internal/ceres/types.cc +++ b/internal/ceres/types.cc
@@ -112,24 +112,6 @@ return false; } -const char* OrderingTypeToString(OrderingType ordering_type) { - switch (ordering_type) { - CASESTR(NATURAL); - CASESTR(USER); - CASESTR(SCHUR); - default: - return "UNKNOWN"; - } -} - -bool StringToOrderingType(string value, OrderingType* type) { - UpperCase(&value); - STRENUM(NATURAL); - STRENUM(USER); - STRENUM(SCHUR); - return false; -} - const char* TrustRegionStrategyTypeToString( TrustRegionStrategyType trust_region_strategy_type) { switch (trust_region_strategy_type) {