dune-pdelab  2.0.0
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 
29  template<typename GFSU, typename GFSV, typename LOP,
30  typename MB, typename DF, typename RF, typename JF,
33  bool nonoverlapping_mode = false>
35  {
36  public:
37 
40 
47 
49  typedef typename MB::template Pattern<Jacobian,GFSV,GFSU> Pattern;
50 
54 
55  typedef typename conditional<
56  nonoverlapping_mode,
60 
63  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,Assembler,LocalAssembler> Traits;
64 
65  template <typename MFT>
67  typedef typename Traits::Jacobian Type;
68  };
69 
71  GridOperator(const GFSU & gfsu_, const CU & cu_, const GFSV & gfsv_, const CV & cv_, LOP & lop_, const MB& mb_ = MB())
72  : global_assembler(gfsu_,gfsv_,cu_,cv_)
73  , dof_exchanger(make_shared<BorderDOFExchanger>(*this))
74  , local_assembler(lop_, cu_, cv_,dof_exchanger)
75  , backend(mb_)
76  {}
77 
79  GridOperator(const GFSU & gfsu_, const GFSV & gfsv_, LOP & lop_, const MB& mb_ = MB())
80  : global_assembler(gfsu_,gfsv_)
81  , dof_exchanger(make_shared<BorderDOFExchanger>(*this))
82  , local_assembler(lop_,dof_exchanger)
83  , backend(mb_)
84  {}
85 
87  const GFSU& trialGridFunctionSpace() const
88  {
89  return global_assembler.trialGridFunctionSpace();
90  }
91 
93  const GFSV& testGridFunctionSpace() const
94  {
95  return global_assembler.testGridFunctionSpace();
96  }
97 
99  typename GFSU::Traits::SizeType globalSizeU () const
100  {
101  return trialGridFunctionSpace().globalSize();
102  }
103 
105  typename GFSV::Traits::SizeType globalSizeV () const
106  {
107  return testGridFunctionSpace().globalSize();
108  }
109 
110  Assembler & assembler() { return global_assembler; }
111 
112  const Assembler & assembler() const { return global_assembler; }
113 
114  LocalAssembler & localAssembler() const { return local_assembler; }
115 
116 
119  template <typename GridOperatorTuple>
122  : index(0), size(Dune::tuple_size<GridOperatorTuple>::value) {}
123 
124  template <typename T>
125  void visit(T& elem) {
126  elem.localAssembler().doPreProcessing = index == 0;
127  elem.localAssembler().doPostProcessing = index == size-1;
128  ++index;
129  }
130 
131  int index;
132  const int size;
133  };
134 
138  template<typename GridOperatorTuple>
139  static void setupGridOperators(GridOperatorTuple tuple)
140  {
141  Dune::ForEachValue<GridOperatorTuple> forEach(tuple);
142  SetupGridOperator<GridOperatorTuple> setup_visitor;
143  forEach.apply(setup_visitor);
144  }
145 
147  template<typename F, typename X>
148  void interpolate (const X& xold, F& f, X& x) const
149  {
150  // Interpolate f into grid function space and set corresponding coefficients
151  Dune::PDELab::interpolate(f,global_assembler.trialGridFunctionSpace(),x);
152 
153  // Copy non-constrained dofs from old time step
155  }
156 
158  void fill_pattern(Pattern & p) const {
159  typedef typename LocalAssembler::LocalPatternAssemblerEngine PatternEngine;
160  PatternEngine & pattern_engine = local_assembler.localPatternAssemblerEngine(p);
161  global_assembler.assemble(pattern_engine);
162  }
163 
165  void residual(const Domain & x, Range & r) const {
166  typedef typename LocalAssembler::LocalResidualAssemblerEngine ResidualEngine;
167  ResidualEngine & residual_engine = local_assembler.localResidualAssemblerEngine(r,x);
168  global_assembler.assemble(residual_engine);
169  }
170 
172  void jacobian(const Domain & x, Jacobian & a) const {
173  typedef typename LocalAssembler::LocalJacobianAssemblerEngine JacobianEngine;
174  JacobianEngine & jacobian_engine = local_assembler.localJacobianAssemblerEngine(a,x);
175  global_assembler.assemble(jacobian_engine);
176  }
177 
179  void jacobian_apply(const Domain & x, Range & r) const {
180  typedef typename LocalAssembler::LocalJacobianApplyAssemblerEngine JacobianApplyEngine;
181  JacobianApplyEngine & jacobian_apply_engine = local_assembler.localJacobianApplyAssemblerEngine(r,x);
182  global_assembler.assemble(jacobian_apply_engine);
183  }
184 
185  void make_consistent(Jacobian& a) const {
186  dof_exchanger->accumulateBorderEntries(*this,a);
187  }
188 
189  void update()
190  {
191  // the DOF exchanger has matrix information, so we need to update it
192  dof_exchanger->update(*this);
193  }
194 
196  const typename Traits::MatrixBackend& matrixBackend() const
197  {
198  return backend;
199  }
200 
201  private:
202  Assembler global_assembler;
203  shared_ptr<BorderDOFExchanger> dof_exchanger;
204 
205  mutable LocalAssembler local_assembler;
206  MB backend;
207 
208  };
209 
210  }
211 }
212 #endif
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: default/assembler.hh:76
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:158
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: gridoperator.hh:172
The assembler for standard DUNE grid.
Definition: default/assembler.hh:23
Definition: borderdofexchanger.hh:566
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: gridoperator.hh:196
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: gridoperator.hh:165
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:66
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
Definition: gridoperator.hh:66
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:71
void update()
Definition: gridoperator.hh:189
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:87
Definition: common/constraintstransformation.hh:127
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: default/assembler.hh:70
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: default/localassembler.hh:162
Assembler & assembler()
Definition: gridoperator.hh:110
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: common/constraints.hh:976
LocalAssembler & localAssembler() const
Definition: gridoperator.hh:114
Standard grid operator implementation.
Definition: gridoperator.hh:34
const int size
Definition: gridoperator.hh:132
void make_consistent(Jacobian &a) const
Definition: gridoperator.hh:185
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: gridoperator.hh:99
DefaultLocalAssembler< GridOperator, LOP, nonoverlapping_mode > LocalAssembler
The local assembler type.
Definition: gridoperator.hh:53
const Assembler & assembler() const
Definition: gridoperator.hh:112
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: gridoperator.hh:63
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: gridoperator.hh:105
conditional< nonoverlapping_mode, NonOverlappingBorderDOFExchanger< GridOperator >, OverlappingBorderDOFExchanger< GridOperator > >::type BorderDOFExchanger
Definition: gridoperator.hh:59
Definition: gridoperator.hh:120
Dune::PDELab::BackendVectorSelector< GFSU, DF >::Type Domain
The type of the domain (solution).
Definition: gridoperator.hh:42
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:226
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: default/localassembler.hh:152
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
static void setupGridOperators(GridOperatorTuple tuple)
Definition: gridoperator.hh:139
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
BackendVectorSelectorHelper< Backend, GridFunctionSpace, FieldType >::Type Type
Definition: backendselector.hh:14
int index
Definition: gridoperator.hh:131
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: default/localassembler.hh:172
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: gridoperator.hh:148
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:93
Traits::Jacobian Type
Definition: gridoperator.hh:67
const F & f
Definition: common/constraints.hh:145
SetupGridOperator()
Definition: gridoperator.hh:121
MBE::template Pattern< Jacobian, GFS, CGGFS > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: gridoperator.hh:49
Dune::PDELab::BackendMatrixSelector< MB, Domain, Range, JF >::Type Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
GridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: gridoperator.hh:79
R p(int i, D x)
Definition: qkdg.hh:62
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: default/localassembler.hh:143
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:191
void jacobian_apply(const Domain &x, Range &r) const
Apply jacobian matrix without explicitly assembling it.
Definition: gridoperator.hh:179
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition: default/assembler.hh:88
Backend::template MatrixHelper< VV, VU, E >::type Type
Definition: backendselector.hh:20
DefaultAssembler< GFSU, GFSV, CU, CV, nonoverlapping_mode > Assembler
The global assembler type.
Definition: gridoperator.hh:39
void visit(T &elem)
Definition: gridoperator.hh:125
Dune::PDELab::BackendVectorSelector< GFSV, RF >::Type Range
The type of the range (residual).
Definition: gridoperator.hh:44
Dune::PDELab::BackendMatrixSelector< MB, Domain, Range, JF >::Type Jacobian
The type of the jacobian.
Definition: gridoperator.hh:46