From: Lela Vukovic (lvukov1_at_gmail.com)
Date: Tue May 18 2010 - 13:08:02 CDT
Hello Jun,
As Felipe pointed out, you need to increase:
1) the number of steps during which the equilibration of the system
(adjustment to the new value of lambda) is performed, defined by
alchEquilSteps,
2) the total number of steps (by which you increase the number of
steps during which the ensemble average is collected), defined by
numSteps.
I don't think that it is a good idea to perform minimization in the
middle of a simulation with finite temperature. After minimization,
your system's temperature may go to zero K - you can check log file to
see if this is the case. So, eliminate minimization, and increase the
equilibration time in order to give enough time to the thermalized
system to readjust to every new value of lambda.
For FEP method, your system needs to reach the equilibrium at every
lambda (equilibration time needs to be long enough for this), and then
your ensemble average collection time needs to be long enough to
result in proper sampling of configurations at the defined lambda.
Lela
On Thu, May 13, 2010 at 8:48 AM, Jun Zhang <coolrainbow_at_yahoo.cn> wrote:
> Hello:
>
> Thank you for your reply. I have just studied FEP for several days but I have encountered many problems. The "minimization", although I know it should appear, appears because without it I will get an infinite energy:
>
> FEP: RESETTING FOR NEW FEP WINDOW LAMBDA SET TO 0 LAMBDA2 0.1
> 188 FEP: WINDOW TO HAVE 100 STEPS OF EQUILIBRATION PRIOR TO FEP DATA COLLECTION.
> FEP: USING CONSTANT TEMPERATURE OF 298 K FOR FEP CALCULATION
> PRESSURE: 0 nan nan nan nan nan nan nan nan nan
> GPRESSURE: 0 nan nan nan nan nan nan nan nan nan
>
> however, with one minimization before each window will eliminate this problem. I don't understand why. And in fact, I have tried very long run but I still can't get a reasonable result.
>
> In my opinion, even if I use a short run(my system is very small), the results of the forward and backward should at least have an opposite sign, but my calculation yielded the same sign. Maybe there exists some stupid problems in my configuration file or system setup. I am still looking for it. Thank you for your ideas.
>
> Yours sincerely
>
>
> Jun Zhang
> Nankai University
> coolrainbow_at_yahoo.cn
>
>
> --- 10年5月13日，周四, Felipe Merino <felmerino_at_uchile.cl> 写道：
>
>> 发件人: Felipe Merino <felmerino_at_uchile.cl>
>> 主题: Re: namd-l: FEP yields the same results for forward and backward computation!
>> 收件人: "Jun Zhang" <coolrainbow_at_yahoo.cn>
>> 抄送: "NAMD MAILLIST" <namd-l_at_ks.uiuc.edu>
>> 日期: 2010年5月13日,周四,下午2:24
>> Hey,
>>
>> What you said in the email does not make much sense, but i
>> suppose that
>> your concern is about the difference between the forward
>> and backward
>> directions.
>>
>> I don't know why you make a minimization before every
>> window.
>>
>> Besides that you are using an extremely short calculation
>> time per
>> window which of course could explain the huge hysteresis
>>
>> best
>>
>> felipe
>>
>> Jun Zhang escribi贸:
>> > Hi Everyone!
>> >
>> > I am performing a free energy calculation using FEP by
>> NAMD2.7b2. When I perform a FEP of a molecular in vacuum, I
>> can get a reasonable result, eithor forward or backward.
>> However, when I solvate the molecular in aqueous, it is
>> surpringly that both direction yields the same results! that
>> is:
>> >
>> > lambda
>> forward
>> backward
>> > 0.10000 9.59495
>> 8.37989
>> > 0.20000 21.96420
>> 17.85860
>> > 0.30000 34.83400
>> 27.51460
>> > 0.40000 47.86530
>> 37.10170
>> > 0.50000 60.88580
>> 46.55320
>> > 0.60000 70.98240
>> 52.97450
>> > 0.70000 80.99340
>> 59.52030
>> > 0.80000 90.86700
>> 67.15570
>> > 0.90000 100.75300
>> 77.51150
>> > 1.00000 110.93300
>> 93.77790
>> >
>> > I think there must be something seriously wrong in my
>> configuration. Yet I cannot find out. My conf file is:
>> >
>> >
>> #############################################################
>> > ## JOB DESCRIPTION
>>
>>
>> ##
>> >
>> #############################################################
>> >
>> > # Alchemical transformation to compute the free erengy
>> of
>> > # A-CH3 <-----> A-H in a Water
>> Box
>> >
>> >
>> >
>> #############################################################
>> > ## ADJUSTABLE PARAMETERS
>>
>> ##
>> >
>> #############################################################
>> >
>> > structure
>> ./a_CH3_H_aq.psf
>> > coordinates
>> ./a_CH3_H_aq.pdb
>> >
>> > set temperature 298
>> > set outputname a_CH3_H_aq_eq
>> >
>> > # Important for restart a simulation
>> > firsttimestep 0
>> >
>> >
>> >
>> #############################################################
>> > ## SIMULATION PARAMETERS
>>
>> ##
>> >
>> #############################################################
>> >
>> > # Input
>> > paraTypeCharmm on
>> > parameters
>> ./par_all27_prot_lipid_MAB.prm
>> > temperature
>> $temperature
>> >
>> >
>> > # Force-Field Parameters
>> > exclude
>> scaled1-4
>> > 1-4scaling 1.0
>> > cutoff
>> 12.
>> > switching
>> on
>> > switchdist 10.
>> > pairlistdist 13.5
>> > LJcorrection yes
>> >
>> > # Integrator Parameters
>> > timestep
>> 1.0 ;# 1fs/step
>> > rigidBonds all
>> ;# needed for 1fs steps
>> > nonbondedFreq 1
>> > fullElectFrequency 2
>> > stepspercycle 10
>> >
>> >
>> > # Constant Temperature Control
>> > langevin
>> on ;# do langevin dynamics
>> > langevinDamping 1
>> ;# damping coefficient (gamma) of 5/ps
>> > langevinTemp $temperature
>> > langevinHydrogen off ;#
>> don't couple langevin bath to hydrogens
>> >
>> >
>> > # Periodic Boundary Conditions
>> > cellBasisVector1 22.
>> 0. 0.
>> > cellBasisVector2
>> 0. 20. 0.
>> > cellBasisVector3 0.
>> 0 22.
>> > cellOrigin 24.95806884765625 -19.506383895874023
>> -13.278444290161133
>> >
>> > wrapAll
>> on
>> >
>> >
>> > # PME (for full-system periodic electrostatics)
>> > PME
>> yes
>> > PMEGridSizeX 25
>> > PMEGridSizeY 25
>> > PMEGridSizeZ 25
>> >
>> > # Constant Pressure Control (variable volume)
>> > useGroupPressure yes ;# needed for
>> rigidBonds
>> > useFlexibleCell no
>> > useConstantArea no
>> >
>> > langevinPiston on
>> > langevinPistonTarget 1.01325 ;# in bar
>> -> 1 atm
>> > langevinPistonPeriod 100.
>> > langevinPistonDecay 50.
>> > langevinPistonTemp $temperature
>> >
>> >
>> > # Output
>> > outputName
>> $outputname
>> >
>> > restartfreq
>> 2000 ;# 1000steps
>> = every 2ps
>> > dcdfreq
>> 2000 ;# 1000steps
>> = every 2ps
>> > xstFreq
>> 2000 ;# 1000steps
>> = every 2ps
>> > outputEnergies 2000
>> ;# 1000steps = every 2ps
>> > outputPressure 2000
>> ;# 1000steps = every 2ps
>> >
>> >
>> >
>> #############################################################
>> > ## EXTRA PARAMETERS
>>
>> ##
>> >
>> #############################################################
>> >
>> >
>> >
>> #############################################################
>> > ## EXECUTION SCRIPT
>>
>> ##
>> >
>> #############################################################
>> >
>> >
>> > alch
>> on
>> > alchType
>> FEP
>> > alchFile
>> a_CH3_H_aq.fep
>> > alchCol
>> B
>> > alchOutFile
>> a_CH3_H_aq_forward.fepout
>> > alchOutFreq
>> 10
>> >
>> > alchVdwLambdaEnd 1.0
>> > alchElecLambdaStart 0.5
>> > alchVdWShiftCoeff 6.0
>> > alchDecouple
>> off
>> >
>> >
>> > alchEquilSteps 100
>> > set numSteps
>> 500
>> > set windowwith 0.1
>> >
>> > for {set i 0.0} { $i < 1.0 } {set i [expr
>> $i+$windowwith]} {
>> >
>> > alchLambda $i
>> > alchLambda2 [expr $i+$windowwith]
>> >
>> > minimize 100
>> > run $numSteps
>> >
>> > }
>> >
>> > Can anyone give me some suggestions? I can provide
>> anything detailed. Thank you!
>> >
>> >
>> > Jun Zhang
>> > Nankai University
>> > coolrainbow_at_yahoo.cn
>> >
>> >
>> >
>> >
>> >
>> >
>>
>>
>
>
>
>
>
This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:55:47 CST