计算化学公社

 找回密码 Forget password
 注册 Register
Views: 3639|回复 Reply: 29

[综合交流] 一篇存在大量计算问题的Angew期刊上的文章的comment

[复制链接 Copy URL]

9

帖子

1

威望

891

eV
积分
920

Level 4 (黑子)

发表于 Post on 2022-11-21 09:46:58 | 显示全部楼层 Show all |阅读模式 Reading model
本帖最后由 Suchent 于 2022-11-21 13:59 编辑

前言
碳硼烷修饰的荧光材料体系具有特殊的光学活性,近年来有许多其相关的研究。然而,其中不乏一些文章在计算上存在较多问题。最近笔者读到今年Angew. Chem. Int. Ed.上新发表的一篇关于碳硼烷修饰分子的文章,其计算过程以及结果分析存在诸多不合理的地方,原文详见:Angew. Chem. Int. Ed. 2022, e20221347(https://onlinelibrary.wiley.com/doi/10.1002/anie.202213470)。从计算级别的选取,激发态的计算与表征,实验光学性质的对应到晶体结构表征等诸多方面都存在严重问题。此外,该文章还存在引用不当的问题,如完全不引用所使用的泛函和基组,在使用了Multiwfn的前提下却完全没有引用任何和Multiwfn相关的工作。而这样一篇文章甚至已经成为了hot paper。其中有许多问题也是相关领域内的常见错误,如果不加以指正势必会对该领域未来的研究工作带来不利影响。因此笔者认为有必要将其中的错误和问题一一指出。



1.     计算级别的选取
首先,作者在几何优化时使用无D3校正的B3LYP泛函。然而,该泛函的选取并不合理,因为该类体系具有明显的大共轭特征(见原文正文Figure 4)。而B3LYP对于大共轭体系的优化存在的问题sob老师在http://sobereva.com/557一文中就已经提到过,此处不再细说。更加离谱的时所使用的基组在正文和SI中存在明显的矛盾之处。在正文的Figure 4处基组描述为6-3+G(d,p),这个基组根本不存在。而在SI Figure S35中的基组又被写成写成了6-31G(d,p)。只能根据猜测作者使用的计算基组可能为6-31+G(d,p)或者6-31G(d,p)中的某一个。这类体系也没有使用弥散函数的必要性。且对于所使用的B3LYP泛函和Pople基组也没有任何引用。文中作者试图研究该类体系在溶液中的吸收和发光,却完全没有使用任何溶剂模型。另外,作者也没有提供计算和优化后的原子坐标,也完全没有说明是否验证结构无虚频。综上所述,该文章中计算所使用的构型和计算级别很难确保是合理的。

2.     激发态计算的问题


2.1  前线轨道
分析跃迁轨道是表征激发态的重要手段。然而这篇文章中的作者却在完全没有确认主要激发跃迁特征是不是HOMO->LUMO的情况下,便直接分析HOMO与LUMO的轨道分布来试图说明电子转移的特征。首先,跃迁轨道和HOMO,LUMO之间并不存在必然联系。其次,即使某个激发态真的是HOMO->LUMO主导的也未必就对体系的吸收光谱存在主要贡献。此外原文图片中还列出了HOMO-LUMO gap来试图说明gap和激发能高低的关联,这一点更是完全没有可比性,即便激发的确由HOMO->LUMO主导,激发能和gap间也存在激子结合能的差异。甚至,作者列举的HOMO和LUMO也并非是从基态得到的,在正文中标注的是S1态的,不具备任何说明性。SI中倒是标注了S0态的HOMO与LUMO,但由于没有给出跃迁轨道贡献信息因此即使强行分析也意义不大。具体的关于前线轨道与基态-激发态gap的关系sob老师也已经讲过,详见http://sobereva.com/543。更为正确的做法是先确定所有激发态中对吸收光谱有主要贡献的激发态,而后对非简并的激发态直接绘制主要的分子轨道等值面图,简并情况较强的可以考虑使用自然跃迁轨道分析或者sob老师提出的空穴-电子分析。另外不难看出,此处的分子轨道等值面图是由Multiwfn绘制,而这篇文章的作者从头到尾都没有说明自己绘制用的软件名称,且无论在正文或是SI中,都完全没有按照Multiwfn主页所要求的那样,引用过任何Multiwfn的工作。原文图片上列举的信息还存在一些其他问题,比如明明仅算了4个化合物,却在Figure S35用compound1-5的字样描述体系等等。总而言之,该类体系的激发态特征计算过程和机理分析存在明显的错误。



2.2  复算验证
针对以上体系的实验性质,我们使用了TD-RI-revDSD-PBEP86-D4/def2-QZVPP级别结合SMD隐式溶剂模型计算Car-[9]CPP激发态性质,试图以更为准确的理论研究解释这个体系的光学性质。使用的结构是以下载到的晶体结构为初猜,基态结构经过ωB97XD/def2-TZVP(-f)级别优化所得,激发态结构在此基础上经过ωB97XD/def2-SVP优化所得。优化后结构经振动分析验证为局域极小点无虚频。所有计算均使用Orca 5.0完成。

首先我们对基态结构利用双杂化泛函进行了垂直激发的计算,对前10个单线态激发态进行了含时密度泛函计算,从而得到了以下结果(以下对HOMO缩写为H,LUMO为L,仅对前五个可能对发光有主要贡献的激发态进行展示,其余激发态f很小):
  
State
  
Excited  Energy/eV
Wavelength/nm
f
Contributions
S1
3.739
331.6
0.5039
  
H  -> L
  
61.0%
  
H-1  -> L+1
  
6.2%
H  -> L+3
5.8%
S2
4.298
288.5
2.3280
  
H-1  -> L
  
27.5%
  
H  -> L+1
  
17.5%
  
H  -> L+4
  
7.7%
H-2  -> L+1
6.0%
S3
4.567
271.5
0.5550
  
H-1  -> L+1
  
20.0%
  
H-2  -> L
  
18.7%
  
H  -> L+3
  
11.2%
  
H-2  -> L+3
  
7.0%
H-3  -> L+1
5.8%
S4
4.873
254.4
0.2297
  
H-2  -> L+1
  
14.8%
  
H-1  -> L+3
  
10.9%
  
H-3  -> L
  
10.7%
H  -> L+8
7.6%
S5
4.931
251.4
0.05056
  
H  -> L+8
  
34.1%
H-5  -> L
7.1%

基态体系与某个激发态之间的振子强度越大,就越容易吸收相应频率的电磁波而跃迁到那个激发态上。根据表格我们可以看出,单线态的第二激发态振子强度明显在这五个态中是最高的,这说明了分子的激发和发射可能直接与第二激发态有关。

随后,我们计算了激发态的垂直发射过程。我们分别优化到势能面S1和S2的极小值点所对应的结构,然后利用双杂化泛函计算了垂直发射的能量。激发态为S1时,吸收和发射数值上结果与垂直激发所得结果相似,而H-1 ->L+1的贡献增大到了20.9%,此时其激发态寿命为τ=3.27×102ps。而激发态为S2时,我们发现S2的垂直发射波长为344.02mm,振子强度高达2.739,此时其主要贡献为H-1 -> L(37.3%)和H-> L+1(31.5%),根据寿命计算公式,此时其激发态寿命为τ=6.48×10-1ps。如果文章报道的实验数据正确且考虑计算误差,该化合物在THF中的最大吸收和激发态寿命都无疑更倾向于S2是它的主要激发态。此外该文章的紫外吸收光谱中,该化合物在主峰右侧明显有一个小肩峰的存在,二者的高度比值和f(S2)/f(S1)较为相似,在scalefactor为0.85时,与计算所得的峰几乎完全一致。这也说明了S1很可能并非对其发光性质是主导的。

此外,该体系不但主要激发态贡献很可能不像原文所描述那样来自S1,而且可以看出,在垂直发射中,S2两个跃迁的贡献是远大于S1的,在这种情况下像文章作者一样讨论HOMO和LUMO布居已经不具备任何意义,其讨论的结果无非是按照预设的结果去进行的。此外该体系明显存在共轭体系大幅弯折,键长拉长的情况,能否使用单参考方法去进行计算得到合理的结果都是未知的,反观原文,抛开分析数据过程的正确性不谈,仅是其使用的B3LYP/6-31G(d,p),本身对这类体系就是很难得到定量合理的结果的。

3.      晶体学问题
原文中还列举了作者得到的三个晶体结构,CCDC编号为2180089,2180090和2180091。本文的计算也是基于这些结构作为初猜的。但是在读取坐标文件时,我们发现了晶体数据的许多问题。2180089和2180090的信噪比非常低,几乎不超过5,而且二者使用了SQUEEZE方法来处理晶胞中的溶剂分子,但是二者根本没有在CIF文件或者SI中说明这一点,晶体报道的分子式与实际给出的也不符。此外我们根据作者写的分子式推测,作者大概率是认为包被的溶剂是DCM或者氯仿,但是去掉溶剂遮罩模型后根本找不到可以指认为Cl的电子残余峰,且溶剂遮罩模型给出的更像是石油醚中的组分。而2180091的分辨率只有0.99Å,这意味着该晶体中长度约为0.96Å的C-H键都是定不准的。此外2180089中C原子的热各向异性明显是异常的,观察其指令发现,原来是加了非常多的限制条件,而且我们不清楚作者为什么在2180089中使用了DISP指令,因为Mo靶的散射因子是已知的,无序自定义。最后作者在SI中说明了单晶测试条件,但是我们根本没有找到所说的Cu靶测试的晶体,相反作者使用的Ga靶根本没有提及。



4.     总结
综上所述,该文章不仅存在引用不当,计算细节前后矛盾等诸多问题,在对于体系发光机理的理论计算研究中更是存在许多误区,从所使用的计算级别到跃迁轨道的分析都存在大量问题。通过在更可靠的计算级别下对Car-[9]CPP体系进行理论计算研究,我们发现该体系的激发特征很可能并非是原文那样简单的从HOMO->LUMO的激发以及从S1->S0的发射过程,可能存在以S2为主导的发射模式。此外,原文中的晶体测量过程也缺乏重要实验和测量细节,所提供的晶体结构文件的准确性同样无法得到保证。总而言之,这篇文章存在较多问题,希望引起同行们的注意。

图3 TD-RI-revDSD-PBEP86-D4/def2-QZVPP级别结合SMD模型计算所得到的吸收光谱,采用Multiwfn绘制。

图3 TD-RI-revDSD-PBEP86-D4/def2-QZVPP级别结合SMD模型计算所得到的吸收光谱,采用Multiwfn绘制。

图2 原文所测得的实验吸收/发射光谱。

图2 原文所测得的实验吸收/发射光谱。

图 1. 原文所列举的HOMO-LUMO分布和能量,其中关于基组的描述前后矛盾。

图 1. 原文所列举的HOMO-LUMO分布和能量,其中关于基组的描述前后矛盾。

评分 Rate

参与人数
Participants 11
威望 +1 eV +29 收起 理由
Reason
jackfish + 5 赞!
zhangzb + 1 赞!
Marc0 + 5 我很赞同
超限制抱怨 + 2 谢谢分享
wuil69 + 2
cyx98 + 5 正解
mizu-bai + 1 とてもいい!
zjxitcc + 2 精品内容
sobereva + 1
CopperL + 5 233333,这都不是拉跨是拉ku了
KangJiaxin + 1

查看全部评分 View all ratings

在理想主义者眼中,什么都是马尔科夫链

4万

帖子

99

威望

4万

eV
积分
89873

管理员

公社社长+计算化学玩家

发表于 Post on 2022-11-21 13:51:35 | 显示全部楼层 Show all
很多实验为主加点计算的文章,而且是用了Multiwfn却不引用的,方法/基组也都该引的不引,我估计大概率都是花钱代算的,那帮毫无责任心的机构算出来什么玩意儿、是否真的靠谱,文章作者自己估计也并不真的清楚。

评分 Rate

参与人数
Participants 4
eV +12 收起 理由
Reason
mj117712 + 3 我很赞同
moritaichi + 3
超限制抱怨 + 3 我很赞同
StormSpirts + 3 我很赞同

查看全部评分 View all ratings

北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办最高质量的各种计算化学类培训:初级量子化学培训班基础(中级)量子化学培训班分子动力学与GROMACS培训班量子化学波函数分析与Multiwfn程序培训班,内容介绍以及往届资料购买请点击链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的最佳途径。培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入“北京科音”微信公众号获取培训最新消息、避免错过网上最有价值的计算化学文章!
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群“思想家公社QQ群”:1号:18616395,2号:466017436,3号:764390338,搜索群号能搜到哪个说明目前哪个能加,合计9000人,讨论范畴相同
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu
Money and papers are rubbish, get a real life!

173

帖子

0

威望

834

eV
积分
1007

Level 4 (黑子)

发表于 Post on 2022-11-21 14:24:40 | 显示全部楼层 Show all
其实我觉得用B3LYP去优化该分子不是最好的选择,但你也不能说它错,因为实验的这些文章里面是不关心具体的C-C键长,跟倾向于算性质,而算个HOMOLUMO什么的用B3LYP效果就很好,至于文章中的激发态确实如您所说的那样有问题

459

帖子

1

威望

3997

eV
积分
4476

Level 6 (一方通行)

发表于 Post on 2022-11-21 14:27:33 | 显示全部楼层 Show all
我遇到过不止一个,把B3LYP/6-31G(d,p)这种方法基组的组合叫机组,注意,是机组

9

帖子

1

威望

891

eV
积分
920

Level 4 (黑子)

 楼主 Author| 发表于 Post on 2022-11-21 14:36:36 | 显示全部楼层 Show all
七尺贱 发表于 2022-11-21 14:24
其实我觉得用B3LYP去优化该分子不是最好的选择,但你也不能说它错,因为实验的这些文章里面是不关心具体的C ...

我也没说B3LYP优化是“错的“,只是说并不合适这类大共轭体系。此外,我认为不关心键长不代表可以用不合适的泛函去优化基态构象,毕竟后续要用S0结构去算吸收光谱。
在理想主义者眼中,什么都是马尔科夫链

173

帖子

0

威望

834

eV
积分
1007

Level 4 (黑子)

发表于 Post on 2022-11-21 14:47:43 | 显示全部楼层 Show all
Suchent 发表于 2022-11-21 14:36
我也没说B3LYP优化是“错的“,只是说并不合适这类大共轭体系。此外,我认为不关心键长不代表可以用不合 ...

确实,后续用到激发态计算用B3LYP不合适

173

帖子

0

威望

834

eV
积分
1007

Level 4 (黑子)

发表于 Post on 2022-11-21 14:52:26 | 显示全部楼层 Show all
本帖最后由 七尺贱 于 2022-11-21 14:55 编辑
Suchent 发表于 2022-11-21 14:36
我也没说B3LYP优化是“错的“,只是说并不合适这类大共轭体系。此外,我认为不关心键长不代表可以用不合 ...

我看了下四个通讯没有一个搞计算的,应该让代算坑了

评分 Rate

参与人数
Participants 1
eV +1 收起 理由
Reason
zjxitcc + 1 也有可能是某个作者瞎算的

查看全部评分 View all ratings

6661

帖子

0

威望

4000

eV
积分
10661

Level 6 (一方通行)

发表于 Post on 2022-11-21 15:26:04 | 显示全部楼层 Show all
我觉得B3LYP优化这个结构未必不适用。B3LYP不适合优化诸如18碳环这样的大共轭体系是事实,但是这不意味着所有比18碳环大的共轭体系都不能用B3LYP优化。像这里的共轭体系主要是由联苯结构组成的,苯环和苯环之间的键级不会比1大太多,相比之下18碳环的全联烯共振式能量并不高,所以泛函稍微有一点SIE问题,就会导致算出来是全联烯结构。可以看一下你优化的wB97XD的结构和作者的结构有没有明显区别。
此外我觉得作者还有一个问题是没有做构象搜索,各个联苯结构的二面角有很多种排列组合方式,会极大地影响苯环之间的共平面性,进而影响光物理光化学性质。如果不做构象搜索的话,很难知道计算出来链长对分子性质的影响,是真的因为链长的变化,还是只是因为某个联苯C-C键的二面角偶然选取得不一致导致的。

评分 Rate

参与人数
Participants 2
eV +8 收起 理由
Reason
Marc0 + 3 我很赞同
七尺贱 + 5 我很赞同

查看全部评分 View all ratings

BDF(https://bdf-manual.readthedocs.io/zh_CN/latest/Introduction.html)、ORCA(https://orcaforum.kofo.mpg.de/index.php)开发团队成员

9

帖子

1

威望

891

eV
积分
920

Level 4 (黑子)

 楼主 Author| 发表于 Post on 2022-11-21 15:35:46 | 显示全部楼层 Show all
wzkchem5 发表于 2022-11-21 15:26
我觉得B3LYP优化这个结构未必不适用。B3LYP不适合优化诸如18碳环这样的大共轭体系是事实,但是这不意味着所 ...

感谢您的建议。对这个体系确实可能确实没有严重的SIE问题。但是原文作者甚至没有给出计算所使用的坐标,所以无法比较。。。此外原作者也没添加任何溶剂模型。构象搜索确实是一个问题,二面角太多了。
在理想主义者眼中,什么都是马尔科夫链

60

帖子

0

威望

1050

eV
积分
1110

Level 4 (黑子)

发表于 Post on 2022-11-21 19:21:50 | 显示全部楼层 Show all
七尺贱 发表于 2022-11-21 14:52
我看了下四个通讯没有一个搞计算的,应该让代算坑了

可以查一下他们有没有瓜丝的版权,让瓜丝来搞他们...

128

帖子

0

威望

968

eV
积分
1096

Level 4 (黑子)

发表于 Post on 2022-11-21 20:46:58 | 显示全部楼层 Show all
fatpig 发表于 2022-11-21 19:21
可以查一下他们有没有瓜丝的版权,让瓜丝来搞他们...

妙啊!!!
一蓑烟雨任平生

459

帖子

1

威望

3997

eV
积分
4476

Level 6 (一方通行)

发表于 Post on 2022-11-21 20:54:34 | 显示全部楼层 Show all
fatpig 发表于 2022-11-21 19:21
可以查一下他们有没有瓜丝的版权,让瓜丝来搞他们...

写的用G16算的,引的G09,滑稽

47

帖子

1

威望

669

eV
积分
736

Level 4 (黑子)

发表于 Post on 2022-11-21 21:05:12 | 显示全部楼层 Show all
这乱搞的风格比较像是找的“*算科技”做的代算 而且他们最近也在到处群发吹代算有发了angew的

评分 Rate

参与人数
Participants 2
eV +5 收起 理由
Reason
ABetaCarw + 3 233333 同感
shuihua + 2 我很赞同

查看全部评分 View all ratings

Sometimes naïve.

388

帖子

1

威望

2028

eV
积分
2436

Level 5 (御坂)

发表于 Post on 2022-11-22 09:25:12 | 显示全部楼层 Show all
fatpig 发表于 2022-11-21 19:21
可以查一下他们有没有瓜丝的版权,让瓜丝来搞他们...

Gaussian 09的site license,南京大学肯定是有的。负责此事的老师曾跟我们说过,他当初花几千块打了官司,抖出了当时瓜丝公司在中国大陆的代理商欺上瞒下的事情。瓜丝公司听闻后,取消了中国大陆的代理商。南京大学趁瓜丝公司没有定下新的中国大陆代理商的时候,直接向瓜丝公司购买了Gaussian 09的site license。据说瓜丝公司打了个折扣。总的来说,那位老师给南京大学省了钱。至于Gaussian 16,我就不清楚了,我只知道2017年夏天瓜丝公司在南京大学举办了workshop,当时用的是Gaussian 16,墨灵格的肖经理当时来过,给机房里装软件。至于瓜丝公司会不会照你说的搞,我更不清楚了。不过黎书华在那个workshop第一天给过致辞,感谢了瓜丝公司。当时黎书华提到国内只有两个课题组买了瓜丝公司的代码,他本人就是其中之一。

74

帖子

0

威望

829

eV
积分
903

Level 4 (黑子)

发表于 Post on 2022-11-22 10:39:17 | 显示全部楼层 Show all
感谢楼主细致的分享。
我觉得楼主可以像sob老师那样,把以上内容整理成一篇 paper 投出去。当然大概率会得罪人,但这么多瑕疵,着实让人看不下去。

评分 Rate

参与人数
Participants 1
eV +1 收起 理由
Reason
七尺贱 + 1 不建议,学术圈就这么小

查看全部评分 View all ratings

本版积分规则 Credits rule

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2023-2-2 09:01 , Processed in 0.375914 second(s), 26 queries .

快速回复 返回顶部 返回列表 Return to list