| 2 |  | * @file Communicator.hpp | 
| 3 |  | * @author Charles Vardeman <cvardema.at.nd.edu> | 
| 4 |  | * @date 08/18/2010 | 
| 5 | – | * @time 11:56am | 
| 5 |  | * @version 1.0 | 
| 6 |  | * | 
| 7 |  | * @section LICENSE | 
| 41 |  | * | 
| 42 |  | * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). | 
| 43 |  | * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). | 
| 44 | < | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). | 
| 45 | < | * [4]  Vardeman & Gezelter, in progress (2009). | 
| 44 | > | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). | 
| 45 | > | * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010). | 
| 46 | > | * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). | 
| 47 |  | */ | 
| 48 |  |  | 
| 49 | < | #ifndef FORCEDECOMPOSITION_COMMUNICATOR_HPP | 
| 50 | < | #define FORCEDECOMPOSITION_COMMUNICATOR_HPP | 
| 49 | > | #ifndef PARALLEL_COMMUNICATOR_HPP | 
| 50 | > | #define PARALLEL_COMMUNICATOR_HPP | 
| 51 |  |  | 
| 52 |  | #include <config.h> | 
| 53 |  | #include <mpi.h> | 
| 54 |  | #include "math/SquareMatrix3.hpp" | 
| 55 |  |  | 
| 56 | + | using namespace std; | 
| 57 |  | namespace OpenMD{ | 
| 58 |  |  | 
| 59 |  | #ifdef IS_MPI | 
| 60 |  |  | 
| 61 | < | enum direction { | 
| 62 | < | Row = 0, | 
| 63 | < | Column = 1 | 
| 61 | > | enum communicatorType { | 
| 62 | > | Global = 0, | 
| 63 | > | Row = 1, | 
| 64 | > | Column = 2 | 
| 65 |  | }; | 
| 66 |  |  | 
| 67 | < | template<typename T> | 
| 68 | < | struct MPITraits | 
| 69 | < | { | 
| 70 | < | static const MPI::Datatype datatype; | 
| 71 | < | static const int dim; | 
| 67 | > | template<class T> | 
| 68 | > | class MPITraits { | 
| 69 | > | public: | 
| 70 | > | static MPI_Datatype Type(); | 
| 71 | > | static int Length() { return 1; }; | 
| 72 |  | }; | 
| 73 |  |  | 
| 74 | < | template<> const MPI::Datatype MPITraits<int>::datatype = MPI_INT; | 
| 75 | < | template<> const int MPITraits<int>::dim = 1; | 
| 76 | < | template<> const MPI::Datatype MPITraits<RealType>::datatype = MPI_REALTYPE; | 
| 77 | < | template<> const int MPITraits<RealType>::dim = 1; | 
| 78 | < | template<> const MPI::Datatype MPITraits<Vector3d>::datatype = MPI_REALTYPE; | 
| 79 | < | template<> const int MPITraits<Vector3d>::dim = 3; | 
| 80 | < | template<> const MPI::Datatype MPITraits<Mat3x3d>::datatype = MPI_REALTYPE; | 
| 81 | < | template<> const int MPITraits<Mat3x3d>::dim = 9; | 
| 74 | > | template<> inline MPI_Datatype MPITraits<int>::Type() { return MPI_INT; } | 
| 75 | > | template<> inline MPI_Datatype MPITraits<RealType>::Type() { return MPI_REALTYPE; } | 
| 76 | > |  | 
| 77 | > | template<class T, unsigned int Dim> | 
| 78 | > | class MPITraits< Vector<T, Dim> > { | 
| 79 | > | public: | 
| 80 | > | static MPI_Datatype Type() { return MPITraits<T>::Type(); } | 
| 81 | > | static int Length() {return Dim;} | 
| 82 | > | }; | 
| 83 | > |  | 
| 84 | > | template<class T> | 
| 85 | > | class MPITraits< Vector3<T> > { | 
| 86 | > | public: | 
| 87 | > | static MPI_Datatype Type() { return MPITraits<T>::Type(); } | 
| 88 | > | static int Length() {return 3;} | 
| 89 | > | }; | 
| 90 | > |  | 
| 91 | > | template<class T, unsigned int R, unsigned int C> | 
| 92 | > | class MPITraits< RectMatrix<T, R, C> > { | 
| 93 | > | public: | 
| 94 | > | static MPI_Datatype Type() { return MPITraits<T>::Type(); } | 
| 95 | > | static int Length() {return R * C;} | 
| 96 | > | }; | 
| 97 | > |  | 
| 98 | > | template<class T> | 
| 99 | > | class MPITraits< SquareMatrix3<T> > { | 
| 100 | > | public: | 
| 101 | > | static MPI_Datatype Type() { return MPITraits<T>::Type(); } | 
| 102 | > | static int Length() {return 9;} | 
| 103 | > | }; | 
| 104 |  |  | 
| 105 | < | template<direction D, typename T> | 
| 105 | > |  | 
| 106 | > | template<communicatorType D> | 
| 107 |  | class Communicator { | 
| 108 |  | public: | 
| 109 |  |  | 
| 110 | < | Communicator<D, T>(int nObjects) { | 
| 110 | > | Communicator<D>() { | 
| 111 |  |  | 
| 112 | < | int nProc = MPI::COMM_WORLD.Get_size(); | 
| 113 | < | int myRank = MPI::COMM_WORLD.Get_rank(); | 
| 112 | > | int nProc; | 
| 113 | > | int myRank; | 
| 114 |  |  | 
| 115 | + | MPI_Comm_size( MPI_COMM_WORLD, &nProc); | 
| 116 | + | MPI_Comm_rank( MPI_COMM_WORLD, &myRank); | 
| 117 | + |  | 
| 118 |  | int nColumnsMax = (int) sqrt(RealType(nProc)); | 
| 119 |  |  | 
| 120 |  | int nColumns; | 
| 122 |  | if (nProc % i == 0) nColumns = i; | 
| 123 |  | } | 
| 124 |  |  | 
| 125 | < | int nRows = nProc / nColumns; | 
| 125 | > | // int nRows = nProc / nColumns; | 
| 126 |  | rowIndex_ = myRank / nColumns; | 
| 127 |  | columnIndex_ = myRank % nColumns; | 
| 128 |  |  | 
| 129 | < | if (D == Row) { | 
| 130 | < | myComm = MPI::COMM_WORLD.Split(rowIndex_, 0); | 
| 131 | < | } else { | 
| 132 | < | myComm = MPI::COMM_WORLD.Split(columnIndex_, 0); | 
| 129 | > | switch(D) { | 
| 130 | > | case Row : | 
| 131 | > | MPI_Comm_split(MPI_COMM_WORLD, rowIndex_, 0, &myComm); | 
| 132 | > | //  myComm = MPI::COMM_WORLD.Split(rowIndex_, 0); | 
| 133 | > | break; | 
| 134 | > | case Column: | 
| 135 | > | MPI_Comm_split(MPI_COMM_WORLD, columnIndex_, 0, &myComm); | 
| 136 | > | //myComm = MPI::COMM_WORLD.Split(columnIndex_, 0); | 
| 137 | > | break; | 
| 138 | > | case Global: | 
| 139 | > | MPI_Comm_split(MPI_COMM_WORLD, myRank, 0, &myComm); | 
| 140 | > | //myComm = MPI::COMM_WORLD.Split(myRank, 0); | 
| 141 |  | } | 
| 106 | – |  | 
| 107 | – | int nCommProcs = myComm.Get_size(); | 
| 142 |  |  | 
| 143 | < | counts.reserve(nCommProcs); | 
| 144 | < | displacements.reserve(nCommProcs); | 
| 143 | > | } | 
| 144 | > |  | 
| 145 | > | MPI_Comm getComm() { return myComm; } | 
| 146 | > |  | 
| 147 | > | private: | 
| 148 | > | int rowIndex_; | 
| 149 | > | int columnIndex_; | 
| 150 | > | MPI_Comm myComm; | 
| 151 | > | }; | 
| 152 | > |  | 
| 153 |  |  | 
| 154 | < | planSize_ = MPITraits<T>::dim * nObjects; | 
| 154 | > | template<typename T> | 
| 155 | > | class Plan { | 
| 156 | > | public: | 
| 157 | > |  | 
| 158 | > | Plan<T>(MPI_Comm comm, int nObjects) : myComm(comm) { | 
| 159 |  |  | 
| 160 | < | myComm.Allgather(&planSize_, 1, MPI::INT, &counts[0], 1, MPI::INT); | 
| 160 | > | int nCommProcs; | 
| 161 | > | MPI_Comm_size( myComm, &nCommProcs ); | 
| 162 |  |  | 
| 163 | < |  | 
| 163 | > | //int nCommProcs = myComm.Get_size(); | 
| 164 | > |  | 
| 165 | > | counts.resize(nCommProcs, 0); | 
| 166 | > | displacements.resize(nCommProcs, 0); | 
| 167 | > |  | 
| 168 | > | planSize_ = MPITraits<T>::Length() * nObjects; | 
| 169 | > |  | 
| 170 | > | MPI_Allgather(&planSize_, 1, MPI_INT, &counts[0], 1, MPI_INT, myComm); | 
| 171 | > | //myComm.Allgather(&planSize_, 1, MPI_INT, &counts[0], 1, MPI_INT); | 
| 172 | > |  | 
| 173 |  | displacements[0] = 0; | 
| 174 |  | for (int i = 1; i < nCommProcs; i++) { | 
| 175 |  | displacements[i] = displacements[i-1] + counts[i-1]; | 
| 120 | – | size_ += count[i-1]; | 
| 176 |  | } | 
| 177 |  |  | 
| 178 |  | size_ = 0; | 
| 181 |  | } | 
| 182 |  | } | 
| 183 |  |  | 
| 184 | < |  | 
| 185 | < | void gather(std::vector<T>& v1, std::vector<T>& v2) { | 
| 184 | > |  | 
| 185 | > | void gather(vector<T>& v1, vector<T>& v2) { | 
| 186 |  |  | 
| 187 | < | myComm.Allgatherv(&v1[0], | 
| 188 | < | planSize_, | 
| 189 | < | MPITraits<T>::datatype, | 
| 190 | < | &v2[0], | 
| 191 | < | &counts[0], | 
| 192 | < | &displacements[0], | 
| 193 | < | MPITraits<T>::datatype); | 
| 194 | < | } | 
| 195 | < |  | 
| 187 | > | // an assert would be helpful here to make sure the vectors are the | 
| 188 | > | // correct geometry | 
| 189 | > |  | 
| 190 | > | MPI_Allgatherv(&v1[0], | 
| 191 | > | planSize_, | 
| 192 | > | MPITraits<T>::Type(), | 
| 193 | > | &v2[0], | 
| 194 | > | &counts[0], | 
| 195 | > | &displacements[0], | 
| 196 | > | MPITraits<T>::Type(), | 
| 197 | > | myComm); | 
| 198 | > | // myComm.Allgatherv(&v1[0], | 
| 199 | > | //                   planSize_, | 
| 200 | > | //                   MPITraits<T>::Type(), | 
| 201 | > | //                   &v2[0], | 
| 202 | > | //                   &counts[0], | 
| 203 | > | //                   &displacements[0], | 
| 204 | > | //                   MPITraits<T>::Type()); | 
| 205 | > | } | 
| 206 |  |  | 
| 207 | < |  | 
| 208 | < | void scatter(std::vector<T>& v1, std::vector<T>& v2) { | 
| 207 | > | void scatter(vector<T>& v1, vector<T>& v2) { | 
| 208 | > | // an assert would be helpful here to make sure the vectors are the | 
| 209 | > | // correct geometry | 
| 210 | > |  | 
| 211 | > | MPI_Reduce_scatter(&v1[0], &v2[0], &counts[0], | 
| 212 | > | MPITraits<T>::Type(), MPI_SUM, myComm); | 
| 213 |  |  | 
| 214 | < | myComm.Reduce_scatter(&v1[0], &v2[0], &counts[0], | 
| 215 | < | MPITraits<T>::datatype, MPI::SUM); | 
| 214 | > | // myComm.Reduce_scatter(&v1[0], &v2[0], &counts[0], | 
| 215 | > | //                       MPITraits<T>::Type(), MPI::SUM); | 
| 216 |  | } | 
| 217 | < |  | 
| 217 | > |  | 
| 218 |  | int getSize() { | 
| 219 |  | return size_; | 
| 220 |  | } | 
| 221 |  |  | 
| 222 |  | private: | 
| 223 |  | int planSize_;     ///< how many are on local proc | 
| 155 | – | int rowIndex_; | 
| 156 | – | int columnIndex_; | 
| 224 |  | int size_; | 
| 225 | < | std::vector<int> counts; | 
| 226 | < | std::vector<int> displacements; | 
| 227 | < | MPI::Intracomm myComm; | 
| 225 | > | vector<int> counts; | 
| 226 | > | vector<int> displacements; | 
| 227 | > | MPI_Comm myComm; | 
| 228 |  | }; | 
| 229 |  |  | 
| 230 |  | #endif |