1#ifndef DUNE_FEM_ISTLMATRIXADAPTER_HH
2#define DUNE_FEM_ISTLMATRIXADAPTER_HH
6#include <dune/common/timer.hh>
7#include <dune/common/version.hh>
9#include <dune/istl/operators.hh>
23 template <
class Matrix>
24 class PreconditionerWrapper ;
27 template <
class MatrixImp>
28 class ISTLParallelMatrixAdapterInterface
29 :
public AssembledLinearOperator< MatrixImp,
30 typename MatrixImp :: RowBlockVectorType,
31 typename MatrixImp :: ColBlockVectorType>
33 typedef ISTLParallelMatrixAdapterInterface< MatrixImp > ThisType;
35 typedef MatrixImp MatrixType;
36 typedef Fem::PreconditionerWrapper< MatrixType > PreconditionAdapterType;
38 typedef typename MatrixType :: RowDiscreteFunctionType RowDiscreteFunctionType;
39 typedef typename MatrixType :: ColDiscreteFunctionType ColumnDiscreteFunctionType;
41 typedef typename RowDiscreteFunctionType :: DiscreteFunctionSpaceType RowSpaceType;
43 typedef typename ColumnDiscreteFunctionType :: DiscreteFunctionSpaceType ColSpaceType;
44 typedef Fem::ParallelScalarProduct<ColumnDiscreteFunctionType> ParallelScalarProductType;
46 typedef typename RowDiscreteFunctionType :: DofStorageType X;
47 typedef typename ColumnDiscreteFunctionType :: DofStorageType Y;
50 typedef MatrixType matrix_type;
51 typedef X domain_type;
53 typedef typename X::field_type field_type;
55#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
57 enum { category=SolverCategory::sequential };
62 ISTLParallelMatrixAdapterInterface (
const ISTLParallelMatrixAdapterInterface &org )
63 : matrix_( org.matrix_ ),
64 rowSpace_( org.rowSpace_ ),
65 colSpace_( org.colSpace_ ),
67 preconditioner_( org.preconditioner_ ),
68 averageCommTime_( org.averageCommTime_ )
73 ISTLParallelMatrixAdapterInterface ( MatrixType &matrix,
74 const RowSpaceType &rowSpace,
75 const ColSpaceType &colSpace,
76 const PreconditionAdapterType& precon )
78 rowSpace_( rowSpace ),
79 colSpace_( colSpace ),
81 preconditioner_( precon ),
86 virtual double averageCommTime()
const {
return averageCommTime_ ; }
89 virtual PreconditionAdapterType &preconditionAdapter() {
return preconditioner_; }
92 virtual const PreconditionAdapterType &preconditionAdapter()
const {
return preconditioner_; }
95 virtual ParallelScalarProductType &scp() {
return scp_; }
98 void apply (
const X &x, Y &y )
const override
100 DUNE_THROW(NotImplemented,
"interface method apply not overloaded!");
104 void applyscaleadd ( field_type alpha,
const X &x, Y &y)
const override
106 DUNE_THROW(NotImplemented,
"interface method applyscaleadd not overloaded!");
110 virtual double residuum(
const Y& rhs, X& x)
const
112 DUNE_THROW(NotImplemented,
"interface method residuum not overloaded!");
117 const MatrixType& getmat ()
const override {
return matrix_; }
119#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
120 SolverCategory::Category category ()
const override {
return SolverCategory::sequential; }
124 void communicate( Y &y )
const
126 if( rowSpace_.grid().comm().size() <= 1 )
129 Dune::Timer commTime;
130 ColumnDiscreteFunctionType tmp(
"LagrangeParallelMatrixAdapter::communicate", colSpace_, y );
131 colSpace_.communicate( tmp );
132 averageCommTime_ += commTime.elapsed();
136 const RowSpaceType &rowSpace_;
137 const ColSpaceType &colSpace_;
138 mutable ParallelScalarProductType scp_;
139 PreconditionAdapterType preconditioner_;
140 mutable double averageCommTime_;
143 template <
class MatrixImp>
144 class LagrangeParallelMatrixAdapter
145 :
public ISTLParallelMatrixAdapterInterface< MatrixImp >
147 typedef LagrangeParallelMatrixAdapter< MatrixImp > ThisType;
148 typedef ISTLParallelMatrixAdapterInterface< MatrixImp > BaseType;
150 typedef MatrixImp MatrixType;
151 typedef Fem::PreconditionerWrapper<MatrixType> PreconditionAdapterType;
153 typedef typename MatrixType :: RowDiscreteFunctionType RowDiscreteFunctionType;
154 typedef typename MatrixType :: ColDiscreteFunctionType ColumnDiscreteFunctionType;
156 typedef typename RowDiscreteFunctionType :: DiscreteFunctionSpaceType RowSpaceType;
158 typedef typename ColumnDiscreteFunctionType :: DiscreteFunctionSpaceType ColSpaceType;
159 typedef Fem::ParallelScalarProduct<ColumnDiscreteFunctionType> ParallelScalarProductType;
161 typedef typename RowDiscreteFunctionType :: DofStorageType X;
162 typedef typename ColumnDiscreteFunctionType :: DofStorageType Y;
165 typedef MatrixType matrix_type;
166 typedef X domain_type;
167 typedef Y range_type;
168 typedef typename X::field_type field_type;
170#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
172 enum { category=SolverCategory::sequential };
176 using BaseType :: matrix_;
177 using BaseType :: scp_ ;
178 using BaseType :: colSpace_ ;
179 using BaseType :: rowSpace_ ;
180 using BaseType :: averageCommTime_;
184 LagrangeParallelMatrixAdapter (
const LagrangeParallelMatrixAdapter &org )
190 LagrangeParallelMatrixAdapter ( MatrixType &matrix,
191 const RowSpaceType &rowSpace,
192 const ColSpaceType &colSpace,
193 const PreconditionAdapterType& precon )
194 : BaseType( matrix, rowSpace, colSpace, precon )
198 void apply (
const X &x, Y &y )
const override
205 void applyscaleadd ( field_type alpha,
const X &x, Y &y)
const override
207 if( rowSpace_.grid().comm().size() <= 1 )
209 matrix_.usmv(alpha,x,y);
216 y.axpy( alpha, tmp );
220 double residuum(
const Y& rhs, X& x)
const override
228 return scp_.norm(tmp);
231#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
232 SolverCategory::Category category ()
const override {
return SolverCategory::sequential; }
236 void communicate( Y &y )
const
238 if( rowSpace_.grid().comm().size() <= 1 )
241 Dune::Timer commTime;
242 ColumnDiscreteFunctionType tmp(
"LagrangeParallelMatrixAdapter::communicate", colSpace_, y );
243 colSpace_.communicate( tmp );
244 averageCommTime_ += commTime.elapsed();
248 template <
class MatrixImp>
249 class DGParallelMatrixAdapter
250 :
public ISTLParallelMatrixAdapterInterface< MatrixImp >
252 typedef DGParallelMatrixAdapter< MatrixImp > ThisType ;
253 typedef ISTLParallelMatrixAdapterInterface< MatrixImp > BaseType;
255 typedef MatrixImp MatrixType;
256 typedef Fem::PreconditionerWrapper<MatrixType> PreconditionAdapterType;
258 typedef typename MatrixType :: RowDiscreteFunctionType RowDiscreteFunctionType;
259 typedef typename MatrixType :: ColDiscreteFunctionType ColumnDiscreteFunctionType;
261 typedef typename RowDiscreteFunctionType :: DiscreteFunctionSpaceType RowSpaceType;
263 typedef typename ColumnDiscreteFunctionType :: DiscreteFunctionSpaceType ColSpaceType;
264 typedef Fem::ParallelScalarProduct<ColumnDiscreteFunctionType> ParallelScalarProductType;
266 typedef typename RowDiscreteFunctionType :: DofStorageType X;
267 typedef typename ColumnDiscreteFunctionType :: DofStorageType Y;
270 typedef MatrixType matrix_type;
271 typedef X domain_type;
272 typedef Y range_type;
273 typedef typename X::field_type field_type;
275#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
277 enum { category=SolverCategory::sequential };
281 using BaseType :: matrix_;
282 using BaseType :: scp_;
283 using BaseType :: rowSpace_;
284 using BaseType :: colSpace_;
285 using BaseType :: averageCommTime_;
289 DGParallelMatrixAdapter (
const DGParallelMatrixAdapter& org)
294 DGParallelMatrixAdapter (MatrixType& A,
295 const RowSpaceType& rowSpace,
296 const ColSpaceType& colSpace,
297 const PreconditionAdapterType& precon )
298 : BaseType( A, rowSpace, colSpace, precon )
302 void apply (
const X& x, Y& y)
const override
311 scp_.deleteNonInterior( y );
315 void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
321 matrix_.usmv(alpha,x,y);
324 scp_.deleteNonInterior( y );
327 double residuum(
const Y& rhs, X& x)
const override
332 typedef typename ParallelScalarProductType :: AuxiliaryDofsType AuxiliaryDofsType;
333 const AuxiliaryDofsType& auxiliaryDofs = scp_.auxiliaryDofs();
335 typedef typename Y :: block_type LittleBlockVectorType;
336 LittleBlockVectorType tmp;
341 const int auxiliarySize = auxiliaryDofs.size();
342 for(
int auxiliary = 0; auxiliary<auxiliarySize; ++auxiliary)
344 const int nextAuxiliary = auxiliaryDofs[auxiliary];
345 for(; i<nextAuxiliary; ++i)
349 typedef typename MatrixType :: row_type row_type;
351 const row_type& row = matrix_[i];
353 typedef typename MatrixType :: ConstColIterator ConstColIterator;
354 ConstColIterator endj = row.end();
355 for (ConstColIterator j = row.begin(); j!=endj; ++j)
357 (*j).umv(x[j.index()], tmp);
364 res += tmp.two_norm2();
369 res = rowSpace_.grid().comm().sum( res );
374#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
375 SolverCategory::Category category ()
const override {
return SolverCategory::sequential; }
379 void communicate(
const X& x)
const
381 if( rowSpace_.grid().comm().size() <= 1 ) return ;
383 Dune::Timer commTime;
386 RowDiscreteFunctionType tmp (
"DGParallelMatrixAdapter::communicate",
390 rowSpace_.communicate( tmp );
393 averageCommTime_ += commTime.elapsed();
double sqrt(const Dune::Fem::Double &v)
Definition: double.hh:977
Definition: bindguard.hh:11