11d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Ceres Solver - A fast non-linear least squares minimizer 21d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Copyright 2013 Google Inc. All rights reserved. 31d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// http://code.google.com/p/ceres-solver/ 41d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 51d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Redistribution and use in source and binary forms, with or without 61d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// modification, are permitted provided that the following conditions are met: 71d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 81d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Redistributions of source code must retain the above copyright notice, 91d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// this list of conditions and the following disclaimer. 101d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Redistributions in binary form must reproduce the above copyright notice, 111d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// this list of conditions and the following disclaimer in the documentation 121d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// and/or other materials provided with the distribution. 131d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// * Neither the name of Google Inc. nor the names of its contributors may be 141d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// used to endorse or promote products derived from this software without 151d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// specific prior written permission. 161d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 171d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 181d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 191d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 201d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 211d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 221d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 231d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS 241d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 251d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 261d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 271d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// POSSIBILITY OF SUCH DAMAGE. 281d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// 291d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Author: sameeragarwal@google.com (Sameer Agarwal) 301d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 311d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#ifndef CERES_INTERNAL_PRECONDITIONER_H_ 321d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#define CERES_INTERNAL_PRECONDITIONER_H_ 331d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 341d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include <vector> 351d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/casts.h" 361d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/compressed_row_sparse_matrix.h" 371d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/linear_operator.h" 381d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#include "ceres/sparse_matrix.h" 3979397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez#include "ceres/types.h" 401d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 411d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingnamespace ceres { 421d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingnamespace internal { 431d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 441d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass BlockSparseMatrix; 451d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass SparseMatrix; 461d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 471d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass Preconditioner : public LinearOperator { 481d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling public: 491d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling struct Options { 501d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling Options() 511d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling : type(JACOBI), 5279397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez visibility_clustering_type(CANONICAL_VIEWS), 53399f7d09e0c45af54b77b4ab9508d6f23759b927Scott Ettinger sparse_linear_algebra_library_type(SUITE_SPARSE), 541d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling num_threads(1), 551d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling row_block_size(Eigen::Dynamic), 561d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling e_block_size(Eigen::Dynamic), 571d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling f_block_size(Eigen::Dynamic) { 581d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 591d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 601d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling PreconditionerType type; 6179397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez VisibilityClusteringType visibility_clustering_type; 62399f7d09e0c45af54b77b4ab9508d6f23759b927Scott Ettinger SparseLinearAlgebraLibraryType sparse_linear_algebra_library_type; 631d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 641d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // If possible, how many threads the preconditioner can use. 651d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int num_threads; 661d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 671d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Hints about the order in which the parameter blocks should be 681d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // eliminated by the linear solver. 691d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 701d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // For example if elimination_groups is a vector of size k, then 711d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // the linear solver is informed that it should eliminate the 721d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // parameter blocks 0 ... elimination_groups[0] - 1 first, and 731d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // then elimination_groups[0] ... elimination_groups[1] - 1 and so 741d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // on. Within each elimination group, the linear solver is free to 751d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // choose how the parameter blocks are ordered. Different linear 761d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // solvers have differing requirements on elimination_groups. 771d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 781d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // The most common use is for Schur type solvers, where there 791d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // should be at least two elimination groups and the first 801d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // elimination group must form an independent set in the normal 811d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // equations. The first elimination group corresponds to the 821d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // num_eliminate_blocks in the Schur type solvers. 831d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling vector<int> elimination_groups; 841d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 851d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // If the block sizes in a BlockSparseMatrix are fixed, then in 861d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // some cases the Schur complement based solvers can detect and 871d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // specialize on them. 881d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 891d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // It is expected that these parameters are set programmatically 901d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // rather than manually. 911d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 921d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Please see schur_complement_solver.h and schur_eliminator.h for 931d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // more details. 941d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int row_block_size; 951d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int e_block_size; 961d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling int f_block_size; 971d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling }; 981d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 9979397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez // If the optimization problem is such that there are no remaining 10079397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez // e-blocks, ITERATIVE_SCHUR with a Schur type preconditioner cannot 10179397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez // be used. This function returns JACOBI if a preconditioner for 10279397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez // ITERATIVE_SCHUR is used. The input preconditioner_type is 10379397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez // returned otherwise. 10479397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez static PreconditionerType PreconditionerForZeroEBlocks( 10579397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez PreconditionerType preconditioner_type); 10679397c21138f54fcff6ec067b44b847f1f7e0e98Carlos Hernandez 1071d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual ~Preconditioner(); 1081d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1091d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Update the numerical value of the preconditioner for the linear 1101d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // system: 1111d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 1121d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // | A | x = |b| 1131d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // |diag(D)| |0| 1141d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 1151d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // for some vector b. It is important that the matrix A have the 1161d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // same block structure as the one used to construct this object. 1171d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // 1181d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // D can be NULL, in which case its interpreted as a diagonal matrix 1191d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // of size zero. 1201d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual bool Update(const LinearOperator& A, const double* D) = 0; 1211d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1221d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // LinearOperator interface. Since the operator is symmetric, 1231d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // LeftMultiply and num_cols are just calls to RightMultiply and 1241d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // num_rows respectively. Update() must be called before 1251d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // RightMultiply can be called. 1261d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual void RightMultiply(const double* x, double* y) const = 0; 1271d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual void LeftMultiply(const double* x, double* y) const { 1281d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return RightMultiply(x, y); 1291d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 1301d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1311d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual int num_rows() const = 0; 1321d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual int num_cols() const { 1331d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return num_rows(); 1341d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 1351d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling}; 1361d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1371d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// This templated subclass of Preconditioner serves as a base class for 1381d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// other preconditioners that depend on the particular matrix layout of 1391d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// the underlying linear operator. 1401d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingtemplate <typename MatrixType> 1411d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass TypedPreconditioner : public Preconditioner { 1421d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling public: 1431d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual ~TypedPreconditioner() {} 1441d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual bool Update(const LinearOperator& A, const double* D) { 1451d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling return UpdateImpl(*down_cast<const MatrixType*>(&A), D); 1461d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling } 1471d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1481d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling private: 1491d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual bool UpdateImpl(const MatrixType& A, const double* D) = 0; 1501d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling}; 1511d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1521d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Preconditioners that depend on acccess to the low level structure 1531d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// of a SparseMatrix. 1541d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingtypedef TypedPreconditioner<SparseMatrix> SparseMatrixPreconditioner; // NOLINT 1551d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingtypedef TypedPreconditioner<BlockSparseMatrix> BlockSparseMatrixPreconditioner; // NOLINT 1561d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingtypedef TypedPreconditioner<CompressedRowSparseMatrix> CompressedRowSparseMatrixPreconditioner; // NOLINT 1571d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1581d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling// Wrap a SparseMatrix object as a preconditioner. 1591d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberlingclass SparseMatrixPreconditionerWrapper : public SparseMatrixPreconditioner { 1601d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling public: 1611d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Wrapper does NOT take ownership of the matrix pointer. 1621d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling explicit SparseMatrixPreconditionerWrapper(const SparseMatrix* matrix); 1631d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual ~SparseMatrixPreconditionerWrapper(); 1641d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1651d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling // Preconditioner interface 1661d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual void RightMultiply(const double* x, double* y) const; 1671d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual int num_rows() const; 1681d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1691d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling private: 1701d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling virtual bool UpdateImpl(const SparseMatrix& A, const double* D); 1711d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling const SparseMatrix* matrix_; 1721d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling}; 1731d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1741d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} // namespace internal 1751d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling} // namespace ceres 1761d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling 1771d2624a10e2c559f8ba9ef89eaa30832c0a83a96Sascha Haeberling#endif // CERES_INTERNAL_PRECONDITIONER_H_ 178