1 #ifndef DUNE_DGFWRITER_HH
2 #define DUNE_DGFWRITER_HH
13 #include <dune/geometry/referenceelements.hh>
43 typedef typename GridView::template Codim< 0 >::Iterator ElementIterator;
45 typedef typename GridView::template Codim< dimGrid >::EntityPointer VertexPointer;
49 typedef GenericReferenceElement< typename Grid::ctype, dimGrid > RefElement;
50 typedef GenericReferenceElements< typename Grid::ctype, dimGrid > RefElements;
58 : gridView_( gridView )
65 void write ( std::ostream &gridout )
const;
71 void write (
const std::string &fileName )
const;
83 gridout.setf( std::ios_base::scientific, std::ios_base::floatfield );
84 gridout.precision( 16 );
86 const IndexSet &indexSet = gridView_.indexSet();
89 gridout <<
"DGF" << std::endl;
91 const Index vxSize = indexSet.size( dimGrid );
92 std::vector< Index > vertexIndex( vxSize, vxSize );
94 gridout <<
"%" <<
" Elements = " << indexSet.size( 0 ) <<
" | Vertices = " << vxSize << std::endl;
97 gridout << std::endl <<
"VERTEX" << std::endl;
98 Index vertexCount = 0;
99 typedef typename ElementIterator :: Entity
Element ;
100 const ElementIterator end = gridView_.template end< 0 >();
101 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
103 const Element& element = *it ;
104 const int numCorners = element.template count< dimGrid > ();
105 for(
int i=0; i<numCorners; ++i )
107 const Index vxIndex = indexSet.subIndex( element, i, dimGrid );
108 assert( vxIndex < vxSize );
109 if( vertexIndex[ vxIndex ] == vxSize )
111 vertexIndex[ vxIndex ] = vertexCount++;
112 gridout << element.geometry().corner( i ) << std::endl;
116 gridout <<
"#" << std::endl;
117 if( vertexCount != vxSize )
118 DUNE_THROW(
GridError,
"Index set reports wrong number of vertices." );
123 gridout << std::endl <<
"SIMPLEX" << std::endl;
124 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
126 const Element& element = *it ;
127 if( !element.type().isSimplex() )
130 std::vector< Index > vertices( dimGrid+1 );
131 for(
size_t i = 0; i < vertices.size(); ++i )
132 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, i, dimGrid ) ];
134 gridout << vertices[ 0 ];
135 for(
size_t i = 1; i < vertices.size(); ++i )
136 gridout <<
" " << vertices[ i ];
137 gridout << std::endl;
139 gridout <<
"#" << std::endl;
143 gridout << std::endl <<
"CUBE" << std::endl;
144 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
146 const Element& element = *it ;
147 if( !element.type().isCube() )
150 std::vector< Index > vertices( 1 << dimGrid );
151 for(
size_t i = 0; i < vertices.size(); ++i )
152 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, i, dimGrid ) ];
154 gridout << vertices[ 0 ];
155 for(
size_t i = 1; i < vertices.size(); ++i )
156 gridout <<
" " << vertices[ i ];
157 gridout << std::endl;
159 gridout <<
"#" << std::endl;
162 gridout << std::endl <<
"BOUNDARYSEGMENTS" << std::endl;
163 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
165 const Element& element = *it ;
166 if( !it->hasBoundaryIntersections() )
169 const RefElement &refElement = RefElements::general( element.type() );
171 const IntersectionIterator iend = gridView_.iend( element ) ;
172 for( IntersectionIterator iit = gridView_.ibegin( element ); iit != iend; ++iit )
174 if( !iit->boundary() )
177 const int boundaryId = iit->boundaryId();
178 if( boundaryId <= 0 )
180 std::cerr <<
"Warning: Ignoring nonpositive boundary id: "
181 << boundaryId <<
"." << std::endl;
185 const int faceNumber = iit->indexInInside();
186 const unsigned int faceSize = refElement.size( faceNumber, 1, dimGrid );
187 std::vector< Index > vertices( faceSize );
188 for(
unsigned int i = 0; i < faceSize; ++i )
190 const int j = refElement.subEntity( faceNumber, 1, i, dimGrid );
191 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, j, dimGrid ) ];
193 gridout << boundaryId <<
" " << vertices[ 0 ];
194 for(
unsigned int i = 1; i < faceSize; ++i )
195 gridout <<
" " << vertices[ i ];
196 gridout << std::endl;
199 gridout <<
"#" << std::endl;
202 gridout << std::endl <<
"GRIDPARAMETER" << std::endl;
204 gridout <<
"#" << std::endl;
206 gridout << std::endl <<
"#" << std::endl;
213 std::ofstream gridout( fileName.c_str() );
219 #endif // #ifndef DUNE_DGFWRITER_HH