MOAB: Mesh Oriented datABase  (version 5.5.0)
h5partial.cpp
Go to the documentation of this file.
1 #include "moab/Core.hpp"
2 #include "moab/Range.hpp"
3 #include "TestRunner.hpp"
4 #include "ReadHDF5.hpp"
5 #include "MBTagConventions.hpp"
6 #include "moab/FileOptions.hpp"
7 
8 #ifdef MOAB_HAVE_MPI
9 #include "moab_mpi.h"
10 #endif
11 
12 #include <vector>
13 #include <cstdlib>
14 #include <iostream>
15 #include <algorithm>
16 #include <limits>
17 
18 using namespace moab;
19 
20 const char TEST_FILE[] = "partial.h5m";
21 #define READ_OPTS "BUFFER_SIZE=256"
22 const char ID_TAG_NAME[] = "test_id_tag";
23 
24 static void test_read_nothing_common( bool non_existant );
25 static void test_read_nodes_common( int num_read_sets, bool blocked_coordinate_io );
26 static void test_read_handle_tag_common( bool var_len );
27 
28 const int MBQUAD_INT = 20;
29 const int NUM_SETS = 10;
30 const int SET_WIDTH = ( MBQUAD_INT + NUM_SETS - 1 ) / NUM_SETS; // ceil(MBQUAD_INT/NUM_SETS)
31 const char LOGICAL_NAME[] = "logical"; // tag storing logical (i,j) coordinates
32 const char CENTROID_NAME[] = "centroid"; // tag storing position of centroid (x,y,0)
33 //! Create a regular MBQUAD_INT^2 element quad mesh with regularly
34 //! spaced coordinates in the range [1,100]. Group elements
35 //! into 10 vertical strips MBQUAD_INT/10 elements wide. Tag elements,
36 //! vertices and/or sets with ID in [1,10] stored in ID_TAG_NAME
37 //! tag. Write new mesh to TEST_FILE.
38 void create_mesh( bool create_element_sets,
39  bool create_vertex_sets,
40  bool tag_elements_with_id,
41  bool tag_vertices_with_id,
42  const char* adj_elem_tag_name = 0,
43  bool var_len_adj_elems = false );
44 // Given a list of vertices adjacent to a quad strip, identify it as one of the
45 // NUM_SETS strips of quads written by create_mesh.
46 int identify_set( Interface& mb, const Range& verts );
47 int identify_set( Interface& mb, EntityHandle set );
48 
49 static Tag check_tag( Interface& mb, const char* name, TagType storage, DataType type, int size );
50 
51 enum GatherTestMode
52 {
53  GATHER_SETS,
54  GATHER_CONTENTS,
55  GATHER_NONE
56 };
57 void test_gather_sets_common( bool contained_sets, GatherTestMode mode, bool no_parent_containing_sets = false );
58 void test_gather_sets_ranged( bool contained_sets, GatherTestMode mode, bool no_parent_containing_sets = false );
59 
60 //! Read a set containing no entities
61 void test_read_empty_set()
62 {
63  test_read_nothing_common( false );
64 }
65 
66 //! Specify ID that doesn't exist in file
67 void test_read_non_existant_set()
68 {
69  test_read_nothing_common( true );
70 }
71 
72 //! Read in the nodes contained in a set.
73 void test_read_one_set_nodes()
74 {
75  test_read_nodes_common( 1, false );
76 }
77 
78 //! Read in the nodes contained in a set.
79 void test_read_one_set_nodes_blocked()
80 {
81  test_read_nodes_common( 1, true );
82 }
83 
84 //! Read in the elems contained in a set
85 void test_read_one_set_elems();
86 
87 //! Read in the polyhedra contained in a set
88 void test_read_one_set_polyhedra();
89 
90 //! Read in the sets contained in a set.
91 //! Should read all sets containing read elements or nodes
92 //! and all sets that are contained the the specified "read"
93 //! set. Test the later here.
94 void test_read_set_sets();
95 
96 //! Read in the nodes contained in a sets.
97 void test_read_two_sets_nodes()
98 {
99  test_read_nodes_common( 2, false );
100 }
101 
102 //! Read in the elems contained in a sets
103 void test_read_two_sets_elems();
104 
105 //! For any set selected to be read by either explicit designation,
106 //! containing read entities, or contained in an explcitly designated
107 //! set, any child sets are also read. Check that here.
108 void test_read_child_sets_only()
109 {
110  test_gather_sets_common( false, GATHER_SETS );
111  test_gather_sets_ranged( false, GATHER_SETS );
112 }
113 void test_read_child_set_contents()
114 {
115  test_gather_sets_common( false, GATHER_CONTENTS );
116  test_gather_sets_ranged( false, GATHER_CONTENTS );
117 }
118 void test_read_no_child_sets()
119 {
120  test_gather_sets_common( false, GATHER_NONE );
121  test_gather_sets_ranged( false, GATHER_NONE );
122 }
123 
124 //! For any set selected to be read by either explicit designation,
125 //! containing read entities, or contained in an explcitly designated
126 //! set, any contained sets are also read. Check that here.
127 void test_read_contained_sets_only()
128 {
129  test_gather_sets_common( true, GATHER_SETS, true );
130  test_gather_sets_ranged( true, GATHER_SETS );
131 }
132 void test_read_contained_set_contents()
133 {
134  test_gather_sets_common( true, GATHER_CONTENTS, true );
135  test_gather_sets_ranged( true, GATHER_CONTENTS );
136 }
137 void test_read_no_contained_sets()
138 {
139  test_gather_sets_common( true, GATHER_NONE, true );
140  test_gather_sets_ranged( true, GATHER_NONE );
141 }
142 
143 //! Read in the sets contained in a set.
144 //! Should read all sets containing read elements or nodes
145 //! and all sets that are contained the the specified "read"
146 //! set. Test the former here.
147 void test_read_containing_sets();
148 
149 //! Test reading of explicit adjacencies
150 void test_read_adjacencies();
151 
152 //! Test reading of sparse double tag data
153 void test_read_double_tag();
154 
155 //! Test reading of sparse opaque tag data
156 void test_read_opaque_tag();
157 
158 //! Test reading of sparse handle tag data
159 void test_read_handle_tag()
160 {
161  test_read_handle_tag_common( false );
162 }
163 
164 //! Test reading of variable-length tag data
165 void test_var_len_tag()
166 {
167  test_read_handle_tag_common( true );
168 }
169 
170 void test_read_tagged_elems();
171 
172 void test_read_tagged_nodes();
173 
174 void test_read_sides();
175 
176 void test_read_ids();
177 
178 void test_read_partial_ids();
179 
180 int main( int argc, char* argv[] )
181 {
182 #ifdef MOAB_HAVE_MPI
183  int fail = MPI_Init( &argc, &argv );
184  if( fail ) return fail;
185 #endif
186 
187  REGISTER_TEST( test_read_empty_set );
188  REGISTER_TEST( test_read_non_existant_set );
189  REGISTER_TEST( test_read_one_set_nodes );
190  REGISTER_TEST( test_read_one_set_nodes_blocked );
191  REGISTER_TEST( test_read_one_set_elems );
192  REGISTER_TEST( test_read_one_set_polyhedra );
193  REGISTER_TEST( test_read_set_sets );
194  REGISTER_TEST( test_read_two_sets_nodes );
195  REGISTER_TEST( test_read_two_sets_elems );
196  REGISTER_TEST( test_read_child_sets_only );
197  REGISTER_TEST( test_read_child_set_contents );
198  REGISTER_TEST( test_read_no_child_sets );
199  REGISTER_TEST( test_read_contained_sets_only );
200  REGISTER_TEST( test_read_contained_set_contents );
201  REGISTER_TEST( test_read_no_contained_sets );
202  REGISTER_TEST( test_read_containing_sets );
203  REGISTER_TEST( test_read_double_tag );
204  REGISTER_TEST( test_read_opaque_tag );
205  REGISTER_TEST( test_read_handle_tag );
206  REGISTER_TEST( test_var_len_tag );
207  REGISTER_TEST( test_read_adjacencies );
208  REGISTER_TEST( test_read_tagged_elems );
209  REGISTER_TEST( test_read_tagged_nodes );
210  REGISTER_TEST( test_read_sides );
211  REGISTER_TEST( test_read_ids );
212  REGISTER_TEST( test_read_partial_ids );
213  int result = RUN_TESTS( argc, argv );
214 
215 #ifdef MOAB_HAVE_MPI
216  fail = MPI_Finalize();
217  if( fail ) return fail;
218 #endif
219 
220  return result;
221 }
222 
223 void test_read_nothing_common( bool non_existant )
224 {
225  ErrorCode rval;
226  Core moab;
227  Interface& mb = moab;
228 
229  // create a few nodes to write to file
230  std::vector< double > coords( 3000 );
231  Range verts;
232  rval = mb.create_vertices( &coords[0], coords.size() / 3, verts );CHECK_ERR( rval );
233 
234  // create three entity sets
235  EntityHandle sets[3];
236  rval = mb.create_meshset( MESHSET_SET, sets[0] );CHECK_ERR( rval );
237  rval = mb.create_meshset( MESHSET_SET, sets[1] );CHECK_ERR( rval );
238  rval = mb.create_meshset( MESHSET_ORDERED, sets[2] );CHECK_ERR( rval );
239 
240  // put all vertices into two of the sets
241  rval = mb.add_entities( sets[0], verts );CHECK_ERR( rval );
242  rval = mb.add_entities( sets[2], verts );CHECK_ERR( rval );
243 
244  // tag all three sets
245  Tag id_tag;
246  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
247  int ids[3] = { 5, 7, 9 };
248  rval = mb.tag_set_data( id_tag, sets, 3, ids );CHECK_ERR( rval );
249 
250  // write mesh
251  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
252  rval = mb.delete_mesh();CHECK_ERR( rval );
253 
254  // now read back in only the empty set
255  EntityHandle file_set;
256  int id = non_existant ? 8 : 7;
257  rval = mb.create_meshset( MESHSET_SET, file_set );CHECK_ERR( rval );
258  rval = mb.load_file( TEST_FILE, &file_set, READ_OPTS, ID_TAG_NAME, &id, 1 );
259  if( non_existant )
260  {
261  CHECK_EQUAL( MB_ENTITY_NOT_FOUND, rval );
262  return;
263  }
264  else
265  CHECK_ERR( rval );
266 
267  // the file should contain exactly two sets (the specified one and the new
268  // file set, and nothing else.)
269  for( EntityType t = MBVERTEX; t < MBENTITYSET; ++t )
270  {
271  int count = -1;
272  rval = mb.get_number_entities_by_type( 0, t, count );CHECK_ERR( rval );
273  CHECK_EQUAL( 0, count );
274  }
275  Range setrange;
276  rval = mb.get_entities_by_type( 0, MBENTITYSET, setrange );CHECK_ERR( rval );
277  CHECK_EQUAL( ( non_existant ? 1 : 2 ), (int)setrange.size() );
278  CHECK( setrange.find( file_set ) != setrange.end() );
279 }
280 
281 static void vtx_coords( int set_id, int j, int num_sets, double coords[3] )
282 {
283  int i = num_sets * j + set_id;
284  coords[0] = i;
285  coords[1] = i + 0.25;
286  coords[2] = i + 0.5;
287 }
288 
289 void test_read_nodes_common( int num_read_sets, bool blocked )
290 {
291  ErrorCode rval;
292  Core moab;
293  Interface& mb = moab;
294 
295  // create 1000 nodes
296  const int num_sets = 2 * num_read_sets;
297  std::vector< EntityHandle > verts( 1000 );
298  std::vector< std::vector< EntityHandle > > set_verts( num_sets );
299  for( size_t i = 0; i < verts.size(); ++i )
300  {
301  double coords[3];
302  int j = i % num_sets;
303  vtx_coords( j + 1, set_verts[j].size(), num_sets, coords );
304  rval = mb.create_vertex( coords, verts[i] );
305  set_verts[j].push_back( verts[i] );CHECK_ERR( rval );
306  }
307 
308  // create two sets, each containing half of the nodes
309  std::vector< EntityHandle > sets( num_sets );
310  for( int i = 0; i < num_sets; ++i )
311  {
312  rval = mb.create_meshset( MESHSET_ORDERED, sets[i] );CHECK_ERR( rval );
313  rval = mb.add_entities( sets[i], &set_verts[i][0], set_verts[i].size() );CHECK_ERR( rval );
314  }
315 
316  // tag both sets
317  Tag id_tag;
318  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
319  std::vector< int > values( num_sets );
320  for( int i = 0; i < num_sets; ++i )
321  values[i] = i + 1;
322  rval = mb.tag_set_data( id_tag, &sets[0], num_sets, &values[0] );CHECK_ERR( rval );
323 
324  // write file
325  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
326  rval = mb.delete_mesh();CHECK_ERR( rval );
327 
328  // now read back in only the specified number of sets
329  std::string opts( READ_OPTS );
330  if( !opts.empty() ) opts += ';';
331  if( blocked )
332  opts += "BLOCKED_COORDINATE_IO=yes";
333  else
334  opts += "BLOCKED_COORDINATE_IO=no";
335 
336  values.resize( num_read_sets );
337  for( int i = 0; i < num_read_sets; ++i )
338  values[i] = 2 * ( i + 1 );
339  EntityHandle file_set;
340  rval = mb.create_meshset( MESHSET_SET, file_set );CHECK_ERR( rval );
341  rval = mb.load_file( TEST_FILE, &file_set, opts.c_str(), ID_TAG_NAME, &values[0], num_read_sets );CHECK_ERR( rval );
342 
343  int count, expected = 0;
344  rval = mb.get_number_entities_by_dimension( 0, 0, count );CHECK_ERR( rval );
345  for( int i = 0; i < num_sets; ++i )
346  if( i % 2 ) expected += set_verts[i].size();
347  CHECK_EQUAL( expected, count );
348 
349  Range sets2;
350  rval = mb.get_entities_by_type( 0, MBENTITYSET, sets2 );CHECK_ERR( rval );
351  CHECK_EQUAL( 1 + num_read_sets, (int)sets2.size() );
352  Range::iterator it = sets2.find( file_set );
353  CHECK( it != sets2.end() );
354  sets2.erase( it );
355 
356  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag );CHECK_ERR( rval );
357  while( !sets2.empty() )
358  {
359  EntityHandle set = sets2.pop_front();
360  int id;
361  rval = mb.tag_get_data( id_tag, &set, 1, &id );CHECK_ERR( rval );
362  CHECK( std::find( values.begin(), values.end(), id ) != values.end() );
363  CHECK( id > 0 );
364  CHECK( (unsigned)id <= set_verts.size() );
365 
366  std::vector< EntityHandle > verts2;
367  rval = mb.get_entities_by_handle( set, verts2 );CHECK_ERR( rval );
368  CHECK_EQUAL( set_verts[id - 1].size(), verts2.size() );
369 
370  for( size_t i = 0; i < verts2.size(); ++i )
371  {
372  double exp_coords[3], coords[3];
373  vtx_coords( id, i, num_sets, exp_coords );
374  rval = mb.get_coords( &verts2[i], 1, coords );CHECK_ERR( rval );
375  CHECK_REAL_EQUAL( exp_coords[0], coords[0], 1e-12 );
376  CHECK_REAL_EQUAL( exp_coords[1], coords[1], 1e-12 );
377  CHECK_REAL_EQUAL( exp_coords[2], coords[2], 1e-12 );
378  }
379  }
380 }
381 
382 //! Create a regular MBQUAD_INT^2 element quad mesh with regularly
383 //! spaced coordinates in the range [1,100]. Group elements
384 //! into 10 vertical strips MBQUAD_INT/10 elements wide. Tag elements,
385 //! vertices and/or sets with ID in [1,10] stored in ID_TAG_NAME
386 //! tag. Write new mesh to TEST_FILE.
387 void create_mesh( bool create_element_sets,
388  bool create_vertex_sets,
389  bool tag_elements_with_id,
390  bool tag_vertices_with_id,
391  const char* adj_elem_tag_name,
392  bool var_len_adj_elems )
393 {
394  Core moab;
395  Interface& mb = moab;
396  ErrorCode rval;
397 
398  // create tags
399  Tag logical_tag, centroid_tag, id_tag;
400  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
401  rval =
402  mb.tag_get_handle( LOGICAL_NAME, 2 * sizeof( int ), MB_TYPE_OPAQUE, logical_tag, MB_TAG_DENSE | MB_TAG_EXCL );CHECK_ERR( rval );
403  rval = mb.tag_get_handle( CENTROID_NAME, 3, MB_TYPE_DOUBLE, centroid_tag, MB_TAG_DENSE | MB_TAG_EXCL );CHECK_ERR( rval );
404 
405  EntityHandle sets[NUM_SETS];
406  if( create_element_sets || create_vertex_sets )
407  {
408  for( int i = 0; i < NUM_SETS; ++i )
409  {
410  rval = mb.create_meshset( MESHSET_ORDERED, sets[i] );CHECK_ERR( rval );
411  int id = i + 1;
412  rval = mb.tag_set_data( id_tag, &sets[i], 1, &id );CHECK_ERR( rval );
413  }
414  }
415 
416  // create elements
417  EntityHandle verts[MBQUAD_INT + 1][MBQUAD_INT + 1], quads[MBQUAD_INT][MBQUAD_INT];
418  for( int i = 0; i <= MBQUAD_INT; ++i )
419  for( int j = 0; j <= MBQUAD_INT; ++j )
420  {
421  double coords[3] = { static_cast< double >( i ), static_cast< double >( j ), 0 };
422  rval = mb.create_vertex( coords, verts[j][i] );CHECK_ERR( rval );
423  int logical[2] = { i, j };
424  rval = mb.tag_set_data( logical_tag, &verts[j][i], 1, logical );CHECK_ERR( rval );
425  rval = mb.tag_set_data( centroid_tag, &verts[j][i], 1, coords );CHECK_ERR( rval );
426  int id = ( i - 1 ) / SET_WIDTH + 1; // Note: assumes SET_WIDTH > 1
427  if( tag_vertices_with_id )
428  {
429  rval = mb.tag_set_data( id_tag, &verts[j][i], 1, &id );CHECK_ERR( rval );
430  }
431  if( create_vertex_sets )
432  {
433  rval = mb.add_entities( sets[id - 1], &verts[j][i], 1 );CHECK_ERR( rval );
434  // Some vertices are shared by quads in different sets.
435  // put such vertices in both sets.
436  int id2 = i / SET_WIDTH + 1;
437  if( id2 != id && id2 <= NUM_SETS )
438  {
439  rval = mb.add_entities( sets[id2 - 1], &verts[j][i], 1 );CHECK_ERR( rval );
440  }
441  }
442  }
443  for( int i = 0; i < MBQUAD_INT; ++i )
444  for( int j = 0; j < MBQUAD_INT; ++j )
445  {
446  EntityHandle conn[4] = { verts[j][i], verts[j][i + 1], verts[j + 1][i + 1], verts[j + 1][i] };
447  rval = mb.create_element( MBQUAD, conn, 4, quads[j][i] );CHECK_ERR( rval );
448  int logical[2] = { i, j };
449  rval = mb.tag_set_data( logical_tag, &quads[j][i], 1, logical );CHECK_ERR( rval );
450  double centroid[3] = { i + 0.5, j + 0.5, 0 };
451  rval = mb.tag_set_data( centroid_tag, &quads[j][i], 1, centroid );CHECK_ERR( rval );
452  int id = i / SET_WIDTH + 1;
453  if( tag_elements_with_id )
454  {
455  rval = mb.tag_set_data( id_tag, &quads[j][i], 1, &id );CHECK_ERR( rval );
456  }
457  if( create_element_sets )
458  {
459  rval = mb.add_entities( sets[id - 1], &quads[j][i], 1 );CHECK_ERR( rval );
460  }
461  }
462 
463  if( adj_elem_tag_name && !var_len_adj_elems )
464  {
465  Tag handle_tag;
466  rval = mb.tag_get_handle( adj_elem_tag_name, 4, MB_TYPE_HANDLE, handle_tag, MB_TAG_DENSE | MB_TAG_EXCL );CHECK_ERR( rval );
467  for( int i = 0; i <= MBQUAD_INT; ++i )
468  for( int j = 0; j <= MBQUAD_INT; ++j )
469  {
470  EntityHandle val[4] = { ( i > 0 && j > 0 ) ? quads[j - 1][i - 1] : 0,
471  ( i > 0 && j < MBQUAD_INT ) ? quads[j][i - 1] : 0,
472  ( i < MBQUAD_INT && j < MBQUAD_INT ) ? quads[j][i] : 0,
473  ( i < MBQUAD_INT && j > 0 ) ? quads[j - 1][i] : 0 };
474  rval = mb.tag_set_data( handle_tag, &verts[j][i], 1, val );CHECK_ERR( rval );
475  }
476  }
477  else if( adj_elem_tag_name && var_len_adj_elems )
478  {
479  Tag handle_tag;
480  rval = mb.tag_get_handle( adj_elem_tag_name, 0, MB_TYPE_HANDLE, handle_tag,
481  MB_TAG_DENSE | MB_TAG_VARLEN | MB_TAG_EXCL );CHECK_ERR( rval );
482  for( int i = 0; i <= MBQUAD_INT; ++i )
483  for( int j = 0; j <= MBQUAD_INT; ++j )
484  {
485  EntityHandle val[4];
486  int num = 0;
487  if( i > 0 && j > 0 ) val[num++] = quads[j - 1][i - 1];
488  if( i > 0 && j < MBQUAD_INT ) val[num++] = quads[j][i - 1];
489  if( i < MBQUAD_INT && j < MBQUAD_INT ) val[num++] = quads[j][i];
490  if( i < MBQUAD_INT && j > 0 ) val[num++] = quads[j - 1][i];
491  const void* ptr = val;
492  rval = mb.tag_set_by_ptr( handle_tag, &verts[j][i], 1, &ptr, &num );CHECK_ERR( rval );
493  }
494  }
495 
496  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
497 }
498 
499 // Given a list of vertices adjacent to a quad strip, identify it as one of the
500 // NUM_SETS strips of quads written by create_mesh.
501 int identify_set( Interface& mb, const Range& verts )
502 {
503  const int COL = SET_WIDTH + 1;
504  CHECK_EQUAL( ( 1 + MBQUAD_INT ) * COL, (int)verts.size() );
505 
506  // Get X range of vertices
507  int min_x = std::numeric_limits< int >::max();
508  int max_x = std::numeric_limits< int >::min();
509  for( Range::const_iterator i = verts.begin(); i != verts.end(); ++i )
510  {
511  double coords[3];
512  ErrorCode rval = mb.get_coords( &*i, 1, coords );CHECK_ERR( rval );
513  // Expect whole-valued coorindates
514  int int_x = (int)coords[0];
515  CHECK( fabs( coords[0] - (double)int_x ) < 1e-12 );
516 
517  if( int_x < min_x ) min_x = int_x;
518  if( int_x > max_x ) max_x = int_x;
519  }
520  CHECK( max_x - min_x < COL );
521 
522  // Calculate ID (return value) from coordinate range)
523  const int ID = min_x / SET_WIDTH + 1;
524 
525  // Now verify that all vertices correctly form a grid
526  EntityHandle grid[MBQUAD_INT + 1][COL];
527  memset( grid, 0, sizeof( grid ) );
528  for( Range::const_iterator i = verts.begin(); i != verts.end(); ++i )
529  {
530  double coords[3];
531  ErrorCode rval = mb.get_coords( &*i, 1, coords );CHECK_ERR( rval );
532  // Expect whole-valued coorindates
533  int x = (int)coords[0] - ( ID - 1 ) * SET_WIDTH, y = (int)coords[1];
534  CHECK( fabs( coords[1] - (double)y ) < 1e-12 );
535  CHECK( fabs( coords[2] ) < 1e-12 );
536  CHECK( y >= 0 && y <= MBQUAD_INT );
537  CHECK_EQUAL( (EntityHandle)0, grid[y][x] );
538  grid[y][x] = *i;
539  }
540 
541  return ID;
542 }
543 int identify_set( Interface& mb, EntityHandle set )
544 {
545  ErrorCode rval;
546  Range verts, elems;
547  rval = mb.get_entities_by_handle( set, elems );CHECK_ERR( rval );
548  Range::iterator it = elems.upper_bound( MBVERTEX );
549  verts.merge( elems.begin(), it );
550  elems.erase( elems.begin(), it );
551  it = elems.lower_bound( MBENTITYSET );
552  elems.erase( it, elems.end() );
553  rval = mb.get_adjacencies( elems, 0, false, verts, Interface::UNION );CHECK_ERR( rval );
554  return identify_set( mb, verts );
555 }
556 
557 //! Read in the elems contained in a set
558 void test_read_one_set_elems()
559 {
560  ErrorCode rval;
561  Core moab;
562  Interface& mb = moab;
563 
564  create_mesh( true, false, false, false );
565 
566  for( int id = 1; id <= NUM_SETS; ++id )
567  {
568  rval = mb.delete_mesh();CHECK_ERR( rval );
569  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, &id, 1 );CHECK_ERR( rval );
570  Range verts;
571  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );
572  int act_id = identify_set( mb, verts );
573  CHECK_EQUAL( id, act_id );
574  }
575 }
576 
577 //! Read in the elems contained in a sets
578 void test_read_two_sets_elems()
579 {
580  ErrorCode rval;
581  Core moab;
582  Interface& mb = moab;
583 
584  create_mesh( true, false, false, false );
585  int ids[2] = { 2, 8 };
586  EntityHandle file_set;
587  rval = mb.create_meshset( MESHSET_SET, file_set );CHECK_ERR( rval );
588  rval = mb.load_file( TEST_FILE, &file_set, READ_OPTS, ID_TAG_NAME, ids, 2 );CHECK_ERR( rval );
589 
590  Range sets;
591  rval = mb.get_entities_by_type( 0, MBENTITYSET, sets );CHECK_ERR( rval );
592  CHECK_EQUAL( 3, (int)sets.size() );
593  Range::iterator it = sets.find( file_set );
594  CHECK( it != sets.end() );
595  sets.erase( it );
596 
597  int id1 = identify_set( mb, sets.front() );
598  int id2 = identify_set( mb, sets.back() );
599  if( id1 == ids[0] )
600  {
601  CHECK_EQUAL( ids[1], id2 );
602  }
603  else
604  {
605  CHECK_EQUAL( ids[1], id1 );
606  CHECK_EQUAL( ids[0], id2 );
607  }
608 }
609 
610 Tag check_tag( Interface& mb, const char* name, TagType storage, DataType type, int size )
611 {
612 
613  Tag tag;
614  ErrorCode rval = mb.tag_get_handle( name, size, type, tag );CHECK_ERR( rval );
615 
616  TagType storage1;
617  rval = mb.tag_get_type( tag, storage1 );CHECK_ERR( rval );
618  CHECK_EQUAL( storage, storage1 );
619 
620  DataType type1;
621  rval = mb.tag_get_data_type( tag, type1 );CHECK_ERR( rval );
622  CHECK_EQUAL( type, type1 );
623 
624  int size1;
625  rval = mb.tag_get_length( tag, size1 );
626  if( size <= 0 )
627  { // variable-length tag
628  CHECK_EQUAL( MB_VARIABLE_DATA_LENGTH, rval );
629  }
630  else
631  {
632  CHECK_ERR( rval );
633  CHECK_EQUAL( size, size1 );
634  }
635 
636  return tag;
637 }
638 
639 //! Test reading of sparse double tag data
640 void test_read_double_tag()
641 {
642  ErrorCode rval;
643  Core moab;
644  Interface& mb = moab;
645 
646  create_mesh( true, false, false, false );
647  int ids[2] = { 1, 4 };
648  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 2 );CHECK_ERR( rval );
649 
650  Tag tag = check_tag( mb, CENTROID_NAME, MB_TAG_DENSE, MB_TYPE_DOUBLE, 3 );
651  Range verts;
652  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
653  CHECK( !verts.empty() );
654  for( Range::iterator i = verts.begin(); i != verts.end(); ++i )
655  {
656  double coords[3], data[3];
657  rval = mb.get_coords( &*i, 1, coords );CHECK_ERR( rval );
658  rval = mb.tag_get_data( tag, &*i, 1, data );CHECK_ERR( rval );
659  CHECK_REAL_EQUAL( coords[0], data[0], 1e-12 );
660  CHECK_REAL_EQUAL( coords[1], data[1], 1e-12 );
661  CHECK_REAL_EQUAL( coords[2], data[2], 1e-12 );
662  }
663 }
664 
665 //! Test reading of sparse opaque tag data
666 void test_read_opaque_tag()
667 {
668  ErrorCode rval;
669  Core moab;
670  Interface& mb = moab;
671 
672  create_mesh( true, false, false, false );
673  int ids[2] = { 1, 4 };
674  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 2 );CHECK_ERR( rval );
675 
676  Tag tag = check_tag( mb, LOGICAL_NAME, MB_TAG_DENSE, MB_TYPE_OPAQUE, 2 * sizeof( int ) );
677  Range verts;
678  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
679  CHECK( !verts.empty() );
680  for( Range::iterator i = verts.begin(); i != verts.end(); ++i )
681  {
682  double coords[3];
683  int data[2];
684  rval = mb.get_coords( &*i, 1, coords );CHECK_ERR( rval );
685  rval = mb.tag_get_data( tag, &*i, 1, data );CHECK_ERR( rval );
686  CHECK_REAL_EQUAL( coords[0], (double)data[0], 1e-12 );
687  CHECK_REAL_EQUAL( coords[1], (double)data[1], 1e-12 );
688  }
689 }
690 
691 static void test_read_handle_tag_common( bool var_len )
692 {
693  ErrorCode rval;
694  Core moab;
695  Interface& mb = moab;
696 
697  const char tag_name[] = "VTX_ADJ";
698  create_mesh( true, false, false, false, tag_name, var_len );
699  int ids[2] = { 7, 10 };
700  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 2 );CHECK_ERR( rval );
701 
702  Tag tag = check_tag( mb, tag_name, MB_TAG_DENSE, MB_TYPE_HANDLE, var_len ? 0 : 4 );
703  Range verts;
704  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
705  CHECK( !verts.empty() );
706  for( Range::iterator i = verts.begin(); i != verts.end(); ++i )
707  {
708  std::vector< EntityHandle > adj, val;
709  rval = mb.get_adjacencies( &*i, 1, 2, false, adj, Interface::UNION );CHECK_ERR( rval );
710  CHECK( !adj.empty() );
711 
712  int num;
713  const void* ptr;
714  rval = mb.tag_get_by_ptr( tag, &*i, 1, &ptr, &num );CHECK_ERR( rval );
715 
716  if( var_len )
717  {
718  CHECK( num > 0 );
719  CHECK( num < 5 );
720  }
721  else
722  {
723  CHECK_EQUAL( 4, num );
724  }
725 
726  val.clear();
727  const EntityHandle* dat = (const EntityHandle*)ptr;
728  for( const EntityHandle* end = dat + num; dat != end; ++dat )
729  if( *dat ) val.push_back( *dat );
730 
731  CHECK_EQUAL( adj.size(), val.size() );
732  std::sort( adj.begin(), adj.end() );
733  std::sort( val.begin(), val.end() );
734  CHECK( adj == val );
735  }
736 }
737 
738 void test_read_tagged_elems()
739 {
740  ErrorCode rval;
741  Core moab;
742  Interface& mb = moab;
743 
744  create_mesh( false, false, true, false );
745  int id = 5;
746  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, &id, 1 );CHECK_ERR( rval );
747 
748  Range verts;
749  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
750  int id2 = identify_set( mb, verts );
751  CHECK_EQUAL( id, id2 );
752 
753  int elems;
754  rval = mb.get_number_entities_by_type( 0, MBQUAD, elems );CHECK_ERR( rval );
755  CHECK_EQUAL( MBQUAD_INT * MBQUAD_INT / NUM_SETS, elems );
756 }
757 
758 void test_read_tagged_nodes()
759 {
760  ErrorCode rval;
761  Core moab;
762  Interface& mb = moab;
763 
764  create_mesh( false, false, false, true );
765  int id = 1; // NOTE: this test will only succeed for ID == 1
766  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, &id, 1 );CHECK_ERR( rval );
767 
768  Range verts;
769  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
770  int id2 = identify_set( mb, verts );
771  CHECK_EQUAL( id, id2 );
772 
773  int elems;
774  rval = mb.get_number_entities_by_type( 0, MBQUAD, elems );CHECK_ERR( rval );
775  CHECK_EQUAL( MBQUAD_INT * MBQUAD_INT / NUM_SETS, elems );
776 }
777 
778 //! Read in the polyhedra contained in a set
779 void test_read_one_set_polyhedra()
780 {
781  ErrorCode rval;
782  Core instance;
783  Interface& mb = instance;
784 
785  // create a 2x2x1 block of hexes, splitting each hex face
786  // into two triangles to form an 12-sided polyhedron
787  EntityHandle verts[18], hexes[4];
788  double coords[18][3] = { { 0, 0, 0 }, { 1, 0, 0 }, { 2, 0, 0 }, { 0, 1, 0 }, { 1, 1, 0 }, { 2, 1, 0 },
789  { 0, 0, 1 }, { 1, 0, 1 }, { 2, 0, 1 }, { 0, 1, 1 }, { 1, 1, 1 }, { 2, 1, 1 },
790  { 0, 0, 2 }, { 1, 0, 2 }, { 2, 0, 2 }, { 0, 1, 2 }, { 1, 1, 2 }, { 2, 1, 2 } };
791  int hexconn[4][8] = { { 0, 1, 4, 3, 6, 7, 10, 9 },
792  { 1, 2, 5, 4, 7, 8, 11, 10 },
793  { 6, 7, 10, 9, 12, 13, 16, 15 },
794  { 7, 8, 11, 10, 13, 14, 17, 16 } };
795  for( int i = 0; i < 18; ++i )
796  {
797  rval = mb.create_vertex( coords[i], verts[i] );CHECK_ERR( rval );
798  }
799  for( int i = 0; i < 4; ++i )
800  {
801  EntityHandle conn[8];
802  for( int j = 0; j < 8; ++j )
803  conn[j] = verts[hexconn[i][j]];
804  rval = mb.create_element( MBHEX, conn, 8, hexes[i] );CHECK_ERR( rval );
805  }
806 
807  Tag tri_tag;
808  rval = mb.tag_get_handle( "tris", 2, MB_TYPE_HANDLE, tri_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
809 
810  std::vector< EntityHandle > quads;
811  EntityHandle tris[12], poly[4];
812  for( int i = 0; i < 4; ++i )
813  {
814  quads.clear();
815  rval = mb.get_adjacencies( &hexes[i], 1, 2, true, quads );CHECK_ERR( rval );
816  CHECK_EQUAL( 6, (int)quads.size() );
817 
818  for( int j = 0; j < 6; ++j )
819  {
820  rval = mb.tag_get_data( tri_tag, &quads[j], 1, tris + 2 * j );
821  if( MB_SUCCESS == rval ) continue;
822  CHECK_EQUAL( MB_TAG_NOT_FOUND, rval );
823  const EntityHandle* conn;
824  int len;
825  rval = mb.get_connectivity( quads[j], conn, len );CHECK_ERR( rval );
826  CHECK_EQUAL( 4, len );
827  EntityHandle tri_conn[2][3] = { { conn[0], conn[1], conn[2] }, { conn[2], conn[3], conn[0] } };
828  rval = mb.create_element( MBTRI, tri_conn[0], 3, tris[2 * j] );CHECK_ERR( rval );
829  rval = mb.create_element( MBTRI, tri_conn[1], 3, tris[2 * j + 1] );CHECK_ERR( rval );
830  rval = mb.tag_set_data( tri_tag, &quads[j], 1, tris + 2 * j );CHECK_ERR( rval );
831  }
832 
833  rval = mb.create_element( MBPOLYHEDRON, tris, 12, poly[i] );CHECK_ERR( rval );
834  }
835 
836  Range all_tri;
837  rval = mb.get_entities_by_type( 0, MBTRI, all_tri );CHECK_ERR( rval );
838  CHECK_EQUAL( 40, (int)all_tri.size() );
839 
840  rval = mb.delete_entities( hexes, 4 );CHECK_ERR( rval );
841  rval = mb.delete_entities( &quads[0], quads.size() );CHECK_ERR( rval );
842 
843  EntityHandle sets[2];
844  rval = mb.create_meshset( 0, sets[0] );CHECK_ERR( rval );
845  rval = mb.add_entities( sets[0], poly, 2 );CHECK_ERR( rval );
846  rval = mb.create_meshset( 0, sets[1] );CHECK_ERR( rval );
847  rval = mb.add_entities( sets[1], poly + 2, 2 );CHECK_ERR( rval );
848 
849  Tag id_tag;
850  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
851  int ids[2] = { 2, 3 };
852  rval = mb.tag_set_data( id_tag, sets, 2, ids );CHECK_ERR( rval );
853 
854  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
855  rval = mb.delete_mesh();CHECK_ERR( rval );
856 
857  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 1 );CHECK_ERR( rval );
858 
859  Range rpoly;
860  rval = mb.get_entities_by_type( 0, MBPOLYHEDRON, rpoly );CHECK_ERR( rval );
861  CHECK_EQUAL( 2, (int)rpoly.size() );
862 
863  Range polyverts;
864  rval = mb.get_adjacencies( rpoly, 0, false, polyverts, Interface::UNION );CHECK_ERR( rval );
865  CHECK_EQUAL( 12, (int)polyverts.size() );
866 
867  for( Range::iterator it = polyverts.begin(); it != polyverts.end(); ++it )
868  {
869  double coords2[3];
870  rval = mb.get_coords( &*it, 1, coords2 );CHECK_ERR( rval );
871  CHECK( coords2[0] > -1e-12 && coords2[0] - 2 < 1e-12 );
872  CHECK( coords2[1] > -1e-12 && coords2[1] - 1 < 1e-12 );
873  CHECK( coords2[2] > -1e-12 && coords2[2] - 1 < 1e-12 );
874  }
875 }
876 
877 //! Read in the sets contained in a set.
878 //! Should read all sets containing read elements or nodes
879 //! and all sets that are contained the the specified "read"
880 //! set. Test the later here.
881 void test_read_set_sets()
882 {
883  ErrorCode rval;
884  Core instance;
885  Interface& mb = instance;
886 
887  Tag id_tag;
888  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
889 
890  // create sets and assign an ID to each
891  const int len = 5;
892  EntityHandle set[2 * len + 2];
893  for( int i = 0; i < 2 * len + 2; ++i )
894  {
895  rval = mb.create_meshset( MESHSET_SET, set[i] );CHECK_ERR( rval );
896  int id = i + 1;
897  rval = mb.tag_set_data( id_tag, set + i, 1, &id );CHECK_ERR( rval );
898  }
899 
900  // make set containment as follows (values are assigned IDs):
901  int cont_ids[2][len] = { { 3, 4, 5, 9, 10 }, { 6, 7, 8, 11, 12 } };
902  for( int i = 0; i < 2; ++i )
903  {
904  EntityHandle contents[len] = { set[cont_ids[i][0] - 1], set[cont_ids[i][1] - 1], set[cont_ids[i][2] - 1],
905  set[cont_ids[i][3] - 1], set[cont_ids[i][4] - 1] };
906  rval = mb.add_entities( set[i], contents, len );CHECK_ERR( rval );
907  }
908 
909  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
910 
911  for( int i = 0; i < 2; ++i )
912  {
913  rval = mb.delete_mesh();CHECK_ERR( rval );
914 
915  EntityHandle file;
916  rval = mb.create_meshset( MESHSET_SET, file );CHECK_ERR( rval );
917  int id = i + 1;
918  rval = mb.load_file( TEST_FILE, &file, READ_OPTS ";SETS=NONE", ID_TAG_NAME, &id, 1 );CHECK_ERR( rval );
919 
920  // check that the total number of sets read is as expected
921  Range sets;
922  rval = mb.get_entities_by_type( 0, MBENTITYSET, sets );CHECK_ERR( rval );
923  Range::iterator it = sets.find( file );
924  if( it != sets.end() ) sets.erase( it );
925  CHECK_EQUAL( len + 1, (int)sets.size() );
926 
927  // check that we read in the set specified by ID to the reader
928  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag );CHECK_ERR( rval );
929  sets.clear();
930  const void* data[] = { &id };
931  rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &id_tag, data, 1, sets );CHECK_ERR( rval );
932  CHECK_EQUAL( 1, (int)sets.size() );
933 
934  // check that it contains the expected sets
935  EntityHandle owner = sets.front();
936  sets.clear();
937  rval = mb.get_entities_by_type( owner, MBENTITYSET, sets );CHECK_ERR( rval );
938  CHECK_EQUAL( len, (int)sets.size() );
939 
940  std::vector< int > expected( cont_ids[i], cont_ids[i] + len );
941  std::vector< int > actual( len );
942  rval = mb.tag_get_data( id_tag, sets, &actual[0] );CHECK_ERR( rval );
943  std::sort( expected.begin(), expected.end() );
944  std::sort( actual.begin(), actual.end() );
945  CHECK( expected == actual );
946  }
947 }
948 
949 static void check_children( bool contents, GatherTestMode mode, Interface& mb, int id, Tag id_tag, EntityHandle file )
950 {
951  // Increase number of expected sets by one if contents is true because
952  // we always read immediately contained (depth 1) sets.
953  const int exp_num_sets = ( mode == GATHER_NONE ) ? 1 + contents : id;
954  const int exp_num_edges = ( mode == GATHER_CONTENTS ) ? id : 1;
955 
956  ErrorCode rval;
957  Range range;
958  rval = mb.get_entities_by_type( 0, MBEDGE, range );CHECK_ERR( rval );
959  CHECK_EQUAL( exp_num_edges, (int)range.size() );
960  range.clear();
961  rval = mb.get_entities_by_type( 0, MBENTITYSET, range );CHECK_ERR( rval );
962  Range::iterator it = range.find( file );
963  CHECK( it != range.end() );
964  range.erase( it );
965  CHECK_EQUAL( exp_num_sets, (int)range.size() );
966 
967  EntityHandle set;
968  const void* val[] = { &id };
969  range.clear();
970  rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &id_tag, val, 1, range );CHECK_ERR( rval );
971  CHECK_EQUAL( 1, (int)range.size() );
972  set = range.front();
973 
974  if( mode == GATHER_NONE )
975  {
976  range.clear();
977  rval = mb.get_entities_by_type( set, MBEDGE, range );CHECK_ERR( rval );
978  CHECK_EQUAL( 1, (int)range.size() );
979  return;
980  }
981 
982  for( int i = id; i > 0; --i )
983  {
984  int act_id;
985  rval = mb.tag_get_data( id_tag, &set, 1, &act_id );CHECK_ERR( rval );
986  CHECK_EQUAL( i, act_id );
987 
988  range.clear();
989  rval = mb.get_entities_by_type( set, MBEDGE, range );CHECK_ERR( rval );
990  if( mode == GATHER_CONTENTS || i == id )
991  {
992  CHECK_EQUAL( 1, (int)range.size() );
993  const EntityHandle* conn;
994  int len;
995  rval = mb.get_connectivity( range.front(), conn, len );CHECK_ERR( rval );
996  CHECK_EQUAL( 2, len );
997  double coords[3];
998  rval = mb.get_coords( conn + 1, 1, coords );CHECK_ERR( rval );
999  CHECK_EQUAL( i, (int)coords[0] );
1000  }
1001  else
1002  {
1003  CHECK( range.empty() );
1004  }
1005 
1006  std::vector< EntityHandle > children;
1007  if( contents )
1008  rval = mb.get_entities_by_type( set, MBENTITYSET, children );
1009  else
1010  rval = mb.get_child_meshsets( set, children );CHECK_ERR( rval );
1011  if( i == 1 )
1012  {
1013  CHECK( children.empty() );
1014  }
1015  else
1016  {
1017  CHECK_EQUAL( 1, (int)children.size() );
1018  set = children[0];
1019  }
1020  }
1021 }
1022 
1023 const char* set_read_opts[] = { "SETS", "CONTENTS", "NONE" };
1024 void test_gather_sets_common( bool contents, GatherTestMode mode, bool no_parent_containing_sets )
1025 {
1026  ErrorCode rval;
1027  Core instance;
1028  Interface& mb = instance;
1029 
1030  Tag id_tag;
1031  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
1032 
1033  // Create a string of edges from [0,INT] along the X axis each 1 unit in length.
1034  // Create a set for edge edge containing the edge and make it the parent of the
1035  // set containing the previous (closer to origin) edge. Assign each set an
1036  // ID that is the X coordinate of the larger of the two vertices of the edge
1037  // contained in the set.
1038  const int INT = 64;
1039  EntityHandle verts[INT + 1], edges[INT], sets[INT];
1040  double coords[] = { 0, 0, 0 };
1041  rval = mb.create_vertex( coords, verts[0] );CHECK_ERR( rval );
1042  for( int i = 0; i < INT; ++i )
1043  {
1044  const int id = i + 1;
1045  coords[0] = id;
1046  rval = mb.create_vertex( coords, verts[id] );CHECK_ERR( rval );
1047  rval = mb.create_element( MBEDGE, verts + i, 2, edges[i] );CHECK_ERR( rval );
1048  rval = mb.create_meshset( MESHSET_SET, sets[i] );CHECK_ERR( rval );
1049  rval = mb.add_entities( sets[i], edges + i, 1 );CHECK_ERR( rval );
1050  rval = mb.tag_set_data( id_tag, sets + i, 1, &id );CHECK_ERR( rval );
1051  if( i > 0 )
1052  {
1053  if( contents )
1054  rval = mb.add_entities( sets[i], sets + ( i - 1 ), 1 );
1055  else
1056  rval = mb.add_child_meshset( sets[i], sets[i - 1] );CHECK_ERR( rval );
1057  }
1058  }
1059 
1060  // Write the data
1061  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
1062 
1063  EntityHandle file;
1064  std::string opt( READ_OPTS );
1065  if( contents )
1066  opt += ";CHILDREN=NONE;SETS=";
1067  else
1068  opt += ";SETS=NONE;CHILDREN=";
1069  opt += set_read_opts[mode];
1070 
1071  if( no_parent_containing_sets ) opt += ";NO_SET_CONTAINING_PARENTS";
1072 
1073  const int test_ids[] = { 2, 7, INT / 3 - 1, INT / 2 + 1, INT - 3 };
1074  const int num_test_ids = sizeof( test_ids ) / sizeof( int );
1075  for( int i = 0; i < num_test_ids; ++i )
1076  {
1077  CHECK( test_ids[i] <= INT );
1078 
1079  rval = mb.delete_mesh();CHECK_ERR( rval );
1080 
1081  rval = mb.create_meshset( MESHSET_SET, file );CHECK_ERR( rval );
1082 
1083  rval = mb.load_file( TEST_FILE, 0, opt.c_str(), ID_TAG_NAME, test_ids + i, 1 );CHECK_ERR( rval );
1084  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag );CHECK_ERR( rval );
1085 
1086  check_children( contents, mode, mb, test_ids[i], id_tag, file );
1087  }
1088 }
1089 
1090 void test_gather_sets_ranged( bool contents, GatherTestMode mode, bool no_parent_containing_sets )
1091 {
1092  ErrorCode rval;
1093  Core instance;
1094  Interface& mb = instance;
1095 
1096  Range verts;
1097  Tag id_tag;
1098  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
1099 
1100  // create four groups of vertices, where all vertices in the same group
1101  // have the same x-coordinate
1102  const int NUM_GRP_VTX = 20;
1103  const int NUM_GRP = 4;
1104  EntityHandle sets[NUM_GRP];
1105  for( int i = 0; i < NUM_GRP; ++i )
1106  {
1107  double coords[3 * NUM_GRP_VTX];
1108  for( int j = 0; j < NUM_GRP_VTX; ++j )
1109  {
1110  coords[3 * j] = i;
1111  coords[3 * j + 1] = j;
1112  coords[3 * j + 2] = 0;
1113  }
1114  rval = mb.create_vertices( coords, NUM_GRP_VTX, verts );CHECK_ERR( rval );
1115 
1116  rval = mb.create_meshset( MESHSET_SET, sets[i] );CHECK_ERR( rval );
1117  rval = mb.add_entities( sets[i], verts );CHECK_ERR( rval );
1118  int id = i + 1;
1119  rval = mb.tag_set_data( id_tag, sets + i, 1, &id );CHECK_ERR( rval );
1120  }
1121 
1122  // place two of the sets inside the others
1123  if( contents )
1124  {
1125  rval = mb.add_entities( sets[0], &sets[1], 1 );CHECK_ERR( rval );
1126  rval = mb.add_entities( sets[2], &sets[3], 1 );CHECK_ERR( rval );
1127  }
1128  else
1129  {
1130  rval = mb.add_child_meshset( sets[0], sets[1] );CHECK_ERR( rval );
1131  rval = mb.add_child_meshset( sets[2], sets[3] );CHECK_ERR( rval );
1132  }
1133 
1134  // Write the data
1135  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
1136 
1137  // Read the data
1138  std::string opt( READ_OPTS );
1139  if( contents )
1140  opt += ";CHILDREN=NONE;SETS=";
1141  else
1142  opt += ";SETS=NONE;CHILDREN=";
1143  opt += set_read_opts[mode];
1144 
1145  if( no_parent_containing_sets ) opt += ";NO_PARENT_CONTAINING_SETS";
1146 
1147  EntityHandle file;
1148  const int read_id = 3;
1149  rval = mb.delete_mesh();CHECK_ERR( rval );
1150  rval = mb.create_meshset( MESHSET_SET, file );CHECK_ERR( rval );
1151  rval = mb.load_file( TEST_FILE, &file, opt.c_str(), ID_TAG_NAME, &read_id, 1 );CHECK_ERR( rval );
1152 
1153  // get any sets that were read it
1154  Range read_sets;
1155  rval = mb.get_entities_by_type( file, MBENTITYSET, read_sets );CHECK_ERR( rval );
1156 
1157  // count number of vertices in each group
1158  int counts[NUM_GRP];
1159  memset( counts, 0, sizeof( counts ) );
1160  verts.clear();
1161  rval = mb.get_entities_by_type( 0, MBVERTEX, verts );CHECK_ERR( rval );
1162  for( Range::iterator it = verts.begin(); it != verts.end(); ++it )
1163  {
1164  double coords[3];
1165  rval = mb.get_coords( &*it, 1, coords );CHECK_ERR( rval );
1166  int i = (int)( coords[0] + 1e-12 );
1167  CHECK( i >= 0 && i < NUM_GRP );
1168  counts[i]++;
1169  }
1170 
1171  // check expected counts
1172  CHECK_EQUAL( 0, counts[0] );
1173  CHECK_EQUAL( 0, counts[1] );
1174  CHECK_EQUAL( NUM_GRP_VTX, counts[2] );
1175  switch( mode )
1176  {
1177  case GATHER_NONE:
1178  CHECK_EQUAL( 0, counts[3] );
1179  CHECK_EQUAL( 1 + contents, (int)read_sets.size() );
1180  break;
1181  case GATHER_SETS:
1182  CHECK_EQUAL( 0, counts[3] );
1183  CHECK_EQUAL( 2, (int)read_sets.size() );
1184  break;
1185  case GATHER_CONTENTS:
1186  CHECK_EQUAL( NUM_GRP_VTX, counts[3] );
1187  CHECK_EQUAL( 2, (int)read_sets.size() );
1188  break;
1189  }
1190 }
1191 
1192 static void check_num_verts( Interface& mb, Tag tag, int id, int num_vtx )
1193 {
1194  ErrorCode rval;
1195  const void* val[] = { &id };
1196  Range range;
1197  rval = mb.get_entities_by_type_and_tag( 0, MBENTITYSET, &tag, val, 1, range );CHECK_ERR( rval );
1198  CHECK_EQUAL( 1, (int)range.size() );
1199 
1200  EntityHandle set = range.front();
1201  range.clear();
1202  rval = mb.get_entities_by_type( set, MBVERTEX, range );CHECK_ERR( rval );
1203  CHECK_EQUAL( num_vtx, (int)range.size() );
1204 }
1205 
1206 //! Read in the sets contained in a set.
1207 //! Should read all sets containing read elements or nodes
1208 //! and all sets that are contained the the specified "read"
1209 //! set. Test the former here.
1210 void test_read_containing_sets()
1211 {
1212  // create mesh decomposed by elements but create
1213  // sets containing all vertices of decomposed elements
1214  // such that adjacent sets share vertices.
1215  create_mesh( false, true, false, false );
1216 
1217  ErrorCode rval;
1218  Core instance;
1219  Interface& mb = instance;
1220 
1221  // read some sets
1222  const int ids[] = { 1, 5, 9 };
1223  const int num_sets = sizeof( ids ) / sizeof( int );
1224  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, num_sets );CHECK_ERR( rval );
1225 
1226  Tag id_tag;
1227  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag );CHECK_ERR( rval );
1228 
1229  // expect all sets adjacent to the specified sets because
1230  // they share vertices.
1231  Range verts;
1232  for( int i = 0; i < num_sets; ++i )
1233  {
1234  if( ids[i] > 1 ) check_num_verts( mb, id_tag, ids[i] - 1, MBQUAD_INT + 1 );
1235  check_num_verts( mb, id_tag, ids[i], ( MBQUAD_INT + 1 ) * ( SET_WIDTH + 1 ) );
1236  if( ids[i] < NUM_SETS ) check_num_verts( mb, id_tag, ids[i] + 1, MBQUAD_INT + 1 );
1237  }
1238 }
1239 
1240 //! Test reading of explicit adjacencies
1241 void test_read_adjacencies()
1242 {
1243  ErrorCode rval;
1244  Core instance;
1245  Interface& mb = instance;
1246 
1247  // create four hexes sharing an edge
1248  EntityHandle verts[3][3][2], hexes[2][2];
1249  for( int k = 0; k < 2; ++k )
1250  {
1251  for( int j = 0; j < 3; ++j )
1252  {
1253  for( int i = 0; i < 3; ++i )
1254  {
1255  double coords[] = { static_cast< double >( i ), static_cast< double >( j ),
1256  static_cast< double >( k ) };
1257  rval = mb.create_vertex( coords, verts[i][j][k] );CHECK_ERR( rval );
1258  }
1259  }
1260  }
1261  for( int j = 0; j < 2; ++j )
1262  {
1263  for( int i = 0; i < 2; ++i )
1264  {
1265  EntityHandle conn[] = { verts[i][j][0], verts[i + 1][j][0], verts[i + 1][j + 1][0], verts[i][j + 1][0],
1266  verts[i][j][1], verts[i + 1][j][1], verts[i + 1][j + 1][1], verts[i][j + 1][1] };
1267  rval = mb.create_element( MBHEX, conn, 8, hexes[i][j] );CHECK_ERR( rval );
1268  }
1269  }
1270 
1271  // create two duplicate edges that connect the vertices common to all four hexes
1272  EntityHandle edge_conn[2] = { verts[1][1][0], verts[1][1][1] };
1273  EntityHandle edges[2];
1274  rval = mb.create_element( MBEDGE, edge_conn, 2, edges[0] );CHECK_ERR( rval );
1275  rval = mb.create_element( MBEDGE, edge_conn, 2, edges[1] );CHECK_ERR( rval );
1276  // mark one edge as adjacent to the left two hexes and the
1277  // other as adjacent to the right two
1278  rval = mb.add_adjacencies( edges[0], hexes[0], 2, true );CHECK_ERR( rval );
1279  rval = mb.add_adjacencies( edges[1], hexes[1], 2, true );CHECK_ERR( rval );
1280  // create two sets containing the front two and the rear two
1281  // hexes, respectively.
1282  EntityHandle sets[2];
1283  rval = mb.create_meshset( MESHSET_SET, sets[0] );CHECK_ERR( rval );
1284  rval = mb.create_meshset( MESHSET_SET, sets[1] );CHECK_ERR( rval );
1285  EntityHandle set1[4] = { hexes[0][0], hexes[1][0], edges[0], edges[1] };
1286  EntityHandle set2[4] = { hexes[0][1], hexes[1][1], edges[0], edges[1] };
1287  rval = mb.add_entities( sets[0], set1, 4 );CHECK_ERR( rval );
1288  rval = mb.add_entities( sets[1], set2, 4 );CHECK_ERR( rval );
1289 
1290  // assign IDs to sets
1291  Tag id_tag;
1292  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
1293  int ids[2] = { 1, 2 };
1294  rval = mb.tag_set_data( id_tag, sets, 2, ids );CHECK_ERR( rval );
1295 
1296  // write mesh
1297  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
1298 
1299  // read mesh
1300  rval = mb.delete_mesh();CHECK_ERR( rval );
1301  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 1 );CHECK_ERR( rval );
1302 
1303  // expect two hexes and two edges
1304  Range range;
1305  rval = mb.get_entities_by_type( 0, MBHEX, range );CHECK_ERR( rval );
1306  CHECK_EQUAL( 2, (int)range.size() );
1307  EntityHandle h1 = range.front(), h2 = range.back();
1308  range.clear();
1309  rval = mb.get_entities_by_type( 0, MBEDGE, range );CHECK_ERR( rval );
1310  CHECK_EQUAL( 2, (int)range.size() );
1311 
1312  // expecte each edge to have one of the hexes
1313  range.clear();
1314  rval = mb.get_adjacencies( &h1, 1, 1, false, range );CHECK_ERR( rval );
1315  CHECK_EQUAL( 1, (int)range.size() );
1316  EntityHandle e1 = range.front();
1317  range.clear();
1318  rval = mb.get_adjacencies( &h2, 1, 1, false, range );CHECK_ERR( rval );
1319  CHECK_EQUAL( 1, (int)range.size() );
1320  EntityHandle e2 = range.front();
1321 
1322  CHECK( e1 != e2 );
1323 }
1324 
1325 void test_read_sides()
1326 {
1327  ErrorCode rval;
1328  Core instance;
1329  Interface& mb = instance;
1330 
1331  // create 4x4 grid of quads with edges
1332  const int INT = 4;
1333  EntityHandle verts[INT + 1][INT + 1];
1334  for( int j = 0; j <= INT; ++j )
1335  {
1336  for( int i = 0; i <= INT; ++i )
1337  {
1338  double coords[3] = { static_cast< double >( i ), static_cast< double >( j ), 0 };
1339  rval = mb.create_vertex( coords, verts[INT - j][i] );CHECK_ERR( rval );
1340  }
1341  }
1342  EntityHandle quads[INT][INT];
1343  for( int j = 0; j < INT; ++j )
1344  {
1345  for( int i = 0; i < INT; ++i )
1346  {
1347  EntityHandle conn[4] = { verts[INT - j][i], verts[INT - j][i + 1], verts[INT - j - 1][i + 1],
1348  verts[INT - j - 1][i] };
1349  rval = mb.create_element( MBQUAD, conn, 4, quads[INT - j - 1][i] );CHECK_ERR( rval );
1350  }
1351  }
1352  Range edges;
1353  rval = mb.get_adjacencies( &quads[0][0], INT * INT, 1, true, edges, Interface::UNION );CHECK_ERR( rval );
1354  CHECK_EQUAL( 40, (int)edges.size() );
1355 
1356  // group quads into two sets
1357  EntityHandle sets[2];
1358  rval = mb.create_meshset( MESHSET_SET, sets[0] );CHECK_ERR( rval );
1359  rval = mb.create_meshset( MESHSET_SET, sets[1] );CHECK_ERR( rval );
1360  rval = mb.add_entities( sets[0], quads[0], INT );CHECK_ERR( rval );
1361  rval = mb.add_entities( sets[1], quads[1], INT );CHECK_ERR( rval );
1362  rval = mb.add_entities( sets[0], quads[2], INT );CHECK_ERR( rval );
1363  rval = mb.add_entities( sets[1], quads[3], INT );CHECK_ERR( rval );
1364 
1365  // assign IDS
1366  Tag id_tag;
1367  rval = mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id_tag, MB_TAG_SPARSE | MB_TAG_EXCL );CHECK_ERR( rval );
1368  int ids[2] = { 4, 5 };
1369  rval = mb.tag_set_data( id_tag, sets, 2, ids );CHECK_ERR( rval );
1370 
1371  // write mesh
1372  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
1373 
1374  // read first set back in
1375  rval = mb.delete_mesh();CHECK_ERR( rval );
1376  rval = mb.load_file( TEST_FILE, 0, READ_OPTS, ID_TAG_NAME, ids, 1 );CHECK_ERR( rval );
1377 
1378  // check expected counts
1379  int count;
1380  rval = mb.get_number_entities_by_type( 0, MBVERTEX, count );CHECK_ERR( rval );
1381  CHECK_EQUAL( ( INT + 1 ) * INT, count );
1382  rval = mb.get_number_entities_by_type( 0, MBQUAD, count );CHECK_ERR( rval );
1383  CHECK_EQUAL( INT * INT / 2, count );
1384  rval = mb.get_number_entities_by_type( 0, MBEDGE, count );CHECK_ERR( rval );
1385  CHECK_EQUAL( 2 * ( INT + 1 ) + INT * INT, count );
1386 
1387  // check edges adjacent to each quad
1388  Range elems;
1389  rval = mb.get_entities_by_type( 0, MBQUAD, elems );CHECK_ERR( rval );
1390  for( Range::iterator it = elems.begin(); it != elems.end(); ++it )
1391  {
1392  edges.clear();
1393  rval = mb.get_adjacencies( &*it, 1, 1, false, edges );CHECK_ERR( rval );
1394  CHECK_EQUAL( 4, (int)edges.size() );
1395  }
1396 }
1397 
1398 const int expected_ids[] = { 2, 4, 6, 8, 10, 12, 14, 16, 18 };
1399 const int expected_vols[] = { 3, 7, 10 };
1400 
1401 void write_id_test_file()
1402 {
1403  Core moab;
1404  Interface& mb = moab;
1405  ErrorCode rval;
1406 
1407  // create 12 entity sets
1408  EntityHandle sets[12];
1409  for( int i = 0; i < 12; ++i )
1410  {
1411  rval = mb.create_meshset( MESHSET_SET, sets[i] );CHECK_ERR( rval );
1412  }
1413 
1414  // create tag handles
1415  Tag id = 0, gid = 0, dim = 0;
1416  mb.tag_get_handle( ID_TAG_NAME, 1, MB_TYPE_INTEGER, id, MB_TAG_SPARSE | MB_TAG_EXCL );
1417  mb.tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, dim, MB_TAG_SPARSE | MB_TAG_EXCL );
1418  gid = mb.globalId_tag();
1419 
1420  // set ID tag on first 10 sets
1421  rval = mb.tag_set_data( id, sets, sizeof( expected_ids ) / sizeof( int ), expected_ids );CHECK_ERR( rval );
1422  // set geom dim on all sets, only three of them have dim == 3
1423  int num_vol = sizeof( expected_vols ) / sizeof( int );
1424  int dims[12], ids[12];
1425  int v = 0;
1426  for( int i = 0; i < 12; ++i )
1427  {
1428  dims[i] = i % 3 + 1;
1429  if( dims[i] == 3 )
1430  {
1431  if( v < num_vol )
1432  ids[i] = expected_vols[v++];
1433  else
1434  ids[i] = expected_vols[0];
1435  }
1436  else
1437  ids[i] = 100;
1438  }
1439  rval = mb.tag_set_data( gid, sets, 12, ids );CHECK_ERR( rval );
1440  rval = mb.tag_set_data( dim, sets, 12, dims );CHECK_ERR( rval );
1441 
1442  rval = mb.write_file( TEST_FILE, "MOAB" );CHECK_ERR( rval );
1443 }
1444 
1445 void test_read_ids()
1446 {
1447  write_id_test_file();
1448 
1449  Core moab;
1450  ReadHDF5 reader( &moab );
1451  FileOptions opts( "" );
1452  ErrorCode rval;
1453  std::vector< int > values;
1454  rval = reader.read_tag_values( TEST_FILE, ID_TAG_NAME, opts, values );
1455  remove( TEST_FILE );CHECK_ERR( rval );
1456 
1457  std::sort( values.begin(), values.end() );
1458  std::vector< int > expected( expected_ids, expected_ids + sizeof( expected_ids ) / sizeof( int ) );
1459  CHECK_EQUAL( expected, values );
1460 }
1461 
1462 void test_read_partial_ids()
1463 {
1464  write_id_test_file();
1465 
1466  const int three = 3;
1467  ReaderIface::IDTag vols = { GEOM_DIMENSION_TAG_NAME, &three, 1 };
1468  ReaderIface::SubsetList subset = { &vols, 1, 0, 0 };
1469 
1470  Core moab;
1471  ReadHDF5 reader( &moab );
1472  FileOptions opts( "" );
1473  ErrorCode rval;
1474  std::vector< int > values;
1475  rval = reader.read_tag_values( TEST_FILE, GLOBAL_ID_TAG_NAME, opts, values, &subset );
1476  remove( TEST_FILE );CHECK_ERR( rval );
1477 
1478  std::sort( values.begin(), values.end() );
1479  std::vector< int > expected( expected_ids, expected_ids + sizeof( expected_ids ) / sizeof( int ) );
1480 }