1 // Ceres Solver - A fast non-linear least squares minimizer
2 // Copyright 2015 Google Inc. All rights reserved.
3 // http://ceres-solver.org/
5 // Redistribution and use in source and binary forms, with or without
6 // modification, are permitted provided that the following conditions are met:
8 // * Redistributions of source code must retain the above copyright notice,
9 // this list of conditions and the following disclaimer.
10 // * Redistributions in binary form must reproduce the above copyright notice,
11 // this list of conditions and the following disclaimer in the documentation
12 // and/or other materials provided with the distribution.
13 // * Neither the name of Google Inc. nor the names of its contributors may be
14 // used to endorse or promote products derived from this software without
15 // specific prior written permission.
17 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 // ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
21 // LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 // SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 // INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 // POSSIBILITY OF SUCH DAMAGE.
29 // Author: sameeragarwal@google.com (Sameer Agarwal)
30 // keir@google.com (Keir Mierle)
32 // The Problem object is used to build and hold least squares problems.
34 #ifndef CERES_PUBLIC_PROBLEM_H_
35 #define CERES_PUBLIC_PROBLEM_H_
42 #include "glog/logging.h"
43 #include "ceres/internal/macros.h"
44 #include "ceres/internal/port.h"
45 #include "ceres/internal/scoped_ptr.h"
46 #include "ceres/types.h"
47 #include "ceres/internal/disable_warnings.h"
54 class LocalParameterization;
63 } // namespace internal
65 // A ResidualBlockId is an opaque handle clients can use to remove residual
66 // blocks from a Problem after adding them.
67 typedef internal::ResidualBlock* ResidualBlockId;
69 // A class to represent non-linear least squares problems. Such
70 // problems have a cost function that is a sum of error terms (known
71 // as "residuals"), where each residual is a function of some subset
72 // of the parameters. The cost function takes the form
75 // SUM --- loss( || r_i1, r_i2,..., r_ik ||^2 ),
80 // r_ij is residual number i, component j; the residual is a
81 // function of some subset of the parameters x1...xk. For
82 // example, in a structure from motion problem a residual
83 // might be the difference between a measured point in an
84 // image and the reprojected position for the matching
85 // camera, point pair. The residual would have two
86 // components, error in x and error in y.
88 // loss(y) is the loss function; for example, squared error or
89 // Huber L1 loss. If loss(y) = y, then the cost function is
90 // non-robustified least squares.
92 // This class is specifically designed to address the important subset
93 // of "sparse" least squares problems, where each component of the
94 // residual depends only on a small number number of parameters, even
95 // though the total number of residuals and parameters may be very
96 // large. This property affords tremendous gains in scale, allowing
97 // efficient solving of large problems that are otherwise
100 // The canonical example of a sparse least squares problem is
101 // "structure-from-motion" (SFM), where the parameters are points and
102 // cameras, and residuals are reprojection errors. Typically a single
103 // residual will depend only on 9 parameters (3 for the point, 6 for
106 // To create a least squares problem, use the AddResidualBlock() and
107 // AddParameterBlock() methods, documented below. Here is an example least
108 // squares problem containing 3 parameter blocks of sizes 3, 4 and 5
109 // respectively and two residual terms of size 2 and 6:
111 // double x1[] = { 1.0, 2.0, 3.0 };
112 // double x2[] = { 1.0, 2.0, 3.0, 5.0 };
113 // double x3[] = { 1.0, 2.0, 3.0, 6.0, 7.0 };
117 // problem.AddResidualBlock(new MyUnaryCostFunction(...), x1);
118 // problem.AddResidualBlock(new MyBinaryCostFunction(...), x2, x3);
120 // Please see cost_function.h for details of the CostFunction object.
121 class CERES_EXPORT Problem {
123 struct CERES_EXPORT Options {
125 : cost_function_ownership(TAKE_OWNERSHIP),
126 loss_function_ownership(TAKE_OWNERSHIP),
127 local_parameterization_ownership(TAKE_OWNERSHIP),
128 enable_fast_removal(false),
129 disable_all_safety_checks(false) {}
131 // These flags control whether the Problem object owns the cost
132 // functions, loss functions, and parameterizations passed into
133 // the Problem. If set to TAKE_OWNERSHIP, then the problem object
134 // will delete the corresponding cost or loss functions on
135 // destruction. The destructor is careful to delete the pointers
136 // only once, since sharing cost/loss/parameterizations is
138 Ownership cost_function_ownership;
139 Ownership loss_function_ownership;
140 Ownership local_parameterization_ownership;
142 // If true, trades memory for faster RemoveResidualBlock() and
143 // RemoveParameterBlock() operations.
145 // By default, RemoveParameterBlock() and RemoveResidualBlock() take time
146 // proportional to the size of the entire problem. If you only ever remove
147 // parameters or residuals from the problem occassionally, this might be
148 // acceptable. However, if you have memory to spare, enable this option to
149 // make RemoveParameterBlock() take time proportional to the number of
150 // residual blocks that depend on it, and RemoveResidualBlock() take (on
151 // average) constant time.
153 // The increase in memory usage is twofold: an additonal hash set per
154 // parameter block containing all the residuals that depend on the parameter
155 // block; and a hash set in the problem containing all residuals.
156 bool enable_fast_removal;
158 // By default, Ceres performs a variety of safety checks when constructing
159 // the problem. There is a small but measurable performance penalty to
160 // these checks, typically around 5% of construction time. If you are sure
161 // your problem construction is correct, and 5% of the problem construction
162 // time is truly an overhead you want to avoid, then you can set
163 // disable_all_safety_checks to true.
165 // WARNING: Do not set this to true, unless you are absolutely sure of what
167 bool disable_all_safety_checks;
170 // The default constructor is equivalent to the
171 // invocation Problem(Problem::Options()).
173 explicit Problem(const Options& options);
177 // Add a residual block to the overall cost function. The cost
178 // function carries with it information about the sizes of the
179 // parameter blocks it expects. The function checks that these match
180 // the sizes of the parameter blocks listed in parameter_blocks. The
181 // program aborts if a mismatch is detected. loss_function can be
182 // NULL, in which case the cost of the term is just the squared norm
185 // The user has the option of explicitly adding the parameter blocks
186 // using AddParameterBlock. This causes additional correctness
187 // checking; however, AddResidualBlock implicitly adds the parameter
188 // blocks if they are not present, so calling AddParameterBlock
189 // explicitly is not required.
191 // The Problem object by default takes ownership of the
192 // cost_function and loss_function pointers. These objects remain
193 // live for the life of the Problem object. If the user wishes to
194 // keep control over the destruction of these objects, then they can
195 // do this by setting the corresponding enums in the Options struct.
197 // Note: Even though the Problem takes ownership of cost_function
198 // and loss_function, it does not preclude the user from re-using
199 // them in another residual block. The destructor takes care to call
200 // delete on each cost_function or loss_function pointer only once,
201 // regardless of how many residual blocks refer to them.
205 // double x1[] = {1.0, 2.0, 3.0};
206 // double x2[] = {1.0, 2.0, 5.0, 6.0};
207 // double x3[] = {3.0, 6.0, 2.0, 5.0, 1.0};
211 // problem.AddResidualBlock(new MyUnaryCostFunction(...), NULL, x1);
212 // problem.AddResidualBlock(new MyBinaryCostFunction(...), NULL, x2, x1);
214 ResidualBlockId AddResidualBlock(
215 CostFunction* cost_function,
216 LossFunction* loss_function,
217 const std::vector<double*>& parameter_blocks);
219 // Convenience methods for adding residuals with a small number of
220 // parameters. This is the common case. Instead of specifying the
221 // parameter block arguments as a vector, list them as pointers.
222 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
223 LossFunction* loss_function,
225 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
226 LossFunction* loss_function,
227 double* x0, double* x1);
228 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
229 LossFunction* loss_function,
230 double* x0, double* x1, double* x2);
231 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
232 LossFunction* loss_function,
233 double* x0, double* x1, double* x2,
235 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
236 LossFunction* loss_function,
237 double* x0, double* x1, double* x2,
238 double* x3, double* x4);
239 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
240 LossFunction* loss_function,
241 double* x0, double* x1, double* x2,
242 double* x3, double* x4, double* x5);
243 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
244 LossFunction* loss_function,
245 double* x0, double* x1, double* x2,
246 double* x3, double* x4, double* x5,
248 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
249 LossFunction* loss_function,
250 double* x0, double* x1, double* x2,
251 double* x3, double* x4, double* x5,
252 double* x6, double* x7);
253 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
254 LossFunction* loss_function,
255 double* x0, double* x1, double* x2,
256 double* x3, double* x4, double* x5,
257 double* x6, double* x7, double* x8);
258 ResidualBlockId AddResidualBlock(CostFunction* cost_function,
259 LossFunction* loss_function,
260 double* x0, double* x1, double* x2,
261 double* x3, double* x4, double* x5,
262 double* x6, double* x7, double* x8,
265 // Add a parameter block with appropriate size to the problem.
266 // Repeated calls with the same arguments are ignored. Repeated
267 // calls with the same double pointer but a different size results
268 // in undefined behaviour.
269 void AddParameterBlock(double* values, int size);
271 // Add a parameter block with appropriate size and parameterization
272 // to the problem. Repeated calls with the same arguments are
273 // ignored. Repeated calls with the same double pointer but a
274 // different size results in undefined behaviour.
275 void AddParameterBlock(double* values,
277 LocalParameterization* local_parameterization);
279 // Remove a parameter block from the problem. The parameterization of the
280 // parameter block, if it exists, will persist until the deletion of the
281 // problem (similar to cost/loss functions in residual block removal). Any
282 // residual blocks that depend on the parameter are also removed, as
283 // described above in RemoveResidualBlock().
285 // If Problem::Options::enable_fast_removal is true, then the
286 // removal is fast (almost constant time). Otherwise, removing a parameter
287 // block will incur a scan of the entire Problem object.
289 // WARNING: Removing a residual or parameter block will destroy the implicit
290 // ordering, rendering the jacobian or residuals returned from the solver
291 // uninterpretable. If you depend on the evaluated jacobian, do not use
292 // remove! This may change in a future release.
293 void RemoveParameterBlock(double* values);
295 // Remove a residual block from the problem. Any parameters that the residual
296 // block depends on are not removed. The cost and loss functions for the
297 // residual block will not get deleted immediately; won't happen until the
298 // problem itself is deleted.
300 // WARNING: Removing a residual or parameter block will destroy the implicit
301 // ordering, rendering the jacobian or residuals returned from the solver
302 // uninterpretable. If you depend on the evaluated jacobian, do not use
303 // remove! This may change in a future release.
304 void RemoveResidualBlock(ResidualBlockId residual_block);
306 // Hold the indicated parameter block constant during optimization.
307 void SetParameterBlockConstant(double* values);
309 // Allow the indicated parameter block to vary during optimization.
310 void SetParameterBlockVariable(double* values);
312 // Returns true if a parameter block is set constant, and false otherwise.
313 bool IsParameterBlockConstant(double* values) const;
315 // Set the local parameterization for one of the parameter blocks.
316 // The local_parameterization is owned by the Problem by default. It
317 // is acceptable to set the same parameterization for multiple
318 // parameters; the destructor is careful to delete local
319 // parameterizations only once. The local parameterization can only
320 // be set once per parameter, and cannot be changed once set.
321 void SetParameterization(double* values,
322 LocalParameterization* local_parameterization);
324 // Get the local parameterization object associated with this
325 // parameter block. If there is no parameterization object
326 // associated then NULL is returned.
327 const LocalParameterization* GetParameterization(double* values) const;
329 // Set the lower/upper bound for the parameter with position "index".
330 void SetParameterLowerBound(double* values, int index, double lower_bound);
331 void SetParameterUpperBound(double* values, int index, double upper_bound);
333 // Number of parameter blocks in the problem. Always equals
334 // parameter_blocks().size() and parameter_block_sizes().size().
335 int NumParameterBlocks() const;
337 // The size of the parameter vector obtained by summing over the
338 // sizes of all the parameter blocks.
339 int NumParameters() const;
341 // Number of residual blocks in the problem. Always equals
342 // residual_blocks().size().
343 int NumResidualBlocks() const;
345 // The size of the residual vector obtained by summing over the
346 // sizes of all of the residual blocks.
347 int NumResiduals() const;
349 // The size of the parameter block.
350 int ParameterBlockSize(const double* values) const;
352 // The size of local parameterization for the parameter block. If
353 // there is no local parameterization associated with this parameter
354 // block, then ParameterBlockLocalSize = ParameterBlockSize.
355 int ParameterBlockLocalSize(const double* values) const;
357 // Is the given parameter block present in this problem or not?
358 bool HasParameterBlock(const double* values) const;
360 // Fills the passed parameter_blocks vector with pointers to the
361 // parameter blocks currently in the problem. After this call,
362 // parameter_block.size() == NumParameterBlocks.
363 void GetParameterBlocks(std::vector<double*>* parameter_blocks) const;
365 // Fills the passed residual_blocks vector with pointers to the
366 // residual blocks currently in the problem. After this call,
367 // residual_blocks.size() == NumResidualBlocks.
368 void GetResidualBlocks(std::vector<ResidualBlockId>* residual_blocks) const;
370 // Get all the parameter blocks that depend on the given residual block.
371 void GetParameterBlocksForResidualBlock(
372 const ResidualBlockId residual_block,
373 std::vector<double*>* parameter_blocks) const;
375 // Get the CostFunction for the given residual block.
376 const CostFunction* GetCostFunctionForResidualBlock(
377 const ResidualBlockId residual_block) const;
379 // Get the LossFunction for the given residual block. Returns NULL
380 // if no loss function is associated with this residual block.
381 const LossFunction* GetLossFunctionForResidualBlock(
382 const ResidualBlockId residual_block) const;
384 // Get all the residual blocks that depend on the given parameter block.
386 // If Problem::Options::enable_fast_removal is true, then
387 // getting the residual blocks is fast and depends only on the number of
388 // residual blocks. Otherwise, getting the residual blocks for a parameter
389 // block will incur a scan of the entire Problem object.
390 void GetResidualBlocksForParameterBlock(
391 const double* values,
392 std::vector<ResidualBlockId>* residual_blocks) const;
394 // Options struct to control Problem::Evaluate.
395 struct EvaluateOptions {
397 : apply_loss_function(true),
401 // The set of parameter blocks for which evaluation should be
402 // performed. This vector determines the order that parameter
403 // blocks occur in the gradient vector and in the columns of the
404 // jacobian matrix. If parameter_blocks is empty, then it is
405 // assumed to be equal to vector containing ALL the parameter
406 // blocks. Generally speaking the parameter blocks will occur in
407 // the order in which they were added to the problem. But, this
408 // may change if the user removes any parameter blocks from the
411 // NOTE: This vector should contain the same pointers as the ones
412 // used to add parameter blocks to the Problem. These parameter
413 // block should NOT point to new memory locations. Bad things will
415 std::vector<double*> parameter_blocks;
417 // The set of residual blocks to evaluate. This vector determines
418 // the order in which the residuals occur, and how the rows of the
419 // jacobian are ordered. If residual_blocks is empty, then it is
420 // assumed to be equal to the vector containing ALL the residual
421 // blocks. Generally speaking the residual blocks will occur in
422 // the order in which they were added to the problem. But, this
423 // may change if the user removes any residual blocks from the
425 std::vector<ResidualBlockId> residual_blocks;
427 // Even though the residual blocks in the problem may contain loss
428 // functions, setting apply_loss_function to false will turn off
429 // the application of the loss function to the output of the cost
430 // function. This is of use for example if the user wishes to
431 // analyse the solution quality by studying the distribution of
432 // residuals before and after the solve.
433 bool apply_loss_function;
438 // Evaluate Problem. Any of the output pointers can be NULL. Which
439 // residual blocks and parameter blocks are used is controlled by
440 // the EvaluateOptions struct above.
442 // Note 1: The evaluation will use the values stored in the memory
443 // locations pointed to by the parameter block pointers used at the
444 // time of the construction of the problem. i.e.,
448 // problem.AddResidualBlock(new MyCostFunction, NULL, &x);
450 // double cost = 0.0;
451 // problem.Evaluate(Problem::EvaluateOptions(), &cost, NULL, NULL, NULL);
453 // The cost is evaluated at x = 1. If you wish to evaluate the
454 // problem at x = 2, then
457 // problem.Evaluate(Problem::EvaluateOptions(), &cost, NULL, NULL, NULL);
459 // is the way to do so.
461 // Note 2: If no local parameterizations are used, then the size of
462 // the gradient vector (and the number of columns in the jacobian)
463 // is the sum of the sizes of all the parameter blocks. If a
464 // parameter block has a local parameterization, then it contributes
465 // "LocalSize" entries to the gradient vector (and the number of
466 // columns in the jacobian).
468 // Note 3: This function cannot be called while the problem is being
469 // solved, for example it cannot be called from an IterationCallback
470 // at the end of an iteration during a solve.
471 bool Evaluate(const EvaluateOptions& options,
473 std::vector<double>* residuals,
474 std::vector<double>* gradient,
475 CRSMatrix* jacobian);
479 friend class Covariance;
480 internal::scoped_ptr<internal::ProblemImpl> problem_impl_;
481 CERES_DISALLOW_COPY_AND_ASSIGN(Problem);
486 #include "ceres/internal/reenable_warnings.h"
488 #endif // CERES_PUBLIC_PROBLEM_H_