dune-pdelab  2.5-dev
interiornode.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH
3 #define DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH
4 
5 #include <array>
6 
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>
11 
12 namespace Dune {
13  namespace PDELab {
14 
18 
22  {
23  std::vector<bool> interior;
24  public:
25  enum{doBoundary=false};
26  enum{doProcessor=false};
27  enum{doSkeleton=false};
28  enum{doVolume=true};
29 
31 
37  template<typename P, typename EG, typename LFS, typename T>
38  void volume (const P& param, const EG& eg, const LFS& lfs, T& trafo) const
39  {
40  typedef typename EG::Entity Entity;
41  enum { dim = Entity::dimension, dimw = Entity::dimensionworld };
42 
43  // update component
44  typename T::RowType empty;
45  typedef typename LFS::Traits::SizeType size_type;
46  typedef FiniteElementInterfaceSwitch<
47  typename LFS::Traits::FiniteElementType
48  > FESwitch;
49  for (size_type i=0; i<lfs.size(); i++){
50  const LocalKey& key = FESwitch::coefficients(lfs.finiteElement()).localKey(i);
51  assert(key.codim() == dim && "InteriorNodeConstraints only work for vertex DOFs");
52  assert(key.index() == 0 && "InteriorNodeConstraints only work for P1 shape functions");
53 
54  // subentity index
55  unsigned int local_idx = key.subEntity();
56 
57  // global idx
58  unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx, dim);
59 
60  // update constraints
61  if (interior[idx])
62  trafo[i] = empty;
63  }
64  }
65 
66  const std::vector<bool> & interiorNodes() const
67  {
68  return interior;
69  }
70 
71  template<typename GV>
72  void updateInteriorNodes(const GV & gv)
73  {
74  // update vector size
75  const int dim = GV::dimension;
76  typedef typename GV::Grid::ctype ctype;
77 
78  interior.resize(gv.indexSet().size(dim));
79  for(int i=0; i< interior.size(); i++)
80  interior[i] = true;
81 
82  // loop over all cells
83  for(const auto& entity : cells(gv))
84  {
85  // find boundary faces & associated vertices
86  for (const auto& intersection : intersection(gv,entity))
87  {
88  if (intersection.boundary())
89  {
90  // boundary face
91  unsigned int f = intersection.indexInInside();
92  // remember associated vertices
93  auto refelem = Dune::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_CONSTRAINTS_INTERIORNODE_HH
Dune::PDELab::InteriorNodeConstraints::doSkeleton
@ doSkeleton
Definition: interiornode.hh:27
dim
static const int dim
Definition: adaptivity.hh:84
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::InteriorNodeConstraints::interiorNodes
const std::vector< bool > & interiorNodes() const
Definition: interiornode.hh:66
Dune::PDELab::InteriorNodeConstraints::doBoundary
@ doBoundary
Definition: interiornode.hh:25
Dune::PDELab::InteriorNodeConstraints::volume
void volume(const P &param, const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:38
Dune::PDELab::InteriorNodeConstraints::doProcessor
@ doProcessor
Definition: interiornode.hh:26
Dune::PDELab::InteriorNodeConstraints
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:21
Dune::PDELab::InteriorNodeConstraints::updateInteriorNodes
void updateInteriorNodes(const GV &gv)
Definition: interiornode.hh:72
Dune::PDELab::InteriorNodeConstraints::doVolume
@ doVolume
Definition: interiornode.hh:28