dune-pdelab  2.7-git
gridoperator.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_GRIDOPERATOR_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_GRIDOPERATOR_HH
3 
4 #include <tuple>
5 
6 #include <dune/common/hybridutilities.hh>
7 
13 
14 namespace Dune{
15  namespace PDELab{
16 
30  template<typename GFSU, typename GFSV, typename LOP,
31  typename MB, typename DF, typename RF, typename JF,
34  >
36  {
37  public:
38 
41 
48 
50  typedef typename MB::template Pattern<Jacobian,GFSV,GFSU> Pattern;
51 
53  typedef DefaultLocalAssembler<
55  LOP,
56  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition
58 
59  // Fix this as soon as the default Partitions are constexpr
60  typedef typename std::conditional<
61  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition,
65 
68  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,Assembler,LocalAssembler> Traits;
69 
70  template <typename MFT>
72  typedef typename Traits::Jacobian Type;
73  };
74 
76  GridOperator(const GFSU & gfsu_, const CU & cu_, const GFSV & gfsv_, const CV & cv_, LOP & lop_, const MB& mb_ = MB())
77  : global_assembler(gfsu_,gfsv_,cu_,cv_)
78  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
79  , local_assembler(lop_, cu_, cv_,dof_exchanger)
80  , backend(mb_)
81  {}
82 
84  GridOperator(const GFSU & gfsu_, const GFSV & gfsv_, LOP & lop_, const MB& mb_ = MB())
85  : global_assembler(gfsu_,gfsv_)
86  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
87  , local_assembler(lop_,dof_exchanger)
88  , backend(mb_)
89  {}
90 
92  const GFSU& trialGridFunctionSpace() const
93  {
94  return global_assembler.trialGridFunctionSpace();
95  }
96 
98  const GFSV& testGridFunctionSpace() const
99  {
100  return global_assembler.testGridFunctionSpace();
101  }
102 
104  typename GFSU::Traits::SizeType globalSizeU () const
105  {
106  return trialGridFunctionSpace().globalSize();
107  }
108 
110  typename GFSV::Traits::SizeType globalSizeV () const
111  {
112  return testGridFunctionSpace().globalSize();
113  }
114 
115  Assembler & assembler() { return global_assembler; }
116 
117  const Assembler & assembler() const { return global_assembler; }
118 
119  LocalAssembler & localAssembler() const { return local_assembler; }
120 
121 
124  template <typename GridOperatorTuple>
126  {
128  : index(0), size(std::tuple_size<GridOperatorTuple>::value) {}
129 
130  template <typename T>
131  void visit(T& elem) {
132  elem.localAssembler().preProcessing(index == 0);
133  elem.localAssembler().postProcessing(index == size-1);
134  ++index;
135  }
136 
137  int index;
138  const int size;
139  };
140 
144  template<typename GridOperatorTuple>
145  static void setupGridOperators(GridOperatorTuple tuple)
146  {
147  SetupGridOperator<GridOperatorTuple> setup_visitor;
148  Hybrid::forEach(tuple,
149  [&](auto &el) { setup_visitor.visit(el); });
150  }
151 
153  template<typename F, typename X>
154  void interpolate (const X& xold, F& f, X& x) const
155  {
156  // Interpolate f into grid function space and set corresponding coefficients
157  Dune::PDELab::interpolate(f,global_assembler.trialGridFunctionSpace(),x);
158 
159  // Copy non-constrained dofs from old time step
161  }
162 
164  void fill_pattern(Pattern & p) const
165  {
166  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
167  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
168  global_assembler.assemble(pattern_engine);
169  }
170 
172  void residual(const Domain & x, Range & r) const
173  {
174  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
175  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
176  global_assembler.assemble(residual_engine);
177  }
178 
180  void jacobian(const Domain & x, Jacobian & a) const
181  {
182  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
183  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
184  global_assembler.assemble(jacobian_engine);
185  }
186 
188  void jacobian_apply(const Domain & update, Range & result) const
189  {
190  if (not local_assembler.localOperator().isLinear)
191  DUNE_THROW(Dune::Exception, "Your trying to use a non linear jacobian apply for a linear problem.");
192  global_assembler.assemble(local_assembler.localJacobianApplyAssemblerEngine(update, result));
193  }
194 
196  void jacobian_apply(const Domain & solution, const Domain & update, Range & result) const
197  {
198  if (local_assembler.localOperator().isLinear)
199  DUNE_THROW(Dune::Exception, "Your trying to use a non linear jacobian apply for a linear problem.");
200  global_assembler.assemble(local_assembler.localJacobianApplyAssemblerEngine(solution, update, result));
201  }
202 
204  void DUNE_DEPRECATED_MSG("nonlinear_jacobian_apply(x,z,r) is deprecated. Please use jacobian_apply(solution, update, result) instead!")
205  nonlinear_jacobian_apply(const Domain & solution, const Domain & update, Range & result) const
206  {
207  if (local_assembler.localOperator().isLinear)
208  DUNE_THROW(Dune::Exception, "Your trying to use a non linear jacobian apply for a linear problem.");
209  global_assembler.assemble(local_assembler.localJacobianApplyAssemblerEngine(solution, update, result));
210  }
211 
212  void make_consistent(Jacobian& a) const
213  {
214  dof_exchanger->accumulateBorderEntries(*this,a);
215  }
216 
217  void update()
218  {
219  // the DOF exchanger has matrix information, so we need to update it
220  dof_exchanger->update(*this);
221  }
222 
224  const typename Traits::MatrixBackend& matrixBackend() const
225  {
226  return backend;
227  }
228 
229  private:
230  Assembler global_assembler;
231  std::shared_ptr<BorderDOFExchanger> dof_exchanger;
232 
233  mutable LocalAssembler local_assembler;
234  MB backend;
235 
236  };
237 
238  }
239 }
240 #endif // DUNE_PDELAB_GRIDOPERATOR_GRIDOPERATOR_HH
Dune::PDELab::GridOperator
Standard grid operator implementation.
Definition: gridoperator.hh:35
Dune::PDELab::GridOperator::result
void const Domain Range &const result
Definition: gridoperator.hh:206
localassembler.hh
Dune::PDELab::Backend::Vector
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
Dune::PDELab::GridOperator::update
void update()
Definition: gridoperator.hh:217
gridoperatorutilities.hh
Dune::PDELab::DefaultLocalPatternAssemblerEngine< DefaultLocalAssembler >
p
const P & p
Definition: constraints.hh:148
Dune::PDELab::GridOperator::residual
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator.hh:172
Dune::PDELab::GridOperator::Assembler
DefaultAssembler< GFSU, GFSV, CU, CV > Assembler
The global assembler type.
Definition: gridoperator.hh:40
Dune::PDELab::GridOperator::GridOperator
GridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: gridoperator.hh:84
Dune::PDELab::copy_nonconstrained_dofs
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:987
borderdofexchanger.hh
Dune::PDELab::GridOperator::SetupGridOperator
Definition: gridoperator.hh:125
Dune::PDELab::GridOperator::jacobian_apply
void jacobian_apply(const Domain &solution, const Domain &update, Range &result) const
Apply jacobian matrix to the vector update without explicitly assembling it.
Definition: gridoperator.hh:196
Dune::PDELab::GridOperatorTraits::MatrixBackend
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
Dune::PDELab::GridOperator::fill_pattern
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:164
Dune::PDELab::GridOperator::DUNE_DEPRECATED_MSG
void DUNE_DEPRECATED_MSG("nonlinear_jacobian_apply(x,z,r) is deprecated. Please use jacobian_apply(solution, update, result) instead!") nonlinear_jacobian_apply(const Domain &solution
Apply jacobian matrix to the vector update without explicitly assembling it.
Dune::PDELab::DefaultLocalAssembler::localResidualAssemblerEngine
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: default/localassembler.hh:159
Dune::PDELab::GridOperator< CGGFS, GFS, CGTODGLOP, MBE, field_type, field_type, field_type, CC, CC >::Jacobian
Dune::PDELab::Backend::Matrix< MBE, Domain, Range, field_type > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:47
Dune::PDELab::GridOperator::GridOperator
GridOperator(const GFSU &gfsu_, const CU &cu_, const GFSV &gfsv_, const CV &cv_, LOP &lop_, const MB &mb_=MB())
Constructor for non trivial constraints.
Definition: gridoperator.hh:76
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::GridOperator::Traits
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator.hh:68
Dune::PDELab::GridOperator::BorderDOFExchanger
std::conditional< GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition, NonOverlappingBorderDOFExchanger< GridOperator >, OverlappingBorderDOFExchanger< GridOperator > >::type BorderDOFExchanger
Definition: gridoperator.hh:64
Dune::PDELab::GridOperator::MatrixContainer::Type
Traits::Jacobian Type
Definition: gridoperator.hh:72
Dune::PDELab::GridOperatorTraits
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
Dune::PDELab::DefaultLocalAssembler::localOperator
LOP & localOperator()
get a reference to the local operator
Definition: default/localassembler.hh:100
Dune::PDELab::GridOperator::SetupGridOperator::index
int index
Definition: gridoperator.hh:137
Dune::PDELab::GridOperator::Range
Dune::PDELab::Backend::Vector< GFSV, RF > Range
The type of the range (residual).
Definition: gridoperator.hh:45
Dune::PDELab::DefaultAssembler::trialGridFunctionSpace
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: default/assembler.hh:67
Dune::PDELab::GridOperator::jacobian
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:180
Dune::PDELab::GridOperator::testGridFunctionSpace
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:98
Dune::PDELab::DefaultLocalAssembler::localPatternAssemblerEngine
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: default/localassembler.hh:150
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::DefaultLocalAssembler::localJacobianApplyAssemblerEngine
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &update, typename Traits::Range &result)
Definition: default/localassembler.hh:179
Dune::PDELab::GridOperator< CGGFS, GFS, CGTODGLOP, MBE, field_type, field_type, field_type, CC, CC >::Pattern
MBE ::template Pattern< Jacobian, GFS, CGGFS > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator.hh:50
Dune::PDELab::GridOperatorTraits::Jacobian
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
Dune::PDELab::LocalAssemblerBase::trialConstraints
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:205
Dune::PDELab::GridOperator::SetupGridOperator::size
const int size
Definition: gridoperator.hh:138
Dune::PDELab::DefaultLocalResidualAssemblerEngine< DefaultLocalAssembler >
Dune::PDELab::interpolate
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:177
assembler.hh
Dune::PDELab::GridOperator::setupGridOperators
static void setupGridOperators(GridOperatorTuple tuple)
Definition: gridoperator.hh:145
Dune::PDELab::Backend::Matrix
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:127
Dune::PDELab::GridOperator::SetupGridOperator::SetupGridOperator
SetupGridOperator()
Definition: gridoperator.hh:127
Dune::PDELab::GridOperator::jacobian_apply
void jacobian_apply(const Domain &update, Range &result) const
Apply jacobian matrix to the vector update without explicitly assembling it.
Definition: gridoperator.hh:188
Dune::PDELab::DefaultLocalAssembler
The local assembler for DUNE grids.
Definition: default/localassembler.hh:31
Dune::PDELab::DefaultAssembler
The assembler for standard DUNE grid.
Definition: default/assembler.hh:22
Dune::PDELab::GridOperator::assembler
const Assembler & assembler() const
Definition: gridoperator.hh:117
Dune::PDELab::GridOperator::make_consistent
void make_consistent(Jacobian &a) const
Definition: gridoperator.hh:212
interpolate.hh
Dune::PDELab::GridOperator::assembler
Assembler & assembler()
Definition: gridoperator.hh:115
Dune::PDELab::EmptyTransformation
Definition: constraintstransformation.hh:111
Dune::PDELab::GridOperator::interpolate
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: gridoperator.hh:154
Dune::PDELab::DefaultLocalJacobianAssemblerEngine< DefaultLocalAssembler >
Dune::PDELab::OverlappingBorderDOFExchanger
Definition: borderdofexchanger.hh:577
Dune::PDELab::GridOperator::MatrixContainer
Definition: gridoperator.hh:71
Dune::PDELab::GridOperator::update
void const Domain & update
Definition: gridoperator.hh:205
Dune::PDELab::GridOperator::matrixBackend
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: gridoperator.hh:224
Dune::PDELab::GridOperator::globalSizeU
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator.hh:104
Dune::PDELab::DefaultAssembler::assemble
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition: default/assembler.hh:85
Dune::PDELab::DefaultLocalAssembler::localJacobianAssemblerEngine
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: default/localassembler.hh:169
Dune::PDELab::GridOperator::Domain
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: gridoperator.hh:43
Dune::PDELab::GridOperator::SetupGridOperator::visit
void visit(T &elem)
Definition: gridoperator.hh:131
Dune::PDELab::NonOverlappingBorderDOFExchanger
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:67
Dune::PDELab::GridOperator::trialGridFunctionSpace
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:92
Dune::PDELab::DefaultAssembler::testGridFunctionSpace
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: default/assembler.hh:73
Dune::PDELab::GridOperator::LocalAssembler
DefaultLocalAssembler< GridOperator, LOP, GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition > LocalAssembler
The local assembler type.
Definition: gridoperator.hh:57
Dune::PDELab::GridOperator::localAssembler
LocalAssembler & localAssembler() const
Definition: gridoperator.hh:119
Dune::PDELab::GridOperator::globalSizeV
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator.hh:110