DOLFIN
DOLFIN C++ interface
EigenMatrix.h
1 // Copyright (C) 2015 Chris Richardson and Garth N. Wells
2 //
3 // This file is part of DOLFIN.
4 //
5 // DOLFIN is free software: you can redistribute it and/or modify
6 // it under the terms of the GNU Lesser General Public License as published by
7 // the Free Software Foundation, either version 3 of the License, or
8 // (at your option) any later version.
9 //
10 // DOLFIN is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU Lesser General Public License for more details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17 
18 #ifndef __DOLFIN_EIGEN_MATRIX_H
19 #define __DOLFIN_EIGEN_MATRIX_H
20 
21 #include <iomanip>
22 #include <memory>
23 #include <string>
24 #include <utility>
25 #include <vector>
26 #include <Eigen/Sparse>
27 
28 #include <dolfin/common/MPI.h>
29 #include <dolfin/common/types.h>
30 #include "EigenVector.h"
31 #include "GenericMatrix.h"
32 #include "GenericVector.h"
33 #include "TensorLayout.h"
34 
35 namespace dolfin
36 {
37 
45 
46  class EigenMatrix : public GenericMatrix
47  {
48  public:
49 
51  typedef Eigen::SparseMatrix<double, Eigen::RowMajor, int> eigen_matrix_type;
52 
54  EigenMatrix();
55 
57  EigenMatrix(std::size_t M, std::size_t N);
58 
60  EigenMatrix(const EigenMatrix& A);
61 
63  virtual ~EigenMatrix();
64 
65  //--- Implementation of the GenericTensor interface ---
66 
68  virtual void init(const TensorLayout& tensor_layout);
69 
71  virtual bool empty() const
72  { return size(0) == 0; }
73 
75  virtual std::size_t size(std::size_t dim) const;
76 
78  virtual std::pair<std::int64_t, std::int64_t>
79  local_range(std::size_t dim) const
80  { return {0, size(dim)}; }
81 
83  std::size_t nnz() const;
84 
86  virtual void zero();
87 
89  virtual void apply(std::string mode);
90 
92  virtual MPI_Comm mpi_comm() const
93  { return _mpi_comm.comm(); }
94 
96  virtual std::string str(bool verbose) const;
97 
98  //--- Implementation of the GenericMatrix interface ---
99 
101  virtual std::shared_ptr<GenericMatrix> copy() const;
102 
104  virtual void resize(std::size_t M, std::size_t N);
105 
112  virtual void init_vector(GenericVector& z, std::size_t dim) const;
113 
115  virtual void get(double* block, std::size_t m, const dolfin::la_index* rows,
116  std::size_t n, const dolfin::la_index* cols) const;
117 
119  virtual void set(const double* block, std::size_t m,
120  const dolfin::la_index* rows, std::size_t n,
121  const dolfin::la_index* cols);
122 
124  virtual void set_local(const double* block, std::size_t m,
125  const dolfin::la_index* rows, std::size_t n,
126  const dolfin::la_index* cols)
127  { set(block, m, rows, n, cols); }
128 
130  virtual void add(const double* block, std::size_t m,
131  const dolfin::la_index* rows, std::size_t n,
132  const dolfin::la_index* cols);
133 
135  virtual void add_local(const double* block, std::size_t m,
136  const dolfin::la_index* rows, std::size_t n,
137  const dolfin::la_index* cols)
138  { add(block, m, rows, n, cols); }
139 
141  virtual void axpy(double a, const GenericMatrix& A,
142  bool same_nonzero_pattern);
143 
145  virtual double norm(std::string norm_type) const;
146 
148  virtual void getrow(std::size_t row, std::vector<std::size_t>& columns,
149  std::vector<double>& values) const;
150 
152  virtual void setrow(std::size_t row_idx,
153  const std::vector<std::size_t>& columns,
154  const std::vector<double>& values);
155 
157  virtual void zero(std::size_t m, const dolfin::la_index* rows);
158 
160  virtual void zero_local(std::size_t m, const dolfin::la_index* rows)
161  { zero(m, rows); }
162 
164  virtual void ident(std::size_t m, const dolfin::la_index* rows);
165 
167  virtual void ident_local(std::size_t m, const dolfin::la_index* rows)
168  { ident(m, rows); }
169 
171  virtual void mult(const GenericVector& x, GenericVector& y) const;
172 
174  virtual void transpmult(const GenericVector& x, GenericVector& y) const;
175 
177  virtual void get_diagonal(GenericVector& x) const;
178 
180  virtual void set_diagonal(const GenericVector& x);
181 
183  virtual const EigenMatrix& operator*= (double a);
184 
186  virtual const EigenMatrix& operator/= (double a);
187 
189  virtual const GenericMatrix& operator= (const GenericMatrix& A);
190 
193  virtual std::tuple<const int*, const int*, const double*, std::size_t>
194  data() const;
195 
196  //--- Special functions ---
197 
199  virtual GenericLinearAlgebraFactory& factory() const;
200 
201  //--- Special Eigen functions ---
202 
204  const eigen_matrix_type& mat() const
205  { return _matA; }
206 
208  eigen_matrix_type& mat()
209  { return _matA; }
210 
212  void compress()
213  { _matA.makeCompressed(); }
214 
217  { return _matA.coeff(i, j); }
218 
220  const EigenMatrix& operator= (const EigenMatrix& A);
221 
222  private:
223 
224  // MPI communicator
225  dolfin::MPI::Comm _mpi_comm;
226 
227  // Eigen matrix object - row major access
228  eigen_matrix_type _matA;
229 
230  };
231 }
232 
233 #endif
virtual MPI_Comm mpi_comm() const
Return MPI communicator.
Definition: EigenMatrix.h:92
virtual void zero_local(std::size_t m, const dolfin::la_index *rows)
Set given rows (local row indices) to zero.
Definition: EigenMatrix.h:160
void compress()
Compress matrix (eliminate all zeros from a sparse matrix)
Definition: EigenMatrix.h:212
double operator()(dolfin::la_index i, dolfin::la_index j) const
Access value of given entry.
Definition: EigenMatrix.h:216
virtual bool empty() const
Return true if empty.
Definition: EigenMatrix.h:71
Definition: adapt.h:29
virtual void set_local(const double *block, std::size_t m, const dolfin::la_index *rows, std::size_t n, const dolfin::la_index *cols)
Set block of values using local indices.
Definition: EigenMatrix.h:124
virtual void ident_local(std::size_t m, const dolfin::la_index *rows)
Set given rows to identity matrix.
Definition: EigenMatrix.h:167
virtual std::shared_ptr< GenericMatrix > copy() const
Return copy of matrix.
Definition: EigenMatrix.cpp:52
virtual std::string str(bool verbose) const
Return informal string representation (pretty-print)
Definition: EigenMatrix.cpp:357
const eigen_matrix_type & mat() const
Return reference to Eigen matrix (const version)
Definition: EigenMatrix.h:204
virtual void getrow(std::size_t row, std::vector< std::size_t > &columns, std::vector< double > &values) const
Get non-zero values of given row.
Definition: EigenMatrix.cpp:107
eigen_matrix_type & mat()
Return reference to Eigen matrix (non-const version)
Definition: EigenMatrix.h:208
Definition: EigenMatrix.h:46
virtual ~EigenMatrix()
Destructor.
Definition: EigenMatrix.cpp:47
virtual void transpmult(const GenericVector &x, GenericVector &y) const
Matrix-vector product, y = A^T x.
Definition: EigenMatrix.cpp:285
virtual void zero()
Set all entries to zero and keep any sparse structure.
Definition: EigenMatrix.cpp:180
Definition: TensorLayout.h:41
virtual void get_diagonal(GenericVector &x) const
Get diagonal of a matrix.
Definition: EigenMatrix.cpp:257
Base class for LinearAlgebra factories.
Definition: GenericLinearAlgebraFactory.h:46
virtual const EigenMatrix & operator*=(double a)
Multiply matrix by given number.
Definition: EigenMatrix.cpp:314
virtual void init(const TensorLayout &tensor_layout)
Initialize zero tensor using tenor layout.
Definition: EigenMatrix.cpp:384
virtual const GenericMatrix & operator=(const GenericMatrix &A)
Assignment operator.
Definition: EigenMatrix.cpp:326
Eigen::SparseMatrix< double, Eigen::RowMajor, int > eigen_matrix_type
Eigen Matrix type.
Definition: EigenMatrix.h:51
MPI_Comm comm() const
Return the underlying MPI_Comm object.
Definition: MPI.cpp:117
virtual const EigenMatrix & operator/=(double a)
Divide matrix by given number.
Definition: EigenMatrix.cpp:320
virtual void init_vector(GenericVector &z, std::size_t dim) const
Definition: EigenMatrix.cpp:139
virtual void ident(std::size_t m, const dolfin::la_index *rows)
Set given rows to identity matrix.
Definition: EigenMatrix.cpp:195
virtual void apply(std::string mode)
Finalize assembly of tensor.
Definition: EigenMatrix.cpp:417
virtual void axpy(double a, const GenericMatrix &A, bool same_nonzero_pattern)
Add multiple of given matrix (AXPY operation)
Definition: EigenMatrix.cpp:422
virtual std::tuple< const int *, const int *, const double *, std::size_t > data() const
Definition: EigenMatrix.cpp:342
virtual void mult(const GenericVector &x, GenericVector &y) const
Matrix-vector product, y = Ax.
Definition: EigenMatrix.cpp:230
virtual void resize(std::size_t M, std::size_t N)
Resize matrix to M x N.
Definition: EigenMatrix.cpp:57
This class defines a common interface for matrices.
Definition: GenericMatrix.h:46
virtual std::pair< std::int64_t, std::int64_t > local_range(std::size_t dim) const
Return local ownership range.
Definition: EigenMatrix.h:79
PetscInt la_index
Index type for compatibility with linear algebra backend(s)
Definition: types.h:32
virtual std::size_t size(std::size_t dim) const
Return size of given dimension.
Definition: EigenMatrix.cpp:65
virtual GenericLinearAlgebraFactory & factory() const
Return linear algebra backend factory.
Definition: EigenMatrix.cpp:25
virtual double norm(std::string norm_type) const
Return norm of matrix.
Definition: EigenMatrix.cpp:77
virtual void add_local(const double *block, std::size_t m, const dolfin::la_index *rows, std::size_t n, const dolfin::la_index *cols)
Add block of values using local indices.
Definition: EigenMatrix.h:135
This class defines a common interface for vectors.
Definition: GenericVector.h:47
EigenMatrix()
Create empty matrix.
Definition: EigenMatrix.cpp:30
Definition: MPI.h:76
virtual void add(const double *block, std::size_t m, const dolfin::la_index *rows, std::size_t n, const dolfin::la_index *cols)
Add block of values using global indices.
Definition: EigenMatrix.cpp:156
virtual void set_diagonal(const GenericVector &x)
Set diagonal of a matrix.
Definition: EigenMatrix.cpp:271
virtual void setrow(std::size_t row_idx, const std::vector< std::size_t > &columns, const std::vector< double > &values)
Set values for given row.
Definition: EigenMatrix.cpp:129
std::size_t nnz() const
Return number of non-zero entries in matrix.
Definition: EigenMatrix.cpp:412