手性性质的量子力学计算已迅速成为归属有机化合物(包括天然产物)绝对构型(absolute configuration, AC)的最流行方法。现在非专业人员很容易以黑盒模式使用含时密度泛函理论(TDDFT)计算电子圆二色性 (ECD)光谱。然而,不加批判的态度很容易给出错误的答案。本文介绍了在进行TDDFT计算ECD图谱时良好计算规范的讨论,通过最近文献中的几个例子来说明几个最重要的关键点。
计算能力的迅速发展[1],彻底改变了化学家和其他科学家归属/指认新分子绝对构型(AC)的基本方式[2]。此前完全依靠Bijvoet方法或化学关联法,随着手性光谱学的出现,尤其是电子圆二色性谱(electronic circular dichroism, ECD)的出现,最近这个问题已成为光谱学的主要目的。然而,长期以来,新化合物ECD光谱的解释完全基于经验或半经验方法,例如光谱比较或扇区规则。
因此,毫不奇怪,诸如激子手性之类的非经验分析方法能够迅速流行并激发了许多研究。然而,只有当手性光学性质的从头计算理论模拟,包括ECD、旋光、OR、振动CD (VCD)、拉曼光学活性(ROA)等等,对现实生活中的分子(例如天然产物和过渡金属络合物)来说成为了一种切实可行的方法的时候,手性光谱的革命或者说复兴才算真正到来。手性光谱的量子力学计算使得指认绝对构型成为可能,而无需任何参考系统或任何化学衍生方法,通常不需要建立对观测性质负责的分子理化机理(所谓的光学活性机理)。
就ECD而言,必须对激发态进行预测才能进行模拟,随着成本-效益非常高的含时密度泛函理论(TDDFT)方法的出现,已经向前迈出了一步[14]。高效的TDDFT可以在合理的时间内在台式PC上模拟中等大小分子的ECD谱,商业或免费软件的TDDFT ECD计算可获取性为非专家的广泛使用开辟了道路[15]。
硬币的另一面是,对计算工具的黑盒使用可能很容易对结果产生意外的负面影响。在AC指认时尤其如此,显而易见,“绝对构型研究中只有两个可能的结论:要么正确,要么错误的”[6]。本文的目的是分享一下我们对TDDFT ECD计算指认有机化合物绝对构型时良好计算规范的一些思考。我们首先介绍一下完成绝对构型指认的典型流程图,然后我们通过一些示例来讨论流程中的每个步骤,以突显每个步骤中最关键的内容。
01结果与讨论典型的计算流程
Figure 1. ECD计算流程图
对给定分子进行ECD计算的研究从输入结构的生成开始。尽管分子构象是在计算过程中确定的,但所有的其他结构特征必须事先知道,包括分子组成(所谓的“平面”结构),非常重要的是其相对构型(RC)分子中的各种手性元素。原则上,后者也可以通过手性光谱或几种方法组合来确定[17,18]。但是,我们强烈建议使用独立的工具,包括X-射线晶体学(只要有)和NMR,来尽可能地阐明RC。在大多数情况下,仔细评估对3D结构敏感的全套NMR数据,例如J-耦合(尤其是3JHH和3JCH)和NOE效应,会得出相对立体化学信息[19]。
必须强调的是,在涉及柔性分子时,应该使用“真实的”分子模型而不仅仅是平面结构来合理解释NOE数据。因此,我们假定:待进行绝对构型确证的分子应该结构已知、相对构型已经确认,并且通过分子建模软件生成了该分子的初始分子模型(具有正确的相对构型!)。ECD计算研究的典型步骤如下(参见图1中的流程图)[15,20]。
(1)在低计算水平下(通常是分子力学,MM),用诸如MMFF(Merck Molecular FF)等比较好的力场进行彻底的构象搜索,以研究构象系综。这可以通过诸如蒙特卡洛或分子动力学之类的构象分析方法来完成[22]。重要的是,在构象搜索过程中,所有可旋转键以及环上的键均应被搜索到,以免错过任何可能的构象[16]。
(2a)在中等计算水平对步骤(1)中发现的所有构象体进行初步的几何优化。通常在B3LYP / 6-31G(d)理论水平完成此操作。虽然如此,我们还是建议用更有效的泛函进行结构和能量计算[23],并尽可能考虑色散校正[24]和双ζ基组。Grimme等人[25]提出的PBEh-3c方法是另一种省时高效的方法,特别是对于大分子,可以得到高质量的几何形状。本步骤的目的是对MM-结构进行筛选,并获得有限数量的低能构象,其落在某一阈值(即5或10kcal/mol)内的能量窗口中。
(2b)在更高的计算水平上对(2a)中选择的结构进行更精细的几何优化。包括:(i)使用更大的基组,例如三重ζ或考虑极化甚至弥散函数。通常,def2-SVP或def2-TZVP[26]之类的Ahlrichs基组比6-31G(d)之类的Pople基组更适合;(ii)使用溶剂模型,例如可极化连续体模型(PCM)[27],COSMO[28]或SMD[29]。
(3)可选地,将步骤(2)中找到的低能构象集在更高的理论水平进行单点能计算(例如,使用比2b步骤更大的基组),或考虑溶剂模型,和/或进行频率计算以验证找到的是真正的能量最小点并获得自由能量或零点校正(ZPC)能量。频率计算必须与几何优化在同样的级别上进行。
(4)用内能、ZPC或自由能,计算步骤(2)或(3)获得的低能构象的将Boltzmann分布[20,30,31]。最后获得一系列构象以及各自的Boltzmann分布。那些Boltzmann分布高于某个阈值(例如1%)的构象在之后的步骤中会用到。
(5)检查所计算的构象集与构象依赖的NMR数据(J-耦合,NOE接触,环电流位移)之间的一致性。在我们看来,这是非常重要的一步,在当前文献中常常被忽略。因为实验数据是几种低能构象的平均值,所以对于非常柔性的分子,此计算/实验一致性检查可能会受阻或非常复杂。这种情况可能会得益于低温NMR实验,该实验在有利的情况下可以直接对构象异构体进行定量[32]。即使在最困难的情况下,也应始终验证分布最大的构象异构体与实验结果没有矛盾,反之亦然。如果出现矛盾,应检查构象分析是否有错并重新计算;另一方面,这也可能暗示了错误的结构指认,比如相对构型指认错误。
(6)计算在步骤(4)中发现的所有构象的吸收光谱和ECD光谱。这是整个过程中对计算要求最高的步骤,因此对计算级别的选择至关重要[9,13,14,20,33,34]。应该谨慎选择TDDFT方法泛函与基组。最重要的是,永远不要只使用单一的泛函/基组组合,而应始终探索不同的组合。至于泛函,至少应该测试两类泛函,即混合泛函(hybrid)和长程校正泛函(range-separated functionals)。越来越多的“精确”或HF交换的常用泛函有B3LYP(20%HF),PBE0(25%),M06(27%),BH&HLYP(50%)和M06-2X(54%)[35]。在ECD计算中,诸如CAM-B3LYP和ωB97X之类的长程校正泛函通常比混合泛函性能更好,因此应始终予以考虑。HF百分比和长程间隔都会影响绝对和相对跃迁能量[36],因此正确地重现能带序列及其能量是决定最终选择的因素之一(参见下面的”TDDFT泛函的选择”部分)。
至于基组,我们建议使用Ahlrichs的def2-SVP或def2-TZVP[26],或者更一般而言,建议使用具有足够宽泛极化函数的双或三-ζ基组。当里德堡态对ECD谱贡献较大时,弥散函数就很有必要[37]。通过使用下面的“TDDFT基础集的选择”一节中讨论的两种方法之一,可以轻松地验证基集收敛。通常,计算(和实验)的ECD谱图越弱,其对计算方法(包括泛函、基组,输入的几何形状和溶剂化模型)的依赖性就越强。吸收/ECD计算的本质是激发态计算:必须以一定方式选择要包含的激发态(root)数量,以使最终计算出的吸收光谱和ECD光谱在波长校正(见下文)后能够在低于观察到的最低波长之下时很好地扩展 。
(7)步骤6的ECD计算会给出波长或能量作为旋转强度函数的列表(棒图),对每个旋转强度进行展峰并在整个光谱范围内求和得到真正的ECD图谱。ECD通常使用高斯带形(Band-shape)展宽[38]。带宽的选择,也就是高斯曲线的标准偏差,通常是在最佳拟合的基础上凭经验完成的。所谓的FWHM(半高全宽)的合理值在0.1–0.3 eV的范围内。这同样适用于根据偶极或振幅强度生成吸收光谱的情况。
(8)用步骤4得到的构象玻尔兹曼权重对步骤7获得的光谱进行加权,然后将所有加权的光谱彼此相加,得到最终的玻尔兹曼平均光谱。
(9)将实验吸收/ECD光谱与步骤8计算的全波长范围光谱进行比较。此时,可以用一个波长偏移量来对校正计算跃迁能量的系统性过高/过低估计(有时称为UV校正,因为通过比较实验值和计算得出的UV吸收光谱更容易确定系统误差)。如果匹配良好,则假定的AC为正确的AC;如果镜像图谱配体良好,则相反的AC是正确的。由于一对对映体的ECD光谱从定义上来说时镜像对称的,因此无需对一对对映体都进行计算。
但是,如果进行计算,一定是精确的镜像关系。尽管计算和实验之间的比较通常是在视觉上进行的,但是我们还是建议使用工具进行定量比较[39,40]。定量比较给出相似因子,正确的AC应该比不正确的AC的相似性因子大得多(理想情况下是1比0)。如果在全波长范围内两种可能的对映异构体的比较结果均较差,则应重新检查计算流程,首先从步骤6中泛函/基组合的选择开始,直到检查第1步的构象生成步骤。
02几何优化与相对能量
Scheme 1. (S)-2-(3-fluoro-(p-tolyl)propyl)naphthalene (1) 以及altersolanol L (2)
与QM计算分子性质有关的主要问题之一是B3LYP杂化泛函的广泛且非关键性的使用[41]。尤其是通常在B3LYP/6-31G(d)理论水平进行几何优化和性质计算(例如相对内能)。这种幼稚的方法有一些缺点,在这里不能完全解决。简而言之,一些基准研究表明,B3LYP不一定是预测标准有机分子的结构和能量的最准确泛函[23]。其他大量可用的混合泛函和长程校正泛函,包括Truhlar的M0x和Head-Gordon的ωB97类泛函在总体上表现更好[42,43]。
此外,当非共价相互作用发挥作用时,对构象能的可靠预测需要包含色散,以作为对正常杂化或双杂化泛函的修正[24,44,45]。如 2-(3-fluoro-(ptolyl)propyl)naphthalene (Scheme 1)可能以展开或折叠形式存在,这取决于分子内π堆积的存在。但是,折叠构象完全被B3LYP计算所忽略,而当B2GPPLYP-D3用于真空几何优化时,折叠构象却是最稳定的构象(图2)。由于两个构象(展开和折叠)的计算得出的ECD谱图非常不同,因此它们的相对分布对计算得出的平均ECD谱图有很大的影响[10]。
Table 1. 在300 K下用五种不同方法估算altersolanol L的玻尔兹曼分布
a排序根据MMFF94能量,仅列出B3LYP/6-31G(d)相对能量在6kcal/mol以内的构象。b小基组:B3LYP/6-31(G);大基组:B3LYP/6-311G (d,p)。c ΔH:分布用内能计算;ΔG:分布用ZPC-自由能计算。dPCM模型的乙腈溶液。
至于基组,尤其是由于分子内基组叠加误差(BSSE)[46-48],例如6-31G(d)之类的双ζ-Pople型基组,其精度不足以准确地计算几何形状和相对能量。至少在单点或频率计算中应考虑使用更大的基组,以估计内能或自由能(流程图中的步骤2b)。诸如标识方法的解析度(resolution of identity method)[49]或球链近似(chain-of-spheres approximation)[50]之类的节省时间方法与大型基组集组合使用或许有用。
为了说明计算方法对构象分布(即使对于相对刚性的有机分子)产生多大的影响,表1报告了在300 K下用五种不同方法估算altersolanol L(2,Scheme 1)[51]的玻尔兹曼分布:
(a) 步骤1MMFF级别计算的能量;
(b)步骤2a在B3LYP / 6- 31G(d)理论水平计算的内能;
(c)在B3LYP /6-311 G(d,p)水平计算的内能;
(d)用PCM溶剂化模型考虑了乙腈的溶剂效应,在与c相同水平下计算的内能。
(e)与d相同水平计算的ZPC校正过自由能。
基组大小和溶剂化模型对构象分布的巨大 影响警告任何人:不要将B3LYP/6-31G(d)几何优化的结果作为输入结构就直接开始ECD(或其他性质)计算。