大豆绝缘油的氧化路径仿真分析研究

时间：2024-07-28

周年荣，何潇，郭新良，何运华，张林山，谭向宇，方正云

（1.云南电网有限责任公司电力科学研究院，云南昆明 650217；2.云南电网有限责任公司，云南昆明 650011）

0 引言

植物绝缘油作为新兴的一种环保型液体电介质，已应用在国内10 kV和35 kV配电变压器中，也开始在110 kV甚至更高电压等级的电力变压器中应用[1-5]。油纸绝缘作为变压器主要的绝缘电介质，直接关系到电力系统的经济性和稳定性。植物绝缘油在运行过程中能够吸收绝缘纸中水分，提高油纸绝缘使用寿命，相比传统矿物油纸绝缘具有更好的经济性[6-7]。植物绝缘油主要由脂肪酸和甘油三酯组成，甘油三酯中脂肪酸基团繁多，因而植物绝缘油的分子结构不同于矿物绝缘油，具有明显不对称性[8]。因此，植物绝缘油的极性高于矿物绝缘油，具有较高的电导率和介质损耗因数；且植物绝缘油的运动黏度较大，会影响变压器的散热性能。

国外有关学者对植物绝缘油的制取方法、理化性能、电气性能等方面进行了系统的研究[9]。1971年，K M KAMATH等[9]对花生油、蓖麻油和椰子油的绝缘性能进行了比较分析，并针对3种油的绝缘情况进行了适当的调整和优化，一年之后将它们投入使用于电力设备中。1987年，巴西报道了蓖麻绝缘油的研究进展[10]。但由于植物绝缘油理化特性的限制，植物绝缘油的工业替代研究暂被搁置。直到90年代，矿物绝缘油的问题日益出现，为了解决这一问题，人们把对植物绝缘油的研究提上日程。近10年，有关植物绝缘油的研究在国内逐渐展开，例如重庆大学、西安交通大学、河南电科院等单位的研究者针对植物绝缘油的各方面性能进行了研究[11]，提升改善了一系列理化和电气性能。目前已经有一部分植物绝缘油成功投入使用。2018年河南省电力公司和国家电网公司电力科学研究院成功实现植物绝缘油的量产，并且设计出与之相匹配的植物绝缘油变压器，通过实地使用测试，证明了其良好的安全性。

然而由于植物绝缘油分子结构上的缺点导致其在电气设备如变压器运行中易氧化老化，致使电气性能下降。且植物绝缘油的介质损耗也较高，运动黏度大，抗氧化能力弱。本文主要对大豆植物绝缘油FR3的氧化安定性进行分子动力学仿真研究。

1 大豆绝缘油的氧化路径

FR3是由美国Cooper公司开发出的一种植物绝缘油。因为FR3的组成成分实质上基本都是脂肪酸，所以FR3的氧化也是油脂的氧化。油脂的氧化过程十分复杂，有学者指出，首先碳碳双键氧化得到氧，接着生成氢过氧化物，通过实验发现，双键的α-亚甲基碰到氧会产生氢过氧化物。随着研究的深入可知，油脂氧化满足自由基链式反应理论。图1为典型的脂质自由基链式反应过程图，RH是不饱和脂，R中有碳碳双键，X是反应后产生的各种自由基[12]。

图1 典型的脂质自由基链式反应Fig.1 Typical chain reaction of lipid free radical

脂肪酸的氧化分为很多种，例如光氧化、自动氧化等。自动氧化的氧化环境为常温下无光照、无催化剂，这时脂肪酸与氧气自主发生氧化反应。周华龙等[13]研究指出，这一自主氧化反应发生在双键连着的α-CH2上。油脂的自动氧化包括链引发、链增长、链终止阶段，氢受到氧分子的攻击产生游离基，此时氧分子同游离基结合为氧游离基，不饱和脂肪酸会越来越少，随着反应时间的逐渐增加，反应逐渐结束，产生的游离基最后会结合在一起形成稳定化合物。

FR3在高温环境下的反应与常温下有着明显的区别，其氧化反应速度会大幅加快，这是因为高温下氧气在反应系统中的溶解度减少了。当氧气处于较高的温度环境时，其有充分的能量去破坏C-C或C-H共价键，形成烃基自由基，紧接着进行自由基链式氧化。当所处温度低于150℃时，烃基自由基的形成大于自由基的降解；当所处温度高于150℃时，烃基自由基的形成小于自由基降解[14]，此时脂质发生氧化反应和热反应。

在高温条件下发生的氧化反应统称为高温氧化反应，根据反应结束后物质的分解和聚合情况可以将高温氧化反应分为氧化分解反应和氧化聚合反应两种类型，高温氧化反应的生成物主要包括氧化甘油三酯、氧化甘油三酯聚合物以及短碳链的甘油三酯及挥发性化合物3种类型，这3种类型氧化产物的生成方式各不相同，但是所用到的化合物种类却较为相似。①当甘油三酯自由基与氢在高温下发生氧化反应或者与羟基在高温下发生氧化反应时就会生成氧化甘油三酯；②当多个甘油三酯自由基R·之间发生相互作用时，其中某一氧化甘油三脂上的就会与另一甘油三酯自由基上的烷氧自由基发生氧化聚合反应，生成氧化甘油三酯聚合物；③短链脂肪酸酯是一种重要的氧化反应物质，可以通过烷氧自由基断裂的方式形成，而甘油脂肪酸链在氧化反应过程中会出现氧化分解现象，形成醛类、醇类、酸类、酮类以及烃类物质等小分子化合产物，统称为短碳链的甘油三酯及挥发性化合物。

2 模拟与计算方法及分子动力学模拟过程

2.1 基于ReaxFF力场的分子模拟

由于键解离能只能反映单个分子键断裂的容易程度，无法反映复杂的热氧化过程的真实情况。而ReaxFF力场使用键序的概念来对键的断裂与形成进行精确描述，可以达到模拟反应过程的目的。因此本研究采用ReaxFF力场构建分子模拟体系进行热氧化过程的模拟，其力场常数基于密度泛函方法获得，使用CHO体系参数[15]。

ReaxFF力场的系统总能量的表达式为式（1）。

式（1）中：EvdWa为范德华力；ECo为非键相互库伦作用力；Eval为键角能量；Epen为键角能量惩罚项；EH-bond为氢键能量项；Elp为孤对电子能量项；Ebond为原子间键能；Eover和Eunder为过配位能量校正；Etors和Econj为四体相互作用能量项。

基于上述ReaxFF力场的理论开展热氧化环境下的分子动力学模拟，其运行步骤如下：构建分子模拟体系，设置模拟时间、模拟温度和平衡时间等参数；依据ti时每个原子对应的坐标（Xti,Yti,Zti），利用ReaxFF力场可以获得任意两个原子（i与j）间的距离rij；获得ti时任意原子间（i与j）的键级，根据文献[16]判断产物是否生成。计算ti时分子动力学模型体系内所有产物的分子个数，计算获得整个系统的能量Esystem；如果ti+1时刻获得的结果与ti时刻相同，则整个动力学反应完成，如若不同，则继续迭代反应直至整个模拟结束。

2.2 键解离能计算方法

本研究使用基于密度泛函理论的量子化学计算方法对构建的三元醇分子模型进行结构优化与频率计算，获得分子结构的热力学信息。为确保计算的准确性,采用B3LYP/6-31G*基组进行分子计算，其中B3LYP方法是目前应用广泛的密度泛函方法[17]。

分子的键解离能（bond-dissociation energy，BDE）可以表示为式（2）。

式（2）中：H(A·)和H(B·)为断键后原子的能量；H(AB)为原子结合的能量。

目前相关交换泛函形式中，B3LYP方法是使用最为广泛的密度泛函方法，对于键解离能的相关计算，由于不需要考虑非键相互作用体系，在计算的过程中不需要添加弥散函数，目前相关基组如6-311G**或TZVP基组能够较好地满足计算的精度要求。由于酯分子原子数较大，导致计算量指数级增长，为保证一定的精度和计算效率，本文的基组设置为6-31G*。

已有文献中尚无甘油三酸酯键解离能相关的试验数据，仅有甘油三油酸酯分子的光谱试验数据。为验证DFT的B3LYP/6-31G*基组的有效性，对甘油三油酸酯中C=O与C=C的振动频率进行了计算，结果如表1所示。

表1 甘油三油酸酯分子不同位置的酯基和碳碳双键振动频率Tab.1 Vibrational frequencies of ester groups and carbon carbon double bonds at different positions of glycerol trioleate

计算得到的振动频率通过美国国家标准与技术研究院（national institute of standards，NIST）提供的修正系数（0.96）对计算的振动频率进行修正。计算值和测量值[18-19]之间的相对误差约为2%，验证了B3LYP/6-31G*机组的可靠性。

2.3 分子动力学模拟过程

为了研究不同酯类油的热氧化过程差异，需要选取不同类型的分子分别进行建立晶胞模型。由于部分分子的密度未知，本研究在COMPASS力场下通过分子模拟获得密度并与试验值比较，以验证模拟结果的可靠性。基于COMPASS力场在1个大气压、293 K的等温等压NPT系综中预平衡10 ns。将预平衡系统在NVT正则系综中在450 K进行5 ns的模拟，然后在1个大气压和293 K下进行500 ps NPT模拟，统计密度稳定后的平均值。三羟甲基丙烷三酸酯在20℃的密度随碳链长度的变化如表2所示。由于目前商用三羟甲基丙烷酸酯种类较少，仅给出了三羟甲基丙烷辛葵酸酯与三羟甲基丙烷三油酸酯的实测值。结果表明计算误差在3%以内，模拟结果具有一定的可靠性。

表2 三羟甲基丙烷酸酯密度（20℃）Tab.2 Density of trimethylolpropane ester(20℃)

分子动力学模拟流程如下：

（1）建立基于计算获得密度的无定型晶胞，控制不同分子体系的体积大小接近。利用基于Becke Lee Yang Parr（BLYP）泛函的广义梯度近似模拟技术对氧分子进行优化，并将自旋多重性设为三重态，然后装入50个氧分子晶胞。无定型晶胞初始结构模型如图2所示，碳原子用灰色表示，氧原子用红色表示，氢原子用白色表示。

图2 无定型晶胞初始结构模型Fig.2 Initial structure model of amorphous cell

（2）基于COMPASS力场的NVT（恒定原子数、体积和温度）正则系综，采用退火动力学方法平衡结构获得在300～600 K温度内的全局最小能量。然后在353 K温度下以0.2 fs的时间步长动态运行1 ns，使退火后的体系进一步平衡。采用Berendsen恒温方法控制体系温度；范德华力的计算采用原子求和法（截止距离为15.5 Å），原子间的长程静电相互作用采用埃瓦尔德法计算。

（3）基于分子模拟的能量最小化后的结果，由于不同力场能量最小化结果不同，采用ReaxFF力场重新进行结构优化。然后进行动态模拟，设置平衡时间为10 ps，时间步长为0.25 fs，分别进行不同温度下的反应动力学模拟，获得热氧化反应动力学结果。

3 大豆油的仿真结果

对FR3和氧气进行分子动力学模拟。先控制温度为常温，在MS中设置环境温度为298 K，然后对优化好的FR3和氧气晶胞进行分子动力学模拟。使用Perl脚本进行产物分析时，发现氧气并没有减少，FR3的成分也没有变化。原因是仿真模拟的时间很短，自然环境条件下不足以让FR3发生自动氧化。因此所有的产物都没有改变。但是可以用该仿真结果计算这个情况下的氧气扩散系数。氧气扩散系数越大，晶胞中的FR3越容易与氧气接触，也就更容易发生氧化反应。计算结果如图3所示。从图3可以得到曲线拟合方程：y=2.587x+14.164，斜率k=2.587，氧气扩散系数为k的1/6，即0.431，远大于氧气在标况下的扩散系数（1.0×10-5）。因此植物绝缘油在常温下也容易自氧化。

图3 氧气扩散的均方根位移曲线Fig.3 Root mean square displacement curve of oxygen diffusion

将分子动力模拟中温度提高，以加速模拟反应速率，得到大豆油分子动力学模拟趋势，因此将分子动力模拟中的环境温度改为1 473 K，其余按照相同的步骤进行仿真模拟。不同的是，FR3晶胞中含有3种不同的长链脂肪酸，因此从结果中分析哪一种脂肪酸最先开始氧化分解。仿真到15 ps时，亚油酸（C18：2）最先开始氧化，然后是只含有一个双键的油酸（C18：1）分子，最后开始分解的是饱和脂肪酸棕榈酸（C16：0）。经过100 ps后，其组成的分子和分子量分布发生了明显的变化。热氧化过程中，甘油三酯中越不饱和的脂肪酸就越易氧化。根据最后仿真结果得出这3种物质的含量变化如图4所示。

图4 模拟过程中FR3成分含量变化Fig.4 The content of FR3 in simulation process

从图4可以看出，棕榈酸（C16：0）的含量表现稳定，只有小幅减少的迹象。油酸（C18：1）的含量发生了较为明显的变化，并没有因为氧化而出现显著的含量降低现象，所占的含量百分比出现了增长的现象。亚油酸（C18：2）含有较多的碳碳双键，随时间的变化含量变化更加强烈。当对油样进行加热升温后，随着加热时间的延长，油酸（C18:1）中包含的高分子量成分含量逐渐升高，虽然植物绝缘油的种类不同，但是经过高温加热后在植物绝缘油中发生的氧化聚合反应和氧化分解反应的类型较为相似，从而使得植物绝缘油中脂肪酸在经过相同时间的加热后其内部包含分子的分子量以及分子量的分布状况都极为相似，说明植物绝缘油中不同种类脂肪酸在高温条件下产生的聚合效应是类似的。

之后使用Perl脚本对FR3的仿真模拟结果进行量化分析，首先统计了FR3碳碳双键C=C以及酯基O=C-H的含量，结果如图5所示。从图5可以看出，脂肪酸中碳碳双键的含量大于酯基。由ReaxFF MD力场模拟结果可知，不饱和碳碳双键中α-C、β-C和γ-C较为活泼，会与氧气发生反应生成氢过氧化物。

图5 FR3仿真中C=C与O=C-H的含量Fig.5 Content of C=C and O=C-H in FR3 simulation

然后统计了引发阶段产物氢过氧化物的含量，如图6所示。从图6可以看出，热氧化初期，氢过氧化物的含量一直处于上升的趋势，但是在40 ps之前增长缓慢。原因为氧化初期，氢过氧化物由不饱和脂肪酸的碳碳双键氧化形成，但是双键在氧化初期活性低，氧化缓慢。模拟时间至70 ps时，氢过氧化物的含量大幅下降，这是由于热氧化至一定时间后，氢过氧化物开始分解生成羟基，导致含量减少。之后生成和分解过程同时进行导致含量波动变化。

图6 FR3仿真中氢过氧化物的含量Fig.6 The content of hydroperoxide in FR3 simulation

利用Origin软件对3种脂肪酸的残余量和加热时间进行多项式拟合，得到相关性高于0.99的回归方程，拟合结果如表3所示。

表3 仿真模拟中棕榈酸、油酸、亚油酸热氧化的动力学参数Tab.3 Kinetic parameters of thermal oxidation of palmitic acid,oleic acid,and linoleic acid in simulation

基于表3中3种脂肪酸在同一温度下的热氧化拟合方程，使用分子动力学分析技术确定高温下脂肪酸的反应级数。为了计算这一部分，将其简化为氧气与脂肪酸发生反应，生成一种氧化产物的过程，其宏观反应速率方程如式（3）所示。

式（3）中：rA为化学反应速率；k表示化学反应速率常数；cA为化学反应过程中脂肪酸的瞬时残余量数值；cB为化学反应过程中的氧气质量浓度；n为化学反应过程中作为反应底物的脂肪酸的反应级数；m为化学反应过程中作为反应底物的O2的反应级数[16]。通过拟合方程求导出的数据来进行线性拟合绘图，得到一条回归方程。通过使用线性拟合曲线的斜率n和截距lnk1得出3种脂肪酸的化学反应级数n以及反应速率常数k1，推导计算结果如表4所示。

表4 FR3中棕榈酸、油酸、亚油酸的热氧化拟合方程Tab.4 Fitting equation of thermal oxidation of palmitic acid,oleic acid,and linoleic acid in FR3

然后计算氧化活化能来比较这3种脂肪酸的抗氧化性，根据Arrhenius方程可以得到式（4）。

式（4）中：T为温度；Ea代表活化能；A为指前因子；R为理想气体常数。通过方程的斜率数值可以发现，呈负数的斜率数值表示lnk与1/T之间的关系是正线性相关的，然后根据计算所得的速率常数值，将lnk对1/T作线性回归分析，通过回归分析结果以及所构回归方程的直线斜率Ea/R以及截距常数数值lnA，联结起来对3种脂肪酸的热氧化降解反应的活化能以及指前因子进行准确计算，计算结果如表5所示。

表5 仿真中棕榈酸、油酸、亚油酸热氧化的活化能(Ea)和指前因子(A)Tab.5 Oxidation energy(Ea)and pre exponential factor(A)of thermal oxidation of palmitic acid,oleic acid,and linoleic acid in simulation

从表5可以看出，棕榈酸、油酸、亚油酸3种物质在氧化反应过程中所需的氧化活化能Ea出现依次降低的趋势，通过所需能量数值的不同可以推导出棕榈酸、油酸、亚油酸3种物质的热氧化安定性存在差异，其稳定性从高到低依次为棕榈酸、油酸、亚油酸。因此可以得出，当相对分子质量差异不大的情况下，不饱和度越高，碳碳双键对温度越敏感，发生反应所需的能量越低。