My Project
co2injectionproblem.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_CO2_INJECTION_PROBLEM_HH
29 #define EWOMS_CO2_INJECTION_PROBLEM_HH
30 
31 #include <opm/models/immiscible/immisciblemodel.hh>
32 #include <opm/simulators/linalg/parallelamgbackend.hh>
33 
34 #include <opm/material/fluidsystems/H2ON2FluidSystem.hpp>
35 #include <opm/material/fluidsystems/BrineCO2FluidSystem.hpp>
36 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
37 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
38 #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
39 #include <opm/material/fluidmatrixinteractions/LinearMaterial.hpp>
40 #include <opm/material/fluidmatrixinteractions/RegularizedBrooksCorey.hpp>
41 #include <opm/material/fluidmatrixinteractions/EffToAbsLaw.hpp>
42 #include <opm/material/fluidmatrixinteractions/MaterialTraits.hpp>
43 #include <opm/material/thermal/SomertonThermalConductionLaw.hpp>
44 #include <opm/material/thermal/ConstantSolidHeatCapLaw.hpp>
45 #include <opm/material/binarycoefficients/Brine_CO2.hpp>
46 #include <opm/material/common/UniformTabulated2DFunction.hpp>
47 
48 #include <dune/grid/yaspgrid.hh>
49 #include <dune/grid/io/file/dgfparser/dgfyasp.hh>
50 
51 #include <dune/common/version.hh>
52 #include <dune/common/fvector.hh>
53 #include <dune/common/fmatrix.hh>
54 
55 #include <sstream>
56 #include <iostream>
57 #include <string>
58 
59 namespace Opm {
61 template <class TypeTag>
62 class Co2InjectionProblem;
63 
64 namespace Co2Injection {
65 #include <opm/material/components/co2tables.inc>
66 }
68 }
69 
70 namespace Opm::Properties {
71 
72 namespace TTag {
74 }
75 
76 // declare the CO2 injection problem specific property tags
77 template<class TypeTag, class MyTypeTag>
78 struct FluidSystemPressureLow { using type = UndefinedProperty; };
79 template<class TypeTag, class MyTypeTag>
80 struct FluidSystemPressureHigh { using type = UndefinedProperty; };
81 template<class TypeTag, class MyTypeTag>
82 struct FluidSystemNumPressure { using type = UndefinedProperty; };
83 template<class TypeTag, class MyTypeTag>
84 struct FluidSystemTemperatureLow { using type = UndefinedProperty; };
85 template<class TypeTag, class MyTypeTag>
86 struct FluidSystemTemperatureHigh { using type = UndefinedProperty; };
87 template<class TypeTag, class MyTypeTag>
88 struct FluidSystemNumTemperature { using type = UndefinedProperty; };
89 
90 template<class TypeTag, class MyTypeTag>
91 struct MaxDepth { using type = UndefinedProperty; };
92 template<class TypeTag, class MyTypeTag>
93 struct Temperature { using type = UndefinedProperty; };
94 template<class TypeTag, class MyTypeTag>
95 struct SimulationName { using type = UndefinedProperty; };
96 
97 // Set the grid type
98 template<class TypeTag>
99 struct Grid<TypeTag, TTag::Co2InjectionBaseProblem> { using type = Dune::YaspGrid<2>; };
100 
101 // Set the problem property
102 template<class TypeTag>
103 struct Problem<TypeTag, TTag::Co2InjectionBaseProblem>
105 
106 // Set fluid configuration
107 template<class TypeTag>
108 struct FluidSystem<TypeTag, TTag::Co2InjectionBaseProblem>
109 {
110 private:
111  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
112  using CO2Tables = Opm::Co2Injection::CO2Tables;
113 
114 public:
115  using type = Opm::BrineCO2FluidSystem<Scalar, CO2Tables>;
116  //using type = Opm::H2ON2FluidSystem<Scalar, /*useComplexRelations=*/false>;
117 };
118 
119 // Set the material Law
120 template<class TypeTag>
121 struct MaterialLaw<TypeTag, TTag::Co2InjectionBaseProblem>
122 {
123 private:
124  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
125  enum { liquidPhaseIdx = FluidSystem::liquidPhaseIdx };
126  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
127 
128  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
129  using Traits = Opm::TwoPhaseMaterialTraits<Scalar,
130  /*wettingPhaseIdx=*/FluidSystem::liquidPhaseIdx,
131  /*nonWettingPhaseIdx=*/FluidSystem::gasPhaseIdx>;
132 
133  // define the material law which is parameterized by effective
134  // saturations
135  using EffMaterialLaw = Opm::RegularizedBrooksCorey<Traits>;
136 
137 public:
138  // define the material law parameterized by absolute saturations
139  using type = Opm::EffToAbsLaw<EffMaterialLaw>;
140 };
141 
142 // Set the thermal conduction law
143 template<class TypeTag>
144 struct ThermalConductionLaw<TypeTag, TTag::Co2InjectionBaseProblem>
145 {
146 private:
147  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
148  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
149 
150 public:
151  // define the material law parameterized by absolute saturations
152  using type = Opm::SomertonThermalConductionLaw<FluidSystem, Scalar>;
153 };
154 
155 // set the energy storage law for the solid phase
156 template<class TypeTag>
157 struct SolidEnergyLaw<TypeTag, TTag::Co2InjectionBaseProblem>
158 { using type = Opm::ConstantSolidHeatCapLaw<GetPropType<TypeTag, Properties::Scalar>>; };
159 
160 // Use the algebraic multi-grid linear solver for this problem
161 template<class TypeTag>
162 struct LinearSolverSplice<TypeTag, TTag::Co2InjectionBaseProblem> { using type = TTag::ParallelAmgLinearSolver; };
163 
164 // Write the Newton convergence behavior to disk?
165 template<class TypeTag>
166 struct NewtonWriteConvergence<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr bool value = false; };
167 
168 // Enable gravity
169 template<class TypeTag>
170 struct EnableGravity<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr bool value = true; };
171 
172 // set the defaults for the problem specific properties
173 template<class TypeTag>
174 struct FluidSystemPressureLow<TypeTag, TTag::Co2InjectionBaseProblem>
175 {
176  using type = GetPropType<TypeTag, Scalar>;
177  static constexpr type value = 3e7;
178 };
179 template<class TypeTag>
180 struct FluidSystemPressureHigh<TypeTag, TTag::Co2InjectionBaseProblem>
181 {
182  using type = GetPropType<TypeTag, Scalar>;
183  static constexpr type value = 4e7;
184 };
185 template<class TypeTag>
186 struct FluidSystemNumPressure<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr int value = 100; };
187 template<class TypeTag>
188 struct FluidSystemTemperatureLow<TypeTag, TTag::Co2InjectionBaseProblem>
189 {
190  using type = GetPropType<TypeTag, Scalar>;
191  static constexpr type value = 290;
192 };
193 template<class TypeTag>
194 struct FluidSystemTemperatureHigh<TypeTag, TTag::Co2InjectionBaseProblem>
195 {
196  using type = GetPropType<TypeTag, Scalar>;
197  static constexpr type value = 500;
198 };
199 template<class TypeTag>
200 struct FluidSystemNumTemperature<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr int value = 100; };
201 
202 template<class TypeTag>
203 struct MaxDepth<TypeTag, TTag::Co2InjectionBaseProblem>
204 {
205  using type = GetPropType<TypeTag, Scalar>;
206  static constexpr type value = 2500;
207 };
208 template<class TypeTag>
209 struct Temperature<TypeTag, TTag::Co2InjectionBaseProblem>
210 {
211  using type = GetPropType<TypeTag, Scalar>;
212  static constexpr type value = 293.15;
213 };
214 template<class TypeTag>
215 struct SimulationName<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr auto value = "co2injection"; };
216 
217 // The default for the end time of the simulation
218 template<class TypeTag>
219 struct EndTime<TypeTag, TTag::Co2InjectionBaseProblem>
220 {
221  using type = GetPropType<TypeTag, Scalar>;
222  static constexpr type value = 1e4;
223 };
224 
225 // The default for the initial time step size of the simulation
226 template<class TypeTag>
227 struct InitialTimeStepSize<TypeTag, TTag::Co2InjectionBaseProblem>
228 {
229  using type = GetPropType<TypeTag, Scalar>;
230  static constexpr type value = 250;
231 };
232 
233 // The default DGF file to load
234 template<class TypeTag>
235 struct GridFile<TypeTag, TTag::Co2InjectionBaseProblem> { static constexpr auto value = "data/co2injection.dgf"; };
236 
237 } // namespace Opm::Properties
238 
239 namespace Opm {
262 template <class TypeTag>
263 class Co2InjectionProblem : public GetPropType<TypeTag, Properties::BaseProblem>
264 {
265  using ParentType = GetPropType<TypeTag, Properties::BaseProblem>;
266 
267  using Scalar = GetPropType<TypeTag, Properties::Scalar>;
268  using Evaluation = GetPropType<TypeTag, Properties::Evaluation>;
269  using GridView = GetPropType<TypeTag, Properties::GridView>;
270  using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
271 
272  enum { dim = GridView::dimension };
273  enum { dimWorld = GridView::dimensionworld };
274 
275  // copy some indices for convenience
276  using Indices = GetPropType<TypeTag, Properties::Indices>;
277  enum { numPhases = FluidSystem::numPhases };
278  enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
279  enum { liquidPhaseIdx = FluidSystem::liquidPhaseIdx };
280  enum { CO2Idx = FluidSystem::CO2Idx };
281  enum { BrineIdx = FluidSystem::BrineIdx };
282  enum { conti0EqIdx = Indices::conti0EqIdx };
283  enum { contiCO2EqIdx = conti0EqIdx + CO2Idx };
284 
285  using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
286  using RateVector = GetPropType<TypeTag, Properties::RateVector>;
287  using BoundaryRateVector = GetPropType<TypeTag, Properties::BoundaryRateVector>;
288  using MaterialLaw = GetPropType<TypeTag, Properties::MaterialLaw>;
289  using Simulator = GetPropType<TypeTag, Properties::Simulator>;
290  using Model = GetPropType<TypeTag, Properties::Model>;
291  using MaterialLawParams = GetPropType<TypeTag, Properties::MaterialLawParams>;
292  using ThermalConductionLaw = GetPropType<TypeTag, Properties::ThermalConductionLaw>;
293  using SolidEnergyLawParams = GetPropType<TypeTag, Properties::SolidEnergyLawParams>;
294  using ThermalConductionLawParams = typename ThermalConductionLaw::Params;
295 
296  using Toolbox = Opm::MathToolbox<Evaluation>;
297  using CoordScalar = typename GridView::ctype;
298  using GlobalPosition = Dune::FieldVector<CoordScalar, dimWorld>;
299  using DimMatrix = Dune::FieldMatrix<Scalar, dimWorld, dimWorld>;
300 
301 public:
305  Co2InjectionProblem(Simulator& simulator)
306  : ParentType(simulator)
307  { }
308 
312  void finishInit()
313  {
314  ParentType::finishInit();
315 
316  eps_ = 1e-6;
317 
318  temperatureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureLow);
319  temperatureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemTemperatureHigh);
320  nTemperature_ = EWOMS_GET_PARAM(TypeTag, unsigned, FluidSystemNumTemperature);
321 
322  pressureLow_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureLow);
323  pressureHigh_ = EWOMS_GET_PARAM(TypeTag, Scalar, FluidSystemPressureHigh);
324  nPressure_ = EWOMS_GET_PARAM(TypeTag, unsigned, FluidSystemNumPressure);
325 
326  maxDepth_ = EWOMS_GET_PARAM(TypeTag, Scalar, MaxDepth);
327  temperature_ = EWOMS_GET_PARAM(TypeTag, Scalar, Temperature);
328 
329  // initialize the tables of the fluid system
330  // FluidSystem::init();
331  FluidSystem::init(/*Tmin=*/temperatureLow_,
332  /*Tmax=*/temperatureHigh_,
333  /*nT=*/nTemperature_,
334  /*pmin=*/pressureLow_,
335  /*pmax=*/pressureHigh_,
336  /*np=*/nPressure_);
337 
338  fineLayerBottom_ = 22.0;
339 
340  // intrinsic permeabilities
341  fineK_ = this->toDimMatrix_(1e-13);
342  coarseK_ = this->toDimMatrix_(1e-12);
343 
344  // porosities
345  finePorosity_ = 0.3;
346  coarsePorosity_ = 0.3;
347 
348  // residual saturations
349  fineMaterialParams_.setResidualSaturation(liquidPhaseIdx, 0.2);
350  fineMaterialParams_.setResidualSaturation(gasPhaseIdx, 0.0);
351  coarseMaterialParams_.setResidualSaturation(liquidPhaseIdx, 0.2);
352  coarseMaterialParams_.setResidualSaturation(gasPhaseIdx, 0.0);
353 
354  // parameters for the Brooks-Corey law
355  fineMaterialParams_.setEntryPressure(1e4);
356  coarseMaterialParams_.setEntryPressure(5e3);
357  fineMaterialParams_.setLambda(2.0);
358  coarseMaterialParams_.setLambda(2.0);
359 
360  fineMaterialParams_.finalize();
361  coarseMaterialParams_.finalize();
362 
363  // parameters for the somerton law thermal conduction
364  computeThermalCondParams_(fineThermalCondParams_, finePorosity_);
365  computeThermalCondParams_(coarseThermalCondParams_, coarsePorosity_);
366 
367  // assume constant heat capacity and granite
368  solidEnergyLawParams_.setSolidHeatCapacity(790.0 // specific heat capacity of granite [J / (kg K)]
369  * 2700.0); // density of granite [kg/m^3]
370  solidEnergyLawParams_.finalize();
371  }
372 
376  static void registerParameters()
377  {
378  ParentType::registerParameters();
379 
380  EWOMS_REGISTER_PARAM(TypeTag, Scalar, FluidSystemTemperatureLow,
381  "The lower temperature [K] for tabulation of the "
382  "fluid system");
383  EWOMS_REGISTER_PARAM(TypeTag, Scalar, FluidSystemTemperatureHigh,
384  "The upper temperature [K] for tabulation of the "
385  "fluid system");
386  EWOMS_REGISTER_PARAM(TypeTag, unsigned, FluidSystemNumTemperature,
387  "The number of intervals between the lower and "
388  "upper temperature");
389 
390  EWOMS_REGISTER_PARAM(TypeTag, Scalar, FluidSystemPressureLow,
391  "The lower pressure [Pa] for tabulation of the "
392  "fluid system");
393  EWOMS_REGISTER_PARAM(TypeTag, Scalar, FluidSystemPressureHigh,
394  "The upper pressure [Pa] for tabulation of the "
395  "fluid system");
396  EWOMS_REGISTER_PARAM(TypeTag, unsigned, FluidSystemNumPressure,
397  "The number of intervals between the lower and "
398  "upper pressure");
399 
400  EWOMS_REGISTER_PARAM(TypeTag, Scalar, Temperature,
401  "The temperature [K] in the reservoir");
402  EWOMS_REGISTER_PARAM(TypeTag, Scalar, MaxDepth,
403  "The maximum depth [m] of the reservoir");
404  EWOMS_REGISTER_PARAM(TypeTag, std::string, SimulationName,
405  "The name of the simulation used for the output "
406  "files");
407  }
408 
413 
417  std::string name() const
418  {
419  std::ostringstream oss;
420  oss << EWOMS_GET_PARAM(TypeTag, std::string, SimulationName)
421  << "_" << Model::name();
422  if (getPropValue<TypeTag, Properties::EnableEnergy>())
423  oss << "_ni";
424  oss << "_" << Model::discretizationName();
425  return oss.str();
426  }
427 
431  void endTimeStep()
432  {
433 #ifndef NDEBUG
434  Scalar tol = this->model().newtonMethod().tolerance()*1e5;
435  this->model().checkConservativeness(tol);
436 
437  // Calculate storage terms
438  PrimaryVariables storageL, storageG;
439  this->model().globalPhaseStorage(storageL, /*phaseIdx=*/0);
440  this->model().globalPhaseStorage(storageG, /*phaseIdx=*/1);
441 
442  // Write mass balance information for rank 0
443  if (this->gridView().comm().rank() == 0) {
444  std::cout << "Storage: liquid=[" << storageL << "]"
445  << " gas=[" << storageG << "]\n" << std::flush;
446  }
447 #endif // NDEBUG
448  }
449 
453  template <class Context>
454  Scalar temperature(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
455  {
456  const auto& pos = context.pos(spaceIdx, timeIdx);
457  if (inHighTemperatureRegion_(pos))
458  return temperature_ + 100;
459  return temperature_;
460  }
461 
465  template <class Context>
466  const DimMatrix& intrinsicPermeability(const Context& context, unsigned spaceIdx,
467  unsigned timeIdx) const
468  {
469  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
470  if (isFineMaterial_(pos))
471  return fineK_;
472  return coarseK_;
473  }
474 
478  template <class Context>
479  Scalar porosity(const Context& context, unsigned spaceIdx, unsigned timeIdx) const
480  {
481  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
482  if (isFineMaterial_(pos))
483  return finePorosity_;
484  return coarsePorosity_;
485  }
486 
490  template <class Context>
491  const MaterialLawParams& materialLawParams(const Context& context,
492  unsigned spaceIdx, unsigned timeIdx) const
493  {
494  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
495  if (isFineMaterial_(pos))
496  return fineMaterialParams_;
497  return coarseMaterialParams_;
498  }
499 
505  template <class Context>
506  const SolidEnergyLawParams&
507  solidEnergyLawParams(const Context& /*context*/,
508  unsigned /*spaceIdx*/,
509  unsigned /*timeIdx*/) const
510  { return solidEnergyLawParams_; }
511 
515  template <class Context>
516  const ThermalConductionLawParams &
517  thermalConductionLawParams(const Context& context,
518  unsigned spaceIdx,
519  unsigned timeIdx) const
520  {
521  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
522  if (isFineMaterial_(pos))
523  return fineThermalCondParams_;
524  return coarseThermalCondParams_;
525  }
526 
528 
533 
537  template <class Context>
538  void boundary(BoundaryRateVector& values, const Context& context,
539  unsigned spaceIdx, unsigned timeIdx) const
540  {
541  const auto& pos = context.pos(spaceIdx, timeIdx);
542  if (onLeftBoundary_(pos)) {
543  Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
544  initialFluidState_(fs, context, spaceIdx, timeIdx);
545  fs.checkDefined();
546 
547  // impose an freeflow boundary condition
548  values.setFreeFlow(context, spaceIdx, timeIdx, fs);
549  }
550  else if (onInlet_(pos)) {
551  RateVector massRate(0.0);
552  massRate[contiCO2EqIdx] = -1e-3; // [kg/(m^3 s)]
553 
554  using FluidState = Opm::ImmiscibleFluidState<Scalar, FluidSystem>;
555  FluidState fs;
556  fs.setSaturation(gasPhaseIdx, 1.0);
557  const auto& pg =
558  context.intensiveQuantities(spaceIdx, timeIdx).fluidState().pressure(gasPhaseIdx);
559  fs.setPressure(gasPhaseIdx, Toolbox::value(pg));
560  fs.setTemperature(temperature(context, spaceIdx, timeIdx));
561 
562  typename FluidSystem::template ParameterCache<Scalar> paramCache;
563  paramCache.updatePhase(fs, gasPhaseIdx);
564  Scalar h = FluidSystem::template enthalpy<FluidState, Scalar>(fs, paramCache, gasPhaseIdx);
565 
566  // impose an forced inflow boundary condition for pure CO2
567  values.setMassRate(massRate);
568  values.setEnthalpyRate(massRate[contiCO2EqIdx] * h);
569  }
570  else
571  // no flow on top and bottom
572  values.setNoFlow();
573  }
574 
575  // \}
576 
581 
585  template <class Context>
586  void initial(PrimaryVariables& values, const Context& context, unsigned spaceIdx,
587  unsigned timeIdx) const
588  {
589  Opm::CompositionalFluidState<Scalar, FluidSystem> fs;
590  initialFluidState_(fs, context, spaceIdx, timeIdx);
591 
592  // const auto& matParams = this->materialLawParams(context, spaceIdx,
593  // timeIdx);
594  // values.assignMassConservative(fs, matParams, /*inEquilibrium=*/true);
595  values.assignNaive(fs);
596  }
597 
604  template <class Context>
605  void source(RateVector& rate,
606  const Context& /*context*/,
607  unsigned /*spaceIdx*/,
608  unsigned /*timeIdx*/) const
609  { rate = Scalar(0.0); }
610 
612 
613 private:
614  template <class Context, class FluidState>
615  void initialFluidState_(FluidState& fs,
616  const Context& context,
617  unsigned spaceIdx,
618  unsigned timeIdx) const
619  {
620  const GlobalPosition& pos = context.pos(spaceIdx, timeIdx);
621 
623  // set temperature
625  fs.setTemperature(temperature(context, spaceIdx, timeIdx));
626 
628  // set saturations
630  fs.setSaturation(FluidSystem::liquidPhaseIdx, 1.0);
631  fs.setSaturation(FluidSystem::gasPhaseIdx, 0.0);
632 
634  // set pressures
636  Scalar densityL = FluidSystem::Brine::liquidDensity(temperature_, Scalar(1e5));
637  Scalar depth = maxDepth_ - pos[dim - 1];
638  Scalar pl = 1e5 - densityL * this->gravity()[dim - 1] * depth;
639 
640  Scalar pC[numPhases];
641  const auto& matParams = this->materialLawParams(context, spaceIdx, timeIdx);
642  MaterialLaw::capillaryPressures(pC, matParams, fs);
643 
644  fs.setPressure(liquidPhaseIdx, pl + (pC[liquidPhaseIdx] - pC[liquidPhaseIdx]));
645  fs.setPressure(gasPhaseIdx, pl + (pC[gasPhaseIdx] - pC[liquidPhaseIdx]));
646 
648  // set composition of the liquid phase
650  fs.setMoleFraction(liquidPhaseIdx, CO2Idx, 0.005);
651  fs.setMoleFraction(liquidPhaseIdx, BrineIdx,
652  1.0 - fs.moleFraction(liquidPhaseIdx, CO2Idx));
653 
654  typename FluidSystem::template ParameterCache<Scalar> paramCache;
655  using CFRP = Opm::ComputeFromReferencePhase<Scalar, FluidSystem>;
656  CFRP::solve(fs, paramCache,
657  /*refPhaseIdx=*/liquidPhaseIdx,
658  /*setViscosity=*/true,
659  /*setEnthalpy=*/true);
660  }
661 
662  bool onLeftBoundary_(const GlobalPosition& pos) const
663  { return pos[0] < eps_; }
664 
665  bool onRightBoundary_(const GlobalPosition& pos) const
666  { return pos[0] > this->boundingBoxMax()[0] - eps_; }
667 
668  bool onInlet_(const GlobalPosition& pos) const
669  { return onRightBoundary_(pos) && (5 < pos[1]) && (pos[1] < 15); }
670 
671  bool inHighTemperatureRegion_(const GlobalPosition& pos) const
672  { return (pos[0] > 20) && (pos[0] < 30) && (pos[1] > 5) && (pos[1] < 35); }
673 
674  void computeThermalCondParams_(ThermalConductionLawParams& params, Scalar poro)
675  {
676  Scalar lambdaWater = 0.6;
677  Scalar lambdaGranite = 2.8;
678 
679  Scalar lambdaWet = std::pow(lambdaGranite, (1 - poro))
680  * std::pow(lambdaWater, poro);
681  Scalar lambdaDry = std::pow(lambdaGranite, (1 - poro));
682 
683  params.setFullySaturatedLambda(gasPhaseIdx, lambdaDry);
684  params.setFullySaturatedLambda(liquidPhaseIdx, lambdaWet);
685  params.setVacuumLambda(lambdaDry);
686  }
687 
688  bool isFineMaterial_(const GlobalPosition& pos) const
689  { return pos[dim - 1] > fineLayerBottom_; }
690 
691  DimMatrix fineK_;
692  DimMatrix coarseK_;
693  Scalar fineLayerBottom_;
694 
695  Scalar finePorosity_;
696  Scalar coarsePorosity_;
697 
698  MaterialLawParams fineMaterialParams_;
699  MaterialLawParams coarseMaterialParams_;
700 
701  ThermalConductionLawParams fineThermalCondParams_;
702  ThermalConductionLawParams coarseThermalCondParams_;
703  SolidEnergyLawParams solidEnergyLawParams_;
704 
705  Scalar temperature_;
706  Scalar maxDepth_;
707  Scalar eps_;
708 
709  unsigned nTemperature_;
710  unsigned nPressure_;
711 
712  Scalar pressureLow_, pressureHigh_;
713  Scalar temperatureLow_, temperatureHigh_;
714 };
715 } // namespace Opm
716 
717 #endif
Problem where is injected under a low permeable layer at a depth of 2700m.
Definition: co2injectionproblem.hh:264
void finishInit()
Definition: co2injectionproblem.hh:312
Scalar porosity(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:479
void boundary(BoundaryRateVector &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:538
void initial(PrimaryVariables &values, const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:586
const MaterialLawParams & materialLawParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:491
void source(RateVector &rate, const Context &, unsigned, unsigned) const
Definition: co2injectionproblem.hh:605
std::string name() const
Definition: co2injectionproblem.hh:417
void endTimeStep()
Definition: co2injectionproblem.hh:431
const SolidEnergyLawParams & solidEnergyLawParams(const Context &, unsigned, unsigned) const
Return the parameters for the heat storage law of the rock.
Definition: co2injectionproblem.hh:507
const ThermalConductionLawParams & thermalConductionLawParams(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:517
Scalar temperature(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:454
Co2InjectionProblem(Simulator &simulator)
Definition: co2injectionproblem.hh:305
static void registerParameters()
Definition: co2injectionproblem.hh:376
const DimMatrix & intrinsicPermeability(const Context &context, unsigned spaceIdx, unsigned timeIdx) const
Definition: co2injectionproblem.hh:466
Definition: co2injectionproblem.hh:82
Definition: co2injectionproblem.hh:88
Definition: co2injectionproblem.hh:80
Definition: co2injectionproblem.hh:78
Definition: co2injectionproblem.hh:86
Definition: co2injectionproblem.hh:84
Definition: co2injectionproblem.hh:91
Definition: co2injectionproblem.hh:95
Definition: co2injectionproblem.hh:73
Definition: co2injectionproblem.hh:93