 |
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
|
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
1 /**
2 * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3 * storing and accessing finite element mesh data.
4 *
5 * Copyright 2004 Sandia Corporation. Under the terms of Contract
6 * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7 * retains certain rights in this software.
8 *
9 * This library is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public
11 * License as published by the Free Software Foundation; either
12 * version 2.1 of the License, or (at your option) any later version.
13 *
14 */
15
16 #ifdef WIN32
17 #ifdef _DEBUG
18 // turn off warnings that say they debugging identifier has been truncated
19 // this warning comes up when using some STL containers
20 #pragma warning( disable : 4786 )
21 #endif
22 #endif
23
24 #include "WriteGMV.hpp"
25
26 #include "moab/Interface.hpp"
27 #include "Internals.hpp"
28 #include "moab/Range.hpp"
29 #include "moab/CN.hpp"
30 #include "MBTagConventions.hpp"
31 #include "moab/WriteUtilIface.hpp"
32 #include <fstream>
33 #include <cassert>
34
35 namespace moab
36 {
37
38 const char* WriteGMV::gmvTypeNames[] = { "", "line", "tri", "quad", "", "tet", "pyramid", "prism", "", "hex", "", "" };
39
40 WriterIface* WriteGMV::factory( Interface* iface )
41 {
42 return new WriteGMV( iface );
43 }
44
45 WriteGMV::WriteGMV( Interface* impl ) : mbImpl( impl )
46 {
47 assert( impl != NULL );
48
49 impl->query_interface( mWriteIface );
50
51 // initialize in case tag_get_handle fails below
52 mMaterialSetTag = 0;
53 mDirichletSetTag = 0;
54 mNeumannSetTag = 0;
55 mHasMidNodesTag = 0;
56 mGeomDimensionTag = 0;
57 mGlobalIdTag = 0;
58
59 //! get and cache predefined tag handles
60 // initialize in case tag_get_handle fails below
61 //! get and cache predefined tag handles
62 int negone = -1;
63 impl->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mMaterialSetTag, MB_TAG_SPARSE | MB_TAG_CREAT,
64 &negone );
65
66 impl->tag_get_handle( DIRICHLET_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mDirichletSetTag, MB_TAG_SPARSE | MB_TAG_CREAT,
67 &negone );
68
69 impl->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mNeumannSetTag, MB_TAG_SPARSE | MB_TAG_CREAT,
70 &negone );
71
72 mGlobalIdTag = impl->globalId_tag();
73
74 int dum_val_array[] = { -1, -1, -1, -1 };
75 impl->tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, mHasMidNodesTag, MB_TAG_SPARSE | MB_TAG_CREAT,
76 dum_val_array );
77 }
78
79 WriteGMV::~WriteGMV()
80 {
81 mbImpl->release_interface( mWriteIface );
82 }
83
84 ErrorCode WriteGMV::write_file( const char* file_name,
85 const EntityHandle output_set,
86 const int user_dimension,
87 const bool mesh,
88 const bool poly_mesh )
89 {
90 // general function for writing a mesh
91
92 ErrorCode result = MB_SUCCESS;
93
94 // initialize file
95
96 if( mesh )
97 {
98 result = local_write_mesh( file_name, output_set, user_dimension, true, false );
99 if( MB_SUCCESS != result ) return result;
100 }
101
102 if( poly_mesh )
103 {
104 result = local_write_mesh( file_name, output_set, user_dimension, false, true );
105 if( MB_SUCCESS != result ) return result;
106 }
107
108 return result;
109 }
110
111 ErrorCode WriteGMV::write_file( const char* filename,
112 const bool,
113 const FileOptions& /*opts*/,
114 const EntityHandle* output_sets,
115 const int num_output_sets,
116 const std::vector< std::string >&,
117 const Tag*,
118 int,
119 int dimension )
120 {
121 EntityHandle output_set = 0;
122 if( output_sets && num_output_sets > 0 )
123 {
124 if( num_output_sets > 1 ) return MB_FAILURE;
125 output_set = output_sets[0];
126 }
127
128 if( dimension == 0 )
129 {
130 mbImpl->get_dimension( dimension );
131 }
132
133 return write_file( filename, output_set, dimension, true, true );
134 }
135
136 ErrorCode WriteGMV::local_write_mesh( const char* file_name,
137 const EntityHandle output_set,
138 const int user_dimension,
139 const bool mesh,
140 const bool poly_mesh )
141 {
142 std::ofstream ofile;
143 ErrorCode result;
144
145 if( mesh )
146 {
147 // need to insert ".gmv"
148 std::string tmp_name( file_name );
149 tmp_name += ".gmv";
150 ofile.open( tmp_name.c_str() );
151 }
152 else if( poly_mesh )
153 {
154 // need to insert ".poly.gmv"
155 std::string tmp_name( file_name );
156 tmp_name += ".poly.gmv";
157 ofile.open( tmp_name.c_str() );
158 }
159
160 ofile << "gmvinput ascii" << std::endl;
161
162 // get elements to be output
163 Range dum_range, elements, all_verts;
164 EntityType otype;
165 if( poly_mesh )
166 {
167 result = mbImpl->get_entities_by_type( output_set, MBPOLYGON, elements, true );
168 if( MB_SUCCESS != result ) return result;
169 }
170 else
171 {
172 for( otype = CN::TypeDimensionMap[user_dimension].first; otype <= CN::TypeDimensionMap[user_dimension].second;
173 otype++ )
174 {
175 if( otype == MBPOLYGON || otype == MBPOLYHEDRON ) continue;
176 dum_range.clear();
177 result = mbImpl->get_entities_by_type( output_set, otype, dum_range, true );
178 if( MB_SUCCESS != result ) return result;
179
180 std::copy( dum_range.begin(), dum_range.end(), range_inserter( elements ) );
181 }
182 }
183
184 // gather the vertices in these elements
185 result = mbImpl->get_adjacencies( elements, 0, false, all_verts, Interface::UNION );
186 if( MB_SUCCESS != result ) return result;
187
188 int num_verts = all_verts.size();
189
190 // allocate coordinate arrays and put pointers to them in a list
191 double* xcoord = new double[num_verts];
192 double* ycoord = new double[num_verts];
193 double* zcoord = new double[num_verts];
194 std::vector< double* > coord_arrays;
195 coord_arrays.push_back( xcoord );
196 coord_arrays.push_back( ycoord );
197 coord_arrays.push_back( zcoord );
198
199 // fill them in, writing id tags at the same time
200 result = mWriteIface->get_node_coords( 3, num_verts, all_verts, mGlobalIdTag, 1, coord_arrays );
201 if( MB_SUCCESS != result ) return result;
202
203 int i, j;
204
205 //========================================
206 // WRITE COORDINATE DATA TO FILE HERE
207
208 ofile << "nodev " << num_verts << std::endl;
209 for( i = 0; i < num_verts; i++ )
210 ofile << xcoord[i] << " " << ycoord[i] << " " << zcoord[i] << std::endl;
211
212 //========================================
213
214 delete[] xcoord;
215 delete[] ycoord;
216 delete[] zcoord;
217
218 // iterate over types in selected dimension
219
220 std::vector< int > connect;
221 std::vector< EntityHandle > connecth;
222
223 if( mesh )
224 {
225 Range sub_range;
226
227 ofile << "cells " << elements.size() << std::endl;
228
229 for( otype = CN::TypeDimensionMap[user_dimension].first; otype <= CN::TypeDimensionMap[user_dimension].second;
230 otype++ )
231 {
232
233 if( otype == MBPOLYGON || otype == MBPOLYHEDRON ) continue;
234
235 // get the first element of this type in the range, and one past the last
236 Range::iterator lower =
237 Range::lower_bound( elements.begin(), elements.end(), CREATE_HANDLE( otype, MB_START_ID, i ) );
238 Range::iterator upper =
239 Range::lower_bound( elements.begin(), elements.end(), CREATE_HANDLE( otype + 1, MB_START_ID, i ) );
240
241 if( lower == upper ) continue;
242
243 // copy these elements into a subrange
244 sub_range.clear();
245 std::copy( lower, upper, range_inserter( sub_range ) );
246
247 // make sure the connectivity array is big enough
248 int verts_per = CN::VerticesPerEntity( otype );
249 if( connect.size() < verts_per * sub_range.size() ) connect.resize( verts_per * sub_range.size() );
250
251 // get the connectivity
252 result = mWriteIface->get_element_connect( sub_range.size(), verts_per, mGlobalIdTag, sub_range,
253 mGlobalIdTag, 1, &connect[0] );
254 if( MB_SUCCESS != result ) return result;
255
256 //========================================
257 // WRITE CONNECTIVITY DATA TO FILE HERE
258
259 for( i = 0; i < (int)sub_range.size(); i++ )
260 {
261 ofile << gmvTypeNames[otype] << " " << verts_per << std::endl;
262 for( j = i * verts_per; j < (int)( i + 1 ) * verts_per; j++ )
263 ofile << connect[j] << " ";
264 ofile << std::endl;
265 }
266
267 //========================================
268 }
269 }
270
271 else if( poly_mesh )
272 {
273
274 // write polygons/hedra, if any
275 Range polygons, polyhedra;
276 result = mbImpl->get_entities_by_type( output_set, MBPOLYGON, polygons, true );
277 if( MB_SUCCESS != result ) return result;
278
279 result = mbImpl->get_entities_by_type( output_set, MBPOLYHEDRON, polyhedra, true );
280 if( MB_SUCCESS != result ) return result;
281
282 if( polygons.size() == 0 ) return result;
283
284 // mark polyhedra with global ids
285 result = mWriteIface->assign_ids( polyhedra, mGlobalIdTag, 1 );
286 if( MB_SUCCESS != result ) return result;
287
288 ofile << "faces " << polygons.size() << " " << polyhedra.size() << std::endl;
289
290 for( Range::iterator rit = polygons.begin(); rit != polygons.end(); ++rit )
291 {
292 // get the vertices
293 connecth.clear();
294 result = mbImpl->get_connectivity( &( *rit ), 1, connecth, true );
295 if( MB_SUCCESS != result ) return result;
296
297 if( 0 == connecth.size() ) continue;
298
299 // get the polyhedra, if any
300 if( user_dimension == 3 )
301 {
302 polyhedra.clear();
303 result = mbImpl->get_adjacencies( Range( *rit, *rit ), 3, false, polyhedra );
304 if( MB_SUCCESS != result ) return result;
305
306 // put them in the connect array
307 connecth.push_back( ( polyhedra.size() > 0 ? *polyhedra.begin() : 0 ) );
308 connecth.push_back( ( polyhedra.size() > 1 ? *polyhedra.rbegin() : 0 ) );
309 }
310
311 // replace handles with ids
312 connect.resize( connecth.size() + 2 );
313
314 // pre-set polyhedra ids in case there aren't any
315 connect[connecth.size()] = 0;
316 connect[connecth.size() + 1] = 0;
317 result =
318 mbImpl->tag_get_data( mGlobalIdTag, &connecth[0], connecth.size() - 2 + polyhedra.size(), &connect[0] );
319 if( MB_SUCCESS != result ) return result;
320
321 // write the data
322 ofile << connecth.size() - 2;
323
324 for( i = 0; i < (int)connecth.size(); i++ )
325 ofile << " " << connect[i];
326
327 ofile << std::endl;
328 }
329 }
330
331 ofile << std::endl << "endgmv" << std::endl;
332
333 ofile.close();
334
335 return MB_SUCCESS;
336 }
337
338 } // namespace moab
1.9.1.