From ab3be6dce408148c516b9203b1c84553d181ea78 Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Fri, 8 Apr 2022 22:31:21 +0200 Subject: [PATCH 1/6] adjustments to account for vaporized water --- ebos/equil/initstateequil.hh | 2 +- opm/simulators/wells/MultisegmentWellEval.cpp | 9 ++-- .../wells/MultisegmentWell_impl.hpp | 1 + opm/simulators/wells/RateConverter.hpp | 6 +-- opm/simulators/wells/SingleWellState.hpp | 1 + opm/simulators/wells/StandardWell.hpp | 3 ++ opm/simulators/wells/StandardWell_impl.hpp | 47 +++++++++++++++++-- opm/simulators/wells/WellState.cpp | 1 + 8 files changed, 58 insertions(+), 12 deletions(-) diff --git a/ebos/equil/initstateequil.hh b/ebos/equil/initstateequil.hh index 83cb541df..d2ed3b212 100644 --- a/ebos/equil/initstateequil.hh +++ b/ebos/equil/initstateequil.hh @@ -273,7 +273,7 @@ private: bGas = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(pvtRegionIdx_, temp_, press); } else { - bGas = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx_, temp_, press, rv); + bGas = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx_, temp_, press, rv, 0.0/*=Rvw*/); } double rho = bGas * FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, pvtRegionIdx_); if (FluidSystem::enableVaporizedOil()) { diff --git a/opm/simulators/wells/MultisegmentWellEval.cpp b/opm/simulators/wells/MultisegmentWellEval.cpp index 72469cfbb..0e034b806 100644 --- a/opm/simulators/wells/MultisegmentWellEval.cpp +++ b/opm/simulators/wells/MultisegmentWellEval.cpp @@ -677,6 +677,7 @@ computeSegmentFluidProperties(const EvalWell& temperature, } EvalWell rv(0.0); + EvalWell rvw(0.0); // gas phase if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); @@ -692,9 +693,9 @@ computeSegmentFluidProperties(const EvalWell& temperature, rv = rvmax; } b[gasCompIdx] = - FluidSystem::gasPvt().inverseFormationVolumeFactor(pvt_region_index, temperature, seg_pressure, rv); + FluidSystem::gasPvt().inverseFormationVolumeFactor(pvt_region_index, temperature, seg_pressure, rv, rvw); visc[gasCompIdx] = - FluidSystem::gasPvt().viscosity(pvt_region_index, temperature, seg_pressure, rv); + FluidSystem::gasPvt().viscosity(pvt_region_index, temperature, seg_pressure, rv, rvw); phase_densities[gasCompIdx] = b[gasCompIdx] * surf_dens[gasCompIdx] + rv * b[gasCompIdx] * surf_dens[oilCompIdx]; } else { // no oil exists @@ -1114,6 +1115,7 @@ getSegmentSurfaceVolume(const EvalWell& temperature, } EvalWell rv(0.0); + EvalWell rvw(0.0); // gas phase if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); @@ -1137,7 +1139,8 @@ getSegmentSurfaceVolume(const EvalWell& temperature, FluidSystem::gasPvt().inverseFormationVolumeFactor(pvt_region_index, temperature, seg_pressure, - rv); + rv, + rvw); } else { // no oil exists b[gasCompIdx] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(pvt_region_index, diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index 112f44e24..e7c7bfc06 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -1509,6 +1509,7 @@ namespace Opm auto& ws = well_state.well(this->index_of_well_); ws.dissolved_gas_rate = 0; ws.vaporized_oil_rate = 0; + ws.vaporized_wat_rate = 0; // for the black oil cases, there will be four equations, // the first three of them are the mass balance equations, the last one is the pressure equations. diff --git a/opm/simulators/wells/RateConverter.hpp b/opm/simulators/wells/RateConverter.hpp index 68e3dec66..2f73f4bde 100644 --- a/opm/simulators/wells/RateConverter.hpp +++ b/opm/simulators/wells/RateConverter.hpp @@ -318,7 +318,7 @@ namespace Opm { if (RegionAttributeHelpers::PhaseUsed::gas(pu)) { // q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s) - const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv); + const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, 0.0 /*=Rvw*/); const double den = bg * detR; coeff[ig] += 1.0 / den; @@ -359,7 +359,7 @@ namespace Opm { } if (RegionAttributeHelpers::PhaseUsed::gas(pu)) { - const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0); + const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0, 0.0); coeff[ig] += 1.0 / bg; } } @@ -447,7 +447,7 @@ namespace Opm { if (RegionAttributeHelpers::PhaseUsed::gas(pu)) { // q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s) - const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv); + const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, 0.0 /*=Rvw*/); const double den = bg * detR; voidage_rates[ig] = surface_rates[ig]; diff --git a/opm/simulators/wells/SingleWellState.hpp b/opm/simulators/wells/SingleWellState.hpp index 3fd99aa77..167efcf09 100644 --- a/opm/simulators/wells/SingleWellState.hpp +++ b/opm/simulators/wells/SingleWellState.hpp @@ -56,6 +56,7 @@ public: double temperature{0}; double dissolved_gas_rate{0}; double vaporized_oil_rate{0}; + double vaporized_wat_rate{0}; std::vector well_potentials; std::vector productivity_index; std::vector surface_rates; diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index 9cefbfdca..561acc4b1 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -280,6 +280,7 @@ namespace Opm std::vector& cq_s, double& perf_dis_gas_rate, double& perf_vap_oil_rate, + double& perf_vap_wat_rate, DeferredLogger& deferred_logger) const; void computePerfRateScalar(const IntensiveQuantities& intQuants, @@ -297,6 +298,7 @@ namespace Opm const Value& bhp, const Value& rs, const Value& rv, + const Value& rvw, std::vector& b_perfcells_dense, const double Tw, const int perf, @@ -306,6 +308,7 @@ namespace Opm std::vector& cq_s, double& perf_dis_gas_rate, double& perf_vap_oil_rate, + double& perf_vap_wat_rate, DeferredLogger& deferred_logger) const; void computeWellRatesWithBhpIterations(const Simulator& ebosSimulator, diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index 362d13813..f317a0e63 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -94,12 +94,15 @@ namespace Opm std::vector& cq_s, double& perf_dis_gas_rate, double& perf_vap_oil_rate, + double& perf_vap_wat_rate, DeferredLogger& deferred_logger) const { const auto& fs = intQuants.fluidState(); const EvalWell pressure = this->extendEval(this->getPerfCellPressure(fs)); const EvalWell rs = this->extendEval(fs.Rs()); const EvalWell rv = this->extendEval(fs.Rv()); + const EvalWell rvw = this->extendEval(fs.Rvw()); + std::vector b_perfcells_dense(this->num_components_, EvalWell{this->numWellEq_ + Indices::numEq, 0.0}); for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) { if (!FluidSystem::phaseIsActive(phaseIdx)) { @@ -139,6 +142,7 @@ namespace Opm bhp, rs, rv, + rvw, b_perfcells_dense, Tw, perf, @@ -148,6 +152,7 @@ namespace Opm cq_s, perf_dis_gas_rate, perf_vap_oil_rate, + perf_vap_wat_rate, deferred_logger); } @@ -167,6 +172,7 @@ namespace Opm const Scalar pressure = this->getPerfCellPressure(fs).value(); const Scalar rs = fs.Rs().value(); const Scalar rv = fs.Rv().value(); + const Scalar rvw = fs.Rvw().value(); std::vector b_perfcells_dense(this->num_components_, 0.0); for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) { if (!FluidSystem::phaseIsActive(phaseIdx)) { @@ -197,6 +203,7 @@ namespace Opm Scalar perf_dis_gas_rate = 0.0; Scalar perf_vap_oil_rate = 0.0; + Scalar perf_vap_wat_rate = 0.0; // surface volume fraction of fluids within wellbore std::vector cmix_s(this->numComponents(), 0.0); @@ -209,6 +216,7 @@ namespace Opm bhp, rs, rv, + rvw, b_perfcells_dense, Tw, perf, @@ -218,6 +226,7 @@ namespace Opm cq_s, perf_dis_gas_rate, perf_vap_oil_rate, + perf_vap_wat_rate, deferred_logger); } @@ -230,6 +239,7 @@ namespace Opm const Value& bhp, const Value& rs, const Value& rv, + const Value& rvw, std::vector& b_perfcells_dense, const double Tw, const int perf, @@ -239,6 +249,7 @@ namespace Opm std::vector& cq_s, double& perf_dis_gas_rate, double& perf_vap_oil_rate, + double& perf_vap_wat_rate, DeferredLogger& deferred_logger) const { // Pressure drawdown (also used to determine direction of flow) @@ -277,6 +288,12 @@ namespace Opm perf_dis_gas_rate = getValue(dis_gas); perf_vap_oil_rate = getValue(vap_oil); } + if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + const Value vap_wat = rvw * cq_sGas; + cq_s[waterCompIdx] += vap_wat; + if (this->isProducer()) + perf_vap_wat_rate = getValue(vap_wat); + } } } else { @@ -353,10 +370,13 @@ namespace Opm // 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 * g_g + rs * q_or * b_o + // q_gs = q_gr * b_g + rs * q_or * b_o + // q_ws = q_wr * b_w + rvw * q_gr * b_g // 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) + // q_wr = 1/b_w * (rvw * q_gr * b_g -q_ws)= 1/b_w * (rvw * 1/d * (q_gs - rs * q_os) - q_ws) = + // 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) const double d = 1.0 - getValue(rv) * getValue(rs); @@ -378,6 +398,18 @@ namespace Opm perf_dis_gas_rate = getValue(rs) * (getValue(cq_s[oilCompIdx]) - getValue(rv) * getValue(cq_s[gasCompIdx])) / d; } } + if (/*has_watVapor && */ FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { + const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + // q_ws = q_wr * b_w + rvw * q_gr * b_g + // d = 1.0 - rs * rv + // q_wr = 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) + + const double d = 1.0 - getValue(rv) * getValue(rs); + // vaporized water in gas + // rvw * q_gr * b_g = q_ws -q_wr *b_w = q_ws - (rvw * (q_gs - rs * q_os) - d * q_ws)/d = rvw * (q_gs -rs *q_os) / d + perf_vap_wat_rate = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; + } } } } @@ -426,6 +458,7 @@ namespace Opm ws.vaporized_oil_rate = 0; ws.dissolved_gas_rate = 0; + ws.vaporized_wat_rate = 0; const int np = this->number_of_phases_; @@ -490,6 +523,7 @@ namespace Opm const auto& comm = this->parallel_well_info_.communication(); ws.dissolved_gas_rate = comm.sum(ws.dissolved_gas_rate); ws.vaporized_oil_rate = comm.sum(ws.vaporized_oil_rate); + ws.vaporized_wat_rate = comm.sum(ws.vaporized_wat_rate); } // accumulate resWell_ and invDuneD_ in parallel to get effects of all perforations (might be distributed) @@ -548,10 +582,11 @@ namespace Opm double perf_dis_gas_rate = 0.; double perf_vap_oil_rate = 0.; + double perf_vap_wat_rate = 0.; double trans_mult = ebosSimulator.problem().template rockCompTransMultiplier(intQuants, cell_idx); const double Tw = this->well_index_[perf] * trans_mult; computePerfRateEval(intQuants, mob, bhp, Tw, perf, allow_cf, - cq_s, perf_dis_gas_rate, perf_vap_oil_rate, deferred_logger); + cq_s, perf_dis_gas_rate, perf_vap_oil_rate, perf_vap_wat_rate, deferred_logger); auto& ws = well_state.well(this->index_of_well_); auto& perf_data = ws.perf_data; @@ -571,6 +606,7 @@ namespace Opm if (this->isProducer()) { ws.dissolved_gas_rate += perf_dis_gas_rate; ws.vaporized_oil_rate += perf_vap_oil_rate; + ws.vaporized_wat_rate += perf_vap_wat_rate; } if constexpr (has_energy) { @@ -696,7 +732,7 @@ namespace Opm if constexpr (has_brine) { // TODO: the application of well efficiency factor has not been tested with an example yet const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); - EvalWell cq_s_sm = cq_s[waterCompIdx]; + EvalWell cq_s_sm = cq_s[waterCompIdx] - perf_vap_wat_rate; if (this->isInjector()) { cq_s_sm *= this->wsalt(); } else { @@ -1303,7 +1339,7 @@ namespace Opm } rv = std::min(rv, rvmax_perf[perf]); - b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg, rv); + b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg, rv, 0.0 /*Rvw*/); } else { b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg); @@ -2017,10 +2053,11 @@ namespace Opm std::vector cq_s(this->num_components_, {this->numWellEq_ + Indices::numEq, 0.}); double perf_dis_gas_rate = 0.; double perf_vap_oil_rate = 0.; + double perf_vap_wat_rate = 0.; double trans_mult = ebos_simulator.problem().template rockCompTransMultiplier(int_quant, cell_idx); const double Tw = this->well_index_[perf] * trans_mult; computePerfRateEval(int_quant, mob, bhp, Tw, perf, allow_cf, - cq_s, perf_dis_gas_rate, perf_vap_oil_rate, deferred_logger); + cq_s, perf_dis_gas_rate, perf_vap_oil_rate, perf_vap_wat_rate, deferred_logger); // TODO: make area a member const double area = 2 * M_PI * this->perf_rep_radius_[perf] * this->perf_length_[perf]; const auto& material_law_manager = ebos_simulator.problem().materialLawManager(); diff --git a/opm/simulators/wells/WellState.cpp b/opm/simulators/wells/WellState.cpp index 16d632e43..919cd5e8e 100644 --- a/opm/simulators/wells/WellState.cpp +++ b/opm/simulators/wells/WellState.cpp @@ -407,6 +407,7 @@ WellState::report(const int* globalCellIdxMap, well.rates.set(rt::dissolved_gas, ws.dissolved_gas_rate); well.rates.set(rt::vaporized_oil, ws.vaporized_oil_rate); + well.rates.set(rt::vaporized_water, ws.vaporized_wat_rate); { auto& curr = well.current_control; From aaf063ed98258e84052d410be1342a81e64706f2 Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Fri, 8 Apr 2022 22:34:54 +0200 Subject: [PATCH 2/6] add three-phase simulator allowing salt precipitation and water vaporization --- flow/flow_brine_precsalt_vapwat.cpp | 24 +++++++ flow/flow_ebos_brine_precsalt_vapwat.cpp | 85 ++++++++++++++++++++++++ flow/flow_ebos_brine_precsalt_vapwat.hpp | 42 ++++++++++++ 3 files changed, 151 insertions(+) create mode 100644 flow/flow_brine_precsalt_vapwat.cpp create mode 100644 flow/flow_ebos_brine_precsalt_vapwat.cpp create mode 100644 flow/flow_ebos_brine_precsalt_vapwat.hpp diff --git a/flow/flow_brine_precsalt_vapwat.cpp b/flow/flow_brine_precsalt_vapwat.cpp new file mode 100644 index 000000000..1a93e3cfb --- /dev/null +++ b/flow/flow_brine_precsalt_vapwat.cpp @@ -0,0 +1,24 @@ +/* + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ +#include "config.h" +#include + + +int main(int argc, char** argv) +{ + return Opm::flowEbosBrinePrecsaltVapwatMainStandalone(argc, argv); +} diff --git a/flow/flow_ebos_brine_precsalt_vapwat.cpp b/flow/flow_ebos_brine_precsalt_vapwat.cpp new file mode 100644 index 000000000..93074a3b3 --- /dev/null +++ b/flow/flow_ebos_brine_precsalt_vapwat.cpp @@ -0,0 +1,85 @@ +/* + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ +#include "config.h" + +#include + +#include +#include +#include +#include + +namespace Opm { +namespace Properties { +namespace TTag { +struct EclFlowBrinePrecsaltVapwatProblem { + using InheritsFrom = std::tuple; +}; +} +template +struct EnableBrine { + static constexpr bool value = true; +}; + +template +struct EnableSaltPrecipitation { + static constexpr bool value = true; +}; + +template +struct EnableEvaporation { + static constexpr bool value = true; +}; +}} + +namespace Opm { +void flowEbosBrinePrecsaltVapwatSetDeck(double setupTime, std::shared_ptr deck, + std::shared_ptr eclState, + std::shared_ptr schedule, + std::shared_ptr summaryConfig) +{ + using TypeTag = Properties::TTag::EclFlowBrinePrecsaltVapwatProblem; + using Vanguard = GetPropType; + + Vanguard::setExternalSetupTime(setupTime); + Vanguard::setExternalDeck(std::move(deck)); + Vanguard::setExternalEclState(std::move(eclState)); + Vanguard::setExternalSchedule(std::move(schedule)); + Vanguard::setExternalSummaryConfig(std::move(summaryConfig)); +} + + +// ----------------- Main program ----------------- +int flowEbosBrinePrecsaltVapwatMain(int argc, char** argv, bool outputCout, bool outputFiles) +{ + // we always want to use the default locale, and thus spare us the trouble + // with incorrect locale settings. + resetLocale(); + + FlowMainEbos + mainfunc {argc, argv, outputCout, outputFiles}; + return mainfunc.execute(); +} + +int flowEbosBrinePrecsaltVapwatMainStandalone(int argc, char** argv) +{ + using TypeTag = Properties::TTag::EclFlowBrinePrecsaltVapwatProblem; + auto mainObject = Opm::Main(argc, argv); + return mainObject.runStatic(); +} + +} diff --git a/flow/flow_ebos_brine_precsalt_vapwat.hpp b/flow/flow_ebos_brine_precsalt_vapwat.hpp new file mode 100644 index 000000000..ffd316525 --- /dev/null +++ b/flow/flow_ebos_brine_precsalt_vapwat.hpp @@ -0,0 +1,42 @@ +/* + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ +#ifndef FLOW_EBOS_BRINE_PRECSALT_VAPWAT_HPP +#define FLOW_EBOS_BRINE_PRECSALT_VAPWAT_HPP + +#include + +namespace Opm { + +class Deck; +class EclipseState; +class Schedule; +class SummaryConfig; + +void flowEbosBrinePrecsaltVapwatSetDeck(double setupTime, std::shared_ptr deck, + std::shared_ptr eclState, + std::shared_ptr schedule, + std::shared_ptr summaryConfig); + +//! \brief Main function used in flow binary. +int flowEbosBrinePrecsaltVapwatMain(int argc, char** argv, bool outputCout, bool outputFiles); + +//! \brief Main function used in flow_brine binary. +int flowEbosBrinePrecsaltVapwatMainStandalone(int argc, char** argv); + +} + +#endif // FLOW_EBOS_BRINE_PRECSALT_VAPWAT_HPP From c2fa5fc5a8d29315fff913f33b7bf8a8b20d8dc0 Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Mon, 11 Apr 2022 21:56:44 +0200 Subject: [PATCH 3/6] adding three-phase simulator including water evaporation and saltpreciopitation --- CMakeLists.txt | 2 +- opm/simulators/flow/Main.hpp | 15 ++++++++++++--- opm/simulators/wells/StandardWell_impl.hpp | 1 + 3 files changed, 14 insertions(+), 4 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index ada75f293..77407f322 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -388,7 +388,7 @@ add_dependencies(moduleVersion opmsimulators) set(FLOW_MODELS blackoil brine energy extbo foam gasoil gaswater oilwater oilwater_brine gaswater_brine oilwater_polymer oilwater_polymer_injectivity micp polymer solvent - gasoil_energy brine_saltprecipitation gaswater_saltprec_vapwat) + gasoil_energy brine_saltprecipitation gaswater_saltprec_vapwat brine_precsalt_vapwat) set(FLOW_VARIANT_MODELS brine_energy onephase onephase_energy) set(FLOW_TGTS) diff --git a/opm/simulators/flow/Main.hpp b/opm/simulators/flow/Main.hpp index 34ebc1f8c..2762a6bce 100644 --- a/opm/simulators/flow/Main.hpp +++ b/opm/simulators/flow/Main.hpp @@ -35,6 +35,7 @@ #include #include #include +#include #include #include #include @@ -685,9 +686,17 @@ private: } } else if (eclipseState_->getSimulationConfig().hasPRECSALT()) { - flowEbosBrineSaltPrecipitationSetDeck( - setupTime_, deck_, eclipseState_, schedule_, summaryConfig_); - return flowEbosBrineSaltPrecipitationMain(argc_, argv_, outputCout_, outputFiles_); + if (eclipseState_->getSimulationConfig().hasVAPWAT()) { + //case with water vaporization into gas phase and salt precipitation + flowEbosBrinePrecsaltVapwatSetDeck( + setupTime_, deck_, eclipseState_, schedule_, summaryConfig_); + return flowEbosBrinePrecsaltVapwatMain(argc_, argv_, outputCout_, outputFiles_); + } + else { + flowEbosBrineSaltPrecipitationSetDeck( + setupTime_, deck_, eclipseState_, schedule_, summaryConfig_); + return flowEbosBrineSaltPrecipitationMain(argc_, argv_, outputCout_, outputFiles_); + } } else { flowEbosBrineSetDeck( diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index f317a0e63..d579e9d54 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -289,6 +289,7 @@ namespace Opm perf_vap_oil_rate = 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()) From 17f236540662248dfed1430056c6b3ff8d3d9587 Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Mon, 11 Apr 2022 22:15:35 +0200 Subject: [PATCH 4/6] clean up comments --- opm/simulators/wells/StandardWell_impl.hpp | 8 ++------ 1 file changed, 2 insertions(+), 6 deletions(-) diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index d579e9d54..2c9e32349 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -372,12 +372,9 @@ namespace Opm // 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 - // q_ws = q_wr * b_w + rvw * q_gr * b_g // 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) - // q_wr = 1/b_w * (rvw * q_gr * b_g -q_ws)= 1/b_w * (rvw * 1/d * (q_gs - rs * q_os) - q_ws) = - // 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) const double d = 1.0 - getValue(rv) * getValue(rs); @@ -403,9 +400,7 @@ namespace Opm const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); // q_ws = q_wr * b_w + rvw * q_gr * b_g - // d = 1.0 - rs * rv - // q_wr = 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) - + // q_wr = 1/b_w * (rvw * 1/d * (q_gs - rs * q_os) - q_ws) = 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) const double d = 1.0 - getValue(rv) * getValue(rs); // vaporized water in gas // rvw * q_gr * b_g = q_ws -q_wr *b_w = q_ws - (rvw * (q_gs - rs * q_os) - d * q_ws)/d = rvw * (q_gs -rs *q_os) / d @@ -733,6 +728,7 @@ namespace Opm if constexpr (has_brine) { // TODO: the application of well efficiency factor has not been tested with an example yet const unsigned waterCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + // Correction salt rate; evaporated water does not contain salt EvalWell cq_s_sm = cq_s[waterCompIdx] - perf_vap_wat_rate; if (this->isInjector()) { cq_s_sm *= this->wsalt(); From 5a23084846575d65885e30d9e3444a2227cb8275 Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Tue, 3 May 2022 11:54:41 +0200 Subject: [PATCH 5/6] restored MultisegmentWellEval.cpp --- opm/simulators/wells/StandardWell.hpp | 2 + opm/simulators/wells/StandardWellEval.cpp | 33 ++++++++- opm/simulators/wells/StandardWellEval.hpp | 1 + opm/simulators/wells/StandardWell_impl.hpp | 78 ++++++++++++++++++---- 4 files changed, 101 insertions(+), 13 deletions(-) diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index 561acc4b1..d3d68db95 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -257,6 +257,7 @@ namespace Opm std::vector& b_perf, std::vector& rsmax_perf, std::vector& rvmax_perf, + std::vector& rvwmax_perf, std::vector& surf_dens_perf) const; void computeWellConnectionDensitesPressures(const Simulator& ebosSimulator, @@ -264,6 +265,7 @@ namespace Opm const std::vector& b_perf, const std::vector& rsmax_perf, const std::vector& rvmax_perf, + const std::vector& rvwmax_perf, const std::vector& surf_dens_perf, DeferredLogger& deferred_logger); diff --git a/opm/simulators/wells/StandardWellEval.cpp b/opm/simulators/wells/StandardWellEval.cpp index c57c6b83f..912b97c45 100644 --- a/opm/simulators/wells/StandardWellEval.cpp +++ b/opm/simulators/wells/StandardWellEval.cpp @@ -831,6 +831,7 @@ computeConnectionDensities(const std::vector& perfComponentRates, const std::vector& b_perf, const std::vector& rsmax_perf, const std::vector& rvmax_perf, + const std::vector& rvwmax_perf, const std::vector& surf_dens_perf, DeferredLogger& deferred_logger) { @@ -939,7 +940,7 @@ computeConnectionDensities(const std::vector& perfComponentRates, // Compute volume ratio. x = mix; - // Subtract dissolved gas from oil phase and vapporized oil from gas phase + // Subtract dissolved gas from oil phase and vapporized oil from gas phase and vaporized water from gas phase if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) { const unsigned gaspos = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); const unsigned oilpos = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); @@ -971,7 +972,37 @@ computeConnectionDensities(const std::vector& perfComponentRates, x[oilpos] = (mix[oilpos] - mix[gaspos]*rv)/d; } } + if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + //matrix system: (mix[oilpos] = q_os, x[oilpos] = bo*q_or, etc...) + //┌ ┐ ┌ ┐ ┌ ┐ + //│mix[oilpos] │ | 1 Rv 0 | |x[oilpos] | + //│mix[gaspos] │ = │ Rs 1 0 │ │x[gaspos] │ + //│mix[waterpos]│ │ 0 Rvw 1 │ │x[waterpos │ + //└ ┘ └ ┘ └ ┘ + const unsigned waterpos = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + double rvw = 0.0; + if (!rvwmax_perf.empty() && mix[gaspos] > 1e-12) { + rvw = std::min(mix[waterpos]/mix[gaspos], rvwmax_perf[perf]); + } + if (rvw > 0.0) { + // Subtract water in gas from water mixture + x[waterpos] = mix[waterpos] - x[gaspos] * rvw; + } + } + } else if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + //no oil + const unsigned gaspos = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); + const unsigned waterpos = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); + double rvw = 0.0; + if (!rvwmax_perf.empty() && mix[gaspos] > 1e-12) { + rvw = std::min(mix[waterpos]/mix[gaspos], rvwmax_perf[perf]); + } + if (rvw > 0.0) { + // Subtract water in gas from water mixture + x[waterpos] = mix[waterpos] - mix[gaspos] * rvw; + } } + double volrat = 0.0; for (int component = 0; component < num_comp; ++component) { volrat += x[component] / b_perf[perf*num_comp+ component]; diff --git a/opm/simulators/wells/StandardWellEval.hpp b/opm/simulators/wells/StandardWellEval.hpp index 27fdeeb43..a460c9b01 100644 --- a/opm/simulators/wells/StandardWellEval.hpp +++ b/opm/simulators/wells/StandardWellEval.hpp @@ -144,6 +144,7 @@ protected: const std::vector& b_perf, const std::vector& rsmax_perf, const std::vector& rvmax_perf, + const std::vector& rvwmax_perf, const std::vector& surf_dens_perf, DeferredLogger& deferred_logger); diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index 2c9e32349..2fa282e54 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -395,16 +395,18 @@ namespace Opm // rs * q_or * b_o = rs * (q_os - rv * q_gs) / d perf_dis_gas_rate = 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_vap_wat_rate = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; + } } - if (/*has_watVapor && */ FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); + else if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { + //no oil const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); - // q_ws = q_wr * b_w + rvw * q_gr * b_g - // q_wr = 1/b_w * (rvw * 1/d * (q_gs - rs * q_os) - q_ws) = 1 / (b_w * d) * (rvw * (q_gs - rs * q_os) - d * q_ws) - const double d = 1.0 - getValue(rv) * getValue(rs); - // vaporized water in gas - // rvw * q_gr * b_g = q_ws -q_wr *b_w = q_ws - (rvw * (q_gs - rs * q_os) - d * q_ws)/d = rvw * (q_gs -rs *q_os) / d - perf_vap_wat_rate = getValue(rvw) * (getValue(cq_s[gasCompIdx]) - getValue(rs) * getValue(cq_s[oilCompIdx])) / d; + perf_vap_wat_rate = getValue(rvw) * getValue(cq_s[gasCompIdx]); } } } @@ -1282,6 +1284,7 @@ namespace Opm std::vector& b_perf, std::vector& rsmax_perf, std::vector& rvmax_perf, + std::vector& rvwmax_perf, std::vector& surf_dens_perf) const { const int nperf = this->number_of_perforations_; @@ -1299,6 +1302,10 @@ namespace Opm rsmax_perf.resize(nperf); rvmax_perf.resize(nperf); } + //rvw is only used if both water and gas is present + if (waterPresent && gasPresent) { + rvwmax_perf.resize(nperf); + } // Compute the average pressure in each well block const auto& perf_press = ws.perf_data.pressure; @@ -1325,7 +1332,35 @@ namespace Opm const unsigned gasCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); const int gaspos = gasCompIdx + perf * this->num_components_; - if (oilPresent) { + if (oilPresent && waterPresent) { + const double oilrate = std::abs(ws.surface_rates[pu.phase_pos[Oil]]); //in order to handle negative rates in producers + const double waterrate = std::abs(ws.surface_rates[pu.phase_pos[Water]]); //in order to handle negative rates in producers + rvmax_perf[perf] = FluidSystem::gasPvt().saturatedOilVaporizationFactor(fs.pvtRegionIndex(), temperature, p_avg); + rvwmax_perf[perf] = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(fs.pvtRegionIndex(), temperature, p_avg); + double rv = 0.0; + double rvw = 0.0; + if (oilrate > 0) { + const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? ws.sum_solvent_rates() : 0.0); + if (gasrate > 0) { + rv = oilrate / gasrate; + } + rv = std::min(rv, rvmax_perf[perf]); + } + if (waterrate > 0) { + const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? ws.sum_solvent_rates() : 0.0); + if (gasrate > 0) { + rvw = waterrate / gasrate; + } + rvw = std::min(rvw, rvwmax_perf[perf]); + } + if (rv > 0.0 || rvw > 0.0){ + b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg, rv, rvw); + } + else { + b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg); + } + } else if (oilPresent) { + //no water const double oilrate = std::abs(ws.surface_rates[pu.phase_pos[Oil]]); //in order to handle negative rates in producers rvmax_perf[perf] = FluidSystem::gasPvt().saturatedOilVaporizationFactor(fs.pvtRegionIndex(), temperature, p_avg); if (oilrate > 0) { @@ -1341,6 +1376,23 @@ namespace Opm else { b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg); } + } else if (waterPresent) { + //no oil + const double waterrate = std::abs(ws.surface_rates[pu.phase_pos[Water]]); //in order to handle negative rates in producers + rvwmax_perf[perf] = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(fs.pvtRegionIndex(), temperature, p_avg); + if (waterrate > 0) { + const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (has_solvent ? ws.sum_solvent_rates() : 0.0); + double rvw = 0.0; + if (gasrate > 0) { + rvw = waterrate / gasrate; + } + rvw = std::min(rvw, rvwmax_perf[perf]); + + b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg, 0.0 /*Rv*/, rvw); + } + else { + b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg); + } } else { b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(fs.pvtRegionIndex(), temperature, p_avg); @@ -1505,6 +1557,7 @@ namespace Opm const std::vector& b_perf, const std::vector& rsmax_perf, const std::vector& rvmax_perf, + const std::vector& rvwmax_perf, const std::vector& surf_dens_perf, DeferredLogger& deferred_logger) { @@ -1561,7 +1614,7 @@ namespace Opm } } - this->computeConnectionDensities(perfRates, b_perf, rsmax_perf, rvmax_perf, surf_dens_perf, deferred_logger); + this->computeConnectionDensities(perfRates, b_perf, rsmax_perf, rvmax_perf, rvwmax_perf, surf_dens_perf, deferred_logger); this->computeConnectionPressureDelta(); } @@ -1583,9 +1636,10 @@ namespace Opm std::vector b_perf; std::vector rsmax_perf; std::vector rvmax_perf; + std::vector rvwmax_perf; std::vector surf_dens_perf; - computePropertiesForWellConnectionPressures(ebosSimulator, well_state, b_perf, rsmax_perf, rvmax_perf, surf_dens_perf); - computeWellConnectionDensitesPressures(ebosSimulator, well_state, b_perf, rsmax_perf, rvmax_perf, surf_dens_perf, deferred_logger); + computePropertiesForWellConnectionPressures(ebosSimulator, well_state, b_perf, rsmax_perf, rvmax_perf, rvwmax_perf, surf_dens_perf); + computeWellConnectionDensitesPressures(ebosSimulator, well_state, b_perf, rsmax_perf, rvmax_perf, rvwmax_perf, surf_dens_perf, deferred_logger); } From 6fe10a3e234b180daeec1ec5968b41a0cdeb054f Mon Sep 17 00:00:00 2001 From: Paul Egberts Date: Wed, 11 May 2022 11:50:54 +0200 Subject: [PATCH 6/6] add throw for case msw is true and vapwat is true --- opm/simulators/wells/MultisegmentWell_impl.hpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index e7c7bfc06..e75fba6c5 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -71,6 +71,10 @@ namespace Opm if constexpr (Base::has_brine) { OPM_THROW(std::runtime_error, "brine is not supported by multisegment well yet"); } + + if constexpr (Base::has_watVapor) { + OPM_THROW(std::runtime_error, "water evaporation is not supported by multisegment well yet"); + } }