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 <iostream>
55#if ! HAVE_MPI
56#include <list>
57#endif
58
59namespace Opm
60{
61struct NNCdata;
62class EclipseGrid;
63class EclipseState;
64}
65
66namespace Dune
67{
68 class CpGrid;
69
70 namespace cpgrid
71 {
72 class CpGridData;
73 template <int> class Entity;
74 template<int,int> class Geometry;
75 class HierarchicIterator;
76 class IntersectionIterator;
77 template<int, PartitionIteratorType> class Iterator;
78 class LevelGlobalIdSet;
79 class GlobalIdSet;
80 class Intersection;
81 class IntersectionIterator;
82 class IndexSet;
83 class IdSet;
84
85 }
86}
87
88namespace Dune
89{
90
92 //
93 // CpGridTraits
94 //
96
98 {
100 typedef CpGrid Grid;
101
110
113
116 template <int cd>
117 struct Codim
118 {
121 typedef cpgrid::Geometry<3-cd, 3> Geometry;
122 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> Geometry;
125 //typedef Dune::Geometry<3-cd, 3, CpGrid, cpgrid::Geometry> LocalGeometry;
128
131
134
137
140 template <PartitionIteratorType pitype>
148 };
149
152 template <PartitionIteratorType pitype>
154 {
156 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> > LevelGridView;
158 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> > LeafGridView;
159
160 };
161
163 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid>> LevelGridView;
165 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid>> LeafGridView;
166
175
177 using Communication = cpgrid::CpGridDataTraits::Communication;
178 using CollectiveCommunication = cpgrid::CpGridDataTraits::CollectiveCommunication;
179 };
180
182 //
183 // CpGridFamily
184 //
186
188 {
189 typedef CpGridTraits Traits;
190 };
191
193 //
194 // CpGrid
195 //
197
199 class CpGrid
200 : public GridDefaultImplementation<3, 3, double, CpGridFamily>
201 {
202 friend class cpgrid::CpGridData;
203 friend class cpgrid::Entity<0>;
204 friend class cpgrid::Entity<1>;
205 friend class cpgrid::Entity<2>;
206 friend class cpgrid::Entity<3>;
207 template<int dim>
208 friend cpgrid::Entity<dim> createEntity(const CpGrid&,int,bool);
209
210 public:
211
212 // --- Typedefs ---
213
214
217
218
219 // --- Methods ---
220
221
223 CpGrid();
224
225 CpGrid(MPIHelper::MPICommunicator comm);
226
227#if HAVE_ECL_INPUT
246 std::vector<std::size_t>
247 processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
248 Opm::EclipseState* ecl_state,
249 bool periodic_extension, bool turn_normals, bool clip_z,
250 bool pinchActive);
251
271 std::vector<std::size_t>
272 processEclipseFormat(const Opm::EclipseGrid* ecl_grid,
273 Opm::EclipseState* ecl_state,
274 bool periodic_extension, bool turn_normals = false, bool clip_z = false);
275
276#endif
277
281 void processEclipseFormat(const grdecl& input_data, bool remove_ij_boundary, bool turn_normals = false);
282
284
290
291
298 void createCartesian(const std::array<int, 3>& dims,
299 const std::array<double, 3>& cellsize,
300 const std::array<int, 3>& shift = {0,0,0});
301
305 const std::array<int, 3>& logicalCartesianSize() const;
306
314 const std::vector<int>& globalCell() const;
315
317 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& currentData() const;
318
320 std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& currentData();
321
329 void getIJK(const int c, std::array<int,3>& ijk) const;
331
335 bool uniqueBoundaryIds() const;
336
339 void setUniqueBoundaryIds(bool uids);
340
341
342 // --- Dune interface below ---
343
345 // \@{
350 std::string name() const;
351
353 int maxLevel() const;
354
356 template<int codim>
357 typename Traits::template Codim<codim>::LevelIterator lbegin (int level) const;
359 template<int codim>
360 typename Traits::template Codim<codim>::LevelIterator lend (int level) const;
361
363 template<int codim, PartitionIteratorType PiType>
364 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lbegin (int level) const;
366 template<int codim, PartitionIteratorType PiType>
367 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator lend (int level) const;
368
370 template<int codim>
371 typename Traits::template Codim<codim>::LeafIterator leafbegin() const;
373 template<int codim>
374 typename Traits::template Codim<codim>::LeafIterator leafend() const;
375
377 template<int codim, PartitionIteratorType PiType>
378 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafbegin() const;
380 template<int codim, PartitionIteratorType PiType>
381 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator leafend() const;
382
384 int size (int level, int codim) const;
385
387 int size (int codim) const;
388
390 int size (int level, GeometryType type) const;
391
393 int size (GeometryType type) const;
394
396 const Traits::GlobalIdSet& globalIdSet() const;
397
399 const Traits::LocalIdSet& localIdSet() const;
400
402 const Traits::LevelIndexSet& levelIndexSet(int level) const;
403
405 const Traits::LeafIndexSet& leafIndexSet() const;
406
410 void globalRefine (int refCount);
411
412 const std::vector<Dune::GeometryType>& geomTypes(const int) const;
413
415 template <int codim>
417
431 void addLgrUpdateLeafView(const std::array<int,3>& cells_per_dim, const std::array<int,3>& startIJK,
432 const std::array<int,3>& endIJK, const std::string& lgr_name);
433
450 void addLgrsUpdateLeafView(const std::vector<std::array<int,3>>& cells_per_dim_vec,
451 const std::vector<std::array<int,3>>& startIJK_vec,
452 const std::vector<std::array<int,3>>& endIJK_vec,
453 const std::vector<std::string>& lgr_name_vec);
454
455 // @brief TO BE DONE
456 const std::map<std::string,int>& getLgrNameToLevel() const;
457
458 // @breif Compute center of an entity/element/cell in the Eclipse way:
459 // - Average of the 4 corners of the bottom face.
460 // - Average of the 4 corners of the top face.
461 // Return average of the previous computations.
462 // @param [in] int Index of a cell.
463 // @return 'eclipse centroid'
464 std::array<double,3> getEclCentroid(const int& idx) const;
465
466 // @breif Compute center of an entity/element/cell in the Eclipse way:
467 // - Average of the 4 corners of the bottom face.
468 // - Average of the 4 corners of the top face.
469 // Return average of the previous computations.
470 // @param [in] Entity<0> Entity
471 // @return 'eclipse centroid'
472 std::array<double,3> getEclCentroid(const cpgrid::Entity<0>& elem) const;
473
474 // @brief Return parent (coarse) intersection (face) of a refined face on the leaf grid view, whose neighboring cells
475 // are two: one coarse cell (equivalent to its origin cell from level 0), and one refined cell
476 // from certain LGR
478
479
496 bool mark(int refCount, const cpgrid::Entity<0>& element);
497
501 int getMark(const cpgrid::Entity<0>& element) const;
502
505 bool preAdapt();
506
508 bool adapt();
509
527 bool adapt(const std::vector<std::array<int,3>>& cells_per_dim_vec,
528 const std::vector<int>& assignRefinedLevel,
529 const std::vector<std::string>& lgr_name_vec,
530 bool isCARFIN = false,
531 const std::vector<std::array<int,3>>& startIJK_vec = std::vector<std::array<int,3>>{},
532 const std::vector<std::array<int,3>>& endIJK_vec = std::vector<std::array<int,3>>{});
533
535 void postAdapt();
537
538 private:
539
541
600 void refineAndProvideMarkedRefinedRelations(/* Marked elements parameters */
601 std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
602 int& markedElem_count,
603 std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
604 std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
605 std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
606 /* Refined cells parameters */
607 std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCell_to_refinedLevelAdRefinedCell,
608 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
609 std::vector<int>& refined_cell_count_vec,
610 const std::vector<int>& assignRefinedLevel,
611 std::vector<std::vector<std::tuple<int,std::vector<int>>>>& preAdapt_parent_to_children_cells_vec,
612 /* Adapted cells parameters */
613 std::map<std::array<int,2>,int>& elemLgrAndElemLgrCell_to_adaptedCell,
614 std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
615 int& cell_count,
616 std::vector<std::vector<int>>& preAdapt_level_to_leaf_cells_vec,
617 /* Additional parameters */
618 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
619
641 std::tuple< std::vector<std::vector<std::array<int,2>>>,
642 std::vector<std::vector<int>>,
643 std::vector<std::array<int,2>>,
644 std::vector<int>> defineChildToParentAndIdxInParentCell( const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
645 const std::vector<int>& refined_cell_count_vec,
646 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
647 const int& cell_count) const;
648
675 std::pair<std::vector<std::vector<int>>, std::vector<std::array<int,2>>>
676 defineLevelToLeafAndLeafToLevelCells(const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCell_to_refinedLevelAndRefinedCell,
677 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
678 const std::vector<int>& refined_cell_count_vec,
679 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCell_to_adaptedCell,
680 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
681 const int& cell_count) const;
682
709 void identifyRefinedCornersPerLevel(std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
710 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner,
711 std::vector<int>& refined_corner_count_vec,
712 std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
713 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
714 const std::vector<int>& assignRefinedLevel,
715 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
716 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
717 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
718
736 void identifyRefinedFacesPerLevel(std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
737 std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
738 std::vector<int>& refined_face_count_vec,
739 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
740 const std::vector<int>& assignRefinedLevel,
741 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
742 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
743
764 void identifyLeafGridCorners(std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
765 std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
766 int& corner_count,
767 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
768 const std::vector<int>& assignRefinedLevel,
769 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
770 std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
771 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
772 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
773
791 void identifyLeafGridFaces(std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
792 std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
793 int& face_count,
794 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
795 const std::vector<int>& assignRefinedLevel,
796 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
797 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
798
800 void populateRefinedCorners(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>>& refined_corners_vec,
801 const std::vector<int>& refined_corner_count_vec,
802 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
803 const int& preAdaptMaxLevel,
804 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner) const;
805
807 void populateRefinedFaces(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>>& refined_faces_vec,
808 std::vector<Dune::cpgrid::EntityVariableBase<enum face_tag>>& mutable_refined_face_tags_vec,
809 std::vector<Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>>& mutable_refine_face_normals_vec,
810 std::vector<Opm::SparseTable<int>>& refined_face_to_point_vec,
811 const std::vector<int>& refined_face_count_vec,
812 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
813 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
814 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
815 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
816 const int& preAdaptMaxLevel,
817 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
818 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner) const;
819
821 void populateRefinedCells(std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>>& refined_cells_vec,
822 std::vector<std::vector<std::array<int,8>>>& refined_cell_to_point_vec,
823 std::vector<std::vector<int>>& refined_global_cell_vec,
824 const std::vector<int>& refined_cell_count_vec,
825 std::vector<cpgrid::OrientedEntityTable<0,1>>& refined_cell_to_face_vec,
826 std::vector<cpgrid::OrientedEntityTable<1,0>>& refined_face_to_cell_vec,
827 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
828 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
829 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
830 const std::vector<Dune::cpgrid::DefaultGeometryPolicy>& refined_geometries_vec,
831 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
832 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
833 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
834 const std::vector<int>& assignRefinedLevel,
835 const int& preAdaptMaxLevel,
836 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
837 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
838 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
839
841 void setRefinedLevelGridsGeometries( /* Refined corner arguments */
842 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>>& refined_corners_vec,
843 const std::vector<int>& refined_corner_count_vec,
844 /* Refined face arguments */
845 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>>& refined_faces_vec,
846 std::vector<Dune::cpgrid::EntityVariableBase<enum face_tag>>& mutable_refined_face_tags_vec,
847 std::vector<Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>>& mutable_refine_face_normals_vec,
848 std::vector<Opm::SparseTable<int>>& refined_face_to_point_vec,
849 const std::vector<int>& refined_face_count_vec,
850 /* Refined cell argumets */
851 std::vector<Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>>& refined_cells_vec,
852 std::vector<std::vector<std::array<int,8>>>& refined_cell_to_point_vec,
853 std::vector<std::vector<int>>& refined_global_cell_vec,
854 std::vector<int>& refined_cell_count_vec,
855 std::vector<cpgrid::OrientedEntityTable<0,1>>& refined_cell_to_face_vec,
856 std::vector<cpgrid::OrientedEntityTable<1,0>>& refined_face_to_cell_vec,
857 /* Auxiliary arguments */
858 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCell_to_elemLgrAndElemLgrCell,
859 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedFace_to_elemLgrAndElemLgrFace,
860 const std::map<std::array<int,2>,std::array<int,2>>& refinedLevelAndRefinedCorner_to_elemLgrAndElemLgrCorner,
861 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
862 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrFace_to_refinedLevelAndRefinedFace,
863 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
864 const std::vector<Dune::cpgrid::DefaultGeometryPolicy>& refined_geometries_vec,
865 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
866 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
867 const std::vector<int>& assignRefinedLevel,
868 const int& preAdaptMaxLevel,
869 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
870 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
871 const std::vector<std::array<int,3>>& cells_per_dim_vec) const;
872
874 void populateLeafGridCorners(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>& adapted_corners,
875 const int& corners_count,
876 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
877 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner) const;
878
880 void populateLeafGridFaces(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>& adapted_faces,
882 Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>& mutable_face_normals,
883 Opm::SparseTable<int>& adapted_face_to_point,
884 const int& face_count,
885 const std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
886 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
887 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
888 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
889 const std::vector<int>& assignRefinedLevel,
890 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
891 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
892 const std::vector<std::array<int,3>>& cells_per_dim_vec,
893 const int& preAdaptMaxLevel) const;
894
896 void populateLeafGridCells(Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>& adapted_cells,
897 std::vector<std::array<int,8>>& adapted_cell_to_point,
898 const int& cell_count,
899 cpgrid::OrientedEntityTable<0,1>& adapted_cell_to_face,
900 cpgrid::OrientedEntityTable<1,0>& adapted_face_to_cell,
901 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
902 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
903 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
904 const Dune::cpgrid::DefaultGeometryPolicy& adapted_geometries,
905 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
906 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
907 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
908 const std::vector<int>& assignRefinedLevel,
909 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
910 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
911 const std::vector<std::array<int,3>>& cells_per_dim_vec,
912 const int& preAdaptMaxLevel) const;
913
915 void updateLeafGridViewGeometries( /* Leaf grid View Corners arguments */
916 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<0,3>>& adapted_corners,
917 const int& corner_count,
918 /* Leaf grid View Faces arguments */
919 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<2,3>>& adapted_faces,
921 Dune::cpgrid::EntityVariableBase<Dune::FieldVector<double,3>>& mutable_face_normals,
922 Opm::SparseTable<int>& adapted_face_to_point,
923 const int& face_count,
924 /* Leaf grid View Cells argumemts */
925 Dune::cpgrid::EntityVariableBase<cpgrid::Geometry<3,3>>& adapted_cells,
926 std::vector<std::array<int,8>>& adapted_cell_to_point,
927 const int& cell_count,
928 cpgrid::OrientedEntityTable<0,1>& adapted_cell_to_face,
929 cpgrid::OrientedEntityTable<1,0>& adapted_face_to_cell,
930 /* Auxiliary arguments */
931 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
932 const std::unordered_map<int,std::array<int,2>>& adaptedFace_to_elemLgrAndElemLgrFace,
933 const std::unordered_map<int,std::array<int,2>>& adaptedCell_to_elemLgrAndElemLgrCell,
934 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrFace_to_adaptedFace,
935 const std::vector<std::vector<std::pair<int, std::vector<int>>>>& faceInMarkedElemAndRefinedFaces,
936 const Dune::cpgrid::DefaultGeometryPolicy& adapted_geometries,
937 const std::map<std::array<int,2>,int>& elemLgrAndElemLgrCorner_to_adaptedCorner,
938 const std::map<std::array<int,2>, std::array<int,2>>& vanishedRefinedCorner_to_itsLastAppearance,
939 const std::vector<std::shared_ptr<Dune::cpgrid::CpGridData>>& markedElem_to_itsLgr,
940 const std::vector<int>& assignRefinedLevel,
941 const std::map<std::array<int,2>,int>& markedElemAndEquivRefinedCorn_to_corner,
942 const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
943 const std::vector<std::array<int,3>>& cells_per_dim_vec,
944 const int& preAdaptMaxLevel) const;
945
946 void updateCornerHistoryLevels(const std::vector<std::vector<std::array<int,2>>>& cornerInMarkedElemWithEquivRefinedCorner,
947 const std::map<std::array<int,2>,std::array<int,2>>& elemLgrAndElemLgrCorner_to_refinedLevelAndRefinedCorner,
948 const std::unordered_map<int,std::array<int,2>>& adaptedCorner_to_elemLgrAndElemLgrCorner,
949 const int& corner_count,
950 const std::vector<std::array<int,2>>& preAdaptGrid_corner_history,
951 const int& preAdaptMaxLevel,
952 const int& newLevels);
953
954
967 void computeGlobalCellLgr(const int& level, const std::array<int,3>& startIJK, std::vector<int>& global_cell_lgr);
968
972 void computeGlobalCellLeafGridViewWithLgrs(std::vector<int>& global_cell_leaf);
973
983 std::array<int,3> getRefinedCornerIJK(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
984
999 std::array<int,3> getRefinedFaceIJK(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1000 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
1001
1006 bool isRefinedCornerInInteriorLgr(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1007
1013 bool isRefinedFaceInInteriorLgr(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1014 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
1015
1021 bool isRefinedNewBornCornerOnLgrBoundary(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1022
1028 bool newRefinedCornerLiesOnEdge(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr) const;
1029
1035 bool isRefinedFaceOnLgrBoundary(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1036 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr) const;
1037
1043 int getParentFaceWhereNewRefinedCornerLiesOn(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr, int elemLgr) const;
1044
1050 std::array<int,2> getParentFacesAssocWithNewRefinedCornLyingOnEdge(const std::array<int,3>& cells_per_dim, int cornerIdxInLgr, int elemLgr) const;
1051
1058 int replaceLgr1CornerIdxByLgr2CornerIdx(const std::array<int,3>& cells_per_dim_lgr1, int cornerIdxLgr1, const std::array<int,3>& cells_per_dim_lgr2) const;
1059
1067 int replaceLgr1CornerIdxByLgr2CornerIdx(const std::array<int,3>& cells_per_dim_lgr1, int cornerIdxLgr1, int elemLgr1, int parentFaceLastAppearanceIdx,
1068 const std::array<int,3>& cells_per_dim_lgr2) const;
1069
1077 int replaceLgr1FaceIdxByLgr2FaceIdx(const std::array<int,3>& cells_per_dim_lgr1, int faceIdxInLgr1,
1078 const std::shared_ptr<cpgrid::CpGridData>& elemLgr1_ptr,
1079 const std::array<int,3>& cells_per_dim_lgr2) const;
1080
1087 int getParentFaceWhereNewRefinedFaceLiesOn(const std::array<int,3>& cells_per_dim, int faceIdxInLgr,
1088 const std::shared_ptr<cpgrid::CpGridData>& elemLgr_ptr,
1089 int elemLgr) const;
1091
1092 public:
1093
1094
1095
1097 unsigned int overlapSize(int) const;
1098
1099
1101 unsigned int ghostSize(int) const;
1102
1104 unsigned int overlapSize(int, int) const;
1105
1107 unsigned int ghostSize(int, int) const;
1108
1110 unsigned int numBoundarySegments() const;
1111
1112 void setPartitioningParams(const std::map<std::string,std::string>& params);
1113
1114 // loadbalance is not part of the grid interface therefore we skip it.
1115
1120 bool loadBalance(int overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltan, double imbalanceTol = 1.1)
1121 {
1122 using std::get;
1123 return get<0>(scatterGrid(defaultTransEdgeWgt, false, nullptr, {}, false, nullptr, true, overlapLayers, partitionMethod, imbalanceTol));
1124 }
1125
1126 // loadbalance is not part of the grid interface therefore we skip it.
1127
1133 bool loadBalanceSerial(int overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltan, int edgeWeightMethod = Dune::EdgeWeightMethod::defaultTransEdgeWgt, double imbalanceTol = 1.1)
1134 {
1135 using std::get;
1136 return get<0>(scatterGrid(EdgeWeightMethod(edgeWeightMethod), false, nullptr, {}, true /*serial partitioning*/, nullptr, true, overlapLayers, partitionMethod, imbalanceTol));
1137 }
1138
1139 // loadbalance is not part of the grid interface therefore we skip it.
1140
1166 std::pair<bool,std::vector<std::pair<std::string,bool>>>
1167 loadBalance(const std::vector<cpgrid::OpmWellType> * wells,
1168 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1169 const double* transmissibilities = nullptr,
1170 int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltan)
1171 {
1172 return scatterGrid(defaultTransEdgeWgt, false, wells, possibleFutureConnections, false, transmissibilities, false, overlapLayers, partitionMethod);
1173 }
1174
1175 // loadbalance is not part of the grid interface therefore we skip it.
1176
1206 std::pair<bool,std::vector<std::pair<std::string,bool>>>
1207 loadBalance(EdgeWeightMethod method, const std::vector<cpgrid::OpmWellType> * wells,
1208 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1209 const double* transmissibilities = nullptr, bool ownersFirst=false,
1210 bool addCornerCells=false, int overlapLayers=1,
1211 int partitionMethod = Dune::PartitionMethod::zoltan,
1212 double imbalanceTol = 1.1)
1213 {
1214 return scatterGrid(method, ownersFirst, wells, possibleFutureConnections, false, transmissibilities, addCornerCells, overlapLayers, partitionMethod, imbalanceTol);
1215 }
1216
1242 template<class DataHandle>
1243 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1244 loadBalance(DataHandle& data,
1245 const std::vector<cpgrid::OpmWellType> * wells,
1246 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1247 const double* transmissibilities = nullptr,
1248 int overlapLayers=1, int partitionMethod = 1)
1249 {
1250 auto ret = loadBalance(wells, possibleFutureConnections, transmissibilities, overlapLayers, partitionMethod);
1251 using std::get;
1252 if (get<0>(ret))
1253 {
1254 scatterData(data);
1255 }
1256 return ret;
1257 }
1258
1293 template<class DataHandle>
1294 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1295 loadBalance(DataHandle& data, EdgeWeightMethod method,
1296 const std::vector<cpgrid::OpmWellType> * wells,
1297 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1298 bool serialPartitioning,
1299 const double* transmissibilities = nullptr, bool ownersFirst=false,
1300 bool addCornerCells=false, int overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltan,
1301 double imbalanceTol = 1.1,
1302 bool allowDistributedWells = false)
1303 {
1304 auto ret = scatterGrid(method, ownersFirst, wells, possibleFutureConnections, serialPartitioning, transmissibilities,
1305 addCornerCells, overlapLayers, partitionMethod, imbalanceTol, allowDistributedWells);
1306 using std::get;
1307 if (get<0>(ret))
1308 {
1309 scatterData(data);
1310 }
1311 return ret;
1312 }
1313
1330 template<class DataHandle>
1331 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1332 loadBalance(DataHandle& data, const std::vector<int>& parts,
1333 const std::vector<cpgrid::OpmWellType> * wells,
1334 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections = {},
1335 bool ownersFirst=false,
1336 bool addCornerCells=false, int overlapLayers=1)
1337 {
1338 using std::get;
1339 auto ret = scatterGrid(defaultTransEdgeWgt, ownersFirst, wells,
1340 possibleFutureConnections,
1341 /* serialPartitioning = */ false,
1342 /* transmissibilities = */ {},
1343 addCornerCells, overlapLayers, /* partitionMethod =*/ Dune::PartitionMethod::simple,
1344 /* imbalanceTol (ignored) = */ 0.0,
1345 /* allowDistributedWells = */ true, parts);
1346 using std::get;
1347 if (get<0>(ret))
1348 {
1349 scatterData(data);
1350 }
1351 return ret;
1352 }
1353
1361 template<class DataHandle>
1362 bool loadBalance(DataHandle& data,
1363 decltype(data.fixedSize(0,0)) overlapLayers=1, int partitionMethod = Dune::PartitionMethod::zoltan)
1364 {
1365 // decltype usage needed to tell the compiler not to use this function if first
1366 // argument is std::vector but rather loadbalance by parts
1367 bool ret = loadBalance(overlapLayers, partitionMethod);
1368 if (ret)
1369 {
1370 scatterData(data);
1371 }
1372 return ret;
1373 }
1374
1386 bool loadBalance(const std::vector<int>& parts, bool ownersFirst=false,
1387 bool addCornerCells=false, int overlapLayers=1)
1388 {
1389 using std::get;
1390 return get<0>(scatterGrid(defaultTransEdgeWgt, ownersFirst, /* wells = */ {},
1391 {},
1392 /* serialPartitioning = */ false,
1393 /* trabsmissibilities = */ {},
1394 addCornerCells, overlapLayers, /* partitionMethod =*/ Dune::PartitionMethod::simple,
1395 /* imbalanceTol (ignored) = */ 0.0,
1396 /* allowDistributedWells = */ true, parts));
1397 }
1398
1411 template<class DataHandle>
1412 bool loadBalance(DataHandle& data, const std::vector<int>& parts, bool ownersFirst=false,
1413 bool addCornerCells=false, int overlapLayers=1)
1414 {
1415 bool ret = loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
1416 if (ret)
1417 {
1418 scatterData(data);
1419 }
1420 return ret;
1421 }
1422
1435 std::vector<int>
1436 zoltanPartitionWithoutScatter(const std::vector<cpgrid::OpmWellType>* wells,
1437 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1438 const double* transmissibilities,
1439 const int numParts,
1440 const double imbalanceTol) const;
1441
1449 template<class DataHandle>
1450 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int /*level*/) const
1451 {
1452 communicate(data, iftype, dir);
1453 }
1454
1462 template<class DataHandle>
1463 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const;
1464
1466 const typename CpGridTraits::Communication& comm () const;
1468
1469 // ------------ End of Dune interface, start of simplified interface --------------
1470
1476
1477 // enum { dimension = 3 }; // already defined
1478
1479 typedef Dune::FieldVector<double, 3> Vector;
1480
1481
1482 const std::vector<double>& zcornData() const;
1483
1484
1485 // Topology
1487 int numCells() const;
1488
1490 int numFaces() const;
1491
1493 int numVertices() const;
1494
1495
1502 int numCellFaces(int cell) const;
1503
1508 int cellFace(int cell, int local_index) const;
1509
1513
1524 int faceCell(int face, int local_index) const;
1525
1532 int numCellFaces() const;
1533
1534 int numFaceVertices(int face) const;
1535
1540 int faceVertex(int face, int local_index) const;
1541
1544 double cellCenterDepth(int cell_index) const;
1545
1546
1547 const Vector faceCenterEcl(int cell_index, int face, const Dune::cpgrid::Intersection& intersection) const;
1548
1549 const Vector faceAreaNormalEcl(int face) const;
1550
1551
1552 // Geometry
1556 const Vector& vertexPosition(int vertex) const;
1557
1560 double faceArea(int face) const;
1561
1564 const Vector& faceCentroid(int face) const;
1565
1569 const Vector& faceNormal(int face) const;
1570
1573 double cellVolume(int cell) const;
1574
1577 const Vector& cellCentroid(int cell) const;
1578
1581 template<int codim>
1583 : public RandomAccessIteratorFacade<CentroidIterator<codim>,
1584 FieldVector<double, 3>,
1585 const FieldVector<double, 3>&, int>
1586 {
1587 public:
1589 typedef typename std::vector<cpgrid::Geometry<3-codim, 3> >::const_iterator
1594 : iter_(iter)
1595 {}
1596
1597 const FieldVector<double,3>& dereference() const
1598 {
1599 return iter_->center();
1600 }
1601 void increment()
1602 {
1603 ++iter_;
1604 }
1605 const FieldVector<double,3>& elementAt(int n)
1606 {
1607 return iter_[n]->center();
1608 }
1609 void advance(int n){
1610 iter_+=n;
1611 }
1612 void decrement()
1613 {
1614 --iter_;
1615 }
1616 int distanceTo(const CentroidIterator& o)
1617 {
1618 return o-iter_;
1619 }
1620 bool equals(const CentroidIterator& o) const{
1621 return o==iter_;
1622 }
1623 private:
1625 GeometryIterator iter_;
1626 };
1627
1629 CentroidIterator<0> beginCellCentroids() const;
1630
1632 CentroidIterator<1> beginFaceCentroids() const;
1633
1634 // Extra
1635 int boundaryId(int face) const;
1636
1643 template<class Cell2FacesRowIterator>
1644 int
1645 faceTag(const Cell2FacesRowIterator& cell_face) const;
1646
1648
1649 // ------------ End of simplified interface --------------
1650
1651 //------------- methods not in the DUNE grid interface.
1652
1657
1658
1667 template<class DataHandle>
1668 void scatterData(DataHandle& handle) const;
1669
1676 template<class DataHandle>
1677 void gatherData(DataHandle& handle) const;
1678
1680 using InterfaceMap = cpgrid::CpGridDataTraits::InterfaceMap;
1681
1711
1715
1717 void switchToGlobalView();
1718
1722
1723#if HAVE_MPI
1725 using ParallelIndexSet = cpgrid::CpGridDataTraits::ParallelIndexSet;
1727 using RemoteIndices = cpgrid::CpGridDataTraits::RemoteIndices;
1728
1730 using CommunicationType = cpgrid::CpGridDataTraits::CommunicationType;
1731
1735 const CommunicationType& cellCommunication() const;
1736
1737 ParallelIndexSet& getCellIndexSet();
1738
1739 RemoteIndices& getCellRemoteIndices();
1740
1741 const ParallelIndexSet& getCellIndexSet() const;
1742
1743 const RemoteIndices& getCellRemoteIndices() const;
1744#endif
1745
1747 const std::vector<int>& sortedNumAquiferCells() const;
1748
1749 private:
1781 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1782 scatterGrid(EdgeWeightMethod method,
1783 bool ownersFirst,
1784 const std::vector<cpgrid::OpmWellType> * wells,
1785 const std::unordered_map<std::string, std::set<int>>& possibleFutureConnections,
1786 bool serialPartitioning,
1787 const double* transmissibilities,
1788 bool addCornerCells,
1789 int overlapLayers,
1790 int partitionMethod = Dune::PartitionMethod::zoltan,
1791 double imbalanceTol = 1.1,
1792 bool allowDistributedWells = true,
1793 const std::vector<int>& input_cell_part = {});
1794
1799 std::vector<std::shared_ptr<cpgrid::CpGridData>> data_;
1801 cpgrid::CpGridData* current_view_data_;
1803 std::vector<std::shared_ptr<cpgrid::CpGridData>> distributed_data_;
1805 std::vector<std::shared_ptr<cpgrid::CpGridData>>* current_data_;
1807 std::map<std::string,int> lgr_names_ = {{"GLOBAL", 0}};
1813 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1814 /*
1815 * @brief Interface for scattering and gathering point data.
1816 *
1817 * @warning Will only update owner cells
1818 */
1819 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1823 std::shared_ptr<cpgrid::GlobalIdSet> global_id_set_ptr_;
1824
1825
1829 std::map<std::string,std::string> partitioningParams;
1830
1831 }; // end Class CpGrid
1832
1833} // end namespace Dune
1834
1835#include <opm/grid/cpgrid/Entity.hpp>
1836#include <opm/grid/cpgrid/Iterators.hpp>
1837#include <opm/grid/cpgrid/CpGridData.hpp>
1838
1839
1840namespace Dune
1841{
1842
1843 namespace Capabilities
1844 {
1846 template <>
1847 struct hasEntity<CpGrid, 0>
1848 {
1849 static const bool v = true;
1850 };
1851
1853 template <>
1854 struct hasEntity<CpGrid, 3>
1855 {
1856 static const bool v = true;
1857 };
1858
1859 template<>
1860 struct canCommunicate<CpGrid,0>
1861 {
1862 static const bool v = true;
1863 };
1864
1865 template<>
1866 struct canCommunicate<CpGrid,3>
1867 {
1868 static const bool v = true;
1869 };
1870
1872 template <>
1873 struct hasBackupRestoreFacilities<CpGrid>
1874 {
1875 static const bool v = false;
1876 };
1877
1878 }
1879
1880 template<class DataHandle>
1881 void CpGrid::communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir) const
1882 {
1883 current_view_data_->communicate(data, iftype, dir);
1884 }
1885
1886
1887 template<class DataHandle>
1888 void CpGrid::scatterData(DataHandle& handle) const
1889 {
1890#if HAVE_MPI
1891 if(distributed_data_.empty())
1892 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balanced grid!");
1893 distributed_data_[0]->scatterData(handle, data_[0].get(), distributed_data_[0].get(), cellScatterGatherInterface(),
1895#else
1896 // Suppress warnings for unused argument.
1897 (void) handle;
1898#endif
1899 }
1900
1901 template<class DataHandle>
1902 void CpGrid::gatherData(DataHandle& handle) const
1903 {
1904#if HAVE_MPI
1905 if(distributed_data_.empty())
1906 OPM_THROW(std::runtime_error, "Moving Data only allowed with a load balance grid!");
1907 distributed_data_[0]->gatherData(handle, data_[0].get(), distributed_data_[0].get());
1908#else
1909 // Suppress warnings for unused argument.
1910 (void) handle;
1911#endif
1912 }
1913
1914
1915 template<class Cell2FacesRowIterator>
1916 int
1917 CpGrid::faceTag(const Cell2FacesRowIterator& cell_face) const
1918 {
1919 // Note that this relies on the following implementation detail:
1920 // The grid is always constructed such that the interior faces constructed
1921 // with orientation set to true are
1922 // oriented along the positive IJK direction. Oriented means that
1923 // the first cell attached to face has the lower index.
1924 // For faces along the boundary (only one cell, always attached at index 0)
1925 // the orientation has to be determined by the orientation of the cell.
1926 // If it is true then in UnstructuredGrid it would be stored at index 0,
1927 // otherwise at index 1.
1928 const int cell = cell_face.getCellIndex();
1929 const int face = *cell_face;
1930 assert (0 <= cell); assert (cell < numCells());
1931 assert (0 <= face); assert (face < numFaces());
1932
1934
1935 const cpgrid::EntityRep<1> f(face, true);
1936 const F2C& f2c = current_view_data_->face_to_cell_[f];
1937 const face_tag tag = current_view_data_->face_tag_[f];
1938
1939 assert ((f2c.size() == 1) || (f2c.size() == 2));
1940
1941 int inside_cell = 0;
1942
1943 if ( f2c.size() == 2 ) // Two cells => interior
1944 {
1945 if ( f2c[1].index() == cell )
1946 {
1947 inside_cell = 1;
1948 }
1949 }
1950 const bool normal_is_in = ! f2c[inside_cell].orientation();
1951
1952 switch (tag) {
1953 case I_FACE:
1954 // LEFT : RIGHT
1955 return normal_is_in ? 0 : 1; // min(I) : max(I)
1956 case J_FACE:
1957 // BACK : FRONT
1958 return normal_is_in ? 2 : 3; // min(J) : max(J)
1959 case K_FACE:
1960 // Note: TOP at min(K) as 'z' measures *depth*.
1961 // TOP : BOTTOM
1962 return normal_is_in ? 4 : 5; // min(K) : max(K)
1963 case NNC_FACE:
1964 // For nnc faces we return the otherwise unused value -1.
1965 return -1;
1966 default:
1967 OPM_THROW(std::logic_error, "Unhandled face tag. This should never happen!");
1968 }
1969 }
1970
1971 template<int dim>
1972 cpgrid::Entity<dim> createEntity(const CpGrid&, int, bool);
1973
1974} // namespace Dune
1975
1976#include <opm/grid/cpgrid/PersistentContainer.hpp>
1977#include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
1978#include "cpgrid/Intersection.hpp"
1979#include "cpgrid/Geometry.hpp"
1980#include "cpgrid/Indexsets.hpp"
1981
1982#endif // OPM_CPGRID_HEADER
An iterator over the centroids of the geometry of the entities.
Definition CpGrid.hpp:1586
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition CpGrid.hpp:1593
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition CpGrid.hpp:1590
[ provides Dune::Grid ]
Definition CpGrid.hpp:201
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:1332
std::string name() const
Get the grid name.
Definition CpGrid.cpp:697
void switchToGlobalView()
Switch to the global view.
Definition CpGrid.cpp:1432
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition CpGrid.cpp:1047
void addLgrUpdateLeafView(const std::array< int, 3 > &cells_per_dim, const std::array< int, 3 > &startIJK, const std::array< int, 3 > &endIJK, const std::string &lgr_name)
Create a grid out of a coarse one and a refinement(LGR) of a selected block-shaped patch of cells fro...
Definition CpGrid.cpp:2030
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition CpGrid.hpp:216
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:1902
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition CpGrid.cpp:954
const std::vector< std::shared_ptr< Dune::cpgrid::CpGridData > > & currentData() const
Returns either data_ or distributed_data_(if non empty).
Definition CpGrid.cpp:619
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition CpGrid.cpp:1366
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition CpGrid.cpp:949
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition CpGrid.hpp:1917
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition CpGrid.cpp:692
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:549
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:1502
unsigned int overlapSize(int) const
------------— Auxiliary methods to support Adaptivity (end) ------------—
Definition CpGrid.cpp:1025
cpgrid::CpGridDataTraits::InterfaceMap InterfaceMap
The type of the map describing communication interfaces.
Definition CpGrid.hpp:1680
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:1133
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:1351
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition CpGrid.cpp:1144
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition CpGrid.cpp:629
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition CpGrid.cpp:1361
int maxLevel() const
Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview),...
Definition CpGrid.cpp:702
const CpGridTraits::Communication & comm() const
Get the collective communication object.
Definition CpGrid.cpp:1072
double faceArea(int face) const
Get the area of a face.
Definition CpGrid.cpp:1341
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition CpGrid.cpp:1346
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition CpGrid.cpp:959
void postAdapt()
Clean up refinement markers - set every element to the mark 0 which represents 'doing nothing'.
Definition CpGrid.cpp:2023
bool loadBalance(int overlapLayers=1, int partitionMethod=Dune::PartitionMethod::zoltan, double imbalanceTol=1.1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1120
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition CpGrid.cpp:925
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition CpGrid.cpp:612
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition CpGrid.cpp:1154
int numVertices() const
Get The number of vertices.
Definition CpGrid.cpp:1093
Dune::cpgrid::Intersection getParentIntersectionFromLgrBoundaryFace(const Dune::cpgrid::Intersection &intersection) const
Definition CpGrid.cpp:1159
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:682
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition CpGrid.cpp:1103
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition CpGrid.cpp:687
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:2037
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(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::zoltan, double imbalanceTol=1.1)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1207
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition CpGrid.cpp:1113
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:1528
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::zoltan, double imbalanceTol=1.1, bool allowDistributedWells=false)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1295
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:1362
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition CpGrid.cpp:1031
cpgrid::Entity< codim > entity(const cpgrid::Entity< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition CpGrid.cpp:1015
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition CpGrid.cpp:966
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:1427
CpGrid()
Default constructor.
Definition CpGrid.cpp:161
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level and PartitionIteratorType.
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::zoltan)
Distributes this grid over the available nodes in a distributed machine.
Definition CpGrid.hpp:1167
const std::vector< int > & sortedNumAquiferCells() const
Get sorted active cell indices of numerical aquifer.
Definition CpGrid.cpp:1376
bool adapt()
Triggers the grid refinement process.
Definition CpGrid.cpp:1558
void switchToDistributedView()
Switch to the distributed view.
Definition CpGrid.cpp:1438
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition CpGrid.cpp:1371
bool preAdapt()
Set mightVanish flags for elements that will be refined in the next adapt() call Need to be called af...
Definition CpGrid.cpp:1545
int getMark(const cpgrid::Entity< 0 > &element) const
Return refinement mark for entity.
Definition CpGrid.cpp:1540
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:1412
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:1108
int numCells() const
Get the number of cells.
Definition CpGrid.cpp:1083
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition CpGrid.cpp:1422
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:1386
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:190
int numFaces() const
Get the number of faces.
Definition CpGrid.cpp:1088
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition CpGrid.cpp:1336
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition CpGrid.cpp:1194
void globalRefine(int refCount)
Refine the grid refCount times using the default refinement rule.
Definition CpGrid.cpp:971
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition CpGrid.hpp:1450
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:1356
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition CpGrid.hpp:1888
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:1244
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:147
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition CpGridData.hpp:981
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:431
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Indexsets.hpp:53
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
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10
@ simple
Use simple approach based on rectangular partitioning the underlying cartesian grid.
Definition GridEnums.hpp:46
@ 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:188
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:142
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition CpGrid.hpp:144
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition CpGrid.hpp:146
Traits associated with a specific codim.
Definition CpGrid.hpp:118
cpgrid::Entity< cd > Entity
The type of the entity.
Definition CpGrid.hpp:127
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition CpGrid.hpp:121
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition CpGrid.hpp:124
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition CpGrid.hpp:133
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition CpGrid.hpp:130
cpgrid::Entity< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition CpGrid.hpp:136
Traits associated with a specific grid partition type.
Definition CpGrid.hpp:154
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:158
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:156
Definition CpGrid.hpp:98
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition CpGrid.hpp:168
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition CpGrid.hpp:107
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition CpGrid.hpp:165
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition CpGrid.hpp:163
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition CpGrid.hpp:172
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition CpGrid.hpp:109
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition CpGrid.hpp:170
GlobalIdSet LocalIdSet
The type of the local id set.
Definition CpGrid.hpp:174
cpgrid::CpGridDataTraits::Communication Communication
The type of the collective communication.
Definition CpGrid.hpp:177
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition CpGrid.hpp:112
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition CpGrid.hpp:105
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition CpGrid.hpp:103
CpGrid Grid
The type that implements the grid.
Definition CpGrid.hpp:100
Raw corner-point specification of a particular geological model.
Definition preprocess.h:56