[gmx-users] dihedral MM scan

[Mohsen Ramezanpour]

Thanks Justin,

I did the scan for the whole range from -180 to 180 in 10 intervals.
I used emtol = 50, and restrain force of 100000 to have a better relaxation
and restraining, respectively.
Unfortunately, my MM profiles do not look similar to GAAMP QM or MM
profiles, although I have used the suggested parameters.

To check the method I am using and the parameters in the .mdp files, I did
the MM scanning on a rotatable bond (CH2-CH2-CH2-CH3) similar to one can be
found in a lipid tail.
This dihedral was detected as a soft dihedral by GAAMP and there was a good
agreement between QM and MM profiles.
Here again, my MM profile from scanning does not seem reasonable.
Most importantly, there is a very large barrier at angle 0 (ca. 300 times
of energies in other angles). Besides, there are sudden increase or
decrease in energy just by 10 degree of rotation.

Any suggestion is appreciated in advance.
Cheers
Mohsen




On Wed, Jan 18, 2017 at 7:55 PM, Justin Lemkul <jalemkul at vt.edu> wrote:

>
>
> On 1/18/17 9:51 PM, Mohsen Ramezanpour wrote:
>
>> Dear Gromacs users,
>>
>> I want to do a MM scan on a dihedral and compare with the available QM
>> profile using Gromacs 5.1.3
>>
>> I used the following (updated compared to my last email on parameters)
>> em.mdp and md.mdp files for relaxation and zero-step md, respectively.
>>
>> *em.mdp:*
>> define                     = -DPOSRESDIH
>> integrator                 = steep
>> emtol                       = 100.0
>> nsteps                     = 50000
>> vdwtype                   = Cut-off
>> vdw-modifier            = none
>> coulombtype             = Cut-off
>> constraints                = h-bonds
>> constraint_algorithm    = LINCS
>> nstlist                        = 0
>> ns-type                      = simple
>> pbc                           = no
>> cutoff-scheme           = group
>> rcoulomb                   = 0
>> rvdw                         = 0
>> rvdw_switch             = 0
>> rlist                            = 0
>>
>> Which, by -DPOSRES I am applying position restraint on specific dihedral
>> of
>> interest.
>> The following lines were included at the end of molecule.itp
>>
>> #ifdef POSRES
>>
>> [ dihedral_restraints ]
>>
>>  1    18  8   9     1   120   0   10000
>>
>> #endif
>>
>>
>>
>> *as for md.mdp:*
>>
>> define                  =
>> integrator              = md
>> dt                        = 0.001
>> nsteps                  = 0
>> nstlog                  = 0
>> nstxout                 = 0
>> nstvout                 = 0
>> nstfout                  = 0
>> nstcalcenergy           = 0
>> nstenergy                 = 0
>> nstxout-compressed      = 100000
>> compressed-x-precision  = 1000
>> coulombtype                  = Cut-off
>> vdwtype                         = Cut-off
>> vdw-modifier                   = none
>> constraints                     = h-bonds
>> constraint_algorithm        = LINCS
>> nstcomm                         = 100
>> comm_mode               = linear
>> comm_grps                 = system
>> refcoord_scaling           = com
>> nstlist                          = 0
>> ns-type                        = simple
>> pbc                            = no
>> cutoff-scheme              = group
>> rcoulomb                = 0
>> rvdw                      = 0
>> rvdw_switch            = 0
>> rlist                         = 0
>> continuation            = yes
>>
>> It seems to work fine, at least because grompp do not complain and it
>> gives
>> me the desired files.
>>
>> There is one problem:
>> Looking at the Potential Energy which zero-step md for angle = 120 gives
>> me
>> a value of 384.682 KJ/mol (ca. 91 Kcal/mol) which is almost 20 times the
>> value from QM profile (4.5 Kcal/mol) by GAAMP server (as I understood, the
>> units for dihedral profiles by GAAMP are in Kcal/mol).
>>
>> I used the parameters GAAMP gave me (partial charges and dihedral
>> parameters).
>>
>> Briefly, the same molecule which I uploaded to GAAMP was used for MM scan
>> in vaccume. No box and infinite cutoff was used. I did an EM step by above
>> em.mdp and let the system to relax. Then, I removed the restraints and did
>> the zero-step md part and got this energy.
>>
>> It is worth mentioning that the molecule is positively charges. BUT, as I
>> understood, this is what GAAMP also uses for the QM scan.
>>
>> What could be the problem?
>>
>> I highly appreciate your comments as expert in this.
>>
>>
> All potential energy scans need to be offset to the global minimum of the
> scan. You can't take the energy of one configuration and compare its value
> to the QM. You have to do the whole scan and offset it.  The relative
> energy is what matters.
>
> -Justin
>
> --
> ==================================================
>
> Justin A. Lemkul, Ph.D.
> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>
> Department of Pharmaceutical Sciences
> School of Pharmacy
> Health Sciences Facility II, Room 629
> University of Maryland, Baltimore
> 20 Penn St.
> Baltimore, MD 21201
>
> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
> http://mackerell.umaryland.edu/~jalemkul
>
> ==================================================
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/Support
> /Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>



-- 
*Rewards work better than punishment ...*