diff --git a/examples/convection_diffusion.cc b/examples/convection_diffusion.cc index 78a43eb1390ab81f8403082e2f8239c645766b51..d5e229c380ac77225b2ea34e66ed9c3422eb3ac0 100644 --- a/examples/convection_diffusion.cc +++ b/examples/convection_diffusion.cc @@ -35,7 +35,7 @@ int main(int argc, char** argv) prob.addDirichletBC(predicate, 0, 0, dbcValues); AdaptInfo adaptInfo("adapt"); - prob.buildAfterAdapt(adaptInfo, Flag(0)); + prob.assemble(adaptInfo); prob.solve(adaptInfo); prob.writeFiles(adaptInfo, true); diff --git a/examples/ellipt.cc b/examples/ellipt.cc index 86b53b9e4a79b041c05ac04d0425febb04f30d7f..d1a07622ca7eee8f9eee49e9e126c3ed3dab8c67 100644 --- a/examples/ellipt.cc +++ b/examples/ellipt.cc @@ -38,14 +38,14 @@ int main(int argc, char** argv) prob.initialize(INIT_ALL); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); - prob.addMatrixOperator(opL, _0, _0); + prob.addMatrixOperator(opL, 0, 0); auto opForce = makeOperator(tag::test{}, f, 6); - prob.addVectorOperator(opForce, _0); + prob.addVectorOperator(opForce, 0); // set boundary condition auto boundary = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8 || x[0] > 1.0-1.e-8 || x[1] > 1.0-1.e-8; }; - prob.addDirichletBC(boundary, _0, _0, g); + prob.addDirichletBC(boundary, 0, 0, g); AdaptInfo adaptInfo("adapt"); @@ -62,17 +62,17 @@ int main(int argc, char** argv) widths.push_back(h); prob.globalBasis().update(gridView); - prob.buildAfterAdapt(adaptInfo, Flag(0)); + prob.assemble(adaptInfo); prob.solve(adaptInfo); - double errorL2 = integrate(sqr(g - prob.getSolution(_0)), gridView, 6); + double errorL2 = integrate(sqr(g - prob.solution(0)), gridView, 6); errL2.push_back(std::sqrt(errorL2)); - double errorH1 = errorL2 + integrate(unary_dot(grad_g - gradientAtQP(prob.getSolution(_0))), gridView, 6); + double errorH1 = errorL2 + integrate(unary_dot(grad_g - gradientAtQP(prob.solution(0))), gridView, 6); errH1.push_back(std::sqrt(errorH1)); #if WRITE_FILES Dune::VTKWriter<typename ElliptProblem::GridView> vtkWriter(gridView); - vtkWriter.addVertexData(prob.getSolution(_0), Dune::VTK::FieldInfo("u", Dune::VTK::FieldInfo::Type::scalar, 1)); + vtkWriter.addVertexData(prob.solution(0), Dune::VTK::FieldInfo("u", Dune::VTK::FieldInfo::Type::scalar, 1)); vtkWriter.write("u_" + std::to_string(i)); #endif } diff --git a/examples/heat.cc b/examples/heat.cc index 8bab7a8fdc841e90bb537a9d249c8b1e875ced5d..3f75d4a4aef060e93362662e66bc00d6040c3445 100644 --- a/examples/heat.cc +++ b/examples/heat.cc @@ -31,16 +31,16 @@ int main(int argc, char** argv) AdaptInfo adaptInfo("adapt"); - auto* invTau = probInstat.getInvTau(); + auto invTau = std::ref(probInstat.invTau()); - auto opTimeLhs = makeOperator(tag::test_trial{}, std::ref(*invTau)); + auto opTimeLhs = makeOperator(tag::test_trial{}, invTau); prob.addMatrixOperator(opTimeLhs, 0, 0); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); prob.addMatrixOperator(opL, 0, 0); auto opTimeRhs = makeOperator(tag::test{}, - invokeAtQP([invTau](double u) { return u * (*invTau); }, prob.getSolution(0)), 2); + invokeAtQP([invTau](double u) { return u * invTau.get(); }, prob.solution(0)), 2); prob.addVectorOperator(opTimeRhs, 0); auto opForce = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0); diff --git a/examples/navier_stokes.cc b/examples/navier_stokes.cc index 762b369deba45fa66ffe002b519cb6bebf03a7ca..c82e313dcae1a52e3b6bf2502570117695ace2be 100644 --- a/examples/navier_stokes.cc +++ b/examples/navier_stokes.cc @@ -48,12 +48,16 @@ int main(int argc, char** argv) auto _v = Dune::Indices::_0; auto _p = Dune::Indices::_1; + auto invTau = std::ref(probInstat.invTau()); + // <1/tau * u, v> - auto opTime = makeOperator(tag::testvec_trialvec{}, density); + auto opTime = makeOperator(tag::testvec_trialvec{}, + density * invTau); prob.addMatrixOperator(opTime, _v, _v); // <1/tau * u^old, v> - auto opTimeOld = makeOperator(tag::testvec{}, density * prob.getSolution(_v)); + auto opTimeOld = makeOperator(tag::testvec{}, + density * invTau * prob.solution(_v)); prob.addVectorOperator(opTimeOld, _v); @@ -62,17 +66,20 @@ int main(int argc, char** argv) prob.addMatrixOperator(opStokes, treepath(), treepath()); // <(u * nabla)u_i^old, v_i> - auto opNonlin1 = makeOperator(tag::testvec_trialvec{}, density * trans(gradientAtQP(prob.getSolution(_v)))); + auto opNonlin1 = makeOperator(tag::testvec_trialvec{}, + density * trans(gradientAtQP(prob.solution(_v)))); prob.addMatrixOperator(opNonlin1, _v, _v); for (std::size_t i = 0; i < AMDIS_DOW; ++i) { // <(u^old * nabla)u_i, v_i> - auto opNonlin2 = makeOperator(tag::test_gradtrial{}, density * prob.getSolution(_v)); + auto opNonlin2 = makeOperator(tag::test_gradtrial{}, + density * prob.solution(_v)); prob.addMatrixOperator(opNonlin2, treepath(_v,i), treepath(_v,i)); } // <(u^old * grad(u_i^old)), v_i> - auto opNonlin3 = makeOperator(tag::testvec{}, trans(gradientAtQP(prob.getSolution(_v))) * prob.getSolution(_v)); + auto opNonlin3 = makeOperator(tag::testvec{}, + trans(gradientAtQP(prob.solution(_v))) * prob.solution(_v)); prob.addVectorOperator(opNonlin3, _v); // define boundary regions @@ -98,8 +105,8 @@ int main(int argc, char** argv) prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0); // set initial conditions - prob.getSolution(_v).interpolate(parabolic_y); - prob.getSolution(_p).interpolate(0.0); + prob.solution(_v).interpolate(parabolic_y); + prob.solution(_p).interpolate(0.0); AdaptInstationary adapt("adapt", prob, adaptInfo, probInstat, adaptInfo); adapt.adapt(); diff --git a/examples/stokes0.cc b/examples/stokes0.cc index 714dae9c7ad403867b5cffddc4ffd0f7d8e36d47..11fe7d2c8cd59208926d1b7fbaa64201d42de9e6 100644 --- a/examples/stokes0.cc +++ b/examples/stokes0.cc @@ -25,73 +25,73 @@ using StokesProblem = ProblemStat<StokesParam>; int main(int argc, char** argv) { - AMDiS::init(argc, argv); + AMDiS::init(argc, argv); - StokesProblem prob("stokes"); - prob.initialize(INIT_ALL); + StokesProblem prob("stokes"); + prob.initialize(INIT_ALL); - double viscosity = 1.0; - Parameters::get("stokes->viscosity", viscosity); + double viscosity = 1.0; + Parameters::get("stokes->viscosity", viscosity); - // tree-paths for components - auto _v = Dune::Indices::_0; - auto _p = Dune::Indices::_1; + // tree-paths for components + auto _v = Dune::Indices::_0; + auto _p = Dune::Indices::_1; - // <viscosity*grad(u_i), grad(v_i)> - for (std::size_t i = 0; i < DOW; ++i) { - auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity); - prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i)); + // <viscosity*grad(u_i), grad(v_i)> + for (std::size_t i = 0; i < DOW; ++i) { + auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity); + prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i)); - // <d_i(v_i), p> - auto opP = makeOperator(tag::partialtest_trial{i}, 1.0); - prob.addMatrixOperator(opP, treepath(_v,i), _p); + // <d_i(v_i), p> + auto opP = makeOperator(tag::partialtest_trial{i}, 1.0); + prob.addMatrixOperator(opP, treepath(_v,i), _p); - // <q, d_i(u_i)> - auto opDiv = makeOperator(tag::test_partialtrial{i}, 1.0); - prob.addMatrixOperator(opDiv, _p, treepath(_v,i)); - } + // <q, d_i(u_i)> + auto opDiv = makeOperator(tag::test_partialtrial{i}, 1.0); + prob.addMatrixOperator(opDiv, _p, treepath(_v,i)); + } - auto opZero = makeOperator(tag::test_trial{}, 0.0); - prob.addMatrixOperator(opZero, _p, _p); + auto opZero = makeOperator(tag::test_trial{}, 0.0); + prob.addMatrixOperator(opZero, _p, _p); - // define boundary regions - auto left = [](auto const& x) { return x[0] < 1.e-8; }; - auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; }; + // define boundary regions + auto left = [](auto const& x) { return x[0] < 1.e-8; }; + auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; }; - // define boundary values - auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW> - { - return {0.0, x[1]*(1.0 - x[1])}; - }; + // define boundary values + auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW> + { + return {0.0, x[1]*(1.0 - x[1])}; + }; - auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW> - { - return {0.0, 0.0}; - }; + auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW> + { + return {0.0, 0.0}; + }; - // set boundary conditions for velocity - prob.addDirichletBC(left, _v, _v, parabolic_y); - prob.addDirichletBC(not_left, _v, _v, zero); + // set boundary conditions for velocity + prob.addDirichletBC(left, _v, _v, parabolic_y); + prob.addDirichletBC(not_left, _v, _v, zero); - prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0); + prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0); - AdaptInfo adaptInfo("adapt"); + AdaptInfo adaptInfo("adapt"); - // assemble and solve system - prob.buildAfterAdapt(adaptInfo, Flag(0)); + // assemble and solve system + prob.assemble(adaptInfo); #ifdef DEBUG_MTL - // write matrix to file - mtl::io::matrix_market_ostream out("matrix_stokes0.mtx"); - out << prob.getSystemMatrix().matrix(); - std::cout << prob.getSystemMatrix().matrix() << '\n'; + // write matrix to file + mtl::io::matrix_market_ostream out("matrix_stokes0.mtx"); + out << prob.systemMatrix().matrix(); + std::cout << prob.systemMatrix().matrix() << '\n'; #endif - prob.solve(adaptInfo); + prob.solve(adaptInfo); - // output solution - prob.writeFiles(adaptInfo); + // output solution + prob.writeFiles(adaptInfo); - AMDiS::finalize(); - return 0; + AMDiS::finalize(); + return 0; } diff --git a/examples/stokes1.cc b/examples/stokes1.cc index a96124061cf265e1f59d89946f49d5d1ba6b190b..b07ef3614f0c364c2f5ccb41c579e7a7ae71b40d 100644 --- a/examples/stokes1.cc +++ b/examples/stokes1.cc @@ -25,74 +25,74 @@ using StokesProblem = ProblemStat<StokesParam>; int main(int argc, char** argv) { - AMDiS::init(argc, argv); + AMDiS::init(argc, argv); - StokesProblem prob("stokes"); - prob.initialize(INIT_ALL); + StokesProblem prob("stokes"); + prob.initialize(INIT_ALL); - double viscosity = 1.0; - Parameters::get("stokes->viscosity", viscosity); + double viscosity = 1.0; + Parameters::get("stokes->viscosity", viscosity); - // tree-paths for components - auto _v = Dune::Indices::_0; - auto _p = Dune::Indices::_1; + // tree-paths for components + auto _v = Dune::Indices::_0; + auto _p = Dune::Indices::_1; - // <viscosity*grad(u_i), grad(v_i)> - for (std::size_t i = 0; i < DOW; ++i) { - auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity); - prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i)); - } + // <viscosity*grad(u_i), grad(v_i)> + for (std::size_t i = 0; i < DOW; ++i) { + auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity); + prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i)); + } - // <d_i(v_i), p> - auto opP = makeOperator(tag::divtestvec_trial{}, 1.0); - prob.addMatrixOperator(opP, _v, _p); + // <d_i(v_i), p> + auto opP = makeOperator(tag::divtestvec_trial{}, 1.0); + prob.addMatrixOperator(opP, _v, _p); - // <q, d_i(u_i)> - auto opDiv = makeOperator(tag::test_divtrialvec{}, 1.0); - prob.addMatrixOperator(opDiv, _p, _v); + // <q, d_i(u_i)> + auto opDiv = makeOperator(tag::test_divtrialvec{}, 1.0); + prob.addMatrixOperator(opDiv, _p, _v); - auto opZero = makeOperator(tag::test_trial{}, 0.0); - prob.addMatrixOperator(opZero, _p, _p); + auto opZero = makeOperator(tag::test_trial{}, 0.0); + prob.addMatrixOperator(opZero, _p, _p); - // define boundary regions - auto left = [](auto const& x) { return x[0] < 1.e-8; }; - auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; }; + // define boundary regions + auto left = [](auto const& x) { return x[0] < 1.e-8; }; + auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; }; - // define boundary values - auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW> - { - return {0.0, x[1]*(1.0 - x[1])}; - }; + // define boundary values + auto parabolic_y = [](auto const& x) -> Dune::FieldVector<double,DOW> + { + return {0.0, x[1]*(1.0 - x[1])}; + }; - auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW> - { - return {0.0, 0.0}; - }; + auto zero = [](auto const& x) -> Dune::FieldVector<double,DOW> + { + return {0.0, 0.0}; + }; - // set boundary conditions for velocity - prob.addDirichletBC(left, _v, _v, parabolic_y); - prob.addDirichletBC(not_left, _v, _v, zero); + // set boundary conditions for velocity + prob.addDirichletBC(left, _v, _v, parabolic_y); + prob.addDirichletBC(not_left, _v, _v, zero); - prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0); + prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0); - AdaptInfo adaptInfo("adapt"); + AdaptInfo adaptInfo("adapt"); - // assemble and solve system - prob.buildAfterAdapt(adaptInfo, Flag(0)); + // assemble and solve system + prob.assemble(adaptInfo); #ifdef DEBUG_MTL - // write matrix to file - mtl::io::matrix_market_ostream out("matrix_stokes1.mtx"); - out << prob.getSystemMatrix().matrix(); - std::cout << prob.getSystemMatrix().matrix() << '\n'; + // write matrix to file + mtl::io::matrix_market_ostream out("matrix_stokes1.mtx"); + out << prob.systemMatrix().matrix(); + std::cout << prob.systemMatrix().matrix() << '\n'; #endif - prob.solve(adaptInfo); + prob.solve(adaptInfo); - // output solution - prob.writeFiles(adaptInfo); + // output solution + prob.writeFiles(adaptInfo); - AMDiS::finalize(); - return 0; + AMDiS::finalize(); + return 0; } diff --git a/examples/stokes3.cc b/examples/stokes3.cc index 21d5b60543f84f7735ab90288648b77d50e75a8e..ce338c56c2cc98ee5d39a8755ba7c5ab593a4571 100644 --- a/examples/stokes3.cc +++ b/examples/stokes3.cc @@ -62,8 +62,7 @@ int main(int argc, char** argv) AdaptInfo adaptInfo("adapt"); // assemble and solve system - prob.buildAfterAdapt(adaptInfo, Flag(0)); - + prob.assemble(adaptInfo); prob.solve(adaptInfo); // output solution diff --git a/examples/vecellipt.cc b/examples/vecellipt.cc index ac53b093a422799a7b23b6fe3233c7699bdf4c3f..5528712cd3e82c4c8a968912d75d2e334ad8692a 100644 --- a/examples/vecellipt.cc +++ b/examples/vecellipt.cc @@ -46,9 +46,9 @@ int main(int argc, char** argv) // write matrix to file if (Parameters::get<int>("elliptMesh->global refinements").value_or(0) < 4) { mtl::io::matrix_market_ostream out("matrix.mtx"); - out << prob.getSystemMatrix().matrix(); + out << prob.systemMatrix().matrix(); - std::cout << prob.getSystemMatrix().matrix() << '\n'; + std::cout << prob.systemMatrix().matrix() << '\n'; } #endif diff --git a/src/amdis/AdaptBase.hpp b/src/amdis/AdaptBase.hpp index 004c893553a0f56f1ca255668f3f4a5159dc123e..fff064bc86efc0841f2357760b5e5aebd77eac91 100644 --- a/src/amdis/AdaptBase.hpp +++ b/src/amdis/AdaptBase.hpp @@ -37,13 +37,13 @@ namespace AMDiS virtual int adapt() = 0; /// Returns \ref name - std::string const& getName() const + std::string const& name() const { return name_; } /// Returns \ref problemIteration - ProblemIterationInterface* getProblemIteration() const + ProblemIterationInterface* problemIteration() const { return problemIteration_; } @@ -55,13 +55,13 @@ namespace AMDiS } /// Returns \ref adaptInfo - AdaptInfo& getAdaptInfo() const + AdaptInfo& adaptInfo() const { return adaptInfo_; } /// Returns \ref problemTime - ProblemTimeInterface* getProblemTime() const + ProblemTimeInterface* problemTime() const { return problemTime_; } @@ -73,7 +73,7 @@ namespace AMDiS } /// Returns \ref initialAdaptInfo - AdaptInfo& getInitialAdaptInfo() const + AdaptInfo& initialAdaptInfo() const { return *initialAdaptInfo_; } diff --git a/src/amdis/AdaptInfo.cpp b/src/amdis/AdaptInfo.cpp index 3db8e71400ef1f9878f2f6197fd9e1350e6754dc..e71e6a6687e44642146ee2eaa6aced8b482d504a 100644 --- a/src/amdis/AdaptInfo.cpp +++ b/src/amdis/AdaptInfo.cpp @@ -29,28 +29,28 @@ namespace AMDiS : name_(name) { // init(); - Parameters::get(name + "->start time", startTime); - time = startTime; + Parameters::get(name + "->start time", startTime_); + time_ = startTime_; - Parameters::get(name + "->timestep", timestep); - Parameters::get(name + "->end time", endTime); - Parameters::get(name + "->max iteration", maxSpaceIteration); - Parameters::get(name + "->max timestep iteration", maxTimestepIteration); - Parameters::get(name + "->max time iteration", maxTimeIteration); - Parameters::get(name + "->min timestep", minTimestep); - Parameters::get(name + "->max timestep", maxTimestep); - Parameters::get(name + "->number of timesteps", nTimesteps); - Parameters::get(name + "->time tolerance", globalTimeTolerance); + Parameters::get(name + "->timestep", timestep_); + Parameters::get(name + "->end time", endTime_); + Parameters::get(name + "->max iteration", maxSpaceIteration_); + Parameters::get(name + "->max timestep iteration", maxTimestepIteration_); + Parameters::get(name + "->max time iteration", maxTimeIteration_); + Parameters::get(name + "->min timestep", minTimestep_); + Parameters::get(name + "->max timestep", maxTimestep_); + Parameters::get(name + "->number of timesteps", nTimesteps_); + Parameters::get(name + "->time tolerance", globalTimeTolerance_); } void AdaptInfo::printTimeErrorLowInfo() const { - for (auto const& scalContent : scalContents) + for (auto const& scalContent : scalContents_) { auto i = scalContent.first; std::cout << " Time error estimate ["<<i<<"] = " - << getTimeEstCombined(i) << "\n" + << timeEstCombined(i) << "\n" << " Time error estimate sum ["<<i<<"] = " << scalContent.second.est_t_sum << "\n" << " Time error estimate max ["<<i<<"] = " @@ -64,17 +64,17 @@ namespace AMDiS void AdaptInfo::reset() { - spaceIteration = -1; - timestepIteration = 0; - timeIteration = 0; - time = 0.0; - timestep = 0.0; - timestepNumber = 0; - solverIterations = 0; - solverResidual = 0.0; + spaceIteration_ = -1; + timestepIteration_ = 0; + timeIteration_ = 0; + time_ = 0.0; + timestep_ = 0.0; + timestepNumber_ = 0; + solverIterations_ = 0; + solverResidual_ = 0.0; - Parameters::get(name_ + "->timestep", timestep); - lastProcessedTimestep=timestep; + Parameters::get(name_ + "->timestep", timestep_); + lastProcessedTimestep_ = timestep_; } } // end namespace AMDiS diff --git a/src/amdis/AdaptInfo.hpp b/src/amdis/AdaptInfo.hpp index 7c2444f4447fbfee3b6f237883d0119b01c9fa18..72e19008774e0712346faf379e65dd9c5b7ce09d 100644 --- a/src/amdis/AdaptInfo.hpp +++ b/src/amdis/AdaptInfo.hpp @@ -86,7 +86,7 @@ namespace AMDiS /// Returns whether space tolerance is reached. virtual bool spaceToleranceReached() const { - for (auto const& scalContent : scalContents) { + for (auto const& scalContent : scalContents_) { if (!(scalContent.second.est_sum < scalContent.second.spaceTolerance)) return false; } @@ -97,7 +97,7 @@ namespace AMDiS /// Returns whether space tolerance of component associated with key is reached. virtual bool spaceToleranceReached(Key key) const { - if (!(getScalContent(key).est_sum < getScalContent(key).spaceTolerance)) + if (!(scalContent(key).est_sum < scalContent(key).spaceTolerance)) return false; else return true; @@ -112,8 +112,8 @@ namespace AMDiS /// Returns whether time tolerance is reached. virtual bool timeToleranceReached() const { - for (auto const& scalContent : scalContents) - if (!(getTimeEstCombined(scalContent.first) < scalContent.second.timeTolerance)) + for (auto const& scalContent : scalContents_) + if (!(timeEstCombined(scalContent.first) < scalContent.second.timeTolerance)) return false; return true; @@ -122,7 +122,7 @@ namespace AMDiS /// Returns whether time tolerance of component associated with key is reached. virtual bool timeToleranceReached(Key key) const { - if (!(getTimeEstCombined(key) < getScalContent(key).timeTolerance)) + if (!(timeEstCombined(key) < scalContent(key).timeTolerance)) return false; else return true; @@ -137,8 +137,8 @@ namespace AMDiS /// Returns whether time error is under its lower bound. virtual bool timeErrorLow() const { - for (auto const& scalContent : scalContents) - if (!(getTimeEstCombined(scalContent.first) < scalContent.second.timeErrLow)) + for (auto const& scalContent : scalContents_) + if (!(timeEstCombined(scalContent.first) < scalContent.second.timeErrLow)) return false; return true; @@ -146,146 +146,146 @@ namespace AMDiS /// Returns the time estimation as a combination /// of maximal and integral time error - double getTimeEstCombined(Key key) const + double timeEstCombined(Key key) const { return - getScalContent(key).est_t_max * getScalContent(key).fac_max + - getScalContent(key).est_t_sum * getScalContent(key).fac_sum; + scalContent(key).est_t_max * scalContent(key).fac_max + + scalContent(key).est_t_sum * scalContent(key).fac_sum; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeEstCombined(const TP& tp) const + double timeEstCombined(const TP& tp) const { - return getTimeEstCombined(to_string(tp)); + return tTimeEstCombined(to_string(tp)); } /// Print debug information about time error and its bound. void printTimeErrorLowInfo() const; - /// Returns \ref spaceIteration. - int getSpaceIteration() const + /// Returns \ref spaceIteration_. + int spaceIteration() const { - return spaceIteration; + return spaceIteration_; } - /// Sets \ref spaceIteration. + /// Sets \ref spaceIteration_. void setSpaceIteration(int it) { - spaceIteration = it; + spaceIteration_ = it; } - /// Returns \ref maxSpaceIteration. - int getMaxSpaceIteration() const + /// Returns \ref maxSpaceIteration_. + int maxSpaceIteration() const { - return maxSpaceIteration; + return maxSpaceIteration_; } - /// Sets \ref maxSpaceIteration. - void setMaxSpaceIteration(int it) + /// Sets \ref maxSpaceIteration_. + void maxSpaceIteration(int it) { - maxSpaceIteration = it; + maxSpaceIteration_ = it; } - /// Increments \ref spaceIteration by 1; + /// Increments \ref spaceIteration_ by 1; void incSpaceIteration() { - spaceIteration++; + spaceIteration_++; } - /// Sets \ref timestepIteration. + /// Sets \ref timestepIteration_. void setTimestepIteration(int it) { - timestepIteration = it; + timestepIteration_ = it; } - /// Returns \ref timestepIteration. - int getTimestepIteration() const + /// Returns \ref timestepIteration_. + int timestepIteration() const { - return timestepIteration; + return timestepIteration_; } - /// Increments \ref timestepIteration by 1; + /// Increments \ref timestepIteration_ by 1; void incTimestepIteration() { - timestepIteration++; + timestepIteration_++; } - /// Returns \ref maxTimestepIteration. - int getMaxTimestepIteration() const + /// Returns \ref maxTimestepIteration_. + int maxTimestepIteration() const { - return maxTimestepIteration; + return maxTimestepIteration_; } /// Sets \ref maxTimestepIteration. void setMaxTimestepIteration(int it) { - maxTimestepIteration = it; + maxTimestepIteration_ = it; } - /// Sets \ref timeIteration. + /// Sets \ref timeIteration_. void setTimeIteration(int it) { - timeIteration = it; + timeIteration_ = it; } - /// Returns \ref timeIteration. - int getTimeIteration() const + /// Returns \ref timeIteration_. + int timeIteration() const { - return timeIteration; + return timeIteration_; } - /// Increments \ref timesIteration by 1; + /// Increments \ref timesIteration_ by 1; void incTimeIteration() { - timeIteration++; + timeIteration_++; } - /// Returns \ref maxTimeIteration. - int getMaxTimeIteration() const + /// Returns \ref maxTimeIteration_. + int maxTimeIteration() const { - return maxTimeIteration; + return maxTimeIteration_; } - /// Sets \ref maxTimeIteration. + /// Sets \ref maxTimeIteration_. void setMaxTimeIteration(int it) { - maxTimeIteration = it; + maxTimeIteration_ = it; } - /// Returns \ref timestepNumber. - int getTimestepNumber() const + /// Returns \ref timestepNumber_. + int timestepNumber() const { - return timestepNumber; + return timestepNumber_; } /// Sets \ref timestepNumber. void setTimestepNumber(int num) { - timestepNumber = std::min(nTimesteps, num); + timestepNumber_ = std::min(nTimesteps_, num); } - /// Returns \ref nTimesteps. - int getNumberOfTimesteps() const + /// Returns \ref nTimesteps_. + int numberOfTimesteps() const { - return nTimesteps; + return nTimesteps_; } /// Sets \ref nTimesteps. void setNumberOfTimesteps(int num) { - nTimesteps = std::max(0, num); + nTimesteps_ = std::max(0, num); } - /// Increments \ref timestepNumber by 1; + /// Increments \ref timestepNumber_ by 1; void incTimestepNumber() { - timestepNumber++; + timestepNumber_++; } /// Sets \ref est_sum. void setEstSum(double e, Key key) { - getScalContent(key).est_sum = e; + scalContent(key).est_sum = e; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -297,7 +297,7 @@ namespace AMDiS /// Sets \ref est_max. void setEstMax(double e, Key key) { - getScalContent(key).est_max = e; + scalContent(key).est_max = e; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -309,7 +309,7 @@ namespace AMDiS /// Sets \ref est_max. void setTimeEstMax(double e, Key key) { - getScalContent(key).est_t_max = e; + scalContent(key).est_t_max = e; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -321,7 +321,7 @@ namespace AMDiS /// Sets \ref est_t_sum. void setTimeEstSum(double e, Key key) { - getScalContent(key).est_t_sum = e; + scalContent(key).est_t_sum = e; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -331,92 +331,92 @@ namespace AMDiS } /// Returns \ref est_sum. - double getEstSum(Key key) const + double estSum(Key key) const { - return getScalContent(key).est_sum; + return scalContent(key).est_sum; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getEstSum(const TP& tp) + double estSum(const TP& tp) { - return getEstSum(to_string(tp)); + return estSum(to_string(tp)); } /// Returns \ref est_t_sum. - double getEstTSum(Key key) const + double estTSum(Key key) const { - return getScalContent(key).est_t_sum; + return scalContent(key).est_t_sum; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getEstTSum(const TP& tp) + double estTSum(const TP& tp) { - return getEstTSum(to_string(tp)); + return estTSum(to_string(tp)); } /// Returns \ref est_max. - double getEstMax(Key key) const + double estMax(Key key) const { - return getScalContent(key).est_max; + return scalContent(key).est_max; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getEstMax(const TP& tp) + double estMax(const TP& tp) { - return getEstMax(to_string(tp)); + return estMax(to_string(tp)); } /// Returns \ref est_max. - double getTimeEstMax(Key key) const + double timeEstMax(Key key) const { - return getScalContent(key).est_t_max; + return scalContent(key).est_t_max; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeEstmax(const TP& tp) + double timeEstmax(const TP& tp) { - return getTimeEstMax(to_string(tp)); + return timeEstMax(to_string(tp)); } /// Returns \ref est_t_sum. - double getTimeEstSum(Key key) const + double timeEstSum(Key key) const { - return getScalContent(key).est_t_sum; + return scalContent(key).est_t_sum; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeEstSum(const TP& tp) + double timeEstSum(const TP& tp) { - return getTimeEstSum(to_string(tp)); + return timeEstSum(to_string(tp)); } - /// Returns \ref est_t the estimated overall time error - double getTimeEst() const + /// Returns \ref timeEst_ the estimated overall time error + double timeEst() const { - return est_t; + return timeEst_; } void setTimeEst(double value) { - est_t = value; + timeEst_ = value; } /// Returns \ref spaceTolerance. - double getSpaceTolerance(Key key) const + double spaceTolerance(Key key) const { - return getScalContent(key).spaceTolerance; + return scalContent(key).spaceTolerance; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getSpaceTolerance(const TP& tp) + double spaceTolerance(const TP& tp) { - return getSpaceTolerance(to_string(tp)); + return spaceTolerance(to_string(tp)); } /// Sets \ref spaceTolerance. void setSpaceTolerance(Key key, double tol) { - getScalContent(key).spaceTolerance = tol; + scalContent(key).spaceTolerance = tol; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -426,163 +426,151 @@ namespace AMDiS } /// Returns \ref timeTolerance. - double getTimeTolerance(Key key) const + double timeTolerance(Key key) const { - return getScalContent(key).timeTolerance; + return scalContent(key).timeTolerance; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeTolerance(const TP& tp) + double timeTolerance(const TP& tp) { - return getTimeTolerance(to_string(tp)); + return timeTolerance(to_string(tp)); } /// Returns \ref timeRelativeTolerance. - double getTimeRelativeTolerance(Key key) const + double timeRelativeTolerance(Key key) const { - return getScalContent(key).timeRelativeTolerance; + return scalContent(key).timeRelativeTolerance; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeRelativeTolerance(const TP& tp) + double timeRelativeTolerance(const TP& tp) { - return getTimeRelativeTolerance(to_string(tp)); + return timeRelativeTolerance(to_string(tp)); } - /// Sets \ref time + /// Sets \ref time_ double setTime(double t) { - time = t; - if (time > endTime) - time = endTime; - if (time < startTime) - time = startTime; + time_ = t; + if (time_ > endTime_) + time_ = endTime_; + if (time_ < startTime_) + time_ = startTime_; - return time; + return time_; } - /// Gets \ref time - double getTime() const + /// Gets \ref time_ + double const& time() const { - return time; + return time_; } - /// Gets \ref &time - double const* getTimePtr() const - { - return &time; - } - - /// Sets \ref timestep + /// Sets \ref timestep_ double setTimestep(double t) { - timestep = t; - if (timestep > maxTimestep) - timestep = maxTimestep; - if (timestep < minTimestep) - timestep = minTimestep; - if (time + timestep > endTime) - timestep = endTime - time; + timestep_ = t; + if (timestep_ > maxTimestep_) + timestep_ = maxTimestep_; + if (timestep_ < minTimestep_) + timestep_ = minTimestep_; + if (time_ + timestep_ > endTime_) + timestep_ = endTime_ - time_; - return timestep; + return timestep_; } - /// Gets \ref timestep - double getTimestep() const + /// Gets \ref timestep_ + double const& timestep() const { - return timestep; + return timestep_; } void setLastProcessedTimestep(double t) { - lastProcessedTimestep = t; + lastProcessedTimestep_ = t; } - double getLastProcessedTimestep() const + double lastProcessedTimestep() const { - return lastProcessedTimestep; + return lastProcessedTimestep_; } /// Returns true, if the end time is reached and no more timestep /// computations must be done. bool reachedEndTime() const { - if (nTimesteps > 0) - return !(timestepNumber < nTimesteps); + if (nTimesteps_ > 0) + return !(timestepNumber_ < nTimesteps_); - return !(std::abs(time - endTime) > std::numeric_limits<double>::epsilon()); + return !(std::abs(time_ - endTime_) > std::numeric_limits<double>::epsilon()); } - /// Sets \ref minTimestep + /// Sets \ref minTimestep_ void setMinTimestep(double t) { - minTimestep = t; + minTimestep_ = t; } - /// Gets \ref minTimestep - double getMinTimestep() const + /// Gets \ref minTimestep_ + double minTimestep() const { - return minTimestep; + return minTimestep_; } - /// Sets \ref maxTimestep - void setMaxTimestep(double t) + /// Sets \ref maxTimestep_ + void maxTimestep(double t) { - maxTimestep = t; + maxTimestep_ = t; } - /// Gets \ref maxTimestep - double getMaxTimestep() const + /// Gets \ref maxTimestep_ + double maxTimestep() const { - return maxTimestep; + return maxTimestep_; } - /// Gets \ref ×tep - double const* getTimestepPtr() const - { - return ×tep; - } - - /// Sets \ref startTime = time + /// Sets \ref startTime_ = time void setStartTime(double time) { - startTime = time; + startTime_ = time; } - /// Sets \ref endTime = time + /// Sets \ref endTime_ = time void setEndTime(double time) { - endTime = time; + endTime_ = time; } - /// Returns \ref startTime - double getStartTime() const + /// Returns \ref startTime_ + double startTime() const { - return startTime; + return startTime_; } - /// Returns \ref endTime - double getEndTime() const + /// Returns \ref endTime_ + double endTime() const { - return endTime; + return endTime_; } /// Returns \ref timeErrLow. - double getTimeErrLow(Key key) const + double timeErrLow(Key key) const { - return getScalContent(key).timeErrLow; + return scalContent(key).timeErrLow; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> - double getTimeErrLow(const TP& tp) + double timeErrLow(const TP& tp) { - return getTimeErrLow(to_string(tp)); + return timeErrLow(to_string(tp)); } /// Returns whether coarsening is allowed or not. bool isCoarseningAllowed(Key key) const { - return (getScalContent(key).coarsenAllowed == 1); + return (scalContent(key).coarsenAllowed == 1); } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -594,7 +582,7 @@ namespace AMDiS /// Returns whether coarsening is allowed or not. bool isRefinementAllowed(Key key) const { - return (getScalContent(key).refinementAllowed == 1); + return (scalContent(key).refinementAllowed == 1); } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -606,7 +594,7 @@ namespace AMDiS /// void allowRefinement(bool allow, Key key) { - getScalContent(key).refinementAllowed = allow; + scalContent(key).refinementAllowed = allow; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -618,7 +606,7 @@ namespace AMDiS /// void allowCoarsening(bool allow, Key key) { - getScalContent(key).coarsenAllowed = allow; + scalContent(key).coarsenAllowed = allow; } template <class TP, REQUIRES( Concepts::PreTreePath<TP> )> @@ -627,72 +615,61 @@ namespace AMDiS return allowCoarsening(allow, to_string(tp)); } - int getSize() const + int size() const { - return int(scalContents.size()); - } - - // TODO: remove from AdaptInfo - bool getRosenbrockMode() const - { - return rosenbrockMode; + return int(scalContents_.size()); } void setSolverIterations(int it) { - solverIterations = it; + solverIterations_ = it; } - int getSolverIterations() const + int solverIterations() const { - return solverIterations; + return solverIterations_; } void setMaxSolverIterations(int it) { - maxSolverIterations = it; + maxSolverIterations_ = it; } - int getMaxSolverIterations() const + int maxSolverIterations() const { - return maxSolverIterations; + return maxSolverIterations_; } void setSolverTolerance(double tol) { - solverTolerance = tol; + solverTolerance_ = tol; } - double getSolverTolerance() const + double solverTolerance() const { - return solverTolerance; + return solverTolerance_; } void setSolverResidual(double res) { - solverResidual = res; + solverResidual_ = res; } - double getSolverResidual() const + double solverResidual() const { - return solverResidual; + return solverResidual_; } void setGlobalTimeTolerance(double tol) { - globalTimeTolerance = tol; + globalTimeTolerance_ = tol; } - double getGlobalTimeTolerance() const + double globalTimeTolerance() const { - return globalTimeTolerance; + return globalTimeTolerance_; } - // TODO: remove from AdaptInfo - void setRosenbrockMode(bool b) - { - rosenbrockMode = b; - } /** \brief * Resets timestep, current time and time boundaries without @@ -702,17 +679,17 @@ namespace AMDiS double newStartTime, double newEndTime) { - time = newStartTime; - startTime = newStartTime; - endTime = newEndTime; - timestep = newTimeStep; - timestepNumber = 0; + time_ = newStartTime; + startTime_ = newStartTime; + endTime_ = newEndTime; + timestep_ = newTimeStep; + timestepNumber_ = 0; } private: - ScalContent& getScalContent(Key key) const + ScalContent& scalContent(Key key) const { - auto result = scalContents.emplace(std::piecewise_construct, std::forward_as_tuple(key), std::forward_as_tuple(name_ + "[" + key + "]") ); + auto result = scalContents_.emplace(std::piecewise_construct, std::forward_as_tuple(key), std::forward_as_tuple(name_ + "[" + key + "]") ); return result.first->second; } @@ -721,83 +698,77 @@ namespace AMDiS std::string name_; /// Current space iteration - int spaceIteration = -1; + int spaceIteration_ = -1; /** \brief * maximal allowed number of iterations of the adaptive procedure; if * maxIteration <= 0, no iteration bound is used */ - int maxSpaceIteration = -1; + int maxSpaceIteration_ = -1; /// Current timestep iteration - int timestepIteration = 0; + int timestepIteration_ = 0; /// Maximal number of iterations for choosing a timestep - int maxTimestepIteration = 30; + int maxTimestepIteration_ = 30; /// Current time iteration - int timeIteration = 0; + int timeIteration_ = 0; /// Maximal number of time iterations - int maxTimeIteration = 30; + int maxTimeIteration_ = 30; /// Actual time, end of time interval for current time step - double time = 0.0; + double time_ = 0.0; /// Initial time - double startTime = 0.0; + double startTime_ = 0.0; /// Final time - double endTime = 1.0; + double endTime_ = 1.0; /// Time step size to be used - double timestep = 0.0; + double timestep_ = 0.0; /// Last processed time step size of finished iteration - double lastProcessedTimestep = 0.0; + double lastProcessedTimestep_ = 0.0; /// Minimal step size - double minTimestep = 0.0; + double minTimestep_ = 0.0; /// Maximal step size - double maxTimestep = 1.0; + double maxTimestep_ = 1.0; /// Number of current time step - int timestepNumber = 0; + int timestepNumber_ = 0; /** \brief * Per default this value is 0 and not used. If it is set to a non-zero value, * the computation of the stationary problem is done nTimesteps times with a * fixed timestep. */ - int nTimesteps = 0; + int nTimesteps_ = 0; /// number of iterations needed of linear or nonlinear solver - int solverIterations = 0; + int solverIterations_ = 0; /// maximal number of iterations needed of linear or nonlinear solver - int maxSolverIterations = 0; + int maxSolverIterations_ = 0; /// - double solverTolerance = 1.e-8; + double solverTolerance_ = 1.e-8; /// - double solverResidual = 0.0; + double solverResidual_ = 0.0; /// tolerance for the overall time error - double globalTimeTolerance = 1.0; + double globalTimeTolerance_ = 1.0; /// Scalar adapt infos - mutable std::map<Key, ScalContent> scalContents; - - /// Is true, if the adaptive procedure was deserialized from a file. TODO: remove deserialization - bool deserialized = false; - - /// Is true, if the time adaption is controlled by a Rosenbrock Method. - bool rosenbrockMode = false; + mutable std::map<Key, ScalContent> scalContents_; /// overall time error estimate - double est_t = 0.0; + double timeEst_ = 0.0; }; } // end namespace AMDiS diff --git a/src/amdis/AdaptInstationary.cpp b/src/amdis/AdaptInstationary.cpp index 671d260bfb528cb78f7a9c3a1a1c4f19c0e4278d..add28e1de18789f05bfad69ac5dcbddc214c33a6 100644 --- a/src/amdis/AdaptInstationary.cpp +++ b/src/amdis/AdaptInstationary.cpp @@ -22,7 +22,7 @@ namespace AMDiS Parameters::get(name_ + "->time delta 2", timeDelta2_); Parameters::get(name_ + "->break when stable", breakWhenStable_); - fixedTimestep_ = (adaptInfo_.getMinTimestep() == adaptInfo_.getMaxTimestep()); + fixedTimestep_ = (adaptInfo_.minTimestep() == adaptInfo_.maxTimestep()); } @@ -31,16 +31,16 @@ namespace AMDiS AMDIS_FUNCNAME("AdaptInstationary::explicitTimeStrategy()"); // estimate before first adaption - if (adaptInfo_.getTime() <= adaptInfo_.getStartTime()) + if (adaptInfo_.time() <= adaptInfo_.startTime()) problemIteration_->oneIteration(adaptInfo_, ESTIMATE); // increment time - adaptInfo_.setTime(adaptInfo_.getTime() + adaptInfo_.getTimestep()); + adaptInfo_.setTime(adaptInfo_.time() + adaptInfo_.timestep()); problemTime_->setTime(adaptInfo_); - msg("time = {}, timestep = {}", adaptInfo_.getTime(), adaptInfo_.getTimestep()); + msg("time = {}, timestep = {}", adaptInfo_.time(), adaptInfo_.timestep()); adaptInfo_.setSpaceIteration(0); @@ -48,7 +48,7 @@ namespace AMDiS problemIteration_->beginIteration(adaptInfo_); problemIteration_->oneIteration(adaptInfo_, FULL_ITERATION); problemIteration_->endIteration(adaptInfo_); - adaptInfo_.setLastProcessedTimestep(adaptInfo_.getTimestep()); + adaptInfo_.setLastProcessedTimestep(adaptInfo_.timestep()); } @@ -58,10 +58,10 @@ namespace AMDiS do { - adaptInfo_.setTime(adaptInfo_.getTime() + adaptInfo_.getTimestep()); + adaptInfo_.setTime(adaptInfo_.time() + adaptInfo_.timestep()); problemTime_->setTime(adaptInfo_); - msg("time = {}, timestep = {}", adaptInfo_.getTime(), adaptInfo_.getTimestep()); + msg("time = {}, timestep = {}", adaptInfo_.time(), adaptInfo_.timestep()); problemIteration_->oneIteration(adaptInfo_, NO_ADAPTION); @@ -69,11 +69,11 @@ namespace AMDiS if (!fixedTimestep_ && !adaptInfo_.timeToleranceReached() && - adaptInfo_.getTimestepIteration() <= adaptInfo_.getMaxTimestepIteration() && - !(adaptInfo_.getTimestep() <= adaptInfo_.getMinTimestep())) + adaptInfo_.timestepIteration() <= adaptInfo_.maxTimestepIteration() && + !(adaptInfo_.timestep() <= adaptInfo_.minTimestep())) { - adaptInfo_.setTime(adaptInfo_.getTime() - adaptInfo_.getTimestep()); - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta1_); + adaptInfo_.setTime(adaptInfo_.time() - adaptInfo_.timestep()); + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta1_); continue; } @@ -83,7 +83,7 @@ namespace AMDiS // === Do space iterations only if the maximum is higher than 0. === - if (adaptInfo_.getMaxSpaceIteration() > 0) + if (adaptInfo_.maxSpaceIteration() > 0) { // === Space iterations. === @@ -96,10 +96,10 @@ namespace AMDiS { if (!fixedTimestep_ && !adaptInfo_.timeToleranceReached() && - !(adaptInfo_.getTimestep() <= adaptInfo_.getMinTimestep())) + !(adaptInfo_.timestep() <= adaptInfo_.minTimestep())) { - adaptInfo_.setTime(adaptInfo_.getTime() - adaptInfo_.getTimestep()); - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta2_); + adaptInfo_.setTime(adaptInfo_.time() - adaptInfo_.timestep()); + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta2_); problemIteration_->endIteration(adaptInfo_); adaptInfo_.incSpaceIteration(); break; @@ -111,7 +111,7 @@ namespace AMDiS } while (!adaptInfo_.spaceToleranceReached() && - adaptInfo_.getSpaceIteration() <= adaptInfo_.getMaxSpaceIteration()); + adaptInfo_.spaceIteration() <= adaptInfo_.maxSpaceIteration()); } else @@ -122,22 +122,22 @@ namespace AMDiS } while(!adaptInfo_.timeToleranceReached() && - !(adaptInfo_.getTimestep() <= adaptInfo_.getMinTimestep()) && - adaptInfo_.getTimestepIteration() <= adaptInfo_.getMaxTimestepIteration()); + !(adaptInfo_.timestep() <= adaptInfo_.minTimestep()) && + adaptInfo_.timestepIteration() <= adaptInfo_.maxTimestepIteration()); - adaptInfo_.setLastProcessedTimestep(adaptInfo_.getTimestep()); + adaptInfo_.setLastProcessedTimestep(adaptInfo_.timestep()); // After successful iteration/timestep the timestep will be changed according // adaption rules for next timestep. // First, check for increase of timestep if (!fixedTimestep_ && adaptInfo_.timeErrorLow()) - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta2_); + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta2_); // Second, check for decrease of timestep if (!fixedTimestep_ && !adaptInfo_.timeToleranceReached() && - !(adaptInfo_.getTimestep() <= adaptInfo_.getMinTimestep())) - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta1_); + !(adaptInfo_.timestep() <= adaptInfo_.minTimestep())) + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta1_); } @@ -146,27 +146,27 @@ namespace AMDiS AMDIS_FUNCNAME("AdaptInstationary::simpleAdaptiveTimeStrategy()"); // estimate before first adaption - if (adaptInfo_.getTime() <= adaptInfo_.getStartTime()) + if (adaptInfo_.time() <= adaptInfo_.startTime()) problemIteration_->oneIteration(adaptInfo_, ESTIMATE); - adaptInfo_.setTime(adaptInfo_.getTime() + adaptInfo_.getTimestep()); + adaptInfo_.setTime(adaptInfo_.time() + adaptInfo_.timestep()); problemTime_->setTime(adaptInfo_); - msg("time = {}, timestep = {}", adaptInfo_.getTime(), adaptInfo_.getTimestep()); + msg("time = {}, timestep = {}", adaptInfo_.time(), adaptInfo_.timestep()); problemIteration_->oneIteration(adaptInfo_, FULL_ITERATION); - adaptInfo_.setLastProcessedTimestep(adaptInfo_.getTimestep()); + adaptInfo_.setLastProcessedTimestep(adaptInfo_.timestep()); // First, check for increase of timestep if (!fixedTimestep_ && adaptInfo_.timeErrorLow()) - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta2_); + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta2_); // Second, check for decrease of timestep if (!fixedTimestep_ && !adaptInfo_.timeToleranceReached() && - !(adaptInfo_.getTimestep() <= adaptInfo_.getMinTimestep())) - adaptInfo_.setTimestep(adaptInfo_.getTimestep() * timeDelta1_); + !(adaptInfo_.timestep() <= adaptInfo_.minTimestep())) + adaptInfo_.setTimestep(adaptInfo_.timestep() * timeDelta1_); } @@ -200,18 +200,18 @@ namespace AMDiS AMDIS_FUNCNAME("AdaptInstationary::adapt()"); int errorCode = 0; - test_exit(adaptInfo_.getTimestep() >= adaptInfo_.getMinTimestep(), + test_exit(adaptInfo_.timestep() >= adaptInfo_.minTimestep(), "timestep < min timestep"); - test_exit(adaptInfo_.getTimestep() <= adaptInfo_.getMaxTimestep(), + test_exit(adaptInfo_.timestep() <= adaptInfo_.maxTimestep(), "timestep > max timestep"); - test_exit(adaptInfo_.getTimestep() > 0, "timestep <= 0!"); + test_exit(adaptInfo_.timestep() > 0, "timestep <= 0!"); - if (adaptInfo_.getTimestepNumber() == 0) + if (adaptInfo_.timestepNumber() == 0) { - adaptInfo_.setTime(adaptInfo_.getStartTime()); - initialAdaptInfo_->setStartTime(adaptInfo_.getStartTime()); - initialAdaptInfo_->setTime(adaptInfo_.getStartTime()); + adaptInfo_.setTime(adaptInfo_.startTime()); + initialAdaptInfo_->setStartTime(adaptInfo_.startTime()); + initialAdaptInfo_->setTime(adaptInfo_.startTime()); problemTime_->setTime(adaptInfo_); @@ -226,7 +226,7 @@ namespace AMDiS oneTimestep(); problemTime_->closeTimestep(adaptInfo_); - if (breakWhenStable_ && (adaptInfo_.getSolverIterations() == 0)) + if (breakWhenStable_ && (adaptInfo_.solverIterations() == 0)) break; } diff --git a/src/amdis/AdaptInstationary.hpp b/src/amdis/AdaptInstationary.hpp index a4f89664d72b012fa7d9b5a2fd8c573b2f09f218..48af810a09f7beee315a59532d7de9e581b79272 100644 --- a/src/amdis/AdaptInstationary.hpp +++ b/src/amdis/AdaptInstationary.hpp @@ -38,7 +38,7 @@ namespace AMDiS } /// Returns \ref strategy - int getStrategy() const + int strategy() const { return strategy_; } diff --git a/src/amdis/AdaptStationary.cpp b/src/amdis/AdaptStationary.cpp index 87964d621762fe8a01318fc9308ee0a3a7c33b8a..eb721da0c6cdd675420c68a6cd8bc4b5c043f0f1 100644 --- a/src/amdis/AdaptStationary.cpp +++ b/src/amdis/AdaptStationary.cpp @@ -18,7 +18,7 @@ AdaptStationary::AdaptStationary(std::string const& name, int AdaptStationary::adapt() { // initial iteration - if (adaptInfo_.getSpaceIteration() == -1) + if (adaptInfo_.spaceIteration() == -1) { problemIteration_->beginIteration(adaptInfo_); problemIteration_->oneIteration(adaptInfo_, NO_ADAPTION); @@ -28,8 +28,8 @@ int AdaptStationary::adapt() // adaption loop while (!adaptInfo_.spaceToleranceReached() && - (adaptInfo_.getSpaceIteration() < adaptInfo_.getMaxSpaceIteration() || - adaptInfo_.getMaxSpaceIteration() < 0) ) + (adaptInfo_.spaceIteration() < adaptInfo_.maxSpaceIteration() || + adaptInfo_.maxSpaceIteration() < 0) ) { problemIteration_->beginIteration(adaptInfo_); diff --git a/src/amdis/Assembler.hpp b/src/amdis/Assembler.hpp index f1d6834ef158a795099657307324bde129c86be4..7452ea985a1cf682ab9c008f25d4eda85da22d7a 100644 --- a/src/amdis/Assembler.hpp +++ b/src/amdis/Assembler.hpp @@ -10,17 +10,17 @@ namespace AMDiS ElementAssembler const& localAssembler) { // assemble element operators - localAssembler(element, operators.element); + localAssembler(element, operators.onElement()); // assemble intersection operators - if (!operators.intersection.empty() - || (!operators.boundary.empty() && element.hasBoundaryIntersections())) + if (!operators.onIntersection().empty() + || (!operators.onBoundary().empty() && element.hasBoundaryIntersections())) { for (auto const& intersection : intersections(gridView, element)) { if (intersection.boundary()) - localAssembler(intersection, operators.boundary); + localAssembler(intersection, operators.onBoundary()); else - localAssembler(intersection, operators.intersection); + localAssembler(intersection, operators.onIntersection()); } } } diff --git a/src/amdis/FileWriter.hpp b/src/amdis/FileWriter.hpp index 6f6488af74291444aa7693d4a0935db3961f05ea..6951a0239d628bb28c1a4738fb80b469f25c7ae8 100644 --- a/src/amdis/FileWriter.hpp +++ b/src/amdis/FileWriter.hpp @@ -99,7 +99,7 @@ namespace AMDiS test_exit(filesystem::exists(dir_), "Output directory '{}' does not exist!",dir_); if (vtkSeqWriter_) - vtkSeqWriter_->write(adaptInfo.getTime(), mode_); + vtkSeqWriter_->write(adaptInfo.time(), mode_); else if (vtkWriter_) vtkWriter_->write(filesystem::path({dir_, filename_}).string(), mode_); } diff --git a/src/amdis/Flag.hpp b/src/amdis/Flag.hpp index 9ac846d4ae2b8c0440840ac0b86308297076af2e..5832b8593053ddfc1a6f84c6e8db4db323efb62e 100644 --- a/src/amdis/Flag.hpp +++ b/src/amdis/Flag.hpp @@ -18,7 +18,7 @@ namespace AMDiS /// Constructs a Flag initialized by f constexpr Flag(const std::uint64_t f) - : flags(f) + : flags_(f) {} /// Copy constructor @@ -30,7 +30,7 @@ namespace AMDiS /// Compares two Flags constexpr bool operator==(Flag const& f) const { - return (flags == f.flags); + return (flags_ == f.flags_); } /// Compares two Flags @@ -43,7 +43,7 @@ namespace AMDiS constexpr Flag& operator=(Flag const& f) { if (this != &f) - flags = f.flags; + flags_ = f.flags_; return *this; } @@ -53,109 +53,109 @@ namespace AMDiS return isAnySet(); } - /// Set \ref flags + /// Set \ref flags_ constexpr void setFlags(const std::uint64_t f) { - flags = f; + flags_ = f; } - /// Set \ref flags + /// Set \ref flags_ constexpr void setFlags(Flag const& f) { - flags = f.flags; + flags_ = f.flags_; } - /// Sets \ref flags to \ref flags | f + /// Sets \ref flags_ to \ref flags_ | f constexpr void setFlag(const std::uint64_t f) { - flags |= f; + flags_ |= f; } - /// Sets \ref flags to \ref flags | f.flags + /// Sets \ref flags_ to \ref flags_ | f.flags_ constexpr void setFlag(Flag const& f) { - flags |= f.flags; + flags_ |= f.flags_; } - /// Sets \ref flags to \ref flags & ~f + /// Sets \ref flags_ to \ref flags_ & ~f constexpr void unsetFlag(const std::uint64_t f) { - flags &= ~f; + flags_ &= ~f; } - /// Sets \ref flags to \ref flags & ~f.flags + /// Sets \ref flags_ to \ref flags_ & ~f.flags_ constexpr void unsetFlag(Flag const& f) { - flags &= ~f.flags; + flags_ &= ~f.flags_; } - constexpr std::uint64_t getFlags() const + constexpr std::uint64_t flags() const { - return flags; + return flags_; } - /// Returns \ref flags | f.flags + /// Returns \ref flags_ | f.flags_ friend Flag operator+(Flag r, Flag const& f) { r.setFlag(f); return r; } - /// Returns \ref flags & ~f.flags + /// Returns \ref flags_ & ~f.flags_ friend Flag operator-(Flag r, Flag const& f) { r.unsetFlag(f); return r; } - /// Returns \ref flags | f.flags + /// Returns \ref flags_ | f.flags_ friend Flag operator|(Flag r, Flag const& f) { r.setFlag(f); return r; } - /// Returns \ref flags & f.flags + /// Returns \ref flags_ & f.flags_ friend Flag operator&(Flag r, Flag const& f) { - r.flags &= f.flags; + r.flags_ &= f.flags_; return r; } - /// Sets \ref flags to \ref flags &= f.flags + /// Sets \ref flags_ to \ref flags_ &= f.flags_ constexpr Flag operator&=(Flag const& f) { - flags &= f.flags; + flags_ &= f.flags_; return *this; } - /// Returns \ref flags ^ f.flags + /// Returns \ref flags_ ^ f.flags_ friend Flag operator^(Flag r, Flag const& f) { - r.flags ^= f.flags; + r.flags_ ^= f.flags_; return r; } - /// Sets \ref flags to \ref flags & f.flags + /// Sets \ref flags_ to \ref flags_ & f.flags_ constexpr Flag& operator|=(Flag const& f) { if (this != &f) - flags |= f.flags; + flags_ |= f.flags_; return *this; } - /// Returns ~\ref flags + /// Returns ~\ref flags_ constexpr Flag operator~() const { Flag r; - r.flags = ~flags; + r.flags_ = ~flags_; return r; } - /// Checks whether all set bits of f.flags are set in \ref flags too. + /// Checks whether all set bits of f.flags_ are set in \ref flags_ too. constexpr bool isSet(Flag const& f) const { - return ((flags&f.flags) == f.flags); + return ((flags_ & f.flags_) == f.flags_); } /// Returns !\ref isSet(f) @@ -164,15 +164,15 @@ namespace AMDiS return !isSet(f); } - /// Returns true if \ref flags != 0 + /// Returns true if \ref flags_ != 0 constexpr bool isAnySet() const { - return (flags != 0); + return (flags_ != 0); } protected: /// Internal flag representation - std::uint64_t flags = 0; + std::uint64_t flags_ = 0; }; } // end namespace AMDiS diff --git a/src/amdis/GridFunctions.hpp b/src/amdis/GridFunctions.hpp index 984f7cf109295ac7d0df727173e5fff2e5330b02..a983cad45aef2ecb9c283dc6767b44cb81a72217 100644 --- a/src/amdis/GridFunctions.hpp +++ b/src/amdis/GridFunctions.hpp @@ -32,7 +32,7 @@ * * 3. Interpolate a GridFunction to a DOFVector: * ``` - * prob.getSolution(_0).interpol(Expression); + * prob.solution(_0).interpol(Expression); * ``` * * 4. Integrate a GridFunction on a GridView: @@ -55,7 +55,7 @@ * polynomial degree of the expression, or a quadrature rule explicitly. * * *Examples:* - * + `auto op1 = makeOperator(B, 1.0 + pow<2>(prob.getSolution(_0)));` + * + `auto op1 = makeOperator(B, 1.0 + pow<2>(prob.solution(_0)));` * + `auto op2 = makeOperator(B, sin(X(0)), 4);` * + `auto op3 = makeOperator(B, sin(X(0)), Dune::QuadratureRules(Dune::GeometryType::simplex, 4));` * + `auto value1 = integrate(sin(X(0)), 4);` diff --git a/src/amdis/LocalAssembler.hpp b/src/amdis/LocalAssembler.hpp index 1d47532dfd273c401c1c5182483ba0e19ecf9ba2..66ca9b1f7d4b984d4aa7f060ebcfd69c023497a5 100644 --- a/src/amdis/LocalAssembler.hpp +++ b/src/amdis/LocalAssembler.hpp @@ -73,7 +73,7 @@ namespace AMDiS Nodes const&... nodes, ElementMatrixVector& elementMatrixVector) final { - ContextGeometry<LocalContext> context{localContext, getElement(), getGeometry()}; + ContextGeometry<LocalContext> context{localContext, element(), geometry()}; assembleImpl(context, nodes..., elementMatrixVector); } @@ -105,14 +105,14 @@ namespace AMDiS public: /// return the bound entity (of codim 0) - Element const& getElement() const + Element const& element() const { assert( element_ ); return *element_; } /// return the geometry of the bound element - Geometry const& getGeometry() const + Geometry const& geometry() const { assert( geometry_ ); return *geometry_; diff --git a/src/amdis/LocalAssemblerList.hpp b/src/amdis/LocalAssemblerList.hpp index 59f4f01bc4c12892d3aa3468e8d4c2eaf57d6410..b9e7fc694542bcdff421069de9050eedbd222729 100644 --- a/src/amdis/LocalAssemblerList.hpp +++ b/src/amdis/LocalAssemblerList.hpp @@ -33,72 +33,91 @@ namespace AMDiS template <class... Nodes> struct Data { - /// The type of local operators associated with grid elements - using ElementOperator = LocalAssemblerBase<Element, Nodes...>; - /// The type of local operators associated with grid intersections - using IntersectionOperator = LocalAssemblerBase<Intersection, Nodes...>; - /// Return whether any operators are stored on the node bool empty() const { - return element.empty() && boundary.empty() && intersection.empty(); + return element_.empty() && boundary_.empty() && intersection_.empty(); } /// Test whether to assemble on the node bool doAssemble() const { - return !assembled || changing; + return !assembled_ || changing_; + } + + operator bool() const + { + return doAssemble() && !empty(); } /// Bind all operators to the grid element and geometry template <class Geo> void bind(Element const& elem, Geo const& geo) { - for (auto& scaled : element) scaled.op->bind(elem,geo); - for (auto& scaled : boundary) scaled.op->bind(elem,geo); - for (auto& scaled : intersection) scaled.op->bind(elem,geo); + for (auto& scaled : element_) scaled.op->bind(elem,geo); + for (auto& scaled : boundary_) scaled.op->bind(elem,geo); + for (auto& scaled : intersection_) scaled.op->bind(elem,geo); } /// Unbind all operators from the element void unbind() { - for (auto& scaled : element) scaled.op->unbind(); - for (auto& scaled : boundary) scaled.op->unbind(); - for (auto& scaled : intersection) scaled.op->unbind(); + for (auto& scaled : element_) scaled.op->unbind(); + for (auto& scaled : boundary_) scaled.op->unbind(); + for (auto& scaled : intersection_) scaled.op->unbind(); + assembled_ = true; } template <class... Args> void push(tag::element_operator<Element>, Args&&... args) { - element.push_back({std::forward<Args>(args)...}); + element_.push_back({std::forward<Args>(args)...}); } template <class... Args> void push(tag::boundary_operator<Intersection> b, Args&&... args) { - boundary.push_back({std::forward<Args>(args)..., b.id}); + boundary_.push_back({std::forward<Args>(args)..., b.id}); } template <class... Args> void push(tag::intersection_operator<Intersection>, Args&&... args) { - intersection.push_back({std::forward<Args>(args)...}); + intersection_.push_back({std::forward<Args>(args)...}); } + auto& onElement() { return element_; } + auto const& onElement() const { return element_; } + + auto& onBoundary() { return boundary_; } + auto const& onBoundary() const { return boundary_; } + + auto& onIntersection() { return intersection_; } + auto const& onIntersection() const { return intersection_; } + + private: + /// The type of local operators associated with grid elements + using ElementOperator = LocalAssemblerBase<Element, Nodes...>; + + /// The type of local operators associated with grid intersections + using IntersectionOperator = LocalAssemblerBase<Intersection, Nodes...>; + /// List of operators to be assembled on grid elements - std::list<DataElement<ElementOperator>> element; + std::vector<DataElement<ElementOperator>> element_; + /// List of operators to be assembled on boundary intersections - std::list<DataElement<IntersectionOperator>> boundary; + std::vector<DataElement<IntersectionOperator>> boundary_; + /// List of operators to be assembled on interior intersections - std::list<DataElement<IntersectionOperator>> intersection; + std::vector<DataElement<IntersectionOperator>> intersection_; /// if false, do reassemble of all operators - bool assembled = false; + bool assembled_ = false; /// if true, or \ref assembled false, do reassemble of all operators - bool changing = true; + bool changing_ = true; }; public: diff --git a/src/amdis/Marker.hpp b/src/amdis/Marker.hpp index fc346c79687c30e1da5a111885cade954cfa3f49..2ec8ed17673b2932048d264b5b3d15ffdb2f2563 100644 --- a/src/amdis/Marker.hpp +++ b/src/amdis/Marker.hpp @@ -65,19 +65,19 @@ namespace AMDiS virtual Flag markGrid(AdaptInfo& adaptInfo); /// Returns \ref elMarkRefine_ of the Marker - int getElMarkRefine() const + int elMarkRefine() const { return elMarkRefine_; } /// Returns \ref elMarkCoarsen_ of the Marker - int getElMarkCoarsen() const + int elMarkCoarsen() const { return elMarkCoarsen_; } /// Returns \ref name_ of the Marker - std::string getName() const + std::string const& name() const { return name_; } @@ -399,6 +399,7 @@ namespace AMDiS GridFct gridFct_; }; + #if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION // Deduction guide for GridFunctionMarker class template <class Grid, class PreGridFct> diff --git a/src/amdis/Marker.inc.hpp b/src/amdis/Marker.inc.hpp index 61bfee4647f5425aeb7cf10a8540b976a3264dcb..cc4687a848be2c4a204a64d1753ed5cd7627f305 100644 --- a/src/amdis/Marker.inc.hpp +++ b/src/amdis/Marker.inc.hpp @@ -74,8 +74,8 @@ template <class Grid> void EstimatorMarker<Grid>::initMarking(AdaptInfo& adaptInfo) { Super::initMarking(adaptInfo); - estSum_ = std::pow(adaptInfo.getEstSum(component_), p_); - estMax_ = adaptInfo.getEstMax(component_); + estSum_ = std::pow(adaptInfo.estSum(component_), p_); + estMax_ = adaptInfo.estMax(component_); this->refineAllowed_ = adaptInfo.isRefinementAllowed(component_); this->coarsenAllowed_ = adaptInfo.isCoarseningAllowed(component_); } @@ -136,11 +136,11 @@ void MSMarker<Grid>::initMarking(AdaptInfo& adaptInfo) double msGammaP = std::pow(msGamma_, this->p_); double msGammaCP = std::pow(msGammaC_, this->p_); - this->markRLimit_ = msGammaP * adaptInfo.getEstMax(this->component_); - this->markCLimit_ = msGammaCP * adaptInfo.getEstMax(this->component_); + this->markRLimit_ = msGammaP * adaptInfo.estMax(this->component_); + this->markCLimit_ = msGammaCP * adaptInfo.estMax(this->component_); msg("start max_est: {}, mark_limits: {}, {}", - adaptInfo.getEstMax(this->component_), this->markRLimit_ , this->markCLimit_); + adaptInfo.estMax(this->component_), this->markRLimit_ , this->markCLimit_); } @@ -151,7 +151,7 @@ void ESMarker<Grid>::initMarking(AdaptInfo& adaptInfo) double esThetaP = std::pow(esTheta_, this->p_); double esThetaCP = std::pow(esThetaC_, this->p_); - double epsP = std::pow(adaptInfo.getSpaceTolerance(this->component_), this->p_); + double epsP = std::pow(adaptInfo.spaceTolerance(this->component_), this->p_); int nLeaves = (this->grid_->leafGridView()).size(0); #if AMDIS_HAS_PARALLEL @@ -176,7 +176,7 @@ Flag GERSMarker<Grid>::markGrid(AdaptInfo& adaptInfo) gersSum_ = 0.0; double LTheta = std::pow(1.0 - gersThetaStar_, this->p_); - double epsP = std::pow(adaptInfo.getSpaceTolerance(this->component_), this->p_); + double epsP = std::pow(adaptInfo.spaceTolerance(this->component_), this->p_); if (this->estSum_ < oldErrSum_) { double improv = this->estSum_ / oldErrSum_; diff --git a/src/amdis/ProblemInstat.hpp b/src/amdis/ProblemInstat.hpp index 9be494845a43915fa7dbc6595be04fb349ad4b26..60cd0d1eb68e77289f324967cab37c8408e581b4 100644 --- a/src/amdis/ProblemInstat.hpp +++ b/src/amdis/ProblemInstat.hpp @@ -52,19 +52,11 @@ namespace AMDiS virtual void closeTimestep(AdaptInfo& adaptInfo) override; /// Returns \ref problemStat. - ProblemType const& getProblemStat() const - { - return problemStat_; - } - - /// Returns \ref problemStat. - ProblemType& getProblemStat() - { - return problemStat_; - } + ProblemType& problemStat() { return problemStat_; } + ProblemType const& problemStat() const { return problemStat_; } /// Returns \ref oldSolution. - std::unique_ptr<SystemVector> const& getOldSolutionVector() const + SystemVector const& oldSolutionVector() const { test_exit_dbg(oldSolution_, "OldSolution need to be created. Call initialize with INIT_UH_OLD."); @@ -73,18 +65,18 @@ namespace AMDiS /// Return a mutable view to a oldSolution component template <class TreePath = RootTreePath> - auto getOldSolution(TreePath const& path = {}) + auto oldSolution(TreePath path = {}) { auto&& tp = makeTreePath(path); - return makeDOFVectorView(*getOldSolutionVector(), tp); + return makeDOFVectorView(oldSolutionVector(), tp); } /// Return a const view to a oldSolution component template <class TreePath = RootTreePath> - auto getOldSolution(TreePath const& path = {}) const + auto oldSolution(TreePath path = {}) const { auto&& tp = makeTreePath(path); - return makeDOFVectorView(*getOldSolutionVector(), tp); + return makeDOFVectorView(oldSolutionVector(), tp); } /// Implementation of \ref ProblemTimeInterface::transferInitialSolution(). diff --git a/src/amdis/ProblemInstat.inc.hpp b/src/amdis/ProblemInstat.inc.hpp index 3ae5156d4bd5ab7f0258a9e2b6cce4ddf505c257..bd1885a8db63a9be853da41b357a2e6a77e60790 100644 --- a/src/amdis/ProblemInstat.inc.hpp +++ b/src/amdis/ProblemInstat.inc.hpp @@ -14,7 +14,7 @@ void ProblemInstat<Traits>::transferInitialSolution(AdaptInfo& adaptInfo) { AMDIS_FUNCNAME("ProblemInstat::transferInitialSolution()"); - test_exit(adaptInfo.getTime() == adaptInfo.getStartTime(), + test_exit(adaptInfo.time() == adaptInfo.startTime(), "after initial solution: time != start time"); problemStat_.writeFiles(adaptInfo, true); } @@ -23,7 +23,7 @@ void ProblemInstat<Traits>::transferInitialSolution(AdaptInfo& adaptInfo) template <class Traits> void ProblemInstat<Traits>::closeTimestep(AdaptInfo& adaptInfo) { - bool force = (adaptInfo.getTime() >= adaptInfo.getEndTime()); + bool force = (adaptInfo.time() >= adaptInfo.endTime()); problemStat_.writeFiles(adaptInfo, force); } @@ -53,7 +53,7 @@ template <class Traits> void ProblemInstat<Traits>::initTimestep(AdaptInfo&) { if (oldSolution_) - *oldSolution_ = problemStat_.getSolutionVector(); + *oldSolution_ = problemStat_.solutionVector(); } } // end namespace AMDiS diff --git a/src/amdis/ProblemInstatBase.cpp b/src/amdis/ProblemInstatBase.cpp index 7fceb390689fad23dad60ec5968dc22572c2b735..6d5db1c4f11a9bfb0739bc08ae35ca889bb5bb53 100644 --- a/src/amdis/ProblemInstatBase.cpp +++ b/src/amdis/ProblemInstatBase.cpp @@ -8,8 +8,8 @@ namespace AMDiS { void ProblemInstatBase::setTime(AdaptInfo& adaptInfo) { - time_ = adaptInfo.getTime(); - tau_ = adaptInfo.getTimestep(); + time_ = adaptInfo.time(); + tau_ = adaptInfo.timestep(); invTau_ = 1.0 / tau_; } diff --git a/src/amdis/ProblemInstatBase.hpp b/src/amdis/ProblemInstatBase.hpp index de055180d5618f3aee83b77c20fd1aeb728c773b..84b644e82e257a6fc57b17b3a6f9a150148ba589 100644 --- a/src/amdis/ProblemInstatBase.hpp +++ b/src/amdis/ProblemInstatBase.hpp @@ -32,7 +32,7 @@ namespace AMDiS /// Implementation of \ref ProblemTimeInterface::setTime(). virtual void setTime(AdaptInfo& adaptInfo) override; - void solve(AdaptInfo&) { /* by default, do nothing */ } + void solve(AdaptInfo&) { /* do nothing */ } /// Implementation of \ref ProblemStatBase::solve(). virtual void solve(AdaptInfo& adaptInfo, bool, bool) override @@ -41,10 +41,10 @@ namespace AMDiS } /// Implementation of \ref ProblemStatBase::estimate(). - virtual void estimate(AdaptInfo&) override { /* by default, do nothing */ } + virtual void estimate(AdaptInfo&) override { /* do nothing */ } /// Implementation of \ref ProblemStatBase::buildAfterAdapt(). - virtual void buildAfterAdapt(AdaptInfo&, Flag, bool, bool) override { /* by default, do nothing */ } + virtual void buildAfterAdapt(AdaptInfo&, Flag, bool, bool) override { /* do nothing */ } /// Implementation of \ref ProblemStatBase::markElements(). virtual Flag markElements(AdaptInfo&) override @@ -59,10 +59,10 @@ namespace AMDiS } /// Implementation of \ref ProblemTimeInterface::closeTimestep(). - virtual void closeTimestep(AdaptInfo&) override { /* by default, do nothing */ } + virtual void closeTimestep(AdaptInfo&) override { /* do nothing */ } - /// Implementation of \ref ProblemStatBase::getName(). - virtual std::string const& getName() const override + /// Implementation of \ref ProblemStatBase::name(). + virtual std::string const& name() const override { return name_; } @@ -70,25 +70,25 @@ namespace AMDiS /// Implementation of \ref ProblemTimeInterface::solveInitialProblem(). virtual void solveInitialProblem(AdaptInfo& adaptInfo) override; - /// Return pointer to current simulation time \ref cTime set in \ref setTime + /// Return reference to current simulation time \ref time_ set in \ref setTime /// from `adaptInfo->getTime()`. - double const* getTime() const + double const& time() const { - return &time_; + return time_; } - /// Return pointer to current simulation timestep \ref tau set in \ref setTime + /// Return reference to current simulation timestep \ref tau_ set in \ref setTime /// from `adaptInfo->getTimestep()`. - double const* getTau() const& + double const& tau() const& { - return &tau_; + return tau_; } - /// Return pointer to current simulation 1.0/timestep \ref invTau set in + /// Return reference to current simulation 1.0/timestep \ref invTau_ set in /// \ref setTime from `1.0 / adaptInfo->getTimestep()`. - double const* getInvTau() const + double const& invTau() const { - return &invTau_; + return invTau_; } protected: diff --git a/src/amdis/ProblemIterationInterface.hpp b/src/amdis/ProblemIterationInterface.hpp index 76f10fe0dbdd5228dd683e22fc4fc2508d2e822e..b5ebd4b1d5c26527695f2f4badaa0c3a1c665678 100644 --- a/src/amdis/ProblemIterationInterface.hpp +++ b/src/amdis/ProblemIterationInterface.hpp @@ -45,19 +45,19 @@ namespace AMDiS virtual void endIteration(AdaptInfo&) { /* by default, do nothing */ } /// Returns number of managed problems - virtual int getNumProblems() const = 0; + virtual int numProblems() const = 0; /** \brief * Returns the problem with the given number. If only one problem * is managed by this master problem, the number hasn't to be given. */ - virtual ProblemStatBase& getProblem(int number = 0) = 0; + virtual ProblemStatBase& problem(int number = 0) = 0; /// Returns the problem with the given name. - virtual ProblemStatBase& getProblem(std::string const& name) = 0; + virtual ProblemStatBase& problem(std::string const& name) = 0; /// Returns the name of the problem. - virtual std::string const& getName() const = 0; + virtual std::string const& name() const = 0; }; } // end namespace AMDiS diff --git a/src/amdis/ProblemStat.cpp b/src/amdis/ProblemStat.cpp index 19e360bb54c01bf2d6767eea60fa077ceee794f4..691db3259d3e93732950d44600f1e8bcd54129ed 100644 --- a/src/amdis/ProblemStat.cpp +++ b/src/amdis/ProblemStat.cpp @@ -5,8 +5,5 @@ namespace AMDiS { - // explicit template instatiation - // template class ProblemStat<YaspGridBasis<2,1>>; - // template class ProblemStat<YaspGridBasis<2,2>>; - // template class ProblemStat<YaspGridBasis<2,1,2>>; + } // end namespace AMDiS diff --git a/src/amdis/ProblemStat.hpp b/src/amdis/ProblemStat.hpp index 5a737a44d2aa2bc6c66bc6cfb8398455417b0aa2..61b8e0f537e31f4ea44bcc3b836a9e875292e099 100644 --- a/src/amdis/ProblemStat.hpp +++ b/src/amdis/ProblemStat.hpp @@ -58,9 +58,9 @@ namespace AMDiS public: // typedefs and static constants using GlobalBasis = typename Traits::GlobalBasis; - using GridView = typename GlobalBasis::GridView; - using Grid = typename GridView::Grid; - using Element = typename GridView::template Codim<0>::Entity; + using GridView = typename GlobalBasis::GridView; + using Grid = typename GridView::Grid; + using Element = typename GridView::template Codim<0>::Entity; using WorldVector = typename Element::Geometry::GlobalCoordinate; /// Dimension of the grid @@ -89,8 +89,7 @@ namespace AMDiS ProblemStat(std::string const& name, Grid& grid) : ProblemStat(name) { - grid_ = Dune::stackobject_to_shared_ptr(grid); - Parameters::get(name + "->mesh", gridName_); + adoptGrid(Dune::stackobject_to_shared_ptr(grid)); } /// \brief Constructor taking a grid reference and a basis reference. @@ -98,8 +97,7 @@ namespace AMDiS ProblemStat(std::string const& name, Grid& grid, GlobalBasis& globalBasis) : ProblemStat(name, grid) { - globalBasis_ = Dune::stackobject_to_shared_ptr(globalBasis); - initGlobalBasis(*globalBasis_); + adoptGlobalBasis(Dune::stackobject_to_shared_ptr(globalBasis)); } @@ -109,50 +107,139 @@ namespace AMDiS * Parameters read in initialize() for problem with name 'PROB' * MESH[0]->global refinements: nr of initial global refinements **/ - void initialize(Flag initFlag, - Self* adoptProblem = nullptr, - Flag adoptFlag = INIT_NOTHING); + void initialize(Flag initFlag, Self* adoptProblem = nullptr, Flag adoptFlag = INIT_NOTHING); - /// Adds an operator to \ref A. + /// Add an operator to \ref A. /** @{ */ + /// Operator evaluated on the whole element + /** + * Adds an operator to the list of element operators to be assembled in + * quadrature points inside the element. + * + * \param op A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator + * \param row TreePath identifying the sub-basis in the global basis tree + * corresponding to the row basis. \see treepath() + * \param col TreePath identifying the sub-basis in the global basis tree + * corresponding to the column basis. \see treepath() + * + * Example: + * ``` + * auto op = makeOperator(tag::test_trial{}, 1.0/tau); + * prob.addMatrixOperator(op, _0, _0); + * ``` + **/ template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath> void addMatrixOperator(Operator const& op, RowTreePath row = {}, ColTreePath col = {}) { systemMatrix_->addOperator(tag::element_operator<Element>{}, op, row, col); } - // operator evaluated on the boundary + /// Operator evaluated on the boundary of the domain with boundary index `b` + /** + * Adds an operator to the list of boundary operators to be assembled in + * quadrature points on the boundary intersections. + * + * \param b Boundary indentifier where to assemble this operator. Can be + * constructed from an integer. \see BoundaryType + * \param op A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator + * \param row TreePath identifying the sub-basis in the global basis tree + * corresponding to the row basis. \see treepath() + * \param col TreePath identifying the sub-basis in the global basis tree + * corresponding to the column basis. \see treepath() + * + * Example: + * ``` + * auto op = makeOperator(tag::test_trial{}, alpha); + * prob.addMatrixOperator(BoundaryType{1}, op, _0, _0); + * ``` + **/ template <class Operator, class RowTreePath = RootTreePath, class ColTreePath = RootTreePath> void addMatrixOperator(BoundaryType b, Operator const& op, RowTreePath row = {}, ColTreePath col = {}) { - systemMatrix_->addOperator(tag::boundary_operator<typename GridView::Intersection>{b}, op, row, col); + using I = typename GridView::Intersection; + systemMatrix_->addOperator(tag::boundary_operator<I>{b}, op, row, col); } /** @} */ - /// Adds an operator to \ref rhs. + /// Add an operator to \ref rhs. /** @{ */ + /// Operator evaluated on the whole element + /** + * Adds an operator to the list of element operators to be assembled in + * quadrature points inside the element. + * + * \param op A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator + * \param path TreePath identifying the sub-basis in the global basis tree + * corresponding to the row basis. \see treepath() + * + * Example: + * ``` + * auto op = makeOperator(tag::test{}, probInstat.getOldSolution(0) / tau); + * prob.addVectorOperator(op, _0); + * ``` + **/ template <class Operator, class TreePath = RootTreePath> void addVectorOperator(Operator const& op, TreePath path = {}) { rhs_->addOperator(tag::element_operator<Element>{}, op, path); } - // operator evaluated on the boundary + /// Operator evaluated on the boundary of the domain with boundary index `b` + /** + * Adds an operator to the list of boundary operators to be assembled in + * quadrature points on the boundary intersections. + * + * \param b Boundary indentifier where to assemble this operator. Can be + * constructed from an integer. \see BoundaryType + * \param op A (pre-) local operator, \see LocalOperator, \see GridFunctionOperator + * \param path TreePath identifying the sub-basis in the global basis tree + * corresponding to the row basis. \see treepath() + * + * Example: + * ``` + * auto op = makeOperator(tag::test{}, [g](auto const& x) { return g(x); }); + * prob.addVectorOperator(BoundaryType{1}, op, _0); + * ``` + **/ template <class Operator, class TreePath = RootTreePath> void addVectorOperator(BoundaryType b, Operator const& op, TreePath path = {}) { - rhs_->addOperator(tag::boundary_operator<typename GridView::Intersection>{b}, op, path); + using I = typename GridView::Intersection; + rhs_->addOperator(tag::boundary_operator<I>{b}, op, path); } /** @} */ - /// Adds a Dirichlet boundary condition + /// Add boundary conditions to the system + /** @{ */ + /// Dirichlet boundary condition + /** + * Enforce Dirichlet boundary values for the solution vector on boundary + * regions identified by the predicate. + * + * \param predicate Functor `bool(WorldVector)` returning true for all + * DOFs on the boundary that should be assigned a value. + * \param row TreePath identifying the sub-basis in the global basis tree + * corresponding to the row basis. \see treepath() + * \param col TreePath identifying the sub-basis in the global basis tree + * corresponding to the column basis. \see treepath() + * \param values Functor `Range(WorldVector)` or any \ref GridFunction + * that is evaluated in the DOFs identified by the predicate. + * + * Example: + * ``` + * prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8; }, 0, 0, + * [](auto const& x) { return 0.0; }); + * ``` + **/ template <class Predicate, class RowTreePath, class ColTreePath, class Values> void addDirichletBC(Predicate const& predicate, RowTreePath row, ColTreePath col, Values const& values); + /** @} */ + public: @@ -167,28 +254,54 @@ namespace AMDiS bool asmMatrix = true, bool asmVector = true) override; + /// \brief Assemble the linear system by calling \ref buildAfterAdapt with + /// `asmMatrix` and `asmVector` set to true. + void assemble(AdaptInfo& adaptInfo) + { + buildAfterAdapt(adaptInfo, Flag{0}, true, true); + } + /// Writes output files. void writeFiles(AdaptInfo& adaptInfo, bool force = false); public: // get-methods + /// Implementation of \ref ProblemStatBase::name + virtual std::string const& name() const override { return name_; } + + + /// Return a reference to the grid, \ref grid + Grid& grid() { return *grid_; } + Grid const& grid() const { return *grid_; } + + /// Return the gridView of the leaf-level + GridView const& gridView() { return globalBasis_->gridView(); } + + /// Return the \ref globalBasis_ + GlobalBasis& globalBasis() { return *globalBasis_; } + GlobalBasis const& globalBasis() const { return *globalBasis_; } + + /// Return a reference to the linear solver, \ref linearSolver + LinearSolverType& solver() { return *linearSolver_; } + LinearSolverType const& solver() const { return *linearSolver_; } + /// Returns a reference to system-matrix, \ref systemMatrix_ - SystemMatrix& getSystemMatrix() { return *systemMatrix_; } - SystemMatrix const& getSystemMatrix() const { return *systemMatrix_; } + SystemMatrix& systemMatrix() { return *systemMatrix_; } + SystemMatrix const& systemMatrix() const { return *systemMatrix_; } /// Returns a reference to the solution vector, \ref solution_ - SystemVector& getSolutionVector() { return *solution_; } - SystemVector const& getSolutionVector() const { return *solution_; } + SystemVector& solutionVector() { return *solution_; } + SystemVector const& solutionVector() const { return *solution_; } /// Return a reference to the rhs system-vector, \ref rhs - SystemVector& getRhsVector() { return *rhs_; } - SystemVector const& getRhsVector() const { return *rhs_; } + SystemVector& rhsVector() { return *rhs_; } + SystemVector const& rhsVector() const { return *rhs_; } /// Return a mutable view to a solution component template <class TreePath = RootTreePath> - auto getSolution(TreePath const& path = {}) + auto solution(TreePath path = {}) { auto&& tp = makeTreePath(path); return makeDOFVectorView(*solution_, tp); @@ -196,16 +309,14 @@ namespace AMDiS /// Return a const view to a solution component template <class TreePath = RootTreePath> - auto getSolution(TreePath const& path = {}) const + auto solution(TreePath path = {}) const { auto&& tp = makeTreePath(path); return makeDiscreteFunction(*solution_, tp); } - /// Return a reference to the linear solver, \ref linearSolver - LinearSolverType& getSolver() { return *linearSolver_; } - LinearSolverType const& getSolver() const { return *linearSolver_; } + public: // set-methods /// Set a new linear solver for the problem void setSolver(std::shared_ptr<LinearSolverType> const& solver) @@ -213,21 +324,29 @@ namespace AMDiS linearSolver_ = solver; } - /// Return a reference to the grid, \ref grid - Grid& grid() { return *grid_; } - Grid const& grid() const { return *grid_; } + void setSolver(LinearSolverType& solver) + { + setSolver(Dune::stackobject_to_shared_ptr(solver)); + } + /// Set the grid. Stores pointer and initializes feSpaces /// matrices and vectors, as well as markers and file-writers. void setGrid(std::shared_ptr<Grid> const& grid) { - grid_ = grid; + adoptGrid(grid); createGlobalBasis(); createMatricesAndVectors(); createMarker(); createFileWriter(); } + void setGrid(Grid& grid) + { + setGrid(Dune::stackobject_to_shared_ptr(grid)); + } + + void addMarker(std::shared_ptr<Marker<Grid>> const& marker) { marker_.push_back(marker); @@ -240,112 +359,46 @@ namespace AMDiS addMarker(Dune::stackobject_to_shared_ptr(marker)); } - /// Return the gridView of the leaf-level - GridView const& gridView() { return globalBasis_->gridView(); } - - /// Return the \ref globalBasis_ - GlobalBasis const& globalBasis() const { return *globalBasis_; } - - /// Return the \ref globalBasis_ - GlobalBasis& globalBasis() { return *globalBasis_; } - - - /// Implementation of \ref ProblemStatBase::getName - virtual std::string const& getName() const override - { - return name_; - } protected: // initialization methods - template <class T, class GV> - using HasCreate = decltype(T::create(std::declval<GV>())); + void createGlobalBasis(); + void createGrid(); + void createMatricesAndVectors(); + void createSolver(); + void createMarker(); + void createFileWriter(); - void createGlobalBasis() + void adoptGlobalBasis(std::shared_ptr<GlobalBasis> const& globalBasis) { - createGlobalBasis(Dune::Std::is_detected<HasCreate,Traits,GridView>{}); + globalBasis_ = globalBasis; initGlobalBasis(*globalBasis_); } - void createGlobalBasis(std::true_type) - { - assert( bool(grid_) ); - static_assert(std::is_same<GridView, typename Grid::LeafGridView>::value, ""); - globalBasis_ = std::make_shared<GlobalBasis>(Traits::create(grid_->leafGridView())); - } - - void createGlobalBasis(std::false_type) - { - error_exit("Cannot create GlobalBasis from type. Pass a BasisCreator instead!"); - } - - void createGrid() + void adoptGrid(std::shared_ptr<Grid> const& grid) { + grid_ = grid; Parameters::get(name_ + "->mesh", gridName_); - grid_ = MeshCreator<Grid>::create(gridName_); - - msg("Create grid:"); - msg("#elements = {}" , grid_->size(0)); - msg("#faces/edges = {}", grid_->size(1)); - msg("#vertices = {}" , grid_->size(dim)); - msg(""); } - void initGlobalBasis(GlobalBasis const& globalBasis) - { - constraints_.init(*globalBasis_, *globalBasis_); - } - - void createMatricesAndVectors() - { - systemMatrix_ = std::make_shared<SystemMatrix>(*globalBasis_, *globalBasis_); - solution_ = std::make_shared<SystemVector>(*globalBasis_); - rhs_ = std::make_shared<SystemVector>(*globalBasis_); - - auto localView = globalBasis_->localView(); - AMDiS::forEachNode_(localView.tree(), [&,this](auto const& node, auto treePath) - { - std::string i = to_string(treePath); - estimates_[i].resize(globalBasis_->gridView().indexSet().size(0)); - for (std::size_t j = 0; j < estimates_[i].size(); j++) - estimates_[i][j] = 0.0; // TODO: Remove when estimate() is implemented - }); - } - - void createSolver() - { - std::string solverName = "cg"; - Parameters::get(name_ + "->solver->name", solverName); - - auto solverCreator - = named(CreatorMap<LinearSolverType>::getCreator(solverName, name_ + "->solver->name")); - - linearSolver_ = solverCreator->create(name_ + "->solver"); - } + private: - void createMarker(); + void createGlobalBasisImpl(std::true_type); + void createGlobalBasisImpl(std::false_type); - void createFileWriter(); + void initGlobalBasis(GlobalBasis const& globalBasis); public: // implementation of iteration interface methods - /** - * \brief Determines the execution order of the single adaption steps. - * - * If adapt is true, mesh adaption will be performed. This allows to avoid - * mesh adaption, e.g. in timestep adaption loops of timestep adaptive - * strategies. - * - * Implementation of \ref StandardProblemIteration::oneIteration. - **/ + /// Implementation of \ref StandardProblemIteration::oneIteration. virtual Flag oneIteration(AdaptInfo& adaptInfo, Flag toDo = FULL_ITERATION) override { return StandardProblemIteration::oneIteration(adaptInfo, toDo); } /// Implementation of \ref ProblemStatBase::estimate. - virtual void estimate(AdaptInfo& adaptInfo) override { /* Does nothing. */ } + virtual void estimate(AdaptInfo& adaptInfo) override { /* do nothing. */ } /// Implementation of \ref ProblemStatBase::refineMesh. virtual Flag adaptGrid(AdaptInfo& adaptInfo) override; @@ -359,10 +412,7 @@ namespace AMDiS std::string name_; /// Grid of this problem. - std::shared_ptr<Grid> grid_; // TODO: generalize to multi-mesh problems - - /// Number of grids - int nGrids = 1; + std::shared_ptr<Grid> grid_; /// Name of the grid std::string gridName_ = "mesh"; @@ -374,7 +424,7 @@ namespace AMDiS std::list<std::shared_ptr<FileWriterInterface>> filewriter_; /// Pointer to the adaptation markers - std::list<std::shared_ptr<Marker<Grid> > > marker_; + std::list<std::shared_ptr<Marker<Grid>>> marker_; /// Pointer to the estimators for this problem // std::vector<Estimator*> estimator; @@ -382,26 +432,27 @@ namespace AMDiS /// An object of the linearSolver Interface std::shared_ptr<LinearSolverType> linearSolver_; - /// A block-matrix that is filled during assembling + /// Matrix that is filled during assembling std::shared_ptr<SystemMatrix> systemMatrix_; - /// A block-vector with the solution components + /// Vector with the solution components std::shared_ptr<SystemVector> solution_; - /// A vector with the local element error estimates - /// for each node in the basis tree, indexed by [to_string(treePath)][element index] - std::map<std::string, std::vector<double> > estimates_; - - /// A block-vector (load-vector) corresponding to the right.hand side + /// Vector (load-vector) corresponding to the right-hand side /// of the equation, filled during assembling std::shared_ptr<SystemVector> rhs_; + /// A vector with the local element error estimates + /// for each node in the basis tree, indexed by [to_string(treePath)][element index] + std::map<std::string, std::vector<double>> estimates_; + private: // some internal data-structures Constraints<GlobalBasis, GlobalBasis> constraints_; }; + #if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION // Deduction rule template <class Grid, class GlobalBasis> @@ -417,13 +468,6 @@ namespace AMDiS return {name, grid, globalBasis}; } - -#ifndef AMDIS_NO_EXTERN_PROBLEMSTAT - // extern template class ProblemStat<YaspGridBasis<2,1>>; - // extern template class ProblemStat<YaspGridBasis<2,2>>; - // extern template class ProblemStat<YaspGridBasis<2,1,2>>; -#endif - } // end namespace AMDiS #include "ProblemStat.inc.hpp" diff --git a/src/amdis/ProblemStat.inc.hpp b/src/amdis/ProblemStat.inc.hpp index e11503babf91dd5751e0f93cd95e0f260b65a75d..37a45b8f98120f9c207ce29b5e560fb3522609a8 100644 --- a/src/amdis/ProblemStat.inc.hpp +++ b/src/amdis/ProblemStat.inc.hpp @@ -36,7 +36,7 @@ void ProblemStat<Traits>::initialize( (adoptFlag.isSet(INIT_MESH) || adoptFlag.isSet(INIT_SYSTEM) || adoptFlag.isSet(INIT_FE_SPACE))) { - grid_ = adoptProblem->grid_; + adoptGrid(adoptProblem->grid_); } } @@ -62,8 +62,7 @@ void ProblemStat<Traits>::initialize( if (adoptProblem && (adoptFlag.isSet(INIT_FE_SPACE) || adoptFlag.isSet(INIT_SYSTEM))) { - globalBasis_ = adoptProblem->globalBasis_; - initGlobalBasis(*globalBasis_); + adoptGlobalBasis(adoptProblem->globalBasis_); } } @@ -76,10 +75,10 @@ void ProblemStat<Traits>::initialize( createMatricesAndVectors(); if (adoptProblem && adoptFlag.isSet(INIT_SYSTEM)) { + systemMatrix_ = adoptProblem->systemMatrix_; solution_ = adoptProblem->solution_; - estimates_ = adoptProblem->estimates_; rhs_ = adoptProblem->rhs_; - systemMatrix_ = adoptProblem->systemMatrix_; + estimates_ = adoptProblem->estimates_; } @@ -117,6 +116,86 @@ void ProblemStat<Traits>::initialize( } +template <class Traits> +void ProblemStat<Traits>::createGrid() +{ + Parameters::get(name_ + "->mesh", gridName_); + grid_ = MeshCreator<Grid>::create(gridName_); + + msg("Create grid:"); + msg("#elements = {}" , grid_->size(0)); + msg("#faces/edges = {}", grid_->size(1)); + msg("#vertices = {}" , grid_->size(dim)); + msg(""); +} + + +template <class T, class GV> +using HasCreate = decltype(T::create(std::declval<GV>())); + + +template <class Traits> +void ProblemStat<Traits>::createGlobalBasis() +{ + createGlobalBasisImpl(Dune::Std::is_detected<HasCreate,Traits,GridView>{}); + initGlobalBasis(*globalBasis_); +} + + +template <class Traits> +void ProblemStat<Traits>::createGlobalBasisImpl(std::true_type) +{ + assert( bool(grid_) ); + static_assert(std::is_same<GridView, typename Grid::LeafGridView>::value, ""); + globalBasis_ = std::make_shared<GlobalBasis>(Traits::create(grid_->leafGridView())); +} + + +template <class Traits> +void ProblemStat<Traits>::createGlobalBasisImpl(std::false_type) +{ + error_exit("Cannot create GlobalBasis from type. Pass a BasisCreator instead!"); +} + + +template <class Traits> +void ProblemStat<Traits>::initGlobalBasis(GlobalBasis const& globalBasis) +{ + constraints_.init(*globalBasis_, *globalBasis_); +} + + +template <class Traits> +void ProblemStat<Traits>::createMatricesAndVectors() +{ + systemMatrix_ = std::make_shared<SystemMatrix>(*globalBasis_, *globalBasis_); + solution_ = std::make_shared<SystemVector>(*globalBasis_); + rhs_ = std::make_shared<SystemVector>(*globalBasis_); + + auto localView = globalBasis_->localView(); + AMDiS::forEachNode_(localView.tree(), [&,this](auto const& node, auto treePath) + { + std::string i = to_string(treePath); + estimates_[i].resize(globalBasis_->gridView().indexSet().size(0)); + for (std::size_t j = 0; j < estimates_[i].size(); j++) + estimates_[i][j] = 0.0; // TODO: Remove when estimate() is implemented + }); +} + + +template <class Traits> +void ProblemStat<Traits>::createSolver() +{ + std::string solverName = "default"; + Parameters::get(name_ + "->solver->name", solverName); + + auto solverCreator + = named(CreatorMap<LinearSolverType>::getCreator(solverName, name_ + "->solver->name")); + + linearSolver_ = solverCreator->create(name_ + "->solver"); +} + + template <class Traits> void ProblemStat<Traits>::createMarker() { @@ -130,7 +209,8 @@ void ProblemStat<Traits>::createMarker() return; std::string tp = to_string(treePath); - auto newMarker = EstimatorMarker<Grid>::createMarker(componentName, tp, estimates_[tp], grid_); + auto newMarker + = EstimatorMarker<Grid>::createMarker(componentName, tp, estimates_[tp], grid_); if (newMarker) marker_.push_back(std::move(newMarker)); @@ -155,7 +235,7 @@ void ProblemStat<Traits>::createFileWriter() if (!Parameters::get<std::string>(componentName + "->filename")) return; - auto writer = makeFileWriterPtr(componentName, this->getSolution(treePath)); + auto writer = makeFileWriterPtr(componentName, this->solution(treePath)); filewriter_.push_back(std::move(writer)); }); } @@ -199,26 +279,26 @@ solve(AdaptInfo& adaptInfo, bool createMatrixData, bool storeMatrixData) linearSolver_->solve(systemMatrix_->matrix(), solution_->vector(), rhs_->vector(), solverInfo); - if (solverInfo.getInfo() > 0) { + if (solverInfo.info() > 0) { msg("solution of discrete system needed {} seconds", t.elapsed()); - if (solverInfo.getAbsResidual() >= 0.0) { - if (solverInfo.getRelResidual() >= 0.0) + if (solverInfo.absResidual() >= 0.0) { + if (solverInfo.relResidual() >= 0.0) msg("Residual norm: ||b-Ax|| = {}, ||b-Ax||/||b|| = {}", - solverInfo.getAbsResidual(), solverInfo.getRelResidual()); + solverInfo.absResidual(), solverInfo.relResidual()); else - msg("Residual norm: ||b-Ax|| = {}", solverInfo.getAbsResidual()); + msg("Residual norm: ||b-Ax|| = {}", solverInfo.absResidual()); } } if (solverInfo.doBreak()) { std::stringstream tol_str; - if (solverInfo.getAbsTolerance() > 0 - && solverInfo.getAbsResidual() > solverInfo.getAbsTolerance()) - tol_str << "absTol = " << solverInfo.getAbsTolerance() << " "; - if (solverInfo.getRelTolerance() > 0 - && solverInfo.getRelResidual() > solverInfo.getRelTolerance()) - tol_str << "relTol = " << solverInfo.getRelTolerance() << " "; + if (solverInfo.absTolerance() > 0 + && solverInfo.absResidual() > solverInfo.absTolerance()) + tol_str << "absTol = " << solverInfo.absTolerance() << " "; + if (solverInfo.relTolerance() > 0 + && solverInfo.relResidual() > solverInfo.relTolerance()) + tol_str << "relTol = " << solverInfo.relTolerance() << " "; error_exit("Tolerance {} could not be reached!", tol_str.str()); } } diff --git a/src/amdis/ProblemStatBase.hpp b/src/amdis/ProblemStatBase.hpp index 9ba0b0b6987ccafd5db725f2c6795ddba45472e5..7c44f1bf7b0b8627f9af91d15af780a5b2bcf70b 100644 --- a/src/amdis/ProblemStatBase.hpp +++ b/src/amdis/ProblemStatBase.hpp @@ -104,7 +104,7 @@ namespace AMDiS virtual void estimate(AdaptInfo& adaptInfo) = 0; /// Returns the name of the problem. - virtual std::string const& getName() const = 0; + virtual std::string const& name() const = 0; }; } // end namespace AMDiS diff --git a/src/amdis/ProblemTimeInterface.hpp b/src/amdis/ProblemTimeInterface.hpp index 59af7c3928cd41da59ba693a00cc14162396a67a..b5dafd1b3ca2494c9ba6f28413016e12ad3afd5f 100644 --- a/src/amdis/ProblemTimeInterface.hpp +++ b/src/amdis/ProblemTimeInterface.hpp @@ -17,7 +17,7 @@ namespace AMDiS virtual ~ProblemTimeInterface() = default; /// Called at the beginning of the adaption loop before any other call - virtual void initTimeInterface() {}; + virtual void initTimeInterface() { /* do nothing */ }; /// Executes all needed operations when the simulation time changes. virtual void setTime(AdaptInfo& adaptInfo) = 0; diff --git a/src/amdis/StandardProblemIteration.cpp b/src/amdis/StandardProblemIteration.cpp index 1c7047d6e53a7fbf6f6dbcf7bda354b6a200aedf..0c6b5c661854a342a33779710ce7d3ef501b9fd2 100644 --- a/src/amdis/StandardProblemIteration.cpp +++ b/src/amdis/StandardProblemIteration.cpp @@ -6,17 +6,17 @@ namespace AMDiS { - ProblemStatBase& StandardProblemIteration::getProblem(int number) + ProblemStatBase& StandardProblemIteration::problem(int number) { AMDIS_FUNCNAME_DBG("StandardProblemIteration::getProblem"); test_exit_dbg(number == 0, "Problem number out of range!"); return problem_; } - ProblemStatBase& StandardProblemIteration::getProblem(std::string const& name) + ProblemStatBase& StandardProblemIteration::problem(std::string const& name) { AMDIS_FUNCNAME_DBG("StandardProblemIteration::getProblem"); - test_exit_dbg(name == problem_.getName(), "Problem name does not match!\n"); + test_exit_dbg(name == problem_.name(), "Problem name does not match!\n"); return problem_; } @@ -25,7 +25,7 @@ namespace AMDiS AMDIS_FUNCNAME("StandardProblemIteration::beginIteration()"); msg(""); - msg("begin of iteration number: {}", (adaptInfo.getSpaceIteration() + 1)); + msg("begin of iteration number: {}", (adaptInfo.spaceIteration() + 1)); msg("============================="); } @@ -52,7 +52,7 @@ namespace AMDiS AMDIS_FUNCNAME("StandardProblemIteration::endIteration()"); msg(""); - msg("end of iteration number: {}", (adaptInfo.getSpaceIteration() + 1)); + msg("end of iteration number: {}", (adaptInfo.spaceIteration() + 1)); msg("============================="); } @@ -77,9 +77,9 @@ namespace AMDiS } /// Returns the name of the problem. - std::string const& StandardProblemIteration::getName() const + std::string const& StandardProblemIteration::name() const { - return problem_.getName(); + return problem_.name(); } } // end namespace AMDiS diff --git a/src/amdis/StandardProblemIteration.hpp b/src/amdis/StandardProblemIteration.hpp index d466b3fd9192fada921b4f86c0c7646e5b81cd48..90c9f584bce14743f43d506d554a2f045e35c90e 100644 --- a/src/amdis/StandardProblemIteration.hpp +++ b/src/amdis/StandardProblemIteration.hpp @@ -31,18 +31,18 @@ namespace AMDiS virtual void endIteration(AdaptInfo& adaptInfo) override; /// Returns the name of the problem. - virtual std::string const& getName() const override; + virtual std::string const& name() const override; - virtual int getNumProblems() const override + virtual int numProblems() const override { return 1; } /// Return the managed ProblemStat \ref problem, by number - virtual ProblemStatBase& getProblem(int number = 0) override; + virtual ProblemStatBase& problem(int number = 0) override; /// Return the managed ProblemStat \ref problem, by name - virtual ProblemStatBase& getProblem(std::string const& name) override; + virtual ProblemStatBase& problem(std::string const& name) override; protected: /// Nested assemblage and mesh adaption. diff --git a/src/amdis/gridfunctions/AnalyticGridFunction.hpp b/src/amdis/gridfunctions/AnalyticGridFunction.hpp index c359b40738e7ae7c1f40973d6eb26b284fb5f7d4..342b2513a3e3de42f960c870fdd22a99b4a5c35b 100644 --- a/src/amdis/gridfunctions/AnalyticGridFunction.hpp +++ b/src/amdis/gridfunctions/AnalyticGridFunction.hpp @@ -139,29 +139,31 @@ namespace AMDiS using LocalFunction = AnalyticLocalFunction<Range(LocalDomain), Element, Function>; public: - /// \brief Constructor. Stores the function `fct` and creates an `EntitySet`. + /// Constructor. Stores the function `fct` and creates an `EntitySet`. AnalyticGridFunction(Function const& fct, GridView const& gridView) : fct_{fct} , entitySet_{gridView} {} - /// \brief Return the evaluated functor at global coordinates + /// Return the evaluated functor at global coordinates Range operator()(Domain const& x) const { return fct_(x); } - /// \brief Return the LocalFunction of the AnalyticGridFunction. + /// Return the \ref AnalyticLocalFunction of the AnalyticGridFunction. LocalFunction localFunction() const { return {fct_}; } + /// Returns \ref entitySet_ EntitySet const& entitySet() const { return entitySet_; } + /// Returns \ref fct_ Function const& fct() const { return fct_; } private: diff --git a/src/amdis/gridfunctions/ConstantGridFunction.hpp b/src/amdis/gridfunctions/ConstantGridFunction.hpp index 0ac6ba3282c878bd37491170eae6ab40a4cf1f63..8d76d2fa7ccded9ae07c371d011edea6b3cdd3fc 100644 --- a/src/amdis/gridfunctions/ConstantGridFunction.hpp +++ b/src/amdis/gridfunctions/ConstantGridFunction.hpp @@ -148,7 +148,6 @@ namespace AMDiS } - namespace Concepts { /** \addtogroup Concepts diff --git a/src/amdis/gridfunctions/DOFVectorView.hpp b/src/amdis/gridfunctions/DOFVectorView.hpp index 734c5e1ebc35037d2a13f252b35f9dae9b972e04..50f12531f34efc840f18f2bc4b696cdbcafc719a 100644 --- a/src/amdis/gridfunctions/DOFVectorView.hpp +++ b/src/amdis/gridfunctions/DOFVectorView.hpp @@ -34,7 +34,7 @@ namespace AMDiS /** * **Example:** * ``` - * auto v = makeDOFVectorView(prob.getSolutionVector(),0); + * auto v = makeDOFVectorView(prob.solutionVector(),0); * v.interpolate([](auto const& x) { return x[0]; }); * ``` **/ @@ -54,7 +54,7 @@ namespace AMDiS /** * **Example:** * ``` - * auto v = makeDOFVectorView(prob.getSolutionVector(),0); + * auto v = makeDOFVectorView(prob.solutionVector(),0); * v.interpolate(v + [](auto const& x) { return x[0]; }); * ``` * Allows to have a reference to the DOFVector in the expression, e.g. as diff --git a/src/amdis/gridfunctions/DerivativeGridFunction.hpp b/src/amdis/gridfunctions/DerivativeGridFunction.hpp index 3b766f7c48a36cf8307084d3ecf2f080db69e206..f1b2abeb9e47cfb8ac379ed5e0b8e9f10d36714a 100644 --- a/src/amdis/gridfunctions/DerivativeGridFunction.hpp +++ b/src/amdis/gridfunctions/DerivativeGridFunction.hpp @@ -9,14 +9,10 @@ namespace AMDiS { - /** - * \addtogroup GridFunctions - * @{ - **/ - /// \class DerivativeGridFunction /// \brief A Gridfunction that returns the derivative when calling localFunction. /** + * \ingroup GridFunctions * Wrapps the GridFunction so that \ref localFunction returns a LocalFunction * representing the derivative of the LocalFunction of the GridFunction. * @@ -146,12 +142,13 @@ namespace AMDiS /// \brief Generator function for DerivativeGridFunction expressions. /// \relates DerivativeGridFunction /** + * \ingroup GridFunctions * Generates a Gridfunction representing the derivative of a GridFunction. * See \ref DerivativeGridFunction. * * **Examples:** - * - `gradientAtQP(prob.getSolution(_0))` - * - `gradientAtQP(X(0) + X(1) + prob.getSolution(_0))` + * - `gradientAtQP(prob.solution(_0))` + * - `gradientAtQP(X(0) + X(1) + prob.solution(_0))` **/ template <class Expr> auto gradientAtQP(Expr const& expr) diff --git a/src/amdis/gridfunctions/FunctorGridFunction.hpp b/src/amdis/gridfunctions/FunctorGridFunction.hpp index 6ee2c1eb5dab061bd66d76b2832ec4d26bd4a4a1..88a6b27d689d807eb1f677fd2f4f193eeaafa17f 100644 --- a/src/amdis/gridfunctions/FunctorGridFunction.hpp +++ b/src/amdis/gridfunctions/FunctorGridFunction.hpp @@ -324,6 +324,4 @@ namespace AMDiS return FunctorPreGridFunction<Functor, std::decay_t<PreGridFcts>...>{f, std::forward<PreGridFcts>(gridFcts)...}; } - /** @} **/ - } // end namespace AMDiS diff --git a/src/amdis/gridfunctions/OperationsGridFunction.hpp b/src/amdis/gridfunctions/OperationsGridFunction.hpp index 3e578d40c83d5bb8c80852de3605ad3ece00d508..8aaad855d619c3f7a0d143960f78925b28ea7c26 100644 --- a/src/amdis/gridfunctions/OperationsGridFunction.hpp +++ b/src/amdis/gridfunctions/OperationsGridFunction.hpp @@ -133,7 +133,7 @@ namespace AMDiS // unary vector operations // @{ - /// \brief Applies a sum-functor to a vector-valued GridFunction. \relates FunctorGridFunction + /// \brief Applies a \ref sum() functor to a vector-valued GridFunction. \relates FunctorGridFunction template <class Vec, REQUIRES(Concepts::AnyGridFunction<Vec>)> auto sum(Vec&& vec) @@ -149,7 +149,7 @@ namespace AMDiS return invokeAtQP(Operation::UnaryDot{}, std::forward<Vec>(vec)); } - /// \brief Applies a one_norm-functor to a vector-valued GridFunction. \relates FunctorGridFunction + /// \brief Applies a \ref one_norm() functor to a vector-valued GridFunction. \relates FunctorGridFunction template <class Vec, REQUIRES(Concepts::AnyGridFunction<Vec>)> auto one_norm(Vec&& vec) @@ -165,7 +165,7 @@ namespace AMDiS return invokeAtQP(Operation::TwoNorm{}, std::forward<Vec>(vec)); } - /// \brief Applies a p_norm-functor to a vector-valued GridFunction. \relates FunctorGridFunction + /// \brief Applies a \ref p_norm() functor to a vector-valued GridFunction. \relates FunctorGridFunction template <int p, class Vec, REQUIRES(Concepts::AnyGridFunction<Vec>)> auto p_norm(Vec&& vec) @@ -173,7 +173,7 @@ namespace AMDiS return invokeAtQP([](auto const& v) { return p_norm<p>(v); }, std::forward<Vec>(vec)); } - /// \brief Applies a infty_norm-functor to a vector-valued GridFunction. \relates FunctorGridFunction + /// \brief Applies a \ref infty_norm() functor to a vector-valued GridFunction. \relates FunctorGridFunction template <class Vec, REQUIRES(Concepts::AnyGridFunction<Vec>)> auto infty_norm(Vec&& vec) @@ -181,7 +181,7 @@ namespace AMDiS return invokeAtQP([](auto const& v) { return infty_norm(v); }, std::forward<Vec>(vec)); } - /// \brief Applies a trans-functor to a matrix-valued GridFunction.\relates FunctorGridFunction + /// \brief Applies \ref Operation::Trans to a matrix-valued GridFunction.\relates FunctorGridFunction template <class Mat, REQUIRES(Concepts::AnyGridFunction<Mat>)> auto trans(Mat&& mat) @@ -214,7 +214,7 @@ namespace AMDiS std::forward<Lhs>(lhs), std::forward<Rhs>(rhs)); } - /// \brief Applies an outer-functor to two vector-valued GridFunctions. \relates FunctorGridFunction + /// \brief Applies an \ref outer() functor to two vector-valued GridFunctions. \relates FunctorGridFunction template <class Lhs, class Rhs, REQUIRES(Concepts::AnyGridFunction<Lhs,Rhs>)> auto outer(Lhs&& lhs, Rhs&& rhs) diff --git a/src/amdis/io/FileWriterInterface.hpp b/src/amdis/io/FileWriterInterface.hpp index 79b8b2eb1677b62e27533b904cadc0a584123a7d..eb2591cc906988e778592db19ecf2ddb61066708 100644 --- a/src/amdis/io/FileWriterInterface.hpp +++ b/src/amdis/io/FileWriterInterface.hpp @@ -31,7 +31,7 @@ namespace AMDiS virtual void writeFiles(AdaptInfo& adaptInfo, bool force) = 0; - std::string getFilename() const + std::string filename() const { return filename_; } @@ -41,7 +41,7 @@ namespace AMDiS filename_ = filename; } - std::string const& getName() const + std::string const& name() const { return name_; } diff --git a/src/amdis/linear_algebra/DOFMatrixBase.inc.hpp b/src/amdis/linear_algebra/DOFMatrixBase.inc.hpp index de0c608b0c99df727db0de9915544611bdd196df..58a8bc124f2a9550619a9279507d91dcf48268bc 100644 --- a/src/amdis/linear_algebra/DOFMatrixBase.inc.hpp +++ b/src/amdis/linear_algebra/DOFMatrixBase.inc.hpp @@ -78,14 +78,12 @@ assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView) forEachNode_(rowLocalView.tree(), [&](auto const& rowNode, auto) { forEachNode_(colLocalView.tree(), [&](auto const& colNode, auto) { auto& matOp = operators_[rowNode][colNode]; - if (matOp.doAssemble() && !matOp.empty()) { + if (matOp) { matOp.bind(element, geometry); - auto matAssembler = [&](auto const& context, auto& operators) { + assembleOperators(gv, element, matOp, [&](auto const& context, auto& operators) { for (auto scaled : operators) scaled.op->assemble(context, rowNode, colNode, elementMatrix_); - }; - assembleOperators(gv, element, matOp, matAssembler); - matOp.assembled = true; + }); matOp.unbind(); } }); diff --git a/src/amdis/linear_algebra/DOFVectorBase.inc.hpp b/src/amdis/linear_algebra/DOFVectorBase.inc.hpp index ce54440e5092e910c802d53a9f9e18f07e790574..980be2ddbcd41fcbd59ed14a5913ad61268cbf42 100644 --- a/src/amdis/linear_algebra/DOFVectorBase.inc.hpp +++ b/src/amdis/linear_algebra/DOFVectorBase.inc.hpp @@ -68,16 +68,12 @@ assemble(LocalView const& localView) forEachNode_(localView.tree(), [&](auto const& node, auto) { auto& rhsOp = operators_[node]; - if (rhsOp.doAssemble() && !rhsOp.empty()) { + if (rhsOp) { rhsOp.bind(element, geometry); - - auto vecAssembler = [&](auto const& context, auto& operators) { + assembleOperators(gv, element, rhsOp, [&](auto const& context, auto& operators) { for (auto scaled : operators) scaled.op->assemble(context, node, elementVector_); - }; - - assembleOperators(gv, element, rhsOp, vecAssembler); - rhsOp.assembled = true; + }); rhsOp.unbind(); } }); diff --git a/src/amdis/linear_algebra/LinearSolver.hpp b/src/amdis/linear_algebra/LinearSolver.hpp index 6dcab9574088fcdb710baa383e713ea743b5c308..c1f5661065514eac5b7a0156b3c06f966b9e6cf6 100644 --- a/src/amdis/linear_algebra/LinearSolver.hpp +++ b/src/amdis/linear_algebra/LinearSolver.hpp @@ -61,7 +61,7 @@ namespace AMDiS runner_->init(A); } - if (solverInfo.getInfo() > 0) + if (solverInfo.info() > 0) msg("fill matrix needed {} seconds", t.elapsed()); int error = runner_->solve(A, x, b, solverInfo); diff --git a/src/amdis/linear_algebra/SolverInfo.hpp b/src/amdis/linear_algebra/SolverInfo.hpp index 36be43d0284fe1158666db81c4d9ec44b7b26f42..23e08f11f7c2c7ccfb24c3874dd6eca805d104f9 100644 --- a/src/amdis/linear_algebra/SolverInfo.hpp +++ b/src/amdis/linear_algebra/SolverInfo.hpp @@ -14,82 +14,83 @@ namespace AMDiS /// The constructor reads needed parameters and sets solver \p prefix. /** * Reads parameters for a solver with name 'NAME': - * NAME->absolute tolerance \ref aTol - * NAME->relative tolerance \ref rTol - * NAME->info \ref info - * NAME->break if tolerance not reached \ref breakTolNotReached + * NAME->absolute tolerance \ref aTol_ + * NAME->relative tolerance \ref rTol_ + * NAME->info \ref info_ + * NAME->break if tolerance not reached \ref breakTolNotReached_ **/ - explicit SolverInfo(std::string prefix) + explicit SolverInfo(std::string const& prefix) + : prefix_(prefix) { - Parameters::get(prefix + "->abolute tolerance", aTol); - Parameters::get(prefix + "->relative tolerance", rTol); - Parameters::get(prefix + "->info", info); - Parameters::get(prefix + "->break if tolerance not reached", breakTolNotReached); + Parameters::get(prefix + "->abolute tolerance", aTol_); + Parameters::get(prefix + "->relative tolerance", rTol_); + Parameters::get(prefix + "->info", info_); + Parameters::get(prefix + "->break if tolerance not reached", breakTolNotReached_); } /** \name getting methods * \{ */ - /// Returns \ref aTol - double getAbsTolerance() const + /// Returns \ref aTol_ + double absTolerance() const { - return aTol; + return aTol_; } - /// Returns \ref rTol - double getRelTolerance() const + /// Returns \ref rTol_ + double relTolerance() const { - return rTol; + return rTol_; } /// Returns error code in last run of an iterative solver - int getErrorCode() const + int errorCode() const { - return error; + return error_; } /// Returns info - int getInfo() const + int info() const { - return info; + return info_; } - /// Returns \ref absResidual - double getAbsResidual() const + /// Returns \ref absResidual_ + double absResidual() const { - return absResidual; + return absResidual_; } - /// Returns \ref relResidual - double getRelResidual() const + /// Returns \ref relResidual_ + double relResidual() const { - return relResidual; + return relResidual_; } - /// Returns the initfile \ref prefix - std::string getPrefix() const + /// Returns the initfile \ref prefix_ + std::string const& prefix() const { - return prefix; + return prefix_; } /// Returns \ref createMatrixData bool doCreateMatrixData() const { - return createMatrixData; + return createMatrixData_; } /// Returns \ref storeMatrixData bool doStoreMatrixData() const { - return storeMatrixData; + return storeMatrixData_; } bool doBreak() const { - return breakTolNotReached && ( - (aTol > 1.e-30 && absResidual > aTol) || - (rTol > 1.e-30 && relResidual > rTol) ); + return breakTolNotReached_ && ( + (aTol_ > 1.e-30 && absResidual_ > aTol_) || + (rTol_ > 1.e-30 && relResidual_ > rTol_) ); } /** \} */ @@ -99,88 +100,88 @@ namespace AMDiS * \{ */ - /// Sets \ref aTol + /// Sets \ref aTol_ void setAbsTolerance(double tol) { - aTol = tol; + aTol_ = tol; } - /// Sets \ref rTol + /// Sets \ref rTol_ void setRelTolerance(double tol) { - rTol = tol; + rTol_ = tol; } - /// Sets \ref aTol + /// Sets \ref aTol_ void setAbsResidual(double r) { - absResidual = r; + absResidual_ = r; } - /// Sets \ref rTol + /// Sets \ref rTol_ void setRelResidual(double r) { - relResidual = r; + relResidual_ = r; } - /// Sets \ref info + /// Sets \ref info_ void setInfo(int i) { - info = i; + info_ = i; } - /// Sets \ref error + /// Sets \ref error_ void setError(int e) { - error = e; + error_ = e; } - /// Sets \ref createMatrixData + /// Sets \ref createMatrixData_ void setCreateMatrixData(bool b) { - createMatrixData = b; + createMatrixData_ = b; } - /// Sets \ref storeMatrixData + /// Sets \ref storeMatrixData_ void setStoreMatrixData(bool b) { - storeMatrixData = b; + storeMatrixData_ = b; } /** \} */ private: /// The initfile prefix to read parameters - std::string prefix; + std::string prefix_; /// The abolute tolerance - double aTol = 0; + double aTol_ = 0; /// The relative tolerance - double rTol = 1.e-6; + double rTol_ = 1.e-6; /// Throw an error if tolerance could not be reached - bool breakTolNotReached = false; + bool breakTolNotReached_ = false; /// The solver verbosity level - int info = 0; + int info_ = 0; /// The absolute residual, default (-1): not set - double absResidual = -1.0; + double absResidual_ = -1.0; /// The relative residual, relative to the two_norm(b), default (-1): not set - double relResidual = -1.0; + double relResidual_ = -1.0; /// The error-code, default (-1): not set - int error = -1; + int error_ = -1; /// If true, the matrix will be initialized and the /// corresponding runner of the system receives the /// matrix in the init() method. - bool createMatrixData = true; + bool createMatrixData_ = true; /// If false, the exit() method of the runner will be called. - bool storeMatrixData = false; + bool storeMatrixData_ = false; }; } // end namespace AMDiS diff --git a/src/amdis/linear_algebra/eigen/DirectRunner.hpp b/src/amdis/linear_algebra/eigen/DirectRunner.hpp index a128154f54fa04182141292cf7575a5e421c422a..3f9c40860331c14b7197ad878757cc056dd41113 100644 --- a/src/amdis/linear_algebra/eigen/DirectRunner.hpp +++ b/src/amdis/linear_algebra/eigen/DirectRunner.hpp @@ -23,7 +23,7 @@ namespace AMDiS public: /// Constructor. - DirectRunner(std::string prefix) + DirectRunner(std::string const& prefix) : solver_{} { SolverConfig<LU>::init(prefix, solver_); diff --git a/src/amdis/linear_algebra/eigen/IterativeRunner.hpp b/src/amdis/linear_algebra/eigen/IterativeRunner.hpp index 59bc1cccf6dbd1d6cf65b3833ae8d50f91b4633d..c3cfff3d93b4cfe448a70f1839e3dcc1570894b9 100644 --- a/src/amdis/linear_algebra/eigen/IterativeRunner.hpp +++ b/src/amdis/linear_algebra/eigen/IterativeRunner.hpp @@ -17,7 +17,7 @@ namespace AMDiS using PreconCfg = PreconConfig<typename IterativeSolver::Preconditioner>; public: - IterativeRunner(std::string prefix) + IterativeRunner(std::string const& prefix) : solver_{} { SolverCfg::init(prefix, solver_); @@ -41,7 +41,7 @@ namespace AMDiS virtual int solve(Matrix const& A, VectorX& x, VectorB const& b, SolverInfo& solverInfo) override { - solver_.setTolerance(solverInfo.getRelTolerance()); + solver_.setTolerance(solverInfo.relTolerance()); x = solver_.solveWithGuess(b, x); auto r = VectorB(b); diff --git a/src/amdis/linear_algebra/eigen/SolverCreator.hpp b/src/amdis/linear_algebra/eigen/SolverCreator.hpp index a1ef5afce34e9e1c11a4b2c17c578ed08158c148..cd629242acbce60e352ee91d086c0592c5c5c651 100644 --- a/src/amdis/linear_algebra/eigen/SolverCreator.hpp +++ b/src/amdis/linear_algebra/eigen/SolverCreator.hpp @@ -30,7 +30,7 @@ namespace AMDiS using SolverBase = LinearSolverInterface<Matrix, VectorX, VectorB>; using Scalar = typename Matrix::Scalar; - virtual std::unique_ptr<SolverBase> create(std::string prefix) override + virtual std::unique_ptr<SolverBase> create(std::string const& prefix) override { // get creator string for preconditioner std::string precon = "no"; @@ -161,7 +161,7 @@ namespace AMDiS #endif // default iterative solver - Map::addCreator("default", bicgstab); + Map::addCreator("default", gmres); // default direct solvers #if HAVE_SUITESPARSE_UMFPACK diff --git a/src/amdis/linear_algebra/istl/DirectRunner.hpp b/src/amdis/linear_algebra/istl/DirectRunner.hpp index 187abadfa5c0767810abf3147a4c28d8ef1e0733..7f71fa864f2fd39b338de182ac725aaae0c68f5e 100644 --- a/src/amdis/linear_algebra/istl/DirectRunner.hpp +++ b/src/amdis/linear_algebra/istl/DirectRunner.hpp @@ -24,7 +24,7 @@ namespace AMDiS public: /// Constructor. - DirectRunner(std::string prefix) + DirectRunner(std::string const& prefix) : solverCreator_(prefix) {} diff --git a/src/amdis/linear_algebra/istl/ISTLRunner.hpp b/src/amdis/linear_algebra/istl/ISTLRunner.hpp index 99c4927e876456f1bfffc9c2586467c8ec6fb69e..61b71e2c15dcc2a8a74ecbe84b5340bf7dd0d581 100644 --- a/src/amdis/linear_algebra/istl/ISTLRunner.hpp +++ b/src/amdis/linear_algebra/istl/ISTLRunner.hpp @@ -19,7 +19,7 @@ namespace AMDiS using ISTLPrecon = ISTLPreconInterface<Matrix, VectorX, VectorB>; public: - ISTLRunner(std::string prefix) + ISTLRunner(std::string const& prefix) : solverCreator_(prefix) { initPrecon(prefix); diff --git a/src/amdis/linear_algebra/istl/ISTL_Preconditioner.hpp b/src/amdis/linear_algebra/istl/ISTL_Preconditioner.hpp index 86b4c0fc0319fbc81f94cea5c0c726f71f68b0c8..32f4f9177cb4e06351cf42096e7d56b2e9bbfef8 100644 --- a/src/amdis/linear_algebra/istl/ISTL_Preconditioner.hpp +++ b/src/amdis/linear_algebra/istl/ISTL_Preconditioner.hpp @@ -36,7 +36,7 @@ namespace AMDiS } }; - ISTLPrecon(std::string prefix) + ISTLPrecon(std::string const& prefix) { Parameters::get(prefix + "->relaxation", w_); } diff --git a/src/amdis/linear_algebra/istl/ISTL_Solver.hpp b/src/amdis/linear_algebra/istl/ISTL_Solver.hpp index e32d52fe3c726a57a951e40ad99819cc59d542a6..fa10d194ee30410f6bcf0ae3217b7340d19187d1 100644 --- a/src/amdis/linear_algebra/istl/ISTL_Solver.hpp +++ b/src/amdis/linear_algebra/istl/ISTL_Solver.hpp @@ -18,7 +18,7 @@ namespace AMDiS template <class ISTLSolver, bool specialized> struct ISTLSolverCreator { - ISTLSolverCreator(std::string prefix) + ISTLSolverCreator(std::string const& prefix) { Parameters::get(prefix + "->info", info_); Parameters::get(prefix + "->max iteration", maxIter_); @@ -43,7 +43,7 @@ namespace AMDiS using ISTLSolver = Dune::RestartedGMResSolver<VectorX>; using Super = ISTLSolverCreator<ISTLSolver, false>; - ISTLSolverCreator(std::string prefix) + ISTLSolverCreator(std::string const& prefix) : Super(prefix) { Parameters::get(prefix + "->restart", restart_); @@ -65,7 +65,7 @@ namespace AMDiS using ISTLSolver = Dune::GeneralizedPCGSolver<VectorX>; using Super = ISTLSolverCreator<ISTLSolver, false>; - ISTLSolverCreator(std::string prefix) + ISTLSolverCreator(std::string const& prefix) : Super(prefix) { Parameters::get(prefix + "->restart", restart_); @@ -86,7 +86,7 @@ namespace AMDiS struct ISTLSolverCreator<Dune::UMFPack<Matrix>, true> { using Solver = Dune::UMFPack<Matrix>; - ISTLSolverCreator(std::string prefix) + ISTLSolverCreator(std::string const& prefix) { Parameters::get(prefix + "->info", verbose_); } @@ -106,7 +106,7 @@ namespace AMDiS struct ISTLSolverCreator<Dune::SuperLU<Matrix>, true> { using Solver = Dune::SuperLU<Matrix>; - ISTLSolverCreator(std::string prefix) + ISTLSolverCreator(std::string const& prefix) { int info = 2; Parameters::get(prefix + "->info", info); @@ -185,11 +185,15 @@ namespace AMDiS Map::addCreator("superlu", superlu); #endif + // default direct solvers #if HAVE_SUITESPARSE_UMFPACK Map::addCreator("direct", umfpack); #elif HAVE_SUPERLU Map::addCreator("direct", superlu); #endif + + // default iterative solver + Map::addCreator("default", gmres); } }; diff --git a/src/amdis/linear_algebra/mtl/ITL_Solver.hpp b/src/amdis/linear_algebra/mtl/ITL_Solver.hpp index c8193bc2bfe73e3982f0b1d6b0292ce097947df0..9356f9c1d0acadb76f6bb360c78529ac532b06ec 100644 --- a/src/amdis/linear_algebra/mtl/ITL_Solver.hpp +++ b/src/amdis/linear_algebra/mtl/ITL_Solver.hpp @@ -315,7 +315,8 @@ namespace AMDiS class fgmres_type { public: - explicit fgmres_type(std::string const& name) : restart(30), orthogonalization(GRAM_SCHMIDT) + explicit fgmres_type(std::string const& name) + : restart(30), orthogonalization(GRAM_SCHMIDT) { Parameters::get(name + "->restart", restart); Parameters::get(name + "->orthogonalization", orthogonalization); @@ -455,9 +456,12 @@ namespace AMDiS #ifdef HAVE_UMFPACK auto umfpack = new UmfpackSolverCreator; Map::addCreator("umfpack", umfpack); + + // default direct solvers Map::addCreator("direct", umfpack); #endif + // default iterative solver Map::addCreator("default", gmres); } }; diff --git a/src/amdis/linear_algebra/mtl/KrylovRunner.hpp b/src/amdis/linear_algebra/mtl/KrylovRunner.hpp index 2927909c207baba7066ab66cba64476561ef54c2..6467e14739ee8de805e01585f1821bb28cf7e2a2 100644 --- a/src/amdis/linear_algebra/mtl/KrylovRunner.hpp +++ b/src/amdis/linear_algebra/mtl/KrylovRunner.hpp @@ -30,7 +30,7 @@ namespace AMDiS public: /// Constructor. - explicit KrylovRunner(std::string prefix) + explicit KrylovRunner(std::string const& prefix) : itlSolver_(prefix) { initPrecon(prefix); @@ -89,11 +89,11 @@ namespace AMDiS r -= A * x; // print information about the solution process - itl::cyclic_iteration<double> iter(r, maxIter_, solverInfo.getRelTolerance(), - solverInfo.getAbsTolerance(), + itl::cyclic_iteration<double> iter(r, maxIter_, solverInfo.relTolerance(), + solverInfo.absTolerance(), printCycle_); - iter.set_quite(solverInfo.getInfo() == 0); - iter.suppress_resume(solverInfo.getInfo() == 0); + iter.set_quite(solverInfo.info() == 0); + iter.suppress_resume(solverInfo.info() == 0); int error = itlSolver_(A, x, b, *L_, *R_, iter); solverInfo.setAbsResidual(iter.resid()); @@ -104,7 +104,7 @@ namespace AMDiS private: /// Create left/right preconditioners from parameters given in the init-file - void initPrecon(std::string prefix) + void initPrecon(std::string const& prefix) { // Creator for the left preconditioner std::string leftPreconName = "", rightPreconName = ""; diff --git a/src/amdis/linear_algebra/mtl/UmfpackRunner.hpp b/src/amdis/linear_algebra/mtl/UmfpackRunner.hpp index 59576c9723621d5f32aed865d8d9e2cbc3777d3b..cb97ff50a83a85eb168fec956f114a438a67e019 100644 --- a/src/amdis/linear_algebra/mtl/UmfpackRunner.hpp +++ b/src/amdis/linear_algebra/mtl/UmfpackRunner.hpp @@ -39,7 +39,7 @@ namespace AMDiS public: /// Constructor. Reads UMFPACK parameters from initfile - UmfpackRunnerBase(std::string prefix) + UmfpackRunnerBase(std::string const& prefix) { Parameters::get(prefix + "->store symbolic", storeSymbolic_); // ? Parameters::get(prefix + "->symmetric strategy", symmetricStrategy_); @@ -93,7 +93,7 @@ namespace AMDiS using SolverType = typename Super::SolverType; public: - UmfpackRunner(std::string prefix) + UmfpackRunner(std::string const& prefix) : Super(prefix) {} diff --git a/test/AdaptInfoTest.cpp b/test/AdaptInfoTest.cpp index 78bc27336d92d9ffad68d2f712ba07b0d8a385e9..0b861c7144d0230aac9fbff0aec8d9ee4505a8b4 100644 --- a/test/AdaptInfoTest.cpp +++ b/test/AdaptInfoTest.cpp @@ -18,21 +18,21 @@ int main() std::string str = "0"; adaptInfo.setEstSum(0.1, tp); - AMDIS_TEST_EQ(adaptInfo.getEstSum(tp), 0.1); + AMDIS_TEST_EQ(adaptInfo.estSum(tp), 0.1); adaptInfo.setEstSum(0.2, root_tp); - AMDIS_TEST_EQ(adaptInfo.getEstSum(root_tp), 0.2); + AMDIS_TEST_EQ(adaptInfo.estSum(root_tp), 0.2); - AMDIS_TEST_EQ(adaptInfo.getSize(), 2); + AMDIS_TEST_EQ(adaptInfo.size(), 2); adaptInfo.setEstSum(0.3, "0"); - AMDIS_TEST_EQ(adaptInfo.getEstSum(tp), 0.3); + AMDIS_TEST_EQ(adaptInfo.estSum(tp), 0.3); adaptInfo.setEstSum(0.4, 0_c); - AMDIS_TEST_EQ(adaptInfo.getEstSum(tp), 0.4); + AMDIS_TEST_EQ(adaptInfo.estSum(tp), 0.4); adaptInfo.setEstSum(0.5, str); - AMDIS_TEST_EQ(adaptInfo.getEstSum(tp), 0.5); + AMDIS_TEST_EQ(adaptInfo.estSum(tp), 0.5); return report_errors(); } \ No newline at end of file diff --git a/test/DiscreteFunctionTest.cpp b/test/DiscreteFunctionTest.cpp index 93314e544df4576d754422f82b433f73e0062cc1..8d483c8f77bc1c391c7d16aea5913ef43598dfe7 100644 --- a/test/DiscreteFunctionTest.cpp +++ b/test/DiscreteFunctionTest.cpp @@ -40,9 +40,9 @@ int main(int argc, char** argv) prob.initialize(INIT_ALL); // create 3 copies of the solution vector - auto U0 = prob.getSolutionVector(); - auto U1 = prob.getSolutionVector(); - auto U2 = prob.getSolutionVector(); + auto U0 = prob.solutionVector(); + auto U1 = prob.solutionVector(); + auto U2 = prob.solutionVector(); auto u0 = makeDOFVectorView(U0); auto u1 = makeDOFVectorView(U1); diff --git a/test/ExpressionsTest.cpp b/test/ExpressionsTest.cpp index 1422351b7f5308589fb67ce532c1bd662fcde311..837ba91945020464d79e5707a9dc2d6010493098 100644 --- a/test/ExpressionsTest.cpp +++ b/test/ExpressionsTest.cpp @@ -25,7 +25,7 @@ int main(int argc, char** argv) ElliptProblem prob("ellipt"); prob.initialize(INIT_ALL); - auto u = prob.getSolution(_0); + auto u = prob.solution(_0); // eval a functor at global coordinates (at quadrature points) auto expr1 = evalAtQP([](Dune::FieldVector<double, 2> const& x) { return x[0] + x[1]; }); @@ -42,7 +42,7 @@ int main(int argc, char** argv) auto expr8 = X(0); // Evaluation of the DOFVector (component) at quadrature points - auto expr9 = prob.getSolution(_0); + auto expr9 = prob.solution(_0); // --------------------------------- diff --git a/test/IntegrateTest.cpp b/test/IntegrateTest.cpp index d348dda35b1c7f71f75c5bbdb792b43603dcbc6b..539d870fce5a3ef28bb17b5d4d417138a84bf93f 100644 --- a/test/IntegrateTest.cpp +++ b/test/IntegrateTest.cpp @@ -23,7 +23,7 @@ int main(int argc, char** argv) ElliptProblem prob("ellipt"); prob.initialize(INIT_ALL); - auto u = prob.getSolution(0); + auto u = prob.solution(0); auto gv = u.basis().gridView(); u << 1.0; diff --git a/test/MarkerTest.cpp b/test/MarkerTest.cpp index f34e0f758aed15022c329a350fa2bf614e1b9acf..18968a7b0f5333200dd00a8accdf17a63bc0ffc5 100644 --- a/test/MarkerTest.cpp +++ b/test/MarkerTest.cpp @@ -50,7 +50,7 @@ int main(int argc, char** argv) break; } - prob.getSolution().interpolate(markerFunc); + prob.solution().interpolate(markerFunc); prob.writeFiles(adaptInfo); AMDiS::finalize();