24#ifndef DUMUX_TEST_1P_PROBLEM_HH
25#define DUMUX_TEST_1P_PROBLEM_HH
27#include <dune/grid/yaspgrid.hh>
41template<
class TypeTag>
55template<
class TypeTag>
74template<
class TypeTag>
78 using TimeManager =
typename GET_PROP_TYPE(TypeTag, TimeManager);
81 using Grid =
typename GridView::Grid;
83 using FluidSystem =
typename GET_PROP_TYPE(TypeTag, FluidSystem);
85 using PrimaryVariables =
typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
86 using BoundaryTypes =
typename GET_PROP_TYPE(TypeTag, BoundaryTypes);
90 dim = GridView::dimension, dimWorld = GridView::dimensionworld
95 using Element =
typename GridView::Traits::template Codim<0>::Entity;
96 using Intersection =
typename GridView::Intersection;
97 using GlobalPosition =
typename Element::Geometry::GlobalCoordinate;
98 using LocalPosition = Dune::FieldVector<Scalar, dim>;
105 delta_ = getParam<Scalar>(
"Problem.Delta", 1e-3);
153 void source(PrimaryVariables &values,
const Element& element)
const
156 values = integratedSource_(element, 4);
165 const Intersection& intersection)
const
167 bcType.setAllDirichlet();
172 const GlobalPosition &globalPos)
const
174 values = exact(globalPos);
179 void neumann(PrimaryVariables &values,
const Intersection& intersection)
const
185 Scalar exact (
const GlobalPosition& globalPos)
const
187 double pi = 4.0*atan(1.0);
189 return (sin(pi*globalPos[0])*sin(pi*globalPos[1]));
192 Dune::FieldVector<Scalar,dim> exactGrad (
const GlobalPosition& globalPos)
const
194 Dune::FieldVector<Scalar,dim> grad(0);
198 double pi = 4.0*atan(1.0);
199 grad[0] = pi*cos(pi*globalPos[0])*sin(pi*globalPos[1]);
200 grad[1] = pi*cos(pi*globalPos[1])*sin(pi*globalPos[0]);
205 Scalar integratedSource_(
const Element& element,
int integrationPoints)
const
208 LocalPosition localPos(0.0);
209 GlobalPosition globalPos(0.0);
210 Scalar halfInterval = 1.0/double(integrationPoints)/2.;
211 for (
int i = 1; i <= integrationPoints; i++)
213 for (
int j = 1; j <= integrationPoints; j++)
215 localPos[0] = double(i)/double(integrationPoints) - halfInterval;
216 localPos[1] = double(j)/double(integrationPoints) - halfInterval;
217 globalPos = element.geometry().global(localPos);
218 source += 1./(integrationPoints*integrationPoints) * evaluateSource_(globalPos);
225 Scalar evaluateSource_(
const GlobalPosition& globalPos)
const
234 Scalar pi = 4.0 * atan(1.0);
235 Scalar x = globalPos[0];
236 Scalar y = globalPos[1];
238 Scalar dpdx = pi * cos(pi * x) * sin(pi * y);
239 Scalar dpdy = pi * sin(pi * x) * cos(pi * y);
240 Scalar dppdxx = -pi * pi * sin(pi * x) * sin(pi * y);
241 Scalar dppdxy = pi * pi * cos(pi * x) * cos(pi * y);
242 Scalar dppdyx = dppdxy;
243 Scalar dppdyy = dppdxx;
244 Scalar kxx = (delta_* x*x + y*y)/(x*x + y*y);
245 Scalar kxy = -(1.0 - delta_) * x * y / (x*x + y*y);
246 Scalar kyy = (x*x + delta_*y*y)/(x*x + y*y);
247 Scalar dkxxdx = 2 * x * y*y * (delta_ - 1.0)/((x*x + y*y) * (x*x + y*y));
248 Scalar dkyydy = 2 * x*x * y * (delta_ - 1.0)/((x*x + y*y) * (x*x + y*y));
249 Scalar dkxydx = (1.0 - delta_) * y * (x*x - y*y) /((x*x + y*y) * (x*x + y*y));
250 Scalar dkxydy = (1.0 - delta_) * x * (y*y - x*x) /((x*x + y*y) * (x*x + y*y));
252 Scalar fx = dkxxdx * dpdx + kxx * dppdxx + dkxydx * dpdy + kxy * dppdyx;
253 Scalar fy = dkxydy * dpdx + kxy * dppdxy + dkyydy * dpdy + kyy * dppdyy;
#define GET_PROP_TYPE(TypeTag, PropTagName)
Definition: propertysystemmacros.hh:283
#define NEW_TYPE_TAG(...)
Definition: propertysystemmacros.hh:130
Setting constant fluid properties via the input file.
A liquid phase consisting of a single component.
spatial parameters for the test problem for diffusion models.
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
SET_TYPE_PROP(FVPressureOneP, Velocity, FVVelocity1P< TypeTag >)
Set velocity reconstruction implementation standard cell centered finite volume schemes as default.
Property tag Velocity
The type velocity reconstruction.
Definition: porousmediumflow/sequential/properties.hh:67
Type tag FVPressureOneP INHERITS_FROM(PressureOneP))
The type tag for the one-phase problems using a standard finite volume model.
std::string viscosity(int phaseIdx) noexcept
I/O name of viscosity for multiphase systems.
Definition: name.hh:74
std::string density(int phaseIdx) noexcept
I/O name of density for multiphase systems.
Definition: name.hh:65
The DUNE grid type.
Definition: common/properties.hh:57
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:69
The type of the spatial parameters object.
Definition: common/properties.hh:221
The type of the fluid system to use.
Definition: common/properties.hh:223
A liquid phase consisting of a single component.
Definition: 1pliquid.hh:46
Base class for all single phase diffusion problems.
Definition: dumux/porousmediumflow/1p/sequential/diffusion/problem.hh:45
SpatialParams & spatialParams()
Returns the spatial parameters object.
Definition: dumux/porousmediumflow/1p/sequential/diffusion/problem.hh:214
void addOutputVtkFields(MultiWriter &writer)
Adds velocity output to the output file.
Definition: sequential/cellcentered/velocity.hh:71
void calculateVelocity()
Function which reconstructs a global velocity field.
Definition: sequential/cellcentered/velocity.hh:96
Base class for definition of an sequential diffusion (pressure) or transport problem.
Definition: onemodelproblem.hh:46
TimeManager & timeManager()
Returns TimeManager object used by the simulation.
Definition: onemodelproblem.hh:551
VtkMultiWriter & resultWriter()
Definition: onemodelproblem.hh:639
test problem for the sequential one-phase model.
Definition: test_1pproblem.hh:76
Scalar referencePressureAtPos(const GlobalPosition &globalPos) const
Returns the reference pressure for evaluation of constitutive relations.
Definition: test_1pproblem.hh:147
void boundaryTypes(BoundaryTypes &bcType, const Intersection &intersection) const
Returns the type of boundary condition.
Definition: test_1pproblem.hh:164
bool shouldWriteRestartFile() const
Definition: test_1pproblem.hh:125
void neumann(PrimaryVariables &values, const Intersection &intersection) const
return neumann condition (flux, [kg/(m^2 s)])
Definition: test_1pproblem.hh:179
std::string name() const
The problem name.
Definition: test_1pproblem.hh:120
void dirichletAtPos(PrimaryVariables &values, const GlobalPosition &globalPos) const
return dirichlet condition (pressure, [Pa])
Definition: test_1pproblem.hh:171
Scalar temperatureAtPos(const GlobalPosition &globalPos) const
Returns the temperature within the domain.
Definition: test_1pproblem.hh:139
void addOutputVtkFields()
Definition: test_1pproblem.hh:128
void source(PrimaryVariables &values, const Element &element) const
source term [kg/(m^3 s)]
Definition: test_1pproblem.hh:153
TestProblemOneP(TimeManager &timeManager, Grid &grid)
Definition: test_1pproblem.hh:102
typename GET_PROP_TYPE(TypeTag, Scalar) Scalar
Definition: test_1pproblem.hh:58
spatial parameters for the test problem for 1-p diffusion models.
Definition: test_1pspatialparams.hh:38
Base class for all single phase diffusion problems.
Defines the properties required for finite volume pressure models.
Finite volume velocity reconstruction.