[gmx-developers] reading velocity in gmx_msd.c

[Mark Abraham]

Did you use a debugger to find out why your segfault is arising?

What GROMACS version are you using? 

A diff against the code with a few lines of context is much easier for us to understand, if we need to look at your code :-)

I think I know what the problem is, but you haven't given enough information for me to be sure :-)

Mark

----- Original Message -----
From: Gaurav Goel <gauravgoeluta at gmail.com>
Date: Saturday, May 22, 2010 6:18
Subject: Re: [gmx-developers] reading velocity in gmx_msd.c
To: Discussion list for GROMACS development <gmx-developers at gromacs.org>

> **************Please ignore the previous mail... I accidentally
> pressed enter before completing the mail***************
> 
> 
> To be able to read velocities in g_msd I decided to do the following:
> 
> 1. Add the following lines to src/gmxlib/trxio.c
> 
> int read_first_xv(int *status,char *fn,
>                  real *t,rvec **x, rvec **v,matrix box)
> {
>  t_trxframe fr;
> 
>  read_first_frame(status,fn,&fr,TRX_NEED_X | TRX_NEED_V);
>  if (*status >= nxframe) {
>    nxframe = *status+1;
>    srenew(xframe,nxframe);
>  }
>  xframe[*status] = fr;
>  *t = xframe[*status].time;
>  *x = xframe[*status].x;
>  *v = xframe[*status].v;
>  copy_mat(xframe[*status].box,box);
> 
>  return xframe[*status].natoms;
> }
> 
> 
> bool read_next_xv(int status,real *t, int natoms, rvec x[], rvec 
> v[], matrix box
> )
> {
>  bool bRet;
> 
>  xframe[status].x = x;
>  xframe[status].v = v;
>  bRet = read_next_frame(status,&xframe[status]);
>  *t = xframe[status].time;
>  copy_mat(xframe[status].box,box);
> 
>  return bRet;
> }
> 
> 2. Add the following lines to /include/statutil.h
> extern int read_first_xv(int *status,char *fn,
>                        real *t,rvec **x, rvec **v,matrix box);
> extern bool read_next_xv(int status,real *t,int natoms,rvec x[], 
> rvec v[],matrix
>  box);
> 
> 3. Make following changes to the gmx_msd.c subroutine:
> /*define the followign pointers*/
> rvec         *v[2],*va[2];
> 
> /*read x and v with these commands*/
> natoms = read_first_xv(&status,fn,&this->t0,&(x[cur]),&(v[cur]),box);
> while (read_next_xv(status,&t,natoms,x[cur],v[cur],box));
> 
> /*repeat the following commands as used for x*/
> snew(v[prev],natoms);
>   va[0] = v[0];
>    va[1] = v[1];
> srenew(this->v0,this->nrestart);
>      snew(this->v0[this->nrestart-1],this-
> >ncoords);
>    if (bFirst) {
> 
>      memcpy(xa[prev],xa[cur],this-
> >ncoords*sizeof(xa[prev][0]));/*error point*/
> memcpy(va[prev],va[cur],this->ncoords*sizeof(va[prev][0]));
> bFirst = False;
> }
> 
> Everything compiled all right. However, when I tried to use this new
> utility I get a
> 'segmentation fault' error at
> ***memcpy(va[prev],va[cur],this->ncoords*sizeof(va[prev][0]));****



> Am I not reading the velocities properly? Is my strategy of repeating
> everything that was done for reading positions for velocities (except
> for doing things like remove pbc) not correct?
> 
> I will greatly appreciate your help.
> 
> Thanks,
> Gaurav
> 
> On Tue, May 18, 2010 at 5:16 PM, David van der Spoel
> <spoel at xray.bmc.uu.se> wrote:
> > On 2010-05-18 23.10, Gaurav Goel wrote:
> >>
> >> I am modifying the 'g_msd' utility to be able to calculate 
> viscosities>> using the corresponding Einstein relationship 
> according to equation
> >> 3.14 in Helfand E., Physical Review, 119, 1, 1960. The 
> equation is
> >> similar to the mean-square displacement equation for self-diffusion
> >> coefficient with the position vector 'r_i ' replaced by 'r_ix p_iy'
> >> where p_iy is the particle momenta in y-direction.
> >>
> >> So basically I need to be able to read both position and velocities
> >> from the trajectory. Presently, gmx_msd.c uses 'read_first_x' and
> >> 'read_next_x' to read the coordinates. Can you please comment 
> on how I
> >> can modify this code to read the corresponding velocities?
> >>
> >> Any help is greatly appreciated.
> >>
> >> thanks,
> >> Gaurav
> >
> > From include/statutil.h:
> >
> > extern bool read_next_frame(const output_env_t oenv,int 
> status,t_trxframe> *fr);
> >  /* Reads the next frame which is in accordance with fr->flags.
> >   * Returns TRUE when succeeded, FALSE otherwise.
> >   */
> >
> > Do also checkout
> > g_tcaf -h
> >
> >
> > --
> > David van der Spoel, Ph.D., Professor of Biology
> > Dept. of Cell & Molec. Biol., Uppsala University.
> > Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
> > spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
> > --
> > gmx-developers mailing list
> > gmx-developers at gromacs.org
> > http://lists.gromacs.org/mailman/listinfo/gmx-developers
> > Please don't post (un)subscribe requests to the list. Use the 
> www interface
> > or send it to gmx-developers-request at gromacs.org.
> >
> -- 
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the 
> www interface or send it to gmx-developers-request at gromacs.org.