AMBER Archive (2008)

Subject: Re: RE: AMBER: from dihe parameters to torsional energy

From: Jeffrey (jeffry20072008_at_yahoo.cn)
Date: Wed Dec 17 2008 - 20:18:29 CST


Dear Ross,

      Thanks very much for your detailed explanation. I get a better understanding of this question now. I have two more questions for this post.

1. I'd like to be sure how to calculate the torsional energy under the condition when multiple lines of parameters for just one dihedral angle exist.

Once more, taking CT-CT-CT-CT in parm99.dat as an example,
CT-CT-CT-CT 1 0.18 0.0 -3. Junmei et al, 1999
CT-CT-CT-CT 1 0.25 180.0 -2. Junmei et al, 1999
CT-CT-CT-CT 1 0.20 180.0 1. Junmei et al, 1999
 
should the total torsional energy for this dihedral be
 E(torsion for CT-CT-CT-CT )
 =
 ( pk1*(1 + Cos[pn1*phi - phase1]) + pk2*(1 + Cos[pn2*phi - phase2]) + pk3*(1 + Cos[pn3*phi - phase3]) )
 
where pk1 = 0.18; pk2 = 0.25; pk3 = 0.20; phase1 = 0; phase2 = Pi; phase3 = Pi; pn1 = 1; pn2 = 2; pn3 = 3;

2. Is there a way to get the energy values in AMBER for each of the terms produced by the following commands
$AMBERHOME/exe/rdparm prmtop
and type 'dihedrals'
?

Thanks in advance.
Have a nice day.

-----
Jeffrey

Dear Jeffrey,
 
My understanding of this is that the two equations are the same thing given certain conditions with respect to phase.
 
The 'real' equation as written in terms of the AMBER force field is:
 
e = pk(1.0+cos(pn*phi-phase))
 
If my trigonometry is correct (it is early in the morning so buyer beware), The cos term can be expanded here as:
 
e = pk(1.0+cos(pn*phi)*cos(-phase))
 
You would not want to do this computationally since cos functions are expensive, however, in the case of phase = 0 then cos(-0) = 1 and if phase = 180 then cos(-180) = -1
 
Hence the above equation becomes:
 
e = pk(1.0(+/-)cos(pn*phi))
 
The net result though is that it is the same function as the original definition above. Hence it should yield the same results. I don't know for sure but I suspect that it is written this way either to make the derivatives easier to express or compute or that the standard way of writing it has some kind of instability / accuracy issue in certain circumstances, i.e. when the dihedral gets close to either 0, 90 or 180.
 
Hence the short answer is that it shouldn't matter which route through the code gets taken here - in both cases you should end up with the same answer.
 
If you want a list of dihedrals from the prmtop then you can run:
 
$AMBERHOME/exe/rdparm prmtop
 
and type 'dihedrals'
 
This will list everything from which given a set of coordinates and the original equation at the top above you should end up with the same energy as is calculated using the expanded terms above.
 
I hope that helps,
All the best
Ross
 
From: owner-amber_at_scripps.edu [mailto:owner-amber_at_scripps.edu] On Behalf Of Jeffrey
Sent: Wednesday, December 17, 2008 1:15 AM
To: amber_at_scripps.edu
Subject: AMBER: from dihe parameters to torsional energy
 
Dear all,
 
         I would like to figure out how dihedral parameters are used to calculate the torsional energy. After checking the subroutine subroutine ephi (src/nmode/ephi.f), I found two expressions were adopted for the dihedral energy:
1>. e = pk(ic) * (1.0+phase*cos(pn(ic)*phi)
  where phase = 1.0 or -1.0, and pn = 1,2,3,4, or 6
2>. e = pk*( 1.0+cos(pn*phi-phase) ----- the old energy form,
if phase angle is other than 0 or pi then assume this angle.
 
Taking CT-CT-CT-CT in parm99.dat as an example,
CT-CT-CT-CT 1 0.18 0.0 -3. Junmei et al, 1999
CT-CT-CT-CT 1 0.25 180.0 -2. Junmei et al, 1999
CT-CT-CT-CT 1 0.20 180.0 1. Junmei et al, 1999
 
should the total torsional energy for this dihedral be one of the two forms:
a. (pk1*(1 + phase*Cos[pn1*phi ]) + pk2*(1 + phase*Cos[pn2*phi ]) + pk3*(1 + phase*Cos[pn3*phi]))
b. (pk1*(1 + Cos[pn1*phi - phase1]) + pk2*(1 + Cos[pn2*phi - phase2]) + pk3*(1 + Cos[pn3*phi - phase3]))
 
where pk1 = 0.18; pk2 = 0.25; pk3 = 0.20; phase1 = 0; phase2 = Pi; phase3 = Pi; pn1 = 1; pn2 = 2; pn3 = 3;
phase=1 or -1
 
 
Another question is how to output each energy term in AMBER to check if I have assigned a correct form of the expression?
 
Any suggestion is greatly appreciated.
Thanks very much for your time.
 
-------
Jeffrey
 
 

__________________________________________________
8O?lW"2aQE;"3,4sH]A?Cb7QSJOd?
http://cn.mail.yahoo.com

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" (in the *body* of the email)
      to majordomo_at_scripps.edu