Re: [gmx-users] Annealing of shell polarizable water model

2011-03-23 Thread Ivan Gladich

Dear Justin,
Dear all
 you guessed right.
I am now using a v-rescale thermostate
and it works fine and it also runs in parallel.
Here below the changed part in my grompp with the v-rescale thermostate.

Thanks a lot
Ivan

###
;OPTIONS FOR ANNELING
annealing = single
annealing_npoints = 2
annealing_time = 0 1000
annealing_temp = 0 160

; OPTIONS FOR WEAK COUPLING ALGORITHMS
; Temperature coupling
Tcoupl   = v-rescale ;Nose-hoover
; Groups to couple separately
tc-grps  = System
; Time constant (ps) and reference temperature (K)
tau_t= 0.1
ref_t= 160.00
##


On 03/22/2011 01:09 PM, Justin A. Lemkul wrote:



Ivan Gladich wrote:

Dear David
Dear all
I did the serial run with the same topology and grompp: even if the 
simulation time is still short (due by the serial run), the 
temperature profiles are the same (see attached file).


As further check, I removed the annealing and the temperature rises 
to 160 K after ~0.5ps without problem.
So I do not think that it is a problem of the shell polarizability 
with the running in parallel.

Could be a problem of the shell polarizability with annealing?



You could test that relatively easily with a simple box of water.  The 
result would be useful.


Other than that, maybe the thermostat itself is causing the problem?  
I always do annealing with a weak coupling method, not Nose-Hoover.  
I've had stability problems with unequilibrated systems with N-H.  In 
principle, it shouldn't matter, but setting tcoupl to Berendsen and/or 
V-rescale would be another very useful diagnostic.


-Justin


Thanks again
Ivan


On 03/22/2011 10:45 AM, Ivan Gladich wrote:

Yes, I am running in parallel...
Now I will try to run in serial to see if the problem persist
Thanks
Ivan

On 03/22/2011 10:00 AM, David van der Spoel wrote:

On 2011-03-22 10.37, Ivan Gladich wrote:

Dear all,
I would like to heat, very slowly, a ice box of 1796 SWM4-NDP water.
This kind of water has 4 sites plus a shell and I am using a small 
time

step (0.1 fms) to heat my ice box from 0 K top 160K in 1 ns.
To do that I used a linear annealing from 0 to 160 K.
Are you running in parallel? Unfortunately polarizable MD is broken 
on more than 1 core. There is a redmine issue for this, and it will 
be fixed soon.


The simulation runs without problem but I cannot reach the desire
temperature.
In other words, if I look my md.log file I can see the ref_t that 
linear
increase from 0 to 160 K in 1 ns but the system temperature seems 
to do

not follow the thermostate temperature.
If I plot the temperature obtained from g_energy, the temperature 
of the

system remains constant at ~36 K.
I attach also my temperature profile up to 600ps. Due to the small 
time

step the simulation takes a bit of time but it is clear that the
temperature remain constant


I have tried to find in the mail list some similar problem without
success...
Here below I report my grompp. Maybe I missed something.
Thank in advance for any suggestions.
Ivan
#3
; VARIOUS PREPROCESSING OPTIONS
title = Ice SWM4-NDP
cpp = /usr/bin/cpp
include =
define =

; RUN CONTROL PARAMETERS
integrator = md
dt = 0.0001
nsteps = 1400

; mode for center of mass motion removal
comm-mode = Linear
; number of steps for center of mass motion removal
nstcomm = 1
; group(s) for center of mass motion removal
comm-grps =

; OUTPUT CONTROL OPTIONS
; Output frequency for coords (x), velocities (v) and forces (f)
nstxout = 0
nstvout = 0
nstfout = 0
; Checkpointing helps you continue after crashes
nstcheckpoint = 1
; Output frequency for energies to log file and energy file
nstlog = 5000
nstenergy = 1000
; Output frequency and precision for xtc file
nstxtcout = 1
xtc-precision = 1000

; NEIGHBORSEARCHING PARAMETERS
; nblist update frequency
nstlist = 5
; ns algorithm (simple or grid)
ns_type = grid
; Periodic boundary conditions: xyz (default), no (vacuum)
; or full (infinite systems only)
pbc = xyz
; nblist cut-off
rlist = 1.1
;domain-decomposition =

; OPTIONS FOR ELECTROSTATICS AND VDW
; Method for doing electrostatics
coulombtype = PME
rcoulomb-switch = 0
rcoulomb = 1.1
; Method for doing Van der Waals
vdw-type = Cut-off
rvdw-switch = 0
rvdw = 1.1

; Apply long range dispersion corrections for Energy and Pressure
DispCorr = EnerPres

; Spacing for the PME/PPPM FFT grid
fourierspacing = 0.12
; FFT grid size, when a value is 0 fourierspacing will be used
fourier_nx = 0
fourier_ny = 0
fourier_nz = 0
; EWALD/PME/PPPM parameters
pme_order = 4
optimize_fft = no

;OPTIONS FOR ANNELING
annealing = single
annealing_npoints = 2
annealing_time = 0 1000
annealing_temp = 0 160

; OPTIONS FOR WEAK COUPLING ALGORITHMS
; Temperature coupling
Tcoupl = Nose-hoover
; Groups to couple separately
tc-grps = System
; Time constant (ps) and reference temperature (K)
tau_t = 0.1
ref_t = 160.00

[gmx-users] Re: Parameterization of a molecule containing a radical, (oxygen) (Justin A. Lemkul)

2011-03-23 Thread Gerrit Groenhof


There is a AMBER FF for this kind of nitroxide spin labels. See PCCP 12 
(2010) 11697.



Hope it helps,

Gerrit

Simone Cirri wrote:

Hi all,
I have a question regarding a new parameterization: actually, I'm
wondering whether or not it is possible to do it.

The molecule is tempol; you cand find its PubChem entry
here http://pubchem.ncbi.nlm.nih.gov/summary/summary.cgi?cid=137994loc=ec_rcs
http://pubchem.ncbi.nlm.nih.gov/summary/summary.cgi?cid=137994loc=ec_rcs
As you can see, the problem of this molecule is its O radical; this, and
the N belonging to the cycle too.
My lab would really benefit from the parameterization (for the AMBER
force field) of tempol, since we commonly use it in our NMR experiments.
We don't know how to treat the oxygen radical and the N bound to it, though.
If this parameterization is impossible to do we will give up, but we
would like to be sure about it.
Thanks in advance for any answer you could provide.
http://pubchem.ncbi.nlm.nih.gov/summary/summary.cgi?cid=137994loc=ec_rcs

This would not be a trivial molecule to parameterize.  You could start by using
antechamber (available free as part of AmberTools), since it allows you to tune
multiplicity (-m flag).  Otherwise, the existing force fields may not provide a
suitable atom type for the O radical or its effects on the rest of the molecule.
   I still don't know if parameterization would be possible for a standard MM
force field (though it very well might), but this might get you started.



-Justin


--

Simone Cirri



--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] how to differentiate crystal water and explicit waters in gromacs ???

2011-03-23 Thread Mark Abraham

On 23/03/2011 4:42 PM, Praveen Kumar Madala wrote:

Hi,

I am new to GROMACS and I am trying to perform a MD run in PBC 
condition with water as solvent.


I want to keep the crystal water molecules different from explicit 
waters,  I believe this helps me in analysis of crystal waters.


You can set up an index group to label them, and then use a 
visualization or other analysis program to see what they do. This 
doesn't require doing anything prior to the simulation. See manual and 
webpage for resources describing use and construction of groups.


Mark
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] single RNA nucleotide in CHARMM with phosphate problem with pdb2gmx

2011-03-23 Thread maria goranovic
Just an update, I decided to create a new residue for an isolated GMP. Well,
not really, because I used pdb2gmx to output a topology of a dimer, and then
modified it carefully to extract the topology of a phosphorylated GTP.

thank you !



On Tue, Mar 22, 2011 at 2:52 PM, maria goranovic
mariagorano...@gmail.comwrote:

 the gro file was the input. I am attaching a pdb as well

 CRYST1100.0100.0100.0  90.00  90.00  90.00 P 1   1
 ATOM  1  P   G   A   4   5.488  -0.842   0.708  1.00  0.00
   P
 ATOM  2  OP1 G   A   4   6.418   0.258   1.047  1.00  0.00
   O
 ATOM  3  OP2 G   A   4   5.367  -1.998   1.624  1.00  0.00
   O
 ATOM  4  O5' G   A   4   4.010  -0.213   0.513  1.00  0.00
   O
 ATOM  5  C5' G   A   4   3.816   1.098   0.014  1.00  0.00
   C
 ATOM  6  C4' G   A   4   2.324   1.432  -0.121  1.00  0.00
   C
 ATOM  7  O4' G   A   4   1.623   0.448  -0.862  1.00  0.00
   O
 ATOM  8  C3' G   A   4   1.584   1.482   1.211  1.00  0.00
   C
 ATOM  9  O3' G   A   4   1.869   2.669   1.934  1.00  0.00
   O
 ATOM 10  C2' G   A   4   0.136   1.413   0.721  1.00  0.00
   C
 ATOM 11  O2' G   A   4  -0.363   2.696   0.397  1.00  0.00
   O
 ATOM 12  C1' G   A   4   0.241   0.597  -0.574  1.00  0.00
   C
 ATOM 13  N9  G   A   4  -0.450  -0.705  -0.438  1.00  0.00
   N
 ATOM 14  C8  G   A   4   0.075  -1.954  -0.237  1.00  0.00
   C
 ATOM 15  N7  G   A   4  -0.808  -2.911  -0.197  1.00  0.00
   N
 ATOM 16  C5  G   A   4  -2.017  -2.251  -0.377  1.00  0.00
   C
 ATOM 17  C6  G   A   4  -3.346  -2.767  -0.450  1.00  0.00
   C
 ATOM 18  O6  G   A   4  -3.705  -3.941  -0.398  1.00  0.00
   O
 ATOM 19  N1  G   A   4  -4.299  -1.762  -0.593  1.00  0.00
   N
 ATOM 20  C2  G   A   4  -4.002  -0.413  -0.682  1.00  0.00
   C
 ATOM 21  N2  G   A   4  -5.023   0.438  -0.804  1.00  0.00
   N
 ATOM 22  N3  G   A   4  -2.751   0.071  -0.652  1.00  0.00
   N
 ATOM 23  C4  G   A   4  -1.809  -0.897  -0.499  1.00  0.00
   C
 ATOM 24  H5' G   A   4   4.293   1.204  -0.960  1.00  0.00
   H
 ATOM 25 H5'' G   A   4   4.269   1.819   0.696  1.00  0.00
   H
 ATOM 26  H4' G   A   4   2.217   2.395  -0.623  1.00  0.00
   H
 ATOM 27  H3' G   A   4   1.819   0.589   1.792  1.00  0.00
   H
 ATOM 28  H2' G   A   4  -0.503   0.940   1.465  1.00  0.00
   H
 ATOM 29 HO2' G   A   4  -1.301   2.614   0.210  1.00  0.00
   H
 ATOM 30  H1' G   A   4  -0.225   1.145  -1.393  1.00  0.00
   H
 ATOM 31  H8  G   A   4   1.129  -2.132  -0.139  1.00  0.00
   H
 ATOM 32  H1  G   A   4  -5.265  -2.065  -0.635  1.00  0.00
   H
 ATOM 33  H21 G   A   4  -5.974   0.099  -0.815  1.00  0.00
   H
 ATOM 34  H22 G   A   4  -4.841   1.428  -0.876  1.00  0.00
   H
 END


 On Tue, Mar 22, 2011 at 2:43 PM, Daniel Adriano Silva M 
 dadri...@gmail.com wrote:

 Hi,

 Is this your input? (a GRO file?). Can you paste the PDB corresponding to
 the RNA nucleotide that you are interested on? I think the file may be small
 enough to paste it here so someone can try to reproduce what is happening to
 you.

 Daniel


 2011/3/22 maria goranovic mariagorano...@gmail.com

 Sorry for the incomplete reply.
 version 4.5.3,
 command: pdb2gmx -f guan.pdb -ignh
 I tried using -ter, but that did not help.


 input coordinate file, pulled out from an RNA molecule.

 Gromacs Runs On Most of All Computer Systems
34
 1GP1   5.549   4.916   5.071
 1G  OP12   5.642   5.026   5.105
 1G  OP23   5.537   4.800   5.162
 1G  O5'4   5.401   4.979   5.051
 1G  C5'5   5.382   5.110   5.001
 1G  C4'6   5.232   5.143   4.988
 1G  O4'7   5.162   5.045   4.914
 1G  C3'8   5.158   5.148   5.121
 1G  O3'9   5.187   5.267   5.193
 1G  C2'   10   5.014   5.141   5.072
 1G  O2'   11   4.964   5.270   5.040
 1G  C1'   12   5.024   5.060   4.943
 1G   N9   13   4.955   4.930   4.956
 1G   C8   14   5.008   4.805   4.976
 1G   N7   15   4.919   4.709   4.980
 1G   C5   16   4.798   4.775   4.962
 1G   C6   17   4.665   4.723   4.955
 1G   O6   18   4.630   4.606   4.960
 1G   N1   19   4.570   4.824   4.941
  1G   C2   20   4.600   4.959   4.932
 1G   N2   21   4.498   5.044   4.920
 1G   N3   22   4.725   5.007   4.935
 1G   C4   23   4.819   4.910   4.950
 1G  H5'   24   5.429   5.120   4.904
 1G H5''   25   5.427   5.182   5.070
 1G  H4'   26   5.222   5.240   4.938
 1G  H3'   27   5.182   5.059   5.179
 1G  H2'   28   4.950   5.094   5.146
 1G HO2'   29   4.870   5.261   5.021
 1G  H1'   30   4.978   5.115   4.861
 1G   H8   31   5.113   4.787   

[gmx-users] angular velocity autocorrelation

2011-03-23 Thread Vigneshwar Ramakrishnan
Dear Users,

I need to calculate the angular velocity correlation function. From the
gromacs discussion forum, I understand that one should be able to use
g_rotacf to calculate the angular velocity correlation function (
http://oldwww.gromacs.org/pipermail/gmx-users/2004-September/012278.html).

However, when I look into the g_rotacf function, if I am not mistaken, the
rotational correlation function gives us the autocorrelation of the normal
vector to the two vectors defined in the input. It does not provide the
angular velocity correlation function.

Could someone give me pointers as to how one can calculate the angular
velocity correlation function from g_rotacf? I need this information to
analyze the density of states distribution.

Any help is greatly appreciated.

Thank you,
Sincerely,
Vignesh



-- 
R.Vigneshwar
Graduate Student,
Dept. of Chemical  Biomolecular Engg,
National University of Singapore,
Singapore

Strive for Excellence, Never be satisfied with the second Best!!

I arise in the morning torn between a desire to improve the world and a
desire to enjoy the world. This makes it hard to plan the day. (E.B. White)
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] Regarding the construction of Dynamics cross correlation map from covar.dat(from g_covar)

2011-03-23 Thread bipin singh
Hello users,

Although there was already much discussion had been done on the procedure to
construct DCCM map from covar.dat(i.e. output from g_covar) but still I am
not clear
with the step wise procedure to calculate the correlation matrix from
covarience matrix(covar.dat).

1) It is stated in the previous posts that to calculate correlation matrix
we have to
divide each i,j-th element of the covariance matrix (covar.dat) by the sqrt
of the ii-th and jj-th diagonal
element, but how to get diagonal of a 865107X3 matrix???
2) As given in the manual the order of elements in covar.dat is is:
x1x1,x1y1, x1z1, x1x2, ...
But I am not clear which are i,j-th ,ii-th and jj-th elements in covar.dat.

My covarience matrix have  865107 rows and 3 columns and the total atoms are
537(sqrt( 865107 X 3) / 3=537 atoms)
Can anyone suggest how to transform the covar.dat into correlation matrix.
I am expecting the stepwise procedure to construct the correlation matrix
and not the theory behind it as it already stated in the
previous posts.

-- 
---
*Regards,*
Bipin Singh
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] Various web issues

2011-03-23 Thread Rossen Apostolov

Hi Justin,

On 3/22/11 6:57 PM, Justin A. Lemkul wrote:
1. manual.gromacs.org no longer has any useful content.  The site 
comes up with syrah.cbr.su.se and nothing else.  Are there changes 
going on with this site? 

I moved the manual to the new server and now it's accessible again.

Cheers,
Rossen
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] %hbond number problem

2011-03-23 Thread gromacs564
Dear Justin
  I download the script plot_hbmap.pl and input perl $0 -s b4md.pdb -map 
hbmap.xpm -index hbond.ndx\n,

but the output info is Unrecognized switch: -bash  (-h will show valid 
options).;

  then,I input perl\n plot_hbmap.pl -s b4md.pdb -map hbmap.xpm -index 
hbond.ndx\n,and the output info is nothing?

 So, could you tell me the correct usage? Very thanks!

There are protein and sol in b4md.pdb file.The hbond.ndx file only contain the 
[ hbonds_Protein ] section:
  1  2603
  1  2619
  1  2634
  1  2635
  1  2   1299
   1942   1943   1917
   1937   1939113
   1937   1939114

-

 xiaodu



My simulation system contains protein, dna and water.
I used your script already for obtaining %exist hydrogen bonds between
protein and dna-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] LINCS warning on galactose molecule

2011-03-23 Thread Anna Marabotti
Dear gmx-users,
I'm experimenting a LINCS warning on my system, an enzyme with ATP,
galactose and Mg. We started from the crystallographic structure, made some
accommodation on the topology of the sugar in order to check for partial
charges (according to Justin Lemkul's suggestions on PRODRG server), then
performed a minimization, NVT and NPT PR-MD and the production MD. After
about 5 ns of simulation the MD stopped with the LINCS error:
 
Step 2562553, time 5125.11 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000102, max 0.001903 (between atoms 7643 and 7642)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 45.3 0.1000 0.1000 0.1000

Step 2562554, time 5125.11 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000104, max 0.001558 (between atoms 7643 and 7642)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 36.5 0.1000 0.1000 0.1000

Step 2562565, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000318, max 0.007712 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 31.5 0.1005 0.1008 0.1000

Step 2562566, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000374, max 0.009231 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 38.2 0.1008 0.1009 0.1000

Step 2562567, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000310, max 0.007658 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 38.6 0.1009 0.1008 0.1000



Step 2563563, time 5127.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.70, max 0.001640 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 57.0 0.0998 0.1002 0.1000

---

Program mdrun, VERSION 4.0.3

Source code file: constr.c, line: 136

Fatal error:

Too many LINCS warnings (1000)

If you know what you are doing you can adjust the lincs warning threshold in
your mdp file

or set the environment variable GMX_MAXCONSTRWARN to -1,

but normally it is better to fix the problem

---

We checked the system to find the problem. The atoms 7641 and 7640
correspond to atom O1 and H12 of the galactose. With VMD we visualized the
trajectory, but we did not see anything anomalous on this bond. Length and
angles are quite normal. This OH group is near the phosphate of ATP but the
distances seem to be compatible with an (expected) H-bond. I don't think the
problem could be an incorrect stabilization because it seems to me that in
that case the problem would have appear during the first step of MD, not
after 5 ns. In any case, stabilization went correctly, with no errors or
warnings. I checked on the gmx-users list and found some messages of Xavier
Periole experimenting some LINCS problems with gmx 4.0.4, that seems to
resemble my situation. I followed the thread, but I did not see anything
that could be of help for me. Some hypotheses were about the gmx version
(mine is 4.0.3), communications between CPUs, virtual sites. Looking at the
trajectory, I don't think it's a problem of periodic images too close. The
thread stops with Xavier saying he would like to test the effects of using
-rdd and -rcon on LINCS warning, but I did not find the rest of the history.

 
Does anybody (and especially Xavier, if possible!) have some hints about my
problem? Or could it be a problem of galactose topology?
 
Thank you very much and best regards
Anna 
__
Anna Marabotti, Ph.D.
Laboratory of Bioinformatics and Computational Biology
Institute of Food Science - CNR
Via Roma, 64
83100 Avellino
Phone: +39 0825 299651
Fax: +39 0825 781585
E-mail: amarabo...@isa.cnr.it
Skype account: annam1972
Web site: http://bioinformatica.isa.cnr.it/anna/anna.htm
 
When a man with a gun meets a man with a pen, the man with the gun is a
dead man
 
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] %hbond number problem

2011-03-23 Thread gromacs564
 
 
 

Dear Justin
  I download the script plot_hbmap.pl and input perl $0 -s b4md.pdb -map 
hbmap.xpm -index hbond.ndx\n,

but the output info is Unrecognized switch: -bash  (-h will show valid 
options).;

  then,I input perl\n plot_hbmap.pl -s b4md.pdb -map hbmap.xpm -index 
hbond.ndx\n,and the output info is nothing?

 So, could you tell me the correct usage? Very thanks!

There are protein and sol in b4md.pdb file.The hbond.ndx file only contain the 
[ hbonds_Protein ] section:
  1  2603
  1  2619
  1  2634
  1  2635
  1  2   1299
   1942   1943   1917
   1937   1939113
   1937   1939114

-

 xiaodu




 -- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] LINCS warning on galactose molecule

2011-03-23 Thread Mark Abraham

On 23/03/2011 9:45 PM, Anna Marabotti wrote:

Dear gmx-users,
I'm experimenting a LINCS warning on my system, an enzyme with ATP, 
galactose and Mg. We started from the crystallographic structure, made 
some accommodation on the topology of the sugar in order to check for 
partial charges (according to Justin Lemkul's suggestions on PRODRG 
server), then performed a minimization, NVT and NPT PR-MD and the 
production MD. After about 5 ns of simulation the MD stopped with the 
LINCS error:


Step 2562553, time 5125.11 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000102, max 0.001903 (between atoms 7643 and 7642)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 45.3 0.1000 0.1000 0.1000

Step 2562554, time 5125.11 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000104, max 0.001558 (between atoms 7643 and 7642)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 36.5 0.1000 0.1000 0.1000

Step 2562565, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000318, max 0.007712 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 31.5 0.1005 0.1008 0.1000

Step 2562566, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000374, max 0.009231 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 38.2 0.1008 0.1009 0.1000

Step 2562567, time 5125.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.000310, max 0.007658 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 38.6 0.1009 0.1008 0.1000



Step 2563563, time 5127.13 (ps) LINCS WARNING

relative constraint deviation after LINCS:

rms 0.70, max 0.001640 (between atoms 7641 and 7640)

bonds that rotated more than 30 degrees:

atom 1 atom 2 angle previous, current, constraint length

7641 7640 57.0 0.0998 0.1002 0.1000

---

Program mdrun, VERSION 4.0.3

Source code file: constr.c, line: 136

Fatal error:

Too many LINCS warnings (1000)

If you know what you are doing you can adjust the lincs warning 
threshold in your mdp file


or set the environment variable GMX_MAXCONSTRWARN to -1,

but normally it is better to fix the problem

---

We checked the system to find the problem. The atoms 7641 and 7640 
correspond to atom O1 and H12 of the galactose. With VMD we visualized 
the trajectory, but we did not see anything anomalous on this bond. 
Length and angles are quite normal. This OH group is near the 
phosphate of ATP but the distances seem to be compatible with an 
(expected) H-bond. I don't think the problem could be an incorrect 
stabilization because it seems to me that in that case the problem 
would have appear during the first step of MD, not after 5 ns. In any 
case, stabilization went correctly, with no errors or warnings. I 
checked on the gmx-users list and found some messages of Xavier 
Periole experimenting some LINCS problems with gmx 4.0.4, that seems 
to resemble my situation. I followed the thread, but I did not see 
anything that could be of help for me. Some hypotheses were about the 
gmx version (mine is 4.0.3), communications between CPUs, virtual 
sites. Looking at the trajectory, I don't think it's a problem of 
periodic images too close. The thread stops with Xavier saying he 
would like to test the effects of using -rdd and -rcon on LINCS 
warning, but I did not find the rest of the history.
Does anybody (and especially Xavier, if possible!) have some hints 
about my problem? Or could it be a problem of galactose topology?


Usually the problem has nothing to do with LINCS or -rdd or -con - 
they're just symptoms of machinery that breaks first.


I'd start by assuming my topologies were not well-formed somehow, and 
divide and conquer to locate the problem. Can you simulate galactose in 
vacuo? In a small water box? Same for ATP. Same for enzyme.


Mark
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] unit of surface tension, KJ/mol?

2011-03-23 Thread Mark Abraham

On 23/03/2011 3:39 PM, Elisabeth wrote:
Thanks Mark for the hints. I am using a pre equilibrated box for both 
runs. Can you again explain how Can I interpret this from fluctuations 
from two runs? Below is the resutls from two versions: Thanks!


IIRC the contents of the .edr file didn't change across these versions, 
but the reporting of computed quantities did, so I'd make my life 
simpler and use only the 4.5.3 g_energy.


Both simulations have a much higher value for the root-mean-squared 
deviation from the average (i.e. standard deviation) than for the 
average. So that means your data set looks something like 1 20 6543 
200 23444 23434 15 500 3444. Look at the .xvg file that is produced. 
You need a heck of a lot of such data to have confidence that your 
sample average reflects the true average. If you haven't got that pile 
of data, then you haven't observed the average with any confidence. This 
is normal for quantities computed from fluctuations in atomic positions. 
They're macroscopic quantities, and have to be observed over a lot of 
microscopic configurations.


Be sure to check visually that your system is doing what you hope it is 
doing.


Mark



4.0.7

Statistics over 733801 steps [ 0. thru 1467.6001 ps ], 1 data sets
All averages are exact over 733801 steps

Energy  Average   RMSD Fluct.  Drift  
Tot-Drift

---
#Surf*SurfTen   203.7783623.993623.98 
-0.00203751   -2.99026



Statistics over 483401 steps [ 500. thru 1466.8000 ps ], 1 data sets
All averages are exact over 483401 steps

Energy  Average   RMSD Fluct.  Drift  
Tot-Drift

---
#Surf*SurfTen   203.1293612.583612.28   
0.168686163.086

-
4.5.3



Statistics over 732501 steps [ 1. through 1466. ps ], 1 data sets
All statistics are over 73251 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   2374.4957021081.3   -3251.38  
(bar nm)


Statistics over 384501 steps [ 1000. through 1769. ps ], 1 
data sets

All statistics are over 38451 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   2610.25   100021109.94638.38  
(bar nm)


Statistics over 483001 steps [ 500. through 1466. ps ], 1 data 
sets

All statistics are over 48301 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   1645.0739020288.3   -174.289  
(bar nm)






On 22 March 2011 23:56, Mark Abraham mark.abra...@anu.edu.au 
mailto:mark.abra...@anu.edu.au wrote:




On 23/03/11, *Elisabeth * katesed...@gmail.com
mailto:katesed...@gmail.com wrote:



On 22 March 2011 22:46, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu wrote:



Elisabeth wrote:



On 22 March 2011 22:31, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu mailto:jalem...@vt.edu
mailto:jalem...@vt.edu wrote:



   Elisabeth wrote:

   Hello,
I did two simulations on the same system using
versions 4.0.7
   and 4.5.3. It seems like the unit of surface
tension is not the
   same in these versions because I am getting ~250
KJ/mol an in
   4.0.7 and ~ 5000bar nm in 4.5.3?! How KJ/mol can
be converted
   into bar nm? Can anyone help please.


   Is this from g_energy output?  In past versions,
everything was
   printed as kJ/mol, even quantities that obviously
weren't, like
   temperature, pressure, etc.


Yes, so why are the results so different. I am using
exactly the same mdp file.!


Any pressure-related quantity is going to be subject to
enormous fluctuations. This has been discussed within the
last few days.  Without seeing the .mdp file and a
description of the system, it's hard for anyone to comment on
what the results might mean.


Thanks. I am working on a pure alkane system, in a box of 3X3X3
which is extended in Z to create liq/air interface. That is 3X3X6
. What could be reason for such a big difference in results from
two version? Thanks alot!


Look at the size of the fluctuations printed by 

Re: [gmx-users] unit of surface tension, KJ/mol?

2011-03-23 Thread Justin A. Lemkul



Mark Abraham wrote:

  On 23/03/2011 3:39 PM, Elisabeth wrote:
Thanks Mark for the hints. I am using a pre equilibrated box for both 
runs. Can you again explain how Can I interpret this from fluctuations 
from two runs? Below is the resutls from two versions: Thanks!


IIRC the contents of the .edr file didn't change across these versions, 
but the reporting of computed quantities did, so I'd make my life 
simpler and use only the 4.5.3 g_energy.




Actually, 4.5.3 may have the bug:

http://redmine.gromacs.org/issues/696

The quantities in the .edr file are correct, but the output from g_energy is 
nonsense.  It was fixed for 4.5.4.


-Justin

Both simulations have a much higher value for the root-mean-squared 
deviation from the average (i.e. standard deviation) than for the 
average. So that means your data set looks something like 1 20 6543 
200 23444 23434 15 500 3444. Look at the .xvg file that is produced. 
You need a heck of a lot of such data to have confidence that your 
sample average reflects the true average. If you haven't got that pile 
of data, then you haven't observed the average with any confidence. This 
is normal for quantities computed from fluctuations in atomic positions. 
They're macroscopic quantities, and have to be observed over a lot of 
microscopic configurations.


Be sure to check visually that your system is doing what you hope it is 
doing.


Mark



4.0.7

Statistics over 733801 steps [ 0. thru 1467.6001 ps ], 1 data sets
All averages are exact over 733801 steps

Energy  Average   RMSD Fluct.  Drift  
Tot-Drift

---
#Surf*SurfTen   203.7783623.993623.98 
-0.00203751   -2.99026



Statistics over 483401 steps [ 500. thru 1466.8000 ps ], 1 data sets
All averages are exact over 483401 steps

Energy  Average   RMSD Fluct.  Drift  
Tot-Drift

---
#Surf*SurfTen   203.1293612.583612.28   
0.168686163.086

-
4.5.3



Statistics over 732501 steps [ 1. through 1466. ps ], 1 data sets
All statistics are over 73251 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   2374.4957021081.3   -3251.38  
(bar nm)


Statistics over 384501 steps [ 1000. through 1769. ps ], 1 
data sets

All statistics are over 38451 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   2610.25   100021109.94638.38  
(bar nm)


Statistics over 483001 steps [ 500. through 1466. ps ], 1 data 
sets

All statistics are over 48301 points

Energy  Average   Err.Est.   RMSD  Tot-Drift
---
#Surf*SurfTen   1645.0739020288.3   -174.289  
(bar nm)






On 22 March 2011 23:56, Mark Abraham mark.abra...@anu.edu.au 
mailto:mark.abra...@anu.edu.au wrote:




On 23/03/11, *Elisabeth * katesed...@gmail.com
mailto:katesed...@gmail.com wrote:



On 22 March 2011 22:46, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu wrote:



Elisabeth wrote:



On 22 March 2011 22:31, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu mailto:jalem...@vt.edu
mailto:jalem...@vt.edu wrote:



   Elisabeth wrote:

   Hello,
I did two simulations on the same system using
versions 4.0.7
   and 4.5.3. It seems like the unit of surface
tension is not the
   same in these versions because I am getting ~250
KJ/mol an in
   4.0.7 and ~ 5000bar nm in 4.5.3?! How KJ/mol can
be converted
   into bar nm? Can anyone help please.


   Is this from g_energy output?  In past versions,
everything was
   printed as kJ/mol, even quantities that obviously
weren't, like
   temperature, pressure, etc.


Yes, so why are the results so different. I am using
exactly the same mdp file.!


Any pressure-related quantity is going to be subject to
enormous fluctuations. This has been discussed within the
last few days.  Without seeing the .mdp file and a
description of the system, it's hard for anyone to comment on
what the results might mean.


Thanks. I am working on a pure alkane system, in a box of 3X3X3

Re: [gmx-users] %hbond number problem

2011-03-23 Thread Justin A. Lemkul



gromacs564 wrote:

Dear Justin
  I download the script plot_hbmap.pl and input perl $0 -s 
b4md.pdb -map hbmap.xpm -index hbond.ndx\n,


but the output info is Unrecognized switch: -bash  (-h will show valid 
options).;


  then,I input perl\n plot_hbmap.pl -s b4md.pdb -map hbmap.xpm 
-index hbond.ndx\n,and the output info is nothing?


 So, could you tell me the correct usage? Very thanks!

There are protein and sol in b4md.pdb file.The hbond.ndx file only 
contain the [ hbonds_Protein ] section:

  1  2603
  1  2619
  1  2634
  1  2635
  1  2   1299
   1942   1943   1917
   1937   1939113
   1937   1939114

-

 xiaodu



My simulation system contains protein, dna and water.
I used your script already for obtaining %exist hydrogen bonds between
protein and dna



Then doesn't this imply that you know how to run it?  Neither of the above 
commands is correct.  The first uses $0, which is a perl construct for the name 
of the program.  The second contains a newline (\n) character, although the 
format of this command is what you should use (without the newline, which is 
printed by the program to the screen as help information).


-Justin

--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] unit of surface tension, KJ/mol?

2011-03-23 Thread Elisabeth
Thanks for your help. :)

On 23 March 2011 07:47, Justin A. Lemkul jalem...@vt.edu wrote:



 Mark Abraham wrote:

  On 23/03/2011 3:39 PM, Elisabeth wrote:

 Thanks Mark for the hints. I am using a pre equilibrated box for both
 runs. Can you again explain how Can I interpret this from fluctuations from
 two runs? Below is the resutls from two versions: Thanks!


 IIRC the contents of the .edr file didn't change across these versions,
 but the reporting of computed quantities did, so I'd make my life simpler
 and use only the 4.5.3 g_energy.


 Actually, 4.5.3 may have the bug:

 http://redmine.gromacs.org/issues/696

 The quantities in the .edr file are correct, but the output from g_energy
 is nonsense.  It was fixed for 4.5.4.

 -Justin

  Both simulations have a much higher value for the root-mean-squared
 deviation from the average (i.e. standard deviation) than for the average.
 So that means your data set looks something like 1 20 6543 200 23444
 23434 15 500 3444. Look at the .xvg file that is produced. You need a heck
 of a lot of such data to have confidence that your sample average reflects
 the true average. If you haven't got that pile of data, then you haven't
 observed the average with any confidence. This is normal for quantities
 computed from fluctuations in atomic positions. They're macroscopic
 quantities, and have to be observed over a lot of microscopic
 configurations.

 Be sure to check visually that your system is doing what you hope it is
 doing.

 Mark


 4.0.7

 Statistics over 733801 steps [ 0. thru 1467.6001 ps ], 1 data sets
 All averages are exact over 733801 steps

 Energy  Average   RMSD Fluct.  Drift
  Tot-Drift

 ---
 #Surf*SurfTen   203.7783623.993623.98 -0.00203751
 -2.99026


 Statistics over 483401 steps [ 500. thru 1466.8000 ps ], 1 data sets
 All averages are exact over 483401 steps

 Energy  Average   RMSD Fluct.  Drift
  Tot-Drift

 ---
 #Surf*SurfTen   203.1293612.583612.28   0.168686
  163.086

 -
 4.5.3



 Statistics over 732501 steps [ 1. through 1466. ps ], 1 data sets
 All statistics are over 73251 points

 Energy  Average   Err.Est.   RMSD  Tot-Drift

 ---
 #Surf*SurfTen   2374.4957021081.3   -3251.38
  (bar nm)

 Statistics over 384501 steps [ 1000. through 1769. ps ], 1 data
 sets
 All statistics are over 38451 points

 Energy  Average   Err.Est.   RMSD  Tot-Drift

 ---
 #Surf*SurfTen   2610.25   100021109.94638.38
  (bar nm)

 Statistics over 483001 steps [ 500. through 1466. ps ], 1 data
 sets
 All statistics are over 48301 points

 Energy  Average   Err.Est.   RMSD  Tot-Drift

 ---
 #Surf*SurfTen   1645.0739020288.3   -174.289
  (bar nm)





 On 22 March 2011 23:56, Mark Abraham mark.abra...@anu.edu.au mailto:
 mark.abra...@anu.edu.au wrote:



On 23/03/11, *Elisabeth * katesed...@gmail.com
mailto:katesed...@gmail.com wrote:



On 22 March 2011 22:46, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu wrote:



Elisabeth wrote:



On 22 March 2011 22:31, Justin A. Lemkul jalem...@vt.edu
mailto:jalem...@vt.edu mailto:jalem...@vt.edu

mailto:jalem...@vt.edu wrote:



   Elisabeth wrote:

   Hello,
I did two simulations on the same system using
versions 4.0.7
   and 4.5.3. It seems like the unit of surface
tension is not the
   same in these versions because I am getting ~250
KJ/mol an in
   4.0.7 and ~ 5000bar nm in 4.5.3?! How KJ/mol can
be converted
   into bar nm? Can anyone help please.


   Is this from g_energy output?  In past versions,
everything was
   printed as kJ/mol, even quantities that obviously
weren't, like
   temperature, pressure, etc.


Yes, so why are the results so different. I am using
exactly the same mdp file.!


Any pressure-related quantity is going to be subject to
enormous fluctuations. This has been discussed within the
last few days.  Without seeing the .mdp file and a
description of the system, it's hard for anyone to 

[gmx-users] no output dgdl file

2011-03-23 Thread Moeed
Hello,

I have a little problem with FE output file. Below is the settings and also
I am including -dgdl in the command I issue but no dgdl (or dhdl) file
generates. I dont figure where the problem lies ! (version 4.5.3).

free_energy  =   yes
init_lambda  =   0
delta_lambda =   0
sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q
couple-lambda1   =   none
couple-intramol  =   yes
nstdhdl  =   10
separate_dhdl_file   =   yes
dh_hist_size =   0
dh_hist_spacing  =   0.1

Best
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Moeed wrote:

Hello,

I have a little problem with FE output file. Below is the settings and 
also I am including -dgdl in the command I issue but no dgdl (or dhdl) 
file generates. I dont figure where the problem lies ! (version 4.5.3).




How long has the simulation been running?  Output is buffered, and I find that 
dhdl.xvg is not updated often.



free_energy  =   yes
init_lambda  =   0 
delta_lambda =   0

sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q   
couple-lambda1   =   none


Note that turning off all interactions simultaneously (especially when using 
soft-core) can be very unstable.  You should (de)couple vdW and Coulombic 
interactions separately.


-Justin


couple-intramol  =   yes
nstdhdl  =   10
separate_dhdl_file   =   yes
dh_hist_size =   0
dh_hist_spacing  =   0.1

Best



--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] compressing a box of water droplets into a homogeneous solution of liquid water

2011-03-23 Thread Patrick Fuchs

Hi Chris,
I experienced the same kind of thing. In the process of building a box 
of liquid (organic solvent), at some point I wanted to get rid of a 
layer of vacuum around my system. So for shrinking the box I used 
similar settings as you and found also that the collapse was going 
slower than I'd have expected.
One solution to accelerate this (if your goal is to shrink the box) is 
to increase the pressure (to say 100 atm). But it's important to stop 
the simulation in time (i.e. once the layer of vacuum has disapeared) 
otherwise the system shrinks too much and density is off.
So to come back to your system which has a very big layer of vacuum 
around, and according to my experience, the volume is probably 
decreasing but too slowly to see anything signigicant (compared to the 
initial value) in 200 ps .

Ciao,

Patrick

Le 21/03/2011 16:53, chris.ne...@utoronto.ca a écrit :

Dear users:

I recently came across a system that was composed of tip4p water vapor
droplets separated by vacuum. This system is what you might get if you
did a NVT simulation of water with a box that was 10 times too large for
the number of water molecules.

I was surprised to see that this system did not collapse to any
significant extent during 200 ps of NPT equilibration at 1 atm using the
Berendsen thermostat with tau_p=1.0 and the sd integrator and a colombic
cut-off. (We also tried a number of other integrator/pressure coupling
combinations with the same results).

I had assumed that such collapse would occur quite rapidly. This does
not seem to be the case (no noticeable contraction within 200 ps).

Has anybody else done anything like this? Can anybody comment on their
expectations/experience of collapse from the gas state to the liquid
state under standard NPT conditions?

Thank you,
Chris.



--
___
 new E-mail address: patrick.fu...@univ-paris-diderot.fr 
Patrick FUCHS
Dynamique des Structures et Interactions des Macromolécules Biologiques
INTS, INSERM UMR-S665, Université Paris Diderot,
6 rue Alexandre Cabanel, 75015 Paris
Tel : +33 (0)1-44-49-30-57 - Fax : +33 (0)1-47-34-74-31
Web Site: http://www.dsimb.inserm.fr/~fuchs
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] compressing a box of water droplets into a homogeneous solution of liquid water

2011-03-23 Thread André Farias de Moura
Hi Chris,

recently one of my students made a mistake during the model system assembling,
setting the initial volume to a value three times larger than expected
by the density.

after 1 microsecond (coarse grained MD, of course) there were droplets
and vacuum
between them, and the volume did not shrink until we applied a 100 bar
pressure to
the external bath (using Berendsen's weak coupling with the default
setting for the
Martini FF).

on the other hand, recently I was trying to model the methanol gas
phase using the
OPLS-AA FF and I found out it was a very tricky kind of model, as
different systems
with different volumes but the same (gas phase) density yielded quite
different results,
ranging from only one droplet (smaller system) to a distribution with
more than 99% of
monomers (larger system). from my standpoint, what happened to the
smaller systems
should be regarded as an artefact due to finite size of the system,
thus only the results
for the larger system got published, as I was not interested in the
discussion of the
artefact itself.

best regards,

André


On Wed, Mar 23, 2011 at 3:48 PM, Patrick Fuchs
patrick.fu...@univ-paris-diderot.fr wrote:
 Hi Chris,
 I experienced the same kind of thing. In the process of building a box of
 liquid (organic solvent), at some point I wanted to get rid of a layer of
 vacuum around my system. So for shrinking the box I used similar settings as
 you and found also that the collapse was going slower than I'd have
 expected.
 One solution to accelerate this (if your goal is to shrink the box) is to
 increase the pressure (to say 100 atm). But it's important to stop the
 simulation in time (i.e. once the layer of vacuum has disapeared) otherwise
 the system shrinks too much and density is off.
 So to come back to your system which has a very big layer of vacuum around,
 and according to my experience, the volume is probably decreasing but too
 slowly to see anything signigicant (compared to the initial value) in 200 ps
 .
 Ciao,

 Patrick

 Le 21/03/2011 16:53, chris.ne...@utoronto.ca a écrit :

 Dear users:

 I recently came across a system that was composed of tip4p water vapor
 droplets separated by vacuum. This system is what you might get if you
 did a NVT simulation of water with a box that was 10 times too large for
 the number of water molecules.

 I was surprised to see that this system did not collapse to any
 significant extent during 200 ps of NPT equilibration at 1 atm using the
 Berendsen thermostat with tau_p=1.0 and the sd integrator and a colombic
 cut-off. (We also tried a number of other integrator/pressure coupling
 combinations with the same results).

 I had assumed that such collapse would occur quite rapidly. This does
 not seem to be the case (no noticeable contraction within 200 ps).

 Has anybody else done anything like this? Can anybody comment on their
 expectations/experience of collapse from the gas state to the liquid
 state under standard NPT conditions?

 Thank you,
 Chris.


 --
 ___
  new E-mail address: patrick.fu...@univ-paris-diderot.fr 
 Patrick FUCHS
 Dynamique des Structures et Interactions des Macromolécules Biologiques
 INTS, INSERM UMR-S665, Université Paris Diderot,
 6 rue Alexandre Cabanel, 75015 Paris
 Tel : +33 (0)1-44-49-30-57 - Fax : +33 (0)1-47-34-74-31
 Web Site: http://www.dsimb.inserm.fr/~fuchs
 --
 gmx-users mailing list    gmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
 Please don't post (un)subscribe requests to the list. Use the www interface
 or send it to gmx-users-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] no output dgdl file

2011-03-23 Thread Moeed
 Hi Justin,


I was looking all the time for a file with dgdl extension throughout the
directory. IIRC in previous versions -dgdl would generate such extension. I
found now .xvg file which has dhdl data. (if I am right dgdl is the very
dhdl in 4.5.3!)



 I have a little problem with FE output file. Below is the settings and also
 I am including -dgdl in the command I issue but no dgdl (or dhdl) file
 generates. I dont figure where the problem lies ! (version 4.5.3).


 How long has the simulation been running?  Output is buffered, and I find
 that dhdl.xvg is not updated often.


  free_energy  =   yes
 init_lambda  =   0 delta_lambda =   0
 sc_alpha =   0.5
 sc-power =   1
 sc_sigma =   0.3
 foreign_lambda   =   0.1
 dhdl_derivatives =   yes
 couple-moltype   =   Polymer
 couple-lambda0   =   vdw-q   couple-lambda1   =   none


 Note that turning off all interactions simultaneously (especially when
 using soft-core) can be very unstable.  You should (de)couple vdW and
 Coulombic interactions separately.


Yes I read your comment about decoupling separately earlier on but I was not
sure if this could depend upon the system so I was thinking of comparing
results from two trials:

trial 1:
couple-lambda0   =   vdw-q
couple-lambda1   =   none

 trail 2 A:
couple-lambda0   =   vdw-q
couple-lambda1   =   vdw
 trail 2 B
couple-lambda0   =   vdw
couple-lambda1   =   none

So:

1- Does this decoupling issue depend on the system or it is a must for all
FE studies?
2- I actually went through many of FE posts and never found how separete
decoupling can be done. from  trail 2 A and B I get dgdl files but how to
get total vdw-q to none from these two? Please shed some light on the
details of procedure..

3- With free_energy turned off my runs take only a few hours but as I shift
to free_energy  =   yes, runs take nearly one day foe each lambda..
Is this normal for s system having 2500 atoms?

 Thanks



 Best


 --
 

 Justin A. Lemkul
 Ph.D. Candidate
 ICTAS Doctoral Scholar
 MILES-IGERT Trainee
 Department of Biochemistry
 Virginia Tech
 Blacksburg, VA
 jalemkul[at]vt.edu | (540) 231-9080
 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

 
 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
 Please don't post (un)subscribe requests to the list. Use the www interface
 or send it to gmx-users-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

[gmx-users] compressing a box of water droplets into a homogeneous solution of liquid water

2011-03-23 Thread chris . neale

Thanks Patrick and Andre!

We repeated this with a few box sizes just to get a quick handle on  
it. The equilibrium volume is about 64 nm^3. If we start with a volume  
of 1000 nm^3 then the overall box does not collapse at all within 200  
ps of NPT Langevin dynamics at 1 atm.If we start with a volume of 200  
nm^3, then it does collapse to approximately 64 nm^3 within 200 ps of  
such simulation.


My best guess is that the rapid collapse is driven by nonbonded  
interactions and thus the rapid collapse does not occur when the  
system is so large with such low density that water forms isolated  
vapour droplets that do not interact with each other by LJ  
interactions. Sure, it is expected to collapse eventually from the 1  
atm pressure coupling, and we have also observed that high pressure  
works, but at 1 atm it might take a very long time to reach equilibrium.


I agree with Andre that none of this matters to regular simulations as  
there is no good reason to go through this type of state when one  
wants to simulate dense liquids. I just found it curious that  
Berendsen pressure coupling at 1 atm was not sufficient to quickly  
equilibrate the volume in a system where the vacuum regions are large  
in comparison to the LJ cutoffs.


Chris.

-- original message --

Hi Chris,
I experienced the same kind of thing. In the process of building a box
of liquid (organic solvent), at some point I wanted to get rid of a
layer of vacuum around my system. So for shrinking the box I used
similar settings as you and found also that the collapse was going
slower than I'd have expected.
One solution to accelerate this (if your goal is to shrink the box) is
to increase the pressure (to say 100 atm). But it's important to stop
the simulation in time (i.e. once the layer of vacuum has disapeared)
otherwise the system shrinks too much and density is off.
So to come back to your system which has a very big layer of vacuum
around, and according to my experience, the volume is probably
decreasing but too slowly to see anything signigicant (compared to the
initial value) in 200 ps .
Ciao,

Patrick

Le 21/03/2011 16:53, chris.ne...@utoronto.ca a écrit :

[Hide Quoted Text]
Dear users:

I recently came across a system that was composed of tip4p water vapor
droplets separated by vacuum. This system is what you might get if you
did a NVT simulation of water with a box that was 10 times too large for
the number of water molecules.

I was surprised to see that this system did not collapse to any
significant extent during 200 ps of NPT equilibration at 1 atm using the
Berendsen thermostat with tau_p=1.0 and the sd integrator and a colombic
cut-off. (We also tried a number of other integrator/pressure coupling
combinations with the same results).

I had assumed that such collapse would occur quite rapidly. This does
not seem to be the case (no noticeable contraction within 200 ps).

Has anybody else done anything like this? Can anybody comment on their
expectations/experience of collapse from the gas state to the liquid
state under standard NPT conditions?

Thank you,
Chris.

--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] account for PBC

2011-03-23 Thread Justin A. Lemkul



ahmet yıldırım wrote:

Dear Justin,

You said You didn't properly account for PBC.  Watch the trajectory - 
you should see distinct jumps across the box. for attached figure.


What should I do?



Exactly as I said.  Watch the trajectory to see if my guess is right.  If the 
molecule is dancing back and forth across PBC, it explains what you're seeing 
with the RMSD.


Fix the PBC issues with trjconv, if that is indeed the case.

http://www.gromacs.org/Documentation/Terminology/Periodic_Boundary_Conditions

-Justin


--
Ahmet YILDIRIM





--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Moeed wrote:


Hi Justin,


I was looking all the time for a file with dgdl extension throughout the 
directory. IIRC in previous versions -dgdl would generate such 
extension. I found now .xvg file which has dhdl data. (if I am right 
dgdl is the very dhdl in 4.5.3!) 



The extension has always been .xvg, but the root filename has changed.  It is 
now perhaps a bit more accurate (delta Hamiltonian).


 


I have a little problem with FE output file. Below is the
settings and also I am including -dgdl in the command I issue
but no dgdl (or dhdl) file generates. I dont figure where the
problem lies ! (version 4.5.3).


How long has the simulation been running?  Output is buffered, and I
find that dhdl.xvg is not updated often.


free_energy  =   yes
init_lambda  =   0 delta_lambda =   0
sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q   couple-lambda1   =   none


Note that turning off all interactions simultaneously (especially
when using soft-core) can be very unstable.  You should (de)couple
vdW and Coulombic interactions separately.


Yes I read your comment about decoupling separately earlier on but I was 
not sure if this could depend upon the system so I was thinking of 
comparing results from two trials:


trial 1:
couple-lambda0   =   vdw-q
couple-lambda1   =   none 


 trail 2 A:
couple-lambda0   =   vdw-q
couple-lambda1   =   vdw

 trail 2 B
couple-lambda0   =   vdw
couple-lambda1   =   none

So:

1- Does this decoupling issue depend on the system or it is a must for 
all FE studies?


No, it pertains to the very stability of all calculations.  Whenever I've 
decoupled charges and LJ simultaneously, either the simulations crash or the 
resulting energies are not reliable.  There are numerous book chapters and 
papers dedicated to these topics.


2- I actually went through many of FE posts and never found how separete 
decoupling can be done. from  trail 2 A and B I get dgdl files but how 
to get total vdw-q to none from these two? Please shed some light on the 
details of procedure..




Add them.  DeltaG is a state function.

3- With free_energy turned off my runs take only a few hours but as I 
shift to free_energy  =   yes, runs take nearly one day foe each 
lambda.. Is this normal for s system having 2500 atoms?
 


Yes.  The free energy code slows down runs quite a bit.

-Justin

--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


RE: [gmx-users] no output dgdl file

2011-03-23 Thread Emanuel Birru
Hi,

 

Your foreign_lambda value is only one, please put all your lambda values
separated by space and you will get the dhdl file. And make sure that
you use -deffnm when you run your mdrun to get all the out put files by
default.

 

Cheers,

 

From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Moeed
Sent: Thursday, 24 March 2011 4:13 AM
To: gmx-users@gromacs.org
Subject: [gmx-users] no output dgdl file

 

Hello,

I have a little problem with FE output file. Below is the settings and
also I am including -dgdl in the command I issue but no dgdl (or dhdl)
file generates. I dont figure where the problem lies ! (version 4.5.3). 

free_energy  =   yes
init_lambda  =   0  
delta_lambda =   0
sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q
couple-lambda1   =   none 
couple-intramol  =   yes
nstdhdl  =   10
separate_dhdl_file   =   yes
dh_hist_size =   0
dh_hist_spacing  =   0.1

Best

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

Hi,

 

Your foreign_lambda value is only one, please put all your lambda values 
separated by space and you will get the dhdl file. And make sure that 


This is not required.  You can specify as many or as few foreign_lambda values 
as you like.


you use –deffnm when you run your mdrun to get all the out put files by 
default.




The -deffnm flag controls the names of the files, not which ones are written.

-Justin

 


Cheers,

 

*From:* gmx-users-boun...@gromacs.org 
[mailto:gmx-users-boun...@gromacs.org] *On Behalf Of *Moeed

*Sent:* Thursday, 24 March 2011 4:13 AM
*To:* gmx-users@gromacs.org
*Subject:* [gmx-users] no output dgdl file

 


Hello,

I have a little problem with FE output file. Below is the settings and 
also I am including -dgdl in the command I issue but no dgdl (or dhdl) 
file generates. I dont figure where the problem lies ! (version 4.5.3).


free_energy  =   yes
init_lambda  =   0 
delta_lambda =   0

sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q   
couple-lambda1   =   none

couple-intramol  =   yes
nstdhdl  =   10
separate_dhdl_file   =   yes
dh_hist_size =   0
dh_hist_spacing  =   0.1

Best



--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

Hi Justin,

Sure what you wrote is correct, what I am trying to tell to him is that
if he has more than one foreign lambda it is better to put all of them
as it is not logical to use only one foreign lambda to calculate FE (0 -
0.1) and he doesn't have delta lambda too. -Deffnm is sure all about


You can do a series of simulations at many values of lambda, each specifying 
their own foreign_lambda.  It seems to me that this method would be more 
reliable, but I have not tested simply attempting a single simulation with all 
values of foreign_lambda.  Using delta_lambda is (from all that I have read) not 
reliable, as there are known issues with slow growth methods.



file names but it also generate all the necessary files including the
xvg's without the need of using -dhdl.



Whether or not one specifies -dhdl or -deffnm is independent of whether or not 
dhdl.xvg (or whatever name) is written.  It is controlled purely by the presence 
of free_energy = yes in the input file.



I would appreciate if you come up with a new solution for him.



As reported by the OP, this issue had been effectively solved already:

http://lists.gromacs.org/pipermail/gmx-users/2011-March/059631.html

-Justin

Cheers 


-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 12:35 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi,

 


Your foreign_lambda value is only one, please put all your lambda
values 
separated by space and you will get the dhdl file. And make sure that 


This is not required.  You can specify as many or as few foreign_lambda
values 
as you like.



you use -deffnm when you run your mdrun to get all the out put files
by 

default.



The -deffnm flag controls the names of the files, not which ones are
written.

-Justin

 


Cheers,

 

*From:* gmx-users-boun...@gromacs.org 
[mailto:gmx-users-boun...@gromacs.org] *On Behalf Of *Moeed

*Sent:* Thursday, 24 March 2011 4:13 AM
*To:* gmx-users@gromacs.org
*Subject:* [gmx-users] no output dgdl file

 


Hello,

I have a little problem with FE output file. Below is the settings and



also I am including -dgdl in the command I issue but no dgdl (or dhdl)



file generates. I dont figure where the problem lies ! (version

4.5.3).

free_energy  =   yes
init_lambda  =   0 
delta_lambda =   0

sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3
foreign_lambda   =   0.1
dhdl_derivatives =   yes
couple-moltype   =   Polymer
couple-lambda0   =   vdw-q   
couple-lambda1   =   none

couple-intramol  =   yes
nstdhdl  =   10
separate_dhdl_file   =   yes
dh_hist_size =   0
dh_hist_spacing  =   0.1

Best





--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


[gmx-users] No cut-off

2011-03-23 Thread Swarnendu Tripathi
Hello everybody,

I want to use the no cut-off option in gromacs for the electrostatic
interactions. In manual it says for this I need to define: pbc=no;
nstlist=0; ns-type=simple and rlist=rcoulomb=rvdw=0 in the .mdp file.

My questions are:

1. If I choose rcoulomb=0 and rlist=rvdw=2.0 is that a problem? Do you
always recommend to use rlist=rcoulomb=rvdw=0 for no cut-off option?

2. How  should I choose the table-extension parameter now? Before, with
pbc=xyz (with cut-off) I used table-extension=rvdw+1/2*length of the longest
diagonal of the box (approximately). I am also using a tabulated potential.

Thank you,

-Swarnendu
-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

RE: [gmx-users] no output dgdl file

2011-03-23 Thread Emanuel Birru
Hi Justin,

It good that the issue is solved. As per my experience if you want to do
a series of simulations it is not necessary to use foreign_lambda, for
each simulation we can give different lambda values starting from 0 to
1( from fully interactive to non-interactive) no need of foreign_lamda
and certainly no need of generating dhdl.xvg to calculate FE. I think
-dhdl is an optional it will not be generated just because free_energy
= yes is present. The problem with the 4.5.3 is not the problem of
generating dhdl data by using single simulation with all foreign_lambda
values. The problem is with g_bar, when I tried to analyse the dhdl.xvg
output (with all the necessary data in it) using g_bar, it doesn't
function properly. The error is related with the source code. I guess
the main advantage of using 4.5.3 to calculated FE using foreign_lambda
was to avoid running series of independent simulations. May be it is
solved on the newer version 4.5.4., didn't try it yet.

Cheers,

  

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 12:54 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:
 Hi Justin,
 
 Sure what you wrote is correct, what I am trying to tell to him is
that
 if he has more than one foreign lambda it is better to put all of them
 as it is not logical to use only one foreign lambda to calculate FE (0
-
 0.1) and he doesn't have delta lambda too. -Deffnm is sure all about

You can do a series of simulations at many values of lambda, each
specifying 
their own foreign_lambda.  It seems to me that this method would be more

reliable, but I have not tested simply attempting a single simulation
with all 
values of foreign_lambda.  Using delta_lambda is (from all that I have
read) not 
reliable, as there are known issues with slow growth methods.

 file names but it also generate all the necessary files including the
 xvg's without the need of using -dhdl.
 

Whether or not one specifies -dhdl or -deffnm is independent of whether
or not 
dhdl.xvg (or whatever name) is written.  It is controlled purely by the
presence 
of free_energy = yes in the input file.

 I would appreciate if you come up with a new solution for him.
 

As reported by the OP, this issue had been effectively solved already:

http://lists.gromacs.org/pipermail/gmx-users/2011-March/059631.html

-Justin

 Cheers 
 
 -Original Message-
 From: gmx-users-boun...@gromacs.org
 [mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
 Sent: Thursday, 24 March 2011 12:35 PM
 To: Discussion list for GROMACS users
 Subject: Re: [gmx-users] no output dgdl file
 
 
 
 Emanuel Birru wrote:
 Hi,

  

 Your foreign_lambda value is only one, please put all your lambda
 values 
 separated by space and you will get the dhdl file. And make sure that

 
 This is not required.  You can specify as many or as few
foreign_lambda
 values 
 as you like.
 
 you use -deffnm when you run your mdrun to get all the out put files
 by 
 default.

 
 The -deffnm flag controls the names of the files, not which ones are
 written.
 
 -Justin
 
  

 Cheers,

  

 *From:* gmx-users-boun...@gromacs.org 
 [mailto:gmx-users-boun...@gromacs.org] *On Behalf Of *Moeed
 *Sent:* Thursday, 24 March 2011 4:13 AM
 *To:* gmx-users@gromacs.org
 *Subject:* [gmx-users] no output dgdl file

  

 Hello,

 I have a little problem with FE output file. Below is the settings
and
 
 also I am including -dgdl in the command I issue but no dgdl (or
dhdl)
 
 file generates. I dont figure where the problem lies ! (version
 4.5.3).
 free_energy  =   yes
 init_lambda  =   0 
 delta_lambda =   0
 sc_alpha =   0.5
 sc-power =   1
 sc_sigma =   0.3
 foreign_lambda   =   0.1
 dhdl_derivatives =   yes
 couple-moltype   =   Polymer
 couple-lambda0   =   vdw-q   
 couple-lambda1   =   none
 couple-intramol  =   yes
 nstdhdl  =   10
 separate_dhdl_file   =   yes
 dh_hist_size =   0
 dh_hist_spacing  =   0.1

 Best

 

-- 


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the 

Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

Hi Justin,

It good that the issue is solved. As per my experience if you want to do
a series of simulations it is not necessary to use foreign_lambda, for
each simulation we can give different lambda values starting from 0 to
1( from fully interactive to non-interactive) no need of foreign_lamda


This is certainly a viable method, but it is thermodynamic integration, not BAR. 
 The purpose of lambda/foreign_lambda is to use the BAR method.



and certainly no need of generating dhdl.xvg to calculate FE. I think
-dhdl is an optional it will not be generated just because free_energy
= yes is present. The problem with the 4.5.3 is not the problem of


It has always been my experience (in versions 3.3.3 and 4.5.3) that if the free 
energy code is activated, this file is written.  I never specify the file names 
independently.



generating dhdl data by using single simulation with all foreign_lambda
values. The problem is with g_bar, when I tried to analyse the dhdl.xvg
output (with all the necessary data in it) using g_bar, it doesn't
function properly. The error is related with the source code. I guess


If you suspect a bug, you should report it.  If you don't use foreign_lambda, 
you can't use the BAR method.  You're doing TI, not BAR.  They are fundamentally 
different, and the analysis for TI is done independently of any Gromacs tool.



the main advantage of using 4.5.3 to calculated FE using foreign_lambda
was to avoid running series of independent simulations. May be it is
solved on the newer version 4.5.4., didn't try it yet.



Version 4.5.3 has worked just fine for all free energy calculations I have done. 
  Perhaps BAR can be done with a single simulation and multiple foreign_lambda, 
but that was not my understanding of the proper procedure (based on some posts 
by developers a few months ago).  Multiple simulations, each at native lambda, 
are conducted with values of foreign_lambda above and below the native value at 
some lambda spacing (except for end points).  Invoking g_bar gives the total 
DeltaG over all lambda intervals according to the BAR algorithm.


-Justin


Cheers,

  


-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 12:54 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi Justin,

Sure what you wrote is correct, what I am trying to tell to him is

that

if he has more than one foreign lambda it is better to put all of them
as it is not logical to use only one foreign lambda to calculate FE (0

-

0.1) and he doesn't have delta lambda too. -Deffnm is sure all about


You can do a series of simulations at many values of lambda, each
specifying 
their own foreign_lambda.  It seems to me that this method would be more


reliable, but I have not tested simply attempting a single simulation
with all 
values of foreign_lambda.  Using delta_lambda is (from all that I have
read) not 
reliable, as there are known issues with slow growth methods.



file names but it also generate all the necessary files including the
xvg's without the need of using -dhdl.



Whether or not one specifies -dhdl or -deffnm is independent of whether
or not 
dhdl.xvg (or whatever name) is written.  It is controlled purely by the
presence 
of free_energy = yes in the input file.



I would appreciate if you come up with a new solution for him.



As reported by the OP, this issue had been effectively solved already:

http://lists.gromacs.org/pipermail/gmx-users/2011-March/059631.html

-Justin

Cheers 


-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 12:35 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi,

 


Your foreign_lambda value is only one, please put all your lambda
values 

separated by space and you will get the dhdl file. And make sure that



This is not required.  You can specify as many or as few

foreign_lambda
values 
as you like.



you use -deffnm when you run your mdrun to get all the out put files
by 

default.


The -deffnm flag controls the names of the files, not which ones are
written.

-Justin

 


Cheers,

 

*From:* gmx-users-boun...@gromacs.org 
[mailto:gmx-users-boun...@gromacs.org] *On Behalf Of *Moeed

*Sent:* Thursday, 24 March 2011 4:13 AM
*To:* gmx-users@gromacs.org
*Subject:* [gmx-users] no output dgdl file

 


Hello,

I have a little problem with FE output file. Below is the settings

and

also I am including -dgdl in the command I issue but no dgdl (or

dhdl)

file generates. I dont figure where the problem lies ! (version

4.5.3).

free_energy  =   yes
init_lambda  =   0 
delta_lambda =   0

sc_alpha =   0.5
sc-power =   1
sc_sigma =   0.3

RE: [gmx-users] no output dgdl file

2011-03-23 Thread Emanuel Birru
Yeah, I am using IT and do analyse the result using another method not
bar. But I used g_bar when I was using the foreign_lambda and simulate
all in a single file. I have already sent my suspect few weeks back. 

I am a bit confused on the g_bar part, when it says -f expects multiple
dhdl files. Do we need to run still multiple independent simulations
using different foreign_lambda values? I do not see why we should run
independent simulations, if we use for couple-lambda0 and couple-lambda1
vdw-q and none respectively.

For IT I am using 4.5.3 and it is working good so far.

Cheers,

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 1:36 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:
 Hi Justin,
 
 It good that the issue is solved. As per my experience if you want to
do
 a series of simulations it is not necessary to use foreign_lambda, for
 each simulation we can give different lambda values starting from 0 to
 1( from fully interactive to non-interactive) no need of foreign_lamda

This is certainly a viable method, but it is thermodynamic integration,
not BAR. 
  The purpose of lambda/foreign_lambda is to use the BAR method.

 and certainly no need of generating dhdl.xvg to calculate FE. I think
 -dhdl is an optional it will not be generated just because
free_energy
 = yes is present. The problem with the 4.5.3 is not the problem of

It has always been my experience (in versions 3.3.3 and 4.5.3) that if
the free 
energy code is activated, this file is written.  I never specify the
file names 
independently.

 generating dhdl data by using single simulation with all
foreign_lambda
 values. The problem is with g_bar, when I tried to analyse the
dhdl.xvg
 output (with all the necessary data in it) using g_bar, it doesn't
 function properly. The error is related with the source code. I guess

If you suspect a bug, you should report it.  If you don't use
foreign_lambda, 
you can't use the BAR method.  You're doing TI, not BAR.  They are
fundamentally 
different, and the analysis for TI is done independently of any Gromacs
tool.

 the main advantage of using 4.5.3 to calculated FE using
foreign_lambda
 was to avoid running series of independent simulations. May be it is
 solved on the newer version 4.5.4., didn't try it yet.
 

Version 4.5.3 has worked just fine for all free energy calculations I
have done. 
   Perhaps BAR can be done with a single simulation and multiple
foreign_lambda, 
but that was not my understanding of the proper procedure (based on some
posts 
by developers a few months ago).  Multiple simulations, each at native
lambda, 
are conducted with values of foreign_lambda above and below the native
value at 
some lambda spacing (except for end points).  Invoking g_bar gives the
total 
DeltaG over all lambda intervals according to the BAR algorithm.

-Justin

 Cheers,
 
   
 
 -Original Message-
 From: gmx-users-boun...@gromacs.org
 [mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
 Sent: Thursday, 24 March 2011 12:54 PM
 To: Gromacs Users' List
 Subject: Re: [gmx-users] no output dgdl file
 
 
 
 Emanuel Birru wrote:
 Hi Justin,

 Sure what you wrote is correct, what I am trying to tell to him is
 that
 if he has more than one foreign lambda it is better to put all of
them
 as it is not logical to use only one foreign lambda to calculate FE
(0
 -
 0.1) and he doesn't have delta lambda too. -Deffnm is sure all about
 
 You can do a series of simulations at many values of lambda, each
 specifying 
 their own foreign_lambda.  It seems to me that this method would be
more
 
 reliable, but I have not tested simply attempting a single simulation
 with all 
 values of foreign_lambda.  Using delta_lambda is (from all that I have
 read) not 
 reliable, as there are known issues with slow growth methods.
 
 file names but it also generate all the necessary files including the
 xvg's without the need of using -dhdl.

 
 Whether or not one specifies -dhdl or -deffnm is independent of
whether
 or not 
 dhdl.xvg (or whatever name) is written.  It is controlled purely by
the
 presence 
 of free_energy = yes in the input file.
 
 I would appreciate if you come up with a new solution for him.

 
 As reported by the OP, this issue had been effectively solved already:
 
 http://lists.gromacs.org/pipermail/gmx-users/2011-March/059631.html
 
 -Justin
 
 Cheers 

 -Original Message-
 From: gmx-users-boun...@gromacs.org
 [mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
 Sent: Thursday, 24 March 2011 12:35 PM
 To: Discussion list for GROMACS users
 Subject: Re: [gmx-users] no output dgdl file



 Emanuel Birru wrote:
 Hi,

  

 Your foreign_lambda value is only one, please put all your lambda
 values 
 separated by space and you will get the dhdl file. And make sure
that
 
 This is not required.  You 

Re: [gmx-users] No cut-off

2011-03-23 Thread Itamar Kass

Hi Swarnendu,

grompp will return errer unless rlist=rcoulomb=rvdw=0, so you should 
stick to it.


Cheers,
Itamar

On 24/03/11 1:13 PM, Swarnendu Tripathi wrote:

Hello everybody,

I want to use the no cut-off option in gromacs for the electrostatic 
interactions. In manual it says for this I need to define: pbc=no; 
nstlist=0; ns-type=simple and rlist=rcoulomb=rvdw=0 in the .mdp file.


My questions are:

1. If I choose rcoulomb=0 and rlist=rvdw=2.0 is that a problem? Do you 
always recommend to use rlist=rcoulomb=rvdw=0 for no cut-off option?


2. How  should I choose the table-extension parameter now? Before, 
with pbc=xyz (with cut-off) I used table-extension=rvdw+1/2*length of 
the longest diagonal of the box (approximately). I am also using a 
tabulated potential.


Thank you,

-Swarnendu


--


In theory, there is no difference between theory and practice. But, in practice, 
there is. - Jan L.A. van de Snepscheut

===
| Itamar Kass, Ph.D.
| Postdoctoral Research Fellow
|
| Department of Biochemistry and Molecular Biology
| Building 77 Clayton Campus
| Wellington Road
| Monash University,
| Victoria 3800
| Australia
|
| Tel: +61 3 9902 9376
| Fax: +61 3 9902 9500
| E-mail: itamar.k...@monash.edu


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

Yeah, I am using IT and do analyse the result using another method not
bar. But I used g_bar when I was using the foreign_lambda and simulate
all in a single file. I have already sent my suspect few weeks back. 



Is there an active redmine issue?  I cannot find any post from you in the list 
archive from the last few months.  What is the issue?



I am a bit confused on the g_bar part, when it says -f expects multiple
dhdl files. Do we need to run still multiple independent simulations
using different foreign_lambda values? I do not see why we should run
independent simulations, if we use for couple-lambda0 and couple-lambda1
vdw-q and none respectively.



As I mentioned to the OP of this thread, simultaneous (de)coupling of vdW and 
Coulombic interactions is not stable.  The output energies are not trustworthy, 
in my experience due to physically unreasonable configurations and the potential 
for numerical singularities.


The multiple files that g_bar expects are not simply from two separate 
processes, however.  It expects dhdl.xvg files from multiple values of native 
lambda that have corresponding foreign_lambda values in it, i.e.:


Interval 1:
init_lambda = 0
foreign_lambda = 0.05

Interval 2:
init_lambda = 0.05
foreign_lambda = 0 0.1

Interval 3:
init_lambda = 0.1
foreign_lambda = 0.05 0.15

etc.

-Justin


For IT I am using 4.5.3 and it is working good so far.

Cheers,

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 1:36 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi Justin,

It good that the issue is solved. As per my experience if you want to

do

a series of simulations it is not necessary to use foreign_lambda, for
each simulation we can give different lambda values starting from 0 to
1( from fully interactive to non-interactive) no need of foreign_lamda


This is certainly a viable method, but it is thermodynamic integration,
not BAR. 
  The purpose of lambda/foreign_lambda is to use the BAR method.



and certainly no need of generating dhdl.xvg to calculate FE. I think
-dhdl is an optional it will not be generated just because

free_energy

= yes is present. The problem with the 4.5.3 is not the problem of


It has always been my experience (in versions 3.3.3 and 4.5.3) that if
the free 
energy code is activated, this file is written.  I never specify the
file names 
independently.



generating dhdl data by using single simulation with all

foreign_lambda

values. The problem is with g_bar, when I tried to analyse the

dhdl.xvg

output (with all the necessary data in it) using g_bar, it doesn't
function properly. The error is related with the source code. I guess


If you suspect a bug, you should report it.  If you don't use
foreign_lambda, 
you can't use the BAR method.  You're doing TI, not BAR.  They are
fundamentally 
different, and the analysis for TI is done independently of any Gromacs

tool.


the main advantage of using 4.5.3 to calculated FE using

foreign_lambda

was to avoid running series of independent simulations. May be it is
solved on the newer version 4.5.4., didn't try it yet.



Version 4.5.3 has worked just fine for all free energy calculations I
have done. 
   Perhaps BAR can be done with a single simulation and multiple
foreign_lambda, 
but that was not my understanding of the proper procedure (based on some
posts 
by developers a few months ago).  Multiple simulations, each at native
lambda, 
are conducted with values of foreign_lambda above and below the native
value at 
some lambda spacing (except for end points).  Invoking g_bar gives the
total 
DeltaG over all lambda intervals according to the BAR algorithm.


-Justin


Cheers,

  


-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 12:54 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi Justin,

Sure what you wrote is correct, what I am trying to tell to him is

that

if he has more than one foreign lambda it is better to put all of

them

as it is not logical to use only one foreign lambda to calculate FE

(0

-

0.1) and he doesn't have delta lambda too. -Deffnm is sure all about

You can do a series of simulations at many values of lambda, each
specifying 
their own foreign_lambda.  It seems to me that this method would be

more

reliable, but I have not tested simply attempting a single simulation
with all 
values of foreign_lambda.  Using delta_lambda is (from all that I have
read) not 
reliable, as there are known issues with slow growth methods.



file names but it also generate all the necessary files including the
xvg's without the need of using -dhdl.


Whether or not one specifies -dhdl or -deffnm is independent of

whether
or not 
dhdl.xvg (or 

RE: [gmx-users] no output dgdl file

2011-03-23 Thread Emanuel Birru
I have tried the free energy calculation using 4.5.3, but when I run
g_bar it still gives me the following and do not generate the output
files.

prod.xvg: 0.0 - 5000.0; lambda = 0.000
foreign lambdas: 0.000 (125001 pts) 0.000 (125001 pts) 0.010 (125001
pts) 0.030 (125001 pts) 0.050 (125001 pts) 0.100 (125001 pts) 0.200
(125001 pts) 0.300 (125001 pts) 0.400 (125001 pts) 0.500 (125001 pts)
0.600 (125001 pts) 0.650 (125001 pts) 0.700 (125001 pts) 0.750 (125001
pts) 0.800 (125001 pts) 0.850 (125001 pts) 0.900 (125001 pts) 0.950
(125001 pts) 1.000 (125001 pts)

No results to calculate.

I put my mdp file for production  run and the output file prod.xvg on
the following links for your consideration.

http://hydra.pharm.monash.edu.au/md_project/production.txt (mdp file)
http://hydra.pharm.monash.edu.au/md_project/prod.txt (xvg file)

Cheers,
Emanuel

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 2:02 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:
 Yeah, I am using IT and do analyse the result using another method not
 bar. But I used g_bar when I was using the foreign_lambda and simulate
 all in a single file. I have already sent my suspect few weeks back. 
 

Is there an active redmine issue?  I cannot find any post from you in
the list 
archive from the last few months.  What is the issue?

 I am a bit confused on the g_bar part, when it says -f expects
multiple
 dhdl files. Do we need to run still multiple independent simulations
 using different foreign_lambda values? I do not see why we should run
 independent simulations, if we use for couple-lambda0 and
couple-lambda1
 vdw-q and none respectively.
 

As I mentioned to the OP of this thread, simultaneous (de)coupling of
vdW and 
Coulombic interactions is not stable.  The output energies are not
trustworthy, 
in my experience due to physically unreasonable configurations and the
potential 
for numerical singularities.

The multiple files that g_bar expects are not simply from two separate 
processes, however.  It expects dhdl.xvg files from multiple values of
native 
lambda that have corresponding foreign_lambda values in it, i.e.:

Interval 1:
init_lambda = 0
foreign_lambda = 0.05

Interval 2:
init_lambda = 0.05
foreign_lambda = 0 0.1

Interval 3:
init_lambda = 0.1
foreign_lambda = 0.05 0.15

etc.

-Justin

 For IT I am using 4.5.3 and it is working good so far.
 
 Cheers,
 
 -Original Message-
 From: gmx-users-boun...@gromacs.org
 [mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
 Sent: Thursday, 24 March 2011 1:36 PM
 To: Gromacs Users' List
 Subject: Re: [gmx-users] no output dgdl file
 
 
 
 Emanuel Birru wrote:
 Hi Justin,

 It good that the issue is solved. As per my experience if you want to
 do
 a series of simulations it is not necessary to use foreign_lambda,
for
 each simulation we can give different lambda values starting from 0
to
 1( from fully interactive to non-interactive) no need of
foreign_lamda
 
 This is certainly a viable method, but it is thermodynamic
integration,
 not BAR. 
   The purpose of lambda/foreign_lambda is to use the BAR method.
 
 and certainly no need of generating dhdl.xvg to calculate FE. I think
 -dhdl is an optional it will not be generated just because
 free_energy
 = yes is present. The problem with the 4.5.3 is not the problem of
 
 It has always been my experience (in versions 3.3.3 and 4.5.3) that if
 the free 
 energy code is activated, this file is written.  I never specify the
 file names 
 independently.
 
 generating dhdl data by using single simulation with all
 foreign_lambda
 values. The problem is with g_bar, when I tried to analyse the
 dhdl.xvg
 output (with all the necessary data in it) using g_bar, it doesn't
 function properly. The error is related with the source code. I guess
 
 If you suspect a bug, you should report it.  If you don't use
 foreign_lambda, 
 you can't use the BAR method.  You're doing TI, not BAR.  They are
 fundamentally 
 different, and the analysis for TI is done independently of any
Gromacs
 tool.
 
 the main advantage of using 4.5.3 to calculated FE using
 foreign_lambda
 was to avoid running series of independent simulations. May be it is
 solved on the newer version 4.5.4., didn't try it yet.

 
 Version 4.5.3 has worked just fine for all free energy calculations I
 have done. 
Perhaps BAR can be done with a single simulation and multiple
 foreign_lambda, 
 but that was not my understanding of the proper procedure (based on
some
 posts 
 by developers a few months ago).  Multiple simulations, each at
native
 lambda, 
 are conducted with values of foreign_lambda above and below the native
 value at 
 some lambda spacing (except for end points).  Invoking g_bar gives the
 total 
 DeltaG over all lambda intervals according to the BAR algorithm.
 
 -Justin
 
 Cheers,

   

 

Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

I have tried the free energy calculation using 4.5.3, but when I run
g_bar it still gives me the following and do not generate the output
files.

prod.xvg: 0.0 - 5000.0; lambda = 0.000
foreign lambdas: 0.000 (125001 pts) 0.000 (125001 pts) 0.010 (125001
pts) 0.030 (125001 pts) 0.050 (125001 pts) 0.100 (125001 pts) 0.200
(125001 pts) 0.300 (125001 pts) 0.400 (125001 pts) 0.500 (125001 pts)
0.600 (125001 pts) 0.650 (125001 pts) 0.700 (125001 pts) 0.750 (125001
pts) 0.800 (125001 pts) 0.850 (125001 pts) 0.900 (125001 pts) 0.950
(125001 pts) 1.000 (125001 pts)

No results to calculate.

I put my mdp file for production  run and the output file prod.xvg on
the following links for your consideration.

http://hydra.pharm.monash.edu.au/md_project/production.txt (mdp file)
http://hydra.pharm.monash.edu.au/md_project/prod.txt (xvg file)



From the above g_bar output and the .xvg file, it seems that lambda=0 is 
considered both native and foreign, which probably causes the algorithm to fail. 
 Note, too, that the results in the first column of the .xvg file (which would 
correspond to some weird difference between lambda=0 and lambda=0, which should 
not exist!) also seem to be nonsensically high.  This is what I mean when I say 
weird results from simultaneous (de)coupling of charges and LJ terms.  The 
systems are not completely stable, yielding weird energies.  The result could 
also be some bizarre result of lambda=0 being considered native and foreign.


Based on the .mdp file, I don't know how this could have happened, but all signs 
point to some sort of problem with the input.


-Justin


Cheers,
Emanuel

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 2:02 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Yeah, I am using IT and do analyse the result using another method not
bar. But I used g_bar when I was using the foreign_lambda and simulate
all in a single file. I have already sent my suspect few weeks back. 



Is there an active redmine issue?  I cannot find any post from you in
the list 
archive from the last few months.  What is the issue?



I am a bit confused on the g_bar part, when it says -f expects

multiple

dhdl files. Do we need to run still multiple independent simulations
using different foreign_lambda values? I do not see why we should run
independent simulations, if we use for couple-lambda0 and

couple-lambda1

vdw-q and none respectively.



As I mentioned to the OP of this thread, simultaneous (de)coupling of
vdW and 
Coulombic interactions is not stable.  The output energies are not
trustworthy, 
in my experience due to physically unreasonable configurations and the
potential 
for numerical singularities.


The multiple files that g_bar expects are not simply from two separate 
processes, however.  It expects dhdl.xvg files from multiple values of
native 
lambda that have corresponding foreign_lambda values in it, i.e.:


Interval 1:
init_lambda = 0
foreign_lambda = 0.05

Interval 2:
init_lambda = 0.05
foreign_lambda = 0 0.1

Interval 3:
init_lambda = 0.1
foreign_lambda = 0.05 0.15

etc.

-Justin


For IT I am using 4.5.3 and it is working good so far.

Cheers,

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 1:36 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file



Emanuel Birru wrote:

Hi Justin,

It good that the issue is solved. As per my experience if you want to

do

a series of simulations it is not necessary to use foreign_lambda,

for

each simulation we can give different lambda values starting from 0

to

1( from fully interactive to non-interactive) no need of

foreign_lamda

This is certainly a viable method, but it is thermodynamic

integration,
not BAR. 
  The purpose of lambda/foreign_lambda is to use the BAR method.



and certainly no need of generating dhdl.xvg to calculate FE. I think
-dhdl is an optional it will not be generated just because

free_energy

= yes is present. The problem with the 4.5.3 is not the problem of

It has always been my experience (in versions 3.3.3 and 4.5.3) that if
the free 
energy code is activated, this file is written.  I never specify the
file names 
independently.



generating dhdl data by using single simulation with all

foreign_lambda

values. The problem is with g_bar, when I tried to analyse the

dhdl.xvg

output (with all the necessary data in it) using g_bar, it doesn't
function properly. The error is related with the source code. I guess

If you suspect a bug, you should report it.  If you don't use
foreign_lambda, 
you can't use the BAR method.  You're doing TI, not BAR.  They are
fundamentally 
different, and the analysis for TI is done independently of any

Gromacs

tool.


the main advantage of 

RE: [gmx-users] no output dgdl file

2011-03-23 Thread Emanuel Birru
Thanks Justin,

 

I will revise my mdp file again and will give it a try once more. Can
you give me an idea on the following too.

 

I want to perform free energy calculation by changing the gd (dihedral
parameters) to fourier dihedral parameters for small molecules
(solutes). But I am getting the following error when I do grompp,  

 

Program grompp, VERSION 4.5.3

Source code file: topsort.c, line: 112

 

Fatal error:

Function type Fourier Dih. not implemented in ip_pert

For more information and tips for troubleshooting, please check the
GROMACS

website at http://www.gromacs.org/Documentation/Errors

 

Is there any simple solution for this or is it a source code error?

 

Here http://hydra.pharm.monash.edu.au/md_project/penteg.top  it is my
topology file.  

 

Cheers,

 

-Original Message-
From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul
Sent: Thursday, 24 March 2011 2:35 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file

 

 

 

Emanuel Birru wrote:

 

 

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] no output dgdl file

2011-03-23 Thread Justin A. Lemkul



Emanuel Birru wrote:

Thanks Justin,

 

I will revise my mdp file again and will give it a try once more. Can 
you give me an idea on the following too.


 

I want to perform free energy calculation by changing the gd (dihedral 
parameters) to fourier dihedral parameters for small molecules 
(solutes). But I am getting the following error when I do grompp,  

 


Program grompp, VERSION 4.5.3

Source code file: topsort.c, line: 112

 


Fatal error:

Function type Fourier Dih. not implemented in ip_pert

For more information and tips for troubleshooting, please check the GROMACS

website at http://www.gromacs.org/Documentation/Errors

 


Is there any simple solution for this or is it a source code error?

 


The manual states that Fourier dihedrals can have any of their coefficients 
perturbed, but this does not appear to actually be implemented.  Please file a 
redmine issue so that it can be fixed.


-Justin



Here http://hydra.pharm.monash.edu.au/md_project/penteg.top it is my 
topology file. 

 


Cheers,

 


-Original Message-
From: gmx-users-boun...@gromacs.org 
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Justin A. Lemkul

Sent: Thursday, 24 March 2011 2:35 PM
To: Gromacs Users' List
Subject: Re: [gmx-users] no output dgdl file

 

 

 


Emanuel Birru wrote:

 

 



--


Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin


--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] No cut-off

2011-03-23 Thread Mark Abraham

On 24/03/2011 1:34 PM, Itamar Kass wrote:

Hi Swarnendu,

grompp will return errer unless rlist=rcoulomb=rvdw=0, so you should 
stick to it.


Cheers,
Itamar

On 24/03/11 1:13 PM, Swarnendu Tripathi wrote:

Hello everybody,

I want to use the no cut-off option in gromacs for the electrostatic 
interactions. In manual it says for this I need to define: pbc=no; 
nstlist=0; ns-type=simple and rlist=rcoulomb=rvdw=0 in the .mdp file.


My questions are:

1. If I choose rcoulomb=0 and rlist=rvdw=2.0 is that a problem? Do 
you always recommend to use rlist=rcoulomb=rvdw=0 for no cut-off option?


2. How  should I choose the table-extension parameter now? Before, 
with pbc=xyz (with cut-off) I used table-extension=rvdw+1/2*length of 
the longest diagonal of the box (approximately). I am also using a 
tabulated potential.


Since there is effectively no box once pbc=no, your tables need only be 
as long as the longest possible interaction, plus some margin for when 
the structure rearranges.


Mark
--
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.

Can't post? Read http://www.gromacs.org/Support/Mailing_Lists


Re: [gmx-users] No cut-off

2011-03-23 Thread Swarnendu Tripathi
Hi all,

Thank you for the reply. The answer I got for the second question about
table-extension, I understand and agree with that.

Regarding my first question I asked, I did not get any error with grompp
even when I tried rlist=rvdw=2.0 and rcoulomb=0. For this I used the no
cut-off conditions which I have mentioned below in my previous e-mail. I am
using Gromacs-4.0.7.  Any further suggestions or comments will be very
helpful.

Thank you,

-Swarnendu

On Thu, Mar 24, 2011 at 12:17 AM, Mark Abraham mark.abra...@anu.edu.auwrote:

 On 24/03/2011 1:34 PM, Itamar Kass wrote:

 Hi Swarnendu,

 grompp will return errer unless rlist=rcoulomb=rvdw=0, so you should stick
 to it.

 Cheers,
 Itamar

 On 24/03/11 1:13 PM, Swarnendu Tripathi wrote:

 Hello everybody,

 I want to use the no cut-off option in gromacs for the electrostatic
 interactions. In manual it says for this I need to define: pbc=no;
 nstlist=0; ns-type=simple and rlist=rcoulomb=rvdw=0 in the .mdp file.

 My questions are:

 1. If I choose rcoulomb=0 and rlist=rvdw=2.0 is that a problem? Do you
 always recommend to use rlist=rcoulomb=rvdw=0 for no cut-off option?

 2. How  should I choose the table-extension parameter now? Before, with
 pbc=xyz (with cut-off) I used table-extension=rvdw+1/2*length of the longest
 diagonal of the box (approximately). I am also using a tabulated potential.


 Since there is effectively no box once pbc=no, your tables need only be as
 long as the longest possible interaction, plus some margin for when the
 structure rearranges.

 Mark

 --
 gmx-users mailing listgmx-users@gromacs.org
 http://lists.gromacs.org/mailman/listinfo/gmx-users
 Please search the archive at
 http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
 Please don't post (un)subscribe requests to the list. Use the www interface
 or send it to gmx-users-requ...@gromacs.org.
 Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists

Re: [gmx-users] No cut-off

2011-03-23 Thread Mark Abraham

On 24/03/2011 4:38 PM, Swarnendu Tripathi wrote:

Hi all,

Thank you for the reply. The answer I got for the second question 
about table-extension, I understand and agree with that.


Regarding my first question I asked, I did not get any error with 
grompp even when I tried rlist=rvdw=2.0 and rcoulomb=0. For this I 
used the no cut-off conditions which I have mentioned below in my 
previous e-mail. I am using Gromacs-4.0.7.  Any further suggestions or 
comments will be very helpful.


VDW interactions normally become insignificant much faster than Coulomb 
interactions, so treating all Coulomb interactions and only some VDW 
interactions is logical. However, performing the search for neighbours 
takes time, and given that you are computing the distances already for 
their Coulomb interaction, it is probably cheaper just to compute all 
the VDW interactions, i.e. rlist=rvdw=0. This might be different if you 
had a lot of atoms that had zero partial charge.


Mark

On Thu, Mar 24, 2011 at 12:17 AM, Mark Abraham 
mark.abra...@anu.edu.au mailto:mark.abra...@anu.edu.au wrote:


On 24/03/2011 1:34 PM, Itamar Kass wrote:

Hi Swarnendu,

grompp will return errer unless rlist=rcoulomb=rvdw=0, so you
should stick to it.

Cheers,
Itamar

On 24/03/11 1:13 PM, Swarnendu Tripathi wrote:

Hello everybody,

I want to use the no cut-off option in gromacs for the
electrostatic interactions. In manual it says for this I
need to define: pbc=no; nstlist=0; ns-type=simple and
rlist=rcoulomb=rvdw=0 in the .mdp file.

My questions are:

1. If I choose rcoulomb=0 and rlist=rvdw=2.0 is that a
problem? Do you always recommend to use
rlist=rcoulomb=rvdw=0 for no cut-off option?

2. How  should I choose the table-extension parameter now?
Before, with pbc=xyz (with cut-off) I used
table-extension=rvdw+1/2*length of the longest diagonal of
the box (approximately). I am also using a tabulated
potential.


Since there is effectively no box once pbc=no, your tables need
only be as long as the longest possible interaction, plus some
margin for when the structure rearranges.

Mark

-- 
gmx-users mailing list gmx-users@gromacs.org

mailto:gmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the www
interface or send it to gmx-users-requ...@gromacs.org
mailto:gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists




-- 
gmx-users mailing listgmx-users@gromacs.org
http://lists.gromacs.org/mailman/listinfo/gmx-users
Please search the archive at 
http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
Please don't post (un)subscribe requests to the list. Use the 
www interface or send it to gmx-users-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists