1 // Ceres Solver - A fast non-linear least squares minimizer 2 // Copyright 2010, 2011, 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: keir@google.com (Keir Mierle) 30 31 #ifndef CERES_INTERNAL_PROGRAM_H_ 32 #define CERES_INTERNAL_PROGRAM_H_ 33 34 #include <set> 35 #include <string> 36 #include <vector> 37 #include "ceres/internal/port.h" 38 39 namespace ceres { 40 namespace internal { 41 42 class ParameterBlock; 43 class ProblemImpl; 44 class ResidualBlock; 45 class TripletSparseMatrix; 46 47 // A nonlinear least squares optimization problem. This is different from the 48 // similarly-named "Problem" object, which offers a mutation interface for 49 // adding and modifying parameters and residuals. The Program contains the core 50 // part of the Problem, which is the parameters and the residuals, stored in a 51 // particular ordering. The ordering is critical, since it defines the mapping 52 // between (residual, parameter) pairs and a position in the jacobian of the 53 // objective function. Various parts of Ceres transform one Program into 54 // another; for example, the first stage of solving involves stripping all 55 // constant parameters and residuals. This is in contrast with Problem, which is 56 // not built for transformation. 57 class Program { 58 public: 59 Program(); 60 explicit Program(const Program& program); 61 62 // The ordered parameter and residual blocks for the program. 63 const vector<ParameterBlock*>& parameter_blocks() const; 64 const vector<ResidualBlock*>& residual_blocks() const; 65 vector<ParameterBlock*>* mutable_parameter_blocks(); 66 vector<ResidualBlock*>* mutable_residual_blocks(); 67 68 // Serialize to/from the program and update states. 69 // 70 // NOTE: Setting the state of a parameter block can trigger the 71 // computation of the Jacobian of its local parameterization. If 72 // this computation fails for some reason, then this method returns 73 // false and the state of the parameter blocks cannot be trusted. 74 bool StateVectorToParameterBlocks(const double *state); 75 void ParameterBlocksToStateVector(double *state) const; 76 77 // Copy internal state to the user's parameters. 78 void CopyParameterBlockStateToUserState(); 79 80 // Set the parameter block pointers to the user pointers. Since this 81 // runs parameter block set state internally, which may call local 82 // parameterizations, this can fail. False is returned on failure. 83 bool SetParameterBlockStatePtrsToUserStatePtrs(); 84 85 // Update a state vector for the program given a delta. 86 bool Plus(const double* state, 87 const double* delta, 88 double* state_plus_delta) const; 89 90 // Set the parameter indices and offsets. This permits mapping backward 91 // from a ParameterBlock* to an index in the parameter_blocks() vector. For 92 // any parameter block p, after calling SetParameterOffsetsAndIndex(), it 93 // is true that 94 // 95 // parameter_blocks()[p->index()] == p 96 // 97 // If a parameter appears in a residual but not in the parameter block, then 98 // it will have an index of -1. 99 // 100 // This also updates p->state_offset() and p->delta_offset(), which are the 101 // position of the parameter in the state and delta vector respectively. 102 void SetParameterOffsetsAndIndex(); 103 104 // Check if the internal state of the program (the indexing and the 105 // offsets) are correct. 106 bool IsValid() const; 107 108 bool ParameterBlocksAreFinite(string* message) const; 109 110 // Returns true if the program has any non-constant parameter blocks 111 // which have non-trivial bounds constraints. 112 bool IsBoundsConstrained() const; 113 114 // Returns false, if the program has any constant parameter blocks 115 // which are not feasible, or any variable parameter blocks which 116 // have a lower bound greater than or equal to the upper bound. 117 bool IsFeasible(string* message) const; 118 119 // Loop over each residual block and ensure that no two parameter 120 // blocks in the same residual block are part of 121 // parameter_blocks as that would violate the assumption that it 122 // is an independent set in the Hessian matrix. 123 bool IsParameterBlockSetIndependent(const set<double*>& independent_set) const; 124 125 // Create a TripletSparseMatrix which contains the zero-one 126 // structure corresponding to the block sparsity of the transpose of 127 // the Jacobian matrix. 128 // 129 // Caller owns the result. 130 TripletSparseMatrix* CreateJacobianBlockSparsityTranspose() const; 131 132 // Create a copy of this program and removes constant parameter 133 // blocks and residual blocks with no varying parameter blocks while 134 // preserving their relative order. 135 // 136 // removed_parameter_blocks on exit will contain the list of 137 // parameter blocks that were removed. 138 // 139 // fixed_cost will be equal to the sum of the costs of the residual 140 // blocks that were removed. 141 // 142 // If there was a problem, then the function will return a NULL 143 // pointer and error will contain a human readable description of 144 // the problem. 145 Program* CreateReducedProgram(vector<double*>* removed_parameter_blocks, 146 double* fixed_cost, 147 string* error) const; 148 149 // See problem.h for what these do. 150 int NumParameterBlocks() const; 151 int NumParameters() const; 152 int NumEffectiveParameters() const; 153 int NumResidualBlocks() const; 154 int NumResiduals() const; 155 156 int MaxScratchDoublesNeededForEvaluate() const; 157 int MaxDerivativesPerResidualBlock() const; 158 int MaxParametersPerResidualBlock() const; 159 int MaxResidualsPerResidualBlock() const; 160 161 // A human-readable dump of the parameter blocks for debugging. 162 // TODO(keir): If necessary, also dump the residual blocks. 163 string ToString() const; 164 165 private: 166 // Remove constant parameter blocks and residual blocks with no 167 // varying parameter blocks while preserving their relative order. 168 // 169 // removed_parameter_blocks on exit will contain the list of 170 // parameter blocks that were removed. 171 // 172 // fixed_cost will be equal to the sum of the costs of the residual 173 // blocks that were removed. 174 // 175 // If there was a problem, then the function will return false and 176 // error will contain a human readable description of the problem. 177 bool RemoveFixedBlocks(vector<double*>* removed_parameter_blocks, 178 double* fixed_cost, 179 string* message); 180 181 // The Program does not own the ParameterBlock or ResidualBlock objects. 182 vector<ParameterBlock*> parameter_blocks_; 183 vector<ResidualBlock*> residual_blocks_; 184 185 friend class ProblemImpl; 186 }; 187 188 } // namespace internal 189 } // namespace ceres 190 191 #endif // CERES_INTERNAL_PROGRAM_H_ 192