|
MOAB: Mesh Oriented datABase
(version 5.4.1)
|
|
MOAB: Mesh Oriented datABase
(version 5.4.1)
|
#include "moab/Core.hpp"#include "moab/AdaptiveKDTree.hpp"#include "moab/Range.hpp"#include "moab/CartVect.hpp"#include "moab/GeomUtil.hpp"#include <iostream>#include <ctime>#include <cstdlib>#include <cassert>#include <sstream>
Include dependency graph for point_location.cpp:Go to the source code of this file.
Defines | |
| #define | CHK(ErrorCode) |
Enumerations | |
| enum | TreeType { UseKDTree, UseNoTree, UseDefaultTree = UseKDTree } |
Functions | |
| static void | fail (ErrorCode error_code, const char *file_name, int line_number) |
| const char * | is_default_tree (TreeType type) |
| static void | usage (char *argv0, bool help=false) |
| static void | generate_random_points (Interface &mesh, size_t num_points, std::vector< CartVect > &points, std::vector< EntityHandle > &point_elems) |
| static void | do_kdtree_test (Interface &mesh, int tree_depth, int elem_per_leaf, long num_test, const std::vector< CartVect > &points, std::vector< EntityHandle > &point_elems, clock_t &build_time, clock_t &test_time, size_t &depth) |
| static void | do_linear_test (Interface &mesh, int tree_depth, int elem_per_leaf, long num_test, const std::vector< CartVect > &points, std::vector< EntityHandle > &point_elems, clock_t &build_time, clock_t &test_time, size_t &depth) |
| int | main (int argc, char *argv[]) |
| static CartVect | random_point_in_hex (Interface &mb, EntityHandle hex) |
Variables | |
| const char * | default_str = "(default)" |
| const char * | empty_str = "" |
| const long | DEFAULT_NUM_TEST = 100000 |
| const long | HARD_MAX_UNIQUE_POINTS = 100000 |
| const long | HARD_MIN_UNIQUE_POINTS = 1000 |
| const long | FRACTION_UNIQUE_POINTS = 100 |
| const int | HexSign [8][3] |
do \ { \ if( MB_SUCCESS != ( ErrorCode ) ) fail( ( ErrorCode ), __FILE__, __LINE__ ); \ } while( false )
Definition at line 12 of file point_location.cpp.
Referenced by do_kdtree_test(), do_linear_test(), generate_random_points(), and random_point_in_hex().
| enum TreeType |
Definition at line 22 of file point_location.cpp.
{
UseKDTree,
UseNoTree,
UseDefaultTree = UseKDTree
};
| void do_kdtree_test | ( | Interface & | mesh, |
| int | tree_depth, | ||
| int | elem_per_leaf, | ||
| long | num_test, | ||
| const std::vector< CartVect > & | points, | ||
| std::vector< EntityHandle > & | point_elems, | ||
| clock_t & | build_time, | ||
| clock_t & | test_time, | ||
| size_t & | depth | ||
| ) | [static] |
Definition at line 310 of file point_location.cpp.
References build_time, moab::AdaptiveKDTree::build_tree(), CHK, moab::Range::clear(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_handle(), moab::Interface::get_entities_by_type(), moab::AdaptiveKDTree::get_info(), init(), leaf, MBHEX, moab::GeomUtil::point_in_trilinear_hex(), and moab::AdaptiveKDTree::point_search().
Referenced by main().
{
ErrorCode rval;
clock_t init = clock();
AdaptiveKDTree tool( &mb );
EntityHandle root;
std::ostringstream options;
if( tree_depth > 0 ) options << "MAX_DEPTH=" << tree_depth << ";";
if( elem_per_leaf > 0 ) options << "MAX_PER_LEAF=" << elem_per_leaf << ";";
Range all_hexes;
rval = mb.get_entities_by_type( 0, MBHEX, all_hexes );
CHK( rval );
FileOptions opt( options.str().c_str() );
rval = tool.build_tree( all_hexes, &root, &opt );
CHK( rval );
all_hexes.clear();
build_time = clock() - init;
EntityHandle leaf;
std::vector< EntityHandle > hexes;
std::vector< EntityHandle >::iterator j;
CartVect coords[8];
const EntityHandle* conn;
int len;
for( long i = 0; i < num_test; ++i )
{
const size_t idx = (size_t)i % points.size();
rval = tool.point_search( points[idx].array(), leaf, 1.0e-10, 1.0e-6, NULL, &root );
CHK( rval );
hexes.clear();
rval = mb.get_entities_by_handle( leaf, hexes );
CHK( rval );
for( j = hexes.begin(); j != hexes.end(); ++j )
{
rval = mb.get_connectivity( *j, conn, len, true );
CHK( rval );
rval = mb.get_coords( conn, 8, reinterpret_cast< double* >( coords ) );
CHK( rval );
if( GeomUtil::point_in_trilinear_hex( coords, points[idx], 1e-12 ) )
{
point_elems[idx] = *j;
break;
}
}
if( j == hexes.end() ) point_elems[idx] = 0;
}
test_time = clock() - build_time;
double tmp_box[3];
unsigned max_d;
tool.get_info( root, tmp_box, tmp_box, max_d );
depth = max_d;
}
| void do_linear_test | ( | Interface & | mesh, |
| int | tree_depth, | ||
| int | elem_per_leaf, | ||
| long | num_test, | ||
| const std::vector< CartVect > & | points, | ||
| std::vector< EntityHandle > & | point_elems, | ||
| clock_t & | build_time, | ||
| clock_t & | test_time, | ||
| size_t & | depth | ||
| ) | [static] |
Definition at line 372 of file point_location.cpp.
References moab::Range::begin(), build_time, CHK, moab::Range::end(), ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), moab::Interface::get_entities_by_type(), init(), MBHEX, and moab::GeomUtil::point_in_trilinear_hex().
Referenced by main().
{
clock_t init = clock();
Range hexes;
ErrorCode rval = mb.get_entities_by_type( 0, MBHEX, hexes );
CHK( rval );
depth = 0;
point_elems.resize( points.size() );
build_time = clock() - init;
CartVect coords[8];
const EntityHandle* conn;
int len;
for( long i = 0; i < num_test; ++i )
{
const size_t idx = (size_t)i % points.size();
for( Range::iterator h = hexes.begin(); h != hexes.end(); ++h )
{
rval = mb.get_connectivity( *h, conn, len, true );
CHK( rval );
rval = mb.get_coords( conn, 8, reinterpret_cast< double* >( coords ) );
CHK( rval );
if( GeomUtil::point_in_trilinear_hex( coords, points[idx], 1e-12 ) )
{
point_elems[idx] = *h;
break;
}
}
}
test_time = clock() - build_time;
}
Definition at line 236 of file point_location.cpp.
{
std::cerr << "Internal error (error code " << error_code << ") at " << file << ":" << line << std::endl;
abort();
}
| void generate_random_points | ( | Interface & | mesh, |
| size_t | num_points, | ||
| std::vector< CartVect > & | points, | ||
| std::vector< EntityHandle > & | point_elems | ||
| ) | [static] |
Definition at line 274 of file point_location.cpp.
References moab::Range::all_of_type(), moab::Range::begin(), CHK, moab::Range::empty(), moab::Range::end(), moab::Range::equal_range(), moab::Range::erase(), ErrorCode, moab::Interface::get_entities_by_dimension(), MBHEX, random_point_in_hex(), and moab::Range::size().
Referenced by main().
{
Range elems;
ErrorCode rval;
rval = mb.get_entities_by_dimension( 0, 3, elems );
CHK( rval );
if( !elems.all_of_type( MBHEX ) )
{
std::cerr << "Warning: ignoring non-hexahedral elements." << std::endl;
std::pair< Range::iterator, Range::iterator > p = elems.equal_range( MBHEX );
elems.erase( p.second, elems.end() );
elems.erase( elems.begin(), p.first );
}
if( elems.empty() )
{
std::cerr << "Input file contains no hexahedral elements." << std::endl;
exit( 1 );
}
points.resize( num_points );
point_elems.resize( num_points );
const size_t num_elem = elems.size();
for( size_t i = 0; i < num_points; ++i )
{
size_t offset = 0;
for( size_t x = num_elem; x > 0; x /= RAND_MAX )
offset += rand();
offset %= num_elem;
point_elems[i] = elems[offset];
points[i] = random_point_in_hex( mb, point_elems[i] );
}
}
| const char* is_default_tree | ( | TreeType | type | ) | [inline] |
Definition at line 31 of file point_location.cpp.
References default_str, empty_str, and UseDefaultTree.
Referenced by usage().
{
return type == UseDefaultTree ? default_str : empty_str;
}
| int main | ( | int | argc, |
| char * | argv[] | ||
| ) |
Definition at line 87 of file point_location.cpp.
References build_time, DEFAULT_NUM_TEST, do_kdtree_test(), do_linear_test(), ErrorCode, moab::Interface::estimated_memory_use(), moab::fail(), FRACTION_UNIQUE_POINTS, generate_random_points(), moab::Interface::get_last_error(), HARD_MAX_UNIQUE_POINTS, HARD_MIN_UNIQUE_POINTS, input_file, moab::Interface::load_file(), mb, MB_SUCCESS, usage, UseDefaultTree, UseKDTree, and UseNoTree.
{
const char* input_file = 0;
long tree_depth = -1;
long elem_per_leaf = -1;
long num_points = DEFAULT_NUM_TEST;
TreeType type = UseDefaultTree;
char* endptr;
// PARSE ARGUMENTS
long* valptr = 0;
for( int i = 1; i < argc; ++i )
{
if( valptr )
{
*valptr = strtol( argv[i], &endptr, 0 );
if( *valptr < 1 || *endptr )
{
std::cerr << "Invalid value for '" << argv[i - 1] << "' flag: " << argv[i] << std::endl;
exit( 1 );
}
valptr = 0;
}
else if( argv[i][0] == '-' && argv[i][1] && !argv[i][2] )
{
switch( argv[i][1] )
{
case 'h':
usage( argv[0], true );
break;
case 'k':
type = UseKDTree;
break;
case 'v':
type = UseNoTree;
break;
case 'd':
valptr = &tree_depth;
break;
case 'e':
valptr = &elem_per_leaf;
break;
case 'n':
valptr = &num_points;
break;
default:
std::cerr << "Invalid flag: " << argv[i] << std::endl;
usage( argv[0] );
break;
}
}
else if( !input_file )
{
input_file = argv[i];
}
else
{
std::cerr << "Unexpected argument: " << argv[i] << std::endl;
usage( argv[0] );
}
}
if( valptr )
{
std::cerr << "Expected value following '" << argv[argc - 1] << "' flag" << std::endl;
usage( argv[0] );
}
if( !input_file )
{
std::cerr << "No input file specified." << std::endl;
usage( argv[0] );
}
// LOAD MESH
Core moab;
Interface& mb = moab;
ErrorCode rval;
std::string init_msg, msg;
mb.get_last_error( init_msg );
rval = mb.load_file( input_file );
if( MB_SUCCESS != rval )
{
std::cerr << input_file << " : failed to read file." << std::endl;
mb.get_last_error( msg );
if( msg != init_msg ) std::cerr << "message : " << msg << std::endl;
}
// GENERATE TEST POINTS
int num_unique = num_points / FRACTION_UNIQUE_POINTS;
if( num_unique > HARD_MAX_UNIQUE_POINTS )
num_unique = HARD_MAX_UNIQUE_POINTS;
else if( num_unique < HARD_MIN_UNIQUE_POINTS )
num_unique = num_points;
std::vector< CartVect > points;
std::vector< EntityHandle > elems;
generate_random_points( mb, num_unique, points, elems );
// GET MEMORY USE BEFORE BUILDING TREE
unsigned long long init_total_storage;
mb.estimated_memory_use( 0, 0, &init_total_storage );
// RUN TIMING TEST
clock_t build_time, test_time;
size_t actual_depth = 0;
std::vector< EntityHandle > results( points.size() );
switch( type )
{
default:
case UseKDTree:
do_kdtree_test( mb, tree_depth, elem_per_leaf, num_points, points, results, build_time, test_time,
actual_depth );
break;
case UseNoTree:
do_linear_test( mb, tree_depth, elem_per_leaf, num_points, points, results, build_time, test_time,
actual_depth );
break;
}
unsigned long long fini_total_storage;
mb.estimated_memory_use( 0, 0, &fini_total_storage );
// VALIDATE RESULTS
int fail = 0;
if( results != elems )
{
std::cout << "WARNING: Test produced invalid results" << std::endl;
fail = 1;
}
// SUMMARIZE RESULTS
std::cout << "Number of test queries: " << num_points << std::endl;
std::cout << "Tree build time: " << ( (double)build_time ) / CLOCKS_PER_SEC << " seconds" << std::endl;
std::cout << "Total query time: " << ( (double)test_time ) / CLOCKS_PER_SEC << " seconds" << std::endl;
std::cout << "Time per query: " << ( (double)test_time ) / CLOCKS_PER_SEC / num_points << " seconds" << std::endl;
std::cout << "Tree depth: " << actual_depth << std::endl;
if( actual_depth )
std::cout << "Total query time/depth: " << ( (double)test_time ) / CLOCKS_PER_SEC / actual_depth << " seconds"
<< std::endl;
std::cout << std::endl;
std::cout << "Bytes before tree construction: " << init_total_storage << std::endl;
std::cout << "Bytes after tree construction: " << fini_total_storage << std::endl;
std::cout << "Difference: " << fini_total_storage - init_total_storage << " bytes" << std::endl;
return fail;
}
| static CartVect random_point_in_hex | ( | Interface & | mb, |
| EntityHandle | hex | ||
| ) | [static] |
Definition at line 245 of file point_location.cpp.
References CHK, ErrorCode, moab::Interface::get_connectivity(), moab::Interface::get_coords(), HexSign, and MB_SUCCESS.
Referenced by generate_random_points().
{
const double f = RAND_MAX / 2;
CartVect xi( ( (double)rand() - f ) / f, ( (double)rand() - f ) / f, ( (double)rand() - f ) / f );
CartVect coords[8];
const EntityHandle* conn;
int len;
ErrorCode rval = mb.get_connectivity( hex, conn, len, true );
if( len != 8 && MB_SUCCESS != rval )
{
std::cerr << "Invalid element" << std::endl;
assert( false );
abort();
}
rval = mb.get_coords( conn, 8, reinterpret_cast< double* >( coords ) );
CHK( rval );
CartVect point( 0, 0, 0 );
for( unsigned i = 0; i < 8; ++i )
{
double coeff = 0.125;
for( unsigned j = 0; j < 3; ++j )
coeff *= 1 + HexSign[i][j] * xi[j];
point += coeff * coords[i];
}
return point;
}
| static void usage | ( | char * | argv0, |
| bool | help = false |
||
| ) | [static] |
Definition at line 41 of file point_location.cpp.
References DEFAULT_NUM_TEST, help(), is_default_tree(), UseKDTree, and UseNoTree.
{
const char* usage_str = "[-k|-v] [-n <N>] [-d <N>] [-e <N>] <input_mesh>";
if( !help )
{
std::cerr << "Usage: " << argv0 << " " << usage_str << std::endl;
std::cerr << " : " << argv0 << " -h" << std::endl;
exit( 1 );
}
std::cout << "Usage: " << argv0 << " " << usage_str << std::endl
<< " -k : Use kD Tree " << is_default_tree( UseKDTree ) << std::endl
<< " -v : Use no tree" << is_default_tree( UseNoTree ) << std::endl
<< " -n : Specify number of test points (default = " << DEFAULT_NUM_TEST << ")" << std::endl
<< " -d : Specify maximum tree depth" << std::endl
<< " -e : Specify maximum elements per leaf" << std::endl
<< " <input_mesh> : Mesh to build and query." << std::endl
<< std::endl;
exit( 0 );
}
| const long DEFAULT_NUM_TEST = 100000 |
Definition at line 36 of file point_location.cpp.
| const char* default_str = "(default)" |
Definition at line 29 of file point_location.cpp.
Referenced by is_default_tree().
| const char* empty_str = "" |
Definition at line 30 of file point_location.cpp.
Referenced by is_default_tree().
| const long FRACTION_UNIQUE_POINTS = 100 |
Definition at line 39 of file point_location.cpp.
Referenced by main().
| const long HARD_MAX_UNIQUE_POINTS = 100000 |
Definition at line 37 of file point_location.cpp.
Referenced by main().
| const long HARD_MIN_UNIQUE_POINTS = 1000 |
Definition at line 38 of file point_location.cpp.
Referenced by main().
| const int HexSign[8][3] |
{ { -1, -1, -1 }, { 1, -1, -1 }, { 1, 1, -1 }, { -1, 1, -1 },
{ -1, -1, 1 }, { 1, -1, 1 }, { 1, 1, 1 }, { -1, 1, 1 } }
Definition at line 242 of file point_location.cpp.
Referenced by random_point_in_hex().
1.7.6.1