本帖最后由 牧生 于 2021-6-15 15:36 编辑
有学生想试试做金属表面和小分子的作用。
我照着培训班的教程也试了一下,按照教程,将金的参数加入到了top文件夹里面相应的itp里面。然后用VMD做了金板,延展,向真空加入水,直接进行模拟,整个过程都没问题。金板仍然保持为板状,原子都没有发生明显的移动。
我初步做的是金板-乙醇水溶液,希望将来可以举一反三,拓展到铁板-酸液和缓蚀剂。。但是在做金-乙醇水溶液就有问题了。我需要请教的如下:
① 使用VMD建立了4×4×2 nm的金板 (第一个方面的问题:如果我要建立VMD没有自带的,比如我需要建立一个铁板,尽管查到了铁的L-J参数 FE 26 0 0 A 0.01514608 0.0000022873928,但是下一步该如何像建立金板那样来建立铁板啊?虽然VMD有new material的选项,看起来应该是可以建立铁板的,但由于没有实例,理解起来有点困难。能否讲解一个实例来帮助理解,铁的unit cell pdb,和top,到底是一种什么样的文件)
经过一天折腾,第一个问题已经解决,使用MS建立晶胞很容易,另存为pdb格式,gromacs也是可以识别的,原子数量,以及盒子尺寸都是可以从MS读取出来的,
② 回到金-小分子水溶液的模拟上面来。
继续延展金板,出现真空层,gmx editconf -f Au.pdb -o Au_plate_box.gro -box 4.078 4.078 6
③然后使用insert命令加入2个ATB得到的乙醇分子, gmx insert-molecules -f Au_plate_box.gro -ci yichun-u.pdb -nmol 2 -o out_box.gro 得到图像;
④ 加水填充盒子的真空部分;
gmx solvate -cp out_box.gro -ptopol.top -o Au_wat.gro
此处有个小问题未解决:gmx solvate是向已经定义好尺寸的盒子里面加满水,但是为什么到盒子外面去了?
⑤ 使用gromos54a7_atb.ff力场,使用教程中的mdp进行npt平衡相模拟,此时就有第三个大的问题了。 始终提示有角度大于30度,不能继续。于是修改mdp中的dt=0.001,仍然提示角度大于30,不能运行。再减小dt=0.0005,可以跑起来了,但是得到的npt.gro图像如下,金原子跑散了,乙醇也跑到了中间位置。
npt如下:
define = integrator = md dt = 0.0005 ; ps (一般这里是0.001或者0.002,但是都跑不起来,改成0.0005可以跑) nsteps = 1000000 ; comm-grps = system energygrps = ; nstxout = 0 nstvout = 0 nstfout = 0 nstlog = 500 nstenergy = 500 nstxout-compressed = 1000 compressed-x-grps = system ; annealing = single single annealing_npoints = 2 2 annealing_time = 0 100 0 100 annealing_temp = 0 298.15 0 298.15 ; pbc = xyz cutoff-scheme = Verlet coulombtype = PME rcoulomb = 1.2 vdwtype = cut-off rvdw = 1.2 DispCorr = EnerPres ; Tcoupl = V-rescale tau_t = 0.2 0.2 tc_grps = SOL non-water ref_t = 298.15 298.15 ; Pcoupl = Berendsen pcoupltype = semiisotropic tau_p = 0.5 ref_p = 1.0 1.0 compressibility = 4.5e-5 4.5e-5 ; gen_vel = no gen_temp = 298.15 gen_seed = -1 ; freezegrps = freezedim = constraints = hbonds
⑥ 由于一直提示角度大于30,应该是没有做能量最小化导致的。如果我进行了能量最小以后再npt,结果会怎样呢。我找了教程中CNT_wat的能量最小化的mdp文件进行能量最小化,跑完以后,发现金原子还是全部都散了。 能量最小化的mdp文件如下: define = -DFLEXIBLE integrator = cg nsteps = 10000 emtol =200.0 emstep = 0.01 ; nstxout = 20 nstlog = 50 nstenergy = 50 ; pbc = xyz cutoff-scheme =Verlet coulombtype = PME rcoulomb = 0.9 vdwtype =Cut-off rvdw = 0.9 DispCorr = EnerPres ; constraints =none
现在我就陷入这个死循环了,能量最小化,或者直接npt,金原子都跑散了。但是和教程中的操作相比,仅仅是多加入了两个乙醇,为何会导致结果差异十万八千里。
经过一天多折腾,金板不会跑散了。
使用相同的方法建立铝板和表面活性剂水溶液,也不会跑散。
使用相同的方法建立铂板和表面活性剂水溶液,也不会跑散。
但是使用相同方法建立铁板和表面活性剂水溶液,铁就会跑散。
培训班好像讲过,只要LJ参数适当,那么,原子就不会跑散。现在,虽然我觉得我自己的操作可能有问题,但又察觉不到问题在哪。
2021.6.15 感觉是铁晶型的问题,用的是α铁的晶体,但用了γ铁的参数。。但是查了一阵,没查到γ铁的晶胞。。
大佬们,能否帮助找一下γ铁晶胞的pdb文件或cif啊,我觉得这个很常规,但是没找到。
|