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