diff --git a/opm/io/eclipse/ERft.hpp b/opm/io/eclipse/ERft.hpp index 14b4ebb57..4c6f2c531 100644 --- a/opm/io/eclipse/ERft.hpp +++ b/opm/io/eclipse/ERft.hpp @@ -44,16 +44,20 @@ public: template const std::vector& getRft(const std::string& name, const std::string& wellName, int year, int month, int day) const; + template + const std::vector& getRft(const std::string& name, int reportIndex) const; std::vector listOfWells() const; std::vector listOfdates() const; - using RftReportList = std::vector>; + using RftReportList = std::vector>; const RftReportList& listOfRftReports() const { return rftReportList; } bool hasRft(const std::string& wellName, const RftDate& date) const; bool hasRft(const std::string& wellName, int year, int month, int day) const; + std::vector listOfRftArrays(int reportIndex ) const; + std::vector listOfRftArrays(const std::string& wellName, const RftDate& date) const; @@ -63,8 +67,10 @@ public: bool hasArray(const std::string& arrayName, const std::string& wellName, const RftDate& date) const; + int numberOfReports() { return numReports; } + private: - std::map> arrIndexRange; + std::map> arrIndexRange; int numReports; std::vector timeList; @@ -72,9 +78,11 @@ private: std::set dateList; RftReportList rftReportList; - std::map,int> reportIndex; // mapping report index to wellName and date (tupe) + std::map,int> reportIndices; // mapping report index to wellName and date (tupe) int getReportIndex(const std::string& wellName, const RftDate& date) const; + + int getArrayIndex(const std::string& name, int reportIndex) const; int getArrayIndex(const std::string& name, const std::string& wellName, const RftDate& date) const; }; diff --git a/src/opm/io/eclipse/ERft.cpp b/src/opm/io/eclipse/ERft.cpp index a79a7b70c..7cab30024 100644 --- a/src/opm/io/eclipse/ERft.cpp +++ b/src/opm/io/eclipse/ERft.cpp @@ -63,13 +63,11 @@ ERft::ERft(const std::string &filename) : EclFile(filename) } for (size_t i = 0; i < first.size(); i++) { - std::pair range; - range.first = first[i]; - + std::tuple range; if (i == first.size() - 1) { - range.second = listOfArrays.size(); + range = std::make_tuple(first[i], listOfArrays.size()); } else { - range.second = first[i+1]; + range = std::make_tuple(first[i], first[i+1]); } arrIndexRange[i] = range; @@ -78,32 +76,33 @@ ERft::ERft(const std::string &filename) : EclFile(filename) numReports = first.size(); for (size_t i = 0; i < wellName.size(); i++) { - std::pair wellDatePair(wellName[i],dates[i]); - reportIndex[wellDatePair] = i; - rftReportList.push_back(wellDatePair); + std::tuple wellDateTuple = std::make_tuple(wellName[i], dates[i]); + std::tuple wellDateTimeTuple = std::make_tuple(wellName[i], dates[i], timeList[i]); + reportIndices[wellDateTuple] = i; + rftReportList.push_back(wellDateTimeTuple); } } bool ERft::hasRft(const std::string& wellName, const RftDate& date) const { - return reportIndex.find({wellName, date}) != reportIndex.end(); + return reportIndices.find({wellName, date}) != reportIndices.end(); } bool ERft::hasRft(const std::string& wellName, int year, int month, int day) const { RftDate date(year, month, day); - return reportIndex.find({wellName,date}) != reportIndex.end(); + return reportIndices.find({wellName,date}) != reportIndices.end(); } int ERft::getReportIndex(const std::string& wellName, const RftDate& date) const { - std::pair> wellDatePair(wellName,date); - auto rIndIt = reportIndex.find(wellDatePair); + std::tuple> wellDatePair(wellName, date); + auto rIndIt = reportIndices.find(wellDatePair); - if (rIndIt == reportIndex.end()) { + if (rIndIt == reportIndices.end()) { int y = std::get<0>(date); int m = std::get<1>(date); int d = std::get<2>(date); @@ -124,8 +123,8 @@ bool ERft::hasArray(const std::string& arrayName, const std::string& wellName, auto searchInd = arrIndexRange.find(reportInd); - int fromInd = searchInd->second.first; - int toInd = searchInd->second.second; + int fromInd = std::get<0>(searchInd->second); + int toInd = std::get<1>(searchInd->second); auto it = std::find(array_name.begin()+fromInd,array_name.begin()+toInd,arrayName); return it != array_name.begin() + toInd; @@ -139,8 +138,8 @@ int ERft::getArrayIndex(const std::string& name, const std::string& wellName, auto searchInd = arrIndexRange.find(rInd); - int fromInd =searchInd->second.first; - int toInd = searchInd->second.second; + int fromInd =std::get<0>(searchInd->second); + int toInd = std::get<1>(searchInd->second); auto it=std::find(array_name.begin()+fromInd,array_name.begin()+toInd,name); if (std::distance(array_name.begin(),it) == toInd) { @@ -157,6 +156,28 @@ int ERft::getArrayIndex(const std::string& name, const std::string& wellName, } +int ERft::getArrayIndex(const std::string& name, int reportIndex) const +{ + if ((reportIndex < 0) || (reportIndex >= numReports)) { + std::string message = "Report index " + std::to_string(reportIndex) + " not found in RFT file."; + OPM_THROW(std::invalid_argument, message); + } + + auto searchInd = arrIndexRange.find(reportIndex); + int fromInd =std::get<0>(searchInd->second); + int toInd = std::get<1>(searchInd->second); + + auto it=std::find(array_name.begin() + fromInd,array_name.begin() + toInd,name); + + if (std::distance(array_name.begin(),it) == toInd) { + std::string message = "Array " + name + " not found for RFT, rft report index: " + std::to_string(reportIndex); + OPM_THROW(std::invalid_argument, message); + } + + return std::distance(array_name.begin(),it); +} + + template<> const std::vector& ERft::getRft(const std::string& name, const std::string &wellName, const RftDate& date) const @@ -277,6 +298,98 @@ ERft::getRft(const std::string& name, const std::string& wellName, } +template<> const std::vector& +ERft::getRft(const std::string& name, int reportIndex) const +{ + int arrInd = getArrayIndex(name, reportIndex); + + if (array_type[arrInd] != REAL) { + std::string message = "Array " + name + " found in RFT file for selected report, but called with wrong type"; + OPM_THROW(std::runtime_error, message); + } + + auto search_array = real_array.find(arrInd); + return search_array->second; +} + + +template<> const std::vector& +ERft::getRft(const std::string& name, int reportIndex) const +{ + int arrInd = getArrayIndex(name, reportIndex); + + if (array_type[arrInd] != DOUB) { + std::string message = "Array " + name + " !!found in RFT file for selected report, but called with wrong type"; + OPM_THROW(std::runtime_error, message); + } + + auto search_array = doub_array.find(arrInd); + return search_array->second; +} + + +template<> const std::vector& +ERft::getRft(const std::string& name, int reportIndex) const +{ + int arrInd = getArrayIndex(name, reportIndex); + + if (array_type[arrInd] != INTE) { + std::string message = "Array " + name + " !!found in RFT file for selected report, but called with wrong type"; + OPM_THROW(std::runtime_error, message); + } + + auto search_array = inte_array.find(arrInd); + return search_array->second; +} + + +template<> const std::vector& +ERft::getRft(const std::string& name, int reportIndex) const +{ + int arrInd = getArrayIndex(name, reportIndex); + + if (array_type[arrInd] != LOGI) { + std::string message = "Array " + name + " !!found in RFT file for selected report, but called with wrong type"; + OPM_THROW(std::runtime_error, message); + } + + auto search_array = logi_array.find(arrInd); + return search_array->second; +} + + +template<> const std::vector& +ERft::getRft(const std::string& name, int reportIndex) const +{ + int arrInd = getArrayIndex(name, reportIndex); + + if (array_type[arrInd] != CHAR) { + std::string message = "Array " + name + " !!found in RFT file for selected report, but called with wrong type"; + OPM_THROW(std::runtime_error, message); + } + + auto search_array = char_array.find(arrInd); + return search_array->second; +} + + +std::vector ERft::listOfRftArrays(int reportIndex) const +{ + if ((reportIndex < 0) || (reportIndex >= numReports)) { + std::string message = "Report index " + std::to_string(reportIndex) + " not found in RFT file."; + OPM_THROW(std::invalid_argument, message); + } + + std::vector list; + auto searchInd = arrIndexRange.find(reportIndex); + + for (int i = std::get<0>(searchInd->second); i < std::get<1>(searchInd->second); i++) { + list.emplace_back(array_name[i], array_type[i], array_size[i]); + } + + return list; +} + std::vector ERft::listOfRftArrays(const std::string& wellName, const RftDate& date) const { @@ -284,7 +397,7 @@ std::vector ERft::listOfRftArrays(const std::string& wellName int rInd = getReportIndex(wellName, date); auto searchInd = arrIndexRange.find(rInd); - for (int i = searchInd->second.first; i < searchInd->second.second; i++) { + for (int i = std::get<0>(searchInd->second); i < std::get<1>(searchInd->second); i++) { list.emplace_back(array_name[i], array_type[i], array_size[i]); } diff --git a/test_util/EclRegressionTest.cpp b/test_util/EclRegressionTest.cpp index 23da4fc73..d2697e277 100644 --- a/test_util/EclRegressionTest.cpp +++ b/test_util/EclRegressionTest.cpp @@ -1006,16 +1006,16 @@ void ECLRegressionTest::results_rft() if (rftReportList1 != rftReportList2) { std::vector rftList1; for (auto& report : rftReportList1) { - std::string well = report.first; - std::tuple date = report.second; + std::string well = std::get<0>(report); + std::tuple date = std::get<1>(report); std::string str1 = well +" (" + std::to_string(std::get<0>(date)) + "/" + std::to_string(std::get<1>(date)) + "/" + std::to_string(std::get<2>(date)) + ")"; rftList1.push_back(str1); } std::vector rftList2; for (auto& report : rftReportList2) { - std::string well = report.first; - std::tuple date = report.second; + std::string well = std::get<0>(report); + std::tuple date = std::get<1>(report); std::string str2 = well +" (" + std::to_string(std::get<0>(date)) + "/" + std::to_string(std::get<1>(date)) + "/" + std::to_string(std::get<2>(date)) + ")"; rftList2.push_back(str2); } @@ -1026,8 +1026,8 @@ void ECLRegressionTest::results_rft() } for (auto& report : rftReportList2) { - std::string well = report.first; - std::tuple date = report.second; + std::string well = std::get<0>(report); + std::tuple date = std::get<1>(report); std::string dateStr = std::to_string(std::get<0>(date)) + "/" + std::to_string(std::get<1>(date)) + "/" + std::to_string(std::get<2>(date)); diff --git a/tests/test_ERft.cpp b/tests/test_ERft.cpp index b3cf2f7d2..f6b44a41e 100644 --- a/tests/test_ERft.cpp +++ b/tests/test_ERft.cpp @@ -90,12 +90,12 @@ BOOST_AUTO_TEST_CASE(TestERft_1) { Date{2017,7,31}, }; - const std::vector> ref_rftList = { - {"PROD", Date{2015,1, 1}}, - {"INJ" , Date{2015,1, 1}}, - {"A-1H", Date{2015,9, 1}}, - {"B-2H", Date{2016,5,31}}, - {"PROD", Date{2017,7,31}} + const std::vector> ref_rftList = { + {"PROD", Date{2015,1, 1}, 0.00000000e+00}, + {"INJ" , Date{2015,1, 1}, 0.00000000e+00}, + {"A-1H", Date{2015,9, 1}, 0.24300000E+03}, + {"B-2H", Date{2016,5,31}, 0.51600000E+03}, + {"PROD", Date{2017,7,31}, 0.94200000E+03} }; std::string testFile="SPE1CASE1.RFT"; @@ -104,7 +104,7 @@ BOOST_AUTO_TEST_CASE(TestERft_1) { std::vector wellList = rft1.listOfWells(); std::vector rftDates = rft1.listOfdates(); - std::vector> rftList = rft1.listOfRftReports(); + std::vector> rftList = rft1.listOfRftReports(); BOOST_CHECK_EQUAL(wellList==ref_wellList, true); BOOST_CHECK_EQUAL(rftDates==ref_dates, true); @@ -129,7 +129,7 @@ BOOST_AUTO_TEST_CASE(TestERft_1) { BOOST_CHECK_THROW(rft1.hasArray("SGAS","C-2H", Date{2016,5,30}),std::invalid_argument); BOOST_CHECK_THROW(rft1.hasArray("XXXX","C-2H", Date{2016,5,30}),std::invalid_argument); - // test member function getRft + // // test member function getRft(name, wellName, date) std::vector vect1=rft1.getRft("CONIPOS","B-2H", Date{2016,5,31}); std::vector vect2=rft1.getRft("PRESSURE","B-2H", Date{2016,5,31}); @@ -139,6 +139,25 @@ BOOST_AUTO_TEST_CASE(TestERft_1) { BOOST_CHECK_EQUAL(vect2.size(), 3); BOOST_CHECK_EQUAL(vect3.size(), 16); + // test member function getRft(name, reportIndex) + + std::vector vect1a=rft1.getRft("CONIPOS", 3); + BOOST_CHECK_EQUAL(vect1.size(), vect1a.size()); + + std::vector vect2a = rft1.getRft("PRESSURE", 3); + BOOST_CHECK_EQUAL(vect2.size(), vect2a.size()); + + for (size_t t = 0; t < vect2.size(); t++){ + BOOST_CHECK_EQUAL(vect2[t], vect2a[t]); + } + + std::vector vect3a = rft1.getRft("WELLETC", 3); + BOOST_CHECK_EQUAL(vect2.size(), vect2a.size()); + + for (size_t t = 0; t < vect3.size(); t++){ + BOOST_CHECK_EQUAL(vect3[t], vect3a[t]); + } + // called with invalid argument, array not existing, wrong well name or wrong date BOOST_CHECK_THROW(std::vector vect11=rft1.getRft("CONIPOS","C-2H", Date{2016,5,31}),std::invalid_argument); BOOST_CHECK_THROW(std::vector vect11=rft1.getRft("CONIPOS","B-2H", Date{2016,5,30}),std::invalid_argument); diff --git a/tests/test_RFT.cpp b/tests/test_RFT.cpp index 1a5802a51..00e7ceb10 100644 --- a/tests/test_RFT.cpp +++ b/tests/test_RFT.cpp @@ -325,7 +325,7 @@ namespace { >{}; for (const auto& wellDate : rft.listOfRftReports()) { - dates[wellDate.first].push_back(wellDate.second); + dates[std::get<0>(wellDate)].push_back(std::get<1>(wellDate)); } // Well OP_1