2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014,2015,2016,2017,2018, by the GROMACS development team, led by
7 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
8 * and including many others, as listed in the AUTHORS file in the
9 * top-level source directory and at http://www.gromacs.org.
11 * GROMACS is free software; you can redistribute it and/or
12 * modify it under the terms of the GNU Lesser General Public License
13 * as published by the Free Software Foundation; either version 2.1
14 * of the License, or (at your option) any later version.
16 * GROMACS is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 * Lesser General Public License for more details.
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with GROMACS; if not, see
23 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
24 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
26 * If you want to redistribute modifications to GROMACS, please
27 * consider that scientific software is very special. Version
28 * control is crucial - bugs must be traceable. We will be happy to
29 * consider code for inclusion in the official distribution, but
30 * derived work must not be called official GROMACS. Details are found
31 * in the README & COPYING files - if they are missing, get the
32 * official version at http://www.gromacs.org.
34 * To help us fund GROMACS development, we humbly ask that you cite
35 * the research papers on the package. Check out http://www.gromacs.org.
46 #include "gromacs/topology/atomprop.h"
47 #include "gromacs/topology/symtab.h"
48 #include "gromacs/utility/compare.h"
49 #include "gromacs/utility/fatalerror.h"
50 #include "gromacs/utility/smalloc.h"
51 #include "gromacs/utility/txtdump.h"
53 const char *ptype_str[eptNR+1] = {
54 "Atom", "Nucleus", "Shell", "Bond", "VSite", nullptr
57 void init_atom(t_atoms *at)
62 at->resinfo = nullptr;
63 at->atomname = nullptr;
64 at->atomtype = nullptr;
65 at->atomtypeB = nullptr;
66 at->pdbinfo = nullptr;
68 at->haveCharge = FALSE;
70 at->haveBState = FALSE;
71 at->havePdbInfo = FALSE;
74 void init_atomtypes(t_atomtypes *at)
77 at->atomnumber = nullptr;
80 void done_atom(t_atoms *at)
91 void done_and_delete_atoms(t_atoms *atoms)
97 void done_atomtypes(t_atomtypes *atype)
100 sfree(atype->atomnumber);
103 void add_t_atoms(t_atoms *atoms, int natom_extra, int nres_extra)
109 srenew(atoms->atomname, atoms->nr+natom_extra);
110 srenew(atoms->atom, atoms->nr+natom_extra);
111 if (nullptr != atoms->pdbinfo)
113 srenew(atoms->pdbinfo, atoms->nr+natom_extra);
115 if (nullptr != atoms->atomtype)
117 srenew(atoms->atomtype, atoms->nr+natom_extra);
119 if (nullptr != atoms->atomtypeB)
121 srenew(atoms->atomtypeB, atoms->nr+natom_extra);
123 for (i = atoms->nr; (i < atoms->nr+natom_extra); i++)
125 atoms->atomname[i] = nullptr;
126 memset(&atoms->atom[i], 0, sizeof(atoms->atom[i]));
127 if (nullptr != atoms->pdbinfo)
129 std::memset(&atoms->pdbinfo[i], 0, sizeof(atoms->pdbinfo[i]));
131 if (nullptr != atoms->atomtype)
133 atoms->atomtype[i] = nullptr;
135 if (nullptr != atoms->atomtypeB)
137 atoms->atomtypeB[i] = nullptr;
140 atoms->nr += natom_extra;
144 srenew(atoms->resinfo, atoms->nres+nres_extra);
145 for (i = atoms->nres; (i < atoms->nres+nres_extra); i++)
147 std::memset(&atoms->resinfo[i], 0, sizeof(atoms->resinfo[i]));
149 atoms->nres += nres_extra;
153 void init_t_atoms(t_atoms *atoms, int natoms, gmx_bool bPdbinfo)
157 snew(atoms->atomname, natoms);
158 atoms->atomtype = nullptr;
159 atoms->atomtypeB = nullptr;
160 snew(atoms->resinfo, natoms);
161 snew(atoms->atom, natoms);
162 atoms->haveMass = FALSE;
163 atoms->haveCharge = FALSE;
164 atoms->haveType = FALSE;
165 atoms->haveBState = FALSE;
166 atoms->havePdbInfo = bPdbinfo;
167 if (atoms->havePdbInfo)
169 snew(atoms->pdbinfo, natoms);
173 atoms->pdbinfo = nullptr;
177 void gmx_pdbinfo_init_default(t_pdbinfo *pdbinfo)
179 pdbinfo->type = epdbATOM;
181 pdbinfo->altloc = ' ';
182 pdbinfo->atomnm[0] = '\0';
183 pdbinfo->occup = 1.0;
185 pdbinfo->bAnisotropic = FALSE;
186 std::fill(pdbinfo->uij, pdbinfo->uij+6, 0.0);
189 t_atoms *copy_t_atoms(const t_atoms *src)
195 init_t_atoms(dst, src->nr, (nullptr != src->pdbinfo));
197 if (nullptr != src->atomname)
199 snew(dst->atomname, src->nr);
201 if (nullptr != src->atomtype)
203 snew(dst->atomtype, src->nr);
205 if (nullptr != src->atomtypeB)
207 snew(dst->atomtypeB, src->nr);
209 for (i = 0; (i < src->nr); i++)
211 dst->atom[i] = src->atom[i];
212 if (nullptr != src->pdbinfo)
214 dst->pdbinfo[i] = src->pdbinfo[i];
216 if (nullptr != src->atomname)
218 dst->atomname[i] = src->atomname[i];
220 if (nullptr != src->atomtype)
222 dst->atomtype[i] = src->atomtype[i];
224 if (nullptr != src->atomtypeB)
226 dst->atomtypeB[i] = src->atomtypeB[i];
229 dst->haveBState = src->haveBState;
230 dst->haveCharge = src->haveCharge;
231 dst->haveMass = src->haveMass;
232 dst->havePdbInfo = src->havePdbInfo;
233 dst->haveType = src->haveType;
234 dst->nres = src->nres;
235 for (i = 0; (i < src->nres); i++)
237 dst->resinfo[i] = src->resinfo[i];
242 void t_atoms_set_resinfo(t_atoms *atoms, int atom_ind, t_symtab *symtab,
243 const char *resname, int resnr, unsigned char ic,
244 int chainnum, char chainid)
248 ri = &atoms->resinfo[atoms->atom[atom_ind].resind];
249 ri->name = put_symtab(symtab, resname);
253 ri->chainnum = chainnum;
254 ri->chainid = chainid;
257 static void pr_atom(FILE *fp, int indent, const char *title, const t_atom *atom, int n)
261 if (available(fp, atom, indent, title))
263 indent = pr_title_n(fp, indent, title, n);
264 for (i = 0; i < n; i++)
266 pr_indent(fp, indent);
267 fprintf(fp, "%s[%6d]={type=%3hu, typeB=%3hu, ptype=%8s, m=%12.5e, "
268 "q=%12.5e, mB=%12.5e, qB=%12.5e, resind=%5d, atomnumber=%3d}\n",
269 title, i, atom[i].type, atom[i].typeB, ptype_str[atom[i].ptype],
270 atom[i].m, atom[i].q, atom[i].mB, atom[i].qB,
271 atom[i].resind, atom[i].atomnumber);
276 static void pr_strings2(FILE *fp, int indent, const char *title,
277 char ***nm, char ***nmB, int n, gmx_bool bShowNumbers)
281 if (available(fp, nm, indent, title))
283 indent = pr_title_n(fp, indent, title, n);
284 for (i = 0; i < n; i++)
286 pr_indent(fp, indent);
287 fprintf(fp, "%s[%d]={name=\"%s\",nameB=\"%s\"}\n",
288 title, bShowNumbers ? i : -1, *(nm[i]), *(nmB[i]));
293 static void pr_resinfo(FILE *fp, int indent, const char *title, const t_resinfo *resinfo, int n,
294 gmx_bool bShowNumbers)
298 if (available(fp, resinfo, indent, title))
300 indent = pr_title_n(fp, indent, title, n);
301 for (i = 0; i < n; i++)
303 pr_indent(fp, indent);
304 fprintf(fp, "%s[%d]={name=\"%s\", nr=%d, ic='%c'}\n",
305 title, bShowNumbers ? i : -1,
306 *(resinfo[i].name), resinfo[i].nr,
307 (resinfo[i].ic == '\0') ? ' ' : resinfo[i].ic);
312 void pr_atoms(FILE *fp, int indent, const char *title, const t_atoms *atoms,
313 gmx_bool bShownumbers)
315 if (available(fp, atoms, indent, title))
317 indent = pr_title(fp, indent, title);
318 pr_atom(fp, indent, "atom", atoms->atom, atoms->nr);
319 pr_strings(fp, indent, "atom", atoms->atomname, atoms->nr, bShownumbers);
320 pr_strings2(fp, indent, "type", atoms->atomtype, atoms->atomtypeB, atoms->nr, bShownumbers);
321 pr_resinfo(fp, indent, "residue", atoms->resinfo, atoms->nres, bShownumbers);
326 void pr_atomtypes(FILE *fp, int indent, const char *title, const t_atomtypes *atomtypes,
327 gmx_bool bShowNumbers)
330 if (available(fp, atomtypes, indent, title))
332 indent = pr_title(fp, indent, title);
333 for (i = 0; i < atomtypes->nr; i++)
335 pr_indent(fp, indent);
337 "atomtype[%3d]={atomnumber=%4d}\n",
338 bShowNumbers ? i : -1, atomtypes->atomnumber[i]);
343 static void cmp_atom(FILE *fp, int index, const t_atom *a1, const t_atom *a2, real ftol, real abstol)
347 cmp_us(fp, "atom.type", index, a1->type, a2->type);
348 cmp_us(fp, "atom.ptype", index, a1->ptype, a2->ptype);
349 cmp_int(fp, "atom.resind", index, a1->resind, a2->resind);
350 cmp_int(fp, "atom.atomnumber", index, a1->atomnumber, a2->atomnumber);
351 cmp_real(fp, "atom.m", index, a1->m, a2->m, ftol, abstol);
352 cmp_real(fp, "atom.q", index, a1->q, a2->q, ftol, abstol);
353 cmp_us(fp, "atom.typeB", index, a1->typeB, a2->typeB);
354 cmp_real(fp, "atom.mB", index, a1->mB, a2->mB, ftol, abstol);
355 cmp_real(fp, "atom.qB", index, a1->qB, a2->qB, ftol, abstol);
359 cmp_us(fp, "atom.type", index, a1->type, a1->typeB);
360 cmp_real(fp, "atom.m", index, a1->m, a1->mB, ftol, abstol);
361 cmp_real(fp, "atom.q", index, a1->q, a1->qB, ftol, abstol);
365 void cmp_atoms(FILE *fp, const t_atoms *a1, const t_atoms *a2, real ftol, real abstol)
369 fprintf(fp, "comparing atoms\n");
373 cmp_int(fp, "atoms->nr", -1, a1->nr, a2->nr);
374 for (i = 0; i < std::min(a1->nr, a2->nr); i++)
376 cmp_atom(fp, i, &(a1->atom[i]), &(a2->atom[i]), ftol, abstol);
381 for (i = 0; (i < a1->nr); i++)
383 cmp_atom(fp, i, &(a1->atom[i]), nullptr, ftol, abstol);
388 void atomsSetMassesBasedOnNames(t_atoms *atoms, gmx_bool printMissingMasses)
392 /* We could decide to anyhow assign then or generate a fatal error,
393 * but it's probably most useful to keep the masses we have.
398 int maxWarn = (printMissingMasses ? 10 : 0);
401 gmx_atomprop_t aps = gmx_atomprop_init();
403 gmx_bool haveMass = TRUE;
404 for (int i = 0; i < atoms->nr; i++)
406 if (!gmx_atomprop_query(aps, epropMass,
407 *atoms->resinfo[atoms->atom[i].resind].name,
413 if (numWarn < maxWarn)
415 fprintf(stderr, "Can not find mass in database for atom %s in residue %d %s\n",
417 atoms->resinfo[atoms->atom[i].resind].nr,
418 *atoms->resinfo[atoms->atom[i].resind].name);
427 atoms->haveMass = haveMass;
429 gmx_atomprop_destroy(aps);