Enable support for dumping trust region minimizer problems.
This support was broken due to the TrustRegionMinimizer refactoring.
It is now enabled again, with the responsibilty for dumping the
problem shifted to the individual TrustRegionStrategy.
There is however one wrinkle, which is perhaps an indication of
poor design to start with. The LinearLeastSquaresProblemProto
carries in it num_eliminate_blocks, something which does not
exist anymore. More importantly, the TrustRegionStrategy does not
have access to this quantity anymore.
Dealing with this will be the subject of a future change.
Change-Id: I358adf6a2e386f4940b617bf950d6c7e87d2635d
diff --git a/internal/ceres/dogleg_strategy.cc b/internal/ceres/dogleg_strategy.cc
index a330ad2..03c0d78 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 @@
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 @@
}
LinearSolver::Summary DoglegStrategy::ComputeGaussNewtonStep(
+ const PerSolveOptions& per_solve_options,
SparseMatrix* jacobian,
const double* residuals) {
const int n = jacobian->num_cols();
@@ -561,6 +563,21 @@
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 7131467..71c785c 100644
--- a/internal/ceres/dogleg_strategy.h
+++ b/internal/ceres/dogleg_strategy.h
@@ -79,8 +79,10 @@
typedef Eigen::Matrix<double, 2, 1, Eigen::DontAlign> Vector2d;
typedef Eigen::Matrix<double, 2, 2, Eigen::DontAlign> 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 9e6a59e..6b29244 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 @@
} 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 6c886a1..d1ee7f0 100644
--- a/internal/ceres/linear_least_squares_problems.cc
+++ b/internal/ceres/linear_least_squares_problems.cc
@@ -574,9 +574,7 @@
}
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 @@
};
#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 @@
}
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 @@
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 @@
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 @@
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 c76ae91..0b550ff 100644
--- a/internal/ceres/linear_least_squares_problems.h
+++ b/internal/ceres/linear_least_squares_problems.h
@@ -73,8 +73,7 @@
// 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 18ce64d..a2ebe2f 100644
--- a/internal/ceres/minimizer.h
+++ b/internal/ceres/minimizer.h
@@ -73,9 +73,12 @@
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 @@
bool jacobi_scaling;
bool use_nonmonotonic_steps;
int max_consecutive_nonmonotonic_steps;
- vector<int> lsqp_iterations_to_dump;
- DumpFormatType lsqp_dump_format_type;
- string lsqp_dump_directory;
+ vector<int> 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 f0ac2f6..3d48700 100644
--- a/internal/ceres/solver_impl.cc
+++ b/internal/ceres/solver_impl.cc
@@ -391,9 +391,13 @@
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 3206dec..d6ae0ab 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::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 @@
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 f150594..f7fd25a 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 <string>
+#include "ceres/internal/port.h"
#include "ceres/types.h"
namespace ceres {
@@ -82,8 +84,22 @@
// 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 {
