3 #ifndef DUNE_ISTL_SUPERMATRIX_HH 4 #define DUNE_ISTL_SUPERMATRIX_HH 10 #include <dune/common/fmatrix.hh> 11 #include <dune/common/fvector.hh> 12 #include <dune/common/typetraits.hh> 34 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
35 float *values,
int *rowindex,
int* colindex,
36 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
38 sCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
46 static void print(
char* name, SuperMatrix*
mat)
48 sPrint_CompCol_Matrix(name, mat);
57 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
58 double *values,
int *rowindex,
int* colindex,
59 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
61 dCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
69 static void print(
char* name, SuperMatrix*
mat)
71 dPrint_CompCol_Matrix(name, mat);
80 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
81 std::complex<float> *values,
int *rowindex,
int* colindex,
82 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
84 cCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast< ::complex*>(values),
85 rowindex, colindex, stype, dtype, mtype);
92 static void print(
char* name, SuperMatrix*
mat)
94 cPrint_CompCol_Matrix(name, mat);
103 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
104 std::complex<double> *values,
int *rowindex,
int* colindex,
105 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
107 zCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast<doublecomplex*>(values),
108 rowindex, colindex, stype, dtype, mtype);
115 static void print(
char* name, SuperMatrix*
mat)
117 zPrint_CompCol_Matrix(name, mat);
134 std::is_same<T,float>::value ? SLU_S :
135 ( std::is_same<T,std::complex<double> >::value ? SLU_Z :
136 ( std::is_same<T,std::complex<float> >::value ? SLU_C : SLU_D ));
185 template<
class B,
class TA,
int n,
int m>
189 template<
class M,
class X,
class TM,
class TD,
class T1>
213 if (this->N_+this->M_*this->Nnz_ != 0)
218 operator SuperMatrix&()
224 operator const SuperMatrix&()
const 233 ::create(&A, this->N_, this->M_, this->colstart[this->N_],
234 this->values,this->rowindex, this->colstart, SLU_NC,
243 ::create(&A, this->N_, this->M_, this->colstart[this->N_],
244 this->values,this->rowindex, this->colstart, SLU_NC,
255 virtual void setMatrix(
const Matrix&
mat,
const std::set<std::size_t>& mrs)
257 if(this->N_+this->M_+this->Nnz_!=0)
259 this->N_=mrs.size()*n;
260 this->M_=mrs.size()*m;
280 SUPERLU_FREE(A.Store);
286 template<
class T,
class A,
int n,
int m>
290 template<
class I,
class S,
class D>
307 ::create(&slumat->A, slumat->N_, slumat->M_, slumat->colstart[this->cols],
308 slumat->values,slumat->rowindex, slumat->colstart, SLU_NC,
312 SuperLUMatrix* slumat;
315 #endif // HAVE_SUPERLU Inititializer for the ColCompMatrix as needed by OverlappingSchwarz.
Definition: colcompmatrix.hh:153
Definition: supermatrix.hh:23
Definition: allocator.hh:7
double float_type
Definition: supermatrix.hh:156
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:422
Definition: supermatrix.hh:123
double float_type
Definition: supermatrix.hh:142
This file implements a vector space as a tensor product of a given vector space. The number of compon...
float float_type
Definition: supermatrix.hh:149
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:42
size_type N() const
number of rows (counted in blocks)
Definition: bcrsmatrix.hh:1894
size_type M() const
number of columns (counted in blocks)
Definition: bcrsmatrix.hh:1900
Definition: colcompmatrix.hh:160
virtual void setMatrix(const Matrix &mat)
Initialize data from given matrix.
Definition: supermatrix.hh:267
Utility class for converting an ISTL Matrix into a column-compressed matrix.
Definition: colcompmatrix.hh:144
static const Dtype_t type
Definition: supermatrix.hh:125
Matrix & mat
Definition: matrixmatrix.hh:345
Definition: matrixutils.hh:25
Implementation of the BCRSMatrix class.
SuperLUMatrix()
Definition: supermatrix.hh:207
SuperMatrixInitializer(SuperLUMatrix &lum)
Definition: supermatrix.hh:296
Sequential overlapping Schwarz preconditioner.
Definition: colcompmatrix.hh:157
Definition: superlu.hh:35
Definition: supermatrix.hh:27
float float_type
Definition: supermatrix.hh:163
SuperLUMatrix(const Matrix &mat)
Constructor that initializes the data.
Definition: supermatrix.hh:204
Definition: supermatrix.hh:176
Definition: supermatrix.hh:129
virtual void setMatrix(const Matrix &mat, const std::set< std::size_t > &mrs)
Initialize data from a given set of matrix rows and columns.
Definition: supermatrix.hh:255
void copyToColCompMatrix(F &initializer, const MRS &mrs)
Definition: colcompmatrix.hh:428
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:457
Matrix::size_type size_type
Definition: supermatrix.hh:198
BCRSMatrix< FieldMatrix< B, n, m >, TA > Matrix
The type of the matrix to convert.
Definition: supermatrix.hh:194
virtual void free()
free allocated space.
Definition: supermatrix.hh:277
Utility class for converting an ISTL Matrix into a SuperLU Matrix.
Definition: supermatrix.hh:172
Provides access to an iterator over an arbitrary subset of matrix rows.
Definition: colcompmatrix.hh:59
virtual ~SuperLUMatrix()
Destructor.
Definition: supermatrix.hh:211
virtual void createMatrix() const
Definition: supermatrix.hh:303
SuperMatrixInitializer()
Definition: supermatrix.hh:300
Dune::SuperLUMatrix< Matrix > SuperLUMatrix
Definition: supermatrix.hh:294
Provides access to an iterator over all matrix rows.
Definition: colcompmatrix.hh:21
SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > & operator=(const SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > &mat)
Definition: supermatrix.hh:239
SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > & operator=(const BCRSMatrix< FieldMatrix< B, n, m >, TA > &mat)
Definition: supermatrix.hh:229