MOAB: Mesh Oriented datABase  (version 5.2.1)
tree_searching_perf.cpp
Go to the documentation of this file.
00001 #include "moab/Core.hpp"
00002 #include "moab/SpatialLocator.hpp"
00003 #include "moab/Tree.hpp"
00004 #include "moab/HomXform.hpp"
00005 #include "moab/ScdInterface.hpp"
00006 #include "moab/CartVect.hpp"
00007 #include "moab/AdaptiveKDTree.hpp"
00008 #include "moab/BVHTree.hpp"
00009 #include "moab/ProgOptions.hpp"
00010 #include "moab/CpuTimer.hpp"
00011 
00012 #ifdef MOAB_HAVE_MPI
00013 #include "moab_mpi.h"
00014 #endif
00015 
00016 #include <cstdlib>
00017 #include <sstream>
00018 
00019 using namespace moab;
00020 
00021 ErrorCode test_locator( SpatialLocator& sl, int npoints, double& cpu_time, double& percent_outside );
00022 ErrorCode create_hex_mesh( Interface& mb, Range& elems, int n, int dim );
00023 
00024 int main( int argc, char** argv )
00025 {
00026 #ifdef MOAB_HAVE_MPI
00027     int fail = MPI_Init( &argc, &argv );
00028     if( fail ) return fail;
00029 #else
00030     // silence the warning of parameters not used, in serial; there should be a smarter way :(
00031     argv[0] = argv[argc - argc];
00032 #endif
00033 
00034     int npoints = 100, dim = 3;
00035     int dints = 1, dleafs = 1, ddeps = 1, csints = 0;
00036 
00037     ProgOptions po;
00038     po.addOpt< int >( "candidateplaneset,c",
00039                       "Candidate plane set (0=SUBDIVISION,1=SUBDIV_SNAP,2=VERTEX_MEDIAN,3=VERTEX_SAMPLE", &csints );
00040     po.addOpt< int >( "ints,i", "Number of doublings of intervals on each side of scd mesh", &dints );
00041     po.addOpt< int >( "leaf,l", "Number of doublings of maximum number of elements per leaf", &dleafs );
00042     po.addOpt< int >( "max_depth,m", "Number of 5-intervals on maximum depth of tree", &ddeps );
00043     po.addOpt< int >( "npoints,n", "Number of query points", &npoints );
00044     po.addOpt< int >( "dim,d", "Dimension of the mesh", &dim );
00045     //  po.addOpt<void>( "print,p", "Print tree details", &print_tree);
00046     po.parseCommandLine( argc, argv );
00047 
00048     std::vector< int > ints, deps, leafs;
00049     ints.push_back( 10 );
00050     for( int i = 1; i < dints; i++ )
00051         ints.push_back( 2 * ints[i - 1] );
00052     deps.push_back( 30 );
00053     for( int i = 1; i < ddeps; i++ )
00054         deps.push_back( deps[i - 1] - 5 );
00055     leafs.push_back( 6 );
00056     for( int i = 1; i < dleafs; i++ )
00057         leafs.push_back( 2 * leafs[i - 1] );
00058 
00059     ErrorCode rval = MB_SUCCESS;
00060     std::cout << "Tree_type"
00061               << " "
00062               << "Elems_per_leaf"
00063               << " "
00064               << "Tree_depth"
00065               << " "
00066               << "Ints_per_side"
00067               << " "
00068               << "N_elements"
00069               << " "
00070               << "search_time"
00071               << " "
00072               << "perc_outside"
00073               << " "
00074               << "initTime"
00075               << " "
00076               << "nodesVisited"
00077               << " "
00078               << "leavesVisited"
00079               << " "
00080               << "numTraversals"
00081               << " "
00082               << "leafObjectTests" << std::endl;
00083 
00084     // outermost iteration: # elements
00085     for( std::vector< int >::iterator int_it = ints.begin(); int_it != ints.end(); ++int_it )
00086     {
00087         Core mb;
00088         Range elems;
00089         rval = create_hex_mesh( mb, elems, *int_it, dim );
00090         if( MB_SUCCESS != rval ) return rval;
00091 
00092         // iteration: tree depth
00093         for( std::vector< int >::iterator dep_it = deps.begin(); dep_it != deps.end(); ++dep_it )
00094         {
00095 
00096             // iteration: tree max elems/leaf
00097             for( std::vector< int >::iterator leafs_it = leafs.begin(); leafs_it != leafs.end(); ++leafs_it )
00098             {
00099 
00100                 // iteration: tree type
00101                 for( int tree_tp = 0; tree_tp < 2; tree_tp++ )
00102                 {
00103                     // create tree
00104                     Tree* tree;
00105                     if( 0 == tree_tp )
00106                         tree = new BVHTree( &mb );
00107                     else
00108                         tree = new AdaptiveKDTree( &mb );
00109 
00110                     std::ostringstream opts;
00111                     opts << "MAX_DEPTH=" << *dep_it << ";MAX_PER_LEAF=" << *leafs_it;
00112                     if( csints )
00113                     {
00114                         if( opts.str().length() > 0 ) opts << ";";
00115                         opts << "PLANE_SET=" << csints;
00116                     }
00117                     FileOptions fo( opts.str().c_str() );
00118                     rval = tree->parse_options( fo );
00119                     if( MB_SUCCESS != rval ) return rval;
00120                     SpatialLocator sl( &mb, elems, tree );
00121 
00122                     // call evaluation
00123                     double cpu_time, perc_outside;
00124                     rval = test_locator( sl, npoints, cpu_time, perc_outside );
00125                     if( MB_SUCCESS != rval ) return rval;
00126 
00127                     std::cout << ( tree_tp == 0 ? "BVH" : "KD" ) << " " << *leafs_it << " " << *dep_it << " " << *int_it
00128                               << " " << ( *int_it ) * ( *int_it ) * ( *int_it ) << " " << cpu_time << " "
00129                               << perc_outside << " ";
00130 
00131                     tree->tree_stats().output_all_stats();
00132 
00133                 }  // tree_tp
00134 
00135             }  // max elems/leaf
00136 
00137         }  // max depth
00138 
00139     }  // # elements
00140 
00141 #ifdef MOAB_HAVE_MPI
00142     fail = MPI_Finalize();
00143     if( fail ) return fail;
00144 #endif
00145 
00146     return 0;
00147 }
00148 
00149 ErrorCode test_locator( SpatialLocator& sl, int npoints, double& cpu_time, double& percent_outside )
00150 {
00151     BoundBox box     = sl.local_box();
00152     CartVect box_del = box.bMax - box.bMin;
00153 
00154     std::vector< CartVect > test_pts( npoints ), test_res( npoints );
00155     std::vector< EntityHandle > ents( npoints );
00156     int* is_in = new int[npoints];
00157 
00158     double denom = 1.0 / (double)RAND_MAX;
00159     for( int i = 0; i < npoints; i++ )
00160     {
00161         // generate a small number of random point to test
00162         double rx = (double)rand() * denom, ry = (double)rand() * denom, rz = (double)rand() * denom;
00163         test_pts[i] = box.bMin + CartVect( rx * box_del[0], ry * box_del[1], rz * box_del[2] );
00164     }
00165 
00166     CpuTimer ct;
00167 
00168     // call spatial locator to locate points
00169     ErrorCode rval = sl.locate_points( test_pts[0].array(), npoints, &ents[0], test_res[0].array(), &is_in[0] );
00170     if( MB_SUCCESS != rval )
00171     {
00172         delete[] is_in;
00173         return rval;
00174     }
00175 
00176     cpu_time = ct.time_elapsed();
00177 
00178     int num_out     = std::count( is_in, is_in + npoints, false );
00179     percent_outside = ( (double)num_out ) / npoints;
00180     delete[] is_in;
00181 
00182     return rval;
00183 }
00184 
00185 ErrorCode create_hex_mesh( Interface& mb, Range& elems, int n, int dim )
00186 {
00187     ScdInterface* scdi;
00188     ErrorCode rval = mb.query_interface( scdi );
00189     if( MB_SUCCESS != rval ) return rval;
00190 
00191     HomCoord high( n - 1, -1, -1 );
00192     if( dim > 1 ) high[1] = n - 1;
00193     if( dim > 2 ) high[2] = n - 1;
00194     ScdBox* new_box;
00195     rval = scdi->construct_box( HomCoord( 0, 0, 0 ), high, NULL, 0, new_box );
00196     if( MB_SUCCESS != rval ) return rval;
00197     rval = mb.release_interface( scdi );
00198     if( MB_SUCCESS != rval ) return rval;
00199 
00200     rval = mb.get_entities_by_dimension( 0, dim, elems );
00201     if( MB_SUCCESS != rval ) return rval;
00202 
00203     return rval;
00204 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines