All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
BlackoilCo2PVT.hpp
1 /*
2  Copyright 2010 SINTEF ICT, Applied Mathematics.
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 3 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 
20 #ifndef OPM_BLACKOILCO2PVT_HEADER_INCLUDED
21 #define OPM_BLACKOILCO2PVT_HEADER_INCLUDED
22 #define OPM_BLACKOILPVT_HEADER_INCLUDED
23 
24 #define OPM_DEPRECATED __attribute__((deprecated))
25 #define OPM_DEPRECATED_MSG(msg) __attribute__((deprecated))
26 
27 #include <opm/material/fluidsystems/BrineCO2FluidSystem.hpp>
28 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
29 #include <opm/material/constraintsolvers/MiscibleMultiPhaseComposition.hpp>
30 #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
31 
32 #include <opm/common/Exceptions.hpp>
33 #include <opm/common/ErrorMacros.hpp>
35 #include <opm/porsol/blackoil/fluid/MiscibilityProps.hpp>
36 #include <opm/porsol/blackoil/fluid/BlackoilDefs.hpp>
37 
38 #include <opm/parser/eclipse/Deck/Deck.hpp>
39 
40 #include <string>
41 #include <iostream>
42 #include <fstream>
43 #include <iomanip>
44 
45 
46 namespace Opm
47 {
49 {
50 public:
51 
52  typedef Opm::FluidSystems::BrineCO2</*Scalar=*/double, Opm::Benchmark3::CO2Tables> FluidSystem;
53  typedef Opm::CompositionalFluidState<double, FluidSystem> CompositionalFluidState;
54 
55  void init(const Opm::Deck& deck);
56 
57  void generateBlackOilTables(double temperature);
58 
59  double getViscosity(double press,
60  const CompVec& surfvol,
61  PhaseIndex phase) const;
62  CompVec surfaceDensities() const;
63  double getSaturation(double press,
64  const CompVec& surfvol,
65  PhaseIndex phase) const;
66  double B (double press,
67  const CompVec& surfvol,
68  PhaseIndex phase) const;
69  double dBdp(double press,
70  const CompVec& surfvol,
71  PhaseIndex phase) const;
72  double R (double press,
73  const CompVec& surfvol,
74  PhaseIndex phase) const;
75  double dRdp(double press,
76  const CompVec& surfvol,
77  PhaseIndex phase) const;
78 
79  void getViscosity(const std::vector<PhaseVec>& pressures,
80  const std::vector<CompVec>& surfvol,
81  std::vector<PhaseVec>& output) const;
82  void B(const std::vector<PhaseVec>& pressures,
83  const std::vector<CompVec>& surfvol,
84  std::vector<PhaseVec>& output) const;
85  void dBdp(const std::vector<PhaseVec>& pressures,
86  const std::vector<CompVec>& surfvol,
87  std::vector<PhaseVec>& output_B,
88  std::vector<PhaseVec>& output_dBdp) const;
89  void R(const std::vector<PhaseVec>& pressures,
90  const std::vector<CompVec>& surfvol,
91  std::vector<PhaseVec>& output) const;
92  void dRdp(const std::vector<PhaseVec>& pressures,
93  const std::vector<CompVec>& surfvol,
94  std::vector<PhaseVec>& output_R,
95  std::vector<PhaseVec>& output_dRdp) const;
96 
97 private:
98  CompVec surfaceDensities_;
99  FluidSystem brineCo2_;
100  double temperature_;
101  struct SubState
102  {
103  Dune::FieldVector<double,2> density;
104  Dune::FieldMatrix<double,2,2> massfrac;
105  Dune::FieldVector<double,2> phaseVolume;
106  Dune::FieldVector<double,2> phaseViscosity;
107  double saturation;
108  };
109  void computeState(SubState& ss, double zBrine, double zCO2, double pressure) const;
110  enum {
111  wPhase = FluidSystem::liquidPhaseIdx,
112  nPhase = FluidSystem::gasPhaseIdx,
113 
114  wComp = FluidSystem::BrineIdx,
115  nComp = FluidSystem::CO2Idx
116  };
117 
118 }; // class BlackoilCo2PVT
119 
120 // ------------ Method implementations --------------
121 
122 void BlackoilCo2PVT::init(const Opm::Deck& /* deck */)
123 {
124  surfaceDensities_[Water] = 1000.;
125  surfaceDensities_[Gas] = 2.0;
126  surfaceDensities_[Oil] = 1000.;
127 
128  temperature_ = 300.;
129 
130  brineCo2_.init();
131 }
132 
133 void BlackoilCo2PVT::generateBlackOilTables(double temperature)
134 {
135  std::cout << "\n Generating pvt tables for the eclipse black oil formulation\n using the oil component as brine and the gas component as co_2." << std::endl;
136  if (std::fabs(temperature-400.) > 100.0) {
137  std::cout << "T=" << temperature << " is outside the allowed range [300,500] Kelvin" << std::endl;
138  exit(-1);
139  }
140 
141  temperature_ = temperature;
142 
143  CompVec z;
144  z[Water] = 0.0;
145  z[Oil] = 1.0;
146  z[Gas] = 1.0e6;
147 
148  std::ofstream black("blackoil_pvt");
149  black.precision(8);
150  black << std::fixed << std::showpoint;
151 
152  const double pMin=150.e5;
153  const double pMax=400.e5;
154  const unsigned int np=11;
155  std::vector<double> pValue(np+1,0.0);
156  std::vector<double> rs(np+1,0.0);
157 
158  pValue[0] = 101330.0;
159  rs[0] = 0.0;
160 
161  // Buble points
162  z[Gas] = 1.0e4;
163  for (unsigned int i=0; i<np; ++i) {
164  pValue[i+1] = pMin + i*(pMax-pMin)/(np-1);
165  rs[i+1] = R(pValue[i+1], z, Liquid);
166  }
167 
168  const unsigned int np_fill_in = 10;
169  const double dr = (rs[1] - rs[0])/(np_fill_in+1);
170  const double dp = (pValue[1] - pValue[0])/(np_fill_in+1);
171  rs.insert(rs.begin()+1,np_fill_in,0.0);
172  pValue.insert(pValue.begin()+1,np_fill_in,0.0);
173  for (unsigned int i=1; i<=np_fill_in; ++i) {
174  rs[i] = rs[i-1] + dr;
175  pValue[i] = pValue[i-1] + dp;
176  }
177 
178  // Brine with dissolved co_2 ("live oil")
179  black << "PVTO\n";
180  black << "-- Rs Pbub Bo Vo\n";
181  black << "-- sm3/sm3 barsa rm3/sm3 cP\n";
182 
183  // Undersaturated
184  for (unsigned int i=0; i<np+np_fill_in; ++i) {
185  z[Gas] = rs[i];
186  black << std::setw(14) << rs[i];
187  for (unsigned int j=i; j<np+1+np_fill_in; ++j) {
188  if (j<=np_fill_in) {
189  if (j==i) black << std::setw(14) << pValue[j]*1.e-5 << std::setw(14) << 1.0-j*0.001 << std::setw(14) << 1.06499;
190  continue;
191  }
192  if (j>i) black << std::endl << std::setw(14) << ' ';
193  black << std::setw(14) << pValue[j]*1.e-5
194  << std::setw(14) << B(pValue[j], z, Liquid)
195  << std::setw(14) << getViscosity(pValue[j], z, Liquid)*1.e3;
196  }
197  black << " /" << std::endl;
198  }
199  black << "/ " << std::endl;
200 
201  // We provide tables for co_2 both with and without dissolved water:
202 
203  // Co_2 neglecting dissolved water ("dry gas")
204  black << "\nPVDG\n";
205  black << "-- Pg Bg Vg\n";
206  black << "-- barsa rm3/sm3 cP\n";
207 
208  for (unsigned int i=0; i<np; ++i) {
209  z[Oil] = 0.0;
210  z[Gas] = 1.0;
211  black << std::setw(14) << pValue[i+np_fill_in+1]*1.e-5
212  << std::setw(14) << B(pValue[i+np_fill_in+1], z, Vapour)
213  << std::setw(14) << getViscosity(pValue[i+np_fill_in+1], z, Vapour)*1.e3
214  << std::endl;
215  }
216  black << "/ " << std::endl;
217 
218  // Co_2 with dissolved water ("wet gas")
219  black << "\nPVTG\n";
220  black << "-- Pg Rv Bg Vg\n";
221  black << "-- barsa sm3/sm3 rm3/sm3 cP\n";
222  for (unsigned int i=0; i<np; ++i) {
223  z[Oil] = 1000.0;
224  z[Gas] = 1.0;
225  black << std::setw(14) << pValue[i+np_fill_in+1]*1.e-5
226  << std::setw(14) << R(pValue[i+np_fill_in+1], z, Vapour)
227  << std::setw(14) << B(pValue[i+np_fill_in+1], z, Vapour)
228  << std::setw(14) << getViscosity(pValue[i+np_fill_in+1], z, Vapour)*1.e3
229  << std::endl;
230  z[Oil] = 0.0;
231  black << std::setw(14) << ' '
232  << std::setw(14) << R(pValue[i+np_fill_in+1], z, Vapour)
233  << std::setw(14) << B(pValue[i+np_fill_in+1], z, Vapour)
234  << std::setw(14) << getViscosity(pValue[i+np_fill_in+1], z, Vapour)*1.e3
235  << " /" << std::endl;
236  }
237  black << "/ " << std::endl;
238  black << std::endl;
239  std::cout << " Pvt tables for temperature=" << temperature << " Kelvin is written to file blackoil_pvt. " << std::endl;
240  std::cout << " NOTE that the file contains tables for both PVDG (dry gas) and PVTG (wet gas)." << std::endl;
241 }
242 
243 BlackoilCo2PVT::CompVec BlackoilCo2PVT::surfaceDensities() const
244 {
245  return surfaceDensities_;
246 }
247 
248 double BlackoilCo2PVT::getViscosity(double press, const CompVec& surfvol, PhaseIndex phase) const
249 {
250  SubState ss;
251  computeState(ss, surfvol[Oil], surfvol[Gas], press);
252  switch(phase) {
253  case Aqua: return 1.0e-10;
254  case Liquid: return ss.phaseViscosity[wPhase];
255  case Vapour: return ss.phaseViscosity[nPhase];
256  };
257  return 1.0e-10;
258 }
259 
260 
261 double BlackoilCo2PVT::getSaturation(double press, const CompVec& surfvol, PhaseIndex phase) const
262 {
263  SubState ss;
264  computeState(ss, surfvol[Oil], surfvol[Gas], press);
265  switch(phase) {
266  case Aqua: return 0.0;
267  case Liquid: return ss.saturation;
268  case Vapour: return 1.0 - ss.saturation;
269  };
270  return 0.0;
271 }
272 
273 
274 double BlackoilCo2PVT::B(double press, const CompVec& surfvol, PhaseIndex phase) const
275 {
276  SubState ss;
277  computeState(ss, surfvol[Oil], surfvol[Gas], press);
278  switch(phase) {
279  case Aqua: return 1.0;
280  case Liquid: return surfaceDensities_[Oil]/(ss.massfrac[wPhase][wComp]*ss.density[wPhase]+1.0e-10);
281  case Vapour: return surfaceDensities_[Gas]/(ss.massfrac[nPhase][nComp]*ss.density[nPhase]+1.0e-10);
282  };
283  return 1.0;
284 }
285 
286 double BlackoilCo2PVT::dBdp(double press, const CompVec& surfvol, PhaseIndex phase) const
287 {
288  const double dp = 100.;
289  return (B(press+dp,surfvol,phase)-B(press,surfvol,phase))/dp;
290 }
291 
292 double BlackoilCo2PVT::R(double press, const CompVec& surfvol, PhaseIndex phase) const
293 {
294  SubState ss;
295  computeState(ss, surfvol[Oil], surfvol[Gas], press);
296  switch(phase) {
297  case Aqua: return 0.0;
298  case Liquid: return (ss.massfrac[wPhase][nComp]*surfaceDensities_[Oil])/(ss.massfrac[wPhase][wComp]*surfaceDensities_[Gas]+1.0e-10);
299  case Vapour: return (ss.massfrac[nPhase][wComp]*surfaceDensities_[Gas])/(ss.massfrac[nPhase][nComp]*surfaceDensities_[Oil]+1.0e-10);
300  };
301  return 0.0;
302 }
303 
304 double BlackoilCo2PVT::dRdp(double press, const CompVec& surfvol, PhaseIndex phase) const
305 {
306  const double dp = 100.;
307  return (R(press+dp,surfvol,phase)-R(press,surfvol,phase))/dp;
308 }
309 
310 void BlackoilCo2PVT::getViscosity(const std::vector<PhaseVec>& pressures,
311  const std::vector<CompVec>& surfvol,
312  std::vector<PhaseVec>& output) const
313 {
314  int num = pressures.size();
315  output.resize(num);
316  SubState ss;
317  for (int i = 0; i < num; ++i) {
318  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]);
319  output[i][Aqua] = 1.0e-10;
320  output[i][Liquid] = ss.phaseViscosity[wPhase];
321  output[i][Vapour] = ss.phaseViscosity[nPhase];
322  }
323 }
324 
325 void BlackoilCo2PVT::B(const std::vector<PhaseVec>& pressures,
326  const std::vector<CompVec>& surfvol,
327  std::vector<PhaseVec>& output) const
328 {
329  int num = pressures.size();
330  output.resize(num);
331  SubState ss;
332  for (int i = 0; i < num; ++i) {
333  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]);
334  output[i][Aqua] = 1.0;
335  output[i][Liquid] = surfaceDensities_[Oil]/(ss.massfrac[wPhase][wComp]*ss.density[wPhase]+1.0e-10);
336  output[i][Vapour] = surfaceDensities_[Gas]/(ss.massfrac[nPhase][nComp]*ss.density[nPhase]+1.0e-10);
337  }
338 }
339 
340 void BlackoilCo2PVT::dBdp(const std::vector<PhaseVec>& pressures,
341  const std::vector<CompVec>& surfvol,
342  std::vector<PhaseVec>& output_B,
343  std::vector<PhaseVec>& output_dBdp) const
344 {
345  int num = pressures.size();
346  output_B.resize(num);
347  output_dBdp.resize(num);
348  SubState ss;
349  const double dp = 100.;
350  for (int i = 0; i < num; ++i) {
351  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]);
352  output_B[i][Aqua] = 1.0;
353  output_B[i][Liquid] = surfaceDensities_[Oil]/(ss.massfrac[wPhase][wComp]*ss.density[wPhase]+1.0e-10);
354  output_B[i][Vapour] = surfaceDensities_[Gas]/(ss.massfrac[nPhase][nComp]*ss.density[nPhase]+1.0e-10);
355  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]+dp);
356  output_dBdp[i][Aqua] = 0.0;
357  output_dBdp[i][Liquid] = (surfaceDensities_[Oil]/(ss.massfrac[wPhase][wComp]*ss.density[wPhase]+1.0e-10) - output_B[i][Liquid])/dp;
358  output_dBdp[i][Vapour] = (surfaceDensities_[Gas]/(ss.massfrac[nPhase][nComp]*ss.density[nPhase]+1.0e-10) - output_B[i][Vapour])/dp;
359  }
360 }
361 
362 void BlackoilCo2PVT::R(const std::vector<PhaseVec>& pressures,
363  const std::vector<CompVec>& surfvol,
364  std::vector<PhaseVec>& output) const
365 {
366  int num = pressures.size();
367  output.resize(num);
368  SubState ss;
369  for (int i = 0; i < num; ++i) {
370  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]);
371  output[i][Aqua] = 0.0;
372  output[i][Liquid] = (ss.massfrac[wPhase][nComp]*surfaceDensities_[Oil])/(ss.massfrac[wPhase][wComp]*surfaceDensities_[Gas]+1.0e-10);
373  output[i][Vapour] = (ss.massfrac[nPhase][wComp]*surfaceDensities_[Gas])/(ss.massfrac[nPhase][nComp]*surfaceDensities_[Oil]+1.0e-10);
374  }
375 }
376 
377 void BlackoilCo2PVT::dRdp(const std::vector<PhaseVec>& pressures,
378  const std::vector<CompVec>& surfvol,
379  std::vector<PhaseVec>& output_R,
380  std::vector<PhaseVec>& output_dRdp) const
381 {
382  int num = pressures.size();
383  output_R.resize(num);
384  output_dRdp.resize(num);
385  SubState ss;
386  const double dp = 100.;
387  for (int i = 0; i < num; ++i) {
388  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]);
389  output_R[i][Aqua] = 0.0;
390  output_R[i][Liquid] = (ss.massfrac[wPhase][nComp]*surfaceDensities_[Oil])/(ss.massfrac[wPhase][wComp]*surfaceDensities_[Gas]+1.0e-10);
391  output_R[i][Vapour] = (ss.massfrac[nPhase][wComp]*surfaceDensities_[Gas])/(ss.massfrac[nPhase][nComp]*surfaceDensities_[Oil]+1.0e-10);
392  computeState(ss, surfvol[i][Oil], surfvol[i][Gas], pressures[i][Liquid]+dp);
393  output_dRdp[i][Aqua] = 0.0;
394  output_dRdp[i][Liquid] = ((ss.massfrac[wPhase][nComp]*surfaceDensities_[Oil])/(ss.massfrac[wPhase][wComp]*surfaceDensities_[Gas]+1.0e-10) - output_R[i][Liquid])/dp;
395  output_dRdp[i][Vapour] = ((ss.massfrac[nPhase][wComp]*surfaceDensities_[Gas])/(ss.massfrac[nPhase][nComp]*surfaceDensities_[Oil]+1.0e-10) - output_R[i][Vapour])/dp;
396  }
397 }
398 
399 void BlackoilCo2PVT::computeState(BlackoilCo2PVT::SubState& ss, double zBrine, double zCO2, double pressure) const
400 {
401 
402  CompositionalFluidState fluidState;
403  fluidState.setTemperature(temperature_);
404  fluidState.setPressure(wPhase, pressure);
405  fluidState.setPressure(nPhase, pressure);
406 
407  double massH20 = surfaceDensities_[Oil]*zBrine;
408  double massCO2 = surfaceDensities_[Gas]*zCO2;
409 
410  // A priori, assume presence of both phases
411  FluidSystem::ParameterCache<double> paramCache;
412  typedef Opm::MiscibleMultiPhaseComposition</*Scalar=*/double, FluidSystem> MMPC;
413  MMPC::solve(fluidState, paramCache, /*setViscosity=*/false, /*setEnthalpy=*/false);
414  ss.density[wPhase] = fluidState.density(wPhase);
415  ss.density[nPhase] = fluidState.density(nPhase);
416  ss.massfrac[wPhase][nComp] = fluidState.massFraction(wPhase, nComp);
417  ss.massfrac[nPhase][wComp] = fluidState.massFraction(nPhase, wComp);
418  ss.massfrac[wPhase][wComp] = 1.0 - ss.massfrac[wPhase][nComp];
419  ss.massfrac[nPhase][nComp] = 1.0 - ss.massfrac[nPhase][wComp];
420 
421  double detX = ss.massfrac[wPhase][wComp]*ss.massfrac[nPhase][nComp]-ss.massfrac[wPhase][nComp]*ss.massfrac[nPhase][wComp];
422  ss.phaseVolume[wPhase] = (massH20*ss.massfrac[nPhase][nComp] - massCO2*ss.massfrac[nPhase][wComp])/(ss.density[wPhase]*detX);
423  ss.phaseVolume[nPhase] = (massCO2*ss.massfrac[wPhase][wComp] - massH20*ss.massfrac[wPhase][nComp])/(ss.density[nPhase]*detX);
424 
425  // Determine number of phase
426  if (ss.phaseVolume[wPhase] > 0.0 && ss.phaseVolume[nPhase] > 0.0) { // Both phases
427  ss.saturation = ss.phaseVolume[wPhase]/(ss.phaseVolume[wPhase]+ss.phaseVolume[nPhase]);
428  fluidState.setSaturation(wPhase, ss.saturation);
429  fluidState.setSaturation(nPhase, 1.0 - ss.saturation);
430  }
431  else if (ss.phaseVolume[wPhase] <= 0.0) { // Wetting phase only
432  ss.saturation = 0.0;
433  // Gas phase:
434  ss.massfrac[nPhase][nComp] = massCO2/(massCO2+massH20);
435  ss.massfrac[nPhase][wComp] = 1.0 - ss.massfrac[nPhase][nComp];
436  double M1 = FluidSystem::molarMass(wComp);
437  double M2 = FluidSystem::molarMass(nComp);
438  double avgMolarMass = M1*M2/(M2 + ss.massfrac[nPhase][nComp]*(M1 - M2));
439  fluidState.setMoleFraction(nPhase, nComp, ss.massfrac[nPhase][nComp]*avgMolarMass/M2);
440  fluidState.setMoleFraction(nPhase, wComp, ss.massfrac[nPhase][wComp]*avgMolarMass/M1);
441  ss.density[nPhase] = brineCo2_.density(fluidState, paramCache, nPhase);
442  fluidState.setDensity(nPhase, ss.density[nPhase]);
443  ss.phaseVolume[nPhase] = (massH20+massCO2)/ss.density[nPhase];
444  fluidState.setSaturation(nPhase, 1.0 - ss.saturation);
445  // Virtual properties of non-existing liquid phase:
446  paramCache.updatePhase(fluidState, /*phaseIdx=*/nPhase);
447  typedef Opm::ComputeFromReferencePhase</*Scalar=*/double, FluidSystem> CFRP;
448  CFRP::solve(fluidState,
449  paramCache,
450  /*refPhaseIdx=*/nPhase,
451  /*setViscosity=*/false,
452  /*setEnthalpy=*/false);
453  ss.massfrac[wPhase][wComp] = fluidState.massFraction(wPhase, wComp);
454  ss.massfrac[wPhase][nComp] = fluidState.massFraction(wPhase, nComp);
455  ss.density[wPhase] = fluidState.density(wPhase);
456  ss.phaseVolume[wPhase] = 0.0;
457  fluidState.setSaturation(wPhase, ss.saturation);
458  }
459  else if (ss.phaseVolume[nPhase] <= 0.0) { // Non-wetting phase only
460  ss.saturation = 1.0;
461  // Liquid phase:
462  ss.massfrac[wPhase][wComp] = massH20/(massCO2+massH20);
463  ss.massfrac[wPhase][nComp] = 1.0 - ss.massfrac[wPhase][wComp];
464  double M1 = FluidSystem::molarMass(wComp);
465  double M2 = FluidSystem::molarMass(nComp);
466  double avgMolarMass = M1*M2/(M2 + ss.massfrac[wPhase][nComp]*(M1 - M2));
467  fluidState.setMoleFraction(wPhase, nComp, ss.massfrac[wPhase][nComp]*avgMolarMass/M2);
468  fluidState.setMoleFraction(wPhase, wComp, ss.massfrac[wPhase][wComp]*avgMolarMass/M1);
469  ss.density[wPhase] = brineCo2_.density(fluidState, paramCache, wPhase);
470  fluidState.setDensity(wPhase, ss.density[wPhase]);
471  ss.phaseVolume[wPhase] = (massH20+massCO2)/ss.density[wPhase];
472  fluidState.setSaturation(wPhase, ss.saturation);
473  // Virtual properties of non-existing gas phase:
474  paramCache.updatePhase(fluidState, /*phaseIdx=*/nPhase);
475  typedef ComputeFromReferencePhase</*Scalar=*/double, FluidSystem> CFRP;
476  CFRP::solve(fluidState,
477  paramCache,
478  /*refPhaseIdx=*/wPhase,
479  /*setViscosity=*/false,
480  /*setEnthalpy=*/false);
481  ss.massfrac[nPhase][nComp] = fluidState.massFraction(nPhase, nComp);
482  ss.massfrac[nPhase][wComp] = fluidState.massFraction(nPhase, wComp);
483  ss.density[nPhase] = fluidState.density(nPhase);
484  ss.phaseVolume[nPhase] = 0.0;
485  fluidState.setSaturation(nPhase, 1.0 - ss.saturation);
486  }
487 
488  ss.phaseViscosity[wPhase] = brineCo2_.viscosity(fluidState, paramCache, wPhase);
489  ss.phaseViscosity[nPhase] = brineCo2_.viscosity(fluidState, paramCache, nPhase);
490 
491 }
492 
493 
494 
495 typedef BlackoilCo2PVT BlackoilPVT;
496 
497 } // Opm
498 
499 
500 #endif // OPM_BLACKOILCO2PVT_HEADER_INCLUDED
Provides the class with the tabulated values of CO2 for the benchmark3 problem.
Definition: BlackoilDefs.hpp:33
Definition: benchmark3co2tables.hh:25261
Definition: BlackoilCo2PVT.hpp:48