MOAB: Mesh Oriented datABase  (version 5.5.0)
h5file_test.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 #include <cstdio>
17 #include <cmath>
18 
19 #include "moab/Core.hpp"
20 #include "moab/Range.hpp"
21 #include "MBTagConventions.hpp"
22 #include "moab/CN.hpp"
23 
24 #define filename "h5test.h5m"
25 
26 #include <iostream>
27 #include <sstream>
28 
29 using namespace moab;
30 
31 // Dense tag name
32 // Pick an ugly name to check special-char escaping in HDF5 file
33 #define tagname "foo/\\/\\"
34 #define bitname "bar\n"
35 #define intname "int tag"
36 #define dblname " dbl "
37 #define handlename "hanlde"
38 
39 const int FACE_SET_ID = 1101;
40 const int VERTEX_SET_ID = 1102;
41 const int REGION_SET_ID = 1103;
42 const int EMPTY_SET_ID = 1100;
43 const int SET_SET_ID = 1105;
44 
45 Interface* iface;
46 
47 void create();
48 
49 bool compare();
50 
51 void moab_error( const char* function );
52 
53 int main( int argc, char* argv[] )
54 {
55  ErrorCode rval;
56  std::string msg;
57 
58  std::string read_opt;
59  std::string write_opt( "DEBUG_BINIO" );
60 
61  for( int i = 1; i < argc; ++i )
62  {
63  long val;
64  char* endptr;
65  if( argv[i][0] == '-' && ( argv[i][1] == 'r' || argv[i][1] == 'w' ) && i + 1 < argc &&
66  ( val = strtol( argv[i + 1], &endptr, 0 ) ) >= 0 && !*endptr )
67  {
68  std::string& s = argv[i][1] == 'r' ? read_opt : write_opt;
69  std::ostringstream str;
70  str << "DEBUG_IO=" << val;
71  if( s.empty() )
72  s = str.str();
73  else
74  {
75  s += ';';
76  s += str.str();
77  }
78  ++i;
79  }
80  else
81  {
82  std::cerr << "Usage: " << argv[0] << " [-r <n>] [-w <n>]" << std::endl;
83  return 1;
84  }
85  }
86 
87  iface = new Core();
88 
89  // create a dodecahedron and inscribed hex
90  fprintf( stderr, "creating... " );
91  create();
92 
93  // write out the dodecahedron
94  fprintf( stderr, "writing... " );
95  rval = iface->write_file( filename, 0, write_opt.c_str() );
96  if( MB_SUCCESS != rval )
97  {
98  fprintf( stderr, "Failed to write \"%s\"\n", filename );
99  if( MB_SUCCESS == iface->get_last_error( msg ) ) fprintf( stderr, "%s\n", msg.c_str() );
100  delete iface;
101  return 1;
102  }
103 
104  // Read back in as a copy of the original
105  fprintf( stderr, "reading... " );
106  rval = iface->load_file( filename, 0, read_opt.c_str() );
107  if( MB_SUCCESS != rval )
108  {
109  fprintf( stderr, "Failed to read \"%s\"\n", filename );
110  if( MB_SUCCESS == iface->get_last_error( msg ) ) fprintf( stderr, "%s\n", msg.c_str() );
111  delete iface;
112  return 1;
113  }
114 
115  // Compare the two.
116  fprintf( stderr, "comparing... " );
117  if( !compare() )
118  {
119  fprintf( stderr, "Comparison failed.\n" );
120  delete iface;
121  return 1;
122  }
123  fprintf( stderr, "success!\n" );
124 
125  // Write both the original and copy to a file
126  fprintf( stderr, "writing... " );
127  rval = iface->write_file( filename, 0, write_opt.c_str() );
128  if( MB_SUCCESS != rval )
129  {
130  fprintf( stderr, "Failed to write \"%s\"\n", filename );
131  if( MB_SUCCESS == iface->get_last_error( msg ) ) fprintf( stderr, "%s\n", msg.c_str() );
132  delete iface;
133  return 1;
134  }
135 
136  // Delete the mesh
137  fprintf( stderr, "clearing db... " );
138  rval = iface->delete_mesh();
139  if( MB_SUCCESS != rval ) moab_error( "delete_mesh" );
140 
141  // Read the two dodecahedrons from the file
142  fprintf( stderr, "reading... " );
143  rval = iface->load_file( filename, 0, read_opt.c_str() );
144  if( MB_SUCCESS != rval )
145  {
146  fprintf( stderr, "Failed to read \"%s\"\n", filename );
147  if( MB_SUCCESS == iface->get_last_error( msg ) ) fprintf( stderr, "%s\n", msg.c_str() );
148  delete iface;
149  return 1;
150  }
151 
152  // Compare them
153  fprintf( stderr, "comparing... " );
154  if( !compare() )
155  {
156  fprintf( stderr, "Comparison failed.\n" );
157  delete iface;
158  return 1;
159  }
160  fprintf( stderr, "success!\n" );
161 
162  // Delete the mesh
163  fprintf( stderr, "cleaning up... " );
164  rval = iface->delete_mesh();
165  if( MB_SUCCESS != rval ) moab_error( "delete_mesh" );
166 
167  // Clean up the file.
168  remove( filename );
169  fprintf( stderr, "done.\n" );
170  delete iface;
171  return 0;
172 }
173 
174 EntityHandle vtx( double x, double y, double z )
175 {
176  const double p[3] = { x, y, z };
177  EntityHandle result;
178  if( MB_SUCCESS != iface->create_vertex( p, result ) ) moab_error( "create_vertex" );
179  return result;
180 }
181 
182 EntityHandle pent( EntityHandle* vtx_list, int i1, int i2, int i3, int i4, int i5 )
183 {
184  const EntityHandle conn[5] = { vtx_list[i1], vtx_list[i2], vtx_list[i3], vtx_list[i4], vtx_list[i5] };
185  EntityHandle result;
186  if( MB_SUCCESS != iface->create_element( MBPOLYGON, conn, 5, result ) ) moab_error( "create_element" );
187  return result;
188 }
189 
190 // Create an entity set containing the passed entities.
191 // options - set options to pass to create_meshset
192 // entities - array of entities to put in set
193 // num_entities - length of 'entities'
194 // reverse - if true, add entities in reverse order
195 // id - value for global id on set
196 EntityHandle make_set( unsigned int options, EntityHandle* entities, size_t num_entities, bool reverse, int id )
197 {
198  EntityHandle handle;
199  if( MB_SUCCESS != iface->create_meshset( options, handle ) ) moab_error( "create_meshset" );
200 
201  if( reverse )
202  {
203  for( int i = (int)num_entities - 1; i >= 0; --i )
204  if( MB_SUCCESS != iface->add_entities( handle, entities + i, 1 ) ) moab_error( "add_entities" );
205  }
206  else
207  {
208  if( MB_SUCCESS != iface->add_entities( handle, entities, num_entities ) ) moab_error( "add_entities" );
209  }
210 
211  Tag id_tag = iface->globalId_tag();
212  if( MB_SUCCESS != iface->tag_set_data( id_tag, &handle, 1, &id ) ) moab_error( "tag_set_data" );
213 
214  return handle;
215 }
216 
217 void create()
218 {
219  // Create dodecahedron
220 
221  // radius
222  const double r = 50.;
223  // center
224  const double x = 0., y = 0., z = 0.;
225  // length of edge of inscribed cube
226  const double cube = r * ( 2.0 / sqrt( 3.0 ) );
227  // length of dodecahedron edge
228  const double edge = cube * ( 2.0 / ( 1.0 + sqrt( 5.0 ) ) );
229  // distance of projection of a dodecahedron vertex to
230  // closest edge of inscribed cube
231  const double p = ( cube - edge ) / 2.0;
232  // distance of projection of a dodecahedron vertex to
233  // closest face of inscribed cube
234  const double d = sqrt( edge * edge - cube * cube / 4 - p );
235  // coordinate values
236  const double c = cube / 2.0;
237  const double a = c + d;
238  const double b = edge / 2.0;
239 
240  // list of vertex handles
241  EntityHandle vertices[20];
242  // list of pentagon handles
243  EntityHandle faces[12];
244  // Dodecahedron handle
245  EntityHandle dodec;
246  // Inscribed Hex handle
247  EntityHandle hex;
248 
249  // Create vertices if inscribed cube
250  vertices[0] = vtx( x - c, y + c, z + c );
251  vertices[1] = vtx( x + c, y + c, z + c );
252  vertices[2] = vtx( x + c, y - c, z + c );
253  vertices[3] = vtx( x - c, y - c, z + c );
254  vertices[4] = vtx( x - c, y + c, z - c );
255  vertices[5] = vtx( x + c, y + c, z - c );
256  vertices[6] = vtx( x + c, y - c, z - c );
257  vertices[7] = vtx( x - c, y - c, z - c );
258 
259  // Create inscribed hex
260  if( MB_SUCCESS != iface->create_element( MBHEX, vertices, 8, hex ) ) moab_error( "create_element" );
261 
262  // Create vertices, 2 "above" each face of inscribed cube
263  // +z face
264  vertices[8] = vtx( x - b, y, z + a );
265  vertices[9] = vtx( x + b, y, z + a );
266  // +x face
267  vertices[10] = vtx( x + a, y + b, z );
268  vertices[11] = vtx( x + a, y - b, z );
269  // -z face
270  vertices[12] = vtx( x - b, y, z - a );
271  vertices[13] = vtx( x + b, y, z - a );
272  // -x face
273  vertices[14] = vtx( x - a, y + b, z );
274  vertices[15] = vtx( x - a, y - b, z );
275  // +y face
276  vertices[16] = vtx( x, y + a, z + b );
277  vertices[17] = vtx( x, y + a, z - b );
278  // -y face
279  vertices[18] = vtx( x, y - a, z + b );
280  vertices[19] = vtx( x, y - a, z - b );
281 
282  // Create petagons
283  faces[0] = pent( vertices, 0, 8, 9, 1, 16 );
284  faces[1] = pent( vertices, 3, 18, 2, 9, 8 );
285  faces[2] = pent( vertices, 2, 11, 10, 1, 9 );
286  faces[3] = pent( vertices, 2, 18, 19, 6, 11 );
287  faces[4] = pent( vertices, 1, 10, 5, 17, 16 );
288  faces[5] = pent( vertices, 5, 10, 11, 6, 13 );
289  faces[6] = pent( vertices, 4, 17, 5, 13, 12 );
290  faces[7] = pent( vertices, 7, 12, 13, 6, 19 );
291  faces[8] = pent( vertices, 4, 12, 7, 15, 14 );
292  faces[9] = pent( vertices, 0, 16, 17, 4, 14 );
293  faces[10] = pent( vertices, 3, 15, 7, 19, 18 );
294  faces[11] = pent( vertices, 0, 14, 15, 3, 8 );
295 
296  // Create dodecahedron
297  if( MB_SUCCESS != iface->create_element( MBPOLYHEDRON, faces, 12, dodec ) ) moab_error( "create_element" );
298 
299  // Create a dense tag
300  int zero = 0;
301  Tag tag;
302  if( MB_SUCCESS != iface->tag_get_handle( tagname, 1, MB_TYPE_INTEGER, tag, MB_TAG_DENSE | MB_TAG_EXCL, &zero ) )
303  moab_error( "tag_get_handle" );
304 
305  // Put dense tag on all vertices.
306  for( int i = 0; i < 20; ++i )
307  if( MB_SUCCESS != iface->tag_set_data( tag, vertices + i, 1, &i ) ) moab_error( "tag_set_data" );
308 
309  // Create bit tag
310  Tag tag2;
311  if( MB_SUCCESS != iface->tag_get_handle( bitname, 2, MB_TYPE_BIT, tag2, MB_TAG_EXCL ) )
312  moab_error( "tag_get_handle" );
313 
314  // Set tag to 0 on Hex
315  char two = '\002';
316  if( MB_SUCCESS != iface->tag_set_data( tag2, &hex, 1, &two ) ) moab_error( "tag_set_data" );
317 
318  // set tag to 1 on dodecahedron
319  char one = '\001';
320  if( MB_SUCCESS != iface->tag_set_data( tag2, &dodec, 1, &one ) ) moab_error( "tag_set_data" );
321 
322  // Create an integer array tag and set some values on the dodecahedron
323  Tag itag;
324  if( MB_SUCCESS != iface->tag_get_handle( intname, 2, MB_TYPE_INTEGER, itag, MB_TAG_SPARSE | MB_TAG_EXCL ) )
325  moab_error( "tag_get_handle(MB_TYPE_INT)" );
326  int idata[] = { static_cast< int >( 0xDEADBEEF ), static_cast< int >( 0xDEFACED ) };
327  if( MB_SUCCESS != iface->tag_set_data( itag, &dodec, 1, idata ) ) moab_error( "tag_set_data(itag)" );
328 
329  // Create a double tag with a non-zero default value, and set on dodecahedron
330  Tag dtag;
331  double ddef = 3.14159;
332  if( MB_SUCCESS != iface->tag_get_handle( dblname, 1, MB_TYPE_DOUBLE, dtag, MB_TAG_DENSE | MB_TAG_EXCL, &ddef ) )
333  moab_error( "tag_get_handle(dtag)" );
334  double dval = 0.333;
335  if( MB_SUCCESS != iface->tag_set_data( dtag, &dodec, 1, &dval ) ) moab_error( "tag_set_data(dtag)" );
336 
337  // Create a tag containing entity handles, with default values
338  Tag htag;
339  EntityHandle hdef[] = { hex, dodec, 0 };
340  if( MB_SUCCESS != iface->tag_get_handle( handlename, 3, MB_TYPE_HANDLE, htag, MB_TAG_SPARSE | MB_TAG_EXCL, hdef ) )
341  moab_error( "tag_get_handle(htag)" );
342  // Set global (mesh) value for tag
343  EntityHandle hgbl[] = { 0, hex, dodec };
344  const EntityHandle root = 0;
345  if( MB_SUCCESS != iface->tag_set_data( htag, &root, 1, hgbl ) ) moab_error( "tag_set_data(hgbl)" );
346  // Store first three entiries of dodec connectivity on dodec
347  EntityHandle hval[] = { faces[0], faces[1], faces[2] };
348  if( MB_SUCCESS != iface->tag_set_data( htag, &dodec, 1, hval ) ) moab_error( "tag_set_data(hgbl)" );
349 
350  // create some sets
351  EntityHandle face_set, vertex_set, region_set, empty_set;
352  EntityHandle regions[] = { dodec, hex };
353  const unsigned empty_flags = MESHSET_ORDERED | MESHSET_TRACK_OWNER;
354  face_set = make_set( MESHSET_SET, faces, 12, false, FACE_SET_ID );
355  vertex_set = make_set( MESHSET_ORDERED, vertices, 20, true, VERTEX_SET_ID );
356  region_set = make_set( MESHSET_SET, regions, 2, false, REGION_SET_ID );
357  empty_set = make_set( empty_flags, 0, 0, true, EMPTY_SET_ID );
358  EntityHandle sets[] = { face_set, vertex_set, region_set, empty_set };
359  make_set( MESHSET_ORDERED, sets, 4, false, SET_SET_ID );
360 
361  // create some set parent-child links
362  if( MB_SUCCESS != iface->add_parent_child( face_set, vertex_set ) ) moab_error( "add_parent_child" );
363  if( MB_SUCCESS != iface->add_child_meshset( region_set, face_set ) ) moab_error( "add_child_meshset" );
364  if( MB_SUCCESS != iface->add_parent_meshset( vertex_set, region_set ) ) moab_error( "add_parent_meshet" );
365 }
366 
367 bool compare_conn( std::vector< EntityHandle >& conn1, std::vector< EntityHandle >& conn2 )
368 {
369  unsigned i;
370 
371  if( conn1.size() != conn2.size() || conn1.size() == 0 )
372  {
373  fprintf( stderr, "Error comparing connectivity: sizes %lu and %lu\n", (unsigned long)conn1.size(),
374  (unsigned long)conn2.size() );
375  return false;
376  }
377 
378  std::vector< double > coords[2];
379  coords[0].resize( 3 * conn1.size() );
380  coords[1].resize( 3 * conn2.size() );
381  if( MB_SUCCESS != iface->get_coords( &conn1[0], conn1.size(), &coords[0][0] ) ||
382  MB_SUCCESS != iface->get_coords( &conn2[0], conn2.size(), &coords[1][0] ) ||
383  coords[0].size() != coords[1].size() )
384  moab_error( "get_coords" );
385 
386  std::vector< double >::iterator citer1 = coords[0].begin(), citer2 = coords[1].begin();
387  for( i = 0; i < conn1.size(); i++ )
388  {
389  double x1 = *( citer1++ ), y1 = *( citer1++ ), z1 = *( citer1++ );
390  double x2 = *( citer2++ ), y2 = *( citer2++ ), z2 = *( citer2++ );
391  if( x1 != x2 || y1 != y2 || z1 != z2 )
392  {
393  fprintf( stderr, "Vertex coords don't match: ( %f, %f, %f ) and ( %f, %f, %f )\n", x1, y1, z1, x2, y2, z2 );
394  return false;
395  }
396  }
397 
398  std::vector< int > tags[2];
399  tags[0].resize( conn1.size() );
400  tags[1].resize( conn2.size() );
401  Tag tag;
402  if( MB_SUCCESS != iface->tag_get_handle( tagname, 1, MB_TYPE_INTEGER, tag ) ) moab_error( "tag_get_handle" );
403  if( MB_SUCCESS != iface->tag_get_data( tag, &conn1[0], conn1.size(), &tags[0][0] ) ||
404  MB_SUCCESS != iface->tag_get_data( tag, &conn2[0], conn2.size(), &tags[1][0] ) )
405  moab_error( "tag_get_data" );
406 
407  std::vector< int >::iterator titer1 = tags[0].begin(), titer2 = tags[1].begin();
408  for( i = 0; i < conn1.size(); i++ )
409  {
410  int t1 = *( titer1++ );
411  int t2 = *( titer2++ );
412  if( t1 != t2 )
413  {
414  fprintf( stderr, "Vertex tags don't match: %d != %d\n", t1, t2 );
415  return false;
416  }
417  }
418 
419  return true;
420 }
421 
422 /* Compare the two sets with the specified global id.
423  If tag_name is not null, compare all entities in the
424  sets using the value of the tag to match entities.
425  Tag must be an integer type.
426  */
427 bool compare_sets( int id, const char* tag_name = 0 )
428 {
429  bool ok;
430  ErrorCode rval;
431 
432  // get sets
433 
434  Tag id_tag = iface->globalId_tag();
435 
436  Range range;
437  const void* tag_data[] = { &id };
438  rval = iface->get_entities_by_type_and_tag( 0, MBENTITYSET, &id_tag, tag_data, 1, range );
439  if( MB_ENTITY_NOT_FOUND == rval || range.size() != 2 )
440  {
441  fprintf( stderr, "Could not find set with id %d pair in file\n", id );
442  return false;
443  }
444  else if( MB_SUCCESS != rval )
445  moab_error( "get_entities_by_type_and_tag" );
446 
447  EntityHandle set1 = *range.begin();
448  EntityHandle set2 = *++range.begin();
449 
450  // Compare set descriptions
451  unsigned opt1 = 0, opt2 = 0;
452  rval = iface->get_meshset_options( set1, opt1 );
453  if( MB_SUCCESS != rval ) moab_error( "get_meshset_options for set1 failed" );
454  rval = iface->get_meshset_options( set2, opt2 );
455  if( MB_SUCCESS != rval ) moab_error( "get_meshset_options for set2 failed" );
456 
457  if( opt1 != opt2 )
458  {
459  fprintf( stderr,
460  "Sets with id %d do not have matching options.\n"
461  "Set 1: track_owner=%s set=%s ordered=%s\n"
462  "Set 2: track_owner=%s set=%s ordered=%s\n",
463  id, opt1 & MESHSET_TRACK_OWNER ? "yes" : "no", opt1 & MESHSET_SET ? "yes" : "no",
464  opt1 & MESHSET_ORDERED ? "yes" : "no", opt2 & MESHSET_TRACK_OWNER ? "yes" : "no",
465  opt2 & MESHSET_SET ? "yes" : "no", opt2 & MESHSET_ORDERED ? "yes" : "no" );
466  return false;
467  }
468 
469  // Compare set contents
470  // First check if same number of entities.
471  // Then select from three possible methods to compare set contents
472  // o If caller gave us a tag to use, compare entities by tag value
473  // o If set is ordered, compare entity types in order
474  // o Otherwise compare counts of entity types in sets
475  std::vector< EntityHandle > list1, list2;
476  rval = iface->get_entities_by_handle( set1, list1 );
477  if( MB_SUCCESS != rval ) moab_error( "get_entities_by_handle for set1 failed" );
478  rval = iface->get_entities_by_handle( set2, list2 );
479  if( MB_SUCCESS != rval ) moab_error( "get_entities_by_handle for set2 failed" );
480 
481  if( list1.size() != list2.size() )
482  {
483  fprintf( stderr,
484  "Sets with id %d do not have the same number of entities.\n"
485  "Set 1 : %u Set 2 : %u\n",
486  id, (unsigned)list1.size(), (unsigned)list2.size() );
487  return false;
488  }
489 
490  if( tag_name ) // compare contents using tag value
491  {
492  Tag tag;
493  if( MB_SUCCESS != iface->tag_get_handle( tag_name, 0, MB_TYPE_OPAQUE, tag, MB_TAG_ANY ) )
494  {
495  fprintf( stderr, "Could not find tag \"%s\" in file.\n", tag_name );
496  return false;
497  }
498 
499  // Make sure tag is integer type
500  DataType type;
501  if( MB_SUCCESS != iface->tag_get_data_type( tag, type ) ) moab_error( "tag_get_data_type" );
502  if( MB_TYPE_INTEGER != type && MB_TYPE_OPAQUE != type ) moab_error( "compare_sets" );
503 
504  std::vector< int > data1( list1.size() ), data2( list2.size() );
505  if( MB_SUCCESS != iface->tag_get_data( tag, &list1[0], list1.size(), &data1[0] ) ||
506  MB_SUCCESS != iface->tag_get_data( tag, &list2[0], list2.size(), &data2[0] ) )
507  moab_error( "tag_get_data" );
508 
509  if( !( opt1 & MESHSET_ORDERED ) )
510  {
511  std::sort( data1.begin(), data1.end() );
512  std::sort( data2.begin(), data2.end() );
513  }
514 
515  for( unsigned i = 0; i < data1.size(); ++i )
516  if( data1[i] != data2[i] )
517  {
518  fprintf( stderr, "Entities in sets with id %d do not have matching tag values.\n", id );
519  return false;
520  }
521  }
522  else if( opt1 & MESHSET_ORDERED ) // compare types of each entity in set
523  {
524  ok = true;
525  for( unsigned i = 0; i < list1.size(); ++i )
526  if( iface->type_from_handle( list1[i] ) != iface->type_from_handle( list2[i] ) )
527  {
528  fprintf( stderr,
529  "Entities at position %u in ordered sets with id %d\n"
530  "have different types.\n",
531  i, id );
532  ok = false;
533  }
534  if( !ok ) return false;
535  }
536  else // compare count of entity types in each set
537  {
538  unsigned counts1[MBMAXTYPE], counts2[MBMAXTYPE];
539  memset( counts1, 0, MBMAXTYPE * sizeof( unsigned ) );
540  memset( counts2, 0, MBMAXTYPE * sizeof( unsigned ) );
541  for( unsigned i = 0; i < list1.size(); ++i )
542  {
543  counts1[iface->type_from_handle( list1[i] )]++;
544  counts2[iface->type_from_handle( list2[i] )]++;
545  }
546 
547  ok = true;
548  for( int j = 0; j < MBMAXTYPE; ++j )
549  if( counts1[j] != counts2[j] )
550  {
551  fprintf( stderr, "Sets with id %d have differing numbers of %s: %u and %u\n", id,
552  CN::EntityTypeName( (EntityType)j ), counts1[j], counts2[j] );
553  ok = false;
554  }
555  if( !ok ) return false;
556  }
557 
558  // Compare set parent/child links using global id
559 
560  ok = true;
561  const char* words[] = { "children", "parents" };
562  std::vector< EntityHandle > adj1, adj2;
563  for( unsigned two = 0; two < 2; ++two )
564  {
565  adj1.clear();
566  adj2.clear();
567  if( two )
568  {
569  if( MB_SUCCESS != iface->get_parent_meshsets( set1, adj1 ) ||
570  MB_SUCCESS != iface->get_parent_meshsets( set2, adj2 ) )
571  moab_error( "get_parent_meshsets" );
572  }
573  else
574  {
575  if( MB_SUCCESS != iface->get_child_meshsets( set1, adj1 ) ||
576  MB_SUCCESS != iface->get_child_meshsets( set2, adj2 ) )
577  moab_error( "get_child_meshsets" );
578  }
579 
580  if( adj1.size() != adj2.size() )
581  {
582  fprintf( stderr, "Sets with id %d have different number of %s: %u and %u\n", id, words[two],
583  (unsigned)adj1.size(), (unsigned)adj2.size() );
584  ok = false;
585  continue;
586  }
587 
588  std::vector< int > ids1( adj1.size() ), ids2( adj2.size() );
589  unsigned i;
590  for( i = 0; i < adj1.size(); ++i )
591  {
592  rval = iface->tag_get_data( id_tag, &adj1[i], 1, &ids1[i] );
593  if( MB_TAG_NOT_FOUND == rval )
594  ids1[i] = 0;
595  else if( MB_SUCCESS != rval )
596  moab_error( "tag_get_data" );
597 
598  rval = iface->tag_get_data( id_tag, &adj2[i], 1, &ids2[i] );
599  if( MB_TAG_NOT_FOUND == rval )
600  ids2[i] = 0;
601  else if( MB_SUCCESS != rval )
602  moab_error( "tag_get_data" );
603  }
604 
605  std::sort( ids1.begin(), ids1.end() );
606  std::sort( ids2.begin(), ids2.end() );
607  for( i = 0; i < ids1.size(); ++i )
608  {
609  if( ids1[i] != ids2[i] )
610  {
611  fprintf( stderr, "Sets with id %d have non-matching %s\n", id, words[two] );
612  ok = false;
613  break;
614  }
615  }
616  }
617  return ok;
618 }
619 
620 bool compare_tags( EntityHandle dod[] )
621 {
622  Tag tag;
623 
624  // Get integer tag handle and characterstics
625  if( MB_SUCCESS != iface->tag_get_handle( intname, 2, MB_TYPE_INTEGER, tag, MB_TAG_SPARSE | MB_TAG_STORE ) )
626  moab_error( "tag_get_handle(intname)" );
627 
628  // integer tag should not have a default value
629  int idata[4];
630  if( MB_ENTITY_NOT_FOUND != iface->tag_get_default_value( tag, idata ) )
631  {
632  fprintf( stderr, "tag_get_default_value() for nonexistant default integer tag value did "
633  "not fail correctly.\n" );
634  return false;
635  }
636 
637  // check data for integer tag on both dodecahedrons
638  if( MB_SUCCESS != iface->tag_get_data( tag, dod, 2, idata ) ) moab_error( "tag_get_data(itag)" );
639  if( idata[0] != (int)0xDEADBEEF || idata[1] != (int)0xDEFACED || idata[2] != idata[0] || idata[3] != idata[1] )
640  {
641  fprintf( stderr, "Incorrect values for integer tag data.\n" );
642  return false;
643  }
644 
645  // Get double tag handle and characterstics
646  if( MB_SUCCESS != iface->tag_get_handle( dblname, 1, MB_TYPE_DOUBLE, tag, MB_TAG_DENSE | MB_TAG_STORE ) )
647  moab_error( "tag_get_handle(dblname)" );
648 
649  // check default value of double tag
650  double ddata[2];
651  if( MB_SUCCESS != iface->tag_get_default_value( tag, ddata ) ) moab_error( "tag_get_default_value" );
652  if( ddata[0] != 3.14159 )
653  {
654  fprintf( stderr, "incorrect default value for double tag.\n" );
655  return false;
656  }
657 
658  // check data for double tag on both dodecahedrons
659  if( MB_SUCCESS != iface->tag_get_data( tag, dod, 2, ddata ) ) moab_error( "tag_get_data()" );
660  if( ddata[0] != .333 || ddata[0] != ddata[1] )
661  {
662  fprintf( stderr, "Incorrect values for double tag data.\n" );
663  return false;
664  }
665 
666  // Get handle tag handle and characterstics
667  if( MB_SUCCESS != iface->tag_get_handle( handlename, 3, MB_TYPE_HANDLE, tag, MB_TAG_SPARSE | MB_TAG_STORE ) )
668  moab_error( "tag_get_handle(handlename)" );
669 
670  // check default value of handle tag
671  // default value will not change after tag is created. As we
672  // do multiple iterations of save/restore, we shouldn't expect
673  // the dodecahedron handles to point to the most recent two.
674  EntityHandle hdata[6];
675  if( MB_SUCCESS != iface->tag_get_default_value( tag, hdata ) ) moab_error( "tag_get_default_value" );
676  if( iface->type_from_handle( hdata[0] ) != MBHEX || hdata[2] != 0 ||
677  iface->type_from_handle( hdata[1] ) != MBPOLYHEDRON )
678  {
679  fprintf( stderr, "incorrect default value for handle tag '%s'\n", handlename );
680  return false;
681  }
682 
683  // check global value for handle tag
684  // global value should be changed each time a new global is read,
685  // so expect one of the two dodecahedrons in the last slot.
686  const EntityHandle root = 0;
687  if( MB_SUCCESS != iface->tag_get_data( tag, &root, 1, hdata ) ) moab_error( "tag_get_data(mesh)" );
688  if( iface->type_from_handle( hdata[1] ) != MBHEX || hdata[0] != 0 || ( hdata[2] != dod[0] && hdata[2] != dod[1] ) )
689  {
690  fprintf( stderr, "incorrect global/mesh value for handle tag.\n" );
691  return false;
692  }
693 
694  // check data on each dodecahedron
695  const EntityHandle* conn;
696  int len;
697  if( MB_SUCCESS != iface->tag_get_data( tag, dod, 2, hdata ) ) moab_error( "tag_get_data()" );
698  if( MB_SUCCESS != iface->get_connectivity( dod[0], conn, len ) ) moab_error( "get_connectivity" );
699  if( memcmp( conn, hdata, 3 * sizeof( EntityHandle ) ) )
700  {
701  fprintf( stderr, "Incorrect values for handle tag data.\n" );
702  return false;
703  }
704  if( MB_SUCCESS != iface->get_connectivity( dod[1], conn, len ) ) moab_error( "get_connectivity" );
705  if( memcmp( conn, hdata + 3, 3 * sizeof( EntityHandle ) ) )
706  {
707  fprintf( stderr, "Incorrect values for handle tag data.\n" );
708  return false;
709  }
710 
711  return true;
712 }
713 
714 bool compare()
715 {
716  Range range;
717  EntityHandle hex[2];
718  EntityHandle dod[2];
719  Tag elemtag;
720  Range::iterator iter;
721 
722  // get tag
723  if( MB_SUCCESS != iface->tag_get_handle( bitname, 2, MB_TYPE_BIT, elemtag ) ) moab_error( "tag_get_handle" );
724 
725  // get two hexes
726  char two = '\002';
727  const void* tarray[] = { &two };
728  if( MB_SUCCESS != iface->get_entities_by_type_and_tag( 0, MBHEX, &elemtag, tarray, 1, range ) )
729  moab_error( "get_entities_by_type_and_tag" );
730  if( range.size() != 2 )
731  {
732  fprintf( stderr, "Expected 2 Hexes. Got %lu\n", (unsigned long)range.size() );
733  exit( 1 );
734  }
735  iter = range.begin();
736  hex[0] = *iter;
737  hex[1] = *++iter;
738 
739  // get two polyhedra
740  range.clear();
741  char one = '\001';
742  const void* oarray[] = { &one };
743  if( MB_SUCCESS != iface->get_entities_by_type_and_tag( 0, MBPOLYHEDRON, &elemtag, oarray, 1, range ) )
744  moab_error( "get_entities_by_type_and_tag" );
745  if( range.size() != 2 )
746  {
747  fprintf( stderr, "Expected 2 Polyhedra. Got %lu\n", (unsigned long)range.size() );
748  exit( 1 );
749  }
750  iter = range.begin();
751  dod[0] = *iter;
752  dod[1] = *++iter;
753 
754  // compare hexes
755  std::vector< EntityHandle > conn[2];
756  if( MB_SUCCESS != iface->get_connectivity( hex, 1, conn[0] ) ||
757  MB_SUCCESS != iface->get_connectivity( hex + 1, 1, conn[1] ) )
758  moab_error( "get_connectivity" );
759  if( !compare_conn( conn[0], conn[1] ) ) return false;
760 
761  // compare polyhedra
762 
763  std::vector< EntityHandle > face[2];
764  conn[0].clear();
765  conn[1].clear();
766  if( MB_SUCCESS != iface->get_connectivity( dod, 1, conn[0], false ) ||
767  MB_SUCCESS != iface->get_connectivity( dod + 1, 1, conn[1], false ) )
768  moab_error( "get_connectivity" );
769  if( conn[0].size() != 12 || conn[1].size() != 12 )
770  {
771  fprintf( stderr,
772  "Expected two dodecahedrons. Got polyhedrons with "
773  "%lu and %lu faces respectively.\n",
774  (unsigned long)conn[0].size(), (unsigned long)conn[1].size() );
775  return false;
776  }
777 
778  for( int i = 0; i < 12; ++i )
779  {
780  face[0].clear();
781  face[1].clear();
782  if( MB_SUCCESS != iface->get_connectivity( &conn[0][i], 1, face[0], false ) ||
783  MB_SUCCESS != iface->get_connectivity( &conn[1][i], 1, face[1], false ) )
784  moab_error( "get_connectivity" );
785  if( !compare_conn( face[0], face[1] ) ) return false;
786  }
787 
788  // compare sets
789 
790  if( !compare_sets( VERTEX_SET_ID, tagname ) || !compare_sets( FACE_SET_ID ) || !compare_sets( REGION_SET_ID ) ||
791  !compare_sets( EMPTY_SET_ID ) || !compare_sets( SET_SET_ID, GLOBAL_ID_TAG_NAME ) )
792  return false;
793 
794  // check tags
795  if( !compare_tags( dod ) ) return false;
796 
797  return true;
798 }
799 
800 void moab_error( const char* str )
801 {
802  std::string msg;
803  fprintf( stderr, "%s() failed.\n", str );
804  if( MB_SUCCESS == iface->get_last_error( msg ) ) fprintf( stderr, "%s\n", msg.c_str() );
805  delete iface;
806  exit( 1 );
807 }