My Project
Loading...
Searching...
No Matches
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// Antonella Ritorto <antonella.ritorto@opm-op.com>
10//
11// $Date$
12//
13// $Revision$
14//
15//===========================================================================
16
17/*
18 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
19 Copyright 2009, 2010, 2014, 2022-2023 Equinor ASA.
20 Copyright 2014, 2015 Dr. Blatt - HPC-Simulartion-Software & Services
21 Copyright 2015 NTNU
22
23 This file is part of The Open Porous Media project (OPM).
24
25 OPM is free software: you can redistribute it and/or modify
26 it under the terms of the GNU General Public License as published by
27 the Free Software Foundation, either version 3 of the License, or
28 (at your option) any later version.
29
30 OPM is distributed in the hope that it will be useful,
31 but WITHOUT ANY WARRANTY; without even the implied warranty of
32 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
33 GNU General Public License for more details.
34
35 You should have received a copy of the GNU General Public License
36 along with OPM. If not, see <http://www.gnu.org/licenses/>.
37*/
38
39#ifndef OPM_CPGRID_HEADER
40#define OPM_CPGRID_HEADER
41
42// Warning suppression for Dune includes.
43#include <opm/grid/utility/platform_dependent/disable_warnings.h> // Not really needed it seems, but alas.
44
45#include <dune/grid/common/grid.hh>
46#include <opm/grid/cpgrid/CpGridDataTraits.hpp>
47#include <opm/grid/cpgrid/OrientedEntityTable.hpp>
49#include <opm/grid/utility/platform_dependent/reenable_warnings.h> // Not really needed it seems, but alas.
50#include "common/GridEnums.hpp"
51#include <opm/grid/utility/OpmWellType.hpp>
52
53#include <iostream>
54#if ! HAVE_MPI
55#include <list>
56#endif
57
58namespace Opm
59{
60class EclipseGrid;
61class EclipseState;
62template<typename Grid, typename GridView> class LookUpData;
63template<typename Grid, typename GridView> class LookUpCartesianData;
64}
65
66namespace Dune
67{
68
69 class CpGrid;
70 template<typename Grid, typename GridView> class LookUpData;
71 template<typename Grid, typename GridView> class LookUpCartesianData;
72
73 namespace cpgrid
74 {
75 class CpGridData;
76 template <int> class Entity;
77 template<int,int> class Geometry;
78 class HierarchicIterator;
79 class IntersectionIterator;
80 template<int, PartitionIteratorType> class Iterator;
81 class LevelGlobalIdSet;
82 class GlobalIdSet;
83 class Intersection;
84 class IntersectionIterator;
85 class IndexSet;
86 class IdSet;
87
88 }
89}
90
91void disjointPatches_check(Dune::CpGrid&,
92 const std::vector<std::array<int,3>>&,
93 const std::vector<std::array<int,3>>&);
94
95void lookup_check(const Dune::CpGrid&);
96
97void refine_and_check(const Dune::cpgrid::Geometry<3, 3>&,
98 const std::array<int, 3>&,
99 bool);
100
101void refinePatch_and_check(Dune::CpGrid&,
102 const std::vector<std::array<int,3>>&,
103 const std::vector<std::array<int,3>>&,
104 const std::vector<std::array<int,3>>&,
105 const std::vector<std::string>&);
106
107void refinePatch_and_check(const std::array<int,3>&,
108 const std::array<int,3>&,
109 const std::array<int,3>&);
110
111void check_global_refine(const Dune::CpGrid&,
112 const Dune::CpGrid&);
113namespace Dune
114{
115
117 //
118 // CpGridTraits
119 //
121
123 {
125 typedef CpGrid Grid;
126
135
138
141 template <int cd>
142 struct Codim
143 {
146 typedef cpgrid::Geometry<3-cd, 3> Geometry;
147 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> Geometry;
150 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> LocalGeometry;
153
156
159
162
165 template <PartitionIteratorType pitype>
173 };
174
177 template <PartitionIteratorType pitype>
179 {
181 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
183 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
184
185 };
186
188 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid>> LevelGridView;
190 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid>> LeafGridView;
191
200
202 using Communication = cpgrid::CpGridDataTraits::Communication;
203 using CollectiveCommunication = cpgrid::CpGridDataTraits::CollectiveCommunication;
204 };
205
207 //
208 // CpGridFamily
209 //
211
213 {
214 typedef CpGridTraits Traits;
215 };
216
218 //
219 // CpGrid
220 //
222
224 class CpGrid
225 : public GridDefaultImplementation<3, 3, double, CpGridFamily>
226 {
227 friend class cpgrid::CpGridData;
228 friend class cpgrid::Entity<0>;
229 friend class cpgrid::Entity<1>;
230 friend class cpgrid::Entity<2>;
231 friend class cpgrid::Entity<3>;
232 template<typename Grid, typename GridView> friend class Opm::LookUpData;
233 template<typename Grid, typename GridView> friend class Opm::LookUpCartesianData;
234 template<int dim>
235 friend cpgrid::Entity<dim> createEntity(const CpGrid&,int,bool);
236 friend void ::disjointPatches_check(Dune::CpGrid&,
237 const std::vector<std::array<int,3>>&,
238 const std::vector<std::array<int,3>>&);
239 friend void ::lookup_check(const Dune::CpGrid&);
240 friend
241 void ::refine_and_check(const Dune::cpgrid::Geometry<3,3>&,
242 const std::array<int,3>&,
243 bool);
244 friend
245 void ::refinePatch_and_check(Dune::CpGrid&,
246 const std::vector<std::array<int,3>>&,
247 const std::vector<std::array<int,3>>&,
248 const std::vector<std::array<int,3>>&,
249 const std::vector<std::string>&);
250 friend
251 void ::refinePatch_and_check(const std::array<int,3>&,
252 const std::array<int,3>&,
253 const std::array<int,3>&);
254 friend
255 void ::check_global_refine(const Dune::CpGrid&,
256 const Dune::CpGrid&);
257
258 public:
259
260 // --- Typedefs ---
261
262
265
266
267 // --- Methods ---
268
269
271 CpGrid();
272
273 CpGrid(MPIHelper::MPICommunicator comm);
274
276
277
280 void readSintefLegacyFormat(const std::string& grid_prefix);
281
282
286 void writeSintefLegacyFormat(const std::string& grid_prefix) const;
287
288
289#if HAVE_ECL_INPUT
308 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
309 Opm::EclipseState* ecl_state,
310 bool periodic_extension, bool turn_normals, bool clip_z,
311 bool pinchActive);
312
332 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
333 Opm::EclipseState* ecl_state,
334 bool periodic_extension, bool turn_normals = false, bool clip_z = false);
335
336#endif
337
341 void processEclipseFormat(const grdecl& input_data, bool remove_ij_boundary, bool turn_normals = false);
342
344
350
351
358 void createCartesian(const std::array<int, 3>& dims,
359 const std::array<double, 3>& cellsize,
360 const std::array<int, 3>& shift = {0,0,0});
361
365 const std::array<int, 3>& logicalCartesianSize() const;
366
374 const std::vector<int>& globalCell() const;
375
383 void getIJK(const int c, std::array<int,3>& ijk) const;
385
389 bool uniqueBoundaryIds() const;
390
393 void setUniqueBoundaryIds(bool uids);
394
395
396 // --- Dune interface below ---
397
399 // \@{
404 std::string name() const;
405
407 int maxLevel() const;
408
410 template<int codim>
411 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
413 template<int codim>
414 typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
415
417 template<int codim, PartitionIteratorType PiType>
418 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
420 template<int codim, PartitionIteratorType PiType>
421 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
422
424 template<int codim>
425 typename Traits::template Codim<codim>::LeafIterator leafbegin() const;
427 template<int codim>
428 typename Traits::template Codim<codim>::LeafIterator leafend() const;
429
431 template<int codim, PartitionIteratorType PiType>
432 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const;
434 template<int codim, PartitionIteratorType PiType>
435 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const;
436
438 int size (int level, int codim) const;
439
441 int size (int codim) const;
442
444 int size (int level, GeometryType type) const;
445
447 int size (GeometryType type) const;
448
450 const Traits::GlobalIdSet& globalIdSet() const;
451
453 const Traits::LocalIdSet& localIdSet() const;
454
456 const Traits::LevelIndexSet& levelIndexSet(int level) const;
457
459 const Traits::LeafIndexSet& leafIndexSet() const;
460
462 void globalRefine (int);
463
464 const std::vector<Dune::GeometryType>& geomTypes(const int) const;
465
467 template <int codim>
469
483 void addLgrUpdateLeafView(const std::array<int,3>& cells_per_dim, const std::array<int,3>& startIJK,
484 const std::array<int,3>& endIJK, const std::string& lgr_name);
485
502 void addLgrsUpdateLeafView(const std::vector<std::array<int,3>>& cells_per_dim_vec,
503 const std::vector<std::array<int,3>>& startIJK_vec,
504 const std::vector<std::array<int,3>>& endIJK_vec,
505 const std::vector<std::string>& lgr_name_vec);
506
507 // @brief TO BE DONE
508 const std::map<std::string,int>& getLgrNameToLevel() const;
509
510 // @breif Compute center of an entity/element/cell in the Eclipse way:
511 // - Average of the 4 corners of the bottom face.
512 // - Average of the 4 corners of the top face.
513 // Return average of the previous computations.
514 // @param [in] int Index of a cell.
515 // @return 'eclipse centroid'
516 std::array<double,3> getEclCentroid(const int& idx) const;
517
518 // @breif Compute center of an entity/element/cell in the Eclipse way:
519 // - Average of the 4 corners of the bottom face.
520 // - Average of the 4 corners of the top face.
521 // Return average of the previous computations.
522 // @param [in] Entity<0> Entity
523 // @return 'eclipse centroid'
524 std::array<double,3> getEclCentroid(const cpgrid::Entity<0>& elem) const;
525
526
527
528 /* No refinement implemented. GridDefaultImplementation's methods will be used.
529
539
540 bool mark(int refCount, const typename Traits::template Codim<0>::EntityPointer & e)
541 {
542 return hostgrid_->mark(refCount, getHostEntity<0>(*e));
543 }
544
548
549 int getMark(const typename Traits::template Codim<0>::EntityPointer & e) const
550 {
551 return hostgrid_->getMark(getHostEntity<0>(*e));
552 }
553
555 bool preAdapt() {
556 return hostgrid_->preAdapt();
557 }
558
559
561 bool adapt()
562 {
563 return hostgrid_->adapt();
564 }
565
567 void postAdapt() {
568 return hostgrid_->postAdapt();
569 }
570
571 end of refinement section */
572
574 unsigned int overlapSize(int) const;
575
576
578 unsigned int ghostSize(int) const;
579
581 unsigned int overlapSize(int, int) const;
582
584 unsigned int ghostSize(int, int) const;
585
587 unsigned int numBoundarySegments() const;
588
589 void setZoltanParams(const std::map<std::string,std::string>& params);
590
591
592 // loadbalance is not part of the grid interface therefore we skip it.
593
599 bool loadBalance(int overlapLayers=1, bool useZoltan=true)
600 {
601 using std::get;
602 return get<0>(scatterGrid(defaultTransEdgeWgt, false, nullptr, false, nullptr, true, overlapLayers, useZoltan ));
603 }
604
605 // loadbalance is not part of the grid interface therefore we skip it.
606
626 std::pair<bool,std::vector<std::pair<std::string,bool>>>
627 loadBalance(const std::vector<cpgrid::OpmWellType> * wells,
628 const double* transmissibilities = nullptr,
629 int overlapLayers=1, bool useZoltan=true)
630 {
631 return scatterGrid(defaultTransEdgeWgt, false, wells, false, transmissibilities, false, overlapLayers, useZoltan);
632 }
633
634 // loadbalance is not part of the grid interface therefore we skip it.
635
659 std::pair<bool,std::vector<std::pair<std::string,bool>>>
660 loadBalance(EdgeWeightMethod method, const std::vector<cpgrid::OpmWellType> * wells,
661 const double* transmissibilities = nullptr, bool ownersFirst=false,
662 bool addCornerCells=false, int overlapLayers=1,
663 bool useZoltan = true)
664 {
665 return scatterGrid(method, ownersFirst, wells, false, transmissibilities, addCornerCells, overlapLayers, useZoltan);
666 }
667
687 template<class DataHandle>
688 std::pair<bool, std::vector<std::pair<std::string,bool> > >
689 loadBalance(DataHandle& data,
690 const std::vector<cpgrid::OpmWellType> * wells,
691 const double* transmissibilities = nullptr,
692 int overlapLayers=1, bool useZoltan = true)
693 {
694 auto ret = loadBalance(wells, transmissibilities, overlapLayers, useZoltan);
695 using std::get;
696 if (get<0>(ret))
697 {
698 scatterData(data);
699 }
700 return ret;
701 }
702
731 template<class DataHandle>
732 std::pair<bool, std::vector<std::pair<std::string,bool> > >
733 loadBalance(DataHandle& data, EdgeWeightMethod method,
734 const std::vector<cpgrid::OpmWellType> * wells,
735 bool serialPartitioning,
736 const double* transmissibilities = nullptr, bool ownersFirst=false,
737 bool addCornerCells=false, int overlapLayers=1, bool useZoltan = true,
738 double zoltanImbalanceTol = 1.1,
739 bool allowDistributedWells = false)
740 {
741 auto ret = scatterGrid(method, ownersFirst, wells, serialPartitioning, transmissibilities,
742 addCornerCells, overlapLayers, useZoltan, zoltanImbalanceTol, allowDistributedWells);
743 using std::get;
744 if (get<0>(ret))
745 {
746 scatterData(data);
747 }
748 return ret;
749 }
750
767 template<class DataHandle>
768 std::pair<bool, std::vector<std::pair<std::string,bool> > >
769 loadBalance(DataHandle& data, const std::vector<int>& parts,
770 const std::vector<cpgrid::OpmWellType> * wells,
771 bool ownersFirst=false,
772 bool addCornerCells=false, int overlapLayers=1)
773 {
774 using std::get;
775 auto ret = scatterGrid(defaultTransEdgeWgt, ownersFirst, wells,
776 /* serialPartitioning = */ false,
777 /* transmissibilities = */ {},
778 addCornerCells, overlapLayers, /* useZoltan =*/ false,
779 /* zoltanImbalanceTol (ignored) = */ 0.0,
780 /* allowDistributedWells = */ true, parts);
781 using std::get;
782 if (get<0>(ret))
783 {
784 scatterData(data);
785 }
786 return ret;
787 }
788
797 template<class DataHandle>
798 bool loadBalance(DataHandle& data,
799 decltype(data.fixedSize(0,0)) overlapLayers=1, bool useZoltan = true)
800 {
801 // decltype usage needed to tell the compiler not to use this function if first
802 // argument is std::vector but rather loadbalance by parts
803 bool ret = loadBalance(overlapLayers, useZoltan);
804 if (ret)
805 {
806 scatterData(data);
807 }
808 return ret;
809 }
810
822 bool loadBalance(const std::vector<int>& parts, bool ownersFirst=false,
823 bool addCornerCells=false, int overlapLayers=1)
824 {
825 using std::get;
826 return get<0>(scatterGrid(defaultTransEdgeWgt, ownersFirst, /* wells = */ {},
827 /* serialPartitioning = */ false,
828 /* trabsmissibilities = */ {},
829 addCornerCells, overlapLayers, /* useZoltan =*/ false,
830 /* zoltanImbalanceTol (ignored) = */ 0.0,
831 /* allowDistributedWells = */ true, parts));
832 }
833
846 template<class DataHandle>
847 bool loadBalance(DataHandle& data, const std::vector<int>& parts, bool ownersFirst=false,
848 bool addCornerCells=false, int overlapLayers=1)
849 {
850 bool ret = loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
851 if (ret)
852 {
853 scatterData(data);
854 }
855 return ret;
856 }
857
863 std::vector<int>
864 zoltanPartitionWithoutScatter(const std::vector<cpgrid::OpmWellType>* wells,
865 const double* transmissibilities,
866 const int numParts,
867 const double zoltanImbalanceTol) const;
868
876 template<class DataHandle>
877 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int /*level*/) const
878 {
879 communicate(data, iftype, dir);
880 }
881
889 template<class DataHandle>
890 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const;
891
893 const typename CpGridTraits::Communication& comm () const;
895
896 // ------------ End of Dune interface, start of simplified interface --------------
897
903
904 // enum { dimension = 3 }; // already defined
905
906 typedef Dune::FieldVector<double, 3> Vector;
907
908
909 const std::vector<double>& zcornData() const;
910
911
912 // Topology
914 int numCells() const;
915
917 int numFaces() const;
918
920 int numVertices() const;
921
922
929 int numCellFaces(int cell) const;
930
935 int cellFace(int cell, int local_index) const;
936
940
951 int faceCell(int face, int local_index) const;
952
959 int numCellFaces() const;
960
961 int numFaceVertices(int face) const;
962
967 int faceVertex(int face, int local_index) const;
968
971 double cellCenterDepth(int cell_index) const;
972
973
974 const Vector faceCenterEcl(int cell_index, int face) const;
975
976 const Vector faceAreaNormalEcl(int face) const;
977
978
979 // Geometry
983 const Vector& vertexPosition(int vertex) const;
984
987 double faceArea(int face) const;
988
991 const Vector& faceCentroid(int face) const;
992
996 const Vector& faceNormal(int face) const;
997
1000 double cellVolume(int cell) const;
1001
1004 const Vector& cellCentroid(int cell) const;
1005
1008 template<int codim>
1010 : public RandomAccessIteratorFacade<CentroidIterator<codim>,
1011 FieldVector<double, 3>,
1012 const FieldVector<double, 3>&, int>
1013 {
1014 public:
1016 typedef typename std::vector<cpgrid::Geometry<3-codim, 3> >::const_iterator
1021 : iter_(iter)
1022 {}
1023
1024 const FieldVector<double,3>& dereference() const
1025 {
1026 return iter_->center();
1027 }
1028 void increment()
1029 {
1030 ++iter_;
1031 }
1032 const FieldVector<double,3>& elementAt(int n)
1033 {
1034 return iter_[n]->center();
1035 }
1036 void advance(int n){
1037 iter_+=n;
1038 }
1039 void decrement()
1040 {
1041 --iter_;
1042 }
1043 int distanceTo(const CentroidIterator& o)
1044 {
1045 return o-iter_;
1046 }
1047 bool equals(const CentroidIterator& o) const{
1048 return o==iter_;
1049 }
1050 private:
1052 GeometryIterator iter_;
1053 };
1054
1056 CentroidIterator<0> beginCellCentroids() const;
1057
1059 CentroidIterator<1> beginFaceCentroids() const;
1060
1061 // Extra
1062 int boundaryId(int face) const;
1063
1070 template<class Cell2FacesRowIterator>
1071 int
1072 faceTag(const Cell2FacesRowIterator& cell_face) const;
1073
1075
1076 // ------------ End of simplified interface --------------
1077
1078 //------------- methods not in the DUNE grid interface.
1079
1084
1085
1094 template<class DataHandle>
1095 void scatterData(DataHandle& handle) const;
1096
1103 template<class DataHandle>
1104 void gatherData(DataHandle& handle) const;
1105
1107 using InterfaceMap = cpgrid::CpGridDataTraits::InterfaceMap;
1108
1138
1142
1144 void switchToGlobalView();
1145
1149
1150#if HAVE_MPI
1152 using ParallelIndexSet = cpgrid::CpGridDataTraits::ParallelIndexSet;
1154 using RemoteIndices = cpgrid::CpGridDataTraits::RemoteIndices;
1155
1157 using CommunicationType = cpgrid::CpGridDataTraits::CommunicationType;
1158
1162 const CommunicationType& cellCommunication() const;
1163
1164 ParallelIndexSet& getCellIndexSet();
1165
1166 RemoteIndices& getCellRemoteIndices();
1167
1168 const ParallelIndexSet& getCellIndexSet() const;
1169
1170 const RemoteIndices& getCellRemoteIndices() const;
1171#endif
1172
1174 const std::vector<int>& sortedNumAquiferCells() const;
1175
1176 private:
1202 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1203 scatterGrid(EdgeWeightMethod method,
1204 bool ownersFirst,
1205 const std::vector<cpgrid::OpmWellType> * wells,
1206 bool serialPartitioning,
1207 const double* transmissibilities,
1208 bool addCornerCells,
1209 int overlapLayers,
1210 bool useZoltan = true,
1211 double zoltanImbalanceTol = 1.1,
1212 bool allowDistributedWells = true,
1213 const std::vector<int>& input_cell_part = {});
1214
1219 std::vector<std::shared_ptr<cpgrid::CpGridData>> data_;
1221 cpgrid::CpGridData* current_view_data_;
1223 std::vector<std::shared_ptr<cpgrid::CpGridData>> distributed_data_;
1225 std::map<std::string,int> lgr_names_ = {{"GLOBAL", 0}};
1231 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1232 /*
1233 * @brief Interface for scattering and gathering point data.
1234 *
1235 * @warning Will only update owner cells
1236 */
1237 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1241 std::shared_ptr<cpgrid::GlobalIdSet> global_id_set_ptr_;
1242
1243
1247 std::map<std::string,std::string> zoltanParams;
1248
1249 }; // end Class CpGrid
1250
1251} // end namespace Dune
1252
1253#include <opm/grid/cpgrid/Entity.hpp>
1254#include <opm/grid/cpgrid/Iterators.hpp>
1255#include <opm/grid/cpgrid/CpGridData.hpp>
1256
1257
1258namespace Dune
1259{
1260
1261 namespace Capabilities
1262 {
1264 template <>
1265 struct hasEntity<CpGrid, 0>
1266 {
1267 static const bool v = true;
1268 };
1269
1271 template <>
1272 struct hasEntity<CpGrid, 3>
1273 {
1274 static const bool v = true;
1275 };
1276
1277 template<>
1278 struct canCommunicate<CpGrid,0>
1279 {
1280 static const bool v = true;
1281 };
1282
1283 template<>
1284 struct canCommunicate<CpGrid,3>
1285 {
1286 static const bool v = true;
1287 };
1288
1290 template <>
1291 struct hasBackupRestoreFacilities<CpGrid>
1292 {
1293 static const bool v = false;
1294 };
1295
1296 }
1297
1298 template<class DataHandle>
1299 void CpGrid::communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
1300 {
1301 current_view_data_->communicate(data, iftype, dir);
1302 }
1303
1304
1305 template<class DataHandle>
1306 void CpGrid::scatterData(DataHandle& handle) const
1307 {
1308#if HAVE_MPI
1309 if(distributed_data_.empty())
1310 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balanced grid!");
1311 distributed_data_[0]->scatterData(handle, data_[0].get(), distributed_data_[0].get(), cellScatterGatherInterface(),
1313#else
1314 // Suppress warnings for unused argument.
1315 (void) handle;
1316#endif
1317 }
1318
1319 template<class DataHandle>
1320 void CpGrid::gatherData(DataHandle& handle) const
1321 {
1322#if HAVE_MPI
1323 if(distributed_data_.empty())
1324 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balance grid!");
1325 distributed_data_[0]->gatherData(handle, data_[0].get(), distributed_data_[0].get());
1326#else
1327 // Suppress warnings for unused argument.
1328 (void) handle;
1329#endif
1330 }
1331
1332
1333 template<class Cell2FacesRowIterator>
1334 int
1335 CpGrid::faceTag(const Cell2FacesRowIterator& cell_face) const
1336 {
1337 // Note that this relies on the following implementation detail:
1338 // The grid is always constructed such that the interior faces constructed
1339 // with orientation set to true are
1340 // oriented along the positive IJK direction. Oriented means that
1341 // the first cell attached to face has the lower index.
1342 // For faces along the boundary (only one cell, always attached at index 0)
1343 // the orientation has to be determined by the orientation of the cell.
1344 // If it is true then in UnstructuredGrid it would be stored at index 0,
1345 // otherwise at index 1.
1346 const int cell = cell_face.getCellIndex();
1347 const int face = *cell_face;
1348 assert (0 <= cell); assert (cell < numCells());
1349 assert (0 <= face); assert (face < numFaces());
1350
1352
1353 const cpgrid::EntityRep<1> f(face, true);
1354 const F2C& f2c = current_view_data_->face_to_cell_[f];
1355 const face_tag tag = current_view_data_->face_tag_[f];
1356
1357 assert ((f2c.size() == 1) || (f2c.size() == 2));
1358
1359 int inside_cell = 0;
1360
1361 if ( f2c.size() == 2 ) // Two cells => interior
1362 {
1363 if ( f2c[1].index() == cell )
1364 {
1365 inside_cell = 1;
1366 }
1367 }
1368 const bool normal_is_in = ! f2c[inside_cell].orientation();
1369
1370 switch (tag) {
1371 case I_FACE:
1372 // LEFT : RIGHT
1373 return normal_is_in ? 0 : 1; // min(I) : max(I)
1374 case J_FACE:
1375 // BACK : FRONT
1376 return normal_is_in ? 2 : 3; // min(J) : max(J)
1377 case K_FACE:
1378 // Note: TOP at min(K) as 'z' measures *depth*.
1379 // TOP : BOTTOM
1380 return normal_is_in ? 4 : 5; // min(K) : max(K)
1381 case NNC_FACE:
1382 // For nnc faces we return the otherwise unused value -1.
1383 return -1;
1384 default:
1385 OPM_THROW(std::logic_error, "Unhandled face tag. This should never happen!");
1386 }
1387 }
1388
1389 template<int dim>
1390 cpgrid::Entity<dim> createEntity(const CpGrid&, int, bool);
1391
1392} // namespace Dune
1393
1394#include <opm/grid/cpgrid/PersistentContainer.hpp>
1395#include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
1396#include "cpgrid/Intersection.hpp"
1397#include "cpgrid/Geometry.hpp"
1398#include "cpgrid/Indexsets.hpp"
1399
1400#endif // OPM_CPGRID_HEADER
An iterator over the centroids of the geometry of the entities.
Definition CpGrid.hpp:1013
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition CpGrid.hpp:1020
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition CpGrid.hpp:1017
[ provides Dune::Grid ]
Definition CpGrid.hpp:226
std::string name() const
Get the grid name.
Definition CpGrid.cpp:633
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:660
void switchToGlobalView()
Switch to the global view.
Definition CpGrid.cpp:1308
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition CpGrid.cpp:1350
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition CpGrid.cpp:980
void addLgrUpdateLeafView(const std::array< int, 3 > &cells_per_dim, const std::array< int, 3 > &startIJK, const std::array< int, 3 > &endIJK, const std::string &lgr_name)
Create a grid out of a coarse one and a refinement(LGR) of a selected block-shaped patch of cells fro...
Definition CpGrid.cpp:1422
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition CpGrid.hpp:264
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:1320
void globalRefine(int)
global refinement
Definition CpGrid.cpp:937
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition CpGrid.cpp:913
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition CpGrid.cpp:1242
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition CpGrid.cpp:908
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition CpGrid.hpp:1335
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition CpGrid.cpp:628
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:534
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:1395
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition CpGrid.cpp:958
cpgrid::CpGridDataTraits::InterfaceMap InterfaceMap
The type of the map describing communication interfaces.
Definition CpGrid.hpp:1107
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:1227
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition CpGrid.cpp:1077
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition CpGrid.cpp:604
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:769
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition CpGrid.cpp:1237
int maxLevel() const
Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview),...
Definition CpGrid.cpp:638
const CpGridTraits::Communication & comm() const
Get the collective communication object.
Definition CpGrid.cpp:1005
double faceArea(int face) const
Get the area of a face.
Definition CpGrid.cpp:1217
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition CpGrid.cpp:1222
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition CpGrid.cpp:918
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition CpGrid.cpp:884
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition CpGrid.cpp:597
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition CpGrid.cpp:1087
bool loadBalance(int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:599
int numVertices() const
Get The number of vertices.
Definition CpGrid.cpp:1026
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:618
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition CpGrid.cpp:1036
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition CpGrid.cpp:623
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:627
void addLgrsUpdateLeafView(const std::vector< std::array< int, 3 > > &cells_per_dim_vec, const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec, const std::vector< std::string > &lgr_name_vec)
Create a grid out of a coarse one and (at most) 2 refinements(LGRs) of selected block-shaped disjoint...
Definition CpGrid.cpp:1429
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:733
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:187
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:798
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition CpGrid.cpp:1046
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:1360
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition CpGrid.cpp:964
cpgrid::Entity< codim > entity(const cpgrid::Entity< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition CpGrid.cpp:948
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:689
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition CpGrid.cpp:932
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:1303
CpGrid()
Default constructor.
Definition CpGrid.cpp:160
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:1252
void switchToDistributedView()
Switch to the distributed view.
Definition CpGrid.cpp:1313
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition CpGrid.cpp:1247
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:847
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:1041
int numCells() const
Get the number of cells.
Definition CpGrid.cpp:1016
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition CpGrid.cpp:1298
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:822
int numFaces() const
Get the number of faces.
Definition CpGrid.cpp:1021
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition CpGrid.cpp:1212
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition CpGrid.cpp:1092
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition CpGrid.hpp:877
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:1232
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition CpGrid.hpp:1306
Definition CpGrid.hpp:71
Definition CpGrid.hpp:70
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:135
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition CpGridData.hpp:832
Represents an entity of a given codim, with positive or negative orientation.
Definition EntityRep.hpp:99
Definition Entity.hpp:65
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
Definition Intersection.hpp:278
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
LookUpCartesianData - To search data via CartesianIndex (cartesianMapper)
Definition LookUpData.hh:155
LookUpData class - To search data via element index.
Definition LookUpData.hh:64
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.cpp:68
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:213
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:167
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition CpGrid.hpp:169
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition CpGrid.hpp:171
Traits associated with a specific codim.
Definition CpGrid.hpp:143
cpgrid::Entity< cd > Entity
The type of the entity.
Definition CpGrid.hpp:152
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition CpGrid.hpp:146
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition CpGrid.hpp:149
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition CpGrid.hpp:158
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition CpGrid.hpp:155
cpgrid::Entity< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition CpGrid.hpp:161
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:179
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:183
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:181
Definition CpGrid.hpp:123
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition CpGrid.hpp:193
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition CpGrid.hpp:132
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:190
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:188
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition CpGrid.hpp:197
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition CpGrid.hpp:134
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition CpGrid.hpp:195
GlobalIdSet LocalIdSet
The type of the local id set.
Definition CpGrid.hpp:199
cpgrid::CpGridDataTraits::Communication Communication
The type of the collective communication.
Definition CpGrid.hpp:202
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition CpGrid.hpp:137
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition CpGrid.hpp:130
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition CpGrid.hpp:128
CpGrid Grid
The type that implements the grid.
Definition CpGrid.hpp:125
Raw corner-point specification of a particular geological model.
Definition preprocess.h:56