diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index ec0f576ff..01225feab 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -762,8 +762,8 @@ namespace Opm const Value drawdown = cell_press_at_perf - perf_press; // producing perforations - if ( drawdown > 0.0) { - // Do nothing is crossflow is not allowed + if (drawdown > 0.0) { + // Do nothing if crossflow is not allowed if (!allow_cf && this->isInjector()) { return; } diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index 043d8dca7..d6afc29b0 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -466,6 +466,54 @@ namespace Opm const double dis_gas_rate, const std::vector& cq_s, const IntensiveQuantities& intQuants) const; + + template + void gasOilPerfRateInj(const std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rv, + const Value& rs, + const Value& pressure, + const Value& rvw, + DeferredLogger& deferred_logger) const; + + template + void gasOilPerfRateProd(std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rv, + const Value& rs, + const Value& rvw) const; + + template + void gasWaterPerfRateProd(std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rvw, + const Value& rsw) const; + + template + void gasWaterPerfRateInj(const std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rvw, + const Value& rsw, + const Value& pressure, + DeferredLogger& deferred_logger) const; + + template + void disOilVapWatVolumeRatio(Value& volumeRatio, + const Value& rvw, + const Value& rsw, + const Value& pressure, + const std::vector& cmix_s, + const std::vector& b_perfcells_dense, + DeferredLogger& deferred_logger) const; + + template + void gasOilVolumeRatio(Value& volumeRatio, + const Value& rv, + const Value& rs, + const Value& pressure, + const std::vector& cmix_s, + const std::vector& b_perfcells_dense, + DeferredLogger& deferred_logger) const; }; } diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index 096caf61f..80172708a 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -37,6 +37,26 @@ #include #include + +namespace { + +template +auto dValueError(const dValue& d, + const std::string& name, + const std::string& methodName, + const Value& Rs, + const Value& Rv, + const Value& pressure) +{ + return fmt::format("Problematic d value {} obtained for well {}" + " during {} calculations with rs {}" + ", rv {} and pressure {}." + " Continue as if no dissolution (rs = 0) and vaporization (rv = 0) " + " for this connection.", d, name, methodName, Rs, Rv, pressure); +} + +} + namespace Opm { @@ -219,8 +239,8 @@ namespace Opm } // producing perforations - if ( drawdown > 0 ) { - //Do nothing if crossflow is not allowed + if (drawdown > 0) { + // Do nothing if crossflow is not allowed if (!allow_cf && this->isInjector()) { return; } @@ -232,46 +252,12 @@ namespace Opm } if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - const Value cq_sOil = cq_s[oilCompIdx]; - const Value cq_sGas = cq_s[gasCompIdx]; - const Value dis_gas = rs * cq_sOil; - const Value vap_oil = rv * cq_sGas; - - cq_s[gasCompIdx] += dis_gas; - cq_s[oilCompIdx] += vap_oil; - - // recording the perforation solution gas rate and solution oil rates - if (this->isProducer()) { - perf_rates.dis_gas = getValue(dis_gas); - perf_rates.vap_oil = getValue(vap_oil); - } - - if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { - const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); - const Value vap_wat = rvw * cq_sGas; - cq_s[waterCompIdx] += vap_wat; - if (this->isProducer()) - perf_rates.vap_wat = getValue(vap_wat); - } + gasOilPerfRateProd(cq_s, perf_rates, rv, rs, rvw); } else if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - const Value cq_sWat = cq_s[waterCompIdx]; - const Value cq_sGas = cq_s[gasCompIdx]; - const Value vap_wat = rvw * cq_sGas; - const Value dis_gas_wat = rsw * cq_sWat; - cq_s[waterCompIdx] += vap_wat; - cq_s[gasCompIdx] += dis_gas_wat; - if (this->isProducer()) { - perf_rates.vap_wat = getValue(vap_wat); - perf_rates.dis_gas_in_water = getValue(dis_gas_wat); - } + gasWaterPerfRateProd(cq_s, perf_rates, rvw, rsw); } - } else { - //Do nothing if crossflow is not allowed + // Do nothing if crossflow is not allowed if (!allow_cf && this->isProducer()) { return; } @@ -289,26 +275,8 @@ namespace Opm Value volumeRatio = bhp * 0.0; // initialize it with the correct type ; if (FluidSystem::enableVaporizedWater() && FluidSystem::enableDissolvedGasInWater()) { - const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - // Incorporate RSW/RVW factors if both water and gas active - const Value d = 1.0 - rvw * rsw; - - if (d <= 0.0) { - std::ostringstream sstr; - sstr << "Problematic d value " << d << " obtained for well " << this->name() - << " during computePerfRate calculations with rsw " << rsw - << ", rvw " << rvw << " and pressure " << pressure - << " obtaining d " << d - << " Continue as if no dissolution (rsw = 0) and vaporization (rvw = 0) " - << " for this connection."; - deferred_logger.debug(sstr.str()); - } - const Value tmp_wat = d > 0.0? (cmix_s[waterCompIdx] - rvw * cmix_s[gasCompIdx]) / d : cmix_s[waterCompIdx]; - volumeRatio += tmp_wat / b_perfcells_dense[waterCompIdx]; - - const Value tmp_gas = d > 0.0? (cmix_s[gasCompIdx] - rsw * cmix_s[waterCompIdx]) / d : cmix_s[gasCompIdx]; - volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx]; + disOilVapWatVolumeRatio(volumeRatio, rvw, rsw, pressure, + cmix_s, b_perfcells_dense, deferred_logger); } else { if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); @@ -321,26 +289,8 @@ namespace Opm } if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - // Incorporate RS/RV factors if both oil and gas active - const Value d = 1.0 - rv * rs; - - if (d <= 0.0) { - std::ostringstream sstr; - sstr << "Problematic d value " << d << " obtained for well " << this->name() - << " during computePerfRate calculations with rs " << rs - << ", rv " << rv << " and pressure " << pressure - << " obtaining d " << d - << " Continue as if no dissolution (rs = 0) and vaporization (rv = 0) " - << " for this connection."; - deferred_logger.debug(sstr.str()); - } - const Value tmp_oil = d > 0.0? (cmix_s[oilCompIdx] - rv * cmix_s[gasCompIdx]) / d : cmix_s[oilCompIdx]; - volumeRatio += tmp_oil / b_perfcells_dense[oilCompIdx]; - - const Value tmp_gas = d > 0.0? (cmix_s[gasCompIdx] - rs * cmix_s[oilCompIdx]) / d : cmix_s[gasCompIdx]; - volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx]; + gasOilVolumeRatio(volumeRatio, rv, rs, pressure, + cmix_s, b_perfcells_dense, deferred_logger); } else { if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) { @@ -354,7 +304,7 @@ namespace Opm } // injecting connections total volumerates at standard conditions - Value cqt_is = cqt_i/volumeRatio; + Value cqt_is = cqt_i / volumeRatio; for (int componentIdx = 0; componentIdx < this->numComponents(); ++componentIdx) { cq_s[componentIdx] = cmix_s[componentIdx] * cqt_is; } @@ -362,68 +312,13 @@ namespace Opm // calculating the perforation solution gas rate and solution oil rates if (this->isProducer()) { if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - // TODO: the formulations here remain to be tested with cases with strong crossflow through production wells - // s means standard condition, r means reservoir condition - // q_os = q_or * b_o + rv * q_gr * b_g - // q_gs = q_gr * b_g + rs * q_or * b_o - // d = 1.0 - rs * rv - // q_or = 1 / (b_o * d) * (q_os - rv * q_gs) - // q_gr = 1 / (b_g * d) * (q_gs - rs * q_os) - - const double d = 1.0 - getValue(rv) * getValue(rs); - - if (d <= 0.0) { - std::ostringstream sstr; - sstr << "Problematic d value " << d << " obtained for well " << this->name() - << " during computePerfRate calculations with rs " << rs - << ", rv " << rv << " and pressure " << pressure - << " obtaining d " << d - << " Continue as if no dissolution (rs = 0) and vaporization (rv = 0) " - << " for this connection."; - deferred_logger.debug(sstr.str()); - } else { - // vaporized oil into gas - // rv * q_gr * b_g = rv * (q_gs - rs * q_os) / d - perf_rates.vap_oil = getValue(rv) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; - // dissolved of gas in oil - // rs * q_or * b_o = rs * (q_os - rv * q_gs) / d - perf_rates.dis_gas = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d; - } - if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { - // q_ws = q_wr * b_w + rvw * q_gr * b_g - // q_wr = 1 / b_w * (q_ws - rvw * q_gr * b_g) = 1 / b_w * (q_ws - rvw * 1 / d (q_gs - rs * q_os)) - // vaporized water in gas - // rvw * q_gr * b_g = q_ws -q_wr *b_w = rvw * (q_gs -rs *q_os) / d - perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; - } + gasOilPerfRateInj(cq_s, perf_rates, + rv, rs, pressure, rvw, deferred_logger); } if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { //no oil - const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); - - const double dw = 1.0 - getValue(rvw) * getValue(rsw); - - if (dw <= 0.0) { - std::ostringstream sstr; - sstr << "Problematic d value " << dw << " obtained for well " << this->name() - << " during computePerfRate calculations with rsw " << rsw - << ", rvw " << rvw << " and pressure " << pressure - << " obtaining d " << dw - << " Continue as if no dissolution (rsw = 0) and vaporization (rvw = 0) " - << " for this connection."; - deferred_logger.debug(sstr.str()); - } else { - // vaporized water into gas - // rvw * q_gr * b_g = rvw * (q_gs - rsw * q_ws) / dw - perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rsw) * getValue(cq_s[waterCompIdx])) / dw; - // dissolved gas in water - // rsw * q_wr * b_w = rsw * (q_ws - rvw * q_gs) / dw - perf_rates.dis_gas_in_water = getValue(rsw) * (getValue(cq_s[waterCompIdx]) - getValue(rvw) * getValue(cq_s[gasCompIdx])) / dw; - } - + gasWaterPerfRateInj(cq_s, perf_rates, rvw, rsw, + pressure, deferred_logger); } } } @@ -2471,4 +2366,206 @@ namespace Opm return {Base::restrictEval(cq_s_zfrac_effective), cq_s_zfrac_effective}; } + + template + template + void + StandardWell:: + gasOilPerfRateInj(const std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rv, + const Value& rs, + const Value& pressure, + const Value& rvw, + DeferredLogger& deferred_logger) const + { + const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + // TODO: the formulations here remain to be tested with cases with strong crossflow through production wells + // s means standard condition, r means reservoir condition + // q_os = q_or * b_o + rv * q_gr * b_g + // q_gs = q_gr * b_g + rs * q_or * b_o + // d = 1.0 - rs * rv + // q_or = 1 / (b_o * d) * (q_os - rv * q_gs) + // q_gr = 1 / (b_g * d) * (q_gs - rs * q_os) + + const double d = 1.0 - getValue(rv) * getValue(rs); + + if (d <= 0.0) { + deferred_logger.debug(dValueError(d, this->name(), + "gasOilPerfRateInj", + rs, rv, pressure)); + } else { + // vaporized oil into gas + // rv * q_gr * b_g = rv * (q_gs - rs * q_os) / d + perf_rates.vap_oil = getValue(rv) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; + // dissolved of gas in oil + // rs * q_or * b_o = rs * (q_os - rv * q_gs) / d + perf_rates.dis_gas = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d; + } + + if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + // q_ws = q_wr * b_w + rvw * q_gr * b_g + // q_wr = 1 / b_w * (q_ws - rvw * q_gr * b_g) = 1 / b_w * (q_ws - rvw * 1 / d (q_gs - rs * q_os)) + // vaporized water in gas + // rvw * q_gr * b_g = q_ws -q_wr *b_w = rvw * (q_gs -rs *q_os) / d + perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; + } + } + + + + template + template + void + StandardWell:: + gasOilPerfRateProd(std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rv, + const Value& rs, + const Value& rvw) const + { + const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + const Value cq_sOil = cq_s[oilCompIdx]; + const Value cq_sGas = cq_s[gasCompIdx]; + const Value dis_gas = rs * cq_sOil; + const Value vap_oil = rv * cq_sGas; + + cq_s[gasCompIdx] += dis_gas; + cq_s[oilCompIdx] += vap_oil; + + // recording the perforation solution gas rate and solution oil rates + if (this->isProducer()) { + perf_rates.dis_gas = getValue(dis_gas); + perf_rates.vap_oil = getValue(vap_oil); + } + + if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + const Value vap_wat = rvw * cq_sGas; + cq_s[waterCompIdx] += vap_wat; + if (this->isProducer()) + perf_rates.vap_wat = getValue(vap_wat); + } + } + + + template + template + void + StandardWell:: + gasWaterPerfRateProd(std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rvw, + const Value& rsw) const + { + const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + const Value cq_sWat = cq_s[waterCompIdx]; + const Value cq_sGas = cq_s[gasCompIdx]; + const Value vap_wat = rvw * cq_sGas; + const Value dis_gas_wat = rsw * cq_sWat; + cq_s[waterCompIdx] += vap_wat; + cq_s[gasCompIdx] += dis_gas_wat; + if (this->isProducer()) { + perf_rates.vap_wat = getValue(vap_wat); + perf_rates.dis_gas_in_water = getValue(dis_gas_wat); + } + } + + + template + template + void + StandardWell:: + gasWaterPerfRateInj(const std::vector& cq_s, + PerforationRates& perf_rates, + const Value& rvw, + const Value& rsw, + const Value& pressure, + DeferredLogger& deferred_logger) const + + { + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + + const double dw = 1.0 - getValue(rvw) * getValue(rsw); + + if (dw <= 0.0) { + deferred_logger.debug(dValueError(dw, this->name(), + "gasWaterPerfRateInj", + rsw, rvw, pressure)); + } else { + // vaporized water into gas + // rvw * q_gr * b_g = rvw * (q_gs - rsw * q_ws) / dw + perf_rates.vap_wat = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rsw) * getValue(cq_s[waterCompIdx])) / dw; + // dissolved gas in water + // rsw * q_wr * b_w = rsw * (q_ws - rvw * q_gs) / dw + perf_rates.dis_gas_in_water = getValue(rsw) * (getValue(cq_s[waterCompIdx]) - getValue(rvw) * getValue(cq_s[gasCompIdx])) / dw; + } + } + + + template + template + void + StandardWell:: + disOilVapWatVolumeRatio(Value& volumeRatio, + const Value& rvw, + const Value& rsw, + const Value& pressure, + const std::vector& cmix_s, + const std::vector& b_perfcells_dense, + DeferredLogger& deferred_logger) const + { + const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + // Incorporate RSW/RVW factors if both water and gas active + const Value d = 1.0 - rvw * rsw; + + if (d <= 0.0) { + deferred_logger.debug(dValueError(d, this->name(), + "disOilVapWatVolumeRatio", + rsw, rvw, pressure)); + } + const Value tmp_wat = d > 0.0 ? (cmix_s[waterCompIdx] - rvw * cmix_s[gasCompIdx]) / d + : cmix_s[waterCompIdx]; + volumeRatio += tmp_wat / b_perfcells_dense[waterCompIdx]; + + const Value tmp_gas = d > 0.0 ? (cmix_s[gasCompIdx] - rsw * cmix_s[waterCompIdx]) / d + : cmix_s[waterCompIdx]; + volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx]; + } + + + template + template + void + StandardWell:: + gasOilVolumeRatio(Value& volumeRatio, + const Value& rv, + const Value& rs, + const Value& pressure, + const std::vector& cmix_s, + const std::vector& b_perfcells_dense, + DeferredLogger& deferred_logger) const + { + const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + // Incorporate RS/RV factors if both oil and gas active + const Value d = 1.0 - rv * rs; + + if (d <= 0.0) { + deferred_logger.debug(dValueError(d, this->name(), + "gasOilVolumeRatio", + rs, rv, pressure)); + } + const Value tmp_oil = d > 0.0? (cmix_s[oilCompIdx] - rv * cmix_s[gasCompIdx]) / d : cmix_s[oilCompIdx]; + volumeRatio += tmp_oil / b_perfcells_dense[oilCompIdx]; + + const Value tmp_gas = d > 0.0? (cmix_s[gasCompIdx] - rs * cmix_s[oilCompIdx]) / d : cmix_s[gasCompIdx]; + volumeRatio += tmp_gas / b_perfcells_dense[gasCompIdx]; + } + } // namespace Opm