Merge pull request #3163 from bska/segment-node-xy

Add XY Coordinates to Segment's Node
This commit is contained in:
Markus Blatt 2022-10-05 09:52:39 +02:00 committed by GitHub
commit 0e341ef6e1
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
5 changed files with 690 additions and 290 deletions

View File

@ -72,11 +72,24 @@ namespace Opm {
Segment();
Segment(const Segment& src, double new_depth, double new_length, double new_volume, double new_x, double new_y);
Segment(const Segment& src, double new_depth, double new_length, double new_x, double new_y);
Segment(const Segment& src, double new_depth, double new_length, double new_volume);
Segment(const Segment& src, double new_depth, double new_length);
Segment(const Segment& src, double new_volume);
Segment(int segment_number_in, int branch_in, int outlet_segment_in, double length_in, double depth_in,
double internal_diameter_in, double roughness_in, double cross_area_in, double volume_in, bool data_ready_in);
Segment(const int segment_number_in,
const int branch_in,
const int outlet_segment_in,
const double length_in,
const double depth_in,
const double internal_diameter_in,
const double roughness_in,
const double cross_area_in,
const double volume_in,
const bool data_ready_in,
const double x_in,
const double y_in);
Segment(const RestartIO::RstSegment& rst_segment);
static Segment serializationTestObject();
@ -86,6 +99,8 @@ namespace Opm {
int outletSegment() const;
double perfLength() const;
double totalLength() const;
double node_X() const;
double node_Y() const;
double depth() const;
double internalDiameter() const;
double roughness() const;
@ -150,6 +165,8 @@ namespace Opm {
serializer(m_cross_area);
serializer(m_volume);
serializer(m_data_ready);
serializer(m_x);
serializer(m_y);
serializer(m_perf_length);
serializer(m_icd);
}
@ -200,14 +217,18 @@ namespace Opm {
// indicate if the data related to 'INC' or 'ABS' is ready
// the volume will be updated at a final step.
bool m_data_ready;
// Length of segment projected onto the X axis. Not used in
// simulations, but needed for the SEG option in WRFTPLT.
double m_x{};
// Length of segment projected onto the Y axis. Not used in
// simulations, but needed for the SEG option in WRFTPLT.
double m_y{};
std::optional<double> m_perf_length;
std::variant<RegularSegment, SICD, AutoICD, Valve> m_icd;
// We are not handling the length of segment projected onto the X-axis and Y-axis.
// They are not used in the simulations and we are not supporting the plotting.
// There are other three properties for segment related to thermal conduction,
// while they are not supported by the keyword at the moment.
// There are three other properties for the segment pertaining to
// thermal conduction. These are not currently supported.
};
}

View File

@ -120,23 +120,22 @@ namespace Opm {
void processINC(double depth_top, double length_top);
void process(LengthDepth length_depth, double depth_top, double length_top);
void addSegment(const Segment& new_segment);
void addSegment(int segment_number,
int branch,
int outlet_segment,
double length,
double depth,
double internal_diameter,
double roughness,
double cross_area,
double volume,
bool data_ready);
void addSegment(const int segment_number,
const int branch,
const int outlet_segment,
const double length,
const double depth,
const double internal_diameter,
const double roughness,
const double cross_area,
const double volume,
const bool data_ready,
const double node_x,
const double node_y);
const Segment& topSegment() const;
// components of the pressure drop to be included
CompPressureDrop m_comp_pressure_drop;
// There are X and Y cooridnate of the nodal point of the top segment
// Since they are not used for simulations and we are not supporting plotting,
// we are not handling them at the moment.
// There are other three properties for segment related to thermal conduction,
// while they are not supported by the keyword at the moment.

View File

@ -17,8 +17,10 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <opm/io/eclipse/rst/segment.hpp>
#include <opm/input/eclipse/Schedule/MSW/Segment.hpp>
#include <opm/io/eclipse/rst/segment.hpp>
#include <opm/input/eclipse/Schedule/MSW/SICD.hpp>
#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
@ -26,52 +28,49 @@
#include <stdexcept>
#include <string>
namespace Opm {
namespace {
static constexpr double invalid_value = -1.e100;
static constexpr double invalid_value = -1.0e100;
}
double if_invalid_value(const double rst_value)
{
return (rst_value == 0.0)
? invalid_value
: rst_value;
}
} // Anonymous namespace
namespace Opm {
Segment::Segment()
: m_segment_number(-1),
m_branch(-1),
m_outlet_segment(-1),
m_total_length(invalid_value),
m_depth(invalid_value),
m_internal_diameter(invalid_value),
m_roughness(invalid_value),
m_cross_area(invalid_value),
m_volume(invalid_value),
m_data_ready(false)
{
}
: m_segment_number(-1)
, m_branch(-1)
, m_outlet_segment(-1)
, m_total_length(invalid_value)
, m_depth(invalid_value)
, m_internal_diameter(invalid_value)
, m_roughness(invalid_value)
, m_cross_area(invalid_value)
, m_volume(invalid_value)
, m_data_ready(false)
, m_x(0.0)
, m_y(0.0)
{}
namespace {
double if_invalid_value(double rst_value) {
if (rst_value == 0)
return invalid_value;
return rst_value;
}
}
Segment::Segment(const RestartIO::RstSegment& rst_segment):
m_segment_number(rst_segment.segment),
m_branch(rst_segment.branch),
m_outlet_segment(rst_segment.outlet_segment),
m_total_length( rst_segment.dist_bhp_ref ),
m_depth(rst_segment.bhp_ref_dz),
m_internal_diameter(if_invalid_value(rst_segment.diameter)),
m_roughness(if_invalid_value(rst_segment.roughness)),
m_cross_area(if_invalid_value(rst_segment.area)),
m_volume(rst_segment.volume),
m_data_ready(true)
Segment::Segment(const RestartIO::RstSegment& rst_segment)
: m_segment_number(rst_segment.segment)
, m_branch(rst_segment.branch)
, m_outlet_segment(rst_segment.outlet_segment)
, m_total_length( rst_segment.dist_bhp_ref )
, m_depth(rst_segment.bhp_ref_dz)
, m_internal_diameter(if_invalid_value(rst_segment.diameter))
, m_roughness(if_invalid_value(rst_segment.roughness))
, m_cross_area(if_invalid_value(rst_segment.area))
, m_volume(rst_segment.volume)
, m_data_ready(true)
, m_x(0.0)
, m_y(0.0)
{
if (rst_segment.segment_type == SegmentType::SICD) {
double scalingFactor = -1; // The scaling factor will be and updated from the simulator.
@ -116,81 +115,119 @@ namespace {
}
}
Segment::Segment(const int segment_number_in,
const int branch_in,
const int outlet_segment_in,
const double length_in,
const double depth_in,
const double internal_diameter_in,
const double roughness_in,
const double cross_area_in,
const double volume_in,
const bool data_ready_in,
const double x_in,
const double y_in)
: m_segment_number(segment_number_in)
, m_branch(branch_in)
, m_outlet_segment(outlet_segment_in)
, m_total_length(length_in)
, m_depth(depth_in)
, m_internal_diameter(internal_diameter_in)
, m_roughness(roughness_in)
, m_cross_area(cross_area_in)
, m_volume(volume_in)
, m_data_ready(data_ready_in)
, m_x(x_in)
, m_y(y_in)
{}
Segment::Segment(int segment_number_in, int branch_in, int outlet_segment_in, double length_in, double depth_in,
double internal_diameter_in, double roughness_in, double cross_area_in,
double volume_in, bool data_ready_in)
: m_segment_number(segment_number_in),
m_branch(branch_in),
m_outlet_segment(outlet_segment_in),
m_total_length(length_in),
m_depth(depth_in),
m_internal_diameter(internal_diameter_in),
m_roughness(roughness_in),
m_cross_area(cross_area_in),
m_volume(volume_in),
m_data_ready(data_ready_in)
Segment::Segment(const Segment& src,
const double new_depth,
const double new_length,
const double new_volume,
const double new_x,
const double new_y)
: Segment { src, new_depth, new_length, new_volume }
{
this->m_x = new_x;
this->m_y = new_y;
}
Segment::Segment(const Segment& src, double new_depth, double new_length):
Segment(src)
{
this->m_depth = new_depth;
this->m_total_length = new_length;
this->m_data_ready = true;
}
Segment::Segment(const Segment& src,
const double new_depth,
const double new_length,
const double new_x,
const double new_y)
: Segment { src, new_depth, new_length }
{
this->m_x = new_x;
this->m_y = new_y;
}
Segment::Segment(const Segment& src, double new_depth, double new_length, double new_volume):
Segment(src, new_depth, new_length)
{
this->m_volume = new_volume;
}
Segment::Segment(const Segment& src,
const double new_depth,
const double new_length,
const double new_volume)
: Segment { src, new_depth, new_length }
{
this->m_volume = new_volume;
}
Segment::Segment(const Segment& src, double new_volume):
Segment(src)
{
this->m_volume = new_volume;
}
Segment::Segment(const Segment& src,
const double new_depth,
const double new_length)
: Segment(src)
{
this->m_depth = new_depth;
this->m_total_length = new_length;
this->m_data_ready = true;
}
Segment Segment::serializationTestObject()
{
Segment result;
result.m_segment_number = 1;
result.m_branch = 2;
result.m_outlet_segment = 3;
result.m_inlet_segments = {4, 5};
result.m_total_length = 6.0;
result.m_depth = 7.0;
result.m_internal_diameter = 8.0;
result.m_roughness = 9.0;
result.m_cross_area = 10.0;
result.m_volume = 11.0;
result.m_data_ready = true;
result.m_icd = SICD::serializationTestObject();
return result;
}
Segment::Segment(const Segment& src, const double new_volume)
: Segment(src)
{
this->m_volume = new_volume;
}
Segment Segment::serializationTestObject()
{
Segment result;
result.m_segment_number = 1;
result.m_branch = 2;
result.m_outlet_segment = 3;
result.m_inlet_segments = {4, 5};
result.m_total_length = 6.0;
result.m_depth = 7.0;
result.m_internal_diameter = 8.0;
result.m_roughness = 9.0;
result.m_cross_area = 10.0;
result.m_volume = 11.0;
result.m_data_ready = true;
result.m_x = 12.0;
result.m_y = 13.0;
result.m_perf_length = 14.0;
result.m_icd = SICD::serializationTestObject();
return result;
}
int Segment::segmentNumber() const {
return m_segment_number;
}
int Segment::branchNumber() const {
return m_branch;
}
int Segment::outletSegment() const {
return m_outlet_segment;
}
double Segment::totalLength() const {
return m_total_length;
}
double Segment::node_X() const { return this->m_x; }
double Segment::node_Y() const { return this->m_y; }
double Segment::depth() const {
return m_depth;
@ -204,12 +241,10 @@ namespace {
return m_internal_diameter;
}
double Segment::roughness() const {
return m_roughness;
}
double Segment::crossArea() const {
return m_cross_area;
}
@ -262,7 +297,10 @@ namespace {
&& this->m_volume == rhs.m_volume
&& this->m_perf_length == rhs.m_perf_length
&& this->m_icd == rhs.m_icd
&& this->m_data_ready == rhs.m_data_ready;
&& this->m_data_ready == rhs.m_data_ready
&& (this->m_x == rhs.m_x)
&& (this->m_y == rhs.m_y)
;
}
bool Segment::operator!=( const Segment& rhs ) const {
@ -285,8 +323,6 @@ namespace {
return std::get<AutoICD>(this->m_icd);
}
void Segment::updateValve(const Valve& input_valve) {
// we need to update some values for the vale
auto valve = input_valve;
@ -318,7 +354,6 @@ namespace {
this->m_icd= valve;
}
void Segment::updateValve__(Valve& valve, const double segment_length) {
if (valve.pipeAdditionalLength() < 0)
valve.setPipeAdditionalLength(segment_length);
@ -331,11 +366,9 @@ namespace {
this->updateValve__(new_valve, segment_length);
}
void Segment::updatePerfLength(double perf_length) {
this->m_perf_length = perf_length;
}
void Segment::updatePerfLength(double perf_length) {
this->m_perf_length = perf_length;
}
const Valve& Segment::valve() const {
return std::get<Valve>(this->m_icd);
@ -370,6 +403,5 @@ namespace {
throw std::invalid_argument("Unhandeled segment type: " + std::to_string(ecl_id));
}
}
}
} // namespace Opm

View File

@ -16,6 +16,18 @@
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <opm/input/eclipse/Schedule/MSW/WellSegments.hpp>
#include <opm/input/eclipse/Schedule/MSW/Segment.hpp>
#include <opm/input/eclipse/Schedule/MSW/SICD.hpp>
#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
#include <opm/input/eclipse/Deck/DeckItem.hpp>
#include <opm/input/eclipse/Deck/DeckKeyword.hpp>
#include <opm/input/eclipse/Deck/DeckRecord.hpp>
#include <cassert>
#include <cmath>
#include <iostream>
@ -30,15 +42,6 @@
#include <fmt/format.h>
#include <opm/input/eclipse/Deck/DeckItem.hpp>
#include <opm/input/eclipse/Deck/DeckKeyword.hpp>
#include <opm/input/eclipse/Deck/DeckRecord.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
#include <opm/input/eclipse/Schedule/MSW/Segment.hpp>
#include <opm/input/eclipse/Schedule/MSW/SICD.hpp>
#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
#include <opm/input/eclipse/Schedule/MSW/WellSegments.hpp>
namespace Opm {
WellSegments::WellSegments(CompPressureDrop compDrop,
@ -112,45 +115,55 @@ namespace Opm {
}
}
void WellSegments::addSegment( const Segment& new_segment ) {
// decide whether to push_back or insert
const int segment_number = new_segment.segmentNumber();
const int segment_index = segmentNumberToIndex(segment_number);
void WellSegments::addSegment(const Segment& new_segment)
{
// Decide whether to push_back or insert
const auto segment_number = new_segment.segmentNumber();
const auto segment_index = segmentNumberToIndex(segment_number);
if (segment_index < 0) {
// New segment object.
const auto new_index = static_cast<int>(this->size());
this->segment_number_to_index.insert_or_assign(segment_number, new_index);
this->m_segments.push_back(new_segment);
}
else {
// Update to existing segment object.
this->m_segments[segment_index] = new_segment;
}
}
if (segment_index < 0) { // it is a new segment
segment_number_to_index[segment_number] = this->size();
m_segments.push_back(new_segment);
} else
m_segments[segment_index] = new_segment;
}
void WellSegments::addSegment(const int segment_number,
const int branch,
const int outlet_segment,
const double length,
const double depth,
const double internal_diameter,
const double roughness,
const double cross_area,
const double volume,
const bool data_ready,
const double node_x,
const double node_y)
{
const auto segment = Segment {
segment_number, branch, outlet_segment,
length, depth, internal_diameter, roughness,
cross_area, volume,
data_ready, node_x, node_y
};
void WellSegments::addSegment(int segment_number,
int branch,
int outlet_segment,
double length,
double depth,
double internal_diameter,
double roughness,
double cross_area,
double volume,
bool data_ready) {
Segment segment(segment_number, branch, outlet_segment, length, depth, internal_diameter, roughness, cross_area, volume, data_ready);
const int segment_index = this->segmentNumberToIndex(segment_number);
if (segment_index < 0) {
segment_number_to_index[segment_number] = this->size();
this->m_segments.push_back( std::move(segment) );
} else
this->m_segments[segment_index] = std::move(segment);
}
this->addSegment(segment);
}
void WellSegments::loadWELSEGS( const DeckKeyword& welsegsKeyword ) {
// for the first record, which provides the information for the top segment
// and information for the whole segment set
void WellSegments::loadWELSEGS(const DeckKeyword& welsegsKeyword)
{
// For the first record, which provides the information for the top
// segment and information for the whole segment set.
const auto& record1 = welsegsKeyword.getRecord(0);
const double invalid_value = Segment::invalidValue(); // meaningless value to indicate unspecified values
// Meaningless value to indicate unspecified values.
const double invalid_value = Segment::invalidValue();
const double depth_top = record1.getItem("DEPTH").getSIDouble(0);
const double length_top = record1.getItem("LENGTH").getSIDouble(0);
@ -158,121 +171,176 @@ namespace Opm {
const LengthDepth length_depth_type = LengthDepthFromString(record1.getItem("INFO_TYPE").getTrimmedString(0));
m_comp_pressure_drop = CompPressureDropFromString(record1.getItem("PRESSURE_COMPONENTS").getTrimmedString(0));
// the main branch is 1 instead of 0
// the segment number for top segment is also 1
if (length_depth_type == LengthDepth::INC)
this->addSegment( 1, 1, 0, 0., 0.,
invalid_value, invalid_value, invalid_value,
volume_top, false);
const auto nodeX_top = record1.getItem("TOP_X").getSIDouble(0);
const auto nodeY_top = record1.getItem("TOP_Y").getSIDouble(0);
else if (length_depth_type == LengthDepth::ABS)
this->addSegment( 1, 1, 0, length_top, depth_top,
invalid_value, invalid_value, invalid_value,
volume_top, true);
// The main branch is 1 instead of 0. The segment number for top
// segment is also 1.
if (length_depth_type == LengthDepth::INC) {
const auto segmentID = 1;
const auto branchID = 1;
const auto outletSegment = 0;
const auto length = 0.0;
const auto depth = 0.0;
const auto internal_diameter = invalid_value;
const auto roughness = invalid_value;
const auto cross_area = invalid_value;
const auto data_ready = false;
// read all the information out from the DECK first then process to get all the required information
this->addSegment(segmentID, branchID, outletSegment, length, depth,
internal_diameter, roughness, cross_area,
volume_top, data_ready, nodeX_top, nodeY_top);
}
else if (length_depth_type == LengthDepth::ABS) {
const auto segmentID = 1;
const auto branchID = 1;
const auto outletSegment = 0;
const auto internal_diameter = invalid_value;
const auto roughness = invalid_value;
const auto cross_area = invalid_value;
const auto data_ready = true;
this->addSegment(segmentID, branchID, outletSegment,
length_top, depth_top,
internal_diameter, roughness, cross_area,
volume_top, data_ready, nodeX_top, nodeY_top);
}
// Read all the information out from the DECK first then process to
// get all the requisite information.
for (size_t recordIndex = 1; recordIndex < welsegsKeyword.size(); ++recordIndex) {
const auto& record = welsegsKeyword.getRecord(recordIndex);
const int segment1 = record.getItem("SEGMENT1").get< int >(0);
const int segment2 = record.getItem("SEGMENT2").get< int >(0);
if ((segment1 < 2) || (segment2 < segment1))
throw std::logic_error("illegal segment number input is found in WELSEGS!\n");
const int segment1 = record.getItem("SEGMENT1").get<int>(0);
const int segment2 = record.getItem("SEGMENT2").get<int>(0);
if (segment1 < 2) {
throw std::logic_error {
fmt::format("Illegal segment 1 number in WELSEGS\n"
"Expected 2..NSEGMX, but got {}", segment1)
};
}
const int branch = record.getItem("BRANCH").get< int >(0);
if ((branch < 1))
throw std::logic_error("illegal branch number input is found in WELSEGS!\n");
if (segment2 < segment1) {
throw std::logic_error {
fmt::format("Illegal segment 2 number in WELSEGS\n"
"Expected {}..NSEGMX, but got {}",
segment1, segment2)
};
}
const int branch = record.getItem("BRANCH").get<int>(0);
if (branch < 1) {
throw std::logic_error {
fmt::format("Illegal branch number input "
"({}) is found in WELSEGS!", branch)
};
}
const double diameter = record.getItem("DIAMETER").getSIDouble(0);
double area = M_PI * diameter * diameter / 4.0;
{
const auto& itemArea = record.getItem("AREA");
if (itemArea.hasValue(0))
if (itemArea.hasValue(0)) {
area = itemArea.getSIDouble(0);
}
}
// if the values are incremental values, then we can just use the values
// if the values are absolute values, then we need to calculate them during the next process
// only the value for the last segment in the range is recorded
// If the values are incremental values, then we can just use
// the values if the values are absolute values, then we need to
// calculate them during the next process only the value for the
// last segment in the range is recorded
const double segment_length = record.getItem("SEGMENT_LENGTH").getSIDouble(0);
// the naming is a little confusing here.
// naming following the definition from the current keyword for the moment
// The naming is a little confusing here. Naming following the
// definition from the current keyword for the moment.
const double depth_change = record.getItem("DEPTH_CHANGE").getSIDouble(0);
double volume = invalid_value;
{
const auto& itemVolume = record.getItem("VOLUME");
if (itemVolume.hasValue(0))
if (itemVolume.hasValue(0)) {
volume = itemVolume.getSIDouble(0);
else if (length_depth_type == LengthDepth::INC)
}
else if (length_depth_type == LengthDepth::INC) {
volume = area * segment_length;
}
}
const double roughness = record.getItem("ROUGHNESS").getSIDouble(0);
for (int segment_number = segment1; segment_number <= segment2; ++segment_number) {
// for the first or the only segment in the range is the one specified in the WELSEGS
// from the second segment in the range, the outlet segment is the previous segment in the range
int outlet_segment;
if (segment_number == segment1)
outlet_segment = record.getItem("JOIN_SEGMENT").get< int >(0);
else
outlet_segment = segment_number - 1;
const auto node_X = record.getItem("LENGTH_X").getSIDouble(0);
const auto node_Y = record.getItem("LENGTH_Y").getSIDouble(0);
if (length_depth_type == LengthDepth::INC) {
this->addSegment( segment_number, branch, outlet_segment, segment_length, depth_change,
diameter, roughness, area, volume, false);
} else if (segment_number == segment2) {
this->addSegment( segment_number, branch, outlet_segment, segment_length, depth_change,
diameter, roughness, area, volume, true);
} else {
this->addSegment( segment_number, branch, outlet_segment, invalid_value, invalid_value,
diameter, roughness, area, volume, false);
}
for (int segment_number = segment1; segment_number <= segment2; ++segment_number) {
// For the first or the only segment in the range is the one
// specified in WELSEGS. From the second segment in the
// range, the outlet segment is the previous segment in the
// range.
const int outlet_segment = (segment_number == segment1)
? record.getItem("JOIN_SEGMENT").get<int>(0)
: segment_number - 1;
const auto data_ready = (length_depth_type != LengthDepth::INC)
&& (segment_number == segment2);
this->addSegment(segment_number, branch, outlet_segment,
segment_length, depth_change, diameter,
roughness, area, volume, data_ready,
node_X, node_Y);
}
}
for (size_t i_segment = 0; i_segment < m_segments.size(); ++i_segment) {
const int segment_number = m_segments[i_segment].segmentNumber();
const int outlet_segment = m_segments[i_segment].outletSegment();
for (const auto& segment : this->m_segments) {
const int outlet_segment = segment.outletSegment();
if (outlet_segment <= 0) { // no outlet segment
continue;
}
const int outlet_segment_index = segment_number_to_index[outlet_segment];
m_segments[outlet_segment_index].addInletSegment(segment_number);
m_segments[outlet_segment_index].addInletSegment(segment.segmentNumber());
}
this->process(length_depth_type, depth_top, length_top);
}
const Segment& WellSegments::getFromSegmentNumber(const int segment_number) const {
// the index of segment in the vector of segments
const int segment_index = segmentNumberToIndex(segment_number);
if (segment_index < 0) {
throw std::runtime_error("Could not indexate the segment " + std::to_string(segment_number)
+ " when trying to get the segment ");
throw std::runtime_error {
fmt::format("Unknown segment number {}", segment_number)
};
}
return m_segments[segment_index];
}
void WellSegments::process(LengthDepth length_depth, double depth_top, double length_top) {
if (length_depth == LengthDepth::ABS)
void WellSegments::process(const LengthDepth length_depth,
const double depth_top,
const double length_top)
{
if (length_depth == LengthDepth::ABS) {
this->processABS();
else if (length_depth == LengthDepth::INC)
}
else if (length_depth == LengthDepth::INC) {
this->processINC(depth_top, length_top);
else
throw std::logic_error("Invalid llength/depth/type in segment data structure");
}
else {
throw std::logic_error {
fmt::format("Invalid length/depth type "
"{} in segment data structure",
static_cast<int>(length_depth))
};
}
}
void WellSegments::processABS() {
const double invalid_value = Segment::invalidValue(); // meaningless value to indicate unspecified/uncompleted values
void WellSegments::processABS()
{
// Meaningless value to indicate unspecified/uncompleted values
const double invalid_value = Segment::invalidValue();
orderSegments();
std::size_t current_index= 1;
while (current_index< size()) {
std::size_t current_index = 1;
while (current_index < size()) {
if (m_segments[current_index].dataReady()) {
current_index++;
continue;
@ -280,7 +348,7 @@ namespace Opm {
const int range_begin = current_index;
const int outlet_segment = m_segments[range_begin].outletSegment();
const int outlet_index= segmentNumberToIndex(outlet_segment);
const int outlet_index = segmentNumberToIndex(outlet_segment);
assert(m_segments[outlet_index].dataReady() == true);
@ -310,30 +378,47 @@ namespace Opm {
const double depth_inc = (depth_last - depth_outlet) / number_segments;
const double volume_segment = m_segments[range_end].crossArea() * length_inc;
const auto x_outlet = m_segments[outlet_index].node_X();
const auto y_outlet = m_segments[outlet_index].node_Y();
const auto dx = (m_segments[range_end].node_X() - x_outlet) / number_segments;
const auto dy = (m_segments[range_end].node_Y() - y_outlet) / number_segments;
for (std::size_t k = range_begin; k <= range_end; ++k) {
const auto& old_segment = this->m_segments[k];
double new_volume, new_length, new_depth;
double new_length, new_depth, new_x, new_y;
if (k == range_end) {
new_length = old_segment.totalLength();
new_depth = old_segment.depth();
} else {
new_length = length_outlet + (k - range_begin + 1) * length_inc;
new_depth = depth_outlet + (k - range_end + 1) * depth_inc;
new_depth = old_segment.depth();
new_x = old_segment.node_X();
new_y = old_segment.node_Y();
}
else {
const auto num_inc = k - range_begin + 1;
new_length = length_outlet + num_inc*length_inc;
new_depth = depth_outlet + num_inc*depth_inc;
new_x = x_outlet + num_inc*dx;
new_y = y_outlet + num_inc*dy;
}
if (old_segment.volume() < 0.5 * invalid_value)
new_volume = volume_segment;
else
new_volume = old_segment.volume();
const auto new_volume = (old_segment.volume() < 0.5 * invalid_value)
? volume_segment
: old_segment.volume();
Segment new_segment(old_segment, new_length, new_depth, new_volume);
addSegment(new_segment);
this->addSegment(Segment {
old_segment, new_length, new_depth,
new_volume, new_x, new_y
});
}
current_index= range_end + 1;
current_index = range_end + 1;
}
// then update the volume for all the segments except the top segment
// this is for the segments specified individually while the volume is not specified.
// Then update the volume for all the segments except the top
// segment. This is for the segments specified individually while
// the volume is not specified.
for (std::size_t i = 1; i < size(); ++i) {
assert(m_segments[i].dataReady());
if (m_segments[i].volume() == invalid_value) {
@ -349,52 +434,65 @@ namespace Opm {
}
}
void WellSegments::processINC(double depth_top, double length_top) {
// update the information inside the WellSegments to be in ABS way
void WellSegments::processINC(double depth_top, double length_top)
{
// Update the information inside the WellSegments to be in ABS way
Segment new_top_segment(this->m_segments[0], depth_top, length_top);
this->addSegment(new_top_segment);
orderSegments();
// begin with the second segment
for (std::size_t i_index= 1; i_index< size(); ++i_index) {
if( m_segments[i_index].dataReady() ) continue;
// Begin with the second segment
for (std::size_t i_index = 1; i_index < size(); ++i_index) {
if (m_segments[i_index].dataReady()) {
continue;
}
// find its outlet segment
// Find its outlet segment
const int outlet_segment = m_segments[i_index].outletSegment();
const int outlet_index= segmentNumberToIndex(outlet_segment);
const int outlet_index = segmentNumberToIndex(outlet_segment);
// assert some information of the outlet_segment
assert(outlet_index>= 0);
assert(outlet_index >= 0);
assert(m_segments[outlet_index].dataReady());
const double outlet_depth = m_segments[outlet_index].depth();
const double outlet_length = m_segments[outlet_index].totalLength();
const double temp_depth = outlet_depth + m_segments[i_index].depth();
const double temp_length = outlet_length + m_segments[i_index].totalLength();
const auto new_x = m_segments[outlet_index].node_X() + m_segments[i_index].node_X();
const auto new_y = m_segments[outlet_index].node_Y() + m_segments[i_index].node_Y();
// applying the calculated length and depth to the current segment
Segment new_segment(this->m_segments[i_index], temp_depth, temp_length);
addSegment(new_segment);
// Applying the calculated length and depth to the current segment
this->addSegment(Segment {
this->m_segments[i_index],
temp_depth, temp_length, new_x, new_y
});
}
}
void WellSegments::orderSegments() {
// re-ordering the segments to make later use easier.
// two principles
// 1. the index of the outlet segment will be stored in the lower index than the segment.
// 2. the segments belong to the same branch will be continuously stored.
void WellSegments::orderSegments()
{
// Re-ordering the segments to make later use easier.
//
// Two principles
//
// 1. A segment's outlet segment is ordered before the segment
// itself.
//
// 2. Segments on the same branch are stored consecutively.
//
// Top segment always at index zero so we only reorder the segments
// in the index range [1 .. size()).
std::size_t current_index = 1;
// top segment will always be the first one
// before this index, the reordering is done.
std::size_t current_index= 1;
// clear the mapping from segment number to store index
// Clear the mapping from segment number to store index.
segment_number_to_index.clear();
// for the top segment
segment_number_to_index[1] = 0;
while (current_index< size()) {
// For the top segment
segment_number_to_index.insert_or_assign(1, 0);
while (current_index < size()) {
// the branch number of the last segment that is done re-ordering
const int last_branch_number = m_segments[current_index-1].branchNumber();
// the one need to be swapped to the current_index.

View File

@ -17,42 +17,38 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <iostream>
#include <memory>
#include <stdexcept>
#include <set>
#define BOOST_TEST_MODULE WellConnectionsTests
#include <boost/test/unit_test.hpp>
#include <string>
#include <opm/common/utility/OpmInputError.hpp>
#include <opm/input/eclipse/Deck/Deck.hpp>
#include <opm/input/eclipse/Deck/DeckItem.hpp>
#include <opm/input/eclipse/Deck/DeckRecord.hpp>
#include <opm/input/eclipse/Deck/DeckKeyword.hpp>
#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
#include <opm/input/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/input/eclipse/EclipseState/Grid/FieldPropsManager.hpp>
#include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
#include <opm/input/eclipse/Schedule/Schedule.hpp>
#include <opm/input/eclipse/Schedule/MSW/SICD.hpp>
#include <opm/input/eclipse/Schedule/MSW/Valve.hpp>
#include <opm/input/eclipse/Schedule/Well/Connection.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
#include "src/opm/input/eclipse/Schedule/MSW/Compsegs.hpp"
#include <opm/input/eclipse/Schedule/Schedule.hpp>
#include <opm/input/eclipse/Schedule/CompletedCells.hpp>
#include <opm/input/eclipse/Schedule/ScheduleGrid.hpp>
#include <opm/input/eclipse/Schedule/Well/Connection.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
#include <opm/input/eclipse/Parser/Parser.hpp>
#include <opm/input/eclipse/Parser/ErrorGuard.hpp>
#include <opm/input/eclipse/Parser/ParseContext.hpp>
#include <opm/input/eclipse/Deck/Deck.hpp>
#include <opm/input/eclipse/Deck/DeckItem.hpp>
#include <opm/input/eclipse/Deck/DeckRecord.hpp>
#include <opm/input/eclipse/Deck/DeckKeyword.hpp>
#include <memory>
#include <set>
#include <stdexcept>
#include <string>
BOOST_AUTO_TEST_CASE(AICDWellTest) {
@ -662,13 +658,15 @@ BOOST_AUTO_TEST_CASE(testwsegvalv) {
BOOST_CHECK_EQUAL(segment2.crossArea(), valv2.pipeCrossArea());
}
Opm::Schedule make_schedule(const std::string& fname) {
Opm::Parser parser;
Opm::Deck deck = parser.parseFile(fname);
Opm::EclipseState st(deck);
return Opm::Schedule(deck, st);
}
namespace {
Opm::Schedule make_schedule(const std::string& fname)
{
Opm::Parser parser;
Opm::Deck deck = parser.parseFile(fname);
Opm::EclipseState st(deck);
return Opm::Schedule(deck, st);
}
} // Anonymous namespace
BOOST_AUTO_TEST_CASE(MSW_SEGMENT_LENGTH) {
@ -737,3 +735,255 @@ BOOST_AUTO_TEST_CASE(MULTIPLE_WELSEGS) {
BOOST_CHECK(segments1 == segments2);
}
BOOST_AUTO_TEST_CASE(Node_XY_ABS_Individual)
{
const auto deck = ::Opm::Parser{}.parseString(R"(RUNSPEC
DIMENS
20 20 20 /
GRID
DXV
20*100 /
DYV
20*100 /
DZV
20*10 /
DEPTHZ
441*2000.0 /
PORO
8000*0.1 /
PERMX
8000*1 /
PERMY
8000*0.1 /
PERMZ
8000*0.01 /
SCHEDULE
WELSPECS
'PROD01' 'P' 20 20 1* OIL /
/
COMPDAT
'PROD01' 20 20 1 20 'OPEN' /
/
WELSEGS
'PROD01' 2512.5 2512.5 1.0e-5 'ABS' 'HF-' 'HO' 123.456 789.012 /
2 2 1 1 2537.5 2537.5 0.3 0.00010 2* 123.456 789.012 /
3 3 1 2 2562.5 2562.5 0.2 0.00010 2* 123.456 789.012 /
4 4 2 2 2737.5 2537.5 0.2 0.00010 2* 123.456 789.012 /
6 6 2 4 3037.5 2539.5 0.2 0.00010 2* 123.456 789.012 /
7 7 2 6 3337.5 2534.5 0.2 0.00010 2* 123.456 789.012 /
8 8 3 7 3337.6 2534.5 0.2 0.00015 2* 123.456 789.012 /
/
)");
const auto es = ::Opm::EclipseState { deck };
const auto sched = ::Opm::Schedule { deck, es, std::make_shared<const ::Opm::Python>() };
const auto& segments = sched[0].wells("PROD01").getSegments();
for (const auto& segment : segments) {
BOOST_CHECK_CLOSE(segment.node_X(), 123.456, 1.0e-8);
BOOST_CHECK_CLOSE(segment.node_Y(), 789.012, 1.0e-8);
}
}
BOOST_AUTO_TEST_CASE(Node_XY_ABS_Range)
{
const auto deck = ::Opm::Parser{}.parseString(R"(RUNSPEC
DIMENS
20 20 20 /
GRID
DXV
20*100 /
DYV
20*100 /
DZV
20*10 /
DEPTHZ
441*2000.0 /
PORO
8000*0.1 /
PERMX
8000*1 /
PERMY
8000*0.1 /
PERMZ
8000*0.01 /
SCHEDULE
WELSPECS
'PROD01' 'P' 20 20 1* OIL /
/
COMPDAT
'PROD01' 20 20 1 20 'OPEN' /
/
WELSEGS
'PROD01' 2512.5 2512.5 1.0e-5 'ABS' 'HF-' 'HO' 123.456 789.012 /
2 2 1 1 2537.5 2537.5 0.3 0.00010 2* 123.456 789.012 /
3 3 1 2 2562.5 2562.5 0.2 0.00010 2* 123.456 789.012 /
4 7 2 2 2737.5 2537.5 0.2 0.00010 2* 123.456 789.012 /
8 8 3 7 3337.6 2534.5 0.2 0.00015 2* 123.456 789.012 /
/
)");
const auto es = ::Opm::EclipseState { deck };
const auto sched = ::Opm::Schedule { deck, es, std::make_shared<const ::Opm::Python>() };
const auto& segments = sched[0].wells("PROD01").getSegments();
for (const auto& segment : segments) {
BOOST_CHECK_CLOSE(segment.node_X(), 123.456, 1.0e-8);
BOOST_CHECK_CLOSE(segment.node_Y(), 789.012, 1.0e-8);
}
}
BOOST_AUTO_TEST_CASE(Node_XY_INC_Individual)
{
const auto deck = ::Opm::Parser{}.parseString(R"(RUNSPEC
DIMENS
20 20 20 /
GRID
DXV
20*100 /
DYV
20*100 /
DZV
20*10 /
DEPTHZ
441*2000.0 /
PORO
8000*0.1 /
PERMX
8000*1 /
PERMY
8000*0.1 /
PERMZ
8000*0.01 /
SCHEDULE
WELSPECS
'PROD01' 'P' 20 20 1* OIL /
/
COMPDAT
'PROD01' 20 20 1 20 'OPEN' /
/
WELSEGS
-- Name Dep 1 Tlen 1 Vol 1 Len&Dep PresDrop
PROD01 2557.18408 0.00000 1* INC 'HF-' 'HO' 12.3 45.6 /
-- First Seg Last Seg Branch Num Outlet Seg Length Depth Change Diam Rough
-- Main Stem Segments
2 2 1 1 5.09434 4.95609 0.15200 0.0000100 2* 10.1 20.2 /
3 3 1 2 10.21718 9.93992 0.15200 0.0000100 2* 10.1 20.2 /
4 4 1 3 10.24573 9.96769 0.15200 0.0000100 2* 10.1 20.2 /
5 5 1 4 10.24574 9.96770 0.15200 0.0000100 2* 10.1 20.2 /
6 6 1 5 6.40355 6.22978 0.15200 0.0000100 2* 10.1 20.2 /
7 7 1 6 6.40355 6.22978 0.15200 0.0000100 2* 10.1 20.2 /
8 8 1 7 10.24567 9.96764 0.15200 0.0000100 2* 10.1 20.2 /
9 9 1 8 10.24571 9.96767 0.15200 0.0000100 2* 10.1 20.2 /
10 10 1 9 10.24570 9.96767 0.15200 0.0000100 2* 10.1 20.2 /
11 11 1 10 10.24571 9.96767 0.15200 0.0000100 2* 10.1 20.2 /
12 12 1 11 5.97902 5.81677 0.15200 0.0000100 2* 10.1 20.2 /
/
)");
const auto es = ::Opm::EclipseState { deck };
const auto sched = ::Opm::Schedule { deck, es, std::make_shared<const ::Opm::Python>() };
auto i = 0;
const auto& segments = sched[0].wells("PROD01").getSegments();
for (const auto& segment : segments) {
BOOST_CHECK_CLOSE(segment.node_X(), 12.3 + i*10.1, 1.0e-8);
BOOST_CHECK_CLOSE(segment.node_Y(), 45.6 + i*20.2, 1.0e-8);
++i;
}
}
BOOST_AUTO_TEST_CASE(Node_XY_INC_Range)
{
const auto deck = ::Opm::Parser{}.parseString(R"(RUNSPEC
DIMENS
20 20 20 /
GRID
DXV
20*100 /
DYV
20*100 /
DZV
20*10 /
DEPTHZ
441*2000.0 /
PORO
8000*0.1 /
PERMX
8000*1 /
PERMY
8000*0.1 /
PERMZ
8000*0.01 /
SCHEDULE
WELSPECS
'PROD01' 'P' 20 20 1* OIL /
/
COMPDAT
'PROD01' 20 20 1 20 'OPEN' /
/
WELSEGS
-- Name Dep 1 Tlen 1 Vol 1 Len&Dep PresDrop
PROD01 2557.18408 0.00000 1* INC 'HF-' 'HO' 12.3 45.6 /
-- First Seg Last Seg Branch Num Outlet Seg Length Depth Change Diam Rough
-- Main Stem Segments
2 12 1 1 5.09434 4.95609 0.15200 0.0000100 2* 10.1 20.2 /
/
)");
const auto es = ::Opm::EclipseState { deck };
const auto sched = ::Opm::Schedule { deck, es, std::make_shared<const ::Opm::Python>() };
auto i = 0;
const auto& segments = sched[0].wells("PROD01").getSegments();
for (const auto& segment : segments) {
BOOST_CHECK_CLOSE(segment.node_X(), 12.3 + i*10.1, 1.0e-8);
BOOST_CHECK_CLOSE(segment.node_Y(), 45.6 + i*20.2, 1.0e-8);
++i;
}
}