欢迎光临散文网 会员登陆 & 注册

吐血整理!ECD计算常见工作流程、知识点归纳资源整合(下)

2022-07-14 10:12 作者:唯理计算  | 我要投稿

在看本次推文前建议大家没看过合集上的先去看一下

吐血整理!ECD计算常见工作流程、知识点归纳资源整合(上),再来看本期内容。

小编在工作过程中,发现有初次接触手性化合物的老师以及在这个领域摸爬滚打多年的实验博士会对基础概念存在误解的情况。这期推文小编就对相关的基础概念和计算化学的计算流程进行一下梳理,以方便对此感兴趣的老师学习,以及相关专业的学生进行查漏补缺。

如果你初次接触手性化合物计算,或者只是在其他人的指导下使用Schrödinger或者等商业软件进行过一些自己也不理解原因的操作。那么阅读本文一定会让你有所收获。

理清概念

首先先说一些需要了解的基础概念。

第一个问题:D-(+)-甘油醛和(S)-2-丁醇分别是啥意思?

如果你对文献中的化学名称中的+-符号或者DLRS不了解是什么意思,那么需要补充一下这部分内容。

该部分内容属于基础有机课必讲的内容。就不去提化学老师讲这部分内容了,对这个内容不了解的老师和同学可以看一下清华大学李艳梅老师的视频,链接如下:

http://content.xuetangx.com/20440333X/video_5_1

http://content.xuetangx.com/20440333X/video_5_2

http://content.xuetangx.com/20440333X/video_5_3

http://content.xuetangx.com/20440333X/video_5_4

http://content.xuetangx.com/20440333X/video_5_5

http://content.xuetangx.com/20440333X/video_5_6

http://content.xuetangx.com/20440333X/video_5_7

http://content.xuetangx.com/20440333X/video_5_8

这些视频都很短小,不会花费大家太多的时间去观看,这里建议大家全看。如果实在没时间,请看678这三个视频。就算不信我,李艳梅老师的实力大家还是要相信的

(此时,小编不得不插播一条课程预告,感兴趣的同学可以自行报名参与,

有机化学课程,重视物理本质,淡化无益细节!)

第二个问题:构型构象是啥意思,有啥区别?

这里从同分异构体开始说起。同分异构体分为构造异构(又称结构异构)和立体异构(又称几何异构)

能够在当前课题中讨论的分子一般是被筛选到,立体异构无法确定的阶段。

立体异构又分为构象异构和构型异构。

以下重点来说说这两者的区别


构型构象的区别

构象异构中,同分异构体之间可以通过化学键的旋转而互相变化。请记住,此处不用打断任何共价键,仅通过旋转化学键就能在各个构象之间相互转换。此时我们大可以把这些就当做是一种物质就行。那些无法自由旋转的双键以及碳环不在讨论之中。

构型异构是原子在大分子中不同空间排列所产生的异构现象。大家平时耳熟能详的括顺反异构、对映异构就是构型异构的内容。简单来说,虽然长得差不多,如果在chemdraw里面画不带楔形键的普通键线式也是没啥区别的。但是一旦3D结构做出来,你会发现,你不去断掉某些共价键重连,你是得不到另外一种异构体的。

如果你平时在文章中看到了有人研究了结构分布的问题,那么大概率你看到的是人家确定了某某构型,然后研究了该构型下不同构象的分布比例。这个分布比例是怎么来的呢?简单来说就是研究不同构象的相对能量,然后带入公式,就可以知道不同构象在某某温度下的分布比例。具体操作方法见卢天老师的博文 http://sobereva.com/165 

大家看博文的时候注意这么一句:下面这张图出自笔者的研究文章Struct. Chem. 25, 1521 (2014),利用Boltzmann分布与相对自由能的关系给出了常温下致癌物benzo[a]pyrene diol epoxide不同构象的出现比例,对anti和syn两种构型是分别算的。对于这句话有不明白的同学,可以回头仔细理解一下构型和构象的区别,以及玻尔兹曼分布的计算前提。


计算EDC流程

接下来说说实际计算ECD时候的流程,这里推荐用惯了相关黑箱软件或者某些公司开发的工作流的同学仔细阅读一下。

各种,当然也可以说全部。全部的ECD计算谱图预测工作流程都是基于如下的流程图进行的(当然此处指的是靠谱的工作软件或者包装后的黑箱工作流)

该图来自文献Pescitelli, G., and Bruhn, T. (2016) Good Computational Practice in the Assignment of Absolute Configurations by TDDFT Calculations of ECD Spectra. Chirality, 28: 466– 474. doi: 10.1002/chir.22600.


第0步


首先要有确定的 relative configuration(RC),如果此时你连自己体系中的物质具体的构型都不确定,那么就不用往下面走了。

如果就是RC确定不清楚,此时应该做表征,比如能够培养出单晶,直接一步到位搞定构型问题。

想用核磁来测定的同学请注意,结构异构体的核磁是不同的,非对映异构体的核磁是不同的,但是对映异构体的核磁是一样的。构象不同也会影响核磁的结果,这无疑给表征带来了不少的影响。

这里能够排除多少可能就尽可能排除,不然这些不确定性都会导致后续的工作的工作量翻倍。

加入已经确定了前文所说的RC,那么开始ECD计算的第一步


第1步

Initial Conformational Analysis

这一步,大家通常会选用各种基于半经验或者力场的工具去做研究。具体的做法是,根据上一步得到的构型去生成足够的构象,然后再去进行优化和能量计算。生成构象这里有很多工具可以选择,但是基本原理无外乎两种

1)穷举或者是基于一定化学规则的穷举。穷举规则和大家说一下,也是某位知名杂志的副主编提到的穷举方法。单键转120-120-120°,双键转180-180°。这种防范很容易就出现天文数字个结构,但也是最全面的。此处说明一下,双键旋转,到底属于构象搜索还是构型异构过程,大家心里清楚是怎么回事就可以。我个人以及手中的教科书,是把双键旋转归类为顺反非对应异构(构型异构)。

2)使用MD方法进行足够次数的退火实验,收集产生的结构。此方法可以较为全面产生所需要研究的构象,有效规避了穷举法一不小心就天文数字个构象的缺点。这种方法缺点也是很明显的。搜索到的结构能不能包含全部低能量结构,这里要考验操刀计算的人的经验和能力。


收集了足够的构象后,此时需要使用半经验工具如MOPAC或者力场如MMFF94去对这些结构进行一定的优化计算。

此时可以基于这些计算得到的能量去进行一个简单地排序,筛选掉明显高能的结构,然后进行下一步的计算。


第2步

上一步的计算机别对于确定构象分布比例还是有些不够。我们应该在DFT计算的精度下去做一个更精确的优化计算来得到能量并进行排序。此处可以选择的工具基本就是Gaussian和ORCA这类的量化老牌工具。也有一些集成工具中会调用老掉牙的泛函和基组去做这个研究,我个人是极度不推荐的。我个人建议是对于第1步得到的结构在第2步这里也遵循一个从粗糙到精确的DFT计算方式。比如先M062X-D3/def2svp计算,算完排个序,然后再用def2TZVP算完排个序。因为还有个第3步,所以此处算到def2svp就可以了。大基组留着下一步用。这一步及以后的DFT计算都需要带上隐式溶剂模型了,建议SMD


第3步

文献中建议大家这一步针对上一步优化后的结构算一次单点能(基于更高级别的基组,比如def2TZVP),然后对上一步优化后的结构以上一步优化时候的泛函基组算一次freq任务。然后将高精度单点能和低精度得到自由能校正量相加,得到一个靠谱的自由能数值。

如果分子较小或者组内计算资源极其丰富,可以直接用M062X-D3/def2TZVP/SMD这种配置算opt+freq。作为计算从业者,我建议高精度单点能加低精度得到的Gcorr已经可以了。当然,关于如何经济且精确算出自由能,这里面门道很多,此处只是举了一个简单易行的方法。


第4步

上一步结束后你能够得到一堆构象的自由能,这个时候借助文初的Boltzmann分布研究工具就可以顺利得到不同构象的分布比例。个人经验:对于这种小分子体系基本上前5个结构就能够95%的结构了。


第5部

这时候我们已经得到了一系列经过计算的结构了,可以根据得到的结构回头再去和NMR数据进行一个对比。这一步状况最多,还有可能直接帮你排除掉一些RC。


第6部

计算在第4步中发现的所有构象的ECD谱。这是最重要的一步,使用TDDFT算过激发态的同学都知道,不同的泛函和基组会得到不同的谱图,这步选择错了方法会导致后面的对比什么的都没有了意义。我们不建议你只算一种泛函,如果有条件应该尝试更多的泛函。长程校正泛函(如CAM-B3LYP和ωB97X)和杂化泛函(如B3LYP(20%HF),PBE0(25%),M06(27%),BH&HLYP(50%)和M06-2X(54%))都试试。不同的泛函对于激发态计算的能量偏差见下图(图中单位eV)

MSE (top) and MAE (bottom) obtained for singlet Excited-states.

*该图来自文献Jacquemin, D., Mennucci, B., & Adamo, C. (2011). Excited-state calculations with TD-DFT: from benchmarks to simulations in complex environments. Physical chemistry chemical physics, 13(38), 16987-16998.

更多关于TDDFT泛函选择的信息,可以看一下本推文对应的原始文献

(10.1002/chir.22600)

至于基组,一般影响不大,别太拉胯就行,常规的6-311G*就可以了,能用到def2-TZVP那就更好了。如果体系中里德堡态对ECD谱贡献较大,此时必须对所用基组加弥散。高斯中对Ahlrichs基组不太好弄,需要自定义基组,这里不多展开说,建议直接用pople基组好弄些。当然,如果你是ORCA用户,可以直接调用加了弥散的Ahlrichs基组,这里必须要夸一下ORCA开发组,很注重用户体验。文献原文也给出了一些基组的建议,不是全部认同,大家可以自行判断选取。溶剂模型呢,可以SMD和PCM都试试,什么好用什么。极端情况下,还可能需要用显式溶剂模型,比如你在苯中研究某个含有大量苯环的分子的ECD,那么这其中溶剂和被研究分子的分子间相互作用是不可以忽视的。为了保证谱图完整性,还需要思考要算几个state,这个就是分子越大要算的越多,如果不清楚,可以先找一个构象尝试算一下,比如算个50个态,看看能不能覆盖住实验谱区间。


第7步

因为计算出来的state是孤立的波长或能量与旋转强度的数值,和我们实验看到的曲线还是有很大区别的,那么这个中间涉及到如何把这些数据转换为谱图。这其中涉及到展宽和加权平均问题。如果对展宽操作的原理感兴趣可以参考http://sobereva.com/224 建议不用过多在意,因为这个过程免费的Multiwfn能够非常智能进行处理,遇到展宽不对的时候,直接换个参数,多试几次就好了。

Multiwfn绘制加权平均谱图的操作,具体可以参考:

http://sobereva.com/383

其中讲解了如果把第6步得到的计算输出文件进行绘图,Multiwfn可以处理不同谱图的加权平均,也可以调整不同的展宽参数。

其实在研究天然化合物的期刊中,使用得更多的软件是SpecDis,用来做谱图拟合什么的,这个软件也是免费的,而且有着简单易懂的手册可供参考,建议大家上手练练。如果这个软件的使用教程,大家有需求,就留言说明一下,这里挖个坑,以后再补

最后再说一句,第12346步都可以使用Molclus联用Gaussian或者ORCA的方式轻松完成,省去很多手动操作的时间,科研就该是该偷懒的地方偷懒,该严肃的地方严肃。你们说是吧

关于Molclus的使用可以参考

http://www.keinsci.com/research/molclus.html

流程算是梳理完了,如果对其中的某些软件使用上有门槛跨不过去,比如高斯,建议报个我们的唯理计算的班学习一下。如果是对整套工作流的时间无从下手,第一推荐的是自己摸索,搞定后那个成就感可以很大的,第二推荐寻求懂这块的人帮你装好软件并铺好道路,以后组里只要把对应文件换个坐标什么的就好了。


吐血整理!ECD计算常见工作流程、知识点归纳资源整合(下)的评论 (共 条)

分享到微博请遵守国家法律