diff --git a/opm/io/eclipse/ERst.hpp b/opm/io/eclipse/ERst.hpp index eeb45ddd6..a0fa0f741 100644 --- a/opm/io/eclipse/ERst.hpp +++ b/opm/io/eclipse/ERst.hpp @@ -51,6 +51,7 @@ public: } int count(const std::string& name, int reportStepNumber) const; + size_t numberOfReportSteps() const { return seqnum.size(); }; const std::vector& listOfReportStepNumbers() const { return seqnum; } diff --git a/python/cxx/eclipse_io.cpp b/python/cxx/eclipse_io.cpp index 221608560..416ac67d4 100644 --- a/python/cxx/eclipse_io.cpp +++ b/python/cxx/eclipse_io.cpp @@ -4,6 +4,7 @@ #include #include +#include #include "export.hpp" #include "converters.hpp" @@ -76,6 +77,52 @@ npArray get_vector_occurrence(Opm::EclIO::EclFile * file_ptr, const std::string& return get_vector_index(file_ptr, array_index); } +bool erst_contains(Opm::EclIO::ERst * file_ptr, std::tuple keyword) +{ + bool hasKeyAtReport = file_ptr->count(std::get<0>(keyword), std::get<1>(keyword)) > 0 ? true : false; + return hasKeyAtReport; +} + +npArray get_erst_by_index(Opm::EclIO::ERst * file_ptr, size_t index, size_t rstep) +{ + auto arrList = file_ptr->listOfRstArrays(rstep); + + if (index >=arrList.size()) + throw std::out_of_range("Array index out of range. "); + + auto array_type = std::get<1>(arrList[index]); + + if (array_type == Opm::EclIO::INTE) + return std::make_tuple (convert::numpy_array( file_ptr->getRst(index, rstep)), array_type); + + if (array_type == Opm::EclIO::REAL) + return std::make_tuple (convert::numpy_array( file_ptr->getRst(index, rstep)), array_type); + + if (array_type == Opm::EclIO::DOUB) + return std::make_tuple (convert::numpy_array( file_ptr->getRst(index, rstep)), array_type); + + if (array_type == Opm::EclIO::LOGI) + return std::make_tuple (convert::numpy_array( file_ptr->getRst(index, rstep)), array_type); + + if (array_type == Opm::EclIO::CHAR) + return std::make_tuple (convert::numpy_string_array( file_ptr->getRst(index, rstep)), array_type); + + throw std::logic_error("Data type not supported"); +} + + +npArray get_erst_vector(Opm::EclIO::ERst * file_ptr, const std::string& key, size_t rstep, size_t occurrence) +{ + if (occurrence >= static_cast(file_ptr->count(key, rstep))) + throw std::out_of_range("file have less than " + std::to_string(occurrence + 1) + " arrays in selected report step"); + + auto array_list = file_ptr->listOfRstArrays(rstep); + + size_t array_index = get_array_index(array_list, key, 0); + + return get_erst_by_index(file_ptr, array_index, rstep); +} + } @@ -102,4 +149,17 @@ void python::common::export_IO(py::module& m) { .def("__get_data", &get_vector_index) .def("__get_data", &get_vector_name) .def("__get_data", &get_vector_occurrence); + + py::class_(m, "ERst") + .def(py::init()) + .def("__has_report_step", &Opm::EclIO::ERst::hasReportStepNumber) + .def("load_report_step", &Opm::EclIO::ERst::loadReportStepNumber) + .def_property_readonly("report_steps", &Opm::EclIO::ERst::listOfReportStepNumbers) + .def("__len__", &Opm::EclIO::ERst::numberOfReportSteps) + .def("count", &Opm::EclIO::ERst::count) + .def("__contains", &erst_contains) + .def("__get_list_of_arrays", &Opm::EclIO::ERst::listOfRstArrays) + .def("__get_data", &get_erst_by_index) + .def("__get_data", &get_erst_vector); + } diff --git a/python/python/opm/_common.py b/python/python/opm/_common.py index 7f5c1c59f..4eb433345 100644 --- a/python/python/opm/_common.py +++ b/python/python/opm/_common.py @@ -21,6 +21,7 @@ from .libopmcommon_python import Schedule from .libopmcommon_python import OpmLog from .libopmcommon_python import SummaryConfig from .libopmcommon_python import EclFile, eclArrType +from .libopmcommon_python import ERst from .libopmcommon_python import SummaryState #from .schedule import Well, Connection, Schedule diff --git a/python/python/opm/io/ecl/__init__.py b/python/python/opm/io/ecl/__init__.py index 4b467a8b1..e863983a6 100644 --- a/python/python/opm/io/ecl/__init__.py +++ b/python/python/opm/io/ecl/__init__.py @@ -1,5 +1,6 @@ from opm._common import eclArrType from opm._common import EclFile +from opm._common import ERst import sys import datetime @@ -33,8 +34,59 @@ def getitem_eclfile(self, arg): return data +def erst_get_list_of_arrays(self, arg): + + if sys.version_info.major==2: + rawData = self.__get_list_of_arrays(arg) + return [ ( x[0].encode("utf-8"), x[1], x[2] ) for x in rawData ] + else: + return self.__get_list_of_arrays(arg) + + +def getitem_erst(self, arg): + + if not isinstance(arg, tuple): + raise ValueError("expecting tuple argument, (index, rstep), (name, rstep) or (name, rstep, occurrence) ") + + if len(arg) == 2: + if isinstance(arg[0], int): + data, array_type = self.__get_data(arg[0], int(arg[1])) + else: + data, array_type = self.__get_data(str(arg[0]), int(arg[1]), 0) # default first occurrence + elif len(arg) == 3: + data, array_type = self.__get_data(str(arg[0]), int(arg[1]), int(arg[2])) + else: + raise ValueError("expecting tuple argument with 2 or 3 argumens: (index, rstep), (name, rstep) or (name, rstep, occurrence) ") + + if array_type == eclArrType.CHAR: + return [ x.decode("utf-8") for x in data ] + + return data + + +def contains_erst(self, arg): + + if isinstance(arg, tuple): + if len(arg) == 2: + return self.__contains((arg[0], arg[1])) + else: + raise ValueError("expecting tuple (array name , report step number) or \ + or report step number") + + elif isinstance(arg, int): + return self.__has_report_step(arg) + + else: + raise ValueError("expecting tuple (array name , report step number) or \ + or report step number") + setattr(EclFile, "__getitem__", getitem_eclfile) setattr(EclFile, "arrays", eclfile_get_list_of_arrays) +setattr(ERst, "__contains__", contains_erst) +setattr(ERst, "arrays", erst_get_list_of_arrays) +setattr(ERst, "__getitem__", getitem_erst) + + diff --git a/python/tests/test_erst.py b/python/tests/test_erst.py new file mode 100755 index 000000000..9d28d78d2 --- /dev/null +++ b/python/tests/test_erst.py @@ -0,0 +1,168 @@ +import unittest +import sys +import numpy as np +import io + +from opm.io.ecl import ERst, eclArrType +try: + from tests.utils import test_path +except ImportError: + from utils import test_path + + +class TestERst(unittest.TestCase): + + def test_reportSteps(self): + + rst1 = ERst(test_path("data/SPE9.UNRST")) + + self.assertTrue( 37 in rst1) + self.assertTrue( 74 in rst1) + + self.assertRaises(ValueError, rst1.load_report_step, 137 ); + + rst1.load_report_step(37); + rst1.load_report_step(74); + + self.assertTrue(37 in rst1.report_steps) + self.assertTrue(74 in rst1.report_steps) + + self.assertEqual(len(rst1), 2) + + + def test_contains(self): + + rst1 = ERst(test_path("data/SPE9.UNRST")) + + self.assertTrue(("PRESSURE", 37) in rst1) + self.assertTrue(("SWAT", 74) in rst1) + + arrayList37=rst1.arrays(37) + + self.assertEqual(len(arrayList37), 21) + + arrName, arrType, arrSize = arrayList37[16] + + self.assertEqual(arrName, "PRESSURE") + self.assertEqual(arrType, eclArrType.REAL) + self.assertEqual(arrSize, 9000) + + + def test_getitem(self): + + rst1 = ERst(test_path("data/SPE9.UNRST")) + + # get first occurrence of ZWEL, report step 37 + zwel1 = rst1[11, 37] + + zwel2 = rst1["ZWEL",37, 0] + zwel3 = rst1["ZWEL",37] + + for v1,v2 in zip (zwel1, zwel2): + self.assertEqual(v1, v2) + + for v1,v2 in zip (zwel1, zwel3): + self.assertEqual(v1, v2) + + self.assertEqual(len(zwel1), 78) + + self.assertEqual(zwel1[0], "INJE1") + self.assertEqual(zwel1[3], "PRODU2") + self.assertEqual(zwel1[6], "PRODU3") + + # get first occurrence of INTEHEAD, report step 37 + inteh = rst1["INTEHEAD",37] + + self.assertEqual(len(inteh), 411) + self.assertTrue(isinstance(inteh, np.ndarray)) + self.assertEqual(inteh.dtype, "int32") + + self.assertEqual(inteh[1], 201702) + self.assertEqual(inteh[9], 25) + self.assertEqual(inteh[64], 6) + self.assertEqual(inteh[65], 1) + self.assertEqual(inteh[66], 2016) + + # get first occurrence of PRESSURE, report step 74 + pres74 = rst1["PRESSURE",74] + + self.assertTrue(isinstance(pres74, np.ndarray)) + self.assertEqual(pres74.dtype, "float32") + self.assertEqual(len(pres74), 9000) + + self.assertAlmostEqual(pres74[0], 2290.9192, 4) + self.assertAlmostEqual(pres74[1], 2254.6619, 4) + self.assertAlmostEqual(pres74[2], 2165.5347, 4) + self.assertAlmostEqual(pres74[3], 1996.2598, 4) + + xcon = rst1["XCON", 74] + self.assertTrue(isinstance(xcon, np.ndarray)) + self.assertEqual(xcon.dtype, "float64") + self.assertEqual(len(xcon), 7540) + + self.assertAlmostEqual(xcon[1], -22.841887080742975, 10) + + logih = rst1["LOGIHEAD", 74] + self.assertTrue(isinstance(logih, np.ndarray)) + self.assertEqual(len(logih), 121) + + for b1, b2 in zip([True, True, False, False, False], logih[0:5]): + self.assertEqual(b1, b2) + + + def test_getby_index(self): + + rst1 = ERst(test_path("data/SPE9.UNRST")) + + self.assertTrue(("SGAS", 74) in rst1) + self.assertTrue(rst1.count("SGAS", 74), 1) + + arrayList74=rst1.arrays(74) + + array_name_list = [item[0] for item in arrayList74] + ind = array_name_list.index("SGAS") + + sgas_a = rst1[ind, 74] + sgas_b = rst1["SGAS", 74] + + self.assertEqual(len(sgas_a), len(sgas_b)) + + for sg1, sg2 in zip(sgas_a, sgas_b): + self.assertEqual(sg1, sg2) + + + def test_list_of_arrays(self): + + refArrList = ["SEQNUM", "INTEHEAD", "LOGIHEAD", "DOUBHEAD", "IGRP", "SGRP", "XGRP", "ZGRP", "IWEL", + "SWEL","XWEL","ZWEL", "ICON", "SCON", "XCON", "STARTSOL","PRESSURE", "RS", "SGAS", "SWAT", "ENDSOL"] + + rst1 = ERst(test_path("data/SPE9.UNRST")) + + array_list_74=rst1.arrays(74) + + self.assertEqual(len(refArrList), len(array_list_74) ) + + for n, (name, arrType, arrSize) in enumerate(array_list_74): + + self.assertEqual(name, refArrList[n]) + + if arrType != eclArrType.MESS: + array = rst1[name, 74] + self.assertEqual(len(array), arrSize) + + if arrType == eclArrType.INTE: + self.assertEqual(array.dtype, "int32") + elif arrType == eclArrType.REAL: + self.assertEqual(array.dtype, "float32") + elif arrType == eclArrType.DOUB: + self.assertEqual(array.dtype, "float64") + elif arrType == eclArrType.LOGI: + self.assertEqual(array.dtype, "bool") + elif arrType == eclArrType.CHAR: + self.assertTrue(isinstance(array, list)) + + +if __name__ == "__main__": + + unittest.main() +