My Project
CpGrid.hpp
1//===========================================================================
2//
3// File: CpGrid.hpp
4//
5// Created: Fri May 29 20:26:36 2009
6//
7// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8// B�rd Skaflestad <bard.skaflestad@sintef.no>
9//
10// $Date$
11//
12// $Revision$
13//
14//===========================================================================
15
16/*
17 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
18 Copyright 2009, 2010, 2014, 2022 Equinor ASA.
19 Copyright 2014, 2015 Dr. Blatt - HPC-Simulartion-Software & Services
20 Copyright 2015 NTNU
21
22 This file is part of The Open Porous Media project (OPM).
23
24 OPM is free software: you can redistribute it and/or modify
25 it under the terms of the GNU General Public License as published by
26 the Free Software Foundation, either version 3 of the License, or
27 (at your option) any later version.
28
29 OPM is distributed in the hope that it will be useful,
30 but WITHOUT ANY WARRANTY; without even the implied warranty of
31 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
32 GNU General Public License for more details.
33
34 You should have received a copy of the GNU General Public License
35 along with OPM. If not, see <http://www.gnu.org/licenses/>.
36*/
37
38#ifndef OPM_CPGRID_HEADER
39#define OPM_CPGRID_HEADER
40
41// Warning suppression for Dune includes.
42#include <opm/grid/utility/platform_dependent/disable_warnings.h> // Not really needed it seems, but alas.
43
44#include <dune/grid/common/grid.hh>
45#include <opm/grid/cpgrid/CpGridDataTraits.hpp>
46#include <opm/grid/cpgrid/OrientedEntityTable.hpp>
48#include <opm/grid/utility/platform_dependent/reenable_warnings.h> // Not really needed it seems, but alas.
49#include "common/GridEnums.hpp"
50#include <opm/grid/utility/OpmWellType.hpp>
51
52#include <iostream>
53#if ! HAVE_MPI
54#include <list>
55#endif
56
57namespace Opm
58{
59class EclipseGrid;
60class EclipseState;
61}
62
63namespace Dune
64{
65
66 class CpGrid;
67
68 namespace cpgrid
69 {
70 class CpGridData;
71 template <int> class Entity;
72 template<int,int> class Geometry;
73 class HierarchicIterator;
74 class IntersectionIterator;
75 template<int, PartitionIteratorType> class Iterator;
76 class LevelGlobalIdSet;
77 class GlobalIdSet;
78 class Intersection;
79 class IntersectionIterator;
80 class IndexSet;
81 class IdSet;
82
83
84 }
85}
86void refine_and_check(const Dune::cpgrid::Geometry<3, 3>&,
87 const std::array<int, 3>&,
88 bool);
89
90void refinePatch_and_check(Dune::CpGrid&,
91 const std::array<int,3>&,
92 const std::array<int,3>&,
93 const std::array<int,3>&);
94
95void refinePatch_and_check(const std::array<int,3>&,
96 const std::array<int,3>&,
97 const std::array<int,3>&);
98
99void check_global_refine(const Dune::CpGrid&,
100 const Dune::CpGrid&);
101namespace Dune
102{
103
105 //
106 // CpGridTraits
107 //
109
111 {
113 typedef CpGrid Grid;
114
123
126
129 template <int cd>
130 struct Codim
131 {
134 typedef cpgrid::Geometry<3-cd, 3> Geometry;
135 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> Geometry;
138 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> LocalGeometry;
141
144
147
150
153 template <PartitionIteratorType pitype>
155 {
160 };
161 };
162
165 template <PartitionIteratorType pitype>
167 {
169 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
171 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
172
173 };
174
176 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid>> LevelGridView;
178 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid>> LeafGridView;
179
188
190 using Communication = cpgrid::CpGridDataTraits::Communication;
191 using CollectiveCommunication = cpgrid::CpGridDataTraits::CollectiveCommunication;
192 };
193
195 //
196 // CpGridFamily
197 //
199
201 {
202 typedef CpGridTraits Traits;
203 };
204
206 //
207 // CpGrid
208 //
210
212 class CpGrid
213 : public GridDefaultImplementation<3, 3, double, CpGridFamily>
214 {
215 friend class cpgrid::CpGridData;
216 friend class cpgrid::Entity<0>;
217 friend class cpgrid::Entity<1>;
218 friend class cpgrid::Entity<2>;
219 friend class cpgrid::Entity<3>;
220 template<int dim>
221 friend cpgrid::Entity<dim> createEntity(const CpGrid&,int,bool);
222 friend
223 void ::refine_and_check(const Dune::cpgrid::Geometry<3,3>&,
224 const std::array<int,3>&,
225 bool);
226 friend
227 void ::refinePatch_and_check(Dune::CpGrid&,
228 const std::array<int,3>&,
229 const std::array<int,3>&,
230 const std::array<int,3>&);
231 friend
232 void ::refinePatch_and_check(const std::array<int,3>&,
233 const std::array<int,3>&,
234 const std::array<int,3>&);
235 friend
236 void ::check_global_refine(const Dune::CpGrid&,
237 const Dune::CpGrid&);
238
239 public:
240
241 // --- Typedefs ---
242
243
246
247
248 // --- Methods ---
249
250
252 CpGrid();
253
254 CpGrid(MPIHelper::MPICommunicator comm);
255
257
258
261 void readSintefLegacyFormat(const std::string& grid_prefix);
262
263
267 void writeSintefLegacyFormat(const std::string& grid_prefix) const;
268
269
270#if HAVE_ECL_INPUT
289 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
290 Opm::EclipseState* ecl_state,
291 bool periodic_extension, bool turn_normals, bool clip_z,
292 bool pinchActive);
293
313 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
314 Opm::EclipseState* ecl_state,
315 bool periodic_extension, bool turn_normals = false, bool clip_z = false);
316
317#endif
318
322 void processEclipseFormat(const grdecl& input_data, bool remove_ij_boundary, bool turn_normals = false);
323
325
331
332
339 void createCartesian(const std::array<int, 3>& dims,
340 const std::array<double, 3>& cellsize,
341 const std::array<int, 3>& shift = {0,0,0});
342
346 const std::array<int, 3>& logicalCartesianSize() const;
347
355 const std::vector<int>& globalCell() const;
356
364 void getIJK(const int c, std::array<int,3>& ijk) const;
366
370 bool uniqueBoundaryIds() const;
371
374 void setUniqueBoundaryIds(bool uids);
375
376
377 // --- Dune interface below ---
378
380 // \@{
385 std::string name() const;
386
388 int maxLevel() const;
389
391 template<int codim>
392 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
394 template<int codim>
395 typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
396
398 template<int codim, PartitionIteratorType PiType>
399 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
401 template<int codim, PartitionIteratorType PiType>
402 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
403
405 template<int codim>
406 typename Traits::template Codim<codim>::LeafIterator leafbegin() const;
408 template<int codim>
409 typename Traits::template Codim<codim>::LeafIterator leafend() const;
410
412 template<int codim, PartitionIteratorType PiType>
413 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const;
415 template<int codim, PartitionIteratorType PiType>
416 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const;
417
419 int size (int level, int codim) const;
420
422 int size (int codim) const;
423
425 int size (int level, GeometryType type) const;
426
428 int size (GeometryType type) const;
429
431 const Traits::GlobalIdSet& globalIdSet() const;
432
434 const Traits::LocalIdSet& localIdSet() const;
435
437 const Traits::LevelIndexSet& levelIndexSet(int level) const;
438
440 const Traits::LeafIndexSet& leafIndexSet() const;
441
443 void globalRefine (int);
444
445 const std::vector<Dune::GeometryType>& geomTypes(const int) const;
446
448 template <int codim>
450
463 void addLgrUpdateLeafView(const std::array<int,3>& cells_per_dim, const std::array<int,3>& startIJK, const std::array<int,3>& endIJK);
464
465 /* No refinement implemented. GridDefaultImplementation's methods will be used.
466
476
477 bool mark(int refCount, const typename Traits::template Codim<0>::EntityPointer & e)
478 {
479 return hostgrid_->mark(refCount, getHostEntity<0>(*e));
480 }
481
485
486 int getMark(const typename Traits::template Codim<0>::EntityPointer & e) const
487 {
488 return hostgrid_->getMark(getHostEntity<0>(*e));
489 }
490
492 bool preAdapt() {
493 return hostgrid_->preAdapt();
494 }
495
496
498 bool adapt()
499 {
500 return hostgrid_->adapt();
501 }
502
504 void postAdapt() {
505 return hostgrid_->postAdapt();
506 }
507
508 end of refinement section */
509
511 unsigned int overlapSize(int) const;
512
513
515 unsigned int ghostSize(int) const;
516
518 unsigned int overlapSize(int, int) const;
519
521 unsigned int ghostSize(int, int) const;
522
524 unsigned int numBoundarySegments() const;
525
526 void setZoltanParams(const std::map<std::string,std::string>& params);
527
528
529 // loadbalance is not part of the grid interface therefore we skip it.
530
536 bool loadBalance(int overlapLayers=1, bool useZoltan=true)
537 {
538 using std::get;
539 return get<0>(scatterGrid(defaultTransEdgeWgt, false, nullptr, false, nullptr, true, overlapLayers, useZoltan ));
540 }
541
542 // loadbalance is not part of the grid interface therefore we skip it.
543
563 std::pair<bool,std::vector<std::pair<std::string,bool>>>
564 loadBalance(const std::vector<cpgrid::OpmWellType> * wells,
565 const double* transmissibilities = nullptr,
566 int overlapLayers=1, bool useZoltan=true)
567 {
568 return scatterGrid(defaultTransEdgeWgt, false, wells, false, transmissibilities, false, overlapLayers, useZoltan);
569 }
570
571 // loadbalance is not part of the grid interface therefore we skip it.
572
596 std::pair<bool,std::vector<std::pair<std::string,bool>>>
597 loadBalance(EdgeWeightMethod method, const std::vector<cpgrid::OpmWellType> * wells,
598 const double* transmissibilities = nullptr, bool ownersFirst=false,
599 bool addCornerCells=false, int overlapLayers=1,
600 bool useZoltan = true)
601 {
602 return scatterGrid(method, ownersFirst, wells, false, transmissibilities, addCornerCells, overlapLayers, useZoltan);
603 }
604
624 template<class DataHandle>
625 std::pair<bool, std::vector<std::pair<std::string,bool> > >
626 loadBalance(DataHandle& data,
627 const std::vector<cpgrid::OpmWellType> * wells,
628 const double* transmissibilities = nullptr,
629 int overlapLayers=1, bool useZoltan = true)
630 {
631 auto ret = loadBalance(wells, transmissibilities, overlapLayers, useZoltan);
632 using std::get;
633 if (get<0>(ret))
634 {
635 scatterData(data);
636 }
637 return ret;
638 }
639
668 template<class DataHandle>
669 std::pair<bool, std::vector<std::pair<std::string,bool> > >
670 loadBalance(DataHandle& data, EdgeWeightMethod method,
671 const std::vector<cpgrid::OpmWellType> * wells,
672 bool serialPartitioning,
673 const double* transmissibilities = nullptr, bool ownersFirst=false,
674 bool addCornerCells=false, int overlapLayers=1, bool useZoltan = true,
675 double zoltanImbalanceTol = 1.1,
676 bool allowDistributedWells = false)
677 {
678 auto ret = scatterGrid(method, ownersFirst, wells, serialPartitioning, transmissibilities,
679 addCornerCells, overlapLayers, useZoltan, zoltanImbalanceTol, allowDistributedWells);
680 using std::get;
681 if (get<0>(ret))
682 {
683 scatterData(data);
684 }
685 return ret;
686 }
687
704 template<class DataHandle>
705 std::pair<bool, std::vector<std::pair<std::string,bool> > >
706 loadBalance(DataHandle& data, const std::vector<int>& parts,
707 const std::vector<cpgrid::OpmWellType> * wells,
708 bool ownersFirst=false,
709 bool addCornerCells=false, int overlapLayers=1)
710 {
711 using std::get;
712 auto ret = scatterGrid(defaultTransEdgeWgt, ownersFirst, wells,
713 /* serialPartitioning = */ false,
714 /* transmissibilities = */ {},
715 addCornerCells, overlapLayers, /* useZoltan =*/ false,
716 /* zoltanImbalanceTol (ignored) = */ 0.0,
717 /* allowDistributedWells = */ true, parts);
718 using std::get;
719 if (get<0>(ret))
720 {
721 scatterData(data);
722 }
723 return ret;
724 }
725
734 template<class DataHandle>
735 bool loadBalance(DataHandle& data,
736 decltype(data.fixedSize(0,0)) overlapLayers=1, bool useZoltan = true)
737 {
738 // decltype usage needed to tell the compiler not to use this function if first
739 // argument is std::vector but rather loadbalance by parts
740 bool ret = loadBalance(overlapLayers, useZoltan);
741 if (ret)
742 {
743 scatterData(data);
744 }
745 return ret;
746 }
747
759 bool loadBalance(const std::vector<int>& parts, bool ownersFirst=false,
760 bool addCornerCells=false, int overlapLayers=1)
761 {
762 using std::get;
763 return get<0>(scatterGrid(defaultTransEdgeWgt, ownersFirst, /* wells = */ {},
764 /* serialPartitioning = */ false,
765 /* trabsmissibilities = */ {},
766 addCornerCells, overlapLayers, /* useZoltan =*/ false,
767 /* zoltanImbalanceTol (ignored) = */ 0.0,
768 /* allowDistributedWells = */ true, parts));
769 }
770
783 template<class DataHandle>
784 bool loadBalance(DataHandle& data, const std::vector<int>& parts, bool ownersFirst=false,
785 bool addCornerCells=false, int overlapLayers=1)
786 {
787 bool ret = loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
788 if (ret)
789 {
790 scatterData(data);
791 }
792 return ret;
793 }
794
800 std::vector<int>
801 zoltanPartitionWithoutScatter(const std::vector<cpgrid::OpmWellType>* wells,
802 const double* transmissibilities,
803 const int numParts,
804 const double zoltanImbalanceTol) const;
805
813 template<class DataHandle>
814 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int /*level*/) const
815 {
816 communicate(data, iftype, dir);
817 }
818
826 template<class DataHandle>
827 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const;
828
830 const typename CpGridTraits::Communication& comm () const;
832
833 // ------------ End of Dune interface, start of simplified interface --------------
834
840
841 // enum { dimension = 3 }; // already defined
842
843 typedef Dune::FieldVector<double, 3> Vector;
844
845
846 const std::vector<double>& zcornData() const;
847
848
849 // Topology
851 int numCells() const;
852
854 int numFaces() const;
855
857 int numVertices() const;
858
859
866 int numCellFaces(int cell) const;
867
872 int cellFace(int cell, int local_index) const;
873
877
888 int faceCell(int face, int local_index) const;
889
896 int numCellFaces() const;
897
898 int numFaceVertices(int face) const;
899
904 int faceVertex(int face, int local_index) const;
905
908 double cellCenterDepth(int cell_index) const;
909
910
911 const Vector faceCenterEcl(int cell_index, int face) const;
912
913 const Vector faceAreaNormalEcl(int face) const;
914
915
916 // Geometry
920 const Vector& vertexPosition(int vertex) const;
921
924 double faceArea(int face) const;
925
928 const Vector& faceCentroid(int face) const;
929
933 const Vector& faceNormal(int face) const;
934
937 double cellVolume(int cell) const;
938
941 const Vector& cellCentroid(int cell) const;
942
945 template<int codim>
947 : public RandomAccessIteratorFacade<CentroidIterator<codim>,
948 FieldVector<double, 3>,
949 const FieldVector<double, 3>&, int>
950 {
951 public:
953 typedef typename std::vector<cpgrid::Geometry<3-codim, 3> >::const_iterator
958 : iter_(iter)
959 {}
960
961 const FieldVector<double,3>& dereference() const
962 {
963 return iter_->center();
964 }
965 void increment()
966 {
967 ++iter_;
968 }
969 const FieldVector<double,3>& elementAt(int n)
970 {
971 return iter_[n]->center();
972 }
973 void advance(int n){
974 iter_+=n;
975 }
976 void decrement()
977 {
978 --iter_;
979 }
980 int distanceTo(const CentroidIterator& o)
981 {
982 return o-iter_;
983 }
984 bool equals(const CentroidIterator& o) const{
985 return o==iter_;
986 }
987 private:
989 GeometryIterator iter_;
990 };
991
993 CentroidIterator<0> beginCellCentroids() const;
994
996 CentroidIterator<1> beginFaceCentroids() const;
997
998 // Extra
999 int boundaryId(int face) const;
1000
1007 template<class Cell2FacesRowIterator>
1008 int
1009 faceTag(const Cell2FacesRowIterator& cell_face) const;
1010
1012
1013 // ------------ End of simplified interface --------------
1014
1015 //------------- methods not in the DUNE grid interface.
1016
1021
1022
1031 template<class DataHandle>
1032 void scatterData(DataHandle& handle) const;
1033
1040 template<class DataHandle>
1041 void gatherData(DataHandle& handle) const;
1042
1044 using InterfaceMap = cpgrid::CpGridDataTraits::InterfaceMap;
1045
1075
1079
1081 void switchToGlobalView();
1082
1086
1087#if HAVE_MPI
1089 using ParallelIndexSet = cpgrid::CpGridDataTraits::ParallelIndexSet;
1091 using RemoteIndices = cpgrid::CpGridDataTraits::RemoteIndices;
1092
1094 using CommunicationType = cpgrid::CpGridDataTraits::CommunicationType;
1095
1099 const CommunicationType& cellCommunication() const;
1100
1101 ParallelIndexSet& getCellIndexSet();
1102
1103 RemoteIndices& getCellRemoteIndices();
1104
1105 const ParallelIndexSet& getCellIndexSet() const;
1106
1107 const RemoteIndices& getCellRemoteIndices() const;
1108#endif
1109
1111 const std::vector<int>& sortedNumAquiferCells() const;
1112
1113 private:
1139 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1140 scatterGrid(EdgeWeightMethod method,
1141 bool ownersFirst,
1142 const std::vector<cpgrid::OpmWellType> * wells,
1143 bool serialPartitioning,
1144 const double* transmissibilities,
1145 bool addCornerCells,
1146 int overlapLayers,
1147 bool useZoltan = true,
1148 double zoltanImbalanceTol = 1.1,
1149 bool allowDistributedWells = true,
1150 const std::vector<int>& input_cell_part = {});
1151
1156 std::vector<std::shared_ptr<cpgrid::CpGridData>> data_;
1158 cpgrid::CpGridData* current_view_data_;
1160 std::vector<std::shared_ptr<cpgrid::CpGridData>> distributed_data_;
1166 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1167 /*
1168 * @brief Interface for scattering and gathering point data.
1169 *
1170 * @warning Will only update owner cells
1171 */
1172 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1176 std::shared_ptr<cpgrid::GlobalIdSet> global_id_set_ptr_;
1177
1178
1182 std::map<std::string,std::string> zoltanParams;
1183
1184 }; // end Class CpGrid
1185
1186} // end namespace Dune
1187
1188#include <opm/grid/cpgrid/Entity.hpp>
1189#include <opm/grid/cpgrid/Iterators.hpp>
1190#include <opm/grid/cpgrid/CpGridData.hpp>
1191
1192
1193namespace Dune
1194{
1195
1196 namespace Capabilities
1197 {
1199 template <>
1200 struct hasEntity<CpGrid, 0>
1201 {
1202 static const bool v = true;
1203 };
1204
1206 template <>
1207 struct hasEntity<CpGrid, 3>
1208 {
1209 static const bool v = true;
1210 };
1211
1212 template<>
1213 struct canCommunicate<CpGrid,0>
1214 {
1215 static const bool v = true;
1216 };
1217
1218 template<>
1219 struct canCommunicate<CpGrid,3>
1220 {
1221 static const bool v = true;
1222 };
1223
1225 template <>
1226 struct hasBackupRestoreFacilities<CpGrid>
1227 {
1228 static const bool v = false;
1229 };
1230
1231 }
1232
1233 template<class DataHandle>
1234 void CpGrid::communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
1235 {
1236 current_view_data_->communicate(data, iftype, dir);
1237 }
1238
1239
1240 template<class DataHandle>
1241 void CpGrid::scatterData(DataHandle& handle) const
1242 {
1243#if HAVE_MPI
1244 if(distributed_data_.empty())
1245 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balanced grid!");
1246 distributed_data_[0]->scatterData(handle, data_[0].get(), distributed_data_[0].get(), cellScatterGatherInterface(),
1248#else
1249 // Suppress warnings for unused argument.
1250 (void) handle;
1251#endif
1252 }
1253
1254 template<class DataHandle>
1255 void CpGrid::gatherData(DataHandle& handle) const
1256 {
1257#if HAVE_MPI
1258 if(distributed_data_.empty())
1259 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balance grid!");
1260 distributed_data_[0]->gatherData(handle, data_[0].get(), distributed_data_[0].get());
1261#else
1262 // Suppress warnings for unused argument.
1263 (void) handle;
1264#endif
1265 }
1266
1267
1268 template<class Cell2FacesRowIterator>
1269 int
1270 CpGrid::faceTag(const Cell2FacesRowIterator& cell_face) const
1271 {
1272 // Note that this relies on the following implementation detail:
1273 // The grid is always constructed such that the interior faces constructed
1274 // with orientation set to true are
1275 // oriented along the positive IJK direction. Oriented means that
1276 // the first cell attached to face has the lower index.
1277 // For faces along the boundary (only one cell, always attached at index 0)
1278 // the orientation has to be determined by the orientation of the cell.
1279 // If it is true then in UnstructuredGrid it would be stored at index 0,
1280 // otherwise at index 1.
1281 const int cell = cell_face.getCellIndex();
1282 const int face = *cell_face;
1283 assert (0 <= cell); assert (cell < numCells());
1284 assert (0 <= face); assert (face < numFaces());
1285
1287
1288 const cpgrid::EntityRep<1> f(face, true);
1289 const F2C& f2c = current_view_data_->face_to_cell_[f];
1290 const face_tag tag = current_view_data_->face_tag_[f];
1291
1292 assert ((f2c.size() == 1) || (f2c.size() == 2));
1293
1294 int inside_cell = 0;
1295
1296 if ( f2c.size() == 2 ) // Two cells => interior
1297 {
1298 if ( f2c[1].index() == cell )
1299 {
1300 inside_cell = 1;
1301 }
1302 }
1303 const bool normal_is_in = ! f2c[inside_cell].orientation();
1304
1305 switch (tag) {
1306 case I_FACE:
1307 // LEFT : RIGHT
1308 return normal_is_in ? 0 : 1; // min(I) : max(I)
1309 case J_FACE:
1310 // BACK : FRONT
1311 return normal_is_in ? 2 : 3; // min(J) : max(J)
1312 case K_FACE:
1313 // Note: TOP at min(K) as 'z' measures *depth*.
1314 // TOP : BOTTOM
1315 return normal_is_in ? 4 : 5; // min(K) : max(K)
1316 case NNC_FACE:
1317 // For nnc faces we return the otherwise unused value -1.
1318 return -1;
1319 default:
1320 OPM_THROW(std::logic_error, "Unhandled face tag. This should never happen!");
1321 }
1322 }
1323
1324 template<int dim>
1325 cpgrid::Entity<dim> createEntity(const CpGrid&, int, bool);
1326
1327} // namespace Dune
1328
1329#include <opm/grid/cpgrid/PersistentContainer.hpp>
1330#include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
1331#include "cpgrid/Intersection.hpp"
1332#include "cpgrid/Geometry.hpp"
1333#include "cpgrid/Indexsets.hpp"
1334
1335#endif // OPM_CPGRID_HEADER
An iterator over the centroids of the geometry of the entities.
Definition: CpGrid.hpp:950
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition: CpGrid.hpp:957
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition: CpGrid.hpp:954
[ provides Dune::Grid ]
Definition: CpGrid.hpp:214
std::string name() const
Get the grid name.
Definition: CpGrid.cpp:584
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:597
void switchToGlobalView()
Switch to the global view.
Definition: CpGrid.cpp:1251
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:1293
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: CpGrid.cpp:923
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition: CpGrid.hpp:245
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
void gatherData(DataHandle &handle) const
Moves data from the distributed view to the global (all data on process) view.
Definition: CpGrid.hpp:1255
void globalRefine(int)
global refinement
Definition: CpGrid.cpp:881
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: CpGrid.cpp:864
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition: CpGrid.cpp:1185
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: CpGrid.cpp:859
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition: CpGrid.hpp:1270
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGrid.cpp:579
void createCartesian(const std::array< int, 3 > &dims, const std::array< double, 3 > &cellsize, const std::array< int, 3 > &shift={0, 0, 0})
Create a cartesian grid.
Definition: CpGrid.cpp:496
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
void processEclipseFormat(const grdecl &input_data, bool remove_ij_boundary, bool turn_normals=false)
Read the Eclipse grid format ('grdecl').
Definition: CpGrid.cpp:1338
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition: CpGrid.cpp:901
cpgrid::CpGridDataTraits::InterfaceMap InterfaceMap
The type of the map describing communication interfaces.
Definition: CpGrid.hpp:1044
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
const Vector & faceNormal(int face) const
Get the unit normal of a face.
Definition: CpGrid.cpp:1170
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition: CpGrid.cpp:1020
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition: CpGrid.cpp:564
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, const std::vector< int > &parts, const std::vector< cpgrid::OpmWellType > *wells, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:706
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition: CpGrid.cpp:1180
int maxLevel() const
Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview),...
Definition: CpGrid.cpp:589
const CpGridTraits::Communication & comm() const
Get the collective communication object.
Definition: CpGrid.cpp:948
double faceArea(int face) const
Get the area of a face.
Definition: CpGrid.cpp:1160
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition: CpGrid.cpp:1165
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: CpGrid.cpp:869
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition: CpGrid.cpp:835
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition: CpGrid.cpp:559
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition: CpGrid.cpp:1030
bool loadBalance(int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:536
int numVertices() const
Get The number of vertices.
Definition: CpGrid.cpp:969
void getIJK(const int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGrid.cpp:569
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition: CpGrid.cpp:979
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGrid.cpp:574
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:564
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lend(int level) const
one past the end on this level
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, bool serialPartitioning, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, bool useZoltan=true, double zoltanImbalanceTol=1.1, bool allowDistributedWells=false)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:670
std::vector< int > zoltanPartitionWithoutScatter(const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities, const int numParts, const double zoltanImbalanceTol) const
Partitions the grid using Zoltan without decomposing and distributing it among processes.
Definition: CpGrid.cpp:149
bool loadBalance(DataHandle &data, decltype(data.fixedSize(0, 0)) overlapLayers=1, bool useZoltan=true)
Distributes this grid and data over the available nodes in a distributed machine.
Definition: CpGrid.hpp:735
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.cpp:989
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim and PartitionIteratorType.
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:1303
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition: CpGrid.cpp:907
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, int overlapLayers=1, bool useZoltan=true)
Distributes this grid and data over the available nodes in a distributed machine.
Definition: CpGrid.hpp:626
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: CpGrid.cpp:876
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
const InterfaceMap & pointScatterGatherInterface() const
Get an interface for gathering/scattering data attached to points with communication.
Definition: CpGrid.cpp:1246
CpGrid()
Default constructor.
Definition: CpGrid.cpp:122
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level and PartitionIteratorType.
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition: CpGrid.cpp:1195
void switchToDistributedView()
Switch to the distributed view.
Definition: CpGrid.cpp:1256
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition: CpGrid.cpp:1190
bool loadBalance(DataHandle &data, const std::vector< int > &parts, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
Distributes this grid and data over the available nodes in a distributed machine.
Definition: CpGrid.hpp:784
const cpgrid::OrientedEntityTable< 0, 1 >::row_type cellFaceRow(int cell) const
Get a list of indices identifying all faces of a cell.
Definition: CpGrid.cpp:984
int numCells() const
Get the number of cells.
Definition: CpGrid.cpp:959
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition: CpGrid.cpp:1241
void addLgrUpdateLeafView(const std::array< int, 3 > &cells_per_dim, const std::array< int, 3 > &startIJK, const std::array< int, 3 > &endIJK)
Create a grid out of a coarse one and a refinement(LGR) of a selected block-shaped patch of cells fro...
Definition: CpGrid.cpp:1365
cpgrid::Entity< codim > entity(const cpgrid::Entity< codim > &) const
given an EntitySeed (or EntityPointer) return an entity object
Definition: CpGrid.cpp:892
bool loadBalance(const std::vector< int > &parts, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:759
int numFaces() const
Get the number of faces.
Definition: CpGrid.cpp:964
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition: CpGrid.cpp:1155
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition: CpGrid.cpp:1035
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition: CpGrid.hpp:814
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
double cellVolume(int cell) const
Get the volume of the cell.
Definition: CpGrid.cpp:1175
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition: CpGrid.hpp:1241
Struct that hods all the data needed to represent a Cpgrid.
Definition: CpGridData.hpp:122
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition: CpGridData.hpp:738
Represents an entity of a given codim, with positive or negative orientation.
Definition: EntityRep.hpp:99
This class encapsulates geometry for vertices, intersections, and cells.
Definition: Geometry.hpp:74
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
Definition: Intersection.hpp:288
Definition: Intersection.hpp:66
Iterator intended to be used as LeafIterator and LevelIterator (no difference due to no adaptivity) f...
Definition: Iterators.hpp:60
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
EdgeWeightMethod
enum for choosing Methods for weighting graph-edges correspoding to cell interfaces in Zoltan's graph...
Definition: GridEnums.hpp:34
@ defaultTransEdgeWgt
Use the transmissibilities as edge weights.
Definition: GridEnums.hpp:38
Holds the implementation of the CpGrid as a pimple.
Definition: CellQuadrature.hpp:29
Low-level corner-point processing routines and supporting data structures.
face_tag
Connection taxonomy.
Definition: preprocess.h:66
@ K_FACE
Connection topologically normal to I-J plane.
Definition: preprocess.h:69
@ J_FACE
Connection topologically normal to I-K plane.
Definition: preprocess.h:68
@ NNC_FACE
Arbitrary non-neighbouring connection.
Definition: preprocess.h:70
@ I_FACE
Connection topologically normal to J-K plane.
Definition: preprocess.h:67
Definition: CpGrid.hpp:201
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:155
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition: CpGrid.hpp:157
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition: CpGrid.hpp:159
Traits associated with a specific codim.
Definition: CpGrid.hpp:131
cpgrid::Entity< cd > Entity
The type of the entity.
Definition: CpGrid.hpp:140
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition: CpGrid.hpp:134
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition: CpGrid.hpp:137
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition: CpGrid.hpp:146
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition: CpGrid.hpp:143
cpgrid::Entity< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:149
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:167
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:171
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:169
Definition: CpGrid.hpp:111
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition: CpGrid.hpp:181
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition: CpGrid.hpp:120
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:178
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:176
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition: CpGrid.hpp:185
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition: CpGrid.hpp:122
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition: CpGrid.hpp:183
GlobalIdSet LocalIdSet
The type of the local id set.
Definition: CpGrid.hpp:187
cpgrid::CpGridDataTraits::Communication Communication
The type of the collective communication.
Definition: CpGrid.hpp:190
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition: CpGrid.hpp:125
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition: CpGrid.hpp:118
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition: CpGrid.hpp:116
CpGrid Grid
The type that implements the grid.
Definition: CpGrid.hpp:113
Raw corner-point specification of a particular geological model.
Definition: preprocess.h:56