diff --git a/test/freeflow/navierstokes/analyticalsolutionvectors.hh b/test/freeflow/navierstokes/analyticalsolutionvectors.hh index fa41298e1965fa447387ed217b3635d540f2cbe8..ec5f4f8d19f851a3f1a5aca2258951ef099ee5c4 100644 --- a/test/freeflow/navierstokes/analyticalsolutionvectors.hh +++ b/test/freeflow/navierstokes/analyticalsolutionvectors.hh @@ -29,71 +29,57 @@ namespace Dumux { * \ingroup NavierStokesTests * \brief Creates and provides the analytical solution in a form that can be written into the vtk output files */ -template<class Scalar, class FaceSolutionVector, class CellCenterSolutionVector, class GridGeometry, class Problem, class Indices> +template<class Problem, class Scalar = double> class NavierStokesAnalyticalSolutionVectors { + using GridGeometry = std::decay_t<decltype(std::declval<Problem>().gridGeometry())>; using GridView = typename GridGeometry::GridView; using FVElementGeometry = typename GridGeometry::LocalView; using SubControlVolume = typename FVElementGeometry::SubControlVolume; using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace; - static constexpr auto dimWorld = GridGeometry::GridView::dimensionworld; + static constexpr int dimWorld = GridGeometry::GridView::dimensionworld; using VelocityVector = Dune::FieldVector<Scalar, dimWorld>; + using Indices = typename Problem::Indices; public: - NavierStokesAnalyticalSolutionVectors(std::shared_ptr<const Problem> problem) + NavierStokesAnalyticalSolutionVectors(std::shared_ptr<const Problem> problem, Scalar tInitial = 0.0) : problem_(problem) - { } + { update(tInitial); } /*! - * \brief Creates the instationary analytical solution in a form that can be written into the vtk output files + * \brief Creates the analytical solution in a form that can be written into the vtk output files */ - void update(Scalar time) + void update(Scalar time = 0.0) { - auto instationaryAnalyticalPressureSolution = [&](const SubControlVolume& scv) - { - return problem_->analyticalSolution(scv.dofPosition(), time)[Indices::pressureIdx]; - }; - - auto instationaryAnalyticalVelocitySolution = [&](const SubControlVolumeFace& scvf) - { - return problem_->analyticalSolution(scvf.center(), time)[Indices::velocity(scvf.directionIndex())]; - }; - - auto instationaryAnalyticalVelocitySolutionAtScvCenter = [&](const SubControlVolume& scv) - { - return problem_->analyticalSolution(scv.center(), time); - }; + analyticalPressure_.resize(problem_->gridGeometry().numCellCenterDofs()); + analyticalVelocity_.resize(problem_->gridGeometry().numCellCenterDofs()); + analyticalVelocityOnFace_.resize(problem_->gridGeometry().numFaceDofs()); - createAnalyticalSolution_(instationaryAnalyticalPressureSolution, instationaryAnalyticalVelocitySolution, instationaryAnalyticalVelocitySolutionAtScvCenter); - } + auto fvGeometry = localView(problem_->gridGeometry()); - /*! - * \brief Creates the stationary analytical solution in a form that can be written into the vtk output files - */ - void update() - { - auto analyticalPressureSolution = [&](const SubControlVolume& scv) - { - return problem_->analyticalSolution(scv.dofPosition())[Indices::pressureIdx]; - }; - - auto analyticalVelocitySolution = [&](const SubControlVolumeFace& scvf) + for (const auto& element : elements((problem_->gridGeometry()).gridView())) { - return problem_->analyticalSolution(scvf.center())[Indices::velocity(scvf.directionIndex())]; - }; + fvGeometry.bindElement(element); + for (const auto& scv : scvs(fvGeometry)) + { + // velocities on faces + for (const auto& scvf : scvfs(fvGeometry)) + { + analyticalVelocityOnFace_[scvf.dofIndex()][scvf.directionIndex()] = problem_->analyticalSolution(scvf.center(), time)[Indices::velocity(scvf.directionIndex())]; + } - auto analyticalVelocitySolutionAtScvCenter = [&](const SubControlVolume& scv) - { - return problem_->analyticalSolution(scv.center()); - }; + analyticalPressure_[scv.dofIndex()] = problem_->analyticalSolution(scv.dofPosition(), time)[Indices::pressureIdx]; - createAnalyticalSolution_(analyticalPressureSolution, analyticalVelocitySolution, analyticalVelocitySolutionAtScvCenter); + for (int dirIdx = 0; dirIdx < dimWorld; ++dirIdx) + analyticalVelocity_[scv.dofIndex()][dirIdx] = problem_->analyticalSolution(scv.center(), time)[Indices::velocity(dirIdx)]; + } + } } /*! * \brief Returns the analytical solution for the pressure */ - auto& getAnalyticalPressureSolution() const + const std::vector<Scalar>& getAnalyticalPressureSolution() const { return analyticalPressure_; } @@ -101,7 +87,7 @@ public: /*! * \brief Returns the analytical solution for the velocity */ - auto& getAnalyticalVelocitySolution() const + const std::vector<VelocityVector>& getAnalyticalVelocitySolution() const { return analyticalVelocity_; } @@ -109,47 +95,15 @@ public: /*! * \brief Returns the analytical solution for the velocity at the faces */ - auto& getAnalyticalVelocitySolutionOnFace() const + const std::vector<VelocityVector>& getAnalyticalVelocitySolutionOnFace() const { return analyticalVelocityOnFace_; } private: - /*! - * \brief Creates the analytical solution in a form that can be written into the vtk output files - */ - template <class LambdaA, class LambdaB, class LambdaC> - void createAnalyticalSolution_(const LambdaA& analyticalPressureSolution, - const LambdaB& analyticalVelocitySolution, - const LambdaC& analyticalVelocitySolutionAtScvCenter) - { - analyticalPressure_.resize((problem_->gridGeometry()).numCellCenterDofs()); - analyticalVelocity_.resize((problem_->gridGeometry()).numCellCenterDofs()); - analyticalVelocityOnFace_.resize((problem_->gridGeometry()).numFaceDofs()); - - for (const auto& element : elements((problem_->gridGeometry()).gridView())) - { - auto fvGeometry = localView((problem_->gridGeometry())); - fvGeometry.bindElement(element); - for (auto&& scv : scvs(fvGeometry)) - { - // velocities on faces - for (auto&& scvf : scvfs(fvGeometry)) - { - analyticalVelocityOnFace_[scvf.dofIndex()][scvf.directionIndex()] = analyticalVelocitySolution(scvf); - } - - analyticalPressure_[scv.dofIndex()] = analyticalPressureSolution(scv); - - for(int dirIdx = 0; dirIdx < dimWorld; ++dirIdx) - analyticalVelocity_[scv.dofIndex()][dirIdx] = analyticalVelocitySolutionAtScvCenter(scv)[Indices::velocity(dirIdx)]; - } - } - } - std::shared_ptr<const Problem> problem_; - CellCenterSolutionVector analyticalPressure_; + std::vector<Scalar> analyticalPressure_; std::vector<VelocityVector> analyticalVelocity_; std::vector<VelocityVector> analyticalVelocityOnFace_; }; diff --git a/test/freeflow/navierstokes/angeli/main.cc b/test/freeflow/navierstokes/angeli/main.cc index d1769df4af96457fdbb058a19cd2908e87a16097..4160b03bee35b0ebdc1b2572c0ff0ead55e3a8da 100644 --- a/test/freeflow/navierstokes/angeli/main.cc +++ b/test/freeflow/navierstokes/angeli/main.cc @@ -46,7 +46,7 @@ #include "properties.hh" -#include "../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> int main(int argc, char** argv) { @@ -109,19 +109,11 @@ int main(int argc, char** argv) auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); gridVariables->init(x); - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - analyticalSolVectors.update(tStart); - // initialize the vtk output module StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); using IOFields = GetPropType<TypeTag, Properties::IOFields>; IOFields::initOutputModule(vtkWriter); // Add model specific output fields + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem, tStart); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); vtkWriter.addFaceField(analyticalSolVectors.getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact"); @@ -159,7 +151,7 @@ int main(int argc, char** argv) using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using L2Error = NavierStokesTestL2Error<Scalar, ModelTraits, PrimaryVariables>; - const auto [l2errorAbs, l2errorRel] = L2Error::calculateL2Error(*problem, x); + const auto [l2errorAbs, l2errorRel] = L2Error::calculateL2Error(*problem, x, timeLoop->time() + timeLoop->timeStepSize()); const int numCellCenterDofs = gridGeometry->numCellCenterDofs(); const int numFaceDofs = gridGeometry->numFaceDofs(); std::cout << std::setprecision(8) << "** L2 error (abs/rel) for " diff --git a/test/freeflow/navierstokes/angeli/problem.hh b/test/freeflow/navierstokes/angeli/problem.hh index b8117d40c4a715a4a8358388a240a3f0409cfbf3..1bb5ef935aba695d6011778ba9ba93101a344fee 100644 --- a/test/freeflow/navierstokes/angeli/problem.hh +++ b/test/freeflow/navierstokes/angeli/problem.hh @@ -164,27 +164,22 @@ public: /*! * \brief Returns the analytical solution of the problem at a given time and position. * \param globalPos The global position + * \param time A parameter for consistent signatures. It is ignored here as this is a stationary test */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar t) const + PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar time) const { const Scalar x = globalPos[0]; const Scalar y = globalPos[1]; PrimaryVariables values; - values[Indices::pressureIdx] = - 0.25 * std::exp(-10.0 * kinematicViscosity_ * M_PI * M_PI * t) * M_PI * M_PI * (4.0 * std::cos(2.0 * M_PI * x) + std::cos(4.0 * M_PI * y))*rho_; - values[Indices::velocityXIdx] = - 2.0 * M_PI * std::exp(- 5.0 * kinematicViscosity_ * M_PI * M_PI * t) * std::cos(M_PI * x) * std::sin(2.0 * M_PI * y); - values[Indices::velocityYIdx] = M_PI * std::exp(- 5.0 * kinematicViscosity_ * M_PI * M_PI * t) * std::sin(M_PI * x) * std::cos(2.0 * M_PI * y); + values[Indices::pressureIdx] = - 0.25 * std::exp(-10.0 * kinematicViscosity_ * M_PI * M_PI * time) * M_PI * M_PI * (4.0 * std::cos(2.0 * M_PI * x) + std::cos(4.0 * M_PI * y))*rho_; + values[Indices::velocityXIdx] = - 2.0 * M_PI * std::exp(- 5.0 * kinematicViscosity_ * M_PI * M_PI * time) * std::cos(M_PI * x) * std::sin(2.0 * M_PI * y); + values[Indices::velocityYIdx] = M_PI * std::exp(- 5.0 * kinematicViscosity_ * M_PI * M_PI * time) * std::sin(M_PI * x) * std::cos(2.0 * M_PI * y); return values; } - //TODO remove as soon as error calculation is able to deal with time-dependent analytical solution - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const - { - return analyticalSolution(globalPos, time_); - } - // \} /*! diff --git a/test/freeflow/navierstokes/channel/1d/main.cc b/test/freeflow/navierstokes/channel/1d/main.cc index 7e8521fac58a14ffc2527ca88694325b6027f5a1..3c3631f149836ecc719d18f164d4dff67e524baf 100644 --- a/test/freeflow/navierstokes/channel/1d/main.cc +++ b/test/freeflow/navierstokes/channel/1d/main.cc @@ -44,7 +44,7 @@ #include "properties.hh" -#include "../../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> int main(int argc, char** argv) { @@ -93,22 +93,16 @@ int main(int argc, char** argv) auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); gridVariables->init(x); - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - analyticalSolVectors.update(); - // intialize the vtk output module using IOFields = GetPropType<TypeTag, Properties::IOFields>; StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); IOFields::initOutputModule(vtkWriter); // Add model specific output fields + + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); vtkWriter.addFaceField(analyticalSolVectors.getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact"); + vtkWriter.write(0.0); // the assembler with time loop for instationary problem diff --git a/test/freeflow/navierstokes/channel/1d/problem.hh b/test/freeflow/navierstokes/channel/1d/problem.hh index 66b9823b1ec39348a8a45145b434f2072187e0ef..7f345c52b9f6def6bc1e65a2b1a13b59716c5a77 100644 --- a/test/freeflow/navierstokes/channel/1d/problem.hh +++ b/test/freeflow/navierstokes/channel/1d/problem.hh @@ -53,7 +53,6 @@ class NavierStokesAnalyticProblem : public NavierStokesProblem<TypeTag> using BoundaryTypes = Dumux::NavierStokesBoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>; using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using NumEqVector = Dumux::NumEqVector<PrimaryVariables>; using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>; @@ -67,6 +66,8 @@ class NavierStokesAnalyticProblem : public NavierStokesProblem<TypeTag> using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>; public: + using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; + NavierStokesAnalyticProblem(std::shared_ptr<const GridGeometry> gridGeometry) : ParentType(gridGeometry) { @@ -207,8 +208,9 @@ public: * \brief Returns the analytical solution of the problem at a given position * * \param globalPos The global position + * \param time A parameter for consistent signatures. It is ignored here as this is a stationary test. */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const + PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar time = 0.0) const { PrimaryVariables values; values[Indices::pressureIdx] = p(globalPos); diff --git a/test/freeflow/navierstokes/channel/2d/main.cc b/test/freeflow/navierstokes/channel/2d/main.cc index f0e9352d90d0f4bc53e1eb70a6a5f63b945b4fa6..1c68c547e0b221097e3ab612b00c9c2bb54a420c 100644 --- a/test/freeflow/navierstokes/channel/2d/main.cc +++ b/test/freeflow/navierstokes/channel/2d/main.cc @@ -45,7 +45,7 @@ #include "properties.hh" -#include "../../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> int main(int argc, char** argv) { @@ -110,12 +110,6 @@ int main(int argc, char** argv) auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); gridVariables->init(x); - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - // initialize the vtk output module using IOFields = GetPropType<TypeTag, Properties::IOFields>; StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); @@ -125,8 +119,7 @@ int main(int argc, char** argv) if (problem->hasAnalyticalSolution()) { - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - analyticalSolVectors.update(); + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); diff --git a/test/freeflow/navierstokes/channel/2d/problem.hh b/test/freeflow/navierstokes/channel/2d/problem.hh index a39a7b46eb96082aff6715bf9a85bfa3606c81c2..89389521fc6164b808d9b0ba42e9f4a129d3f6d2 100644 --- a/test/freeflow/navierstokes/channel/2d/problem.hh +++ b/test/freeflow/navierstokes/channel/2d/problem.hh @@ -57,7 +57,6 @@ class ChannelTestProblem : public NavierStokesProblem<TypeTag> using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; using FVElementGeometry = typename GridGeometry::LocalView; using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using NumEqVector = Dumux::NumEqVector<PrimaryVariables>; using Scalar = GetPropType<TypeTag, Properties::Scalar>; @@ -76,6 +75,8 @@ class ChannelTestProblem : public NavierStokesProblem<TypeTag> }; public: + using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; + ChannelTestProblem(std::shared_ptr<const GridGeometry> gridGeometry) : ParentType(gridGeometry) { @@ -269,8 +270,9 @@ public: * \brief Return the analytical solution of the problem at a given position * * \param globalPos The global position + * \param time A parameter for consistent signatures. It is ignored here as this analytical solution is for a stationary version of the test only. */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const + PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar time = 0.0) const { Scalar x = globalPos[0]; Scalar y = globalPos[1]; diff --git a/test/freeflow/navierstokes/donea/main.cc b/test/freeflow/navierstokes/donea/main.cc index 73bddb0732682e87c8c672a8e3c24f56646be02f..03ca64680bbe4622678c8e94d877ca621f825ba7 100644 --- a/test/freeflow/navierstokes/donea/main.cc +++ b/test/freeflow/navierstokes/donea/main.cc @@ -45,7 +45,7 @@ #include "properties.hh" -#include "../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> int main(int argc, char** argv) { @@ -95,22 +95,16 @@ int main(int argc, char** argv) auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); gridVariables->init(x); - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - analyticalSolVectors.update(); - // intialize the vtk output module using IOFields = GetPropType<TypeTag, Properties::IOFields>; StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); IOFields::initOutputModule(vtkWriter); // Add model specific output fields + + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); vtkWriter.addFaceField(analyticalSolVectors.getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact"); + vtkWriter.write(0.0); // use the staggered FV assembler diff --git a/test/freeflow/navierstokes/donea/problem.hh b/test/freeflow/navierstokes/donea/problem.hh index 7d8fb0ce0efb7194453ee97635eff2e4656e5957..73631c676b7b431dd54a0a9b4948ae84ad2663f9 100644 --- a/test/freeflow/navierstokes/donea/problem.hh +++ b/test/freeflow/navierstokes/donea/problem.hh @@ -50,7 +50,6 @@ class DoneaTestProblem : public NavierStokesProblem<TypeTag> using BoundaryTypes = Dumux::NavierStokesBoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>; using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using NumEqVector = Dumux::NumEqVector<PrimaryVariables>; using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>; @@ -61,6 +60,8 @@ class DoneaTestProblem : public NavierStokesProblem<TypeTag> using GlobalPosition = typename Element::Geometry::GlobalCoordinate; public: + using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; + DoneaTestProblem(std::shared_ptr<const GridGeometry> gridGeometry) : ParentType(gridGeometry) { @@ -171,8 +172,9 @@ public: * \brief Return the analytical solution of the problem at a given position * * \param globalPos The global position + * \param time A parameter for consistent signatures. It is ignored here as this is a stationary test */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const + PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar time = 0.0) const { Scalar x = globalPos[0]; Scalar y = globalPos[1]; diff --git a/test/freeflow/navierstokes/kovasznay/main.cc b/test/freeflow/navierstokes/kovasznay/main.cc index 0b3f316cb1fa8771791645564df133555cb8f9e1..9e2659443a7696d171e7890d996464d5f1615d69 100644 --- a/test/freeflow/navierstokes/kovasznay/main.cc +++ b/test/freeflow/navierstokes/kovasznay/main.cc @@ -46,7 +46,7 @@ #include "properties.hh" -#include "../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> int main(int argc, char** argv) { @@ -123,22 +123,16 @@ int main(int argc, char** argv) auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); gridVariables->init(x); - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - analyticalSolVectors.update(); - // intialize the vtk output module using IOFields = GetPropType<TypeTag, Properties::IOFields>; StaggeredVtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); IOFields::initOutputModule(vtkWriter); // Add model specific output fields + + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); vtkWriter.addFaceField(analyticalSolVectors.getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact"); + vtkWriter.write(0.0); // the assembler with time loop for instationary problem diff --git a/test/freeflow/navierstokes/kovasznay/problem.hh b/test/freeflow/navierstokes/kovasznay/problem.hh index a8bc9b170c425bac07b4d77dbac459ee79f7511c..a6bdb73ec2a93cc9cda534ffd4bff04b4eb40e58 100644 --- a/test/freeflow/navierstokes/kovasznay/problem.hh +++ b/test/freeflow/navierstokes/kovasznay/problem.hh @@ -53,7 +53,6 @@ class KovasznayTestProblem : public NavierStokesProblem<TypeTag> using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; using FVElementGeometry = typename GridGeometry::LocalView; using SubControlVolume = typename GridGeometry::SubControlVolume; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>; using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using NumEqVector = Dumux::NumEqVector<PrimaryVariables>; @@ -66,6 +65,8 @@ class KovasznayTestProblem : public NavierStokesProblem<TypeTag> static constexpr auto upwindSchemeOrder = getPropValue<TypeTag, Properties::UpwindSchemeOrder>(); public: + using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; + KovasznayTestProblem(std::shared_ptr<const GridGeometry> gridGeometry) : ParentType(gridGeometry) { @@ -174,8 +175,9 @@ public: * \brief Returns the analytical solution of the problem at a given position. * * \param globalPos The global position + * \param time A parameter for consistent signatures. It is ignored here as this is a stationary test */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const + PrimaryVariables analyticalSolution(const GlobalPosition& globalPos, Scalar time = 0.0) const { Scalar x = globalPos[0]; Scalar y = globalPos[1]; diff --git a/test/freeflow/navierstokes/l2error.hh b/test/freeflow/navierstokes/l2error.hh index 9a5f866ebd7c73887f82bfb9019ce0410b26fb97..631df5a88f86f23ad5300e0df90290ad3c8a21e1 100644 --- a/test/freeflow/navierstokes/l2error.hh +++ b/test/freeflow/navierstokes/l2error.hh @@ -48,7 +48,7 @@ public: * \param curSol Vector containing the current solution */ template<class Problem, class SolutionVector> - static auto calculateL2Error(const Problem& problem, const SolutionVector& curSol) + static auto calculateL2Error(const Problem& problem, const SolutionVector& curSol, Scalar time = 0.0) { using GridGeometry = std::decay_t<decltype(problem.gridGeometry())>; using Extrusion = Extrusion_t<GridGeometry>; @@ -72,7 +72,7 @@ public: // treat cell-center dofs const auto dofIdxCellCenter = scv.dofIndex(); const auto& posCellCenter = scv.dofPosition(); - const auto analyticalSolutionCellCenter = problem.analyticalSolution(posCellCenter)[Indices::pressureIdx]; + const auto analyticalSolutionCellCenter = problem.analyticalSolution(posCellCenter, time)[Indices::pressureIdx]; const auto numericalSolutionCellCenter = curSol[GridGeometry::cellCenterIdx()][dofIdxCellCenter][Indices::pressureIdx - ModelTraits::dim()]; sumError[Indices::pressureIdx] += squaredDiff_(analyticalSolutionCellCenter, numericalSolutionCellCenter) * Extrusion::volume(scv); sumReference[Indices::pressureIdx] += analyticalSolutionCellCenter * analyticalSolutionCellCenter * Extrusion::volume(scv); @@ -83,7 +83,7 @@ public: { const int dofIdxFace = scvf.dofIndex(); const int dirIdx = scvf.directionIndex(); - const auto analyticalSolutionFace = problem.analyticalSolution(scvf.center())[Indices::velocity(dirIdx)]; + const auto analyticalSolutionFace = problem.analyticalSolution(scvf.center(), time)[Indices::velocity(dirIdx)]; const auto numericalSolutionFace = curSol[GridGeometry::faceIdx()][dofIdxFace][0]; directionIndex[dofIdxFace] = dirIdx; errorVelocity[dofIdxFace] = squaredDiff_(analyticalSolutionFace, numericalSolutionFace); diff --git a/test/freeflow/navierstokes/sincos/main.cc b/test/freeflow/navierstokes/sincos/main.cc index c965dc6c2550adef9a2151d24ca9c85135fe5371..230d6c5393e623570d471c544115e4d00fd0282b 100644 --- a/test/freeflow/navierstokes/sincos/main.cc +++ b/test/freeflow/navierstokes/sincos/main.cc @@ -48,7 +48,7 @@ #include "properties.hh" -#include "../analyticalsolutionvectors.hh" +#include <test/freeflow/navierstokes/analyticalsolutionvectors.hh> /*! * \brief Creates analytical solution. @@ -89,18 +89,17 @@ auto createSource(const Problem& problem) return source; } -template<class Problem, class SolutionVector, class GridGeometry> -void printL2Error(const Problem& problem, const SolutionVector& x, const GridGeometry& gridGeometry) +template<class Problem, class SolutionVector, class GridGeometry, class Scalar = double> +void printL2Error(const Problem& problem, const SolutionVector& x, const GridGeometry& gridGeometry, Scalar time = 0.0) { using namespace Dumux; - using Scalar = double; using TypeTag = Properties::TTag::SincosTest; using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>; using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; using L2Error = NavierStokesTestL2Error<Scalar, ModelTraits, PrimaryVariables>; - const auto l2error = L2Error::calculateL2Error(*problem, x); + const auto l2error = L2Error::calculateL2Error(*problem, x, time); const int numCellCenterDofs = gridGeometry->numCellCenterDofs(); const int numFaceDofs = gridGeometry->numFaceDofs(); std::cout << std::setprecision(8) << "** L2 error (abs/rel) for " @@ -189,26 +188,14 @@ int main(int argc, char** argv) vtkWriter.addField(sourceX, "sourceX"); vtkWriter.addField(sourceY, "sourceY"); - const bool isStationary = getParam<bool>("Problem.IsStationary"); - - // create analytical solution vectors - using Scalar = GetPropType<TypeTag, Properties::Scalar>; - using FaceSolutionVector = GetPropType<TypeTag, Properties::FaceSolutionVector>; - using CellCenterSolutionVector = GetPropType<TypeTag, Properties::CellCenterSolutionVector>; - using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices; - - Dumux::NavierStokesAnalyticalSolutionVectors<Scalar, FaceSolutionVector, CellCenterSolutionVector, GridGeometry, Problem, Indices> analyticalSolVectors(problem); - if (isStationary) - analyticalSolVectors.update(); - else - analyticalSolVectors.update(0.0); - + Dumux::NavierStokesAnalyticalSolutionVectors<Problem> analyticalSolVectors(problem, 0.0); vtkWriter.addField(analyticalSolVectors.getAnalyticalPressureSolution(), "pressureExact"); vtkWriter.addField(analyticalSolVectors.getAnalyticalVelocitySolution(), "velocityExact"); vtkWriter.addFaceField(analyticalSolVectors.getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact"); vtkWriter.write(0.0); + const bool isStationary = getParam<bool>("Problem.IsStationary"); // the assembler with time loop for instationary problem using Assembler = StaggeredFVAssembler<TypeTag, DiffMethod::numeric>; @@ -261,7 +248,7 @@ int main(int argc, char** argv) if (shouldPrintL2Error) { - printL2Error(problem, x, gridGeometry); + printL2Error(problem, x, gridGeometry, timeLoop->time()+timeLoop->timeStepSize()); } // advance to the time loop to the next step diff --git a/test/freeflow/navierstokes/sincos/problem.hh b/test/freeflow/navierstokes/sincos/problem.hh index 3c43fbea8fe35cf54d818480ea984c5642445116..f78eb27e1732f44539f6806197ff73a3ef56baa7 100644 --- a/test/freeflow/navierstokes/sincos/problem.hh +++ b/test/freeflow/navierstokes/sincos/problem.hh @@ -185,16 +185,6 @@ public: return values; } - /*! - * \brief Returns the analytical solution of the problem at a given position. - * - * \param globalPos The global position - */ - PrimaryVariables analyticalSolution(const GlobalPosition& globalPos) const - { - return analyticalSolution(globalPos, time_+timeStepSize_); - } - // \} /*!