dune-pdelab  2.4.1
onestepparameter.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil -*-
2 // vi: set et ts=2 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
5 #define DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
6 
7 #include <dune/common/exceptions.hh>
8 #include <dune/common/fmatrix.hh>
9 #include <dune/common/fvector.hh>
10 
11 namespace Dune {
12  namespace PDELab {
13 
18 
22  template<class R>
24  {
25  public:
26  typedef R RealType;
27 
30  virtual bool implicit () const = 0;
31 
34  virtual unsigned s () const = 0;
35 
39  virtual R a (int r, int i) const = 0;
40 
44  virtual R b (int r, int i) const = 0;
45 
49  virtual R d (int r) const = 0;
50 
53  virtual std::string name () const = 0;
54 
57  };
58 
59 
68  template<class R>
70  {
73 
74  public:
77  : theta(theta_)
78  {
79  D[0] = 0.0; D[1] = 1.0;
80  A[0][0] = -1.0; A[0][1] = 1.0;
81  B[0][0] = 1.0-theta; B[0][1] = theta;
82  }
83 
86  virtual bool implicit () const
87  {
88  return (theta > 0.0);
89  }
90 
93  virtual unsigned s () const
94  {
95  return 1;
96  }
97 
101  virtual R a (int r, int i) const
102  {
103  return A[r-1][i];
104  }
105 
109  virtual R b (int r, int i) const
110  {
111  return B[r-1][i];
112  }
113 
117  virtual R d (int i) const
118  {
119  return D[i];
120  }
121 
124  virtual std::string name () const
125  {
126  return std::string("one step theta");
127  }
128 
129  private:
130  R theta;
131  Dune::FieldVector<R,2> D;
132  Dune::FieldMatrix<R,1,2> A;
133  Dune::FieldMatrix<R,1,2> B;
134  };
135 
136 
143  template<class R>
145  {
146  public:
147 
149  : OneStepThetaParameter<R>(0.0)
150  {}
151 
154  virtual std::string name () const
155  {
156  return std::string("explicit Euler");
157  }
158 
159  };
160 
161 
168  template<class R>
170  {
171  public:
172 
174  : OneStepThetaParameter<R>(1.0)
175  {}
176 
179  virtual std::string name () const
180  {
181  return std::string("implicit Euler");
182  }
183 
184  };
185 
186 
193  template<class R>
195  {
196  public:
197 
199  {
200  D[0] = 0.0; D[1] = 1.0; D[2] = 1.0;
201 
202  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
203  A[1][0] = -0.5; A[1][1] = -0.5; A[1][2] = 1.0;
204 
205  B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0;
206  B[1][0] = 0.0; B[1][1] = 0.5; B[1][2] = 0.0;
207  }
208 
211  virtual bool implicit () const
212  {
213  return false;
214  }
215 
218  virtual unsigned s () const
219  {
220  return 2;
221  }
222 
226  virtual R a (int r, int i) const
227  {
228  return A[r-1][i];
229  }
230 
234  virtual R b (int r, int i) const
235  {
236  return B[r-1][i];
237  }
238 
242  virtual R d (int i) const
243  {
244  return D[i];
245  }
246 
249  virtual std::string name () const
250  {
251  return std::string("Heun");
252  }
253 
254  private:
255  Dune::FieldVector<R,3> D;
256  Dune::FieldMatrix<R,2,3> A;
257  Dune::FieldMatrix<R,2,3> B;
258  };
259 
266  template<class R>
268  {
269  public:
270 
272  {
273  D[0] = 0.0; D[1] = 1.0; D[2] = 0.5; D[3] = 1.0;
274 
275  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
276  A[1][0] = -0.75; A[1][1] = -0.25; A[1][2] = 1.0; A[1][3] = 0.0;
277  A[2][0] = -1.0/3.0; A[2][1] = 0.0; A[2][2] = -2.0/3.0; A[2][3] = 1.0;
278 
279  B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0; B[0][3] = 0.0;
280  B[1][0] = 0.0; B[1][1] = 0.25; B[1][2] = 0.0; B[1][3] = 0.0;
281  B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = 2.0/3.0; B[2][3] = 0.0;
282 
283  }
284 
287  virtual bool implicit () const
288  {
289  return false;
290  }
291 
294  virtual unsigned s () const
295  {
296  return 3;
297  }
298 
302  virtual R a (int r, int i) const
303  {
304  return A[r-1][i];
305  }
306 
310  virtual R b (int r, int i) const
311  {
312  return B[r-1][i];
313  }
314 
318  virtual R d (int i) const
319  {
320  return D[i];
321  }
322 
325  virtual std::string name () const
326  {
327  return std::string("Shu's third order method");
328  }
329 
330  private:
331  Dune::FieldVector<R,4> D;
332  Dune::FieldMatrix<R,3,4> A;
333  Dune::FieldMatrix<R,3,4> B;
334  };
335 
336 
343  template<class R>
345  {
346  public:
347 
349  {
350  D[0] = 0.0; D[1] = 0.5; D[2] = 0.5; D[3] = 1.0; D[4] = 1.0;
351 
352  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0; A[0][4] = 0.0;
353  A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0; A[1][3] = 0.0; A[1][4] = 0.0;
354  A[2][0] = -1.0; A[2][1] = 0.0; A[2][2] = 0.0; A[2][3] = 1.0; A[2][4] = 0.0;
355  A[3][0] = -1.0; A[3][1] = 0.0; A[3][2] = 0.0; A[3][3] = 0.0; A[3][4] = 1.0;
356 
357  B[0][0] = 0.5; B[0][1] = 0.0; B[0][2] = 0.0; B[0][3] = 0.0; B[0][4] = 0.0;
358  B[1][0] = 0.0; B[1][1] = 0.5; B[1][2] = 0.0; B[1][3] = 0.0; B[1][4] = 0.0;
359  B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = 1.0; B[2][3] = 0.0; B[2][4] = 0.0;
360  B[3][0] = 1.0/6.0; B[3][1] = 1.0/3.0; B[3][2] = 1.0/3.0; B[3][3] = 1.0/6.0; B[3][4] = 0.0;
361 
362  }
363 
366  virtual bool implicit () const
367  {
368  return false;
369  }
370 
373  virtual unsigned s () const
374  {
375  return 4;
376  }
377 
381  virtual R a (int r, int i) const
382  {
383  return A[r-1][i];
384  }
385 
389  virtual R b (int r, int i) const
390  {
391  return B[r-1][i];
392  }
393 
397  virtual R d (int i) const
398  {
399  return D[i];
400  }
401 
404  virtual std::string name () const
405  {
406  return std::string("RK4");
407  }
408 
409  private:
410  Dune::FieldVector<R,5> D;
411  Dune::FieldMatrix<R,4,5> A;
412  Dune::FieldMatrix<R,4,5> B;
413  };
414 
415 
416 
417 
424  template<class R>
426  {
427  public:
428 
430  {
431  alpha = 1.0 - 0.5*sqrt(2.0);
432 
433  D[0] = 0.0; D[1] = alpha; D[2] = 1.0;
434 
435  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
436  A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0;
437 
438  B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0;
439  B[1][0] = 0.0; B[1][1] = 1.0-alpha; B[1][2] = alpha;
440  }
441 
444  virtual bool implicit () const
445  {
446  return true;
447  }
448 
451  virtual unsigned s () const
452  {
453  return 2;
454  }
455 
459  virtual R a (int r, int i) const
460  {
461  return A[r-1][i];
462  }
463 
467  virtual R b (int r, int i) const
468  {
469  return B[r-1][i];
470  }
471 
475  virtual R d (int i) const
476  {
477  return D[i];
478  }
479 
482  virtual std::string name () const
483  {
484  return std::string("Alexander (order 2)");
485  }
486 
487  private:
488  R alpha;
489  Dune::FieldVector<R,3> D;
490  Dune::FieldMatrix<R,2,3> A;
491  Dune::FieldMatrix<R,2,3> B;
492  };
493 
500  template<class R>
502  {
503  public:
504 
506  {
507  R alpha, theta, thetap, beta;
508  theta = 1.0 - 0.5*sqrt(2.0);
509  thetap = 1.0-2.0*theta;
510  alpha = 2.0-sqrt(2.0);
511  beta = 1.0-alpha;
512 
513  D[0] = 0.0; D[1] = theta; D[2] = 1.0-theta; D[3] = 1.0;
514 
515  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
516  A[1][0] = 0.0; A[1][1] = -1.0; A[1][2] = 1.0; A[1][3] = 0.0;
517  A[2][0] = 0.0; A[2][1] = 0.0; A[2][2] = -1.0; A[2][3] = 1.0;
518 
519  B[0][0] = beta*theta; B[0][1] = alpha*theta; B[0][2] = 0.0; B[0][3] = 0.0;
520  B[1][0] = 0.0; B[1][1] = alpha*thetap; B[1][2] = alpha*theta; B[1][3] = 0.0;
521  B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = beta*theta; B[2][3] = alpha*theta;
522  }
523 
526  virtual bool implicit () const
527  {
528  return true;
529  }
530 
533  virtual unsigned s () const
534  {
535  return 3;
536  }
537 
541  virtual R a (int r, int i) const
542  {
543  return A[r-1][i];
544  }
545 
549  virtual R b (int r, int i) const
550  {
551  return B[r-1][i];
552  }
553 
557  virtual R d (int i) const
558  {
559  return D[i];
560  }
561 
564  virtual std::string name () const
565  {
566  return std::string("Fractional step theta");
567  }
568 
569  private:
570  Dune::FieldVector<R,4> D;
571  Dune::FieldMatrix<R,3,4> A;
572  Dune::FieldMatrix<R,3,4> B;
573  };
574 
582  template<class R>
584  {
585  public:
586 
588  {
589  R alpha = 0.4358665215;
590 
591  // Newton iteration for alpha
592  for (int i=1; i<=10; i++)
593  {
594  alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
595  // std::cout.precision(16);
596  // std::cout << "alpha "
597  // << std::setw(8) << i << " "
598  // << std::scientific << alpha << std::endl;
599  }
600 
601  R tau2 = (1.0+alpha)*0.5;
602  R b1 = -(6.0*alpha*alpha -16.0*alpha + 1.0)*0.25;
603  R b2 = (6*alpha*alpha - 20.0*alpha + 5.0)*0.25;
604 
605  // std::cout.precision(16);
606  // std::cout << "alpha " << std::scientific << alpha << std::endl;
607  // std::cout << "tau2 " << std::scientific << tau2 << std::endl;
608  // std::cout << "b1 " << std::scientific << b1 << std::endl;
609  // std::cout << "b2 " << std::scientific << b2 << std::endl;
610 
611  D[0] = 0.0; D[1] = alpha; D[2] = tau2; D[3] = 1.0;
612 
613  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
614  A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0; A[1][3] = 0.0;
615  A[2][0] = -1.0; A[2][1] = 0.0; A[2][2] = 0.0; A[2][3] = 1.0;
616 
617  B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0; B[0][3] = 0.0;
618  B[1][0] = 0.0; B[1][1] = tau2-alpha; B[1][2] = alpha; B[1][3] = 0.0;
619  B[2][0] = 0.0; B[2][1] = b1; B[2][2] = b2; B[2][3] = alpha;
620  }
621 
624  virtual bool implicit () const
625  {
626  return true;
627  }
628 
631  virtual unsigned s () const
632  {
633  return 3;
634  }
635 
639  virtual R a (int r, int i) const
640  {
641  return A[r-1][i];
642  }
643 
647  virtual R b (int r, int i) const
648  {
649  return B[r-1][i];
650  }
651 
655  virtual R d (int i) const
656  {
657  return D[i];
658  }
659 
662  virtual std::string name () const
663  {
664  return std::string("Alexander (claims order 3)");
665  }
666 
667  private:
668  R alpha, theta, thetap, beta;
669  Dune::FieldVector<R,4> D;
670  Dune::FieldMatrix<R,3,4> A;
671  Dune::FieldMatrix<R,3,4> B;
672  };
673 
674  } // end namespace PDELab
675 } // end namespace Dune
676 #endif
virtual R d(int r) const =0
Return entries of the d Vector.
Parameters to turn the OneStepMethod into an Alexander scheme.
Definition: onestepparameter.hh:425
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:482
Alexander2Parameter()
Definition: onestepparameter.hh:429
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:234
OneStepThetaParameter(R theta_)
construct OneStepThetaParameter class
Definition: onestepparameter.hh:76
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:211
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:475
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:109
Parameters to turn the ExplicitOneStepMethod into a Heun scheme.
Definition: onestepparameter.hh:194
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:242
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:154
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:249
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:662
Alexander3Parameter()
Definition: onestepparameter.hh:587
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:564
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:444
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:467
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:226
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:381
Definition: adaptivity.hh:27
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:318
virtual unsigned s() const =0
Return number of stages of the method.
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:218
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:631
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:389
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:451
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:302
virtual R b(int r, int i) const =0
Return entries of the B matrix.
Parameters to turn the ExplicitOneStepMethod into a third order strong stability preserving (SSP) sch...
Definition: onestepparameter.hh:267
Parameters to turn the OneStepMethod into an one step theta method.
Definition: onestepparameter.hh:69
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:639
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:404
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:101
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:179
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:549
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:541
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:655
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:557
RK4Parameter()
Definition: onestepparameter.hh:348
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:624
R RealType
Definition: onestepparameter.hh:26
HeunParameter()
Definition: onestepparameter.hh:198
Parameters to turn the OneStepMethod into an implicit Euler method.
Definition: onestepparameter.hh:169
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:647
ImplicitEulerParameter()
Definition: onestepparameter.hh:173
virtual bool implicit() const =0
Return true if method is implicit.
ExplicitEulerParameter()
Definition: onestepparameter.hh:148
Parameters to turn the OneStepMethod into an Alexander3 scheme.
Definition: onestepparameter.hh:583
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:459
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:310
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:93
virtual std::string name() const =0
Return name of the scheme.
virtual R a(int r, int i) const =0
Return entries of the A matrix.
virtual ~TimeSteppingParameterInterface()
every abstract base class has a virtual destructor
Definition: onestepparameter.hh:56
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:325
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:86
Shu3Parameter()
Definition: onestepparameter.hh:271
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:124
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:533
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:294
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:287
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:373
Parameters to turn the ExplicitOneStepMethod into a classical fourth order Runge-Kutta method...
Definition: onestepparameter.hh:344
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:526
FractionalStepParameter()
Definition: onestepparameter.hh:505
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:366
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:117
Parameters to turn the OneStepMethod into a fractional step theta scheme.
Definition: onestepparameter.hh:501
Parameters to turn the ExplicitOneStepMethod into an explicit Euler method.
Definition: onestepparameter.hh:144
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:397