1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
4 * This source code is part of
8 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2008, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
36 /*! \file groupcoord.h
38 * @brief Assemble atom positions for comparison with a reference set.
40 * This file contains functions to assemble the positions of a subset of the
41 * atoms and to do operations on it like determining the center of mass, or
42 * doing translations and rotations. These functions are useful when
43 * a subset of the positions needs to be compared to some set of reference
44 * positions, as e.g. done for essential dynamics.
54 #include "types/commrec.h"
57 /*! \brief Select local atoms of a group.
59 * Selects the indices of local atoms of a group and stores them in anrs_loc[0..nr_loc].
60 * If you need the positions of the group's atoms on all nodes, provide a coll_ind[0..nr]
61 * array and pass it on to communicate_group_positions. Thus the collective array
62 * will always have the same atom order (ascending indices).
64 * \param[in] ga2la Global to local atom index conversion data.
65 * \param[in] nr The total number of atoms that the group contains.
66 * \param[in] anrs The global atom number of the group's atoms.
67 * \param[out] nr_loc The number of group atoms present on the local node.
68 * \param[out] anrs_loc The local atom numbers of the group.
69 * \param[in,out] nalloc_loc Local allocation size of anrs_loc array.
70 * \param[out] coll_ind If not NULL this array must be of size nr. It stores
71 * for each local atom where it belongs in the global
72 * (collective) array such that it can be gmx_summed
73 * in the communicate_group_positions routine.
75 extern void dd_make_local_group_indices(gmx_ga2la_t ga2la,
76 const int nr, int anrs[], int *nr_loc,
77 int *anrs_loc[], int *nalloc_loc,
81 /*! \brief Assemble local positions into a collective array present on all nodes.
83 * Communicate the positions of the group's atoms such that every node has all of
84 * them. Unless running on huge number of cores, this is not a big performance impact
85 * as long as the collective subset [0..nr] is kept small. The atom indices are
86 * retrieved from anrs_loc[0..nr_loc]. If you call the routine for the serial case,
87 * provide an array coll_ind[i] = i for i in 1..nr.
89 * \param[in] cr Pointer to MPI communication data.
90 * \param[out] xcoll Collective array of positions, idential on all nodes
91 * after this routine has been called.
92 * \param[in,out] shifts Collective array of shifts for xcoll, needed to make
93 * the group whole. This array remembers the shifts
94 * since the start of the simulation (where the group
95 * is whole) and must therefore not be changed outside
97 * \param[out] extra_shifts Extra shifts since last time step, only needed as
98 * buffer variable [0..nr].
99 * \param[in] bNS Neighborsearching/domain redecomposition has been
100 * performed at the begin of this time step such that
101 * the shifts have changed and need to be updated.
102 * \param[in] x_loc Pointer to the local atom positions this node has.
103 * \param[in] nr Total number of atoms in the group.
104 * \param[in] nr_loc Number of group atoms on the local node.
105 * \param[in] anrs_loc Array of the local atom indices.
106 * \param[in] coll_ind This array of size nr stores for each local atom where
107 * it belongs in the collective array so that the local
108 * contributions can be gmx_summed. It is provided by
109 * dd_make_local_group_indices.
110 * \param[in,out] xcoll_old Positions from the last time step, used to make the
112 * \param[in] box Simulation box matrix, needed to shift xcoll such that
113 * the group becomes whole.
115 extern void communicate_group_positions(t_commrec *cr, rvec *xcoll, ivec *shifts,
116 ivec *extra_shifts, const gmx_bool bNS,
117 rvec *x_loc, const int nr, const int nr_loc,
118 int *anrs_loc, int *coll_ind, rvec *xcoll_old,
122 /*! \brief Calculates the center of the positions x locally.
124 * Calculates the center of mass (if masses are given in the weight array) or
125 * the geometrical center (if NULL is passed as weight).
127 * \param[in] x Positions.
128 * \param[in] weight Can be NULL or an array of weights. If masses are
129 * given as weights, the COM is calculated.
130 * \param[in] nr Number of positions and weights if present.
131 * \param[out] center The (weighted) center of the positions.
134 extern void get_center(rvec x[], real weight[], const int nr, rvec center);
137 /*! \brief Calculates the sum of the positions x locally.
139 * Calculates the (weighted) sum of position vectors and returns the sum of
140 * weights, which is needed when local contributions shall be summed to a
141 * global weighted center.
143 * \param[in] x Array of positions.
144 * \param[in] weight Can be NULL or an array of weights.
145 * \param[in] nr Number of positions and weights if present.
146 * \param[out] dsumvec The (weighted) sum of the positions.
147 * \return Sum of weights.
150 extern double get_sum_of_positions(rvec x[], real weight[], const int nr, dvec dsumvec);
153 /*! \brief Calculates the global center of all local arrays x.
155 * Get the center from local positions [0..nr_loc], this involves communication.
156 * Not that the positions must already have the correct PBC representation. Use
157 * this routine if no collective coordinates are assembled from which the center
158 * could be calculated without communication.
160 * \param[in] cr Pointer to MPI communication data.
161 * \param[in] x_loc Array of local positions [0..nr_loc].
162 * \param[in] weight_loc Array of local weights, these are the masses if the
163 * center of mass is to be calculated.
164 * \param[in] nr_loc The number of positions on the local node.
165 * \param[in] nr_group The number of positions in the whole group. Since
166 * this is known anyway, we do not need to communicate
167 * and sum nr_loc if we pass it over.
168 * \param[out] center The (weighted) center of all x_loc from all the
171 extern void get_center_comm(t_commrec *cr, rvec x_loc[], real weight_loc[],
172 int nr_loc, int nr_group, rvec center);
175 /*! \brief Translate positions.
177 * Add a translation vector to the positions x.
179 * \param[in,out] x Array of positions.
180 * \param[in] nr Number of entries in the position array.
181 * \param[in] transvec Translation vector to be added to all positions.
184 extern void translate_x(rvec x[], const int nr, const rvec transvec);
187 /*! \brief Rotate positions.
189 * Rotate the positions with the rotation matrix.
191 * \param[in,out] x Array of positions.
192 * \param[in] nr Number of entries in the position array.
193 * \param[in] rmat Rotation matrix to operate on all positions.
196 extern void rotate_x(rvec x[], const int nr, matrix rmat);