added: ERst class for reading restart simulation data from ECL files

CMakeLists.txt

@@ -166,6 +166,7 @@ if(ENABLE_ECL_INPUT)
               examples/test_util/EclUtil.cpp
               examples/test_util/EGrid.cpp
               examples/test_util/ERft.cpp
+              examples/test_util/ERst.cpp
               examples/test_util/summaryComparator.cpp
               examples/test_util/summaryIntegrationTest.cpp
               examples/test_util/summaryRegressionTest.cpp)
@@ -177,7 +178,7 @@ if(ENABLE_ECL_INPUT)
   set(_libs testutil opmcommon
             ${Boost_UNIT_TEST_FRAMEWORK_LIBRARY})
   foreach(test test_compareSummary test_EclFilesComparator test_EclIO
-               test_EGrid test_ERft)
+               test_EGrid test_ERft test_ERst)
     opm_add_test(${test} CONDITION ENABLE_ECL_INPUT
                          LIBRARIES ${_libs}
                          WORKING_DIRECTORY ${PROJECT_BINARY_DIR}/tests)

CMakeLists_files.cmake

@@ -323,6 +323,8 @@ if(ENABLE_ECL_INPUT)
     tests/ECLFILE.FINIT
     tests/SPE1CASE1.EGRID
     tests/SPE1CASE1.RFT
+    tests/SPE1_TESTCASE.UNRST
+    tests/SPE1_TESTCASE.FUNRST
   )
   list (APPEND EXAMPLE_SOURCE_FILES
     examples/opmi.cpp

examples/test_util/ERst.cpp

@@ -0,0 +1,178 @@
+/*
+   Copyright 2019 Equinor ASA.
+
+   This file is part of the Open Porous Media project (OPM).
+
+   OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by
+   the Free Software Foundation, either version 3 of the License, or
+   (at your option) any later version.
+
+   OPM is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+   */
+
+#include "ERst.hpp"
+
+#include <string>
+#include <string.h>
+#include <sstream>
+#include <iterator>
+#include <iomanip>
+#include <algorithm>
+
+
+ERst::ERst(const std::string& filename) : EclFile(filename)
+{
+    loadData("SEQNUM");
+
+    std::vector<int> firstIndex;
+
+    for (size_t i = 0;  i < array_name.size(); i++) {
+        if (array_name[i] == "SEQNUM") {
+            auto seqn = get<int>(i);
+            seqnum.push_back(seqn[0]);
+            firstIndex.push_back(i);
+        }
+    }
+
+
+    for (size_t i = 0; i < seqnum.size(); i++) {
+        std::pair<int,int> range;
+        range.first = firstIndex[i];
+
+        if (i != seqnum.size() - 1) {
+            range.second = firstIndex[i+1];
+        } else {
+            range.second = array_name.size();
+        }
+
+        arrIndexRange[seqnum[i]] = range;
+    }
+
+    nReports = seqnum.size();
+
+    for (int i = 0; i < nReports; i++) {
+        reportLoaded[seqnum[i]] = false;
+    }
+};
+
+
+bool ERst::hasReportStepNumber(int number) const
+{
+    auto search = arrIndexRange.find(number);
+    return search != arrIndexRange.end();
+}
+
+
+void ERst::loadReportStepNumber(int number)
+{
+    if (!hasReportStepNumber(number)) {
+        std::string message="Trying to load non existing report step number " + std::to_string(number);
+        OPM_THROW(std::invalid_argument, message);
+    }
+
+    std::vector<int> arrayIndexList;
+    arrayIndexList.reserve(arrIndexRange[number].second - arrIndexRange[number].first + 1);
+
+    for (int i = arrIndexRange[number].first; i < arrIndexRange[number].second; i++) {
+        arrayIndexList.push_back(i);
+    }
+
+    loadData(arrayIndexList);
+
+    reportLoaded[number] = true;
+}
+
+
+std::vector<EclFile::EclEntry> ERst::listOfRstArrays(int reportStepNumber)
+{
+    std::vector<EclEntry> list;
+
+    if (!hasReportStepNumber(reportStepNumber)) {
+        std::string message = "Trying to get list of arrays from non existing report step number  " + std::to_string(reportStepNumber);
+        OPM_THROW(std::invalid_argument, message);
+    }
+
+    const auto& rng = this->arrIndexRange[reportStepNumber];
+    list.reserve(rng.second - rng.first);
+
+    for (int i = rng.first;  i < rng.second; i++) {
+        list.emplace_back(array_name[i], array_type[i], array_size[i]);
+    }
+
+    return list;
+}
+
+
+int ERst::getArrayIndex(const std::string& name, int number) const
+{
+    if (!hasReportStepNumber(number)) {
+        std::string message = "Trying to get vector " + name + " from non existing sequence " + std::to_string(number);
+        OPM_THROW(std::invalid_argument, message);
+    }
+
+    auto search = reportLoaded.find(number);
+
+    if (!search->second) {
+        std::string message = "Data not loaded for sequence " + std::to_string(number);
+        OPM_THROW(std::runtime_error, message);
+    }
+
+    auto range_it = arrIndexRange.find(number);
+
+    std::pair<int,int> indexRange = range_it->second;
+
+    auto it = std::find(array_name.begin() + indexRange.first,
+                        array_name.begin() + indexRange.second, name);
+
+    if (std::distance(array_name.begin(),it) == indexRange.second) {
+        std::string message = "Array " + name + " not found in sequence " + std::to_string(number);
+        OPM_THROW(std::runtime_error, message);
+    }
+
+    return std::distance(array_name.begin(), it);
+}
+
+
+template<>
+const std::vector<int>& ERst::getRst<int>(const std::string& name, int reportStepNumber)
+{
+    int ind = getArrayIndex(name, reportStepNumber);
+    return getImpl(ind, EIOD::INTE, inte_array, "integer");
+}
+
+
+template<>
+const std::vector<float>& ERst::getRst<float>(const std::string& name, int reportStepNumber)
+{
+    int ind = getArrayIndex(name, reportStepNumber);
+    return getImpl(ind, EIOD::REAL, real_array, "float");
+}
+
+
+template<>
+const std::vector<double>& ERst::getRst<double>(const std::string& name, int reportStepNumber)
+{
+    int ind = getArrayIndex(name, reportStepNumber);
+    return getImpl(ind, EIOD::DOUB, doub_array, "double");
+}
+
+
+template<>
+const std::vector<bool>& ERst::getRst<bool>(const std::string& name, int reportStepNumber)
+{
+    int ind = getArrayIndex(name, reportStepNumber);
+    return getImpl(ind, EIOD::LOGI, logi_array, "bool");
+}
+
+template<>
+const std::vector<std::string>& ERst::getRst<std::string>(const std::string& name, int reportStepNumber)
+{
+    int ind = getArrayIndex(name, reportStepNumber);
+    return getImpl(ind, EIOD::CHAR, char_array, "string");
+}

examples/test_util/ERst.hpp

@@ -0,0 +1,57 @@
+/*
+   Copyright 2019 Equinor ASA.
+
+   This file is part of the Open Porous Media project (OPM).
+
+   OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by
+   the Free Software Foundation, either version 3 of the License, or
+   (at your option) any later version.
+
+   OPM is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+   */
+
+#ifndef ERST_HPP
+#define ERST_HPP
+
+
+#include "EclFile.hpp"
+
+#include <iostream>
+#include <string>
+#include <fstream>
+#include <vector>
+#include <ctime>
+#include <map>
+
+
+class ERst : public EclFile
+{
+public:
+    explicit ERst(const std::string& filename);
+    bool hasReportStepNumber(int number) const;
+
+    void loadReportStepNumber(int number);
+
+    template <typename T>
+    const std::vector<T>& getRst(const std::string& name, int reportStepNumber);
+
+    const std::vector<int>& listOfReportStepNumbers() const { return seqnum; }
+
+    std::vector<EclEntry> listOfRstArrays(int reportStepNumber);
+
+private:
+    int nReports;
+    std::vector<int> seqnum;                           // report step numbers, from SEQNUM array in restart file
+    std::unordered_map<int,bool> reportLoaded;
+    std::map<int, std::pair<int,int>> arrIndexRange;   // mapping report step number to array indeces (start and end)
+
+    int getArrayIndex(const std::string& name, int seqnum) const;
+};
+
+#endif

tests/test_ERst.cpp

@@ -0,0 +1,219 @@
+/*
+   +   Copyright 2019 Equinor ASA.
+   +
+   +   This file is part of the Open Porous Media project (OPM).
+   +
+   +   OPM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by
+   +   the Free Software Foundation, either version 3 of the License, or
+   +   (at your option) any later version.
+   +
+   +   OPM is distributed in the hope that it will be useful,
+   +   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   +   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   +   GNU General Public License for more details.
+   +
+   +   You should have received a copy of the GNU General Public License
+   +   along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+   +   */
+
+
+#include "config.h"
+#include <iostream>
+#include <stdio.h>
+#include <iomanip>
+#include <math.h>
+#include <algorithm>
+#include <tuple>
+
+#include <examples/test_util/ERst.hpp>
+#include <examples/test_util/EclOutput.hpp>
+
+#define BOOST_TEST_MODULE Test EclIO
+#include <boost/test/unit_test.hpp>
+
+template<typename InputIterator1, typename InputIterator2>
+bool
+range_equal(InputIterator1 first1, InputIterator1 last1,
+            InputIterator2 first2, InputIterator2 last2)
+{
+    while(first1 != last1 && first2 != last2)
+    {
+        if(*first1 != *first2) return false;
+        ++first1;
+        ++first2;
+    }
+    return (first1 == last1) && (first2 == last2);
+}
+
+bool compare_files(const std::string& filename1, const std::string& filename2)
+{
+    std::ifstream file1(filename1);
+    std::ifstream file2(filename2);
+
+    std::istreambuf_iterator<char> begin1(file1);
+    std::istreambuf_iterator<char> begin2(file2);
+
+    std::istreambuf_iterator<char> end;
+
+    return range_equal(begin1, end, begin2, end);
+}
+
+
+template <typename T>
+bool operator==(const std::vector<T> & t1, const std::vector<T> & t2)
+{
+    return std::equal(t1.begin(), t1.end(), t2.begin(), t2.end());
+}
+
+
+BOOST_AUTO_TEST_CASE(TestERst_1) {
+
+    std::string testFile="SPE1_TESTCASE.UNRST";
+    std::vector<int> refReportStepNumbers= {1,2,5,10,15,25,50,100,120};
+
+    ERst rst1(testFile);
+    rst1.loadReportStepNumber(5);
+    
+    std::vector<int> reportStepNumbers = rst1.listOfReportStepNumbers();
+
+    BOOST_CHECK_EQUAL(reportStepNumbers==refReportStepNumbers, true);
+
+    BOOST_CHECK_EQUAL(rst1.hasReportStepNumber(4), false);
+    BOOST_CHECK_EQUAL(rst1.hasReportStepNumber(5), true);
+
+    // try loading non-existing report step, should throw exception
+    BOOST_CHECK_THROW(rst1.loadReportStepNumber(4) , std::invalid_argument );
+
+
+    // try to get a list of vectors from non-existing report step, should throw exception
+    std::vector<std::tuple<std::string, EIOD::eclArrType, int>> rstArrays; // = rst1.listOfRstArrays(4);
+    BOOST_CHECK_THROW(rstArrays = rst1.listOfRstArrays(4), std::invalid_argument);
+
+    // non exising report step number, should throw exception 
+    
+    BOOST_CHECK_THROW(std::vector<int> vect1=rst1.getRst<int>("ICON",0) , std::invalid_argument );
+    BOOST_CHECK_THROW(std::vector<float> vect2=rst1.getRst<float>("PRESSURE",0) , std::invalid_argument );
+    BOOST_CHECK_THROW(std::vector<double> vect3=rst1.getRst<double>("XGRP",0) , std::invalid_argument );
+    BOOST_CHECK_THROW(std::vector<bool> vect4=rst1.getRst<bool>("LOGIHEAD",0) , std::invalid_argument );
+    BOOST_CHECK_THROW(std::vector<std::string> vect4=rst1.getRst<std::string>("ZWEL",0) , std::invalid_argument );
+
+    // report step number exists, but data is not loaded. Should throw exception
+    BOOST_CHECK_THROW(std::vector<int> vect1=rst1.getRst<int>("ICON",10) , std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<float> vect2=rst1.getRst<float>("PRESSURE",10) , std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<double> vect3=rst1.getRst<double>("XGRP",10) , std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<bool> vect4=rst1.getRst<bool>("LOGIHEAD",10) , std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<std::string> vect4=rst1.getRst<std::string>("ZWEL",10) , std::runtime_error );
+
+    // calling getRst<T> member function with wrong type, should throw exception
+
+    BOOST_CHECK_THROW(std::vector<float> vect1=rst1.getRst<float>("ICON",5) , std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<int> vect2=rst1.getRst<int>("PRESSURE",5), std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<float> vect3=rst1.getRst<float>("XGRP",5), std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<double> vect4=rst1.getRst<double>("LOGIHEAD",5), std::runtime_error );
+    BOOST_CHECK_THROW(std::vector<bool> vect5=rst1.getRst<bool>("ZWEL",5), std::runtime_error );
+
+    rst1.loadReportStepNumber(25);
+
+    std::vector<int> vect1 = rst1.getRst<int>("ICON",25);
+    std::vector<float> vect2 = rst1.getRst<float>("PRESSURE",25);
+    std::vector<double> vect3 = rst1.getRst<double>("XGRP",25);
+    std::vector<bool> vect4 = rst1.getRst<bool>("LOGIHEAD",25);
+    std::vector<std::string> vect5 = rst1.getRst<std::string>("ZWEL",25);
+}
+
+
+static void readAndWrite(EclOutput& eclTest, ERst& rst1,
+                         const std::string& name, int seqnum,
+                         EIOD::eclArrType arrType)
+{
+    if (arrType == EIOD::INTE) {
+        std::vector<int> vect = rst1.getRst<int>(name, seqnum);
+        eclTest.write(name, vect);
+    } else if (arrType == EIOD::REAL) {
+        std::vector<float> vect = rst1.getRst<float>(name, seqnum);
+        eclTest.write(name, vect);
+    } else if (arrType == EIOD::DOUB) {
+        std::vector<double> vect = rst1.getRst<double>(name, seqnum);
+        eclTest.write(name, vect);
+    } else if (arrType == EIOD::LOGI) {
+        std::vector<bool> vect = rst1.getRst<bool>(name, seqnum);
+        eclTest.write(name, vect);
+    } else if (arrType == EIOD::CHAR) {
+        std::vector<std::string> vect = rst1.getRst<std::string>(name, seqnum);
+        eclTest.write(name, vect);
+    } else if (arrType == EIOD::MESS) {
+        eclTest.write(name, std::vector<char>());
+    } else {
+        std::cout << "unknown type " << std::endl;
+        exit(1);
+    }
+}
+
+
+BOOST_AUTO_TEST_CASE(TestERst_2) {
+    
+    std::string testFile="SPE1_TESTCASE.UNRST";
+    std::string outFile="TEST.UNRST";
+
+    // using API for ERst to read all array from a binary unified restart file1
+    // Then write the data back to a new file and check that new file is identical with input file 
+    
+    ERst rst1(testFile);
+
+    {
+        EclOutput eclTest(outFile, false);
+
+        std::vector<int> seqnums = rst1.listOfReportStepNumbers();
+
+        for (size_t i = 0; i < seqnums.size(); i++) {
+            rst1.loadReportStepNumber(seqnums[i]);
+
+            auto rstArrays = rst1.listOfRstArrays(seqnums[i]);
+
+            for (auto& array : rstArrays) {
+                std::string name = std::get<0>(array);
+                EIOD::eclArrType arrType = std::get<1>(array);
+                readAndWrite(eclTest, rst1, name, seqnums[i], arrType);
+            }
+        }
+    }
+
+    BOOST_CHECK_EQUAL(compare_files(testFile, outFile), true);
+
+    if (remove(outFile.c_str()) == -1) {
+        std::cout << " > Warning! temporary file was not deleted" << std::endl;
+    };
+}
+
+BOOST_AUTO_TEST_CASE(TestERst_3) {
+
+    std::string testFile="SPE1_TESTCASE.FUNRST";
+    std::string outFile="TEST.FUNRST";
+
+    // using API for ERst to read all array from a formatted unified restart file1
+    // Then write the data back to a new file and check that new file is identical with input file 
+
+    ERst rst1(testFile);
+
+    EclOutput eclTest(outFile, true);
+
+    std::vector<int> seqnums = rst1.listOfReportStepNumbers();
+    for (unsigned int i=0; i<seqnums.size(); i++) {
+        rst1.loadReportStepNumber(seqnums[i]);
+
+        auto rstArrays = rst1.listOfRstArrays(seqnums[i]);
+
+        for (auto& array : rstArrays) {
+            std::string name = std::get<0>(array);
+            EIOD::eclArrType arrType = std::get<1>(array);
+            readAndWrite(eclTest, rst1, name, seqnums[i], arrType);
+        }
+    }
+
+    BOOST_CHECK_EQUAL(compare_files(testFile, outFile), true);
+
+    if (remove(outFile.c_str() )== -1) {
+        std::cout << " > Warning! temporary file was not deleted" << std::endl;
+    };
+}
+