Skip to content
Snippets Groups Projects
Commit 2788e058 authored by Michael Sinsbeck's avatar Michael Sinsbeck
Browse files

Copied test/decoupled/2p to test/decoupled/2padaptive as a trunc for extension...

Copied test/decoupled/2p to test/decoupled/2padaptive as a trunc for extension to grid adaption. Not modified yet

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@6551 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent 0cca4c10
No related branches found
No related tags found
No related merge requests found
...@@ -73,6 +73,7 @@ AC_CONFIG_FILES([dumux.pc ...@@ -73,6 +73,7 @@ AC_CONFIG_FILES([dumux.pc
test/decoupled/Makefile test/decoupled/Makefile
test/decoupled/1p/Makefile test/decoupled/1p/Makefile
test/decoupled/2p/Makefile test/decoupled/2p/Makefile
test/decoupled/2padaptive/Makefile
test/decoupled/2p2c/Makefile test/decoupled/2p2c/Makefile
test/material/Makefile test/material/Makefile
test/material/tabulation/Makefile test/material/tabulation/Makefile
......
add_definitions(-DYASPGRID -DGRIDDIM=2 -DENABLE_UG)
# add build targets
ADD_EXECUTABLE("test_impes_adaptive" test_impes.cc)
TARGET_LINK_LIBRARIES("test_impes_adaptive" ${DumuxLinkLibraries})
# add required libraries and includes to the build flags
LINK_DIRECTORIES(${DumuxLinkDirectories})
INCLUDE_DIRECTORIES(${DumuxIncludeDirectories})
# make sure the grids are present in the build directory
add_custom_command(TARGET "test_transport"
POST_BUILD
COMMAND ${CMAKE_COMMAND} -E
copy_directory
"${CMAKE_CURRENT_SOURCE_DIR}/grids"
"${CMAKE_CURRENT_BINARY_DIR}/grids")
check_PROGRAMS = test_impes_adaptive
noinst_HEADERS = *.hh
EXTRA_DIST = grids/*.dgf CMakeLists.txt
test_impes_adaptive_SOURCES = test_impes_adaptive.cc
include $(top_srcdir)/am/global-rules
DGF
Interval
0 0 % first corner
1 1 % second corner
10 10 % cells in x and y direction
#
BOUNDARYDOMAIN
default 1 % all boundaries have id 1
#BOUNDARYDOMAIN
# unitcube.dgf
/*****************************************************************************
* Copyright (C) 20010 by Markus Wolff *
* Copyright (C) 2007-2008 by Bernd Flemisch *
* Copyright (C) 2008-2009 by Andreas Lauser *
* Institute of Hydraulic Engineering *
* University of Stuttgart, Germany *
* email: <givenname>.<name>@iws.uni-stuttgart.de *
* *
* 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
*
* \ingroup IMPETtests
* \brief test for the sequential 2p model
*/
#include "config.h"
#include "test_impes_adaptive_problem.hh"
#include <dune/grid/common/gridinfo.hh>
#include <dune/common/exceptions.hh>
#include <dune/common/mpihelper.hh>
#include <iostream>
#include <boost/format.hpp>
////////////////////////
// the main function
////////////////////////
void usage(const char *progname)
{
std::cout << boost::format("usage: %s [--restart restartTime] tEnd\n")%progname;
exit(1);
}
int main(int argc, char** argv)
{
try {
typedef TTAG(IMPESTestProblem) TypeTag;
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
typedef GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
typedef Dune::FieldVector<Scalar, Grid::dimensionworld> GlobalPosition;
static const int dim = Grid::dimension;
// initialize MPI, finalize is done automatically on exit
Dune::MPIHelper::instance(argc, argv);
////////////////////////////////////////////////////////////
// parse the command line arguments
////////////////////////////////////////////////////////////
if (argc < 2)
usage(argv[0]);
// deal with the restart stuff
int argPos = 1;
bool restart = false;
double restartTime = 0;
if (std::string("--restart") == argv[argPos]) {
restart = true;
++argPos;
std::istringstream(argv[argPos++]) >> restartTime;
}
if (argc - argPos != 1) {
usage(argv[0]);
}
// read the initial time step and the end time
double tEnd, dt;
std::istringstream(argv[argPos++]) >> tEnd;
dt = tEnd;
////////////////////////////////////////////////////////////
// create the grid
////////////////////////////////////////////////////////////
Dune::FieldVector<int,dim> numCells;
numCells[0] = 30;
numCells[1] = 6;
Dune::FieldVector<double,dim> lowerLeft;
lowerLeft[0] = 0;
lowerLeft[1] = 0;
Dune::FieldVector<double,dim> upperRight;
upperRight[0] = 300;
upperRight[1] = 60;
Dune::FieldVector<bool,dim> periodic(false);
Grid grid(Dune::MPIHelper::getCommunicator(), upperRight, numCells, periodic, /*overlap=*/1);
grid.loadBalance();
////////////////////////////////////////////////////////////
// instantiate and run the concrete problem
////////////////////////////////////////////////////////////
TimeManager timeManager;
Problem problem(timeManager, grid.leafView());
// load restart file if necessarry
if (restart)
problem.deserialize(restartTime);
timeManager.init(problem, 0, dt, tEnd, !restart);
timeManager.run();
return 0;
}
catch (Dune::Exception &e) {
std::cerr << "Dune reported error: " << e << std::endl;
}
catch (...) {
std::cerr << "Unknown exception thrown!\n";
throw;
}
return 3;
}
/*****************************************************************************
* Copyright (C) 2007-2008 by Klaus Mosthaf *
* Copyright (C) 2007-2008 by Bernd Flemisch *
* Copyright (C) 2008-2009 by Andreas Lauser *
* Institute of Hydraulic Engineering *
* University of Stuttgart, Germany *
* email: <givenname>.<name>@iws.uni-stuttgart.de *
* *
* 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 problem for the sequential 2p model
*/
#ifndef DUMUX_TEST_IMPES_PROBLEM_HH
#define DUMUX_TEST_IMPES_PROBLEM_HH
#if HAVE_UG
#include <dune/grid/uggrid.hh>
#endif
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/sgrid.hh>
#include <dumux/material/fluidsystems/liquidphase.hh>
#include <dumux/material/components/simpleh2o.hh>
#include <dumux/material/components/oil.hh>
#include <dumux/decoupled/2p/impes/impesproblem2p.hh>
#include <dumux/decoupled/2p/diffusion/fv/fvvelocity2p.hh>
#include <dumux/decoupled/2p/transport/fv/fvsaturation2p.hh>
#include <dumux/decoupled/2p/transport/fv/capillarydiffusion.hh>
#include <dumux/decoupled/2p/transport/fv/gravitypart.hh>
#include "test_impes_adaptive_spatialparams.hh"
#include<dumux/decoupled/2p/transport/fv/evalcflflux_coats.hh>
namespace Dumux
{
template<class TypeTag>
class TestIMPESProblem;
//////////
// Specify the properties
//////////
namespace Properties
{
NEW_TYPE_TAG(IMPESTestProblem, INHERITS_FROM(DecoupledTwoP, Transport, TestIMPESSpatialParams));
// Set the grid type
SET_PROP(IMPESTestProblem, Grid)
{
typedef Dune::YaspGrid<2> type;
//typedef Dune::SGrid<2, 2> type;
};
// Set the problem property
SET_TYPE_PROP(IMPESTestProblem, Problem, Dumux::TestIMPESProblem<TTAG(IMPESTestProblem)>);
// Set the model properties
SET_TYPE_PROP(IMPESTestProblem, TransportModel, Dumux::FVSaturation2P<TTAG(IMPESTestProblem)>);
SET_TYPE_PROP(IMPESTestProblem, DiffusivePart, Dumux::CapillaryDiffusion<TypeTag>);
SET_TYPE_PROP(IMPESTestProblem, ConvectivePart, Dumux::GravityPart<TypeTag>);
SET_PROP(IMPESTestProblem, PressureModel)
{
typedef Dumux::FVVelocity2P<TTAG(IMPESTestProblem)> type;
};
//SET_INT_PROP(IMPESTestProblem, Formulation,
// DecoupledTwoPCommonIndices::pnSn);
// Set the wetting phase
SET_PROP(IMPESTestProblem, WettingPhase)
{
private:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
public:
typedef Dumux::LiquidPhase<Scalar, Dumux::SimpleH2O<Scalar> > type;
};
// Set the non-wetting phase
SET_PROP(IMPESTestProblem, NonwettingPhase)
{
private:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
public:
typedef Dumux::LiquidPhase<Scalar, Dumux::SimpleH2O<Scalar> > type;
};
// Enable gravity
SET_BOOL_PROP(IMPESTestProblem, EnableGravity, false);
SET_TYPE_PROP(IMPESTestProblem, EvalCflFluxFunction, Dumux::EvalCflFluxCoats<TypeTag>);
SET_SCALAR_PROP(IMPESTestProblem, CFLFactor, 0.95);
}
/*!
* \ingroup IMPETtests
*
* \brief test problem for the sequential 2p model
*
* Water is injected from the left side into a rectangular 2D domain also
* filled with water. Upper and lower boundary is closed (Neumann = 0),
* and there is free outflow on the right side.
*
* To run the simulation execute the following line in shell:
* <tt>./test_impes 1e8</tt>,
* where the argument defines the simulation endtime.
*/
template<class TypeTag = TTAG(IMPESTestProblem)>
class TestIMPESProblem: public IMPESProblem2P<TypeTag>
{
typedef IMPESProblem2P<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TwoPIndices)) Indices;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidState)) FluidState;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(TimeManager)) TimeManager;
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
};
enum
{
wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
pWIdx = Indices::pwIdx,
SwIdx = Indices::SwIdx,
eqIdxPress = Indices::pressEqIdx,
eqIdxSat = Indices::satEqIdx
};
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef typename GridView::Traits::template Codim<0>::Entity Element;
typedef typename GridView::Intersection Intersection;
typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(BoundaryTypes)) BoundaryTypes;
typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes))::PrimaryVariables PrimaryVariables;
public:
TestIMPESProblem(TimeManager &timeManager, const GridView &gridView) :
ParentType(timeManager, gridView)
{
}
/*!
* \name Problem parameters
*/
// \{
/*!
* \brief The problem name.
*
* This is used as a prefix for files generated by the simulation.
*/
const char *name() const
{
return "test2p";
}
bool shouldWriteRestartFile() const
{
return false;
}
/*!
* \brief Returns the temperature within the domain.
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
Scalar temperature(const Element& element) const
{
return 273.15 + 10; // -> 10°C
}
// \}
//! Returns the reference pressure for evaluation of constitutive relations
Scalar referencePressure(const Element& element) const
{
return 1e5; // -> 10°C
}
void source(PrimaryVariables &values,const Element& element) const
{
values = 0;
}
/*!
* \brief Returns the type of boundary condition.
*
* BC for pressure equation can be dirichlet (pressure) or neumann (flux).
*
* BC for saturation equation can be dirichlet (saturation), neumann (flux), or outflow.
*/
void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos) const
{
if (globalPos[0] < eps_)
{
bcTypes.setAllDirichlet();
}
else if (globalPos[0] > this->bboxMax()[0] - eps_)
{
bcTypes.setNeumann(eqIdxPress);
bcTypes.setOutflow(eqIdxSat);
}
// all other boundaries
else
{
bcTypes.setAllNeumann();
}
}
//! set dirichlet condition (pressure [Pa], saturation [-])
void dirichletAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
values = 0;
if (globalPos[0] < eps_)
{
if (GET_PARAM(TypeTag, bool, EnableGravity))
{
Scalar pRef = referencePressureAtPos(globalPos);
Scalar temp = temperatureAtPos(globalPos);
Scalar sat = 1;
FluidState fluidState;
fluidState.update(sat, pRef, pRef, temp);
values[pWIdx] = (2e5 + (this->bboxMax()[dim-1] - globalPos[dim-1]) * FluidSystem::phaseDensity(wPhaseIdx, temp, pRef, fluidState) * this->gravity().two_norm());
}
else
{
values[pWIdx] = 2e5;
}
values[SwIdx] = 0.8;
}
else
{
values[pWIdx] = 2e5;
values[SwIdx] = 0.2;
}
}
//! set neumann condition for phases (flux, [kg/(m^2 s)])
void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) const
{
values = 0;
if (globalPos[0] > this->bboxMax()[0] - eps_)
{
values[nPhaseIdx] = 3e-4;
}
}
//! return initial solution -> only saturation values have to be given!
void initial(PrimaryVariables &values,
const Element& element) const
{
values[pWIdx] = 0;
values[SwIdx] = 0.2;
}
private:
static constexpr Scalar eps_ = 1e-6;
};
} //end namespace
#endif
/*****************************************************************************
* Copyright (C) 2008-2009 by Markus Wolff *
* Institute of Hydraulic Engineering *
* University of Stuttgart, Germany *
* email: <givenname>.<name>@iws.uni-stuttgart.de *
* *
* 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 spatial parameters for the sequential 2p test
*/
#ifndef TEST_IMPES_SPATIALPARAMETERS_HH
#define TEST_IMPES_SPATIALPARAMETERS_HH
#include <dumux/material/spatialparameters/fvspatialparameters.hh>
#include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
namespace Dumux
{
//forward declaration
template<class TypeTag>
class TestIMPESSpatialParams;
namespace Properties
{
// The spatial parameters TypeTag
NEW_TYPE_TAG(TestIMPESSpatialParams);
// Set the spatial parameters
SET_TYPE_PROP(TestIMPESSpatialParams, SpatialParameters, Dumux::TestIMPESSpatialParams<TypeTag>);
// Set the material law
SET_PROP(TestIMPESSpatialParams, MaterialLaw)
{
private:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef RegularizedBrooksCorey<Scalar> RawMaterialLaw;
public:
typedef EffToAbsLaw<RawMaterialLaw> type;
};
}
/*!
*
* \ingroup IMPETtests
* \brief spatial parameters for the sequential 2p test
*/
template<class TypeTag>
class TestIMPESSpatialParams: public FVSpatialParameters<TypeTag>
{
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef FVSpatialParameters<TypeTag> ParentType;
typedef typename Grid::ctype CoordScalar;
enum
{dim=Grid::dimension, dimWorld=Grid::dimensionworld};
typedef typename Grid::Traits::template Codim<0>::Entity Element;
typedef Dune::FieldVector<CoordScalar, dimWorld> GlobalPosition;
typedef Dune::FieldMatrix<Scalar,dim,dim> FieldMatrix;
public:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(MaterialLaw)) MaterialLaw;
typedef typename MaterialLaw::Params MaterialLawParams;
Scalar intrinsicPermeability (const Element& element) const
{
return 1.0e-7;
}
double porosity(const Element& element) const
{
return 0.2;
}
// return the parameter object for the Brooks-Corey material law which depends on the position
// const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
const MaterialLawParams& materialLawParams(const Element& element) const
{
return materialLawParams_;
}
TestIMPESSpatialParams(const GridView& gridView)
: ParentType(gridView)
{
// residual saturations
materialLawParams_.setSwr(0.2);
materialLawParams_.setSnr(0.2);
// // parameters for the Brooks-Corey Law
// // entry pressures
materialLawParams_.setPe(0);
// // Brooks-Corey shape parameters
materialLawParams_.setLambda(2);
// parameters for the linear
// entry pressures function
// materialLawParams_.setEntryPC(0);
// materialLawParams_.setMaxPC(0);
}
private:
MaterialLawParams materialLawParams_;
};
} // end namespace
#endif
SUBDIRS = 1p \ SUBDIRS = 1p \
2p \ 2p \
2padaptive \
2p2c 2p2c
EXTRA_DIST = CMakeLists.txt EXTRA_DIST = CMakeLists.txt
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment