*
* Copyright (c) 1991-2000, University of Groningen, The Netherlands.
* Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2011,2014, by the GROMACS development team, led by
+ * Copyright (c) 2011,2014,2015, by the GROMACS development team, led by
* Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
* and including many others, as listed in the AUTHORS file in the
* top-level source directory and at http://www.gromacs.org.
#include "vec.h"
#include "macros.h"
#include "names.h"
+#include "gmx_fatal.h"
/* these MUST correspond to the enum in hackblock.h */
const char *btsNames[ebtsNR] = { "bonds", "angles", "dihedrals", "impropers", "exclusions", "cmap" };
return bRet;
}
-gmx_bool rbonded_atoms_exist_in_list(t_rbonded *b, t_rbonded blist[], int nlist, int natoms)
+int
+rbonded_find_atoms_in_list(t_rbonded *b, t_rbonded blist[], int nlist, int natoms)
{
int i, k;
- gmx_bool matchFound = FALSE;
+ int foundPos = -1;
gmx_bool atomsMatch;
- for (i = 0; i < nlist && !matchFound; i++)
+ for (i = 0; i < nlist && foundPos < 0; i++)
{
atomsMatch = TRUE;
for (k = 0; k < natoms && atomsMatch; k++)
atomsMatch = atomsMatch && !strcmp(b->a[k], blist[i].a[natoms-1-k]);
}
}
- matchFound = atomsMatch;
+ if (atomsMatch)
+ {
+ foundPos = i;
+ /* If all the atoms AND all the parameters match, it is likely that
+ * the user made a copy-and-paste mistake (since it would be much cheaper
+ * to just bump the force constant 2x if you really want it twice).
+ * Since we only have the unparsed string here we can only detect
+ * EXACT matches (including identical whitespace).
+ */
+ if (!strcmp(b->s, blist[i].s))
+ {
+ gmx_warning("Duplicate line found in or between hackblock and rtp entries");
+ }
+ }
}
- return matchFound;
+ return foundPos;
}
gmx_bool merge_t_bondeds(t_rbondeds s[], t_rbondeds d[], gmx_bool bMin, gmx_bool bPlus)
{
int i, j;
gmx_bool bBondsRemoved;
+ int nbHackblockStart;
+ int index;
bBondsRemoved = FALSE;
for (i = 0; i < ebtsNR; i++)
{
if (s[i].nb > 0)
{
+ /* Record how many bonds we have in the destination when we start.
+ *
+ * If an entry is present in the hackblock (destination), we will
+ * not add the one from the main rtp, since the point is for hackblocks
+ * to overwrite it. However, if there is no hackblock entry we do
+ * allow multiple main rtp entries since some forcefield insist on that.
+ *
+ * We accomplish this by checking the position we find an entry in,
+ * rather than merely checking whether it exists at all.
+ * If that index is larger than the original (hackblock) destination
+ * size, it was added from the main rtp, and then we will allow more
+ * such entries. In contrast, if the entry found has a lower index
+ * it is a hackblock entry meant to override the main rtp, and then
+ * we don't add the main rtp one.
+ */
+ nbHackblockStart = d[i].nb;
+
/* make space */
srenew(d[i].b, d[i].nb + s[i].nb);
for (j = 0; j < s[i].nb; j++)
{
/* Check if this bonded string already exists before adding.
- * We are merging from the main rtp to the hackblocks, so this
- * will mean the hackblocks overwrite the man rtp, as intended.
+ * We are merging from the main RTP to the hackblocks, so this
+ * will mean the hackblocks overwrite the man RTP, as intended.
+ */
+ index = rbonded_find_atoms_in_list(&s[i].b[j], d[i].b, d[i].nb, btsNiatoms[i]);
+ /* - If we did not find this interaction at all, the index will be -1,
+ * and then we should definitely add it to the merged hackblock and rtp.
+ *
+ * Alternatively, if it was found, index will be >=0.
+ * - In case this index is lower than the original number of entries,
+ * it is already present as a *hackblock* entry, and those should
+ * always override whatever we have listed in the RTP. Thus, we
+ * should just keep that one and not add anything from the RTP.
+ * - Finally, if it was found, but with an index higher than
+ * the original number of entries, it comes from the RTP rather
+ * than hackblock, and then we must have added it ourselves
+ * in a previous iteration. In that case it is a matter of
+ * several entries for the same sequence of atoms, and we allow
+ * that in the RTP. In this case we should simply copy all of
+ * them, including this one.
*/
- if (!rbonded_atoms_exist_in_list(&s[i].b[j], d[i].b, d[i].nb, btsNiatoms[i]))
+ if (index < 0 || index >= nbHackblockStart)
{
if (!(bMin && contains_char(&s[i].b[j], '-'))
&& !(bPlus && contains_char(&s[i].b[j], '+')))
bBondsRemoved = TRUE;
}
}
+ else
+ {
+ /* This is the common case where a hackblock entry simply
+ * overrides the RTP, so we cannot warn here.
+ */
+ }
}
}
}