dune-pdelab  2.7-git
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 
42  template<class R>
44  {
45  public:
46  typedef R RealType;
47 
50  virtual bool implicit () const = 0;
51 
54  virtual unsigned s () const = 0;
55 
59  virtual R a (int r, int i) const = 0;
60 
64  virtual R b (int r, int i) const = 0;
65 
69  virtual R d (int r) const = 0;
70 
73  virtual std::string name () const = 0;
74 
77  };
78 
79 
88  template<class R>
90  {
93 
94  public:
97  : theta(theta_)
98  {
99  D[0] = 0.0; D[1] = 1.0;
100  A[0][0] = -1.0; A[0][1] = 1.0;
101  B[0][0] = 1.0-theta; B[0][1] = theta;
102  }
103 
106  virtual bool implicit () const override
107  {
108  return (theta > 0.0);
109  }
110 
113  virtual unsigned s () const override
114  {
115  return 1;
116  }
117 
121  virtual R a (int r, int i) const override
122  {
123  return A[r-1][i];
124  }
125 
129  virtual R b (int r, int i) const override
130  {
131  return B[r-1][i];
132  }
133 
137  virtual R d (int i) const override
138  {
139  return D[i];
140  }
141 
144  virtual std::string name () const override
145  {
146  return std::string("one step theta");
147  }
148 
149  private:
150  R theta;
151  Dune::FieldVector<R,2> D;
152  Dune::FieldMatrix<R,1,2> A;
153  Dune::FieldMatrix<R,1,2> B;
154  };
155 
156 
163  template<class R>
165  {
166  public:
167 
169  : OneStepThetaParameter<R>(0.0)
170  {}
171 
174  virtual std::string name () const override
175  {
176  return std::string("explicit Euler");
177  }
178 
179  };
180 
181 
188  template<class R>
190  {
191  public:
192 
194  : OneStepThetaParameter<R>(1.0)
195  {}
196 
199  virtual std::string name () const override
200  {
201  return std::string("implicit Euler");
202  }
203 
204  };
205 
206 
213  template<class R>
215  {
216  public:
217 
219  {
220  D[0] = 0.0; D[1] = 1.0; D[2] = 1.0;
221 
222  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
223  A[1][0] = -0.5; A[1][1] = -0.5; A[1][2] = 1.0;
224 
225  B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0;
226  B[1][0] = 0.0; B[1][1] = 0.5; B[1][2] = 0.0;
227  }
228 
231  virtual bool implicit () const override
232  {
233  return false;
234  }
235 
238  virtual unsigned s () const override
239  {
240  return 2;
241  }
242 
246  virtual R a (int r, int i) const override
247  {
248  return A[r-1][i];
249  }
250 
254  virtual R b (int r, int i) const override
255  {
256  return B[r-1][i];
257  }
258 
262  virtual R d (int i) const override
263  {
264  return D[i];
265  }
266 
269  virtual std::string name () const override
270  {
271  return std::string("Heun");
272  }
273 
274  private:
275  Dune::FieldVector<R,3> D;
276  Dune::FieldMatrix<R,2,3> A;
277  Dune::FieldMatrix<R,2,3> B;
278  };
279 
286  template<class R>
288  {
289  public:
290 
292  {
293  D[0] = 0.0; D[1] = 1.0; D[2] = 0.5; D[3] = 1.0;
294 
295  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
296  A[1][0] = -0.75; A[1][1] = -0.25; A[1][2] = 1.0; A[1][3] = 0.0;
297  A[2][0] = -1.0/3.0; A[2][1] = 0.0; A[2][2] = -2.0/3.0; A[2][3] = 1.0;
298 
299  B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0; B[0][3] = 0.0;
300  B[1][0] = 0.0; B[1][1] = 0.25; B[1][2] = 0.0; B[1][3] = 0.0;
301  B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = 2.0/3.0; B[2][3] = 0.0;
302 
303  }
304 
307  virtual bool implicit () const override
308  {
309  return false;
310  }
311 
314  virtual unsigned s () const override
315  {
316  return 3;
317  }
318 
322  virtual R a (int r, int i) const override
323  {
324  return A[r-1][i];
325  }
326 
330  virtual R b (int r, int i) const override
331  {
332  return B[r-1][i];
333  }
334 
338  virtual R d (int i) const override
339  {
340  return D[i];
341  }
342 
345  virtual std::string name () const override
346  {
347  return std::string("Shu's third order method");
348  }
349 
350  private:
351  Dune::FieldVector<R,4> D;
352  Dune::FieldMatrix<R,3,4> A;
353  Dune::FieldMatrix<R,3,4> B;
354  };
355 
356 
363  template<class R>
365  {
366  public:
367 
369  {
370  D[0] = 0.0; D[1] = 0.5; D[2] = 0.5; D[3] = 1.0; D[4] = 1.0;
371 
372  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;
373  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;
374  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;
375  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;
376 
377  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;
378  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;
379  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;
380  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;
381 
382  }
383 
386  virtual bool implicit () const override
387  {
388  return false;
389  }
390 
393  virtual unsigned s () const override
394  {
395  return 4;
396  }
397 
401  virtual R a (int r, int i) const override
402  {
403  return A[r-1][i];
404  }
405 
409  virtual R b (int r, int i) const override
410  {
411  return B[r-1][i];
412  }
413 
417  virtual R d (int i) const override
418  {
419  return D[i];
420  }
421 
424  virtual std::string name () const override
425  {
426  return std::string("RK4");
427  }
428 
429  private:
430  Dune::FieldVector<R,5> D;
431  Dune::FieldMatrix<R,4,5> A;
432  Dune::FieldMatrix<R,4,5> B;
433  };
434 
435 
436 
437 
444  template<class R>
446  {
447  public:
448 
450  {
451  alpha = 1.0 - 0.5*sqrt(2.0);
452 
453  D[0] = 0.0; D[1] = alpha; D[2] = 1.0;
454 
455  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
456  A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0;
457 
458  B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0;
459  B[1][0] = 0.0; B[1][1] = 1.0-alpha; B[1][2] = alpha;
460  }
461 
464  virtual bool implicit () const override
465  {
466  return true;
467  }
468 
471  virtual unsigned s () const override
472  {
473  return 2;
474  }
475 
479  virtual R a (int r, int i) const override
480  {
481  return A[r-1][i];
482  }
483 
487  virtual R b (int r, int i) const override
488  {
489  return B[r-1][i];
490  }
491 
495  virtual R d (int i) const override
496  {
497  return D[i];
498  }
499 
502  virtual std::string name () const override
503  {
504  return std::string("Alexander (order 2)");
505  }
506 
507  private:
508  R alpha;
509  Dune::FieldVector<R,3> D;
510  Dune::FieldMatrix<R,2,3> A;
511  Dune::FieldMatrix<R,2,3> B;
512  };
513 
520  template<class R>
522  {
523  public:
524 
526  {
527  R alpha, theta, thetap, beta;
528  theta = 1.0 - 0.5*sqrt(2.0);
529  thetap = 1.0-2.0*theta;
530  alpha = 2.0-sqrt(2.0);
531  beta = 1.0-alpha;
532 
533  D[0] = 0.0; D[1] = theta; D[2] = 1.0-theta; D[3] = 1.0;
534 
535  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
536  A[1][0] = 0.0; A[1][1] = -1.0; A[1][2] = 1.0; A[1][3] = 0.0;
537  A[2][0] = 0.0; A[2][1] = 0.0; A[2][2] = -1.0; A[2][3] = 1.0;
538 
539  B[0][0] = beta*theta; B[0][1] = alpha*theta; B[0][2] = 0.0; B[0][3] = 0.0;
540  B[1][0] = 0.0; B[1][1] = alpha*thetap; B[1][2] = alpha*theta; B[1][3] = 0.0;
541  B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = beta*theta; B[2][3] = alpha*theta;
542  }
543 
546  virtual bool implicit () const override
547  {
548  return true;
549  }
550 
553  virtual unsigned s () const override
554  {
555  return 3;
556  }
557 
561  virtual R a (int r, int i) const override
562  {
563  return A[r-1][i];
564  }
565 
569  virtual R b (int r, int i) const override
570  {
571  return B[r-1][i];
572  }
573 
577  virtual R d (int i) const override
578  {
579  return D[i];
580  }
581 
584  virtual std::string name () const override
585  {
586  return std::string("Fractional step theta");
587  }
588 
589  private:
590  Dune::FieldVector<R,4> D;
591  Dune::FieldMatrix<R,3,4> A;
592  Dune::FieldMatrix<R,3,4> B;
593  };
594 
602  template<class R>
604  {
605  public:
606 
608  {
609  R alpha = 0.4358665215;
610 
611  // Newton iteration for alpha
612  for (int i=1; i<=10; i++)
613  {
614  alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
615  // std::cout.precision(16);
616  // std::cout << "alpha "
617  // << std::setw(8) << i << " "
618  // << std::scientific << alpha << std::endl;
619  }
620 
621  R tau2 = (1.0+alpha)*0.5;
622  R b1 = -(6.0*alpha*alpha -16.0*alpha + 1.0)*0.25;
623  R b2 = (6*alpha*alpha - 20.0*alpha + 5.0)*0.25;
624 
625  // std::cout.precision(16);
626  // std::cout << "alpha " << std::scientific << alpha << std::endl;
627  // std::cout << "tau2 " << std::scientific << tau2 << std::endl;
628  // std::cout << "b1 " << std::scientific << b1 << std::endl;
629  // std::cout << "b2 " << std::scientific << b2 << std::endl;
630 
631  D[0] = 0.0; D[1] = alpha; D[2] = tau2; D[3] = 1.0;
632 
633  A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
634  A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0; A[1][3] = 0.0;
635  A[2][0] = -1.0; A[2][1] = 0.0; A[2][2] = 0.0; A[2][3] = 1.0;
636 
637  B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0; B[0][3] = 0.0;
638  B[1][0] = 0.0; B[1][1] = tau2-alpha; B[1][2] = alpha; B[1][3] = 0.0;
639  B[2][0] = 0.0; B[2][1] = b1; B[2][2] = b2; B[2][3] = alpha;
640  }
641 
644  virtual bool implicit () const override
645  {
646  return true;
647  }
648 
651  virtual unsigned s () const override
652  {
653  return 3;
654  }
655 
659  virtual R a (int r, int i) const override
660  {
661  return A[r-1][i];
662  }
663 
667  virtual R b (int r, int i) const override
668  {
669  return B[r-1][i];
670  }
671 
675  virtual R d (int i) const override
676  {
677  return D[i];
678  }
679 
682  virtual std::string name () const override
683  {
684  return std::string("Alexander (claims order 3)");
685  }
686 
687  private:
688  R alpha, theta, thetap, beta;
689  Dune::FieldVector<R,4> D;
690  Dune::FieldMatrix<R,3,4> A;
691  Dune::FieldMatrix<R,3,4> B;
692  };
693 
694  } // end namespace PDELab
695 } // end namespace Dune
696 #endif // DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
Dune::PDELab::FractionalStepParameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:584
Dune::PDELab::Alexander2Parameter
Parameters to turn the OneStepMethod into an Alexander scheme.
Definition: onestepparameter.hh:445
Dune::PDELab::TimeSteppingParameterInterface::RealType
R RealType
Definition: onestepparameter.hh:46
Dune::PDELab::TimeSteppingParameterInterface::a
virtual R a(int r, int i) const =0
Return entries of the A matrix.
Dune::PDELab::HeunParameter
Parameters to turn the ExplicitOneStepMethod into a Heun scheme.
Definition: onestepparameter.hh:214
Dune::PDELab::OneStepThetaParameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:106
Dune::PDELab::Shu3Parameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:330
Dune::PDELab::OneStepThetaParameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:129
Dune::PDELab::Alexander2Parameter::Alexander2Parameter
Alexander2Parameter()
Definition: onestepparameter.hh:449
Dune::PDELab::Alexander3Parameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:682
Dune::PDELab::ExplicitEulerParameter
Parameters to turn the ExplicitOneStepMethod into an explicit Euler method.
Definition: onestepparameter.hh:164
Dune::PDELab::RK4Parameter
Parameters to turn the ExplicitOneStepMethod into a classical fourth order Runge-Kutta method.
Definition: onestepparameter.hh:364
Dune::PDELab::Shu3Parameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:307
Dune::PDELab::FractionalStepParameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:553
Dune::PDELab::TimeSteppingParameterInterface::implicit
virtual bool implicit() const =0
Return true if method is implicit.
Dune::PDELab::RK4Parameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:386
Dune::PDELab::Alexander2Parameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:471
Dune::PDELab::FractionalStepParameter::FractionalStepParameter
FractionalStepParameter()
Definition: onestepparameter.hh:525
Dune::PDELab::FractionalStepParameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:577
Dune::PDELab::HeunParameter::HeunParameter
HeunParameter()
Definition: onestepparameter.hh:218
Dune::PDELab::OneStepThetaParameter
Parameters to turn the OneStepMethod into an one step theta method.
Definition: onestepparameter.hh:89
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::RK4Parameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:424
Dune::PDELab::HeunParameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:246
Dune::PDELab::HeunParameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:254
Dune::PDELab::OneStepThetaParameter::OneStepThetaParameter
OneStepThetaParameter(R theta_)
construct OneStepThetaParameter class
Definition: onestepparameter.hh:96
Dune::PDELab::FractionalStepParameter
Parameters to turn the OneStepMethod into a fractional step theta scheme.
Definition: onestepparameter.hh:521
Dune::PDELab::ImplicitEulerParameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:199
Dune::PDELab::Alexander2Parameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:464
Dune::PDELab::FractionalStepParameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:546
Dune::PDELab::Alexander3Parameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:644
Dune::PDELab::Alexander3Parameter::Alexander3Parameter
Alexander3Parameter()
Definition: onestepparameter.hh:607
Dune::PDELab::Alexander3Parameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:659
Dune::PDELab::TimeSteppingParameterInterface::b
virtual R b(int r, int i) const =0
Return entries of the B matrix.
Dune::PDELab::HeunParameter::implicit
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:231
Dune::PDELab::ExplicitEulerParameter::ExplicitEulerParameter
ExplicitEulerParameter()
Definition: onestepparameter.hh:168
Dune::PDELab::OneStepThetaParameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:113
Dune::PDELab::ImplicitEulerParameter::ImplicitEulerParameter
ImplicitEulerParameter()
Definition: onestepparameter.hh:193
Dune::PDELab::Alexander2Parameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:487
Dune::PDELab::ExplicitEulerParameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:174
Dune::PDELab::Alexander2Parameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:502
Dune::PDELab::Alexander3Parameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:667
Dune::PDELab::TimeSteppingParameterInterface::d
virtual R d(int r) const =0
Return entries of the d Vector.
Dune::PDELab::RK4Parameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:401
Dune::PDELab::OneStepThetaParameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:144
Dune::PDELab::HeunParameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:269
Dune::PDELab::OneStepThetaParameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:137
Dune::PDELab::RK4Parameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:409
Dune::PDELab::Shu3Parameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:322
Dune::PDELab::RK4Parameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:393
Dune::PDELab::ImplicitEulerParameter
Parameters to turn the OneStepMethod into an implicit Euler method.
Definition: onestepparameter.hh:189
Dune::PDELab::RK4Parameter::RK4Parameter
RK4Parameter()
Definition: onestepparameter.hh:368
Dune::PDELab::TimeSteppingParameterInterface::s
virtual unsigned s() const =0
Return number of stages of the method.
Dune::PDELab::Shu3Parameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:314
Dune::PDELab::HeunParameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:262
Dune::PDELab::Shu3Parameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:338
Dune::PDELab::Alexander3Parameter
Parameters to turn the OneStepMethod into an Alexander3 scheme.
Definition: onestepparameter.hh:603
Dune::PDELab::Alexander3Parameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:675
Dune::PDELab::FractionalStepParameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:561
Dune::PDELab::TimeSteppingParameterInterface::name
virtual std::string name() const =0
Return name of the scheme.
Dune::PDELab::Shu3Parameter
Parameters to turn the ExplicitOneStepMethod into a third order strong stability preserving (SSP) sch...
Definition: onestepparameter.hh:287
Dune::PDELab::HeunParameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:238
Dune::PDELab::TimeSteppingParameterInterface
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:43
Dune::PDELab::Alexander2Parameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:479
Dune::PDELab::Alexander2Parameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:495
Dune::PDELab::FractionalStepParameter::b
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:569
Dune::PDELab::OneStepThetaParameter::a
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:121
Dune::PDELab::Alexander3Parameter::s
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:651
Dune::PDELab::Shu3Parameter::name
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:345
Dune::PDELab::RK4Parameter::d
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:417
Dune::PDELab::Shu3Parameter::Shu3Parameter
Shu3Parameter()
Definition: onestepparameter.hh:291
Dune::PDELab::TimeSteppingParameterInterface::~TimeSteppingParameterInterface
virtual ~TimeSteppingParameterInterface()
every abstract base class has a virtual destructor
Definition: onestepparameter.hh:76