dune-pdelab  2.5-dev
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_GRIDOPERATOR_COMMON_ASSEMBLERUTILITIES_HH
5 #define DUNE_PDELAB_GRIDOPERATOR_COMMON_ASSEMBLERUTILITIES_HH
6 
7 #include <algorithm>
8 #include <tuple>
9 
12 
13 namespace Dune{
14  namespace PDELab{
15 
21  template<typename GO>
23  {
24 
26  typedef typename GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace;
27 
29  typedef typename GO::Traits::TestGridFunctionSpace TestGridFunctionSpace;
30 
31 
33  typedef typename GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints;
34 
36  typedef typename GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints;
37 
38 
40  typedef typename GO::Traits::MatrixBackend MatrixBackend;
41 
42 
44  typedef typename GO::Traits::DomainField DomainField;
45 
47  typedef typename GO::Traits::Domain Solution;
48 
49 
51  typedef typename GO::Traits::RangeField RangeField;
52 
54  typedef typename GO::Traits::Range Residual;
55 
56 
58  typedef typename GO::Traits::JacobianField JacobianField;
59 
61  typedef typename GO::Traits::Jacobian Jacobian;
62 
64  typedef typename MatrixBackend::template Pattern<
65  Jacobian,
69 
71  typedef typename GO::BorderDOFExchanger BorderDOFExchanger;
72 
73  };
74 
75  // ********************************************************************************
76  // default local pattern implementation
77  // ********************************************************************************
78 
91  class SparsityLink : public std::tuple<int,int>
92  {
93  public:
96  {}
97 
99  SparsityLink (int i, int j)
100  : std::tuple<int,int>(i,j)
101  {}
102 
104  inline int i () const
105  {
106  return std::get<0>(*this);
107  }
108 
110  inline int j () const
111  {
112  return std::get<1>(*this);
113  }
114 
116  void set (int i, int j)
117  {
118  std::get<0>(*this) = i;
119  std::get<1>(*this) = j;
120  }
121  };
122 
130  : public std::vector<SparsityLink>
131  {
132 
133  // make push_back() inaccessible
134  using std::vector<SparsityLink>::push_back;
135 
136  public:
137 
139 
148  template<typename LFSV, typename LFSU>
149  void addLink(const LFSV& lfsv, std::size_t i, const LFSU& lfsu, std::size_t j)
150  {
151  std::vector<SparsityLink>::push_back(
152  SparsityLink(
153  lfsv.localIndex(i),
154  lfsu.localIndex(j)
155  )
156  );
157  }
158 
159  };
160 
161 
162 
163  // ********************************************************************************
164  // Assembler base class
165  // ********************************************************************************
166 
179  template<typename B,
180  typename CU=EmptyTransformation,
181  typename CV=EmptyTransformation>
183  {
184  public:
185 
186  typedef typename B::size_type SizeType;
187 
191  {}
192 
194  LocalAssemblerBase (const CU& cu, const CV& cv)
195  :pconstraintsu(&cu), pconstraintsv(&cv)
196  {}
197 
199  const CU& trialConstraints() const
200  {
201  return *pconstraintsu;
202  }
203 
205  const CV& testConstraints() const
206  {
207  return *pconstraintsv;
208  }
209 
210 
216  template<typename X>
217  typename std::enable_if<
218  AlwaysTrue<X>::value && !std::is_same<
219  CV,
221  >::value
222  >::type
223  forwardtransform(X & x, const bool postrestrict = false) const
224  {
225  typedef typename CV::const_iterator global_col_iterator;
226  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit){
227  typedef typename global_col_iterator::value_type::first_type GlobalIndex;
228  const GlobalIndex & contributor = cit->first;
229 
230  typedef typename global_col_iterator::value_type::second_type ContributedMap;
231  typedef typename ContributedMap::const_iterator global_row_iterator;
232  const ContributedMap & contributed = cit->second;
233  global_row_iterator it = contributed.begin();
234  global_row_iterator eit = contributed.end();
235 
236  for(;it!=eit;++it)
237  {
238  // typename X::block_type block(x[contributor]);
239  // block *= it->second;
240  // x[it->first] += block;
241  x[it->first] += it->second * x[contributor];
242  }
243  }
244 
245  if(postrestrict)
246  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit)
247  x[cit->first]=0.;
248  }
249 
250 
251  // Disable forwardtransform for EmptyTransformation
252  template<typename X>
253  typename std::enable_if<
254  AlwaysTrue<X>::value && std::is_same<
255  CV,
257  >::value
258  >::type
259  forwardtransform(X & x, const bool postrestrict = false) const
260  {}
261 
262 
267  template<typename X>
268  typename std::enable_if<
269  AlwaysTrue<X>::value && !std::is_same<
270  CV,
272  >::value
273  >::type
274  backtransform(X & x, const bool prerestrict = false) const
275  {
276  typedef typename CV::const_iterator global_col_iterator;
277  for (global_col_iterator cit=pconstraintsv->begin(); cit!=pconstraintsv->end(); ++cit){
278  typedef typename global_col_iterator::value_type::first_type GlobalIndex;
279  const GlobalIndex & contributor = cit->first;
280 
281  typedef typename global_col_iterator::value_type::second_type ContributedMap;
282  typedef typename ContributedMap::const_iterator global_row_iterator;
283  const ContributedMap & contributed = cit->second;
284  global_row_iterator it = contributed.begin();
285  global_row_iterator eit = contributed.end();
286 
287  if(prerestrict)
288  x[contributor] = 0.;
289 
290  for(;it!=eit;++it)
291  {
292  // typename X::block_type block(x[it->first]);
293  // block *= it->second;
294  // x[contributor] += block;
295  x[contributor] += it->second * x[it->first]; // PB: 27 Sep 12 this was the old version
296  }
297  }
298  }
299 
300  // disable backtransform for empty transformation
301  template<typename X>
302  typename std::enable_if<
303  AlwaysTrue<X>::value && std::is_same<
304  CV,
306  >::value
307  >::type
308  backtransform(X & x, const bool prerestrict = false) const
309  {}
310 
311 
312  protected:
313 
315  template<typename GCView, typename T>
316  void eread (const GCView& globalcontainer_view, LocalMatrix<T>& localcontainer) const
317  {
318  for (int i = 0; i < localcontainer.N(); ++i)
319  for (int j = 0; j < localcontainer.M(); ++j)
320  localcontainer(i,j) = globalcontainer_view(i,j);
321  }
322 
324  template<typename T, typename GCView>
325  void ewrite (const LocalMatrix<T>& localcontainer, GCView& globalcontainer_view) const
326  {
327  for (int i = 0; i < localcontainer.N(); ++i)
328  for (int j = 0; j < localcontainer.M(); ++j)
329  globalcontainer_view(i,j) = localcontainer(i,j);
330  }
331 
333  template<typename T, typename GCView>
334  void eadd (const LocalMatrix<T>& localcontainer, GCView& globalcontainer_view) const
335  {
336  for (int i = 0; i < localcontainer.N(); ++i)
337  for (int j = 0; j < localcontainer.M(); ++j)
338  globalcontainer_view.add(i,j,localcontainer(i,j));
339  }
340 
342  template<typename M, typename GCView>
343  typename std::enable_if<
344  AlwaysTrue<M>::value && !std::is_same<
345  CV,
347  >::value
348  >::type
349  scatter_jacobian(M& local_container, GCView& global_container_view, bool symmetric_mode) const
350  {
351  typedef typename GCView::RowIndexCache LFSVIndexCache;
352  typedef typename GCView::ColIndexCache LFSUIndexCache;
353 
354  const LFSVIndexCache& lfsv_indices = global_container_view.rowIndexCache();
355  const LFSUIndexCache& lfsu_indices = global_container_view.colIndexCache();
356 
357  if (lfsv_indices.constraintsCachingEnabled() && lfsu_indices.constraintsCachingEnabled())
358  if (symmetric_mode)
359  etadd_symmetric(local_container,global_container_view);
360  else
361  etadd(local_container,global_container_view);
362  else
363  {
364 
365  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
366  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
367 
368  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
369  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
370 
371  // optionally clear out columns that belong to Dirichlet-constrained DOFs to keep matrix symmetric
372  if (symmetric_mode)
373  {
374  typedef typename LFSUIndexCache::ContainerIndex CI;
375 
376  for (size_t j = 0; j < lfsu_indices.size(); ++j)
377  {
378  const CI& container_index = lfsu_indices.containerIndex(j);
379  const typename CU::const_iterator cit = pconstraintsu->find(container_index);
380  if (cit != pconstraintsu->end())
381  {
382  // make sure we only have Dirichlet constraints
383  assert(cit->second.empty());
384  // clear out the current column
385  for (size_t i = 0; i < lfsv_indices.size(); ++i)
386  {
387  // we do not need to update the residual, since the solution
388  // (i.e. the correction) for the Dirichlet DOF is 0 by definition
389  local_container(lfsv,i,lfsu,j) = 0.0;
390  }
391  }
392  }
393  }
394 
395  // write entries without considering constraints.
396  // Dirichlet-constrained rows will be fixed in a postprocessing step.
397  for (auto it = local_container.begin(); it != local_container.end(); ++it)
398  {
399  // skip 0 entries because they might not be present in the pattern
400  if (*it == 0.0)
401  continue;
402  global_container_view.add(it.row(),it.col(),*it);
403  }
404  }
405  }
406 
407  // specialization for empty constraints container
408  template<typename M, typename GCView>
409  typename std::enable_if<
410  AlwaysTrue<M>::value && std::is_same<
411  CV,
413  >::value
414  >::type
415  scatter_jacobian(M& local_container, GCView& global_container_view, bool symmetric_mode) const
416  {
417  // write entries without considering constraints.
418  // Dirichlet-constrained rows will be fixed in a postprocessing step.
419  for (auto it = local_container.begin(); it != local_container.end(); ++it)
420  {
421  // skip 0 entries because they might not be present in the pattern
422  if (*it == 0.0)
423  continue;
424  global_container_view.add(it.row(),it.col(),*it);
425  }
426  }
427 
431  template<typename M, typename GCView>
432  void etadd_symmetric (M& localcontainer, GCView& globalcontainer_view) const
433  {
434 
435  typedef typename GCView::RowIndexCache LFSVIndexCache;
436  typedef typename GCView::ColIndexCache LFSUIndexCache;
437 
438  const LFSVIndexCache& lfsv_indices = globalcontainer_view.rowIndexCache();
439  const LFSUIndexCache& lfsu_indices = globalcontainer_view.colIndexCache();
440 
441  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
442  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
443 
444  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
445  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
446 
447  for (size_t j = 0; j < lfsu_indices.size(); ++j)
448  {
449  if (lfsu_indices.isConstrained(j) && lfsu_indices.isDirichletConstraint(j))
450  {
451  // clear out the current column
452  for (size_t i = 0; i < lfsv_indices.size(); ++i)
453  {
454  // we do not need to update the residual, since the solution
455  // (i.e. the correction) for the Dirichlet DOF is 0 by definition
456  localcontainer(lfsv,i,lfsu,j) = 0.0;
457  }
458  }
459  }
460 
461  // hand off to standard etadd() method
462  etadd(localcontainer,globalcontainer_view);
463  }
464 
465 
466  template<typename M, typename GCView>
467  void etadd (const M& localcontainer, GCView& globalcontainer_view) const
468  {
469 
470  typedef typename M::value_type T;
471 
472  typedef typename GCView::RowIndexCache LFSVIndexCache;
473  typedef typename GCView::ColIndexCache LFSUIndexCache;
474 
475  const LFSVIndexCache& lfsv_indices = globalcontainer_view.rowIndexCache();
476  const LFSUIndexCache& lfsu_indices = globalcontainer_view.colIndexCache();
477 
478  typedef typename LFSVIndexCache::LocalFunctionSpace LFSV;
479  const LFSV& lfsv = lfsv_indices.localFunctionSpace();
480 
481  typedef typename LFSUIndexCache::LocalFunctionSpace LFSU;
482  const LFSU& lfsu = lfsu_indices.localFunctionSpace();
483 
484  for (size_t i = 0; i<lfsv_indices.size(); ++i)
485  for (size_t j = 0; j<lfsu_indices.size(); ++j)
486  {
487 
488  if (localcontainer(lfsv,i,lfsu,j) == 0.0)
489  continue;
490 
491  const bool constrained_v = lfsv_indices.isConstrained(i);
492  const bool constrained_u = lfsu_indices.isConstrained(j);
493 
494  typedef typename LFSVIndexCache::ConstraintsIterator VConstraintsIterator;
495  typedef typename LFSUIndexCache::ConstraintsIterator UConstraintsIterator;
496 
497  if (constrained_v)
498  {
499  if (lfsv_indices.isDirichletConstraint(i))
500  continue;
501 
502  for (VConstraintsIterator vcit = lfsv_indices.constraintsBegin(i); vcit != lfsv_indices.constraintsEnd(i); ++vcit)
503  if (constrained_u)
504  {
505  if (lfsu_indices.isDirichletConstraint(j))
506  {
507  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight();
508  if (value != 0.0)
509  globalcontainer_view.add(vcit->containerIndex(),j,value);
510  }
511  else
512  {
513  for (UConstraintsIterator ucit = lfsu_indices.constraintsBegin(j); ucit != lfsu_indices.constraintsEnd(j); ++ucit)
514  {
515  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight() * ucit->weight();
516  if (value != 0.0)
517  globalcontainer_view.add(vcit->containerIndex(),ucit->containerIndex(),value);
518  }
519  }
520  }
521  else
522  {
523  T value = localcontainer(lfsv,i,lfsu,j) * vcit->weight();
524  if (value != 0.0)
525  globalcontainer_view.add(vcit->containerIndex(),j,value);
526  }
527  }
528  else
529  {
530  if (constrained_u)
531  {
532  if (lfsu_indices.isDirichletConstraint(j))
533  {
534  T value = localcontainer(lfsv,i,lfsu,j);
535  if (value != 0.0)
536  globalcontainer_view.add(i,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) * ucit->weight();
543  if (value != 0.0)
544  globalcontainer_view.add(i,ucit->containerIndex(),value);
545  }
546  }
547  }
548  else
549  globalcontainer_view.add(i,j,localcontainer(lfsv,i,lfsu,j));
550  }
551  }
552  }
553 
554 
555  template<typename Pattern, typename RI, typename CI>
556  typename std::enable_if<
558  >::type
559  add_diagonal_entry(Pattern& pattern, const RI& ri, const CI& ci) const
560  {
561  if (ri == ci)
562  pattern.add_link(ri,ci);
563  }
564 
565  template<typename Pattern, typename RI, typename CI>
566  typename std::enable_if<
568  >::type
569  add_diagonal_entry(Pattern& pattern, const RI& ri, const CI& ci) const
570  {}
571 
576  template<typename P, typename LFSVIndices, typename LFSUIndices, typename Index>
577  void add_entry(P & globalpattern,
578  const LFSVIndices& lfsv_indices, Index i,
579  const LFSUIndices& lfsu_indices, Index j) const
580  {
581  typedef typename LFSVIndices::ConstraintsIterator VConstraintsIterator;
582  typedef typename LFSUIndices::ConstraintsIterator UConstraintsIterator;
583 
584  const bool constrained_v = lfsv_indices.isConstrained(i);
585  const bool constrained_u = lfsu_indices.isConstrained(j);
586 
587  add_diagonal_entry(globalpattern,
588  lfsv_indices.containerIndex(i),
589  lfsu_indices.containerIndex(j)
590  );
591 
592  if(!constrained_v)
593  {
594  if (!constrained_u || lfsu_indices.isDirichletConstraint(j))
595  {
596  globalpattern.add_link(lfsv_indices.containerIndex(i),lfsu_indices.containerIndex(j));
597  }
598  else
599  {
600  for (UConstraintsIterator gurit = lfsu_indices.constraintsBegin(j); gurit != lfsu_indices.constraintsEnd(j); ++gurit)
601  globalpattern.add_link(lfsv_indices.containerIndex(i),gurit->containerIndex());
602  }
603  }
604  else
605  {
606  if (lfsv_indices.isDirichletConstraint(i))
607  {
608  globalpattern.add_link(lfsv_indices.containerIndex(i),lfsu_indices.containerIndex(j));
609  }
610  else
611  {
612  for(VConstraintsIterator gvrit = lfsv_indices.constraintsBegin(i); gvrit != lfsv_indices.constraintsEnd(i); ++gvrit)
613  {
614  if (!constrained_u || lfsu_indices.isDirichletConstraint(j))
615  {
616  globalpattern.add_link(gvrit->containerIndex(),lfsu_indices.containerIndex(j));
617  }
618  else
619  {
620  for (UConstraintsIterator gurit = lfsu_indices.constraintsBegin(j); gurit != lfsu_indices.constraintsEnd(j); ++gurit)
621  globalpattern.add_link(gvrit->containerIndex(),gurit->containerIndex());
622  }
623  }
624  }
625  }
626  }
627 
631  template<typename GFSV, typename GC, typename C>
632  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const C& c) const
633  {
634  typedef typename C::const_iterator global_row_iterator;
635  for (global_row_iterator cit = c.begin(); cit != c.end(); ++cit)
636  globalcontainer.clear_row(cit->first,1);
637  }
638 
639  template<typename GFSV, typename GC>
640  void set_trivial_rows(const GFSV& gfsv, GC& globalcontainer, const EmptyTransformation& c) const
641  {
642  }
643 
644  template<typename GFSV, typename GC>
645  void handle_dirichlet_constraints(const GFSV& gfsv, GC& globalcontainer) const
646  {
647  globalcontainer.flush();
648  set_trivial_rows(gfsv,globalcontainer,*pconstraintsv);
649  globalcontainer.finalize();
650  }
651 
652  /* constraints */
653  const CU* pconstraintsu;
654  const CV* pconstraintsv;
655  static CU emptyconstraintsu;
656  static CV emptyconstraintsv;
657  };
658 
659  template<typename B, typename CU, typename CV>
661  template<typename B, typename CU, typename CV>
663 
664  } // namespace PDELab
665 } // namespace Dune
666 
667 #endif // DUNE_PDELAB_GRIDOPERATOR_COMMON_ASSEMBLERUTILITIES_HH
Dune::PDELab::SparsityLink::j
int j() const
Return second component.
Definition: assemblerutilities.hh:110
Dune::PDELab::LocalAssemblerBase::LocalAssemblerBase
LocalAssemblerBase(const CU &cu, const CV &cv)
construct AssemblerSpace, with constraints
Definition: assemblerutilities.hh:194
Dune::PDELab::LocalAssemblerBase::etadd
void etadd(const M &localcontainer, GCView &globalcontainer_view) const
Definition: assemblerutilities.hh:467
Dune::PDELab::LocalAssemblerTraits::BorderDOFExchanger
GO::BorderDOFExchanger BorderDOFExchanger
The helper class to exchange data on the processor boundary.
Definition: assemblerutilities.hh:71
Dune::PDELab::SparsityLink::i
int i() const
Return first component.
Definition: assemblerutilities.hh:104
Dune::PDELab::LocalAssemblerTraits::DomainField
GO::Traits::DomainField DomainField
The field type of the domain (solution).
Definition: assemblerutilities.hh:44
Dune::PDELab::LocalAssemblerBase::scatter_jacobian
std::enable_if< AlwaysTrue< M >::value &&std::is_same< CV, EmptyTransformation >::value >::type scatter_jacobian(M &local_container, GCView &global_container_view, bool symmetric_mode) const
Definition: assemblerutilities.hh:415
Dune::PDELab::LocalAssemblerBase
Base class for local assembler.
Definition: assemblerutilities.hh:182
Dune::PDELab::LocalAssemblerTraits::Solution
GO::Traits::Domain Solution
The type of the domain (solution).
Definition: assemblerutilities.hh:47
Dune::PDELab::LocalAssemblerTraits::JacobianField
GO::Traits::JacobianField JacobianField
The field type of the jacobian.
Definition: assemblerutilities.hh:58
Dune::PDELab::SparsityLink::SparsityLink
SparsityLink(int i, int j)
Initialize all components.
Definition: assemblerutilities.hh:99
Dune::PDELab::SparsityLink
Entry in sparsity pattern.
Definition: assemblerutilities.hh:91
Dune::PDELab::LocalAssemblerTraits::RangeField
GO::Traits::RangeField RangeField
The field type of the range (residual).
Definition: assemblerutilities.hh:51
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::LocalAssemblerBase::pconstraintsv
const CV * pconstraintsv
Definition: assemblerutilities.hh:654
Dune::PDELab::LocalAssemblerBase::ewrite
void ewrite(const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
write local stiffness matrix for entity
Definition: assemblerutilities.hh:325
Dune::PDELab::LocalAssemblerBase::set_trivial_rows
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:632
Dune::PDELab::LocalAssemblerTraits
Definition: assemblerutilities.hh:22
Dune::PDELab::LocalAssemblerBase::eadd
void eadd(const LocalMatrix< T > &localcontainer, GCView &globalcontainer_view) const
write local stiffness matrix for entity
Definition: assemblerutilities.hh:334
Dune::PDELab::LocalAssemblerTraits::TrialGridFunctionSpaceConstraints
GO::Traits::TrialGridFunctionSpaceConstraints TrialGridFunctionSpaceConstraints
The type of the trial grid function space constraints.
Definition: assemblerutilities.hh:33
Dune::PDELab::LocalAssemblerBase::SizeType
B::size_type SizeType
Definition: assemblerutilities.hh:186
Dune::PDELab::LocalAssemblerTraits::MatrixBackend
GO::Traits::MatrixBackend MatrixBackend
The matrix backend of the grid operator.
Definition: assemblerutilities.hh:40
Dune::PDELab::LocalAssemblerBase::add_entry
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:577
Dune::PDELab::LocalAssemblerBase::backtransform
std::enable_if< AlwaysTrue< X >::value &&std::is_same< CV, EmptyTransformation >::value >::type backtransform(X &x, const bool prerestrict=false) const
Definition: assemblerutilities.hh:308
Dune::PDELab::LocalAssemblerBase::pconstraintsu
const CU * pconstraintsu
Definition: assemblerutilities.hh:653
Dune::PDELab::LocalAssemblerTraits::TrialGridFunctionSpace
GO::Traits::TrialGridFunctionSpace TrialGridFunctionSpace
The trial grid function space.
Definition: assemblerutilities.hh:26
Dune::PDELab::LocalAssemblerBase::LocalAssemblerBase
LocalAssemblerBase()
construct AssemblerSpace
Definition: assemblerutilities.hh:189
Dune::PDELab::LocalAssemblerBase::forwardtransform
std::enable_if< AlwaysTrue< X >::value &&std::is_same< CV, EmptyTransformation >::value >::type forwardtransform(X &x, const bool postrestrict=false) const
Definition: assemblerutilities.hh:259
Dune::PDELab::SparsityLink::SparsityLink
SparsityLink()
Standard constructor for uninitialized local index.
Definition: assemblerutilities.hh:95
Dune::PDELab::EmptyTransformation
Definition: constraintstransformation.hh:111
Dune::PDELab::LocalAssemblerTraits::MatrixPattern
MatrixBackend::template Pattern< Jacobian, TestGridFunctionSpace, TrialGridFunctionSpace > MatrixPattern
The matrix pattern.
Definition: assemblerutilities.hh:68
Dune::PDELab::LocalAssemblerBase::eread
void eread(const GCView &globalcontainer_view, LocalMatrix< T > &localcontainer) const
read local stiffness matrix for entity
Definition: assemblerutilities.hh:316
Dune::PDELab::LocalAssemblerTraits::Residual
GO::Traits::Range Residual
The type of the range (residual).
Definition: assemblerutilities.hh:54
constraintstransformation.hh
Dune::PDELab::SparsityLink::set
void set(int i, int j)
Set both components.
Definition: assemblerutilities.hh:116
Dune::PDELab::LocalAssemblerBase::forwardtransform
std::enable_if< AlwaysTrue< X >::value &&!std::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:223
Dune::PDELab::LocalAssemblerTraits::TestGridFunctionSpaceConstraints
GO::Traits::TestGridFunctionSpaceConstraints TestGridFunctionSpaceConstraints
The type of the test grid function space constraints.
Definition: assemblerutilities.hh:36
Dune::PDELab::LocalAssemblerBase::emptyconstraintsv
static CV emptyconstraintsv
Definition: assemblerutilities.hh:656
Dune::PDELab::LocalMatrix
A dense matrix for storing data associated with the degrees of freedom of a pair of LocalFunctionSpac...
Definition: localmatrix.hh:184
Dune::PDELab::LocalAssemblerTraits::TestGridFunctionSpace
GO::Traits::TestGridFunctionSpace TestGridFunctionSpace
The test grid function space.
Definition: assemblerutilities.hh:29
Dune::PDELab::LocalAssemblerBase::add_diagonal_entry
std::enable_if< !std::is_same< RI, CI >::value >::type add_diagonal_entry(Pattern &pattern, const RI &ri, const CI &ci) const
Definition: assemblerutilities.hh:569
Dune::PDELab::LocalAssemblerBase::emptyconstraintsu
static CU emptyconstraintsu
Definition: assemblerutilities.hh:655
Dune::PDELab::LocalSparsityPattern
Layout description for a sparse linear operator.
Definition: assemblerutilities.hh:129
Dune::PDELab::LocalSparsityPattern::addLink
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:149
Dune::PDELab::LocalAssemblerBase::trialConstraints
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:199
Dune::PDELab::LocalAssemblerBase::backtransform
std::enable_if< AlwaysTrue< X >::value &&!std::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:274
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::LocalAssemblerBase::scatter_jacobian
std::enable_if< AlwaysTrue< M >::value &&!std::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:349
localmatrix.hh
Dune::PDELab::LocalAssemblerBase::handle_dirichlet_constraints
void handle_dirichlet_constraints(const GFSV &gfsv, GC &globalcontainer) const
Definition: assemblerutilities.hh:645
Dune::PDELab::LocalAssemblerBase::etadd_symmetric
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:432
Dune::PDELab::LocalAssemblerBase::set_trivial_rows
void set_trivial_rows(const GFSV &gfsv, GC &globalcontainer, const EmptyTransformation &c) const
Definition: assemblerutilities.hh:640
Dune::PDELab::LocalAssemblerTraits::Jacobian
GO::Traits::Jacobian Jacobian
The type of the jacobian.
Definition: assemblerutilities.hh:61
Dune::PDELab::LocalAssemblerBase::add_diagonal_entry
std::enable_if< std::is_same< RI, CI >::value >::type add_diagonal_entry(Pattern &pattern, const RI &ri, const CI &ci) const
Definition: assemblerutilities.hh:559
Dune::PDELab::LocalAssemblerBase::testConstraints
const CV & testConstraints() const
get the constraints on the test grid function space
Definition: assemblerutilities.hh:205