19 #ifndef OPM_PARALLELDEBUGOUTPUT_HEADER_INCLUDED
20 #define OPM_PARALLELDEBUGOUTPUT_HEADER_INCLUDED
22 #include <opm/core/grid.h>
23 #include <opm/core/simulator/SimulatorState.hpp>
24 #include <opm/core/simulator/WellState.hpp>
25 #include <opm/core/wells/WellsManager.hpp>
30 #if HAVE_DUNE_CORNERPOINT
31 #include <dune/grid/common/p2pcommunicator.hh>
45 virtual bool collectToIORank(
const SimulatorState& localReservoirState,
46 const WellState& localWellState ) = 0;
56 template <
class Gr
idImpl>
67 Opm::EclipseStateConstPtr ,
74 const WellState& localWellState )
76 globalState_ = &localReservoirState;
77 wellState_ = &localWellState;
83 virtual bool isIORank ()
const {
return true; }
85 virtual int numCells()
const {
return grid_.number_of_cells; }
86 virtual const int*
globalCell()
const {
return grid_.global_cell; }
89 #if HAVE_DUNE_CORNERPOINT
91 class ParallelDebugOutput< Dune::CpGrid> :
public ParallelDebugOutputInterface
94 typedef Dune::CpGrid Grid;
95 typedef typename Grid :: CollectiveCommunication CollectiveCommunication;
104 GlobalCellIndex() : globalId_(-1), localIndex_(-1), isInterior_(true) {}
105 void setGhost() { isInterior_ =
false; }
107 void setId(
const int globalId ) { globalId_ = globalId; }
108 void setIndex(
const int localIndex ) { localIndex_ = localIndex; }
110 int localIndex ()
const {
return localIndex_; }
111 int id ()
const {
return globalId_; }
112 bool isInterior()
const {
return isInterior_; }
115 typedef typename Dune::PersistentContainer< Grid, GlobalCellIndex > GlobalIndexContainer;
117 static const int dimension = Grid :: dimension ;
119 typedef typename Grid :: LeafGridView GridView;
120 typedef GridView AllGridView;
122 typedef Dune :: Point2PointCommunicator< Dune :: SimpleMessageBuffer > P2PCommunicatorType;
123 typedef typename P2PCommunicatorType :: MessageBufferType MessageBufferType;
125 typedef std::vector< GlobalCellIndex > LocalIndexMapType;
127 typedef std::vector<int> IndexMapType;
128 typedef std::vector< IndexMapType > IndexMapStorageType;
130 class DistributeIndexMapping :
public P2PCommunicatorType::DataHandleInterface
133 const std::vector<int>& distributedGlobalIndex_;
134 IndexMapType& localIndexMap_;
135 IndexMapStorageType& indexMaps_;
136 std::map< const int, const int > globalPosition_;
138 std::set< int > checkPosition_;
142 DistributeIndexMapping(
const std::vector<int>& globalIndex,
143 const std::vector<int>& distributedGlobalIndex,
144 IndexMapType& localIndexMap,
145 IndexMapStorageType& indexMaps )
146 : distributedGlobalIndex_( distributedGlobalIndex ),
147 localIndexMap_( localIndexMap ),
148 indexMaps_( indexMaps ),
151 const size_t size = globalIndex.size();
153 for (
size_t index = 0; index < size; ++index )
155 globalPosition_.insert( std::make_pair( globalIndex[ index ], index ) );
159 if( ! indexMaps_.empty() )
162 IndexMapType& indexMap = indexMaps_.back();
163 const size_t localSize = localIndexMap_.size();
164 indexMap.resize( localSize );
165 for(
size_t i=0; i<localSize; ++i )
167 const int id = distributedGlobalIndex_[ localIndexMap_[ i ] ];
168 indexMap[ i ] = globalPosition_[ id ] ;
170 assert( checkPosition_.find(
id ) == checkPosition_.end() );
171 checkPosition_.insert(
id );
177 void pack(
const int link, MessageBufferType& buffer )
181 OPM_THROW(std::logic_error,
"link in method pack is not 0 as execpted");
185 const int size = localIndexMap_.size();
186 buffer.write( size );
188 for(
int index = 0; index < size; ++index )
190 const int globalIdx = distributedGlobalIndex_[ localIndexMap_[ index ] ];
191 buffer.write( globalIdx );
195 void unpack(
const int link, MessageBufferType& buffer )
198 IndexMapType& indexMap = indexMaps_[ link ];
199 assert( ! globalPosition_.empty() );
203 buffer.read( numCells );
204 indexMap.resize( numCells );
205 for(
int index = 0; index <
numCells; ++index )
208 buffer.read( globalId );
209 assert( globalPosition_.find( globalId ) != globalPosition_.end() );
210 indexMap[ index ] = globalPosition_[ globalId ];
212 assert( checkPosition_.find( globalId ) == checkPosition_.end() );
213 checkPosition_.insert( globalId );
222 Opm::EclipseStateConstPtr eclipseState,
224 const double* permeability )
225 : toIORankComm_( otherGrid.comm() ),
226 isIORank_( otherGrid.comm().rank() == ioRank )
228 const CollectiveCommunication& comm = otherGrid.comm();
229 if( comm.size() > 1 )
231 std::set< int > send, recv;
235 Dune::CpGrid globalGrid( otherGrid );
236 globalGrid.switchToGlobalView();
239 globalReservoirState_.init( globalGrid.numCells(), globalGrid.numFaces(), numPhases );
243 WellsManager wells_manager(eclipseState,
245 Opm::UgGridHelpers::numCells( globalGrid ),
246 Opm::UgGridHelpers::globalCell( globalGrid ),
247 Opm::UgGridHelpers::cartDims( globalGrid ),
248 Opm::UgGridHelpers::dimensions( globalGrid ),
249 Opm::UgGridHelpers::cell2Faces( globalGrid ),
250 Opm::UgGridHelpers::beginFaceCentroids( globalGrid ),
254 const Wells* wells = wells_manager.c_wells();
255 globalWellState_.init(wells, globalReservoirState_, globalWellState_ );
258 globalIndex_ = globalGrid.globalCell();
260 unsigned int count = 0;
261 auto gridView = globalGrid.leafGridView();
262 for(
auto it = gridView.begin< 0 >(),
263 end = gridView.end< 0 >(); it != end; ++it, ++count )
266 assert( count == globalIndex_.size() );
268 for(
int i=0; i<comm.size(); ++i)
278 send.insert( ioRank );
281 localIndexMap_.clear();
282 localIndexMap_.reserve( otherGrid.size( 0 ) );
284 unsigned int index = 0;
285 auto localView = otherGrid.leafGridView();
286 for(
auto it = localView.begin< 0 >(),
287 end = localView.end< 0 >(); it != end; ++it, ++index )
289 const auto element = *it ;
291 if( element.partitionType() == Dune :: InteriorEntity )
293 localIndexMap_.push_back( index );
298 toIORankComm_.insertRequest( send, recv );
304 indexMaps_.resize( comm.size() );
308 DistributeIndexMapping distIndexMapping( globalIndex_, otherGrid.globalCell(), localIndexMap_, indexMaps_ );
309 toIORankComm_.exchange( distIndexMapping );
314 globalIndex_ = otherGrid.globalCell();
318 class PackUnPackSimulatorState :
public P2PCommunicatorType::DataHandleInterface
320 const SimulatorState& localState_;
322 const WellState& localWellState_;
323 WellState& globalWellState_;
324 const IndexMapType& localIndexMap_;
325 const IndexMapStorageType& indexMaps_;
328 PackUnPackSimulatorState(
const SimulatorState& localState,
329 SimulatorState& globalState,
330 const WellState& localWellState,
332 const IndexMapType& localIndexMap,
333 const IndexMapStorageType& indexMaps,
335 : localState_( localState ),
336 globalState_( globalState ),
337 localWellState_( localWellState ),
338 globalWellState_( globalWellState ),
339 localIndexMap_( localIndexMap ),
340 indexMaps_( indexMaps )
345 for(
size_t i=globalState_.cellData().size();
346 i<localState.cellData().size(); ++i )
348 const size_t components = localState.cellData()[ i ].size() / localState.numCells();
349 assert( components * localState.numCells() == localState.cellData()[ i ].size() );
350 globalState_.registerCellData( localState.cellDataNames()[ i ], components );
353 MessageBufferType buffer;
356 doUnpack( indexMaps.back(), buffer );
361 void pack(
const int link, MessageBufferType& buffer )
365 OPM_THROW(std::logic_error,
"link in method pack is not 0 as execpted");
369 const size_t numCells = localState_.numCells();
370 const size_t numCellData = localState_.cellData().size();
371 for(
size_t d=0; d<numCellData; ++d )
373 const std::vector< double >& data = localState_.cellData()[ d ];
374 const size_t stride = data.size() /
numCells ;
375 assert( numCells * stride == data.size() );
377 for(
size_t i=0; i<stride; ++i )
380 write( buffer, localIndexMap_, data, i, stride );
385 writeWells( buffer );
388 void doUnpack(
const IndexMapType& indexMap, MessageBufferType& buffer )
391 const size_t numCells = globalState_.numCells();
392 const size_t numCellData = globalState_.cellData().size();
393 for(
size_t d=0; d<numCellData; ++d )
395 std::vector< double >& data = globalState_.cellData()[ d ];
396 const size_t stride = data.size() /
numCells ;
397 assert( numCells * stride == data.size() );
399 for(
size_t i=0; i<stride; ++i )
402 read( buffer, indexMap, data, i, stride );
411 void unpack(
const int link, MessageBufferType& buffer )
413 doUnpack( indexMaps_[ link ], buffer );
417 template <
class Vector>
418 void write( MessageBufferType& buffer,
const IndexMapType& localIndexMap,
419 const Vector& vector,
420 const unsigned int offset = 0,
const unsigned int stride = 1 )
const
422 unsigned int size = localIndexMap.size();
423 buffer.write( size );
424 assert( vector.size() >= stride * size );
425 for(
unsigned int i=0; i<size; ++i )
427 const unsigned int index = localIndexMap[ i ] * stride + offset;
428 assert( index < vector.size() );
429 buffer.write( vector[ index ] );
433 template <
class Vector>
434 void read( MessageBufferType& buffer,
435 const IndexMapType& indexMap,
437 const unsigned int offset = 0,
const unsigned int stride = 1 )
const
439 unsigned int size = 0;
441 assert( size == indexMap.size() );
442 for(
unsigned int i=0; i<size; ++i )
444 const unsigned int index = indexMap[ i ] * stride + offset;
445 assert( index < vector.size() );
446 buffer.read( vector[ index ] );
450 void writeString( MessageBufferType& buffer,
const std::string& s)
const
452 const int size = s.size();
453 buffer.write( size );
454 for(
int i=0; i<size; ++i )
456 buffer.write( s[ i ] );
460 void readString( MessageBufferType& buffer, std::string& s)
const
465 for(
int i=0; i<size; ++i )
467 buffer.read( s[ i ] );
471 void writeWells( MessageBufferType& buffer )
const
473 int nWells = localWellState_.wellMap().size();
474 buffer.write( nWells );
475 auto end = localWellState_.wellMap().end();
476 for(
auto it = localWellState_.wellMap().begin(); it != end; ++it )
478 const std::string& name = it->first;
479 const int wellIdx = it->second[ 0 ];
482 writeString( buffer, name );
485 buffer.write( localWellState_.bhp()[ wellIdx ] );
486 buffer.write( localWellState_.thp()[ wellIdx ] );
487 const int wellRateIdx = wellIdx * localWellState_.numPhases();
488 for(
int np=0; np<localWellState_.numPhases(); ++np )
489 buffer.write( localWellState_.wellRates()[ wellRateIdx + np ] );
496 void readWells( MessageBufferType& buffer )
499 buffer.read( nWells );
502 for(
int well = 0; well < nWells ; ++well )
505 readString( buffer, name );
508 auto it = globalWellState_.wellMap().find( name );
509 if( it == globalWellState_.wellMap().end() )
511 OPM_THROW(std::logic_error,
"global state does not contain well " << name );
513 const int wellIdx = it->second[ 0 ];
515 buffer.read( globalWellState_.bhp()[ wellIdx ] );
516 buffer.read( globalWellState_.thp()[ wellIdx ] );
517 const int wellRateIdx = wellIdx * globalWellState_.numPhases();
518 for(
int np=0; np<globalWellState_.numPhases(); ++np )
519 buffer.read( globalWellState_.wellRates()[ wellRateIdx + np ] );
529 const WellState& localWellState )
531 PackUnPackSimulatorState packUnpack( localReservoirState, globalReservoirState_,
532 localWellState, globalWellState_,
533 localIndexMap_, indexMaps_,
isIORank() );
536 toIORankComm_.exchange( packUnpack );
539 toIORankComm_.barrier();
554 return toIORankComm_.size() > 1;
557 int numCells ()
const {
return globalIndex_.size(); }
560 assert( ! globalIndex_.empty() );
561 return globalIndex_.data();
565 P2PCommunicatorType toIORankComm_;
566 IndexMapType globalIndex_;
567 IndexMapType localIndexMap_;
568 IndexMapStorageType indexMaps_;
569 SimulatorState globalReservoirState_;
571 WellStateFullyImplicitBlackoil globalWellState_;
573 const bool isIORank_;
575 #endif // #if HAVE_DUNE_CORNERPOINT
virtual bool isIORank() const
Definition: ParallelDebugOutput.hpp:83
virtual bool isParallel() const =0
Definition: AdditionalObjectDeleter.hpp:22
virtual const WellState & globalWellState() const =0
Definition: ParallelDebugOutput.hpp:37
virtual const int * globalCell() const =0
virtual const SimulatorState & globalReservoirState() const =0
virtual bool isIORank() const =0
virtual const SimulatorState & globalReservoirState() const
Definition: ParallelDebugOutput.hpp:81
virtual int numCells() const
Definition: ParallelDebugOutput.hpp:85
virtual int numCells() const =0
const WellState * wellState_
Definition: ParallelDebugOutput.hpp:63
virtual ~ParallelDebugOutputInterface()
Definition: ParallelDebugOutput.hpp:42
virtual bool collectToIORank(const SimulatorState &localReservoirState, const WellState &localWellState)
Definition: ParallelDebugOutput.hpp:73
Definition: ParallelDebugOutput.hpp:57
virtual bool isParallel() const
Definition: ParallelDebugOutput.hpp:84
const GridImpl & grid_
Definition: ParallelDebugOutput.hpp:60
virtual const WellState & globalWellState() const
Definition: ParallelDebugOutput.hpp:82
virtual bool collectToIORank(const SimulatorState &localReservoirState, const WellState &localWellState)=0
const SimulatorState * globalState_
Definition: ParallelDebugOutput.hpp:62
ParallelDebugOutput(const GridImpl &grid, Opm::EclipseStateConstPtr, const int, const double *)
Definition: ParallelDebugOutput.hpp:66
ParallelDebugOutputInterface()
Definition: ParallelDebugOutput.hpp:40
virtual const int * globalCell() const
Definition: ParallelDebugOutput.hpp:86