dune-pdelab  2.7-git
fastdg/assembler.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil -*-
2 // vi: set et ts=2 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_ASSEMBLER_HH
4 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_ASSEMBLER_HH
5 
6 #include <dune/common/typetraits.hh>
12 
13 namespace Dune{
14  namespace PDELab{
15 
24  template<typename GFSU, typename GFSV, typename CU, typename CV, bool nonoverlapping_mode=false>
26  public:
27 
30  using EntitySet = typename GFSU::Traits::EntitySet;
31  using Element = typename EntitySet::Element;
32  using Intersection = typename EntitySet::Intersection;
34 
37  typedef GFSU TrialGridFunctionSpace;
38  typedef GFSV TestGridFunctionSpace;
40 
42  typedef typename GFSU::Traits::SizeType SizeType;
43 
46 
47  FastDGAssembler (const GFSU& gfsu_, const GFSV& gfsv_, const CU& cu_, const CV& cv_)
48  : gfsu(gfsu_)
49  , gfsv(gfsv_)
50  , cu(cu_)
51  , cv(cv_)
52  , lfsu(gfsu_)
53  , lfsv(gfsv_)
54  , lfsun(gfsu_)
55  , lfsvn(gfsv_)
56  { }
57 
58  FastDGAssembler (const GFSU& gfsu_, const GFSV& gfsv_)
59  : gfsu(gfsu_)
60  , gfsv(gfsv_)
61  , cu()
62  , cv()
63  , lfsu(gfsu_)
64  , lfsv(gfsv_)
65  , lfsun(gfsu_)
66  , lfsvn(gfsv_)
67  { }
68 
70  const GFSU& trialGridFunctionSpace() const
71  {
72  return gfsu;
73  }
74 
76  const GFSV& testGridFunctionSpace() const
77  {
78  return gfsv;
79  }
80 
81  // Assembler (const GFSU& gfsu_, const GFSV& gfsv_)
82  // : gfsu(gfsu_), gfsv(gfsv_), lfsu(gfsu_), lfsv(gfsv_),
83  // lfsun(gfsu_), lfsvn(gfsv_),
84  // sub_triangulation(ST(gfsu_.gridview(),Dune::PDELab::NoSubTriangulationImp()))
85  // { }
86 
87  template<class LocalAssemblerEngine>
88  void assemble(LocalAssemblerEngine & assembler_engine) const
89  {
90  const bool fast = true;
91  typedef LFSIndexCache<LFSU,CU,fast> LFSUCache;
92 
93  typedef LFSIndexCache<LFSV,CV,fast> LFSVCache;
94 
95  const bool needs_constraints_caching = assembler_engine.needsConstraintsCaching(cu,cv);
96 
97  LFSUCache lfsu_cache(lfsu,cu,needs_constraints_caching);
98  LFSVCache lfsv_cache(lfsv,cv,needs_constraints_caching);
99  LFSUCache lfsun_cache(lfsun,cu,needs_constraints_caching);
100  LFSVCache lfsvn_cache(lfsvn,cv,needs_constraints_caching);
101 
102  // Notify assembler engine about oncoming assembly
103  assembler_engine.preAssembly();
104 
105  // Extract integration requirements from the local assembler
106  const bool require_uv_skeleton = assembler_engine.requireUVSkeleton();
107  const bool require_v_skeleton = assembler_engine.requireVSkeleton();
108  const bool require_uv_boundary = assembler_engine.requireUVBoundary();
109  const bool require_v_boundary = assembler_engine.requireVBoundary();
110  const bool require_uv_processor = assembler_engine.requireUVBoundary();
111  const bool require_v_processor = assembler_engine.requireVBoundary();
112  const bool require_uv_post_skeleton = assembler_engine.requireUVVolumePostSkeleton();
113  const bool require_v_post_skeleton = assembler_engine.requireVVolumePostSkeleton();
114  const bool require_skeleton_two_sided = assembler_engine.requireSkeletonTwoSided();
115 
116  auto entity_set = gfsu.entitySet();
117  auto& index_set = entity_set.indexSet();
118 
119  // Traverse grid view
120  for (const auto& element : elements(entity_set,assembler_engine.partition()))
121  {
122  // Compute unique id
123  auto ids = index_set.uniqueIndex(element);
124 
125  ElementGeometry<Element> eg(element);
126 
127  if(assembler_engine.assembleCell(eg))
128  continue;
129 
130  // Bind local test function space to element
131  lfsv.bind(element,std::integral_constant<bool,fast>{});
132  lfsv_cache.update();
133 
134  // Notify assembler engine about bind
135  assembler_engine.onBindLFSV(eg,lfsv_cache);
136 
137  // Volume integration
138  assembler_engine.assembleVVolume(eg,lfsv_cache);
139 
140  // Bind local trial function space to element
141  lfsu.bind(element,std::integral_constant<bool,fast>{});
142  lfsu_cache.update();
143 
144  // Notify assembler engine about bind
145  assembler_engine.onBindLFSUV(eg,lfsu_cache,lfsv_cache);
146 
147  // Load coefficients of local functions
148  assembler_engine.loadCoefficientsLFSUInside(lfsu_cache);
149 
150  // Volume integration
151  assembler_engine.assembleUVVolume(eg,lfsu_cache,lfsv_cache);
152 
153  // Skip if no intersection iterator is needed
154  if (require_uv_skeleton || require_v_skeleton ||
155  require_uv_boundary || require_v_boundary ||
156  require_uv_processor || require_v_processor)
157  {
158  // Traverse intersections
159  unsigned int intersection_index = 0;
160  for(const auto& intersection : intersections(entity_set,element))
161  {
162 
163  IntersectionGeometry<Intersection> ig(intersection,intersection_index);
164 
165  auto intersection_data = classifyIntersection(entity_set,intersection);
166  auto intersection_type = std::get<0>(intersection_data);
167  auto& outside_element = std::get<1>(intersection_data);
168 
169  switch (intersection_type)
170  {
172  // the specific ordering of the if-statements in the old code caused periodic
173  // boundary intersection to be handled the same as skeleton intersections
175  if (require_uv_skeleton || require_v_skeleton)
176  {
177  // compute unique id for neighbor
178 
179  auto idn = index_set.uniqueIndex(outside_element);
180 
181  // Visit face if id is bigger
182  bool visit_face = ids > idn
183  or require_skeleton_two_sided
184  or (assembler_engine.partition() == Partitions::interiorBorder
185  and outside_element.partitionType() != InteriorEntity);
186 
187  // unique vist of intersection
188  if (visit_face)
189  {
190  // Bind local test space to neighbor element
191  lfsvn.bind(outside_element,std::integral_constant<bool,fast>{});
192  lfsvn_cache.update();
193 
194  // Notify assembler engine about binds
195  assembler_engine.onBindLFSVOutside(ig,lfsv_cache,lfsvn_cache);
196 
197  // Skeleton integration
198  assembler_engine.assembleVSkeleton(ig,lfsv_cache,lfsvn_cache);
199 
200  if(require_uv_skeleton){
201 
202  // Bind local trial space to neighbor element
203  lfsun.bind(outside_element,std::integral_constant<bool,fast>{});
204  lfsun_cache.update();
205 
206  // Notify assembler engine about binds
207  assembler_engine.onBindLFSUVOutside(ig,
208  lfsu_cache,lfsv_cache,
209  lfsun_cache,lfsvn_cache);
210 
211  // Load coefficients of local functions
212  assembler_engine.loadCoefficientsLFSUOutside(lfsun_cache);
213 
214  // Skeleton integration
215  assembler_engine.assembleUVSkeleton(ig,lfsu_cache,lfsv_cache,lfsun_cache,lfsvn_cache);
216 
217  // Notify assembler engine about unbinds
218  assembler_engine.onUnbindLFSUVOutside(ig,
219  lfsu_cache,lfsv_cache,
220  lfsun_cache,lfsvn_cache);
221  }
222 
223  // Notify assembler engine about unbinds
224  assembler_engine.onUnbindLFSVOutside(ig,lfsv_cache,lfsvn_cache);
225  }
226  }
227  break;
228 
230  if(require_uv_boundary || require_v_boundary )
231  {
232 
233  // Boundary integration
234  assembler_engine.assembleVBoundary(ig,lfsv_cache);
235 
236  if(require_uv_boundary){
237  // Boundary integration
238  assembler_engine.assembleUVBoundary(ig,lfsu_cache,lfsv_cache);
239  }
240  }
241  break;
242 
244  if(require_uv_processor || require_v_processor )
245  {
246 
247  // Processor integration
248  assembler_engine.assembleVProcessor(ig,lfsv_cache);
249 
250  if(require_uv_processor){
251  // Processor integration
252  assembler_engine.assembleUVProcessor(ig,lfsu_cache,lfsv_cache);
253  }
254  }
255  break;
256  } // switch
257 
258  ++intersection_index;
259  } // iit
260  } // do skeleton
261 
262  if(require_uv_post_skeleton || require_v_post_skeleton){
263  // Volume integration
264  assembler_engine.assembleVVolumePostSkeleton(eg,lfsv_cache);
265 
266  if(require_uv_post_skeleton){
267  // Volume integration
268  assembler_engine.assembleUVVolumePostSkeleton(eg,lfsu_cache,lfsv_cache);
269  }
270  }
271 
272  // Notify assembler engine about unbinds
273  assembler_engine.onUnbindLFSUV(eg,lfsu_cache,lfsv_cache);
274 
275  // Notify assembler engine about unbinds
276  assembler_engine.onUnbindLFSV(eg,lfsv_cache);
277 
278  } // it
279 
280  // Notify assembler engine that assembly is finished
281  assembler_engine.postAssembly(gfsu,gfsv);
282 
283  }
284 
285  private:
286 
287  /* global function spaces */
288  const GFSU& gfsu;
289  const GFSV& gfsv;
290 
291  typename std::conditional<
293  const CU,
294  const CU&
295  >::type cu;
296  typename std::conditional<
298  const CV,
299  const CV&
300  >::type cv;
301 
302  /* local function spaces */
305  // local function spaces in local cell
306  mutable LFSU lfsu;
307  mutable LFSV lfsv;
308  // local function spaces in neighbor
309  mutable LFSU lfsun;
310  mutable LFSV lfsvn;
311 
312  }; // end class FastDGAssembler
313  } // end namespace PDELab
314 } // end namespace Dune
315 #endif
Dune::PDELab::LocalAssemblerEngine::requireVBoundary
bool requireVBoundary() const
Dune::PDELab::FastDGAssembler< GFSU, GFSV, Dune::PDELab::EmptyTransformation, Dune::PDELab::EmptyTransformation >::Element
typename EntitySet::Element Element
Definition: fastdg/assembler.hh:31
Dune::PDELab::IntersectionType::processor
@ processor
processor boundary intersection (neighbor() == false && boundary() == false) or outside entity not in...
Dune::PDELab::LocalAssemblerEngine::postAssembly
void postAssembly()
Called last thing after assembling.
Dune::PDELab::LocalAssemblerEngine::assembleCell
bool assembleCell(const EG &eg)
Dune::PDELab::FastDGAssembler::FastDGAssembler
FastDGAssembler(const GFSU &gfsu_, const GFSV &gfsv_, const CU &cu_, const CV &cv_)
Definition: fastdg/assembler.hh:47
Dune::PDELab::LocalAssemblerEngine::onUnbindLFSUVOutside
void onUnbindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
Dune::PDELab::IntersectionType::skeleton
@ skeleton
skeleton intersection (neighbor() == true && boundary() == false)
Dune::PDELab::LocalAssemblerEngine::onUnbindLFSUV
void onUnbindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
Dune::PDELab::FastDGAssembler::TrialGridFunctionSpace
GFSU TrialGridFunctionSpace
Definition: fastdg/assembler.hh:37
Dune::PDELab::LocalAssemblerEngine::assembleVSkeleton
void assembleVSkeleton(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
Dune::PDELab::LocalFunctionSpace< GFSU, TrialSpaceTag >
Dune::PDELab::LocalAssemblerEngine::assembleVBoundary
void assembleVBoundary(const IG &ig, const LFSV_S &lfsv_s)
Dune::PDELab::LFSIndexCache
Definition: lfsindexcache.hh:977
Dune::PDELab::FastDGAssembler::SizeType
GFSU::Traits::SizeType SizeType
Size type as used in grid function space.
Definition: fastdg/assembler.hh:42
Dune::PDELab::classifyIntersection
std::tuple< IntersectionType, typename EntitySet::Element > classifyIntersection(const EntitySet &entity_set, const Intersection &is)
Classifies the type of an intersection wrt to the passed EntitySet.
Definition: intersectiontype.hh:37
Dune::PDELab::IntersectionGeometry
Wrap intersection.
Definition: geometrywrapper.hh:56
Dune::PDELab::LocalAssemblerEngine::onBindLFSV
void onBindLFSV(const EG &eg, const LFSV_S &lfsv_s)
Dune::PDELab::FastDGAssembler::TestGridFunctionSpace
GFSV TestGridFunctionSpace
Definition: fastdg/assembler.hh:38
Dune::PDELab::LocalAssemblerEngine::requireUVVolumePostSkeleton
bool requireUVVolumePostSkeleton() const
Dune::PDELab::FastDGAssembler< GFSU, GFSV, Dune::PDELab::EmptyTransformation, Dune::PDELab::EmptyTransformation >::Intersection
typename EntitySet::Intersection Intersection
Definition: fastdg/assembler.hh:32
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::LocalAssemblerEngine::assembleVVolume
void assembleVVolume(const EG &eg, const LFSV &lfsv)
Dune::PDELab::FastDGAssembler::isGalerkinMethod
static const bool isGalerkinMethod
Static check on whether this is a Galerkin method.
Definition: fastdg/assembler.hh:45
Dune::PDELab::LocalAssemblerEngine::preAssembly
void preAssembly()
Called directly before assembling.
Dune::PDELab::IntersectionType::periodic
@ periodic
periodic boundary intersection (neighbor() == true && boundary() == true)
Dune::PDELab::LocalAssemblerEngine::assembleVVolumePostSkeleton
void assembleVVolumePostSkeleton(const EG &eg, const LFSV &lfsv)
Dune::PDELab::IntersectionType::boundary
@ boundary
domain boundary intersection (neighbor() == false && boundary() == true)
Dune::PDELab::FastDGAssembler::trialGridFunctionSpace
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: fastdg/assembler.hh:70
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::LocalAssemblerEngine::loadCoefficientsLFSUOutside
void loadCoefficientsLFSUOutside(const LFSU_N &lfsu_n)
Dune::PDELab::FastDGAssembler::FastDGAssembler
FastDGAssembler(const GFSU &gfsu_, const GFSV &gfsv_)
Definition: fastdg/assembler.hh:58
Dune::PDELab::LocalAssemblerEngine::onBindLFSUVOutside
void onBindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
Dune::PDELab::LocalAssemblerEngine::assembleUVProcessor
void assembleUVProcessor(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
Dune::PDELab::LocalAssemblerEngine::onUnbindLFSV
void onUnbindLFSV(const EG &eg, const LFSV_S &lfsv_s)
Dune::PDELab::FastDGAssembler::testGridFunctionSpace
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: fastdg/assembler.hh:76
Dune::PDELab::LocalAssemblerEngine::requireVVolumePostSkeleton
bool requireVVolumePostSkeleton() const
Dune::PDELab::LocalAssemblerEngine
The local assembler engine which handles the integration parts as provided by the global assemblers.
Definition: common/assembler.hh:33
geometrywrapper.hh
Dune::PDELab::LocalAssemblerEngine::loadCoefficientsLFSUInside
void loadCoefficientsLFSUInside(const LFSU_S &lfsu_s)
Dune::PDELab::ElementGeometry
Wrap element.
Definition: geometrywrapper.hh:15
Dune::PDELab::LocalAssemblerEngine::assembleUVVolume
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
Dune::PDELab::FastDGAssembler::assemble
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition: fastdg/assembler.hh:88
Dune::PDELab::LocalAssemblerEngine::onBindLFSVOutside
void onBindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
localfunctionspace.hh
intersectiontype.hh
Dune::PDELab::LocalAssemblerEngine::requireSkeletonTwoSided
bool requireSkeletonTwoSided() const
lfsindexcache.hh
Dune::PDELab::LocalAssemblerEngine::onBindLFSUV
void onBindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
ig
const IG & ig
Definition: constraints.hh:149
Dune::PDELab::LocalAssemblerEngine::onUnbindLFSVOutside
void onUnbindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
Dune::PDELab::LocalAssemblerEngine::assembleUVSkeleton
void assembleUVSkeleton(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
Dune::PDELab::LocalAssemblerEngine::assembleUVBoundary
void assembleUVBoundary(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
Dune::PDELab::FastDGAssembler
The fast DG assembler for standard DUNE grid.
Definition: fastdg/assembler.hh:25
Dune::PDELab::FastDGAssembler< GFSU, GFSV, Dune::PDELab::EmptyTransformation, Dune::PDELab::EmptyTransformation >::EntitySet
typename GFSU::Traits::EntitySet EntitySet
Definition: fastdg/assembler.hh:30
assemblerutilities.hh
Dune::PDELab::LocalAssemblerEngine::assembleUVVolumePostSkeleton
void assembleUVVolumePostSkeleton(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
Dune::PDELab::LocalAssemblerEngine::assembleVProcessor
void assembleVProcessor(const IG &ig, const LFSV_S &lfsv_s)
Dune::PDELab::LocalAssemblerEngine::requireVSkeleton
bool requireVSkeleton() const
Dune::PDELab::LocalAssemblerEngine::requireUVSkeleton
bool requireUVSkeleton() const
Dune::PDELab::LocalAssemblerEngine::requireUVBoundary
bool requireUVBoundary() const