#include "gromacs/domdec/domdec.h"
#include "gromacs/domdec/domdec_network.h"
#include "gromacs/domdec/ga2la.h"
-#include "gromacs/gmxlib/chargegroup.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/forcerec.h"
bool bSettle = false;
//! \brief All bonded interactions have to be assigned?
bool bBCheck = false;
- //! \brief Are there bondeds/exclusions between charge-groups?
- bool bInterCGInteractions = false;
+ //! \brief Are there bondeds/exclusions between atoms?
+ bool bInterAtomicInteractions = false;
//! \brief Reverse ilist for all moltypes
std::vector<reverse_ilist_t> ril_mt;
//! \brief The size of ril_mt[?].index summed over all entries
rt.bSettle = bSettle;
rt.bBCheck = bBCheck;
- rt.bInterCGInteractions = mtop->bIntermolecularInteractions;
+ rt.bInterAtomicInteractions = mtop->bIntermolecularInteractions;
rt.ril_mt.resize(mtop->moltype.size());
rt.ril_mt_tot_size = 0;
std::vector<int> nint_mt;
for (size_t mt = 0; mt < mtop->moltype.size(); mt++)
{
const gmx_moltype_t &molt = mtop->moltype[mt];
- if (molt.cgs.nr > 1)
+ if (molt.atoms.nr > 1)
{
- rt.bInterCGInteractions = true;
+ rt.bInterAtomicInteractions = true;
}
/* Make the atom to interaction list for this molecule type */
int nbonded_local;
gmx_reverse_top_t *rt;
- if (dd->reverse_top->bInterCGInteractions)
+ if (dd->reverse_top->bInterAtomicInteractions)
{
nzone_bondeds = zones->n;
}
bRCheckMB = FALSE;
bRCheck2B = FALSE;
- if (dd->reverse_top->bInterCGInteractions)
+ if (dd->reverse_top->bInterAtomicInteractions)
{
/* We need to check to which cell bondeds should be assigned */
rc = dd_cutoff_twobody(dd);
}
}
-/*! \brief Return a vector of the charge group index for all atoms */
-static std::vector<int> make_at2cg(const t_block &cgs)
-{
- std::vector<int> at2cg(cgs.index[cgs.nr]);
- for (int cg = 0; cg < cgs.nr; cg++)
- {
- for (int a = cgs.index[cg]; a < cgs.index[cg + 1]; a++)
- {
- at2cg[a] = cg;
- }
- }
-
- return at2cg;
-}
-
t_blocka *makeBondedLinks(const gmx_mtop_t *mtop,
cginfo_mb_t *cginfo_mb)
{
t_blocka *link;
cginfo_mb_t *cgi_mb;
- /* For each charge group make a list of other charge groups
- * in the system that a linked to it via bonded interactions
+ /* For each atom make a list of other atoms in the system
+ * that a linked to it via bonded interactions
* which are also stored in reverse_top.
*/
reverse_ilist_t ril_intermol;
if (mtop->bIntermolecularInteractions)
{
- if (ncg_mtop(mtop) < mtop->natoms)
- {
- gmx_fatal(FARGS, "The combination of intermolecular interactions, charge groups and domain decomposition is not supported. Use cutoff-scheme=Verlet (which removes the charge groups) or run without domain decomposition.");
- }
-
t_atoms atoms;
atoms.nr = mtop->natoms;
}
snew(link, 1);
- snew(link->index, ncg_mtop(mtop)+1);
+ snew(link->index, mtop->natoms + 1);
link->nalloc_a = 0;
link->a = nullptr;
continue;
}
const gmx_moltype_t &molt = mtop->moltype[molb.type];
- const t_block &cgs = molt.cgs;
- std::vector<int> a2c = make_at2cg(cgs);
/* Make a reverse ilist in which the interactions are linked
* to all atoms, not only the first atom as in gmx_reverse_top.
* The constraints are discarded here.
int mol;
for (mol = 0; mol < (mtop->bIntermolecularInteractions ? molb.nmol : 1); mol++)
{
- for (int cg = 0; cg < cgs.nr; cg++)
+ for (int a = 0; a < molt.atoms.nr; a++)
{
- int cg_gl = cg_offset + cg;
+ int cg_gl = cg_offset + a;
link->index[cg_gl+1] = link->index[cg_gl];
- for (int a = cgs.index[cg]; a < cgs.index[cg + 1]; a++)
+ int i = ril.index[a];
+ while (i < ril.index[a+1])
{
- int i = ril.index[a];
- while (i < ril.index[a+1])
+ int ftype = ril.il[i++];
+ int nral = NRAL(ftype);
+ /* Skip the ifunc index */
+ i++;
+ for (int j = 0; j < nral; j++)
{
- int ftype = ril.il[i++];
- int nral = NRAL(ftype);
- /* Skip the ifunc index */
- i++;
- for (int j = 0; j < nral; j++)
+ int aj = ril.il[i + j];
+ if (aj != a)
{
- int aj = ril.il[i + j];
- if (a2c[aj] != cg)
- {
- check_link(link, cg_gl, cg_offset+a2c[aj]);
- }
+ check_link(link, cg_gl, cg_offset + aj);
}
- i += nral_rt(ftype);
}
+ i += nral_rt(ftype);
+ }
- if (mtop->bIntermolecularInteractions)
+ if (mtop->bIntermolecularInteractions)
+ {
+ int i = ril_intermol.index[a];
+ while (i < ril_intermol.index[a+1])
{
- int i = ril_intermol.index[a];
- while (i < ril_intermol.index[a+1])
+ int ftype = ril_intermol.il[i++];
+ int nral = NRAL(ftype);
+ /* Skip the ifunc index */
+ i++;
+ for (int j = 0; j < nral; j++)
{
- int ftype = ril_intermol.il[i++];
- int nral = NRAL(ftype);
- /* Skip the ifunc index */
- i++;
- for (int j = 0; j < nral; j++)
- {
- /* Here we assume we have no charge groups;
- * this has been checked above.
- */
- int aj = ril_intermol.il[i + j];
- check_link(link, cg_gl, aj);
- }
- i += nral_rt(ftype);
+ /* Here we assume we have no charge groups;
+ * this has been checked above.
+ */
+ int aj = ril_intermol.il[i + j];
+ check_link(link, cg_gl, aj);
}
+ i += nral_rt(ftype);
}
}
if (link->index[cg_gl+1] - link->index[cg_gl] > 0)
{
- SET_CGINFO_BOND_INTER(cgi_mb->cginfo[cg]);
+ SET_CGINFO_BOND_INTER(cgi_mb->cginfo[a]);
ncgi++;
}
}
- cg_offset += cgs.nr;
+ cg_offset += molt.atoms.nr;
}
- int nlink_mol = link->index[cg_offset] - link->index[cg_offset - cgs.nr];
+ int nlink_mol = link->index[cg_offset] - link->index[cg_offset - molt.atoms.nr];
if (debug)
{
- fprintf(debug, "molecule type '%s' %d cgs has %d cg links through bonded interac.\n", *molt.name, cgs.nr, nlink_mol);
+ fprintf(debug, "molecule type '%s' %d atoms has %d atom links through bonded interac.\n", *molt.name, molt.atoms.nr, nlink_mol);
}
if (molb.nmol > mol)
srenew(link->a, link->nalloc_a);
for (; mol < molb.nmol; mol++)
{
- for (int cg = 0; cg < cgs.nr; cg++)
+ for (int a = 0; a < molt.atoms.nr; a++)
{
- int cg_gl = cg_offset + cg;
+ int cg_gl = cg_offset + a;
link->index[cg_gl + 1] =
- link->index[cg_gl + 1 - cgs.nr] + nlink_mol;
+ link->index[cg_gl + 1 - molt.atoms.nr] + nlink_mol;
for (int j = link->index[cg_gl]; j < link->index[cg_gl+1]; j++)
{
- link->a[j] = link->a[j - nlink_mol] + cgs.nr;
+ link->a[j] = link->a[j - nlink_mol] + molt.atoms.nr;
}
if (link->index[cg_gl+1] - link->index[cg_gl] > 0 &&
cg_gl - cgi_mb->cg_start < cgi_mb->cg_mod)
ncgi++;
}
}
- cg_offset += cgs.nr;
+ cg_offset += molt.atoms.nr;
}
}
}
if (debug)
{
- fprintf(debug, "Of the %d charge groups %d are linked via bonded interactions\n", ncg_mtop(mtop), ncgi);
+ fprintf(debug, "Of the %d atoms %d are linked via bonded interactions\n",
+ mtop->natoms, ncgi);
}
return link;
/*! \brief Set the distance, function type and atom indices for the longest distance between charge-groups of molecule type \p molt for two-body and multi-body bonded interactions */
static void bonded_cg_distance_mol(const gmx_moltype_t *molt,
- const std::vector<int> &at2cg,
gmx_bool bBCheck, gmx_bool bExcl, rvec *cg_cm,
bonded_distance_t *bd_2b,
bonded_distance_t *bd_mb)
{
for (int ai = 0; ai < nral; ai++)
{
- int cgi = at2cg[il.iatoms[i+1+ai]];
+ int atomI = il.iatoms[i + 1 + ai];
for (int aj = ai + 1; aj < nral; aj++)
{
- int cgj = at2cg[il.iatoms[i+1+aj]];
- if (cgi != cgj)
+ int atomJ = il.iatoms[i + 1 + aj];
+ if (atomI != atomJ)
{
- real rij2 = distance2(cg_cm[cgi], cg_cm[cgj]);
+ real rij2 = distance2(cg_cm[atomI], cg_cm[atomJ]);
update_max_bonded_distance(rij2, ftype,
- il.iatoms[i+1+ai],
- il.iatoms[i+1+aj],
+ atomI, atomJ,
(nral == 2) ? bd_2b : bd_mb);
}
}
const t_blocka *excls = &molt->excls;
for (int ai = 0; ai < excls->nr; ai++)
{
- int cgi = at2cg[ai];
for (int j = excls->index[ai]; j < excls->index[ai+1]; j++)
{
- int cgj = at2cg[excls->a[j]];
- if (cgi != cgj)
+ int aj = excls->a[j];
+ if (ai != aj)
{
- real rij2 = distance2(cg_cm[cgi], cg_cm[cgj]);
+ real rij2 = distance2(cg_cm[ai], cg_cm[aj]);
/* There is no function type for exclusions, use -1 */
- update_max_bonded_distance(rij2, -1, ai, excls->a[j], bd_2b);
+ update_max_bonded_distance(rij2, -1, ai, aj, bd_2b);
}
}
}
return hasVsite;
}
-//! Compute charge group centers of mass for molecule \p molt
-static void get_cgcm_mol(const gmx_moltype_t *molt,
- const gmx_ffparams_t *ffparams,
- int ePBC, t_graph *graph, const matrix box,
- const rvec *x, rvec *xs, rvec *cg_cm)
+//! Returns coordinates not broken over PBC for a molecule
+static void getWholeMoleculeCoordinates(const gmx_moltype_t *molt,
+ const gmx_ffparams_t *ffparams,
+ int ePBC, t_graph *graph, const matrix box,
+ const rvec *x, rvec *xs)
{
int n, i;
* unchanged, just to be 100% sure that we do not affect
* binary reproducibility of simulations.
*/
- n = molt->cgs.index[molt->cgs.nr];
+ n = molt->atoms.nr;
for (i = 0; i < n; i++)
{
copy_rvec(x[i], xs[i]);
ffparams->iparams.data(), ilist,
epbcNONE, TRUE, nullptr, nullptr);
}
-
- calc_cgcm(nullptr, 0, molt->cgs.nr, &molt->cgs, xs, cg_cm);
}
void dd_bonded_cg_distance(const gmx::MDLogger &mdlog,
gmx_bool bExclRequired;
int at_offset;
t_graph graph;
- rvec *xs, *cg_cm;
+ rvec *xs;
bonded_distance_t bd_2b = { 0, -1, -1, -1 };
bonded_distance_t bd_mb = { 0, -1, -1, -1 };
for (const gmx_molblock_t &molb : mtop->molblock)
{
const gmx_moltype_t &molt = mtop->moltype[molb.type];
- if (molt.cgs.nr == 1 || molb.nmol == 0)
+ if (molt.atoms.nr == 1 || molb.nmol == 0)
{
at_offset += molb.nmol*molt.atoms.nr;
}
mk_graph_moltype(molt, &graph);
}
- std::vector<int> at2cg = make_at2cg(molt.cgs);
snew(xs, molt.atoms.nr);
- snew(cg_cm, molt.cgs.nr);
for (int mol = 0; mol < molb.nmol; mol++)
{
- get_cgcm_mol(&molt, &mtop->ffparams, ir->ePBC, &graph, box,
- x+at_offset, xs, cg_cm);
+ getWholeMoleculeCoordinates(&molt, &mtop->ffparams, ir->ePBC, &graph, box,
+ x+at_offset, xs);
bonded_distance_t bd_mol_2b = { 0, -1, -1, -1 };
bonded_distance_t bd_mol_mb = { 0, -1, -1, -1 };
- bonded_cg_distance_mol(&molt, at2cg, bBCheck, bExclRequired, cg_cm,
+ bonded_cg_distance_mol(&molt, bBCheck, bExclRequired, xs,
&bd_mol_2b, &bd_mol_mb);
/* Process the mol data adding the atom index offset */
at_offset += molt.atoms.nr;
}
- sfree(cg_cm);
sfree(xs);
if (ir->ePBC != epbcNONE)
{
if (mtop->bIntermolecularInteractions)
{
- if (ncg_mtop(mtop) < mtop->natoms)
- {
- gmx_fatal(FARGS, "The combination of intermolecular interactions, charge groups and domain decomposition is not supported. Use cutoff-scheme=Verlet (which removes the charge groups) or run without domain decomposition.");
- }
-
GMX_RELEASE_ASSERT(mtop->intermolecular_ilist, "We should have an ilist when intermolecular interactions are on");
bonded_distance_intermol(*mtop->intermolecular_ilist,