[gmx-users] Re: How to set the sigma and epsilon for Cu2+ in OPLS-AA/L

[Dr. Vitaly Chaban]

On Mon, Apr 8, 2013 at 3:36 PM, fantasticqhl <fantasticqhl at gmail.com> wrote:

>  Dear Dr. Vitaly Chaban,
>
> Thanks very much for your patient explanation. Yeah, you are right, that
> is what I want to know: how you tuned this parameter?
>
> Since then, if I want to set a new atom type and I know its vdw radius, so
> how should I set the sigma for it based on the vdw radius,
>


You cannot set the sigma based ONLY on the VDW radius.



>  which should be in agreement with OPLS-AA/L force filed? Could you give
> me some suggestions?
> I guess that I have to tune it by myself this time, right? Thanks in
> advance!
>


I would do the following:

1) Optimize ion-ligand complex using ab initio. Write down binding energy
and optimal distance;
2) Construct topology for classical MD using approximate sigma;
3) Calculate energy and distance from classical MD;
4) Compare them to distance and energy from ab initio;
5) If you are not satisfied, adjust your sigma;
6) Repeat classical MD until the difference between ion-ligand distance in
classical MD becomes reasonably similar to that in ab initio.

To preserve compatibility with OPLS, use the same level of theory in ab
initio, which they used when derived OPLS. Keep in mind that their original
level of theory is not so perfect...


Dr. Vitaly Chaban



>
> All the best,
> Qinghua Liao
>
>
>
>
> On 04/08/2013 03:21 PM, Dr. Vitaly Chaban wrote:
>
> I think your misunderstanding comes from the belief that sigma (as they
> are tabulated in the force field files) should *exactly correspond* to the
> VDW diameter, as in encyclopedia.
>
>  This is simply not the case. In reality, sigmas in the force fields are
> "tuned" in order to give right interatomic distances AFTER you turn on all
> the necessary potentials (Coulombic attraction in case of OPLS).
>
>  Dr. Vitaly Chaban
>
>
>
>
>
> On Mon, Apr 8, 2013 at 3:14 PM, fantasticqhl <fantasticqhl at gmail.com>wrote:
>
>>  Dear Dr. Vitaly Chaban,
>>
>> Thanks very much for concern on my research! We are going to the use the
>> bonded model together with Coulomb and LJ potentials.
>>
>> My problem is that vdw radius and its sigma do not follow the equation
>> of  Rvdw = pow(2, 1/6)*sigma in the OPLS force field files,
>> not just for copper. That's why I sent these e-mails for suggestions. I
>> am sorry for the unclear.
>>
>> All the best,
>> Qinghua Liao
>>
>>
>> On 04/08/2013 01:22 PM, Dr. Vitaly Chaban wrote:
>>
>> Dear Qinghua Liao -
>>
>>  In that case, I am just wishing you luck with the copper containing
>> systems.
>>
>>  Are you going to simulate copper-ligand interactions using Coulomb+LJ
>> potential only? I would guess it is a chemical bonding case. Maybe the
>> Morse potential (additionally) can be of better service?
>>
>>
>>  Dr. Vitaly Chaban
>>
>>
>>
>>
>> On Mon, Apr 8, 2013 at 1:09 PM, fantasticqhl <fantasticqhl at gmail.com>wrote:
>>
>>>  Dear Dr. Vitaly Chaban,
>>>
>>> Thanks very much for your explanation. I guess that I get what you mean
>>> now! Thanks!
>>>
>>> All the best,
>>> Qinghua Liao
>>>
>>> On 04/07/2013 11:35 AM, Dr. Vitaly Chaban wrote:
>>>
>>> The equation is a direct consequence of LJ-12-6 equation. This equation
>>> is used in OPLS and most other force fields.
>>>
>>>  The difference you found originate from the fact that, besides LJ
>>> potential, there is much stronger Coulomb potential in the copper-ion case.
>>> If you run simulations, you will see that copper-ligand distance is smaller
>>> than the sum of their sigmas multiplied by pow (2, 1/6).
>>>
>>>
>>>  Dr. Vitaly Chaban
>>>
>>>
>>>
>>>
>>>
>>>
>>>  On Sun, Apr 7, 2013 at 11:28 AM, fantasticqhl <fantasticqhl at gmail.com>wrote:
>>>
>>>>  Dear Dr. Vitaly Chaban,
>>>>
>>>> Thanks for the explanation. I know this equation. However, the van der
>>>> Waals radius and its counterpart sigma in OPLS-AA/L force field files do
>>>> not follow this equation.
>>>>
>>>> For example, the vdw radius of copper ion is 1.4 angstrom, and its
>>>> sigma is  2.08470e-01 (I guess the unit is nm). pow(2, 1/6) is more than 1,
>>>> so obviously this equation
>>>> does not work with copper. So do other atoms. I guess that there might
>>>> be an additional coefficient for this equation in gromacs. That's the
>>>> purpose for asking. Thanks very much!
>>>>
>>>>
>>>> All the best,
>>>> Qinghua
>>>>
>>>> On 04/07/2013 10:48 AM, Dr. Vitaly Chaban wrote:
>>>>
>>>> Dear Qinghua -
>>>>
>>>>  The formal relation is diameter = pow (2, 1/6) * sigma, provided that
>>>> you have only LJ potential in your interacting subsystem.
>>>>
>>>>  If this is not the case, an optimal sigma can only be found
>>>> iteratively.
>>>>
>>>>
>>>>  Dr. Vitaly Chaban
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> On Sun, Apr 7, 2013 at 10:36 AM, fantasticqhl <fantasticqhl at gmail.com>wrote:
>>>>
>>>>> Dear Dr. Vitaly Chaban,
>>>>>
>>>>> Thanks very much for your reply. My question is the relationship
>>>>> between van der Waals radius and sigma in the OPLS-AA/L force filed files
>>>>> of Gromacs.
>>>>>
>>>>> Of course I did ab initio optimizations of my system, but I do not
>>>>> know there is some relation between the optimal bond length (copper--atom
>>>>> of the ligand) and sigma.
>>>>> Could you please be more clear and give a little detailed explanation?
>>>>> Thanks very much!
>>>>>
>>>>> All the best,
>>>>> Qinghua
>>>>>
>>>>> On 04/06/2013 06:07 PM, Dr. Vitaly Chaban wrote:
>>>>>
>>>>>> In systems of such kind, everything will depend on the atom of the
>>>>>> ligand,
>>>>>> which coordinated by copper ion.
>>>>>>
>>>>>> Perform ab initio geometry optimization and find the optimal
>>>>>> distance. Then
>>>>>> adjust sigma(s).
>>>>>>
>>>>>> Dr. Vitaly Chaban
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> There is a copper ion with four ligands in my system. I am going to
>>>>>>
>>>>>>> study this system using MD simulations.
>>>>>>> For the vdW parameters, R*=1.74 angstrom and epsilon=1.14 kcal.mol
>>>>>>> from
>>>>>>> one paper will be used in our
>>>>>>> simulations. I already found the parameters of copper ion (Cu2+) in
>>>>>>> the
>>>>>>> OPLS-AA/L force field files:
>>>>>>> sigma= 2.08470e-01, epsilon=4.76976e+00, which are for Cu2+ without
>>>>>>> ligands. The two epsilon are the same,
>>>>>>> just with different units.
>>>>>>>
>>>>>>> My question is that I do not know how to convert the vdW radius to
>>>>>>> sigma. I found that the vdw radius of copper is
>>>>>>> 1.4 angstrom, and the sigma in the force field file is 2.08470e-01.
>>>>>>> Could someone tell me how to do the converting?
>>>>>>>
>>>>>>> Thanks very much!
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>
>>>>
>>>>
>>>
>>>
>>
>>
>
>