陶瓷基复合材料氧化行为与剩余强度数值分析

刘宝瑞 李尧 侯传涛 王建民 张东

陶瓷基复合材料氧化行为与剩余强度数值分析

刘宝瑞1李尧1侯传涛1王建民1张东2

（1 北京强度环境研究所可靠性与环境工程技术重点实验室，北京 100076；2 北京宇航系统工程研究所，北京 100076）

陶瓷基复合材料是当前高速飞行器热结构最广泛使用的材料之一，但在高温环境下长时间使用中面临着不可避免的氧化损伤问题。本文以C/SiC复合材料为对象，基于失重率-时间关系曲线，建立了包含温度和应力影响的氧化扩展速率模型，并提出了一种氧化损伤的刚度退化准则和剩余强度模型。通过编写UMAT子程序，在Abaqus软件框架下实现了陶瓷基复合材料氧化行为和剩余强度的数值分析。通过5个典型算例的计算分析以及与试验数据对比，验证了本文提出的模型和方法的有效性。

陶瓷基复合材料；氧化行为；剩余强度；数值分析

0 引言

以C/SiC为代表的陶瓷基复合材料是一类新型的热结构材料，具有耐高温、高比强、高比模等特点，其长时间使用温度高达1600℃，在先进飞行器等热结构设计中得到大量应用。陶瓷基复合材料具有优良的高温力学性能和化学稳定性，但是在400℃以上，随着O2、H2O等氧化性气体通过微裂纹和孔洞侵入结构内部，将导致C纤维、PyC界面层等发生氧化损伤，使得材料力学性能下降。外部严酷复杂的力热载荷环境会引起表面抗氧化涂层和SiC基体内部的微裂纹、孔洞等缺陷发生扩展演化，导致局部氧化性气体的进入量增加、加快氧化反应过程，严重时可能导致复合材料结构的整体灾难性破坏，已成为陶瓷基复合材料结构长时间使用面临的关键问题。

国内外在相关飞行器研制的牵引下，针对复合材料结构的氧化问题开展了大量研究工作。Lamouroux[1-3]和Naslain[4]等研究发现，复合材料内部基体裂纹、纤维与热解碳之间的分离间隙、层间间隙等会形成相互连通的网状孔隙，使得氧气、水蒸气等氧化性气体能够进入到复合材料内部，与C/SiC复合材料表面和内部的C元素发生氧化反应，导致材料性能下降。Medford[5]发展了最早的氧扩散分析模型之一，用于预测航天飞机RCC翼前缘面板的氧化行为，引入涂层厚度方向氧气浓度梯度驱动的稳态扩散过程假设，模拟氧气通过SiC涂层裂纹进入到RCC内部的扩散过程。Eckel[6]发展了一个类似的模型，采用气体稳态扩散方程模拟氧气在管状通道中的扩散过程，随后Halbig[7]采用了Eckel的模型对C/SiC复合材料碳纤维的氧化退化过程进行了模拟。Bacos等[8]在引入动量守恒，考虑材料内部氧化性气体扩散过程和化学反应过程，发展了C/C复合材料氧化行为分析模型[9]，应用于飞行器再入过程中C/C复合材料结构氧化行为计算与分析[10]。Sullivan[11-12]基于理想气体在多孔固体中的流动机理，建立了C/SiC复合材料结构氧化行为预测的解析理论和数值分析方法。国内成来飞、殷小玮、卢国锋等人[13-18]对C/SiC复合材料的氧化问题进行了长期的试验与分析研究，积累了大量有价值的试验数据，获得了各种环境下C/SiC的氧化机理；高魁垠等人[19]对C/SiC复合材料的力学性能分散性、非线性等问题进行了研究。以上这些学者在复合材料氧化失效机理、细观尺度氧化行为建模分析等方面进行了研究，但是普遍没有和复合材料宏观的力学性能建立关联，不适用于宏观C/SiC复合材料结构氧化损伤和剩余强度分析。本文基于试验获得不同温度下的失重率随时间变化规律曲线，建立了一种便捷的、适用于试验数据拟合的氧化扩展速率数学模型并通过引入应力加速因子模拟应力的影响，提出了一种氧化损伤的刚度退化准则和剩余强度模型，最后通过编写UMAT子程序，在ABAQUS®框架下实现了氧化扩展和氧化损伤的数值分析，并对典型算例进行计算与分析。

1 氧化损伤分析数学模型

复合材料氧化试验中，失重率是少数可以获取的、与氧化程度相关的物理量，因此氧化扩展分析大多基于失重率的变化。卢国锋等[16]和曹素等[17]针对C/SiC复合材料通过开展静态空气中的氧化试验研究获得材料氧化失重率变化规律，对于不含涂层的C/SiC复合材料，随着温度的升高，氧化速率显著升高，氧化失重率随时间的变化关系总体呈线性规律，且完全氧化后的失重量约为37%；对于含SiC涂层的C/SiC复合材料，氧化失重率存在较长的线性段，且在700℃左右氧化速率最大。

在上述学者研究的基础上，以CVI工艺制备的C/SiC复合材料薄壁结构为研究对象，引入如下假设：1）氧化扩展速率正比于失重率；2）结构表面的涂层细观裂纹分布均匀，结构氧化扩展和氧化反应从表面向中心均匀进行；3）C/SiC复合材料薄壁结构的材料力学性能正交各向异性。

对于C/SiC复合材料薄壁结构，根据氧化扩展速率和氧化失重率成正比的假设，某特定温度的空气环境下C/SiC复合材料氧化扩展深度为

式中，为氧化失重率，为比例系数，为氧化时间。根据氧化反应从表面向中心均匀进行的假设，氧化扩展的最短路径垂直于厚度方向，通过测量某时刻的氧化深度和失重率，即可获得比例系数，根据试验数据，对于3mm厚的复合材料平板，取值为4.05mm-1，由此可以获得氧化深度随时间的变化关系，进一步可以获得厚度方向氧化扩展速率

采用上述方法，获得各个试验温度状态下的氧化扩展速率，其他温度状态采用线性插值的方法计算相应温度点的氧化速率。

殷小玮[15]、卢国锋[16]、曹素[17]等人的研究表明，复合材料在不同温度下的氧化失重率随时间的变化规律存在显著的线性段。基于文献[15]和[17]的试验数据，通过公式(1)、(2)计算得到含SiC涂层C/SiC复合材料氧化扩展速率的结果如图1所示。从图1中可以看出，复合材料在400℃以下氧化速率为0，在700℃时氧化速率达到峰值，随着温度的进一步升高，表面抗氧化涂层以及基体内部的裂纹发生闭合，在1100℃时氧化速率相对降到较低的水平。随着温度进一步升高，复合材料结构内部缺陷逐渐增多，氧化速率缓慢上升，超过1400℃以后，由于表面抗氧化涂层的逐渐趋于熔融状态，出现较多大尺寸的孔洞等缺陷，导致复合材料本体裸露，使得氧化速率急剧上升[20]。

图1 不同温度下C/SiC的氧化扩展速率

针对应力的影响，引入应力氧化加速因子s，无应力状态下应力加速因子s=1，拉伸应力下s>1，压缩应力下s<1。由于变形和应力水平呈正比，因此假设应力氧化加速因子与应力水平呈正比，对于拉伸应力，应力氧化加速因子为

式中，st表示拉伸氧化加速因子，t表示比例系数，pmax表示面内最大主应力，t表示材料拉伸强度。对于压缩应力，应力氧化加速因子为

式中，sc表示压缩氧化加速因子，c为比例系数，pmin表示面内最小主应力，c表示材料压缩强度，0表示裂纹闭合应力，根据文献[21]的研究，对于C/SiC复合材料而言，0=128.5MPa。如果面内最小主应力超过裂纹闭合应力，那么如果应力继续增大，由压缩应力引起的裂纹闭合导致的复合材料氧化速率减缓效应将不再增加。

对于每个积分点，应力氧化加速因子按照如下准则进行选取

由此便可以根据结构表面应力场获得结构表面应力氧化加速因子分布规律，从而获得结构有限元模型每个积分点的氧化扩展速率模型

式中，0表示初始氧化速率，n表示考虑局部应力分布影响的氧化速率。

建立同时考虑力学性能非线性和氧化损伤的C/SiC复合材料正交各向异性本构模型，如式（7）所示

通过改变式(7)中复合材料的模量i和ij来表征材料力学性能非线性和氧化后力学性能的衰减。在面内方向上，采用多项式函数建立无氧化条件的本构关系，氧化影响通过氧化模量损伤因子e进行考虑，如式(8)和(9)所示

式中，i和i（=1,2,…）为系数，由无氧化条件下试验应力-应变曲线拟合得到，具体取值参考文献[20]，i为面内主方向的应变，12为面内剪切应变，e为氧化模量损伤因子，取值为e=10-4。

在面外方向上，近似采用线性本构，考虑氧化损伤的模量采用式(10)和(11)表征

式中，i、ij为氧化后材料模量，i,0、ij,0为氧化前材料初始模量。

氧化后材料的剩余强度如式（12）所示

式中，b为氧化后材料的剩余强度，b,0为氧化前材料初始强度，s为氧化强度损伤因子，取值为：s=0.02。

2 算例分析

2.1 平板结构氧化损伤行为与剩余强度分析

针对典型C/SiC平板进行力热氧耦合环境下力学行为分析。平板一端固支，另一端施加拉伸载荷。为了降低计算规模，对模型进行简化，从C/SiC平板中切取一小片试样，试样上下两个表面为氧化起始界面，侧面的一端固支，另一端施加位移载荷，如图2所示。分两个载荷步进行计算：1)输入氧化时间、氧化温度，计算氧化过程；2)施加载荷，计算氧化后结构的剩余强度。

图2 典型C/SiC试样有限元模型

针对1300℃空气环境下的C/SiC复合材料的氧化过程和氧化后力学行为进行计算分析。计算获得C/SiC平板沿厚度方向的氧化扩展深度随时间变化如图3所示，曲线图如图4所示。

图3 1300℃空气环境下C/SiC氧化扩展深度

在1300℃空气环境下不同氧化时间作用后对平板施加拉伸载荷，获得载荷位移曲线如图5所示，由于平板的破坏位移主要取决于中间未氧化层的断裂应变，因此不同氧化时间后的破坏位移大小基本一致。剩余拉伸强度随氧化时间的变化如图6所示，近似呈现出线性递减的趋势。

图4 1300℃空气环境下C/SiC氧化扩展深度随时间变化

图5 1300℃空气环境下不同氧化时间后的载荷位移曲线

图6 1300℃空气环境下不同氧化时间后剩余拉伸强度

针对典型C/SiC复合材料试样，计算不同温度下氧化10小时后的氧化扩展情况，计算结果如图7所示。可以看出，500℃氧化10小时后材料近似无氧化，700℃时氧化最为严重。

进一步计算不同温度氧化10小时后材料的剩余力学性能和强度，获得氧化后材料的拉伸载荷与位移关系曲线，如图8所示。采用式（13）计算材料试样的等效拉伸强度，获得不同温度氧化10小时后材料的剩余拉伸强度如图9所示，计算与试验结果[18]的规律基本一致，在700℃时强度最低，500℃以下及1100℃的强度相对较高

式中，σT为拉伸破坏强度，FT为拉伸破坏载荷，A为试验件横截面积。

图8 空气环境不同温度氧化10h后拉伸载荷与位移曲线

图9 空气环境不同温度氧化10小时后C/SiC剩余强度与温度的关系

之所以计算剩余拉伸强度，是因为目前的失效判据在计算拉伸强度时精度最高，弯曲载荷材料试样的应力分布相对比较复杂，失效判据的预测精度较低。而试验中考虑到试验实施的可行性和数据获取的有效性，氧化后剩余强度通常采用弯曲强度来表征。虽然计算和试验采用的强度表征量不同，但获得的氧化后剩余强度变化规律在理论上应是大体一致的。

2.2 温度场分布不均匀的复合材料平板氧化扩展行为分析

对于温度场分布不均匀的复合材料结构，由于不同位置温度不同导致不同点氧化扩展速率不一样，一定时间后的氧化损伤也呈现出不同的分布规律。针对某C/SiC复合材料平板，假设其受到单面加热，温度场分布呈中间高、四周低的规律，加热面温度高、背热面温度低，如图10所示，中心峰值温度1000℃，四周最低温度382℃。对该温度场分布状态下复合材料平板在经过10小时后的氧化扩展情况进行计算分析，计算获得平板的氧化损伤扩散如图11所示。从图中可以看出，在10小时以内，低于500℃的区域氧化损伤基本没有扩展；中心1000℃的区域，由于该温度下复合材料裂纹闭合，氧化损伤扩展相对较小；在700℃左右的区域，氧化损伤最严重。模型的计算结果较好地模拟了由于温度场不均匀导致的氧化损伤不均匀现象，氧化后的平板的不同区域呈现出不同程度的氧化损伤情况。

图10 复合材料平板温度场分布（1/4模型）

图11 平板氧化扩展计算结果（1/4模型）

2.3 含表面涂层裂纹的平板氧化扩展行为分析

在实际工程应用过程中，经常出现复合材料结构表面涂层发生损伤的情况，为了评估含涂层损伤的复合材料结构在热力氧耦合环境下的力学行为，建立了考虑涂层损伤的复合材料结构氧化扩展模拟方法。

对于涂层损伤部位，如划伤、片状剥落等，可以近似简化为不含涂层的复合材料结构表面进行处理。根据文献[16]中试验结果，可以获得不含涂层的C/SiC复合材料氧化速率，如图12所示，从曲线中可以看出，不含涂层的C/SiC复合材料氧化速率随着温度的增加呈现出指数级增大的规律。

图12 不含抗氧化涂层的C/SiC复合材料氧化速率

为了计算含抗氧化涂层裂纹损伤的复合材料结构氧化问题，需要判断结构中每一个积分点在任意时刻的损伤状态，其中需要计算结构中任意点P距离涂层损伤位置的距离，假设结构中任意点距离线状损伤两端点的距离分别为a和b，线段长度为c，则采用海伦公式计算三角形面积S

由此获得结构中任意点P距离线状损伤的最短距离

对于线段两端点为锐角的情况：即为氧化扩展最短路径（|2-2|<2）；对于线段两端点存在钝角的情况(|2-2|>2)，扩展最短路径取Min(,)。从表面扩散到该点的时间为0，从裂纹处扩散到该点的时间为1，比较0和1，选用时间较小的值。

假设某3mm厚的C/SiC平板表面存在约5mm长的抗氧化涂层划伤，通过计算研究该涂层损伤导致的局部氧化扩展行为。截取涂层周边一小块样品为分析域建立有限元分析模型，如图13所示。假设试样暴露于1300℃的空气环境中10小时，计算试样的氧化扩展过程。计算获得试样不同时刻的氧化扩展情况如图14所示。浅色为未氧化部分，深色为发生氧化损伤的部分。从图中可以看出，通过涂层划伤部分的氧化扩展速率明显比未损伤部分要快，并且形成了明显凹坑。

图13 含表面涂层裂纹损伤的C/SiC平板有限元模型

Fig.13 Finite element model of C/SiC plate with surface coating crack damage

图14 氧化扩展过程计算结果(1/4模型)

2.4 考虑结构应力的复合材料结构氧化行为分析

针对某C/SiC复合材料开孔板，计算结构在承受轴向拉伸载荷产生的应力场分布情况下的氧化行为。建立其1/4有限元分析模型，在对称面上施加对称边界条件，在结构的一端施加0.05mm位移模拟轴向拉伸载荷，如图15所示。

图15 C/SiC开孔板有限元模型（1/4模型）

计算获得结构的主应力场分布，如图16所示，可看出开孔的边缘在垂直轴向的方向承受拉应力，在沿轴向的方向承受压应力；计算获得结构在800℃环境下氧化10小时后应力氧化加速因子的分布如图17所示，从图17中可以看出，在面内拉应力较大的区域应力氧化加速因子相对较大，在面内压应力较大的区域应力氧化加速因子相对较小。图18所示为开孔板结构氧化扩展情况分布云图，从图中可以看出，在承受拉伸应力较大的区域，氧化扩展程度相比无应力状态明显更深，而承受压缩应力较大的区域，氧化扩展程度较浅。

图16 开孔板结构主应力场分布

图17 开孔板结构应力氧化加速因子分布

图18 开孔板结构氧化扩展分布

从该算例可以看出，本文发展的模型很好地模拟了应力场分布对氧化扩展过程的影响，通过修正应力氧化加速因子的方程，能够与试验结果取得更好的一致性。

2.5 考虑结构应力、温度分布的复合材料结构氧化行为分析

针对2.4节所示的C/SiC复合材料开孔板，在结构轴向拉伸载荷下产生的应力场分布的基础上，施加不均匀温度场分布（如图19所示），开展力热氧耦合计算，获得结构氧化10小时后的损伤扩展情况如图20所示，从图中可以看出，模型由于应力场、温度场分布的不均匀，氧化损伤扩散的程度呈现出特殊的分布规律，在拉应力较大、温度接近700℃的区域，氧化扩展程度较深，压应力较大、温度较低或接近温度1000℃的区域，氧化扩展程度较浅。该算例验证了模型具备力热氧耦合计算分析的能力。

图19 开孔板结构温度场分布

图20 开孔板结构氧化扩展分布

3 总结

建立了同时考虑力学性能非线性和氧化损伤的C/SiC复合材料正交各向异性本构模型，预测了典型复合材料平板的氧化扩展过程与剩余强度，通过试验数据对比分析验证了方法的有效性；针对温度场分布不均匀平板、表面涂层裂纹损伤平板、同时考虑结构应力和温度分布的复合材料结构等进行计算分析，预测了结构的氧化损伤情况；研究方法可为复杂环境下飞行器热结构强度评估和优化设计提供技术支撑。

[1] Lamouroux F, Bourrat X, Nasalain R, et al.Structure/ oxidation behavior relationship in the carbonaceous constituents of 2D–C/PyC/SiC composites[J].Carbon, 1993, 31(8):1273-1288.

[2] Lamouroux F , Camus G , J Thébault.Kinetics and mechanisms of oxidation of 2D Woven C/SiC composites: I, experimental approach[J].Journal of the American Ceramic Society, 1994, 77(8):2049-2057.

[3] Lamouroux F, Naslain R, Jouin J.Kinetics and mechanisms of oxidation of 2D woven C/SiC composites: II theoretical approach[J].Journal of the American Ceramic Society,1994,77(8): 2058-2068．

[4] Naslain R , Guette A , F Rebillat, et al.Oxidation mechanisms and kinetics of sic-matrix composites and their constituents[J].Journal of Materials Science, 2004, 39(24):7303-7316.

[5] Medford J E.Prediction of oxidation performance of reinforced Carbon-Carbon material for space shuttle leading edges[C].AIAA 10th Thermophysics Conference, Denver, CO, 1975.

[6] Eckel A J , Cawley J D , Parthasarathy T A.Oxidation kinetics of a continuous Carbon phase in a nonreactive matrix[J].J.Am.Ceram.Soc., 1995,78(4): 972-980.

[7] Halbig M C.The oxidation kinetics of continuous Carbon fibers in a cracked ceramic matrix composite[R].NASA/TM-2001-210520, 2001.

[8] Bacos M P ,Dorvaux J M ,Lavigne O, et al.C/C composite oxidation model: I.morphological experimental investigations[J].Carbon, 2000, 38(1):77-92.

[9] Bacos M P, Cochon J L, Dorvaux J M, et al.C/C composite oxidation model: II.Oxidation experimental investigations[J].Carbon, 2013, 38(1):93-103.

[10] Bacos M P , Dorvaux J M , Lavigne O , et al.C/C composite oxidation model: III.Physical basis, limitations and applications[J].Carbon, 2000, 38(1):105-117.

[11] Sullivan R M .A model for the oxidation of carbon silicon carbide composite structures[J].Carbon, 2005, 43(2):275-285.

[12] Sullivan R M.A model for the oxidation of C/SiC composite structures[R].NASA/TM-2003-212720, 2003.

[13] Cheng L, Xu Y , Zhang L, et al.Oxidation behavior of three dimensional C/SiC composites in air and combustion gas environments[J].Carbon, 2000, 38(15):2103-2108.

[14] Cheng L, Xu Y, Zhang L, et al.Effect of carbon interlayer on oxidation behavior of C/SiC composites with a coating from room temperature to 1500°C- ScienceDirect[J].Materials Science and Engineering: A, 2001 (1-2):219-225.

[15] 殷小玮.3D C/SiC复合材料的环境氧化行为[D].西安:西北工业大学, 2001:86-87.

[16] 卢国锋, 乔生儒, 弓满锋, 等.C/Si-C-N复合材料的制备及其氧化行为研究[J].材料工程, 2010(3):13-17.[Lu Guofeng, Qiao Shengru, Gong Manfeng, et al.Fabrication and oxidation behavior of C/Si-C-N composite[J].Journal of Materials Engineering, 2010 (3):13-17.]

[17] 曹素, 刘永胜, 左新章, 等.涂敷含硼硅玻璃SiC涂层的C/SiC复合材料空气氧化行为[J].复合材料学报, 2011, 28(2):142-148.[Cao Su, Liu Yongsheng, Zuo Xinzhang, et al.Oxidation behaviors of C/SiC composites coated with SiC coatings containing borosilicate glass[J].Acta Materiae Compositae Sinica, 2011, 28(2):142-148.]

[18] 李建章, 张立同, 成来飞, 等.高温氧化气氛下3D C/SiC质量变化率与剩余强度的相关性[J].复合材料学报, 2007,24(4):101-105.[Li Jianzhang, Zhang Litong, Cheng Laifei, Relationship between mass variation and residual strength change with temperature on 3D C/SiC composites in oxidizing atmosphere [J].Acta Materiae Compositae Sinica, 2007,24(4):101-105.]

[19] 高魁垠, 李海波, 吴建国, 等.2D-C/SiC 复合材料螺栓连接结构可靠性分析[J].强度与环境, 2020, 47(1): 33-40.[Gao Kuiyin, Li Haibo, Wu Jianguo, et al Reliability analysis of 2D-C/SiC composite bolted joints[J].Structure & Environment Engineering, 2020, 47(1): 33-40.]

[20] 张立同.纤维增韧碳化硅陶瓷复合材料—模拟、表征与设计[M].北京,化学工业出版社, 2009:152-154.

[21] 李俊, 矫桂琼, 王波, 等.二维编织C/SiC复合材料非线性损伤本构模型与应用[J].复合材料学报, 2013, 30(1):165-171.[Li Jun, Jiao Guiqiong, Wang Bo, A non-linear damage constitutive model for 2D woven C/SiC composite material and its application[J].Acta Materiae Compositae Sinica, 2013, 30(1):165-171

Numerical Analysis Method of Oxidation Behavior and Residual Strength for Ceramic Matrix Composite

LIU Bao-rui1LI Yao1HOU Chuan-tao1WANG Jian-min1ZHANG Dong2

（1 Science and Technology on Reliability and Environmental Engineering Laboratory, Beijing Institute of Structure and Environment Engineering, Beijing 100076, China；2 Beijing Institute of Aerospace System Engineering, Beijing 100076, China）

Ceramic matrix composites (CMCs), one of the most widely used materials in the hot structures of high-speed aerial vehicles; suffer from inevitable oxidation problem in the long-time service under high temperature.Taking the C/SiC composites as a prototype, this work build an oxidation diffusion speed model including the effects of temperatures and stresses based on the curves of the weight loss versus service time, and proposed a rule of stiffness degradation and a residual strength model for oxidation damage.Through UMAT subroutine, this work has realized the numerical analysis of oxidation behaviors and residual strength for CMCs in the framework of ABAQUS software.The proposed models and method are validated through analysis of the calculation results of five typical examples and the comparisons with experimental data.

Ceramic matrix composite; Oxidation behavior; Residual strength; Numerical analysis

TB332, V25, V23

A

1006-3919(2021)06-0001-08

10.19447/j.cnki.11-1773/v.2021.06.001

2021-07-15；

2021-10-22

科技部国家重点研发计划（2021YFB3801700）;重点实验室基金（6142911180512）

刘宝瑞（1985—），男，高工，研究方向：飞行器结构强度评估；（100076）北京市9200信箱72分箱.