diff --git a/docs/source/solving.rst b/docs/source/solving.rst index 80aee8636..d3f2d0faf 100644 --- a/docs/source/solving.rst +++ b/docs/source/solving.rst @@ -1172,29 +1172,32 @@ elimination group [LiSaad]_. #. ``it`` is the time take by the current iteration. #. ``tt`` is the the total time taken by the minimizer. -.. member:: vector Solver::Options::lsqp_iterations_to_dump +.. member:: vector Solver::Options::trust_region_minimizer_iterations_to_dump Default: ``empty`` - List of iterations at which the optimizer should dump the linear - least squares problem to disk. Useful for testing and - benchmarking. If ``empty``, no problems are dumped. + List of iterations at which the trust region minimizer should dump + the trust region problem. Useful for testing and benchmarking. If + ``empty``, no problems are dumped. -.. member:: string Solver::Options::lsqp_dump_directory +.. member:: string Solver::Options::trust_region_problem_dump_directory Default: ``/tmp`` - If :member:`Solver::Options::lsqp_iterations_to_dump` is non-empty, then - this setting determines the directory to which the files containing - the linear least squares problems are written to. + Directory to which the problems should be written to. Should be + non-empty if + :member:`Solver::Options::trust_region_minimizer_iterations_to_dump` is + non-empty and + :member:`Solver::Options::trust_region_problem_dump_format_type` is not + ``CONSOLE``. -.. member:: DumpFormatType Solver::Options::lsqp_dump_format +.. member:: DumpFormatType Solver::Options::trust_region_problem_dump_format Default: ``TEXTFILE`` - The format in which linear least squares problems should be logged - when :member:`Solver::Options::lsqp_iterations_to_dump` is non-empty. - There are three options: + The format in which trust region problems should be logged when + :member:`Solver::Options::trust_region_minimizer_iterations_to_dump` + is non-empty. There are three options: * ``CONSOLE`` prints the linear least squares problem in a human readable format to ``stderr``. The Jacobian is printed as a @@ -1203,22 +1206,25 @@ elimination group [LiSaad]_. problems. * ``PROTOBUF`` Write out the linear least squares problem to the - directory pointed to by :member:`Solver::Options::lsqp_dump_directory` as + directory pointed to by + :member:`Solver::Options::trust_region_problem_dump_directory` as a protocol buffer. ``linear_least_squares_problems.h/cc`` contains routines for loading these problems. For details on the on disk format used, see ``matrix.proto``. The files are named - ``lm_iteration_???.lsqp``. This requires that ``protobuf`` be + ``ceres_solver_iteration_???.bin``. This requires that ``protobuf`` be linked into Ceres Solver. * ``TEXTFILE`` Write out the linear least squares problem to the - directory pointed to by member::`Solver::Options::lsqp_dump_directory` as + directory pointed to by + :member:`Solver::Options::trust_region_problem_dump_directory` as text files which can be read into ``MATLAB/Octave``. The Jacobian is dumped as a text file containing :math:`(i,j,s)` triplets, the vectors :math:`D`, `x` and `f` are dumped as text files containing a list of their values. - A ``MATLAB/Octave`` script called ``lm_iteration_???.m`` is also - output, which can be used to parse and load the problem into memory. + A ``MATLAB/Octave`` script called + ``ceres_solver_iteration_???.m`` is also output, which can be + used to parse and load the problem into memory. .. member:: bool Solver::Options::check_gradients diff --git a/include/ceres/solver.h b/include/ceres/solver.h index 8bf87f9d2..a86a58629 100644 --- a/include/ceres/solver.h +++ b/include/ceres/solver.h @@ -106,8 +106,8 @@ class Solver { inner_iteration_ordering = NULL; logging_type = PER_MINIMIZER_ITERATION; minimizer_progress_to_stdout = false; - lsqp_dump_directory = "/tmp"; - lsqp_dump_format_type = TEXTFILE; + trust_region_problem_dump_directory = "/tmp"; + trust_region_problem_dump_format_type = TEXTFILE; check_gradients = false; gradient_check_relative_precision = 1e-8; numeric_derivative_relative_step_size = 1e-6; @@ -490,14 +490,17 @@ class Solver { // is sent to STDOUT. bool minimizer_progress_to_stdout; - // List of iterations at which the optimizer should dump the - // linear least squares problem to disk. Useful for testing and - // benchmarking. If empty (default), no problems are dumped. - // - // This is ignored if protocol buffers are disabled. - vector lsqp_iterations_to_dump; - string lsqp_dump_directory; - DumpFormatType lsqp_dump_format_type; + // List of iterations at which the minimizer should dump the trust + // region problem. Useful for testing and benchmarking. If empty + // (default), no problems are dumped. + vector trust_region_minimizer_iterations_to_dump; + + // Directory to which the problems should be written to. Should be + // non-empty if trust_region_minimizer_iterations_to_dump is + // non-empty and trust_region_problem_dump_format_type is not + // CONSOLE. + string trust_region_problem_dump_directory; + DumpFormatType trust_region_problem_dump_format_type; // Finite differences options ---------------------------------------------- diff --git a/internal/ceres/dogleg_strategy.cc b/internal/ceres/dogleg_strategy.cc index a330ad2c7..03c0d78c9 100644 --- a/internal/ceres/dogleg_strategy.cc +++ b/internal/ceres/dogleg_strategy.cc @@ -34,6 +34,7 @@ #include "Eigen/Dense" #include "ceres/array_utils.h" #include "ceres/internal/eigen.h" +#include "ceres/linear_least_squares_problems.h" #include "ceres/linear_solver.h" #include "ceres/polynomial.h" #include "ceres/sparse_matrix.h" @@ -127,7 +128,7 @@ TrustRegionStrategy::Summary DoglegStrategy::ComputeStep( ComputeCauchyPoint(jacobian); LinearSolver::Summary linear_solver_summary = - ComputeGaussNewtonStep(jacobian, residuals); + ComputeGaussNewtonStep(per_solve_options, jacobian, residuals); TrustRegionStrategy::Summary summary; summary.residual_norm = linear_solver_summary.residual_norm; @@ -507,6 +508,7 @@ bool DoglegStrategy::FindMinimumOnTrustRegionBoundary(Vector2d* minimum) const { } LinearSolver::Summary DoglegStrategy::ComputeGaussNewtonStep( + const PerSolveOptions& per_solve_options, SparseMatrix* jacobian, const double* residuals) { const int n = jacobian->num_cols(); @@ -561,6 +563,21 @@ LinearSolver::Summary DoglegStrategy::ComputeGaussNewtonStep( solve_options, gauss_newton_step_.data()); + if (per_solve_options.dump_format_type == CONSOLE || + (per_solve_options.dump_format_type != CONSOLE && + !per_solve_options.dump_filename_base.empty())) { + if (!DumpLinearLeastSquaresProblem(per_solve_options.dump_filename_base, + per_solve_options.dump_format_type, + jacobian, + solve_options.D, + residuals, + gauss_newton_step_.data(), + 0)) { + LOG(ERROR) << "Unable to dump trust region problem." + << " Filename base: " << per_solve_options.dump_filename_base; + } + } + if (linear_solver_summary.termination_type == FAILURE || !IsArrayValid(n, gauss_newton_step_.data())) { mu_ *= mu_increase_factor_; diff --git a/internal/ceres/dogleg_strategy.h b/internal/ceres/dogleg_strategy.h index 7131467d6..71c785cc3 100644 --- a/internal/ceres/dogleg_strategy.h +++ b/internal/ceres/dogleg_strategy.h @@ -79,8 +79,10 @@ class DoglegStrategy : public TrustRegionStrategy { typedef Eigen::Matrix Vector2d; typedef Eigen::Matrix Matrix2d; - LinearSolver::Summary ComputeGaussNewtonStep(SparseMatrix* jacobian, - const double* residuals); + LinearSolver::Summary ComputeGaussNewtonStep( + const PerSolveOptions& per_solve_options, + SparseMatrix* jacobian, + const double* residuals); void ComputeCauchyPoint(SparseMatrix* jacobian); void ComputeGradient(SparseMatrix* jacobian, const double* residuals); void ComputeTraditionalDoglegStep(double* step); diff --git a/internal/ceres/levenberg_marquardt_strategy.cc b/internal/ceres/levenberg_marquardt_strategy.cc index 9e6a59e38..6b29244e9 100644 --- a/internal/ceres/levenberg_marquardt_strategy.cc +++ b/internal/ceres/levenberg_marquardt_strategy.cc @@ -34,6 +34,7 @@ #include "Eigen/Core" #include "ceres/array_utils.h" #include "ceres/internal/eigen.h" +#include "ceres/linear_least_squares_problems.h" #include "ceres/linear_solver.h" #include "ceres/sparse_matrix.h" #include "ceres/trust_region_strategy.h" @@ -111,9 +112,24 @@ TrustRegionStrategy::Summary LevenbergMarquardtStrategy::ComputeStep( } else { VectorRef(step, num_parameters) *= -1.0; } - reuse_diagonal_ = true; + if (per_solve_options.dump_format_type == CONSOLE || + (per_solve_options.dump_format_type != CONSOLE && + !per_solve_options.dump_filename_base.empty())) { + if (!DumpLinearLeastSquaresProblem(per_solve_options.dump_filename_base, + per_solve_options.dump_format_type, + jacobian, + solve_options.D, + residuals, + step, + 0)) { + LOG(ERROR) << "Unable to dump trust region problem." + << " Filename base: " << per_solve_options.dump_filename_base; + } + } + + TrustRegionStrategy::Summary summary; summary.residual_norm = linear_solver_summary.residual_norm; summary.num_iterations = linear_solver_summary.num_iterations; diff --git a/internal/ceres/linear_least_squares_problems.cc b/internal/ceres/linear_least_squares_problems.cc index 6c886a1be..d1ee7f0c3 100644 --- a/internal/ceres/linear_least_squares_problems.cc +++ b/internal/ceres/linear_least_squares_problems.cc @@ -574,9 +574,7 @@ LinearLeastSquaresProblem* LinearLeastSquaresProblem3() { } namespace { -bool DumpLinearLeastSquaresProblemToConsole(const string& directory, - int iteration, - const SparseMatrix* A, +bool DumpLinearLeastSquaresProblemToConsole(const SparseMatrix* A, const double* D, const double* b, const double* x, @@ -602,8 +600,7 @@ bool DumpLinearLeastSquaresProblemToConsole(const string& directory, }; #ifndef CERES_NO_PROTOCOL_BUFFERS -bool DumpLinearLeastSquaresProblemToProtocolBuffer(const string& directory, - int iteration, +bool DumpLinearLeastSquaresProblemToProtocolBuffer(const string& filename_base, const SparseMatrix* A, const double* D, const double* b, @@ -632,18 +629,14 @@ bool DumpLinearLeastSquaresProblemToProtocolBuffer(const string& directory, } lsqp.set_num_eliminate_blocks(num_eliminate_blocks); - string format_string = JoinPath(directory, - "lm_iteration_%03d.lsqp"); - string filename = - StringPrintf(format_string.c_str(), iteration); - LOG(INFO) << "Dumping least squares problem for iteration " << iteration - << " to disk. File: " << filename; + + const string filename = filename_base + ".bin"; + LOG(INFO) << "Dumping least squares problem to disk. File: " << filename; WriteStringToFileOrDie(lsqp.SerializeAsString(), filename); return true; } #else -bool DumpLinearLeastSquaresProblemToProtocolBuffer(const string& directory, - int iteration, +bool DumpLinearLeastSquaresProblemToProtocolBuffer(const string& filename_base, const SparseMatrix* A, const double* D, const double* b, @@ -669,31 +662,25 @@ void WriteArrayToFileOrDie(const string& filename, fclose(fptr); } -bool DumpLinearLeastSquaresProblemToTextFile(const string& directory, - int iteration, +bool DumpLinearLeastSquaresProblemToTextFile(const string& filename_base, const SparseMatrix* A, const double* D, const double* b, const double* x, int num_eliminate_blocks) { CHECK_NOTNULL(A); - string format_string = JoinPath(directory, - "lm_iteration_%03d"); - string filename_prefix = - StringPrintf(format_string.c_str(), iteration); - - LOG(INFO) << "writing to: " << filename_prefix << "*"; + LOG(INFO) << "writing to: " << filename_base << "*"; string matlab_script; StringAppendF(&matlab_script, - "function lsqp = lm_iteration_%03d()\n", iteration); + "function lsqp = load_trust_region_problem()\n"); StringAppendF(&matlab_script, "lsqp.num_rows = %d;\n", A->num_rows()); StringAppendF(&matlab_script, "lsqp.num_cols = %d;\n", A->num_cols()); { - string filename = filename_prefix + "_A.txt"; + string filename = filename_base + "_A.txt"; FILE* fptr = fopen(filename.c_str(), "w"); CHECK_NOTNULL(fptr); A->ToTextFile(fptr); @@ -709,34 +696,33 @@ bool DumpLinearLeastSquaresProblemToTextFile(const string& directory, if (D != NULL) { - string filename = filename_prefix + "_D.txt"; + string filename = filename_base + "_D.txt"; WriteArrayToFileOrDie(filename, D, A->num_cols()); StringAppendF(&matlab_script, "lsqp.D = load('%s', '-ascii');\n", filename.c_str()); } if (b != NULL) { - string filename = filename_prefix + "_b.txt"; + string filename = filename_base + "_b.txt"; WriteArrayToFileOrDie(filename, b, A->num_rows()); StringAppendF(&matlab_script, "lsqp.b = load('%s', '-ascii');\n", filename.c_str()); } if (x != NULL) { - string filename = filename_prefix + "_x.txt"; + string filename = filename_base + "_x.txt"; WriteArrayToFileOrDie(filename, x, A->num_cols()); StringAppendF(&matlab_script, "lsqp.x = load('%s', '-ascii');\n", filename.c_str()); } - string matlab_filename = filename_prefix + ".m"; + string matlab_filename = filename_base + ".m"; WriteStringToFileOrDie(matlab_script, matlab_filename); return true; } } // namespace -bool DumpLinearLeastSquaresProblem(const string& directory, - int iteration, +bool DumpLinearLeastSquaresProblem(const string& filename_base, DumpFormatType dump_format_type, const SparseMatrix* A, const double* D, @@ -745,19 +731,16 @@ bool DumpLinearLeastSquaresProblem(const string& directory, int num_eliminate_blocks) { switch (dump_format_type) { case CONSOLE: - return DumpLinearLeastSquaresProblemToConsole(directory, - iteration, - A, D, b, x, + return DumpLinearLeastSquaresProblemToConsole(A, D, b, x, num_eliminate_blocks); case PROTOBUF: return DumpLinearLeastSquaresProblemToProtocolBuffer( - directory, - iteration, + filename_base, A, D, b, x, num_eliminate_blocks); + case TEXTFILE: - return DumpLinearLeastSquaresProblemToTextFile(directory, - iteration, + return DumpLinearLeastSquaresProblemToTextFile(filename_base, A, D, b, x, num_eliminate_blocks); default: diff --git a/internal/ceres/linear_least_squares_problems.h b/internal/ceres/linear_least_squares_problems.h index c76ae91c7..0b550ff9d 100644 --- a/internal/ceres/linear_least_squares_problems.h +++ b/internal/ceres/linear_least_squares_problems.h @@ -73,8 +73,7 @@ LinearLeastSquaresProblem* LinearLeastSquaresProblem3(); // Write the linear least squares problem to disk. The exact format // depends on dump_format_type. -bool DumpLinearLeastSquaresProblem(const string& directory, - int iteration, +bool DumpLinearLeastSquaresProblem(const string& filename_base, DumpFormatType dump_format_type, const SparseMatrix* A, const double* D, diff --git a/internal/ceres/minimizer.h b/internal/ceres/minimizer.h index 18ce64d2d..a2ebe2f1b 100644 --- a/internal/ceres/minimizer.h +++ b/internal/ceres/minimizer.h @@ -73,9 +73,12 @@ class Minimizer { use_nonmonotonic_steps = options.use_nonmonotonic_steps; max_consecutive_nonmonotonic_steps = options.max_consecutive_nonmonotonic_steps; - 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; + trust_region_problem_dump_directory = + options.trust_region_problem_dump_directory; + trust_region_minimizer_iterations_to_dump = + options.trust_region_minimizer_iterations_to_dump; + trust_region_problem_dump_format_type = + options.trust_region_problem_dump_format_type; max_num_consecutive_invalid_steps = options.max_num_consecutive_invalid_steps; min_trust_region_radius = options.min_trust_region_radius; @@ -110,9 +113,9 @@ class Minimizer { bool jacobi_scaling; bool use_nonmonotonic_steps; int max_consecutive_nonmonotonic_steps; - vector lsqp_iterations_to_dump; - DumpFormatType lsqp_dump_format_type; - string lsqp_dump_directory; + vector trust_region_minimizer_iterations_to_dump; + DumpFormatType trust_region_problem_dump_format_type; + string trust_region_problem_dump_directory; int max_num_consecutive_invalid_steps; double min_trust_region_radius; LineSearchDirectionType line_search_direction_type; diff --git a/internal/ceres/solver_impl.cc b/internal/ceres/solver_impl.cc index f0ac2f6c1..3d48700c6 100644 --- a/internal/ceres/solver_impl.cc +++ b/internal/ceres/solver_impl.cc @@ -391,9 +391,13 @@ void SolverImpl::TrustRegionSolve(const Solver::Options& original_options, summary->num_threads_given = original_options.num_threads; summary->num_threads_used = options.num_threads; - if (options.lsqp_iterations_to_dump.size() > 0) { - LOG(WARNING) << "Dumping linear least squares problems to disk is" - " currently broken. Ignoring Solver::Options::lsqp_iterations_to_dump"; + if (options.trust_region_minimizer_iterations_to_dump.size() > 0 && + options.trust_region_problem_dump_format_type != CONSOLE && + options.trust_region_problem_dump_directory.empty()) { + summary->error = + "Solver::Options::trust_region_problem_dump_directory is empty."; + LOG(ERROR) << summary->error; + return; } event_logger.AddEvent("Init"); diff --git a/internal/ceres/trust_region_minimizer.cc b/internal/ceres/trust_region_minimizer.cc index 3206decb2..d6ae0ab1a 100644 --- a/internal/ceres/trust_region_minimizer.cc +++ b/internal/ceres/trust_region_minimizer.cc @@ -41,6 +41,7 @@ #include "Eigen/Core" #include "ceres/array_utils.h" #include "ceres/evaluator.h" +#include "ceres/file.h" #include "ceres/internal/eigen.h" #include "ceres/internal/scoped_ptr.h" #include "ceres/linear_least_squares_problems.h" @@ -70,25 +71,8 @@ void TrustRegionMinimizer::EstimateScale(const SparseMatrix& jacobian, void TrustRegionMinimizer::Init(const Minimizer::Options& options) { options_ = options; - sort(options_.lsqp_iterations_to_dump.begin(), - options_.lsqp_iterations_to_dump.end()); -} - -bool TrustRegionMinimizer::MaybeDumpLinearLeastSquaresProblem( - const int iteration, - const SparseMatrix* jacobian, - const double* residuals, - const double* step) const { - // TODO(sameeragarwal): Since the use of trust_region_radius has - // moved inside TrustRegionStrategy, its not clear how we dump the - // regularization vector/matrix anymore. - // - // Also num_eliminate_blocks is not visible to the trust region - // minimizer either. - // - // Both of these indicate that this is the wrong place for this - // code, and going forward this should needs fixing/refactoring. - return true; + sort(options_.trust_region_minimizer_iterations_to_dump.begin(), + options_.trust_region_minimizer_iterations_to_dump.end()); } void TrustRegionMinimizer::Minimize(const Minimizer::Options& options, @@ -212,33 +196,38 @@ void TrustRegionMinimizer::Minimize(const Minimizer::Options& options, break; } - iteration_summary = IterationSummary(); - iteration_summary = summary->iterations.back(); - iteration_summary.iteration = summary->iterations.back().iteration + 1; - iteration_summary.step_is_valid = false; - iteration_summary.step_is_successful = false; - const double strategy_start_time = WallTimeInSeconds(); TrustRegionStrategy::PerSolveOptions per_solve_options; per_solve_options.eta = options_.eta; + if (find(options_.trust_region_minimizer_iterations_to_dump.begin(), + options_.trust_region_minimizer_iterations_to_dump.end(), + iteration_summary.iteration) != + options_.trust_region_minimizer_iterations_to_dump.end()) { + per_solve_options.dump_format_type = + options_.trust_region_problem_dump_format_type; + per_solve_options.dump_filename_base = + JoinPath(options_.trust_region_problem_dump_directory, + StringPrintf("ceres_solver_iteration_%03d", + iteration_summary.iteration)); + } else { + per_solve_options.dump_format_type = TEXTFILE; + per_solve_options.dump_filename_base.clear(); + } + TrustRegionStrategy::Summary strategy_summary = strategy->ComputeStep(per_solve_options, jacobian, residuals.data(), trust_region_step.data()); + iteration_summary = IterationSummary(); + iteration_summary.iteration = summary->iterations.back().iteration + 1; iteration_summary.step_solver_time_in_seconds = WallTimeInSeconds() - strategy_start_time; iteration_summary.linear_solver_iterations = strategy_summary.num_iterations; - - if (!MaybeDumpLinearLeastSquaresProblem(iteration_summary.iteration, - jacobian, - residuals.data(), - trust_region_step.data())) { - LOG(FATAL) << "Tried writing linear least squares problem: " - << options.lsqp_dump_directory << "but failed."; - } + iteration_summary.step_is_valid = false; + iteration_summary.step_is_successful = false; double model_cost_change = 0.0; if (strategy_summary.termination_type != FAILURE) { diff --git a/internal/ceres/trust_region_strategy.h b/internal/ceres/trust_region_strategy.h index f150594bb..f7fd25a81 100644 --- a/internal/ceres/trust_region_strategy.h +++ b/internal/ceres/trust_region_strategy.h @@ -31,6 +31,8 @@ #ifndef CERES_INTERNAL_TRUST_REGION_STRATEGY_H_ #define CERES_INTERNAL_TRUST_REGION_STRATEGY_H_ +#include +#include "ceres/internal/port.h" #include "ceres/types.h" namespace ceres { @@ -82,8 +84,22 @@ class TrustRegionStrategy { // Per solve options. struct PerSolveOptions { + PerSolveOptions() + : eta(0), + dump_filename_base(""), + dump_format_type(TEXTFILE) { + } + // Forcing sequence for inexact solves. double eta; + + // If non-empty and dump_format_type is not CONSOLE, the trust + // regions strategy will write the linear system to file(s) with + // name starting with dump_filename_base. If dump_format_type is + // CONSOLE then dump_filename_base will be ignored and the linear + // system will be written to the standard error. + string dump_filename_base; + DumpFormatType dump_format_type; }; struct Summary {