Branch data Line data Source code
1 : : /**
2 : : * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3 : : * storing and accessing finite element mesh data.
4 : : *
5 : : * Copyright 2004 Sandia Corporation. Under the terms of Contract
6 : : * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7 : : * retains certain rights in this software.
8 : : *
9 : : * This library is free software; you can redistribute it and/or
10 : : * modify it under the terms of the GNU Lesser General Public
11 : : * License as published by the Free Software Foundation; either
12 : : * version 2.1 of the License, or (at your option) any later version.
13 : : *
14 : : */
15 : :
16 : : #ifdef __MFC_VER
17 : : #pragma warning( disable : 4786 )
18 : : #endif
19 : :
20 : : #ifdef WIN32 /* windows */
21 : : #define _USE_MATH_DEFINES // For M_PI
22 : : #endif
23 : : #include "moab/Skinner.hpp"
24 : : #include "moab/Range.hpp"
25 : : #include "moab/CN.hpp"
26 : : #include <vector>
27 : : #include <set>
28 : : #include <algorithm>
29 : : #include <math.h>
30 : : #include <assert.h>
31 : : #include <iostream>
32 : : #include "moab/Util.hpp"
33 : : #include "Internals.hpp"
34 : : #include "MBTagConventions.hpp"
35 : : #include "moab/Core.hpp"
36 : : #include "AEntityFactory.hpp"
37 : : #include "moab/ScdInterface.hpp"
38 : :
39 : : #ifdef M_PI
40 : : #define SKINNER_PI M_PI
41 : : #else
42 : : #define SKINNER_PI 3.1415926535897932384626
43 : : #endif
44 : :
45 : : namespace moab
46 : : {
47 : :
48 : 38 : Skinner::~Skinner()
49 : : {
50 : : // delete the adjacency tag
51 : 38 : }
52 : :
53 : 0 : ErrorCode Skinner::initialize()
54 : : {
55 : : // go through and mark all the target dimension entities
56 : : // that already exist as not deleteable
57 : : // also get the connectivity tags for each type
58 : : // also populate adjacency information
59 : : EntityType type;
60 : 0 : DimensionPair target_ent_types = CN::TypeDimensionMap[mTargetDim];
61 : :
62 : 0 : void* null_ptr = NULL;
63 : :
64 : : ErrorCode result = thisMB->tag_get_handle( "skinner adj", sizeof( void* ), MB_TYPE_OPAQUE, mAdjTag,
65 [ # # ][ # # ]: 0 : MB_TAG_DENSE | MB_TAG_CREAT, &null_ptr );MB_CHK_ERR( result );
[ # # ][ # # ]
66 : :
67 [ # # ]: 0 : if( mDeletableMBTag == 0 )
68 : : {
69 : : result =
70 [ # # ][ # # ]: 0 : thisMB->tag_get_handle( "skinner deletable", 1, MB_TYPE_BIT, mDeletableMBTag, MB_TAG_BIT | MB_TAG_CREAT );MB_CHK_ERR( result );
[ # # ][ # # ]
71 : : }
72 : :
73 [ # # ]: 0 : Range entities;
74 : :
75 : : // go through each type at this dimension
76 [ # # ][ # # ]: 0 : for( type = target_ent_types.first; type <= target_ent_types.second; ++type )
77 : : {
78 : : // get the entities of this type in the MB
79 [ # # ]: 0 : thisMB->get_entities_by_type( 0, type, entities );
80 : :
81 : : // go through each entity of this type in the MB
82 : : // and set its deletable tag to NO
83 [ # # ][ # # ]: 0 : Range::iterator iter, end_iter;
84 [ # # ]: 0 : end_iter = entities.end();
85 [ # # ][ # # ]: 0 : for( iter = entities.begin(); iter != end_iter; ++iter )
[ # # ][ # # ]
86 : : {
87 : 0 : unsigned char bit = 0x1;
88 [ # # ][ # # ]: 0 : result = thisMB->tag_set_data( mDeletableMBTag, &( *iter ), 1, &bit );
89 [ # # ]: 0 : assert( MB_SUCCESS == result );
90 : : // add adjacency information too
91 [ # # ][ # # ]: 0 : if( TYPE_FROM_HANDLE( *iter ) != MBVERTEX ) add_adjacency( *iter );
[ # # ][ # # ]
[ # # ]
92 : : }
93 : : }
94 : :
95 : 0 : return MB_SUCCESS;
96 : : }
97 : :
98 : 0 : ErrorCode Skinner::deinitialize()
99 : : {
100 : : ErrorCode result;
101 [ # # ]: 0 : if( 0 != mDeletableMBTag )
102 : : {
103 [ # # ]: 0 : result = thisMB->tag_delete( mDeletableMBTag );
104 [ # # ][ # # ]: 0 : mDeletableMBTag = 0;MB_CHK_ERR( result );
[ # # ]
105 : : }
106 : :
107 : : // remove the adjacency tag
108 [ # # ]: 0 : std::vector< std::vector< EntityHandle >* > adj_arr;
109 : 0 : std::vector< std::vector< EntityHandle >* >::iterator i;
110 [ # # ]: 0 : if( 0 != mAdjTag )
111 : : {
112 [ # # ][ # # ]: 0 : for( EntityType t = MBVERTEX; t != MBMAXTYPE; ++t )
113 : : {
114 [ # # ]: 0 : Range entities;
115 [ # # ][ # # ]: 0 : result = thisMB->get_entities_by_type_and_tag( 0, t, &mAdjTag, 0, 1, entities );MB_CHK_ERR( result );
[ # # ][ # # ]
116 [ # # ][ # # ]: 0 : adj_arr.resize( entities.size() );
117 [ # # ][ # # ]: 0 : result = thisMB->tag_get_data( mAdjTag, entities, &adj_arr[0] );MB_CHK_ERR( result );
[ # # ][ # # ]
[ # # ]
118 [ # # ][ # # ]: 0 : for( i = adj_arr.begin(); i != adj_arr.end(); ++i )
[ # # ][ # # ]
119 [ # # ][ # # ]: 0 : delete *i;
120 : 0 : }
121 : :
122 [ # # ]: 0 : result = thisMB->tag_delete( mAdjTag );
123 [ # # ][ # # ]: 0 : mAdjTag = 0;MB_CHK_ERR( result );
[ # # ]
124 : : }
125 : :
126 : 0 : return MB_SUCCESS;
127 : : }
128 : :
129 : 0 : ErrorCode Skinner::add_adjacency( EntityHandle entity )
130 : : {
131 : 0 : std::vector< EntityHandle >* adj = NULL;
132 : : const EntityHandle* nodes;
133 : : int num_nodes;
134 [ # # ][ # # ]: 0 : ErrorCode result = thisMB->get_connectivity( entity, nodes, num_nodes, true );MB_CHK_ERR( result );
[ # # ][ # # ]
135 [ # # ]: 0 : const EntityHandle* iter = std::min_element( nodes, nodes + num_nodes );
136 : :
137 [ # # ]: 0 : if( iter == nodes + num_nodes ) return MB_SUCCESS;
138 : :
139 : : // add this entity to the node
140 [ # # ][ # # ]: 0 : if( thisMB->tag_get_data( mAdjTag, iter, 1, &adj ) == MB_SUCCESS && adj != NULL ) { adj->push_back( entity ); }
[ # # ][ # # ]
[ # # ]
141 : : // create a new vector and add it
142 : : else
143 : : {
144 [ # # ][ # # ]: 0 : adj = new std::vector< EntityHandle >;
145 [ # # ]: 0 : adj->push_back( entity );
146 [ # # ][ # # ]: 0 : result = thisMB->tag_set_data( mAdjTag, iter, 1, &adj );MB_CHK_ERR( result );
[ # # ][ # # ]
147 : : }
148 : :
149 : 0 : return MB_SUCCESS;
150 : : }
151 : :
152 : 0 : void Skinner::add_adjacency( EntityHandle entity, const EntityHandle* nodes, const int num_nodes )
153 : : {
154 : 0 : std::vector< EntityHandle >* adj = NULL;
155 [ # # ]: 0 : const EntityHandle* iter = std::min_element( nodes, nodes + num_nodes );
156 : :
157 [ # # ]: 0 : if( iter == nodes + num_nodes ) return;
158 : :
159 : : // should not be setting adjacency lists in ho-nodes
160 [ # # ][ # # ]: 0 : assert( TYPE_FROM_HANDLE( entity ) == MBPOLYGON ||
[ # # ][ # # ]
161 [ # # ]: 0 : num_nodes == CN::VerticesPerEntity( TYPE_FROM_HANDLE( entity ) ) );
162 : :
163 : : // add this entity to the node
164 [ # # ][ # # ]: 0 : if( thisMB->tag_get_data( mAdjTag, iter, 1, &adj ) == MB_SUCCESS && adj != NULL ) { adj->push_back( entity ); }
[ # # ][ # # ]
[ # # ]
165 : : // create a new vector and add it
166 : : else
167 : : {
168 [ # # ][ # # ]: 0 : adj = new std::vector< EntityHandle >;
169 [ # # ]: 0 : adj->push_back( entity );
170 [ # # ]: 0 : thisMB->tag_set_data( mAdjTag, iter, 1, &adj );
171 : : }
172 : : }
173 : :
174 : 0 : ErrorCode Skinner::find_geometric_skin( const EntityHandle meshset, Range& forward_target_entities )
175 : : {
176 : : // attempts to find whole model skin, using geom topo sets first then
177 : : // normal find_skin function
178 : 0 : bool debug = true;
179 : :
180 : : // look for geom topo sets
181 : : Tag geom_tag;
182 : : ErrorCode result =
183 [ # # ]: 0 : thisMB->tag_get_handle( GEOM_DIMENSION_TAG_NAME, 1, MB_TYPE_INTEGER, geom_tag, MB_TAG_SPARSE | MB_TAG_CREAT );
184 : :
185 [ # # ]: 0 : if( MB_SUCCESS != result ) return result;
186 : :
187 : : // get face sets (dimension = 2)
188 [ # # ]: 0 : Range face_sets;
189 : 0 : int two = 2;
190 : 0 : const void* two_ptr = &two;
191 [ # # ]: 0 : result = thisMB->get_entities_by_type_and_tag( meshset, MBENTITYSET, &geom_tag, &two_ptr, 1, face_sets );
192 : :
193 [ # # ]: 0 : Range::iterator it;
194 [ # # ]: 0 : if( MB_SUCCESS != result )
195 : 0 : return result;
196 [ # # ][ # # ]: 0 : else if( face_sets.empty() )
197 : 0 : return MB_ENTITY_NOT_FOUND;
198 : :
199 : : // ok, we have face sets; use those to determine skin
200 [ # # ]: 0 : Range skin_sets;
201 [ # # ][ # # ]: 0 : if( debug ) std::cout << "Found " << face_sets.size() << " face sets total..." << std::endl;
[ # # ][ # # ]
[ # # ][ # # ]
202 : :
203 [ # # ][ # # ]: 0 : for( it = face_sets.begin(); it != face_sets.end(); ++it )
[ # # ][ # # ]
[ # # ]
204 : : {
205 : : int num_parents;
206 [ # # ][ # # ]: 0 : result = thisMB->num_parent_meshsets( *it, &num_parents );
207 [ # # ]: 0 : if( MB_SUCCESS != result )
208 : 0 : return result;
209 [ # # ]: 0 : else if( num_parents == 1 )
210 [ # # ][ # # ]: 0 : skin_sets.insert( *it );
211 : : }
212 : :
213 [ # # ][ # # ]: 0 : if( debug ) std::cout << "Found " << skin_sets.size() << " 1-parent face sets..." << std::endl;
[ # # ][ # # ]
[ # # ][ # # ]
214 : :
215 [ # # ][ # # ]: 0 : if( skin_sets.empty() ) return MB_FAILURE;
216 : :
217 : : // ok, we have the shell; gather up the elements, putting them all in forward for now
218 [ # # ][ # # ]: 0 : for( it = skin_sets.begin(); it != skin_sets.end(); ++it )
[ # # ][ # # ]
[ # # ]
219 : : {
220 [ # # ][ # # ]: 0 : result = thisMB->get_entities_by_handle( *it, forward_target_entities, true );
221 [ # # ]: 0 : if( MB_SUCCESS != result ) return result;
222 : : }
223 : :
224 : 0 : return result;
225 : : }
226 : :
227 : 72 : ErrorCode Skinner::find_skin( const EntityHandle meshset, const Range& source_entities, bool get_vertices,
228 : : Range& output_handles, Range* output_reverse_handles, bool create_vert_elem_adjs,
229 : : bool create_skin_elements, bool look_for_scd )
230 : : {
231 [ - + ]: 72 : if( source_entities.empty() ) return MB_SUCCESS;
232 : :
233 [ + + ]: 72 : if( look_for_scd )
234 : : {
235 : 2 : ErrorCode rval = find_skin_scd( source_entities, get_vertices, output_handles, create_skin_elements );
236 : : // if it returns success, it's all scd, and we don't need to do anything more
237 [ + - ]: 2 : if( MB_SUCCESS == rval ) return rval;
238 : : }
239 : :
240 [ - + ]: 70 : Core* this_core = dynamic_cast< Core* >( thisMB );
241 [ + - ][ + + ]: 70 : if( this_core && create_vert_elem_adjs && !this_core->a_entity_factory()->vert_elem_adjacencies() )
[ + + ][ + + ]
242 : 11 : this_core->a_entity_factory()->create_vert_elem_adjacencies();
243 : :
244 : : return find_skin_vertices( meshset, source_entities, get_vertices ? &output_handles : 0,
245 [ + + ][ + + ]: 70 : get_vertices ? 0 : &output_handles, output_reverse_handles, create_skin_elements );
246 : : }
247 : :
248 : 2 : ErrorCode Skinner::find_skin_scd( const Range& source_entities, bool get_vertices, Range& output_handles,
249 : : bool create_skin_elements )
250 : : {
251 : : // get the scd interface and check if it's been initialized
252 : 2 : ScdInterface* scdi = NULL;
253 [ + - ]: 2 : ErrorCode rval = thisMB->query_interface( scdi );
254 [ - + ]: 2 : if( !scdi ) return MB_FAILURE;
255 : :
256 : : // ok, there's scd mesh; see if the entities passed in are all in a scd box
257 : : // a box needs to be wholly included in entities for this to work
258 [ + - ][ + - ]: 4 : std::vector< ScdBox* > boxes, myboxes;
259 [ + - ]: 4 : Range myrange;
260 [ + - ]: 2 : rval = scdi->find_boxes( boxes );
261 [ - + ]: 2 : if( MB_SUCCESS != rval ) return rval;
262 [ + - ][ + - ]: 4 : for( std::vector< ScdBox* >::iterator bit = boxes.begin(); bit != boxes.end(); ++bit )
[ + + ]
263 : : {
264 [ + - ][ + - ]: 2 : Range belems( ( *bit )->start_element(), ( *bit )->start_element() + ( *bit )->num_elements() - 1 );
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ]
265 [ + - ][ + - ]: 2 : if( source_entities.contains( belems ) )
266 : : {
267 [ + - ][ + - ]: 2 : myboxes.push_back( *bit );
268 [ + - ]: 2 : myrange.merge( belems );
269 : : }
270 : 2 : }
271 [ + - ][ + - ]: 2 : if( myboxes.empty() || myrange.size() != source_entities.size() ) return MB_FAILURE;
[ + - ][ - + ]
[ - + ]
272 : :
273 : : // ok, we're all structured; get the skin for each box
274 [ + - ][ + - ]: 4 : for( std::vector< ScdBox* >::iterator bit = boxes.begin(); bit != boxes.end(); ++bit )
[ + + ]
275 : : {
276 [ + - ][ + - ]: 2 : rval = skin_box( *bit, get_vertices, output_handles, create_skin_elements );
277 [ - + ]: 2 : if( MB_SUCCESS != rval ) return rval;
278 : : }
279 : :
280 : 4 : return MB_SUCCESS;
281 : : }
282 : :
283 : 2 : ErrorCode Skinner::skin_box( ScdBox* box, bool get_vertices, Range& output_handles, bool create_skin_elements )
284 : : {
285 [ + - ][ + - ]: 2 : HomCoord bmin = box->box_min(), bmax = box->box_max();
286 : :
287 : : // don't support 1d boxes
288 [ + - ][ + - ]: 2 : if( bmin.j() == bmax.j() && bmin.k() == bmax.k() ) return MB_FAILURE;
[ - + ][ # # ]
[ # # ][ # # ]
[ - + ]
289 : :
290 [ + - ][ + - ]: 2 : int dim = ( bmin.k() == bmax.k() ? 1 : 2 );
[ - + ]
291 : :
292 : : ErrorCode rval;
293 : : EntityHandle ent;
294 : :
295 : : // i=min
296 [ + - ][ + - ]: 12 : for( int k = bmin.k(); k < bmax.k(); k++ )
[ + + ]
297 : : {
298 [ + - ][ + - ]: 110 : for( int j = bmin.j(); j < bmax.j(); j++ )
[ + + ]
299 : : {
300 : 100 : ent = 0;
301 [ + - ][ + - ]: 100 : rval = box->get_adj_edge_or_face( dim, bmin.i(), j, k, 0, ent, create_skin_elements );
302 [ - + ]: 100 : if( MB_SUCCESS != rval ) return rval;
303 [ + - ][ + - ]: 100 : if( ent ) output_handles.insert( ent );
304 : : }
305 : : }
306 : : // i=max
307 [ + - ][ + - ]: 12 : for( int k = bmin.k(); k < bmax.k(); k++ )
[ + + ]
308 : : {
309 [ + - ][ + - ]: 110 : for( int j = bmin.j(); j < bmax.j(); j++ )
[ + + ]
310 : : {
311 : 100 : ent = 0;
312 [ + - ][ + - ]: 100 : rval = box->get_adj_edge_or_face( dim, bmax.i(), j, k, 0, ent, create_skin_elements );
313 [ - + ]: 100 : if( MB_SUCCESS != rval ) return rval;
314 [ + - ][ + - ]: 100 : if( ent ) output_handles.insert( ent );
315 : : }
316 : : }
317 : : // j=min
318 [ + - ][ + - ]: 12 : for( int k = bmin.k(); k < bmax.k(); k++ )
[ + + ]
319 : : {
320 [ + - ][ + - ]: 110 : for( int i = bmin.i(); i < bmax.i(); i++ )
[ + + ]
321 : : {
322 : 100 : ent = 0;
323 [ + - ][ + - ]: 100 : rval = box->get_adj_edge_or_face( dim, i, bmin.j(), k, 1, ent, create_skin_elements );
324 [ - + ]: 100 : if( MB_SUCCESS != rval ) return rval;
325 [ + - ][ + - ]: 100 : if( ent ) output_handles.insert( ent );
326 : : }
327 : : }
328 : : // j=max
329 [ + - ][ + - ]: 12 : for( int k = bmin.k(); k < bmax.k(); k++ )
[ + + ]
330 : : {
331 [ + - ][ + - ]: 110 : for( int i = bmin.i(); i < bmax.i(); i++ )
[ + + ]
332 : : {
333 : 100 : ent = 0;
334 [ + - ][ + - ]: 100 : rval = box->get_adj_edge_or_face( dim, i, bmax.j(), k, 1, ent, create_skin_elements );
335 [ - + ]: 100 : if( MB_SUCCESS != rval ) return rval;
336 [ + - ][ + - ]: 100 : if( ent ) output_handles.insert( ent );
337 : : }
338 : : }
339 : : // k=min
340 [ + - ][ + - ]: 22 : for( int j = bmin.j(); j < bmax.j(); j++ )
[ + + ]
341 : : {
342 [ + - ][ + - ]: 220 : for( int i = bmin.i(); i < bmax.i(); i++ )
[ + + ]
343 : : {
344 : 200 : ent = 0;
345 [ + - ][ + - ]: 200 : rval = box->get_adj_edge_or_face( dim, i, j, bmin.k(), 2, ent, create_skin_elements );
346 [ - + ]: 200 : if( MB_SUCCESS != rval ) return rval;
347 [ + - ][ + - ]: 200 : if( ent ) output_handles.insert( ent );
348 : : }
349 : : }
350 : : // k=max
351 [ + - ][ + - ]: 22 : for( int j = bmin.j(); j < bmax.j(); j++ )
[ + + ]
352 : : {
353 [ + - ][ + - ]: 220 : for( int i = bmin.i(); i < bmax.i(); i++ )
[ + + ]
354 : : {
355 : 200 : ent = 0;
356 [ + - ][ + - ]: 200 : rval = box->get_adj_edge_or_face( dim, i, j, bmax.k(), 2, ent, create_skin_elements );
357 [ - + ]: 200 : if( MB_SUCCESS != rval ) return rval;
358 [ + - ][ + - ]: 200 : if( ent ) output_handles.insert( ent );
359 : : }
360 : : }
361 : :
362 [ + + ]: 2 : if( get_vertices )
363 : : {
364 [ + - ]: 1 : Range verts;
365 [ + - ]: 1 : rval = thisMB->get_adjacencies( output_handles, 0, true, verts, Interface::UNION );
366 [ - + ]: 1 : if( MB_SUCCESS != rval ) return rval;
367 [ + - ][ + - ]: 1 : output_handles.merge( verts );
368 : : }
369 : :
370 : 2 : return MB_SUCCESS;
371 : : }
372 : :
373 : 0 : ErrorCode Skinner::find_skin_noadj(const Range &source_entities,
374 : : Range &forward_target_entities,
375 : : Range &reverse_target_entities/*,
376 : : bool create_vert_elem_adjs*/)
377 : : {
378 [ # # ][ # # ]: 0 : if( source_entities.empty() ) return MB_FAILURE;
379 : :
380 : : // get our working dimensions
381 [ # # ][ # # ]: 0 : EntityType type = thisMB->type_from_handle( *( source_entities.begin() ) );
[ # # ]
382 [ # # ]: 0 : const int source_dim = CN::Dimension( type );
383 : 0 : mTargetDim = source_dim - 1;
384 : :
385 : : // make sure we can handle the working dimensions
386 [ # # ][ # # ]: 0 : if( mTargetDim < 0 || source_dim > 3 ) return MB_FAILURE;
387 : :
388 [ # # ]: 0 : initialize();
389 : :
390 [ # # ][ # # ]: 0 : Range::const_iterator iter, end_iter;
391 [ # # ]: 0 : end_iter = source_entities.end();
392 : : const EntityHandle* conn;
393 : : EntityHandle match;
394 : :
395 : : direction direct;
396 : : ErrorCode result;
397 : : // assume we'll never have more than 32 vertices on a facet (checked
398 : : // with assert later)
399 : : EntityHandle sub_conn[32];
400 [ # # ]: 0 : std::vector< EntityHandle > tmp_conn_vec;
401 : : int num_nodes, num_sub_nodes, num_sides;
402 : : int sub_indices[32]; // Also, assume that no polygon has more than 32 nodes
403 : : // we could increase that, but we will not display it right in visit moab h5m , anyway
404 : : EntityType sub_type;
405 : :
406 : : // for each source entity
407 [ # # ][ # # ]: 0 : for( iter = source_entities.begin(); iter != end_iter; ++iter )
[ # # ][ # # ]
408 : : {
409 : : // get the connectivity of this entity
410 : 0 : int actual_num_nodes_polygon = 0;
411 [ # # ][ # # ]: 0 : result = thisMB->get_connectivity( *iter, conn, num_nodes, false, &tmp_conn_vec );
412 [ # # ]: 0 : if( MB_SUCCESS != result ) return result;
413 : :
414 [ # # ][ # # ]: 0 : type = thisMB->type_from_handle( *iter );
415 [ # # ]: 0 : Range::iterator seek_iter;
416 : :
417 : : // treat separately polygons (also, polyhedra will need special handling)
418 [ # # ]: 0 : if( MBPOLYGON == type )
419 : : {
420 : : // treat padded polygons, if existing; count backwards, see how many of the last nodes
421 : : // are repeated assume connectivity is fine, otherwise we could be in trouble
422 : 0 : actual_num_nodes_polygon = num_nodes;
423 [ # # ][ # # ]: 0 : while( actual_num_nodes_polygon >= 3 &&
424 : 0 : conn[actual_num_nodes_polygon - 1] == conn[actual_num_nodes_polygon - 2] )
425 : 0 : actual_num_nodes_polygon--;
426 : 0 : num_sides = actual_num_nodes_polygon;
427 : 0 : sub_type = MBEDGE;
428 : 0 : num_sub_nodes = 2;
429 : : }
430 : : else // get connectivity of each n-1 dimension entity
431 [ # # ]: 0 : num_sides = CN::NumSubEntities( type, mTargetDim );
432 [ # # ]: 0 : for( int i = 0; i < num_sides; i++ )
433 : : {
434 [ # # ]: 0 : if( MBPOLYGON == type )
435 : : {
436 : 0 : sub_conn[0] = conn[i];
437 : 0 : sub_conn[1] = conn[i + 1];
438 [ # # ]: 0 : if( i + 1 == actual_num_nodes_polygon ) sub_conn[1] = conn[0];
439 : : }
440 : : else
441 : : {
442 [ # # ]: 0 : CN::SubEntityNodeIndices( type, num_nodes, mTargetDim, i, sub_type, num_sub_nodes, sub_indices );
443 [ # # ]: 0 : assert( (size_t)num_sub_nodes <= sizeof( sub_indices ) / sizeof( sub_indices[0] ) );
444 [ # # ]: 0 : for( int j = 0; j < num_sub_nodes; j++ )
445 : 0 : sub_conn[j] = conn[sub_indices[j]];
446 : : }
447 : :
448 : : // see if we can match this connectivity with
449 : : // an existing entity
450 [ # # ]: 0 : find_match( sub_type, sub_conn, num_sub_nodes, match, direct );
451 : :
452 : : // if there is no match, create a new entity
453 [ # # ]: 0 : if( match == 0 )
454 : : {
455 : 0 : EntityHandle tmphndl = 0;
456 : : int indices[MAX_SUB_ENTITY_VERTICES];
457 : : EntityType new_type;
458 : : int num_new_nodes;
459 [ # # ]: 0 : if( MBPOLYGON == type )
460 : : {
461 : 0 : new_type = MBEDGE;
462 : 0 : num_new_nodes = 2;
463 : : }
464 : : else
465 : : {
466 [ # # ]: 0 : CN::SubEntityNodeIndices( type, num_nodes, mTargetDim, i, new_type, num_new_nodes, indices );
467 [ # # ]: 0 : for( int j = 0; j < num_new_nodes; j++ )
468 : 0 : sub_conn[j] = conn[indices[j]];
469 : : }
470 [ # # ]: 0 : result = thisMB->create_element( new_type, sub_conn, num_new_nodes, tmphndl );
471 [ # # ]: 0 : assert( MB_SUCCESS == result );
472 [ # # ][ # # ]: 0 : add_adjacency( tmphndl, sub_conn, CN::VerticesPerEntity( new_type ) );
473 [ # # ]: 0 : forward_target_entities.insert( tmphndl );
474 : : }
475 : : // if there is a match, delete the matching entity
476 : : // if we can.
477 : : else
478 : : {
479 [ # # ][ # # ]: 0 : if( ( seek_iter = forward_target_entities.find( match ) ) != forward_target_entities.end() )
[ # # ][ # # ]
480 : : {
481 [ # # ]: 0 : forward_target_entities.erase( seek_iter );
482 [ # # ]: 0 : remove_adjacency( match );
483 [ # # ][ # # ]: 0 : if( /*!use_adjs &&*/ entity_deletable( match ) )
484 : : {
485 [ # # ]: 0 : result = thisMB->delete_entities( &match, 1 );
486 [ # # ]: 0 : assert( MB_SUCCESS == result );
487 : : }
488 : : }
489 [ # # ][ # # ]: 0 : else if( ( seek_iter = reverse_target_entities.find( match ) ) != reverse_target_entities.end() )
[ # # ][ # # ]
490 : : {
491 [ # # ]: 0 : reverse_target_entities.erase( seek_iter );
492 [ # # ]: 0 : remove_adjacency( match );
493 [ # # ][ # # ]: 0 : if( /*!use_adjs &&*/ entity_deletable( match ) )
494 : : {
495 [ # # ]: 0 : result = thisMB->delete_entities( &match, 1 );
496 [ # # ]: 0 : assert( MB_SUCCESS == result );
497 : : }
498 : : }
499 : : else
500 : : {
501 [ # # ][ # # ]: 0 : if( direct == FORWARD ) { forward_target_entities.insert( match ); }
502 : : else
503 : : {
504 [ # # ]: 0 : reverse_target_entities.insert( match );
505 : : }
506 : : }
507 : : }
508 : : }
509 : : }
510 : :
511 [ # # ]: 0 : deinitialize();
512 : :
513 : 0 : return MB_SUCCESS;
514 : : }
515 : :
516 : 0 : void Skinner::find_match( EntityType type, const EntityHandle* conn, const int num_nodes, EntityHandle& match,
517 : : Skinner::direction& direct )
518 : : {
519 : 0 : match = 0;
520 : :
521 [ # # ]: 0 : if( type == MBVERTEX )
522 : : {
523 : 0 : match = *conn;
524 : 0 : direct = FORWARD;
525 : 0 : return;
526 : : }
527 : :
528 [ # # ]: 0 : const EntityHandle* iter = std::min_element( conn, conn + num_nodes );
529 : :
530 : 0 : std::vector< EntityHandle >* adj = NULL;
531 : :
532 [ # # ]: 0 : ErrorCode result = thisMB->tag_get_data( mAdjTag, iter, 1, &adj );
533 [ # # ][ # # ]: 0 : if( result == MB_FAILURE || adj == NULL ) { return; }
534 : :
535 : 0 : std::vector< EntityHandle >::iterator jter, end_jter;
536 : 0 : end_jter = adj->end();
537 : :
538 : : const EntityHandle* tmp;
539 : : int num_verts;
540 : :
541 [ # # ][ # # ]: 0 : for( jter = adj->begin(); jter != end_jter; ++jter )
[ # # ]
542 : : {
543 : : EntityType tmp_type;
544 [ # # ][ # # ]: 0 : tmp_type = thisMB->type_from_handle( *jter );
545 : :
546 [ # # ]: 0 : if( type != tmp_type ) continue;
547 : :
548 [ # # ][ # # ]: 0 : result = thisMB->get_connectivity( *jter, tmp, num_verts, false );
549 [ # # ][ # # ]: 0 : assert( MB_SUCCESS == result && num_verts >= CN::VerticesPerEntity( type ) );
[ # # ]
550 : : // FIXME: connectivity_match appears to work only for linear elements,
551 : : // so ignore higher-order nodes.
552 [ # # ][ # # ]: 0 : if( connectivity_match( conn, tmp, CN::VerticesPerEntity( type ), direct ) )
[ # # ]
553 : : {
554 [ # # ]: 0 : match = *jter;
555 : 0 : break;
556 : : }
557 : : }
558 : : }
559 : :
560 : 0 : bool Skinner::connectivity_match( const EntityHandle* conn1, const EntityHandle* conn2, const int num_verts,
561 : : Skinner::direction& direct )
562 : : {
563 : 0 : const EntityHandle* iter = std::find( conn2, conn2 + num_verts, conn1[0] );
564 [ # # ]: 0 : if( iter == conn2 + num_verts ) return false;
565 : :
566 : 0 : bool they_match = true;
567 : :
568 : : int i;
569 : 0 : unsigned int j = iter - conn2;
570 : :
571 : : // first compare forward
572 [ # # ]: 0 : for( i = 1; i < num_verts; ++i )
573 : : {
574 [ # # ]: 0 : if( conn1[i] != conn2[( j + i ) % num_verts] )
575 : : {
576 : 0 : they_match = false;
577 : 0 : break;
578 : : }
579 : : }
580 : :
581 [ # # ]: 0 : if( they_match == true )
582 : : {
583 : : // need to check for reversed edges here
584 [ # # ][ # # ]: 0 : direct = ( num_verts == 2 && j ) ? REVERSE : FORWARD;
585 : 0 : return true;
586 : : }
587 : :
588 : 0 : they_match = true;
589 : :
590 : : // then compare reverse
591 : 0 : j += num_verts;
592 [ # # ]: 0 : for( i = 1; i < num_verts; )
593 : : {
594 [ # # ]: 0 : if( conn1[i] != conn2[( j - i ) % num_verts] )
595 : : {
596 : 0 : they_match = false;
597 : 0 : break;
598 : : }
599 : 0 : ++i;
600 : : }
601 [ # # ]: 0 : if( they_match ) { direct = REVERSE; }
602 : 0 : return they_match;
603 : : }
604 : :
605 : 0 : ErrorCode Skinner::remove_adjacency( EntityHandle entity )
606 : : {
607 [ # # ]: 0 : std::vector< EntityHandle > nodes, *adj = NULL;
608 [ # # ][ # # ]: 0 : ErrorCode result = thisMB->get_connectivity( &entity, 1, nodes );MB_CHK_ERR( result );
[ # # ][ # # ]
609 [ # # ]: 0 : std::vector< EntityHandle >::iterator iter = std::min_element( nodes.begin(), nodes.end() );
610 : :
611 [ # # ][ # # ]: 0 : if( iter == nodes.end() ) return MB_FAILURE;
612 : :
613 : : // remove this entity from the node
614 [ # # ][ # # ]: 0 : if( thisMB->tag_get_data( mAdjTag, &( *iter ), 1, &adj ) == MB_SUCCESS && adj != NULL )
[ # # ][ # # ]
[ # # ]
615 : : {
616 [ # # ]: 0 : iter = std::find( adj->begin(), adj->end(), entity );
617 [ # # ][ # # ]: 0 : if( iter != adj->end() ) adj->erase( iter );
[ # # ]
618 : : }
619 : :
620 : 0 : return result;
621 : : }
622 : :
623 : 0 : bool Skinner::entity_deletable( EntityHandle entity )
624 : : {
625 : 0 : unsigned char deletable = 0;
626 [ # # ]: 0 : ErrorCode result = thisMB->tag_get_data( mDeletableMBTag, &entity, 1, &deletable );
627 [ # # ]: 0 : assert( MB_SUCCESS == result );
628 [ # # ][ # # ]: 0 : if( MB_SUCCESS == result && deletable == 1 ) return false;
629 : 0 : return true;
630 : : }
631 : :
632 : 0 : ErrorCode Skinner::classify_2d_boundary( const Range& boundary, const Range& bar_elements, EntityHandle boundary_edges,
633 : : EntityHandle inferred_edges, EntityHandle non_manifold_edges,
634 : : EntityHandle other_edges, int& number_boundary_nodes )
635 : : {
636 [ # # ][ # # ]: 0 : Range bedges, iedges, nmedges, oedges;
[ # # ][ # # ]
637 : : ErrorCode result =
638 [ # # ][ # # ]: 0 : classify_2d_boundary( boundary, bar_elements, bedges, iedges, nmedges, oedges, number_boundary_nodes );MB_CHK_ERR( result );
[ # # ][ # # ]
639 : :
640 : : // now set the input meshsets to the output ranges
641 [ # # ][ # # ]: 0 : result = thisMB->clear_meshset( &boundary_edges, 1 );MB_CHK_ERR( result );
[ # # ][ # # ]
642 [ # # ][ # # ]: 0 : result = thisMB->add_entities( boundary_edges, bedges );MB_CHK_ERR( result );
[ # # ][ # # ]
643 : :
644 [ # # ][ # # ]: 0 : result = thisMB->clear_meshset( &inferred_edges, 1 );MB_CHK_ERR( result );
[ # # ][ # # ]
645 [ # # ][ # # ]: 0 : result = thisMB->add_entities( inferred_edges, iedges );MB_CHK_ERR( result );
[ # # ][ # # ]
646 : :
647 [ # # ][ # # ]: 0 : result = thisMB->clear_meshset( &non_manifold_edges, 1 );MB_CHK_ERR( result );
[ # # ][ # # ]
648 [ # # ][ # # ]: 0 : result = thisMB->add_entities( non_manifold_edges, nmedges );MB_CHK_ERR( result );
[ # # ][ # # ]
649 : :
650 [ # # ][ # # ]: 0 : result = thisMB->clear_meshset( &other_edges, 1 );MB_CHK_ERR( result );
[ # # ][ # # ]
651 [ # # ][ # # ]: 0 : result = thisMB->add_entities( other_edges, oedges );MB_CHK_ERR( result );
[ # # ][ # # ]
652 : :
653 : 0 : return MB_SUCCESS;
654 : : }
655 : :
656 : 0 : ErrorCode Skinner::classify_2d_boundary( const Range& boundary, const Range& bar_elements, Range& boundary_edges,
657 : : Range& inferred_edges, Range& non_manifold_edges, Range& other_edges,
658 : : int& number_boundary_nodes )
659 : : {
660 : :
661 : : // clear out the edge lists
662 : :
663 [ # # ]: 0 : boundary_edges.clear();
664 [ # # ]: 0 : inferred_edges.clear();
665 [ # # ]: 0 : non_manifold_edges.clear();
666 [ # # ]: 0 : other_edges.clear();
667 : :
668 : 0 : number_boundary_nodes = 0;
669 : :
670 : : // make sure we have something to work with
671 [ # # ][ # # ]: 0 : if( boundary.empty() ) { return MB_FAILURE; }
672 : :
673 : : // get our working dimensions
674 [ # # ][ # # ]: 0 : EntityType type = thisMB->type_from_handle( *( boundary.begin() ) );
[ # # ]
675 [ # # ]: 0 : const int source_dim = CN::Dimension( type );
676 : :
677 : : // make sure we can handle the working dimensions
678 [ # # ]: 0 : if( source_dim != 2 ) { return MB_FAILURE; }
679 : 0 : mTargetDim = source_dim - 1;
680 : :
681 : : // initialize
682 [ # # ]: 0 : initialize();
683 : :
684 : : // additional initialization for this routine
685 : : // define a tag for MBEDGE which counts the occurrences of the edge below
686 : : // default should be 0 for existing edges, if any
687 : :
688 : : Tag count_tag;
689 : 0 : int default_count = 0;
690 : : ErrorCode result =
691 [ # # ][ # # ]: 0 : thisMB->tag_get_handle( 0, 1, MB_TYPE_INTEGER, count_tag, MB_TAG_DENSE | MB_TAG_CREAT, &default_count );MB_CHK_ERR( result );
[ # # ][ # # ]
692 : :
693 [ # # ][ # # ]: 0 : Range::const_iterator iter, end_iter;
694 [ # # ]: 0 : end_iter = boundary.end();
695 : :
696 [ # # ]: 0 : std::vector< EntityHandle > conn;
697 : : EntityHandle sub_conn[2];
698 : : EntityHandle match;
699 : :
700 [ # # ]: 0 : Range edge_list;
701 [ # # ]: 0 : Range boundary_nodes;
702 : : Skinner::direction direct;
703 : :
704 : : EntityType sub_type;
705 : : int num_edge, num_sub_ent_vert;
706 : : const short* edge_verts;
707 : :
708 : : // now, process each entity in the boundary
709 : :
710 [ # # ][ # # ]: 0 : for( iter = boundary.begin(); iter != end_iter; ++iter )
[ # # ][ # # ]
711 : : {
712 : : // get the connectivity of this entity
713 : 0 : conn.clear();
714 [ # # ][ # # ]: 0 : result = thisMB->get_connectivity( &( *iter ), 1, conn, false );
715 [ # # ]: 0 : assert( MB_SUCCESS == result );
716 : :
717 : : // add node handles to boundary_node range
718 [ # # ][ # # ]: 0 : std::copy( conn.begin(), conn.begin() + CN::VerticesPerEntity( type ), range_inserter( boundary_nodes ) );
[ # # ][ # # ]
719 : :
720 [ # # ][ # # ]: 0 : type = thisMB->type_from_handle( *iter );
721 : :
722 : : // get connectivity of each n-1 dimension entity (edge in this case)
723 : 0 : const struct CN::ConnMap* conn_map = &( CN::mConnectivityMap[type][0] );
724 [ # # ]: 0 : num_edge = CN::NumSubEntities( type, 1 );
725 [ # # ]: 0 : for( int i = 0; i < num_edge; i++ )
726 : : {
727 [ # # ]: 0 : edge_verts = CN::SubEntityVertexIndices( type, 1, i, sub_type, num_sub_ent_vert );
728 [ # # ][ # # ]: 0 : assert( sub_type == MBEDGE && num_sub_ent_vert == 2 );
729 [ # # ]: 0 : sub_conn[0] = conn[edge_verts[0]];
730 [ # # ]: 0 : sub_conn[1] = conn[edge_verts[1]];
731 : 0 : int num_sub_nodes = conn_map->num_corners_per_sub_element[i];
732 : :
733 : : // see if we can match this connectivity with
734 : : // an existing entity
735 [ # # ]: 0 : find_match( MBEDGE, sub_conn, num_sub_nodes, match, direct );
736 : :
737 : : // if there is no match, create a new entity
738 [ # # ]: 0 : if( match == 0 )
739 : : {
740 : 0 : EntityHandle tmphndl = 0;
741 : : int indices[MAX_SUB_ENTITY_VERTICES];
742 : : EntityType new_type;
743 : : int num_new_nodes;
744 [ # # ]: 0 : CN::SubEntityNodeIndices( type, conn.size(), 1, i, new_type, num_new_nodes, indices );
745 [ # # ]: 0 : for( int j = 0; j < num_new_nodes; j++ )
746 [ # # ]: 0 : sub_conn[j] = conn[indices[j]];
747 : :
748 [ # # ]: 0 : result = thisMB->create_element( new_type, sub_conn, num_new_nodes, tmphndl );
749 [ # # ]: 0 : assert( MB_SUCCESS == result );
750 [ # # ]: 0 : add_adjacency( tmphndl, sub_conn, num_sub_nodes );
751 : : // target_entities.insert(tmphndl);
752 [ # # ]: 0 : edge_list.insert( tmphndl );
753 : : int count;
754 [ # # ]: 0 : result = thisMB->tag_get_data( count_tag, &tmphndl, 1, &count );
755 [ # # ]: 0 : assert( MB_SUCCESS == result );
756 : 0 : count++;
757 [ # # ]: 0 : result = thisMB->tag_set_data( count_tag, &tmphndl, 1, &count );
758 [ # # ]: 0 : assert( MB_SUCCESS == result );
759 : : }
760 : : else
761 : : {
762 : : // We found a match, we must increment the count on the match
763 : : int count;
764 [ # # ]: 0 : result = thisMB->tag_get_data( count_tag, &match, 1, &count );
765 [ # # ]: 0 : assert( MB_SUCCESS == result );
766 : 0 : count++;
767 [ # # ]: 0 : result = thisMB->tag_set_data( count_tag, &match, 1, &count );
768 [ # # ]: 0 : assert( MB_SUCCESS == result );
769 : :
770 : : // if the entity is not deletable, it was pre-existing in
771 : : // the database. We therefore may need to add it to the
772 : : // edge_list. Since it will not hurt the range, we add
773 : : // whether it was added before or not
774 [ # # ][ # # ]: 0 : if( !entity_deletable( match ) ) { edge_list.insert( match ); }
[ # # ]
775 : : }
776 : : }
777 : : }
778 : :
779 : : // Any bar elements in the model should be classified separately
780 : : // If the element is in the skin edge_list, then it should be put in
781 : : // the non-manifold edge list. Edges not in the edge_list are stand-alone
782 : : // bars, and we make them simply boundary elements
783 : :
784 [ # # ][ # # ]: 0 : if( !bar_elements.empty() )
785 : : {
786 [ # # ]: 0 : Range::iterator bar_iter;
787 [ # # ][ # # ]: 0 : for( iter = bar_elements.begin(); iter != bar_elements.end(); ++iter )
[ # # ][ # # ]
[ # # ]
788 : : {
789 [ # # ]: 0 : EntityHandle handle = *iter;
790 [ # # ]: 0 : bar_iter = edge_list.find( handle );
791 [ # # ][ # # ]: 0 : if( bar_iter != edge_list.end() )
[ # # ]
792 : : {
793 : : // it is in the list, erase it and put in non-manifold list
794 [ # # ]: 0 : edge_list.erase( bar_iter );
795 [ # # ]: 0 : non_manifold_edges.insert( handle );
796 : : }
797 : : else
798 : : {
799 : : // not in the edge list, make it a boundary edge
800 [ # # ]: 0 : boundary_edges.insert( handle );
801 : : }
802 : : }
803 : : }
804 : :
805 : : // now all edges should be classified. Go through the edge_list,
806 : : // and put all in the appropriate lists
807 : :
808 [ # # ][ # # ]: 0 : Range::iterator edge_iter, edge_end_iter;
809 [ # # ]: 0 : edge_end_iter = edge_list.end();
810 : : int count;
811 [ # # ][ # # ]: 0 : for( edge_iter = edge_list.begin(); edge_iter != edge_end_iter; ++edge_iter )
[ # # ][ # # ]
812 : : {
813 : : // check the count_tag
814 [ # # ][ # # ]: 0 : result = thisMB->tag_get_data( count_tag, &( *edge_iter ), 1, &count );
815 [ # # ]: 0 : assert( MB_SUCCESS == result );
816 [ # # ][ # # ]: 0 : if( count == 1 ) { boundary_edges.insert( *edge_iter ); }
[ # # ]
817 [ # # ]: 0 : else if( count == 2 )
818 : : {
819 [ # # ][ # # ]: 0 : other_edges.insert( *edge_iter );
820 : : }
821 : : else
822 : : {
823 [ # # ][ # # ]: 0 : non_manifold_edges.insert( *edge_iter );
824 : : }
825 : : }
826 : :
827 : : // find the inferred edges from the other_edge_list
828 : :
829 : 0 : double min_angle_degrees = 20.0;
830 [ # # ]: 0 : find_inferred_edges( const_cast< Range& >( boundary ), other_edges, inferred_edges, min_angle_degrees );
831 : :
832 : : // we now want to remove the inferred_edges from the other_edges
833 : :
834 [ # # ]: 0 : Range temp_range;
835 : :
836 : : std::set_difference( other_edges.begin(), other_edges.end(), inferred_edges.begin(), inferred_edges.end(),
837 [ # # ][ # # ]: 0 : range_inserter( temp_range ), std::less< EntityHandle >() );
[ # # ][ # # ]
[ # # ][ # # ]
838 : :
839 [ # # ]: 0 : other_edges = temp_range;
840 : :
841 : : // get rid of count tag and deinitialize
842 : :
843 [ # # ]: 0 : result = thisMB->tag_delete( count_tag );
844 [ # # ]: 0 : assert( MB_SUCCESS == result );
845 [ # # ]: 0 : deinitialize();
846 : :
847 : : // set the node count
848 [ # # ]: 0 : number_boundary_nodes = boundary_nodes.size();
849 : :
850 : 0 : return MB_SUCCESS;
851 : : }
852 : :
853 : 0 : void Skinner::find_inferred_edges( Range& skin_boundary, Range& candidate_edges, Range& inferred_edges,
854 : : double reference_angle_degrees )
855 : : {
856 : :
857 : : // mark all the entities in the skin boundary
858 : : Tag mark_tag;
859 [ # # ]: 0 : ErrorCode result = thisMB->tag_get_handle( 0, 1, MB_TYPE_BIT, mark_tag, MB_TAG_CREAT );
860 [ # # ]: 0 : assert( MB_SUCCESS == result );
861 : 0 : unsigned char bit = true;
862 [ # # ]: 0 : result = thisMB->tag_clear_data( mark_tag, skin_boundary, &bit );
863 [ # # ]: 0 : assert( MB_SUCCESS == result );
864 : :
865 : : // find the cosine of the reference angle
866 : :
867 : 0 : double reference_cosine = cos( reference_angle_degrees * SKINNER_PI / 180.0 );
868 : :
869 : : // check all candidate edges for an angle greater than the minimum
870 : :
871 [ # # ][ # # ]: 0 : Range::iterator iter, end_iter = candidate_edges.end();
872 [ # # ]: 0 : std::vector< EntityHandle > adjacencies;
873 : 0 : std::vector< EntityHandle >::iterator adj_iter;
874 : : EntityHandle face[2];
875 : :
876 [ # # ][ # # ]: 0 : for( iter = candidate_edges.begin(); iter != end_iter; ++iter )
[ # # ][ # # ]
877 : : {
878 : :
879 : : // get the 2D elements connected to this edge
880 : 0 : adjacencies.clear();
881 [ # # ][ # # ]: 0 : result = thisMB->get_adjacencies( &( *iter ), 1, 2, false, adjacencies );
882 [ # # ]: 0 : if( MB_SUCCESS != result ) continue;
883 : :
884 : : // there should be exactly two, that is why the edge is classified as nonBoundary
885 : : // and manifold
886 : :
887 : 0 : int faces_found = 0;
888 [ # # ][ # # ]: 0 : for( adj_iter = adjacencies.begin(); adj_iter != adjacencies.end() && faces_found < 2; ++adj_iter )
[ # # ][ # # ]
[ # # ]
[ # # # # ]
889 : : {
890 : : // we need to find two of these which are in the skin
891 : 0 : unsigned char is_marked = 0;
892 [ # # ][ # # ]: 0 : result = thisMB->tag_get_data( mark_tag, &( *adj_iter ), 1, &is_marked );
893 [ # # ]: 0 : assert( MB_SUCCESS == result );
894 [ # # ]: 0 : if( is_marked )
895 : : {
896 [ # # ]: 0 : face[faces_found] = *adj_iter;
897 : 0 : faces_found++;
898 : : }
899 : : }
900 : :
901 : : // assert(faces_found == 2 || faces_found == 0);
902 [ # # ]: 0 : if( 2 != faces_found ) continue;
903 : :
904 : : // see if the two entities have a sufficient angle
905 : :
906 [ # # ][ # # ]: 0 : if( has_larger_angle( face[0], face[1], reference_cosine ) ) { inferred_edges.insert( *iter ); }
[ # # ][ # # ]
907 : : }
908 : :
909 [ # # ]: 0 : result = thisMB->tag_delete( mark_tag );
910 [ # # ]: 0 : assert( MB_SUCCESS == result );
911 : 0 : }
912 : :
913 : 0 : bool Skinner::has_larger_angle( EntityHandle& entity1, EntityHandle& entity2, double reference_angle_cosine )
914 : : {
915 : : // compare normals to get angle. We assume that the surface quads
916 : : // which we test here will be approximately planar
917 : :
918 : : double norm[2][3];
919 [ # # ]: 0 : Util::normal( thisMB, entity1, norm[0][0], norm[0][1], norm[0][2] );
920 [ # # ]: 0 : Util::normal( thisMB, entity2, norm[1][0], norm[1][1], norm[1][2] );
921 : :
922 : 0 : double cosine = norm[0][0] * norm[1][0] + norm[0][1] * norm[1][1] + norm[0][2] * norm[1][2];
923 : :
924 [ # # ]: 0 : if( cosine < reference_angle_cosine ) { return true; }
925 : :
926 : 0 : return false;
927 : : }
928 : :
929 : : // get skin entities of prescribed dimension
930 : 33 : ErrorCode Skinner::find_skin( const EntityHandle this_set, const Range& entities, int dim, Range& skin_entities,
931 : : bool create_vert_elem_adjs, bool create_skin_elements )
932 : : {
933 [ + - ]: 33 : Range tmp_skin;
934 : : ErrorCode result =
935 [ + - ]: 33 : find_skin( this_set, entities, ( dim == 0 ), tmp_skin, 0, create_vert_elem_adjs, create_skin_elements );
936 [ + - ][ + - ]: 33 : if( MB_SUCCESS != result || tmp_skin.empty() ) return result;
[ - + ][ - + ]
937 : :
938 [ + - ][ + + ]: 33 : if( tmp_skin.all_of_dimension( dim ) )
939 : : {
940 [ + - ][ + - ]: 32 : if( skin_entities.empty() )
941 [ + - ]: 32 : skin_entities.swap( tmp_skin );
942 : : else
943 [ # # ]: 32 : skin_entities.merge( tmp_skin );
944 : : }
945 : : else
946 : : {
947 [ + - ][ - + ]: 1 : result = thisMB->get_adjacencies( tmp_skin, dim, create_skin_elements, skin_entities, Interface::UNION );MB_CHK_ERR( result );
[ # # ][ # # ]
948 [ + - ][ + - ]: 1 : if( this_set ) result = thisMB->add_entities( this_set, skin_entities );
949 : : }
950 : :
951 : 33 : return result;
952 : : }
953 : :
954 : 70 : ErrorCode Skinner::find_skin_vertices( const EntityHandle this_set, const Range& entities, Range* skin_verts,
955 : : Range* skin_elems, Range* skin_rev_elems, bool create_skin_elems,
956 : : bool corners_only )
957 : : {
958 : : ErrorCode rval;
959 [ + - ][ - + ]: 70 : if( entities.empty() ) return MB_SUCCESS;
960 : :
961 [ + - ][ + - ]: 70 : const int dim = CN::Dimension( TYPE_FROM_HANDLE( entities.front() ) );
[ + - ]
962 [ + - ][ + - ]: 70 : if( dim < 1 || dim > 3 || !entities.all_of_dimension( dim ) ) return MB_TYPE_OUT_OF_RANGE;
[ + - ][ - + ]
[ - + ]
963 : :
964 : : // are we skinning all entities
965 [ + - ]: 70 : size_t count = entities.size();
966 : : int num_total;
967 [ + - ]: 70 : rval = thisMB->get_number_entities_by_dimension( this_set, dim, num_total );
968 [ - + ]: 70 : if( MB_SUCCESS != rval ) return rval;
969 : 70 : bool all = ( count == (size_t)num_total );
970 : :
971 : : // Create a bit tag for fast intersection with input entities range.
972 : : // If we're skinning all the entities in the mesh, we really don't
973 : : // need the tag. To save memory, just create it with a default value
974 : : // of one and don't set it. That way MOAB will return 1 for all
975 : : // entities.
976 : : Tag tag;
977 [ + + ]: 70 : char bit = all ? 1 : 0;
978 [ + - ]: 70 : rval = thisMB->tag_get_handle( NULL, 1, MB_TYPE_BIT, tag, MB_TAG_CREAT, &bit );
979 [ - + ]: 70 : if( MB_SUCCESS != rval ) return rval;
980 : :
981 : : // tag all entities in input range
982 [ + + ]: 70 : if( !all )
983 : : {
984 [ + - ]: 30 : std::vector< unsigned char > vect( count, 1 );
985 [ + - ][ + - ]: 30 : rval = thisMB->tag_set_data( tag, entities, &vect[0] );
986 [ - + ]: 30 : if( MB_SUCCESS != rval )
987 : : {
988 [ # # ]: 0 : thisMB->tag_delete( tag );
989 [ + - ]: 30 : return rval;
990 : 30 : }
991 : : }
992 : :
993 [ + + + - ]: 70 : switch( dim )
994 : : {
995 : : case 1:
996 [ + - ]: 4 : if( skin_verts )
997 [ + - ]: 4 : rval = find_skin_vertices_1D( tag, entities, *skin_verts );
998 [ # # ]: 0 : else if( skin_elems )
999 [ # # ]: 0 : rval = find_skin_vertices_1D( tag, entities, *skin_elems );
1000 : : else
1001 : 0 : rval = MB_SUCCESS;
1002 : 4 : break;
1003 : : case 2:
1004 : : rval = find_skin_vertices_2D( this_set, tag, entities, skin_verts, skin_elems, skin_rev_elems,
1005 [ + - ]: 40 : create_skin_elems, corners_only );
1006 : 40 : break;
1007 : : case 3:
1008 : : rval = find_skin_vertices_3D( this_set, tag, entities, skin_verts, skin_elems, skin_rev_elems,
1009 [ + - ]: 26 : create_skin_elems, corners_only );
1010 : 26 : break;
1011 : : default:
1012 : 0 : rval = MB_TYPE_OUT_OF_RANGE;
1013 : 0 : break;
1014 : : }
1015 : :
1016 [ + - ]: 70 : thisMB->tag_delete( tag );
1017 : 70 : return rval;
1018 : : }
1019 : :
1020 : 4 : ErrorCode Skinner::find_skin_vertices_1D( Tag tag, const Range& edges, Range& skin_verts )
1021 : : {
1022 : : // This rather simple algorithm is provided for completeness
1023 : : // (not sure how often one really wants the 'skin' of a chain
1024 : : // or tangle of edges.)
1025 : : //
1026 : : // A vertex is on the skin of the edges if it is contained in exactly
1027 : : // one of the edges *in the input range*.
1028 : : //
1029 : : // This function expects the caller to have tagged all edges in the
1030 : : // input range with a value of one for the passed bit tag, and all
1031 : : // other edges with a value of zero. This allows us to do a faster
1032 : : // intersection with the input range and the edges adjacent to a vertex.
1033 : :
1034 : : ErrorCode rval;
1035 [ + - ]: 4 : Range::iterator hint = skin_verts.begin();
1036 : :
1037 : : // All input entities must be edges.
1038 [ + - ][ - + ]: 4 : if( !edges.all_of_dimension( 1 ) ) return MB_TYPE_OUT_OF_RANGE;
1039 : :
1040 : : // get all the vertices
1041 [ + - ]: 4 : Range verts;
1042 [ + - ]: 4 : rval = thisMB->get_connectivity( edges, verts, true );
1043 [ - + ]: 4 : if( MB_SUCCESS != rval ) return rval;
1044 : :
1045 : : // Test how many edges each input vertex is adjacent to.
1046 [ + - ]: 8 : std::vector< char > tag_vals;
1047 [ + - ]: 8 : std::vector< EntityHandle > adj;
1048 : : int n;
1049 [ + - ][ + - ]: 28 : for( Range::const_iterator it = verts.begin(); it != verts.end(); ++it )
[ + - ][ + - ]
[ + + ]
1050 : : {
1051 : : // get edges adjacent to vertex
1052 : 24 : adj.clear();
1053 [ + - ][ + - ]: 24 : rval = thisMB->get_adjacencies( &*it, 1, 1, false, adj );
1054 [ - + ]: 24 : if( MB_SUCCESS != rval ) return rval;
1055 [ - + ]: 24 : if( adj.empty() ) continue;
1056 : :
1057 : : // Intersect adjacent edges with the input list of edges
1058 [ + - ]: 24 : tag_vals.resize( adj.size() );
1059 [ + - ][ + - ]: 24 : rval = thisMB->tag_get_data( tag, &adj[0], adj.size(), &tag_vals[0] );
[ + - ]
1060 [ - + ]: 24 : if( MB_SUCCESS != rval ) return rval;
1061 : : #ifdef MOAB_OLD_STD_COUNT
1062 : : n = 0;
1063 : : std::count( tag_vals.begin(), tag_vals.end(), '\001', n );
1064 : : #else
1065 [ + - ]: 24 : n = std::count( tag_vals.begin(), tag_vals.end(), '\001' );
1066 : : #endif
1067 : : // If adjacent to only one input edge, then vertex is on skin
1068 [ + + ][ + - ]: 24 : if( n == 1 ) { hint = skin_verts.insert( hint, *it ); }
[ + - ]
1069 : : }
1070 : :
1071 : 8 : return MB_SUCCESS;
1072 : : }
1073 : :
1074 : : // A Container for storing a list of element sides adjacent
1075 : : // to a vertex. The template parameter is the number of
1076 : : // corners for the side.
1077 : : template < unsigned CORNERS >
1078 : 236 : class AdjSides
1079 : : {
1080 : : public:
1081 : : /**
1082 : : * This struct is used to for a reduced representation of element
1083 : : * "sides" adjacent to a give vertex. As such, it
1084 : : * a) does not store the initial vertex that all sides are adjacent to
1085 : : * b) orders the remaining vertices in a specific way for fast comparison.
1086 : : *
1087 : : * For edge elements, only the opposite vertex is stored.
1088 : : * For triangle elements, only the other two vertices are stored,
1089 : : * and they are stored with the smaller of those two handles first.
1090 : : * For quad elements, only the other three vertices are stored.
1091 : : * They are stored such that the vertex opposite the implicit (not
1092 : : * stored) vertex is always in slot 1. The other two vertices
1093 : : * (in slots 0 and 2) are ordered such that the handle of the one in
1094 : : * slot 0 is smaller than the handle in slot 2.
1095 : : *
1096 : : * For each side, the adj_elem field is used to store the element that
1097 : : * it is a side of as long as the element is considered to be on the skin.
1098 : : * The adj_elem field is cleared (set to zero) to indicate that this
1099 : : * side is no longer considered to be on the skin (and is the side of
1100 : : * more than one element.)
1101 : : */
1102 : : struct Side
1103 : : {
1104 : : EntityHandle handles[CORNERS - 1]; //!< side vertices, except for implicit one
1105 : : EntityHandle adj_elem; //!< element that this is a side of, or zero
1106 : 8720 : bool skin() const
1107 : : {
1108 : 8720 : return 0 != adj_elem;
1109 : : }
1110 : :
1111 : : /** construct from connectivity of side
1112 : : *\param array The connectivity of the element side.
1113 : : *\param idx The index of the implicit vertex (contained
1114 : : * in all sides in the list.)
1115 : : *\param adj The element that this is a side of.
1116 : : */
1117 : 73922 : Side( const EntityHandle* array, int idx, EntityHandle adj, unsigned short ) : adj_elem( adj )
1118 : : {
1119 : : switch( CORNERS )
1120 : : {
1121 : : case 3:
1122 : 50 : handles[1] = array[( idx + 2 ) % CORNERS];
1123 : : case 2:
1124 : 73922 : handles[0] = array[( idx + 1 ) % CORNERS];
1125 : 73922 : break;
1126 : : default:
1127 : : assert( false );
1128 : : break;
1129 : : }
1130 [ + + ]: 50 : if( CORNERS == 3 && handles[1] > handles[0] ) std::swap( handles[0], handles[1] );
1131 : 73922 : }
1132 : :
1133 : : /** construct from connectivity of parent element
1134 : : *\param array The connectivity of the parent element
1135 : : *\param idx The index of the implicit vertex (contained
1136 : : * in all sides in the list.) This is an index
1137 : : * into 'indices', not 'array'.
1138 : : *\param adj The element that this is a side of.
1139 : : *\param indices The indices into 'array' at which the vertices
1140 : : * representing the side occur.
1141 : : */
1142 : 120 : Side( const EntityHandle* array, int idx, EntityHandle adj, unsigned short, const short* indices )
1143 : 120 : : adj_elem( adj )
1144 : : {
1145 : : switch( CORNERS )
1146 : : {
1147 : : case 3:
1148 : 120 : handles[1] = array[indices[( idx + 2 ) % CORNERS]];
1149 : : case 2:
1150 : 120 : handles[0] = array[indices[( idx + 1 ) % CORNERS]];
1151 : 120 : break;
1152 : : default:
1153 : : assert( false );
1154 : : break;
1155 : : }
1156 [ + + ]: 120 : if( CORNERS == 3 && handles[1] > handles[0] ) std::swap( handles[0], handles[1] );
1157 : 120 : }
1158 : :
1159 : : // Compare two side instances. Relies in the ordering of
1160 : : // vertex handles as described above.
1161 : 213405 : bool operator==( const Side& other ) const
1162 : : {
1163 : : switch( CORNERS )
1164 : : {
1165 : : case 3:
1166 [ + + ][ + + ]: 232 : return handles[0] == other.handles[0] && handles[1] == other.handles[1];
1167 : : case 2:
1168 : 213173 : return handles[0] == other.handles[0];
1169 : : default:
1170 : : assert( false );
1171 : : return false;
1172 : : }
1173 : : }
1174 : : };
1175 : :
1176 : : private:
1177 : : std::vector< Side > data; //!< List of sides
1178 : : size_t skin_count; //!< Cached count of sides that are skin
1179 : :
1180 : : public:
1181 : : typedef typename std::vector< Side >::iterator iterator;
1182 : : typedef typename std::vector< Side >::const_iterator const_iterator;
1183 : 1469 : const_iterator begin() const
1184 : : {
1185 : 1469 : return data.begin();
1186 : : }
1187 : 10189 : const_iterator end() const
1188 : : {
1189 : 10189 : return data.end();
1190 : : }
1191 : :
1192 : 21312 : void clear()
1193 : : {
1194 : 21312 : data.clear();
1195 : 21312 : skin_count = 0;
1196 : 21312 : }
1197 : : bool empty() const
1198 : : {
1199 : : return data.empty();
1200 : : }
1201 : :
1202 : 236 : AdjSides() : skin_count( 0 ) {}
1203 : :
1204 : 28486 : size_t num_skin() const
1205 : : {
1206 : 28486 : return skin_count;
1207 : : }
1208 : :
1209 : : /** \brief insert side, specifying side connectivity
1210 : : *
1211 : : * Either insert a new side, or if the side is already in the
1212 : : * list, mark it as not on the skin.
1213 : : *
1214 : : *\param handles The connectivity of the element side.
1215 : : *\param skip_idx The index of the implicit vertex (contained
1216 : : * in all sides in the list.)
1217 : : *\param adj_elem The element that this is a side of.
1218 : : *\param elem_side Which side of adj_elem are we storing
1219 : : * (CN side number.)
1220 : : */
1221 : 71808 : void insert( const EntityHandle* handles, int skip_idx, EntityHandle adj_elem, unsigned short elem_side )
1222 : : {
1223 [ + - ][ # # ]: 71808 : Side side( handles, skip_idx, adj_elem, elem_side );
[ + - ]
1224 [ + - ][ # # ]: 71808 : iterator p = std::find( data.begin(), data.end(), side );
[ + - ]
1225 [ + - ][ + + ]: 71808 : if( p == data.end() )
[ # # ][ # # ]
[ + - ][ + + ]
1226 : : {
1227 [ + - ][ # # ]: 37164 : data.push_back( side );
[ + - ]
1228 : 37164 : ++skin_count; // not in list yet, so skin side (so far)
1229 : : }
1230 [ + - ][ + - ]: 34644 : else if( p->adj_elem )
[ # # ][ # # ]
[ + - ][ + - ]
1231 : : {
1232 [ + - ][ # # ]: 34644 : p->adj_elem = 0; // mark as not on skin
[ + - ]
1233 : 34644 : --skin_count; // decrement cached count of skin elements
1234 : : }
1235 : 71808 : }
1236 : :
1237 : : /** \brief insert side, specifying list of indices into parent element
1238 : : * connectivity.
1239 : : *
1240 : : * Either insert a new side, or if the side is already in the
1241 : : * list, mark it as not on the skin.
1242 : : *
1243 : : *\param handles The connectivity of the parent element
1244 : : *\param skip_idx The index of the implicit vertex (contained
1245 : : * in all sides in the list.) This is an index
1246 : : * into 'indices', not 'handles'.
1247 : : *\param adj_elem The element that this is a side of (parent handle).
1248 : : *\param indices The indices into 'handles' at which the vertices
1249 : : * representing the side occur.
1250 : : *\param elem_side Which side of adj_elem are we storing
1251 : : * (CN side number.)
1252 : : */
1253 : 73320 : void insert( const EntityHandle* handles, int skip_idx, EntityHandle adj_elem, unsigned short elem_side,
1254 : : const short* indices )
1255 : : {
1256 [ + - ][ + - ]: 73320 : Side side( handles, skip_idx, adj_elem, elem_side, indices );
1257 [ + - ][ + - ]: 73320 : iterator p = std::find( data.begin(), data.end(), side );
1258 [ + - ][ + + ]: 73320 : if( p == data.end() )
[ + - ][ + + ]
1259 : : {
1260 [ + - ][ + - ]: 42350 : data.push_back( side );
1261 : 42350 : ++skin_count; // not in list yet, so skin side (so far)
1262 : : }
1263 [ + - ][ + - ]: 30970 : else if( p->adj_elem )
[ + - ][ + - ]
1264 : : {
1265 [ + - ][ + - ]: 30970 : p->adj_elem = 0; // mark as not on skin
1266 : 30970 : --skin_count; // decrement cached count of skin elements
1267 : : }
1268 : 73320 : }
1269 : :
1270 : : /**\brief Search list for a given side, and if found, mark as not skin.
1271 : : *
1272 : : *\param other Connectivity of side
1273 : : *\param skip_index Index in 'other' at which implicit vertex occurs.
1274 : : *\param elem_out If return value is true, the element that the side is a
1275 : : * side of. If return value is false, not set.
1276 : : *\return true if found and marked not-skin, false if not found.
1277 : : *
1278 : : * Given the connectivity of some existing element, check if it occurs
1279 : : * in the list. If it does, clear the "is skin" state of the side so
1280 : : * that we know that we don't need to later create the side element.
1281 : : */
1282 : 9666 : bool find_and_unmark( const EntityHandle* other, int skip_index, EntityHandle& elem_out )
1283 : : {
1284 [ + - ][ + - ]: 9666 : Side s( other, skip_index, 0, 0 );
[ + - ]
1285 [ + - ][ + - ]: 9666 : iterator p = std::find( data.begin(), data.end(), s );
[ + - ]
1286 [ + - ][ + - ]: 9666 : if( p == data.end() || !p->adj_elem )
[ + - ][ + + ]
[ + - ]
[ + + # # ]
[ + - ][ + - ]
[ + - ][ - + ]
[ + - ]
[ - + # # ]
[ + - ][ + + ]
[ + - ][ + + ]
[ + - ]
[ + + # # ]
1287 : 938 : return false;
1288 : : else
1289 : : {
1290 [ + - ][ + - ]: 8728 : elem_out = p->adj_elem;
[ + - ]
1291 [ + - ][ + - ]: 8728 : p->adj_elem = 0; // clear "is skin" state for side
[ + - ]
1292 : 8728 : --skin_count; // decrement cached count of skin sides
1293 : 9666 : return true;
1294 : : }
1295 : : }
1296 : : };
1297 : :
1298 : : /** construct from connectivity of side
1299 : : *\param array The connectivity of the element side.
1300 : : *\param idx The index of the implicit vertex (contained
1301 : : * in all sides in the list.)
1302 : : *\param adj The element that this is a side of.
1303 : : */
1304 : : template <>
1305 : 7552 : AdjSides< 4 >::Side::Side( const EntityHandle* array, int idx, EntityHandle adj, unsigned short ) : adj_elem( adj )
1306 : : {
1307 : 7552 : const unsigned int CORNERS = 4;
1308 : 7552 : handles[2] = array[( idx + 3 ) % CORNERS];
1309 : 7552 : handles[1] = array[( idx + 2 ) % CORNERS];
1310 : 7552 : handles[0] = array[( idx + 1 ) % CORNERS];
1311 [ + + ]: 7552 : if( handles[2] > handles[0] ) std::swap( handles[0], handles[2] );
1312 : 7552 : }
1313 : :
1314 : : /** construct from connectivity of parent element
1315 : : *\param array The connectivity of the parent element
1316 : : *\param idx The index of the implicit vertex (contained
1317 : : * in all sides in the list.) This is an index
1318 : : * into 'indices', not 'array'.
1319 : : *\param adj The element that this is a side of.
1320 : : *\param indices The indices into 'array' at which the vertices
1321 : : * representing the side occur.
1322 : : */
1323 : : template <>
1324 : 73200 : AdjSides< 4 >::Side::Side( const EntityHandle* array, int idx, EntityHandle adj, unsigned short, const short* indices )
1325 : 73200 : : adj_elem( adj )
1326 : : {
1327 : 73200 : const unsigned int CORNERS = 4;
1328 : 73200 : handles[2] = array[indices[( idx + 3 ) % CORNERS]];
1329 : 73200 : handles[1] = array[indices[( idx + 2 ) % CORNERS]];
1330 : 73200 : handles[0] = array[indices[( idx + 1 ) % CORNERS]];
1331 [ + + ]: 73200 : if( handles[2] > handles[0] ) std::swap( handles[0], handles[2] );
1332 : 73200 : }
1333 : :
1334 : : // Compare two side instances. Relies in the ordering of
1335 : : // vertex handles as described above.
1336 : : template <>
1337 : 393423 : bool AdjSides< 4 >::Side::operator==( const Side& other ) const
1338 : : {
1339 [ + + ][ + + ]: 393423 : return handles[0] == other.handles[0] && handles[1] == other.handles[1] && handles[2] == other.handles[2];
[ + - ]
1340 : : }
1341 : :
1342 : : // Utility function used by find_skin_vertices_2D and
1343 : : // find_skin_vertices_3D to create elements representing
1344 : : // the skin side of a higher-dimension element if one
1345 : : // does not already exist.
1346 : : //
1347 : : // Some arguments may seem redundant, but they are used
1348 : : // to create the correct order of element when the input
1349 : : // element contains higher-order nodes.
1350 : : //
1351 : : // This function always creates elements that have a "forward"
1352 : : // orientation with respect to the parent element (have
1353 : : // nodes ordered the same as CN returns for the "side").
1354 : : //
1355 : : // elem - The higher-dimension element for which to create
1356 : : // a lower-dim element representing the side.
1357 : : // side_type - The EntityType of the desired side.
1358 : : // side_conn - The connectivity of the new side.
1359 : 1552 : ErrorCode Skinner::create_side( const EntityHandle this_set, EntityHandle elem, EntityType side_type,
1360 : : const EntityHandle* side_conn, EntityHandle& side_elem )
1361 : : {
1362 : 1552 : const int max_side = 9;
1363 : : const EntityHandle* conn;
1364 : : int len, side_len, side, sense, offset, indices[max_side];
1365 : : ErrorCode rval;
1366 [ + - ]: 1552 : EntityType type = TYPE_FROM_HANDLE( elem ), tmp_type;
1367 [ + - ]: 1552 : const int ncorner = CN::VerticesPerEntity( side_type );
1368 [ + - ]: 1552 : const int d = CN::Dimension( side_type );
1369 [ + - ]: 1552 : std::vector< EntityHandle > storage;
1370 : :
1371 : : // Get the connectivity of the parent element
1372 [ + - ]: 1552 : rval = thisMB->get_connectivity( elem, conn, len, false, &storage );
1373 [ - + ]: 1552 : if( MB_SUCCESS != rval ) return rval;
1374 : :
1375 : : // treat separately MBPOLYGON; we want to create the edge in the
1376 : : // forward sense always ; so figure out the sense first, then get out
1377 [ - + ][ # # ]: 1552 : if( MBPOLYGON == type && 1 == d && MBEDGE == side_type )
[ # # ]
1378 : : {
1379 : : // first find the first vertex in the conn list
1380 : 0 : int i = 0;
1381 [ # # ]: 0 : for( i = 0; i < len; i++ )
1382 : : {
1383 [ # # ]: 0 : if( conn[i] == side_conn[0] ) break;
1384 : : }
1385 [ # # ]: 0 : if( len == i ) return MB_FAILURE; // not found, big error
1386 : : // now, what if the polygon is padded?
1387 : : // the previous index is fine always. but the next one could be trouble :(
1388 : 0 : int prevIndex = ( i + len - 1 ) % len;
1389 : 0 : int nextIndex = ( i + 1 ) % len;
1390 : : // if the next index actually point to the same node, as current, it means it is padded
1391 [ # # ]: 0 : if( conn[nextIndex] == conn[i] )
1392 : : {
1393 : : // it really means we are at the end of proper nodes, the last nodes are repeated, so it
1394 : : // should be the first node
1395 : 0 : nextIndex = 0; // this is the first node!
1396 : : }
1397 : 0 : EntityHandle conn2[2] = { side_conn[0], side_conn[1] };
1398 [ # # ]: 0 : if( conn[prevIndex] == side_conn[1] )
1399 : : {
1400 : : // reverse, so the edge will be forward
1401 : 0 : conn2[0] = side_conn[1];
1402 : 0 : conn2[1] = side_conn[0];
1403 : : }
1404 [ # # ]: 0 : else if( conn[nextIndex] != side_conn[1] )
1405 : 0 : return MB_FAILURE; // it is not adjacent to the polygon
1406 : :
1407 [ # # ][ # # ]: 0 : rval = thisMB->create_element( MBEDGE, conn2, 2, side_elem );MB_CHK_ERR( rval );
[ # # ][ # # ]
1408 [ # # ]: 0 : if( this_set )
1409 : : {
1410 [ # # ][ # # ]: 0 : rval = thisMB->add_entities( this_set, &side_elem, 1 );MB_CHK_ERR( rval );
[ # # ][ # # ]
1411 : : }
1412 : 0 : return MB_SUCCESS;
1413 : : }
1414 : : // Find which side we are creating and get indices of all nodes
1415 : : // (including higher-order, if any.)
1416 [ + - ]: 1552 : CN::SideNumber( type, conn, side_conn, ncorner, d, side, sense, offset );
1417 [ + - ]: 1552 : CN::SubEntityNodeIndices( type, len, d, side, tmp_type, side_len, indices );
1418 [ - + ]: 1552 : assert( side_len <= max_side );
1419 [ - + ]: 1552 : assert( side_type == tmp_type );
1420 : :
1421 : : // NOTE: re-create conn array even when no higher-order nodes
1422 : : // because we want it to always be forward with respect
1423 : : // to the side ordering.
1424 : : EntityHandle side_conn_full[max_side];
1425 [ + + ]: 6964 : for( int i = 0; i < side_len; ++i )
1426 : 5412 : side_conn_full[i] = conn[indices[i]];
1427 : :
1428 [ + - ][ - + ]: 1552 : rval = thisMB->create_element( side_type, side_conn_full, side_len, side_elem );MB_CHK_ERR( rval );
[ # # ][ # # ]
1429 [ + + ]: 1552 : if( this_set )
1430 : : {
1431 [ + - ][ - + ]: 22 : rval = thisMB->add_entities( this_set, &side_elem, 1 );MB_CHK_ERR( rval );
[ # # ][ # # ]
1432 : : }
1433 : 1552 : return MB_SUCCESS;
1434 : 1552 : ;
1435 : : }
1436 : :
1437 : : // Test if an edge is reversed with respect CN's ordering
1438 : : // for the "side" of a face.
1439 : 10 : bool Skinner::edge_reversed( EntityHandle face, const EntityHandle* edge_ends )
1440 : : {
1441 : : const EntityHandle* conn;
1442 : : int len, idx;
1443 [ + - ]: 10 : ErrorCode rval = thisMB->get_connectivity( face, conn, len, true );
1444 [ - + ]: 10 : if( MB_SUCCESS != rval )
1445 : : {
1446 : 0 : assert( false );
1447 : : return false;
1448 : : }
1449 [ + - ]: 10 : idx = std::find( conn, conn + len, edge_ends[0] ) - conn;
1450 [ - + ]: 10 : if( idx == len )
1451 : : {
1452 : 0 : assert( false );
1453 : : return false;
1454 : : }
1455 : 10 : return ( edge_ends[1] == conn[( idx + len - 1 ) % len] );
1456 : : }
1457 : :
1458 : : // Test if a 2D element representing the side or face of a
1459 : : // volume element is reversed with respect to the CN node
1460 : : // ordering for the corresponding region element side.
1461 : 230 : bool Skinner::face_reversed( EntityHandle region, const EntityHandle* face_corners, EntityType face_type )
1462 : : {
1463 : : const EntityHandle* conn;
1464 : : int len, side, sense, offset;
1465 [ + - ]: 230 : ErrorCode rval = thisMB->get_connectivity( region, conn, len, true );
1466 [ - + ]: 230 : if( MB_SUCCESS != rval )
1467 : : {
1468 : 0 : assert( false );
1469 : : return false;
1470 : : }
1471 [ + - ]: 230 : short r = CN::SideNumber( TYPE_FROM_HANDLE( region ), conn, face_corners, CN::VerticesPerEntity( face_type ),
1472 [ + - ][ + - ]: 460 : CN::Dimension( face_type ), side, sense, offset );
[ + - ]
1473 [ - + ]: 230 : assert( 0 == r );
1474 [ + - ][ + + ]: 230 : return ( !r && sense == -1 );
1475 : : }
1476 : :
1477 : 40 : ErrorCode Skinner::find_skin_vertices_2D( const EntityHandle this_set, Tag tag, const Range& faces, Range* skin_verts,
1478 : : Range* skin_edges, Range* reversed_edges, bool create_edges,
1479 : : bool corners_only )
1480 : : {
1481 : : // This function iterates over all the vertices contained in the
1482 : : // input face list. For each such vertex, it then iterates over
1483 : : // all of the sides of the face elements adjacent to the vertex.
1484 : : // If an adjacent side is the side of only one of the input
1485 : : // faces, then that side is on the skin.
1486 : : //
1487 : : // This algorithm will visit each skin vertex exactly once. It
1488 : : // will visit each skin side once for each vertex in the side.
1489 : : //
1490 : : // This function expects the caller to have created the passed bit
1491 : : // tag and set it to one only for the faces in the passed range. This
1492 : : // tag is used to do a fast intersection of the faces adjacent to a
1493 : : // vertex with the faces in the input range (discard any for which the
1494 : : // tag is not set to one.)
1495 : :
1496 : : ErrorCode rval;
1497 : 40 : std::vector< EntityHandle >::iterator i, j;
1498 [ + - ]: 40 : Range::iterator hint;
1499 [ + + ][ + - ]: 40 : if( skin_verts ) hint = skin_verts->begin();
1500 [ + - ]: 40 : std::vector< EntityHandle > storage;
1501 : : const EntityHandle* conn;
1502 : : int len;
1503 [ + + ][ + + ]: 40 : bool find_edges = skin_edges || create_edges;
1504 : 40 : bool printed_nonconformal_ho_warning = false;
1505 : : EntityHandle face;
1506 : :
1507 [ + - ][ - + ]: 40 : if( !faces.all_of_dimension( 2 ) ) return MB_TYPE_OUT_OF_RANGE;
1508 : :
1509 : : // get all the vertices
1510 [ + - ]: 80 : Range verts;
1511 [ + - ]: 40 : rval = thisMB->get_connectivity( faces, verts, true );
1512 [ - + ]: 40 : if( MB_SUCCESS != rval ) return rval;
1513 : :
1514 [ + - ]: 80 : std::vector< char > tag_vals;
1515 [ + - ]: 80 : std::vector< EntityHandle > adj;
1516 [ + - ]: 80 : AdjSides< 2 > adj_edges;
1517 [ + - ][ + - ]: 6622 : for( Range::const_iterator it = verts.begin(); it != verts.end(); ++it )
[ + - ][ + - ]
[ + + ]
1518 : : {
1519 : 6582 : bool higher_order = false;
1520 : :
1521 : : // get all adjacent faces
1522 : 6582 : adj.clear();
1523 [ + - ][ + - ]: 6582 : rval = thisMB->get_adjacencies( &*it, 1, 2, false, adj );
1524 [ - + ]: 6582 : if( MB_SUCCESS != rval ) return rval;
1525 [ - + ]: 6582 : if( adj.empty() ) continue;
1526 : :
1527 : : // remove those not in the input list (intersect with input list)
1528 : 6582 : i = j = adj.begin();
1529 [ + - ]: 6582 : tag_vals.resize( adj.size() );
1530 [ + - ][ + - ]: 6582 : rval = thisMB->tag_get_data( tag, &adj[0], adj.size(), &tag_vals[0] );
[ + - ]
1531 [ - + ]: 6582 : if( MB_SUCCESS != rval ) return rval;
1532 : : // remove non-tagged entries
1533 : 6582 : i = j = adj.begin();
1534 [ + - ][ + - ]: 42844 : for( ; i != adj.end(); ++i )
[ + + ]
1535 [ + - ][ + - ]: 36262 : if( tag_vals[i - adj.begin()] ) *( j++ ) = *i;
[ + + ][ + - ]
[ + - ][ + - ]
1536 [ + - ]: 6582 : adj.erase( j, adj.end() );
1537 : :
1538 : : // For each adjacent face, check the edges adjacent to the current vertex
1539 [ + - ]: 6582 : adj_edges.clear(); // other vertex for adjacent edges
1540 [ + - ][ + - ]: 42054 : for( i = adj.begin(); i != adj.end(); ++i )
[ + + ]
1541 : : {
1542 [ + - ][ + - ]: 35472 : rval = thisMB->get_connectivity( *i, conn, len, false, &storage );
1543 [ - + ]: 35472 : if( MB_SUCCESS != rval ) return rval;
1544 : :
1545 : : // For a single face element adjacent to this vertex, there
1546 : : // will be exactly two sides (edges) adjacent to the vertex.
1547 : : // Find the other vertex for each of the two edges.
1548 : :
1549 : : EntityHandle prev, next; // vertices of two adjacent edge-sides
1550 [ + - ][ + - ]: 35472 : const int idx = std::find( conn, conn + len, *it ) - conn;
1551 [ - + ]: 35472 : assert( idx != len );
1552 : :
1553 [ + - ][ + - ]: 35472 : if( TYPE_FROM_HANDLE( *i ) == MBTRI && len > 3 )
[ + + ][ + + ]
[ + + ]
1554 : : {
1555 : 12 : len = 3;
1556 : 12 : higher_order = true;
1557 [ - + ]: 12 : if( idx > 2 )
1558 : : { // skip higher-order nodes for now
1559 [ # # ]: 0 : if( !printed_nonconformal_ho_warning )
1560 : : {
1561 : 0 : printed_nonconformal_ho_warning = true;
1562 [ # # ]: 0 : std::cerr << "Non-conformal higher-order mesh detected in skinner: "
1563 [ # # ][ # # ]: 0 : << "vertex " << ID_FROM_HANDLE( *it ) << " is a corner in "
[ # # ][ # # ]
[ # # ]
1564 [ # # ][ # # ]: 0 : << "some elements and a higher-order node in others" << std::endl;
1565 : : }
1566 : 0 : continue;
1567 : : }
1568 : : }
1569 [ + - ][ + - ]: 35460 : else if( TYPE_FROM_HANDLE( *i ) == MBQUAD && len > 4 )
[ + + ][ + + ]
[ + + ]
1570 : : {
1571 : 16 : len = 4;
1572 : 16 : higher_order = true;
1573 [ - + ]: 16 : if( idx > 3 )
1574 : : { // skip higher-order nodes for now
1575 [ # # ]: 0 : if( !printed_nonconformal_ho_warning )
1576 : : {
1577 : 0 : printed_nonconformal_ho_warning = true;
1578 [ # # ]: 0 : std::cerr << "Non-conformal higher-order mesh detected in skinner: "
1579 [ # # ][ # # ]: 0 : << "vertex " << ID_FROM_HANDLE( *it ) << " is a corner in "
[ # # ][ # # ]
[ # # ]
1580 [ # # ][ # # ]: 0 : << "some elements and a higher-order node in others" << std::endl;
1581 : : }
1582 : 0 : continue;
1583 : : }
1584 : : }
1585 : :
1586 : : // so it must be a MBPOLYGON
1587 : 35472 : const int prev_idx = ( idx + len - 1 ) % len; // this should be fine, always, even for padded case
1588 : 35472 : prev = conn[prev_idx];
1589 : 35472 : next = conn[( idx + 1 ) % len];
1590 [ - + ]: 35472 : if( next == conn[idx] ) // it must be the padded case, so roll to the beginning
1591 : 0 : next = conn[0];
1592 : :
1593 : : // Insert sides (edges) in our list of candidate skin sides
1594 [ + - ][ + - ]: 35472 : adj_edges.insert( &prev, 1, *i, prev_idx );
1595 [ + - ][ + - ]: 35472 : adj_edges.insert( &next, 1, *i, idx );
1596 : : }
1597 : :
1598 : : // If vertex is not on skin, advance to next vertex.
1599 : : // adj_edges handled checking for duplicates on insertion.
1600 : : // If every candidate skin edge occurred more than once (was
1601 : : // not in fact on the skin), then we're done with this vertex.
1602 [ + - ][ + + ]: 6582 : if( 0 == adj_edges.num_skin() ) continue;
1603 : :
1604 : : // If user requested Range of *vertices* on the skin...
1605 [ + + ]: 1092 : if( skin_verts )
1606 : : {
1607 : : // Put skin vertex in output list
1608 [ + - ][ + - ]: 56 : hint = skin_verts->insert( hint, *it );
1609 : :
1610 : : // Add mid edge nodes to vertex list
1611 [ + - ][ + + ]: 56 : if( !corners_only && higher_order )
1612 : : {
1613 [ + - ][ + - ]: 80 : for( AdjSides< 2 >::const_iterator p = adj_edges.begin(); p != adj_edges.end(); ++p )
[ + - ][ + - ]
[ + + ]
1614 : : {
1615 [ + - ][ + - ]: 24 : if( p->skin() )
[ + + ]
1616 : : {
1617 [ + - ]: 20 : face = p->adj_elem;
1618 [ + - ]: 20 : EntityType type = TYPE_FROM_HANDLE( face );
1619 : :
1620 [ + - ]: 20 : rval = thisMB->get_connectivity( face, conn, len, false );
1621 [ - + ]: 20 : if( MB_SUCCESS != rval ) return rval;
1622 [ + - ][ - + ]: 20 : if( !CN::HasMidEdgeNodes( type, len ) ) continue;
1623 : :
1624 [ + - ][ + - ]: 20 : EntityHandle ec[2] = { *it, p->handles[0] };
1625 : : int side, sense, offset;
1626 [ + - ]: 20 : CN::SideNumber( type, conn, ec, 2, 1, side, sense, offset );
1627 [ + - ]: 20 : offset = CN::HONodeIndex( type, len, 1, side );
1628 [ + - ][ - + ]: 20 : assert( offset >= 0 && offset < len );
1629 [ + - ]: 20 : skin_verts->insert( conn[offset] );
1630 : : }
1631 : : }
1632 : : }
1633 : : }
1634 : :
1635 : : // If user requested Range of *edges* on the skin...
1636 [ + + ]: 1092 : if( find_edges )
1637 : : {
1638 : : // Search list of existing adjacent edges for any that are on the skin
1639 : 1072 : adj.clear();
1640 [ + - ][ + - ]: 1072 : rval = thisMB->get_adjacencies( &*it, 1, 1, false, adj );
1641 [ - + ]: 1072 : if( MB_SUCCESS != rval ) return rval;
1642 [ + - ][ + - ]: 3584 : for( i = adj.begin(); i != adj.end(); ++i )
[ + + ]
1643 : : {
1644 [ + - ][ + - ]: 2512 : rval = thisMB->get_connectivity( *i, conn, len, true );
1645 [ - + ]: 2512 : if( MB_SUCCESS != rval ) return rval;
1646 : :
1647 : : // bool equality expression within find_and_unmark call
1648 : : // will be evaluate to the index of *it in the conn array.
1649 : : //
1650 : : // Note that the order of the terms in the if statement is important.
1651 : : // We want to unmark any existing skin edges even if we aren't
1652 : : // returning them. Otherwise we'll end up creating duplicates
1653 : : // if create_edges is true and skin_edges is not.
1654 [ + - ][ + - ]: 2512 : if( adj_edges.find_and_unmark( conn, ( conn[1] == *it ), face ) && skin_edges )
[ + + ][ + + ]
[ + + ]
1655 : : {
1656 [ + + ][ + - ]: 1666 : if( reversed_edges && edge_reversed( face, conn ) )
[ + + ][ + + ]
1657 [ + - ][ + - ]: 4 : reversed_edges->insert( *i );
1658 : : else
1659 [ + - ][ + - ]: 1666 : skin_edges->insert( *i );
1660 : : }
1661 : : }
1662 : : }
1663 : :
1664 : : // If the user requested that we create new edges for sides
1665 : : // on the skin for which there is no existing edge, and there
1666 : : // are still skin sides for which no corresponding edge was found...
1667 [ + + ][ + - ]: 1092 : if( create_edges && adj_edges.num_skin() )
[ + + ][ + + ]
1668 : : {
1669 : : // Create any skin edges that don't exist
1670 [ + - ][ + - ]: 1831 : for( AdjSides< 2 >::const_iterator p = adj_edges.begin(); p != adj_edges.end(); ++p )
[ + - ][ + - ]
[ + + ]
1671 : : {
1672 [ + - ][ + - ]: 1422 : if( p->skin() )
[ + + ]
1673 : : {
1674 [ + - ][ + - ]: 442 : EntityHandle edge, ec[] = { *it, p->handles[0] };
1675 [ + - ][ + - ]: 442 : rval = create_side( this_set, p->adj_elem, MBEDGE, ec, edge );
1676 [ - + ]: 442 : if( MB_SUCCESS != rval ) return rval;
1677 [ + + ][ + - ]: 442 : if( skin_edges ) skin_edges->insert( edge );
1678 : : }
1679 : : }
1680 : : }
1681 : :
1682 : : } // end for each vertex
1683 : :
1684 : 80 : return MB_SUCCESS;
1685 : : }
1686 : :
1687 : 26 : ErrorCode Skinner::find_skin_vertices_3D( const EntityHandle this_set, Tag tag, const Range& entities,
1688 : : Range* skin_verts, Range* skin_faces, Range* reversed_faces,
1689 : : bool create_faces, bool corners_only )
1690 : : {
1691 : : // This function iterates over all the vertices contained in the
1692 : : // input vol elem list. For each such vertex, it then iterates over
1693 : : // all of the sides of the vol elements adjacent to the vertex.
1694 : : // If an adjacent side is the side of only one of the input
1695 : : // elements, then that side is on the skin.
1696 : : //
1697 : : // This algorithm will visit each skin vertex exactly once. It
1698 : : // will visit each skin side once for each vertex in the side.
1699 : : //
1700 : : // This function expects the caller to have created the passed bit
1701 : : // tag and set it to one only for the elements in the passed range. This
1702 : : // tag is used to do a fast intersection of the elements adjacent to a
1703 : : // vertex with the elements in the input range (discard any for which the
1704 : : // tag is not set to one.)
1705 : : //
1706 : : // For each vertex, iterate over each adjacent element. Construct
1707 : : // lists of the sides of each adjacent element that contain the vertex.
1708 : : //
1709 : : // A list of three-vertex sides is kept for all triangular faces,
1710 : : // included three-vertex faces of type MBPOLYGON. Putting polygons
1711 : : // in the same list ensures that we find polyhedron and non-polyhedron
1712 : : // elements that are adjacent.
1713 : : //
1714 : : // A list of four-vertex sides is kept for all quadrilateral faces,
1715 : : // including four-vertex faces of type MBPOLYGON.
1716 : : //
1717 : : // Sides with more than four vertices must have an explicit MBPOLYGON
1718 : : // element representing them because MBPOLYHEDRON connectivity is a
1719 : : // list of faces rather than vertices. So the third list (vertices>=5),
1720 : : // need contain only the handle of the face rather than the vertex handles.
1721 : :
1722 : : ErrorCode rval;
1723 : 26 : std::vector< EntityHandle >::iterator i, j;
1724 [ + - ]: 26 : Range::iterator hint;
1725 [ + + ][ + - ]: 26 : if( skin_verts ) hint = skin_verts->begin();
1726 [ + - ][ + - ]: 52 : std::vector< EntityHandle > storage, storage2; // temp storage for conn lists
1727 : : const EntityHandle *conn, *conn2;
1728 : : int len, len2;
1729 [ + + ][ + + ]: 26 : bool find_faces = skin_faces || create_faces;
1730 : : int clen, side, sense, offset, indices[9];
1731 : : EntityType face_type;
1732 : : EntityHandle elem;
1733 : 26 : bool printed_nonconformal_ho_warning = false;
1734 : :
1735 [ + - ][ - + ]: 26 : if( !entities.all_of_dimension( 3 ) ) return MB_TYPE_OUT_OF_RANGE;
1736 : :
1737 : : // get all the vertices
1738 [ + - ]: 52 : Range verts;
1739 [ + - ]: 26 : rval = thisMB->get_connectivity( entities, verts, true );
1740 [ - + ]: 26 : if( MB_SUCCESS != rval ) return rval;
1741 : : // if there are polyhedra in the input list, need to make another
1742 : : // call to get vertices from faces
1743 [ + - ][ + + ]: 26 : if( !verts.all_of_dimension( 0 ) )
1744 : : {
1745 [ + - ]: 3 : Range::iterator it = verts.upper_bound( MBVERTEX );
1746 [ + - ]: 3 : Range pfaces;
1747 [ + - ][ + - ]: 3 : pfaces.merge( it, verts.end() );
1748 [ + - ][ + - ]: 3 : verts.erase( it, verts.end() );
1749 [ + - ]: 3 : rval = thisMB->get_connectivity( pfaces, verts, true );
1750 [ - + ]: 3 : if( MB_SUCCESS != rval ) return rval;
1751 [ + - ][ - + ]: 3 : assert( verts.all_of_dimension( 0 ) );
[ + - ]
1752 : : }
1753 : :
1754 [ + - ]: 52 : AdjSides< 4 > adj_quads; // 4-node sides adjacent to a vertex
1755 [ + - ]: 52 : AdjSides< 3 > adj_tris; // 3-node sides adjacent to a vertex
1756 [ + - ]: 52 : AdjSides< 2 > adj_poly; // n-node sides (n>5) adjacent to vertex
1757 : : // (must have an explicit polygon, so store
1758 : : // polygon handle rather than vertices.)
1759 [ + - ]: 52 : std::vector< char > tag_vals;
1760 [ + - ]: 52 : std::vector< EntityHandle > adj;
1761 [ + - ][ + - ]: 4936 : for( Range::const_iterator it = verts.begin(); it != verts.end(); ++it )
[ + - ][ + - ]
[ + + ]
1762 : : {
1763 : 4910 : bool higher_order = false;
1764 : :
1765 : : // get all adjacent elements
1766 : 4910 : adj.clear();
1767 [ + - ][ + - ]: 4910 : rval = thisMB->get_adjacencies( &*it, 1, 3, false, adj );
1768 [ - + ]: 4910 : if( MB_SUCCESS != rval ) return rval;
1769 [ - + ]: 4910 : if( adj.empty() ) continue;
1770 : :
1771 : : // remove those not tagged (intersect with input range)
1772 : 4910 : i = j = adj.begin();
1773 [ + - ]: 4910 : tag_vals.resize( adj.size() );
1774 [ + - ][ + - ]: 4910 : rval = thisMB->tag_get_data( tag, &adj[0], adj.size(), &tag_vals[0] );
[ + - ]
1775 [ - + ]: 4910 : if( MB_SUCCESS != rval ) return rval;
1776 [ + - ][ + - ]: 29678 : for( ; i != adj.end(); ++i )
[ + + ]
1777 [ + - ][ + - ]: 24768 : if( tag_vals[i - adj.begin()] ) *( j++ ) = *i;
[ + + ][ + - ]
[ + - ][ + - ]
1778 [ + - ]: 4910 : adj.erase( j, adj.end() );
1779 : :
1780 : : // Build lists of sides of 3D element adjacent to the current vertex
1781 [ + - ]: 4910 : adj_quads.clear(); // store three other vertices for each adjacent quad face
1782 [ + - ]: 4910 : adj_tris.clear(); // store two other vertices for each adjacent tri face
1783 [ + - ]: 4910 : adj_poly.clear(); // store handle of each adjacent polygonal face
1784 : : int idx;
1785 [ + - ][ + - ]: 29638 : for( i = adj.begin(); i != adj.end(); ++i )
[ + + ]
1786 : : {
1787 [ + - ][ + - ]: 24728 : const EntityType type = TYPE_FROM_HANDLE( *i );
1788 : :
1789 : : // Special case for POLYHEDRA
1790 [ + + ]: 24728 : if( type == MBPOLYHEDRON )
1791 : : {
1792 [ + - ][ + - ]: 288 : rval = thisMB->get_connectivity( *i, conn, len );
1793 [ - + ]: 288 : if( MB_SUCCESS != rval ) return rval;
1794 [ + + ]: 2592 : for( int k = 0; k < len; ++k )
1795 : : {
1796 [ + - ]: 2304 : rval = thisMB->get_connectivity( conn[k], conn2, len2, true, &storage2 );
1797 [ - + ]: 2304 : if( MB_SUCCESS != rval ) return rval;
1798 [ + - ][ + - ]: 2304 : idx = std::find( conn2, conn2 + len2, *it ) - conn2;
1799 [ + + ]: 2304 : if( idx == len2 ) // vertex not in this face
1800 : 1440 : continue;
1801 : :
1802 : : // Treat 3- and 4-vertex faces specially, so that
1803 : : // if the mesh contains both elements and polyhedra,
1804 : : // we don't miss one type adjacent to the other.
1805 [ - + + ]: 864 : switch( len2 )
1806 : : {
1807 : : case 3:
1808 [ # # ][ # # ]: 0 : adj_tris.insert( conn2, idx, *i, k );
1809 : 0 : break;
1810 : : case 4:
1811 [ + - ][ + - ]: 576 : adj_quads.insert( conn2, idx, *i, k );
1812 : 576 : break;
1813 : : default:
1814 [ + - ][ + - ]: 288 : adj_poly.insert( conn + k, 1, *i, k );
1815 : 288 : break;
1816 : : }
1817 : : }
1818 : : }
1819 : : else
1820 : : {
1821 [ + - ][ + - ]: 24440 : rval = thisMB->get_connectivity( *i, conn, len, false, &storage );
1822 [ - + ]: 24440 : if( MB_SUCCESS != rval ) return rval;
1823 : :
1824 [ + - ][ + - ]: 24440 : idx = std::find( conn, conn + len, *it ) - conn;
1825 [ - + ]: 24440 : assert( idx != len );
1826 : :
1827 [ + - ][ + + ]: 24440 : if( len > CN::VerticesPerEntity( type ) )
1828 : : {
1829 : 64 : higher_order = true;
1830 : : // skip higher-order nodes for now
1831 [ + - ][ - + ]: 64 : if( idx >= CN::VerticesPerEntity( type ) )
1832 : : {
1833 [ # # ]: 0 : if( !printed_nonconformal_ho_warning )
1834 : : {
1835 : 0 : printed_nonconformal_ho_warning = true;
1836 [ # # ]: 0 : std::cerr << "Non-conformal higher-order mesh detected in skinner: "
1837 [ # # ][ # # ]: 0 : << "vertex " << ID_FROM_HANDLE( *it ) << " is a corner in "
[ # # ][ # # ]
[ # # ]
1838 [ # # ][ # # ]: 0 : << "some elements and a higher-order node in others" << std::endl;
1839 : : }
1840 : 0 : continue;
1841 : : }
1842 : : }
1843 : :
1844 : : // For each side of the element...
1845 [ + - ]: 24440 : const int num_faces = CN::NumSubEntities( type, 2 );
1846 [ + + ]: 170976 : for( int f = 0; f < num_faces; ++f )
1847 : : {
1848 : : int num_vtx;
1849 [ + - ]: 146536 : const short* face_indices = CN::SubEntityVertexIndices( type, 2, f, face_type, num_vtx );
1850 [ + - ]: 146536 : const short face_idx = std::find( face_indices, face_indices + num_vtx, (short)idx ) - face_indices;
1851 : : // skip sides that do not contain vertex from outer loop
1852 [ + + ]: 146536 : if( face_idx == num_vtx ) continue; // current vertex not in this face
1853 : :
1854 [ - + ]: 73320 : assert( num_vtx <= 4 ); // polyhedra handled above
1855 [ + + - ]: 73320 : switch( face_type )
1856 : : {
1857 : : case MBTRI:
1858 [ + - ][ + - ]: 120 : adj_tris.insert( conn, face_idx, *i, f, face_indices );
1859 : 120 : break;
1860 : : case MBQUAD:
1861 [ + - ][ + - ]: 73200 : adj_quads.insert( conn, face_idx, *i, f, face_indices );
1862 : 73200 : break;
1863 : : default:
1864 : 73320 : return MB_TYPE_OUT_OF_RANGE;
1865 : : }
1866 : : }
1867 : : }
1868 : : } // end for (adj[3])
1869 : :
1870 : : // If vertex is not on skin, advance to next vertex
1871 [ + - ][ + - ]: 4910 : if( 0 == ( adj_tris.num_skin() + adj_quads.num_skin() + adj_poly.num_skin() ) ) continue;
[ + - ][ + + ]
1872 : :
1873 : : // If user requested that skin *vertices* be passed back...
1874 [ + + ]: 3002 : if( skin_verts )
1875 : : {
1876 : : // Put skin vertex in output list
1877 [ + - ][ + - ]: 1814 : hint = skin_verts->insert( hint, *it );
1878 : :
1879 : : // Add mid-edge and mid-face nodes to vertex list
1880 [ + - ][ + + ]: 1814 : if( !corners_only && higher_order )
1881 : : {
1882 [ + - ][ # # ]: 24 : for( AdjSides< 3 >::const_iterator t = adj_tris.begin(); t != adj_tris.end(); ++t )
[ + - ][ + - ]
[ - + ]
1883 : : {
1884 [ # # ][ # # ]: 0 : if( t->skin() )
[ # # ]
1885 : : {
1886 [ # # ]: 0 : elem = t->adj_elem;
1887 [ # # ]: 0 : EntityType type = TYPE_FROM_HANDLE( elem );
1888 : :
1889 [ # # ]: 0 : rval = thisMB->get_connectivity( elem, conn, len, false );
1890 [ # # ]: 0 : if( MB_SUCCESS != rval ) return rval;
1891 [ # # ][ # # ]: 0 : if( !CN::HasMidNodes( type, len ) ) continue;
1892 : :
1893 [ # # ][ # # ]: 0 : EntityHandle ec[3] = { *it, t->handles[0], t->handles[1] };
[ # # ]
1894 [ # # ]: 0 : CN::SideNumber( type, conn, ec, 3, 2, side, sense, offset );
1895 [ # # ]: 0 : CN::SubEntityNodeIndices( type, len, 2, side, face_type, clen, indices );
1896 [ # # ]: 0 : assert( MBTRI == face_type );
1897 [ # # ]: 0 : for( int k = 3; k < clen; ++k )
1898 [ # # ]: 0 : skin_verts->insert( conn[indices[k]] );
1899 : : }
1900 : : }
1901 [ + - ][ + - ]: 1902 : for( AdjSides< 4 >::const_iterator q = adj_quads.begin(); q != adj_quads.end(); ++q )
[ + - ][ + - ]
[ + + ]
1902 : : {
1903 [ + - ][ + - ]: 88 : if( q->skin() )
[ + + ]
1904 : : {
1905 [ + - ]: 80 : elem = q->adj_elem;
1906 [ + - ]: 80 : EntityType type = TYPE_FROM_HANDLE( elem );
1907 : :
1908 [ + - ]: 80 : rval = thisMB->get_connectivity( elem, conn, len, false );
1909 [ - + ]: 80 : if( MB_SUCCESS != rval ) return rval;
1910 [ + - ][ - + ]: 80 : if( !CN::HasMidNodes( type, len ) ) continue;
1911 : :
1912 [ + - ][ + - ]: 80 : EntityHandle ec[4] = { *it, q->handles[0], q->handles[1], q->handles[2] };
[ + - ][ + - ]
1913 [ + - ]: 80 : CN::SideNumber( type, conn, ec, 4, 2, side, sense, offset );
1914 [ + - ]: 80 : CN::SubEntityNodeIndices( type, len, 2, side, face_type, clen, indices );
1915 [ - + ]: 80 : assert( MBQUAD == face_type );
1916 [ + + ]: 480 : for( int k = 4; k < clen; ++k )
1917 [ + - ]: 400 : skin_verts->insert( conn[indices[k]] );
1918 : : }
1919 : : }
1920 : : }
1921 : : }
1922 : :
1923 : : // If user requested that we pass back the list of 2D elements
1924 : : // representing the skin of the mesh...
1925 [ + + ]: 3002 : if( find_faces )
1926 : : {
1927 : : // Search list of adjacent faces for any that are on the skin
1928 : 2078 : adj.clear();
1929 [ + - ][ + - ]: 2078 : rval = thisMB->get_adjacencies( &*it, 1, 2, false, adj );
1930 [ - + ]: 2078 : if( MB_SUCCESS != rval ) return rval;
1931 : :
1932 [ + - ][ + - ]: 9232 : for( i = adj.begin(); i != adj.end(); ++i )
[ + + ]
1933 : : {
1934 [ + - ][ + - ]: 7154 : rval = thisMB->get_connectivity( *i, conn, len, true );
1935 [ - + ]: 7154 : if( MB_SUCCESS != rval ) return rval;
1936 [ + - ][ + - ]: 7154 : const int idx2 = std::find( conn, conn + len, *it ) - conn;
1937 [ - + ]: 7154 : if( idx2 >= len )
1938 : : {
1939 [ # # ]: 0 : assert( printed_nonconformal_ho_warning );
1940 : 0 : continue;
1941 : : }
1942 : :
1943 : : // Note that the order of the terms in the if statements below
1944 : : // is important. We want to unmark any existing skin faces even
1945 : : // if we aren't returning them. Otherwise we'll end up creating
1946 : : // duplicates if create_faces is true.
1947 [ + + ]: 7154 : if( 3 == len )
1948 : : {
1949 [ + - ][ + - ]: 50 : if( adj_tris.find_and_unmark( conn, idx2, elem ) && skin_faces )
[ + - ][ + - ]
1950 : : {
1951 [ + + ][ + - ]: 50 : if( reversed_faces && face_reversed( elem, conn, MBTRI ) )
[ + + ][ + + ]
1952 [ + - ][ + - ]: 6 : reversed_faces->insert( *i );
1953 : : else
1954 [ + - ][ + - ]: 50 : skin_faces->insert( *i );
1955 : : }
1956 : : }
1957 [ + + ]: 7104 : else if( 4 == len )
1958 : : {
1959 [ + - ][ + + ]: 6976 : if( adj_quads.find_and_unmark( conn, idx2, elem ) && skin_faces )
[ + + ][ + + ]
1960 : : {
1961 [ + + ][ + - ]: 3822 : if( reversed_faces && face_reversed( elem, conn, MBQUAD ) )
[ - + ][ - + ]
1962 [ # # ][ # # ]: 0 : reversed_faces->insert( *i );
1963 : : else
1964 [ + - ][ + - ]: 3822 : skin_faces->insert( *i );
1965 : : }
1966 : : }
1967 : : else
1968 : : {
1969 [ + - ][ + - ]: 128 : if( adj_poly.find_and_unmark( &*i, 1, elem ) && skin_faces ) skin_faces->insert( *i );
[ + + ][ + - ]
[ + + ][ + - ]
[ + - ]
1970 : : }
1971 : : }
1972 : : }
1973 : :
1974 : : // If user does not want use to create new faces representing
1975 : : // sides for which there is currently no explicit element,
1976 : : // skip the remaining code and advance the outer loop to the
1977 : : // next vertex.
1978 [ + + ]: 3002 : if( !create_faces ) continue;
1979 : :
1980 : : // Polyhedra always have explicitly defined faces, so
1981 : : // there is no way we could need to create such a face.
1982 [ + - ][ - + ]: 2034 : assert( 0 == adj_poly.num_skin() );
1983 : :
1984 : : // Create any skin tris that don't exist
1985 [ + - ][ + + ]: 2034 : if( adj_tris.num_skin() )
1986 : : {
1987 [ + - ][ + - ]: 60 : for( AdjSides< 3 >::const_iterator t = adj_tris.begin(); t != adj_tris.end(); ++t )
[ + - ][ + - ]
[ + + ]
1988 : : {
1989 [ + - ][ + - ]: 42 : if( t->skin() )
[ + + ]
1990 : : {
1991 [ + - ][ + - ]: 22 : EntityHandle tri, c[3] = { *it, t->handles[0], t->handles[1] };
[ + - ]
1992 [ + - ][ + - ]: 22 : rval = create_side( this_set, t->adj_elem, MBTRI, c, tri );
1993 [ - + ]: 22 : if( MB_SUCCESS != rval ) return rval;
1994 [ + - ][ + - ]: 22 : if( skin_faces ) skin_faces->insert( tri );
1995 : : }
1996 : : }
1997 : : }
1998 : :
1999 : : // Create any skin quads that don't exist
2000 [ + - ][ + + ]: 2034 : if( adj_quads.num_skin() )
2001 : : {
2002 [ + - ][ + - ]: 8128 : for( AdjSides< 4 >::const_iterator q = adj_quads.begin(); q != adj_quads.end(); ++q )
[ + - ][ + - ]
[ + + ]
2003 : : {
2004 [ + - ][ + - ]: 7144 : if( q->skin() )
[ + + ]
2005 : : {
2006 [ + - ][ + - ]: 1088 : EntityHandle quad, c[4] = { *it, q->handles[0], q->handles[1], q->handles[2] };
[ + - ][ + - ]
2007 [ + - ][ + - ]: 1088 : rval = create_side( this_set, q->adj_elem, MBQUAD, c, quad );
2008 [ - + ]: 1088 : if( MB_SUCCESS != rval ) return rval;
2009 [ + + ][ + - ]: 1088 : if( skin_faces ) skin_faces->insert( quad );
2010 : : }
2011 : : }
2012 : : }
2013 : : } // end for each vertex
2014 : :
2015 : 52 : return MB_SUCCESS;
2016 : : }
2017 : :
2018 [ + - ][ + - ]: 228 : } // namespace moab
|
