[gmx-users] problems with g_wham
Gavin Melaugh
gmelaugh01 at qub.ac.uk
Wed Jul 7 19:39:06 CEST 2010
Hi Chris
Thanks very much for your help. The reason I used pull_init1 = 0 is
because this is the value that is used in the tutorial. Maybe I am
misinterpreting its meaning. In the tutorial it says that you do not
have to specify a value separately for each configuration. I am now
thinking that in my case this does not apply ? Should I set pull_init to
the distance between the molecules at the first configuration?
Many thanks
Gavin
chris.neale at utoronto.ca wrote:
> Gavin, let's take g_wham out of the picture for now to simplify
> things. First, you umbrella is pretty strange. You want them to have a
> distance of 0.0 nm (pull_init1 = 0)???
>
> So let's look at making a histogram from your raw data.
>
> 1. cat pullx.xvg |awk '{print sqrt($1,$5*$5+$6*$6+$7*$7)}' > my.dist
> 2. ./histo my.dist 0.01 > my.histo
> ** select the bin width 0.01 so that the curve is still smooth but
> also keep the number small for precision.
> ** compile the histo.c code below to get ./histo as follows (gcc
> histo.c -o histo -lm)
> 3. Make a plot of the histogram and post it on the internet (e.g.
> using photobucket website)
> 4. reply with a link to your histogram.
>
> #### PS: i make no claim about the accuracy of the histo.c program so
> if you use it for other activities then its at your own risk. In fact,
> I think that there are cases where the very fist value is reported
> incorrectly, but that should not matter here.
>
> #include <stdio.h>
> #include <stdlib.h>
> #include "math.h"
>
> const int LINE_SIZE=1000;
>
> void showUsage(const char *c){
> printf("Usage: %s <input_file> <binwidth> [<min> <max>]\n",c);
> }
>
> int main(int argn, char *args[]){
> char linein[LINE_SIZE];
> FILE *f;
> int *data;
> int nbin,tot,i;
> double min,max,bin,c;
>
> if ((argn!=3 &&argn!=5)|| sscanf(args[2],"%lf",&bin)!=1){
> showUsage(args[0]);
> exit(1);
> }
>
> if((f=fopen(args[1],"r"))==NULL){
> printf("Error: unable to open file\n");
> showUsage(args[0]);
> exit(1);
> }
> while(fgets(linein,LINE_SIZE,f)!=NULL){
> if((sscanf(linein,"%lf",&c))!=1)continue;
> if(c<min)min=c;
> if(c>max)max=c;
> }
> fclose(f);
>
> if (argn==5){
> if(sscanf(args[3],"%lf",&min)!=1 || sscanf(args[4],"%lf",&max)!=1){
> showUsage(args[0]);
> exit(1);
> }
> }
>
> min=floor(min/bin)*bin;
> max=ceil(max/bin)*bin;
>
> nbin=(int)ceil((max-min)/bin);
> data=(int *)malloc((nbin+1)*sizeof(int));
> if(data==NULL){
> printf("Memory allocation error\n");
> exit(1);
> }
> for(i=0; i<=nbin; ++i){
> data[i]=0;
> }
>
> tot=0;
> if((f=fopen(args[1],"r"))==NULL){
> printf("Error: unable to open file on second attempt ...
> programming error?\n");
> exit(1);
> }
> while(fgets(linein,LINE_SIZE,f)!=NULL){
> if((sscanf(linein,"%lf",&c))!=1)continue;
> ++data[(int)ceil((c-min)/bin)];
> ++tot;
> }
> fclose(f);
>
> //Bin 0 only holds the absolute minimum value
> i=1;
>
> printf("%lf\t%lf\n",min+bin*((double)i-0.5),(double)(data[0]+data[1])/(double)tot);
>
> for(i=2; i<=nbin; ++i){
>
> printf("%lf\t%lf\n",min+bin*((double)i-0.5),(double)data[i]/(double)tot);
> }
> }
>
>
>
> -- original message --
>
> Hi Guys
>
> I am having great difficulty generating a histogram of the correct shape
> from g_wham. I am running umbrella sampling on one configuration fo two
> molecules at a fixed distance using the following
> pull parameters. ( I am just running one configuration to begin with so
> that can make sure that the shape of my distribution is correct).
>
> pull = umbrella
> pull_geometry = distance
> pull_dim = Y Y Y
> pull_start = yes
> pull_ngroups = 1
> pull_group0 = cage_1
> pull_group1 = cage_2
> pull_init1 = 0
> pull_rate1 = 0.0
> pull_k1 = 2000
> pull_nstxout = 100
> pull_nstfout = 100
>
> I have set the force constant quite high in the hope that the system
> will sample the equilibrium distance much more often than the distances
> accessible form the harmonic potential. In the plot of count vs distance
> I get an upside down bell curve in which there is a larger distribution
> at the extremeties of the curve. Does anyone know why this may be
> happening? If so please let me know. Am I right in assuming that the x
> axis in the histo.xvg file is the distance between the centres of masses
> of the two molecules?
>
> Cheers
>
> Gavin
>
>
More information about the gromacs.org_gmx-users
mailing list