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