1 /* 2 * This file is part of the GROMACS molecular simulation package. 3 * 4 * Copyright (c) 2008-2018, The GROMACS development team. 5 * Copyright (c) 2019,2020, by the GROMACS development team, led by 6 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, 7 * and including many others, as listed in the AUTHORS file in the 8 * top-level source directory and at http://www.gromacs.org. 9 * 10 * GROMACS is free software; you can redistribute it and/or 11 * modify it under the terms of the GNU Lesser General Public License 12 * as published by the Free Software Foundation; either version 2.1 13 * of the License, or (at your option) any later version. 14 * 15 * GROMACS is distributed in the hope that it will be useful, 16 * but WITHOUT ANY WARRANTY; without even the implied warranty of 17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 18 * Lesser General Public License for more details. 19 * 20 * You should have received a copy of the GNU Lesser General Public 21 * License along with GROMACS; if not, see 22 * http://www.gnu.org/licenses, or write to the Free Software Foundation, 23 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. 24 * 25 * If you want to redistribute modifications to GROMACS, please 26 * consider that scientific software is very special. Version 27 * control is crucial - bugs must be traceable. We will be happy to 28 * consider code for inclusion in the official distribution, but 29 * derived work must not be called official GROMACS. Details are found 30 * in the README & COPYING files - if they are missing, get the 31 * official version at http://www.gromacs.org. 32 * 33 * To help us fund GROMACS development, we humbly ask that you cite 34 * the research papers on the package. Check out http://www.gromacs.org. 35 */ 36 37 /*! \libinternal \file 38 * 39 * \brief This file declares functions for (mostly) the domdec module 40 * to use MPI functionality 41 * 42 * \todo Wrap the raw dd_bcast in md.cpp into a higher-level function 43 * in the domdec module, then this file can be module-internal. 44 * 45 * \author Berk Hess <hess@kth.se> 46 * \inlibraryapi 47 * \ingroup module_domdec 48 */ 49 #ifndef GMX_DOMDEC_DOMDEC_NETWORK_H 50 #define GMX_DOMDEC_DOMDEC_NETWORK_H 51 52 #include "gromacs/math/vectypes.h" 53 54 struct gmx_domdec_t; 55 56 namespace gmx 57 { 58 template<typename> 59 class ArrayRef; 60 } 61 /* \brief */ 62 enum 63 { 64 dddirForward, 65 dddirBackward 66 }; 67 68 /*! \brief Move T values in the communication region one cell along 69 * the domain decomposition 70 * 71 * Moves in the dimension indexed by ddDimensionIndex, either forward 72 * (direction=dddirFoward) or backward (direction=dddirBackward). 73 * 74 * \todo This function template is deprecated, new calls should be 75 * made to the version taking ArrayRef parameters and this function 76 * template removed when unused. 77 */ 78 template<typename T> 79 void ddSendrecv(const gmx_domdec_t* dd, 80 int ddDimensionIndex, 81 int direction, 82 T* sendBuffer, 83 int numElementsToSend, 84 T* receiveBuffer, 85 int numElementsToReceive); 86 87 //! Extern declaration for int specialization 88 extern template void ddSendrecv<int>(const gmx_domdec_t* dd, 89 int ddDimensionIndex, 90 int direction, 91 int* buf_s, 92 int n_s, 93 int* buf_r, 94 int n_r); 95 96 //! Extern declaration for real specialization 97 extern template void ddSendrecv<real>(const gmx_domdec_t* dd, 98 int ddDimensionIndex, 99 int direction, 100 real* buf_s, 101 int n_s, 102 real* buf_r, 103 int n_r); 104 105 //! Extern declaration for rvec specialization 106 extern template void ddSendrecv<rvec>(const gmx_domdec_t* dd, 107 int ddDimensionIndex, 108 int direction, 109 rvec* buf_s, 110 int n_s, 111 rvec* buf_r, 112 int n_r); 113 114 /*! \brief Move a view of T values in the communication region one 115 * cell along the domain decomposition 116 * 117 * Moves in the dimension indexed by ddDimensionIndex, either forward 118 * (direction=dddirFoward) or backward (direction=dddirBackward). 119 */ 120 template<typename T> 121 void ddSendrecv(const gmx_domdec_t* dd, 122 int ddDimensionIndex, 123 int direction, 124 gmx::ArrayRef<T> sendBuffer, 125 gmx::ArrayRef<T> receiveBuffer); 126 127 //! Extern declaration for int specialization 128 extern template void ddSendrecv<int>(const gmx_domdec_t* dd, 129 int ddDimensionIndex, 130 int direction, 131 gmx::ArrayRef<int> sendBuffer, 132 gmx::ArrayRef<int> receiveBuffer); 133 134 //! Extern declaration for real specialization 135 extern template void ddSendrecv<real>(const gmx_domdec_t* dd, 136 int ddDimensionIndex, 137 int direction, 138 gmx::ArrayRef<real> sendBuffer, 139 gmx::ArrayRef<real> receiveBuffer); 140 141 //! Extern declaration for gmx::RVec specialization 142 extern template void ddSendrecv<gmx::RVec>(const gmx_domdec_t* dd, 143 int ddDimensionIndex, 144 int direction, 145 gmx::ArrayRef<gmx::RVec> sendBuffer, 146 gmx::ArrayRef<gmx::RVec> receiveBuffer); 147 148 /*! \brief Move revc's in the comm. region one cell along the domain decomposition 149 * 150 * Moves in dimension indexed by ddimind, simultaneously in the forward 151 * and backward directions. 152 */ 153 void dd_sendrecv2_rvec(const struct gmx_domdec_t* dd, 154 int ddimind, 155 rvec* buf_s_fw, 156 int n_s_fw, 157 rvec* buf_r_fw, 158 int n_r_fw, 159 rvec* buf_s_bw, 160 int n_s_bw, 161 rvec* buf_r_bw, 162 int n_r_bw); 163 164 165 /* The functions below perform the same operations as the MPI functions 166 * with the same name appendices, but over the domain decomposition 167 * nodes only. 168 * The DD master node is the master for these operations. 169 */ 170 171 /*! \brief Broadcasts \p nbytes from \p data on \p DDMASTERRANK to all PP ranks */ 172 void dd_bcast(const gmx_domdec_t* dd, int nbytes, void* data); 173 174 /*! \brief Copies \p nbytes from \p src to \p dest on \p DDMASTERRANK 175 * and then broadcasts to \p dest on all PP ranks */ 176 void dd_bcastc(const gmx_domdec_t* dd, int nbytes, void* src, void* dest); 177 178 /*! \brief Scatters \p nbytes from \p src on \p DDMASTERRANK to all PP ranks, received in \p dest */ 179 void dd_scatter(const gmx_domdec_t* dd, int nbytes, const void* src, void* dest); 180 181 /*! \brief Gathers \p nbytes from \p src on all PP ranks, received in \p dest on \p DDMASTERRANK */ 182 void dd_gather(const gmx_domdec_t* dd, int nbytes, const void* src, void* dest); 183 184 /*! \brief Scatters \p scounts bytes from \p src on \p DDMASTERRANK to all PP ranks, receiving \p rcount bytes in \p dest. 185 * 186 * See man MPI_Scatterv for details of how to construct scounts and disps. 187 * If rcount==0, rbuf is allowed to be NULL */ 188 void dd_scatterv(const gmx_domdec_t* dd, int* scounts, int* disps, const void* sbuf, int rcount, void* rbuf); 189 190 /*! \brief Gathers \p rcount bytes from \p src on all PP ranks, received in \p scounts bytes in \p dest on \p DDMASTERRANK. 191 * 192 * See man MPI_Gatherv for details of how to construct scounts and disps. 193 * 194 * If scount==0, sbuf is allowed to be NULL */ 195 void dd_gatherv(const gmx_domdec_t* dd, int scount, const void* sbuf, int* rcounts, int* disps, void* rbuf); 196 197 #endif 198