Newer
Older
/*----------------------------------------------------------------------
PuReMD - Purdue ReaxFF Molecular Dynamics Program
Copyright (2010) Purdue University
Hasan Metin Aktulga, haktulga@cs.purdue.edu
Joseph Fogarty, jcfogart@mail.usf.edu
Sagar Pandit, pandit@usf.edu
Ananth Y Grama, ayg@cs.purdue.edu
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
published by the Free Software Foundation; either version 2 of
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
See the GNU General Public License for more details:
<http://www.gnu.org/licenses/>.
----------------------------------------------------------------------*/
Kurt A. O'Hearn
committed
#if (defined(HAVE_CONFIG_H) && !defined(__CONFIG_H_))
#define __CONFIG_H_
#include "../../common/include/config.h"
#endif
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
#include "comm_tools.h"
Kurt A. O'Hearn
committed
#include "tool_box.h"
#include "reax_basic_comm.h"
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
#include "reax_comm_tools.h"
Kurt A. O'Hearn
committed
#include "reax_tool_box.h"
Kurt A. O'Hearn
committed
typedef void (*dist_packer)( void const * const, mpi_out_data * const );
typedef void (*coll_unpacker)( void const * const, void * const,
mpi_out_data * const );
Kurt A. O'Hearn
committed
static void int_packer( void const * const dummy, mpi_out_data * const out_buf )
Kurt A. O'Hearn
committed
int *buf, *out;
buf = (int *) dummy;
out = (int *) out_buf->out_atoms;
for ( i = 0; i < out_buf->cnt; ++i )
{
out[i] = buf[ out_buf->index[i] ];
}
}
Kurt A. O'Hearn
committed
static void real_packer( void const * const dummy, mpi_out_data * const out_buf )
{
int i;
Kurt A. O'Hearn
committed
real *buf, *out;
buf = (real *) dummy;
out = (real *) out_buf->out_atoms;
for ( i = 0; i < out_buf->cnt; ++i )
{
out[i] = buf[ out_buf->index[i] ];
}
}
Kurt A. O'Hearn
committed
static void rvec_packer( void const * const dummy, mpi_out_data * const out_buf )
{
int i;
rvec *buf, *out;
Kurt A. O'Hearn
committed
buf = (rvec *) dummy;
out = (rvec *) out_buf->out_atoms;
for ( i = 0; i < out_buf->cnt; ++i )
{
Kurt A. O'Hearn
committed
memcpy( &out[i], &buf[out_buf->index[i]], sizeof(rvec) );
}
}
Kurt A. O'Hearn
committed
static void rvec2_packer( void const * const dummy, mpi_out_data * const out_buf )
{
int i;
rvec2 *buf, *out;
Kurt A. O'Hearn
committed
buf = (rvec2 *) dummy;
out = (rvec2 *) out_buf->out_atoms;
for ( i = 0; i < out_buf->cnt; ++i )
{
Kurt A. O'Hearn
committed
memcpy( &out[i], &buf[out_buf->index[i]], sizeof(rvec2) );
}
}
Kurt A. O'Hearn
committed
static void int_unpacker( void const * const dummy_in, void * const dummy_buf,
mpi_out_data * const out_buf )
Kurt A. O'Hearn
committed
int i, *in, *buf;
Kurt A. O'Hearn
committed
in = (int*) dummy_in;
buf = (int*) dummy_buf;
Kurt A. O'Hearn
committed
for ( i = 0; i < out_buf->cnt; ++i )
{
//TODO: used in SAI, purpose?
if ( buf[ out_buf->index[i] ] == -1 && in[i] != -1 )
Kurt A. O'Hearn
committed
buf[ out_buf->index[i] ] = in[i];
Kurt A. O'Hearn
committed
}
Kurt A. O'Hearn
committed
static void real_unpacker( void const * const dummy_in, void * const dummy_buf,
mpi_out_data * const out_buf )
{
int i;
real *in, *buf;
in = (real*) dummy_in;
buf = (real*) dummy_buf;
for ( i = 0; i < out_buf->cnt; ++i )
{
buf[ out_buf->index[i] ] += in[i];
}
}
Kurt A. O'Hearn
committed
static void rvec_unpacker( void const * const dummy_in, void * const dummy_buf,
mpi_out_data * const out_buf )
{
int i;
rvec *in, *buf;
in = (rvec*) dummy_in;
buf = (rvec*) dummy_buf;
for ( i = 0; i < out_buf->cnt; ++i )
{
rvec_Add( buf[ out_buf->index[i] ], in[i] );
}
}
Kurt A. O'Hearn
committed
static void rvec2_unpacker( void const * const dummy_in, void * const dummy_buf,
mpi_out_data * const out_buf )
{
int i;
rvec2 *in, *buf;
in = (rvec2*) dummy_in;
buf = (rvec2*) dummy_buf;
for ( i = 0; i < out_buf->cnt; ++i )
{
buf[ out_buf->index[i] ][0] += in[i][0];
buf[ out_buf->index[i] ][1] += in[i][1];
}
}
Kurt A. O'Hearn
committed
static void * Get_Buffer_Offset( void const * const buffer,
int offset, int type )
Kurt A. O'Hearn
committed
{
void * ptr;
switch ( type )
{
Kurt A. O'Hearn
committed
ptr = &((int *) buffer)[offset];
Kurt A. O'Hearn
committed
case REAL_PTR_TYPE:
Kurt A. O'Hearn
committed
ptr = &((real *) buffer)[offset];
Kurt A. O'Hearn
committed
break;
case RVEC_PTR_TYPE:
Kurt A. O'Hearn
committed
ptr = &((rvec *) buffer)[offset];
Kurt A. O'Hearn
committed
break;
case RVEC2_PTR_TYPE:
Kurt A. O'Hearn
committed
ptr = &((rvec2 *) buffer)[offset];
Kurt A. O'Hearn
committed
break;
default:
fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
exit( UNKNOWN_OPTION );
break;
}
return ptr;
}
Kurt A. O'Hearn
committed
static dist_packer Get_Packer( int type )
Kurt A. O'Hearn
committed
dist_packer func_ptr;
switch ( type )
Kurt A. O'Hearn
committed
case INT_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &int_packer;
Kurt A. O'Hearn
committed
break;
case REAL_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &real_packer;
break;
case RVEC_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &rvec_packer;
break;
case RVEC2_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &rvec2_packer;
break;
default:
fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
exit( UNKNOWN_OPTION );
break;
Kurt A. O'Hearn
committed
return func_ptr;
Kurt A. O'Hearn
committed
static coll_unpacker Get_Unpacker( int type )
Kurt A. O'Hearn
committed
coll_unpacker func_ptr;
switch ( type )
Kurt A. O'Hearn
committed
case INT_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &int_unpacker;
Kurt A. O'Hearn
committed
break;
case REAL_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &real_unpacker;
break;
case RVEC_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &rvec_unpacker;
break;
case RVEC2_PTR_TYPE:
Kurt A. O'Hearn
committed
func_ptr = &rvec2_unpacker;
break;
default:
fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
exit( UNKNOWN_OPTION );
break;
Kurt A. O'Hearn
committed
return func_ptr;
Kurt A. O'Hearn
committed
void Dist( reax_system const * const system, mpi_datatypes * const mpi_data,
void const * const buf, int buf_type, MPI_Datatype type )
Kurt A. O'Hearn
committed
#if defined(NEUTRAL_TERRITORY)
Kurt A. O'Hearn
committed
int d, index, ret;
Kurt A. O'Hearn
committed
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req[6];
MPI_Status stat[6];
dist_packer pack;
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_nt_buffers;
pack = Get_Packer( buf_type );
/* initiate recvs */
for ( d = 0; d < 6; ++d )
{
Kurt A. O'Hearn
committed
ret = MPI_Irecv( Get_Buffer_Offset( buf, system->my_nt_nbrs[d].atoms_str, buf_type ),
system->my_nt_nbrs[d].atoms_cnt, type,
system->my_nt_nbrs[d].receive_rank, d, comm, &req[d] );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
Kurt A. O'Hearn
committed
for ( d = 0; d < 6; ++d )
Kurt A. O'Hearn
committed
{
/* send both messages in dimension d */
Kurt A. O'Hearn
committed
pack( buf, &out_bufs[d] );
ret = MPI_Send( out_bufs[d].out_atoms, out_bufs[d].cnt, type,
system->my_nt_nbrs[d].rank, d, comm );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
Kurt A. O'Hearn
committed
for ( d = 0; d < 6; ++d )
Kurt A. O'Hearn
committed
{
Kurt A. O'Hearn
committed
ret = MPI_Waitany( MAX_NT_NBRS, req, &index, stat );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
#else
Kurt A. O'Hearn
committed
int d, cnt1, cnt2, ret;
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req1, req2;
MPI_Status stat1, stat2;
Kurt A. O'Hearn
committed
const neighbor_proc *nbr1, *nbr2;
dist_packer pack;
Kurt A. O'Hearn
committed
MPI_Aint extent, lower_bound;
size_t type_size;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Type_get_extent( type, &lower_bound, &extent );
Check_MPI_Error( ret, __FILE__, __LINE__ );
type_size = MPI_Aint_add( lower_bound, extent );
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_buffers;
pack = Get_Packer( buf_type );
Kurt A. O'Hearn
committed
nbr1 = &system->my_nbrs[2 * d];
nbr2 = &system->my_nbrs[2 * d + 1];
Kurt A. O'Hearn
committed
/* pack MPI buffers and initiate sends */
check_smalloc( &out_bufs[2 * d].out_atoms,
&out_bufs[2 * d].out_atoms_size,
Kurt A. O'Hearn
committed
type_size * out_bufs[2 * d].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist::mpi_data->out_atoms" );
Kurt A. O'Hearn
committed
check_srealloc( (void **) &out_bufs[2 * d].index,
&out_bufs[2 * d].index_size,
Kurt A. O'Hearn
committed
sizeof(int) * out_bufs[2 * d].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist::mpi_data->index" );
Kurt A. O'Hearn
committed
pack( buf, &out_bufs[2 * d] );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Isend( out_bufs[2 * d].out_atoms, out_bufs[2 * d].cnt,
type, nbr1->rank, 2 * d, comm, &req1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
check_smalloc( &out_bufs[2 * d + 1].out_atoms,
&out_bufs[2 * d + 1].out_atoms_size,
Kurt A. O'Hearn
committed
type_size * out_bufs[2 * d + 1].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist::mpi_data->out_atoms" );
Kurt A. O'Hearn
committed
check_srealloc( (void **) &out_bufs[2 * d + 1].index,
&out_bufs[2 * d + 1].index_size,
Kurt A. O'Hearn
committed
sizeof(int) * out_bufs[2 * d + 1].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist::mpi_data->index" );
Kurt A. O'Hearn
committed
pack( buf, &out_bufs[2 * d + 1] );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Isend( out_bufs[2 * d + 1].out_atoms, out_bufs[2 * d + 1].cnt,
type, nbr2->rank, 2 * d + 1, comm, &req2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
/* recv both messages in dimension d */
ret = MPI_Probe( nbr1->rank, 2 * d + 1, comm, &stat1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat1, type, &cnt1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt1 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
#if defined(DEBUG)
if ( cnt1 + nbr1->atoms_str > system->total_cap )
{
fprintf( stderr, "[ERROR] Dist: not enough space in recv buffer for nbr1 (dim = %d)\n", d );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
#endif
Kurt A. O'Hearn
committed
ret = MPI_Recv( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
Kurt A. O'Hearn
committed
cnt1, type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Probe( nbr2->rank, 2 * d, comm, &stat2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat2, type, &cnt2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt2 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
#if defined(DEBUG)
if ( cnt2 + nbr2->atoms_str > system->total_cap )
{
fprintf( stderr, "[ERROR] Dist: not enough space in recv buffer for nbr2 (dim = %d)\n", d );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
#endif
Kurt A. O'Hearn
committed
ret = MPI_Recv( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
Kurt A. O'Hearn
committed
cnt2, type, nbr2->rank, 2 * d, comm, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req1, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req2, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
#endif
Kurt A. O'Hearn
committed
void Dist_FS( reax_system const * const system, mpi_datatypes * const mpi_data,
void const * const buf, int buf_type, MPI_Datatype type )
Kurt A. O'Hearn
committed
int d, cnt1, cnt2, ret;
Kurt A. O'Hearn
committed
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req1, req2;
MPI_Status stat1, stat2;
const neighbor_proc *nbr1, *nbr2;
dist_packer pack;
Kurt A. O'Hearn
committed
MPI_Aint extent, lower_bound;
size_t type_size;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Type_get_extent( type, &lower_bound, &extent );
Check_MPI_Error( ret, __FILE__, __LINE__ );
type_size = MPI_Aint_add( lower_bound, extent );
Kurt A. O'Hearn
committed
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_buffers;
pack = Get_Packer( buf_type );
for ( d = 0; d < 3; ++d )
{
nbr1 = &system->my_nbrs[2 * d];
nbr2 = &system->my_nbrs[2 * d + 1];
Kurt A. O'Hearn
committed
/* pack MPI buffers and initiate sends */
check_smalloc( &out_bufs[2 * d].out_atoms,
&out_bufs[2 * d].out_atoms_size,
Kurt A. O'Hearn
committed
type_size * out_bufs[2 * d].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist_FS::mpi_data->out_atoms" );
Kurt A. O'Hearn
committed
check_srealloc( (void **) &out_bufs[2 * d].index,
&out_bufs[2 * d].index_size,
Kurt A. O'Hearn
committed
sizeof(int) * out_bufs[2 * d].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist_FS::mpi_data->index" );
Kurt A. O'Hearn
committed
pack( buf, &out_bufs[2 * d] );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Isend( out_bufs[2 * d].out_atoms, out_bufs[2 * d].cnt,
type, nbr1->rank, 2 * d, comm, &req1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
check_smalloc( &out_bufs[2 * d + 1].out_atoms,
&out_bufs[2 * d + 1].out_atoms_size,
Kurt A. O'Hearn
committed
type_size * out_bufs[2 * d + 1].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist_FS::mpi_data->out_atoms" );
Kurt A. O'Hearn
committed
check_srealloc( (void **) &out_bufs[2 * d + 1].index,
Kurt A. O'Hearn
committed
&out_bufs[2 * d + 1].index_size,
Kurt A. O'Hearn
committed
sizeof(int) * out_bufs[2 * d + 1].cnt,
Kurt A. O'Hearn
committed
TRUE, SAFE_ZONE, "Dist_FS::mpi_data->index" );
Kurt A. O'Hearn
committed
pack( buf, &out_bufs[2 * d + 1] );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Isend( out_bufs[2 * d + 1].out_atoms, out_bufs[2 * d + 1].cnt,
type, nbr2->rank, 2 * d + 1, comm, &req2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
/* recv both messages in dimension d */
ret = MPI_Probe( nbr1->rank, 2 * d + 1, comm, &stat1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat1, type, &cnt1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt1 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
ret = MPI_Recv( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
Kurt A. O'Hearn
committed
cnt1, type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Probe( nbr2->rank, 2 * d, comm, &stat2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat2, type, &cnt2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt2 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
ret = MPI_Recv( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
Kurt A. O'Hearn
committed
cnt2, type, nbr2->rank, 2 * d, comm, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req1, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req2, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
}
void Coll( reax_system const * const system, mpi_datatypes * const mpi_data,
void * const buf, int buf_type, MPI_Datatype type )
{
#if defined(NEUTRAL_TERRITORY)
Kurt A. O'Hearn
committed
int d, index, ret;
Kurt A. O'Hearn
committed
void *in[6];
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req[6];
MPI_Status stat[6];
coll_unpacker unpack;
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_nt_buffers;
unpack = Get_Unpacker( buf_type );
for ( d = 0; d < 6; ++d )
{
in[d] = mpi_data->in_nt_buffer[d];
Kurt A. O'Hearn
committed
ret = MPI_Irecv( in[d], out_bufs[d].cnt, type,
system->my_nt_nbrs[d].rank, d, comm, &req[d] );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
for ( d = 0; d < 6; ++d )
{
/* send both messages in direction d */
Kurt A. O'Hearn
committed
ret = MPI_Send( Get_Buffer_Offset( buf, system->my_nt_nbrs[d].atoms_str, buf_type ),
system->my_nt_nbrs[d].atoms_cnt, type,
system->my_nt_nbrs[d].receive_rank, d, comm );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
}
Kurt A. O'Hearn
committed
for ( d = 0; d < 6; ++d )
Kurt A. O'Hearn
committed
{
Kurt A. O'Hearn
committed
ret = MPI_Waitany( MAX_NT_NBRS, req, &index, stat);
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
unpack( in[index], buf, &out_bufs[index] );
}
#else
Kurt A. O'Hearn
committed
int d, cnt1, cnt2, ret;
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req1, req2;
MPI_Status stat1, stat2;
Kurt A. O'Hearn
committed
const neighbor_proc *nbr1, *nbr2;
coll_unpacker unpack;
Kurt A. O'Hearn
committed
MPI_Aint extent, lower_bound;
size_t type_size;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Type_get_extent( type, &lower_bound, &extent );
Check_MPI_Error( ret, __FILE__, __LINE__ );
type_size = MPI_Aint_add( lower_bound, extent );
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_buffers;
unpack = Get_Unpacker( buf_type );
Kurt A. O'Hearn
committed
nbr1 = &system->my_nbrs[2 * d];
Kurt A. O'Hearn
committed
nbr2 = &system->my_nbrs[2 * d + 1];
/* send both messages in dimension d */
ret = MPI_Isend( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
nbr1->atoms_cnt, type, nbr1->rank, 2 * d, comm, &req1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Isend( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
nbr2->atoms_cnt, type, nbr2->rank, 2 * d + 1, comm, &req2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
/* recvs and unpack messages */
ret = MPI_Probe( nbr1->rank, 2 * d + 1, comm, &stat1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat1, type, &cnt1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt1 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
check_smalloc( &mpi_data->in1_buffer, &mpi_data->in1_buffer_size,
Kurt A. O'Hearn
committed
type_size * cnt1, TRUE, SAFE_ZONE, "Coll::mpi_data->in1_buffer" );
Kurt A. O'Hearn
committed
ret = MPI_Recv( mpi_data->in1_buffer, cnt1,
Kurt A. O'Hearn
committed
type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Probe( nbr2->rank, 2 * d, comm, &stat2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat2, type, &cnt2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
if ( cnt2 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
check_smalloc( &mpi_data->in2_buffer, &mpi_data->in2_buffer_size,
Kurt A. O'Hearn
committed
type_size * cnt2, TRUE, SAFE_ZONE, "Coll::mpi_data->in2_buffer" );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Recv( mpi_data->in2_buffer, cnt2,
Kurt A. O'Hearn
committed
type, nbr2->rank, 2 * d, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
ret = MPI_Wait( &req1, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req2, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
unpack( mpi_data->in1_buffer, buf, &out_bufs[2 * d] );
unpack( mpi_data->in2_buffer, buf, &out_bufs[2 * d + 1] );
Kurt A. O'Hearn
committed
}
#endif
}
void Coll_FS( reax_system const * const system, mpi_datatypes * const mpi_data,
void * const buf, int buf_type, MPI_Datatype type )
{
Kurt A. O'Hearn
committed
int d, cnt1, cnt2, ret;
Kurt A. O'Hearn
committed
mpi_out_data *out_bufs;
MPI_Comm comm;
MPI_Request req1, req2;
MPI_Status stat1, stat2;
const neighbor_proc *nbr1, *nbr2;
coll_unpacker unpack;
Kurt A. O'Hearn
committed
MPI_Aint extent, lower_bound;
size_t type_size;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Type_get_extent( type, &lower_bound, &extent );
Check_MPI_Error( ret, __FILE__, __LINE__ );
type_size = MPI_Aint_add( lower_bound, extent );
Kurt A. O'Hearn
committed
comm = mpi_data->comm_mesh3D;
out_bufs = mpi_data->out_buffers;
unpack = Get_Unpacker( buf_type );
for ( d = 2; d >= 0; --d )
{
nbr1 = &system->my_nbrs[2 * d];
Kurt A. O'Hearn
committed
nbr2 = &system->my_nbrs[2 * d + 1];
/* send both messages in dimension d */
ret = MPI_Isend( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
nbr1->atoms_cnt, type, nbr1->rank, 2 * d, comm, &req1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Isend( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
nbr2->atoms_cnt, type, nbr2->rank, 2 * d + 1, comm, &req2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
/* recvs and unpack messages */
ret = MPI_Probe( nbr1->rank, 2 * d + 1, comm, &stat1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat1, type, &cnt1 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
if ( cnt1 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
check_smalloc( &mpi_data->in1_buffer, &mpi_data->in1_buffer_size,
Kurt A. O'Hearn
committed
type_size * cnt1, TRUE, SAFE_ZONE, "Coll_FS::mpi_data->in1_buffer" );
Kurt A. O'Hearn
committed
ret = MPI_Recv( mpi_data->in1_buffer, cnt1,
Kurt A. O'Hearn
committed
type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
ret = MPI_Probe( nbr2->rank, 2 * d, comm, &stat2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Get_count( &stat2, type, &cnt2 );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
if ( cnt2 == MPI_UNDEFINED )
{
fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
}
Kurt A. O'Hearn
committed
check_smalloc( &mpi_data->in2_buffer, &mpi_data->in2_buffer_size,
Kurt A. O'Hearn
committed
type_size * cnt2, TRUE, SAFE_ZONE, "Coll_FS::mpi_data->in2_buffer" );
Kurt A. O'Hearn
committed
ret = MPI_Recv( mpi_data->in2_buffer, cnt2,
Kurt A. O'Hearn
committed
type, nbr2->rank, 2 * d, comm, MPI_STATUS_IGNORE );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
ret = MPI_Wait( &req1, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
ret = MPI_Wait( &req2, MPI_STATUS_IGNORE );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
unpack( mpi_data->in1_buffer, buf, &out_bufs[2 * d] );
unpack( mpi_data->in2_buffer, buf, &out_bufs[2 * d + 1] );
Kurt A. O'Hearn
committed
real Parallel_Norm( real const * const v, const int n, MPI_Comm comm )
Kurt A. O'Hearn
committed
int i, ret;
Kurt A. O'Hearn
committed
real sum_l, norm_sqr;
Kurt A. O'Hearn
committed
sum_l = 0.0;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
/* compute local part of vector 2-norm */
Kurt A. O'Hearn
committed
sum_l += SQR( v[i] );
Kurt A. O'Hearn
committed
ret = MPI_Allreduce( &sum_l, &norm_sqr, 1, MPI_DOUBLE, MPI_SUM, comm );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
real Parallel_Dot( real const * const v1, real const * const v2,
const int n, MPI_Comm comm )
Kurt A. O'Hearn
committed
int i, ret;
Kurt A. O'Hearn
committed
real dot_l, dot;
Kurt A. O'Hearn
committed
dot_l = 0.0;
Kurt A. O'Hearn
committed
Kurt A. O'Hearn
committed
/* compute local part of inner product */
Kurt A. O'Hearn
committed
dot_l += v1[i] * v2[i];
Kurt A. O'Hearn
committed
ret = MPI_Allreduce( &dot_l, &dot, 1, MPI_DOUBLE, MPI_SUM, comm );
Kurt A. O'Hearn
committed
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
return dot;
Kurt A. O'Hearn
committed
real Parallel_Vector_Acc( real const * const v, const int n,
MPI_Comm comm )
Kurt A. O'Hearn
committed
int i, ret;
Kurt A. O'Hearn
committed
/* compute local part of vector element-wise sum */
Kurt A. O'Hearn
committed
my_acc = 0.0;
Kurt A. O'Hearn
committed
ret = MPI_Allreduce( &my_acc, &res, 1, MPI_DOUBLE, MPI_SUM, comm );
Check_MPI_Error( ret, __FILE__, __LINE__ );
}
/*****************************************************************************/
#if defined(TEST_FORCES)
Kurt A. O'Hearn
committed
void Coll_ids_at_Master( reax_system *system, storage *workspace,
mpi_datatypes *mpi_data )
Kurt A. O'Hearn
committed
int i, *id_list, ret;
Kurt A. O'Hearn
committed
ret = MPI_Gather( &system->n, 1, MPI_INT, workspace->rcounts, 1, MPI_INT,
MASTER_NODE, MPI_COMM_WORLD );
Check_MPI_Error( ret, __FILE__, __LINE__ );
if ( system->my_rank == MASTER_NODE )
{
workspace->displs[0] = 0;
for ( i = 1; i < system->nprocs; ++i )
workspace->displs[i] = workspace->displs[i - 1] + workspace->rcounts[i - 1];
Kurt A. O'Hearn
committed
id_list = (int*) smalloc( system->n * sizeof(int), "Coll_ids_at_Master::id_list" );
Kurt A. O'Hearn
committed
ret = MPI_Gatherv( id_list, system->n, MPI_INT, workspace->id_all,
workspace->rcounts, workspace->displs, MPI_INT,
MASTER_NODE, MPI_COMM_WORLD );
Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn
committed
sfree( id_list, "Coll_ids_at_Master::id_list" );
Kurt A. O'Hearn
committed
#if defined(DEBUG_FOCUS)
if ( system->my_rank == MASTER_NODE )
{
for ( i = 0 ; i < system->bigN; ++i )
fprintf( stderr, "id_all[%d]: %d\n", i, workspace->id_all[i] );
void Coll_rvecs_at_Master( reax_system *system, storage *workspace,
Kurt A. O'Hearn
committed
int ret;
ret = MPI_Gatherv( v, system->n, mpi_data->mpi_rvec, workspace->f_all,
workspace->rcounts, workspace->displs, mpi_data->mpi_rvec,
Kurt A. O'Hearn
committed
MASTER_NODE, MPI_COMM_WORLD );
Check_MPI_Error( ret, __FILE__, __LINE__ );