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 
33  template<typename GFSU, typename GFSV, typename LOP,
34  typename MB, typename DF, typename RF, typename JF,
37  >
39  {
40  public:
41 
44 
51 
53  typedef typename MB::template Pattern<Jacobian,GFSV,GFSU> Pattern;
54 
56  typedef DefaultLocalAssembler<
58  LOP,
59  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition
61 
62  // Fix this as soon as the default Partitions are constexpr
63  typedef typename std::conditional<
64  GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition,
68 
71  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,Assembler,LocalAssembler> Traits;
72 
73  template <typename MFT>
75  typedef typename Traits::Jacobian Type;
76  };
77 
79  GridOperator(const GFSU & gfsu_, const CU & cu_, const GFSV & gfsv_, const CV & cv_, LOP & lop_, const MB& mb_ = MB())
80  : global_assembler(gfsu_,gfsv_,cu_,cv_)
81  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
82  , local_assembler(lop_, cu_, cv_,dof_exchanger)
83  , backend(mb_)
84  {}
85 
87  GridOperator(const GFSU & gfsu_, const GFSV & gfsv_, LOP & lop_, const MB& mb_ = MB())
88  : global_assembler(gfsu_,gfsv_)
89  , dof_exchanger(std::make_shared<BorderDOFExchanger>(*this))
90  , local_assembler(lop_,dof_exchanger)
91  , backend(mb_)
92  {}
93 
95  const GFSU& trialGridFunctionSpace() const
96  {
97  return global_assembler.trialGridFunctionSpace();
98  }
99 
101  const GFSV& testGridFunctionSpace() const
102  {
103  return global_assembler.testGridFunctionSpace();
104  }
105 
107  typename GFSU::Traits::SizeType globalSizeU () const
108  {
109  return trialGridFunctionSpace().globalSize();
110  }
111 
113  typename GFSV::Traits::SizeType globalSizeV () const
114  {
115  return testGridFunctionSpace().globalSize();
116  }
117 
118  Assembler & assembler() { return global_assembler; }
119 
120  const Assembler & assembler() const { return global_assembler; }
121 
122  LocalAssembler & localAssembler() const { return local_assembler; }
123 
124 
127  template <typename GridOperatorTuple>
129  {
131  : index(0), size(std::tuple_size<GridOperatorTuple>::value) {}
132 
133  template <typename T>
134  void visit(T& elem) {
135  elem.localAssembler().preProcessing(index == 0);
136  elem.localAssembler().postProcessing(index == size-1);
137  ++index;
138  }
139 
140  int index;
141  const int size;
142  };
143 
147  template<typename GridOperatorTuple>
148  static void setupGridOperators(GridOperatorTuple tuple)
149  {
150  SetupGridOperator<GridOperatorTuple> setup_visitor;
151  Hybrid::forEach(tuple,
152  [&](auto &el) { setup_visitor.visit(el); });
153  }
154 
156  template<typename F, typename X>
157  void interpolate (const X& xold, F& f, X& x) const
158  {
159  // Interpolate f into grid function space and set corresponding coefficients
160  Dune::PDELab::interpolate(f,global_assembler.trialGridFunctionSpace(),x);
161 
162  // Copy non-constrained dofs from old time step
163  Dune::PDELab::copy_nonconstrained_dofs(local_assembler.trialConstraints(),xold,x);
164  }
165 
167  void fill_pattern(Pattern & p) const
168  {
169  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
170  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
171  global_assembler.assemble(pattern_engine);
172  }
173 
175  void residual(const Domain & x, Range & r) const
176  {
177  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
178  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
179  global_assembler.assemble(residual_engine);
180  }
181 
183  void jacobian(const Domain & x, Jacobian & a) const
184  {
185  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
186  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
187  global_assembler.assemble(jacobian_engine);
188  }
189 
191  void jacobian_apply(const Domain & z, Range & r) const
192  {
193  typedef typename LocalAssembler::LocalJacobianApplyAssemblerEngine JacobianApplyEngine;
194  JacobianApplyEngine & jacobian_apply_engine = local_assembler.localJacobianApplyAssemblerEngine(r,z);
195  global_assembler.assemble(jacobian_apply_engine);
196  }
197 
199  void nonlinear_jacobian_apply(const Domain & x, const Domain & z, Range & r) const
200  {
201  global_assembler.assemble(local_assembler.localNonlinearJacobianApplyAssemblerEngine(r,x,z));
202  }
203 
204 
205  void make_consistent(Jacobian& a) const
206  {
207  dof_exchanger->accumulateBorderEntries(*this,a);
208  }
209 
210  void update()
211  {
212  // the DOF exchanger has matrix information, so we need to update it
213  dof_exchanger->update(*this);
214  }
215 
217  const typename Traits::MatrixBackend& matrixBackend() const
218  {
219  return backend;
220  }
221 
222  private:
223  Assembler global_assembler;
224  shared_ptr<BorderDOFExchanger> dof_exchanger;
225 
226  mutable LocalAssembler local_assembler;
227  MB backend;
228 
229  };
230 
231  }
232 }
233 #endif // DUNE_PDELAB_GRIDOPERATOR_GRIDOPERATOR_HH
Dune::PDELab::Backend::Vector< CGGFS, field_type > Domain
The type of the domain (solution).
Definition: gridoperator.hh:46
Definition: borderdofexchanger.hh:577
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator.hh:71
Dune::PDELab::Backend::Matrix< MBE, Domain, Range, field_type > Jacobian
The type of the jacobian.
Definition: gridoperator.hh:50
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:127
Dune::PDELab::Backend::Vector< GFS, field_type > Range
The type of the range (residual).
Definition: gridoperator.hh:48
void update()
Definition: gridoperator.hh:210
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:167
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:989
MBE ::template Pattern< Jacobian, GFS, CGGFS > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator.hh:53
void jacobian_apply(const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:191
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:204
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: gridoperator.hh:217
Traits::Jacobian Type
Definition: gridoperator.hh:75
void nonlinear_jacobian_apply(const Domain &x, const Domain &z, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:199
Assembler & assembler()
Definition: gridoperator.hh:118
Definition: constraintstransformation.hh:111
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const int size
Definition: gridoperator.hh:141
void make_consistent(Jacobian &a) const
Definition: gridoperator.hh:205
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: gridoperator.hh:157
const P & p
Definition: constraints.hh:147
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
std::conditional< GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition, NonOverlappingBorderDOFExchanger< GridOperator >, OverlappingBorderDOFExchanger< GridOperator > >::type BorderDOFExchanger
Definition: gridoperator.hh:67
DefaultAssembler< GFSU, GFSV, CU, CV > Assembler
The global assembler type.
Definition: gridoperator.hh:43
Standard grid operator implementation.
Definition: gridoperator.hh:38
int index
Definition: gridoperator.hh:140
static void setupGridOperators(GridOperatorTuple tuple)
Definition: gridoperator.hh:148
void visit(T &elem)
Definition: gridoperator.hh:134
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:101
const Assembler & assembler() const
Definition: gridoperator.hh:120
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator.hh:175
STL namespace.
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator.hh:107
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:79
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:95
GridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: gridoperator.hh:87
Definition: gridoperator.hh:128
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
The local assembler for DUNE grids.
Definition: default/localassembler.hh:32
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
LocalAssembler & localAssembler() const
Definition: gridoperator.hh:122
The assembler for standard DUNE grid.
Definition: default/assembler.hh:22
std::size_t index
Definition: interpolate.hh:118
SetupGridOperator()
Definition: gridoperator.hh:130
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:67
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator.hh:113
DefaultLocalAssembler< GridOperator, LOP, GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition > LocalAssembler
The local assembler type.
Definition: gridoperator.hh:60
Definition: gridoperator.hh:74
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:183