Move Static Dake Model Calculation to Connection Class

Certain parts of the Dake model D-factor correlation relation depend
only on static parameters at the cell level.  Calculate the static
product once, and store it as a member of Connection::CTFProperties.
In particular, the new member

    double CTFProperties::static_dfac_corr_coeff

stores the product

    A * (K/K0)**B * poro**C * Ke/(h*rw)

with 'A', 'B', and 'C', being the Dake model parameters.  We compute
new values for the 'static_dfac_corr_coeff' as part of processing
the WDFACCOR keyword.

This static product will also be output to the restart file in
follow-up commits, and will be used to simplify evaluating the full
dynamic Dake model correlation.
This commit is contained in:
Bård Skaflestad 2023-11-10 15:47:46 +01:00
parent 736eb9f7db
commit e0fe776e2e
6 changed files with 152 additions and 38 deletions

View File

@ -119,6 +119,10 @@ namespace Opm {
/// for gas.
double d_factor{};
/// Product of certain static elements of D-factor correlation
/// law (WDFACCOR keyword).
double static_dfac_corr_coeff{};
/// Denominator in peaceman's formula-i.e., log(r0/rw) + skin.
double peaceman_denom{};
@ -157,6 +161,7 @@ namespace Opm {
serializer(this->connection_length);
serializer(this->skin_factor);
serializer(this->d_factor);
serializer(this->static_dfac_corr_coeff);
serializer(this->peaceman_denom);
}
};
@ -241,6 +246,7 @@ namespace Opm {
void setKe(double Ke);
void setCF(double CF);
void setDefaultSatTabId(bool id);
void setStaticDFacCorrCoeff(const double c);
void scaleWellPi(double wellPi);
bool prepareWellPIScaling();

View File

@ -39,6 +39,7 @@ namespace Opm {
class FieldPropsManager;
class KeywordLocation;
class ScheduleGrid;
class WDFAC;
} // namespace Opm
namespace Opm {
@ -88,6 +89,7 @@ namespace Opm {
void loadCOMPDAT(const DeckRecord& record,
const ScheduleGrid& grid,
const std::string& wname,
const WDFAC& wdfac,
const KeywordLocation& location);
void loadCOMPTRAJ(const DeckRecord& record,
@ -101,6 +103,9 @@ namespace Opm {
const std::string& wname,
const KeywordLocation& location);
void applyDFactorCorrelation(const ScheduleGrid& grid,
const WDFAC& wdfac);
int getHeadI() const;
int getHeadJ() const;
const std::vector<double>& getMD() const;

View File

@ -196,37 +196,51 @@ namespace {
this->snapshots.back().network.update( std::move( ext_network ));
}
void Schedule::handleCOMPDAT(HandlerContext& handlerContext) {
void Schedule::handleCOMPDAT(HandlerContext& handlerContext) {
std::unordered_set<std::string> wells;
for (const auto& record : handlerContext.keyword) {
const std::string& wellNamePattern = record.getItem("WELL").getTrimmedString(0);
auto wellnames = this->wellNames(wellNamePattern, handlerContext,
isWList(handlerContext.currentStep,
wellNamePattern));
const auto wellNamePattern = record.getItem("WELL").getTrimmedString(0);
const auto wellnames =
this->wellNames(wellNamePattern, handlerContext,
isWList(handlerContext.currentStep, wellNamePattern));
for (const auto& name : wellnames) {
auto well2 = this->snapshots.back().wells.get(name);
auto connections = std::shared_ptr<WellConnections>( new WellConnections( well2.getConnections()));
connections->loadCOMPDAT(record, handlerContext.grid, name, handlerContext.keyword.location());
if (well2.updateConnections(connections, handlerContext.grid)) {
auto connections = std::make_shared<WellConnections>(well2.getConnections());
const auto origWellConnSetIsEmpty = connections->empty();
connections->loadCOMPDAT(record, handlerContext.grid, name,
well2.getWDFAC(), handlerContext.keyword.location());
const auto newWellConnSetIsEmpty = connections->empty();
if (well2.updateConnections(std::move(connections), handlerContext.grid)) {
auto wdfac = std::make_shared<WDFAC>(well2.getWDFAC());
wdfac->updateWDFACType(*connections);
wdfac->updateWDFACType(well2.getConnections());
well2.updateWDFAC(std::move(wdfac));
this->snapshots.back().wells.update( well2 );
wells.insert( name );
this->snapshots.back().wells.update(well2);
wells.insert(name);
}
if (connections->empty() && well2.getConnections().empty()) {
if (origWellConnSetIsEmpty && newWellConnSetIsEmpty) {
const auto& location = handlerContext.keyword.location();
auto msg = fmt::format("Problem with COMPDAT/{}\n"
"In {} line {}\n"
"Well {} is not connected to grid - will remain SHUT", name, location.filename, location.lineno, name);
const auto msg = fmt::format(R"(Problem with COMPDAT/{}
In {} line {}
Well {} is not connected to grid - will remain SHUT)",
name, location.filename,
location.lineno, name);
OpmLog::warning(msg);
}
this->snapshots.back().wellgroup_events().addEvent( name, ScheduleEvents::COMPLETION_CHANGE);
this->snapshots.back().wellgroup_events()
.addEvent(name, ScheduleEvents::COMPLETION_CHANGE);
}
}
this->snapshots.back().events().addEvent(ScheduleEvents::COMPLETION_CHANGE);
// In the case the wells reference depth has been defaulted in the
@ -234,9 +248,10 @@ namespace {
// reference depth exactly when the COMPDAT keyword has been completely
// processed.
for (const auto& wname : wells) {
auto& well = this->snapshots.back().wells.get( wname );
auto well = this->snapshots.back().wells.get(wname);
well.updateRefDepth();
this->snapshots.back().wells.update( std::move(well));
this->snapshots.back().wells.update(std::move(well));
}
if (! wells.empty()) {
@ -268,20 +283,30 @@ namespace {
this->snapshots.back().events().addEvent(ScheduleEvents::COMPLETION_CHANGE);
}
void Schedule::handleCOMPTRAJ(HandlerContext& handlerContext) {
void Schedule::handleCOMPTRAJ(HandlerContext& handlerContext)
{
// Keyword WELTRAJ must be read first
std::unordered_set<std::string> wells;
external::cvf::ref<external::cvf::BoundingBoxTree> cellSearchTree = nullptr;
for (const auto& record : handlerContext.keyword) {
const std::string& wellNamePattern = record.getItem("WELL").getTrimmedString(0);
auto wellnames = this->wellNames(wellNamePattern, handlerContext );
const auto wellNamePattern = record.getItem("WELL").getTrimmedString(0);
const auto wellnames = this->wellNames(wellNamePattern, handlerContext);
for (const auto& name : wellnames) {
auto well2 = this->snapshots.back().wells.get(name);
auto connections = std::make_shared<WellConnections>(WellConnections(well2.getConnections()));
// cellsearchTree is calculated only once and is used to calculated cell intersections of the perforations specified in COMPTRAJ
connections->loadCOMPTRAJ(record, handlerContext.grid, name, handlerContext.keyword.location(), cellSearchTree);
// In the case that defaults are used in WELSPECS for headI/J the headI/J are calculated based on the well trajectory data
auto connections = std::make_shared<WellConnections>(well2.getConnections());
// cellsearchTree is calculated only once and is used to
// calculated cell intersections of the perforations
// specified in COMPTRAJ
connections->loadCOMPTRAJ(record, handlerContext.grid, name,
handlerContext.keyword.location(),
cellSearchTree);
// In the case that defaults are used in WELSPECS for
// headI/J the headI/J are calculated based on the well
// trajectory data
well2.updateHead(connections->getHeadI(), connections->getHeadJ());
if (well2.updateConnections(connections, handlerContext.grid)) {
this->snapshots.back().wells.update( well2 );
@ -295,19 +320,23 @@ namespace {
"Well {} is not connected to grid - will remain SHUT", name, location.filename, location.lineno, name);
OpmLog::warning(msg);
}
this->snapshots.back().wellgroup_events().addEvent( name, ScheduleEvents::COMPLETION_CHANGE);
this->snapshots.back().wellgroup_events()
.addEvent(name, ScheduleEvents::COMPLETION_CHANGE);
}
}
this->snapshots.back().events().addEvent(ScheduleEvents::COMPLETION_CHANGE);
// In the case the wells reference depth has been defaulted in the
// WELSPECS keyword we need to force a calculation of the wells
// reference depth exactly when the WELCOML keyword has been completely
// processed.
// reference depth exactly when the WELCOML keyword has been
// completely processed.
for (const auto& wname : wells) {
auto& well = this->snapshots.back().wells.get( wname );
auto& well = this->snapshots.back().wells.get(wname);
well.updateRefDepth();
this->snapshots.back().wells.update( std::move(well));
this->snapshots.back().wells.update(std::move(well));
}
if (! wells.empty()) {
@ -1641,25 +1670,55 @@ File {} line {}.)", wname, location.keyword, location.filename, location.lineno)
auto wdfac = std::make_shared<WDFAC>(well.getWDFAC());
wdfac->updateWDFAC( record );
wdfac->updateTotalCF(well.getConnections());
if (well.updateWDFAC(std::move(wdfac)))
this->snapshots.back().wells.update( std::move(well) );
if (well.updateWDFAC(std::move(wdfac))) {
this->snapshots.back().wells.update(std::move(well));
handlerContext.affected_well(well_name);
handlerContext.record_well_structure_change();
this->snapshots.back().events()
.addEvent(ScheduleEvents::COMPLETION_CHANGE);
this->snapshots.back().wellgroup_events()
.addEvent(well_name, ScheduleEvents::COMPLETION_CHANGE);
}
}
}
}
void Schedule::handleWDFACCOR(HandlerContext& handlerContext) {
for (const auto& record : handlerContext.keyword) {
const std::string& wellNamePattern = record.getItem("WELLNAME").getTrimmedString(0);
const std::string wellNamePattern = record.getItem("WELLNAME").getTrimmedString(0);
const auto well_names = wellNames(wellNamePattern, handlerContext.currentStep);
if (well_names.empty())
this->invalidNamePattern(wellNamePattern, handlerContext);
for (const auto& well_name : well_names) {
auto well = this->snapshots.back().wells.get(well_name);
auto conns = std::make_shared<WellConnections>(well.getConnections());
auto wdfac = std::make_shared<WDFAC>(well.getWDFAC());
wdfac->updateWDFACCOR( record );
if (well.updateWDFAC(std::move(wdfac)))
this->snapshots.back().wells.update( std::move(well) );
wdfac->updateWDFACCOR(record);
conns->applyDFactorCorrelation(handlerContext.grid, *wdfac);
const auto updateConns =
well.updateConnections(std::move(conns), handlerContext.grid);
if (well.updateWDFAC(std::move(wdfac)) || updateConns) {
this->snapshots.back().wells.update(std::move(well));
handlerContext.affected_well(well_name);
handlerContext.record_well_structure_change();
if (updateConns) {
this->snapshots.back().events()
.addEvent(ScheduleEvents::COMPLETION_CHANGE);
this->snapshots.back().wellgroup_events()
.addEvent(well_name, ScheduleEvents::COMPLETION_CHANGE);
}
}
}
}
}

View File

@ -78,7 +78,8 @@ namespace Opm
props.connection_length = 7.0;
props.skin_factor = 8.0;
props.d_factor = 9.0;
props.peaceman_denom = 10.0;
props.static_dfac_corr_coeff = 10.0;
props.peaceman_denom = 11.0;
return props;
}
@ -94,6 +95,7 @@ namespace Opm
&& (this->connection_length == that.connection_length)
&& (this->skin_factor == that.skin_factor)
&& (this->d_factor == that.d_factor)
&& (this->static_dfac_corr_coeff == that.static_dfac_corr_coeff)
&& (this->peaceman_denom == that.peaceman_denom)
;
}
@ -394,6 +396,11 @@ namespace Opm
return true;
}
void Connection::setStaticDFacCorrCoeff(const double c)
{
this->ctf_properties_.static_dfac_corr_coeff = c;
}
std::string Connection::str() const
{
std::stringstream ss;

View File

@ -153,6 +153,22 @@ namespace {
return peacemanDenominator(ctf_props.r0, ctf_props.rw, ctf_props.skin_factor);
}
double staticForchheimerCoeffcient(const Opm::Connection::CTFProperties& ctf_props,
const double porosity,
const Opm::WDFAC& wdfac)
{
// Reference/background permeability against which to scale cell
// permeability for model evaluation.
constexpr auto Kref = 1.0*Opm::prefix::milli*Opm::unit::darcy;
const auto& corr_coeff = wdfac.getDFactorCorrelationCoefficients();
return corr_coeff.coeff_a
* std::pow(ctf_props.Ke / Kref, corr_coeff.exponent_b)
* std::pow(porosity , corr_coeff.exponent_c)
* ctf_props.Ke / (ctf_props.connection_length * ctf_props.rw);
}
// Calculate permeability thickness Kh for line segment in a cell for x,y,z directions
std::array<double, 3>
permThickness(const external::cvf::Vec3d& effective_connection,
@ -339,6 +355,7 @@ namespace Opm {
void WellConnections::loadCOMPDAT(const DeckRecord& record,
const ScheduleGrid& grid,
const std::string& wname,
const WDFAC& wdfac,
const KeywordLocation& location)
{
const auto& itemI = record.getItem("I");
@ -492,6 +509,9 @@ The cell ({},{},{}) in well {} is not active and the connection will be ignored)
// PolymerMW module.
ctf_props.re = std::sqrt(D[0] * D[1] / angle * 2);
ctf_props.static_dfac_corr_coeff =
staticForchheimerCoeffcient(ctf_props, props->poro, wdfac);
auto prev = std::find_if(this->m_connections.begin(),
this->m_connections.end(),
[I, J, k](const Connection& c)
@ -734,6 +754,21 @@ CF and Kh items for well {} must both be specified or both defaulted/negative)",
this->md.push_back(record.getItem("MD").getSIDouble(0));
}
void WellConnections::applyDFactorCorrelation(const ScheduleGrid& grid,
const WDFAC& wdfac)
{
for (auto& conn : this->m_connections) {
const auto& complCell = grid.get_cell(conn.getI(), conn.getJ(), conn.getK());
if (! complCell.is_active()) {
continue;
}
conn.setStaticDFacCorrCoeff
(staticForchheimerCoeffcient(conn.ctfProperties(),
complCell.props->poro, wdfac));
}
}
std::size_t WellConnections::size() const
{
return m_connections.size();

View File

@ -35,6 +35,7 @@
#include <opm/input/eclipse/Schedule/Schedule.hpp>
#include <opm/input/eclipse/Schedule/ScheduleGrid.hpp>
#include <opm/input/eclipse/Schedule/Well/Well.hpp>
#include <opm/input/eclipse/Schedule/Well/WDFAC.hpp>
#include <opm/input/eclipse/Schedule/Well/WellConnections.hpp>
#include <opm/common/OpmLog/KeywordLocation.hpp>
@ -65,6 +66,7 @@ namespace {
};
const auto deck = Opm::Parser{}.parseString(compdat_keyword);
const auto wdfac = Opm::WDFAC{};
const auto loc = Opm::KeywordLocation{};
const Opm::EclipseGrid grid { 10, 10, 10 };
@ -77,7 +79,7 @@ namespace {
const auto sg = Opm::ScheduleGrid { grid, field_props, completed_cells };
for (const auto& rec : deck["COMPDAT"][0]) {
connections.loadCOMPDAT(rec, sg, "WELL", loc);
connections.loadCOMPDAT(rec, sg, "WELL", wdfac, loc);
}
return connections;