4 #ifndef DUNE_PDELAB_CONSTRAINTS_CONFORMING_HH
5 #define DUNE_PDELAB_CONSTRAINTS_CONFORMING_HH
10 #include <dune/common/exceptions.hh>
12 #include <dune/geometry/referenceelements.hh>
13 #include <dune/geometry/type.hh>
15 #include <dune/grid/common/grid.hh>
17 #include <dune/localfunctions/common/interfaceswitch.hh>
19 #include <dune/typetree/typetree.hh>
37 :
public TypeTree::LeafNode
52 template<
typename P,
typename IG,
typename LFS,
typename T>
53 void boundary (
const P& param,
const IG&
ig,
const LFS& lfs, T& trafo)
const
55 typedef FiniteElementInterfaceSwitch<
56 typename LFS::Traits::FiniteElementType
58 typedef FieldVector<typename IG::ctype, IG::mydimension> FaceCoord;
60 const int face =
ig.indexInInside();
63 auto refelem = referenceElement(
ig.inside().geometry());
64 auto face_refelem = referenceElement(
ig.geometry());
67 typename T::RowType empty;
69 const FaceCoord testpoint = face_refelem.position(0,0);
72 if (!param.isDirichlet(
ig,testpoint))
76 i<std::size_t(FESwitch::coefficients(lfs.finiteElement()).size());
81 FESwitch::coefficients(lfs.finiteElement()).localKey(i).codim();
83 if (codim==0)
continue;
85 for (
int j=0; j<refelem.size(face,1,codim); j++){
87 if (
static_cast<int>(FESwitch::coefficients(lfs.finiteElement()).
88 localKey(i).subEntity())
89 == refelem.subEntity(face,1,j,codim))
90 trafo[lfs.dofIndex(i)] = empty;
108 template<
typename IG,
typename LFS,
typename T>
111 typedef FiniteElementInterfaceSwitch<
112 typename LFS::Traits::FiniteElementType
116 const int face =
ig.indexInInside();
118 auto refelem = referenceElement(
ig.inside().geometry());
121 typename T::RowType empty;
124 for (
size_t i=0; i<FESwitch::coefficients(lfs.finiteElement()).size();
129 FESwitch::coefficients(lfs.finiteElement()).localKey(i).codim();
131 if (codim==0)
continue;
133 for (
int j=0; j<refelem.size(face,1,codim); j++)
134 if (FESwitch::coefficients(lfs.finiteElement()).localKey(i).
135 subEntity() == std::size_t(refelem.subEntity(face,1,j,codim)))
136 trafo[lfs.dofIndex(i)] = empty;
142 template<
typename GV>
154 template<
typename P,
typename EG,
typename LFS,
typename T>
155 void volume (
const P& param,
const EG& eg,
const LFS& lfs, T& trafo)
const
157 typedef FiniteElementInterfaceSwitch<
158 typename LFS::Traits::FiniteElementType
161 auto& entity = eg.entity();
164 if (entity.partitionType()==Dune::InteriorEntity)
167 typedef typename FESwitch::Coefficients Coefficients;
168 const Coefficients& coeffs = FESwitch::coefficients(lfs.finiteElement());
171 typename T::RowType empty;
173 auto ref_el = referenceElement(entity.geometry());
176 for (
size_t i = 0; i < coeffs.size(); ++i)
178 size_t codim = coeffs.localKey(i).codim();
179 size_t sub_entity = coeffs.localKey(i).subEntity();
181 size_t entity_index = _gv.indexSet().subIndex(entity,sub_entity,codim);
184 size_t index = _gt_offsets[gt_index] + entity_index;
188 trafo[lfs.dofIndex(i)] = empty;
196 std::fill(_gt_offsets.begin(),_gt_offsets.end(),0);
198 for (
size_t codim = 0; codim <= GV::dimension; ++codim)
200 if (gfs.ordering().contains(codim))
202 for (
auto gt : _gv.indexSet().types(codim))
207 std::partial_sum(_gt_offsets.begin(),_gt_offsets.end(),_gt_offsets.begin());
209 _ghosts.assign(_gt_offsets.back(),
true);
211 typedef typename GV::template Codim<0>::
212 template Partition<Interior_Partition>::Iterator Iterator;
214 for(Iterator it = _gv.template begin<0, Interior_Partition>(),
215 end = _gv.template end<0, Interior_Partition>();
221 auto ref_el = referenceElement(entity.geometry());
223 for (
size_t codim = 0; codim <= GV::dimension; ++codim)
224 if (gfs.ordering().contains(codim))
226 for (
int i = 0; i < ref_el.size(codim); ++i)
228 size_t entity_index = _gv.indexSet().subIndex(entity,i,codim);
230 size_t index = _gt_offsets[gt_index] + entity_index;
232 _ghosts[
index] =
false;
241 , _rank(gv.comm().rank())
242 , _gt_offsets(GlobalGeometryTypeIndex::size(GV::dimension) + 1)
249 std::vector<bool> _ghosts;
250 std::vector<size_t> _gt_offsets;
257 #endif // DUNE_PDELAB_CONSTRAINTS_CONFORMING_HH