Mesh Oriented datABase  (version 5.5.1)
An array-based unstructured mesh library
WriteGMV.cpp
Go to the documentation of this file.
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 
41 {
42  return new WriteGMV( iface );
43 }
44 
45 WriteGMV::WriteGMV( Interface* impl ) : mbImpl( impl )
46 {
47  assert( impl != NULL );
48 
50 
51  // initialize in case tag_get_handle fails below
52  mMaterialSetTag = 0;
53  mDirichletSetTag = 0;
54  mNeumannSetTag = 0;
55  mHasMidNodesTag = 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;
64  &negone );
65 
67  &negone );
68 
70  &negone );
71 
72  mGlobalIdTag = impl->globalId_tag();
73 
74  int dum_val_array[] = { -1, -1, -1, -1 };
76  dum_val_array );
77 }
78 
80 {
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