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: sameeragarwal@google.com (Sameer Agarwal) 30 31 #ifndef CERES_INTERNAL_COMPRESSED_ROW_SPARSE_MATRIX_H_ 32 #define CERES_INTERNAL_COMPRESSED_ROW_SPARSE_MATRIX_H_ 33 34 #include <vector> 35 #include "ceres/internal/macros.h" 36 #include "ceres/internal/port.h" 37 #include "ceres/sparse_matrix.h" 38 #include "ceres/types.h" 39 #include "glog/logging.h" 40 41 namespace ceres { 42 43 struct CRSMatrix; 44 45 namespace internal { 46 47 class TripletSparseMatrix; 48 49 class CompressedRowSparseMatrix : public SparseMatrix { 50 public: 51 // Build a matrix with the same content as the TripletSparseMatrix 52 // m. TripletSparseMatrix objects are easier to construct 53 // incrementally, so we use them to initialize SparseMatrix 54 // objects. 55 // 56 // We assume that m does not have any repeated entries. 57 explicit CompressedRowSparseMatrix(const TripletSparseMatrix& m); 58 59 // Use this constructor only if you know what you are doing. This 60 // creates a "blank" matrix with the appropriate amount of memory 61 // allocated. However, the object itself is in an inconsistent state 62 // as the rows and cols matrices do not match the values of 63 // num_rows, num_cols and max_num_nonzeros. 64 // 65 // The use case for this constructor is that when the user knows the 66 // size of the matrix to begin with and wants to update the layout 67 // manually, instead of going via the indirect route of first 68 // constructing a TripletSparseMatrix, which leads to more than 69 // double the peak memory usage. 70 CompressedRowSparseMatrix(int num_rows, 71 int num_cols, 72 int max_num_nonzeros); 73 74 // Build a square sparse diagonal matrix with num_rows rows and 75 // columns. The diagonal m(i,i) = diagonal(i); 76 CompressedRowSparseMatrix(const double* diagonal, int num_rows); 77 78 virtual ~CompressedRowSparseMatrix(); 79 80 // SparseMatrix interface. 81 virtual void SetZero(); 82 virtual void RightMultiply(const double* x, double* y) const; 83 virtual void LeftMultiply(const double* x, double* y) const; 84 virtual void SquaredColumnNorm(double* x) const; 85 virtual void ScaleColumns(const double* scale); 86 87 virtual void ToDenseMatrix(Matrix* dense_matrix) const; 88 virtual void ToTextFile(FILE* file) const; num_rows()89 virtual int num_rows() const { return num_rows_; } num_cols()90 virtual int num_cols() const { return num_cols_; } num_nonzeros()91 virtual int num_nonzeros() const { return rows_[num_rows_]; } values()92 virtual const double* values() const { return &values_[0]; } mutable_values()93 virtual double* mutable_values() { return &values_[0]; } 94 95 // Delete the bottom delta_rows. 96 // num_rows -= delta_rows 97 void DeleteRows(int delta_rows); 98 99 // Append the contents of m to the bottom of this matrix. m must 100 // have the same number of columns as this matrix. 101 void AppendRows(const CompressedRowSparseMatrix& m); 102 103 void ToCRSMatrix(CRSMatrix* matrix) const; 104 105 // Low level access methods that expose the structure of the matrix. cols()106 const int* cols() const { return &cols_[0]; } mutable_cols()107 int* mutable_cols() { return &cols_[0]; } 108 rows()109 const int* rows() const { return &rows_[0]; } mutable_rows()110 int* mutable_rows() { return &rows_[0]; } 111 row_blocks()112 const vector<int>& row_blocks() const { return row_blocks_; } mutable_row_blocks()113 vector<int>* mutable_row_blocks() { return &row_blocks_; } 114 col_blocks()115 const vector<int>& col_blocks() const { return col_blocks_; } mutable_col_blocks()116 vector<int>* mutable_col_blocks() { return &col_blocks_; } 117 118 // Destructive array resizing method. 119 void SetMaxNumNonZeros(int num_nonzeros); 120 121 // Non-destructive array resizing method. set_num_rows(const int num_rows)122 void set_num_rows(const int num_rows) { num_rows_ = num_rows; } set_num_cols(const int num_cols)123 void set_num_cols(const int num_cols) { num_cols_ = num_cols; } 124 125 void SolveLowerTriangularInPlace(double* solution) const; 126 void SolveLowerTriangularTransposeInPlace(double* solution) const; 127 128 CompressedRowSparseMatrix* Transpose() const; 129 130 static CompressedRowSparseMatrix* CreateBlockDiagonalMatrix( 131 const double* diagonal, 132 const vector<int>& blocks); 133 134 // Compute the sparsity structure of the product m.transpose() * m 135 // and create a CompressedRowSparseMatrix corresponding to it. 136 // 137 // Also compute a "program" vector, which for every term in the 138 // outer product points to the entry in the values array of the 139 // result matrix where it should be accumulated. 140 // 141 // This program is used by the ComputeOuterProduct function below to 142 // compute the outer product. 143 // 144 // Since the entries of the program are the same for rows with the 145 // same sparsity structure, the program only stores the result for 146 // one row per row block. The ComputeOuterProduct function reuses 147 // this information for each row in the row block. 148 static CompressedRowSparseMatrix* CreateOuterProductMatrixAndProgram( 149 const CompressedRowSparseMatrix& m, 150 vector<int>* program); 151 152 // Compute the values array for the expression m.transpose() * m, 153 // where the matrix used to store the result and a program have been 154 // created using the CreateOuterProductMatrixAndProgram function 155 // above. 156 static void ComputeOuterProduct(const CompressedRowSparseMatrix& m, 157 const vector<int>& program, 158 CompressedRowSparseMatrix* result); 159 160 private: 161 int num_rows_; 162 int num_cols_; 163 vector<int> rows_; 164 vector<int> cols_; 165 vector<double> values_; 166 167 // If the matrix has an underlying block structure, then it can also 168 // carry with it row and column block sizes. This is auxilliary and 169 // optional information for use by algorithms operating on the 170 // matrix. The class itself does not make use of this information in 171 // any way. 172 vector<int> row_blocks_; 173 vector<int> col_blocks_; 174 175 CERES_DISALLOW_COPY_AND_ASSIGN(CompressedRowSparseMatrix); 176 }; 177 178 } // namespace internal 179 } // namespace ceres 180 181 #endif // CERES_INTERNAL_COMPRESSED_ROW_SPARSE_MATRIX_H_ 182