OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

support ZMF explicit initialization

Browse Source
This commit is contained in:
Kai Bao 2024-09-23 11:51:53 +02:00
parent 3acf0cbd14
commit 562da9ef70
1 changed files with 69 additions and 28 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

97
opm/simulators/flow/FlowProblemComp.hpp
View File

@ -350,55 +350,66 @@ protected:
const auto& vanguard = simulator.vanguard();
const auto& eclState = vanguard.eclState();
const auto& fp = eclState.fieldProps();
// const bool has_xmf = fp.has_double("XMF");
assert(fp.has_double("XMF"));
// const bool has_ymf = fp.has_double("YMF");
assert(fp.has_double("YMF"));
const bool has_temp = fp.has_double("TEMPI");
const bool has_pressure = fp.has_double("PRESSURE");
if (!has_pressure)
throw std::runtime_error("The ECL input file requires the presence of the PRESSURE "
"keyword if the model is initialized explicitly");
const bool has_xmf = fp.has_double("XMF");
const bool has_ymf = fp.has_double("YMF");
const bool has_zmf = fp.has_double("ZMF");
if (! (has_zmf || (has_xmf && has_ymf))) {
throw std::runtime_error("The ECL input file requires the presence of ZMF or XMF and YMF "
"keyword if the model is initialized explicitly");
}
if (has_zmf && (has_xmf || has_ymf)) {
throw std::runtime_error("The ECL input file can not handle explicit initialization "
"with both ZMF and XMF or YMF");
}
if (has_xmf != has_ymf) {
throw std::runtime_error("The ECL input file needs XMF and YMF combined to do the explicit "
"initializtion when using XMF or YMF");
}
const bool has_temp = fp.has_double("TEMPI");
// const bool has_gas = fp.has_double("SGAS");
assert(fp.has_double("SGAS"));
if (!has_pressure)
throw std::runtime_error("The ECL input file requires the presence of the PRESSURE "
"keyword if the model is initialized explicitly");
std::size_t numDof = this->model().numGridDof();
initialFluidStates_.resize(numDof);
std::vector<double> xmfData;
std::vector<double> ymfData;
std::vector<double> waterSaturationData;
std::vector<double> gasSaturationData;
std::vector<double> soilData;
std::vector<double> pressureData;
std::vector<double> tempiData;
if (waterPhaseIdx > 0 && Indices::numPhases > 1)
const bool water_active = FluidSystem::phaseIsActive(waterPhaseIdx);
const bool gas_active = FluidSystem::phaseIsActive(gasPhaseIdx);
const bool oil_active = FluidSystem::phaseIsActive(oilPhaseIdx);
if (water_active && Indices::numPhases > 1)
waterSaturationData = fp.get_double("SWAT");
else
waterSaturationData.resize(numDof);
pressureData = fp.get_double("PRESSURE");
assert(waterPhaseIdx < 0);
xmfData = fp.get_double("XMF");
ymfData = fp.get_double("YMF");
if (has_temp) {
tempiData = fp.get_double("TEMPI");
} else {
; // TODO: throw?
}
if (gasPhaseIdx > 0) // && FluidSystem::phaseIsActive(oilPhaseIdx))
if (gas_active) // && FluidSystem::phaseIsActive(oilPhaseIdx))
gasSaturationData = fp.get_double("SGAS");
else
gasSaturationData.resize(numDof);
// constexpr std::size_t num_comps = 3;
for (std::size_t dofIdx = 0; dofIdx < numDof; ++dofIdx) {
auto& dofFluidState = initialFluidStates_[dofIdx];
// dofFluidState.setPvtRegionIndex(pvtRegionIndex(dofIdx));
@ -407,11 +418,11 @@ protected:
assert(std::isfinite(temperatureLoc) && temperatureLoc > 0);
dofFluidState.setTemperature(temperatureLoc);
if (gasPhaseIdx > 0) {
dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
gasSaturationData[dofIdx]);
if (gas_active) {
dofFluidState.setSaturation(FluidSystem::gasPhaseIdx,
gasSaturationData[dofIdx]);
}
if (oilPhaseIdx > 0) {
if (oil_active) {
dofFluidState.setSaturation(FluidSystem::oilPhaseIdx,
1.0
- waterSaturationData[dofIdx]
@ -438,13 +449,43 @@ protected:
dofFluidState.setPressure(phaseIdx, pressure);
}
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
const std::size_t data_idx = compIdx * numDof + dofIdx;
const Scalar xmf = xmfData[data_idx];
const Scalar ymf = ymfData[data_idx];
if (has_xmf && has_ymf) {
const auto& xmfData = fp.get_double("XMF");
const auto& ymfData = fp.get_double("YMF");
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
const std::size_t data_idx = compIdx * numDof + dofIdx;
const Scalar xmf = xmfData[data_idx];
const Scalar ymf = ymfData[data_idx];
dofFluidState.setMoleFraction(FluidSystem::oilPhaseIdx, compIdx, xmf);
dofFluidState.setMoleFraction(FluidSystem::gasPhaseIdx, compIdx, ymf);
dofFluidState.setMoleFraction(FluidSystem::oilPhaseIdx, compIdx, xmf);
dofFluidState.setMoleFraction(FluidSystem::gasPhaseIdx, compIdx, ymf);
}
}
if (has_zmf) {
const auto& zmfData = fp.get_double("ZMF");
for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
const std::size_t data_idx = compIdx * numDof + dofIdx;
const Scalar zmf = zmfData[data_idx];
// TODO: I know we do this for co2_ptflash example, but maybe we should not do it this way?
if (gas_active) {
const auto& sat_gas = dofFluidState.saturation(gasPhaseIdx);
if (sat_gas > 0.) {
dofFluidState.setMoleFraction(FluidSystem::gasPhaseIdx, compIdx, zmf);
} else {
dofFluidState.setMoleFraction(FluidSystem::gasPhaseIdx, compIdx, 0.);
}
}
if (oil_active) {
dofFluidState.setMoleFraction(FluidSystem::oilPhaseIdx, compIdx, zmf);
const auto& sat_oil = dofFluidState.saturation(oilPhaseIdx);
if (sat_oil > 0.) {
dofFluidState.setMoleFraction(FluidSystem::oilPhaseIdx, compIdx, zmf);
} else {
dofFluidState.setMoleFraction(FluidSystem::oilPhaseIdx, compIdx, 0.);
}
}
}
}
}
}