dune-alugrid 3.0.0
fromtogridfactory.hh
Go to the documentation of this file.
1#ifndef DUNE_ALUGRID_FROMTOGRIDFACTORY_HH
2#define DUNE_ALUGRID_FROMTOGRIDFACTORY_HH
3
4#include <map>
5#include <vector>
6
7#include <dune/common/version.hh>
8
9#include <dune/grid/common/gridfactory.hh>
10#include <dune/grid/common/exceptions.hh>
11
14
15namespace Dune
16{
17
18 // External Forward Declarations
19 // -----------------------------
20
21 template< class ToGrid >
23
24
25 // FromToGridFactory for ALUGrid
26 // -----------------------------
27
28 template< int dim, int dimworld, ALUGridElementType eltype, ALUGridRefinementType refineType, class Comm >
29 class FromToGridFactory< ALUGrid< dim, dimworld, eltype, refineType, Comm > >
30 {
31 public:
32 // type of grid that is converted to
34 protected:
36
37 std::vector< unsigned int > ordering_ ;
38
39 public:
40 template <class FromGrid, class Vector>
41 Grid* convert( const FromGrid& grid, Vector& cellData, std::vector<unsigned int>& ordering )
42 {
43 int rank = 0;
44#if HAVE_MPI
45 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
46#endif
47
48 // create ALUGrid GridFactory
49 GridFactory< Grid > factory;
50
51 // only insert elements on rank 0
52 if( rank == 0 )
53 {
54 typedef typename FromGrid :: LeafGridView GridView ;
55 typedef typename GridView :: IndexSet IndexSet ;
56 typedef typename IndexSet :: IndexType IndexType ;
57 typedef typename GridView :: template Codim< 0 > :: Iterator ElementIterator ;
58 typedef typename ElementIterator::Entity Entity ;
59 typedef typename GridView :: IntersectionIterator IntersectionIterator ;
60 typedef typename IntersectionIterator :: Intersection Intersection ;
61
62 GridView gridView = grid.leafGridView();
63 const IndexSet &indexSet = gridView.indexSet();
64
65 // map global vertex ids to local ones
66 std::map< IndexType, unsigned int > vtxMap;
67
68 const int numVertices = (1 << dim);
69 std::vector< unsigned int > vertices( numVertices );
70 typedef std::pair< Dune::GeometryType, std::vector< unsigned int > > ElementPair;
71 std::vector< ElementPair > elements;
72 if( ! ordering.empty() )
73 elements.resize( ordering.size() );
74
75 int nextElementIndex = 0;
76 const ElementIterator end = gridView.template end< 0 >();
77 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
78 {
79 const Entity &entity = *it;
80
81 // insert vertices and element
82 const typename Entity::Geometry geo = entity.geometry();
83 alugrid_assert( numVertices == geo.corners() );
84 for( int i = 0; i < numVertices; ++i )
85 {
86 const IndexType vtxId = indexSet.subIndex( entity, i, dim );
87 std::pair< typename std::map< IndexType, unsigned int >::iterator, bool > result
88 = vtxMap.insert( std::make_pair( vtxId, vtxMap.size() ) );
89 if( result.second )
90 factory.insertVertex( geo.corner( i ), vtxId );
91 vertices[ i ] = result.first->second;
92 }
93 if( ordering.empty() )
94 {
95 factory.insertElement( entity.type(), vertices );
96 }
97 else
98 {
99 // store element applying the reordering
100 elements[ ordering[ nextElementIndex++ ] ] = ElementPair( entity.type(), vertices ) ;
101 }
102 }
103
104 if( ! ordering.empty() )
105 {
106 // insert elements using reordered list
107 for( auto it = elements.begin(), end = elements.end(); it != end; ++it )
108 {
109 factory.insertElement( (*it).first, (*it).second );
110 }
111 }
112
113 nextElementIndex = 0;
114 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
115 {
116 const Entity &entity = *it;
117
118 const int elementIndex = ordering.empty() ? nextElementIndex++ : ordering[ nextElementIndex++ ];
119 const IntersectionIterator iend = gridView.iend( entity );
120 for( IntersectionIterator iit = gridView.ibegin( entity ); iit != iend; ++iit )
121 {
122 const Intersection &intersection = *iit;
123 const int faceNumber = intersection.indexInInside();
124 // insert boundary face in case of domain boundary
125 if( intersection.boundary() )
126 factory.insertBoundary( elementIndex, faceNumber );
127
128 // for parallel grids we can check if we are at a process border
129#if DUNE_VERSION_NEWER(DUNE_GRID,3,0)
130 const double isParallel = true;
131#else
132 const double isParallel = Capabilities::isParallel< FromGrid > :: v;
133#endif
134 if( isParallel &&
135 intersection.neighbor() &&
136 intersection.outside().partitionType() != InteriorEntity )
137 {
138 // insert process boundary if the neighboring element has a different rank
139 factory.insertProcessBorder( elementIndex, faceNumber );
140 }
141 }
142 }
143 }
144
145 // create grid pointer (behaving like a shared_ptr)
146 Grid* newgrid = factory.createGrid( true, true, std::string("FromToGrid") );
147
148 if( ! cellData.empty() )
149 {
150 Vector oldCellData( cellData );
151 auto macroView = newgrid->levelGridView( 0 );
152 size_t idx = 0;
153 for( auto it = macroView.template begin<0>(), end = macroView.template end<0>();
154 it != end; ++it, ++idx )
155 {
156 const int insertionIndex = ordering.empty() ?
157 factory.insertionIndex( *it ) : ordering[ factory.insertionIndex( *it ) ]; ;
158 cellData[ idx ] = oldCellData[ insertionIndex ] ;
159 }
160 }
161
162 // store the ordering from the factory, if it was not provided
163 if( ordering.empty() )
164 ordering = factory.ordering();
165
166#if HAVE_MPI
167 MPI_Barrier( MPI_COMM_WORLD );
168#endif
169
170 return newgrid;
171 }
172
173 template <class FromGrid, class Vector>
174 Grid* convert( const FromGrid& fromGrid, Vector& cellData )
175 {
176 return convert( fromGrid, cellData, ordering_ );
177 }
178
179 template <class FromGrid>
180 Grid* convert( const FromGrid& fromGrid )
181 {
182 std::vector<int> dummy(0);
183 return convert( fromGrid, dummy, ordering_ );
184 }
185 protected:
186 template <int codim, class Entity>
187 int subEntities ( const Entity& entity ) const
188 {
189#if DUNE_VERSION_NEWER(DUNE_GRID,2,4)
190 return entity.subEntities( codim );
191#else
192 return entity.template count< codim > ();
193#endif
194 }
195 };
196
197} // namespace Dune
198
199#endif // #ifndef DUNE_ALUGRID_STRUCTUREDGRIDFACTORY_HH
#define alugrid_assert(EX)
Definition alugrid_assert.hh:20
Definition alu3dinclude.hh:50
Definition alu3dinclude.hh:80
unstructured parallel implementation of the DUNE grid interface
Definition alugrid.hh:31
Definition fromtogridfactory.hh:22
Grid * convert(const FromGrid &grid, Vector &cellData, std::vector< unsigned int > &ordering)
Definition fromtogridfactory.hh:41
FromToGridFactory< Grid > This
Definition fromtogridfactory.hh:35
int subEntities(const Entity &entity) const
Definition fromtogridfactory.hh:187
std::vector< unsigned int > ordering_
Definition fromtogridfactory.hh:37
ALUGrid< dim, dimworld, eltype, refineType, Comm > Grid
Definition fromtogridfactory.hh:33
Grid * convert(const FromGrid &fromGrid)
Definition fromtogridfactory.hh:180
Grid * convert(const FromGrid &fromGrid, Vector &cellData)
Definition fromtogridfactory.hh:174