OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-11 08:51:55 -06:00
Code Issues Packages Projects Releases Wiki Activity

Uses more speaking names for template parameter and adds more documentation.

Browse Source
This commit is contained in:
Markus Blatt 2015-12-01 14:49:07 +01:00
parent e05e3fa478
commit 4adf8487ea
2 changed files with 137 additions and 68 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

63
opm/autodiff/ParallelOverlappingILU0.hpp
View File

@ -37,20 +37,26 @@ namespace Dune
namespace Amg
{
template<class M, class X, class Y, class C>
struct ConstructionTraits<Opm::ParallelOverlappingILU0<M,X,Y,C> >
/// \brief Tells AMG how to construct the Opm::ParallelOverlappingILU0 smoother
/// \tparam Matrix The type of the Matrix.
/// \tparam Domain The type of the Vector representing the domain.
/// \tparam Range The type of the Vector representing the range.
/// \tparam ParallelInfo The type of the parallel information object
/// used, e.g. Dune::OwnerOverlapCommunication
template<class Matrix, class Domain, class Range, class ParallelInfo>
struct ConstructionTraits<Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> >
{
typedef Dune::SeqILU0<M,X,Y> T;
typedef DefaultParallelConstructionArgs<T,C> Arguments;
typedef Dune::SeqILU0<Matrix,Domain,Range> T;
typedef DefaultParallelConstructionArgs<T,ParallelInfo> Arguments;
typedef ConstructionTraits<T> SeqConstructionTraits;
static inline Opm::ParallelOverlappingILU0<M,X,Y,C>* construct(Arguments& args)
static inline Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo>* construct(Arguments& args)
{
return new Opm::ParallelOverlappingILU0<M,X,Y,C>(args.getMatrix(),
return new Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo>(args.getMatrix(),
args.getComm(),
args.getArgs().relaxationFactor);
}
static inline void deconstruct(Opm::ParallelOverlappingILU0<M,X,Y,C>* bp)
static inline void deconstruct(Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo>* bp)
{
delete bp;
}
@ -64,17 +70,33 @@ struct ConstructionTraits<Opm::ParallelOverlappingILU0<M,X,Y,C> >
namespace Opm
{
template<class M, class X, class Y, typename C>
class ParallelOverlappingILU0 : public Dune::Preconditioner<X,Y> {
/// \brief A two-step version of an overlapping Schwarz preconditioner using one step ILU0 as
///
/// This preconditioner differs from a ParallelRestrictedOverlappingSchwarz with
/// Dune:SeqILU0 in the follwing way:
/// During apply we make sure that the current residual is consistent (i.e.
/// each process knows the same value for each index. The we solve
/// Ly= d for y and make y consistent again. Last we solve Ux = y and
/// make sure that x is consistent consistent.
/// In contrast for ParallelRestrictedOverlappingSchwarz we solve (LU)x = d for x
/// without forcing consistency between the two steps.
/// \tparam Matrix The type of the Matrix.
/// \tparam Domain The type of the Vector representing the domain.
/// \tparam Range The type of the Vector representing the range.
/// \tparam ParallelInfo The type of the parallel information object
/// used, e.g. Dune::OwnerOverlapCommunication
template<class Matrix, class Domain, class Range, class ParallelInfo>
class ParallelOverlappingILU0
: public Dune::Preconditioner<Domain,Range> {
public:
//! \brief The matrix type the preconditioner is for.
typedef typename Dune::remove_const<M>::type matrix_type;
typedef typename Dune::remove_const<Matrix>::type matrix_type;
//! \brief The domain type of the preconditioner.
typedef X domain_type;
typedef Domain domain_type;
//! \brief The range type of the preconditioner.
typedef Y range_type;
typedef Range range_type;
//! \brief The field type of the preconditioner.
typedef typename X::field_type field_type;
typedef typename Domain::field_type field_type;
// define the category
enum {
@ -88,7 +110,8 @@ public:
\param A The matrix to operate on.
\param w The relaxation factor.
*/
ParallelOverlappingILU0 (const M& A, const C& comm, field_type w)
ParallelOverlappingILU0 (const Matrix& A, const ParallelInfo& comm,
field_type w)
: ilu_(A), comm_(comm), w_(w)
{
int ilu_setup_successful = 1;
@ -118,7 +141,7 @@ public:
\copydoc Preconditioner::pre(X&,Y&)
*/
virtual void pre (X& x, Y& b)
virtual void pre (Domain& x, Range& b)
{
DUNE_UNUSED_PARAMETER(x);
DUNE_UNUSED_PARAMETER(b);
@ -129,9 +152,9 @@ public:
\copydoc Preconditioner::apply(X&,const Y&)
*/
virtual void apply (X& v, const Y& d)
virtual void apply (Domain& v, const Range& d)
{
Y& md = const_cast<Y&>(d);
Range& md = const_cast<Range&>(d);
comm_.copyOwnerToAll(md,md);
auto endrow=ilu_.end();
for ( auto row = ilu_.begin(); row != endrow; ++row )
@ -163,15 +186,15 @@ public:
\copydoc Preconditioner::post(X&)
*/
virtual void post (X& x)
virtual void post (Range& x)
{
DUNE_UNUSED_PARAMETER(x);
}
private:
//! \brief The ILU0 decomposition of the matrix.
matrix_type ilu_;
const C& comm_;
Matrix ilu_;
const ParallelInfo& comm_;
//! \brief The relaxation factor to use.
field_type w_;

142
opm/autodiff/ParallelRestrictedAdditiveSchwarz.hpp
View File

@ -36,28 +36,59 @@ namespace Dune
namespace Amg
{
template<class X, class Y, class C, class T>
struct ConstructionTraits<Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T> >
/// \brief Tells AMG how to construct the Opm::ParallelOverlappingILU0 smoother
/// \tparam Domain The type of the Vector representing the domain.
/// \tparam Range The type of the Vector representing the range.
/// \tparam ParallelInfo The type of the parallel information object
/// used, e.g. Dune::OwnerOverlapCommunication
/// \tparam SeqPreconditioner The underlying sequential preconditioner to use.
template<class Range, class Domain, class ParallelInfo, class SeqPreconditioner>
struct ConstructionTraits<Opm::ParallelRestrictedOverlappingSchwarz<Range,
Domain,
ParallelInfo,
SeqPreconditioner> >
{
typedef DefaultParallelConstructionArgs<T,C> Arguments;
typedef ConstructionTraits<T> SeqConstructionTraits;
static inline Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T>* construct(Arguments& args)
typedef DefaultParallelConstructionArgs<SeqPreconditioner,ParallelInfo> Arguments;
typedef ConstructionTraits<SeqPreconditioner> SeqConstructionTraits;
/// \brief Construct a parallel restricted overlapping schwarz preconditioner.
static inline Opm::ParallelRestrictedOverlappingSchwarz<Range,
Domain,
ParallelInfo,
SeqPreconditioner>*
construct(Arguments& args)
{
return new Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T>(*SeqConstructionTraits::construct(args),
args.getComm());
return new Opm::ParallelRestrictedOverlappingSchwarz
<Range,Domain,ParallelInfo,SeqPreconditioner>(*SeqConstructionTraits
::construct(args),
args.getComm());
}
static inline void deconstruct(Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T>* bp)
/// \brief Deconstruct and free a parallel restricted overlapping schwarz preconditioner.
static inline void deconstruct(Opm::ParallelRestrictedOverlappingSchwarz
<Range,Domain,ParallelInfo,SeqPreconditioner>* bp)
{
SeqConstructionTraits::deconstruct(static_cast<T*>(&bp->preconditioner));
SeqConstructionTraits
::deconstruct(static_cast<SeqPreconditioner*>(&bp->preconditioner));
delete bp;
}
};
template<class X, class Y, class C, class T>
struct SmootherTraits<Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T> >
/// \brief Tells AMG how to use Opm::ParallelOverlappingILU0 smoother
/// \tparam Domain The type of the Vector representing the domain.
/// \tparam Range The type of the Vector representing the range.
/// \tparam ParallelInfo The type of the parallel information object
/// used, e.g. Dune::OwnerOverlapCommunication
/// \tparam SeqPreconditioner The underlying sequential preconditioner to use.
template<class Range, class Domain, class ParallelInfo, class SeqPreconditioner>
struct SmootherTraits<Opm::ParallelRestrictedOverlappingSchwarz<Range,
Domain,
ParallelInfo,
SeqPreconditioner> >
{
typedef DefaultSmootherArgs<typename T::matrix_type::field_type> Arguments;
typedef DefaultSmootherArgs<typename SeqPreconditioner::matrix_type::field_type> Arguments;
};
@ -67,26 +98,43 @@ struct SmootherTraits<Opm::ParallelRestrictedOverlappingSchwarz<X,Y,C,T> >
namespace Opm{
/**
* @brief Block parallel preconditioner.
*
* This is essentially a wrapper that takes a sequential
* preconditioner. In each step the sequential preconditioner
* is applied and then all owner data points are updated on
* all other processes.
*/
template<class X, class Y, class C, class T=Dune::Preconditioner<X,Y> >
class ParallelRestrictedOverlappingSchwarz : public Dune::Preconditioner<X,Y> {
friend class Dune::Amg::ConstructionTraits<ParallelRestrictedOverlappingSchwarz<X,Y,C,T> >;
/// \brief Block parallel preconditioner.
///
/// This is essentially a wrapper that takes a sequential
/// preconditioner. In each step the sequential preconditioner
/// is applied to the whole subdomain and then all owner data
/// points are updated on all other processes from the processor
/// that knows the complete matrix row for this data point (in dune-istl
/// speak that is the one that owns the data).
///
/// Note that this is different from the usual approach in dune-istl where
/// the application of the sequential preconditioner only takes place on
/// the (owner) partition of the process disregarding any overlap/ghost region.
///
/// For more information see https://www.cs.colorado.edu/~cai/papers/rash.pdf
///
/// \tparam Domain The type of the Vector representing the domain.
/// \tparam Range The type of the Vector representing the range.
/// \tparam ParallelInfo The type of the parallel information object
/// used, e.g. Dune::OwnerOverlapCommunication
/// \tparam SeqPreconditioner The underlying sequential preconditioner to use.
template<class Range, class Domain, class ParallelInfo, class SeqPreconditioner=Dune::Preconditioner<Range,Domain> >
class ParallelRestrictedOverlappingSchwarz
: public Dune::Preconditioner<Range,Domain> {
friend class Dune::Amg
::ConstructionTraits<ParallelRestrictedOverlappingSchwarz<Range,
Domain,
ParallelInfo,
SeqPreconditioner> >;
public:
//! \brief The domain type of the preconditioner.
typedef X domain_type;
typedef Domain domain_type;
//! \brief The range type of the preconditioner.
typedef Y range_type;
typedef Range range_type;
//! \brief The field type of the preconditioner.
typedef typename X::field_type field_type;
typedef typename Domain::field_type field_type;
//! \brief The type of the communication object.
typedef C communication_type;
typedef ParallelInfo communication_type;
// define the category
enum {
@ -101,8 +149,8 @@ public:
\param c The communication object for syncing overlap and copy
data points. (E.~g. OwnerOverlapCommunication )
*/
ParallelRestrictedOverlappingSchwarz (T& p, const communication_type& c)
: preconditioner(p), communication(c)
ParallelRestrictedOverlappingSchwarz (SeqPreconditioner& p, const communication_type& c)
: preconditioner_(p), communication_(c)
{ }
/*!
@ -110,10 +158,10 @@ public:
\copydoc Preconditioner::pre(X&,Y&)
*/
virtual void pre (X& x, Y& b)
virtual void pre (Domain& x, Range& b)
{
communication.copyOwnerToAll(x,x); // make dirichlet values consistent
preconditioner.pre(x,b);
communication_.copyOwnerToAll(x,x); // make dirichlet values consistent
preconditioner_.pre(x,b);
}
/*!
@ -121,23 +169,21 @@ public:
\copydoc Preconditioner::apply(X&,const Y&)
*/
virtual void apply (X& v, const Y& d)
virtual void apply (Domain& v, const Range& d)
{
Y& md = const_cast<Y&>(d);
communication.copyOwnerToAll(md,md);
preconditioner.apply(v,d);
communication.copyOwnerToAll(v,v);
communication.project(md);
apply<true>(v, d);
}
template<bool forward>
void apply (X& v, const Y& d)
void apply (Domain& v, const Range& d)
{
Y& md = const_cast<Y&>(d);
communication.copyOwnerToAll(md,md);
preconditioner.template apply<forward>(v,d);
communication.copyOwnerToAll(v,v);
communication.project(md);
// hack us a mutable d to prevent copying.
Range& md = const_cast<Range&>(d);
communication_.copyOwnerToAll(md,md);
preconditioner_.template apply<forward>(v,d);
communication_.copyOwnerToAll(v,v);
// Make sure that d is the same as at the beginning of apply.
communication_.project(md);
}
/*!
@ -145,17 +191,17 @@ public:
\copydoc Preconditioner::post(X&)
*/
virtual void post (X& x)
virtual void post (Range& x)
{
preconditioner.post(x);
preconditioner_.post(x);
}
private:
//! \brief a sequential preconditioner
T& preconditioner;
SeqPreconditioner& preconditioner_;
//! \brief the communication object
const communication_type& communication;
const communication_type& communication_;
};