国产大型宽体客机发动机中的高温合金涡轮盘轴、压气机盘、挡板等关键热端承力构件,在服役过程中须长时间承受极端高温和较大载荷。在制造环节,受新型大尺寸锻件摩擦焊接、一体化机加工和表面强化等技术影响,容易引入较大梯度的残余应力分布,导致盘环件在自由状态下发生变形,显著影响构件的尺寸稳定性[1-4],并可能引发装配失效甚至干涉破坏。在服役阶段,内部高水平残余应力场与热、力载荷场相耦合,直接改变构件内部的受力状态,对其疲劳性能造成显著影响[5-7]。因此,及时、无损地测量大型高温合金锻件中的残余应力及其分布,有助于实现关键热端承力构件的寿命延长和轻量化设计,提高构件抗损伤容限能力,满足航空发动机日益苛刻的安全性和低成本要求[8-9]。
残余应力主要产生于热处理环节,尤其是在淬火与冷却阶段[10]。受温度场分布不均和相变引起的体积变化等因素影响,锻件内部极易形成梯度较大的三向残余应力场。这类应力场在随后的时效、机加工或表面处理过程中并不能完全消除,而是通过局部再分配持续在构件内部保留[11-12]。与塑性加工或表面强化工艺所引入的表层残余应力相比,热处理引起的深部残余应力更难预测和调控,且在盘锻件中普遍达到数百兆帕的水平,成为影响其长期稳定性与寿命的关键因素。因此,对发动机高温合金锻件深部残余应力水平的检测与表征具有重要意义。
目前,国内外针对残余应力的检测已经形成了多种手段,主要可分为有损和无损两类[13-15]。常规的有损检测如钻孔法、层剥离法、切割法等,虽然能够获得局部应力水平,但由于对零件造成不可逆的破坏,不适用于航空发动机关键构件的实际检测[16-19]。无损方法主要包括 X 射线衍射(X-ray diffraction,XRD)、中子衍射和同步辐射等方法[20-22]。X 射线衍射方法成熟可靠,适用于高精度测量表面及近表面残余应力,但受限于 X 射线穿透深度,该方法仅能检测盘锻件表面1~2 mm范围[23-24]。中子衍射法虽有厘米级穿透能力,但需依赖大型中子源装置,设备昂贵且测试周期长,难以在工程中推广[25-29]。同步辐射技术同样受限于设施条件和检测成本[30]。总体而言,目前衍射类无损检测方法在深部残余应力工程检测领域难以满足航空发动机关键构件的需求。
超声检测作为一种高效、低成本且适合复杂结构的无损检测方法,逐渐受到研究者的关注[31-32]。超声波在材料中传播时不仅具有强穿透性,其传播速度还对应力状态变化高度敏感,因此有望用于材料内部残余应力的无损表征。根据超弹性本构理论,材料的应力状态会使等效弹性常数发生变化,进而导致超声波速改变。因此,通过测量超声波不同传播模式和不同方向声束的波速变化,有望实现材料内部三向残余应力状态的反演和表征[33-34]。目前,已有文献报道利用超声临界折射纵波[35-36](longitudinal critically refracted wave, LCR)、瑞利波[37-38]、多模式超声导波[33]或超声体波[39-40]波速变化进行残余应力反演的方法。受制于波传播模式的限制,前3种方法主要适用于薄壁零件或厚度较小的结构件中,用于表征近表面残余应力梯度。而超声体波具备穿透大厚度结构的能力,适用于发动机金属盘锻件的三向残余应力场无损表征。
尽管如此,超声波与残余应力的耦合机制依然较为复杂。现有声弹性理论尚未建立完善的任意三向应力对任意传播方向超声体波波速变化的完整解析模型。因此,仍需进一步开展声弹性理论模型的深入研究,建立超声参量与残余应力变化的定量表征关系,为反演方法提供理论基础。
与此同时,如何结合工程实际需求,建立具有工程应用价值的高温合金深部残余应力场无损表征与成像方法,也是目前亟须研究的重点。声弹性理论模型中引入了非线性材料本构,基于该模型的反演方法存在未知变量维度高、超声透射路径设计困难等问题。需针对发动机关键构件,开展基于超声测量的残余应力场反演与成像方法研究,设计覆盖待测区域的声束路径及超声信号采集方法,实现三向残余应力场的无损成像。
本文基于超声测量方法,以航空发动机金属盘锻件所用典型高温合金淬火棒材为对象,提出了两种用于高温合金锻件残余应力场无损表征的方法,实现了任意位置三向应力的无损成像,并基于波传播有限元仿真,完成了两种方法的对比分析。
1 考虑任意三向应力作用的声弹性理论模型
弹性波在内部存在残余应力的金属零件中传播会受到声弹性效应影响,其波速变化与三向残余应力存在定量关系。目前的声弹性理论研究主要关注波传播方向与应力方向平行的单向应力情形,较少涉及三向应力共同作用下任意方向波速变化的声弹性理论模型。
1.1 超弹性本构下的应力-应变关系
为研究弹性波在内部存在残余应力物体中的传播规律,可采用自然状态、初始状态和最终状态描述该过程,如图1a所示。自然状态下物体处于无应力、无应变的未变形状态,其位置用矢量ξ或其分量ζα表示。初始状态下,物体受残余应力作用发生变形,位置矢量用X或分量XJ表示。最终状态是在已变形的物体上叠加超声波小扰动后进一步发生的动态微小变形,位置矢量用x或分量xj表示。三种状态下的所有物理量和材料属性分别用上标o(自然状态)、i(初始状态)、f(最终状态)表示,三者的相应分量分别以希腊字母、大写罗马字母和小写罗马字母为下标来表示。从自然状态到初始状态、自然状态到最终状态以及初始状态到最终状态的变形分别记为ui=X-ξ,uf=x-ξ和u=x-X。
采用自然状态坐标描述初始状态和最终状态的拉格朗日应变,可表示为
$\begin{array}{l}E_{\alpha \beta}^{\mathrm{i}}=\frac{1}{2}\left(\frac{\partial u_{\alpha}^{\mathrm{i}}}{\partial \xi_{\beta}}+\frac{\partial u_{\beta}^{\mathrm{i}}}{\partial \xi_{\alpha}}+\frac{\partial u_{\gamma}^{\mathrm{i}}}{\partial \xi_{\alpha}} \frac{\partial u_{\gamma}^{\mathrm{i}}}{\partial \xi_{\beta}}\right), \\E_{\alpha \beta}^{\mathrm{f}}=\frac{1}{2}\left(\frac{\partial u_{\alpha}^{\mathrm{f}}}{\partial \xi_{\beta}}+\frac{\partial u_{\beta}^{\mathrm{f}}}{\partial \xi_{\alpha}}+\frac{\partial u_{\gamma}^{\mathrm{f}}}{\partial \xi_{\alpha}} \frac{\partial u_{\gamma}^{\mathrm{f}}}{\partial \xi_{\beta}}\right).\end{array}$
式中:$E_{\alpha \beta}^{\mathrm{i}}$为初始状态下的拉格朗日应变;$E_{\alpha \beta}^{f}$为最终状态下的拉格朗日应变。由于超声扰动相较初始状态变形很小,即$u_{\alpha} \ll u_{\alpha}^{\mathrm{i}}, E_{\alpha \beta}^{\mathrm{f}}-E_{\alpha \beta}^{\mathrm{i}} \ll E_{\alpha \beta}^{\mathrm{i}}$,因此可以推导超声小扰动引起的应变近似为
$\begin{aligned}E_{\alpha \beta}=E_{\alpha \beta}^{\mathrm{f}}- & E_{\alpha \beta}^{\mathrm{i}}= \\& \frac{1}{2}\left(\frac{\partial u_{\alpha}}{\partial \xi_{\beta}}+\frac{\partial u_{\beta}}{\partial \xi_{\alpha}}+\frac{\partial u_{\gamma}^{\mathrm{i}}}{\partial \xi_{\alpha}} \frac{\partial u_{\gamma}}{\partial \xi_{\beta}}+\frac{\partial u_{\gamma}^{\mathrm{i}}}{\partial \xi_{\beta}} \frac{\partial u_{\gamma}}{\partial \xi_{\alpha}}\right) 。\end{aligned}$
某点的应力既可由柯西应力张量t表示,也可用第二比奥拉-克希霍夫应力张量T表示。根据物体的连续性和质量守恒定理,可以得到三种状态下密度的关系,从而推导初始状态下柯西应力张量和自然坐标下第二比奥拉-克希霍夫应力张量之间的关系
$\dot{t}_{J K}^{\mathrm{i}}=\frac{\rho^{\mathrm{i}}}{\rho^{\circ}} \frac{\partial X_{K}}{\partial \xi_{\beta}} \frac{\partial X_{J}}{\partial \xi_{\alpha}} T_{\alpha \beta}^{\mathrm{i}}$
式中,ρi、ρo分别为初始状态、自然状态下物体的密度。在初始状态物体上叠加超声小扰动后,会形成应力增量T=Tf-Ti,该应力增量在自然坐标系下可表示为$T_{\alpha \beta}=T_{\alpha \beta}^{f}-t_{\alpha \beta}^{i}$。物体从初始状态变换到最终状态为动态变形,其运动方程用自然坐标可表示为
$\frac{\partial}{\partial \xi_{\alpha}}\left(T_{\alpha \beta}^{f} \frac{\partial x_{j}}{\partial \xi_{\alpha}}\right)=\rho^{o} \frac{\partial^{2} x_{j}}{\partial t^{2}} 。$
超弹性指材料应力-应变存在非线性关系,则比奥拉-克希霍夫应力Tαβ与应变能函数W(Eαβ)满足Taβ=∂W/∂Eαβ。材料的应变能函数可表示为
$W(E)=\frac{1}{2} c_{\alpha \beta \gamma \delta} E_{\alpha \beta} E_{\gamma \delta}+\frac{1}{6} c_{\alpha \beta \gamma \delta \zeta \eta} E_{\alpha \beta} E_{\gamma \delta} E_{\zeta \eta}+\cdots 。$
式中,cαβγδ和cαβγδζη分别为材料的二阶和三阶弹性常数。
将预应变和超声小扰动引起的应力微小增量代入运动方程(4),且自然状态到初始状态为静态过程(即ti或Ti须平衡),可得声弹性波在自然坐标下的控制方程为
$\frac{\partial}{\partial \xi_{\beta}}\left[T_{\gamma \beta}^{\mathrm{i}} \frac{\partial u_{\alpha}}{\partial \xi_{\gamma}}+\Gamma_{\alpha \beta \gamma \delta} \frac{\partial u_{\gamma}}{\partial \xi_{\delta}}\right]=\rho^{\circ} \frac{\partial^{2} u_{\alpha}}{\partial t^{2}} 。$
式中,$\Gamma_{\alpha \beta \gamma \delta}=c_{\alpha \beta \gamma \delta}+c_{\alpha \beta \lambda \delta} \frac{\partial u_{\gamma}^{\mathrm{i}}}{\partial \zeta_{\lambda}}+c_{\lambda \beta \gamma \delta} \frac{\partial u_{\alpha}^{\mathrm{i}}}{\partial \zeta_{\lambda}}+c_{\alpha \beta \gamma \delta \zeta \eta} e_{\zeta \eta}^{\mathrm{i}}$,$e_{\zeta \eta}^{\mathrm{i}}$为初始柯西应力张量的分量。
同理可得到以初始坐标表示的声弹性控制方程为
$\frac{\partial}{\partial X_{J}}\left[\left(\delta_{I K} t_{J L}^{\mathrm{i}}+C_{I J K L}\right) \frac{\partial u_{K}}{\partial X_{L}}\right]=\rho^{\mathrm{i}} \frac{\partial^{2} u_{I}}{\partial t^{2}} 。$
式中,$C_{I J K L}=c_{I J K L}+c_{I J K L M N} e_{M N}^{\mathrm{i}}+c_{M J K L} \frac{\partial u_{I}^{\mathrm{i}}}{\partial X_{M}}+c_{I M K L} • \frac{\partial u_{J}^{\mathrm{i}}}{\partial X_{M}}+c_{I J M I} \frac{\partial u_{K}^{\mathrm{i}}}{\partial X_{M}}+c_{I J K M} \frac{\partial u_{L}^{\mathrm{i}}}{\partial X_{M}}$为初始坐标下的等效刚度。
上述以初始坐标表示的任意方向声弹性理论模型,适用于计算任意已知二阶及三阶弹性常数材料中超声传播的情形。本文研究的材料采用Murnaghan模型,其二阶弹性常数(Lamé常数)为λ和μ,三阶弹性常数为l、m和n。
1.2 任意单向应力下的声弹性理论模型
以单向应力σ为例,具体坐标系和符号如图1b。将σ所在轴设为x1轴,x1轴与超声传播方向构成x1y1平面,且x1轴与y1轴垂直。按右手法则选定z1轴方向,建立局部坐标系。由于超声在该平面内传播,假定超声传播方向为(cosθ,sinθ,0),其中θ为Ox1y1z1坐标系下纵波传播方向与y1轴的夹角,θ0为空间坐标系Oxyz下超声波与x轴的夹角。
为解析波速与应力的线性关系,将超声传播的解uI=UIeiK(NJXJ-Vt)代入式(7)中可得特征方程的展开式为
$\left(\begin{array}{ccc}D_{11}-\rho^{\circ} V^{2} & D_{12} & 0 \\D_{21} & D_{22}-\rho^{\circ} V^{2} & 0 \\0 & 0 & D_{33}-\rho^{\circ} V^{2}\end{array}\right)=0 。$
式中,D11、D22、D33为等效刚度系数,其为σ、θ的函数,并可由二阶、三阶弹性常数计算得到。
该方程有三组解,纵波波速为方程(D11-ρoV2)·(D22-ρoV2)-D12D21=0的最大根。假定纵波波速为VL,解得
ρoVL2= $\frac{({D}_{11}+{D}_{22})+\sqrt{({D}_{11}-{D}_{22}{)}^{2}+4{D}_{12}^{2}}}{2}$ 。
省略式(9)根号下内部关于σ的一次项,则可以解析地得到ρoVL2与应力之间的关系式,其形如一次函数。由于应力导致的声速变化相对于无应力状态时的声速较小,因此在之后的研究中主要关注应力变化量引起的声速变化量。对式(9)两边同时求导,则可得出在单向应力下纵波波速和应力之间的因变关系
$\mathrm{d} V_{\mathrm{L}}=\frac{\left[K_{1}+g\left(K_{1}\right) \sqrt{S_{1}}\right]}{4 \rho^{\circ} V_{L}} \mathrm{~d} \sigma 。$
式中,K1、S1为θ的函数,可由二阶、三阶弹性常数计算。当K1≥0时,g(K1)=1,当K1<0时,g(K1)=-1。由于应力导致的声速变化很小,故可近似地将VL≈VL0,VL0表示零应力下的纵波波速。
基于该单向应力声弹性理论模型,可通过拉伸实验标定各向同性材料的高阶弹性常数。在拉伸方向上施加不同幅值应力,测量超声以不同角度或模式穿过应力区域的波速变化,即可解得二阶、三阶弹性常数。与文献[33]理论模型对比,设置与波传播方向一致的应力σ=108、217、325 MPa,纵波波速变化量分别为ΔVL=8.4、16.1、24.5 m/s,两者之间误差分别为2.38%、2.48%、0.82%,验证了本文理论模型的准确性。
1.3 任意三向应力下的声弹性理论模型
进一步推导更一般情形,即三向应力作用下超声波速的理论模型,期望获得类似单向应力的线性关系。推导以空间坐标系为基础,应力用张量表达,由于应力张量包含较多未知量,空间坐标系下的控制方程在残余应力表征中较为复杂。为简化后续计算及残余应力表征,用主应力等效表达各点的应力张量,并在局部坐标系下进行表示。
超声在截面内传播,对应空间坐标系下的传播方向为 $\overrightarrow{{n}_{0}}$ =(cosθ0,sinθ0,0)。三向应力情形的坐标系设置及方向说明见图1c。三向应力矢量在局部坐标系x1y1z1下为 $\overrightarrow{\sigma }$ =(aσ,bσ,cσ),其中σ为三向应力大小,a为其在x1轴的投影系数,b为在y1轴的投影系数,c为在z1轴的投影系数。因截面在空间坐标系下的三向应力局部坐标系z1轴通常与z轴一致,将空间坐标系下的超声传播方向转换至应力主轴坐标系Ox1y1z1,对超声传播方向的坐标变换公式为
将超声传播解uI=UIeiK(NJXJ-Vt)代入式(8)中,可得到特征方程,同理可解得纵波波速变化量和应力变化量之间的线性解析表达式
$\mathrm{d} V_{\mathrm{L}}=\frac{[K+g(K) \sqrt{S}]}{4 \rho^{\circ} V_{\mathrm{L}}} \mathrm{~d} \sigma 。$
式中,K、S为θ的函数,可由二阶、三阶弹性常数计算。该式建立了任意三向应力和任意传播方向纵波波速变化量的显式定量关系,因此通过测量超声纵波波速变化量,即可反映棒材内部残余应力大小。通过测量超声沿不同角度与棒材交互的信号,即可反演三向残余应力场的分布情形。
进一步基于上述推导分析不同方向残余应力对超声波速变化的影响。棒材经淬火后可能产生300~600 MPa的轴向、周向、径向残余应力,最大残余应力通常分布在芯部,呈现“内拉外压”特点。基于理论模型计算,超声波以不同夹角穿透100 mm、残余应力为300 MPa区域后,纵波波速变化如图2所示。可见,纵波对与传播方向平行的残余应力最为敏感,波速变化达22.7 m/s;对与传播方向垂直的残余应力则不敏感,波速变化仅为11.2 m/s,在实际检测中难以准确测量。因此,需要采集多角度超声体纵波透射金属盘锻件的信号,以增强对各方向残余应力的敏感性。同时,三向残余应力场在棒材内部连续分布,物理模型反演面临未知量多、约束方程少等问题,难以直接实现精确反演,可考虑数据驱动的神经网络反演方法,通过多角度超声与残余应力交互数据,实现残余应力场的高精度表征。
2 含残余应力棒材中的超声波有限元仿真
镍基合金锻造棒材内部的残余应力主要来自淬火工艺的冷却阶段。采用有限元仿真方法,对淬火棒材中的残余应力产生过程及超声波传播过程进行模拟,提取超声参量与残余应力场的相关规律。
首先利用热力耦合有限元方法,数值模拟航空发动机典型高温合金棒材淬火过程中的温度场和应力场,获得淬火工艺后棒材内部的残余应力分布。进一步基于时域显式算法,对淬火后超声波在含残余应力金属棒材中的传播进行仿真,分析超声波在不同残余应力场下的传播速度变化,建立“残余应力场-超声测量信号”数据库。
2.1 基于热力耦合的棒材淬火过程有限元仿真方法
高温合金固溶后的淬火属于典型的瞬态热传导问题,其瞬态热传导控制方程为
κ $\left(\frac{{\partial }^{2}T}{\partial {x}^{2}}\right)$ +κ $\left(\frac{{\partial }^{2}T}{\partial {y}^{2}}\right)$ +κ $\left(\frac{{\partial }^{2}T}{\partial {z}^{2}}\right)$ +q=cρ $\frac{\partial T}{\partial t}$ 。
式中:κ为热导率;T为瞬态温度;q为内热源强度;t为时间;c为比热容;ρ为材料密度。
淬火过程中,零件由初始温度逐渐冷却至设定水温,其材料初始温度和终止温度均为已知,因此模型外表面采用第三类边界条件,数学表达式为
$-\left.\kappa \frac{\partial T}{\partial n}\right|_{\Gamma}=\left.h\left(T-T_{f}\right)\right|_{\Gamma}$。
式中:Tf为流体介质温度;h为对流换热系数;Γ为域的全部边界。
固体力学的基本控制方程为
$\rho \frac{\partial^{2} \boldsymbol{u}(x, t)}{\partial t^{2}}=\nabla \cdot \boldsymbol{\sigma}(x, t)+\boldsymbol{b}(x, t) 。$
式中:u(x,t)为位移;σ(x,t)为结构应力;b(x,t)为结构受到的体力。三者均随位置和时间变化。
有限元模拟过程中,材料采用热弹塑性本构模型进行计算。该模型中总应变分为弹性应变、塑性应变和热应变3部分,表达式为
ε=εe+εp+εTh。
式中:ε为总应变张量;εe为弹性应变张量;εp为塑性应变张量;εTh为热应变张量。
除上述控制方程和本构关系外,求解热力耦合问题还需给出应变与结构位移之间的几何方程:
$\boldsymbol{\varepsilon}(x, t)=\frac{1}{2}\left[\nabla \boldsymbol{u}(x, t)+(\nabla \boldsymbol{u}(x, t))^{\mathrm{T}}\right] 。$
热传导分析中所涉及的材料热力学参数包括热膨胀系数、热导率、比热容、密度以及杨氏模量等,具体见表1。
表1 GH4169高温合金热力学参数Tab.1 Thermodynamic parameters of GH4169 high-temperature alloy |
| 温度T/℃ | 热膨胀系数α/ (10-6·K-1) | 热传导系数κ/ (W·m-1·K-1) | 热容Cp/ (J·kg-1·K-1) | 密度ρ/ (kg·m-3) | 杨氏模量 E/GPa |
|---|---|---|---|---|---|
| 0 | 12.8 | 11.02 | 424 | 8 226 | 197 |
| 100 | 13.1 | 12.75 | 434 | 8 190 | 197 |
| 200 | 13.4 | 14.36 | 448 | 8 160 | 197 |
| 300 | 13.8 | 15.96 | 463 | 8 130 | 197 |
| 400 | 14.2 | 17.60 | 480 | 8 090 | 197 |
| 500 | 14.0 | 19.18 | 500 | 8 050 | 196 |
| 600 | 15.1 | 20.77 | 525 | 8 010 | 194 |
| 700 | 15.7 | 22.36 | 560 | 7 960 | 187 |
| 800 | 16.4 | 23.95 | 605 | 7 910 | 165 |
| 850 | 16.8 | 24.53 | 625 | 7 890 | 145 |
| 900 | 17.1 | 25.10 | 636 | 7 860 | 130 |
| 1 000 | 17.5 | 26.83 | 645 | 7 810 | 105 |
热-力耦合仿真模拟了980 ℃棒材从固溶处理炉中取出后在20 ℃水中进行淬火的过程。仿真获得的芯部温度随时间变化曲线如图4所示,棒材芯部温度在淬火过程中呈现出先快后慢的下降趋势。初始淬火时,材料表面与水快速接触,温度迅速下降,从而形成较大的表面与芯部温度梯度,导致热量由芯部向表面传导,芯部降温较快。随着冷却的持续,表面温度已大幅下降,温度梯度减小,热量传导的驱动力减弱,此时芯部降温速度下降。最终在1 500 s时水淬过程结束,棒材整体温度降至水温(20 ℃)。
棒材淬火后,中心截面的三向残余应力分布如图5所示。可见,周向与轴向方向在内部为拉应力、靠近表面为压应力,呈现“内拉外压”特征。自由表面处的径向应力为零,芯部的环向和径向应力相等,符合物理规律。仿真结果表明,在当前热处理工艺下,径向残余应力幅值范围为0~300 MPa,周向残余应力幅值范围为-350~250 MPa,而轴向残余应力幅值范围为-400~550 MPa。
2.2 基于固体力学的含残余应力棒材超声波传播仿真
进一步采用超声波时域显式仿真方法,对淬火后含有残余应力的棒材内超声波传播行为进行数值模拟,旨在建立“残余应力场-超声测量信号”数据库。
以1号监测点激励所获得的波场为例,分析超声波在含残余应力棒材中的传播情况,如图7所示。棒材内可同时激发出纵波、横波及Rayleigh波。由于样件边界约束,在边界处会出现与纵波波速相同、沿边界传播的头波(head wave)。其中,纵波传播速度最快,且对残余应力的变化最为敏感。因此,通过提取各监测点最早到达的纵波信号,分析纵波波速的变化即可反演棒材内部残余应力的分布。
提取与激励位置存在夹角的监测点所接收到的超声时域信号。4号监测点(间隔15°)与激励点距离较近,不同模式的波尚未完全分离,可能影响纵波波速的计算。而在8号监测点(间隔35°),纵波信号已与其他模式波有效分开;在21号监测点(间隔100°),纵波信号受到头波干扰,波速计算存在误差。因此,选取8~18、27~47、56~66共43个监测点,形成数据格式为72×43的超声时域信号,用于残余应力场反演。
3 基于超声测量的淬火棒材残余应力反演成像
3.1 基于迭代重建算法的淬火棒材三向残余应力层析成像
根据1.3节理论分析,三向残余应力会影响超声波在任意传播方向上的波速,且该影响难以解耦。考虑实际工程中的淬火棒材坯料生产过程,采用相同热处理工艺生产的同一批次棒材,其内部残余应力分布方向通常具有较高一致性。因此,可将热处理环节后棒材内部的三向残余应力方向分布规律作为先验信息,将三向应力未知变量简化为应力幅值问题,从而避免直接求解高维物理场带来的计算难题,提高方法的效率与工程适用性。
基于迭代重建算法的淬火棒材三向残余应力层析成像方法,如图8所示。将待重建区域离散为N=n×n个正方形网格,每个网格的应力幅值用二维函数f(x,y)表示。纵波在一段声程Lin内的传播时间为t。将dVL=(-Lin/t2)dt代入式(12),有
dσ=- $\frac{4{L}_{in}{V}_{L0}{\rho }^{o}}{(K+g(K\left)\sqrt{S}\right){t}^{2}}$ ·dt。
式中:K、S为θ的函数;VL0为零应力下的纵波波速;ρo为自然状态下物体的密度。由于应力引起的声时变化较小,可近似认为t≈t0。因此,对于一条完整的声束路径,残余应力导致的声时变化可表示为
Δt=∑dt=∑ $\left[-\frac{(K+g(K\left)\sqrt{S}\right){L}_{in}}{4{V}_{L0}^{3}{\rho }^{o}}\right]$ dσ。
当三向残余应力方向已知时,式(19)中括号内的内容即已完全确定。由此,原问题转化为线性方程组求解,其未知量为3N个,方程数为所采集全部超声声束的数量。采用Kaczmarz交替投影法求解此线性方程组。以2.1节棒材淬火过程热仿真结果中的残余应力方向及2.2节淬火棒材超声波传播仿真结果中波传播信号作为输入,反演获得三向残余应力场的成像结果(如图9所示)。
图8 基于迭代重建算法的淬火棒材三向残余应力层析成像方法流程图Fig.8 Flowchart of the triaxial residual stress tomographic imaging method for the quenched bar based on an iterative reconstruction algorithm |
三向残余应力层析成像结果呈现“内拉外压”特征,符合通常规律。其中,径向残余应力最大值约325 MPa,并在表面处为0,满足自由边界条件。周向残余应力最大值约320 MPa,轴向残余应力最大约700 MPa。反演结果中,高残余应力区域与热仿真结果相符,应力幅值偏大。
由于高温合金盘锻件内部残余应力多呈区域分布,工程关注点主要在高应力区域,以预防疲劳裂纹萌生和薄壁部位变形。为适应实际工程需求,以芯部残余应力高值区域的平均值作为评价指标,定量比较两种方法的反演精度。采用迭代重建方法的反演结果显示,径向残余应力幅值为0~325 MPa,周向为-400~320 MPa,轴向为-420~700 MPa。反演结果与设定的残余应力真值处于同一应力水平。反演结果中三向残余应力幅值和分布范围较好反映了实际锻件的应力状态。
3.2 基于神经网络方法的淬火棒材三向残余应力反演
在已有大量超声测量残余应力数据集的前提下,基于神经网络的数据驱动反演方法可以在无需先验信息的条件下,建立超声测量数据与残余应力场之间的映射关系,避免直接求解式(12)描述的复杂物理规律。而基于2.1与2.2节建立的有限元仿真方法,可快速获得棒材在常用热处理工艺后的内部残余应力分布及相应超声纵波结果。通过设置不同热处理工艺参数,可批量生成“残余应力场-超声测量信号”数据,为基于神经网络的残余应力场反演方法提供数据集,数据集来源与分配数量见表3。
表3 数据集来源与分配数量Tab.3 Sources of the dataset and allocation quantities |
| 数据来源 | 训练集 | 验证集 | 测试集 | 数据总量 |
|---|---|---|---|---|
| 有限元仿真 | 4 000 | 1 000 | 100 | 5 100 |
| 理论计算 | 4 000 | 1 000 | 100 | 5 100 |
不同数据集间无重复样本。其中,2.2节的仿真结果被归入测试集。
为降低神经网络接受冗余信息的负担,输入数据需预处理。为排除棒材几何尺寸变化对神经网络精度的干扰,选取不同声束路径对应的声时差作为神经网络输入,数据形式为二维矩阵[Δtij]i×j,其中,i、j分别表示激励与接收超声压电单元的编号,Δt为当前路径纵波到达时间减去同路径零应力时纵波到达时间。
搭建的BP神经网络架构如图10所示。经前处理后的声时差数据作为输入,经过两层全连接后映射为三向残余应力大小。网络中激活函数选用ReLU,反向传播采用随机梯度下降法,最大训练轮数设为50轮,批量大小为20样本,学习率为0.000 1,损失函数为均方误差。输出数据为目标反演区域内每一像素点的三向残余应力,形式为σn×3,其中n为目标反演区域内的像素点数,此处n取10 201。
以第2.2节所述样本为例,分析训练后的神经网络性能。基于神经网络法的残余应力反演成像结果如图11所示。
反演结果表现出与图9所示相同的规律。其中,径向残余应力最大值约315 MPa,且表面处为0;周向残余应力最大值约350 MPa,轴向残余应力最大值约600 MPa。神经网络反演结果表现出关于棒材轴线的对称性,且精度高于基于迭代重建算法的反演结果。
神经网络反演结果显示,径向残余应力幅值为0~315 MPa,周向为-400~350 MPa,轴向为-400~600 MPa。两种方法的反演结果均与设定的残余应力真值处于同一应力水平。其中,迭代重建法对锻件内部残余应力特殊分布具有更高敏感性。例如,若棒材淬火受夹具影响导致传热不均,迭代重建法能有效反映内部残余应力的非对称性。而神经网络方法则充分学习了典型热处理工艺下的残余应力分布规律,具备对复杂模式的非线性拟合能力,因此其芯部高应力区平均误差最小,反映了在典型工艺条件下的高效反演能力。
4 结论
本文针对航空发动机高温合金锻件中残余应力场的无损表征需求,提出并建立了两类基于超声测量的三向残余应力成像方法。第一种为基于声弹性理论和迭代重建算法的三向残余应力层析成像方法,在应力主方向已知且工艺条件稳定时,可实现三向残余应力分布的有效重建,适用于热处理工艺较固定的棒材或锻件。第二种为基于神经网络的三向残余应力反演方法,依托高质量残余应力数据库,可适应工艺复杂、结构形貌多样的锻件与零部件,具备较强泛化能力与快速反演特性。两种方法互为补充,构成了适用于高温合金构件全生命周期的分阶段、可追踪残余应力无损表征技术体系。
本文研究表明,超声测量法在发动机关键构件残余应力检测中具有良好的应用前景,为传统检测方法难以实现的深层应力无损获取提供了新思路。未来工作需更充分考虑锻件表面加工硬化、表面强化处理及组织不均对超声响应的影响,建立统一的误差标定和补偿体系。在检测环节上,对于几何外形复杂的坯料或盘锻件,可开展超声体波透射路径设计,搭建快速信号采集系统,实现高信噪比和高覆盖率的数据采集。同时,结合多模态无损检测方法和更高精度的数值模拟,有望进一步提升残余应力反演的准确性与鲁棒性。总体来看,本文提出的方法为航空发动机关键热端构件全生命周期残余应力监测提供了重要技术支撑,对提升构件服役安全性和可靠性具有重要工程意义。