Bug Summary

File:gromacs/gmxana/gmx_density.c
Location:line 525, column 9
Description:Array access results in a null pointer dereference

Annotated Source Code

1/*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2004, The GROMACS development team.
6 * Copyright (c) 2013,2014, 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.
10 *
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.
15 *
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.
20 *
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.
25 *
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.
33 *
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.
36 */
37#ifdef HAVE_CONFIG_H1
38#include <config.h>
39#endif
40
41#include <ctype.h>
42#include <math.h>
43#include <stdlib.h>
44#include <string.h>
45
46#include "gromacs/utility/cstringutil.h"
47#include "typedefs.h"
48#include "gromacs/utility/smalloc.h"
49#include "macros.h"
50#include "gstat.h"
51#include "gromacs/math/vec.h"
52#include "gromacs/fileio/xvgr.h"
53#include "viewit.h"
54#include "pbc.h"
55#include "gromacs/utility/futil.h"
56#include "gromacs/commandline/pargs.h"
57#include "index.h"
58#include "gromacs/fileio/tpxio.h"
59#include "gromacs/fileio/trxio.h"
60#include "physics.h"
61#include "gmx_ana.h"
62#include "macros.h"
63
64#include "gromacs/utility/fatalerror.h"
65
66typedef struct {
67 char *atomname;
68 int nr_el;
69} t_electron;
70
71/****************************************************************************/
72/* This program calculates the partial density across the box. */
73/* Peter Tieleman, Mei 1995 */
74/****************************************************************************/
75
76/* used for sorting the list */
77int compare(void *a, void *b)
78{
79 t_electron *tmp1, *tmp2;
80 tmp1 = (t_electron *)a; tmp2 = (t_electron *)b;
81
82 return strcmp(tmp1->atomname, tmp2->atomname)__extension__ ({ size_t __s1_len, __s2_len; (__builtin_constant_p
(tmp1->atomname) && __builtin_constant_p (tmp2->
atomname) && (__s1_len = strlen (tmp1->atomname), __s2_len
= strlen (tmp2->atomname), (!((size_t)(const void *)((tmp1
->atomname) + 1) - (size_t)(const void *)(tmp1->atomname
) == 1) || __s1_len >= 4) && (!((size_t)(const void
*)((tmp2->atomname) + 1) - (size_t)(const void *)(tmp2->
atomname) == 1) || __s2_len >= 4)) ? __builtin_strcmp (tmp1
->atomname, tmp2->atomname) : (__builtin_constant_p (tmp1
->atomname) && ((size_t)(const void *)((tmp1->atomname
) + 1) - (size_t)(const void *)(tmp1->atomname) == 1) &&
(__s1_len = strlen (tmp1->atomname), __s1_len < 4) ? (
__builtin_constant_p (tmp2->atomname) && ((size_t)
(const void *)((tmp2->atomname) + 1) - (size_t)(const void
*)(tmp2->atomname) == 1) ? __builtin_strcmp (tmp1->atomname
, tmp2->atomname) : (__extension__ ({ const unsigned char *
__s2 = (const unsigned char *) (const char *) (tmp2->atomname
); int __result = (((const unsigned char *) (const char *) (tmp1
->atomname))[0] - __s2[0]); if (__s1_len > 0 &&
__result == 0) { __result = (((const unsigned char *) (const
char *) (tmp1->atomname))[1] - __s2[1]); if (__s1_len >
1 && __result == 0) { __result = (((const unsigned char
*) (const char *) (tmp1->atomname))[2] - __s2[2]); if (__s1_len
> 2 && __result == 0) __result = (((const unsigned
char *) (const char *) (tmp1->atomname))[3] - __s2[3]); }
} __result; }))) : (__builtin_constant_p (tmp2->atomname)
&& ((size_t)(const void *)((tmp2->atomname) + 1) -
(size_t)(const void *)(tmp2->atomname) == 1) && (
__s2_len = strlen (tmp2->atomname), __s2_len < 4) ? (__builtin_constant_p
(tmp1->atomname) && ((size_t)(const void *)((tmp1
->atomname) + 1) - (size_t)(const void *)(tmp1->atomname
) == 1) ? __builtin_strcmp (tmp1->atomname, tmp2->atomname
) : (- (__extension__ ({ const unsigned char *__s2 = (const unsigned
char *) (const char *) (tmp1->atomname); int __result = (
((const unsigned char *) (const char *) (tmp2->atomname))[
0] - __s2[0]); if (__s2_len > 0 && __result == 0) {
__result = (((const unsigned char *) (const char *) (tmp2->
atomname))[1] - __s2[1]); if (__s2_len > 1 && __result
== 0) { __result = (((const unsigned char *) (const char *) (
tmp2->atomname))[2] - __s2[2]); if (__s2_len > 2 &&
__result == 0) __result = (((const unsigned char *) (const char
*) (tmp2->atomname))[3] - __s2[3]); } } __result; })))) :
__builtin_strcmp (tmp1->atomname, tmp2->atomname)))); }
)
;
83}
84
85int get_electrons(t_electron **eltab, const char *fn)
86{
87 char buffer[256]; /* to read in a line */
88 char tempname[80]; /* buffer to hold name */
89 int tempnr;
90
91 FILE *in;
92 int nr; /* number of atomstypes to read */
93 int i;
94
95 if (!(in = gmx_ffopen(fn, "r")))
96 {
97 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 97
, "Couldn't open %s. Exiting.\n", fn);
98 }
99
100 if (NULL((void*)0) == fgets(buffer, 255, in))
101 {
102 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 102
, "Error reading from file %s", fn);
103 }
104
105 if (sscanf(buffer, "%d", &nr) != 1)
106 {
107 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 107
, "Invalid number of atomtypes in datafile\n");
108 }
109
110 snew(*eltab, nr)(*eltab) = save_calloc("*eltab", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 110, (nr), sizeof(*(*eltab)))
;
111
112 for (i = 0; i < nr; i++)
113 {
114 if (fgets(buffer, 255, in) == NULL((void*)0))
115 {
116 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 116
, "reading datafile. Check your datafile.\n");
117 }
118 if (sscanf(buffer, "%s = %d", tempname, &tempnr) != 2)
119 {
120 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 120
, "Invalid line in datafile at line %d\n", i+1);
121 }
122 (*eltab)[i].nr_el = tempnr;
123 (*eltab)[i].atomname = strdup(tempname)(__extension__ (__builtin_constant_p (tempname) && ((
size_t)(const void *)((tempname) + 1) - (size_t)(const void *
)(tempname) == 1) ? (((const char *) (tempname))[0] == '\0' ?
(char *) calloc ((size_t) 1, (size_t) 1) : ({ size_t __len =
strlen (tempname) + 1; char *__retval = (char *) malloc (__len
); if (__retval != ((void*)0)) __retval = (char *) memcpy (__retval
, tempname, __len); __retval; })) : __strdup (tempname)))
;
124 }
125 gmx_ffclose(in);
126
127 /* sort the list */
128 fprintf(stderrstderr, "Sorting list..\n");
129 qsort ((void*)*eltab, nr, sizeof(t_electron),
130 (int(*)(const void*, const void*))compare);
131
132 return nr;
133}
134
135void center_coords(t_atoms *atoms, matrix box, rvec x0[], int axis)
136{
137 int i, m;
138 real tmass, mm;
139 rvec com, shift, box_center;
140
141 tmass = 0;
142 clear_rvec(com);
143 for (i = 0; (i < atoms->nr); i++)
144 {
145 mm = atoms->atom[i].m;
146 tmass += mm;
147 for (m = 0; (m < DIM3); m++)
148 {
149 com[m] += mm*x0[i][m];
150 }
151 }
152 for (m = 0; (m < DIM3); m++)
153 {
154 com[m] /= tmass;
155 }
156 calc_box_center(ecenterDEF, box, box_center);
157 rvec_sub(box_center, com, shift);
158 shift[axis] -= box_center[axis];
159
160 for (i = 0; (i < atoms->nr); i++)
161 {
162 rvec_dec(x0[i], shift);
163 }
164}
165
166void calc_electron_density(const char *fn, atom_id **index, int gnx[],
167 double ***slDensity, int *nslices, t_topology *top,
168 int ePBC,
169 int axis, int nr_grps, real *slWidth,
170 t_electron eltab[], int nr, gmx_bool bCenter,
171 const output_env_t oenv)
172{
173 rvec *x0; /* coordinates without pbc */
174 matrix box; /* box (3x3) */
175 double invvol;
176 int natoms; /* nr. atoms in trj */
177 t_trxstatus *status;
178 int i, n, /* loop indices */
179 nr_frames = 0, /* number of frames */
180 slice; /* current slice */
181 t_electron *found; /* found by bsearch */
182 t_electron sought; /* thingie thought by bsearch */
183 gmx_rmpbc_t gpbc = NULL((void*)0);
184
185 real t,
186 z;
187
188 if (axis < 0 || axis >= DIM3)
189 {
190 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 190
, "Invalid axes. Terminating\n");
191 }
192
193 if ((natoms = read_first_x(oenv, &status, fn, &t, &x0, box)) == 0)
194 {
195 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 195
, "Could not read coordinates from statusfile\n");
196 }
197
198 if (!*nslices)
199 {
200 *nslices = (int)(box[axis][axis] * 10); /* default value */
201 }
202 fprintf(stderrstderr, "\nDividing the box in %d slices\n", *nslices);
203
204 snew(*slDensity, nr_grps)(*slDensity) = save_calloc("*slDensity", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 204, (nr_grps), sizeof(*(*slDensity)))
;
205 for (i = 0; i < nr_grps; i++)
206 {
207 snew((*slDensity)[i], *nslices)((*slDensity)[i]) = save_calloc("(*slDensity)[i]", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 207, (*nslices), sizeof(*((*slDensity)[i])))
;
208 }
209
210 gpbc = gmx_rmpbc_init(&top->idef, ePBC, top->atoms.nr);
211 /*********** Start processing trajectory ***********/
212 do
213 {
214 gmx_rmpbc(gpbc, natoms, box, x0);
215
216 if (bCenter)
217 {
218 center_coords(&top->atoms, box, x0, axis);
219 }
220
221 *slWidth = box[axis][axis]/(*nslices);
222 invvol = *nslices/(box[XX0][XX0]*box[YY1][YY1]*box[ZZ2][ZZ2]);
223
224 for (n = 0; n < nr_grps; n++)
225 {
226 for (i = 0; i < gnx[n]; i++) /* loop over all atoms in index file */
227 {
228 z = x0[index[n][i]][axis];
229 while (z < 0)
230 {
231 z += box[axis][axis];
232 }
233 while (z > box[axis][axis])
234 {
235 z -= box[axis][axis];
236 }
237
238 /* determine which slice atom is in */
239 slice = (z / (*slWidth));
240 sought.nr_el = 0;
241 sought.atomname = strdup(*(top->atoms.atomname[index[n][i]]))(__extension__ (__builtin_constant_p (*(top->atoms.atomname
[index[n][i]])) && ((size_t)(const void *)((*(top->
atoms.atomname[index[n][i]])) + 1) - (size_t)(const void *)(*
(top->atoms.atomname[index[n][i]])) == 1) ? (((const char *
) (*(top->atoms.atomname[index[n][i]])))[0] == '\0' ? (char
*) calloc ((size_t) 1, (size_t) 1) : ({ size_t __len = strlen
(*(top->atoms.atomname[index[n][i]])) + 1; char *__retval
= (char *) malloc (__len); if (__retval != ((void*)0)) __retval
= (char *) memcpy (__retval, *(top->atoms.atomname[index[
n][i]]), __len); __retval; })) : __strdup (*(top->atoms.atomname
[index[n][i]]))))
;
242
243 /* now find the number of electrons. This is not efficient. */
244 found = (t_electron *)
245 bsearch((const void *)&sought,
246 (const void *)eltab, nr, sizeof(t_electron),
247 (int(*)(const void*, const void*))compare);
248
249 if (found == NULL((void*)0))
250 {
251 fprintf(stderrstderr, "Couldn't find %s. Add it to the .dat file\n",
252 *(top->atoms.atomname[index[n][i]]));
253 }
254 else
255 {
256 (*slDensity)[n][slice] += (found->nr_el -
257 top->atoms.atom[index[n][i]].q)*invvol;
258 }
259 free(sought.atomname);
260 }
261 }
262 nr_frames++;
263 }
264 while (read_next_x(oenv, status, &t, x0, box));
265 gmx_rmpbc_done(gpbc);
266
267 /*********** done with status file **********/
268 close_trj(status);
269
270/* slDensity now contains the total number of electrons per slice, summed
271 over all frames. Now divide by nr_frames and volume of slice
272 */
273
274 fprintf(stderrstderr, "\nRead %d frames from trajectory. Counting electrons\n",
275 nr_frames);
276
277 for (n = 0; n < nr_grps; n++)
278 {
279 for (i = 0; i < *nslices; i++)
280 {
281 (*slDensity)[n][i] /= nr_frames;
282 }
283 }
284
285 sfree(x0)save_free("x0", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 285, (x0))
; /* free memory used by coordinate array */
286}
287
288void calc_density(const char *fn, atom_id **index, int gnx[],
289 double ***slDensity, int *nslices, t_topology *top, int ePBC,
290 int axis, int nr_grps, real *slWidth, gmx_bool bCenter,
291 const output_env_t oenv)
292{
293 rvec *x0; /* coordinates without pbc */
294 matrix box; /* box (3x3) */
295 double invvol;
296 int natoms; /* nr. atoms in trj */
297 t_trxstatus *status;
298 int **slCount, /* nr. of atoms in one slice for a group */
299 i, j, n, /* loop indices */
300 teller = 0,
301 ax1 = 0, ax2 = 0,
302 nr_frames = 0, /* number of frames */
303 slice; /* current slice */
304 real t,
305 z;
306 char *buf; /* for tmp. keeping atomname */
307 gmx_rmpbc_t gpbc = NULL((void*)0);
308
309 if (axis < 0 || axis >= DIM3)
310 {
311 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 311
, "Invalid axes. Terminating\n");
312 }
313
314 if ((natoms = read_first_x(oenv, &status, fn, &t, &x0, box)) == 0)
315 {
316 gmx_fatal(FARGS0, "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 316
, "Could not read coordinates from statusfile\n");
317 }
318
319 if (!*nslices)
320 {
321 *nslices = (int)(box[axis][axis] * 10); /* default value */
322 fprintf(stderrstderr, "\nDividing the box in %d slices\n", *nslices);
323 }
324
325 snew(*slDensity, nr_grps)(*slDensity) = save_calloc("*slDensity", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 325, (nr_grps), sizeof(*(*slDensity)))
;
326 for (i = 0; i < nr_grps; i++)
327 {
328 snew((*slDensity)[i], *nslices)((*slDensity)[i]) = save_calloc("(*slDensity)[i]", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 328, (*nslices), sizeof(*((*slDensity)[i])))
;
329 }
330
331 gpbc = gmx_rmpbc_init(&top->idef, ePBC, top->atoms.nr);
332 /*********** Start processing trajectory ***********/
333 do
334 {
335 gmx_rmpbc(gpbc, natoms, box, x0);
336
337 if (bCenter)
338 {
339 center_coords(&top->atoms, box, x0, axis);
340 }
341
342 *slWidth = box[axis][axis]/(*nslices);
343 invvol = *nslices/(box[XX0][XX0]*box[YY1][YY1]*box[ZZ2][ZZ2]);
344 teller++;
345
346 for (n = 0; n < nr_grps; n++)
347 {
348 for (i = 0; i < gnx[n]; i++) /* loop over all atoms in index file */
349 {
350 z = x0[index[n][i]][axis];
351 while (z < 0)
352 {
353 z += box[axis][axis];
354 }
355 while (z > box[axis][axis])
356 {
357 z -= box[axis][axis];
358 }
359
360 /* determine which slice atom is in */
361 slice = (int)(z / (*slWidth));
362 (*slDensity)[n][slice] += top->atoms.atom[index[n][i]].m*invvol;
363 }
364 }
365 nr_frames++;
366 }
367 while (read_next_x(oenv, status, &t, x0, box));
368 gmx_rmpbc_done(gpbc);
369
370 /*********** done with status file **********/
371 close_trj(status);
372
373 /* slDensity now contains the total mass per slice, summed over all
374 frames. Now divide by nr_frames and volume of slice
375 */
376
377 fprintf(stderrstderr, "\nRead %d frames from trajectory. Calculating density\n",
378 nr_frames);
379
380 for (n = 0; n < nr_grps; n++)
381 {
382 for (i = 0; i < *nslices; i++)
383 {
384 (*slDensity)[n][i] /= nr_frames;
385 }
386 }
387
388 sfree(x0)save_free("x0", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 388, (x0))
; /* free memory used by coordinate array */
389}
390
391void plot_density(double *slDensity[], const char *afile, int nslices,
392 int nr_grps, char *grpname[], real slWidth,
393 const char **dens_opt,
394 gmx_bool bSymmetrize, const output_env_t oenv)
395{
396 FILE *den;
397 const char *ylabel = NULL((void*)0);
398 int slice, n;
399 real ddd;
400
401 switch (dens_opt[0][0])
402 {
403 case 'm': ylabel = "Density (kg m\\S-3\\N)"; break;
404 case 'n': ylabel = "Number density (nm\\S-3\\N)"; break;
405 case 'c': ylabel = "Charge density (e nm\\S-3\\N)"; break;
406 case 'e': ylabel = "Electron density (e nm\\S-3\\N)"; break;
407 }
408
409 den = xvgropen(afile, "Partial densities", "Box (nm)", ylabel, oenv);
410
411 xvgr_legend(den, nr_grps, (const char**)grpname, oenv);
412
413 for (slice = 0; (slice < nslices); slice++)
414 {
415 fprintf(den, "%12g ", slice * slWidth);
416 for (n = 0; (n < nr_grps); n++)
417 {
418 if (bSymmetrize)
419 {
420 ddd = (slDensity[n][slice]+slDensity[n][nslices-slice-1])*0.5;
421 }
422 else
423 {
424 ddd = slDensity[n][slice];
425 }
426 if (dens_opt[0][0] == 'm')
427 {
428 fprintf(den, " %12g", ddd*AMU(1.6605402e-27)/(NANO(1e-9)*NANO(1e-9)*NANO(1e-9)));
429 }
430 else
431 {
432 fprintf(den, " %12g", ddd);
433 }
434 }
435 fprintf(den, "\n");
436 }
437
438 gmx_ffclose(den);
439}
440
441int gmx_density(int argc, char *argv[])
442{
443 const char *desc[] = {
444 "[THISMODULE] computes partial densities across the box, using an index file.[PAR]",
445 "For the total density of NPT simulations, use [gmx-energy] instead.",
446 "[PAR]",
447 "Densities are in kg/m^3, and number densities or electron densities can also be",
448 "calculated. For electron densities, a file describing the number of",
449 "electrons for each type of atom should be provided using [TT]-ei[tt].",
450 "It should look like:[BR]",
451 " [TT]2[tt][BR]",
452 " [TT]atomname = nrelectrons[tt][BR]",
453 " [TT]atomname = nrelectrons[tt][BR]",
454 "The first line contains the number of lines to read from the file.",
455 "There should be one line for each unique atom name in your system.",
456 "The number of electrons for each atom is modified by its atomic",
457 "partial charge."
458 };
459
460 output_env_t oenv;
461 static const char *dens_opt[] =
462 { NULL((void*)0), "mass", "number", "charge", "electron", NULL((void*)0) };
463 static int axis = 2; /* normal to memb. default z */
464 static const char *axtitle = "Z";
465 static int nslices = 50; /* nr of slices defined */
466 static int ngrps = 1; /* nr. of groups */
467 static gmx_bool bSymmetrize = FALSE0;
468 static gmx_bool bCenter = FALSE0;
469 t_pargs pa[] = {
470 { "-d", FALSE0, etSTR, {&axtitle},
471 "Take the normal on the membrane in direction X, Y or Z." },
472 { "-sl", FALSE0, etINT, {&nslices},
473 "Divide the box in this number of slices." },
474 { "-dens", FALSE0, etENUM, {dens_opt},
475 "Density"},
476 { "-ng", FALSE0, etINT, {&ngrps},
477 "Number of groups of which to compute densities." },
478 { "-symm", FALSE0, etBOOL, {&bSymmetrize},
479 "Symmetrize the density along the axis, with respect to the center. Useful for bilayers." },
480 { "-center", FALSE0, etBOOL, {&bCenter},
481 "Shift the center of mass along the axis to zero. This means if your axis is Z and your box is bX, bY, bZ, the center of mass will be at bX/2, bY/2, 0."}
482 };
483
484 const char *bugs[] = {
485 "When calculating electron densities, atomnames are used instead of types. This is bad.",
486 };
487
488 double **density; /* density per slice */
489 real slWidth; /* width of one slice */
490 char **grpname; /* groupnames */
491 int nr_electrons; /* nr. electrons */
492 int *ngx; /* sizes of groups */
493 t_electron *el_tab; /* tabel with nr. of electrons*/
494 t_topology *top; /* topology */
495 int ePBC;
496 atom_id **index; /* indices for all groups */
497 int i;
498
499 t_filenm fnm[] = { /* files for g_density */
500 { efTRX, "-f", NULL((void*)0), ffREAD1<<1 },
501 { efNDX, NULL((void*)0), NULL((void*)0), ffOPTRD(1<<1 | 1<<3) },
502 { efTPX, NULL((void*)0), NULL((void*)0), ffREAD1<<1 },
503 { efDAT, "-ei", "electrons", ffOPTRD(1<<1 | 1<<3) }, /* file with nr. of electrons */
504 { efXVG, "-o", "density", ffWRITE1<<2 },
505 };
506
507#define NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))) asize(fnm)((int)(sizeof(fnm)/sizeof((fnm)[0])))
508
509 if (!parse_common_args(&argc, argv, PCA_CAN_VIEW(1<<5) | PCA_CAN_TIME((1<<6) | (1<<7) | (1<<14)) | PCA_BE_NICE(1<<13),
1
Taking false branch
510 NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm, asize(pa)((int)(sizeof(pa)/sizeof((pa)[0]))), pa, asize(desc)((int)(sizeof(desc)/sizeof((desc)[0]))), desc, asize(bugs)((int)(sizeof(bugs)/sizeof((bugs)[0]))), bugs,
511 &oenv))
512 {
513 return 0;
514 }
515
516 if (bSymmetrize && !bCenter)
517 {
518 fprintf(stderrstderr, "Can not symmetrize without centering. Turning on -center\n");
519 bCenter = TRUE1;
520 }
521 /* Calculate axis */
522 axis = toupper(axtitle[0])(__extension__ ({ int __res; if (sizeof (axtitle[0]) > 1) {
if (__builtin_constant_p (axtitle[0])) { int __c = (axtitle[
0]); __res = __c < -128 || __c > 255 ? __c : (*__ctype_toupper_loc
())[__c]; } else __res = toupper (axtitle[0]); } else __res =
(*__ctype_toupper_loc ())[(int) (axtitle[0])]; __res; }))
- 'X';
523
524 top = read_top(ftp2fn(efTPX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), &ePBC); /* read topology file */
525 if (dens_opt[0][0] == 'n')
2
Array access results in a null pointer dereference
526 {
527 for (i = 0; (i < top->atoms.nr); i++)
528 {
529 top->atoms.atom[i].m = 1;
530 }
531 }
532 else if (dens_opt[0][0] == 'c')
533 {
534 for (i = 0; (i < top->atoms.nr); i++)
535 {
536 top->atoms.atom[i].m = top->atoms.atom[i].q;
537 }
538 }
539
540 snew(grpname, ngrps)(grpname) = save_calloc("grpname", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 540, (ngrps), sizeof(*(grpname)))
;
541 snew(index, ngrps)(index) = save_calloc("index", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 541, (ngrps), sizeof(*(index)))
;
542 snew(ngx, ngrps)(ngx) = save_calloc("ngx", "/home/alexxy/Develop/gromacs/src/gromacs/gmxana/gmx_density.c"
, 542, (ngrps), sizeof(*(ngx)))
;
543
544 get_index(&top->atoms, ftp2fn_null(efNDX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), ngrps, ngx, index, grpname);
545
546 if (dens_opt[0][0] == 'e')
547 {
548 nr_electrons = get_electrons(&el_tab, ftp2fn(efDAT, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm));
549 fprintf(stderrstderr, "Read %d atomtypes from datafile\n", nr_electrons);
550
551 calc_electron_density(ftp2fn(efTRX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), index, ngx, &density,
552 &nslices, top, ePBC, axis, ngrps, &slWidth, el_tab,
553 nr_electrons, bCenter, oenv);
554 }
555 else
556 {
557 calc_density(ftp2fn(efTRX, NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), index, ngx, &density, &nslices, top,
558 ePBC, axis, ngrps, &slWidth, bCenter, oenv);
559 }
560
561 plot_density(density, opt2fn("-o", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm),
562 nslices, ngrps, grpname, slWidth, dens_opt,
563 bSymmetrize, oenv);
564
565 do_view(oenv, opt2fn("-o", NFILE((int)(sizeof(fnm)/sizeof((fnm)[0]))), fnm), "-nxy"); /* view xvgr file */
566 return 0;
567}