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/DefaultGeometryPolicy.hpp>
48#include <opm/grid/cpgrid/OrientedEntityTable.hpp>
50#include <opm/grid/utility/platform_dependent/reenable_warnings.h> // Not really needed it seems, but alas.
51#include "common/GridEnums.hpp"
52#include <opm/grid/utility/OpmWellType.hpp>
53
54#include <set>
55
56namespace Opm
57{
58struct NNCdata;
59class EclipseGrid;
60class EclipseState;
61}
62
63namespace Dune
64{
65 class CpGrid;
66
67 namespace cpgrid
68 {
69 class CpGridData;
70 template <int> class Entity;
71 template<int,int> class Geometry;
72 class HierarchicIterator;
73 class IntersectionIterator;
74 template<int, PartitionIteratorType> class Iterator;
75 class LevelGlobalIdSet;
76 class GlobalIdSet;
77 class Intersection;
78 class IntersectionIterator;
79 class IndexSet;
80 class IdSet;
81
82 }
83}
84
85namespace Dune
86{
87
89 //
90 // CpGridTraits
91 //
93
95 {
97 typedef CpGrid Grid;
98
107
110
113 template <int cd>
114 struct Codim
115 {
118 typedef cpgrid::Geometry<3-cd, 3> Geometry;
119 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> Geometry;
122 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> LocalGeometry;
125
128
131
134
137 template <PartitionIteratorType pitype>
145 };
146
149 template <PartitionIteratorType pitype>
151 {
153 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
155 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
156
157 };
158
160 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid>> LevelGridView;
162 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid>> LeafGridView;
163
172
174 using Communication = cpgrid::CpGridDataTraits::Communication;
175 using CollectiveCommunication = cpgrid::CpGridDataTraits::CollectiveCommunication;
176 };
177
179 //
180 // CpGridFamily
181 //
183
185 {
186 typedef CpGridTraits Traits;
187 };
188
190 //
191 // CpGrid
192 //
194
196 class CpGrid
197 : public GridDefaultImplementation<3, 3, double, CpGridFamily>
198 {
199 friend class cpgrid::CpGridData;
200 friend class cpgrid::Entity<0>;
201 friend class cpgrid::Entity<1>;
202 friend class cpgrid::Entity<2>;
203 friend class cpgrid::Entity<3>;
204 template<int dim>
205 friend cpgrid::Entity<dim> createEntity(const CpGrid&,int,bool);
206
207 public:
208
209 // --- Typedefs ---
210
211
214
215
216 // --- Methods ---
217
218
220 CpGrid();
221
222 explicit CpGrid(MPIHelper::MPICommunicator comm);
223
224#if HAVE_ECL_INPUT
243 std::vector<std::size_t>
244 processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
245 Opm::EclipseState* ecl_state,
246 bool periodic_extension, bool turn_normals, bool clip_z,
247 bool pinchActive);
248
268 std::vector<std::size_t>
269 processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
270 Opm::EclipseState* ecl_state,
271 bool periodic_extension, bool turn_normals = false, bool clip_z = false);
272
273#endif
274
278 void processEclipseFormat(const grdecl& input_data, bool remove_ij_boundary, bool turn_normals = false);
279
281
287
288
295 void createCartesian(const std::array<int, 3>& dims,
296 const std::array<double, 3>& cellsize,
297 const std::array<int, 3>& shift = {0,0,0});
298
302 const std::array<int, 3>& logicalCartesianSize() const;
303
311 const std::vector<int>& globalCell() const;
312
314 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& currentData() const;
315
317 std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& currentData();
318
326 void getIJK(const int c, std::array<int,3>& ijk) const;
328
332 bool uniqueBoundaryIds() const;
333
336 void setUniqueBoundaryIds(bool uids);
337
338
339 // --- Dune interface below ---
340
342 // \@{
347 std::string name() const;
348
350 int maxLevel() const;
351
353 template<int codim>
354 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
356 template<int codim>
357 typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
358
360 template<int codim, PartitionIteratorType PiType>
361 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
363 template<int codim, PartitionIteratorType PiType>
364 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
365
367 template<int codim>
368 typename Traits::template Codim<codim>::LeafIterator leafbegin() const;
370 template<int codim>
371 typename Traits::template Codim<codim>::LeafIterator leafend() const;
372
374 template<int codim, PartitionIteratorType PiType>
375 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const;
377 template<int codim, PartitionIteratorType PiType>
378 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const;
379
381 int size (int level, int codim) const;
382
384 int size (int codim) const;
385
387 int size (int level, GeometryType type) const;
388
390 int size (GeometryType type) const;
391
393 const Traits::GlobalIdSet& globalIdSet() const;
394
396 const Traits::LocalIdSet& localIdSet() const;
397
399 const Traits::LevelIndexSet& levelIndexSet(int level) const;
400
402 const Traits::LeafIndexSet& leafIndexSet() const;
403
407 void globalRefine (int refCount);
408
409 const std::vector<Dune::GeometryType>& geomTypes(const int) const;
410
412 template <int codim>
414
431 void addLgrsUpdateLeafView(const std::vector<std::array<int,3>>& cells_per_dim_vec,
432 const std::vector<std::array<int,3>>& startIJK_vec,
433 const std::vector<std::array<int,3>>& endIJK_vec,
434 const std::vector<std::string>& lgr_name_vec);
435
436 // @brief TO BE DONE
437 const std::map<std::string,int>& getLgrNameToLevel() const;
438
439 // @breif Compute center of an entity/element/cell in the Eclipse way:
440 // - Average of the 4 corners of the bottom face.
441 // - Average of the 4 corners of the top face.
442 // Return average of the previous computations.
443 // @param [in] int Index of a cell.
444 // @return 'eclipse centroid'
445 std::array<double,3> getEclCentroid(const int& idx) const;
446
447 // @breif Compute center of an entity/element/cell in the Eclipse way:
448 // - Average of the 4 corners of the bottom face.
449 // - Average of the 4 corners of the top face.
450 // Return average of the previous computations.
451 // @param [in] Entity<0> Entity
452 // @return 'eclipse centroid'
453 std::array<double,3> getEclCentroid(const cpgrid::Entity<0>& elem) const;
454
455 // @brief Return parent (coarse) intersection (face) of a refined face on the leaf grid view, whose neighboring cells
456 // are two: one coarse cell (equivalent to its origin cell from level 0), and one refined cell
457 // from certain LGR
459
476 bool mark(int refCount, const cpgrid::Entity<0>& element);
477
481 int getMark(const cpgrid::Entity<0>& element) const;
482
485 bool preAdapt();
486
489 bool adapt();
490
504 bool adapt(const std::vector<std::array<int,3>>& cells_per_dim_vec,
505 const std::vector<int>& assignRefinedLevel,
506 const std::vector<std::string>& lgr_name_vec,
507 const std::vector<std::array<int,3>>& startIJK_vec = std::vector<std::array<int,3>>{},
508 const std::vector<std::array<int,3>>& endIJK_vec = std::vector<std::array<int,3>>{});
509
511 void postAdapt();
513
514 private:
515
517
576 void refineAndProvideMarkedRefinedRelations(/* Marked elements parameters */
577 std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
578 int& markedElem_count,
579 std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
580 std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
581 std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
582 /* Refined cells parameters */
583 std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCell_to_refinedLevelAdRefinedCell,
584 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
585 std::vector<int>& refined_cell_count_vec,
586 const std::vector<int>& assignRefinedLevel,
587 std::vector<std::vector<std::tuple<int,std::vector<int>>>>& preAdapt_parent_to_children_cells_vec,
588 /* Adapted cells parameters */
589 std::map<std::array<int,2>,int>& elemLgrAndElemLgrCell_to_adaptedCell,
590 std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
591 int& cell_count,
592 std::vector<std::vector<int>>& preAdapt_level_to_leaf_cells_vec,
593 /* Additional parameters */
594 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
595
617 std::tuple< std::vector<std::vector<std::array<int,2>>>,
618 std::vector<std::vector<int>>,
619 std::vector<std::array<int,2>>,
620 std::vector<int>> defineChildToParentAndIdxInParentCell( const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
621 const std::vector<int>& refined_cell_count_vec,
622 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
623 const int& cell_count) const;
624
651 std::pair<std::vector<std::vector<int>>, std::vector<std::array<int,2>>>
652 defineLevelToLeafAndLeafToLevelCells(const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCell_to_refinedLevelAndRefinedCell,
653 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
654 const std::vector<int>& refined_cell_count_vec,
655 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCell_to_adaptedCell,
656 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
657 const int& cell_count) const;
658
685 void identifyRefinedCornersPerLevel(std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
686 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner,
687 std::vector<int>& refined_corner_count_vec,
688 std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
689 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
690 const std::vector<int>& assignRefinedLevel,
691 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
692 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
693 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
694
712 void identifyRefinedFacesPerLevel(std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
713 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
714 std::vector<int>& refined_face_count_vec,
715 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
716 const std::vector<int>& assignRefinedLevel,
717 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
718 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
719
740 void identifyLeafGridCorners(std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
741 std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
742 int& corner_count,
743 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
744 const std::vector<int>& assignRefinedLevel,
745 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
746 std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
747 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
748 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
749
767 void identifyLeafGridFaces(std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
768 std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
769 int& face_count,
770 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
771 const std::vector<int>& assignRefinedLevel,
772 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
773 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
774
776 void populateRefinedCorners(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>>& refined_corners_vec,
777 const std::vector<int>& refined_corner_count_vec,
778 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
779 const int& preAdaptMaxLevel,
780 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner) const;
781
783 void populateRefinedFaces(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>>& refined_faces_vec,
784 std::vector<Dune::cpgrid::EntityVariableBase<enum face_tag>>& mutable_refined_face_tags_vec,
785 std::vector<Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>>& mutable_refine_face_normals_vec,
786 std::vector<Opm::SparseTable<int>>& refined_face_to_point_vec,
787 const std::vector<int>& refined_face_count_vec,
788 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
789 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
790 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
791 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
792 const int& preAdaptMaxLevel,
793 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
794 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner) const;
795
797 void populateRefinedCells(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>>& refined_cells_vec,
798 std::vector<std::vector<std::array<int,8>>>& refined_cell_to_point_vec,
799 std::vector<std::vector<int>>& refined_global_cell_vec,
800 const std::vector<int>& refined_cell_count_vec,
801 std::vector<cpgrid::OrientedEntityTable<0,1>>& refined_cell_to_face_vec,
802 std::vector<cpgrid::OrientedEntityTable<1,0>>& refined_face_to_cell_vec,
803 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
804 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
805 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
806 const std::vector<Dune::cpgrid::DefaultGeometryPolicy>& refined_geometries_vec,
807 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
808 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
809 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
810 const std::vector<int>& assignRefinedLevel,
811 const int& preAdaptMaxLevel,
812 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
813 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
814 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
815
817 void setRefinedLevelGridsGeometries( /* Refined corner arguments */
818 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>>& refined_corners_vec,
819 const std::vector<int>& refined_corner_count_vec,
820 /* Refined face arguments */
821 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>>& refined_faces_vec,
822 std::vector<Dune::cpgrid::EntityVariableBase<enum face_tag>>& mutable_refined_face_tags_vec,
823 std::vector<Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>>& mutable_refine_face_normals_vec,
824 std::vector<Opm::SparseTable<int>>& refined_face_to_point_vec,
825 const std::vector<int>& refined_face_count_vec,
826 /* Refined cell argumets */
827 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>>& refined_cells_vec,
828 std::vector<std::vector<std::array<int,8>>>& refined_cell_to_point_vec,
829 std::vector<std::vector<int>>& refined_global_cell_vec,
830 const std::vector<int>& refined_cell_count_vec,
831 std::vector<cpgrid::OrientedEntityTable<0,1>>& refined_cell_to_face_vec,
832 std::vector<cpgrid::OrientedEntityTable<1,0>>& refined_face_to_cell_vec,
833 /* Auxiliary arguments */
834 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
835 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
836 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner,
837 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
838 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
839 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
840 const std::vector<Dune::cpgrid::DefaultGeometryPolicy>& refined_geometries_vec,
841 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
842 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
843 const std::vector<int>& assignRefinedLevel,
844 const int& preAdaptMaxLevel,
845 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
846 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
847 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
848
850 void populateLeafGridCorners(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>& adapted_corners,
851 const int& corners_count,
852 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
853 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner) const;
854
856 void populateLeafGridFaces(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>& adapted_faces,
858 Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>& mutable_face_normals,
859 Opm::SparseTable<int>& adapted_face_to_point,
860 const int& face_count,
861 const std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
862 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
863 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
864 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
865 const std::vector<int>& assignRefinedLevel,
866 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
867 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
868 const std::vector<std::array<int,3>>& cells_per_dim_vec,
869 const int& preAdaptMaxLevel) const;
870
872 void populateLeafGridCells(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>& adapted_cells,
873 std::vector<std::array<int,8>>& adapted_cell_to_point,
874 const int& cell_count,
875 cpgrid::OrientedEntityTable<0,1>& adapted_cell_to_face,
876 cpgrid::OrientedEntityTable<1,0>& adapted_face_to_cell,
877 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
878 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
879 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
880 const Dune::cpgrid::DefaultGeometryPolicy& adapted_geometries,
881 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
882 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
883 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
884 const std::vector<int>& assignRefinedLevel,
885 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
886 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
887 const std::vector<std::array<int,3>>& cells_per_dim_vec,
888 const int& preAdaptMaxLevel) const;
889
891 void updateLeafGridViewGeometries( /* Leaf grid View Corners arguments */
892 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>& adapted_corners,
893 const int& corner_count,
894 /* Leaf grid View Faces arguments */
895 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>& adapted_faces,
897 Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>& mutable_face_normals,
898 Opm::SparseTable<int>& adapted_face_to_point,
899 const int& face_count,
900 /* Leaf grid View Cells argumemts */
901 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>& adapted_cells,
902 std::vector<std::array<int,8>>& adapted_cell_to_point,
903 const int& cell_count,
904 cpgrid::OrientedEntityTable<0,1>& adapted_cell_to_face,
905 cpgrid::OrientedEntityTable<1,0>& adapted_face_to_cell,
906 /* Auxiliary arguments */
907 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
908 const std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
909 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
910 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
911 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
912 const Dune::cpgrid::DefaultGeometryPolicy& adapted_geometries,
913 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
914 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
915 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
916 const std::vector<int>& assignRefinedLevel,
917 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
918 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
919 const std::vector<std::array<int,3>>& cells_per_dim_vec,
920 const int& preAdaptMaxLevel) const;
921
922 void updateCornerHistoryLevels(const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
923 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
924 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
925 const int& corner_count,
926 const std::vector<std::array<int,2>>& preAdaptGrid_corner_history,
927 const int& preAdaptMaxLevel,
928 const int& newLevels);
929
930 void globalIdsPartitionTypesLgrAndLeafGrids(const std::vector<int>& assignRefinedLevel,
931 const std::vector<std::array<int,3>>& cells_per_dim_vec,
932 const std::vector<int>& lgr_with_at_least_one_active_cell);
933
940 void getFirstChildGlobalIds([[maybe_unused]] std::vector<int>& parentToFirstChildGlobalIds);
941 public:
948 private:
949
962 void computeGlobalCellLgr(const int& level, const std::array<int,3>& startIJK, std::vector<int>& global_cell_lgr);
963
967 void computeGlobalCellLeafGridViewWithLgrs(std::vector<int>& global_cell_leaf);
968
978 std::array<int,3> getRefinedCornerIJK(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
979
994 std::array<int,3> getRefinedFaceIJK(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
995 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
996
1001 bool isRefinedCornerInInteriorLgr(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1002
1008 bool isRefinedFaceInInteriorLgr(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1009 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
1010
1016 bool isRefinedNewBornCornerOnLgrBoundary(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1017
1023 bool newRefinedCornerLiesOnEdge(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1024
1030 bool isRefinedFaceOnLgrBoundary(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1031 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
1032
1038 int getParentFaceWhereNewRefinedCornerLiesOn(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr, int elemLgr) const;
1039
1045 std::array<int,2> getParentFacesAssocWithNewRefinedCornLyingOnEdge(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr, int elemLgr) const;
1046
1047 protected:
1054 int replaceLgr1CornerIdxByLgr2CornerIdx(const std::array<int,3>& cells_per_dim_lgr1, int cornerIdxLgr1, const std::array<int,3>& cells_per_dim_lgr2) const;
1055
1063 int replaceLgr1CornerIdxByLgr2CornerIdx(const std::array<int,3>& cells_per_dim_lgr1, int cornerIdxLgr1, int elemLgr1, int parentFaceLastAppearanceIdx,
1064 const std::array<int,3>& cells_per_dim_lgr2) const;
1065
1073 int replaceLgr1FaceIdxByLgr2FaceIdx(const std::array<int,3>& cells_per_dim_lgr1, int faceIdxInLgr1,
1074 const std::shared_ptr<cpgrid::CpGridData>& elemLgr1_ptr,
1075 const std::array<int,3>& cells_per_dim_lgr2) const;
1076 private:
1077
1084 int getParentFaceWhereNewRefinedFaceLiesOn(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1085 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr,
1086 int elemLgr) const;
1087
1089
1097 bool nonNNCsSelectedCellsLGR( const std::vector<std::array<int,3>>& startIJK_vec,
1098 const std::vector<std::array<int,3>>& endIJK_vec) const;
1099
1110 void detectActiveLgrs(const std::vector<std::array<int,3>>& startIJK_vec,
1111 const std::vector<std::array<int,3>>& endIJK_vec,
1112 std::vector<int>& lgr_with_at_least_one_active_cell);
1113
1127 void markElemAssignLevelDetectActiveLgrs(const std::vector<std::array<int,3>>& startIJK_vec,
1128 const std::vector<std::array<int,3>>& endIJK_vec,
1129 std::vector<int>& assignRefinedLevel,
1130 std::vector<int>& lgr_with_at_least_one_active_cell);
1131
1137 template<class T>
1138 void computeOnLgrParents(const std::vector<std::array<int,3>>& startIJK_vec,
1139 const std::vector<std::array<int,3>>& endIJK_vec,
1140 T func);
1141
1155 void predictMinCellAndPointGlobalIdPerProcess(const std::vector<int>& assignRefinedLevel,
1156 const std::vector<std::array<int,3>>& cells_per_dim_vec,
1157 const std::vector<int>& lgr_with_at_least_one_active_cell,
1158 int& min_globalId_cell_in_proc,
1159 int& min_globalId_point_in_proc) const;
1160
1169 void assignCellIdsAndCandidatePointIds( std::vector<std::vector<int>>& localToGlobal_cells_per_level,
1170 std::vector<std::vector<int>>& localToGlobal_points_per_level,
1171 int min_globalId_cell_in_proc,
1172 int min_globalId_point_in_proc,
1173 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
1174
1188 void selectWinnerPointIds(std::vector<std::vector<int>>& localToGlobal_points_per_level,
1189 const std::vector<std::tuple<int,std::vector<int>>>& parent_to_children,
1190 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
1191
1193 void populateCellIndexSetRefinedGrid(int level);
1194
1196 void populateCellIndexSetLeafGridView();
1197
1199 void populateLeafGlobalIdSet();
1200
1201 public:
1202
1207 std::vector<std::unordered_map<std::size_t, std::size_t>> mapLocalCartesianIndexSetsToLeafIndexSet() const;
1208
1210 std::vector<std::array<int,2>> mapLeafIndexSetToLocalCartesianIndexSets() const;
1211
1213 unsigned int overlapSize(int) const;
1214
1215
1217 unsigned int ghostSize(int) const;
1218
1220 unsigned int overlapSize(int, int) const;
1221
1223 unsigned int ghostSize(int, int) const;
1224
1226 unsigned int numBoundarySegments() const;
1227
1228 void setPartitioningParams(const std::map<std::string,std::string>& params);
1229
1230 // loadbalance is not part of the grid interface therefore we skip it.
1231
1236 bool loadBalance(int overlapLayers=1,
1237 int partitionMethod = Dune::PartitionMethod::zoltan,
1238 double imbalanceTol = 1.1,
1239 int level =-1)
1240 {
1241 using std::get;
1242 return get<0>(scatterGrid(defaultTransEdgeWgt /*method*/,
1243 false /*ownersFirst*/,
1244 nullptr /*wells*/,
1245 {} /*possibleFutureConnections*/,
1246 false /*serialPartitioning*/,
1247 nullptr /*transmissibilities*/,
1248 true /*addCornerCells*/,
1249 overlapLayers,
1250 partitionMethod,
1251 imbalanceTol,
1252 false /*allowDistributedWells*/,
1253 {} /*input_cell_part*/,
1254 level));
1255 }
1256
1257 // loadbalance is not part of the grid interface therefore we skip it.
1258
1264 bool loadBalanceSerial(int overlapLayers=1,
1265 int partitionMethod = Dune::PartitionMethod::zoltan,
1266 int edgeWeightMethod = Dune::EdgeWeightMethod::defaultTransEdgeWgt,
1267 double imbalanceTol = 1.1)
1268 {
1269 using std::get;
1270 return get<0>(scatterGrid(EdgeWeightMethod(edgeWeightMethod), false, nullptr, {}, true /*serial partitioning*/, nullptr, true, overlapLayers, partitionMethod, imbalanceTol));
1271 }
1272
1273 // loadbalance is not part of the grid interface therefore we skip it.
1274
1300 std::pair<bool,std::vector<std::pair<std::string,bool>>>
1301 loadBalance(const std::vector<cpgrid::OpmWellType> * wells,
1302 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1303 const double* transmissibilities = nullptr,
1304 int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltanGoG)
1305 {
1306 return scatterGrid(defaultTransEdgeWgt, false, wells, possibleFutureConnections, false, transmissibilities, false, overlapLayers, partitionMethod);
1307 }
1308
1309 // loadbalance is not part of the grid interface therefore we skip it.
1310
1340 std::pair<bool,std::vector<std::pair<std::string,bool>>>
1341 loadBalance(EdgeWeightMethod method, const std::vector<cpgrid::OpmWellType> * wells,
1342 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1343 const double* transmissibilities = nullptr, bool ownersFirst=false,
1344 bool addCornerCells=false, int overlapLayers=1,
1345 int partitionMethod = Dune::PartitionMethod::zoltanGoG,
1346 double imbalanceTol = 1.1)
1347 {
1348 return scatterGrid(method, ownersFirst, wells, possibleFutureConnections, false, transmissibilities, addCornerCells, overlapLayers, partitionMethod, imbalanceTol);
1349 }
1350
1376 template<class DataHandle>
1377 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1378 loadBalance(DataHandle& data,
1379 const std::vector<cpgrid::OpmWellType> * wells,
1380 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1381 const double* transmissibilities = nullptr,
1382 int overlapLayers=1, int partitionMethod = 1)
1383 {
1384 auto ret = loadBalance(wells, possibleFutureConnections, transmissibilities, overlapLayers, partitionMethod);
1385 using std::get;
1386 if (get<0>(ret))
1387 {
1388 scatterData(data);
1389 }
1390 return ret;
1391 }
1392
1427 template<class DataHandle>
1428 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1429 loadBalance(DataHandle& data, EdgeWeightMethod method,
1430 const std::vector<cpgrid::OpmWellType> * wells,
1431 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1432 bool serialPartitioning,
1433 const double* transmissibilities = nullptr, bool ownersFirst=false,
1434 bool addCornerCells=false, int overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltanGoG,
1435 double imbalanceTol = 1.1,
1436 bool allowDistributedWells = false)
1437 {
1438 auto ret = scatterGrid(method, ownersFirst, wells, possibleFutureConnections, serialPartitioning, transmissibilities,
1439 addCornerCells, overlapLayers, partitionMethod, imbalanceTol, allowDistributedWells);
1440 using std::get;
1441 if (get<0>(ret))
1442 {
1443 scatterData(data);
1444 }
1445 return ret;
1446 }
1447
1464 template<class DataHandle>
1465 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1466 loadBalance(DataHandle& data, const std::vector<int>& parts,
1467 const std::vector<cpgrid::OpmWellType> * wells,
1468 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1469 bool ownersFirst=false,
1470 bool addCornerCells=false, int overlapLayers=1)
1471 {
1472 using std::get;
1473 auto ret = scatterGrid(defaultTransEdgeWgt, ownersFirst, wells,
1474 possibleFutureConnections,
1475 /* serialPartitioning = */ false,
1476 /* transmissibilities = */ {},
1477 addCornerCells, overlapLayers, /* partitionMethod =*/ Dune::PartitionMethod::simple,
1478 /* imbalanceTol (ignored) = */ 0.0,
1479 /* allowDistributedWells = */ true, parts);
1480 using std::get;
1481 if (get<0>(ret))
1482 {
1483 scatterData(data);
1484 }
1485 return ret;
1486 }
1487
1495 template<class DataHandle>
1496 bool loadBalance(DataHandle& data,
1497 decltype(data.fixedSize(0,0)) overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltan)
1498 {
1499 // decltype usage needed to tell the compiler not to use this function if first
1500 // argument is std::vector but rather loadbalance by parts
1501 bool ret = loadBalance(overlapLayers, partitionMethod);
1502 if (ret)
1503 {
1504 scatterData(data);
1505 }
1506 return ret;
1507 }
1508
1520 bool loadBalance(const std::vector<int>& parts, bool ownersFirst=false,
1521 bool addCornerCells=false, int overlapLayers=1)
1522 {
1523 using std::get;
1524 return get<0>(scatterGrid(defaultTransEdgeWgt, ownersFirst, /* wells = */ {},
1525 {},
1526 /* serialPartitioning = */ false,
1527 /* trabsmissibilities = */ {},
1528 addCornerCells, overlapLayers, /* partitionMethod =*/ Dune::PartitionMethod::simple,
1529 /* imbalanceTol (ignored) = */ 0.0,
1530 /* allowDistributedWells = */ true, parts));
1531 }
1532
1545 template<class DataHandle>
1546 bool loadBalance(DataHandle& data, const std::vector<int>& parts, bool ownersFirst=false,
1547 bool addCornerCells=false, int overlapLayers=1)
1548 {
1549 bool ret = loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
1550 if (ret)
1551 {
1552 scatterData(data);
1553 }
1554 return ret;
1555 }
1556
1569 std::vector<int>
1570 zoltanPartitionWithoutScatter(const std::vector<cpgrid::OpmWellType>* wells,
1571 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1572 const double* transmissibilities,
1573 const int numParts,
1574 const double imbalanceTol) const;
1575
1583 template<class DataHandle>
1584 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int /*level*/) const
1585 {
1586 communicate(data, iftype, dir);
1587 }
1588
1596 template<class DataHandle>
1597 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const;
1598
1600 const typename CpGridTraits::Communication& comm () const;
1602
1603 // ------------ End of Dune interface, start of simplified interface --------------
1604
1610
1611 // enum { dimension = 3 }; // already defined
1612
1613 typedef Dune::FieldVector<double, 3> Vector;
1614
1615
1616 const std::vector<double>& zcornData() const;
1617
1618
1619 // Topology
1624 int numCells(int level = -1) const;
1625
1630 int numFaces(int level = -1) const;
1631
1633 int numVertices() const;
1634
1635
1644 int numCellFaces(int cell, int level = -1) const;
1645
1652 int cellFace(int cell, int local_index, int level = -1) const;
1653
1657
1670 int faceCell(int face, int local_index, int level = -1) const;
1671
1678 int numCellFaces() const;
1679
1680 int numFaceVertices(int face) const;
1681
1686 int faceVertex(int face, int local_index) const;
1687
1690 double cellCenterDepth(int cell_index) const;
1691
1692
1693 const Vector faceCenterEcl(int cell_index, int face, const Dune::cpgrid::Intersection& intersection) const;
1694
1695 const Vector faceAreaNormalEcl(int face) const;
1696
1697
1698 // Geometry
1702 const Vector& vertexPosition(int vertex) const;
1703
1706 double faceArea(int face) const;
1707
1710 const Vector& faceCentroid(int face) const;
1711
1715 const Vector& faceNormal(int face) const;
1716
1719 double cellVolume(int cell) const;
1720
1723 const Vector& cellCentroid(int cell) const;
1724
1727 template<int codim>
1729 : public RandomAccessIteratorFacade<CentroidIterator<codim>,
1730 FieldVector<double, 3>,
1731 const FieldVector<double, 3>&, int>
1732 {
1733 public:
1735 typedef typename std::vector<cpgrid::Geometry<3-codim, 3> >::const_iterator
1740 : iter_(iter)
1741 {}
1742
1743 const FieldVector<double,3>& dereference() const
1744 {
1745 return iter_->center();
1746 }
1747 void increment()
1748 {
1749 ++iter_;
1750 }
1751 const FieldVector<double,3>& elementAt(int n)
1752 {
1753 return iter_[n]->center();
1754 }
1755 void advance(int n){
1756 iter_+=n;
1757 }
1758 void decrement()
1759 {
1760 --iter_;
1761 }
1762 int distanceTo(const CentroidIterator& o)
1763 {
1764 return o-iter_;
1765 }
1766 bool equals(const CentroidIterator& o) const{
1767 return o==iter_;
1768 }
1769 private:
1771 GeometryIterator iter_;
1772 };
1773
1775 CentroidIterator<0> beginCellCentroids() const;
1776
1778 CentroidIterator<1> beginFaceCentroids() const;
1779
1780 // Extra
1781 int boundaryId(int face) const;
1782
1789 template<class Cell2FacesRowIterator>
1790 int
1791 faceTag(const Cell2FacesRowIterator& cell_face) const;
1792
1794
1795 // ------------ End of simplified interface --------------
1796
1797 //------------- methods not in the DUNE grid interface.
1798
1803
1804
1813 template<class DataHandle>
1814 void scatterData(DataHandle& handle) const;
1815
1822 template<class DataHandle>
1823 void gatherData(DataHandle& handle) const;
1824
1826 using InterfaceMap = cpgrid::CpGridDataTraits::InterfaceMap;
1827
1857
1861
1863 void switchToGlobalView();
1864
1868
1869#if HAVE_MPI
1871 using ParallelIndexSet = cpgrid::CpGridDataTraits::ParallelIndexSet;
1873 using RemoteIndices = cpgrid::CpGridDataTraits::RemoteIndices;
1874
1876 using CommunicationType = cpgrid::CpGridDataTraits::CommunicationType;
1877
1881 const CommunicationType& cellCommunication() const;
1882
1883 ParallelIndexSet& getCellIndexSet();
1884
1885 RemoteIndices& getCellRemoteIndices();
1886
1887 const ParallelIndexSet& getCellIndexSet() const;
1888
1889 const RemoteIndices& getCellRemoteIndices() const;
1890#endif
1891
1893 const std::vector<int>& sortedNumAquiferCells() const;
1894
1895 private:
1927 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1928 scatterGrid(EdgeWeightMethod method,
1929 bool ownersFirst,
1930 const std::vector<cpgrid::OpmWellType> * wells,
1931 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1932 bool serialPartitioning,
1933 const double* transmissibilities,
1934 bool addCornerCells,
1935 int overlapLayers,
1936 int partitionMethod = Dune::PartitionMethod::zoltanGoG,
1937 double imbalanceTol = 1.1,
1938 bool allowDistributedWells = true,
1939 const std::vector<int>& input_cell_part = {},
1940 int level = -1);
1941
1946 std::vector<std::shared_ptr<cpgrid::CpGridData>> data_;
1948 cpgrid::CpGridData* current_view_data_;
1950 std::vector<std::shared_ptr<cpgrid::CpGridData>> distributed_data_;
1952 std::vector<std::shared_ptr<cpgrid::CpGridData>>* current_data_;
1954 std::map<std::string,int> lgr_names_ = {{"GLOBAL", 0}};
1960 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1961 /*
1962 * @brief Interface for scattering and gathering point data.
1963 *
1964 * @warning Will only update owner cells
1965 */
1966 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1970 std::shared_ptr<cpgrid::GlobalIdSet> global_id_set_ptr_;
1971
1972
1976 std::map<std::string,std::string> partitioningParams;
1977
1978 }; // end Class CpGrid
1979
1980} // end namespace Dune
1981
1982#include <opm/grid/cpgrid/Entity.hpp>
1983#include <opm/grid/cpgrid/Iterators.hpp>
1984#include <opm/grid/cpgrid/CpGridData.hpp>
1985
1986
1987namespace Dune
1988{
1989
1990 namespace Capabilities
1991 {
1993 template <>
1994 struct hasEntity<CpGrid, 0>
1995 {
1996 static const bool v = true;
1997 };
1998
2000 template <>
2001 struct hasEntity<CpGrid, 3>
2002 {
2003 static const bool v = true;
2004 };
2005
2006 template<>
2007 struct canCommunicate<CpGrid,0>
2008 {
2009 static const bool v = true;
2010 };
2011
2012 template<>
2013 struct canCommunicate<CpGrid,3>
2014 {
2015 static const bool v = true;
2016 };
2017
2019 template <>
2020 struct hasBackupRestoreFacilities<CpGrid>
2021 {
2022 static const bool v = false;
2023 };
2024
2025 }
2026
2027 template<class DataHandle>
2028 void CpGrid::communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
2029 {
2030 current_view_data_->communicate(data, iftype, dir);
2031 }
2032
2033
2034 template<class DataHandle>
2035 void CpGrid::scatterData(DataHandle& handle) const
2036 {
2037#if HAVE_MPI
2038 if(distributed_data_.empty())
2039 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balanced grid!");
2040 distributed_data_[0]->scatterData(handle, data_[0].get(), distributed_data_[0].get(), cellScatterGatherInterface(),
2042#else
2043 // Suppress warnings for unused argument.
2044 (void) handle;
2045#endif
2046 }
2047
2048 template<class DataHandle>
2049 void CpGrid::gatherData(DataHandle& handle) const
2050 {
2051#if HAVE_MPI
2052 if(distributed_data_.empty())
2053 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balance grid!");
2054 distributed_data_[0]->gatherData(handle, data_[0].get(), distributed_data_[0].get());
2055#else
2056 // Suppress warnings for unused argument.
2057 (void) handle;
2058#endif
2059 }
2060
2061
2062 template<class Cell2FacesRowIterator>
2063 int
2064 CpGrid::faceTag(const Cell2FacesRowIterator& cell_face) const
2065 {
2066 // Note that this relies on the following implementation detail:
2067 // The grid is always constructed such that the interior faces constructed
2068 // with orientation set to true are
2069 // oriented along the positive IJK direction. Oriented means that
2070 // the first cell attached to face has the lower index.
2071 // For faces along the boundary (only one cell, always attached at index 0)
2072 // the orientation has to be determined by the orientation of the cell.
2073 // If it is true then in UnstructuredGrid it would be stored at index 0,
2074 // otherwise at index 1.
2075 const int cell = cell_face.getCellIndex();
2076 const int face = *cell_face;
2077 assert (0 <= cell); assert (cell < numCells());
2078 assert (0 <= face); assert (face < numFaces());
2079
2081
2082 const cpgrid::EntityRep<1> f(face, true);
2083 const F2C& f2c = current_view_data_->face_to_cell_[f];
2084 const face_tag tag = current_view_data_->face_tag_[f];
2085
2086 assert ((f2c.size() == 1) || (f2c.size() == 2));
2087
2088 int inside_cell = 0;
2089
2090 if ( f2c.size() == 2 ) // Two cells => interior
2091 {
2092 if ( f2c[1].index() == cell )
2093 {
2094 inside_cell = 1;
2095 }
2096 }
2097 const bool normal_is_in = ! f2c[inside_cell].orientation();
2098
2099 switch (tag) {
2100 case I_FACE:
2101 // LEFT : RIGHT
2102 return normal_is_in ? 0 : 1; // min(I) : max(I)
2103 case J_FACE:
2104 // BACK : FRONT
2105 return normal_is_in ? 2 : 3; // min(J) : max(J)
2106 case K_FACE:
2107 // Note: TOP at min(K) as 'z' measures *depth*.
2108 // TOP : BOTTOM
2109 return normal_is_in ? 4 : 5; // min(K) : max(K)
2110 case NNC_FACE:
2111 // For nnc faces we return the otherwise unused value -1.
2112 return -1;
2113 default:
2114 OPM_THROW(std::logic_error, "Unhandled face tag. This should never happen!");
2115 }
2116 }
2117
2118 template<int dim>
2119 cpgrid::Entity<dim> createEntity(const CpGrid&, int, bool);
2120
2121} // namespace Dune
2122
2123#include <opm/grid/cpgrid/PersistentContainer.hpp>
2124#include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
2125#include "cpgrid/Intersection.hpp"
2126#include "cpgrid/Geometry.hpp"
2127#include "cpgrid/Indexsets.hpp"
2128
2129#endif // OPM_CPGRID_HEADER
An iterator over the centroids of the geometry of the entities.
Definition CpGrid.hpp:1732
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition CpGrid.hpp:1739
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition CpGrid.hpp:1736
[ provides Dune::Grid ]
Definition CpGrid.hpp:198
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, const std::vector< int > &parts, const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections={}, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1466
std::string name() const
Get the grid name.
Definition CpGrid.cpp:759
std::vector< std::unordered_map< std::size_t, std::size_t > > mapLocalCartesianIndexSetsToLeafIndexSet() const
Compute for each level grid, a map from the global_cell_[ cell index in level grid ] to the leaf inde...
Definition CpGrid.cpp:724
void switchToGlobalView()
Switch to the global view.
Definition CpGrid.cpp:1860
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition CpGrid.cpp:1112
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition CpGrid.hpp:213
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:2049
int numFaces(int level=-1) const
Get the number of faces.
Definition CpGrid.cpp:1154
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition CpGrid.cpp:1016
std::vector< std::array< int, 2 > > mapLeafIndexSetToLocalCartesianIndexSets() const
Reverse map: from leaf index cell to { level, local/level Cartesian index of the cell }.
Definition CpGrid.cpp:734
const std::vector< std::shared_ptr< Dune::cpgrid::CpGridData > > & currentData() const
Returns either data_ or distributed_data_(if non empty).
Definition CpGrid.cpp:662
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition CpGrid.cpp:1794
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition CpGrid.cpp:1011
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections, bool serialPartitioning, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltanGoG, double imbalanceTol=1.1, bool allowDistributedWells=false)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1429
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition CpGrid.hpp:2064
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition CpGrid.cpp:754
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:592
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:1930
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition CpGrid.cpp:1090
cpgrid::CpGridDataTraits::InterfaceMap InterfaceMap
The type of the map describing communication interfaces.
Definition CpGrid.hpp:1826
bool loadBalanceSerial(int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltan, int edgeWeightMethod=Dune::EdgeWeightMethod::defaultTransEdgeWgt, double imbalanceTol=1.1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1264
int replaceLgr1CornerIdxByLgr2CornerIdx(const std::array< int, 3 > &cells_per_dim_lgr1, int cornerIdxLgr1, const std::array< int, 3 > &cells_per_dim_lgr2) const
A refined corner appears in two single-cell-refinements.
Definition CpGrid.cpp:4522
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:1779
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition CpGrid.cpp:1217
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition CpGrid.cpp:672
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition CpGrid.cpp:1789
bool loadBalance(int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltan, double imbalanceTol=1.1, int level=-1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1236
int maxLevel() const
Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview),...
Definition CpGrid.cpp:764
const CpGridTraits::Communication & comm() const
Get the collective communication object.
Definition CpGrid.cpp:1137
double faceArea(int face) const
Get the area of a face.
Definition CpGrid.cpp:1769
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition CpGrid.cpp:1774
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition CpGrid.cpp:1021
void postAdapt()
Clean up refinement markers - set every element to the mark 0 which represents 'doing nothing'.
Definition CpGrid.cpp:2542
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition CpGrid.cpp:987
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition CpGrid.cpp:655
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition CpGrid.cpp:1227
int numVertices() const
Get The number of vertices.
Definition CpGrid.cpp:1160
Dune::cpgrid::Intersection getParentIntersectionFromLgrBoundaryFace(const Dune::cpgrid::Intersection &intersection) const
Definition CpGrid.cpp:1232
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:744
void syncDistributedGlobalCellIds()
Synchronizes cell global ids across processes after load balancing.
Definition CpGrid.cpp:2702
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition CpGrid.cpp:749
int faceCell(int face, int local_index, int level=-1) const
Get the index identifying a cell attached to a face.
Definition CpGrid.cpp:1184
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:2760
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(const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections={}, const double *transmissibilities=nullptr, int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltanGoG)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1301
int replaceLgr1FaceIdxByLgr2FaceIdx(const std::array< int, 3 > &cells_per_dim_lgr1, int faceIdxInLgr1, const std::shared_ptr< cpgrid::CpGridData > &elemLgr1_ptr, const std::array< int, 3 > &cells_per_dim_lgr2) const
A new refined face lays on the boudndary of a single-cell-refinement appears in at most two single-ce...
Definition CpGrid.cpp:4630
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim and PartitionIteratorType.
bool mark(int refCount, const cpgrid::Entity< 0 > &element)
------------— Adaptivity (begin) ------------—
Definition CpGrid.cpp:1956
bool loadBalance(DataHandle &data, decltype(data.fixedSize(0, 0)) overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltan)
Distributes this grid and data over the available nodes in a distributed machine.
Definition CpGrid.hpp:1496
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition CpGrid.cpp:1096
cpgrid::Entity< codim > entity(const cpgrid::Entity< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition CpGrid.cpp:1080
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition CpGrid.cpp:1028
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:1855
CpGrid()
Default constructor.
Definition CpGrid.cpp:166
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:1804
bool adapt()
Triggers the grid refinement process.
Definition CpGrid.cpp:1995
void switchToDistributedView()
Switch to the distributed view.
Definition CpGrid.cpp:1866
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition CpGrid.cpp:1799
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections={}, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltanGoG, double imbalanceTol=1.1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1341
bool preAdapt()
Set mightVanish flags for elements that will be refined in the next adapt() call Need to be called af...
Definition CpGrid.cpp:1982
int getMark(const cpgrid::Entity< 0 > &element) const
Return refinement mark for entity.
Definition CpGrid.cpp:1977
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:1546
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:1179
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition CpGrid.cpp:1850
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:1520
std::vector< int > zoltanPartitionWithoutScatter(const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections, const double *transmissibilities, const int numParts, const double imbalanceTol) const
Partitions the grid using Zoltan without decomposing and distributing it among processes.
Definition CpGrid.cpp:195
int numCells(int level=-1) const
Get the number of cells.
Definition CpGrid.cpp:1148
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition CpGrid.cpp:1764
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition CpGrid.cpp:1622
void globalRefine(int refCount)
Refine the grid refCount times using the default refinement rule.
Definition CpGrid.cpp:1033
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition CpGrid.hpp:1584
int cellFace(int cell, int local_index, int level=-1) const
Get a specific face of a cell.
Definition CpGrid.cpp:1172
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:1784
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition CpGrid.hpp:2035
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance(DataHandle &data, const std::vector< cpgrid::OpmWellType > *wells, const std::unordered_map< std::string, std::set< int > > &possibleFutureConnections={}, const double *transmissibilities=nullptr, int overlapLayers=1, int partitionMethod=1)
Distributes this grid and data over the available nodes in a distributed machine.
Definition CpGrid.hpp:1378
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:138
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition CpGridData.hpp:1021
Definition DefaultGeometryPolicy.hpp:53
Represents an entity of a given codim, with positive or negative orientation.
Definition EntityRep.hpp:99
Base class for EntityVariable and SignedEntityVariable.
Definition EntityRep.hpp:219
Definition Entity.hpp:71
This class encapsulates geometry for vertices, intersections, and cells.
Definition Geometry.hpp:75
The global id set for Dune.
Definition Indexsets.hpp:450
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Indexsets.hpp:56
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
A SparseTable stores a table with rows of varying size as efficiently as possible.
Definition SparseTable.hpp:55
The namespace Dune is the main namespace for all Dune code.
Definition CartesianIndexMapper.hpp:10
@ simple
Use simple approach based on rectangular partitioning the underlying cartesian grid.
Definition GridEnums.hpp:46
@ zoltanGoG
use Zoltan on GraphOfGrid for partitioning
Definition GridEnums.hpp:52
@ zoltan
Use Zoltan for partitioning.
Definition GridEnums.hpp:48
EdgeWeightMethod
enum for choosing Methods for weighting graph-edges correspoding to cell interfaces in Zoltan's or Me...
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:185
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:139
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition CpGrid.hpp:141
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition CpGrid.hpp:143
Traits associated with a specific codim.
Definition CpGrid.hpp:115
cpgrid::Entity< cd > Entity
The type of the entity.
Definition CpGrid.hpp:124
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition CpGrid.hpp:118
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition CpGrid.hpp:121
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition CpGrid.hpp:130
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition CpGrid.hpp:127
cpgrid::Entity< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition CpGrid.hpp:133
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:151
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:155
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:153
Definition CpGrid.hpp:95
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition CpGrid.hpp:165
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition CpGrid.hpp:104
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:162
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:160
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition CpGrid.hpp:169
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition CpGrid.hpp:106
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition CpGrid.hpp:167
GlobalIdSet LocalIdSet
The type of the local id set.
Definition CpGrid.hpp:171
cpgrid::CpGridDataTraits::Communication Communication
The type of the collective communication.
Definition CpGrid.hpp:174
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition CpGrid.hpp:109
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition CpGrid.hpp:102
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition CpGrid.hpp:100
CpGrid Grid
The type that implements the grid.
Definition CpGrid.hpp:97
Raw corner-point specification of a particular geological model.
Definition preprocess.h:56