152 using ParentType = GetPropType<TypeTag, Properties::BaseProblem>;
154 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
155 using GridView = GetPropType<TypeTag, Properties::GridView>;
156 using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
157 using MaterialLawParams = GetPropType<TypeTag, Properties::MaterialLawParams>;
158 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
159 using EqVector = GetPropType<TypeTag, Properties::EqVector>;
160 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
161 using BoundaryRateVector = GetPropType<TypeTag, Properties::BoundaryRateVector>;
162 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
163 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
164 using Model = GetPropType<TypeTag, Properties::Model>;
167 using Indices = GetPropType<TypeTag, Properties::Indices>;
170 conti0EqIdx = Indices::conti0EqIdx,
173 numPhases = FluidSystem::numPhases,
176 NAPLIdx = FluidSystem::NAPLIdx,
177 H2OIdx = FluidSystem::H2OIdx,
178 airIdx = FluidSystem::airIdx,
181 waterPhaseIdx = FluidSystem::waterPhaseIdx,
182 gasPhaseIdx = FluidSystem::gasPhaseIdx,
183 naplPhaseIdx = FluidSystem::naplPhaseIdx,
186 dim = GridView::dimension,
187 dimWorld = GridView::dimensionworld
190 using CoordScalar =
typename GridView::ctype;
191 using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
192 using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
199 : ParentType(simulator)
208 ParentType::finishInit();
210 temperature_ = 273.15 + 10.0;
211 FluidSystem::init(temperature_ - 1,
219 fineK_ = this->toDimMatrix_(1e-11);
220 coarseK_ = this->toDimMatrix_(1e-11);
226 materialParams_.setSwr(0.12);
227 materialParams_.setSwrx(0.12);
228 materialParams_.setSnr(0.07);
229 materialParams_.setSgr(0.03);
232 materialParams_.setVgAlpha(0.0005);
233 materialParams_.setVgN(4.);
234 materialParams_.setkrRegardsSnr(
false);
236 materialParams_.finalize();
237 materialParams_.checkDefined();
258 std::ostringstream oss;
259 oss <<
"infiltration_" << Model::name();
269 this->model().checkConservativeness();
273 this->model().globalStorage(storage);
276 if (this->gridView().comm().rank() == 0) {
277 std::cout <<
"Storage: " << storage << std::endl << std::flush;
285 template <
class Context>
289 {
return temperature_; }
294 template <
class Context>
298 unsigned timeIdx)
const
300 const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
301 if (isFineMaterial_(pos))
309 template <
class Context>
313 {
return porosity_; }
318 template <
class Context>
319 const MaterialLawParams&
323 {
return materialParams_; }
335 template <
class Context>
337 const Context& context,
339 unsigned timeIdx)
const
341 const auto& pos = context.pos(spaceIdx, timeIdx);
343 if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
344 Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
346 initialFluidState_(fs, context, spaceIdx, timeIdx);
348 values.setFreeFlow(context, spaceIdx, timeIdx, fs);
350 else if (onInlet_(pos)) {
351 RateVector molarRate(0.0);
352 molarRate[conti0EqIdx + NAPLIdx] = -0.001;
354 values.setMolarRate(molarRate);
355 Opm::Valgrind::CheckDefined(values);
371 template <
class Context>
373 const Context& context,
375 unsigned timeIdx)
const
377 Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
379 initialFluidState_(fs, context, spaceIdx, timeIdx);
382 values.assignMassConservative(fs, matParams,
true);
383 Opm::Valgrind::CheckDefined(values);
392 template <
class Context>
397 { rate = Scalar(0.0); }
402 bool onLeftBoundary_(
const GlobalPosition& pos)
const
403 {
return pos[0] < eps_; }
405 bool onRightBoundary_(
const GlobalPosition& pos)
const
406 {
return pos[0] > this->boundingBoxMax()[0] - eps_; }
408 bool onLowerBoundary_(
const GlobalPosition& pos)
const
409 {
return pos[1] < eps_; }
411 bool onUpperBoundary_(
const GlobalPosition& pos)
const
412 {
return pos[1] > this->boundingBoxMax()[1] - eps_; }
414 bool onInlet_(
const GlobalPosition& pos)
const
415 {
return onUpperBoundary_(pos) && 50 < pos[0] && pos[0] < 75; }
417 template <
class Flu
idState,
class Context>
418 void initialFluidState_(FluidState& fs,
const Context& context,
419 unsigned spaceIdx,
unsigned timeIdx)
const
421 const GlobalPosition pos = context.pos(spaceIdx, timeIdx);
425 Scalar densityW = 1000.0;
426 Scalar pc = 9.81 * densityW * (y - (5 - 5e-4 * x));
432 Scalar Sw = matParams.Swr();
433 Scalar Swr = matParams.Swr();
434 Scalar Sgr = matParams.Sgr();
441 Opm::Valgrind::CheckDefined(Sw);
442 Opm::Valgrind::CheckDefined(Sg);
444 fs.setSaturation(waterPhaseIdx, Sw);
445 fs.setSaturation(gasPhaseIdx, Sg);
446 fs.setSaturation(naplPhaseIdx, 0);
449 fs.setTemperature(temperature_);
452 Scalar pcAll[numPhases];
454 if (onLeftBoundary_(pos))
456 MaterialLaw::capillaryPressures(pcAll, matParams, fs);
457 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
458 fs.setPressure(phaseIdx, pg + (pcAll[phaseIdx] - pcAll[gasPhaseIdx]));
461 fs.setMoleFraction(gasPhaseIdx, H2OIdx, 1e-6);
462 fs.setMoleFraction(gasPhaseIdx, airIdx,
463 1 - fs.moleFraction(gasPhaseIdx, H2OIdx));
464 fs.setMoleFraction(gasPhaseIdx, NAPLIdx, 0);
466 using CFRP = Opm::ComputeFromReferencePhase<Scalar, FluidSystem>;
467 typename FluidSystem::template ParameterCache<Scalar> paramCache;
468 CFRP::solve(fs, paramCache, gasPhaseIdx,
472 fs.setMoleFraction(waterPhaseIdx, H2OIdx,
473 1 - fs.moleFraction(waterPhaseIdx, H2OIdx));
476 bool isFineMaterial_(
const GlobalPosition& pos)
const
477 {
return 70. <= pos[0] && pos[0] <= 85. && 7.0 <= pos[1] && pos[1] <= 7.50; }
484 MaterialLawParams materialParams_;