Skip to content
Snippets Groups Projects
Commit 13513a25 authored by Jonathan Youett's avatar Jonathan Youett
Browse files

Change interface such that it can inherit from the multigridtransfer base class:

The critical bitfield, determining which dofs are truncated, is now a member of the class, so the standard transferoperator methods can be used
parent d06940f7
Branches
No related tags found
No related merge requests found
...@@ -22,7 +22,7 @@ template< ...@@ -22,7 +22,7 @@ template<
class BitVectorType = Dune::BitSetVector<VectorType::block_type::dimension>, class BitVectorType = Dune::BitSetVector<VectorType::block_type::dimension>,
class MatrixType = Dune::BCRSMatrix< typename Dune::FieldMatrix< class MatrixType = Dune::BCRSMatrix< typename Dune::FieldMatrix<
typename VectorType::field_type, VectorType::block_type::dimension, VectorType::block_type::dimension> > > typename VectorType::field_type, VectorType::block_type::dimension, VectorType::block_type::dimension> > >
class TruncatedMGTransfer class TruncatedMGTransfer : virtual public MultigridTransfer<VectorType, BitVectorType, MatrixType>
{ {
enum {blocksize = VectorType::block_type::dimension}; enum {blocksize = VectorType::block_type::dimension};
...@@ -32,31 +32,19 @@ class TruncatedMGTransfer ...@@ -32,31 +32,19 @@ class TruncatedMGTransfer
public: public:
/** \brief Default constructor */ /** \brief Default constructor */
TruncatedMGTransfer() : recompute_(NULL) TruncatedMGTransfer() : recompute_(nullptr), critical_(nullptr) {}
{}
/** \brief Restrict level fL of f and store the result in level cL of t /** \brief Set recompute bitfield. */
* void setRecomputeBitField(const Dune::BitSetVector<1>* recompute)
* \param critical Has to contain an entry for each degree of freedom. {
* Those dofs with a set bit are treated as critical. recompute_ = recompute;
*/ }
virtual void restrict(const VectorType& f, VectorType &t, const BitVectorType& critical) const = 0;
/** \brief Prolong level cL of f and store the result in level fL of t
*
* \param critical Has to contain an entry for each degree of freedom.
* Those dofs with a set bit are treated as critical.
*/
virtual void prolong(const VectorType& f, VectorType &t, const BitVectorType& critical) const = 0;
/** \brief Galerkin assemble a coarse stiffness matrix /** \brief Set critical bitfield. */
* void setCriticalBitField(const BitVectorType* critical)
* \param critical Has to contain an entry for each degree of freedom. {
* Those dofs with a set bit are treated as critical. critical_ = critical;
*/ }
virtual void galerkinRestrict(const MatrixType& fineMat,
MatrixType& coarseMat,
const BitVectorType& critical) const = 0;
/** \brief Bitfield specifying a subset of dofs which need to be recomputed /** \brief Bitfield specifying a subset of dofs which need to be recomputed
* when doing Galerkin restriction * when doing Galerkin restriction
...@@ -68,6 +56,12 @@ public: ...@@ -68,6 +56,12 @@ public:
* on the problem this can lead to considerable time savings. * on the problem this can lead to considerable time savings.
*/ */
const Dune::BitSetVector<1>* recompute_; const Dune::BitSetVector<1>* recompute_;
/**
* \brief Has to contain an entry for each degree of freedom.
* Those dofs with a set bit are treated as critical.
*/
const BitVectorType* critical_;
}; };
#endif #endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment