Gromacs  2026.0-dev-20241106-9ba7f4d
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Functions
#include "gmxpre.h"
#include "gromacs/domdec/localtopology.h"
#include <cstdio>
#include <algorithm>
#include <array>
#include <iterator>
#include <memory>
#include <vector>
#include "gromacs/domdec/domdec_internal.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/domdec/domdec_zones.h"
#include "gromacs/domdec/ga2la.h"
#include "gromacs/domdec/options.h"
#include "gromacs/domdec/reversetopology.h"
#include "gromacs/math/functions.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdtypes/atominfo.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/mdatom.h"
#include "gromacs/pbcutil/pbc.h"
#include "gromacs/topology/idef.h"
#include "gromacs/topology/ifunc.h"
#include "gromacs/topology/mtop_util.h"
#include "gromacs/topology/topology.h"
#include "gromacs/topology/topsort.h"
#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/fixedcapacityvector.h"
#include "gromacs/utility/gmxassert.h"
#include "gromacs/utility/listoflists.h"
#include "gromacs/utility/range.h"
#include "gromacs/utility/real.h"
#include "gromacs/utility/strconvert.h"
+ Include dependency graph for localtopology.cpp:

Description

This file defines functions used by the domdec module while building topologies local to a domain.

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se

Classes

struct  AtomIndexSet
 Container for returning molecule type and index information for an atom. More...
 
struct  AtomInMolblock
 Container for returning molecule type and index information for an atom. More...
 

Functions

static AtomInMolblock atomInMolblockFromGlobalAtomnr (ArrayRef< const MolblockIndices > molblockIndices, const int globalAtomIndex)
 Return global topology molecule information for global atom index i_gl.
 
static ArrayRef< const int > add_ifunc_for_vsites (const gmx_ga2la_t &ga2la, const int nral, const bool isLocalVsite, const AtomIndexSet &atomIndexSet, ArrayRef< const int > iatoms, InteractionList *il)
 Store a vsite at the end of il, returns an array with parameter and atom indices. More...
 
static void add_posres (int mol, int a_mol, int numAtomsInMolecule, const gmx_molblock_t &molb, gmx::ArrayRef< int > iatoms, const t_iparams *ip_in, InteractionDefinitions *idef)
 Store a position restraint in idef and iatoms, complex because the parameters are different for each entry.
 
static void add_fbposres (int mol, int a_mol, int numAtomsInMolecule, const gmx_molblock_t &molb, gmx::ArrayRef< int > iatoms, const t_iparams *ip_in, InteractionDefinitions *idef)
 Store a flat-bottomed position restraint in idef and iatoms, complex because the parameters are different for each entry.
 
static void add_vsite (const gmx_ga2la_t &ga2la, const reverse_ilist_t &reverseIlist, const int ftype, const int nral, const bool isLocalVsite, const AtomIndexSet &atomIndexSet, ArrayRef< const int > iatoms, InteractionDefinitions *idef)
 Store a virtual site interaction, complex because of PBC and recursion.
 
static real dd_dist2 (const t_pbc *pbc_null, ArrayRef< const RVec > coordinates, const int i, const int j)
 Returns the squared distance between atoms i and j.
 
static void combine_idef (InteractionDefinitions *dest, gmx::ArrayRef< const thread_work_t > src)
 Append t_idef structures 1 to nsrc in src to *dest.
 
template<bool haveSingleDomain>
static int assignInteractionsForAtom (const AtomIndexSet &atomIndexSet, const reverse_ilist_t &reverseIlist, const gmx_ga2la_t &ga2la, const gmx::DomdecZones &zones, const bool gmx_unused checkDistanceMultiBody, const ivec gmx_unused rcheck, const bool gmx_unused checkDistanceTwoBody, const real gmx_unused cutoffSquared, const t_pbc gmx_unused *pbc_null, ArrayRef< const RVec > gmx_unused coordinates, InteractionDefinitions *idef, const int iz, const DDBondedChecking ddBondedChecking)
 Determine whether the local domain has responsibility for any of the bonded interactions for local atom atomIndex and assign those to the local domain. More...
 
static int assignPositionRestraintsForAtom (const AtomIndexSet &atomIndexSet, const int moleculeIndex, const int numAtomsInMolecule, const reverse_ilist_t &reverseIlist, const gmx_molblock_t &molb, const t_iparams *ip_in, InteractionDefinitions *idef)
 Determine whether the local domain has responsibility for any of the position restraints for local atom atomIndex and assign those to the local domain. More...
 
template<bool haveSingleDomain>
static int make_bondeds_zone (const gmx_reverse_top_t &rt, ArrayRef< const int > globalAtomIndices, const gmx_ga2la_t &ga2la, const gmx::DomdecZones &zones, const std::vector< gmx_molblock_t > &molb, const bool checkDistanceMultiBody, const ivec rcheck, const bool checkDistanceTwoBody, const real cutoffSquared, const t_pbc *pbc_null, ArrayRef< const RVec > coordinates, const t_iparams *ip_in, InteractionDefinitions *idef, int izone, const gmx::Range< int > &atomRange)
 This function looks up and assigns bonded interactions for zone iz. More...
 
template<bool haveSingleDomain>
static void make_exclusions_zone (ArrayRef< const int > globalAtomIndices, const gmx_ga2la_t &ga2la, const gmx::DomdecZones &zones, ArrayRef< const MolblockIndices > molblockIndices, const std::vector< gmx_moltype_t > &moltype, gmx::ArrayRef< const int32_t > atomInfo, ListOfLists< int > *lexcls, int iz, int at_start, int at_end, const gmx::ArrayRef< const int > intermolecularExclusionGroup)
 Set the exclusion data for i-zone iz.
 
static int make_local_bondeds_excls (const gmx_domdec_t &dd, const gmx::DomdecZones &zones, const gmx_mtop_t &mtop, ArrayRef< const int32_t > atomInfo, const bool checkDistanceMultiBody, const ivec rcheck, const gmx_bool checkDistanceTwoBody, const real cutoff, const t_pbc *pbc_null, ArrayRef< const RVec > coordinates, InteractionDefinitions *idef, ListOfLists< int > *lexcls)
 Generate and store all required local bonded interactions in idef and local exclusions in lexcls. More...
 
int dd_make_local_top (const gmx_domdec_t &dd, const gmx::DomdecZones &zones, int npbcdim, matrix box, rvec cellsize_min, const ivec npulse, t_forcerec *fr, ArrayRef< const RVec > coordinates, const gmx_mtop_t &mtop, gmx::ArrayRef< const int32_t > atomInfo, gmx_localtop_t *ltop)
 

Function Documentation

static ArrayRef<const int> add_ifunc_for_vsites ( const gmx_ga2la_t ga2la,
const int  nral,
const bool  isLocalVsite,
const AtomIndexSet atomIndexSet,
ArrayRef< const int >  iatoms,
InteractionList il 
)
inlinestatic

Store a vsite at the end of il, returns an array with parameter and atom indices.

This routine is very similar to add_ifunc, but vsites interactions have more work to do than other kinds of interactions, and the complex way nral (and thus vector contents) depends on ftype confuses static analysis tools unless we fuse the vsite atom-indexing organization code with the ifunc-adding code, so that they can see that nral is the same value.

Parameters
[in]ga2laThe global to local atom index mapping
[in]nralThe number of atoms involved in this vsite
[in]isLocalVsiteWhether all atoms involved in this vsite are local
[in]atomIndexSetThe atom index set for the virtual site
[in]iatomsIndices, 0: interaction type, 1: vsite (unused), 2 ...: constructing atoms
[in,out]ilThe interaction list to add this vsite to
Returns
an array with the parameter index and the NRAL atom indices
template<bool haveSingleDomain>
static int assignInteractionsForAtom ( const AtomIndexSet atomIndexSet,
const reverse_ilist_t &  reverseIlist,
const gmx_ga2la_t ga2la,
const gmx::DomdecZones zones,
const bool gmx_unused  checkDistanceMultiBody,
const ivec gmx_unused  rcheck,
const bool gmx_unused  checkDistanceTwoBody,
const real gmx_unused  cutoffSquared,
const t_pbc gmx_unused pbc_null,
ArrayRef< const RVec > gmx_unused  coordinates,
InteractionDefinitions idef,
const int  iz,
const DDBondedChecking  ddBondedChecking 
)
inlinestatic

Determine whether the local domain has responsibility for any of the bonded interactions for local atom atomIndex and assign those to the local domain.

Returns
The total number of bonded interactions for this atom for which this domain is responsible.
static int assignPositionRestraintsForAtom ( const AtomIndexSet atomIndexSet,
const int  moleculeIndex,
const int  numAtomsInMolecule,
const reverse_ilist_t &  reverseIlist,
const gmx_molblock_t molb,
const t_iparams ip_in,
InteractionDefinitions idef 
)
inlinestatic

Determine whether the local domain has responsibility for any of the position restraints for local atom atomIndex and assign those to the local domain.

Returns
The total number of bonded interactions for this atom for which this domain is responsible.
template<bool haveSingleDomain>
static int make_bondeds_zone ( const gmx_reverse_top_t rt,
ArrayRef< const int >  globalAtomIndices,
const gmx_ga2la_t ga2la,
const gmx::DomdecZones zones,
const std::vector< gmx_molblock_t > &  molb,
const bool  checkDistanceMultiBody,
const ivec  rcheck,
const bool  checkDistanceTwoBody,
const real  cutoffSquared,
const t_pbc pbc_null,
ArrayRef< const RVec >  coordinates,
const t_iparams ip_in,
InteractionDefinitions idef,
int  izone,
const gmx::Range< int > &  atomRange 
)
static

This function looks up and assigns bonded interactions for zone iz.

With thread parallelizing each thread acts on a different atom range: at_start to at_end.

static int make_local_bondeds_excls ( const gmx_domdec_t &  dd,
const gmx::DomdecZones zones,
const gmx_mtop_t &  mtop,
ArrayRef< const int32_t >  atomInfo,
const bool  checkDistanceMultiBody,
const ivec  rcheck,
const gmx_bool  checkDistanceTwoBody,
const real  cutoff,
const t_pbc pbc_null,
ArrayRef< const RVec >  coordinates,
InteractionDefinitions idef,
ListOfLists< int > *  lexcls 
)
static

Generate and store all required local bonded interactions in idef and local exclusions in lexcls.

Returns
Total count of bonded interactions in the local topology on this domain