From 0bff57af338ec58d4e3f07c349c9e8553c7b5b12 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Wed, 18 Dec 2024 15:12:28 +0100 Subject: [PATCH] Collect Statistics on Connection Level Fracturing Process This commit adds a new member function Fracture::assignGeomechWellState(ConnFracStatistics&) which computes statistical measures (minimum value, maximum value, average value/mean, and standard deviation) of fracture pressures, injection flow rate, and fracture width for a single fracture. We also add a new member function FractureModel::assignGeomechWellState(WellState& wellState) which calls the above for each fracture in the model. Collectively, these enable emitting the new summary vectors CFR{MAX,MIN,AVG,STD} for the quantity as 'P' (pressure), 'IR' (injection rate), and 'WD' (width) through the regular mechanism. --- opm/geomech/Fracture.cpp | 39 +++++++++++++++++++++-- opm/geomech/Fracture.hpp | 12 ++++++- opm/geomech/FractureModel.cpp | 54 +++++++++++++++++++++++++++++++- opm/geomech/FractureModel.hpp | 13 +++++++- opm/geomech/eclproblemgeomech.hh | 3 ++ 5 files changed, 115 insertions(+), 6 deletions(-) diff --git a/opm/geomech/Fracture.cpp b/opm/geomech/Fracture.cpp index 75b85d2..c7d1ca0 100644 --- a/opm/geomech/Fracture.cpp +++ b/opm/geomech/Fracture.cpp @@ -5,6 +5,8 @@ #include #include #include +#include + #include #include // needed for printSparseMatrix?? #include // needed for printSparseMatrix?? @@ -1040,8 +1042,36 @@ std::vector> Fracture::wellIndices() const{ return wellIndices; } -void -Fracture::writePressureSystem() const{ +template +void Fracture::assignGeomechWellState(ConnFracStatistics& stats) const +{ + using Quantity = typename ConnFracStatistics::Quantity; + + constexpr auto pressIx = static_cast>(Quantity::Pressure); + constexpr auto rateIx = static_cast>(Quantity::FlowRate); + constexpr auto widthIx = static_cast>(Quantity::Width); + + const auto nCells = this->reservoir_cells_.size(); + + stats.reset(); + + for (auto cellIx = 0*nCells; cellIx < nCells; ++cellIx) { + auto samplePoint = typename ConnFracStatistics::SamplePoint{}; + + samplePoint[pressIx] = this->fracture_pressure_[cellIx][0]; + + samplePoint[rateIx] = this->leakof_[cellIx] + * (this->fracture_pressure_[cellIx][0] - + this->reservoir_pressure_[cellIx]); + + samplePoint[widthIx] = this->fracture_width_[cellIx][0]; + + stats.addSamplePoint(samplePoint); + } +} + +void Fracture::writePressureSystem() const +{ if(prm_.get("write_pressure_system")){ Dune::storeMatrixMarket(*pressure_matrix_, "pressure_matrix"); Dune::storeMatrixMarket(rhs_pressure_, "pressure_rhs"); @@ -1372,11 +1402,14 @@ void Fracture::printMechMatrix() const // debug purposes // } - // template void // Fracture::updateReservoirCells(const external::cvf::ref& cellSearchTree, // const Dune::CpGrid& grid3D); // template void // Fracture::updateReservoirCells>(const external::cvf::ref& cellSearchTree, // const Dune::PolyhedralGrid<3,3,double>& grid3D); + +template void Fracture::assignGeomechWellState(ConnFracStatistics&) const; +template void Fracture::assignGeomechWellState(ConnFracStatistics&) const; + } // namespace Opm diff --git a/opm/geomech/Fracture.hpp b/opm/geomech/Fracture.hpp index 877d096..a4d561e 100644 --- a/opm/geomech/Fracture.hpp +++ b/opm/geomech/Fracture.hpp @@ -53,6 +53,11 @@ #include +namespace Opm { + template + class ConnFracStatistics; +} + namespace Opm { struct WellInfo { @@ -163,12 +168,17 @@ class Fracture void setFractureGrid(std::unique_ptr gptr = nullptr); // a hack to allow use of another grid std::vector> wellIndices() const; WellInfo& wellInfo(){return wellinfo_;} + const WellInfo& wellInfo() const { return wellinfo_; } std::vector leakOfRate() const; double injectionPressure() const; void setPerfPressure(double perfpressure){perf_pressure_ = perfpressure;} Dune::FieldVector stress(Dune::FieldVector obs) const; Dune::FieldVector strain(Dune::FieldVector obs) const; - Dune::FieldVector disp(Dune::FieldVector obs) const; + Dune::FieldVector disp(Dune::FieldVector obs) const; + + template + void assignGeomechWellState(ConnFracStatistics& stats) const; + private: size_t numFractureCells() const { return grid_->leafGridView().size(0); } diff --git a/opm/geomech/FractureModel.cpp b/opm/geomech/FractureModel.cpp index ab71de0..872f330 100644 --- a/opm/geomech/FractureModel.cpp +++ b/opm/geomech/FractureModel.cpp @@ -4,6 +4,21 @@ #include #include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + namespace Opm{ void FractureModel::addWell(std::string name, const std::vector& points, @@ -180,8 +195,10 @@ namespace Opm{ } } } + std::vector> - FractureModel::getExtraWellIndices(std::string wellname) const{ + FractureModel::getExtraWellIndices(const std::string& wellname) const + { // for now just do a search bool addconnections = prm_.get("addconnections"); if (addconnections) { @@ -203,4 +220,39 @@ namespace Opm{ std::vector> empty; return empty; } + + template + void FractureModel::assignGeomechWellState(WellState& wellState) const + { + const auto nWells = this->wells_.size(); + for (auto i = 0*nWells; i < nWells; ++i) { + const auto wsIx = wellState.index(this->wells_[i].name()); + if (! wsIx.has_value()) { continue; } + + auto& perfData = wellState[*wsIx].perf_data; + + if (perfData.connFracStatistics.size() != perfData.cell_index.size()) { + perfData.connFracStatistics.resize(perfData.cell_index.size()); + } + + for (const auto& fracture : this->well_fractures_[i]) { + auto perfPos = std::find(perfData.cell_index.begin(), + perfData.cell_index.end(), + fracture.wellInfo().well_cell); + if (perfPos == perfData.cell_index.end()) { continue; } + + // Possibly just "fracture.wellInfo().perf" instead. + const auto perfIx = std::distance(perfData.cell_index.begin(), perfPos); + + fracture.assignGeomechWellState(perfData.connFracStatistics[perfIx]); + } + } + } } + +// =========================================================================== +// Explicit specialisations. No other code below separator. +// =========================================================================== + +template void Opm::FractureModel::assignGeomechWellState(WellState&) const; +template void Opm::FractureModel::assignGeomechWellState(WellState&) const; diff --git a/opm/geomech/FractureModel.hpp b/opm/geomech/FractureModel.hpp index 15d205d..7306989 100644 --- a/opm/geomech/FractureModel.hpp +++ b/opm/geomech/FractureModel.hpp @@ -23,6 +23,11 @@ std::vector as_vector(ptree const& pt, ptree::key_type const& key) return r; } +namespace Opm { + template + class WellState; +} + // template // std::vector opm_as_vector(const Opm::PropertyTree& pt, const std::string& key) // { @@ -135,7 +140,13 @@ class FractureModel{ well.updateReservoirProperties(simulator); } } - std::vector> getExtraWellIndices(std::string wellname) const; + + std::vector> + getExtraWellIndices(const std::string& wellname) const; + + template + void assignGeomechWellState(WellState& wellState) const; + bool addPertsToSchedule(){return prm_.get("addperfs_to_schedule");} // probably this should be collected in one loop since all do full loop over fracture ++ well Dune::FieldVector stress(Dune::FieldVector obs) const; diff --git a/opm/geomech/eclproblemgeomech.hh b/opm/geomech/eclproblemgeomech.hh index 84584ba..4074527 100644 --- a/opm/geomech/eclproblemgeomech.hh +++ b/opm/geomech/eclproblemgeomech.hh @@ -263,6 +263,9 @@ namespace Opm{ assert(false); this->addConnectionsToWell(); } + + this->geomechModel_.fractureModel() + .assignGeomechWellState(this->wellModel_.wellState()); } }