1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
34 * GROwing Monsters And Cloning Shrimps
47 #include "gmx_fatal.h"
52 #include "sighandler.h"
61 #define PP_PME_CHARGE (1<<0)
62 #define PP_PME_CHARGEB (1<<1)
63 #define PP_PME_COORD (1<<2)
64 #define PP_PME_FEP (1<<3)
65 #define PP_PME_ENER_VIR (1<<4)
66 #define PP_PME_FINISH (1<<5)
67 #define PP_PME_SWITCH (1<<6)
69 #define PME_PP_SIGSTOP (1<<0)
70 #define PME_PP_SIGSTOPNSS (1<<1)
72 typedef struct gmx_pme_pp {
74 MPI_Comm mpi_comm_mysim;
76 int nnode; /* The number of PP node to communicate with */
77 int *node; /* The PP node ranks */
78 int node_peer; /* The peer PP node rank */
79 int *nat; /* The number of atom for each PP node */
80 int flags_charge; /* The flags sent along with the last charges */
92 typedef struct gmx_pme_comm_n_box {
100 ivec grid_size; /* For PME grid tuning */
101 real ewaldcoeff; /* For PME grid tuning */
102 } gmx_pme_comm_n_box_t;
109 gmx_stop_cond_t stop_cond;
110 } gmx_pme_comm_vir_ene_t;
115 gmx_pme_pp_t gmx_pme_pp_init(t_commrec *cr)
117 struct gmx_pme_pp *pme_pp;
123 pme_pp->mpi_comm_mysim = cr->mpi_comm_mysim;
124 MPI_Comm_rank(cr->mpi_comm_mygroup,&rank);
125 get_pme_ddnodes(cr,rank,&pme_pp->nnode,&pme_pp->node,&pme_pp->node_peer);
126 snew(pme_pp->nat,pme_pp->nnode);
127 snew(pme_pp->req,2*pme_pp->nnode);
128 snew(pme_pp->stat,2*pme_pp->nnode);
130 pme_pp->flags_charge = 0;
136 /* This should be faster with a real non-blocking MPI implementation */
137 /* #define GMX_PME_DELAYED_WAIT */
139 static void gmx_pme_send_q_x_wait(gmx_domdec_t *dd)
143 MPI_Waitall(dd->nreq_pme,dd->req_pme,MPI_STATUSES_IGNORE);
149 static void gmx_pme_send_q_x(t_commrec *cr, int flags,
150 real *chargeA, real *chargeB,
153 int maxshift_x, int maxshift_y,
154 gmx_large_int_t step)
157 gmx_pme_comm_n_box_t *cnb;
164 fprintf(debug,"PP node %d sending to PME node %d: %d%s%s\n",
165 cr->sim_nodeid,dd->pme_nodeid,n,
166 flags & PP_PME_CHARGE ? " charges" : "",
167 flags & PP_PME_COORD ? " coordinates" : "");
169 #ifdef GMX_PME_DELAYED_WAIT
170 /* When can not use cnb until pending communication has finished */
171 gmx_pme_send_x_q_wait(dd);
174 if (dd->pme_receive_vir_ener) {
175 /* Peer PP node: communicate all data */
182 cnb->maxshift_x = maxshift_x;
183 cnb->maxshift_y = maxshift_y;
184 cnb->lambda = lambda;
186 if (flags & PP_PME_COORD)
187 copy_mat(box,cnb->box);
189 MPI_Isend(cnb,sizeof(*cnb),MPI_BYTE,
190 dd->pme_nodeid,0,cr->mpi_comm_mysim,
191 &dd->req_pme[dd->nreq_pme++]);
193 } else if (flags & PP_PME_CHARGE) {
195 /* Communicate only the number of atoms */
196 MPI_Isend(&n,sizeof(n),MPI_BYTE,
197 dd->pme_nodeid,0,cr->mpi_comm_mysim,
198 &dd->req_pme[dd->nreq_pme++]);
204 if (flags & PP_PME_CHARGE) {
205 MPI_Isend(chargeA,n*sizeof(real),MPI_BYTE,
206 dd->pme_nodeid,1,cr->mpi_comm_mysim,
207 &dd->req_pme[dd->nreq_pme++]);
209 if (flags & PP_PME_CHARGEB) {
210 MPI_Isend(chargeB,n*sizeof(real),MPI_BYTE,
211 dd->pme_nodeid,2,cr->mpi_comm_mysim,
212 &dd->req_pme[dd->nreq_pme++]);
214 if (flags & PP_PME_COORD) {
215 MPI_Isend(x[0],n*sizeof(rvec),MPI_BYTE,
216 dd->pme_nodeid,3,cr->mpi_comm_mysim,
217 &dd->req_pme[dd->nreq_pme++]);
221 #ifndef GMX_PME_DELAYED_WAIT
222 /* Wait for the data to arrive */
223 /* We can skip this wait as we are sure x and q will not be modified
224 * before the next call to gmx_pme_send_x_q or gmx_pme_receive_f.
226 gmx_pme_send_q_x_wait(dd);
231 void gmx_pme_send_q(t_commrec *cr,
232 gmx_bool bFreeEnergy, real *chargeA, real *chargeB,
233 int maxshift_x, int maxshift_y)
237 flags = PP_PME_CHARGE;
239 flags |= PP_PME_CHARGEB;
241 gmx_pme_send_q_x(cr,flags,
242 chargeA,chargeB,NULL,NULL,0,maxshift_x,maxshift_y,-1);
245 void gmx_pme_send_x(t_commrec *cr, matrix box, rvec *x,
246 gmx_bool bFreeEnergy, real lambda,
248 gmx_large_int_t step)
252 flags = PP_PME_COORD;
256 flags |= PP_PME_ENER_VIR;
258 gmx_pme_send_q_x(cr,flags,NULL,NULL,box,x,lambda,0,0,step);
261 void gmx_pme_send_finish(t_commrec *cr)
265 flags = PP_PME_FINISH;
267 gmx_pme_send_q_x(cr,flags,NULL,NULL,NULL,NULL,0,0,0,-1);
270 void gmx_pme_send_switch(t_commrec *cr, ivec grid_size, real ewaldcoeff)
273 gmx_pme_comm_n_box_t cnb;
275 if (cr->dd->pme_receive_vir_ener)
277 cnb.flags = PP_PME_SWITCH;
278 copy_ivec(grid_size,cnb.grid_size);
279 cnb.ewaldcoeff = ewaldcoeff;
281 /* We send this, uncommon, message blocking to simplify the code */
282 MPI_Send(&cnb,sizeof(cnb),MPI_BYTE,
283 cr->dd->pme_nodeid,0,cr->mpi_comm_mysim);
288 int gmx_pme_recv_q_x(struct gmx_pme_pp *pme_pp,
289 real **chargeA, real **chargeB,
290 matrix box, rvec **x,rvec **f,
291 int *maxshift_x, int *maxshift_y,
292 gmx_bool *bFreeEnergy,real *lambda,
294 gmx_large_int_t *step,
295 ivec grid_size, real *ewaldcoeff)
297 gmx_pme_comm_n_box_t cnb;
298 int nat=0,q,messages,sender;
303 /* avoid compiler warning about unused variable without MPI support */
307 /* Receive the send count, box and time step from the peer PP node */
308 MPI_Recv(&cnb,sizeof(cnb),MPI_BYTE,
310 pme_pp->mpi_comm_mysim,MPI_STATUS_IGNORE);
314 fprintf(debug,"PME only node receiving:%s%s%s%s\n",
315 (cnb.flags & PP_PME_CHARGE) ? " charges" : "",
316 (cnb.flags & PP_PME_COORD ) ? " coordinates" : "",
317 (cnb.flags & PP_PME_FINISH) ? " finish" : "",
318 (cnb.flags & PP_PME_SWITCH) ? " switch" : "");
321 if (cnb.flags & PP_PME_SWITCH)
323 /* Special case, receive the new parameters and return */
324 copy_ivec(cnb.grid_size,grid_size);
325 *ewaldcoeff = cnb.ewaldcoeff;
330 if (cnb.flags & PP_PME_CHARGE) {
331 /* Receive the send counts from the other PP nodes */
332 for(sender=0; sender<pme_pp->nnode; sender++) {
333 if (pme_pp->node[sender] == pme_pp->node_peer) {
334 pme_pp->nat[sender] = cnb.natoms;
336 MPI_Irecv(&(pme_pp->nat[sender]),sizeof(pme_pp->nat[0]),
338 pme_pp->node[sender],0,
339 pme_pp->mpi_comm_mysim,&pme_pp->req[messages++]);
342 MPI_Waitall(messages, pme_pp->req, pme_pp->stat);
346 for(sender=0; sender<pme_pp->nnode; sender++)
347 nat += pme_pp->nat[sender];
349 if (nat > pme_pp->nalloc) {
350 pme_pp->nalloc = over_alloc_dd(nat);
351 srenew(pme_pp->chargeA,pme_pp->nalloc);
352 if (cnb.flags & PP_PME_CHARGEB)
353 srenew(pme_pp->chargeB,pme_pp->nalloc);
354 srenew(pme_pp->x,pme_pp->nalloc);
355 srenew(pme_pp->f,pme_pp->nalloc);
358 /* maxshift is sent when the charges are sent */
359 *maxshift_x = cnb.maxshift_x;
360 *maxshift_y = cnb.maxshift_y;
362 /* Receive the charges in place */
363 for(q=0; q<((cnb.flags & PP_PME_CHARGEB) ? 2 : 1); q++) {
365 charge_pp = pme_pp->chargeA;
367 charge_pp = pme_pp->chargeB;
369 for(sender=0; sender<pme_pp->nnode; sender++) {
370 if (pme_pp->nat[sender] > 0) {
371 MPI_Irecv(charge_pp+nat,
372 pme_pp->nat[sender]*sizeof(real),
374 pme_pp->node[sender],1+q,
375 pme_pp->mpi_comm_mysim,
376 &pme_pp->req[messages++]);
377 nat += pme_pp->nat[sender];
379 fprintf(debug,"Received from PP node %d: %d "
381 pme_pp->node[sender],pme_pp->nat[sender]);
386 pme_pp->flags_charge = cnb.flags;
389 if (cnb.flags & PP_PME_COORD) {
390 if (!(pme_pp->flags_charge & PP_PME_CHARGE))
391 gmx_incons("PME-only node received coordinates before charges"
394 /* The box, FE flag and lambda are sent along with the coordinates
396 copy_mat(cnb.box,box);
397 *bFreeEnergy = (cnb.flags & PP_PME_FEP);
398 *lambda = cnb.lambda;
399 *bEnerVir = (cnb.flags & PP_PME_ENER_VIR);
401 if (*bFreeEnergy && !(pme_pp->flags_charge & PP_PME_CHARGEB))
402 gmx_incons("PME-only node received free energy request, but "
403 "did not receive B-state charges");
405 /* Receive the coordinates in place */
407 for(sender=0; sender<pme_pp->nnode; sender++) {
408 if (pme_pp->nat[sender] > 0) {
409 MPI_Irecv(pme_pp->x[nat],pme_pp->nat[sender]*sizeof(rvec),
411 pme_pp->node[sender],3,
412 pme_pp->mpi_comm_mysim,&pme_pp->req[messages++]);
413 nat += pme_pp->nat[sender];
415 fprintf(debug,"Received from PP node %d: %d "
417 pme_pp->node[sender],pme_pp->nat[sender]);
422 /* Wait for the coordinates and/or charges to arrive */
423 MPI_Waitall(messages, pme_pp->req, pme_pp->stat);
425 } while (!(cnb.flags & (PP_PME_COORD | PP_PME_FINISH)));
430 *chargeA = pme_pp->chargeA;
431 *chargeB = pme_pp->chargeB;
436 return ((cnb.flags & PP_PME_FINISH) ? -1 : nat);
439 static void receive_virial_energy(t_commrec *cr,
440 matrix vir,real *energy,real *dvdlambda,
443 gmx_pme_comm_vir_ene_t cve;
445 if (cr->dd->pme_receive_vir_ener) {
448 "PP node %d receiving from PME node %d: virial and energy\n",
449 cr->sim_nodeid,cr->dd->pme_nodeid);
451 MPI_Recv(&cve,sizeof(cve),MPI_BYTE,cr->dd->pme_nodeid,1,cr->mpi_comm_mysim,
454 memset(&cve,0,sizeof(cve));
457 m_add(vir,cve.vir,vir);
458 *energy = cve.energy;
459 *dvdlambda += cve.dvdlambda;
460 *pme_cycles = cve.cycles;
462 if ( cve.stop_cond != gmx_stop_cond_none )
464 gmx_set_stop_condition(cve.stop_cond);
472 void gmx_pme_receive_f(t_commrec *cr,
473 rvec f[], matrix vir,
474 real *energy, real *dvdlambda,
479 #ifdef GMX_PME_DELAYED_WAIT
480 /* Wait for the x request to finish */
481 gmx_pme_send_q_x_wait(cr->dd);
484 natoms = cr->dd->nat_home;
486 if (natoms > cr->dd->pme_recv_f_alloc)
488 cr->dd->pme_recv_f_alloc = over_alloc_dd(natoms);
489 srenew(cr->dd->pme_recv_f_buf, cr->dd->pme_recv_f_alloc);
493 MPI_Recv(cr->dd->pme_recv_f_buf[0],
494 natoms*sizeof(rvec),MPI_BYTE,
495 cr->dd->pme_nodeid,0,cr->mpi_comm_mysim,
499 for(i=0; i<natoms; i++)
500 rvec_inc(f[i],cr->dd->pme_recv_f_buf[i]);
503 receive_virial_energy(cr,vir,energy,dvdlambda,pme_cycles);
506 void gmx_pme_send_force_vir_ener(struct gmx_pme_pp *pme_pp,
508 real energy, real dvdlambda,
511 gmx_pme_comm_vir_ene_t cve;
512 int messages,ind_start,ind_end,receiver;
516 /* Now the evaluated forces have to be transferred to the PP nodes */
519 for (receiver=0; receiver<pme_pp->nnode; receiver++) {
521 ind_end = ind_start + pme_pp->nat[receiver];
523 if (MPI_Isend(f[ind_start],(ind_end-ind_start)*sizeof(rvec),MPI_BYTE,
524 pme_pp->node[receiver],0,
525 pme_pp->mpi_comm_mysim,&pme_pp->req[messages++]) != 0)
526 gmx_comm("MPI_Isend failed in do_pmeonly");
530 /* send virial and energy to our last PP node */
531 copy_mat(vir,cve.vir);
533 cve.dvdlambda = dvdlambda;
534 /* check for the signals to send back to a PP node */
535 cve.stop_cond = gmx_get_stop_condition();
540 fprintf(debug,"PME node sending to PP node %d: virial and energy\n",
543 MPI_Isend(&cve,sizeof(cve),MPI_BYTE,
545 pme_pp->mpi_comm_mysim,&pme_pp->req[messages++]);
547 /* Wait for the forces to arrive */
548 MPI_Waitall(messages, pme_pp->req, pme_pp->stat);