MOAB: Mesh Oriented datABase  (version 5.4.1)
TopologyInfo.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2004 Lawrence Livermore National Laboratory.  Under
00005     the terms of Contract B545069 with the University of Wisconsin --
00006     Madison, Lawrence Livermore National Laboratory retains certain
00007     rights in this software.
00008 
00009     This library is free software; you can redistribute it and/or
00010     modify it under the terms of the GNU Lesser General Public
00011     License as published by the Free Software Foundation; either
00012     version 2.1 of the License, or (at your option) any later version.
00013 
00014     This library is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017     Lesser General Public License for more details.
00018 
00019     You should have received a copy of the GNU Lesser General Public License
00020     (lgpl.txt) along with this library; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 
00023     [email protected]
00024 
00025   ***************************************************************** */
00026 
00027 #include "MsqError.hpp"
00028 #include "TopologyInfo.hpp"
00029 
00030 #include <cstring>
00031 #include <cassert>
00032 
00033 namespace MBMesquite
00034 {
00035 
00036 TopologyInfo TopologyInfo::instance;
00037 
00038 const char long_polygon_name[]        = "Polygon";
00039 const char long_triangle_name[]       = "Triangle";
00040 const char long_quadrilateral_name[]  = "Quadrilateral";
00041 const char long_polyhedron_name[]     = "Polyhedron";
00042 const char long_tetrahedron_name[]    = "Tetrahedron";
00043 const char long_hexahedron_name[]     = "Hexahedron";
00044 const char long_prism_name[]          = "Prism";
00045 const char long_pyramid_name[]        = "Pyramd";
00046 const char long_septahedron_name[]    = "Septahedron";
00047 const char short_polygon_name[]       = "Polygon";
00048 const char short_triangle_name[]      = "Tri";
00049 const char short_quadrilateral_name[] = "Quad";
00050 const char short_polyhedron_name[]    = "Polyhedron";
00051 const char short_tetrahedron_name[]   = "Tet";
00052 const char short_hexahedron_name[]    = "Hex";
00053 const char short_prism_name[]         = "Pri";
00054 const char short_pyramid_name[]       = "Pyr";
00055 const char short_septahedron_name[]   = "Sept";
00056 
00057 TopologyInfo::TopologyInfo()
00058 {
00059     memset( dimMap, 0, sizeof( dimMap ) );
00060     memset( adjMap, 0, sizeof( adjMap ) );
00061     memset( edgeMap, 0, sizeof( edgeMap ) );
00062     memset( faceMap, 0, sizeof( faceMap ) );
00063     memset( vertAdjMap, 0, sizeof( vertAdjMap ) );
00064     memset( shortNames, 0, sizeof( shortNames ) );
00065     memset( longNames, 0, sizeof( longNames ) );
00066 
00067     longNames[POLYGON]       = long_polygon_name;
00068     longNames[TRIANGLE]      = long_triangle_name;
00069     longNames[QUADRILATERAL] = long_quadrilateral_name;
00070     longNames[POLYHEDRON]    = long_polyhedron_name;
00071     longNames[TETRAHEDRON]   = long_tetrahedron_name;
00072     longNames[HEXAHEDRON]    = long_hexahedron_name;
00073     longNames[PRISM]         = long_prism_name;
00074     longNames[PYRAMID]       = long_pyramid_name;
00075     longNames[SEPTAHEDRON]   = long_septahedron_name;
00076 
00077     shortNames[POLYGON]       = short_polygon_name;
00078     shortNames[TRIANGLE]      = short_triangle_name;
00079     shortNames[QUADRILATERAL] = short_quadrilateral_name;
00080     shortNames[POLYHEDRON]    = short_polyhedron_name;
00081     shortNames[TETRAHEDRON]   = short_tetrahedron_name;
00082     shortNames[HEXAHEDRON]    = short_hexahedron_name;
00083     shortNames[PRISM]         = short_prism_name;
00084     shortNames[PYRAMID]       = short_pyramid_name;
00085     shortNames[SEPTAHEDRON]   = short_septahedron_name;
00086 
00087     dimMap[POLYGON]       = 2;
00088     dimMap[TRIANGLE]      = 2;
00089     dimMap[QUADRILATERAL] = 2;
00090     dimMap[POLYHEDRON]    = 3;
00091     dimMap[TETRAHEDRON]   = 3;
00092     dimMap[HEXAHEDRON]    = 3;
00093     dimMap[PRISM]         = 3;
00094     dimMap[PYRAMID]       = 3;
00095     dimMap[SEPTAHEDRON]   = 3;
00096 
00097     adjMap[TRIANGLE][0] = 3;
00098     adjMap[TRIANGLE][1] = 3;
00099     adjMap[TRIANGLE][2] = 1;
00100     adjMap[TRIANGLE][3] = 0;
00101 
00102     adjMap[QUADRILATERAL][0] = 4;
00103     adjMap[QUADRILATERAL][1] = 4;
00104     adjMap[QUADRILATERAL][2] = 1;
00105     adjMap[QUADRILATERAL][3] = 0;
00106 
00107     adjMap[TETRAHEDRON][0] = 4;
00108     adjMap[TETRAHEDRON][1] = 6;
00109     adjMap[TETRAHEDRON][2] = 4;
00110     adjMap[TETRAHEDRON][3] = 1;
00111 
00112     adjMap[HEXAHEDRON][0] = 8;
00113     adjMap[HEXAHEDRON][1] = 12;
00114     adjMap[HEXAHEDRON][2] = 6;
00115     adjMap[HEXAHEDRON][3] = 1;
00116 
00117     adjMap[PRISM][0] = 6;
00118     adjMap[PRISM][1] = 9;
00119     adjMap[PRISM][2] = 5;
00120     adjMap[PRISM][3] = 1;
00121 
00122     adjMap[PYRAMID][0] = 5;
00123     adjMap[PYRAMID][1] = 8;
00124     adjMap[PYRAMID][2] = 5;
00125     adjMap[PYRAMID][3] = 1;
00126 
00127     adjMap[SEPTAHEDRON][0] = 7;
00128     adjMap[SEPTAHEDRON][1] = 11;
00129     adjMap[SEPTAHEDRON][2] = 6; /* See description in TSTT mesh interface doc */
00130     adjMap[SEPTAHEDRON][3] = 1;
00131 
00132     int side;
00133     for( side = 0; side < 3; ++side )
00134     {
00135         edgeMap[TRIANGLE - FIRST_FACE][side][0] = side;
00136         edgeMap[TRIANGLE - FIRST_FACE][side][1] = ( side + 1 ) % 3;
00137     }
00138     for( side = 0; side < 4; ++side )
00139     {
00140         edgeMap[QUADRILATERAL - FIRST_FACE][side][0] = side;
00141         edgeMap[QUADRILATERAL - FIRST_FACE][side][1] = ( side + 1 ) % 4;
00142     }
00143     for( side = 0; side < 3; ++side )
00144     {
00145         edgeMap[TETRAHEDRON - FIRST_FACE][side][0] = side;
00146         edgeMap[TETRAHEDRON - FIRST_FACE][side][1] = ( side + 1 ) % 3;
00147     }
00148     for( side = 3; side < 6; ++side )
00149     {
00150         edgeMap[TETRAHEDRON - FIRST_FACE][side][0] = side - 3;
00151         edgeMap[TETRAHEDRON - FIRST_FACE][side][1] = 3;
00152     }
00153     for( side = 0; side < 4; ++side )
00154     {
00155         edgeMap[HEXAHEDRON - FIRST_FACE][side][0] = side;
00156         edgeMap[HEXAHEDRON - FIRST_FACE][side][1] = ( side + 1 ) % 4;
00157     }
00158     for( side = 4; side < 8; ++side )
00159     {
00160         edgeMap[HEXAHEDRON - FIRST_FACE][side][0] = side - 4;
00161         edgeMap[HEXAHEDRON - FIRST_FACE][side][1] = side;
00162     }
00163     for( side = 8; side < 12; ++side )
00164     {
00165         edgeMap[HEXAHEDRON - FIRST_FACE][side][0] = side - 4;
00166         edgeMap[HEXAHEDRON - FIRST_FACE][side][1] = 4 + ( side + 1 ) % 4;
00167     }
00168     for( side = 0; side < 3; ++side )
00169     {
00170         edgeMap[PRISM - FIRST_FACE][side][0] = side;
00171         edgeMap[PRISM - FIRST_FACE][side][1] = ( side + 1 ) % 3;
00172     }
00173     for( side = 3; side < 6; ++side )
00174     {
00175         edgeMap[PRISM - FIRST_FACE][side][0] = side - 3;
00176         edgeMap[PRISM - FIRST_FACE][side][1] = side;
00177     }
00178     for( side = 6; side < 9; ++side )
00179     {
00180         edgeMap[PRISM - FIRST_FACE][side][0] = side - 3;
00181         edgeMap[PRISM - FIRST_FACE][side][1] = 3 + ( side + 1 ) % 3;
00182     }
00183     for( side = 0; side < 4; ++side )
00184     {
00185         edgeMap[PYRAMID - FIRST_FACE][side][0] = side;
00186         edgeMap[PYRAMID - FIRST_FACE][side][1] = ( side + 1 ) % 4;
00187     }
00188     for( side = 4; side < 8; ++side )
00189     {
00190         edgeMap[PYRAMID - FIRST_FACE][side][0] = side - 4;
00191         edgeMap[PYRAMID - FIRST_FACE][side][1] = 4;
00192     }
00193 
00194     for( side = 0; side < 3; ++side )
00195     {
00196         faceMap[TETRAHEDRON - FIRST_VOL][side][0] = 3;
00197         faceMap[TETRAHEDRON - FIRST_VOL][side][1] = side;
00198         faceMap[TETRAHEDRON - FIRST_VOL][side][2] = ( side + 1 ) % 3;
00199         faceMap[TETRAHEDRON - FIRST_VOL][side][3] = 3;
00200     }
00201     faceMap[TETRAHEDRON - FIRST_VOL][3][0] = 3;
00202     faceMap[TETRAHEDRON - FIRST_VOL][3][1] = 2;
00203     faceMap[TETRAHEDRON - FIRST_VOL][3][2] = 1;
00204     faceMap[TETRAHEDRON - FIRST_VOL][3][3] = 0;
00205 
00206     for( side = 0; side < 4; ++side )
00207     {
00208         faceMap[HEXAHEDRON - FIRST_VOL][side][0] = 4;
00209         faceMap[HEXAHEDRON - FIRST_VOL][side][1] = side;
00210         faceMap[HEXAHEDRON - FIRST_VOL][side][2] = ( side + 1 ) % 4;
00211         faceMap[HEXAHEDRON - FIRST_VOL][side][3] = 4 + ( side + 1 ) % 4;
00212         faceMap[HEXAHEDRON - FIRST_VOL][side][4] = side + 4;
00213     }
00214     faceMap[HEXAHEDRON - FIRST_VOL][4][0] = 4;
00215     faceMap[HEXAHEDRON - FIRST_VOL][4][1] = 3;
00216     faceMap[HEXAHEDRON - FIRST_VOL][4][2] = 2;
00217     faceMap[HEXAHEDRON - FIRST_VOL][4][3] = 1;
00218     faceMap[HEXAHEDRON - FIRST_VOL][4][4] = 0;
00219     faceMap[HEXAHEDRON - FIRST_VOL][5][0] = 4;
00220     faceMap[HEXAHEDRON - FIRST_VOL][5][1] = 4;
00221     faceMap[HEXAHEDRON - FIRST_VOL][5][2] = 5;
00222     faceMap[HEXAHEDRON - FIRST_VOL][5][3] = 6;
00223     faceMap[HEXAHEDRON - FIRST_VOL][5][4] = 7;
00224 
00225     for( side = 0; side < 4; ++side )
00226     {
00227         faceMap[PYRAMID - FIRST_VOL][side][0] = 3;
00228         faceMap[PYRAMID - FIRST_VOL][side][1] = side;
00229         faceMap[PYRAMID - FIRST_VOL][side][2] = ( side + 1 ) % 4;
00230         faceMap[PYRAMID - FIRST_VOL][side][3] = 4;
00231     }
00232     faceMap[PYRAMID - FIRST_VOL][4][0] = 4;
00233     faceMap[PYRAMID - FIRST_VOL][4][1] = 3;
00234     faceMap[PYRAMID - FIRST_VOL][4][2] = 2;
00235     faceMap[PYRAMID - FIRST_VOL][4][3] = 1;
00236     faceMap[PYRAMID - FIRST_VOL][4][4] = 0;
00237 
00238     for( side = 0; side < 3; ++side )
00239     {
00240         faceMap[PRISM - FIRST_VOL][side][0] = 4;
00241         faceMap[PRISM - FIRST_VOL][side][1] = side;
00242         faceMap[PRISM - FIRST_VOL][side][2] = ( side + 1 ) % 3;
00243         faceMap[PRISM - FIRST_VOL][side][3] = 3 + ( side + 1 ) % 3;
00244         faceMap[PRISM - FIRST_VOL][side][4] = side + 3;
00245     }
00246     faceMap[PRISM - FIRST_VOL][3][0] = 3;
00247     faceMap[PRISM - FIRST_VOL][3][1] = 2;
00248     faceMap[PRISM - FIRST_VOL][3][2] = 1;
00249     faceMap[PRISM - FIRST_VOL][3][3] = 0;
00250     faceMap[PRISM - FIRST_VOL][4][0] = 3;
00251     faceMap[PRISM - FIRST_VOL][4][1] = 3;
00252     faceMap[PRISM - FIRST_VOL][4][2] = 4;
00253     faceMap[PRISM - FIRST_VOL][4][3] = 5;
00254 
00255     int i;
00256     for( i = 0; i < 3; ++i )
00257     {
00258         vertAdjMap[TRIANGLE - FIRST_FACE][i][0] = 2;
00259         vertAdjMap[TRIANGLE - FIRST_FACE][i][1] = ( i + 1 ) % 3;
00260         vertAdjMap[TRIANGLE - FIRST_FACE][i][2] = ( i + 2 ) % 3;
00261     }
00262 
00263     for( i = 0; i < 4; ++i )
00264     {
00265         vertAdjMap[QUADRILATERAL - FIRST_FACE][i][0] = 2;
00266         vertAdjMap[QUADRILATERAL - FIRST_FACE][i][1] = ( i + 1 ) % 4;
00267         vertAdjMap[QUADRILATERAL - FIRST_FACE][i][2] = ( i + 3 ) % 4;
00268     }
00269 
00270     unsigned tet_corner_data[] = { 1, 2, 3, 0, 3, 2, 3, 0, 1, 2, 1, 0 };
00271     for( i = 0; i < 4; ++i )
00272     {
00273         vertAdjMap[TETRAHEDRON - FIRST_FACE][i][0] = 3;
00274         for( unsigned j = 0; j < 3; ++j )
00275             vertAdjMap[TETRAHEDRON - FIRST_FACE][i][j + 1] = tet_corner_data[3 * i + j];
00276     }
00277 
00278     for( i = 0; i < 4; ++i )
00279     {
00280         vertAdjMap[PYRAMID - FIRST_FACE][i][0] = 3;
00281         vertAdjMap[PYRAMID - FIRST_FACE][i][1] = ( i + 1 ) % 4;
00282         vertAdjMap[PYRAMID - FIRST_FACE][i][2] = ( i + 3 ) % 4;
00283         vertAdjMap[PYRAMID - FIRST_FACE][i][3] = 4;
00284     }
00285     vertAdjMap[PYRAMID - FIRST_FACE][4][0] = 4;
00286     for( i = 0; i < 4; i++ )
00287         vertAdjMap[PYRAMID - FIRST_FACE][4][i + 1] = 3 - i;
00288 
00289     for( i = 0; i < 4; ++i )
00290     {
00291         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][0] = 3;
00292         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][1] = ( i + 1 ) % 4;
00293         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][2] = ( i + 3 ) % 4;
00294         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][3] = i + 4;
00295     }
00296     for( i = 4; i < 8; ++i )
00297     {
00298         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][0] = 3;
00299         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][1] = ( i + 3 ) % 4 + 4;
00300         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][2] = ( i + 1 ) % 4 + 4;
00301         vertAdjMap[HEXAHEDRON - FIRST_FACE][i][3] = i - 4;
00302     }
00303 
00304     for( i = 0; i < 3; ++i )
00305     {
00306         vertAdjMap[PRISM - FIRST_FACE][i][0] = 3;
00307         vertAdjMap[PRISM - FIRST_FACE][i][1] = ( i + 1 ) % 3;
00308         vertAdjMap[PRISM - FIRST_FACE][i][2] = ( i + 2 ) % 3;
00309         vertAdjMap[PRISM - FIRST_FACE][i][3] = i + 3;
00310     }
00311     for( i = 3; i < 6; ++i )
00312     {
00313         vertAdjMap[PRISM - FIRST_FACE][i][0] = 3;
00314         vertAdjMap[PRISM - FIRST_FACE][i][1] = ( i + 2 ) % 3 + 3;
00315         vertAdjMap[PRISM - FIRST_FACE][i][2] = ( i + 1 ) % 3 + 3;
00316         vertAdjMap[PRISM - FIRST_FACE][i][3] = i - 3;
00317     }
00318 
00319     // Build reverse vertex-vertex adjacency index map
00320     const EntityTopology types[] = { TRIANGLE, QUADRILATERAL, TETRAHEDRON, PYRAMID, PRISM, HEXAHEDRON };
00321     const int num_types          = sizeof( types ) / sizeof( types[0] );
00322     for( i = 0; i < num_types; ++i )
00323     {
00324         const unsigned num_vert = corners( types[i] );
00325         for( unsigned v = 0; v < num_vert; ++v )
00326         {
00327             unsigned num_v_adj;
00328             const unsigned* v_adj = adjacent_vertices( types[i], v, num_v_adj );
00329             unsigned* reverse     = revVertAdjIdx[types[i] - FIRST_FACE][v];
00330             reverse[0]            = num_v_adj;
00331 
00332             for( unsigned j = 0; j < num_v_adj; ++j )
00333             {
00334                 unsigned num_j_adj, k;
00335                 const unsigned* j_adj = adjacent_vertices( types[i], v_adj[j], num_j_adj );
00336                 for( k = 0; k < num_j_adj && j_adj[k] != v; ++k )
00337                     ;
00338                 assert( k < num_j_adj );  // If this fails, vertAdjMap is corrupt!
00339                 reverse[j + 1] = k;
00340             }
00341         }
00342     }
00343 }
00344 
00345 void TopologyInfo::higher_order( EntityTopology topo,
00346                                  unsigned num_nodes,
00347                                  bool& midedge,
00348                                  bool& midface,
00349                                  bool& midvol,
00350                                  MsqError& err )
00351 {
00352     int ho  = higher_order( topo, num_nodes, err );
00353     midedge = (bool)( ( ho & ( 1 << 1 ) ) >> 1 );
00354     midface = (bool)( ( ho & ( 1 << 2 ) ) >> 2 );
00355     midvol  = (bool)( ( ho & ( 1 << 3 ) ) >> 3 );
00356 }
00357 
00358 int TopologyInfo::higher_order( EntityTopology topo, unsigned num_nodes, MsqError& err )
00359 {
00360     int result = 0;
00361     if( topo == POLYGON )  // polygons currently do not have higher order elements
00362         return 0;
00363 
00364     if( topo >= MIXED || num_nodes < instance.adjMap[topo][0] )
00365     {
00366         MSQ_SETERR( err )( "Invalid element topology", MsqError::INVALID_ARG );
00367         return 0;
00368     }
00369 
00370     unsigned dim = instance.dimMap[topo];
00371     assert( num_nodes >= instance.adjMap[topo][0] );
00372     unsigned nodes = num_nodes - instance.adjMap[topo][0];
00373     unsigned edges = instance.adjMap[topo][1];
00374     unsigned faces = instance.adjMap[topo][2];
00375     if( edges && nodes >= edges )
00376     {
00377         nodes -= edges;
00378         result |= 1 << 1;
00379     }
00380     if( faces && nodes >= faces )
00381     {
00382         nodes -= faces;
00383         result |= 1 << 2;
00384     }
00385     if( 1 == nodes )
00386     {
00387         if( 2 == dim )
00388         {
00389             nodes -= 1;
00390             result |= 1 << 2;
00391         }
00392         else if( 3 == dim )
00393         {
00394             nodes -= 1;
00395             result |= 1 << 3;
00396         }
00397     }
00398 
00399     if( nodes )
00400     {
00401         MSQ_SETERR( err )( "Invalid element topology", MsqError::INVALID_STATE );
00402     }
00403 
00404     return result;
00405 }
00406 
00407 int TopologyInfo::higher_order_from_side( EntityTopology topo,
00408                                           unsigned num_nodes,
00409                                           unsigned side_dimension,
00410                                           unsigned side_number,
00411                                           MsqError& err )
00412 {
00413     bool mids[4] = { true };
00414     higher_order( topo, num_nodes, mids[1], mids[2], mids[3], err );
00415     MSQ_ERRZERO( err );
00416 
00417     if( side_dimension > dimension( topo ) || side_number > adjacent( topo, side_dimension ) )
00418     {
00419         MSQ_SETERR( err )( MsqError::INVALID_ARG, "Invalid side number: %u\n", side_number );
00420         return 0;
00421     }
00422 
00423     if( !mids[side_dimension] ) return -1;
00424 
00425     int result = side_number;
00426     switch( side_dimension )
00427     {
00428         case 3:
00429             if( mids[2] ) result += faces( topo );
00430         case 2:
00431             if( mids[1] ) result += edges( topo );
00432         case 1:
00433             result += corners( topo );
00434         case 0:
00435             break;
00436         default:
00437             MSQ_SETERR( err )( MsqError::INVALID_ARG, "Invalid dimension: %u\n", side_dimension );
00438             return 0;
00439     }
00440     return result;
00441 }
00442 
00443 void TopologyInfo::side_from_higher_order( EntityTopology topo,
00444                                            unsigned num_nodes,
00445                                            unsigned node_number,
00446                                            unsigned& side_dim_out,
00447                                            unsigned& side_num_out,
00448                                            MsqError& err )
00449 {
00450     bool midedge, midface, midvol;
00451     higher_order( topo, num_nodes, midedge, midface, midvol, err );MSQ_ERRRTN( err );
00452     side_num_out = node_number;
00453 
00454     if( side_num_out < corners( topo ) )
00455     {
00456         side_dim_out = 0;
00457         return;
00458     }
00459     side_num_out -= corners( topo );
00460 
00461     if( midedge )
00462     {
00463         if( side_num_out < edges( topo ) )
00464         {
00465             side_dim_out = 1;
00466             return;
00467         }
00468         side_num_out -= edges( topo );
00469     }
00470 
00471     if( midface )
00472     {
00473         if( side_num_out < faces( topo ) )
00474         {
00475             side_dim_out = 2;
00476             return;
00477         }
00478         side_num_out -= faces( topo );
00479     }
00480 
00481     if( midvol && side_num_out == 0 )
00482     {
00483         side_dim_out = 3;
00484         return;
00485     }
00486 
00487     MSQ_SETERR( err )( MsqError::INVALID_ARG, "Invalid node index\n" );
00488 }
00489 
00490 const unsigned* TopologyInfo::edge_vertices( EntityTopology topo, unsigned edge, MsqError& err )
00491 {
00492     if( topo < (EntityTopology)FIRST_FACE || topo > (EntityTopology)LAST_VOL || edge >= edges( topo ) )
00493     {
00494         MSQ_SETERR( err )( MsqError::INVALID_ARG );
00495         topo = (EntityTopology)FIRST_FACE;
00496         edge = 0;
00497     }
00498 
00499     return instance.edgeMap[topo - FIRST_FACE][edge];
00500 }
00501 
00502 const unsigned* TopologyInfo::edge_vertices( EntityTopology topo, unsigned edge )
00503 {
00504     if( topo < (EntityTopology)FIRST_FACE || topo > (EntityTopology)LAST_VOL || edge >= edges( topo ) )
00505     {
00506         return 0;
00507     }
00508     return instance.edgeMap[topo - FIRST_FACE][edge];
00509 }
00510 
00511 const unsigned* TopologyInfo::face_vertices( EntityTopology topo, unsigned face, unsigned& length, MsqError& err )
00512 {
00513     if( topo < (EntityTopology)FIRST_VOL || topo > (EntityTopology)LAST_VOL || face >= faces( topo ) )
00514     {
00515         MSQ_SETERR( err )( MsqError::INVALID_ARG );
00516         topo = (EntityTopology)FIRST_VOL;
00517         face = 0;
00518     }
00519 
00520     length = instance.faceMap[topo - FIRST_VOL][face][0];
00521     return instance.faceMap[topo - FIRST_VOL][face] + 1;
00522 }
00523 const unsigned* TopologyInfo::face_vertices( EntityTopology topo, unsigned face, unsigned& length )
00524 {
00525     if( topo < (EntityTopology)FIRST_VOL || topo > (EntityTopology)LAST_VOL || face >= faces( topo ) )
00526     {
00527         return 0;
00528     }
00529 
00530     length = instance.faceMap[topo - FIRST_VOL][face][0];
00531     return instance.faceMap[topo - FIRST_VOL][face] + 1;
00532 }
00533 
00534 const unsigned* TopologyInfo::side_vertices( EntityTopology topo,
00535                                              unsigned dim,
00536                                              unsigned side,
00537                                              unsigned& count_out,
00538                                              MsqError& err )
00539 {
00540     static const unsigned all[] = { 0, 1, 2, 3, 4, 5, 6, 7 };
00541     const unsigned* result;
00542 
00543     if( dim != 0 && dim == dimension( topo ) )
00544     {
00545         count_out = corners( topo );
00546         result    = all;
00547     }
00548     else if( dim == 1 )
00549     {
00550         count_out = 2;
00551         result    = edge_vertices( topo, side, err );
00552     }
00553     else if( dim == 2 )
00554     {
00555         result = face_vertices( topo, side, count_out, err );
00556     }
00557     else
00558     {
00559         MSQ_SETERR( err )( MsqError::INVALID_ARG );
00560         count_out = 0;
00561         result    = 0;
00562     }
00563     return result;
00564 }
00565 const unsigned* TopologyInfo::side_vertices( EntityTopology topo, unsigned dim, unsigned side, unsigned& count_out )
00566 {
00567     static const unsigned all[] = { 0, 1, 2, 3, 4, 5, 6, 7 };
00568     const unsigned* result;
00569 
00570     if( dim != 0 && dim == dimension( topo ) )
00571     {
00572         count_out = corners( topo );
00573         result    = all;
00574     }
00575     else if( dim == 1 )
00576     {
00577         count_out = 2;
00578         result    = edge_vertices( topo, side );
00579     }
00580     else if( dim == 2 )
00581     {
00582         result = face_vertices( topo, side, count_out );
00583     }
00584     else
00585     {
00586         result = 0;
00587     }
00588     return result;
00589 }
00590 
00591 void TopologyInfo::side_number( EntityTopology topo,
00592                                 unsigned num_nodes,
00593                                 unsigned node_index,
00594                                 unsigned& side_dim_out,
00595                                 unsigned& side_num_out,
00596                                 MsqError& err )
00597 {
00598     if( topo >= (EntityTopology)MIXED || num_nodes < instance.adjMap[topo][0] )
00599     {
00600         MSQ_SETERR( err )( "Invalid element topology", MsqError::INVALID_ARG );
00601         return;
00602     }
00603 
00604     unsigned nodes = instance.adjMap[topo][0];
00605     unsigned edges = instance.adjMap[topo][1];
00606     unsigned faces = instance.adjMap[topo][2];
00607     side_num_out   = node_index;
00608 
00609     if( side_num_out < nodes )
00610     {
00611         side_dim_out = 0;
00612         return;
00613     }
00614     num_nodes -= nodes;
00615     side_num_out -= nodes;
00616 
00617     if( edges && num_nodes >= edges )
00618     {
00619         if( side_num_out < edges )
00620         {
00621             side_dim_out = 1;
00622             return;
00623         }
00624         num_nodes -= edges;
00625         side_num_out -= edges;
00626     }
00627     if( faces && num_nodes >= faces )
00628     {
00629         if( side_num_out < faces )
00630         {
00631             side_dim_out = 2;
00632             return;
00633         }
00634         num_nodes -= faces;
00635         side_num_out -= faces;
00636     }
00637     if( side_num_out == 0 )
00638     {
00639         side_dim_out = instance.dimMap[topo];
00640         side_num_out = 0;
00641         return;
00642     }
00643 
00644     MSQ_SETERR( err )( MsqError::INVALID_ARG );
00645 }
00646 
00647 const unsigned* TopologyInfo::adjacent_vertices( EntityTopology topo, unsigned index, unsigned& num_adj_out )
00648 {
00649     const unsigned count = corners( topo );
00650     if( !count || index >= count )
00651     {
00652         num_adj_out = 0;
00653         return 0;
00654     }
00655 
00656     const unsigned* vect = instance.vertAdjMap[topo - FIRST_FACE][index];
00657     num_adj_out          = vect[0];
00658     return vect + 1;
00659 }
00660 
00661 const unsigned* TopologyInfo::reverse_vertex_adjacency_offsets( EntityTopology topo,
00662                                                                 unsigned index,
00663                                                                 unsigned& num_adj_out )
00664 {
00665     const unsigned count = corners( topo );
00666     if( !count || index >= count )
00667     {
00668         num_adj_out = 0;
00669         return 0;
00670     }
00671 
00672     const unsigned* vect = instance.revVertAdjIdx[topo - FIRST_FACE][index];
00673     num_adj_out          = vect[0];
00674     return vect + 1;
00675 }
00676 
00677 bool TopologyInfo::compare_sides( const size_t* verts1,
00678                                   EntityTopology type1,
00679                                   unsigned side1,
00680                                   const size_t* verts2,
00681                                   EntityTopology type2,
00682                                   unsigned side2,
00683                                   unsigned side_dim,
00684                                   MsqError& err )
00685 {
00686     const unsigned *conn1, *conn2;
00687     unsigned len1, len2;
00688 
00689     conn1 = side_vertices( type1, side_dim, side1, len1, err );
00690     MSQ_ERRZERO( err );
00691     conn2 = side_vertices( type2, side_dim, side2, len2, err );
00692     MSQ_ERRZERO( err );
00693 
00694     // obviously not the same if different number of vertices
00695     // (triangular face cannot match quadrilateral face)
00696     if( len1 != len2 ) return false;
00697 
00698     // Find location (i) in vertices of side of second element
00699     // that matches the first vertex in the side of the first
00700     // element.
00701     unsigned i, j;
00702     for( i = 0; i < len2; ++i )
00703         if( verts1[conn1[0]] == verts2[conn2[i]] ) break;
00704     // If not found, then no match
00705     if( i == len2 ) return false;
00706 
00707     // Try comparing side connectivity in forward order
00708     for( j = 1; j < len1; ++j )
00709         if( verts1[conn1[j]] != verts2[conn2[( i + j ) % len2]] ) break;
00710     // If they match, we're done
00711     if( j == len1 ) return true;
00712 
00713     // Try comparing in reverse order
00714     for( j = 1; j < len1; ++j )
00715         if( verts1[conn1[j]] != verts2[conn2[( i + len2 - j ) % len2]] ) return false;
00716     // If here, matched in reverse order
00717     return true;
00718 }
00719 
00720 unsigned TopologyInfo::find_edge( EntityTopology topo,
00721                                   const unsigned* side_vertices,
00722                                   bool& reversed_out,
00723                                   MsqError& err )
00724 {
00725     if( dimension( topo ) <= 1 )
00726     {
00727         MSQ_SETERR( err )( MsqError::INVALID_ARG, "Invalid element dimension" );
00728         return (unsigned)-1;
00729     }
00730 
00731     for( unsigned i = 0; i < edges( topo ); ++i )
00732     {
00733         const unsigned* edge = edge_vertices( topo, i, err );
00734         MSQ_ERRZERO( err );
00735 
00736         if( edge[0] == side_vertices[0] && edge[1] == side_vertices[1] )
00737         {
00738             reversed_out = false;
00739             return i;
00740         }
00741 
00742         if( edge[0] == side_vertices[1] && edge[1] == side_vertices[0] )
00743         {
00744             reversed_out = true;
00745             return i;
00746         }
00747     }
00748 
00749     MSQ_SETERR( err )( MsqError::INVALID_ARG, "No such edge" );
00750     return (unsigned)-1;
00751 }
00752 
00753 unsigned TopologyInfo::find_face( EntityTopology topo,
00754                                   const unsigned* side_vertices,
00755                                   unsigned num_vertices,
00756                                   bool& reversed_out,
00757                                   MsqError& err )
00758 {
00759     if( dimension( topo ) <= 2 )
00760     {
00761         MSQ_SETERR( err )( MsqError::INVALID_ARG, "Invalid element dimension" );
00762         return (unsigned)-1;
00763     }
00764 
00765     for( unsigned i = 0; i < faces( topo ); ++i )
00766     {
00767         unsigned j, n, offset;
00768         const unsigned* face = face_vertices( topo, i, n, err );
00769         MSQ_ERRZERO( err );
00770         if( n != num_vertices ) continue;
00771 
00772         for( offset = 0; offset < num_vertices; ++offset )
00773             if( face[offset] == side_vertices[0] ) break;
00774         if( offset == num_vertices ) continue;
00775 
00776         for( j = 1; j < num_vertices; ++j )
00777             if( side_vertices[j] != face[( offset + j ) % num_vertices] ) break;
00778         if( j == num_vertices )
00779         {
00780             reversed_out = false;
00781             return i;
00782         }
00783 
00784         for( j = 1; j < num_vertices; ++j )
00785             if( side_vertices[j] != face[( offset + num_vertices - j ) % num_vertices] ) break;
00786         if( j == num_vertices )
00787         {
00788             reversed_out = true;
00789             return i;
00790         }
00791     }
00792 
00793     MSQ_SETERR( err )( MsqError::INVALID_ARG, "No such face" );
00794     return (unsigned)-1;
00795 }
00796 
00797 void TopologyInfo::find_side( EntityTopology topo,
00798                               const unsigned* side_vertices,
00799                               unsigned num_vertices,
00800                               unsigned& dimension_out,
00801                               unsigned& number_out,
00802                               bool& reversed_out,
00803                               MsqError& err )
00804 {
00805     switch( num_vertices )
00806     {
00807         case 1:
00808             dimension_out = 0;
00809             number_out    = *side_vertices;
00810             reversed_out  = false;
00811             if( *side_vertices >= corners( topo ) ) MSQ_SETERR( err )
00812             ( MsqError::INVALID_ARG, "Invalid corner number: %u\n", *side_vertices );
00813             break;
00814         case 2:
00815             dimension_out = 1;
00816             number_out    = find_edge( topo, side_vertices, reversed_out, err );MSQ_CHKERR( err );
00817             break;
00818         case 3:
00819         case 4:
00820             dimension_out = 2;
00821             number_out    = find_face( topo, side_vertices, num_vertices, reversed_out, err );MSQ_CHKERR( err );
00822             break;
00823         default:
00824             MSQ_SETERR( err )
00825             ( MsqError::UNSUPPORTED_ELEMENT, "Invalid number of side vertices: %u\n", num_vertices );
00826             break;
00827     }
00828 }
00829 
00830 }  // namespace MBMesquite
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines