My Project
lensproblem.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
28 #ifndef EWOMS_LENS_PROBLEM_HH
29 #define EWOMS_LENS_PROBLEM_HH
30 
31 #include <opm/models/io/structuredgridvanguard.hh>
32 #include <opm/models/immiscible/immiscibleproperties.hh>
33 #include <opm/models/discretization/common/fvbaseadlocallinearizer.hh>
34 #include <opm/models/discretization/ecfv/ecfvdiscretization.hh>
35 #include <opm/models/common/transfluxmodule.hh>
36 #include <opm/material/fluidmatrixinteractions/RegularizedVanGenuchten.hpp>
37 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
38 #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
39 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
40 #include <opm/material/fluidsystems/TwoPhaseImmiscibleFluidSystem.hpp>
41 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
42 #include <opm/material/components/SimpleH2O.hpp>
43 #include <opm/material/components/Dnapl.hpp>
44 
45 #include <dune/common/version.hh>
46 #include <dune/common/fvector.hh>
47 #include <dune/common/fmatrix.hh>
48 
49 #include <sstream>
50 #include <string>
51 #include <iostream>
52 
53 namespace Opm {
54 template <class TypeTag>
55 class LensProblem;
56 }
57 
58 namespace Opm::Properties {
59 
60 // Create new type tags
61 namespace TTag {
62 struct LensBaseProblem { using InheritsFrom = std::tuple<StructuredGridVanguard>; };
63 } // end namespace TTag
64 
65 // declare the properties specific for the lens problem
66 template<class TypeTag, class MyTypeTag>
67 struct LensLowerLeftX { using type = UndefinedProperty; };
68 template<class TypeTag, class MyTypeTag>
69 struct LensLowerLeftY { using type = UndefinedProperty; };
70 template<class TypeTag, class MyTypeTag>
71 struct LensLowerLeftZ { using type = UndefinedProperty; };
72 template<class TypeTag, class MyTypeTag>
73 struct LensUpperRightX { using type = UndefinedProperty; };
74 template<class TypeTag, class MyTypeTag>
75 struct LensUpperRightY { using type = UndefinedProperty; };
76 template<class TypeTag, class MyTypeTag>
77 struct LensUpperRightZ { using type = UndefinedProperty; };
78 
79 // Set the problem property
80 template<class TypeTag>
81 struct Problem<TypeTag, TTag::LensBaseProblem> { using type = Opm::LensProblem<TypeTag>; };
82 
83 // Use Dune-grid's YaspGrid
84 template<class TypeTag>
85 struct Grid<TypeTag, TTag::LensBaseProblem> { using type = Dune::YaspGrid<2>; };
86 
87 // Set the wetting phase
88 template<class TypeTag>
89 struct WettingPhase<TypeTag, TTag::LensBaseProblem>
90 {
91 private:
92  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
93 
94 public:
95  using type = Opm::LiquidPhase<Scalar, Opm::SimpleH2O<Scalar> >;
96 };
97 
98 // Set the non-wetting phase
99 template<class TypeTag>
100 struct NonwettingPhase<TypeTag, TTag::LensBaseProblem>
101 {
102 private:
103  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
104 
105 public:
106  using type = Opm::LiquidPhase<Scalar, Opm::DNAPL<Scalar> >;
107 };
108 
109 // Set the material Law
110 template<class TypeTag>
111 struct MaterialLaw<TypeTag, TTag::LensBaseProblem>
112 {
113 private:
114  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
115  enum { wettingPhaseIdx = FluidSystem::wettingPhaseIdx };
116  enum { nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx };
117 
118  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
119  using Traits = Opm::TwoPhaseMaterialTraits<Scalar,
120  /*wettingPhaseIdx=*/FluidSystem::wettingPhaseIdx,
121  /*nonWettingPhaseIdx=*/FluidSystem::nonWettingPhaseIdx>;
122 
123  // define the material law which is parameterized by effective
124  // saturations
125  using EffectiveLaw = Opm::RegularizedVanGenuchten<Traits>;
126 
127 public:
128  // define the material law parameterized by absolute saturations
129  using type = Opm::EffToAbsLaw<EffectiveLaw>;
130 };
131 
132 // Write the solutions of individual newton iterations?
133 template<class TypeTag>
134 struct NewtonWriteConvergence<TypeTag, TTag::LensBaseProblem> { static constexpr bool value = false; };
135 
136 // Use forward differences instead of central differences
137 template<class TypeTag>
138 struct NumericDifferenceMethod<TypeTag, TTag::LensBaseProblem> { static constexpr int value = +1; };
139 
140 // Enable gravity
141 template<class TypeTag>
142 struct EnableGravity<TypeTag, TTag::LensBaseProblem> { static constexpr bool value = true; };
143 
144 // define the properties specific for the lens problem
145 template<class TypeTag>
146 struct LensLowerLeftX<TypeTag, TTag::LensBaseProblem>
147 {
148  using type = GetPropType<TypeTag, Scalar>;
149  static constexpr type value = 1.0;
150 };
151 template<class TypeTag>
152 struct LensLowerLeftY<TypeTag, TTag::LensBaseProblem>
153 {
154  using type = GetPropType<TypeTag, Scalar>;
155  static constexpr type value = 2.0;
156 };
157 template<class TypeTag>
158 struct LensLowerLeftZ<TypeTag, TTag::LensBaseProblem>
159 {
160  using type = GetPropType<TypeTag, Scalar>;
161  static constexpr type value = 0.0;
162 };
163 template<class TypeTag>
164 struct LensUpperRightX<TypeTag, TTag::LensBaseProblem>
165 {
166  using type = GetPropType<TypeTag, Scalar>;
167  static constexpr type value = 4.0;
168 };
169 template<class TypeTag>
170 struct LensUpperRightY<TypeTag, TTag::LensBaseProblem>
171 {
172  using type = GetPropType<TypeTag, Scalar>;
173  static constexpr type value = 3.0;
174 };
175 template<class TypeTag>
176 struct LensUpperRightZ<TypeTag, TTag::LensBaseProblem>
177 {
178  using type = GetPropType<TypeTag, Scalar>;
179  static constexpr type value = 1.0;
180 };
181 
182 template<class TypeTag>
183 struct DomainSizeX<TypeTag, TTag::LensBaseProblem>
184 {
185  using type = GetPropType<TypeTag, Scalar>;
186  static constexpr type value = 6.0;
187 };
188 template<class TypeTag>
189 struct DomainSizeY<TypeTag, TTag::LensBaseProblem>
190 {
191  using type = GetPropType<TypeTag, Scalar>;
192  static constexpr type value = 4.0;
193 };
194 template<class TypeTag>
195 struct DomainSizeZ<TypeTag, TTag::LensBaseProblem>
196 {
197  using type = GetPropType<TypeTag, Scalar>;
198  static constexpr type value = 1.0;
199 };
200 
201 template<class TypeTag>
202 struct CellsX<TypeTag, TTag::LensBaseProblem> { static constexpr int value = 48; };
203 template<class TypeTag>
204 struct CellsY<TypeTag, TTag::LensBaseProblem> { static constexpr int value = 32; };
205 template<class TypeTag>
206 struct CellsZ<TypeTag, TTag::LensBaseProblem> { static constexpr int value = 16; };
207 
208 // The default for the end time of the simulation
209 template<class TypeTag>
210 struct EndTime<TypeTag, TTag::LensBaseProblem>
211 {
212  using type = GetPropType<TypeTag, Scalar>;
213  static constexpr type value = 30e3;
214 };
215 
216 // The default for the initial time step size of the simulation
217 template<class TypeTag>
218 struct InitialTimeStepSize<TypeTag, TTag::LensBaseProblem>
219 {
220  using type = GetPropType<TypeTag, Scalar>;
221  static constexpr type value = 250;
222 };
223 
224 // By default, include the intrinsic permeability tensor to the VTK output files
225 template<class TypeTag>
226 struct VtkWriteIntrinsicPermeabilities<TypeTag, TTag::LensBaseProblem> { static constexpr bool value = true; };
227 
228 // enable the storage cache by default for this problem
229 template<class TypeTag>
230 struct EnableStorageCache<TypeTag, TTag::LensBaseProblem> { static constexpr bool value = true; };
231 
232 // enable the cache for intensive quantities by default for this problem
233 template<class TypeTag>
234 struct EnableIntensiveQuantityCache<TypeTag, TTag::LensBaseProblem> { static constexpr bool value = true; };
235 
236 } // namespace Opm::Properties
237 
238 namespace Opm {
239 
263 template <class TypeTag>
264 class LensProblem : public GetPropType<TypeTag, Properties::BaseProblem>
265 {
266  using ParentType = GetPropType<TypeTag, Properties::BaseProblem>;
267 
268  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
269  using GridView = GetPropType<TypeTag, Properties::GridView>;
270  using Indices = GetPropType<TypeTag, Properties::Indices>;
271  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
272  using WettingPhase = GetPropType<TypeTag, Properties::WettingPhase>;
273  using NonwettingPhase = GetPropType<TypeTag, Properties::NonwettingPhase>;
274  using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
275  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
276  using Model = GetPropType<TypeTag, Properties::Model>;
277 
278  enum {
279  // number of phases
280  numPhases = FluidSystem::numPhases,
281 
282  // phase indices
283  wettingPhaseIdx = FluidSystem::wettingPhaseIdx,
284  nonWettingPhaseIdx = FluidSystem::nonWettingPhaseIdx,
285 
286  // equation indices
287  contiNEqIdx = Indices::conti0EqIdx + nonWettingPhaseIdx,
288 
289  // Grid and world dimension
290  dim = GridView::dimension,
291  dimWorld = GridView::dimensionworld
292  };
293 
294  using EqVector = GetPropType<TypeTag, Properties::EqVector>;
295  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
296  using BoundaryRateVector = GetPropType<TypeTag, Properties::BoundaryRateVector>;
297  using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
298  using MaterialLawParams = GetPropType<TypeTag, Properties::MaterialLawParams>;
299 
300  using CoordScalar = typename GridView::ctype;
301  using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
302 
303  using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
304 
305 public:
309  LensProblem(Simulator& simulator)
310  : ParentType(simulator)
311  { }
312 
316  void finishInit()
317  {
318  ParentType::finishInit();
319 
320  eps_ = 3e-6;
321  FluidSystem::init();
322 
323  temperature_ = 273.15 + 20; // -> 20°C
324  lensLowerLeft_[0] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftX);
325  lensLowerLeft_[1] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftY);
326  lensUpperRight_[0] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightX);
327  lensUpperRight_[1] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightY);
328 
329  if (dimWorld == 3) {
330  lensLowerLeft_[2] = EWOMS_GET_PARAM(TypeTag, Scalar, LensLowerLeftZ);
331  lensUpperRight_[2] = EWOMS_GET_PARAM(TypeTag, Scalar, LensUpperRightZ);
332  }
333 
334  // residual saturations
335  lensMaterialParams_.setResidualSaturation(wettingPhaseIdx, 0.18);
336  lensMaterialParams_.setResidualSaturation(nonWettingPhaseIdx, 0.0);
337  outerMaterialParams_.setResidualSaturation(wettingPhaseIdx, 0.05);
338  outerMaterialParams_.setResidualSaturation(nonWettingPhaseIdx, 0.0);
339 
340  // parameters for the Van Genuchten law: alpha and n
341  lensMaterialParams_.setVgAlpha(0.00045);
342  lensMaterialParams_.setVgN(7.3);
343  outerMaterialParams_.setVgAlpha(0.0037);
344  outerMaterialParams_.setVgN(4.7);
345 
346  lensMaterialParams_.finalize();
347  outerMaterialParams_.finalize();
348 
349  lensK_ = this->toDimMatrix_(9.05e-12);
350  outerK_ = this->toDimMatrix_(4.6e-10);
351 
352  if (dimWorld == 3) {
353  this->gravity_ = 0;
354  this->gravity_[1] = -9.81;
355  }
356  }
357 
361  static void registerParameters()
362  {
363  ParentType::registerParameters();
364 
365  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftX,
366  "The x-coordinate of the lens' lower-left corner "
367  "[m].");
368  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftY,
369  "The y-coordinate of the lens' lower-left corner "
370  "[m].");
371  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightX,
372  "The x-coordinate of the lens' upper-right corner "
373  "[m].");
374  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightY,
375  "The y-coordinate of the lens' upper-right corner "
376  "[m].");
377 
378  if (dimWorld == 3) {
379  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensLowerLeftZ,
380  "The z-coordinate of the lens' lower-left "
381  "corner [m].");
382  EWOMS_REGISTER_PARAM(TypeTag, Scalar, LensUpperRightZ,
383  "The z-coordinate of the lens' upper-right "
384  "corner [m].");
385  }
386  }
387 
391  static std::string briefDescription()
392  {
393  std::string thermal = "isothermal";
394  bool enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>();
395  if (enableEnergy)
396  thermal = "non-isothermal";
397 
398  std::string deriv = "finite difference";
399  using LLS = GetPropType<TypeTag, Properties::LocalLinearizerSplice>;
400  bool useAutoDiff = std::is_same<LLS, Properties::TTag::AutoDiffLocalLinearizer>::value;
401  if (useAutoDiff)
402  deriv = "automatic differentiation";
403 
404  std::string disc = "vertex centered finite volume";
405  using D = GetPropType<TypeTag, Properties::Discretization>;
406  bool useEcfv = std::is_same<D, Opm::EcfvDiscretization<TypeTag>>::value;
407  if (useEcfv)
408  disc = "element centered finite volume";
409 
410  return std::string("")+
411  "Ground remediation problem where a dense oil infiltrates "+
412  "an aquifer with an embedded low-permability lens. " +
413  "This is the binary for the "+thermal+" variant using "+deriv+
414  "and the "+disc+" discretization";
415  }
416 
421 
425  template <class Context>
426  const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
427  unsigned timeIdx) const
428  {
429  const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
430 
431  if (isInLens_(globalPos))
432  return lensK_;
433  return outerK_;
434  }
435 
439  template <class Context>
440  Scalar porosity(const Context& /*context*/,
441  unsigned /*spaceIdx*/,
442  unsigned /*timeIdx*/) const
443  { return 0.4; }
444 
448  template <class Context>
449  const MaterialLawParams& materialLawParams(const Context& context,
450  unsigned spaceIdx, unsigned timeIdx) const
451  {
452  const GlobalPosition& globalPos = context.pos(spaceIdx, timeIdx);
453 
454  if (isInLens_(globalPos))
455  return lensMaterialParams_;
456  return outerMaterialParams_;
457  }
458 
462  template <class Context>
463  Scalar temperature(const Context& /*context*/,
464  unsigned /*spaceIdx*/,
465  unsigned /*timeIdx*/) const
466  { return temperature_; }
467 
469 
474 
478  std::string name() const
479  {
480  using LLS = GetPropType<TypeTag, Properties::LocalLinearizerSplice>;
481 
482  bool useAutoDiff = std::is_same<LLS, Properties::TTag::AutoDiffLocalLinearizer>::value;
483 
484  using FM = GetPropType<TypeTag, Properties::FluxModule>;
485  bool useTrans = std::is_same<FM, Opm::TransFluxModule<TypeTag>>::value;
486 
487  std::ostringstream oss;
488  oss << "lens_" << Model::name()
489  << "_" << Model::discretizationName()
490  << "_" << (useAutoDiff?"ad":"fd");
491  if (useTrans)
492  oss << "_trans";
493 
494  return oss.str();
495  }
496 
501  { }
502 
507  { }
508 
512  void endTimeStep()
513  {
514 #ifndef NDEBUG
515  //this->model().checkConservativeness();
516 
517  // Calculate storage terms
518  EqVector storage;
519  this->model().globalStorage(storage);
520 
521  // Write mass balance information for rank 0
522  if (this->gridView().comm().rank() == 0) {
523  std::cout << "Storage: " << storage << std::endl << std::flush;
524  }
525 #endif // NDEBUG
526  }
527 
529 
534 
538  template <class Context>
539  void boundary(BoundaryRateVector& values,
540  const Context& context,
541  unsigned spaceIdx,
542  unsigned timeIdx) const
543  {
544  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
545 
546  if (onLeftBoundary_(pos) || onRightBoundary_(pos)) {
547  // free flow boundary. we assume incompressible fluids
548  Scalar densityW = WettingPhase::density(temperature_, /*pressure=*/Scalar(1e5));
549  Scalar densityN = NonwettingPhase::density(temperature_, /*pressure=*/Scalar(1e5));
550 
551  Scalar T = temperature(context, spaceIdx, timeIdx);
552  Scalar pw, Sw;
553 
554  // set wetting phase pressure and saturation
555  if (onLeftBoundary_(pos)) {
556  Scalar height = this->boundingBoxMax()[1] - this->boundingBoxMin()[1];
557  Scalar depth = this->boundingBoxMax()[1] - pos[1];
558  Scalar alpha = (1 + 1.5 / height);
559 
560  // hydrostatic pressure scaled by alpha
561  pw = 1e5 - alpha * densityW * this->gravity()[1] * depth;
562  Sw = 1.0;
563  }
564  else {
565  Scalar depth = this->boundingBoxMax()[1] - pos[1];
566 
567  // hydrostatic pressure
568  pw = 1e5 - densityW * this->gravity()[1] * depth;
569  Sw = 1.0;
570  }
571 
572  // specify a full fluid state using pw and Sw
573  const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
574 
575  Opm::ImmiscibleFluidState<Scalar, FluidSystem,
576  /*storeEnthalpy=*/false> fs;
577  fs.setSaturation(wettingPhaseIdx, Sw);
578  fs.setSaturation(nonWettingPhaseIdx, 1 - Sw);
579  fs.setTemperature(T);
580 
581  Scalar pC[numPhases];
582  MaterialLaw::capillaryPressures(pC, matParams, fs);
583  fs.setPressure(wettingPhaseIdx, pw);
584  fs.setPressure(nonWettingPhaseIdx, pw + pC[nonWettingPhaseIdx] - pC[wettingPhaseIdx]);
585 
586  fs.setDensity(wettingPhaseIdx, densityW);
587  fs.setDensity(nonWettingPhaseIdx, densityN);
588 
589  fs.setViscosity(wettingPhaseIdx, WettingPhase::viscosity(temperature_, fs.pressure(wettingPhaseIdx)));
590  fs.setViscosity(nonWettingPhaseIdx, NonwettingPhase::viscosity(temperature_, fs.pressure(nonWettingPhaseIdx)));
591 
592  // impose an freeflow boundary condition
593  values.setFreeFlow(context, spaceIdx, timeIdx, fs);
594  }
595  else if (onInlet_(pos)) {
596  RateVector massRate(0.0);
597  massRate = 0.0;
598  massRate[contiNEqIdx] = -0.04; // kg / (m^2 * s)
599 
600  // impose a forced flow boundary
601  values.setMassRate(massRate);
602  }
603  else {
604  // no flow boundary
605  values.setNoFlow();
606  }
607  }
608 
610 
615 
619  template <class Context>
620  void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx, unsigned timeIdx) const
621  {
622  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
623  Scalar depth = this->boundingBoxMax()[1] - pos[1];
624 
625  Opm::ImmiscibleFluidState<Scalar, FluidSystem> fs;
626  fs.setPressure(wettingPhaseIdx, /*pressure=*/1e5);
627 
628  Scalar Sw = 1.0;
629  fs.setSaturation(wettingPhaseIdx, Sw);
630  fs.setSaturation(nonWettingPhaseIdx, 1 - Sw);
631 
632  fs.setTemperature(temperature_);
633 
634  typename FluidSystem::template ParameterCache<Scalar> paramCache;
635  paramCache.updatePhase(fs, wettingPhaseIdx);
636  Scalar densityW = FluidSystem::density(fs, paramCache, wettingPhaseIdx);
637 
638  // hydrostatic pressure (assuming incompressibility)
639  Scalar pw = 1e5 - densityW * this->gravity()[1] * depth;
640 
641  // calculate the capillary pressure
642  const MaterialLawParams& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
643  Scalar pC[numPhases];
644  MaterialLaw::capillaryPressures(pC, matParams, fs);
645 
646  // make a full fluid state
647  fs.setPressure(wettingPhaseIdx, pw);
648  fs.setPressure(nonWettingPhaseIdx, pw + (pC[wettingPhaseIdx] - pC[nonWettingPhaseIdx]));
649 
650  // assign the primary variables
651  values.assignNaive(fs);
652  }
653 
660  template <class Context>
661  void source(RateVector& rate,
662  const Context& /*context*/,
663  unsigned /*spaceIdx*/,
664  unsigned /*timeIdx*/) const
665  { rate = Scalar(0.0); }
666 
668 
669 private:
670  bool isInLens_(const GlobalPosition& pos) const
671  {
672  for (unsigned i = 0; i < dim; ++i) {
673  if (pos[i] < lensLowerLeft_[i] - eps_ || pos[i] > lensUpperRight_[i]
674  + eps_)
675  return false;
676  }
677  return true;
678  }
679 
680  bool onLeftBoundary_(const GlobalPosition& pos) const
681  { return pos[0] < this->boundingBoxMin()[0] + eps_; }
682 
683  bool onRightBoundary_(const GlobalPosition& pos) const
684  { return pos[0] > this->boundingBoxMax()[0] - eps_; }
685 
686  bool onLowerBoundary_(const GlobalPosition& pos) const
687  { return pos[1] < this->boundingBoxMin()[1] + eps_; }
688 
689  bool onUpperBoundary_(const GlobalPosition& pos) const
690  { return pos[1] > this->boundingBoxMax()[1] - eps_; }
691 
692  bool onInlet_(const GlobalPosition& pos) const
693  {
694  Scalar width = this->boundingBoxMax()[0] - this->boundingBoxMin()[0];
695  Scalar lambda = (this->boundingBoxMax()[0] - pos[0]) / width;
696  return onUpperBoundary_(pos) && 0.5 < lambda && lambda < 2.0 / 3.0;
697  }
698 
699  GlobalPosition lensLowerLeft_;
700  GlobalPosition lensUpperRight_;
701 
702  DimMatrix lensK_;
703  DimMatrix outerK_;
704  MaterialLawParams lensMaterialParams_;
705  MaterialLawParams outerMaterialParams_;
706 
707  Scalar temperature_;
708  Scalar eps_;
709 };
710 
711 } // namespace Opm
712 
713 #endif
Soil contamination problem where DNAPL infiltrates a fully water saturated medium.
Definition: lensproblem.hh:265
Scalar temperature(const Context &, unsigned, unsigned) const
Definition: lensproblem.hh:463
static void registerParameters()
Definition: lensproblem.hh:361
void beginTimeStep()
Definition: lensproblem.hh:500
static std::string briefDescription()
Definition: lensproblem.hh:391
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:620
const MaterialLawParams & materialLawParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:449
Scalar porosity(const Context &, unsigned, unsigned) const
Definition: lensproblem.hh:440
LensProblem(Simulator &simulator)
Definition: lensproblem.hh:309
void finishInit()
Definition: lensproblem.hh:316
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:539
const DimMatrix & intrinsicPermeability(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: lensproblem.hh:426
std::string name() const
Definition: lensproblem.hh:478
void endTimeStep()
Definition: lensproblem.hh:512
void beginIteration()
Definition: lensproblem.hh:506
void source(RateVector &rate, const Context &, unsigned, unsigned) const
Definition: lensproblem.hh:661
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: lensproblem.hh:62