dune-pdelab 2.7-git
Loading...
Searching...
No Matches
common/assembler.hh
Go to the documentation of this file.
1#ifndef DUNE_PDELAB_GRIDOPERATOR_COMMON_ASSEMBLER_HH
2#define DUNE_PDELAB_GRIDOPERATOR_COMMON_ASSEMBLER_HH
4#include <assemblerutilties.hh>
5
6namespace Dune {
7 namespace PDELab {
8
12
23 public:
24 template<class LocalAssemblerEngine>
25 void assemble(LocalAssemblerEngine & local_assembler_engine);
26 };
27
28
34 public:
37
40
48 bool requireSkeleton() const;
50 bool requireUVVolume() const;
51 bool requireVVolume() const;
52 bool requireUVSkeleton() const;
53 bool requireVSkeleton() const;
54 bool requireUVBoundary() const;
55 bool requireVBoundary() const;
56 bool requireUVProcessor() const;
57 bool requireVProcessor() const;
63
83 template<typename EG>
84 bool assembleCell(const EG & eg);
85
89 template<typename EG, typename LFSU, typename LFSV>
90 void assembleUVVolume(const EG & eg, const LFSU & lfsu, const LFSV & lfsv);
91
95 template<typename EG, typename LFSV>
96 void assembleVVolume(const EG & eg, const LFSV & lfsv);
97
101 template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N>
102 void assembleUVSkeleton(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
103 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
104
108 template<typename IG, typename LFSV_S, typename LFSV_N>
109 void assembleVSkeleton(const IG & ig, const LFSV_S & lfsv_s, const LFSV_N & lfsv_n);
110
114 template<typename IG, typename LFSU_S, typename LFSV_S>
115 void assembleUVBoundary(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
116
120 template<typename IG, typename LFSV_S>
121 void assembleVBoundary(const IG & ig, const LFSV_S & lfsv_s);
122
128 template<typename IG, typename LFSU_S, typename LFSV_S>
129 void assembleUVProcessor(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
130
136 template<typename IG, typename LFSV_S>
137 void assembleVProcessor(const IG & ig, const LFSV_S & lfsv_s);
138
139 template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N,
140 typename LFSU_C, typename LFSV_C>
142 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
143 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n,
144 const LFSU_C & lfsu_c, const LFSV_C & lfsv_c);
145
146 template<typename IG, typename LFSV_S, typename LFSV_N, typename LFSV_C>
148 const LFSV_S & lfsv_s,
149 const LFSV_N & lfsv_n,
150 const LFSV_C & lfsv_c);
151
156 template<typename EG, typename LFSU, typename LFSV>
157 void assembleUVVolumePostSkeleton(const EG & eg, const LFSU & lfsu, const LFSV & lfsv);
158
163 template<typename EG, typename LFSV>
164 void assembleVVolumePostSkeleton(const EG & eg, const LFSV & lfsv);
165
168
171
184 void onBindLFSUV(const EG & eg,
185 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
186 void onBindLFSV(const EG & eg,
187 const LFSV_S & lfsv_s);
188 void onBindLFSUVInside(const IG & ig,
189 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
190 void onBindLFSVInside(const IG & ig,
191 const LFSV_S & lfsv_s);
192 void onBindLFSUVOutside(const IG & ig,
193 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
194 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
195 void onBindLFSVOutside(const IG & ig,
196 const LFSV_S & lfsv_s,
197 const LFSV_N & lfsv_n);
198 void onBindLFSUVCoupling(const IG & ig,
199 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
200 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n
201 const LFSU_Coupling & lfsu_coupling, const LFSV_Coupling & lfsv_coupling);
202 void onBindLFSVCoupling(const IG & ig,
203 const LFSV_S & lfsv_s,
204 const LFSV_N & lfsv_n,
205 const LFSV_Coupling & lfsv_coupling);
206
207 void onUnbindLFSUV(const EG & eg,
208 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
209 void onUnbindLFSV(const EG & eg,
210 const LFSV_S & lfsv_s);
211 void onUnbindLFSUVInside(const IG & ig,
212 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s);
213 void onUnbindLFSVInside(const IG & ig,
214 const LFSV_S & lfsv_s);
215 void onUnbindLFSUVOutside(const IG & ig,
216 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
217 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n);
218 void onUnbindLFSVOutside(const IG & ig,
219 const LFSV_S & lfsv_s,
220 const LFSV_N & lfsv_n);
221 void onUnbindLFSUVCoupling(const IG & ig,
222 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
223 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n,
224 const LFSU_Coupling & lfsu_coupling, const LFSV_Coupling & lfsv_coupling);
225 void onUnbindLFSVCoupling(const IG & ig,
226 const LFSV_S & lfsv_s,
227 const LFSV_N & lfsv_n,
228 const LFSV_Coupling & lfsv_coupling);
229
240 void loadCoefficientsLFSUInside(const LFSU_S & lfsu_s);
241 void loadCoefficientsLFSUOutside(const LFSU_N & lfsu_n);
242 void loadCoefficientsLFSUCoupling(const LFSU_Coupling & lfsu_coupling);
255 void setSolution(const X& x);
256 void setPattern(const P& p);
257 void setJacobian(const J & j);
258 void setResidual(const R& r);
261 };
262
274 template<typename B, typename CU, typename CV>
276 {
277 public:
278
284 template<class TT>
285 void setTime(TT time);
286
288 template<typename TT>
289 void preStep (TT time, TT dt, std::size_t stages);
290
292 void postStep ();
293
295 template<typename TT>
296 void preStage (TT time, std::size_t stage);
297
299 void postStage ();
300
302 template<typename TT>
303 TT suggestTimestep (TT dt) const;
304
308 template<class RF>
309 void setWeight(RF weight);
310
329 };
330
344 template<typename GFSU, typename GFSV,
345 typename MB, typename DF, typename RF, typename JF>
347 public:
348
350 typedef GridOperatorTraits
351 <GFSU,GFSV,MB,DF,RF,JF,CU,CV,AssemblerInterface,LocalAssemblerInterface> Traits;
352
354 template<typename P>
355 void fill_pattern (P& globalpattern) const;
356
359 template<typename X, typename R>
360 void residual (const X& x, R& r) const;
361
364 template<typename X, typename A>
365 void jacobian (const X& x, A& a) const;
366
369 Assembler & assembler();
372
375 const GFSU& trialGridFunctionSpace() const;
376 const GFSV& testGridFunctionSpace() const;
377 typename GFSU::Traits::SizeType globalSizeU () const;
378 typename GFSV::Traits::SizeType globalSizeV () const;
380
382
388 template<typename F>
389 void interpolate(const typename Traits::Domain& xold,
390 const F& f,
391 typename Traits::Domain& xnew);
392
397
412 template<typename GridOperatorTuple>
413 static void setupGridOperators(GridOperatorTuple& tuple);
414
415 };
417 };
418};
419
522#endif // DUNE_PDELAB_GRIDOPERATOR_COMMON_ASSEMBLER_HH
const P & p
Definition constraints.hh:148
const IG & ig
Definition constraints.hh:149
For backward compatibility – Do not use this!
Definition adaptivity.hh:28
The global assembler which performs the traversing of the integration parts.
Definition common/assembler.hh:22
void assemble(LocalAssemblerEngine &local_assembler_engine)
The local assembler engine which handles the integration parts as provided by the global assemblers.
Definition common/assembler.hh:33
void onUnbindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
void onBindLFSUV(const EG &eg, 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 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 assembleUVProcessor(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
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 loadCoefficientsLFSUCoupling(const LFSU_Coupling &lfsu_coupling)
void loadCoefficientsLFSUInside(const LFSU_S &lfsu_s)
void assembleVProcessor(const IG &ig, const LFSV_S &lfsv_s)
void assembleVEnrichedCoupling(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n, const LFSV_C &lfsv_c)
const LocalAssembler & localAssembler()
Access to the superior local assembler object.
void onUnbindLFSUVInside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void onUnbindLFSVInside(const IG &ig, 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 preAssembly()
Called directly before assembling.
void assembleVVolumePostSkeleton(const EG &eg, const LFSV &lfsv)
void onUnbindLFSVCoupling(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n, const LFSV_Coupling &lfsv_coupling)
void assembleUVVolumePostSkeleton(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
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)
void assembleVBoundary(const IG &ig, const LFSV_S &lfsv_s)
void onBindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void onBindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
LocalAssemblerInterface LocalAssembler
The type of the local assembler.
Definition common/assembler.hh:36
void onUnbindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assembleVVolume(const EG &eg, const LFSV &lfsv)
void onBindLFSVInside(const IG &ig, 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 onUnbindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void postAssembly()
Called last thing after assembling.
void assembleUVBoundary(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
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)
void loadCoefficientsLFSUOutside(const LFSU_N &lfsu_n)
void onBindLFSUVInside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
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)
The local assembler which provides the engines that drive the global assembler.
Definition common/assembler.hh:276
LocalResidualAssemblerEngine & localResidualAssemblerEngine(R &r, const X &x)
void preStep(TT time, TT dt, std::size_t stages)
Notify local assembler about upcoming time step.
void setWeight(RF weight)
Set current weight of assembling.
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(A &a, const X &x)
LocalPatternAssemblerEngine & localPatternAssemblerEngine(P &p)
void postStage()
Notify local assembler about completion of time step stage.
LocalResidualJacobianAssemblerEngine & localResidualJacobianAssemblerEngine(R &r, A &a, const X &x)
void postStep()
Notify local assembler about completion of time step.
void preStage(TT time, std::size_t stage)
Notify local assembler about upcoming time step stage.
void setTime(TT time)
Set current time of assembling.
TT suggestTimestep(TT dt) const
Suggest a valid time step size.
The grid operator represents an operator mapping which corresponds to the (possibly nonlinear) algebr...
Definition common/assembler.hh:346
GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, AssemblerInterface, LocalAssemblerInterface > Traits
The traits class.
Definition common/assembler.hh:351
void interpolate(const typename Traits::Domain &xold, const F &f, typename Traits::Domain &xnew)
Interpolate xnew from f, taking unconstrained values from xold.
GFSV::Traits::SizeType globalSizeV() const
GFSU::Traits::SizeType globalSizeU() const
void residual(const X &x, R &r) const
static void setupGridOperators(GridOperatorTuple &tuple)
LocalAssemblerInterface & localAssembler()
const GFSU & trialGridFunctionSpace() const
void jacobian(const X &x, A &a) const
const GFSV & testGridFunctionSpace() const
void fill_pattern(P &globalpattern) const
Determines the sparsity pattern of the jacobian matrix.
Base class for local assembler.
Definition assemblerutilities.hh:189
Traits class for the grid operator.
Definition gridoperatorutilities.hh:34