gsMpiComm.h

Go to the documentation of this file.

#pragma once

#include <gsCore/gsForwardDeclarations.h>

namespace gismo
{

#ifndef GISMO_WITH_MPI
typedef int MPI_Comm;
typedef int MPI_Group;
typedef int MPI_Request;
#endif

class gsSerialGroup
{
public:
    static int compare (const gsSerialGroup&)
    {
        return 0;
    }

    static int compare (const int&)
    {
        return 0;
    }

    static gsSerialGroup diff (const gsSerialGroup&)
    {
        return gsSerialGroup();
    }

    static gsSerialGroup diff (const int&)
    {
        return gsSerialGroup();
    }

    static gsSerialGroup intersect (const gsSerialGroup&)
    {
        return gsSerialGroup();
    }

    static gsSerialGroup intersect (const int&)
    {
        return gsSerialGroup();
    }

    static gsSerialGroup unite (const gsSerialGroup&)
    {
        return gsSerialGroup();
    }

    static gsSerialGroup unite (const int&)
    {
        return gsSerialGroup();
    }

    static int rank ()
    {
        return 0;
    }

    static int size ()
    {
        return 1;
    }

    static std::ostream &print(std::ostream &os)
    {
        os << "gsSerialGroup : rank = " << rank() << ", size = " << size();
        return os;
    }

    const MPI_Group* operator& () const
    {
        static MPI_Group m_group(0);
        return &m_group;
    }

};

inline std::ostream& operator<<(std::ostream& os, const gsSerialGroup& obj)
{
    obj.print(os);
    return os;
}

class GISMO_EXPORT gsSerialStatus
{
public:
    struct MPI_Status {};
    static int rank () { return 0; }

    static int tag () { return 0; }

    template<typename T>
        static int size ()
    {
        return 1;
    }

    const MPI_Status* operator& () const
    {
        static MPI_Status status;
        return &status;
    }

    static std::ostream &print(std::ostream &os)
    {
        os << "gsSerialStatus : rank = " << gsSerialStatus::rank()
           << ", tag = " << gsSerialStatus::tag();
        return os;
    }
};

inline std::ostream& operator<<(std::ostream& os, const gsSerialStatus& obj)
{
    obj.print(os);
    return os;
}

class gsSerialRequest
{
public:
    static int cancel ()
    {
        return 0;
    }

    static gsSerialStatus status ()
    {
        return gsSerialStatus();
    }

    static gsSerialStatus test ()
    {
        return gsSerialStatus();
    }

    static gsSerialStatus wait ()
    {
        return gsSerialStatus();
    }

    static gsSerialStatus waitAny (int, gsSerialRequest [], int*)
    {
        return gsSerialStatus();
    }

    static gsSerialRequest getNullRequest()
    {
        gsSerialRequest request;
        return request;
    }

    const MPI_Request* operator& () const
    {
        static MPI_Request req(0);
        return &req;
    }

    static std::ostream &print(std::ostream &os)
    {
        os << "gsSerialRequest";
        return os;
    }
};

inline std::ostream& operator<<(std::ostream& os, const gsSerialRequest& obj)
{
    obj.print(os);
    return os;
}

class gsSerialComm
{
public:
    int rank () const { return 0; }

    static int size () { return 1; }

    static std::string name() { return "gsSerialComm"; }

    static int compare ( gsSerialComm)
    {
        return 0;
    }

    gsSerialComm duplicate () const
    {
        return gsSerialComm(*this);
    }

    static int group ( const MPI_Group*)
    {
        return 0;
    }

    gsSerialComm split (int, int) const
    {
        return gsSerialComm(*this);
    }

#ifdef GISMO_WITH_MPI
    operator MPI_Comm () const  { return MPI_COMM_SELF;}
#else
    operator MPI_Comm () const  { return 0;}
    //    typedef int MPI_Group;
    //    typedef int MPI_Request;
    //    struct MPI_Status {};
#endif

public:

    template<typename T>
    static T sum (T& in)
    {
        return in;
    }

    template<typename T>
    static int sum (T* inout, int len)
    {
        return 0;
    }

    template<typename T>
    static T prod (T& in)
    {
        return in;
    }

    template<typename T>
    static int prod (T* inout, int len)
    {
        return 0;
    }

    template<typename T>
    static T min (T& in)
    {
        return in;
    }

    template<typename T>
    static int min (T* inout, int len)
    {
        return 0;
    }

    template<typename T>
    static T max (T& in)
    {
        return in;
    }

    template<typename T>
    static int max (T* inout, int len)
    {
        return 0;
    }

    static int barrier ()
    {
        return 0;
    }

    static gsSerialStatus probe(int, int = 0)
    {
        return gsSerialStatus();
    }

    static gsSerialStatus iprobe (int, int*, int = 0)
    {
        return gsSerialStatus();
    }

    template<typename T>
    static int send (T*, int, int, int = 0)
    {
        return 0;
    }

    template<typename T>
    static int isend (T*, int, int, MPI_Request, int = 0)
    {
        return 0;
    }

    template<typename T>
    static int recv (T*, int, int, int = 0, const void* = NULL)
    {
        return 0;
    }

    template<typename T>
    static int irecv (T*, int, int, MPI_Request, int = 0)
    {
        return 0;
    }

    template<typename T>
    static int broadcast (T* , int , int )
    {
        return 0;
    }

    template<typename T>
    static int gather (T* in, T* out, int len, int) // note out must have same size as in
    {
        // copy_n(in, len, out);
        for (int i=0; i<len; i++)
            out[i] = in[i];
        return 0;
    }

    template<typename T>
    static int gatherv (T* in, int sendlen, T* out, int* recvlen, int* displ, int root)
    {
        for (int i=*displ; i<sendlen; i++)
            out[i] = in[i];
        return 0;
    }

    template<typename T>
    static int scatter (T* send, T* recv, int len, int root) // note out must have same size as in
    {
        for (int i=0; i<len; i++)
            recv[i] = send[i];
        return 0;
    }

    template<typename T>
    static int scatterv (T* send, int* sendlen, int* displ, T* recv, int recvlen, int root)
    {
        for (int i=*displ; i<*sendlen; i++)
            recv[i] = send[i];
        return 0;
    }

    template<typename T>
    static int allgather(T* sbuf, int count, T* rbuf)
    {
        for(T* end=sbuf+count; sbuf < end; ++sbuf, ++rbuf)
            *rbuf=*sbuf;
        return 0;
    }

    template<typename T>
    static int allgatherv (T* in, int sendlen, T* out, int* recvlen, int* displ)
    {
        for (int i=*displ; i<sendlen; i++)
            out[i] = in[i];
        return 0;
    }

    template<typename BinaryFunction, typename Type>
    static int allreduce(Type* inout, int len)
    {
        return 0;
    }

    template<typename BinaryFunction, typename Type>
    static void allreduce(Type* in, Type* out, int len)
    {
        std::copy(in, in+len, out);
        return;
    }
};

#ifdef GISMO_WITH_MPI

class gsMpiGroup
{
public:
    gsMpiGroup()
    {}

    gsMpiGroup(MPI_Group& group) : m_group(group)
    {}

    ~gsMpiGroup()
    {
        MPI_Group_free(&m_group);
    }

    int compare (const gsMpiGroup& other) const
    {
        int result;
        MPI_Group_compare(m_group, *(&other), &result);
        return result;
    }

    int compare (const MPI_Group& other) const
    {
        int result;
        MPI_Group_compare(m_group, other, &result);
        return result;
    }

    gsMpiGroup diff (const gsMpiGroup& other) const
    {
        gsMpiGroup diff_;
        MPI_Group_difference(m_group, *(&other), &diff_);
        return diff_;
    }

    gsMpiGroup diff (const MPI_Group& other) const
    {
        gsMpiGroup diff_;
        MPI_Group_difference(m_group, other, &diff_);
        return diff_;
    }

    gsMpiGroup intersect (const gsMpiGroup& other) const
    {
        gsMpiGroup intersect_;
        MPI_Group_intersection(m_group, *(&other), &intersect_);
        return intersect_;
    }

    gsMpiGroup intersect (const MPI_Group& other) const
    {
        gsMpiGroup intersect_;
        MPI_Group_intersection(m_group, other, &intersect_);
        return intersect_;
    }

    gsMpiGroup unite (const gsMpiGroup& other) const
    {
        gsMpiGroup union_;
        MPI_Group_union(m_group, *(&other), &union_);
        return union_;
    }

    gsMpiGroup unite (const MPI_Group& other) const
    {
        gsMpiGroup union_;
        MPI_Group_union(m_group, other, &union_);
        return union_;
    }

    int rank () const
    {
        int rank_;
        MPI_Group_rank(m_group, &rank_);
        return rank_;
    }

    int size () const
    {
        int size_;
        MPI_Group_size(m_group, &size_);
        return size_;
    }

    std::ostream &print(std::ostream &os) const
    {
        os << "gsMpiGroup : rank = " << rank() << ", size = " << size();
        return os;
    }

    MPI_Group* operator& ()
    {
        return &m_group;
    }

    const MPI_Group* operator& () const
    {
        return &m_group;
    }

private:
    MPI_Group m_group;
};

inline std::ostream& operator<<(std::ostream& os, const gsMpiGroup& obj)
{
    obj.print(os);
    return os;
}

class gsMpiStatus : public MPI_Status
{
public:
    int rank () const { return MPI_SOURCE; }

    int tag () const { return MPI_TAG; }

    template<typename T>
    int size () const
    {
        int count;
        MPI_Get_count(this, MPITraits<T>::getType(), &count);
        return count;
    }

    std::ostream &print(std::ostream &os) const
    {
        os << "gsMpiStatus : rank = " << rank() << ", tag = " << tag();
        return os;
    }
};

inline std::ostream& operator<<(std::ostream& os, const gsMpiStatus& obj)
{
    obj.print(os);
    return os;
}

class gsMpiRequest
{
public:
    int cancel ()
    {
        return MPI_Cancel(&m_request);
    }

    int free ()
    {
        return MPI_Request_free(&m_request);
    }

    gsMpiStatus status () const
    {
        gsMpiStatus status;
        int flag;
        MPI_Request_get_status(m_request, &flag, &status);
        return status;
    }

    gsMpiStatus test ()
    {
        gsMpiStatus status;
        int flag;
        MPI_Test(&m_request, &flag, &status);
        return status;
    }

    gsMpiStatus wait ()
    {
        gsMpiStatus status;
        MPI_Wait(&m_request, &status);
        return status;
    }

    std::ostream &print(std::ostream &os) const
    {
        os << "gsMpiRequest";
        return os;
    }

    MPI_Request* operator& ()
    {
        return &m_request;
    }

    const MPI_Request* operator& () const
    {
        return &m_request;
    }

    static gsMpiStatus waitAny (int numberRequests, gsMpiRequest requests[], int* outIndex)
    {
        gsMpiStatus status;
        MPI_Request mpiRequests[numberRequests];
        for(int i = 0; i < numberRequests; i++)
        {
            mpiRequests[i] = requests[i].m_request;
        }

        MPI_Waitany(numberRequests, mpiRequests, outIndex, &status);
        return status;
    }

    static gsMpiRequest getNullRequest()
    {
        gsMpiRequest request;
        request.m_request = MPI_REQUEST_NULL;
        return request;
    }

private:
    MPI_Request m_request;
};

inline std::ostream& operator<<(std::ostream& os, const gsMpiRequest& obj)
{
    obj.print(os);
    return os;
}

class GISMO_EXPORT gsMpiComm
{
    friend class gsMpi;

public:

    gsMpiComm() : rank_(-1), size_(0) { }

    gsMpiComm(const MPI_Comm & _comm)
    : m_comm(_comm)
    {
        if(_comm != MPI_COMM_NULL)
        {
#           ifndef NDEBUG
            int initialized = 0;
            MPI_Initialized(&initialized);
            GISMO_ENSURE(1==initialized,
                         "You must call gsMpi::init(..) in your main() function"
                         " before using gsMpiComm");
            MPI_Comm_set_errhandler(m_comm, ErrHandler);
#           endif
            MPI_Comm_rank(m_comm, &rank_);
            MPI_Comm_size(m_comm, &size_);
        }
        else
        {
            size_ = 0;
            rank_ =-1;
        }
    }

    gsMpiComm(const gsSerialComm &) : m_comm(MPI_COMM_SELF) { }

    typedef MPI_Comm Communicator;

    int rank () const { return rank_; }

    int size () const { return size_; }

    std::string name() const
    {
        char str[MPI_MAX_OBJECT_NAME];
        int len;
        MPI_Comm_get_name(m_comm, str, &len);
        return std::string(str, len);
    }

    int compare ( MPI_Comm comm ) const
    {
        int result;
        MPI_Comm_compare(m_comm, comm, &result);
        return result;
    }

    MPI_Comm duplicate () const
    {
        MPI_Comm comm;
        MPI_Comm_dup(m_comm, &comm);
        return comm;
    }

    int group (MPI_Group* group_) const
    {
        return MPI_Comm_group(m_comm, group_);
    }

    MPI_Comm split (int color, int key) const
    {
        MPI_Comm comm;
        MPI_Comm_split(m_comm, color, key, &comm);
        return comm;
    }

    operator MPI_Comm () const { return m_comm; }

private:
    int rank_;
    int size_;

    MPI_Comm m_comm;

#   ifndef NDEBUG
protected:

    // Mpi error handling
    static void ErrCallBack(MPI_Comm *comm, int *err_code, ...)
    {
        char err_string[MPI_MAX_ERROR_STRING];
        int err_length, err_class;

        int rank;
        MPI_Comm_get_name(*comm, err_string, &err_length);
        MPI_Comm_rank(*comm, &rank);
        gsWarn << "MPI error ("<<*err_code<<") at process "<< rank
               <<" of "<< err_string <<"\n";
        MPI_Error_class(*err_code, &err_class);
        MPI_Error_string(err_class, err_string, &err_length);
        gsWarn <<"gsMpi error class: "<<err_class <<" ("<< err_string <<")\n";
        MPI_Error_string(*err_code, err_string, &err_length);
        gsWarn <<"gsMpi error      : "<<*err_code <<" ("<< err_string <<")\n";
        throw std::runtime_error("GISMO_ERROR: " + std::string(err_string, err_length));
        //MPI_Abort(*comm, *err_code);
    }

    static MPI_Errhandler ErrHandler;
#   endif

public:

    template<typename T>
    T sum (T& in) const
    {
        T out;
        allreduce<std::plus<T> >(&in,&out,1);
        return out;
    }

    template<typename T>
    int sum (T* inout, int len) const
    {
        return allreduce<std::plus<T> >(inout,len);
    }

    template<typename T>
    int sum (T* inout, int len, int root) const
    {
        return reduce<std::plus<T> >(inout,len,root);
    }

    template<typename T>
    int sum (T* in, T* out, int len, int root) const
    {
        return reduce<std::plus<T> >(in,out,len,root);
    }

    template<typename T>
    int isum (T* in,T* out, int len, int root, MPI_Request* req) const
    {
        return iallreduce<std::plus<T> >(in, out,len,req);
    }

    template<typename T>
    int isum (T* inout, int len, int root, MPI_Request* req) const
    {
        return ireduce<std::plus<T> >(inout,len,root,req);
    }

    template<typename T>
    int isum (T* inout, int len, MPI_Request* req) const
    {
        return iallreduce<std::plus<T> >(inout,len,req);
    }


    template<typename T>
    T prod (T& in) const
    {
        T out;
        allreduce<std::multiplies<T> >(&in,&out,1);
        return out;
    }

    template<typename T>
    int prod (T* inout, int len) const
    {
        return allreduce<std::multiplies<T> >(inout,len);
    }

    template<typename T>
    T min (T& in) const
    {
        T out;
        allreduce<Min<T> >(&in,&out,1);
        return out;
    }

    template<typename T>
    int min (T* inout, int len) const
    {
        return allreduce<Min<T> >(inout,len);
    }


    template<typename T>
    T max (T& in) const
    {
        T out;
        allreduce<Max<T> >(&in,&out,1);
        return out;
    }

    template<typename T>
    int max (T* inout, int len) const
    {
        return allreduce<Max<T> >(inout,len);
    }

    int barrier () const
    {
        return MPI_Barrier(m_comm);
    }

    gsMpiStatus probe (int source, int tag = 0) const
    {
        gsMpiStatus status;
        MPI_Probe(source,tag,m_comm,&status);
        return status;
    }

    gsMpiStatus iprobe (int source, int* flag, int tag = 0) const
    {
        gsMpiStatus status;
        MPI_Iprobe(source,tag,m_comm,flag,&status);
        return status;
    }

    template<typename T>
    int send (T* in, int len, int dest, int tag = 0) const
    {
        return MPI_Send(in,len,MPITraits<T>::getType(),
                        dest,tag,m_comm);
    }

    template<typename T>
    int isend (T* in, int len, int dest, MPI_Request* req, int tag = 0) const
    {
        return MPI_Isend(in,len,MPITraits<T>::getType(),
                         dest,tag,m_comm,req);
    }

    template<typename T>
    int recv (T* out, int len, int source, int tag = 0, MPI_Status* status = NULL) const
    {
        return MPI_Recv(out,len,MPITraits<T>::getType(),
                        source,tag,m_comm,(status == NULL ? MPI_STATUS_IGNORE : status));
    }

    template<typename T>
    int irecv (T* out, int len, int source, MPI_Request* req, int tag = 0) const
    {
        return MPI_Irecv(out,len,MPITraits<T>::getType(),
                         source,tag,m_comm,req);
    }

    template<typename T>
    int broadcast (T* inout, int len, int root) const
    {
        return MPI_Bcast(inout,len,MPITraits<T>::getType(),root,m_comm);
    }

    template<typename T>
    int gather (T* in, T* out, int len, int root) const
    {
        return MPI_Gather(in,len,MPITraits<T>::getType(),
                          out,len,MPITraits<T>::getType(),
                          root,m_comm);
    }

    template<typename T>
    int gatherv (T* in, int sendlen, T* out, int* recvlen, int* displ, int root) const
    {
        return MPI_Gatherv(in,sendlen,MPITraits<T>::getType(),
                           out,recvlen,displ,MPITraits<T>::getType(),
                           root,m_comm);
    }

    template<typename T>
    int scatter (T* send, T* recv, int len, int root) const
    {
        return MPI_Scatter(send,len,MPITraits<T>::getType(),
                           recv,len,MPITraits<T>::getType(),
                           root,m_comm);
    }

    template<typename T>
    int scatterv (T* send, int* sendlen, int* displ, T* recv, int recvlen, int root) const
    {
        return MPI_Scatterv(send,sendlen,displ,MPITraits<T>::getType(),
                            recv,recvlen,MPITraits<T>::getType(),
                            root,m_comm);
    }

    template<typename T>
    int alltoall (T* send, T* recv, int sendcount, int recvcount) const
    {
        return MPI_Alltoall(send,sendcount,MPITraits<T>::getType(),
                            recv,recvcount,MPITraits<T>::getType(),
                            m_comm);
    }

    template<typename T>
    int alltoallv (T* send, int* sendcount, int* senddispl, T* recv, int* recvcount, int* recvdispl) const
    {
        return MPI_Alltoallv(send,sendcount,senddispl,MPITraits<T>::getType(),
                             recv,recvcount,recvdispl,MPITraits<T>::getType(),
                             m_comm);
    }

    template<typename T, typename T1>
    int allgather(T* sbuf, int count, T1* rbuf) const
    {
        return MPI_Allgather(sbuf, count, MPITraits<T>::getType(),
                             rbuf, count, MPITraits<T1>::getType(),
                             m_comm);
    }

    template<typename T>
    int allgatherv (T* in, int sendlen, T* out, int* recvlen, int* displ) const
    {
        return MPI_Allgatherv(in,sendlen,MPITraits<T>::getType(),
                              out,recvlen,displ,MPITraits<T>::getType(),
                              m_comm);
    }

#ifndef MPI_IN_PLACE
#define MPI_IN_PLACE inout
#define MASK_MPI_IN_PLACE
/*
  #  ifdef _MSC_VER
  #    pragma message ("Masking MPI_IN_PLACE (not found in MPI version used).")
  #  else
  #    warning Masking MPI_IN_PLACE (not found in MPI version used).
  #  endif
*/
#endif

    template<typename BinaryFunction, typename Type>
    int allreduce(Type* inout, int len) const
    {
        // Type* out = new Type[len];
        // int ret = allreduce<BinaryFunction>(inout,out,len);
        // std::copy(out, out+len, inout);
        // delete[] out;
        // return ret;
        return MPI_Allreduce(MPI_IN_PLACE, inout, len, MPITraits<Type>::getType(),
                             (Generic_MPI_Op<Type, BinaryFunction>::get()),m_comm);
    }


    template<typename BinaryFunction, typename Type>
    int allreduce(Type* in, Type* out, int len) const
    {
        return MPI_Allreduce(in, out, len, MPITraits<Type>::getType(),
                             (Generic_MPI_Op<Type, BinaryFunction>::get()),m_comm);
    }



    template<typename BinaryFunction, typename Type>
    int iallreduce(Type* in, Type* out, int len, MPI_Request* req) const
    {
        return MPI_Iallreduce(in, out, len, MPITraits<Type>::getType(),
                              (Generic_MPI_Op<Type, BinaryFunction>::get()),m_comm,req);
    }

    template<typename BinaryFunction, typename Type>
    int iallreduce(Type* inout, int len, MPI_Request* req) const
    {
        return MPI_Iallreduce(MPI_IN_PLACE, inout, len, MPITraits<Type>::getType(),
                              (Generic_MPI_Op<Type, BinaryFunction>::get()),m_comm,req);
    }

    template<typename BinaryFunction, typename Type>
    int reduce(Type* inout, int len,int root) const
    {
        int ret;
        if(root == rank())
            ret = MPI_Reduce(MPI_IN_PLACE, inout, len, MPITraits<Type>::getType(),
                             (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm);
        else
            ret = MPI_Reduce(inout, NULL, len, MPITraits<Type>::getType(),
                             (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm);
        return ret;
    }

    template<typename BinaryFunction, typename Type>
    int reduce(Type* in,Type* out, int len,int root) const
    {
        return MPI_Reduce(in, out, len, MPITraits<Type>::getType(),
                          (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm);
    }

    template<typename BinaryFunction, typename Type>
    int ireduce(Type* inout, int len,int root, MPI_Request* req) const
    {
        int ret;
        if(root == rank())
            ret = MPI_Ireduce(MPI_IN_PLACE, inout, len, MPITraits<Type>::getType(),
                              (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm,req);
        else
            ret = MPI_Ireduce(inout, inout, len, MPITraits<Type>::getType(),
                              (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm,req);
        return ret;
    }

#ifdef MASK_MPI_IN_PLACE
#undef MPI_IN_PLACE
#undef MASK_MPI_IN_PLACE
#endif

    template<typename BinaryFunction, typename Type>
    int ireduce(Type* in, Type* out, int len, int root, MPI_Request* req) const
    {
        return MPI_Ireduce(in, out, len, MPITraits<Type>::getType(),
                           (Generic_MPI_Op<Type, BinaryFunction>::get()),root,m_comm,req);
    }
};

#else
// If we compile without MPI, then all we have is the gsSerialComm
typedef gsSerialComm gsMpiComm;
typedef gsSerialGroup gsMpiGroup;
typedef gsSerialStatus gsMpiStatus;
typedef gsSerialRequest gsMpiRequest;
#endif

}