Home | History | Annotate | Download | only in ceres
      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 (at) 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;
     89   virtual int num_rows() const { return num_rows_; }
     90   virtual int num_cols() const { return num_cols_; }
     91   virtual int num_nonzeros() const { return rows_[num_rows_]; }
     92   virtual const double* values() const { return &values_[0]; }
     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.
    106   const int* cols() const { return &cols_[0]; }
    107   int* mutable_cols() { return &cols_[0]; }
    108 
    109   const int* rows() const { return &rows_[0]; }
    110   int* mutable_rows() { return &rows_[0]; }
    111 
    112   const vector<int>& row_blocks() const { return row_blocks_; }
    113   vector<int>* mutable_row_blocks() { return &row_blocks_; }
    114 
    115   const vector<int>& col_blocks() const { return 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.
    122   void set_num_rows(const int num_rows) { num_rows_ = num_rows; }
    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