diff --git a/opm/parser/eclipse/EclipseState/Aquancon.hpp b/opm/parser/eclipse/EclipseState/Aquancon.hpp index a7a9074f6..c7ab9389f 100755 --- a/opm/parser/eclipse/EclipseState/Aquancon.hpp +++ b/opm/parser/eclipse/EclipseState/Aquancon.hpp @@ -47,6 +47,7 @@ namespace Opm { std::vector> influx_coeff; // Size = size(global_index) std::vector influx_multiplier; // Size = size(global_index) std::vector reservoir_face_dir; // Size = size(global_index) + std::vector record_index; }; Aquancon(const EclipseGrid& grid, const Deck& deck); @@ -55,7 +56,7 @@ namespace Opm { private: - void logic_application(std::vector& output_vector); + std::vector logic_application(const std::vector original_vector); void collate_function(std::vector& output_vector, std::vector& m_aqurecord, diff --git a/src/opm/parser/eclipse/EclipseState/Aquancon.cpp b/src/opm/parser/eclipse/EclipseState/Aquancon.cpp index fd8be5070..edb6c3681 100755 --- a/src/opm/parser/eclipse/EclipseState/Aquancon.cpp +++ b/src/opm/parser/eclipse/EclipseState/Aquancon.cpp @@ -18,8 +18,10 @@ */ #include #include +#include #include #include +#include namespace Opm { namespace{ @@ -33,7 +35,8 @@ namespace Opm { std::vector> influx_coeff_per_record; //Aquifer influx coefficient std::vector influx_mult_per_record; //Aquifer influx coefficient Multiplier // Cell face to connect aquifer to - std::vector face_per_record; + std::vector face_per_record; + std::vector record_index_per_record; }; } @@ -53,8 +56,8 @@ namespace Opm { m_aquiferID_per_record.resize(aquanconKeyword.size()); // We now do a loop over each record entry in aquancon - for (size_t aquanconRecordIdx = 0; aquanconRecordIdx < aquanconKeyword.size(); ++ aquanconRecordIdx) - { + for (size_t aquanconRecordIdx = 0; aquanconRecordIdx < aquanconKeyword.size(); ++aquanconRecordIdx) + { const auto& aquanconRecord = aquanconKeyword.getRecord(aquanconRecordIdx); m_aquiferID_per_record.at(aquanconRecordIdx) = aquanconRecord.getItem("AQUIFER_ID").template get(0); @@ -98,13 +101,14 @@ namespace Opm { m_aqurecord.at(aquanconRecordIdx).influx_mult_per_record.resize(global_index_per_record_size,influx_mult); m_aqurecord.at(aquanconRecordIdx).face_per_record.resize(global_index_per_record_size,faceDir); + m_aqurecord.at(aquanconRecordIdx).record_index_per_record.resize(global_index_per_record_size,aquanconRecordIdx); } // Collate_function collate_function(m_aquoutput, m_aqurecord, m_aquiferID_per_record, m_maxAquID); // Logic for grid connection applied here - logic_application(m_aquoutput); + m_aquoutput = logic_application(m_aquoutput); } @@ -149,22 +153,102 @@ namespace Opm { m_aqurecord.at(record_index_matching_id).face_per_record.begin(), m_aqurecord.at(record_index_matching_id).face_per_record.end() ); + // This is for the record index in order for us to know which one is updated + output_vector.at(i - 1).record_index.insert( + output_vector.at(i - 1).record_index.end(), + m_aqurecord.at(record_index_matching_id).record_index_per_record.begin(), + m_aqurecord.at(record_index_matching_id).record_index_per_record.end() + ); } } } - void Aquancon::logic_application(std::vector& output_vector) + std::vector Aquancon::logic_application(const std::vector original_vector) { - // Find if Global index is repeated for each aquifer, if so, select only the first one + std::vector output_vector = original_vector; + + // Create a local struct to couple each element for easy sorting + struct pair_elements + { + size_t global_index; + std::shared_ptr influx_coeff; + double influx_multiplier; + int reservoir_face_dir; + int record_index; + }; + + // Create a working buffer + std::vector working_buffer; + + // Iterate through each aquifer IDs (This is because each element in the original vector represents an aquifer ID) for (auto aquconvec = output_vector.begin(); aquconvec != output_vector.end(); ++aquconvec) { - std::sort(aquconvec->global_index.begin(), aquconvec->global_index.end()); - auto it = std::unique ( aquconvec->global_index.begin(), aquconvec->global_index.end() ); - aquconvec->global_index.resize( std::distance(aquconvec->global_index.begin(),it) ); - } + //Begin to fill the working buffer + working_buffer.clear(); + working_buffer.resize(aquconvec->global_index.size()); + for (size_t i = 0; i < aquconvec->global_index.size(); ++i ) + { + working_buffer.at(i).global_index = aquconvec->global_index.at(i); + working_buffer.at(i).influx_coeff = aquconvec->influx_coeff.at(i); + working_buffer.at(i).influx_multiplier = aquconvec->influx_multiplier.at(i); + working_buffer.at(i).reservoir_face_dir = aquconvec->reservoir_face_dir.at(i); + working_buffer.at(i).record_index = aquconvec->record_index.at(i); + } - //TODO: Find if face on outside of reservoir or adjoins an inactive cell - //TODO: Total number of grid blocks connected to aquifer must not exceed item 6 of AQUDIMS + // Sort by ascending order the working_buffer vector in order of priority: + // 1) global_index, then 2) record_index + + std::sort( working_buffer.begin(), + working_buffer.end(), + [](pair_elements& element1, pair_elements& element2) -> bool + { + if (element1.global_index == element2.global_index) + { + return element1.record_index < element2.record_index; + } + else + return element1.global_index < element2.global_index; + } + ); + + // We then proceed to obtain unique elements of the global_index, and we apply the + // following behaviour (as mentioned in the Eclipse 2014.1 Reference Manual p.345): + // If a reservoir cell is defined more than once, its previous value for the + // aquifer influx coefficient is added to the present value. + + auto i2 = std::unique( working_buffer.begin(), + working_buffer.end(), + [](pair_elements& element1, pair_elements& element2) -> bool + { + if (element1.global_index == element2.global_index) + { + *(element1.influx_coeff) += *(element2.influx_coeff); + return true; + } + + return false; + } + ); + working_buffer.resize( std::distance(working_buffer.begin(),i2) ); + + // We now resize and fill the output_vector elements + aquconvec->global_index.resize(working_buffer.size()); + aquconvec->influx_coeff.resize(working_buffer.size()); + aquconvec->influx_multiplier.resize(working_buffer.size()); + aquconvec->reservoir_face_dir.resize(working_buffer.size()); + aquconvec->record_index.resize(working_buffer.size()); + for (size_t i = 0; i < working_buffer.size(); ++i) + { + aquconvec->global_index.at(i) = working_buffer.at(i).global_index; + aquconvec->influx_coeff.at(i) = working_buffer.at(i).influx_coeff; + aquconvec->influx_multiplier.at(i) = working_buffer.at(i).influx_multiplier; + aquconvec->reservoir_face_dir.at(i) = working_buffer.at(i).reservoir_face_dir; + aquconvec->record_index.at(i) = working_buffer.at(i).record_index; + } + + } + + return output_vector; } void Aquancon::convert_record_id_to_aquifer_id(std::vector& record_indices_matching_id, diff --git a/tests/parser/AquanconTests.cpp b/tests/parser/AquanconTests.cpp index 7105f1b87..eae1bce75 100755 --- a/tests/parser/AquanconTests.cpp +++ b/tests/parser/AquanconTests.cpp @@ -50,6 +50,10 @@ inline Deck createAQUANCONDeck() { "\n" "AQUANCON\n" " 1 1 1 1 1 1 1 J- 1.0 1.0 NO /\n" + " 1 1 3 1 3 3 3 I+ 0.5 1.0 NO /\n" + " 1 1 3 1 3 3 3 J+ 0.75 1.0 NO /\n" + " 1 1 3 1 3 3 3 J- 2.75 1.0 NO /\n" + " 1 2 3 2 3 1 1 I+ 2.75 1.0 NO /\n" "/ \n"; Parser parser; @@ -65,13 +69,31 @@ inline std::vector init_aquancon(){ return aquifers; } +inline std::vector fill_result(){ + auto deck = createAQUANCONDeck(); + EclipseState eclState( deck ); + Aquancon aqucon( eclState.getInputGrid(), deck); + std::vector aquifers = aqucon.getAquOutput(); + + return aquifers; +} + BOOST_AUTO_TEST_CASE(AquanconTest){ std::vector< Aquancon::AquanconOutput > aquifers = init_aquancon(); - for (const auto& it : aquifers){ - for (size_t i = 0; i < it.global_index.size(); ++i){ - BOOST_CHECK_EQUAL(it.aquiferID , 1); - BOOST_CHECK_EQUAL(it.global_index.at(i) , 0); - BOOST_CHECK_EQUAL(it.reservoir_face_dir.at(i) , 8); - } - } -} + std::vector< Aquancon::AquanconOutput > expected_output = fill_result(); + + BOOST_CHECK_EQUAL(aquifers.size(), expected_output.size()); + for (size_t i = 0; i < aquifers.size(); ++i) + { + BOOST_CHECK_EQUAL_COLLECTIONS( aquifers.at(i).global_index.begin(), aquifers.at(i).global_index.end(), + expected_output.at(i).global_index.begin(), expected_output.at(i).global_index.end() ); + BOOST_CHECK_EQUAL_COLLECTIONS( aquifers.at(i).influx_coeff.begin(), aquifers.at(i).influx_coeff.end(), + expected_output.at(i).influx_coeff.begin(), expected_output.at(i).influx_coeff.end() ); + BOOST_CHECK_EQUAL_COLLECTIONS( aquifers.at(i).influx_multiplier.begin(), aquifers.at(i).influx_multiplier.end(), + expected_output.at(i).influx_multiplier.begin(), expected_output.at(i).influx_multiplier.end() ); + BOOST_CHECK_EQUAL_COLLECTIONS( aquifers.at(i).reservoir_face_dir.begin(), aquifers.at(i).reservoir_face_dir.end(), + expected_output.at(i).reservoir_face_dir.begin(), expected_output.at(i).reservoir_face_dir.end() ); + BOOST_CHECK_EQUAL_COLLECTIONS( aquifers.at(i).record_index.begin(), aquifers.at(i).record_index.end(), + expected_output.at(i).record_index.begin(), expected_output.at(i).record_index.end() ); + } +}