diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 41056ea..15f55c7 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -8,6 +8,9 @@ list( APPEND MAIN_SOURCE_FILES opm/output/eclipse/Summary.cpp opm/output/eclipse/writeECLData.cpp opm/output/vtk/writeVtkData.cpp + opm/test_util/EclFilesComparator.cpp + opm/test_util/summaryRegressionTest.cpp + opm/test_util/summaryComparator.cpp ) list (APPEND PUBLIC_HEADER_FILES @@ -23,22 +26,27 @@ list (APPEND PUBLIC_HEADER_FILES opm/output/eclipse/Summary.hpp opm/output/eclipse/writeECLData.hpp opm/output/vtk/writeVtkData.hpp + opm/test_util/EclFilesComparator.hpp + opm/test_util/summaryRegressionTest.hpp + opm/test_util/summaryComparator.hpp ) list (APPEND TEST_SOURCE_FILES + tests/test_compareSummary.cpp + tests/test_EclFilesComparator.cpp tests/test_EclipseWriter.cpp - tests/test_RFT.cpp - tests/test_writenumwells.cpp tests/test_Restart.cpp - tests/test_Wells.cpp + tests/test_RFT.cpp tests/test_Summary.cpp -) + tests/test_Wells.cpp + tests/test_writenumwells.cpp + ) # originally generated with the command: # find tests -name '*.xml' -a ! -wholename '*/not-unit/*' -printf '\t%p\n' | sort list (APPEND TEST_DATA_FILES - tests/testBlackoilState3.DATA - tests/summary_deck.DATA tests/FIRST_SIM.DATA + tests/summary_deck.DATA + tests/testBlackoilState3.DATA tests/testRFT.DATA -) + ) diff --git a/opm/test_util/EclFilesComparator.cpp b/opm/test_util/EclFilesComparator.cpp new file mode 100644 index 0000000..c2905ed --- /dev/null +++ b/opm/test_util/EclFilesComparator.cpp @@ -0,0 +1,305 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include +#include + +#include +#include +#include +#include + +bool ECLFilesComparator::keywordValidForComparing(const std::string& keyword) const { + auto it = std::find(keywords1.begin(), keywords1.end(), keyword); + if (it == keywords1.end()) { + std::cout << "Could not find keyword: " << keyword << std::endl; + OPM_THROW(std::runtime_error, "Keyword does not exist in first file."); + } + it = find(keywords2.begin(), keywords2.end(), keyword); + if (it == keywords2.end()) { + std::cout << "Could not find keyword: " << keyword << std::endl; + OPM_THROW(std::runtime_error, "Keyword does not exist in second file."); + } + ecl_kw_type* ecl_kw = ecl_file_iget_named_kw(ecl_file1, keyword.c_str(), 0); + ecl_type_enum ecl_kw_type = ecl_kw_get_type(ecl_kw); + switch (ecl_kw_type) { + case ECL_DOUBLE_TYPE: + return true; + case ECL_FLOAT_TYPE: + return true; + case ECL_CHAR_TYPE: + std::cout << "\nKeyword " << keyword << " is of type ECL_CHAR_TYPE, " + << "which is not supported for comparison." << std::endl; + return false; + case ECL_INT_TYPE: + std::cout << "\nKeyword " << keyword << " is of type ECL_INT_TYPE, " + << "which is not supported for comparison." << std::endl; + return false; + case ECL_BOOL_TYPE: + std::cout << "\nKeyword " << keyword << " is of type ECL_BOOL_TYPE, " + << "which is not supported for comparison." << std::endl; + return false; + case ECL_MESS_TYPE: + std::cout << "\nKeyword " << keyword << " is of type ECL_MESS_TYPE, " + << "which is not supported for comparison." << std::endl; + return false; + default: + std::cout << "\nKeyword " << keyword << "has undefined type." << std::endl; + return false; + } +} + + + +void ECLFilesComparator::printResultsForKeyword(const std::string& keyword) const { + const double absDeviationAverage = average(absDeviation); + const double relDeviationAverage = average(relDeviation); + std::cout << "Average absolute deviation = " << absDeviationAverage << std::endl; + std::cout << "Median absolute deviation = " << median(absDeviation) << std::endl; + std::cout << "Average relative deviation = " << relDeviationAverage << std::endl; + std::cout << "Median relative deviation = " << median(relDeviation) << "\n\n"; +} + + + +void ECLFilesComparator::deviationsForOccurence(const std::string& keyword, int occurence) { + ecl_kw_type* ecl_kw1 = ecl_file_iget_named_kw(ecl_file1, keyword.c_str(), occurence); + ecl_kw_type* ecl_kw2 = ecl_file_iget_named_kw(ecl_file2, keyword.c_str(), occurence); + const unsigned int numCells1 = ecl_kw_get_size(ecl_kw1); + const unsigned int numCells2 = ecl_kw_get_size(ecl_kw2); + if (numCells1 != numCells2) { + std::cout << "For keyword " << keyword << " and occurence " << occurence << ":\n"; + OPM_THROW(std::runtime_error, "Number of active cells are different."); + } + std::vector values1(numCells1), values2(numCells2); + ecl_kw_get_data_as_double(ecl_kw1, values1.data()); + ecl_kw_get_data_as_double(ecl_kw2, values2.data()); + + auto it = std::find(keywordWhitelist.begin(), keywordWhitelist.end(), keyword); + for (unsigned int cell = 0; cell < values1.size(); cell++) { + deviationsForCell(values1[cell], values2[cell], keyword, occurence, cell, it == keywordWhitelist.end()); + } +} + + + +void ECLFilesComparator::deviationsForCell(double val1, double val2, const std::string& keyword, int occurence, int cell, bool allowNegativeValues) { + int i, j, k; + ecl_grid_get_ijk1(ecl_grid1, cell, &i, &j, &k); + // Coordinates from this function are zero-based, hence incrementing + i++, j++, k++; + + if (showValues) { + std::cout << "Occurence = " << occurence << std::endl; + std::cout << "Grid coordinate = (" << i << ", " << j << ", " << k << ")\n"; + std::cout << "(value1, value2) = (" << val1 << ", " << val2 << ")\n\n"; + } + if (!allowNegativeValues) { + if (val1 < 0) { + if (std::abs(val1) > absTolerance) { + std::cout << "For keyword " << keyword << " (which is listed for positive values only), occurence " << occurence + << ", and grid coordinate (" << i << ", " << j << ", " << k << "):\n" + << "the vector has a value of " << val1 << ", which exceeds the absolute tolerance of " << absTolerance << ".\n"; + OPM_THROW(std::runtime_error, "Negative value in first file."); + } + val1 = 0; + } + if (val2 < 0) { + if (std::abs(val2) > absTolerance) { + std::cout << "For keyword " << keyword << " (which is listed for positive values only), occurence " << occurence + << ", and grid coordinate (" << i << ", " << j << ", " << k << "):\n" + << "the vector has a value of " << val2 << ", which exceeds the absolute tolerance of " << absTolerance << ".\n"; + OPM_THROW(std::runtime_error, "Negative value in second file."); + } + val2 = 0; + } + } + Deviation dev = calculateDeviations(val1, val2); + if (dev.abs > absTolerance && dev.rel > relTolerance) { + std::cout << "For keyword " << keyword << ", occurence " << occurence + << ", and grid coordinate (" << i << ", " << j << ", " << k << "):\n" + << "(first value, second value) = (" << val1 << ", " << val2 << ")\n"; + std::cout << "The absolute deviation is " << dev.abs << ", and the tolerance limit is " << absTolerance << ".\n"; + std::cout << "The relative deviation is " << dev.rel << ", and the tolerance limit is " << relTolerance << ".\n"; + OPM_THROW(std::runtime_error, "Absolute deviation exceeds tolerance."); + } + if (dev.abs != -1) { + absDeviation.push_back(dev.abs); + } + if (dev.rel != -1) { + relDeviation.push_back(dev.rel); + } +} + + + +ECLFilesComparator::ECLFilesComparator(eclFileEnum fileType, const std::string& basename1, const std::string& basename2, double absTolerance, double relTolerance): + fileType(fileType), absTolerance(absTolerance), relTolerance(relTolerance) { + + std::string file1, file2; + if (fileType == RESTARTFILE) { + file1 = basename1 + ".UNRST"; + file2 = basename2 + ".UNRST"; + } + else if (fileType == INITFILE) { + file1 = basename1 + ".INIT"; + file2 = basename2 + ".INIT"; + } + + ecl_file1 = ecl_file_open(file1.c_str(), 0); + ecl_file2 = ecl_file_open(file2.c_str(), 0); + ecl_grid1 = ecl_grid_load_case(basename1.c_str()); + ecl_grid2 = ecl_grid_load_case(basename2.c_str()); + + if (ecl_file1 == nullptr) { + std::cout << "Tried to open: " << file1 << std::endl; + OPM_THROW(std::runtime_error, "Error opening first file."); + } + if (ecl_file2 == nullptr) { + std::cout << "Tried to open: " << file2 << std::endl; + OPM_THROW(std::runtime_error, "Error opening second file."); + } + if (ecl_grid1 == nullptr) { + std::cout << "Tried to open: " << basename1 + ".EGRID" << std::endl; + OPM_THROW(std::runtime_error, "Error opening first file."); + } + if (ecl_grid2 == nullptr) { + std::cout << "Tried to open: " << basename2 + ".EGRID" << std::endl; + OPM_THROW(std::runtime_error, "Error opening second file."); + } + + unsigned int numKeywords1 = ecl_file_get_num_distinct_kw(ecl_file1); + unsigned int numKeywords2 = ecl_file_get_num_distinct_kw(ecl_file2); + for (unsigned int i = 0; i < numKeywords1; ++i) { + std::string keyword(ecl_file_iget_distinct_kw(ecl_file1, i)); + keywords1.push_back(keyword); + } + for (unsigned int i = 0; i < numKeywords2; ++i) { + std::string keyword(ecl_file_iget_distinct_kw(ecl_file2, i)); + keywords2.push_back(keyword); + } +} + + + +ECLFilesComparator::~ECLFilesComparator() { + ecl_file_close(ecl_file1); + ecl_file_close(ecl_file2); + ecl_grid_free(ecl_grid1); + ecl_grid_free(ecl_grid2); +} + + + +void ECLFilesComparator::printKeywords() const { + std::cout << "\nKeywords for first file:\n"; + for (const auto& it : keywords1) std::cout << it << std::endl; + std::cout << "\nKeywords for second file:\n"; + for (const auto& it : keywords2) std::cout << it << std::endl; +} + + + +bool ECLFilesComparator::gridCompare() const { + const int globalGridCount1 = ecl_grid_get_global_size(ecl_grid1); + const int activeGridCount1 = ecl_grid_get_active_size(ecl_grid1); + const int globalGridCount2 = ecl_grid_get_global_size(ecl_grid2); + const int activeGridCount2 = ecl_grid_get_active_size(ecl_grid2); + if (globalGridCount1 != globalGridCount2) { + OPM_THROW(std::runtime_error, "Grid file: The total number of cells are different."); + return false; + } + if (activeGridCount1 != activeGridCount2) { + OPM_THROW(std::runtime_error, "Grid file: The number of active cells are different."); + return false; + } + return true; +} + + + +void ECLFilesComparator::results() { + if (keywords1.size() != keywords2.size()) { + OPM_THROW(std::runtime_error, "Number of keywords are not equal."); + } + for (const auto& it : keywords1) resultsForKeyword(it); +} + + + +void ECLFilesComparator::resultsForKeyword(const std::string keyword) { + std::cout << "\nKeyword " << keyword << ":\n\n"; + const unsigned int occurrences1 = ecl_file_get_num_named_kw(ecl_file1, keyword.c_str()); + const unsigned int occurrences2 = ecl_file_get_num_named_kw(ecl_file2, keyword.c_str()); + if (occurrences1 != occurrences2) { + std::cout << "For keyword " << keyword << ":\n"; + OPM_THROW(std::runtime_error, "Number of keyword occurrences are not equal."); + } + if (!keywordValidForComparing(keyword)) { + return; + } + for (unsigned int occurence = 0; occurence < occurrences1; ++occurence) { + deviationsForOccurence(keyword, occurence); + } + printResultsForKeyword(keyword); + absDeviation.clear(); + relDeviation.clear(); +} + + + +Deviation ECLFilesComparator::calculateDeviations(double val1, double val2) { + val1 = std::abs(val1); + val2 = std::abs(val2); + Deviation deviation; + if (val1 != 0 || val2 != 0) { + deviation.abs = std::abs(val1 - val2); + if (val1 != 0 && val2 != 0) { + deviation.rel = deviation.abs/(std::max(val1, val2)); + } + } + return deviation; +} + + + +double ECLFilesComparator::median(std::vector vec) { + if (vec.empty()) { + return 0; + } + else { + size_t n = vec.size()/2; + nth_element(vec.begin(), vec.begin() + n, vec.end()); + if (vec.size() % 2 == 0) { + return 0.5*(vec[n-1]+vec[n]); + } + else { + return vec[n]; + } + } +} + + + +double ECLFilesComparator::average(const std::vector& vec) { + if (vec.empty()) { + return 0; + } + double sum = std::accumulate(vec.begin(), vec.end(), 0.0); + return sum/vec.size(); +} diff --git a/opm/test_util/EclFilesComparator.hpp b/opm/test_util/EclFilesComparator.hpp new file mode 100644 index 0000000..75c75e6 --- /dev/null +++ b/opm/test_util/EclFilesComparator.hpp @@ -0,0 +1,111 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#ifndef ECLFILESCOMPARATOR_HPP +#define ECLFILESCOMPARATOR_HPP + +#include +#include + +/// Prototype for eclipse file struct, defined in the ERT library +struct ecl_file_struct; +typedef struct ecl_file_struct ecl_file_type; +/// Prototype for eclipse grid struct, defined in the ERT library +struct ecl_grid_struct; +typedef struct ecl_grid_struct ecl_grid_type; + +enum eclFileEnum {RESTARTFILE, INITFILE}; + + + +//! Deviation struct. +/*! The member variables are default initialized to -1, which is an invalid deviation value.*/ +struct Deviation { + double abs = -1; + double rel = -1; +}; + + + +//! A class for handeling ECLIPSE files. +/*! ECLFilesComparator opens ECLIPSE files (.UNRST or .INIT in addition to .EGRID or .GRID) from two simulations. + * Then the functions gridCompare, results and resultsForKeyword can be invoked to compare the files from the simulations. + * Note that when either absolute or relative deviation is larger than the member variables (double) absTolerance and (double) relTolerance respectively, + * an exception is thrown. This is also the case if a value is negative and its absolute value is larger than absTolerance. */ +class ECLFilesComparator { + private: + eclFileEnum fileType; + ecl_file_type* ecl_file1 = nullptr; + ecl_grid_type* ecl_grid1 = nullptr; + ecl_file_type* ecl_file2 = nullptr; + ecl_grid_type* ecl_grid2 = nullptr; + double absTolerance = 0; + double relTolerance = 0; + bool showValues = false; + std::vector keywords1, keywords2; + std::vector absDeviation, relDeviation; + + // Keywords which should not contain negative values, i.e. uses allowNegativeValues = false in deviationsForCell(): + const std::vector keywordWhitelist = {"SGAS", "SWAT", "PRESSURE"}; + + bool keywordValidForComparing(const std::string& keyword) const; + void printResultsForKeyword(const std::string& keyword) const; + void deviationsForOccurence(const std::string& keyword, int occurence); + void deviationsForCell(double val1, double val2, const std::string& keyword, int occurence, int cell, bool allowNegativeValues = true); + public: + //! \brief Open ECLIPSE files and set tolerances and keywords. + //! \param[in] fileType Specifies which filetype to be compared, possible inputs are "RESTARTFILE" and "INITFILE". + //! \param[in] basename1 Full path without file extension to the first case. + //! \param[in] basename2 Full path without file extension to the second case. + //! \param[in] absTolerance Tolerance for absolute deviation. + //! \param[in] relTolerance Tolerance for relative deviation. + //! \details The content of the ECLIPSE files specified in the input is stored in the ecl_file_type and ecl_grid_type member variables. In addition the keywords and absolute and relative tolerances (member variables) are set. If the constructor is unable to open one of the ECLIPSE files, an exception will be thrown. + ECLFilesComparator(eclFileEnum fileType, const std::string& basename1, const std::string& basename2, double absTolerance, double relTolerance); + //! \brief Closing the ECLIPSE files. + ~ECLFilesComparator(); + + //! \brief Setting the member variable showValues. + //! \param[in] showValues Variable which is copied to #showValues + //! \details When showValues is set to true, resultsForKeyword() will print all values side by side from the two files. + void setShowValues(bool showValues) {this->showValues = showValues;} + //! \brief Print all keywords stored in #keywords1 and #keywords2 + void printKeywords() const; + //! \brief Comparing grids from the to gridfiles. + //! \details This functions compares the number of active and global cells of the two grids. If the number of cells differ, an exception is thrown. + bool gridCompare() const; + //! \brief Calculating deviations for all keywords. + //! \details results() loops over all keywords, checks if they are supported for comparison, and calls resultsForKeyword() for each keyword. + void results(); + //! \brief Calculating deviations for a specific keyword. + //! \param[in] keyword Keyword which should be compared. + //! \details This function loops through every report step and every cell and compares the values for the given keyword from the two input cases. If the absolute or relative deviation between the two values for each step exceeds #absTolerance or #relTolerance respectively, or a value is both negative and has a larger absolute value than #absTolerance, an exception is thrown. The function also temporarily stores the absolute and relative deviation in #absDeviation and #relDeviation, and prints average and median values. + void resultsForKeyword(const std::string keyword); + + //! \brief Calculate deviations for two values. + //! \details Using absolute values of the input arguments: If one of the values are non-zero, the Deviation::abs returned is the difference between the two input values. In addition, if both values are non-zero, the Deviation::rel returned is the absolute deviation divided by the largest value. + static Deviation calculateDeviations(double val1, double val2); + //! \brief Calculate median of a vector. + //! \details Returning the median of the input vector, i.e. the middle value of the sorted vector if the number of elements is odd or the mean of the two middle values if the number of elements are even. + static double median(std::vector vec); + //! \brief Calculate average of a vector. + //! \details Returning the average of the input vector, i.e. the sum of all values divided by the number of elements. + static double average(const std::vector& vec); +}; + +#endif diff --git a/opm/test_util/summaryComparator.cpp b/opm/test_util/summaryComparator.cpp new file mode 100644 index 0000000..27672b2 --- /dev/null +++ b/opm/test_util/summaryComparator.cpp @@ -0,0 +1,210 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include +#include +#include +#include +#include +#include + +SummaryComparator::SummaryComparator(const char* basename1, const char* basename2, double absoluteTolerance, double relativeTolerance){ + ecl_sum1 = ecl_sum_fread_alloc_case(basename1, ":"); + ecl_sum2 = ecl_sum_fread_alloc_case(basename2, ":"); + if (ecl_sum1 == nullptr || ecl_sum2 == nullptr) { + OPM_THROW(std::runtime_error, "Not able to open files"); + } + this->absoluteTolerance = absoluteTolerance; + this->relativeTolerance = relativeTolerance; + keys1 = stringlist_alloc_new(); + keys2 = stringlist_alloc_new(); + ecl_sum_select_matching_general_var_list( ecl_sum1 , "*" , this->keys1); + stringlist_sort(this->keys1 , nullptr ); + ecl_sum_select_matching_general_var_list( ecl_sum2 , "*" , this->keys2); + stringlist_sort(this->keys2 , nullptr ); + + if(stringlist_get_size(keys1) <= stringlist_get_size(keys2)){ + this->keysShort = this->keys1; + this->keysLong = this->keys2; + }else{ + this->keysShort = this->keys2; + this->keysLong = this->keys1; + } +} + + + +SummaryComparator::~SummaryComparator(){ + ecl_sum_free(ecl_sum1); + ecl_sum_free(ecl_sum2); + stringlist_free(keys1); + stringlist_free(keys2); +} + + + +Deviation SummaryComparator::calculateDeviations(double val1, double val2){ + double absDev, relDev; + Deviation deviation; + absDev = std::abs(val1 - val2); + deviation.abs = absDev; + if (val1 != 0 || val2 != 0) { + relDev = absDev/double(std::max(std::abs(val1), std::abs(val2))); + deviation.rel = relDev; + } + return deviation; +} + + + +void SummaryComparator::setTimeVecs(std::vector &timeVec1,std::vector &timeVec2){ + timeVec1.reserve(ecl_sum_get_data_length(ecl_sum1)); + for (int time_index = 0; time_index < ecl_sum_get_data_length(ecl_sum1); time_index++){ + timeVec1.push_back(ecl_sum_iget_sim_days(ecl_sum1 , time_index )); + } + timeVec2.reserve(ecl_sum_get_data_length(ecl_sum2)); + for (int time_index = 0; time_index < ecl_sum_get_data_length(ecl_sum2); time_index++){ + timeVec2.push_back(ecl_sum_iget_sim_days(ecl_sum2 , time_index )); + } +} + + + +void SummaryComparator::getDataVecs(std::vector &dataVec1,std::vector &dataVec2,const char* keyword){ + dataVec1.reserve(ecl_sum_get_data_length(ecl_sum1)); + for (int time_index = 0; time_index < ecl_sum_get_data_length(ecl_sum1); time_index++){ + dataVec1.push_back(ecl_sum_iget(ecl_sum1, time_index, ecl_sum_get_general_var_params_index( ecl_sum1 , keyword ))); + } + dataVec2.reserve(ecl_sum_get_data_length(ecl_sum2)); + for (int time_index = 0; time_index < ecl_sum_get_data_length(ecl_sum2); time_index++){ + + dataVec2.push_back(ecl_sum_iget(ecl_sum2, time_index, ecl_sum_get_general_var_params_index( ecl_sum2 , keyword ))); + } +} + + + +void SummaryComparator::setDataSets(std::vector &timeVec1,std::vector &timeVec2){ + if(timeVec1.size()< timeVec2.size()){ + ecl_sum_fileShort = this->ecl_sum1; + ecl_sum_fileLong = this->ecl_sum2; + } + else{ + ecl_sum_fileShort = this->ecl_sum2; + ecl_sum_fileLong = this->ecl_sum1; + } +} + + + +void SummaryComparator::chooseReference(std::vector &timeVec1,std::vector &timeVec2,std::vector &dataVec1,std::vector &dataVec2){ + if(timeVec1.size() <= timeVec2.size()){ + referenceVec = &timeVec1; // time vector + referenceDataVec = &dataVec1; //data vector + checkVec = &timeVec2; + checkDataVec = &dataVec2; + } + else{ + referenceVec = &timeVec2; + referenceDataVec = &dataVec2; + checkVec = &timeVec1; + checkDataVec = &dataVec1; + } +} + + + +void SummaryComparator::getDeviation(size_t refIndex, size_t &checkIndex, Deviation &dev){ + if((*referenceVec)[refIndex] == (*checkVec)[checkIndex]){ + dev = SummaryComparator::calculateDeviations((*referenceDataVec)[refIndex], (*checkDataVec)[checkIndex]); + checkIndex++; + return; + } + else if((*referenceVec)[refIndex]<(*checkVec)[checkIndex]){ + double value = SummaryComparator::unitStep((*checkDataVec)[checkIndex-1]); + dev = SummaryComparator::calculateDeviations((*referenceDataVec)[refIndex], value); + checkIndex++; + return; + } + else{ + checkIndex++; + getDeviation(refIndex, checkIndex , dev); + } + if(checkIndex == checkVec->size() -1 ){ + return; + } +} + + + +double SummaryComparator::average(std::vector &vec){ + return std::accumulate(vec.begin(), vec.end(), 0.0) / vec.size(); +} + + + +double SummaryComparator::interpolation(double check_value, double check_value_prev , double time_array[3]){ + //does a Linear Polation + double time_check = time_array[0]; //From the check time vector + double time_check_prev = time_array[1];//From the check time vector + double time_reference = time_array[2];// The time of interest, from the reference vector + double sloap, factor, lp_value; + sloap = (check_value - check_value_prev)/double(time_check - time_check_prev); + factor = (time_reference - time_check_prev)/double(time_check - time_check_prev); + lp_value = check_value_prev + factor*sloap*(time_check - time_check_prev); + return lp_value; +} + + + +double SummaryComparator::median(std::vector vec) { + // Sorts and returns the median in a std::vector + if(vec.empty()){ + return 0; + } + else { + std::sort(vec.begin(), vec.end()); + if(vec.size() % 2 == 0) + return (vec[vec.size()/2 - 1] + vec[vec.size()/2]) / 2; + else + return vec[vec.size()/2]; + } +} + + + +void SummaryComparator::printUnits(){ + std::vector timeVec1, timeVec2; + setTimeVecs(timeVec1, timeVec2); // Sets the time vectors, they are equal for all keywords (WPOR:PROD01 etc) + setDataSets(timeVec1, timeVec2); + for (int jvar = 0; jvar < stringlist_get_size(keysLong); jvar++){ + std::cout << stringlist_iget(keysLong, jvar) << " unit: " << ecl_sum_get_unit(ecl_sum_fileShort, stringlist_iget(keysLong, jvar)) << std::endl; + } +} + + + +const char* SummaryComparator::getUnit(const char* keyword){ + return ecl_sum_get_unit(ecl_sum_fileShort, keyword); //Called only when the keywords are equal in the getDeviations()-function +} + + + +double SummaryComparator::unitStep(double value){ + return value; +} diff --git a/opm/test_util/summaryComparator.hpp b/opm/test_util/summaryComparator.hpp new file mode 100644 index 0000000..ac16402 --- /dev/null +++ b/opm/test_util/summaryComparator.hpp @@ -0,0 +1,140 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include +#include +#include +#include + +//! \brief Prototyping struct, encapsuling the stringlist libraries. +struct stringlist_struct; +typedef struct stringlist_struct stringlist_type; +//! \brief Prototyping struct, encapsuling the ert libraries. +struct ecl_sum_struct; +typedef struct ecl_sum_struct ecl_sum_type; + + + + +//! \brief Struct for storing the deviation between two values. +struct Deviation { + double abs = 0;//!< The absolute deviation + double rel = 0; //!< The relative deviation +}; + + + + +class SummaryComparator { + private: + double absoluteTolerance = 0;//!< The maximum absolute deviation that is allowed between two values. + double relativeTolerance = 0;//!< The maximum relative deviation that is allowed between twi values. + protected: + ecl_sum_type * ecl_sum1 = nullptr; //!< Struct that contains file1 + ecl_sum_type * ecl_sum2 = nullptr; //!< Struct that contains file2 + ecl_sum_type * ecl_sum_fileShort = nullptr; //!< For keeping track of the file with most/fewest timesteps + ecl_sum_type * ecl_sum_fileLong = nullptr; //!< For keeping track of the file with most/fewest timesteps + stringlist_type* keys1 = nullptr; //!< For storing all the keywords of file1 + stringlist_type* keys2 = nullptr; //!< For storing all the keywords of file2 + stringlist_type * keysShort = nullptr; //!< For keeping track of the file with most/fewest keywords + stringlist_type * keysLong = nullptr; //!< For keeping track of the file with most/fewest keywords + std::vector * referenceVec = nullptr; //!< For storing the values of each time step for the file containing the fewer time steps. + std::vector * referenceDataVec = nullptr; //!< For storing the data corresponding to each time step for the file containing the fewer time steps. + std::vector * checkVec = nullptr; //!< For storing the values of each time step for the file containing the more time steps. + std::vector * checkDataVec = nullptr; //!< For storing the data values corresponding to each time step for the file containing the more time steps. + + //! \brief Calculate deviation between two data values and stores it in a Deviation struct. + //! \param[in] refIndex Index in reference data + //! \param[in] checkindex Index in data to be checked. + //! \param[out] dev Holds the result from the comparison on return. + //! \details Uses the #referenceVec as basis, and checks its values against the values in #checkDataVec. The function is reccursive, and will update the iterative index j of the #checkVec until #checkVec[j] >= #referenceVec[i]. \n When #referenceVec and #checkVec have the same time value (i.e. #referenceVec[i] == #checkVec[j]) a direct comparison is used, \n when this is not the case, when #referenceVec[i] do not excist as an element in #checkVec, a value is generated, either by the principle of unit step or by interpolation. + void getDeviation(size_t refIndex, size_t &checkIndex, Deviation &dev); + + //! \brief Figure out which data file contains the most / less timesteps and assign member variable pointers accordingly. + //! \param[in] timeVec1 Data from first file + //! \param[in] timeVec2 Data from second file + //! \details Figure out which data file that cointains the more/fewer time steps and assigns the private member variable pointers #ecl_sum_fileShort / #ecl_sum_fileLong to the correct data sets #ecl_sum1 / #ecl_sum2. + void setDataSets(std::vector &timeVec1,std::vector &timeVec2); + + //! \brief Reads in the time values of each time step. + //! \param[in] timeVec1 Vector for storing the time steps from file1. + //! \param[in] timeVec2 Vector for storing the time steps from file2. + void setTimeVecs(std::vector &timeVec1,std::vector &timeVec2); + + //! \brief Read the data for one specific keyword into two separate vectors. + //! \param[in] dataVec1 Vector for storing the data for one specific keyword from file1. + //! \param[in] dataVec2 Vector for storing the data for one specific keyword from file2. + //! \details The two data files do not necessarily have the same amount of data values, but the values must correspond to the same interval in time. Thus possible to interpolate values. + void getDataVecs(std::vector &dataVec1,std::vector &dataVec2, const char* keyword); + + //! \brief Sets one data set as a basis and the other as values to check against. + //! \param[in] timeVec1 Used to figure out which dataset that have the more/fewer time steps. + //! \param[in] timeVec2 Used to figure out which dataset that have the more/fewer time steps. + //! \param[in] dataVec1 For assiging the the correct pointer to the data vector. + //! \param[in] dataVec2 For assiging the the correct pointer to the data vector. + //! \details Figures out which time vector that contains the fewer elements. Sets this as #referenceVec and its corresponding data as #referenceDataVec. \n The remaining data set is set as #checkVec (the time vector) and #checkDataVec. + void chooseReference(std::vector &timeVec1,std::vector &timeVec2,std::vector &dataVec1,std::vector &dataVec2); + + //! \brief Returns the relative tolerance. + double getRelTolerance(){return this->relativeTolerance;} + + //! \brief Returns the absolute tolerance. + double getAbsTolerance(){return this->absoluteTolerance;} + + //! \brief Returns the unit of the values of a keyword + //! \param[in] keyword The keyword of interest. + //! \param[out] ret The unit of the keyword as a const char*. + const char* getUnit(const char* keyword); + void printUnits(); + + public: + //! \brief Creates an SummaryComparator class object + //! \param[in] basename1 Path to file1 without extension. + //! \param[in] basename1 Path to file2 without extension. + //! \param[in] absoluteTolerance The absolute tolerance which is to be used in the test. + //! \param[in] relativeTolerance The relative tolerance which is to be used in the test. + //! \details The constructor creates an object of the class, and openes the files, an exception is thrown if the opening of the files fails. \n It creates stringlists, in which keywords are to be stored, and figures out which keylist that contains the more/less keywords. \n Also the private member variables aboluteTolerance and relativeTolerance are set. + SummaryComparator(const char* basename1, const char* basename2, double absoluteTolerance, double relativeTolerance); + //! \brief Destructor + //! \details The destructor takes care of the allocated memory in which data has been stored. + ~SummaryComparator(); + + //! \brief Calculates the deviation between two values + //! \param[in] val1 The first value of interest. + //! \param[in] val2 The second value if interest. + //! \param[out] ret Returns a Deviation struct. + //! \details The function takes two values, calculates the absolute and relative deviation and returns the result as a Deviation struct. + static Deviation calculateDeviations( double val1, double val2); + static double median(std::vector vec); + + //! \brief Function which linearly interpolates a value + //! \param[in] checkValue + //! \param[in] checkValuePrev + //! \param[in] timeArray[3] Array of doubles which contains the times of interest. \n timeArray[0]: The value of the time step in which checkVec exceeds the reference time, i.e. the first j which satisfies checkVec[j] > referenceVec[i].\n timeArray[1]: The value of the time step before checkVec exceeds referenceVec[i], i.e. checkVec[j-1].\n timeArray[2]: The value of the time step at which the comparison is to be comitted, i.e. referenceVec[i]. + //! \param[out] ret Returns the interpolated value. + //! \details In this case: The interpolation generates a value from the principle of linear polation. + static double interpolation(double checkValue, double checkValuePrev, double timeArray[3]); + + //! \brief Unit step function + //! \param[in] value The input value should be the last know value + //! \param[out] ret Return the unit-step-function value. + //! \details In this case: The unit step function is used when the data from the two data set, which is to be compared, don't match in time. \n The unit step function is then to be called on the #checkDataVec 's value at the last time step which is before the time of comparison. + static double unitStep(double value);//!< Returns a value based on the unit step principle. + static double average(std::vector &vec); +}; diff --git a/opm/test_util/summaryRegressionTest.cpp b/opm/test_util/summaryRegressionTest.cpp new file mode 100644 index 0000000..0d2b4d0 --- /dev/null +++ b/opm/test_util/summaryRegressionTest.cpp @@ -0,0 +1,101 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include +#include +#include +#include + +void RegressionTest::getRegressionTest(){ + std::vector timeVec1, timeVec2; + setTimeVecs(timeVec1, timeVec2); // Sets the time vectors, they are equal for all keywords (WPOR:PROD01 etc) + setDataSets(timeVec1, timeVec2); //Figures which dataset that contains more/less values pr keyword vector. + std::vector dataVec1, dataVec2, absDevVec, relDevVec; + int ivar = 0; + if(stringlist_get_size(keysShort) != stringlist_get_size(keysLong)){ + OPM_THROW(std::invalid_argument, "Different ammont of keywords in the two summary files."); + } + + //Iterates over all keywords from the restricted file, use iterator "ivar". Searches for a match in the file with more keywords, use the iterator "jvar". + while(ivar < stringlist_get_size(keysShort)){ + const char* keyword = stringlist_iget(keysShort, ivar); + for (int jvar = 0; jvar < stringlist_get_size(keysLong); jvar++){ + if (strcmp(keyword, stringlist_iget(keysLong, jvar)) == 0){ //When the keywords are equal, proceed in comparing summary files. + checkForKeyword(timeVec1, timeVec2, keyword); + break; + } + //will only enter here if no keyword match + if(jvar == stringlist_get_size(keysLong)-1){ + std::cout << "Could not find keyword: " << stringlist_iget(keysShort, ivar) << std::endl; + OPM_THROW(std::invalid_argument, "No match on keyword"); + } + } + ivar++; + } +} + + + +void RegressionTest::getRegressionTest(const char* keyword){ + std::vector timeVec1, timeVec2; + setTimeVecs(timeVec1, timeVec2); // Sets the time vectors, they are equal for all keywords (WPOR:PROD01 etc) + setDataSets(timeVec1, timeVec2); //Figures which dataset that contains more/less values pr keyword vector. + + if(stringlist_contains(keysShort,keyword) && stringlist_contains(keysLong, keyword)){ + checkForKeyword(timeVec1, timeVec2, keyword); + return; + } + std::cout << "The keyword suggested, " << keyword << ", is not supported by one or both of the summary files. Please use a different keyword." << std::endl; + OPM_THROW(std::invalid_argument, "Input keyword from user does not exist in/is not common for the two summary files."); +} + + + +void RegressionTest::checkDeviation(Deviation deviation, const char* keyword, int reportStep){ + double absTol = getAbsTolerance(); + double relTol = getRelTolerance(); + + if (deviation.rel > relTol && deviation.abs > absTol){ + std::cout << "For keyword " << keyword << " at reportstep " << reportStep <& timeVec1, std::vector& timeVec2, const char* keyword){ + std::vector dataVec1, dataVec2; + getDataVecs(dataVec1,dataVec2,keyword); + chooseReference(timeVec1, timeVec2,dataVec1,dataVec2); + startTest(keyword); +} + + + +void RegressionTest::startTest(const char* keyword){ + size_t jvar = 0; + Deviation deviation; + for (size_t ivar = 0; ivar < referenceVec->size(); ivar++){ + getDeviation(ivar, jvar, deviation);//Reads from the protected member variables in the super class. + checkDeviation(deviation, keyword,ivar); + } +} diff --git a/opm/test_util/summaryRegressionTest.hpp b/opm/test_util/summaryRegressionTest.hpp new file mode 100644 index 0000000..6dd41dc --- /dev/null +++ b/opm/test_util/summaryRegressionTest.hpp @@ -0,0 +1,59 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include + +//! \details The class inherits from the SummaryComparator class, which takes care of all file reading. \n The RegressionTest class compares the values from the two different files and throws exceptions when the deviation is unsatisfying. +class RegressionTest: public SummaryComparator { + private: + //! \brief Gathers the correct data for comparison for a specific keyword + //! \param[in] timeVec1 The time steps of file 1. + //! \param[in] timeVec2 The time steps of file 2. + //! \param[in] keyword The keyword of interest + //! \details The function prepares a regression test by gathering the data, stroing it into two vectors, \n deciding which is to be used as a reference/basis and calling the test function. + void checkForKeyword(std::vector& timeVec1, std::vector& timeVec2, const char* keyword); + + //! \brief The regression test + //! \param[in] keyword The keyword common for both the files. The vectors associated with the keyword are used for comparison. + //! \details Start test uses the private member variables, pointers of std::vector type, which are set to point to the correct vectors in SummaryComparison::chooseReference(...). \n The function iterates over the referenceVev/basis and for each iteration it calculates the deviation with SummaryComparison::getDeviation(..) and stors it in a Deviation struct. \n SummaryComparison::getDeviation takes the int jvar as an reference input, and using it as an iterative index for the values which are to be compared with the basis. \n Thus, by updating the jvar variable every time a deviation is calculated, one keep track jvar and do not have to iterate over already checked values. + void startTest(const char* keyword); + + //! \brief Caluculates a deviation, throws exceptions and writes and error message. + //! \param[in] deviation Deviation struct + //! \param[in] keyword The keyword that the data that are being compared belongs to. + //! \param[in] reportStep The report step of which the deviation originates from. + //! \details The function checks the values of the Deviation struct against the absolute and relative tolerance, which are private member values of the super class. \n When comparing against the relative tolerance an additional term is added, the absolute deviation has to be greater than 1e-6 for the function to throw an exception. \n When the deviations are too great, the function writes out which keyword, and at what report step the deviation is too great before throwing an exception. + void checkDeviation(Deviation deviation, const char* keyword, int reportStep); + + public: + //! \brief Constructor, creates an object of RefressionTest class. + //! \param[in] basename1 Path to file1 without extension. + //! \param[in] basename1 Path to file2 without extension. + //! \param[in] absoluteTolerance The absolute tolerance which is to be used in the test. + //! \param[in] relativeTolerance The relative tolerance which is to be used in the test. + //! \details The constructor calls the constructor of the super class. + RegressionTest(const char* basename1, const char* basename2, double relativeTolerance, double absoluteTolerance): + SummaryComparator(basename1, basename2,relativeTolerance, absoluteTolerance) {} + + //! \details The function executes a regression test for all the keywords. If the two files do not match in amount of keywords, an exception is thrown. + void getRegressionTest(); + + //! \details The function executes a regression test for one specific keyword. If one or both of the files do not have the keyword, an exception is thrown. + void getRegressionTest(const char* keyword);///< Regression test for a certain keyword of the files. +}; diff --git a/tests/test_EclFilesComparator.cpp b/tests/test_EclFilesComparator.cpp new file mode 100644 index 0000000..9a56f6b --- /dev/null +++ b/tests/test_EclFilesComparator.cpp @@ -0,0 +1,72 @@ +/* + Copyright 2016 Statoil 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 . + */ + +#include "config.h" + +#if HAVE_DYNAMIC_BOOST_TEST +#define BOOST_TEST_DYN_LINK +#endif + +#define BOOST_TEST_MODULE EclFilesComparatorTest + +#include +#include + +BOOST_AUTO_TEST_CASE(deviation) { + double a = 1; + double b = 3; + const double tol = 1.0e-14; + + Deviation dev = ECLFilesComparator::calculateDeviations(a,b); + + BOOST_CHECK_EQUAL(dev.abs, 2.0); + BOOST_CHECK_CLOSE(dev.rel, 2.0/3, tol); + + a = 0; + + dev = ECLFilesComparator::calculateDeviations(a,b); + + BOOST_CHECK_EQUAL(dev.abs, 3.0); + BOOST_CHECK_EQUAL(dev.rel, -1); +} + + + +BOOST_AUTO_TEST_CASE(median) { + std::vector vec = {1,3,4,5}; + + double median = ECLFilesComparator::median(vec); + + BOOST_CHECK_EQUAL(median, 3.5); + + vec = {1,4,5}; + median = ECLFilesComparator::median(vec); + + BOOST_CHECK_EQUAL(median, 4); +} + + + +BOOST_AUTO_TEST_CASE(average) { + std::vector vec = {1,3,4,5}; + const double tol = 1.0e-14; + + double average = ECLFilesComparator::average(vec); + + BOOST_CHECK_CLOSE(average, 13.0/4, tol); +} diff --git a/tests/test_compareSummary.cpp b/tests/test_compareSummary.cpp new file mode 100644 index 0000000..59d3785 --- /dev/null +++ b/tests/test_compareSummary.cpp @@ -0,0 +1,74 @@ +/* + Copyright 2016 Statoil 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 . + */ + + +#include "config.h" +#include + + +#if HAVE_DYNAMIC_BOOST_TEST +#define BOOST_TEST_DYN_LINK +#endif + +#define BOOST_TEST_MODULE CalculationTest +#include + + +BOOST_AUTO_TEST_CASE(deviation){ + double a = 5; + double b = 10; + const double tol = 1.0e-14; + + Deviation dev = SummaryComparator::calculateDeviations(a,b); + BOOST_CHECK_EQUAL(dev.abs, 5); + BOOST_CHECK_CLOSE(dev.rel, 5.0/10.0, tol); +} + + +BOOST_AUTO_TEST_CASE(median){ + std::vector vec = {8.6, 0.6 ,0 , 3.0, 7.2}; + BOOST_CHECK_EQUAL(3,SummaryComparator::median(vec)); + vec.pop_back(); + BOOST_CHECK_EQUAL(1.8, SummaryComparator::median(vec)); + +} + +BOOST_AUTO_TEST_CASE(interpolation){ + double timeArray[3] = {6,1,2}; + double linearCheckValues[2] = {6,11}; + double linearCheckValues_2[2] = {2,12}; + double constCheckValues[2] = {3,3}; + + double linearLP = SummaryComparator::interpolation(linearCheckValues[1],linearCheckValues[0], timeArray); + double constLP = SummaryComparator::interpolation(constCheckValues[1], constCheckValues[0], timeArray); + double linearLP_2 = SummaryComparator::interpolation(linearCheckValues_2[1], linearCheckValues_2[0], timeArray); + BOOST_CHECK_EQUAL(linearLP,7); + BOOST_CHECK_EQUAL(constLP,3); + BOOST_CHECK_EQUAL(linearLP_2,4); +} + + +BOOST_AUTO_TEST_CASE(average) { + + std::vector vec = {1,2,3,4,5,6}; + const double tol = 1.0e-14; + + double average = SummaryComparator::average(vec); + + BOOST_CHECK_CLOSE(average, 21.0/6, tol); +}