From 6db925189cfe0a1b2904fc5046063e47f6fb3f2c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Ove=20S=C3=A6vareid?= Date: Tue, 6 Jul 2021 12:54:16 +0200 Subject: [PATCH] Some regularisation of tracer rate reporting. --- ebos/ecltracermodel.hh | 29 ++++++++++++++++++++++++++++- 1 file changed, 28 insertions(+), 1 deletion(-) diff --git a/ebos/ecltracermodel.hh b/ebos/ecltracermodel.hh index c8cf58a60..0887ac0a8 100644 --- a/ebos/ecltracermodel.hh +++ b/ebos/ecltracermodel.hh @@ -390,6 +390,11 @@ protected: if (well.getStatus() == Well::Status::SHUT) continue; + if (!well.isProducer()) //Injection rates already reported during assembly + continue; + + Scalar rateWellPos = 0.0; + Scalar rateWellNeg = 0.0; std::array cartesianCoordinate; for (auto& connection : well.getConnections()) { @@ -402,11 +407,33 @@ protected: const size_t cartIdx = simulator_.vanguard().cartesianIndex(cartesianCoordinate); const int I = this->cartToGlobal_[cartIdx]; Scalar rate = simulator_.problem().wellModel().well(well.name())->volumetricSurfaceRateForConnection(I, tr.phaseIdx_); - if (rate < 0 && well.isProducer()) { //Injection rates already reported during assembly + if (rate < 0) { + rateWellNeg += rate; for (int tIdx =0; tIdx < tr.numTracer(); ++tIdx) { this->wellTracerRate_.at(std::make_pair(well.name(),this->tracerNames_[tr.idx_[tIdx]])) += rate*tr.concentration_[tIdx][I]; } } + else { + rateWellPos += rate; + } + } + + Scalar rateWellTotal = rateWellNeg + rateWellPos; + + // TODO: Some inconsistencies here that perhaps should be clarified. The "offical" rate as reported below is + // occasionally significant different from the sum over connections (as calculated above). Only observed + // for small values, neglible for the rate itself, but matters when used to calculate tracer concentrations. + std::size_t well_index = simulator_.problem().wellModel().wellState().wellIndex(well.name()); + Scalar official_well_rate_total = simulator_.problem().wellModel().wellState().wellRates(well_index)[tr.phaseIdx_]; + + rateWellTotal = official_well_rate_total; + + if (rateWellTotal > rateWellNeg) { // Cross flow + const Scalar bucketPrDay = 10.0/(1000.*3600.*24.); // ... keeps (some) trouble away + const Scalar factor = (rateWellTotal < -bucketPrDay) ? rateWellTotal/rateWellNeg : 0.0; + for (int tIdx =0; tIdx < tr.numTracer(); ++tIdx) { + this->wellTracerRate_.at(std::make_pair(well.name(),this->tracerNames_[tr.idx_[tIdx]])) *= factor; + } } } }