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;
72 D[0] = 0.0; D[1] = 1.0;
73 A[0][0] = -1.0; A[0][1] = 1.0;
74 B[0][0] = 1.0; B[0][1] = 0.0;
86 virtual unsigned s ()
const
94 virtual R
a (
int r,
int i)
const
102 virtual R
b (
int r,
int i)
const
110 virtual R
d (
int i)
const
117 virtual std::string
name ()
const
119 return std::string(
"explicit Euler");
123 Dune::FieldVector<R,2> D;
124 Dune::FieldMatrix<R,1,2> A;
125 Dune::FieldMatrix<R,1,2> B;
147 D[0] = 0.0; D[1] = 1.0;
148 A[0][0] = -1.0; A[0][1] = 1.0;
149 B[0][0] = 1.0-theta; B[0][1] = theta;
164 virtual unsigned s ()
const
172 virtual R
a (
int r,
int i)
const
180 virtual R
b (
int r,
int i)
const
188 virtual R
d (
int i)
const
195 virtual std::string
name ()
const
197 return std::string(
"one step theta");
202 Dune::FieldVector<R,2> D;
203 Dune::FieldMatrix<R,1,2> A;
204 Dune::FieldMatrix<R,1,2> B;
220 D[0] = 0.0; D[1] = 1.0; D[2] = 1.0;
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;
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;
238 virtual unsigned s ()
const
246 virtual R
a (
int r,
int i)
const
254 virtual R
b (
int r,
int i)
const
262 virtual R
d (
int i)
const
269 virtual std::string
name ()
const
271 return std::string(
"Heun");
275 Dune::FieldVector<R,3> D;
276 Dune::FieldMatrix<R,2,3> A;
277 Dune::FieldMatrix<R,2,3> B;
293 D[0] = 0.0; D[1] = 1.0; D[2] = 0.5; D[3] = 1.0;
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;
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;
314 virtual unsigned s ()
const
322 virtual R
a (
int r,
int i)
const
330 virtual R
b (
int r,
int i)
const
338 virtual R
d (
int i)
const
345 virtual std::string
name ()
const
347 return std::string(
"Shu's third order method");
351 Dune::FieldVector<R,4> D;
352 Dune::FieldMatrix<R,3,4> A;
353 Dune::FieldMatrix<R,3,4> B;
370 D[0] = 0.0; D[1] = 0.5; D[2] = 0.5; D[3] = 1.0; D[4] = 1.0;
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;
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;
393 virtual unsigned s ()
const
401 virtual R
a (
int r,
int i)
const
409 virtual R
b (
int r,
int i)
const
417 virtual R
d (
int i)
const
424 virtual std::string
name ()
const
426 return std::string(
"RK4");
430 Dune::FieldVector<R,5> D;
431 Dune::FieldMatrix<R,4,5> A;
432 Dune::FieldMatrix<R,4,5> B;
451 alpha = 1.0 - 0.5*sqrt(2.0);
453 D[0] = 0.0; D[1] = alpha; D[2] = 1.0;
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;
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;
471 virtual unsigned s ()
const
479 virtual R
a (
int r,
int i)
const
487 virtual R
b (
int r,
int i)
const
495 virtual R
d (
int i)
const
502 virtual std::string
name ()
const
504 return std::string(
"Alexander (order 2)");
509 Dune::FieldVector<R,3> D;
510 Dune::FieldMatrix<R,2,3> A;
511 Dune::FieldMatrix<R,2,3> B;
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);
533 D[0] = 0.0; D[1] = theta; D[2] = 1.0-theta; D[3] = 1.0;
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;
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;
553 virtual unsigned s ()
const
561 virtual R
a (
int r,
int i)
const
569 virtual R
b (
int r,
int i)
const
577 virtual R
d (
int i)
const
584 virtual std::string
name ()
const
586 return std::string(
"Fractional step theta");
590 Dune::FieldVector<R,4> D;
591 Dune::FieldMatrix<R,3,4> A;
592 Dune::FieldMatrix<R,3,4> B;
609 R alpha = 0.4358665215;
612 for (
int i=1; i<=10; i++)
614 alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
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;
631 D[0] = 0.0; D[1] = alpha; D[2] = tau2; D[3] = 1.0;
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;
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;
651 virtual unsigned s ()
const
659 virtual R
a (
int r,
int i)
const
667 virtual R
b (
int r,
int i)
const
675 virtual R
d (
int i)
const
682 virtual std::string
name ()
const
684 return std::string(
"Alexander (claims order 3)");
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;
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:79
Parameters to turn the OneStepMethod into an Alexander3 scheme.
Definition: onestepparameter.hh:603
Alexander3Parameter()
Definition: onestepparameter.hh:607
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:659
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:110
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:409
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:188
Parameters to turn the OneStepMethod into an one step theta method.
Definition: onestepparameter.hh:137
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:502
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:417
FractionalStepParameter()
Definition: onestepparameter.hh:525
virtual std::string name() const =0
Return name of the scheme.
Parameters to turn the OneStepMethod into an Alexander scheme.
Definition: onestepparameter.hh:445
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:94
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:577
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:262
Base parameter class for time stepping scheme parameters.
Definition: onestepparameter.hh:23
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:651
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:424
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:314
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:495
virtual ~TimeSteppingParameterInterface()
every abstract base class has a virtual destructor
Definition: onestepparameter.hh:56
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:561
virtual unsigned s() const =0
Return number of stages of the method.
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:479
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:471
HeunParameter()
Definition: onestepparameter.hh:218
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:246
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:393
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:553
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:338
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:667
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:401
virtual R d(int i) const
Return entries of the d Vector.
Definition: onestepparameter.hh:675
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:86
Parameters to turn the ExplicitOneStepMethod into a classical fourth order Runge-Kutta method...
Definition: onestepparameter.hh:364
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:231
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:180
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:569
Parameters to turn the OneStepMethod into a fractional step theta scheme.
Definition: onestepparameter.hh:521
virtual R a(int r, int i) const =0
Return entries of the A matrix.
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:682
OneStepThetaParameter(R theta_)
construct OneStepThetaParameter class
Definition: onestepparameter.hh:144
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:172
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:307
Parameters to turn the ExplicitOneStepMethod into a third order strong stability preserving (SSP) sch...
Definition: onestepparameter.hh:287
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: onestepparameter.hh:322
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:238
Definition: adaptivity.hh:27
virtual bool implicit() const =0
Return true if method is implicit.
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:117
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:195
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:345
virtual unsigned s() const
Return number of stages s of the method.
Definition: onestepparameter.hh:164
RK4Parameter()
Definition: onestepparameter.hh:368
Parameters to turn the ExplicitOneStepMethod into an explicite Euler method.
Definition: onestepparameter.hh:66
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:269
virtual std::string name() const
Return name of the scheme.
Definition: onestepparameter.hh:584
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:386
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:330
virtual R d(int r) const =0
Return entries of the d Vector.
Parameters to turn the ExplicitOneStepMethod into a Heun scheme.
Definition: onestepparameter.hh:214
virtual R b(int r, int i) const =0
Return entries of the B matrix.
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:254
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:546
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:464
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:154
R RealType
Definition: onestepparameter.hh:26
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:487
ExplicitEulerParameter()
Definition: onestepparameter.hh:70
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: onestepparameter.hh:102
virtual bool implicit() const
Return true if method is implicit.
Definition: onestepparameter.hh:644
Shu3Parameter()
Definition: onestepparameter.hh:291
Alexander2Parameter()
Definition: onestepparameter.hh:449