48#ifndef OPM_CPGRIDDATA_HEADER
49#define OPM_CPGRIDDATA_HEADER
52#include <dune/common/parallel/mpihelper.hh>
54#include <dune/istl/owneroverlapcopy.hh>
57#include <dune/common/parallel/communication.hh>
58#include <dune/common/parallel/variablesizecommunicator.hh>
59#include <dune/grid/common/gridenums.hh>
62#include <opm/input/eclipse/EclipseState/Grid/NNC.hpp>
67#include "Entity2IndexDataHandle.hpp"
68#include "CpGridDataTraits.hpp"
71#include "Geometry.hpp"
90class LevelGlobalIdSet;
91class PartitionTypeIndicator;
92template<
int,
int>
class Geometry;
93template<
int>
class Entity;
94template<
int>
class EntityRep;
99 const std::vector<std::array<int,3>>&,
100 const std::vector<std::array<int,3>>&);
105 const std::array<int, 3>&,
108void refinePatch_and_check(
const std::array<int,3>&,
109 const std::array<int,3>&,
110 const std::array<int,3>&);
113 const std::vector<std::array<int,3>>&,
114 const std::vector<std::array<int,3>>&,
115 const std::vector<std::array<int,3>>&,
116 const std::vector<std::string>&);
127template<
class T,
int i>
struct Mover;
144 const std::vector<std::array<int,3>>&,
145 const std::vector<std::array<int,3>>&);
152 const std::array<int, 3>&,
155 void ::refinePatch_and_check(
const std::array<int,3>&,
156 const std::array<int,3>&,
157 const std::array<int,3>&);
161 const std::vector<std::array<int,3>>&,
162 const std::vector<std::array<int,3>>&,
163 const std::vector<std::array<int,3>>&,
164 const std::vector<std::string>&);
175#ifndef MAX_DATA_COMMUNICATED_PER_ENTITY
195 explicit CpGridData(MPIHelper::MPICommunicator comm, std::vector<std::shared_ptr<CpGridData>>& data);
200 explicit CpGridData(std::vector<std::shared_ptr<CpGridData>>& data);
208 int size(
int codim)
const;
211 int size (GeometryType type)
const
214 return size(3 - type.dim());
234 void readEclipseFormat(
const std::string& filename,
bool periodic_extension,
bool turn_normals =
false);
245 void processEclipseFormat(
const Opm::Deck& deck,
bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false,
246 const std::vector<double>& poreVolume = std::vector<double>());
260 std::vector<std::size_t>
processEclipseFormat(
const Opm::EclipseGrid* ecl_grid, Opm::EclipseState* ecl_state,
261 bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false,
262 bool pinchActive =
true);
278 Opm::EclipseState* ecl_state,
280 std::array<std::set<std::pair<int, int>>, 2>& nnc,
281 bool remove_ij_boundary,
bool turn_normals,
bool pinchActive,
282 double tolerance_unique_points);
291 void getIJK(
int c, std::array<int,3>& ijk)
const
293 int gc = global_cell_[c];
294 ijk[0] = gc % logical_cartesian_size_[0]; gc /= logical_cartesian_size_[0];
295 ijk[1] = gc % logical_cartesian_size_[1];
296 ijk[2] = gc / logical_cartesian_size_[1];
307 std::array<int,3> getPatchDim(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
316 std::vector<int> getPatchCorners(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
325 std::vector<int> getPatchFaces(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
334 std::vector<int> getPatchCells(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
343 std::vector<int> getPatchBoundaryCorners(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
350 bool disjointPatches(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
360 getPatchesCells(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
376 Geometry<3,3> cellifyPatch(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK,
378 std::array<int,8>& cellifiedPatch_to_point,
379 std::array<int,8>& allcorners_cellifiedPatch)
const;
385 std::array<double,3> getAverageArr(
const std::vector<std::array<double,3>>& vec)
const;
409 std::tuple< const std::shared_ptr<CpGridData>,
410 const std::vector<std::array<int,2>>,
411 const std::vector<std::tuple<int,std::vector<int>>>,
412 const std::tuple<int, std::vector<int>>,
413 const std::vector<std::array<int,2>>,
414 const std::vector<std::array<int,2>>>
415 refineSingleCell(
const std::array<int,3>& cells_per_dim,
const int& parent_idx)
const;
429 std::tuple< std::shared_ptr<CpGridData>,
430 const std::vector<std::array<int,2>>,
431 const std::vector<std::tuple<int,std::vector<int>>>,
432 const std::vector<std::tuple<int,std::vector<int>>>,
433 const std::vector<std::tuple<int,std::vector<int>>>,
434 const std::vector<std::array<int,2>>,
435 const std::vector<std::array<int,2>>>
436 refinePatch(
const std::array<int,3>& cells_per_dim,
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
444 std::array<double,3> computeEclCentroid(
const int idx)
const;
452 std::array<double,3> computeEclCentroid(
const Entity<0>& elem)
const;
455 void computeUniqueBoundaryIds();
462 return use_unique_boundary_ids_;
469 use_unique_boundary_ids_ = uids;
470 if (use_unique_boundary_ids_ && unique_boundary_ids_.empty()) {
471 computeUniqueBoundaryIds();
494 return logical_cartesian_size_;
502 const std::vector<int>& cell_part);
509 template<
class DataHandle>
510 void communicate(DataHandle& data, InterfaceType iftype, CommunicationDirection dir);
516 using CollectiveCommunication = CpGridDataTraits::CollectiveCommunication;
522 using Communicator = CpGridDataTraits::Communicator;
525 using InterfaceMap = CpGridDataTraits::InterfaceMap;
528 using CommunicationType = CpGridDataTraits::CommunicationType;
531 using ParallelIndexSet = CpGridDataTraits::ParallelIndexSet;
534 using RemoteIndices = CpGridDataTraits::RemoteIndices;
539 CommunicationType& cellCommunication()
547 const CommunicationType& cellCommunication()
const
552 ParallelIndexSet& cellIndexSet()
554 return cellCommunication().indexSet();
557 const ParallelIndexSet& cellIndexSet()
const
559 return cellCommunication().indexSet();
562 RemoteIndices& cellRemoteIndices()
564 return cellCommunication().remoteIndices();
567 const RemoteIndices& cellRemoteIndices()
const
569 return cellCommunication().remoteIndices();
576 return aquifer_cells_;
582 void populateGlobalCellIndexSet();
591 template<
class DataHandle>
592 void gatherData(DataHandle& data,
CpGridData* global_view,
602 template<
int codim,
class DataHandle>
603 void gatherCodimData(DataHandle& data,
CpGridData* global_data,
612 template<
class DataHandle>
613 void scatterData(DataHandle& data,
CpGridData* global_data,
614 CpGridData* distributed_data,
const InterfaceMap& cell_inf,
615 const InterfaceMap& point_inf);
624 template<
int codim,
class DataHandle>
625 void scatterCodimData(DataHandle& data,
CpGridData* global_data,
636 template<
int codim,
class DataHandle>
638 const Interface& interface);
648 template<
int codim,
class DataHandle>
650 const InterfaceMap& interface);
654 void computeGeometry(
CpGrid& grid,
656 const std::vector<int>& globalAquiferCells,
659 std::vector<int>& aquiferCells,
661 const std::vector< std::array<int,8> >& cell2Points);
677 std::vector< std::array<int,8> > cell_to_point_;
684 std::array<int, 3> logical_cartesian_size_;
691 std::vector<int> global_cell_;
697 typedef FieldVector<double, 3> PointType;
703 std::unique_ptr<cpgrid::IndexSet> index_set_;
705 std::shared_ptr<const cpgrid::IdSet> local_id_set_;
707 std::shared_ptr<LevelGlobalIdSet> global_id_set_;
709 std::shared_ptr<PartitionTypeIndicator> partition_type_indicator_;
713 std::vector<std::shared_ptr<CpGridData>>* level_data_ptr_;
716 std::vector<int> level_to_leaf_cells_;
// {level LGR, {child0, child1, ...}}
718 std::vector<std::tuple<int,std::vector<int>>> parent_to_children_cells_;
720 std::array<int,3> cells_per_dim_;
// {level, cell index in that level}
723 std::vector<std::array<int,2>> leaf_to_level_cells_;
// {level parent cell, parent cell index}
726 std::vector<std::array<int,2>> child_to_parent_cells_;
733 bool use_unique_boundary_ids_;
740 std::vector<double> zcorn;
743 std::vector<int> aquifer_cells_;
748 CommunicationType cell_comm_;
751 std::tuple<Interface,Interface,Interface,Interface,Interface> cell_interfaces_;
756 std::tuple<InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap>
760 std::tuple<InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap>
767 const EntityVariable<Geometry<3 - codim, 3>, codim>& geomVector()
const
773 template<
int>
friend class Entity;
774 template<
int>
friend class EntityRep;
775 friend class Intersection;
776 friend class PartitionTypeIndicator;
790T& getInterface(InterfaceType iftype,
791 std::tuple<T,T,T,T,T>& interfaces)
796 return std::get<0>(interfaces);
798 return std::get<1>(interfaces);
800 return std::get<2>(interfaces);
802 return std::get<3>(interfaces);
804 return std::get<4>(interfaces);
806 OPM_THROW(std::runtime_error,
"Invalid Interface type was used during communication");
811template<
int codim,
class DataHandle>
812void CpGridData::communicateCodim(Entity2IndexDataHandle<DataHandle, codim>& data, CommunicationDirection dir,
813 const Interface& interface)
815 this->
template communicateCodim<codim>(data, dir, interface.interfaces());
818template<
int codim,
class DataHandle>
819void CpGridData::communicateCodim(Entity2IndexDataHandle<DataHandle, codim>& data_wrapper, CommunicationDirection dir,
820 const InterfaceMap& interface)
822 Communicator comm(ccobj_, interface);
824 if(dir==ForwardCommunication)
825 comm.forward(data_wrapper);
827 comm.backward(data_wrapper);
831template<
class DataHandle>
833 CommunicationDirection dir)
836 if(data.contains(3,0))
839 communicateCodim<0>(data_wrapper, dir, getInterface(iftype, cell_interfaces_));
841 if(data.contains(3,3))
844 communicateCodim<3>(data_wrapper, dir, getInterface(iftype, point_interfaces_));
856#include <opm/grid/cpgrid/Entity.hpp>
857#include <opm/grid/cpgrid/Indexsets.hpp>
871 data=buffer_[index_++];
873 void write(
const T& data)
875 buffer_[index_++]=data;
881 void resize(std::size_t size)
883 buffer_.resize(size);
887 std::vector<T> buffer_;
888 typename std::vector<T>::size_type index_;
890template<
class DataHandle,
int codim>
895template<
class DataHandle>
898 BaseMover(DataHandle& data)
902 void moveData(
const E& from,
const E& to)
904 std::size_t size=data_.size(from);
906 data_.gather(buffer, from);
908 data_.scatter(buffer, to, size);
911 MoveBuffer<typename DataHandle::DataType> buffer;
915template<
class DataHandle>
916struct Mover<DataHandle,0> :
public BaseMover<DataHandle>
918 Mover<DataHandle,0>(DataHandle& data, CpGridData* gatherView,
919 CpGridData* scatterView)
920 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
923 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
925 Entity<0> from_entity=Entity<0>(*gatherView_, from_cell_index,
true);
926 Entity<0> to_entity=Entity<0>(*scatterView_, to_cell_index,
true);
927 this->moveData(from_entity, to_entity);
929 CpGridData* gatherView_;
930 CpGridData* scatterView_;
933template<
class DataHandle>
934struct Mover<DataHandle,1> :
public BaseMover<DataHandle>
936 Mover<DataHandle,1>(DataHandle& data, CpGridData* gatherView,
937 CpGridData* scatterView)
938 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
941 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
943 typedef typename OrientedEntityTable<0,1>::row_type row_type;
944 EntityRep<0> from_cell=EntityRep<0>(from_cell_index,
true);
945 EntityRep<0> to_cell=EntityRep<0>(to_cell_index,
true);
946 OrientedEntityTable<0,1>& table = gatherView_->cell_to_face_;
947 row_type from_faces=table.operator[](from_cell);
948 row_type to_faces=scatterView_->cell_to_face_[to_cell];
950 for(
int i=0; i<from_faces.size(); ++i)
951 this->moveData(from_faces[i], to_faces[i]);
953 CpGridData *gatherView_;
954 CpGridData *scatterView_;
957template<
class DataHandle>
958struct Mover<DataHandle,3> :
public BaseMover<DataHandle>
960 Mover<DataHandle,3>(DataHandle& data, CpGridData* gatherView,
961 CpGridData* scatterView)
962 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
964 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
966 const std::array<int,8>& from_cell_points=
967 gatherView_->cell_to_point_[from_cell_index];
968 const std::array<int,8>& to_cell_points=
969 scatterView_->cell_to_point_[to_cell_index];
970 for(std::size_t i=0; i<8; ++i)
972 this->moveData(Entity<3>(*gatherView_, from_cell_points[i],
true),
973 Entity<3>(*scatterView_, to_cell_points[i],
true));
976 CpGridData* gatherView_;
977 CpGridData* scatterView_;
982template<
class DataHandle>
983void CpGridData::scatterData(DataHandle& data, CpGridData* global_data,
984 CpGridData* distributed_data,
const InterfaceMap& cell_inf,
985 const InterfaceMap& point_inf)
988 if(data.contains(3,0))
990 Entity2IndexDataHandle<DataHandle, 0> data_wrapper(*global_data, *distributed_data, data);
991 communicateCodim<0>(data_wrapper, ForwardCommunication, cell_inf);
993 if(data.contains(3,3))
995 Entity2IndexDataHandle<DataHandle, 3> data_wrapper(*global_data, *distributed_data, data);
996 communicateCodim<3>(data_wrapper, ForwardCommunication, point_inf);
1001template<
int codim,
class DataHandle>
1002void CpGridData::scatterCodimData(DataHandle& data, CpGridData* global_data,
1003 CpGridData* distributed_data)
1005 CpGridData *gather_view, *scatter_view;
1006 gather_view=global_data;
1007 scatter_view=distributed_data;
1009 mover::Mover<DataHandle,codim> mover(data, gather_view, scatter_view);
1012 for(
auto index=distributed_data->cellIndexSet().begin(),
1013 end = distributed_data->cellIndexSet().end();
1014 index!=end; ++index)
1016 std::size_t from=index->global();
1017 std::size_t to=index->local();
1025template<
int codim,
class T,
class F>
1026void visitInterior(CpGridData& distributed_data, T begin, T endit, F& func)
1028 for(T it=begin; it!=endit; ++it)
1030 Entity<codim> entity(distributed_data, it-begin,
true);
1031 PartitionType pt = entity.partitionType();
1032 if(pt==Dune::InteriorEntity)
1039template<
class DataHandle>
1040struct GlobalIndexSizeGatherer
1042 GlobalIndexSizeGatherer(DataHandle& data_,
1043 std::vector<int>& ownedGlobalIndices_,
1044 std::vector<int>& ownedSizes_)
1045 : data(data_), ownedGlobalIndices(ownedGlobalIndices_), ownedSizes(ownedSizes_)
1048 template<
class T,
class E>
1049 void operator()(T& i, E& entity)
1051 ownedGlobalIndices.push_back(i);
1052 ownedSizes.push_back(data.size(entity));
1055 std::vector<int>& ownedGlobalIndices;
1056 std::vector<int>& ownedSizes;
1059template<
class DataHandle>
1062 DataGatherer(mover::MoveBuffer<typename DataHandle::DataType>& buffer_,
1064 : buffer(buffer_), data(data_)
1067 template<
class T,
class E>
1068 void operator()(T& , E& entity)
1070 data.gather(buffer, entity);
1072 mover::MoveBuffer<typename DataHandle::DataType>& buffer;
1078template<
class DataHandle>
1079void CpGridData::gatherData(DataHandle& data, CpGridData* global_data,
1080 CpGridData* distributed_data)
1083 if(data.contains(3,0))
1084 gatherCodimData<0>(data, global_data, distributed_data);
1085 if(data.contains(3,3))
1086 gatherCodimData<3>(data, global_data, distributed_data);
1090template<
int codim,
class DataHandle>
1091void CpGridData::gatherCodimData(DataHandle& data, CpGridData* global_data,
1092 CpGridData* distributed_data)
1096 const std::vector<int>& mapping =
1097 distributed_data->global_id_set_->getMapping<codim>();
1101 std::vector<int> owned_global_indices;
1102 std::vector<int> owned_sizes;
1103 owned_global_indices.reserve(mapping.size());
1104 owned_sizes.reserve(mapping.size());
1106 GlobalIndexSizeGatherer<DataHandle> gisg(data, owned_global_indices, owned_sizes);
1107 visitInterior<codim>(*distributed_data, mapping.begin(), mapping.end(), gisg);
1110 int no_indices=owned_sizes.size();
1113 if ( owned_global_indices.empty() )
1114 owned_global_indices.resize(1);
1115 if ( owned_sizes.empty() )
1116 owned_sizes.resize(1);
1117 std::vector<int> no_indices_to_recv(distributed_data->ccobj_.size());
1118 distributed_data->ccobj_.allgather(&no_indices, 1, &(no_indices_to_recv[0]));
1122 std::vector<int> displ(distributed_data->ccobj_.size()+1, 0);
1123 std::transform(displ.begin(), displ.end()-1, no_indices_to_recv.begin(), displ.begin()+1,
1125 int global_size=displ[displ.size()-1];
1126 std::vector<int> global_indices(global_size);
1127 std::vector<int> global_sizes(global_size);
1128 MPI_Allgatherv(&(owned_global_indices[0]), no_indices, MPITraits<int>::getType(),
1129 &(global_indices[0]), &(no_indices_to_recv[0]), &(displ[0]),
1130 MPITraits<int>::getType(),
1131 distributed_data->ccobj_);
1132 MPI_Allgatherv(&(owned_sizes[0]), no_indices, MPITraits<int>::getType(),
1133 &(global_sizes[0]), &(no_indices_to_recv[0]), &(displ[0]),
1134 MPITraits<int>::getType(),
1135 distributed_data->ccobj_);
1136 std::vector<int>().swap(owned_global_indices);
1138 std::vector<int> no_data_send(distributed_data->ccobj_.size());
1139 for(
typename std::vector<int>::iterator begin=no_data_send.begin(),
1140 i=begin, end=no_data_send.end(); i!=end; ++i)
1141 *i = std::accumulate(global_sizes.begin()+displ[i-begin],
1142 global_sizes.begin()+displ[i-begin+1], std::size_t());
1144 std::vector<int>().swap(owned_sizes);
1147 std::transform(displ.begin(), displ.end()-1, no_data_send.begin(), displ.begin()+1,
1148 std::plus<std::size_t>());
1150 int no_data_recv = displ[displ.size()-1];
1153 mover::MoveBuffer<typename DataHandle::DataType> local_data_buffer, global_data_buffer;
1154 if ( no_data_send[distributed_data->ccobj_.rank()] )
1156 local_data_buffer.resize(no_data_send[distributed_data->ccobj_.rank()]);
1160 local_data_buffer.resize(1);
1162 global_data_buffer.resize(no_data_recv);
1164 DataGatherer<DataHandle> gatherer(local_data_buffer, data);
1165 visitInterior<codim>(*distributed_data, mapping.begin(), mapping.end(), gatherer);
1166 MPI_Allgatherv(&(local_data_buffer.buffer_[0]), no_data_send[distributed_data->ccobj_.rank()],
1167 MPITraits<typename DataHandle::DataType>::getType(),
1168 &(global_data_buffer.buffer_[0]), &(no_data_send[0]), &(displ[0]),
1169 MPITraits<typename DataHandle::DataType>::getType(),
1170 distributed_data->ccobj_);
1171 Entity2IndexDataHandle<DataHandle, codim> edata(*global_data, data);
1173 for(
int i=0; i< codim; ++i)
1174 offset+=global_data->size(i);
1176 typename std::vector<int>::const_iterator s=global_sizes.begin();
1177 for(
typename std::vector<int>::const_iterator i=global_indices.begin(),
1178 end=global_indices.end();
1181 edata.scatter(global_data_buffer, *i-offset, *s);
[ provides Dune::Grid ]
Definition CpGrid.hpp:226
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:135
@ MAX_DATA_PER_CELL
The maximum data items allowed per cell (DUNE < 2.5.2)
Definition CpGridData.hpp:182
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition CpGridData.hpp:211
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the grid.
Definition CpGridData.hpp:492
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition CpGridData.hpp:832
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition CpGridData.hpp:460
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition readSintefLegacyFormat.cpp:67
int size(int codim) const
number of leaf entities per codim in this process
Definition CpGridData.cpp:144
void readEclipseFormat(const std::string &filename, bool periodic_extension, bool turn_normals=false)
Read the Eclipse grid format ('grdecl').
CpGridDataTraits::Communication Communication
The type of the collective communication.
Definition CpGridData.hpp:515
~CpGridData()
Destructor.
Definition CpGridData.cpp:97
void getIJK(int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition CpGridData.hpp:291
const IndexSet & indexSet() const
Get the index set.
Definition CpGridData.hpp:485
std::tuple< std::shared_ptr< CpGridData >, const std::vector< std::array< int, 2 > >, const std::vector< std::tuple< int, std::vector< int > > >, const std::vector< std::tuple< int, std::vector< int > > >, const std::vector< std::tuple< int, std::vector< int > > >, const std::vector< std::array< int, 2 > >, const std::vector< std::array< int, 2 > > > refinePatch(const std::array< int, 3 > &cells_per_dim, const std::array< int, 3 > &startIJK, const std::array< int, 3 > &endIJK) const
Refine a (connected block-shaped) patch of cells.
Definition CpGridData.cpp:2062
CpGridDataTraits::MPICommunicator MPICommunicator
The type of the mpi communicator.
Definition CpGridData.hpp:513
void processEclipseFormat(const grdecl &input_data, std::array< std::set< std::pair< int, int > >, 2 > &nnc, bool remove_ij_boundary, bool turn_normals, bool pinchActive, double tolerance_unique_points)
Read the Eclipse grid format ('grdecl').
Definition processEclipseFormat.cpp:268
std::tuple< const std::shared_ptr< CpGridData >, const std::vector< std::array< int, 2 > >, const std::vector< std::tuple< int, std::vector< int > > >, const std::tuple< int, std::vector< int > >, const std::vector< std::array< int, 2 > >, const std::vector< std::array< int, 2 > > > refineSingleCell(const std::array< int, 3 > &cells_per_dim, const int &parent_idx) const
Refine a single cell and return a shared pointer of CpGridData type.
Definition CpGridData.cpp:1933
void distributeGlobalGrid(CpGrid &grid, const CpGridData &view_data, const std::vector< int > &cell_part)
Redistribute a global grid.
Definition CpGridData.cpp:1513
const std::vector< double > & zcornData() const
Return the internalized zcorn copy from the grid processing, if no cells were adjusted during the min...
Definition CpGridData.hpp:478
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition CpGridData.hpp:467
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition writeSintefLegacyFormat.cpp:73
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition CpGridData.hpp:574
Definition DefaultGeometryPolicy.hpp:53
const EntityVariable< cpgrid::Geometry< 3 - codim, 3 >, codim > & geomVector() const
Definition DefaultGeometryPolicy.hpp:86
Wrapper that turns a data handle suitable for dune-grid into one based on integers instead of entitie...
Definition Entity2IndexDataHandle.hpp:56
A class design to hold a variable with a value for each entity of the given codimension,...
Definition EntityRep.hpp:267
This class encapsulates geometry for vertices, intersections, and cells.
Definition Geometry.hpp:75
The global id set for Dune.
Definition Indexsets.hpp:304
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Indexsets.hpp:55
Represents the topological relationships between sets of entities, for example cells and faces.
Definition OrientedEntityTable.hpp:139
A class design to hold a variable with a value for each entity of the given codimension,...
Definition EntityRep.hpp:302
A SparseTable stores a table with rows of varying size as efficiently as possible.
Definition SparseTable.hpp:55
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10
Holds the implementation of the CpGrid as a pimple.
Definition CellQuadrature.cpp:68
Low-level corner-point processing routines and supporting data structures.
MPIHelper::MPICommunicator MPICommunicator
The type of the collective communication.
Definition CpGridDataTraits.hpp:55
AttributeSet
The type of the set of the attributes.
Definition CpGridDataTraits.hpp:65
Definition CpGridData.hpp:127
Raw corner-point specification of a particular geological model.
Definition preprocess.h:56