dune-pdelab  2.0.0
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  typedef typename GFSU::Traits::GridViewType GV;
29  typedef typename GV::Traits::template Codim<0>::Iterator ElementIterator;
30  typedef typename GV::Traits::template Codim<0>::Entity Element;
31  typedef typename GV::IntersectionIterator IntersectionIterator;
32  typedef typename IntersectionIterator::Intersection Intersection;
34 
37  typedef GFSU TrialGridFunctionSpace;
38  typedef GFSV TestGridFunctionSpace;
40 
42  typedef typename GFSU::Traits::SizeType SizeType;
43 
46 
47  DefaultAssembler (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  DefaultAssembler (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  typedef typename GV::Traits::template Codim<0>::Entity Element;
91 
92  typedef LFSIndexCache<LFSU,CU> LFSUCache;
93 
94  typedef LFSIndexCache<LFSV,CV> LFSVCache;
95 
96  const bool needs_constraints_caching = assembler_engine.needsConstraintsCaching(cu,cv);
97 
98  LFSUCache lfsu_cache(lfsu,cu,needs_constraints_caching);
99  LFSVCache lfsv_cache(lfsv,cv,needs_constraints_caching);
100  LFSUCache lfsun_cache(lfsun,cu,needs_constraints_caching);
101  LFSVCache lfsvn_cache(lfsvn,cv,needs_constraints_caching);
102 
103  // Notify assembler engine about oncoming assembly
104  assembler_engine.preAssembly();
105 
106  // Map each cell to unique id
107  ElementMapper<GV> cell_mapper(gfsu.gridView());
108 
109  // Extract integration requirements from the local assembler
110  const bool require_uv_skeleton = assembler_engine.requireUVSkeleton();
111  const bool require_v_skeleton = assembler_engine.requireVSkeleton();
112  const bool require_uv_boundary = assembler_engine.requireUVBoundary();
113  const bool require_v_boundary = assembler_engine.requireVBoundary();
114  const bool require_uv_processor = assembler_engine.requireUVBoundary();
115  const bool require_v_processor = assembler_engine.requireVBoundary();
116  const bool require_uv_post_skeleton = assembler_engine.requireUVVolumePostSkeleton();
117  const bool require_v_post_skeleton = assembler_engine.requireVVolumePostSkeleton();
118  const bool require_skeleton_two_sided = assembler_engine.requireSkeletonTwoSided();
119 
120  // Traverse grid view
121  for (ElementIterator it = gfsu.gridView().template begin<0>();
122  it!=gfsu.gridView().template end<0>(); ++it)
123  {
124  // Compute unique id
125  const typename GV::IndexSet::IndexType ids = cell_mapper.map(*it);
126 
128 
129  if(assembler_engine.assembleCell(eg))
130  continue;
131 
132  // Bind local test function space to element
133  lfsv.bind( *it );
134  lfsv_cache.update();
135 
136  // Notify assembler engine about bind
137  assembler_engine.onBindLFSV(eg,lfsv_cache);
138 
139  // Volume integration
140  assembler_engine.assembleVVolume(eg,lfsv_cache);
141 
142  // Bind local trial function space to element
143  lfsu.bind( *it );
144  lfsu_cache.update();
145 
146  // Notify assembler engine about bind
147  assembler_engine.onBindLFSUV(eg,lfsu_cache,lfsv_cache);
148 
149  // Load coefficients of local functions
150  assembler_engine.loadCoefficientsLFSUInside(lfsu_cache);
151 
152  // Volume integration
153  assembler_engine.assembleUVVolume(eg,lfsu_cache,lfsv_cache);
154 
155  // Skip if no intersection iterator is needed
156  if (require_uv_skeleton || require_v_skeleton ||
157  require_uv_boundary || require_v_boundary ||
158  require_uv_processor || require_v_processor)
159  {
160  // Traverse intersections
161  unsigned int intersection_index = 0;
162  IntersectionIterator endit = gfsu.gridView().iend(*it);
163  IntersectionIterator iit = gfsu.gridView().ibegin(*it);
164  for(; iit!=endit; ++iit, ++intersection_index)
165  {
166 
167  IntersectionGeometry<Intersection> ig(*iit,intersection_index);
168 
169  switch (IntersectionType::get(*iit))
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  const typename GV::IndexSet::IndexType idn = cell_mapper.map(*(iit->outside()));
180 
181  // Visit face if id is bigger
182  bool visit_face = ids > idn || require_skeleton_two_sided;
183 
184  // unique vist of intersection
185  if (visit_face)
186  {
187  // Bind local test space to neighbor element
188  lfsvn.bind(*(iit->outside()));
189  lfsvn_cache.update();
190 
191  // Notify assembler engine about binds
192  assembler_engine.onBindLFSVOutside(ig,lfsv_cache,lfsvn_cache);
193 
194  // Skeleton integration
195  assembler_engine.assembleVSkeleton(ig,lfsv_cache,lfsvn_cache);
196 
197  if(require_uv_skeleton){
198 
199  // Bind local trial space to neighbor element
200  lfsun.bind(*(iit->outside()));
201  lfsun_cache.update();
202 
203  // Notify assembler engine about binds
204  assembler_engine.onBindLFSUVOutside(ig,
205  lfsu_cache,lfsv_cache,
206  lfsun_cache,lfsvn_cache);
207 
208  // Load coefficients of local functions
209  assembler_engine.loadCoefficientsLFSUOutside(lfsun_cache);
210 
211  // Skeleton integration
212  assembler_engine.assembleUVSkeleton(ig,lfsu_cache,lfsv_cache,lfsun_cache,lfsvn_cache);
213 
214  // Notify assembler engine about unbinds
215  assembler_engine.onUnbindLFSUVOutside(ig,
216  lfsu_cache,lfsv_cache,
217  lfsun_cache,lfsvn_cache);
218  }
219 
220  // Notify assembler engine about unbinds
221  assembler_engine.onUnbindLFSVOutside(ig,lfsv_cache,lfsvn_cache);
222  }
223  }
224  break;
225 
227  if(require_uv_boundary || require_v_boundary )
228  {
229 
230  // Boundary integration
231  assembler_engine.assembleVBoundary(ig,lfsv_cache);
232 
233  if(require_uv_boundary){
234  // Boundary integration
235  assembler_engine.assembleUVBoundary(ig,lfsu_cache,lfsv_cache);
236  }
237  }
238  break;
239 
241  if(require_uv_processor || require_v_processor )
242  {
243 
244  // Processor integration
245  assembler_engine.assembleVProcessor(ig,lfsv_cache);
246 
247  if(require_uv_processor){
248  // Processor integration
249  assembler_engine.assembleUVProcessor(ig,lfsu_cache,lfsv_cache);
250  }
251  }
252  break;
253  } // switch
254 
255  } // iit
256  } // do skeleton
257 
258  if(require_uv_post_skeleton || require_v_post_skeleton){
259  // Volume integration
260  assembler_engine.assembleVVolumePostSkeleton(eg,lfsv_cache);
261 
262  if(require_uv_post_skeleton){
263  // Volume integration
264  assembler_engine.assembleUVVolumePostSkeleton(eg,lfsu_cache,lfsv_cache);
265  }
266  }
267 
268  // Notify assembler engine about unbinds
269  assembler_engine.onUnbindLFSUV(eg,lfsu_cache,lfsv_cache);
270 
271  // Notify assembler engine about unbinds
272  assembler_engine.onUnbindLFSV(eg,lfsv_cache);
273 
274  } // it
275 
276  // Notify assembler engine that assembly is finished
277  assembler_engine.postAssembly(gfsu,gfsv);
278 
279  }
280 
281  private:
282 
283  /* global function spaces */
284  const GFSU& gfsu;
285  const GFSV& gfsv;
286 
287  typename conditional<
289  const CU,
290  const CU&
291  >::type cu;
292  typename conditional<
294  const CV,
295  const CV&
296  >::type cv;
297 
298  /* local function spaces */
301  // local function spaces in local cell
302  mutable LFSU lfsu;
303  mutable LFSV lfsv;
304  // local function spaces in neighbor
305  mutable LFSU lfsun;
306  mutable LFSV lfsvn;
307 
308  };
309 
310  }
311 }
312 #endif
Definition: lfsindexcache.hh:948
void assembleVProcessor(const IG &ig, const LFSV_S &lfsv_s)
void assembleVVolume(const EG &eg, const LFSV &lfsv)
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: default/assembler.hh:76
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)
The assembler for standard DUNE grid.
Definition: default/assembler.hh:23
void assembleVSkeleton(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void onBindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void loadCoefficientsLFSUOutside(const LFSU_N &lfsu_n)
static const bool isGalerkinMethod
Static check on whether this is a Galerkin method.
Definition: default/assembler.hh:45
DefaultAssembler(const GFSU &gfsu_, const GFSV &gfsv_, const CU &cu_, const CV &cv_)
Definition: default/assembler.hh:47
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
domain boundary intersection (neighbor() == false && boundary() == true)
Definition: assemblerutilities.hh:85
GFSU TrialGridFunctionSpace
Definition: default/assembler.hh:37
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: default/assembler.hh:70
void onUnbindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void onBindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void assembleVBoundary(const IG &ig, const LFSV_S &lfsv_s)
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)
void postAssembly()
Called last thing after assembling.
void assembleUVVolumePostSkeleton(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
skeleton intersection (neighbor() == true && boundary() == false)
Definition: assemblerutilities.hh:84
void onUnbindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
const IG & ig
Definition: common/constraints.hh:146
void preAssembly()
Called directly before assembling.
Class providing a consecutive index for codim 0 entities of a GridView.
Definition: elementmapper.hh:117
GV::Traits::template Codim< 0 >::Iterator ElementIterator
Definition: default/assembler.hh:29
IntersectionIterator::Intersection Intersection
Definition: default/assembler.hh:32
periodic boundary intersection (neighbor() == true && boundary() == true)
Definition: assemblerutilities.hh:86
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
processor boundary intersection (neighbor() == false && boundary() == false)
Definition: assemblerutilities.hh:83
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)
The local assembler engine which handles the integration parts as provided by the global assemblers...
Definition: common/assembler.hh:31
void onUnbindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
GFSU::Traits::GridViewType GV
Definition: default/assembler.hh:28
const EG & eg
Definition: common/constraints.hh:277
void onBindLFSV(const EG &eg, const LFSV_S &lfsv_s)
GV::Traits::template Codim< 0 >::Entity Element
Definition: default/assembler.hh:30
void loadCoefficientsLFSUInside(const LFSU_S &lfsu_s)
GV::IntersectionIterator IntersectionIterator
Definition: default/assembler.hh:31
void assembleUVProcessor(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition: default/assembler.hh:88
GFSV TestGridFunctionSpace
Definition: default/assembler.hh:38
void assembleVVolumePostSkeleton(const EG &eg, const LFSV &lfsv)
void assembleUVBoundary(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
DefaultAssembler(const GFSU &gfsu_, const GFSV &gfsv_)
Definition: default/assembler.hh:58
static Type get(const Intersection &is)
Returns the type of the intersection.
Definition: assemblerutilities.hh:91
GFSU::Traits::SizeType SizeType
Size type as used in grid function space.
Definition: default/assembler.hh:42
Wrap intersection.
Definition: geometrywrapper.hh:56
Wrap element.
Definition: geometrywrapper.hh:15