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