烯丙基正離子旋轉異構反應的計算化學實驗設計

王亞妮，張學鵬

陜西師范大學化學化工學院，西安 710119

隨著計算機技術的不斷突破和量子化學理論的不斷完善，理論計算在材料、催化、合成、生物醫藥等不同領域都得到了廣泛應用[1]。計算化學已受到廣大科研工作者的高度重視，使學生學習并掌握一定計算化學實驗技能也愈發必要。

目前，國內很多高校都已將計算化學引入本科或低年級研究生教學[2-6]。自新冠肺炎疫情爆發以來，為響應國家抗擊疫情的號召，全國各大高校按照教育部提出的“停課不停教、停課不停學”的要求，積極開展了多種形式的線上教學。在化學專業課程的線上教學過程中，一些知識點比較抽象，學生較難快速地理解和掌握，同時教師也較難及時了解學生的學習進度。特別是在化學實驗課的線上教學過程中，學生無法通過線上課堂進行實驗操作，教師無法對儀器設備、藥品進行實物講解，使得學生對實驗現象的了解以及實驗整體的把握都不好。相對地，計算化學可以較好地適應線上教學，實驗安全性高，教學成本低，并且可以幫助學生更好地了解和掌握微觀抽象的化學知識點。在實驗過程中，教師與學生只需通過計算機，就可以共同完成實驗，教師可以根據計算結果，實時掌握學生的實驗進展，對存在的問題進行第一時間的溝通與處理。并且計算化學可以將抽象的反應勢能面、過渡態、分子軌道、原子電荷、靜電勢等概念進行圖形化處理，讓學生可以直觀、深入地掌握這些知識[7]。計算化學實驗的開設不僅可以幫助學生掌握計算化學的基本概念與實驗技能，還可以與無機化學、物理化學、有機化學等專業課程的理論知識結合，幫助學生加深對一些復雜的化學概念的理解[8-12]。

本文以烯丙基正離子的旋轉異構反應為例，設計了一個完整的計算化學實驗，該實驗設計可以較好地進行線下與線上教學。烯丙基正離子的分子結構簡單，計算耗時較少，可以在規定課程時間內完成教學。學生通過該實驗設計能夠進一步加深旋轉異構反應機理的掌握，包含反應物結構的確定、反應勢能面的構建、活化能的計算。我們希望能夠通過此計算化學實驗課程設計，幫助學生掌握基本的計算化學概念與實驗技能，從而對其后面的課程學習以及科研有所幫助。

1 實驗目的

1)了解計算化學基本概念，了解結構優化、單點能計算、過渡態尋找、布局分析等計算化學常用研究手段。

2)了解和掌握烯丙基正離子旋轉異構反應機理。

2 實驗原理

隨著計算機性能的不斷提升和計算化學理論與程序的不斷完善，計算化學已經成為化學、材料、環境等領域不可或缺的研究手段之一。目前，常用的電子結構理論基于量子力學的基本原理和方法，通過求解薛定諤方程得到分子的能量和相關性質。薛定諤方程可以寫為：Ψ = EΨ，式中為哈密頓(Hamiltonian)算符，Ψ為波函數，E為能量。分子體系的算符由以下幾項組成：

基于不同的處理手段，現今常用的電子結構理論有從頭算方法(ab initio)、半經驗方法以及密度泛函理論(Density Functional Theory，DFT)。從頭算方法在求解薛定諤方程時對哈密頓算符不做過多的簡化，不引入經驗參數來擬合雙電子積分[13]。目前，從頭算方法包括基于Hartree-Fock方程的HF方法[14]，以及在HF方法基礎上引入電子相關作用而發展起來的后HF (Post-HF)，主要有多級微擾方法(Many-body Perturbation method，MPn)[15]、耦合簇方法(Coupled Cluster method，CC)[16]、組態相關作用法(Configuration Interaction method，CI)[17-19]等計算方法。HF方法引入了波恩-奧本海默絕熱近似、非相對論近似和單電子近似，考慮了自旋電子之間的交換作用，但忽略了電子間的庫倫相關作用，對于弱相互作用(比如氫鍵、鹵鍵、π-相互作用等)比較明顯的體系描述不好。半經驗方法在求解薛定諤方程時通過省略多中心積分或者直接使用擬合的經驗參數來近似求解自洽場分子軌道方程，極大提高了計算效率。PM6是目前使用比較廣泛的半經驗方法[20]，其支持過渡金屬(支持上限到鉍元素)，在大分子量體系計算、結構初猜與篩選等方面有較廣泛的應用。密度泛函理論使用電子密度函數取代多電子波函數來描述相互作用，極大簡化了薛定諤方程的求解。密度泛函理論[21-23]依然引入了絕熱近似、非相對論近似以及單電子近似，與HF方法相比，DFT在單電子近似的基礎上進一步考慮了交換能與關聯能，對電子相關作用描述得更好一些，但依然不能較好地描述弱相互作用。DFT方法的計算速度比HF方法快或者相當，但計算精度大大增加。本文中使用的B3LYP泛函[24]是最常使用的DFT泛函之一，經常用于體系的結構優化以及能量計算。

量化計算中用基組來構建體系的波函數，常用的基組由多個高斯型函數(GTO)組成。一般基組越大，計算的精度也越高。例如6-311++G(d,p)，其表示用6個GTO來構建一個內層軌道，外層價軌道劈裂為三組，分別用3個、1個和1個GTO來構建，“d”表示給非氫原子添加d軌道，“p”表示給氫原子添加p軌道，第一個和第二個“+”分別表示給非氫原子和氫原子添加彌散函數。通過給軌道添加d或p角動量可以改變軌道形狀，從而易與其他軌道成鍵。另一方面，通過給軌道添加彌散函數，可以使軌道占據更大的空間，從而有利于弱相互作用的描述。

Gaussian是目前使用最廣泛的計算化學軟件之一，它功能強大，簡單易學，支持多種計算方法。GaussView是輔助Gaussian程序使用的圖形界面軟件，主要用于分子建模和計算結果分析。考慮到Gaussian和GaussView為收費軟件，本實驗也可以通過Chem3D、HyperChem等免費軟件進行分子建模，使用開源軟件ORCA[25]、GAMESS-US[26]、Dalton[27]等完成量化計算。

烯丙基正離子常在SN1親核取代反應中產生，是有機反應重要的活潑中間體，其可以發生氧化、偶聯、親電取代等多類化學反應。著名的Friedel-Crafts烷基化反應[28]是典型的芳烴烯丙基化的方法。Suji-Trost烯丙基取代反應[29]可以合成出多種用于反應的過渡金屬和親核試劑。深入了解烯丙基正離子的反應活性和化學性質有利于幫助我們更好地掌握相關有機反應的機理，進一步設計合成更多天然化合物。目前，烯丙基正離子旋轉異構反應的機理一般認為通過以下三個步驟進行：① 烯丙基正離子中間碳原子上的氫原子遷移到α位碳原子上，形成甲基；② 發生甲基與烯烴之間C―C鍵旋轉；③ 氫原子從α位碳原子遷移回中間碳原子上。

根據分子結構的大小，以及計算耗時的考慮，本實驗主要使用DFT的B3LYP泛函和6-31G(d,p)基組對反應物結構進行優化和頻率計算，并在此基礎上使用更高精度的B3LYP/6-311+G(d,p)方法計算反應物的熱力學能。另外，我們還將B3LYP的結果與HF、PM6計算的結果進行了對比，從而讓學生更加直觀地了解計算化學方法的差異。需要注意的是，本實驗設計的教學內容較多，常規4小時的實驗時間可能無法完成全部的實驗內容，指導教師可以基于學情酌情選擇適當的教學內容。

3 實驗設備

臺式機或者筆記本電腦(安裝Windows XP，Vista 8, 10或11)，Gaussian 09[30]軟件，GaussView 5.0或6.0軟件，Origin軟件，Chemdraw軟件。

4 實驗步驟

4.1 構建分子初始結構

使用GaussView軟件構建反應物(reactant，RC)、反應過渡態(transition state，TS1、TS2和TS3)、反應中間體(intermediate，IM1和IM2)以及產物(product，PC)的初始分子結構，其結構如補充材料中圖S1所示。具體操作如下：點擊Create New Molecule Group按鈕建立工作區，然后點擊Element Fragment按鈕在彈出的元素周期表中選擇需要的元素以及構型，在工作區構建各反應分子。烯烴基團也可以通過點擊R-Group Fragment按鈕在彈出的窗口中直接選取烯烴基團。分子建模提示：由于烯丙基正離子的α位碳原子失去了一對電子，帶有部分正電荷，其p軌道和C＝C雙鍵的π鍵之間具有p-π共軛效應，所以反應物RC應呈平面結構。過渡態TS1涉及氫遷移，具有即將斷裂和即將形成的C―H鍵，所以應有一個氫原子位于中間碳原子和α位碳原子的上方。氫遷移完成后，產生一個中間體IM1。IM1通過甲基旋轉過渡態TS2產生中間體IM2。IM2再發生氫遷移，得到過渡態TS3。最后氫原子重新回到中間碳原子上后，得到仍是平面結構的產物PC。

4.2 編輯輸入文件

使用GaussView完成分子建模后，點擊Calculate → Gaussian Calculation Setup進行計算命令設置，如圖1(a)。

圖1 (a) GaussView軟件的計算命令設置界面；(b) Gaussian 09軟件執行量化計算界面

(1) Job Type部分：其中Energy表示能量計算，Optimization表示結構優化計算，Frequency表示分子頻率計算，Opt + Freq表示分子結構優化后繼續計算分子的頻率，IRC表示計算過渡態分子的內稟反應坐標，Scan表示對分子的鍵長、鍵角、二面角等參數做勢能面掃描，Stability和NMR分別表示穩定性計算和核磁譜圖計算。本實驗中依次選擇Opt + Freq和Energy完成結構優化和熱力學能的計算。本界面其他選項采用默認設置。

(2) Method部分：用來選擇實驗所用的計算方法。當Job Type選擇Opt + Freq時，依次選擇DFT →B3LYP泛函和6-31G(d,p)基組，設置結構優化的命令；當Job Type選擇Energy時，依次選擇DFT →B3LYP泛函和6-311+G(d,p)基組，設置能量計算的命令。方法對比時可以選擇Semi-empirical → PM6或Hartree-Fock設置半經驗或從頭算計算。本界面Charge統一設為1，表示整個體系帶一個正電荷，Spin設為Singlet，表示體系為單線態。

(3) Title部分：可以自由設置計算任務的標題。

(4) Link 0部分：計算控制設置。內存(Memory)和CPU核數(Shared Processors)可以采用默認設置，或者自定義(內存最大可設為2 GB，CPU核數建議為電腦核數的一半)。Chkpoint File指定chk文件的命名和保存路徑(chk文件在本實驗中可以用來生成分子軌道、靜電勢圖)。其他選項不設置。

(5) General、Guess部分本實驗不設置。

(6) Pop部分：設置原子電荷的計算命令，也可以在Additional Keyword部分中直接輸入原子電荷計算的關鍵詞，比如Pop = NPA，表示計算NPA電荷。

(7) Preview部分：查看編輯完成的輸入文件的內容，RC結構優化的輸入文件可參考表1。為了簡化學生建模與輸入文件編輯過程，本實驗使用的所有輸入文件都提供在補充材料中(Supporting Information)。

表1 RC分子結構優化計算輸入文件的編輯和解釋

輸入文件設置好以后點擊左下角的Submit按鈕運行Gaussian 09進行量化計算。或者也可以先打開Gaussian 09，選擇File，點擊Open，打開編輯完成的輸入文件，檢查各命令行的關鍵詞無誤，點擊Run命令，開始量化計算，如圖1(b)。

4.3 尋找過渡態和IRC計算

反應起始物RC、中間體IM1和IM2，以及產物PC可以通過4.2小節中計算步驟的描述直接進行結構優化計算得到它們穩定的分子結構。對于反應過渡態TS1、TS2和TS3，本實驗使用TS法尋找其穩定的結構[31]，在編輯過渡態計算的輸入文件時，Job Type部分選擇Opt + Freq，設定Optimize to a為TS(Berny)，Force Constants為Calculate at First Point。得到優化后的結構后，點擊Result → Vibrations查看其是否有且只有一個虛頻(過渡態結構有且只有一個數值為負的振動頻率，即虛頻(imaginary frequency))。所有分子優化結束后，按4.2小節所描述的步驟對它們進行能量計算。

點擊File → Save將優化后的過渡態結構保存為新的輸入文件，在Job Type中選擇IRC計算，設定Force Constants為Calculate Once，點擊Compute more points，設定N = 200，完成IRC計算的輸入文件編輯。最后可以得到由過渡態結構出發，連接反應物和產物兩個極小點的勢能面曲線，進一步確認該過渡態結構是否正確。需要注意的是，如果IRC計算報錯，可以一開始在Additional Keyword部分輸入IRC = LQA。IRC計算比較復雜，可以作為教師演示實驗或者給有興趣的同學額外計算。

4.4 原子電荷計算

本實驗使用Mulliken、NPA、Merz-Kollman (MK)原子電荷計算方法[32-37]，計算過渡態TS1分子的原子電荷。原子電荷計算可以通過Calculate → Gaussian Calculation Setup → Pop設置，也可以直接在Additional Keyword里添加Pop = NPA或Pop = MK關鍵詞。Mulliken電荷由Gaussian 09軟件自動計算。上述三種不同的原子電荷計算方法在下面的章節會有簡單介紹。

5 計算結果分析與討論

5.1 分子的結構優化

圖2為各反應分子優化后的結構示意圖。由此圖可知，此旋轉異構反應的RC分子的C1―C2和C1―C3鍵長大小接近，分子呈平面構型，證實了此分子具有p-π共軛效應。在TS1分子的結構中，C1―C2鍵長減小更接近C＝C，C1―H6鍵長增加，逐漸生成C3―H6鍵，鍵角A (C1―H6―C3)為56.6°，表示H6原子開始向C3遷移。IM1分子的結構中，C1―C2鍵長為1.27 ? (1 ? = 0.1 nm)，為C＝C雙鍵；C1―C3鍵長為1.41 ?，為C―C單鍵，此時H6原子遷移到C3原子上，H6與H5同側，C3原子由一開始的sp2雜化轉換為sp3雜化。甲基旋轉過渡態TS2中發生C1―C3單鍵旋轉，二面角D (H5―C2―C3―H6)由IM1中的-0.2°旋轉為TS2中的90.2°。IM2相比IM1，發生了甲基旋轉，H6與H5對側。H6原子由甲基遷移回到C1原子上的TS3與TS1類似。烯丙基正離子旋轉異構反應的PC回到與RC類似的結構，但相比于RC，PC中H6與H4同側，H6的位置發生了改變。

圖2 B3LYP/6-31G(d,p)方法下優化的分子結構以及原子序號、關鍵鍵長(?)和二面角數據

表2為B3LYP/6-31G(d,p)、B3LYP/6-31G(d)、HF/6-31G(d)和PM6方法下，優化TS1分子的結構參數和耗時數據。由表2可知：PM6方法優化的TS1分子中將要斷裂的C1―H6鍵長明顯比HF和B3LYP方法大，HF和B3LYP方法差異較小。耗時方面，PM6、HF、B3LYP方法耗時依次增加。在B3LYP泛函下，6-31G(d)和6-31G(d,p)基組得到的結構參數差別較小，但后者的耗時較大。該計算方法對比內容可以作為教師演示實驗或者給有興趣的同學額外計算。

表2 B3LYP/6-31G(d,p)、B3LYP/6-31G(d)、HF/6-31G(d)和PM6方法優化TS1分子的結構參數和耗時數據＊

5.2 熱力學數據和反應能壘圖

在GaussView軟件中，點擊Result → Summary → Thermo查看各分子的熱力學數據。其中Electronic Energy (EE)是分子的電子能，Thermal Correction to Energy/Enthalpy/Free Energy是分子的內能/焓/吉布斯自由能的校正因子，EE + Thermal Energy/Enthalpy/Free Energy Correction是分子的內能/焓/吉布斯自由能的大小。高精度的焓/吉布斯自由能=高精度單點能+結構優化的焓/吉布斯自由能校正因子。在此旋轉異構反應中，各分子的單點能、焓值、吉布斯自由能值以及相對于反應物的能量差值如表3所示。

由表3中的各分子熱力學數據和圖3的旋轉異構反應吉布斯自由能圖可知，RC和PC分子的吉布斯自由能大小一致，互為旋轉異構體，而且TS1和TS3分子以及IM1和IM2分子的吉布斯自由能大小也分別一致，顯示了前后兩次氫遷移過程的對稱性。甲基旋轉過渡態的反應能壘幾乎可忽略不計，此旋轉異構反應的決速步為氫遷移過程，其吉布斯自由能壘108.8 kJ·mol-1(26.0 kcal·mol-1)，表明此旋轉異構反應在室溫下較難發生。同時，由表3可知，該烯丙基正離子旋轉異構反應各步驟的焓變與自由能變的程度接近，說明反應過程的熵變對反應影響不大。

表3 B3LYP/6-311+G(d,p)//B3LYP/6-31G(d,p)方法下各反應分子的熱力學數據＊

圖3 B3LYP/6-311+G(d,p)//B3LYP/6-31G(d,p)方法下烯丙基正離子旋轉異構反應的吉布斯自由能圖

5.3 過渡態和IRC分析

使用GaussView軟件查看TS1、TS2和TS3分子的輸出文件，點擊Result → Vibrations查看頻率計算結果，發現這三個分子都只有一個虛頻，分別為-835.44 cm-1、-143.62 cm-1和-835.38 cm-1，點擊Start Animation可以看出TS1和TS3的虛頻對應的分子振動模式都是氫遷移過程，TS2的虛頻對應的分子振動模式是甲基旋轉過程。點擊Result → IRC/Path分析IRC計算結果并導出IRC數據，通過Origin軟件作圖得到圖4所示的以過渡態分子結構為坐標零點，連接反應物和產物的IRC曲線。IRC計算可以進一步確認過渡態TS1、TS2和TS3分子結構的正確性。

圖4 B3LYP/6-31G(d,p)方法下TS1、TS2和TS3分子的內稟反應坐標(IRC)

5.4 原子電荷

本實驗選取了三種原子電荷計算方法。其中Mulliken方法是Gaussian系統默認的計算原子電荷的方法，其將成鍵電子平均分配給成鍵原子；NPA方法按照體系的波函數劃分電子云空間分布，得到原子電荷，常用于分子體系的電荷計算；MK方法是典型的擬合靜電勢電荷的計算方法，常用于周期性體系的電荷計算。通過GaussView軟件查看原子電荷計算的輸出文件，點擊Result → Charge Distribution得到在不同原子電荷計算方法下TS1分子上各原子的電荷值。由表4數據知，不同原子電荷計算方法對于TS1分子中不同碳原子和氫原子上電荷值大小的分配存在較大差異。因此，在進行原子電荷計算時，要根據體系的性質和原子電荷計算方法的特點，選擇合適的電荷計算方法。

表4 B3LYP/6-31G(d,p)方法下TS1分子在不同電荷計算方法下的原子電荷值＊

5.5 HOMO、LUMO軌道和靜電勢圖

打開Gaussian 09軟件，點擊Utilities → formchk，將在B3LYP/6-311+G(d,p)下得到的chk文件轉換為fchk文件。使用GaussView軟件打開fchk文件，點擊Tools → Mos查看分子的軌道信息。選中目標分子軌道后，軌道會被標黃，設置Isovalue值等于0.07，點擊Visualize → Update，生成相應的分子軌道圖并保存，即得到圖5所示的各分子HOMO、LUMO軌道圖。RC分子的HOMO軌道是一個大共軛π鍵，表明其3個C原子間存在p-π共軛效應。隨著此旋轉異構反應的進行，TS1、TS3分子和IM1、IM2分子的HOMO軌道大π鍵在消失，表明分子的p-π共軛效應也在消失。

圖5 B3LYP/6-311+G(d,p)方法下各分子的HOMO、LUMO軌道

當打開fchk文件時，依次點擊Result → Surfaces/Contours → Cube Actions → New Cube，在Tpye選項中選擇Total Density，生成分子的Cube文件，調整Density值為0.11，點擊Surface Actions → New Mapped Surface，在Type選項中選擇ESP，得到圖6所示的各分子靜電勢圖。由圖6可知，RC分子的正電荷被兩側的碳原子平分，這兩個碳原子具有相同的反應活性。隨著旋轉異構反應的進行，TS1分子的兩側碳原子上的正電荷變小，IM1分子上正電荷已被轉移到中間碳原子上，之后隨著反應的進行正電荷又重新回到兩側碳原子上。

圖6 B3LYP/6-311+G(d,p)方法下各分子的靜電勢圖

HOMO、LUMO軌道和靜電勢圖的計算操作較為復雜，可以作為教師演示實驗或者給有興趣的同學額外計算。

6 實驗時間安排

本實驗面向高年級本科生或低年級研究生開設，若為線上教學需要在課前指導學生下載并安裝所需軟件，計劃使用4課時完成實驗教學。具體安排如下：1課時介紹計算化學基本概念和烯丙基正離子旋轉異構反應機理；0.5課時完成分子建模，編輯輸入文件；2.5課時完成各部分量化計算。課后安排學生獨立撰寫實驗報告，完成思考題。

7 思考題

(1)根據以上計算出的反應吉布斯自由能壘大小，判斷烯丙基正離子的旋轉異構反應在常溫常壓下能否發生？

(2)根據本實驗得到的熱力學值和反應勢能面圖，結合物理化學中的相關知識，思考反應熱力學與動力學的區別。

8 教學反思

此烯丙基正離子旋轉異構反應的計算化學實驗課程，包括了計算化學實驗的基本概念與實驗步驟。整個實驗教學內容比較基礎，方便沒有量化計算基礎的學生學習。通過本實驗的學習，可以幫助學生認識計算化學，激發學生對計算化學的學習興趣，還可以更好的理解課本中一些較為抽象的概念。

筆者在三年間多次將此實驗應用于大三本科生創新實驗、大四本科生開展畢業設計前以及低年級研究生入門計算化學的教學實驗中，均取得較好的效果。每個學生單獨進行實驗，并撰寫實驗報告，通過較短的時間就可以讓學生了解并掌握一定的計算化學實驗基本操作。另外，學生只需準備個人電腦就可以完成此實驗，十分適合疫情期間的線上教學以及線下課堂的教學。

補充材料：可通過鏈接http://www.dxhx.pku.edu.cn免費下載。