addressing review comments for OPM/opm-common#2246

This commit is contained in:
Kai Bao
2021-01-27 22:50:50 +01:00
parent af882dc483
commit abbb56c173
8 changed files with 66 additions and 65 deletions
@@ -30,21 +30,21 @@ namespace Opm {
struct NumericalAquiferCell {
NumericalAquiferCell(const DeckRecord&, const EclipseGrid&, const FieldPropsManager&);
NumericalAquiferCell() = default;
size_t aquifer_id; // aquifer id
size_t I, J, K; // indices for the grid block
double area; // cross-sectional area
size_t aquifer_id;
size_t I, J, K;
double area;
double length;
double porosity;
double permeability;
double depth; // by default the grid block depth will be used
std::optional<double> init_pressure; // by default, the grid pressure from equilibration will be used
int pvttable; // by default, the block PVTNUM
int sattable; // saturation table number, by default, the block value
double transmissibility;
double depth;
std::optional<double> init_pressure;
int pvttable;
int sattable;
size_t global_index;
double cellVolume() const;
double poreVolume() const;
double transmissiblity() const;
bool operator == (const NumericalAquiferCell& other) const;
template<class Serializer>
@@ -60,7 +60,6 @@ namespace Opm {
serializer(this->depth);
serializer(this->pvttable);
serializer(this->sattable);
serializer(this->transmissibility);
serializer(this->global_index);
}
};
@@ -23,6 +23,7 @@
#include <opm/parser/eclipse/EclipseState/Grid/FaceDir.hpp>
#include <map>
#include <array>
#include <vector>
namespace Opm {
@@ -28,7 +28,7 @@ namespace Opm {
class Deck;
class EclipseGrid;
class FieldPropsManager;
class SingleNumericalAquifer;
struct NumericalAquiferCell;
class NumericalAquifers {
public:
@@ -37,24 +37,20 @@ namespace Opm {
size_t numAquifer() const;
bool hasAquifer(size_t aquifer_id) const;
bool hasCell(size_t global_index) const;
const SingleNumericalAquifer& getAquifer(size_t aquifer_id) const;
const NumericalAquiferCell& getCell(size_t global_index) const;
bool operator==(const NumericalAquifers& other) const;
std::unordered_map<size_t, const NumericalAquiferCell*> allAquiferCells() const;
static NumericalAquifers serializeObject();
template <class Serializer>
void serializeOp(Serializer& serializer)
{
serializer.map(this->aquifers_);
// TODO: serialize pointers is problematic, for pointers should we regenerate the pointers?
// TODO: if yes, we need to write the function to generate this->aquifer_cells_
// serializer.map(this->aquifer_cells_);
serializer.map(this->m_aquifers);
}
private:
std::unordered_map <size_t, SingleNumericalAquifer> aquifers_;
std::unordered_map<size_t, const NumericalAquiferCell*> aquifer_cells_;
std::unordered_map <size_t, SingleNumericalAquifer> m_aquifers;
void addAquiferCell(const NumericalAquiferCell& aqu_cell);
@@ -70,8 +70,6 @@ namespace Opm {
} else {
this->sattable = satnum[active_index];
}
this->transmissibility = 2. * this->permeability * this->area / this->length;
}
double NumericalAquiferCell::cellVolume() const {
@@ -91,13 +89,16 @@ namespace Opm {
this->init_pressure == other.init_pressure &&
this->pvttable == other.pvttable &&
this->sattable == other.sattable &&
this->transmissibility == other.transmissibility &&
this->global_index == other.global_index;
}
double NumericalAquiferCell::poreVolume() const {
return this->porosity * this->cellVolume();
}
double NumericalAquiferCell::transmissiblity() const {
return 2. * this->permeability * this->area / this->length;;
}
}
@@ -34,9 +34,10 @@ namespace Opm {
std::map<size_t, std::map<size_t, NumericalAquiferConnection>>
NumericalAquiferConnection::generateConnections(const Deck &deck, const EclipseGrid &grid)
{
std::map<size_t, std::map<size_t, NumericalAquiferConnection>> connections;
using AQUCON=ParserKeywords::AQUCON;
if ( !deck.hasKeyword<AQUCON>() ) return connections;
if ( !deck.hasKeyword<AQUCON>() ) return {};
std::map<size_t, std::map<size_t, NumericalAquiferConnection>> connections;
const auto& aqucon_keywords = deck.getKeywordList<AQUCON>();
for (const auto& keyword : aqucon_keywords) {
@@ -83,16 +84,15 @@ namespace Opm {
std::vector<NumericalAquiferConnection> cons;
const size_t i1 = record.getItem<AQUCON::I1>().get<int>(0) - 1;
const size_t j1 = record.getItem<AQUCON::J1>().get<int>(0) -1;
const size_t j1 = record.getItem<AQUCON::J1>().get<int>(0) - 1;
const size_t k1 = record.getItem<AQUCON::K1>().get<int>(0) - 1;
const size_t i2 = record.getItem<AQUCON::I2>().get<int>(0) - 1;
const size_t j2 = record.getItem<AQUCON::J2>().get<int>(0) -1;
const size_t j2 = record.getItem<AQUCON::J2>().get<int>(0) - 1;
const size_t k2 = record.getItem<AQUCON::K2>().get<int>(0) - 1;
const std::string str_allow_internal_cells = record.getItem<AQUCON::ALLOW_INTERNAL_CELLS>().getTrimmedString(0);
// whether the connection face can connect to active/internal cells
// by default NO, which means basically the aquifer should be outside of the reservoir
const bool allow_internal_cells = DeckItem::to_bool(str_allow_internal_cells);
const bool allow_internal_cells = DeckItem::to_bool( record.getItem<AQUCON::ALLOW_INTERNAL_CELLS>().getTrimmedString(0) );
const FaceDir::DirEnum face_dir
= FaceDir::FromString(record.getItem<AQUCON::CONNECT_FACE>().getTrimmedString(0));
@@ -103,13 +103,10 @@ namespace Opm {
if (!grid.cellActive(i, j, k)) {
continue;
}
/* if (!actnum[grid.getGlobalIndex(i, j, k)]) {
continue;
}*/
if (allow_internal_cells ||
!AquiferHelpers::neighborCellInsideReservoirAndActive(grid, i, j, k, face_dir)) {
const size_t global_index = grid.getGlobalIndex(i, j, k);
cons.emplace_back(NumericalAquiferConnection{i, j, k, global_index, allow_internal_cells, record});
cons.emplace_back(i, j, k, global_index, allow_internal_cells, record);
}
}
}
@@ -25,10 +25,13 @@
#include <opm/common/OpmLog/KeywordLocation.hpp>
#include <opm/parser/eclipse/Deck/Deck.hpp>
#include <opm/parser/eclipse/Deck/DeckRecord.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferCell.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquifers.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/SingleNumericalAquifer.hpp>
#include <set>
namespace Opm {
@@ -37,18 +40,20 @@ namespace Opm {
using AQUNUM=ParserKeywords::AQUNUM;
if ( !deck.hasKeyword<AQUNUM>() ) return;
std::set<size_t> cells;
// there might be multiple keywords of keyword AQUNUM, it is not totally
// clear about the rules here. For now, we take care of all the keywords
const auto& aqunum_keywords = deck.getKeywordList<AQUNUM>();
for (const auto& keyword : aqunum_keywords) {
for (const auto& record : *keyword) {
const NumericalAquiferCell aqu_cell(record, grid, field_props);
if (this->hasCell(aqu_cell.global_index)) {
if (cells.count(aqu_cell.global_index) > 0) {
auto error = fmt::format("Numerical aquifer cell at ({}, {}, {}) is declared more than once",
aqu_cell.I + 1, aqu_cell.J + 1, aqu_cell.K + 1);
throw OpmInputError(error, keyword->location());
} else {
this->addAquiferCell(aqu_cell);
cells.insert(aqu_cell.global_index);
}
}
}
@@ -60,34 +65,21 @@ namespace Opm {
void NumericalAquifers::addAquiferCell(const NumericalAquiferCell& aqu_cell) {
const size_t id = aqu_cell.aquifer_id;
if (!this->hasAquifer(id)) {
this->aquifers_.insert(std::make_pair(id, SingleNumericalAquifer{id}));
this->m_aquifers.insert(std::make_pair(id, SingleNumericalAquifer{id}));
}
auto& aquifer = this->aquifers_.at(id);
auto& aquifer = this->m_aquifers.at(id);
aquifer.addAquiferCell(aqu_cell);
this->aquifer_cells_.insert(std::pair{aqu_cell.global_index, aquifer.getCellPrt(aquifer.numCells())});
}
bool NumericalAquifers::hasCell(const size_t cell_global_index) const {
const auto& cells = this->aquifer_cells_;
return (cells.find(cell_global_index) != cells.end());
}
const NumericalAquiferCell& NumericalAquifers::getCell(const size_t cell_global_index) const {
assert(this->hasCell(cell_global_index));
return *(this->aquifer_cells_.at(cell_global_index));
}
bool NumericalAquifers::hasAquifer(const size_t aquifer_id) const {
return (this->aquifers_.find(aquifer_id) != this->aquifers_.end());
return (this->m_aquifers.find(aquifer_id) != this->m_aquifers.end());
}
void NumericalAquifers::addAquiferConnections(const Deck& deck, const EclipseGrid& grid) {
const auto aquifer_connections = NumericalAquiferConnection::generateConnections(deck, grid);
for (auto& pair : this->aquifers_) {
for (auto& pair : this->m_aquifers) {
const size_t aqu_id = pair.first;
const auto& aqu_cons = aquifer_connections.find(aqu_id);
if (aqu_cons == aquifer_connections.end()) {
@@ -97,13 +89,15 @@ namespace Opm {
auto& aquifer = pair.second;
const auto& cons = aqu_cons->second;
const auto all_aquifer_cells = this->allAquiferCells();
// For now, there is no two aquifers can be connected to one cell
// aquifer can not connect to aquifer cells
for (const auto& con : cons) {
const auto& aqu_con = con.second;
const size_t con_global_index = aqu_con.global_index;
if (this->hasCell(con_global_index)) {
const size_t cell_aquifer_id = this->getCell(con_global_index).aquifer_id;
const auto cell_iter = all_aquifer_cells.find(con_global_index);
if (cell_iter != all_aquifer_cells.end()) {
const size_t cell_aquifer_id = cell_iter->second->aquifer_id;
auto msg = fmt::format("Problem with keyword AQUCON \n"
"Aquifer connection declared at grid cell ({}, {}, {}), is a aquifer cell "
"of Aquifer {}, and will be removed",
@@ -118,24 +112,37 @@ namespace Opm {
}
bool NumericalAquifers::operator==(const NumericalAquifers& other) const {
return this->aquifers_ == other.aquifers_;
return this->m_aquifers == other.m_aquifers;
}
size_t NumericalAquifers::numAquifer() const {
return this->aquifers_.size();
return this->m_aquifers.size();
}
NumericalAquifers NumericalAquifers::serializeObject() {
NumericalAquifers result;
result.aquifers_ = {{1, SingleNumericalAquifer{1}}};
result.m_aquifers = {{1, SingleNumericalAquifer{1}}};
return result;
}
const SingleNumericalAquifer& NumericalAquifers::getAquifer(const size_t aquifer_id) const {
if ( !this->hasAquifer(aquifer_id)) {
const auto msg = fmt::format(" There is no numerical aquifer {}", aquifer_id);
throw std::runtime_error(msg);
const auto iter = this->m_aquifers.find(aquifer_id);
if ( iter != this->m_aquifers.end() ) {
return iter->second;
} else {
const auto msg = fmt::format(" There is no numerical aquifer {}", aquifer_id);
throw std::runtime_error(msg);
}
return this->aquifers_.at(aquifer_id);
}
std::unordered_map<size_t, const NumericalAquiferCell*> NumericalAquifers::allAquiferCells() const {
std::unordered_map<size_t, const NumericalAquiferCell*> cells;
for (const auto& [id, aquifer] : this->m_aquifers) {
for (size_t i = 0; i < aquifer.numCells(); ++i) {
const NumericalAquiferCell* cell_ptr = aquifer.getCellPrt(i);
cells.insert(std::make_pair(cell_ptr->global_index, cell_ptr));
}
}
return cells;
}
}
@@ -17,9 +17,9 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/SingleNumericalAquifer.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferConnection.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferCell.hpp>
#include <opm/parser/eclipse/EclipseState/Aquifer/NumericalAquifer/SingleNumericalAquifer.hpp>
namespace Opm {
SingleNumericalAquifer::SingleNumericalAquifer(const size_t aqu_id)
+5 -5
View File
@@ -618,11 +618,11 @@ BOOST_AUTO_TEST_CASE(NumericalAquiferTest){
const Opm::NumericalAquifers num_aqu{numaquifer_deck, grid, ecl_state.fieldProps()};
BOOST_CHECK(num_aqu.hasAquifer(1));
BOOST_CHECK(num_aqu.numAquifer() == 1);
BOOST_CHECK(num_aqu.hasCell(0));
BOOST_CHECK(!num_aqu.hasCell(1));
BOOST_CHECK(num_aqu.hasCell(2));
BOOST_CHECK(num_aqu.hasCell(3));
BOOST_CHECK(!num_aqu.hasAquifer(2));
const auto all_aquifer_cells = num_aqu.allAquiferCells();
BOOST_CHECK(all_aquifer_cells.count(0) > 0);
BOOST_CHECK(all_aquifer_cells.count(2) > 0);
BOOST_CHECK(all_aquifer_cells.count(3) > 0);
BOOST_CHECK(all_aquifer_cells.count(1) == 0);
const auto& aquifer = num_aqu.getAquifer(1);
BOOST_CHECK(aquifer.numCells() == 3);