Skip to content

Commit

Permalink
Save and Restart ICD Scaling Factor for SICDs Too
Browse files Browse the repository at this point in the history 
While here, extract common item handing for AICDs and SICDs out
to a separate helper routine to avoid repeating the logic on the
output side.
  • Loading branch information
@bska
bska committed Mar 22, 2024
1 parent 7aa17bc commit 9091496
Show file tree
Hide file tree
Showing 4 changed files with 94 additions and 69 deletions.
2 changes: 1 addition & 1 deletion opm/input/eclipse/Schedule/MSW/SICD.cpp
View file
Expand Up @@ -74,7 +74,7 @@ namespace Opm {
, m_method_flow_scaling (rstSegment.icd_scaling_mode)
, m_max_absolute_rate (rstSegment.max_valid_flow_rate)
, m_status (from_int<ICDStatus>(rstSegment.icd_status))
, m_scaling_factor (-1.0)
, m_scaling_factor (rstSegment.icd_scaling_factor)
{}

SICD::SICD(const double strength,
Expand Down
17 changes: 17 additions & 0 deletions opm/io/eclipse/rst/segment.cpp
View file
Expand Up @@ -33,7 +33,23 @@ double area_to_si(const Opm::UnitSystem& unit_system, const double raw_value)
return unit_system.to_si(M::length, unit_system.to_si(M::length, raw_value));
}

double
load_icd_scaling_factor(const Opm::UnitSystem& unit_system,
const int* iseg,
const double* rseg)
{
const auto scalingFactor = rseg[VI::RSeg::ScalingFactor];
const auto scalingMethod = iseg[VI::ISeg::ICDScalingMode];

if ((scalingMethod == 1) ||
((scalingMethod < 0) && (rseg[VI::RSeg::ICDLength] < 0.0)))
{
// Scaling factor is a length. Convert to SI.
return unit_system.to_si(M::length, scalingFactor);
}

// Scaling factor is a relative measure. No unit conversion needed.
return scalingFactor;
}

} // Anonymous namespace
Expand Down Expand Up @@ -72,6 +88,7 @@ Opm::RestartIO::RstSegment::RstSegment(const UnitSystem& unit_system,
, max_emulsion_ratio( rseg[VI::RSeg::MaxEmulsionRatio])
, max_valid_flow_rate( unit_system.to_si(M::rate, rseg[VI::RSeg::MaxValidFlowRate]))
, icd_length( unit_system.to_si(M::length, rseg[VI::RSeg::ICDLength]))
, icd_scaling_factor(load_icd_scaling_factor(unit_system, iseg, rseg))
, valve_area_fraction( rseg[VI::RSeg::ValveAreaFraction])
, aicd_flowrate_exponent( rseg[VI::RSeg::FlowRateExponent])
, aicd_viscosity_exponent( rseg[VI::RSeg::ViscFuncExponent])
Expand Down
1 change: 1 addition & 0 deletions opm/io/eclipse/rst/segment.hpp
View file
Expand Up @@ -68,6 +68,7 @@ struct RstSegment
double max_emulsion_ratio{};
double max_valid_flow_rate{};
double icd_length{};
double icd_scaling_factor{};
double valve_area_fraction{};

double aicd_flowrate_exponent{};
Expand Down
143 changes: 75 additions & 68 deletions opm/output/eclipse/AggregateMSWData.cpp
View file
Expand Up @@ -549,20 +549,20 @@ namespace {
// about their origins, other than the fact that they depend on
// the active unit conventions.
switch (uType) {
case UType::UNIT_TYPE_METRIC:
return 1.340e-15f;
case UType::UNIT_TYPE_METRIC:
return 1.340e-15f;

case UType::UNIT_TYPE_FIELD:
return 2.892e-14f;
case UType::UNIT_TYPE_FIELD:
return 2.892e-14f;

case UType::UNIT_TYPE_LAB:
return 7.615e-14f;
case UType::UNIT_TYPE_LAB:
return 7.615e-14f;

case UType::UNIT_TYPE_PVT_M:
return 1.322e-15f;
case UType::UNIT_TYPE_PVT_M:
return 1.322e-15f;

default:
break;
default:
break;
}

throw std::invalid_argument {
Expand Down Expand Up @@ -603,16 +603,14 @@ namespace {
}

template <class RSegArray>
void assignSpiralICDCharacteristics(const ::Opm::Segment& segment,
const ::Opm::UnitSystem& usys,
const int baseIndex,
RSegArray& rSeg)
void assignICDBaseCharacteristics(const ::Opm::SICD& sicd,
const ::Opm::UnitSystem& usys,
const int baseIndex,
RSegArray& rSeg)
{
using Ix = ::Opm::RestartIO::Helpers::VectorItems::RSeg::index;
using M = ::Opm::UnitSystem::measure;

const auto& sicd = segment.spiralICD();

rSeg[baseIndex + Ix::DeviceBaseStrength] =
usys.from_si(M::icd_strength, sicd.strength());

Expand All @@ -622,74 +620,85 @@ namespace {
rSeg[baseIndex + Ix::CalibrFluidViscosity] =
usys.from_si(M::viscosity, sicd.viscosityCalibration());

rSeg[baseIndex + Ix::CriticalWaterFraction] = sicd.criticalValue();
rSeg[baseIndex + Ix::CriticalWaterFraction] =
sicd.criticalValue();

rSeg[baseIndex + Ix::TransitionRegWidth] =
sicd.widthTransitionRegion();

rSeg[baseIndex + Ix::MaxEmulsionRatio] =
sicd.maxViscosityRatio();

const auto& max_rate = sicd.maxAbsoluteRate();
rSeg[baseIndex + Ix::MaxValidFlowRate] = max_rate.has_value() ? usys.from_si(M::geometric_volume_rate, max_rate.value()) : -1;

rSeg[baseIndex + Ix::ICDLength] =
usys.from_si(M::length, sicd.length());

// The scalingFactor's output value depends on the scaling
// method (item 11 of WSEGAICD/WSEGSICD). If either
//
// 1. Item 11 is 1, or
// 2. Item 11 is negative and the scalingFactor is negative
//
// then the scalingFactor is a length and needs unit conversion.
// Otherwise the scalingFactor is a dimensionless relative
// measure and does not need unit conversion.
const auto convertScalingFactor =
(sicd.methodFlowScaling() == 1) ||
((sicd.methodFlowScaling() < 0) && (sicd.length() < 0.0));

rSeg[baseIndex + Ix::ScalingFactor] = convertScalingFactor
? usys.from_si(M::length, sicd.scalingFactor())
: sicd.scalingFactor();
}

template <class RSegArray>
void assignAICDCharacteristics(const ::Opm::Segment& segment,
void assignSpiralICDCharacteristics(const ::Opm::Segment& segment,
const ::Opm::UnitSystem& usys,
const int baseIndex,
RSegArray& rSeg)
{
using Ix = ::Opm::RestartIO::Helpers::VectorItems::RSeg::index;
using M = ::Opm::UnitSystem::measure;
const double infinity = 1.0e+20;

const auto& aicd = segment.autoICD();

rSeg[baseIndex + Ix::DeviceBaseStrength] =
usys.from_si(M::aicd_strength, aicd.strength());
const auto& sicd = segment.spiralICD();

// The value of the scalingFactor depends on the option used:
// if 1. item 11 on the WSEGAICD keyword is 1, or
// 2. item 11 is negative plus the value of the scalingFactor is negative then
// the scalingFactor is an absolute length and need unit conversion
// In other cases the scalingFactor is a relative length - and is unitless and do not need unit conversion
//
rSeg[baseIndex + Ix::ScalingFactor] = ((aicd.methodFlowScaling() == 1) ||
((aicd.methodFlowScaling() < 0) && (aicd.length() < 0))) ?
usys.from_si(M::length, aicd.scalingFactor()) : aicd.scalingFactor();
assignICDBaseCharacteristics(sicd, usys, baseIndex, rSeg);

rSeg[baseIndex + Ix::CalibrFluidDensity] =
usys.from_si(M::density, aicd.densityCalibration());
{
const auto& max_rate = sicd.maxAbsoluteRate();

rSeg[baseIndex + Ix::CalibrFluidViscosity] =
usys.from_si(M::viscosity, aicd.viscosityCalibration());
rSeg[baseIndex + Ix::MaxValidFlowRate] = max_rate.has_value()
? usys.from_si(M::geometric_volume_rate, max_rate.value())
: -1.0;
}
}

rSeg[baseIndex + Ix::CriticalWaterFraction] = aicd.criticalValue();
template <class RSegArray>
void assignAICDCharacteristics(const ::Opm::Segment& segment,
const ::Opm::UnitSystem& usys,
const int baseIndex,
RSegArray& rSeg)
{
using Ix = ::Opm::RestartIO::Helpers::VectorItems::RSeg::index;
using M = ::Opm::UnitSystem::measure;
const double infinity = 1.0e+20;

rSeg[baseIndex + Ix::TransitionRegWidth] =
aicd.widthTransitionRegion();
const auto& aicd = segment.autoICD();

rSeg[baseIndex + Ix::MaxEmulsionRatio] =
aicd.maxViscosityRatio();
assignICDBaseCharacteristics(aicd, usys, baseIndex, rSeg);

rSeg[baseIndex + Ix::FlowRateExponent] =
aicd.flowRateExponent();

rSeg[baseIndex + Ix::ViscFuncExponent] =
aicd.viscExponent();

const auto& max_rate = aicd.maxAbsoluteRate();
if (max_rate.has_value())
rSeg[baseIndex + Ix::MaxValidFlowRate] = usys.from_si(M::geometric_volume_rate, max_rate.value());
else
rSeg[baseIndex + Ix::MaxValidFlowRate] = -2 * infinity;
{
const auto& max_rate = aicd.maxAbsoluteRate();

rSeg[baseIndex + Ix::ICDLength] =
usys.from_si(M::length, aicd.length());
rSeg[baseIndex + Ix::MaxValidFlowRate] = max_rate.has_value()
? usys.from_si(M::geometric_volume_rate, max_rate.value())
: -2 * infinity;
}

rSeg[baseIndex + Ix::flowFractionOilDensityExponent] =
aicd.oilDensityExponent();
Expand All @@ -708,10 +717,8 @@ namespace {

rSeg[baseIndex + Ix::flowFractionGasViscosityExponent] =
aicd.gasViscExponent();

}


template <class RSegArray>
void assignSegmentTypeCharacteristics(const ::Opm::Segment& segment,
const ::Opm::UnitSystem& usys,
Expand All @@ -732,9 +739,9 @@ namespace {
}

template <class RSegArray>
void assignTracerData(std::size_t segment_offset,
void assignTracerData(const std::size_t segment_offset,
const Opm::Runspec& runspec,
RSegArray& rSeg)
RSegArray& rSeg)
{
auto tracer_offset = segment_offset + Opm::RestartIO::InteHEAD::numRsegElem(runspec.phases());
auto tracer_end = tracer_offset + runspec.tracers().water_tracers() * 8;
Expand All @@ -743,14 +750,14 @@ namespace {
}

template <class RSegArray>
void staticContrib_useMSW(const Opm::Runspec& runspec,
const Opm::Well& well,
const std::vector<int>& inteHead,
const Opm::EclipseGrid& grid,
const Opm::UnitSystem& units,
const ::Opm::SummaryState& smry,
const Opm::data::Wells& wr,
RSegArray& rSeg)
void staticContrib(const Opm::Runspec& runspec,
const Opm::Well& well,
const std::vector<int>& inteHead,
const Opm::EclipseGrid& grid,
const Opm::UnitSystem& units,
const Opm::SummaryState& smry,
const Opm::data::Wells& wr,
RSegArray& rSeg)
{
using Ix = ::Opm::RestartIO::Helpers::VectorItems::RSeg::index;
if (well.isMultiSegment()) {
Expand Down Expand Up @@ -1333,12 +1340,12 @@ captureDeclaredMSWData(const Schedule& sched,
MSWLoop(msw, [&units, &inteHead, &sched, &grid, &smry, &wr, this]
(const Well& well, const std::size_t mswID)
{
auto imsw = this->iSeg_[mswID];
auto rmsw = this->rSeg_[mswID];
auto iseg = this->iSeg_[mswID];
auto rseg = this->rSeg_[mswID];

ISeg::staticContrib(well, inteHead, imsw);
RSeg::staticContrib_useMSW(sched.runspec(), well, inteHead,
grid, units, smry, wr, rmsw);
ISeg::staticContrib(well, inteHead, iseg);
RSeg::staticContrib(sched.runspec(), well, inteHead,
grid, units, smry, wr, rseg);
});

// Extract contributions to the ILBS and ILBR arrays.
Expand Down

0 comments on commit 9091496

Please sign in to comment.