diff --git a/src/Makefile.am b/src/Makefile.am index 23d83a43ffc5b26aadf274c63e0df6099f9d9802..0f184e850fc1e2eeaa7a86e73db317216da8b436 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -3,6 +3,7 @@ bin_PROGRAMS = sand-wedge-2D sand-wedge-3D common_sources = \ assemblers.cc \ boundary_writer.cc \ + coupledtimestepper.cc \ enumparser.cc \ fixedpointiterator.cc \ friction_writer.cc \ diff --git a/src/coupledtimestepper.cc b/src/coupledtimestepper.cc new file mode 100644 index 0000000000000000000000000000000000000000..17432aa14e5410baefa744b15b5f4ced29c35873 --- /dev/null +++ b/src/coupledtimestepper.cc @@ -0,0 +1,53 @@ +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + +#include <dune/solvers/norms/energynorm.hh> + +#include "coupledtimestepper.hh" +#include "fixedpointiterator.hh" + +template <class Factory, class StateUpdater, class VelocityUpdater> +CoupledTimeStepper<Factory, StateUpdater, VelocityUpdater>::CoupledTimeStepper( + double finalTime, Factory &factory, Dune::ParameterTree const &parset, + std::shared_ptr<Nonlinearity> globalFriction, + std::shared_ptr<StateUpdater> stateUpdater, + std::shared_ptr<VelocityUpdater> velocityUpdater, + std::function<void(double, Vector &)> externalForces) + : finalTime_(finalTime), + factory_(factory), + parset_(parset), + globalFriction_(globalFriction), + stateUpdater_(stateUpdater), + velocityUpdater_(velocityUpdater), + externalForces_(externalForces) {} + +template <class Factory, class StateUpdater, class VelocityUpdater> +int CoupledTimeStepper<Factory, StateUpdater, VelocityUpdater>::step( + double relativeTime, double relativeTau) { + stateUpdater_->nextTimeStep(); + velocityUpdater_->nextTimeStep(); + + auto const newRelativeTime = relativeTime + relativeTau; + Vector ell; + externalForces_(newRelativeTime, ell); + + Matrix velocityMatrix; + Vector velocityRHS; + Vector velocityIterate; + + auto const tau = relativeTau * finalTime_; + stateUpdater_->setup(tau); + velocityUpdater_->setup(ell, tau, newRelativeTime, velocityRHS, + velocityIterate, velocityMatrix); + EnergyNorm<Matrix, Vector> const velocityMatrixNorm(velocityMatrix); + + FixedPointIterator<Factory, StateUpdater, VelocityUpdater> fixedPointIterator( + factory_, parset_, globalFriction_); + auto const iterations = + fixedPointIterator.run(stateUpdater_, velocityUpdater_, velocityMatrix, + velocityMatrixNorm, velocityRHS, velocityIterate); + return iterations; +} + +#include "coupledtimestepper_tmpl.cc" diff --git a/src/coupledtimestepper.hh b/src/coupledtimestepper.hh new file mode 100644 index 0000000000000000000000000000000000000000..52ece50ce098f81648901eef69fb299f436131c3 --- /dev/null +++ b/src/coupledtimestepper.hh @@ -0,0 +1,35 @@ +#ifndef SRC_COUPLEDTIMESTEPPER_HH +#define SRC_COUPLEDTIMESTEPPER_HH + +#include <functional> +#include <memory> + +#include <dune/common/parametertree.hh> + +template <class Factory, class StateUpdater, class VelocityUpdater> +class CoupledTimeStepper { + using Vector = typename Factory::Vector; + using Matrix = typename Factory::Matrix; + using ConvexProblem = typename Factory::ConvexProblem; + using Nonlinearity = typename ConvexProblem::NonlinearityType; + +public: + CoupledTimeStepper(double finalTime, Factory &factory, + Dune::ParameterTree const &parset, + std::shared_ptr<Nonlinearity> globalFriction, + std::shared_ptr<StateUpdater> stateUpdater, + std::shared_ptr<VelocityUpdater> velocityUpdater, + std::function<void(double, Vector &)> externalForces); + + int step(double relativeTime, double relativeTau); + +private: + double finalTime_; + Factory &factory_; + Dune::ParameterTree const &parset_; + std::shared_ptr<Nonlinearity> globalFriction_; + std::shared_ptr<StateUpdater> stateUpdater_; + std::shared_ptr<VelocityUpdater> velocityUpdater_; + std::function<void(double, Vector &)> externalForces_; +}; +#endif diff --git a/src/coupledtimestepper_tmpl.cc b/src/coupledtimestepper_tmpl.cc new file mode 100644 index 0000000000000000000000000000000000000000..d72f58449333d68c48b23b80265a0985f05351e4 --- /dev/null +++ b/src/coupledtimestepper_tmpl.cc @@ -0,0 +1,26 @@ +#ifndef DIM +#error DIM unset +#endif + +#include <dune/common/function.hh> + +#include <dune/tnnmg/problem-classes/convexproblem.hh> + +#include <dune/tectonic/globalfriction.hh> +#include <dune/tectonic/myblockproblem.hh> + +#include "explicitgrid.hh" +#include "explicitvectors.hh" + +#include "solverfactory.hh" +#include "state/stateupdater.hh" +#include "timestepping.hh" + +using Function = Dune::VirtualFunction<double, double>; +using Factory = SolverFactory< + DIM, MyBlockProblem<ConvexProblem<GlobalFriction<Matrix, Vector>, Matrix>>, + Grid>; +using MyStateUpdater = StateUpdater<ScalarVector, Vector>; +using MyVelocityUpdater = TimeSteppingScheme<Vector, Matrix, Function, DIM>; + +template class CoupledTimeStepper<Factory, MyStateUpdater, MyVelocityUpdater>; diff --git a/src/sand-wedge.cc b/src/sand-wedge.cc index 1f35a207b671688a9880ab01d706e80ba868bb39..520cfc4c3d19e7012fba34f095759ce3b07fc2b3 100644 --- a/src/sand-wedge.cc +++ b/src/sand-wedge.cc @@ -68,6 +68,7 @@ #include "assemblers.hh" #include "tobool.hh" +#include "coupledtimestepper.hh" #include "enumparser.hh" #include "enums.hh" #include "fixedpointiterator.hh" @@ -93,64 +94,6 @@ void initPython() { Python::run("sys.path.append('" datadir "')"); } -template <class Factory, class StateUpdater, class VelocityUpdater> -class CoupledTimeStepper { - using Vector = typename Factory::Vector; - using Matrix = typename Factory::Matrix; - using ConvexProblem = typename Factory::ConvexProblem; - using Nonlinearity = typename ConvexProblem::NonlinearityType; - -public: - CoupledTimeStepper(double finalTime, Factory &factory, - Dune::ParameterTree const &parset, - std::shared_ptr<Nonlinearity> globalFriction, - std::shared_ptr<StateUpdater> stateUpdater, - std::shared_ptr<VelocityUpdater> velocityUpdater, - std::function<void(double, Vector &)> externalForces) - : finalTime_(finalTime), - factory_(factory), - parset_(parset), - globalFriction_(globalFriction), - stateUpdater_(stateUpdater), - velocityUpdater_(velocityUpdater), - externalForces_(externalForces) {} - - int step(double relativeTime, double relativeTau) { - stateUpdater_->nextTimeStep(); - velocityUpdater_->nextTimeStep(); - - auto const newRelativeTime = relativeTime + relativeTau; - Vector ell; - externalForces_(newRelativeTime, ell); - - Matrix velocityMatrix; - Vector velocityRHS; - Vector velocityIterate; - - auto const tau = relativeTau * finalTime_; - stateUpdater_->setup(tau); - velocityUpdater_->setup(ell, tau, newRelativeTime, velocityRHS, - velocityIterate, velocityMatrix); - EnergyNorm<Matrix, Vector> const velocityMatrixNorm(velocityMatrix); - - FixedPointIterator<Factory, StateUpdater, VelocityUpdater> - fixedPointIterator(factory_, parset_, globalFriction_); - auto const iterations = fixedPointIterator.run( - stateUpdater_, velocityUpdater_, velocityMatrix, velocityMatrixNorm, - velocityRHS, velocityIterate); - return iterations; - } - -private: - double finalTime_; - Factory &factory_; - Dune::ParameterTree const &parset_; - std::shared_ptr<Nonlinearity> globalFriction_; - std::shared_ptr<StateUpdater> stateUpdater_; - std::shared_ptr<VelocityUpdater> velocityUpdater_; - std::function<void(double, Vector &)> externalForces_; -}; - int main(int argc, char *argv[]) { try { Dune::ParameterTree parset;