compiled and running now.

This commit is contained in:
Kai Bao 2021-11-12 13:07:01 +01:00 committed by Trine Mykkeltvedt
parent 52f23ec882
commit df50b67ed3
3 changed files with 56 additions and 25 deletions

View File

@ -191,10 +191,13 @@ public:
// the EOS has three intersections with the pressure, // the EOS has three intersections with the pressure,
// i.e. the molar volume of gas is the largest one and the // i.e. the molar volume of gas is the largest one and the
// molar volume of liquid is the smallest one // molar volume of liquid is the smallest one
#warning HACK, should investigate why
if (isGasPhase) if (isGasPhase)
Vm = Z[2]*RT/p; // Vm = Z[2]*RT/p;
Vm = max(1e-7, Z[2]*RT/p);
else else
Vm = Z[0]*RT/p; // Vm = Z[2]*RT/p;
Vm = max(1e-7, Z[0]*RT/p);
} }
else if (numSol == 1) { else if (numSol == 1) {
// the EOS only has one intersection with the pressure, // the EOS only has one intersection with the pressure,
@ -203,27 +206,28 @@ public:
Evaluation VmCubic = Z[0]*RT/p; Evaluation VmCubic = Z[0]*RT/p;
Vm = VmCubic; Vm = VmCubic;
#warning, should investigate why here
// find the extrema (if they are present) // find the extrema (if they are present)
Evaluation Vmin, Vmax, pmin, pmax; // Evaluation Vmin, Vmax, pmin, pmax;
if (findExtrema_(Vmin, Vmax, // if (findExtrema_(Vmin, Vmax,
pmin, pmax, // pmin, pmax,
a, b, T)) // a, b, T))
{ // {
if (isGasPhase) // if (isGasPhase)
Vm = std::max(Vmax, VmCubic); // Vm = std::max(Vmax, VmCubic);
else { // else {
if (Vmin > 0) // if (Vmin > 0)
Vm = std::min(Vmin, VmCubic); // Vm = std::min(Vmin, VmCubic);
else // else
Vm = VmCubic; // Vm = VmCubic;
} // }
} // }
else { // else {
// the EOS does not exhibit any physically meaningful // // the EOS does not exhibit any physically meaningful
// extrema, and the fluid is critical... // // extrema, and the fluid is critical...
Vm = VmCubic; // Vm = VmCubic;
handleCriticalFluid_(Vm, fs, params, phaseIdx, isGasPhase); // handleCriticalFluid_(Vm, fs, params, phaseIdx, isGasPhase);
} // }
} }
Valgrind::CheckDefined(Vm); Valgrind::CheckDefined(Vm);

View File

@ -87,12 +87,28 @@ public:
* R. Reid, et al.: The Properties of Gases and Liquids, * R. Reid, et al.: The Properties of Gases and Liquids,
* 4th edition, McGraw-Hill, 1987, pp. 42-44, 143-145 * 4th edition, McGraw-Hill, 1987, pp. 42-44, 143-145
*/ */
#warning should check why this function is changed
template <class FluidState, class Params, class LhsEval = typename FluidState::Scalar> template <class FluidState, class Params, class LhsEval = typename FluidState::Scalar>
static LhsEval computeFugacityCoefficient(const FluidState& fs, static LhsEval computeFugacityCoefficient(const FluidState& fs2,
const Params& params, const Params& params2,
unsigned phaseIdx, unsigned phaseIdx,
unsigned compIdx) unsigned compIdx)
{ {
#warning also a HACK, should investigate why
auto fs = fs2;
double sumx = 0.0;
for (int i = 0; i < FluidState::numComponents; ++i)
sumx += Opm::scalarValue(fs.moleFraction(phaseIdx, i));
if (sumx < 0.95) {
double alpha = 0.95/sumx;
std::cerr << "normalize: " << sumx
<< " alpha: " << alpha << "\n";
for (int i = 0; i < FluidState::numComponents; ++i)
fs.setMoleFraction(phaseIdx, i, alpha*fs.moleFraction(phaseIdx, i));
}
auto params = params2;
params.updatePhase(fs, phaseIdx);
// note that we normalize the component mole fractions, so // note that we normalize the component mole fractions, so
// that their sum is 100%. This increases numerical stability // that their sum is 100%. This increases numerical stability
// considerably if the fluid state is not physical. // considerably if the fluid state is not physical.

View File

@ -78,7 +78,18 @@ void testChiFlash()
fs.setDensity(FluidSystem::gasPhaseIdx, FluidSystem::density(fs, paramCache, FluidSystem::gasPhaseIdx)); fs.setDensity(FluidSystem::gasPhaseIdx, FluidSystem::density(fs, paramCache, FluidSystem::gasPhaseIdx));
ComponentVector zInit(0.); // TODO; zInit needs to be normalized. ComponentVector zInit(0.); // TODO; zInit needs to be normalized.
const double flash_tolerance = 1.e-8; {
Scalar sumMoles = 0.0;
for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
Scalar tmp = Opm::getValue(fs.molarity(phaseIdx, compIdx) * fs.saturation(phaseIdx));
zInit[compIdx] += Opm::max(tmp, 1e-8);
sumMoles += tmp;
}
}
zInit /= sumMoles;
}
const double flash_tolerance = -1.; // just to test the setup in co2-compositional
const int flash_verbosity = 1; const int flash_verbosity = 1;
const std::string flash_twophase_method = "ssi"; const std::string flash_twophase_method = "ssi";