dune-istl 2.9.0
Loading...
Searching...
No Matches
supermatrix.hh
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
2// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4// vi: set et ts=4 sw=2 sts=2:
5#ifndef DUNE_ISTL_SUPERMATRIX_HH
6#define DUNE_ISTL_SUPERMATRIX_HH
7
8#if HAVE_SUPERLU
9
10#include "bcrsmatrix.hh"
11#include "bvector.hh"
12#include <dune/common/fmatrix.hh>
13#include <dune/common/fvector.hh>
14#include <dune/common/typetraits.hh>
15#include <limits>
16
18
19#include "superlufunctions.hh"
20
21namespace Dune
22{
23
24 template<class T>
27
28 template<class T>
30 {};
31
32#if __has_include("slu_sdefs.h")
33 template<>
35 {
36 static void create(SuperMatrix *mat, int n, int m, int offset,
37 float *values, int *rowindex, int* colindex,
38 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
39 {
40 sCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
41 stype, dtype, mtype);
42 }
43 };
44
45 template<>
46 struct SuperMatrixPrinter<float>
47 {
48 static void print(char* name, SuperMatrix* mat)
49 {
50 sPrint_CompCol_Matrix(name, mat);
51 }
52 };
53#endif
54
55#if __has_include("slu_ddefs.h")
56 template<>
57 struct SuperMatrixCreateSparseChooser<double>
58 {
59 static void create(SuperMatrix *mat, int n, int m, int offset,
60 double *values, int *rowindex, int* colindex,
61 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
62 {
63 dCreate_CompCol_Matrix(mat, n, m, offset, values, rowindex, colindex,
64 stype, dtype, mtype);
65 }
66 };
67
68 template<>
69 struct SuperMatrixPrinter<double>
70 {
71 static void print(char* name, SuperMatrix* mat)
72 {
73 dPrint_CompCol_Matrix(name, mat);
74 }
75 };
76#endif
77
78#if __has_include("slu_cdefs.h")
79 template<>
80 struct SuperMatrixCreateSparseChooser<std::complex<float> >
81 {
82 static void create(SuperMatrix *mat, int n, int m, int offset,
83 std::complex<float> *values, int *rowindex, int* colindex,
84 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
85 {
86 cCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast< ::complex*>(values),
87 rowindex, colindex, stype, dtype, mtype);
88 }
89 };
90
91 template<>
92 struct SuperMatrixPrinter<std::complex<float> >
93 {
94 static void print(char* name, SuperMatrix* mat)
95 {
96 cPrint_CompCol_Matrix(name, mat);
97 }
98 };
99#endif
100
101#if __has_include("slu_zdefs.h")
102 template<>
103 struct SuperMatrixCreateSparseChooser<std::complex<double> >
104 {
105 static void create(SuperMatrix *mat, int n, int m, int offset,
106 std::complex<double> *values, int *rowindex, int* colindex,
107 Stype_t stype, Dtype_t dtype, Mtype_t mtype)
108 {
109 zCreate_CompCol_Matrix(mat, n, m, offset, reinterpret_cast<doublecomplex*>(values),
110 rowindex, colindex, stype, dtype, mtype);
111 }
112 };
113
114 template<>
115 struct SuperMatrixPrinter<std::complex<double> >
116 {
117 static void print(char* name, SuperMatrix* mat)
118 {
119 zPrint_CompCol_Matrix(name, mat);
120 }
121 };
122#endif
123
124 template<class T>
126 {
127 static const Dtype_t type;
128 };
129
130 template<class T>
132 {};
133
134 template<class T>
135 const Dtype_t BaseGetSuperLUType<T>::type =
136 std::is_same<T,float>::value ? SLU_S :
137 ( std::is_same<T,std::complex<double> >::value ? SLU_Z :
138 ( std::is_same<T,std::complex<float> >::value ? SLU_C : SLU_D ));
139
140 template<>
141 struct GetSuperLUType<double>
142 : public BaseGetSuperLUType<double>
143 {
144 typedef double float_type;
145 };
146
147 template<>
148 struct GetSuperLUType<float>
149 : public BaseGetSuperLUType<float>
150 {
151 typedef float float_type;
152 };
153
154 template<>
155 struct GetSuperLUType<std::complex<double> >
156 : public BaseGetSuperLUType<std::complex<double> >
157 {
158 typedef double float_type;
159 };
160
161 template<>
162 struct GetSuperLUType<std::complex<float> >
163 : public BaseGetSuperLUType<std::complex<float> >
164 {
165 typedef float float_type;
166
167 };
168
173 template<class M>
175 {};
176
177 template<class M>
180
181 template<class T>
182 class SuperLU;
183
184 template<class M, class X, class TM, class TD, class T1>
186
187 template<class T, bool flag>
189
193 template<class B, class TA>
195 : public ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>
196 {
197 template<class M, class X, class TM, class TD, class T1>
199 friend struct SuperMatrixInitializer<BCRSMatrix<B,TA> >;
200 public:
203
205
206 typedef typename Matrix::size_type size_type;
207
212 explicit SuperLUMatrix(const Matrix& mat) : ISTL::Impl::BCCSMatrix<BCRSMatrix<B,TA>, int>(mat)
213 {}
214
215 SuperLUMatrix() : ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>()
216 {}
217
220 {
221 if (this->N_+this->M_*this->Nnz_ != 0)
222 free();
223 }
224
226 operator SuperMatrix&()
227 {
228 return A;
229 }
230
232 operator const SuperMatrix&() const
233 {
234 return A;
235 }
236
238 {
239 if (this->N_ + this->M_ + this->Nnz_!=0)
240 free();
241
242 using Matrix = BCRSMatrix<B,TA>;
245 ISTL::Impl::BCCSMatrixInitializer<Matrix, int> initializer(*this);
246
247 copyToBCCSMatrix(initializer, mat);
248
250 ::create(&A, this->N_, this->M_, this->colstart[this->N_],
251 this->values,this->rowindex, this->colstart, SLU_NC,
252 static_cast<Dtype_t>(GetSuperLUType<typename Matrix::field_type>::type), SLU_GE);
253 return *this;
254 }
255
257 {
258 if (this->N_ + this->M_ + this->Nnz_!=0)
259 free();
260
261 using Matrix = BCRSMatrix<B,TA>;
264 ISTL::Impl::BCCSMatrixInitializer<Matrix, int> initializer(*this);
265
266 copyToBCCSMatrix(initializer, mat);
267
269 ::create(&A, this->N_, this->M_, this->colstart[this->N_],
270 this->values,this->rowindex, this->colstart, SLU_NC,
271 static_cast<Dtype_t>(GetSuperLUType<B>::type), SLU_GE);
272 return *this;
273 }
274
281 virtual void setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs)
282 {
283 if(this->N_+this->M_+this->Nnz_!=0)
284 free();
285 this->N_=mrs.size()*MatrixDimension<typename Matrix::block_type>::rowdim(*(mat[0].begin()));
286 this->M_=mrs.size()*MatrixDimension<typename Matrix::block_type>::coldim(*(mat[0].begin()));
287 SuperMatrixInitializer<Matrix> initializer(*this);
288
289 copyToBCCSMatrix(initializer, ISTL::Impl::MatrixRowSubset<Matrix,std::set<std::size_t> >(mat,mrs));
290 }
291
293 virtual void setMatrix(const Matrix& mat)
294 {
297 SuperMatrixInitializer<Matrix> initializer(*this);
298
299 copyToBCCSMatrix(initializer, mat);
300 }
301
303 virtual void free()
304 {
305 ISTL::Impl::BCCSMatrix<typename BCRSMatrix<B,TA>::field_type, int>::free();
306 SUPERLU_FREE(A.Store);
307 }
308 private:
309 SuperMatrix A;
310 };
311
312 template<class B, class A>
314 : public ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>
315 {
316 template<class I, class S, class D>
318 public:
321
322 SuperMatrixInitializer(SuperLUMatrix& lum) : ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>(lum)
323 ,slumat(&lum)
324 {}
325
326 SuperMatrixInitializer() : ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>()
327 {}
328
329 virtual void createMatrix() const
330 {
331 ISTL::Impl::BCCSMatrixInitializer<BCRSMatrix<B,A>, int>::createMatrix();
333 ::create(&slumat->A, slumat->N_, slumat->M_, slumat->colstart[this->cols],
334 slumat->values,slumat->rowindex, slumat->colstart, SLU_NC,
335 static_cast<Dtype_t>(GetSuperLUType<typename Matrix::field_type>::type), SLU_GE);
336 }
337 private:
338 SuperLUMatrix* slumat;
339 };
340}
341#endif // HAVE_SUPERLU
342#endif
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Implementation of the BCRSMatrix class.
Matrix & mat
Definition matrixmatrix.hh:347
STL namespace.
Definition allocator.hh:11
Initializer for SuperLU Matrices representing the subdomains.
Definition overlappingschwarz.hh:47
static auto coldim(const M &A)
Definition matrixutils.hh:219
static auto rowdim(const M &A)
Definition matrixutils.hh:214
A sparse block matrix with compressed row storage.
Definition bcrsmatrix.hh:466
A::size_type size_type
The type for the index access and the size.
Definition bcrsmatrix.hh:500
Sequential overlapping Schwarz preconditioner.
Definition overlappingschwarz.hh:755
Definition overlappingschwarz.hh:694
SuperLu Solver.
Definition superlu.hh:271
Definition supermatrix.hh:26
Definition supermatrix.hh:30
Definition supermatrix.hh:126
static const Dtype_t type
Definition supermatrix.hh:127
Definition supermatrix.hh:132
double float_type
Definition supermatrix.hh:144
float float_type
Definition supermatrix.hh:151
double float_type
Definition supermatrix.hh:158
float float_type
Definition supermatrix.hh:165
Utility class for converting an ISTL Matrix into a SuperLU Matrix.
Definition supermatrix.hh:175
Definition supermatrix.hh:179
virtual void free()
free allocated space.
Definition supermatrix.hh:303
SuperLUMatrix< BCRSMatrix< B, TA > > & operator=(const SuperLUMatrix< BCRSMatrix< B, TA > > &mat)
Definition supermatrix.hh:256
SuperLUMatrix< BCRSMatrix< B, TA > > & operator=(const BCRSMatrix< B, TA > &mat)
Definition supermatrix.hh:237
SuperLUMatrix(const Matrix &mat)
Constructor that initializes the data.
Definition supermatrix.hh:212
virtual void setMatrix(const Matrix &mat)
Initialize data from given matrix.
Definition supermatrix.hh:293
SuperLUMatrix()
Definition supermatrix.hh:215
Matrix::size_type size_type
Definition supermatrix.hh:206
BCRSMatrix< B, TA > Matrix
The type of the matrix to convert.
Definition supermatrix.hh:202
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:281
virtual ~SuperLUMatrix()
Destructor.
Definition supermatrix.hh:219
Dune::SuperLUMatrix< Matrix > SuperLUMatrix
Definition supermatrix.hh:320
BCRSMatrix< B, A > Matrix
Definition supermatrix.hh:319
SuperMatrixInitializer()
Definition supermatrix.hh:326
virtual void createMatrix() const
Definition supermatrix.hh:329
SuperMatrixInitializer(SuperLUMatrix &lum)
Definition supermatrix.hh:322