diff --git a/CHANGELOG.md b/CHANGELOG.md index 4a5d1c7066fcdd55447d17a9fe00bbe3944154d4..759e0a8e9dbf5c4f220e9c4bdb8d094a32c86249 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,6 +1,8 @@ # Master (will become release 2.8) -- ... +- constructBoundaryDofs: + - Small interface change in the template parameters: drop `blocksize` and replace it by `BitSetVector` + - The method can now handle generic `dune-functions` basis types, as long as we have consistency in the data types ## Deprecations and removals diff --git a/dune/fufem/functiontools/boundarydofs.hh b/dune/fufem/functiontools/boundarydofs.hh index ab0db3bef17fb5492480f93a388d0fd6079e4251..d70b34d4438354327d51dab5bd24009e501d651f 100644 --- a/dune/fufem/functiontools/boundarydofs.hh +++ b/dune/fufem/functiontools/boundarydofs.hh @@ -10,59 +10,56 @@ #include <dune/fufem/boundarypatch.hh> +#include <dune/functions/backends/concepts.hh> +#include <dune/functions/backends/istlvectorbackend.hh> +#include <dune/functions/functionspacebases/subentitydofs.hh> + /** \brief For a given basis and boundary patch, determine all degrees of freedom on the patch. - \param boundaryDofs Bit is set if corresponding dofs belongs to the boundary patch + A bit of boundaryDofs is set if corresponding dof belongs to the boundary patch. + + \param boundaryDofs BitSetVector corresponding to the basis DOFs \param sfinae Dummy parameter to instantiate this method only for dune-functions bases */ -template <class GridView, class Basis, int blocksize> +template <class GridView, class Basis, class BitSetVector> void constructBoundaryDofs(const BoundaryPatch<GridView>& boundaryPatch, const Basis& basis, - Dune::BitSetVector<blocksize>& boundaryDofs, + BitSetVector& boundaryDofs, typename Basis::LocalView* sfinae = nullptr) { // Check consistency of the input static_assert((std::is_same<GridView, typename Basis::GridView>::value), "BoundaryPatch and global basis must be for the same grid view!"); + // Small helper function to wrap vectors using istlVectorBackend + // if they do not already satisfy the VectorBackend interface. + auto toVectorBackend = [&](auto& v) -> decltype(auto) { + if constexpr (Dune::models<Dune::Functions::Concept::VectorBackend<Basis>, decltype(v)>()) { + return v; + } else { + return Dune::Functions::istlVectorBackend(v); + } + }; + ////////////////////////////////////////////////////////// // Init output bitfield ////////////////////////////////////////////////////////// - boundaryDofs.resize(basis.size({})); - boundaryDofs.unsetAll(); + // ensure IstlVectorBackend interface + auto boundaryDofsBackend = toVectorBackend(boundaryDofs); - auto localView = basis.localView(); -#if DUNE_VERSION_LT(DUNE_FUNCTIONS,2,6) - auto localIndexSet = basis.localIndexSet(); -#endif + boundaryDofsBackend.resize(basis); + boundaryDofsBackend = false; - for (auto it = boundaryPatch.begin(); it != boundaryPatch.end(); ++it) + auto localView = basis.localView(); + auto seDOFs = Dune::Functions::subEntityDOFs(basis); + for(const auto& intersection : boundaryPatch) { - localView.bind(it->inside()); -#if DUNE_VERSION_LT(DUNE_FUNCTIONS,2,6) - localIndexSet.bind(localView); -#endif - const auto& localCoefficients = localView.tree().finiteElement().localCoefficients(); - - for (size_t i=0; i<localCoefficients.size(); i++) - { - //////////////////////////////////////////////////// - // Test whether dof is on the boundary face - //////////////////////////////////////////////////// - - unsigned int entity = localCoefficients.localKey(i).subEntity(); - unsigned int codim = localCoefficients.localKey(i).codim(); - - if (it.containsInsideSubentity(entity, codim)) -#if DUNE_VERSION_LT(DUNE_FUNCTIONS,2,6) - boundaryDofs[localIndexSet.index(i)] = true; -#else - boundaryDofs[localView.index(i)] = true; -#endif - } + localView.bind(intersection.inside()); + for(auto localIndex: seDOFs.bind(localView, intersection)) + boundaryDofsBackend[localView.index(localIndex)] = true; } } diff --git a/dune/fufem/test/CMakeLists.txt b/dune/fufem/test/CMakeLists.txt index 593af5df70a0e7604c79992702c60490a9a94a52..455c6c79edc9efbeefeeb2b0747e25278e65c74f 100644 --- a/dune/fufem/test/CMakeLists.txt +++ b/dune/fufem/test/CMakeLists.txt @@ -8,6 +8,7 @@ set(GRID_BASED_TESTS boundarypatchprolongatortest coarsegridfunctionwrappertest composedfunctiontest + constructboundarydofstest dunefunctionsassemblertest dunefunctionsipdgassemblertest functionintegratortest diff --git a/dune/fufem/test/constructboundarydofstest.cc b/dune/fufem/test/constructboundarydofstest.cc new file mode 100644 index 0000000000000000000000000000000000000000..8c3facfa6e7a77afc1f0e974257622daed0957ed --- /dev/null +++ b/dune/fufem/test/constructboundarydofstest.cc @@ -0,0 +1,86 @@ +#include <config.h> + +#include <dune/common/parallel/mpihelper.hh> +#include <dune/common/test/testsuite.hh> +#include <dune/common/bitsetvector.hh> +#include <dune/common/fvector.hh> + +#include <dune/fufem/functiontools/boundarydofs.hh> +#include <dune/fufem/boundarypatch.hh> +#include <dune/fufem/test/common.hh> + +#include <dune/functions/functionspacebases/defaultglobalbasis.hh> +#include <dune/functions/functionspacebases/lagrangebasis.hh> +#include <dune/functions/functionspacebases/powerbasis.hh> + +#include <dune/istl/bvector.hh> + +using namespace Dune; + +template<int order, class GW> +bool checkConstructBoundaryDofs(const GW& gw) +{ + // create active vertices + const auto dim = GW::dimension; + const auto size = gw.size(dim); + BitSetVector<dim> vertices(size); + // set some values + vertices[0].set(); + vertices[size / 2].set(); + vertices[size-1].set(); + + // create BoundaryPatch + BoundaryPatch bp(gw,vertices); + + // create scalar and powerBasis + using namespace Functions::BasisFactory; + auto lagrangeBasis = makeBasis(gw, lagrange<order>()); + auto powerBasis = makeBasis(gw, power<3>(lagrange<order>())); + + // create common bitsetvector + BitSetVector<1> lagrangeDOFs; + BitSetVector<3> powerDOFs; + + // ... and functions-compatible bitVectors + Dune::BlockVector<Dune::FieldVector<char,1>> lagrangeBits; + Dune::BlockVector<Dune::FieldVector<char,3>> powerBits; + auto lagrangeBackend = Dune::Functions::istlVectorBackend(lagrangeBits); + auto powerBackend = Dune::Functions::istlVectorBackend(powerBits); + + + // check whether this compiles + constructBoundaryDofs(bp,lagrangeBasis,lagrangeDOFs); + constructBoundaryDofs(bp,powerBasis,powerDOFs); + constructBoundaryDofs(bp,lagrangeBasis,lagrangeBackend); + constructBoundaryDofs(bp,powerBasis,powerBackend); + + return true; +} + + +struct Suite +{ + template<class GridType> + bool check(const GridType& grid) + { + auto maxLevel = grid.maxLevel(); + auto gridView = grid.levelGridView(maxLevel); + + bool passed = checkConstructBoundaryDofs<0>(gridView); + passed = passed and checkConstructBoundaryDofs<1>(gridView); + passed = passed and checkConstructBoundaryDofs<2>(gridView); + return passed; + } +}; + +int main(int argc, char** argv) +{ + Dune::MPIHelper::instance(argc, argv); + + Suite tests; + + bool passed = checkWithStructuredGrid<2>(tests, 3); + passed = passed and checkWithStructuredGrid<3>(tests, 3); + + return passed ? 0 : 1; +}