diff --git a/CMakeLists.txt b/CMakeLists.txt index acbc8b450..aa2d4dd6c 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,18 +1,11 @@ -cmake_minimum_required( - VERSION 2.8.12) # If you are creating Python wrappers for Windows, the actual - # version requirement is 3.4 +cmake_minimum_required(VERSION 3.12) project(ecl C CXX) include(GNUInstallDirs) include(TestBigEndian) -if(NOT DEFINED CMAKE_MACOSX_RPATH) - # There is some weirdness around this variable, the default value is different - # depending on the cmake version, see policy CMP0042. A more explicit way to - # treat this would be `cmake_policy(SET CMP0042 NEW)` but that would fail on - # CMake 2.8.12 (the variable exists but the policy doesn't) - set(CMAKE_MACOSX_RPATH ON) -endif() +set(CMAKE_CXX_STANDARD 17) +set(CMAKE_CXX_STANDARD_REQUIRED TRUE) # ----------------------------------------------------------------- @@ -174,7 +167,6 @@ if(MSVC) endif() list(APPEND CMAKE_MODULE_PATH ${PROJECT_SOURCE_DIR}/cmake/Modules) -find_package(CXX11Features) # ----------------------------------------------------------------- @@ -402,21 +394,37 @@ if(ENABLE_PYTHON AND ERT_WINDOWS) endif() endif() -if(ENABLE_PYTHON) - # If finding the Python interpreter and required packages fails in the - # python/CMakeLists.txt file the ENABLE_PYTHON will be set to OFF. - add_subdirectory(python) -endif() - if(ENABLE_PYTHON) if(NOT ${BUILD_SHARED_LIBS}) message(FATAL_ERROR "The Python wrappers require shared libraries") endif() endif() -if(NOT SKBUILD) - # Avoid installing when calling from python setup.py +if(SKBUILD) + # Building from Python + install(TARGETS ecl LIBRARY DESTINATION python/ecl/.libs) + install(TARGETS _native DESTINATION python/ecl) + install(DIRECTORY lib/include/ DESTINATION python/ecl/.include) + install(DIRECTORY ${CMAKE_BINARY_DIR}/lib/include/ + DESTINATION python/ecl/.include) + if(APPLICATIONS) + install(TARGETS ${APPLICATIONS} RUNTIME DESTINATION python/ecl/.bin) + endif() +else() install(EXPORT ecl-config DESTINATION share/cmake/ecl) export(TARGETS ecl FILE eclConfig.cmake) export(PACKAGE ecl) + + install( + TARGETS ecl + EXPORT ecl-config + ARCHIVE DESTINATION ${CMAKE_INSTALL_LIBDIR} + LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} + RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) + install(DIRECTORY include/ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}) + install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/include/ + DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}) + if(APPLICATIONS) + install(TARGETS ${APPLICATIONS} RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) + endif() endif() diff --git a/applications/CMakeLists.txt b/applications/CMakeLists.txt index 3b23ba6d3..b2ca2fa62 100644 --- a/applications/CMakeLists.txt +++ b/applications/CMakeLists.txt @@ -1,10 +1,7 @@ project(libecl-applications) function(target_link_libecl target) - target_link_libraries(${target} ecl) - if(SKBUILD) - set_target_properties(${target} PROPERTIES INSTALL_RPATH "$ORIGIN/../.libs") - endif() + target_link_libraries(${target} ecl-static) endfunction() if(BUILD_APPLICATIONS) @@ -94,9 +91,8 @@ if(apps) PROPERTY LINK_FLAGS ${OpenMP_C_FLAGS}) endforeach() endif() - install( - TARGETS ${apps} - ARCHIVE DESTINATION ${CMAKE_INSTALL_LIBDIR} - LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} - RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) + + set(APPLICATIONS + ${apps} + PARENT_SCOPE) endif() diff --git a/cmake/Modules/FindCXX11Features.cmake b/cmake/Modules/FindCXX11Features.cmake deleted file mode 100644 index 03090334d..000000000 --- a/cmake/Modules/FindCXX11Features.cmake +++ /dev/null @@ -1,25 +0,0 @@ -# -# Module that checks for supported C++11 (former C++0x) features. -# - -include(CheckCXXCompilerFlag) -if(NOT MSVC) - check_cxx_compiler_flag("-std=c++11" COMPILER_SUPPORTS_CXX11) - check_cxx_compiler_flag("-std=c++0x" COMPILER_SUPPORTS_CXX0X) - - if(COMPILER_SUPPORTS_CXX11) - set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11") - elseif(COMPILER_SUPPORTS_CXX0X) - set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++0x") - else() - message( - SEND_ERROR - "The compiler ${CMAKE_CXX_COMPILER} has no C++11 support. Please use a different C++ compiler." - ) - endif() -endif() - -set(CMAKE_CXX_STANDARD_REQUIRED ON) -set(CMAKE_CXX_EXTENSIONS OFF) - -set(CXX11FEATURES_FOUND TRUE) diff --git a/lib/CMakeLists.txt b/lib/CMakeLists.txt index 3bf0ad518..356ceee8d 100644 --- a/lib/CMakeLists.txt +++ b/lib/CMakeLists.txt @@ -57,93 +57,97 @@ configure_file(ert_api_config.h.in include/ert/util/ert_api_config.h) configure_file(build_config.hpp.in include/ert/util/build_config.hpp) configure_file(ert_api_config.hpp.in include/ert/util/ert_api_config.hpp) -add_library( - ecl - util/rng.cpp - util/lookup_table.cpp - util/statistics.cpp - util/mzran.cpp - util/hash_node.cpp - util/hash_sll.cpp - util/path.cpp - util/hash.cpp - util/node_data.cpp - util/node_ctype.cpp - util/util.cpp - util/util_symlink.cpp - util/util_lfs.c - util/util_unlink.cpp - util/vector.cpp - util/parser.cpp - util/stringlist.cpp - util/buffer.cpp - util/string_util.cpp - util/type_vector_functions.cpp - util/ecl_version.cpp - util/perm_vector.cpp - util/test_util.cpp - util/cxx_string_util.cpp - ${opt_srcs} - ecl/ecl_rsthead.cpp - ecl/ecl_sum_tstep.cpp - ecl/ecl_rst_file.cpp - ecl/ecl_init_file.cpp - ecl/ecl_grid_cache.cpp - ecl/smspec_node.cpp - ecl/ecl_kw_grdecl.cpp - ecl/ecl_file_kw.cpp - ecl/ecl_file_view.cpp - ecl/ecl_grav.cpp - ecl/ecl_grav_calc.cpp - ecl/ecl_smspec.cpp - ecl/ecl_unsmry_loader.cpp - ecl/ecl_sum_data.cpp - ecl/ecl_sum_file_data.cpp - ecl/ecl_util.cpp - ecl/ecl_kw.cpp - ecl/ecl_sum.cpp - ecl/ecl_sum_vector.cpp - ecl/fortio.c - ecl/ecl_rft_file.cpp - ecl/ecl_rft_node.cpp - ecl/ecl_rft_cell.cpp - ecl/ecl_grid.cpp - ecl/ecl_coarse_cell.cpp - ecl/ecl_box.cpp - ecl/ecl_io_config.cpp - ecl/ecl_file.cpp - ecl/ecl_region.cpp - ecl/ecl_subsidence.cpp - ecl/ecl_grid_dims.cpp - ecl/grid_dims.cpp - ecl/nnc_info.cpp - ecl/ecl_grav_common.cpp - ecl/nnc_vector.cpp - ecl/ecl_nnc_export.cpp - ecl/ecl_nnc_data.cpp - ecl/ecl_nnc_geometry.cpp - ecl/layer.cpp - ecl/fault_block.cpp - ecl/fault_block_layer.cpp - ecl/ecl_type.cpp - ecl/ecl_type_python.cpp - ecl/well_state.cpp - ecl/well_conn.cpp - ecl/well_info.cpp - ecl/well_ts.cpp - ecl/well_conn_collection.cpp - ecl/well_segment.cpp - ecl/well_segment_collection.cpp - ecl/well_branch_collection.cpp - ecl/well_rseg_loader.cpp - ecl/FortIO.cpp - ecl/EclFilename.cpp - geometry/geo_surface.cpp - geometry/geo_util.cpp - geometry/geo_pointset.cpp - geometry/geo_region.cpp - geometry/geo_polygon.cpp - geometry/geo_polygon_collection.cpp) +set(ECL_SOURCES + util/rng.cpp + util/lookup_table.cpp + util/statistics.cpp + util/mzran.cpp + util/hash_node.cpp + util/hash_sll.cpp + util/path.cpp + util/hash.cpp + util/node_data.cpp + util/node_ctype.cpp + util/util.cpp + util/util_symlink.cpp + util/util_lfs.c + util/util_unlink.cpp + util/vector.cpp + util/parser.cpp + util/stringlist.cpp + util/buffer.cpp + util/string_util.cpp + util/type_vector_functions.cpp + util/ecl_version.cpp + util/perm_vector.cpp + util/test_util.cpp + util/cxx_string_util.cpp + ${opt_srcs} + ecl/ecl_rsthead.cpp + ecl/ecl_sum_tstep.cpp + ecl/ecl_rst_file.cpp + ecl/ecl_init_file.cpp + ecl/ecl_grid_cache.cpp + ecl/smspec_node.cpp + ecl/ecl_kw_grdecl.cpp + ecl/ecl_file_kw.cpp + ecl/ecl_file_view.cpp + ecl/ecl_grav.cpp + ecl/ecl_grav_calc.cpp + ecl/ecl_smspec.cpp + ecl/ecl_unsmry_loader.cpp + ecl/ecl_sum_data.cpp + ecl/ecl_sum_file_data.cpp + ecl/ecl_util.cpp + ecl/ecl_kw.cpp + ecl/ecl_sum.cpp + ecl/ecl_sum_vector.cpp + ecl/fortio.c + ecl/ecl_rft_file.cpp + ecl/ecl_rft_node.cpp + ecl/ecl_rft_cell.cpp + ecl/ecl_grid.cpp + ecl/ecl_coarse_cell.cpp + ecl/ecl_box.cpp + ecl/ecl_io_config.cpp + ecl/ecl_file.cpp + ecl/ecl_region.cpp + ecl/ecl_subsidence.cpp + ecl/ecl_grid_dims.cpp + ecl/grid_dims.cpp + ecl/nnc_info.cpp + ecl/ecl_grav_common.cpp + ecl/nnc_vector.cpp + ecl/ecl_nnc_export.cpp + ecl/ecl_nnc_data.cpp + ecl/ecl_nnc_geometry.cpp + ecl/layer.cpp + ecl/fault_block.cpp + ecl/fault_block_layer.cpp + ecl/ecl_type.cpp + ecl/ecl_type_python.cpp + ecl/well_state.cpp + ecl/well_conn.cpp + ecl/well_info.cpp + ecl/well_ts.cpp + ecl/well_conn_collection.cpp + ecl/well_segment.cpp + ecl/well_segment_collection.cpp + ecl/well_branch_collection.cpp + ecl/well_rseg_loader.cpp + ecl/FortIO.cpp + ecl/EclFilename.cpp + geometry/geo_surface.cpp + geometry/geo_util.cpp + geometry/geo_pointset.cpp + geometry/geo_region.cpp + geometry/geo_polygon.cpp + geometry/geo_polygon_collection.cpp) + +add_library(ecl-static STATIC ${ECL_SOURCES}) +set_property(TARGET ecl-static PROPERTY POSITION_INDEPENDENT_CODE ON) + +add_library(ecl SHARED ${ECL_SOURCES}) if(ERT_WINDOWS) set_target_properties(ecl PROPERTIES PREFIX "lib") @@ -152,9 +156,20 @@ if(ERT_WINDOWS) endif() endif() +target_link_libraries(ecl-static PUBLIC ${m} ${dl} ${pthread} ${blas} ${zlib} + ${shlwapi}) target_link_libraries(ecl PUBLIC ${m} ${dl} ${pthread} ${blas} ${zlib} ${shlwapi}) +target_include_directories( + ecl-static + PUBLIC $ + $ + $ + PRIVATE ${ZLIB_INCLUDE_DIRS} util include + ${CMAKE_CURRENT_SOURCE_DIR}/private-include + ${CMAKE_CURRENT_BINARY_DIR}/include) + target_include_directories( ecl PUBLIC $ @@ -164,6 +179,15 @@ target_include_directories( ${CMAKE_CURRENT_SOURCE_DIR}/private-include ${CMAKE_CURRENT_BINARY_DIR}/include) +target_compile_definitions( + ecl-static + PRIVATE -DGIT_COMMIT=${GIT_COMMIT} + -DGIT_COMMIT_SHORT=${GIT_COMMIT_SHORT} + -DECL_VERSION_MAJOR=${ECL_VERSION_MAJOR} + -DECL_VERSION_MINOR=${ECL_VERSION_MINOR} + -DECL_VERSION_MICRO=${ECL_VERSION_MICRO} + $<$:HOST_BIG_ENDIAN>) + target_compile_definitions( ecl PRIVATE -DGIT_COMMIT=${GIT_COMMIT} @@ -184,19 +208,15 @@ if(ERT_USE_OPENMP) target_link_libraries(ecl PUBLIC ${OpenMP_EXE_LINKER_FLAGS}) endif() -set_target_properties( - ecl PROPERTIES VERSION ${ECL_VERSION_MAJOR}.${ECL_VERSION_MINOR} - SOVERSION ${ECL_VERSION_MAJOR}) +if(NOT SKBUILD) + set_target_properties( + ecl PROPERTIES VERSION ${ECL_VERSION_MAJOR}.${ECL_VERSION_MINOR} + SOVERSION ${ECL_VERSION_MAJOR}) +endif() -install( - TARGETS ecl - EXPORT ecl-config - ARCHIVE DESTINATION ${CMAKE_INSTALL_LIBDIR} - LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR} - RUNTIME DESTINATION ${CMAKE_INSTALL_BINDIR}) -install(DIRECTORY include/ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}) -install(DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/include/ - DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}) +if(SKBUILD) + add_subdirectory(python-binding) +endif() if(NOT BUILD_TESTS) return() @@ -363,7 +383,7 @@ foreach( well_segment_collection ecl_file) add_executable(${name} ecl/tests/${name}.cpp) - target_link_libraries(${name} ecl) + target_link_libraries(${name} ecl-static) target_include_directories( ${name} PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/private-include) add_test(NAME ${name} COMMAND ${name}) diff --git a/lib/ecl/ecl_smspec.cpp b/lib/ecl/ecl_smspec.cpp index 37b1f2280..df063a09e 100644 --- a/lib/ecl/ecl_smspec.cpp +++ b/lib/ecl/ecl_smspec.cpp @@ -46,6 +46,8 @@ #include #include +#include "detail/ecl/ecl_smspec.hpp" + #ifdef HAVE_FNMATCH #include #endif @@ -97,68 +99,8 @@ */ -#define ECL_SMSPEC_ID 806647 #define PARAMS_GLOBAL_DEFAULT -99 -typedef std::map node_map; - -struct ecl_smspec_struct { - UTIL_TYPE_ID_DECLARATION; - /* - All the hash tables listed below here are different ways to access - smspec_node instances. The actual smspec_node instances are - owned by the smspec_nodes vector; - */ - node_map field_var_index; - node_map misc_var_index; /* Variables like 'TCPU' and 'NEWTON'. */ - node_map - gen_var_index /* This is "everything" - things can either be found as gen_var("WWCT:OP_X") or as well_var("WWCT" , "OP_X") */ - ; - - std::map - well_var_index; /* Indexes for all well variables: - {well1: {var1: index1 , var2: index2} , well2: {var1: index1 , var2: index2}} */ - std::map - group_var_index; /* Indexes for group variables.*/ - std::map - region_var_index; /* The stored index is an offset. */ - std::map block_var_index; /* Block variables like BPR */ - std::map> - well_completion_var_index; /* Indexes for completion indexes .*/ - - std::vector> smspec_nodes; - bool write_mode; - bool need_nums; - std::vector index_map; - std::map inv_index_map; - int params_size; - - int time_seconds; - int grid_dims[3]; /* Grid dimensions - in DIMENS[1,2,3] */ - int num_regions; - int Nwells, param_offset; - std::string - key_join_string; /* The string used to join keys when building gen_key keys - typically ":" - - but arbitrary - NOT necessary to be able to invert the joining. */ - std::string - header_file; /* FULL path to the currenbtly loaded header_file. */ - - bool - formatted; /* Has this summary instance been loaded from a formatted (i.e. FSMSPEC file) or unformatted (i.e. SMSPEC) file. */ - time_t sim_start_time; /* When did the simulation start - worldtime. */ - - int time_index; /* The fields time_index, day_index, month_index and year_index */ - int day_index; /* are used by the ecl_sum_data object to locate per. timestep */ - int month_index; /* time information. */ - int year_index; - bool has_lgr; - std::vector params_default; - - std::string restart_case; - ert_ecl_unit_enum unit_system; - int restart_step; -}; - /** About indexing: --------------- diff --git a/lib/ecl/ecl_sum.cpp b/lib/ecl/ecl_sum.cpp index 414b265ae..400d36eaa 100644 --- a/lib/ecl/ecl_sum.cpp +++ b/lib/ecl/ecl_sum.cpp @@ -42,6 +42,7 @@ #include #include "detail/util/path.hpp" +#include "detail/ecl/ecl_sum.hpp" /** The ECLIPSE summary data is organised in a header file (.SMSPEC) @@ -87,27 +88,6 @@ ecl_sum_data.c files. */ -#define ECL_SUM_ID 89067 - -struct ecl_sum_struct { - UTIL_TYPE_ID_DECLARATION; - ecl_smspec_type *smspec; /* Internalized version of the SMSPEC file. */ - ecl_sum_data_type *data; /* The data - can be NULL. */ - ecl_sum_type *restart_case; - - bool fmt_case; - bool unified; - char *key_join_string; - char * - path; /* The path - as given for the case input. Can be NULL for cwd. */ - char *abs_path; /* Absolute path. */ - char *base; /* Only the basename. */ - char * - ecl_case; /* This is the current case, with optional path component. == path + base*/ - char * - ext; /* Only to support selective loading of formatted|unformatted and unified|multiple. (can be NULL) */ -}; - UTIL_SAFE_CAST_FUNCTION(ecl_sum, ECL_SUM_ID); UTIL_IS_INSTANCE_FUNCTION(ecl_sum, ECL_SUM_ID); diff --git a/lib/ecl/ecl_sum_data.cpp b/lib/ecl/ecl_sum_data.cpp index 76e6bfbdb..a2ec041bc 100644 --- a/lib/ecl/ecl_sum_data.cpp +++ b/lib/ecl/ecl_sum_data.cpp @@ -23,6 +23,7 @@ #include #include #include +#include #include #include @@ -42,6 +43,7 @@ #include #include +#include "detail/ecl/ecl_sum_data.hpp" #include "detail/ecl/ecl_sum_file_data.hpp" /* @@ -57,183 +59,6 @@ ecl_sum_file_data instance. */ -namespace { - -/* - The class CaseIndex and the struct IndexNode are used to maintain a list of - the ecl_sum_file_data instances, and lookup the correct one based one various - time related arguments. -*/ - -struct IndexNode { - IndexNode(int d, int o, int l) { - this->data_index = d; - this->offset = o; - this->length = l; - } - - int end() const { return this->offset + this->length; } - - int data_index; - int offset; - int length; - int report1; - int report2; - time_t time1; - time_t time2; - double days1; - double days2; - std::vector params_map; -}; - -class CaseIndex { -public: - IndexNode &add(int length) { - int offset = 0; - int data_index = this->index.size(); - - if (!this->index.empty()) - offset = this->index.back().end(); - - this->index.emplace_back(data_index, offset, length); - return this->index.back(); - } - - /* - The lookup_time() and lookup_report() methods will lookup which file_data - instance corresponds to the time/report argument. The methods will return two - pointers to file_data instances, if the argument is inside one file_data - instance the pointers will be equal - otherwise they will point to the - file_data instance before and after the argument: - - File 1 File 2 - |------|-----|------| |----|----------|---| - /|\ /|\ - | | - | | - A B - - For time A the lookup_time function will return whereas for time - B the function will return . - */ - - std::pair - lookup_time(time_t sim_time) const { - auto iter = this->index.begin(); - auto next = this->index.begin(); - if (sim_time < iter->time1) - throw std::invalid_argument("Simulation time out of range"); - - ++next; - while (true) { - double t1 = iter->time1; - double t2 = iter->time2; - - if (sim_time >= t1) { - if (sim_time <= t2) - return std::make_pair( - &(*iter), &(*iter)); - - if (next == this->index.end()) - throw std::invalid_argument("Simulation days out of range"); - - if (sim_time < next->time1) - return std::make_pair( - &(*iter), &(*next)); - } - ++next; - ++iter; - } - } - - std::pair - lookup_days(double days) const { - auto iter = this->index.begin(); - auto next = this->index.begin(); - if (days < iter->days1) - throw std::invalid_argument("Simulation days out of range"); - - ++next; - while (true) { - double d1 = iter->days1; - double d2 = iter->days2; - - if (days >= d1) { - if (days <= d2) - return std::make_pair( - &(*iter), &(*iter)); - - if (next == this->index.end()) - throw std::invalid_argument("Simulation days out of range"); - - if (days < next->days1) - return std::make_pair( - &(*iter), &(*next)); - } - ++next; - ++iter; - } - } - - const IndexNode &lookup(int internal_index) const { - for (const auto &node : this->index) - if (internal_index >= node.offset && internal_index < node.end()) - return node; - - throw std::invalid_argument("Internal error when looking up index: " + - std::to_string(internal_index)); - } - - const IndexNode &lookup_report(int report) const { - for (const auto &node : this->index) - if (node.report1 <= report && node.report2 >= report) - return node; - - throw std::invalid_argument("Internal error when looking up report: " + - std::to_string(report)); - } - - /* - This will check that we have a datafile which report range covers the - report argument, in adition there can be 'holes' in the series - that must - be checked by actually querying the data_file object. - */ - - bool has_report(int report) const { - for (const auto &node : this->index) - if (node.report1 <= report && node.report2 >= report) - return true; - - return false; - } - - IndexNode &back() { return this->index.back(); } - - void clear() { this->index.clear(); } - - int length() const { return this->index.back().end(); } - - std::vector::const_iterator begin() const { - return this->index.begin(); - } - - std::vector::const_iterator end() const { - return this->index.end(); - } - -private: - std::vector index; -}; - -} // namespace - -struct ecl_sum_data_struct { - const ecl_smspec_type *smspec; - std::vector - data_files; // List of ecl_sum_file_data instances - CaseIndex index; -}; - static void ecl_sum_data_build_index(ecl_sum_data_type *self); static double ecl_sum_data_iget_sim_seconds(const ecl_sum_data_type *data, int internal_index); @@ -824,11 +649,11 @@ static double ecl_sum_data_interp_get(const ecl_sum_data_type *data, ecl_sum_data_iget(data, time_index2, params_index) * weight2; } -static double ecl_sum_data_vector_iget(const ecl_sum_data_type *data, - time_t sim_time, int params_index, - bool is_rate, int time_index1, - int time_index2, double weight1, - double weight2) { +extern "C" double ecl_sum_data_vector_iget(const ecl_sum_data_type *data, + time_t sim_time, int params_index, + bool is_rate, int time_index1, + int time_index2, double weight1, + double weight2) { double value = 0.0; if (is_rate) { @@ -1071,10 +896,10 @@ ecl_sum_data_alloc_time_vector(const ecl_sum_data_type *data, return time_vector; } -static void ecl_sum_data_init_double_vector__(const ecl_sum_data_type *data, - int main_params_index, - double *output_data, - bool report_only) { +static void ecl_sum_data_init_float_vector__(const ecl_sum_data_type *data, + int main_params_index, + float *output_data, + bool report_only) { int offset = 0; for (const auto &index_node : data->index) { const auto &data_file = data->data_files[index_node.data_index]; @@ -1085,7 +910,7 @@ static void ecl_sum_data_init_double_vector__(const ecl_sum_data_type *data, const ecl::smspec_node &smspec_node = ecl_smspec_iget_node_w_params_index(data->smspec, main_params_index); - double default_value = smspec_node.get_default(); + float default_value = smspec_node.get_default(); offset += data_file->get_data_report(params_index, index_node.length, &output_data[offset], default_value); @@ -1114,13 +939,17 @@ void ecl_sum_data_init_datetime64_vector(const ecl_sum_data_type *data, void ecl_sum_data_init_double_vector(const ecl_sum_data_type *data, int params_index, double *output_data) { - ecl_sum_data_init_double_vector__(data, params_index, output_data, false); + auto len = data->index.length(); + auto output_data_f32 = std::make_unique(len); + ecl_sum_data_init_float_vector__(data, params_index, output_data_f32.get(), + false); + std::copy_n(output_data_f32.get(), len, output_data); } double_vector_type * ecl_sum_data_alloc_data_vector(const ecl_sum_data_type *data, int params_index, bool report_only) { - std::vector output_data; + std::vector output_data; if (report_only) output_data.resize(1 + ecl_sum_data_get_last_report_step(data) - ecl_sum_data_get_first_report_step(data)); @@ -1130,15 +959,13 @@ ecl_sum_data_alloc_data_vector(const ecl_sum_data_type *data, int params_index, if (params_index >= ecl_smspec_get_params_size(data->smspec)) throw std::out_of_range("Out of range"); - ecl_sum_data_init_double_vector__(data, params_index, output_data.data(), - report_only); + ecl_sum_data_init_float_vector__(data, params_index, output_data.data(), + report_only); double_vector_type *data_vector = double_vector_alloc(output_data.size(), 0); - { - double *tmp_data = double_vector_get_ptr(data_vector); - memcpy(tmp_data, output_data.data(), - output_data.size() * sizeof(double)); - } + + double *tmp_data = double_vector_get_ptr(data_vector); + std::copy(output_data.cbegin(), output_data.cend(), tmp_data); return data_vector; } diff --git a/lib/ecl/ecl_sum_file_data.cpp b/lib/ecl/ecl_sum_file_data.cpp index 8c816c06f..85781dca3 100644 --- a/lib/ecl/ecl_sum_file_data.cpp +++ b/lib/ecl/ecl_sum_file_data.cpp @@ -7,6 +7,7 @@ #include #include #include +#include #include "detail/ecl/ecl_sum_file_data.hpp" #include "detail/ecl/ecl_unsmry_loader.hpp" @@ -356,7 +357,7 @@ void ecl_sum_file_data::build_index() { int offset = ecl_smspec_get_first_step(this->ecl_smspec) - 1; std::vector report_steps = this->loader->report_steps(offset); std::vector sim_time = this->loader->sim_time(); - std::vector sim_seconds = this->loader->sim_seconds(); + std::vector sim_seconds = this->loader->sim_seconds(); for (int i = 0; i < this->loader->length(); i++) { this->index.add(sim_time[i], sim_seconds[i], report_steps[i]); @@ -395,10 +396,10 @@ int ecl_sum_file_data::get_time_report(int end_index, time_t *data) { return offset; } -void ecl_sum_file_data::get_data(int params_index, int length, double *data) { +void ecl_sum_file_data::get_data(int params_index, int length, float *data) { if (this->loader) { const auto tmp_data = loader->get_vector(params_index); - memcpy(data, tmp_data.data(), length * sizeof data); + std::copy_n(tmp_data.cbegin(), length, data); } else { for (int time_index = 0; time_index < length; time_index++) data[time_index] = this->iget(time_index, params_index); @@ -406,7 +407,7 @@ void ecl_sum_file_data::get_data(int params_index, int length, double *data) { } int ecl_sum_file_data::get_data_report(int params_index, int end_index, - double *data, double default_value) { + float *data, float default_value) { int offset = 0; for (int report_step = this->first_report(); diff --git a/lib/ecl/ecl_unsmry_loader.cpp b/lib/ecl/ecl_unsmry_loader.cpp index 5ccc99e86..29cf893e0 100644 --- a/lib/ecl/ecl_unsmry_loader.cpp +++ b/lib/ecl/ecl_unsmry_loader.cpp @@ -44,13 +44,13 @@ unsmry_loader::~unsmry_loader() { ecl_file_close(file); } int unsmry_loader::length() const { return this->m_length; } -std::vector unsmry_loader::get_vector(int pos) const { +std::vector unsmry_loader::get_vector(int pos) const { if (pos >= size) throw std::out_of_range( "unsmry_loader::get_vector pos: " + std::to_string(pos) + " PARAMS_SIZE: " + std::to_string(size)); - std::vector data(this->length()); + std::vector data(this->length()); int_vector_type *index_map = int_vector_alloc(1, pos); char buffer[4]; @@ -69,7 +69,7 @@ std::vector unsmry_loader::get_vector(int pos) const { } // This is horribly inefficient -double unsmry_loader::iget(int time_index, int params_index) const { +float unsmry_loader::iget(int time_index, int params_index) const { int_vector_type *index_map = int_vector_alloc(1, params_index); float value; ecl_file_view_index_fload_kw(this->file_view, PARAMS_KW, time_index, @@ -128,7 +128,7 @@ std::vector unsmry_loader::report_steps(int offset) const { std::vector unsmry_loader::sim_time() const { if (this->time_index >= 0) { - const std::vector sim_seconds = this->sim_seconds(); + const std::vector sim_seconds = this->sim_seconds(); std::vector st(this->length(), this->sim_start); for (size_t i = 0; i < st.size(); i++) @@ -150,16 +150,16 @@ std::vector unsmry_loader::sim_time() const { } } -std::vector unsmry_loader::sim_seconds() const { +std::vector unsmry_loader::sim_seconds() const { if (this->time_index >= 0) { - std::vector seconds = this->get_vector(this->time_index); + std::vector seconds = this->get_vector(this->time_index); for (size_t i = 0; i < seconds.size(); i++) seconds[i] *= this->time_seconds; return seconds; } else { std::vector st = this->sim_time(); - std::vector seconds(st.size()); + std::vector seconds(st.size()); for (size_t i = 0; i < st.size(); i++) seconds[i] = util_difftime_seconds(this->sim_start, st[i]); diff --git a/lib/ecl/tests/ecl_unsmry_loader_test.cpp b/lib/ecl/tests/ecl_unsmry_loader_test.cpp index 88b7a75f2..c2b60bc87 100644 --- a/lib/ecl/tests/ecl_unsmry_loader_test.cpp +++ b/lib/ecl/tests/ecl_unsmry_loader_test.cpp @@ -45,9 +45,9 @@ void test_load() { ecl::unsmry_loader *loader = new ecl::unsmry_loader(ecl_sum_get_smspec(ecl_sum), "CASE.UNSMRY", 0); - const std::vector FOPT_value = loader->get_vector(1); - const std::vector BPR_value = loader->get_vector(2); - const std::vector WWCT_value = loader->get_vector(3); + const std::vector FOPT_value = loader->get_vector(1); + const std::vector BPR_value = loader->get_vector(2); + const std::vector WWCT_value = loader->get_vector(3); test_assert_int_equal(FOPT_value.size(), 4); test_assert_double_equal(FOPT_value[3], 6.0); test_assert_double_equal(BPR_value[2], 10.0); diff --git a/lib/include/ert/ecl/ecl_util.hpp b/lib/include/ert/ecl/ecl_util.hpp index 9c08c652b..39b6358c6 100644 --- a/lib/include/ert/ecl/ecl_util.hpp +++ b/lib/include/ert/ecl/ecl_util.hpp @@ -28,6 +28,7 @@ extern "C" { #include #include #include +#include typedef enum { ECL_OTHER_FILE = 0, @@ -100,13 +101,6 @@ typedef enum { } #define ECL_PHASE_ENUM_SIZE 3 -typedef enum { - ECL_METRIC_UNITS = 1, - ECL_FIELD_UNITS = 2, - ECL_LAB_UNITS = 3, - ECL_PVT_M_UNITS = 4 -} ert_ecl_unit_enum; - // For unformatted files: #define ECL_BOOL_TRUE_INT \ -1 // Binary representation: 11111111 11111111 11111111 1111111 diff --git a/lib/include/ert/ecl/enums.h b/lib/include/ert/ecl/enums.h new file mode 100644 index 000000000..b8cea06b7 --- /dev/null +++ b/lib/include/ert/ecl/enums.h @@ -0,0 +1,8 @@ +#pragma once + +typedef enum { + ECL_METRIC_UNITS = 1, + ECL_FIELD_UNITS = 2, + ECL_LAB_UNITS = 3, + ECL_PVT_M_UNITS = 4 +} ert_ecl_unit_enum; diff --git a/lib/private-include/detail/ecl/ecl_smspec.hpp b/lib/private-include/detail/ecl/ecl_smspec.hpp new file mode 100644 index 000000000..9abfec5ab --- /dev/null +++ b/lib/private-include/detail/ecl/ecl_smspec.hpp @@ -0,0 +1,72 @@ +#pragma once + +#include +#include +#include +#include +#include + +#define ECL_SMSPEC_ID 806647 + +namespace ecl { +class smspec_node; +} + +typedef std::map node_map; + +typedef struct ecl_smspec_struct { + int __type_id; + /* + All the hash tables listed below here are different ways to access + smspec_node instances. The actual smspec_node instances are + owned by the smspec_nodes vector; + */ + node_map field_var_index; + node_map misc_var_index; /* Variables like 'TCPU' and 'NEWTON'. */ + node_map + gen_var_index /* This is "everything" - things can either be found as gen_var("WWCT:OP_X") or as well_var("WWCT" , "OP_X") */ + ; + + std::map + well_var_index; /* Indexes for all well variables: + {well1: {var1: index1 , var2: index2} , well2: {var1: index1 , var2: index2}} */ + std::map + group_var_index; /* Indexes for group variables.*/ + std::map + region_var_index; /* The stored index is an offset. */ + std::map block_var_index; /* Block variables like BPR */ + std::map> + well_completion_var_index; /* Indexes for completion indexes .*/ + + std::vector> smspec_nodes; + bool write_mode; + bool need_nums; + std::vector index_map; + std::map inv_index_map; + int params_size; + + int time_seconds; + int grid_dims[3]; /* Grid dimensions - in DIMENS[1,2,3] */ + int num_regions; + int Nwells, param_offset; + std::string + key_join_string; /* The string used to join keys when building gen_key keys - typically ":" - + but arbitrary - NOT necessary to be able to invert the joining. */ + std::string + header_file; /* FULL path to the currenbtly loaded header_file. */ + + bool + formatted; /* Has this summary instance been loaded from a formatted (i.e. FSMSPEC file) or unformatted (i.e. SMSPEC) file. */ + time_t sim_start_time; /* When did the simulation start - worldtime. */ + + int time_index; /* The fields time_index, day_index, month_index and year_index */ + int day_index; /* are used by the ecl_sum_data object to locate per. timestep */ + int month_index; /* time information. */ + int year_index; + bool has_lgr; + std::vector params_default; + + std::string restart_case; + ert_ecl_unit_enum unit_system; + int restart_step; +} ecl_smspec_type; diff --git a/lib/private-include/detail/ecl/ecl_sum.hpp b/lib/private-include/detail/ecl/ecl_sum.hpp new file mode 100644 index 000000000..fbb9c15a1 --- /dev/null +++ b/lib/private-include/detail/ecl/ecl_sum.hpp @@ -0,0 +1,26 @@ +#pragma + +typedef struct ecl_smspec_struct ecl_smspec_type; +typedef struct ecl_sum_data_struct ecl_sum_data_type; +typedef struct ecl_sum_struct ecl_sum_type; + +#define ECL_SUM_ID 89067 + +struct ecl_sum_struct { + int __type_id; + ecl_smspec_type *smspec; /* Internalized version of the SMSPEC file. */ + ecl_sum_data_type *data; /* The data - can be NULL. */ + ecl_sum_type *restart_case; + + bool fmt_case; + bool unified; + char *key_join_string; + char * + path; /* The path - as given for the case input. Can be NULL for cwd. */ + char *abs_path; /* Absolute path. */ + char *base; /* Only the basename. */ + char * + ecl_case; /* This is the current case, with optional path component. == path + base*/ + char * + ext; /* Only to support selective loading of formatted|unformatted and unified|multiple. (can be NULL) */ +}; diff --git a/lib/private-include/detail/ecl/ecl_sum_data.hpp b/lib/private-include/detail/ecl/ecl_sum_data.hpp new file mode 100644 index 000000000..444b53a92 --- /dev/null +++ b/lib/private-include/detail/ecl/ecl_sum_data.hpp @@ -0,0 +1,179 @@ +#pragma once + +#include + +typedef struct ecl_smspec_struct ecl_smspec_type; + +namespace detail { + +struct IndexNode { + IndexNode(int d, int o, int l) { + this->data_index = d; + this->offset = o; + this->length = l; + } + + int end() const { return this->offset + this->length; } + + int data_index; + int offset; + int length; + int report1; + int report2; + time_t time1; + time_t time2; + double days1; + double days2; + std::vector params_map; +}; + +class CaseIndex { +public: + IndexNode &add(int length) { + int offset = 0; + int data_index = this->index.size(); + + if (!this->index.empty()) + offset = this->index.back().end(); + + this->index.emplace_back(data_index, offset, length); + return this->index.back(); + } + + /* + The lookup_time() and lookup_report() methods will lookup which file_data + instance corresponds to the time/report argument. The methods will return two + pointers to file_data instances, if the argument is inside one file_data + instance the pointers will be equal - otherwise they will point to the + file_data instance before and after the argument: + + File 1 File 2 + |------|-----|------| |----|----------|---| + /|\ /|\ + | | + | | + A B + + For time A the lookup_time function will return whereas for time + B the function will return . + */ + + std::pair + lookup_time(time_t sim_time) const { + auto iter = this->index.begin(); + auto next = this->index.begin(); + if (sim_time < iter->time1) + throw std::invalid_argument("Simulation time out of range"); + + ++next; + while (true) { + double t1 = iter->time1; + double t2 = iter->time2; + + if (sim_time >= t1) { + if (sim_time <= t2) + return std::make_pair( + &(*iter), &(*iter)); + + if (next == this->index.end()) + throw std::invalid_argument("Simulation days out of range"); + + if (sim_time < next->time1) + return std::make_pair( + &(*iter), &(*next)); + } + ++next; + ++iter; + } + } + + std::pair + lookup_days(double days) const { + auto iter = this->index.begin(); + auto next = this->index.begin(); + if (days < iter->days1) + throw std::invalid_argument("Simulation days out of range"); + + ++next; + while (true) { + double d1 = iter->days1; + double d2 = iter->days2; + + if (days >= d1) { + if (days <= d2) + return std::make_pair( + &(*iter), &(*iter)); + + if (next == this->index.end()) + throw std::invalid_argument("Simulation days out of range"); + + if (days < next->days1) + return std::make_pair( + &(*iter), &(*next)); + } + ++next; + ++iter; + } + } + + const IndexNode &lookup(int internal_index) const { + for (const auto &node : this->index) + if (internal_index >= node.offset && internal_index < node.end()) + return node; + + throw std::invalid_argument("Internal error when looking up index: " + + std::to_string(internal_index)); + } + + const IndexNode &lookup_report(int report) const { + for (const auto &node : this->index) + if (node.report1 <= report && node.report2 >= report) + return node; + + throw std::invalid_argument("Internal error when looking up report: " + + std::to_string(report)); + } + + /* + This will check that we have a datafile which report range covers the + report argument, in adition there can be 'holes' in the series - that must + be checked by actually querying the data_file object. + */ + + bool has_report(int report) const { + for (const auto &node : this->index) + if (node.report1 <= report && node.report2 >= report) + return true; + + return false; + } + + IndexNode &back() { return this->index.back(); } + + void clear() { this->index.clear(); } + + size_t length() const { return this->index.back().end(); } + + std::vector::const_iterator begin() const { + return this->index.begin(); + } + + std::vector::const_iterator end() const { + return this->index.end(); + } + +private: + std::vector index; +}; + +} // namespace detail + +namespace ecl { +class ecl_sum_file_data; +} + +typedef struct ecl_sum_data_struct { + const ecl_smspec_type *smspec; + std::vector data_files; + detail::CaseIndex index; +} ecl_sum_data_type; diff --git a/lib/private-include/detail/ecl/ecl_sum_file_data.hpp b/lib/private-include/detail/ecl/ecl_sum_file_data.hpp index dfce4e41b..2eabe9cd6 100644 --- a/lib/private-include/detail/ecl/ecl_sum_file_data.hpp +++ b/lib/private-include/detail/ecl/ecl_sum_file_data.hpp @@ -2,16 +2,20 @@ #include #include -#include +#define INVALID_MINISTEP_NR -1 +#define INVALID_TIME_T 0 -#include -#include -#include +typedef struct ecl_file_struct ecl_file_type; +typedef struct ecl_file_view_struct ecl_file_view_type; +typedef struct ecl_smspec_struct ecl_smspec_type; +typedef struct fortio_struct fortio_type; +typedef struct vector_struct vector_type; +typedef struct stringlist_struct stringlist_type; +typedef struct ecl_sum_tstep_struct ecl_sum_tstep_type; namespace ecl { -#define INVALID_MINISTEP_NR -1 -#define INVALID_TIME_T 0 +class unsmry_loader; struct IndexNode { @@ -89,7 +93,7 @@ class ecl_sum_file_data { int length_before(time_t end_time) const; void get_time(int length, time_t *data); - void get_data(int params_index, int length, double *data); + void get_data(int params_index, int length, float *data); int length() const; time_t get_data_start() const; time_t get_sim_end() const; @@ -105,8 +109,8 @@ class ecl_sum_file_data { bool report_step_equal(const ecl_sum_file_data &other, bool strict) const; int report_before(time_t end_time) const; int get_time_report(int max_internal_index, time_t *data); - int get_data_report(int params_index, int max_internal_index, double *data, - double default_value); + int get_data_report(int params_index, int max_internal_index, float *data, + float default_value); int first_report() const; int last_report() const; int iget_report(int time_index) const; diff --git a/lib/private-include/detail/ecl/ecl_unsmry_loader.hpp b/lib/private-include/detail/ecl/ecl_unsmry_loader.hpp index b65f3e9c2..11f770ce2 100644 --- a/lib/private-include/detail/ecl/ecl_unsmry_loader.hpp +++ b/lib/private-include/detail/ecl/ecl_unsmry_loader.hpp @@ -14,15 +14,15 @@ class unsmry_loader { int file_options); ~unsmry_loader(); - std::vector get_vector(int pos) const; - std::vector sim_seconds() const; + std::vector get_vector(int pos) const; + std::vector sim_seconds() const; std::vector sim_time() const; int length() const; time_t iget_sim_time(int time_index) const; double iget_sim_seconds(int time_index) const; std::vector report_steps(int offset) const; - double iget(int time_index, int params_index) const; + float iget(int time_index, int params_index) const; private: int size; //Number of entries in the smspec index diff --git a/lib/python-binding/CMakeLists.txt b/lib/python-binding/CMakeLists.txt new file mode 100644 index 000000000..e4988c791 --- /dev/null +++ b/lib/python-binding/CMakeLists.txt @@ -0,0 +1,17 @@ +execute_process( + COMMAND "${PYTHON_EXECUTABLE}" -c + "import pybind11; print(pybind11.get_cmake_dir())" + OUTPUT_VARIABLE _tmp_dir + OUTPUT_STRIP_TRAILING_WHITESPACE COMMAND_ECHO STDOUT) +list(APPEND CMAKE_PREFIX_PATH "${_tmp_dir}") + +find_package(pybind11 CONFIG REQUIRED) + +pybind11_add_module(_native MODULE init.cpp summary.cpp) +target_link_libraries(_native PRIVATE ecl-static) +target_include_directories( + _native + PRIVATE "${CMAKE_SOURCE_DIR}/lib/include" + "${CMAKE_SOURCE_DIR}/lib/private-include" + "${CMAKE_BINARY_DIR}/lib/include") +target_compile_definitions(_native PRIVATE VERSION_INFO=${PROJECT_VERSION}) diff --git a/lib/python-binding/init.cpp b/lib/python-binding/init.cpp new file mode 100644 index 000000000..25b0515f6 --- /dev/null +++ b/lib/python-binding/init.cpp @@ -0,0 +1,6 @@ +#include +namespace py = pybind11; + +void init_native_summary(py::module_); + +PYBIND11_MODULE(_native, m) { init_native_summary(m.def_submodule("summary")); } diff --git a/lib/python-binding/summary.cpp b/lib/python-binding/summary.cpp new file mode 100644 index 000000000..f3b83021a --- /dev/null +++ b/lib/python-binding/summary.cpp @@ -0,0 +1,275 @@ +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include "ert/ecl/ecl_sum.hpp" +#include "ert/ecl/ecl_sum_data.hpp" +#include "ert/ecl/smspec_node.hpp" + +#include "detail/ecl/ecl_sum.hpp" +#include "detail/ecl/ecl_sum_data.hpp" +#include "detail/ecl/ecl_sum_file_data.hpp" +#include "detail/ecl/ecl_smspec.hpp" + +extern "C" double ecl_sum_data_vector_iget(const ecl_sum_data_type *data, + time_t sim_time, int params_index, + bool is_rate, int time_index1, + int time_index2, double weight1, + double weight2); + +void ecl_sum_data_init_interp_from_sim_time(const ecl_sum_data_type *data, + time_t sim_time, int *index1, + int *index2, double *weight1, + double *weight2); + +namespace py = pybind11; +using std::chrono::system_clock; + +using namespace pybind11::literals; + +using std::optional; +using std::string; +using std::vector; + +namespace ecl { +namespace data { + +void assert_ecl_type(const void *obj, int type_id) { + if (!(*reinterpret_cast(obj) == type_id)) { + throw py::type_error("Unexpected ecl type"); + } +} + +} // namespace data +} // namespace ecl + +namespace { +template const T *_get_object_from_cwrap(py::object &obj) { + static auto basecclass = py::module_::import("cwrap").attr("BaseCClass"); + if (!py::isinstance(obj, basecclass)) { + throw py::type_error("Expected a BaseCClass type"); + } + + auto address = py::cast(obj.attr("_address")()); + auto object = reinterpret_cast(address); + ecl::data::assert_ecl_type(object, ECL_SUM_ID); + return object; +} + +const ecl::smspec_node &_get_smspec_by_name(const ecl_sum_type *sum, + const std::string &key) { + const auto &index = sum->smspec->gen_var_index; + auto it = index.find(key.c_str()); + if (it == index.end()) + throw py::key_error("No such key"); + return *it->second; +} + +const ecl::smspec_node &_get_node_by_params_index(const ecl_smspec_type *smspec, + int params_index) { + auto node_index = smspec->inv_index_map.at(params_index); + return *smspec->smspec_nodes[node_index]; +} + +/** + * @brief Get an array of ecl::smspec_node by names + * + * @note This function runs in at least O(n^2) when keys is not none, which is + * unnecessarily slow. Wildcards are not supported. + * + * @param smspec[in] The main smspec object, containing nodes + * @param keys[in] An array of keywords which we match + * + * @return Array of ecl::smspec_nodes + */ +vector +_get_nodes_by_keys(const ecl_smspec_type *smspec, + const optional> &keys) { + std::vector nodes{}; + if (keys.has_value()) { + for (const auto &key : *keys) { + try { + nodes.push_back(smspec->gen_var_index.at(key)); + } catch (std::out_of_range &e) { + throw py::key_error(e.what()); + } + } + } else { + for (const auto &node : smspec->smspec_nodes) { + nodes.push_back(node.get()); + } + } + return nodes; +} + +} // namespace + +/** + * @brief Generate a numpy 1-d array with report steps only + * + * @param self ecl.summary.EclSum object + * @param keyword + */ +py::array report_to_numpy(py::object self, const string &keyword) { + const auto ecl_sum = _get_object_from_cwrap(self); + const auto data = ecl_sum->data; + const auto &smspec = _get_smspec_by_name(ecl_sum, keyword); + + size_t size = 1 + ecl_sum_data_get_last_report_step(data) - + ecl_sum_data_get_first_report_step(data); + + py::array_t array(size); + int offset = 0; + auto output_data = array.mutable_data(); + for (const auto &index_node : data->index) { + const auto &data_file = data->data_files[index_node.data_index]; + const auto ¶ms_map = index_node.params_map; + int params_index = params_map[smspec.get_params_index()]; + + const ecl::smspec_node &smspec_node = + _get_node_by_params_index(data->smspec, smspec.get_params_index()); + float default_value = smspec_node.get_default(); + offset += + data_file->get_data_report(params_index, index_node.length, + &output_data[offset], default_value); + } + + return std::move(array); +} + +/** + * @brief Generate a numpy array for the given list of keywords + * + * @param self ecl.summary.EclSum object + * @param keywords Array of keywords (wildcards not supported) + */ +py::array to_numpy(py::object self, const optional> &keywords) { + const auto ecl_sum = _get_object_from_cwrap(self); + const auto data = ecl_sum->data; + const auto smspec = ecl_sum->smspec; + + // Find all relevant smspec_nodes + const auto nodes = _get_nodes_by_keys(smspec, keywords); + + const ssize_t num_rows = ecl_sum_data_get_length(data); + const ssize_t num_cols = nodes.size(); + + py::array_t array{{num_rows, num_cols}}; + for (ssize_t time_index{}; time_index < num_rows; ++time_index) { + for (ssize_t key_index{}; key_index < num_cols; ++key_index) { + const auto node = nodes[key_index]; + const auto param_index = node->get_params_index(); + + array.mutable_at(time_index, key_index) = + ecl_sum_data_iget(data, time_index, param_index); + } + } + + return std::move(array); +} + +/** + * @brief ecl_sum_type to numpy matrix by interpolating + * + * Read data from ECLIPSE summary files into a numpy array. This array has + * keys.size() columns and time_points rows. If @p keys is + * None, then the columns are all keywords defined in the summary file. If a + * time point is not present in the summary file, linearly interpolate between + * its two nearest time points. + * + * @param self ecl.summary.EclSum object + * @param keywords Array of keywords (wildcards not supported) + * @param time_points Array of timestamps to which we'll interpolate + * + * @return @p keys x @p time_points numpy matrix + */ +py::array to_numpy_interp(py::object self, const optional> &keys, + const vector &time_points) { + const auto ecl_sum = _get_object_from_cwrap(self); + const auto data = ecl_sum->data; + const auto smspec = ecl_sum->smspec; + + // Find all relevant smspec_nodes + const auto nodes = _get_nodes_by_keys(smspec, keys); + + const auto start_time = ecl_sum_data_get_data_start(data); + const auto end_time = ecl_sum_data_get_sim_end(data); + + const ssize_t num_rows = time_points.size(); + const ssize_t num_cols = nodes.size(); + + py::array_t array{{num_rows, num_cols}}; + + auto inner_loop_clamp = [&](auto time_index, auto fn) { + for (ssize_t key_index{}; key_index < num_cols; ++key_index) { + const auto node = nodes[key_index]; + const auto param_index = node->get_params_index(); + + float value = 0.f; + if (!node->is_rate()) + value = fn(data, param_index); + array.mutable_at(time_index, key_index) = value; + } + }; + + auto inner_loop = [&](auto time_index, auto sim_time) { + double weight1, weight2; + int time_index1, time_index2; + + ecl_sum_data_init_interp_from_sim_time( + data, sim_time, &time_index1, &time_index2, &weight1, &weight2); + + for (ssize_t key_index{}; key_index < num_cols; ++key_index) { + const auto node = nodes[key_index]; + const auto param_index = node->get_params_index(); + + float value = ecl_sum_data_vector_iget( + data, sim_time, param_index, node->is_rate(), time_index1, + time_index2, weight1, weight2); + array.mutable_at(time_index, key_index) = value; + } + }; + + for (ssize_t time_index{}; time_index < num_rows; ++time_index) { + auto sim_time = system_clock::to_time_t(time_points[time_index]); + if (sim_time < start_time) { + inner_loop_clamp(time_index, &ecl_sum_data_iget_first_value); + } else if (sim_time > end_time) { + inner_loop_clamp(time_index, &ecl_sum_data_iget_last_value); + } else { + inner_loop(time_index, sim_time); + } + } + return std::move(array); +} + +/** + * @brief ecl._native.summary module set-up + * + * @param m Submodule object for ecl._native.summary + */ +void init_native_summary(py::module_ m) { + using keywords_t = optional>; + using time_points_t = optional>; + + m.def( + "to_numpy", + [](py::object self, const keywords_t &keywords, + const time_points_t &time_points) -> py::array { + if (time_points) + return to_numpy_interp(self, keywords, *time_points); + else + return to_numpy(self, keywords); + }, + py::arg{"self"}, py::arg{"keywords"}, py::arg{"time_indices"}); + m.def("report_to_numpy", &report_to_numpy, py::arg{"self"}, + py::arg{"keyword"}); +} diff --git a/lib/util/tests/ert_util_vector_test.cpp b/lib/util/tests/ert_util_vector_test.cpp index 43145fbdb..af5e63d51 100644 --- a/lib/util/tests/ert_util_vector_test.cpp +++ b/lib/util/tests/ert_util_vector_test.cpp @@ -18,6 +18,7 @@ #include #include +#include #include #include #include diff --git a/lib/util/type_vector_functions.cpp b/lib/util/type_vector_functions.cpp index 7d1d99e1a..b7219750f 100644 --- a/lib/util/type_vector_functions.cpp +++ b/lib/util/type_vector_functions.cpp @@ -18,6 +18,7 @@ #include +#include #include #include #include diff --git a/pyproject.toml b/pyproject.toml index 8280b3eee..66d21e009 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,2 +1,2 @@ [build-system] -requires = ["setuptools", "setuptools_scm", "wheel", "scikit-build", "cmake", "ninja"] +requires = ["setuptools", "setuptools_scm", "wheel", "scikit-build", "cmake", "ninja", "pybind11"] diff --git a/python/ecl/summary/ecl_sum.py b/python/ecl/summary/ecl_sum.py index b98649fff..dc6d0f6b8 100644 --- a/python/ecl/summary/ecl_sum.py +++ b/python/ecl/summary/ecl_sum.py @@ -21,6 +21,8 @@ libecl/src directory. """ +from typing import List, Optional +import sys import warnings import numpy import datetime @@ -43,6 +45,7 @@ from ecl.summary.ecl_sum_vector import EclSumVector from ecl.summary.ecl_smspec_node import EclSMSPECNode from ecl import EclPrototype, EclUnitTypeEnum +from ecl._native import summary as _native # , EclSumKeyWordVector @@ -91,12 +94,14 @@ def date2num(dt): class EclSum(BaseCClass): TYPE_NAME = "ecl_sum" + _fread_alloc_case2 = EclPrototype( "void* ecl_sum_fread_alloc_case2__(char*, char*, bool, bool, int)", bind=False, ) _fread_alloc = EclPrototype( - "void* ecl_sum_fread_alloc(char*, stringlist, char*, bool)", bind=False + "void* ecl_sum_fread_alloc(char*, stringlist, char*, bool, bool, int)", + bind=False, ) _create_restart_writer = EclPrototype( "ecl_sum_obj ecl_sum_alloc_restart_writer2(char*, char*, int, bool, bool, char*, time_t, bool, int, int, int)", @@ -276,7 +281,7 @@ def load(cls, smspec_file, unsmry_file, key_join_string=":", include_restart=Tru data_files = StringList() data_files.append(unsmry_file) c_ptr = cls._fread_alloc( - smspec_file, data_files, key_join_string, include_restart + smspec_file, data_files, key_join_string, include_restart, False, 0 ) if c_ptr is None: raise IOError("Failed to create summary instance") @@ -474,6 +479,17 @@ def _make_time_vector(self, time_index): time_points.append(t) return time_points + def to_numpy( + self, + keyword: str, + time_index: Optional[List[datetime.datetime]] = None, + report_only: bool = False, + ) -> numpy.ndarray: + if report_only: + return _native.report_to_numpy(self, keyword) + else: + return _native.to_numpy(self, [keyword], time_index) + def numpy_vector(self, key, time_index=None, report_only=False): """Will return numpy vector of all the values corresponding to @key. @@ -497,30 +513,10 @@ def numpy_vector(self, key, time_index=None, report_only=False): ValueError exception. """ - if key not in self: - raise KeyError("No such key:%s" % key) - - if report_only: - if time_index is None: - time_index = self.report_dates - else: - raise ValueError("Can not suuply both time_index and report_only=True") + if time_index is not None and report_only is True: + raise ValueError - if time_index is None: - np_vector = numpy.zeros(len(self)) - self._init_numpy_vector( - key, np_vector.ctypes.data_as(ctypes.POINTER(ctypes.c_double)) - ) - return np_vector - else: - time_vector = self._make_time_vector(time_index) - np_vector = numpy.zeros(len(time_vector)) - self._init_numpy_vector_interp( - key, - time_vector, - np_vector.ctypes.data_as(ctypes.POINTER(ctypes.c_double)), - ) - return np_vector + return self.to_numpy(key, time_index, report_only) # .astype(numpy.float64) @property def numpy_dates(self): @@ -583,35 +579,26 @@ def pandas_frame(self, time_index=None, column_keys=None): 2010-04-01 672.7 620.4 78.7 .... """ - from ecl.summary import EclSumKeyWordVector - if column_keys is None: - keywords = EclSumKeyWordVector(self, add_keywords=True) - else: - keywords = EclSumKeyWordVector(self) - for key in column_keys: - keywords.add_keywords(key) + column_keys = self.keys() - if len(keywords) == 0: + if len(column_keys) == 0: raise ValueError("No valid key") - if time_index is None: - time_index = self.numpy_dates - data = numpy.zeros([len(time_index), len(keywords)]) - EclSum._init_pandas_frame( - self, keywords, data.ctypes.data_as(ctypes.POINTER(ctypes.c_double)) - ) - else: - time_points = self._make_time_vector(time_index) - data = numpy.zeros([len(time_points), len(keywords)]) - EclSum._init_pandas_frame_interp( - self, - keywords, - time_points, - data.ctypes.data_as(ctypes.POINTER(ctypes.c_double)), - ) - - frame = pandas.DataFrame(index=time_index, columns=list(keywords), data=data) + try: + if time_index is None: + time_index = self.numpy_dates + data = _native.to_numpy(self, column_keys, None) + else: + data = _native.to_numpy(self, column_keys, time_index) + except KeyError as exc: + # For backwards-compatibility reasons, rethrow the KeyError as a + # ValueError + raise ValueError from exc + + frame = pandas.DataFrame( + index=time_index, columns=column_keys or self.keys(), data=data + ) return frame @staticmethod diff --git a/python/tests/ecl_tests/test_ecl_util.py b/python/tests/ecl_tests/test_ecl_util.py index 2f7cff40f..be671750c 100644 --- a/python/tests/ecl_tests/test_ecl_util.py +++ b/python/tests/ecl_tests/test_ecl_util.py @@ -24,6 +24,8 @@ def test_enums(self): source_file_path = "lib/include/ert/ecl/ecl_util.hpp" self.assertEnumIsFullyDefined(EclFileEnum, "ecl_file_enum", source_file_path) self.assertEnumIsFullyDefined(EclPhaseEnum, "ecl_phase_enum", source_file_path) + + source_file_path = "lib/include/ert/ecl/enums.h" self.assertEnumIsFullyDefined( EclUnitTypeEnum, "ert_ecl_unit_enum", source_file_path ) diff --git a/setup.py b/setup.py index c12880151..b2c24a26a 100644 --- a/setup.py +++ b/setup.py @@ -73,13 +73,10 @@ def utility_wrappers(): "-DECL_VERSION=" + version, "-DBUILD_APPLICATIONS=" + ("ON" if sys.platform == "linux" else "OFF"), "-DCMAKE_EXPORT_COMPILE_COMMANDS=ON", - "-DCMAKE_INSTALL_BINDIR=python/ecl/.bin", - "-DCMAKE_INSTALL_LIBDIR=python/ecl/.libs", - "-DCMAKE_INSTALL_INCLUDEDIR=python/ecl/.include", # we can safely pass OSX_DEPLOYMENT_TARGET as it's ignored on # everything not OS X. We depend on C++11, which makes our minimum # supported OS X release 10.9 - "-DCMAKE_OSX_DEPLOYMENT_TARGET=10.9", + "-DCMAKE_OSX_DEPLOYMENT_TARGET=10.13", ], # skbuild's test imples develop, which is pretty obnoxious instead, use a # manually integrated pytest.