Mesh Oriented datABase  (version 5.4.1)
Array-based unstructured mesh datastructure
point_locater.hpp
Go to the documentation of this file.
00001 /**
00002  * point_locater.hpp
00003  * Ryan H. Lewis
00004  * Copyright 2012
00005  */
00006 #include <vector>
00007 #include <iostream>
00008 #include <ctime>
00009 
00010 #include <sstream>
00011 #include <fstream>
00012 
00013 #ifndef POINT_LOCATER_HPP
00014 #define POINT_LOCATER_HPP
00015 
00016 namespace moab
00017 {
00018 
00019 template < typename _Tree, typename _Boxes >
00020 class Point_search
00021 {
00022 
00023     // public types
00024   public:
00025     typedef _Tree Tree;
00026     typedef _Boxes Boxes;
00027     //  typedef typename Tree::Elements::value_type Element;
00028     typedef int Error;
00029 
00030     // private types
00031   private:
00032     typedef Point_search< _Tree, _Boxes > Self;
00033     // public methods
00034   public:
00035     // Constructor
00036     Point_search( Tree& _tree, Boxes& _boxes ) : tree_( _tree ), boxes( _boxes ) {}
00037 
00038     // Copy constructor
00039     Point_search( Self& s ) : tree_( s.tree_ ), boxes( s.boxes ) {}
00040 
00041     // private functionality
00042   private:
00043     // TODO: deprecate this
00044     template < typename Point_map, typename List >
00045     void resolve_boxes( const Point_map& query_points, List& list )
00046     {
00047         /*
00048      typedef typename std::vector< bool> Bitmask;
00049      typedef typename Points::const_iterator Point;
00050      typedef typename Tree::const_element_iterator Element;
00051      typedef typename std::vector< std::size_t> Processor_list;
00052      typedef typename List::value_type Tuple;
00053      const Element & end = tree_.element_end();
00054      Bitmask located_points( query_points.size(), 0);
00055      std::size_t index=0;
00056      for( Point i = query_points.begin();
00057             i != query_points.end(); ++i,++index){
00058          const Element  & element = tree_.find( *i);
00059          if(element != end){
00060              located_points[ index] = 1;
00061          }
00062      }
00063      for(int i = 0; i < located_points.size(); ++i){
00064          if(!located_points[ i]){
00065              Processor_list processors;
00066              const Point & point = query_point.begin()+i;
00067              resolve_box_for_point( point, processors);
00068              for( std::size_t p = processors.begin();
00069                       p != processors.end(); ++p){
00070                  list.push_back( Tuple( *point, *p) );
00071              }
00072          }
00073      }
00074      */
00075     }
00076 
00077     // public functionality
00078   public:
00079     template < typename Point_map, typename Entities, typename Communicator >
00080     Error locate_points( Point_map& query_points, Entities& entities, Communicator& comm, double tol )
00081     {
00082         /*TODO: implement a parallel location algorithm here
00083         original algorithm:
00084         all-scatter/gather bounding boxes
00085         for each point perform box checks
00086         scatter/gather points
00087         perform location and solution interpolation
00088         //send results back
00089         ------------
00090         new algorithm
00091         Assume power num_proc = n^2;
00092         all-reduce bounding boxes
00093         divide box up into n x n grid
00094         each source processor can communicate its elements to grid proc
00095         each target processor then communicates its points to grid proc
00096         grid proc sends points to correct source proc.
00097         */
00098         return 0;
00099     }
00100 
00101     template < typename Points, typename Entities >
00102     Error locate_points( const Points& query_points, Entities& entities, double tol ) const
00103     {
00104         typedef typename Points::const_iterator Point_iterator;
00105         typedef typename Entities::value_type Result;
00106         Entities result;
00107         result.reserve( query_points.size() );
00108         for( Point_iterator i = query_points.begin(); i != query_points.end(); ++i )
00109         {
00110             Result h;
00111             tree_.find( *i, tol, h );
00112             result.push_back( h );
00113         }
00114         entities = result;
00115         return 0;
00116     }
00117 
00118     template < typename Points, typename Entities >
00119     Error bruteforce_locate_points( const Points& query_points, Entities& entities, double tol ) const
00120     {
00121         // TODO: this could be faster with caching, but of course this is
00122         // really just for testing
00123         typedef typename Points::const_iterator Point_iterator;
00124         typedef typename Entities::value_type::first_type Entity_handle;
00125         Entities result;
00126         result.reserve( query_points.size() );
00127         std::size_t count = 0;
00128         std::stringstream ss;
00129         typename Entities::iterator j = entities.begin();
00130         for( Point_iterator i = query_points.begin(); i != query_points.end(); ++i, ++j )
00131         {
00132             if( j->first == 0 )
00133             {
00134                 const Entity_handle h = tree_.bruteforce_find( *i, tol );
00135                 if( h == 0 )
00136                 {
00137                     ++count;
00138                     for( int k = 0; k < 3; ++k )
00139                     {
00140                         ss << ( *i )[k];
00141                         if( k < 2 )
00142                         {
00143                             ss << ", ";
00144                         }
00145                         else
00146                         {
00147                             ss << std::endl;
00148                         }
00149                     }
00150                 }
00151             }
00152         }
00153         std::ofstream out( "unlocated_pts" );
00154         out << ss.str();
00155         out.close();
00156         std::cout << count << " vertices are not contained in _any_ elements!" << std::endl;
00157         return 0;
00158     }
00159 
00160     // public accessor methods
00161   public:
00162     Tree& tree() const
00163     {
00164         return tree_;
00165     }
00166 
00167     // private data members
00168   private:
00169     const Tree& tree_;
00170     const Boxes& boxes;
00171 };  // class Point_search
00172 
00173 }  // namespace moab
00174 
00175 #endif  // POINT_LOCATER_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines