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 
32 namespace Opm
33 {
34 
35 //template<class M, class X, class Y, class C>
36 //class ParallelOverlappingILU0;
37 template<class Matrix, class Domain, class Range, class ParallelInfo>
38 class ParallelOverlappingILU0;
39 
40 template<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 
70 namespace Dune
71 {
72 
73 namespace Amg
74 {
75 
76 
77 template<class M, class X, class Y, class C>
78 struct SmootherTraits<Opm::ParallelOverlappingILU0<M,X,Y,C> >
79 {
81 };
82 
89 template<class Matrix, class Domain, class Range, class ParallelInfo>
90 struct 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 
125 namespace Opm
126 {
142 template<class Matrix, class Domain, class Range, class ParallelInfoT>
144  : public Dune::PreconditionerWithUpdate<Domain,Range>
145 {
146  using ParallelInfo = ParallelInfoT;
147 
148 public:
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 
161 protected:
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 
214 public:
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 
335 protected:
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_;
351  Range reorderedD_;
353  Domain reorderedV_;
354 
355  const ParallelInfo* comm_;
357  const field_type w_;
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