28#ifndef EWOMS_POWER_INJECTION_PROBLEM_HH
29#define EWOMS_POWER_INJECTION_PROBLEM_HH
31#include <opm/models/immiscible/immisciblemodel.hh>
32#include <opm/models/io/cubegridvanguard.hh>
34#include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
35#include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
36#include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
37#include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
38#include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp>
39#include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
40#include <opm/material/components/SimpleH2O.hpp>
41#include <opm/material/components/Air.hpp>
43#include <dune/grid/yaspgrid.hh>
45#include <dune/common/version.hh>
46#include <dune/common/fvector.hh>
47#include <dune/common/fmatrix.hh>
55template <
class TypeTag>
56class PowerInjectionProblem;
59namespace Opm::Properties {
66template<
class TypeTag>
67struct Grid<TypeTag, TTag::PowerInjectionBaseProblem> {
using type = Dune::YaspGrid<1>; };
70template<
class TypeTag>
71struct Vanguard<TypeTag, TTag::PowerInjectionBaseProblem> {
using type = Opm::CubeGridVanguard<TypeTag>; };
74template<
class TypeTag>
78template<
class TypeTag>
79struct WettingPhase<TypeTag, TTag::PowerInjectionBaseProblem>
82 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
85 using type = Opm::LiquidPhase<Scalar, Opm::SimpleH2O<Scalar> >;
89template<
class TypeTag>
90struct NonwettingPhase<TypeTag, TTag::PowerInjectionBaseProblem>
93 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
96 using type = Opm::GasPhase<Scalar, Opm::Air<Scalar> >;
100template<
class TypeTag>
101struct MaterialLaw<TypeTag, TTag::PowerInjectionBaseProblem>
104 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
105 enum { wettingPhaseIdx = FluidSystem::wettingPhaseIdx };
106 enum { nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx };
108 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
109 using Traits = Opm::TwoPhaseMaterialTraits<Scalar,
110 FluidSystem::wettingPhaseIdx,
111 FluidSystem::nonWettingPhaseIdx>;
115 using EffectiveLaw = Opm::RegularizedVanGenuchten<Traits>;
119 using type = Opm::EffToAbsLaw<EffectiveLaw>;
123template<
class TypeTag>
124struct VtkWriteFilterVelocities<TypeTag, TTag::PowerInjectionBaseProblem> {
static constexpr bool value =
true; };
127template<
class TypeTag>
128struct EnableGravity<TypeTag, TTag::PowerInjectionBaseProblem> {
static constexpr bool value =
false; };
131template<
class TypeTag>
132struct DomainSizeX<TypeTag, TTag::PowerInjectionBaseProblem>
134 using type = GetPropType<TypeTag, Scalar>;
135 static constexpr type value = 100.0;
137template<
class TypeTag>
138struct DomainSizeY<TypeTag, TTag::PowerInjectionBaseProblem>
140 using type = GetPropType<TypeTag, Scalar>;
141 static constexpr type value = 1.0;
143template<
class TypeTag>
144struct DomainSizeZ<TypeTag, TTag::PowerInjectionBaseProblem>
146 using type = GetPropType<TypeTag, Scalar>;
147 static constexpr type value = 1.0;
150template<
class TypeTag>
151struct CellsX<TypeTag, TTag::PowerInjectionBaseProblem> {
static constexpr int value = 250; };
152template<
class TypeTag>
153struct CellsY<TypeTag, TTag::PowerInjectionBaseProblem> {
static constexpr int value = 1; };
154template<
class TypeTag>
155struct CellsZ<TypeTag, TTag::PowerInjectionBaseProblem> {
static constexpr int value = 1; };
158template<
class TypeTag>
159struct EndTime<TypeTag, TTag::PowerInjectionBaseProblem>
161 using type = GetPropType<TypeTag, Scalar>;
162 static constexpr type value = 100;
166template<
class TypeTag>
167struct InitialTimeStepSize<TypeTag, TTag::PowerInjectionBaseProblem>
169 using type = GetPropType<TypeTag, Scalar>;
170 static constexpr type value = 1e-3;
188template <
class TypeTag>
191 using ParentType = GetPropType<TypeTag, Properties::BaseProblem>;
193 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
194 using GridView = GetPropType<TypeTag, Properties::GridView>;
195 using Indices = GetPropType<TypeTag, Properties::Indices>;
196 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
197 using WettingPhase = GetPropType<TypeTag, Properties::WettingPhase>;
198 using NonwettingPhase = GetPropType<TypeTag, Properties::NonwettingPhase>;
199 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
200 using EqVector = GetPropType<TypeTag, Properties::EqVector>;
201 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
202 using BoundaryRateVector = GetPropType<TypeTag, Properties::BoundaryRateVector>;
203 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
207 numPhases = FluidSystem::numPhases,
210 wettingPhaseIdx = FluidSystem::wettingPhaseIdx,
211 nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx,
214 contiNEqIdx = Indices::conti0EqIdx + nonWettingPhaseIdx,
217 dim = GridView::dimension,
218 dimWorld = GridView::dimensionworld
221 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
222 using MaterialLawParams = GetPropType<TypeTag, Properties::MaterialLawParams>;
224 using CoordScalar =
typename GridView::ctype;
225 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
227 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
234 : ParentType(simulator)
242 ParentType::finishInit();
247 temperature_ = 273.15 + 26.6;
251 materialParams_.setVgAlpha(0.00045);
252 materialParams_.setVgN(7.3);
253 materialParams_.finalize();
255 K_ = this->toDimMatrix_(5.73e-08);
257 setupInitialFluidState_();
270 std::ostringstream oss;
271 oss <<
"powerinjection_";
272 if (std::is_same<GetPropType<TypeTag, Properties::FluxModule>,
273 Opm::DarcyFluxModule<TypeTag> >::value)
276 oss <<
"forchheimer";
278 if (std::is_same<GetPropType<TypeTag, Properties::LocalLinearizerSplice>,
279 Properties::TTag::AutoDiffLocalLinearizer>::value)
293 this->model().checkConservativeness();
297 this->model().globalStorage(storage);
300 if (this->gridView().comm().rank() == 0) {
301 std::cout <<
"Storage: " << storage << std::endl << std::flush;
315 template <
class Context>
324 template <
class Context>
333 template <
class Context>
342 template <
class Context>
343 const MaterialLawParams&
347 {
return materialParams_; }
352 template <
class Context>
356 {
return temperature_; }
371 template <
class Context>
373 const Context& context,
375 unsigned timeIdx)
const
377 const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
379 if (onLeftBoundary_(pos)) {
380 RateVector massRate(0.0);
382 massRate[contiNEqIdx] = -1.00;
385 values.setMassRate(massRate);
387 else if (onRightBoundary_(pos))
389 values.setFreeFlow(context, spaceIdx, timeIdx, initialFluidState_);
404 template <
class Context>
411 values.assignNaive(initialFluidState_);
420 template <
class Context>
425 { rate = Scalar(0.0); }
430 bool onLeftBoundary_(
const GlobalPosition& pos)
const
431 {
return pos[0] < this->boundingBoxMin()[0] + eps_; }
433 bool onRightBoundary_(
const GlobalPosition& pos)
const
434 {
return pos[0] > this->boundingBoxMax()[0] - eps_; }
436 void setupInitialFluidState_()
438 initialFluidState_.setTemperature(temperature_);
441 initialFluidState_.setSaturation(wettingPhaseIdx, Sw);
442 initialFluidState_.setSaturation(nonWettingPhaseIdx, 1 - Sw);
445 initialFluidState_.setPressure(wettingPhaseIdx, p);
446 initialFluidState_.setPressure(nonWettingPhaseIdx, p);
448 typename FluidSystem::template ParameterCache<Scalar> paramCache;
449 paramCache.updateAll(initialFluidState_);
450 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
451 initialFluidState_.setDensity(phaseIdx,
452 FluidSystem::density(initialFluidState_, paramCache, phaseIdx));
453 initialFluidState_.setViscosity(phaseIdx,
454 FluidSystem::viscosity(initialFluidState_, paramCache, phaseIdx));
459 MaterialLawParams materialParams_;
461 Opm::ImmiscibleFluidState<Scalar, FluidSystem> initialFluidState_;
1D Problem with very fast injection of gas on the left.
Definition: powerinjectionproblem.hh:190
void finishInit()
Definition: powerinjectionproblem.hh:240
Scalar porosity(const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:334
Scalar temperature(const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:353
void endTimeStep()
Definition: powerinjectionproblem.hh:290
void initial(PrimaryVariables &values, const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:405
Scalar ergunCoefficient(const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:325
const MaterialLawParams & materialLawParams(const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:344
std::string name() const
Definition: powerinjectionproblem.hh:268
void source(RateVector &rate, const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:421
PowerInjectionProblem(Simulator &simulator)
Definition: powerinjectionproblem.hh:233
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: powerinjectionproblem.hh:372
const DimMatrix & intrinsicPermeability(const Context &, unsigned, unsigned) const
Definition: powerinjectionproblem.hh:316
Definition: powerinjectionproblem.hh:62