From 9514e44d139c655ae2e40ece52fa4ec53654360e Mon Sep 17 00:00:00 2001 From: Kai Bao Date: Fri, 23 Jun 2017 15:54:05 +0200 Subject: [PATCH] adding function computeConnectionDensities to StandardWell which is copied from WellDensitySegmented, while avoid using of the Wells struct. TODO: find a better place to put these long functions. --- opm/autodiff/StandardWell.hpp | 8 +++ opm/autodiff/StandardWell_impl.hpp | 98 ++++++++++++++++++++++++++++++ 2 files changed, 106 insertions(+) diff --git a/opm/autodiff/StandardWell.hpp b/opm/autodiff/StandardWell.hpp index 2c0f44258..61be4d166 100644 --- a/opm/autodiff/StandardWell.hpp +++ b/opm/autodiff/StandardWell.hpp @@ -219,6 +219,14 @@ namespace Opm std::vector& rsmax_perf, std::vector& rvmax_perf, std::vector& surf_dens_perf) const; + + // TODO: not total sure whether it is a good idea to put here + // the major reason to put here is to avoid the usage of Wells struct + void computeConnectionDensities(const std::vector& perfComponentRates, + const std::vector& b_perf, + const std::vector& rsmax_perf, + const std::vector& rvmax_perf, + const std::vector& surf_dens_perf); }; } diff --git a/opm/autodiff/StandardWell_impl.hpp b/opm/autodiff/StandardWell_impl.hpp index 6a2bfb8a5..2f6159b03 100644 --- a/opm/autodiff/StandardWell_impl.hpp +++ b/opm/autodiff/StandardWell_impl.hpp @@ -1524,4 +1524,102 @@ namespace Opm } } + + + + + template + void + StandardWell:: + computeConnectionDensities(const std::vector& perfComponentRates, + const std::vector& b_perf, + const std::vector& rsmax_perf, + const std::vector& rvmax_perf, + const std::vector& surf_dens_perf) + { + // Verify that we have consistent input. + const int np = numberOfPhases(); + const int nperf = numberOfPerforations(); + const int num_comp = numComponents(); + const PhaseUsage* phase_usage = phase_usage_; + + // 1. Compute the flow (in surface volume units for each + // component) exiting up the wellbore from each perforation, + // taking into account flow from lower in the well, and + // in/out-flow at each perforation. + std::vector q_out_perf(nperf*num_comp); + + // TODO: investigate whether we should use the following techniques to calcuate the composition of flows in the wellbore + // Iterate over well perforations from bottom to top. + for (int perf = nperf - 1; perf >= 0; --perf) { + for (int component = 0; component < num_comp; ++component) { + if (perf == nperf - 1) { + // This is the bottom perforation. No flow from below. + q_out_perf[perf*num_comp+ component] = 0.0; + } else { + // Set equal to flow from below. + q_out_perf[perf*num_comp + component] = q_out_perf[(perf+1)*num_comp + component]; + } + // Subtract outflow through perforation. + q_out_perf[perf*num_comp + component] -= perfComponentRates[perf*num_comp + component]; + } + } + + // 2. Compute the component mix at each perforation as the + // absolute values of the surface rates divided by their sum. + // Then compute volume ratios (formation factors) for each perforation. + // Finally compute densities for the segments associated with each perforation. + const int gaspos = phase_usage->phase_pos[BlackoilPhases::Vapour]; + const int oilpos = phase_usage->phase_pos[BlackoilPhases::Liquid]; + std::vector mix(num_comp,0.0); + std::vector x(num_comp); + std::vector surf_dens(num_comp); + std::vector dens(nperf); + + for (int perf = 0; perf < nperf; ++perf) { + // Find component mix. + const double tot_surf_rate = std::accumulate(q_out_perf.begin() + num_comp*perf, + q_out_perf.begin() + num_comp*(perf+1), 0.0); + if (tot_surf_rate != 0.0) { + for (int component = 0; component < num_comp; ++component) { + mix[component] = std::fabs(q_out_perf[perf*num_comp + component]/tot_surf_rate); + } + } else { + // No flow => use well specified fractions for mix. + for (int phase = 0; phase < np; ++phase) { + mix[phase] = compFrac()[phase]; + } + // intialize 0.0 for comIdx >= np; + } + // Compute volume ratio. + x = mix; + double rs = 0.0; + double rv = 0.0; + if (!rsmax_perf.empty() && mix[oilpos] > 0.0) { + rs = std::min(mix[gaspos]/mix[oilpos], rsmax_perf[perf]); + } + if (!rvmax_perf.empty() && mix[gaspos] > 0.0) { + rv = std::min(mix[oilpos]/mix[gaspos], rvmax_perf[perf]); + } + if (rs != 0.0) { + // Subtract gas in oil from gas mixture + x[gaspos] = (mix[gaspos] - mix[oilpos]*rs)/(1.0 - rs*rv); + } + if (rv != 0.0) { + // Subtract oil in gas from oil mixture + x[oilpos] = (mix[oilpos] - mix[gaspos]*rv)/(1.0 - rs*rv);; + } + double volrat = 0.0; + for (int component = 0; component < num_comp; ++component) { + volrat += x[component] / b_perf[perf*num_comp+ component]; + } + for (int component = 0; component < num_comp; ++component) { + surf_dens[component] = surf_dens_perf[perf*num_comp+ component]; + } + + // Compute segment density. + perf_densities_[perf] = std::inner_product(surf_dens.begin(), surf_dens.end(), mix.begin(), 0.0) / volrat; + } + } + }