Skip to content
Snippets Groups Projects
basic_comm.c 11.5 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/>.
  ----------------------------------------------------------------------*/

#include "reax_types.h"
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*, mpi_out_data* );
typedef void (*coll_unpacker)( void*, void*, mpi_out_data* );


static void int_packer( void *dummy, mpi_out_data *out_buf )
{
    int i;
    int *buf = (int*) dummy;
    int *out = (int*) out_buf->out_atoms;

    for ( i = 0; i < out_buf->cnt; ++i )
    {
        out[i] = buf[ out_buf->index[i] ];
    }
}


static void real_packer( void *dummy, mpi_out_data *out_buf )
{
    int i;
    real *buf = (real*) dummy;
    real *out = (real*) out_buf->out_atoms;

    for ( i = 0; i < out_buf->cnt; ++i )
    {
        out[i] = buf[ out_buf->index[i] ];
    }
}


static void rvec_packer( void *dummy, mpi_out_data *out_buf )
{
    int i;
    rvec *buf, *out;

    buf = (rvec*) dummy;
    out = (rvec*) out_buf->out_atoms;

    for ( i = 0; i < out_buf->cnt; ++i )
    {
        memcpy( out + i, buf + out_buf->index[i], sizeof(rvec) );
    }
}


static void rvec2_packer( void *dummy, mpi_out_data *out_buf )
{
    int i;
    rvec2 *buf, *out;

    buf = (rvec2*) dummy;
    out = (rvec2*) out_buf->out_atoms;

    for ( i = 0; i < out_buf->cnt; ++i )
    {
        memcpy( out + i, buf + out_buf->index[i], sizeof(rvec2) );
    }
}


static void int_unpacker( void *dummy_in, void *dummy_buf, mpi_out_data *out_buf )
{
        int i;
        int *in, *buf;

        in = (int*) dummy_in;
        buf = (int*) dummy_buf;

        for ( i = 0; i < out_buf->cnt; ++i )
        {
            buf[ out_buf->index[i] ] += in[i];
        }
}


static void real_unpacker( void *dummy_in, void *dummy_buf, mpi_out_data *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];
    }
}


static void rvec_unpacker( void *dummy_in, void *dummy_buf, mpi_out_data *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] );

        fprintf( stderr, "rvec_unpacker: cnt=%d  i =%d  index[i]=%d\n",
                out_buf->cnt, i, out_buf->index[i] );
#endif
    }
}


static void rvec2_unpacker( void *dummy_in, void *dummy_buf, mpi_out_data *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];
    }
}


static void * Get_Buffer_Offset( const void * const buffer,
        const int offset, const int type )
{
    void * ptr;

    switch ( type )
    {
        case INT_PTR_TYPE:
            ptr = (int *) buffer + offset;
            break;

        case REAL_PTR_TYPE:
            ptr = (real *) buffer + offset;
            break;

        case RVEC_PTR_TYPE:
            ptr = (rvec *) buffer + offset;
            break;

        case RVEC2_PTR_TYPE:
            ptr = (rvec2 *) buffer + offset;
            break;

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

    return ptr;
}


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

        case REAL_PTR_TYPE:
            ptr = real_packer;
            break;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        case RVEC_PTR_TYPE:
            ptr = rvec_packer;
            break;
        case RVEC2_PTR_TYPE:
            ptr = rvec2_packer;
            break;
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( const int type )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        case REAL_PTR_TYPE:
            ptr = real_unpacker;
            break;

        case RVEC_PTR_TYPE:
            ptr = rvec_unpacker;
            break;

        case RVEC2_PTR_TYPE:
            ptr = rvec2_unpacker;
            break;

        default:
            fprintf( stderr, "[ERROR] unknown pointer type. Terminating...\n" );
            exit( UNKNOWN_OPTION );
            break;
void Dist( const reax_system * const system, mpi_datatypes * const mpi_data,
        void *buf, int buf_type, MPI_Datatype type )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    int d;
    mpi_out_data *out_bufs;
    MPI_Comm comm;
    MPI_Request req1, req2;
    MPI_Status stat1, stat2;
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 )
    {
        /* initiate recvs */
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        if ( nbr1->atoms_cnt )
            MPI_Irecv( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
                    nbr1->atoms_cnt, type, nbr1->rank, 2 * d + 1, comm, &req1 );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        if ( nbr2->atoms_cnt )
            MPI_Irecv( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
                    nbr2->atoms_cnt, type, nbr2->rank, 2 * d, comm, &req2 );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        /* send both messages in dimension d */
        if ( out_bufs[2 * d].cnt )
        {
            pack( buf, out_bufs + (2 * d) );
            MPI_Send( out_bufs[2 * d].out_atoms, out_bufs[2 * d].cnt,
                    type, nbr1->rank, 2 * d, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        }

        if ( out_bufs[2 * d + 1].cnt )
        {
            pack( buf, out_bufs + (2 * d + 1) );
            MPI_Send( out_bufs[2 * d + 1].out_atoms, out_bufs[2 * d + 1].cnt,
                    type, nbr2->rank, 2 * d + 1, comm );
        if( nbr1->atoms_cnt )
        {
            MPI_Wait( &req1, &stat1 );
        }
        if( nbr2->atoms_cnt )
        {
            MPI_Wait( &req2, &stat2 );
        }
void Coll( const reax_system * const system, mpi_datatypes * const mpi_data,
        void *buf, int buf_type, MPI_Datatype type )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    int d;
    mpi_out_data *out_bufs;
    MPI_Comm comm;
    MPI_Request req1, req2;
    MPI_Status stat1, stat2;
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 )
    {
        /* initiate recvs */
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        if ( out_bufs[2 * d].cnt )
            MPI_Irecv( mpi_data->in1_buffer, out_bufs[2 * d].cnt,
                    type, nbr1->rank, 2 * d + 1, comm, &req1 );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        if ( out_bufs[2 * d + 1].cnt )
            MPI_Irecv( mpi_data->in2_buffer, out_bufs[2 * d + 1].cnt,
                    type, nbr2->rank, 2 * d, comm, &req2 );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        /* send both messages in dimension d */
        if ( nbr1->atoms_cnt )
            MPI_Send( Get_Buffer_Offset( buf, nbr1->atoms_str, buf_type ),
                    nbr1->atoms_cnt, type, nbr1->rank, 2 * d, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

        if ( nbr2->atoms_cnt )
            MPI_Send( Get_Buffer_Offset( buf, nbr2->atoms_str, buf_type ),
                    nbr2->atoms_cnt, type, nbr2->rank, 2 * d + 1, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        fprintf( stderr, "p%d coll[%d] nbr1: str=%d cnt=%d recv=%d\n",
                system->my_rank, d, nbr1->atoms_str, nbr1->atoms_cnt,
                out_bufs[2 * d].cnt );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        fprintf( stderr, "p%d coll[%d] nbr2: str=%d cnt=%d recv=%d\n",
                system->my_rank, d, nbr2->atoms_str, nbr2->atoms_cnt,
                out_bufs[2 * d + 1].cnt );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
#endif

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        if ( out_bufs[2 * d].cnt )
        {
            MPI_Wait( &req1, &stat1 );
            unpack( mpi_data->in1_buffer, buf, &out_bufs[2 * d] );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        }

        if ( out_bufs[2 * d + 1].cnt )
        {
            MPI_Wait( &req2, &stat2 );
            unpack( mpi_data->in2_buffer, buf, &out_bufs[2 * d + 1] );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        }
real Parallel_Norm( const real * const v, const int n, MPI_Comm comm )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    int  i;
    real my_sum, norm_sqr;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    /* compute local part of vector 2-norm */
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    my_sum = 0;
    for ( i = 0; i < n; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        my_sum += SQR( v[i] );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    MPI_Allreduce( &my_sum, &norm_sqr, 1, MPI_DOUBLE, MPI_SUM, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    return SQRT( norm_sqr );
real Parallel_Dot( const real * const v1, const real * const v2,
        const int n, MPI_Comm comm )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    int  i;
    real my_dot, res;
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
    my_dot = 0;
    for ( i = 0; i < n; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        my_dot += v1[i] * v2[i];
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    MPI_Allreduce( &my_dot, &res, 1, MPI_DOUBLE, MPI_SUM, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    return res;
real Parallel_Vector_Acc( const real * const v, const int n,
        MPI_Comm comm )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
{
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    int  i;
    real my_acc, res;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

    /* compute local part of vector element-wise sum */
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    my_acc = 0;
    for ( i = 0; i < n; ++i )
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
        my_acc += v[i];
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    MPI_Allreduce( &my_acc, &res, 1, MPI_DOUBLE, MPI_SUM, comm );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    return res;
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
    int i;
    int *id_list;
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed

Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
    MPI_Gather( &system->n, 1, MPI_INT, workspace->rcounts, 1, MPI_INT,
            MASTER_NODE, MPI_COMM_WORLD );
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 = (int*) smalloc( system->n * sizeof(int), "Coll_ids_at_Master::id_list" );
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

    MPI_Gatherv( id_list, system->n, MPI_INT, workspace->id_all,
            workspace->rcounts, workspace->displs, MPI_INT, MASTER_NODE,
            MPI_COMM_WORLD );
    sfree( id_list, "Coll_ids_at_Master::id_list" );
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
{
    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 );
Kurt A. O'Hearn's avatar
Kurt A. O'Hearn committed
}

#endif