My Project
Loading...
Searching...
No Matches
PAvgCalculator.hpp
1/*
2 Copyright 2020, 2023 Equinor ASA.
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef PAVG_CALCULATOR_HPP
21#define PAVG_CALCULATOR_HPP
22
23#include <array>
24#include <cstddef>
25#include <functional>
26#include <memory>
27#include <optional>
28#include <string>
29#include <unordered_map>
30#include <utility>
31#include <vector>
32
33namespace Opm {
34
35class Connection;
36class GridDims;
37class PAvg;
38class PAvgDynamicSourceData;
39class WellConnections;
40
41} // namespace Opm
42
43namespace Opm {
44
49{
50protected:
51 class Accumulator;
52
53public:
55 class Result
56 {
57 private:
60 friend class Accumulator;
61
63 friend Result
64 linearCombination(const double alpha, Result x,
65 const double beta , const Result& y);
66
67 public:
69 enum class WBPMode
70 {
71 WBP, //< Connecting cells
72 WBP4, //< Immediate neighbours
73 WBP5, //< Connecting cells and immediate neighbours
74 WBP9, //< Connecting cells, immediate, and diagonal neighbours
75 };
76
81 double value(const WBPMode type) const
82 {
83 return this->wbp_[this->index(type)];
84 }
85
86 private:
88 static constexpr auto NumModes =
89 static_cast<std::size_t>(WBPMode::WBP9) + 1;
90
92 using WBPStore = std::array<double, NumModes>;
93
95 WBPStore wbp_{};
96
102 Result& set(const WBPMode type, const double wbp)
103 {
104 this->wbp_[this->index(type)] = wbp;
105
106 return *this;
107 }
108
113 constexpr WBPStore::size_type index(const WBPMode mode) const
114 {
115 return static_cast<WBPStore::size_type>(mode);
116 }
117 };
118
121 {
122 public:
129 {
130 this->wb_ = &wbSrc;
131 return *this;
132 }
133
140 {
141 this->wc_ = &wcSrc;
142 return *this;
143 }
144
147 {
148 return *this->wb_;
149 }
150
153 {
154 return *this->wc_;
155 }
156
157 private:
159 const PAvgDynamicSourceData* wb_{nullptr};
160
162 const PAvgDynamicSourceData* wc_{nullptr};
163 };
164
171 PAvgCalculator(const GridDims& cellIndexMap,
172 const WellConnections& connections);
173
176
187 void pruneInactiveWBPCells(const std::vector<bool>& isActive);
188
202 const PAvg& controls,
203 const double gravity,
204 const double refDepth);
205
208 const std::vector<std::size_t>& allWBPCells() const
209 {
210 return this->contributingCells_;
211 }
212
223 std::vector<std::size_t> allWellConnections() const;
224
231 {
232 return this->averagePressures_;
233 }
234
235protected:
238 {
239 public:
244 using LocalRunningAverages = std::array<double, 8>;
245
248
251
256
261
266
271
277 Accumulator& addCentre(const double weight,
278 const double press);
279
286 Accumulator& addRectangular(const double weight,
287 const double press);
288
294 Accumulator& addDiagonal(const double weight,
295 const double press);
296
305 Accumulator& add(const double weight,
306 const Accumulator& other);
307
310
313
323 void commitContribution(const double innerWeight = -1.0);
324
325 // Please note that member functions \c getRunningAverages() and \c
326 // assignRunningAverages() are concessions to parallel/MPI runs, and
327 // especially for simulation runs with distributed wells. In this
328 // situation we need a way to access, communicate/collect/sum, and
329 // assign partial results. Moreover, the \c LocalRunningAverages
330 // should be treated opaquely apart from applying a global reduction
331 // operation. In other words, the intended/expected use case is
332 //
333 // Accumulator a{}
334 // ...
335 // auto avg = a.getRunningAverages()
336 // MPI_Allreduce(avg, MPI_SUM)
337 // a.assignRunningAverages(avg)
338 //
339 // Any other use is probably a bug and the above is the canonical
340 // implementation of member function collectGlobalContributions() in
341 // MPI aware sub classes of PAvgCalculator.
342
345
350
355
356 private:
358 class Impl;
359
361 std::unique_ptr<Impl> pImpl_;
362 };
363
366
369
370private:
372 using ContrIndexType = std::vector<std::size_t>::size_type;
373
378 using SetupMap = std::unordered_map<std::size_t, ContrIndexType>;
379
381 enum class NeighbourKind
382 {
384 Rectangular,
385
387 Diagonal,
388 };
389
392 struct PAvgConnection
393 {
402 PAvgConnection(const double ctf_arg,
403 const double depth_arg,
404 const ContrIndexType cell_arg)
405 : ctf (ctf_arg)
406 , depth(depth_arg)
407 , cell (cell_arg)
408 {}
409
411 double ctf{};
412
414 double depth{};
415
417 ContrIndexType cell{};
418
422 std::vector<ContrIndexType> rectNeighbours{};
423
427 std::vector<ContrIndexType> diagNeighbours{};
428 };
429
432 std::vector<PAvgConnection> connections_{};
433
435 std::vector<std::vector<PAvgConnection>::size_type> openConns_{};
436
439 std::vector<std::size_t> contributingCells_{};
440
444 Result averagePressures_{};
445
459 void addConnection(const GridDims& cellIndexMap,
460 const Connection& conn,
461 SetupMap& setupHelperMap);
462
478 void accumulateLocalContributions(const Sources& sources,
479 const PAvg& controls,
480 const double gravity,
481 const double refDepth);
482
490 virtual void collectGlobalContributions();
491
499 void assignResults(const PAvg& controls);
500
517 void addNeighbour(std::optional<std::size_t> neighbour,
518 NeighbourKind neighbourKind,
519 SetupMap& setupHelperMap);
520
524 std::size_t lastConnsCell() const;
525
537 void addNeighbours_X(const GridDims& grid, SetupMap& setupHelperMap);
538
550 void addNeighbours_Y(const GridDims& grid, SetupMap& setupHelperMap);
551
563 void addNeighbours_Z(const GridDims& grid, SetupMap& setupHelperMap);
564
604 template <typename ConnIndexMap, typename CTFPressureWeightFunction>
605 void accumulateLocalContributions(const Sources& sources,
606 const PAvg& controls,
607 const std::vector<double>& connDP,
608 ConnIndexMap connIndex,
609 CTFPressureWeightFunction ctfPressWeight);
610
642 template <typename ConnIndexMap>
643 void accumulateLocalContributions(const Sources& sources,
644 const PAvg& controls,
645 const std::vector<double>& connDP,
646 ConnIndexMap&& connIndex);
647
659 void accumulateLocalContribOpen(const Sources& sources,
660 const PAvg& controls,
661 const std::vector<double>& connDP);
662
674 void accumulateLocalContribAll(const Sources& sources,
675 const PAvg& controls,
676 const std::vector<double>& connDP);
677
708 template <typename ConnIndexMap>
709 std::vector<double>
710 connectionPressureOffsetWell(const std::size_t nconn,
711 const Sources& sources,
712 const double gravity,
713 const double refDepth,
714 ConnIndexMap connIndex) const;
715
747 template <typename ConnIndexMap>
748 std::vector<double>
749 connectionPressureOffsetRes(const std::size_t nconn,
750 const Sources& sources,
751 const double gravity,
752 const double refDepth,
753 ConnIndexMap connIndex) const;
754
773 std::vector<double>
774 connectionPressureOffset(const Sources& sources,
775 const PAvg& controls,
776 const double gravity,
777 const double refDepth) const;
778};
779
797PAvgCalculator::Result
799 const double beta , const PAvgCalculator::Result& y);
800
801} // namespace Opm
802
803#endif // PAVG_CALCULATOR_HPP
Definition GridDims.hpp:31
Accumulate weighted running averages of cell contributions to WBP.
Definition PAvgCalculator.hpp:238
LocalRunningAverages getRunningAverages() const
Get buffer of intermediate, local results.
void prepareContribution()
Zero out/clear WBP term buffer.
Accumulator(const Accumulator &rhs)
Copy constructor.
Accumulator & add(const double weight, const Accumulator &other)
Add contribution from other accumulator.
void prepareAccumulation()
Zero out/clear WBP result buffer.
std::array< double, 8 > LocalRunningAverages
Collection of running averages and their associate weights.
Definition PAvgCalculator.hpp:244
Accumulator & operator=(Accumulator &&rhs)
Move assignment operator.
Accumulator & addDiagonal(const double weight, const double press)
Add contribution from diagonal, level 2 neighbouring cell.
Accumulator(Accumulator &&rhs)
Move constructor.
Accumulator & addCentre(const double weight, const double press)
Add contribution from centre/connecting cell.
void commitContribution(const double innerWeight=-1.0)
Accumulate current source term into result buffer whilst applying any user-prescribed term weighting.
Accumulator & operator=(const Accumulator &rhs)
Assignment operator.
Accumulator & addRectangular(const double weight, const double press)
Add contribution from direct, rectangular, level 1 neighbouring cell.
void assignRunningAverages(const LocalRunningAverages &avg)
Assign coalesced/global contributions.
Result getFinalResult() const
Calculate final WBP results from individual contributions.
Result of block-averaging well pressure procedure.
Definition PAvgCalculator.hpp:56
friend Result linearCombination(const double alpha, Result x, const double beta, const Result &y)
Grant internal data member access to combination function.
double value(const WBPMode type) const
Retrieve numerical value of specific block-averaged well pressure.
Definition PAvgCalculator.hpp:81
WBPMode
Kind of block-averaged well pressure.
Definition PAvgCalculator.hpp:70
References to source contributions owned by other party.
Definition PAvgCalculator.hpp:121
const PAvgDynamicSourceData & wellConns() const
Get read-only access to connection-level contributions.
Definition PAvgCalculator.hpp:152
Sources & wellBlocks(const PAvgDynamicSourceData &wbSrc)
Provide reference to cell-level contributions (pressure, pore-volume, mixture density) owned by other...
Definition PAvgCalculator.hpp:128
const PAvgDynamicSourceData & wellBlocks() const
Get read-only access to cell-level contributions.
Definition PAvgCalculator.hpp:146
Sources & wellConns(const PAvgDynamicSourceData &wcSrc)
Provide reference to connection-level contributions (pressure, pore-volume, mixture density) owned by...
Definition PAvgCalculator.hpp:139
Facility for deriving well-level pressure values from selected block-averaging procedures.
Definition PAvgCalculator.hpp:49
const std::vector< std::size_t > & allWBPCells() const
List of all cells, global indices in natural ordering, that contribute to the block-average pressures...
Definition PAvgCalculator.hpp:208
virtual ~PAvgCalculator()
Destructor.
const Result & averagePressures() const
Block-average pressure derived from selection of source cells.
Definition PAvgCalculator.hpp:230
Accumulator accumPV_
Average pressures weighted by pore-volume.
Definition PAvgCalculator.hpp:368
void pruneInactiveWBPCells(const std::vector< bool > &isActive)
Finish construction by pruning inactive cells.
Accumulator accumCTF_
Average pressures weighted by connection transmissibility factor.
Definition PAvgCalculator.hpp:365
PAvgCalculator(const GridDims &cellIndexMap, const WellConnections &connections)
Constructor.
std::vector< std::size_t > allWellConnections() const
List all reservoir connections that potentially contribute to this block-averaging pressure calculati...
void inferBlockAveragePressures(const Sources &sources, const PAvg &controls, const double gravity, const double refDepth)
Compute block-average well-level pressure values from collection of source contributions and user-def...
Dynamic source data for block-average pressure calculations.
Definition PAvgDynamicSourceData.hpp:35
Definition PAvg.hpp:30
Definition WellConnections.hpp:45
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
PAvgCalculator::Result linearCombination(const double alpha, PAvgCalculator::Result x, const double beta, const PAvgCalculator::Result &y)
Form linear combination of WBP result objects.