3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
9 * VERSION 3.3.99_development_20071104
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-2006, 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 * Groningen Machine for Chemical Simulation
55 int main(int argc,char *argv[])
57 static char *desc[] = {
58 "[TT]g_anavel[tt] computes temperature profiles in a sample. The sample",
59 "can be analysed radial, i.e. the temperature as a function of",
60 "distance from the center, cylindrical, i.e. as a function of distance",
61 "from the vector (0,0,1) through the center of the box, or otherwise",
62 "(will be specified later)"
65 { efTRN, "-f", NULL, ffREAD },
66 { efTPX, "-s", NULL, ffREAD },
67 { efXPM, "-o", "xcm", ffWRITE }
69 #define NFILE asize(fnm)
71 static int mode = 0, nlevels = 10;
72 static real tmax = 300, xmax = -1;
74 { "-mode", FALSE, etINT, {&mode}, "mode" },
75 { "-nlevels", FALSE, etINT, {&nlevels}, "number of levels" },
76 { "-tmax", FALSE, etREAL, {&tmax}, "max temperature in output" },
77 { "-xmax", FALSE, etREAL, {&xmax}, "max distance from center" }
83 int i,j,idum,step,nframe=0,index;
84 real temp,rdum,hboxx,hboxy,scale,xnorm=0;
86 real *t_x=NULL,*t_y,hi=0;
91 gmx_bool bHaveV,bReadV;
92 t_rgb rgblo = { 0, 0, 1 },rgbhi = { 1, 0, 0 };
93 int flags = TRX_READ_X | TRX_READ_V;
97 CopyRight(stderr,argv[0]);
98 parse_common_args(&argc,argv,PCA_CAN_TIME | PCA_BE_NICE ,NFILE,fnm,
99 asize(pa),pa,asize(desc),desc,0,NULL);
101 top = read_top(ftp2fn(efTPX,NFILE,fnm));
103 read_first_frame(&status,ftp2fn(efTRX,NFILE,fnm),&fr,flags);
111 nmax = (0.5*sqrt(sqr(box[XX][XX])+sqr(box[YY][YY])))*scale;
115 for(i=0; (i<=nmax); i++) {
120 srenew(profile,++nframe);
121 snew(profile[nframe-1],nmax+1);
123 t_x[nframe-1] = fr.time*1000;
124 hboxx = box[XX][XX]/2;
125 hboxy = box[YY][YY]/2;
126 for(i=0; (i<fr.natoms); i++) {
127 /* determine position dependent on mode */
130 xnorm = sqrt(sqr(fr.x[i][XX]-hboxx) + sqr(fr.x[i][YY]-hboxy));
133 gmx_fatal(FARGS,"Unknown mode %d",mode);
137 temp = top->atoms.atom[i].m*iprod(fr.v[i],fr.v[i])/(2*BOLTZ);
141 profile[nframe-1][index] += temp;
144 for(i=0; (i<=nmax); i++) {
146 profile[nframe-1][i] /= npts[i];
149 } while (read_next_frame(status,&fr));
152 fp = ftp2FILE(efXPM,NFILE,fnm,"w");
153 write_xpm(fp,0,"Temp. profile","T (a.u.)",
155 nframe,nmax+1,t_x,t_y,profile,0,tmax,
156 rgblo,rgbhi,&nlevels);