Add GradientProblem and GradientProblemSolver. The line search minimizer in Ceres does not require that the problems that is solving is a sum of squares. Over the past year there have been multiple requests to expose this algorithm on its own so that it can be used to solve unconstrained non-linear minimization problems on its own. With this change, a new optimization problem called GradientProblem is introduced which is basically a thin wrapper around a user defined functor that evaluates cost and gradients (FirstOrderFunction) and an optional LocalParameterization. Corresponding to it, a GradientProblemSolver and its associated options and summary structs are introduced too. An example that uses the new API to find the minimum of Rosenbrock's function is also added. Change-Id: I42bf687540da25de991e9bdb00e321239244e8b4
diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index dbbcb81..e26dc9c 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt
@@ -47,6 +47,9 @@ ADD_EXECUTABLE(curve_fitting curve_fitting.cc) TARGET_LINK_LIBRARIES(curve_fitting ceres) +ADD_EXECUTABLE(rosenbrock rosenbrock.cc) +TARGET_LINK_LIBRARIES(rosenbrock ceres) + ADD_EXECUTABLE(curve_fitting_c curve_fitting.c) TARGET_LINK_LIBRARIES(curve_fitting_c ceres) # As this is a C file #including <math.h> we have to explicitly add the math
diff --git a/examples/rosenbrock.cc b/examples/rosenbrock.cc new file mode 100644 index 0000000..da4ee63 --- /dev/null +++ b/examples/rosenbrock.cc
@@ -0,0 +1,74 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2014 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#include "ceres/ceres.h" +#include "glog/logging.h" + +// f(x,y) = (1-x)^2 + 100(y - x^2)^2; +class Rosenbrock : public ceres::FirstOrderFunction { + public: + virtual ~Rosenbrock() {} + + virtual bool Evaluate(const double* parameters, + double* cost, + double* gradient) const { + const double x = parameters[0]; + const double y = parameters[1]; + + cost[0] = (1.0 - x) * (1.0 - x) + 100.0 * (y - x * x) * (y - x * x); + if (gradient != NULL) { + gradient[0] = -2.0 * (1.0 - x) - 200.0 * (y - x * x) * 2.0 * x; + gradient[1] = 200.0 * (y - x * x); + } + return true; + } + + virtual int NumParameters() const { return 2; } +}; + + +int main(int argc, char** argv) { + google::InitGoogleLogging(argv[0]); + + double parameters[2] = {-1.2, 1.0}; + + ceres::GradientProblemSolver::Options options; + options.minimizer_progress_to_stdout = true; + + ceres::GradientProblemSolver::Summary summary; + ceres::GradientProblem problem(new Rosenbrock()); + ceres::Solve(options, problem, parameters, &summary); + + std::cout << summary.FullReport() << "\n"; + std::cout << "Initial x: " << -1.2 << " y: " << 1.0 << "\n"; + std::cout << "Final x: " << parameters[0] + << " y: " << parameters[1] << "\n"; + return 0; +}
diff --git a/include/ceres/ceres.h b/include/ceres/ceres.h index 0b3ae0d..7c8981e 100644 --- a/include/ceres/ceres.h +++ b/include/ceres/ceres.h
@@ -42,6 +42,8 @@ #include "ceres/crs_matrix.h" #include "ceres/dynamic_autodiff_cost_function.h" #include "ceres/dynamic_numeric_diff_cost_function.h" +#include "ceres/gradient_problem.h" +#include "ceres/gradient_problem_solver.h" #include "ceres/iteration_callback.h" #include "ceres/jet.h" #include "ceres/local_parameterization.h"
diff --git a/include/ceres/gradient_problem.h b/include/ceres/gradient_problem.h new file mode 100644 index 0000000..ae3ac4f --- /dev/null +++ b/include/ceres/gradient_problem.h
@@ -0,0 +1,127 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2014 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#ifndef CERES_PUBLIC_GRADIENT_PROBLEM_H_ +#define CERES_PUBLIC_GRADIENT_PROBLEM_H_ + +#include "ceres/internal/macros.h" +#include "ceres/internal/port.h" +#include "ceres/internal/scoped_ptr.h" + +namespace ceres { + +class FirstOrderFunction; +class LocalParameterization; + +// Instances of GradientProblem represent general non-linear +// optimization problems that must be solved using just the value of +// the objective function and its gradient. Unlike the Problem class, +// which can only be used to model non-linear least squares problems, +// instances of GradientProblem not restricted in the form of the +// objective function. +// +// Structurally GradientProblem is a composition of a +// FirstOrderFunction and optionally a LocalParameterization. +// +// The FirstOrderFunction is responsible for evaluating the cost and +// gradient of the objective function. +// +// The LocalParameterization is responsible for going back and forth +// between the ambient space and the local tangent space. (See +// local_parameterization.h for more details). When a +// LocalParameterization is not provided, then the tangent space is +// assumed to coincide with the ambient Euclidean space that the +// gradient vector lives in. +// +// Example usage: +// +// The following demonstrate the problem construction for Rosenbrock's function +// +// f(x,y) = (1-x)^2 + 100(y - x^2)^2; +// +// class Rosenbrock : public ceres::FirstOrderFunction { +// public: +// virtual ~Rosenbrock() {} +// +// virtual bool Evaluate(const double* parameters, +// double* cost, +// double* gradient) const { +// const double x = parameters[0]; +// const double y = parameters[1]; +// +// cost[0] = (1.0 - x) * (1.0 - x) + 100.0 * (y - x * x) * (y - x * x); +// if (gradient != NULL) { +// gradient[0] = -2.0 * (1.0 - x) - 200.0 * (y - x * x) * 2.0 * x; +// gradient[1] = 200.0 * (y - x * x); +// } +// return true; +// }; +// +// virtual int NumParameters() const { return 2; }; +// }; +// +// ceres::GradientProblem problem(new Rosenbrock()); +class CERES_EXPORT GradientProblem { + public: + // Takes ownership of the function. + explicit GradientProblem(FirstOrderFunction* function); + + // Takes ownership of the function and the parameterization. + GradientProblem(FirstOrderFunction* function, + LocalParameterization* parameterization); + + int NumParameters() const; + int NumLocalParameters() const; + + // This call is not thread safe. + bool Evaluate(const double* parameters, double* cost, double* gradient) const; + bool Plus(const double* x, const double* delta, double* x_plus_delta) const; + + private: + internal::scoped_ptr<FirstOrderFunction> function_; + internal::scoped_ptr<LocalParameterization> parameterization_; + internal::scoped_array<double> scratch_; +}; + +// A FirstOrderFunction object implements the evaluation of a function +// and its gradient. +class CERES_EXPORT FirstOrderFunction { + public: + virtual ~FirstOrderFunction() {} + // cost is never NULL. gradient may be null. + virtual bool Evaluate(const double* const parameters, + double* cost, + double* gradient) const = 0; + virtual int NumParameters() const = 0; +}; + +} // namespace ceres + +#endif // CERES_PUBLIC_GRADIENT_PROBLEM_H_
diff --git a/include/ceres/gradient_problem_solver.h b/include/ceres/gradient_problem_solver.h new file mode 100644 index 0000000..484d88e --- /dev/null +++ b/include/ceres/gradient_problem_solver.h
@@ -0,0 +1,365 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2014 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#ifndef CERES_PUBLIC_GRADIENT_PROBLEM_SOLVER_H_ +#define CERES_PUBLIC_GRADIENT_PROBLEM_SOLVER_H_ + +#include <cmath> +#include <string> +#include <vector> +#include "ceres/internal/macros.h" +#include "ceres/internal/port.h" +#include "ceres/iteration_callback.h" +#include "ceres/types.h" +#include "ceres/internal/disable_warnings.h" + +namespace ceres { + +class GradientProblem; + +class CERES_EXPORT GradientProblemSolver { + public: + virtual ~GradientProblemSolver(); + + // The options structure contains, not surprisingly, options that control how + // the solver operates. The defaults should be suitable for a wide range of + // problems; however, better performance is often obtainable with tweaking. + // + // The constants are defined inside types.h + struct CERES_EXPORT Options { + // Default constructor that sets up a generic sparse problem. + Options() { + line_search_direction_type = LBFGS; + line_search_type = WOLFE; + nonlinear_conjugate_gradient_type = FLETCHER_REEVES; + max_lbfgs_rank = 20; + use_approximate_eigenvalue_bfgs_scaling = false; + line_search_interpolation_type = CUBIC; + min_line_search_step_size = 1e-9; + line_search_sufficient_function_decrease = 1e-4; + max_line_search_step_contraction = 1e-3; + min_line_search_step_contraction = 0.6; + max_num_line_search_step_size_iterations = 20; + max_num_line_search_direction_restarts = 5; + line_search_sufficient_curvature_decrease = 0.9; + max_line_search_step_expansion = 10.0; + max_num_iterations = 50; + max_solver_time_in_seconds = 1e9; + num_threads = 1; + function_tolerance = 1e-6; + gradient_tolerance = 1e-10; + logging_type = PER_MINIMIZER_ITERATION; + minimizer_progress_to_stdout = false; + } + + // Returns true if the options struct has a valid + // configuration. Returns false otherwise, and fills in *error + // with a message describing the problem. + bool IsValid(string* error) const; + + // Minimizer options ---------------------------------------- + LineSearchDirectionType line_search_direction_type; + LineSearchType line_search_type; + NonlinearConjugateGradientType nonlinear_conjugate_gradient_type; + + // The LBFGS hessian approximation is a low rank approximation to + // the inverse of the Hessian matrix. The rank of the + // approximation determines (linearly) the space and time + // complexity of using the approximation. Higher the rank, the + // better is the quality of the approximation. The increase in + // quality is however is bounded for a number of reasons. + // + // 1. The method only uses secant information and not actual + // derivatives. + // + // 2. The Hessian approximation is constrained to be positive + // definite. + // + // So increasing this rank to a large number will cost time and + // space complexity without the corresponding increase in solution + // quality. There are no hard and fast rules for choosing the + // maximum rank. The best choice usually requires some problem + // specific experimentation. + // + // For more theoretical and implementation details of the LBFGS + // method, please see: + // + // Nocedal, J. (1980). "Updating Quasi-Newton Matrices with + // Limited Storage". Mathematics of Computation 35 (151): 773–782. + int max_lbfgs_rank; + + // As part of the (L)BFGS update step (BFGS) / right-multiply step (L-BFGS), + // the initial inverse Hessian approximation is taken to be the Identity. + // However, Oren showed that using instead I * \gamma, where \gamma is + // chosen to approximate an eigenvalue of the true inverse Hessian can + // result in improved convergence in a wide variety of cases. Setting + // use_approximate_eigenvalue_bfgs_scaling to true enables this scaling. + // + // It is important to note that approximate eigenvalue scaling does not + // always improve convergence, and that it can in fact significantly degrade + // performance for certain classes of problem, which is why it is disabled + // by default. In particular it can degrade performance when the + // sensitivity of the problem to different parameters varies significantly, + // as in this case a single scalar factor fails to capture this variation + // and detrimentally downscales parts of the jacobian approximation which + // correspond to low-sensitivity parameters. It can also reduce the + // robustness of the solution to errors in the jacobians. + // + // Oren S.S., Self-scaling variable metric (SSVM) algorithms + // Part II: Implementation and experiments, Management Science, + // 20(5), 863-874, 1974. + bool use_approximate_eigenvalue_bfgs_scaling; + + // Degree of the polynomial used to approximate the objective + // function. Valid values are BISECTION, QUADRATIC and CUBIC. + // + // BISECTION corresponds to pure backtracking search with no + // interpolation. + LineSearchInterpolationType line_search_interpolation_type; + + // If during the line search, the step_size falls below this + // value, it is truncated to zero. + double min_line_search_step_size; + + // Line search parameters. + + // Solving the line search problem exactly is computationally + // prohibitive. Fortunately, line search based optimization + // algorithms can still guarantee convergence if instead of an + // exact solution, the line search algorithm returns a solution + // which decreases the value of the objective function + // sufficiently. More precisely, we are looking for a step_size + // s.t. + // + // f(step_size) <= f(0) + sufficient_decrease * f'(0) * step_size + // + double line_search_sufficient_function_decrease; + + // In each iteration of the line search, + // + // new_step_size >= max_line_search_step_contraction * step_size + // + // Note that by definition, for contraction: + // + // 0 < max_step_contraction < min_step_contraction < 1 + // + double max_line_search_step_contraction; + + // In each iteration of the line search, + // + // new_step_size <= min_line_search_step_contraction * step_size + // + // Note that by definition, for contraction: + // + // 0 < max_step_contraction < min_step_contraction < 1 + // + double min_line_search_step_contraction; + + // Maximum number of trial step size iterations during each line search, + // if a step size satisfying the search conditions cannot be found within + // this number of trials, the line search will terminate. + int max_num_line_search_step_size_iterations; + + // Maximum number of restarts of the line search direction algorithm before + // terminating the optimization. Restarts of the line search direction + // algorithm occur when the current algorithm fails to produce a new descent + // direction. This typically indicates a numerical failure, or a breakdown + // in the validity of the approximations used. + int max_num_line_search_direction_restarts; + + // The strong Wolfe conditions consist of the Armijo sufficient + // decrease condition, and an additional requirement that the + // step-size be chosen s.t. the _magnitude_ ('strong' Wolfe + // conditions) of the gradient along the search direction + // decreases sufficiently. Precisely, this second condition + // is that we seek a step_size s.t. + // + // |f'(step_size)| <= sufficient_curvature_decrease * |f'(0)| + // + // Where f() is the line search objective and f'() is the derivative + // of f w.r.t step_size (d f / d step_size). + double line_search_sufficient_curvature_decrease; + + // During the bracketing phase of the Wolfe search, the step size is + // increased until either a point satisfying the Wolfe conditions is + // found, or an upper bound for a bracket containing a point satisfying + // the conditions is found. Precisely, at each iteration of the + // expansion: + // + // new_step_size <= max_step_expansion * step_size. + // + // By definition for expansion, max_step_expansion > 1.0. + double max_line_search_step_expansion; + + // Maximum number of iterations for the minimizer to run for. + int max_num_iterations; + + // Maximum time for which the minimizer should run for. + double max_solver_time_in_seconds; + + // Number of threads used by Ceres for evaluating the cost and + // jacobians. + int num_threads; + + // Minimizer terminates when + // + // (new_cost - old_cost) < function_tolerance * old_cost; + // + double function_tolerance; + + // Minimizer terminates when + // + // max_i |x - Project(Plus(x, -g(x))| < gradient_tolerance + // + // This value should typically be 1e-4 * function_tolerance. + double gradient_tolerance; + + // Logging options --------------------------------------------------------- + + LoggingType logging_type; + + // By default the Minimizer progress is logged to VLOG(1), which + // is sent to STDERR depending on the vlog level. If this flag is + // set to true, and logging_type is not SILENT, the logging output + // is sent to STDOUT. + bool minimizer_progress_to_stdout; + + // If true, the user's parameter blocks are updated at the end of + // every Minimizer iteration, otherwise they are updated when the + // Minimizer terminates. This is useful if, for example, the user + // wishes to visualize the state of the optimization every + // iteration. + bool update_state_every_iteration; + + // Callbacks that are executed at the end of each iteration of the + // Minimizer. An iteration may terminate midway, either due to + // numerical failures or because one of the convergence tests has + // been satisfied. In this case none of the callbacks are + // executed. + + // Callbacks are executed in the order that they are specified in + // this vector. By default, parameter blocks are updated only at + // the end of the optimization, i.e when the Minimizer + // terminates. This behaviour is controlled by + // update_state_every_variable. If the user wishes to have access + // to the update parameter blocks when his/her callbacks are + // executed, then set update_state_every_iteration to true. + // + // The solver does NOT take ownership of these pointers. + vector<IterationCallback*> callbacks; + }; + + struct CERES_EXPORT Summary { + Summary(); + + // A brief one line description of the state of the solver after + // termination. + string BriefReport() const; + + // A full multiline description of the state of the solver after + // termination. + string FullReport() const; + + bool IsSolutionUsable() const; + + // Minimizer summary ------------------------------------------------- + TerminationType termination_type; + + // Reason why the solver terminated. + string message; + + // Cost of the problem (value of the objective function) before + // the optimization. + double initial_cost; + + // Cost of the problem (value of the objective function) after the + // optimization. + double final_cost; + + // IterationSummary for each minimizer iteration in order. + vector<IterationSummary> iterations; + + // Sum total of all time spent inside Ceres when Solve is called. + double total_time_in_seconds; + + // Time (in seconds) spent evaluating the residual vector. + double cost_evaluation_time_in_seconds; + + // Time (in seconds) spent evaluating the jacobian matrix. + double gradient_evaluation_time_in_seconds; + + // Number of parameters in the probem. + int num_parameters; + + // Dimension of the tangent space of the problem. + int num_local_parameters; + + // Type of line search direction used. + LineSearchDirectionType line_search_direction_type; + + // Type of the line search algorithm used. + LineSearchType line_search_type; + + // When performing line search, the degree of the polynomial used + // to approximate the objective function. + LineSearchInterpolationType line_search_interpolation_type; + + // If the line search direction is NONLINEAR_CONJUGATE_GRADIENT, + // then this indicates the particular variant of non-linear + // conjugate gradient used. + NonlinearConjugateGradientType nonlinear_conjugate_gradient_type; + + // If the type of the line search direction is LBFGS, then this + // indicates the rank of the Hessian approximation. + int max_lbfgs_rank; + }; + + // Once a least squares problem has been built, this function takes + // the problem and optimizes it based on the values of the options + // parameters. Upon return, a detailed summary of the work performed + // by the preprocessor, the non-linear minmizer and the linear + // solver are reported in the summary object. + virtual void Solve(const GradientProblemSolver::Options& options, + const GradientProblem& problem, + double* parameters, + GradientProblemSolver::Summary* summary); +}; + +// Helper function which avoids going through the interface. +CERES_EXPORT void Solve(const GradientProblemSolver::Options& options, + const GradientProblem& problem, + double* parameters, + GradientProblemSolver::Summary* summary); + +} // namespace ceres + +#include "ceres/internal/reenable_warnings.h" + +#endif // CERES_PUBLIC_GRADIENT_PROBLEM_SOLVER_H_
diff --git a/internal/ceres/CMakeLists.txt b/internal/ceres/CMakeLists.txt index 9c49503..e4747ea 100644 --- a/internal/ceres/CMakeLists.txt +++ b/internal/ceres/CMakeLists.txt
@@ -64,6 +64,8 @@ evaluator.cc file.cc gradient_checking_cost_function.cc + gradient_problem.cc + gradient_problem_solver.cc implicit_schur_complement.cc incomplete_lq_factorization.cc iterative_schur_complement_solver.cc
diff --git a/internal/ceres/gradient_problem.cc b/internal/ceres/gradient_problem.cc new file mode 100644 index 0000000..8f9a842 --- /dev/null +++ b/internal/ceres/gradient_problem.cc
@@ -0,0 +1,81 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2014 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#include "ceres/gradient_problem.h" +#include "ceres/local_parameterization.h" +#include "glog/logging.h" + +namespace ceres { + +GradientProblem::GradientProblem(FirstOrderFunction* function) + : function_(function), + parameterization_( + new IdentityParameterization(function_->NumParameters())), + scratch_(new double[function_->NumParameters()]) { +} + +GradientProblem::GradientProblem(FirstOrderFunction* function, + LocalParameterization* parameterization) + : function_(function), + parameterization_(parameterization), + scratch_(new double[function_->NumParameters()]) { + CHECK_EQ(function_->NumParameters(), parameterization_->GlobalSize()); +} + +int GradientProblem::NumParameters() const { + return function_->NumParameters(); +} + +int GradientProblem::NumLocalParameters() const { + return parameterization_->LocalSize(); +} + + +bool GradientProblem::Evaluate(const double* parameters, + double* cost, + double* gradient) const { + if (gradient == NULL) { + return function_->Evaluate(parameters, cost, NULL); + } + + return (function_->Evaluate(parameters, cost, scratch_.get()) && + parameterization_->MultiplyByJacobian(parameters, + 1, + scratch_.get(), + gradient)); +} + +bool GradientProblem::Plus(const double* x, + const double* delta, + double* x_plus_delta) const { + return parameterization_->Plus(x, delta, x_plus_delta); +} + +} // namespace ceres
diff --git a/internal/ceres/gradient_problem_evaluator.h b/internal/ceres/gradient_problem_evaluator.h new file mode 100644 index 0000000..20053de --- /dev/null +++ b/internal/ceres/gradient_problem_evaluator.h
@@ -0,0 +1,98 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2013 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#ifndef CERES_INTERNAL_GRADIENT_PROBLEM_EVALUATOR_H_ +#define CERES_INTERNAL_GRADIENT_PROBLEM_EVALUATOR_H_ + +#include <map> +#include <string> + +#include "ceres/evaluator.h" +#include "ceres/execution_summary.h" +#include "ceres/gradient_problem.h" +#include "ceres/internal/port.h" +#include "ceres/wall_time.h" + +namespace ceres { +namespace internal { + +class GradientProblemEvaluator : public Evaluator { + public: + explicit GradientProblemEvaluator(const GradientProblem& problem) + : problem_(problem) {} + virtual ~GradientProblemEvaluator() {} + virtual SparseMatrix* CreateJacobian() const { return NULL; } + virtual bool Evaluate(const EvaluateOptions& evaluate_options, + const double* state, + double* cost, + double* residuals, + double* gradient, + SparseMatrix* jacobian) { + CHECK(jacobian == NULL); + ScopedExecutionTimer total_timer("Evaluator::Total", &execution_summary_); + ScopedExecutionTimer call_type_timer( + gradient == NULL ? "Evaluator::Cost" : "Evaluator::Gradient", + &execution_summary_); + return problem_.Evaluate(state, cost, gradient); + } + + virtual bool Plus(const double* state, + const double* delta, + double* state_plus_delta) const { + return problem_.Plus(state, delta, state_plus_delta); + } + + virtual int NumParameters() const { + return problem_.NumParameters(); + } + + virtual int NumEffectiveParameters() const { + return problem_.NumLocalParameters(); + } + + virtual int NumResiduals() const { return 1; } + + virtual map<string, int> CallStatistics() const { + return execution_summary_.calls(); + } + + virtual map<string, double> TimeStatistics() const { + return execution_summary_.times(); + } + + private: + const GradientProblem& problem_; + ::ceres::internal::ExecutionSummary execution_summary_; +}; + +} // namespace internal +} // namespace ceres + +#endif // CERES_INTERNAL_GRADIENT_PROBLEM_EVALUATOR_H_
diff --git a/internal/ceres/gradient_problem_solver.cc b/internal/ceres/gradient_problem_solver.cc new file mode 100644 index 0000000..4024f4c --- /dev/null +++ b/internal/ceres/gradient_problem_solver.cc
@@ -0,0 +1,270 @@ +// Ceres Solver - A fast non-linear least squares minimizer +// Copyright 2014 Google Inc. All rights reserved. +// http://code.google.com/p/ceres-solver/ +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are met: +// +// * Redistributions of source code must retain the above copyright notice, +// this list of conditions and the following disclaimer. +// * Redistributions in binary form must reproduce the above copyright notice, +// this list of conditions and the following disclaimer in the documentation +// and/or other materials provided with the distribution. +// * Neither the name of Google Inc. nor the names of its contributors may be +// used to endorse or promote products derived from this software without +// specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +// POSSIBILITY OF SUCH DAMAGE. +// +// Author: sameeragarwal@google.com (Sameer Agarwal) + +#include "ceres/gradient_problem_solver.h" + +#include "ceres/callbacks.h" +#include "ceres/gradient_problem.h" +#include "ceres/gradient_problem_evaluator.h" +#include "ceres/internal/eigen.h" +#include "ceres/internal/port.h" +#include "ceres/map_util.h" +#include "ceres/minimizer.h" +#include "ceres/solver.h" +#include "ceres/solver_utils.h" +#include "ceres/stringprintf.h" +#include "ceres/types.h" +#include "ceres/wall_time.h" + +namespace ceres { +using internal::StringPrintf; +using internal::StringAppendF; + +namespace { + +Solver::Options GradientProblemSolverOptionsToSolverOptions( + const GradientProblemSolver::Options& options) { +#define COPY_OPTION(x) solver_options.x = options.x + + Solver::Options solver_options; + solver_options.minimizer_type = LINE_SEARCH; + COPY_OPTION(line_search_direction_type); + COPY_OPTION(line_search_type); + COPY_OPTION(nonlinear_conjugate_gradient_type); + COPY_OPTION(max_lbfgs_rank); + COPY_OPTION(use_approximate_eigenvalue_bfgs_scaling); + COPY_OPTION(line_search_interpolation_type); + COPY_OPTION(min_line_search_step_size); + COPY_OPTION(line_search_sufficient_function_decrease); + COPY_OPTION(max_line_search_step_contraction); + COPY_OPTION(min_line_search_step_contraction); + COPY_OPTION(max_num_line_search_step_size_iterations); + COPY_OPTION(max_num_line_search_direction_restarts); + COPY_OPTION(line_search_sufficient_curvature_decrease); + COPY_OPTION(max_line_search_step_expansion); + COPY_OPTION(max_num_iterations); + COPY_OPTION(max_solver_time_in_seconds); + COPY_OPTION(function_tolerance); + COPY_OPTION(gradient_tolerance); + COPY_OPTION(logging_type); + COPY_OPTION(minimizer_progress_to_stdout); + COPY_OPTION(callbacks); + return solver_options; +#undef COPY_OPTION +} + + +} // namespace + +GradientProblemSolver::~GradientProblemSolver() { +} + +void GradientProblemSolver::Solve(const GradientProblemSolver::Options& options, + const GradientProblem& problem, + double* parameters_ptr, + GradientProblemSolver::Summary* summary) { + using internal::scoped_ptr; + using internal::WallTimeInSeconds; + using internal::Minimizer; + using internal::GradientProblemEvaluator; + using internal::LoggingCallback; + using internal::SetSummaryFinalCost; + + double start_time = WallTimeInSeconds(); + Solver::Options solver_options = + GradientProblemSolverOptionsToSolverOptions(options); + + *CHECK_NOTNULL(summary) = Summary(); + summary->num_parameters = problem.NumParameters(); + summary->num_local_parameters = problem.NumLocalParameters(); + summary->line_search_direction_type = options.line_search_direction_type; // NOLINT + summary->line_search_interpolation_type = options.line_search_interpolation_type; // NOLINT + summary->line_search_type = options.line_search_type; + summary->max_lbfgs_rank = options.max_lbfgs_rank; + summary->nonlinear_conjugate_gradient_type = options.nonlinear_conjugate_gradient_type; // NOLINT + + // Check validity + if (!solver_options.IsValid(&summary->message)) { + LOG(ERROR) << "Terminating: " << summary->message; + return; + } + + // Assuming that the parameter blocks in the program have been + Minimizer::Options minimizer_options; + minimizer_options = Minimizer::Options(solver_options); + minimizer_options.evaluator.reset(new GradientProblemEvaluator(problem)); + + scoped_ptr<IterationCallback> logging_callback; + if (options.logging_type != SILENT) { + logging_callback.reset( + new LoggingCallback(LINE_SEARCH, options.minimizer_progress_to_stdout)); + minimizer_options.callbacks.insert(minimizer_options.callbacks.begin(), + logging_callback.get()); + } + + scoped_ptr<Minimizer> minimizer(Minimizer::Create(LINE_SEARCH)); + Vector solution(problem.NumParameters()); + VectorRef parameters(parameters_ptr, problem.NumParameters()); + solution = parameters; + + Solver::Summary solver_summary; + solver_summary.fixed_cost = 0.0; + solver_summary.preprocessor_time_in_seconds = 0.0; + solver_summary.postprocessor_time_in_seconds = 0.0; + + minimizer->Minimize(minimizer_options, solution.data(), &solver_summary); + + summary->termination_type = solver_summary.termination_type; + summary->message = solver_summary.message; + summary->initial_cost = solver_summary.initial_cost; + summary->final_cost = solver_summary.final_cost; + summary->iterations = solver_summary.iterations; + + if (summary->IsSolutionUsable()) { + parameters = solution; + SetSummaryFinalCost(summary); + } + + const map<string, double>& evaluator_time_statistics = + minimizer_options.evaluator->TimeStatistics(); + summary->cost_evaluation_time_in_seconds = + FindWithDefault(evaluator_time_statistics, "Evaluator::Residual", 0.0); + summary->gradient_evaluation_time_in_seconds = + FindWithDefault(evaluator_time_statistics, "Evaluator::Jacobian", 0.0); + + summary->total_time_in_seconds = WallTimeInSeconds() - start_time; +} + +// Invalid values for most fields, to ensure that we are not +// accidentally reporting default values. +GradientProblemSolver::Summary::Summary() + : termination_type(FAILURE), + message("ceres::GradientProblemSolve was not called."), + initial_cost(-1.0), + final_cost(-1.0), + total_time_in_seconds(-1.0), + cost_evaluation_time_in_seconds(-1.0), + gradient_evaluation_time_in_seconds(-1.0), + num_parameters(-1), + num_local_parameters(-1), + line_search_direction_type(LBFGS), + line_search_type(ARMIJO), + line_search_interpolation_type(BISECTION), + nonlinear_conjugate_gradient_type(FLETCHER_REEVES), + max_lbfgs_rank(-1) { +} + +bool GradientProblemSolver::Summary::IsSolutionUsable() const { + return internal::IsSolutionUsable(*this); +} + +string GradientProblemSolver::Summary::BriefReport() const { + return StringPrintf("Ceres GradientProblemSolver Report: " + "Iterations: %d, " + "Initial cost: %e, " + "Final cost: %e, " + "Termination: %s", + static_cast<int>(iterations.size()), + initial_cost, + final_cost, + TerminationTypeToString(termination_type)); +}; + +string GradientProblemSolver::Summary::FullReport() const { + using internal::VersionString; + + string report = string("\nSolver Summary (v " + VersionString() + ")\n\n"); + + StringAppendF(&report, "Parameters % 25d\n", num_parameters); + if (num_local_parameters != num_parameters) { + StringAppendF(&report, "Local parameters % 25d\n", + num_local_parameters); + } + + string line_search_direction_string; + if (line_search_direction_type == LBFGS) { + line_search_direction_string = StringPrintf("LBFGS (%d)", max_lbfgs_rank); + } else if (line_search_direction_type == NONLINEAR_CONJUGATE_GRADIENT) { + line_search_direction_string = + NonlinearConjugateGradientTypeToString( + nonlinear_conjugate_gradient_type); + } else { + line_search_direction_string = + LineSearchDirectionTypeToString(line_search_direction_type); + } + + StringAppendF(&report, "Line search direction %19s\n", + line_search_direction_string.c_str()); + + const string line_search_type_string = + StringPrintf("%s %s", + LineSearchInterpolationTypeToString( + line_search_interpolation_type), + LineSearchTypeToString(line_search_type)); + StringAppendF(&report, "Line search type %19s\n", + line_search_type_string.c_str()); + StringAppendF(&report, "\n"); + + StringAppendF(&report, "\nCost:\n"); + StringAppendF(&report, "Initial % 30e\n", initial_cost); + if (termination_type != FAILURE && + termination_type != USER_FAILURE) { + StringAppendF(&report, "Final % 30e\n", final_cost); + StringAppendF(&report, "Change % 30e\n", + initial_cost - final_cost); + } + + StringAppendF(&report, "\nMinimizer iterations % 16d\n", + static_cast<int>(iterations.size())); + + StringAppendF(&report, "\nTime (in seconds):\n"); + + StringAppendF(&report, "\n Cost evaluation %23.3f\n", + cost_evaluation_time_in_seconds); + StringAppendF(&report, " Gradient evaluation %23.3f\n", + gradient_evaluation_time_in_seconds); + + StringAppendF(&report, "Total %25.3f\n\n", + total_time_in_seconds); + + StringAppendF(&report, "Termination: %25s (%s)\n", + TerminationTypeToString(termination_type), message.c_str()); + return report; +}; + +void Solve(const GradientProblemSolver::Options& options, + const GradientProblem& problem, + double* parameters, + GradientProblemSolver::Summary* summary) { + GradientProblemSolver solver; + solver.Solve(options, problem, parameters, summary); +} + +} // namespace ceres
diff --git a/jni/Android.mk b/jni/Android.mk index 496956c..a11833f 100644 --- a/jni/Android.mk +++ b/jni/Android.mk
@@ -136,6 +136,8 @@ $(CERES_SRC_PATH)/evaluator.cc \ $(CERES_SRC_PATH)/file.cc \ $(CERES_SRC_PATH)/gradient_checking_cost_function.cc \ + $(CERES_SRC_PATH)/gradient_problem.cc \ + $(CERES_SRC_PATH)/gradient_problem_solver.cc \ $(CERES_SRC_PATH)/implicit_schur_complement.cc \ $(CERES_SRC_PATH)/iterative_schur_complement_solver.cc \ $(CERES_SRC_PATH)/lapack.cc \