[gmx-users] OPLS parametrizaton

[Marcelo Silva]

Em 02-04-2011 0:33, Justin A. Lemkul escreveu:
>
>
> Marcelo Silva wrote:
>> Hi,
>>
>> I am trying to simulate mixtures of ethanol and trifluoroethanol and 
>> I would like to use the same parameters used on OPLS parametrization.
>>
>> On the OPLS1996 paper the authors starts by stating that they've used 
>> 267 molecules and "The box size varied from approximately 26 26  26 
>> Å for methanol to 37  37  37 Å for cyclohexane". In my simulations 
>> I am considering using 500 molecules (each box vector would be around 
>> 37 Å).
>>
>> "In most cases, the intermolecular nonbonded interactions were 
>> truncated at 11Å based on roughly the center-of-mass separations with 
>> quadratic smoothing of the interaction energy to zero over the last 
>> 0.5 Å (...) A standard correction was made for the Lennard-Jones 
>> interactions neglected beyond the cutoff."
>>
>> Trying to interpret what is written I think the correct way to treat 
>> cut-offs and energy corrections would be like this:
>>
>> ---------------------------------------
>> Cut-offs
>>
>> rlist                    = 1.4
>> coulombtype              = Reaction-field-zero
>> rcoulomb                 = 1.1
>> rcuolomb_shift      = 0.6
>> vdwtype               = shift
>> rvdw                     = 1.1
>> rvdw_switch        = 0.6
>>
>> Dispersion Correction
>>
>> DispCorr             = EnerPress
>>
>> ---------------------------------------
>>
>> However I have some doubts because I'm using a box with a different 
>> number of molecules/ different size (should I change the cutoffs to 
>> 15Å like it was for cyclohexane?) and I'm not sure about using 
>> Reaction-field-zero instead of the shift function for coulombic 
>> interactions.
>>
>> Could you please help me confirm this parameters so I can run a 
>> simulation in the same conditions as OPLS was parameterized (although 
>> using MC and not MD).
>>
>
> I don't think there is a perfectly clear answer to this.  OPLS was 
> initially designed for MC with some special functional forms and is 
> now applied to MD using all manner of settings.  Commonly, one sees 
> rcoulomb=rvdw=rlist=1.0 or 1.4 in conjunction with PME.  There have 
> been a variety of systematic force field comparisons in recent years.  
> You should consult those to see what their specific recommendations 
> are with respect to proper settings, although many of them focus on 
> proteins.
>
> -Justin
>
>> Best regards,
>>
>> Marcelo
>

Just to test it, I used this options (nstlist = 5 with timestep 1fs):

; Cut-offs

rlist                    = 1.4
coulombtype              = Switch
rcoulomb                 = 1.1
rcuolomb_shift      = 1.05
vdwtype               = Switch
rvdw                     = 1.1
rvdw_switch        = 1.05

; Dispersion Correction

DispCorr             = EnerPres

and compared to vdw plain cutoff at 1.1 and PME with cutoff 1.1 nm 
(rlist = 1.1 nm also). The big difference in the output seems to be with 
the Coulomb (SR) term:

-1.25490e+05 kj/mol with Switch
-5.38895e+03  kj/mol   with PME

With switch I get a low density for ethanol (about 0.3 g/ml). What could 
be the problem with the simulation using switch? By the way does gromacs 
use a different switch function depending on the forcefield?

Best regards,

Marcelo