首页> 中国专利> 一种基于单纯形法的半无限域基体涂层结构材料弹性性质获取方法

一种基于单纯形法的半无限域基体涂层结构材料弹性性质获取方法

摘要

一种基于单纯形法的半无限域基体涂层结构材料弹性性质获取方法,材料的弹性系数与材料制备工艺密切相关;利用超声波对材料弹性性质进行测量是无损检测领域很有前景的测量方法之一。基于声学显微镜技术,自行开发的涂层材料弹性系数的超声测量系统,采用线聚焦PVDF探头,通过纵波和表面波波速的同时测量,可实现材料的弹性系数无损检测。本方法具有如下有益效果,不必再将整幅理论与实际曲线的波速平方的残差总和最小值都计算出来;不必把所有的理论频散曲线都确定,这种方式耗时且数据量冗长,单纯形的方法可以加速并化简程序。

著录项

  • 公开/公告号CN103926329A

    专利类型发明专利

  • 公开/公告日2014-07-16

    原文格式PDF

  • 申请/专利权人 北京工业大学;

    申请/专利号CN201410135724.7

  • 申请日2014-04-04

  • 分类号G01N29/44(20060101);G01N29/04(20060101);G06F19/00(20110101);

  • 代理机构11203 北京思海天达知识产权代理有限公司;

  • 代理人沈波

  • 地址 100124 北京市朝阳区平乐园100号

  • 入库时间 2023-12-17 00:20:51

法律信息

  • 法律状态公告日

    法律状态信息

    法律状态

  • 2017-08-25

    授权

    授权

  • 2014-08-13

    实质审查的生效 IPC(主分类):G01N29/44 申请日:20140404

    实质审查的生效

  • 2014-07-16

    公开

    公开

说明书

技术领域

本发明涉及一种基于单纯形法的半无限域基体涂层结构材料弹 性性质获取方法,属于超声导波无损检测与评估领域。

背景技术

表面工程是再造技术的重要组成部分,同时又为先进制造技术的 发展提供了很好的技术支撑,表面工程是指利用物理、化学或其他方 法在材料表面制备一层不同于基体材料且有一定厚度的强化、防护覆 盖层。涂层材料广泛应用于航空航天、机械、石油、化工、核电及微 器件与微制造等领域,其主要功能是防腐、耐摩擦、抗氧化等表面防 护。鉴于涂层在工业领域应用越来越广泛,现代工业对装备可靠性、 安全性的要求不断提高,涂层强度和失效分析越来越重要,同时涂层 的力学性能的表征可以及时对工业镀覆技术进行改进提供参考。

声速和材料显微结构参数中一部分参数相关,超声波声速在材料 中传播速度的变化通常与显微结构特性和力学性能的显著变化相关, 同时由于超声波具有穿透能力强,频带宽和可实现对材料的无损检测 等优势,超声波速法在材料质量检测和评价中得到广泛应用。

当超声波波长小于基底厚度时,视为半无限域基体涂层结构材 料,通过小尺寸材料弹性常数测量系统,根据V(f,z)频域分析法, 采用线聚焦PVDF探头,进行频散分析。将纵波和表面波波速的同时 测量,可实现材料的弹性系数无损检测。经由一种新的反演算法,利 用该线聚焦超声探头来获取半无限域基体涂层结构材料的表面波波 速。该方法基于单纯形法诱导目标函数于频散特征方程的系数矩阵行 列式中,半无限域基体涂层结构材料的弹性性质与试样密度都可以通 过声学性质及所测密度来获得。

发明内容

本发明的目的是为了解决半无限域基体涂层结构材料波速提取 的问题,提出一种先进的材料波速提取方法。

为实现上述目的,本发明采用的技术方案为一种基于单纯形法的 半无限域基体涂层结构材料弹性性质获取方法,该方法的实现步骤如 下,

步骤1确立获得表面波波速反演的仿真目标函数

通过实验获得的频率f与表面波波速c与改变纵波波速CL,横波 波速CT,涂层密度ρ,涂层厚度h后的残差值进行叠加,经过叠加后 的目标函数在某一组CL,CT,ρ,h值时最小,此时可反演出其横、纵波波 速、镀层密度、镀层厚度。

Πs=Σj=1Ns[K(fjs,cjs,CLg,CTg,ρg,hg)]

步骤2搭建表面波波速测试系统

为了方便散焦步进测量,搭建了一套进行散焦步进测量的测试系 统,如图1所示。该测试系统包括试样1、水槽与水2、换能器3、 移动平台4、脉冲激励/接收仪5、示波器6、GPIB总线7、PXI总控 制系统8、移动伺服马达9、旋转轴10;其中,在移动平台4下面安 装换能器3,换能器3与脉冲激励/接收仪5相连,脉冲激励/接收仪5 与示波器6相连,示波器6通过GPIB总线7与PXI总控制系统8相 连,PXI总控制系统8与移动伺服马达9相连,同时PXI总控制系统 8与旋转轴10相连。

步骤3聚焦面数据采集。

将被测试样置于换能器的聚焦面,脉冲激励/接收仪5在发出一个 带宽为10-200MHz的脉冲后转换为接收状态,当接收到反射信号后, 将信号传输进示波器6,示波器的采样频率为fS,fS为0.5-5GHz,采 样点数为Ns,Ns的取值范围为10000-100000点;经过示波器的低通 滤波后,通过GPIB总线7存储进PXI总控制系统8。

步骤4换能器散焦测量。

将换能器垂直向下移动一个距离Δz0,Δz0的取值范围为1-50μm, 待移动完成后进行数据采集,采样频率为fS,采样点数为Ns;采集结 束后再将换能器垂直向下移动Δz0进行数据采集,如此循环往复,共 移动距离z,z的取值范围为2-20mm,因此将得到M组电压数据, M由z与Δz0共同决定,为40-20000组。

步骤5散焦数据时域傅里叶变换。

将所有数据沿散焦距离排列好,对测得的数据进行时域傅里叶变 换:

Ai[k]=Σn=0Ns-1xi[n]e-j2πnk/Ns

其中:Ai为时域傅里叶变换后的频谱值,xi代表一组电压数据, i=0,1,2…M-1,k=0,1,2…Ns-1,j代表虚部。

步骤6时域傅里叶变换后数据空间傅里叶变换。

为了得到精确的振荡周期Δz,需要对时域傅里叶变换的结果再 进行沿散焦距离方向的空间傅里叶变换,将散焦距离z变换至z-1域:

Bi[k]=Σm=0M-1Am[k]e-j2πmi/M

其中:Bi为空间傅里叶变换后的频谱值,Am代表沿散焦方向的时域傅 里叶变换的频谱值,i=0,1,2…M-1,k=0,1,2…Ns-1,j代表虚部。沿 z-1域的曲线峰值即为振荡周期Δz的倒数。

步骤7振荡周期多模态追踪。

对多模态中散焦测量起始处及终止处对每个模态的极大值通过 最优路径进行连续追踪,即可得到连续的z-1值,其倒数即为Δz。

步骤8表面波波速提取。

V(z)曲线理论,可根据如下公式进行波速的计算:

vLamb=vw·[1-(1-vw2·f·Δz)2]-1/2

将水的超声波波速vW、每一个极大值对应的频率f与Δz代入其中: Δz为V(z)曲线振荡周期,vw为水中的超声波波速,f为换能器的激 励频率,vLamb为材料的兰姆波波速;测量被测材料的V(z)曲线振荡周 期是波速提取的关键。

通过实验数据获得其频散波速曲线。

步骤9单纯形法确定波速值

将试验测得数据带入目标函数中,通过单纯形法改变CL,CT,ρ,h四 个变量值,在其保证目标函数值最小即接近0时,即确定各值。

本发明具有以下优点:(1)不必再将整幅理论与实际曲线的波速 平方的残差总和最小值都计算出来。(2)不必把所有的理论频散曲线 都确定,这种方式耗时且数据量冗长,单纯形的方法可以加速并化简 程序。

附图说明

图1为为本方法的实施方法流程图。

图2为散焦测量系统示意图。

图3为表面波传播示意图。

图4为Lamb波在半无限域基体涂层结构材料中传播示意图。

图5为聚焦面时域波形图。

图6为7.5MHz频率下V(z)振荡曲线图。

图7为7.5MHz频率下z-1域曲线图。

图中:1、试样,2、水槽与水,3、换能器,4、移动平台,5、脉冲 激励/接收仪,6、示波器,7、GPIB总线,8、PXI总控制系统,9、 移动伺服马达,10、旋转轴。

具体实施方式

以下结合具体实例对本发明的内容做进一步的详细说明:

步骤1确立仿真目标函数

通过对各项同性材料试件内部部分波线性组合,将其定义在一个 关于边界条件的特征方程中,从而得到频散关系。薄膜厚度为d的半 无限域基体试样置于笛卡尔坐标系下,如图4所示,x1-x3平面是 Lamb波传播的矢面,而坐标x1代表波的传播方向。全局坐标 (x1,x2,x3)为下面的理论服务。

步骤1.1:根据运动方程:σij,i+ρfj=ρüj(可忽略重力ρfj影响); 及本构方程:σij=cijklεkl;几何方程:推知控制波的传播 方程ci·,=ρ..

σij,i=ρu..jσij=cijklϵklϵij=12(ui,j+uj,i)σij=cijklϵkl=cijkl12(uk,l+ul,k)=cijkl·uk,lcijkl·uk,li=ρu..j.

其中σij为应力张量,ρ是密度,uj是质点位移向量,t为时间, 圆点代表时间的不同,指针从逗号后开始代表空间坐标的不同,重复 指针的合计就如张量标志自动的被设定在自此之后。cijkl是弹性系数 在之后可以将其收缩在CIJ系数中,εkl是应变张量,其与质点位移有 关。

根据弹性动力学理论,现有一个平面谐波沿着x1传播,其角频率 为ω,该波的相应物理量在x2方向上独立,波的形式uk如下:

(u1,u2,u3)=(U,V,W)·e(x1+αx3)·e-jωt

其中ξ为波数,而α则表示一个未知系数,(ξ,0,αξ)为该波的传播 向量,(U,V,W)为该平面波相应的振幅。将uk带入控制波的传播方程 cijkl·uk,li=ρüj得到:[Kij]3×3·UVW=000,该线性代数形式[Kij]3×3其 中:

K11=C11+2·C15·α+C55·α20·c2

K22=C66+2·C46·α+C44·α20·c2

K33=C55+2·C35·α+C33·α20·c2

K12=K21=C16+(C14+C56)·α+C45·α2

K13=K31=C15+(C13+C55)·α+C35·α2

K23=K32=C56+(C36+C45)·α+C34·α2

为了使ω,ξ存在值,α的值要使K矩阵行列式为零,也就是说α 可以作为K矩阵的特征值,向量(U,V,W)是其相对相应的特征向量。

在半无限域基体涂层结构材料的波场中可以利用以下的弹性系 数关系。首先将CIJ定义在一个半无限域基体涂层结构材料中,得到 的弹性系数为:

CIJ(1)=λ(1)+2μ(1)λ(1)λ(1)000λ(1)+2μ(1)λ(1)000λ(1)+2μ(1)000μ(1)00Sym.μ(1)0μ(1)

CIJ(2)=λ(2)+2μ(2)λ(2)λ(2)000λ(2)+2μ(2)λ(2)000λ(2)+2μ(2)000μ(2)00Sym.μ(2)0μ(2)

λ与μ为拉梅系数,我们定义CIJ=cijkl,其中ij→I或J,定义 11→1,22→2,33→3,23或32→4,31或13→5,12 或21→6。

获得α的表示式,从而可知该波的传播方式。通常,在各向同性 材料中,α的值有四个。

α1,α3=±(ω/cpξ)2-1α2,α4=±(ω/ctξ)2-1

cp为纵波波速,ct为横波波速。将不同的α值下对应的波的形式 累加,得到以下位移、应力公式来表示波场。其中Aq为未知系数。

可解得特征值,这些特征值是给定ω,ξ情况下波运动的基本解, 波在半无限域基体涂层结构材料中的位移场可以通过特征值的线性 组合表示,不同是,基底中没有通过半无限域边界反射回来的波。可 得到以下位移、应力公式来表示波场:

uk(1)=Ukq(1)·Aq(1)·e(x1+αq(1)x3)·e-jωtq=1,2,3,4

uk(2)=Ukq(2)·Aq(2)·e(x1+αq(2)x3)·e-jωtq=1,2

σij(1)=cijkl(1)·uk,l(1)

q=1,2,3,4

σij(2)=cijkl(2)·uk,l(2)

σ31(2)σ32(2)σ33(2)=D1q(2)D2q(2)D3q(2)·Aq(2)·e(x1+αq(2)x3)·e-jωtq=1,2

获得边界条件的系数矩阵,

在x3=d为自由载荷状态(T33=T13=0),根据位移与应力在x3=0处 的边界条件,得到以下关系式:

σ31(1)|x3=d=0

σ33(1)|x3=d=0

u1(1)|x3=0=u1(2)|x3=0

u3(1)|x3=0=u3(2)|x3=0

σ31(1)|x2=0=σ31(2)|x2=0

σ33(1)|x2=0=σ33(1)|x2=0

并最终可得:D11(1)D12(1)D13(1)D14(1)00D31(1)D32(1)D33(1)D34(1)00U11(1)U12(1)U13(1)U14(1)-U11(2)-U12(2)U31(1)U32(1)U33(1)U34(1)-U31(2)-U32(2)D11(1)D12(1)D13(1)D14(1)-D11(2)-D12(2)D31(1)D32(1)D33(1)D34(1)-D31(2)-D32(2)A1(1)A2(1)A3(1)A4(1)A1(2)A2(2)=0

Ukq(1)=11110000α1(1)1-α2(1)α3(1)1-α4(1)Uq(1)Vq(1)Wq(1)

此时:Ukq(2)=1100α1(2)1-α2(2)Uq(1)Vq(1)Wq(1)

根据:Wq=ρc2-C11-C55αq2(C13+C55)αq

可得:

W1=ρc2-C11-C55α12(C13+C55)α1=(α12+1)ρcp2-(λ+2μ)-μα12(λ+μ)α1=(α12+1)(λ+2μ)-(λ+2μ)-μα12(λ+μ)α1=α1

W2=ρc2-C11-C55α22(C13+C55)α2=(α22+1)ρcs2-(λ+2μ)-μα22(λ+μ)α2=(α22+1)μ-(λ+2μ)-μα22(λ+μ)α2=-1α2

W3=α3=-α1

W4=-1α4=1α2

经由Dpq表达式,整理、去除为零部分后可得:

D1q=C55(Wqq)

D3q=C31+(C33Wqq

由边界条件[M]6×6行列式可得六个方向的频散曲线,再得到其频散 特征函数K,为了得到非零解,函数K的值应为0,CL,CT,ρ,h都为已 知,改变f,c的值使目标函数值最小保证其尽量趋近于0,频散曲线, 频散特征函数,仿真的目标函数Πs分别如下:

σ31|x3=dσ33(1)|x3=du1(1)|x3=0-u1(2)|x3=0u3(1)|x3=0=u3(2)|x3=0σ31(1)|x3=0=σ31(2)|x3=0σ33(1)|x3=0=σ33(1)|x3=0=[M]6×6A1A2A3A4A5A6={0}

K(f,c;CL,CT,ρ,h)=log10([M]6×6)

Πs=Σj=1Ns[K(fjs,cjs,CLg,CTg,ρg,hg)]

步骤1.2:确立波速提取的公式。

在单频激励/接收的情况下,图3所示的漏表面波传播示意图中, 上表面的直接反射回波I传播的时间与漏表面波L的传播时间分别 为:

t1=2(R-Δz)vw

t2=2(R-ΔzcosθSAW)vw+2·Δz·tanθSAWvSAW

其中R为聚焦半径,Δz为散焦距离,vw为水的超声波波速,θSAW为产生表面波的瑞利角,vSAW为材料的表面波波速。因此两者的时间 差为:

Δt=t2-t1=2(1-cosθSAW)vw·Δz

即:

cosθSAW=1-vw·Δt2·Δt

将Snell定律:

sinθSAW=vwvSAWθSAW=sin-1(vwvSAW)

代入(4)后,可得:

vwvSAW=1-(1-vw2.ΔtΔz)2

此时如果Δz恰为一个V(z)曲线的振荡周期时,1/Δt则为换能器的 激励频率f。如果Δz能够确定,便可使用如下公式进行表面波波速的 计算:

vSAW=vw·[1-(1-vw2·f·Δz)2]-1/2

因此,测量被测材料的V(z)曲线振荡周期成为波速提取的重点。

步骤2:搭建测试系统。

为了方便散焦步进测量,搭建了一套进行散焦步进测量的测试系 统,如图2所示。该测试系统主要包括:试样1、水槽与水2、换能 器3、移动平台4、脉冲激励/接收仪5、示波器6、GPIB总线7、PXI 总控制系统8、移动伺服马达9、旋转轴10。其中,在移动平台4下 面安装换能器3,换能器3与脉冲激励/接收仪5相连,脉冲激励/接 收仪5与示波器6相连,示波器6通过GPIB总线7与PXI总控制系 统8相连,PXI总控制系统8与移动伺服马达9相连,同时PXI总控 制系统8与旋转轴10相连。

步骤3:聚焦面数据采集。

以长方体碳化钨为被测试样,其尺寸为40mm×40mm×10mm,将 换能器3聚焦到试样的上表面,通过脉冲激励/接收仪5在发出一个 带宽为10-200MHz的脉冲后转换为接收状态,当接收到反射信号后, 将信号传输进示波器6,示波器的采样频率fS=2.5GHz,采样点数 Ns=10000。经过示波器的低通滤波后,通过GPIB总线7存储进PXI 总控制系统,聚焦面的时域波形如图5所示。

步骤4:散焦测量。

将换能器朝试样方向移动Δz0=10μm,待移动完成后进行电压数据 采集,采集结束后再将换能器朝试样方向移动Δz0=10μm进行数据采 集,采样频率fS=2.5GHz,采样点数Ns=10000,如此循环往复,共移 动4mm,因此将得到400组电压数据,将聚焦面的电压数据包含在 内共得到M=401组电压数据。

步骤5:时域傅里叶变换。

将测得的数据进行时域傅里叶变换。

Ai[k]=Σn=0Ns-1xi[n]e-j2πnk/Ns

其中:Ai为时域傅里叶变换后的频谱值,xi代表一组电压数据, i=0,1,2…M-1,k=0,1,2…Ns-1,j代表虚部,Ns=10000。例如,所得 Ai[k],i=0,1,2…M-1,k=0,1,2…Ns-1。7.5MHz频率下的振荡曲线如 图6所示。

步骤6:空间傅里叶变换。

为了得到精确的振荡周期Δz,需要对时域傅里叶变换的结果再进 行沿散焦距离方向的空间傅里叶变换,将散焦距离z变换至z-1域:

Bi[k]=Σm=0M-1Am[k]e-j2πmi/M

其中:Bi为空间傅里叶变换后的频谱值,Am代表沿散焦方向的时 域傅里叶变换的频谱值,i=0,1,2…M-1,k=0,1,2Ns-1,M=401,j 代表虚部。所得Bi[k],i=0,1,2…M-1,k=0,1,2…Ns-1。例,7.5MHz 频率下z-1域的曲线如图7所示。

步骤7:模态追踪。

对2.5-22.5MHz范围内的峰值进行追踪,即可找出该频率段连续 的Δz值。

步骤8:波速提取。

将水中的超声波波速vW=1500m/s,每一个峰值对应的频率与Δz带 入公式,即可得到该频率段内连续的表面 波波速。

步骤9:单纯形法获得横、纵波波速

通过单纯形法使目标函数残差绝对值最接近零,单纯形法改变 CL,CT,ρ,h使目标函数值的达到所需范围(单纯形法简介)。

单纯形搜索法通过构造单纯形来逼近极小点,每构造一个单纯 形,确定其最高点和最低点,然后通过扩展或压缩、反射构造新的单 纯形,目地是使得极小点能包含于单纯形内。

本方法中未知量有四个(CL,CT,ρ,h),为四维变量问题。

用单纯形搜索法求无约束问题minF(x),x∈Rn的算法步骤如下:

①选取初始单纯形{x0,x1,…,xn},反映系数α>1,紧缩系数 θ∈(0,1),扩展系数γ>1,收缩系数β∈(0,1)及精度ε>0,置k=0;

②将单纯形的n+1个顶点按目标函数值的大小重新编号,使顶点 的编号满足F(x0)≤F(x1)≤…≤F(xn-1)≤F(xn);

③令若{1n+1Σj=0n[F(xj)-F(xn+1)]2}2/1ϵ停止迭代输出 x0,否则转入④;

④计算xn+2=xn+1+α(xn+1-xn),若F(xn+2)<F(x0),转⑤,否则当 F(xn+2)<F(xn-1)时转⑥,当F(xn+2)≥F(xn-1)转⑦;

⑤计算xn+3=xn+1+γ(xn+2-xn+1),若F(xn+3)<F(x0),令xn=xn+3,转②, 否则转⑥;

⑥令xn=xn+2,转②

⑦令xn={xiF(xi)=min(F(xn),F(xn+2))},计算xn+4=xn+1+β(xn-xn+1), 若F(xn+4)<F(xn),令xn=xn+4,转②,否则转⑧;

⑧令xj=x0+θ(xj-x0),j=0,1,…,n,转②

其中,Ns为通过的仿真频散曲线得到的数据点数,这样只 有那些相应残差范围内的K被累加计算起来,就不必再将整幅 理论与实际曲线的波速平方的残差总和最小值都计算出来。并且以前 的方法,必须把所有的理论频散曲线都确定,这种方式耗时且数据量 冗长,因此,单纯形的方法可以加速并化简程序。

去获取专利,查看全文>

相似文献

  • 专利
  • 中文文献
  • 外文文献
获取专利

客服邮箱:kefu@zhangqiaokeyan.com

京公网安备:11010802029741号 ICP备案号:京ICP备15016152号-6 六维联合信息科技 (北京) 有限公司©版权所有
  • 客服微信

  • 服务号