[gmx-users] incomparable results with files produced with gromacs3.2.1 and 3.3.1
Sampo Karkola
sampo.karkola at helsinki.fi
Fri Mar 23 09:38:43 CET 2007
Hi,
I would also like to add that when simulating the protein alone, I
observe rising energies in md simulation. Only minor differences in
energy minimisation energy values are observed, if the protein is first
energy minimised in vacuum or not. The energy minimisation in a water
box goes succesfully with decreasing potential energy but when the md
simulation starts, the energy rises immediately, not depending on the
initial vacuum minimisation. gmxcheck for the md tpr file produces this
output:
**********************
Checking coordinate file 1fdt_md.tpr
Reading file 1fdt_md.tpr, VERSION 3.3.1 (single precision)
Reading file 1fdt_md.tpr, VERSION 3.3.1 (single precision)
29863 atoms in file
coordinates found
box found
velocities found
Kinetic energy: 115455 (kJ/mol)
Assuming the number of degrees of freedom to be Natoms * 3 or Natoms * 2,
the velocities correspond to a temperature of the system
of 309.994 K or 464.991 K respectively.
Checking for atoms closer than 0.8 and not between 0.4 and 0.7,
relative to sum of Van der Waals distance:
WARNING: masses will be determined based on residue and atom names,
this can deviate from the real mass of the atom type
In case you use free energy of solvation predictions:
++++ PLEASE READ AND CITE THE FOLLOWING REFERENCE ++++
D. Eisenberg and A. D. McLachlan
Solvation energy in protein folding and binding
Nature 319 (1986) pp. 199-203
-------- -------- --- Thank You --- -------- --------
Opening library file /usr/share/gromacs/top/aminoacids.dat
Opening library file /usr/share/gromacs/top/atommass.dat
Opening library file /usr/share/gromacs/top/vdwradii.dat
Opening library file /usr/share/gromacs/top/dgsolv.dat
#Entries in atommass.dat: 82 vdwradii.dat: 26 dgsolv.dat: 7
*************************
and does not end. gmxcheck -c 1fdt_md.tpr gets stuck running but does
not produce more output. This happens also for a successfull simulation
of my aromatase model.
gmxcheck for md trr gives the following
*************************
Checking file 1fdt_md.trr
trn version: GMX_trn_file (single precision)
Reading frame 0 time 0.000
# Atoms 29863
Last frame 50 time 500.000
Item #frames Timestep (ps)
Step 51 10
Time 51 10
Lambda 51 10
Coords 51 10
Velocities 51 10
Forces 0
Box 51 10
gcq#159: "With a Lead Filled Snowshoe" (F. Zappa)
*************************
Looking forward to your suggestions,
Sampo
Sampo Karkola wrote:
> Hi Tsjerk and Mark,
>
> thanks for the answers.
>
> The topologies differ in one atom which is changed from NR in the old
> topology to N in new topology and also, there is one dihedral angle not
> defined in the old topology and then zero is used by grompp. In the new
> topology, the dihedral is defined as an improper one to keep the atoms
> in plane. I did a test run where I changed that particular dihedral to
> correspond to the old working topology but the results are the same,
> i.e. energies rising.
>
> About the energy minimisations in vacuum, I have not performed them when
> using gromacs 3.2.1, and I have upgraded to 3.3.1 so I can't compare. I
> really do not feel like installing 3.2.1 again :)
>
> It's not only the ligand. Also the protein itself (17beta-HSD, pdb code
> 1fdt) and the cofactor (NADP+, NDPP topology from ffG43a1.rtp) give
> these rising energies when simulated alone, so I do not think that the
> error is in the top files, as mentioned by Mark. But you are welcome to
> take a look at http://www.helsinki.fi/~karkola/topologies.tar
>
> Thanks again,
>
> Sampo
>
>
> Tsjerk Wassenaar wrote:
>> Hi Sampo,
>>
>> Thank you for the elaborate account of your problem :)
>>
>> So, the problems seem not to be due to grompp, since, when you use the
>> old .gro .mdp .top files and run grompp on them, you get the behaviour
>> expected from the simulation. Therefore, something happens before. It
>> would be very helpful if you could dig a bit deeper. Starting with the
>> original pdb files, try to track where the divergent behaviour comes
>> from. Are the topologies generated different? Does energy minimization
>> in vacuum behave differently between the versions and do these yield
>> different results? Etc., etc.
>> Besides, does it have to do with your ligand, i.e. do you get the same
>> results when you use the protein only?
>>
>> Best,
>>
>> Tsjerk
>>
>>
>> On 3/22/07, Sampo Karkola <sampo.karkola at helsinki.fi> wrote:
>>> Dear list,
>>>
>>> I am trying to simulate a crystal structure of a protein with a cofactor
>>> (NADP, I've used NDPP topology) and a (self-made) ligand in a box of
>>> water. The preparation of the simulation went normally: pdb2gmx,
>>> editconf, genbox, grompp (for ions), genion, grompp (for em) mdrun (em),
>>> grompp (for md), mdrun (md). Now, in the md simulation, the potential
>>> and total energies of the system rise immediately from the beginning,
>>> stabilize to a certain level and later the run crashed with LINCS
>>> warnings. Earlier I have succesfully simulated a model of an another
>>> enzyme with similar procedure and in that run the potential energy
>>> decreased immediately and stabilized with no crashes. I tried to
>>> simulate the protein alone, the cofactor alone and the ligand alone and
>>> the energies rise in all cases.
>>>
>>> Luckily, I had an older simulation (with gromacs 3.2.1, now I'm using
>>> 3.3.1) with the same ligand (except one atom type is changed from NR in
>>> the old topology to N in new topology) alone and there the energies
>>> decreased as they should. To study the differences in the ligand-alone
>>> simulations, I performed an mdrun with the old tpr file (mdrun -s
>>> old.tpr -o new_run -c new_run...) and the energies decreased again, as
>>> they should. Next, I grompped a new tpr-file from the old gro, top and
>>> mdp files, got a new tpr file and simulated it resulting in decreasing
>>> energies again. But when I try to do the whole process from the
>>> beginning using the same pdb file as in the older run something goes
>>> differently and the mdrun produces rising energies. The energy
>>> minimisation with the new run takes over 3000 steps while the old run
>>> was minimised in 34 steps. I've checked that the gro and top files
>>> between the old (working) and the new (non-working) runs are identical.
>>> The only differences I find are in the tpr files of the energy
>>> minimisation and md simulation and they look like this (compared with
>>> diff):
>>>
>>> Energy minimisation (hm2 is the new bad behaving simulation and hm5 is
>>> the old, good one)
>>>
>>> 1c1
>>> < hm2_em.tpr:
>>> ---
>>> > hm5_em.tpr:
>>> 78c78
>>> < sc_power = 0
>>> ---
>>> > sc_power = 2
>>> 89c89
>>> < dihre-fc = 1000
>>> ---
>>> > dihre-fc = 1.4013e-42
>>> 142c142
>>> < scalefactor = 1
>>> ---
>>> > scalefactor = 0
>>> 22983c22983
>>> < atom[ 4]={type= 2, typeB= 2, ptype= Atom, m=
>>> 1.40067e+01, q= 4.80000e-02, mB= 1.40067e+01, qB= 4.80000e-02, resnr=
>>> 0} grpnrs=[ 0 0 0 0 0 0 0
>>> 0 0 0 ]}
>>> ---
>>> > atom[ 4]={type= 1, typeB= 1, ptype= Atom, m=
>>> 1.40067e+01, q= 4.80000e-02, mB= 1.40067e+01, qB= 4.80000e-02, resnr=
>>> 0} grpnrs=[ 0 0 0 0 0 0 0
>>> 0 0 0 ]}
>>> 22985c22985
>>> < atom[ 6]={type= 3, typeB= 3, ptype= Atom, m=
>>> 1.59994e+01, q=-5.56000e-01, mB= 1.59994e+01, qB=-5.56000e-01, resnr=
>>> 0} grpnrs=[ 0 0 0 0 0 0 0
>>> 0 0 0 ]}
>>> ---
>>> > atom[ 6]={type= 2, typeB= 2, ptype= Atom, m=
>>> 1.59994e+01, q=-5.56000e-01, mB= 1.59994e+01, qB=-5.56000e-01, resnr=
>>> 0} grpnrs=[ 0 0 0 0 0 0 0
>>> 0 0 0 ]}
>>>
>>> plus a lot of these atom definitions.
>>>
>>>
>>> MD simulation
>>>
>>> 1c1
>>> < hm2_md.tpr:
>>> ---
>>> > hm5_md.tpr:
>>> 78c78
>>> < sc_power = 0
>>> ---
>>> > sc_power = 2
>>> 89c89
>>> < dihre-fc = 1000
>>> ---
>>> > dihre-fc = 1.4013e-42
>>> 143c143
>>> < scalefactor = 1
>>> ---
>>> > scalefactor = 0
>>> 160,11566c160,11566
>>> < x[ 0]={ 2.47200e+00, 2.49800e+00, 2.51500e+00}
>>> < x[ 1]={ 2.37100e+00, 2.47900e+00, 2.61100e+00}
>>> < x[ 2]={ 2.39700e+00, 2.50800e+00, 2.74600e+00}
>>> < x[ 3]={ 2.52300e+00, 2.55600e+00, 2.78600e+00}
>>> < x[ 4]={ 2.62400e+00, 2.58300e+00, 2.69000e+00}
>>> < x[ 5]={ 2.59800e+00, 2.54700e+00, 2.55500e+00}
>>> < x[ 6]={ 2.69600e+00, 2.51500e+00, 2.48800e+00}
>>> < x[ 7]={ 2.42800e+00, 2.46600e+00, 2.38500e+00}
>>> < x[ 8]={ 2.29600e+00, 2.42000e+00, 2.38200e+00}
>>> < x[ 9]={ 2.21300e+00, 2.41600e+00, 2.54400e+00}
>>> < x[ 10]={ 2.49800e+00, 2.49100e+00, 2.26500e+00}
>>> < x[ 11]={ 2.45600e+00, 2.39700e+00, 2.14900e+00}
>>> < x[ 12]={ 2.31700e+00, 2.37400e+00, 2.14600e+00}
>>> < x[ 13]={ 2.24100e+00, 2.36700e+00, 2.26400e+00}
>>> < x[ 14]={ 2.73500e+00, 2.67900e+00, 2.71500e+00}
>>> < x[ 15]={ 2.82900e+00, 2.63400e+00, 2.82800e+00}
>>>
>>> plus a lot more coordination, velocity, atom type and interaction
>>> definitions.
>>>
>>> I have no dihedral restraints applied and those sc_power, dihre-fc and
>>> scalepower definitions are different even in the old, working simulation
>>> and in the one which was grompped again and producing correct results
>>> with decreasing energies.
>>>
>>> So the problem is that the old gro, top and mdp files give correct
>>> results even with grompped again, but a new procedure starting from the
>>> pdb file give bad results.
>>>
>>> I was hoping that someone could direct me to the correct direction and
>>> find out what's wrong. All suggestions are highly appreciated and
>>> tested!
>>>
>>>
>>> Kind regards,
>>>
>>> Sampo Karkola
>>>
>>>
>>> __________________________________
>>>
>>> Sampo Karkola
>>> Doctoral Student
>>>
>>> Laboratory of Organic Chemistry
>>> Department of Chemistry
>>> Faculty of Science
>>> PO Box 55 (A.I. Virtasen aukio 1)
>>> FIN-00014 University of Helsinki
>>> Finland
>>>
>>> email: sampo.karkola at helsinki.fi
>>>
>>> tel +358 9 19150369
>>> fax +358 9 19150366
>>>
>>>
>>> _______________________________________________
>>> gmx-users mailing list gmx-users at gromacs.org
>>> http://www.gromacs.org/mailman/listinfo/gmx-users
>>> Please don't post (un)subscribe requests to the list. Use the
>>> www interface or send it to gmx-users-request at gromacs.org.
>>> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>>>
>>
>>
> _______________________________________________
> gmx-users mailing list gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please don't post (un)subscribe requests to the list. Use the www
> interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>
More information about the gromacs.org_gmx-users
mailing list