Re: [gmx-users] g_dipole ? => dipole moment ?

[David van der Spoel]

On 2010-10-19 02.01, Chih-Ying Lin wrote:


Hi

 From David =>
"dipole is 48.0 sum of q_i x_i, therefore if you
have large charge separation you will get a large dipole."


1. What does the coefficient 48.0 represent?

Unit conversion from e nm to D.


2. Is q_i = partial charge or charge on the atoms / united atoms?
3. Is x_i = the distance between the mass center of atom from the mass
center of the molecule?
4. What does the large charge separation mean?  Do you mean the charged
molecule?  Or, Do you mean the molecule with a long carbon chain?



Please do a calculation on paper and you will understand.



Thank you
Lin






On 2010-10-18 03.30, Chih-Ying Lin wrote:

 HI
 I confined one molecule in the center of box and issue the g_dipole

command.

 The average dipole moment is still around 32.
 It is the molecule with 33 atoms / united atoms of most carbon groups,
 isn't the dipole moment around 32 too high?
 How can I test next  and  know that the dipole moment around 32 is
 acceptable?

By calculating on paper: dipole is 48.0 sum of q_i x_i, therefore if you
have large charge separation you will get a large dipole.


 Thank you
 Lin
 On 2010-10-16 21.36, Chih-Ying Lin wrote:
 >
 > Hi
 > I issue the g_dipole command on Gromacs => And, the following
 > information is shown.
 > There are 10 molecules in the selection,
 > Does the Average =32.1611 refer to the average for a single over the
 > simulation time?
 > Or, the Average = 32.1611 summing for all the 10 molecules over the
 > simulation time?
 > If the average = 32.1611 for a single molecule with 33 atoms / united
 > atoms of most carbon groups, isn't the dipole moment too high?
 I think this is the average per molecule. You may need to run trjconv
 -pbc whole, because mdrun may break molecules in two parts, meaning that
 the molecule becomes as big as the box.

 >
 >
 >
 >
 > What does "will subtract their charge at their center of mass"  this
 mean?
 > Why "will subtract their charge at their center of mass"  ?
 >
 >
 >





--
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se
--
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] A very strang e problem about version of pdb2gmx

[英雄不再寂寞]

Dear gmxers,
   I have reinstalled the gmx-4.5.1 and this time the pdb2gmx is found to work 
using the same input and parameter files. So I consider that the "very strange" 
problem is due to some modified file which is unknown yet for me. Thanks again 
for your attention.
 Sincerely, yours
 Chaofu Wu, Dr.
  -- 原始邮件 --
  发件人: "Justin A. Lemkul";
 发送时间: 2010年10月17日(星期天) 晚上8:21
 收件人: "Discussion list for GROMACS users"; 
 
 主题: Re: [gmx-users] A very strange problem about version of pdb2gmx

英雄不再寂寞 wrote:
> Dear gmxers,
>   Yesterday, I posted a mail titled "pdb2gmx stop at "AtomType 1" " 
> where I mentioned that I used gmx-4.5.1. The resulting problem is again 
> posted below. Today, I turn to the gmx-4.0.7, and use the same pdb 
> and rtp files with the name modified to compatible to the version. 
> Strangely, it is found that the pdb2gmx of this version can generate the 
> top file without any errors. Could you give me clues to why the latest 
> version can not whereas the older version can? To get the clear help, my 
> pdb and rtp files are also posted below. Thanks a lot for any reply.

Your system works for me under version 4.5.1, after fixing the naming mismatch 
you have between your .pdb file ("MOL") and .rtp file ("DET").  I suppose your 
actual input files have the correct names, though, otherwise this combination 
would have generated a fatal error, not a hang.

-Justin

> Chaofu Wu, Dr.
>  
> My command line and the screen output:
> xiaowu...@linux-s38y:~/workshop 
> > pdb2gmx -f deta.pdb -o 
> deta.gro -p deta.top -i deta.itp -ter
>  :-)  G  R  O  M  A  C  S  (-:
>Gromacs Runs One Microsecond At Cannonball Speeds
> :-)  VERSION 4.5.1  (-:
> ..
> Opening force field file oplsaa.ff/aminoacids.r2b
> Reading deta.pdb...
> Read 20 atoms
> Analyzing pdb file
> Splitting PDB chains based on TER records or changing chain id.
> There are 1 chains and 0 blocks of water and 1 residues with 20 atoms
>   chain  #res #atoms
>   1 ' ' 1 20 
> All occupancies are one
> Opening force field file oplsaa.ff/atomtypes.atp
> Atomtype 1
> 
> My pdb file:
> ATOM  1  N1  MOL 1  -5.199   1.343   0.334  1.00  
> 0.00   N   
> ATOM  2  C2  MOL 1  -4.582   0.639   1.470  1.00  
> 0.00   C   
> ATOM  3  C3  MOL 1  -3.897   1.462   2.593  1.00  
> 0.00   C   
> ATOM  4  N4  MOL 1  -3.194   0.615   3.574  1.00  
> 0.00   N   
> ATOM  5  C5  MOL 1  -2.373   1.331   4.566  1.00  
> 0.00   C   
> ATOM  6  C6  MOL 1  -1.520   0.397   5.464  1.00  
> 0.00   C   
> ATOM  7  N7  MOL 1  -0.684   1.133   6.432  1.00  
> 0.00   N   
> ATOM  8  H12 MOL 1  -4.614   2.067  -0.035  1.00  
> 0.00   H   
> ATOM  9  H11 MOL 1  -6.128   1.694   0.483  1.00  
> 0.00   H   
> ATOM 10  H22 MOL 1  -3.826  -0.053   1.056  1.00  
> 0.00   H   
> ATOM 11  H21 MOL 1  -5.345  -0.016   1.926  1.00  
> 0.00   H   
> ATOM 12  H31 MOL 1  -4.628   2.144   3.065  1.00  
> 0.00   H   
> ATOM 13  H32 MOL 1  -3.141   2.116   2.124  1.00  
> 0.00   H   
> ATOM 14  H41 MOL 1  -3.854   0.018   4.043  1.00  
> 0.00   H   
> ATOM 15  H51 MOL 1  -2.980   2.003   5.204  1.00  
> 0.00   H   
> ATOM 16  H52 MOL 1  -1.679   1.992   4.016  1.00  
> 0.00   H   
> ATOM 17  H62 MOL 1  -0.857  -0.180   4.795  1.00  
> 0.00   H   
> ATOM 18  H61 MOL 1  -2.188  -0.346   5.939  1.00  
> 0.00   H   
> ATOM 19  H72 MOL 1  -0.023   0.541   6.916  1.00  
> 0.00   H   
> ATOM 20  H71 MOL 1  -1.253   1.554   7.143  1.00  
> 0.00   H   
> TER
> My rtp file:
> [ DET ]
>  [ atoms ]
> N1opls_900-0.900  1
>H11opls_909 0.360  1
>H12opls_909 0.360  1
> C2opls_906 0.060  1
>H21opls_911 0.060  1
>H22opls_911 0.060  1
> C3opls_907 0.080  2
>H31opls_911 0.060  2
>H32opls_911 0.060  2
> N4opls_901-0.780  2
>H41opls_910 0.380  2
> C5opls_907 0.080  2
>H51opls_911 0.060  2
>H52opls_911 0.060  2
> C6opls_906 0.060  3
>H61opls_911 0.060  3
>H62opls_911 0.060  3
> N7opls_900-0.900  3
>H71opls_909 0.360  3
>H72opls_909 0.360  3
>  [ bonds ]
> N1   H11
> N1   H12
> N1C2
> C2   H21
> C2   H22
> C2C3
> C3   H31
> C3   H32
> C3N4
> N4   H41
> N4C5
> C5   H51
> C5   H52
> C5C6

[gmx-users] g_dipole ? => Calculate bond dipole moment for a molecule of multiple atoms by hand?

[Chih-Ying Lin]

Hi
According to the following website,

http://en.wikipedia.org/wiki/Bond_dipole_moment


[image: \mu = \delta \, d].
The bond dipole is modeled as +δ -- δ- with a distance *d* between the partial
charges  +δ and δ-.
For a complete molecule the total molecular dipole moment may be
approximated as the vector sum of individual bond dipole moments.


However, for a molecule of multiple atoms,
There may be more than one bond connected on one atom.
 E
|
B - A - C
partial charge of atom_A = -0.5
partial charge of atom_B = 0.2
partial charge of atom_C = 0.35
partial charge of atom_E = 0.4



Which partial charges should I use when I calculate bond-dipole-moment of
A-B ?
Which partial charges should I use when I calculate bond-dipole-moment of
A-C ?
Which partial charges should I use when I calculate bond-dipole-moment of
A-E ?

Thank you
Lin








On 2010-10-18 03.30, Chih-Ying Lin wrote:
> HI
> I confined one molecule in the center of box and issue the g_dipole
command.
> The average dipole moment is still around 32.
> It is the molecule with 33 atoms / united atoms of most carbon groups,
> isn't the dipole moment around 32 too high?
> How can I test next  and  know that the dipole moment around 32 is
> acceptable?
By calculating on paper: dipole is 48.0 sum of q_i x_i, therefore if you
have large charge separation you will get a large dipole.

> Thank you
> Lin
> On 2010-10-16 21.36, Chih-Ying Lin wrote:
>  >
>  > Hi
>  > I issue the g_dipole command on Gromacs => And, the following
>  > information is shown.
>  > There are 10 molecules in the selection,
>  > Does the Average =32.1611 refer to the average for a single over the
>  > simulation time?
>  > Or, the Average = 32.1611 summing for all the 10 molecules over the
>  > simulation time?
>  > If the average = 32.1611 for a single molecule with 33 atoms / united
>  > atoms of most carbon groups, isn't the dipole moment too high?
> I think this is the average per molecule. You may need to run trjconv
> -pbc whole, because mdrun may break molecules in two parts, meaning that
> the molecule becomes as big as the box.
>
>  >
>  >
>  >
>  >
>  > What does "will subtract their charge at their center of mass"  this
> mean?
>  > Why "will subtract their charge at their center of mass"  ?
>  >
>  >
>  >
-- 
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] Energy Minimization

[Dallas Warren]

Look at atom 1149, exactly as the error message has noted.

 

Is it, or any of it's neighbours, not quite where it "should" be?

 

Catch ya,

Dr. Dallas Warren

Medicinal Chemistry and Drug Action

Monash Institute of Pharmaceutical Sciences, Monash University
381 Royal Parade, Parkville VIC 3010
dallas.war...@monash.edu

+61 3 9903 9304
-
When the only tool you own is a hammer, every problem begins to resemble
a nail. 

 

From: gmx-users-boun...@gromacs.org
[mailto:gmx-users-boun...@gromacs.org] On Behalf Of Sai Pooja
Sent: Tuesday, 19 October 2010 7:26 AM
To: Discussion list for GROMACS users
Subject: [gmx-users] Energy Minimization

 

Hi all,

 

I am using CHARMM forcefield with tables for Protein SOL interactions
for alanine dipeptide. 

Tables supplied: table.xvg, table_Protein_SOL.xvg, tablep.xvg (Standard
tables for 6-12 interactions used)

Combination rule changed from '2' to '1' in forcefield.itp file.

 

In the energy minimization step, using mdrun the following problem is
encountered:

 

Polak-Ribiere Conjugate Gradients:
   Tolerance (Fmax)   =  1.0e+00
   Number of steps=1
   F-max =  2.98523e+10 on atom 4
   F-Norm=  1.32072e+09


step -1: Water molecule starting at atom 1149 can not be settled.
Check for bad contacts and/or reduce the timestep if appropriate.
Wrote pdb files with previous and current coordinates
title   = Energy Minimization   ; Title of run

; The following line tell the program the standard locations where to
find certain files
cpp = /lib/cpp  ; Preprocessor

; Define can be used to control processes
;define  = -DFLEXIBLE
define  = -DPOSRES

; Parameters describing what to do, when to stop and what to save
integrator  = cg; Algorithm (steep = steepest descent
minimization)
emtol   = 1.0   ; Stop minimization when the maximum
force < 1.0 kJ/mol
nsteps  = 1 ; Maximum number of (minimization) steps
to perform
nstenergy   = 1 ; Write energies to disk every nstenergy
steps
energygrps  = Protein SOL
energygrp_table = Protein SOL

; Parameters describing how to find the neighbors of each atom and how
to calculate the interactions
ns_type = grid  ; Method to determine neighbor list
(simple, grid)
coulombtype = User  ; Treatment of long range electrostatic
interactions
rcoulomb= 1.0   ; long range electrostatic cut-off
rvdw= 1.0   ; long range Van der Waals cut-off
constraints = none  ; Bond types to replace by constraints
pbc = xyz   ; Periodic Boundary Conditions (yes/no)



 

Any suggestions?

 

Pooja

 


-- 
Quaerendo Invenietis-Seek and you shall discover.

-- 
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] editconf

[Payman Pirzadeh]

Dear Mark and Justin,
Thanks for very informative comments. I'll try to rotate and then solvate.

Paymon

-Original Message-
From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org]
On Behalf Of Mark Abraham
Sent: October 18, 2010 7:26 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] editconf

  On 19/10/2010 12:22 PM, Payman Pirzadeh wrote:
> After rotation, wouldn't a minimization with very very small steps solve
the
> problem?

Only if the solute surface was very close to being symmetric about the 
rotation axis. Since that never happens, editconf doesn't even have the 
ability to try.

Mark
> Paymon
>
> -Original Message-
> From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org]
> On Behalf Of Mark Abraham
> Sent: October 18, 2010 7:14 PM
> To: Discussion list for GROMACS users
> Subject: Re: [gmx-users] editconf
>
>On 19/10/2010 8:11 AM, Paymon Pirzadeh wrote:
>> Hello,
>> I am trying to rotate my protein in my simulation box (solvent molecules
>> are there as well). I issue the following command:
>>editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
>> 180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328
>>
>> I select my group (protein) in the first group selection prompt and
>> whole system as output in the second prompt. However, in the output file
>> the whole system is rotated rather than protein only. Also, despite the
>> -box switch used, the dimensions of the box changes.
>> Is there any ideas on what I am missing in my procedure?
> editconf can't rotate the solute with solvent present - what would you
> like it to do about the voids and clashes that would be created? Do your
> rotation before you solvate.
>
> 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] editconf

[Mark Abraham]

 On 19/10/2010 12:22 PM, Payman Pirzadeh wrote:

After rotation, wouldn't a minimization with very very small steps solve the
problem?


Only if the solute surface was very close to being symmetric about the 
rotation axis. Since that never happens, editconf doesn't even have the 
ability to try.


Mark

Paymon

-Original Message-
From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org]
On Behalf Of Mark Abraham
Sent: October 18, 2010 7:14 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] editconf

   On 19/10/2010 8:11 AM, Paymon Pirzadeh wrote:

Hello,
I am trying to rotate my protein in my simulation box (solvent molecules
are there as well). I issue the following command:
   editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328

I select my group (protein) in the first group selection prompt and
whole system as output in the second prompt. However, in the output file
the whole system is rotated rather than protein only. Also, despite the
-box switch used, the dimensions of the box changes.
Is there any ideas on what I am missing in my procedure?

editconf can't rotate the solute with solvent present - what would you
like it to do about the voids and clashes that would be created? Do your
rotation before you solvate.

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] editconf

[Justin A. Lemkul]

Payman Pirzadeh wrote:

After rotation, wouldn't a minimization with very very small steps solve the
problem?



If any atoms overlap, you can't count on this working.  It would be hit or miss. 
 The proper way to build the system would be to orient the solute, then fill 
the box with solvent.  Far more reliable.


-Justin


Paymon

-Original Message-
From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org]
On Behalf Of Mark Abraham
Sent: October 18, 2010 7:14 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] editconf

  On 19/10/2010 8:11 AM, Paymon Pirzadeh wrote:

Hello,
I am trying to rotate my protein in my simulation box (solvent molecules
are there as well). I issue the following command:
  editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328

I select my group (protein) in the first group selection prompt and
whole system as output in the second prompt. However, in the output file
the whole system is rotated rather than protein only. Also, despite the
-box switch used, the dimensions of the box changes.
Is there any ideas on what I am missing in my procedure?


editconf can't rotate the solute with solvent present - what would you 
like it to do about the voids and clashes that would be created? Do your 
rotation before you solvate.


Mark


--


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] editconf

[Payman Pirzadeh]

After rotation, wouldn't a minimization with very very small steps solve the
problem?

Paymon

-Original Message-
From: gmx-users-boun...@gromacs.org [mailto:gmx-users-boun...@gromacs.org]
On Behalf Of Mark Abraham
Sent: October 18, 2010 7:14 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] editconf

  On 19/10/2010 8:11 AM, Paymon Pirzadeh wrote:
> Hello,
> I am trying to rotate my protein in my simulation box (solvent molecules
> are there as well). I issue the following command:
>   editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
> 180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328
>
> I select my group (protein) in the first group selection prompt and
> whole system as output in the second prompt. However, in the output file
> the whole system is rotated rather than protein only. Also, despite the
> -box switch used, the dimensions of the box changes.
> Is there any ideas on what I am missing in my procedure?

editconf can't rotate the solute with solvent present - what would you 
like it to do about the voids and clashes that would be created? Do your 
rotation before you solvate.

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] editconf

[Mark Abraham]

 On 19/10/2010 8:11 AM, Paymon Pirzadeh wrote:

Hello,
I am trying to rotate my protein in my simulation box (solvent molecules
are there as well). I issue the following command:
  editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328

I select my group (protein) in the first group selection prompt and
whole system as output in the second prompt. However, in the output file
the whole system is rotated rather than protein only. Also, despite the
-box switch used, the dimensions of the box changes.
Is there any ideas on what I am missing in my procedure?


editconf can't rotate the solute with solvent present - what would you 
like it to do about the voids and clashes that would be created? Do your 
rotation before you solvate.


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] g_dipole ? => dipole moment ?

[Chih-Ying Lin]

Hi

>From David =>
"dipole is 48.0 sum of q_i x_i, therefore if you
have large charge separation you will get a large dipole."


1. What does the coefficient 48.0 represent?
2. Is q_i = partial charge or charge on the atoms / united atoms?
3. Is x_i = the distance between the mass center of atom from the mass
center of the molecule?
4. What does the large charge separation mean?  Do you mean the charged
molecule?  Or, Do you mean the molecule with a long carbon chain?


Thank you
Lin






On 2010-10-18 03.30, Chih-Ying Lin wrote:
> HI
> I confined one molecule in the center of box and issue the g_dipole
command.
> The average dipole moment is still around 32.
> It is the molecule with 33 atoms / united atoms of most carbon groups,
> isn't the dipole moment around 32 too high?
> How can I test next  and  know that the dipole moment around 32 is
> acceptable?
By calculating on paper: dipole is 48.0 sum of q_i x_i, therefore if you
have large charge separation you will get a large dipole.

> Thank you
> Lin
> On 2010-10-16 21.36, Chih-Ying Lin wrote:
>  >
>  > Hi
>  > I issue the g_dipole command on Gromacs => And, the following
>  > information is shown.
>  > There are 10 molecules in the selection,
>  > Does the Average =32.1611 refer to the average for a single over the
>  > simulation time?
>  > Or, the Average = 32.1611 summing for all the 10 molecules over the
>  > simulation time?
>  > If the average = 32.1611 for a single molecule with 33 atoms / united
>  > atoms of most carbon groups, isn't the dipole moment too high?
> I think this is the average per molecule. You may need to run trjconv
> -pbc whole, because mdrun may break molecules in two parts, meaning that
> the molecule becomes as big as the box.
>
>  >
>  >
>  >
>  >
>  > What does "will subtract their charge at their center of mass"  this
> mean?
>  > Why "will subtract their charge at their center of mass"  ?
>  >
>  >
>  >
-- 
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] editconf

[Paymon Pirzadeh]

Hello,
I am trying to rotate my protein in my simulation box (solvent molecules
are there as well). I issue the following command:
 editconf_d_mpi -f AFPIII_Ih0001_54_26_em.gro -n protein.ndx -rotate 15
180 0 -o AFPIII_Ih0001_54_26_rotated.gro -box 5.40022 6.23564 20.59328

I select my group (protein) in the first group selection prompt and
whole system as output in the second prompt. However, in the output file
the whole system is rotated rather than protein only. Also, despite the
-box switch used, the dimensions of the box changes.
Is there any ideas on what I am missing in my procedure?
Regards,

Paymon  

-- 
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] gromacs installation problem

[MoJie Duan]

Hi everyone,  I met a problem in the installation of gromacs-4.5. There isn’t any mistake in the configuration step, but in the “make” step, it exit with the error:  make[4]: *** No rule to make target `/usr/incluDe/gnu/stubs-64.h', needed by `sgetri.lo'.  Stop. make[4]: Leaving directory `/home/mduan/soft/gromacs-4.5.1/src/gmxlib/gmx_lapack'   I can find the ‘stubs-64.h’ file in the directory ‘/usr/include/gnu on my computer’. Anyone know how to fix this?  Thanks, Mojie  

-- 
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] Energy Minimization

[Sai Pooja]

Hi all,

I am using CHARMM forcefield with tables for Protein SOL interactions for
alanine dipeptide.
Tables supplied: table.xvg, table_Protein_SOL.xvg, tablep.xvg (Standard
tables for 6-12 interactions used)
Combination rule changed from '2' to '1' in forcefield.itp file.

In the energy minimization step, using mdrun the following problem is
encountered:

Polak-Ribiere Conjugate Gradients:
   Tolerance (Fmax)   =  1.0e+00
   Number of steps=1
   F-max =  2.98523e+10 on atom 4
   F-Norm=  1.32072e+09

step -1: Water molecule starting at atom 1149 can not be settled.
Check for bad contacts and/or reduce the timestep if appropriate.
Wrote pdb files with previous and current coordinates
title   = Energy Minimization   ; Title of run
; The following line tell the program the standard locations where to find
certain files
cpp = /lib/cpp  ; Preprocessor
; Define can be used to control processes
;define  = -DFLEXIBLE
define  = -DPOSRES
; Parameters describing what to do, when to stop and what to save
integrator  = cg; Algorithm (steep = steepest descent
minimization)
emtol   = 1.0   ; Stop minimization when the maximum force <
1.0 kJ/mol
nsteps  = 1 ; Maximum number of (minimization) steps to
perform
nstenergy   = 1 ; Write energies to disk every nstenergy
steps
energygrps  = Protein SOL
energygrp_table = Protein SOL
; Parameters describing how to find the neighbors of each atom and how to
calculate the interactions
ns_type = grid  ; Method to determine neighbor list (simple,
grid)
coulombtype = User  ; Treatment of long range electrostatic interactions
rcoulomb= 1.0   ; long range electrostatic cut-off
rvdw= 1.0   ; long range Van der Waals cut-off
constraints = none  ; Bond types to replace by constraints
pbc = xyz   ; Periodic Boundary Conditions (yes/no)


Any suggestions?

Pooja


-- 
Quaerendo Invenietis-Seek and you shall discover.
-- 
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] Re: PRODRG server

[Lucio Ricardo Montero Valenzuela]

The message says that you have atoms with the same name. Give each atom an
unique name and try running it again. If you have still trouble, try specifying
explicitly the protonation states.
Hope it helps.
Lucio.
Mensaje citado por "Justin A. Lemkul" :

>
>
> mohsen ramezanpour wrote:
> > Dear justin
> > no,I had not done it.
> > But this two warning are not present in other my molecules,for example I
> > tested another molecule but the massage was:
> >
> > defaults to zero instead of generating an error
> >
> > PRODRG>
> > PRODRG>
> > PRODRG> Starting up PRODRG version 071121.0636
> > PRODRG> PRODRG written/copyrighted by Daan van Aalten
> > PRODRG> and Alexander Schuettelkopf
> > PRODRG>
> > PRODRG> Questions/comments to d...@davapc1.bioch.dundee.ac.uk
> 
> >
> > PRODRG>
> > PRODRG> When using this software in a publication, cite:
> > PRODRG> A. W. Schuettelkopf and D. M. F. van Aalten (2004).
> > PRODRG> PRODRG - a tool for high-throughput crystallography
> > PRODRG> of protein-ligand complexes.
> >
> > PRODRG> Acta Crystallogr. D60, 1355--1363.
> > PRODRG>
> > PRODRG>
> > PRODRG> PDB mode detected.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> >
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> >
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> >
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> >
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> > PRODRG> WARNING: deleted hydrogen  H   from your input.
> >
> > PRODRG> WARNING: duplicate atom name F   .
> > PRODRG> WARNING: duplicate atom name C   .
> > PRODRG> WARNING: atoms with same name found. Auto-renaming.
> > PRODRG> Molecule complexity index: 2.09.
> > PRODRG>   2 explicit hydrogen(s) added.
> >
> > PRODRG>  25 bonds4 ambiguous
> > PRODRG>  36 bond angles 20 ambiguous
> > PRODRG>  18 improper dihedrals   0 ambiguous
> > PRODRG>   7 dihedrals2 ambiguous
> > PRODRG>   7 partial charges  4 ambiguous
> >
> > PRODRG> Net charge on molecule:   0.000
> > PRODRG> Using charge groups.
> > PRODRG> Writing GROMACS topology.
> > PRODRG> GROMACS topology quality on 0-10 scale:  5.7
> > GENDRG> Best structure was iteration  931 with   7.66034174
> >
> > PRODRG> Spawning GROMACS version 3.3.3...
> > PRODRG> RMSD from GROMOS bond ideality (Angstrom) :   0.019
> > PRODRG> RMSD from GROMOS angle ideality (degrees) :   2.332
> > PRODRG> RMSD from GROMOS plane ideality (degrees) :   1.169
> >
> > PRODRG> Number of improper improper dihedrals :   0
> > PRODRG> RMSD from starting bonds (Angstrom)  :   0.036
> > PRODRG> RMSD from starting angles (degrees)  :   3.155
> > PRODRG> RMSD from starting planes (degrees)  :   1.494
> >
> > PRODRG> RMSD from starting coords (Angstrom) :   1.930
> > PRODRG> Writing: SCRHWMMPG
> > Unfortunately PRODRG crashed
> >
> > beside I have done these 6 month ago.exactly with this molecules and I
> could recieve zip files from server.
> >
> > but I can not do the same action with the same files?what do you think?
> >
>
> Sounds like the PRODRG server is having issues.  Try contacting the PRODRG
> developers directly.  You won't be able to get sufficient support on this
> forum.
>
> -Justin
>
> >
> > PRODRG> Normal program end
> >
> >
> >
> > On Sat, Oct 16, 2010 at 4:04 PM, Justin A. Lemkul  > > wrote:
> >
> >
> > These error messages are likely more important:
> >
> >
> > PRODRG> WARNING: multiplicity of generated molecule is not 1.
> > PRODRG> WARNING: bond type assignment failed at C13
> >
> > There is something about the chemical structure involving this atom
> > that PRODRG can't handle.  Based on the information you've provided,
> > however, there's no real way anyone can give you much advice.  Have
> > you tried drawing the molecule with the JME editor provided by
> > PRODRG instead of using a .pdb file?
> >
> > -Justin
> >
> > mohsen ramezanpour wrote:
> >
> > actually the final massage was :Unfortunately PRODRG crashed
> >
> > On Sat, Oct 16, 2010 at 11:36 AM, mohsen ramezanpour
> >  > 
> > 

Re: [gmx-users] position restraints and out of bounds atom index

[Alexandre Suman de Araujo]

Since anybody helped me with this issue and I solved the problem, I'll 
report the solution here to contribute with the comunity.


As I pointed in the original e-mail, the problem was the difference 
between the numbering scheme of .gro file and .itp file.


I created 3 .pdb files, one for each chain and with pdb2gmx I created 3 
.gro files. With these 3 .gro files I created 3 .ndx files with the 
groups I'd like to apply PR and from these .ndx files I created the PR 
.itp files, now with the numbering agreeing with the one of the topology 
files.
The include of these PR.itp files in .top file is the same I mentioned 
in the original e-mail.


With this, everything worked fine.

Cheers

**
Alexandre Suman de Araujo*
Faculdade de Ciências Farmacêuticas de Ribeirão Preto*
Universidade de São Paulo*
Dep. de Física e Química *
Grupo de Física Biológica * e-mail: asara...@fcfrp.usp.br*
Av. do Café, s/n° * e-mail: ale.su...@gmail.com  *
CEP: 14040-903* Phone: +55 (16) 3602-4172*
Ribeirão Preto, SP, Brasil* Phone: +55 (16) 3602-4222*
**


Em 16-10-2010 15:41, Alexandre Suman de Araujo escreveu:

Hi all

I'm trying to simulate a system with 3 chains (A, H and L).

Pb2gmx gave me 3 .itp files (one for each chain), each one with atom 
number starting from 1, and a .top file with these 3 .itp included. It 
gave me also a .gro file with the 3 chains numbered consecutively (the 
atom number does not restart from 1 for each chain).


With this .gro file I made an index file using make_ndx with 3 groups 
containing the backbone of each chain. With this index file I built 3 
position restraint files with genrestr, one for each chain.


Following the instructions from website, I removed the #ifdef 
POSRES/#include/#endif lines from the end of the .itp file of each 
chain and added them directly in the .top file as shown below:



;
;File 'proteina.gro.top' was generated
;By user: onbekend (0)
;On host: onbekend
;At date: Thu Oct 14 14:22:01 2010
;
;This is a standalone topology file
;
;It was generated using program:
;pdb2gmx - VERSION 4.5.1
;
;Command line was:
;pdb2gmx -f 2R29CL_original.pdb -o protein.gro -p proteina.gro 
-his -ignh

;
;Force field was read from the standard Gromacs share directory.
;

; Include forcefield parameters
#include "oplsaa.ff/forcefield.itp"

; Include chain topologies
#include "proteina.gro_Protein_chain_A.itp"
; Include Position restraint file for chain A
#ifdef POSRES
#include "posre_Protein_chain_A.itp"
#endif

#include "proteina.gro_Protein_chain_H.itp"
; Include Position restraint file for chain H
#ifdef POSRES
#include "posre_Protein_chain_H.itp"
#endif

#include "proteina.gro_Protein_chain_L.itp"
; Include Position restraint file for chain L
#ifdef POSRES
#include "posre_Protein_chain_L.itp"
#endif

; Include water topology
#include "oplsaa.ff/tip3p.itp"

#ifdef POSRES_WATER
; Position restraint for each water oxygen
[ position_restraints ]
;  i funct   fcxfcyfcz
   11   1000   1000   1000
#endif

; Include topology for ions
#include "oplsaa.ff/ions.itp"

[ system ]
; Name
ENVELOPE PROTEIN E; HEAVY CHAIN OF FAB 1A1D-2; LIGHT CHAIN OF FAB 
1A1D-2 in water


[ molecules ]
; Compound#mols
Protein_chain_A 1
Protein_chain_H 1
Protein_chain_L 1
SOL59
SOL97
SOL95
SOL 30308
NA  91
CL  97


When I run a simulation without the "define = -DPOSRES" definition in 
the .mdp file, grompp runs ok. However, when I say in .mdp file to use 
position restraints grompp give me the following error:


Fatal error:
[ file posre_Protein_chain_H.itp, line 334 ]:
Atom index (3221) in position_restraints out of bounds (1-3212).
This probably means that you have inserted topology section 
"position_restraints"

in a part belonging to a different molecule than you intended to.
In that case move the "position_restraints" section to the right 
molecule.
For more information and tips for troubleshooting, please check the 
GROMACS

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

I know that this was already discussed in list but nothing described 
there worked.


It is clear that the problem is due to the different numbering scheme 
of the .gro file and .itp files. The .gro file does not restart the 
numbering scheme when a new chain starts, on the other hand, the .itp 
file of each chain starts from 1 in atom number field. Since the 
position restraint file is generated based on numbering scheme of the 
.gro file, in some cases (my case is one of them) the PR file refers 
to a atom number that does not exists in .itp file.


I think the way to correct this was generate

RE: [gmx-users] Intrinsically disordered proteins

[#ZHAO LINA#]

Hi, 

1. Thanks for the webpage you provided, I noticed there is a avogadro there, 
and  right now coincidently I have a question relevant to the avogadro,

When I tried the avogadro-packmol-extensions, I had installed the avogadro and 
packmol well, 
but I do not know

where the avogadro plugins/modules can be put ?
I tried the put the packmolextension.so under usr/lib/avogadro/1_0/extensions/ 
but seems not work,

Might someone has tried and knew the answer, actually  I had asked on the other 
place, but none feedback.

2.  Wow, it's going to be lots of works. 

lina 

From: gmx-users-boun...@gromacs.org [gmx-users-boun...@gromacs.org] on behalf 
of Justin A. Lemkul [jalem...@vt.edu]
Sent: Monday, October 18, 2010 9:37 PM
To: Discussion list for GROMACS users
Subject: Re: [gmx-users] Intrinsically disordered proteins

#ZHAO LINA# wrote:
> Hi,
>
> 1. For those intrinsically disordered proteins, the sequence is known,
> how the simulations will be set up, I mean, the first PDB will be
> needed, how to get this one? (Ideally, not necessarily to be practical)
>

http://www.gromacs.org/Documentation/File_Formats/Coordinate_File#Sources

> 2. suppose I got a PDB, there were several models there, let's say 16,
> is it acceptable to just take one, say model 1?
>

Making an arbitrary choice is not a good idea.  You have to justify that the
model you selected is in some way representative of the ensemble.  You could do
that by demonstrating that your simulation overlapped the configuration space of
the ensemble in the PDB file, i.e. sampled all the expected configurations, thus
the starting structure was not an inappropriate bias.

-Justin

> Thanks for any answering,
>
> lina
>

--


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] Simulation parameter problem about protein unfolding

[Chen]

This information is helpful!

I think I need to deal with the simulation result really carefully. Thanks for 
the help! And also thanks Massimo for the information about expected gyration 
radius of a purely non-interacting random chain!


At 2010-10-18 21:01:23，chris.ne...@utoronto.ca wrote:

>Chen,
>
>to answer this, please see David's post and my reply. To answer if  
>your unfolded state ensemble is accurate, you will need to know some  
>experimentally evaluated quantities of the state that you are trying  
>to reproduce: 498K, 26 atm, a certain amount of salt, etc, (but not  
>including the forcefield model or cutoffs, which are obviously  
>intrinsically part of the model).
>
>Then, if your data does not match this, then you can conclude that the  
>simulation time was too short or the parameters were not valid. Also,  
>I guarantee it was... you'll probably need something more like 10 ms -  
>10 s (what is the unfolding rate under these conditions?), which is  
>rather difficult to obtain).
>
>Comparing a temperature denatured state to an AFM denatured state, one  
>does not expect them to match. The point here was to ask: are you sure  
>that your ensemble of conformations is incorrect given the conditions?
>
>Chris.
>
>-- original message --
>
>> At 2010-10-18 12:56:31£¬chris.neale at utoronto.ca wrote:
>> Generally, forcefields are not parameterized for temperatures other   
>> than 298K, so simulations are not expected to reproduce the expected  
>>  properties (like boiling water and the correct temperature   
>> denaturation of proteins).
>>
>> There's almost certainly other issues here (including the fact that   
>> I'm entirely sure that you can get a lot more than 24 ns of  
>> simulation  on a 54 aa protein; and 26 atom of pressure seems pretty  
>> arbitrary)  but it will come down to this eventually.
>>
>> Just because you found a paper in which they get a denatured state   
>> does not imply that they got the correct denatured state.
>>
>> Chris.
>>
>
>Hi, Chris! Thanks for the reply! I have not conducted unfolding
>before, so I compared my result with other's to make sure it is reasonable in
>some extend. The 26 atm pressure comes from experiment result (Haar et al.,
>1984) mentioned in some MD related papers (e.g. 'J. Mol. Biol. (2000) 296,
>1257-1282').
>
>  I've searched the maillist
>and noticed some issues about High temperature simulations. I'm not sure
>whether the 'ilmm' force field has been optimized for high temperature
>simulation. I also noticed some users asked about MD parameters in  
>'unfolding a
>protein'. And the parameters they used are different from ours' (e.g. the
>'rlist', 'rcoulomb' and 'temperature or pressure couple algorithm').
>
>  I just want to make sure I
>didn't make mistakes in these parameters which maybe cause the protein keeping
>in a compact state! The radius of gyration of the protein fluctuated around
>1.1nm (never bigger than 1.4nm) during our 30ns simulation now. If the MD
>parameters I chose have no problem, then I could keep going. Any comment?
>
>
>
>--
>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] Intrinsically disordered proteins

[Justin A. Lemkul]

#ZHAO LINA# wrote:

Hi,

1. For those intrinsically disordered proteins, the sequence is known, 
how the simulations will be set up, I mean, the first PDB will be 
needed, how to get this one? (Ideally, not necessarily to be practical)




http://www.gromacs.org/Documentation/File_Formats/Coordinate_File#Sources

2. suppose I got a PDB, there were several models there, let's say 16, 
is it acceptable to just take one, say model 1?




Making an arbitrary choice is not a good idea.  You have to justify that the 
model you selected is in some way representative of the ensemble.  You could do 
that by demonstrating that your simulation overlapped the configuration space of 
the ensemble in the PDB file, i.e. sampled all the expected configurations, thus 
the starting structure was not an inappropriate bias.


-Justin


Thanks for any answering,

lina



--


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] Intrinsically disordered proteins

[ms]

On 18/10/10 14:26, #ZHAO LINA# wrote:

Hi,

1. For those intrinsically disordered proteins, the sequence is known, how the 
simulations will be set up, I mean, the first PDB will be needed, how to get 
this one? (Ideally, not necessarily to be practical)

2. suppose I got a PDB, there were several models there, let's say 16, is it 
acceptable to just take one, say model 1?

Thanks for any answering,


That's a theoretical question which has no real "true" answer, given 
that there is no "true" structure for an IDP, but only an ensemble of 
structures.


I've seen different approaches, from using ROSETTA predictions to simply 
using random coil structures, to carefully picking ensembles from 
previous experimentally-restrained simulations. Pick up your favourite 
approach in the literature, also taking into account what you want to do.


m.

--
Massimo Sandal, Ph.D.
http://devicerandom.org
--
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] Intrinsically disordered proteins

[#ZHAO LINA#]

Hi,

1. For those intrinsically disordered proteins, the sequence is known, how the 
simulations will be set up, I mean, the first PDB will be needed, how to get 
this one? (Ideally, not necessarily to be practical)

2. suppose I got a PDB, there were several models there, let's say 16, is it 
acceptable to just take one, say model 1?

Thanks for any answering,

lina

-- 
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] Simulation parameter problem about protein unfolding

[ms]

On 18/10/10 13:12, Chen wrote:



At 2010-10-18 16:38:30??"David van der Spoel"  wrote:



On 2010-10-18 06.56, chris.ne...@utoronto.ca wrote:
Generally, forcefields are not parameterized for temperatures other than
298K, so simulations are not expected to reproduce the expected
properties (like boiling water and the correct temperature denaturation
of proteins).

There's almost certainly other issues here (including the fact that I'm
entirely sure that you can get a lot more than 24 ns of simulation on a
54 aa protein; and 26 atom of pressure seems pretty arbitrary) but it
will come down to this eventually.

Just because you found a paper in which they get a denatured state does
not imply that they got the correct denatured state.


There is no correct denatured state. There are infinitely many. Check
out recent work on NMR of "unfolded" proteins.



I thought about the unfolding state space is huge. I just wonder
whether the relative low radius of gyration value space sampled by us is caused
by some error setting in MD parameters. If no one here finds there's problem in
my MD parameters, then I can keep going. Thanks!


In theory, you should be able to calculate the expected gyration radius 
of a purely non-interacting random chain (which one can take as the 
limit of an unfolded protein). given the distance between monomeric 
units (i.e. C-alpha for example here). It is a three-dimensional random 
walk.

http://en.wikipedia.org/wiki/Radius_of_gyration#Molecular_applications

Doing the calculation for the equation of the three-dimensional random 
walk of your 54 aa protein, with N=54 and a=0.35 nm, it turns out that 
the expected gyration radius is:


>>> (np.sqrt(54) * 0.35) * (1/np.sqrt(6))
1.0503 nm

So, yes, your average protein gyration radius is consistent with the 
protein being a random coil -that is, unfolded.


In truth your random coil ball should be a bit larger because in the 
calculation we ignore the chain self-avoidance, but on the other hand we 
have that a protein chain is not completely non-self-interacting, for 
obvious reasons. Therefore it seems you can conclude that your protein 
is reasonably close to being unfolded.


Hope it helps,
Massimo

--
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] Backing up trajectories

[Esztermann, Ansgar]

Hi everyone,

I am just wondering how backups of Gromacs trajectories are handled "out 
there". 
We are using TSM. The system is very reliable, but it does not scale well in 
conjunction with Gromacs: trajectories can get very large when -append is used 
(dozens of Gigabytes), but TSM cannot do block-incremental backups (at least on 
Linux).
So every time a trajectory is backed up, the whole xx GB are copied to the 
backup server. There are mechanisms in place to limit .trr/.xtc backups to once 
a week or so, but for runs lasting several months, this still means a huge 
amount of wasted network bandwidth and server disk space.
Since computational power increases faster than network bandwidth does, this is 
bound to become a bottleneck sooner or later.

Is anyone willing to share how backups are handled at your site?


Thanks,

A.

-- 
Ansgar Esztermann
DV-Systemadministration
Max-Planck-Institut für biophysikalische Chemie, Abteilung 105

--
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] Simulation parameter problem about protein unfolding

[chris . neale]

Chen,

to answer this, please see David's post and my reply. To answer if  
your unfolded state ensemble is accurate, you will need to know some  
experimentally evaluated quantities of the state that you are trying  
to reproduce: 498K, 26 atm, a certain amount of salt, etc, (but not  
including the forcefield model or cutoffs, which are obviously  
intrinsically part of the model).


Then, if your data does not match this, then you can conclude that the  
simulation time was too short or the parameters were not valid. Also,  
I guarantee it was... you'll probably need something more like 10 ms -  
10 s (what is the unfolding rate under these conditions?), which is  
rather difficult to obtain).


Comparing a temperature denatured state to an AFM denatured state, one  
does not expect them to match. The point here was to ask: are you sure  
that your ensemble of conformations is incorrect given the conditions?


Chris.

-- original message --


At 2010-10-18 12:56:31£¬chris.neale at utoronto.ca wrote:
Generally, forcefields are not parameterized for temperatures other   
than 298K, so simulations are not expected to reproduce the expected  
 properties (like boiling water and the correct temperature   
denaturation of proteins).


There's almost certainly other issues here (including the fact that   
I'm entirely sure that you can get a lot more than 24 ns of  
simulation  on a 54 aa protein; and 26 atom of pressure seems pretty  
arbitrary)  but it will come down to this eventually.


Just because you found a paper in which they get a denatured state   
does not imply that they got the correct denatured state.


Chris.



Hi, Chris! Thanks for the reply! I have not conducted unfolding
before, so I compared my result with other's to make sure it is reasonable in
some extend. The 26 atm pressure comes from experiment result (Haar et al.,
1984) mentioned in some MD related papers (e.g. 'J. Mol. Biol. (2000) 296,
1257-1282').

 I've searched the maillist
and noticed some issues about High temperature simulations. I'm not sure
whether the 'ilmm' force field has been optimized for high temperature
simulation. I also noticed some users asked about MD parameters in  
'unfolding a

protein'. And the parameters they used are different from ours' (e.g. the
'rlist', 'rcoulomb' and 'temperature or pressure couple algorithm').

 I just want to make sure I
didn't make mistakes in these parameters which maybe cause the protein keeping
in a compact state! The radius of gyration of the protein fluctuated around
1.1nm (never bigger than 1.4nm) during our 30ns simulation now. If the MD
parameters I chose have no problem, then I could keep going. Any comment?



--
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] Simulation parameter problem about protein unfolding

[chris . neale]

Thanks David. I think that we're saying the same thing, but there is  
indeed a 'correct' denatured state at a given temperature and pressure  
and protein concentration and salt concentration, etc. It was this to  
which I was referring.


Chris.

On 2010-10-18 06.56, chris.neale at utoronto.ca wrote:

Generally, forcefields are not parameterized for temperatures other than
298K, so simulations are not expected to reproduce the expected
properties (like boiling water and the correct temperature denaturation
of proteins).

There's almost certainly other issues here (including the fact that I'm
entirely sure that you can get a lot more than 24 ns of simulation on a
54 aa protein; and 26 atom of pressure seems pretty arbitrary) but it
will come down to this eventually.

Just because you found a paper in which they get a denatured state does
not imply that they got the correct denatured state.


There is no correct denatured state. There are infinitely many. Check
out recent work on NMR of "unfolded" proteins.


Chris.

-- original message --

Hi All,

I met a problem when I try to unfold a protein using Gromacs, It seemed
the protein cannot be totally unfolded!

The simulated system has one Engrailed Homeodomain (En) protein (a three
helix bundle protein with 54 residues, 629 atoms), total 4848 spce
waters, and 7 Cl- used to neutralize the system in a 5.752(nm)^3 water
BOX. I use the NTP ensemble with T=498K and P=26atm. The system has 1nm
thick water in each side of the En protein, and the density of the
system has been adjusted to0.829 g/cm3.

The simulation lasted 24ns. The helixes disappeared at about 4ns. And
after that some beta sheet formed in the N terminal of the protein.
However, the protein kept in a compact state during the 24ns simulation.
The radius of gyration of the protein fluctuated around 1.1nm during the
simulation.

I've also noticed similar simulation done by others. For example, David
Becka and Valerie Daggett reviewed their En protein unfolding in paper
"Methods for molecular dynamics simulations of protein folding/unfolding
in solution. Methods 34 (2004) 112¨C120". The simulations were performed
with ENCAD and ilmm, and used an 0.8nm cutoff range. And the ensemble
they used is NVE as I know. A stretched unfolding state occurred at
about 5ns in their 60ns simulation in 498K, with little helix structure.

I was wondering whether the difference is caused by using different MD
software and force field, or by some wrong parameters in my .mdp file.
I've also conducted another 18ns simulation, and the result is almost
the same. I list he mdp file below. Any comment is appreciated!





--


--
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] Simulation parameter problem about protein unfolding

[Chen]

>At 2010-10-18 16:38:30，"David van der Spoel"  wrote:

>>On 2010-10-18 06.56, chris.ne...@utoronto.ca wrote:
>> Generally, forcefields are not parameterized for temperatures other than
>> 298K, so simulations are not expected to reproduce the expected
>> properties (like boiling water and the correct temperature denaturation
>> of proteins).
>>
>> There's almost certainly other issues here (including the fact that I'm
>> entirely sure that you can get a lot more than 24 ns of simulation on a
>> 54 aa protein; and 26 atom of pressure seems pretty arbitrary) but it
>> will come down to this eventually.
>>
>> Just because you found a paper in which they get a denatured state does
>> not imply that they got the correct denatured state.
>>
>There is no correct denatured state. There are infinitely many. Check 
>out recent work on NMR of "unfolded" proteins.


I thought about the unfolding state space is huge. I just wonder
whether the relative low radius of gyration value space sampled by us is caused
by some error setting in MD parameters. If no one here finds there's problem in
my MD parameters, then I can keep going. Thanks!

At 2010-10-18 16:19:02，"Gerrit Groenhof"  wrote:Hi,

You write it yourself: In paper you mention, they have used a 0.8 nm cutoff 
range for both electrostatics and cutoff. You are doing something different by 
using PME for the electrostatics. Also you are using a much longer cutoff for 
the VDW interactions. 
If you want to reproduce their results, you need to stick to to the parameters 
mentioned in that paper and use cut-offs as well.

BTW I do not believe this to be a good idea though.

Gerrit



I didn't want to simply reproduce their results I've mentioned. I
guess even if I set the same cutoff, the result wouldn't change too much or
would be different for some other reason (e.g. the force field or the MD
software :) ). Anyway I would try it! 

I chose PME for the electrostatics and set longer cutoff for the VDW
in order to get more precise interaction calculation. I've noticed the force
field bias at high temperature and even the couple algorithm would flaw the
protein unfolding. However, that's what I could do by now. I would take these
parameters if no one here finds there's problem. Thanks!



;-
integrator   = md
dt   = 0.002
nsteps   = 300
nstxout  = 500
nstvout  = 500
nstlog   = 500
nstenergy= 500
nstxtcout= 1000

energygrps   = protein non-protein

nstlist  = 10
ns_type  = grid
pbc  = xyz
rlist= 1.0

coulombtype  = PME
rcoulomb = 1.0
vdwtype  = cut-off
rvdw = 1.4
fourierspacing   = 0.12
fourier_nx   = 0
fourier_ny   = 0
fourier_nz   = 0
pme_order= 4
ewald_rtol   = 1e-5
optimize_fft = yes

tcoupl   = v-rescale
tc_grps  = protein non-protein
tau_t= 0.1  0.1
ref_t= 498  498
Pcoupl   = Parrinello-Rahman
pcoupltype   = isotropic
tau_p= 0.5
compressibility  = 4.5e-5
ref_p= 26.0

gen_vel  = yes
gen_temp = 498
gen_seed = 173529

constraints  = hbonds
lincs_order  = 10
;--- 
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] about sa_surface_tension

[Per Larsson]

Hi!

This is unfortunately a bit confusing, but prior to the release of 4.5.2 
gromacs did not use the value of
the sa_surface_tension that was specified in the mdp-file. It was always set to 
2.092, and then hardcode to take a certain
value (other than 2.092) depending on the GB-model. 

2.25936 is the value used together with OBC/HCT models, it is slightly lower 
with the Still model. 
In 4.5.2 the default will be -1, which, if left unchanged, will be changed to 
proper values (and a note written about this).
The effect on the result will be very small, since the 
sa_surface_tension-parameter is only used for the non-polar part
of solvation, which in turn is a very minor part of the total force. 

Also in 4.5.2, setting sa_surface_tension=0 means that no non-polar solvation 
is being calculated. 

Sorry for the confusion
/Per

18 okt 2010 kl. 14.07 skrev Christian Mücksch:

> Dear all,
> 
> when using implicit solvent models, the sa_surface_tension parameter is set 
> to 2.092 kJ/mol/nm2 by default.
> 
> However in Gromacs benchmarks with dihydrofolate reductase this parameter is 
> set to 2.25936.
> 
> Can someone please explain this parameter and in how far this would change 
> the simulation outcome.
> 
> Kind regards,
> 
> -- 
> Christian Mücksch
> Department of Physics
> TU Kaiserslautern
> Erwin Schrödinger Straße
> 67663 Kaiserslautern
> Germany
> 
> Phone: +49 (0)631 205 4287
> Email: mueck...@rhrk.uni-kl.de
> 
> -- 
> 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] Reg: Running Energy minimisation for dichloroethane (DCE)

[Justin A. Lemkul]

vinothkumar mohanakrishnan wrote:

Dear Justin

what corrections i should make to ffoplsaan.itp to make it correct. what 
i should give instead of CAB and CLAA?.




You must use atom types, not names.  Unfortunately, you've chosen to use atom 
names, which are also types, which makes all of this quite confusing if you're 
not sure what you're doing.


You defined two new atom types - opls_966 (CAB) and opls_967 (CLAA), which are 
the only indicators you are allowed to use if introducing new types.  Thus, 
references to atom names (CAC, CLAD) will generate fatal errors.


i copied the .rtp .atp .itp files from 
usr/local/gromacs/share/gromacs/top to my working directory and added 
these parameters to the corresponding files. what you mean is should i 
need to add these parameters to the source directory ( 
usr/local/gromacs/share/gromacs/top)?




Your topology needs to be consistent with whatever files need to be included. 
By default, Gromacs checks the working directory first, but if you've moved to a 
new (sub)directory to carry out further steps, the grompp will not find your 
modified files, but will instead locate only the default force field files in 
$GMXLIB.  Either keep all your work in one directory (which can get messy), or 
make use of the "include" keyword in the .mdp file.  Any directory specified 
there will be searched after the working directory, but before $GMXLIB.


-Justin


Regards
Vinoth

On Mon, Oct 18, 2010 at 5:27 PM, Justin A. Lemkul > wrote:




vinothkumar mohanakrishnan wrote:

Hi all

i added new atomtype for dichloroethane (DCE) and added the
corresponding parameters in the .rtp, .atp, .bon.itp, .nb.itp
respectively. given below are my additions to the corresponding
files respectively.

*ffoplsaa.rtp*

[ DCE ]
 [ atoms ]
 CLAA opls_967 -0.2270  1
 CAB  opls_966  0.2270  1
 CAC  opls_966  0.2270  2
 CLAD opls_967 -0.2270  2

[ bonds ]
 CLAACAB
 CABCAC
 CACCLAD

[ angles ]
 CLAACABCAC
 CABCACCLAD

[ dihedrals ]
 CLAACABCACCLAD

*ffoplsaa.atp*

 opls_966   14.02700  ; CH2 for DCE
 opls_967   35.45300  ; CL for DCE

*ffoplsaabon.itp*

[bondtypes]
CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE

[angletypes]
CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE

[dihedraltypes]
CLAA   CAB CAC   CLAD 3 20.76096  -0.4184  
27.011904  0.0   0.0   0.0 ; for DCE



You shouldn't be using atom names in the [*types] directives.  What
you should be using are the interpolated types from ffoplsaanb.itp,
thus you have only utilized opls_966 (CAB) and opls_967 (CLAA).

*ffoplsaanb.itp*


opls_966   CAB614.027000.227   A3.98000e-01
 4.76976e-01 ; CH2 of DCE
opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01
 0.20920e+01 ; Cl of DCE

*dce.pdb*

HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00
20.00CL
HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00
20.00 C
HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00
20.00 C
HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00
20.00CL

After adding all this, when i run grompp i get the error as
*fatal error Unknown bond_atomtype CLAA*. can any one tell me
why this happens?.


Have you been sure to #include the correct (modified) force field
files?  That is, if you made a local copy and adjusted them, these
won't be the files that pdb2gmx will #include by default.

-Justin

Regards
Vinoth


-- 



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 pos

[gmx-users] about sa_surface_tension

[Christian Mücksch]

 Dear all,

when using implicit solvent models, the sa_surface_tension parameter is 
set to 2.092 kJ/mol/nm2 by default.


However in Gromacs benchmarks with dihydrofolate reductase this 
parameter is set to 2.25936.


Can someone please explain this parameter and in how far this would 
change the simulation outcome.


Kind regards,

--
Christian Mücksch
Department of Physics
TU Kaiserslautern
Erwin Schrödinger Straße
67663 Kaiserslautern
Germany

Phone: +49 (0)631 205 4287
Email: mueck...@rhrk.uni-kl.de

--
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] The trouble with dihedral restraints: frozen peptide backbones

[ms]

On 18/10/10 13:06, Mark Abraham wrote:

On 18/10/2010 10:53 PM, ms wrote:

On 18/10/10 03:30, Mark Abraham wrote:

To enforce chirality in such a toy, I thought that a simple
(naive?) but functional idea could be that of enforcing a not-
too-hard and wide-bottomed dihedral restrain on the phi angle,
like that:

[ dihedral_restraints ]
; ai aj ak
al type label phi dphi kfac power
; phi C'(n-1) - N - CA - C'
3
5 7 8 1 1 -90 75 2 1


Try mdrun -debug on a short run. What does it have to say about the
parameter values for this dihedral restraint? Likewise gmxdump -s on
the .tpr (though you will have to search around for F_DIHRES among
other things, probably).



I have run a 500 ps run with -debug on, but I am unsure of where
should I look for what it says of relevance -I see the ctab,rtab,dtab
xvg's but I doubt they are of interest.


There should be a much larger .log file, probably named differently.
Check what it reports about dihedral restraints - text searching will be
your friend, here.

Mark


Thanks a lot. As far as I can see I see an endless amount of the following:

dihres[0]: 2 4 6 7 : phi=-2.239583, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[1]: 7 9 11 12 : phi=-2.264564, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[2]: 12 14 16 17 : phi=-2.286253, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[3]: 17 19 21 22 : phi=-2.285502, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[4]: 22 24 26 27 : phi=-2.325288, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[5]: 27 29 31 32 : phi=-2.333271, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[6]: 32 34 36 37 : phi=-2.344587, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[7]: 37 39 41 42 : phi=-2.359801, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[8]: 42 44 46 47 : phi=-2.347612, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[9]: 47 49 51 52 : phi=-2.375774, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[10]: 52 54 56 57 : phi=-2.370417, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[11]: 57 59 61 62 : phi=-2.385128, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[12]: 62 64 66 67 : phi=-2.411767, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[13]: 67 69 71 72 : phi=-2.383690, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[14]: 72 74 76 77 : phi=-2.444273, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[15]: 77 79 81 82 : phi=-2.391789, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[16]: 82 84 86 87 : phi=-2.416251, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[17]: 87 89 91 92 : phi=-2.419053, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[18]: 92 94 96 97 : phi=-2.441001, dphi=1.308997, kfac=20.00, 
power=1, label=1
dihres[19]: 97 99 101 102 : phi=-2.437035, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[20]: 102 104 106 107 : phi=-2.405856, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[21]: 107 109 111 112 : phi=-2.430783, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[22]: 112 114 116 117 : phi=-2.405508, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[23]: 117 119 121 122 : phi=-2.386797, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[24]: 122 124 126 127 : phi=-2.388632, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[25]: 127 129 131 132 : phi=-2.419342, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[26]: 132 134 136 137 : phi=-2.396203, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[27]: 137 139 141 142 : phi=-2.409687, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[28]: 142 144 146 147 : phi=-2.372591, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[29]: 147 149 151 152 : phi=-2.376767, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[30]: 152 154 156 157 : phi=-2.355252, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[31]: 157 159 161 162 : phi=-2.356343, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[32]: 162 164 166 167 : phi=-2.344373, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[33]: 167 169 171 172 : phi=-2.337336, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[34]: 172 174 176 177 : phi=-2.299121, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[35]: 177 179 181 182 : phi=-2.294940, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[36]: 182 184 186 187 : phi=-2.297313, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[37]: 187 189 191 192 : phi=-2.287532, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[38]: 192 194 196 197 : phi=-2.255120, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[39]: 197 199 201 202 : phi=-2.221712, dphi=1.308997, 
kfac=20.00, power=1, label=1
dihres[40]: 202 204 206 207 : phi=-2.115131, dphi=1.308997, 
kfac=20.00, power=1, label=1


Now, the 1.308997 value stays fixed and seems to be the radiants value 
of the 75 degrees dphi parameter I've put in the topology. The phi 
value,

Re: [gmx-users] Reg: Running Energy minimisation for dichloroethane (DCE)

[vinothkumar mohanakrishnan]

Dear Justin

what corrections i should make to ffoplsaan.itp to make it correct. what i
should give instead of CAB and CLAA?.

i copied the .rtp .atp .itp files from usr/local/gromacs/share/gromacs/top
to my working directory and added these parameters to the corresponding
files. what you mean is should i need to add these parameters to the source
directory ( usr/local/gromacs/share/gromacs/top)?

Regards
Vinoth

On Mon, Oct 18, 2010 at 5:27 PM, Justin A. Lemkul  wrote:

>
>
> vinothkumar mohanakrishnan wrote:
>
>> Hi all
>>
>> i added new atomtype for dichloroethane (DCE) and added the corresponding
>> parameters in the .rtp, .atp, .bon.itp, .nb.itp respectively. given below
>> are my additions to the corresponding files respectively.
>>
>> *ffoplsaa.rtp*
>>
>> [ DCE ]
>>  [ atoms ]
>>  CLAA opls_967 -0.2270  1
>>  CAB  opls_966  0.2270  1
>>  CAC  opls_966  0.2270  2
>>  CLAD opls_967 -0.2270  2
>>
>> [ bonds ]
>>  CLAACAB
>>  CABCAC
>>  CACCLAD
>>
>> [ angles ]
>>  CLAACABCAC
>>  CABCACCLAD
>>
>> [ dihedrals ]
>>  CLAACABCACCLAD
>>
>> *ffoplsaa.atp*
>>
>>  opls_966   14.02700  ; CH2 for DCE
>>  opls_967   35.45300  ; CL for DCE
>>
>> *ffoplsaabon.itp*
>>
>> [bondtypes]
>> CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
>> CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
>> CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE
>>
>> [angletypes]
>> CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
>> CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE
>>
>> [dihedraltypes]
>> CLAA   CAB CAC   CLAD 3 20.76096  -0.4184   27.011904  0.0
>>   0.0   0.0 ; for DCE
>>
>>
> You shouldn't be using atom names in the [*types] directives.  What you
> should be using are the interpolated types from ffoplsaanb.itp, thus you
> have only utilized opls_966 (CAB) and opls_967 (CLAA).
>
>  *ffoplsaanb.itp*
>>
>>
>> opls_966   CAB614.027000.227   A3.98000e-01
>>  4.76976e-01 ; CH2 of DCE
>> opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01
>>  0.20920e+01 ; Cl of DCE
>>
>> *dce.pdb*
>>
>> HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00 20.00
>>CL
>> HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00 20.00
>> C
>> HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00 20.00
>> C
>> HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00 20.00
>>CL
>>
>> After adding all this, when i run grompp i get the error as *fatal error
>> Unknown bond_atomtype CLAA*. can any one tell me why this happens?.
>>
>>
> Have you been sure to #include the correct (modified) force field files?
>  That is, if you made a local copy and adjusted them, these won't be the
> files that pdb2gmx will #include by default.
>
> -Justin
>
>  Regards
>> Vinoth
>>
>>
> --
> 
>
> 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

Re: [gmx-users] The trouble with dihedral restraints: frozen peptide backbones

[ms]

On 18/10/10 13:06, Mark Abraham wrote:

On 18/10/2010 10:53 PM, ms wrote:

On 18/10/10 03:30, Mark Abraham wrote:

To enforce chirality in such a toy, I thought that a simple
(naive?) but functional idea could be that of enforcing a not-
too-hard and wide-bottomed dihedral restrain on the phi angle,
like that:

[ dihedral_restraints ]
; ai aj ak
al type label phi dphi kfac power
; phi C'(n-1) - N - CA - C'
3
5 7 8 1 1 -90 75 2 1


Try mdrun -debug on a short run. What does it have to say about the
parameter values for this dihedral restraint? Likewise gmxdump -s on
the .tpr (though you will have to search around for F_DIHRES among
other things, probably).



I have run a 500 ps run with -debug on, but I am unsure of where
should I look for what it says of relevance -I see the ctab,rtab,dtab
xvg's but I doubt they are of interest.


There should be a much larger .log file, probably named differently.
Check what it reports about dihedral restraints - text searching will be
your friend, here.


Found - mdrun0.log . It's a 4.2 Gb file (argh), I am now trying to grep 
it and see.


m.

--
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] The trouble with dihedral restraints: frozen peptide backbones

[Mark Abraham]

 On 18/10/2010 10:53 PM, ms wrote:

On 18/10/10 03:30, Mark Abraham wrote:

To enforce chirality in such a toy, I thought that a simple
(naive?) but functional idea could be that of enforcing a not-
too-hard and wide-bottomed dihedral restrain on the phi angle,
like that:

[ dihedral_restraints ]
; ai   ajak
al  type  label  phi  dphi  kfac  power
; phi C'(n-1) - N - CA - C'
   3
5   7   8   11  -90 752  1


Try mdrun -debug on a short run. What does it have to say about the 
parameter values for this dihedral restraint? Likewise gmxdump -s on 
the .tpr (though you will have to search around for F_DIHRES among 
other things, probably).




I have run a 500 ps run with -debug on, but I am unsure of where 
should I look for what it says of relevance -I see the ctab,rtab,dtab 
xvg's but I doubt they are of interest.


There should be a much larger .log file, probably named differently. 
Check what it reports about dihedral restraints - text searching will be 
your friend, here.


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] Re: Amber to gromacs topology conversion for non-amino acid

[Justin A. Lemkul]

Please keep all Gromacs-related correspondence on the gmx-users list.  As I have 
said many times, I am not a private tutor.  I am CC'ing this message to the 
list; please continue all further discussion there.


There are a number of topology-generating scripts on the Gromacs website.  Have 
a look there:


http://www.gromacs.org/Downloads/User_contributions/Other_software

-Justin

swapnil chavan wrote:

Dear Dr. Lemkul,
I want to run MD simulation for ligand-protein complex. I can't use
dundee server for topology files preparation for this ligand. Could
you please tell me another way by which I can generate that topology
and forcefield files on my desk (without submission to online server).
Thank you in advance.



--


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] Reg: Running Energy minimisation for dichloroethane (DCE)

[ms]

On 18/10/10 12:49, vinothkumar mohanakrishnan wrote:

Dear Mark

I generated the topology using pdb2gmx using dce.pdb. i checked the .top
file and i contains necessary .itp files. the error message that i get is
given below. any help is highly appreciated.


What Mark meant is that you should copy-and-paste here the actual 
commands you typed, and the relevant file lines.


m.
--
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] The trouble with dihedral restraints: frozen peptide backbones

[ms]

On 18/10/10 03:30, Mark Abraham wrote:

To enforce chirality in such a toy, I thought that a simple
(naive?) but functional idea could be that of enforcing a not-
too-hard and wide-bottomed dihedral restrain on the phi angle,
like that:

[ dihedral_restraints ]
; ai   ajak
al  type  label  phi  dphi  kfac  power
; phi C'(n-1) - N - CA - C'
   3
5   7   8   11  -90 752  1


Try mdrun -debug on a short run. What does it have to say about the parameter 
values for this dihedral restraint? Likewise gmxdump -s on the .tpr (though you 
will have to search around for F_DIHRES among other things, probably).


In gmxdump -s output I find:
functype[32]=DIHRES, label=1, power=   1 phi=-9.e+01, dphi= 
7.5000e+01, kfac= 2.e+00)


which seems consistent with what I put in the topology. I looked a bit 
in the dump but it's quite difficult for me to understand what it is 
dumping and how it can be relevant.


If you need other info, the full dump, etc., let me know.

thanks!!
m.
--
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] Reg: Running Energy minimisation for dichloroethane (DCE)

[Justin A. Lemkul]

vinothkumar mohanakrishnan wrote:

Hi all

i added new atomtype for dichloroethane (DCE) and added the 
corresponding parameters in the .rtp, .atp, .bon.itp, .nb.itp 
respectively. given below are my additions to the corresponding files 
respectively.


*ffoplsaa.rtp*

[ DCE ]
 [ atoms ]
 CLAA opls_967 -0.2270  1
 CAB  opls_966  0.2270  1
 CAC  opls_966  0.2270  2
 CLAD opls_967 -0.2270  2

[ bonds ]
 CLAACAB
 CABCAC
 CACCLAD

[ angles ]
 CLAACABCAC
 CABCACCLAD

[ dihedrals ]
 CLAACABCACCLAD

*ffoplsaa.atp*

 opls_966   14.02700  ; CH2 for DCE
 opls_967   35.45300  ; CL for DCE

*ffoplsaabon.itp*

[bondtypes]
CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE

[angletypes]
CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE

[dihedraltypes]
CLAA   CAB CAC   CLAD 3 20.76096  -0.4184   27.011904  
0.0   0.0   0.0 ; for DCE




You shouldn't be using atom names in the [*types] directives.  What you should 
be using are the interpolated types from ffoplsaanb.itp, thus you have only 
utilized opls_966 (CAB) and opls_967 (CLAA).



*ffoplsaanb.itp*

opls_966   CAB614.027000.227   A3.98000e-01  
4.76976e-01 ; CH2 of DCE
opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01  
0.20920e+01 ; Cl of DCE


*dce.pdb*

HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00 
20.00CL
HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00 
20.00 C
HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00 
20.00 C
HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00 
20.00CL


After adding all this, when i run grompp i get the error as *fatal error 
Unknown bond_atomtype CLAA*. can any one tell me why this happens?.




Have you been sure to #include the correct (modified) force field files?  That 
is, if you made a local copy and adjusted them, these won't be the files that 
pdb2gmx will #include by default.


-Justin


Regards
Vinoth



--


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] The trouble with dihedral restraints: frozen peptide backbones

[ms]

On 18/10/10 03:30, Mark Abraham wrote:

To enforce chirality in such a toy, I thought that a simple
(naive?) but functional idea could be that of enforcing a not-
too-hard and wide-bottomed dihedral restrain on the phi angle,
like that:

[ dihedral_restraints ]
; ai   ajak
al  type  label  phi  dphi  kfac  power
; phi C'(n-1) - N - CA - C'
   3
5   7   8   11  -90 752  1


Try mdrun -debug on a short run. What does it have to say about the parameter 
values for this dihedral restraint? Likewise gmxdump -s on the .tpr (though you 
will have to search around for F_DIHRES among other things, probably).



I have run a 500 ps run with -debug on, but I am unsure of where should 
I look for what it says of relevance -I see the ctab,rtab,dtab xvg's but 
I doubt they are of interest.


m.
--
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] Reg: Running Energy minimisation for dichloroethane (DCE)

[vinothkumar mohanakrishnan]

Dear Mark

I generated the topology using pdb2gmx using dce.pdb. i checked the .top
file and i contains necessary .itp files. the error message that i get is
given below. any help is highly appreciated.

Program grompp, VERSION 4.0.7
Source code file: toppush.c, line: 620

Fatal error:
Unknown bond_atomtype CLAA

Regards
Vinoth


On Mon, Oct 18, 2010 at 5:12 PM, Mark Abraham wrote:

>  On 18/10/2010 10:31 PM, vinothkumar mohanakrishnan wrote:
>
> Hi all
>
> i added new atomtype for dichloroethane (DCE) and added the corresponding
> parameters in the .rtp, .atp, .bon.itp, .nb.itp respectively. given below
> are my additions to the corresponding files respectively.
>
> *ffoplsaa.rtp*
>
> [ DCE ]
>  [ atoms ]
>  CLAA opls_967 -0.2270  1
>  CAB  opls_966  0.2270  1
>  CAC  opls_966  0.2270  2
>  CLAD opls_967 -0.2270  2
>
> [ bonds ]
>  CLAACAB
>  CABCAC
>  CACCLAD
>
> [ angles ]
>  CLAACABCAC
>  CABCACCLAD
>
> [ dihedrals ]
>  CLAACABCACCLAD
>
> *ffoplsaa.atp*
>
>  opls_966   14.02700  ; CH2 for DCE
>  opls_967   35.45300  ; CL for DCE
>
> *ffoplsaabon.itp*
>
> [bondtypes]
> CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
> CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
> CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE
>
> [angletypes]
> CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
> CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE
>
> [dihedraltypes]
> CLAA   CAB CAC   CLAD 3 20.76096  -0.4184   27.011904
> 0.0   0.0   0.0 ; for DCE
>
> *ffoplsaanb.itp*
>
> opls_966   CAB614.027000.227   A3.98000e-01
> 4.76976e-01 ; CH2 of DCE
> opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01
> 0.20920e+01 ; Cl of DCE
>
> *dce.pdb*
>
> HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00
> 20.00CL
> HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00
> 20.00 C
> HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00
> 20.00 C
> HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00
> 20.00CL
>
> After adding all this, when i run grompp i get the error as *fatal error
> Unknown bond_atomtype CLAA*. can any one tell me why this happens?.
>
>
> I can't tell until I know how you've invoked grompp, and how you generated
> the .top file it needs, and which line of that .top grompp doesn't like.
> When asking for help, please copy and paste all your relevant commands, and
> any error messages in full, and any relevant chunk of the .top. Check
> whether that .top file is #including the correct .itp files.
>
> 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] Reg: Running Energy minimisation for dichloroethane (DCE)

[Mark Abraham]

 On 18/10/2010 10:31 PM, vinothkumar mohanakrishnan wrote:

Hi all

i added new atomtype for dichloroethane (DCE) and added the 
corresponding parameters in the .rtp, .atp, .bon.itp, .nb.itp 
respectively. given below are my additions to the corresponding files 
respectively.


*ffoplsaa.rtp*

[ DCE ]
 [ atoms ]
 CLAA opls_967 -0.2270  1
 CAB  opls_966  0.2270  1
 CAC  opls_966  0.2270  2
 CLAD opls_967 -0.2270  2

[ bonds ]
 CLAACAB
 CABCAC
 CACCLAD

[ angles ]
 CLAACABCAC
 CABCACCLAD

[ dihedrals ]
 CLAACABCACCLAD

*ffoplsaa.atp*

 opls_966   14.02700  ; CH2 for DCE
 opls_967   35.45300  ; CL for DCE

*ffoplsaabon.itp*

[bondtypes]
CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE

[angletypes]
CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE

[dihedraltypes]
CLAA   CAB CAC   CLAD 3 20.76096  -0.4184   27.011904  
0.0   0.0   0.0 ; for DCE


*ffoplsaanb.itp*

opls_966   CAB614.027000.227   A3.98000e-01  
4.76976e-01 ; CH2 of DCE
opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01  
0.20920e+01 ; Cl of DCE


*dce.pdb*

HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00 
20.00CL
HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00 
20.00 C
HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00 
20.00 C
HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00 
20.00CL


After adding all this, when i run grompp i get the error as *fatal 
error Unknown bond_atomtype CLAA*. can any one tell me why this happens?.


I can't tell until I know how you've invoked grompp, and how you 
generated the .top file it needs, and which line of that .top grompp 
doesn't like. When asking for help, please copy and paste all your 
relevant commands, and any error messages in full, and any relevant 
chunk of the .top. Check whether that .top file is #including the 
correct .itp files.


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] Reg: Running Energy minimisation for dichloroethane (DCE)

[vinothkumar mohanakrishnan]

 Hi all

i added new atomtype for dichloroethane (DCE) and added the corresponding
parameters in the .rtp, .atp, .bon.itp, .nb.itp respectively. given below
are my additions to the corresponding files respectively.

*ffoplsaa.rtp*

[ DCE ]
 [ atoms ]
 CLAA opls_967 -0.2270  1
 CAB  opls_966  0.2270  1
 CAC  opls_966  0.2270  2
 CLAD opls_967 -0.2270  2

[ bonds ]
 CLAACAB
 CABCAC
 CACCLAD

[ angles ]
 CLAACABCAC
 CABCACCLAD

[ dihedrals ]
 CLAACABCACCLAD

*ffoplsaa.atp*

 opls_966   14.02700  ; CH2 for DCE
 opls_967   35.45300  ; CL for DCE

*ffoplsaabon.itp*

[bondtypes]
CLAA  CAB 10.17870   194137.6   ; CL-CH2 for DCE
CAB   CAC 10.15300   259408.0   ; CH2-CH2 for DCE
CAC   CLAD10.17870   194137.6   ; CL-CH2 for DCE

[angletypes]
CLAA  CAB   CAC   1   108.200368.192   ; C-C-CL for DCE
CAB   CAC   CLAD  1   108.200368.192   ; C-C-CL for DCE

[dihedraltypes]
CLAA   CAB CAC   CLAD 3 20.76096  -0.4184   27.011904  0.0
0.0   0.0 ; for DCE

*ffoplsaanb.itp*

opls_966   CAB614.027000.227   A3.98000e-01  4.76976e-01
; CH2 of DCE
opls_967   CLAA   17   35.45300   -0.227   A3.16000e-01  0.20920e+01
; Cl of DCE

*dce.pdb*

HETATM1 CLAA DCE 1   6.300  -2.280   1.360  1.00
20.00CL
HETATM2  CAB DCE 1   5.060  -3.540   1.500  1.00
20.00 C
HETATM3  CAC DCE 1   4.320  -3.740   0.170  1.00
20.00 C
HETATM4 CLAD DCE 1   3.270  -2.350  -0.210  1.00
20.00CL

After adding all this, when i run grompp i get the error as *fatal error
Unknown bond_atomtype CLAA*. can any one tell me why this happens?.

Regards
Vinoth
-- 
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] simulation of a protein including calcium ion

[Justin A. Lemkul]

leila karami wrote:

Dear Justin

Is there any way for this problem?



Please keep Tsjerk's comments 
(http://lists.gromacs.org/pipermail/gmx-users/2010-October/054938.html) in mind. 
 If the calcium coordination complex plays only a minor structural role, you 
can probably justify using distance restraints (added to the topology after 
pdb2gmx).  If it has some other functionally-significant role, then you have 
harder decisions to make.


In either case, you should turn to the literature to understand how others 
typically address metal coordination complexes and make an educated choice that 
suits your needs best.  Posting to the mailing list and getting a few tips 
should not take place of doing your homework thoroughly.


-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] residence time of water molecule and life time of hydrogen

[Mark Abraham]

 On 18/10/2010 10:15 PM, atila petrosian wrote:

*Dear Justin*
**
*I read manual and specially g_hbond. but manual doesn't gime me 
information about *



  residence time of water molecule and life time of hydrogen.

In your reading, did you see that g_hbond has a -life option? Did you 
try this?


Regarding "residence time of a water molecule", you should make yourself 
aware of how one might calculate the quantity for your simulation, and 
then see if a GROMACS tool does something useful to that end. If you 
need help, please be very specific about what you're trying to calculate 
- people who might give free help are not interested in spending their 
valuable time guessing what it is you actually want. :-)


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] simulation of a protein including calcium ion

[leila karami]

Dear Justin

Is there any way for this problem?
-- 
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] residence time of water molecule and life time of hydrogen

[Justin A. Lemkul]

atila petrosian wrote:

*Dear Justin*
** 
*I read manual and specially g_hbond. but manual doesn't gime me 
information about *



  residence time of water molecule and life time of hydrogen.



The more important piece of advice I gave you was to search the mailing list 
archive.  This same question has been asked dozens of times before, and there 
are plenty of helpful posts on this topic.


-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

[gmx-users] residence time of water molecule and life time of hydrogen

[atila petrosian]

*Dear Justin*
**
*I read manual and specially g_hbond. but manual doesn't gime me information
about * residence time of water molecule and life time of hydrogen.
-- 
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] simulation of a protein including calcium ion

[Justin A. Lemkul]

leila karami wrote:

Dear Justin

thanks for your attention

in my case calcium ion is bonded to 4 oxygen atoms of C=O group of 
protein. yes, Calcium ions are present in nearly all the force fields in 
Gromacs. but there are only parameters in * nb.itp and not in * bon.itp. 
how to obtain these new parameters for my case?





Deriving new bonded parameters is not a simple task, and may or may not be 
appropriate.  A coordination complex may be held in place by electrostatic 
interactions, or you may choose to use distance restraints to keep the bonded 
geometry in place.  If you think that actual (chemical) bonds should be present, 
then there is not a suitable force field for you to use, since the charge 
transfer effects will invalidate the charges assigned to C=O groups by the force 
field.


-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

[gmx-users] simulation of a protein including calcium ion

[leila karami]

Dear Justin

thanks for your attention

in my case calcium ion is bonded to 4 oxygen atoms of C=O group of protein.
yes, Calcium ions are present in nearly all the force fields in Gromacs. but
there are only parameters in * nb.itp and not in * bon.itp. how to obtain
these new parameters for my case?
-- 
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] residence time of water molecule and life time of hydrogen

[Justin A. Lemkul]

atila petrosian wrote:

Hi gromacs users

I am beginner in gromacs. I did md simulation of a protein by gromacs 
and now

I want to obtain residence time of water molecule and life time of hydrogen
bonds. Can I obtain both of them using gromacs?
please guide me by detail.



Start by searching the gmx-users list archive.  This question gets asked 
frequently, so there are plenty of helpful posts.  Also, read the manual about 
g_hbond.


-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

[gmx-users] residence time of water molecule and life time of hydrogen

[atila petrosian]

Hi gromacs users

I am beginner in gromacs. I did md simulation of a protein by gromacs and
now
I want to obtain residence time of water molecule and life time of hydrogen
bonds. Can I obtain both of them using gromacs?
please guide me by detail.
-- 
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] Simulation parameter problem about protein unfolding

[David van der Spoel]

On 2010-10-18 06.56, chris.ne...@utoronto.ca wrote:

Generally, forcefields are not parameterized for temperatures other than
298K, so simulations are not expected to reproduce the expected
properties (like boiling water and the correct temperature denaturation
of proteins).

There's almost certainly other issues here (including the fact that I'm
entirely sure that you can get a lot more than 24 ns of simulation on a
54 aa protein; and 26 atom of pressure seems pretty arbitrary) but it
will come down to this eventually.

Just because you found a paper in which they get a denatured state does
not imply that they got the correct denatured state.

There is no correct denatured state. There are infinitely many. Check 
out recent work on NMR of "unfolded" proteins.




Chris.

-- original message --

Hi All,

I met a problem when I try to unfold a protein using Gromacs, It seemed
the protein cannot be totally unfolded!

The simulated system has one Engrailed Homeodomain (En) protein (a three
helix bundle protein with 54 residues, 629 atoms), total 4848 spce
waters, and 7 Cl- used to neutralize the system in a 5.752(nm)^3 water
BOX. I use the NTP ensemble with T=498K and P=26atm. The system has 1nm
thick water in each side of the En protein, and the density of the
system has been adjusted to0.829 g/cm3.

The simulation lasted 24ns. The helixes disappeared at about 4ns. And
after that some beta sheet formed in the N terminal of the protein.
However, the protein kept in a compact state during the 24ns simulation.
The radius of gyration of the protein fluctuated around 1.1nm during the
simulation.

I've also noticed similar simulation done by others. For example, David
Becka and Valerie Daggett reviewed their En protein unfolding in paper
"Methods for molecular dynamics simulations of protein folding/unfolding
in solution. Methods 34 (2004) 112¨C120". The simulations were performed
with ENCAD and ilmm, and used an 0.8nm cutoff range. And the ensemble
they used is NVE as I know. A stretched unfolding state occurred at
about 5ns in their 60ns simulation in 498K, with little helix structure.

I was wondering whether the difference is caused by using different MD
software and force field, or by some wrong parameters in my .mdp file.
I've also conducted another 18ns simulation, and the result is almost
the same. I list he mdp file below. Any comment is appreciated!





--
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
sp...@xray.bmc.uu.sehttp://folding.bmc.uu.se
--
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] Re: gmx-users Digest, Vol 78, Issue 123

[Gerrit Groenhof]

3. Simulation parameter problem about protein unfolding (Chen)



Hi,

You write it yourself: In paper you mention, they have used a 0.8 nm 
cutoff range for both electrostatics and cutoff. You are doing something 
different by using PME for the electrostatics. Also you are using a much 
longer cutoff for the VDW interactions.


If you want to reproduce their results, you need to stick to to the 
parameters mentioned in that paper and use cut-offs as well.


BTW I do not believe this to be a good idea though.

Gerrit



Hi All,

I met a problem when I try to unfold a protein using Gromacs, It seemed the 
protein cannot be totally unfolded!

The simulated system has one Engrailed Homeodomain (En) protein (a three helix 
bundle protein with 54 residues, 629 atoms), total 4848 spce waters, and 7 Cl- 
used to neutralize the system in a 5.752(nm)3  water BOX. I use the NTP 
ensemble with T=498K and P=26atm. The system has 1nm thick water in each side 
of the En protein, and the density of the system has been adjusted to0.829 
g/cm3.

The simulation lasted 24ns. The helixes disappeared at about 4ns. And after 
that some beta sheet formed in the N terminal of the protein. However, the 
protein kept in a compact state during the 24ns simulation. The radius of 
gyration of the protein fluctuated around 1.1nm during the simulation.

I've also noticed similar simulation done by others. For example, David Becka and Valerie 
Daggett reviewed their En protein unfolding in paper "Methods for molecular dynamics 
simulations of protein folding/unfolding in solution. Methods 34 (2004) 112¨C120". 
The simulations were performed with ENCAD and ilmm, and used an 0.8nm cutoff range. And 
the ensemble they used is NVE as I know. A stretched unfolding state occurred at about 
5ns in their 60ns simulation in 498K, with little helix structure.

I was wondering whether the difference is caused by using different MD software 
and force field, or by some wrong parameters in my .mdp file. I've also 
conducted another 18ns simulation, and the result is almost the same. I list he 
mdp file below. Any comment is appreciated!

-- 
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] Simulation parameter problem about protein unfolding

[Chen]

>At 2010-10-18 12:56:31，chris.ne...@utoronto.ca wrote:
>Generally, forcefields are not parameterized for temperatures other  
>than 298K, so simulations are not expected to reproduce the expected  
>properties (like boiling water and the correct temperature  
>denaturation of proteins).
>
>There's almost certainly other issues here (including the fact that  
>I'm entirely sure that you can get a lot more than 24 ns of simulation  
>on a 54 aa protein; and 26 atom of pressure seems pretty arbitrary)  
>but it will come down to this eventually.
>
>Just because you found a paper in which they get a denatured state  
>does not imply that they got the correct denatured state.
>
>Chris.
>

Hi, Chris! Thanks for the reply! I have not conducted unfolding
before, so I compared my result with other's to make sure it is reasonable in
some extend. The 26 atm pressure comes from experiment result (Haar et al.,
1984) mentioned in some MD related papers (e.g. 'J. Mol. Biol. (2000) 296,
1257-1282').

 I've searched the maillist
and noticed some issues about High temperature simulations. I'm not sure
whether the 'ilmm' force field has been optimized for high temperature
simulation. I also noticed some users asked about MD parameters in 'unfolding a
protein'. And the parameters they used are different from ours' (e.g. the
'rlist', 'rcoulomb' and 'temperature or pressure couple algorithm').

 I just want to make sure I
didn't make mistakes in these parameters which maybe cause the protein keeping
in a compact state! The radius of gyration of the protein fluctuated around
1.1nm (never bigger than 1.4nm) during our 30ns simulation now. If the MD
parameters I chose have no problem, then I could keep going. Any comment?




>Hi,
>If I can add to this discussion, I think that your results are very 
>reasonable. Proteins in solution are not straight lines, but fold to 
>some extent. If you wish to have en elongated protein you have to pull 
>it like it is being done in AFM experiments.
>
>Itamar

Hi, Itamar! Thanks for the reply! I didn't intend to make the
unfolding protein like a straight line. I just 'feel' the protein in 498K
should fluctuate a little larger. So I want to make sure I didn't make mistake
in MD parameters setting or water BOX preparing.


>-- original message --
>
>Hi All,
>
>I met a problem when I try to unfold a protein using Gromacs, It  
>seemed the protein cannot be totally unfolded!
>
>The simulated system has one Engrailed Homeodomain (En) protein (a  
>three helix bundle protein with 54 residues, 629 atoms), total 4848  
>spce waters, and 7 Cl- used to neutralize the system in a 5.752(nm)^3  
>water BOX. I use the NTP ensemble with T=498K and P=26atm. The system  
>has 1nm thick water in each side of the En protein, and the density of  
>the system has been adjusted to0.829 g/cm3.
>
>The simulation lasted 24ns. The helixes disappeared at about 4ns. And  
>after that some beta sheet formed in the N terminal of the protein.  
>However, the protein kept in a compact state during the 24ns  
>simulation. The radius of gyration of the protein fluctuated around  
>1.1nm during the simulation.
>
>I've also noticed similar simulation done by others. For example,  
>David Becka and Valerie Daggett reviewed their En protein unfolding in  
>paper "Methods for molecular dynamics simulations of protein  
>folding/unfolding in solution. Methods 34 (2004) 112¨C120". The  
>simulations were performed with ENCAD and ilmm, and used an 0.8nm  
>cutoff range. And the ensemble they used is NVE as I know. A stretched  
>unfolding state occurred at about 5ns in their 60ns simulation in  
>498K, with little helix structure.
>
>I was wondering whether the difference is caused by using different MD  
>software and force field, or by some wrong parameters in my .mdp file.  
>I've also conducted another 18ns simulation, and the result is almost  
>the same. I list he mdp file below. Any comment is appreciated!
>
-- 
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 does g_msd calculates MSD?

[Javier Cerezo]

Hello Julian.

I've not gone though the C-code neither. However, from what I have read 
and understood about MSD calculation, both calculations you proposed 
should not give the same results. This is the way it is calculated.


1. You take different starting times over tyour trajectory, every 
-trestart picoseconds, from -b (let's take 0ps in this explanation) to -e.


2. From every restarting point (tr) you calculate the MSD, by 
calculating the diplacement of atomic positions compared to the position 
at the corresponding tr.


3. With the previous result you'll have a lot of MSD curves, 
corresponding to each restating point, e.g. one where t0=0ps, other 
t0=0+trestart, and so on. Actually, you are refering all curves to its 
corresponding t0, so you have all of them starting at 0.


4. The average of all of these curves is done. With that you obtain 
properly averaged results (which is essential for a good MSD 
calculation). However, note that only the restarting point at 0ps of 
your original trayectory (or the -b time you selected) will contain all 
the points of the curve since, for example, the last restating point wil 
conntribute with the few points remaining to the end of the simulation 
(-e time). For that reason, the average is made with less points as the 
time from the reference is increased, and at the end of the curve, the 
results are not well averaged any more.


Hope that helps.

Javier

El 15/10/10 23:58, Julian escribió:

Hi everybody,

We're trying to understand how gromacs g_msd calculates MSD (without 
reading through the C code, I really don't know much about it, and 
there are hundreds of lines).
(What we want to do is to calculate MSD properly in our supercooled 
water simulations, i.e., to choose correctly the values for -b, -e  
and -trestart in order to get the longest, reliable MSD data)



We found a mail from Gaurav Goel (quoted below) which gives a 
reasonable explanation on the topic.


If we understand this correctly, then the output of

g_msd -b 50 -e 100 -trestart 50

should be the same as the second half (with the proper shifting) of
the output of

g_msd -b 0 -e 100 -trestart 50

But it is not, as anyone can verify with any simulation.
.. what are we missing?

Thanks in advance,
Julian

--
Julian Gelman Constantin
Department of Inorganic, Analytic and Chemical Physics (DQIAQF)
School of Exact and Natural Sciences
University of Buenos Aires, Argentina


---
Gaurav Goel gauravgoeluta at gmail.com
Tue Jul 13 00:56:02 CEST 2010

On Mon, Jul 12, 2010 at 5:22 PM, Ricardo Cuya Guizado
 wrote:
> Dear gromacs users
> I make a MD of 20 ns of a solute in water
> With the g_msd program the msd vs the time was obtained
> In the plot, I observed a linear behaviour of the MSD from 0 to 15 
ns and a

> plateau with no linear tendence at the last 5 ns arpoximately.
> In order to know if the observed plateau was due to the data or is 
due to

> the way as the algorithm process the data, I divided the MD in two
> trajectories and obtained the msd for each one.
> From 0-10ns, the plot observed shows a linear tendence en the begining
and a
> plateau with no linear tendence from 9 to 10 ns.
> From 10-20 ns the plot observed was linear from 10 to 18 ns and not 
linear

> at the last, the same plateau was observed.
> Comparing the plots there are not equivalent,.
> Why g_msd produces a non linear plot at the last of the calculation 
and the

> plateau is ever produces.
> Somebody will explain the way as the g_msd algorithm work? and why 
the plot

> are no equivalent or why there must be equivalent?

I will explain how the g_msd algorithm works and hopefully that will
answer all your questions above. What you see in the output file is
average-MSD versus time. This average is done over all the particles
in the group you selected and over multiple time origins (this last
option can be selected with the -trestart parameter). Also, time in
column 1 is time difference from the start of your trajectory to
current time.

E.g., let's say you collected a trajectory over 5 time units and
choose -trestart=1 time unit and -dt=1 time unit.

dt=1 means you'll have 6 configurations for your analysis (including
the configuration at t=0).

trestart=1 means you'll have 5 distinct trajectories for your analysis:
Trajectory 1: 0-5
T2: 1-5
T3: 2-5
T4: 3-5
T5: 4-5

Now you can notice that all 5 trajectories contribute to the average
MSD after 1 time unit (T1-T5), 4 trajectories contribute to the
average MSD after 2 time units (T1-T4),  3 trajectories to the average
MSD after 3 time units (T1-T3), , and only one trajectory to the
MSD after 5 time units (T1). Of course, this assumes that trestart is
large enough that all all these trajectories are uncorrelated.

So, it's clear that longer the time interval at which you want to
evaluate the MSD lesser the number of trajectories used to evaluate
it...and hence, higher error in MSD values at longer times. That might
explain deviation from linear behaviour at