3 #ifndef DUNE_ISTL_SUPERMATRIX_HH
4 #define DUNE_ISTL_SUPERMATRIX_HH
7 #ifdef SUPERLU_POST_2005_VERSION
10 #define SUPERLU_NTYPE 1
14 #include "slu_sdefs.h"
18 #include "slu_ddefs.h"
22 #include "slu_cdefs.h"
26 #include "slu_zdefs.h"
50 #include <dune/common/fmatrix.hh>
51 #include <dune/common/fvector.hh>
52 #include <dune/common/typetraits.hh>
72 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
73 float *values,
int *rowindex,
int* colindex,
74 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
76 sCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
84 static void print(
char* name, SuperMatrix*
mat)
86 sPrint_CompCol_Matrix(name, mat);
93 struct SuperMatrixCreateSparseChooser<double>
95 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
96 double *values,
int *rowindex,
int* colindex,
97 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
99 dCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
100 stype, dtype, mtype);
105 struct SuperMatrixPrinter<double>
107 static void print(
char* name, SuperMatrix*
mat)
109 dPrint_CompCol_Matrix(name, mat);
116 struct SuperMatrixCreateSparseChooser<std::complex<float> >
118 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
119 std::complex<float> *values,
int *rowindex,
int* colindex,
120 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
122 cCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast< ::complex*>(values),
123 rowindex, colindex, stype, dtype, mtype);
128 struct SuperMatrixPrinter<std::complex<float> >
130 static void print(
char* name, SuperMatrix*
mat)
132 cPrint_CompCol_Matrix(name, mat);
139 struct SuperMatrixCreateSparseChooser<std::complex<double> >
141 static void create(SuperMatrix *
mat,
int n,
int m,
int offset,
142 std::complex<double> *values,
int *rowindex,
int* colindex,
143 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
145 zCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast<doublecomplex*>(values),
146 rowindex, colindex, stype, dtype, mtype);
151 struct SuperMatrixPrinter<std::complex<double> >
153 static void print(
char* name, SuperMatrix*
mat)
155 zPrint_CompCol_Matrix(name, mat);
172 Dune::is_same<T,float>::value ? SLU_S :
173 ( Dune::is_same<T,std::complex<double> >::value ? SLU_Z :
174 ( Dune::is_same<T,std::complex<float> >::value ? SLU_C : SLU_D ));
223 template<
class B,
class TA,
int n,
int m>
227 template<
class M,
class X,
class TM,
class TD,
class T1>
251 if (this->N_+this->M_*this->Nnz_ != 0)
256 operator SuperMatrix&()
262 operator const SuperMatrix&()
const
271 ::create(&
A, this->N_, this->M_, this->colstart[this->N_],
272 this->values,this->rowindex, this->colstart, SLU_NC,
281 ::create(&
A, this->N_, this->M_, this->colstart[this->N_],
282 this->values,this->rowindex, this->colstart, SLU_NC,
295 if(this->N_+this->M_+this->Nnz_!=0)
297 this->N_=mrs.size()*n;
298 this->M_=mrs.size()*m;
318 SUPERLU_FREE(
A.Store);
324 template<
class T,
class A,
int n,
int m>
328 template<
class I,
class S,
class D>
345 ::create(&slumat->A, slumat->N_, slumat->M_, slumat->colstart[this->cols],
346 slumat->values,slumat->rowindex, slumat->colstart, SLU_NC,
Definition: supermatrix.hh:214
Definition: supermatrix.hh:61
virtual void free()
free allocated space.
Definition: supermatrix.hh:315
BCRSMatrix< FieldMatrix< B, n, m >, TA > Matrix
The type of the matrix to convert.
Definition: supermatrix.hh:232
size_type M() const
number of columns (counted in blocks)
Definition: bcrsmatrix.hh:1855
Definition: supermatrix.hh:161
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:447
SuperLUMatrix()
Definition: supermatrix.hh:245
Inititializer for the ColCompMatrix as needed by OverlappingSchwarz.
Definition: colcompmatrix.hh:153
SuperMatrixInitializer(SuperLUMatrix &lum)
Definition: supermatrix.hh:334
static void print(char *name, SuperMatrix *mat)
Definition: supermatrix.hh:84
Sequential overlapping Schwarz preconditioner.
Definition: colcompmatrix.hh:157
Matrix & A
Definition: matrixmatrix.hh:216
void copyToColCompMatrix(F &initializer, const MRS &mrs)
Definition: colcompmatrix.hh:431
SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > & operator=(const BCRSMatrix< FieldMatrix< B, n, m >, TA > &mat)
Definition: supermatrix.hh:267
virtual void createMatrix() const
Definition: supermatrix.hh:341
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:412
Definition: colcompmatrix.hh:160
static void create(SuperMatrix *mat, int n, int m, int offset, float *values, int *rowindex, int *colindex, Stype_t stype, Dtype_t dtype, Mtype_t mtype)
Definition: supermatrix.hh:72
float float_type
Definition: supermatrix.hh:201
size_type N() const
number of rows (counted in blocks)
Definition: bcrsmatrix.hh:1849
Provides access to an iterator over all matrix rows.
Definition: colcompmatrix.hh:21
Definition: supermatrix.hh:65
static const Dtype_t type
Definition: supermatrix.hh:163
virtual void setMatrix(const Matrix &mat)
Initialize data from given matrix.
Definition: supermatrix.hh:305
Matrix & mat
Definition: matrixmatrix.hh:343
virtual ~SuperLUMatrix()
Destructor.
Definition: supermatrix.hh:249
SuperMatrixInitializer()
Definition: supermatrix.hh:338
Definition: supermatrix.hh:167
double float_type
Definition: supermatrix.hh:180
Utility class for converting an ISTL Matrix into a SuperLU Matrix.
Definition: supermatrix.hh:210
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:42
Definition: matrixutils.hh:25
SuperLUMatrix(const Matrix &mat)
Constructor that initializes the data.
Definition: supermatrix.hh:242
Definition: superlu.hh:80
Implementation of the BCRSMatrix class.
Utility class for converting an ISTL Matrix into a column-compressed matrix.
Definition: colcompmatrix.hh:144
BCRSMatrix< FieldMatrix< T, n, m >, A > Matrix
Definition: supermatrix.hh:331
double float_type
Definition: supermatrix.hh:194
This file implements a vector space as a tensor product of a given vector space. The number of compon...
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:293
Provides access to an iterator over an arbitrary subset of matrix rows.
Definition: colcompmatrix.hh:59
SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > & operator=(const SuperLUMatrix< BCRSMatrix< FieldMatrix< B, n, m >, TA > > &mat)
Definition: supermatrix.hh:277
float float_type
Definition: supermatrix.hh:187
Dune::SuperLUMatrix< Matrix > SuperLUMatrix
Definition: supermatrix.hh:332
Matrix::size_type size_type
Definition: supermatrix.hh:236