2 #ifndef DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
3 #define DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
5 #include <dune/common/array.hh>
6 #include <dune/grid/common/gridenums.hh>
7 #include <dune/geometry/referenceelements.hh>
8 #include <dune/localfunctions/common/interfaceswitch.hh>
9 #include <dune/localfunctions/common/localkey.hh>
22 std::vector<bool> interior;
36 template<
typename EG,
typename LFS,
typename T>
37 void volume (
const EG&
eg,
const LFS& lfs, T& trafo)
const
39 typedef typename EG::Entity Entity;
40 enum {
dim = Entity::dimension, dimw = Entity::dimensionworld };
43 typename T::RowType empty;
44 typedef typename LFS::Traits::SizeType size_type;
45 typedef FiniteElementInterfaceSwitch<
46 typename LFS::Traits::FiniteElementType
48 for (size_type i=0; i<lfs.size(); i++){
49 const LocalKey& key = FESwitch::coefficients(lfs.finiteElement()).localKey(i);
50 assert(key.codim() ==
dim &&
"InteriorNodeConstraints only work for vertex DOFs");
51 assert(key.index() == 0 &&
"InteriorNodeConstraints only work for P1 shape functions");
53 unsigned int local_idx = key.subEntity();
56 unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx,
dim);
73 const int dim = GV::dimension;
74 typedef typename GV::Grid::ctype ctype;
77 interior.resize(gv.indexSet().size(dim));
78 for(
int i=0; i< interior.size(); i++)
81 typedef typename GV::template Codim<0>::Iterator ElementIterator;
82 typedef typename ElementIterator::Entity Entity;
85 for(ElementIterator it = gv.template begin<0>(); it != gv.template end<0>(); ++it)
87 const Entity &entity = *it;
90 typedef typename GV::IntersectionIterator IntersectionIterator;
91 IntersectionIterator iit = gv.ibegin(entity);
92 const IntersectionIterator iend = gv.iend(entity);
93 for (; iit != iend; ++iit)
98 unsigned int f = iit->indexInInside();
100 const ReferenceElement<ctype,dim> & refelem =
101 ReferenceElements<ctype,dim>::simplex();
102 assert(entity.geometry().type().isSimplex() &&
"InteriorNodeConstraints only work for simplicial meshes");
103 unsigned int sz = refelem.size(f,1, dim);
105 for (
unsigned int v = 0; v < sz; ++v)
107 unsigned int local_idx = refelem.subEntity (f,1, v,dim);
108 unsigned int idx = gv.indexSet().subIndex(entity, local_idx, dim);
109 interior[idx] =
false;
121 #endif // DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
void updateInteriorNodes(const GV &gv)
Definition: interiornode.hh:70
static const int dim
Definition: adaptivity.hh:82
const std::vector< bool > & interiorNodes() const
Definition: interiornode.hh:64
void volume(const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:37
Definition: interiornode.hh:24
Definition: interiornode.hh:25
Definition: interiornode.hh:27
Definition: interiornode.hh:26
const F & f
Definition: common/constraints.hh:145
const EG & eg
Definition: common/constraints.hh:277
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:20