Commit 6f24254d authored by Timo Koch's avatar Timo Koch
Browse files

[amg] Rewrite internal structure + deprecate TypeTag version

parent f7efc3ee
......@@ -52,29 +52,21 @@ namespace Dumux {
* \brief A linear solver based on the ISTL AMG preconditioner
* and the ISTL BiCGSTAB solver.
*/
template <class GridView, class AmgTraits>
class ParallelAMGBackend : public LinearSolver
template <class LinearSolverTraits>
class AMGBiCGSTABBackend : public LinearSolver
{
using Grid = typename GridView::Grid;
using LinearOperator = typename AmgTraits::LinearOperator;
using ScalarProduct = typename AmgTraits::ScalarProduct;
using VType = typename AmgTraits::VType;
using Comm = typename AmgTraits::Comm;
using Smoother = typename AmgTraits::Smoother;
using AMGType = Dune::Amg::AMG<typename AmgTraits::LinearOperator, VType, Smoother,Comm>;
using BCRSMat = typename AmgTraits::LinearOperator::matrix_type;
using DofMapper = typename AmgTraits::DofMapper;
public:
/*!
* \brief Construct the backend for the sequential case only
*
* \param paramGroup the parameter group for parameter lookup
*/
ParallelAMGBackend(const std::string& paramGroup = "")
AMGBiCGSTABBackend(const std::string& paramGroup = "")
: LinearSolver(paramGroup)
, isParallel_(Dune::MPIHelper::getCollectiveCommunication().size() > 1)
, firstCall_(true)
{
if (Dune::MPIHelper::getCollectiveCommunication().size() > 1)
if (isParallel_)
DUNE_THROW(Dune::InvalidStateException, "Using sequential constructor for parallel run. Use signature with gridView and dofMapper!");
}
......@@ -85,11 +77,12 @@ public:
* \param dofMapper an index mapper for dof entities
* \param paramGroup the parameter group for parameter lookup
*/
ParallelAMGBackend(const GridView& gridView,
const DofMapper& dofMapper,
AMGBiCGSTABBackend(const typename LinearSolverTraits::GridView& gridView,
const typename LinearSolverTraits::DofMapper& dofMapper,
const std::string& paramGroup = "")
: LinearSolver(paramGroup)
, phelper_(std::make_shared<ParallelISTLHelper<GridView, AmgTraits>>(gridView, dofMapper))
, phelper_(std::make_shared<ParallelISTLHelper<LinearSolverTraits>>(gridView, dofMapper))
, isParallel_(Dune::MPIHelper::getCollectiveCommunication().size() > 1)
, firstCall_(true)
{}
......@@ -103,37 +96,11 @@ public:
template<class Matrix, class Vector>
bool solve(Matrix& A, Vector& x, Vector& b)
{
std::shared_ptr<Comm> comm;
std::shared_ptr<LinearOperator> fop;
std::shared_ptr<ScalarProduct> sp;
#if HAVE_MPI
if constexpr (AmgTraits::isParallel)
prepareLinearAlgebraParallel<AmgTraits>(A, b, comm, fop, sp, *phelper_, firstCall_);
else
prepareLinearAlgebraSequential<AmgTraits>(A, comm, fop, sp);
solveSequentialOrParallel_(A, x, b);
#else
prepareLinearAlgebraSequential<AmgTraits>(A, comm, fop, sp);
solveSequential_(A, x, b);
#endif
using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
using Criterion = Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<BCRSMat, Dune::Amg::FirstDiagonal>>;
//! \todo Check whether the default accumulation mode atOnceAccu is needed.
//! \todo make parameters changeable at runtime from input file / parameter tree
Dune::Amg::Parameters params(15,2000,1.2,1.6,Dune::Amg::atOnceAccu);
params.setDefaultValuesIsotropic(Grid::dimension);
params.setDebugLevel(this->verbosity());
Criterion criterion(params);
SmootherArgs smootherArgs;
smootherArgs.iterations = 1;
smootherArgs.relaxationFactor = 1;
AMGType amg(*fop, criterion, smootherArgs, *comm);
Dune::BiCGSTABSolver<VType> solver(*fop, *sp, amg, this->residReduction(), this->maxIter(),
comm->communicator().rank() == 0 ? this->verbosity() : 0);
solver.apply(x, b, result_);
firstCall_ = false;
return result_.converged;
}
......@@ -143,7 +110,7 @@ public:
*/
std::string name() const
{
return "AMG preconditioned BiCGSTAB solver";
return "AMG-preconditioned BiCGSTAB solver";
}
/*!
......@@ -156,11 +123,97 @@ public:
private:
std::shared_ptr<ParallelISTLHelper<GridView, AmgTraits>> phelper_;
#if HAVE_MPI
template<class Matrix, class Vector>
void solveSequentialOrParallel_(Matrix& A, Vector& x, Vector& b)
{
if constexpr (LinearSolverTraits::canCommunicate)
{
if (isParallel_)
solveParallel_(A, x, b);
else
solveSequential_(A, x, b);
}
else
{
solveSequential_(A, x, b);
}
}
template<class Matrix, class Vector>
void solveParallel_(Matrix& A, Vector& x, Vector& b)
{
using Traits = typename LinearSolverTraits::template Parallel<Matrix, Vector>;
using Comm = typename Traits::Comm;
using LinearOperator = typename Traits::LinearOperator;
using ScalarProduct = typename Traits::ScalarProduct;
std::shared_ptr<Comm> comm;
std::shared_ptr<LinearOperator> linearOperator;
std::shared_ptr<ScalarProduct> scalarProduct;
prepareLinearAlgebraParallel<LinearSolverTraits>(A, b, comm, linearOperator, scalarProduct, *phelper_, firstCall_);
using SeqSmoother = Dune::SeqSSOR<Matrix, Vector, Vector>;
using Smoother = typename Traits::template Preconditioner<SeqSmoother>;
solveWithAmg_<Smoother>(A, x, b, linearOperator, comm, scalarProduct);
}
#endif // HAVE_MPI
template<class Matrix, class Vector>
void solveSequential_(Matrix& A, Vector& x, Vector& b)
{
using Comm = Dune::Amg::SequentialInformation;
using Traits = typename LinearSolverTraits::template Sequential<Matrix, Vector>;
using LinearOperator = typename Traits::LinearOperator;
using ScalarProduct = typename Traits::ScalarProduct;
auto comm = std::make_shared<Comm>();
auto linearOperator = std::make_shared<LinearOperator>(A);
auto scalarProduct = std::make_shared<ScalarProduct>();
using Smoother = Dune::SeqSSOR<Matrix, Vector, Vector>;
solveWithAmg_<Smoother>(A, x, b, linearOperator, comm, scalarProduct);
}
template<class Smoother, class Matrix, class Vector, class LinearOperator, class Comm, class ScalarProduct>
void solveWithAmg_(Matrix& A, Vector& x, Vector& b,
std::shared_ptr<LinearOperator>& linearOperator,
std::shared_ptr<Comm>& comm,
std::shared_ptr<ScalarProduct>& scalarProduct)
{
using SmootherArgs = typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
using Criterion = Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<Matrix, Dune::Amg::FirstDiagonal>>;
//! \todo Check whether the default accumulation mode atOnceAccu is needed.
//! \todo make parameters changeable at runtime from input file / parameter tree
Dune::Amg::Parameters params(15, 2000, 1.2, 1.6, Dune::Amg::atOnceAccu);
params.setDefaultValuesIsotropic(LinearSolverTraits::GridView::dimension);
params.setDebugLevel(this->verbosity());
Criterion criterion(params);
SmootherArgs smootherArgs;
smootherArgs.iterations = 1;
smootherArgs.relaxationFactor = 1;
using Amg = Dune::Amg::AMG<LinearOperator, Vector, Smoother, Comm>;
auto amg = std::make_shared<Amg>(*linearOperator, criterion, smootherArgs, *comm);
Dune::BiCGSTABSolver<Vector> solver(*linearOperator, *scalarProduct, *amg, this->residReduction(), this->maxIter(),
comm->communicator().rank() == 0 ? this->verbosity() : 0);
solver.apply(x, b, result_);
}
std::shared_ptr<ParallelISTLHelper<LinearSolverTraits>> phelper_;
Dune::InverseOperatorResult result_;
bool isParallel_;
bool firstCall_;
};
// deprecation helper -> will be removed after 3.2
template<class GridView, class Traits>
using ParallelAMGBackend
= AMGBiCGSTABBackend<Traits>;
} // namespace Dumux
#include <dumux/common/properties.hh>
......@@ -175,9 +228,9 @@ namespace Dumux {
* \note This is an adaptor using a TypeTag
*/
template<class TypeTag>
using AMGBackend = ParallelAMGBackend<GetPropType<TypeTag, Properties::GridView>, LinearSolverTraits<GetPropType<TypeTag, Properties::JacobianMatrix>,
GetPropType<TypeTag, Properties::SolutionVector>,
GetPropType<TypeTag, Properties::GridGeometry>>>;
using AMGBackend [[deprecated("Use AMGBiCGSTABBackend instead. Will be removed after 3.2!")]]
= ParallelAMGBackend<GetPropType<TypeTag, Properties::GridView>,
LinearSolverTraits<GetPropType<TypeTag, Properties::GridGeometry>>>;
} // namespace Dumux
......
......@@ -33,7 +33,7 @@ namespace Dumux {
//! The type traits required for using the AMG backend
template<class MType, class VType, class GridGeometry>
using AmgTraits [[deprecated("Use LinearSolverTraits<MType, VType, GridGeometry> instead. AmgTraits will be removed after 3.2!")]] = LinearSolverTraits<MType, VType, GridGeometry>;
using AmgTraits [[deprecated("Use LinearSolverTraits<GridGeometry> instead. AmgTraits will be removed after 3.2!")]] = LinearSolverTraits<GridGeometry>;
} // end namespace Dumux
......
......@@ -64,114 +64,108 @@ struct canCommunicate<Dune::UGGrid<dim>, codim>
namespace Dumux {
//! The implementation is specialized for the different discretizations
template<class MType, class VType, class GridGeometry, DiscretizationMethod discMethod> struct LinearSolverTraitsImpl;
template<class GridGeometry, DiscretizationMethod discMethod>
struct LinearSolverTraitsImpl;
//! The type traits required for using the IstlFactoryBackend
template<class MType, class VType, class GridGeometry>
using LinearSolverTraits = LinearSolverTraitsImpl<MType, VType, GridGeometry, GridGeometry::discMethod>;
template<class GridGeometry>
using LinearSolverTraits = LinearSolverTraitsImpl<GridGeometry, GridGeometry::discMethod>;
//! NonoverlappingSolverTraits used by discretization with non-overlapping parallel model
template <class MType, class VType, bool isParallel>
class NonoverlappingSolverTraits
//! sequential solver traits
template<class MType, class VType>
struct SequentialSolverTraits
{
public:
using Comm = Dune::Amg::SequentialInformation;
using LinearOperator = Dune::MatrixAdapter<MType,VType,VType>;
using Matrix = MType;
using Vector = VType;
using LinearOperator = Dune::MatrixAdapter<MType, VType, VType>;
using ScalarProduct = Dune::SeqScalarProduct<VType>;
using Smoother = Dune::SeqSSOR<MType,VType, VType>;
template<class SeqPreconditioner>
using Preconditioner = SeqPreconditioner;
};
#if HAVE_MPI
template <class MType, class VType>
class NonoverlappingSolverTraits<MType, VType, true>
struct NonoverlappingSolverTraits
{
public:
using Matrix = MType;
using Vector = VType;
using Comm = Dune::OwnerOverlapCopyCommunication<Dune::bigunsignedint<96>, int>;
using LinearOperator = Dune::NonoverlappingSchwarzOperator<MType, VType, VType, Comm>;
using ScalarProduct = Dune::NonoverlappingSchwarzScalarProduct<VType, Comm>;
template<class SeqPreconditioner>
using Preconditioner = Dune::NonoverlappingBlockPreconditioner<Comm, SeqPreconditioner>;
};
template <class MType, class VType>
struct OverlappingSolverTraits
{
public:
using Comm = Dune::OwnerOverlapCopyCommunication<Dune::bigunsignedint<96>,int>;
using LinearOperator = Dune::NonoverlappingSchwarzOperator<MType,VType, VType,Comm>;
using ScalarProduct = Dune::NonoverlappingSchwarzScalarProduct<VType,Comm>;
using Smoother = Dune::NonoverlappingBlockPreconditioner<Comm,Dune::SeqSSOR<MType,VType, VType> >;
using Matrix = MType;
using Vector = VType;
using Comm = Dune::OwnerOverlapCopyCommunication<Dune::bigunsignedint<96>, int>;
using LinearOperator = Dune::OverlappingSchwarzOperator<MType, VType, VType, Comm>;
using ScalarProduct = Dune::OverlappingSchwarzScalarProduct<VType, Comm>;
template<class SeqPreconditioner>
using Preconditioner = Dune::BlockPreconditioner<VType, VType, Comm, SeqPreconditioner>;
};
#endif
//! Box: use non-overlapping model
template<class Matrix, class Vector, class GridGeometry>
struct LinearSolverTraitsImpl<Matrix, Vector, GridGeometry, DiscretizationMethod::box>
template<class GridGeometry>
struct LinearSolverTraitsImpl<GridGeometry, DiscretizationMethod::box>
{
using GridView = typename GridGeometry::GridView;
using Grid = typename GridGeometry::GridView::Traits::Grid;
enum {
dofCodim = Grid::dimension,
isNonOverlapping = true,
// TODO: see above for description of this workaround, remove second line if fixed upstream
isParallel = Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v
};
static constexpr bool canCommunicate = Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v;
using MType = Matrix;
using VType = Dune::BlockVector<Dune::FieldVector<typename Vector::block_type::value_type, Vector::block_type::dimension>>;
using SolverTraits = NonoverlappingSolverTraits<MType, VType, isParallel>;
using Comm = typename SolverTraits::Comm;
using LinearOperator = typename SolverTraits::LinearOperator;
using ScalarProduct = typename SolverTraits::ScalarProduct;
using Smoother = typename SolverTraits::Smoother;
using DofMapper = typename GridGeometry::VertexMapper;
};
//! OverlappingSolverTraits used by discretization with overlapping parallel model
template <class MType, class VType, bool isParallel>
class OverlappingSolverTraits
{
public:
using Comm = Dune::Amg::SequentialInformation;
using LinearOperator = Dune::MatrixAdapter<MType,VType,VType>;
using ScalarProduct = Dune::SeqScalarProduct<VType>;
using Smoother = Dune::SeqSSOR<MType,VType, VType>;
};
static constexpr int dofCodim = Grid::dimension;
static constexpr bool isNonOverlapping = true;
// TODO: see above for description of this workaround, remove second line if fixed upstream
static constexpr bool canCommunicate =
Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v;
#if HAVE_MPI
template <class MType, class VType>
class OverlappingSolverTraits<MType, VType, true>
{
public:
using Comm = Dune::OwnerOverlapCopyCommunication<Dune::bigunsignedint<96>,int>;
using LinearOperator = Dune::OverlappingSchwarzOperator<MType,VType, VType,Comm>;
using ScalarProduct = Dune::OverlappingSchwarzScalarProduct<VType,Comm>;
using Smoother = Dune::BlockPreconditioner<VType,VType,Comm,Dune::SeqSSOR<MType,VType, VType> >;
template<class Matrix, class Vector>
using Sequential = SequentialSolverTraits<Matrix, Vector>;
template<class Matrix, class Vector>
using Parallel = std::conditional_t<isNonOverlapping,
NonoverlappingSolverTraits<Matrix, Vector>,
OverlappingSolverTraits<Matrix, Vector>>;
};
#endif
//! Cell-centered tpfa: use overlapping model
template<class Matrix, class Vector, class GridGeometry>
struct LinearSolverTraitsImpl<Matrix, Vector, GridGeometry, DiscretizationMethod::cctpfa>
template<class GridGeometry>
struct LinearSolverTraitsImpl<GridGeometry, DiscretizationMethod::cctpfa>
{
using GridView = typename GridGeometry::GridView;
using Grid = typename GridGeometry::GridView::Traits::Grid;
enum {
dofCodim = 0,
isNonOverlapping = false,
// TODO: see above for description of this workaround, remove second line if fixed upstream
isParallel = Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v
};
static constexpr bool canCommunicate = Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v;
using MType = Matrix;
using VType = Dune::BlockVector<Dune::FieldVector<typename Vector::block_type::value_type, Vector::block_type::dimension>>;
using SolverTraits = OverlappingSolverTraits<MType, VType, isParallel>;
using Comm = typename SolverTraits::Comm;
using LinearOperator = typename SolverTraits::LinearOperator;
using ScalarProduct = typename SolverTraits::ScalarProduct;
using Smoother = typename SolverTraits::Smoother;
using DofMapper = typename GridGeometry::ElementMapper;
using DofMapper = typename GridGeometry::VertexMapper;
static constexpr int dofCodim = 0;
static constexpr bool isNonOverlapping = false;
// TODO: see above for description of this workaround, remove second line if fixed upstream
static constexpr bool canCommunicate =
Dune::Capabilities::canCommunicate<Grid, dofCodim>::v
|| Dumux::Temp::Capabilities::canCommunicate<Grid, dofCodim>::v;
template<class Matrix, class Vector>
using Sequential = SequentialSolverTraits<Matrix, Vector>;
template<class Matrix, class Vector>
using Parallel = std::conditional_t<isNonOverlapping,
NonoverlappingSolverTraits<Matrix, Vector>,
OverlappingSolverTraits<Matrix, Vector>>;
};
template<class Matrix, class Vector, class GridGeometry>
struct LinearSolverTraitsImpl<Matrix, Vector, GridGeometry, DiscretizationMethod::ccmpfa>
: public LinearSolverTraitsImpl<Matrix, Vector, GridGeometry, DiscretizationMethod::cctpfa> {};
//! Cell-centered mpfa: use overlapping model
template<class GridGeometry>
struct LinearSolverTraitsImpl<GridGeometry, DiscretizationMethod::ccmpfa>
: public LinearSolverTraitsImpl<GridGeometry, DiscretizationMethod::cctpfa> {};
} // end namespace Dumux
......
......@@ -22,8 +22,10 @@
* \brief Provides a helper class for nonoverlapping
* decomposition.
*/
#ifndef DUMUX_PARALLELHELPERS_HH
#define DUMUX_PARALLELHELPERS_HH
#ifndef DUMUX_LINEAR_PARALLELHELPERS_HH
#define DUMUX_LINEAR_PARALLELHELPERS_HH
#if HAVE_MPI
#include <dune/common/version.hh>
#include <dune/geometry/dimension.hh>
......@@ -40,9 +42,10 @@ namespace Dumux {
* decomposition of all degrees of freedom
*/
// operator that resets result to zero at constrained DOFS
template<class GridView, class LinearSolverTraits>
template<class LinearSolverTraits>
class ParallelISTLHelper
{
using GridView = typename LinearSolverTraits::GridView;
using DofMapper = typename LinearSolverTraits::DofMapper;
enum { dofCodim = LinearSolverTraits::dofCodim };
......@@ -312,11 +315,6 @@ public:
: gridView_(gridView), mapper_(mapper), initialized_(false)
{}
[[deprecated("The verbose argument has no effect. Use ParallelISTLHelper(gridView, mapper) instead. Will be removed after 3.2!")]]
ParallelISTLHelper(const GridView& gridView, const DofMapper& mapper, int verbose)
: gridView_(gridView), mapper_(mapper), initialized_(false)
{}
// \brief Initializes the markers for ghosts and owners with the correct size and values.
//
void initGhostsAndOwners()
......@@ -344,22 +342,6 @@ public:
initialized_ = true;
}
// keep only DOFs assigned to this processor
template<class W>
[[deprecated("This function has no effect. Will be removed after 3.2!")]]
void mask(W& w) const
{}
// access to mask vector
[[deprecated("Will be removed after 3.2!")]]
std::size_t mask(std::size_t i) const
{ return isOwned_[i]; }
// access to ghost vector
[[deprecated("Will be removed after 3.2!")]]
std::size_t ghost(std::size_t i) const
{ return isGhost_[i]; }
bool isGhost(std::size_t i) const
{ return isGhost_[i] == ghostMarker_; }
......@@ -367,22 +349,12 @@ public:
template<class B, class A>
void makeNonOverlappingConsistent(Dune::BlockVector<B,A>& v)
{
#if HAVE_MPI
ConsistencyBoxGatherScatter<Dune::BlockVector<B,A> > gs(v, mapper_);
if (gridView_.comm().size() > 1)
gridView_.communicate(gs, Dune::InteriorBorder_InteriorBorder_Interface,
Dune::ForwardCommunication);
#endif
}
#if HAVE_MPI
template<typename MatrixType, typename Comm>
[[deprecated("Use createParallelIndexSet(comm) instead. Will be removed after 3.2!")]]
void createIndexSetAndProjectForAMG(MatrixType& m, Comm& c)
{ createParallelIndexSet(c); }
/*!
* \brief Creates a parallel index set
*
......@@ -452,8 +424,6 @@ public:
comm.remoteIndices().template rebuild<false>();
}
#endif
//! Return the dofMapper
const DofMapper& dofMapper() const
{ return mapper_; }
......@@ -462,24 +432,6 @@ public:
const GridView& gridView() const
{ return gridView_; }
template<class Comm>
Dune::OwnerOverlapCopyAttributeSet::AttributeSet
getAttribute_(const Comm& comm, const bool isOwned, const bool isGhost) const
{
if (isOwned)
return Dune::OwnerOverlapCopyAttributeSet::owner;
#if DUNE_VERSION_GTE(DUNE_ISTL, 2, 7)
else if (isGhost && ( comm.category() ==
static_cast<int>(Dune::SolverCategory::nonoverlapping)) )
#else
else if (isGhost && ( comm.getSolverCategory() ==
static_cast<int>(Dune::SolverCategory::nonoverlapping)) )
#endif
return Dune::OwnerOverlapCopyAttributeSet::overlap;
else
return Dune::OwnerOverlapCopyAttributeSet::copy;
}
template<class Comm, class GlobalIndices>
void resizeIndexSet_(Comm& comm, const GlobalIndices& globalIndices) const
{
......@@ -490,16 +442,32 @@ public:
const auto globalIndex = globalIndices[localIndex];
if (globalIndex != std::numeric_limits<typename GlobalIndices::value_type>::max())
{
const bool isOwned = isOwned_[localIndex] > 0;
const auto attr = getAttribute_(comm, isOwned, isGhost(localIndex));
using LocalIndex = typename Comm::ParallelIndexSet::LocalIndex;
comm.indexSet().add(globalIndex, LocalIndex{localIndex, attr});
const bool isOwned = isOwned_[localIndex] > 0;
const auto attr = getAttribute_(comm, isOwned, isGhost(localIndex));
using LocalIndex = typename Comm::ParallelIndexSet::LocalIndex;
comm.indexSet().add(globalIndex, LocalIndex{localIndex, attr});
}
}
comm.indexSet().endResize();
}
template<class Comm>
Dune::OwnerOverlapCopyAttributeSet::AttributeSet
getAttribute_(const Comm& comm, const bool isOwned, const bool isGhost) const
{
if (isOwned)
return Dune::OwnerOverlapCopyAttributeSet::owner;
#if DUNE_VERSION_GTE(DUNE_ISTL, 2, 7)
else if (isGhost && (comm.category() == static_cast<int>(Dune::SolverCategory::nonoverlapping)) )
#else
else if (isGhost && (comm.getSolverCategory() == static_cast<int>(Dune::SolverCategory::nonoverlapping)) )
#endif
return Dune::OwnerOverlapCopyAttributeSet::overlap;
else
return Dune::OwnerOverlapCopyAttributeSet::copy;
}
private:
const GridView gridView_; //!< the grid view
const DofMapper& mapper_; //!< the dof mapper
......@@ -512,20 +480,15 @@ private:
/*!
* \ingroup Linear
* \brief Helper class for adding up matrix entries on border.
* \tparam GridView The grid view to work on
* \tparam LinearSolverTraits traits class
*/
template<class GridView, class LinearSolverTraits>
class EntityExchanger
template<class Matrix, class GridView, class DofMapper, int dofCodim>
class ParallelMatrixHelper
{
using Matrix = typename LinearSolverTraits::MType;
enum { dim = GridView::dimension };
enum { dofCodim = LinearSolverTraits::dofCodim };
static constexpr int dim = GridView::dimension;
using Grid = typename GridView::Traits::Grid;
using BlockType = typename Matrix::block_type;
using IDS = typename Grid::Traits::GlobalIdSet;
using IdType = typename IDS::IdType;
using DofMapper = typename LinearSolverTraits::DofMapper;
/*!
* \brief A DataHandle class to exchange matrix sparsity patterns.
......@@ -545,22 +508,27 @@ class EntityExchanger
* </pre>
* If we look at vertex 7 and the corresponding entries in the matrix for P0,
* there will be entries for (7,4) and (7,8), but not for (7,2).
* The MatrixPatternExchange class will find these entries and returns a vector "sparsity",
* The MatrixPatternExchange class will find these entries and fills a vector "sparsity",
* that contains all missing connections.
<