Dear Dr. Case:
I have read your new great tutorial on Sander TI, but I am still
confused by a (maybe stupid) question on vdw part: why did you
integrate the dV/d_lam from 0.05 to 0.98 instead of 0 to 1? I ask
this question because I performed an almost the same perturbation
simulation with yours used 80 windows (both in water and vacuum),
i.e. clambda = 0, 0.0125, ... 0.9875, 1. I totally got 80 values of
dV/d_lam ( because when clambda = 1, Sander did not print out the
value of dV/d_lam). Then I integrated them with trapezoidal rule and
got the result of -5.37kal/mol in water and -3.80kal/mol in vacuum.
Clearly, my final result (-3.80-(-5.37) = 1.57kal/mol )of vdw
contribution of toluene's solvation free enery is much smaller then
yours(2.70kal/mol) and Shirts et al (2.45kal/mol). But the question
is, when I picked out the dV/d_lams of the following clambda (0.05,
0.15,...0.85,0.9, 0.9125... 0.9625, 0.975) and integrated them with
trapezoidal rule, I got the result of (-3.80-(-7.02) ) 3.22kal/mol .
So, could you give me some advice on my method?
My system including an uncharged toluene molecule and 11A TIP3P
waters. Before TI, 4ns equilibration was performed. Following are my
TI input files and 80 dV/d_lams. By the way, I also did this with
Gibbs dynamic windows method, the result is 2.94kal/mol. And my
result of charging free energy of toluene is :
-1.62(vacuum)-(+0.62)(water) = -2.24kal/mol. (50 windows and 40ps
eqili. and 40ps produ. ) I think the difference between our results
of this part is coming from the number of windows we used.
Best regards,
Chunhu
40ps equilibration.
&cntrl
imin = 0 ,
irest = 1 , ntx = 5,
ntb = 2 , ntp = 1,
ntt = 1 , iwrap = 1,
ntf = 2 , ntc = 2, tol = 0.000001,
temp0 = 300.0,
nstlim = 40000, dt = 0.001,
ntpr = 500, ntwr = 500, ntwx = 0,
icfe = 1, clambda = $clambda,klambda =6,
cut = 9., ntr = 1
&ewald
skinnb=1.
/
Hold the molecule fixed
50
RES 1 1
END
END
40ps data production
&cntrl
imin = 0 ,
irest = 1 , ntx = 5,
ntb = 2 , ntp = 1,
ntt = 1 , iwrap = 1,
ntf = 2 , ntc = 2, tol = 0.000001,
temp0 = 300. ,
nstlim = 40000, dt = 0.001,
ntpr = 500, ntwr = 500, ntwx = 0,
icfe = 1, clambda = $clambda,klambda = 6,
cut = 9., ntr = 1
&ewald
skinnb=1.
/
Hold the molecule fixed
50
RES 1 1
END
END
dV/d_lam of 80 clambdas :
47.0559
45.0056
42.1897
37.0196
33.0724
29.1933
26.7656
26.2567
22.4820
19.6780
17.4702
15.1580
13.5812
11.5357
9.3842
8.4599
7.4047
6.3753
3.5732
2.9784
1.6973
-0.1802
-0.5715
-1.6963
-1.4059
-2.5515
-5.3735
-4.4047
-5.3718
-7.5483
-9.2175
-8.6574
-9.6055
-8.8855
-10.1884
-11.3757
-13.1140
-10.6789
-12.7776
-13.0568
-12.2951
-14.2990
-12.7457
-15.0621
-14.7040
-14.3622
-16.3769
-13.0978
-16.9654
-17.0538
-18.8635
-18.1810
-19.0062
-19.7110
-20.8599
-17.6701
-19.1727
-20.1025
-19.3362
-22.2121
-20.1404
-19.1419
-27.0331
-18.8095
-23.5859
-22.4232
-21.4871
-18.9847
-19.2507
-20.5511
-20.4623
-18.4538
-17.2204
-15.4475
-15.1130
-13.8613
-16.4275
-12.8686
-9.6216
-5.8242
-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber.scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu
Received on Sun Jan 23 2005 - 10:53:00 PST