dune-pdelab  2.5-dev
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 // TODO: Remove deprecated nonoverlapping parameter. This will break
5 // boilerplate code and needs to be done with care.
6 
7 #include <tuple>
8 
9 #include <dune/common/hybridutilities.hh>
10 
16 
17 namespace Dune{
18  namespace PDELab{
19 
20  namespace impl {
21 
22  template<int i>
24  {
25 
27 "the nonoverlapping_mode parameter on the GridOperator has been deprecated and will be removed after PDELab 2.4."
28 "The correct mode is now automatically deduced from the EntitySet of the function space.")
29  {}
30 
31  };
32 
33  template<>
35  {};
36 
37  }
38 
52  template<typename GFSU, typename GFSV, typename LOP,
53  typename MB, typename DF, typename RF, typename JF,
56  int nonoverlapping_mode = -1>
58  : public impl::warn_on_deprecated_nonoverlapping_mode_parameter<nonoverlapping_mode>
59  {
60  public:
61 
62  static_assert(nonoverlapping_mode == -1 ||
63  nonoverlapping_mode == 0 ||
64  nonoverlapping_mode == 1,
65  "invalid value for nonoverlapping_mode! This parameter is also deprecated in PDELab 2.4, so please remove it from your typedefs!");
66 
69 
76 
78  typedef typename MB::template Pattern<Jacobian,GFSV,GFSU> Pattern;
79 
81  typedef DefaultLocalAssembler<
83  LOP,
84  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition
86 
87  // Fix this as soon as the default Partitions are constexpr
88  typedef typename std::conditional<
89  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition,
93 
96  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,Assembler,LocalAssembler> Traits;
97 
98  template <typename MFT>
100  typedef typename Traits::Jacobian Type;
101  };
102 
104  GridOperator(const GFSU & gfsu_, const CU & cu_, const GFSV & gfsv_, const CV & cv_, LOP & lop_, const MB& mb_ = MB())
105  : global_assembler(gfsu_,gfsv_,cu_,cv_)
106  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
107  , local_assembler(lop_, cu_, cv_,dof_exchanger)
108  , backend(mb_)
109  {}
110 
112  GridOperator(const GFSU & gfsu_, const GFSV & gfsv_, LOP & lop_, const MB& mb_ = MB())
113  : global_assembler(gfsu_,gfsv_)
114  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
115  , local_assembler(lop_,dof_exchanger)
116  , backend(mb_)
117  {}
118 
120  const GFSU& trialGridFunctionSpace() const
121  {
122  return global_assembler.trialGridFunctionSpace();
123  }
124 
126  const GFSV& testGridFunctionSpace() const
127  {
128  return global_assembler.testGridFunctionSpace();
129  }
130 
132  typename GFSU::Traits::SizeType globalSizeU () const
133  {
134  return trialGridFunctionSpace().globalSize();
135  }
136 
138  typename GFSV::Traits::SizeType globalSizeV () const
139  {
140  return testGridFunctionSpace().globalSize();
141  }
142 
143  Assembler & assembler() { return global_assembler; }
144 
145  const Assembler & assembler() const { return global_assembler; }
146 
147  LocalAssembler & localAssembler() const { return local_assembler; }
148 
149 
152  template <typename GridOperatorTuple>
155  : index(0), size(std::tuple_size<GridOperatorTuple>::value) {}
156 
157  template <typename T>
158  void visit(T& elem) {
159  elem.localAssembler().doPreProcessing = index == 0;
160  elem.localAssembler().doPostProcessing = index == size-1;
161  ++index;
162  }
163 
164  int index;
165  const int size;
166  };
167 
171  template<typename GridOperatorTuple>
172  static void setupGridOperators(GridOperatorTuple tuple)
173  {
174  SetupGridOperator<GridOperatorTuple> setup_visitor;
175  Hybrid::forEach(tuple,
176  [&](auto &el) { setup_visitor.visit(el); });
177  }
178 
180  template<typename F, typename X>
181  void interpolate (const X& xold, F& f, X& x) const
182  {
183  // Interpolate f into grid function space and set corresponding coefficients
184  Dune::PDELab::interpolate(f,global_assembler.trialGridFunctionSpace(),x);
185 
186  // Copy non-constrained dofs from old time step
187  Dune::PDELab::copy_nonconstrained_dofs(local_assembler.trialConstraints(),xold,x);
188  }
189 
191  void fill_pattern(Pattern & p) const {
192  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
193  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
194  global_assembler.assemble(pattern_engine);
195  }
196 
198  void residual(const Domain & x, Range & r) const {
199  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
200  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
201  global_assembler.assemble(residual_engine);
202  }
203 
205  void jacobian(const Domain & x, Jacobian & a) const {
206  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
207  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
208  global_assembler.assemble(jacobian_engine);
209  }
210 
212  void jacobian_apply(const Domain & z, Range & r) const {
213  typedef typename LocalAssembler::LocalJacobianApplyAssemblerEngine JacobianApplyEngine;
214  JacobianApplyEngine & jacobian_apply_engine = local_assembler.localJacobianApplyAssemblerEngine(r,z);
215  global_assembler.assemble(jacobian_apply_engine);
216  }
217 
219  void nonlinear_jacobian_apply(const Domain & x, const Domain & z, Range & r) const {
220  global_assembler.assemble(local_assembler.localNonlinearJacobianApplyAssemblerEngine(r,x,z));
221  }
222 
223 
224  void make_consistent(Jacobian& a) const {
225  dof_exchanger->accumulateBorderEntries(*this,a);
226  }
227 
228  void update()
229  {
230  // the DOF exchanger has matrix information, so we need to update it
231  dof_exchanger->update(*this);
232  }
233 
235  const typename Traits::MatrixBackend& matrixBackend() const
236  {
237  return backend;
238  }
239 
240  private:
241  Assembler global_assembler;
242  shared_ptr<BorderDOFExchanger> dof_exchanger;
243 
244  mutable LocalAssembler local_assembler;
245  MB backend;
246 
247  };
248 
249  }
250 }
251 #endif // DUNE_PDELAB_GRIDOPERATOR_GRIDOPERATOR_HH
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
int index
Definition: gridoperator.hh:164
STL namespace.
Definition: constraintstransformation.hh:111
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
void jacobian_apply(const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:212
std::conditional< GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition, NonOverlappingBorderDOFExchanger< GridOperator >, OverlappingBorderDOFExchanger< GridOperator > >::type BorderDOFExchanger
Definition: gridoperator.hh:92
Definition: gridoperator.hh:99
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:191
void visit(T &elem)
Definition: gridoperator.hh:158
void nonlinear_jacobian_apply(const Domain &x, const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:219
const int size
Definition: gridoperator.hh:165
void update()
Definition: gridoperator.hh:228
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:120
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator.hh:132
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:191
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
DefaultLocalAssembler< GridOperator, LOP, GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition > LocalAssembler
The local assembler type.
Definition: gridoperator.hh:85
SetupGridOperator()
Definition: gridoperator.hh:154
Standard grid operator implementation.
Definition: gridoperator.hh:57
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:104
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: gridoperator.hh:181
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:126
GridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: gridoperator.hh:112
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:127
MBE ::template Pattern< Jacobian, GFS, CGGFS > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator.hh:78
const Assembler & assembler() const
Definition: gridoperator.hh:145
Dune::PDELab::Backend::Vector< CGGFS, field_type > Domain
The type of the domain (solution).
Definition: gridoperator.hh:71
Dune::PDELab::Backend::Vector< GFS, field_type > Range
The type of the range (residual).
Definition: gridoperator.hh:73
const P & p
Definition: constraints.hh:147
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator.hh:96
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
void make_consistent(Jacobian &a) const
Definition: gridoperator.hh:224
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: gridoperator.hh:235
LocalAssembler & localAssembler() const
Definition: gridoperator.hh:147
warn_on_deprecated_nonoverlapping_mode_parameter()
Definition: gridoperator.hh:26
The local assembler for DUNE grids.
Definition: default/localassembler.hh:32
Assembler & assembler()
Definition: gridoperator.hh:143
static void setupGridOperators(GridOperatorTuple tuple)
Definition: gridoperator.hh:172
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:67
Definition: gridoperator.hh:153
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:205
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator.hh:138
Traits::Jacobian Type
Definition: gridoperator.hh:100
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator.hh:198
Dune::PDELab::Backend::Matrix< MBE, Domain, Range, field_type > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:75
Definition: borderdofexchanger.hh:577
The assembler for standard DUNE grid.
Definition: default/assembler.hh:23