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 OPM_TIMEBLOCK(assemble);
44 // We assume that calculateExplicitQuantities() and
45 // prepareTimeStep() have been called already for the entire
46 // well model, so we do not need to do it here (when
47 // iterationIdx is 0).
48
49 DeferredLogger local_deferredLogger;
50 this->updateWellControls(local_deferredLogger, domain);
51 this->assembleWellEq(dt, domain, local_deferredLogger);
52}
53
54template<typename TypeTag>
55void
56BlackoilWellModelNldd<TypeTag>::
57assembleWellEq(const double dt,
58 const Domain& domain,
59 DeferredLogger& deferred_logger)
60{
61 OPM_TIMEBLOCK(assembleWellEq);
62 for (const auto& well : wellModel_.localNonshutWells()) {
63 if (this->well_domain().at(well->name()) == domain.index) {
64 well->assembleWellEq(wellModel_.simulator(),
65 dt,
66 wellModel_.wellState(),
67 wellModel_.groupState(),
69 }
70 }
71}
72
73template<typename TypeTag>
74void
75BlackoilWellModelNldd<TypeTag>::
76addWellPressureEquations(PressureMatrix& /*jacobian*/,
77 const BVector& /*weights*/,
78 const bool /*use_well_weights*/,
79 const int /*domainIndex*/) const
80{
81 throw std::logic_error("CPRW is not yet implemented for NLDD subdomains");
82 // To fix this function, rdofs should be the size of the domain, and the nw should be the number of wells in the domain
83 // int nw = this->numLocalWellsEnd(); // should number of wells in the domain
84 // int rdofs = local_num_cells_; // should be the size of the domain
85 // for ( int i = 0; i < nw; i++ ) {
86 // int wdof = rdofs + i;
87 // jacobian[wdof][wdof] = 1.0;// better scaling ?
88 // }
89
90 // for ( const auto& well : well_container_ ) {
91 // if (well_domain_.at(well->name()) == domainIndex) {
92 // weights should be the size of the domain
93 // well->addWellPressureEquations(jacobian, weights, pressureVarIndex, use_well_weights, this->wellState());
94 // }
95 // }
96}
97
98template<typename TypeTag>
99void
100BlackoilWellModelNldd<TypeTag>::
101recoverWellSolutionAndUpdateWellState(const BVector& x,
102 const int domainIdx)
103{
104 // Note: no point in trying to do a parallel gathering
105 // try/catch here, as this function is not called in
106 // parallel but for each individual domain of each rank.
107 DeferredLogger local_deferredLogger;
108 for (const auto& well : wellModel_.localNonshutWells()) {
109 if (this->well_domain().at(well->name()) == domainIdx) {
110 const auto& cells = well->cells();
111 x_local_.resize(cells.size());
112
113 for (size_t i = 0; i < cells.size(); ++i) {
114 x_local_[i] = x[cells[i]];
115 }
116 well->recoverWellSolutionAndUpdateWellState(wellModel_.simulator(),
117 x_local_,
118 wellModel_.wellState(),
120 }
121 }
122 // TODO: avoid losing the logging information that could
123 // be stored in the local_deferredlogger in a parallel case.
124 if (wellModel_.terminalOutput()) {
125 local_deferredLogger.logMessages();
126 }
127}
128
129template<typename TypeTag>
130ConvergenceReport
131BlackoilWellModelNldd<TypeTag>::
132getWellConvergence(const Domain& domain,
133 const std::vector<Scalar>& B_avg,
134 DeferredLogger& local_deferredLogger) const
135{
136 const int iterationIdx = wellModel_.simulator().model().newtonMethod().numIterations();
137 const bool relax_tolerance = iterationIdx > wellModel_.numStrictIterations();
138
139 ConvergenceReport report;
140 for (const auto& well : wellModel_.localNonshutWells()) {
141 if ((this->well_domain().at(well->name()) == domain.index)) {
142 if (well->isOperableAndSolvable() || well->wellIsStopped()) {
143 report += well->getWellConvergence(wellModel_.simulator(),
144 wellModel_.wellState(),
145 B_avg,
148 } else {
149 ConvergenceReport xreport;
150 using CR = ConvergenceReport;
151 xreport.setWellFailed({CR::WellFailure::Type::Unsolvable,
152 CR::Severity::Normal, -1, well->name()});
153 report += xreport;
154 }
155 }
156 }
157
158 // Log debug messages for NaN or too large residuals.
159 if (wellModel_.terminalOutput()) {
160 for (const auto& f : report.wellFailures()) {
161 if (f.severity() == ConvergenceReport::Severity::NotANumber) {
162 local_deferredLogger.debug("NaN residual found with phase " +
163 std::to_string(f.phase()) +
164 " for well " + f.wellName());
165 } else if (f.severity() == ConvergenceReport::Severity::TooLarge) {
166 local_deferredLogger.debug("Too large residual found with phase " +
167 std::to_string(f.phase()) +
168 " for well " + f.wellName());
169 }
170 }
171 }
172 return report;
173}
174
175template<typename TypeTag>
176void
177BlackoilWellModelNldd<TypeTag>::
178updateWellControls(DeferredLogger& deferred_logger,
179 const Domain& domain)
180{
181 OPM_TIMEBLOCK(updateWellControls);
182 if (!wellModel_.wellsActive()) {
183 return;
184 }
185
186 // TODO: decide on and implement an approach to handling of
187 // group controls, network and similar for domain solves.
188
189 // Check only individual well constraints and communicate.
190 for (const auto& well : wellModel_.localNonshutWells()) {
191 if (this->well_domain().at(well->name()) == domain.index) {
192 constexpr auto mode = WellInterface<TypeTag>::IndividualOrGroup::Individual;
193 well->updateWellControl(wellModel_.simulator(),
194 mode,
195 wellModel_.wellState(),
196 wellModel_.groupState(),
198 }
199 }
200}
201
202template <typename TypeTag>
203void
204BlackoilWellModelNldd<TypeTag>::
205setupDomains(const std::vector<Domain>& domains)
206{
207 std::vector<const SubDomainIndices*> genDomains;
208 std::transform(domains.begin(), domains.end(),
209 std::back_inserter(genDomains),
210 [](const auto& domain)
211 { return static_cast<const SubDomainIndices*>(&domain); });
212 this->calcDomains(genDomains);
213}
214
215} // namespace Opm
216
217#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