Re: FEP calculations: optimizing the dlambda, number of steps

From: Shyno Mathew (sm3334_at_columbia.edu)
Date: Thu Dec 12 2013 - 09:37:22 CST

Hello Jerome,

1. Yes, I did. here are the contents of my .conf file:

readexclusions no
exclude scaled1-4
1-4scaling 0.833333333
scnb 2.0

temperature $temperature

# OUTPUT FREQUENCIES

outputenergies 10
outputtiming 10
outputpressure 10
restartfreq 10
XSTFreq 10
dcdfreq 10

# OUTPUT AND RESTART

outputname 05collab
restartname 05collab

binaryoutput yes
binaryrestart yes

# CONSTANT-T
langevin on
langevinTemp $temperature
langevinDamping 5.0
langevinHydrogen off

# Periodic Boundary Conditions
cellBasisVector1 124.16300010681152 0. 0.
cellBasisVector2 0. 82.71599864959717 0.
cellBasisVector3 0. 0. 128.3490011692047
cellOrigin -20.06950092315674 -22.470999240875244 46.29550063610077

# PME
PME yes
PMEGridSizeX 128
PMEGridSizeY 96
PMEGridSizeZ 160

# WRAP WATER FOR OUTPUT

wrapAll on

# CONSTANT-P

langevinPiston on
langevinPistonTarget 1.01325 ;# in bar -> 1 atm
langevinPistonPeriod 100.
langevinPistonDecay 50.
langevinPistonTemp $temperature

#StrainRate 0.0 0.0 0.0
useGroupPressure yes

useFlexibleCell no
useConstantArea no

# SPACE PARTITIONING

stepspercycle 10
margin 3.0

# CUT-OFFS

switching on
switchdist 10.0
cutoff 12.0
pairlistdist 13.5

# RESPA PROPAGATOR

timestep 1.0
fullElectFrequency 2
nonbondedFreq 1

# SHAKE

rigidbonds all

# FEP PARAMETERS

alch on
alchType FEP
alchFile 04collab_wet.fep
alchCol B
alchOutFile 05collab.fepout
alchOutFreq 10
alchVdwLambdaEnd 1.0
alchElecLambdaStart 0.5
alchVdWShiftCoeff 5.0
alchDecouple off
alchEquilSteps 10000

set Lambda0 0.0
set dLambda 0.05
while {$Lambda0 < 1.0} {
        alchLambda $Lambda0
        set Lambda0 [expr \$Lambda0 + \$dLambda]
        alchLambda2 $Lambda0
        run 20000
}

2. Since I am replacing hydrogen with a bulky group at few locations, can
it cause the water molecules to be pushed away?

thanks,
Shyno

On Thu, Dec 12, 2013 at 4:05 AM, Jérôme Hénin <jerome.henin_at_ibpc.fr> wrote:

> Did you enable a barostat?
>
> Jerome
>
> ------------------------------
>
> to make my question more clear:
> I have a cubic water box, but as the FEP progresses, shape of the water
> box becomes more irregular.
> All the water molecular are still there, but boundaries look more
> distorted!!
>
>
>
> On Wed, Dec 11, 2013 at 4:00 PM, Shyno Mathew <sm3334_at_columbia.edu> wrote:
>
>> Hi all,
>> Initially while doing FEP, I wasn't outputting the dcd file!!
>> Now, I did a short run with dcdfreq specification. When I look at the
>> trajectory, I see the water box is chopped little on each corner.
>> In the beginning the water box looks fine (rectangular shape). As FEP
>> progresses, some water molecules are disappearing.
>> Does this mean something is wrong with the FEP calculations I did?
>>
>> Look forward to hearing an answer!!
>>
>> thanks,
>> Shyno
>>
>>
>>
>> On Tue, Dec 10, 2013 at 11:28 AM, Shyno Mathew <sm3334_at_columbia.edu>wrote:
>>
>>> Dear all,
>>> I would like to confirm the FEP part in my .conf file is correct?
>>> 1. Here is how it looks like for the forward FEP :
>>> alch on
>>> alchType FEP
>>> alchFile 04collab_wet.fep
>>> alchCol B
>>> alchOutFile 05collab.fepout
>>> alchOutFreq 100
>>> alchVdwLambdaEnd 1.0
>>> alchElecLambdaStart 0.5
>>> alchVdWShiftCoeff 5.0
>>> alchDecouple off
>>> alchEquilSteps 5000
>>>
>>> set Lambda0 0.0
>>> set dLambda 0.05
>>> while {$Lambda0 < 1.0} {
>>> alchLambda $Lambda0
>>> set Lambda0 [expr \$Lambda0 + \$dLambda]
>>> alchLambda2 $Lambda0
>>> run 15000
>>> }
>>>
>>> 2. Now for the reverse FEP everything is same except the following:
>>> set Lambda0 1.0
>>> set dLambda -0.05
>>> while {$Lambda0 > 0.0}
>>>
>>> thanks in advance for your help,
>>>
>>> best,
>>> Shyno
>>>
>>>
>>> On Tue, Dec 3, 2013 at 6:12 PM, Shyno Mathew <sm3334_at_columbia.edu>wrote:
>>>
>>>> Hello Jerome,
>>>> Is there a way to make sure my system is equilibrated at each strata?
>>>> Currently, I do 10000 steps of equilibration followed by 10000 steps of
>>>> data collection. But I am not sure 10,000 steps of equilibration is enough?
>>>> Also, what method/algorithm namd uses to do FEP?
>>>> I read BAR is another method to calculate free energy.
>>>> Since the ParseFEPplugin was showing us some errors, we are trying to
>>>> write a script that calculates BAR.
>>>>
>>>> thanks,
>>>> Shyno
>>>>
>>>>
>>>>
>>>> On Fri, Nov 22, 2013 at 1:09 PM, Shyno Mathew <sm3334_at_columbia.edu>wrote:
>>>>
>>>>> Hi Jerome,
>>>>> I was using "minimize" incorrectly.
>>>>> I meant for each window I have 5000 steps of equilibration followed by
>>>>> 10,000 steps of data collection.
>>>>> So I should try increasing the equilibration steps.
>>>>>
>>>>> thanks,
>>>>> Shyno
>>>>>
>>>>>
>>>>>
>>>>> On Fri, Nov 22, 2013 at 12:01 PM, Shyno Mathew <sm3334_at_columbia.edu>wrote:
>>>>>
>>>>>> thanks again Jerome. Sorry to bother with multiple questions.
>>>>>> Here is how the FEP section of my .conf file looks like (I saw it in
>>>>>> one of the FEP tutorials).
>>>>>> So for each window it does 5000 step minimization followed by 10,000
>>>>>> step run. What you are suggesting is just do minimization at the very first
>>>>>> window only?
>>>>>>
>>>>>> alch on
>>>>>> alchType FEP
>>>>>> alchFile 04collab_wet.fep
>>>>>> alchCol B
>>>>>> alchOutFile 05collab.fepout
>>>>>> alchOutFreq 100
>>>>>> alchVdwLambdaEnd 1.0
>>>>>> alchElecLambdaStart 0.5
>>>>>> alchVdWShiftCoeff 5.0
>>>>>> alchDecouple off
>>>>>> alchEquilSteps 5000
>>>>>>
>>>>>> set Lambda0 0.0
>>>>>> set dLambda 0.05
>>>>>> while {$Lambda0 < 1.0} {
>>>>>> alchLambda $Lambda0
>>>>>> set Lambda0 [expr \$Lambda0 + \$dLambda]
>>>>>> alchLambda2 $Lambda0
>>>>>> run 15000
>>>>>> }
>>>>>>
>>>>>> thanks,
>>>>>> Shyno
>>>>>>
>>>>>>
>>>>>>
>>>>>> On Fri, Nov 22, 2013 at 11:49 AM, Jérôme Hénin <jerome.henin_at_ibpc.fr>wrote:
>>>>>>
>>>>>>> ----- Original Message -----
>>>>>>> > Thanks for the reply Jerome.
>>>>>>> > Now, for each window I have 5000 minimization followed by 10000
>>>>>>> steps. I
>>>>>>> > did try going to 10^6 steps, but didn't give significant
>>>>>>> difference in free
>>>>>>> > energy.
>>>>>>>
>>>>>>> In general I wouldn't recommend minimizing between windows. You have
>>>>>>> to thermalize to room temperature again every time, and discard the
>>>>>>> additional equilibration time.
>>>>>>>
>>>>>>> > So just to be clear, in the previous message by saying "increase
>>>>>>> the
>>>>>>> > sampling time for just those two by running extra simulations" you
>>>>>>> meant
>>>>>>> > going to smaller dlambda values between 0.45 and 0.55?
>>>>>>>
>>>>>>> I meant running more MD for the same lambda and dlambda values.
>>>>>>>
>>>>>>> > Also, does it make sense to do FEP just between 0.45 and 0.55 and
>>>>>>> then use
>>>>>>> > it with data for all other windows?
>>>>>>>
>>>>>>> Yes, you can merge any amount of data from separate runs. How to do
>>>>>>> it properly depends on what free energy estimator you use. For the FEP
>>>>>>> estimator, it's just a matter of taking the log of a weighted average of
>>>>>>> exponenials of deltaGs.
>>>>>>>
>>>>>>> Jerome
>>>>>>>
>>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> Shyno Mathew
>>>>>> PhD student
>>>>>> Department of Chemical Engineering
>>>>>> Columbia University
>>>>>>
>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> Shyno Mathew
>>>>> PhD student
>>>>> Department of Chemical Engineering
>>>>> Columbia University
>>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> Shyno Mathew
>>>> PhD student
>>>> Department of Chemical Engineering
>>>> Columbia University
>>>>
>>>
>>>
>>>
>>> --
>>> Shyno Mathew
>>> PhD student
>>> Department of Chemical Engineering
>>> Columbia University
>>>
>>
>>
>>
>> --
>> Shyno Mathew
>> PhD student
>> Department of Chemical Engineering
>> Columbia University
>>
>
>
>
> --
> Shyno Mathew
> PhD student
> Department of Chemical Engineering
> Columbia University
>
>
>

-- 
Shyno Mathew
PhD student
Department of Chemical Engineering
Columbia University

This archive was generated by hypermail 2.1.6 : Tue Dec 31 2013 - 23:24:05 CST