Commit d29f8a59 authored by Markus Wolff's avatar Markus Wolff
Browse files

added test for decoupled 1p model



git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@4166 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent dc29a0ab
# tests where program to build and program to run are equal
check_PROGRAMS = test_diffusion
check_PROGRAMS = test_diffusion test_1p
test_diffusion_SOURCES = test_diffusion.cc
test_1p_SOURCES = test_1p.cc
include $(top_srcdir)/am/global-rules
// $Id: test_diffusion.cc 4144 2010-08-24 10:10:47Z bernd $
/*****************************************************************************
* 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, as long as this copyright notice *
* is included in its original form. *
* *
* This program is distributed WITHOUT ANY WARRANTY. *
*****************************************************************************/
#include "config.h"
#include <iostream>
#include <boost/format.hpp>
#include <dune/common/exceptions.hh>
#include <dune/common/mpihelper.hh>
#include <dune/grid/common/gridinfo.hh>
#include "test_1p_problem.hh"
#include "benchmarkresult.hh"
////////////////////////
// the main function
////////////////////////
void usage(const char *progname)
{
std::cout << boost::format("usage: %s #refine [delta]\n")%progname;
exit(1);
}
int main(int argc, char** argv)
{
try {
typedef TTAG(TestProblemOneP) TypeTag;
typedef GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
typedef GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
static const int dim = Grid::dimension;
typedef Dune::FieldVector<Scalar, dim> GlobalPosition;
// initialize MPI, finalize is done automatically on exit
Dune::MPIHelper::instance(argc, argv);
////////////////////////////////////////////////////////////
// parse the command line arguments
////////////////////////////////////////////////////////////
if (argc != 2 && argc != 3)
usage(argv[0]);
int numRefine;
std::istringstream(argv[1]) >> numRefine;
double delta = 1e-3;
if (argc == 3)
std::istringstream(argv[2]) >> delta;
////////////////////////////////////////////////////////////
// create the grid
////////////////////////////////////////////////////////////
Dune::FieldVector<int,dim> N(1);
GlobalPosition L(0.0);
GlobalPosition H(1.0);
Grid grid(N,L,H);
grid.globalRefine(numRefine);
////////////////////////////////////////////////////////////
// instantiate and run the concrete problem
////////////////////////////////////////////////////////////
// Dune::Timer timer;
// bool consecutiveNumbering = true;
typedef GET_PROP_TYPE(TypeTag, PTAG(Problem)) Problem;
Problem problem(grid.leafView(), delta);
// timer.reset();
problem.init();
problem.writeOutput();
// double time = timer.elapsed();
// Dumux::ResultEvaluation result;
// result.evaluate(grid.leafView(), problem, problem.variables().pressure(), problem.variables().velocity(), consecutiveNumbering);
//
// std::cout.setf(std::ios_base::scientific, std::ios_base::floatfield);
// std::cout.precision(2);
// std::cout << "\t error press \t error grad\t sumflux\t erflm\t\t uMin\t\t uMax\t\t time" << std::endl;
// std::cout << "2pfa\t " << result.relativeL2Error << "\t " << result.ergrad << "\t " << result.sumflux
// << "\t " << result.erflm << "\t " << result.uMin << "\t " << result.uMax << "\t " << time << std::endl;
return 0;
}
catch (Dune::Exception &e) {
std::cerr << "Dune reported error: " << e << std::endl;
}
catch (...) {
std::cerr << "Unknown exception thrown!\n";
throw;
}
return 3;
}
// $Id: test_diffusion_problem.hh 3655 2010-05-26 17:13:50Z bernd $
/*****************************************************************************
* 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, as long as this copyright notice *
* is included in its original form. *
* *
* This program is distributed WITHOUT ANY WARRANTY. *
*****************************************************************************/
#ifndef DUMUX_TEST_2P_PROBLEM_HH
#define DUMUX_TEST_2P_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/unit.hh>
#include <dumux/decoupled/1p/diffusion/diffusionproblem1p.hh>
#include <dumux/decoupled/1p/diffusion/fv/fvvelocity1p.hh>
#include "test_diffusion_spatialparams.hh"
namespace Dumux
{
template<class TypeTag>
class TestProblemOneP;
//////////
// Specify the properties
//////////
namespace Properties
{
NEW_TYPE_TAG(TestProblemOneP, INHERITS_FROM(DecoupledOneP))
;
// Set the grid type
SET_PROP(TestProblemOneP, Grid)
{// typedef Dune::YaspGrid<2> type;
typedef Dune::SGrid<2, 2> type;
};
// Set the wetting phase
SET_PROP(TestProblemOneP, Fluid)
{
private:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
public:
typedef Dumux::LiquidPhase<Scalar, Dumux::Unit<Scalar> > type;
};
// Set the spatial parameters
SET_PROP(TestProblemOneP, SpatialParameters)
{
private:
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Grid)) Grid;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar;
public:
typedef Dumux::TestDiffusionSpatialParams<TypeTag> type;
};
// Enable gravity
SET_BOOL_PROP(TestProblemOneP, EnableGravity, false);
// Set the model
SET_TYPE_PROP(TestProblemOneP, Model, Dumux::FVVelocity1P<TypeTag>);
//Set the problem
SET_TYPE_PROP(TestProblemOneP, Problem, Dumux::TestProblemOneP<TTAG(TestProblemOneP)>);
}
/*!
* \ingroup DecoupledProblems
*/
template<class TypeTag = TTAG(TestProblemOneP)>
class TestProblemOneP: public DiffusionProblem1P<TypeTag, TestProblemOneP<TypeTag> >
{
typedef TestProblemOneP<TypeTag> ThisType;
typedef DiffusionProblem1P<TypeTag, ThisType> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(GridView)) GridView;
typedef typename GET_PROP_TYPE(TypeTag, PTAG(Fluid)) Fluid;
enum
{
dim = GridView::dimension, dimWorld = GridView::dimensionworld
};
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;
public:
TestProblemOneP(const GridView &gridView, const double delta = 1.0) :
ParentType(gridView), delta_(delta)
{
this->spatialParameters().setDelta(delta_);
}
/*!
* \name Problem parameters
*/
// \{
/*!
* \brief The problem name.
*
* This is used as a prefix for files generated by the simulation.
*/
const char *name() const
{
return "test_1p";
}
bool shouldWriteRestartFile() const
{ return false; }
/*!
* \brief Returns the temperature within the domain.
*
* This problem assumes a temperature of 10 degrees Celsius.
*/
Scalar temperature(const GlobalPosition& globalPos, const Element& element) const
{
return 273.15 + 10; // -> 10°C
}
// \}
Scalar referencePressure(const GlobalPosition& globalPos, const Element& element) const
{
return 1e5; // -> 10°C
}
Scalar source(const GlobalPosition& globalPos, const Element& element)
{
double pi = 4.0*atan(1.0);
double rt = globalPos[0]*globalPos[0]+globalPos[1]*globalPos[1];
double ux = pi*cos(pi*globalPos[0])*sin(pi*globalPos[1]);
double uy = pi*cos(pi*globalPos[1])*sin(pi*globalPos[0]);
double kxx = (delta_*globalPos[0]*globalPos[0] + globalPos[1]*globalPos[1])/rt;
double kxy = -(1.0 - delta_)*globalPos[0]*globalPos[1]/rt;
double kyy = (globalPos[0]*globalPos[0] + delta_*globalPos[1]*globalPos[1])/rt;
double f0 = sin(pi*globalPos[0])*sin(pi*globalPos[1])*pi*pi*(1.0 + delta_)*(globalPos[0]*globalPos[0] + globalPos[1]*globalPos[1])
+ cos(pi*globalPos[0])*sin(pi*globalPos[1])*pi*(1.0 - 3.0*delta_)*globalPos[0]
+ cos(pi*globalPos[1])*sin(pi*globalPos[0])*pi*(1.0 - 3.0*delta_)*globalPos[1]
+ cos(pi*globalPos[1])*cos(pi*globalPos[0])*2.0*pi*pi*(1.0 - delta_)*globalPos[0]*globalPos[1];
Scalar result=(f0 + 2.0*(globalPos[0]*(kxx*ux + kxy*uy) + globalPos[1]*(kxy*ux + kyy*uy)))/rt;
return (result);
}
typename BoundaryConditions::Flags bctype(const GlobalPosition& globalPos, const Intersection& intersection) const
{
return BoundaryConditions::dirichlet;
}
Scalar dirichlet(const GlobalPosition& globalPos, const Intersection& intersection) const
{
return (exact(globalPos));
}
Scalar neumann(const GlobalPosition& globalPos, const Intersection& intersection) const
{
return 0.0;
}
Scalar exact (const GlobalPosition& globalPos) const
{
double pi = 4.0*atan(1.0);
return (sin(pi*globalPos[0])*sin(pi*globalPos[1]));
}
Dune::FieldVector<Scalar,dim> exactGrad (const GlobalPosition& globalPos) const
{
Dune::FieldVector<Scalar,dim> grad(0);
double pi = 4.0*atan(1.0);
grad[0] = pi*cos(pi*globalPos[0])*sin(pi*globalPos[1]);
grad[1] = pi*cos(pi*globalPos[1])*sin(pi*globalPos[0]);
return grad;
}
private:
double delta_;
};
} //end namespace
#endif
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment