dune-istl  2.8.0
supermatrix.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_ISTL_SUPERMATRIX_HH
4 #define DUNE_ISTL_SUPERMATRIX_HH
5 
6 #if HAVE_SUPERLU
7 
8 #include "bcrsmatrix.hh"
9 #include "bvector.hh"
10 #include <dune/common/fmatrix.hh>
11 #include <dune/common/fvector.hh>
12 #include <dune/common/typetraits.hh>
13 #include <limits>
14 
16 
17 #include "superlufunctions.hh"
18 
19 namespace Dune
20 {
21 
22  template<class T>
24  {};
25 
26  template<class T>
28  {};
29 
30 #if __has_include("slu_sdefs.h")
31  template<>
32  struct SuperMatrixCreateSparseChooser<float>
33  {
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)
37  {
38  sCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
39  stype, dtype, mtype);
40  }
41  };
42 
43  template<>
44  struct SuperMatrixPrinter<float>
45  {
46  static void print(char* name, SuperMatrix* mat)
47  {
48  sPrint_CompCol_Matrix(name, mat);
49  }
50  };
51 #endif
52 
53 #if __has_include("slu_ddefs.h")
54  template<>
55  struct SuperMatrixCreateSparseChooser<double>
56  {
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)
60  {
61  dCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
62  stype, dtype, mtype);
63  }
64  };
65 
66  template<>
67  struct SuperMatrixPrinter<double>
68  {
69  static void print(char* name, SuperMatrix* mat)
70  {
71  dPrint_CompCol_Matrix(name, mat);
72  }
73  };
74 #endif
75 
76 #if __has_include("slu_cdefs.h")
77  template<>
78  struct SuperMatrixCreateSparseChooser<std::complex<float> >
79  {
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)
83  {
84  cCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast< ::complex*>(values),
85  rowindex, colindex, stype, dtype, mtype);
86  }
87  };
88 
89  template<>
90  struct SuperMatrixPrinter<std::complex<float> >
91  {
92  static void print(char* name, SuperMatrix* mat)
93  {
94  cPrint_CompCol_Matrix(name, mat);
95  }
96  };
97 #endif
98 
99 #if __has_include("slu_zdefs.h")
100  template<>
101  struct SuperMatrixCreateSparseChooser<std::complex<double> >
102  {
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)
106  {
107  zCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast<doublecomplex*>(values),
108  rowindex, colindex, stype, dtype, mtype);
109  }
110  };
111 
112  template<>
113  struct SuperMatrixPrinter<std::complex<double> >
114  {
115  static void print(char* name, SuperMatrix* mat)
116  {
117  zPrint_CompCol_Matrix(name, mat);
118  }
119  };
120 #endif
121 
122  template<class T>
124  {
125  static const Dtype_t type;
126  };
127 
128  template<class T>
130  {};
131 
132  template<class T>
133  const Dtype_t BaseGetSuperLUType<T>::type =
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 ));
137 
138  template<>
139  struct GetSuperLUType<double>
140  : public BaseGetSuperLUType<double>
141  {
142  typedef double float_type;
143  };
144 
145  template<>
146  struct GetSuperLUType<float>
147  : public BaseGetSuperLUType<float>
148  {
149  typedef float float_type;
150  };
151 
152  template<>
153  struct GetSuperLUType<std::complex<double> >
154  : public BaseGetSuperLUType<std::complex<double> >
155  {
156  typedef double float_type;
157  };
158 
159  template<>
160  struct GetSuperLUType<std::complex<float> >
161  : public BaseGetSuperLUType<std::complex<float> >
162  {
163  typedef float float_type;
164 
165  };
166 
171  template<class M>
173  {};
174 
175  template<class M>
177  {};
178 
179  template<class T>
180  class SuperLU;
181 
182  template<class M, class X, class TM, class TD, class T1>
183  class SeqOverlappingSchwarz;
184 
185  template<class T, bool flag>
187 
191  template<class B, class TA>
192  class SuperLUMatrix<BCRSMatrix<B,TA> >
193  : public ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>
194  {
195  template<class M, class X, class TM, class TD, class T1>
196  friend class SeqOverlappingSchwarz;
197  friend struct SuperMatrixInitializer<BCRSMatrix<B,TA> >;
198  public:
201 
203 
204  typedef typename Matrix::size_type size_type;
205 
210  explicit SuperLUMatrix(const Matrix& mat) : ISTL::Impl::BCCSMatrix<BCRSMatrix<B,TA>, int>(mat)
211  {}
212 
213  SuperLUMatrix() : ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>()
214  {}
215 
217  virtual ~SuperLUMatrix()
218  {
219  if (this->N_+this->M_*this->Nnz_ != 0)
220  free();
221  }
222 
224  operator SuperMatrix&()
225  {
226  return A;
227  }
228 
230  operator const SuperMatrix&() const
231  {
232  return A;
233  }
234 
236  {
237  if (this->N_ + this->M_ + this->Nnz_!=0)
238  free();
239 
240  using Matrix = BCRSMatrix<B,TA>;
243  ISTL::Impl::BCCSMatrixInitializer<Matrix, int> initializer(*this);
244 
245  copyToBCCSMatrix(initializer, mat);
246 
248  ::create(&A, this->N_, this->M_, this->colstart[this->N_],
249  this->values,this->rowindex, this->colstart, SLU_NC,
250  static_cast<Dtype_t>(GetSuperLUType<typename Matrix::field_type>::type), SLU_GE);
251  return *this;
252  }
253 
255  {
256  if (this->N_ + this->M_ + this->Nnz_!=0)
257  free();
258 
259  using Matrix = BCRSMatrix<B,TA>;
262  ISTL::Impl::BCCSMatrixInitializer<Matrix, int> initializer(*this);
263 
264  copyToBCCSMatrix(initializer, mat);
265 
267  ::create(&A, this->N_, this->M_, this->colstart[this->N_],
268  this->values,this->rowindex, this->colstart, SLU_NC,
269  static_cast<Dtype_t>(GetSuperLUType<B>::type), SLU_GE);
270  return *this;
271  }
272 
279  virtual void setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs)
280  {
281  if(this->N_+this->M_+this->Nnz_!=0)
282  free();
283  this->N_=mrs.size()*MatrixDimension<typename Matrix::block_type>::rowdim(*(mat[0].begin()));
284  this->M_=mrs.size()*MatrixDimension<typename Matrix::block_type>::coldim(*(mat[0].begin()));
285  SuperMatrixInitializer<Matrix> initializer(*this);
286 
287  copyToBCCSMatrix(initializer, ISTL::Impl::MatrixRowSubset<Matrix,std::set<std::size_t> >(mat,mrs));
288  }
289 
291  virtual void setMatrix(const Matrix& mat)
292  {
295  SuperMatrixInitializer<Matrix> initializer(*this);
296 
297  copyToBCCSMatrix(initializer, mat);
298  }
299 
301  virtual void free()
302  {
303  ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>::free();
304  SUPERLU_FREE(A.Store);
305  }
306  private:
307  SuperMatrix A;
308  };
309 
310  template<class B, class A>
312  : public ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>
313  {
314  template<class I, class S, class D>
316  public:
319 
320  SuperMatrixInitializer(SuperLUMatrix& lum) : ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>(lum)
321  ,slumat(&lum)
322  {}
323 
324  SuperMatrixInitializer() : ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>()
325  {}
326 
327  virtual void createMatrix() const
328  {
329  ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>::createMatrix();
331  ::create(&slumat->A, slumat->N_, slumat->M_, slumat->colstart[this->cols],
332  slumat->values,slumat->rowindex, slumat->colstart, SLU_NC,
333  static_cast<Dtype_t>(GetSuperLUType<typename Matrix::field_type>::type), SLU_GE);
334  }
335  private:
336  SuperLUMatrix* slumat;
337  };
338 }
339 #endif // HAVE_SUPERLU
340 #endif
Implementation of the BCRSMatrix class.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Matrix & mat
Definition: matrixmatrix.hh:345
Definition: allocator.hh:9
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:45
static auto coldim(const M &A)
Definition: matrixutils.hh:217
static auto rowdim(const M &A)
Definition: matrixutils.hh:212
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:464
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:498
Sequential overlapping Schwarz preconditioner.
Definition: overlappingschwarz.hh:753
Definition: overlappingschwarz.hh:692
SuperLu Solver.
Definition: superlu.hh:269
Definition: supermatrix.hh:24
Definition: supermatrix.hh:28
Definition: supermatrix.hh:124
static const Dtype_t type
Definition: supermatrix.hh:125
Definition: supermatrix.hh:130
double float_type
Definition: supermatrix.hh:142
float float_type
Definition: supermatrix.hh:149
double float_type
Definition: supermatrix.hh:156
float float_type
Definition: supermatrix.hh:163
Utility class for converting an ISTL Matrix into a SuperLU Matrix.
Definition: supermatrix.hh:173
Definition: supermatrix.hh:177
virtual void free()
free allocated space.
Definition: supermatrix.hh:301
SuperLUMatrix(const Matrix &mat)
Constructor that initializes the data.
Definition: supermatrix.hh:210
SuperLUMatrix< BCRSMatrix< B, TA > > & operator=(const BCRSMatrix< B, TA > &mat)
Definition: supermatrix.hh:235
virtual void setMatrix(const Matrix &mat)
Initialize data from given matrix.
Definition: supermatrix.hh:291
SuperLUMatrix()
Definition: supermatrix.hh:213
Matrix::size_type size_type
Definition: supermatrix.hh:204
BCRSMatrix< B, TA > Matrix
The type of the matrix to convert.
Definition: supermatrix.hh:200
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:279
virtual ~SuperLUMatrix()
Destructor.
Definition: supermatrix.hh:217
SuperLUMatrix< BCRSMatrix< B, TA > > & operator=(const SuperLUMatrix< BCRSMatrix< B, TA > > &mat)
Definition: supermatrix.hh:254
Dune::SuperLUMatrix< Matrix > SuperLUMatrix
Definition: supermatrix.hh:318
BCRSMatrix< B, A > Matrix
Definition: supermatrix.hh:317
SuperMatrixInitializer()
Definition: supermatrix.hh:324
virtual void createMatrix() const
Definition: supermatrix.hh:327
SuperMatrixInitializer(SuperLUMatrix &lum)
Definition: supermatrix.hh:320