Commit bd9cc7c0 authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

Matrix-Vector facade class

parent 6804ddd9
...@@ -10,6 +10,7 @@ HIDE_SCOPE_NAMES = YES ...@@ -10,6 +10,7 @@ HIDE_SCOPE_NAMES = YES
HIDE_UNDOC_CLASSES = NO HIDE_UNDOC_CLASSES = NO
INTERNAL_DOCS = NO INTERNAL_DOCS = NO
MARKDOWN_SUPPORT = YES MARKDOWN_SUPPORT = YES
GENERATE_XML = YES
EXCLUDE_SYMBOLS = AMDiS::Impl \ EXCLUDE_SYMBOLS = AMDiS::Impl \
AMDiS::Math::Impl_ \ AMDiS::Math::Impl_ \
...@@ -30,6 +31,7 @@ PREDEFINED += HAVE_UMFPACK \ ...@@ -30,6 +31,7 @@ PREDEFINED += HAVE_UMFPACK \
INPUT += @top_srcdir@/src/amdis \ INPUT += @top_srcdir@/src/amdis \
@top_srcdir@/src/amdis/common \ @top_srcdir@/src/amdis/common \
@top_srcdir@/src/amdis/functions \
@top_srcdir@/src/amdis/gridfunctions \ @top_srcdir@/src/amdis/gridfunctions \
@top_srcdir@/src/amdis/linearalgebra \ @top_srcdir@/src/amdis/linearalgebra \
@top_srcdir@/src/amdis/linearalgebra/mtl \ @top_srcdir@/src/amdis/linearalgebra/mtl \
...@@ -47,7 +49,8 @@ INPUT += @top_srcdir@/src/amdis \ ...@@ -47,7 +49,8 @@ INPUT += @top_srcdir@/src/amdis \
# subdirectory from a directory tree whose root is specified with the INPUT tag. # subdirectory from a directory tree whose root is specified with the INPUT tag.
EXCLUDE += @top_srcdir@/src/amdis/linearalgebra/eigen \ EXCLUDE += @top_srcdir@/src/amdis/linearalgebra/eigen \
@top_srcdir@/src/amdis/linearalgebra/istl @top_srcdir@/src/amdis/linearalgebra/istl \
@top_srcdir@/src/amdis/linearalgebra/petsc
# The EXAMPLE_PATH tag can be used to specify one or more files or # The EXAMPLE_PATH tag can be used to specify one or more files or
# directories that contain example code fragments that are included (see # directories that contain example code fragments that are included (see
......
...@@ -20,14 +20,12 @@ namespace AMDiS ...@@ -20,14 +20,12 @@ namespace AMDiS
* \tparam CB Basis of matrix columns * \tparam CB Basis of matrix columns
* \tparam T Coefficient type to store in the matrix * \tparam T Coefficient type to store in the matrix
**/ **/
template <class RB, class CB, class T = double> template <class RB, class CB, class T = double, class Traits = BackendTraits<RB,T>>
class BiLinearForm class BiLinearForm
: public MatrixBase<RB,CB,MatrixBackend<BackendTraits<RB,T>>> : public MatrixFacade<T, typename Traits::SparsityPattern, Traits::template MatrixImpl>
{ {
using Self = BiLinearForm; using Self = BiLinearForm;
using Traits = BackendTraits<RB,T>; using Super = MatrixFacade<T, typename Traits::SparsityPattern, Traits::template MatrixImpl>;
using Backend = MatrixBackend<Traits>;
using Super = MatrixBase<RB,CB,Backend>;
public: public:
/// The type of the finite element space / basis of the row /// The type of the finite element space / basis of the row
...@@ -47,23 +45,21 @@ namespace AMDiS ...@@ -47,23 +45,21 @@ namespace AMDiS
public: public:
/// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into /// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into
/// a new shared_ptr. /// a new shared_ptr.
template <class RB_, class CB_> BiLinearForm(std::shared_ptr<RB> const& rowBasis, std::shared_ptr<CB> const& colBasis)
BiLinearForm(RB_&& rowBasis, CB_&& colBasis) : Super(*rowBasis, *colBasis)
: Super(FWD(rowBasis), FWD(colBasis)) , rowBasis_(rowBasis)
, colBasis_(colBasis)
{ {
operators_.init(*Super::rowBasis(), *Super::colBasis()); operators_.init(*rowBasis, *colBasis);
}
/// Prepare the matrix for insertion. Clears all the entries.
void init(SymmetryStructure symmetry = SymmetryStructure::unknown)
{
Super::init(symmetry);
auto const rowSize = this->rowBasis()->localView().maxSize(); auto const rowSize = rowBasis->localView().maxSize();
auto const colSize = this->colBasis()->localView().maxSize(); auto const colSize = colBasis->localView().maxSize();
elementMatrix_.resize(rowSize, colSize); elementMatrix_.resize(rowSize, colSize);
} }
std::shared_ptr<RowBasis const> const& rowBasis() const { return rowBasis_; }
std::shared_ptr<ColBasis const> const& colBasis() const { return colBasis_; }
/// \brief Associate a local operator with this BiLinearForm /// \brief Associate a local operator with this BiLinearForm
/** /**
* Stores an operator in a list that gets assembled during a call to \ref assemble(). * Stores an operator in a list that gets assembled during a call to \ref assemble().
...@@ -118,6 +114,9 @@ namespace AMDiS ...@@ -118,6 +114,9 @@ namespace AMDiS
/// List of operators associated to row/col node /// List of operators associated to row/col node
MatrixOperators<RowBasis,ColBasis,ElementMatrix> operators_; MatrixOperators<RowBasis,ColBasis,ElementMatrix> operators_;
std::shared_ptr<RowBasis const> rowBasis_;
std::shared_ptr<ColBasis const> colBasis_;
}; };
......
...@@ -9,9 +9,9 @@ ...@@ -9,9 +9,9 @@
namespace AMDiS { namespace AMDiS {
template <class RB, class CB, class T> template <class RB, class CB, class T, class Traits>
template <class ContextTag, class Expr, class RowTreePath, class ColTreePath> template <class ContextTag, class Expr, class RowTreePath, class ColTreePath>
void BiLinearForm<RB,CB,T>:: void BiLinearForm<RB,CB,T,Traits>::
addOperator(ContextTag contextTag, Expr const& expr, addOperator(ContextTag contextTag, Expr const& expr,
RowTreePath row, ColTreePath col) RowTreePath row, ColTreePath col)
{ {
...@@ -24,16 +24,16 @@ addOperator(ContextTag contextTag, Expr const& expr, ...@@ -24,16 +24,16 @@ addOperator(ContextTag contextTag, Expr const& expr,
auto j = child(this->colBasis()->localView().tree(), makeTreePath(col)); auto j = child(this->colBasis()->localView().tree(), makeTreePath(col));
using LocalContext = typename ContextTag::type; using LocalContext = typename ContextTag::type;
using Traits = DefaultAssemblerTraits<LocalContext, ElementMatrix>; using Tr = DefaultAssemblerTraits<LocalContext, ElementMatrix>;
auto op = makeLocalOperator<LocalContext>(expr, this->rowBasis()->gridView()); auto op = makeLocalOperator<LocalContext>(expr, this->rowBasis()->gridView());
auto localAssembler = makeUniquePtr(makeAssembler<Traits>(std::move(op), i, j)); auto localAssembler = makeUniquePtr(makeAssembler<Tr>(std::move(op), i, j));
operators_[i][j].push(contextTag, std::move(localAssembler)); operators_[i][j].push(contextTag, std::move(localAssembler));
} }
template <class RB, class CB, class T> template <class RB, class CB, class T, class Traits>
void BiLinearForm<RB,CB,T>:: void BiLinearForm<RB,CB,T,Traits>::
assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView) assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView)
{ {
elementMatrix_.resize(rowLocalView.size(), colLocalView.size()); elementMatrix_.resize(rowLocalView.size(), colLocalView.size());
...@@ -58,8 +58,8 @@ assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView) ...@@ -58,8 +58,8 @@ assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView)
} }
template <class RB, class CB, class T> template <class RB, class CB, class T, class Traits>
void BiLinearForm<RB,CB,T>:: void BiLinearForm<RB,CB,T,Traits>::
assemble(SymmetryStructure symmetry) assemble(SymmetryStructure symmetry)
{ {
auto rowLocalView = this->rowBasis()->localView(); auto rowLocalView = this->rowBasis()->localView();
......
...@@ -33,38 +33,40 @@ namespace AMDiS ...@@ -33,38 +33,40 @@ namespace AMDiS
/** /**
* \tparam GB Basis of the vector * \tparam GB Basis of the vector
* \tparam T Type of the coefficients * \tparam T Type of the coefficients
* \tparam TraitsType Collection of parameter for the linear-algebra backend
**/ **/
template <class GB, class T = double> template <class GB, class T = double, class TraitsType = BackendTraits<GB,T>>
class DOFVector class DOFVector
: public VectorBase<GB, VectorBackend<BackendTraits<GB,T>>> : public VectorFacade<T, TraitsType::template VectorImpl>
, private Observer<event::preAdapt> , private Observer<event::preAdapt>
, private Observer<event::adapt> , private Observer<event::adapt>
, private Observer<event::postAdapt> , private Observer<event::postAdapt>
{ {
using Self = DOFVector; using Self = DOFVector;
using Super = VectorBase<GB, VectorBackend<BackendTraits<GB,T>>>; using Coefficients = VectorFacade<T, TraitsType::template VectorImpl>;
public: public:
using Backend = VectorBackend<BackendTraits<GB,T>>;
using GlobalBasis = GB; using GlobalBasis = GB;
/// The index/size - type /// The index/size - type
using size_type = typename GlobalBasis::size_type; using size_type = typename GlobalBasis::size_type;
/// The type of the elements of the DOFVector /// The type of the elements of the DOFVector
using value_type = typename Backend::value_type; using value_type = T;
/// Wrapper for the implementation of the transfer of data during grid adaption /// Wrapper for the implementation of the transfer of data during grid adaption
using DataTransfer = DataTransferWrapper<Self>; using DataTransfer = DataTransferWrapper<Self>;
/// Collection of parameter for the linear-algebra backend
using Traits = TraitsType;
public: public:
/// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer. /// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
DOFVector(std::shared_ptr<GB> const& basis, DOFVector(std::shared_ptr<GB> const& basis,
DataTransferOperation op = DataTransferOperation::INTERPOLATE) DataTransferOperation op = DataTransferOperation::INTERPOLATE)
: Super(basis) : Coefficients(*basis)
, Observer<event::preAdapt>(basis->gridView().grid()) , Observer<event::preAdapt>(basis->gridView().grid())
, Observer<event::adapt>(*basis) , Observer<event::adapt>(basis)
, Observer<event::postAdapt>(basis->gridView().grid()) , Observer<event::postAdapt>(basis->gridView().grid())
, dataTransfer_(op, basis) , dataTransfer_(op, basis)
, basis_(basis) , basis_(basis)
...@@ -83,19 +85,20 @@ namespace AMDiS ...@@ -83,19 +85,20 @@ namespace AMDiS
: DOFVector(Dune::wrap_or_move(FWD(basis)), op) : DOFVector(Dune::wrap_or_move(FWD(basis)), op)
{} {}
std::shared_ptr<GlobalBasis const> const& basis() const { return basis_; }
template <class TreePath = RootTreePath> template <class TreePath = RootTreePath>
auto child(TreePath const& path = {}) auto child(TreePath const& path = {})
{ {
auto&& tp = makeTreePath(path); auto&& tp = makeTreePath(path);
return makeDiscreteFunction(*this, tp); return makeDiscreteFunction(coefficients(), *basis_, tp);
} }
template <class TreePath = RootTreePath> template <class TreePath = RootTreePath>
auto child(TreePath const& path = {}) const auto child(TreePath const& path = {}) const
{ {
auto&& tp = makeTreePath(path); auto&& tp = makeTreePath(path);
return makeDiscreteFunction(*this, tp); return makeDiscreteFunction(coefficients(), *basis_, tp);
} }
...@@ -103,7 +106,7 @@ namespace AMDiS ...@@ -103,7 +106,7 @@ namespace AMDiS
/// reference to this DOFVector in the expression. /// reference to this DOFVector in the expression.
/// See \ref DiscreteFunction::interpolate_noalias /// See \ref DiscreteFunction::interpolate_noalias
template <class Expr, class Tag = tag::average> template <class Expr, class Tag = tag::average>
void interpolate_noalias(Expr&& expr, Tag strategy) void interpolate_noalias(Expr&& expr, Tag strategy = {})
{ {
child().interpolate_noalias(FWD(expr), strategy); child().interpolate_noalias(FWD(expr), strategy);
} }
...@@ -111,7 +114,7 @@ namespace AMDiS ...@@ -111,7 +114,7 @@ namespace AMDiS
/// Interpolation of GridFunction to DOFVector. /// Interpolation of GridFunction to DOFVector.
/// See \ref DiscreteFunction::interpolate /// See \ref DiscreteFunction::interpolate
template <class Expr, class Tag = tag::average> template <class Expr, class Tag = tag::average>
void interpolate(Expr&& expr, Tag strategy) void interpolate(Expr&& expr, Tag strategy = {})
{ {
child().interpolate(FWD(expr), strategy); child().interpolate(FWD(expr), strategy);
} }
...@@ -126,6 +129,16 @@ namespace AMDiS ...@@ -126,6 +129,16 @@ namespace AMDiS
} }
Coefficients const& coefficients() const
{
return static_cast<Coefficients const&>(*this);
}
Coefficients& coefficients()
{
return static_cast<Coefficients&>(*this);
}
/// Write DOFVector to file /// Write DOFVector to file
void backup(std::string const& filename); void backup(std::string const& filename);
...@@ -192,7 +205,7 @@ namespace AMDiS ...@@ -192,7 +205,7 @@ namespace AMDiS
/// to \ref DataTransferOperation::INTERPOLATE. /// to \ref DataTransferOperation::INTERPOLATE.
DataTransfer dataTransfer_; DataTransfer dataTransfer_;
std::shared_ptr<GB> basis_; std::shared_ptr<GlobalBasis const> basis_;
}; };
#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION #if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
...@@ -218,6 +231,20 @@ namespace AMDiS ...@@ -218,6 +231,20 @@ namespace AMDiS
return {FWD(basis), op}; return {FWD(basis), op};
} }
/// A Generator for a mutable \ref DiscreteFunction
template <class GB, class T, class Path = RootTreePath>
auto makeDiscreteFunction(DOFVector<GB,T>& dofVec, Path const& path = {})
{
return dofVec.child(path);
}
/// A Generator for a mutable \ref DiscreteFunction
template <class GB, class T, class Path = RootTreePath>
auto makeDiscreteFunction(DOFVector<GB,T> const& dofVec, Path const& path = {})
{
return dofVec.child(path);
}
} // end namespace AMDiS } // end namespace AMDiS
#include <amdis/DOFVector.inc.hpp> #include <amdis/DOFVector.inc.hpp>
...@@ -10,8 +10,8 @@ ...@@ -10,8 +10,8 @@
namespace AMDiS { namespace AMDiS {
template <class GB, class B> template <class GB, class T, class Traits>
void DOFVector<GB,B>:: void DOFVector<GB,T,Traits>::
backup(std::string const& filename) backup(std::string const& filename)
{ {
std::ofstream out(filename, std::ios::binary); std::ofstream out(filename, std::ios::binary);
...@@ -35,8 +35,8 @@ backup(std::string const& filename) ...@@ -35,8 +35,8 @@ backup(std::string const& filename)
} }
template <class GB, class B> template <class GB, class T, class Traits>
void DOFVector<GB,B>:: void DOFVector<GB,T,Traits>::
restore(std::string const& filename) restore(std::string const& filename)
{ {
std::ifstream in(filename, std::ios::binary); std::ifstream in(filename, std::ios::binary);
......
...@@ -82,7 +82,7 @@ preAdapt(C const& coeff, bool mightCoarsen) ...@@ -82,7 +82,7 @@ preAdapt(C const& coeff, bool mightCoarsen)
auto maxLevel = gv.grid().maxLevel(); auto maxLevel = gv.grid().maxLevel();
using std::sqrt; using std::sqrt;
typename Grid::ctype const checkInsideTolerance = sqrt(std::numeric_limits<typename Grid::ctype>::epsilon()); typename Grid::ctype const checkInsideTolerance = sqrt(std::numeric_limits<typename Grid::ctype>::epsilon());
for (auto const& e : elements(gv, typename C::Backend::Traits::PartitionSet{})) for (auto const& e : elements(gv, typename C::Traits::PartitionSet{}))
{ {
auto father = e; auto father = e;
while (father.mightVanish() && father.hasFather()) while (father.mightVanish() && father.hasFather())
...@@ -162,7 +162,7 @@ void DataTransfer<C,B>::adapt(C& coeff) ...@@ -162,7 +162,7 @@ void DataTransfer<C,B>::adapt(C& coeff)
mapper_.update(); mapper_.update();
std::vector<bool> finished(mapper_.size(), false); std::vector<bool> finished(mapper_.size(), false);
for (const auto& e : elements(gv, typename C::Backend::Traits::PartitionSet{})) for (const auto& e : elements(gv, typename C::Traits::PartitionSet{}))
{ {
auto index = mapper_.index(e); auto index = mapper_.index(e);
if (finished[index]) if (finished[index])
......
...@@ -37,7 +37,7 @@ ...@@ -37,7 +37,7 @@
#endif #endif
#include <amdis/linearalgebra/MatrixBase.hpp> #include <amdis/linearalgebra/MatrixFacade.hpp>
#include <amdis/linearalgebra/VectorBase.hpp> #include <amdis/linearalgebra/VectorFacade.hpp>
#include <amdis/linearalgebra/LinearSolver.hpp> #include <amdis/linearalgebra/LinearSolver.hpp>
#include <amdis/linearalgebra/SolverInfo.hpp> #include <amdis/linearalgebra/SolverInfo.hpp>
...@@ -17,24 +17,20 @@ namespace AMDiS ...@@ -17,24 +17,20 @@ namespace AMDiS
* *
* \tparam GB Basis of the vector * \tparam GB Basis of the vector
* \tparam T Coefficient type to store in the vector * \tparam T Coefficient type to store in the vector
* \tparam Traits Collection of parameter for the linear-algebra backend
**/ **/
template <class GB, class T = double> template <class GB, class T = double, class Traits = BackendTraits<GB,T>>
class LinearForm class LinearForm
: public VectorBase<GB, VectorBackend<BackendTraits<GB,T>>> : public VectorFacade<T, Traits::template VectorImpl>
{ {
using Super = VectorBase<GB, VectorBackend<BackendTraits<GB,T>>>; using Super = VectorFacade<T, Traits::template VectorImpl>;
using Self = LinearForm; using Self = LinearForm;
public: public:
using Backend = VectorBackend<BackendTraits<GB,T>>;
/// The type of the functionspace basis associated to this linearform /// The type of the functionspace basis associated to this linearform
using GlobalBasis = GB; using GlobalBasis = GB;
using LocalView = typename GlobalBasis::LocalView; using LocalView = typename GlobalBasis::LocalView;
/// A traits class collecting several parameters of the linear algebra backend
using Traits = typename Backend::Traits;
/// The type of the elements of the DOFVector /// The type of the elements of the DOFVector
using CoefficientType = typename Traits::CoefficientType; using CoefficientType = typename Traits::CoefficientType;
...@@ -43,24 +39,18 @@ namespace AMDiS ...@@ -43,24 +39,18 @@ namespace AMDiS
public: public:
/// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer. /// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
template <class GB_, explicit LinearForm(std::shared_ptr<GB> const& basis)
Dune::disableCopyMove<Self, GB_> = 0> : Super(*basis)
explicit LinearForm(GB_&& basis) , basis_(basis)
: Super(FWD(basis))
{ {
operators_.init(*Super::basis()); operators_.init(*basis);
}
/// Prepare the LinearForm for insertion of values, finish the insertion with auto const localSize = basis->localView().maxSize();
/// \ref finish().
void init(bool clear)
{
Super::init(clear);
auto const localSize = this->basis()->localView().maxSize();
elementVector_.resize(localSize); elementVector_.resize(localSize);
} }
std::shared_ptr<GlobalBasis const> const& basis() const { return basis_; }
/// \brief Associate a local operator with this LinearForm /// \brief Associate a local operator with this LinearForm
/** /**
* Stores an operator in a list that gets assembled during a call to \ref assemble(). * Stores an operator in a list that gets assembled during a call to \ref assemble().
...@@ -112,6 +102,8 @@ namespace AMDiS ...@@ -112,6 +102,8 @@ namespace AMDiS
/// List of operators associated to nodes, filled in \ref addOperator(). /// List of operators associated to nodes, filled in \ref addOperator().
VectorOperators<GlobalBasis,ElementVector> operators_; VectorOperators<GlobalBasis,ElementVector> operators_;
std::shared_ptr<GlobalBasis const> basis_;
}; };
......
...@@ -10,9 +10,9 @@ ...@@ -10,9 +10,9 @@
namespace AMDiS { namespace AMDiS {
template <class GB, class B> template <class GB, class T, class Traits>
template <class ContextTag, class Expr, class TreePath> template <class ContextTag, class Expr, class TreePath>
void LinearForm<GB,B>:: void LinearForm<GB,T,Traits>::
addOperator(ContextTag contextTag, Expr const& expr, TreePath path) addOperator(ContextTag contextTag, Expr const& expr, TreePath path)
{ {
static_assert( Concepts::PreTreePath<TreePath>, static_assert( Concepts::PreTreePath<TreePath>,
...@@ -21,16 +21,16 @@ addOperator(ContextTag contextTag, Expr const& expr, TreePath path) ...@@ -21,16 +21,16 @@ addOperator(ContextTag contextTag, Expr const& expr, TreePath path)
auto i = child(this->basis()->localView().tree(), makeTreePath(path)); auto i = child(this->basis()->localView().tree(), makeTreePath(path));
using LocalContext = typename ContextTag::type; using LocalContext = typename ContextTag::type;
using Traits = DefaultAssemblerTraits<LocalContext, ElementVector>; using Tr = DefaultAssemblerTraits<LocalContext, ElementVector>;
auto op = makeLocalOperator<LocalContext>(expr, this->basis()->gridView()); auto op = makeLocalOperator<LocalContext>(expr, this->basis()->gridView());
auto localAssembler = makeUniquePtr(makeAssembler<Traits>(std::move(op), i)); auto localAssembler = makeUniquePtr(makeAssembler<Tr>(std::move(op), i));
operators_[i].push(contextTag, std::move(localAssembler)); operators_[i].push(contextTag, std::move(localAssembler));
} }
template <class GB, class B> template <class GB, class T, class Traits>
void LinearForm<GB,B>:: void LinearForm<GB,T,Traits>::
assemble(typename GB::LocalView const& localView) assemble(typename GB::LocalView const& localView)
{ {
elementVector_.resize(localView.size()); elementVector_.resize(localView.size());
...@@ -53,8 +53,8 @@ assemble(typename GB::LocalView const& localView) ...@@ -53,8 +53,8 @@ assemble(typename GB::LocalView const& localView)
} }
template <class GB, class B> template <class GB, class T, class Traits>
void LinearForm<GB,B>:: void LinearForm<GB,T,Traits>::
assemble() assemble()
{ {
auto localView = this->basis()->localView(); auto localView = this->basis()->localView();
......