计算化学公社

标题: 求助：TiO2和Bi2O3的结构发生严重变形，无法维持刚性 [打印本页]

作者
Author: confidence 时间: 2022-10-11 15:24
标题: 求助：TiO2和Bi2O3的结构发生严重变形，无法维持刚性
本帖最后由 confidence 于 2022-10-12 10:43 编辑

各位老师好，我想研究蛋白质在TiO2和Bi2O3上的作用，所以首先建立了TiO2和Bi2O3的模型，在进行能量最小化时就出现了如图所示问题， (, 下载次数 Times of downloads: 47) (, 下载次数 Times of downloads: 42) ，由图可以看到，在进行能量最小化后（实际上emtol达不到1000，说明此时也是不够稳定），TiO2的结构严重变形。

本人按照sob老师http://sobereva.com/soft/Sobtop#FAQ8中提及的可能原因进行了排查，（1）mdp文件设置：并未发现明显问题（该文件跑蛋白质时能顺利跑到emtol=100）；（2）grompp时只出现warning原子名不匹配(如Ti20-Ti)，这个应当是可以忽略的；（3）结构文件：确认了gro文件和top文件里记录的分子信息一致，且和itp文件里的原子也是顺序一致，LJ参数由文献获取，确保了没有单位转换错误；（4）拓扑文件：利用sobtop生成的itp文件中包含了对应的bond 和angle参数（dihedrals由于问过Sob老师我全部删除了），并未发现明显的参数缺失。

基于以上检查我觉得是力常数不够大，不足以维持TiO2的刚性。（但是我对比了由sobtop自动猜的力常数和我利用Hessian生成的力常数，发现我获取的力常数是更大的，详细见附件）
因此我想请教一下各位老师，看看我的做法或者文件中究竟是哪里出现了错误导致这种情况呢？谢谢各位老师。

PS：涉及的文件都在附件中 (, 下载次数 Times of downloads: 18) ，Bi2O3也是一样的情况。
我尝试着将bond和angle的力常数都增大了十倍，发现跑起来的情况还是一样的。目前不知道究竟是哪里有问题，还请各位老师有空的时候看看呢。

作者
Author: sobereva 时间: 2022-10-12 12:51
把原子电荷都设为0再试试
并且[angles]里随机找一些angle项看看当前设的平衡键角和gro里的实际键角是否对应

作者
Author: confidence 时间: 2022-10-12 15:03
本帖最后由 confidence 于 2022-10-12 16:07 编辑

sobereva 发表于 2022-10-12 12:51
把原子电荷都设为0再试试
并且[angles]里随机找一些angle项看看当前设的平衡键角和gro里的实际键角是否对 ...

谢谢sob老师的回答，我试了一下。如果把原子电荷都设为0的话，能量最小化就能到emtol=100，但是同样的还是会和上图一样发生严重形变。
[angles]中我查找了发现，的确是有一部分平衡键角和gro中的实际键角不对应（偏差在几度到几十度都有），老师这种情况下怎么办呢？（考虑过自己修改，但是数目太多了，差的数值没有规律）
我尝试了一下TiO2不存在真空区的情况，发现是可以正常运行em且结构没有变形（但直接做md会报错pressure scaling more ....），那sob老师我能否将em后的结构进行冻结/位置限制来与其他物质作用呢？

作者
Author: sobereva 时间: 2022-10-12 19:22

confidence 发表于 2022-10-12 15:03
谢谢sob老师的回答，我试了一下。如果把原子电荷都设为0的话，能量最小化就能到emtol=100，但是同样的还 ...

让sobtop直接猜参数，这样平衡键角肯定跟结构文件对应

作者
Author: confidence 时间: 2022-10-12 22:21
本帖最后由 confidence 于 2022-10-13 09:21 编辑

sobereva 发表于 2022-10-12 19:22
让sobtop直接猜参数，这样平衡键角肯定跟结构文件对应

是的呢sob老师，我就是让sobtop自动猜参数后再去修改对应的力常数（此时发现都是对应的了）。老师我直接用sobtop猜参数的去跑也是不行，em和md都会报错，难道是因为我加真空后（比如z方向加真空），在z方向上原本的周期性结构就不存在了吗？（没加氢饱和）
不知道老师您还有什么建议吗？不知道从何查找错误了。感谢老师

作者
Author: sobereva 时间: 2022-10-13 10:54

confidence 发表于 2022-10-12 22:21
是的呢sob老师，我就是让sobtop自动猜参数后再去修改对应的力常数（此时发现都是对应的了）。老师我直接 ...

若给sobtop用的mol2本身就是带真空的，就不用考虑z方向因为真空区导致该有的键没了的问题

作者
Author: confidence 时间: 2022-10-13 11:17

sobereva 发表于 2022-10-13 10:54
若给sobtop用的mol2本身就是带真空的，就不用考虑z方向因为真空区导致该有的键没了的问题

老师您说的这个我不太理解了。如果我要用sobtop生成itp时给的mol2文件是带有真空的话，就与我做振动分析获取力常数时提供的mol2文件不完全一致了（会导致有些键的键结情况与晶体固有结构不同，如Ti-O是六个键，但是加真空后的mol2文件有些Ti-O是两个键）

作者
Author: sobereva 时间: 2022-10-13 11:23

confidence 发表于 2022-10-13 11:17
老师您说的这个我不太理解了。如果我要用sobtop生成itp时给的mol2文件是带有真空的话，就与我做振动分析 ...

前面我说了，用sobtop自动猜的参数，没振动分析的事

作者
Author: confidence 时间: 2022-10-13 11:32
本帖最后由 confidence 于 2022-10-13 17:07 编辑

sobereva 发表于 2022-10-13 11:23
前面我说了，用sobtop自动猜的参数，没振动分析的事

嗯嗯那我明白老师您的意思了。可是老师这样不就有点矛盾了么？我一开始是想利用振动分析来精确获取力常数从而将k值写入到我研究体系之中（之前问过sob老师您说对于TiO2这类体系做振动分析不能存在真空区），那如果我将处在真空下的TiO2用sobtop自动猜获取参数，大部分的数值和我之前精确算的数值应该是不能匹配的吧？（甚至原子类型都无法对应上，之前确定的是TI-O有6个键，O-Ti有3个键）

难道老师我要全部都用sobtop猜的参数来跑吗？这样会不会精度不够呢？
老师我又试了在TiO2在真空下用sobtop获取参数，的确这时候在真空下em和md都能运行，但是加上水运行就不行了。
此外老师我发现TiO2的itp文件用sobtop生成时fudgeLJ和fudgeQQ是0.5和0.833，而我主top中采用的charmm36对应值为1.0和1.0，老师我这是错误的力场混用了么？

作者
Author: sobereva 时间: 2022-10-14 11:28

confidence 发表于 2022-10-13 11:32
嗯嗯那我明白老师您的意思了。可是老师这样不就有点矛盾了么？我一开始是想利用振动分析来精确获取力常数 ...

本来这个slab就是刚性的，动力学过程中能维持住就够了，当前也不是刻意计算这个slab的弹性之类的对bonded参数依赖显著的问题，所以用sobtop直接猜的便可。

若没特殊情况，[defaults]一律结合AMBER的

作者
Author: confidence 时间: 2022-10-14 14:33

sobereva 发表于 2022-10-14 11:28
本来这个slab就是刚性的，动力学过程中能维持住就够了，当前也不是刻意计算这个slab的弹性之类的对bonded ...

好的谢谢sob老师的回答，那我的力常数就全都用sobtop猜的便行。

老师因为我的体系中还有蛋白质，打算用的是charmm36力场，那如果我这个TiO2的[defaults]用的AMBER的，最后我是保留AMBER的[defaults]吗？

作者
Author: sobereva 时间: 2022-10-14 14:35

confidence 发表于 2022-10-14 14:33
好的谢谢sob老师的回答，那我的力常数就全都用sobtop猜的便行。

老师因为我的体系中还有蛋白质，打算 ...

蛋白质用CHARMM比起AMBER没任何特殊好处

作者
Author: confidence 时间: 2022-10-14 14:57
本帖最后由 confidence 于 2022-10-14 15:04 编辑

sobereva 发表于 2022-10-14 14:35
蛋白质用CHARMM比起AMBER没任何特殊好处

嗯嗯，是的sob老师，之所以用charmm是因为TiO2的LJ参数在文献中找到说是适用于charmm的 (, 下载次数 Times of downloads: 39) 那老师如果我蛋白质换成amber力场的，我TiO2的LJ参数还能用找到的这个数值吗？

作者
Author: sobereva 时间: 2022-10-14 20:05

confidence 发表于 2022-10-14 14:57
嗯嗯，是的sob老师，之所以用charmm是因为TiO2的LJ参数在文献中找到说是适用于charmm的那老师如果我蛋白 ...

用于结合AMBER力场原理上不会有问题

作者
Author: confidence 时间: 2022-10-19 17:24

sobereva 发表于 2022-10-14 20:05
用于结合AMBER力场原理上不会有问题

嗯嗯sob老师，我按照您说的将蛋白的力场换成了AMBER，之后跑的过程中又出现了新的问题。
1、 (, 下载次数 Times of downloads: 40) 如图所示，水分子的H原子基本上全都朝向了TiO2，我看过sob老师回复过一些帖子说是水被极化了，我这个TiO2与水接触表面暴露的是O原子，是这个原因么？（但是我整体TiO2呈现电中性，不带电荷）2、如果是这个原因的话，我是不是应该对TiO2表面进行加氢饱和呢？对于这些晶体结构的用来跑分子动力学是否都应该饱和而不是简单扩胞后放真空就行呢？
3、 (, 下载次数 Times of downloads: 39) 最后我想研究的是异质结对蛋白的作用，测试后发现如图TiO2与Bi2O3在水环境下跑em的话立马就segmentation fault，而我将两者距离增加到0.5nm就可以正常跑em。那我想请教一下对于研究这种异质结的，在md中两者的距离应该多少比较合适呢？

作者
Author: sobereva 时间: 2022-10-19 18:22

confidence 发表于 2022-10-19 17:24
嗯嗯sob老师，我按照您说的将蛋白的力场换成了AMBER，之后跑的过程中又出现了新的问题。
1、如图所示， ...

1 TiO2的氧带明显负电荷，水的氢自发冲着氧是极其正常的。但目前看起来水有晶态特征，所以要注意考虑下控温的问题，说不定表面上那层水的温度很低。建议把水单独作为一个控温组

2 TiO2表面是什么样和晶型、晶面有很大关系，很多文章都做了讨论，注意恰当注意表面结构的选取。看比如Langmuir 2014, 30, 14832−14840

3 尽可能和实际情况一致。距离5埃从常理上来说太远。结合处的结构应当恰当考虑，可能需要用bond项描述，还可能需要利用第一性原理程序优化异质结接触部分

作者
Author: confidence 时间: 2022-10-20 12:25

sobereva 发表于 2022-10-19 18:22
1 TiO2的氧带明显负电荷，水的氢自发冲着氧是极其正常的。但目前看起来水有晶态特征，所以要注意考虑下控 ...

谢谢sob老师的回答
1、我尝试着把水单独作为控温组，简单跑了一下0.5ns的npt后发现水还是存在一样的情况
2、我这个TiO2就只是直接扩胞，并未进行切面，不知道这样可以不？还是说老师对于这种晶体结构必然要切面呢？
3、老师如果要对结合处结构bond项进行描述，是指将这整个异质结（TiO2/Bi2O3）的itp文件用sobtop生成么？

作者
Author: sobereva 时间: 2022-10-20 16:15

confidence 发表于 2022-10-20 12:25
谢谢sob老师的回答
1、我尝试着把水单独作为控温组，简单跑了一下0.5ns的npt后发现水还是存在一样的情况 ...

2 如果你是指直接对原胞平移复制，得到表面，那有没有意义就看你要研究什么晶面与水的作用了。

3 是

作者
Author: confidence 时间: 2022-10-24 17:46

sobereva 发表于 2022-10-20 16:15
2 如果你是指直接对原胞平移复制，得到表面，那有没有意义就看你要研究什么晶面与水的作用了。

3 是

谢谢sob老师的回答，经过这几天的尝试，又遇到了一些新的问题
1、跑TiO2与蛋白时，发现md结果如图 (, 下载次数 Times of downloads: 47) ，其pressure达到了-10000bar左右，与设定的1.013bar严重不符，看了相关的控压设置感觉也没发现哪里有问题，因此想麻烦老师看看是哪里的错误？（或者说这结果能用么？虽然个人觉得大概率不行）2、按照sob老师说的，我将TiO2与Bi2O3构建成异质结后 (, 下载次数 Times of downloads: 39) 用sobtop生成了itp文件，将其在真空中跑固定后跑em能正常运行（不固定的话该结构变形会比较严重，但还是能保持，不会太散）。但是我在盒子中加水后跑em或者md的时候立马就会报segmentation fault。我检查了对应的gro和top顺序是没错误的，单独跑也是能行的，但就不明白为啥加水了就不行，麻烦老师看看怎么回事呢？

ps：附件包含了两个模拟涉及到的mdp、gro文件 (, 下载次数 Times of downloads: 2) (, 下载次数 Times of downloads: 1)

麻烦各位老师帮我看看是怎么回事呢，谢谢各位老师

作者
Author: sobereva 时间: 2022-10-24 23:55

confidence 发表于 2022-10-24 17:46
谢谢sob老师的回答，经过这几天的尝试，又遇到了一些新的问题
1、跑TiO2与蛋白时，发现md结果如图，其pr ...

我没时间具体看你的文件。如果不是各向异性控压，用各向异性控压。如果已经是的话，调调压浴参数，或者尝试其它压浴

Bi2O3的情况，一点点加水，反复测试

作者
Author: confidence 时间: 2022-10-25 19:42
本帖最后由 confidence 于 2022-10-26 09:40 编辑

sobereva 发表于 2022-10-24 23:55
我没时间具体看你的文件。如果不是各向异性控压，用各向异性控压。如果已经是的话，调调压浴参数，或者尝 ...

谢谢sob老师的回复
1、我按照您说的改为各向异性控压后的确能顺利跑完了，而且不会存在上图那些洞。但是老师为什么对我这个体系还需要用各项异性呢？看您回复过的是对于膜体系需要各向异性（因为膜在垂直和水平侧差别很多）
2、我尝试把水删掉一部分如图， (, 下载次数 Times of downloads: 36) 就可以跑em，但是连emtol5000都到不了，拿去跑npt也是立即就segmentation fault。之后我又试着换了另一个晶面的Bi2O3来构建异质结， (, 下载次数 Times of downloads: 42) 此时仅仅加水的话是能顺利跑完em/npt/md的（虽然在em的过程中出现了一些pdb文件，这应该是不太合理的情况）。但是当我把蛋白质加进来后，想跑em就出现了许多的core文件 (, 下载次数 Times of downloads: 43) 。请问老师这又是怎么回事呢？
老师我又测试了一下，最终的TiO2/Bi2O3异质结在水环境中能正常跑完所有流程，TiO2/Bi2O3/蛋白质在真空下也能跑完所有流程。但是TiO2/Bi2O3/蛋白质在水环境中就会出现上图所示报错，这不科学呀，请问老师究竟是哪里有问题呢？

作者
Author: sobereva 时间: 2022-10-26 14:27

confidence 发表于 2022-10-25 19:42
谢谢sob老师的回复
1、我按照您说的改为各向异性控压后的确能顺利跑完了，而且不会存在上图那些洞。但是 ...

1 各项同性控压时盒子三个方向变化的幅度是相同的，显然这会导致由于slab在XY方向上撑着盒子，令Z方向没法自如地变化尺寸来消掉真空区

2 只能一点点调试

作者
Author: confidence 时间: 2022-10-26 16:10
本帖最后由 confidence 于 2022-10-26 16:22 编辑

sobereva 发表于 2022-10-26 14:27
1 各项同性控压时盒子三个方向变化的幅度是相同的，显然这会导致由于slab在XY方向上撑着盒子，令Z方向没 ...

嗯嗯谢谢sob老师，经过不断将水分子删减后我发现当我构建如图所示的初始模型时 (, 下载次数 Times of downloads: 36) ，就能顺利的把em/npt/md跑完。那老师我这样的模型是能够使用的把？
不过为啥会需要删减水分子成这样子我也还是云里雾里，不知道老师您这边能不能给个解释呢？

作者
Author: sobereva 时间: 2022-10-27 05:05

confidence 发表于 2022-10-26 16:10
嗯嗯谢谢sob老师，经过不断将水分子删减后我发现当我构建如图所示的初始模型时，就能顺利的把em/npt/md跑 ...

可以模拟

没有明确的解释，无非就是初始结构有不合理接触之类问题，或者mdp设置等等影响模拟稳定性问题

作者
Author: zhanging 时间: 2023-5-16 20:45
TIO2的参数是那篇文献中给出的，可以给一下吗？

作者
Author: confidence 时间: 2023-5-17 09:45

zhanging 发表于 2023-5-16 20:45
TIO2的参数是那篇文献中给出的，可以给一下吗？

Wu X, Hao P, He F, et al. Molecular dynamics simulations of BSA absorptions on pure and formate-contaminated rutile (1 1 0) surface[J]. Applied Surface Science, 2020, 533: 147574.

作者
Author: jim 时间: 2023-8-4 22:30

confidence 发表于 2023-5-17 09:45
Wu X, Hao P, He F, et al. Molecular dynamics simulations of BSA absorptions on pure and formate-co ...

太感谢了，也在做TiO2，请问您是用CHARMM生成TiO2拓扑文件的是吧？

作者
Author: confidence 时间: 2023-8-5 10:55

jim 发表于 2023-8-4 22:30
太感谢了，也在做TiO2，请问您是用CHARMM生成TiO2拓扑文件的是吧？

用sobtop生成itp的

作者
Author: jim 时间: 2023-8-5 11:54

confidence 发表于 2023-8-5 10:55
用sobtop生成itp的

谢谢！我见您最初选择的是CHARMM，但是sobtop中没有，请问您选择下面哪一项呢？
1 Assign UFF atom types for all undetermined atoms
2 Assign GAFF atom types first and then UFF, then go to next step
3 Assign GAFF atom types as much as possible
4 Assign AMBER atom types as much as possible
5 Assign atom types according to the rule in assign_AT.dat file
6 Manually set atom type for selected atoms
7 Make all atom types to undetermined state
9 Show all current atom types on screen
能否加您的QQ详细请教一下？谢谢

作者
Author: confidence 时间: 2023-8-5 16:24

jim 发表于 2023-8-5 11:54
谢谢！我见您最初选择的是CHARMM，但是sobtop中没有，请问您选择下面哪一项呢？
1 Assign UFF atom type ...

我最初选择CHARMM力场是想给蛋白生成拓扑文件，最后在sob老师指导下改用了AMBER力场。

至于TiO2的话，就是用了这里的第二项

作者
Author: jim 时间: 2023-8-8 10:33

confidence 发表于 2023-8-5 16:24
我最初选择CHARMM力场是想给蛋白生成拓扑文件，最后在sob老师指导下改用了AMBER力场。

至于TiO2的话， ...

请问：你的复合物TiO2的电荷怎么算的？求指点

欢迎光临计算化学公社 (http://bbs.keinsci.com/)