[1pnc][test] add the 1p2c problem

to test the 1pnc model in next

test/porousmediumflow/1pnc/implicit/1p2ctestproblem.hh

+// -*- 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 2 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
+ * \brief Definition of a problem, for the 1p2c problem:
+ * Component transport of nitrogen dissolved in the water phase.
+ */
+#ifndef DUMUX_1P2C_TEST_PROBLEM_HH
+#define DUMUX_1P2C_TEST_PROBLEM_HH
+
+#include <dumux/implicit/box/properties.hh>
+#include <dumux/implicit/cellcentered/tpfa/properties.hh>
+#include <dumux/implicit/cellcentered/mpfa/properties.hh>
+#include <dumux/porousmediumflow/1p2c/implicit/model.hh>
+#include <dumux/porousmediumflow/implicit/problem.hh>
+
+#include <dumux/material/fluidsystems/h2on2.hh>
+#include "1p2ctestspatialparams.hh"
+
+#define NONISOTHERMAL 0
+
+namespace Dumux
+{
+
+template <class TypeTag>
+class OnePTwoCTestProblem;
+
+namespace Properties
+{
+#if NONISOTHERMAL
+NEW_TYPE_TAG(OnePTwoCOutflowProblem, INHERITS_FROM(OnePTwoCNI));
+NEW_TYPE_TAG(OnePTwoCOutflowBoxProblem, INHERITS_FROM(BoxModel, OnePTwoCOutflowProblem));
+NEW_TYPE_TAG(OnePTwoCOutflowCCProblem, INHERITS_FROM(CCTpfaModel, OnePTwoCOutflowProblem));
+#else
+NEW_TYPE_TAG(OnePTwoCTestProblem, INHERITS_FROM(OnePTwoC));
+NEW_TYPE_TAG(OnePTwoCTestBoxProblem, INHERITS_FROM(BoxModel, OnePTwoCTestProblem));
+NEW_TYPE_TAG(OnePTwoCTestCCProblem, INHERITS_FROM(CCTpfaModel, OnePTwoCTestProblem));
+NEW_TYPE_TAG(OnePTwoCTestCCMpfaProblem, INHERITS_FROM(CCMpfaModel, OnePTwoCTestProblem));
+#endif
+
+// Set the grid type
+#if HAVE_UG
+SET_TYPE_PROP(OnePTwoCTestProblem, Grid, Dune::UGGrid<2>);
+#else
+SET_TYPE_PROP(OnePTwoCTestProblem, Grid, Dune::YaspGrid<2>);
+#endif
+
+// Set the problem property
+SET_TYPE_PROP(OnePTwoCTestProblem, Problem, OnePTwoCTestProblem<TypeTag>);
+
+// Set fluid configuration
+SET_TYPE_PROP(OnePTwoCTestProblem,
+              FluidSystem,
+              FluidSystems::H2ON2<typename GET_PROP_TYPE(TypeTag, Scalar), false>);
+
+// Set the spatial parameters
+SET_TYPE_PROP(OnePTwoCTestProblem, SpatialParams, OnePTwoCTestSpatialParams<TypeTag>);
+
+// Define whether mole(true) or mass (false) fractions are used
+SET_BOOL_PROP(OnePTwoCTestProblem, UseMoles, true);
+
+// Enable velocity output
+SET_BOOL_PROP(OnePTwoCTestProblem, VtkAddVelocity, true);
+
+// Disable gravity
+SET_BOOL_PROP(OnePTwoCTestProblem, ProblemEnableGravity, false);
+}
+
+
+/*!
+ * \ingroup OnePTwoCModel
+ * \ingroup ImplicitTestProblems
+ *
+ * \brief Definition of a problem, for the 1p2c problem:
+ * Nitrogen is dissolved in the water phase and
+ * is transported with the water flow from the left side to the right.
+ *
+ * The model domain is 1 m times 1 m with a discretization length of 0.05 m
+ * and homogeneous soil properties (\f$ \mathrm{K=10e-10, \Phi=0.4, \tau=0.28}\f$).
+ * Initially the domain is filled with pure water.
+ *
+ * At the left side, a Dirichlet condition defines a nitrogen mole fraction
+ * of 0.3 mol/mol.
+ * The water phase flows from the left side to the right due to the applied pressure
+ * gradient of 1e5 Pa/m. The nitrogen is transported with the water flow
+ * and leaves the domain at the right boundary
+ * where again Dirichlet boundary conditions are applied. Here, the nitrogen mole
+ * fraction is set to 0.0 mol/mol.
+ *
+ * The model is able to use either mole or mass fractions. The property useMoles can be set to either true or false in the
+ * problem file. Make sure that the according units are used in the problem setup. The default setting for useMoles is true.
+ *
+ * This problem uses the \ref OnePTwoCModel model.
+ *
+ * To run the simulation execute the following line in shell:
+ * <tt>./test_box1p2c -parameterFile ./test_box1p2c.input</tt> or
+ * <tt>./test_cc1p2c -parameterFile ./test_cc1p2c.input</tt>
+ */
+template <class TypeTag>
+class OnePTwoCTestProblem : public ImplicitPorousMediaProblem<TypeTag>
+{
+    using ParentType = ImplicitPorousMediaProblem<TypeTag>;
+
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+    using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
+    using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+    using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
+    using TimeManager = typename GET_PROP_TYPE(TypeTag, TimeManager);
+    using BoundaryTypes = typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
+    using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
+
+    // copy some indices for convenience
+    enum
+    {
+        // indices of the primary variables
+        pressureIdx = Indices::pressureIdx,
+        massOrMoleFracIdx = Indices::massOrMoleFracIdx,
+
+        // indices for the non-isothermal case
+#if NONISOTHERMAL
+        temperatureIdx = Indices::temperatureIdx,
+        energyEqIdx = Indices::energyEqIdx,
+#endif
+
+        // indices of the equations
+        conti0EqIdx = Indices::conti0EqIdx,
+        transportEqIdx = Indices::transportEqIdx
+    };
+
+    //! property that defines whether mole or mass fractions are used
+    static const bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
+
+    static const int dimWorld = GridView::dimensionworld;
+    using GlobalPosition = Dune::FieldVector<Scalar, dimWorld>;
+
+public:
+    OnePTwoCTestProblem(TimeManager &timeManager, const GridView &gridView)
+        : ParentType(timeManager, gridView)
+    {
+        //initialize fluid system
+        FluidSystem::init();
+
+        name_ = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, std::string, Problem, Name);
+
+        // stating in the console whether mole or mass fractions are used
+        if(useMoles)
+            std::cout<<"problem uses mole fractions"<<std::endl;
+        else
+            std::cout<<"problem uses mass fractions"<<std::endl;
+    }
+
+    /*!
+     * \name Problem parameters
+     */
+    // \{
+
+    /*!
+     * \brief The problem name.
+     *
+     * This is used as a prefix for files generated by the simulation.
+     */
+    const std::string& name() const
+    {
+        return name_;
+    }
+
+#if !NONISOTHERMAL
+    /*!
+     * \brief Returns the temperature within the domain [K].
+     *
+     * This problem assumes a temperature of 20 degrees Celsius.
+     */
+    Scalar temperature() const
+    { return 273.15 + 20; } // in [K]
+
+#endif
+
+    // \}
+
+    /*!
+     * \name Boundary conditions
+     */
+    // \{
+
+    /*!
+     * \brief Specifies which kind of boundary condition should be
+     *        used for which equation on a given boundary segment.
+     *
+     * \param globalPos The position for which the bc type should be evaluated
+     */
+    BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
+    {
+        BoundaryTypes values;
+
+        if(globalPos[0] < eps_ || globalPos[0] > this->bBoxMax()[0] - eps_)
+            values.setAllDirichlet();
+        else
+            values.setAllNeumann();
+
+        return values;
+    }
+
+    /*!
+     * \brief Evaluate the boundary conditions for a dirichlet
+     *        boundary segment.
+     *
+     * \param globalPos The position for which the bc type should be evaluated
+     */
+    PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
+    {
+        PrimaryVariables values = initial_(globalPos);
+
+        // condition for the N2 molefraction at left boundary
+        if (globalPos[0] < eps_)
+        {
+            values[massOrMoleFracIdx] = 2.0e-5;
+#if NONISOTHERMAL
+            values[temperatureIdx] = 300.;
+#endif
+        }
+
+        return values;
+    }
+
+    /*!
+     * \brief Evaluate the boundary conditions for a Neumann
+     *        boundary segment.
+     *
+     * For this method, the \a priVars parameter stores the mass flux
+     * in normal direction of each component. Negative values mean
+     * influx.
+     *
+     * The units must be according to either using mole or mass fractions. (mole/(m^2*s) or kg/(m^2*s))
+     */
+    PrimaryVariables neumannAtPos(const GlobalPosition& globalPos) const
+    { return PrimaryVariables(0.0); }
+
+    // \}
+
+    /*!
+     * \name Volume terms
+     */
+    // \{
+
+    /*!
+     * \brief Evaluate 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 annihilate 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))
+     */
+    PrimaryVariables sourceAtPos(const GlobalPosition &globalPos) const
+    { return PrimaryVariables(0.0); }
+
+    /*!
+     * \brief Evaluate the initial value for a control volume.
+     *
+     * \param globalPos The position for which the initial condition should be evaluated
+     *
+     * For this method, the \a values parameter stores primary
+     * variables.
+     */
+    PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
+    { return initial_(globalPos); }
+
+    // \}
+
+private:
+    // the internal method for the initial condition
+    PrimaryVariables initial_(const GlobalPosition &globalPos) const
+    {
+        PrimaryVariables priVars;
+        priVars[pressureIdx] = 2e5 - 1e5*globalPos[0]; // initial condition for the pressure
+        priVars[massOrMoleFracIdx] = 0.0;  // initial condition for the N2 molefraction
+#if NONISOTHERMAL
+        priVars[temperatureIdx] = 290.;
+#endif
+        return priVars;
+    }
+
+    static constexpr Scalar eps_ = 1e-6;
+    std::string name_;
+};
+
+} //end namespace Dumux
+
+#endif

test/porousmediumflow/1pnc/implicit/1p2ctestspatialparams.hh

+// -*- 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 2 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
+ *
+ * \brief Definition of the spatial parameters for the 1p2c
+ *        outlfow problem.
+ */
+#ifndef DUMUX_1P2C_TEST_SPATIAL_PARAMS_HH
+#define DUMUX_1P2C_TEST_SPATIAL_PARAMS_HH
+
+#include <dumux/material/spatialparams/implicit1p.hh>
+
+namespace Dumux
+{
+
+/*!
+ * \ingroup OnePTwoCModel
+ * \ingroup ImplicitTestProblems
+ *
+ * \brief Definition of the spatial parameters for the 1p2c
+ *        outflow problem.
+ */
+template<class TypeTag>
+class OnePTwoCTestSpatialParams : public ImplicitSpatialParamsOneP<TypeTag>
+{
+    using ParentType = ImplicitSpatialParamsOneP<TypeTag>;
+    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
+    using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
+    using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
+    using Element = typename GridView::template Codim<0>::Entity;
+    using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
+    using ElementSolutionVector = typename GET_PROP_TYPE(TypeTag, ElementSolutionVector);
+
+    static const int dimWorld = GridView::dimensionworld;
+    using GlobalPosition = typename Dune::FieldVector<Scalar, dimWorld>;
+
+public:
+    // export permeability type
+    using PermeabilityType = Scalar;
+
+    OnePTwoCTestSpatialParams(const Problem& problem, const GridView &gridView)
+        : ParentType(problem, gridView)
+    {
+        permeability_ = 1e-10;
+        porosity_ = 0.4;
+
+        // heat conductivity of granite
+        lambdaSolid_ = 2.8;
+    }
+
+    /*!
+     * \brief Define the intrinsic permeability \f$\mathrm{[m^2]}\f$.
+     *
+     * \param globalPos The global position
+     */
+    PermeabilityType permeabilityAtPos(const GlobalPosition& globalPos) const
+    { return permeability_; }
+
+    /*!
+     * \brief Define the porosity \f$\mathrm{[-]}\f$.
+     *
+     * \param globalPos The global position
+     */
+    Scalar porosityAtPos(const GlobalPosition& globalPos) const
+    { return porosity_; }
+
+    /*!
+     * \brief Define the dispersivity.
+     *
+     * \param element The finite element
+     * \param scv The sub-control volume
+     * \param elemSol The solution for all dofs of the element
+     */
+    Scalar dispersivity(const Element &element,
+                        const SubControlVolume& scv,
+                        const ElementSolutionVector& elemSol) const
+    { return 0; }
+
+    /*!
+     * \brief Returns the heat capacity \f$[J / (kg K)]\f$ of the rock matrix.
+     *
+     * This is only required for non-isothermal models.
+     *
+     * \param element The element
+     * \param scv The sub control volume
+     * \param elemSol The element solution vector
+     */
+    Scalar solidHeatCapacity(const Element &element,
+                             const SubControlVolume& scv,
+                             const ElementSolutionVector& elemSol) const
+    { return 790; /*specific heat capacity of granite [J / (kg K)]*/ }
+
+    /*!
+     * \brief Returns the mass density \f$[kg / m^3]\f$ of the rock matrix.
+     *
+     * This is only required for non-isothermal models.
+     *
+     * \param element The element
+     * \param scv The sub control volume
+     * \param elemSol The element solution vector
+     */
+    Scalar solidDensity(const Element &element,
+                        const SubControlVolume& scv,
+                        const ElementSolutionVector& elemSol) const
+    { return 2700; /*density of granite [kg/m^3]*/ }
+
+    /*!
+     * \brief Returns the thermal conductivity \f$\mathrm{[W/(m K)]}\f$ of the porous material.
+     *
+     * \param element The element
+     * \param scv The sub control volume
+     * \param elemSol The element solution vector
+     */
+    Scalar solidThermalConductivity(const Element &element,
+                                    const SubControlVolume& scv,
+                                    const ElementSolutionVector& elemSol) const
+    { return lambdaSolid_; }
+
+
+
+private:
+    Scalar permeability_;
+    Scalar porosity_;
+    Scalar lambdaSolid_;
+};
+
+} // end namespace Dumux
+
+#endif

test/porousmediumflow/1pnc/implicit/CMakeLists.txt

+add_input_file_links()
+
+# isothermal tests
+add_executable(test test_box1p2c.cc)
+
+# add_executable(test_cc1pnc test_cc1pnc.cc)
+
+
+
+
+# add_dumux_test(test_box1p2c test_box1p2c test_box1p2c.cc
+#                python ${CMAKE_SOURCE_DIR}/bin/testing/runtest.py
+#                  --script fuzzy
+#                  --files ${CMAKE_SOURCE_DIR}/test/references/outflowbox-reference.vtu
+#                          ${CMAKE_CURRENT_BINARY_DIR}/outflowbox-00009.vtu
+#                  --command "${CMAKE_CURRENT_BINARY_DIR}/test_box1p2c"
+#                  --zeroThreshold {"velocity":5e-16})
+#
+# add_dumux_test(test_cc1p2c test_cc1p2c test_cc1p2c.cc
+#                python ${CMAKE_SOURCE_DIR}/bin/testing/runtest.py
+#                  --script fuzzy
+#                  --files ${CMAKE_SOURCE_DIR}/test/references/1p2ctestcc-reference.vtu
+#                          ${CMAKE_CURRENT_BINARY_DIR}/1p2ctestcc-00009.vtu
+#                  --command "${CMAKE_CURRENT_BINARY_DIR}/test_cc1p2c"
+#                  --zeroThreshold {"velocity":5e-16})
+#
+# #install sources
+# install(FILES
+#
+# 1p2cnispatialparams.hh
+# 1p2coutflowproblem.hh
+# 1p2coutflowspatialparams.hh
+# test_box1p2c.cc
+# test_box1p2cniconduction.cc
+# test_box1p2cniconvection.cc
+# test_cc1p2c.cc
+# test_cc1p2cniconduction.cc
+# test_cc1p2cniconvection.cc
+# DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/test/implicit/1pnc)

test/porousmediumflow/1pnc/implicit/test_box1p2c.cc

+// -*- 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 2 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
+ *
+ * \brief test for the 1p2c box model
+ */
+#include <config.h>
+#include "1p2ctestproblem.hh"
+#include <dumux/common/start.hh>
+
+/*!
+ * \brief Provides an interface for customizing error messages associated with
+ *        reading in parameters.
+ *
+ * \param progName  The name of the program, that was tried to be started.
+ * \param errorMsg  The error message that was issued by the start function.
+ *                  Comprises the thing that went wrong and a general help message.
+ */
+void usage(const char *progName, const std::string &errorMsg)
+{
+    if (errorMsg.size() > 0) {
+        std::string errorMessageOut = "\nUsage: ";
+                    errorMessageOut += progName;
+                    errorMessageOut += " [options]\n";
+                    errorMessageOut += errorMsg;
+                    errorMessageOut += "\n\nThe list of mandatory options for this program is:\n"
+                                        "\t-TimeManager.TEnd      End of the simulation [s] \n"
+                                        "\t-TimeManager.DtInitial Initial timestep size [s] \n"
+                                        "\t-Grid.File             Name of the file containing the grid \n"
+                                        "\t                       definition in DGF format\n";
+        std::cout << errorMessageOut
+                  << "\n";
+    }
+}
+
+int main(int argc, char** argv)
+{
+    typedef TTAG(OnePTwoCTestBoxProblem) ProblemTypeTag;
+    return Dumux::start<ProblemTypeTag>(argc, argv, usage);
+}

test/porousmediumflow/1pnc/implicit/test_box1p2c.input

+[TimeManager]
+DtInitial = 1 # [s]
+TEnd = 100 # [s]
+
+[Grid]
+LowerLeft = 0 0
+UpperRight = 1 1
+Cells = 20 2
+
+[Problem]
+Name = outflowbox # name passed to the output routines
+EnableGravity = 0 # disable gravity
+
+[Vtk]
+AddVelocity = 1 # enable velocity output

test/porousmediumflow/CMakeLists.txt

 add_subdirectory("1p")
 add_subdirectory("1p2c")
+add_subdirectory("1pnc")
 add_subdirectory("2p")
 add_subdirectory("2p1c")
 add_subdirectory("2p2c")