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#include <string>
42#include <map>
43#include <array>
44#include <unordered_set>
45#include <opm/grid/utility/ErrorMacros.hpp>
46
47// Warning suppression for Dune includes.
48#include <opm/grid/utility/platform_dependent/disable_warnings.h>
49
50#include <dune/common/version.hh>
51
52#if HAVE_MPI
53#if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
54#include <dune/common/parallel/variablesizecommunicator.hh>
55#else
56#include <opm/grid/utility/VariableSizeCommunicator.hpp>
57#endif
58#endif
59
60#include <dune/grid/common/capabilities.hh>
61#include <dune/grid/common/grid.hh>
62#include <dune/grid/common/gridenums.hh>
63
64#include <opm/grid/utility/platform_dependent/reenable_warnings.h>
65
66#include "cpgrid/Intersection.hpp"
67#include "cpgrid/Entity.hpp"
68#include "cpgrid/Geometry.hpp"
69#include "cpgrid/Iterators.hpp"
70#include "cpgrid/Indexsets.hpp"
71#include "cpgrid/DefaultGeometryPolicy.hpp"
72#include "common/GridEnums.hpp"
73#include "common/Volumes.hpp"
75
76#include <opm/grid/utility/OpmParserIncludes.hpp>
77
78#include <iostream>
79#if ! HAVE_MPI
80#include <list>
81#endif
82
83namespace Dune
84{
85
86 class CpGrid;
87
88 namespace cpgrid
89 {
90 class CpGridData;
91 }
92
94 //
95 // CpGridTraits
96 //
98
100 {
102 typedef CpGrid Grid;
103
112
115
118 template <int cd>
119 struct Codim
120 {
123 typedef cpgrid::Geometry<3-cd, 3> Geometry;
124 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> Geometry;
127 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> LocalGeometry;
130
133
136
139
142 template <PartitionIteratorType pitype>
144 {
149 };
150 };
151
154 template <PartitionIteratorType pitype>
156 {
158 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
160 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
161
162 };
163
165 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
167 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
168
177
179
180 typedef Dune::MPIHelper::MPICommunicator MPICommunicator;
181#if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
182 using Communication = Dune::Communication<MPICommunicator>;
183 using CollectiveCommunication = Communication;
184#else
185 using CollectiveCommunication = Dune::CollectiveCommunication<MPICommunicator>;
186 using Communication = Dune::CollectiveCommunication<MPICommunicator>;
187#endif
188 };
189
191 //
192 // CpGridFamily
193 //
195
197 {
198 typedef CpGridTraits Traits;
199 };
200
202 //
203 // CpGrid
204 //
206
208 class CpGrid
209 : public GridDefaultImplementation<3, 3, double, CpGridFamily >
210 {
211 friend class cpgrid::CpGridData;
212 template<int dim>
213 friend cpgrid::Entity<dim> createEntity(const CpGrid&,int,bool);
214
215 public:
216
217 // --- Typedefs ---
218
219
222
223
224 // --- Methods ---
225
226
228 CpGrid();
229
230 CpGrid(MPIHelper::MPICommunicator comm);
231
233
234
237 void readSintefLegacyFormat(const std::string& grid_prefix);
238
239
243 void writeSintefLegacyFormat(const std::string& grid_prefix) const;
244
245
246#if HAVE_ECL_INPUT
265 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
266 Opm::EclipseState* ecl_state,
267 bool periodic_extension, bool turn_normals, bool clip_z,
268 bool pinchActive);
269
289 std::vector<std::size_t> processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
290 Opm::EclipseState* ecl_state,
291 bool periodic_extension, bool turn_normals = false, bool clip_z = false);
292
293#endif
294
298 void processEclipseFormat(const grdecl& input_data, bool remove_ij_boundary, bool turn_normals = false);
299
301
307
308
311 void createCartesian(const std::array<int, 3>& dims,
312 const std::array<double, 3>& cellsize);
313
317 const std::array<int, 3>& logicalCartesianSize() const
318 {
319 return current_view_data_->logical_cartesian_size_;
320 }
321
329 const std::vector<int>& globalCell() const
330 {
331 return current_view_data_->global_cell_;
332 }
333
341 void getIJK(const int c, std::array<int,3>& ijk) const
342 {
343 current_view_data_->getIJK(c, ijk);
344 }
346
350 bool uniqueBoundaryIds() const
351 {
352 return current_view_data_->uniqueBoundaryIds();
353 }
354
357 void setUniqueBoundaryIds(bool uids)
358 {
359 current_view_data_->setUniqueBoundaryIds(uids);
360 }
361
362 // --- Dune interface below ---
363
365 // \@{
370 std::string name() const
371 {
372 return "CpGrid";
373 }
374
375
378 int maxLevel() const
379 {
380 return 0;
381 }
382
383
385 template<int codim>
386 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const
387 {
388 if (level<0 || level>maxLevel())
389 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
390 return cpgrid::Iterator<codim, All_Partition>(*current_view_data_, 0, true);
391 }
392
393
395 template<int codim>
396 typename Traits::template Codim<codim>::LevelIterator lend (int level) const
397 {
398 if (level<0 || level>maxLevel())
399 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
400 return cpgrid::Iterator<codim,All_Partition>(*current_view_data_, size(codim), true );
401
402 }
403
404
406 template<int codim, PartitionIteratorType PiType>
407 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const
408 {
409 if (level<0 || level>maxLevel())
410 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
411 return cpgrid::Iterator<codim,PiType>(*current_view_data_, 0, true );
412 }
413
414
416 template<int codim, PartitionIteratorType PiType>
417 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const
418 {
419 if (level<0 || level>maxLevel())
420 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
421 return cpgrid::Iterator<codim,PiType>(*current_view_data_, size(codim), true);
422 }
423
424
426 template<int codim>
427 typename Traits::template Codim<codim>::LeafIterator leafbegin() const
428 {
429 return cpgrid::Iterator<codim, All_Partition>(*current_view_data_, 0, true);
430 }
431
432
434 template<int codim>
435 typename Traits::template Codim<codim>::LeafIterator leafend() const
436 {
437 return cpgrid::Iterator<codim, All_Partition>(*current_view_data_, size(codim), true);
438 }
439
440
442 template<int codim, PartitionIteratorType PiType>
443 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const
444 {
445 return cpgrid::Iterator<codim, PiType>(*current_view_data_, 0, true);
446 }
447
448
450 template<int codim, PartitionIteratorType PiType>
451 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const
452 {
453 return cpgrid::Iterator<codim, PiType>(*current_view_data_, size(codim), true);
454 }
455
456
458 int size (int level, int codim) const
459 {
460 if (level<0 || level>maxLevel())
461 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
462 return size(codim);
463 }
464
465
467 int size (int codim) const
468 {
469 return current_view_data_->size(codim);
470 }
471
472
474 int size (int level, GeometryType type) const
475 {
476 if (level<0 || level>maxLevel())
477 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
478 return size(type);
479 }
480
481
483 int size (GeometryType type) const
484 {
485 return current_view_data_->size(type);
486 }
487
488
490 const Traits::GlobalIdSet& globalIdSet() const
491 {
492 return global_id_set_;
493 }
494
495
497 const Traits::LocalIdSet& localIdSet() const
498 {
499 return global_id_set_;
500 }
501
502
504 const Traits::LevelIndexSet& levelIndexSet(int level) const
505 {
506 if (level<0 || level>maxLevel())
507 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!");
508 return *current_view_data_->index_set_;
509 }
510
511
513 const Traits::LeafIndexSet& leafIndexSet() const
514 {
515 return *current_view_data_->index_set_;
516 }
517
518
520 void globalRefine (int)
521 {
522 std::cout << "Warning: Global refinement not implemented, yet." << std::endl;
523 }
524
525 const std::vector< Dune :: GeometryType >& geomTypes( const int codim ) const
526 {
527 return leafIndexSet().geomTypes( codim );
528 }
529
531 template <int codim>
533 {
534 return seed;
535 }
536
537 /* No refinement implemented. GridDefaultImplementation's methods will be used.
538
548
549 bool mark(int refCount, const typename Traits::template Codim<0>::EntityPointer & e)
550 {
551 return hostgrid_->mark(refCount, getHostEntity<0>(*e));
552 }
553
557
558 int getMark(const typename Traits::template Codim<0>::EntityPointer & e) const
559 {
560 return hostgrid_->getMark(getHostEntity<0>(*e));
561 }
562
564 bool preAdapt() {
565 return hostgrid_->preAdapt();
566 }
567
568
570 bool adapt()
571 {
572 return hostgrid_->adapt();
573 }
574
576 void postAdapt() {
577 return hostgrid_->postAdapt();
578 }
579
580 end of refinement section */
581
583 unsigned int overlapSize(int) const {
584 return 1;
585 }
586
587
589 unsigned int ghostSize(int) const {
590 return 0;
591 }
592
593
595 unsigned int overlapSize(int, int) const {
596 return 1;
597 }
598
599
601 unsigned int ghostSize(int, int) const {
602 return 0;
603 }
604
606 unsigned int numBoundarySegments() const
607 {
608 if( uniqueBoundaryIds() )
609 {
610 return current_view_data_->unique_boundary_ids_.size();
611 }
612 else
613 {
614 unsigned int numBndSegs = 0;
615 const int num_faces = numFaces();
616 for (int i = 0; i < num_faces; ++i) {
617 cpgrid::EntityRep<1> face(i, true);
618 if (current_view_data_->face_to_cell_[face].size() == 1) {
619 ++numBndSegs;
620 }
621 }
622 return numBndSegs;
623 }
624 }
625
626 void setZoltanParams(const std::map<std::string,std::string>& params)
627 {
628 zoltanParams = params;
629 }
630
631 // loadbalance is not part of the grid interface therefore we skip it.
632
638 bool loadBalance(int overlapLayers=1, bool useZoltan=true)
639 {
640 using std::get;
641 return get<0>(scatterGrid(defaultTransEdgeWgt, false, nullptr, false, nullptr, true, overlapLayers, useZoltan ));
642 }
643
644 // loadbalance is not part of the grid interface therefore we skip it.
645
665 std::pair<bool, std::vector<std::pair<std::string,bool> > >
666 loadBalance(const std::vector<cpgrid::OpmWellType> * wells,
667 const double* transmissibilities = nullptr,
668 int overlapLayers=1, bool useZoltan=true)
669 {
670 return scatterGrid(defaultTransEdgeWgt, false, wells, false, transmissibilities, false, overlapLayers, useZoltan);
671 }
672
673 // loadbalance is not part of the grid interface therefore we skip it.
674
698 std::pair<bool, std::vector<std::pair<std::string,bool> > >
699 loadBalance(EdgeWeightMethod method, const std::vector<cpgrid::OpmWellType> * wells,
700 const double* transmissibilities = nullptr, bool ownersFirst=false,
701 bool addCornerCells=false, int overlapLayers=1,
702 bool useZoltan = true)
703 {
704 return scatterGrid(method, ownersFirst, wells, false, transmissibilities, addCornerCells, overlapLayers, useZoltan);
705 }
706
726 template<class DataHandle>
727 std::pair<bool, std::vector<std::pair<std::string,bool> > >
728 loadBalance(DataHandle& data,
729 const std::vector<cpgrid::OpmWellType> * wells,
730 const double* transmissibilities = nullptr,
731 int overlapLayers=1, bool useZoltan = true)
732 {
733 auto ret = loadBalance(wells, transmissibilities, overlapLayers, useZoltan);
734 using std::get;
735 if (get<0>(ret))
736 {
737 scatterData(data);
738 }
739 return ret;
740 }
741
770 template<class DataHandle>
771 std::pair<bool, std::vector<std::pair<std::string,bool> > >
772 loadBalance(DataHandle& data, EdgeWeightMethod method,
773 const std::vector<cpgrid::OpmWellType> * wells,
774 bool serialPartitioning,
775 const double* transmissibilities = nullptr, bool ownersFirst=false,
776 bool addCornerCells=false, int overlapLayers=1, bool useZoltan = true,
777 double zoltanImbalanceTol = 1.1,
778 bool allowDistributedWells = false)
779 {
780 auto ret = scatterGrid(method, ownersFirst, wells, serialPartitioning, transmissibilities,
781 addCornerCells, overlapLayers, useZoltan, zoltanImbalanceTol, allowDistributedWells);
782 using std::get;
783 if (get<0>(ret))
784 {
785 scatterData(data);
786 }
787 return ret;
788 }
789
806 template<class DataHandle>
807 std::pair<bool, std::vector<std::pair<std::string,bool> > >
808 loadBalance(DataHandle& data, const std::vector<int>& parts,
809 const std::vector<cpgrid::OpmWellType> * wells,
810 bool ownersFirst=false,
811 bool addCornerCells=false, int overlapLayers=1)
812 {
813 using std::get;
814 auto ret = scatterGrid(defaultTransEdgeWgt, ownersFirst, wells,
815 /* serialPartitioning = */ false,
816 /* transmissibilities = */ {},
817 addCornerCells, overlapLayers, /* useZoltan =*/ false,
818 /* zoltanImbalanceTol (ignored) = */ 0.0,
819 /* allowDistributedWells = */ true, parts);
820 using std::get;
821 if (get<0>(ret))
822 {
823 scatterData(data);
824 }
825 return ret;
826 }
835 template<class DataHandle>
836 bool loadBalance(DataHandle& data,
837#if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7)
838 decltype(data.fixedSize(0,0)) overlapLayers=1, bool useZoltan = true)
839#else
840 decltype(data.fixedsize(0,0)) overlapLayers=1, bool useZoltan = true)
841#endif
842 {
843 // decltype usage needed to tell the compiler not to use this function if first
844 // argument is std::vector but rather loadbalance by parts
845 bool ret = loadBalance(overlapLayers, useZoltan);
846 if (ret)
847 {
848 scatterData(data);
849 }
850 return ret;
851 }
852
864 bool loadBalance(const std::vector<int>& parts, bool ownersFirst=false,
865 bool addCornerCells=false, int overlapLayers=1)
866 {
867 using std::get;
868 return get<0>(scatterGrid(defaultTransEdgeWgt, ownersFirst, /* wells = */ {},
869 /* serialPartitioning = */ false,
870 /* trabsmissibilities = */ {},
871 addCornerCells, overlapLayers, /* useZoltan =*/ false,
872 /* zoltanImbalanceTol (ignored) = */ 0.0,
873 /* allowDistributedWells = */ true, parts));
874 }
875
888 template<class DataHandle>
889 bool loadBalance(DataHandle& data, const std::vector<int>& parts, bool ownersFirst=false,
890 bool addCornerCells=false, int overlapLayers=1)
891 {
892 bool ret = loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
893 if (ret)
894 {
895 scatterData(data);
896 }
897 return ret;
898 }
899
905 std::vector<int> zoltanPartitionWithoutScatter(const std::vector<cpgrid::OpmWellType> * wells,
906 const double* transmissibilities, int numParts,
907 const double zoltanImbalanceTol);
908
916 template<class DataHandle>
917 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int /*level*/) const
918 {
919 communicate(data, iftype, dir);
920 }
921
929 template<class DataHandle>
930 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
931 {
932 current_view_data_->communicate(data, iftype, dir);
933 }
934
936 const typename CpGridTraits::Communication& comm () const
937 {
938 return current_view_data_->ccobj_;
939 }
941
942 // ------------ End of Dune interface, start of simplified interface --------------
943
949
950 // enum { dimension = 3 }; // already defined
951
952 typedef Dune::FieldVector<double, 3> Vector;
953
954
955 const std::vector<double>& zcornData() const {
956 return current_view_data_->zcornData();
957 }
958
959
960 // Topology
962 int numCells() const
963 {
964 return current_view_data_->cell_to_face_.size();
965 }
967 int numFaces() const
968 {
969 return current_view_data_->face_to_cell_.size();
970 }
972 int numVertices() const
973 {
974 return current_view_data_->geomVector<3>().size();
975 }
982 int numCellFaces(int cell) const
983 {
984 return current_view_data_->cell_to_face_[cpgrid::EntityRep<0>(cell, true)].size();
985 }
990 int cellFace(int cell, int local_index) const
991 {
992 return current_view_data_->cell_to_face_[cpgrid::EntityRep<0>(cell, true)][local_index].index();
993 }
994
998 {
999 return current_view_data_->cell_to_face_[cpgrid::EntityRep<0>(cell, true)];
1000 }
1011 int faceCell(int face, int local_index) const
1012 {
1013 // In the parallel case we store non-existent cells for faces along
1014 // the front region. Theses marked with index std::numeric_limits<int>::max(),
1015 // orientation might be arbitrary, though.
1017 = current_view_data_->face_to_cell_[cpgrid::EntityRep<1>(face, true)];
1018 bool a = (local_index == 0);
1019 bool b = r[0].orientation();
1020 bool use_first = a ? b : !b;
1021 // The number of valid cells.
1022 int r_size = r.size();
1023 // In the case of only one valid cell, this is the index of it.
1024 int index = 0;
1025 if(r[0].index()==std::numeric_limits<int>::max()){
1026 assert(r_size==2);
1027 --r_size;
1028 index=1;
1029 }
1030 if(r.size()>1 && r[1].index()==std::numeric_limits<int>::max())
1031 {
1032 assert(r_size==2);
1033 --r_size;
1034 }
1035 if (r_size == 2) {
1036 return use_first ? r[0].index() : r[1].index();
1037 } else {
1038 return use_first ? r[index].index() : -1;
1039 }
1040 }
1047 int numCellFaces() const
1048 {
1049 return current_view_data_->cell_to_face_.dataSize();
1050 }
1051 int numFaceVertices(int face) const
1052 {
1053 return current_view_data_->face_to_point_[face].size();
1054 }
1059 int faceVertex(int face, int local_index) const
1060 {
1061 return current_view_data_->face_to_point_[face][local_index];
1062 }
1065 double cellCenterDepth(int cell_index) const
1066 {
1067 // Here cell center depth is computed as a raw average of cell corner depths.
1068 // This generally gives slightly different results than using the cell centroid.
1069 double zz = 0.0;
1070 const int nv = current_view_data_->cell_to_point_[cell_index].size();
1071 const int nd = 3;
1072 for (int i=0; i<nv; ++i) {
1073 zz += vertexPosition(current_view_data_->cell_to_point_[cell_index][i])[nd-1];
1074 }
1075 return zz/nv;
1076 }
1077
1078 const Vector faceCenterEcl(int cell_index, int face) const
1079 {
1080 // This method is an alternative to the method faceCentroid(...).
1081 // The face center is computed as a raw average of cell corners.
1082 // For faulted cells this gives different results then average of face nodes
1083 // that seems to agree more with eclipse.
1084 // This assumes the cell nodes are ordered
1085 // 6---7
1086 // | T |
1087 // 4---5
1088 // 2---3
1089 // | B |
1090 // 0---1
1091
1092 // this follows the DUNE reference cube
1093 static const int faceVxMap[ 6 ][ 4 ] = { {0, 2, 4, 6}, // face 0
1094 {1, 3, 5, 7}, // face 1
1095 {0, 1, 4, 5}, // face 2
1096 {2, 3, 6, 7}, // face 3
1097 {0, 1, 2, 3}, // face 4
1098 {4, 5, 6, 7} // face 5
1099 };
1100
1101
1102 assert (current_view_data_->cell_to_point_[cell_index].size() == 8);
1103 Vector center(0.0);
1104 for( int i=0; i<4; ++i )
1105 {
1106 center += vertexPosition(current_view_data_->cell_to_point_[cell_index][ faceVxMap[ face ][ i ] ]);
1107 }
1108
1109 for (int i=0; i<3; ++i) {
1110 center[i] /= 4;
1111 }
1112 return center;
1113
1114 }
1115
1116 const Vector faceAreaNormalEcl(int face) const
1117 {
1118 // same implementation as ResInsight
1119 const int nd = Vector::dimension;
1120 const int nv = numFaceVertices(face);
1121 switch (nv)
1122 {
1123 case 0:
1124 case 1:
1125 case 2:
1126 {
1127 return Vector(0.0);
1128 }
1129 break;
1130 case 3:
1131 {
1132 Vector a = vertexPosition(current_view_data_->face_to_point_[face][0]) - vertexPosition(current_view_data_->face_to_point_[face][2]);
1133 Vector b = vertexPosition(current_view_data_->face_to_point_[face][1]) - vertexPosition(current_view_data_->face_to_point_[face][2]);
1134 Vector areaNormal = cross(a,b);
1135 for (int i=0; i<nd; ++i) {
1136 areaNormal[i] /= 2;
1137 }
1138 return areaNormal;
1139 }
1140 break;
1141 case 4:
1142 {
1143 Vector a = vertexPosition(current_view_data_->face_to_point_[face][0]) - vertexPosition(current_view_data_->face_to_point_[face][2]);
1144 Vector b = vertexPosition(current_view_data_->face_to_point_[face][1]) - vertexPosition(current_view_data_->face_to_point_[face][3]);
1145 Vector areaNormal = cross(a,b);
1146 areaNormal *= 0.5;
1147 return areaNormal;
1148 }
1149 break;
1150 default:
1151 {
1152 int h = (nv - 1)/2;
1153 int k = (nv % 2) ? 0 : nv - 1;
1154
1155 Vector areaNormal(0.0);
1156 // First quads
1157 for (int i = 1; i < h; ++i)
1158 {
1159 Vector a = vertexPosition(current_view_data_->face_to_point_[face][2*i]) - vertexPosition(current_view_data_->face_to_point_[face][0]);
1160 Vector b = vertexPosition(current_view_data_->face_to_point_[face][2*i+1]) - vertexPosition(current_view_data_->face_to_point_[face][2*i-1]);
1161 areaNormal += cross(a,b);
1162 }
1163
1164 // Last triangle or quad
1165 Vector a = vertexPosition(current_view_data_->face_to_point_[face][2*h]) - vertexPosition(current_view_data_->face_to_point_[face][0]);
1166 Vector b = vertexPosition(current_view_data_->face_to_point_[face][k]) - vertexPosition(current_view_data_->face_to_point_[face][2*h-1]);
1167 areaNormal += cross(a,b);
1168
1169 areaNormal *= 0.5;
1170
1171 return areaNormal;
1172 }
1173
1174 }
1175 }
1176
1177 // Geometry
1181 const Vector& vertexPosition(int vertex) const
1182 {
1183 return current_view_data_->geomVector<3>()[cpgrid::EntityRep<3>(vertex, true)].center();
1184 }
1187 double faceArea(int face) const
1188 {
1189 return current_view_data_->geomVector<1>()[cpgrid::EntityRep<1>(face, true)].volume();
1190 }
1193 const Vector& faceCentroid(int face) const
1194 {
1195 return current_view_data_->geomVector<1>()[cpgrid::EntityRep<1>(face, true)].center();
1196 }
1200 const Vector& faceNormal(int face) const
1201 {
1202 return current_view_data_->face_normals_.get(face);
1203 }
1206 double cellVolume(int cell) const
1207 {
1208 return current_view_data_->geomVector<0>()[cpgrid::EntityRep<0>(cell, true)].volume();
1209 }
1212 const Vector& cellCentroid(int cell) const
1213 {
1214 return current_view_data_->geomVector<0>()[cpgrid::EntityRep<0>(cell, true)].center();
1215 }
1216
1219 template<int codim>
1221 : public RandomAccessIteratorFacade<CentroidIterator<codim>,
1222 FieldVector<double, 3>,
1223 const FieldVector<double, 3>&, int>
1224 {
1225 public:
1227 typedef typename std::vector<cpgrid::Geometry<3-codim, 3> >::const_iterator
1232 : iter_(iter)
1233 {}
1234
1235 const FieldVector<double, 3>& dereference() const
1236 {
1237 return iter_->center();
1238 }
1239 void increment()
1240 {
1241 ++iter_;
1242 }
1243 const FieldVector<double, 3>& elementAt(int n)
1244 {
1245 return iter_[n]->center();
1246 }
1247 void advance(int n)
1248 {
1249 iter_+=n;
1250 }
1251 void decrement()
1252 {
1253 --iter_;
1254 }
1255 int distanceTo(const CentroidIterator& o)
1256 {
1257 return o-iter_;
1258 }
1259 bool equals(const CentroidIterator& o) const
1260 {
1261 return o==iter_;
1262 }
1263 private:
1265 GeometryIterator iter_;
1266 };
1267
1270 {
1271 return CentroidIterator<0>(current_view_data_->geomVector<0>().begin());
1272 }
1273
1276 {
1277 return CentroidIterator<1>(current_view_data_->geomVector<1>().begin());
1278 }
1279
1280 // Extra
1281 int boundaryId(int face) const
1282 {
1283 // Note that this relies on the following implementation detail:
1284 // The grid is always construct such that the faces where
1285 // orientation() returns true are oriented along the positive IJK
1286 // direction. Oriented means that the first cell attached to face
1287 // has the lower index.
1288 int ret = 0;
1289 cpgrid::EntityRep<1> f(face, true);
1290 if (current_view_data_->face_to_cell_[f].size() == 1) {
1291 if (current_view_data_->uniqueBoundaryIds()) {
1292 // Use the unique boundary ids.
1293 ret = current_view_data_->unique_boundary_ids_[f];
1294 } else {
1295 // Use the face tag based ids, i.e. 1-6 for i-, i+, j-, j+, k-, k+.
1296 const bool normal_is_in =
1297 !(current_view_data_->face_to_cell_[f][0].orientation());
1298 enum face_tag tag = current_view_data_->face_tag_[f];
1299 switch (tag) {
1300 case I_FACE:
1301 // LEFT : RIGHT
1302 ret = normal_is_in ? 1 : 2; // min(I) : max(I)
1303 break;
1304 case J_FACE:
1305 // BACK : FRONT
1306 ret = normal_is_in ? 3 : 4; // min(J) : max(J)
1307 break;
1308 case K_FACE:
1309 // Note: TOP at min(K) as 'z' measures *depth*.
1310 // TOP : BOTTOM
1311 ret = normal_is_in ? 5 : 6; // min(K) : max(K)
1312 break;
1313 case NNC_FACE:
1314 // This should not be possible, as NNC "faces" always
1315 // have two cell neighbours.
1316 OPM_THROW(std::logic_error, "NNC face at boundary. This should never happen!");
1317 }
1318 }
1319 }
1320 return ret;
1321 }
1322
1329 template<class Cell2FacesRowIterator>
1330 int
1331 faceTag(const Cell2FacesRowIterator& cell_face) const
1332 {
1333 // Note that this relies on the following implementation detail:
1334 // The grid is always constructed such that the interior faces constructed
1335 // with orientation set to true are
1336 // oriented along the positive IJK direction. Oriented means that
1337 // the first cell attached to face has the lower index.
1338 // For faces along the boundary (only one cell, always attached at index 0)
1339 // the orientation has to be determined by the orientation of the cell.
1340 // If it is true then in UnstructuredGrid it would be stored at index 0,
1341 // otherwise at index 1.
1342 const int cell = cell_face.getCellIndex();
1343 const int face = *cell_face;
1344 assert (0 <= cell); assert (cell < numCells());
1345 assert (0 <= face); assert (face < numFaces());
1346
1348
1349 const cpgrid::EntityRep<1> f(face, true);
1350 const F2C& f2c = current_view_data_->face_to_cell_[f];
1351 const face_tag tag = current_view_data_->face_tag_[f];
1352
1353 assert ((f2c.size() == 1) || (f2c.size() == 2));
1354
1355 int inside_cell = 0;
1356
1357 if ( f2c.size() == 2 ) // Two cells => interior
1358 {
1359 if ( f2c[1].index() == cell )
1360 {
1361 inside_cell = 1;
1362 }
1363 }
1364 const bool normal_is_in = ! f2c[inside_cell].orientation();
1365
1366 switch (tag) {
1367 case I_FACE:
1368 // LEFT : RIGHT
1369 return normal_is_in ? 0 : 1; // min(I) : max(I)
1370 case J_FACE:
1371 // BACK : FRONT
1372 return normal_is_in ? 2 : 3; // min(J) : max(J)
1373 case K_FACE:
1374 // Note: TOP at min(K) as 'z' measures *depth*.
1375 // TOP : BOTTOM
1376 return normal_is_in ? 4 : 5; // min(K) : max(K)
1377 case NNC_FACE:
1378 // For nnc faces we return the otherwise unused value -1.
1379 return -1;
1380 default:
1381 OPM_THROW(std::logic_error, "Unhandled face tag. This should never happen!");
1382 }
1383 }
1384
1386
1387 // ------------ End of simplified interface --------------
1388
1389 //------------- methods not in the DUNE grid interface.
1390
1395
1396
1405 template<class DataHandle>
1406 void scatterData(DataHandle& handle) const
1407 {
1408#if HAVE_MPI
1409 if(distributed_data_.empty())
1410 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balanced grid!");
1411 distributed_data_[0]->scatterData(handle, data_[0].get(), distributed_data_[0].get(), cellScatterGatherInterface(),
1413#else
1414 // Suppress warnings for unused argument.
1415 (void) handle;
1416#endif
1417 }
1418
1425 template<class DataHandle>
1426 void gatherData(DataHandle& handle) const
1427 {
1428#if HAVE_MPI
1429 if(distributed_data_.empty())
1430 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balance grid!");
1431 distributed_data_[0]->gatherData(handle, data_[0].get(), distributed_data_[0].get());
1432#else
1433 // Suppress warnings for unused argument.
1434 (void) handle;
1435#endif
1436 }
1437#if HAVE_MPI
1438#if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
1440 using InterfaceMap = VariableSizeCommunicator<>::InterfaceMap;
1441#else
1443 using InterfaceMap = Opm::VariableSizeCommunicator<>::InterfaceMap;
1444#endif
1445#else
1446 // bogus definition for the non parallel type. VariableSizeCommunicator not
1447 // availabe
1448
1450 typedef std::map<int, std::list<int> > InterfaceMap;
1451#endif
1452
1482 {
1483 return *cell_scatter_gather_interfaces_;
1484 }
1485
1489 {
1490 return *point_scatter_gather_interfaces_;
1491 }
1492
1495 {
1496 current_view_data_=data_[0].get();
1497 }
1498
1501 {
1502 if (distributed_data_.empty())
1503 OPM_THROW(std::logic_error, "No distributed view available in grid");
1504 current_view_data_=distributed_data_[0].get();
1505 }
1507
1508#if HAVE_MPI
1510 typedef cpgrid::CpGridData::ParallelIndexSet ParallelIndexSet;
1512 typedef cpgrid::CpGridData::RemoteIndices RemoteIndices;
1513
1515 using CommunicationType = cpgrid::CpGridData::CommunicationType;
1516
1520 const CommunicationType& cellCommunication() const
1521 {
1522 return current_view_data_->cellCommunication();
1523 }
1524
1525 ParallelIndexSet& getCellIndexSet()
1526 {
1527 return current_view_data_->cellIndexSet();
1528 }
1529
1530 RemoteIndices& getCellRemoteIndices()
1531 {
1532 return current_view_data_->cellRemoteIndices();
1533 }
1534
1535 const ParallelIndexSet& getCellIndexSet() const
1536 {
1537 return current_view_data_->cellIndexSet();
1538 }
1539
1540 const RemoteIndices& getCellRemoteIndices() const
1541 {
1542 return current_view_data_->cellRemoteIndices();
1543 }
1544
1545#endif
1546
1548 const std::vector<int>& sortedNumAquiferCells() const
1549 {
1550 return current_view_data_->sortedNumAquiferCells();
1551 }
1552
1553 private:
1579 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1580 scatterGrid(EdgeWeightMethod method,
1581 bool ownersFirst,
1582 const std::vector<cpgrid::OpmWellType> * wells,
1583 bool serialPartitioning,
1584 const double* transmissibilities,
1585 bool addCornerCells,
1586 int overlapLayers,
1587 bool useZoltan = true,
1588 double zoltanImbalanceTol = 1.1,
1589 bool allowDistributedWells = true,
1590 const std::vector<int>& input_cell_part = {});
1591
1596 std::vector<std::shared_ptr<cpgrid::CpGridData>> data_;
1598 cpgrid::CpGridData* current_view_data_;
1600 std::vector<std::shared_ptr<cpgrid::CpGridData>> distributed_data_;
1606 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1607 /*
1608 * @brief Interface for scattering and gathering point data.
1609 *
1610 * @warning Will only update owner cells
1611 */
1612 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1616 cpgrid::GlobalIdSet global_id_set_;
1617
1621 std::map<std::string,std::string> zoltanParams;
1622
1623 }; // end Class CpGrid
1624
1625
1626
1627 namespace Capabilities
1628 {
1630 template <>
1631 struct hasEntity<CpGrid, 0>
1632 {
1633 static const bool v = true;
1634 };
1635
1637 template <>
1638 struct hasEntity<CpGrid, 3>
1639 {
1640 static const bool v = true;
1641 };
1642
1643 template<>
1644 struct canCommunicate<CpGrid,0>
1645 {
1646 static const bool v = true;
1647 };
1648
1649 template<>
1650 struct canCommunicate<CpGrid,3>
1651 {
1652 static const bool v = true;
1653 };
1654
1656 template <>
1657 struct hasBackupRestoreFacilities<CpGrid>
1658 {
1659 static const bool v = false;
1660 };
1661
1662 }
1663
1664
1665 template<int dim>
1666 cpgrid::Entity<dim> createEntity(const CpGrid& grid,int index,bool orientation)
1667 {
1668 return cpgrid::Entity<dim>(*grid.current_view_data_, index, orientation);
1669 }
1670
1671} // namespace Dune
1672
1673#include <opm/grid/cpgrid/PersistentContainer.hpp>
1674#include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
1675#endif // OPM_CPGRID_HEADER
An iterator over the centroids of the geometry of the entities.
Definition: CpGrid.hpp:1224
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition: CpGrid.hpp:1231
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition: CpGrid.hpp:1228
[ provides Dune::Grid ]
Definition: CpGrid.hpp:210
std::string name() const
Get the grid name.
Definition: CpGrid.hpp:370
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition: CpGrid.hpp:1269
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:699
void switchToGlobalView()
Switch to the global view.
Definition: CpGrid.hpp:1494
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:527
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: CpGrid.hpp:606
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition: CpGrid.hpp:221
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:386
void gatherData(DataHandle &handle) const
Moves data from the distributed view to the global (all data on process) view.
Definition: CpGrid.hpp:1426
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition: CpGrid.hpp:1275
void globalRefine(int)
global refinement
Definition: CpGrid.hpp:520
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: CpGrid.hpp:483
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition: CpGrid.hpp:1331
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGrid.hpp:357
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:435
void processEclipseFormat(const grdecl &input_data, bool remove_ij_boundary, bool turn_normals=false)
Read the Eclipse grid format ('grdecl').
Definition: CpGrid.cpp:572
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition: CpGrid.hpp:583
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:396
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition: CpGrid.hpp:1047
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition: CpGrid.hpp:329
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:808
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: CpGrid.hpp:497
int maxLevel() const
Return maximum level defined in this grid.
Definition: CpGrid.hpp:378
const CpGridTraits::Communication & comm() const
Get the collective communication object.
Definition: CpGrid.hpp:936
double faceArea(int face) const
Get the area of a face.
Definition: CpGrid.hpp:1187
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition: CpGrid.hpp:1181
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition: CpGrid.hpp:458
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition: CpGrid.hpp:317
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition: CpGrid.hpp:1059
bool loadBalance(int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:638
std::vector< int > zoltanPartitionWithoutScatter(const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities, int numParts, const double zoltanImbalanceTol)
Partitions the grid using Zoltan without decomposing and distributing it among processes.
Definition: CpGrid.cpp:139
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:836
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition: CpGrid.hpp:1193
unsigned int ghostSize(int, int) const
Size of the ghost cell layer on a given level.
Definition: CpGrid.hpp:601
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition: CpGrid.hpp:1481
int numVertices() const
Get The number of vertices.
Definition: CpGrid.hpp:972
void getIJK(const int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGrid.hpp:341
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition: CpGrid.hpp:990
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGrid.hpp:350
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:666
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:417
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:772
void createCartesian(const std::array< int, 3 > &dims, const std::array< double, 3 > &cellsize)
Create a cartesian grid.
Definition: CpGrid.cpp:469
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: CpGrid.hpp:490
cpgrid::Entity< codim > entity(const cpgrid::Entity< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition: CpGrid.hpp:532
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.hpp:1011
int numCellFaces(int cell) const
Get the number of faces of a cell.
Definition: CpGrid.hpp:982
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:443
unsigned int overlapSize(int, int) const
Size of the overlap on a given level.
Definition: CpGrid.hpp:595
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:537
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition: CpGrid.hpp:589
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir) const
The new communication interface.
Definition: CpGrid.hpp:930
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition: CpGrid.hpp:1212
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:728
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:427
const Vector & faceNormal(int face) const
Get the unit normal of a face.
Definition: CpGrid.hpp:1200
CpGrid()
Default constructor.
Definition: CpGrid.cpp:119
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:407
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition: CpGrid.hpp:1548
void switchToDistributedView()
Switch to the distributed view.
Definition: CpGrid.hpp:1500
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:889
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: CpGrid.hpp:513
const cpgrid::OrientedEntityTable< 0, 1 >::row_type cellFaceRow(int cell) const
Get a list of indices identifying all faces of a cell.
Definition: CpGrid.hpp:997
int numCells() const
Get the number of cells.
Definition: CpGrid.hpp:962
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:864
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: CpGrid.hpp:474
int numFaces() const
Get the number of faces.
Definition: CpGrid.hpp:967
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition: CpGrid.hpp:1065
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition: CpGrid.hpp:917
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: CpGrid.hpp:504
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGrid.hpp:467
const InterfaceMap & pointScatterGatherInterface() const
Get an interface for gathering/scattering data attached to points with communication.
Definition: CpGrid.hpp:1488
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:451
double cellVolume(int cell) const
Get the volume of the cell.
Definition: CpGrid.hpp:1206
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition: CpGrid.hpp:1406
std::map< int, std::list< int > > InterfaceMap
The type of the map describing communication interfaces.
Definition: CpGrid.hpp:1450
Struct that hods all the data needed to represent a Cpgrid.
Definition: CpGridData.hpp:123
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition: CpGridData.hpp:614
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGridData.hpp:249
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGridData.cpp:144
void getIJK(int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGridData.hpp:235
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:267
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGridData.hpp:256
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition: CpGridData.hpp:368
Represents an entity of a given codim, with positive or negative orientation.
Definition: EntityRep.hpp:98
This class encapsulates geometry for both vertices, intersections and cells.
Definition: Geometry.hpp:71
The global id set for Dune.
Definition: Indexsets.hpp:325
Only needs to provide interface for doing nothing.
Definition: Iterators.hpp:104
Definition: Indexsets.hpp:54
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:56
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
int size() const
Returns the number of rows in the table.
Definition: SparseTable.hpp:121
int dataSize() const
Returns the number of data elements.
Definition: SparseTable.hpp:141
int size() const
Returns the number of rows in the table.
Definition: SparseTable.hpp:121
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
FieldVector< T, 3 > cross(const FieldVector< T, 3 > &a, const FieldVector< T, 3 > &b)
Definition: Volumes.hpp:58
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
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:197
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:144
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition: CpGrid.hpp:146
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition: CpGrid.hpp:148
Traits associated with a specific codim.
Definition: CpGrid.hpp:120
cpgrid::Entity< cd > Entity
The type of the entity.
Definition: CpGrid.hpp:129
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition: CpGrid.hpp:123
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition: CpGrid.hpp:126
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition: CpGrid.hpp:135
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition: CpGrid.hpp:132
cpgrid::Entity< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:138
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:156
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:160
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:158
Definition: CpGrid.hpp:100
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition: CpGrid.hpp:170
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition: CpGrid.hpp:109
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:167
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:165
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition: CpGrid.hpp:174
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition: CpGrid.hpp:111
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition: CpGrid.hpp:172
GlobalIdSet LocalIdSet
The type of the local id set.
Definition: CpGrid.hpp:176
Dune::MPIHelper::MPICommunicator MPICommunicator
The type of the collective communication.
Definition: CpGrid.hpp:180
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition: CpGrid.hpp:114
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition: CpGrid.hpp:107
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition: CpGrid.hpp:105
CpGrid Grid
The type that implements the grid.
Definition: CpGrid.hpp:102
Raw corner-point specification of a particular geological model.
Definition: preprocess.h:56