[MOAB-dev] commit/MOAB: 3 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Mon Jun 30 22:28:47 CDT 2014
3 new commits in MOAB:
https://bitbucket.org/fathomteam/moab/commits/dbf0a5098d31/
Changeset: dbf0a5098d31
Branch: None
User: nray
Date: 2014-07-01 05:28:15
Summary: Added the HalfFacetRep class and one test code.
Affected #: 7 files
diff --git a/MeshFiles/unittest/Makefile.am b/MeshFiles/unittest/Makefile.am
index bc4f522..2bd2a8a 100644
--- a/MeshFiles/unittest/Makefile.am
+++ b/MeshFiles/unittest/Makefile.am
@@ -34,5 +34,6 @@ EXTRA_DIST = 125hex.g \
BedCrop2.h5m \
mpas_p8.h5m \
Homme_2pt.h5m \
- surfrandomtris-4part.h5m
+ surfrandomtris-4part.h5m \
+ hexes_mixed.vtk
diff --git a/MeshFiles/unittest/hexes_mixed.vtk b/MeshFiles/unittest/hexes_mixed.vtk
new file mode 100755
index 0000000..8ee5c11
--- /dev/null
+++ b/MeshFiles/unittest/hexes_mixed.vtk
@@ -0,0 +1,45 @@
+# vtk DataFile Version 3.0
+This file was written using writevtk_unstr.m
+ASCII
+DATASET UNSTRUCTURED_GRID
+POINTS 18 double
+0 -1 0
+0 1 0
+1 -1 0
+1 1 0
+-1 1 0
+-1 -1 0
+0 -1 1
+0 1 1
+1 -1 1
+1 1 1
+-1 1 1
+-1 -1 1
+0 -1 -1
+0 1 -1
+1 -1 -1
+1 1 -1
+-1 1 -1
+-1 -1 -1
+
+CELLS 9 59
+2 0 1
+4 0 1 7 6
+4 0 2 3 1
+4 0 12 13 1
+4 0 1 4 5
+8 0 2 3 1 6 8 9 7
+8 5 0 1 4 11 6 7 10
+8 12 14 15 13 0 2 3 1
+8 17 12 13 16 5 0 1 4
+
+CELL_TYPES 9
+3
+9
+9
+9
+9
+12
+12
+12
+12
This diff is so big that we needed to truncate the remainder.
https://bitbucket.org/fathomteam/moab/commits/8fa155080b49/
Changeset: 8fa155080b49
Branch: None
User: nray
Date: 2014-07-01 05:28:15
Summary: Added interface to AHF functionalities through the get_adjacencies
functions.
Affected #: 10 files
diff --git a/configure.ac b/configure.ac
index 5de30e3..87ad9ee 100644
--- a/configure.ac
+++ b/configure.ac
@@ -657,6 +657,13 @@ MB_OPTIONAL_TOOL([mbdepth], [yes])
MB_OPTIONAL_TOOL([refiner], [no])
MB_OPTIONAL_TOOL([h5mtools], [yes])
MB_OPTIONAL_TOOL([mbcslam], [no])
+MB_OPTIONAL_TOOL([ahf], [no])
+
+
+if test "xyes" = "x$ENABLE_ahf"; then
+ AM_CPPFLAGS="$AM_CPPFLAGS -DUSE_AHF"
+fi
+
if test "xyes" = "x$ENABLE_refiner"; then
if test "xyes" != "x$WITH_MPI"; then
@@ -1333,3 +1340,27 @@ if test "x$WARN_PARALLEL_HDF5_NO_COMPLEX" = "xyes"; then
*************************************************************************])
fi
+################################################################################
+# AHF OPTIONS
+################################################################################
+#AC_ARG_ENABLE([ahf],
+#[AC_HELP_STRING([--enable-ahf],[use the half facet representation])],
+#[ case "${enableval}" in
+# yes) ahf1=true ;;
+# no) ahf1=false ;;
+# *) AC_MSG_ERROR([bad value]) ;;
+#esac],[ahf1=false])
+ #ENABLE_AHF=$enableval],[ENABLE_AHF=no] )
+
+#if test "x$ENABLE_AHF" = "xyes"; then
+# AM_CPPFLAGS="$AM_CPPFLAGS -DUSE_AHF"
+#fi
+#AM_CONDITIONAL([ENABLE_AHF],[test x$ahf1 = xtrue])
+#AC_SUBST(ENABLE_AHF)
+
+
+#test "xyes" != "x$WITH_AHF" || AM_CPPFLAGS="$AM_CPPFLAGS -DUSE_AHF"
+#AM_CONDITIONAL(USE_AHF, [test "xno" != "x$WITH_AHF"])
+
+
+
diff --git a/src/Core.cpp b/src/Core.cpp
index 0bd882f..3641fd5 100644
--- a/src/Core.cpp
+++ b/src/Core.cpp
@@ -56,6 +56,10 @@
#include <errno.h>
#include <string.h>
+#ifdef USE_AHF
+#include "moab/HalfFacetRep.hpp"
+#endif
+
#ifdef USE_MPI
/* Leave ParallelComm.hpp before mpi.h or MPICH2 will fail
* because its C++ headers do not like SEEK_* macros.
@@ -272,6 +276,12 @@ ErrorCode Core::initialize()
geom_dimension_tag();
globalId_tag();
+#ifdef USE_AHF
+ ahfRep = new HalfFacetRep(this);
+ if (!ahfRep)
+ return MB_MEMORY_ALLOCATION_FAILED;
+#endif
+
return MB_SUCCESS;
}
@@ -322,6 +332,11 @@ void Core::deinitialize()
if (mpiFinalize)
MPI_Finalize();
#endif
+
+#ifdef USE_AHF
+ delete ahfRep;
+ ahfRep = 0;
+#endif
}
ErrorCode Core::query_interface_type( const std::type_info& type, void*& ptr )
@@ -1338,6 +1353,80 @@ ErrorCode get_adjacencies_union( Core* gMB,
return result;
}
+
+///////////////////////////
+///
+#ifdef USE_AHF
+template <typename ITER> static inline
+ErrorCode get_adjacencies_intersection_ahf(Core *mb,
+ ITER begin, ITER end,
+ const int to_dimension,
+ std::vector<EntityHandle>& adj_entities )
+{
+ const size_t SORT_THRESHOLD = 200;
+ std::vector<EntityHandle> temp_vec;
+ std::vector<EntityHandle>::iterator adj_it, w_it;
+ ErrorCode result = MB_SUCCESS;
+
+ if (begin == end) {
+ adj_entities.clear(); // intersection
+ return MB_SUCCESS;
+ }
+
+ // First iteration is a special case if input list is empty.
+ // Rather than returning nothing (intersecting with empty
+ // input list), we begin with the adjacencies for the first entity.
+ if (adj_entities.empty()) {
+ EntityType type = TYPE_FROM_HANDLE(*begin);
+
+ if(to_dimension == 0 && type != MBPOLYHEDRON)
+ result = mb->get_connectivity(&(*begin), 1, adj_entities);
+ else
+ result = mb->a_half_facet_rep()->get_adjacencies(*begin, to_dimension, adj_entities);
+ if (MB_SUCCESS != result)
+ return result;
+ ++begin;
+ }
+
+ for (ITER from_it = begin; from_it != end; from_it++)
+ {
+ // running results kept in adj_entities; clear temp_vec, which is working space
+ temp_vec.clear();
+
+ // get the next set of adjacencies
+ EntityType type = TYPE_FROM_HANDLE(*from_it);
+ if(to_dimension == 0 && type != MBPOLYHEDRON)
+ result = mb->get_connectivity(&(*from_it), 1, temp_vec);
+ else
+ result = mb->a_half_facet_rep()->get_adjacencies(*from_it, to_dimension, temp_vec);
+ if (MB_SUCCESS != result)
+ return result;
+
+ // otherwise intersect with the current set of results
+ w_it = adj_it = adj_entities.begin();
+ if (temp_vec.size()*adj_entities.size() < SORT_THRESHOLD) {
+ for (; adj_it != adj_entities.end(); ++adj_it)
+ if (std::find(temp_vec.begin(), temp_vec.end(), *adj_it) != temp_vec.end())
+ { *w_it = *adj_it; ++w_it; }
+ }
+ else {
+ std::sort( temp_vec.begin(), temp_vec.end() );
+ for (; adj_it != adj_entities.end(); ++adj_it)
+ if (std::binary_search(temp_vec.begin(), temp_vec.end(), *adj_it))
+ { *w_it = *adj_it; ++w_it; }
+ }
+ adj_entities.erase( w_it, adj_entities.end() );
+
+ // we're intersecting, so if there are no more results, we're done
+ if (adj_entities.empty())
+ break;
+ }
+
+ return MB_SUCCESS;
+}
+#endif
+///////////////////////////////////////////
+
template <typename ITER> static inline
ErrorCode get_adjacencies_intersection( Core* mb,
ITER begin, ITER end,
@@ -1445,37 +1534,116 @@ ErrorCode Core::get_adjacencies( const EntityHandle *from_entities,
const int to_dimension,
const bool create_if_missing,
std::vector<EntityHandle> &adj_entities,
- const int operation_type )
+ const int operation_type
+#ifdef USE_AHF
+ , const bool use_ahf)
+#else
+ )
+#endif
{
- if (operation_type == Interface::INTERSECT)
- return get_adjacencies_intersection( this, from_entities, from_entities+num_entities,
- to_dimension, create_if_missing, adj_entities );
- else if (operation_type != Interface::UNION)
- return MB_FAILURE;
- // do union
- ErrorCode result;
- std::vector<EntityHandle> tmp_storage;
- const EntityHandle* conn;
- int len;
- for (int i = 0; i < num_entities; ++i) {
- if(to_dimension == 0 && TYPE_FROM_HANDLE(from_entities[0]) != MBPOLYHEDRON) {
- result = get_connectivity(from_entities[i], conn, len, false, &tmp_storage);
- adj_entities.insert( adj_entities.end(), conn, conn+len );
- if (MB_SUCCESS != result)
- return result;
+#ifdef USE_AHF
+ bool can_handle = true;
+ int source_dim = this->dimension_from_handle(from_entities[0]);
+
+ if ((source_dim > to_dimension) && (to_dimension != 0))
+ {
+ std::cout<<"Currently Not Supported by MOAB_AHF: Down Adjacencies"<<std::endl;
+ std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
+ can_handle = false;
}
- else {
- result = aEntityFactory->get_adjacencies(from_entities[i], to_dimension,
- create_if_missing, adj_entities);
- if (MB_SUCCESS != result)
- return result;
+ else if (((source_dim == 0) && (to_dimension == 2))||((source_dim == 0) && (to_dimension == 3)))
+ {
+ std::cout<<"Currently Not Supported by MOAB_AHF: Vertex to face/cell"<<std::endl;
+ std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
+ can_handle = false;
+ }
+ else if (TYPE_FROM_HANDLE(from_entities[0]) == MBPOLYHEDRON)
+ {
+ std::cout<<"Currently Not Supported by MOAB_AHF: Polyhedron Meshes"<<std::endl;
+ std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
+ can_handle = false;
+ }
+ else if (create_if_missing)
+ {
+ std::cout<<"Currently Not Supporteded by MOAB_AHF: create_if_missing capability "<<std::endl;
+ std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
+ can_handle = false;
}
- }
- std::sort( adj_entities.begin(), adj_entities.end() );
- adj_entities.erase( std::unique( adj_entities.begin(), adj_entities.end() ), adj_entities.end() );
- return MB_SUCCESS;
+ if ((use_ahf) && (can_handle))
+ {
+
+ if (operation_type == Interface::INTERSECT)
+ return get_adjacencies_intersection_ahf(this, from_entities, from_entities+num_entities,
+ to_dimension, adj_entities );
+ else if (operation_type != Interface::UNION)
+ return MB_FAILURE;
+
+ // do union
+ ErrorCode result;
+ std::vector<EntityHandle> tmp_storage;
+ const EntityHandle* conn;
+ int len;
+ for (int i = 0; i < num_entities; ++i) {
+ if(to_dimension == 0 && TYPE_FROM_HANDLE(from_entities[0]) != MBPOLYHEDRON) {
+ result = get_connectivity(from_entities[i], conn, len, false, &tmp_storage);
+ adj_entities.insert( adj_entities.end(), conn, conn+len );
+ if (MB_SUCCESS != result)
+ return result;
+ }
+ else {
+ result = ahfRep->get_adjacencies(from_entities[i], to_dimension, adj_entities);
+ if (MB_SUCCESS != result)
+ return result;
+ }
+ }
+ std::sort( adj_entities.begin(), adj_entities.end() );
+ adj_entities.erase( std::unique( adj_entities.begin(), adj_entities.end() ), adj_entities.end() );
+
+ }
+ else
+ {
+
+#endif
+
+ if (operation_type == Interface::INTERSECT)
+ return get_adjacencies_intersection( this, from_entities, from_entities+num_entities,
+ to_dimension, create_if_missing, adj_entities );
+ else if (operation_type != Interface::UNION)
+ return MB_FAILURE;
+
+ // do union
+ ErrorCode result;
+ std::vector<EntityHandle> tmp_storage;
+ const EntityHandle* conn;
+ int len;
+ for (int i = 0; i < num_entities; ++i) {
+ if(to_dimension == 0 && TYPE_FROM_HANDLE(from_entities[0]) != MBPOLYHEDRON) {
+ result = get_connectivity(from_entities[i], conn, len, false, &tmp_storage);
+ adj_entities.insert( adj_entities.end(), conn, conn+len );
+ if (MB_SUCCESS != result)
+ return result;
+ }
+ else {
+ result = aEntityFactory->get_adjacencies(from_entities[i], to_dimension,
+ create_if_missing, adj_entities);
+ if (MB_SUCCESS != result)
+ return result;
+ }
+ }
+ std::sort( adj_entities.begin(), adj_entities.end() );
+ adj_entities.erase( std::unique( adj_entities.begin(), adj_entities.end() ), adj_entities.end() );
+
+ //return MB_SUCCESS;
+
+
+#ifdef USE_AHF
+ }
+#endif
+
+return MB_SUCCESS;
+
}
diff --git a/src/HalfFacetRep.cpp b/src/HalfFacetRep.cpp
index 95efb7e..f3818a6 100755
--- a/src/HalfFacetRep.cpp
+++ b/src/HalfFacetRep.cpp
@@ -17,21 +17,32 @@
#pragma warning (disable : 4786)
#endif
+#include "moab/HalfFacetRep.hpp"
#include <iostream>
#include <assert.h>
-#include <time.h>
#include <vector>
-#include <queue>
-#include <stack>
#include "moab/Core.hpp"
#include "moab/Range.hpp"
#include "moab/CN.hpp"
-#include "moab/HalfFacetRep.hpp"
namespace moab {
const int MAXSIZE = 500;
+ HalfFacetRep::HalfFacetRep(Core *impl)
+ {
+ assert(NULL != impl);
+ mb = impl;
+ mInitAHFmaps = false;
+ }
+
+ HalfFacetRep::~HalfFacetRep()
+ {
+ ErrorCode result;
+ result = deinitialize();
+ //if (MB_SUCCESS != result) return result;
+ }
+
MESHTYPE HalfFacetRep::get_mesh_type(int nverts, int nedges, int nfaces, int ncells)
{
MESHTYPE mesh_type;
@@ -58,7 +69,9 @@ namespace moab {
* initialize *
******************************************************/
- ErrorCode HalfFacetRep::initialize(){
+ ErrorCode HalfFacetRep::initialize()
+ {
+ mInitAHFmaps = true;
ErrorCode error;
@@ -436,36 +449,100 @@ namespace moab {
* User interface for adjacency functions *
********************************************************/
- ErrorCode HalfFacetRep::get_upward_incidences(EntityHandle ent, int out_dim, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * lids)
+ ErrorCode HalfFacetRep::get_adjacencies(const EntityHandle source_entity,
+ const unsigned int target_dimension,
+ std::vector<EntityHandle> &target_entities)
+ {
+
+ ErrorCode error;
+
+ unsigned int source_dimension = mb->dimension_from_handle(source_entity);
+
+ if (mInitAHFmaps == false)
+ {
+ error = initialize();
+ if (MB_SUCCESS != error) return error;
+ }
+
+ if ((source_dimension == 0) && (target_dimension == 1))
+ {
+ error = get_up_adjacencies_1d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+
+ else if ((source_dimension == 1) && (target_dimension == 2))
+ {
+ error = get_up_adjacencies_2d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+ else if ((source_dimension == 1) && (target_dimension == 3))
+ {
+ error = get_up_adjacencies_edg_3d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+ else if ((source_dimension == 2) && (target_dimension ==3))
+ {
+ error = get_up_adjacencies_face_3d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+ else if (source_dimension == target_dimension)
+ {
+ if (target_dimension == 1)
+ {
+ error = get_neighbor_adjacencies_1d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+
+ else if (target_dimension == 2)
+ {
+ error = get_neighbor_adjacencies_2d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+ else if (target_dimension == 3)
+ {
+ error = get_neighbor_adjacencies_3d(source_entity, target_entities);
+ if (MB_SUCCESS != error) return error;
+ }
+ }
+ return MB_SUCCESS;
+ }
+
+
+ ErrorCode HalfFacetRep::get_up_adjacencies(EntityHandle ent,
+ int out_dim,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * lids)
{
ErrorCode error;
int in_dim = mb->dimension_from_handle(ent);
if ((in_dim == 0) && (out_dim == 1))
{
- error = get_upward_incidences_1d(ent, adjents, local_id, lids);
+ error = get_up_adjacencies_1d(ent, adjents, local_id, lids);
if (MB_SUCCESS != error) return error;
}
else if ((in_dim == 1) && (out_dim == 2))
{
- error = get_upward_incidences_2d(ent, adjents, local_id, lids);
+ error = get_up_adjacencies_2d(ent, adjents, local_id, lids);
if (MB_SUCCESS != error) return error;
}
else if ((in_dim == 1) && (out_dim == 3))
{
- error = get_upward_incidences_edg_3d(ent, adjents, local_id, lids);
+ error = get_up_adjacencies_edg_3d(ent, adjents, local_id, lids);
if (MB_SUCCESS != error) return error;
}
else if ((in_dim == 2) && (out_dim ==3))
{
- error = get_upward_incidences_face_3d(ent, adjents, local_id, lids);
+ error = get_up_adjacencies_face_3d(ent, adjents, local_id, lids);
if (MB_SUCCESS != error) return error;
}
return MB_SUCCESS;
}
- ErrorCode HalfFacetRep::get_neighbor_adjacencies(EntityHandle ent, std::vector<EntityHandle> &adjents)
+ ErrorCode HalfFacetRep::get_neighbor_adjacencies(EntityHandle ent,
+ std::vector<EntityHandle> &adjents)
{
ErrorCode error;
int in_dim = mb->dimension_from_handle(ent);
@@ -619,8 +696,11 @@ namespace moab {
return MB_SUCCESS;
}
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_1d( EntityHandle vid, std::vector< EntityHandle > &adjents, bool local_id, std::vector<int> * lvids){
- ErrorCode error;
+ ErrorCode HalfFacetRep::get_up_adjacencies_1d( EntityHandle vid,
+ std::vector< EntityHandle > &adjents,
+ bool local_id,
+ std::vector<int> * lvids){
+ ErrorCode error;
EntityHandle start_eid, eid, sibeid[2];
int start_lid, lid, siblid[2];
@@ -655,7 +735,8 @@ namespace moab {
return MB_SUCCESS;
}
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_neighbor_adjacencies_1d( EntityHandle eid, std::vector<EntityHandle> &adjents){
+ ErrorCode HalfFacetRep::get_neighbor_adjacencies_1d( EntityHandle eid,
+ std::vector<EntityHandle> &adjents){
ErrorCode error;
@@ -912,9 +993,65 @@ namespace moab {
return MB_SUCCESS;
}
+ ///////////////////////////////////////////////////////////////////
+ /* ErrorCode HalfFacetRep::get_up_adjacencies_2d(EntityHandle vid, std::vector<EntityHandle> &adjents)
+ {
+ ErrorCode error;
+
+ EntityHandle fid; int lid;
+ error = mb->tag_get_data(v2he_fid, &vid, 1, &fid);
+ if (MB_SUCCESS != error) return error;
+ error = mb->tag_get_data(v2he_leid, &vid, 1, &lid);
+ if (MB_SUCCESS != error) return error;
+
+ if (fid != 0)
+ {
+ adjents.push_back(fid);
+
+ EntityHandle queue_fid[MAXSIZE], trackfaces[MAXSIZE];
+ int queue_lid[MAXSIZE];
+ for (int i = 0; i< MAXSIZE; i++)
+ {
+ queue_fid[i] = 0;
+ queue_lid[i] = 0;
+ trackfaces[i] = 0;
+ }
+ int qsize = 0, tcount = -1;
+ int num_qvals = 0;
+ error = gather_halfedges(vid, fid, lid, queue_fid, queue_lid, &qsize, trackfaces, &tcount);
+ if (MB_SUCCESS != error) return error;
+
+ while (num_qvals < *qsize)
+ {
+ EntityHandle curfid = queue_fid[num_qvals];
+ int curlid = queue_lid[num_qvals];
+ num_qvals += 1;
+
+ EntityHandle he2_fid; int he2_lid;
+ error = another_halfedge(vid, curfid, curlid, &he2_fid, &he2_lid);
+ if (MB_SUCCESS != error) return error;
+ bool found_ent = find_match_in_array(he2_fid, trackfaces, tcount[0]);
+
+ if (found_ent)
+ continue;
+ tcount[0] += 1;
+ trackfaces[tcount[0]] = he2_fid;
+
+ error = get_up_adjacencies_2d(he2_fid, he2_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
+ if (MB_SUCCESS != error) return error;
+
+ adjents.push_back(he2_fid);
+
+ }
+ }
+ }
+*/
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_2d( EntityHandle eid, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * leids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_2d( EntityHandle eid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * leids)
{
// Given an explicit edge eid, find the incident faces.
@@ -928,7 +1065,7 @@ namespace moab {
// Step 2: If there is a corresponding half-edge, collect all sibling half-edges and store the incident faces.
if (found)
{
- error = get_upward_incidences_2d(he_fid, he_lid, true, adjents, local_id, leids);
+ error = get_up_adjacencies_2d(he_fid, he_lid, true, adjents, local_id, leids);
if (MB_SUCCESS != error) return error;
}
@@ -936,7 +1073,12 @@ namespace moab {
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_2d( EntityHandle fid, int leid, bool add_inent, std::vector<EntityHandle> &fids, bool local_id, std::vector<int> * leids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_2d( EntityHandle fid,
+ int leid,
+ bool add_inent,
+ std::vector<EntityHandle> &fids,
+ bool local_id,
+ std::vector<int> * leids)
{
// Given an implicit half-edge <fid, leid>, find the incident half-edges.
ErrorCode error;
@@ -982,7 +1124,13 @@ namespace moab {
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_2d(EntityHandle fid, int lid, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount)
+ ErrorCode HalfFacetRep::get_up_adjacencies_2d(EntityHandle fid,
+ int lid,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount)
{
ErrorCode error;
@@ -1035,7 +1183,9 @@ namespace moab {
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- bool HalfFacetRep::find_matching_halfedge( EntityHandle eid, EntityHandle *hefid, int *helid){
+ bool HalfFacetRep::find_matching_halfedge( EntityHandle eid,
+ EntityHandle *hefid,
+ int *helid){
ErrorCode error;
std::vector<EntityHandle> conn(2);
error = mb->get_connectivity(&eid, 1, conn);
@@ -1071,7 +1221,14 @@ namespace moab {
return found;
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::gather_halfedges( EntityHandle vid, EntityHandle he_fid, int he_lid, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount)
+ ErrorCode HalfFacetRep::gather_halfedges( EntityHandle vid,
+ EntityHandle he_fid,
+ int he_lid,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount)
{
ErrorCode error;
EntityHandle he2_fid = 0; int he2_lid = 0;
@@ -1086,16 +1243,20 @@ namespace moab {
tcount[0] += 1;
trackfaces[tcount[0]] = he_fid;
- error = get_upward_incidences_2d(he_fid, he_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
+ error = get_up_adjacencies_2d(he_fid, he_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
if (MB_SUCCESS != error) return error;
- error = get_upward_incidences_2d(he2_fid, he2_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
+ error = get_up_adjacencies_2d(he2_fid, he2_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
if (MB_SUCCESS != error) return error;
return MB_SUCCESS;
}
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::another_halfedge( EntityHandle vid, EntityHandle he_fid, int he_lid, EntityHandle *he2_fid, int *he2_lid)
+ ErrorCode HalfFacetRep::another_halfedge( EntityHandle vid,
+ EntityHandle he_fid,
+ int he_lid,
+ EntityHandle *he2_fid,
+ int *he2_lid)
{
ErrorCode error;
int nepf;
@@ -1121,7 +1282,14 @@ namespace moab {
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- bool HalfFacetRep::collect_and_compare(std::vector<EntityHandle> &edg_vert, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount, EntityHandle *he_fid, int *he_lid)
+ bool HalfFacetRep::collect_and_compare(std::vector<EntityHandle> &edg_vert,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount,
+ EntityHandle *he_fid,
+ int *he_lid)
{
ErrorCode error;
int nepf = local_maps_2d(*_faces.begin());
@@ -1162,7 +1330,7 @@ namespace moab {
error = another_halfedge(edg_vert[0], curfid, curlid, &he2_fid, &he2_lid);
if (MB_SUCCESS != error) return error;
- error = get_upward_incidences_2d(he2_fid, he2_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
+ error = get_up_adjacencies_2d(he2_fid, he2_lid, queue_fid, queue_lid, qsize, trackfaces, tcount);
if (MB_SUCCESS != error) return error;
counter += 1;
@@ -1174,14 +1342,15 @@ namespace moab {
}
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_neighbor_adjacencies_2d( EntityHandle fid, std::vector<EntityHandle> &adjents)
+ ErrorCode HalfFacetRep::get_neighbor_adjacencies_2d( EntityHandle fid,
+ std::vector<EntityHandle> &adjents)
{
ErrorCode error;
if (fid != 0){
int nepf = local_maps_2d(fid);
for (int lid = 0; lid < nepf; ++lid){
- error = get_upward_incidences_2d(fid, lid, false, adjents);
+ error = get_up_adjacencies_2d(fid, lid, false, adjents);
if (MB_SUCCESS != error) return error;
}
}
@@ -1212,7 +1381,7 @@ namespace moab {
int id = nepf*(*f-firstF)+l;
if (!trackF[id])
{
- error = get_upward_incidences_2d(*f,l, false, adj_fids, true, &adj_lids);
+ error = get_up_adjacencies_2d(*f,l, false, adj_fids, true, &adj_lids);
if (MB_SUCCESS != error) return error;
total_edges -= adj_fids.size();
@@ -1530,7 +1699,10 @@ namespace moab {
}
//////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_edg_3d( EntityHandle eid, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * leids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_edg_3d( EntityHandle eid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * leids)
{
ErrorCode error;
@@ -1542,7 +1714,7 @@ namespace moab {
//Find all incident cells
if (found)
{
- error =get_upward_incidences_edg_3d(cid, leid, adjents, local_id, leids);
+ error =get_up_adjacencies_edg_3d(cid, leid, adjents, local_id, leids);
if (MB_SUCCESS != error) return error;
}
@@ -1550,7 +1722,11 @@ namespace moab {
}
//////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_edg_3d( EntityHandle cid, int leid, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * leids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_edg_3d( EntityHandle cid,
+ int leid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * leids)
{
ErrorCode error;
@@ -1648,7 +1824,10 @@ namespace moab {
}
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_face_3d( EntityHandle fid, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * lfids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_face_3d( EntityHandle fid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * lfids)
{
ErrorCode error;
@@ -1657,14 +1836,18 @@ namespace moab {
bool found = find_matching_halfface(fid, &cid, &lid);
if (found){
- error = get_upward_incidences_face_3d(cid, lid, adjents, local_id, lfids);
+ error = get_up_adjacencies_face_3d(cid, lid, adjents, local_id, lfids);
if (MB_SUCCESS != error) return error;
}
return MB_SUCCESS;
}
///////////////////////////////////////////
- ErrorCode HalfFacetRep::get_upward_incidences_face_3d( EntityHandle cid, int lfid, std::vector<EntityHandle> &adjents, bool local_id, std::vector<int> * lfids)
+ ErrorCode HalfFacetRep::get_up_adjacencies_face_3d( EntityHandle cid,
+ int lfid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id,
+ std::vector<int> * lfids)
{
ErrorCode error;
@@ -1701,7 +1884,9 @@ namespace moab {
return MB_SUCCESS;
}
/////////////////////////////////////////////////////////////////
- bool HalfFacetRep::find_matching_implicit_edge_in_cell( EntityHandle eid, EntityHandle *cid, int *leid)
+ bool HalfFacetRep::find_matching_implicit_edge_in_cell( EntityHandle eid,
+ EntityHandle *cid,
+ int *leid)
{
ErrorCode error;
diff --git a/src/Makefile.am b/src/Makefile.am
index e2bf713..6fac039 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -226,6 +226,12 @@ if HAVE_CGM
AM_CPPFLAGS += @CGM_CPPFLAGS@ -DCGM @MOAB_CGM_DEFINES@
endif
+# if ENABLE_ahf
+ # libMOAB_la_SOURCES += HalfFacetRep.cpp
+ # nobase_libMOAB_la_include_HEADERS += moab/HalfFacetRep.hpp
+# AM_CPPFLAGS += -DUSE_AHF
+# endif
+
MBCN_protos.h: MBCN.h $(top_srcdir)/itaps/mkprotos.sh
$(AM_V_GEN)$(top_srcdir)/itaps/mkprotos.sh MBCN MOAB $< $@ MOAB_FCDefs.h
diff --git a/src/moab/Core.hpp b/src/moab/Core.hpp
index aa7d9da..941129a 100644
--- a/src/moab/Core.hpp
+++ b/src/moab/Core.hpp
@@ -35,6 +35,10 @@ class EntitySequence;
class FileOptions;
class SetIterator;
+#ifdef USE_AHF
+class HalfFacetRep;
+#endif
+
#ifdef XPCOM_MB
#define MBCORE_CID \
@@ -311,12 +315,18 @@ public:
get_adjacencies( from_entities, MB_1D_ENTITY, adjacencies );
\endcode */
- virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
- const int num_entities,
- const int to_dimension,
- const bool create_if_missing,
- std::vector<EntityHandle>& adj_entities,
- const int operation_type = Interface::INTERSECT);
+ virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
+ const int num_entities,
+ const int to_dimension,
+ const bool create_if_missing,
+ std::vector<EntityHandle>& adj_entities,
+ const int operation_type = Interface::INTERSECT
+ #ifdef USE_AHF
+ , const bool use_ahf = false);
+ #else
+ );
+ #endif
+
virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
const int num_entities,
@@ -1126,7 +1136,12 @@ public:
//! return the a_entity_factory pointer
AEntityFactory *a_entity_factory() { return aEntityFactory; }
const AEntityFactory *a_entity_factory() const { return aEntityFactory; }
-
+
+#ifdef USE_AHF
+ HalfFacetRep *a_half_facet_rep() { return ahfRep; }
+ const HalfFacetRep *a_half_facet_rep() const {return ahfRep; }
+#endif
+
//! return set of registered IO tools
ReaderWriterSet* reader_writer_set() { return readerWriterSet; }
@@ -1347,6 +1362,10 @@ private:
//! list of iterators
std::vector<SetIterator*> setIterators;
+
+#ifdef USE_AHF
+ HalfFacetRep *ahfRep;
+#endif
};
diff --git a/src/moab/HalfFacetRep.hpp b/src/moab/HalfFacetRep.hpp
index eca5e91..3d1069a 100755
--- a/src/moab/HalfFacetRep.hpp
+++ b/src/moab/HalfFacetRep.hpp
@@ -52,14 +52,15 @@ namespace moab {
VOLUME_MIXED //Volume mesh with embedded curves and surfaces
};
-
+ class Core;
+
class HalfFacetRep{
public:
- HalfFacetRep(Interface *impl) : mb(impl) {}
+ HalfFacetRep(Core *impl);
- ~HalfFacetRep() {}
+ ~HalfFacetRep();
// User interface functions
@@ -74,6 +75,12 @@ namespace moab {
//! Prints the tag values.
ErrorCode print_tags();
+
+ ErrorCode get_adjacencies(const EntityHandle source_entity,
+ const unsigned int target_dimension,
+ std::vector<EntityHandle> &target_entities);
+
+
//! Get the upward incidences associated with an entity.
/** Given an entity of dimension <em>d</em>, gather all the incident <em>D(>d)</em> dimensional entities.
* Parameters:
@@ -84,7 +91,11 @@ namespace moab {
* \param lids Vector in which the local id's are returned.
*/
- ErrorCode get_upward_incidences(EntityHandle ent, int out_dim, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * lids = NULL );
+ ErrorCode get_up_adjacencies(EntityHandle ent,
+ int out_dim,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * lids = NULL );
//! Get the same-dimensional entities connected with an entity.
/** Given an entity of dimension <em>d</em>, gather all the entities connected via <em>d-1</em> dimensional entities.
@@ -94,7 +105,8 @@ namespace moab {
* \param adjents Vector in which the neighbor entities are returned.
*/
- ErrorCode get_neighbor_adjacencies(EntityHandle ent, std::vector<EntityHandle> &adjents);
+ ErrorCode get_neighbor_adjacencies(EntityHandle ent,
+ std::vector<EntityHandle> &adjents);
// 1D Maps and queries
@@ -127,7 +139,10 @@ namespace moab {
* \param lvids Vector returning the local vertex id's
*/
- ErrorCode get_upward_incidences_1d(EntityHandle vid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * lvids = NULL);
+ ErrorCode get_up_adjacencies_1d(EntityHandle vid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * lvids = NULL);
//! Given an edge, finds vertex-connected neighbor edges
/** Given an edge, it gathers all the incident edges of each vertex of the edge.
@@ -136,7 +151,8 @@ namespace moab {
* \param adjents Vector returning neighbor edges
*/
- ErrorCode get_neighbor_adjacencies_1d(EntityHandle eid, std::vector<EntityHandle> &adjents);
+ ErrorCode get_neighbor_adjacencies_1d(EntityHandle eid,
+ std::vector<EntityHandle> &adjents);
// 2D Maps and queries
@@ -161,6 +177,9 @@ namespace moab {
ErrorCode determine_incident_halfedges(Range &faces);
+ /* ErrorCode get_up_adjacencies_2d(EntityHandle vid,
+ std::vector<EntityHandle> &adjents);*/
+
//! Given an edge, finds the faces incident on it.
/** Given an edge, it first finds a matching half-edge corresponding to eid, and then
* collects all the incident half-edges/faces via the sibhes map.
@@ -171,7 +190,10 @@ namespace moab {
* \param leids Vector returning local edge ids
*/
- ErrorCode get_upward_incidences_2d(EntityHandle eid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * leids = NULL);
+ ErrorCode get_up_adjacencies_2d(EntityHandle eid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * leids = NULL);
//! Given a half-edge <fid, leid>, finds the faces incident on it.
/**
@@ -184,7 +206,12 @@ namespace moab {
* \param leids Vector returning local edge ids
*/
- ErrorCode get_upward_incidences_2d(EntityHandle fid, int leid, bool add_inent, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * leids = NULL);
+ ErrorCode get_up_adjacencies_2d(EntityHandle fid,
+ int leid,
+ bool add_inent,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * leids = NULL);
//! Given an edge, finds edge-connected neighbor face
/** Given an face, it gathers all the neighbor faces of each local edge of the face.
@@ -193,7 +220,8 @@ namespace moab {
* \param adjents Vector returning neighbor faces
*/
- ErrorCode get_neighbor_adjacencies_2d(EntityHandle fid, std::vector<EntityHandle> &adjents);
+ ErrorCode get_neighbor_adjacencies_2d(EntityHandle fid,
+ std::vector<EntityHandle> &adjents);
//! Given a range of faces, finds the total number of edges.
@@ -229,7 +257,8 @@ namespace moab {
* \param isborder: A dense tag over all vertices of size 1. Value is true for a border vertex, otherwise is false.
*/
- ErrorCode determine_border_vertices( Range &cells, Tag isborder);
+ ErrorCode determine_border_vertices( Range &cells,
+ Tag isborder);
//! Given an edge, finds the cells incident on it.
/** Given an edge, it first finds a matching local edge in a cell corresponding to eid, and then
@@ -241,7 +270,10 @@ namespace moab {
* \param leids Vector returning local edge ids
*/
- ErrorCode get_upward_incidences_edg_3d(EntityHandle eid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * leids = NULL);
+ ErrorCode get_up_adjacencies_edg_3d(EntityHandle eid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * leids = NULL);
//! Given a local edge <cid, leid>, finds the cells incident on it.
/** Given a local edge, it gathers all the incident cells via the sibhfs map.
@@ -253,7 +285,10 @@ namespace moab {
* \param leids Vector returning local edge ids
*/
- ErrorCode get_upward_incidences_edg_3d(EntityHandle cid, int leid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * leids = NULL);
+ ErrorCode get_up_adjacencies_edg_3d(EntityHandle cid,
+ int leid, std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * leids = NULL);
//! Given an face, finds the cells incident on it.
/** Given an face, it first finds a matching half-face in a cell corresponding to face, and then
@@ -265,7 +300,10 @@ namespace moab {
* \param leids Vector returning local face ids
*/
- ErrorCode get_upward_incidences_face_3d(EntityHandle fid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * lfids = NULL);
+ ErrorCode get_up_adjacencies_face_3d(EntityHandle fid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * lfids = NULL);
//! Given a local face <cid, lfid>, finds the cells incident on it.
/** Given a local face, it gathers all the incident cells via the sibhfs map.
@@ -277,7 +315,11 @@ namespace moab {
* \param lfids Vector returning local face ids
*/
- ErrorCode get_upward_incidences_face_3d(EntityHandle cid, int lfid, std::vector<EntityHandle> &adjents, bool local_id = false, std::vector<int> * lfids = NULL);
+ ErrorCode get_up_adjacencies_face_3d(EntityHandle cid,
+ int lfid,
+ std::vector<EntityHandle> &adjents,
+ bool local_id = false,
+ std::vector<int> * lfids = NULL);
//! Given a cell, finds face-connected neighbor cells
/** Given a cell, it gathers all the neighbor cells of each local face of the cell.
@@ -286,12 +328,17 @@ namespace moab {
* \param adjents Vector returning neighbor cells
*/
- ErrorCode get_neighbor_adjacencies_3d(EntityHandle cid, std::vector<EntityHandle> &adjents);
+ ErrorCode get_neighbor_adjacencies_3d(EntityHandle cid,
+ std::vector<EntityHandle> &adjents);
protected:
- Interface * mb;
+ Core * mb;
+
+ HalfFacetRep();
+
+ bool mInitAHFmaps;
enum {
MAX_VERTICES = 8,
@@ -351,14 +398,23 @@ namespace moab {
* queue of half-edges, if they do not already exist in the queue. This function is used to increment the
* search space for finding a matching half-edge.
* Parameters:
- * \param <EntityHandle he_fid, int he_lid> The query half-edge
- * \param queue_fid, queue_lid, qsize Array of faces and local edge ids. qsize is the current size of the queue_fid.
+ * \param he_fid EntityHandle of query half-edge
+ * \param he_lid Local id of query half-edge
+ * \param queue_fid
+ * \param queue_lid
+ * \param qsize Array of faces and local edge ids. qsize is the current size of the queue_fid.
* \param trackfaces Array containing faces. If fid of an incident half-edge doesn't belong to trackfaces,
* the half-edge is added to the queue.
* \param tcount Current size of trackfaces
*/
- ErrorCode get_upward_incidences_2d(EntityHandle he_fid, int he_lid, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount);
+ ErrorCode get_up_adjacencies_2d(EntityHandle he_fid,
+ int he_lid,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount);
//! Given an edge, finds a matching half-edge in the surface.
/** Given an edge eid, it first collects few half-edges belonging to one-ring neighborhood of
@@ -369,7 +425,9 @@ namespace moab {
* \param hefid, helid: Returns the matching half-edge corresponding to the query edge.
*/
- bool find_matching_halfedge( EntityHandle eid, EntityHandle *hefid, int *helid);
+ bool find_matching_halfedge( EntityHandle eid,
+ EntityHandle *hefid,
+ int *helid);
//! Gather half-edges to a queue of half-edges.
/** Given a vertex vid, and a half-edge <he_fid,he_lid>, add another half-edge in the same face sharing the vertex
@@ -381,7 +439,14 @@ namespace moab {
* \param trackfaces, tcount
*/
- ErrorCode gather_halfedges( EntityHandle vid, EntityHandle he_fid, int he_lid, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount);
+ ErrorCode gather_halfedges( EntityHandle vid,
+ EntityHandle he_fid,
+ int he_lid,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount);
//! Obtains another half-edge belonging to the same face as the input half-edge
/** It uses the local maps to find another half-edge that is either incident or outgoing depending
@@ -392,7 +457,11 @@ namespace moab {
* \param <EntityHandle he2_fid, int he2_lid>: Returns another half-edge in the same he_fid sharing vid.
*/
- ErrorCode another_halfedge( EntityHandle vid, EntityHandle he_fid, int he_lid, EntityHandle *he2_fid, int *he2_lid);
+ ErrorCode another_halfedge( EntityHandle vid,
+ EntityHandle he_fid,
+ int he_lid,
+ EntityHandle *he2_fid,
+ int *he2_lid);
//! Collect and compare to find a matching half-edge with the given edge connectivity.
/** Given edge connectivity, compare to an input list of half-edges to find a matching half-edge
@@ -404,7 +473,14 @@ namespace moab {
* \param <EntityHandle he_fid, int he_lid>: Returns matching half-edge
*/
- bool collect_and_compare(std::vector<EntityHandle> &edg_vert, EntityHandle *queue_fid, int *queue_lid, int *qsize, EntityHandle *trackfaces, int *tcount, EntityHandle *he_fid, int *he_lid);
+ bool collect_and_compare(std::vector<EntityHandle> &edg_vert,
+ EntityHandle *queue_fid,
+ int *queue_lid,
+ int *qsize,
+ EntityHandle *trackfaces,
+ int *tcount,
+ EntityHandle *he_fid,
+ int *he_lid);
//! The local maps for 3D entities.
@@ -469,7 +545,9 @@ namespace moab {
* \param leid Returns the local id of the edge corresponding to the input edge w.r.t the incident cell.
*/
- bool find_matching_implicit_edge_in_cell( EntityHandle eid, EntityHandle *cid, int *leid);
+ bool find_matching_implicit_edge_in_cell( EntityHandle eid,
+ EntityHandle *cid,
+ int *leid);
//! Given a face, finds a matching local face in an incident cell.
/** Find a local face with the same connectivity as the input face, belonging to an incident cell.
@@ -479,9 +557,15 @@ namespace moab {
* \param lfid Returns the local id of the face corresponding to the input face w.r.t the incident cell.
*/
- bool find_matching_halfface(EntityHandle fid, EntityHandle *cid, int *leid);
+ bool find_matching_halfface(EntityHandle fid,
+ EntityHandle *cid,
+ int *leid);
- bool find_match_in_array(EntityHandle ent, EntityHandle *ent_list, int count, bool get_index = false, int *index = NULL);
+ bool find_match_in_array(EntityHandle ent,
+ EntityHandle *ent_list,
+ int count,
+ bool get_index = false,
+ int *index = NULL);
};
diff --git a/src/moab/Interface.hpp b/src/moab/Interface.hpp
index 09f4219..bfd5fb0 100644
--- a/src/moab/Interface.hpp
+++ b/src/moab/Interface.hpp
@@ -624,7 +624,12 @@ public:
const int to_dimension,
const bool create_if_missing,
std::vector<EntityHandle>& adj_entities,
- const int operation_type = Interface::INTERSECT) = 0;
+ const int operation_type = Interface::INTERSECT
+ #ifdef USE_AHF
+ , const bool use_ahf = false) = 0;
+ #else
+ ) = 0;
+ #endif
//! Get the adjacencies associated with a vector of entities to entities of a specfied dimension.
/** Identical to vector-based get_adjacencies function, except results are returned in a
diff --git a/test/Makefile.am b/test/Makefile.am
index fe7edb6..34432d0 100644
--- a/test/Makefile.am
+++ b/test/Makefile.am
@@ -61,6 +61,10 @@ if HDF5_FILE
lloyd_smoother_test
endif
+if ENABLE_ahf
+ TESTS += ahf_mbintf_test
+endif
+
# merge_test \ # input files no longer exist?
# test_tag_server \ # fails
@@ -154,6 +158,11 @@ lloyd_smoother_test_SOURCES = lloyd_smoother_test.cpp
ahf_test_SOURCES = test_ahf_moab.cpp
+if ENABLE_ahf
+ahf_mbintf_test_SOURCES = test_ahf_mb_interface.cpp
+endif
+
+
if PARALLEL
moab_test_CPPFLAGS += -I$(top_srcdir)/src/parallel
kd_tree_test_CPPFLAGS += -I$(top_srcdir)/src/parallel
diff --git a/test/test_ahf_mb_interface.cpp b/test/test_ahf_mb_interface.cpp
new file mode 100755
index 0000000..9bf76ef
--- /dev/null
+++ b/test/test_ahf_mb_interface.cpp
@@ -0,0 +1,174 @@
+/*This function tests the AHF datastructures on CST meshes*/
+#include <iostream>
+#include <vector>
+#include <algorithm>
+#include "moab/Core.hpp"
+#include "moab/Range.hpp"
+#include "moab/MeshTopoUtil.hpp"
+#include "moab/HalfFacetRep.hpp"
+#include "TestUtil.hpp"
+
+using namespace moab;
+
+#ifdef MESHDIR
+static const char example[] = STRINGIFY(MESHDIR) "/hexes_mixed.vtk";
+#else
+static const char example[] = "/hexes_mixed.vtk";
+#endif
+
+void ahf_mbintf_test()
+{
+
+ Core moab;
+ Interface* mbImpl = &moab;
+ MeshTopoUtil mtu(mbImpl);
+
+ ErrorCode error = mbImpl->load_file(example);
+ CHECK_ERR(error);
+
+ /*Create ranges for handles of explicit elements of the mixed mesh*/
+ Range verts, edges, faces, cells;
+ error = mbImpl->get_entities_by_dimension( 0, 0, verts);
+ error = mbImpl->get_entities_by_dimension( 0, 1, edges);
+ error = mbImpl->get_entities_by_dimension( 0, 2, faces);
+ error = mbImpl->get_entities_by_dimension( 0, 3, cells);
+
+ //Perform queries
+ std::vector<EntityHandle> adjents;
+ Range mbents, ahfents;
+
+ //1D Queries //
+ //IQ1: For every vertex, obtain incident edges
+ for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 1, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 1, false, mbents );
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(),mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ1: For every edge, obtain neighbor edges
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 1, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 0, 1, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ // 2D Queries
+ //IQ2: For every edge, obtain incident faces
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 2, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 2, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ2: For every face, obtain neighbor faces
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 2, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 1, 2, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ // 3D Queries
+ // IQ 31: For every edge, obtain incident cells
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 3, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //IQ32: For every face, obtain incident cells
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 3, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ3: For every cell, obtain neighbor cells
+ for (Range::iterator i = cells.begin(); i != cells.end(); ++i) {
+ adjents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 3, false, adjents, Interface::INTERSECT, true);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 2, 3, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+}
+
+int main(int argc, char *argv[])
+{
+ int result = 0;
+
+ argv[0] = argv[argc - argc]; // Followed read_mpas_nc.cpp test for removing warnings in serial mode about unused variables.
+
+ result += RUN_TEST(ahf_mbintf_test);
+
+ return result;
+}
+
diff --git a/test/test_ahf_moab.cpp b/test/test_ahf_moab.cpp
index ff97fb1..93bd947 100755
--- a/test/test_ahf_moab.cpp
+++ b/test/test_ahf_moab.cpp
@@ -34,7 +34,7 @@ void ahf_test()
error = mbImpl->get_entities_by_dimension( 0, 3, cells);
// Create an ahf instance
- HalfFacetRep ahf(mbImpl);
+ HalfFacetRep ahf(&moab);
// Call the initialize function which creates the maps for each dimension
ahf.initialize();
@@ -47,7 +47,7 @@ void ahf_test()
//IQ1: For every vertex, obtain incident edges
for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
adjents.clear();
- error = ahf.get_upward_incidences( *i, 1, adjents);
+ error = ahf.get_up_adjacencies( *i, 1, adjents);
CHECK_ERR(error);
mbents.clear();
error = mbImpl->get_adjacencies( &*i, 1, 1, false, mbents );
@@ -82,7 +82,7 @@ void ahf_test()
//IQ2: For every edge, obtain incident faces
for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
adjents.clear();
- error = ahf.get_upward_incidences( *i, 2, adjents);
+ error = ahf.get_up_adjacencies( *i, 2, adjents);
CHECK_ERR(error);
mbents.clear();
error = mbImpl->get_adjacencies( &*i, 1, 2, false, mbents);
@@ -117,7 +117,7 @@ void ahf_test()
// IQ 31: For every edge, obtain incident cells
for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
adjents.clear();
- error = ahf.get_upward_incidences( *i, 3, adjents);
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
CHECK_ERR(error);
mbents.clear();
error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
@@ -134,7 +134,7 @@ void ahf_test()
//IQ32: For every face, obtain incident cells
for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
adjents.clear();
- error = ahf.get_upward_incidences( *i, 3, adjents);
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
CHECK_ERR(error);
mbents.clear();
error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
@@ -165,7 +165,7 @@ void ahf_test()
CHECK(!mbents.size());
}
- ahf.deinitialize();
+ //ahf.deinitialize();
}
https://bitbucket.org/fathomteam/moab/commits/58f615d2bcf2/
Changeset: 58f615d2bcf2
Branch: nray/feature_ahf
User: nray
Date: 2014-07-01 05:28:15
Summary: Removed tag_iterate from sibling half-facet maps. Changed input to
test_ahf_moab to read in meshfile from command line. Added a test code
for timing and memory requirements.
Affected #: 9 files
diff --git a/src/Core.cpp b/src/Core.cpp
index 3641fd5..692495b 100644
--- a/src/Core.cpp
+++ b/src/Core.cpp
@@ -1540,6 +1540,12 @@ ErrorCode Core::get_adjacencies( const EntityHandle *from_entities,
#else
)
#endif
+/*ErrorCode Core::get_adjacencies( const EntityHandle *from_entities,
+ const int num_entities,
+ const int to_dimension,
+ const bool create_if_missing,
+ std::vector<EntityHandle> &adj_entities,
+ const int operation_type)*/
{
#ifdef USE_AHF
@@ -1558,6 +1564,12 @@ ErrorCode Core::get_adjacencies( const EntityHandle *from_entities,
std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
can_handle = false;
}
+ else if (to_dimension == 4)
+ {
+ std::cout<<"Currently Not Supported by MOAB_AHF: meshsets"<<std::endl;
+ std::cout<<"Reverting to MOAB adjacency functionality"<<std::endl;
+ can_handle = false;
+ }
else if (TYPE_FROM_HANDLE(from_entities[0]) == MBPOLYHEDRON)
{
std::cout<<"Currently Not Supported by MOAB_AHF: Polyhedron Meshes"<<std::endl;
@@ -1571,7 +1583,9 @@ ErrorCode Core::get_adjacencies( const EntityHandle *from_entities,
can_handle = false;
}
+
if ((use_ahf) && (can_handle))
+ //if (can_handle)
{
if (operation_type == Interface::INTERSECT)
diff --git a/src/HalfFacetRep.cpp b/src/HalfFacetRep.cpp
index f3818a6..449a48a 100755
--- a/src/HalfFacetRep.cpp
+++ b/src/HalfFacetRep.cpp
@@ -27,7 +27,7 @@
namespace moab {
- const int MAXSIZE = 500;
+ const int MAXSIZE = 150;
HalfFacetRep::HalfFacetRep(Core *impl)
{
@@ -573,15 +573,6 @@ namespace moab {
{
ErrorCode error;
- EntityHandle start_edge = *edges.begin();
- int count, *sibhvs_lvid_ptr; EntityHandle *sibhvs_eid_ptr;
- error = mb->tag_iterate(sibhvs_eid, edges.begin(), edges.end(), count, (void*&)sibhvs_eid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) edges.size());
- error = mb->tag_iterate(sibhvs_lvid, edges.begin(), edges.end(), count, (void*&)sibhvs_lvid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) edges.size());
-
//Step 1: Create an index list storing the starting position for each vertex
int nv = _verts.size();
int *is_index = new int[nv+1];
@@ -595,8 +586,10 @@ namespace moab {
error = mb->get_connectivity(&*eid, 1, conn);
if (MB_SUCCESS != error) return error;
- is_index[conn[0]-*_verts.begin()+1] += 1;
- is_index[conn[1]-*_verts.begin()+1] += 1;
+ int index = _verts.index(conn[0]);
+ is_index[index+1] += 1;
+ index = _verts.index(conn[1]);
+ is_index[index+1] += 1;
}
is_index[0] = 0;
@@ -615,7 +608,7 @@ namespace moab {
for (int j = 0; j< 2; j++)
{
- int v = conn[j] - *_verts.begin();
+ int v = _verts.index(conn[j]);
v2hv_map_eid[is_index[v]] = *eid;
v2hv_map_lvid[is_index[v]] = j;
is_index[v] += 1;
@@ -633,26 +626,59 @@ namespace moab {
error = mb->get_connectivity(&*eid, 1, conn);
if (MB_SUCCESS != error) return error;
+ EntityHandle sibeid[2];
+ int siblvid[2];
+
+ error = mb->tag_get_data(sibhvs_eid, &*eid, 1, sibeid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhvs_lvid, &*eid, 1, siblvid);
+ if (MB_SUCCESS != error) return error;
+
for (int k =0; k<2; k++)
{
- if (sibhvs_eid_ptr[2*(*eid-start_edge)+k] != 0)
+
+ if (sibeid[k] != 0)
continue;
- int v = conn[k] - *_verts.begin();
+ int v = _verts.index(conn[k]);
int last = is_index[v+1] - 1;
if (last > is_index[v])
{
EntityHandle prev_eid = v2hv_map_eid[last];
int prev_lvid = v2hv_map_lvid[last];
+
+ EntityHandle seteid[2];
+ int setlvid[2];
+ error = mb->tag_get_data(sibhvs_eid, &prev_eid, 1, seteid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhvs_lvid, &prev_eid, 1, setlvid);
+ if (MB_SUCCESS != error) return error;
+
for (int i=is_index[v]; i<=last; i++)
{
EntityHandle cur_eid = v2hv_map_eid[i];
int cur_lvid = v2hv_map_lvid[i];
- sibhvs_eid_ptr[2*(prev_eid-start_edge) + prev_lvid] = cur_eid;
- sibhvs_lvid_ptr[2*(prev_eid-start_edge) + prev_lvid] = cur_lvid;
+
+
+ seteid[prev_lvid] = cur_eid;
+ setlvid[prev_lvid] = cur_lvid;
+
+ error = mb->tag_set_data(sibhvs_eid, &prev_eid, 1, seteid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_set_data(sibhvs_lvid, &prev_eid, 1, setlvid);
+ if (MB_SUCCESS != error) return error;
prev_eid = cur_eid;
prev_lvid = cur_lvid;
+ error = mb->tag_get_data(sibhvs_eid, &prev_eid, 1, seteid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhvs_lvid, &prev_eid, 1, setlvid);
+ if (MB_SUCCESS != error) return error;
+
}
}
}
@@ -663,6 +689,7 @@ namespace moab {
delete [] is_index;
delete [] v2hv_map_eid;
delete [] v2hv_map_lvid;
+
return MB_SUCCESS;
}
@@ -797,7 +824,7 @@ namespace moab {
int nepf;
if (type == MBTRI) nepf = 3;
else if (type ==MBQUAD) nepf = 4;
-
+
return nepf;
}
/////////////////////////////////////////////////////////////////////////////////
@@ -807,16 +834,13 @@ namespace moab {
// next: Local ids of next edges
// prev: Local ids of prev edges
- if (nepf == 3)
- {
- next[0] = 1; next[1] = 2; next[2] = 0;
- prev[0] = 2; prev[1] = 0; prev[2] = 1;
- }
- else if (nepf == 4)
- {
- next[0] = 1; next[1] = 2; next[2] = 3; next[3] = 0;
- prev[0] = 3; prev[1] = 0; prev[2] = 1; prev[3] = 2;
- }
+ for (int k=0; k<nepf-1; k++)
+ {
+ next[k]=k+1;
+ prev[k+1]=k;
+ }
+ next[nepf-1]=0;
+ prev[0] = nepf-1;
return MB_SUCCESS;
}
@@ -826,15 +850,6 @@ namespace moab {
ErrorCode error;
EntityHandle start_face = *faces.begin();
- // Create pointers to the tags for direct access to memory
- int count, *sibhes_leid_ptr; EntityHandle *sibhes_fid_ptr;
- error = mb->tag_iterate(sibhes_fid, faces.begin(), faces.end(), count, (void*&)sibhes_fid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) faces.size());
- error = mb->tag_iterate(sibhes_leid, faces.begin(), faces.end(), count, (void*&)sibhes_leid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) faces.size());
-
int nepf = local_maps_2d(start_face);
int * next = new int[nepf];
int * prev = new int[nepf];
@@ -847,6 +862,7 @@ namespace moab {
for (int i =0; i<nv+1; i++)
is_index[i] = 0;
+ int index;
std::vector<EntityHandle> conn(nepf);
for (Range::iterator fid = faces.begin(); fid != faces.end(); ++fid)
{
@@ -856,7 +872,8 @@ namespace moab {
for (int i = 0; i<nepf; i++)
{
- is_index[conn[i]-*_verts.begin()+1] += 1;
+ index = _verts.index(conn[i]);
+ is_index[index+1] += 1;
}
}
is_index[0] = 0;
@@ -877,7 +894,7 @@ namespace moab {
for (int j = 0; j< nepf; j++)
{
- int v = conn[j] - *_verts.begin();
+ int v = _verts.index(conn[j]);
v2nv[is_index[v]] = conn[next[j]];
v2he_map_fid[is_index[v]] = *fid;
v2he_map_leid[is_index[v]] = j;
@@ -889,6 +906,7 @@ namespace moab {
is_index[i+1] = is_index[i];
is_index[0] = 0;
+
//Step 3: Fill up sibling half-verts map
for (Range::iterator fid = faces.begin(); fid != faces.end(); ++fid)
{
@@ -896,13 +914,22 @@ namespace moab {
error = mb->get_connectivity(&*fid, 1, conn);
if (MB_SUCCESS != error) return error;
+ EntityHandle *sibfid = new EntityHandle[nepf];
+ int *sibleid = new int[nepf];
+
+ error = mb->tag_get_data(sibhes_fid, &*fid, 1, sibfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhes_leid, &*fid, 1, sibleid);
+ if (MB_SUCCESS != error) return error;
+
for (int k =0; k<nepf; k++)
{
- if (sibhes_fid_ptr[nepf*(*fid-start_face)+k] != 0)
+ if (sibfid[k] != 0)
continue;
- int v = conn[k] - *_verts.begin();
- int vn = conn[next[k]] - *_verts.begin();
+ int v = _verts.index(conn[k]);
+ int vn = _verts.index(conn[next[k]]);
EntityHandle first_fid = *fid;
int first_leid = k;
@@ -910,38 +937,87 @@ namespace moab {
EntityHandle prev_fid = *fid;
int prev_leid = k;
- for (int index = is_index[vn]; index <= is_index[vn+1]-1; index++)
+ EntityHandle *setfid = new EntityHandle[nepf];
+ int *setleid = new int[nepf];
+
+ for (index = is_index[vn]; index <= is_index[vn+1]-1; index++)
{
if (v2nv[index] == conn[k])
{
EntityHandle cur_fid = v2he_map_fid[index];
int cur_leid = v2he_map_leid[index];
- sibhes_fid_ptr[nepf*(prev_fid-start_face)+prev_leid] = cur_fid;
- sibhes_leid_ptr[nepf*(prev_fid-start_face)+prev_leid] = cur_leid;
+
+ error = mb->tag_get_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+
+ setfid[prev_leid] = cur_fid;
+ setleid[prev_leid] = cur_leid;
+
+ error = mb->tag_set_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_set_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+
prev_fid = cur_fid;
prev_leid = cur_leid;
+
}
}
- for (int index = is_index[v]; index <= is_index[v+1]-1; index++)
+ for (index = is_index[v]; index <= is_index[v+1]-1; index++)
{
if ((v2nv[index] == conn[next[k]])&&(v2he_map_fid[index] != *fid))
{
+
EntityHandle cur_fid = v2he_map_fid[index];
int cur_leid = v2he_map_leid[index];
- sibhes_fid_ptr[nepf*(prev_fid-start_face)+prev_leid] = cur_fid;
- sibhes_leid_ptr[nepf*(prev_fid-start_face)+prev_leid] = cur_leid;
+
+ error = mb->tag_get_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+
+ setfid[prev_leid] = cur_fid;
+ setleid[prev_leid] = cur_leid;
+
+ error = mb->tag_set_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_set_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+
prev_fid = cur_fid;
prev_leid = cur_leid;
}
}
if (prev_fid != first_fid){
- sibhes_fid_ptr[nepf*(prev_fid-start_face) + prev_leid] = first_fid;
- sibhes_leid_ptr[nepf*(prev_fid-start_face) + prev_leid] = first_leid;
- }
+ error = mb->tag_get_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+
+ setfid[prev_leid] = first_fid;
+ setleid[prev_leid] = first_leid;
+
+ error = mb->tag_set_data(sibhes_fid, &prev_fid, 1, setfid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_set_data(sibhes_leid, &prev_fid, 1, setleid);
+ if (MB_SUCCESS != error) return error;
+ }
+ delete [] setfid;
+ delete [] setleid;
}
+ delete [] sibfid;
+ delete [] sibleid;
}
delete [] next;
@@ -1311,7 +1387,7 @@ namespace moab {
error = mb->get_connectivity(&curfid, 1, conn);
if (MB_SUCCESS != error) return error;
- int id = next[curlid];
+ int id = next[curlid];
if (((conn[curlid]==edg_vert[0])&&(conn[id]==edg_vert[1]))||((conn[curlid]==edg_vert[1])&&(conn[id]==edg_vert[0]))){
*he_fid = curfid;
*he_lid = curlid;
@@ -1336,7 +1412,7 @@ namespace moab {
counter += 1;
}
- delete [] next;
+ delete [] next;
delete [] prev;
return found;
}
@@ -1417,7 +1493,7 @@ namespace moab {
}
/////////////////////////////////////////////
- HalfFacetRep::LocalMaps3D HalfFacetRep::lConnMap3D[4] =
+ const HalfFacetRep::LocalMaps3D HalfFacetRep::lConnMap3D[4] =
{
// Tet
{4, 6, 4, {3,3,3,3}, {{0,1,3},{1,2,3},{2,0,3},{0,2,1}}, {3,3,3,3}, {{0,2,3},{0,1,3},{1,2,3},{0,1,2}}, {{0,1},{1,2},{2,0},{0,3},{1,3},{2,3}}, {{3,0},{3,1},{3,2},{0,2},{0,1},{1,2}}, {{0,4,3},{1,5,4},{2,3,5},{2,1,0}}, {{-1,0,2,3},{0,-1,1,4},{2,1,-1,5},{3,4,5,-1}}},
@@ -1439,16 +1515,6 @@ namespace moab {
{
ErrorCode error;
EntityHandle start_cell = *cells.begin();
-
- int count, *sibhfs_lfid_ptr;
- EntityHandle *sibhfs_cid_ptr;
- error = mb->tag_iterate(sibhfs_cid, cells.begin(), cells.end(), count, (void*&)sibhfs_cid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) cells.size());
- error = mb->tag_iterate(sibhfs_lfid, cells.begin(), cells.end(), count, (void*&)sibhfs_lfid_ptr);
- if (MB_SUCCESS != error) return error;
- assert(count == (int) cells.size());
-
int index = get_index_from_type(start_cell);
int nvpc = lConnMap3D[index].num_verts_in_cell;
int nfpc = lConnMap3D[index].num_faces_in_cell;
@@ -1459,6 +1525,7 @@ namespace moab {
for (int i =0; i<nv+1; i++)
is_index[i] = 0;
+ int vindex;
std::vector<EntityHandle> conn(nvpc);
for (Range::iterator cid = cells.begin(); cid != cells.end(); ++cid)
{
@@ -1476,7 +1543,8 @@ namespace moab {
if (v <= conn[id])
v = conn[id];
}
- is_index[v-*_verts.begin()+1] += 1;
+ vindex = _verts.index(v);
+ is_index[vindex+1] += 1;
}
}
is_index[0] = 0;
@@ -1518,7 +1586,7 @@ namespace moab {
error = local_maps_2d(nvF, next, prev);
if (MB_SUCCESS != error) return error;
- int v = vmax - *_verts.begin();
+ int v = _verts.index(vmax);
v2oe_v1[is_index[v]] = vs[next[lv]];
v2oe_v2[is_index[v]] = vs[prev[lv]];
v2hf_map_cid[is_index[v]] = *cid;
@@ -1542,9 +1610,18 @@ namespace moab {
error = mb->get_connectivity(&*cid, 1, conn);
if (MB_SUCCESS != error) return error;
+ EntityHandle *sibcid = new EntityHandle[nfpc];
+ int *siblfid = new int[nfpc];
+
+ error = mb->tag_get_data(sibhfs_cid, &*cid, 1, sibcid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhfs_lfid, &*cid, 1, siblfid);
+ if (MB_SUCCESS != error) return error;
+
for (int i =0; i<nfpc; i++)
{
- if (sibhfs_cid_ptr[nfpc*(*cid-start_cell)+i] != 0)
+ if (sibcid[i] != 0)
continue;
@@ -1568,31 +1645,60 @@ namespace moab {
error = local_maps_2d(nvF, next, prev);
if (MB_SUCCESS != error) return error;
- int v = vmax - *_verts.begin();
+ int v = _verts.index(vmax);
EntityHandle v1 = vs[prev[lv]];
EntityHandle v2 = vs[next[lv]];
+ EntityHandle *setcid = new EntityHandle[nfpc];
+ int *setlfid = new int[nfpc];
+
for (int ind = is_index[v]; ind <= is_index[v+1]-1; ind++)
{
if ((v2oe_v1[ind] == v1)&&(v2oe_v2[ind] == v2))
{
-
+ // Map to opposite hf
EntityHandle cur_cid = v2hf_map_cid[ind];
int cur_lfid = v2hf_map_lfid[ind];
- sibhfs_cid_ptr[nfpc*(*cid-start_cell)+i] = cur_cid;
- sibhfs_lfid_ptr[nfpc*(*cid-start_cell)+i] = cur_lfid;
+ error = mb->tag_get_data(sibhfs_cid, &*cid, 1, setcid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhfs_lfid, &*cid, 1, setlfid);
+ if (MB_SUCCESS != error) return error;
+
+ setcid[i] = cur_cid;
+ setlfid[i] = cur_lfid;
+
+ error = mb->tag_set_data(sibhfs_cid, &*cid, 1, setcid);
+ if (MB_SUCCESS != error) return error;
+ error = mb->tag_set_data(sibhfs_lfid, &*cid, 1, setlfid);
+ if (MB_SUCCESS != error) return error;
+
+ //Map opposite hf to current cell
+ error = mb->tag_get_data(sibhfs_cid, &cur_cid, 1, setcid);
+ if (MB_SUCCESS != error) return error;
+
+ error = mb->tag_get_data(sibhfs_lfid, &cur_cid, 1, setlfid);
+ if (MB_SUCCESS != error) return error;
- sibhfs_cid_ptr[nfpc*(cur_cid-start_cell)+cur_lfid] = *cid;
- sibhfs_lfid_ptr[nfpc*(cur_cid-start_cell)+cur_lfid] = i;
+ setcid[cur_lfid] = *cid;
+ setlfid[cur_lfid] = i;
+ error = mb->tag_set_data(sibhfs_cid, &cur_cid, 1, setcid);
+ if (MB_SUCCESS != error) return error;
+ error = mb->tag_set_data(sibhfs_lfid, &cur_cid, 1, setlfid);
+ if (MB_SUCCESS != error) return error;
}
}
delete [] next;
delete [] prev;
delete [] vs;
+ delete [] setcid;
+ delete [] setlfid;
}
+ delete sibcid;
+ delete siblfid;
}
@@ -2115,10 +2221,8 @@ namespace moab {
bool HalfFacetRep::find_match_in_array(EntityHandle ent, EntityHandle *ent_list, int count, bool get_index, int *index)
{
bool found = false;
- // std::cout<<"count = "<<count<<std::endl;
for (int i = 0; i<= count; i++)
{
- // std::cout<<"For i = "<<i<<": ent-ent_list[i] = "<<(int)(ent-ent_list[i])<<std::endl;
if (!((int)(ent - ent_list[i])))
{
found = true;
diff --git a/src/moab/Core.hpp b/src/moab/Core.hpp
index 941129a..20aafdc 100644
--- a/src/moab/Core.hpp
+++ b/src/moab/Core.hpp
@@ -328,6 +328,16 @@ public:
#endif
+ /*virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
+ const int num_entities,
+ const int to_dimension,
+ const bool create_if_missing,
+ std::vector<EntityHandle>& adj_entities,
+ const int operation_type = Interface::INTERSECT);*/
+
+
+
+
virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
const int num_entities,
const int to_dimension,
diff --git a/src/moab/HalfFacetRep.hpp b/src/moab/HalfFacetRep.hpp
index 3d1069a..344dbd1 100755
--- a/src/moab/HalfFacetRep.hpp
+++ b/src/moab/HalfFacetRep.hpp
@@ -533,7 +533,7 @@ namespace moab {
int lookup_leids[MAX_VERTICES][MAX_VERTICES];
};
- static LocalMaps3D lConnMap3D[4];
+ static const LocalMaps3D lConnMap3D[4];
int get_index_from_type(EntityHandle cid);
diff --git a/src/moab/Interface.hpp b/src/moab/Interface.hpp
index bfd5fb0..57252e9 100644
--- a/src/moab/Interface.hpp
+++ b/src/moab/Interface.hpp
@@ -630,6 +630,12 @@ public:
#else
) = 0;
#endif
+ /*virtual ErrorCode get_adjacencies(const EntityHandle *from_entities,
+ const int num_entities,
+ const int to_dimension,
+ const bool create_if_missing,
+ std::vector<EntityHandle>& adj_entities,
+ const int operation_type = Interface::INTERSECT) = 0; */
//! Get the adjacencies associated with a vector of entities to entities of a specfied dimension.
/** Identical to vector-based get_adjacencies function, except results are returned in a
diff --git a/test/Makefile.am b/test/Makefile.am
index 34432d0..1354e3c 100644
--- a/test/Makefile.am
+++ b/test/Makefile.am
@@ -61,9 +61,9 @@ if HDF5_FILE
lloyd_smoother_test
endif
-if ENABLE_ahf
- TESTS += ahf_mbintf_test
-endif
+# if ENABLE_ahf
+# TESTS += ahf_mbintf_test
+# endif
# merge_test \ # input files no longer exist?
# test_tag_server \ # fails
@@ -158,9 +158,9 @@ lloyd_smoother_test_SOURCES = lloyd_smoother_test.cpp
ahf_test_SOURCES = test_ahf_moab.cpp
-if ENABLE_ahf
-ahf_mbintf_test_SOURCES = test_ahf_mb_interface.cpp
-endif
+# if ENABLE_ahf
+# ahf_mbintf_test_SOURCES = test_ahf_mb_interface.cpp
+# endif
if PARALLEL
diff --git a/test/perf/Makefile.am b/test/perf/Makefile.am
index f946f94..32e2a57 100644
--- a/test/perf/Makefile.am
+++ b/test/perf/Makefile.am
@@ -7,13 +7,14 @@ AM_CPPFLAGS += -I$(top_srcdir)/src \
LDADD = $(top_builddir)/src/libMOAB.la
-check_PROGRAMS = perf seqperf adj_time perftool
+check_PROGRAMS = perf seqperf adj_time perftool ahf_mem_time
noinst_PROGRAMS =
perf_SOURCES = perf.cpp
seqperf_SOURCES = seqperf.cpp
adj_time_SOURCES = adj_time.cpp
perftool_SOURCES = perftool.cpp
+ahf_mem_time_SOURCES = test_ahf_mem_time.cpp
if ENABLE_imesh
LDADD += $(top_builddir)/itaps/imesh/libiMesh.la
diff --git a/test/perf/test_ahf_mem_time.cpp b/test/perf/test_ahf_mem_time.cpp
new file mode 100755
index 0000000..9c72ae7
--- /dev/null
+++ b/test/perf/test_ahf_mem_time.cpp
@@ -0,0 +1,281 @@
+/*This function tests the AHF datastructures on CST meshes*/
+#include <iostream>
+#include <assert.h>
+#include <time.h>
+#include <vector>
+#include "moab/Core.hpp"
+#include "moab/Range.hpp"
+#include "moab/MeshTopoUtil.hpp"
+#include "moab/HalfFacetRep.hpp"
+#include <sys/time.h>
+
+using namespace moab;
+
+#ifdef MESHDIR
+std::string TestDir(STRINGIFY(MESHDIR));
+#else
+std::string TestDir(".");
+#endif
+
+std::string filename;
+
+double wtime() {
+ double y = -1;
+ struct timeval cur_time;
+ gettimeofday(&cur_time, NULL);
+ y = (double)(cur_time.tv_sec) + (double)(cur_time.tv_usec)*1.e-6;
+ return (y);
+}
+
+int main(int argc, char **argv)
+{
+ // Read the input mesh
+ filename = TestDir + "/hexes_mixed.vtk";
+
+ if (argc==1)
+ std::cout<<"Using default input file:"<<filename<<std::endl;
+ else if (argc==2)
+ filename = argv[1];
+ else {
+ std::cerr << "Usage: " << argv[0] << " [filename]" << std::endl;
+ return 1;
+ }
+
+ ErrorCode error;
+ Core moab;
+ Interface* mbImpl = &moab;
+ MeshTopoUtil mtu(mbImpl);
+
+ error = mbImpl->load_file( filename.c_str());
+ if (MB_SUCCESS != error) {
+ std::cerr << filename <<": failed to load file." << std::endl;
+ return error;
+ }
+
+ //Create ranges for handles of explicit elements of the mixed mesh
+ Range verts, edges, faces, cells;
+ error = mbImpl->get_entities_by_dimension( 0, 0, verts);
+ error = mbImpl->get_entities_by_dimension( 0, 1, edges);
+ error = mbImpl->get_entities_by_dimension( 0, 2, faces);
+ error = mbImpl->get_entities_by_dimension( 0, 3, cells);
+
+ int nverts = verts.size();
+ int nedges = edges.size();
+ int nfaces = faces.size();
+ int ncells = cells.size();
+
+ std::cout<<"nverts = "<<nverts<<", nedges = "<<nedges<<", nfaces = "<<nfaces<<", ncells = "<<ncells<<std::endl;
+
+
+ //Storage Costs before calling ahf functionalities
+ unsigned long sTotS, sTAS, sES, sAES, sAS, sAAS, sTS, sATS;
+ sTotS = sTAS = sES = sAES = sAS = sAAS = sTS = sATS = 0;
+ mbImpl->estimated_memory_use(NULL, 0, &sTotS, &sTAS, &sES, &sAES, &sAS, &sAAS, NULL, 0, &sTS, &sATS);
+ std::cout<<std::endl;
+ std::cout<<"Total storage = "<<sTotS<<std::endl;
+ std::cout<<"Total amortized storage = "<<sTAS<<std::endl;
+ std::cout<<"Entity storage = "<<sES<<std::endl;
+ std::cout<<"Amortized entity storage = "<<sAES<<std::endl;
+ std::cout<<"Adjacency storage = "<<sAS<<std::endl;
+ std::cout<<"Amortized adjacency storage = "<<sAAS<<std::endl;
+ std::cout<<"Tag storage = "<<sTS<<std::endl;
+ std::cout<<"Amortized tag storage = "<<sATS<<std::endl;
+ std::cout<<std::endl;
+
+
+ double time_start, time_elapsed, time_avg;
+
+ // Create an ahf instance
+ HalfFacetRep ahf(&moab);
+
+ // Call the initialize function which creates the maps for each dimension
+ time_start = wtime();
+
+ ahf.initialize();
+
+ time_elapsed = wtime() - time_start;
+ std::cout << "Time taken to construct the MDS = "<<time_elapsed<<" secs"<<std::endl;
+
+ //Storage Costs after calling ahf initialize
+ unsigned long TotS, TAS, ES, AES, AS, AAS, TS, ATS;
+ TotS = TAS = ES = AES = AS = AAS = TS = ATS = 0;
+ mbImpl->estimated_memory_use(NULL, 0, &TotS, &TAS, &ES, &AES, &AS, &AAS, NULL, 0, &TS, &ATS);
+ std::cout<<std::endl;
+ std::cout<<"Total storage = "<<TotS<<std::endl;
+ std::cout<<"Total amortized storage = "<<TAS<<std::endl;
+ std::cout<<"Entity storage = "<<ES<<std::endl;
+ std::cout<<"Amortized entity storage = "<<AES<<std::endl;
+ std::cout<<"Adjacency storage = "<<AS<<std::endl;
+ std::cout<<"Amortized adjacency storage = "<<AAS<<std::endl;
+ std::cout<<"Tag storage = "<<TS<<std::endl;
+ std::cout<<"Amortized tag storage = "<<ATS<<std::endl;
+ std::cout<<std::endl;
+
+ //Perform queries
+ std::vector<EntityHandle> adjents;
+ std::vector<int> lids;
+ Range mbents;
+
+ //1D Queries //
+ //IQ1: For every vertex, obtain incident edges
+ std::cout<<"1D QUERIES"<<std::endl;
+ time_start = wtime();
+ for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 1, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)verts.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute incident edges to a vertex = "<< time_avg<<" secs" <<std::endl;
+
+ error = mbImpl->get_adjacencies( &*verts.begin(), 1, 1, false, mbents );
+ time_start = wtime();
+ for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 1, false, mbents );
+ }
+ time_avg = (wtime()-time_start)/(double)verts.size();
+ std::cout<<"MOAB: Average time taken to compute incident edges to a vertex = "<< time_avg<<" secs" <<std::endl;
+ std::cout<<std::endl;
+
+
+ //NQ1: For every edge, obtain neighbor edges
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute neighbor edges of an edge = "<<time_avg<<" secs" << std::endl;
+
+ error = mtu.get_bridge_adjacencies( *edges.begin(), 0, 1, mbents);
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 0, 1, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB: Average time taken to compute neighbor edges of an edge = "<<time_avg<<" secs" << std::endl;
+ std::cout<<std::endl;
+
+
+ // 2D Queries
+ //IQ21: For every edge, obtain incident faces
+ std::cout<<"2D QUERIES"<<std::endl;
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 2, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute incident faces on an edge = "<<time_avg<<" secs" <<std::endl;
+
+ error = mbImpl->get_adjacencies( &*edges.begin(), 1, 2, false, mbents);
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 2, false, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB: Average time taken to compute incident faces on an edge = "<<time_avg<<" secs" <<std::endl;
+ std::cout<<std::endl;
+
+ //NQ2: For every face, obtain neighbor faces
+ time_start = wtime();
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)faces.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute neighbor faces of a face = "<< time_avg<<" secs" <<std::endl;
+
+ error = mtu.get_bridge_adjacencies( *faces.begin(), 1, 2, mbents);
+ time_start = wtime();
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 1, 2, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)faces.size();
+ std::cout<<"MOAB: Average time taken to compute neighbor faces of a face = "<< time_avg<<" secs" <<std::endl;
+ std::cout<<std::endl;
+
+
+ // 3D Queries
+ // IQ 31: For every edge, obtain incident cells
+ std::cout<<"3D QUERIES"<<std::endl;
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute incident cells on an edge = "<<time_avg <<" secs"<<std::endl;
+
+ error = mbImpl->get_adjacencies(&*edges.begin(), 1, 3, false, mbents);
+ time_start = wtime();
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)edges.size();
+ std::cout<<"MOAB: Average time taken to compute incident cells on an edge = "<<time_avg <<" secs"<<std::endl;
+ std::cout<<std::endl;
+
+ //IQ32: For every face, obtain incident cells
+ time_start = wtime();
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
+
+ }
+ time_avg = (wtime()-time_start)/(double)faces.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute incident cells on a face = "<<time_avg <<" secs"<<std::endl;
+
+ error = mbImpl->get_adjacencies(&*faces.begin(), 1, 3, false, mbents);
+ time_start = wtime();
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)faces.size();
+ std::cout<<"MOAB: Average time taken to compute incident cells on a face = "<<time_avg <<" secs"<<std::endl;
+ std::cout<<std::endl;
+
+ //NQ3: For every cell, obtain neighbor cells
+ time_start = wtime();
+ for (Range::iterator i = cells.begin(); i != cells.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ }
+ time_avg = (wtime()-time_start)/(double)cells.size();
+ std::cout<<"MOAB_AHF: Average time taken to compute neighbor cells of a cell = "<< time_avg <<" secs" << std::endl;
+
+ error = mtu.get_bridge_adjacencies( *cells.begin(), 2, 3, mbents);
+ time_start = wtime();
+ for (Range::iterator i = cells.begin(); i != cells.end(); ++i) {
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 2, 3, mbents);
+ }
+ time_avg = (wtime()-time_start)/(double)cells.size();
+ std::cout<<"MOAB: Average time taken to compute neighbor cells of a cell = "<< time_avg <<" secs" << std::endl;
+ std::cout<<std::endl;
+
+ ahf.deinitialize();
+
+ //Storage Costs after calling ahf deinitialize
+ unsigned long eTotS, eTAS, eES, eAES, eAS, eAAS, eTS, eATS;
+ eTotS = eTAS = eES = eAES = eAS = eAAS = eTS = eATS = 0;
+ mbImpl->estimated_memory_use(NULL, 0, &eTotS, &eTAS, &eES, &eAES, &eAS, &eAAS, NULL, 0, &eTS, &eATS);
+ std::cout<<std::endl;
+ std::cout<<"Total storage = "<<eTotS<<std::endl;
+ std::cout<<"Total amortized storage = "<<eTAS<<std::endl;
+ std::cout<<"Entity storage = "<<eES<<std::endl;
+ std::cout<<"Amortized entity storage = "<<eAES<<std::endl;
+ std::cout<<"Adjacency storage = "<<eAS<<std::endl;
+ std::cout<<"Amortized adjacency storage = "<<eAAS<<std::endl;
+ std::cout<<"Tag storage = "<<eTS<<std::endl;
+ std::cout<<"Amortized tag storage = "<<eATS<<std::endl;
+ std::cout<<std::endl;
+
+ return 0;
+}
+
diff --git a/test/test_ahf_moab.cpp b/test/test_ahf_moab.cpp
index 93bd947..e415186 100755
--- a/test/test_ahf_moab.cpp
+++ b/test/test_ahf_moab.cpp
@@ -11,172 +11,203 @@
using namespace moab;
#ifdef MESHDIR
-static const char example[] = STRINGIFY(MESHDIR) "/hexes_mixed.vtk";
+std::string TestDir(STRINGIFY(MESHDIR));
#else
-static const char example[] = "/hexes_mixed.vtk";
+std::string TestDir(".");
#endif
-void ahf_test()
-{
-
- Core moab;
- Interface* mbImpl = &moab;
- MeshTopoUtil mtu(mbImpl);
-
- ErrorCode error = mbImpl->load_file(example);
- CHECK_ERR(error);
-
- /*Create ranges for handles of explicit elements of the mixed mesh*/
- Range verts, edges, faces, cells;
- error = mbImpl->get_entities_by_dimension( 0, 0, verts);
- error = mbImpl->get_entities_by_dimension( 0, 1, edges);
- error = mbImpl->get_entities_by_dimension( 0, 2, faces);
- error = mbImpl->get_entities_by_dimension( 0, 3, cells);
-
- // Create an ahf instance
- HalfFacetRep ahf(&moab);
-
- // Call the initialize function which creates the maps for each dimension
- ahf.initialize();
-
- //Perform queries
- std::vector<EntityHandle> adjents;
- Range mbents, ahfents;
-
- //1D Queries //
- //IQ1: For every vertex, obtain incident edges
- for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
- adjents.clear();
- error = ahf.get_up_adjacencies( *i, 1, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mbImpl->get_adjacencies( &*i, 1, 1, false, mbents );
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(),mbents.size());
-
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
-
- //NQ1: For every edge, obtain neighbor edges
- for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
- adjents.clear();
- error = ahf.get_neighbor_adjacencies( *i, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mtu.get_bridge_adjacencies( *i, 0, 1, mbents);
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(), mbents.size());
-
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
-
- // 2D Queries
- //IQ2: For every edge, obtain incident faces
- for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
- adjents.clear();
- error = ahf.get_up_adjacencies( *i, 2, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mbImpl->get_adjacencies( &*i, 1, 2, false, mbents);
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(), mbents.size());
+std::string filename;
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
+int number_tests_successful = 0;
+int number_tests_failed = 0;
- //NQ2: For every face, obtain neighbor faces
- for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
- adjents.clear();
- error = ahf.get_neighbor_adjacencies( *i, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mtu.get_bridge_adjacencies( *i, 1, 2, mbents);
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(), mbents.size());
-
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
+void handle_error_code(ErrorCode rv, int &number_failed, int &number_successful)
+{
+ if (rv == MB_SUCCESS) {
+ std::cout << "Success";
+ number_successful++;
+ } else {
+ std::cout << "Failure";
+ number_failed++;
}
+}
- // 3D Queries
- // IQ 31: For every edge, obtain incident cells
- for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
- adjents.clear();
- error = ahf.get_up_adjacencies( *i, 3, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(), mbents.size());
-
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
- //IQ32: For every face, obtain incident cells
- for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
- adjents.clear();
- error = ahf.get_up_adjacencies( *i, 3, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
- CHECK_ERR(error);
-
- CHECK_EQUAL(adjents.size(), mbents.size());
+ErrorCode ahf_test(Core *moab)
+{
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
+ Interface* mbImpl = &*moab;
+ MeshTopoUtil mtu(mbImpl);
- //NQ3: For every cell, obtain neighbor cells
- for (Range::iterator i = cells.begin(); i != cells.end(); ++i) {
- adjents.clear();
- error = ahf.get_neighbor_adjacencies( *i, adjents);
- CHECK_ERR(error);
- mbents.clear();
- error = mtu.get_bridge_adjacencies( *i, 2, 3, mbents);
+ ErrorCode error = mbImpl->load_file(filename.c_str());
CHECK_ERR(error);
- CHECK_EQUAL(adjents.size(), mbents.size());
-
- std::sort(adjents.begin(), adjents.end());
- std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
- mbents = subtract(mbents, ahfents);
- CHECK(!mbents.size());
- }
-
- //ahf.deinitialize();
+ /*Create ranges for handles of explicit elements of the mixed mesh*/
+ Range verts, edges, faces, cells;
+ error = mbImpl->get_entities_by_dimension( 0, 0, verts);
+ error = mbImpl->get_entities_by_dimension( 0, 1, edges);
+ error = mbImpl->get_entities_by_dimension( 0, 2, faces);
+ error = mbImpl->get_entities_by_dimension( 0, 3, cells);
+
+ // Create an ahf instance
+ HalfFacetRep ahf(&*moab);
+
+ // Call the initialize function which creates the maps for each dimension
+ ahf.initialize();
+
+ //ahf.print_tags();
+
+ //Perform queries
+ std::vector<EntityHandle> adjents;
+ Range mbents, ahfents;
+
+ //1D Queries //
+ //IQ1: For every vertex, obtain incident edges
+ for (Range::iterator i = verts.begin(); i != verts.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 1, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 1, false, mbents );
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(),mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ1: For every edge, obtain neighbor edges
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 0, 1, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ // 2D Queries
+ //IQ2: For every edge, obtain incident faces
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 2, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies( &*i, 1, 2, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ2: For every face, obtain neighbor faces
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 1, 2, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ // 3D Queries
+ // IQ 31: For every edge, obtain incident cells
+ for (Range::iterator i = edges.begin(); i != edges.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //IQ32: For every face, obtain incident cells
+ for (Range::iterator i = faces.begin(); i != faces.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_up_adjacencies( *i, 3, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mbImpl->get_adjacencies(&*i, 1, 3, false, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ //NQ3: For every cell, obtain neighbor cells
+ for (Range::iterator i = cells.begin(); i != cells.end(); ++i) {
+ adjents.clear();
+ error = ahf.get_neighbor_adjacencies( *i, adjents);
+ CHECK_ERR(error);
+ mbents.clear();
+ error = mtu.get_bridge_adjacencies( *i, 2, 3, mbents);
+ CHECK_ERR(error);
+
+ CHECK_EQUAL(adjents.size(), mbents.size());
+
+ std::sort(adjents.begin(), adjents.end());
+ std::copy(adjents.begin(), adjents.end(), range_inserter(ahfents));
+ mbents = subtract(mbents, ahfents);
+ CHECK(!mbents.size());
+ }
+
+ return MB_SUCCESS;
}
int main(int argc, char *argv[])
{
- int result = 0;
-
- argv[0] = argv[argc - argc]; // Followed read_mpas_nc.cpp test for removing warnings in serial mode about unused variables.
-
- result += RUN_TEST(ahf_test);
-
- return result;
+ filename = TestDir + "/hexes_mixed.vtk";
+
+ if (argc==1)
+ std::cout<<"Using default input file:"<<filename<<std::endl;
+ else if (argc==2)
+ filename = argv[1];
+ else {
+ std::cerr << "Usage: " << argv[0] << " [filename]" << std::endl;
+ return 1;
+ }
+
+ Core moab;
+ ErrorCode result;
+
+ std::cout<<" ahf_test: ";
+ result = ahf_test(&moab);
+ handle_error_code(result, number_tests_failed, number_tests_successful);
+ std::cout<<"\n";
+
+ return number_tests_failed;
}
Repository URL: https://bitbucket.org/fathomteam/moab/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the moab-dev
mailing list