src/gromacs/mdlib/pme_pp.c

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2004, The GROMACS development team.
   6  * Copyright (c) 2013,2014, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37
  38 #include "gmxpre.h"
  39
  40 #include "config.h"
  41
  42
  43 #include <stdio.h>
  44 #include <string.h>
  45 #include <math.h>
  46 #include "gromacs/legacyheaders/typedefs.h"
  47 #include "gromacs/legacyheaders/types/commrec.h"
  48 #include "gromacs/utility/smalloc.h"
  49 #include "gromacs/utility/fatalerror.h"
  50 #include "gromacs/math/vec.h"
  51 #include "gromacs/legacyheaders/pme.h"
  52 #include "gromacs/legacyheaders/network.h"
  53 #include "gromacs/legacyheaders/domdec.h"
  54 #include "gromacs/legacyheaders/sighandler.h"
  55
  56 #include "gromacs/utility/gmxmpi.h"
  57
  58 enum {
  59     eCommType_ChargeA, eCommType_ChargeB, eCommType_SQRTC6A, eCommType_SQRTC6B,
  60     eCommType_SigmaA, eCommType_SigmaB, eCommType_NR, eCommType_COORD,
  61     eCommType_CNB
  62 };
  63
  64 /* Some parts of the code(gmx_pme_send_q, gmx_pme_recv_q_x) assume
  65  * that the six first flags are exactly in this order.
  66  * If more PP_PME_...-flags are to be introduced be aware of some of
  67  * the PME-specific flags in pme.h. Currently, they are also passed
  68  * through here.
  69  */
  70
  71 #define PP_PME_CHARGE         (1<<0)
  72 #define PP_PME_CHARGEB        (1<<1)
  73 #define PP_PME_SQRTC6         (1<<2)
  74 #define PP_PME_SQRTC6B        (1<<3)
  75 #define PP_PME_SIGMA          (1<<4)
  76 #define PP_PME_SIGMAB         (1<<5)
  77 #define PP_PME_COORD          (1<<6)
  78 #define PP_PME_FEP_Q          (1<<7)
  79 #define PP_PME_FEP_LJ         (1<<8)
  80 #define PP_PME_ENER_VIR       (1<<9)
  81 #define PP_PME_FINISH         (1<<10)
  82 #define PP_PME_SWITCHGRID     (1<<11)
  83 #define PP_PME_RESETCOUNTERS  (1<<12)
  84
  85 #define PME_PP_SIGSTOP        (1<<0)
  86 #define PME_PP_SIGSTOPNSS     (1<<1)
  87
  88 typedef struct gmx_pme_pp {
  89 #ifdef GMX_MPI
  90     MPI_Comm     mpi_comm_mysim;
  91 #endif
  92     int          nnode;        /* The number of PP node to communicate with  */
  93     int         *node;         /* The PP node ranks                          */
  94     int          node_peer;    /* The peer PP node rank                      */
  95     int         *nat;          /* The number of atom for each PP node        */
  96     int          flags_charge; /* The flags sent along with the last charges */
  97     real        *chargeA;
  98     real        *chargeB;
  99     real        *sqrt_c6A;
 100     real        *sqrt_c6B;
 101     real        *sigmaA;
 102     real        *sigmaB;
 103     rvec        *x;
 104     rvec        *f;
 105     int          nalloc;
 106 #ifdef GMX_MPI
 107     MPI_Request *req;
 108     MPI_Status  *stat;
 109 #endif
 110 } t_gmx_pme_pp;
 111
 112 typedef struct gmx_pme_comm_n_box {
 113     int             natoms;
 114     matrix          box;
 115     int             maxshift_x;
 116     int             maxshift_y;
 117     real            lambda_q;
 118     real            lambda_lj;
 119     int             flags;
 120     gmx_int64_t     step;
 121     ivec            grid_size;    /* For PME grid tuning */
 122     real            ewaldcoeff_q; /* For PME grid tuning */
 123     real            ewaldcoeff_lj;
 124 } gmx_pme_comm_n_box_t;
 125
 126 typedef struct {
 127     matrix          vir_q;
 128     matrix          vir_lj;
 129     real            energy_q;
 130     real            energy_lj;
 131     real            dvdlambda_q;
 132     real            dvdlambda_lj;
 133     float           cycles;
 134     gmx_stop_cond_t stop_cond;
 135 } gmx_pme_comm_vir_ene_t;
 136
 137
 138
 139
 140 gmx_pme_pp_t gmx_pme_pp_init(t_commrec gmx_unused *cr)
 141 {
 142     struct gmx_pme_pp *pme_pp;
 143     int                rank;
 144
 145     snew(pme_pp, 1);
 146
 147 #ifdef GMX_MPI
 148     pme_pp->mpi_comm_mysim = cr->mpi_comm_mysim;
 149     MPI_Comm_rank(cr->mpi_comm_mygroup, &rank);
 150     get_pme_ddnodes(cr, rank, &pme_pp->nnode, &pme_pp->node, &pme_pp->node_peer);
 151     snew(pme_pp->nat, pme_pp->nnode);
 152     snew(pme_pp->req, eCommType_NR*pme_pp->nnode);
 153     snew(pme_pp->stat, eCommType_NR*pme_pp->nnode);
 154     pme_pp->nalloc       = 0;
 155     pme_pp->flags_charge = 0;
 156 #endif
 157
 158     return pme_pp;
 159 }
 160
 161 /* This should be faster with a real non-blocking MPI implementation */
 162 /* #define GMX_PME_DELAYED_WAIT */
 163
 164 static void gmx_pme_send_coeffs_coords_wait(gmx_domdec_t gmx_unused *dd)
 165 {
 166 #ifdef GMX_MPI
 167     if (dd->nreq_pme)
 168     {
 169         MPI_Waitall(dd->nreq_pme, dd->req_pme, MPI_STATUSES_IGNORE);
 170         dd->nreq_pme = 0;
 171     }
 172 #endif
 173 }
 174
 175 static void gmx_pme_send_coeffs_coords(t_commrec *cr, int flags,
 176                                        real gmx_unused *chargeA, real gmx_unused *chargeB,
 177                                        real gmx_unused *c6A, real gmx_unused *c6B,
 178                                        real gmx_unused *sigmaA, real gmx_unused *sigmaB,
 179                                        matrix box, rvec gmx_unused *x,
 180                                        real lambda_q, real lambda_lj,
 181                                        int maxshift_x, int maxshift_y,
 182                                        gmx_int64_t step)
 183 {
 184     gmx_domdec_t         *dd;
 185     gmx_pme_comm_n_box_t *cnb;
 186     int                   n;
 187
 188     dd = cr->dd;
 189     n  = dd->nat_home;
 190
 191     if (debug)
 192     {
 193         fprintf(debug, "PP rank %d sending to PME rank %d: %d%s%s%s%s\n",
 194                 cr->sim_nodeid, dd->pme_nodeid, n,
 195                 flags & PP_PME_CHARGE ? " charges" : "",
 196                 flags & PP_PME_SQRTC6 ? " sqrtC6" : "",
 197                 flags & PP_PME_SIGMA  ? " sigma" : "",
 198                 flags & PP_PME_COORD  ? " coordinates" : "");
 199     }
 200
 201 #ifdef GMX_PME_DELAYED_WAIT
 202     /* When can not use cnb until pending communication has finished */
 203     gmx_pme_send_coeffs_coords_wait(dd);
 204 #endif
 205
 206     if (dd->pme_receive_vir_ener)
 207     {
 208         /* Peer PP node: communicate all data */
 209         if (dd->cnb == NULL)
 210         {
 211             snew(dd->cnb, 1);
 212         }
 213         cnb = dd->cnb;
 214
 215         cnb->flags      = flags;
 216         cnb->natoms     = n;
 217         cnb->maxshift_x = maxshift_x;
 218         cnb->maxshift_y = maxshift_y;
 219         cnb->lambda_q   = lambda_q;
 220         cnb->lambda_lj  = lambda_lj;
 221         cnb->step       = step;
 222         if (flags & PP_PME_COORD)
 223         {
 224             copy_mat(box, cnb->box);
 225         }
 226 #ifdef GMX_MPI
 227         MPI_Isend(cnb, sizeof(*cnb), MPI_BYTE,
 228                   dd->pme_nodeid, eCommType_CNB, cr->mpi_comm_mysim,
 229                   &dd->req_pme[dd->nreq_pme++]);
 230 #endif
 231     }
 232     else if (flags & (PP_PME_CHARGE | PP_PME_SQRTC6 | PP_PME_SIGMA))
 233     {
 234 #ifdef GMX_MPI
 235         /* Communicate only the number of atoms */
 236         MPI_Isend(&n, sizeof(n), MPI_BYTE,
 237                   dd->pme_nodeid, eCommType_CNB, cr->mpi_comm_mysim,
 238                   &dd->req_pme[dd->nreq_pme++]);
 239 #endif
 240     }
 241
 242 #ifdef GMX_MPI
 243     if (n > 0)
 244     {
 245         if (flags & PP_PME_CHARGE)
 246         {
 247             MPI_Isend(chargeA, n*sizeof(real), MPI_BYTE,
 248                       dd->pme_nodeid, eCommType_ChargeA, cr->mpi_comm_mysim,
 249                       &dd->req_pme[dd->nreq_pme++]);
 250         }
 251         if (flags & PP_PME_CHARGEB)
 252         {
 253             MPI_Isend(chargeB, n*sizeof(real), MPI_BYTE,
 254                       dd->pme_nodeid, eCommType_ChargeB, cr->mpi_comm_mysim,
 255                       &dd->req_pme[dd->nreq_pme++]);
 256         }
 257         if (flags & PP_PME_SQRTC6)
 258         {
 259             MPI_Isend(c6A, n*sizeof(real), MPI_BYTE,
 260                       dd->pme_nodeid, eCommType_SQRTC6A, cr->mpi_comm_mysim,
 261                       &dd->req_pme[dd->nreq_pme++]);
 262         }
 263         if (flags & PP_PME_SQRTC6B)
 264         {
 265             MPI_Isend(c6B, n*sizeof(real), MPI_BYTE,
 266                       dd->pme_nodeid, eCommType_SQRTC6B, cr->mpi_comm_mysim,
 267                       &dd->req_pme[dd->nreq_pme++]);
 268         }
 269         if (flags & PP_PME_SIGMA)
 270         {
 271             MPI_Isend(sigmaA, n*sizeof(real), MPI_BYTE,
 272                       dd->pme_nodeid, eCommType_SigmaA, cr->mpi_comm_mysim,
 273                       &dd->req_pme[dd->nreq_pme++]);
 274         }
 275         if (flags & PP_PME_SIGMAB)
 276         {
 277             MPI_Isend(sigmaB, n*sizeof(real), MPI_BYTE,
 278                       dd->pme_nodeid, eCommType_SigmaB, cr->mpi_comm_mysim,
 279                       &dd->req_pme[dd->nreq_pme++]);
 280         }
 281         if (flags & PP_PME_COORD)
 282         {
 283             MPI_Isend(x[0], n*sizeof(rvec), MPI_BYTE,
 284                       dd->pme_nodeid, eCommType_COORD, cr->mpi_comm_mysim,
 285                       &dd->req_pme[dd->nreq_pme++]);
 286         }
 287     }
 288
 289 #ifndef GMX_PME_DELAYED_WAIT
 290     /* Wait for the data to arrive */
 291     /* We can skip this wait as we are sure x and q will not be modified
 292      * before the next call to gmx_pme_send_x_q or gmx_pme_receive_f.
 293      */
 294     gmx_pme_send_coeffs_coords_wait(dd);
 295 #endif
 296 #endif
 297 }
 298
 299 void gmx_pme_send_parameters(t_commrec *cr,
 300                              const interaction_const_t *ic,
 301                              gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj,
 302                              real *chargeA, real *chargeB,
 303                              real *sqrt_c6A, real *sqrt_c6B,
 304                              real *sigmaA, real *sigmaB,
 305                              int maxshift_x, int maxshift_y)
 306 {
 307     int flags;
 308
 309     flags = 0;
 310     if (EEL_PME(ic->eeltype))
 311     {
 312         flags |= PP_PME_CHARGE;
 313     }
 314     if (EVDW_PME(ic->vdwtype))
 315     {
 316         flags |= (PP_PME_SQRTC6 | PP_PME_SIGMA);
 317     }
 318     if (bFreeEnergy_q || bFreeEnergy_lj)
 319     {
 320         /* Assumes that the B state flags are in the bits just above
 321          * the ones for the A state. */
 322         flags |= (flags << 1);
 323     }
 324
 325     gmx_pme_send_coeffs_coords(cr, flags,
 326                                chargeA, chargeB,
 327                                sqrt_c6A, sqrt_c6B, sigmaA, sigmaB,
 328                                NULL, NULL, 0, 0, maxshift_x, maxshift_y, -1);
 329 }
 330
 331 void gmx_pme_send_coordinates(t_commrec *cr, matrix box, rvec *x,
 332                               gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj,
 333                               real lambda_q, real lambda_lj,
 334                               gmx_bool bEnerVir, int pme_flags,
 335                               gmx_int64_t step)
 336 {
 337     int flags;
 338
 339     flags = pme_flags | PP_PME_COORD;
 340     if (bFreeEnergy_q)
 341     {
 342         flags |= PP_PME_FEP_Q;
 343     }
 344     if (bFreeEnergy_lj)
 345     {
 346         flags |= PP_PME_FEP_LJ;
 347     }
 348     if (bEnerVir)
 349     {
 350         flags |= PP_PME_ENER_VIR;
 351     }
 352     gmx_pme_send_coeffs_coords(cr, flags, NULL, NULL, NULL, NULL, NULL, NULL,
 353                                box, x, lambda_q, lambda_lj, 0, 0, step);
 354 }
 355
 356 void gmx_pme_send_finish(t_commrec *cr)
 357 {
 358     int flags;
 359
 360     flags = PP_PME_FINISH;
 361
 362     gmx_pme_send_coeffs_coords(cr, flags, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, 0, 0, 0, 0, -1);
 363 }
 364
 365 void gmx_pme_send_switchgrid(t_commrec gmx_unused *cr,
 366                              ivec gmx_unused       grid_size,
 367                              real gmx_unused       ewaldcoeff_q,
 368                              real gmx_unused       ewaldcoeff_lj)
 369 {
 370 #ifdef GMX_MPI
 371     gmx_pme_comm_n_box_t cnb;
 372
 373     /* Only let one PP node signal each PME node */
 374     if (cr->dd->pme_receive_vir_ener)
 375     {
 376         cnb.flags = PP_PME_SWITCHGRID;
 377         copy_ivec(grid_size, cnb.grid_size);
 378         cnb.ewaldcoeff_q  = ewaldcoeff_q;
 379         cnb.ewaldcoeff_lj = ewaldcoeff_lj;
 380
 381         /* We send this, uncommon, message blocking to simplify the code */
 382         MPI_Send(&cnb, sizeof(cnb), MPI_BYTE,
 383                  cr->dd->pme_nodeid, eCommType_CNB, cr->mpi_comm_mysim);
 384     }
 385 #endif
 386 }
 387
 388 void gmx_pme_send_resetcounters(t_commrec gmx_unused *cr, gmx_int64_t gmx_unused step)
 389 {
 390 #ifdef GMX_MPI
 391     gmx_pme_comm_n_box_t cnb;
 392
 393     /* Only let one PP node signal each PME node */
 394     if (cr->dd->pme_receive_vir_ener)
 395     {
 396         cnb.flags = PP_PME_RESETCOUNTERS;
 397         cnb.step  = step;
 398
 399         /* We send this, uncommon, message blocking to simplify the code */
 400         MPI_Send(&cnb, sizeof(cnb), MPI_BYTE,
 401                  cr->dd->pme_nodeid, eCommType_CNB, cr->mpi_comm_mysim);
 402     }
 403 #endif
 404 }
 405
 406 int gmx_pme_recv_coeffs_coords(struct gmx_pme_pp          *pme_pp,
 407                                int                        *natoms,
 408                                real                      **chargeA,
 409                                real                      **chargeB,
 410                                real                      **sqrt_c6A,
 411                                real                      **sqrt_c6B,
 412                                real                      **sigmaA,
 413                                real                      **sigmaB,
 414                                matrix gmx_unused           box,
 415                                rvec                      **x,
 416                                rvec                      **f,
 417                                int gmx_unused             *maxshift_x,
 418                                int gmx_unused             *maxshift_y,
 419                                gmx_bool gmx_unused        *bFreeEnergy_q,
 420                                gmx_bool gmx_unused        *bFreeEnergy_lj,
 421                                real gmx_unused            *lambda_q,
 422                                real gmx_unused            *lambda_lj,
 423                                gmx_bool gmx_unused        *bEnerVir,
 424                                int                        *pme_flags,
 425                                gmx_int64_t gmx_unused     *step,
 426                                ivec gmx_unused             grid_size,
 427                                real gmx_unused            *ewaldcoeff_q,
 428                                real gmx_unused            *ewaldcoeff_lj)
 429 {
 430     gmx_pme_comm_n_box_t cnb;
 431     int                  nat = 0, q, messages, sender;
 432     real                *charge_pp;
 433
 434     messages = 0;
 435
 436     /* avoid compiler warning about unused variable without MPI support */
 437     cnb.flags  = 0;
 438     *pme_flags = 0;
 439 #ifdef GMX_MPI
 440     do
 441     {
 442         /* Receive the send count, box and time step from the peer PP node */
 443         MPI_Recv(&cnb, sizeof(cnb), MPI_BYTE,
 444                  pme_pp->node_peer, eCommType_CNB,
 445                  pme_pp->mpi_comm_mysim, MPI_STATUS_IGNORE);
 446
 447         if (debug)
 448         {
 449             fprintf(debug, "PME only rank receiving:%s%s%s%s%s\n",
 450                     (cnb.flags & PP_PME_CHARGE)        ? " charges" : "",
 451                     (cnb.flags & PP_PME_COORD )        ? " coordinates" : "",
 452                     (cnb.flags & PP_PME_FINISH)        ? " finish" : "",
 453                     (cnb.flags & PP_PME_SWITCHGRID)    ? " switch grid" : "",
 454                     (cnb.flags & PP_PME_RESETCOUNTERS) ? " reset counters" : "");
 455         }
 456
 457         if (cnb.flags & PP_PME_SWITCHGRID)
 458         {
 459             /* Special case, receive the new parameters and return */
 460             copy_ivec(cnb.grid_size, grid_size);
 461             *ewaldcoeff_q  = cnb.ewaldcoeff_q;
 462             *ewaldcoeff_lj = cnb.ewaldcoeff_lj;
 463             return pmerecvqxSWITCHGRID;
 464         }
 465
 466         if (cnb.flags & PP_PME_RESETCOUNTERS)
 467         {
 468             /* Special case, receive the step and return */
 469             *step = cnb.step;
 470
 471             return pmerecvqxRESETCOUNTERS;
 472         }
 473
 474         if (cnb.flags & (PP_PME_CHARGE | PP_PME_SQRTC6 | PP_PME_SIGMA))
 475         {
 476             /* Receive the send counts from the other PP nodes */
 477             for (sender = 0; sender < pme_pp->nnode; sender++)
 478             {
 479                 if (pme_pp->node[sender] == pme_pp->node_peer)
 480                 {
 481                     pme_pp->nat[sender] = cnb.natoms;
 482                 }
 483                 else
 484                 {
 485                     MPI_Irecv(&(pme_pp->nat[sender]), sizeof(pme_pp->nat[0]),
 486                               MPI_BYTE,
 487                               pme_pp->node[sender], eCommType_CNB,
 488                               pme_pp->mpi_comm_mysim, &pme_pp->req[messages++]);
 489                 }
 490             }
 491             MPI_Waitall(messages, pme_pp->req, pme_pp->stat);
 492             messages = 0;
 493
 494             nat = 0;
 495             for (sender = 0; sender < pme_pp->nnode; sender++)
 496             {
 497                 nat += pme_pp->nat[sender];
 498             }
 499
 500             if (nat > pme_pp->nalloc)
 501             {
 502                 pme_pp->nalloc = over_alloc_dd(nat);
 503                 if (cnb.flags & PP_PME_CHARGE)
 504                 {
 505                     srenew(pme_pp->chargeA, pme_pp->nalloc);
 506                 }
 507                 if (cnb.flags & PP_PME_CHARGEB)
 508                 {
 509                     srenew(pme_pp->chargeB, pme_pp->nalloc);
 510                 }
 511                 if (cnb.flags & PP_PME_SQRTC6)
 512                 {
 513                     srenew(pme_pp->sqrt_c6A, pme_pp->nalloc);
 514                 }
 515                 if (cnb.flags & PP_PME_SQRTC6B)
 516                 {
 517                     srenew(pme_pp->sqrt_c6B, pme_pp->nalloc);
 518                 }
 519                 if (cnb.flags & PP_PME_SIGMA)
 520                 {
 521                     srenew(pme_pp->sigmaA, pme_pp->nalloc);
 522                 }
 523                 if (cnb.flags & PP_PME_SIGMAB)
 524                 {
 525                     srenew(pme_pp->sigmaB, pme_pp->nalloc);
 526                 }
 527                 srenew(pme_pp->x, pme_pp->nalloc);
 528                 srenew(pme_pp->f, pme_pp->nalloc);
 529             }
 530
 531             /* maxshift is sent when the charges are sent */
 532             *maxshift_x = cnb.maxshift_x;
 533             *maxshift_y = cnb.maxshift_y;
 534
 535             /* Receive the charges in place */
 536             for (q = 0; q < eCommType_NR; q++)
 537             {
 538                 if (!(cnb.flags & (PP_PME_CHARGE<<q)))
 539                 {
 540                     continue;
 541                 }
 542                 switch (q)
 543                 {
 544                     case eCommType_ChargeA: charge_pp = pme_pp->chargeA;  break;
 545                     case eCommType_ChargeB: charge_pp = pme_pp->chargeB;  break;
 546                     case eCommType_SQRTC6A: charge_pp = pme_pp->sqrt_c6A; break;
 547                     case eCommType_SQRTC6B: charge_pp = pme_pp->sqrt_c6B; break;
 548                     case eCommType_SigmaA:  charge_pp = pme_pp->sigmaA;   break;
 549                     case eCommType_SigmaB:  charge_pp = pme_pp->sigmaB;   break;
 550                     default: gmx_incons("Wrong eCommType");
 551                 }
 552                 nat = 0;
 553                 for (sender = 0; sender < pme_pp->nnode; sender++)
 554                 {
 555                     if (pme_pp->nat[sender] > 0)
 556                     {
 557                         MPI_Irecv(charge_pp+nat,
 558                                   pme_pp->nat[sender]*sizeof(real),
 559                                   MPI_BYTE,
 560                                   pme_pp->node[sender], q,
 561                                   pme_pp->mpi_comm_mysim,
 562                                   &pme_pp->req[messages++]);
 563                         nat += pme_pp->nat[sender];
 564                         if (debug)
 565                         {
 566                             fprintf(debug, "Received from PP rank %d: %d %s\n",
 567                                     pme_pp->node[sender], pme_pp->nat[sender],
 568                                     (q == eCommType_ChargeA ||
 569                                      q == eCommType_ChargeB) ? "charges" : "params");
 570                         }
 571                     }
 572                 }
 573             }
 574
 575             pme_pp->flags_charge = cnb.flags;
 576         }
 577
 578         if (cnb.flags & PP_PME_COORD)
 579         {
 580             if (!(pme_pp->flags_charge & (PP_PME_CHARGE | PP_PME_SQRTC6)))
 581             {
 582                 gmx_incons("PME-only rank received coordinates before charges and/or C6-values"
 583                            );
 584             }
 585
 586             /* The box, FE flag and lambda are sent along with the coordinates
 587              *  */
 588             copy_mat(cnb.box, box);
 589             *bFreeEnergy_q  = ((cnb.flags & GMX_PME_DO_COULOMB) &&
 590                                (cnb.flags & PP_PME_FEP_Q));
 591             *bFreeEnergy_lj = ((cnb.flags & GMX_PME_DO_LJ) &&
 592                                (cnb.flags & PP_PME_FEP_LJ));
 593             *lambda_q       = cnb.lambda_q;
 594             *lambda_lj      = cnb.lambda_lj;
 595             *bEnerVir       = (cnb.flags & PP_PME_ENER_VIR);
 596             *pme_flags      = cnb.flags;
 597
 598             if (*bFreeEnergy_q && !(pme_pp->flags_charge & PP_PME_CHARGEB))
 599             {
 600                 gmx_incons("PME-only rank received free energy request, but "
 601                            "did not receive B-state charges");
 602             }
 603
 604             if (*bFreeEnergy_lj && !(pme_pp->flags_charge & PP_PME_SQRTC6B))
 605             {
 606                 gmx_incons("PME-only rank received free energy request, but "
 607                            "did not receive B-state C6-values");
 608             }
 609
 610             /* Receive the coordinates in place */
 611             nat = 0;
 612             for (sender = 0; sender < pme_pp->nnode; sender++)
 613             {
 614                 if (pme_pp->nat[sender] > 0)
 615                 {
 616                     MPI_Irecv(pme_pp->x[nat], pme_pp->nat[sender]*sizeof(rvec),
 617                               MPI_BYTE,
 618                               pme_pp->node[sender], eCommType_COORD,
 619                               pme_pp->mpi_comm_mysim, &pme_pp->req[messages++]);
 620                     nat += pme_pp->nat[sender];
 621                     if (debug)
 622                     {
 623                         fprintf(debug, "Received from PP rank %d: %d "
 624                                 "coordinates\n",
 625                                 pme_pp->node[sender], pme_pp->nat[sender]);
 626                     }
 627                 }
 628             }
 629         }
 630
 631         /* Wait for the coordinates and/or charges to arrive */
 632         MPI_Waitall(messages, pme_pp->req, pme_pp->stat);
 633         messages = 0;
 634     }
 635     while (!(cnb.flags & (PP_PME_COORD | PP_PME_FINISH)));
 636
 637     *step = cnb.step;
 638 #endif
 639
 640     *natoms   = nat;
 641     *chargeA  = pme_pp->chargeA;
 642     *chargeB  = pme_pp->chargeB;
 643     *sqrt_c6A = pme_pp->sqrt_c6A;
 644     *sqrt_c6B = pme_pp->sqrt_c6B;
 645     *sigmaA   = pme_pp->sigmaA;
 646     *sigmaB   = pme_pp->sigmaB;
 647     *x        = pme_pp->x;
 648     *f        = pme_pp->f;
 649
 650     return ((cnb.flags & PP_PME_FINISH) ? pmerecvqxFINISH : pmerecvqxX);
 651 }
 652
 653 static void receive_virial_energy(t_commrec *cr,
 654                                   matrix vir_q, real *energy_q,
 655                                   matrix vir_lj, real *energy_lj,
 656                                   real *dvdlambda_q, real *dvdlambda_lj,
 657                                   float *pme_cycles)
 658 {
 659     gmx_pme_comm_vir_ene_t cve;
 660
 661     if (cr->dd->pme_receive_vir_ener)
 662     {
 663         if (debug)
 664         {
 665             fprintf(debug,
 666                     "PP rank %d receiving from PME rank %d: virial and energy\n",
 667                     cr->sim_nodeid, cr->dd->pme_nodeid);
 668         }
 669 #ifdef GMX_MPI
 670         MPI_Recv(&cve, sizeof(cve), MPI_BYTE, cr->dd->pme_nodeid, 1, cr->mpi_comm_mysim,
 671                  MPI_STATUS_IGNORE);
 672 #else
 673         memset(&cve, 0, sizeof(cve));
 674 #endif
 675
 676         m_add(vir_q, cve.vir_q, vir_q);
 677         m_add(vir_lj, cve.vir_lj, vir_lj);
 678         *energy_q      = cve.energy_q;
 679         *energy_lj     = cve.energy_lj;
 680         *dvdlambda_q  += cve.dvdlambda_q;
 681         *dvdlambda_lj += cve.dvdlambda_lj;
 682         *pme_cycles    = cve.cycles;
 683
 684         if (cve.stop_cond != gmx_stop_cond_none)
 685         {
 686             gmx_set_stop_condition(cve.stop_cond);
 687         }
 688     }
 689     else
 690     {
 691         *energy_q   = 0;
 692         *energy_lj  = 0;
 693         *pme_cycles = 0;
 694     }
 695 }
 696
 697 void gmx_pme_receive_f(t_commrec *cr,
 698                        rvec f[], matrix vir_q, real *energy_q,
 699                        matrix vir_lj, real *energy_lj,
 700                        real *dvdlambda_q, real *dvdlambda_lj,
 701                        float *pme_cycles)
 702 {
 703     int natoms, i;
 704
 705 #ifdef GMX_PME_DELAYED_WAIT
 706     /* Wait for the x request to finish */
 707     gmx_pme_send_coeffs_coords_wait(cr->dd);
 708 #endif
 709
 710     natoms = cr->dd->nat_home;
 711
 712     if (natoms > cr->dd->pme_recv_f_alloc)
 713     {
 714         cr->dd->pme_recv_f_alloc = over_alloc_dd(natoms);
 715         srenew(cr->dd->pme_recv_f_buf, cr->dd->pme_recv_f_alloc);
 716     }
 717
 718 #ifdef GMX_MPI
 719     MPI_Recv(cr->dd->pme_recv_f_buf[0],
 720              natoms*sizeof(rvec), MPI_BYTE,
 721              cr->dd->pme_nodeid, 0, cr->mpi_comm_mysim,
 722              MPI_STATUS_IGNORE);
 723 #endif
 724
 725     for (i = 0; i < natoms; i++)
 726     {
 727         rvec_inc(f[i], cr->dd->pme_recv_f_buf[i]);
 728     }
 729
 730
 731     receive_virial_energy(cr, vir_q, energy_q, vir_lj, energy_lj, dvdlambda_q, dvdlambda_lj, pme_cycles);
 732 }
 733
 734 void gmx_pme_send_force_vir_ener(struct gmx_pme_pp *pme_pp,
 735                                  rvec gmx_unused *f,
 736                                  matrix vir_q, real energy_q,
 737                                  matrix vir_lj, real energy_lj,
 738                                  real dvdlambda_q, real dvdlambda_lj,
 739                                  float cycles)
 740 {
 741     gmx_pme_comm_vir_ene_t cve;
 742     int                    messages, ind_start, ind_end, receiver;
 743
 744     cve.cycles = cycles;
 745
 746     /* Now the evaluated forces have to be transferred to the PP nodes */
 747     messages = 0;
 748     ind_end  = 0;
 749     for (receiver = 0; receiver < pme_pp->nnode; receiver++)
 750     {
 751         ind_start = ind_end;
 752         ind_end   = ind_start + pme_pp->nat[receiver];
 753 #ifdef GMX_MPI
 754         if (MPI_Isend(f[ind_start], (ind_end-ind_start)*sizeof(rvec), MPI_BYTE,
 755                       pme_pp->node[receiver], 0,
 756                       pme_pp->mpi_comm_mysim, &pme_pp->req[messages++]) != 0)
 757         {
 758             gmx_comm("MPI_Isend failed in do_pmeonly");
 759         }
 760 #endif
 761     }
 762
 763     /* send virial and energy to our last PP node */
 764     copy_mat(vir_q, cve.vir_q);
 765     copy_mat(vir_lj, cve.vir_lj);
 766     cve.energy_q     = energy_q;
 767     cve.energy_lj    = energy_lj;
 768     cve.dvdlambda_q  = dvdlambda_q;
 769     cve.dvdlambda_lj = dvdlambda_lj;
 770     /* check for the signals to send back to a PP node */
 771     cve.stop_cond = gmx_get_stop_condition();
 772
 773     cve.cycles = cycles;
 774
 775     if (debug)
 776     {
 777         fprintf(debug, "PME rank sending to PP rank %d: virial and energy\n",
 778                 pme_pp->node_peer);
 779     }
 780 #ifdef GMX_MPI
 781     MPI_Isend(&cve, sizeof(cve), MPI_BYTE,
 782               pme_pp->node_peer, 1,
 783               pme_pp->mpi_comm_mysim, &pme_pp->req[messages++]);
 784
 785     /* Wait for the forces to arrive */
 786     MPI_Waitall(messages, pme_pp->req, pme_pp->stat);
 787 #endif
 788 }