Merge pull request #3125 from bska/aquct-dimless-time-and-press-smry

Report Dimensionless Time and Pressure for CT Aquifers
This commit is contained in:
Joakim Hove 2021-03-25 07:17:43 +01:00 committed by GitHub
commit d4a715aba5
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23

View File

@ -26,7 +26,9 @@
#include <opm/output/data/Aquifer.hpp> #include <opm/output/data/Aquifer.hpp>
#include <exception> #include <exception>
#include <memory>
#include <stdexcept> #include <stdexcept>
#include <utility>
namespace Opm namespace Opm
{ {
@ -80,7 +82,12 @@ public:
} }
data.volume = this->W_flux_.value(); data.volume = this->W_flux_.value();
data.initPressure = this->pa0_; data.initPressure = this->pa0_;
data.type = Opm::data::AquiferType::CarterTracey; data.type = Opm::data::AquiferType::CarterTracy;
data.aquCT = std::make_shared<data::CarterTracyData>();
data.aquCT->dimensionless_time = this->dimensionless_time_;
data.aquCT->dimensionless_pressure = this->dimensionless_pressure_;
return data; return data;
} }
@ -92,20 +99,36 @@ protected:
Scalar mu_w_{1}; // water viscosity Scalar mu_w_{1}; // water viscosity
Scalar fluxValue_{0}; // value of flux Scalar fluxValue_{0}; // value of flux
Scalar dimensionless_time_{0};
Scalar dimensionless_pressure_{0};
void assignRestartData(const data::AquiferData& /* xaq */) override void assignRestartData(const data::AquiferData& /* xaq */) override
{ {
throw std::runtime_error {"Restart-based initialization not currently supported " throw std::runtime_error {"Restart-based initialization not currently supported "
"for Carter-Tracey analytic aquifers"}; "for Carter-Tracey analytic aquifers"};
} }
inline void getInfluenceTableValues(Scalar& pitd, Scalar& pitd_prime, const Scalar& td) std::pair<Scalar, Scalar>
getInfluenceTableValues(const Scalar td_plus_dt)
{ {
// We use the opm-common numeric linear interpolator // We use the opm-common numeric linear interpolator
pitd = Opm::linearInterpolation(aquct_data_.td, aquct_data_.pi, td); this->dimensionless_pressure_ =
pitd_prime = Opm::linearInterpolationDerivative(aquct_data_.td, aquct_data_.pi, td); Opm::linearInterpolation(this->aquct_data_.td,
this->aquct_data_.pi,
this->dimensionless_time_);
const auto PItd =
Opm::linearInterpolation(this->aquct_data_.td,
this->aquct_data_.pi, td_plus_dt);
const auto PItdprime =
Opm::linearInterpolationDerivative(this->aquct_data_.td,
this->aquct_data_.pi, td_plus_dt);
return std::make_pair(PItd, PItdprime);
} }
inline Scalar dpai(int idx) Scalar dpai(const int idx) const
{ {
Scalar dp = this->pa0_ Scalar dp = this->pa0_
+ this->rhow_.at(idx).value() * this->gravity_() * (this->cell_depth_.at(idx) - this->aquiferDepth()) + this->rhow_.at(idx).value() * this->gravity_() * (this->cell_depth_.at(idx) - this->aquiferDepth())
@ -114,24 +137,28 @@ protected:
} }
// This function implements Eqs 5.8 and 5.9 of the EclipseTechnicalDescription // This function implements Eqs 5.8 and 5.9 of the EclipseTechnicalDescription
inline void calculateEqnConstants(Scalar& a, Scalar& b, const int idx, const Simulator& simulator) std::pair<Scalar, Scalar>
calculateEqnConstants(const int idx, const Simulator& simulator)
{ {
const Scalar td_plus_dt = (simulator.timeStepSize() + simulator.time()) / this->Tc_; const Scalar td_plus_dt = (simulator.timeStepSize() + simulator.time()) / this->Tc_;
const Scalar td = simulator.time() / this->Tc_; this->dimensionless_time_ = simulator.time() / this->Tc_;
Scalar PItdprime = 0.;
Scalar PItd = 0.; const auto [PItd, PItdprime] = this->getInfluenceTableValues(td_plus_dt);
getInfluenceTableValues(PItd, PItdprime, td_plus_dt);
a = 1.0 / this->Tc_ * ((beta_ * dpai(idx)) - (this->fluxValue_ * PItdprime)) / (PItd - td * PItdprime); const auto denom = this->Tc_ * (PItd - this->dimensionless_time_*PItdprime);
b = beta_ / (this->Tc_ * (PItd - td * PItdprime)); const auto a = (this->beta_*dpai(idx) - this->fluxValue_*PItdprime) / denom;
const auto b = this->beta_ / denom;
return std::make_pair(a, b);
} }
// This function implements Eq 5.7 of the EclipseTechnicalDescription // This function implements Eq 5.7 of the EclipseTechnicalDescription
inline void calculateInflowRate(int idx, const Simulator& simulator) override inline void calculateInflowRate(int idx, const Simulator& simulator) override
{ {
Scalar a, b; const auto [a, b] = this->calculateEqnConstants(idx, simulator);
calculateEqnConstants(a, b, idx, simulator);
this->Qai_.at(idx) this->Qai_.at(idx) = this->alphai_.at(idx) *
= this->alphai_.at(idx) * (a - b * (this->pressure_current_.at(idx) - this->pressure_previous_.at(idx))); (a - b*(this->pressure_current_.at(idx) - this->pressure_previous_.at(idx)));
} }
inline void calculateAquiferConstants() override inline void calculateAquiferConstants() override