基于ReaxFF力场的对硝基苯酚臭氧氧化分子动力学模拟

王子民 郑 默 谢勇冰 李晓霞,4,* 曾 鸣曹宏斌,* 郭 力,4

(1中国矿业大学(北京)化学与环境工程学院，北京100083；2中国科学院过程工程研究所复杂系统国家重点实验室，北京100190；3中国科学院过程工程研究所环境技术与工程研究部，北京100190；4中国科学院大学化学与化工学院，北京100049)

基于ReaxFF力场的对硝基苯酚臭氧氧化分子动力学模拟

王子民1,2郑 默2谢勇冰3李晓霞2,4,*曾 鸣1曹宏斌3,*郭 力2,4

(1中国矿业大学(北京)化学与环境工程学院，北京100083；2中国科学院过程工程研究所复杂系统国家重点实验室，北京100190；3中国科学院过程工程研究所环境技术与工程研究部，北京100190；4中国科学院大学化学与化工学院，北京100049)

本文采用基于ReaxFF反应力场的分子动力学方法(ReaxFF MD)，利用自主研发的国际首个基于GPU加速的ReaxFF MD程序系统GMD-Reax和独特的化学反应分析工具VARxMD，探索臭氧氧化对硝基苯酚的反应机理。通过模拟考察了300 K恒温条件下臭氧氧化水中对硝基苯酚的行为，获得了酚结构开环、CO2生成、主要自由基(·OH、·O2、·O)及团簇型自由基的数量演变趋势，并可定性描述六元环开环和CO2生成均遵循伪一级反应动力学规律。反应机理分析表明酚类分子在水溶液中被臭氧氧化的路径主要经过攫氢、六元环开环、碳链的氧化分解三个阶段，也揭示了自由基和团簇型自由基在臭氧降解对硝基苯酚时所发挥的重要作用。本工作是应用ReaxFF MD分子模拟方法对常温水环境下臭氧降解酚类污染物反应机理研究的一个尝试，可为深入认识该机理及相关的实验、理论研究提供一定的参考。

对硝基苯酚；臭氧；反应分子动力学模拟；自由基行为；反应机理

1 引 言

煤焦废水处理是水环境治理中的重大问题之一，一直备受关注1,2。焦化废水主要来源于原煤高温干馏、焦炉煤气初冷和化工产品精制等工业过程，含大量难降解有机污染物，主要包括酚类及其衍生物、喹啉类化合物、苯、吡啶、萘类化合物等。这些污染物可长期存在于自然环境中，危害人类及其它生物的生命健康3,4。随着人类工业的发展，焦化废水的处理方法在不断优化改进，从最早多段曝气的活性污泥法5到应用广泛的生物法6,7，逐步提高了废水中有机物的去除效率8。但是，这些方法并不能保证水质的绝对安全9，生物法处理的出水中仍含有少量苯系物、杂环化合物等有毒有机物质，因此针对生物段出水的深度处理非常重要。臭氧氧化法作为深度处理废水的一种有效方法，具有无二次污染、反应迅速、适用范围广、氧化彻底等优势10，逐渐成为人们关注的热点。目前，臭氧氧化法已成功应用于纺织、皮革、石油化工厂污水的处理11，并逐步应用于焦化废水的深度处理。

臭氧氧化包括自由基间接氧化和臭氧分子直接氧化两种路径。间接氧化反应比直接氧化更为迅速，可氧化的污染物种类也更多11。尽管自由基对臭氧氧化过程非常重要10,12，但目前的实验手段还难以对自由基进行连续原位检测，加上焦化废水中有毒有机污染物组成复杂，人们对臭氧降解有机物的化学反应机理的认识并不全面。Utsumi等12利用电子自旋共振技术(ESR)检测臭氧在水环境中产生的活泼·OH的浓度，获得·OH在臭氧间接氧化过程中起主导作用的结论。常温下臭氧在水中溶解度低13，臭氧氧化过程产生的·OH 在水中存活时间极短14，因此，充分认识·OH自由基在臭氧氧化中的作用，提高其数量和反应活性对提高臭氧的氧化效率非常重要。将臭氧化与其它方法联用，如 O3/H2O215、O3/H2O2/UV16、催化臭氧氧化17,18都可能提高·OH自由基的浓度，强化臭氧化作用。利用高效液相色谱(HPLC)、离子色谱(ion chromatograph)或气质联用仪(GC-MS)可检测臭氧降解酚类污染物中产物的化学结构，并大致推断出酚-醌-羧酸-二氧化碳降解历程19,20。尽管如此，由于不稳定的降解产物出现时间极短，常规实验手段难以捕捉瞬时出现的中间体结构，也难以连续检测降解反应过程中相关小分子自由基的数量和演变趋势，无法获得详细的臭氧降解有机污染物的反应机理。近几年，分子模拟方法的快速发展为探索臭氧降解焦化废水中有机污染物的反应机理提供了另一条途径。

反应分子动力学方法(ReaxFF MD，Reactive Force Field Molecular Dynamics)是由 van Duin和Goddard等21共同发展的化学反应力场ReaxFF与分子动力学相结合的一种模拟方法。ReaxFF21基于键级描述化学体系中的键相互作用，利用 DFT拟合与优化原子与原子间的相互作用参数，并通过电负性平衡方法(EEM)22动态更新模拟体系中每个原子的电荷。与量子化学方法相比，ReaxFF MD无需预先定义化学反应路径，并可应用于较大规模体系(＞ 103个原子)的模拟。

近年来，ReaxFF力场得到快速发展，针对不同元素的力场参数不断被提出23,24，扩展了ReaxFF的应用范围，包括燃烧25、热解26,27、催化28等过程的反应机理研究。其中，Rahaman等29发展了可应用于水环境下有机化学反应的ReaxFF力场参数，并用于描述氨基酸的异构化过程，获得了甘氨酸在水环境中从中性到两性离子状态的转变机理，并发现其中的质子转移反应可能是经由单一水分子完成。姜丹丹等30利用ReaxFF MD模拟了以苯并芘为代表的焦化废水在超临界水环境下降解过程中的结构变化，得到了三种主要的降解机理。张金利等31通过ReaxFF MD研究在超临界水环境下降解印染废水中的分散橙 25(Disperse Orange 25，DO25)，发现温度、分子比例(DO25/H2O/O2)及反应时间对降解过程均有影响，其中自由基 HO·、HO2·等可加速 DO25的降解。这些应用显示了ReaxFF MD在应用于水环境下降解有机污染物的反应研究的潜力。

为了进一步提高ReaxFF MD的计算效率和模拟规模，郑默等32开发了国际上首个基于图形处理器(GPU)加速的ReaxFF MD程序GMD-Reax。ReaxFF MD的GPU并行化使得大规模反应体系(～10000个原子)在桌面计算机上的快速模拟成为可能。刘健33和韩君易34等开发的 VARxMD(Visualization and Analysis of ReaxFF Molecular Dynamics)是国际上首个实现对 ReaxFF MD模拟轨迹蕴含的化学反应进行自动化系统分析和可视化的程序系统。VARxMD可识别两个时刻之间的化学反应，自动获得模拟过程中发生的完整化学反应列表，进而对物种进行分析统计，显示物种和化学反应的 2D、3D结构，并基于全局物种结构特征筛选化学反应，为复杂化学反应的分析提供了便利。GMD-Reax与VARxMD的成功研发及应用25−27,35，显著提高了模拟效率，减少了分析ReaxFF MD模拟轨迹所消耗的人力和时间，有助于获得完整的化学反应机理信息。

本文采用基于 GPU 高性能计算的GMD-Reax，利用VARxMD独特的化学反应分析能力，以对硝基苯酚作为焦化废水中有机污染物的代表性组分36,37，对大规模有机污染物体系进行了臭氧氧化的ReaxFF MD模拟，从分子层面获得了有机物降解规律以及降解过程中重要自由基的演化规律和影响，期望为进一步认识臭氧降解对硝基苯酚的机理提供帮助。

2 模拟策略

2.1 模型构建

利用 ReaxFF MD从分子水平探索臭氧氧化有机物的反应机理，构建合理的模型体系是模拟的基础，大规模模型对描述反应多样性至关重要，但是会受到后续模拟计算代价的制约。酚类污染物是焦化废水中的特征有机污染物，生物法降解的效果不好38，本文选择对硝基苯酚作为臭氧氧化的目标污染物。焦化废水生化出水酚的实际浓度多为 100−300 mg·L−1，甚至更低39,40。从模拟的角度看，依照实际浓度构建的体系不利于观察臭氧与对硝基苯酚的反应情况，也不利于得到降解过程中的规律。考虑到实际实验中臭氧分子是连续通入，为了在可接受的计算代价内观察到反应趋势，本文假设酚分子附近局部反应环境中的臭氧过量。此外，体系中的水环境会带来大量的分子间氢键，易形成水团簇结构，不仅消耗宝贵的计算资源，也使化学反应分析更为复杂。因此本文进一步假设酚类有机物周围的反应环境为臭氧分子的数量多于水分子的数量。基于上述假设，同时为了考察羟基与臭氧的协同作用，构建了不同组成的臭氧氧化水中对硝基苯酚的模型体系，表1列出了体系的组成、性质及模拟温度(为简略表示，对硝基苯酚缩写为p-NP)。

模型构建过程中，首先基于 Materials Studio(MS)41中的Modules模块，使用Amorphous Cell中Construction功能构建臭氧-对硝基苯酚-羟基-水的初始结构模型，密度均设为 1.0 g·cm−3，各体系构建时均采用 Dreiding42力场。接着利用 MS中的Forcite模块对初始体系进行几何优化(Geometry Optimization)，降低初始模型的能量，同时防止对硝基苯酚体系中各分子和重要官能团相互重叠。对优化后的初始体系在300 K下进行NVT(恒温恒容)分子动力学驰豫，同样采用Dreiding力场，模拟时长200 ps，模拟步长为1.0 fs，采用Berendsen43方法控温。最后再次进行几何构型最优化。图1为含300个对硝基苯酚、500个臭氧、500个水分子的初始体系3D结构。

2.2 ReaxFF MD模拟细节及反应分析

ReaxFF力场通过键级来描述体系中的成键和断键，ReaxFF与MD结合可模拟较大的体系，且不需要对反应路径进行预先定义，有潜力研究臭氧降解酚分子的相关自由基和反应机理。本文的ReaxFF MD模拟研究采用NVT系综，为了保持一致，模拟温度采用与实验温度接近的300 K。模拟采用周期性边界条件，Berendsen方法控温，控温系数为0.1 ps，同时采用velocity-Verlet算法动态地更新原子在每个时间步的坐标和速度，模拟步长为0.4 fs，共模拟200 ps。模拟程序采用基于GPU并行的GMD-Reax32，模拟在CentOS 5.4服务器(处理器型号为英特尔-至强E5620、频率为2.4 GHz、内存为2 GB、GPU卡型号为C2050)上进行，并每隔0.4 ps输出一次反应轨迹和键级轨迹。每两帧轨迹对应的两个时刻之间有大量化学反应发生；且由于该模拟体系含有大量水，在ReaxFF MD模拟体系中会出现大量的水团簇分子，因此模拟结果的反应机理分析是一项有挑战性的工作。本文采用VARxMD程序系统33,34对ReaxFF MD模拟结果进行化学反应的自动分析和物种自动检索，提升了化学反应机理探索的效率。

表1 利用Materials Studio构建的臭氧-有机物-水模型体系的组成和其它细节Table 1 Details of ozone-organic-water models constructed using Materials Studio.

图1 利用Materials Studio构建的含7500个原子的对硝基苯酚体系I(包括300个对硝基苯酚、500个臭氧、500个水)在300 K下几何优化后的3D初始构象Fig.1 Optimized ozone-p-nitrophenol-water model containing 7500 atoms (model I, 300 p-nitrophenol,500 ozone and 500 H2O) at 300 K.

3 结果与讨论

3.1 主要分子数量的演变趋势

臭氧降解对硝基苯酚的终产物为水和二氧化碳，而芳香环开环是降解对硝基苯酚过程中的一个重要步骤，因此本文将芳香环和终产物 CO2的数量作为评价臭氧降解对硝基苯酚的两个指标。在模拟开始的极短时间内臭氧与水快速反应生成大量自由基，这些自由基会促使酚环活化，使其形成携带自由基的芳香环(为简略描述，下文统称其为六元环)。图2(a)中列出了采用ReaxFF MD模拟臭氧降解四个对硝基苯酚体系过程中六元环结构的数量随时间的演化趋势。对比不含羟基、不同臭氧含量的体系(I和II)发现，初始体系中臭氧数量的增加，可加快六元环开环速率。对比含1000个臭氧的两个体系(II和 III)发现，在初始体系中额外加入一定数量的羟基自由基，也可加快六元环开环速率。为了考察开环速率，利用一级反应动力学对各体系六元环的数量演变趋势进行拟合，其中含臭氧的体系相关系数均为0.98−0.99，遵循伪一级反应动力学规律。图 2(b)是模拟四个体系过程中得到的 CO2数量随时间的演变趋势，可看出未加入臭氧、仅含 250个羟基自由基的体系(IV)自始至终均没有生成CO2，而加入臭氧的三个体系在24−36 ps的阶段开始生成CO2分子，其数量在36−200 ps的阶段内随时间稳步上升，且开环数量最多的体系中生成的 CO2分子数量也最多。利用一级反应动力学对此阶段 CO2的数量演变趋势进行拟合，得到相关系数为0.98−0.99。

通过图2中六元环和CO2的数量演变可知，两者数量均呈现较好的线性关系，说明利用ReaxFF MD模拟臭氧降解水中对硝基苯酚体系可定性描述其遵循伪一级反应动力学规律，该现象与杨德敏等44通过实验得出的臭氧降解焦化废水的生化出水化学需氧量(COD)符合表观一级反应动力学的结论一致。

图2 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚(p-NP)体系在300 K条件下主要分子的数量随时间的演变Fig.2 Evolution tendency of dominant molecules in ozonation of p-NP in water obtained from NVT ReaxFF MD simulation at 300 K by using VARxMD.(a) aromatic-ring; (b) CO2.

图3 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚(p-NP)体系在300 K条件下初始物质及团簇的数量随时间的演变Fig.3 Evolution tendency of original components and clusters in ozonation of p-NP in water obtained from NVT ReaxFF MD simulation at 300 K by using VARxMD.(a) O3; (b) H2O; (c) ·OH; (d) cluster.

图 3给出了模拟获得的四个对硝基苯酚体系中初始物质臭氧、水、羟基及水团簇数量随时间的演化趋势。如图3(a)所示，在仅含500个和1000个臭氧的体系(I和II)中，臭氧数量均在36 ps内快速下降了超过 70%，在 0−12 ps的下降速度最快；其后一直保持下降趋势，虽然速度趋缓，但在200 ps时基本消耗殆尽。与此同时，如图3(b)所示，体系I和II中的H2O的数量在0−12 ps阶段快速下降至 108个(I)或 42个(含更多臭氧的体系II)；其后，体系I、II中H2O的数量均在36 ps时下降到最低点(I中79个、II中22个)，然后开始缓慢上升。图 3(c)中，体系 I中的·OH数量在0−12 ps阶段迅速上升至150个，体系II的·OH数量也上升至 182个的现象表明，臭氧的快速减少一部分生成了羟基自由基。图 3(d)所给出的各体系中观察到的团簇(包括 H2O-·OH、·OH-·OH 以及少量的H2O-H2O)数量随时间演变趋势表明，体系I的团簇数量在0−12 ps内快速上升至175个，体系II团簇数量也同时快速上升至194个，说明臭氧与大部分水快速反应，除了生成·OH自由基，一部分水可与·OH、或水与水分子结合成为团簇状态，解释了图3 (b)中水快速下降的现象。当初始体系中臭氧数量增加，水分子数量下降幅度更大，所生成的游离·OH以及团簇数量增多。说明臭氧浓度的提升增加了体系中自由基及团簇的数量，可能以此促进对硝基苯酚的降解。

体系III因在含1000个臭氧的体系II中引入了250个·OH，其早期快速生成的·OH及团簇数量均比其他体系多。在不含臭氧的体系IV中，H2O及·OH的数量在模拟开始时均迅速下降，而团簇数量迅速上升，接着各物质数量在40−200 ps的阶段内保持稳定。这一观察说明，水体系中的 H2O及·OH均会相互结合形成团簇。此外，通过分析四个体系中水分子趋势图还可发现体系I、II、III中H2O数量在36−200 ps阶段呈现缓慢上升状态，而体系IV中H2O数量在此阶段并未出现整体上升的现象。这是因为体系I、II、III在36 ps时已经有不同数量的对硝基苯酚分子被完全降解，生成CO2(图2(b))和H2O，而体系IV中并未出现对硝基苯酚完全降解的现象。

通过各体系中臭氧、水分子、·OH及团簇数量的演变可知，臭氧与水分子反应会消耗大量水分子；此外，体系中不仅存在游离的·OH，还有大量·OH与水分子以团簇状结构存在于水体系中。对硝基苯酚体系中六元环开环数量和·OH、团簇数量呈正相关关系，臭氧的存在可显著加快开环速率以及 CO2生成速率，说明在构建体系时加入的臭氧在降解对硝基苯酚时起主要作用，·OH自由基和团簇可进一步提高中间产物降解为水和二氧化碳的效率。

3.2 臭氧降解对硝基苯酚中自由基的作用

臭氧降解对硝基苯酚时，臭氧遇水易分解生成大量自由基，因此自由基驱动的反应占主导45。通过 ReaxFF MD的模拟发现，绝大多数自由基为·OH、·O2、·O。考虑到体系可能呈现出多样化的反应途径、且有利于得到规律性结论，下文选取体系III，即含300个对硝基苯酚、1000个臭氧、250个羟基、500个水的体系作为考察重点。

3.2.1 小分子自由基

图4 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K条件下与·OH有关的化学反应数量随时间的演变Fig.4 Evolution trend of ·OH radical involved reactions in ozonation of p-nitrophenol in water obtained from NVT ReaxFF MD simulation at 300 K by using VARxMD (model III).

为了考察臭氧与·OH的协同作用，在构建p-NP体系III时除了1000个臭氧，额外加入了250个·OH自由基。模拟结果表明，该体系中六元环开环速率及 CO2生成速率均最快。并且从 3.1节的分析可知，臭氧与水反应可生成大量·OH自由基，体系III中·OH数量虽然多一些，但其演化并未发生太大变化(图3(c))，在模拟过程中维持一个较稳定的状态((220 ± 30)个)。为了解释这种现象，对体系中·OH相关的化学反应演化进行了观察，发现在不同的时刻均有大量·OH生成和消耗。为了大致得出·OH的生成情况，每隔6 ps提取一次反应列表，并对总反应列表中生成·OH的反应进行筛选，共得到2936个反应，其中由O3作为反应物参与的反应有264个。图4给出不同时刻相关反应的数量，可看出在反应的起始时刻，即0−6 ps内几乎所有·OH的生成反应都有O3参与；之后，由 O3参与生成·OH的化学反应迅速变少。由 O3参与生成·OH的化学反应主要集中在 0−54 ps阶段内；在 54 ps后的阶段内大部分·OH的生成与O3关系不大。借助VARxMD分析54 ps后的化学反应，发现体系中存在大量团簇型自由基，这些自由基主要由·OH、·O 等与水分子通过氢键结合形成，在所模拟的反应过程中，这些团簇型自由基也会发生内部断键从而生成小分子自由基。该现象说明在反应起始阶段，臭氧与水反应生成大量·OH，同时有几乎同等数量的游离·OH参与形成团簇，这些团簇在之后模拟过程中会断键形成新的游离·OH自由基，从而使体系中·OH数量维持在一个较稳定的范围，见图3(c)。

臭氧分解会生成大量的自由基，这些自由基在体系模拟中会相互转化，如表 2所示。臭氧分解生成的·O2自由基演化如图5(a)所示，在0−96 ps的阶段内，体系中·O2自由基的数量生成趋势与图3(a)中臭氧消耗规律基本一致：在 0−12 ps内·O2的数量从0迅速升至572个，在12−96 ps内缓慢上升至815个；在96 ps后开始呈现出缓慢下降的趋势。利用 VARxMD详细观察体系 p-NP中·O2的行为后发现：·O2主要参与 CO2的生成，通过与·CO反应生成·CO3自由基，接着O―O断裂生成 CO2分子。如图 5(a)所示，·CO3在 36 ps时首次出现，几乎与图 2(b)中 CO2分子在同一时间生成；而在 43−200 ps阶段，·CO3自由基的数量起伏不定，这是因为·CO3自由基是CO2的主要前驱体，且·CO3→ CO2反应迅速，致使·CO3存在时间极短，因此各时间步出现的·CO3数量并没有呈现出太强的规律性。

表2 利用VARxMD观察ReaxFF MD模拟结果获得的在300 K条件下臭氧降解水中对硝基苯酚过程中主要自由基的生成与消耗反应Table 2 Generation and consumption reactions of major radicals in ozonation of p-nitrophenol obtained by VARxMD in the 300 K-NVT ReaxFF MD simulations.

图5 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K条件下主要自由基的数量随时间的演变Fig.5 Evolution of dominant radicals in ozonation of p-nitrophenol in water (model III) obtained with aid of VARxMD in the 300 K-NVT ReaxFF MD simulation.(a) ·O2, CO3 radical; (b) ·O.

臭氧在水体系中除了生成·O2，还会生成少量的·O，同样有助于其降解过程，如图 5(b)所示。在前40 ps的模拟阶段，体系p-NP中的·O的数量迅速上升，随后·O的数量稳定为约 80个。借助VARxMD追踪·O的消耗途径，发现·O主要有四种消耗途径：(1) ·O通过氢键作用与水分子团簇结合，使水团簇变为具有自由基性质的团簇，对六元环的开环有一定的帮助作用；(2) ·O相互结合，形成稳定的氧气分子存在于密闭的模拟体系中；(3) ·O 之间相互结合形成具有进攻性的·O2自由基，一定程度上增加了臭氧直接分解的·O2数量，该反应途径可能是·O2与·O 数量上有所差异的原因之一；(4) ·O 是生成·OH 自由基的一种中间体，·O与水分子结合形成·HOHO，接着断键生成两分子·OH自由基。

从上述分析可知，臭氧降解对硝基苯酚的过程与自由基·OH、·O2、·O的生成和消耗都有一定的关系，·OH与其它自由基或水分子形成团簇型自由基活化了整个体系，·O2主要作用于CO2的生成，·O通过上述四个途径对帮助降解酚分子发挥一定的作用。

3.2.2 团簇型自由基

·OH的氧化能力强于臭氧46，在臭氧深度氧化酚类分子时起到至关重要的作用14,45。但是，通过本文对小分子自由基(·OH、·O等)的反应轨迹分析发现，在300 K恒温条件下模拟臭氧降解对硝基苯酚体系中，·OH、·O 之间或者与水分子间更易以氢键相连，缔结成大小不一的团簇型自由基，这些团簇型自由基的形成同样需要消耗一定量的水分子，这应该是体系中水分子数量在模拟开始就迅速下降的一个原因(如图3)。团簇型自由基具有和·OH、·O 等小分子自由基相似的反应活性，会攻击酚分子的侧链或六元环，发生加成反应或攫氢反应，从而促进酚的开环和降解。这一反应过程会加速六元环的开环和 CO2的生成，因此，团簇型自由基在各个体系数量的差异有可能是造成六元环的开环和 CO2的生成趋势有所差别的原因之一(如图2)。

为了进行区分，本文将含 7个及以上氢原子的HmOn团簇定义为大团簇，将含6个及以下氢原子的HmOn团簇定为小团簇，图6是ReaxFF MD模拟得出的对硝基苯酚体系(p-NP)中大团簇和小团簇数量随时间的演化趋势。从图6中可以看出，反应开始时的0−10 ps阶段，体系p-NP中大团簇、小团簇数量均迅速升高，小团簇的数量更多，利用VRAxMD观察可知，体系中最大团簇含H原子达到19个；在10−36 ps内大团簇数量继续上升，而小团簇数量有少量下降。说明在模拟开始时，小分子自由基与水会随机碰撞结合形成大小不一的自由基团簇，小团簇同样会相互碰撞形成大团簇。在48−84 ps阶段内大团簇数量出现大幅度波动，对应着O3参与生成·OH及其他自由基的化学反应的结束，说明大团簇的形成与O3密切相关。与此同时小团簇保持数量增多，猜测可能是O3生成的各类自由基促进大团簇断键形成小团簇导致。84 ps之后大团簇数量缓慢下降，小团簇数量继续缓慢上升。

图6 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K条件下的团簇分子数量随时间的演变Fig.6 Evolution trend of water-clusters in p-nitrophenol model obtained from NVT ReaxFF MD simulation at 300 K by using VARxMD (model III).(a) large water-clusters (containing more than or equal to seven H atoms); (b) small water-clusters (containing less than or equal to six H atoms).

通过对比图6中大、小团簇数量的演化猜测，在模拟过程中团簇的数量发生变化与对硝基苯酚的降解过程有一定的关系。利用VARxMD详细观察体系中团簇的行为，发现在臭氧降解酚的过程中，·OH直接进攻六元环使其开环的现象并不常见，更多的是具有自由基性质的团簇对六元环的进攻。这些团簇分子中大部分包裹·OH或·O自由基，·OH (·O)与水分子、·OH (·O)与·OH (·O)间通过氢键结合成大小不一的(·OH)x·(H2O)y形式的团簇分子，具有与自由基相似的攻击性。与小团簇相比，大团簇所占的空间范围更大，相对有更多的机会与酚分子发生碰撞。大团簇攻击六元环使其支链增长形成大分子，在密闭盒子中可能会因空间构型的原因使大团簇支链断裂成小团簇，因此会出现大团簇和小团簇数量交替上升或下降的现象。该模拟结果说明水体系中可能会存在团簇型自由基，更易进攻酚类污染物发生开环反应。

3.3 对硝基苯酚在ReaxFF MD模拟中降解机理分析

为了在可接受的模拟运行时间内让大部分对硝基苯酚的降解反应发生，模拟中采用的时间步长为0.4 fs，一定程度上加速了氧化反应过程，也有可能掩盖了一些反应的细节。借助VARxMD发现，臭氧降解对硝基苯酚的过程中，六元环开环后大多直接断裂为小分子片段，并没有经过多次开环闭环的振荡反应。图 7为模拟体系中含不同碳原子数的分子片断的数量演化。

图7 利用VARxMD得到的NVT ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K条件下含C片断数量(C1指片断分子中含1个碳)随时间的演变Fig.7 Evolution trend of carbon fragments in ozonation of p-nitrophenol in water (model III)obtained by VARxMD from the 300 K-NVT ReaxFF MD simulation (C1 represents fragment containing only one carbon).

表3 通过VARxMD得到ReaxFF MD模拟臭氧降解对硝基苯酚体系在300 K-NVT的三个主要反应阶段示例Table 3 Stages of p-nitrophenol degradation observed in 300 K-NVT ReaxFF MD simulation with the aid of VARxMD.

水环境中的酚分子间可能会因氢键作用与水分子等瞬时连结形成更大的分子团，由于未发生开环反应，图7中将此与C6分子归为一类，表示为“＞ C6”。如图7所示，＞ C6物质数量随着模拟时间的增加从300个缓慢减少至200左右，表示越来越多的六元环被臭氧降解为小分子片段。体系中 C1物质主要以·CO、·CO3、CO 及 CO2的形式存在，其分子数量从降解过程的30 ps开始迅速上升，生成速率远快于其它的 C2−C5分子片断。通过VARxMD分析发现，C2物质主要由六元环开环断裂生成，它的数量同样随着时间缓慢上升。C2物质在降解过程中迅速转化为·CO自由基和 CO2分子，使得 C1物质的数量快速上升。C4和 C5片段数量都比较少，仅有不足10个，生成趋势不明显。不同含碳片段数量的演化说明六元环开环形成链状后更易于发生C−C均裂，形成C2、C3中间体，很少出现六元碳链中伯碳与仲碳断键形成C1、C5中间体的现象。利用VARxMD可显示两个输出时间步之间的反应物及生成物 2D、3D结构这一功能，可进一步跟踪观察到有机物在每一时间步的氧化情况，表 3给出对硝基苯酚在水环境下被臭氧降解过程的反应示例，可分为羟基攫氢、酚环开环、链状结构断裂分解三个阶段：(1)酚环或者酚羟基上的攫氢反应，大多是由羟基或者含羟基的自由基团簇进攻酚环引发；自由基团簇中的氧与酚羟基中的氢形成氢键，导致酚羟基脱氢形成自由基；接着另一分子羟基或自由基团簇加成至羟基的α碳位，并继续发生六元环的脱氢反应，形成苯醌结构。模拟观察到的攫氢反应机理与通常认为的臭氧分子亲电加成到苯环的机理有所差别，可能的原因是模拟时间步长较大(0.4 fs)，缺失了部分的反应细节，也有可能是ReaxFF力场对该机理的描述不够准确，需要进一步的工作加以证实。(2)羟基或自由基团簇对六元环的进攻，活化了六元环，其中，这些自由基更倾向于先进攻酚羟基的邻位，这可由图 8所示的开环过程中邻位中间体分子数量远大于间位中间体分子的观察所佐证；当环上 C＝O双键或 C―O·单键自由基等结构越来越多，六元环的空间构型已发生很大形变，接着发生开环反应，大部分会形成链状分子；小部分开环后会缩为五元环，经较长时间后转化为小分子。(3)开环后的链状分子含羰基、硝基等官能团，经较长时间的振荡反应后断键形成小分子结构。

图8 利用VARxMD得到的ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K条件下酚羟基邻、间位被进攻的中间体分子数量随时间的演变Fig.8 Evolution trend of intermediates representing the p-nitrophenol’ reaction sites attacked by radicals in ozonation of p-nitrophenol in water obtained by VARxMD from the 300 K-NVT ReaxFF MD simulation of model III.

图9 利用VARxMD得到ReaxFF MD模拟臭氧降解水中对硝基苯酚体系(III)在300 K-NVT降解条件下重要中间体片段随时间的演变Fig.9 Evolution trend of dominant intermediates in the ozonation of p-nitrophenol in water obtained by using VARxMD from the 300 K-NVT ReaxFF MD simulation of model III.(a) key intermediates; (b) C＝O moieties.

上述所观察到的降解机理与图 9中主要中间体的数量演化趋势一致。如图 9(a)所示，六元碳链的数量在0−70 ps处于上升阶段，在70−150 ps处于动态稳定阶段，在这一阶段内既有链状结构的生成也有碳链的分解，但须经过较长时间震荡才会最终分解为小分子片断，接着C6碳链分子有缓慢下降的趋势；而五元环结构在40 ps时首次出现，其数量仅不足 5个，在随后的模拟时间内处于稳定状态。该现象说明大多数六元环开环后形成了碳链。C6碳链断裂后形成的中间体多数转化为 OCCO自由基，并分解为·CO，模拟体系中存在大量的·O2提高了·CO 与其碰撞的概率，加快了·CO 向 CO2转变的过程。CO2的前驱体主要为·CO和含羧酸结构片段，并且·CO居多(体系III中发现仅有2个含羧酸结构片段)。据文献报道在臭氧降解酚类分子的过程中会出现较多羧酸类中间体，但在模拟过程中并未发现太多此类中间体，这可能与体系中过多臭氧及羟基自由基有关。臭氧及羟基自由基的大量存在很容易使―COOH官能团中的氢或羟基脱落，形成H2O或·HO2、·HO3等自由基，导致―C＝O结构数量增加，体系 III中―C＝O结构(CO、CO2分子除外)数量演变如图9(b)所示，其数量在0−12 ps阶段迅速上升到264个，这是因酚羟基被体系中活性自由基攫氢生成酮所致；在24−200 ps阶段缓慢上升，说明自由基进攻六元环后最终形成大量―C＝O结构。

4 结 论

本文选取对硝基苯酚作为焦化废水生物段出水的目标污染物，通过构建较大规模的模型，利用基于反应力场的ReaxFF MD模拟了对硝基苯酚在水环境下的降解反应。结合VARxMD分析工具，得到了体系中对硝基苯酚降解初始过程中重要自由基、中间体和产物随时间的演变趋势。

分析ReaxFF MD模拟结果的轨迹文件表明：在臭氧及羟基数量充足的水环境中，模拟可定性描述臭氧降解对硝基苯酚的初始反应行为，对硝基苯酚的降解与二氧化碳的生成遵循伪一级反应动力学规律。自由基的演化趋势表明，O3在水环境中可快速分解形成·O2、·OH、·O等自由基，·O2在降解末期可加快中间产物向CO2转化，提高降解速率，·O在对硝基苯酚的降解过程中以四种转化路径发挥其作用。此外，·OH、·O2、·O在水环境中易形成水-自由基的团簇，隐藏在团簇状结构中的自由基使团簇有了攻击性，可帮助六元环中C形成―C＝O的结构，利于开环。借助VARxMD分析对硝基苯酚的降解路径发现，酚类有机物在水环境下被臭氧降解的过程主要分为三步，攫氢、六元环开环、碳链的氧化分解。

值得注意的是，在水体系中可检测出活性氧(·O2)的主要存在形式有两种：单线态分子氧(1O2)和超氧自由基(·O2−)，由于两者电子排布不同使1O2与·O2−具有不同的化学活性，在臭氧催化氧化水中酚类有毒有机污染物时发挥了不同的作用47,48，由于模拟方法的限制，本文尚不能从·O2中严格区分出1O2和·O2−。

(1) Zhou, S. Y.; Watanabe, H.; Wei, C.; Wang, D. Z.; Zhou, J. T.;Tatarazako, N.; Masunaga, S.; Zhang, Y. Ecotox. Ecotoxicol.Environ. Safe. 2015, 115, 217. doi: 10.1016/j.ecoenv.2015.02.017

(2) Jin, X. W.; Li, E. C.; Lu, S. G.; Qiu, Z. F.; Sui, Q. J. Environ.Sci-China 2013, 25 (8), 1565. doi: 10.1016/S1001-0742(12)60212-5

(3) Zhu, Q.; Liu, Z. W.; Wang, J. R. Nonferrous Metal Science and Engineering 2011, 2 (5), 93. [朱 强, 刘祖文, 王建如. 有色金属科学与工程, 2011, 2 (5), 93.] doi:10.13264/j.cnki.ysjskx.2011.05.003

(4) Wang, F.; Hu, Y. R.; Guo, C.; Huang, W.; Liu, C. Z. Bioresour.Technol. 2012, 110, 120. doi: 10.1016/j.biortech.2012.01.184

(5) Hall, D. A.; Nellist, G. R. J. Chem. Technol. Biotechnol. 1959, 9 (11),565. doi: 10.1002/jctb.5010091102

(6) Yang, P.; Wang, B. Environ. Protec. Chem. Ind. 2001, 21 (3), 144.[杨 平, 王 彬. 化工环保, 2001, 21 (3), 144.]doi: 10.3969/j.issn.1006-1878.2001.03.005

(7) Morper, M.; Jell, A. Linde Reports Sci. Technol. 2000, 62, 20.

(8) Han, T.; Chen, Z. S.; Lin, C.; Wen, Z. W.; Wei, C. H. Acta Sci.Circumst. 2016, 36 (1), 149. [韩 涛, 陈梓晟, 林 冲, 文泽伟,韦朝海. 环境科学学报, 2016, 36 (1), 149.]doi: 10.13671/j.hjkxxb.2015.0480

(9) Ren, Y.; Wei, C. H.; Wu, C. F.; Li, G. B. Acta Sci. Circumst. 2007,27 (7), 1094. [任 源, 韦朝海, 吴超飞, 李国保. 环境科学学报,2007, 27 (7), 1094.] doi: 10.13671/j.hjkxxb.2007.07.004

(10) Zhang, H.; Yang, Z. R.; Chen, H. Q. Environmental Protection 1995,No. 6, 11. [张 晖, 杨卓如, 陈焕钦. 环境保护, 1995, No. 6, 11.]doi: 10.14026/j.cnki.0253-9705.1995.06.004

(11) Chang, E. -E.; Hsing, H. J.; Chiang, P. C.; Chen, M. Y.; Shyng, J. Y.J. Hazard. Mater. 2008, 156 (1−3), 560.doi: 10.1016/j.jhazmat.2007.12.106

(12) Utsumi, H.; Hakoda, M.; Shimbara, S.; Nagaoka, H.; Chung, Y. S.;Hamada, A. Water Sci. Technol. 1994, 30 (9), 91.

(13) Wang, H. R.; Wang, S.; Li, Y. Q.; Wang, F. Y.; Li, Y. K.; Yin, J.Chinese Journal of Disinfection. 2009, 26 (5), 481. [王华然, 王 尚,李昀桥, 王福玉, 李迎凯, 尹 静. 中国消毒学杂志, 2009, 26 (5),481.]

(14) Sun, H. F. The treatment of Organic Wastewater by Hydroxyl Radical Activated Oxygen. Ph. D. Dissertation, Beijing University of Chemical Technology, Beijing, 2010. [孙华峰. 羟基自由基活性氧处理有机废水的研究[D]. 北京: 北京化工大学, 2010.]

(15) Roche, P.; Prados, M. Ozone Sci. Eng. 1995, 17 (6), 657.doi: 10.1080/01919512.1995.10555777

(16) Rosenfeldt, E. J.; Linden, K. G.; Canonica, S.; Gunten, V. U. Water Res. 2006, 40 (20), 3695. doi: 10.1016/j.watres.2006.09.008

(17) Wei, Q.; Qiao, S. F.; Sun, B. C.; Zou, H. K.; Chen, J. F.; Shao, L.RSC Adv. 2015, 5 (113), 93386. doi: 10.1039/C5RA14198B

(18) Tanaka, K.; Abe, K.; Hisanaga, T. J. Photochem. Photobiol. A 1996,101 (1), 85. doi: 10.1016/S1010-6030(96)04393-6

(19) Zazo, J. A.; Casas, J. A.; Mohedano, A. F.; Gilarranz, M. A.;Rodriguez, J. J. Environ. Sci. Technol. 2005, 39 (23), 9295. doi:10.1021/es050452h

(20) Valsania, M. C.; Fasano, F.; Richardson, S. D.; Vincenti, M. Water Res. 2012, 46 (8), 2795. doi: 10.1016/j.watres.2012.02.040

(21) van Duin, A. C. T.; Dasgupta, S.; Lorant, F.; Goddard, W. A. J. Phys.Chem. A 2001, 105 (41), 9396. doi: 10.1021/jp004368u

(22) Mortier, W. J.; Ghosh, S. K.; Shankar, S. J. Am. Chem. Soc. 1986,108 (15), 4315. doi: 10.1021/ja00275a013

(23) Mattsson, T. R.; Lane, J. M. D.; Cochrane, K. R.; Desjarlais, M. P.;Thompson, A. P.; Pierce, F.; Grest, G. S. Phys. Rev. 2010, 81 (5),054103. doi: 10.1103/PhysRevB.81.054103

(24) Chenoweth, K.; van Duin, A. C. T.; Goddard, W. A. J. Phys. Chem.A 2008, 112 (5), 1040. doi: 10.1021/jp709896w

(25) Liu, X. L.; Li, X. X.; Han, S.; Qiao, X. J.; Zhong, B. J.; Guo, L. Acta Phys. -Chim. Sin. 2016, 32 (6), 1424. [刘晓龙, 李晓霞, 韩 嵩, 乔显杰, 钟北京, 郭 力. 物理化学学报, 2016, 32 (6), 1424.]doi: 10.3866/PKU.WHXB201603233

(26) Zheng, M., Wang, Z., Li, X. X.; Qiao, X. J.; Song, W. L.; Guo, L.Fuel 2016, 177, 130. doi: 10.1016/j.fuel.2016.03.008

(27) Zhang, T. T.; Li, X. X.; Qiao, X. J.; Zheng, M.; Guo, L.; Song, W. L.;Lin, W. G. Energy and Fuel 2016, 30 (4), 3140.doi: 10.1021/acs.energyfuels.6b00247

(28) Castro-Marcano, F.; van Duin, A. C. T. Combustion and Flame 2013,160 (4), 766. doi: 10.1016/j.combustflame.2012.12.007

(29) Rahaman, O.; van Duin, A. C. T.; Goddard, W. A., III.; Doren, D. J.J. Phys. Chem. B 2011, 115 (2), 249. doi: 10.1021/jp108642r

(30) Jiang, D. D.; Zhang, J. L.; Wei, L.; Zhang, M.; Han, Y. Coking Wastewater Treatment via Supercritical Water: Reaxff Reactive Molecular Dynamics Simulation. In Engineering Science and Technology, Annual Conference of Chemical Industry and Engineering Society of China, Beijing, Oct 16–18, 2015; Chemical Industry and Engineering Society of China, Ed.; Chemical Industry and Engineering Society of China: Beijing, 2015; 414. [姜丹丹, 张金利, 李 韡, 张 铭, 韩 优. 基于ReaxFF力场对超临界水处理焦化废水反应机理的研究. In 工程科技, 2015年中国化工学会年会, 北京, Oct 16–18, 2015; 中国化工学会编; 北京: 中国化工学会, 2015: 414.]

(31) Zhang, J. L.; Gu, J. T.; Han, Y.; Li, W.; Gan, Z. X.; Gu, J. J. J. Mol.Model. 2015, 21 (3), 54. doi: 10.1007/s00894-015-2603-7

(32) Zheng, M.; Li, X. X.; Guo, L. J. Mol. Graph. Model. 2013, 41, 1.doi: 10.1016/j.jmgm.2013.02.001

(33) Liu, J.; Li, X. X.; Guo, L.; Zheng, M.; Han, J. Y.; Yuan, X. L.; Nie, F.G.; Liu, X. L. J. Mol. Graph. Model. 2014, 53 (9), 13.doi: 10.1016/j.jmgm.2014.07.002

(34) Han, J. Y.; Li, X. X.; Guo, L.; Zheng, M.; Qiao, X. J.; Liu, X. L.;Gao, M. J.; Zhang, T. T.; Han, S. Comp. & App. Chem. 2015, 32 (5),519. [韩君易, 李晓霞, 郭 力, 郑 默, 乔显杰, 刘晓龙, 高明杰,张婷婷, 韩 嵩. 计算机与应用化学, 2015, 32 (5), 519.]doi: 10.11719/com.app.chem20150502

(35) Li, X. X.; Zheng, M.; Han, J. Y.; Guo, L.; Liu, X. L.; Qiao, X. J.;Yuan, X. L. Scientia Sinica Chimica 2015, 45 (4), 373. [李晓霞, 郑默, 韩君易, 郭 力, 刘晓龙, 乔显杰, 袁小龙. 中国科学: 化学,2015, 45 (4), 373.] doi: 10.1360/N032014-00281

(36) Shao, N.; Wang, B. H.; Yuan, D. D. Contemp. Chem. Ind. 2015, 44(11), 2615. [邵 楠, 王宝辉, 苑丹丹. 当代化工, 2015, 44 (11),2615.] doi: 10.13840/j.cnki.cn21-1457/tq.2015.11.033

(37) Zhang, W. H.; Wei, C. H. Environmental Protection of Chemical Industry 2015, 35 (3), 273. [张万辉, 韦朝海. 化工环保, 2015, 35(3), 273.] doi: 10.3969/j.issn.1006-1878.2015.03.010

(38) Liu, X. Q.; Li, G. B.; Wu, H. Z.; Zhang, W. H.; Guan, Q. Q.; Feng, C.H.; Wu, C. F.; Hu, Y.; Wei, C. H. Environ. Chem. 2012, 31 (10),1487. [刘显清, 李国保, 吴海珍, 张万辉, 关清卿, 冯春华, 吴超飞, 胡 芸, 韦朝海. 环境化学, 2012, 31 (10), 1487.]

(39) Li, J. F.; Wu, J.; Sun, H. F.; Cheng, F. Q.; Liu, Y. Desalination 2016,380 (15), 43. doi: 10.1016/j.desal.2015.11.020

(40) Zhang, S. H.; Zheng, J.; Chen, Z. Q. Sep. Purif. Technol. 2014, 130(20), 610. doi: 10.1016/j.seppur.2014.06.019

(41) http://accelrys.com/products/collaborative-science/bioviamaterials-st udio (accessed Jan 18, 2017).

(42) Mayo, S. L.; Olafson, B. D.; Goddard, W. A. J. Phys. Chem. 1990,94 (26), 8897. doi: 10.1021/j100389a010

(43) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola,A.; Haak, J. R. J. Chem. Phys. 1984, 81 (8), 3684. doi:http://dx.doi.org/10.1063/1.448118

(44) Yang, D. M.; Wang, B.; Yuan, J. M. Acta Petrolei Sinica (Petroleum Processing Section) 2012, 28 (4), 683. [杨德敏, 王 兵, 袁建梅.石油学报(石油化工), 2012, 28 (4), 683.] doi:10.3969/j.issn.1001-8719.2012.04.026

(45) Gharbani, P.; Khosravi, M.; Tabatabaii, S. M.; Zare, K.; Dastmalchi,S.; Mehrizad, A. Int. J. Environ. Sci. Tech. 2010, 7 (2), 377.doi: 10.1007/BF03326147

(46) Tang, Q. L. Ozonation of Various Organic Substances. Ph. D.Dissertation, Shanghai Jiao Tong University, Shanghai, 2013. [唐庆丽. 臭氧氧化降解多种有机污染物的研究[D]. 上海: 上海交通大学，2013.]

(47) Wang, Y. X.; Xie, Y. B.; Sun, H. Q.; Xiao, J. D.; Cao, H. B.; Wang,S. B. Catal. Sci. Technol. 2016, 6 (9), 2918.doi: 10.1039/C5CY01967B

(48) Wang, Y. X.; Xie, Y. B.; Sun, H. Q.; Xiao, J. D.; Cao, H. B.;Wang, S. B. ACS Appl. Mater. Interfaces 2016, 8 (15), 9710.doi: 10.1021/acsami.6b011

Molecular Dynamics Simulation of Ozonation of p-Nitrophenol at Room Temperature with ReaxFF Force Field

WANG Zi-Min1,2ZHENG Mo2XIE Yong-Bing3LI Xiao-Xia2,4,*ZENG Ming1CAO Hong-Bin3,*GUO Li2,4
(1School of Chemical and Environmental Engineering, China University of Mining and Technology, Beijing, Beijing 100083, P. R.China;2State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Science, Beijing 100190, P. R. China;3Research Center for Process Pollution Control, National Engineering Laboratory for Hydrometallurgical Cleaner Production Technology, Institute of Process Engineering, Chinese Academy of Science, Beijing 100190,P. R. China;4School of Chemistry and Chemical Engineering, University of Chinese Academy of Science,Beijing 100049, P. R. China)

Understanding the reaction mechanism of phenol ozonation in coking wastewater is very important for industrial applications of the ozonation process. Ozonation of p-nitrophenol in water at 300 K was simulated by ReaxFF force field molecular dynamics (ReaxFF MD) employing the GPU-enabled high-performance code of GMD-Reax and a unique code of VARxMD developed in authors’ group.Evolution trends of aromatic ring opening, CO2generation, dominant radicals (·OH, ·O2, ·O), and H2O clusters were obtained. The simulated CO2generation and reduction of aromatic ring could be described with pseudo-first-order kinetics. Moreover, the reaction pathways for the ozonation of p-nitrophenol can bedivided into three stages: hydrogen abstraction, opening of the six-membered-ring structure, and the breaking of C−C bonds. The simulations revealed the important role of radicals and water clusters in the ozonation of p-nitrophenol. This work is an attempt to investigate the ozonation mechanism of phenols in aqueous solutions at room temperature using ReaxFF MD, which should be helpful for further experimental or theoretical investigation of the mechanism.

p-Nitrophenol; Ozonation; ReaxFF MD; Radical behavior; Reaction mechanism

January 26, 2017; Revised: April 6, 2017; Published online: April 13, 2017.

O643

10.3866/PKU.WHXB201704132 www.whxb.pku.edu.cn

*Corresponding authors. LI Xiao-Xia, Email: xxia@ipe.ac.cn; Tel: +86-10-82544944. CAO Hong-Bin, Email: hbcao@ipe.ac.cn; Tel: +86-10-82544844.

The project was supported by the National Natural Science Foundation of China (21373227), the National Science Fund for Distinguished Young Scholars,China (51425405), and China’s State Key Laboratory of Multiphase Complex Systems, China (COM2015A004).

国家自然科学基金(21373227), 国家杰出青年科学基金(51425405)和多相复杂系统国家重点实验室基金(COM2015A004)资助项目

© Editorial office of Acta Physico-Chimica Sinica