diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index fe3c5f3b3..90b3a0d01 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -265,23 +265,14 @@ namespace Opm // calculate the properties for the well connections // to calulate the pressure difference between well connections. + using WellConnectionProps = typename StdWellEval::StdWellConnections::Properties; void computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator, const WellState& well_state, - std::vector& b_perf, - std::vector& rsmax_perf, - std::vector& rvmax_perf, - std::vector& rvwmax_perf, - std::vector& rswmax_perf, - std::vector& surf_dens_perf) const; + WellConnectionProps& props) const; void computeWellConnectionDensitesPressures(const Simulator& ebosSimulator, const WellState& well_state, - const std::vector& b_perf, - const std::vector& rsmax_perf, - const std::vector& rvmax_perf, - const std::vector& rvwmax_perf, - const std::vector& rswmax_perf, - const std::vector& surf_dens_perf, + const WellConnectionProps& props, DeferredLogger& deferred_logger); void computeWellConnectionPressures(const Simulator& ebosSimulator, diff --git a/opm/simulators/wells/StandardWellConnections.cpp b/opm/simulators/wells/StandardWellConnections.cpp index 28e787a14..f9bda316d 100644 --- a/opm/simulators/wells/StandardWellConnections.cpp +++ b/opm/simulators/wells/StandardWellConnections.cpp @@ -89,12 +89,7 @@ computePressureDelta() template void StandardWellConnections:: computeDensities(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& rswmax_perf, - const std::vector& surf_dens_perf, + const Properties& props, DeferredLogger& deferred_logger) { // Verify that we have consistent input. @@ -208,11 +203,11 @@ computeDensities(const std::vector& perfComponentRates, const unsigned oilpos = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); Scalar rs = 0.0; Scalar rv = 0.0; - if (!rsmax_perf.empty() && mix[oilpos] > 1e-12) { - rs = std::min(mix[gaspos]/mix[oilpos], rsmax_perf[perf]); + if (!props.rsmax_perf.empty() && mix[oilpos] > 1e-12) { + rs = std::min(mix[gaspos] / mix[oilpos], props.rsmax_perf[perf]); } - if (!rvmax_perf.empty() && mix[gaspos] > 1e-12) { - rv = std::min(mix[oilpos]/mix[gaspos], rvmax_perf[perf]); + if (!props.rvmax_perf.empty() && mix[gaspos] > 1e-12) { + rv = std::min(mix[oilpos] / mix[gaspos], props.rvmax_perf[perf]); } const Scalar d = 1.0 - rs*rv; if (d <= 0.0) { @@ -246,12 +241,12 @@ computeDensities(const std::vector& perfComponentRates, const unsigned waterpos = Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx); const unsigned gaspos = Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx); Scalar rvw = 0.0; - if (!rvwmax_perf.empty() && mix[gaspos] > 1e-12) { - rvw = std::min(mix[waterpos]/mix[gaspos], rvwmax_perf[perf]); + if (!props.rvwmax_perf.empty() && mix[gaspos] > 1e-12) { + rvw = std::min(mix[waterpos] / mix[gaspos], props.rvwmax_perf[perf]); } Scalar rsw = 0.0; - if (!rswmax_perf.empty() && mix[waterpos] > 1e-12) { - rsw = std::min(mix[gaspos]/mix[waterpos], rswmax_perf[perf]); + if (!props.rswmax_perf.empty() && mix[waterpos] > 1e-12) { + rsw = std::min(mix[gaspos] / mix[waterpos], props.rswmax_perf[perf]); } const Scalar d = 1.0 - rsw*rvw; if (d <= 0.0) { @@ -277,10 +272,10 @@ computeDensities(const std::vector& perfComponentRates, Scalar volrat = 0.0; for (int component = 0; component < num_comp; ++component) { - volrat += x[component] / b_perf[perf*num_comp+ component]; + volrat += x[component] / props.b_perf[perf * num_comp + component]; } for (int component = 0; component < num_comp; ++component) { - surf_dens[component] = surf_dens_perf[perf*num_comp+ component]; + surf_dens[component] = props.surf_dens_perf[perf * num_comp + component]; } // Compute segment density. @@ -296,17 +291,12 @@ computePropertiesForPressures(const WellState& well_state, const std::function& pvtRegionIdx, const std::function& solventInverseFormationVolumeFactor, const std::function& solventRefDensity, - std::vector& b_perf, - std::vector& rsmax_perf, - std::vector& rvmax_perf, - std::vector& rvwmax_perf, - std::vector& rswmax_perf, - std::vector& surf_dens_perf) const + Properties& props) const { const int nperf = well_.numPerfs(); const PhaseUsage& pu = well_.phaseUsage(); - b_perf.resize(nperf * well_.numComponents()); - surf_dens_perf.resize(nperf * well_.numComponents()); + props.b_perf.resize(nperf * well_.numComponents()); + props.surf_dens_perf.resize(nperf * well_.numComponents()); const auto& ws = well_state.well(well_.indexOfWell()); static constexpr int Water = BlackoilPhases::Aqua; @@ -319,13 +309,13 @@ computePropertiesForPressures(const WellState& well_state, //rs and rv are only used if both oil and gas is present if (oilPresent && gasPresent) { - rsmax_perf.resize(nperf); - rvmax_perf.resize(nperf); + props.rsmax_perf.resize(nperf); + props.rvmax_perf.resize(nperf); } //rvw is only used if both water and gas is present if (waterPresent && gasPresent) { - rvwmax_perf.resize(nperf); - rswmax_perf.resize(nperf); + props.rvwmax_perf.resize(nperf); + props.rswmax_perf.resize(nperf); } // Compute the average pressure in each well block @@ -349,16 +339,16 @@ computePropertiesForPressures(const WellState& well_state, // TODO support mutual solubility in water and oil assert(!FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); const double waterrate = std::abs(ws.surface_rates[pu.phase_pos[Water]]); - rswmax_perf[perf] = FluidSystem::waterPvt().saturatedGasDissolutionFactor(region_idx, temperature, p_avg, saltConcentration); + props.rswmax_perf[perf] = FluidSystem::waterPvt().saturatedGasDissolutionFactor(region_idx, temperature, p_avg, saltConcentration); if (waterrate > 0) { const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0); if (gasrate > 0) { rsw = waterrate / gasrate; } - rsw = std::min(rsw, rswmax_perf[perf]); + rsw = std::min(rsw, props.rswmax_perf[perf]); } } - b_perf[ waterCompIdx + perf * well_.numComponents()] = FluidSystem::waterPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rsw, saltConcentration); + props.b_perf[waterCompIdx + perf * well_.numComponents()] = FluidSystem::waterPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rsw, saltConcentration); } if (gasPresent) { @@ -368,8 +358,8 @@ computePropertiesForPressures(const WellState& well_state, 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(region_idx, temperature, p_avg); - rvwmax_perf[perf] = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(region_idx, temperature, p_avg); + props.rvmax_perf[perf] = FluidSystem::gasPvt().saturatedOilVaporizationFactor(region_idx, temperature, p_avg); + props.rvwmax_perf[perf] = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(region_idx, temperature, p_avg); double rv = 0.0; double rvw = 0.0; if (oilrate > 0) { @@ -377,58 +367,58 @@ computePropertiesForPressures(const WellState& well_state, if (gasrate > 0) { rv = oilrate / gasrate; } - rv = std::min(rv, rvmax_perf[perf]); + rv = std::min(rv, props.rvmax_perf[perf]); } if (waterrate > 0) { const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0); if (gasrate > 0) { rvw = waterrate / gasrate; } - rvw = std::min(rvw, rvwmax_perf[perf]); + rvw = std::min(rvw, props.rvwmax_perf[perf]); } if (rv > 0.0 || rvw > 0.0){ - b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rv, rvw); + props.b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rv, rvw); } else { - b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, 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(region_idx, temperature, p_avg); + props.rvmax_perf[perf] = FluidSystem::gasPvt().saturatedOilVaporizationFactor(region_idx, temperature, p_avg); if (oilrate > 0) { const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0); double rv = 0.0; if (gasrate > 0) { rv = oilrate / gasrate; } - rv = std::min(rv, rvmax_perf[perf]); + rv = std::min(rv, props.rvmax_perf[perf]); - b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rv, 0.0 /*Rvw*/); + props.b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rv, 0.0 /*Rvw*/); } else { - b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, 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(region_idx, temperature, p_avg); + props.rvwmax_perf[perf] = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(region_idx, temperature, p_avg); if (waterrate > 0) { const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0); double rvw = 0.0; if (gasrate > 0) { rvw = waterrate / gasrate; } - rvw = std::min(rvw, rvwmax_perf[perf]); + rvw = std::min(rvw, props.rvwmax_perf[perf]); - b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, 0.0 /*Rv*/, rvw); + props.b_perf[gaspos] = FluidSystem::gasPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, 0.0 /*Rv*/, rvw); } else { - b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); } } else { - b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[gaspos] = FluidSystem::gasPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); } } @@ -436,7 +426,7 @@ computePropertiesForPressures(const WellState& well_state, const unsigned oilCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx); const int oilpos = oilCompIdx + perf * well_.numComponents(); if (gasPresent) { - rsmax_perf[perf] = FluidSystem::oilPvt().saturatedGasDissolutionFactor(region_idx, temperature, p_avg); + props.rsmax_perf[perf] = FluidSystem::oilPvt().saturatedGasDissolutionFactor(region_idx, temperature, p_avg); const double gasrate = std::abs(ws.surface_rates[pu.phase_pos[Gas]]) - (Indices::enableSolvent ? ws.sum_solvent_rates() : 0.0); if (gasrate > 0) { const double oilrate = std::abs(ws.surface_rates[pu.phase_pos[Oil]]); @@ -444,13 +434,13 @@ computePropertiesForPressures(const WellState& well_state, if (oilrate > 0) { rs = gasrate / oilrate; } - rs = std::min(rs, rsmax_perf[perf]); - b_perf[oilpos] = FluidSystem::oilPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rs); + rs = std::min(rs, props.rsmax_perf[perf]); + props.b_perf[oilpos] = FluidSystem::oilPvt().inverseFormationVolumeFactor(region_idx, temperature, p_avg, rs); } else { - b_perf[oilpos] = FluidSystem::oilPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[oilpos] = FluidSystem::oilPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); } } else { - b_perf[oilpos] = FluidSystem::oilPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); + props.b_perf[oilpos] = FluidSystem::oilPvt().saturatedInverseFormationVolumeFactor(region_idx, temperature, p_avg); } } @@ -461,13 +451,13 @@ computePropertiesForPressures(const WellState& well_state, } const unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx)); - surf_dens_perf[well_.numComponents() * perf + compIdx] = FluidSystem::referenceDensity( phaseIdx, region_idx ); + props.surf_dens_perf[well_.numComponents() * perf + compIdx] = FluidSystem::referenceDensity( phaseIdx, region_idx ); } // We use cell values for solvent injector if constexpr (Indices::enableSolvent) { - b_perf[well_.numComponents() * perf + Indices::contiSolventEqIdx] = solventInverseFormationVolumeFactor(cell_idx); - surf_dens_perf[well_.numComponents() * perf + Indices::contiSolventEqIdx] = solventRefDensity(cell_idx); + props.b_perf[well_.numComponents() * perf + Indices::contiSolventEqIdx] = solventInverseFormationVolumeFactor(cell_idx); + props.surf_dens_perf[well_.numComponents() * perf + Indices::contiSolventEqIdx] = solventRefDensity(cell_idx); } } } @@ -479,18 +469,13 @@ computeProperties(const WellState& well_state, const std::function& mobility, const std::function& solventInverseFormationVolumeFactor, const std::function& solventMobility, - const std::vector& b_perf, - const std::vector& rsmax_perf, - const std::vector& rvmax_perf, - const std::vector& rvwmax_perf, - const std::vector& rswmax_perf, - const std::vector& surf_dens_perf, + const Properties& props, DeferredLogger& deferred_logger) { // Compute densities const int nperf = well_.numPerfs(); const int np = well_.numPhases(); - std::vector perfRates(b_perf.size(),0.0); + std::vector perfRates(props.b_perf.size(),0.0); const auto& ws = well_state.well(well_.indexOfWell()); const auto& perf_data = ws.perf_data; const auto& perf_rates_state = perf_data.phase_rates; @@ -549,7 +534,7 @@ computeProperties(const WellState& well_state, } } - this->computeDensities(perfRates, b_perf, rsmax_perf, rvmax_perf, rvwmax_perf, rswmax_perf, surf_dens_perf, deferred_logger); + this->computeDensities(perfRates, props, deferred_logger); this->computePressureDelta(); } diff --git a/opm/simulators/wells/StandardWellConnections.hpp b/opm/simulators/wells/StandardWellConnections.hpp index 79a338bdf..45b6eb02c 100644 --- a/opm/simulators/wells/StandardWellConnections.hpp +++ b/opm/simulators/wells/StandardWellConnections.hpp @@ -55,12 +55,7 @@ public: const std::function& pvtRegionIdx, const std::function& solventInverseFormationVolumeFactor, const std::function& solventRefDensity, - std::vector& b_perf, - std::vector& rsmax_perf, - std::vector& rvmax_perf, - std::vector& rvwmax_perf, - std::vector& rswmax_perf, - std::vector& surf_dens_perf) const; + Properties& props) const; //! \brief Compute connection properties (densities, pressure drop, ...) void computeProperties(const WellState& well_state, @@ -68,12 +63,7 @@ public: const std::function& mobility, const std::function& solventInverseFormationVolumeFactor, const std::function& solventMobility, - const std::vector& b_perf, - const std::vector& rsmax_perf, - const std::vector& rvmax_perf, - const std::vector& rvwmax_perf, - const std::vector& rswmax_perf, - const std::vector& surf_dens_perf, + const Properties& props, DeferredLogger& deferred_logger); //! \brief Returns density for first perforation. @@ -92,12 +82,7 @@ private: // TODO: not total sure whether it is a good idea to put this function here // the major reason to put here is to avoid the usage of Wells struct void computeDensities(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& rswmax_perf, - const std::vector& surf_dens_perf, + const Properties& props, DeferredLogger& deferred_logger); const WellInterfaceIndices& well_; //!< Reference to well interface diff --git a/opm/simulators/wells/StandardWellEval.hpp b/opm/simulators/wells/StandardWellEval.hpp index ecd6b0e3f..2ca901864 100644 --- a/opm/simulators/wells/StandardWellEval.hpp +++ b/opm/simulators/wells/StandardWellEval.hpp @@ -49,6 +49,7 @@ class StandardWellEval { protected: using PrimaryVariables = StandardWellPrimaryVariables; + using StdWellConnections = StandardWellConnections; static constexpr int Bhp = PrimaryVariables::Bhp; static constexpr int WQTotal= PrimaryVariables::WQTotal; static constexpr int numWellConservationEq = PrimaryVariables::numWellConservationEq; @@ -102,7 +103,7 @@ protected: std::vector F0_; StandardWellEquations linSys_; //!< Linear equation system - StandardWellConnections connections_; //!< Connection level values + StdWellConnections connections_; //!< Connection level values }; } diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index e57c9f125..b3b628df3 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -1350,12 +1350,7 @@ namespace Opm StandardWell:: computePropertiesForWellConnectionPressures(const Simulator& ebosSimulator, const WellState& well_state, - std::vector& b_perf, - std::vector& rsmax_perf, - std::vector& rvmax_perf, - std::vector& rvwmax_perf, - std::vector& rswmax_perf, - std::vector& surf_dens_perf) const + WellConnectionProps& props) const { std::function getTemperature = [&ebosSimulator](int cell_idx, int phase_idx) @@ -1389,12 +1384,7 @@ namespace Opm getPvtRegionIdx, getInvFac, getSolventDensity, - b_perf, - rsmax_perf, - rvmax_perf, - rvwmax_perf, - rswmax_perf, - surf_dens_perf); + props); } @@ -1519,12 +1509,7 @@ namespace Opm StandardWell:: computeWellConnectionDensitesPressures(const Simulator& ebosSimulator, const WellState& well_state, - const std::vector& b_perf, - const std::vector& rsmax_perf, - const std::vector& rvmax_perf, - const std::vector& rvwmax_perf, - const std::vector& rswmax_perf, - const std::vector& surf_dens_perf, + const WellConnectionProps& props, DeferredLogger& deferred_logger) { std::function invB = @@ -1553,12 +1538,7 @@ namespace Opm mobility, invFac, solventMobility, - b_perf, - rsmax_perf, - rvmax_perf, - rvwmax_perf, - rswmax_perf, - surf_dens_perf, + props, deferred_logger); } @@ -1576,21 +1556,10 @@ namespace Opm // 1. Compute properties required by computePressureDelta(). // Note that some of the complexity of this part is due to the function // taking std::vector arguments, and not Eigen objects. - StdWellEval::StdWellConnections::Properties props; - computePropertiesForWellConnectionPressures(ebosSimulator, well_state, - props.b_perf, - props.rsmax_perf, - props.rvmax_perf, - props.rvwmax_perf, - props.rswmax_perf, - props.surf_dens_perf); + WellConnectionProps props; + computePropertiesForWellConnectionPressures(ebosSimulator, well_state, props); computeWellConnectionDensitesPressures(ebosSimulator, well_state, - props.b_perf, - props.rsmax_perf, - props.rvmax_perf, - props.rvwmax_perf, - props.rswmax_perf, - props.surf_dens_perf, deferred_logger); + props, deferred_logger); }