1 #ifndef DUNE_ALUGRID_FROMTOGRIDFACTORY_HH
2 #define DUNE_ALUGRID_FROMTOGRIDFACTORY_HH
7 #include <dune/common/version.hh>
9 #include <dune/grid/common/gridfactory.hh>
10 #include <dune/grid/common/exceptions.hh>
21 template<
class ToGr
id >
28 template<
int dim,
int dimworld, ALUGr
idElementType eltype, ALUGr
idRefinementType refineType,
class Comm >
40 template <
class FromGr
id,
class Vector>
41 Grid*
convert(
const FromGrid& grid, Vector& cellData, std::vector<unsigned int>& ordering )
45 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
49 GridFactory< Grid > factory;
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 ;
62 GridView gridView = grid.leafGridView();
63 const IndexSet &indexSet = gridView.indexSet();
66 std::map< IndexType, unsigned int > vtxMap;
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() );
75 int nextElementIndex = 0;
76 const ElementIterator end = gridView.template end< 0 >();
77 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
79 const Entity &entity = *it;
82 const typename Entity::Geometry geo = entity.geometry();
84 for(
int i = 0; i < numVertices; ++i )
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() ) );
90 factory.insertVertex( geo.corner( i ), vtxId );
91 vertices[ i ] = result.first->second;
93 if( ordering.empty() )
95 factory.insertElement( entity.type(), vertices );
100 elements[ ordering[ nextElementIndex++ ] ] = ElementPair( entity.type(), vertices ) ;
104 if( ! ordering.empty() )
107 for(
auto it = elements.begin(), end = elements.end(); it != end; ++it )
109 factory.insertElement( (*it).first, (*it).second );
113 nextElementIndex = 0;
114 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
116 const Entity &entity = *it;
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 )
122 const Intersection &intersection = *iit;
123 const int faceNumber = intersection.indexInInside();
125 if( intersection.boundary() )
126 factory.insertBoundary( elementIndex, faceNumber );
129 if( intersection.neighbor() &&
130 intersection.outside().partitionType() != InteriorEntity )
133 factory.insertProcessBorder( elementIndex, faceNumber );
140 Grid* newgrid = factory.createGrid(
true,
true, std::string(
"FromToGrid") );
142 if( ! cellData.empty() )
144 Vector oldCellData( cellData );
147 for(
auto it = macroView.template begin<0>(), end = macroView.template end<0>();
148 it != end; ++it, ++idx )
150 const int insertionIndex = ordering.empty() ?
151 factory.insertionIndex( *it ) : ordering[ factory.insertionIndex( *it ) ]; ;
152 cellData[ idx ] = oldCellData[ insertionIndex ] ;
157 if( ordering.empty() )
158 ordering = factory.ordering();
161 MPI_Barrier( MPI_COMM_WORLD );
167 template <
class FromGr
id,
class Vector>
170 return convert( fromGrid, cellData, ordering_ );
173 template <
class FromGr
id>
176 std::vector<int> dummy(0);
177 return convert( fromGrid, dummy, ordering_ );
180 template <
int codim,
class Entity>
183 return entity.subEntities( codim );
#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
Partition< pitype >::LevelGridView levelGridView(int level) const
Definition: alugrid.hh:165
Definition: fromtogridfactory.hh:22
FromToGridFactory< Grid > This
Definition: fromtogridfactory.hh:35
int subEntities(const Entity &entity) const
Definition: fromtogridfactory.hh:181
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, Vector &cellData)
Definition: fromtogridfactory.hh:168
Grid * convert(const FromGrid &grid, Vector &cellData, std::vector< unsigned int > &ordering)
Definition: fromtogridfactory.hh:41
Grid * convert(const FromGrid &fromGrid)
Definition: fromtogridfactory.hh:174