Go to the documentation of this file.
4 #ifndef DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
5 #define DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
7 #include <dune/common/exceptions.hh>
8 #include <dune/common/fmatrix.hh>
9 #include <dune/common/fvector.hh>
34 virtual unsigned s ()
const = 0;
39 virtual R
a (
int r,
int i)
const = 0;
44 virtual R
b (
int r,
int i)
const = 0;
49 virtual R
d (
int r)
const = 0;
53 virtual std::string
name ()
const = 0;
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;
93 virtual unsigned s ()
const override
101 virtual R
a (
int r,
int i)
const override
109 virtual R
b (
int r,
int i)
const override
117 virtual R
d (
int i)
const override
124 virtual std::string
name ()
const override
126 return std::string(
"one step theta");
131 Dune::FieldVector<R,2> D;
132 Dune::FieldMatrix<R,1,2> A;
133 Dune::FieldMatrix<R,1,2> B;
154 virtual std::string
name ()
const override
156 return std::string(
"explicit Euler");
179 virtual std::string
name ()
const override
181 return std::string(
"implicit Euler");
200 D[0] = 0.0; D[1] = 1.0; D[2] = 1.0;
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;
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;
218 virtual unsigned s ()
const override
226 virtual R
a (
int r,
int i)
const override
234 virtual R
b (
int r,
int i)
const override
242 virtual R
d (
int i)
const override
249 virtual std::string
name ()
const override
251 return std::string(
"Heun");
255 Dune::FieldVector<R,3> D;
256 Dune::FieldMatrix<R,2,3> A;
257 Dune::FieldMatrix<R,2,3> B;
273 D[0] = 0.0; D[1] = 1.0; D[2] = 0.5; D[3] = 1.0;
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;
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;
294 virtual unsigned s ()
const override
302 virtual R
a (
int r,
int i)
const override
310 virtual R
b (
int r,
int i)
const override
318 virtual R
d (
int i)
const override
325 virtual std::string
name ()
const override
327 return std::string(
"Shu's third order method");
331 Dune::FieldVector<R,4> D;
332 Dune::FieldMatrix<R,3,4> A;
333 Dune::FieldMatrix<R,3,4> B;
350 D[0] = 0.0; D[1] = 0.5; D[2] = 0.5; D[3] = 1.0; D[4] = 1.0;
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;
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;
373 virtual unsigned s ()
const override
381 virtual R
a (
int r,
int i)
const override
389 virtual R
b (
int r,
int i)
const override
397 virtual R
d (
int i)
const override
404 virtual std::string
name ()
const override
406 return std::string(
"RK4");
410 Dune::FieldVector<R,5> D;
411 Dune::FieldMatrix<R,4,5> A;
412 Dune::FieldMatrix<R,4,5> B;
431 alpha = 1.0 - 0.5*sqrt(2.0);
433 D[0] = 0.0; D[1] = alpha; D[2] = 1.0;
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;
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;
451 virtual unsigned s ()
const override
459 virtual R
a (
int r,
int i)
const override
467 virtual R
b (
int r,
int i)
const override
475 virtual R
d (
int i)
const override
482 virtual std::string
name ()
const override
484 return std::string(
"Alexander (order 2)");
489 Dune::FieldVector<R,3> D;
490 Dune::FieldMatrix<R,2,3> A;
491 Dune::FieldMatrix<R,2,3> B;
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);
513 D[0] = 0.0; D[1] = theta; D[2] = 1.0-theta; D[3] = 1.0;
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;
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;
533 virtual unsigned s ()
const override
541 virtual R
a (
int r,
int i)
const override
549 virtual R
b (
int r,
int i)
const override
557 virtual R
d (
int i)
const override
564 virtual std::string
name ()
const override
566 return std::string(
"Fractional step theta");
570 Dune::FieldVector<R,4> D;
571 Dune::FieldMatrix<R,3,4> A;
572 Dune::FieldMatrix<R,3,4> B;
589 R alpha = 0.4358665215;
592 for (
int i=1; i<=10; i++)
594 alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
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;
611 D[0] = 0.0; D[1] = alpha; D[2] = tau2; D[3] = 1.0;
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;
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;
631 virtual unsigned s ()
const override
639 virtual R
a (
int r,
int i)
const override
647 virtual R
b (
int r,
int i)
const override
655 virtual R
d (
int i)
const override
662 virtual std::string
name ()
const override
664 return std::string(
"Alexander (claims order 3)");
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;
676 #endif // DUNE_PDELAB_INSTATIONARY_ONESTEPPARAMETER_HH
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:564
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:467
Parameters to turn the ExplicitOneStepMethod into a classical fourth order Runge-Kutta method.
Definition: onestepparameter.hh:344
virtual unsigned s() const =0
Return number of stages of the method.
Parameters to turn the OneStepMethod into an implicit Euler method.
Definition: onestepparameter.hh:169
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:211
virtual bool implicit() const =0
Return true if method is implicit.
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:310
RK4Parameter()
Definition: onestepparameter.hh:348
Parameters to turn the OneStepMethod into an one step theta method.
Definition: onestepparameter.hh:69
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:242
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:639
Shu3Parameter()
Definition: onestepparameter.hh:271
Parameters to turn the ExplicitOneStepMethod into a third order strong stability preserving (SSP) sch...
Definition: onestepparameter.hh:267
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:444
virtual ~TimeSteppingParameterInterface()
every abstract base class has a virtual destructor
Definition: onestepparameter.hh:56
virtual std::string name() const =0
Return name of the scheme.
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:397
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:557
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:93
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:218
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:154
HeunParameter()
Definition: onestepparameter.hh:198
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:234
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:302
Parameters to turn the OneStepMethod into an Alexander scheme.
Definition: onestepparameter.hh:425
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:662
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:647
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:287
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:533
virtual R b(int r, int i) const =0
Return entries of the B matrix.
Parameters to turn the OneStepMethod into an Alexander3 scheme.
Definition: onestepparameter.hh:583
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:86
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:655
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:631
FractionalStepParameter()
Definition: onestepparameter.hh:505
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:404
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:459
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:294
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:482
Alexander3Parameter()
Definition: onestepparameter.hh:587
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:541
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:624
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:381
virtual R a(int r, int i) const =0
Return entries of the A matrix.
Parameters to turn the ExplicitOneStepMethod into an explicit Euler method.
Definition: onestepparameter.hh:144
R RealType
Definition: onestepparameter.hh:26
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:373
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:366
virtual unsigned s() const override
Return number of stages s of the method.
Definition: onestepparameter.hh:451
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:101
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:318
Parameters to turn the OneStepMethod into a fractional step theta scheme.
Definition: onestepparameter.hh:501
Alexander2Parameter()
Definition: onestepparameter.hh:429
OneStepThetaParameter(R theta_)
construct OneStepThetaParameter class
Definition: onestepparameter.hh:76
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:249
Parameters to turn the ExplicitOneStepMethod into a Heun scheme.
Definition: onestepparameter.hh:194
ImplicitEulerParameter()
Definition: onestepparameter.hh:173
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:124
virtual R a(int r, int i) const override
Return entries of the A matrix.
Definition: onestepparameter.hh:226
ExplicitEulerParameter()
Definition: onestepparameter.hh:148
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:549
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:475
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:109
virtual R d(int r) const =0
Return entries of the d Vector.
virtual bool implicit() const override
Return true if method is implicit.
Definition: onestepparameter.hh:526
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:179
virtual std::string name() const override
Return name of the scheme.
Definition: onestepparameter.hh:325
virtual R d(int i) const override
Return entries of the d Vector.
Definition: onestepparameter.hh:117
virtual R b(int r, int i) const override
Return entries of the B matrix.
Definition: onestepparameter.hh:389