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

Calculate Segment Phase and Mixture Densities for Summary Output

Browse Source 
This commit adds logic and backing storage in the SegmentState to
provide the segment level summary vectors

  - SDENx -- Phase density of phase 'x' (O, G, W)
  - SDENM -- Mixture density without flowing fraction exponents
  - SMDEN -- Mixture density with flowing fraction exponents

We defer the calculation of SDENM and, especially, SMDEN, to the
MultisegmentWellSegments class since this class maintains the
current flowing fractions.
This commit is contained in:
Bård Skaflestad 2023-07-04 18:40:13 +02:00
parent 4ac439475b
commit 77fe28979e
6 changed files with 212 additions and 16 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

147
opm/simulators/wells/MultisegmentWellSegments.cpp
View File

@ -23,6 +23,7 @@
#include <opm/common/ErrorMacros.hpp>
#include <opm/input/eclipse/Schedule/MSW/AICD.hpp>
#include <opm/input/eclipse/Schedule/MSW/WellSegments.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
@ -37,10 +38,23 @@
#include <opm/simulators/utils/DeferredLogger.hpp>
#include <opm/simulators/wells/MSWellHelpers.hpp>
#include <opm/simulators/wells/SegmentState.hpp>
#include <opm/simulators/wells/WellInterfaceGeneric.hpp>
#include <opm/core/props/BlackoilPhases.hpp>
#include <fmt/format.h>
#include <algorithm>
#include <array>
#include <cmath>
#include <cstddef>
#include <functional>
#include <stdexcept>
#include <string>
#include <utility>
#include <vector>
namespace Opm
{
@ -317,17 +331,15 @@ updateUpwindingSegments(const PrimaryVariables& primary_variables)
}
}
template<class FluidSystem, class Indices, class Scalar>
typename MultisegmentWellSegments<FluidSystem,Indices,Scalar>::EvalWell
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
getHydroPressureLoss(const int seg,
const int seg_density) const
{
const int seg_density) const
{
return densities_[seg_density] * well_.gravity() * depth_diffs_[seg];
}
template<class FluidSystem, class Indices, class Scalar>
Scalar MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
getPressureDiffSegPerf(const int seg,
@ -779,6 +791,133 @@ accelerationPressureLoss(const int seg) const
return accelerationPressureLoss;
}
template <class FluidSystem, class Indices, class Scalar>
void
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
copyPhaseDensities(const PhaseUsage& pu, SegmentState& segSol) const
{
auto* rho = segSol.phase_density.data();
const auto phaseMap = std::vector {
std::pair { BlackoilPhases::Liquid, FluidSystem::oilPhaseIdx },
std::pair { BlackoilPhases::Vapour, FluidSystem::gasPhaseIdx },
std::pair { BlackoilPhases::Aqua , FluidSystem::waterPhaseIdx },
};
// Densities stored in 'rho' as
// [{ p0, p1, ..., (np - 1), mixture, mixture_with_exponents },
// { p0, p1, ..., (np - 1), mixture, mixture_with_exponents },
// ...
// { p0, p1, ..., (np - 1), mixture, mixture_with_exponents }]
// Stride is np + 2.
for (const auto& [boPhase, fsPhaseIdx] : phaseMap) {
if (pu.phase_used[boPhase]) {
this->copyPhaseDensities(fsPhaseIdx, pu.num_phases + 2,
rho + pu.phase_pos[boPhase]);
}
}
// Mixture densities.
for (auto seg = 0*this->densities_.size(); seg < this->densities_.size(); ++seg) {
const auto mixOffset = seg*(pu.num_phases + 2) + pu.num_phases;
rho[mixOffset + 0] = this->mixtureDensity(seg);
rho[mixOffset + 1] = this->mixtureDensityWithExponents(seg);
}
}
template <class FluidSystem, class Indices, class Scalar>
void
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
copyPhaseDensities(const unsigned phaseIdx,
const std::size_t stride,
double* dens) const
{
const auto compIdx = Indices::canonicalToActiveComponentIndex
(FluidSystem::solventComponentIndex(phaseIdx));
for (const auto& phase_density : this->phase_densities_) {
*dens = phase_density[compIdx].value();
dens += stride;
}
}
template <class FluidSystem, class Indices, class Scalar>
double
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
mixtureDensity(const int seg) const
{
auto mixDens = 0.0;
const auto& rho = this->phase_densities_[seg];
const auto& q = this->phase_fractions_[seg];
for (const auto& phIdx : {
FluidSystem::oilPhaseIdx,
FluidSystem::gasPhaseIdx,
FluidSystem::waterPhaseIdx
})
{
if (! FluidSystem::phaseIsActive(phIdx)) {
continue;
}
const auto compIdx = Indices::
canonicalToActiveComponentIndex(phIdx);
mixDens += q[compIdx].value() * rho[compIdx].value();
}
return mixDens;
}
template <class FluidSystem, class Indices, class Scalar>
double
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
mixtureDensityWithExponents(const int seg) const
{
if (const auto& segment = this->well_.wellEcl().getSegments()[seg];
segment.isAICD())
{
return this->mixtureDensityWithExponents(segment.autoICD(), seg);
}
// No other segment type includes exponents of flowing fractions when
// calculating the mixture/emulsion density.
return this->mixtureDensity(seg);
}
template <class FluidSystem, class Indices, class Scalar>
double
MultisegmentWellSegments<FluidSystem,Indices,Scalar>::
mixtureDensityWithExponents(const AutoICD& aicd, const int seg) const
{
auto mixDens = 0.0;
const auto& rho = this->phase_densities_[seg];
const auto& q = this->phase_fractions_[seg];
constexpr auto densityExponents = std::array {
std::pair { FluidSystem::oilPhaseIdx , &AutoICD::oilDensityExponent },
std::pair { FluidSystem::gasPhaseIdx , &AutoICD::gasDensityExponent },
std::pair { FluidSystem::waterPhaseIdx, &AutoICD::waterDensityExponent },
};
for (const auto& [fsPhaseIdx, densityExponent] : densityExponents) {
if (FluidSystem::phaseIsActive(fsPhaseIdx)) {
const auto compIdx = Indices::
canonicalToActiveComponentIndex(fsPhaseIdx);
// exp = (aicd.*densityExponent)() in native syntax.
const auto exp = std::invoke(densityExponent, aicd);
mixDens += std::pow(q[compIdx].value(), exp) * rho[compIdx].value();
}
}
return mixDens;
}
#define INSTANCE(...) \
template class MultisegmentWellSegments<BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,__VA_ARGS__,double>;

30
opm/simulators/wells/MultisegmentWellSegments.hpp
View File

@ -24,13 +24,20 @@
#include <opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp>
#include <cstddef>
#include <vector>
namespace Opm
{
namespace Opm {
class UnitSystem;
class WellInterfaceGeneric;
class AutoICD;
struct PhaseUsage;
class SegmentState;
class UnitSystem;
class WellInterfaceGeneric;
} // namespace Opm
namespace Opm {
template<typename FluidSystem, typename Indices, typename Scalar>
class MultisegmentWellSegments
@ -52,7 +59,7 @@ public:
void updateUpwindingSegments(const PrimaryVariables& primary_variables);
EvalWell getHydroPressureLoss(const int seg,
const int seg_side) const;
const int seg_side) const;
//! Pressure difference between segment and perforation.
Scalar getPressureDiffSegPerf(const int seg,
@ -109,6 +116,9 @@ public:
return perforation_depth_diffs_[perf];
}
void copyPhaseDensities(const PhaseUsage& pu,
SegmentState& segSol) const;
private:
// TODO: trying to use the information from the Well opm-parser as much
// as possible, it will possibly be re-implemented later for efficiency reason.
@ -149,8 +159,16 @@ private:
std::vector<std::vector<EvalWell>> phase_viscosities_;
const WellInterfaceGeneric& well_;
void copyPhaseDensities(const unsigned phaseIdx,
const std::size_t stride,
double* dens) const;
double mixtureDensity(const int seg) const;
double mixtureDensityWithExponents(const int seg) const;
double mixtureDensityWithExponents(const AutoICD& aicd, const int seg) const;
};
}
} // namespace Opm
#endif // OPM_MULTISEGMENTWELL_SEGMENTS_HEADER_INCLUDED

6
opm/simulators/wells/MultisegmentWell_impl.hpp
View File

@ -600,6 +600,12 @@ namespace Opm
this->primary_variables_.copyToWellState(*this, getRefDensity(), stop_or_zero_rate_target,
well_state, summary_state, deferred_logger);
{
auto& ws = well_state.well(this->index_of_well_);
this->segments_.copyPhaseDensities(ws.pu, ws.segments);
}
Base::calculateReservoirRates(well_state.well(this->index_of_well_));
}

15
opm/simulators/wells/SegmentState.cpp
View File

@ -61,6 +61,7 @@ SegmentState::SegmentState(int num_phases, const WellSegments& segments)
, phase_velocity (segments.size() * num_phases)
, phase_holdup (segments.size() * num_phases)
, phase_viscosity (segments.size() * num_phases)
, phase_density (segments.size() * (num_phases + 2)) // +2 for mixture with and without exponents
, pressure (segments.size())
, pressure_drop_friction (segments.size())
, pressure_drop_hydrostatic(segments.size())
@ -77,12 +78,13 @@ SegmentState SegmentState::serializationTestObject()
result.phase_resv_rates = {7.0, 8.0};
result.phase_velocity = {9.0};
result.phase_holdup = {10.0, 11.0};
result.phase_viscosity = {12.0};
result.pressure = {13.0, 14.0};
result.pressure_drop_friction = {15.0};
result.pressure_drop_hydrostatic = {16.0, 17.0};
result.pressure_drop_accel = {18.0};
result.m_segment_number = {19, 20};
result.phase_viscosity = {12.0, 12.5};
result.phase_density = {13.0, 13.5, 13.6, 13.75};
result.pressure = {14.0, 15.0};
result.pressure_drop_friction = {16.0};
result.pressure_drop_hydrostatic = {17.0, 18.0};
result.pressure_drop_accel = {19.0};
result.m_segment_number = {20, 21};
return result;
}
@ -124,6 +126,7 @@ bool SegmentState::operator==(const SegmentState& rhs) const
this->phase_velocity == rhs.phase_velocity &&
this->phase_holdup == rhs.phase_holdup &&
this->phase_viscosity == rhs.phase_viscosity &&
this->phase_density == rhs.phase_density &&
this->pressure == rhs.pressure &&
this->pressure_drop_friction == rhs.pressure_drop_friction &&
this->pressure_drop_hydrostatic == rhs.pressure_drop_hydrostatic &&

19
opm/simulators/wells/SegmentState.hpp
View File

@ -56,6 +56,7 @@ public:
serializer(phase_velocity);
serializer(phase_holdup);
serializer(phase_viscosity);
serializer(phase_density);
serializer(pressure);
serializer(pressure_drop_friction);
serializer(pressure_drop_hydrostatic);
@ -81,6 +82,24 @@ public:
/// Segment condition phase viscosities.
std::vector<double> phase_viscosity;
/// Segment condition phase densities.
///
/// \code num_phases + 2 \endcode elements per segment, with the two
/// additional elements being fluid mixture densities. Fluid mixture
/// densities are calculated with and without ICD-related flowing phase
/// fraction exponents, respectively. Mixture density *without* flowing
/// phase fraction exponents is stored at offset \c num_phases, and
/// mixture density *with* flowing phase fraction exponents is stored at
/// offset \code num_phases + 1 \endcode.
///
/// In other words, the phase and mixture densities are stored as
///
/// [{ p0, p1, ..., (np - 1), mixture, mixture_with_exponents },
/// { p0, p1, ..., (np - 1), mixture, mixture_with_exponents },
/// ...
/// { p0, p1, ..., (np - 1), mixture, mixture_with_exponents }]
std::vector<double> phase_density;
std::vector<double> pressure;
std::vector<double> pressure_drop_friction;
std::vector<double> pressure_drop_hydrostatic;

11
opm/simulators/wells/WellState.cpp
View File

@ -887,6 +887,7 @@ WellState::reportSegmentResults(const int well_id,
const int seg_no) const
{
using PhaseQuant = data::SegmentPhaseQuantity::Item;
using PhaseDensity = data::SegmentPhaseDensity::Item;
const auto& segments = this->well(well_id).segments;
if (segments.empty()) {
@ -911,6 +912,7 @@ WellState::reportSegmentResults(const int well_id,
const auto* velocity = &segments.phase_velocity[seg_ix * pu.num_phases];
const auto* holdup = &segments.phase_holdup[seg_ix * pu.num_phases];
const auto* viscosity = &segments.phase_viscosity[seg_ix * pu.num_phases];
const auto* density = &segments.phase_density[seg_ix * (pu.num_phases + 2)]; // +2 for mixture densities
if (pu.phase_used[Water]) {
const auto iw = pu.phase_pos[Water];
@ -920,6 +922,7 @@ WellState::reportSegmentResults(const int well_id,
seg_res.velocity.set(PhaseQuant::Water, velocity[iw]);
seg_res.holdup.set(PhaseQuant::Water, holdup[iw]);
seg_res.viscosity.set(PhaseQuant::Water, viscosity[iw]);
seg_res.density.set(PhaseDensity::Water, density[iw]);
}
if (pu.phase_used[Oil]) {
@ -931,6 +934,7 @@ WellState::reportSegmentResults(const int well_id,
seg_res.velocity.set(PhaseQuant::Oil, velocity[io]);
seg_res.holdup.set(PhaseQuant::Oil, holdup[io]);
seg_res.viscosity.set(PhaseQuant::Oil, viscosity[io]);
seg_res.density.set(PhaseDensity::Oil, density[io]);
}
if (pu.phase_used[Gas]) {
@ -942,10 +946,17 @@ WellState::reportSegmentResults(const int well_id,
seg_res.velocity.set(PhaseQuant::Gas, velocity[ig]);
seg_res.holdup.set(PhaseQuant::Gas, holdup[ig]);
seg_res.viscosity.set(PhaseQuant::Gas, viscosity[ig]);
seg_res.density.set(PhaseDensity::Gas, density[ig]);
}
seg_res.segNumber = seg_no;
// Recall: The mixture density *without* exponents is stored at offset
// 'num_phases' and the mixture density *with* exponents is stored at
// offset 'num_phases + 1'.
seg_res.density.set(PhaseDensity::Mixture, density[pu.num_phases]);
seg_res.density.set(PhaseDensity::MixtureWithExponents, density[pu.num_phases + 1]);
return seg_res;
}