dune-pdelab  2.4.1
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 P, typename EG, typename LFS, typename T>
37  void volume (const P& param, 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 
53  // subentity index
54  unsigned int local_idx = key.subEntity();
55 
56  // global idx
57  unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx, dim);
58 
59  // update constraints
60  if (interior[idx])
61  trafo[i] = empty;
62  }
63  }
64 
65  const std::vector<bool> & interiorNodes() const
66  {
67  return interior;
68  }
69 
70  template<typename GV>
71  void updateInteriorNodes(const GV & gv)
72  {
73  // update vector size
74  const int dim = GV::dimension;
75  typedef typename GV::Grid::ctype ctype;
76 
77  interior.resize(gv.indexSet().size(dim));
78  for(int i=0; i< interior.size(); i++)
79  interior[i] = true;
80 
81  // loop over all cells
82  for(const auto& entity : cells(gv))
83  {
84  // find boundary faces & associated vertices
85  for (const auto& intersection : intersection(gv,entity))
86  {
87  if (intersection.boundary())
88  {
89  // boundary face
90  unsigned int f = intersection.indexInInside();
91  // remember associated vertices
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);
96  assert(sz == dim);
97  for (unsigned int v = 0; v < sz; ++v)
98  {
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;
102  }
103  }
104  }
105  }
106  }
107  };
109 
110  } // end namespace PDELab
111 } // end namespace Dune
112 
113 #endif // DUNE_PDELAB_INTERIORNODECONSTRAINTS_HH
static const int dim
Definition: adaptivity.hh:83
Definition: adaptivity.hh:27
const std::vector< bool > & interiorNodes() const
Definition: interiornode.hh:65
void volume(const P &param, const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:37
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:20
void updateInteriorNodes(const GV &gv)
Definition: interiornode.hh:71