Merge branch release-2016
[alexxy/gromacs.git] / src / gromacs / pulling / pull.cpp
index b77a7fcacf8d2357d3b7d1fd132f4408a7644692..a944b712ebe1d87e583d07cc31ce232526c1caae 100644 (file)
@@ -66,7 +66,7 @@
 #include "gromacs/mdtypes/mdatom.h"
 #include "gromacs/pbcutil/pbc.h"
 #include "gromacs/pulling/pull_internal.h"
-#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/mtop_lookup.h"
 #include "gromacs/topology/topology.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/exceptions.h"
@@ -1328,7 +1328,7 @@ static void do_constraint(struct pull_t *pull, t_pbc *pbc,
         /* Accumulate the forces, in case we have multiple constraint steps */
         pcrd->f_scal += dr_tot[c]/((pull->group[pcrd->params.group[0]].invtm + pull->group[pcrd->params.group[1]].invtm)*dt*dt);
 
-        if (vir != NULL && pcrd->params.eGeom != epullgDIRPBC && bMaster)
+        if (vir != nullptr && pcrd->params.eGeom != epullgDIRPBC && bMaster)
         {
             double f_invr;
 
@@ -1366,7 +1366,7 @@ static void add_virial_coord_dr(tensor vir, const dvec dr, const dvec f)
 /* Adds the pull contribution to the virial */
 static void add_virial_coord(tensor vir, const pull_coord_work_t *pcrd)
 {
-    if (vir != NULL && pcrd->params.eGeom != epullgDIRPBC)
+    if (vir != nullptr && pcrd->params.eGeom != epullgDIRPBC)
     {
         /* Add the pull contribution for each distance vector to the virial. */
         add_virial_coord_dr(vir, pcrd->dr01, pcrd->f01);
@@ -1565,8 +1565,8 @@ void register_external_pull_potential(struct pull_t *pull,
                                       int            coord_index,
                                       const char    *provider)
 {
-    GMX_RELEASE_ASSERT(pull != NULL, "register_external_pull_potential called before init_pull");
-    GMX_RELEASE_ASSERT(provider != NULL, "register_external_pull_potential called with NULL as provider name");
+    GMX_RELEASE_ASSERT(pull != nullptr, "register_external_pull_potential called before init_pull");
+    GMX_RELEASE_ASSERT(provider != nullptr, "register_external_pull_potential called with NULL as provider name");
 
     if (coord_index < 0 || coord_index > pull->ncoord - 1)
     {
@@ -1582,7 +1582,7 @@ void register_external_pull_potential(struct pull_t *pull,
                   provider, coord_index + 1, epull_names[pcrd->params.eType], epull_names[epullEXTERNAL]);
     }
 
-    GMX_RELEASE_ASSERT(pcrd->params.externalPotentialProvider != NULL, "The external potential provider string for a pull coordinate is NULL");
+    GMX_RELEASE_ASSERT(pcrd->params.externalPotentialProvider != nullptr, "The external potential provider string for a pull coordinate is NULL");
 
     if (gmx_strcasecmp(provider, pcrd->params.externalPotentialProvider) != 0)
     {
@@ -1703,7 +1703,7 @@ real pull_potential(struct pull_t *pull, t_mdatoms *md, t_pbc *pbc,
     {
         real dVdl = 0;
 
-        pull_calc_coms(cr, pull, md, pbc, t, x, NULL);
+        pull_calc_coms(cr, pull, md, pbc, t, x, nullptr);
 
         for (int c = 0; c < pull->ncoord; c++)
         {
@@ -1715,7 +1715,7 @@ real pull_potential(struct pull_t *pull, t_mdatoms *md, t_pbc *pbc,
             }
 
             do_pull_pot_coord(pull, c, pbc, t, lambda,
-                              &V, MASTER(cr) ? vir : NULL, &dVdl);
+                              &V, MASTER(cr) ? vir : nullptr, &dVdl);
 
             /* Distribute the force over the atoms in the pulled groups */
             apply_forces_coord(pull, c, md, f);
@@ -1772,13 +1772,13 @@ static void make_local_pull_group(gmx_ga2la_t *ga2la,
             {
                 pg->nalloc_loc = over_alloc_dd(pg->nat_loc+1);
                 srenew(pg->ind_loc, pg->nalloc_loc);
-                if (pg->epgrppbc == epgrppbcCOS || pg->params.weight != NULL)
+                if (pg->epgrppbc == epgrppbcCOS || pg->params.weight != nullptr)
                 {
                     srenew(pg->weight_loc, pg->nalloc_loc);
                 }
             }
             pg->ind_loc[pg->nat_loc] = ii;
-            if (pg->params.weight != NULL)
+            if (pg->params.weight != nullptr)
             {
                 pg->weight_loc[pg->nat_loc] = pg->params.weight[i];
             }
@@ -1805,13 +1805,13 @@ void dd_make_local_pull_groups(t_commrec *cr, struct pull_t *pull, t_mdatoms *md
     }
     else
     {
-        ga2la = NULL;
+        ga2la = nullptr;
     }
 
     /* We always make the master node participate, such that it can do i/o
      * and to simplify MC type extensions people might have.
      */
-    bMustParticipate = (comm->bParticipateAll || dd == NULL || DDMASTER(dd));
+    bMustParticipate = (comm->bParticipateAll || dd == nullptr || DDMASTER(dd));
 
     for (g = 0; g < pull->ngroup; g++)
     {
@@ -1853,7 +1853,7 @@ void dd_make_local_pull_groups(t_commrec *cr, struct pull_t *pull, t_mdatoms *md
             (comm->bParticipate &&
              comm->must_count >= comm->setup_count - history_count);
 
-        if (debug && dd != NULL)
+        if (debug && dd != nullptr)
         {
             fprintf(debug, "Our DD rank (%3d) pull #atoms>0 or master: %d, will be part %d\n",
                     dd->rank, bMustParticipate, bWillParticipate);
@@ -1917,13 +1917,6 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
                                   const gmx_mtop_t *mtop,
                                   const t_inputrec *ir, real lambda)
 {
-    int                   i, ii, d, nfrozen, ndim;
-    real                  m, w, mbd;
-    double                tmass, wmass, wwmass;
-    const gmx_groups_t   *groups;
-    gmx_mtop_atomlookup_t alook;
-    t_atom               *atom;
-
     if (EI_ENERGY_MINIMIZATION(ir->eI) || ir->eI == eiBD)
     {
         /* There are no masses in the integrator.
@@ -1941,8 +1934,8 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
     {
         pg->nat_loc    = 0;
         pg->nalloc_loc = 0;
-        pg->ind_loc    = NULL;
-        pg->weight_loc = NULL;
+        pg->ind_loc    = nullptr;
+        pg->weight_loc = nullptr;
     }
     else
     {
@@ -1958,21 +1951,20 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
         }
     }
 
-    groups = &mtop->groups;
-
-    alook = gmx_mtop_atomlookup_init(mtop);
+    const gmx_groups_t *groups = &mtop->groups;
 
-    nfrozen = 0;
-    tmass   = 0;
-    wmass   = 0;
-    wwmass  = 0;
-    for (i = 0; i < pg->params.nat; i++)
+    /* Count frozen dimensions and (weighted) mass */
+    int    nfrozen = 0;
+    double tmass   = 0;
+    double wmass   = 0;
+    double wwmass  = 0;
+    int    molb    = 0;
+    for (int i = 0; i < pg->params.nat; i++)
     {
-        ii = pg->params.ind[i];
-        gmx_mtop_atomnr_to_atom(alook, ii, &atom);
+        int ii = pg->params.ind[i];
         if (bConstraint && ir->opts.nFreeze)
         {
-            for (d = 0; d < DIM; d++)
+            for (int d = 0; d < DIM; d++)
             {
                 if (pulldim_con[d] == 1 &&
                     ir->opts.nFreeze[ggrpnr(groups, egcFREEZE, ii)][d])
@@ -1981,14 +1973,17 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
                 }
             }
         }
+        const t_atom &atom = mtopGetAtomParameters(mtop, ii, &molb);
+        real          m;
         if (ir->efep == efepNO)
         {
-            m = atom->m;
+            m = atom.m;
         }
         else
         {
-            m = (1 - lambda)*atom->m + lambda*atom->mB;
+            m = (1 - lambda)*atom.m + lambda*atom.mB;
         }
+        real w;
         if (pg->params.nweight > 0)
         {
             w = pg->params.weight[i];
@@ -2006,13 +2001,14 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
         }
         else if (ir->eI == eiBD)
         {
+            real mbd;
             if (ir->bd_fric)
             {
                 mbd = ir->bd_fric*ir->delta_t;
             }
             else
             {
-                if (groups->grpnr[egcTC] == NULL)
+                if (groups->grpnr[egcTC] == nullptr)
                 {
                     mbd = ir->delta_t/ir->opts.tau_t[0];
                 }
@@ -2030,8 +2026,6 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
         wwmass += m*w*w;
     }
 
-    gmx_mtop_atomlookup_destroy(alook);
-
     if (wmass == 0)
     {
         /* We can have single atom groups with zero mass with potential pulling
@@ -2072,8 +2066,8 @@ static void init_pull_group_index(FILE *fplog, t_commrec *cr,
     }
     else
     {
-        ndim = 0;
-        for (d = 0; d < DIM; d++)
+        int ndim = 0;
+        for (int d = 0; d < DIM; d++)
         {
             ndim += pulldim_con[d]*pg->params.nat;
         }
@@ -2106,8 +2100,8 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
     /* Copy the pull parameters */
     pull->params       = *pull_params;
     /* Avoid pointer copies */
-    pull->params.group = NULL;
-    pull->params.coord = NULL;
+    pull->params.group = nullptr;
+    pull->params.coord = nullptr;
 
     pull->ncoord       = pull_params->ncoord;
     pull->ngroup       = pull_params->ngroup;
@@ -2396,7 +2390,7 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
                     }
                     else
                     {
-                        if (pgrp->params.weight != NULL)
+                        if (pgrp->params.weight != nullptr)
                         {
                             gmx_fatal(FARGS, "Pull groups can not have relative weights and cosine weighting at same time");
                         }
@@ -2453,9 +2447,9 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
 
 #if GMX_MPI
     /* Use a sub-communicator when we have more than 32 ranks */
-    comm->bParticipateAll = (cr == NULL || !DOMAINDECOMP(cr) ||
+    comm->bParticipateAll = (cr == nullptr || !DOMAINDECOMP(cr) ||
                              cr->dd->nnodes <= 32 ||
-                             getenv("GMX_PULL_PARTICIPATE_ALL") != NULL);
+                             getenv("GMX_PULL_PARTICIPATE_ALL") != nullptr);
     /* This sub-commicator is not used with comm->bParticipateAll,
      * so we can always initialize it to NULL.
      */
@@ -2469,21 +2463,21 @@ init_pull(FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir,
     comm->setup_count     = 0;
     comm->must_count      = 0;
 
-    if (!comm->bParticipateAll && fplog != NULL)
+    if (!comm->bParticipateAll && fplog != nullptr)
     {
         fprintf(fplog, "Will use a sub-communicator for pull communication\n");
     }
 
-    comm->rbuf     = NULL;
-    comm->dbuf     = NULL;
-    comm->dbuf_cyl = NULL;
+    comm->rbuf     = nullptr;
+    comm->dbuf     = nullptr;
+    comm->dbuf_cyl = nullptr;
 
     /* We still need to initialize the PBC reference coordinates */
     pull->bSetPBCatoms = TRUE;
 
     /* Only do I/O when we are doing dynamics and if we are the MASTER */
-    pull->out_x = NULL;
-    pull->out_f = NULL;
+    pull->out_x = nullptr;
+    pull->out_f = nullptr;
     if (bOutFile)
     {
         /* Check for px and pf filename collision, if we are writing