diff --git a/test/porousmediumflow/1p/implicit/convergence/CMakeLists.txt b/test/porousmediumflow/1p/implicit/convergence/CMakeLists.txt index 578ab409c81a46eeba1347cfa4e7bf7fc6f84dda..fbc093277d8ef3b376b8199f65d01938c821a903 100644 --- a/test/porousmediumflow/1p/implicit/convergence/CMakeLists.txt +++ b/test/porousmediumflow/1p/implicit/convergence/CMakeLists.txt @@ -1,41 +1,2 @@ -dune_symlink_to_source_files(FILES "params.input" "ploterrors.py") - -# executable for tpfa tests -add_executable(test_1p_convergence_tpfa EXCLUDE_FROM_ALL main.cc) -target_compile_definitions(test_1p_convergence_tpfa PUBLIC "TYPETAG=OnePIncompressibleTpfa") - -# executable for box tests -add_executable(test_1p_convergence_box EXCLUDE_FROM_ALL main.cc) -target_compile_definitions(test_1p_convergence_box PUBLIC "TYPETAG=OnePIncompressibleBox") - -# using tpfa and conforming refinement -dumux_add_test(NAME test_1p_convergence_tpfa_conforming - TARGET test_1p_convergence_tpfa - CMAKE_GUARD "( dune-functions_FOUND )" - LABELS porousmediumflow 1p - COMMAND ./test_1p_convergence_tpfa - CMD_ARGS params.input -Problem.Name test_1p_convergence_tpfa_conforming) - -# using tpfa and nonconforming refinement -dumux_add_test(NAME test_1p_convergence_tpfa_nonconforming - TARGET test_1p_convergence_tpfa - CMAKE_GUARD "( dune-functions_FOUND )" - LABELS porousmediumflow 1p - COMMAND ./test_1p_convergence_tpfa - CMD_ARGS params.input -Problem.Name test_1p_convergence_tpfa_conforming) - -# using box and conforming refinement -dumux_add_test(NAME test_1p_convergence_box_conforming - TARGET test_1p_convergence_box - CMAKE_GUARD "( dune-functions_FOUND )" - LABELS porousmediumflow 1p - COMMAND ./test_1p_convergence_box - CMD_ARGS params.input -Problem.Name test_1p_convergence_box_conforming) - -# using box and nonconforming refinement -dumux_add_test(NAME test_1p_convergence_box_nonconforming - TARGET test_1p_convergence_box - CMAKE_GUARD "( dune-functions_FOUND )" - LABELS porousmediumflow 1p - COMMAND ./test_1p_convergence_box - CMD_ARGS params.input -Problem.Name test_1p_convergence_box_conforming) +add_subdirectory(analyticsolution) +add_subdirectory(discretesolution) diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/CMakeLists.txt b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..a50498793eb468554f9b708aa9af7dfb01baf218 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/CMakeLists.txt @@ -0,0 +1,72 @@ +add_input_file_links() +dune_symlink_to_source_files(FILES "convergencetest.py") + +# if UG is found, we use UG for all tests. Otherwise, we use +# YaspGrid for the tests on unstructured grids and enforce +# the tests on unstructured grids to be skipped +add_executable(test_1p_convergence_analytic_tpfa EXCLUDE_FROM_ALL main.cc) +add_executable(test_1p_convergence_analytic_mpfa EXCLUDE_FROM_ALL main.cc) +add_executable(test_1p_convergence_analytic_box EXCLUDE_FROM_ALL main.cc) + +if (dune-uggrid_FOUND) + target_compile_definitions(test_1p_convergence_analytic_tpfa PUBLIC "TYPETAG=OnePConvergenceTpfa" PUBLIC "GRIDTYPE=Dune::UGGrid<2>") + target_compile_definitions(test_1p_convergence_analytic_mpfa PUBLIC "TYPETAG=OnePConvergenceMpfa" PUBLIC "GRIDTYPE=Dune::UGGrid<2>") + target_compile_definitions(test_1p_convergence_analytic_box PUBLIC "TYPETAG=OnePConvergenceBox" PUBLIC "GRIDTYPE=Dune::UGGrid<2>") +else() + target_compile_definitions(test_1p_convergence_analytic_tpfa PUBLIC "TYPETAG=OnePConvergenceTpfa" PUBLIC "GRIDTYPE=Dune::YaspGrid<2>") + target_compile_definitions(test_1p_convergence_analytic_mpfa PUBLIC "TYPETAG=OnePConvergenceMpfa" PUBLIC "GRIDTYPE=Dune::YaspGrid<2>") + target_compile_definitions(test_1p_convergence_analytic_box PUBLIC "TYPETAG=OnePConvergenceBox" PUBLIC "GRIDTYPE=Dune::YaspGrid<2>") +endif() + +# using tpfa and structured grid with zero off-diagonal permeability entries +dumux_add_test(NAME test_1p_convergence_analytic_tpfa_structured + TARGET test_1p_convergence_analytic_tpfa + LABELS porousmediumflow 1p + COMMAND ./convergencetest.py + CMD_ARGS test_1p_convergence_analytic_tpfa + tpfa_structured params.input -Problem.C 0.0) + +# using mpfa and structured grid +dumux_add_test(NAME test_1p_convergence_analytic_mpfa_structured + TARGET test_1p_convergence_analytic_mpfa + LABELS porousmediumflow 1p + COMMAND ./convergencetest.py + CMD_ARGS test_1p_convergence_analytic_mpfa mpfa_structured) + +# using box and structured grid +dumux_add_test(NAME test_1p_convergence_analytic_box_structured + TARGET test_1p_convergence_analytic_box + LABELS porousmediumflow 1p + COMMAND ./convergencetest.py + CMD_ARGS test_1p_convergence_analytic_box box_structured) + +# tests on ustructured grids in case UG is found +if (dune-uggrid_FOUND) + # using mpfa and unstructured grid + dumux_add_test(NAME test_1p_convergence_analytic_mpfa_unstructured + TARGET test_1p_convergence_analytic_mpfa + LABELS porousmediumflow 1p + COMMAND ./convergencetest.py + CMD_ARGS test_1p_convergence_analytic_mpfa mpfa_unstructured + params.input -Grid.File ../../incompressible/grids/randomlydistorted.dgf) + + # using box and unstructured grid + dumux_add_test(NAME test_1p_convergence_analytic_box_unstructured + TARGET test_1p_convergence_analytic_box + LABELS porousmediumflow 1p + COMMAND ./convergencetest.py + CMD_ARGS test_1p_convergence_analytic_box box_unstructured + params.input -Grid.File ../../incompressible/grids/randomlydistorted.dgf) + +# otherwise define test with the same name but enforce skip +else() + dumux_add_test(NAME test_1p_convergence_analytic_mpfa_unstructured + SOURCES main.cc + CMAKE_GUARD dune-uggrid_FOUND + LABELS porousmediumflow 1p) + + dumux_add_test(NAME test_1p_convergence_analytic_box_unstructured + SOURCES main.cc + CMAKE_GUARD dune-uggrid_FOUND + LABELS porousmediumflow 1p) +endif() diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/convergencetest.py b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/convergencetest.py new file mode 100755 index 0000000000000000000000000000000000000000..b2991031c263021f1b2b0f337bd20073ab4e4dfe --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/convergencetest.py @@ -0,0 +1,74 @@ +#!/usr/bin/env python3 + +from math import * +import subprocess +import sys +import os + +if len(sys.argv) < 3: + sys.stderr.write("Please provide the following arguments:\n" \ + " - the name of the executable\n" \ + " - the name to be used for generated output files\n" \ + " - (optional) runtime arguments to be passed to the executable\n") + sys.exit(1) + +executableName = sys.argv[1] +testName = sys.argv[2] +testArgs = [str(i) for i in sys.argv][3:] if len(sys.argv) > 3 else ['params.input'] + +# remove the old log files +if os.path.exists(testName + '.log'): + subprocess.call(['rm', testName + '.log']) + print("Removed old log file ({})!".format(testName + '.log')) + +# do the runs with different refinement +for i in [0, 1, 2, 3]: + subprocess.call(['./' + executableName] + + testArgs + + ['-Problem.Name', testName] + + ['-Grid.Refinement', str(i)]) + +def computeRates(): + # check the rates and append them to the log file + logFile = open(testName + '.log', "r+") + + errorP = [] + for line in logFile: + line = line.strip("\n") + line = line.strip("\[ConvergenceTest\]") + line = line.split() + errorP.append(float(line[2])) + + resultsP = [] + logFile.truncate(0) + logFile.write("\nConvergence rates computed for {}:\n\n".format(testName)) + logFile.write("n\terrorP\t\trateP\n") + logFile.write("-"*50 + "\n") + for i in range(len(errorP)-1): + if isnan(errorP[i]) or isinf(errorP[i]): + continue + if not ((errorP[i] < 1e-12 or errorP[i+1] < 1e-12)): + rateP = (log(errorP[i])-log(errorP[i+1]))/log(2) + message = "{}\t{:0.4e}\t{:0.4e}\n".format(i, errorP[i], rateP) + logFile.write(message) + resultsP.append(rateP) + else: + logFile.write("error: exact solution!?") + i = len(errorP)-1 + message = "{}\t{:0.4e}\n".format(i, errorP[i], "") + logFile.write(message) + logFile.close() + + subprocess.call(['cat', testName + '.log']) + + return {"p" : resultsP} + +rates = computeRates() + +def mean(numbers): + return float(sum(numbers)) / len(numbers) + +# check the rates, we expect rates around 2 +if mean(rates["p"]) < 1.8: + sys.stderr.write("*"*70 + "\n" + "The convergence rates for pressure were not close enough to 2! Test failed.\n" + "*"*70 + "\n") + sys.exit(1) diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/main.cc b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/main.cc new file mode 100644 index 0000000000000000000000000000000000000000..126598f7ae154e49902f9d26cdca3943e1d160d9 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/main.cc @@ -0,0 +1,182 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \ingroup OnePTests + * \brief Convergence test with analytic solution + */ +#include <config.h> +#include <iostream> +#include <iomanip> + +#include <dune/common/parallel/mpihelper.hh> +#include <dune/grid/io/file/dgfparser/dgfexception.hh> + +#include <dumux/nonlinear/newtonsolver.hh> +#include <dumux/linear/seqsolverbackend.hh> + +#include <dumux/common/properties.hh> +#include <dumux/common/parameters.hh> +#include <dumux/common/dumuxmessage.hh> + +#include <dumux/io/vtkoutputmodule.hh> +#include <dumux/io/grid/gridmanager.hh> + +#include <dumux/assembly/fvassembler.hh> + +#include "properties.hh" + +template<class Problem, class SolutionVector> +void printL2Error(const Problem& problem, const SolutionVector& x) +{ + using namespace Dumux; + using Scalar = double; + + Scalar l2error = 0.0; + + for (const auto& element : elements(problem.gridGeometry().gridView())) + { + auto fvGeometry = localView(problem.gridGeometry()); + fvGeometry.bindElement(element); + + for (auto&& scv : scvs(fvGeometry)) + { + const auto dofIdx = scv.dofIndex(); + const Scalar delta = x[dofIdx] - problem.analyticalSolution(scv.dofPosition())[2/*pressureIdx*/]; + l2error += scv.volume()*(delta*delta); + } + } + using std::sqrt; + l2error = sqrt(l2error); + + const auto numDofs = problem.gridGeometry().numDofs(); + std::ostream tmp(std::cout.rdbuf()); + tmp << std::setprecision(8) << "** L2 error (abs) for " + << std::setw(6) << numDofs << " dofs " + << std::scientific + << "L2 error = " << l2error + << std::endl; + + // write the norm into a log file + std::ofstream logFile; + logFile.open(problem.name() + ".log", std::ios::app); + logFile << "[ConvergenceTest] L2(p) = " << l2error << std::endl; + logFile.close(); +} + +int main(int argc, char** argv) try +{ + using namespace Dumux; + + using TypeTag = Properties::TTag::TYPETAG; + + // initialize MPI, finalize is done automatically on exit + const auto& mpiHelper = Dune::MPIHelper::instance(argc, argv); + + // print dumux start message + if (mpiHelper.rank() == 0) + DumuxMessage::print(/*firstCall=*/true); + + //////////////////////////////////////////////////////////// + // parse the command line arguments and input file + //////////////////////////////////////////////////////////// + + // parse command line arguments + Parameters::init(argc, argv); + + ////////////////////////////////////////////////////////////////////// + // try to create a grid (from the given grid file or the input file) + ///////////////////////////////////////////////////////////////////// + using Grid = GetPropType<TypeTag, Properties::Grid>; + GridManager<Grid> gridManager; + gridManager.init(); + + // we compute on the leaf grid view + const auto& leafGridView = gridManager.grid().leafGridView(); + + // create the finite volume grid geometry + using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; + auto gridGeometry = std::make_shared<GridGeometry>(leafGridView); + gridGeometry->update(); + + // the problem (boundary conditions) + using Problem = GetPropType<TypeTag, Properties::Problem>; + auto problem = std::make_shared<Problem>(gridGeometry); + + // the solution vector + using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>; + SolutionVector x(gridGeometry->numDofs()); + + // the grid variables + using GridVariables = GetPropType<TypeTag, Properties::GridVariables>; + auto gridVariables = std::make_shared<GridVariables>(problem, gridGeometry); + gridVariables->init(x); + + // intialize the vtk output module + VtkOutputModule<GridVariables, SolutionVector> vtkWriter(*gridVariables, x, problem->name()); + using VelocityOutput = GetPropType<TypeTag, Properties::VelocityOutput>; + vtkWriter.addVelocityOutput(std::make_shared<VelocityOutput>(*gridVariables)); + using IOFields = GetPropType<TypeTag, Properties::IOFields>; + IOFields::initOutputModule(vtkWriter); // Add model specific output fields + vtkWriter.write(0.0); + + // create assembler & linear solver + using Assembler = FVAssembler<TypeTag, DiffMethod::numeric>; + auto assembler = std::make_shared<Assembler>(problem, gridGeometry, gridVariables); + + using LinearSolver = ILU0BiCGSTABBackend; + auto linearSolver = std::make_shared<LinearSolver>(); + + // the non-linear solver + using NewtonSolver = Dumux::NewtonSolver<Assembler, LinearSolver>; + NewtonSolver nonLinearSolver(assembler, linearSolver); + + // linearize & solve + nonLinearSolver.solve(x); + + // output result to vtk + vtkWriter.write(1.0); + + printL2Error(*problem, x); + + if (mpiHelper.rank() == 0) + Parameters::print(); + + return 0; + +} +catch (Dumux::ParameterException &e) +{ + std::cerr << std::endl << e << " ---> Abort!" << std::endl; + return 1; +} +catch (Dune::DGFException & e) +{ + std::cerr << "DGF exception thrown (" << e << + "). Most likely, the DGF file name is wrong " + "or the DGF file is corrupted, " + "e.g. missing hash at end of file or wrong number (dimensions) of entries." + << " ---> Abort!" << std::endl; + return 2; +} +catch (Dune::Exception &e) +{ + std::cerr << "Dune reported error: " << e << " ---> Abort!" << std::endl; + return 3; +} diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/params.input b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/params.input new file mode 100644 index 0000000000000000000000000000000000000000..3e7b5c20e6908e6476d71fb45737b0920f109297 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/params.input @@ -0,0 +1,13 @@ +[Grid] +UpperRight = 1 1 +Cells = 10 10 +Refinement = 0 + +[Problem] +EnableGravity = false +Name = convergence +C = 0.9 + +[Component] +LiquidDensity = 1.0 +LiquidKinematicViscosity = 1.0 diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/problem.hh b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/problem.hh new file mode 100644 index 0000000000000000000000000000000000000000..a9b6ce136dca3781828987b674b452e11c16f0a1 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/problem.hh @@ -0,0 +1,185 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \ingroup OnePTests + * \brief The problem setup for the convergence test with analytic solution + */ +#ifndef DUMUX_CONVERGENCE_TEST_ONEP_PROBLEM_HH +#define DUMUX_CONVERGENCE_TEST_ONEP_PROBLEM_HH + +#include <cmath> +#include <dune/common/fvector.hh> + +#include <dumux/common/properties.hh> +#include <dumux/common/boundarytypes.hh> +#include <dumux/porousmediumflow/problem.hh> + +namespace Dumux { + +/*! + * \ingroup OnePTests + * \brief The problem setup for the convergence test with analytic solution + */ +template <class TypeTag> +class ConvergenceProblem : public PorousMediumFlowProblem<TypeTag> +{ + using ParentType = PorousMediumFlowProblem<TypeTag>; + using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView; + using Scalar = GetPropType<TypeTag, Properties::Scalar>; + using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>; + using NumEqVector = GetPropType<TypeTag, Properties::NumEqVector>; + using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>; + using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; + using Element = typename GridView::template Codim<0>::Entity; + using GlobalPosition = typename Element::Geometry::GlobalCoordinate; + + static constexpr auto velocityXIdx = 0; + static constexpr auto velocityYIdx = 1; + static constexpr auto pressureIdx = 2; + +public: + /*! + * \brief The constructor. + * \param gridGeometry The finite-volume grid geometry + */ + ConvergenceProblem(std::shared_ptr<const GridGeometry> gridGeometry) + : ParentType(gridGeometry) + , c_(getParam<Scalar>("Problem.C")) + {} + + /*! + * \name Problem parameters + */ + // \{ + + /*! + * \brief Returns the temperature within the domain in [K]. + */ + Scalar temperature() const + { return 273.15 + 10; } // 10°C + // \} + + /*! + * \name Boundary conditions + */ + // \{ + + /*! + * \brief Specifies which kind of boundary condition should be + * used for which equation on a given boundary control volume. + * + * \param globalPos The position of the center of the finite volume + */ + BoundaryTypes boundaryTypesAtPos(const GlobalPosition& globalPos) const + { + BoundaryTypes values; + values.setAllDirichlet(); + return values; + } + + /*! + * \brief Evaluates Dirichlet boundary conditions. + * \param globalPos The center of the finite volume which ought to be set. + */ + PrimaryVariables dirichletAtPos(const GlobalPosition& globalPos) const + { + const auto p = analyticalSolution(globalPos)[pressureIdx]; + return PrimaryVariables(p); + } + + // \} + + /*! + * \name Volume terms + */ + // \{ + + /*! + * \brief Evaluates the source term for all phases within a given + * sub-control volume. + * + * For this method, the \a priVars parameter stores the rate mass + * of a component is generated or annihilated per volume + * unit. Positive values mean that mass is created, negative ones + * mean that it vanishes. + * + * The units must be according to either using mole or mass fractions (mole/(m^3*s) or kg/(m^3*s)). + */ + NumEqVector sourceAtPos(const GlobalPosition& globalPos) const + { + const Scalar x = globalPos[0]; + const Scalar y = globalPos[1]; + using std::exp; + using std::sin; + using std::cos; + const Scalar cosOmegaX = cos(omega_*x); + static const Scalar expTwo = exp(2); + const Scalar expYPlusOne = exp(y+1); + + const Scalar result = ( -(c_*cosOmegaX + 1)*exp(y - 1) + + 1.5*c_*expYPlusOne*cosOmegaX + + omega_*omega_*(expYPlusOne - expTwo + 2)) + * sin(omega_*x); + + return NumEqVector(result); + } + + // \} + + /*! + * \brief Evaluates the initial value for a control volume. + * \param globalPos The position for which the initial condition should be evaluated + */ + PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const + { return PrimaryVariables(0.0); } + + /*! + * \brief Returns the analytical solution of the problem at a given position. + * \param globalPos The global position + */ + auto analyticalSolution(const GlobalPosition& globalPos) const + { + Dune::FieldVector<Scalar, 3> sol(0.0); + const Scalar x = globalPos[0]; + const Scalar y = globalPos[1]; + using std::exp; using std::sin; using std::cos; + const Scalar sinOmegaX = sin(omega_*x); + const Scalar cosOmegaX = cos(omega_*x); + static const Scalar expTwo = exp(2); + const Scalar expYPlusOne = exp(y+1); + + sol[pressureIdx] = (expYPlusOne + 2 - expTwo)*sinOmegaX + 10.0; + sol[velocityXIdx] = c_/(2*omega_)*expYPlusOne*sinOmegaX*sinOmegaX + -omega_*(expYPlusOne + 2 - expTwo)*cosOmegaX; + sol[velocityYIdx] = (0.5*c_*(expYPlusOne + 2 - expTwo)*cosOmegaX + -(c_*cosOmegaX + 1)*exp(y-1))*sinOmegaX; + + return sol; + } + +private: + static constexpr Scalar eps_ = 1e-7; + static constexpr Scalar omega_ = M_PI; + Scalar c_; +}; + +} // end namespace Dumux + +#endif diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/properties.hh b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/properties.hh new file mode 100644 index 0000000000000000000000000000000000000000..c72ed927486667cf3a5e6202726566c7d3185210 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/properties.hh @@ -0,0 +1,87 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \ingroup OnePTests + * \brief The properties for the convergence test with analytic solution + */ +#ifndef DUMUX_CONVERGENCE_TEST_ONEP_PROPERTIES_HH +#define DUMUX_CONVERGENCE_TEST_ONEP_PROPERTIES_HH + +#include <dune/grid/yaspgrid.hh> +#if HAVE_UG +#include <dune/grid/uggrid.hh> +#endif + +#include <dumux/discretization/cctpfa.hh> +#include <dumux/discretization/ccmpfa.hh> +#include <dumux/discretization/box.hh> + +#include <dumux/material/components/constant.hh> +#include <dumux/material/fluidsystems/1pliquid.hh> +#include <dumux/porousmediumflow/1p/model.hh> + +#include "spatialparams.hh" +#include "problem.hh" + +namespace Dumux::Properties { + +// Create new type tags +namespace TTag { +struct OnePConvergence { using InheritsFrom = std::tuple<OneP>; }; +struct OnePConvergenceTpfa { using InheritsFrom = std::tuple<OnePConvergence, CCTpfaModel>; }; +struct OnePConvergenceMpfa { using InheritsFrom = std::tuple<OnePConvergence, CCMpfaModel>; }; +struct OnePConvergenceBox { using InheritsFrom = std::tuple<OnePConvergence, BoxModel>; }; +} // end namespace TTag + +// Set the problem property +template<class TypeTag> +struct Problem<TypeTag, TTag::OnePConvergence> { using type = Dumux::ConvergenceProblem<TypeTag>; }; + +// the fluid system +template<class TypeTag> +struct FluidSystem<TypeTag, TTag::OnePConvergence> +{ + using Scalar = GetPropType<TypeTag, Properties::Scalar>; + using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::Constant<1, Scalar> > ; +}; + +// Set the grid type +template<class TypeTag> +struct Grid<TypeTag, TTag::OnePConvergence> { using type = GRIDTYPE; }; + +template<class TypeTag> +struct SpatialParams<TypeTag, TTag::OnePConvergence> +{ + using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>; + using Scalar = GetPropType<TypeTag, Properties::Scalar>; + using type = ConvergenceTestSpatialParams<GridGeometry, Scalar>; +}; + +// Enable caching +template<class TypeTag> +struct EnableGridVolumeVariablesCache<TypeTag, TTag::OnePConvergence> { static constexpr bool value = true; }; +template<class TypeTag> +struct EnableGridFluxVariablesCache<TypeTag, TTag::OnePConvergence> { static constexpr bool value = true; }; +template<class TypeTag> +struct EnableGridGeometryCache<TypeTag, TTag::OnePConvergence> { static constexpr bool value = true; }; + +} // end namespace Dumux::Properties + +#endif diff --git a/test/porousmediumflow/1p/implicit/convergence/analyticsolution/spatialparams.hh b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/spatialparams.hh new file mode 100644 index 0000000000000000000000000000000000000000..d57ca0acdf8a954fb30f7512d9641a45f5d8cf5d --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/analyticsolution/spatialparams.hh @@ -0,0 +1,103 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 3 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \ingroup OnePTests + * \brief The spatial params of the incompressible single-phase convergence test with analytic solution + */ + +#ifndef DUMUX_CONVERGENCE_TEST_ONEP_SPATIALPARAMS_HH +#define DUMUX_CONVERGENCE_TEST_ONEP_SPATIALPARAMS_HH + +#include <cmath> +#include <dumux/common/math.hh> +#include <dumux/common/parameters.hh> +#include <dune/common/fmatrix.hh> +#include <dumux/material/spatialparams/fv1p.hh> + +namespace Dumux { + +/*! + * \ingroup OnePTests + * \brief The spatial params of the incompressible single-phase convergence test with analytic solution + */ +template<class GridGeometry, class Scalar> +class ConvergenceTestSpatialParams +: public FVSpatialParamsOneP<GridGeometry, Scalar, + ConvergenceTestSpatialParams<GridGeometry, Scalar>> +{ + using GridView = typename GridGeometry::GridView; + using ParentType = FVSpatialParamsOneP<GridGeometry, Scalar, + ConvergenceTestSpatialParams<GridGeometry, Scalar>>; + + using Element = typename GridView::template Codim<0>::Entity; + using GlobalPosition = typename Element::Geometry::GlobalCoordinate; + + static constexpr int dimWorld = GridView::dimensionworld; + using DimWorldMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>; + + +public: + // export permeability type + using PermeabilityType = DimWorldMatrix; + + ConvergenceTestSpatialParams(std::shared_ptr<const GridGeometry> gridGeometry) + : ParentType(gridGeometry) + { + c_ = getParam<Scalar>("Problem.C"); + } + + /*! + * \brief Defines the intrinsic permeability \f$\mathrm{[m^2]}\f$. + * + * \param globalPos The global position where we evaluate + */ + PermeabilityType permeabilityAtPos(const GlobalPosition& globalPos) const + { + PermeabilityType K(0.0); + + using std::cos; + using std::sin; + using std::exp; + + const Scalar x = globalPos[0]; + K[0][0] = 1.0; + K[0][1] = -c_/(2*omega_) * sin(omega_*x); + K[1][0] = K[0][1]; + K[1][1] = exp(-2)*(1 + c_*cos(omega_*x)); + + return K; + } + + /*! \brief Defines the porosity in [-]. + * + * \param globalPos The global position + */ + Scalar porosityAtPos(const GlobalPosition& globalPos) const + { return 0.4; } + +private: + static constexpr Scalar omega_ = M_PI; + Scalar permeability_; + Scalar c_; +}; + +} // end namespace Dumux + +#endif // DUMUX_CONVERGENCE_TEST_ONEP_SPATIALPARAMS_HH diff --git a/test/porousmediumflow/1p/implicit/convergence/discretesolution/CMakeLists.txt b/test/porousmediumflow/1p/implicit/convergence/discretesolution/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..65422ee2f420158a8ff0398d4b674e90327c5478 --- /dev/null +++ b/test/porousmediumflow/1p/implicit/convergence/discretesolution/CMakeLists.txt @@ -0,0 +1,41 @@ +dune_symlink_to_source_files(FILES "params.input") + +# executable for tpfa tests +add_executable(test_1p_convergence_tpfa EXCLUDE_FROM_ALL main.cc) +target_compile_definitions(test_1p_convergence_tpfa PUBLIC "TYPETAG=OnePIncompressibleTpfa") + +# executable for box tests +add_executable(test_1p_convergence_box EXCLUDE_FROM_ALL main.cc) +target_compile_definitions(test_1p_convergence_box PUBLIC "TYPETAG=OnePIncompressibleBox") + +# using tpfa and conforming refinement +dumux_add_test(NAME test_1p_convergence_tpfa_conforming + TARGET test_1p_convergence_tpfa + CMAKE_GUARD "( dune-functions_FOUND )" + LABELS porousmediumflow 1p + COMMAND ./test_1p_convergence_tpfa + CMD_ARGS params.input -Problem.Name test_1p_convergence_tpfa_conforming) + +# using tpfa and nonconforming refinement +dumux_add_test(NAME test_1p_convergence_tpfa_nonconforming + TARGET test_1p_convergence_tpfa + CMAKE_GUARD "( dune-functions_FOUND )" + LABELS porousmediumflow 1p + COMMAND ./test_1p_convergence_tpfa + CMD_ARGS params.input -Problem.Name test_1p_convergence_tpfa_conforming) + +# using box and conforming refinement +dumux_add_test(NAME test_1p_convergence_box_conforming + TARGET test_1p_convergence_box + CMAKE_GUARD "( dune-functions_FOUND )" + LABELS porousmediumflow 1p + COMMAND ./test_1p_convergence_box + CMD_ARGS params.input -Problem.Name test_1p_convergence_box_conforming) + +# using box and nonconforming refinement +dumux_add_test(NAME test_1p_convergence_box_nonconforming + TARGET test_1p_convergence_box + CMAKE_GUARD "( dune-functions_FOUND )" + LABELS porousmediumflow 1p + COMMAND ./test_1p_convergence_box + CMD_ARGS params.input -Problem.Name test_1p_convergence_box_conforming) diff --git a/test/porousmediumflow/1p/implicit/convergence/main.cc b/test/porousmediumflow/1p/implicit/convergence/discretesolution/main.cc similarity index 100% rename from test/porousmediumflow/1p/implicit/convergence/main.cc rename to test/porousmediumflow/1p/implicit/convergence/discretesolution/main.cc diff --git a/test/porousmediumflow/1p/implicit/convergence/params.input b/test/porousmediumflow/1p/implicit/convergence/discretesolution/params.input similarity index 100% rename from test/porousmediumflow/1p/implicit/convergence/params.input rename to test/porousmediumflow/1p/implicit/convergence/discretesolution/params.input diff --git a/test/porousmediumflow/1p/implicit/convergence/problem.hh b/test/porousmediumflow/1p/implicit/convergence/discretesolution/problem.hh similarity index 100% rename from test/porousmediumflow/1p/implicit/convergence/problem.hh rename to test/porousmediumflow/1p/implicit/convergence/discretesolution/problem.hh diff --git a/test/porousmediumflow/1p/implicit/convergence/solver.hh b/test/porousmediumflow/1p/implicit/convergence/discretesolution/solver.hh similarity index 100% rename from test/porousmediumflow/1p/implicit/convergence/solver.hh rename to test/porousmediumflow/1p/implicit/convergence/discretesolution/solver.hh diff --git a/test/porousmediumflow/1p/implicit/convergence/spatialparams.hh b/test/porousmediumflow/1p/implicit/convergence/discretesolution/spatialparams.hh similarity index 100% rename from test/porousmediumflow/1p/implicit/convergence/spatialparams.hh rename to test/porousmediumflow/1p/implicit/convergence/discretesolution/spatialparams.hh