|
|
使用Amber的ABMD计算自由能。使用二维来推演一个C-O键的断裂与另一个的形成，colva如下
&colvar
cv_type = 'DISTANCE'
cv_ni = 2, cv_i = 8335, 8317
cv_min = 1.5, cv_max = 4 ! min is not needed for DISTANCE
resolution = 0.2 ! required for mode = FLOODING
&end
&colvar
cv_type = 'DISTANCE'
cv_ni = 2, cv_i = 8335, 8322
cv_min = 1.5, cv_max = 4 ! min is not needed for DISTANCE
resolution = 0.2 ! required for mode = FLOODING
&end
脚本控制文件如下
&cntrl
imin=0, ! Perform MD
ntx=5, ! Read coordinates and velocities from ASCII formatted rst7 coordinate file
irest=1 ! restart simulation
nstlim=50000, ! Number of MD steps
dt=0.002, ! Timestep (ps)
ntp=1, ! Use the Berendsen barostat for constant pressure simulation
ntc=2, ! Enable SHAKE to constrain all bonds involving hydrogen
ntf=2, ! Setting to not calculate force for SHAKE constrained bonds
ntt=3, ! Langevin thermostat
gamma_ln=2.0 ! Collision Frequency for thermostat
ig=-1, ! Random seed for thermostat
temp0=300 ! Simulation temperature (K)
............................略中间一些读写参数
&qmmm
qmmask=':509-513', ! Residue 509-512 should be treated using QM
qm_theory='DFTB3', ! Use the DFTB3 semi-empirical Hamiltonian
qmcut=14.0 ! Use 14 angstrom cut off for QM region
qmcharge= 3 ! The charge of the quantum section
qmshake=1 ! Shake QM H atoms
&end
&abmd
mode = 'FLOODING' ! sets the execution mode
monitor_file= 'FW_1step.out' ! sets the name of the file to which value(s) of reaction coordinate(s)
monitor_freq=50 ! the frequency of the output to the monitor_file
timescale=1 ! τF, the flooding timescale in picoseconds
umbrella_file='basic_1step.nc' ! biasing potential file name
selection_freq=50 ! positive integer number that sets the frequency of the resampling algorithm (in MD steps)
selection_constant=0.001 ! Parameter C is to determine how strong the selection mechanism is. If C is too large, all the walkers will be replaced with the most dominant one.If C is too small, there will be no killing/duplicating of walkers
selection_epsilon=0 ! positive real number (typically less than unity) that sets the stopping criterion parameter ε
wt_temperature=10000 ! The smaller the T’; the smoother/slower the convergence
wt_umbrella_file='wt_basic_1step.nc' !the file name of true biasing potential after modification by 1 + ( T /T’ ) in which T is thereference temperature of the system (temp0)
cv_file = 'cv.rst'
&end
运行完后我发现在我的初始结构还能得到能量值
cv#1 cv#2 E
......................
1.500000e+00 2.950000e+00 -4.077634e+00
1.500000e+00 3.000000e+00 -3.342485e+00
1.500000e+00 3.050000e+00 -2.698778e+00
1.500000e+00 3.100000e+00 -2.157272e+00
1.500000e+00 3.150000e+00 -1.817898e+00
1.500000e+00 3.200000e+00 -1.734833e+00
--------------------
但是随着cv#1 增大及 cv#2减少所有的值都是0
...........................
1.950000e+00 3.750000e+00 -9.039146e-49
1.950000e+00 3.800000e+00 -3.118966e-48
1.950000e+00 3.850000e+00 -3.550809e-48
1.950000e+00 3.900000e+00 -1.548028e-48
1.950000e+00 3.950000e+00 -2.153079e-49
1.950000e+00 4.000000e+00 0.000000e+00
2.000000e+00 1.500000e+00 0.000000e+00
2.000000e+00 1.550000e+00 0.000000e+00
2.000000e+00 1.600000e+00 0.000000e+00
2.000000e+00 1.650000e+00 0.000000e+00
2.000000e+00 1.700000e+00 0.000000e+00
2.000000e+00 1.750000e+00 0.000000e+00
请问是采样不充分吗？
还是我的参数设置有错误？
谢谢
|
|
发表于 Post on 2019-5-27 16:46:54
收藏 Add to favorites