20#ifndef OPM_PARALLELOVERLAPPINGILU0_HEADER_INCLUDED
21#define OPM_PARALLELOVERLAPPINGILU0_HEADER_INCLUDED
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>
37template<
class Matrix,
class Domain,
class Range,
class ParallelInfo>
38class ParallelOverlappingILU0;
42 :
public Dune::Amg::DefaultSmootherArgs<F>
77template<
class M,
class X,
class Y,
class C>
78struct SmootherTraits<
Opm::ParallelOverlappingILU0<M,X,Y,C> >
89template<
class Matrix,
class Domain,
class Range,
class ParallelInfo>
90struct ConstructionTraits<
Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> >
93 using Arguments = DefaultParallelConstructionArgs<T,ParallelInfo>;
95#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
104 new T(args.getMatrix(),
106 args.getArgs().getN(),
107 args.getArgs().relaxationFactor,
108 args.getArgs().getMilu()) );
111#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
113 static inline void deconstruct(
T* bp)
142template<
class Matrix,
class Domain,
class Range,
class ParallelInfoT>
146 using ParallelInfo = ParallelInfoT;
158 using block_type =
typename matrix_type::block_type;
159 using size_type =
typename matrix_type::size_type;
164 CRS() : nRows_( 0 ) {}
166 size_type rows()
const {
return nRows_; }
168 size_type nonZeros()
const
170 assert( rows_[ rows() ] != size_type(-1) );
171 return rows_[ rows() ];
174 void resize(
const size_type nRows )
176 if( nRows_ != nRows )
179 rows_.resize( nRows_+1, size_type(-1) );
183 void reserveAdditional(
const size_type nonZeros )
185 const size_type needed = values_.size() + nonZeros ;
186 if( values_.capacity() < needed )
188 const size_type estimate = needed * 1.1;
189 values_.reserve( estimate );
190 cols_.reserve( estimate );
194 void push_back(
const block_type& value,
const size_type index )
196 values_.push_back( value );
197 cols_.push_back( index );
208 std::vector<size_type> rows_;
209 std::vector<block_type> values_;
210 std::vector<size_type> cols_;
215 Dune::SolverCategory::Category category()
const override;
233 bool reorder_sphere =
true);
248 const ParallelInfo& comm,
const int n,
const field_type w,
250 bool reorder_sphere =
true);
266 bool redblack =
false,
267 bool reorder_sphere =
true);
285 bool reorder_sphere =
true);
302 const ParallelInfo& comm,
304 size_type interiorSize,
bool redblack =
false,
305 bool reorder_sphere =
true);
312 void pre (Domain&, Range&)
override
320 void apply (Domain& v,
const Range& d)
override;
323 void copyOwnerToAll( V& v )
const;
333 void update()
override;
342 void reorderBack(
const Range& reorderedV, Range& v);
347 std::vector< block_type > inv_;
355 const ParallelInfo* comm_;
358 const bool relaxation_;
359 size_type interiorSize_;
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