dune-pdelab  2.0.0
assemblerutilities.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_ASSEMBLERUTILITIES_HH
5 #define DUNE_PDELAB_ASSEMBLERUTILITIES_HH
6 
9 
10 namespace Dune{
11  namespace PDELab{
12 
18  template<typename GO>
20  {
21 
23  typedef typename GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace;
24 
26  typedef typename GO::Traits::TestGridFunctionSpace TestGridFunctionSpace;
27 
28 
30  typedef typename GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints;
31 
33  typedef typename GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints;
34 
35 
37  typedef typename GO::Traits::MatrixBackend MatrixBackend;
38 
39 
41  typedef typename GO::Traits::DomainField DomainField;
42 
44  typedef typename GO::Traits::Domain Solution;
45 
46 
48  typedef typename GO::Traits::RangeField RangeField;
49 
51  typedef typename GO::Traits::Range Residual;
52 
53 
55  typedef typename GO::Traits::JacobianField JacobianField;
56 
58  typedef typename GO::Traits::Jacobian Jacobian;
59 
61  typedef typename MatrixBackend::template Pattern<
62  Jacobian,
66 
68  typedef typename GO::BorderDOFExchanger BorderDOFExchanger;
69 
70  };
71 
73 
80  {
81 
82  enum Type {
83  processor = 0,
84  skeleton = 1,
85  boundary = 2,
86  periodic = 3
87  };
88 
90  template<typename Intersection>
91  static Type get(const Intersection& is)
92  {
93  return static_cast<Type>(1*is.neighbor() + 2*is.boundary());
94  }
95 
96  };
97 
98 
99 
100  // ********************************************************************************
101  // default local pattern implementation
102  // ********************************************************************************
103 
116  class SparsityLink : public Dune::tuple<int,int>
117  {
118  public:
121  {}
122 
124  SparsityLink (int i, int j)
125  : Dune::tuple<int,int>(i,j)
126  {}
127 
129  inline int i () const
130  {
131  return Dune::get<0>(*this);
132  }
133 
135  inline int j () const
136  {
137  return Dune::get<1>(*this);
138  }
139 
141  void set (int i, int j)
142  {
143  Dune::get<0>(*this) = i;
144  Dune::get<1>(*this) = j;
145  }
146  };
147 
155  : public std::vector<SparsityLink>
156  {
157 
158  public:
159 
160  void push_back(const SparsityLink& link)
161  DUNE_DEPRECATED_MSG("The std::vector-like interface to LocalSparsityPattern is deprecated, use addLink() instead.")
162  {
163  std::vector<SparsityLink>::push_back(link);
164  }
165 
167 
176  template<typename LFSV, typename LFSU>
177  void addLink(const LFSV& lfsv, std::size_t i, const LFSU& lfsu, std::size_t j)
178  {
179  std::vector<SparsityLink>::push_back(
180  SparsityLink(
181  lfsv.localIndex(i),
182  lfsu.localIndex(j)
183  )
184  );
185  }
186 
187  };
188 
189 
190 
191  // ********************************************************************************
192  // Assembler base class
193  // ********************************************************************************
194 
207  template<typename B,
208  typename CU=EmptyTransformation,
209  typename CV=EmptyTransformation>
211  public:
212 
213  typedef typename B::size_type SizeType;
214 
218  {}
219 
221  LocalAssemblerBase (const CU& cu, const CV& cv)
222  :pconstraintsu(&cu), pconstraintsv(&cv)
223  {}
224 
226  const CU& trialConstraints() const
227  {
228  return *pconstraintsu;
229  }
230 
232  const CV& testConstraints() const
233  {
234  return *pconstraintsv;
235  }
236 
237 
243  template<typename X>
244  typename enable_if<
245  AlwaysTrue<X>::value && !is_same<
246  CV,
248  >::value
249  >::type
250  forwardtransform(X & x, const bool postrestrict = false) const
251  {
252  typedef typename CV::const_iterator global_col_iterator;
253  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit){
254  typedef typename global_col_iterator::value_type::first_type GlobalIndex;
255  const GlobalIndex & contributor = cit->first;
256 
257  typedef typename global_col_iterator::value_type::second_type ContributedMap;
258  typedef typename ContributedMap::const_iterator global_row_iterator;
259  const ContributedMap & contributed = cit->second;
260  global_row_iterator it = contributed.begin();
261  global_row_iterator eit = contributed.end();
262 
263  for(;it!=eit;++it)
264  {
265  // typename X::block_type block(x[contributor]);
266  // block *= it->second;
267  // x[it->first] += block;
268  x[it->first] += it->second * x[contributor];
269  }
270  }
271 
272  if(postrestrict)
273  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit)
274  x[cit->first]=0.;
275  }
276 
277 
278  // Disable forwardtransform for EmptyTransformation
279  template<typename X>
280  typename enable_if<
281  AlwaysTrue<X>::value && is_same<
282  CV,
284  >::value
285  >::type
286  forwardtransform(X & x, const bool postrestrict = false) const
287  {}
288 
289 
294  template<typename X>
295  typename enable_if<
296  AlwaysTrue<X>::value && !is_same<
297  CV,
299  >::value
300  >::type
301  backtransform(X & x, const bool prerestrict = false) const
302  {
303  typedef typename CV::const_iterator global_col_iterator;
304  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit){
305  typedef typename global_col_iterator::value_type::first_type GlobalIndex;
306  const GlobalIndex & contributor = cit->first;
307 
308  typedef typename global_col_iterator::value_type::second_type ContributedMap;
309  typedef typename ContributedMap::const_iterator global_row_iterator;
310  const ContributedMap & contributed = cit->second;
311  global_row_iterator it = contributed.begin();
312  global_row_iterator eit = contributed.end();
313 
314  if(prerestrict)
315  x[contributor] = 0.;
316 
317  for(;it!=eit;++it)
318  {
319  // typename X::block_type block(x[it->first]);
320  // block *= it->second;
321  // x[contributor] += block;
322  x[contributor] += it->second * x[it->first]; // PB: 27 Sep 12 this was the old version
323  }
324  }
325  }
326 
327  // disable backtransform for empty transformation
328  template<typename X>
329  typename enable_if<
330  AlwaysTrue<X>::value && is_same<
331  CV,
333  >::value
334  >::type
335  backtransform(X & x, const bool prerestrict = false) const
336  {}
337 
338 
339  protected:
340 
342  template<typename GCView, typename T>
343  void eread (const GCView& globalcontainer_view, LocalMatrix<T>& localcontainer) const
344  {
345  for (int i = 0; i < localcontainer.N(); ++i)
346  for (int j = 0; j < localcontainer.M(); ++j)
347  localcontainer(i,j) = globalcontainer_view(i,j);
348  }
349 
351  template<typename T, typename GCView>
352  void ewrite (const LocalMatrix<T>& localcontainer, GCView& globalcontainer_view) const
353  {
354  for (int i = 0; i < localcontainer.N(); ++i)
355  for (int j = 0; j < localcontainer.M(); ++j)
356  globalcontainer_view(i,j) = localcontainer(i,j);
357  }
358 
360  template<typename T, typename GCView>
361  void eadd (const LocalMatrix<T>& localcontainer, GCView& globalcontainer_view) const
362  {
363  for (int i = 0; i < localcontainer.N(); ++i)
364  for (int j = 0; j < localcontainer.M(); ++j)
365  globalcontainer_view.add(i,j,localcontainer(i,j));
366  }
367 
369  template<typename M, typename GCView>
370  typename enable_if<
371  AlwaysTrue<M>::value && !is_same<
372  CV,
374  >::value
375  >::type
376  scatter_jacobian(M& local_container, GCView& global_container_view, bool symmetric_mode) const
377  {
378  typedef typename GCView::RowIndexCache LFSVIndexCache;
379  typedef typename GCView::ColIndexCache LFSUIndexCache;
380 
381  const LFSVIndexCache& lfsv_indices = global_container_view.rowIndexCache();
382  const LFSUIndexCache& lfsu_indices = global_container_view.colIndexCache();
383 
384  if (lfsv_indices.constraintsCachingEnabled() && lfsu_indices.constraintsCachingEnabled())
385  if (symmetric_mode)
386  etadd_symmetric(local_container,global_container_view);
387  else
388  etadd(local_container,global_container_view);
389  else
390  {
391 
392  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
393  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
394 
395  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
396  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
397 
398  // optionally clear out columns that belong to Dirichlet-constrained DOFs to keep matrix symmetric
399  if (symmetric_mode)
400  {
401  typedef typename LFSUIndexCache::ContainerIndex CI;
402 
403  for (size_t j = 0; j < lfsu_indices.size(); ++j)
404  {
405  const CI& container_index = lfsu_indices.containerIndex(j);
406  const typename CU::const_iterator cit = pconstraintsu->find(container_index);
407  if (cit != pconstraintsu->end())
408  {
409  // make sure we only have Dirichlet constraints
410  assert(cit->second.empty());
411  // clear out the current column
412  for (size_t i = 0; i < lfsv_indices.size(); ++i)
413  {
414  // we do not need to update the residual, since the solution
415  // (i.e. the correction) for the Dirichlet DOF is 0 by definition
416  local_container(lfsv,i,lfsu,j) = 0.0;
417  }
418  }
419  }
420  }
421 
422  // write entries without considering constraints.
423  // Dirichlet-constrained rows will be fixed in a postprocessing step.
424  for (auto it = local_container.begin(); it != local_container.end(); ++it)
425  {
426  // skip 0 entries because they might not be present in the pattern
427  if (*it == 0.0)
428  continue;
429  global_container_view.add(it.row(),it.col(),*it);
430  }
431  }
432  }
433 
434  // specialization for empty constraints container
435  template<typename M, typename GCView>
436  typename enable_if<
437  AlwaysTrue<M>::value && is_same<
438  CV,
440  >::value
441  >::type
442  scatter_jacobian(M& local_container, GCView& global_container_view, bool symmetric_mode) const
443  {
444  // write entries without considering constraints.
445  // Dirichlet-constrained rows will be fixed in a postprocessing step.
446  for (auto it = local_container.begin(); it != local_container.end(); ++it)
447  {
448  // skip 0 entries because they might not be present in the pattern
449  if (*it == 0.0)
450  continue;
451  global_container_view.add(it.row(),it.col(),*it);
452  }
453  }
454 
458  template<typename M, typename GCView>
459  void etadd_symmetric (M& localcontainer, GCView& globalcontainer_view) const
460  {
461 
462  typedef typename GCView::RowIndexCache LFSVIndexCache;
463  typedef typename GCView::ColIndexCache LFSUIndexCache;
464 
465  const LFSVIndexCache& lfsv_indices = globalcontainer_view.rowIndexCache();
466  const LFSUIndexCache& lfsu_indices = globalcontainer_view.colIndexCache();
467 
468  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
469  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
470 
471  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
472  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
473 
474  for (size_t j = 0; j < lfsu_indices.size(); ++j)
475  {
476  if (lfsu_indices.isConstrained(j) && lfsu_indices.isDirichletConstraint(j))
477  {
478  // clear out the current column
479  for (size_t i = 0; i < lfsv_indices.size(); ++i)
480  {
481  // we do not need to update the residual, since the solution
482  // (i.e. the correction) for the Dirichlet DOF is 0 by definition
483  localcontainer(lfsv,i,lfsu,j) = 0.0;
484  }
485  }
486  }
487 
488  // hand off to standard etadd() method
489  etadd(localcontainer,globalcontainer_view);
490  }
491 
492 
493  template<typename M, typename GCView>
494  void etadd (const M& localcontainer, GCView& globalcontainer_view) const
495  {
496 
497  typedef typename M::value_type T;
498 
499  typedef typename GCView::RowIndexCache LFSVIndexCache;
500  typedef typename GCView::ColIndexCache LFSUIndexCache;
501 
502  const LFSVIndexCache& lfsv_indices = globalcontainer_view.rowIndexCache();
503  const LFSUIndexCache& lfsu_indices = globalcontainer_view.colIndexCache();
504 
505  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
506  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
507 
508  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
509  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
510 
511  for (size_t i = 0; i<lfsv_indices.size(); ++i)
512  for (size_t j = 0; j<lfsu_indices.size(); ++j)
513  {
514 
515  if (localcontainer(lfsv,i,lfsu,j) == 0.0)
516  continue;
517 
518  const bool constrained_v = lfsv_indices.isConstrained(i);
519  const bool constrained_u = lfsu_indices.isConstrained(j);
520 
521  typedef typename LFSVIndexCache::ConstraintsIterator VConstraintsIterator;
522  typedef typename LFSUIndexCache::ConstraintsIterator UConstraintsIterator;
523 
524  if (constrained_v)
525  {
526  if (lfsv_indices.isDirichletConstraint(i))
527  continue;
528 
529  for (VConstraintsIterator vcit = lfsv_indices.constraintsBegin(i); vcit != lfsv_indices.constraintsEnd(i); ++vcit)
530  if (constrained_u)
531  {
532  if (lfsu_indices.isDirichletConstraint(j))
533  {
534  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight();
535  if (value != 0.0)
536  globalcontainer_view.add(vcit->containerIndex(),j,value);
537  }
538  else
539  {
540  for (UConstraintsIterator ucit = lfsu_indices.constraintsBegin(j); ucit != lfsu_indices.constraintsEnd(j); ++ucit)
541  {
542  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight() * ucit->weight();
543  if (value != 0.0)
544  globalcontainer_view.add(vcit->containerIndex(),ucit->containerIndex(),value);
545  }
546  }
547  }
548  else
549  {
550  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight();
551  if (value != 0.0)
552  globalcontainer_view.add(vcit->containerIndex(),j,value);
553  }
554  }
555  else
556  {
557  if (constrained_u)
558  {
559  if (lfsu_indices.isDirichletConstraint(j))
560  {
561  T value = localcontainer(lfsv,i,lfsu,j);
562  if (value != 0.0)
563  globalcontainer_view.add(i,j,value);
564  }
565  else
566  {
567  for (UConstraintsIterator ucit = lfsu_indices.constraintsBegin(j); ucit != lfsu_indices.constraintsEnd(j); ++ucit)
568  {
569  T value = localcontainer(lfsv,i,lfsu,j) * ucit->weight();
570  if (value != 0.0)
571  globalcontainer_view.add(i,ucit->containerIndex(),value);
572  }
573  }
574  }
575  else
576  globalcontainer_view.add(i,j,localcontainer(lfsv,i,lfsu,j));
577  }
578  }
579  }
580 
581 
582  template<typename Pattern, typename RI, typename CI>
583  typename enable_if<
585  >::type
586  add_diagonal_entry(Pattern& pattern, const RI& ri, const CI& ci) const
587  {
588  if (ri == ci)
589  pattern.add_link(ri,ci);
590  }
591 
592  template<typename Pattern, typename RI, typename CI>
593  typename enable_if<
595  >::type
596  add_diagonal_entry(Pattern& pattern, const RI& ri, const CI& ci) const
597  {}
598 
603  template<typename P, typename LFSVIndices, typename LFSUIndices, typename Index>
604  void add_entry(P & globalpattern,
605  const LFSVIndices& lfsv_indices, Index i,
606  const LFSUIndices& lfsu_indices, Index j) const
607  {
608  typedef typename LFSVIndices::ConstraintsIterator VConstraintsIterator;
609  typedef typename LFSUIndices::ConstraintsIterator UConstraintsIterator;
610 
611  const bool constrained_v = lfsv_indices.isConstrained(i);
612  const bool constrained_u = lfsu_indices.isConstrained(j);
613 
614  add_diagonal_entry(globalpattern,
615  lfsv_indices.containerIndex(i),
616  lfsu_indices.containerIndex(j)
617  );
618 
619  if(!constrained_v)
620  {
621  if (!constrained_u || lfsu_indices.isDirichletConstraint(j))
622  {
623  globalpattern.add_link(lfsv_indices.containerIndex(i),lfsu_indices.containerIndex(j));
624  }
625  else
626  {
627  for (UConstraintsIterator gurit = lfsu_indices.constraintsBegin(j); gurit != lfsu_indices.constraintsEnd(j); ++gurit)
628  globalpattern.add_link(lfsv_indices.containerIndex(i),gurit->containerIndex());
629  }
630  }
631  else
632  {
633  if (lfsv_indices.isDirichletConstraint(i))
634  {
635  globalpattern.add_link(lfsv_indices.containerIndex(i),lfsu_indices.containerIndex(j));
636  }
637  else
638  {
639  for(VConstraintsIterator gvrit = lfsv_indices.constraintsBegin(i); gvrit != lfsv_indices.constraintsEnd(i); ++gvrit)
640  {
641  if (!constrained_u || lfsu_indices.isDirichletConstraint(j))
642  {
643  globalpattern.add_link(gvrit->containerIndex(),lfsu_indices.containerIndex(j));
644  }
645  else
646  {
647  for (UConstraintsIterator gurit = lfsu_indices.constraintsBegin(j); gurit != lfsu_indices.constraintsEnd(j); ++gurit)
648  globalpattern.add_link(gvrit->containerIndex(),gurit->containerIndex());
649  }
650  }
651  }
652  }
653  }
654 
658  template<typename GFSV, typename GC, typename C>
659  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const C& c) const
660  {
661  typedef typename C::const_iterator global_row_iterator;
662  for (global_row_iterator cit = c.begin(); cit != c.end(); ++cit)
663  globalcontainer.clear_row(cit->first,1);
664  }
665 
666  template<typename GFSV, typename GC>
667  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const EmptyTransformation& c) const
668  {
669  }
670 
671  template<typename GFSV, typename GC>
672  void handle_dirichlet_constraints(const GFSV& gfsv, GC& globalcontainer) const
673  {
674  globalcontainer.flush();
675  set_trivial_rows(gfsv,globalcontainer,*pconstraintsv);
676  globalcontainer.finalize();
677  }
678 
679  /* constraints */
680  const CU* pconstraintsu;
681  const CV* pconstraintsv;
682  static CU emptyconstraintsu;
683  static CV emptyconstraintsv;
684  };
685 
686  template<typename B, typename CU, typename CV>
688  template<typename B, typename CU, typename CV>
690 
691  } // namespace PDELab
692 } // namespace Dune
693 
694 #endif //DUNE_PDELAB_ASSEMBLERUTILITIES_HH
const CU * pconstraintsu
Definition: assemblerutilities.hh:680
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:51
void etadd(const M &localcontainer, GCView &globalcontainer_view) const
Definition: assemblerutilities.hh:494
Base class for local assembler.
Definition: assemblerutilities.hh:210
SparsityLink()
Standard constructor for uninitialized local index.
Definition: assemblerutilities.hh:120
enable_if< AlwaysTrue< X >::value &&is_same< CV, EmptyTransformation >::value >::type backtransform(X &x, const bool prerestrict=false) const
Definition: assemblerutilities.hh:335
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const EmptyTransformation &c) const
Definition: assemblerutilities.hh:667
int i() const
Return first component.
Definition: assemblerutilities.hh:129
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:44
const CV & testConstraints() const
get the constraints on the test grid function space
Definition: assemblerutilities.hh:232
domain boundary intersection (neighbor() == false && boundary() == true)
Definition: assemblerutilities.hh:85
GO::Traits::RangeField RangeField
The field type of the range (residual).
Definition: assemblerutilities.hh:48
Translation helper from intersection method return values to intersection type.
Definition: assemblerutilities.hh:79
Definition: common/constraintstransformation.hh:127
GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace
The trial grid function space.
Definition: assemblerutilities.hh:23
enable_if< AlwaysTrue< M >::value &&is_same< CV, EmptyTransformation >::value >::type scatter_jacobian(M &local_container, GCView &global_container_view, bool symmetric_mode) const
Definition: assemblerutilities.hh:442
SparsityLink(int i, int j)
Initialize all components.
Definition: assemblerutilities.hh:124
void set(int i, int j)
Set both components.
Definition: assemblerutilities.hh:141
enable_if< AlwaysTrue< X >::value &&is_same< CV, EmptyTransformation >::value >::type forwardtransform(X &x, const bool postrestrict=false) const
Definition: assemblerutilities.hh:286
GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints
The type of the test grid function space constraints.
Definition: assemblerutilities.hh:33
GO::Traits::DomainField DomainField
The field type of the domain (solution).
Definition: assemblerutilities.hh:41
GO::BorderDOFExchanger BorderDOFExchanger
The helper class to exchange data on the processor boundary.
Definition: assemblerutilities.hh:68
const CV * pconstraintsv
Definition: assemblerutilities.hh:681
skeleton intersection (neighbor() == true && boundary() == false)
Definition: assemblerutilities.hh:84
void push_back(const SparsityLink &link)
Definition: assemblerutilities.hh:160
A dense matrix for storing data associated with the degrees of freedom of a pair of LocalFunctionSpac...
Definition: common/localmatrix.hh:184
void add_entry(P &globalpattern, const LFSVIndices &lfsv_indices, Index i, const LFSUIndices &lfsu_indices, Index j) const
Adding matrix entry to pattern with respect to the constraints contributions. This assembles the entr...
Definition: assemblerutilities.hh:604
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:30
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:65
void etadd_symmetric(M &localcontainer, GCView &globalcontainer_view) const
Add local matrix to global matrix, and apply Dirichlet constraints in a symmetric fashion...
Definition: assemblerutilities.hh:459
enable_if< AlwaysTrue< X >::value &&!is_same< CV, EmptyTransformation >::value >::type backtransform(X &x, const bool prerestrict=false) const
Transforms a vector from to . If prerestrict == true then is applied instead of the full transform...
Definition: assemblerutilities.hh:301
enable_if< is_same< RI, CI >::value >::type add_diagonal_entry(Pattern &pattern, const RI &ri, const CI &ci) const
Definition: assemblerutilities.hh:586
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:226
enable_if< AlwaysTrue< M >::value &&!is_same< CV, EmptyTransformation >::value >::type scatter_jacobian(M &local_container, GCView &global_container_view, bool symmetric_mode) const
Scatter local jacobian to global container.
Definition: assemblerutilities.hh:376
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const C &c) const
insert dirichlet constraints for row and assemble T^T_U in constrained rows
Definition: assemblerutilities.hh:659
LocalAssemblerBase(const CU &cu, const CV &cv)
construct AssemblerSpace, with constraints
Definition: assemblerutilities.hh:221
void handle_dirichlet_constraints(const GFSV &gfsv, GC &globalcontainer) const
Definition: assemblerutilities.hh:672
GO::Traits::TestGridFunctionSpace TestGridFunctionSpace
The test grid function space.
Definition: assemblerutilities.hh:26
Layout description for a sparse linear operator.
Definition: assemblerutilities.hh:154
periodic boundary intersection (neighbor() == true && boundary() == true)
Definition: assemblerutilities.hh:86
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
Type
Definition: assemblerutilities.hh:82
GO::Traits::JacobianField JacobianField
The field type of the jacobian.
Definition: assemblerutilities.hh:55
GO::Traits::MatrixBackend MatrixBackend
The matrix backend of the grid operator.
Definition: assemblerutilities.hh:37
Entry in sparsity pattern.
Definition: assemblerutilities.hh:116
processor boundary intersection (neighbor() == false && boundary() == false)
Definition: assemblerutilities.hh:83
void addLink(const LFSV &lfsv, std::size_t i, const LFSU &lfsu, std::size_t j)
Adds a link between DOF i of lfsv and DOF j of lfsu.
Definition: assemblerutilities.hh:177
void eread(const GCView &globalcontainer_view, LocalMatrix< T > &localcontainer) const
read local stiffness matrix for entity
Definition: assemblerutilities.hh:343
enable_if< !is_same< RI, CI >::value >::type add_diagonal_entry(Pattern &pattern, const RI &ri, const CI &ci) const
Definition: assemblerutilities.hh:596
void ewrite(const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
write local stiffness matrix for entity
Definition: assemblerutilities.hh:352
static CV emptyconstraintsv
Definition: assemblerutilities.hh:683
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:58
static CU emptyconstraintsu
Definition: assemblerutilities.hh:682
Definition: assemblerutilities.hh:19
void eadd(const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
write local stiffness matrix for entity
Definition: assemblerutilities.hh:361
int j() const
Return second component.
Definition: assemblerutilities.hh:135
enable_if< AlwaysTrue< X >::value &&!is_same< CV, EmptyTransformation >::value >::type forwardtransform(X &x, const bool postrestrict=false) const
Transforms a vector from to . If postrestrict == true then is applied instead of the full transfor...
Definition: assemblerutilities.hh:250
LocalAssemblerBase()
construct AssemblerSpace
Definition: assemblerutilities.hh:216
B::size_type SizeType
Definition: assemblerutilities.hh:213