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." );
121 gridout << std::endl <<
"SIMPLEX" << std::endl;
122 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
124 const Element& element = *it ;
125 if( ! element.type().isSimplex() )
128 std::vector< Index > vertices( dimGrid+1 );
129 for(
size_t i = 0; i < vertices.size(); ++i )
130 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, i, dimGrid ) ];
132 gridout << vertices[ 0 ];
133 for(
size_t i = 1; i < vertices.size(); ++i )
134 gridout <<
" " << vertices[ i ];
135 gridout << std::endl;
137 gridout <<
"#" << std::endl;
140 gridout << std::endl <<
"CUBE" << std::endl;
141 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
143 const Element& element = *it ;
144 if( !element.type().isCube() )
147 std::vector< Index > vertices( 1 << dimGrid );
148 for(
size_t i = 0; i < vertices.size(); ++i )
149 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, i, dimGrid ) ];
151 gridout << vertices[ 0 ];
152 for(
size_t i = 1; i < vertices.size(); ++i )
153 gridout <<
" " << vertices[ i ];
154 gridout << std::endl;
156 gridout <<
"#" << std::endl;
159 gridout << std::endl <<
"BOUNDARYSEGMENTS" << std::endl;
160 for( ElementIterator it = gridView_.template begin< 0 >(); it != end; ++it )
162 const Element& element = *it ;
163 if( !it->hasBoundaryIntersections() )
166 const RefElement &refElement = RefElements::general( element.type() );
168 const IntersectionIterator iend = gridView_.iend( element ) ;
169 for( IntersectionIterator iit = gridView_.ibegin( element ); iit != iend; ++iit )
171 if( !iit->boundary() )
174 const int boundaryId = iit->boundaryId();
175 if( boundaryId <= 0 )
177 std::cerr <<
"Warning: Ignoring nonpositive boundary id: "
178 << boundaryId <<
"." << std::endl;
182 const int faceNumber = iit->indexInInside();
183 const unsigned int faceSize = refElement.size( faceNumber, 1, dimGrid );
184 std::vector< Index > vertices( faceSize );
185 for(
unsigned int i = 0; i < faceSize; ++i )
187 const int j = refElement.subEntity( faceNumber, 1, i, dimGrid );
188 vertices[ i ] = vertexIndex[ indexSet.subIndex( element, j, dimGrid ) ];
190 gridout << boundaryId <<
" " << vertices[ 0 ];
191 for(
unsigned int i = 1; i < faceSize; ++i )
192 gridout <<
" " << vertices[ i ];
193 gridout << std::endl;
196 gridout <<
"#" << std::endl;
199 gridout << std::endl <<
"GRIDPARAMETER" << std::endl;
201 gridout <<
"#" << std::endl;
203 gridout << std::endl <<
"#" << std::endl;
210 std::ofstream gridout( fileName.c_str() );
216 #endif // #ifndef DUNE_DGFWRITER_HH