My Project
ParallelOverlappingILU0.hpp
1/*
2 Copyright 2015 Dr. Blatt - HPC-Simulation-Software & Services
3 Copyright 2015 Statoil AS
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20#ifndef OPM_PARALLELOVERLAPPINGILU0_HEADER_INCLUDED
21#define OPM_PARALLELOVERLAPPINGILU0_HEADER_INCLUDED
22
23#include <opm/simulators/linalg/MILU.hpp>
24#include <opm/simulators/linalg/PreconditionerWithUpdate.hpp>
25#include <dune/common/version.hh>
26#include <dune/istl/paamg/smoother.hh>
27
28#include <cstddef>
29#include <vector>
30#include <type_traits>
31
32namespace Opm
33{
34
35//template<class M, class X, class Y, class C>
36//class ParallelOverlappingILU0;
37template<class Matrix, class Domain, class Range, class ParallelInfo>
38class ParallelOverlappingILU0;
39
40template<class F>
42 : public Dune::Amg::DefaultSmootherArgs<F>
43{
44 public:
46 : milu_(milu), n_(0)
47 {}
48 void setMilu(MILU_VARIANT milu)
49 {
50 milu_ = milu;
51 }
52 MILU_VARIANT getMilu() const
53 {
54 return milu_;
55 }
56 void setN(int n)
57 {
58 n_ = n;
59 }
60 int getN() const
61 {
62 return n_;
63 }
64 private:
65 MILU_VARIANT milu_;
66 int n_;
67};
68} // end namespace Opm
69
70namespace Dune
71{
72
73namespace Amg
74{
75
76
77template<class M, class X, class Y, class C>
78struct SmootherTraits<Opm::ParallelOverlappingILU0<M,X,Y,C> >
79{
81};
82
89template<class Matrix, class Domain, class Range, class ParallelInfo>
90struct ConstructionTraits<Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> >
91{
93 using Arguments = DefaultParallelConstructionArgs<T,ParallelInfo>;
94
95#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
96 using ParallelOverlappingILU0Pointer = std::shared_ptr<T>;
97#else
99#endif
100
101 static inline ParallelOverlappingILU0Pointer construct(Arguments& args)
102 {
104 new T(args.getMatrix(),
105 args.getComm(),
106 args.getArgs().getN(),
107 args.getArgs().relaxationFactor,
108 args.getArgs().getMilu()) );
109 }
110
111#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
112 // this method is not needed anymore in 2.7 since std::shared_ptr is used
113 static inline void deconstruct(T* bp)
114 {
115 delete bp;
116 }
117#endif
118
119};
120
121} // end namespace Amg
122
123} // end namespace Dune
124
125namespace Opm
126{
142template<class Matrix, class Domain, class Range, class ParallelInfoT>
144 : public Dune::PreconditionerWithUpdate<Domain,Range>
145{
146 using ParallelInfo = ParallelInfoT;
147
148public:
150 using matrix_type = typename std::remove_const<Matrix>::type;
152 using domain_type = Domain;
154 using range_type = Range;
156 using field_type = typename Domain::field_type;
157
158 using block_type = typename matrix_type::block_type;
159 using size_type = typename matrix_type::size_type;
160
161protected:
162 struct CRS
163 {
164 CRS() : nRows_( 0 ) {}
165
166 size_type rows() const { return nRows_; }
167
168 size_type nonZeros() const
169 {
170 assert( rows_[ rows() ] != size_type(-1) );
171 return rows_[ rows() ];
172 }
173
174 void resize( const size_type nRows )
175 {
176 if( nRows_ != nRows )
177 {
178 nRows_ = nRows ;
179 rows_.resize( nRows_+1, size_type(-1) );
180 }
181 }
182
183 void reserveAdditional( const size_type nonZeros )
184 {
185 const size_type needed = values_.size() + nonZeros ;
186 if( values_.capacity() < needed )
187 {
188 const size_type estimate = needed * 1.1;
189 values_.reserve( estimate );
190 cols_.reserve( estimate );
191 }
192 }
193
194 void push_back( const block_type& value, const size_type index )
195 {
196 values_.push_back( value );
197 cols_.push_back( index );
198 }
199
200 void clear()
201 {
202 rows_.clear();
203 values_.clear();
204 cols_.clear();
205 nRows_= 0;
206 }
207
208 std::vector<size_type> rows_;
209 std::vector<block_type> values_;
210 std::vector<size_type> cols_;
211 size_type nRows_;
212 };
213
214public:
215 Dune::SolverCategory::Category category() const override;
216
230 ParallelOverlappingILU0 (const Matrix& A,
231 const int n, const field_type w,
232 MILU_VARIANT milu, bool redblack = false,
233 bool reorder_sphere = true);
234
247 ParallelOverlappingILU0 (const Matrix& A,
248 const ParallelInfo& comm, const int n, const field_type w,
249 MILU_VARIANT milu, bool redblack = false,
250 bool reorder_sphere = true);
251
264 ParallelOverlappingILU0 (const Matrix& A,
265 const field_type w, MILU_VARIANT milu,
266 bool redblack = false,
267 bool reorder_sphere = true);
268
282 ParallelOverlappingILU0 (const Matrix& A,
283 const ParallelInfo& comm, const field_type w,
284 MILU_VARIANT milu, bool redblack = false,
285 bool reorder_sphere = true);
286
301 ParallelOverlappingILU0 (const Matrix& A,
302 const ParallelInfo& comm,
303 const field_type w, MILU_VARIANT milu,
304 size_type interiorSize, bool redblack = false,
305 bool reorder_sphere = true);
306
312 void pre (Domain&, Range&) override
313 {}
314
320 void apply (Domain& v, const Range& d) override;
321
322 template <class V>
323 void copyOwnerToAll( V& v ) const;
324
330 void post (Range&) override
331 {}
332
333 void update() override;
334
335protected:
337 Range& reorderD(const Range& d);
338
340 Domain& reorderV(Domain& v);
341
342 void reorderBack(const Range& reorderedV, Range& v);
343
346 CRS upper_;
347 std::vector< block_type > inv_;
349 std::vector< std::size_t > ordering_;
354
355 const ParallelInfo* comm_;
358 const bool relaxation_;
359 size_type interiorSize_;
360 const Matrix* A_;
361 int iluIteration_;
362 MILU_VARIANT milu_;
363 bool redBlack_;
364 bool reorderSphere_;
365};
366
367} // end namespace Opm
368#endif
Interface class adding the update() method to the preconditioner interface.
Definition: PreconditionerWithUpdate.hpp:32
Definition: ParallelOverlappingILU0.hpp:43
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
Definition: ParallelOverlappingILU0.hpp:145
ParallelOverlappingILU0(const Matrix &A, const int n, const field_type w, MILU_VARIANT milu, bool redblack=false, bool reorder_sphere=true)
Constructor.
Definition: ParallelOverlappingILU0_impl.hpp:197
Domain & reorderV(Domain &v)
Reorder V if needed and return a reference to it.
Definition: ParallelOverlappingILU0_impl.hpp:533
typename std::remove_const< Matrix >::type matrix_type
The matrix type the preconditioner is for.
Definition: ParallelOverlappingILU0.hpp:150
Domain domain_type
The domain type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:152
std::vector< std::size_t > ordering_
the reordering of the unknowns
Definition: ParallelOverlappingILU0.hpp:349
Domain reorderedV_
The reordered left hand side.
Definition: ParallelOverlappingILU0.hpp:353
void post(Range &) override
Clean up.
Definition: ParallelOverlappingILU0.hpp:330
Range reorderedD_
The reordered right hand side.
Definition: ParallelOverlappingILU0.hpp:351
const field_type w_
The relaxation factor to use.
Definition: ParallelOverlappingILU0.hpp:357
Range & reorderD(const Range &d)
Reorder D if needed and return a reference to it.
Definition: ParallelOverlappingILU0_impl.hpp:509
Range range_type
The range type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:154
CRS lower_
The ILU0 decomposition of the matrix.
Definition: ParallelOverlappingILU0.hpp:345
void apply(Domain &v, const Range &d) override
Apply the preconditoner.
Definition: ParallelOverlappingILU0_impl.hpp:286
void pre(Domain &, Range &) override
Prepare the preconditioner.
Definition: ParallelOverlappingILU0.hpp:312
typename Domain::field_type field_type
The field type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:156
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
MILU_VARIANT
Definition: MILU.hpp:34
@ ILU
Do not perform modified ILU.
Definition: ParallelOverlappingILU0.hpp:163