2 #ifndef DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH
3 #define DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH
7 #include <dune/grid/common/gridenums.hh>
8 #include <dune/geometry/referenceelements.hh>
9 #include <dune/localfunctions/common/interfaceswitch.hh>
10 #include <dune/localfunctions/common/localkey.hh>
12 #include <dune/typetree/typetree.hh>
24 :
public TypeTree::LeafNode
26 std::vector<bool> interior;
40 template<
typename P,
typename EG,
typename LFS,
typename T>
41 void volume (
const P& param,
const EG& eg,
const LFS& lfs, T& trafo)
const
43 typedef typename EG::Entity Entity;
44 enum {
dim = Entity::dimension };
47 typename T::RowType empty;
48 typedef typename LFS::Traits::SizeType size_type;
49 typedef FiniteElementInterfaceSwitch<
50 typename LFS::Traits::FiniteElementType
52 for (size_type i=0; i<lfs.size(); i++){
53 const LocalKey& key = FESwitch::coefficients(lfs.finiteElement()).localKey(i);
54 assert(key.codim() ==
dim &&
"InteriorNodeConstraints only work for vertex DOFs");
55 assert(key.index() == 0 &&
"InteriorNodeConstraints only work for P1 shape functions");
58 unsigned int local_idx = key.subEntity();
61 unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx,
dim);
65 trafo[lfs.dofIndex(i)] = empty;
79 const int dim = GV::dimension;
80 typedef typename GV::Grid::ctype ctype;
82 interior.resize(gv.indexSet().size(
dim));
83 for(
int i=0; i< interior.size(); i++)
87 for(
const auto& entity : elements(gv))
90 for (
const auto& intersection : intersections(gv,entity))
92 if (intersection.boundary())
95 unsigned int f = intersection.indexInInside();
97 auto refelem = Dune::ReferenceElements<ctype,dim>::simplex();
98 assert(entity.geometry().type().isSimplex() &&
"InteriorNodeConstraints only work for simplicial meshes");
99 unsigned int sz = refelem.size(f,1,
dim);
101 for (
unsigned int v = 0; v < sz; ++v)
103 unsigned int local_idx = refelem.subEntity (f,1, v,
dim);
104 unsigned int idx = gv.indexSet().subIndex(entity, local_idx,
dim);
105 interior[idx] =
false;
static const int dim
Definition: adaptivity.hh:84
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:25
@ doSkeleton
Definition: interiornode.hh:30
void updateInteriorNodes(const GV &gv)
Definition: interiornode.hh:76
@ doProcessor
Definition: interiornode.hh:29
@ doVolume
Definition: interiornode.hh:31
const std::vector< bool > & interiorNodes() const
Definition: interiornode.hh:70
void volume(const P ¶m, const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:41
@ doBoundary
Definition: interiornode.hh:28