dune-pdelab  2.0.0
interiornode.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
3 #define DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
4 
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>
10 
11 namespace Dune {
12  namespace PDELab {
13 
17 
21  {
22  std::vector<bool> interior;
23  public:
24  enum{doBoundary=false};
25  enum{doProcessor=false};
26  enum{doSkeleton=false};
27  enum{doVolume=true};
28 
30 
36  template<typename EG, typename LFS, typename T>
37  void volume (const EG& eg, const LFS& lfs, T& trafo) const
38  {
39  typedef typename EG::Entity Entity;
40  enum { dim = Entity::dimension, dimw = Entity::dimensionworld };
41 
42  // update component
43  typename T::RowType empty;
44  typedef typename LFS::Traits::SizeType size_type;
45  typedef FiniteElementInterfaceSwitch<
46  typename LFS::Traits::FiniteElementType
47  > FESwitch;
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");
52  // subentity index
53  unsigned int local_idx = key.subEntity();
54  // global idx
55 
56  unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx, dim);
57 
58  // update constraints
59  if (interior[idx])
60  trafo[i] = empty;
61  }
62  }
63 
64  const std::vector<bool> & interiorNodes() const
65  {
66  return interior;
67  }
68 
69  template<typename GV>
70  void updateInteriorNodes(const GV & gv)
71  {
72  // update vector size
73  const int dim = GV::dimension;
74  typedef typename GV::Grid::ctype ctype;
75 
76 
77  interior.resize(gv.indexSet().size(dim));
78  for(int i=0; i< interior.size(); i++)
79  interior[i] = true;
80 
81  typedef typename GV::template Codim<0>::Iterator ElementIterator;
82  typedef typename ElementIterator::Entity Entity;
83 
84  // loop over all cells
85  for(ElementIterator it = gv.template begin<0>(); it != gv.template end<0>(); ++it)
86  {
87  const Entity &entity = *it;
88 
89  // find boundary faces & associated vertices
90  typedef typename GV::IntersectionIterator IntersectionIterator;
91  IntersectionIterator iit = gv.ibegin(entity);
92  const IntersectionIterator iend = gv.iend(entity);
93  for (; iit != iend; ++iit)
94  {
95  if (iit->boundary())
96  {
97  // boundary face
98  unsigned int f = iit->indexInInside();
99  // remember associated vertices
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);
104  assert(sz == dim);
105  for (unsigned int v = 0; v < sz; ++v)
106  {
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;
110  }
111  }
112  }
113  }
114  }
115  };
117 
118  } // end namespace PDELab
119 } // end namespace Dune
120 
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
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