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 >
38 typedef decltype(std::declval< Dune::GridFactoryInterface< Grid >* >()->createGrid()) GridPtrType;
43 template <
class FromGr
id,
class Vector>
44 GridPtrType
convert(
const FromGrid& grid, Vector& cellData, std::vector<unsigned int>& ordering )
48 MPI_Comm_rank( MPI_COMM_WORLD, &rank );
52 GridFactory< Grid > factory;
57 typedef typename FromGrid :: LeafGridView GridView ;
58 typedef typename GridView :: IndexSet IndexSet ;
59 typedef typename IndexSet :: IndexType IndexType ;
60 typedef typename GridView :: template Codim< 0 > :: Iterator ElementIterator ;
61 typedef typename ElementIterator::Entity Entity ;
62 typedef typename GridView :: IntersectionIterator IntersectionIterator ;
63 typedef typename IntersectionIterator :: Intersection Intersection ;
65 GridView gridView = grid.leafGridView();
66 const IndexSet &indexSet = gridView.indexSet();
69 std::map< IndexType, unsigned int > vtxMap;
71 const int numVertices = (1 << dim);
72 std::vector< unsigned int > vertices( numVertices );
73 typedef std::pair< Dune::GeometryType, std::vector< unsigned int > > ElementPair;
74 std::vector< ElementPair > elements;
75 if( ! ordering.empty() )
76 elements.resize( ordering.size() );
78 int nextElementIndex = 0;
79 const ElementIterator end = gridView.template end< 0 >();
80 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
82 const Entity &entity = *it;
85 const typename Entity::Geometry geo = entity.geometry();
87 for(
int i = 0; i < numVertices; ++i )
89 const IndexType vtxId = indexSet.subIndex( entity, i, dim );
90 std::pair< typename std::map< IndexType, unsigned int >::iterator,
bool > result
91 = vtxMap.insert( std::make_pair( vtxId, vtxMap.size() ) );
93 factory.insertVertex( geo.corner( i ), vtxId );
94 vertices[ i ] = result.first->second;
96 if( ordering.empty() )
98 factory.insertElement( entity.type(), vertices );
103 elements[ ordering[ nextElementIndex++ ] ] = ElementPair( entity.type(), vertices ) ;
107 if( ! ordering.empty() )
110 for(
auto it = elements.begin(), end = elements.end(); it != end; ++it )
112 factory.insertElement( (*it).first, (*it).second );
116 nextElementIndex = 0;
117 for( ElementIterator it = gridView.template begin< 0 >(); it != end; ++it )
119 const Entity &entity = *it;
121 const int elementIndex = ordering.empty() ? nextElementIndex++ : ordering[ nextElementIndex++ ];
122 const IntersectionIterator iend = gridView.iend( entity );
123 for( IntersectionIterator iit = gridView.ibegin( entity ); iit != iend; ++iit )
125 const Intersection &intersection = *iit;
126 const int faceNumber = intersection.indexInInside();
128 if( intersection.boundary() )
129 factory.insertBoundary( elementIndex, faceNumber );
132 if( intersection.neighbor() &&
133 intersection.outside().partitionType() != InteriorEntity )
136 factory.insertProcessBorder( elementIndex, faceNumber );
143 GridPtrType newgrid = factory.createGrid(
true,
true, std::string(
"FromToGrid") );
145 if( ! cellData.empty() )
147 Vector oldCellData( cellData );
148 auto macroView = newgrid->levelGridView( 0 );
150 for(
auto it = macroView.template begin<0>(), end = macroView.template end<0>();
151 it != end; ++it, ++idx )
153 const int insertionIndex = ordering.empty() ?
154 factory.insertionIndex( *it ) : ordering[ factory.insertionIndex( *it ) ]; ;
155 cellData[ idx ] = oldCellData[ insertionIndex ] ;
160 if( ordering.empty() )
161 ordering = factory.ordering();
164 MPI_Barrier( MPI_COMM_WORLD );
170 template <
class FromGr
id,
class Vector>
171 GridPtrType
convert(
const FromGrid& fromGrid, Vector& cellData )
173 return convert( fromGrid, cellData, ordering_ );
176 template <
class FromGr
id>
177 GridPtrType
convert(
const FromGrid& fromGrid )
179 std::vector<int> dummy(0);
180 return convert( fromGrid, dummy, ordering_ );
183 template <
int codim,
class Entity>
186 return entity.subEntities( codim );
#define alugrid_assert(EX)
Definition: alugrid_assert.hh:20
Definition: alu3dinclude.hh:33
Definition: alu3dinclude.hh:63
unstructured parallel implementation of the DUNE grid interface
Definition: alugrid.hh:31
Definition: fromtogridfactory.hh:22
FromToGridFactory< Grid > This
Definition: fromtogridfactory.hh:35
int subEntities(const Entity &entity) const
Definition: fromtogridfactory.hh:184
std::vector< unsigned int > ordering_
Definition: fromtogridfactory.hh:40
ALUGrid< dim, dimworld, eltype, refineType, Comm > Grid
Definition: fromtogridfactory.hh:33
GridPtrType convert(const FromGrid &grid, Vector &cellData, std::vector< unsigned int > &ordering)
Definition: fromtogridfactory.hh:44
GridPtrType convert(const FromGrid &fromGrid)
Definition: fromtogridfactory.hh:177
GridPtrType convert(const FromGrid &fromGrid, Vector &cellData)
Definition: fromtogridfactory.hh:171