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::array<int, 3>&,
102void refinePatch_and_check(
const std::array<int,3>&,
103 const std::array<int,3>&,
104 const std::array<int,3>&);
107 const std::vector<std::array<int,3>>&,
108 const std::vector<std::array<int,3>>&,
109 const std::vector<std::array<int,3>>&,
110 const std::vector<std::string>&);
115void fieldProp_check(
const Dune::CpGrid& grid, Opm::EclipseGrid eclGrid, std::string deck_string);
123template<
class T,
int i>
struct Mover;
139 const std::array<int, 3>&,
142 void ::refinePatch_and_check(
const std::array<int,3>&,
143 const std::array<int,3>&,
144 const std::array<int,3>&);
148 const std::vector<std::array<int,3>>&,
149 const std::vector<std::array<int,3>>&,
150 const std::vector<std::array<int,3>>&,
151 const std::vector<std::string>&);
158 void ::fieldProp_check(
const Dune::CpGrid& grid, Opm::EclipseGrid eclGrid, std::string deck_string);
165#ifndef MAX_DATA_COMMUNICATED_PER_ENTITY
185 explicit CpGridData(MPIHelper::MPICommunicator comm, std::vector<std::shared_ptr<CpGridData>>& data);
190 explicit CpGridData(std::vector<std::shared_ptr<CpGridData>>& data);
198 int size(
int codim)
const;
201 int size (GeometryType type)
const
204 return size(3 - type.dim());
224 void readEclipseFormat(
const std::string& filename,
bool periodic_extension,
bool turn_normals =
false);
235 void processEclipseFormat(
const Opm::Deck& deck,
bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false,
236 const std::vector<double>& poreVolume = std::vector<double>());
250 std::vector<std::size_t>
processEclipseFormat(
const Opm::EclipseGrid* ecl_grid, Opm::EclipseState* ecl_state,
251 bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false,
252 bool pinchActive =
true);
268 Opm::EclipseState* ecl_state,
270 std::array<std::set<std::pair<int, int>>, 2>& nnc,
271 bool remove_ij_boundary,
bool turn_normals,
bool pinchActive,
272 double tolerance_unique_points);
281 void getIJK(
int c, std::array<int,3>& ijk)
const
283 int gc = global_cell_[c];
284 ijk[0] = gc % logical_cartesian_size_[0]; gc /= logical_cartesian_size_[0];
285 ijk[1] = gc % logical_cartesian_size_[1];
286 ijk[2] = gc / logical_cartesian_size_[1];
294 bool disjointPatches(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
304 getPatchesCells(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
308 bool hasNNCs(
const std::vector<int>& cellIndices)
const;
316 void validStartEndIJKs(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
326 bool patchesShareFace(
const std::vector<std::array<int,3>>& startIJK_vec,
const std::vector<std::array<int,3>>& endIJK_vec)
const;
336 std::array<int,3> getPatchDim(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
345 std::vector<int> getPatchCorners(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
354 std::vector<int> getPatchFaces(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
363 std::vector<int> getPatchCells(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
372 std::vector<int> getPatchBoundaryCorners(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
381 std::array<std::vector<int>,6> getBoundaryPatchFaces(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
384 std::array<std::vector<double>,3> getWidthsLengthsHeights(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
400 Geometry<3,3> cellifyPatch(
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK,
402 std::array<int,8>& cellifiedPatch_to_point,
403 std::array<int,8>& allcorners_cellifiedPatch)
const;
409 std::array<double,3> getAverageArr(
const std::vector<std::array<double,3>>& vec)
const;
433 std::tuple< const std::shared_ptr<CpGridData>,
434 const std::vector<std::array<int,2>>,
435 const std::vector<std::tuple<int,std::vector<int>>>,
436 const std::tuple<int, std::vector<int>>,
437 const std::vector<std::array<int,2>>,
438 const std::vector<std::array<int,2>>>
439 refineSingleCell(
const std::array<int,3>& cells_per_dim,
const int& parent_idx)
const;
453 std::tuple< std::shared_ptr<CpGridData>,
454 const std::vector<std::array<int,2>>,
455 const std::vector<std::tuple<int,std::vector<int>>>,
456 const std::vector<std::tuple<int,std::vector<int>>>,
457 const std::vector<std::tuple<int,std::vector<int>>>,
458 const std::vector<std::array<int,2>>,
459 const std::vector<std::array<int,2>>>
460 refinePatch(
const std::array<int,3>& cells_per_dim,
const std::array<int,3>& startIJK,
const std::array<int,3>& endIJK)
const;
468 std::array<double,3> computeEclCentroid(
const int idx)
const;
476 std::array<double,3> computeEclCentroid(
const Entity<0>& elem)
const;
479 void computeUniqueBoundaryIds();
486 return use_unique_boundary_ids_;
493 use_unique_boundary_ids_ = uids;
494 if (use_unique_boundary_ids_ && unique_boundary_ids_.empty()) {
495 computeUniqueBoundaryIds();
517 return *local_id_set_;
525 return logical_cartesian_size_;
533 const std::vector<int>& cell_part);
540 template<
class DataHandle>
541 void communicate(DataHandle& data, InterfaceType iftype, CommunicationDirection dir);
547 using CollectiveCommunication = CpGridDataTraits::CollectiveCommunication;
553 using Communicator = CpGridDataTraits::Communicator;
556 using InterfaceMap = CpGridDataTraits::InterfaceMap;
559 using CommunicationType = CpGridDataTraits::CommunicationType;
562 using ParallelIndexSet = CpGridDataTraits::ParallelIndexSet;
565 using RemoteIndices = CpGridDataTraits::RemoteIndices;
570 CommunicationType& cellCommunication()
578 const CommunicationType& cellCommunication()
const
583 ParallelIndexSet& cellIndexSet()
585 return cellCommunication().indexSet();
588 const ParallelIndexSet& cellIndexSet()
const
590 return cellCommunication().indexSet();
593 RemoteIndices& cellRemoteIndices()
595 return cellCommunication().remoteIndices();
598 const RemoteIndices& cellRemoteIndices()
const
600 return cellCommunication().remoteIndices();
607 return aquifer_cells_;
613 void populateGlobalCellIndexSet();
622 template<
class DataHandle>
623 void gatherData(DataHandle& data,
CpGridData* global_view,
633 template<
int codim,
class DataHandle>
634 void gatherCodimData(DataHandle& data,
CpGridData* global_data,
643 template<
class DataHandle>
644 void scatterData(DataHandle& data,
CpGridData* global_data,
645 CpGridData* distributed_data,
const InterfaceMap& cell_inf,
646 const InterfaceMap& point_inf);
655 template<
int codim,
class DataHandle>
656 void scatterCodimData(DataHandle& data,
CpGridData* global_data,
667 template<
int codim,
class DataHandle>
669 const Interface& interface);
679 template<
int codim,
class DataHandle>
681 const InterfaceMap& interface);
685 void computeGeometry(
CpGrid& grid,
687 const std::vector<int>& globalAquiferCells,
690 std::vector<int>& aquiferCells,
692 const std::vector< std::array<int,8> >& cell2Points);
708 std::vector< std::array<int,8> > cell_to_point_;
715 std::array<int, 3> logical_cartesian_size_;
722 std::vector<int> global_cell_;
728 typedef FieldVector<double, 3> PointType;
734 std::unique_ptr<cpgrid::IndexSet> index_set_;
736 std::shared_ptr<const cpgrid::IdSet> local_id_set_;
738 std::shared_ptr<LevelGlobalIdSet> global_id_set_;
740 std::shared_ptr<PartitionTypeIndicator> partition_type_indicator_;
744 std::vector<std::shared_ptr<CpGridData>>* level_data_ptr_;
747 std::vector<int> level_to_leaf_cells_;
// {level LGR, {child0, child1, ...}}
749 std::vector<std::tuple<int,std::vector<int>>> parent_to_children_cells_;
751 std::array<int,3> cells_per_dim_;
// {level, cell index in that level}
754 std::vector<std::array<int,2>> leaf_to_level_cells_;
// {level parent cell, parent cell index}
757 std::vector<std::array<int,2>> child_to_parent_cells_;
764 bool use_unique_boundary_ids_;
771 std::vector<double> zcorn;
774 std::vector<int> aquifer_cells_;
779 CommunicationType cell_comm_;
782 std::tuple<Interface,Interface,Interface,Interface,Interface> cell_interfaces_;
787 std::tuple<InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap>
791 std::tuple<InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap,InterfaceMap>
798 const EntityVariable<Geometry<3 - codim, 3>, codim>& geomVector()
const
804 template<
int>
friend class Entity;
805 template<
int>
friend class EntityRep;
806 friend class Intersection;
807 friend class PartitionTypeIndicator;
821T& getInterface(InterfaceType iftype,
822 std::tuple<T,T,T,T,T>& interfaces)
827 return std::get<0>(interfaces);
829 return std::get<1>(interfaces);
831 return std::get<2>(interfaces);
833 return std::get<3>(interfaces);
835 return std::get<4>(interfaces);
837 OPM_THROW(std::runtime_error,
"Invalid Interface type was used during communication");
842template<
int codim,
class DataHandle>
843void CpGridData::communicateCodim(Entity2IndexDataHandle<DataHandle, codim>& data, CommunicationDirection dir,
844 const Interface& interface)
846 this->
template communicateCodim<codim>(data, dir, interface.interfaces());
849template<
int codim,
class DataHandle>
850void CpGridData::communicateCodim(Entity2IndexDataHandle<DataHandle, codim>& data_wrapper, CommunicationDirection dir,
851 const InterfaceMap& interface)
853 Communicator comm(ccobj_, interface);
855 if(dir==ForwardCommunication)
856 comm.forward(data_wrapper);
858 comm.backward(data_wrapper);
862template<
class DataHandle>
864 CommunicationDirection dir)
867 if(data.contains(3,0))
870 communicateCodim<0>(data_wrapper, dir, getInterface(iftype, cell_interfaces_));
872 if(data.contains(3,3))
875 communicateCodim<3>(data_wrapper, dir, getInterface(iftype, point_interfaces_));
887#include <opm/grid/cpgrid/Entity.hpp>
888#include <opm/grid/cpgrid/Indexsets.hpp>
902 data=buffer_[index_++];
904 void write(
const T& data)
906 buffer_[index_++]=data;
912 void resize(std::size_t size)
914 buffer_.resize(size);
918 std::vector<T> buffer_;
919 typename std::vector<T>::size_type index_;
921template<
class DataHandle,
int codim>
926template<
class DataHandle>
929 BaseMover(DataHandle& data)
933 void moveData(
const E& from,
const E& to)
935 std::size_t size=data_.size(from);
937 data_.gather(buffer, from);
939 data_.scatter(buffer, to, size);
942 MoveBuffer<typename DataHandle::DataType> buffer;
946template<
class DataHandle>
947struct Mover<DataHandle,0> :
public BaseMover<DataHandle>
949 Mover<DataHandle,0>(DataHandle& data, CpGridData* gatherView,
950 CpGridData* scatterView)
951 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
954 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
956 Entity<0> from_entity=Entity<0>(*gatherView_, from_cell_index,
true);
957 Entity<0> to_entity=Entity<0>(*scatterView_, to_cell_index,
true);
958 this->moveData(from_entity, to_entity);
960 CpGridData* gatherView_;
961 CpGridData* scatterView_;
964template<
class DataHandle>
965struct Mover<DataHandle,1> :
public BaseMover<DataHandle>
967 Mover<DataHandle,1>(DataHandle& data, CpGridData* gatherView,
968 CpGridData* scatterView)
969 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
972 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
974 typedef typename OrientedEntityTable<0,1>::row_type row_type;
975 EntityRep<0> from_cell=EntityRep<0>(from_cell_index,
true);
976 EntityRep<0> to_cell=EntityRep<0>(to_cell_index,
true);
977 OrientedEntityTable<0,1>& table = gatherView_->cell_to_face_;
978 row_type from_faces=table.operator[](from_cell);
979 row_type to_faces=scatterView_->cell_to_face_[to_cell];
981 for(
int i=0; i<from_faces.size(); ++i)
982 this->moveData(from_faces[i], to_faces[i]);
984 CpGridData *gatherView_;
985 CpGridData *scatterView_;
988template<
class DataHandle>
989struct Mover<DataHandle,3> :
public BaseMover<DataHandle>
991 Mover<DataHandle,3>(DataHandle& data, CpGridData* gatherView,
992 CpGridData* scatterView)
993 : BaseMover<DataHandle>(data), gatherView_(gatherView), scatterView_(scatterView)
995 void operator()(std::size_t from_cell_index,std::size_t to_cell_index)
997 const std::array<int,8>& from_cell_points=
998 gatherView_->cell_to_point_[from_cell_index];
999 const std::array<int,8>& to_cell_points=
1000 scatterView_->cell_to_point_[to_cell_index];
1001 for(std::size_t i=0; i<8; ++i)
1003 this->moveData(Entity<3>(*gatherView_, from_cell_points[i],
true),
1004 Entity<3>(*scatterView_, to_cell_points[i],
true));
1007 CpGridData* gatherView_;
1008 CpGridData* scatterView_;
1013template<
class DataHandle>
1014void CpGridData::scatterData(DataHandle& data, CpGridData* global_data,
1015 CpGridData* distributed_data,
const InterfaceMap& cell_inf,
1016 const InterfaceMap& point_inf)
1019 if(data.contains(3,0))
1021 Entity2IndexDataHandle<DataHandle, 0> data_wrapper(*global_data, *distributed_data, data);
1022 communicateCodim<0>(data_wrapper, ForwardCommunication, cell_inf);
1024 if(data.contains(3,3))
1026 Entity2IndexDataHandle<DataHandle, 3> data_wrapper(*global_data, *distributed_data, data);
1027 communicateCodim<3>(data_wrapper, ForwardCommunication, point_inf);
1032template<
int codim,
class DataHandle>
1033void CpGridData::scatterCodimData(DataHandle& data, CpGridData* global_data,
1034 CpGridData* distributed_data)
1036 CpGridData *gather_view, *scatter_view;
1037 gather_view=global_data;
1038 scatter_view=distributed_data;
1040 mover::Mover<DataHandle,codim> mover(data, gather_view, scatter_view);
1043 for(
auto index=distributed_data->cellIndexSet().begin(),
1044 end = distributed_data->cellIndexSet().end();
1045 index!=end; ++index)
1047 std::size_t from=index->global();
1048 std::size_t to=index->local();
1056template<
int codim,
class T,
class F>
1057void visitInterior(CpGridData& distributed_data, T begin, T endit, F& func)
1059 for(T it=begin; it!=endit; ++it)
1061 Entity<codim> entity(distributed_data, it-begin,
true);
1062 PartitionType pt = entity.partitionType();
1063 if(pt==Dune::InteriorEntity)
1070template<
class DataHandle>
1071struct GlobalIndexSizeGatherer
1073 GlobalIndexSizeGatherer(DataHandle& data_,
1074 std::vector<int>& ownedGlobalIndices_,
1075 std::vector<int>& ownedSizes_)
1076 : data(data_), ownedGlobalIndices(ownedGlobalIndices_), ownedSizes(ownedSizes_)
1079 template<
class T,
class E>
1080 void operator()(T& i, E& entity)
1082 ownedGlobalIndices.push_back(i);
1083 ownedSizes.push_back(data.size(entity));
1086 std::vector<int>& ownedGlobalIndices;
1087 std::vector<int>& ownedSizes;
1090template<
class DataHandle>
1093 DataGatherer(mover::MoveBuffer<typename DataHandle::DataType>& buffer_,
1095 : buffer(buffer_), data(data_)
1098 template<
class T,
class E>
1099 void operator()(T& , E& entity)
1101 data.gather(buffer, entity);
1103 mover::MoveBuffer<typename DataHandle::DataType>& buffer;
1109template<
class DataHandle>
1110void CpGridData::gatherData(DataHandle& data, CpGridData* global_data,
1111 CpGridData* distributed_data)
1114 if(data.contains(3,0))
1115 gatherCodimData<0>(data, global_data, distributed_data);
1116 if(data.contains(3,3))
1117 gatherCodimData<3>(data, global_data, distributed_data);
1121template<
int codim,
class DataHandle>
1122void CpGridData::gatherCodimData(DataHandle& data, CpGridData* global_data,
1123 CpGridData* distributed_data)
1127 const std::vector<int>& mapping =
1128 distributed_data->global_id_set_->getMapping<codim>();
1132 std::vector<int> owned_global_indices;
1133 std::vector<int> owned_sizes;
1134 owned_global_indices.reserve(mapping.size());
1135 owned_sizes.reserve(mapping.size());
1137 GlobalIndexSizeGatherer<DataHandle> gisg(data, owned_global_indices, owned_sizes);
1138 visitInterior<codim>(*distributed_data, mapping.begin(), mapping.end(), gisg);
1141 int no_indices=owned_sizes.size();
1144 if ( owned_global_indices.empty() )
1145 owned_global_indices.resize(1);
1146 if ( owned_sizes.empty() )
1147 owned_sizes.resize(1);
1148 std::vector<int> no_indices_to_recv(distributed_data->ccobj_.size());
1149 distributed_data->ccobj_.allgather(&no_indices, 1, &(no_indices_to_recv[0]));
1153 std::vector<int> displ(distributed_data->ccobj_.size()+1, 0);
1154 std::transform(displ.begin(), displ.end()-1, no_indices_to_recv.begin(), displ.begin()+1,
1156 int global_size=displ[displ.size()-1];
1157 std::vector<int> global_indices(global_size);
1158 std::vector<int> global_sizes(global_size);
1159 MPI_Allgatherv(&(owned_global_indices[0]), no_indices, MPITraits<int>::getType(),
1160 &(global_indices[0]), &(no_indices_to_recv[0]), &(displ[0]),
1161 MPITraits<int>::getType(),
1162 distributed_data->ccobj_);
1163 MPI_Allgatherv(&(owned_sizes[0]), no_indices, MPITraits<int>::getType(),
1164 &(global_sizes[0]), &(no_indices_to_recv[0]), &(displ[0]),
1165 MPITraits<int>::getType(),
1166 distributed_data->ccobj_);
1167 std::vector<int>().swap(owned_global_indices);
1169 std::vector<int> no_data_send(distributed_data->ccobj_.size());
1170 for(
typename std::vector<int>::iterator begin=no_data_send.begin(),
1171 i=begin, end=no_data_send.end(); i!=end; ++i)
1172 *i = std::accumulate(global_sizes.begin()+displ[i-begin],
1173 global_sizes.begin()+displ[i-begin+1], std::size_t());
1175 std::vector<int>().swap(owned_sizes);
1178 std::transform(displ.begin(), displ.end()-1, no_data_send.begin(), displ.begin()+1,
1179 std::plus<std::size_t>());
1181 int no_data_recv = displ[displ.size()-1];
1184 mover::MoveBuffer<typename DataHandle::DataType> local_data_buffer, global_data_buffer;
1185 if ( no_data_send[distributed_data->ccobj_.rank()] )
1187 local_data_buffer.resize(no_data_send[distributed_data->ccobj_.rank()]);
1191 local_data_buffer.resize(1);
1193 global_data_buffer.resize(no_data_recv);
1195 DataGatherer<DataHandle> gatherer(local_data_buffer, data);
1196 visitInterior<codim>(*distributed_data, mapping.begin(), mapping.end(), gatherer);
1197 MPI_Allgatherv(&(local_data_buffer.buffer_[0]), no_data_send[distributed_data->ccobj_.rank()],
1198 MPITraits<typename DataHandle::DataType>::getType(),
1199 &(global_data_buffer.buffer_[0]), &(no_data_send[0]), &(displ[0]),
1200 MPITraits<typename DataHandle::DataType>::getType(),
1201 distributed_data->ccobj_);
1202 Entity2IndexDataHandle<DataHandle, codim> edata(*global_data, data);
1204 for(
int i=0; i< codim; ++i)
1205 offset+=global_data->size(i);
1207 typename std::vector<int>::const_iterator s=global_sizes.begin();
1208 for(
typename std::vector<int>::const_iterator i=global_indices.begin(),
1209 end=global_indices.end();
1212 edata.scatter(global_data_buffer, *i-offset, *s);
[ provides Dune::Grid ]
Definition CpGrid.hpp:238
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:131
@ MAX_DATA_PER_CELL
The maximum data items allowed per cell (DUNE < 2.5.2)
Definition CpGridData.hpp:172
std::vector< int > getPatchesCells(const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec) const
Compute cell indices of selected patches of cells (Cartesian grid required).
Definition CpGridData.cpp:2041
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition CpGridData.hpp:201
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the grid.
Definition CpGridData.hpp:523
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition CpGridData.hpp:863
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition CpGridData.hpp:484
bool patchesShareFace(const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec) const
Determine if a finite amount of patches (of cells) share a face.
Definition CpGridData.cpp:1974
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:146
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:546
~CpGridData()
Destructor.
Definition CpGridData.cpp:99
void getIJK(int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition CpGridData.hpp:281
const IndexSet & indexSet() const
Get the index set.
Definition CpGridData.hpp:509
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:2295
CpGridDataTraits::MPICommunicator MPICommunicator
The type of the mpi communicator.
Definition CpGridData.hpp:544
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:288
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:2159
void distributeGlobalGrid(CpGrid &grid, const CpGridData &view_data, const std::vector< int > &cell_part)
Redistribute a global grid.
Definition CpGridData.cpp:1515
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:502
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition CpGridData.hpp:491
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition writeSintefLegacyFormat.cpp:73
bool hasNNCs(const std::vector< int > &cellIndices) const
Check all cells selected for refinement have no NNCs (no neighbor connections).
Definition CpGridData.cpp:2052
const cpgrid::IdSet & localIdSet() const
Get the local index set.
Definition CpGridData.hpp:515
void checkCuboidShape(const std::vector< int > &cellIdx_vec) const
Check that every cell to be refined has cuboid shape.
Definition CpGridData.cpp:1778
bool disjointPatches(const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec) const
Determine if a finite amount of patches (of cells) are disjoint, namely, they do not share any corner...
Definition CpGridData.cpp:1925
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition CpGridData.hpp:605
void validStartEndIJKs(const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec) const
Check startIJK and endIJK of each patch of cells to be refined are valid, i.e.
Definition CpGridData.cpp:2076
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:199
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:123
Raw corner-point specification of a particular geological model.
Definition preprocess.h:56