Merge pull request #4691 from bska/wbp-comm-api-change

Prepare for Revised Implementation of WBPn
This commit is contained in:
Bård Skaflestad 2023-06-09 15:15:32 +02:00 committed by GitHub
commit 9fde77079e
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
8 changed files with 307 additions and 257 deletions

View File

@ -531,55 +531,48 @@ public:
class PackUnPackWBPData : public P2PCommunicatorType::DataHandleInterface
{
const std::map<std::size_t, double>& localWBPData_;
std::map<std::size_t, double>& globalWBPValues_;
const data::WellBlockAveragePressures& localWBPData_;
data::WellBlockAveragePressures& globalWBPValues_;
public:
PackUnPackWBPData(const std::map<std::size_t, double>& localWBPData,
std::map<std::size_t, double>& globalWBPValues,
bool isIORank)
: localWBPData_(localWBPData)
PackUnPackWBPData(const data::WellBlockAveragePressures& localWBPData,
data::WellBlockAveragePressures& globalWBPValues,
const bool isIORank)
: localWBPData_ (localWBPData)
, globalWBPValues_(globalWBPValues)
{
if (isIORank) {
MessageBufferType buffer;
pack(0, buffer);
// pass a dummyLink to satisfy virtual class
int dummyLink = -1;
unpack(dummyLink, buffer);
if (! isIORank) {
return;
}
MessageBufferType buffer;
pack(0, buffer);
// Pass a dummy link to satisfy base class API requirement
const int dummyLink = -1;
unpack(dummyLink, buffer);
}
// pack all data associated with link
// Pack all data associated with link
void pack(int link, MessageBufferType& buffer)
{
// we should only get one link
if (link != 0)
throw std::logic_error("link in method pack is not 0 as expected");
// write all block data
unsigned int size = localWBPData_.size();
buffer.write(size);
for (const auto& [global_index, wbp_value] : localWBPData_) {
buffer.write(global_index);
buffer.write(wbp_value);
// We should only get one link
if (link != 0) {
throw std::logic_error {
"link (" + std::to_string(link) +
") in method pack() is not 0 as expected"
};
}
// Write all local, per-well, WBP data
this->localWBPData_.write(buffer);
}
// unpack all data associated with link
void unpack(int /*link*/, MessageBufferType& buffer)
// Unpack all data associated with link
void unpack([[maybe_unused]] const int link,
MessageBufferType& buffer)
{
// read all block data
unsigned int size = 0;
buffer.read(size);
for (size_t i = 0; i < size; ++i) {
std::size_t idx;
double data;
buffer.read(idx);
buffer.read(data);
globalWBPValues_[idx] = data;
}
this->globalWBPValues_.read(buffer);
}
};
@ -984,21 +977,21 @@ CollectDataToIORank(const Grid& grid, const EquilGrid* equilGrid,
template <class Grid, class EquilGrid, class GridView>
void CollectDataToIORank<Grid,EquilGrid,GridView>::
collect(const data::Solution& localCellData,
collect(const data::Solution& localCellData,
const std::map<std::pair<std::string, int>, double>& localBlockData,
const std::map<std::size_t, double>& localWBPData,
const data::Wells& localWellData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const data::Aquifers& localAquiferData,
const WellTestState& localWellTestState,
const EclInterRegFlowMap& localInterRegFlows,
const data::Wells& localWellData,
const data::WellBlockAveragePressures& localWBPData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const data::Aquifers& localAquiferData,
const WellTestState& localWellTestState,
const EclInterRegFlowMap& localInterRegFlows,
const std::array<std::pair<std::string, std::pair<std::vector<int>, std::vector<double>>>, 3>& localFlowsn,
const std::array<std::pair<std::string, std::pair<std::vector<int>, std::vector<double>>>, 3>& localFloresn)
{
globalCellData_ = {};
globalBlockData_.clear();
globalWBPData_.clear();
globalWellData_.clear();
globalWBPData_.values.clear();
globalGroupAndNetworkData_.clear();
globalAquiferData_.clear();
globalWellTestState_.clear();
@ -1055,8 +1048,8 @@ collect(const data::Solution& localCellData,
PackUnPackWBPData packUnpackWBPData {
localWBPData,
this->globalWBPData_,
this->isIORank()
this->globalWBPData_,
this->isIORank()
};
PackUnPackAquiferData packUnpackAquiferData {
@ -1110,28 +1103,37 @@ template <class Grid, class EquilGrid, class GridView>
int CollectDataToIORank<Grid,EquilGrid,GridView>::
localIdxToGlobalIdx(unsigned localIdx) const
{
if (!isParallel())
if (!isParallel()) {
return localIdx;
}
if (localIdxToGlobalIdx_.empty())
if (this->localIdxToGlobalIdx_.empty()) {
throw std::logic_error("index map is not created on this rank");
}
if (localIdx > localIdxToGlobalIdx_.size())
if (localIdx > this->localIdxToGlobalIdx_.size()) {
throw std::logic_error("local index is outside map range");
}
return localIdxToGlobalIdx_[localIdx];
return this->localIdxToGlobalIdx_[localIdx];
}
template <class Grid, class EquilGrid, class GridView>
bool CollectDataToIORank<Grid,EquilGrid,GridView>::
isCartIdxOnThisRank(int cartIdx) const
{
if (!isParallel())
if (! this->isParallel()) {
return true;
}
assert(!needsReordering);
auto candidate = std::lower_bound(sortedCartesianIdx_.begin(), sortedCartesianIdx_.end(), cartIdx);
return (candidate != sortedCartesianIdx_.end() && *candidate == cartIdx);
assert (! needsReordering);
auto candidate = std::lower_bound(this->sortedCartesianIdx_.begin(),
this->sortedCartesianIdx_.end(),
cartIdx);
return (candidate != sortedCartesianIdx_.end())
&& (*candidate == cartIdx);
}

View File

@ -25,15 +25,21 @@
#include <opm/output/data/Aquifer.hpp>
#include <opm/output/data/Cells.hpp>
#include <opm/output/data/Groups.hpp>
#include <opm/output/data/Solution.hpp>
#include <opm/output/data/Wells.hpp>
#include <opm/output/data/Groups.hpp>
#include <opm/input/eclipse/Schedule/Well/WellTestState.hpp>
#include <opm/grid/common/p2pcommunicator.hh>
#include <ebos/eclinterregflows.hh>
#include <array>
#include <map>
#include <set>
#include <string>
#include <type_traits>
#include <utility>
#include <vector>
@ -67,20 +73,17 @@ public:
const std::set<std::string>& fipRegionsInterregFlow = {});
// gather solution to rank 0 for EclipseWriter
void collect(const data::Solution& localCellData,
void collect(const data::Solution& localCellData,
const std::map<std::pair<std::string, int>, double>& localBlockData,
const std::map<std::size_t, double>& localWBPData,
const data::Wells& localWellData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const data::Aquifers& localAquiferData,
const WellTestState& localWellTestState,
const EclInterRegFlowMap& interRegFlows,
const data::Wells& localWellData,
const data::WellBlockAveragePressures& localWBPData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const data::Aquifers& localAquiferData,
const WellTestState& localWellTestState,
const EclInterRegFlowMap& interRegFlows,
const std::array<std::pair<std::string, std::pair<std::vector<int>, std::vector<double>>>, 3>& localFlowsn,
const std::array<std::pair<std::string, std::pair<std::vector<int>, std::vector<double>>>, 3>& localFloresn);
const std::map<std::size_t, double>& globalWBPData() const
{ return this->globalWBPData_; }
const std::map<std::pair<std::string, int>, double>& globalBlockData() const
{ return globalBlockData_; }
@ -90,8 +93,8 @@ public:
const data::Wells& globalWellData() const
{ return globalWellData_; }
const WellTestState& globalWellTestState() const
{ return this->globalWellTestState_; }
const data::WellBlockAveragePressures& globalWBPData() const
{ return this->globalWBPData_; }
const data::GroupAndNetworkValues& globalGroupAndNetworkData() const
{ return globalGroupAndNetworkData_; }
@ -99,6 +102,9 @@ public:
const data::Aquifers& globalAquiferData() const
{ return globalAquiferData_; }
const WellTestState& globalWellTestState() const
{ return this->globalWellTestState_; }
EclInterRegFlowMap& globalInterRegFlows()
{ return this->globalInterRegFlows_; }
@ -144,8 +150,8 @@ protected:
std::vector<int> globalRanks_;
data::Solution globalCellData_;
std::map<std::pair<std::string, int>, double> globalBlockData_;
std::map<std::size_t, double> globalWBPData_;
data::Wells globalWellData_;
data::WellBlockAveragePressures globalWBPData_;
data::GroupAndNetworkValues globalGroupAndNetworkData_;
data::Aquifers globalAquiferData_;
WellTestState globalWellTestState_;

View File

@ -234,11 +234,6 @@ public:
return this->interRegionFlows_.wantInterRegflowSummary();
}
const std::map<std::size_t, double>& getWBPData() const
{
return this->wbpData_;
}
const std::map<std::pair<std::string, int>, double>& getBlockData()
{
return blockData_;
@ -552,7 +547,6 @@ protected:
std::map<size_t, Scalar> waterConnectionSaturations_;
std::map<size_t, Scalar> gasConnectionSaturations_;
std::map<std::pair<std::string, int>, double> blockData_;
std::map<std::size_t , double> wbpData_;
std::optional<Inplace> initialInplace_;
bool local_data_valid_;

View File

@ -22,12 +22,14 @@
*/
#include <config.h>
#include <ebos/eclgenericwriter.hh>
#include <opm/grid/CpGrid.hpp>
#include <opm/grid/cpgrid/GridHelpers.hpp>
#include <opm/grid/polyhedralgrid.hh>
#include <opm/grid/utility/cartesianToCompressed.hpp>
#if HAVE_DUNE_ALUGRID
#include "eclalugridvanguard.hh"
#include <dune/alugrid/grid.hh>
@ -40,13 +42,14 @@
#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
#include <opm/input/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
#include <opm/input/eclipse/Schedule/Action/State.hpp>
#include <opm/input/eclipse/Schedule/Schedule.hpp>
#include <opm/input/eclipse/Schedule/SummaryState.hpp>
#include <opm/input/eclipse/Schedule/UDQ/UDQConfig.hpp>
#include <opm/input/eclipse/Schedule/UDQ/UDQState.hpp>
#include <opm/input/eclipse/Schedule/Well/PAvgCalculatorCollection.hpp>
#include <opm/input/eclipse/Schedule/Well/WellMatcher.hpp>
#include <opm/input/eclipse/Units/UnitSystem.hpp>
#include <dune/grid/common/mcmgmapper.hh>
@ -65,7 +68,12 @@
#include <mpi.h>
#endif
#include <algorithm>
#include <array>
#include <cassert>
#include <functional>
#include <map>
#include <memory>
#include <string>
#include <unordered_map>
#include <utility>
@ -169,22 +177,22 @@ struct EclWriteTasklet : public Opm::TaskletInterface
, reportStepNum_(reportStepNum)
, isSubStep_(isSubStep)
, secondsElapsed_(secondsElapsed)
, restartValue_(restartValue)
, restartValue_(std::move(restartValue))
, writeDoublePrecision_(writeDoublePrecision)
{ }
{}
// callback to eclIO serial writeTimeStep method
void run()
{
eclIO_.writeTimeStep(actionState_,
wtestState_,
summaryState_,
udqState_,
reportStepNum_,
isSubStep_,
secondsElapsed_,
restartValue_,
writeDoublePrecision_);
this->eclIO_.writeTimeStep(this->actionState_,
this->wtestState_,
this->summaryState_,
this->udqState_,
this->reportStepNum_,
this->isSubStep_,
this->secondsElapsed_,
std::move(this->restartValue_),
this->writeDoublePrecision_);
}
};
@ -210,38 +218,30 @@ EclGenericWriter(const Schedule& schedule,
cartMapper,
equilCartMapper,
summaryConfig.fip_regions_interreg_flow())
, grid_(grid)
, gridView_(gridView)
, schedule_(schedule)
, eclState_(eclState)
, summaryConfig_(summaryConfig)
, cartMapper_(cartMapper)
, grid_ (grid)
, gridView_ (gridView)
, schedule_ (schedule)
, eclState_ (eclState)
, summaryConfig_ (summaryConfig)
, cartMapper_ (cartMapper)
, equilCartMapper_(equilCartMapper)
, equilGrid_(equilGrid)
, equilGrid_ (equilGrid)
{
if (collectToIORank_.isIORank()) {
eclIO_.reset(new EclipseIO(eclState_,
UgGridHelpers::createEclipseGrid(*equilGrid, eclState_.getInputGrid()),
schedule_,
summaryConfig_, "", enableEsmry));
if (this->collectToIORank_.isIORank()) {
this->eclIO_ = std::make_unique<EclipseIO>
(this->eclState_,
UgGridHelpers::createEclipseGrid(*equilGrid, eclState_.getInputGrid()),
this->schedule_, this->summaryConfig_, "", enableEsmry);
}
const auto& wbp_calculators = eclIO_->summary().wbp_calculators( schedule.size() - 1 );
wbp_index_list_ = wbp_calculators.index_list();
}
if (collectToIORank_.isParallel()) {
const auto& comm = grid_.comm();
unsigned long size = wbp_index_list_.size();
comm.broadcast(&size, 1, collectToIORank_.ioRank);
if (!collectToIORank_.isIORank())
wbp_index_list_.resize( size );
comm.broadcast(wbp_index_list_.data(), size, collectToIORank_.ioRank);
}
// create output thread if enabled and rank is I/O rank
// async output is enabled by default if pthread are enabled
int numWorkerThreads = 0;
if (enableAsyncOutput && collectToIORank_.isIORank())
if (enableAsyncOutput && collectToIORank_.isIORank()) {
numWorkerThreads = 1;
taskletRunner_.reset(new TaskletRunner(numWorkerThreads));
}
this->taskletRunner_.reset(new TaskletRunner(numWorkerThreads));
}
template<class Grid, class EquilGrid, class GridView, class ElementMapper, class Scalar>
@ -257,12 +257,18 @@ void EclGenericWriter<Grid,EquilGrid,GridView,ElementMapper,Scalar>::
writeInit(const std::function<unsigned int(unsigned int)>& map)
{
if (collectToIORank_.isIORank()) {
std::map<std::string, std::vector<int> > integerVectors;
if (collectToIORank_.isParallel())
std::map<std::string, std::vector<int>> integerVectors;
if (collectToIORank_.isParallel()) {
integerVectors.emplace("MPI_RANK", collectToIORank_.globalRanks());
}
auto cartMap = cartesianToCompressed(equilGrid_->size(0), UgGridHelpers::globalCell(*equilGrid_));
eclIO_->writeInitial(computeTrans_(cartMap, map), integerVectors, exportNncStructure_(cartMap, map));
eclIO_->writeInitial(computeTrans_(cartMap, map),
integerVectors,
exportNncStructure_(cartMap, map));
}
#if HAVE_MPI
if (collectToIORank_.isParallel()) {
const auto& comm = grid_.comm();
@ -450,8 +456,9 @@ doWriteOutput(const int reportStepNum,
const bool needsReordering = this->collectToIORank_.doesNeedReordering();
RestartValue restartValue {
(isParallel || needsReordering) ? this->collectToIORank_.globalCellData()
: std::move(localCellData),
(isParallel || needsReordering)
? this->collectToIORank_.globalCellData()
: std::move(localCellData),
isParallel ? this->collectToIORank_.globalWellData()
: std::move(localWellData),
@ -515,32 +522,32 @@ doWriteOutput(const int reportStepNum,
template<class Grid, class EquilGrid, class GridView, class ElementMapper, class Scalar>
void EclGenericWriter<Grid,EquilGrid,GridView,ElementMapper,Scalar>::
evalSummary(const int reportStepNum,
const Scalar curTime,
const std::map<std::size_t, double>& wbpData,
const data::Wells& localWellData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const std::map<int,data::AquiferData>& localAquiferData,
evalSummary(const int reportStepNum,
const Scalar curTime,
const data::Wells& localWellData,
const data::WellBlockAveragePressures& localWBPData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const std::map<int,data::AquiferData>& localAquiferData,
const std::map<std::pair<std::string, int>, double>& blockData,
const std::map<std::string, double>& miscSummaryData,
const std::map<std::string, std::vector<double>>& regionData,
const Inplace& inplace,
const Inplace& initialInPlace,
const EclInterRegFlowMap& interRegionFlowMap,
SummaryState& summaryState,
UDQState& udqState)
const std::map<std::string, double>& miscSummaryData,
const std::map<std::string, std::vector<double>>& regionData,
const Inplace& inplace,
const Inplace& initialInPlace,
const EclInterRegFlowMap& interRegFlows,
SummaryState& summaryState,
UDQState& udqState)
{
if (collectToIORank_.isIORank()) {
const auto& summary = eclIO_->summary();
auto wbp_calculators = summary.wbp_calculators(reportStepNum);
for (const auto& [global_index, pressure] : wbpData)
wbp_calculators.add_pressure( global_index, pressure );
const auto& wellData = this->collectToIORank_.isParallel()
? this->collectToIORank_.globalWellData()
: localWellData;
const auto& wbpData = this->collectToIORank_.isParallel()
? this->collectToIORank_.globalWBPData()
: localWBPData;
const auto& groupAndNetworkData = this->collectToIORank_.isParallel()
? this->collectToIORank_.globalGroupAndNetworkData()
: localGroupAndNetworkData;
@ -553,24 +560,29 @@ evalSummary(const int reportStepNum,
reportStepNum,
curTime,
wellData,
wbpData,
groupAndNetworkData,
miscSummaryData,
initialInPlace,
inplace,
wbp_calculators,
regionData,
blockData,
aquiferData,
getInterRegFlowsAsMap(interRegionFlowMap));
getInterRegFlowsAsMap(interRegFlows));
// Off-by-one-fun: The reportStepNum argument corresponds to the
// report step these results will be written to, whereas the
// argument to UDQ function evaluation corresponds to the report
// step we are currently on.
auto udq_step = reportStepNum - 1;
const auto& udq_config = schedule_.getUDQConfig(udq_step);
udq_config.eval( udq_step, schedule_.wellMatcher(udq_step), summaryState, udqState);
const auto udq_step = reportStepNum - 1;
this->schedule_.getUDQConfig(udq_step)
.eval(udq_step,
this->schedule_.wellMatcher(udq_step),
summaryState,
udqState);
}
#if HAVE_MPI
if (collectToIORank_.isParallel()) {
EclMpiSerializer ser(grid_.comm());

View File

@ -116,10 +116,10 @@ protected:
data::GroupAndNetworkValues&& localGroupAndNetworkData,
data::Aquifers&& localAquiferData,
WellTestState&& localWTestState,
const Action::State& actionState,
const UDQState& udqState,
const SummaryState& summaryState,
const std::vector<Scalar>& thresholdPressure,
const Action::State& actionState,
const UDQState& udqState,
const SummaryState& summaryState,
const std::vector<Scalar>& thresholdPressure,
Scalar curTime,
Scalar nextStepSize,
bool doublePrecision,
@ -128,20 +128,20 @@ protected:
bool isFloresn,
std::array<std::pair<std::string, std::pair<std::vector<int>, std::vector<double>>>, 3>&& floresn);
void evalSummary(int reportStepNum,
Scalar curTime,
const std::map<std::size_t, double>& wbpData,
const data::Wells& localWellData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const std::map<int,data::AquiferData>& localAquiferData,
void evalSummary(int reportStepNum,
Scalar curTime,
const data::Wells& localWellData,
const data::WellBlockAveragePressures& localWBPData,
const data::GroupAndNetworkValues& localGroupAndNetworkData,
const std::map<int,data::AquiferData>& localAquiferData,
const std::map<std::pair<std::string, int>, double>& blockData,
const std::map<std::string, double>& miscSummaryData,
const std::map<std::string, std::vector<double>>& regionData,
const Inplace& inplace,
const Inplace& initialInPlace,
const EclInterRegFlowMap& interRegionFlowMap,
SummaryState& summaryState,
UDQState& udqState);
const std::map<std::string, double>& miscSummaryData,
const std::map<std::string, std::vector<double>>& regionData,
const Inplace& inplace,
const Inplace& initialInPlace,
const EclInterRegFlowMap& interRegFlows,
SummaryState& summaryState,
UDQState& udqState);
CollectDataToIORankType collectToIORank_;
const Grid& grid_;
@ -156,7 +156,6 @@ protected:
const Dune::CartesianIndexMapper<Grid>& cartMapper_;
const Dune::CartesianIndexMapper<EquilGrid>* equilCartMapper_;
const EquilGrid* equilGrid_;
std::vector<std::size_t> wbp_index_list_;
SimulatorReportSingle sub_step_report_;
SimulatorReport simulation_report_;
mutable std::vector<NNCdata> outputNnc_;

View File

@ -138,7 +138,6 @@ class EclOutputBlackOilModule : public EclGenericOutputBlackoilModule<GetPropTyp
public:
template <class CollectDataToIORankType>
EclOutputBlackOilModule(const Simulator& simulator,
const std::vector<std::size_t>& wbp_index_list,
const CollectDataToIORankType& collectToIORank)
: BaseType(simulator.vanguard().eclState(),
simulator.vanguard().schedule(),
@ -159,22 +158,17 @@ public:
this->createLocalRegion_(region_pair.second);
}
auto isCartIdxOnThisRank = [&collectToIORank](int idx)
{
return collectToIORank.isCartIdxOnThisRank(idx);
};
auto isCartIdxOnThisRank = [&collectToIORank](const int idx) {
return collectToIORank.isCartIdxOnThisRank(idx);
};
this->setupBlockData(isCartIdxOnThisRank);
for (const auto& global_index : wbp_index_list) {
if (collectToIORank.isCartIdxOnThisRank(global_index - 1)) {
this->wbpData_.emplace(global_index, 0.0);
}
}
this->forceDisableFipOutput_ =
EWOMS_GET_PARAM(TypeTag, bool, ForceDisableFluidInPlaceOutput);
this->forceDisableFipOutput_ = EWOMS_GET_PARAM(TypeTag, bool, ForceDisableFluidInPlaceOutput);
this->forceDisableFipresvOutput_ = EWOMS_GET_PARAM(TypeTag, bool, ForceDisableResvFluidInPlaceOutput);
this->forceDisableFipresvOutput_ =
EWOMS_GET_PARAM(TypeTag, bool, ForceDisableResvFluidInPlaceOutput);
}
/*!
@ -199,25 +193,32 @@ public:
* write to ECL output files
*/
void
allocBuffers(unsigned bufferSize, unsigned reportStepNum, const bool substep, const bool log, const bool isRestart)
allocBuffers(const unsigned bufferSize,
const unsigned reportStepNum,
const bool substep,
const bool log,
const bool isRestart)
{
if (!std::is_same<Discretization, EcfvDiscretization<TypeTag>>::value)
if (! std::is_same<Discretization, EcfvDiscretization<TypeTag>>::value) {
return;
}
const auto& problem = this->simulator_.problem();
this->doAllocBuffers(bufferSize,
reportStepNum,
substep,
log,
isRestart,
simulator_.problem().vapparsActive(std::max(simulator_.episodeIndex(), 0)),
simulator_.problem().materialLawManager()->enableHysteresis(),
simulator_.problem().tracerModel().numTracers(),
simulator_.problem().eclWriter()->getOutputNnc().size());
problem.vapparsActive(std::max(simulator_.episodeIndex(), 0)),
problem.materialLawManager()->enableHysteresis(),
problem.tracerModel().numTracers(),
problem.eclWriter()->getOutputNnc().size());
}
/*!
* \brief Modify the internal buffers according to the intensive quanties relevant
* for an element
* \brief Modify the internal buffers according to the intensive
* quanties relevant for an element
*/
void processElement(const ElementContext& elemCtx)
{
@ -550,7 +551,6 @@ public:
= FluidSystem::viscosity(fsInitial, gasPhaseIdx, intQuants.pvtRegionIndex());
}
// Adding Well RFT data
const auto cartesianIdx = elemCtx.simulator().vanguard().cartesianIndex(globalDofIdx);
if (this->oilConnectionPressures_.count(cartesianIdx) > 0) {
@ -562,25 +562,17 @@ public:
if (this->gasConnectionSaturations_.count(cartesianIdx) > 0) {
this->gasConnectionSaturations_[cartesianIdx] = getValue(fs.saturation(gasPhaseIdx));
}
if (this->wbpData_.count(cartesianIdx) > 0) {
if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
this->wbpData_[cartesianIdx] = getValue(fs.pressure(oilPhaseIdx));
} else if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
this->wbpData_[cartesianIdx] = getValue(fs.pressure(gasPhaseIdx));
} else if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
this->wbpData_[cartesianIdx] = getValue(fs.pressure(waterPhaseIdx));
}
}
// tracers
const auto& tracerModel = simulator_.problem().tracerModel();
if (!this->tracerConcentrations_.empty()) {
for (int tracerIdx = 0; tracerIdx < tracerModel.numTracers(); tracerIdx++) {
if (this->tracerConcentrations_[tracerIdx].empty())
if (! this->tracerConcentrations_.empty()) {
for (int tracerIdx = 0; tracerIdx < tracerModel.numTracers(); ++tracerIdx) {
if (this->tracerConcentrations_[tracerIdx].empty()) {
continue;
}
this->tracerConcentrations_[tracerIdx][globalDofIdx]
= tracerModel.tracerConcentration(tracerIdx, globalDofIdx);
this->tracerConcentrations_[tracerIdx][globalDofIdx] =
tracerModel.tracerConcentration(tracerIdx, globalDofIdx);
}
}
}

View File

@ -148,22 +148,28 @@ public:
simulator.vanguard().eclState(),
simulator.vanguard().summaryConfig(),
simulator.vanguard().grid(),
simulator.vanguard().grid().comm().rank() == 0 ? &simulator.vanguard().equilGrid() : nullptr,
((simulator.vanguard().grid().comm().rank() == 0)
? &simulator.vanguard().equilGrid()
: nullptr),
simulator.vanguard().gridView(),
simulator.vanguard().cartesianIndexMapper(),
simulator.vanguard().grid().comm().rank() == 0 ? &simulator.vanguard().equilCartesianIndexMapper() : nullptr,
EWOMS_GET_PARAM(TypeTag, bool, EnableAsyncEclOutput), EWOMS_GET_PARAM(TypeTag, bool, EnableEsmry))
((simulator.vanguard().grid().comm().rank() == 0)
? &simulator.vanguard().equilCartesianIndexMapper()
: nullptr),
EWOMS_GET_PARAM(TypeTag, bool, EnableAsyncEclOutput),
EWOMS_GET_PARAM(TypeTag, bool, EnableEsmry))
, simulator_(simulator)
{
#ifdef HAVE_DAMARIS
this->damarisUpdate_ = enableDamarisOutput_();
#endif
this->eclOutputModule_ = std::make_unique<EclOutputBlackOilModule<TypeTag>>(simulator, this->wbp_index_list_, this->collectToIORank_);
this->wbp_index_list_.clear();
this->eclOutputModule_ = std::make_unique<EclOutputBlackOilModule<TypeTag>>
(simulator, this->collectToIORank_);
}
~EclWriter()
{ }
{}
const EquilGrid& globalGrid() const
{
@ -204,6 +210,7 @@ public:
simulator_.vanguard().setupTime();
const auto localWellData = simulator_.problem().wellModel().wellData();
const auto localWBP = data::WellBlockAveragePressures{};
const auto localGroupAndNetworkData = simulator_.problem().wellModel()
.groupAndNetworkData(reportStepNum);
@ -220,8 +227,8 @@ public:
this->collectToIORank_.collect({},
eclOutputModule_->getBlockData(),
eclOutputModule_->getWBPData(),
localWellData,
localWBP,
localGroupAndNetworkData,
localAquiferData,
localWellTestState,
@ -283,26 +290,34 @@ public:
if (this->simulation_report_.success.total_newton_iterations != 0) {
miscSummaryData["MSUMNEWT"] = this->simulation_report_.success.total_newton_iterations;
}
{
OPM_TIMEBLOCK(evalSummary);
this->evalSummary(reportStepNum, curTime,
this->collectToIORank_.isParallel() ?
this->collectToIORank_.globalWBPData() :
this->eclOutputModule_->getWBPData(),
localWellData,
localGroupAndNetworkData,
localAquiferData,
this->collectToIORank_.isParallel() ?
this->collectToIORank_.globalBlockData() :
this->eclOutputModule_->getBlockData(),
miscSummaryData, regionData,
inplace,
eclOutputModule_->initialInplace(),
this->collectToIORank_.isParallel()
? this->collectToIORank_.globalInterRegFlows()
: this->eclOutputModule_->getInterRegFlows(),
summaryState(), udqState());
OPM_TIMEBLOCK(evalSummary);
const auto& blockData = this->collectToIORank_.isParallel()
? this->collectToIORank_.globalBlockData()
: this->eclOutputModule_->getBlockData();
const auto& interRegFlows = this->collectToIORank_.isParallel()
? this->collectToIORank_.globalInterRegFlows()
: this->eclOutputModule_->getInterRegFlows();
this->evalSummary(reportStepNum,
curTime,
localWellData,
localWBP,
localGroupAndNetworkData,
localAquiferData,
blockData,
miscSummaryData,
regionData,
inplace,
this->eclOutputModule_->initialInplace(),
interRegFlows,
this->summaryState(),
this->udqState());
}
{
OPM_TIMEBLOCK(outputXXX);
eclOutputModule_->outputProdLog(reportStepNum, isSubStep, forceDisableProdOutput);
@ -314,9 +329,11 @@ public:
void writeOutput(bool isSubStep)
{
OPM_TIMEBLOCK(writeOutput);
const int reportStepNum = simulator_.episodeIndex() + 1;
this->prepareLocalCellData(isSubStep, reportStepNum);
this->eclOutputModule_->outputErrorLog(simulator_.gridView().comm());
#ifdef HAVE_DAMARIS
if (EWOMS_GET_PARAM(TypeTag, bool, EnableDamarisOutput)) {
// N.B. damarisUpdate_ should be set to true if at any time the model geometry changes
@ -340,35 +357,44 @@ public:
}
}
#endif
// output using eclWriter if enabled
auto localWellData = simulator_.problem().wellModel().wellData();
auto localGroupAndNetworkData = simulator_.problem().wellModel()
.groupAndNetworkData(reportStepNum);
auto localAquiferData = simulator_.problem().aquiferModel().aquiferData();
auto localWellTestState = simulator_.problem().wellModel().wellTestState();
auto flowsn = this->eclOutputModule_->getFlowsn();
const bool isFlowsn = this->eclOutputModule_->hasFlowsn();
auto floresn = this->eclOutputModule_->getFloresn();
auto flowsn = this->eclOutputModule_->getFlowsn();
const bool isFloresn = this->eclOutputModule_->hasFloresn();
auto floresn = this->eclOutputModule_->getFloresn();
data::Solution localCellData = {};
if (! isSubStep) {
this->eclOutputModule_->assignToSolution(localCellData);
// add cell data to perforations for Rft output
// Add cell data to perforations for RFT output
this->eclOutputModule_->addRftDataToWells(localWellData, reportStepNum);
}
if (this->collectToIORank_.isParallel()|| this->collectToIORank_.doesNeedReordering()) {
if (this->collectToIORank_.isParallel() ||
this->collectToIORank_.doesNeedReordering())
{
// Note: We don't need WBP (well-block averaged pressures) or
// inter-region flow rate values in order to create restart file
// output. There's consequently no need to collect those
// properties on the I/O rank.
this->collectToIORank_.collect(localCellData,
eclOutputModule_->getBlockData(),
eclOutputModule_->getWBPData(),
this->eclOutputModule_->getBlockData(),
localWellData,
/* wbpData = */ {},
localGroupAndNetworkData,
localAquiferData,
localWellTestState,
{},
/* interRegFlows = */ {},
flowsn,
floresn);
}
@ -376,6 +402,7 @@ public:
if (this->collectToIORank_.isIORank()) {
const Scalar curTime = simulator_.time() + simulator_.timeStepSize();
const Scalar nextStepSize = simulator_.problem().nextTimeStepSize();
this->doWriteOutput(reportStepNum, isSubStep,
std::move(localCellData),
std::move(localWellData),
@ -385,13 +412,11 @@ public:
this->actionState(),
this->udqState(),
this->summaryState(),
simulator_.problem().thresholdPressure().getRestartVector(),
this->simulator_.problem().thresholdPressure().getRestartVector(),
curTime, nextStepSize,
EWOMS_GET_PARAM(TypeTag, bool, EclOutputDoublePrecision),
isFlowsn,
std::move(flowsn),
isFloresn,
std::move(floresn));
isFlowsn, std::move(flowsn),
isFloresn, std::move(floresn));
}
}
@ -487,7 +512,7 @@ public:
Scalar restartTimeStepSize() const
{ return restartTimeStepSize_; }
template<class Serializer>
template <class Serializer>
void serializeOp(Serializer& serializer)
{
serializer(*eclOutputModule_);
@ -496,10 +521,12 @@ public:
private:
static bool enableEclOutput_()
{ return EWOMS_GET_PARAM(TypeTag, bool, EnableEclOutput); }
#ifdef HAVE_DAMARIS
static bool enableDamarisOutput_()
{ return EWOMS_GET_PARAM(TypeTag, bool, EnableDamarisOutput); }
#endif
const EclipseState& eclState() const
{ return simulator_.vanguard().eclState(); }
@ -519,7 +546,8 @@ private:
const int reportStepNum)
{
OPM_TIMEBLOCK(prepareLocalCellData);
if (eclOutputModule_->localDataValid()) {
if (this->eclOutputModule_->localDataValid()) {
return;
}
@ -527,54 +555,68 @@ private:
const int numElements = gridView.size(/*codim=*/0);
const bool log = this->collectToIORank_.isIORank();
eclOutputModule_->allocBuffers(numElements, reportStepNum,
isSubStep, log, /*isRestart*/ false);
this->eclOutputModule_->
allocBuffers(numElements, reportStepNum,
isSubStep, log, /*isRestart*/ false);
ElementContext elemCtx(simulator_);
OPM_BEGIN_PARALLEL_TRY_CATCH();
{
OPM_TIMEBLOCK(prepareCellBasedData);
for (const auto& elem : elements(gridView)) {
elemCtx.updatePrimaryStencil(elem);
elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
eclOutputModule_->processElement(elemCtx);
OPM_BEGIN_PARALLEL_TRY_CATCH();
{
OPM_TIMEBLOCK(prepareCellBasedData);
for (const auto& elem : elements(gridView)) {
elemCtx.updatePrimaryStencil(elem);
elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
this->eclOutputModule_->processElement(elemCtx);
}
}
}
if(!simulator_.model().linearizer().getFlowsInfo().empty()){
if (! this->simulator_.model().linearizer().getFlowsInfo().empty()) {
OPM_TIMEBLOCK(prepareFlowsData);
for (const auto& elem : elements(gridView)) {
elemCtx.updatePrimaryStencil(elem);
elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
eclOutputModule_->processElementFlows(elemCtx);
this->eclOutputModule_->processElementFlows(elemCtx);
}
}
{
OPM_TIMEBLOCK(prepareBlockData);
for (const auto& elem : elements(gridView)) {
elemCtx.updatePrimaryStencil(elem);
elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
eclOutputModule_->processElementBlockData(elemCtx);
}
OPM_TIMEBLOCK(prepareBlockData);
for (const auto& elem : elements(gridView)) {
elemCtx.updatePrimaryStencil(elem);
elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
this->eclOutputModule_->processElementBlockData(elemCtx);
}
}
{
OPM_TIMEBLOCK(prepareFluidInPlace);
OPM_TIMEBLOCK(prepareFluidInPlace);
#ifdef _OPENMP
#pragma omp parallel for
#endif
for (int dofIdx=0; dofIdx < numElements; ++dofIdx){
const auto& intQuants = *(simulator_.model().cachedIntensiveQuantities(dofIdx, /*timeIdx=*/0));
for (int dofIdx = 0; dofIdx < numElements; ++dofIdx) {
const auto& intQuants = *simulator_.model().cachedIntensiveQuantities(dofIdx, /*timeIdx=*/0);
const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
eclOutputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
this->eclOutputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
}
}
}
eclOutputModule_->validateLocalData();
OPM_END_PARALLEL_TRY_CATCH("EclWriter::prepareLocalCellData() failed: ", simulator_.vanguard().grid().comm());
this->eclOutputModule_->validateLocalData();
OPM_END_PARALLEL_TRY_CATCH("EclWriter::prepareLocalCellData() failed: ",
this->simulator_.vanguard().grid().comm());
}
void captureLocalFluxData()
{
OPM_TIMEBLOCK(captureLocalData);
const auto& gridView = this->simulator_.vanguard().gridView();
const auto timeIdx = 0u;
@ -612,6 +654,7 @@ private:
Simulator& simulator_;
std::unique_ptr<EclOutputBlackOilModule<TypeTag>> eclOutputModule_;
Scalar restartTimeStepSize_;
#ifdef HAVE_DAMARIS
bool damarisUpdate_ = false; ///< Whenever this is true writeOutput() will set up Damaris offsets of model fields
#endif

View File

@ -211,6 +211,8 @@ TEST_FOR_TYPE_NAMED(data::SegmentPressures, SegmentPressures)
TEST_FOR_TYPE_NAMED(data::Solution, Solution)
TEST_FOR_TYPE_NAMED(data::Well, dataWell)
TEST_FOR_TYPE_NAMED(data::Wells, Wells)
TEST_FOR_TYPE_NAMED(data::WellBlockAvgPress, dataWBPObject)
TEST_FOR_TYPE_NAMED(data::WellBlockAveragePressures, dataWBPCollection)
TEST_FOR_TYPE(Deck)
TEST_FOR_TYPE(DeckItem)
TEST_FOR_TYPE(DeckKeyword)