OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-17 15:32:57 -06:00
Code Issues Packages Projects Releases Wiki Activity

Add trivial well equation when there is no flow

Browse Source 
Fixes case when all perforations in a well are crossflowing while
crossflow is not allowed
This commit is contained in:
Tor Harald Sandve 2016-09-29 15:29:55 +02:00
parent 07339be80e
commit 84c33e0e96
1 changed files with 18 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

26
opm/autodiff/StandardWellsDense.hpp
View File

@ -226,6 +226,7 @@ namespace Opm {
const double volume = 0.002831684659200; // 0.1 cu ft;
for (int w = 0; w < nw; ++w) {
double total_flux = 0.0;
for (int perf = wells().well_connpos[w] ; perf < wells().well_connpos[w+1]; ++perf) {
@ -243,6 +244,7 @@ namespace Opm {
// subtract sum of phase fluxes in the well equations.
resWell_[w][flowPhaseToEbosCompIdx(p1)] -= cq_s[p1].value;
total_flux += cq_s[p1].value;
// assemble the jacobians
for (int p2 = 0; p2 < np; ++p2) {
@ -260,15 +262,23 @@ namespace Opm {
// Store the perforation pressure for later usage.
well_state.perfPress()[perf] = well_state.bhp()[w] + wellPerforationPressureDiffs()[perf];
}
// add vol * dF/dt + Q to the well equations;
for (int p1 = 0; p1 < np; ++p1) {
EvalWell resWell_loc = (wellVolumeFraction(w, p1) - F0_[w + nw*p1]) * volume / dt;
resWell_loc += getQs(w, p1);
for (int p2 = 0; p2 < np; ++p2) {
invDuneD_[w][w][flowPhaseToEbosCompIdx(p1)][flowToEbosPvIdx(p2)] += resWell_loc.derivatives[p2+3];
if (std::abs(total_flux) == 0) { // add a trivial equation for the case when there is no flow
for (int p1 = 0; p1 < np; ++p1) {
invDuneD_[w][w][flowPhaseToEbosCompIdx(p1)][flowToEbosPvIdx(p1)] += 1;
resWell_[w][flowPhaseToEbosCompIdx(p1)] += 0;
}
}
else {
// add vol * dF/dt + Q to the well equations;
for (int p1 = 0; p1 < np; ++p1) {
EvalWell resWell_loc = (wellVolumeFraction(w, p1) - F0_[w + nw*p1]) * volume / dt;
resWell_loc += getQs(w, p1);
for (int p2 = 0; p2 < np; ++p2) {
invDuneD_[w][w][flowPhaseToEbosCompIdx(p1)][flowToEbosPvIdx(p2)] += resWell_loc.derivatives[p2+3];
}
resWell_[w][flowPhaseToEbosCompIdx(p1)] += resWell_loc.value;
}
resWell_[w][flowPhaseToEbosCompIdx(p1)] += resWell_loc.value;
}
}