Merge pull request #720 from GitPaean/throw_for_problem_in_compsegs
throw when end distance not bigger than start distance in COMPSEGS
This commit is contained in:
@@ -27,8 +27,9 @@
|
||||
#include <opm/parser/eclipse/EclipseState/Schedule/MSW/WellSegments.hpp>
|
||||
|
||||
namespace Opm {
|
||||
WellConnections * newConnectionsWithSegments(const DeckKeyword& compsegs, const WellConnections& input_connections,
|
||||
const WellSegments& segments, const EclipseGrid& grid);
|
||||
WellConnections * newConnectionsWithSegments(const DeckKeyword& compsegs, const WellConnections& input_connections,
|
||||
const WellSegments& segments, const EclipseGrid& grid,
|
||||
const ParseContext& parseContext, ErrorGuard& errors);
|
||||
}
|
||||
|
||||
#endif
|
||||
|
||||
@@ -207,7 +207,7 @@ namespace Opm
|
||||
void handleCOMPDAT( const DeckKeyword& keyword, size_t currentStep, const EclipseGrid& grid, const Eclipse3DProperties& eclipseProperties, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
void handleCOMPLUMP( const DeckKeyword& keyword, size_t currentStep );
|
||||
void handleWELSEGS( const DeckKeyword& keyword, size_t currentStep);
|
||||
void handleCOMPSEGS( const DeckKeyword& keyword, size_t currentStep, const EclipseGrid& grid);
|
||||
void handleCOMPSEGS( const DeckKeyword& keyword, size_t currentStep, const EclipseGrid& grid, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
void handleWCONINJE( const SCHEDULESection&, const DeckKeyword& keyword, size_t currentStep, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
void handleWPOLYMER( const DeckKeyword& keyword, size_t currentStep, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
void handleWSOLVENT( const DeckKeyword& keyword, size_t currentStep, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
|
||||
@@ -142,7 +142,7 @@ public:
|
||||
bool updateInjection(std::shared_ptr<WellInjectionProperties> injection);
|
||||
|
||||
bool handleWELSEGS(const DeckKeyword& keyword);
|
||||
bool handleCOMPSEGS(const DeckKeyword& keyword, const EclipseGrid& grid);
|
||||
bool handleCOMPSEGS(const DeckKeyword& keyword, const EclipseGrid& grid, const ParseContext& parseContext, ErrorGuard& errors);
|
||||
bool handleWELOPEN(const DeckRecord& record, WellCompletion::StateEnum status);
|
||||
bool handleCOMPLUMP(const DeckRecord& record);
|
||||
bool handleWPIMULT(const DeckRecord& record);
|
||||
|
||||
@@ -333,6 +333,9 @@ namespace Opm {
|
||||
|
||||
const static std::string SCHEDULE_WELL_ERROR;
|
||||
|
||||
const static std::string SCHEDULE_COMPSEGS_INVALID;
|
||||
const static std::string SCHEDULE_COMPSEGS_NOT_SUPPORTED;
|
||||
|
||||
/*
|
||||
The SIMULATOR_KEYWORD_ errormodes are for the situation where the
|
||||
parser recognizes, and correctly parses a keyword, but we know that
|
||||
|
||||
@@ -19,6 +19,7 @@
|
||||
|
||||
#include <cmath>
|
||||
#include <iostream>
|
||||
#include <sstream>
|
||||
|
||||
#include <opm/parser/eclipse/Parser/ParserKeywords/C.hpp>
|
||||
|
||||
@@ -50,12 +51,13 @@ namespace Opm {
|
||||
{
|
||||
}
|
||||
|
||||
std::vector< Compsegs > Compsegs::compsegsFromCOMPSEGSKeyword( const DeckKeyword& compsegsKeyword, const EclipseGrid& grid) {
|
||||
std::vector< Compsegs > Compsegs::compsegsFromCOMPSEGSKeyword(const DeckKeyword& compsegsKeyword, const EclipseGrid& grid,
|
||||
const ParseContext& parseContext, ErrorGuard& errors) {
|
||||
|
||||
// only handle the second record here
|
||||
// The first record here only contains the well name
|
||||
std::vector< Compsegs > compsegs;
|
||||
|
||||
// The first record in the keyword only contains the well name
|
||||
// looping from the second record in the keyword
|
||||
for (size_t recordIndex = 1; recordIndex < compsegsKeyword.size(); ++recordIndex) {
|
||||
const auto& record = compsegsKeyword.getRecord(recordIndex);
|
||||
// following the coordinate rule for connections
|
||||
@@ -64,6 +66,8 @@ namespace Opm {
|
||||
const int K = record.getItem<ParserKeywords::COMPSEGS::K>().get< int >(0) - 1;
|
||||
const int branch = record.getItem<ParserKeywords::COMPSEGS::BRANCH>().get< int >(0);
|
||||
|
||||
const std::string& well_name = compsegsKeyword.getRecord(0).getItem("WELL").getTrimmedString(0);
|
||||
|
||||
double distance_start;
|
||||
double distance_end;
|
||||
if (record.getItem<ParserKeywords::COMPSEGS::DISTANCE_START>().hasValue(0)) {
|
||||
@@ -74,23 +78,38 @@ namespace Opm {
|
||||
// TODO: the end of the previous connection or range
|
||||
// 'previous' should be in term of the input order
|
||||
// since basically no specific order for the connections
|
||||
throw std::runtime_error("this way to obtain DISTANCE_START not implemented yet!");
|
||||
const std::string msg = "This way to obtain DISTANCE_START in keyword COMPSEGS "
|
||||
"is not implemented yet for well " + well_name;
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, msg, errors);
|
||||
}
|
||||
if (record.getItem<ParserKeywords::COMPSEGS::DISTANCE_END>().hasValue(0)) {
|
||||
distance_end = record.getItem<ParserKeywords::COMPSEGS::DISTANCE_END>().getSIDouble(0);
|
||||
} else {
|
||||
// TODO: the distance_start plus the thickness of the grid block
|
||||
throw std::runtime_error("this way to obtain DISTANCE_END not implemented yet!");
|
||||
const std::string msg = "This way to obtain DISTANCE_END in keyword COMPSEGS "
|
||||
"is not implemented yet for well " + well_name;
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, msg, errors);
|
||||
}
|
||||
|
||||
if (distance_end <= distance_start) {
|
||||
std::ostringstream sstr;
|
||||
sstr << " The end of the perforations need be to further down than the start of the perforations\n "
|
||||
<< " well " << well_name << " " << I + 1 << " " << J + 1 << " " << K + 1 << " in keyword COMPSEGS\n";
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_INVALID, sstr.str(), errors);
|
||||
}
|
||||
|
||||
if( !record.getItem< ParserKeywords::COMPSEGS::DIRECTION >().hasValue( 0 ) &&
|
||||
!record.getItem< ParserKeywords::COMPSEGS::DISTANCE_END >().hasValue( 0 ) ) {
|
||||
throw std::runtime_error("the direction has to be specified when DISTANCE_END in the record is not specified");
|
||||
const std::string msg = "The direction has to be specified when DISTANCE_END "
|
||||
"is not specified in keyword COMPSEGS for well " + well_name;
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_INVALID, msg, errors);
|
||||
}
|
||||
|
||||
if( record.getItem< ParserKeywords::COMPSEGS::END_IJK >().hasValue( 0 ) &&
|
||||
!record.getItem< ParserKeywords::COMPSEGS::DIRECTION >().hasValue( 0 ) ) {
|
||||
throw std::runtime_error("the direction has to be specified when END_IJK in the record is specified");
|
||||
const std::string msg = "The direction has to be specified when END_IJK "
|
||||
"is specified in keyword COMPSEGS for well " + well_name;
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_INVALID, msg, errors);
|
||||
}
|
||||
|
||||
/*
|
||||
@@ -113,7 +132,9 @@ namespace Opm {
|
||||
|
||||
if (center_depth < 0.) {
|
||||
//TODO: get the depth from COMPDAT data.
|
||||
throw std::runtime_error("this way to obtain CENTER_DISTANCE not implemented yet either!");
|
||||
const std::string msg = "This way to obtain CENTER_DISTANCE in keyword COMPSEGS "
|
||||
"is not implemented yet for well " + well_name;
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, msg, errors);
|
||||
}
|
||||
|
||||
int segment_number;
|
||||
@@ -136,7 +157,10 @@ namespace Opm {
|
||||
seqIndex);
|
||||
}
|
||||
} else { // a range is defined. genrate a range of Compsegs
|
||||
throw std::runtime_error("entering COMPSEGS entries with a range is not supported yet!");
|
||||
std::ostringstream sstr;
|
||||
sstr << "COMPSEGS entries can only be input for single connection, not supporting COMPSEGS entries specified with a range yet.\n"
|
||||
<< " well " << well_name << " " << I + 1 << " " << J + 1 << " " << K + 1 << " in keyword COMPSEGS\n";
|
||||
parseContext.handleError(ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, sstr.str(), errors);
|
||||
}
|
||||
}
|
||||
|
||||
@@ -169,11 +193,10 @@ namespace Opm {
|
||||
}
|
||||
|
||||
if (segment_number == 0) {
|
||||
throw std::runtime_error("The perforation failed in finding a related segment \n");
|
||||
}
|
||||
|
||||
if (compseg.center_depth < 0.) {
|
||||
throw std::runtime_error("Obtaining perforation depth from COMPDAT data is not supported yet");
|
||||
std::ostringstream sstr;
|
||||
sstr << "The connection specified in COMPSEGS with index of " << compseg.m_i + 1 << " "
|
||||
<< compseg.m_j + 1 << " " << compseg.m_k + 1 << " failed in finding a related segment";
|
||||
throw std::runtime_error(sstr.str());
|
||||
}
|
||||
|
||||
compseg.segment_number = segment_number;
|
||||
|
||||
@@ -58,14 +58,15 @@ namespace Opm {
|
||||
|
||||
void calculateCenterDepthWithSegments(const WellSegments& segment_set);
|
||||
|
||||
static std::vector< Compsegs > compsegsFromCOMPSEGSKeyword( const DeckKeyword& compsegsKeyword, const EclipseGrid& grid);
|
||||
static std::vector< Compsegs > compsegsFromCOMPSEGSKeyword(const DeckKeyword& compsegsKeyword, const EclipseGrid& grid,
|
||||
const ParseContext& parseContext, ErrorGuard& errors);
|
||||
|
||||
// get the segment number information and depth information based on the information from WellSegments
|
||||
static void processCOMPSEGS(std::vector< Compsegs >& compsegs, const WellSegments& segment_set );
|
||||
static void processCOMPSEGS(std::vector< Compsegs >& compsegs, const WellSegments& segment_set);
|
||||
|
||||
// update the segment related information for Connections
|
||||
static void updateConnectionsWithSegment(const std::vector< Compsegs >& compsegs,
|
||||
const EclipseGrid& grid,
|
||||
const EclipseGrid& grid,
|
||||
WellConnections& connection_set);
|
||||
|
||||
};
|
||||
|
||||
@@ -27,10 +27,12 @@ namespace Opm {
|
||||
WellConnections * newConnectionsWithSegments(const DeckKeyword& compsegs,
|
||||
const WellConnections& input_connections,
|
||||
const WellSegments& segment_set,
|
||||
const EclipseGrid& grid)
|
||||
const EclipseGrid& grid,
|
||||
const ParseContext& parseContext,
|
||||
ErrorGuard& errors)
|
||||
{
|
||||
WellConnections * new_connection_set = new WellConnections(input_connections);
|
||||
std::vector<Compsegs> compsegs_vector = Compsegs::compsegsFromCOMPSEGSKeyword( compsegs, grid );
|
||||
std::vector<Compsegs> compsegs_vector = Compsegs::compsegsFromCOMPSEGSKeyword( compsegs, grid, parseContext, errors);
|
||||
Compsegs::processCOMPSEGS(compsegs_vector, segment_set);
|
||||
Compsegs::updateConnectionsWithSegment(compsegs_vector, grid, *new_connection_set);
|
||||
return new_connection_set;
|
||||
|
||||
@@ -275,7 +275,7 @@ namespace Opm {
|
||||
handleWELSEGS(keyword, currentStep);
|
||||
|
||||
else if (keyword.name() == "COMPSEGS")
|
||||
handleCOMPSEGS(keyword, currentStep, grid);
|
||||
handleCOMPSEGS(keyword, currentStep, grid, parseContext, errors);
|
||||
|
||||
else if (keyword.name() == "WELOPEN")
|
||||
handleWELOPEN(keyword, currentStep, parseContext, errors);
|
||||
@@ -1669,13 +1669,14 @@ namespace Opm {
|
||||
}
|
||||
}
|
||||
|
||||
void Schedule::handleCOMPSEGS( const DeckKeyword& keyword, size_t currentStep, const EclipseGrid& grid) {
|
||||
void Schedule::handleCOMPSEGS( const DeckKeyword& keyword, size_t currentStep, const EclipseGrid& grid,
|
||||
const ParseContext& parseContext, ErrorGuard& errors) {
|
||||
const auto& record1 = keyword.getRecord(0);
|
||||
const std::string& well_name = record1.getItem("WELL").getTrimmedString(0);
|
||||
{
|
||||
auto& dynamic_state = this->wells_static.at(well_name);
|
||||
auto well_ptr = std::make_shared<Well2>( *dynamic_state[currentStep] );
|
||||
if (well_ptr->handleCOMPSEGS(keyword, grid))
|
||||
if (well_ptr->handleCOMPSEGS(keyword, grid, parseContext, errors))
|
||||
this->updateWell(well_ptr, currentStep);
|
||||
}
|
||||
}
|
||||
|
||||
@@ -334,8 +334,10 @@ bool Well2::updateSolventFraction(double solvent_fraction) {
|
||||
}
|
||||
|
||||
|
||||
bool Well2::handleCOMPSEGS(const DeckKeyword& keyword, const EclipseGrid& grid) {
|
||||
std::shared_ptr<WellConnections> new_connection_set( newConnectionsWithSegments(keyword, *this->connections, *this->segments , grid) );
|
||||
bool Well2::handleCOMPSEGS(const DeckKeyword& keyword, const EclipseGrid& grid,
|
||||
const ParseContext& parseContext, ErrorGuard& errors) {
|
||||
std::shared_ptr<WellConnections> new_connection_set( newConnectionsWithSegments(keyword, *this->connections, *this->segments , grid,
|
||||
parseContext, errors) );
|
||||
return this->updateConnections(new_connection_set);
|
||||
}
|
||||
|
||||
|
||||
@@ -110,6 +110,8 @@ namespace Opm {
|
||||
|
||||
addKey(UDQ_PARSE_ERROR, InputError::THROW_EXCEPTION);
|
||||
this->addKey(SCHEDULE_WELL_ERROR, InputError::THROW_EXCEPTION);
|
||||
this->addKey(SCHEDULE_COMPSEGS_INVALID, InputError::THROW_EXCEPTION);
|
||||
this->addKey(SCHEDULE_COMPSEGS_NOT_SUPPORTED, InputError::THROW_EXCEPTION);
|
||||
}
|
||||
|
||||
void ParseContext::initEnv() {
|
||||
@@ -344,4 +346,7 @@ namespace Opm {
|
||||
|
||||
const std::string ParseContext::UDQ_PARSE_ERROR = "UDQ_PARSE_ERROR";
|
||||
const std::string ParseContext::SCHEDULE_WELL_ERROR = "SCHEDULE_WELL_ERROR";
|
||||
|
||||
const std::string ParseContext::SCHEDULE_COMPSEGS_INVALID = "SCHEDULE_COMPSEG_INVALID";
|
||||
const std::string ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED = "SCHEDULE_COMPSEGS_NOT_SUPPORTED";
|
||||
}
|
||||
|
||||
@@ -87,7 +87,13 @@ BOOST_AUTO_TEST_CASE(MultisegmentWellTest) {
|
||||
segment_set.loadWELSEGS(welsegs);
|
||||
|
||||
BOOST_CHECK_EQUAL(6U, segment_set.size());
|
||||
const Opm::WellConnections * new_connection_set = Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid);
|
||||
|
||||
Opm::ErrorGuard errorGuard;
|
||||
Opm::ParseContext parseContext;
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_INVALID, Opm::InputError::THROW_EXCEPTION);
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, Opm::InputError::THROW_EXCEPTION);
|
||||
Opm::WellConnections* new_connection_set = nullptr;
|
||||
BOOST_CHECK_NO_THROW(new_connection_set = Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid, parseContext, errorGuard));
|
||||
|
||||
BOOST_CHECK_EQUAL(7U, new_connection_set->size());
|
||||
|
||||
@@ -121,3 +127,117 @@ BOOST_AUTO_TEST_CASE(MultisegmentWellTest) {
|
||||
BOOST_CHECK_EQUAL(segment_number_connection7, 7);
|
||||
BOOST_CHECK_EQUAL(center_depth_connection7, 2534.5);
|
||||
}
|
||||
|
||||
BOOST_AUTO_TEST_CASE(WrongDistanceCOMPSEGS) {
|
||||
Opm::WellCompletion::DirectionEnum dir = Opm::WellCompletion::DirectionEnum::Z;
|
||||
Opm::WellConnections connection_set(10,10);
|
||||
Opm::EclipseGrid grid(20,20,20);
|
||||
connection_set.add(Opm::Connection( 19, 0, 0, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 19, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 19, 0, 2, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
|
||||
connection_set.add(Opm::Connection( 18, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 17, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 16, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 15, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
|
||||
BOOST_CHECK_EQUAL( 7U , connection_set.size() );
|
||||
|
||||
const std::string compsegs_string =
|
||||
"WELSEGS \n"
|
||||
"'PROD01' 2512.5 2512.5 1.0e-5 'ABS' 'H--' 'HO' /\n"
|
||||
"2 2 1 1 2537.5 2537.5 0.3 0.00010 /\n"
|
||||
"3 3 1 2 2562.5 2562.5 0.2 0.00010 /\n"
|
||||
"4 4 2 2 2737.5 2537.5 0.2 0.00010 /\n"
|
||||
"6 6 2 4 3037.5 2539.5 0.2 0.00010 /\n"
|
||||
"7 7 2 6 3337.5 2534.5 0.2 0.00010 /\n"
|
||||
"/\n"
|
||||
"\n"
|
||||
"COMPSEGS\n"
|
||||
"PROD01 / \n"
|
||||
"20 1 1 1 2512.5 2525.0 /\n"
|
||||
"20 1 2 1 2525.0 2550.0 /\n"
|
||||
"20 1 3 1 2550.0 2545.0 /\n"
|
||||
"19 1 2 2 2637.5 2837.5 /\n"
|
||||
"18 1 2 2 2837.5 3037.5 /\n"
|
||||
"17 1 2 2 3037.5 3237.5 /\n"
|
||||
"16 1 2 2 3237.5 3437.5 /\n"
|
||||
"/\n";
|
||||
|
||||
Opm::Parser parser;
|
||||
Opm::Deck deck = parser.parseString(compsegs_string);
|
||||
|
||||
const Opm::DeckKeyword compsegs = deck.getKeyword("COMPSEGS");
|
||||
BOOST_CHECK_EQUAL( 8U, compsegs.size() );
|
||||
|
||||
Opm::WellSegments segment_set;
|
||||
const Opm::DeckKeyword welsegs = deck.getKeyword("WELSEGS");
|
||||
segment_set.loadWELSEGS(welsegs);
|
||||
|
||||
BOOST_CHECK_EQUAL(6U, segment_set.size());
|
||||
|
||||
Opm::ErrorGuard errorGuard;
|
||||
Opm::ParseContext parseContext;
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_INVALID, Opm::InputError::THROW_EXCEPTION);
|
||||
BOOST_CHECK_THROW(Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid, parseContext, errorGuard), std::invalid_argument);
|
||||
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_INVALID, Opm::InputError::IGNORE);
|
||||
BOOST_CHECK_NO_THROW(Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid, parseContext, errorGuard));
|
||||
}
|
||||
|
||||
BOOST_AUTO_TEST_CASE(NegativeDepthCOMPSEGS) {
|
||||
Opm::WellCompletion::DirectionEnum dir = Opm::WellCompletion::DirectionEnum::Z;
|
||||
Opm::WellConnections connection_set(10,10);
|
||||
Opm::EclipseGrid grid(20,20,20);
|
||||
connection_set.add(Opm::Connection( 19, 0, 0, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 19, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 19, 0, 2, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, dir,0, 0., 0., true) );
|
||||
|
||||
connection_set.add(Opm::Connection( 18, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 17, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 16, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
connection_set.add(Opm::Connection( 15, 0, 1, 1, 0.0, Opm::WellCompletion::OPEN , 200, 17.29, 0.25, 0.0, 0.0, 0, Opm::WellCompletion::DirectionEnum::X,0, 0., 0., true) );
|
||||
|
||||
BOOST_CHECK_EQUAL( 7U , connection_set.size() );
|
||||
|
||||
const std::string compsegs_string =
|
||||
"WELSEGS \n"
|
||||
"'PROD01' 2512.5 2512.5 1.0e-5 'ABS' 'H--' 'HO' /\n"
|
||||
"2 2 1 1 2537.5 2537.5 0.3 0.00010 /\n"
|
||||
"3 3 1 2 2562.5 2562.5 0.2 0.00010 /\n"
|
||||
"4 4 2 2 2737.5 2537.5 0.2 0.00010 /\n"
|
||||
"6 6 2 4 3037.5 2539.5 0.2 0.00010 /\n"
|
||||
"7 7 2 6 3337.5 2534.5 0.2 0.00010 /\n"
|
||||
"/\n"
|
||||
"\n"
|
||||
"COMPSEGS\n"
|
||||
"PROD01 / \n"
|
||||
"20 1 1 1 2512.5 2525.0 /\n"
|
||||
"20 1 2 1 2525.0 2550.0 /\n"
|
||||
"20 1 3 1 2550.0 2575.0 /\n"
|
||||
"19 1 2 2 2637.5 2837.5 2* -8./\n"
|
||||
"18 1 2 2 2837.5 3037.5 /\n"
|
||||
"17 1 2 2 3037.5 3237.5 /\n"
|
||||
"16 1 2 2 3237.5 3437.5 /\n"
|
||||
"/\n";
|
||||
|
||||
Opm::Parser parser;
|
||||
Opm::Deck deck = parser.parseString(compsegs_string);
|
||||
|
||||
const Opm::DeckKeyword compsegs = deck.getKeyword("COMPSEGS");
|
||||
BOOST_CHECK_EQUAL( 8U, compsegs.size() );
|
||||
|
||||
Opm::WellSegments segment_set;
|
||||
const Opm::DeckKeyword welsegs = deck.getKeyword("WELSEGS");
|
||||
segment_set.loadWELSEGS(welsegs);
|
||||
|
||||
BOOST_CHECK_EQUAL(6U, segment_set.size());
|
||||
|
||||
Opm::ErrorGuard errorGuard;
|
||||
Opm::ParseContext parseContext;
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, Opm::InputError::THROW_EXCEPTION);
|
||||
BOOST_CHECK_THROW(Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid, parseContext, errorGuard), std::invalid_argument);
|
||||
|
||||
parseContext.update(Opm::ParseContext::SCHEDULE_COMPSEGS_NOT_SUPPORTED, Opm::InputError::IGNORE);
|
||||
BOOST_CHECK_NO_THROW(Opm::newConnectionsWithSegments(compsegs, connection_set, segment_set, grid, parseContext, errorGuard));
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user