19 #ifndef OPM_PARALLELDEBUGOUTPUT_HEADER_INCLUDED 20 #define OPM_PARALLELDEBUGOUTPUT_HEADER_INCLUDED 22 #include <unordered_set> 24 #include <opm/common/data/SimulationDataContainer.hpp> 27 #include <opm/core/grid.h> 28 #include <opm/core/simulator/WellState.hpp> 29 #include <opm/core/wells/WellsManager.hpp> 31 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp> 32 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp> 33 #include <opm/core/wells/DynamicListEconLimited.hpp> 36 #include <dune/grid/common/p2pcommunicator.hh> 57 virtual bool collectToIORank(
const SimulationDataContainer& localReservoirState,
59 const data::Solution& localCellData,
60 const int wellStateStepNumber ) = 0;
62 virtual const SimulationDataContainer& globalReservoirState()
const = 0 ;
63 virtual const data::Solution& globalCellData()
const = 0 ;
65 virtual bool isIORank()
const = 0;
66 virtual bool isParallel()
const = 0;
67 virtual int numCells()
const = 0 ;
68 virtual const int* globalCell()
const = 0;
71 template <
class Gr
idImpl>
75 const GridImpl& grid_;
77 const SimulationDataContainer* globalState_;
79 const data::Solution* globalCellData_;
85 const Opm::PhaseUsage& )
91 const data::Solution& localCellData,
94 globalState_ = &localReservoirState;
95 wellState_ = &localWellState;
96 globalCellData_ = &localCellData;
100 virtual const SimulationDataContainer& globalReservoirState()
const {
return *globalState_; }
101 virtual const data::Solution& globalCellData()
const 103 return *globalCellData_;
105 virtual const WellStateFullyImplicitBlackoil& globalWellState()
const {
return *wellState_; }
106 virtual bool isIORank ()
const {
return true; }
107 virtual bool isParallel ()
const {
return false; }
108 virtual int numCells()
const {
return Opm::AutoDiffGrid::numCells(grid_); }
109 virtual const int* globalCell()
const {
return Opm::AutoDiffGrid::globalCell(grid_); }
114 class ParallelDebugOutput<
Dune::CpGrid> :
public ParallelDebugOutputInterface
117 typedef Dune::CpGrid Grid;
118 typedef typename Grid :: CollectiveCommunication CollectiveCommunication;
121 class GlobalCellIndex
127 GlobalCellIndex() : globalId_(-1), localIndex_(-1), isInterior_(true) {}
128 void setGhost() { isInterior_ =
false; }
130 void setId(
const int globalId ) { globalId_ = globalId; }
131 void setIndex(
const int localIndex ) { localIndex_ = localIndex; }
133 int localIndex ()
const {
return localIndex_; }
134 int id ()
const {
return globalId_; }
135 bool isInterior()
const {
return isInterior_; }
138 typedef typename Dune::PersistentContainer< Grid, GlobalCellIndex > GlobalIndexContainer;
140 static const int dimension = Grid :: dimension ;
142 typedef typename Grid :: LeafGridView GridView;
143 typedef GridView AllGridView;
145 typedef Dune :: Point2PointCommunicator< Dune :: SimpleMessageBuffer > P2PCommunicatorType;
146 typedef typename P2PCommunicatorType :: MessageBufferType MessageBufferType;
148 typedef std::vector< GlobalCellIndex > LocalIndexMapType;
150 typedef std::vector<int> IndexMapType;
151 typedef std::vector< IndexMapType > IndexMapStorageType;
153 class DistributeIndexMapping :
public P2PCommunicatorType::DataHandleInterface
156 const std::vector<int>& distributedGlobalIndex_;
157 IndexMapType& localIndexMap_;
158 IndexMapStorageType& indexMaps_;
159 std::map< const int, const int > globalPosition_;
161 std::set< int > checkPosition_;
165 DistributeIndexMapping(
const std::vector<int>& globalIndex,
166 const std::vector<int>& distributedGlobalIndex,
167 IndexMapType& localIndexMap,
168 IndexMapStorageType& indexMaps )
169 : distributedGlobalIndex_( distributedGlobalIndex ),
170 localIndexMap_( localIndexMap ),
171 indexMaps_( indexMaps ),
174 const size_t size = globalIndex.size();
176 for (
size_t index = 0; index < size; ++index )
178 globalPosition_.insert( std::make_pair( globalIndex[ index ], index ) );
182 if( ! indexMaps_.empty() )
185 IndexMapType& indexMap = indexMaps_.back();
186 const size_t localSize = localIndexMap_.size();
187 indexMap.resize( localSize );
188 for(
size_t i=0; i<localSize; ++i )
190 const int id = distributedGlobalIndex_[ localIndexMap_[ i ] ];
191 indexMap[ i ] = globalPosition_[ id ] ;
193 assert( checkPosition_.find(
id ) == checkPosition_.end() );
194 checkPosition_.insert(
id );
200 void pack(
const int link, MessageBufferType& buffer )
204 OPM_THROW(std::logic_error,
"link in method pack is not 0 as execpted");
208 const int size = localIndexMap_.size();
209 buffer.write( size );
211 for(
int index = 0; index < size; ++index )
213 const int globalIdx = distributedGlobalIndex_[ localIndexMap_[ index ] ];
214 buffer.write( globalIdx );
218 void unpack(
const int link, MessageBufferType& buffer )
221 IndexMapType& indexMap = indexMaps_[ link ];
222 assert( ! globalPosition_.empty() );
226 buffer.read( numCells );
227 indexMap.resize( numCells );
228 for(
int index = 0; index < numCells; ++index )
231 buffer.read( globalId );
232 assert( globalPosition_.find( globalId ) != globalPosition_.end() );
233 indexMap[ index ] = globalPosition_[ globalId ];
235 assert( checkPosition_.find( globalId ) == checkPosition_.end() );
236 checkPosition_.insert( globalId );
249 ParallelDebugOutput(
const Dune::CpGrid& otherGrid,
250 const EclipseState& eclipseState,
252 const Opm::PhaseUsage& phaseUsage)
254 eclipseState_( eclipseState ),
255 globalCellData_(new data::Solution),
257 phaseUsage_(phaseUsage)
261 Dune::CpGrid distributed_grid = otherGrid;
263 const CollectiveCommunication& comm = otherGrid.comm();
264 if( comm.size() > 1 )
266 std::set< int > send, recv;
267 distributed_grid.switchToDistributedView();
268 toIORankComm_ = distributed_grid.comm();
269 isIORank_ = (distributed_grid.comm().rank() == ioRank);
275 grid_.reset(
new Dune::CpGrid(otherGrid ) );
276 grid_->switchToGlobalView();
277 Dune::CpGrid& globalGrid = *grid_;
280 globalReservoirState_.reset(
new SimulationDataContainer( globalGrid.numCells(), globalGrid.numFaces(), numPhases ));
283 globalIndex_ = globalGrid.globalCell();
285 unsigned int count = 0;
286 auto gridView = globalGrid.leafGridView();
287 for(
auto it = gridView.begin< 0 >(),
288 end = gridView.end< 0 >(); it != end; ++it, ++count )
291 assert( count == globalIndex_.size() );
293 for(
int i=0; i<comm.size(); ++i)
305 globalReservoirState_.reset(
new SimulationDataContainer( 0, 0, 0));
306 send.insert( ioRank );
309 localIndexMap_.clear();
310 localIndexMap_.reserve( distributed_grid.size( 0 ) );
312 unsigned int index = 0;
313 auto localView = distributed_grid.leafGridView();
314 for(
auto it = localView.begin< 0 >(),
315 end = localView.end< 0 >(); it != end; ++it, ++index )
317 const auto element = *it ;
319 if( element.partitionType() == Dune :: InteriorEntity )
321 localIndexMap_.push_back( index );
326 toIORankComm_.insertRequest( send, recv );
332 indexMaps_.resize( comm.size() );
336 DistributeIndexMapping distIndexMapping( globalIndex_, distributed_grid.globalCell(), localIndexMap_, indexMaps_ );
337 toIORankComm_.exchange( distIndexMapping );
342 globalIndex_ = distributed_grid.globalCell();
346 class PackUnPackSimulationDataContainer :
public P2PCommunicatorType::DataHandleInterface
348 const data::Solution& localCellData_;
349 data::Solution& globalCellData_;
350 const WellStateFullyImplicitBlackoil& localWellState_;
351 WellStateFullyImplicitBlackoil& globalWellState_;
352 const IndexMapType& localIndexMap_;
353 const IndexMapStorageType& indexMaps_;
356 PackUnPackSimulationDataContainer( std::size_t numGlobalCells,
357 const data::Solution& localCellData,
358 data::Solution& globalCellData,
359 const WellStateFullyImplicitBlackoil& localWellState,
360 WellStateFullyImplicitBlackoil& globalWellState,
361 const IndexMapType& localIndexMap,
362 const IndexMapStorageType& indexMaps,
363 const bool isIORank )
364 : localCellData_( localCellData ),
365 globalCellData_( globalCellData ),
366 localWellState_( localWellState ),
367 globalWellState_( globalWellState ),
368 localIndexMap_( localIndexMap ),
369 indexMaps_( indexMaps )
375 for (
const auto& pair : localCellData_) {
376 const std::string& key = pair.first;
377 std::size_t container_size = numGlobalCells;
378 auto ret = globalCellData_.insert(key, pair.second.dim,
379 std::vector<double>(container_size),
382 DUNE_UNUSED_PARAMETER(ret.second);
385 MessageBufferType buffer;
388 doUnpack( indexMaps.back(), buffer );
393 void pack(
const int link, MessageBufferType& buffer )
397 OPM_THROW(std::logic_error,
"link in method pack is not 0 as execpted");
401 for (
const auto& pair : localCellData_) {
402 const auto& data = pair.second.data;
405 write( buffer, localIndexMap_, data);
409 writeWells( buffer );
412 void doUnpack(
const IndexMapType& indexMap, MessageBufferType& buffer )
416 for (
auto& pair : localCellData_) {
417 const std::string& key = pair.first;
418 auto& data = globalCellData_.data(key);
421 read( buffer, indexMap, data);
429 void unpack(
const int link, MessageBufferType& buffer )
431 doUnpack( indexMaps_[ link ], buffer );
435 template <
class Vector>
436 void write( MessageBufferType& buffer,
const IndexMapType& localIndexMap,
437 const Vector& vector,
438 const unsigned int offset = 0,
const unsigned int stride = 1 )
const 440 unsigned int size = localIndexMap.size();
441 buffer.write( size );
442 assert( vector.size() >= stride * size );
443 for(
unsigned int i=0; i<size; ++i )
445 const unsigned int index = localIndexMap[ i ] * stride + offset;
446 assert( index < vector.size() );
447 buffer.write( vector[ index ] );
451 template <
class Vector>
452 void read( MessageBufferType& buffer,
453 const IndexMapType& indexMap,
455 const unsigned int offset = 0,
const unsigned int stride = 1 )
const 457 unsigned int size = 0;
459 assert( size == indexMap.size() );
460 for(
unsigned int i=0; i<size; ++i )
462 const unsigned int index = indexMap[ i ] * stride + offset;
463 assert( index < vector.size() );
464 buffer.read( vector[ index ] );
468 void writeString( MessageBufferType& buffer,
const std::string& s)
const 470 const int size = s.size();
471 buffer.write( size );
472 for(
int i=0; i<size; ++i )
474 buffer.write( s[ i ] );
478 void readString( MessageBufferType& buffer, std::string& s)
const 483 for(
int i=0; i<size; ++i )
485 buffer.read( s[ i ] );
489 void writeWells( MessageBufferType& buffer )
const 491 int nWells = localWellState_.wellMap().size();
492 buffer.write( nWells );
493 auto end = localWellState_.wellMap().end();
494 for(
auto it = localWellState_.wellMap().begin(); it != end; ++it )
496 const std::string& name = it->first;
497 const int wellIdx = it->second[ 0 ];
500 writeString( buffer, name );
503 buffer.write( localWellState_.bhp()[ wellIdx ] );
504 buffer.write( localWellState_.thp()[ wellIdx ] );
505 const int wellRateIdx = wellIdx * localWellState_.numPhases();
506 for(
int np=0; np<localWellState_.numPhases(); ++np )
507 buffer.write( localWellState_.wellRates()[ wellRateIdx + np ] );
510 buffer.write(localWellState_.currentControls()[ wellIdx ]);
515 const int end_con = it->second[1] + it->second[2];
517 for(
int con = it->second[1]; con < end_con; ++con )
519 buffer.write( localWellState_.perfRates()[ con ] );
522 for(
int con = it->second[1]; con < end_con; ++con )
524 buffer.write( localWellState_.perfPress()[ con ] );
528 const int np = localWellState_.perfPhaseRates().size() /
529 localWellState_.perfRates().size();
531 for(
int con = it->second[1]*np; con < end_con*np; ++con )
533 buffer.write( localWellState_.perfPhaseRates()[ con ] );
538 void readWells( MessageBufferType& buffer )
541 buffer.read( nWells );
544 for(
int well = 0; well < nWells ; ++well )
547 readString( buffer, name );
550 auto it = globalWellState_.wellMap().find( name );
551 if( it == globalWellState_.wellMap().end() )
553 OPM_THROW(std::logic_error,
"global state does not contain well " << name );
555 const int wellIdx = it->second[ 0 ];
557 buffer.read( globalWellState_.bhp()[ wellIdx ] );
558 buffer.read( globalWellState_.thp()[ wellIdx ] );
559 const int wellRateIdx = wellIdx * globalWellState_.numPhases();
560 for(
int np=0; np<globalWellState_.numPhases(); ++np )
561 buffer.read( globalWellState_.wellRates()[ wellRateIdx + np ] );
564 buffer.read(globalWellState_.currentControls()[ wellIdx ]);
569 const int end_con = it->second[1] + it->second[2];
571 for(
int con = it->second[1]; con < end_con; ++con )
573 buffer.read( globalWellState_.perfRates()[ con ] );
576 for(
int con = it->second[1]; con < end_con; ++con )
578 buffer.read( globalWellState_.perfPress()[ con ] );
582 const int np = globalWellState_.perfPhaseRates().size() /
583 globalWellState_.perfRates().size();
585 for(
int con = it->second[1]*np; con < end_con*np; ++con )
587 buffer.read( globalWellState_.perfPhaseRates()[ con ] );
595 const WellStateFullyImplicitBlackoil& localWellState,
596 const data::Solution& localCellData,
597 const int wellStateStepNumber )
601 Dune::CpGrid& globalGrid = *grid_;
605 const DynamicListEconLimited dynamic_list_econ_limited;
607 WellsManager wells_manager(eclipseState_,
609 Opm::UgGridHelpers::numCells( globalGrid ),
610 Opm::UgGridHelpers::globalCell( globalGrid ),
611 Opm::UgGridHelpers::cartDims( globalGrid ),
612 Opm::UgGridHelpers::dimensions( globalGrid ),
613 Opm::UgGridHelpers::cell2Faces( globalGrid ),
614 Opm::UgGridHelpers::beginFaceCentroids( globalGrid ),
615 dynamic_list_econ_limited,
621 std::unordered_set<std::string>());
623 const Wells* wells = wells_manager.c_wells();
624 globalWellState_.init(wells, *globalReservoirState_, globalWellState_, phaseUsage_ );
625 globalCellData_->clear();
628 PackUnPackSimulationDataContainer packUnpack( numCells(),
629 localCellData, *globalCellData_,
630 localWellState, globalWellState_,
631 localIndexMap_, indexMaps_,
635 toIORankComm_.exchange( packUnpack );
638 toIORankComm_.barrier();
643 const std::map<std::string, std::vector<double> > no_extra_data;
644 solutionToSim(*globalCellData_, no_extra_data, phaseUsage_, *globalReservoirState_);
649 const SimulationDataContainer& globalReservoirState()
const {
return *globalReservoirState_; }
651 const data::Solution& globalCellData()
const 653 return *globalCellData_;
656 const WellStateFullyImplicitBlackoil& globalWellState()
const {
return globalWellState_; }
663 bool isParallel()
const 665 return toIORankComm_.size() > 1;
668 int numCells ()
const {
return globalIndex_.size(); }
669 const int* globalCell ()
const 671 assert( ! globalIndex_.empty() );
672 return globalIndex_.data();
676 std::unique_ptr< Dune::CpGrid > grid_;
677 const EclipseState& eclipseState_;
678 P2PCommunicatorType toIORankComm_;
679 IndexMapType globalIndex_;
680 IndexMapType localIndexMap_;
681 IndexMapStorageType indexMaps_;
682 std::unique_ptr<SimulationDataContainer> globalReservoirState_;
683 std::unique_ptr<data::Solution> globalCellData_;
685 WellStateFullyImplicitBlackoil globalWellState_;
689 Opm::PhaseUsage phaseUsage_;
691 #endif // #if HAVE_OPM_GRID Definition: ParallelDebugOutput.hpp:72
Definition: ISTLSolver.hpp:44
virtual bool collectToIORank(const SimulationDataContainer &localReservoirState, const WellStateFullyImplicitBlackoil &localWellState, const data::Solution &localCellData, const int)
gather solution to rank 0 for EclipseWriter
Definition: ParallelDebugOutput.hpp:89
virtual bool collectToIORank(const SimulationDataContainer &localReservoirState, const WellStateFullyImplicitBlackoil &localWellState, const data::Solution &localCellData, const int wellStateStepNumber)=0
gather solution to rank 0 for EclipseWriter
Definition: ParallelDebugOutput.hpp:42
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
bool isIORank(const boost::any ¶llel_info)
Return true if this is a serial run, or rank zero on an MPI run.
Definition: NewtonIterationUtilities.cpp:294
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
void solutionToSim(const data::Solution &sol, const std::map< std::string, std::vector< double > > &, PhaseUsage phases, SimulationDataContainer &state)
Copies the following fields from sol into state (all conditionally): PRESSURE, TEMP, SWAT, SGAS, RS, RV, SSOL Also handles extra data such as hysteresis parameters, SOMAX, etc.
Definition: Compat.cpp:134