21 #ifndef OPM_DUNEMATRIX_HEADER_INCLUDED
22 #define OPM_DUNEMATRIX_HEADER_INCLUDED
24 #ifdef DUNE_BCRSMATRIX_HH
25 #error This header must be included before any bcrsmatrix.hh is included (directly or indirectly)
28 #include <opm/common/utility/platform_dependent/disable_warnings.h>
30 #include <Eigen/Eigen>
31 #include <Eigen/Sparse>
33 #include <dune/common/fmatrix.hh>
34 #include <dune/common/version.hh>
36 #if DUNE_VERSION_NEWER(DUNE_ISTL,2,4)
37 #include <dune/istl/bcrsmatrix.hh>
40 #define private protected
41 #include <dune/istl/bcrsmatrix.hh>
45 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
50 class DuneMatrix :
public Dune::BCRSMatrix< Dune::FieldMatrix<double, 1, 1> >
53 DuneMatrix(
const int rows,
const int cols,
const int* ia,
const int* ja,
const double* sa)
56 init( rows, cols, ia, ja, sa );
60 DuneMatrix(
const Eigen::SparseMatrix<double, Eigen::RowMajor>& matrix )
63 const int rows = matrix.rows();
64 const int cols = matrix.cols();
65 const int* ia = matrix.outerIndexPtr();
66 const int* ja = matrix.innerIndexPtr();
67 const double* sa = matrix.valuePtr();
69 init( rows, cols, ia, ja, sa );
73 void init(
const int rows,
const int cols,
const int* ia,
const int* ja,
const double* sa)
75 typedef Dune::BCRSMatrix< Dune::FieldMatrix<double, 1, 1> > Super;
76 typedef Super::block_type block_type;
77 this->build_mode = Super::unknown;
78 this->ready = Super::built;
83 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 3)
84 this->allocationSize = this->nnz;
86 this->overflowsize = -1.0;
91 this->a = this->allocator_.allocate(this->nnz);
92 static_assert(
sizeof(block_type) ==
sizeof(
double),
"This constructor requires a block type that is the same as a double.");
93 std::copy(sa, sa + this->nnz, reinterpret_cast<double*>(this->a));
94 this->j.reset(this->sizeAllocator_.allocate(this->nnz));
95 std::copy(ja, ja +this-> nnz, this->j.get());
96 this->r = rowAllocator_.allocate(rows);
97 for (
int row = 0; row < rows; ++row) {
98 this->r[row].set(ia[row+1] - ia[row], this->a + ia[row], this->j.get() + ia[row]);
105 #endif // OPM_DUNEMATRIX_HEADER_INCLUDED
Definition: AdditionalObjectDeleter.hpp:22
void init(const int rows, const int cols, const int *ia, const int *ja, const double *sa)
Definition: DuneMatrix.hpp:73
DuneMatrix(const int rows, const int cols, const int *ia, const int *ja, const double *sa)
Definition: DuneMatrix.hpp:53
Definition: DuneMatrix.hpp:50
DuneMatrix(const Eigen::SparseMatrix< double, Eigen::RowMajor > &matrix)
create an ISTL BCRSMatrix from a Eigen::SparseMatrix
Definition: DuneMatrix.hpp:60