14 #include <stk_util/parallel/ParallelComm.hpp> 15 #include <stk_util/parallel/ParallelReduce.hpp> 21 #if defined( STK_HAS_MPI ) 23 enum { STK_MPI_TAG_SIZING = 0 , STK_MPI_TAG_DATA = 1 };
30 const CommBuffer *
const send ,
31 const CommBuffer *
const recv ,
34 typedef unsigned char * ucharp ;
36 static const char method[] =
"stk_classic::CommAll::communicate" ;
43 std::vector<int> tmp( p_size * 4 );
45 int *
const send_counts = (tmp.empty() ? NULL : & tmp[0]) ;
46 int *
const send_displs = send_counts + p_size ;
47 int *
const recv_counts = send_displs + p_size ;
48 int *
const recv_displs = recv_counts + p_size ;
50 unsigned char *
const ps =
static_cast<ucharp
>(send[0].buffer());
51 unsigned char *
const pr =
static_cast<ucharp
>(recv[0].buffer());
53 for (
unsigned i = 0 ; i < p_size ; ++i ) {
54 const CommBuffer & send_buf = send[i] ;
55 const CommBuffer & recv_buf = recv[i] ;
57 send_counts[i] = send_buf.capacity();
58 recv_counts[i] = recv_buf.capacity();
60 send_displs[i] =
static_cast<ucharp
>(send_buf.buffer()) - ps ;
61 recv_displs[i] =
static_cast<ucharp
>(recv_buf.buffer()) - pr ;
64 result = MPI_Alltoallv( ps , send_counts , send_displs , MPI_BYTE ,
65 pr , recv_counts , recv_displs , MPI_BYTE ,
68 if ( MPI_SUCCESS != result ) {
69 msg << method <<
" GLOBAL ERROR: " << result <<
" == MPI_Alltoallv" ;
73 return MPI_SUCCESS == result ;
77 const CommBuffer *
const send ,
78 const CommBuffer *
const recv ,
81 static const char method[] =
"stk_classic::CommAll::communicate" ;
82 static const int mpi_tag = STK_MPI_TAG_DATA ;
84 int result = MPI_SUCCESS ;
93 unsigned num_recv = 0 ;
95 for (
unsigned i = 0 ; i < p_size ; ++i ) {
96 if ( recv[i].capacity() ) { ++num_recv ; }
102 MPI_Request request_null = MPI_REQUEST_NULL ;
103 std::vector<MPI_Request> request( num_recv , request_null );
104 std::vector<MPI_Status> status( num_recv );
108 for (
unsigned i = 0 ; result == MPI_SUCCESS && i < p_size ; ++i ) {
109 const unsigned recv_size = recv[i].capacity();
110 void *
const recv_buf = recv[i].buffer();
112 result = MPI_Irecv( recv_buf , recv_size , MPI_BYTE ,
113 i , mpi_tag , p_comm , & request[count] );
118 if ( MPI_SUCCESS != result ) {
119 msg << method <<
" LOCAL[" << p_rank <<
"] ERROR: " 120 << result <<
" == MPI_Irecv , " ;
126 int local_error = MPI_SUCCESS == result ? 0 : 1 ;
127 int global_error = 0 ;
129 result = MPI_Allreduce( & local_error , & global_error ,
130 1 , MPI_INT , MPI_SUM , p_comm );
132 if ( MPI_SUCCESS != result ) {
133 msg << method <<
" GLOBAL ERROR: " << result <<
" == MPI_Allreduce" ;
135 else if ( global_error ) {
136 result = MPI_ERR_UNKNOWN ;
145 for (
unsigned i = 0 ; MPI_SUCCESS == result && i < p_size ; ++i ) {
146 const int dst = ( i + p_rank ) % p_size ;
147 const unsigned send_size = send[dst].capacity();
148 void *
const send_buf = send[dst].buffer();
150 result = MPI_Rsend( send_buf , send_size , MPI_BYTE ,
151 dst , mpi_tag , p_comm );
155 if ( MPI_SUCCESS != result ) {
156 msg << method <<
" LOCAL ERROR: " << result <<
" == MPI_Rsend , " ;
159 MPI_Request *
const p_request = (request.empty() ? NULL : & request[0]) ;
160 MPI_Status *
const p_status = (status.empty() ? NULL : & status[0]) ;
162 result = MPI_Waitall( num_recv , p_request , p_status );
165 if ( MPI_SUCCESS != result ) {
166 msg << method <<
" LOCAL[" << p_rank <<
"] ERROR: " 167 << result <<
" == MPI_Waitall , " ;
171 for (
unsigned i = 0 ; i < num_recv ; ++i ) {
172 MPI_Status *
const recv_status = & status[i] ;
173 const int recv_proc = recv_status->MPI_SOURCE ;
174 const int recv_tag = recv_status->MPI_TAG ;
175 const int recv_plan = recv[recv_proc].capacity();
178 MPI_Get_count( recv_status , MPI_BYTE , & recv_count );
180 if ( recv_tag != mpi_tag || recv_count != recv_plan ) {
181 msg << method <<
" LOCAL[" << p_rank <<
"] ERROR: Recv[" 182 << recv_proc <<
"] Size( " 183 << recv_count <<
" != " << recv_plan <<
" ) , " ;
184 result = MPI_ERR_UNKNOWN ;
191 return MPI_SUCCESS == result ;
203 const CommBuffer *
const send ,
204 const CommBuffer *
const recv ,
206 {
return send == recv ; }
209 const CommBuffer *
const send ,
210 const CommBuffer *
const recv ,
212 {
return send == recv ; }
223 size_t align_quad(
size_t n )
225 enum { Size = 4 *
sizeof(int) };
226 return n + CommBufferAlign<Size>::align(n);
233 void CommBuffer::pack_overflow()
const 235 std::ostringstream os ;
236 os <<
"stk_classic::CommBuffer::pack<T>(...){ overflow by " ;
239 throw std::overflow_error( os.str() );
242 void CommBuffer::unpack_overflow()
const 244 std::ostringstream os ;
245 os <<
"stk_classic::CommBuffer::unpack<T>(...){ overflow by " ;
248 throw std::overflow_error( os.str() );
251 void CommAll::rank_error(
const char * method ,
unsigned p )
const 253 std::ostringstream os ;
254 os <<
"stk_classic::CommAll::" << method
255 <<
"(" << p <<
") ERROR: Not in [0:" << m_size <<
")" ;
256 throw std::range_error( os.str() );
261 CommBuffer::CommBuffer()
262 : m_beg(NULL), m_ptr(NULL), m_end(NULL)
265 CommBuffer::~CommBuffer()
268 void CommBuffer::deallocate(
const unsigned number , CommBuffer * buffers )
270 if ( NULL != buffers ) {
271 for (
unsigned i = 0 ; i < number ; ++i ) {
272 ( buffers + i )->~CommBuffer();
278 CommBuffer * CommBuffer::allocate(
279 const unsigned number ,
const unsigned *
const size )
281 const size_t n_base = align_quad( number *
sizeof(CommBuffer) );
282 size_t n_size = n_base ;
284 if ( NULL != size ) {
285 for (
unsigned i = 0 ; i < number ; ++i ) {
286 n_size += align_quad( size[i] );
292 void *
const p_malloc = malloc( n_size );
294 CommBuffer *
const b_base =
295 p_malloc != NULL ?
reinterpret_cast<CommBuffer*
>(p_malloc)
296 : reinterpret_cast<CommBuffer*>( NULL );
298 if ( p_malloc != NULL ) {
300 for (
unsigned i = 0 ; i < number ; ++i ) {
301 new( b_base + i ) CommBuffer();
304 if ( NULL != size ) {
306 ucharp ptr =
reinterpret_cast<ucharp
>( p_malloc );
310 for (
unsigned i = 0 ; i < number ; ++i ) {
311 CommBuffer & b = b_base[i] ;
314 b.m_end = ptr + size[i] ;
315 ptr += align_quad( size[i] );
329 CommBuffer::deallocate( m_size , m_send );
330 if ( 1 < m_size ) { CommBuffer::deallocate( m_size , m_recv ); }
341 m_size( 0 ), m_rank( 0 ),
357 m_send = CommBuffer::allocate( m_size , NULL );
359 if ( NULL == m_send ) {
360 std::string msg(
"stk_classic::CommAll::CommAll FAILED malloc");
361 throw std::runtime_error(msg);
365 bool CommAll::allocate_buffers(
const unsigned num_msg_bounds ,
366 const bool symmetric ,
367 const bool local_flag )
369 const unsigned zero = 0 ;
370 std::vector<unsigned> tmp( m_size , zero );
372 for (
unsigned i = 0 ; i < m_size ; ++i ) {
373 tmp[i] = m_send[i].size();
376 const unsigned *
const send_size = (tmp.empty() ? NULL : & tmp[0]) ;
377 const unsigned *
const recv_size = symmetric ? (tmp.empty() ? NULL : & tmp[0]) : NULL ;
379 return allocate_buffers( m_comm, num_msg_bounds,
380 send_size, recv_size, local_flag );
385 void CommAll::reset_buffers()
388 CommBuffer * m = m_send ;
389 CommBuffer *
const me = m + m_size ;
390 for ( ; m != me ; ++m ) { m->reset(); }
392 if ( m_recv && 1 < m_size ) {
393 CommBuffer * m = m_recv ;
394 CommBuffer *
const me = m + m_size ;
395 for ( ; m != me ; ++m ) { m->reset(); }
401 void CommAll::swap_send_recv()
403 if ( m_recv == NULL ) {
406 msg(
"stk_classic::CommAll::swap_send_recv(){ NULL recv buffers }" );
407 throw std::logic_error( msg );
410 CommBuffer * tmp_msg = m_send ;
418 const unsigned num_msg_bounds ,
419 const unsigned *
const send_size ,
420 const unsigned *
const recv_size ,
421 const bool local_flag )
423 static const char method[] =
"stk_classic::CommAll::allocate_buffers" ;
424 const unsigned uzero = 0 ;
426 CommBuffer::deallocate( m_size , m_send );
427 CommBuffer::deallocate( m_size , m_recv );
432 m_bound = num_msg_bounds ;
434 std::ostringstream msg ;
440 const bool send_none = NULL == send_size ;
442 std::vector<unsigned> tmp_send ;
444 if ( send_none ) { tmp_send.resize( m_size , uzero ); }
446 const unsigned *
const send = send_none ? (tmp_send.empty() ? NULL : & tmp_send[0]) : send_size ;
448 m_send = CommBuffer::allocate( m_size , send );
452 std::vector<unsigned> tmp_recv ;
454 const bool recv_tbd = NULL == recv_size ;
458 tmp_recv.resize( m_size , uzero );
460 unsigned *
const r = (tmp_recv.empty() ? NULL : & tmp_recv[0]) ;
462 comm_sizes( m_comm , m_bound , m_max , send , r );
465 const unsigned *
const recv = recv_tbd ? (tmp_recv.empty() ? NULL : & tmp_recv[0]) : recv_size ;
467 m_recv = CommBuffer::allocate( m_size , recv );
474 bool error_alloc = m_send == NULL || m_recv == NULL ;
485 enum { Length = 2 + 2 * NPSum };
487 int local_result[ Length ];
488 int global_result[ Length ];
490 Copy<Length>( local_result , 0 );
492 local_result[ Length - 2 ] = error_alloc ;
493 local_result[ Length - 1 ] = local_flag ;
495 if ( ! error_alloc ) {
497 const unsigned r = 2 * ( m_rank % NPSum );
499 for (
unsigned i = 0 ; i < m_size ; ++i ) {
500 const unsigned n_send = m_send[i].capacity();
501 const unsigned n_recv = m_recv[i].capacity();
503 const unsigned s = 2 * ( i % NPSum );
505 local_result[s] += n_send ? 1 : 0 ;
506 local_result[s+1] += n_send ;
508 local_result[r] -= n_recv ? 1 : 0 ;
509 local_result[r+1] -= n_recv ;
517 Copy<Length>(global_result, local_result);
522 error_alloc = global_result[ Length - 2 ] ;
523 global_flag = global_result[ Length - 1 ] ;
527 for (
unsigned i = 0 ; ok && i < 2 * NPSum ; ++i ) {
528 ok = 0 == global_result[i] ;
531 if ( error_alloc || ! ok ) {
532 msg << method <<
" ERROR:" ;
533 if ( error_alloc ) { msg <<
" Failed memory allocation ," ; }
534 if ( ! ok ) { msg <<
" Parallel inconsistent send/receive ," ; }
535 throw std::runtime_error( msg.str() );
543 void CommAll::communicate()
545 static const char method[] =
"stk_classic::CommAll::communicate" ;
547 std::ostringstream msg ;
551 for (
unsigned i = 0 ; i < m_size ; ++i ) {
553 if ( m_send[i].remaining() ) {
554 msg << method <<
" LOCAL[" << m_rank <<
"] ERROR: Send[" << i
555 <<
"] Buffer not filled." ;
556 throw std::underflow_error( msg.str() );
567 if ( m_bound < m_max ) {
568 ok = all_to_all_dense( m_comm , m_send , m_recv , msg );
571 ok = all_to_all_sparse( m_comm , m_send , m_recv , msg );
574 if ( ! ok ) {
throw std::runtime_error( msg.str() ); }
581 CommBroadcast::CommBroadcast(
ParallelMachine comm ,
unsigned root_rank )
585 m_root_rank( root_rank ),
589 bool CommBroadcast::allocate_buffer(
const bool local_flag )
591 static const char method[] =
"stk_classic::CommBroadcast::allocate_buffer" ;
593 unsigned root_rank_min = m_root_rank ;
594 unsigned root_rank_max = m_root_rank ;
595 unsigned root_send_size = m_root_rank == m_rank ? m_buffer.size() : 0 ;
596 unsigned flag = local_flag ;
598 all_reduce( m_comm , ReduceMin<1>( & root_rank_min ) &
599 ReduceMax<1>( & root_rank_max ) &
600 ReduceMax<1>( & root_send_size ) &
601 ReduceBitOr<1>( & flag ) );
603 if ( root_rank_min != root_rank_max ) {
605 msg.append( method );
606 msg.append(
" FAILED: inconsistent root processor" );
607 throw std::runtime_error( msg );
610 m_buffer.m_beg =
static_cast<CommBuffer::ucharp
>( malloc( root_send_size ) );
611 m_buffer.m_ptr = m_buffer.m_beg ;
612 m_buffer.m_end = m_buffer.m_beg + root_send_size ;
617 CommBroadcast::~CommBroadcast()
620 if ( m_buffer.m_beg ) { free( static_cast<void*>( m_buffer.m_beg ) ); }
622 m_buffer.m_beg = NULL ;
623 m_buffer.m_ptr = NULL ;
624 m_buffer.m_end = NULL ;
627 CommBuffer & CommBroadcast::recv_buffer()
632 CommBuffer & CommBroadcast::send_buffer()
634 static const char method[] =
"stk_classic::CommBroadcast::send_buffer" ;
636 if ( m_root_rank != m_rank ) {
638 msg.append( method );
639 msg.append(
" FAILED: is not root processor" );
640 throw std::runtime_error( msg );
646 void CommBroadcast::communicate()
648 #if defined( STK_HAS_MPI ) 650 const int count = m_buffer.capacity();
651 void *
const buf = m_buffer.buffer();
653 const int result = MPI_Bcast( buf, count, MPI_BYTE, m_root_rank, m_comm);
655 if ( MPI_SUCCESS != result ) {
656 std::ostringstream msg ;
657 msg <<
"stk_classic::CommBroadcast::communicate ERROR : " 658 << result <<
" == MPI_Bcast" ;
659 throw std::runtime_error( msg.str() );
670 CommGather::~CommGather()
673 free( static_cast<void*>( m_send.m_beg ) );
675 if ( NULL != m_recv_count ) { free( static_cast<void*>( m_recv_count ) ); }
677 if ( NULL != m_recv ) { CommBuffer::deallocate( m_size , m_recv ); }
685 if ( NULL != m_recv ) {
686 for (
unsigned i = 0 ; i < m_size ; ++i ) { m_recv[i].reset(); }
690 CommBuffer & CommGather::recv_buffer(
unsigned p )
692 static CommBuffer empty ;
694 return m_size <= p ? empty : (
695 m_size <= 1 ? m_send : m_recv[p] );
701 unsigned root_rank ,
unsigned send_size )
705 m_root_rank( root_rank ),
711 m_send.m_beg =
static_cast<CommBuffer::ucharp
>( malloc( send_size ) );
712 m_send.m_ptr = m_send.m_beg ;
713 m_send.m_end = m_send.m_beg + send_size ;
715 #if defined( STK_HAS_MPI ) 719 const bool is_root = m_rank == m_root_rank ;
722 m_recv_count =
static_cast<int*
>( malloc(2*m_size*
sizeof(
int)) );
723 m_recv_displ = m_recv_count + m_size ;
726 MPI_Gather( & send_size , 1 , MPI_INT ,
727 m_recv_count , 1 , MPI_INT ,
728 m_root_rank , m_comm );
731 m_recv = CommBuffer::allocate( m_size ,
732 reinterpret_cast<unsigned*>( m_recv_count ) );
734 for (
unsigned i = 0 ; i < m_size ; ++i ) {
735 m_recv_displ[i] = m_recv[i].m_beg - m_recv[0].m_beg ;
745 void CommGather::communicate()
747 #if defined( STK_HAS_MPI ) 751 const int send_count = m_send.capacity();
753 void *
const send_buf = m_send.buffer();
754 void *
const recv_buf = m_rank == m_root_rank ? m_recv->buffer() : NULL ;
756 MPI_Gatherv( send_buf , send_count , MPI_BYTE ,
757 recv_buf , m_recv_count , m_recv_displ , MPI_BYTE ,
758 m_root_rank , m_comm );
769 #if defined( STK_HAS_MPI ) 772 const unsigned *
const send_size ,
773 unsigned *
const recv_size ,
776 static const char method[] =
"stk_classic::comm_dense_sizes" ;
778 const unsigned zero = 0 ;
781 std::vector<unsigned> send_buf( p_size * 2 , zero );
782 std::vector<unsigned> recv_buf( p_size * 2 , zero );
784 for (
unsigned i = 0 ; i < p_size ; ++i ) {
785 const unsigned i2 = i * 2 ;
786 send_buf[i2] = send_size[i] ;
787 send_buf[i2+1] = local_flag ;
791 unsigned *
const ps = (send_buf.empty() ? NULL : & send_buf[0]) ;
792 unsigned *
const pr = (recv_buf.empty() ? NULL : & recv_buf[0]) ;
794 MPI_Alltoall( ps , 2 , MPI_UNSIGNED , pr , 2 , MPI_UNSIGNED , comm );
796 if ( MPI_SUCCESS != result ) {
798 msg.append( method );
799 msg.append(
" FAILED: MPI_SUCCESS != MPI_Alltoall" );
800 throw std::runtime_error( msg );
804 bool global_flag = false ;
806 for (
unsigned i = 0 ; i < p_size ; ++i ) {
807 const unsigned i2 = i * 2 ;
808 recv_size[i] = recv_buf[i2] ;
809 if ( recv_buf[i2+1] ) { global_flag = true ; }
821 void sum_np_max_2_op(
824 const int np = *len - 2 ;
825 unsigned * ind = (
unsigned *) inv ;
826 unsigned * outd = (
unsigned *) outv ;
831 for (
int i = 0 ; i < np ; ++i ) {
836 if ( outd[0] < ind[0] ) { outd[0] = ind[0] ; }
837 if ( outd[1] < ind[1] ) { outd[1] = ind[1] ; }
845 const unsigned num_msg_bound ,
846 unsigned & num_msg_maximum ,
847 const unsigned *
const send_size ,
848 unsigned *
const recv_size ,
851 static const char method[] =
"stk_classic::comm_unknown_sizes" ;
852 const unsigned uzero = 0 ;
854 static MPI_Op mpi_op = MPI_OP_NULL ;
856 if ( mpi_op == MPI_OP_NULL ) {
858 MPI_Op_create( sum_np_max_2_op , 1 , & mpi_op );
866 std::ostringstream msg ;
868 num_msg_maximum = 0 ;
870 unsigned num_recv = 0 ;
871 unsigned max_msg = 0 ;
872 bool global_flag = false ;
875 std::vector<unsigned> send_buf( p_size + 2 , uzero );
876 std::vector<unsigned> recv_buf( p_size + 2 , uzero );
878 unsigned *
const p_send = (send_buf.empty() ? NULL : & send_buf[0]) ;
879 unsigned *
const p_recv = (recv_buf.empty() ? NULL : & recv_buf[0]) ;
881 for (
unsigned i = 0 ; i < p_size ; ++i ) {
883 if ( send_size[i] ) {
888 send_buf[p_size] = max_msg ;
889 send_buf[p_size+1] = local_flag ;
891 result = MPI_Allreduce(p_send,p_recv,p_size+2,MPI_UNSIGNED,mpi_op,comm);
893 if ( result != MPI_SUCCESS ) {
895 msg << method <<
" ERROR: " << result <<
" == MPI_AllReduce" ;
896 throw std::runtime_error( msg.str() );
899 num_recv = recv_buf[ p_rank ] ;
900 max_msg = recv_buf[ p_size ] ;
901 global_flag = recv_buf[ p_size + 1 ] ;
907 for (
unsigned i = 0 ; i < p_size ; ++i ) {
908 if ( max_msg < recv_buf[i] ) { max_msg = recv_buf[i] ; }
912 num_msg_maximum = max_msg ;
914 if ( num_msg_bound < max_msg ) {
918 MPI_Alltoall( (
void*) send_size , 1 , MPI_UNSIGNED ,
919 recv_size , 1 , MPI_UNSIGNED , comm );
921 if ( MPI_SUCCESS != result ) {
923 msg << method <<
" ERROR: " << result <<
" == MPI_Alltoall" ;
924 throw std::runtime_error( msg.str() );
927 else if ( max_msg ) {
930 const int mpi_tag = STK_MPI_TAG_SIZING ;
932 MPI_Request request_null = MPI_REQUEST_NULL ;
933 std::vector<MPI_Request> request( num_recv , request_null );
934 std::vector<MPI_Status> status( num_recv );
935 std::vector<unsigned> buf( num_recv );
939 for (
unsigned i = 0 ; i < num_recv ; ++i ) {
940 unsigned *
const p_buf = & buf[i] ;
941 MPI_Request *
const p_request = & request[i] ;
942 result = MPI_Irecv( p_buf , 1 , MPI_UNSIGNED ,
943 MPI_ANY_SOURCE , mpi_tag , comm , p_request );
944 if ( MPI_SUCCESS != result ) {
946 msg << method <<
" ERROR: " << result <<
" == MPI_Irecv" ;
947 throw std::runtime_error( msg.str() );
954 for (
unsigned i = 0 ; i < p_size ; ++i ) {
955 int dst = ( i + p_rank ) % p_size ;
956 unsigned value = send_size[dst] ;
958 result = MPI_Send( & value , 1 , MPI_UNSIGNED , dst , mpi_tag , comm );
959 if ( MPI_SUCCESS != result ) {
961 msg << method <<
" ERROR: " << result <<
" == MPI_Send" ;
962 throw std::runtime_error( msg.str() );
970 MPI_Request *
const p_request = (request.empty() ? NULL : & request[0]) ;
971 MPI_Status *
const p_status = (status.empty() ? NULL : & status[0]) ;
972 result = MPI_Waitall( num_recv , p_request , p_status );
974 if ( MPI_SUCCESS != result ) {
976 msg << method <<
" ERROR: " << result <<
" == MPI_Waitall" ;
977 throw std::runtime_error( msg.str() );
982 for (
unsigned i = 0 ; i < num_recv ; ++i ) {
983 MPI_Status *
const recv_status = & status[i] ;
984 const int recv_proc = recv_status->MPI_SOURCE ;
985 const int recv_tag = recv_status->MPI_TAG ;
988 MPI_Get_count( recv_status , MPI_UNSIGNED , & recv_count );
990 if ( recv_tag != mpi_tag || recv_count != 1 ) {
991 msg << method <<
" ERROR: Received buffer mismatch " ;
992 msg <<
"P" << p_rank <<
" <- P" << recv_proc ;
993 msg <<
" " << 1 <<
" != " << recv_count ;
994 throw std::runtime_error( msg.str() );
997 const unsigned r_size = buf[i] ;
998 recv_size[ recv_proc ] = r_size ;
1002 return global_flag ;
1013 unsigned & num_msg_maximum ,
1014 const unsigned *
const send_size ,
1015 unsigned *
const recv_size ,
1018 num_msg_maximum = send_size[0] ? 1 : 0 ;
1020 recv_size[0] = send_size[0] ;
1026 const unsigned *
const send_size ,
1027 unsigned *
const recv_size ,
1030 recv_size[0] = send_size[0] ;
bool comm_dense_sizes(ParallelMachine comm, const unsigned *const send_size, unsigned *const recv_size, bool local_flag)
void all_reduce_sum(ParallelMachine comm, const double *local, double *global, unsigned count)
Parallel summation to all processors.
unsigned parallel_machine_rank(ParallelMachine parallel_machine)
Member function parallel_machine_rank ...
ParallelMachine parallel_machine_null()
parallel_machine_null returns MPI_COMM_NULL if MPI is enabled.
unsigned parallel_machine_size(ParallelMachine parallel_machine)
Member function parallel_machine_size ...
bool comm_sizes(ParallelMachine comm, const unsigned num_msg_bound, unsigned &num_msg_maximum, const unsigned *const send_size, unsigned *const recv_size, bool local_flag)
MPI_Datatype ParallelDatatype
void reset(MPI_Comm new_comm)
Function reset determines new parallel_size and parallel_rank. Flushes, closes, and reopens log files...
void all_reduce(ParallelMachine, const ReduceOp &)
eastl::iterator_traits< InputIterator >::difference_type count(InputIterator first, InputIterator last, const T &value)