28#ifndef EWOMS_GROUND_WATER_PROBLEM_HH
29#define EWOMS_GROUND_WATER_PROBLEM_HH
31#include <opm/models/immiscible/immiscibleproperties.hh>
32#include <opm/simulators/linalg/parallelistlbackend.hh>
34#include <opm/material/components/SimpleH2O.hpp>
35#include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
36#include <opm/material/fluidsystems/LiquidPhase.hpp>
38#include <dune/grid/yaspgrid.hh>
39#include <dune/grid/io/file/dgfparser/dgfyasp.hh>
41#include <dune/common/version.hh>
42#include <dune/common/fmatrix.hh>
43#include <dune/common/fvector.hh>
49template <
class TypeTag>
50class GroundWaterProblem;
53namespace Opm::Properties {
59template<
class TypeTag,
class MyTypeTag>
61template<
class TypeTag,
class MyTypeTag>
63template<
class TypeTag,
class MyTypeTag>
65template<
class TypeTag,
class MyTypeTag>
67template<
class TypeTag,
class MyTypeTag>
69template<
class TypeTag,
class MyTypeTag>
71template<
class TypeTag,
class MyTypeTag>
73template<
class TypeTag,
class MyTypeTag>
76template<
class TypeTag>
77struct Fluid<TypeTag, TTag::GroundWaterBaseProblem>
80 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
83 using type = Opm::LiquidPhase<Scalar, Opm::SimpleH2O<Scalar> >;
87template<
class TypeTag>
88struct Grid<TypeTag, TTag::GroundWaterBaseProblem> {
using type = Dune::YaspGrid<2>; };
91template<
class TypeTag>
92struct Problem<TypeTag, TTag::GroundWaterBaseProblem>
95template<
class TypeTag>
98 using type = GetPropType<TypeTag, Scalar>;
99 static constexpr type value = 0.25;
101template<
class TypeTag>
104 using type = GetPropType<TypeTag, Scalar>;
105 static constexpr type value = 0.25;
107template<
class TypeTag>
110 using type = GetPropType<TypeTag, Scalar>;
111 static constexpr type value = 0.25;
113template<
class TypeTag>
116 using type = GetPropType<TypeTag, Scalar>;
117 static constexpr type value = 0.75;
119template<
class TypeTag>
122 using type = GetPropType<TypeTag, Scalar>;
123 static constexpr type value = 0.75;
125template<
class TypeTag>
128 using type = GetPropType<TypeTag, Scalar>;
129 static constexpr type value = 0.75;
131template<
class TypeTag>
134 using type = GetPropType<TypeTag, Scalar>;
135 static constexpr type value = 1e-10;
137template<
class TypeTag>
140 using type = GetPropType<TypeTag, Scalar>;
141 static constexpr type value = 1e-12;
145template<
class TypeTag>
146struct EnableGravity<TypeTag, TTag::GroundWaterBaseProblem> {
static constexpr bool value =
true; };
149template<
class TypeTag>
150struct EndTime<TypeTag, TTag::GroundWaterBaseProblem>
152 using type = GetPropType<TypeTag, Scalar>;
153 static constexpr type value = 1;
157template<
class TypeTag>
158struct InitialTimeStepSize<TypeTag, TTag::GroundWaterBaseProblem>
160 using type = GetPropType<TypeTag, Scalar>;
161 static constexpr type value = 1;
165template<
class TypeTag>
166struct GridFile<TypeTag, TTag::GroundWaterBaseProblem> {
static constexpr auto value =
"./data/groundwater_2d.dgf"; };
170template<
class TypeTag>
171struct LinearSolverSplice<TypeTag, TTag::GroundWaterBaseProblem> {
using type = TTag::ParallelIstlLinearSolver; };
173template<
class TypeTag>
174struct LinearSolverWrapper<TypeTag, TTag::GroundWaterBaseProblem>
175{
using type = Opm::Linear::SolverWrapperConjugatedGradients<TypeTag>; };
192template <
class TypeTag>
195 using ParentType = GetPropType<TypeTag, Properties::BaseProblem>;
197 using GridView = GetPropType<TypeTag, Properties::GridView>;
198 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
199 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
202 using Indices = GetPropType<TypeTag, Properties::Indices>;
204 numPhases = FluidSystem::numPhases,
207 dim = GridView::dimension,
208 dimWorld = GridView::dimensionworld,
211 pressure0Idx = Indices::pressure0Idx
214 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
215 using EqVector = GetPropType<TypeTag, Properties::EqVector>;
216 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
217 using BoundaryRateVector = GetPropType<TypeTag, Properties::BoundaryRateVector>;
218 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
219 using Model = GetPropType<TypeTag, Properties::Model>;
221 using CoordScalar =
typename GridView::ctype;
222 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
224 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
231 : ParentType(simulator)
239 ParentType::finishInit();
243 lensLowerLeft_[0] = Parameters::get<TypeTag, Properties::LensLowerLeftX>();
245 lensLowerLeft_[1] = Parameters::get<TypeTag, Properties::LensLowerLeftY>();
247 lensLowerLeft_[2] = Parameters::get<TypeTag, Properties::LensLowerLeftY>();
249 lensUpperRight_[0] = Parameters::get<TypeTag, Properties::LensUpperRightX>();
251 lensUpperRight_[1] = Parameters::get<TypeTag, Properties::LensUpperRightY>();
253 lensUpperRight_[2] = Parameters::get<TypeTag, Properties::LensUpperRightY>();
255 intrinsicPerm_ = this->toDimMatrix_(Parameters::get<TypeTag, Properties::Permeability>());
256 intrinsicPermLens_ = this->toDimMatrix_(Parameters::get<TypeTag, Properties::PermeabilityLens>());
264 ParentType::registerParameters();
266 Parameters::registerParam<TypeTag, Properties::LensLowerLeftX>
267 (
"The x-coordinate of the lens' lower-left corner [m].");
268 Parameters::registerParam<TypeTag, Properties::LensUpperRightX>
269 (
"The x-coordinate of the lens' upper-right corner [m].");
272 Parameters::registerParam<TypeTag, Properties::LensLowerLeftY>
273 (
"The y-coordinate of the lens' lower-left corner [m].");
274 Parameters::registerParam<TypeTag, Properties::LensUpperRightY>
275 (
"The y-coordinate of the lens' upper-right corner [m].");
279 Parameters::registerParam<TypeTag, Properties::LensLowerLeftZ>
280 (
"The z-coordinate of the lens' lower-left corner [m].");
281 Parameters::registerParam<TypeTag, Properties::LensUpperRightZ>
282 (
"The z-coordinate of the lens' upper-right corner [m].");
285 Parameters::registerParam<TypeTag, Properties::Permeability>
286 (
"The intrinsic permeability [m^2] of the ambient material.");
287 Parameters::registerParam<TypeTag, Properties::PermeabilityLens>
288 (
"The intrinsic permeability [m^2] of the lens.");
301 std::ostringstream oss;
302 oss <<
"groundwater_" << Model::name();
312 this->model().checkConservativeness();
316 this->model().globalStorage(storage);
319 if (this->gridView().comm().rank() == 0) {
320 std::cout <<
"Storage: " << storage << std::endl << std::flush;
328 template <
class Context>
332 {
return 273.15 + 10; }
337 template <
class Context>
346 template <
class Context>
349 unsigned timeIdx)
const
351 if (isInLens_(context.pos(spaceIdx, timeIdx)))
352 return intrinsicPermLens_;
354 return intrinsicPerm_;
366 template <
class Context>
367 void boundary(BoundaryRateVector& values,
const Context& context,
368 unsigned spaceIdx,
unsigned timeIdx)
const
370 const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
372 if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos)) {
374 Scalar T =
temperature(context, spaceIdx, timeIdx);
375 if (onLowerBoundary_(globalPos))
380 Opm::ImmiscibleFluidState<Scalar, FluidSystem,
382 fs.setSaturation(0, 1.0);
383 fs.setPressure(0, pressure);
384 fs.setTemperature(T);
386 typename FluidSystem::template ParameterCache<Scalar> paramCache;
387 paramCache.updateAll(fs);
388 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
389 fs.setDensity(phaseIdx, FluidSystem::density(fs, paramCache, phaseIdx));
390 fs.setViscosity(phaseIdx, FluidSystem::viscosity(fs, paramCache, phaseIdx));
394 values.setFreeFlow(context, spaceIdx, timeIdx, fs);
412 template <
class Context>
419 values[pressure0Idx] = 1.0e+5;
425 template <
class Context>
430 { rate = Scalar(0.0); }
435 bool onLowerBoundary_(
const GlobalPosition& pos)
const
436 {
return pos[dim - 1] < eps_; }
438 bool onUpperBoundary_(
const GlobalPosition& pos)
const
439 {
return pos[dim - 1] > this->boundingBoxMax()[dim - 1] - eps_; }
441 bool isInLens_(
const GlobalPosition& pos)
const
443 return lensLowerLeft_[0] <= pos[0] && pos[0] <= lensUpperRight_[0]
444 && lensLowerLeft_[1] <= pos[1] && pos[1] <= lensUpperRight_[1];
447 GlobalPosition lensLowerLeft_;
448 GlobalPosition lensUpperRight_;
450 DimMatrix intrinsicPerm_;
451 DimMatrix intrinsicPermLens_;
Test for the immisicible VCVF discretization with only a single phase.
Definition groundwaterproblem.hh:194
Scalar porosity(const Context &, unsigned, unsigned) const
Definition groundwaterproblem.hh:338
void endTimeStep()
Definition groundwaterproblem.hh:309
GroundWaterProblem(Simulator &simulator)
Definition groundwaterproblem.hh:230
static void registerParameters()
Definition groundwaterproblem.hh:262
std::string name() const
Definition groundwaterproblem.hh:299
Scalar temperature(const Context &, unsigned, unsigned) const
Definition groundwaterproblem.hh:329
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition groundwaterproblem.hh:367
void source(RateVector &rate, const Context &, unsigned, unsigned) const
Definition groundwaterproblem.hh:426
const DimMatrix & intrinsicPermeability(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition groundwaterproblem.hh:347
void finishInit()
Definition groundwaterproblem.hh:237
void initial(PrimaryVariables &values, const Context &, unsigned, unsigned) const
Definition groundwaterproblem.hh:413
Definition groundwaterproblem.hh:60
Definition groundwaterproblem.hh:62
Definition groundwaterproblem.hh:64
Definition groundwaterproblem.hh:66
Definition groundwaterproblem.hh:68
Definition groundwaterproblem.hh:70
Definition groundwaterproblem.hh:74
Definition groundwaterproblem.hh:72
Definition groundwaterproblem.hh:56