internal/ceres/compressed_row_sparse_matrix.h

   1 // Ceres Solver - A fast non-linear least squares minimizer
   2 // Copyright 2015 Google Inc. All rights reserved.
   3 // http://ceres-solver.org/
   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   enum StorageType {
  52     UNSYMMETRIC,
  53     // Matrix is assumed to be symmetric but only the lower triangular
  54     // part of the matrix is stored.
  55     LOWER_TRIANGULAR,
  56     // Matrix is assumed to be symmetric but only the upper triangular
  57     // part of the matrix is stored.
  58     UPPER_TRIANGULAR
  59   };
  60
  61   // Create a matrix with the same content as the TripletSparseMatrix
  62   // input. We assume that input does not have any repeated
  63   // entries.
  64   //
  65   // The storage type of the matrix is set to UNSYMMETRIC.
  66   //
  67   // Caller owns the result.
  68   static CompressedRowSparseMatrix* FromTripletSparseMatrix(
  69       const TripletSparseMatrix& input);
  70
  71   // Create a matrix with the same content as the TripletSparseMatrix
  72   // input transposed. We assume that input does not have any repeated
  73   // entries.
  74   //
  75   // The storage type of the matrix is set to UNSYMMETRIC.
  76   //
  77   // Caller owns the result.
  78   static CompressedRowSparseMatrix* FromTripletSparseMatrixTransposed(
  79       const TripletSparseMatrix& input);
  80
  81   // Use this constructor only if you know what you are doing. This
  82   // creates a "blank" matrix with the appropriate amount of memory
  83   // allocated. However, the object itself is in an inconsistent state
  84   // as the rows and cols matrices do not match the values of
  85   // num_rows, num_cols and max_num_nonzeros.
  86   //
  87   // The use case for this constructor is that when the user knows the
  88   // size of the matrix to begin with and wants to update the layout
  89   // manually, instead of going via the indirect route of first
  90   // constructing a TripletSparseMatrix, which leads to more than
  91   // double the peak memory usage.
  92   //
  93   // The storage type is set to UNSYMMETRIC.
  94   CompressedRowSparseMatrix(int num_rows,
  95                             int num_cols,
  96                             int max_num_nonzeros);
  97
  98   // Build a square sparse diagonal matrix with num_rows rows and
  99   // columns. The diagonal m(i,i) = diagonal(i);
 100   //
 101   // The storage type is set to UNSYMMETRIC
 102   CompressedRowSparseMatrix(const double* diagonal, int num_rows);
 103
 104   // SparseMatrix interface.
 105   virtual ~CompressedRowSparseMatrix();
 106   virtual void SetZero();
 107   virtual void RightMultiply(const double* x, double* y) const;
 108   virtual void LeftMultiply(const double* x, double* y) const;
 109   virtual void SquaredColumnNorm(double* x) const;
 110   virtual void ScaleColumns(const double* scale);
 111
 112   virtual void ToDenseMatrix(Matrix* dense_matrix) const;
 113   virtual void ToTextFile(FILE* file) const;
 114   virtual int num_rows() const { return num_rows_; }
 115   virtual int num_cols() const { return num_cols_; }
 116   virtual int num_nonzeros() const { return rows_[num_rows_]; }
 117   virtual const double* values() const { return &values_[0]; }
 118   virtual double* mutable_values() { return &values_[0]; }
 119
 120   // Delete the bottom delta_rows.
 121   // num_rows -= delta_rows
 122   void DeleteRows(int delta_rows);
 123
 124   // Append the contents of m to the bottom of this matrix. m must
 125   // have the same number of columns as this matrix.
 126   void AppendRows(const CompressedRowSparseMatrix& m);
 127
 128   void ToCRSMatrix(CRSMatrix* matrix) const;
 129
 130   CompressedRowSparseMatrix* Transpose() const;
 131
 132   // Destructive array resizing method.
 133   void SetMaxNumNonZeros(int num_nonzeros);
 134
 135   // Non-destructive array resizing method.
 136   void set_num_rows(const int num_rows) { num_rows_ = num_rows; }
 137   void set_num_cols(const int num_cols) { num_cols_ = num_cols; }
 138
 139   // Low level access methods that expose the structure of the matrix.
 140   const int* cols() const { return &cols_[0]; }
 141   int* mutable_cols() { return &cols_[0]; }
 142
 143   const int* rows() const { return &rows_[0]; }
 144   int* mutable_rows() { return &rows_[0]; }
 145
 146   const StorageType storage_type() const { return storage_type_; }
 147   void set_storage_type(const StorageType storage_type) {
 148     storage_type_ = storage_type;
 149   }
 150
 151   const std::vector<int>& row_blocks() const { return row_blocks_; }
 152   std::vector<int>* mutable_row_blocks() { return &row_blocks_; }
 153
 154   const std::vector<int>& col_blocks() const { return col_blocks_; }
 155   std::vector<int>* mutable_col_blocks() { return &col_blocks_; }
 156
 157   // Create a block diagonal CompressedRowSparseMatrix with the given
 158   // block structure. The individual blocks are assumed to be laid out
 159   // contiguously in the diagonal array, one block at a time.
 160   //
 161   // Caller owns the result.
 162   static CompressedRowSparseMatrix* CreateBlockDiagonalMatrix(
 163       const double* diagonal,
 164       const std::vector<int>& blocks);
 165
 166   // Options struct to control the generation of random block sparse
 167   // matrices in compressed row sparse format.
 168   //
 169   // The random matrix generation proceeds as follows.
 170   //
 171   // First the row and column block structure is determined by
 172   // generating random row and column block sizes that lie within the
 173   // given bounds.
 174   //
 175   // Then we walk the block structure of the resulting matrix, and with
 176   // probability block_density detemine whether they are structurally
 177   // zero or not. If the answer is no, then we generate entries for the
 178   // block which are distributed normally.
 179   struct RandomMatrixOptions {
 180     RandomMatrixOptions()
 181         : num_row_blocks(0),
 182           min_row_block_size(0),
 183           max_row_block_size(0),
 184           num_col_blocks(0),
 185           min_col_block_size(0),
 186           max_col_block_size(0),
 187           block_density(0.0) {
 188     }
 189
 190     int num_row_blocks;
 191     int min_row_block_size;
 192     int max_row_block_size;
 193     int num_col_blocks;
 194     int min_col_block_size;
 195     int max_col_block_size;
 196
 197     // 0 < block_density <= 1 is the probability of a block being
 198     // present in the matrix. A given random matrix will not have
 199     // precisely this density.
 200     double block_density;
 201   };
 202
 203   // Create a random CompressedRowSparseMatrix whose entries are
 204   // normally distributed and whose structure is determined by
 205   // RandomMatrixOptions.
 206   //
 207   // Caller owns the result.
 208   static CompressedRowSparseMatrix* CreateRandomMatrix(
 209       const RandomMatrixOptions& options);
 210
 211  private:
 212   static CompressedRowSparseMatrix* FromTripletSparseMatrix(
 213       const TripletSparseMatrix& input, bool transpose);
 214
 215   int num_rows_;
 216   int num_cols_;
 217   std::vector<int> rows_;
 218   std::vector<int> cols_;
 219   std::vector<double> values_;
 220   StorageType storage_type_;
 221
 222   // If the matrix has an underlying block structure, then it can also
 223   // carry with it row and column block sizes. This is auxilliary and
 224   // optional information for use by algorithms operating on the
 225   // matrix. The class itself does not make use of this information in
 226   // any way.
 227   std::vector<int> row_blocks_;
 228   std::vector<int> col_blocks_;
 229 };
 230
 231 }  // namespace internal
 232 }  // namespace ceres
 233
 234 #endif  // CERES_INTERNAL_COMPRESSED_ROW_SPARSE_MATRIX_H_