dune-pdelab  2.4.1
gridoperator.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_GRIDOPERATOR_HH
2 #define DUNE_PDELAB_GRIDOPERATOR_HH
3 
4 #include <dune/common/tupleutility.hh>
5 
11 
12 namespace Dune{
13  namespace PDELab{
14 
15  namespace impl {
16 
17  template<int i>
19  {
20 
22 "the nonoverlapping_mode parameter on the GridOperator has been deprecated and will be removed after PDELab 2.4."
23 "The correct mode is now automatically deduced from the EntitySet of the function space.")
24  {}
25 
26  };
27 
28  template<>
30  {};
31 
32  }
33 
47  template<typename GFSU, typename GFSV, typename LOP,
48  typename MB, typename DF, typename RF, typename JF,
51  int nonoverlapping_mode = -1>
53  : public impl::warn_on_deprecated_nonoverlapping_mode_parameter<nonoverlapping_mode>
54  {
55  public:
56 
57  static_assert(nonoverlapping_mode == -1 ||
58  nonoverlapping_mode == 0 ||
59  nonoverlapping_mode == 1,
60  "invalid value for nonoverlapping_mode! This parameter is also deprecated in PDELab 2.4, so please remove it from your typedefs!");
61 
64 
71 
73  typedef typename MB::template Pattern<Jacobian,GFSV,GFSU> Pattern;
74 
76  typedef DefaultLocalAssembler<
78  LOP,
79  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition
81 
82  // Fix this as soon as the default Partitions are constexpr
83  typedef typename conditional<
84  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition,
88 
91  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,Assembler,LocalAssembler> Traits;
92 
93  template <typename MFT>
95  typedef typename Traits::Jacobian Type;
96  };
97 
99  GridOperator(const GFSU & gfsu_, const CU & cu_, const GFSV & gfsv_, const CV & cv_, LOP & lop_, const MB& mb_ = MB())
100  : global_assembler(gfsu_,gfsv_,cu_,cv_)
101  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
102  , local_assembler(lop_, cu_, cv_,dof_exchanger)
103  , backend(mb_)
104  {}
105 
107  GridOperator(const GFSU & gfsu_, const GFSV & gfsv_, LOP & lop_, const MB& mb_ = MB())
108  : global_assembler(gfsu_,gfsv_)
109  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
110  , local_assembler(lop_,dof_exchanger)
111  , backend(mb_)
112  {}
113 
115  const GFSU& trialGridFunctionSpace() const
116  {
117  return global_assembler.trialGridFunctionSpace();
118  }
119 
121  const GFSV& testGridFunctionSpace() const
122  {
123  return global_assembler.testGridFunctionSpace();
124  }
125 
127  typename GFSU::Traits::SizeType globalSizeU () const
128  {
129  return trialGridFunctionSpace().globalSize();
130  }
131 
133  typename GFSV::Traits::SizeType globalSizeV () const
134  {
135  return testGridFunctionSpace().globalSize();
136  }
137 
138  Assembler & assembler() { return global_assembler; }
139 
140  const Assembler & assembler() const { return global_assembler; }
141 
142  LocalAssembler & localAssembler() const { return local_assembler; }
143 
144 
147  template <typename GridOperatorTuple>
150  : index(0), size(Dune::tuple_size<GridOperatorTuple>::value) {}
151 
152  template <typename T>
153  void visit(T& elem) {
154  elem.localAssembler().doPreProcessing = index == 0;
155  elem.localAssembler().doPostProcessing = index == size-1;
156  ++index;
157  }
158 
159  int index;
160  const int size;
161  };
162 
166  template<typename GridOperatorTuple>
167  static void setupGridOperators(GridOperatorTuple tuple)
168  {
169  Dune::ForEachValue<GridOperatorTuple> forEach(tuple);
170  SetupGridOperator<GridOperatorTuple> setup_visitor;
171  forEach.apply(setup_visitor);
172  }
173 
175  template<typename F, typename X>
176  void interpolate (const X& xold, F& f, X& x) const
177  {
178  // Interpolate f into grid function space and set corresponding coefficients
179  Dune::PDELab::interpolate(f,global_assembler.trialGridFunctionSpace(),x);
180 
181  // Copy non-constrained dofs from old time step
182  Dune::PDELab::copy_nonconstrained_dofs(local_assembler.trialConstraints(),xold,x);
183  }
184 
186  void fill_pattern(Pattern & p) const {
187  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
188  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
189  global_assembler.assemble(pattern_engine);
190  }
191 
193  void residual(const Domain & x, Range & r) const {
194  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
195  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
196  global_assembler.assemble(residual_engine);
197  }
198 
200  void jacobian(const Domain & x, Jacobian & a) const {
201  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
202  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
203  global_assembler.assemble(jacobian_engine);
204  }
205 
207  void jacobian_apply(const Domain & z, Range & r) const {
208  typedef typename LocalAssembler::LocalJacobianApplyAssemblerEngine JacobianApplyEngine;
209  JacobianApplyEngine & jacobian_apply_engine = local_assembler.localJacobianApplyAssemblerEngine(r,z);
210  global_assembler.assemble(jacobian_apply_engine);
211  }
212 
214  void nonlinear_jacobian_apply(const Domain & x, const Domain & z, Range & r) const {
215  global_assembler.assemble(local_assembler.localNonlinearJacobianApplyAssemblerEngine(r,x,z));
216  }
217 
218 
219  void make_consistent(Jacobian& a) const {
220  dof_exchanger->accumulateBorderEntries(*this,a);
221  }
222 
223  void update()
224  {
225  // the DOF exchanger has matrix information, so we need to update it
226  dof_exchanger->update(*this);
227  }
228 
230  const typename Traits::MatrixBackend& matrixBackend() const
231  {
232  return backend;
233  }
234 
235  private:
236  Assembler global_assembler;
237  shared_ptr<BorderDOFExchanger> dof_exchanger;
238 
239  mutable LocalAssembler local_assembler;
240  MB backend;
241 
242  };
243 
244  }
245 }
246 #endif
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
STL namespace.
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator.hh:91
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:134
The local assembler for DUNE grids.
Definition: default/localassembler.hh:32
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
LocalAssembler & localAssembler() const
Definition: gridoperator.hh:142
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:121
Definition: borderdofexchanger.hh:576
static const unsigned int value
Definition: gridfunctionspace/tags.hh:177
void make_consistent(Jacobian &a) const
Definition: gridoperator.hh:219
Dune::PDELab::Backend::Vector< GFS, field_type > Range
The type of the range (residual).
Definition: gridoperator.hh:68
Dune::PDELab::Backend::Matrix< MBE, Domain, Range, field_type > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:70
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator.hh:133
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:99
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:115
const int size
Definition: gridoperator.hh:160
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
const Assembler & assembler() const
Definition: gridoperator.hh:140
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:191
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: gridoperator.hh:176
Definition: gridoperator.hh:148
Definition: adaptivity.hh:27
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator.hh:127
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:186
int index
Definition: gridoperator.hh:159
The assembler for standard DUNE grid.
Definition: default/assembler.hh:23
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
Traits::Jacobian Type
Definition: gridoperator.hh:95
static void setupGridOperators(GridOperatorTuple tuple)
Definition: gridoperator.hh:167
void update()
Definition: gridoperator.hh:223
void jacobian_apply(const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:207
Definition: gridoperator.hh:94
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator.hh:193
Standard grid operator implementation.
Definition: gridoperator.hh:52
Dune::PDELab::Backend::Vector< CGGFS, field_type > Domain
The type of the domain (solution).
Definition: gridoperator.hh:66
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:66
SetupGridOperator()
Definition: gridoperator.hh:149
GridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: gridoperator.hh:107
void visit(T &elem)
Definition: gridoperator.hh:153
Definition: constraintstransformation.hh:111
const P & p
Definition: constraints.hh:147
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:113
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: gridoperator.hh:230
warn_on_deprecated_nonoverlapping_mode_parameter()
Definition: gridoperator.hh:21
conditional< GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition, NonOverlappingBorderDOFExchanger< GridOperator >, OverlappingBorderDOFExchanger< GridOperator > >::type BorderDOFExchanger
Definition: gridoperator.hh:87
Assembler & assembler()
Definition: gridoperator.hh:138
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:200
MBE::template Pattern< Jacobian, GFS, CGGFS > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator.hh:73
DefaultLocalAssembler< GridOperator, LOP, GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition > LocalAssembler
The local assembler type.
Definition: gridoperator.hh:80
void nonlinear_jacobian_apply(const Domain &x, const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:214