3 * This source code is part of
7 * 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-2004, 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 * Gromacs Runs On Most of All Computer Systems
35 #include "visibility.h"
42 /* Should be called after generating or reading mtop,
43 * to set some compute intesive variables to avoid
44 * N^2 operations later on.
48 gmx_mtop_finalize(gmx_mtop_t *mtop);
51 /* Returns the total number of charge groups in mtop */
54 ncg_mtop(const gmx_mtop_t *mtop);
56 /* Removes the charge groups, i.e. makes single atom charge groups, in mtop */
58 void gmx_mtop_remove_chargegroups(gmx_mtop_t *mtop);
61 /* Abstract data type for looking up atoms by global atom number */
62 typedef struct gmx_mtop_atomlookup *gmx_mtop_atomlookup_t;
64 /* Initialize atom lookup by global atom number */
67 gmx_mtop_atomlookup_init(const gmx_mtop_t *mtop);
69 /* As gmx_mtop_atomlookup_init, but optimized for atoms involved in settle */
72 gmx_mtop_atomlookup_settle_init(const gmx_mtop_t *mtop);
74 /* Destroy a gmx_mtop_atomlookup_t data structure */
77 gmx_mtop_atomlookup_destroy(gmx_mtop_atomlookup_t alook);
80 /* Returns a pointer to the t_atom struct belonging to atnr_global.
81 * This can be an expensive operation, so if possible use
82 * one of the atom loop constructs below.
86 gmx_mtop_atomnr_to_atom(const gmx_mtop_atomlookup_t alook,
91 /* Returns a pointer to the molecule interaction array ilist_mol[F_NRE]
92 * and the local atom number in the molecule belonging to atnr_global.
96 gmx_mtop_atomnr_to_ilist(const gmx_mtop_atomlookup_t alook,
98 t_ilist **ilist_mol,int *atnr_offset);
101 /* Returns the molecule block index
102 * and the molecule number in the block
103 * and the atom number offset for the atom indices in moltype
104 * belonging to atnr_global.
108 gmx_mtop_atomnr_to_molblock_ind(const gmx_mtop_atomlookup_t alook,
110 int *molb,int *molnr,int *atnr_mol);
113 /* Returns atom name, global resnr and residue name of atom atnr_global */
116 gmx_mtop_atominfo_global(const gmx_mtop_t *mtop,int atnr_global,
117 char **atomname,int *resnr,char **resname);
120 /* Abstract type for atom loop over all atoms */
121 typedef struct gmx_mtop_atomloop_all *gmx_mtop_atomloop_all_t;
123 /* Initialize an atom loop over all atoms in the system.
124 * The order of the atoms will be as in the state struct.
125 * Only use this when you really need to loop over all atoms,
126 * i.e. when you use groups which might differ per molecule,
127 * otherwise use gmx_mtop_atomloop_block.
130 gmx_mtop_atomloop_all_t
131 gmx_mtop_atomloop_all_init(const gmx_mtop_t *mtop);
133 /* Loop to the next atom.
134 * When not at the end:
135 * returns TRUE and at_global,
136 * writes the global atom number in *at_global
137 * and sets the pointer atom to the t_atom struct of that atom.
138 * When at the end, destroys aloop and returns FALSE.
140 * gmx_mtop_atomloop_all_t aloop;
141 * aloop = gmx_mtop_atomloop_all_init(mtop)
142 * while (gmx_mtop_atomloop_all_next(aloop,&at_global,&atom)) {
148 gmx_mtop_atomloop_all_next(gmx_mtop_atomloop_all_t aloop,
149 int *at_global,t_atom **atom);
151 /* Return the atomname, the residue number and residue name
152 * of the current atom in the loop.
155 gmx_mtop_atomloop_all_names(gmx_mtop_atomloop_all_t aloop,
156 char **atomname,int *resnr,char **resname);
158 /* Return the a pointer to the moltype struct of the current atom
159 * in the loop and the atom number in the molecule.
162 gmx_mtop_atomloop_all_moltype(gmx_mtop_atomloop_all_t aloop,
163 gmx_moltype_t **moltype,int *at_mol);
166 /* Abstract type for atom loop over atoms in all molecule blocks */
167 typedef struct gmx_mtop_atomloop_block *gmx_mtop_atomloop_block_t;
169 /* Initialize an atom loop over atoms in all molecule blocks the system.
172 gmx_mtop_atomloop_block_t
173 gmx_mtop_atomloop_block_init(const gmx_mtop_t *mtop);
175 /* Loop to the next atom.
176 * When not at the end:
178 * sets the pointer atom to the t_atom struct of that atom
179 * and return the number of molecules corresponding to this atom.
180 * When at the end, destroys aloop and returns FALSE.
182 * gmx_mtop_atomloop_block_t aloop;
183 * aloop = gmx_mtop_atomloop_block_init(mtop)
184 * while (gmx_mtop_atomloop_block_next(aloop,&atom,&nmol)) {
190 gmx_mtop_atomloop_block_next(gmx_mtop_atomloop_block_t aloop,
191 t_atom **atom,int *nmol);
194 /* Abstract type for ilist loop over all ilists */
195 typedef struct gmx_mtop_ilistloop *gmx_mtop_ilistloop_t;
197 /* Initialize an ilist loop over all molecule types in the system. */
200 gmx_mtop_ilistloop_init(const gmx_mtop_t *mtop);
203 /* Loop to the next molecule,
204 * When not at the end:
205 * returns TRUE and a pointer to the next array ilist_mol[F_NRE],
206 * writes the number of molecules for this ilist in *nmol.
207 * When at the end, destroys iloop and returns FALSE.
211 gmx_mtop_ilistloop_next(gmx_mtop_ilistloop_t iloop,
212 t_ilist **ilist_mol,int *nmol);
215 /* Abstract type for ilist loop over all ilists of all molecules */
216 typedef struct gmx_mtop_ilistloop_all *gmx_mtop_ilistloop_all_t;
218 /* Initialize an ilist loop over all molecule types in the system.
219 * Only use this when you really need to loop over all molecules,
220 * i.e. when you use groups which might differ per molecule,
221 * otherwise use gmx_mtop_ilistloop.
224 gmx_mtop_ilistloop_all_t
225 gmx_mtop_ilistloop_all_init(const gmx_mtop_t *mtop);
227 /* Loop to the next molecule,
228 * When not at the end:
229 * returns TRUE and a pointer to the next array ilist_mol[F_NRE],
230 * writes the atom offset which should be added to iatoms in atnr_offset.
231 * When at the end, destroys iloop and returns FALSE.
235 gmx_mtop_ilistloop_all_next(gmx_mtop_ilistloop_all_t iloop,
236 t_ilist **ilist_mol,int *atnr_offset);
239 /* Returns the total number of interactions in the system of type ftype */
242 gmx_mtop_ftype_count(const gmx_mtop_t *mtop,int ftype);
245 /* Returns a charge group index for the whole system */
248 gmx_mtop_global_cgs(const gmx_mtop_t *mtop);
251 /* Returns a single t_atoms struct for the whole system */
254 gmx_mtop_global_atoms(const gmx_mtop_t *mtop);
257 /* Make all charge groups the size of one atom.
258 * When bKeepSingleMolCG==TRUE keep charge groups for molecules
259 * that consist of a single charge group.
262 gmx_mtop_make_atomic_charge_groups(gmx_mtop_t *mtop,gmx_bool bKeepSingleMolCG);
265 /* Generate a 'local' topology for the whole system.
266 * When ir!=NULL the free energy interactions will be sorted to the end.
270 gmx_mtop_generate_local_top(const gmx_mtop_t *mtop,const t_inputrec *ir);
273 /* Converts a gmx_mtop_t struct to t_topology.
274 * All memory relating only to mtop will be freed.
278 gmx_mtop_t_to_t_topology(gmx_mtop_t *mtop);