dune-fem 2.8.0
Loading...
Searching...
No Matches
istladapter.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_OPERATOR_LINEAR_ISTLADAPTER_HH
2#define DUNE_FEM_OPERATOR_LINEAR_ISTLADAPTER_HH
3
4#if HAVE_DUNE_ISTL
5#include <dune/common/version.hh>
6
7#include <dune/istl/operators.hh>
8
11
12namespace Dune
13{
14
15 namespace Fem
16 {
17
18 // ISTLLinearOperatorAdapter
19 // -------------------------
20
21 template< class Operator >
22 class ISTLLinearOperatorAdapter
23 : public Dune::LinearOperator< typename Operator::DomainFunctionType::DofStorageType, typename Operator::RangeFunctionType::DofStorageType >
24 {
25 typedef ISTLLinearOperatorAdapter< Operator > ThisType;
26 typedef Dune::LinearOperator< typename Operator::DomainFunctionType::DofStorageType, typename Operator::RangeFunctionType::DofStorageType > BaseType;
27
28 typedef typename Operator::DomainFunctionType DomainFunctionType;
29 typedef typename Operator::RangeFunctionType RangeFunctionType;
30
31 public:
32#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
33 enum {category=SolverCategory::sequential};
34#endif // #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
35
36 typedef Operator OperatorType;
37
38 typedef typename DomainFunctionType::DiscreteFunctionSpaceType DomainFunctionSpaceType;
39 typedef typename RangeFunctionType::DiscreteFunctionSpaceType RangeFunctionSpaceType;
40
41 typedef typename BaseType::domain_type domain_type;
42 typedef typename BaseType::range_type range_type;
43 typedef typename BaseType::field_type field_type;
44
45 ISTLLinearOperatorAdapter ( const OperatorType &op,
46 const DomainFunctionSpaceType &domainSpace,
47 const RangeFunctionSpaceType &rangeSpace )
48 : op_( op ),
49 domainSpace_( domainSpace ),
50 rangeSpace_( rangeSpace )
51 {}
52
53 virtual void apply ( const domain_type &x, range_type &y ) const override
54 {
55 const DomainFunctionType fx( "ISTLLinearOperatorAdapter::apply::x", domainSpace_, x );
56 RangeFunctionType fy( "ISTLLinearOperatorAdapter::apply::y", rangeSpace_, y );
57 op_( fx, fy );
58 }
59
60 virtual void applyscaleadd ( field_type alpha, const domain_type &x, range_type &y ) const override
61 {
62 const DomainFunctionType fx( "ISTLLinearOperatorAdapter::applyscaleadd::x", domainSpace_, x );
63 RangeFunctionType fy( "ISTLLinearOperatorAdapter::applyscaleadd::y", rangeSpace_ );
64 op_( fx, fy );
65 y.axpy( alpha, fy.blockVector() );
66 }
67
68#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
69 SolverCategory::Category category () const override { return SolverCategory::sequential; }
70#endif // #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
71
72 protected:
73 const OperatorType &op_;
74 const DomainFunctionSpaceType &domainSpace_;
75 const RangeFunctionSpaceType &rangeSpace_;
76 };
77
78 template< class Operator >
79 class ISTLMatrixFreeOperatorAdapter : public ISTLLinearOperatorAdapter< Operator >
80 {
81 typedef ISTLMatrixFreeOperatorAdapter< Operator > ThisType;
82 typedef ISTLLinearOperatorAdapter< Operator > BaseType;
83
84 typedef typename Operator::DomainFunctionType DomainFunctionType;
85 typedef typename Operator::RangeFunctionType RangeFunctionType;
86
87 public:
88#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
89 enum {category=SolverCategory::sequential};
90#endif // #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
91
92 typedef Operator OperatorType;
93
94 typedef typename DomainFunctionType::DiscreteFunctionSpaceType DomainFunctionSpaceType;
95 typedef typename RangeFunctionType::DiscreteFunctionSpaceType RangeFunctionSpaceType;
96
97 template <class DF>
98 struct IsISTLBlockVectorDiscreteFunction
99 {
100 static const bool value = false ;
101 };
102
103 template <class Space, class Block>
104 struct IsISTLBlockVectorDiscreteFunction< ISTLBlockVectorDiscreteFunction< Space, Block > >
105 {
106 static const bool value = true ;
107 };
108
109
110 static_assert( IsISTLBlockVectorDiscreteFunction< DomainFunctionType > :: value,
111 "ISTLMatrixFreeOperatorAdapter only works with ISTLBlockVectorDiscreteFunction" );
112 static_assert( IsISTLBlockVectorDiscreteFunction< RangeFunctionType > :: value,
113 "ISTLMatrixFreeOperatorAdapter only works with ISTLBlockVectorDiscreteFunction" );
114
115 typedef typename BaseType::domain_type domain_type;
116 typedef typename BaseType::range_type range_type;
117 typedef typename BaseType::field_type field_type;
118
119 typedef Fem::ParallelScalarProduct< RangeFunctionType > ParallelScalarProductType;
120 typedef Fem::IdentityPreconditionerWrapper< domain_type, range_type > PreconditionAdapterType;
121
123 ISTLMatrixFreeOperatorAdapter ( const OperatorType &op,
124 const DomainFunctionSpaceType &domainSpace,
125 const RangeFunctionSpaceType &rangeSpace )
126 : BaseType( op, domainSpace, rangeSpace ),
127 scp_( rangeSpace ),
128 preconditioner_()
129 {}
130
132 ISTLMatrixFreeOperatorAdapter ( const ISTLMatrixFreeOperatorAdapter& other )
133 : BaseType( other ),
134 scp_( other.rangeSpace_ ),
135 preconditioner_()
136 {}
137
139 PreconditionAdapterType& preconditionAdapter() { return preconditioner_; }
141 const PreconditionAdapterType& preconditionAdapter() const { return preconditioner_; }
142
143 virtual double residuum(const range_type& rhs, domain_type& x) const
144 {
145 range_type tmp( rhs );
146
147 this->apply(x,tmp);
148 tmp -= rhs;
149
150 // return global sum of residuum
151 return scp_.norm(tmp);
152 }
153
154#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
155 SolverCategory::Category category () const override { return SolverCategory::sequential; }
156#endif // #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
157
159 ParallelScalarProductType& scp() const { return scp_; }
160
162 double averageCommTime () const { return 0; }
163
164 protected:
165 mutable ParallelScalarProductType scp_;
166 PreconditionAdapterType preconditioner_;
167 };
168
169 } // namespace Fem
170
171} // namespace Dune
172
173#endif // #if HAVE_DUNE_ISTL
174
175#endif // #ifndef DUNE_FEM_OPERATOR_LINEAR_ISTLADAPTER_HH
Definition: bindguard.hh:11