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