[MOAB-dev] r2658 - MOAB/trunk
kraftche at cae.wisc.edu
kraftche at cae.wisc.edu
Thu Feb 26 14:13:48 CST 2009
Author: kraftche
Date: 2009-02-26 14:13:48 -0600 (Thu, 26 Feb 2009)
New Revision: 2658
Modified:
MOAB/trunk/AEntityFactory.cpp
MOAB/trunk/HigherOrderFactory.cpp
MOAB/trunk/MBCN.cpp
MOAB/trunk/MBCN.hpp
MOAB/trunk/MBCNArrays.hpp
MOAB/trunk/MBSkinner.cpp
MOAB/trunk/mbcn_test.cc
Log:
o Fix bug getting number of edges, edge indices, etc. for MBKNIFE.
o Add new form of MBCN::HasMidNodes that returns an integer with bits set
for mid nodes at the corresponding dimension.
o Re-implement all forms of MBCN::HasMid*Nodes as table lookups.
o Rename MBCN::num_nodes_per_sub_element to MBCN::num_corners_per_sub_element
o Add new function: SubEntityNodeIndices. Similar to SubEntityVertexIndices,
but also returns any higher order nodes on the side.
o Add unit tests for SubEntityNodeIndices
Modified: MOAB/trunk/AEntityFactory.cpp
===================================================================
--- MOAB/trunk/AEntityFactory.cpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/AEntityFactory.cpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -806,7 +806,7 @@
const MBCN::ConnMap &cmap =
MBCN::mConnectivityMap[source_type][target_dimension-1];
- int verts_per_sub = cmap.num_nodes_per_sub_element[j];
+ int verts_per_sub = cmap.num_corners_per_sub_element[j];
// get the corner vertices
for (int i = 0; i < verts_per_sub; i++)
@@ -992,24 +992,27 @@
// shouldn't be here if source entity is not an element
assert(conn_len > 1);
- // create necessary entities
- if (create_if_missing) {
- rval = get_adjacency_ptr( conn[0], vtx_adj );
- if (MB_SUCCESS != rval)
- return rval;
- assert(vtx_adj != NULL); // should contain at least source_entity
+ // create necessary entities. this only makes sense if there exists of a
+ // dimension greater than the target dimension.
+ if (create_if_missing && target_dimension < 3 && MBCN::Dimension(src_type) < 2) {
+ for (size_t i = 0; i < conn_len; ++i) {
+ rval = get_adjacency_ptr( conn[i], vtx_adj );
+ if (MB_SUCCESS != rval)
+ return rval;
+ assert(vtx_adj != NULL); // should contain at least source_entity
- std::vector<MBEntityHandle> tmp2, tmp(*vtx_adj); // copy in case adjacency vector is changed
- for (size_t j = 0; j < tmp.size(); ++j) {
- if (MBCN::Dimension(TYPE_FROM_HANDLE(tmp[j])) <= (int)target_dimension)
- continue;
- if (TYPE_FROM_HANDLE(tmp[j]) == MBENTITYSET)
- break;
+ std::vector<MBEntityHandle> tmp2, tmp(*vtx_adj); // copy in case adjacency vector is changed
+ for (size_t j = 0; j < tmp.size(); ++j) {
+ if (MBCN::Dimension(TYPE_FROM_HANDLE(tmp[j])) <= (int)target_dimension)
+ continue;
+ if (TYPE_FROM_HANDLE(tmp[j]) == MBENTITYSET)
+ break;
- tmp2.clear();
- rval = get_down_adjacency_elements( tmp[j], target_dimension, tmp2, true, option );
- if (MB_SUCCESS != rval)
- return rval;
+ tmp2.clear();
+ rval = get_down_adjacency_elements( tmp[j], target_dimension, tmp2, true, option );
+ if (MB_SUCCESS != rval)
+ return rval;
+ }
}
}
Modified: MOAB/trunk/HigherOrderFactory.cpp
===================================================================
--- MOAB/trunk/HigherOrderFactory.cpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/HigherOrderFactory.cpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -343,7 +343,7 @@
tmp_face_conn[0] = element[entity_faces.conn[i][0]];
tmp_face_conn[1] = element[entity_faces.conn[i][1]];
tmp_face_conn[2] = element[entity_faces.conn[i][2]];
- if(entity_faces.num_nodes_per_sub_element[i] == 4)
+ if(entity_faces.num_corners_per_sub_element[i] == 4)
tmp_face_conn[3] = element[entity_faces.conn[i][3]];
else
tmp_face_conn[3] = 0;
@@ -359,7 +359,7 @@
{
EntitySequence* tmp_sequence = NULL;
double sum_coords[3] = {0,0,0};
- int max_nodes = entity_faces.num_nodes_per_sub_element[i];
+ int max_nodes = entity_faces.num_corners_per_sub_element[i];
for(int k=0; k<max_nodes; k++)
{
seq_manager->find(tmp_face_conn[k], tmp_sequence);
Modified: MOAB/trunk/MBCN.cpp
===================================================================
--- MOAB/trunk/MBCN.cpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/MBCN.cpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -65,6 +65,50 @@
return MBMAXTYPE;
}
+void MBCN::SubEntityNodeIndices( const MBEntityType this_topo,
+ const int num_nodes,
+ const int sub_dimension,
+ const int sub_index,
+ MBEntityType& subentity_topo,
+ int& num_sub_entity_nodes,
+ int sub_entity_conn[] )
+{
+ // If asked for a node, the special case...
+ if (sub_dimension == 0) {
+ assert( sub_index < num_nodes );
+ subentity_topo = MBVERTEX;
+ num_sub_entity_nodes = 1;
+ sub_entity_conn[0] = sub_index;
+ return;
+ }
+
+ const int ho_bits = HasMidNodes( this_topo, num_nodes );
+ subentity_topo = SubEntityType(this_topo, sub_dimension, sub_index);
+ num_sub_entity_nodes = VerticesPerEntity(subentity_topo);
+ const short* corners = mConnectivityMap[this_topo][sub_dimension-1].conn[sub_index];
+ std::copy( corners, corners+num_sub_entity_nodes, sub_entity_conn );
+
+ int sub_sub_corners[MB_MAX_SUB_ENTITY_VERTICES];
+ int side, sense, offset;
+ for (int dim = 1; dim <= sub_dimension; ++dim) {
+ if (!(ho_bits & (1<<dim)))
+ continue;
+
+ const short num_mid = NumSubEntities( subentity_topo, dim );
+ for (int i = 0; i < num_mid; ++i) {
+ const MBEntityType sub_sub_topo = SubEntityType( subentity_topo, dim, i );
+ const int sub_sub_num_vert = VerticesPerEntity( sub_sub_topo );
+ SubEntityVertexIndices( subentity_topo, dim, i, sub_sub_corners );
+
+ for (int j = 0; j < sub_sub_num_vert; ++j)
+ sub_sub_corners[j] = corners[sub_sub_corners[j]];
+ SideNumber( this_topo, sub_sub_corners, sub_sub_num_vert, dim, side, sense, offset );
+ sub_entity_conn[num_sub_entity_nodes++] = HONodeIndex( this_topo, num_nodes, dim, side );
+ }
+ }
+}
+
+
//! return the vertices of the specified sub entity
//! \param parent_conn Connectivity of parent entity
//! \param parent_type Entity type of parent entity
@@ -108,7 +152,7 @@
((source_dim > 0 &&
*it1 < mConnectivityMap[this_type][source_dim-1].num_sub_elements) ||
(source_dim == 0 &&
- *it1 < mConnectivityMap[this_type][Dimension(this_type)-1].num_nodes_per_sub_element[0])) &&
+ *it1 < mConnectivityMap[this_type][Dimension(this_type)-1].num_corners_per_sub_element[0])) &&
*it1 >= 0);
@@ -122,7 +166,7 @@
// less than source_dim, which should give the connectivity of that sub element
const ConnMap &cm = mConnectivityMap[this_type][source_dim-1];
std::copy(cm.conn[source_indices[0]],
- cm.conn[source_indices[0]]+cm.num_nodes_per_sub_element[source_indices[0]],
+ cm.conn[source_indices[0]]+cm.num_corners_per_sub_element[source_indices[0]],
std::back_inserter(index_list));
return 0;
}
Modified: MOAB/trunk/MBCN.hpp
===================================================================
--- MOAB/trunk/MBCN.hpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/MBCN.hpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -38,6 +38,7 @@
#include <vector>
#include <algorithm>
+#include <assert.h>
#include "MBEntityType.h"
@@ -67,6 +68,11 @@
public:
+ enum { MAX_NODES_PER_ELEMENT = 27 };
+ enum { MID_EDGE_BIT = 1<<1,
+ MID_FACE_BIT = 1<<2,
+ MID_REGION_BIT = 1<<3 };
+
//! enum used to specify operation type
enum {INTERSECT = 0, UNION};
@@ -82,7 +88,7 @@
short int num_sub_elements;
// Number of nodes in each sub-element of this dimension
- short int num_nodes_per_sub_element[MB_MAX_SUB_ENTITIES];
+ short int num_corners_per_sub_element[MB_MAX_SUB_ENTITIES];
// Type of each sub-element
MBEntityType target_type[MB_MAX_SUB_ENTITIES];
@@ -93,7 +99,7 @@
// mConnectivityMap[i=entity type][j=0,1,2]:
// num_sub_elements = # bounding edges(j=0) or faces(j=1) for entity type i, or self (j=2)
- // num_nodes_per_sub_element[k] (k=0..num_sub_elements-1) = number of nodes in sub-facet k
+ // num_corners_per_sub_element[k] (k=0..num_sub_elements-1) = number of nodes in sub-facet k
// (can vary over sub-facets, e.g. faces bounding a pyramid) or self (j=2)
// target_type[k] = entity type of sub-facet k (e.g. MBTRI or MBQUAD bounding a pyramid) or self (j=2)
// conn[k][l] (l=0..MBCN::VerticesPerEntity[target_type[k]]) = vertex connectivity of sub-facet k,
@@ -117,6 +123,9 @@
// efficient sorting)
// (not documented with Doxygen)
static const UpConnMap mUpConnMap[MBMAXTYPE][4][4];
+
+ // Mid-node bits indexed by number of nodes in element
+ static const unsigned char midNodesPerType[MBMAXTYPE][MAX_NODES_PER_ELEMENT+1];
//! Permutation and reverse permutation vectors
static short int permuteVec[MBMAXTYPE][3][MB_MAX_SUB_ENTITIES+1];
@@ -166,6 +175,22 @@
const int sub_dimension,
const int sub_index,
int sub_entity_conn[]);
+
+ //! return the node indices of the specified sub-entity.
+ //! \param this_topo The topology of the queried element type
+ //! \param num_nodes The number of nodes in the queried element type.
+ //! \param sub_dimension Dimension of sub-entity
+ //! \param sub_index Index of sub-entity
+ //! \param sub_entity_topo (Output) Topology of requested sub-entity.
+ //! \param num_sub_entity_nodes (Output) Number of nodes in the requested sub-entity.
+ //! \param sub_entity_conn (Output) Connectivity of sub-entity
+ static void SubEntityNodeIndices(const MBEntityType this_topo,
+ const int num_nodes,
+ const int sub_dimension,
+ const int sub_index,
+ MBEntityType& sub_entity_topo,
+ int& num_sub_entity_nodes,
+ int sub_entity_conn[]);
//! return the vertices of the specified sub entity
//! \param parent_conn Connectivity of parent entity
@@ -370,11 +395,16 @@
//! are present.
//! \param this_type Type of entity for which sub-entity connectivity is being queried
//! \param num_verts Number of nodes defining entity
- //! \param mid_nodes If <em>mid_nodes[i], i=1..2</em> is true, indicates that mid-edge
+ //! \param mid_nodes If <em>mid_nodes[i], i=1..2</em> is non-zero, indicates that mid-edge
//! (i=1), mid-face (i=2), and/or mid-region (i=3) nodes are likely
static void HasMidNodes(const MBEntityType this_type,
const int num_verts,
int mid_nodes[4]);
+ //! Same as above, except returns a single integer withthe bits, from
+ //! least significant to most significant set to one if the cooresponding
+ //! mid nodes on sub entities of the least dimension (0) to the highest
+ //! dimension (3) are present in the elment type.
+ static int HasMidNodes( const MBEntityType this_type, const int num_verts );
//! given data about an element and a vertex in that element, return the dimension
//! and index of the sub-entity that the vertex resolves. If it does not resolve a
@@ -418,7 +448,7 @@
inline short int MBCN::VerticesPerEntity(const MBEntityType t)
{
- return (MBVERTEX == t ? 1 : mConnectivityMap[t][mConnectivityMap[t][0].topo_dimension-1].num_nodes_per_sub_element[0]);
+ return (MBVERTEX == t ? 1 : mConnectivityMap[t][mConnectivityMap[t][0].topo_dimension-1].num_corners_per_sub_element[0]);
}
inline short int MBCN::NumSubEntities(const MBEntityType t, const int d)
@@ -454,75 +484,39 @@
inline bool MBCN::HasMidEdgeNodes(const MBEntityType this_type,
const int num_nodes)
{
- // poly elements never have mid nodes as far as canonical ordering is concerned
- if (MBPOLYGON == this_type || MBPOLYHEDRON == this_type) return false;
-
- if (num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1)) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) +
- NumSubEntities(this_type, 2)) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) +
- NumSubEntities(this_type, 2) + 1) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) + 1) )
- return true;
-
- else
- return false;
+ const int bits = HasMidNodes( this_type, num_nodes );
+ return static_cast<bool>( (bits & (1<<1)) >> 1 );
}
inline bool MBCN::HasMidFaceNodes(const MBEntityType this_type,
const int num_nodes)
{
- // poly elements never have mid nodes as far as canonical ordering is concerned
- if (MBPOLYGON == this_type || MBPOLYHEDRON == this_type) return false;
- // edges cannot have mid-face nodes
- if (Dimension(this_type) < 2) return false;
-
- if (num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 2)) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) +
- NumSubEntities(this_type, 2)) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) +
- NumSubEntities(this_type, 2) + 1) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 2) + 1) )
- return true;
-
- else
- return false;
+ const int bits = HasMidNodes( this_type, num_nodes );
+ return static_cast<bool>( (bits & (1<<2)) >> 2 );
}
inline bool MBCN::HasMidRegionNodes(const MBEntityType this_type,
const int num_nodes)
{
- // poly elements never have mid nodes as far as canonical ordering is concerned
- if (Dimension(this_type) != 3 || MBPOLYHEDRON == this_type) return false;
+ const int bits = HasMidNodes( this_type, num_nodes );
+ return static_cast<bool>( (bits & (1<<3)) >> 3 );
+}
- if (num_nodes == (VerticesPerEntity(this_type) + 1) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) + 1) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 1) +
- NumSubEntities(this_type, 2) + 1) ||
- num_nodes == (VerticesPerEntity(this_type) + NumSubEntities(this_type, 2) + 1) )
- return true;
-
- else
- return false;
+inline int MBCN::HasMidNodes( const MBEntityType this_type, const int num_nodes )
+{
+ assert( (unsigned)num_nodes <= (unsigned)MAX_NODES_PER_ELEMENT );
+ return midNodesPerType[this_type][num_nodes];
}
+
inline void MBCN::HasMidNodes(const MBEntityType this_type, const int num_nodes,
int mid_nodes[4])
{
- if (num_nodes <= VerticesPerEntity(this_type) ||
- // poly elements never have mid nodes as far as canonical ordering is concerned
- MBPOLYGON == this_type || MBPOLYHEDRON == this_type) {
- mid_nodes[0] = 0;
- mid_nodes[1] = 0;
- mid_nodes[2] = 0;
- mid_nodes[3] = 0;
- return;
- }
-
+ const int bits = HasMidNodes( this_type, num_nodes );
mid_nodes[0] = 0;
- mid_nodes[1] = HasMidEdgeNodes(this_type, num_nodes);
- mid_nodes[2] = HasMidFaceNodes(this_type, num_nodes);
- mid_nodes[3] = HasMidRegionNodes(this_type, num_nodes);
+ mid_nodes[1] = bits & MID_EDGE_BIT;
+ mid_nodes[2] = bits & MID_FACE_BIT;
+ mid_nodes[3] = bits & MID_REGION_BIT;
}
//! Set permutation or reverse permutation vector
Modified: MOAB/trunk/MBCNArrays.hpp
===================================================================
--- MOAB/trunk/MBCNArrays.hpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/MBCNArrays.hpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -18,8 +18,9 @@
const MBCN::ConnMap MBCN::mConnectivityMap[MBMAXTYPE][3] =
{
// vertex-edge
- {{ 0, 0 , {0}, {MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE,
- MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE}, {{0}} },
+ {{ 0, 0 , {0},
+ {MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE,
+ MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE}, {{0}} },
// vertex-face
{ 0, 0 , {0}, {MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE,
MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE}, {{0}} },
@@ -103,9 +104,9 @@
{{0,1,2,3,4,5}} }},
// knife-edge
- {{ 3, 8, {2,2,2,2,2,2,2,2}, {MBEDGE,MBEDGE,MBEDGE,MBEDGE,MBEDGE,MBEDGE,
- MBEDGE,MBEDGE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE},
- { {0,1},{1,2},{2,3},{3,0},{0,4},{1,5},{2,6},{3,5} } },
+ {{ 3, 10, {2,2,2,2,2,2,2,2,2,2}, {MBEDGE,MBEDGE,MBEDGE,MBEDGE,MBEDGE,MBEDGE,
+ MBEDGE,MBEDGE,MBEDGE,MBEDGE, MBMAXTYPE, MBMAXTYPE},
+ { {0,1},{1,2},{2,3},{3,0},{0,4},{1,5},{2,6},{3,5},{4,5},{5,6} } },
// knife-face
{ 3, 5, {4,4,4,4,4}, {MBQUAD,MBQUAD,MBQUAD,MBQUAD,MBQUAD, MBMAXTYPE,
MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE, MBMAXTYPE},
@@ -687,4 +688,34 @@
{{{0}, {{0}} }, {{0}, {{0}} }, {{0}, {{0}} }} // source dim 3
} // end type MBENTITYSET
};
+
+const unsigned char E = MBCN::MID_EDGE_BIT;
+const unsigned char F = MBCN::MID_FACE_BIT;
+const unsigned char R = MBCN::MID_REGION_BIT;
+const unsigned char MBCN::midNodesPerType[MBMAXTYPE][MAX_NODES_PER_ELEMENT+1] = {
+// vertex
+ { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// edge
+ { 0, 0, 0, E, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// tri
+ { 0, 0, 0, 0, F, 0, E, E|F, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// quad
+ { 0, 0, 0, 0, 0, F, 0, 0, E, E|F, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// polygon
+ { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// tet 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12,13, 14, 15
+ { 0, 0, 0, 0, 0, R, 0, 0, F, F|R, E, E|R, 0, 0, E|F, E|F|R, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// pyramid 3, 4, 5, 6, 7, 8, 9,10, 11, 12,13, 14, 15,16,17, 18, 19
+ { 0, 0, 0, 0, 0, 0, R, 0, 0, 0, F, F|R, 0, E, E|R, 0, 0, 0, E|F, E|F|R, 0, 0, 0, 0, 0, 0, 0, 0 },
+// prism 2, 3, 4, 5, 6, 7, 8, 9,10,11, 12, 13,14,15, 16, 17,18,19, 20, 21
+ { 0, 0, 0, 0, 0, 0, 0, R, 0, 0, 0, F, F|R, 0, 0, E, E|R, 0, 0, 0, E|F, E|F|R, 0, 0, 0, 0, 0, 0 },
+// knife 2, 3, 4, 5, 6, 7, 8, 9,10,11, 12,13, 14,15,16,17, 18, 19,20,21, 22, 23
+ { 0, 0, 0, 0, 0, 0, 0, 0, R, 0, 0, 0, F, F|R, 0, 0, 0, E, E|R, 0, 0, 0, E|F, E|F|R, 0, 0, 0, 0 },
+// hex 1, 2, 3, 4, 5, 6, 7, 8, 9, 10,11,12,13,14,15, 16,17,18,19,20, 21, 22,23,24,25, 26, 27
+ { 0, 0, 0, 0, 0, 0, 0, 0, 0, R, 0, 0, 0, 0, F, F|R, 0, 0, 0, 0, E, E|R, 0, 0, 0, 0, E|F, E|F|R },
+// polyhedron
+ { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+// set
+ { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 },
+ };
#endif
Modified: MOAB/trunk/MBSkinner.cpp
===================================================================
--- MOAB/trunk/MBSkinner.cpp 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/MBSkinner.cpp 2009-02-26 20:13:48 UTC (rev 2658)
@@ -302,7 +302,7 @@
const struct MBCN::ConnMap* conn_map = &(MBCN::mConnectivityMap[type][mTargetDim-1]);
for(int i=0; i<conn_map->num_sub_elements; i++)
{
- int num_sub_nodes = conn_map->num_nodes_per_sub_element[i];
+ int num_sub_nodes = conn_map->num_corners_per_sub_element[i];
assert(num_sub_nodes <= 32);
for(int j=0; j<num_sub_nodes; j++)
sub_conn[j] = conn[conn_map->conn[i][j]];
@@ -652,7 +652,7 @@
const struct MBCN::ConnMap* conn_map = &(MBCN::mConnectivityMap[type][0]);
for(int i=0; i<conn_map->num_sub_elements; i++)
{
- int num_sub_nodes = conn_map->num_nodes_per_sub_element[i];
+ int num_sub_nodes = conn_map->num_corners_per_sub_element[i];
assert(num_sub_nodes <= 32);
for(int j=0; j<num_sub_nodes; j++)
sub_conn[j] = conn[conn_map->conn[i][j]];
Modified: MOAB/trunk/mbcn_test.cc
===================================================================
--- MOAB/trunk/mbcn_test.cc 2009-02-25 19:39:48 UTC (rev 2657)
+++ MOAB/trunk/mbcn_test.cc 2009-02-26 20:13:48 UTC (rev 2658)
@@ -38,10 +38,34 @@
void test_opposite_side_tet();
void test_opposite_side_hex();
-void test_has_mid_nodes();
+void test_has_mid_nodes( MBEntityType type );
+void test_has_mid_nodes_edge() { test_has_mid_nodes(MBEDGE); }
+void test_has_mid_nodes_tri() { test_has_mid_nodes(MBTRI); }
+void test_has_mid_nodes_quad() { test_has_mid_nodes(MBQUAD); }
+void test_has_mid_nodes_tet() { test_has_mid_nodes(MBTET); }
+void test_has_mid_nodes_pyr() { test_has_mid_nodes(MBPYRAMID); }
+void test_has_mid_nodes_pri() { test_has_mid_nodes(MBPRISM); }
+void test_has_mid_nodes_knife(){ test_has_mid_nodes(MBKNIFE); }
+void test_has_mid_nodes_hex() { test_has_mid_nodes(MBHEX); }
+
void test_ho_node_parent();
void test_ho_node_index();
+void test_sub_entity_nodes( MBEntityType parent, int sub_dimension );
+void test_sub_entity_nodes( MBEntityType parent, int num_nodes, int sub_dimension );
+void test_sub_entity_nodes_tri_edges() { test_sub_entity_nodes(MBTRI, 1 ); }
+void test_sub_entity_nodes_quad_edges() { test_sub_entity_nodes(MBQUAD, 1 ); }
+void test_sub_entity_nodes_tet_edges() { test_sub_entity_nodes(MBTET, 1 ); }
+void test_sub_entity_nodes_tet_faces() { test_sub_entity_nodes(MBTET, 2 ); }
+void test_sub_entity_nodes_pyr_edges() { test_sub_entity_nodes(MBPYRAMID, 1 ); }
+void test_sub_entity_nodes_pyr_faces() { test_sub_entity_nodes(MBPYRAMID, 2 ); }
+void test_sub_entity_nodes_pri_edges() { test_sub_entity_nodes(MBPRISM, 1 ); }
+void test_sub_entity_nodes_pri_faces() { test_sub_entity_nodes(MBPRISM, 2 ); }
+void test_sub_entity_nodes_kni_edges() { test_sub_entity_nodes(MBKNIFE, 1 ); }
+void test_sub_entity_nodes_kni_faces() { test_sub_entity_nodes(MBKNIFE, 2 ); }
+void test_sub_entity_nodes_hex_edges() { test_sub_entity_nodes(MBHEX, 1 ); }
+void test_sub_entity_nodes_hex_faces() { test_sub_entity_nodes(MBHEX, 2 ); }
+
int main()
{
int result = 0;
@@ -82,7 +106,28 @@
result += RUN_TEST(test_opposite_side_tet);
result += RUN_TEST(test_opposite_side_hex);
- result += RUN_TEST(test_has_mid_nodes);
+ result += RUN_TEST(test_has_mid_nodes_edge);
+ result += RUN_TEST(test_has_mid_nodes_tri);
+ result += RUN_TEST(test_has_mid_nodes_quad);
+ result += RUN_TEST(test_has_mid_nodes_tet);
+ result += RUN_TEST(test_has_mid_nodes_pyr);
+ result += RUN_TEST(test_has_mid_nodes_pri);
+ result += RUN_TEST(test_has_mid_nodes_knife);
+ result += RUN_TEST(test_has_mid_nodes_hex);
+
+ result += RUN_TEST(test_sub_entity_nodes_tri_edges);
+ result += RUN_TEST(test_sub_entity_nodes_quad_edges);
+ result += RUN_TEST(test_sub_entity_nodes_tet_edges);
+ result += RUN_TEST(test_sub_entity_nodes_tet_faces);
+ result += RUN_TEST(test_sub_entity_nodes_pyr_edges);
+ result += RUN_TEST(test_sub_entity_nodes_pyr_faces);
+ result += RUN_TEST(test_sub_entity_nodes_pri_edges);
+ result += RUN_TEST(test_sub_entity_nodes_pri_faces);
+ result += RUN_TEST(test_sub_entity_nodes_kni_edges);
+ result += RUN_TEST(test_sub_entity_nodes_kni_faces);
+ result += RUN_TEST(test_sub_entity_nodes_hex_edges);
+ result += RUN_TEST(test_sub_entity_nodes_hex_faces);
+
result += RUN_TEST(test_ho_node_parent);
result += RUN_TEST(test_ho_node_index);
return result;
@@ -182,7 +227,7 @@
CHECK_EQUAL( 5, MBCN::NumSubEntities(MBPRISM, 2));
CHECK_EQUAL( 7, MBCN::NumSubEntities(MBKNIFE, 0));
- CHECK_EQUAL( 8, MBCN::NumSubEntities(MBKNIFE, 1));
+ CHECK_EQUAL(10, MBCN::NumSubEntities(MBKNIFE, 1));
CHECK_EQUAL( 5, MBCN::NumSubEntities(MBKNIFE, 2));
CHECK_EQUAL( 8, MBCN::NumSubEntities(MBHEX, 0));
@@ -822,7 +867,7 @@
CHECK_EQUAL( 4, idx );
}
-void test_has_mid_nodes()
+void test_has_mid_nodes(MBEntityType type)
{
const int combinations[][4] = { { 0, 0, 0, 0 },
{ 0, 1, 0, 0 },
@@ -833,8 +878,6 @@
{ 0, 0, 1, 1 },
{ 0, 1, 1, 1 } };
- for (const MBEntityType* t = elem_types; *t != MBMAXTYPE; ++t) {
- const MBEntityType type = *t;
const int dim = MBCN::Dimension(type);
// calculate number of valid combinations of ho node flags
int num_comb = 1;
@@ -857,12 +900,11 @@
int results[4] = { 0, -1, -1, -1 };
MBCN::HasMidNodes( type, num_vtx, results );
- CHECK_EQUAL( 0, results[0] );
- CHECK_EQUAL( ho_nodes[1], results[1] );
- CHECK_EQUAL( ho_nodes[2], results[2] );
- CHECK_EQUAL( ho_nodes[3], results[3] );
+ CHECK_EQUAL( 0, !!results[0] );
+ CHECK_EQUAL( ho_nodes[1], !!results[1] );
+ CHECK_EQUAL( ho_nodes[2], !!results[2] );
+ CHECK_EQUAL( ho_nodes[3], !!results[3] );
}
- }
}
void test_ho_node_parent()
@@ -985,3 +1027,100 @@
} // for each type
}
+void test_sub_entity_nodes( MBEntityType parent, int sub_dimension )
+{
+ const int num_corner = MBCN::VerticesPerEntity( parent );
+ const int num_edge = MBCN::NumSubEntities( parent, 1 );
+ const int num_face = MBCN::NumSubEntities( parent, 2 );
+
+ switch (MBCN::Dimension(parent)) {
+ case 3:
+ test_sub_entity_nodes( parent, num_corner+num_face, sub_dimension );
+ test_sub_entity_nodes( parent, num_corner+num_edge+num_face, sub_dimension );
+ test_sub_entity_nodes( parent, num_corner+num_face+1, sub_dimension );
+ test_sub_entity_nodes( parent, num_corner+num_edge+num_face+1, sub_dimension );
+ case 2:
+ test_sub_entity_nodes( parent, num_corner+num_edge, sub_dimension );
+ test_sub_entity_nodes( parent, num_corner+num_edge+1, sub_dimension );
+ case 1:
+ test_sub_entity_nodes( parent, num_corner, sub_dimension );
+ test_sub_entity_nodes( parent, num_corner+1, sub_dimension );
+ break;
+ default:
+ CHECK(false);
+ }
+}
+
+void test_sub_entity_nodes( MBEntityType parent, int num_nodes, int sub_dimension )
+{
+ const int num_sub = MBCN::NumSubEntities( parent, sub_dimension );
+ const int parent_ho = MBCN::HasMidNodes( parent, num_nodes );
+ int child_ho = 0;
+ for (int d = 1; d <= sub_dimension; ++d)
+ child_ho |= (parent_ho & (1<<d));
+
+ // first test the types
+ for (int i = 0; i < num_sub; ++i) {
+ int num, conn[MB_MAX_SUB_ENTITY_VERTICES];
+ MBEntityType type;
+ MBCN::SubEntityNodeIndices( parent, num_nodes, sub_dimension, i, type, num, conn );
+ CHECK_EQUAL( MBCN::SubEntityType(parent, sub_dimension, i), type );
+ }
+
+ // now test that they have the correct number of higher-order node
+ for (int i = 0; i < num_sub; ++i) {
+ int num, conn[MB_MAX_SUB_ENTITY_VERTICES];
+ MBEntityType type;
+ MBCN::SubEntityNodeIndices( parent, num_nodes, sub_dimension, i, type, num, conn );
+ const int ho = MBCN::HasMidNodes( type, num );
+ CHECK_EQUAL( child_ho, ho );
+ }
+
+ // now test the actual indices
+ for (int i = 0; i < num_sub; ++i) {
+ int num, conn[MB_MAX_SUB_ENTITY_VERTICES], corners[MB_MAX_SUB_ENTITY_VERTICES];
+ MBEntityType type;
+ MBCN::SubEntityNodeIndices( parent, num_nodes, sub_dimension, i, type, num, conn );
+
+ // check corner indices against SubEntityVertexIndices
+ const int num_corner = MBCN::VerticesPerEntity(type);
+ CHECK( num >= num_corner );
+ MBCN::SubEntityVertexIndices( parent, sub_dimension, i, corners );
+ for (int j = 0; j < num_corner; ++j)
+ CHECK_EQUAL( corners[j], conn[j] );
+
+ // check mid-edge indices, if present
+ int idx = num_corner;
+ if (child_ho & MBCN::MID_EDGE_BIT) {
+ // for each edge in the sub-entity type
+ const int num_edge = MBCN::NumSubEntities( type, 1 );
+ for (int j = 0; j < num_edge; ++j) {
+ // get edge indices for sub-entity connectivity
+ int edge_ends[2];
+ MBCN::SubEntityVertexIndices( type, 1, j, edge_ends );
+ // convert to indices into parent type's connectivity
+ CHECK( edge_ends[0] < num_corner );
+ edge_ends[0] = corners[edge_ends[0]];
+ CHECK( edge_ends[1] < num_corner );
+ edge_ends[1] = corners[edge_ends[1]];
+ // find edge index in parent element
+ int side, sense, off;
+ int result = MBCN::SideNumber( parent, edge_ends, 2, 1, side, sense, off );
+ CHECK_EQUAL( 0, result );
+ // get location in parent entity connectivity for mid-edge node
+ int loc = MBCN::HONodeIndex( parent, num_nodes, 1, side );
+ CHECK_EQUAL( loc, conn[idx++] );
+ }
+ }
+
+ // check mid-face indices, if present
+ if (child_ho & MBCN::MID_FACE_BIT) {
+ CHECK_EQUAL( 2, MBCN::Dimension(type) );
+ int loc = MBCN::HONodeIndex( parent, num_nodes, 2, i );
+ CHECK_EQUAL( loc, conn[idx++] );
+ }
+
+ // make sure there were no extra node indices returned
+ CHECK_EQUAL( idx, num );
+ }
+}
More information about the moab-dev
mailing list