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<class Scalar> class PAvgCalculator;
46
64template<class Scalar>
65typename PAvgCalculator<Scalar>::Result
66linearCombination(const Scalar alpha, typename PAvgCalculator<Scalar>::Result x,
67 const Scalar beta , const typename PAvgCalculator<Scalar>::Result& y);
68
72template<class Scalar>
74{
75protected:
76 class Accumulator;
77
78public:
80 class Result
81 {
82 private:
85 friend class Accumulator;
86
88 friend Result
89 linearCombination<Scalar>(const Scalar alpha, Result x,
90 const Scalar beta , const Result& y);
91
92 public:
94 enum class WBPMode
95 {
96 WBP, //< Connecting cells
97 WBP4, //< Immediate neighbours
98 WBP5, //< Connecting cells and immediate neighbours
99 WBP9, //< Connecting cells, immediate, and diagonal neighbours
100 };
101
106 Scalar value(const WBPMode type) const
107 {
108 return this->wbp_[this->index(type)];
109 }
110
111 private:
113 static constexpr auto NumModes =
114 static_cast<std::size_t>(WBPMode::WBP9) + 1;
115
117 using WBPStore = std::array<Scalar, NumModes>;
118
120 WBPStore wbp_{};
121
127 Result& set(const WBPMode type, const Scalar wbp)
128 {
129 this->wbp_[this->index(type)] = wbp;
130
131 return *this;
132 }
133
138 constexpr typename WBPStore::size_type index(const WBPMode mode) const
139 {
140 return static_cast<typename WBPStore::size_type>(mode);
141 }
142 };
143
146 {
147 public:
154 {
155 this->wb_ = &wbSrc;
156 return *this;
157 }
158
165 {
166 this->wc_ = &wcSrc;
167 return *this;
168 }
169
172 {
173 return *this->wb_;
174 }
175
178 {
179 return *this->wc_;
180 }
181
182 private:
184 const PAvgDynamicSourceData<Scalar>* wb_{nullptr};
185
187 const PAvgDynamicSourceData<Scalar>* wc_{nullptr};
188 };
189
196 PAvgCalculator(const GridDims& cellIndexMap,
197 const WellConnections& connections);
198
201
212 void pruneInactiveWBPCells(const std::vector<bool>& isActive);
213
226 void inferBlockAveragePressures(const Sources& sources,
227 const PAvg& controls,
228 const Scalar gravity,
229 const Scalar refDepth);
230
233 const std::vector<std::size_t>& allWBPCells() const
234 {
235 return this->contributingCells_;
236 }
237
248 std::vector<std::size_t> allWellConnections() const;
249
256 {
257 return this->averagePressures_;
258 }
259
260protected:
263 {
264 public:
269 using LocalRunningAverages = std::array<Scalar, 8>;
270
272 Accumulator();
273
276
280 Accumulator(const Accumulator& rhs);
281
286
290 Accumulator& operator=(const Accumulator& rhs);
291
296
302 Accumulator& addCentre(const Scalar weight,
303 const Scalar press);
304
311 Accumulator& addRectangular(const Scalar weight,
312 const Scalar press);
313
319 Accumulator& addDiagonal(const Scalar weight,
320 const Scalar press);
321
330 Accumulator& add(const Scalar weight,
331 const Accumulator& other);
332
334 void prepareAccumulation();
335
337 void prepareContribution();
338
348 void commitContribution(const Scalar innerWeight = -1.0);
349
350 // Please note that member functions \c getRunningAverages() and \c
351 // assignRunningAverages() are concessions to parallel/MPI runs, and
352 // especially for simulation runs with distributed wells. In this
353 // situation we need a way to access, communicate/collect/sum, and
354 // assign partial results. Moreover, the \c LocalRunningAverages
355 // should be treated opaquely apart from applying a global reduction
356 // operation. In other words, the intended/expected use case is
357 //
358 // Accumulator a{}
359 // ...
360 // auto avg = a.getRunningAverages()
361 // MPI_Allreduce(avg, MPI_SUM)
362 // a.assignRunningAverages(avg)
363 //
364 // Any other use is probably a bug and the above is the canonical
365 // implementation of member function collectGlobalContributions() in
366 // MPI aware sub classes of PAvgCalculator.
367
370
375
379 Result getFinalResult() const;
380
381 private:
383 class Impl;
384
386 std::unique_ptr<Impl> pImpl_;
387 };
388
391
394
395private:
397 using ContrIndexType = std::vector<std::size_t>::size_type;
398
403 using SetupMap = std::unordered_map<std::size_t, ContrIndexType>;
404
406 enum class NeighbourKind
407 {
409 Rectangular,
410
412 Diagonal,
413 };
414
417 struct PAvgConnection
418 {
427 PAvgConnection(const Scalar ctf_arg,
428 const Scalar depth_arg,
429 const ContrIndexType cell_arg)
430 : ctf (ctf_arg)
431 , depth(depth_arg)
432 , cell (cell_arg)
433 {}
434
436 Scalar ctf{};
437
439 Scalar depth{};
440
442 ContrIndexType cell{};
443
447 std::vector<ContrIndexType> rectNeighbours{};
448
452 std::vector<ContrIndexType> diagNeighbours{};
453 };
454
459 typename std::vector<PAvgConnection>::size_type numInputConns_{};
460
463 std::vector<PAvgConnection> connections_{};
464
466 std::vector<typename std::vector<PAvgConnection>::size_type> openConns_{};
467
476 std::vector<typename std::vector<PAvgConnection>::size_type> inputConn_{};
477
480 std::vector<std::size_t> contributingCells_{};
481
485 Result averagePressures_{};
486
500 void addConnection(const GridDims& cellIndexMap,
501 const Connection& conn,
502 SetupMap& setupHelperMap);
503
514 void pruneInactiveConnections(const std::vector<bool>& isActive);
515
531 void accumulateLocalContributions(const Sources& sources,
532 const PAvg& controls,
533 const Scalar gravity,
534 const Scalar refDepth);
535
543 virtual void collectGlobalContributions();
544
552 void assignResults(const PAvg& controls);
553
570 void addNeighbour(std::optional<std::size_t> neighbour,
571 NeighbourKind neighbourKind,
572 SetupMap& setupHelperMap);
573
577 std::size_t lastConnsCell() const;
578
590 void addNeighbours_X(const GridDims& grid, SetupMap& setupHelperMap);
591
603 void addNeighbours_Y(const GridDims& grid, SetupMap& setupHelperMap);
604
616 void addNeighbours_Z(const GridDims& grid, SetupMap& setupHelperMap);
617
657 template <typename ConnIndexMap, typename CTFPressureWeightFunction>
658 void accumulateLocalContributions(const Sources& sources,
659 const PAvg& controls,
660 const std::vector<Scalar>& connDP,
661 ConnIndexMap connIndex,
662 CTFPressureWeightFunction ctfPressWeight);
663
695 template <typename ConnIndexMap>
696 void accumulateLocalContributions(const Sources& sources,
697 const PAvg& controls,
698 const std::vector<Scalar>& connDP,
699 ConnIndexMap&& connIndex);
700
712 void accumulateLocalContribOpen(const Sources& sources,
713 const PAvg& controls,
714 const std::vector<Scalar>& connDP);
715
727 void accumulateLocalContribAll(const Sources& sources,
728 const PAvg& controls,
729 const std::vector<Scalar>& connDP);
730
761 template <typename ConnIndexMap>
762 std::vector<Scalar>
763 connectionPressureOffsetWell(const std::size_t nconn,
764 const Sources& sources,
765 const Scalar gravity,
766 const Scalar refDepth,
767 ConnIndexMap connIndex) const;
768
800 template <typename ConnIndexMap>
801 std::vector<Scalar>
802 connectionPressureOffsetRes(const std::size_t nconn,
803 const Sources& sources,
804 const Scalar gravity,
805 const Scalar refDepth,
806 ConnIndexMap connIndex) const;
807
826 std::vector<Scalar>
827 connectionPressureOffset(const Sources& sources,
828 const PAvg& controls,
829 const Scalar gravity,
830 const Scalar refDepth) const;
831};
832
833} // namespace Opm
834
835#endif // PAVG_CALCULATOR_HPP
Accumulate weighted running averages of cell contributions to WBP.
Definition PAvgCalculator.hpp:263
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
Result getFinalResult() const
Calculate final WBP results from individual contributions.
Definition PAvgCalculator.cpp:598
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
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:269
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
Result of block-averaging well pressure procedure.
Definition PAvgCalculator.hpp:81
Scalar value(const WBPMode type) const
Retrieve numerical value of specific block-averaged well pressure.
Definition PAvgCalculator.hpp:106
WBPMode
Kind of block-averaged well pressure.
Definition PAvgCalculator.hpp:95
References to source contributions owned by other party.
Definition PAvgCalculator.hpp:146
const PAvgDynamicSourceData< Scalar > & wellConns() const
Get read-only access to connection-level contributions.
Definition PAvgCalculator.hpp:177
Sources & wellConns(const PAvgDynamicSourceData< Scalar > &wcSrc)
Provide reference to connection-level contributions (pressure, pore-volume, mixture density) owned by...
Definition PAvgCalculator.hpp:164
Sources & wellBlocks(const PAvgDynamicSourceData< Scalar > &wbSrc)
Provide reference to cell-level contributions (pressure, pore-volume, mixture density) owned by other...
Definition PAvgCalculator.hpp:153
const PAvgDynamicSourceData< Scalar > & wellBlocks() const
Get read-only access to cell-level contributions.
Definition PAvgCalculator.hpp:171
Facility for deriving well-level pressure values from selected block-averaging procedures.
Definition PAvgCalculator.hpp:74
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:721
Accumulator accumPV_
Average pressures weighted by pore-volume.
Definition PAvgCalculator.hpp:393
Accumulator accumCTF_
Average pressures weighted by connection transmissibility factor.
Definition PAvgCalculator.hpp:390
const Result & averagePressures() const
Block-average pressure derived from selection of source cells.
Definition PAvgCalculator.hpp:255
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:233
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:707
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