DUNE PDELab (git)

interiornode.hh
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#include <dune/typetree/typetree.hh>
13
14namespace Dune {
15 namespace PDELab {
16
20
24 : public TypeTree::LeafNode
25 {
26 std::vector<bool> interior;
27 public:
28 enum{doBoundary=false};
29 enum{doProcessor=false};
30 enum{doSkeleton=false};
31 enum{doVolume=true};
32
34
40 template<typename P, typename EG, typename LFS, typename T>
41 void volume (const P& param, const EG& eg, const LFS& lfs, T& trafo) const
42 {
43 typedef typename EG::Entity Entity;
44 enum { dim = Entity::dimension };
45
46 // update component
47 typename T::RowType empty;
48 typedef typename LFS::Traits::SizeType size_type;
50 typename LFS::Traits::FiniteElementType
51 > FESwitch;
52 for (size_type i=0; i<lfs.size(); i++){
53 const LocalKey& key = FESwitch::coefficients(lfs.finiteElement()).localKey(i);
54 assert(key.codim() == dim && "InteriorNodeConstraints only work for vertex DOFs");
55 assert(key.index() == 0 && "InteriorNodeConstraints only work for P1 shape functions");
56
57 // subentity index
58 unsigned int local_idx = key.subEntity();
59
60 // global idx
61 unsigned int idx = lfs.gridFunctionSpace().gridView().indexSet().subIndex(eg.entity(), local_idx, dim);
62
63 // update constraints
64 if (interior[idx])
65 trafo[lfs.dofIndex(i)] = empty;
66 }
67
68 }
69
70 const std::vector<bool> & interiorNodes() const
71 {
72 return interior;
73 }
74
75 template<typename GV>
76 void updateInteriorNodes(const GV & gv)
77 {
78 // update vector size
79 const int dim = GV::dimension;
80 typedef typename GV::Grid::ctype ctype;
81
82 interior.resize(gv.indexSet().size(dim));
83 for(std::size_t i=0; i< interior.size(); i++)
84 interior[i] = true;
85
86 // loop over all cells
87 for(const auto& entity : elements(gv))
88 {
89 // find boundary faces & associated vertices
90 for (const auto& intersection : intersections(gv,entity))
91 {
92 if (intersection.boundary())
93 {
94 // boundary face
95 unsigned int f = intersection.indexInInside();
96 // remember associated vertices
98 assert(entity.geometry().type().isSimplex() && "InteriorNodeConstraints only work for simplicial meshes");
99 unsigned int sz = refelem.size(f,1, dim);
100 assert(sz == dim);
101 for (unsigned int v = 0; v < sz; ++v)
102 {
103 unsigned int local_idx = refelem.subEntity (f,1, v,dim);
104 unsigned int idx = gv.indexSet().subIndex(entity, local_idx, dim);
105 interior[idx] = false;
106 }
107 }
108 }
109 }
110 }
111 };
113
114 } // end namespace PDELab
115} // end namespace Dune
116
117#endif // DUNE_PDELAB_CONSTRAINTS_INTERIORNODE_HH
Wrapper class for entities.
Definition: entity.hh:66
static constexpr int dimension
Know the grid dimension.
Definition: entity.hh:109
Describe position of one degree of freedom.
Definition: localkey.hh:24
constexpr unsigned int index() const noexcept
Return offset within subentity.
Definition: localkey.hh:70
constexpr unsigned int codim() const noexcept
Return codim of associated entity.
Definition: localkey.hh:63
constexpr unsigned int subEntity() const noexcept
Return number of associated subentity.
Definition: localkey.hh:56
constraints all DOFs associated with interior vertices This allows to implement surface FEM using sta...
Definition: interiornode.hh:25
void volume(const P &param, const EG &eg, const LFS &lfs, T &trafo) const
volume constraints
Definition: interiornode.hh:41
Base class for leaf nodes in a dune-typetree.
Definition: leafnode.hh:28
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::bool_constant<(sizeof...(II)==0)> empty(std::integer_sequence< T, II... >)
Checks whether the sequence is empty.
Definition: integersequence.hh:80
Switch for uniform treatment of finite element with either the local or the global interface.
Definition: interfaceswitch.hh:30
static const ReferenceElement & simplex()
get simplex reference elements
Definition: referenceelements.hh:162
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)