2 #ifndef DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH 3 #define DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_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 P,
typename EG,
typename LFS,
typename T>
37 void volume (
const P& param,
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");
54 unsigned int local_idx = key.subEntity();
57 unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx,
dim);
74 const int dim = GV::dimension;
75 typedef typename GV::Grid::ctype ctype;
77 interior.resize(gv.indexSet().size(dim));
78 for(
int i=0; i< interior.size(); i++)
82 for(
const auto& entity : cells(gv))
85 for (
const auto& intersection : intersection(gv,entity))
87 if (intersection.boundary())
90 unsigned int f = intersection.indexInInside();
92 const ReferenceElement<ctype,dim> & refelem =
93 ReferenceElements<ctype,dim>::simplex();
94 assert(entity.geometry().type().isSimplex() &&
"InteriorNodeConstraints only work for simplicial meshes");
95 unsigned int sz = refelem.size(f,1, dim);
97 for (
unsigned int v = 0; v < sz; ++v)
99 unsigned int local_idx = refelem.subEntity (f,1, v,dim);
100 unsigned int idx = gv.indexSet().subIndex(entity, local_idx, dim);
101 interior[idx] =
false;
113 #endif // DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH static const int dim
Definition: adaptivity.hh:83
const std::vector< bool > & interiorNodes() const
Definition: interiornode.hh:65
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:20
Definition: interiornode.hh:25
void volume(const P ¶m, const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:37
void updateInteriorNodes(const GV &gv)
Definition: interiornode.hh:71
Definition: interiornode.hh:24
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
Definition: interiornode.hh:26
Definition: interiornode.hh:27