 |
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
|
Mesh Oriented datABase
(version 5.5.1)
An array-based unstructured mesh library
|
1 /**
2 * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3 * storing and accessing finite element mesh data.
4 *
5 * Copyright 2004 Sandia Corporation. Under the terms of Contract
6 * DE-AC04-94AL85000 with Sandia Coroporation, 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 /**
17 * \class ReadRTT
18 * \brief ReadRTT based on ReadNASTRAN
19 *
20 * See:
21 *
22 * \author Andrew Davis
23 */
24
25 #include "ReadRTT.hpp"
26
27 #include <iostream>
28 #include <sstream>
29 #include <fstream>
30 #include <vector>
31 #include <cstdlib>
32 #include <map>
33 #include <cassert>
34 #include <cmath>
35
36 #include "moab/Interface.hpp"
37 #include "moab/ReadUtilIface.hpp"
38 #include "Internals.hpp" // for MB_START_ID
39 #include "moab/Range.hpp"
40 #include "moab/FileOptions.hpp"
41 #include "FileTokenizer.hpp"
42 #include "MBTagConventions.hpp"
43 #include "moab/CN.hpp"
44 #include "moab/ErrorHandler.hpp"
45 #include "moab/GeomTopoTool.hpp"
46
47 namespace moab
48 {
49
50 ReaderIface* ReadRTT::factory( Interface* iface )
51 {
52 return new ReadRTT( iface );
53 }
54
55 // constructor
56 ReadRTT::ReadRTT( Interface* impl )
57 : MBI( impl ), geom_tag( 0 ), id_tag( 0 ), name_tag( 0 ), category_tag( 0 ), faceting_tol_tag( 0 )
58 {
59 assert( NULL != impl );
60 myGeomTool = new GeomTopoTool( impl );
61 MBI->query_interface( readMeshIface );
62 assert( NULL != readMeshIface );
63
64 // this section copied from ReadCGM initalisation
65 int negone = -1;
66 double zero = 0.;
67 ErrorCode rval;
68 rval = MBI->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, geom_tag, MB_TAG_SPARSE | MB_TAG_CREAT,
69 &negone );MB_CHK_ERR_CONT(rval);
70 id_tag = MBI->globalId_tag();
71 rval = MBI->tag_get_handle( NAME_TAG_NAME, NAME_TAG_SIZE, MB_TYPE_OPAQUE, name_tag, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_ERR_CONT(rval);
72 rval = MBI->tag_get_handle( CATEGORY_TAG_NAME, CATEGORY_TAG_SIZE, MB_TYPE_OPAQUE, category_tag,
73 MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_ERR_CONT(rval);
74 rval =
75 MBI->tag_get_handle( "FACETING_TOL", 1, MB_TYPE_DOUBLE, faceting_tol_tag, MB_TAG_SPARSE | MB_TAG_CREAT, &zero );MB_CHK_ERR_CONT(rval);
76 }
77
78 // destructor
79 ReadRTT::~ReadRTT()
80 {
81 if( readMeshIface )
82 {
83 MBI->release_interface( readMeshIface );
84 readMeshIface = 0;
85 }
86
87 delete myGeomTool;
88 }
89
90 ErrorCode ReadRTT::read_tag_values( const char* /*file_name*/,
91 const char* /*tag_name*/,
92 const FileOptions& /*opts*/,
93 std::vector< int >& /*tag_values_out*/,
94 const SubsetList* /*subset_list*/ )
95 {
96 return MB_NOT_IMPLEMENTED;
97 }
98
99 // load the file as called by the Interface function
100 ErrorCode ReadRTT::load_file( const char* filename,
101 const EntityHandle*,
102 const FileOptions&,
103 const ReaderIface::SubsetList* subset_list,
104 const Tag* /*file_id_tag*/ )
105 {
106 ErrorCode rval;
107
108 // at this time there is no support for reading a subset of the file
109 if( subset_list )
110 {
111 std::cout << "Subset reading not supported for RTT meshes" << std::endl;
112 return MB_UNSUPPORTED_OPERATION;
113 }
114
115 // test to see if file exists
116 FILE* file = NULL;
117 file = fopen( filename, "r" );
118 if( file == NULL ) return MB_FILE_DOES_NOT_EXIST;
119 // otherwise close the file
120 fclose( file );
121
122 // read the header
123 rval = ReadRTT::read_header( filename );
124 if( rval != MB_SUCCESS ) return rval;
125
126 // read the side_flag data
127 std::vector< side > side_data;
128 rval = ReadRTT::read_sides( filename, side_data );
129 if( rval != MB_SUCCESS ) return rval;
130
131 // read the cell data
132 std::vector< cell > cell_data;
133 rval = ReadRTT::read_cells( filename, cell_data );
134 if( rval != MB_SUCCESS ) return rval;
135
136 // read the node data
137 std::vector< node > node_data;
138 rval = ReadRTT::read_nodes( filename, node_data );
139 if( rval != MB_SUCCESS ) return rval;
140
141 // read the facet data
142 std::vector< facet > facet_data;
143 rval = ReadRTT::read_facets( filename, facet_data );
144 if( rval != MB_SUCCESS ) return rval;
145
146 // read the tetrahedra data
147 std::vector< tet > tet_data;
148 rval = ReadRTT::read_tets( filename, tet_data );
149 if( rval != MB_SUCCESS ) return rval;
150
151 // make the map of surface number in the rttmesh to the surface meshset
152 std::map< int, EntityHandle > surface_map; // corrsespondance of surface number to entity handle
153 rval = ReadRTT::generate_topology( side_data, cell_data, surface_map );
154 if( rval != MB_SUCCESS ) return rval;
155
156 // generate the rest of the database, triangles to surface meshsets etc
157 rval = ReadRTT::build_moab( node_data, facet_data, tet_data, surface_map );
158 if( rval != MB_SUCCESS ) return rval;
159
160 return MB_SUCCESS;
161 }
162
163 /*
164 * builds the topology of the problem
165 */
166 ErrorCode ReadRTT::generate_topology( std::vector< side > side_data,
167 std::vector< cell > cell_data,
168 std::map< int, EntityHandle >& surface_map )
169 {
170
171 ErrorCode rval;
172 std::vector< EntityHandle > entmap[4];
173 int num_ents[4]; // number of entities in each dimension
174
175 const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
176
177 std::vector< int > surface_numbers; // the surface numbers in the problem
178
179 // corresponds to number of cad like surfaces and cad like volumes
180 num_ents[2] = side_data.size();
181 num_ents[3] = cell_data.size();
182
183 // loop over surfaces & volumes
184 for( int dim = 2; dim <= 3; dim++ )
185 {
186 for( int i = 0; i != num_ents[dim]; i++ )
187 {
188 EntityHandle handle;
189 // create a meshset for each entity surface/volume
190 rval = MBI->create_meshset( dim == 1 ? MESHSET_ORDERED : MESHSET_SET, handle );
191 // if failure
192 if( rval != MB_SUCCESS ) return rval;
193
194 // collect the entity handles into an
195 entmap[dim].push_back( handle );
196
197 // set the dimension tag
198 rval = MBI->tag_set_data( geom_tag, &handle, 1, &dim );
199 // if fail
200 if( MB_SUCCESS != rval ) return rval;
201 // if we are a surface
202 if( dim == 2 )
203 {
204 // tag the id onto the surface meshset
205 rval = MBI->tag_set_data( id_tag, &handle, 1, &side_data[i].id );
206 // inesert entity into the map
207 surface_map[side_data[i].id] = handle;
208 }
209 else
210 {
211 // otherwise we set the volume tag data, loop is only 2 & 3 dim
212 rval = MBI->tag_set_data( id_tag, &handle, 1, &cell_data[i].id );
213 }
214 // if fail
215 if( MB_SUCCESS != rval ) return rval;
216 // set the category tag
217 rval = MBI->tag_set_data( category_tag, &handle, 1, &geom_categories[dim] );
218 if( MB_SUCCESS != rval ) return rval;
219 }
220 }
221
222 // generate parent child links
223 // best to loop over the surfaces and assign them to volumes, we can then assign facets to
224 // to each surface
225 generate_parent_child_links( num_ents, entmap, side_data, cell_data );
226
227 // set the surface senses
228 set_surface_senses( num_ents, entmap, side_data, cell_data );
229
230 // set the group data
231 rval = setup_group_data( entmap );
232
233 return MB_SUCCESS;
234 }
235
236 /*
237 * builds the moab representation of the mesh
238 */
239 ErrorCode ReadRTT::build_moab( std::vector< node > node_data,
240 std::vector< facet > facet_data,
241 std::vector< tet > tet_data,
242 std::map< int, EntityHandle > surface_map )
243 {
244
245 ErrorCode rval; // reusable return value
246 EntityHandle file_set; // the file handle
247 // create the file set
248 rval = MBI->create_meshset( MESHSET_SET, file_set );
249 if( MB_SUCCESS != rval ) return rval;
250
251 // create the vertices
252 EntityHandle handle;
253 std::vector< node >::iterator it; // iterate over the nodes
254 Range mb_coords; // range of coordinates
255 for( it = node_data.begin(); it != node_data.end(); ++it )
256 {
257 node tmp = *it;
258 double coords[3] = { tmp.x, tmp.y, tmp.z };
259 rval = MBI->create_vertex( coords, handle );
260 if( MB_SUCCESS != rval ) return rval;
261 mb_coords.insert( handle ); // inesert handle into the coordinate range
262 }
263
264 // add verts to set
265 rval = MBI->add_entities( file_set, mb_coords );
266
267 // create sense tag
268 Tag side_id_tag, surface_number_tag;
269 // int zero = 0;
270 rval = MBI->tag_get_handle( "SIDEID_TAG", 1, MB_TYPE_INTEGER, side_id_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
271 rval =
272 MBI->tag_get_handle( "SURFACE_NUMBER", 1, MB_TYPE_INTEGER, surface_number_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
273
274 // create the facets
275 EntityHandle triangle;
276 std::vector< facet >::iterator it_f;
277 // range of triangles
278 Range mb_tris;
279 // loop over the facet data
280 for( it_f = facet_data.begin(); it_f != facet_data.end(); ++it_f )
281 {
282 facet tmp = *it_f;
283 // get the nodes for the triangle
284 EntityHandle tri_nodes[3] = { mb_coords[tmp.connectivity[0] - 1], mb_coords[tmp.connectivity[1] - 1],
285 mb_coords[tmp.connectivity[2] - 1] };
286 // create a triangle element
287 rval = MBI->create_element( MBTRI, tri_nodes, 3, triangle );
288 // tag in side id on the triangle
289 rval = MBI->tag_set_data( side_id_tag, &triangle, 1, &tmp.side_id );
290 // tag the surface number on the triangle
291 rval = MBI->tag_set_data( surface_number_tag, &triangle, 1, &tmp.surface_number );
292 // insert vertices and triangles into the appropriate surface meshset
293 EntityHandle meshset_handle = surface_map[tmp.surface_number];
294 // also set surface tag
295 rval = MBI->tag_set_data( side_id_tag, &meshset_handle, 1, &tmp.side_id );
296 rval = MBI->tag_set_data( surface_number_tag, &meshset_handle, 1, &tmp.surface_number );
297 // add vertices to the mesh
298 rval = MBI->add_entities( meshset_handle, &( *tri_nodes ), 3 );
299 // add triangles to the meshset
300 rval = MBI->add_entities( meshset_handle, &triangle, 1 );
301 // ineter triangles into large run
302 mb_tris.insert( triangle );
303 }
304 // add tris to set to fileset
305 rval = MBI->add_entities( file_set, mb_tris );
306
307 // create material number tag
308 Tag mat_num_tag;
309 // int zero = 0;
310 rval = MBI->tag_get_handle( "MATERIAL_NUMBER", 1, MB_TYPE_INTEGER, mat_num_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
311
312 // create the tets
313 EntityHandle tetra; // handle for a specific tet
314 std::vector< tet >::iterator it_t;
315 Range mb_tets;
316 // loop over all tets
317 for( it_t = tet_data.begin(); it_t != tet_data.end(); ++it_t )
318 {
319 tet tmp = *it_t;
320 // get the handles for the tet
321 EntityHandle tet_nodes[4] = { mb_coords[tmp.connectivity[0] - 1], mb_coords[tmp.connectivity[1] - 1],
322 mb_coords[tmp.connectivity[2] - 1], mb_coords[tmp.connectivity[3] - 1] };
323 // create the tet
324 rval = MBI->create_element( MBTET, tet_nodes, 4, tetra );
325 int mat_number = tmp.material_number;
326 // tag the tet with the material number
327 rval = MBI->tag_set_data( mat_num_tag, &tetra, 1, &mat_number );
328 // set the tag data
329 mb_tets.insert( tetra );
330 }
331 // add tris to set
332 rval = MBI->add_entities( file_set, mb_tets );
333
334 return MB_SUCCESS;
335 }
336
337 /*
338 * read the header data from the filename pointed to
339 */
340 ErrorCode ReadRTT::read_header( const char* filename )
341 {
342 std::ifstream input_file( filename ); // filename for rtt file
343 // file ok?
344 if( !input_file.good() )
345 {
346 std::cout << "Problems reading file = " << filename << std::endl;
347 return MB_FAILURE;
348 }
349
350 // if it works
351 std::string line;
352 moab::ErrorCode rval = MB_FAILURE;
353 if( input_file.is_open() )
354 {
355 while( std::getline( input_file, line ) )
356 {
357 if( line.compare( "header" ) == 0 )
358 {
359 rval = get_header_data( input_file );
360 }
361 }
362 input_file.close();
363 }
364 return rval;
365 }
366
367 /*
368 * reads the side data from the filename pointed to
369 */
370 ErrorCode ReadRTT::read_sides( const char* filename, std::vector< side >& side_data )
371 {
372 std::string line; // the current line being read
373 std::ifstream input_file( filename ); // filestream for rttfile
374 // file ok?
375 if( !input_file.good() )
376 {
377 std::cout << "Problems reading file = " << filename << std::endl;
378 return MB_FAILURE;
379 }
380 // if it works
381 if( input_file.is_open() )
382 {
383 while( std::getline( input_file, line ) )
384 {
385 if( line.compare( " 2 FACES\0" ) == 0 )
386 {
387 // read lines until find end nodes
388 while( std::getline( input_file, line ) )
389 {
390 if( line.compare( "end_side_flags\0" ) == 0 ) break;
391 side data = ReadRTT::get_side_data( line );
392 side_data.push_back( data );
393 }
394 }
395 }
396 input_file.close();
397 }
398 if( side_data.size() == 0 ) return MB_FAILURE;
399 return MB_SUCCESS;
400 }
401
402 /*
403 * reads the cell data from the filename pointed to
404 */
405 ErrorCode ReadRTT::read_cells( const char* filename, std::vector< cell >& cell_data )
406 {
407 std::string line; // the current line being read
408 std::ifstream input_file( filename ); // filestream for rttfile
409 // file ok?
410 if( !input_file.good() )
411 {
412 std::cout << "Problems reading file = " << filename << std::endl;
413 return MB_FAILURE;
414 }
415 // if it works
416 if( input_file.is_open() )
417 {
418 while( std::getline( input_file, line ) )
419 {
420 if( line.compare( " 1 REGIONS\0" ) == 0 )
421 {
422 // read lines until find end nodes
423 while( std::getline( input_file, line ) )
424 {
425 if( line.compare( "end_cell_flags\0" ) == 0 ) break;
426 cell data = ReadRTT::get_cell_data( line );
427 cell_data.push_back( data );
428 }
429 }
430 }
431 input_file.close();
432 }
433 if( cell_data.size() == 0 ) return MB_FAILURE;
434 return MB_SUCCESS;
435 }
436
437 /*
438 * Reads the node data fromt the filename pointed to
439 */
440 ErrorCode ReadRTT::read_nodes( const char* filename, std::vector< node >& node_data )
441 {
442 std::string line; // the current line being read
443 std::ifstream input_file( filename ); // filestream for rttfile
444 // file ok?
445 if( !input_file.good() )
446 {
447 std::cout << "Problems reading file = " << filename << std::endl;
448 return MB_FAILURE;
449 }
450
451 // if it works
452 if( input_file.is_open() )
453 {
454 while( std::getline( input_file, line ) )
455 {
456 if( line.compare( "nodes\0" ) == 0 )
457 {
458 // read lines until find end nodes
459 while( std::getline( input_file, line ) )
460 {
461 if( line.compare( "end_nodes\0" ) == 0 ) break;
462 node data = ReadRTT::get_node_data( line );
463 node_data.push_back( data );
464 }
465 }
466 }
467 input_file.close();
468 }
469 if( node_data.size() == 0 ) return MB_FAILURE;
470 return MB_SUCCESS;
471 }
472
473 /*
474 * Reads the facet data fromt the filename pointed to
475 */
476 ErrorCode ReadRTT::read_facets( const char* filename, std::vector< facet >& facet_data )
477 {
478 std::string line; // the current line being read
479 std::ifstream input_file( filename ); // filestream for rttfile
480 // file ok?
481 if( !input_file.good() )
482 {
483 std::cout << "Problems reading file = " << filename << std::endl;
484 return MB_FAILURE;
485 }
486
487 // if it works
488 if( input_file.is_open() )
489 {
490 while( std::getline( input_file, line ) )
491 {
492 if( line.compare( "sides\0" ) == 0 )
493 {
494 // read lines until find end nodes
495 while( std::getline( input_file, line ) )
496 {
497 if( line.compare( "end_sides\0" ) == 0 ) break;
498 facet data = ReadRTT::get_facet_data( line );
499 facet_data.push_back( data );
500 }
501 }
502 }
503 input_file.close();
504 }
505 if( facet_data.size() == 0 ) return MB_FAILURE;
506 return MB_SUCCESS;
507 }
508
509 /*
510 * Reads the facet data fromt the filename pointed to
511 */
512 ErrorCode ReadRTT::read_tets( const char* filename, std::vector< tet >& tet_data )
513 {
514 std::string line; // the current line being read
515 std::ifstream input_file( filename ); // filestream for rttfile
516 // file ok?
517 if( !input_file.good() )
518 {
519 std::cout << "Problems reading file = " << filename << std::endl;
520 return MB_FAILURE;
521 }
522 // if it works
523 if( input_file.is_open() )
524 {
525 while( std::getline( input_file, line ) )
526 {
527 if( line.compare( "cells\0" ) == 0 )
528 {
529 // read lines until find end nodes
530 while( std::getline( input_file, line ) )
531 {
532 if( line.compare( "end_cells\0" ) == 0 ) break;
533 tet data = ReadRTT::get_tet_data( line );
534 tet_data.push_back( data );
535 }
536 }
537 }
538 input_file.close();
539 }
540 if( tet_data.size() == 0 ) return MB_FAILURE;
541 return MB_SUCCESS;
542 }
543
544 /*
545 * given the open file handle read until we find
546 */
547 ErrorCode ReadRTT::get_header_data( std::ifstream& input_file )
548 {
549 std::string line;
550 while( std::getline( input_file, line ) )
551 {
552
553 // tokenize the line
554 std::istringstream iss( line );
555 std::vector< std::string > split_string;
556 do
557 {
558 std::string sub_string;
559 iss >> sub_string;
560 split_string.push_back( sub_string );
561 } while( iss );
562
563 // if we find version
564 if( line.find( "version" ) != std::string::npos )
565 {
566 if( split_string[1].find( "v" ) != std::string::npos &&
567 split_string[0].find( "version" ) != std::string::npos )
568 {
569 header_data.version = split_string[1];
570 }
571 }
572
573 if( line.find( "title" ) != std::string::npos )
574 {
575 header_data.title = split_string[1];
576 }
577 if( line.find( "date" ) != std::string::npos )
578 {
579 header_data.date = split_string[1];
580 }
581 if( line.find( "end_header" ) != std::string::npos )
582 {
583 return MB_SUCCESS;
584 }
585 }
586
587 // otherwise we never found the end_header keyword
588 return MB_FAILURE;
589 }
590
591 /*
592 * given the string sidedata, get the id number, senses and names of the sides
593 */
594 ReadRTT::side ReadRTT::get_side_data( std::string sidedata )
595 {
596 side new_side;
597 std::vector< std::string > tokens;
598 tokens = ReadRTT::split_string( sidedata, ' ' );
599
600 // set the side id
601 if( tokens.size() != 2 )
602 {
603 MB_SET_ERR_RET_VAL( "Error, too many tokens found from side_data", new_side );
604 }
605 // create the new side
606 new_side.id = std::atoi( tokens[0].c_str() );
607
608 std::vector< std::string > cell_names = ReadRTT::split_string( tokens[1], '/' );
609 // get the boundary
610 boundary new_bnd = ReadRTT::split_name( cell_names[0] );
611 // set the surface sense and name
612 new_side.senses[0] = new_bnd.sense;
613 new_side.names[0] = new_bnd.name;
614 //
615 if( cell_names.size() > 1 )
616 {
617 boundary bnd = ReadRTT::split_name( cell_names[1] );
618 new_side.senses[1] = bnd.sense;
619 new_side.names[1] = bnd.name;
620 }
621 else
622 {
623 new_side.senses[1] = 0;
624 new_side.names[1] = "\0";
625 }
626
627 return new_side;
628 }
629
630 /*
631 * given the string celldata, get the id number and name of each cell
632 */
633 ReadRTT::cell ReadRTT::get_cell_data( std::string celldata )
634 {
635 cell new_cell;
636 std::vector< std::string > tokens;
637 tokens = ReadRTT::split_string( celldata, ' ' );
638
639 // set the side id
640 if( tokens.size() != 2 )
641 {
642 MB_SET_ERR_RET_VAL( "Error, too many tokens found from cell_data", new_cell );
643 }
644 // create the new side
645 new_cell.id = std::atoi( tokens[0].c_str() );
646 new_cell.name = tokens[1];
647
648 return new_cell;
649 }
650
651 /*
652 * given the string nodedata, get the id number and coordinates of the node
653 */
654 ReadRTT::node ReadRTT::get_node_data( std::string nodedata )
655 {
656 node new_node;
657 std::vector< std::string > tokens;
658 tokens = ReadRTT::split_string( nodedata, ' ' );
659
660 // set the side id
661 if( tokens.size() != 5 )
662 {
663 MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_node_data", new_node );
664 }
665 new_node.id = std::atoi( tokens[0].c_str() );
666 new_node.x = std::atof( tokens[1].c_str() );
667 new_node.y = std::atof( tokens[2].c_str() );
668 new_node.z = std::atof( tokens[3].c_str() );
669 return new_node;
670 }
671
672 /*
673 * given the string nodedata, get the id number, connectivity and sense data
674 */
675 ReadRTT::facet ReadRTT::get_facet_data( std::string facetdata )
676 {
677 facet new_facet;
678 std::vector< std::string > tokens;
679 tokens = ReadRTT::split_string( facetdata, ' ' );
680
681 // set the side id
682 if( tokens.size() != 7 )
683 {
684 MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_facet_data", new_facet );
685 }
686
687 new_facet.id = std::atoi( tokens[0].c_str() );
688 // branch on the rtt version number
689 if( header_data.version == "v1.0.0" )
690 {
691 new_facet.connectivity[0] = std::atoi( tokens[1].c_str() );
692 new_facet.connectivity[1] = std::atoi( tokens[2].c_str() );
693 new_facet.connectivity[2] = std::atoi( tokens[3].c_str() );
694 new_facet.side_id = std::atoi( tokens[4].c_str() );
695 new_facet.surface_number = std::atoi( tokens[5].c_str() );
696 }
697 else if( header_data.version == "v1.0.1" )
698 {
699 new_facet.connectivity[0] = std::atoi( tokens[2].c_str() );
700 new_facet.connectivity[1] = std::atoi( tokens[3].c_str() );
701 new_facet.connectivity[2] = std::atoi( tokens[4].c_str() );
702 new_facet.side_id = std::atoi( tokens[5].c_str() );
703 new_facet.surface_number = std::atoi( tokens[6].c_str() );
704 }
705 else
706 {
707 MB_SET_ERR_RET_VAL( "Error, version number not understood", new_facet );
708 }
709
710 return new_facet;
711 }
712
713 /*
714 * given the string tetdata, get the id number, connectivity and mat num of the tet
715 */
716 ReadRTT::tet ReadRTT::get_tet_data( std::string tetdata )
717 {
718 tet new_tet;
719 std::vector< std::string > tokens;
720 tokens = ReadRTT::split_string( tetdata, ' ' );
721
722 // set the side id
723 if( tokens.size() != 7 )
724 {
725 MB_SET_ERR_RET_VAL( "Error, too many tokens found from get_tet_data", new_tet );
726 }
727 new_tet.id = std::atoi( tokens[0].c_str() );
728 // branch on the version number
729 if( header_data.version == "v1.0.0" )
730 {
731 new_tet.connectivity[0] = std::atoi( tokens[1].c_str() );
732 new_tet.connectivity[1] = std::atoi( tokens[2].c_str() );
733 new_tet.connectivity[2] = std::atoi( tokens[3].c_str() );
734 new_tet.connectivity[3] = std::atoi( tokens[4].c_str() );
735 new_tet.material_number = std::atoi( tokens[5].c_str() );
736 }
737 else if( header_data.version == "v1.0.1" )
738 {
739 new_tet.connectivity[0] = std::atoi( tokens[2].c_str() );
740 new_tet.connectivity[1] = std::atoi( tokens[3].c_str() );
741 new_tet.connectivity[2] = std::atoi( tokens[4].c_str() );
742 new_tet.connectivity[3] = std::atoi( tokens[5].c_str() );
743 new_tet.material_number = std::atoi( tokens[6].c_str() );
744 }
745 else
746 {
747 MB_SET_ERR_RET_VAL( "Error, version number not supported", new_tet );
748 }
749
750 return new_tet;
751 }
752
753 /*
754 * splits string into sense and name, to later facilitate the building
755 * of sense data, strips off the tailing @ if it exists
756 */
757 ReadRTT::boundary ReadRTT::split_name( std::string atilla_cellname )
758 {
759 boundary new_boundary;
760 // default initialisation
761 new_boundary.sense = 0;
762 new_boundary.name = "\0";
763 // +ve sense
764 if( atilla_cellname.find( "+" ) != std::string::npos )
765 {
766 new_boundary.sense = 1;
767 // look for the @# we do not want it
768 std::size_t found = atilla_cellname.find( "@" );
769 if( found != std::string::npos )
770 new_boundary.name = atilla_cellname.substr( 3, found );
771 else
772 new_boundary.name = atilla_cellname.substr( 3, atilla_cellname.length() );
773 }
774 else if( atilla_cellname.find( "-" ) != std::string::npos )
775 {
776 // negative sense
777 new_boundary.sense = -1;
778 new_boundary.name = atilla_cellname.substr( 3, atilla_cellname.length() );
779 }
780 return new_boundary;
781 }
782
783 /*
784 * splits a string in a vector of strings split by spaces
785 */
786 std::vector< std::string > ReadRTT::split_string( std::string string_to_split, char split_char )
787 {
788 std::istringstream ss( string_to_split );
789 std::vector< std::string > tokens;
790 while( !ss.eof() )
791 {
792 std::string x; // here's a nice, empty string
793 std::getline( ss, x, split_char ); // try to read the next field into it
794 tokens.push_back( x );
795 }
796
797 // remove empty tokens
798 std::vector< std::string >::iterator it;
799 for( it = tokens.begin(); it != tokens.end(); )
800 {
801 std::string string = *it;
802 if( string.compare( "\0" ) == 0 )
803 it = tokens.erase( it );
804 else
805 ++it;
806 }
807 return tokens;
808 }
809
810 /*
811 * Generate the parent-child links bwetween the cell and surface meshsets
812 */
813 void ReadRTT::generate_parent_child_links( int num_ents[4],
814 std::vector< EntityHandle > entity_map[4],
815 std::vector< side > side_data,
816 std::vector< cell > cell_data )
817 {
818 ErrorCode rval; // return value
819 // loop over the number of surfaces
820 for( int i = 0; i < num_ents[2]; i++ )
821 {
822 // get the surface handle
823 EntityHandle surf_handle = entity_map[2][i];
824 // there are volumes that share this face
825 for( unsigned int shared = 0; shared <= 1; shared++ )
826 {
827 std::string parent_name = side_data[i].names[shared];
828 // find the @ sign
829 unsigned pos = parent_name.find( "@" );
830 parent_name = parent_name.substr( 0, pos );
831
832 // loop over tets looking for matching name
833 for( int j = 0; j < num_ents[3]; j++ )
834 {
835 // if match found
836 if( cell_data[j].name.compare( parent_name ) == 0 )
837 {
838 EntityHandle cell_handle = entity_map[3][j];
839 // parent
840 rval = MBI->add_parent_child( cell_handle, surf_handle );
841 if( rval != MB_SUCCESS )
842 {
843 std::cerr << "Failed to add parent child relationship" << std::endl;
844 }
845 }
846 }
847 }
848 }
849 return;
850 }
851
852 /*
853 * sets the sense of the surfaces wrt to volumes using geom topo tool
854 */
855 void ReadRTT::set_surface_senses( int num_ents[4],
856 std::vector< EntityHandle > entity_map[4],
857 std::vector< side > side_data,
858 std::vector< cell > cell_data )
859 {
860
861 ErrorCode rval; // return value
862 // loop over the number of surfaces
863 for( int i = 0; i < num_ents[2]; i++ )
864 {
865 EntityHandle surf_handle = entity_map[2][i];
866 // there are 2 volumes that share this face
867 for( unsigned int shared = 0; shared <= 1; shared++ )
868 {
869 std::string parent_name = side_data[i].names[shared];
870 unsigned pos = parent_name.find( "@" );
871 parent_name = parent_name.substr( 0, pos );
872 // loop over tets looking for matching name
873 for( int j = 0; j < num_ents[3]; j++ )
874 {
875 // if match found
876 if( cell_data[j].name.compare( parent_name ) == 0 )
877 {
878 EntityHandle cell_handle = entity_map[3][j];
879 // in rtt mesh +represents the inside and -represents outside
880 // in moab reverse is outside and forward is inside
881 if( side_data[i].senses[shared] == 1 )
882 rval = myGeomTool->set_sense( surf_handle, cell_handle, SENSE_FORWARD );
883 else if( side_data[i].senses[shared] == -1 )
884 rval = myGeomTool->set_sense( surf_handle, cell_handle, SENSE_REVERSE );
885 else
886 rval = myGeomTool->set_sense( surf_handle, 0, SENSE_REVERSE );
887
888 if( rval != MB_SUCCESS )
889 {
890 std::cerr << "Failed to set sense appropriately" << std::endl;
891 }
892 }
893 }
894 }
895 }
896 return;
897 }
898
899 /*
900 * Add all entities that are to be part of the graveyard
901 */
902 ErrorCode ReadRTT::setup_group_data( std::vector< EntityHandle > entity_map[4] )
903 {
904 ErrorCode rval; // error codes
905 EntityHandle handle;
906 handle = create_group( "graveyard_comp", 1 );
907
908 // add any volume to group graveyard, it is ignored by dag
909 EntityHandle vol_handle = entity_map[3][0];
910 rval = MBI->add_entities( handle, &vol_handle, 1 );
911 return rval;
912 }
913
914 /*
915 * create a new group of with the name group name and id
916 */
917 EntityHandle ReadRTT::create_group( std::string group_name, int id )
918 {
919 ErrorCode rval;
920 // category tags
921 const char geom_categories[][CATEGORY_TAG_SIZE] = { "Vertex\0", "Curve\0", "Surface\0", "Volume\0", "Group\0" };
922
923 EntityHandle handle;
924 rval = MBI->create_meshset( MESHSET_SET, handle );
925 if( MB_SUCCESS != rval ) return rval;
926
927 rval = MBI->tag_set_data( name_tag, &handle, 1, group_name.c_str() );
928 if( MB_SUCCESS != rval ) return MB_FAILURE;
929
930 rval = MBI->tag_set_data( id_tag, &handle, 1, &id );
931 if( MB_SUCCESS != rval ) return MB_FAILURE;
932
933 rval = MBI->tag_set_data( category_tag, &handle, 1, &geom_categories[4] );
934 if( MB_SUCCESS != rval ) return MB_FAILURE;
935
936 return handle;
937 }
938
939 } // namespace moab
1.9.1.