[gmx-developers] reading velocity in gmx_msd.c
Gaurav Goel
gauravgoeluta at gmail.com
Mon May 24 17:07:01 CEST 2010
(gdb) where
#0 0x0000003e2b482b6b in memcpy () from /lib64/libc.so.6
#1 0x0000000000403992 in corr_loop ()
#2 0x00000000004066df in do_corr ()
#3 0x0000000000407f84 in gmx_visco ()
#4 0x000000000040226b in main ()
-Gaurav
On Mon, May 24, 2010 at 11:01 AM, David van der Spoel
<spoel at xray.bmc.uu.se> wrote:
> On 2010-05-24 16.47, Gaurav Goel wrote:
>>
>> I had allocated memory for both v and va using
>> snew(v[prev],natoms);
>>
>> and
>> snew(va[0],this->ncoords);
>> snew(va[1],this->ncoords);
>>
>> Results from the DEBUG script:
>>
>> $: more DEBUG
>> #!/bin/bash -f
>> xterm -e gdb -x gdb_cmds
>> /home/gxg133/Download/gromacs-4.0.5/src/tools/g_visco
>>
>> $: more gdb_cmds
>> break _gmx_error
>> break gmx_visco
>> run
>>
>> $: ./DEBUG
>>
>> Breakpoint 1 at 0x411670
>> Breakpoint 2 at 0x4079f0
>>
>> Breakpoint 2, 0x00000000004079f0 in gmx_visco ()
>> Missing separate debuginfos, use: debuginfo-install
>> e2fsprogs-libs-1.41.4-12.fc11.x86_64 glibc-2.10.2-1.x86_64
>> libICE-1.0.4-7.fc11.x86_64 libSM-1.1.0-4.fc11.x86_64
>> libX11-1.2.2-1.fc11.x86_64 libXau-1.0.4-5.fc11.x86_64
>> libxcb-1.2-5.fc11.x86_64
>>
>> (gdb) n
>> Single stepping until exit from function gmx_visco,
>> which has no line number information.
>> :-) G R O M A C S (-:
>>
>> GRoups of Organic Molecules in ACtion for Science
>>
>> :-) VERSION 4.0.5 (-
>> .....
>> .....
>>
>> Group 0 ( System) has 2160 elements
>> Group 1 ( Si1) has 720 elements
>> Group 2 ( O1) has 1440 elements
>> Select a group: 1
>> Selected 1: 'Si1'
>> Reading frame 0 time 0.000
>> Program received signal SIGSEGV, Segmentation fault.
>> 0x0000003e2b482b6b in memcpy () from /lib64/libc.so.6
>> ----------------------
>>
>> So it seems the error is happening in "memcpy()". Any suggestions on
>> how I proceed from here?
>>
> In an interactive debugger session you go up up up until you can see where
> it is going wrong (or type "where").
>
>
>> -Gaurav
>>
>> On Sat, May 22, 2010 at 10:44 AM, Mark Abraham<mark.abraham at anu.edu.au>
>> wrote:
>>>
>>> ----- Original Message -----
>>> From: David van der Spoel<spoel at xray.bmc.uu.se>
>>> Date: Sunday, May 23, 2010 0:08
>>> Subject: Re: [gmx-developers] reading velocity in gmx_msd.c
>>> To: Discussion list for GROMACS development<gmx-developers at gromacs.org>
>>>
>>>> On 2010-05-22 15.56, Gaurav Goel wrote:
>>>>>
>>>>> I've not used the debugger. I located the error point by
>>>>
>>>> writing to
>>>>>
>>>>> stdout from different blocks in the code.
>>>
>>> That works for a while, but what you really want to know is whether
>>> there's memory allocated for v, which a debugger will tell you right away.
>>> If you look back at the code for xa, there's an if statement regarding bMol
>>> that allocates xa or not. You probably didn't do that for va.
>>>
>>> Mark
>>>
>>>> did you allocate memory for the velocities?
>>>>>
>>>>> I'm using v. 4.0.5.
>>>>>
>>>>> Goal: Calculate viscosity with an analogous mean-squared displacement
>>>>> method since it converges much faster. So, I need to practically
>>>>> repeat the diffusivity calculation except that [r_i(t) -
>>>>
>>>> r_i(0)]**2 is
>>>>>
>>>>> replaced by [r_ix(t)*p_iy(t) - r_ix(0)p_iy(0)]**2. This method
>>>>> converges much faster... advantages and constraints on
>>>>
>>>> application of
>>>>>
>>>>> this method are in Phys. Rev. A, 43(8), 4289 (1991).
>>>>>
>>>>> Diff against the code:
>>>>> -----------------
>>>>> 1d0
>>>>> < /*
>>>>> 75a75
>>>>>>
>>>>>> rvec **v0;
>>>>>
>>>>> 83c83
>>>>> < typedef real t_calc_func(t_corr
>>>>
>>>> *,int,atom_id[],int,rvec[],rvec,bool,matrix);>---
>>>>>>
>>>>>> typedef real t_calc_func(t_corr
>>>>
>>>> *,int,atom_id[],int,rvec[],rvec[],rvec,bool,matrix);>111a112
>>>>>>
>>>>>> this-
>>>>>
>>>>> v0 = NULL;
>>>>> 153c154
>>>>> < for(i=0; (i<this->nrestart); i++)
>>>>> ---
>>>>>>
>>>>>> for(i=0; (i<this->nrestart); i++) {
>>>>>
>>>>> 154a156,157
>>>>>>
>>>>>> sfree(this->v0[i]);
>>>>>> }
>>>>>
>>>>> 155a159
>>>>>>
>>>>>> sfree(this->v0);
>>>>>
>>>>> 195c199
>>>>> < static void calc_corr(t_corr *this,int nr,int
>>>>
>>>> nx,atom_id index[],rvec xc[],
>>>>>
>>>>> ---
>>>>>>
>>>>>> static void calc_corr(t_corr *this,int nr,int nx,atom_id
>>>>
>>>> index[],rvec xc[],rvec vc[],
>>>>>
>>>>> 206a211
>>>>>>
>>>>>> memcpy(this->v0[this-
>>>>>
>>>>> nlast],vc,this->ncoords*sizeof(vc[0]));
>>>>> 222c227
>>>>> < g =
>>>>
>>>> calc1(this,nx,index,nx0,xc,dcom,bTen,mat);>---
>>>>>>
>>>>>> g =
>>>>
>>>> calc1(this,nx,index,nx0,xc,vc,dcom,bTen,mat);>235c240
>>>>>
>>>>> < static real calc1_norm(t_corr *this,int nx,atom_id
>>>>
>>>> index[],int nx0,rvec xc[],
>>>>>
>>>>> ---
>>>>>>
>>>>>> static real calc1_norm(t_corr *this,int nx,atom_id index[],int
>>>>
>>>> nx0,rvec xc[],rvec vc[],
>>>>>
>>>>> 306c311
>>>>> < static real calc_one_mw(t_corr *this,int ix,int
>>>>
>>>> nx0,rvec xc[],real *tm,
>>>>>
>>>>> ---
>>>>>>
>>>>>> static real calc_one_mw(t_corr *this,int ix,int nx0,rvec
>>>>
>>>> xc[],rvec vc[],real *tm,
>>>>>
>>>>> 350c355
>>>>> < static real calc1_mw(t_corr *this,int nx,atom_id
>>>>
>>>> index[],int nx0,rvec xc[],
>>>>>
>>>>> ---
>>>>>>
>>>>>> static real calc1_mw(t_corr *this,int nx,atom_id index[],int
>>>>
>>>> nx0,rvec xc[],rvec vc[],
>>>>>
>>>>> 359c364
>>>>> < g +=
>>>>
>>>> calc_one_mw(this,index[i],nx0,xc,&tm,dcom,bTen,mat);>---
>>>>>>
>>>>>> g +=
>>>>
>>>> calc_one_mw(this,index[i],nx0,xc,vc,&tm,dcom,bTen,mat);>429c434
>>>>>
>>>>> < static real calc1_mol(t_corr *this,int nx,atom_id
>>>>
>>>> index[],int nx0,rvec xc[],
>>>>>
>>>>> ---
>>>>>>
>>>>>> static real calc1_mol(t_corr *this,int nx,atom_id index[],int
>>>>
>>>> nx0,rvec xc[],rvec vc[],
>>>>>
>>>>> 440c445
>>>>> < g =
>>>>
>>>> calc_one_mw(this,i,nx0,xc,&mm,dcom,bTen,mat);>---
>>>>>>
>>>>>> g =
>>>>
>>>> calc_one_mw(this,i,nx0,xc,vc,&mm,dcom,bTen,mat);>529a535
>>>>>>
>>>> rvec *v[2],*va[2];
>>>>>
>>>>> 536c542,543
>>>>> < natoms = read_first_x(&status,fn,&this-
>>>>> t0,&(x[cur]),box);>---
>>>>>>
>>>>>> natoms = read_first_xv(&status,fn,&this-
>>>>>
>>>>> t0,&(x[cur]),&(v[cur]),box);>>
>>>>> 545c552
>>>>> <
>>>>> ---
>>>>>>
>>>>>> snew(v[prev],natoms);
>>>>>
>>>>> 549a557,558
>>>>>>
>>>>>> snew(va[0],this->ncoords);
>>>>>> snew(va[1],this->ncoords);
>>>>>
>>>>> 553a563,564
>>>>>>
>>>>>> va[0] = v[0];
>>>>>> va[1] = v[1];
>>>>>
>>>>> 575d585
>>>>> <
>>>>> 577a588,591
>>>>>>
>>>>>> srenew(this->v0,this-
>>>>>
>>>>> nrestart);>> snew(this-
>>>>> v0[this->nrestart-1],this->ncoords);
>>>>>>
>>>>> 615a630
>>>>>>
>>>> memcpy(va[prev],va[cur],this->ncoords*sizeof(va[prev][0]));
>>>>>
>>>>> 628c643,644
>>>>> <
>>>>> ---
>>>>>>
>>>>>>
>>>>> 631a648
>>>>>>
>>>>> 635c652,653
>>>>> <
>>>>
>>>> calc_corr(this,i,gnx[i],index[i],xa[cur],bRmCOMM,com,calc1,bTen);>---
>>>>>>
>>>>>>
>>>>
>>>> calc_corr(this,i,gnx[i],index[i],xa[cur],va[cur],bRmCOMM,com,calc1,bTen);>641c659
>>>>>
>>>>> < } while
>>>>
>>>> (read_next_x(status,&t,natoms,x[cur],box));>---
>>>>>>
>>>>>> } while
>>>>
>>>> (read_next_xv(status,&t,natoms,x[cur],v[cur],box));>905a924,925
>>>>>>
>>>>>>
>>>>> 946c966
>>>>> <
>>>>> ---
>>>>>>
>>>>> --------------------------------------------------------
>>>>>
>>>>> -Gaurav
>>>>>
>>>>> On Fri, May 21, 2010 at 8:30 PM, Mark
>>>>
>>>> Abraham<mark.abraham at anu.edu.au> wrote:
>>>>>>
>>>>>> 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.
>>>>>>
>>>>>> --
>>>>>> 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.
>>>>>>
>>>>
>>>>
>>>> --
>>>> 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.
>>>
>
>
> --
> 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.
>
More information about the gromacs.org_gmx-developers
mailing list