反应密度泛函理论的构建与初步应用

唐伟强，谢鹏，徐小飞，赵双良，

（1 华东理工大学化学工程联合国家重点实验室，上海200237； 2 广西大学化学化工学院，广西石化资源加工及过程强化技术重点实验室，广西南宁530004）

引 言

能源、资源和环境是人类社会赖以生存和发展的基础。随着人口的增加和社会的发展，人类消耗的能源和资源急剧增加。能源利用效率低、资源浪费严重、生产过程不够清洁，不仅引起一些资源和化石能源的快速枯竭，污染物的排放也快速增长，环境污染加重，大量的温室气体排放更是引起全球气候变化、极端灾害性天气频发，成为全球政治的焦点。在这样的大形势下，发展绿色化工、提高过程效率，降低能耗和污染具有重要的社会意义和经济价值。

根据“绿色化学十二原则”，绿色化工的理念体现在反应路径高效环保、能耗低、反应速率快、反应转化率和选择性高[1]。大多数化学反应都是在溶剂中发生的，溶剂能较大程度上调控化学反应的效率。例如，Han 等[2]发现以超临界CO2为溶剂可使苯酚加氢制备环己酮反应的产率和反应速率大大提高，解决了传统反应路线选择性差、效率低这一挑战性难题。最近，Kong 等[3]发现溶剂效应对脱羧/羧基化反应的实现是至关重要的。此外，如图1所示，实验发现混合溶剂的组成及其配比对化学反应有着显著的影响[4−5]。选择合适的反应溶剂已成为反应工程中调控化学反应选择性和转化率的控制步骤。

图1 混合溶剂对化学反应的影响：（a）量子点在TiO2膜上的吸附量严重依赖分散的溶剂[4]；（b）1,2−丙二醇在1,4−二氧六环（DIO）水溶液中脱水生成丙醛，在二甲基亚砜（DMSO）水溶液中生成丙酮，反应速率在不同混合比例下相差较大[5]Fig.1 The effect of mixed solvents on chemical reactions:(a)The absorbance of quantum dots adsorbed on the TiO2 film is heavily dependent on the dispersed solvent[4];(b)1,2−propanediol was dehydrated in 1,4−dioxane(DIO)aqueous solution to form propanal and in dimethyl sulfoxide(DMSO)aqueous solution to form acetone,and the reaction rates differed greatly under different mixed solvent ratios[5]

溶剂特性（包括种类、结构和流动特征等）对反应的影响可从两个角度来剖析。从反应自由能来看，反应自由能分布受到溶剂的分子作用力影响。当溶质（反应物或产物）在溶剂中溶解时，溶剂和溶质之间产生分子作用力[6]。这些作用力的存在使得反应物分子在反应过程中改变原来的结构而转化为另外一种物质时，必然受到溶剂和其他溶质的影响。溶剂特性、溶剂的微观结构等强烈地影响反应物和过渡态的稳定性，从而影响反应过程和反应速度，影响反应的活化能[6−8]。同一反应使用不同的溶剂，反应效果相差甚大[9]。这是因为不同溶剂中反应体系微环境不同。在化学反应工程中，体系总的反应速率常数（kT）可以表示为：

式中，kӨ是本征反应速率常数，Δk 是反应微环境及其结构演化对反应速率常数的影响。因此可通过改变微环境调控Δk，从而提高反应速率。如表1 所示，在溴丁烷和氰化钠的取代反应C4H9Br +NaCNC4H9CN + NaBr 中，溶剂发挥着重要的调控作用[10]：溴丁烷和氰化钠可以发生取代反应，但是如果反应溶剂为甲醇和水时，反应虽可以进行，但反应速率很慢，产率低。若采用二甲基亚砜、二甲基甲酰胺和乙腈作反应溶剂时，其反应速率比以甲醇作溶剂时分别快1300、2800 和5000 倍。从反应−传递耦合的观点来看，化学反应速率还与反应底物浓度、分子扩散速率和催化剂种类、结构和形貌等因素密切相关。而反应溶剂特性既决定了反应底物分子的饱和溶解浓度，也影响了其扩散速率，进而影响反应−分子传递耦合的匹配程度，此外，反应溶剂流动速率还影响产物分子的界面脱附、传递及反应平衡走向等。以双氧水的绿色制备为例，氢气和氧气通过催化反应生成双氧水，被视为现行蒽醌法合成技术的绿色替代路线，具有重要战略意义和应用价值，受到学术界和工业界的广泛关注。人们开发了多种高活性双金属催化剂，但受制于双氧水反应强放热特性，该项技术研发陷入了“活性越高越易爆炸”怪圈。分析表明，采用纳微反应器并选择合适的溶剂，可降低催化位点附近反应底物浓度并调控反应－传递匹配程度，从而有效规避爆炸风险，并稳定提高反应的选择性和转化率。

表1 溴丁烷和氰化钠的取代反应C4H9Br + NaCNC4H9CN + NaBr在不同溶剂中的相对反应速率［10］Table 1 Relative rate of SN2 displacement of 1-bromobutane by azide C4H9Br + NaCN C4H9CN +NaBr in various solvents［10］

表1 溴丁烷和氰化钠的取代反应C4H9Br + NaCNC4H9CN + NaBr在不同溶剂中的相对反应速率［10］Table 1 Relative rate of SN2 displacement of 1-bromobutane by azide C4H9Br + NaCN C4H9CN +NaBr in various solvents［10］

溶剂methanol water dimethyl sulfoxide N,N−dimethylformamide acetonitrile介电常数32.6 78.5 48.9 36.7 37.5相对反应速率1 7 1300 2800 5000

迄今为止，国内外研究人员对反应溶剂的影响和调控主要基于经验或试错法，对化学反应的溶剂效应严重缺乏机理上的理解。困难源自于实验和理论两个方面。一方面，实验上可通过红外光谱、拉曼光谱、核磁共振、高分辨粉末X 射线衍射（HRPD）等高精度分析手段检测反应体系中某些关键位置的动态原位信息[11−15]，但这些信息往往都是局部的，管中窥豹式的高精度观测能为机理性理解提供实验佐证，但难以给出溶剂在热力学和动力学上对反应影响的完整内在规律。另一方面，液相化学反应中往往涉及电子的转移和溶剂的结构变化，是典型的介尺度问题。在以往绝大部分研究中，电子、分子等不同尺度上的物理化学过程都是单独处理或解耦研究的，很少涉及不同尺度之间的耦合与关联，耦合量子力学和统计力学的介尺度模型缺失。实验和理论上的困难，严重阻滞了复杂化学反应过程中溶剂效应的深入研究，也不利于绿色化工新技术的发展。

反应活化能和反应自由能计算是研究反应机理最重要的途径之一。从物理化学的角度来看，反应自由能为判断反应路径提供了重要的判据，活化自由能的大小则决定了反应发生的概率。本文从溶剂对化学反应自由能分布影响出发，综述了目前发展的理论模型，并聚焦于本课题组发展的反应密度 泛 函 理 论（reaction density functional theory，RxDFT），分别介绍了RxDFT 的构建及其在水相、有机相、界面体系和限域体系中的应用，分析了不同反应环境对化学反应自由能分布的影响，总结了溶剂效应的影响机制，最后展望了自洽反应密度泛函理论的构建、反应−扩散耦合研究、聚合物反应密度泛函理论及反应密度泛函理论在反应溶剂筛选、界面反应和电解液设计中的应用。

1 溶剂对反应自由能分布影响的研究进展

发展对反应活化能与反应自由能的准确、快速计算一直是计算化学中重要的研究课题之一，到目前为止，人们也发展了多种研究方法[6−8]。总地来说，这些方法可以概括为三大类。第一大类是采用溶 剂 效 应 的 零 阶 近 似（zeroth−order approximation）[7−8]，即将溶剂视为连续介质，在计算过程中通过介电常数的参数设置来体现溶剂的性质。目前，大部分的量子化学计算都是采用这种方法[16−17]。这种方法相对简单且计算成本低，而且对于部分反应，也可以做出定性的解释。但事实上，在很多化学反应中，特别是溶剂为极性或带电溶剂，且反应物与产物分子有极性转变时，考虑溶剂效应后反应自由能有很大的甚至有定性的区别[6,18]。这种区别对于研究反应路径与反应机理极为关键，但难以通过连续介质模型来描述。

第二大类方法是采用经典力学与量子力学耦合的办法，即最早由Car 等[19]提出的量子力学/分子力学（quantum mechanics/molecular mechanics，QM/MM）模拟方法，把溶剂的影响严格地考虑到反应计算过程中[20]。QM/MM 方法原则上可以对反应动力学包括反应路径等进行严格描述。其核心思想是把反应体系分为两个区域，分而治之并合理考虑两个区域之间的耦合作用[19−20]。核心区域发生本征反应，该区域（量子系统）一般只包含参与反应的分子、离子（有时也包括部分催化剂）。外面的区域是溶剂系统，不参与到反应中来，但是与量子系统有相互作用，从而影响反应的发生及反应路径的选择。对于尚未发生反应的反应底物，也包含在溶剂系统中。溶剂系统可以用分子力学来描述，其对量子系统的影响体现在交叉Hamilton 算子上[20]。对于量子系统来说，Hamilton算子由两部分组成：

式中，EQM和EQM/MM是对应的能量本征值。很显然，如果忽略ĤQM/MM（或EQM/MM），式(3)就回归到反应介质效应的零阶近似。

QM/MM 计算较为复杂且计算资源需求非常大[21]，其应用受到了较多的限制，而零阶近似可能有定性的误差。因此人们提出了第三类研究方法，即把式(2)中的耦合项ĤQM/MM用更高级的统计力学理论来描述[20,22−25]。到目前为止，发展了一系列研究方法，其中包括多构型自洽场理论（MCSCF）[26−27]、溶剂平均静电势理论（ASEP）[28−29]、平均场近似（MFA）[30]和RISM/SCF[31]等。在这些方法中，比较引人注意的是RISM/SCF 理论。RISM（reference interacting site−site model）本身是现代经典统计力学中一个比较高级的理论，最初由美国加州大学伯克利分校Chandler 等[32]提出，后经日本Hirata 等[33−35]发展壮大的基于分子水平的积分方程方法。它可以通过自洽的方法快速预测出溶质周围的溶剂结构，继而能预测出溶质的溶剂化自由能（solvation free energy，SFE）[36−37]。RISM 与量化计算用自洽的方法结合就发展成为RISM/SCF 理论，最初由Ten−no 等[38−39]提出，到目前为止已经应用到很多反应系统，如研究水的自电离现象[40]、甲酰胺的互变异构现象[41]等。

值得说明的是，上述这些与量子计算耦合的统计力学理论中，自变量是溶剂系统内组元的密度分布，而不是每个溶剂分子的坐标和动量。因此与QM/MM方法相比，这种与统计力学结合的理论有先天的计算效率上的优势；同时，它又基于较为严格的理论基础，比起溶剂效应的零阶近似来说更为合理，适用范围也更广泛。在具体计算的过程中，溶剂的结构往往与溶质（反应体系）耦合在一起，需要通过多次迭代后，两者才会达到最后的自洽[42]。一旦溶质的构型因为反应发生变化，溶剂的微观结构也会跟着弛豫。这个过程虽然涉及到动力学，但RISM/SCF理论还不能用来描述系统的动力学性质。这可能一方面是因为非平衡态统计理论还在发展之中，另外一方面RISM 理论作为一种积分方程方法，并不适合用来描述诸如扩散等动力学性质。因此，在反应和扩散耦合比较紧密的情况下，需要回归到QM/MM 方法，或者用更好的统计理论来代替RISM与量子力学结合。在这个方向上，经典密度泛函理论是一个很好的选择。

经典密度泛函理论与RISM 理论一样，也是基于分子水平的现代统计力学理论，它是密度泛函理论的一个分支[43]。与量子密度泛函理论（QDFT）依赖于电子密度分布作为泛函变量不同，经典密度泛函理论（CDFT）以系统组分的空间密度分布作为泛函变量，可以预测不涉及化学反应的流体系统的结构和热力学性质，它的处理结果可以与宏观热力学直接联系起来[44]。基于系统组分的几何复杂度，已经发展了三种不同类型CDFT：原子DFT、分子DFT和聚合物DFT[43]。顾名思义，原子DFT 可用于研究由球形粒子组成的简单流体系统[30]。分子DFT 主要是考虑非球形刚性分子的取向贡献来研究分子流体系统[45]。聚合物DFT 则考虑了聚合物分子的链接贡献[46−47]。尽管这三种不同类型的CDFT 的泛函形式明显不同，但它们有着相同的变分框架[43]，即系统自由能或者巨势都可以表达为系统组分密度的泛函。基于此，本课题组将不同尺度下的统计密度泛函和量子密度泛函统一在同一个变分框架下（图2），提出了多尺度平衡态密度泛函理论的统一框架。该统一框架为界面体系多尺度研究提供了高效模型。

近几十年来，CDFT 已发展成为化学、化工及相关领域中最强大最多面的计算工具之一，并在计算多孔介质中气体吸附、受限空间流体等高度非均相流体的结构和热力学性质等方面展现了快速、精确的计算能力[48−49]。其中，分子CDFT 被用于准确预测溶剂化自由能，最近又被用于快速筛选水合自由能[50−51]。近些年，平衡态密度泛函又推广到非平衡态，在扩散、吸附、输运等动力学性质描述方面与动力学模拟相比表现了很好的一致性[52−53]。原子DFT 的动态版本，通常称为动态DFT，在高性能超级电容器的计算辅助设计中得到了广泛的发展和应用[54−57]。

Petrosyan 等[58]首 次 提 出 将CDFT 和QDFT 相 结合来考察溶剂对化学反应的影响，且发展了Joint DFT（JDFT）。在JDFT 框架中，自由能泛函可以写成：

式中，AHK是Hohenberg−Kohn(HK)定理中的电子能量泛函，ACDFT是液体的自由能泛函，ΔA[n,{Nα},V(r)]是溶质与溶剂相互作用有关的耦合自由能泛函。JDFT 通常用极化连续介质模型来处理溶剂，溶剂密度直接从溶质电子密度计算出来，而不是一个独立的变量[59−60]。溶质和溶剂之间的联系采用所谓的电子密度泛函方法来描述[60−61]，从而使电子密度成为JDFT 中唯一的变量[61]，但预测的氢气和氧气密度分布和实验结果有明显的偏差[59]。

本课题组自2014 年以来探索QDFT 与CDFT 结合[43]，发展了反应密度泛函理论，并研究了几类水相、有机相和限域流体中化学反应的溶剂效应，阐述了溶剂对反应路径影响的微观机理，从而为良性反应溶剂选择和溶剂化效应机理研究提供了可行模型。

2 反应密度泛函理论

2.1 反应密度泛函理论框架

图2 密度泛函理论统一框架示意图[43]Fig.2 Schematic graph for the unified framework of density functional theories[43]

在反应密度泛函理论框架下，溶液中的反应自由能FR描述了整个反应的自由能变化，决定了反应能垒高低和体系自由能的变化，它由两部分构成，即反应分子发生构型变化时的本征能量，以及这种本征反应导致的周围液相环境微观结构变化从而引起的液相介质的自由能变化[62−65]。理论上可表示为[6,66]：

式中，ER是本征自由能，也就是忽略溶剂环境发生反应所需要的能量，可采用量子密度泛函计算。ΔFSol表示化学反应发生前后反应体系的溶剂化自由能之差，也就是反应前后溶剂环境的自由能变化，可采用分子密度泛函理论来计算。如图3 所示，考虑溶剂中的一个典型简单化学反应A +，则ΔFSol=。当ΔFSol相对较小时，式(5)回归到溶剂效应的零阶近似。

图3 典型化学反应在气相和溶液中的热力学循环（Fsol为反应物/产物分子的溶剂化自由能，ER表示在气相中的反应自由能，FR表示在液相中的反应自由能）[66]Fig.3 Thermodynamic circle for a representative chemical reaction()occurs in gas phase and in solution(Fsol is the solvation free energies of reagent/product molecules，ER denotes the reaction free energy in gas phase,and FR is the reaction free energy in solution)[66]

2.2 反应自由能和活化自由能

根据过渡态理论（transition state theory，TST）[67]，如图4（a）中所示，气相中量子系统在反应前后处于平衡状态，其自由能分别设为F0和F2。两者之差即为反应自由能，即ER= F2−F0。假设过渡态系统自由能是F1，则反应自由能垒或活化自由能为EA=F1−F0。相似地，在溶液中反应自由能可以表示为FR= F′2−F0。如图4（b）所示，反应自由能垒或活化自由能为FA= F′1−F0。由式(5)可知，由于本征反应与溶剂效应的解耦，FR和FA可以通过将相应的本征反应自由能ER和EA与溶剂化自由能差相加得到。

图4 在气相中（a）和在溶剂中（b）的反应的典型自由能分布（过渡态与反应物之间的能量差为反应的自由能垒或活化自由能，其代表了一个反应系统发生所必须获得的最小能量。生成物和反应物的能量之差为反应的自由能）[66]Fig.4 Typical energy profiles for a reaction in gas phase(a)and in a solvent(b)(The difference between the energies of transition state and reactant is the energy barrier or activation free energy of reaction;it represents the minimum energy that a reacting system must acquire for the transformation to take place.The difference between the energies of the product and the reactant is the free energy of reaction)[66]

在反应密度泛函理论框架下，可假设溶剂的弛豫通常比反应过程快得多[9,68−69]，即溶剂化自由能的计算可以与本征反应自由能的计算解耦。因此，溶剂化自由能可以表示为溶质周围溶剂密度分布函数的泛函：

式中，ρ(x)是溶质分子的密度分布，x=(r,Ω)表示溶剂分子的坐标和空间取向。ρ0和μ 分别是溶剂的体相密度和化学势，V 和T 分别是系统体积和热力学温度。体相化学式由体相密度通过状态方程确定，系统的体积选取需足够大，其数值并不影响计算结果。在下面的泛函变量符号中省略μ、V 和T。Vext(x)表示溶剂系统的外势，它来源于溶剂和溶质的相互作用。随着化学反应的发生，量子体系的分子结构和电荷分布会发生变化，Vext(x)也会随之发生变化。Vext(x)是由量子系统和溶剂系统之间的分子原子相互作用造成的，和其他多尺度方法一样，它通常包括库仑相互作用和Lennard−Jones（L−J）相互作用[31,58]，其表达式为：

对于给定的Vext(x)，溶剂系统的巨势可以由式（8）计算得到：

式中，kB是Boltzmann 常数，Λ 是溶剂分子的有效热波长，Fexc[ ρ(x)]是由溶剂分子间相互作用引起的过剩自由能泛函。到目前为止，已经开发了密度展开的方法和自由能分解法等来构建过剩自由能泛函。

在热力学平衡条件下，溶剂密度分布使溶剂体系的自由能泛函达到最小。换言之，通过对式(6)中溶剂化自由能泛函最小化，可以得到平衡时的溶剂密度分布和溶质的溶剂化自由能[50,71]，即：

因为CDFT 直接将自由能与系统组分的密度分布关联起来，从而绕过了冗长的热力学积分，大大减少了计算成本。这点与分子模拟相比具有明显的优势[71]。利用CDFT 可定量预测不同中性溶质和离子的水合自由能[50]。另一方面，与RISM 理论相比，CDFT 的理论框架更加灵活，而且它可拓展到非平衡系统用来描述扩散动力学[54−57]。

3 反应密度泛函理论初步应用

3.1 水溶液对化学反应的影响

3.1.1 甘氨酸水合反应路径研究 互变异构化现象在化学和生物领域普遍存在[72]。在溶液环境中研究互变异构化已引起了越来越多的理论研究兴趣[73−82]。这是因为溶液中的互变异构化反应不仅可作为模型系统来研究，而且这些反应具有实际应用的重要性[83−85]。互变异构化反应也同时可成为生物质能转化中酮化和醛醇缩合反应甚至分子炸药分解的决速步[86−87]。

甘氨酸（glycine）是最简单的氨基酸，白色或浅黄色晶体，易溶于水，有甜味。甘氨酸在气相或真空中是以中性形式（N）存在的；而在水溶液中，由于与水分子的静电相互作用，以两性离子形式（Z）存在。尽管实验测定从中性形式到两性离子形式的自由能变化约为−7.3 kcal/mol[73,88]（1 kal=4.18 kJ），但甘氨酸在水溶液中的结构转变路径有两种可能，分别对应两种反应机制：非辅助协同机制（nonassisted concerted mechanism）和辅助协同机制（assisted concerted mechanism）[82]，两种反应机制的中间过渡态分别为TS1和TS2。如图5所示，在非协同反应机制中，中性结构中的氧原子上的氢原子朝向氮原子方向，在过渡态结构（TS1）中氢原子与氮原子形成了很弱的键，由两个碳原子、一个氧原子、一个氮原子和一个氢原子组成一个五元环的结构，之后氢原子与氧原子之间的键发生了断裂，形成了两性离子结构，完成了质子转移。在协同反应机制中，水分子作为传递氢原子的媒介参加了反应。在中性结构中，氧原子上的氢原子朝向水分子中氧原子，同时氧原子中的一个氢原子朝向氨基乙酸中的氮原子，在过渡态结构（TS2）中由于氧原子上的氢原子与水分子中的氧原子，水分子中的氢原子和氨基乙酸中的氮原子成键，由两个碳原子、两个氧原子、一个氮原子和两个氢原子形成了一个七元环结构，之后甘氨酸中的氧原子与氢原子断键，水分子中的氧原子与连接在氮原子上的氢原子发生断键，形成了两性离子结构，完成了质子转移[82]。

图6是采用反应密度泛函理论预测甘氨酸互变异构化反应分别在气相、水溶液环境中的反应路径。在MP2、B3LYP、M06−2X 水平上计算本征自由能，得到的自由能分布表明在气相中的化学反应是吸热的。RxDFT 预测的非辅助协同机制和辅助协同机制的自由能分布图表明甘氨酸异构化反应在水溶液中是放热的。RxDFT 理论预测的结果与文献报道采用MD/MP2方法的计算结果和实验测量结果吻合度较好[82]，这表明RxDFT 方法的准确性和可靠性。

图5 非辅助协同机制（a）和辅助协同机制（b）中甘氨酸中性、过渡态和两性离子结构的示意图[66]Fig.5 Schematic representations and atomic numbering for the neutral,transition,and zwitterionic structure of glycine with nonassisted concerted mechanism(a)and assisted concerted mechanism(b)[66]

图6 非辅助协同机制（a）和辅助协同机制（b）的甘氨酸自由能分布图[66]（G表示气相反应，S表示水溶液反应。MD/MP2方法的自由能分布来自文献[82]）Fig.6 Free energy profiles for glycine using a nonassisted mechanism (a)and an assisted mechanism(b)[66](The label G denotes the reaction in gas phase and S denotes the reaction in aqueous solution.The free energy profiles of MD/MP2 method are taken from Ref.[82])

图7 四个互变异构体（A、B、C、D）与相应的过渡态（TS1、TS2、TS3、TS4）之间的反应方案（AOO的分子结构是在没有水分子参与的情况下展现的）[96]Fig.7 Reaction scheme among the four tautomers(A,B,C,and D)with the transition states(TS1,TS2,TS3,and TS4)(The molecular structures of AOO are depicted in the absence of water)[96]

对于简单的互变异构化反应，气相和水溶液中几何优化和能量分析表明，异构体A 和C 分别是最稳定和最不稳定的。四个异构体之间的稳定性顺序是A > B > D > C。如图8 所示，溶剂效应使所有简单互变异构化反应的自由能垒都增加，使B →C和C →D 异构化反应自由能增加，使A →B 和C →D 反应自由能降低。如图9 所示，对于水协助的互变异构化反应，水分子的协助降低了B →C 和C→D 异构化反应的反应自由能，增加了A →B 和D →A 异构化反应的反应自由能。在水分子的协助下，所有反应路径中的自由能垒都大大减少。而且这些趋势与水分子的数量无关。水分子的协助作用可以归结为两个相反的效应，即过渡态环的大小和氢位移步效应，且二者竞争协调，当两个水分子参与协助时，能垒最低。其次，分析了溶剂对水分子协助互变异构化反应的影响，发现溶剂可以进一步降低反应的自由能垒，虽然不是所有反应的自由能都降低。一般来说，水溶液的存在主要是通过水分子的协助作用和溶剂化效应来降低自由能垒。溶剂化效应对自由能垒几乎没有贡献，溶剂分子的协助作用对互变异构化反应的溶剂效应起主导作用。

图8 气相和水溶液中非辅助协同机制互变异构化反应的自由能分布（青色曲线与绿色曲线互相重叠。温度298.15 K）[96]Fig.8 Free energy profiles for the nonassisted mechanism tautomerization reactions in the gas phase and aqueous solution(It should be noted that the cyan curve overlaps and the greencurve.The temperature is 298.15 K)[96]

3.1.3 水溶液对SN2 反应的影响 双分子亲核取代反应是构建官能团最常见、最有价值的反应之一，因此在化学和生物化学中具有重要意义[97−98]。在X−++Y−反应中，根据亲核试剂X−和离去基团Y−的不同，可以将SN2 反应分为两大类：对称（Type Ⅰ）和非对称（Type Ⅱ）SN2 反应。顾名思义，对称SN2 反应中的亲核试剂和离去基团是一样的，而非对称SN2 反应中亲核试剂与离去基团不同[99]。SN2反应的反应机理最早由Ingold[100]在1953年提出，随后Brauman 等[101−102]提出了SN2 反应的详细原子模型。这个模型假设在气相中反应首先形成X−…CH3Y 离子−偶极复合物（ion−dipole complex，IDC），然后X−…CH3Y离子−偶极复合物越过能垒形成XCH3…Y−离子−偶极复合物或者分解成反应物。

图9 水溶液中0～3个水分子参与的互变异构化反应的自由能分布（温度为298.15 K）[96]Fig.9 Free energy profiles for the assisted mechanism tautomerization reaction in the absence and presence of 1—3 water molecules in aqueous solution(The temperature is 298.15 K)[96]

自由能分布是研究反应机理和揭示反应动力学的重要性质之一[97,103−104]。由于SN2 反应对周围溶剂的存在很敏感，所以这些反应在气相和在溶液中的自由能分布有很大的差别[105]。如图10所示，在气相中SN2 反应通常沿反应坐标具有双阱自由能分布。自由能分布图有一个单一的过渡状态（TS）和两个能量最小的稳定态，这两个稳定态对应于两个能量稳定的IDC[106−107]。然而，在溶液中，这两个能量最小的稳定态完全消失不见了，自由能分布图的两端是平整的，并且在反应物和生成物之间存在一个巨大的自由能垒[107]。此外，据报道，在质子溶剂中SN2 反应的速率常数比在非质子溶剂中要小几个数量级（表1）[10,108]。

图10 碘离子与碘化甲烷反应的中间产物的微水化效应[107]Fig.10 Microhydration effects on the intermediates of the SN2 reaction of iodide anion with methyl iodide[107]

在气相和在溶液中调查对称和非对称SN2 反应引起了广泛的实验和理论研究兴趣[99,107]。Chandrasekhar 等[105]证明了水溶液对典型的对称SN2反应Cl−+ CH3ClClCH3+ Cl−的自由能分布有显著影响。Doi 等[107]采用红外光解光谱（IR photodissociation，IRPD）阐 明 了 微 溶 剂 化 对I−+CH3IICH3+ I−SN2 反应的作用。对于非对称SN2 反 应，Gao 等[109−110]采 用QM/MM 方 法 得 到 了NH3+ CH3ClNH3CH+3+ Cl−SN2 反应的二维自由能面。他们的研究结果表明，溶剂效应对产物和过渡态具有很强的稳定性作用。这是由于在反应过程中发生了大量的电荷分离。此外，还有很多实验采用流动余辉（flowing afterglow）技术研究了溶剂对SN2 不对称反应的影响[111−116]。尽管化学反应在溶液中的实验技术和预测模型取得了令人满意的进展，但目前溶剂效应对化学反应影响的潜在分子机制仍不清楚。这个挑战主要来自于反应子系统与溶剂环境之间复杂和动态的相互作用。

如图11 所示，采用RxDFT 考察了对称SN2 反应Cl−+ CH3ClClCH3+ Cl−和 非 对 称SN2 反 应NH3+ CH3ClNH3CH3++ Cl−的自由能分布。对于对称SN2反应，量子DFT在M06−2X理论水平上获得的气相的本征自由能分布图具有一个稳定的IDC 的双阱自由能分布。在水溶液中，活化自由能增大，而双阱的能量分布则由于TS 和IDC 在溶剂中的不稳定而消失。对于非对称SN2 反应，它在气相中是吸热反应。然而，它在水溶液中变成了放热反应，伴随着活化自由能的减少以及过渡态位置的位移。对于这些SN2 反应，RxDFT 预测的活化自由能和反应自由能与QM/MM和实验数据的结果一致。

3.2 乙腈溶液对SN2反应的影响

图12 和图13 为采用RxDFT 方法预测的八个对称SN2 反应和四个不对称SN2 反应在乙腈溶液中的自由能分布。如图12 所示，在所选的八个对称SN2反应中，不同的烷基氯化物在气相中的活化自由能差异很大。每增加一个甲基的取代都会增加活化自由能，这些都是空间位阻效应引起的。在乙腈溶液中，极化效应降低了约6 kcal/mol 的自由能垒，溶剂化效应提高了约18 kcal/mol 的自由能垒。因此，溶剂化效应对自由能垒起主导作用。如图13所示，所选择的四种不对称SN2 反应，由于环境温度下它们具有较高的自由能垒和吸热的反应自由能，因此在没有催化剂或溶剂辅助下，它们显然是不容易发生的。在乙腈溶液中，由于溶剂作用使过渡态和产物更稳定，活化能垒显著降低，反应变为放热反应。极化效应和溶剂化效应都使自由能垒和反应自由能降低，而且溶剂化效应起主导作用。对于这些SN2 反应，可以看出RxDFT 预测的活化自由能和反应自由能与用直接与热力学循环方法结合SMD−M062X 溶剂化模型、Monte Carlo 模拟和实验得到的结果比较吻合[119−121]。

3.3 固液界面对亲核加成反应自由能分布的影响

固液界面处发生的化学反应往往会展现出不同于均相化学反应的性质。为了探究固液界面处的复杂溶剂效应，本课题组将RxDFT 方法拓展到氢氧根与甲醛在类石墨烯基底附近反应的研究中。通过RxDFT 方法计算了反应物、过渡态和产物在距离固体基底不同距离时的溶剂化自由能，发现溶质距离固体基底越近溶剂化自由能越小；计算了溶质在距离固体基底不同距离时周围的溶剂密度分布，发现当距离较小时基底和溶质附近水的波动发生重叠，这种重叠关系导致了溶剂化自由能的降低；计算了溶质在距离固体基底不同距离时的反应自由能和反应活化能，发现当反应核心区域贴近基底时，反应活化能相较于体相环境中的能量值有所降低（图14）。这意味着该反应更容易发生在石墨烯−水的界面附近[70]。

图11 对称SN2反应Cl−+ CH3Cl ClCH3 + Cl（−a）和非对称SN2反应NH3 + CH3Cl NH3CH3+ + Cl（−b）的自由能分布(1 Å=0.1 nm)[117]Fig.11 Free energy profile for the symmetric Cl−+ CH3Cl ClCH3 + Cl−(a)and asymmetric NH3 + CH3Cl NH3CH3+ + Cl−(b)SN2 reaction in aqueous solution or in the gas phase[117]

图12 八个对称SN2反应的相对自由能分布（温度为298.15 K）[118]Fig.12 Relative free energy profiles for the eight selected symmetric SN2 reactions in gas phase and in the acetonitrile solution(The temperature is 298.15 K)[118]

图13 四个非对称SN2反应的相对自由能分布（温度为298.15 K）[118]Fig.13 Relative free energy profiles for the four selected asymmetric SN2 reactions in the gas phase and acetonitrile solution(The temperature is 298.15 K)[118]

图14 反应自由能（a）和活化自由能（b）随团簇与界面之间距离的变化[70]Fig.14 Reaction free energy(a)and activation free energy(b)as a function of the cluster−wall distance[70]

3.4 表面润湿对受限空间流体溶剂化效应的影响及其应用

当反应发生在受限空间溶剂环境中时，在受限空间流体溶剂化效应的影响下，会呈现与均相溶剂中不一致的反应路径或反应机理[122−123]。烷烃分子构象转变是一种特殊的反应。据报道，烷烃分子在不同的溶液体系中可呈现出直链（extended）−螺旋（helical）构象转变的现象[124−125]，这种构象的转变可以在实际中有多种用途，比如选择性分子识别[124,126]以及分子自组装等[127]，研究烷烃分子构象转变的微观机理以及表面润湿调控机制对实际应用中对分子构象精准调控有重要指导意义。为了给出烷烃由直链构象转变为螺旋构象的微观机理，采用RxDFT 方法的框架研究不同孔宽度、壁面润湿性等不同受限效应对不同链长的烷烃分子构象转变的影响。如图15所示，烷烃分子在孔径较小的纳米圆管水溶液中更易形成螺旋构象。圆管管壁润湿性越强，即管壁越亲水，烷烃分子越倾向于由直链构象转变为螺旋构象，且烷烃分子链越长，这种转变倾向越明显。当孔径较小时，螺旋烷烃与直链烷烃系统的巨势之差随孔径变化而“振荡”变化，即烷烃分子构象随孔径变化呈现“直链−螺旋”的交替变化现象，随着孔径增大，烷烃分子逐步恢复稳定直链构象。本节中所做的理论预测结果能够与前人的实验结果很好地吻合。

4 未来展望

4.1 自洽反应密度泛函理论的构建

图15 单个烷烃分子在充满水的空腔中的构象变化（a）；充满水的圆柱孔中庚烷的巨势差随孔内表面润湿参数的变化曲线（b）；在充满水的圆柱孔中的巨势差随烷烃分子碳链长度NC的变化曲线（c）；充满水的圆柱孔中庚烷的巨势差随孔径大小的变化曲线（d）[128]Fig.15 Conformational change of a single alkane molecule confined in water−filled cavitands (a);Grand potential difference of heptane in water−filled cylindrical pore varying with the wetting parameter of the pore inner surface (b);Grand potential difference of alkane molecule with varying NC in water−filled cylindrical pore(c);Grand potential difference of heptane in water−filled cylindrical pore varying with the pore size(d)[128]

目前，RxDFT框架下的QDFT和CDFT的耦合模式仅限QDFT 对CDFT 的单向参数传递。在RxDFT具体的计算过程中，仅将QDFT 方法计算得到反应体系的结构参数、各原子的偏电荷结合各原子的L−J 参数，作为外势参数输入CDFT 进而计算得到反应体系周围溶剂分子的密度分布ρ(x)和溶剂化自由能Fsol。此时的溶剂分子密度分布相对于初始状态已经发生了变化，该变化并没有反馈到QDFT 的计算中，导致溶质分子的波函数没有得到溶剂的充分极化，CDFT与QDFT在RxDFT理论框架下并没有达到最终的自洽。如果溶剂分子为较小刚性分子，溶剂分子产生的电场对溶剂分子的密度分布影响并不大，此时没有QDFT与CDFT间的迭代耦合也可以实现相当不错的计算精度。当溶质分子为较大的柔性分子时，溶剂分子产生的电场对溶剂分子密度分布的影响就不能忽略，此时QDFT和CDFT间的迭代自洽对于RxDFT 的计算精度就显得尤为重要。在下一阶段，拟将基于RxDFT 方法发展为一种自洽反应密度泛函理论(self−consistent reaction density functional theory，sc−RxDFT)方法。

要实现sc−RxDFT 框架下的QDFT 和CDFT的迭代耦合，需要实现两方面的自洽。首先是溶剂对溶质分子轨道充分的极化。在QDFT 计算中，溶质分子的溶剂化可以通过在其Fock 算符后加上一个溶剂化修正项来实现，该修正项可表示为溶剂密度分布的函数，溶剂密度分布则通过CDFT 求解得出；另一方面是溶质分子形成的电场对溶剂分子密度分布产生影响。在CDFT计算中，把QDFT计算得到的溶质分子的原子偏电荷及原子坐标结合原子的L−J参数传递给CDFT，进而计算出溶剂化自由能和溶剂新的密度分布。然后再将新的溶剂密度分布传递到QDFT 中进行新一轮的计算，这样经过多次迭代，直至溶质分子自由能和溶剂化自由能实现自洽。具体计算流程请参考图16。sc−RxDFT 方法的建立将为RxDFT 在化学、化工诸多领域的应用打下坚实的基础。

4.2 聚合物反应密度泛函理论

经典催化系统往往具有循环利用率低、毒性高和催化效率难以调控等缺点。在催化剂表面铆接高分子可以克服这些缺点，达到定向调控催化活性的目的。铆接高分子一般对控制条件（如铆接密度、高分子链长、温度、压强和pH 等）的变化具有高度可逆的响应。通过改变控制条件可以精确控制高分子的体积、形状、表面积、溶解度和力学性能等性质。这些高度响应的性质可以用来改性催化剂表面，达到调控催化反应的目的。目前人们对表界面性质与纳微结构的关系、界面作用机理和反应动力学过程等仍然缺乏足够的认识，研究的难点主要在于原位实验观测手段（特别是微观或介观的手段）的缺乏或不够精确，而常规的实验试错法难以获得高分子表面改性与反应效率之间的完整构效关系。聚合物密度泛函理论适用于聚合物体系，衔接了微观结构性质和宏观热力学性质，是研究聚合物表界面调控的重要理论工具。其以聚合物的空间密度分布为变量，把自由能和系统组元密度分布直接联系起来，通过求泛函极值，可同时获得宏观热力学和系统组元的空间密度泛函。然而，相对于平衡态体系，描述非平衡态条件下化学反应的聚合物反应密度泛函理论亟待完善和发展。这一方面是因为聚合物具有复杂的链构型分布，另一方面是由于聚合物体系具有多尺度性，从理论上描述聚合物体系的化学反应具有一定的难度。从分子尺度建立聚合物反应密度泛函理论，研究催化反应中聚合物改性对表面催化反应速率的影响机制，考虑铆接聚合物与反应物、产物分子之间的相互作用，总结控制条件对反应效率的影响规律，具有重要的理论和应用价值。

4.3 反应-扩散耦合研究

图16 RxDFT和sc−RxDFT的计算流程图Fig.16 The scheme of RxDFT and sc−RxDFT

反应−扩散（reaction−diffusion，RD）过程在化工、催化、材料加工和生物领域普遍存在[129]，是化工过程、分离、催化反应等动力学过程中最重要和最基本的过程之一，在化学工程的核心“三传一反”中占据重要地位，对化工过程中反应分子的微观传递行为及催化剂内部发生的化学反应的深入研究有助于化工反应器和催化剂及化工过程的设计。成功设计新的催化剂和反应器需要掌握催化剂反应微孔空间内分子的催化反应与分子扩散的耦合效应，可以说，由限定在催化剂表界面及孔道内分子发生的化学反应及与扩散协同构成的复杂过程，涉及“三传一反”理论框架在具有非均匀结构体系中的拓展应用，是当前化工学科中最富挑战性的研究方向之一。然而，目前反应−扩散耦合主要是采用宏观扩散理论与反应理论来耦合研究。这种耦合方法虽能解决一些体相的反应问题，但对于界面反应，理论预测与实验结果相差较大[130]。根本原因在于这些宏观理论对于界面体系（特别是组元密度分布高度非均相的界面体系）的性质描述不够准确。在界面区域，流体密度高度不均匀，需要采用分子层面的理论才能对反应−传递耦合进行合理描述。而大部分反应过程为界面反应，如界面催化等，界面的属性影响反应与传递的耦合，改变界面性质调控反应−传递耦合是提高化工反应效率的重要途径之一。因此，发展面向表界面体系的反应−扩散耦合理论研究具有十分重要的意义，具有极大的潜在应用价值。

4.4 RxDFT在反应溶剂筛选及界面反应中的应用

催化反应是化学工程的核心组成，是实现物质转化、工业生产的重要步骤。迄今为止研究人员针对气−固反应中催化剂的开发和负载，发展了各种性能优越的催化剂和巧妙的催化颗粒/基团负载方式；而对于液−固反应或气−液−固三相反应，反应的溶剂调控就显得尤为重要。在这些催化反应过程中，溶剂效应既影响反应体系组元分子的空间密度分布，也影响了其反应−传递耦合的匹配程度。具体体现在对催化位点周围反应物局部密度（溶解度）的调节，对反应物、产物分子扩散动力学特性的改变，以及对传热效率的影响等。

基于研究人员的经验积累和实验试错，人们发展出一系列绿色反应溶剂，如超临界二氧化碳、离子液体等，但目前仍缺乏对溶剂效应的深入理解。因此针对纳微界面反应，溶剂调控反应效率目前没有相关理论依据支撑，找到合适溶剂是一个高成本、高耗时的过程。理论计算是这一问题的解决方案之一。然而，目前针对界面反应，常用的软件要么不考虑溶剂的作用要么采用极化连续模型处理溶剂效应。例如，在VASPsol 中，通过向Poisson 方程中添加与电子密度分布相关的介电函数来体现溶剂对化学反应的影响，且需要预先设定大量的合理参数来提高计算的准确性[131]。因此，RxDFT 在界面催化领域具有巨大的潜在应用价值。同时，将RxDFT 并入现有软件有利于其推广和应用，加速反应溶剂的筛选，为高选择性和高转化率的实现和突破提供思路。

4.5 RxDFT在电解液设计中的应用

电解液是电化学电池和超级电容器的重要组成部分。在电解液中，由于离子与溶剂分子的强相互作用，离子总是以溶剂化的形式存在，且溶剂化结构直接影响电解液的稳定性。近年来，溶剂化层中溶剂的作用被不断认识，并成为电解液设计的一种有效手段。理解电解液溶剂化结构及其构效关系，开发设计更加稳定、高效的电解液体系，是实现电化学电池和超级电容器实用化的必然要求。最近，清华大学张强教授团队Zhang 等[132]通过在溶剂化层中引入硝酸根阴离子NO−3，形成了更大的溶剂化团簇，并促进双（氟磺酰基）亚胺阴离子FSI—的分解，形成富含LiF 的界面层，拓宽了电解液的稳定窗口。这些溶剂化过程可以看作离子与溶剂分子之间发生的缔合反应。由于涉及到的溶剂种类是多元的，因此需将RxDFT 拓展到混合溶剂体系，发展面向混合溶剂的RxDFT，为高效电解液体系的设计提供理论指导。

5 总 结

液相反应中的溶剂调控是发展绿色化学的重要内容。溶剂对液相反应的能垒、速率甚至反应机理都有重要影响，选择合适的溶剂已成为调控反应选择性和转化率的控制步骤。目前对溶剂的遴选主要基于经验或试错，能够深入描述溶剂对反应自由能影响的定量模型依然缺乏。本文综述了反应密度泛函理论的构建及其在水相、有机相、界面体系和限域体系中的应用，成功阐述了溶剂对几类化学反应的影响机理。RxDFT 能够定量描述溶剂对反应自由能的影响，从而为液相反应中的溶剂化效应机理研究提供了可行模型，为绿色化学发展中的溶剂定向调控提供了理论依据。