Documentation update.

Update modeling.rst and solving.rst to reflect
changes to the API.

Change-Id: Id1a8adfed1486f08e5fd67c5af2d29708a26490c
This commit is contained in:
Sameer Agarwal
2014-04-29 21:40:57 -07:00
parent d48e705022
commit 8f4dcb25f1
6 changed files with 179 additions and 111 deletions
+6 -4
View File
@@ -48,6 +48,12 @@ Bibliography
preconditioning for bundle adjustment**, *In Proceedings of the
IEEE Conference on Computer Vision and Pattern Recognition*, 2012.
.. [Kanzow] C. Kanzow, N. Yamashita and M. Fukushima,
**LevenbergMarquardt methods with strong local convergence
properties for solving nonlinear equations with convex
constraints**, *Journal of Computational and Applied Mathematics*,
177(2):375397, 2005.
.. [Levenberg] K. Levenberg, **A method for the solution of certain
nonlinear problems in least squares**, *Quart. Appl. Math*,
2(2):164168, 1944.
@@ -113,7 +119,3 @@ Bibliography
Levenberg Marquardt Method for Large Sparse Nonlinear Least
Squares**, *Journal of the Australian Mathematical Society Series
B*, 26(4):387403, 1985.
+2 -1
View File
@@ -61,7 +61,8 @@ Getting started
git clone https://ceres-solver.googlesource.com/ceres-solver
* Read the :ref:`chapter-tutorial`, browse :ref:`chapter-modeling` and :ref:`chapter-solving`.
* Read the :ref:`chapter-tutorial`, browse the chapters on the
:ref:`chapter-modeling` API and the :ref:`chapter-solving` API.
* Join the `mailing list
<https://groups.google.com/forum/?fromgroups#!forum/ceres-solver>`_
and ask questions.
+107 -80
View File
@@ -4,45 +4,67 @@
.. _`chapter-modeling`:
============
Modeling API
============
========
Modeling
========
Recall that Ceres solves robustified non-linear least squares problems
of the form
Ceres solver consists of two distinct parts. A modeling API which
provides a rich set of tools to construct an optimization problem one
term at a time and a solver API that controls the minimization
algorithm. This chapter is devoted to the task of modeling
optimization problems using Ceres. :ref:`chapter-solving` discusses
the various ways in which an optimization problem can be solved using
Ceres.
.. math:: \frac{1}{2}\sum_{i=1} \rho_i\left(\left\|f_i\left(x_{i_1}, ... ,x_{i_k}\right)\right\|^2\right).
:label: ceresproblem
Ceres solves robustified bounds constrained non-linear least squares
problems of the form:
The expression
.. math:: :label: ceresproblem
\min_{\mathbf{x}} &\quad \frac{1}{2}\sum_{i}
\rho_i\left(\left\|f_i\left(x_{i_1},
... ,x_{i_k}\right)\right\|^2\right) \\
\text{s.t.} &\quad l_j \le x_j \le u_j
In Ceres parlance, the expression
:math:`\rho_i\left(\left\|f_i\left(x_{i_1},...,x_{i_k}\right)\right\|^2\right)`
is known as a ``ResidualBlock``, where :math:`f_i(\cdot)` is a
:class:`CostFunction` that depends on the parameter blocks
:math:`\left[x_{i_1},... , x_{i_k}\right]`. In most optimization
problems small groups of scalars occur together. For example the three
components of a translation vector and the four components of the
quaternion that define the pose of a camera. We refer to such a group
of small scalars as a ``ParameterBlock``. Of course a
``ParameterBlock`` can just be a single parameter. :math:`\rho_i` is a
:class:`LossFunction`. A :class:`LossFunction` is a scalar function
that is used to reduce the influence of outliers on the solution of
non-linear least squares problems.
is known as a **residual block**, where :math:`f_i(\cdot)` is a
:class:`CostFunction` that depends on the **parameter blocks**
:math:`\left\{x_{i_1},... , x_{i_k}\right\}`.
In this chapter we will describe the various classes that are part of
Ceres Solver's modeling API, and how they can be used to construct an
optimization problem. Once a problem has been constructed, various
methods for solving them will be discussed in
:ref:`chapter-solving`. It is by design that the modeling and the
solving APIs are orthogonal to each other. This enables
switching/tweaking of various solver parameters without having to
touch the problem once it has been successfully modeled.
In most optimization problems small groups of scalars occur
together. For example the three components of a translation vector and
the four components of the quaternion that define the pose of a
camera. We refer to such a group of scalars as a **parameter block**. Of
course a parameter block can be just a single scalar too.
:math:`\rho_i` is a :class:`LossFunction`. A :class:`LossFunction` is
a scalar valued function that is used to reduce the influence of
outliers on the solution of non-linear least squares problems.
:math:`l_j` and :math:`u_j` are lower and upper bounds on the
:math:parameter block `x_j`.
As a special case, when :math:`\rho_i(x) = x`, i.e., the identity
function, and :math:`l_j = -\infty` and :math:`u_j = \infty` we get
the more familiar unconstrained `non-linear least squares problem
<http://en.wikipedia.org/wiki/Non-linear_least_squares>`_.
.. math:: :label: ceresproblemunconstrained
\frac{1}{2}\sum_{i} \left\|f_i\left(x_{i_1}, ... ,x_{i_k}\right)\right\|^2.
:class:`CostFunction`
---------------------
The single biggest task when modeling a problem is specifying the
residuals and their derivatives. This is done using
:class:`CostFunction` objects.
For each term in the objective function, a :class:`CostFunction` is
responsible for computing a vector of residuals and if asked a vector
of Jacobian matrices, i.e., given :math:`\left[x_{i_1}, ... ,
x_{i_k}\right]`, compute the vector
:math:`f_i\left(x_{i_1},...,x_{i_k}\right)` and the matrices
.. math:: J_{ij} = \frac{\partial}{\partial x_{i_j}}f_i\left(x_{i_1},...,x_{i_k}\right),\quad \forall j \in \{i_1,..., i_k\}
.. class:: CostFunction
@@ -61,22 +83,14 @@ residuals and their derivatives. This is done using
void set_num_residuals(int num_residuals);
};
Given parameter blocks :math:`\left[x_{i_1}, ... , x_{i_k}\right]`,
a :class:`CostFunction` is responsible for computing a vector of
residuals and if asked a vector of Jacobian matrices, i.e., given
:math:`\left[x_{i_1}, ... , x_{i_k}\right]`, compute the vector
:math:`f_i\left(x_{i_1},...,x_{i_k}\right)` and the matrices
.. math:: J_{ij} = \frac{\partial}{\partial x_{i_j}}f_i\left(x_{i_1},...,x_{i_k}\right),\quad \forall j \in \{i_1,..., i_k\}
The signature of the :class:`CostFunction` (number and sizes of
input parameter blocks and number of outputs) is stored in
:member:`CostFunction::parameter_block_sizes_` and
:member:`CostFunction::num_residuals_` respectively. User code
inheriting from this class is expected to set these two members
with the corresponding accessors. This information will be verified
by the :class:`Problem` when added with
:func:`Problem::AddResidualBlock`.
The signature of the :class:`CostFunction` (number and sizes of input
parameter blocks and number of outputs) is stored in
:member:`CostFunction::parameter_block_sizes_` and
:member:`CostFunction::num_residuals_` respectively. User code
inheriting from this class is expected to set these two members with
the corresponding accessors. This information will be verified by the
:class:`Problem` when added with :func:`Problem::AddResidualBlock`.
.. function:: bool CostFunction::Evaluate(double const* const* parameters, double* residuals, double** jacobians)
@@ -114,19 +128,6 @@ residuals and their derivatives. This is done using
This can be used to communicate numerical failures in Jacobian
computations for instance.
A more interesting and common use is to impose constraints on the
parameters. If the initial values of the parameter blocks satisfy
the constraints, then returning false whenever the constraints are
not satisfied will prevent the solver from moving into the
infeasible region. This is not a very sophisticated mechanism for
enforcing constraints, but is often good enough for things like
non-negativity constraints.
Note that it is important that the initial values of the parameter
block must be feasible, otherwise the solver will declare a
numerical problem at iteration 0.
:class:`SizedCostFunction`
--------------------------
@@ -177,7 +178,12 @@ residuals and their derivatives. This is done using
int N9 = 0> // Number of parameters in block 9.
class AutoDiffCostFunction : public
SizedCostFunction<kNumResiduals, N0, N1, N2, N3, N4, N5, N6, N7, N8, N9> {
};
public:
explicit AutoDiffCostFunction(CostFunctor* functor);
// Ignore the template parameter kNumResiduals and use
// num_residuals instead.
AutoDiffCostFunction(CostFunctor* functor, int num_residuals);
}
To get an auto differentiated cost function, you must define a
class with a templated ``operator()`` (a functor) that computes the
@@ -189,9 +195,6 @@ residuals and their derivatives. This is done using
The function must write the computed value in the last argument
(the only non-``const`` one) and return true to indicate success.
Please see :class:`CostFunction` for details on how the return
value may be used to impose simple constraints on the parameter
block.
For example, consider a scalar error :math:`e = k - x^\top y`,
where both :math:`x` and :math:`y` are two-dimensional vector
@@ -802,6 +805,8 @@ residuals and their derivatives. This is done using
.. _`section-loss_function`:
:class:`LossFunction`
---------------------
@@ -1228,10 +1233,10 @@ Instances
.. class:: Problem
:class:`Problem` holds the robustified non-linear least squares
problem :eq:`ceresproblem`. To create a least squares problem, use
the :func:`Problem::AddResidualBlock` and
:func:`Problem::AddParameterBlock` methods.
:class:`Problem` holds the robustified bounds constrained
non-linear least squares problem :eq:`ceresproblem`. To create a
least squares problem, use the :func:`Problem::AddResidualBlock`
and :func:`Problem::AddParameterBlock` methods.
For example a problem containing 3 parameter blocks of sizes 3, 4
and 5 respectively and two residual blocks of size 2 and 6:
@@ -1407,54 +1412,76 @@ Instances
parameterizations only once. The local parameterization can only be
set once per parameter, and cannot be changed once set.
.. function:: int Problem::NumParameterBlocks() const
.. function LocalParameterization* Problem::GetParameterization(double* values) const;
Get the local parameterization object associated with this
parameter block. If there is no parameterization object associated
then `NULL` is returned
.. function void Problem::SetParameterLowerBound(double* values, int index, double lower_bound);
Set the lower bound for the parameter at position `index` in the
parameter block corresponding to `values`. By default the lower
bound is :math:`-\infty`.
.. function void Problem::SetParameterUpperBound(double* values, int index, double upper_bound);
Set the upper bound for the parameter at position `index` in the
parameter block corresponding to `values`. By default the value is
:math:`\infty`.
.. function int Problem::NumParameterBlocks() const
Number of parameter blocks in the problem. Always equals
parameter_blocks().size() and parameter_block_sizes().size().
.. function:: int Problem::NumParameters() const
.. function int Problem::NumParameters() const
The size of the parameter vector obtained by summing over the sizes
of all the parameter blocks.
.. function:: int Problem::NumResidualBlocks() const
.. function int Problem::NumResidualBlocks() const
Number of residual blocks in the problem. Always equals
residual_blocks().size().
.. function:: int Problem::NumResiduals() const
.. function int Problem::NumResiduals() const
The size of the residual vector obtained by summing over the sizes
of all of the residual blocks.
.. function:: int Problem::ParameterBlockSize(const double* values) const
.. function int Problem::ParameterBlockSize(const double* values) const
The size of the parameter block.
.. function:: int Problem::ParameterBlockLocalSize(const double* values) const
.. function int Problem::ParameterBlockLocalSize(const double* values) const
The size of local parameterization for the parameter block. If
there is no local parameterization associated with this parameter
block, then ``ParameterBlockLocalSize`` = ``ParameterBlockSize``.
The size of local parameterization for the parameter block. If
there is no local parameterization associated with this parameter
block, then ``ParameterBlockLocalSize`` = ``ParameterBlockSize``.
.. function:: void Problem::GetParameterBlocks(vector<double*>* parameter_blocks) const
.. function bool Problem::HasParameterBlock(const double* values) const;
Is the give parameter block present in the problem or not?
.. function void Problem::GetParameterBlocks(vector<double*>* parameter_blocks) const
Fills the passed ``parameter_blocks`` vector with pointers to the
parameter blocks currently in the problem. After this call,
``parameter_block.size() == NumParameterBlocks``.
.. function:: void Problem::GetResidualBlocks(vector<ResidualBlockId>* residual_blocks) const
.. function void Problem::GetResidualBlocks(vector<ResidualBlockId>* residual_blocks) const
Fills the passed `residual_blocks` vector with pointers to the
residual blocks currently in the problem. After this call,
`residual_blocks.size() == NumResidualBlocks`.
.. function:: void Problem::GetParameterBlocksForResidualBlock(const ResidualBlockId residual_block, vector<double*>* parameter_blocks) const
.. function void Problem::GetParameterBlocksForResidualBlock(const ResidualBlockId residual_block, vector<double*>* parameter_blocks) const
Get all the parameter blocks that depend on the given residual
block.
.. function:: void Problem::GetResidualBlocksForParameterBlock(const double* values, vector<ResidualBlockId>* residual_blocks) const
.. function void Problem::GetResidualBlocksForParameterBlock(const double* values, vector<ResidualBlockId>* residual_blocks) const
Get all the residual blocks that depend on the given parameter
block.
@@ -1465,7 +1492,7 @@ Instances
blocks for a parameter block will incur a scan of the entire
:class:`Problem` object.
.. function:: bool Problem::Evaluate(const Problem::EvaluateOptions& options, double* cost, vector<double>* residuals, vector<double>* gradient, CRSMatrix* jacobian)
.. function bool Problem::Evaluate(const Problem::EvaluateOptions& options, double* cost, vector<double>* residuals, vector<double>* gradient, CRSMatrix* jacobian)
Evaluate a :class:`Problem`. Any of the output pointers can be
`NULL`. Which residual blocks and parameter blocks are used is
+57 -21
View File
@@ -5,9 +5,9 @@
.. _chapter-solving:
===========
Solving API
===========
=======
Solving
=======
Introduction
============
@@ -23,7 +23,8 @@ variables, and
:math:`m`-dimensional function of :math:`x`. We are interested in
solving the following optimization problem [#f1]_ .
.. math:: \arg \min_x \frac{1}{2}\|F(x)\|^2\ .
.. math:: \arg \min_x \frac{1}{2}\|F(x)\|^2\ . \\
L \le x \le U
:label: nonlinsq
Here, the Jacobian :math:`J(x)` of :math:`F(x)` is an :math:`m\times
@@ -80,15 +81,20 @@ Trust Region Methods
The basic trust region algorithm looks something like this.
1. Given an initial point :math:`x` and a trust region radius :math:`\mu`.
2. :math:`\arg \min_{\Delta x} \frac{1}{2}\|J(x)\Delta
x + F(x)\|^2` s.t. :math:`\|D(x)\Delta x\|^2 \le \mu`
2. Solve
.. math::
\arg \min_{\Delta x}& \frac{1}{2}\|J(x)\Delta x + F(x)\|^2 \\
\text{such that} &\|D(x)\Delta x\|^2 \le \mu\\
&L \le x + \Delta x \le U.
3. :math:`\rho = \frac{\displaystyle \|F(x + \Delta x)\|^2 -
\|F(x)\|^2}{\displaystyle \|J(x)\Delta x + F(x)\|^2 -
\|F(x)\|^2}`
4. if :math:`\rho > \epsilon` then :math:`x = x + \Delta x`.
5. if :math:`\rho > \eta_1` then :math:`\rho = 2 \rho`
6. else if :math:`\rho < \eta_2` then :math:`\rho = 0.5 * \rho`
7. Goto 2.
7. Go to 2.
Here, :math:`\mu` is the trust region radius, :math:`D(x)` is some
matrix used to define a metric on the domain of :math:`F(x)` and
@@ -102,21 +108,27 @@ in the value of :math:`\rho`.
The key computational step in a trust-region algorithm is the solution
of the constrained optimization problem
.. math:: \arg\min_{\Delta x} \frac{1}{2}\|J(x)\Delta x + F(x)\|^2\quad \text{such that}\quad \|D(x)\Delta x\|^2 \le \mu
.. math::
\arg \min_{\Delta x}& \frac{1}{2}\|J(x)\Delta x + F(x)\|^2 \\
\text{such that} &\|D(x)\Delta x\|^2 \le \mu\\
&L \le x + \Delta x \le U.
:label: trp
There are a number of different ways of solving this problem, each
giving rise to a different concrete trust-region algorithm. Currently
Ceres, implements two trust-region algorithms - Levenberg-Marquardt
and Dogleg. The user can choose between them by setting
:member:`Solver::Options::trust_region_strategy_type`.
and Dogleg, each of which is augmented with a line search if bounds
constrained are present [Kanzow]_. The user can choose between them by
setting :member:`Solver::Options::trust_region_strategy_type`.
.. rubric:: Footnotes
.. [#f1] At the level of the non-linear solver, the block
structure is not relevant, therefore our discussion here is
in terms of an optimization problem defined over a state
vector of size :math:`n`.
.. [#f1] At the level of the non-linear solver, the block structure is
not relevant, therefore our discussion here is in terms of an
optimization problem defined over a state vector of size
:math:`n`. Similarly the presence of loss functions is also
ignored as the problem is internally converted into a pure
non-linear least squares problem.
.. _section-levenberg-marquardt:
@@ -346,10 +358,9 @@ region strategies.
Line Search Methods
===================
**The implementation of line search algorithms in Ceres Solver is
fairly new and not very well tested, so for now this part of the
solver should be considered beta quality. We welcome reports of your
experiences both good and bad on the mailinglist.**
The line search method in Ceres Solver cannot handle bounds
constraints right now, so it can only be used for solving
unconstrained problems.
Line search algorithms
@@ -1208,7 +1219,7 @@ elimination group [LiSaad]_.
Number of threads used by the linear solver.
.. member:: ParameterBlockOrdering* Solver::Options::linear_solver_ordering
.. member:: shared_ptr<ParameterBlockOrdering> Solver::Options::linear_solver_ordering
Default: ``NULL``
@@ -1245,6 +1256,22 @@ elimination group [LiSaad]_.
expense of an extra copy of the Jacobian matrix. Setting
``use_postordering`` to ``true`` enables this tradeoff.
.. member:: bool Solver::Options::dynamic_sparsity
Some non-linear least squares problems are symbolically dense but
numerically sparse. i.e. at any given state only a small number of
Jacobian entries are non-zero, but the position and number of
non-zeros is different depending on the state. For these problems
it can be useful to factorize the sparse jacobian at each solver
iteration instead of including all of the zero entries in a single
general factorization.
If your problem does not have this property (or you do not know),
then it is probably best to keep this false, otherwise it will
likely lead to worse performance.
This settings affects the `SPARSE_NORMAL_CHOLESKY` solver.
.. member:: int Solver::Options::min_linear_solver_iterations
Default: ``1``
@@ -1304,7 +1331,7 @@ elimination group [LiSaad]_.
inner iterations in subsequent trust region minimizer iterations is
disabled.
.. member:: ParameterBlockOrdering* Solver::Options::inner_iteration_ordering
.. member:: shared_ptr<ParameterBlockOrdering> Solver::Options::inner_iteration_ordering
Default: ``NULL``
@@ -1829,6 +1856,16 @@ The three arrays will be:
A full multiline description of the state of the solver after
termination.
.. function:: bool Solver::Summary::IsSolutionUsable() const
Whether the solution returned by the optimization algorithm can be
relied on to be numerically sane. This will be the case if
`Solver::Summary:termination_type` is set to `CONVERGENCE`,
`USER_SUCCESS` or `NO_CONVERGENCE`, i.e., either the solver
converged by meeting one of the convergence tolerances or because
the user indicated that it had converged or it ran to the maximum
number of iterations or time.
.. member:: MinimizerType Solver::Summary::minimizer_type
Type of minimization algorithm used.
@@ -1841,7 +1878,6 @@ The three arrays will be:
Reason why the solver terminated.
.. member:: double Solver::Summary::initial_cost
Cost of the problem (value of the objective function) before the
+6 -4
View File
@@ -10,7 +10,7 @@ HEAD
New Features
------------
#. Parameters can now have upper and/or lower bounds when using the
#. Support for bounds constrained optimization problems when using the
trust region minimizer.
#. Problems in which the sparsity structure of the Jacobian changes
over the course of the optimization can now be solved much more
@@ -34,8 +34,8 @@ New Features
parameter block is already present in the problem.
#. ``Problem::GetParameterization`` can be used to access the
parameterization associated with a parameter block.
#. Added the (2,4,9) template specialization for PartitionedMatrixView
and SchurEliminator.
#. Added the (2,4,9) and (2,4,8) template specialization for
PartitionedMatrixView and SchurEliminator.
#. An example demonstrating the use of
DynamicAutoDiffCostFunction. (Joydeep Biswas)
#. An example demonstrating the use of dynamic sparsity (Richard
@@ -92,6 +92,7 @@ Backward Incompatible API Changes
Bug Fixes
---------
#. Fixed errant verbose levels (Bjorn Piltz)
#. Variety of code cleanups, optimizations and bug fixes to the line
search minimizer code (Alex Stewart)
#. Fixed ``BlockSparseMatrix::Transpose`` when the matrix has row and
@@ -121,7 +122,8 @@ Bug Fixes
Schönberger)
#. METIS_FOUND is never set. Changed the commit to fit the setting of
the other #._FOUND definitions. (Andreas Franek)
#. Variety of bug fixes to the ``CMake`` build system (Alex Stewart)
#. Variety of bug fixes and cleanups to the ``CMake`` build system
(Alex Stewart)
#. Removed fictious shared library target from the NDK build.
#. Solver::Options now uses ``shared_ptr`` to handle ownership of
``Solver::Options::linear_solver_ordering`` and
+1 -1
View File
@@ -365,7 +365,7 @@ class CERES_EXPORT Solver {
// Minimizer terminates when
//
// max_i |gradient_i| < gradient_tolerance * max_i|initial_gradient_i|
// max_i |x - Project(Plus(x, -g(x))| < gradient_tolerance
//
// This value should typically be 1e-4 * function_tolerance.
double gradient_tolerance;