DUNE PDELab (git)

prestageengine.hh
1#ifndef DUNE_PDELAB_GRIDOPERATOR_ONESTEP_PRESTAGEENGINE_HH
2#define DUNE_PDELAB_GRIDOPERATOR_ONESTEP_PRESTAGEENGINE_HH
3
4#include <dune/pdelab/gridoperator/onestep/enginebase.hh>
5#include <cmath>
6#include <vector>
7
8namespace Dune{
9 namespace PDELab{
10
18 template<typename OSLA>
19 class OneStepLocalPreStageAssemblerEngine
20 : public OneStepLocalAssemblerEngineBase<OSLA,
21 typename OSLA::LocalAssemblerDT0::LocalResidualAssemblerEngine,
22 typename OSLA::LocalAssemblerDT1::LocalResidualAssemblerEngine
23 >
24 {
25
26 typedef OneStepLocalAssemblerEngineBase<OSLA,
27 typename OSLA::LocalAssemblerDT0::LocalResidualAssemblerEngine,
28 typename OSLA::LocalAssemblerDT1::LocalResidualAssemblerEngine
29 > BaseT;
30
31 using BaseT::la;
32 using BaseT::lae0;
33 using BaseT::lae1;
34 using BaseT::implicit;
35 using BaseT::setLocalAssemblerEngineDT0;
36 using BaseT::setLocalAssemblerEngineDT1;
37
38 public:
40 typedef OSLA LocalAssembler;
41
44 typedef typename OSLA::LocalAssemblerDT0 LocalAssemblerDT0;
45 typedef typename OSLA::LocalAssemblerDT1 LocalAssemblerDT1;
46
47 typedef typename LocalAssemblerDT0::LocalResidualAssemblerEngine ResidualEngineDT0;
48 typedef typename LocalAssemblerDT1::LocalResidualAssemblerEngine ResidualEngineDT1;
50
52 typedef typename OSLA::Traits::Residual Residual;
53 typedef typename Residual::ElementType ResidualElement;
54
56 typedef typename OSLA::Traits::Solution Solution;
57 typedef typename Solution::ElementType SolutionElement;
58
60 typedef typename OSLA::Real Real;
61
63 typedef std::vector<Solution*> Solutions;
64
71 OneStepLocalPreStageAssemblerEngine(LocalAssembler & la_)
72 : BaseT(la_)
73 , invalid_residual(nullptr)
74 , invalid_solutions(nullptr)
75 , const_residual_0(invalid_residual)
76 , const_residual_1(invalid_residual)
77 , solutions(invalid_solutions)
78 {}
79
82 bool requireSkeleton() const
83 { return lae0->requireSkeleton() || lae1->requireSkeleton(); }
85
88 void setSolutions(const Solutions & solutions_)
89 {
90 solutions = &solutions_;
91 }
92
95 void setConstResiduals(Residual & const_residual_0_, Residual & const_residual_1_)
96 {
97 const_residual_0 = &const_residual_0_;
98 const_residual_1 = &const_residual_1_;
99
100 // Initialize the engines of the two wrapped local assemblers
101 assert(solutions != invalid_solutions);
102 setLocalAssemblerEngineDT0(la.la0.localResidualAssemblerEngine(*const_residual_0,*((*solutions)[0])));
103 setLocalAssemblerEngineDT1(la.la1.localResidualAssemblerEngine(*const_residual_1,*((*solutions)[0])));
104 }
105
108 void setConstResidual(Residual & const_residual_)
109 {
110 const_residual_0 = &const_residual_;
111 const_residual_1 = &const_residual_;
112
113 // Initialize the engines of the two wrapped local assemblers
114 assert(solutions != invalid_solutions);
115 setLocalAssemblerEngineDT0(la.la0.localResidualAssemblerEngine(*const_residual_0,*((*solutions)[0])));
116 setLocalAssemblerEngineDT1(la.la1.localResidualAssemblerEngine(*const_residual_1,*((*solutions)[0])));
117 }
118
125 template<typename EG, typename LFSU, typename LFSV>
126 void onBindLFSUV(const EG& eg, const LFSU& lfsu, const LFSV& lfsv)
127 {}
128 template<typename IG, typename LFSU_S, typename LFSV_S>
129 void onBindLFSUVInside(const IG & ig,
130 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s)
131 {}
132 template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N>
133 void onBindLFSUVOutside(const IG & ig,
134 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
135 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n)
136 {}
138
145 template<typename LFSU>
146 void loadCoefficientsLFSUInside(const LFSU & lfsu_s)
147 {}
148 template<typename LFSU>
149 void loadCoefficientsLFSUOutside(const LFSU & lfsu_n)
150 {}
151 template<typename LFSU>
152 void loadCoefficientsLFSUCoupling(const LFSU & lfsu_c)
153 {}
155
159 void setTimeInLastStage()
160 {
161 la.la1.setTime(la.time+la.osp_method->d(la.stage)*la.dt);
162 }
163
166 void preAssembly()
167 {
168 lae0->preAssembly();
169 lae1->preAssembly();
170
171 *const_residual_0 = 0.0;
172 *const_residual_1 = 0.0;
173
174 // Extract the coefficients of the time step scheme
175 a.resize(la.stage);
176 b.resize(la.stage);
177 d.resize(la.stage);
178 do0.resize(la.stage);
179 do1.resize(la.stage);
180 for (int i=0; i<la.stage; ++i){
181 a[i] = la.osp_method->a(la.stage,i);
182 b[i] = la.osp_method->b(la.stage,i);
183 d[i] = la.osp_method->d(i);
184 do0[i] = ( std::abs(b[i]) > 1E-6 );
185 do1[i] = ( std::abs(a[i]) > 1E-6 );
186 }
187
188 // prepare local operators for stage
189 la.la0.preStage(la.time+la.osp_method->d(la.stage)*la.dt,la.stage);
190 la.la1.preStage(la.time+la.osp_method->d(la.stage)*la.dt,la.stage);
191 }
192
193 template<typename GFSU, typename GFSV>
194 void postAssembly(const GFSU& gfsu, const GFSV& gfsv)
195 {
196 lae0->postAssembly(gfsu,gfsv);
197 lae1->postAssembly(gfsu,gfsv);
198 }
200
203
204 template<typename EG, typename LFSU, typename LFSV>
205 void assembleUVVolume(const EG & eg, const LFSU & lfsu, const LFSV & lfsv)
206 {
207 for (int s=0; s<la.stage; ++s)
208 {
209 // Reset the time in the local assembler
210 la.la0.setTime(la.time+d[s]*la.dt);
211 la.la1.setTime(la.time+d[s]*la.dt);
212
213 if(do0[s])
214 {
215 lae0->setSolution(*((*solutions)[s]));
216 lae0->onBindLFSUV(eg,lfsu,lfsv);
217 lae0->loadCoefficientsLFSUInside(lfsu);
218 la.la0.setWeight(b[s]*la.dt_factor0);
219 lae0->assembleUVVolume(eg,lfsu,lfsv);
220 }
221
222 if(do1[s])
223 {
224 lae1->setSolution(*((*solutions)[s]));
225 lae1->onBindLFSUV(eg,lfsu,lfsv);
226 lae1->loadCoefficientsLFSUInside(lfsu);
227 la.la1.setWeight(a[s]*la.dt_factor1);
228 lae1->assembleUVVolume(eg,lfsu,lfsv);
229 }
230 }
231 }
232
233 template<typename EG, typename LFSV>
234 void assembleVVolume(const EG & eg, const LFSV & lfsv)
235 {
236 for (int s=0; s<la.stage; ++s)
237 {
238 // Reset the time in the local assembler
239 la.la0.setTime(la.time+d[s]*la.dt);
240 la.la1.setTime(la.time+d[s]*la.dt);
241
242 if(do0[s])
243 {
244 la.la0.setWeight(b[s]*la.dt_factor0);
245 lae0->assembleVVolume(eg,lfsv);
246 }
247
248 if(do1[s])
249 {
250 la.la1.setWeight(a[s]*la.dt_factor1);
251 lae1->assembleVVolume(eg,lfsv);
252 }
253
254 }
255 }
256
257 template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N>
258 void assembleUVSkeleton(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
259 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n)
260 {
261 for (int s=0; s<la.stage; ++s)
262 {
263 // Reset the time in the local assembler
264 la.la0.setTime(la.time+d[s]*la.dt);
265 la.la1.setTime(la.time+d[s]*la.dt);
266
267 if(do0[s])
268 {
269 la.la0.setWeight(b[s]*la.dt_factor0);
270 lae0->setSolution(*((*solutions)[s]));
271 lae0->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
272 lae0->onBindLFSUVOutside(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n);
273 lae0->loadCoefficientsLFSUInside(lfsu_s);
274 lae0->loadCoefficientsLFSUOutside(lfsu_n);
275 la.la0.setWeight(b[s]*la.dt_factor0);
276 lae0->assembleUVSkeleton(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n);
277 }
278
279 if(do1[s])
280 {
281 la.la1.setWeight(a[s]*la.dt_factor1);
282 lae1->setSolution(*((*solutions)[s]));
283 lae1->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
284 lae1->onBindLFSUVOutside(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n);
285 lae1->loadCoefficientsLFSUInside(lfsu_s);
286 lae1->loadCoefficientsLFSUOutside(lfsu_n);
287 la.la1.setWeight(a[s]*la.dt_factor1);
288 lae1->assembleUVSkeleton(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n);
289 }
290 }
291 }
292
293 template<typename IG, typename LFSV_S, typename LFSV_N>
294 void assembleVSkeleton(const IG & ig, const LFSV_S & lfsv_s, const LFSV_N & lfsv_n)
295 {
296 for (int s=0; s<la.stage; ++s)
297 {
298 // Reset the time in the local assembler
299 la.la0.setTime(la.time+d[s]*la.dt);
300 la.la1.setTime(la.time+d[s]*la.dt);
301
302 if(do0[s])
303 {
304 la.la0.setWeight(b[s]*la.dt_factor0);
305 lae0->assembleVSkeleton(ig,lfsv_s,lfsv_n);
306 }
307
308 if(do1[s])
309 {
310 la.la1.setWeight(a[s]*la.dt_factor1);
311 lae1->assembleVSkeleton(ig,lfsv_s,lfsv_n);
312 }
313 }
314 }
315
316 template<typename IG, typename LFSU_S, typename LFSV_S>
317 void assembleUVBoundary(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s)
318 {
319 for (int s=0; s<la.stage; ++s)
320 {
321 // Reset the time in the local assembler
322 la.la0.setTime(la.time+d[s]*la.dt);
323 la.la1.setTime(la.time+d[s]*la.dt);
324
325 if(do0[s])
326 {
327 lae0->setSolution(*((*solutions)[s]));
328 lae0->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
329 lae0->loadCoefficientsLFSUInside(lfsu_s);
330 la.la0.setWeight(b[s]*la.dt_factor0);
331 lae0->assembleUVBoundary(ig,lfsu_s,lfsv_s);
332 }
333
334 if(do1[s])
335 {
336 lae1->setSolution(*((*solutions)[s]));
337 lae1->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
338 lae1->loadCoefficientsLFSUInside(lfsu_s);
339 la.la1.setWeight(a[s]*la.dt_factor1);
340 lae1->assembleUVBoundary(ig,lfsu_s,lfsv_s);
341 }
342 }
343 }
344
345 template<typename IG, typename LFSV_S>
346 void assembleVBoundary(const IG & ig, const LFSV_S & lfsv_s)
347 {
348 for (int s=0; s<la.stage; ++s)
349 {
350 // Reset the time in the local assembler
351 la.la0.setTime(la.time+d[s]*la.dt);
352 la.la1.setTime(la.time+d[s]*la.dt);
353
354 if(do0[s])
355 {
356 la.la0.setWeight(b[s]*la.dt_factor0);
357 lae0->assembleVBoundary(ig,lfsv_s);
358 }
359
360 if(do1[s])
361 {
362 la.la1.setWeight(a[s]*la.dt_factor1);
363 lae1->assembleVBoundary(ig,lfsv_s);
364 }
365 }
366 }
367
368 template<typename IG, typename LFSU_S, typename LFSV_S>
369 void assembleUVProcessor(const IG & ig, const LFSU_S & lfsu_s, const LFSV_S & lfsv_s)
370 {
371 for (int s=0; s<la.stage; ++s)
372 {
373 // Reset the time in the local assembler
374 la.la0.setTime(la.time+d[s]*la.dt);
375 la.la1.setTime(la.time+d[s]*la.dt);
376
377 if(do0[s])
378 {
379 lae0->setSolution(*((*solutions)[s]));
380 lae0->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
381 lae0->loadCoefficientsLFSUInside(lfsu_s);
382 la.la0.setWeight(b[s]*la.dt_factor0);
383 lae0->assembleUVProcessor(ig,lfsu_s,lfsv_s);
384 }
385
386 if(do1[s])
387 {
388 lae1->setSolution(*((*solutions)[s]));
389 lae1->onBindLFSUVInside(ig,lfsu_s,lfsv_s);
390 lae1->loadCoefficientsLFSUInside(lfsu_s);
391 la.la1.setWeight(a[s]*la.dt_factor1);
392 lae1->assembleUVProcessor(ig,lfsu_s,lfsv_s);
393 }
394 }
395 }
396
397 template<typename IG, typename LFSV_S>
398 void assembleVProcessor(const IG & ig, const LFSV_S & lfsv_s)
399 {
400 for (int s=0; s<la.stage; ++s)
401 {
402 // Reset the time in the local assembler
403 la.la0.setTime(la.time+d[s]*la.dt);
404 la.la1.setTime(la.time+d[s]*la.dt);
405
406 if(do0[s])
407 {
408 la.la0.setWeight(b[s]*la.dt_factor0);
409 lae0->assembleVProcessor(ig,lfsv_s);
410 }
411
412 if(do1[s])
413 {
414 la.la1.setWeight(a[s]*la.dt_factor1);
415 lae1->assembleVProcessor(ig,lfsv_s);
416 }
417 }
418 }
419
420 template<typename IG, typename LFSU_S, typename LFSV_S, typename LFSU_N, typename LFSV_N,
421 typename LFSU_C, typename LFSV_C>
422 void assembleUVEnrichedCoupling(const IG & ig,
423 const LFSU_S & lfsu_s, const LFSV_S & lfsv_s,
424 const LFSU_N & lfsu_n, const LFSV_N & lfsv_n,
425 const LFSU_C & lfsu_c, const LFSV_C & lfsv_c)
426 {
427 for (int s=0; s<la.stage; ++s)
428 {
429 // Reset the time in the local assembler
430 la.la0.setTime(la.time+d[s]*la.dt);
431 la.la1.setTime(la.time+d[s]*la.dt);
432
433 if(do0[s])
434 {
435 lae0->setSolution(*((*solutions)[s]));
436 lae0->loadCoefficientsLFSUInside(lfsu_s);
437 lae0->loadCoefficientsLFSUOutside(lfsu_n);
438 lae0->loadCoefficientsLFSUCoupling(lfsu_c);
439 la.la0.setWeight(b[s]*la.dt_factor0);
440 lae0->assembleUVEnrichedCoupling(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n,lfsu_c,lfsv_c);
441 }
442
443 if(do1[s])
444 {
445 lae1->setSolution(*((*solutions)[s]));
446 lae1->loadCoefficientsLFSUInside(lfsu_s);
447 lae1->loadCoefficientsLFSUOutside(lfsu_n);
448 lae1->loadCoefficientsLFSUCoupling(lfsu_c);
449 la.la1.setWeight(a[s]*la.dt_factor1);
450 lae1->assembleUVEnrichedCoupling(ig,lfsu_s,lfsv_s,lfsu_n,lfsv_n,lfsu_c,lfsv_c);
451 }
452 }
453 }
454
455 template<typename IG, typename LFSV_S, typename LFSV_N, typename LFSV_C>
456 void assembleVEnrichedCoupling(const IG & ig,
457 const LFSV_S & lfsv_s,
458 const LFSV_N & lfsv_n,
459 const LFSV_C & lfsv_c)
460 {
461 for (int s=0; s<la.stage; ++s)
462 {
463 // Reset the time in the local assembler
464 la.la0.setTime(la.time+d[s]*la.dt);
465 la.la1.setTime(la.time+d[s]*la.dt);
466
467 if(do0[s])
468 {
469 la.la0.setWeight(b[s]*la.dt_factor0);
470 lae0->assembleVEnrichedCoupling(ig,lfsv_s,lfsv_n,lfsv_c);
471 }
472
473 if(do1[s])
474 {
475 la.la1.setWeight(a[s]*la.dt_factor1);
476 lae1->assembleVEnrichedCoupling(ig,lfsv_s,lfsv_n,lfsv_c);
477 }
478
479 }
480 }
481
482 template<typename EG, typename LFSU, typename LFSV>
483 void assembleUVVolumePostSkeleton(const EG & eg, const LFSU & lfsu, const LFSV & lfsv)
484 {
485 for (int s=0; s<la.stage; ++s)
486 {
487 // Reset the time in the local assembler
488 la.la0.setTime(la.time+d[s]*la.dt);
489 la.la1.setTime(la.time+d[s]*la.dt);
490
491 if(do0[s])
492 {
493 lae0->setSolution(*((*solutions)[s]));
494 lae0->onBindLFSUV(eg,lfsu,lfsv);
495 lae0->loadCoefficientsLFSUInside(lfsu);
496 la.la0.setWeight(b[s]*la.dt_factor0);
497 lae0->assembleUVVolumePostSkeleton(eg,lfsu,lfsv);
498 }
499
500 if(do1[s])
501 {
502 lae1->setSolution(*((*solutions)[s]));
503 lae1->onBindLFSUV(eg,lfsu,lfsv);
504 lae1->loadCoefficientsLFSUInside(lfsu);
505 la.la1.setWeight(a[s]*la.dt_factor1);
506 lae1->assembleUVVolumePostSkeleton(eg,lfsu,lfsv);
507 }
508
509 }
510 }
511
512 template<typename EG, typename LFSV>
513 void assembleVVolumePostSkeleton(const EG & eg, const LFSV & lfsv)
514 {
515 for (int s=0; s<la.stage; ++s)
516 {
517 // Reset the time in the local assembler
518 la.la0.setTime(la.time+d[s]*la.dt);
519 la.la1.setTime(la.time+d[s]*la.dt);
520
521 if(do0[s])
522 {
523 la.la0.setWeight(b[s]*la.dt_factor0);
524 lae0->assembleVVolumePostSkeleton(eg,lfsv);
525 }
526
527 if(do1[s])
528 {
529 la.la1.setWeight(a[s]*la.dt_factor1);
530 lae1->assembleVVolumePostSkeleton(eg,lfsv);
531 }
532 }
533 }
535
536 private:
537
539 Residual * const invalid_residual;
540
542 Solutions * const invalid_solutions;
543
548 Residual * const_residual_0;
549 Residual * const_residual_1;
551
553 const Solutions * solutions;
554
556 std::vector<Real> a;
557 std::vector<Real> b;
558 std::vector<Real> d;
559 std::vector<bool> do0;
560 std::vector<bool> do1;
561
562 }; // End of class OneStepLocalPreStageAssemblerEngine
563
564 }
565}
566#endif // DUNE_PDELAB_GRIDOPERATOR_ONESTEP_PRESTAGEENGINE_HH
Dune::PDELab::OneStepLocalAssemblerEngineBase
The local assembler engine for UDG sub triangulations which assembles the residual vector.
Definition: enginebase.hh:16
Dune::PDELab::OneStepLocalAssemblerEngineBase::la
const LocalAssembler & la
Definition: enginebase.hh:472
Dune::PDELab::OneStepLocalPreStageAssemblerEngine
The local assembler engine for one step methods which assembles the constant part of the residual vec...
Definition: prestageengine.hh:24
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::setSolutions
void setSolutions(const Solutions &solutions_)
Definition: prestageengine.hh:88
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::Solutions
std::vector< Solution * > Solutions
The type of the solution container.
Definition: prestageengine.hh:63
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::Solution
OSLA::Traits::Solution Solution
The type of the solution vector.
Definition: prestageengine.hh:56
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::loadCoefficientsLFSUInside
void loadCoefficientsLFSUInside(const LFSU &lfsu_s)
Definition: prestageengine.hh:146
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::assembleUVVolume
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
Definition: prestageengine.hh:205
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::OneStepLocalPreStageAssemblerEngine
OneStepLocalPreStageAssemblerEngine(LocalAssembler &la_)
Constructor.
Definition: prestageengine.hh:71
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::Real
OSLA::Real Real
The type for real numbers.
Definition: prestageengine.hh:60
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::onBindLFSUV
void onBindLFSUV(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
Definition: prestageengine.hh:126
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::requireSkeleton
bool requireSkeleton() const
Definition: prestageengine.hh:82
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::Residual
OSLA::Traits::Residual Residual
The type of the residual vector.
Definition: prestageengine.hh:52
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::preAssembly
void preAssembly()
Definition: prestageengine.hh:166
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::setConstResiduals
void setConstResiduals(Residual &const_residual_0_, Residual &const_residual_1_)
Definition: prestageengine.hh:95
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::LocalAssembler
OSLA LocalAssembler
The type of the wrapping local assembler.
Definition: prestageengine.hh:40
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::LocalAssemblerDT0
OSLA::LocalAssemblerDT0 LocalAssemblerDT0
Definition: prestageengine.hh:44
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::setTimeInLastStage
void setTimeInLastStage()
Definition: prestageengine.hh:159
Dune::PDELab::OneStepLocalPreStageAssemblerEngine::setConstResidual
void setConstResidual(Residual &const_residual_)
Definition: prestageengine.hh:108
Dune
Dune namespace.
Definition: alignedallocator.hh:13
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 23, 23:29, 2024)