20#ifndef PAVG_CALCULATOR_HPP
21#define PAVG_CALCULATOR_HPP
29#include <unordered_map>
38class PAvgDynamicSourceData;
65 const double beta ,
const Result& y);
83 return this->wbp_[this->index(type)];
88 static constexpr auto NumModes =
89 static_cast<std::size_t
>(WBPMode::WBP9) + 1;
92 using WBPStore = std::array<double, NumModes>;
102 Result& set(
const WBPMode type,
const double wbp)
104 this->wbp_[this->index(type)] = wbp;
113 constexpr WBPStore::size_type index(
const WBPMode mode)
const
115 return static_cast<WBPStore::size_type
>(mode);
172 const WellConnections& connections);
202 const PAvg& controls,
203 const double gravity,
204 const double refDepth);
210 return this->contributingCells_;
232 return this->averagePressures_;
361 std::unique_ptr<Impl> pImpl_;
372 using ContrIndexType = std::vector<std::size_t>::size_type;
378 using SetupMap = std::unordered_map<std::size_t, ContrIndexType>;
381 enum class NeighbourKind
392 struct PAvgConnection
402 PAvgConnection(
const double ctf_arg,
403 const double depth_arg,
404 const ContrIndexType cell_arg)
417 ContrIndexType cell{};
422 std::vector<ContrIndexType> rectNeighbours{};
427 std::vector<ContrIndexType> diagNeighbours{};
434 std::vector<PAvgConnection>::size_type numInputConns_{};
438 std::vector<PAvgConnection> connections_{};
441 std::vector<std::vector<PAvgConnection>::size_type> openConns_{};
451 std::vector<std::vector<PAvgConnection>::size_type> inputConn_{};
455 std::vector<std::size_t> contributingCells_{};
460 Result averagePressures_{};
475 void addConnection(
const GridDims& cellIndexMap,
476 const Connection& conn,
477 SetupMap& setupHelperMap);
489 void pruneInactiveConnections(
const std::vector<bool>& isActive);
506 void accumulateLocalContributions(
const Sources& sources,
507 const PAvg& controls,
508 const double gravity,
509 const double refDepth);
518 virtual void collectGlobalContributions();
527 void assignResults(
const PAvg& controls);
545 void addNeighbour(std::optional<std::size_t> neighbour,
546 NeighbourKind neighbourKind,
547 SetupMap& setupHelperMap);
552 std::size_t lastConnsCell()
const;
565 void addNeighbours_X(
const GridDims& grid, SetupMap& setupHelperMap);
578 void addNeighbours_Y(
const GridDims& grid, SetupMap& setupHelperMap);
591 void addNeighbours_Z(
const GridDims& grid, SetupMap& setupHelperMap);
632 template <
typename ConnIndexMap,
typename CTFPressureWeightFunction>
633 void accumulateLocalContributions(
const Sources& sources,
634 const PAvg& controls,
635 const std::vector<double>& connDP,
636 ConnIndexMap connIndex,
637 CTFPressureWeightFunction ctfPressWeight);
670 template <
typename ConnIndexMap>
671 void accumulateLocalContributions(
const Sources& sources,
672 const PAvg& controls,
673 const std::vector<double>& connDP,
674 ConnIndexMap&& connIndex);
687 void accumulateLocalContribOpen(
const Sources& sources,
688 const PAvg& controls,
689 const std::vector<double>& connDP);
702 void accumulateLocalContribAll(
const Sources& sources,
703 const PAvg& controls,
704 const std::vector<double>& connDP);
736 template <
typename ConnIndexMap>
738 connectionPressureOffsetWell(
const std::size_t nconn,
739 const Sources& sources,
740 const double gravity,
741 const double refDepth,
742 ConnIndexMap connIndex)
const;
775 template <
typename ConnIndexMap>
777 connectionPressureOffsetRes(
const std::size_t nconn,
778 const Sources& sources,
779 const double gravity,
780 const double refDepth,
781 ConnIndexMap connIndex)
const;
802 connectionPressureOffset(
const Sources& sources,
803 const PAvg& controls,
804 const double gravity,
805 const double refDepth)
const;
825PAvgCalculator::Result
826linearCombination(
const double alpha, PAvgCalculator::Result x,
827 const double beta ,
const PAvgCalculator::Result& y);
Accumulate weighted running averages of cell contributions to WBP.
Definition PAvgCalculator.hpp:238
Accumulator & operator=(const Accumulator &rhs)
Assignment operator.
Definition PAvgCalculator.cpp:507
void prepareContribution()
Zero out/clear WBP term buffer.
Definition PAvgCalculator.cpp:557
Accumulator()
Constructor.
Definition PAvgCalculator.cpp:492
LocalRunningAverages getRunningAverages() const
Get buffer of intermediate, local results.
Definition PAvgCalculator.cpp:568
Result getFinalResult() const
Calculate final WBP results from individual contributions.
Definition PAvgCalculator.cpp:581
Accumulator & addDiagonal(const double weight, const double press)
Add contribution from diagonal, level 2 neighbouring cell.
Definition PAvgCalculator.cpp:537
Accumulator & addCentre(const double weight, const double press)
Add contribution from centre/connecting cell.
Definition PAvgCalculator.cpp:521
void prepareAccumulation()
Zero out/clear WBP result buffer.
Definition PAvgCalculator.cpp:552
std::array< double, 8 > LocalRunningAverages
Collection of running averages and their associate weights.
Definition PAvgCalculator.hpp:244
void commitContribution(const double innerWeight=-1.0)
Accumulate current source term into result buffer whilst applying any user-prescribed term weighting.
Definition PAvgCalculator.cpp:562
Accumulator & addRectangular(const double weight, const double press)
Add contribution from direct, rectangular, level 1 neighbouring cell.
Definition PAvgCalculator.cpp:529
void assignRunningAverages(const LocalRunningAverages &avg)
Assign coalesced/global contributions.
Definition PAvgCalculator.cpp:575
~Accumulator()
Destructor.
Accumulator & add(const double weight, const Accumulator &other)
Add contribution from other accumulator.
Definition PAvgCalculator.cpp:545
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.
Definition PAvgCalculator.cpp:1146
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
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.
Definition PAvgCalculator.cpp:616
Accumulator accumCTF_
Average pressures weighted by connection transmissibility factor.
Definition PAvgCalculator.hpp:365
PAvgCalculator(const GridDims &cellIndexMap, const WellConnections &connections)
Constructor.
Definition PAvgCalculator.cpp:600
std::vector< std::size_t > allWellConnections() const
List all reservoir connections that potentially contribute to this block-averaging pressure calculati...
Definition PAvgCalculator.cpp:694
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...
Definition PAvgCalculator.cpp:682
virtual ~PAvgCalculator()
Destructor.
Dynamic source data for block-average pressure calculations.
Definition PAvgDynamicSourceData.hpp:35
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30