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/istl/paamg/smoother.hh>
36template<
class Matrix,
class Domain,
class Range,
class ParallelInfo>
37class ParallelOverlappingILU0;
41 :
public Dune::Amg::DefaultSmootherArgs<F>
76template<
class M,
class X,
class Y,
class C>
77struct SmootherTraits<
Opm::ParallelOverlappingILU0<M,X,Y,C> >
88template<
class Matrix,
class Domain,
class Range,
class ParallelInfo>
89struct ConstructionTraits<
Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> >
92 using Arguments = DefaultParallelConstructionArgs<T,ParallelInfo>;
94 using ParallelOverlappingILU0Pointer = std::shared_ptr<T>;
96 static inline ParallelOverlappingILU0Pointer construct(Arguments& args)
98 return ParallelOverlappingILU0Pointer(
99 new T(args.getMatrix(),
101 args.getArgs().getN(),
102 args.getArgs().relaxationFactor,
103 args.getArgs().getMilu()) );
128template<
class Matrix,
class Domain,
class Range,
class ParallelInfoT>
132 using ParallelInfo = ParallelInfoT;
144 using block_type =
typename matrix_type::block_type;
145 using size_type =
typename matrix_type::size_type;
150 CRS() : nRows_( 0 ) {}
152 size_type rows()
const {
return nRows_; }
154 size_type nonZeros()
const
156 assert( rows_[ rows() ] != size_type(-1) );
157 return rows_[ rows() ];
160 void resize(
const size_type nRows )
162 if( nRows_ != nRows )
165 rows_.resize( nRows_+1, size_type(-1) );
169 void reserveAdditional(
const size_type nonZeros )
171 const size_type needed = values_.size() + nonZeros ;
172 if( values_.capacity() < needed )
174 const size_type estimate = needed * 1.1;
175 values_.reserve( estimate );
176 cols_.reserve( estimate );
180 void push_back(
const block_type& value,
const size_type index )
182 values_.push_back( value );
183 cols_.push_back( index );
194 std::vector<size_type> rows_;
195 std::vector<block_type> values_;
196 std::vector<size_type> cols_;
201 Dune::SolverCategory::Category category()
const override;
219 bool reorder_sphere =
true);
234 const ParallelInfo& comm,
const int n,
const field_type w,
236 bool reorder_sphere =
true);
252 bool redblack =
false,
253 bool reorder_sphere =
true);
271 bool reorder_sphere =
true);
288 const ParallelInfo& comm,
290 size_type interiorSize,
bool redblack =
false,
291 bool reorder_sphere =
true);
298 void pre (Domain&, Range&)
override
306 void apply (Domain& v,
const Range& d)
override;
309 void copyOwnerToAll( V& v )
const;
319 void update()
override;
328 void reorderBack(
const Range& reorderedV, Range& v);
334 std::vector< block_type > inv_;
342 const ParallelInfo* comm_;
345 const bool relaxation_;
346 size_type interiorSize_;
Interface class adding the update() method to the preconditioner interface.
Definition: PreconditionerWithUpdate.hpp:32
Definition: ParallelOverlappingILU0.hpp:42
A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as.
Definition: ParallelOverlappingILU0.hpp:131
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:201
std::unique_ptr< Matrix > ILU_
The ILU0 decomposition of the matrix.
Definition: ParallelOverlappingILU0.hpp:331
Domain & reorderV(Domain &v)
Reorder V if needed and return a reference to it.
Definition: ParallelOverlappingILU0_impl.hpp:555
typename std::remove_const< Matrix >::type matrix_type
The matrix type the preconditioner is for.
Definition: ParallelOverlappingILU0.hpp:136
Domain domain_type
The domain type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:138
std::vector< std::size_t > ordering_
the reordering of the unknowns
Definition: ParallelOverlappingILU0.hpp:336
Domain reorderedV_
The reordered left hand side.
Definition: ParallelOverlappingILU0.hpp:340
void post(Range &) override
Clean up.
Definition: ParallelOverlappingILU0.hpp:316
Range reorderedD_
The reordered right hand side.
Definition: ParallelOverlappingILU0.hpp:338
const field_type w_
The relaxation factor to use.
Definition: ParallelOverlappingILU0.hpp:344
Range & reorderD(const Range &d)
Reorder D if needed and return a reference to it.
Definition: ParallelOverlappingILU0_impl.hpp:531
Range range_type
The range type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:140
void apply(Domain &v, const Range &d) override
Apply the preconditoner.
Definition: ParallelOverlappingILU0_impl.hpp:290
void pre(Domain &, Range &) override
Prepare the preconditioner.
Definition: ParallelOverlappingILU0.hpp:298
typename Domain::field_type field_type
The field type of the preconditioner.
Definition: ParallelOverlappingILU0.hpp:142
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:149