Skip to content
Snippets Groups Projects
basic_comm.c 20.4 KiB
Newer Older
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
/*----------------------------------------------------------------------
  PuReMD - Purdue ReaxFF Molecular Dynamics Program
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  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
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  published by the Free Software Foundation; either version 2 of
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  the License, or (at your option) any later version.
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
  See the GNU General Public License for more details:
  <http://www.gnu.org/licenses/>.
  ----------------------------------------------------------------------*/

#if (defined(HAVE_CONFIG_H) && !defined(__CONFIG_H_))
  #define __CONFIG_H_
  #include "../../common/include/config.h"
#endif
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
#if defined(PURE_REAX)
  #include "basic_comm.h"
  #include "vector.h"
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
#elif defined(LAMMPS_REAX)
  #include "reax_basic_comm.h"
  #include "reax_vector.h"
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
#endif


typedef void (*dist_packer)( void const * const, void * const, int const * const, int );
typedef void (*coll_unpacker)( void const * const, void * const, int const * const, int );
static void int_packer( void const * const src, void * const dest,
        int const * const index, int cnt )
static void real_packer( void const * const src, void * const dest,
        int const * const index, int cnt )
static void rvec_packer( void const * const src, void * const dest,
        int const * const index, int cnt )
//        memcpy( &dest_[i], &src_[index[i]], sizeof(rvec) );
        dest_[i][0] = src_[index[i]][0];
        dest_[i][1] = src_[index[i]][1];
        dest_[i][2] = src_[index[i]][2];
static void rvec2_packer( void const * const src, void * const dest,
        int const * const index, int cnt )
    src_ = (rvec2 *) src;
    dest_ = (rvec2 *) dest;
//        memcpy( &dest_[i], &src_[index[i]], sizeof(rvec2) );
        dest_[i][0] = src_[index[i]][0];
        dest_[i][1] = src_[index[i]][1];
static void int_unpacker( void const * const src, void * const dest,
        int const * const index, int cnt )
static void real_unpacker( void const * const src, void * const dest,
        int const * const index, int cnt )
static void rvec_unpacker( void const * const src, void * const dest,
        int const * const index, int cnt )
static void rvec2_unpacker( void const * const src, void * const dest,
        int const * const index, int cnt )
        dest_[index[i]][0] += src_[i][0];
        dest_[index[i]][1] += src_[i][1];
static void * Get_Buffer_Offset( void const * const buffer,
        int offset, int type )
            break;

        default:
            fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
            exit( UNKNOWN_OPTION );
            break;
    }

    return ptr;
}


static dist_packer Get_Packer( int type )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        default:
            fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
            exit( UNKNOWN_OPTION );
            break;
static coll_unpacker Get_Unpacker( int type )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

            break;

        default:
            fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
            exit( UNKNOWN_OPTION );
            break;
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's avatar
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 )
    {
        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__ );
        pack( buf, out_bufs[d].out_atoms, out_bufs[d].index, out_bufs[d].cnt );
        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__ );
        ret = MPI_Waitany( MAX_NT_NBRS, req, &index, stat );
        Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    mpi_out_data *out_bufs;
    MPI_Comm comm;
    MPI_Request req1, req2;
    MPI_Status stat1, stat2;
    MPI_Aint extent, lower_bound;
    size_t type_size;
    ret = MPI_Type_get_extent( type, &lower_bound, &extent );
    Check_MPI_Error( ret, __FILE__, __LINE__ );
//    type_size = MPI_Aint_add( lower_bound, extent );
    type_size = extent;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    comm = mpi_data->comm_mesh3D;
    out_bufs = mpi_data->out_buffers;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    for ( d = 0; d < 3; ++d )
    {
        nbr1 = &system->my_nbrs[2 * d];
        nbr2 = &system->my_nbrs[2 * d + 1];
        srealloc_check( (void **) &out_bufs[2 * d].index,
                TRUE, SAFE_ZONE, __FILE__, __LINE__ );
        pack( buf, out_bufs[2 * d].out_atoms,
                out_bufs[2 * d].index, out_bufs[2 * d].cnt );
        if ( nbr1->rank != system->my_rank )
        {
            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__ );
        }
        smalloc_check( &out_bufs[2 * d + 1].out_atoms,
        srealloc_check( (void **) &out_bufs[2 * d + 1].index,
                TRUE, SAFE_ZONE, __FILE__, __LINE__ );
        pack( buf, out_bufs[2 * d + 1].out_atoms,
                out_bufs[2 * d + 1].index, out_bufs[2 * d + 1].cnt );
        if ( nbr2->rank != system->my_rank )
        {
            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__ );
        }
            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__ );

            if ( cnt1 == MPI_UNDEFINED )
            {
                fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }
            else 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 );
            }

            ret = MPI_Recv( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
                    cnt1, type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
            Check_MPI_Error( ret, __FILE__, __LINE__ );
            pack( buf, Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
                    out_bufs[2 * d + 1].index, out_bufs[2 * d + 1].cnt );
            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__ );

            if ( cnt2 == MPI_UNDEFINED )
            {
                fprintf( stderr, "[ERROR] MPI_Get_count returned MPI_UNDEFINED\n" );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }
            else 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 );
            }

            ret = MPI_Recv( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
                    cnt2, type, nbr2->rank, 2 * d, comm, MPI_STATUS_IGNORE );
            Check_MPI_Error( ret, __FILE__, __LINE__ );
            pack( buf, Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
                    out_bufs[2 * d].index, out_bufs[2 * d].cnt );
            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__ );
}


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)
    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];

        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__ );
    }

    for ( d = 0; d < 6; ++d )
    {
        /* send both messages in direction d */
        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__ );
        ret = MPI_Waitany( MAX_NT_NBRS, req, &index, stat);
        Check_MPI_Error( ret, __FILE__, __LINE__ );
        unpack( in[index], buf, out_bufs[index].index, out_bufs[index].cnt );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    mpi_out_data *out_bufs;
    MPI_Comm comm;
    MPI_Request req1, req2;
    MPI_Status stat1, stat2;
    MPI_Aint extent, lower_bound;
    size_t type_size;
    ret = MPI_Type_get_extent( type, &lower_bound, &extent );
    Check_MPI_Error( ret, __FILE__, __LINE__ );
//    type_size = MPI_Aint_add( lower_bound, extent );
    type_size = extent;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    comm = mpi_data->comm_mesh3D;
    out_bufs = mpi_data->out_buffers;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    for ( d = 2; d >= 0; --d )
    {
        nbr2 = &system->my_nbrs[2 * d + 1];
        
        /* send both messages in dimension d */
        if ( nbr1->rank != system->my_rank )
        {
            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__ );
        }
        if ( nbr2->rank != system->my_rank )
        {
            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__ );
        }
            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__ );

            if ( cnt1 == MPI_UNDEFINED )
            {
                fprintf( stderr, "[ERROR] Coll: MPI_Get_count returned MPI_UNDEFINED\n" );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }
            else if ( cnt1 != out_bufs[2 * d].cnt )
            {
                fprintf( stderr, "[ERROR] Coll: counts mismatch (MPI_Get_count = %d, out_bufs[%d].cnt = %d\n",
                      cnt1, 2 * d, out_bufs[2 * d].cnt );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }

            smalloc_check( &mpi_data->in1_buffer, &mpi_data->in1_buffer_size,
                    type_size * cnt1, TRUE, SAFE_ZONE, __FILE__, __LINE__ );

            ret = MPI_Recv( mpi_data->in1_buffer, cnt1,
                    type, nbr1->rank, 2 * d + 1, comm, MPI_STATUS_IGNORE );
            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__ );

            if ( cnt2 == MPI_UNDEFINED )
            {
                fprintf( stderr, "[ERROR] Coll: MPI_Get_count returned MPI_UNDEFINED\n" );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }
            else if ( cnt2 != out_bufs[2 * d + 1].cnt )
            {
                fprintf( stderr, "[ERROR] Coll: counts mismatch (MPI_Get_count = %d, out_bufs[%d].cnt = %d\n",
                      cnt2, 2 * d + 1, out_bufs[2 * d + 1].cnt );
                MPI_Abort( MPI_COMM_WORLD, RUNTIME_ERROR );
            }

            smalloc_check( &mpi_data->in2_buffer, &mpi_data->in2_buffer_size,
                    type_size * cnt2, TRUE, SAFE_ZONE, __FILE__, __LINE__ );

            ret = MPI_Recv( mpi_data->in2_buffer, 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 );
            Check_MPI_Error( ret, __FILE__, __LINE__ );
            unpack( mpi_data->in1_buffer, buf, out_bufs[2 * d].index,
                    out_bufs[2 * d].cnt );
            unpack( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
                    buf, out_bufs[2 * d + 1].index, out_bufs[2 * d + 1].cnt );
        }
        if ( nbr2->rank != system->my_rank )
        {
            unpack( mpi_data->in2_buffer, buf, out_bufs[2 * d + 1].index,
                    out_bufs[2 * d + 1].cnt );
        }
        else
        {
            unpack( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
                    buf, out_bufs[2 * d].index, out_bufs[2 * d].cnt );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    }
real Parallel_Dot( real const * const v1, real const * const v2,
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    /* compute local part of inner product */
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    for ( i = 0; i < n; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    ret = MPI_Allreduce( &dot_l, &dot, 1, MPI_DOUBLE, MPI_SUM, comm );
    Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
}


/*****************************************************************************/
#if defined(TEST_FORCES)
void Coll_ids_at_Master( reax_system *system, storage *workspace,
        mpi_datatypes *mpi_data )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
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__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    if ( system->my_rank == MASTER_NODE )
    {
        workspace->displs[0] = 0;
        for ( i = 1; i < system->nprocs; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
            workspace->displs[i] = workspace->displs[i - 1] + workspace->rcounts[i - 1];
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    }
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    id_list = smalloc( system->n * sizeof(int), __FILE__, __LINE__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    for ( i = 0; i < system->n; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        id_list[i] = system->my_atoms[i].orig_id;
Kurt A. O'Hearn's avatar
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__ );
    sfree( id_list, __FILE__, __LINE__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    if ( system->my_rank == MASTER_NODE )
    {
        for ( i = 0 ; i < system->bigN; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
            fprintf( stderr, "id_all[%d]: %d\n", i, workspace->id_all[i] );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    }
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
void Coll_rvecs_at_Master( reax_system *system, storage *workspace,
        mpi_datatypes *mpi_data, rvec* v )
Kurt A. O'Hearn's avatar
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,
            MASTER_NODE, MPI_COMM_WORLD );
    Check_MPI_Error( ret, __FILE__, __LINE__ );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
}

#endif