MOAB: Mesh Oriented datABase
(version 5.2.1)
|
00001 /* ***************************************************************** 00002 MESQUITE -- The Mesh Quality Improvement Toolkit 00003 00004 Copyright 2004 Sandia Corporation and Argonne National 00005 Laboratory. Under the terms of Contract DE-AC04-94AL85000 00006 with Sandia Corporation, the U.S. Government retains certain 00007 rights in this software. 00008 00009 This library is free software; you can redistribute it and/or 00010 modify it under the terms of the GNU Lesser General Public 00011 License as published by the Free Software Foundation; either 00012 version 2.1 of the License, or (at your option) any later version. 00013 00014 This library is distributed in the hope that it will be useful, 00015 but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00017 Lesser General Public License for more details. 00018 00019 You should have received a copy of the GNU Lesser General Public License 00020 (lgpl.txt) along with this library; if not, write to the Free Software 00021 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00022 00023 diachin2@llnl.gov, djmelan@sandia.gov, mbrewer@sandia.gov, 00024 pknupp@sandia.gov, tleurent@mcs.anl.gov, tmunson@mcs.anl.gov 00025 00026 ***************************************************************** */ 00027 /*! 00028 \file ParallelHelper.hpp 00029 \brief 00030 00031 Implements ParallelHelper Class 00032 00033 \author Martin Isenburg 00034 \date 2008-03-04 00035 */ 00036 00037 #ifndef Mesquite_ParallelHelper_hpp 00038 #define Mesquite_ParallelHelper_hpp 00039 00040 #include "ParallelHelperInterface.hpp" 00041 00042 #include <vector> 00043 #include <map> 00044 00045 namespace MBMesquite 00046 { 00047 typedef struct VertexIdMapKey 00048 { 00049 size_t glob_id; 00050 int proc_id; 00051 } VertexIdMapKey; 00052 00053 struct VertexIdLessFunc 00054 { 00055 bool operator()( const VertexIdMapKey& that1, const VertexIdMapKey& that2 ) const 00056 { 00057 return ( ( that1.proc_id < that2.proc_id ) || 00058 ( ( that1.proc_id == that2.proc_id ) && ( that1.glob_id < that2.glob_id ) ) ); 00059 } 00060 }; 00061 00062 typedef std::map< VertexIdMapKey, int, VertexIdLessFunc > VertexIdMap; 00063 00064 int get_parallel_rank(); 00065 int get_parallel_size(); 00066 double reduce_parallel_max( double value ); 00067 void parallel_barrier(); 00068 00069 class ParallelHelperImpl : public ParallelHelper 00070 { 00071 public: 00072 enum CommunicationModel 00073 { 00074 TrulyNonBlocking = 0, 00075 TrulyNonBlockingAvoidAllReduce, 00076 NonBlocking, 00077 NonBlockingAvoidAllReduce, 00078 Blocking, 00079 BlockingAvoidAllReduce 00080 }; 00081 00082 ParallelHelperImpl(); 00083 ~ParallelHelperImpl(); 00084 00085 // function called by application during set-up 00086 void set_parallel_mesh( ParallelMesh* mesh ); 00087 void set_communicator( size_t comm ); 00088 void set_communicator( const void* comm ) 00089 { 00090 set_communicator( reinterpret_cast< size_t >( comm ) ); 00091 } 00092 void set_communication_model( int model, MsqError& err ); 00093 void set_generate_random_numbers( int grn, MsqError& err ); 00094 00095 protected: 00096 friend class VertexMover; 00097 // functions called by VertexMover::loop_over_mesh() 00098 void smoothing_init( MsqError& err ); 00099 void compute_first_independent_set( std::vector< Mesh::VertexHandle >& fixed_vertices ); 00100 void communicate_first_independent_set( MsqError& err ); 00101 bool compute_next_independent_set(); 00102 bool get_next_partition_boundary_vertex( MBMesquite::Mesh::VertexHandle& vertex_handle ); 00103 void communicate_next_independent_set( MsqError& err ); 00104 void smoothing_close( MsqError& err ); 00105 00106 protected: 00107 friend class QualityAssessor; 00108 // functions called by QualityAssessor::loop_over_mesh() 00109 int get_rank() const; 00110 int get_nprocs() const; 00111 bool is_our_element( Mesh::ElementHandle element_handle, MsqError& err ) const; 00112 bool is_our_vertex( Mesh::VertexHandle vertex_handle, MsqError& err ) const; 00113 void communicate_min_max_to_all( double* minimum, double* maximum, MsqError& ) const; 00114 void communicate_min_max_to_zero( double* minimum, double* maximum, MsqError& ) const; 00115 void communicate_sums_to_zero( size_t* freeElementCount, int* invertedElementCount, size_t* elementCount, 00116 int* invertedSampleCount, size_t* sampleCount, long unsigned int* count, 00117 long unsigned int* invalid, double* sum, double* sqrSum, MsqError& ) const; 00118 void communicate_power_sum_to_zero( double* pMean, MsqError& ) const; 00119 void communicate_histogram_to_zero( std::vector< int >& histogram, MsqError& ) const; 00120 void communicate_any_true( bool& value, MsqError& err ) const; 00121 void communicate_all_true( bool& value, MsqError& err ) const; 00122 00123 private: 00124 ParallelMesh* mesh; 00125 00126 size_t communicator; 00127 int communication_model; 00128 00129 int rank; 00130 int nprocs; 00131 00132 // variables for VertexMover::loop_over_mesh() 00133 int generate_random_numbers; 00134 std::vector< MBMesquite::Mesh::VertexHandle > vertices; 00135 int num_vertex; 00136 std::vector< char > vtx_in_partition_boundary; 00137 int num_vtx_partition_boundary; 00138 int num_vtx_partition_boundary_local; 00139 int num_vtx_partition_boundary_remote; 00140 std::vector< MBMesquite::Mesh::VertexHandle > part_vertices; 00141 std::vector< int > part_proc_owner; 00142 std::vector< size_t > part_gid; 00143 std::vector< int > part_smoothed_flag; 00144 std::vector< double > part_rand_number; 00145 int num_exportVtx; 00146 std::vector< size_t > exportVtxGIDs; 00147 std::vector< int > exportVtxLIDs; 00148 std::vector< int > exportProc; 00149 std::vector< bool > in_independent_set; 00150 VertexIdMap vid_map; 00151 int total_num_vertices_to_smooth; 00152 int total_num_vertices_to_recv; 00153 std::vector< int > neighbourProcSend; 00154 std::vector< int > neighbourProcRecv; 00155 std::vector< int > neighbourProcSendRemain; 00156 std::vector< int > neighbourProcRecvRemain; 00157 int num_already_smoothed_vertices; 00158 int num_already_recv_vertices; 00159 std::vector< std::vector< int > > vtx_off_proc_list; 00160 std::vector< int > neighbourProc; 00161 int iteration; 00162 int global_work_remains; 00163 int next_vtx_partition_boundary; 00164 /* for exchanging unused ghost node information */ 00165 int unghost_num_vtx; 00166 std::vector< MBMesquite::Mesh::VertexHandle > unghost_vertices; 00167 int unghost_num_procs; 00168 std::vector< int > unghost_procs; 00169 std::vector< int > unghost_procs_num_vtx; 00170 std::vector< int > unghost_procs_offset; 00171 int update_num_vtx; 00172 std::vector< size_t > update_gid; 00173 int update_num_procs; 00174 std::vector< int > update_procs; 00175 std::vector< int > update_procs_num_vtx; 00176 std::vector< int > update_procs_offset; 00177 00178 // functions for VertexMover::loop_over_mesh() 00179 void compute_independent_set(); 00180 int comm_smoothed_vtx_b( MsqError& err ); 00181 int comm_smoothed_vtx_b_no_all( MsqError& err ); 00182 int comm_smoothed_vtx_nb( MsqError& err ); 00183 int comm_smoothed_vtx_nb_no_all( MsqError& err ); 00184 int comm_smoothed_vtx_tnb( MsqError& err ); 00185 int comm_smoothed_vtx_tnb_no_all( MsqError& err ); 00186 }; 00187 00188 } // namespace MBMesquite 00189 #endif // Mesquite_ParallelHelper_hpp