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;
38template<class Scalar> class PAvgDynamicSourceData;
39class WellConnections;
40
41} // namespace Opm
42
43namespace Opm {
44
45template <typename Scalar>
46class PAvgCalculatorResult;
47
65template <typename Scalar>
66PAvgCalculatorResult<Scalar>
67linearCombination(const Scalar alpha, PAvgCalculatorResult<Scalar> x,
68 const Scalar beta , const PAvgCalculatorResult<Scalar>& y);
69
71template <typename Scalar>
73{
74private:
76 template <typename T>
79 const T beta , const PAvgCalculatorResult<T>& y);
80
81public:
83 enum class WBPMode
84 {
85 WBP, //< Connecting cells
86 WBP4, //< Immediate neighbours
87 WBP5, //< Connecting cells and immediate neighbours
88 WBP9, //< Connecting cells, immediate, and diagonal neighbours
89 };
90
96 PAvgCalculatorResult& set(const WBPMode type, const Scalar wbp)
97 {
98 this->wbp_[this->index(type)] = wbp;
99
100 return *this;
101 }
102
107 Scalar value(const WBPMode type) const
108 {
109 return this->wbp_[this->index(type)];
110 }
111
112private:
114 static constexpr auto NumModes =
115 static_cast<std::size_t>(WBPMode::WBP9) + 1;
116
118 using WBPStore = std::array<Scalar, NumModes>;
119
121 WBPStore wbp_{};
122
127 constexpr typename WBPStore::size_type index(const WBPMode mode) const
128 {
129 return static_cast<typename WBPStore::size_type>(mode);
130 }
131};
132
136template<class Scalar>
138{
139protected:
140 class Accumulator;
141
142public:
145 {
146 public:
153 {
154 this->wb_ = &wbSrc;
155 return *this;
156 }
157
164 {
165 this->wc_ = &wcSrc;
166 return *this;
167 }
168
171 {
172 return *this->wb_;
173 }
174
177 {
178 return *this->wc_;
179 }
180
181 private:
183 const PAvgDynamicSourceData<Scalar>* wb_{nullptr};
184
186 const PAvgDynamicSourceData<Scalar>* wc_{nullptr};
187 };
188
195 PAvgCalculator(const GridDims& cellIndexMap,
196 const WellConnections& connections);
197
200
211 void pruneInactiveWBPCells(const std::vector<bool>& isActive);
212
225 void inferBlockAveragePressures(const Sources& sources,
226 const PAvg& controls,
227 const Scalar gravity,
228 const Scalar refDepth);
229
232 const std::vector<std::size_t>& allWBPCells() const
233 {
234 return this->contributingCells_;
235 }
236
247 std::vector<std::size_t> allWellConnections() const;
248
255 {
256 return this->averagePressures_;
257 }
258
259protected:
262 {
263 public:
268 using LocalRunningAverages = std::array<Scalar, 8>;
269
271 Accumulator();
272
275
279 Accumulator(const Accumulator& rhs);
280
285
289 Accumulator& operator=(const Accumulator& rhs);
290
295
301 Accumulator& addCentre(const Scalar weight,
302 const Scalar press);
303
310 Accumulator& addRectangular(const Scalar weight,
311 const Scalar press);
312
318 Accumulator& addDiagonal(const Scalar weight,
319 const Scalar press);
320
329 Accumulator& add(const Scalar weight,
330 const Accumulator& other);
331
333 void prepareAccumulation();
334
336 void prepareContribution();
337
347 void commitContribution(const Scalar innerWeight = -1.0);
348
349 // Please note that member functions \c getRunningAverages() and \c
350 // assignRunningAverages() are concessions to parallel/MPI runs, and
351 // especially for simulation runs with distributed wells. In this
352 // situation we need a way to access, communicate/collect/sum, and
353 // assign partial results. Moreover, the \c LocalRunningAverages
354 // should be treated opaquely apart from applying a global reduction
355 // operation. In other words, the intended/expected use case is
356 //
357 // Accumulator a{}
358 // ...
359 // auto avg = a.getRunningAverages()
360 // MPI_Allreduce(avg, MPI_SUM)
361 // a.assignRunningAverages(avg)
362 //
363 // Any other use is probably a bug and the above is the canonical
364 // implementation of member function collectGlobalContributions() in
365 // MPI aware sub classes of PAvgCalculator.
366
369
374
379
380 private:
382 class Impl;
383
385 std::unique_ptr<Impl> pImpl_;
386 };
387
390
393
394private:
396 using ContrIndexType = std::vector<std::size_t>::size_type;
397
402 using SetupMap = std::unordered_map<std::size_t, ContrIndexType>;
403
405 enum class NeighbourKind
406 {
408 Rectangular,
409
411 Diagonal,
412 };
413
416 struct PAvgConnection
417 {
426 PAvgConnection(const Scalar ctf_arg,
427 const ContrIndexType cell_arg)
428 : ctf { ctf_arg }
429 , cell { cell_arg }
430 {}
431
433 Scalar ctf{};
434
436 ContrIndexType cell{};
437
441 std::vector<ContrIndexType> rectNeighbours{};
442
446 std::vector<ContrIndexType> diagNeighbours{};
447 };
448
453 typename std::vector<PAvgConnection>::size_type numInputConns_{};
454
457 std::vector<PAvgConnection> connections_{};
458
460 std::vector<typename std::vector<PAvgConnection>::size_type> openConns_{};
461
469 std::vector<typename std::vector<PAvgConnection>::size_type> inputConn_{};
470
473 std::vector<std::size_t> contributingCells_{};
474
478 PAvgCalculatorResult<Scalar> averagePressures_{};
479
493 void addConnection(const GridDims& cellIndexMap,
494 const Connection& conn,
495 SetupMap& setupHelperMap);
496
507 void pruneInactiveConnections(const std::vector<bool>& isActive);
508
524 void accumulateLocalContributions(const Sources& sources,
525 const PAvg& controls,
526 const Scalar gravity,
527 const Scalar refDepth);
528
536 virtual void collectGlobalContributions();
537
545 void assignResults(const PAvg& controls);
546
563 void addNeighbour(std::optional<std::size_t> neighbour,
564 NeighbourKind neighbourKind,
565 SetupMap& setupHelperMap);
566
570 std::size_t lastConnsCell() const;
571
583 void addNeighbours_X(const GridDims& grid, SetupMap& setupHelperMap);
584
596 void addNeighbours_Y(const GridDims& grid, SetupMap& setupHelperMap);
597
609 void addNeighbours_Z(const GridDims& grid, SetupMap& setupHelperMap);
610
656 template <typename ConnIndexMap, typename CTFPressureWeightFunction>
657 void accumulateLocalContributions(const Sources& sources,
658 const PAvg& controls,
659 const Scalar gravity,
660 const Scalar refDepth,
661 const std::vector<Scalar>& connDensity,
662 ConnIndexMap connIndex,
663 CTFPressureWeightFunction ctfPressWeight);
664
702 template <typename ConnIndexMap>
703 void accumulateLocalContributions(const Sources& sources,
704 const PAvg& controls,
705 const Scalar gravity,
706 const Scalar refDepth,
707 const std::vector<Scalar>& connDensity,
708 ConnIndexMap&& connIndex);
709
727 void accumulateLocalContribOpen(const Sources& sources,
728 const PAvg& controls,
729 const Scalar gravity,
730 const Scalar refDepth,
731 const std::vector<Scalar>& connDP);
732
750 void accumulateLocalContribAll(const Sources& sources,
751 const PAvg& controls,
752 const Scalar gravity,
753 const Scalar refDepth,
754 const std::vector<Scalar>& connDensity);
755
781 template <typename ConnIndexMap>
782 std::vector<Scalar>
783 connectionDensityWell(const std::size_t nconn,
784 const Sources& sources,
785 ConnIndexMap connIndex) const;
786
813 template <typename ConnIndexMap>
814 std::vector<Scalar>
815 connectionDensityRes(const std::size_t nconn,
816 const Sources& sources,
817 ConnIndexMap connIndex) const;
818
833 std::vector<Scalar>
834 connectionDensity(const Sources& sources,
835 const PAvg& controls,
836 const Scalar gravity) const;
837};
838
839} // namespace Opm
840
841#endif // PAVG_CALCULATOR_HPP
Result of block-averaging well pressure procedure.
Definition PAvgCalculator.hpp:73
WBPMode
Kind of block-averaged well pressure.
Definition PAvgCalculator.hpp:84
PAvgCalculatorResult & set(const WBPMode type, const Scalar wbp)
Assign single block-averaged pressure result.
Definition PAvgCalculator.hpp:96
friend PAvgCalculatorResult< T > linearCombination(const T alpha, PAvgCalculatorResult< T > x, const T beta, const PAvgCalculatorResult< T > &y)
Grant internal data member access to combination function.
Scalar value(const WBPMode type) const
Retrieve numerical value of specific block-averaged well pressure.
Definition PAvgCalculator.hpp:107
Accumulate weighted running averages of cell contributions to WBP.
Definition PAvgCalculator.hpp:262
Accumulator()
Constructor.
Definition PAvgCalculator.cpp:495
void prepareAccumulation()
Zero out/clear WBP result buffer.
Definition PAvgCalculator.cpp:565
Accumulator & addCentre(const Scalar weight, const Scalar press)
Add contribution from centre/connecting cell.
Definition PAvgCalculator.cpp:530
Accumulator & addDiagonal(const Scalar weight, const Scalar press)
Add contribution from diagonal, level 2 neighbouring cell.
Definition PAvgCalculator.cpp:548
void prepareContribution()
Zero out/clear WBP term buffer.
Definition PAvgCalculator.cpp:571
void assignRunningAverages(const LocalRunningAverages &avg)
Assign coalesced/global contributions.
Definition PAvgCalculator.cpp:591
LocalRunningAverages getRunningAverages() const
Get buffer of intermediate, local results.
Definition PAvgCalculator.cpp:584
Accumulator & addRectangular(const Scalar weight, const Scalar press)
Add contribution from direct, rectangular, level 1 neighbouring cell.
Definition PAvgCalculator.cpp:539
PAvgCalculatorResult< Scalar > getFinalResult() const
Calculate final WBP results from individual contributions.
Definition PAvgCalculator.cpp:598
Accumulator & add(const Scalar weight, const Accumulator &other)
Add contribution from other accumulator.
Definition PAvgCalculator.cpp:557
std::array< Scalar, 8 > LocalRunningAverages
Collection of running averages and their associate weights.
Definition PAvgCalculator.hpp:268
void commitContribution(const Scalar innerWeight=-1.0)
Accumulate current source term into result buffer whilst applying any user-prescribed term weighting.
Definition PAvgCalculator.cpp:577
Accumulator & operator=(const Accumulator &rhs)
Assignment operator.
Definition PAvgCalculator.cpp:514
References to source contributions owned by other party.
Definition PAvgCalculator.hpp:145
const PAvgDynamicSourceData< Scalar > & wellConns() const
Get read-only access to connection-level contributions.
Definition PAvgCalculator.hpp:176
Sources & wellConns(const PAvgDynamicSourceData< Scalar > &wcSrc)
Provide reference to connection-level contributions (pressure, pore-volume, mixture density) owned by...
Definition PAvgCalculator.hpp:163
Sources & wellBlocks(const PAvgDynamicSourceData< Scalar > &wbSrc)
Provide reference to cell-level contributions (pressure, pore-volume, mixture density) owned by other...
Definition PAvgCalculator.hpp:152
const PAvgDynamicSourceData< Scalar > & wellBlocks() const
Get read-only access to cell-level contributions.
Definition PAvgCalculator.hpp:170
Facility for deriving well-level pressure values from selected block-averaging procedures.
Definition PAvgCalculator.hpp:138
virtual ~PAvgCalculator()
Destructor.
void pruneInactiveWBPCells(const std::vector< bool > &isActive)
Finish construction by pruning inactive cells.
Definition PAvgCalculator.cpp:636
std::vector< std::size_t > allWellConnections() const
List all reservoir connections that potentially contribute to this block-averaging pressure calculati...
Definition PAvgCalculator.cpp:729
const PAvgCalculatorResult< Scalar > & averagePressures() const
Block-average pressure derived from selection of source cells.
Definition PAvgCalculator.hpp:254
Accumulator accumPV_
Average pressures weighted by pore-volume.
Definition PAvgCalculator.hpp:392
Accumulator accumCTF_
Average pressures weighted by connection transmissibility factor.
Definition PAvgCalculator.hpp:389
PAvgCalculator(const GridDims &cellIndexMap, const WellConnections &connections)
Constructor.
Definition PAvgCalculator.cpp:618
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:232
void inferBlockAveragePressures(const Sources &sources, const PAvg &controls, const Scalar gravity, const Scalar refDepth)
Compute block-average well-level pressure values from collection of source contributions and user-def...
Definition PAvgCalculator.cpp:715
Dynamic source data for block-average pressure calculations.
Definition PAvgDynamicSourceData.hpp:36
Definition PAvg.hpp:30
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30