22 #ifndef OPM_REDISTRIBUTEDATAHANDLES_HEADER
23 #define OPM_REDISTRIBUTEDATAHANDLES_HEADER
25 #include <unordered_set>
27 #include <type_traits>
30 #include <opm/core/simulator/BlackoilState.hpp>
32 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
33 #include <opm/autodiff/ExtractParallelGridInformationToISTL.hpp>
34 #include <opm/autodiff/createGlobalCellArray.hpp>
36 #include<boost/any.hpp>
42 inline std::unordered_set<std::string>
43 distributeGridAndData( Grid& ,
47 BlackoilPropsAdFromDeck& ,
49 std::shared_ptr<BlackoilPropsAdFromDeck::MaterialLawManager>&,
54 return std::unordered_set<std::string>();
62 template<
class Iter1,
class Iter2=Iter1>
65 typedef typename std::iterator_traits<Iter1>::value_type DataType2;
68 typedef typename std::iterator_traits<Iter1>::value_type DataType;
74 const Iter2& receive_begin,
76 : send_(send_begin), receive_(receive_begin), size_(size)
78 static_assert(std::is_same<DataType,DataType2>::value,
79 "Iter1 and Iter2 have to have the same value_type!");
82 template<
class Buffer>
83 void gather(Buffer& buffer, std::size_t i)
85 for(
auto index = i*size(i), end = (i+1)*size(i);
88 buffer.write(send_[index]);
92 template<
class Buffer>
93 void scatter(Buffer& buffer, std::size_t i, std::size_t s OPM_OPTIM_UNUSED)
98 for(
auto index = i*size(i), end = (i+1)*size(i);
101 buffer.read(receive_[index]);
110 std::size_t size(std::size_t)
121 #if HAVE_OPM_GRID && HAVE_MPI
122 class ThresholdPressureDataHandle
127 typedef double DataType;
133 ThresholdPressureDataHandle(
const Dune::CpGrid& sendGrid,
134 const Dune::CpGrid& recvGrid,
135 const std::vector<double>& sendPressures,
136 std::vector<double>& recvPressures)
137 : sendGrid_(sendGrid), recvGrid_(recvGrid), sendPressures_(sendPressures),
138 recvPressures_(recvPressures)
141 bool fixedsize(
int ,
int )
146 std::size_t size(
const T& e)
148 if ( T::codimension == 0)
150 return sendGrid_.numCellFaces(e.index());
154 OPM_THROW(std::logic_error,
"Data handle can only be used for elements");
157 template<
class B,
class T>
158 void gather(B& buffer,
const T& e)
160 assert( T::codimension == 0);
161 for (
int i=0; i< sendGrid_.numCellFaces(e.index()); ++i )
163 buffer.write(sendPressures_[sendGrid_.cellFace(e.index(), i)]);
166 template<
class B,
class T>
167 void scatter(B& buffer,
const T& e, std::size_t )
169 assert( T::codimension == 0);
170 for (
int i=0; i< recvGrid_.numCellFaces(e.index()); ++i )
174 recvPressures_[recvGrid_.cellFace(e.index(), i)]=val;
177 bool contains(
int dim,
int codim)
179 return dim==3 && codim==0;
184 const Dune::CpGrid& sendGrid_;
186 const Dune::CpGrid& recvGrid_;
188 const std::vector<double>& sendPressures_;
190 std::vector<double>& recvPressures_;
194 class GeologyDataHandle
198 typedef double DataType;
204 GeologyDataHandle(
const Dune::CpGrid& sendGrid,
205 const Dune::CpGrid& recvGrid,
206 const DerivedGeology& sendGeology,
207 DerivedGeology& recvGeology)
208 : sendGrid_(sendGrid), recvGrid_(recvGrid), sendGeology_(sendGeology),
209 recvGeology_(recvGeology)
212 bool fixedsize(
int ,
int )
217 std::size_t size(
const T& e)
219 if ( T::codimension == 0)
221 return 1 + sendGrid_.numCellFaces(e.index());
225 OPM_THROW(std::logic_error,
"Data handle can only be used for elements");
228 template<
class B,
class T>
229 void gather(B& buffer,
const T& e)
231 assert( T::codimension == 0);
232 buffer.write(sendGeology_.poreVolume()[e.index()]);
233 for (
int i=0; i< sendGrid_.numCellFaces(e.index()); ++i )
235 buffer.write(sendGeology_.transmissibility()[sendGrid_.cellFace(e.index(), i)]);
238 template<
class B,
class T>
239 void scatter(B& buffer,
const T& e, std::size_t )
241 assert( T::codimension == 0);
244 recvGeology_.poreVolume()[e.index()]=val;
245 for (
int i=0; i< recvGrid_.numCellFaces(e.index()); ++i )
248 recvGeology_.transmissibility()[recvGrid_.cellFace(e.index(), i)]=val;
251 bool contains(
int dim,
int codim)
253 return dim==3 && codim==0;
257 const Dune::CpGrid& sendGrid_;
259 const Dune::CpGrid& recvGrid_;
261 const DerivedGeology& sendGeology_;
263 DerivedGeology& recvGeology_;
267 class BlackoilStateDataHandle
271 typedef double DataType;
277 BlackoilStateDataHandle(
const Dune::CpGrid& sendGrid,
278 const Dune::CpGrid& recvGrid,
279 const BlackoilState& sendState,
280 BlackoilState& recvState)
281 : sendGrid_(sendGrid), recvGrid_(recvGrid), sendState_(sendState), recvState_(recvState)
284 recvState.surfacevol().resize(recvGrid.numCells()*sendState.numPhases(),
285 std::numeric_limits<double>::max());
286 recvState.hydroCarbonState().resize(recvGrid.numCells());
289 bool fixedsize(
int ,
int )
295 std::size_t size(
const T& e)
297 if ( T::codimension == 0)
299 return 2 * sendState_.numPhases() + 5 + 2*sendGrid_.numCellFaces(e.index());
303 OPM_THROW(std::logic_error,
"Data handle can only be used for elements");
307 template<
class B,
class T>
308 void gather(B& buffer,
const T& e)
310 assert( T::codimension == 0);
312 for (
size_t i=0; i<sendState_.numPhases(); ++i )
314 buffer.write(sendState_.surfacevol()[e.index()*sendState_.numPhases()+i]);
316 buffer.write(sendState_.gasoilratio()[e.index()]);
317 buffer.write(sendState_.rv()[e.index()]);
318 buffer.write(sendState_.pressure()[e.index()]);
319 buffer.write(sendState_.temperature()[e.index()]);
321 double hydroCarbonState_ = sendState_.hydroCarbonState()[e.index()];
322 buffer.write(hydroCarbonState_);
324 for (
size_t i=0; i<sendState_.numPhases(); ++i )
326 buffer.write(sendState_.saturation()[e.index()*sendState_.numPhases()+i]);
328 for (
int i=0; i<sendGrid_.numCellFaces(e.index()); ++i )
330 buffer.write(sendState_.facepressure()[sendGrid_.cellFace(e.index(), i)]);
332 for (
int i=0; i<sendGrid_.numCellFaces(e.index()); ++i )
334 buffer.write(sendState_.faceflux()[sendGrid_.cellFace(e.index(), i)]);
337 template<
class B,
class T>
338 void scatter(B& buffer,
const T& e, std::size_t size_arg)
340 assert( T::codimension == 0);
341 assert( size_arg == 2 * recvState_.numPhases() + 5 +2*recvGrid_.numCellFaces(e.index()));
342 static_cast<void>(size_arg);
345 for (
size_t i=0; i<recvState_.numPhases(); ++i )
348 recvState_.surfacevol()[e.index()*sendState_.numPhases()+i]=val;
351 recvState_.gasoilratio()[e.index()]=val;
353 recvState_.rv()[e.index()]=val;
355 recvState_.pressure()[e.index()]=val;
357 recvState_.temperature()[e.index()]=val;
360 recvState_.hydroCarbonState()[e.index()]=
static_cast<HydroCarbonState
>(val);
362 for (
size_t i=0; i<recvState_.numPhases(); ++i )
365 recvState_.saturation()[e.index()*sendState_.numPhases()+i]=val;
367 for (
int i=0; i<recvGrid_.numCellFaces(e.index()); ++i )
370 recvState_.facepressure()[recvGrid_.cellFace(e.index(), i)]=val;
372 for (
int i=0; i<recvGrid_.numCellFaces(e.index()); ++i )
375 recvState_.faceflux()[recvGrid_.cellFace(e.index(), i)]=val;
378 bool contains(
int dim,
int codim)
380 return dim==3 && codim==0;
384 const Dune::CpGrid& sendGrid_;
386 const Dune::CpGrid& recvGrid_;
388 const BlackoilState& sendState_;
390 BlackoilState& recvState_;
394 class BlackoilPropsDataHandle
398 typedef double DataType;
402 BlackoilPropsDataHandle(
const BlackoilPropsAdFromDeck& sendProps,
403 BlackoilPropsAdFromDeck& recvProps)
404 : sendProps_(sendProps), recvProps_(recvProps),
408 if ( sendProps.satOilMax_.size()>0 )
412 recvProps_.satOilMax_.resize(recvProps_.cellPvtRegionIdx_.size(),
413 -std::numeric_limits<double>::max());
418 bool fixedsize(
int ,
int )
424 std::size_t size(
const T&)
426 if ( T::codimension == 0)
432 OPM_THROW(std::logic_error,
"Data handle can only be used for elements");
435 template<
class B,
class T>
436 void gather(B& buffer,
const T& e)
438 assert( T::codimension == 0);
440 buffer.write(sendProps_.cellPvtRegionIndex()[e.index()]);
441 for( std::size_t i = 0; i < 9; ++i )
443 buffer.write(sendProps_.rock_.permeability_[e.index()*9+i]);
445 buffer.write(sendProps_.rock_.porosity_[e.index()]);
447 buffer.write(sendProps_.satOilMax_[e.index()]);
450 template<
class B,
class T>
451 void scatter(B& buffer,
const T& e, std::size_t size_arg)
453 assert( T::codimension == 0);
454 assert( size_arg==size_ ); (void) size_arg;
457 recvProps_.cellPvtRegionIdx_[e.index()]=val;
458 for( std::size_t i = 0; i < 9; ++i )
461 recvProps_.rock_.permeability_[e.index()*9+i]
465 recvProps_.rock_.porosity_[e.index()]=val;
468 recvProps_.satOilMax_[e.index()]=val;
471 bool contains(
int dim,
int codim)
473 return dim==3 && codim==0;
477 const BlackoilPropsAdFromDeck& sendProps_;
479 BlackoilPropsAdFromDeck& recvProps_;
488 std::unordered_set<std::string>
489 distributeGridAndData( Dune::CpGrid& grid,
490 const Opm::Deck& deck,
491 const EclipseState& eclipseState,
492 BlackoilState& state,
493 BlackoilPropsAdFromDeck& properties,
494 DerivedGeology& geology,
495 std::shared_ptr<BlackoilPropsAdFromDeck::MaterialLawManager>& material_law_manager,
496 std::vector<double>& threshold_pressures,
497 boost::any& parallelInformation,
498 const bool useLocalPerm)
500 Dune::CpGrid global_grid ( grid );
501 global_grid.switchToGlobalView();
505 auto my_defunct_wells = get<1>(grid.loadBalance(&eclipseState,
506 geology.transmissibility().data()));
507 grid.switchToDistributedView();
508 std::vector<int> compressedToCartesianIdx;
510 typedef BlackoilPropsAdFromDeck::MaterialLawManager MaterialLawManager;
511 auto distributed_material_law_manager = std::make_shared<MaterialLawManager>();
512 distributed_material_law_manager->initFromDeck(deck, eclipseState, compressedToCartesianIdx);
517 typedef Dune::CpGrid::ParallelIndexSet IndexSet;
518 const IndexSet& local_indices = grid.getCellIndexSet();
519 for (
auto index : local_indices )
521 distributed_material_law_manager->materialLawParamsPointerReferenceHack(index.local()) =
522 material_law_manager->materialLawParamsPointerReferenceHack(index.global());
524 distributed_material_law_manager->oilWaterScaledEpsInfoDrainagePointerReferenceHack(index.local()) =
525 material_law_manager->oilWaterScaledEpsInfoDrainagePointerReferenceHack(index.global());
527 BlackoilPropsAdFromDeck distributed_props(properties,
528 distributed_material_law_manager,
530 BlackoilState distributed_state(grid.numCells(), grid.numFaces(), state.numPhases());
531 BlackoilStateDataHandle state_handle(global_grid, grid,
532 state, distributed_state);
533 BlackoilPropsDataHandle props_handle(properties,
535 grid.scatterData(state_handle);
536 grid.scatterData(props_handle);
539 DerivedGeology distributed_geology(grid,
540 distributed_props, eclipseState,
541 useLocalPerm, geology.gravity());
542 GeologyDataHandle geo_handle(global_grid, grid,
543 geology, distributed_geology);
544 grid.scatterData(geo_handle);
546 std::vector<double> distributed_pressures;
548 if( !threshold_pressures.empty() )
550 if( threshold_pressures.size() !=
551 static_cast<std::size_t
>(UgGridHelpers::numFaces(global_grid)) )
553 OPM_THROW(std::runtime_error,
"NNCs not yet supported for parallel runs. "
554 << UgGridHelpers::numFaces(grid) <<
" faces but " <<
555 threshold_pressures.size()<<
" threshold pressure values");
557 distributed_pressures.resize(UgGridHelpers::numFaces(grid));
558 ThresholdPressureDataHandle press_handle(global_grid, grid,
560 distributed_pressures);
561 grid.scatterData(press_handle);
565 properties = distributed_props;
566 geology = distributed_geology;
567 state = distributed_state;
568 material_law_manager = distributed_material_law_manager;
569 threshold_pressures = distributed_pressures;
570 extractParallelGridInformationToISTL(grid, parallelInformation);
572 return my_defunct_wells;
A handle that copies a fixed number data per index.
Definition: RedistributeDataHandles.hpp:63
void createGlobalCellArray(const Grid &grid, std::vector< int > &dest)
Create a mapping from a global cell index of a grid to the logically Cartesian index of the ECL deck...
Definition: createGlobalCellArray.hpp:31
FixedSizeIterCopyHandle(const Iter1 &send_begin, const Iter2 &receive_begin, std::size_t size=1)
Constructor.
Definition: RedistributeDataHandles.hpp:73