4 #ifndef DUNE_PDELAB_ONESTEP_HH
5 #define DUNE_PDELAB_ONESTEP_HH
14 #include <dune/common/exceptions.hh>
15 #include <dune/common/fmatrix.hh>
16 #include <dune/common/fvector.hh>
17 #include <dune/common/ios_state.hh>
42 D[0] = 0.0; D[1] = 1.0;
43 A[0][0] = -1.0; A[0][1] = 1.0;
44 B[0][0] = 1.0; B[0][1] = 0.0;
56 virtual unsigned s ()
const
64 virtual R
a (
int r,
int i)
const
72 virtual R
b (
int r,
int i)
const
80 virtual R
d (
int i)
const
87 virtual std::string
name ()
const
89 return std::string(
"explicit Euler");
93 Dune::FieldVector<R,2> D;
94 Dune::FieldMatrix<R,1,2> A;
95 Dune::FieldMatrix<R,1,2> B;
117 D[0] = 0.0; D[1] = 1.0;
118 A[0][0] = -1.0; A[0][1] = 1.0;
119 B[0][0] = 1.0-theta; B[0][1] = theta;
134 virtual unsigned s ()
const
142 virtual R
a (
int r,
int i)
const
150 virtual R
b (
int r,
int i)
const
158 virtual R
d (
int i)
const
165 virtual std::string
name ()
const
167 return std::string(
"one step theta");
172 Dune::FieldVector<R,2> D;
173 Dune::FieldMatrix<R,1,2> A;
174 Dune::FieldMatrix<R,1,2> B;
190 D[0] = 0.0; D[1] = 1.0; D[2] = 1.0;
192 A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
193 A[1][0] = -0.5; A[1][1] = -0.5; A[1][2] = 1.0;
195 B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0;
196 B[1][0] = 0.0; B[1][1] = 0.5; B[1][2] = 0.0;
208 virtual unsigned s ()
const
216 virtual R
a (
int r,
int i)
const
224 virtual R
b (
int r,
int i)
const
232 virtual R
d (
int i)
const
239 virtual std::string
name ()
const
241 return std::string(
"Heun");
245 Dune::FieldVector<R,3> D;
246 Dune::FieldMatrix<R,2,3> A;
247 Dune::FieldMatrix<R,2,3> B;
263 D[0] = 0.0; D[1] = 1.0; D[2] = 0.5; D[3] = 1.0;
265 A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
266 A[1][0] = -0.75; A[1][1] = -0.25; A[1][2] = 1.0; A[1][3] = 0.0;
267 A[2][0] = -1.0/3.0; A[2][1] = 0.0; A[2][2] = -2.0/3.0; A[2][3] = 1.0;
269 B[0][0] = 1.0; B[0][1] = 0.0; B[0][2] = 0.0; B[0][3] = 0.0;
270 B[1][0] = 0.0; B[1][1] = 0.25; B[1][2] = 0.0; B[1][3] = 0.0;
271 B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = 2.0/3.0; B[2][3] = 0.0;
284 virtual unsigned s ()
const
292 virtual R
a (
int r,
int i)
const
300 virtual R
b (
int r,
int i)
const
308 virtual R
d (
int i)
const
315 virtual std::string
name ()
const
317 return std::string(
"Shu's third order method");
321 Dune::FieldVector<R,4> D;
322 Dune::FieldMatrix<R,3,4> A;
323 Dune::FieldMatrix<R,3,4> B;
340 D[0] = 0.0; D[1] = 0.5; D[2] = 0.5; D[3] = 1.0; D[4] = 1.0;
342 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;
343 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;
344 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;
345 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;
347 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;
348 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;
349 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;
350 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;
363 virtual unsigned s ()
const
371 virtual R
a (
int r,
int i)
const
379 virtual R
b (
int r,
int i)
const
387 virtual R
d (
int i)
const
394 virtual std::string
name ()
const
396 return std::string(
"RK4");
400 Dune::FieldVector<R,5> D;
401 Dune::FieldMatrix<R,4,5> A;
402 Dune::FieldMatrix<R,4,5> B;
421 alpha = 1.0 - 0.5*sqrt(2.0);
423 D[0] = 0.0; D[1] = alpha; D[2] = 1.0;
425 A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0;
426 A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0;
428 B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0;
429 B[1][0] = 0.0; B[1][1] = 1.0-alpha; B[1][2] = alpha;
441 virtual unsigned s ()
const
449 virtual R
a (
int r,
int i)
const
457 virtual R
b (
int r,
int i)
const
465 virtual R
d (
int i)
const
472 virtual std::string
name ()
const
474 return std::string(
"Alexander (order 2)");
479 Dune::FieldVector<R,3> D;
480 Dune::FieldMatrix<R,2,3> A;
481 Dune::FieldMatrix<R,2,3> B;
497 R alpha, theta, thetap, beta;
498 theta = 1.0 - 0.5*sqrt(2.0);
499 thetap = 1.0-2.0*theta;
500 alpha = 2.0-sqrt(2.0);
503 D[0] = 0.0; D[1] = theta; D[2] = 1.0-theta; D[3] = 1.0;
505 A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
506 A[1][0] = 0.0; A[1][1] = -1.0; A[1][2] = 1.0; A[1][3] = 0.0;
507 A[2][0] = 0.0; A[2][1] = 0.0; A[2][2] = -1.0; A[2][3] = 1.0;
509 B[0][0] = beta*theta; B[0][1] = alpha*theta; B[0][2] = 0.0; B[0][3] = 0.0;
510 B[1][0] = 0.0; B[1][1] = alpha*thetap; B[1][2] = alpha*theta; B[1][3] = 0.0;
511 B[2][0] = 0.0; B[2][1] = 0.0; B[2][2] = beta*theta; B[2][3] = alpha*theta;
523 virtual unsigned s ()
const
531 virtual R
a (
int r,
int i)
const
539 virtual R
b (
int r,
int i)
const
547 virtual R
d (
int i)
const
554 virtual std::string
name ()
const
556 return std::string(
"Fractional step theta");
560 Dune::FieldVector<R,4> D;
561 Dune::FieldMatrix<R,3,4> A;
562 Dune::FieldMatrix<R,3,4> B;
579 R alpha = 0.4358665215;
582 for (
int i=1; i<=10; i++)
584 alpha = alpha - (alpha*(alpha*alpha-3.0*(alpha-0.5))-1.0/6.0)/(3.0*alpha*(alpha-2.0)+1.5);
591 R tau2 = (1.0+alpha)*0.5;
592 R b1 = -(6.0*alpha*alpha -16.0*alpha + 1.0)*0.25;
593 R b2 = (6*alpha*alpha - 20.0*alpha + 5.0)*0.25;
601 D[0] = 0.0; D[1] = alpha; D[2] = tau2; D[3] = 1.0;
603 A[0][0] = -1.0; A[0][1] = 1.0; A[0][2] = 0.0; A[0][3] = 0.0;
604 A[1][0] = -1.0; A[1][1] = 0.0; A[1][2] = 1.0; A[1][3] = 0.0;
605 A[2][0] = -1.0; A[2][1] = 0.0; A[2][2] = 0.0; A[2][3] = 1.0;
607 B[0][0] = 0.0; B[0][1] = alpha; B[0][2] = 0.0; B[0][3] = 0.0;
608 B[1][0] = 0.0; B[1][1] = tau2-alpha; B[1][2] = alpha; B[1][3] = 0.0;
609 B[2][0] = 0.0; B[2][1] = b1; B[2][2] = b2; B[2][3] = alpha;
621 virtual unsigned s ()
const
629 virtual R
a (
int r,
int i)
const
637 virtual R
b (
int r,
int i)
const
645 virtual R
d (
int i)
const
652 virtual std::string
name ()
const
654 return std::string(
"Alexander (claims order 3)");
658 R alpha, theta, thetap, beta;
659 Dune::FieldVector<R,4> D;
660 Dune::FieldMatrix<R,3,4> A;
661 Dune::FieldMatrix<R,3,4> B;
707 template<
class R,
class IGOS>
716 CFLTimeController (R cfl_, R target_,
const IGOS& igos_) : cfl(cfl_), target(target_), igos(igos_)
728 RealType suggested = cfl*igos.suggestTimestep(givendt);
729 if (time+2.0*suggested<target)
731 if (time+suggested<target)
732 return 0.5*(target-time);
777 template<
class T,
class IGOS,
class PDESOLVER,
class TrlV,
class TstV = TrlV>
780 typedef typename PDESOLVER::Result PDESolverResult;
798 IGOS& igos_, PDESOLVER& pdesolver_)
799 : method(&method_), igos(igos_), pdesolver(pdesolver_), verbosityLevel(1), step(1), res()
801 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
808 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
811 verbosityLevel = level;
854 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
857 ios_base_all_saver format_attribute_saver(std::cout);
862 std::vector<TrlV*> x(1);
865 if (verbosityLevel>=1){
866 std::ios_base::fmtflags oldflags = std::cout.flags();
867 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
868 << std::setw(6) << step
870 << std::setw(12) << std::setprecision(4) << std::scientific
873 << std::setw(12) << std::setprecision(4) << std::scientific
876 << std::setw(12) << std::setprecision(4) << std::scientific
879 std::cout.flags(oldflags);
883 igos.preStep(*method,time,dt);
886 for (
unsigned r=1; r<=method->
s(); ++r)
888 if (verbosityLevel>=2){
889 std::ios_base::fmtflags oldflags = std::cout.flags();
890 std::cout <<
"STAGE "
893 << std::setw(12) << std::setprecision(4) << std::scientific
894 << time+method->
d(r)*dt
896 std::cout.flags(oldflags);
907 if (r>1) xnew = *(x[r-1]);
912 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
921 pdesolver.apply(*x[r]);
926 PDESolverResult pderes = pdesolver.result();
938 PDESolverResult pderes = pdesolver.result();
949 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
965 if (verbosityLevel>=1){
966 std::ios_base::fmtflags oldflags = std::cout.flags();
973 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
975 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
977 std::cout.flags(oldflags);
995 T
apply (T time, T dt, TrlV& xold, F&
f, TrlV& xnew)
1001 ios_base_all_saver format_attribute_saver(std::cout);
1003 std::vector<TrlV*> x(1);
1006 if (verbosityLevel>=1){
1007 std::ios_base::fmtflags oldflags = std::cout.flags();
1008 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
1009 << std::setw(6) << step
1011 << std::setw(12) << std::setprecision(4) << std::scientific
1014 << std::setw(12) << std::setprecision(4) << std::scientific
1017 << std::setw(12) << std::setprecision(4) << std::scientific
1020 std::cout.flags(oldflags);
1024 igos.preStep(*method,time,dt);
1027 for (
unsigned r=1; r<=method->
s(); ++r)
1029 if (verbosityLevel>=2){
1030 std::ios_base::fmtflags oldflags = std::cout.flags();
1031 std::cout <<
"STAGE "
1034 << std::setw(12) << std::setprecision(4) << std::scientific
1035 << time+method->
d(r)*dt
1036 <<
"." << std::endl;
1037 std::cout.flags(oldflags);
1052 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
1056 igos.interpolate(r,*x[r-1],f,*x[r]);
1060 pdesolver.apply(*x[r]);
1065 PDESolverResult pderes = pdesolver.result();
1077 PDESolverResult pderes = pdesolver.result();
1088 for (
unsigned i=1; i<method->
s(); ++i)
delete x[i];
1104 if (verbosityLevel>=1){
1105 std::ios_base::fmtflags oldflags = std::cout.flags();
1112 std::cout <<
"::: assemble time " << std::setw(12) << std::setprecision(4) << std::scientific
1114 std::cout <<
"::: lin solve time " << std::setw(12) << std::setprecision(4) << std::scientific
1116 std::cout.flags(oldflags);
1126 PDESOLVER& pdesolver;
1141 template<
class T,
class IGOS,
class LS,
class TrlV,
class TstV = TrlV,
class TC = SimpleTimeController<T> >
1144 typedef typename TrlV::ElementType Real;
1145 typedef typename IGOS::template MatrixContainer<Real>::Type M;
1161 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
1165 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
1166 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
1183 : method(&method_), igos(igos_), ls(ls_), verbosityLevel(1), step(1), D(igos),
1184 tc(&tc_), allocated(false)
1187 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
1192 if (allocated)
delete tc;
1198 if (igos.trialGridFunctionSpace().gridView().comm().rank()>0)
1201 verbosityLevel = level;
1219 DUNE_THROW(
Exception,
"explicit one step method called with implicit scheme");
1230 T
apply (T time, T dt, TrlV& xold, TrlV& xnew)
1232 DefaultLimiter limiter;
1233 return apply(time,dt,xold,xnew,limiter);
1236 template<
typename Limiter>
1237 T
apply (T time, T dt, TrlV& xold, TrlV& xnew, Limiter& limiter)
1240 ios_base_all_saver format_attribute_saver(std::cout);
1242 mytag <<
"ExplicitOneStepMethod::apply(): ";
1244 std::vector<TrlV*> x(1);
1246 if(verbosityLevel>=4)
1247 std::cout << mytag <<
"Creating residual vectors alpha and beta..."
1249 TstV alpha(igos.testGridFunctionSpace()), beta(igos.testGridFunctionSpace());
1250 if(verbosityLevel>=4)
1252 <<
"Creating residual vectors alpha and beta... done."
1255 if (verbosityLevel>=1){
1256 std::ios_base::fmtflags oldflags = std::cout.flags();
1257 std::cout <<
"TIME STEP [" << method->
name() <<
"] "
1258 << std::setw(6) << step
1260 << std::setw(12) << std::setprecision(4) << std::scientific
1263 << std::setw(12) << std::setprecision(4) << std::scientific
1266 << std::setw(12) << std::setprecision(4) << std::scientific
1269 std::cout.flags(oldflags);
1273 if(verbosityLevel>=4)
1274 std::cout << mytag <<
"Preparing assembler..." << std::endl;
1275 igos.preStep(*method,time,dt);
1276 if(verbosityLevel>=4)
1277 std::cout << mytag <<
"Preparing assembler... done." << std::endl;
1280 for(
unsigned r=1; r<=method->
s(); ++r)
1283 stagetag <<
"stage " << r <<
": ";
1284 if (verbosityLevel>=4)
1285 std::cout << stagetag <<
"Start." << std::endl;
1287 if (verbosityLevel>=2){
1288 std::ios_base::fmtflags oldflags = std::cout.flags();
1289 std::cout <<
"STAGE "
1292 << std::setw(12) << std::setprecision(4) << std::scientific
1293 << time+method->
d(r)*dt
1294 <<
"." << std::endl;
1295 std::cout.flags(oldflags);
1303 if (r>1) xnew = *(x[r-1]);
1308 x.push_back(
new TrlV(igos.trialGridFunctionSpace()));
1310 *(x[r]) = *(x[r-1]);
1316 if (verbosityLevel>=4) std::cout <<
"assembling D, alpha, beta ..." << std::endl;
1322 limiter.prestage(*x[r-1]);
1324 if(verbosityLevel>=4)
1325 std::cout << stagetag <<
"Assembling residual..." << std::endl;
1326 igos.explicit_jacobian_residual(r,x,D,alpha,beta);
1327 if(verbosityLevel>=4)
1328 std::cout << stagetag <<
"Assembling residual... done."
1335 newdt = std::min(newdt, dt);
1337 if (verbosityLevel>=4){
1338 std::ios_base::fmtflags oldflags = std::cout.flags();
1339 std::cout << stagetag
1341 << std::setw(12) << std::setprecision(4) << std::scientific
1343 <<
" suggested dt: "
1344 << std::setw(12) << std::setprecision(4) << std::scientific
1347 std::cout.flags(oldflags);
1350 if (verbosityLevel>=2 && newdt!=dt)
1352 std::ios_base::fmtflags oldflags = std::cout.flags();
1353 std::cout <<
"changed dt to "
1354 << std::setw(12) << std::setprecision(4) << std::scientific
1357 std::cout.flags(oldflags);
1363 if (verbosityLevel>=4)
1364 std::cout << stagetag
1365 <<
"Combining residuals with selected dt..."
1367 alpha.axpy(dt,beta);
1368 if (verbosityLevel>=4)
1369 std::cout << stagetag
1370 <<
"Combining residuals with selected dt... done."
1374 if (verbosityLevel>=4)
1375 std::cout << stagetag <<
"Solving diagonal system..."
1377 ls.apply(D,*x[r],alpha,0.99);
1378 if (verbosityLevel>=4)
1379 std::cout << stagetag <<
"Solving diagonal system... done."
1383 limiter.poststage(*x[r]);
1386 if (verbosityLevel>=4)
1387 std::cout << stagetag <<
"Cleanup..." << std::endl;
1389 if (verbosityLevel>=4)
1390 std::cout << stagetag <<
"Cleanup... done" << std::endl;
1392 if (verbosityLevel>=4)
1393 std::cout << stagetag <<
"Finished." << std::endl;
1397 for(
unsigned i=1; i<method->
s(); ++i)
delete x[i];
1400 if (verbosityLevel>=4)
1401 std::cout << mytag <<
"Cleanup..." << std::endl;
1403 if (verbosityLevel>=4)
1404 std::cout << mytag <<
"Cleanup... done." << std::endl;
1413 class DefaultLimiter
1416 template<
typename V>
1420 template<
typename V>
1421 void poststage(V& v)
1425 const TimeSteppingParameterInterface<T> *method;
1431 TimeControllerInterface<T> *tc;
1441 sprintf(basename,
"%s",basename_);
1447 sprintf(basename,
"%s",basename_.c_str());
1452 sprintf(fname,
"%s-%05d",basename,i_);
1458 sprintf(fname,
"%s-%05d",basename,i);
CFLTimeController(R cfl_, R target_, const IGOS &igos_)
Definition: instationary/onestep.hh:716
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:371
Parameters to turn the OneStepMethod into an Alexander3 scheme.
Definition: instationary/onestep.hh:573
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:239
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:64
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:394
virtual std::string name() const =0
Return name of the scheme.
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: instationary/onestep.hh:842
Alexander2Parameter()
Definition: instationary/onestep.hh:419
OneStepThetaParameter(R theta_)
construct OneStepThetaParameter class
Definition: instationary/onestep.hh:114
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:150
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:72
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:208
OneStepMethodResult()
Definition: instationary/onestep.hh:765
OneStepMethodResult Result
Definition: instationary/onestep.hh:783
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:472
void increment()
Definition: instationary/onestep.hh:1462
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:224
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:284
Base parameter class for time stepping scheme parameters.
Definition: timesteppingparameterinterface.hh:19
OneStepMethodPartialResult()
Definition: instationary/onestep.hh:752
Controller interface for adaptive time stepping.
Definition: instationary/onestep.hh:670
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:531
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:142
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:165
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:554
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:49
Do one step of an explicit time-stepping scheme.
Definition: instationary/onestep.hh:1142
virtual RealType suggestTimestep(RealType time, RealType givendt)
Return name of the scheme.
Definition: instationary/onestep.hh:726
unsigned int timesteps
Definition: instationary/onestep.hh:746
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:201
Default time controller; just returns given dt.
Definition: instationary/onestep.hh:688
~ExplicitOneStepMethod()
Definition: instationary/onestep.hh:1190
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:621
void setStepNumber(int newstep)
change number of current step
Definition: instationary/onestep.hh:1205
Parameters to turn the OneStepMethod into an Alexander scheme.
Definition: instationary/onestep.hh:415
ExplicitEulerParameter()
Definition: instationary/onestep.hh:40
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: instationary/onestep.hh:1196
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:277
Shu3Parameter()
Definition: instationary/onestep.hh:261
R RealType
Definition: instationary/onestep.hh:711
HeunParameter()
Definition: instationary/onestep.hh:188
virtual unsigned s() const =0
Return number of stages of the method.
Insert standard boilerplate into log messages.
Definition: logtag.hh:181
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:379
double linear_solver_time
Definition: instationary/onestep.hh:748
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:516
Alexander3Parameter()
Definition: instationary/onestep.hh:577
T apply(T time, T dt, TrlV &xold, TrlV &xnew, Limiter &limiter)
Definition: instationary/onestep.hh:1237
void setStepNumber(int newstep)
change number of current step
Definition: instationary/onestep.hh:815
OneStepMethodPartialResult total
Definition: instationary/onestep.hh:763
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:523
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: instationary/onestep.hh:854
const PDESOLVER & getPDESolver() const
Access to the (non) linear solver.
Definition: instationary/onestep.hh:818
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:363
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:158
T apply(T time, T dt, TrlV &xold, TrlV &xnew)
do one step;
Definition: instationary/onestep.hh:1230
T apply(T time, T dt, TrlV &xold, F &f, TrlV &xnew)
do one step; This is a version which interpolates constraints at the start of each stage ...
Definition: instationary/onestep.hh:995
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:315
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:457
Parameters to turn the OneStepMethod into an one step theta method.
Definition: instationary/onestep.hh:107
FractionalStepParameter()
Definition: instationary/onestep.hh:495
void setMethod(const TimeSteppingParameterInterface< T > &method_)
redefine the method to be used; can be done before every step
Definition: instationary/onestep.hh:1215
Parameters to turn the ExplicitOneStepMethod into a Heun scheme.
Definition: instationary/onestep.hh:184
OneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, PDESOLVER &pdesolver_)
construct a new one step scheme
Definition: instationary/onestep.hh:797
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:387
R RealType
Definition: instationary/onestep.hh:673
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_)
construct a new one step scheme
Definition: instationary/onestep.hh:1160
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:637
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:645
Parameters to turn the ExplicitOneStepMethod into an explicite Euler method.
Definition: instationary/onestep.hh:36
PDESOLVER & getPDESolver()
Access to the (non) linear solver.
Definition: instationary/onestep.hh:824
ExplicitOneStepMethod(const TimeSteppingParameterInterface< T > &method_, IGOS &igos_, LS &ls_, TC &tc_)
construct a new one step scheme
Definition: instationary/onestep.hh:1182
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:80
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:300
int linear_solver_iterations
Definition: instationary/onestep.hh:749
virtual ~TimeControllerInterface()
every abstract base class has a virtual destructor
Definition: instationary/onestep.hh:680
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:652
virtual bool implicit() const =0
Return true if method is implicit.
OneStepMethodPartialResult successful
Definition: instationary/onestep.hh:764
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:124
virtual RealType suggestTimestep(RealType time, RealType givendt)=0
Return name of the scheme.
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:465
const char * getName()
Definition: instationary/onestep.hh:1456
virtual R d(int r) const =0
Return entries of the d Vector.
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:441
Parameters to turn the ExplicitOneStepMethod into a third order strong stability preserving (SSP) sch...
Definition: instationary/onestep.hh:257
limit time step to maximum dt * CFL number
Definition: instationary/onestep.hh:708
R RealType
Definition: instationary/onestep.hh:691
Parameters to turn the OneStepMethod into a fractional step theta scheme.
Definition: instationary/onestep.hh:491
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:216
int nonlinear_solver_iterations
Definition: instationary/onestep.hh:750
const F & f
Definition: common/constraints.hh:145
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:292
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:547
CFLTimeController(R cfl_, const IGOS &igos_)
Definition: instationary/onestep.hh:713
Definition: instationary/onestep.hh:761
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:308
Do one step of a time-stepping scheme.
Definition: instationary/onestep.hh:778
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:614
Parameters to turn the ExplicitOneStepMethod into a classical fourth order Runge-Kutta method...
Definition: instationary/onestep.hh:334
void setVerbosityLevel(int level)
change verbosity level; 0 means completely quiet
Definition: instationary/onestep.hh:806
void setTarget(R target_)
Definition: instationary/onestep.hh:719
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:56
virtual R d(int i) const
Return entries of the d Vector.
Definition: instationary/onestep.hh:232
Base class for all PDELab exceptions.
Definition: exceptions.hh:17
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:356
FilenameHelper(const char *basename_, int i_=0)
Definition: instationary/onestep.hh:1438
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:449
virtual R a(int r, int i) const
Return entries of the A matrix.
Definition: instationary/onestep.hh:629
virtual bool implicit() const
Return true if method is implicit.
Definition: instationary/onestep.hh:434
virtual RealType suggestTimestep(RealType time, RealType givendt)
Return name of the scheme.
Definition: instationary/onestep.hh:695
const char * getName(int i_)
Definition: instationary/onestep.hh:1450
RK4Parameter()
Definition: instationary/onestep.hh:338
virtual R b(int r, int i) const
Return entries of the B matrix.
Definition: instationary/onestep.hh:539
Definition: instationary/onestep.hh:1435
double assembler_time
Definition: instationary/onestep.hh:747
FilenameHelper(const std::string &basename_, int i_=0)
Definition: instationary/onestep.hh:1444
Definition: instationary/onestep.hh:744
const Result & result() const
Definition: instationary/onestep.hh:829
virtual unsigned s() const
Return number of stages s of the method.
Definition: instationary/onestep.hh:134
virtual std::string name() const
Return name of the scheme.
Definition: instationary/onestep.hh:87