 |
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 #include <cstdlib> // For exit()
2 #include <vector>
3 #include <map>
4 #include <iostream>
5 #include <string>
6 #include <algorithm>
7
8 #include "moab/CN.hpp"
9 #include "moab/Range.hpp"
10 #include "moab/Interface.hpp"
11 #include "MBTagConventions.hpp"
12 #include "Internals.hpp"
13 #include "moab/ReadUtilIface.hpp"
14 #include "moab/FileOptions.hpp"
15 #include "ReadCCMIO.hpp"
16 #include "moab/MeshTopoUtil.hpp"
17
18 #include "ccmio.h"
19
20 /*
21 * CCMIO file structure
22 *
23 * Root
24 * State(kCCMIOState)
25 * Processor*
26 * Vertices
27 * ->ReadVerticesx, ReadMap
28 * Topology
29 * Boundary faces*(kCCMIOBoundaryFaces)
30 * ->ReadFaces, ReadFaceCells, ReadMap
31 * Internal faces(kCCMIOInternalFaces)
32 * Cells (kCCMIOCells)
33 * ->ReadCells (mapID), ReadMap, ReadCells (cellTypes)
34 * Solution
35 * Phase
36 * Field
37 * FieldData
38 * Problem(kCCMIOProblemDescription)
39 * CellType* (kCCMIOCellType)
40 * Index (GetEntityIndex), MaterialId(ReadOpti), MaterialType(ReadOptstr),
41 * PorosityId(ReadOpti), SpinId(ReadOpti), GroupId(ReadOpti)
42 *
43 * MaterialType (CCMIOReadOptstr in readexample)
44 * constants (see readexample)
45 * lagrangian data (CCMIOReadLagrangianData)
46 * vertices label (CCMIOEntityDescription)
47 * restart info: char solver[], iteratoins, time, char timeUnits[], angle
48 * (CCMIOReadRestartInfo, kCCMIORestartData), reference data?
49 * phase:
50 * field: char name[], dims, CCMIODataType datatype, char units[]
51 * dims = kCCMIOScalar (CCMIOReadFieldDataf),
52 * kCCMIOVector (CCMIOReadMultiDimensionalFieldData),
53 * kCCMIOTensor
54 * MonitoringSets: num, name (CellSet, VertexSet, BoundarySet, BlockSet, SplineSet, CoupleSet)
55 * CCMIOGetProstarSet, CCMIOReadOpt1i,
56 */
57
58 enum DataType
59 {
60 kScalar,
61 kVector,
62 kVertex,
63 kCell,
64 kInternalFace,
65 kBoundaryFace,
66 kBoundaryData,
67 kBoundaryFaceData,
68 kCellType
69 };
70
71 namespace moab
72 {
73
74 static int const kNValues = 10; // Number of values of each element to print
75 static char const kDefaultState[] = "default";
76 static char const kUnitsName[] = "Units";
77 static int const kVertOffset = 2;
78 static int const kCellInc = 4;
79
80 #define CHK_SET_CCMERR( ccm_err_code, ccm_err_msg ) \
81 { \
82 if( kCCMIONoErr != ( ccm_err_code ) && kCCMIONoFileErr != ( ccm_err_code ) && \
83 kCCMIONoNodeErr != ( ccm_err_code ) ) \
84 MB_SET_ERR( MB_FAILURE, ccm_err_msg ); \
85 }
86
87 ReaderIface* ReadCCMIO::factory( Interface* iface )
88 {
89 return new ReadCCMIO( iface );
90 }
91
92 ReadCCMIO::ReadCCMIO( Interface* impl )
93 : mMaterialIdTag( 0 ), mMaterialTypeTag( 0 ), mRadiationTag( 0 ), mPorosityIdTag( 0 ), mSpinIdTag( 0 ),
94 mGroupIdTag( 0 ), mColorIdxTag( 0 ), mProcessorIdTag( 0 ), mLightMaterialTag( 0 ), mFreeSurfaceMaterialTag( 0 ),
95 mThicknessTag( 0 ), mProstarRegionNumberTag( 0 ), mBoundaryTypeTag( 0 ), mCreatingProgramTag( 0 ), mbImpl( impl ),
96 hasSolution( false )
97 {
98 assert( impl != NULL );
99
100 impl->query_interface( readMeshIface );
101
102 // Initialize in case tag_get_handle fails below
103 mMaterialSetTag = 0;
104 mDirichletSetTag = 0;
105 mNeumannSetTag = 0;
106 mHasMidNodesTag = 0;
107 mGlobalIdTag = 0;
108 mNameTag = 0;
109
110 //! Get and cache predefined tag handles
111 const int negone = -1;
112 ErrorCode result = impl->tag_get_handle( MATERIAL_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mMaterialSetTag,
113 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get MATERIAL_SET tag" );
114
115 result = impl->tag_get_handle( DIRICHLET_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mDirichletSetTag,
116 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get DIRICHLET_SET tag" );
117
118 result = impl->tag_get_handle( NEUMANN_SET_TAG_NAME, 1, MB_TYPE_INTEGER, mNeumannSetTag,
119 MB_TAG_CREAT | MB_TAG_SPARSE, &negone );MB_CHK_SET_ERR_RET( result, "Failed to get NEUMANN_SET tag" );
120
121 const int negonearr[] = { -1, -1, -1, -1 };
122 result = impl->tag_get_handle( HAS_MID_NODES_TAG_NAME, 4, MB_TYPE_INTEGER, mHasMidNodesTag,
123 MB_TAG_CREAT | MB_TAG_SPARSE, negonearr );MB_CHK_SET_ERR_RET( result, "Failed to get HAS_MID_NODES tag" );
124
125 mGlobalIdTag = impl->globalId_tag();
126
127 result =
128 impl->tag_get_handle( NAME_TAG_NAME, NAME_TAG_SIZE, MB_TYPE_OPAQUE, mNameTag, MB_TAG_CREAT | MB_TAG_SPARSE );MB_CHK_SET_ERR_RET( result, "Failed to get NAME tag" );
129 }
130
131 ReadCCMIO::~ReadCCMIO()
132 {
133 mbImpl->release_interface( readMeshIface );
134 }
135
136 ErrorCode ReadCCMIO::load_file( const char* file_name,
137 const EntityHandle* file_set,
138 const FileOptions& /* opts */,
139 const ReaderIface::SubsetList* subset_list,
140 const Tag* /* file_id_tag */ )
141 {
142 CCMIOID rootID, problemID, stateID, processorID, verticesID, topologyID, solutionID;
143 CCMIOError error = kCCMIONoErr;
144
145 if( subset_list )
146 {
147 MB_SET_ERR( MB_UNSUPPORTED_OPERATION, "Reading subset of files not supported for CCMOI data" );
148 }
149
150 CCMIOOpenFile( &error, file_name, kCCMIORead, &rootID );
151 CHK_SET_CCMERR( error, "Problem opening file" );
152
153 // Get the file state
154 ErrorCode rval = get_state( rootID, problemID, stateID );MB_CHK_SET_ERR( rval, "Failed to get state" );
155
156 // Get processors
157 std::vector< CCMIOSize_t > procs;
158 bool has_solution = false;
159 rval = get_processors( stateID, processorID, verticesID, topologyID, solutionID, procs, has_solution );MB_CHK_SET_ERR( rval, "Failed to get processors" );
160
161 std::vector< CCMIOSize_t >::iterator vit;
162 Range new_ents, *new_ents_ptr = NULL;
163 if( file_set ) new_ents_ptr = &new_ents;
164
165 for( vit = procs.begin(); vit != procs.end(); ++vit )
166 {
167 rval = read_processor( stateID, problemID, processorID, verticesID, topologyID, *vit, new_ents_ptr );MB_CHK_SET_ERR( rval, "Failed to read processors" );
168 }
169
170 // Load some meta-data
171 rval = load_metadata( rootID, problemID, stateID, processorID, file_set );MB_CHK_SET_ERR( rval, "Failed to load some meta-data" );
172
173 // Now, put all this into the file set, if there is one
174 if( file_set )
175 {
176 rval = mbImpl->add_entities( *file_set, new_ents );MB_CHK_SET_ERR( rval, "Failed to add new entities to file set" );
177 }
178
179 return rval;
180 }
181
182 ErrorCode ReadCCMIO::get_state( CCMIOID rootID, CCMIOID& problemID, CCMIOID& stateID )
183 {
184 CCMIOError error = kCCMIONoErr;
185
186 // First try default
187 CCMIOGetState( &error, rootID, "default", &problemID, &stateID );
188 if( kCCMIONoErr != error )
189 {
190 CCMIOSize_t i = CCMIOSIZEC( 0 );
191 CCMIOError tmp_error = kCCMIONoErr;
192 CCMIONextEntity( &tmp_error, rootID, kCCMIOState, &i, &stateID );
193 if( kCCMIONoErr == tmp_error ) CCMIONextEntity( &error, rootID, kCCMIOProblemDescription, &i, &problemID );
194 }
195 CHK_SET_CCMERR( error, "Couldn't find state" );
196
197 return MB_SUCCESS;
198 }
199
200 ErrorCode ReadCCMIO::load_metadata( CCMIOID rootID,
201 CCMIOID problemID,
202 CCMIOID /* stateID */,
203 CCMIOID processorID,
204 const EntityHandle* file_set )
205 {
206 // Read the simulation title.
207 CCMIOError error = kCCMIONoErr;
208 ErrorCode rval = MB_SUCCESS;
209 CCMIONode rootNode;
210 if( kCCMIONoErr == CCMIOGetEntityNode( &error, rootID, &rootNode ) )
211 {
212 char* name = NULL;
213 CCMIOGetTitle( &error, rootNode, &name );
214
215 if( NULL != name && strlen( name ) != 0 )
216 {
217 // Make a tag for it and tag the read set
218 Tag simname;
219 rval = mbImpl->tag_get_handle( "Title", strlen( name ), MB_TYPE_OPAQUE, simname,
220 MB_TAG_CREAT | MB_TAG_SPARSE );MB_CHK_SET_ERR( rval, "Simulation name tag not found or created" );
221 EntityHandle set = file_set ? *file_set : 0;
222 rval = mbImpl->tag_set_data( simname, &set, 1, name );MB_CHK_SET_ERR( rval, "Problem setting simulation name tag" );
223 }
224 if( name ) free( name );
225 }
226
227 // Creating program
228 EntityHandle dumh = ( file_set ? *file_set : 0 );
229 rval = get_str_option( "CreatingProgram", dumh, mCreatingProgramTag, processorID );MB_CHK_SET_ERR( rval, "Trouble getting CreatingProgram tag" );
230
231 rval = load_matset_data( problemID );MB_CHK_SET_ERR( rval, "Failure loading matset data" );
232
233 rval = load_neuset_data( problemID );MB_CHK_SET_ERR( rval, "Failure loading neuset data" );
234
235 return rval;
236 }
237
238 ErrorCode ReadCCMIO::load_matset_data( CCMIOID problemID )
239 {
240 // Make sure there are matsets
241 if( newMatsets.empty() ) return MB_SUCCESS;
242
243 // ... walk through each cell type
244 CCMIOSize_t i = CCMIOSIZEC( 0 );
245 CCMIOID next;
246 CCMIOError error = kCCMIONoErr;
247
248 while( CCMIONextEntity( NULL, problemID, kCCMIOCellType, &i, &next ) == kCCMIONoErr )
249 {
250 // Get index, corresponding set, and label with material set tag
251 int mindex;
252 CCMIOGetEntityIndex( &error, next, &mindex );
253 std::map< int, EntityHandle >::iterator mit = newMatsets.find( mindex );
254 if( mit == newMatsets.end() )
255 // No actual faces for this matset; continue to next
256 continue;
257
258 EntityHandle dum_ent = mit->second;
259 ErrorCode rval = mbImpl->tag_set_data( mMaterialSetTag, &dum_ent, 1, &mindex );MB_CHK_SET_ERR( rval, "Trouble setting material set tag" );
260
261 // Set name
262 CCMIOSize_t len;
263 CCMIOEntityLabel( &error, next, &len, NULL );
264 std::vector< char > opt_string2( GETINT32( len ) + 1, '\0' );
265 CCMIOEntityLabel( &error, next, NULL, &opt_string2[0] );
266 if( opt_string2.size() >= NAME_TAG_SIZE )
267 opt_string2[NAME_TAG_SIZE - 1] = '\0';
268 else
269 ( opt_string2.resize( NAME_TAG_SIZE, '\0' ) );
270 rval = mbImpl->tag_set_data( mNameTag, &dum_ent, 1, &opt_string2[0] );MB_CHK_SET_ERR( rval, "Trouble setting name tag for material set" );
271
272 // Material id
273 rval = get_int_option( "MaterialId", dum_ent, mMaterialIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting MaterialId tag" );
274
275 rval = get_str_option( "MaterialType", dum_ent, mMaterialTypeTag, next );MB_CHK_SET_ERR( rval, "Trouble getting MaterialType tag" );
276
277 rval = get_int_option( "Radiation", dum_ent, mRadiationTag, next );MB_CHK_SET_ERR( rval, "Trouble getting Radiation option" );
278
279 rval = get_int_option( "PorosityId", dum_ent, mPorosityIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting PorosityId option" );
280
281 rval = get_int_option( "SpinId", dum_ent, mSpinIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting SpinId option" );
282
283 rval = get_int_option( "GroupId", dum_ent, mGroupIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting GroupId option" );
284
285 rval = get_int_option( "ColorIdx", dum_ent, mColorIdxTag, next );MB_CHK_SET_ERR( rval, "Trouble getting ColorIdx option" );
286
287 rval = get_int_option( "ProcessorId", dum_ent, mProcessorIdTag, next );MB_CHK_SET_ERR( rval, "Trouble getting ProcessorId option" );
288
289 rval = get_int_option( "LightMaterial", dum_ent, mLightMaterialTag, next );MB_CHK_SET_ERR( rval, "Trouble getting LightMaterial option" );
290
291 rval = get_int_option( "FreeSurfaceMaterial", dum_ent, mFreeSurfaceMaterialTag, next );MB_CHK_SET_ERR( rval, "Trouble getting FreeSurfaceMaterial option" );
292
293 rval = get_dbl_option( "Thickness", dum_ent, mThicknessTag, next );MB_CHK_SET_ERR( rval, "Trouble getting Thickness option" );
294 }
295
296 return MB_SUCCESS;
297 }
298
299 ErrorCode ReadCCMIO::get_int_option( const char* opt_str, EntityHandle seth, Tag& tag, CCMIOID node )
300 {
301 int idum;
302 ErrorCode rval;
303 if( kCCMIONoErr == CCMIOReadOpti( NULL, node, opt_str, &idum ) )
304 {
305 if( !tag )
306 {
307 rval = mbImpl->tag_get_handle( opt_str, 1, MB_TYPE_INTEGER, tag, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" );
308 }
309
310 rval = mbImpl->tag_set_data( tag, &seth, 1, &idum );MB_CHK_SET_ERR( rval, "Failed to set tag data" );
311 }
312
313 return MB_SUCCESS;
314 }
315
316 ErrorCode ReadCCMIO::get_dbl_option( const char* opt_str, EntityHandle seth, Tag& tag, CCMIOID node )
317 {
318 float fdum;
319 if( kCCMIONoErr == CCMIOReadOptf( NULL, node, opt_str, &fdum ) )
320 {
321 ErrorCode rval;
322 if( !tag )
323 {
324 rval = mbImpl->tag_get_handle( opt_str, 1, MB_TYPE_DOUBLE, tag, MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" );
325 }
326
327 double dum_dbl = fdum;
328 rval = mbImpl->tag_set_data( tag, &seth, 1, &dum_dbl );MB_CHK_SET_ERR( rval, "Failed to set tag data" );
329 }
330
331 return MB_SUCCESS;
332 }
333
334 ErrorCode ReadCCMIO::get_str_option( const char* opt_str,
335 EntityHandle seth,
336 Tag& tag,
337 CCMIOID node,
338 const char* other_tag_name )
339 {
340 int len;
341 CCMIOError error = kCCMIONoErr;
342 std::vector< char > opt_string;
343 if( kCCMIONoErr != CCMIOReadOptstr( NULL, node, opt_str, &len, NULL ) ) return MB_SUCCESS;
344
345 opt_string.resize( len );
346 CCMIOReadOptstr( &error, node, opt_str, &len, &opt_string[0] );
347 ErrorCode rval = MB_SUCCESS;
348 if( !tag )
349 {
350 rval = mbImpl->tag_get_handle( other_tag_name ? other_tag_name : opt_str, NAME_TAG_SIZE, MB_TYPE_OPAQUE, tag,
351 MB_TAG_SPARSE | MB_TAG_CREAT );MB_CHK_SET_ERR( rval, "Failed to get tag handle" );
352 }
353
354 if( opt_string.size() > NAME_TAG_SIZE )
355 opt_string[NAME_TAG_SIZE - 1] = '\0';
356 else
357 ( opt_string.resize( NAME_TAG_SIZE, '\0' ) );
358
359 rval = mbImpl->tag_set_data( tag, &seth, 1, &opt_string[0] );MB_CHK_SET_ERR( rval, "Failed to set tag data" );
360
361 return MB_SUCCESS;
362 }
363
364 ErrorCode ReadCCMIO::load_neuset_data( CCMIOID problemID )
365 {
366 CCMIOSize_t i = CCMIOSIZEC( 0 );
367 CCMIOID next;
368
369 // Make sure there are matsets
370 if( newNeusets.empty() ) return MB_SUCCESS;
371
372 while( CCMIONextEntity( NULL, problemID, kCCMIOBoundaryRegion, &i, &next ) == kCCMIONoErr )
373 {
374 // Get index, corresponding set, and label with neumann set tag
375 int mindex;
376 CCMIOError error = kCCMIONoErr;
377 CCMIOGetEntityIndex( &error, next, &mindex );
378 std::map< int, EntityHandle >::iterator mit = newNeusets.find( mindex );
379 if( mit == newNeusets.end() )
380 // No actual faces for this neuset; continue to next
381 continue;
382
383 EntityHandle dum_ent = mit->second;
384 ErrorCode rval = mbImpl->tag_set_data( mNeumannSetTag, &dum_ent, 1, &mindex );MB_CHK_SET_ERR( rval, "Trouble setting neumann set tag" );
385
386 // Set name
387 rval = get_str_option( "BoundaryName", dum_ent, mNameTag, next, NAME_TAG_NAME );MB_CHK_SET_ERR( rval, "Trouble creating BoundaryName tag" );
388
389 // BoundaryType
390 rval = get_str_option( "BoundaryType", dum_ent, mBoundaryTypeTag, next );MB_CHK_SET_ERR( rval, "Trouble creating BoundaryType tag" );
391
392 // ProstarRegionNumber
393 rval = get_int_option( "ProstarRegionNumber", dum_ent, mProstarRegionNumberTag, next );MB_CHK_SET_ERR( rval, "Trouble creating ProstarRegionNumber tag" );
394 }
395
396 return MB_SUCCESS;
397 }
398
399 ErrorCode ReadCCMIO::read_processor( CCMIOID /* stateID */,
400 CCMIOID problemID,
401 CCMIOID processorID,
402 CCMIOID verticesID,
403 CCMIOID topologyID,
404 CCMIOSize_t proc,
405 Range* new_ents )
406 {
407 ErrorCode rval;
408
409 // vert_map fields: s: none, i: gid, ul: vert handle, r: none
410 // TupleList vert_map(0, 1, 1, 0, 0);
411 TupleList vert_map;
412 rval = read_vertices( proc, processorID, verticesID, topologyID, new_ents, vert_map );MB_CHK_SET_ERR( rval, "Failed to read vertices" );
413
414 rval = read_cells( proc, problemID, verticesID, topologyID, vert_map, new_ents );MB_CHK_SET_ERR( rval, "Failed to read cells" );
415
416 return rval;
417 }
418
419 ErrorCode ReadCCMIO::read_cells( CCMIOSize_t /* proc */,
420 CCMIOID problemID,
421 CCMIOID /* verticesID */,
422 CCMIOID topologyID,
423 TupleList& vert_map,
424 Range* new_ents )
425 {
426 // Read the faces.
427 // face_map fields: s:forward/reverse, i: cell id, ul: face handle, r: none
428 ErrorCode rval;
429 #ifdef TUPLE_LIST
430 TupleList face_map( 1, 1, 1, 0, 0 );
431 #else
432 TupleList face_map;
433 SenseList sense_map;
434 #endif
435 rval = read_all_faces( topologyID, vert_map, face_map,
436 #ifndef TUPLE_LIST
437 sense_map,
438 #endif
439 new_ents );MB_CHK_SET_ERR( rval, "Failed to read all cells" );
440
441 // Read the cell topology types, if any exist in the file
442 std::map< int, int > cell_topo_types;
443 rval = read_topology_types( topologyID, cell_topo_types );MB_CHK_SET_ERR( rval, "Problem reading cell topo types" );
444
445 // Now construct the cells; sort the face map by cell ids first
446 #ifdef TUPLE_LIST
447 rval = face_map.sort( 1 );MB_CHK_SET_ERR( rval, "Couldn't sort face map by cell id" );
448 #endif
449 std::vector< EntityHandle > new_cells;
450 rval = construct_cells( face_map,
451 #ifndef TUPLE_LIST
452 sense_map,
453 #endif
454 vert_map, cell_topo_types, new_cells );MB_CHK_SET_ERR( rval, "Failed to construct cells" );
455 if( new_ents )
456 {
457 Range::iterator rit = new_ents->end();
458 std::vector< EntityHandle >::reverse_iterator vit;
459 for( vit = new_cells.rbegin(); vit != new_cells.rend(); ++vit )
460 rit = new_ents->insert( rit, *vit );
461 }
462
463 rval = read_gids_and_types( problemID, topologyID, new_cells );MB_CHK_SET_ERR( rval, "Failed to read gids and types" );
464
465 return MB_SUCCESS;
466 }
467
468 ErrorCode ReadCCMIO::read_topology_types( CCMIOID& topologyID, std::map< int, int >& cell_topo_types )
469 {
470 CCMIOError error = kCCMIONoErr;
471 CCMIOID cellID, mapID;
472 CCMIOSize_t ncells;
473 CCMIOGetEntity( &error, topologyID, kCCMIOCells, 0, &cellID );
474 CCMIOEntitySize( &error, cellID, &ncells, NULL );
475 int num_cells = GETINT32( ncells );
476
477 // First, do a dummy read to see if we even have topo types in this mesh
478 int dum_int;
479 CCMIOReadOpt1i( &error, cellID, "CellTopologyType", &dum_int, CCMIOINDEXC( kCCMIOStart ),
480 CCMIOINDEXC( kCCMIOStart ) + 1 );
481 if( kCCMIONoErr != error ) return MB_SUCCESS;
482
483 // OK, we have topo types; first get the map node
484 std::vector< int > dum_ints( num_cells );
485 CCMIOReadCells( &error, cellID, &mapID, &dum_ints[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOStart ) + 1 );
486 CHK_SET_CCMERR( error, "Failed to get the map node" );
487
488 // Now read the map
489 CCMIOReadMap( &error, mapID, &dum_ints[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
490 CHK_SET_CCMERR( error, "Failed to get cell ids" );
491 int i;
492 for( i = 0; i < num_cells; i++ )
493 cell_topo_types[dum_ints[i]] = 0;
494
495 // Now read the cell topo types for real, reusing cell_topo_types
496 std::vector< int > topo_types( num_cells );
497 CCMIOReadOpt1i( &error, cellID, "CellTopologyType", &topo_types[0], CCMIOINDEXC( kCCMIOStart ),
498 CCMIOINDEXC( kCCMIOEnd ) );
499 CHK_SET_CCMERR( error, "Failed to get cell topo types" );
500 for( i = 0; i < num_cells; i++ )
501 cell_topo_types[dum_ints[i]] = topo_types[i];
502
503 return MB_SUCCESS;
504 }
505
506 ErrorCode ReadCCMIO::read_gids_and_types( CCMIOID /* problemID */,
507 CCMIOID topologyID,
508 std::vector< EntityHandle >& cells )
509 {
510 // Get the cells entity and number of cells
511 CCMIOSize_t dum_cells;
512 int num_cells;
513 CCMIOError error = kCCMIONoErr;
514 CCMIOID cellsID, mapID;
515 CCMIOGetEntity( &error, topologyID, kCCMIOCells, 0, &cellsID );
516 CCMIOEntitySize( &error, cellsID, &dum_cells, NULL );
517 num_cells = GETINT32( dum_cells );
518
519 // Check the number of cells against how many are in the cell array
520 if( num_cells != (int)cells.size() ) MB_SET_ERR( MB_FAILURE, "Number of cells doesn't agree" );
521
522 // Read the gid map and set global ids
523 std::vector< int > cell_gids( num_cells );
524 CCMIOReadCells( &error, cellsID, &mapID, NULL, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
525 CHK_SET_CCMERR( error, "Couldn't read cells" );
526 CCMIOReadMap( &error, mapID, &cell_gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
527 CHK_SET_CCMERR( error, "Couldn't read cell id map" );
528
529 ErrorCode rval = mbImpl->tag_set_data( mGlobalIdTag, &cells[0], cells.size(), &cell_gids[0] );MB_CHK_SET_ERR( rval, "Couldn't set gids tag" );
530
531 // Now read cell material types; reuse cell_gids
532 CCMIOReadCells( &error, cellsID, NULL, &cell_gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
533 CHK_SET_CCMERR( error, "Trouble reading cell types" );
534
535 // Create the matsets
536 std::map< int, Range > matset_ents;
537 for( int i = 0; i < num_cells; i++ )
538 matset_ents[cell_gids[i]].insert( cells[i] );
539
540 for( std::map< int, Range >::iterator mit = matset_ents.begin(); mit != matset_ents.end(); ++mit )
541 {
542 EntityHandle matset;
543 rval = mbImpl->create_meshset( MESHSET_SET, matset );MB_CHK_SET_ERR( rval, "Couldn't create material set" );
544 newMatsets[mit->first] = matset;
545
546 rval = mbImpl->add_entities( matset, mit->second );MB_CHK_SET_ERR( rval, "Couldn't add entities to material set" );
547 }
548
549 return MB_SUCCESS;
550 }
551
552 ErrorCode ReadCCMIO::construct_cells( TupleList& face_map,
553 #ifndef TUPLE_LIST
554 SenseList& sense_map,
555 #endif
556 TupleList& /* vert_map */,
557 std::map< int, int >& cell_topo_types,
558 std::vector< EntityHandle >& new_cells )
559 {
560 std::vector< EntityHandle > facehs;
561 std::vector< int > senses;
562 EntityHandle cell;
563 ErrorCode tmp_rval, rval = MB_SUCCESS;
564 EntityType this_type = MBMAXTYPE;
565 bool has_mid_nodes = false;
566 #ifdef TUPLE_LIST
567 unsigned int i = 0;
568 while( i < face_map.n )
569 {
570 // Pull out face handles bounding the same cell
571 facehs.clear();
572 int this_id = face_map.get_int( i );
573 unsigned int inext = i;
574 while( face_map.get_int( inext ) == this_id && inext <= face_map.n )
575 {
576 inext++;
577 EntityHandle face = face_map.get_ulong( inext );
578 facehs.push_back( face );
579 senses.push_back( face_map.get_short( inext ) );
580 }
581 this_type = MBMAXTYPE;
582 has_mid_nodes = false;
583 #else
584 std::map< int, std::vector< EntityHandle > >::iterator fmit;
585 std::map< int, std::vector< int > >::iterator smit;
586 std::map< int, int >::iterator typeit;
587 for( fmit = face_map.begin(), smit = sense_map.begin(); fmit != face_map.end(); ++fmit, ++smit )
588 {
589 // Pull out face handles bounding the same cell
590 facehs.clear();
591 int this_id = ( *fmit ).first;
592 facehs = ( *fmit ).second;
593 senses.clear();
594 senses = ( *smit ).second;
595 typeit = cell_topo_types.find( this_id );
596 if( typeit != cell_topo_types.end() )
597 {
598 rval = ccmio_to_moab_type( typeit->second, this_type, has_mid_nodes );
599 }
600 else
601 {
602 this_type = MBMAXTYPE;
603 has_mid_nodes = false;
604 }
605 #endif
606 tmp_rval = create_cell_from_faces( facehs, senses, this_type, has_mid_nodes, cell );
607 if( MB_SUCCESS != tmp_rval )
608 rval = tmp_rval;
609 else
610 {
611 new_cells.push_back( cell );
612 // Tag cell with global id
613 tmp_rval = mbImpl->tag_set_data( mGlobalIdTag, &cell, 1, &this_id );
614 if( MB_SUCCESS != tmp_rval ) rval = tmp_rval;
615 }
616 }
617
618 return rval;
619 }
620
621 ErrorCode ReadCCMIO::ccmio_to_moab_type( int ccm_type, EntityType& moab_type, bool& has_mid_nodes )
622 {
623 switch( ccm_type )
624 {
625 case 1:
626 moab_type = MBVERTEX;
627 break;
628 case 2:
629 case 28:
630 moab_type = MBEDGE;
631 break;
632 case 29:
633 moab_type = MBMAXTYPE;
634 break;
635 case 3:
636 case 4:
637 moab_type = MBQUAD;
638 break;
639 case 11:
640 case 21:
641 moab_type = MBHEX;
642 break;
643 case 12:
644 case 22:
645 moab_type = MBPRISM;
646 break;
647 case 13:
648 case 23:
649 moab_type = MBTET;
650 break;
651 case 14:
652 case 24:
653 moab_type = MBPYRAMID;
654 break;
655 case 255:
656 moab_type = MBPOLYHEDRON;
657 break;
658 default:
659 moab_type = MBMAXTYPE;
660 }
661
662 switch( ccm_type )
663 {
664 case 28:
665 case 4:
666 case 21:
667 case 22:
668 case 23:
669 case 24:
670 has_mid_nodes = true;
671 break;
672 default:
673 break;
674 }
675
676 return MB_SUCCESS;
677 }
678
679 ErrorCode ReadCCMIO::create_cell_from_faces( std::vector< EntityHandle >& facehs,
680 std::vector< int >& senses,
681 EntityType this_type,
682 bool /* has_mid_nodes */,
683 EntityHandle& cell )
684 {
685 ErrorCode rval;
686
687 // Test up front to see if they're one type
688 EntityType face_type = mbImpl->type_from_handle( facehs[0] );
689 bool same_type = true;
690 for( std::vector< EntityHandle >::iterator vit = facehs.begin(); vit != facehs.end(); ++vit )
691 {
692 if( face_type != mbImpl->type_from_handle( *vit ) )
693 {
694 same_type = false;
695 break;
696 }
697 }
698
699 std::vector< EntityHandle > verts;
700 EntityType input_type = this_type;
701 std::vector< EntityHandle > storage;
702 MeshTopoUtil mtu( mbImpl );
703
704 // Preset this to maxtype, so we get an affirmative choice in loop
705 this_type = MBMAXTYPE;
706
707 if( ( MBTET == input_type || MBMAXTYPE == input_type ) && same_type && face_type == MBTRI && facehs.size() == 4 )
708 {
709 // Try to get proper connectivity for tet
710
711 // Get connectivity of first face, and reverse it if sense is forward, since
712 // base face always points into entity
713 rval = mbImpl->get_connectivity( &facehs[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
714 if( senses[0] > 0 ) std::reverse( verts.begin(), verts.end() );
715
716 // Get the 4th vertex through the next tri
717 const EntityHandle* conn;
718 int conn_size;
719 rval = mbImpl->get_connectivity( facehs[1], conn, conn_size, true, &storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
720 int i = 0;
721 while( std::find( verts.begin(), verts.end(), conn[i] ) != verts.end() && i < conn_size )
722 i++;
723
724 // If i is not at the end of the verts, found the apex; otherwise fall back to polyhedron
725 if( conn_size != i )
726 {
727 this_type = MBTET;
728 verts.push_back( conn[i] );
729 }
730 }
731 else if( ( MBHEX == input_type || MBMAXTYPE == input_type ) && same_type && MBQUAD == face_type &&
732 facehs.size() == 6 )
733 {
734 // Build hex from quads
735 // Algorithm:
736 // - verts = vertices from 1st quad
737 // - Find quad q1 sharing verts[0] and verts[1]
738 // - Find quad q2 sharing other 2 verts in q1
739 // - Find v1 = opposite vert from verts[1] in q1 , v2 = opposite from verts[0]
740 // - Get i = offset of v1 in verts2 of q2, rotate verts2 by i
741 // - If verts2[(i + 1) % 4] != v2, flip verts2 by switching verts2[1] and verts2[3]
742 // - append verts2 to verts
743
744 // Get the other vertices for this hex; need to find the quad with no common vertices
745 Range tmp_faces, tmp_verts;
746 // Get connectivity of first face, and reverse it if sense is forward, since
747 // base face always points into entity
748 rval = mbImpl->get_connectivity( &facehs[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
749 if( senses[0] > 0 ) std::reverse( verts.begin(), verts.end() );
750
751 // Get q1, which shares 2 vertices with q0
752 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) );
753 rval = mbImpl->get_adjacencies( &verts[0], 2, 2, false, tmp_faces );
754 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" );
755 tmp_faces.erase( facehs[0] );
756 EntityHandle q1 = *tmp_faces.begin();
757 // Get other 2 verts of q1
758 rval = mbImpl->get_connectivity( &q1, 1, tmp_verts );MB_CHK_SET_ERR( rval, "Couldn't get adj verts" );
759 tmp_verts.erase( verts[0] );
760 tmp_verts.erase( verts[1] );
761 // Get q2
762 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) );
763 rval = mbImpl->get_adjacencies( tmp_verts, 2, false, tmp_faces );
764 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" );
765 tmp_faces.erase( q1 );
766 EntityHandle q2 = *tmp_faces.begin();
767 // Get verts in q2
768 rval = mbImpl->get_connectivity( &q2, 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get adj vertices" );
769
770 // Get verts in q1 opposite from v[1] and v[0] in q0
771 EntityHandle v0 = 0, v1 = 0;
772 rval = mtu.opposite_entity( q1, verts[1], v0 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" );
773 rval = mtu.opposite_entity( q1, verts[0], v1 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" );
774 if( v0 && v1 )
775 {
776 // Offset of v0 in q2, then rotate and flip
777 unsigned int ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin();
778 if( 4 == ioff ) MB_SET_ERR( MB_FAILURE, "Trouble finding offset" );
779
780 if( storage[( ioff + 1 ) % 4] != v1 )
781 {
782 std::reverse( storage.begin(), storage.end() );
783 ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin();
784 }
785 if( 0 != ioff ) std::rotate( storage.begin(), storage.begin() + ioff, storage.end() );
786
787 // Copy into verts, and make hex
788 std::copy( storage.begin(), storage.end(), std::back_inserter( verts ) );
789 this_type = MBHEX;
790 }
791 }
792
793 if( MBMAXTYPE == this_type && facehs.size() == 5 )
794 {
795 // Some preliminaries
796 std::vector< EntityHandle > tris, quads;
797 for( unsigned int i = 0; i < 5; i++ )
798 {
799 if( MBTRI == mbImpl->type_from_handle( facehs[i] ) )
800 tris.push_back( facehs[i] );
801 else if( MBQUAD == mbImpl->type_from_handle( facehs[i] ) )
802 quads.push_back( facehs[i] );
803 }
804
805 // Check for prisms
806 if( 2 == tris.size() && 3 == quads.size() )
807 {
808 // OK, we have the right number of tris and quads; try to find the proper verts
809
810 // Get connectivity of first tri, and reverse if necessary
811 int index = std::find( facehs.begin(), facehs.end(), tris[0] ) - facehs.begin();
812 rval = mbImpl->get_connectivity( &tris[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
813 if( senses[index] > 0 ) std::reverse( verts.begin(), verts.end() );
814
815 // Now align vertices of other tri, through a quad, similar to how we did hexes
816 // Get q1, which shares 2 vertices with t0
817 Range tmp_faces, tmp_verts;
818 std::copy( facehs.begin(), facehs.end(), range_inserter( tmp_faces ) );
819 rval = mbImpl->get_adjacencies( &verts[0], 2, 2, false, tmp_faces );
820 if( MB_SUCCESS != rval || tmp_faces.size() != 2 ) MB_SET_ERR( MB_FAILURE, "Couldn't get adj face" );
821 tmp_faces.erase( tris[0] );
822 EntityHandle q1 = *tmp_faces.begin();
823 // Get verts in q1
824 rval = mbImpl->get_connectivity( &q1, 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get adj vertices" );
825
826 // Get verts in q1 opposite from v[1] and v[0] in q0
827 EntityHandle v0 = 0, v1 = 0;
828 rval = mtu.opposite_entity( q1, verts[1], v0 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" );
829 rval = mtu.opposite_entity( q1, verts[0], v1 );MB_CHK_SET_ERR( rval, "Couldn't get the opposite side entity" );
830 if( v0 && v1 )
831 {
832 // Offset of v0 in t2, then rotate and flip
833 storage.clear();
834 rval = mbImpl->get_connectivity( &tris[1], 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
835
836 index = std::find( facehs.begin(), facehs.end(), tris[1] ) - facehs.begin();
837 if( senses[index] < 0 ) std::reverse( storage.begin(), storage.end() );
838 unsigned int ioff = std::find( storage.begin(), storage.end(), v0 ) - storage.begin();
839 if( 3 == ioff ) MB_SET_ERR( MB_FAILURE, "Trouble finding offset" );
840 for( unsigned int i = 0; i < 3; i++ )
841 verts.push_back( storage[( ioff + i ) % 3] );
842
843 this_type = MBPRISM;
844 }
845 }
846 else if( tris.size() == 4 && quads.size() == 1 )
847 {
848 // Check for pyramid
849 // Get connectivity of first tri, and reverse if necessary
850 int index = std::find( facehs.begin(), facehs.end(), quads[0] ) - facehs.begin();
851 rval = mbImpl->get_connectivity( &quads[0], 1, verts );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
852 if( senses[index] > 0 ) std::reverse( verts.begin(), verts.end() );
853
854 // Get apex node
855 rval = mbImpl->get_connectivity( &tris[0], 1, storage );MB_CHK_SET_ERR( rval, "Couldn't get connectivity" );
856 for( unsigned int i = 0; i < 3; i++ )
857 {
858 if( std::find( verts.begin(), verts.end(), storage[i] ) == verts.end() )
859 {
860 verts.push_back( storage[i] );
861 break;
862 }
863 }
864
865 if( 5 == verts.size() ) this_type = MBPYRAMID;
866 }
867 else
868 {
869 // Dummy else clause to stop in debugger
870 this_type = MBMAXTYPE;
871 }
872 }
873
874 if( MBMAXTYPE != input_type && input_type != this_type && this_type != MBMAXTYPE )
875 std::cerr << "Warning: types disagree (cell_topo_type = " << CN::EntityTypeName( input_type )
876 << ", faces indicate type " << CN::EntityTypeName( this_type ) << std::endl;
877
878 if( MBMAXTYPE != input_type && this_type == MBMAXTYPE && input_type != MBPOLYHEDRON )
879 std::cerr << "Warning: couldn't find proper connectivity for specified topo_type = "
880 << CN::EntityTypeName( input_type ) << std::endl;
881
882 // Now make the element; if we fell back to polyhedron, use faces, otherwise use verts
883 if( MBPOLYHEDRON == this_type || MBMAXTYPE == this_type )
884 {
885 rval = mbImpl->create_element( MBPOLYHEDRON, &facehs[0], facehs.size(), cell );MB_CHK_SET_ERR( rval, "create_element failed" );
886 }
887 else
888 {
889 rval = mbImpl->create_element( this_type, &verts[0], verts.size(), cell );MB_CHK_SET_ERR( rval, "create_element failed" );
890 }
891
892 return MB_SUCCESS;
893 }
894
895 ErrorCode ReadCCMIO::read_all_faces( CCMIOID topologyID,
896 TupleList& vert_map,
897 TupleList& face_map,
898 #ifndef TUPLE_LIST
899 SenseList& sense_map,
900 #endif
901 Range* new_faces )
902 {
903 CCMIOSize_t index;
904 CCMIOID faceID;
905 ErrorCode rval;
906 CCMIOError error = kCCMIONoErr;
907
908 // Get total # internal/bdy faces, size the face map accordingly
909 #ifdef TUPLE_LIST
910 index = CCMIOSIZEC( 0 );
911 int nbdy_faces = 0;
912 CCMIOSize_t nf;
913 error = kCCMIONoErr;
914 while( kCCMIONoErr == CCMIONextEntity( NULL, topologyID, kCCMIOBoundaryFaces, &index, &faceID ) )
915 {
916 CCMIOEntitySize( &error, faceID, &nf, NULL );
917 nbdy_faces += nf;
918 }
919 CCMIOGetEntity( &error, topologyID, kCCMIOInternalFaces, 0, &faceID );
920 CCMIOEntitySize( &error, faceID, &nf, NULL );
921
922 int nint_faces = nf;
923 face_map.resize( 2 * nint_faces + nbdy_faces );
924 #endif
925
926 // Get multiple blocks of bdy faces
927 index = CCMIOSIZEC( 0 );
928 while( kCCMIONoErr == CCMIONextEntity( NULL, topologyID, kCCMIOBoundaryFaces, &index, &faceID ) )
929 {
930 rval = read_faces( faceID, kCCMIOBoundaryFaces, vert_map, face_map,
931 #ifndef TUPLE_LIST
932 sense_map,
933 #endif
934 new_faces );MB_CHK_SET_ERR( rval, "Trouble reading boundary faces" );
935 }
936
937 // Now get internal faces
938 CCMIOGetEntity( &error, topologyID, kCCMIOInternalFaces, 0, &faceID );
939 CHK_SET_CCMERR( error, "Couldn't get internal faces" );
940
941 rval = read_faces( faceID, kCCMIOInternalFaces, vert_map, face_map,
942 #ifndef TUPLE_LIST
943 sense_map,
944 #endif
945 new_faces );MB_CHK_SET_ERR( rval, "Trouble reading internal faces" );
946
947 return rval;
948 }
949
950 ErrorCode ReadCCMIO::read_faces( CCMIOID faceID,
951 CCMIOEntity bdy_or_int,
952 TupleList& vert_map,
953 TupleList& face_map,
954 #ifndef TUPLE_LIST
955 SenseList& sense_map,
956 #endif
957 Range* new_faces )
958 {
959 if( kCCMIOInternalFaces != bdy_or_int && kCCMIOBoundaryFaces != bdy_or_int )
960 MB_SET_ERR( MB_FAILURE, "Face type isn't boundary or internal" );
961
962 CCMIOSize_t dum_faces;
963 CCMIOError error = kCCMIONoErr;
964 CCMIOEntitySize( &error, faceID, &dum_faces, NULL );
965 int num_faces = GETINT32( dum_faces );
966
967 // Get the size of the face connectivity array (not really a straight connect
968 // array, has n, connect(n), ...)
969 CCMIOSize_t farray_size = CCMIOSIZEC( 0 );
970 CCMIOReadFaces( &error, faceID, bdy_or_int, NULL, &farray_size, NULL, CCMIOINDEXC( kCCMIOStart ),
971 CCMIOINDEXC( kCCMIOEnd ) );
972 CHK_SET_CCMERR( error, "Trouble reading face connectivity length" );
973
974 // Allocate vectors for holding farray and cells for each face; use new for finer
975 // control of de-allocation
976 int num_sides = ( kCCMIOInternalFaces == bdy_or_int ? 2 : 1 );
977 int* farray = new int[GETINT32( farray_size )];
978
979 // Read farray and make the faces
980 CCMIOID mapID;
981 CCMIOReadFaces( &error, faceID, bdy_or_int, &mapID, NULL, farray, CCMIOINDEXC( kCCMIOStart ),
982 CCMIOINDEXC( kCCMIOEnd ) );
983 CHK_SET_CCMERR( error, "Trouble reading face connectivity" );
984
985 std::vector< EntityHandle > face_handles;
986 ErrorCode rval = make_faces( farray, vert_map, face_handles, num_faces );MB_CHK_SET_ERR( rval, "Failed to make the faces" );
987
988 // Read face cells and make tuples
989 int* face_cells;
990 if( num_sides * num_faces < farray_size )
991 face_cells = new int[num_sides * num_faces];
992 else
993 face_cells = farray;
994 CCMIOReadFaceCells( &error, faceID, bdy_or_int, face_cells, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
995 CHK_SET_CCMERR( error, "Trouble reading face cells" );
996
997 int* tmp_ptr = face_cells;
998 for( unsigned int i = 0; i < face_handles.size(); i++ )
999 {
1000 #ifdef TUPLE_LIST
1001 short forward = 1, reverse = -1;
1002 face_map.push_back( &forward, tmp_ptr++, &face_handles[i], NULL );
1003 if( 2 == num_sides ) face_map.push_back( &reverse, tmp_ptr++, &face_handles[i], NULL );
1004 #else
1005 face_map[*tmp_ptr].push_back( face_handles[i] );
1006 sense_map[*tmp_ptr++].push_back( 1 );
1007 if( 2 == num_sides )
1008 {
1009 face_map[*tmp_ptr].push_back( face_handles[i] );
1010 sense_map[*tmp_ptr++].push_back( -1 );
1011 }
1012 #endif
1013 }
1014
1015 // Now read & set face gids, reuse face_cells 'cuz we know it's big enough
1016 CCMIOReadMap( &error, mapID, face_cells, CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
1017 CHK_SET_CCMERR( error, "Trouble reading face gids" );
1018
1019 rval = mbImpl->tag_set_data( mGlobalIdTag, &face_handles[0], face_handles.size(), face_cells );MB_CHK_SET_ERR( rval, "Couldn't set face global ids" );
1020
1021 // Make a neumann set for these faces if they're all in a boundary face set
1022 if( kCCMIOBoundaryFaces == bdy_or_int )
1023 {
1024 EntityHandle neuset;
1025 rval = mbImpl->create_meshset( MESHSET_SET, neuset );MB_CHK_SET_ERR( rval, "Failed to create neumann set" );
1026
1027 // Don't trust entity index passed in
1028 int index;
1029 CCMIOGetEntityIndex( &error, faceID, &index );
1030 newNeusets[index] = neuset;
1031
1032 rval = mbImpl->add_entities( neuset, &face_handles[0], face_handles.size() );MB_CHK_SET_ERR( rval, "Failed to add faces to neumann set" );
1033
1034 // Now tag as neumann set; will add id later
1035 int dum_val = 0;
1036 rval = mbImpl->tag_set_data( mNeumannSetTag, &neuset, 1, &dum_val );MB_CHK_SET_ERR( rval, "Failed to tag neumann set" );
1037 }
1038
1039 if( new_faces )
1040 {
1041 std::sort( face_handles.begin(), face_handles.end() );
1042 std::copy( face_handles.rbegin(), face_handles.rend(), range_inserter( *new_faces ) );
1043 }
1044
1045 return MB_SUCCESS;
1046 }
1047
1048 ErrorCode ReadCCMIO::make_faces( int* farray,
1049 TupleList& vert_map,
1050 std::vector< EntityHandle >& new_faces,
1051 int num_faces )
1052 {
1053 std::vector< EntityHandle > verts;
1054 ErrorCode tmp_rval = MB_SUCCESS, rval = MB_SUCCESS;
1055
1056 for( int i = 0; i < num_faces; i++ )
1057 {
1058 int num_verts = *farray++;
1059 verts.resize( num_verts );
1060
1061 // Fill in connectivity by looking up by gid in vert tuple_list
1062 for( int j = 0; j < num_verts; j++ )
1063 {
1064 #ifdef TUPLE_LIST
1065 int tindex = vert_map.find( 1, farray[j] );
1066 if( -1 == tindex )
1067 {
1068 tmp_rval = MB_FAILURE;
1069 break;
1070 }
1071 verts[j] = vert_map.get_ulong( tindex, 0 );
1072 #else
1073 verts[j] = ( vert_map[farray[j]] )[0];
1074 #endif
1075 }
1076 farray += num_verts;
1077
1078 if( MB_SUCCESS == tmp_rval )
1079 {
1080 // Make face
1081 EntityType ftype = ( 3 == num_verts ? MBTRI : ( 4 == num_verts ? MBQUAD : MBPOLYGON ) );
1082 EntityHandle faceh;
1083 tmp_rval = mbImpl->create_element( ftype, &verts[0], num_verts, faceh );
1084 if( faceh ) new_faces.push_back( faceh );
1085 }
1086
1087 if( MB_SUCCESS != tmp_rval ) rval = tmp_rval;
1088 }
1089
1090 return rval;
1091 }
1092
1093 ErrorCode ReadCCMIO::read_vertices( CCMIOSize_t /* proc */,
1094 CCMIOID /* processorID */,
1095 CCMIOID verticesID,
1096 CCMIOID /* topologyID */,
1097 Range* verts,
1098 TupleList& vert_map )
1099 {
1100 CCMIOError error = kCCMIONoErr;
1101
1102 // Pre-read the number of vertices, so we can pre-allocate & read directly in
1103 CCMIOSize_t nverts = CCMIOSIZEC( 0 );
1104 CCMIOEntitySize( &error, verticesID, &nverts, NULL );
1105 CHK_SET_CCMERR( error, "Couldn't get number of vertices" );
1106
1107 // Get # dimensions
1108 CCMIOSize_t dims;
1109 float scale;
1110 CCMIOReadVerticesf( &error, verticesID, &dims, NULL, NULL, NULL, CCMIOINDEXC( 0 ), CCMIOINDEXC( 1 ) );
1111 CHK_SET_CCMERR( error, "Couldn't get number of dimensions" );
1112
1113 // Allocate vertex space
1114 EntityHandle node_handle = 0;
1115 std::vector< double* > arrays;
1116 readMeshIface->get_node_coords( 3, GETINT32( nverts ), MB_START_ID, node_handle, arrays );
1117
1118 // Read vertex coords
1119 CCMIOID mapID;
1120 std::vector< double > tmp_coords( GETINT32( dims ) * GETINT32( nverts ) );
1121 CCMIOReadVerticesd( &error, verticesID, &dims, &scale, &mapID, &tmp_coords[0], CCMIOINDEXC( 0 ),
1122 CCMIOINDEXC( 0 + nverts ) );
1123 CHK_SET_CCMERR( error, "Trouble reading vertex coordinates" );
1124
1125 // Copy interleaved coords into moab blocked coordinate space
1126 int i = 0, threei = 0;
1127 for( ; i < nverts; i++ )
1128 {
1129 arrays[0][i] = tmp_coords[threei++];
1130 arrays[1][i] = tmp_coords[threei++];
1131 if( 3 == GETINT32( dims ) )
1132 arrays[2][i] = tmp_coords[threei++];
1133 else
1134 arrays[2][i] = 0.0;
1135 }
1136
1137 // Scale, if necessary
1138 if( 1.0 != scale )
1139 {
1140 for( i = 0; i < nverts; i++ )
1141 {
1142 arrays[0][i] *= scale;
1143 arrays[1][i] *= scale;
1144 if( 3 == GETINT32( dims ) ) arrays[2][i] *= scale;
1145 }
1146 }
1147
1148 // Read gids for vertices
1149 std::vector< int > gids( GETINT32( nverts ) );
1150 CCMIOReadMap( &error, mapID, &gids[0], CCMIOINDEXC( kCCMIOStart ), CCMIOINDEXC( kCCMIOEnd ) );
1151 CHK_SET_CCMERR( error, "Trouble reading vertex global ids" );
1152
1153 // Put new vertex handles into range, and set gids for them
1154 Range new_verts( node_handle, node_handle + nverts - 1 );
1155 ErrorCode rval = mbImpl->tag_set_data( mGlobalIdTag, new_verts, &gids[0] );MB_CHK_SET_ERR( rval, "Couldn't set gids on vertices" );
1156
1157 // Pack vert_map with global ids and handles for these vertices
1158 #ifdef TUPLE_LIST
1159 vert_map.resize( GETINT32( nverts ) );
1160 for( i = 0; i < GETINT32( nverts ); i++ )
1161 {
1162 vert_map.push_back( NULL, &gids[i], &node_handle, NULL );
1163 #else
1164 for( i = 0; i < GETINT32( nverts ); i++ )
1165 {
1166 ( vert_map[gids[i]] ).push_back( node_handle );
1167 #endif
1168 node_handle += 1;
1169 }
1170
1171 if( verts ) verts->merge( new_verts );
1172
1173 return MB_SUCCESS;
1174 }
1175
1176 ErrorCode ReadCCMIO::get_processors( CCMIOID stateID,
1177 CCMIOID& processorID,
1178 CCMIOID& verticesID,
1179 CCMIOID& topologyID,
1180 CCMIOID& solutionID,
1181 std::vector< CCMIOSize_t >& procs,
1182 bool& /* has_solution */ )
1183 {
1184 CCMIOSize_t proc = CCMIOSIZEC( 0 );
1185 CCMIOError error = kCCMIONoErr;
1186
1187 CCMIONextEntity( &error, stateID, kCCMIOProcessor, &proc, &processorID );
1188 CHK_SET_CCMERR( error, "CCMIONextEntity() failed" );
1189 if( CCMIOReadProcessor( NULL, processorID, &verticesID, &topologyID, NULL, &solutionID ) != kCCMIONoErr )
1190 {
1191 // Maybe no solution; try again
1192 CCMIOReadProcessor( &error, processorID, &verticesID, &topologyID, NULL, NULL );
1193 CHK_SET_CCMERR( error, "CCMIOReadProcessor() failed" );
1194 hasSolution = false;
1195 }
1196
1197 procs.push_back( proc );
1198
1199 return MB_SUCCESS;
1200 }
1201
1202 ErrorCode ReadCCMIO::read_tag_values( const char* /* file_name */,
1203 const char* /* tag_name */,
1204 const FileOptions& /* opts */,
1205 std::vector< int >& /* tag_values_out */,
1206 const SubsetList* /* subset_list */ )
1207 {
1208 return MB_FAILURE;
1209 }
1210
1211 } // namespace moab
1.9.1.