Mesh Oriented datABase  (version 5.5.0)
An array-based unstructured mesh library
CN.cpp
Go to the documentation of this file.
1 /**
2  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3  * storing and accessing finite element mesh data.
4  *
5  * Copyright 2004 Sandia Corporation. Under the terms of Contract
6  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7  * retains certain rights in this software.
8  *
9  * This library is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 2.1 of the License, or (at your option) any later version.
13  *
14  */
15 
16 #include "moab/CN.hpp"
17 #include "MBCNArrays.hpp"
18 #include "MBCN.h"
19 #include <cassert>
20 #include <cstring>
21 #include <iterator>
22 
23 namespace moab
24 {
25 
26 const char* CN::entityTypeNames[] = { "Vertex", "Edge", "Tri", "Quad", "Polygon", "Tet", "Pyramid",
27  "Prism", "Knife", "Hex", "Polyhedron", "EntitySet", "MaxType" };
28 
29 short int CN::numberBasis = 0;
30 
31 short int CN::permuteVec[MBMAXTYPE][3][MAX_SUB_ENTITIES + 1];
33 
38 
39 short CN::increasingInts[] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
40  20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39 };
41 
42 //! set the basis of the numbering system; may or may not do things besides setting the
43 //! member variable
44 void CN::SetBasis( const int in_basis )
45 {
46  numberBasis = in_basis;
47 }
48 
49 /// Get the dimension pair corresponding to a dimension
50 DimensionPair CN::getDimPair( int entity_type )
51 {
52  return TypeDimensionMap[entity_type];
53 }
54 
55 //! return a type for the given name
56 EntityType CN::EntityTypeFromName( const char* name )
57 {
58  for( EntityType i = MBVERTEX; i < MBMAXTYPE; i++ )
59  {
60  if( 0 == strcmp( name, entityTypeNames[i] ) ) return i;
61  }
62 
63  return MBMAXTYPE;
64 }
65 
66 void CN::SubEntityNodeIndices( const EntityType this_topo,
67  const int num_nodes,
68  const int sub_dimension,
69  const int sub_index,
70  EntityType& subentity_topo,
71  int& num_sub_entity_nodes,
72  int sub_entity_conn[] )
73 {
74  // If asked for a node, the special case...
75  if( sub_dimension == 0 )
76  {
77  assert( sub_index < num_nodes );
78  subentity_topo = MBVERTEX;
79  num_sub_entity_nodes = 1;
80  sub_entity_conn[0] = sub_index;
81  return;
82  }
83 
84  const int ho_bits = HasMidNodes( this_topo, num_nodes );
85  subentity_topo = SubEntityType( this_topo, sub_dimension, sub_index );
86  num_sub_entity_nodes = VerticesPerEntity( subentity_topo );
87  const short* corners = mConnectivityMap[this_topo][sub_dimension - 1].conn[sub_index];
88  std::copy( corners, corners + num_sub_entity_nodes, sub_entity_conn );
89 
90  int sub_sub_corners[MAX_SUB_ENTITY_VERTICES];
91  int side, sense, offset;
92  for( int dim = 1; dim <= sub_dimension; ++dim )
93  {
94  if( !( ho_bits & ( 1 << dim ) ) ) continue;
95 
96  const short num_mid = NumSubEntities( subentity_topo, dim );
97  for( int i = 0; i < num_mid; ++i )
98  {
99  const EntityType sub_sub_topo = SubEntityType( subentity_topo, dim, i );
100  const int sub_sub_num_vert = VerticesPerEntity( sub_sub_topo );
101  SubEntityVertexIndices( subentity_topo, dim, i, sub_sub_corners );
102 
103  for( int j = 0; j < sub_sub_num_vert; ++j )
104  sub_sub_corners[j] = corners[sub_sub_corners[j]];
105  SideNumber( this_topo, sub_sub_corners, sub_sub_num_vert, dim, side, sense, offset );
106  sub_entity_conn[num_sub_entity_nodes++] = HONodeIndex( this_topo, num_nodes, dim, side );
107  }
108  }
109 }
110 
111 //! return the vertices of the specified sub entity
112 //! \param parent_conn Connectivity of parent entity
113 //! \param parent_type Entity type of parent entity
114 //! \param sub_dimension Dimension of sub-entity being queried
115 //! \param sub_index Index of sub-entity being queried
116 //! \param sub_entity_conn Connectivity of sub-entity, based on parent_conn and canonical
117 //! ordering for parent_type
118 //! \param num_sub_vertices Number of vertices in sub-entity
119 void CN::SubEntityConn( const void* parent_conn,
120  const EntityType parent_type,
121  const int sub_dimension,
122  const int sub_index,
123  void* sub_entity_conn,
124  int& num_sub_vertices )
125 {
126  static int sub_indices[MAX_SUB_ENTITY_VERTICES];
127 
128  SubEntityVertexIndices( parent_type, sub_dimension, sub_index, sub_indices );
129 
130  num_sub_vertices = VerticesPerEntity( SubEntityType( parent_type, sub_dimension, sub_index ) );
131  void** parent_conn_ptr = static_cast< void** >( const_cast< void* >( parent_conn ) );
132  void** sub_conn_ptr = static_cast< void** >( sub_entity_conn );
133  for( int i = 0; i < num_sub_vertices; i++ )
134  sub_conn_ptr[i] = parent_conn_ptr[sub_indices[i]];
135 }
136 
137 //! given an entity and a target dimension & side number, get that entity
138 short int CN::AdjacentSubEntities( const EntityType this_type,
139  const int* source_indices,
140  const int num_source_indices,
141  const int source_dim,
142  const int target_dim,
143  std::vector< int >& index_list,
144  const int operation_type )
145 {
146  // first get all the vertex indices
147  std::vector< int > tmp_indices;
148  const int* it1 = source_indices;
149 
150  assert( source_dim <= 3 && target_dim >= 0 && target_dim <= 3 &&
151  // make sure we're not stepping off the end of the array;
152  ( ( source_dim > 0 && *it1 < mConnectivityMap[this_type][source_dim - 1].num_sub_elements ) ||
153  ( source_dim == 0 &&
154  *it1 < mConnectivityMap[this_type][Dimension( this_type ) - 1].num_corners_per_sub_element[0] ) ) &&
155  *it1 >= 0 );
156 
157 #define MUC CN::mUpConnMap[this_type][source_dim][target_dim]
158 
159  // if we're looking for the vertices of a single side, return them in
160  // the canonical ordering; otherwise, return them in sorted order
161  if( num_source_indices == 1 && 0 == target_dim && source_dim != target_dim )
162  {
163 
164  // element of mConnectivityMap should be for this type and for one
165  // less than source_dim, which should give the connectivity of that sub element
166  const ConnMap& cm = mConnectivityMap[this_type][source_dim - 1];
167  std::copy( cm.conn[source_indices[0]],
168  cm.conn[source_indices[0]] + cm.num_corners_per_sub_element[source_indices[0]],
169  std::back_inserter( index_list ) );
170  return 0;
171  }
172 
173  // now go through source indices, folding adjacencies into target list
174  for( it1 = source_indices; it1 != source_indices + num_source_indices; it1++ )
175  {
176  // *it1 is the side index
177  // at start of iteration, index_list has the target list
178 
179  // if a union, or first iteration and index list was empty, copy the list
180  if( operation_type == CN::UNION || ( it1 == source_indices && index_list.empty() ) )
181  {
182  std::copy( MUC.targets_per_source_element[*it1],
183  MUC.targets_per_source_element[*it1] + MUC.num_targets_per_source_element[*it1],
184  std::back_inserter( index_list ) );
185  }
186  else
187  {
188  // else we're intersecting, and have a non-empty list; intersect with this target list
189  tmp_indices.clear();
190  for( int i = MUC.num_targets_per_source_element[*it1] - 1; i >= 0; i-- )
191  if( std::find( index_list.begin(), index_list.end(), MUC.targets_per_source_element[*it1][i] ) !=
192  index_list.end() )
193  tmp_indices.push_back( MUC.targets_per_source_element[*it1][i] );
194  // std::set_intersection(MUC.targets_per_source_element[*it1],
195  // MUC.targets_per_source_element[*it1]+
196  // MUC.num_targets_per_source_element[*it1],
197  // index_list.begin(), index_list.end(),
198  // std::back_inserter(tmp_indices));
199  index_list.swap( tmp_indices );
200 
201  // if we're at this point and the list is empty, the intersection will be NULL;
202  // return if so
203  if( index_list.empty() ) return 0;
204  }
205  }
206 
207  if( operation_type == CN::UNION && num_source_indices != 1 )
208  {
209  // need to sort then unique the list
210  std::sort( index_list.begin(), index_list.end() );
211  index_list.erase( std::unique( index_list.begin(), index_list.end() ), index_list.end() );
212  }
213 
214  return 0;
215 }
216 
217 template < typename T >
218 static short int side_number( const T* parent_conn,
219  const EntityType parent_type,
220  const T* child_conn,
221  const int child_num_verts,
222  const int child_dim,
223  int& side_no,
224  int& sense,
225  int& offset )
226 {
227  int parent_num_verts = CN::VerticesPerEntity( parent_type );
228  int side_indices[8];
229  assert( sizeof( side_indices ) / sizeof( side_indices[0] ) >= (size_t)child_num_verts );
230 
231  for( int i = 0; i < child_num_verts; i++ )
232  {
233  side_indices[i] = std::find( parent_conn, parent_conn + parent_num_verts, child_conn[i] ) - parent_conn;
234  if( side_indices[i] == parent_num_verts ) return -1;
235  }
236 
237  return CN::SideNumber( parent_type, &side_indices[0], child_num_verts, child_dim, side_no, sense, offset );
238 }
239 
240 short int CN::SideNumber( const EntityType parent_type,
241  const int* parent_conn,
242  const int* child_conn,
243  const int child_num_verts,
244  const int child_dim,
245  int& side_no,
246  int& sense,
247  int& offset )
248 {
249  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
250 }
251 
252 short int CN::SideNumber( const EntityType parent_type,
253  const unsigned int* parent_conn,
254  const unsigned int* child_conn,
255  const int child_num_verts,
256  const int child_dim,
257  int& side_no,
258  int& sense,
259  int& offset )
260 {
261  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
262 }
263 short int CN::SideNumber( const EntityType parent_type,
264  const long* parent_conn,
265  const long* child_conn,
266  const int child_num_verts,
267  const int child_dim,
268  int& side_no,
269  int& sense,
270  int& offset )
271 {
272  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
273 }
274 short int CN::SideNumber( const EntityType parent_type,
275  const unsigned long* parent_conn,
276  const unsigned long* child_conn,
277  const int child_num_verts,
278  const int child_dim,
279  int& side_no,
280  int& sense,
281  int& offset )
282 {
283  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
284 }
285 
286 short int CN::SideNumber( const EntityType parent_type,
287  const unsigned long long* parent_conn,
288  const unsigned long long* child_conn,
289  const int child_num_verts,
290  const int child_dim,
291  int& side_no,
292  int& sense,
293  int& offset )
294 
295 {
296  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
297 }
298 
299 short int CN::SideNumber( const EntityType parent_type,
300  void* const* parent_conn,
301  void* const* child_conn,
302  const int child_num_verts,
303  const int child_dim,
304  int& side_no,
305  int& sense,
306  int& offset )
307 {
308  return side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, side_no, sense, offset );
309 }
310 
311 short int CN::SideNumber( const EntityType parent_type,
312  const int* child_conn_indices,
313  const int child_num_verts,
314  const int child_dim,
315  int& side_no,
316  int& sense,
317  int& offset )
318 {
319  int parent_dim = Dimension( parent_type );
320  int parent_num_verts = VerticesPerEntity( parent_type );
321 
322  // degenerate case (vertex), output == input
323  if( child_dim == 0 )
324  {
325  if( child_num_verts != 1 ) return -1;
326  side_no = *child_conn_indices;
327  sense = offset = 0;
328  }
329 
330  // given a parent and child element, find the corresponding side number
331 
332  // dim_diff should be -1, 0 or 1 (same dimension, one less dimension, two less, resp.)
333  if( child_dim > parent_dim || child_dim < 0 ) return -1;
334 
335  // different types of same dimension won't be the same
336  if( parent_dim == child_dim && parent_num_verts != child_num_verts )
337  {
338  side_no = -1;
339  sense = 0;
340  return 0;
341  }
342 
343  // loop over the sub-elements, comparing to child connectivity
344  int sub_conn_indices[10];
345  for( int i = 0; i < NumSubEntities( parent_type, child_dim ); i++ )
346  {
347  int sub_size = VerticesPerEntity( SubEntityType( parent_type, child_dim, i ) );
348  if( sub_size != child_num_verts ) continue;
349 
350  SubEntityVertexIndices( parent_type, child_dim, i, sub_conn_indices );
351  bool they_match = ConnectivityMatch( child_conn_indices, sub_conn_indices, sub_size, sense, offset );
352  if( they_match )
353  {
354  side_no = i;
355  return 0;
356  }
357  }
358 
359  // if we've gotten here, we don't match
360  side_no = -1;
361 
362  // return value is no success
363  return 1;
364 }
365 
366 //! return the dimension and index of the opposite side, given parent entity type and child
367 //! dimension and index. This function is only defined for certain types of parent/child types:
368 //! (Parent, Child dim->Opposite dim):
369 //! (Tri, 1->0), (Tri, 0->1), (Quad, 1->1), (Quad, 0->0),
370 //! (Tet, 2->0), (Tet, 1->1), (Tet, 0->2),
371 //! (Hex, 2->2), (Hex, 1->1)(diagonally across element), (Hex, 0->0) (diagonally across element)
372 //! All other parent types and child dimensions return an error.
373 //!
374 //! \param parent_type The type of parent element
375 //! \param child_type The type of child element
376 //! \param child_index The index of the child element
377 //! \param opposite_index The index of the opposite element
378 //! \return status Returns 0 if successful, -1 if not
379 short int CN::OppositeSide( const EntityType parent_type,
380  const int child_index,
381  const int child_dim,
382  int& opposite_index,
383  int& opposite_dim )
384 {
385  switch( parent_type )
386  {
387  case MBEDGE:
388  if( 0 != child_dim )
389  return -1;
390  else
391  opposite_index = 1 - child_index;
392  opposite_dim = 0;
393  break;
394 
395  case MBTRI:
396  switch( child_dim )
397  {
398  case 0:
399  opposite_dim = 1;
400  opposite_index = ( child_index + 1 ) % 3;
401  break;
402  case 1:
403  opposite_dim = 0;
404  opposite_index = ( child_index + 2 ) % 3;
405  break;
406  default:
407  return -1;
408  }
409  break;
410 
411  case MBQUAD:
412  switch( child_dim )
413  {
414  case 0:
415  case 1:
416  opposite_dim = child_dim;
417  opposite_index = ( child_index + 2 ) % 4;
418  break;
419  default:
420  return -1;
421  }
422  break;
423 
424  case MBTET:
425  switch( child_dim )
426  {
427  case 0:
428  opposite_dim = 2;
429  opposite_index = ( child_index + 1 ) % 3 + 2 * ( child_index / 3 );
430  break;
431  case 1:
432  opposite_dim = 1;
433  opposite_index = child_index < 3 ? 3 + ( child_index + 2 ) % 3 : ( child_index + 1 ) % 3;
434  break;
435  case 2:
436  opposite_dim = 0;
437  opposite_index = ( child_index + 2 ) % 3 + child_index / 3;
438  break;
439  default:
440  return -1;
441  }
442  break;
443  case MBHEX:
444  opposite_dim = child_dim;
445  switch( child_dim )
446  {
447  case 0:
448  opposite_index = child_index < 4 ? 4 + ( child_index + 2 ) % 4 : ( child_index - 2 ) % 4;
449  break;
450  case 1:
451  opposite_index = 4 * ( 2 - child_index / 4 ) + ( child_index + 2 ) % 4;
452  break;
453  case 2:
454  opposite_index = child_index < 4 ? ( child_index + 2 ) % 4 : 9 - child_index;
455  break;
456  default:
457  return -1;
458  }
459  break;
460 
461  default:
462  return -1;
463  }
464 
465  return 0;
466 }
467 
468 template < typename T >
469 inline bool connectivity_match( const T* conn1_i, const T* conn2_i, const int num_vertices, int& direct, int& offset )
470 {
471 
472  bool they_match;
473 
474  // special test for 2 handles, since we don't want to wrap the list in this
475  // case
476  if( num_vertices == 2 )
477  {
478  they_match = false;
479  if( conn1_i[0] == conn2_i[0] && conn1_i[1] == conn2_i[1] )
480  {
481  direct = 1;
482  they_match = true;
483  offset = 0;
484  }
485  else if( conn1_i[0] == conn2_i[1] && conn1_i[1] == conn2_i[0] )
486  {
487  they_match = true;
488  direct = -1;
489  offset = 1;
490  }
491  }
492 
493  else
494  {
495  const T* iter;
496  iter = std::find( &conn2_i[0], &conn2_i[num_vertices], conn1_i[0] );
497  if( iter == &conn2_i[num_vertices] ) return false;
498 
499  they_match = true;
500 
501  offset = iter - conn2_i;
502  int i;
503 
504  // first compare forward
505  for( i = 1; i < num_vertices; ++i )
506  {
507  if( conn1_i[i] != conn2_i[( offset + i ) % num_vertices] )
508  {
509  they_match = false;
510  break;
511  }
512  }
513 
514  if( they_match == true )
515  {
516  direct = 1;
517  return they_match;
518  }
519 
520  they_match = true;
521 
522  // then compare reverse
523  for( i = 1; i < num_vertices; i++ )
524  {
525  if( conn1_i[i] != conn2_i[( offset + num_vertices - i ) % num_vertices] )
526  {
527  they_match = false;
528  break;
529  }
530  }
531  if( they_match )
532  {
533  direct = -1;
534  }
535  }
536 
537  return they_match;
538 }
539 
540 bool CN::ConnectivityMatch( const int* conn1_i, const int* conn2_i, const int num_vertices, int& direct, int& offset )
541 {
542  return connectivity_match< int >( conn1_i, conn2_i, num_vertices, direct, offset );
543 }
544 
545 bool CN::ConnectivityMatch( const unsigned int* conn1_i,
546  const unsigned int* conn2_i,
547  const int num_vertices,
548  int& direct,
549  int& offset )
550 {
551  return connectivity_match< unsigned int >( conn1_i, conn2_i, num_vertices, direct, offset );
552 }
553 
554 bool CN::ConnectivityMatch( const long* conn1_i, const long* conn2_i, const int num_vertices, int& direct, int& offset )
555 {
556  return connectivity_match< long >( conn1_i, conn2_i, num_vertices, direct, offset );
557 }
558 
559 bool CN::ConnectivityMatch( const unsigned long* conn1_i,
560  const unsigned long* conn2_i,
561  const int num_vertices,
562  int& direct,
563  int& offset )
564 {
565  return connectivity_match< unsigned long >( conn1_i, conn2_i, num_vertices, direct, offset );
566 }
567 
568 bool CN::ConnectivityMatch( const unsigned long long* conn1_i,
569  const unsigned long long* conn2_i,
570  const int num_vertices,
571  int& direct,
572  int& offset )
573 {
574  return connectivity_match< unsigned long long >( conn1_i, conn2_i, num_vertices, direct, offset );
575 }
576 
577 bool CN::ConnectivityMatch( void* const* conn1_i,
578  void* const* conn2_i,
579  const int num_vertices,
580  int& direct,
581  int& offset )
582 {
583  return connectivity_match< void* >( conn1_i, conn2_i, num_vertices, direct, offset );
584 }
585 
586 //! for an entity of this type and a specified subfacet (dimension and index), return
587 //! the index of the higher order node for that entity in this entity's connectivity array
588 short int CN::HONodeIndex( const EntityType this_type,
589  const int num_verts,
590  const int subfacet_dim,
591  const int subfacet_index )
592 {
593  int i;
594  int has_mids[4];
595  HasMidNodes( this_type, num_verts, has_mids );
596 
597  // if we have no mid nodes on the subfacet_dim, we have no index
598  if( subfacet_index != -1 && !has_mids[subfacet_dim] ) return -1;
599 
600  // put start index at last index (one less than the number of vertices
601  // plus the index basis)
602  int index = VerticesPerEntity( this_type ) - 1 + numberBasis;
603 
604  // for each subfacet dimension less than the target subfacet dim which has mid nodes,
605  // add the number of subfacets of that dimension to the index
606  for( i = 1; i < subfacet_dim; i++ )
607  if( has_mids[i] ) index += NumSubEntities( this_type, i );
608 
609  // now add the index of this subfacet, or one if we're asking about the entity as a whole
610  if( subfacet_index == -1 && has_mids[subfacet_dim] )
611  // want the index of the last ho node on this subfacet
612  index += NumSubEntities( this_type, subfacet_dim );
613 
614  else if( subfacet_index != -1 && has_mids[subfacet_dim] )
615  index += subfacet_index + 1 - numberBasis;
616 
617  // that's it
618  return index;
619 }
620 
621 //! given data about an element and a vertex in that element, return the dimension
622 //! and index of the sub-entity that the vertex resolves. If it does not resolve a
623 //! sub-entity, either because it's a corner node or it's not in the element, -1 is
624 //! returned in both return values
625 void CN::HONodeParent( EntityType elem_type, int num_verts, int ho_index, int& parent_dim, int& parent_index )
626 {
627  // begin with error values
628  parent_dim = parent_index = -1;
629 
630  // given the number of verts and the element type, get the hasmidnodes solution
631  int has_mids[4];
632  HasMidNodes( elem_type, num_verts, has_mids );
633 
634  int index = VerticesPerEntity( elem_type ) - 1;
635  const int dim = Dimension( elem_type );
636 
637  // keep a running sum of the ho node indices for this type of element, and stop
638  // when you get to the dimension which has the ho node
639  for( int i = 1; i < dim; i++ )
640  {
641  if( has_mids[i] )
642  {
643  if( ho_index <= index + NumSubEntities( elem_type, i ) )
644  {
645  // the ho_index resolves an entity of dimension i, so set the return values
646  // and break out of the loop
647  parent_dim = i;
648  parent_index = ho_index - index - 1;
649  return;
650  }
651  else
652  {
653  index += NumSubEntities( elem_type, i );
654  }
655  }
656  }
657 
658  // mid region node case
659  if( has_mids[dim] && ho_index == index + 1 )
660  {
661  parent_dim = dim;
662  parent_index = 0;
663  }
664 }
665 
666 const char* CN::EntityTypeName( const EntityType this_type )
667 {
668  return entityTypeNames[this_type];
669 }
670 
671 } // namespace moab
672 
673 using moab::CN;
674 using moab::EntityType;
675 
676 //! get the basis of the numbering system
677 void MBCN_GetBasis( int* rval )
678 {
679  *rval = CN::GetBasis();
680 }
681 
682 //! set the basis of the numbering system
683 void MBCN_SetBasis( const int in_basis )
684 {
685  CN::SetBasis( in_basis );
686 }
687 
688 //! return the string type name for this type
689 void MBCN_EntityTypeName( const int this_type, char* rval, int rval_len )
690 {
691  const char* rval_tmp = CN::EntityTypeName( (EntityType)this_type );
692  int rval_len_tmp = strlen( rval_tmp );
693  rval_len_tmp = ( rval_len_tmp < rval_len ? rval_len_tmp : rval_len );
694  strncpy( rval, rval_tmp, rval_len_tmp );
695 }
696 
697 //! given a name, find the corresponding entity type
698 void MBCN_EntityTypeFromName( const char* name, int* rval )
699 {
700  *rval = CN::EntityTypeFromName( name );
701 }
702 
703 //! return the topological entity dimension
704 void MBCN_Dimension( const int t, int* rval )
705 {
706  *rval = CN::Dimension( (EntityType)t );
707 }
708 
709 //! return the number of (corner) vertices contained in the specified type.
710 void MBCN_VerticesPerEntity( const int t, int* rval )
711 {
712  *rval = CN::VerticesPerEntity( (EntityType)t );
713 }
714 
715 //! return the number of sub-entities bounding the entity.
716 void MBCN_NumSubEntities( const int t, const int d, int* rval )
717 {
718  *rval = CN::NumSubEntities( (EntityType)t, d );
719 }
720 
721 //! return the type of a particular sub-entity.
722 //! \param this_type Type of entity for which sub-entity type is being queried
723 //! \param sub_dimension Topological dimension of sub-entity whose type is being queried
724 //! \param index Index of sub-entity whose type is being queried
725 //! \return type Entity type of sub-entity with specified dimension and index
726 void MBCN_SubEntityType( const int this_type, const int sub_dimension, const int index, int* rval )
727 
728 {
729 
730  *rval = CN::SubEntityType( (EntityType)this_type, sub_dimension, index );
731 }
732 
733 //! return the vertex indices of the specified sub-entity.
734 //! \param this_type Type of entity for which sub-entity connectivity is being queried
735 //! \param sub_dimension Dimension of sub-entity
736 //! \param sub_index Index of sub-entity
737 //! \param sub_entity_conn Connectivity of sub-entity (returned to calling function)
738 void MBCN_SubEntityVertexIndices( const int this_type,
739  const int sub_dimension,
740  const int sub_index,
741  int sub_entity_conn[] )
742 {
743  CN::SubEntityVertexIndices( (EntityType)this_type, sub_dimension, sub_index, sub_entity_conn );
744 }
745 
746 //! return the vertices of the specified sub entity
747 //! \param parent_conn Connectivity of parent entity
748 //! \param parent_type Entity type of parent entity
749 //! \param sub_dimension Dimension of sub-entity being queried
750 //! \param sub_index Index of sub-entity being queried
751 //! \param sub_entity_conn Connectivity of sub-entity, based on parent_conn and canonical
752 //! ordering for parent_type
753 //! \param num_sub_vertices Number of vertices in sub-entity
754 // void MBCN_SubEntityConn(const void *parent_conn, const int parent_type,
755 // const int sub_dimension,
756 // const int sub_index,
757 // void *sub_entity_conn, int &num_sub_vertices) {return
758 // CN::SubEntityConn();}
759 
760 //! For a specified set of sides of given dimension, return the intersection
761 //! or union of all sides of specified target dimension adjacent to those sides.
762 //! \param this_type Type of entity for which sub-entity connectivity is being queried
763 //! \param source_indices Indices of sides being queried
764 //! \param num_source_indices Number of entries in <em>source_indices</em>
765 //! \param source_dim Dimension of source entity
766 //! \param target_dim Dimension of target entity
767 //! \param index_list Indices of target entities (returned)
768 //! \param num_indices Number of indices of target entities (returned)
769 //! \param operation_type Specify either CN::INTERSECT (0) or CN::UNION (1) to get intersection
770 //! or union of target entity lists over source entities
771 //! \param rval Error code indicating success or failure (returned)
772 void MBCN_AdjacentSubEntities( const int this_type,
773  const int* source_indices,
774  const int num_source_indices,
775  const int source_dim,
776  const int target_dim,
777  int* index_list,
778  int* num_indices,
779  const int operation_type,
780  int* rval )
781 {
782  std::vector< int > tmp_index_list;
783  *rval = CN::AdjacentSubEntities( (EntityType)this_type, source_indices, num_source_indices, source_dim, target_dim,
784  tmp_index_list, operation_type );
785  std::copy( tmp_index_list.begin(), tmp_index_list.end(), index_list );
786  *num_indices = tmp_index_list.size();
787 }
788 
789 //! return the side index represented in the input sub-entity connectivity
790 //! \param parent_type Entity type of parent entity
791 //! \param child_conn_indices Child connectivity to query, specified as indices
792 //! into the connectivity list of the parent.
793 //! \param child_num_verts Number of values in <em>child_conn_indices</em>
794 //! \param child_dim Dimension of child entity being queried
795 //! \param side_no Side number of child entity (returned)
796 //! \param sense Sense of child entity with respect to order in <em>child_conn</em> (returned)
797 //! \param offset Offset of <em>child_conn</em> with respect to canonical ordering data (returned)
798 //! \return status Returns zero if successful, -1 if not
799 void MBCN_SideNumber( const int parent_type,
800  const int* child_conn_indices,
801  const int child_num_verts,
802  const int child_dim,
803  int* side_no,
804  int* sense,
805  int* offset )
806 {
807  CN::SideNumber( (EntityType)parent_type, child_conn_indices, child_num_verts, child_dim, *side_no, *sense,
808  *offset );
809 }
810 
811 void MBCN_SideNumberInt( const int* parent_conn,
812  const EntityType parent_type,
813  const int* child_conn,
814  const int child_num_verts,
815  const int child_dim,
816  int* side_no,
817  int* sense,
818  int* offset )
819 {
820  moab::side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, *side_no, *sense, *offset );
821 }
822 
823 void MBCN_SideNumberUint( const unsigned int* parent_conn,
824  const EntityType parent_type,
825  const unsigned int* child_conn,
826  const int child_num_verts,
827  const int child_dim,
828  int* side_no,
829  int* sense,
830  int* offset )
831 {
832  moab::side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, *side_no, *sense, *offset );
833 }
834 
835 void MBCN_SideNumberLong( const long* parent_conn,
836  const EntityType parent_type,
837  const long* child_conn,
838  const int child_num_verts,
839  const int child_dim,
840  int* side_no,
841  int* sense,
842  int* offset )
843 {
844  moab::side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, *side_no, *sense, *offset );
845 }
846 
847 void MBCN_SideNumberUlong( const unsigned long* parent_conn,
848  const EntityType parent_type,
849  const unsigned long* child_conn,
850  const int child_num_verts,
851  const int child_dim,
852  int* side_no,
853  int* sense,
854  int* offset )
855 {
856  moab::side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, *side_no, *sense, *offset );
857 }
858 
859 void MBCN_SideNumberVoid( void* const* parent_conn,
860  const EntityType parent_type,
861  void* const* child_conn,
862  const int child_num_verts,
863  const int child_dim,
864  int* side_no,
865  int* sense,
866  int* offset )
867 {
868  moab::side_number( parent_conn, parent_type, child_conn, child_num_verts, child_dim, *side_no, *sense, *offset );
869 }
870 
871 //! return the dimension and index of the opposite side, given parent entity type and child
872 //! dimension and index. This function is only defined for certain types of parent/child types:
873 //! (Parent, Child dim->Opposite dim):
874 //! (Tri, 1->0), (Tri, 0->1), (Quad, 1->1), (Quad, 0->0),
875 //! (Tet, 2->0), (Tet, 1->1), (Tet, 0->2),
876 //! (Hex, 2->2), (Hex, 1->1)(diagonally across element), (Hex, 0->0) (diagonally across element)
877 //! All other parent types and child dimensions return an error.
878 //!
879 //! \param parent_type The type of parent element
880 //! \param child_type The type of child element
881 //! \param child_index The index of the child element
882 //! \param opposite_index The index of the opposite element
883 //! \return status Returns 0 if successful, -1 if not
884 void MBCN_OppositeSide( const int parent_type,
885  const int child_index,
886  const int child_dim,
887  int* opposite_index,
888  int* opposite_dim,
889  int* rval )
890 {
891  *rval = CN::OppositeSide( (EntityType)parent_type, child_index, child_dim, *opposite_index, *opposite_dim );
892 }
893 
894 //! given two connectivity arrays, determine whether or not they represent the same entity.
895 //! \param conn1 Connectivity array of first entity
896 //! \param conn2 Connectivity array of second entity
897 //! \param num_vertices Number of entries in <em>conn1</em> and <em>conn2</em>
898 //! \param direct If positive, entities have the same sense (returned)
899 //! \param offset Offset of <em>conn2</em>'s first vertex in <em>conn1</em>
900 //! \return rval Returns true if <em>conn1</em> and <em>conn2</em> match
901 void MBCN_ConnectivityMatchInt( const int* conn1,
902  const int* conn2,
903  const int num_vertices,
904  int* direct,
905  int* offset,
906  int* rval )
907 {
908  *rval = CN::ConnectivityMatch( conn1, conn2, num_vertices, *direct, *offset );
909 }
910 
911 void MBCN_ConnectivityMatchUint( const unsigned int* conn1,
912  const unsigned int* conn2,
913  const int num_vertices,
914  int* direct,
915  int* offset,
916  int* rval )
917 {
918  *rval = CN::ConnectivityMatch( conn1, conn2, num_vertices, *direct, *offset );
919 }
920 
921 void MBCN_ConnectivityMatchLong( const long* conn1,
922  const long* conn2,
923  const int num_vertices,
924  int* direct,
925  int* offset,
926  int* rval )
927 {
928  *rval = CN::ConnectivityMatch( conn1, conn2, num_vertices, *direct, *offset );
929 }
930 
931 void MBCN_ConnectivityMatchUlong( const unsigned long* conn1,
932  const unsigned long* conn2,
933  const int num_vertices,
934  int* direct,
935  int* offset,
936  int* rval )
937 {
938  *rval = CN::ConnectivityMatch( conn1, conn2, num_vertices, *direct, *offset );
939 }
940 
941 void MBCN_ConnectivityMatchVoid( void* const* conn1,
942  void* const* conn2,
943  const int num_vertices,
944  int* direct,
945  int* offset,
946  int* rval )
947 {
948  *rval = CN::ConnectivityMatch( conn1, conn2, num_vertices, *direct, *offset );
949 }
950 
951 //! true if entities of a given type and number of nodes indicates mid edge nodes are present.
952 //! \param this_type Type of entity for which sub-entity connectivity is being queried
953 //! \param num_verts Number of nodes defining entity
954 //! \return int Returns true if <em>this_type</em> combined with <em>num_nodes</em> indicates
955 //! mid-edge nodes are likely
956 void MBCN_HasMidEdgeNodes( const int this_type, const int num_verts, int* rval )
957 {
958  *rval = CN::HasMidEdgeNodes( (EntityType)this_type, num_verts );
959 }
960 
961 //! true if entities of a given type and number of nodes indicates mid face nodes are present.
962 //! \param this_type Type of entity for which sub-entity connectivity is being queried
963 //! \param num_verts Number of nodes defining entity
964 //! \return int Returns true if <em>this_type</em> combined with <em>num_nodes</em> indicates
965 //! mid-face nodes are likely
966 void MBCN_HasMidFaceNodes( const int this_type, const int num_verts, int* rval )
967 {
968  *rval = CN::HasMidFaceNodes( (EntityType)this_type, num_verts );
969 }
970 
971 //! true if entities of a given type and number of nodes indicates mid region nodes are present.
972 //! \param this_type Type of entity for which sub-entity connectivity is being queried
973 //! \param num_verts Number of nodes defining entity
974 //! \return int Returns true if <em>this_type</em> combined with <em>num_nodes</em> indicates
975 //! mid-region nodes are likely
976 void MBCN_HasMidRegionNodes( const int this_type, const int num_verts, int* rval )
977 {
978  *rval = CN::HasMidRegionNodes( (EntityType)this_type, num_verts );
979 }
980 
981 //! true if entities of a given type and number of nodes indicates mid edge/face/region nodes
982 //! are present.
983 //! \param this_type Type of entity for which sub-entity connectivity is being queried
984 //! \param num_verts Number of nodes defining entity
985 //! \param mid_nodes If <em>mid_nodes[i], i=1..3</em> is true, indicates that mid-edge
986 //! (i=1), mid-face (i=2), and/or mid-region (i=3) nodes are likely
987 void MBCN_HasMidNodes( const int this_type, const int num_verts, int mid_nodes[4] )
988 {
989  return CN::HasMidNodes( (EntityType)this_type, num_verts, mid_nodes );
990 }
991 
992 //! given data about an element and a vertex in that element, return the dimension
993 //! and index of the sub-entity that the vertex resolves. If it does not resolve a
994 //! sub-entity, either because it's a corner node or it's not in the element, -1 is
995 //! returned in both return values.
996 //! \param elem_type Type of entity being queried
997 //! \param num_nodes The number of nodes in the element connectivity
998 //! \param ho_node_index The position of the HO node in the connectivity list (zero based)
999 //! \param parent_dim Dimension of sub-entity high-order node resolves (returned)
1000 //! \param parent_index Index of sub-entity high-order node resolves (returned)
1001 void MBCN_HONodeParent( int elem_type, int num_nodes, int ho_node_index, int* parent_dim, int* parent_index )
1002 {
1003  return CN::HONodeParent( (EntityType)elem_type, num_nodes, ho_node_index, *parent_dim, *parent_index );
1004 }
1005 
1006 //! for an entity of this type with num_verts vertices, and a specified subfacet
1007 //! (dimension and index), return the index of the higher order node for that entity
1008 //! in this entity's connectivity array
1009 //! \param this_type Type of entity being queried
1010 //! \param num_verts Number of vertices for the entity being queried
1011 //! \param subfacet_dim Dimension of sub-entity being queried
1012 //! \param subfacet_index Index of sub-entity being queried
1013 //! \return index Index of sub-entity's higher-order node
1014 void MBCN_HONodeIndex( const int this_type,
1015  const int num_verts,
1016  const int subfacet_dim,
1017  const int subfacet_index,
1018  int* rval )
1019 
1020 {
1021 
1022  *rval = CN::HONodeIndex( (EntityType)this_type, num_verts, subfacet_dim, subfacet_index );
1023 }
1024 
1025 namespace moab
1026 {
1027 
1028 template < typename T >
1029 inline int permute_this( EntityType t, const int dim, T* conn, const int indices_per_ent, const int num_entries )
1030 {
1031  T tmp_conn[MAX_SUB_ENTITIES];
1032  assert( indices_per_ent <= CN::permuteVec[t][dim][MAX_SUB_ENTITIES] );
1033  if( indices_per_ent > CN::permuteVec[t][dim][MAX_SUB_ENTITIES] ) return 1;
1034  short int* tvec = CN::permuteVec[t][dim];
1035  T* pvec = conn;
1036  for( int j = 0; j < num_entries; j++ )
1037  {
1038  for( int i = 0; i < indices_per_ent; i++ )
1039  tmp_conn[tvec[i]] = pvec[i];
1040  memcpy( pvec, tmp_conn, indices_per_ent * sizeof( T ) );
1041  pvec += indices_per_ent;
1042  }
1043 
1044  return 0;
1045 }
1046 
1047 template < typename T >
1048 inline int rev_permute_this( EntityType t, const int dim, T* conn, const int indices_per_ent, const int num_entries )
1049 {
1050  T tmp_conn[MAX_SUB_ENTITIES];
1051  assert( indices_per_ent <= CN::revPermuteVec[t][dim][MAX_SUB_ENTITIES] );
1052  if( indices_per_ent > CN::revPermuteVec[t][dim][MAX_SUB_ENTITIES] ) return 1;
1053  short int* tvec = CN::revPermuteVec[t][dim];
1054  T* pvec = conn;
1055  for( int j = 0; j < num_entries; j++ )
1056  {
1057  for( int i = 0; i < indices_per_ent; i++ )
1058  tmp_conn[i] = pvec[tvec[i]];
1059  memcpy( pvec, tmp_conn, indices_per_ent * sizeof( T ) );
1060  pvec += indices_per_ent;
1061  }
1062 
1063  return 0;
1064 }
1065 
1066 short int CN::Dimension( const EntityType t )
1067 {
1068  return mConnectivityMap[t][0].topo_dimension;
1069 }
1070 
1071 short int CN::VerticesPerEntity( const EntityType t )
1072 {
1073  return ( MBVERTEX == t
1074  ? (short int)1
1075  : mConnectivityMap[t][mConnectivityMap[t][0].topo_dimension - 1].num_corners_per_sub_element[0] );
1076 }
1077 
1078 short int CN::NumSubEntities( const EntityType t, const int d )
1079 {
1080  return ( t != MBVERTEX && d > 0 ? mConnectivityMap[t][d - 1].num_sub_elements
1081  : ( d ? (short int)-1 : VerticesPerEntity( t ) ) );
1082 }
1083 
1084 //! return the type of a particular sub-entity.
1085 EntityType CN::SubEntityType( const EntityType this_type, const int sub_dimension, const int index )
1086 {
1087  return ( !sub_dimension ? MBVERTEX
1088  : ( Dimension( this_type ) == sub_dimension && 0 == index
1089  ? this_type
1090  : mConnectivityMap[this_type][sub_dimension - 1].target_type[index] ) );
1091 }
1092 
1093 const short* CN::SubEntityVertexIndices( const EntityType this_type,
1094  const int sub_dimension,
1095  const int index,
1096  EntityType& sub_type,
1097  int& n )
1098 {
1099  if( sub_dimension == 0 )
1100  {
1101  n = 1;
1102  sub_type = MBVERTEX;
1103  return increasingInts + index;
1104  }
1105  else
1106  {
1107  const CN::ConnMap& map = mConnectivityMap[this_type][sub_dimension - 1];
1108  sub_type = map.target_type[index];
1109  n = map.num_corners_per_sub_element[index];
1110  return map.conn[index];
1111  }
1112 }
1113 
1114 //! Permute this vector
1115 inline int CN::permuteThis( const EntityType t, const int dim, int* pvec, const int num_indices, const int num_entries )
1116 {
1117  return permute_this( t, dim, pvec, num_indices, num_entries );
1118 }
1119 inline int CN::permuteThis( const EntityType t,
1120  const int dim,
1121  unsigned int* pvec,
1122  const int num_indices,
1123  const int num_entries )
1124 {
1125  return permute_this( t, dim, pvec, num_indices, num_entries );
1126 }
1127 inline int CN::permuteThis( const EntityType t,
1128  const int dim,
1129  long* pvec,
1130  const int num_indices,
1131  const int num_entries )
1132 {
1133  return permute_this( t, dim, pvec, num_indices, num_entries );
1134 }
1135 inline int CN::permuteThis( const EntityType t,
1136  const int dim,
1137  void** pvec,
1138  const int num_indices,
1139  const int num_entries )
1140 {
1141  return permute_this( t, dim, pvec, num_indices, num_entries );
1142 }
1143 
1144 //! Reverse permute this vector
1145 inline int CN::revPermuteThis( const EntityType t,
1146  const int dim,
1147  int* pvec,
1148  const int num_indices,
1149  const int num_entries )
1150 {
1151  return rev_permute_this( t, dim, pvec, num_indices, num_entries );
1152 }
1153 inline int CN::revPermuteThis( const EntityType t,
1154  const int dim,
1155  unsigned int* pvec,
1156  const int num_indices,
1157  const int num_entries )
1158 {
1159  return rev_permute_this( t, dim, pvec, num_indices, num_entries );
1160 }
1161 inline int CN::revPermuteThis( const EntityType t,
1162  const int dim,
1163  long* pvec,
1164  const int num_indices,
1165  const int num_entries )
1166 {
1167  return rev_permute_this( t, dim, pvec, num_indices, num_entries );
1168 }
1169 inline int CN::revPermuteThis( const EntityType t,
1170  const int dim,
1171  void** pvec,
1172  const int num_indices,
1173  const int num_entries )
1174 {
1175  return rev_permute_this( t, dim, pvec, num_indices, num_entries );
1176 }
1177 
1178 } // namespace moab