张颖 王兴 徐忠锋† 任洁茹 张艳宁 周贤明 梁昌慧 张小安3)

1) (咸阳师范学院物理与电子工程学院,咸阳 712000)

2) (西安交通大学物理学院,西安 710049)

3) (中国科学院近代物理研究所,兰州 730000)

1 引言

荷能粒子束(重离子束、β、X 和γ射线等)在癌症治疗、超高分辨生命体成像和放射性靶向药物研发等方面一直是研究热点[1–3],候氢[4]等利用分子动力学以靶受到辐射剂量为目标函数,首次对治疗过程中的关键参量束流强度和束流线型进行优化,为放射治疗提出重要的理论基础.研究表明[5],粒子束入射到生命体时在10–15—10–12s 的时间尺度内、在其径迹每沉积1 MeV 能量会产生5×104个次级离子或电子,这些次级离子或电子会引起生命体分子的诱变、重组和其他潜在的对DNA分子的致命损伤,如DNA 的单键或双键断裂.Boudaïffa 等[6]研究发现,能量在1—20 eV 的电子,与分子发生共振相互作用,即共振电子吸附形成瞬态分子(其寿命在飞秒量级),然后发生分子解离 (dissociative electron attachment,DEA),这在DNA 辐射分解的初期阶段起着至关重要的作用,并且可能在它们热化分解之前很久就已经引起了实质性的损害,这个结果对传统的观点(只有离子能量大于生命分子的离化阈值或解离能时才会对生命体有毒害)提出根本性挑战.进一步的研究发现[7,8],在次级离子中,能量低于4 eV 的自由电子会进入生命体分子形状共振形成的反键(π*)分子轨道,形成瞬态分子,它发生二次解离时,会产生具有很高反应活性的自由基和离子,进而更多地破坏生物分子的结构和功能,对生物分子造成不可逆的辐射损伤.

低能电子对生物体系的辐射损伤并不都是有害的,比如进行化疗和放疗时,低能电子吸附解离治疗药物,其产物可以破坏肿瘤细胞的DNA,使其无法进行DNA 复制,从而达到减少或去除肿瘤组织的目的[9].研究表明[10],在构成核糖核酸(ribonucleic acid,RNA) 4 种碱基中,脲嘧啶(uracil) 取代了DNA 中的胸腺嘧啶而成为RNA的特征碱基,脲嘧啶的特别点位容易被铁、硫等元素发生化学修饰,硫脲嘧啶(2-thiouracil)是尿嘧啶中碳原子位置上的氧原子被硫原子所取代,属于核苷类似物.近年的研究表明[11,12],反应性F 或S 修饰后的尿嘧啶自由基产生的碱基(氟脲嘧啶和硫脲嘧啶等) 成为DNA 损伤的前兆物,有望成为生命体演化、病毒变异过程中分子诊断的标志物.硫脲嘧啶由于其具有特殊的光毒性,是治疗白血病和银屑病的首选药物.因此,硫脲嘧啶及其结构类似物在近年来受到相当大的关注,并且由于其显著的生物学和化学治疗特性,已经从实验和理论的角度对其进行了广泛的研究[13–16].另一方面,Plekan等[17]的研究表明硫脲嘧啶及其互变异构体成生物传感器、生物芯片等相关研究的热点之一.目前,对其与低能电子的吸附解离过程研究的较少.就仅有的实验报道[18,19]来看,不同的实验条件和离子探测手段造成实验结果存在明显的差异,其中脱氢过程是源于碳氢键的断裂还是氮氢键断裂,尚不能确定.理论上,研究多集中于对吸附解离过程中负离子的束缚态和解离过程所涉及的共振类型的研究,Kopyra 等[18]利用含有赝势的Schwinger multichannel 变分法,Dolgounitcheva 等[20]采用二阶微扰理论方法和耦合簇方法,Costa 等[21,22]分别利用单中心展开的量子散射方法,这些研究得出相同的结论,即硫脲嘧啶分子与低能电子吸附解离过程中发生了形状共振.由于硫脲嘧啶分子在气相中存在多种互变异构体,针对硫脲嘧啶分子及其互变异构体的动态解离过程的理论计算未曾报道过.因此,我们无法深刻认识低能电子与硫脲嘧啶分子吸附解离中的脱氢过程.虽然上述3 种方法很好地描述了吸附解离过程,然而硫脲嘧啶分子的吸附解离过程中存在解离通道的多样性和复杂性的问题,这些方法在计算量上也遇到巨大的挑战.此外,这些量子计算的方法在处理化学键的离散态和电子连续性耦合的问题并不理想[23–26].生物分子解离动力学过程主要是由离化和电子吸附两个过程引发的,而从头算分子动力学模拟恰好可以合理地描述这些过程.基于以上考虑,本文将采用从头算动力学模拟中的玻恩-奥本海默分子动力学方法来描述硫脲嘧啶分子及其互变异构体在吸附解离过程中的动态变化情况,可以更为深刻地认识硫脲嘧啶分子解离过程中脱氢的原因.

2 计算模型与方法

2.1 计算模型

在玻恩-奥本海默分子动力学 (Born-Oppenheimer molecular dynamics,BOMD) 方法中,电子运动的问题用Schrödindger 方程来进行计算;对于离子运动的问题,可以忽略其受到量子效应的影响,利用经典动力学法(牛顿力学)来处理,其运动方程具体表示为

其中,MI表示第I个离子的质量,RI表示第I个离子的位置,FI表示第I个离子受到的力.

根据费曼-海尔曼定理,即可计算离子的受力情况:

式中,ψ0为电子的基态波函数,He表示为哈密顿算子.由BOMD 分子动力学中拉格朗日函数LBO可知:

进行每一步动力学模拟时,依据原子核的空间位置对体系的电子结构进行重新计算,使得电子达到基态.正是由于计算量大,过程比较耗时,并且计算时需要考虑原子中电子结构,BOMD 模拟目前只适合用于小分子体系.López-Tarifa 等[27]运用BOMD 模型结合含时密度泛函理论模拟研究了气相下质子碰撞尿嘧啶分子并形成2 价离子的解离过程.利用分子动力学模拟结果成功分辨了离子飞行时间符合谱中多种复杂的解离通道.这不仅表明BOMD 方法结合密度泛函理论可以很好地计算并还原气相条件下分子动态解离过程,还近一步表明结合密度泛函理论可以大大降低计算量.基于此,本文将通过使用BOMD 方法结合密度泛函理论以及6-311++G(2d,2p) 基组模拟硫脲嘧啶及其互变异构体与低能电子吸附解离动态过程,从而更为深刻地认识硫脲嘧啶分子解离过程中脱氢的原因.

2.2 计算方法

采用密度泛函理论中的B3LYP 杂化泛函方法对硫脲嘧啶分子在气相条件下存在的几何构型进行结构优化和频率计算,以获得最稳定的4 种构型.这4 种几何构型均未出现虚频.为了得到更为准确的能量,分别采用多体微扰理论中的(secondorder Moller-Plesset perturbational,MP2)方法和耦合簇 (coupled cluster singles and doubles,CCSD)方法在优化后的构型基础上进行能量计算.所有的能量都进行了零点振动能修正.以上计算采用的是6-311++G(2d,2p) 基组,该基组是在劈裂价键基组的基础上先对碳原子分别增加弥散s,p 高斯函数和极化函数,对氢原子增加弥散s 高斯函数和极化函数;再对碳原子增加两组d 轨道函数,对氢原子增加两组p 轨道函数.

采用BOMD 模拟研究硫脲嘧啶及其互变异构体动态解离过程.BOMD 动力学模拟是从电子垂直吸附到硫脲嘧啶分子及其互变异构体的Franck-Condon 区域开始计算.对于硫脲嘧啶解离的准经典轨迹取样是从初始内能为1.0 eV 开始采样,即电子吸附能在1.0 eV 附近开始模拟吸附解离过程,动力学模拟的时间步长是0.1 fs,模拟的时间尺度是400 fs.充分考虑到计算量和准确描述电子相关,所以BOMD 动力学模拟是采取密度泛函理论中的B3LYP 杂化泛函方法和6-311++G(2d,2p) 基组相结合的方式完成的.所有研究的解离路径在模拟过程中没有设置任何约束条件.由于受到计算资源的限制,所以无法在有限的解离路径数目中展示出吸附解离过程中所有可能存在的解离过程.本文所有计算均在Gaussian 09 软件包下完成.

3 结果与讨论

3.1 硫脲嘧啶分子互变异构体基态几何结构的相对稳定性

根据实验提供硫脲嘧啶分子气相下可能存在的互变异构体几何构型的基础上,采用密度泛函理论中的B3LYP 杂化泛函的方法在6-311++G(2d,2p)基组水平上进行无限制全梯度优化获得4 种几何构型,同时进行频率分析,所获得的振动频率并未出现负值,这表明硫脲嘧啶的3 种互变异构体在势能面上的局域最小值点,是稳定的互变异构体,其几何结构如图1 所示.优化后硫脲嘧啶及其3 个互变异构体的几何结构均为平面结构,具有C1对称性.表1 列出了采用CCSD 方法和MP2 方法计算分子的相对能量,和相关文献报道的计算结果,所有计算结果均进行了零点振动能的修正.在此基础上结合能的计算还采用了Boys 和Benardi[28]的完全均衡校正(counterposie procedure,CP)消除基组重叠误差(basis set superposition error,BSSE).根据表1 可以看出,采用CCSD 方法计算获得硫脲嘧啶互变异构体的相对能量分别要比硫脲嘧啶高46.31,93.59 和 70.43 kJ/mol,同时采用MP2计算获得硫脲嘧啶互变异构体的能量同样分别要比硫脲嘧啶高30.88,61.75 和50.17 kJ/mol.硫脲嘧啶及其互变异构体相对稳定性从大到小的表述顺序为: 2-thiouracil ＞ Tautomer 1 ＞ Tautomer 3＞Tautomer 2.通过比较不同方法计算结果可以发现,CCSD 和MP2 计算方法获得硫脲嘧啶以及互变异构体的相对稳定性描述与文献[18]中报道的B3LYP 和CCSD(T)计算方法获得的结果是一致的.另外采用CCSD 方法对硫脲嘧啶分子及其互变异构体的偶极矩也进行了计算,发现获得的偶极矩结果均大于2.5 deb (1 deb=3.33564×10–30C·m–1).现有的研究结果表明,较大的偶极矩可以促进吸附解离过程中偶极束缚负离子的形成.较大的偶极矩同时也意味着存在较强的静电偶极场,其会率先俘获入射的低能电子.但是在静电偶极场作用下电荷发生转移的概率极低,不利于采用BOMD 模拟方法观察硫脲嘧啶及其互变异构体的吸附解离过程.经过CP 校正后硫脲嘧啶的结合能最大,互变异构体2 的结合能最小.

表1 硫脲嘧啶分子及其互变异构体稳定构型的相对能量,基组重叠误差BSSE,CP 校正后的结合能 ΔECP 和偶极矩μTable 1.Relative energies,BSSE,CP-corrected binding energies ΔECP and dipole moment of 2-thiouracil and its tautomers.

图1 采用B3LYP/6-311++G(2d,2p)方法优化后基态硫脲嘧啶分子及其互变异构体的几何构型 (图中红色为氧原子,白色为氢原子,灰色为碳原子,蓝色为氮原子,黄色为硫原子)Fig.1.Optimized geometries of the ground state 2-thiouracil and its tautomers.The geometries have been optimized using the B3LYP method and the 6-311++G(2d,2p) basis set (Colors in graphical representations correspond to atom types as follows,red: oxygen;white: hydrogen;gray: carbon;blue: nitrogen;yellow: sulfur).

3.2 BOMD 模拟

图2 所示的是室温条件下BOMD 模拟硫脲嘧啶及其互变异构体1 解离动态过程.在硫脲嘧啶在吸附解离过程中化学键N1—H 键在67.2 fs 发生断裂.当化学键N3—H 键的键长被拉伸超过3 Å时发生断裂(107.4 fs).当BOMD 模拟进行到373.6 fs时,化学键 C6—H 键发生断裂.在以上的解离过程中,发现在N1,N3和C6的位点均失去了氢原子,进而导致闭壳层脱氢负离子(TU-H)–的形成.对于硫脲嘧啶在解离过程所形成的脱氢负离子,在BOMD 模拟过程中发现其解离通道为

图2 硫脲嘧啶及其互变异构体1 在解离过程中分子结构随时间演化图Fig.2.The evolution of the geometrical structures of 2-thiouracil and tautomer 1 during dissociation.

由于充分考虑到计算时间和动态解离计算过程中的自洽场计算收敛的问题,对脱氢负离子形成所需的解离通道进行简单的统计分析发现,大约90.3%解离事例均出现脱氢解离通道.从图2 还发现,在BOMD 模拟硫脲嘧啶解离过程中,化学键C—S键在384 fs 发生断裂,即在C2的位点失去了硫原子,从而形成了负离子碎片(TU-S)–,其解离通道可具体表示为

同样对负离子碎片(TU−S)−形成所需的解离通道进行简单的统计分析后发现,大约83.1%解离事例均出现(5)式的解离过程.BOMD 模拟硫脲嘧啶动态解离过程的结果与文献[18]报道的实验结果符合得非常好.

图3 是在室温条件下BOMD 模拟硫脲嘧啶互变异构体2 和互变异构体3 的解离动态过程.由于低能电子吸附解离过程高度依赖于分子构型,另外互变异构体2,互变异构体3 与硫脲嘧啶分子构型相比存在明显的差异,所以BOMD 模拟互变异构体2 以及互变异构体3 的动态解离过程与硫脲嘧啶动态解离过程存在显著的差异.在硫脲嘧啶互变异构体2 解离过程中,化学键O—H 键最先发生断裂(158.7 fs),随后是化学键C2—H 键在284.1 fs发生断裂,最后是化学键S—H 键在388.4 fs 发生断裂.由于化学键C2—H 键断裂,C2位点也就失去了氢原子,所以在互变异构体2 解离过程中也形成了闭壳层脱氢负离子(TU-H)–,其解离通道与(1)式相同.另外,O 位点和S 位点同样失去了氢原子,所以在解离过程中形成了另一种负离子碎片(TU-2H)–,其解离通道可以按照如下表示:

图3 硫脲嘧啶互变异构体2 和其互变异构体3 在解离过程中分子结构随时间演化图Fig.3.The evolution of the geometrical structures of tautomer 2 and tautomer 3 during dissociation.

或者

在进行BOMD 模拟解离过程中对形成的通道(6)式和(7)式进行样本抽样并统计分析获得,约43.1% 解离事例出现(6)式解离过程,约31.5%解离事例出现(7)式解离过程.目前开展的硫脲嘧啶吸附解离实验研究中已经探测到负离子碎片(TU-2H)–的存在,这就意味着实验过程中解离通道(6)式和(7)式可能是同时存在的.

在图3 所示BOMD 动力学模拟硫脲嘧啶互变异构体3 解离动态过程中,当模拟进行到159.6 fs,化学键N1—H 键断裂.当化学键C6—H 键的键长被拉伸超过4 Å时发生断裂(227.4 fs).由于在解离过程中N1和C6位点均失去了氢原子,所以在互变异构体3 的解离过程中同样也形成了脱氢负离子(TU-H)–,其解离通道与(4)式相同.从图3 中还可以看出,化学键S—H 键也发生断裂,其时间在393.6 fs.这与硫脲嘧啶互变异构体2 的BOMD模拟解离过程十分相似.

3.3 解离动力学分析

采用密度泛函的B3LYP 方法结合6-311++G(2d,2p)基组对硫脲嘧啶在吸附解离过程中所出现过渡态(transition state,TS)和产物(product,P)的几何构型进行了全优化.通过振动频率分析,确认了过渡态的真实性,其中各个过渡态均有唯一虚频,并且利用内禀反应坐标理论对过渡态进行了验证.为了获得更为精确的能量信息,再进一步使用CCSD/6-311++G(2d,2p)方法对所有过渡态的能量进行单点能校正.从图4 可以看出,硫脲嘧啶需要克服势垒83.65 kJ/mol 才能发生解离.整个解离过程出现3 个过渡态TS1,TS2 和TS3,这3 个过渡态涉及了解离过程中不同化学键断裂的情况,其中过渡态TS 1 的几何构型涉及化学键N—H 键的断裂,过渡态TS2 的几何构型涉及化学键C—H 键的断裂,而过渡态TS3 的几何构型涉及化学键C—S 键的断裂.过渡态TS2 能量比过渡态TS1 高出14.34 kJ/mol,过渡态TS3 的能量则比过渡态TS1 高出83.47 kJ/mol,所有的能量值均为正值,并且解离过程为吸热过程.从热力学角度综合考虑,在硫脲嘧啶解离过程中化学键N—H 键比化学键C—H 键更容易断裂.

图4 硫脲嘧啶解离过程势能剖面示意图Fig.4.The potential energy profile of the 2-thiouracil dissociation.

硫脲嘧啶在BOMD 模拟吸附解离过程中化学键N—H 键和C—H 键的势能曲线被绘制在图5中.在构建势能曲线时,分别对化学键N—H 键和C—H 键进行固定步长的逐点扫描,扫描范围从0.01—3.00 Å.在扫描势能曲线时使用的是DFT/B3LYP 的理论方法,这个方法虽然在研究生物分子吸附解离时会有一定的误差,但已经证实了其是一种可靠的定性研究吸附解离过程能量路径的方法.从势能曲线图可以知道,化学键N1—H 键和N3—H 键解离势垒分别为73.63 kJ/mol和97.94 kJ/mol.随着化学键N1—H 键和N3—H键的键长不断增大,N1—H 键和N3—H 键的势能曲线也随即上升.当化学键N1—H 键和N3—H 键的键长分别为0.96 Å和0.97 Å,势能曲线达到其最大值,随后N1—H 键和N3—H 键的势能曲线逐步稳定.在硫脲嘧啶吸附解离过程中化学键N1—H 键比N3—H 键更容易发生断裂.与化学键N—H 键势能曲线相比,化学键C6—H 键的解离势垒明显高化学键于N—H 键.这就表明在吸附解离动力学过程中,化学键N—H 键要比C—H 键更容易发生断裂.

图5 B3LYP/6-311++G(2d,2p)方法计算获得硫脲嘧啶解离过程中化学键N—H 键和C—H 键的势能曲线Fig.5.The potential energy curves for the N—H and the C—H bond in 2-thiouracil dissociation obtained from B3LYP/6-311++G(2d,2p) calculations.

采用相同的方法获得硫脲嘧啶互变异构体2 解离过程的势能剖面见图6.计算发现2 个中间体(intermediate,IM),3 个过渡态以及1 个解离产物.结果表明,互变异构体2 与硫脲嘧啶的解离过程类似,互变异构体2 需要克服势垒123.28 kJ/mol 才能解离,并且吸收的热量显著增大.解离过程主要包括化学键C—H 键、化学键O—H 键和化学键S—H 键的断裂.互变异构体2的解离产物能量为57.32 kJ/mol,而硫脲嘧啶的解离产物能量为12.15 kJ/mol.二者相比之下,硫脲嘧啶可以解离出热力学更为稳定的产物(结构见图4),其解离过程也更易进行一些.

图6 硫脲嘧啶互变异构体2 解离过程势能剖面示意图Fig.6.The potential energy profile of the tautomer 2 dissociation.

使用B3LYP/6-311++G(2d,2p) 的方法构建了硫脲嘧啶互变异构体2 在BOMD 模拟吸附解离过程中化学键C—H 键、O—H 键和S—H 键的势能曲线,如图7 所示,在构建曲线时,分别对化学键C—H 键、O—H 键和S—H 键的键进行固定步长的逐点扫描,扫描范围从0.01—3.00 Å.从势能曲线图可以看出,硫脲嘧啶互变异构体2 中的化学键O—H 键的解离势垒(115.05 kJ/mol)低于其C—H 键和S—H 键解离势垒(164.01 kJ/mol和185.35 kJ/mol),硫脲嘧啶互变异构体2 的化学键C—H 键、O—H 键和S—H 键的势能曲线变化趋势相同,首先化学键C—H 键、O—H 键和S—H键的键长在逐渐拉伸的过程中,其势能曲线迅速上升.当化学键C—H 键、O—H 键和S—H 键的键长被逐渐拉伸至2.25,1.72 和2.01 Å时,势能迅速达到最大值,势能曲线逐步平稳.这就表明在硫脲嘧啶互变异构体2 吸附解离过程中,化学键O—H键最容易发生断裂,其次是C—H 键发生断裂,最后是S—H 键发生断裂.这从能量角度解释了BOMD模拟吸附解离过程中化学键C—H 键、O—H 键和S—H 键在不同时间断裂的原因.

图7 B3LYP/6-311++G(2d,2p)方法计算获得硫脲嘧啶互变异构体2 解离过程中化学键C—H 键、O—H 键和S—H键的势能曲线Fig.7.The potential energy curves for the C—H,the O—H and the S—H bond in tautomer 2 dissociation obtained from B3LYP/6-311++G(2d,2p) calculations.

图8 是通过B3LYP/6-311++G(2d,2p)的计算方法获得硫脲嘧啶互变异构体3 在BOMD 模拟吸附解离过程中化学键N1—H 键的势能曲线.在构建曲线时,对化学键N—H 键进行了固定步长的逐点扫描,扫描范围从0.01—3.00 Å.从图8 可以看出,硫脲嘧啶互变异构体3 的化学键N1—H 键的解离势垒是102.63 kJ/mol.随着互变异构体3中化学键N1—H 键的键长逐渐拉伸,其势能曲线也随即上升,直至互变异构体3 的化学键N1—H键的键长被拉伸至1.01 Å时,势能达到最大值,随后势能曲线逐步稳定.与硫脲嘧啶的化学键N1—H 键的势能曲线相比,互变异构体3 的化学键N1—H 键解离势垒明显高于硫脲嘧啶N1—H 键的解离势垒.这就说明在吸附解离动力学过程中,硫脲嘧啶的化学键N1—H 键要比互变异构体3的化学键N1—H 键更容易发生断裂.这就解释了BOMD 模拟解离过程中硫脲嘧啶的化学键N1—H键断裂时间快于互变异构体3 的化学键N1—H 键断裂时间的现象.势能曲线为研究硫脲嘧啶吸附解离中脱氢过程提供了必要的证据.

图8 B3LYP/6-311++G(2d,2p)方法计算获得硫脲嘧啶互变异构体3 解离过程化学键N1—H 键的势能曲线Fig.8.The potential energy curves for the N1—H bond in tautomer 3 dissociation obtained from B3LYP/6-311++G(2d,2p) calculations.

在硫脲嘧啶及其互变异构体与低能电子吸附解离过程中形成的主要负离子碎片是闭壳层脱氢负离子(TU-H)–.在BOMD 模拟过程中,发现解离通道(4)会发生在不同的位置,例如,在硫脲嘧啶吸附解离过程中,氢自由基是通过N1—H 键,N3—H 键和C—H 键化学键断裂而形成的;在互变异构体2 吸附解离过程中,氢自由基的形成是源自于C—H 键,O—H 键和S—H 键化学键的断裂;而在互变异构体3 吸附解离过程中,由于N1—H键,C—H 键和S—H 键化学键的断裂,使得氢自由基的产生,对于所有可能形成的(4)式解离通道,采用CCSD 方法计算了不同位置脱氢所需键解离能,硫脲嘧啶互变异构体的电子亲和能以及脱氢负离子(TU-H)–的出现势.表2 列出了硫脲嘧啶及其3 个互变异构体分子不同位置脱氢所需的键解离能(bond dissociation energy,BE),绝热电子亲和能(adiabatic electron affinity,EA) 和(TU-H)–负离子出现势(appearance energy,AE) .计算结果从能量的角度表明,在硫脲嘧啶分子及其互变异构体中,N1位置最容易失去氢原子而形成负离子碎片(TU-H)–.在硫脲嘧啶和互变异构体1 中,C 位置失去氢原子而形成负离子碎片(TU-H)–相对较难,其值分别为378.22 kJ/mol 和382.08 kJ/mol.值得注意的是,在互变异构体2 吸附解离过程中,负离子碎片(TU-H)–的形成也相对较难.Kopyra等[18]采用G2(MP2)计算方法也得到相似的结果.N1位置最容易失去氢原子而形成负离子碎片(TU-H)–,其解离能为53.07 kJ/mol,其次是N3位置,解离能为96.49 kJ/mol.值得注意的是,硫脲嘧啶与3 个互变异构体解离过程中,都具有较小的绝热电子亲和能(小于289.46 kJ/mol),尤其是硫脲嘧啶绝热电子亲和能小于213.23 kJ/mol.这就表明,入射电子会大量吸附在化学键N1—H 键,N3—H 键和C—H 键的周围,对其断裂有明显的促进作用.同时,CCSD 计算的结果也得到了硫脲嘧啶分子与低能电子(289.46 kJ/mol) 吸附解离实验定性的验证.

表2 由CCSD 计算方法获得硫脲嘧啶及其互变异构体分子不同位置脱氢所需的键解离能,绝热电子亲和能和(TU-H)–负离子出现势(能量单位为kJ/mol)Table 2.The bond dissociation energy,the adiabatic electron affinity and the appearance potential of the (TU-H)– obtained from CCSD calculation method of dehydrogenation of 2-thiouracil and its tautomer at different positions (all values in kJ/mol).

在BOMD 模拟硫脲嘧啶分子及其互变异构体吸附解离过程中,获得了负离子碎片(TU-2H)–.该负离子碎片可能是分别通过(5)式和(6)式解离通道形成.基于此,采用CCSD 方法对硫脲嘧啶分子及其互变异构体所有可能形成(5)式和(6)式的解离通道进行最小解离能计算,确认失去两个H 原子所需的最小解离和失去一个H2分子所需的最小解离能,其结果列于表3.从表3 可以看出,解离通道(6)式所需的解离能范围为89.10—362.78 kJ/mol,而解离通道(5)式所需的解离能在446.73—651.27 kJ/mol 范围内.在比较两个解离通道的最小解离能之后,得出在硫脲嘧啶及其互变异构体分子吸附解离过程中更容易失去一个H2分子.在解离过程中,硫脲嘧啶分子最容易失去氢分子的位置是N1和N3位置,其最小解离能量是88.77 kJ/mol;其次是互变异构体1 中N1和N3位置,最小解离能为91.66 kJ/mol.

表3 由CCSD 计算方法获得硫脲嘧啶及其互变异构体解离通道(6)式和(7)式所需的最小解离能(能量单位为kJ/mol).Table 3.The minimum dissociation energy required to obtain 2-thiouracil and its tautomer dissociation channels (6) and (7) by CCSD calculation method (all values in kJ/mol).

4 结论

通过运用BOMD 模型结合密度泛函理论的B3LYP 杂化泛函,以及6-311++G(2d,2p) 基组,探究了气相下硫脲嘧啶与低能电子吸附解离动力学过程.通过高精度的理论计算(CCSD/6-311++G(2d,2p)),表明硫脲嘧啶及其3 个互变异构体可以在气相下稳定共存.BOMD 模拟结果发现硫脲嘧啶及其互变异构体在吸附解离过程中一个重要的脱氢现象,即特定位置的化学键N—H键,C—H 键断裂,使其在N 位点和C 位点失去氢原子,形成解离过程中最主要的负离子碎片——闭壳层脱氢负离子(TU-H)–;同时,还观察到了负离子(TU-S)–和(TU-2H)–的形成.为了更为深刻认识脱氢现象,分别对断裂的化学键进行势能曲线扫描,键解离能,电子亲和能以及负离子(TU-H)–出现势进行计算对比,发现在所有断裂的化学键中N—H 键从能量的角度上最容易发生断裂,负离子(TU-H)–的形成主要源于N—H 键的断裂.本研究的理论计算结果与现有实验研究结果符合得较好,表明所选取的计算方法是完全可靠的.BOMD 模拟不仅可以动态还原硫脲嘧啶与低能电子吸附解离过程,更重要的是可以深刻地认识硫脲嘧啶分子吸附解离过程中脱氢的机理和解离通道.