1 // Ceres Solver - A fast non-linear least squares minimizer 2 // Copyright 2012 Google Inc. All rights reserved. 3 // http://code.google.com/p/ceres-solver/ 4 // 5 // Redistribution and use in source and binary forms, with or without 6 // modification, are permitted provided that the following conditions are met: 7 // 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. 16 // 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. 28 // 29 // Author: moll.markus (at) arcor.de (Markus Moll) 30 // sameeragarwal (at) google.com (Sameer Agarwal) 31 32 #ifndef CERES_INTERNAL_POLYNOMIAL_SOLVER_H_ 33 #define CERES_INTERNAL_POLYNOMIAL_SOLVER_H_ 34 35 #include <vector> 36 #include "ceres/internal/eigen.h" 37 #include "ceres/internal/port.h" 38 39 namespace ceres { 40 namespace internal { 41 42 // All polynomials are assumed to be the form 43 // 44 // sum_{i=0}^N polynomial(i) x^{N-i}. 45 // 46 // and are given by a vector of coefficients of size N + 1. 47 48 // Evaluate the polynomial at x using the Horner scheme. 49 inline double EvaluatePolynomial(const Vector& polynomial, double x) { 50 double v = 0.0; 51 for (int i = 0; i < polynomial.size(); ++i) { 52 v = v * x + polynomial(i); 53 } 54 return v; 55 } 56 57 // Use the companion matrix eigenvalues to determine the roots of the 58 // polynomial. 59 // 60 // This function returns true on success, false otherwise. 61 // Failure indicates that the polynomial is invalid (of size 0) or 62 // that the eigenvalues of the companion matrix could not be computed. 63 // On failure, a more detailed message will be written to LOG(ERROR). 64 // If real is not NULL, the real parts of the roots will be returned in it. 65 // Likewise, if imaginary is not NULL, imaginary parts will be returned in it. 66 bool FindPolynomialRoots(const Vector& polynomial, 67 Vector* real, 68 Vector* imaginary); 69 70 // Return the derivative of the given polynomial. It is assumed that 71 // the input polynomial is at least of degree zero. 72 Vector DifferentiatePolynomial(const Vector& polynomial); 73 74 // Find the minimum value of the polynomial in the interval [x_min, 75 // x_max]. The minimum is obtained by computing all the roots of the 76 // derivative of the input polynomial. All real roots within the 77 // interval [x_min, x_max] are considered as well as the end points 78 // x_min and x_max. Since polynomials are differentiable functions, 79 // this ensures that the true minimum is found. 80 void MinimizePolynomial(const Vector& polynomial, 81 double x_min, 82 double x_max, 83 double* optimal_x, 84 double* optimal_value); 85 86 // Structure for storing sample values of a function. 87 // 88 // Clients can use this struct to communicate the value of the 89 // function and or its gradient at a given point x. 90 struct FunctionSample { 91 FunctionSample() 92 : x(0.0), 93 value(0.0), 94 value_is_valid(false), 95 gradient(0.0), 96 gradient_is_valid(false) { 97 } 98 99 double x; 100 double value; // value = f(x) 101 bool value_is_valid; 102 double gradient; // gradient = f'(x) 103 bool gradient_is_valid; 104 }; 105 106 // Given a set of function value and/or gradient samples, find a 107 // polynomial whose value and gradients are exactly equal to the ones 108 // in samples. 109 // 110 // Generally speaking, 111 // 112 // degree = # values + # gradients - 1 113 // 114 // Of course its possible to sample a polynomial any number of times, 115 // in which case, generally speaking the spurious higher order 116 // coefficients will be zero. 117 Vector FindInterpolatingPolynomial(const vector<FunctionSample>& samples); 118 119 // Interpolate the function described by samples with a polynomial, 120 // and minimize it on the interval [x_min, x_max]. Depending on the 121 // input samples, it is possible that the interpolation or the root 122 // finding algorithms may fail due to numerical difficulties. But the 123 // function is guaranteed to return its best guess of an answer, by 124 // considering the samples and the end points as possible solutions. 125 void MinimizeInterpolatingPolynomial(const vector<FunctionSample>& samples, 126 double x_min, 127 double x_max, 128 double* optimal_x, 129 double* optimal_value); 130 131 } // namespace internal 132 } // namespace ceres 133 134 #endif // CERES_INTERNAL_POLYNOMIAL_SOLVER_H_ 135
