dune-pdelab  2.4-dev
common/assembler.hh
Go to the documentation of this file.
2 #include <assemblerutilties.hh>
3 
4 namespace Dune {
5  namespace PDELab {
6 
10 
21  public:
22  template<class LocalAssemblerEngine>
23  void assemble(LocalAssemblerEngine & local_assembler_engine);
24  };
25 
26 
32  public:
35 
37  const LocalAssembler & localAssembler();
38 
46  bool requireSkeleton() const;
47  bool requireSkeletonTwoSided() const;
48  bool requireUVVolume() const;
49  bool requireVVolume() const;
50  bool requireUVSkeleton() const;
51  bool requireVSkeleton() const;
52  bool requireUVBoundary() const;
53  bool requireVBoundary() const;
54  bool requireUVProcessor() const;
55  bool requireVProcessor() const;
56  bool requireUVEnrichedCoupling() const;
57  bool requireVEnrichedCoupling() const;
58  bool requireUVVolumePostSkeleton() const;
59  bool requireVVolumePostSkeleton() const;
61 
81  template<typename EG>
82  bool assembleCell(const EG & eg);
83 
87  template<typename EG, typename LFSU, typename LFSV>
88  void assembleUVVolume(const EG & eg, const LFSU & lfsu, const LFSV & lfsv);
89 
93  template<typename EG, typename LFSV>
94  void assembleVVolume(const EG & eg, const LFSV & lfsv);
95 
99  template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N>
100  void assembleUVSkeleton(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
101  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
102 
106  template<typename IG, typename LFSV_S, typename LFSV_N>
107  void assembleVSkeleton(const IG & ig, const LFSV_S & lfsv_s, const LFSV_N & lfsv_n);
108 
112  template<typename IG, typename LFSU_S, typename LFSV_S>
113  void assembleUVBoundary(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
114 
118  template<typename IG, typename LFSV_S>
119  void assembleVBoundary(const IG & ig, const LFSV_S & lfsv_s);
120 
126  template<typename IG, typename LFSU_S, typename LFSV_S>
127  void assembleUVProcessor(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
128 
134  template<typename IG, typename LFSV_S>
135  void assembleVProcessor(const IG & ig, const LFSV_S & lfsv_s);
136 
137  template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N,
138  typename LFSU_C, typename LFSV_C>
139  void assembleUVEnrichedCoupling(const IG & ig,
140  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
141  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n,
142  const LFSU_C & lfsu_c, const LFSV_C & lfsv_c);
143 
144  template<typename IG, typename LFSV_S, typename LFSV_N, typename LFSV_C>
145  void assembleVEnrichedCoupling(const IG & ig,
146  const LFSV_S & lfsv_s,
147  const LFSV_N & lfsv_n,
148  const LFSV_C & lfsv_c);
149 
154  template<typename EG, typename LFSU, typename LFSV>
155  void assembleUVVolumePostSkeleton(const EG & eg, const LFSU & lfsu, const LFSV & lfsv);
156 
161  template<typename EG, typename LFSV>
162  void assembleVVolumePostSkeleton(const EG & eg, const LFSV & lfsv);
163 
165  void preAssembly();
166 
168  void postAssembly();
169 
182  void onBindLFSUV(const EG & eg,
183  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
184  void onBindLFSV(const EG & eg,
185  const LFSV_S & lfsv_s);
186  void onBindLFSUVInside(const IG & ig,
187  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
188  void onBindLFSVInside(const IG & ig,
189  const LFSV_S & lfsv_s);
190  void onBindLFSUVOutside(const IG & ig,
191  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
192  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
193  void onBindLFSVOutside(const IG & ig,
194  const LFSV_S & lfsv_s,
195  const LFSV_N & lfsv_n);
196  void onBindLFSUVCoupling(const IG & ig,
197  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
198  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n
199  const LFSU_Coupling & lfsu_coupling, const LFSV_Coupling & lfsv_coupling);
200  void onBindLFSVCoupling(const IG & ig,
201  const LFSV_S & lfsv_s,
202  const LFSV_N & lfsv_n,
203  const LFSV_Coupling & lfsv_coupling);
204 
205  void onUnbindLFSUV(const EG & eg,
206  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
207  void onUnbindLFSV(const EG & eg,
208  const LFSV_S & lfsv_s);
209  void onUnbindLFSUVInside(const IG & ig,
210  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
211  void onUnbindLFSVInside(const IG & ig,
212  const LFSV_S & lfsv_s);
213  void onUnbindLFSUVOutside(const IG & ig,
214  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
215  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
216  void onUnbindLFSVOutside(const IG & ig,
217  const LFSV_S & lfsv_s,
218  const LFSV_N & lfsv_n);
219  void onUnbindLFSUVCoupling(const IG & ig,
220  const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
221  const LFSU_N & lfsu_n, const LFSV_N & lfsv_n,
222  const LFSU_Coupling & lfsu_coupling, const LFSV_Coupling & lfsv_coupling);
223  void onUnbindLFSVCoupling(const IG & ig,
224  const LFSV_S & lfsv_s,
225  const LFSV_N & lfsv_n,
226  const LFSV_Coupling & lfsv_coupling);
227 
238  void loadCoefficientsLFSUInside(const LFSU_S & lfsu_s);
239  void loadCoefficientsLFSUOutside(const LFSU_N & lfsu_n);
240  void loadCoefficientsLFSUCoupling(const LFSU_Coupling & lfsu_coupling);
253  void setSolution(const X& x);
254  void setPattern(const P& p);
255  void setJacobian(const J & j);
256  void setResidual(const R& r);
259  };
260 
272  template<typename B, typename CU, typename CV>
274  {
275  public:
276 
281  template<class TT>
283  void setTime(TT time);
284 
286  template<typename TT>
287  void preStep (TT time, TT dt, std::size_t stages);
288 
290  void postStep ();
291 
293  template<typename TT>
294  void preStage (TT time, std::size_t stage);
295 
297  void postStage ();
298 
300  template<typename TT>
301  TT suggestTimestep (TT dt) const;
302 
305  template<class RF>
307  void setWeight(RF weight);
308 
327  };
328 
342  template<typename GFSU, typename GFSV,
343  typename MB, typename DF, typename RF, typename JF>
345  public:
346 
348  typedef GridOperatorTraits
349  <GFSU,GFSV,MB,DF,RF,JF,CU,CV,AssemblerInterface,LocalAssemblerInterface> Traits;
350 
352  template<typename P>
353  void fill_pattern (P& globalpattern) const;
354 
357  template<typename X, typename R>
358  void residual (const X& x, R& r) const;
359 
362  template<typename X, typename A>
363  void jacobian (const X& x, A& a) const;
364 
367  Assembler & assembler();
370 
373  const GFSU& trialGridFunctionSpace() const;
374  const GFSV& testGridFunctionSpace() const;
375  typename GFSU::Traits::SizeType globalSizeU () const;
376  typename GFSV::Traits::SizeType globalSizeV () const;
378 
380 
386  template<typename F>
387  void interpolate(const typename Traits::Domain& xold,
388  const F& f,
389  typename Traits::Domain& xnew);
390 
395 
410  template<typename GridOperatorTuple>
411  static void setupGridOperators(GridOperatorTuple& tuple);
412 
413  };
415  };
416 };
417 
void onBindLFSVInside(const IG &ig, const LFSV_S &lfsv_s)
LocalAssemblerInterface & localAssembler()
void onUnbindLFSUVCoupling(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n, const LFSU_Coupling &lfsu_coupling, const LFSV_Coupling &lfsv_coupling)
void jacobian(const X &x, A &a) const
void assemble(LocalAssemblerEngine &local_assembler_engine)
void assembleUVSkeleton(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
void postAssembly()
Called last thing after assembling.
void preStage(TT time, std::size_t stage)
Notify local assembler about upcoming time step stage.
static void setupGridOperators(GridOperatorTuple &tuple)
void assembleUVVolumePostSkeleton(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
void onBindLFSUVCoupling(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n const LFSU_Coupling &lfsu_coupling, const LFSV_Coupling &lfsv_coupling)
void onUnbindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: gridoperatorutilities.hh:58
void loadCoefficientsLFSUCoupling(const LFSU_Coupling &lfsu_coupling)
LocalPatternAssemblerEngine & localPatternAssemblerEngine(P &p)
GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, AssemblerInterface, LocalAssemblerInterface > Traits
The traits class.
Definition: common/assembler.hh:349
void postStep()
Notify local assembler about completion of time step.
void assembleUVProcessor(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void onUnbindLFSVInside(const IG &ig, const LFSV_S &lfsv_s)
const GFSV & testGridFunctionSpace() const
void assembleUVBoundary(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void onUnbindLFSVCoupling(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n, const LFSV_Coupling &lfsv_coupling)
void assembleVEnrichedCoupling(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n, const LFSV_C &lfsv_c)
const IG & ig
Definition: constraints.hh:147
void postStage()
Notify local assembler about completion of time step stage.
void residual(const X &x, R &r) const
void onUnbindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void assembleVVolume(const EG &eg, const LFSV &lfsv)
void onBindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void assembleUVEnrichedCoupling(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n, const LFSU_C &lfsu_c, const LFSV_C &lfsv_c)
void onBindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
GFSV::Traits::SizeType globalSizeV() const
void loadCoefficientsLFSUInside(const LFSU_S &lfsu_s)
void onUnbindLFSUVInside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
Base class for local assembler.
Definition: assemblerutilities.hh:210
void assembleVVolumePostSkeleton(const EG &eg, const LFSV &lfsv)
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(A &a, const X &x)
const GFSU & trialGridFunctionSpace() const
The local assembler engine which handles the integration parts as provided by the global assemblers...
Definition: common/assembler.hh:31
Definition: adaptivity.hh:27
void onUnbindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void onBindLFSVCoupling(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n, const LFSV_Coupling &lfsv_coupling)
void assembleVProcessor(const IG &ig, const LFSV_S &lfsv_s)
void setTime(TT time)
Set current time of assembling.
void onBindLFSUVInside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assembleVSkeleton(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void fill_pattern(P &globalpattern) const
Determines the sparsity pattern of the jacobian matrix.
void interpolate(const typename Traits::Domain &xold, const F &f, typename Traits::Domain &xnew)
Interpolate xnew from f, taking unconstrained values from xold.
void onBindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void loadCoefficientsLFSUOutside(const LFSU_N &lfsu_n)
LocalResidualAssemblerEngine & localResidualAssemblerEngine(R &r, const X &x)
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
const P & p
Definition: constraints.hh:146
void assembleVBoundary(const IG &ig, const LFSV_S &lfsv_s)
void onUnbindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
TT suggestTimestep(TT dt) const
Suggest a valid time step size.
GFSU::Traits::SizeType globalSizeU() const
The global assembler which performs the traversing of the integration parts.
Definition: common/assembler.hh:20
The local assembler which provides the engines that drive the global assembler.
Definition: common/assembler.hh:273
const LocalAssembler & localAssembler()
Access to the superior local assembler object.
void preAssembly()
Called directly before assembling.
void preStep(TT time, TT dt, std::size_t stages)
Notify local assembler about upcoming time step.
LocalAssemblerInterface LocalAssembler
The type of the local assembler.
Definition: common/assembler.hh:34
LocalResidualJacobianAssemblerEngine & localResidualJacobianAssemblerEngine(R &r, A &a, const X &x)
const EG & eg
Definition: constraints.hh:280
void setWeight(RF weight)
Set current weight of assembling.
void onBindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
The grid operator represents an operator mapping which corresponds to the (possibly nonlinear) algebr...
Definition: common/assembler.hh:344