From 5785d088ea5ac2a708ce0ec5cbedb38ec00f6f8c Mon Sep 17 00:00:00 2001 From: Joakim Hove Date: Thu, 27 Mar 2014 15:41:19 +0100 Subject: [PATCH] Added functionality in Completion class to adjust I and J to wellHead value if they have been defaulted (i.e. have value -1). This is unfortunate because we loose quite a lot of const correctness. --- .../EclipseState/Schedule/Completion.cpp | 30 ++++++++++++------- .../EclipseState/Schedule/Completion.hpp | 5 ++-- .../EclipseState/Schedule/Schedule.cpp | 4 +-- .../eclipse/EclipseState/Schedule/Well.cpp | 8 +++-- .../eclipse/EclipseState/Schedule/Well.hpp | 4 +-- .../EclipseState/Schedule/tests/WellTests.cpp | 14 ++++----- .../IntegrationTests/CompletionsFromDeck.cpp | 10 +++---- .../ScheduleCreateFromDeck.cpp | 25 ++++++++++++++++ 8 files changed, 68 insertions(+), 32 deletions(-) diff --git a/opm/parser/eclipse/EclipseState/Schedule/Completion.cpp b/opm/parser/eclipse/EclipseState/Schedule/Completion.cpp index 39258c5c9..3dcbd9ea2 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/Completion.cpp +++ b/opm/parser/eclipse/EclipseState/Schedule/Completion.cpp @@ -50,8 +50,8 @@ namespace Opm { disentangled, and each completion is returned separately. */ - std::pair > Completion::completionsFromCOMPDATRecord( DeckRecordConstPtr compdatRecord ) { - std::vector completions; + std::pair > Completion::completionsFromCOMPDATRecord( DeckRecordConstPtr compdatRecord ) { + std::vector completions; std::string well = compdatRecord->getItem("WELL")->getString(0); // We change from eclipse's 1 - n, to a 0 - n-1 solution int I = compdatRecord->getItem("I")->getInt(0) - 1; @@ -66,16 +66,16 @@ namespace Opm { throw std::invalid_argument("The connection factor item can not be defaulted"); } double CF = compdatRecord->getItem("CF")->getSIDouble(0); - double diameter = compdatRecord->getItem("DIAMETER")->getSIDouble(0); double skinFactor = compdatRecord->getItem("SKIN")->getRawDouble(0); + for (int k = K1; k <= K2; k++) { - CompletionConstPtr completion(new Completion(I , J , k , state , CF, diameter, skinFactor )); + CompletionPtr completion(new Completion(I , J , k , state , CF, diameter, skinFactor )); completions.push_back( completion ); } - return std::pair >( well , completions ); + return std::pair >( well , completions ); } /* @@ -89,18 +89,18 @@ namespace Opm { */ - std::map > Completion::completionsFromCOMPDATKeyword( DeckKeywordConstPtr compdatKeyword ) { - std::map > completionMapList; + std::map > Completion::completionsFromCOMPDATKeyword( DeckKeywordConstPtr compdatKeyword ) { + std::map > completionMapList; for (size_t recordIndex = 0; recordIndex < compdatKeyword->size(); recordIndex++) { - std::pair > wellCompletionsPair = completionsFromCOMPDATRecord( compdatKeyword->getRecord( recordIndex )); + std::pair > wellCompletionsPair = completionsFromCOMPDATRecord( compdatKeyword->getRecord( recordIndex )); std::string well = wellCompletionsPair.first; - std::vector& newCompletions = wellCompletionsPair.second; + std::vector& newCompletions = wellCompletionsPair.second; if (completionMapList.find(well) == completionMapList.end()) - completionMapList[well] = std::vector(); + completionMapList[well] = std::vector(); { - std::vector& currentCompletions = completionMapList.find(well)->second; + std::vector& currentCompletions = completionMapList.find(well)->second; for (size_t ic = 0; ic < newCompletions.size(); ic++) currentCompletions.push_back( newCompletions[ic] ); @@ -109,6 +109,14 @@ namespace Opm { return completionMapList; } + void Completion::fixDefaultIJ(int wellHeadI , int wellHeadJ) { + if (m_i < 0) + m_i = wellHeadI; + + if (m_j < 0) + m_j = wellHeadJ; + } + int Completion::getI() const { return m_i; diff --git a/opm/parser/eclipse/EclipseState/Schedule/Completion.hpp b/opm/parser/eclipse/EclipseState/Schedule/Completion.hpp index 5b25c57a3..ddb271e91 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/Completion.hpp +++ b/opm/parser/eclipse/EclipseState/Schedule/Completion.hpp @@ -42,9 +42,10 @@ namespace Opm { double getCF() const; double getDiameter() const; double getSkinFactor() const; + void fixDefaultIJ(int wellHeadI , int wellHeadJ); - static std::map > > completionsFromCOMPDATKeyword( DeckKeywordConstPtr compdatKeyword ); - static std::pair > > completionsFromCOMPDATRecord( DeckRecordConstPtr compdatRecord ); + static std::map > > completionsFromCOMPDATKeyword( DeckKeywordConstPtr compdatKeyword ); + static std::pair > > completionsFromCOMPDATRecord( DeckRecordConstPtr compdatRecord ); private: int m_i, m_j, m_k; diff --git a/opm/parser/eclipse/EclipseState/Schedule/Schedule.cpp b/opm/parser/eclipse/EclipseState/Schedule/Schedule.cpp index 8c6fcb1aa..03aba60a1 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/Schedule.cpp +++ b/opm/parser/eclipse/EclipseState/Schedule/Schedule.cpp @@ -439,8 +439,8 @@ namespace Opm { } void Schedule::handleCOMPDAT(DeckKeywordConstPtr keyword , size_t currentStep) { - std::map > completionMapList = Completion::completionsFromCOMPDATKeyword( keyword ); - std::map >::iterator iter; + std::map > completionMapList = Completion::completionsFromCOMPDATKeyword( keyword ); + std::map >::iterator iter; for( iter= completionMapList.begin(); iter != completionMapList.end(); iter++) { const std::string wellName = iter->first; diff --git a/opm/parser/eclipse/EclipseState/Schedule/Well.cpp b/opm/parser/eclipse/EclipseState/Schedule/Well.cpp index 2296e5359..81dc04349 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/Well.cpp +++ b/opm/parser/eclipse/EclipseState/Schedule/Well.cpp @@ -153,13 +153,15 @@ namespace Opm { return m_completions->get( timeStep ); } - void Well::addCompletions(size_t time_step , const std::vector& newCompletions) { + void Well::addCompletions(size_t time_step , const std::vector& newCompletions) { CompletionSetConstPtr currentCompletionSet = m_completions->get(time_step); CompletionSetPtr newCompletionSet = CompletionSetPtr( currentCompletionSet->shallowCopy() ); - for (size_t ic = 0; ic < newCompletions.size(); ic++) + for (size_t ic = 0; ic < newCompletions.size(); ic++) { + newCompletions[ic]->fixDefaultIJ( m_headI , m_headJ ); newCompletionSet->add( newCompletions[ic] ); - + } + m_completions->add( time_step , newCompletionSet); } diff --git a/opm/parser/eclipse/EclipseState/Schedule/Well.hpp b/opm/parser/eclipse/EclipseState/Schedule/Well.hpp index b7e336838..bc465d7b7 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/Well.hpp +++ b/opm/parser/eclipse/EclipseState/Schedule/Well.hpp @@ -150,8 +150,8 @@ namespace Opm { bool isProducer(size_t timeStep) const; bool isInjector(size_t timeStep) const; void addWELSPECS(DeckRecordConstPtr deckRecord); - void addCompletions(size_t time_step , const std::vector& newCompletions); - CompletionSetConstPtr getCompletions(size_t timeStep) const; + void addCompletions(size_t time_step , const std::vector& newCompletions); + CompletionSetConstPtr getCompletions(size_t timeStep) const; void setProductionProperties(size_t timeStep , const WellProductionProperties properties); WellProductionProperties getProductionPropertiesCopy(size_t timeStep) const; diff --git a/opm/parser/eclipse/EclipseState/Schedule/tests/WellTests.cpp b/opm/parser/eclipse/EclipseState/Schedule/tests/WellTests.cpp index da517877c..ab7cd47b0 100644 --- a/opm/parser/eclipse/EclipseState/Schedule/tests/WellTests.cpp +++ b/opm/parser/eclipse/EclipseState/Schedule/tests/WellTests.cpp @@ -168,13 +168,13 @@ BOOST_AUTO_TEST_CASE(UpdateCompletions) { Opm::CompletionSetConstPtr completions = well.getCompletions( 0 ); BOOST_CHECK_EQUAL( 0U , completions->size()); - std::vector newCompletions; - std::vector newCompletions2; - Opm::CompletionConstPtr comp1(new Opm::Completion( 10 , 10 , 10 , Opm::AUTO , 99.0, 22.3, 33.2)); - Opm::CompletionConstPtr comp2(new Opm::Completion( 10 , 11 , 10 , Opm::SHUT , 99.0, 22.3, 33.2)); - Opm::CompletionConstPtr comp3(new Opm::Completion( 10 , 10 , 12 , Opm::OPEN , 99.0, 22.3, 33.2)); - Opm::CompletionConstPtr comp4(new Opm::Completion( 10 , 10 , 12 , Opm::SHUT , 99.0, 22.3, 33.2)); - Opm::CompletionConstPtr comp5(new Opm::Completion( 10 , 10 , 13 , Opm::OPEN , 99.0, 22.3, 33.2)); + std::vector newCompletions; + std::vector newCompletions2; + Opm::CompletionPtr comp1(new Opm::Completion( 10 , 10 , 10 , Opm::AUTO , 99.0, 22.3, 33.2)); + Opm::CompletionPtr comp2(new Opm::Completion( 10 , 11 , 10 , Opm::SHUT , 99.0, 22.3, 33.2)); + Opm::CompletionPtr comp3(new Opm::Completion( 10 , 10 , 12 , Opm::OPEN , 99.0, 22.3, 33.2)); + Opm::CompletionPtr comp4(new Opm::Completion( 10 , 10 , 12 , Opm::SHUT , 99.0, 22.3, 33.2)); + Opm::CompletionPtr comp5(new Opm::Completion( 10 , 10 , 13 , Opm::OPEN , 99.0, 22.3, 33.2)); //std::vector newCompletions2{ comp4 , comp5}; Newer c++ diff --git a/opm/parser/eclipse/IntegrationTests/CompletionsFromDeck.cpp b/opm/parser/eclipse/IntegrationTests/CompletionsFromDeck.cpp index 0a0c159fe..f9d056ff4 100644 --- a/opm/parser/eclipse/IntegrationTests/CompletionsFromDeck.cpp +++ b/opm/parser/eclipse/IntegrationTests/CompletionsFromDeck.cpp @@ -43,13 +43,13 @@ BOOST_AUTO_TEST_CASE( CreateCompletionsFromRecord ) { DeckRecordConstPtr line0 = COMPDAT1->getRecord(0); DeckRecordConstPtr line1 = COMPDAT1->getRecord(1); - std::pair< std::string , std::vector > completionsList = Completion::completionsFromCOMPDATRecord( line0 ); + std::pair< std::string , std::vector > completionsList = Completion::completionsFromCOMPDATRecord( line0 ); BOOST_CHECK_EQUAL( "W_1" , completionsList.first ); BOOST_CHECK_EQUAL( 3U , completionsList.second.size() ); { - CompletionConstPtr completion0 = completionsList.second[0]; - CompletionConstPtr completion2 = completionsList.second[2]; + CompletionPtr completion0 = completionsList.second[0]; + CompletionPtr completion2 = completionsList.second[2]; BOOST_CHECK_EQUAL( 29 , completion0->getI() ); BOOST_CHECK_EQUAL( 36 , completion0->getJ() ); @@ -77,7 +77,7 @@ BOOST_AUTO_TEST_CASE( CreateCompletionsFromKeyword ) { DeckPtr deck = parser->parseFile(scheduleFile.string()); DeckKeywordConstPtr COMPDAT1 = deck->getKeyword("COMPDAT" , 1); - std::map< std::string , std::vector > completions = Completion::completionsFromCOMPDATKeyword( COMPDAT1 ); + std::map< std::string , std::vector > completions = Completion::completionsFromCOMPDATKeyword( COMPDAT1 ); BOOST_CHECK_EQUAL( 3U , completions.size() ); @@ -89,7 +89,7 @@ BOOST_AUTO_TEST_CASE( CreateCompletionsFromKeyword ) { BOOST_CHECK_EQUAL( 5U , completions.find("W_2")->second.size() ); BOOST_CHECK_EQUAL( 5U , completions.find("W_3")->second.size() ); - std::vector W_3Completions = completions.find("W_3")->second; + std::vector W_3Completions = completions.find("W_3")->second; CompletionConstPtr completion0 = W_3Completions[0]; CompletionConstPtr completion4 = W_3Completions[4]; diff --git a/opm/parser/eclipse/IntegrationTests/ScheduleCreateFromDeck.cpp b/opm/parser/eclipse/IntegrationTests/ScheduleCreateFromDeck.cpp index a854319c3..ebb69b9f9 100644 --- a/opm/parser/eclipse/IntegrationTests/ScheduleCreateFromDeck.cpp +++ b/opm/parser/eclipse/IntegrationTests/ScheduleCreateFromDeck.cpp @@ -537,3 +537,28 @@ BOOST_AUTO_TEST_CASE(WellTestWGRUPCONWellPropertiesSet) { BOOST_CHECK_EQUAL(0.5, well3->getGuideRateScalingFactor(0)); } + +BOOST_AUTO_TEST_CASE(TestDefaultedCOMPDATIJ) { + ParserPtr parser(new Parser()); + boost::filesystem::path scheduleFile("testdata/integration_tests/SCHEDULE/SCHEDULE_COMPDAT_DEFAULT_IJ"); + const char * deckString = "\n\ +START\n\ +\n\ +10 MAI 2007 /\n\ +\n\ +SCHEDULE\n\ +WELSPECS \n\ + 'W1' 'OP' 11 21 3.33 'OIL' 7* / \n\ +/\n\ +COMPDAT \n\ + 'W1' 2* 1 1 'OPEN' 1* 32.948 0.311 3047.839 2* 'X' 22.100 /\n\ +/\n"; + DeckPtr deck = parser->parseString(deckString); + ScheduleConstPtr sched(new Schedule(deck)); + WellConstPtr well = sched->getWell("W1"); + CompletionSetConstPtr completions = well->getCompletions(0); + BOOST_CHECK_EQUAL( 10 , completions->get(0)->getI() ); + BOOST_CHECK_EQUAL( 20 , completions->get(0)->getJ() ); +} + +