My Project
Loading...
Searching...
No Matches
BlackoilWellModelNldd_impl.hpp
1/*
2 Copyright 2016 - 2019 SINTEF Digital, Mathematics & Cybernetics.
3 Copyright 2016 - 2018 Equinor ASA.
4 Copyright 2017 Dr. Blatt - HPC-Simulation-Software & Services
5 Copyright 2016 - 2018 Norce AS
6
7 This file is part of the Open Porous Media project (OPM).
8
9 OPM is free software: you can redistribute it and/or modify
10 it under the terms of the GNU General Public License as published by
11 the Free Software Foundation, either version 3 of the License, or
12 (at your option) any later version.
13
14 OPM is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 GNU General Public License for more details.
18
19 You should have received a copy of the GNU General Public License
20 along with OPM. If not, see <http://www.gnu.org/licenses/>.
21*/
22
23#ifndef OPM_BLACKOILWELLMODEL_NLDD_IMPL_HEADER_INCLUDED
24#define OPM_BLACKOILWELLMODEL_NLDD_IMPL_HEADER_INCLUDED
25
26// Improve IDE experience
27#ifndef OPM_BLACKOILWELLMODEL_NLDD_HEADER_INCLUDED
28#include <config.h>
29#include <opm/simulators/wells/BlackoilWellModelNldd.hpp>
30#endif
31
32#include <algorithm>
33
34namespace Opm {
35
36template<typename TypeTag>
37void
38BlackoilWellModelNldd<TypeTag>::
39assemble(const int /*iterationIdx*/,
40 const double dt,
41 const Domain& domain)
42{
43 // We assume that calculateExplicitQuantities() and
44 // prepareTimeStep() have been called already for the entire
45 // well model, so we do not need to do it here (when
46 // iterationIdx is 0).
47
48 DeferredLogger local_deferredLogger;
49 this->updateWellControls(local_deferredLogger, domain);
50 this->assembleWellEq(dt, domain, local_deferredLogger);
51}
52
53template<typename TypeTag>
54void
55BlackoilWellModelNldd<TypeTag>::
56assembleWellEq(const double dt,
57 const Domain& domain,
58 DeferredLogger& deferred_logger)
59{
60 for (const auto& well : wellModel_.localNonshutWells()) {
61 if (this->well_domain().at(well->name()) == domain.index) {
62 well->assembleWellEq(wellModel_.simulator(),
63 dt,
64 wellModel_.wellState(),
65 wellModel_.groupState(),
67 }
68 }
69}
70
71template<typename TypeTag>
72void
73BlackoilWellModelNldd<TypeTag>::
74addWellPressureEquations(PressureMatrix& /*jacobian*/,
75 const BVector& /*weights*/,
76 const bool /*use_well_weights*/,
77 const int /*domainIndex*/) const
78{
79 throw std::logic_error("CPRW is not yet implemented for NLDD subdomains");
80 // To fix this function, rdofs should be the size of the domain, and the nw should be the number of wells in the domain
81 // int nw = this->numLocalWellsEnd(); // should number of wells in the domain
82 // int rdofs = local_num_cells_; // should be the size of the domain
83 // for ( int i = 0; i < nw; i++ ) {
84 // int wdof = rdofs + i;
85 // jacobian[wdof][wdof] = 1.0;// better scaling ?
86 // }
87
88 // for ( const auto& well : well_container_ ) {
89 // if (well_domain_.at(well->name()) == domainIndex) {
90 // weights should be the size of the domain
91 // well->addWellPressureEquations(jacobian, weights, pressureVarIndex, use_well_weights, this->wellState());
92 // }
93 // }
94}
95
96template<typename TypeTag>
97void
98BlackoilWellModelNldd<TypeTag>::
99recoverWellSolutionAndUpdateWellState(const BVector& x,
100 const int domainIdx)
101{
102 // Note: no point in trying to do a parallel gathering
103 // try/catch here, as this function is not called in
104 // parallel but for each individual domain of each rank.
105 DeferredLogger local_deferredLogger;
106 for (const auto& well : wellModel_.localNonshutWells()) {
107 if (this->well_domain().at(well->name()) == domainIdx) {
108 const auto& cells = well->cells();
109 x_local_.resize(cells.size());
110
111 for (size_t i = 0; i < cells.size(); ++i) {
112 x_local_[i] = x[cells[i]];
113 }
114 well->recoverWellSolutionAndUpdateWellState(wellModel_.simulator(),
115 x_local_,
116 wellModel_.wellState(),
118 }
119 }
120 // TODO: avoid losing the logging information that could
121 // be stored in the local_deferredlogger in a parallel case.
122 if (wellModel_.terminalOutput()) {
123 local_deferredLogger.logMessages();
124 }
125}
126
127template<typename TypeTag>
128ConvergenceReport
129BlackoilWellModelNldd<TypeTag>::
130getWellConvergence(const Domain& domain,
131 const std::vector<Scalar>& B_avg,
132 DeferredLogger& local_deferredLogger) const
133{
134 const int iterationIdx = wellModel_.simulator().model().newtonMethod().numIterations();
135 const bool relax_tolerance = iterationIdx > wellModel_.numStrictIterations();
136
137 ConvergenceReport report;
138 for (const auto& well : wellModel_.localNonshutWells()) {
139 if ((this->well_domain().at(well->name()) == domain.index)) {
140 if (well->isOperableAndSolvable() || well->wellIsStopped()) {
141 report += well->getWellConvergence(wellModel_.simulator(),
142 wellModel_.wellState(),
143 B_avg,
146 } else {
147 ConvergenceReport xreport;
148 using CR = ConvergenceReport;
149 xreport.setWellFailed({CR::WellFailure::Type::Unsolvable,
150 CR::Severity::Normal, -1, well->name()});
151 report += xreport;
152 }
153 }
154 }
155
156 // Log debug messages for NaN or too large residuals.
157 if (wellModel_.terminalOutput()) {
158 for (const auto& f : report.wellFailures()) {
159 if (f.severity() == ConvergenceReport::Severity::NotANumber) {
160 local_deferredLogger.debug("NaN residual found with phase " +
161 std::to_string(f.phase()) +
162 " for well " + f.wellName());
163 } else if (f.severity() == ConvergenceReport::Severity::TooLarge) {
164 local_deferredLogger.debug("Too large residual found with phase " +
165 std::to_string(f.phase()) +
166 " for well " + f.wellName());
167 }
168 }
169 }
170 return report;
171}
172
173template<typename TypeTag>
174void
175BlackoilWellModelNldd<TypeTag>::
176updateWellControls(DeferredLogger& deferred_logger,
177 const Domain& domain)
178{
179 if (!wellModel_.wellsActive()) {
180 return;
181 }
182
183 // TODO: decide on and implement an approach to handling of
184 // group controls, network and similar for domain solves.
185
186 // Check only individual well constraints and communicate.
187 for (const auto& well : wellModel_.localNonshutWells()) {
188 if (this->well_domain().at(well->name()) == domain.index) {
189 constexpr auto mode = WellInterface<TypeTag>::IndividualOrGroup::Individual;
190 well->updateWellControl(wellModel_.simulator(),
191 mode,
192 wellModel_.wellState(),
193 wellModel_.groupState(),
195 }
196 }
197}
198
199template <typename TypeTag>
200void
201BlackoilWellModelNldd<TypeTag>::
202setupDomains(const std::vector<Domain>& domains)
203{
204 std::vector<const SubDomainIndices*> genDomains;
205 std::transform(domains.begin(), domains.end(),
206 std::back_inserter(genDomains),
207 [](const auto& domain)
208 { return static_cast<const SubDomainIndices*>(&domain); });
209 this->calcDomains(genDomains);
210}
211
212} // namespace Opm
213
214#endif
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242