Re: [AMBER] Calculation of Urea_Urea Interaction energy

From: MOHD HOMAIDUR RAHMAN <rahmanhpu.gmail.com>
Date: Mon, 8 Sep 2014 13:26:30 +0530

Dear Jason, Daniel and All Amber user's

Thank you for your help and guidance

I updated AmberTool 13 to AmberTool !4 and I try to understand how code is
working.
For this I check with few combination, but the result I got is not clear to
me.
I calculated Interaction Energy by pairwise command in cpptraj 14 like

Input file

trajin 6_md1_ureatest.crd
reference 6_md1_ureatest.rst
autoimage
pairwise uu :MOL out pair_uu.txt cuteelec 12.0 cutevdw 12.0

and in other seperately

trajin ura_leap.crd
reference ura_leap.crd
pairwise uu :MOL out pair_uu.txt cuteelec 12.0 cutevdw 12.0

Case 1
I given input for trajin and reference both leap generated crd file and
with (electrostatic, vander waal cutoff 12.0) or without cutoff and I am
expecting zero kcal/mol interaction energy. But I got very high value and
the value is not effected by cutoff distance.

#Frame uu_cut[EVDW] uu_cut[EELEC] uu_ncut[EVDW]
uu_ncut[EELEC]
       1 45847167.8970 29301.6094 45847167.8970
29301.6094

The same value I am also getting when I am using the leap generated crd and
md gerated rst file.

Case 2
I given input to calculate Interaction energy between two pairs atom. For
this I taken urea oxygen and all urea hydrogen. Again in this case I am
getting high positive value. But I am expecting some negative value because
oxygen have negative charge (-0.641711) and hydrogen have positive charge
(0.338638). In this case I used reference file as leap generated crd file.

#Frame uu_cut[EVDW] uu_cut[EELEC]
       1 16926599.0449 1284114.2386

Case 3
In this case I just used residue wise pairwise Interaction Energy
calculation. For reference input, I used both leap generated crd and md
generated restart file. Both cases giving same result and the result is
given below.

#Frame uu_cut[EVDW] uu_cut[EELEC] uu_ncut[EVDW]
uu_ncut[EELEC]
       1 -2711.7978 24536.6033 -2711.7978
        24536.6033

All this result not making sense to me, So I need guidance and help from
you all expert.

Other query how pairwise deal with 1-4 electrostatic and this value how
contributing in total interaction energy. Because my mdout show total EELEC
positive and and 1-4 EEL negative and this value also not making sense to
me.

------------------------------------------------------------------------------
      A V E R A G E S O V E R 500000 S T E P S

 NSTEP = 500000 TIME(PS) = 1840.000 TEMP(K) = 297.87 PRESS =
1.0
 Etot = -87724.5924 EKtot = 4260.9239 EPtot =
 -91985.5163
 BOND = 2486.7436 ANGLE = 1560.5391 DIHED =
 6173.2503
 1-4 NB = -76.4015 1-4 EEL = -121540.5102 VDWAALS =
-3548.8691
 EELEC = 22959.7315 EHBOND = 0.0000 RESTRAINT =
0.0000
 EKCMT = 564.2516 VIRIAL = 563.3281 VOLUME =
41876.3880
                                                    Density =
1.4290
 Ewald error estimate: 0.2502E-03
 ------------------------------------------------------------------------------

All this I am doing to reproduce the result already publish by Sandip Paul
and Grenfell N. Patey in J. Am. Chem. Soc., 2007, 129 (14), pp 4476–4482.
For more information I am attaching his paper here. Please see the table: 3
in this paper where he calculated Average Interaction Energies for the urea
, tmao and water Systems. In this paper he used OPLS parameter for urea and
I also simulated with same parameter. I am able to reproduce water-water,
urea-water, tmao-water but unable to reproduce tmao-tmao and urea-urea
average interaction energy.

Sorry for very long paragraph

Thank you in advance for your guidance and help

Regards
Rahman

'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
*Md Homaidur Rahman*
PhD Research Scholar
Computational Biophysics Lab.
Department of Biotechnology
Indian Institute of Technology-Madras
Chennai-600 036
*Mobile No = +91- 7845991785*
'''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''

On Sat, Sep 6, 2014 at 8:11 AM, Daniel Roe <daniel.r.roe.gmail.com> wrote:

> Hi,
>
> Just some general advice - you are better off upgrading to the latest
> version of AmberTools/Cpptraj (14,
> http://ambermd.org/AmberTools14-get.html) before attempting any
> calculation. Version 13 is quite out of date at this point.
>
> -Dan
>
> On Fri, Sep 5, 2014 at 12:06 PM, MOHD HOMAIDUR RAHMAN
> <rahmanhpu.gmail.com> wrote:
> > Dear Jason
> >
> > Thank you for help
> >
> >
> >> Dear Amber Users
> >> I want to calculate interaction energy between two urea molecules in
> >> homogeneous urea solution (only urea molecules are there). For that I
> used
> >> antechamber modules to generate urea "prep" and "frc" file. The
> generated
> >> file contain GAFF force field parameter. After production MD run, I
> > process
> >> the trajectory to calculate pairwise interaction energy. I used the
> >> pairwise command but I am getting + 36.34 kcal/mol while the reported
> >> theoretical value for dimer is -10.33 kcal/mol (MP2/aug-cc-pVDZ).
> >>
> >> I am using "CPPTRAJ: Trajectory Analysis. V13.0".
> >>
> >
> > What's the actual command you're using? (Not just the command name, the
> > _full_ command).
> >
> > I am using the following command
> >
> > cpptraj -p urea.top -i pairwise.in
> >
> > In pairwise.in
> >
> > trajin ../md_file/
> > md1
> > .crd
> >
> > reference ../md_file/urea
> > _leap
> > .crd
> >
> > pairwise U_U :URA out pairwise_
> > IE
> > .txt cuteelec 12.0 cutevdw 12.0
> >
> >
> > Is there are bugs in this command??
> >>
> >
> > cpptraj is part of Amber. Amber has no bugs. So cpptraj has no bugs.
> >
> > Kidding, of course, but in cases like this I've found that it's best to
> > assume I don't fully understand what I'm computing or exactly what I
> asked
> > the program to compute. This is why knowing the exact command is
> important.
> >
> > Also, there's no reason to expect that two urea molecules in bulk urea
> have
> > the same interaction as the urea dimer... You're also comparing a
> minimum
> > energy structure (MP2/aug-cc-pVDZ) with what is likely an ensemble
> property
> > if you're computing this value over a trajectory. When doing comparisons
> > like this, it's very important that you make every effort to make sure
> you
> > are making a fair comparison.
> >
> > The main point is that, How urea interaction energy showing positive
> value
> > ?
> > It means that the system is thermodynamically very unstable due to
> positive
> > interaction energy (unfavorable interaction). The value should be equal
> to
> > zero
> > or below zero in any range (never exactly match with quantum calculation)
> > because we are calculating over trajectory.
> >
> > Or may be the high positive Interaction value arise due to wrong
> non-bonded
> > parameter.
> >
> > Thanks
> > Rahman
> > _______________________________________________
> > AMBER mailing list
> > AMBER.ambermd.org
> > http://lists.ambermd.org/mailman/listinfo/amber
>
>
>
> --
> -------------------------
> Daniel R. Roe, PhD
> Department of Medicinal Chemistry
> University of Utah
> 30 South 2000 East, Room 307
> Salt Lake City, UT 84112-5820
> http://home.chpc.utah.edu/~cheatham/
> (801) 587-9652
> (801) 585-6208 (Fax)
>
> _______________________________________________
> AMBER mailing list
> AMBER.ambermd.org
> http://lists.ambermd.org/mailman/listinfo/amber
>


_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber

Received on Mon Sep 08 2014 - 01:00:02 PDT
Custom Search