diff --git a/examples/CMakeLists.txt b/examples/CMakeLists.txt index 239b1ad4d250898f07c6cf27f1d2c042c056a4c5..740e0444e564905231bbbcd0d5ffad913746138a 100644 --- a/examples/CMakeLists.txt +++ b/examples/CMakeLists.txt @@ -35,4 +35,4 @@ endif () if (dune-foamgrid_FOUND) add_amdis_executable(NAME surface.2d SOURCES surface.cc DIM 2 DOW 3 ALBERTA_GRID) add_dependencies(examples surface.2d) -endif () \ No newline at end of file +endif () diff --git a/examples/boundary.cc b/examples/boundary.cc index d08280c7cdc551473b84b1f11400bbf026261f5f..f54924aee43462f96b8f1796e8d525582d9545ec 100644 --- a/examples/boundary.cc +++ b/examples/boundary.cc @@ -11,10 +11,8 @@ #include <amdis/Integrate.hpp> #include <amdis/ProblemStat.hpp> #include <amdis/LocalOperators.hpp> -#include <amdis/common/Literals.hpp> using namespace AMDiS; -using namespace Dune::Indices; // 1 component with polynomial degree 2 using Param = YaspGridBasis<GRIDDIM, 2>; @@ -28,7 +26,7 @@ void run(SetBoundary setBoundary) setBoundary(prob.boundaryManager()); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); - prob.addMatrixOperator(opL, 0, 0); + prob.addMatrixOperator(opL); auto f = [](auto const& x) { double r2 = dot(x,x); @@ -36,18 +34,18 @@ void run(SetBoundary setBoundary) return -(400.0 * r2 - 20.0 * x.size()) * ux; }; auto opForce = makeOperator(tag::test{}, f, 6); - prob.addVectorOperator(opForce, 0); + prob.addVectorOperator(opForce); // set boundary condition auto g = [](auto const& x){ return std::exp(-10.0 * dot(x,x)); }; - prob.addDirichletBC(BoundaryType{1}, 0, 0, g); + prob.addDirichletBC(BoundaryType{1}, g); AdaptInfo adaptInfo("adapt"); prob.assemble(adaptInfo); prob.solve(adaptInfo); - double errorL2 = integrate(sqr(g - prob.solution(0)), prob.gridView(), 6); + double errorL2 = integrate(sqr(g - prob.solution()), prob.gridView(), 6); msg("error_L2 = {}", errorL2); } @@ -69,17 +67,17 @@ void run_periodic() prob.boundaryManager().setBoxBoundary({-1,-1,1,1}); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); - prob.addMatrixOperator(opL, 0, 0); + prob.addMatrixOperator(opL); auto opForce = makeOperator(tag::test{}, 1.0, 6); - prob.addVectorOperator(opForce, 0); + prob.addVectorOperator(opForce); // set boundary condition auto g = [](auto const& x) { return std::sin(0.5*M_PI*x[1])*std::sin(2*M_PI * (0.25 + x[0])) + std::cos(0.5*M_PI*x[1])*std::sin(2*M_PI * (-0.25 + x[0])); }; - prob.addDirichletBC(BoundaryType{1}, 0, 0, g); + prob.addDirichletBC(BoundaryType{1}, g); typename ElliptProblem::WorldMatrix A{{1.0,0.0}, {0.0,1.0}}; typename ElliptProblem::WorldVector b{1.0, 0.0}; diff --git a/examples/convection_diffusion.cc b/examples/convection_diffusion.cc index 0ba82e9a98dd213590372d85d9e6b505464ea488..e613e13753e30ee0a7f8551cc86b6498fffcdbb6 100644 --- a/examples/convection_diffusion.cc +++ b/examples/convection_diffusion.cc @@ -6,7 +6,6 @@ #include <amdis/AMDiS.hpp> #include <amdis/ProblemStat.hpp> #include <amdis/localoperators/ConvectionDiffusionOperator.hpp> -#include <amdis/common/Literals.hpp> using namespace AMDiS; @@ -19,20 +18,18 @@ int main(int argc, char** argv) { Environment env(argc, argv); - using namespace Dune::Indices; - ElliptProblem prob("ellipt"); prob.initialize(INIT_ALL); // -div(A*grad(u)) + div(b*u) + c*u = f auto opCD = convectionDiffusion(/*A=*/1.0, /*b=*/0.0, /*c=*/1.0, /*f=*/1.0); - prob.addMatrixOperator(opCD, 0, 0); - prob.addVectorOperator(opCD, 0); + prob.addMatrixOperator(opCD); + prob.addVectorOperator(opCD); // set boundary condition auto predicate = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8; }; // define boundary auto dbcValues = [](auto const& x){ return 0.0; }; // set value - prob.addDirichletBC(predicate, 0, 0, dbcValues); + prob.addDirichletBC(predicate, dbcValues); AdaptInfo adaptInfo("adapt"); prob.assemble(adaptInfo); diff --git a/examples/ellipt.cc b/examples/ellipt.cc index 147f6d2ec6836baf75890f8c2cfe4a28a7d5d396..419242e94bc91f93279757062d3d8427d96c2308 100644 --- a/examples/ellipt.cc +++ b/examples/ellipt.cc @@ -8,7 +8,6 @@ using Grid = Dune::YaspGrid<GRIDDIM>; using namespace AMDiS; -using namespace Dune::Indices; int main(int argc, char** argv) { @@ -34,13 +33,13 @@ 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); auto opForce = makeOperator(tag::test{}, f, 6); - prob.addVectorOperator(opForce, 0); + prob.addVectorOperator(opForce); // set boundary condition - prob.addDirichletBC(1, 0, 0, g); + prob.addDirichletBC(1, g); AdaptInfo adaptInfo("adapt"); @@ -65,9 +64,9 @@ int main(int argc, char** argv) prob.assemble(adaptInfo); prob.solve(adaptInfo); - double errorL2 = integrate(sqr(g - prob.solution(0)), prob.gridView(), 6); + double errorL2 = integrate(sqr(g - prob.solution()), prob.gridView(), 6); errL2.push_back(std::sqrt(errorL2)); - double errorH1 = errorL2 + integrate(unary_dot(grad_g - gradientAtQP(prob.solution(0))), prob.gridView(), 6); + double errorH1 = errorL2 + integrate(unary_dot(grad_g - gradientAtQP(prob.solution())), prob.gridView(), 6); errH1.push_back(std::sqrt(errorH1)); } diff --git a/examples/heat.cc b/examples/heat.cc index 9fe769a479bf6c8b75a8fc7dd9b7d34de6cbfd1c..d419d0ddb0c739f754ddbb1c97cbaeaac3b46997 100644 --- a/examples/heat.cc +++ b/examples/heat.cc @@ -9,7 +9,6 @@ #include <amdis/ProblemInstat.hpp> #include <amdis/ProblemStat.hpp> #include <amdis/GridFunctions.hpp> -#include <amdis/common/Literals.hpp> using namespace AMDiS; @@ -34,23 +33,23 @@ int main(int argc, char** argv) auto invTau = std::ref(probInstat.invTau()); auto opTimeLhs = makeOperator(tag::test_trial{}, invTau); - prob.addMatrixOperator(opTimeLhs, 0, 0); + prob.addMatrixOperator(opTimeLhs); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); - prob.addMatrixOperator(opL, 0, 0); + prob.addMatrixOperator(opL); auto opTimeRhs = makeOperator(tag::test{}, - invokeAtQP([invTau](double u) { return u * invTau.get(); }, prob.solution(0)), 2); - prob.addVectorOperator(opTimeRhs, 0); + invokeAtQP([invTau](double u) { return u * invTau.get(); }, prob.solution()), 2); + prob.addVectorOperator(opTimeRhs); auto opForce = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0); - prob.addVectorOperator(opForce, 0); + prob.addVectorOperator(opForce); // set boundary condition auto predicate = [](auto const& p){ return p[0] < 1.e-8 || p[1] < 1.e-8; }; auto dbcValues = [](auto const& p){ return 0.0; }; - prob.addDirichletBC(predicate, 0, 0, dbcValues); + prob.addDirichletBC(predicate, dbcValues); AdaptInstationary adapt("adapt", prob, adaptInfo, probInstat, adaptInfo); adapt.adapt(); diff --git a/examples/init/boundary.dat.2d b/examples/init/boundary.dat.2d index d4e033a7e5ffd03b146ad1b891e4bc4b12ebe735..411c5c9a3ee43242ab84c4d78ab047bd9947def2 100644 --- a/examples/init/boundary.dat.2d +++ b/examples/init/boundary.dat.2d @@ -8,7 +8,7 @@ ellipt->solver->relative tolerance: 1.e-10 ellipt->solver->info: 1 ellipt->solver->precon: ilu -ellipt->output[0]->format: vtk -ellipt->output[0]->filename: boundary.2d -ellipt->output[0]->name: u -ellipt->output[0]->output directory: ./output +ellipt->output->format: vtk +ellipt->output->filename: boundary.2d +ellipt->output->name: u +ellipt->output->output directory: ./output diff --git a/examples/init/cahn_hilliard.dat.2d b/examples/init/cahn_hilliard.dat.2d index 7af9db935e5e680770a887bd46ee2bf1dbf9e9de..2013aead79285dbdfface47ecadde799da2e4e61 100644 --- a/examples/init/cahn_hilliard.dat.2d +++ b/examples/init/cahn_hilliard.dat.2d @@ -16,7 +16,7 @@ ch->solver: direct ch->solver->max iteration: 1000 ch->solver->relative tolerance: 1e-6 ch->solver->break if tolerance not reached: 1 -ch->solver->info: 2 +ch->solver->info: 0 ch->output[0]->format: vtk ch->output[0]->filename: phi.2d @@ -34,4 +34,4 @@ ch->output[1]->animation: 1 adapt->timestep: 0.001 adapt->start time: 0.0 -adapt->end time: 1.0 +adapt->end time: 0.5 diff --git a/examples/init/ellipt.dat.2d b/examples/init/ellipt.dat.2d index 0ef31af1f0abbaadffad556eb61ea3aac807830c..ac3b7e1d9bc93e161d9b80c3ee7c45813ac5c162 100644 --- a/examples/init/ellipt.dat.2d +++ b/examples/init/ellipt.dat.2d @@ -13,7 +13,7 @@ ellipt->solver->max iteration: 10000 ellipt->solver->relative tolerance: 1.e-10 ellipt->solver->precon: ilu -ellipt->output[0]->format: vtk backup -ellipt->output[0]->filename: ellipt.2d -ellipt->output[0]->name: u -ellipt->output[0]->output directory: ./output +ellipt->output->format: vtk backup +ellipt->output->filename: ellipt.2d +ellipt->output->name: u +ellipt->output->output directory: ./output diff --git a/examples/init/ellipt.dat.3d b/examples/init/ellipt.dat.3d index baef5567dc843c2cac02f28f5503726d8effe9ad..edbe4bf38c09409c8be660041f1e51ae015c890f 100644 --- a/examples/init/ellipt.dat.3d +++ b/examples/init/ellipt.dat.3d @@ -13,7 +13,7 @@ ellipt->solver->max iteration: 10000 ellipt->solver->relative tolerance: 1.e-10 ellipt->solver->precon: ilu -ellipt->output[0]->format: vtk -ellipt->output[0]->filename: ellipt.3d -ellipt->output[0]->name: u -ellipt->output[0]->output directory: ./output +ellipt->output->format: vtk +ellipt->output->filename: ellipt.3d +ellipt->output->name: u +ellipt->output->output directory: ./output diff --git a/examples/init/heat.dat.2d b/examples/init/heat.dat.2d index c23a811c95db9563f6ed48bfbf43682dade1266d..fef62797813e2f1fa5f5bbb2b711a0328f543a20 100644 --- a/examples/init/heat.dat.2d +++ b/examples/init/heat.dat.2d @@ -12,12 +12,12 @@ heat->solver->absolute tolerance: 1e-6 heat->solver->info: 1 heat->solver->precon: ilu -heat->output[0]->format: vtk backup -heat->output[0]->name: u -heat->output[0]->filename: heat.2d -heat->output[0]->output directory: output -heat->output[0]->mode: 1 -heat->output[0]->animation: 1 +heat->output->format: vtk backup +heat->output->name: u +heat->output->filename: heat.2d +heat->output->output directory: output +heat->output->mode: 1 +heat->output->animation: 1 adapt->timestep: 0.1 adapt->start time: 0.0 diff --git a/examples/init/surface.dat.2d b/examples/init/surface.dat.2d index 282ea4d0c865db5294402a95335563888cd27f4e..971637a66301075531a907481b0eaef9a7cfb828 100644 --- a/examples/init/surface.dat.2d +++ b/examples/init/surface.dat.2d @@ -9,7 +9,7 @@ surface->solver->relative tolerance: 1.e-6 surface->solver->info: 1 surface->solver->precon: ilu -surface->output[0]->format: vtk -surface->output[0]->filename: surface.2d -surface->output[0]->name: u -surface->output[0]->output directory: ./output +surface->output->format: vtk +surface->output->filename: surface.2d +surface->output->name: u +surface->output->output directory: ./output diff --git a/examples/neumann.cc b/examples/neumann.cc index 1037b8137d3dee812c58ef727ef7b41ae9aab4f1..7e6e6bfc7342330aabf1ef686bc779bde3f5a289 100644 --- a/examples/neumann.cc +++ b/examples/neumann.cc @@ -15,10 +15,8 @@ #include <amdis/Integrate.hpp> #include <amdis/ProblemStat.hpp> #include <amdis/LocalOperators.hpp> -#include <amdis/common/Literals.hpp> using namespace AMDiS; -using namespace Dune::Indices; template <class Grid> void run(Grid& grid) @@ -28,18 +26,18 @@ void run(Grid& grid) using Traits = LagrangeBasis<Grid, 2>; ProblemStat<Traits> prob("ellipt", grid); prob.initialize(INIT_ALL); - prob.boundaryManager().setBoxBoundary({1,2,2,1}); + prob.boundaryManager()->setBoxBoundary({1,2,2,1}); auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0); - prob.addMatrixOperator(opL, 0, 0); + prob.addMatrixOperator(opL); // set dirichlet boundary condition auto g = [](auto const& x) { return 0.0; }; - prob.addDirichletBC(BoundaryType{1}, 0, 0, g); + prob.addDirichletBC(BoundaryType{1}, g); // set neumann boundary condition auto opNeumann = makeOperator(tag::test{}, 1.0); - prob.addVectorOperator(BoundaryType{2}, opNeumann, 0); + prob.addVectorOperator(BoundaryType{2}, opNeumann); AdaptInfo adaptInfo("adapt"); diff --git a/examples/periodic.cc b/examples/periodic.cc index 8236d99c8bf05b7df63d1141fb3aaee34359360c..35f7a81d81ff4809af15ebf10ef526ce4e50ad99 100644 --- a/examples/periodic.cc +++ b/examples/periodic.cc @@ -11,7 +11,6 @@ #include <amdis/Integrate.hpp> #include <amdis/ProblemStat.hpp> #include <amdis/LocalOperators.hpp> -#include <amdis/common/Literals.hpp> using namespace AMDiS; diff --git a/examples/surface.cc b/examples/surface.cc index 7de309bcd43b39c4d818e47302047b27926a1a1f..f0bed4c8871b379fc6870a2761a366b85f037be7 100644 --- a/examples/surface.cc +++ b/examples/surface.cc @@ -3,7 +3,6 @@ #include <amdis/Integrate.hpp> #include <amdis/LocalOperators.hpp> #include <amdis/ProblemStat.hpp> -#include <amdis/common/Literals.hpp> #include <dune/grid/albertagrid/albertareader.hh> #include <dune/foamgrid/foamgrid.hh> @@ -12,7 +11,6 @@ #include "SphereMapping.hpp" using namespace AMDiS; -using namespace Dune::Indices; using Grid = Dune::FoamGrid<2,3>; using Basis = LagrangeBasis<Grid, 1>; @@ -44,14 +42,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); double kappa = Parameters::get<double>("surface->kappa").value_or(1.0); auto opM = makeOperator(tag::test_trial{}, -kappa*kappa); - prob.addMatrixOperator(opM, 0, 0); + prob.addMatrixOperator(opM); auto opForce = makeOperator(tag::test{}, X(0) + X(1) + X(2)); - prob.addVectorOperator(opForce, 0); + prob.addVectorOperator(opForce); AdaptInfo adaptInfo("adapt"); prob.assemble(adaptInfo); diff --git a/src/amdis/ProblemStat.hpp b/src/amdis/ProblemStat.hpp index 383d0b8894d754b7d071f6a79447b4b40f4df8e2..13bffdc68e9579cfab3cd0f69918d7118418cfc5 100644 --- a/src/amdis/ProblemStat.hpp +++ b/src/amdis/ProblemStat.hpp @@ -266,6 +266,12 @@ namespace AMDiS RowTreePath row, ColTreePath col, Values const& values); + template <class Identifier, class Values> + void addDirichletBC(Identifier&& id, Values&& values) + { + addDirichletBC(FWD(id), RootTreePath{}, RootTreePath{}, FWD(values)); + } + /// Add a periodic boundary conditions to the system, by specifying a face transformation /// y = A*x + b of coordinates. We assume, that A is orthonormal. void addPeriodicBC(BoundaryType id, WorldMatrix const& A, WorldVector const& b); diff --git a/src/amdis/ProblemStatTraits.hpp b/src/amdis/ProblemStatTraits.hpp index 1dfbd6d8a6e2f79ebb269c3856fa01ea6b71bca7..742581e7cd6961b52aa5a1706256f50429b522b3 100644 --- a/src/amdis/ProblemStatTraits.hpp +++ b/src/amdis/ProblemStatTraits.hpp @@ -18,6 +18,17 @@ namespace AMDiS template <bool same, int... degs> struct LagrangePreBasisCreatorImpl; + // specialization for single node basis + template <int deg> + struct LagrangePreBasisCreatorImpl<true, deg> + { + static auto create() + { + using namespace Dune::Functions::BasisFactory; + return lagrange<deg>(); + } + }; + // specialization for composite basis template <int... degs> struct LagrangePreBasisCreatorImpl<false, degs...> diff --git a/test/ExpressionsTest.cpp b/test/ExpressionsTest.cpp index aa1ecae2532ef3811b80c5e3009d3d5b70470f6f..ae18fb438ce63809ba1bb69aff51363bd18eb9ad 100644 --- a/test/ExpressionsTest.cpp +++ b/test/ExpressionsTest.cpp @@ -7,7 +7,6 @@ #include <amdis/Integrate.hpp> #include <amdis/LocalOperators.hpp> #include <amdis/ProblemStat.hpp> -#include <amdis/common/Literals.hpp> #include <amdis/common/FieldMatVec.hpp> @@ -20,12 +19,10 @@ int main(int argc, char** argv) { Environment env(argc, argv); - using namespace Dune::Indices; - ElliptProblem prob("ellipt"); prob.initialize(INIT_ALL); - auto u = prob.solution(_0); + auto u = prob.solution(); // 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 +39,7 @@ int main(int argc, char** argv) auto expr8 = X(0); // Evaluation of the DOFVector (component) at quadrature points - auto expr9 = prob.solution(_0); + auto expr9 = prob.solution(); // --------------------------------- diff --git a/test/IntegrateTest.cpp b/test/IntegrateTest.cpp index adb28a08472de8a7f56dd96179519532c51c6ca7..c4a97a1a05c7146ae4f9c41673cd496ad4825287 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.solution(0); + auto u = prob.solution(); auto gv = u.basis().gridView(); u << 1.0;