CN109933922A - 一种矩形金属薄板直线穿透裂纹的识别方法 - Google Patents

Abstract

本发明为一种矩形金属薄板直线穿透裂纹的识别方法，应用于机械装备***的金属薄板裂纹损伤识别，本方法首先用传感器测量待测裂纹金属板受力后的侧边上某些点的位移信息，根据实测板的实际尺寸建立仿真板坐标系，并用拉丁超立方采样法生成近似模型位置值；然后用水平集法确定仿真板等参单元的类型，并用扩展有限元法修改不同等参单元的加强函数以计算仿真板的不同等参单元的位移模式。接着依次计算各位置值对应点的位移，并以此为基础建立径向基近似模型(RBF)。最后根据实测板受力后的位移信息使用粒子群算法在近似模型中识别实测板裂纹的首尾坐标。本发明所用近似模型，降低了仿真问题中的计算量，减少计算成本，耗时短结果精确。

Description

一种矩形金属薄板直线穿透裂纹的识别方法

技术领域

本发明属于机械损伤监测、诊断技术领域，具体涉及一种矩形金属薄板直线穿透裂纹的识别方法。

背景技术

板类结构是机械单元主要构成之一，由于板材重量轻、耗材少，被广泛地用于多种实际结构中，如火箭尾翼、太阳能板等。如果这些组合结构中的板出现裂纹，这些裂纹将导致整体结构的强度与振动特性等方面的变化，从而有可能产生结构破坏以至产生重大事故，造成巨大损失。但早期的初始微小裂纹不易被发现，容易被忽略。脆断发生均是构件内部存在的宏观尺寸(肉眼可见的、0.1mm以上)的微小裂纹扩展造成的。金属薄板的裂纹损伤识别，从原理上讲,与其他传统检测方式类似。首先由传感器测量出缺陷对外界激励的响应,再使用适当的方法进行识别分析。不同于传统的检测方式，现有的识别方式是根据金属薄板的实测结果，通过数值分析方法建立缺陷的反演模型，进而识别出结构内的缺陷分布。这是一种行之有效的方式。但是已有的方法往往只能寻优平行于板边界的裂纹，算法效率低，在反问题寻优中容易陷入局部最优解，算法改进复杂，且寻优耗时长。

发明内容

本发明的目的在于提供一种矩形金属薄板直线穿透裂纹的识别方法，应用于机械装备***的金属薄板损伤裂纹识别，本方法首先用传感器测量待测裂纹金属板受力后侧边某些点的位移信息，根据实测板的实际尺寸建立仿真板坐标系分布，并用拉丁超立方采样法生成近似模型位置值；然后用水平集法确定仿真板等参单元的类型，并用扩展有限元法修改不同等参单元的加强函数以计算仿真板的不同等参单元的位移模式。接着依次计算各位置值对应点的位移，并以此为基础建立径向基近似模型(RBF)。最后根据实测板受力后的位移信息使用粒子群算法在近似模型(RBF)中识别出实测板裂纹的首尾坐标位置。

本发明设计的一种矩形金属薄板直线穿透裂纹的识别方法主要包括如下步骤：

步骤1、实测裂纹金属板受力变形后的位移信息

实测板为厚度均匀的金属薄板，其有一条直线穿透裂纹。实测板的长为2L，宽为2W，L/W为1/1至2/1。板厚小于L或W的1％。因为板的厚度小于小于板的边长，且受到的力为垂直于厚度方向的静态拉伸力或切向力，分析过程中忽略厚度方向剪切力的影响。

检查是否有与实测板相同的仿真板的近似模型存在，若无，进入步骤2；若已有，直接进入步骤6；

步骤2、建立仿真板的直角坐标系，并生成位置值；

步骤3、使用水平集法确定仿真板的等参单元类型；

步骤4、使用扩展有限元法修改等参单元的加强函数并计算仿真板的总***移模式，按步骤2获取的位置值一一计算各位置值对应的仿真板的总***移；

步骤5、构建径向基近似模型；

步骤6、使用粒子群算法在近似模型中识别出实测板上裂纹的位置。

所述步骤1实测裂纹金属板受力变形后的位移信息的具体执行步骤为：

所述实测板的下边固定，在板的两条侧边上安装N个位移传感器，N为4、6或者8。两侧边安装的传感器型号规格和数量相同，且位置对称。在板的上边施加水平切向力或者垂直拉伸力，记录各传感器所测得的位移r^*为各传感器的编号，r^*∈(1,N)。该位移用此点的直角坐标位置表示，

即

所述步骤2建立仿真板的直角坐标系，并生成位置值的具体执行步骤为：

步骤21根据实测板的实际尺寸建立仿真板直角坐标系并划分网格；

以与步骤1的实测板长、宽和厚相同的矩形板为仿真板，直角坐标系原点设置在仿真板的中心，二坐标轴平行于矩形两边。仿真板矩形4顶点的坐标依次为X_tip1(-L，W)，Y_tip1(-L，-W)，X_tip2(L，W)，Y_tip2(L，-W)。

将横坐标分成nx份，纵坐标分成ny份，nx，ny的取值范围是30～50。

此时，将仿真板4顶点的坐标和nx，ny带入matlab(矩阵实验室计算软件)中生成仿真板的网格。仿真板得到nx×ny个网格单元，即矩形等参单元，每个矩形等参单元的4个顶点称为结点，每个结点在所述仿真板的直角坐标系中的坐标为u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)。

步骤22用拉丁超立方采样法生成位置值；

在裂纹识别过程中，要在仿真板众多模拟裂纹中找到与实测板上裂纹首尾最接近的位置，需涉及到多条不同的裂纹，计算很多的位置。且每生成一个裂纹位置，都需要进行一次扩展有限元建模分析。为了节约计算成本，减少计算量，本发明引入近似模型。

生成近似模型首先要确定位置值，过多的位置值会造成浪费，不利于节省成本，过少的位置值又不能准确的反映源模型的特性，本发明采用拉丁超立方采样获取位置值。

拉丁超立方采样如下：

221、采样区间为[-L，L]，位置值的组数M为20～50。

222、在X轴上取两个坐标值X_n,，X_n+1；同样在Y轴上取两个坐标值Y_n，Y_n+1，所述n取值为0～(M-1)，(X_n,Y_n；X_n+1,Y_n+1)组合成位置值每个包含两个坐标，这两个坐标就是一条裂纹的首尾坐标，即：

位置值表示某条裂纹的一组首尾坐标；

223、重复步骤222M次，然后在matlab中筛选数据，以保所选取的M组位置值没有重复。

224、在步骤223所取的M组位置值中随机取一组值作为首坐标(X₁,Y₁)和尾坐标(X₂,Y₂)，代入仿真板直角坐标系生成一条直线，该直线模拟一条裂纹。

以下步骤3和4模拟单条直线裂纹，求仿真板等参单元的位移模式。

所述步骤3使用水平集法确定仿真板的等参单元类型的具体执行步骤为：

步骤31水平集法划分仿真板区域

水平集法包括裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)。在含有一直线裂纹的仿真板上，裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)相互垂直，Ψ(x)＝0和φ^k(x)＝0的3条直线延长到仿真板的边界将仿真板划分为4个不同区域：

Ψ(x)>0,φ^k(x)>0；

Ψ(x)<0,φ^k(x)>0；

Ψ(x)>0，φ^k(x)<0；

Ψ(x)<0,φ^k(x)<0；

Ψ(x)>0、φ^k(x)>0和Ψ(x)<0、φ^k(x)>0均分别为2个子区域，对称地位于裂纹首尾两端外侧。

步骤32水平集法判断裂纹单元；

计算仿真板上各个等参单元的4结点uS到裂纹面Ψ(x)的距离Г和4个结点uS到波前φ^k(x)的距离γ，等参单元的类型判断如下：

321、计算所得Г的最大值Г_max与Г的最小值Г_min之积，当所得结果小于0，且φ^k(x)大于0，此时该等参单元为裂纹贯穿单元；

322、计算所得γ的最大值与计算所得γ的最小值之积，当所得结果小于0，且φ^k(x)大于0，该等参单元为裂尖单元；

323、不符合321和322条件的是常规等参单元。

所述步骤4、使用扩展有限元法修改等参单元的加强函数并计算仿真板的总***移模式的具体执行步骤为：

用扩展有限元法对等参单元进行修改，保持仿真板整***移函数的连续性。对于不同的等参单元，用不同的加强函数进行修改。

多数情况只有边界条件或者等参单元结点的坐标，为求解等参单元内任意一点的位移，用一个函数对等参单元内任一点近似插值，模拟表示其位移。用插值形函数求位移u(x)的公式如下：

其中，u_S为等参单元四个结点坐标值u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)；N_S(x)为所选择的插值形函数；N^m为等参单元的结点总集。

所述步骤4使用扩展有限元法修改等参单元的加强函数并计算仿真板的总***移模式的具体执行步骤如下：

步骤41加强裂纹贯穿单元的阶跃函数及裂纹贯穿单元的位移模式u(x1)

裂纹贯穿单元的特征是两侧位移不连续，选取阶跃函数H(X)进行加强：

式中，X^*是考察点X在裂纹面上的投影点；n^*为裂纹面的单位外法向。

加强后的裂纹贯穿单元的位移模式u(x1)如下：

式中，N^disc为裂纹贯穿单元的结点集合，j^{^}表示贯穿单元的第j^{^}个结点；N_j^(x)为裂纹贯穿单元的插值形函数；H(X)为阶跃函数；a_j^为结点加强变量，本法a_j^＝2。

步骤42加强裂尖单元的裂尖改进函数及裂尖单元的位移模式u(x2)

裂尖单元通过裂尖改进函数进行加强，裂尖改进函数F_k(X)表示如下：

式中，r和θ为裂纹的两个端点，即首尾坐标(X_n,Y_n,)、(X_n+1,Y_n+1)的局部极坐标，r即裂尖加强区域圆的半径，2r大于单个网格的长度、小于裂纹的长度；θ为30°～45°。

加强后的裂尖单元位移模式u(x2)如下：

式中，N^asy为裂尖加强单元k_o的集合，k表示裂尖单元k_o的第k个结点，为裂尖单元对应的插值形函数，F_k(X)为裂尖单元的改进函数，b_k ^α对应结点加强变量，当为一条裂纹时设定b_k ^α＝8。

步骤43常规等参单元的位移模式u(x3)

常规等参单元的位移模式如下：

式中，N^c为常规等参单元结点；N_S(x)为常规等参单元对应的插值形函数；

步骤44仿真板的总***移模式

综上仿真板所包含的三种等参单元的位移模式，仿真板的总***移模式为：

f＝u(x1)+u(x2)+u(x3)即

其中，u(x1)表示加强后的裂纹贯穿单元位移模式，u(x2)表示加强后的裂尖单元的位移模式，u(x3)代表常规等参单元的位移模式；

最后将M组位置值分别代入上式(8)中，计算输出M组位置值的仿真板总***移f(f₁，f₂…f_M)。

所述步骤5构建径向基(RBF)近似模型的的具体执行步骤为：

步骤2至4已经计算了一条贯穿直线裂纹的位移模式，若逐条计算每一单直线裂纹的位移模式，很耗时且复杂；

本发明用近似模型通过位置值的输入和输出信息，模拟源模型中所求裂纹位置与设计变量(位移)之间的关系，替代复杂的源模型。

径向基构建的裂纹近似模型用矩阵简化表示，即线性方程组：

A_ijλ＝f_ij (9)

其中， 表示裂纹首尾端的位置坐标，表示和裂纹首尾端的位置坐标有关的向量，是matlab为方便存储赋予的值，与维度相同。f_ij为该裂纹对应的仿真板总***移。i*,j*∈(0,M)。

由步骤223选取的M组位置值，即M条直线裂纹的首尾端位置 作为模拟裂纹的位置坐标将由步骤4.4计算所得的各组位置值对应的仿真板总***移f_ij，代入公式(9)中，通过matlab计算得到一组插值系数λ＝[λ₁,λ₂,…λ_M]。

求得λ后，对于新的模拟裂纹，即由其首尾坐标由公式(9)直接求得该模拟裂纹的位移模式f_ij。

初步完成径向基近似模型的构建。

所述步骤6使用粒子群算法在近似模型中识别出实测板上裂纹位置的具体执行步骤为：

步骤61、建立目标函数

仿真板模拟的任意单条裂纹，用如下向量表示：

其中，(X_n,Y_n)是裂纹首端坐标；(X_n+1,Y_n+1)是裂纹尾端坐标。在步骤5建立的近似模型内，每一个向量都对应该模拟直线裂纹的仿真位移计算值，用表示：

对于不同的向量所得到的不同。

建立目标函数将步骤1传感器测量的位移值与仿真板上模拟裂纹对应点上的仿真位移值进行比较，寻找二者之间误差最小，即目标函数值最小的模拟裂纹，即识别裂纹位置。

所述目标函数设为：

在二范数中，是步骤1实测板侧边的位移传感器实际测量结果，为仿真板上对与步骤1位移传感器所处位置对应点位移的仿真计算结果。当目标函数最小，即表示仿真板上的此条模拟裂纹对应的位移与实测板上实际测量结果最接近。

步骤62、设置粒子群算法的参数

本法的粒子群算法(PSO)是：由Λ个粒子组成的一个粒子群在d＝4维的设计空间中运动，每个粒子在运动的过程中，结合自己经历的最优位置个体极值Pbest和粒子群体中其他成员经历的最优位置全局极值Gbest来变化自己当前所处的位置。设定所有粒子都根据目标函数来判断其目前位置的好坏，且在运动过程中，每一个粒子都是具有记忆性的，记得其所搜寻到的最佳位置。

粒子群算法首先初始化Λ个粒子组成的一个随机粒子群，然后通过找到目标函数最优解更新迭代个体极值和全局极值。粒子的速度决定下次运动的方向和距离。粒子速度值由以下公式计算:

i∈(1，Λ)；由公式(12)得知，目标函数有四个变量，即每个粒子有四个维度，d＝4。

本法Z_id＝(Z_i1,Z_i2…Z_id-1,Z_id)，表示粒子i在d维的位置，表示粒子i迭代了K^*代；V_id＝(V_i1,V_i2…V_id-1,V_id)表示粒子i在第d维的速度，个体极值是第i个粒子每个维度迭代K^*次经历的最优位置，全局极值表示此时其他全部粒子每个维度迭代K^*次经历的最优位置，ω是惯性权重，取值为从0.9到0.4的线性下降，其大小决定了对当前粒子速度所占的比重；K^*是迭代次数；c1和c2是学习因子，代表粒子具备自我归纳和向群体中优秀粒子学习的能力；r₁和r₂是均衡分布在[0,1]范围内的matlab随机选择的值。

粒子群算法(PSO)的参数设置如下：

621、初始粒子数Λ设为25、50或100；

622、粒子的寻优范围即仿真板的长(-L，L)或宽(-W，W)；

623、粒子单次运动最远距离设置为(0.1L～0.2L)或(0.1W～0.2W)，即最大寻优范围的(0.05～0.1)倍；

624、加速系数C₁＝C₂＝2；

625、随机数r₁和r₂为matlab于[0,1]区间分别随机取值，r₁≠r₂；

626、本法迭代次数K^*设置最高为100，即1≤K^*≤100；

627、惯性权重ω，取值为从0.9到0.4的线性下降。

步骤63、粒子群算法结合近似模型的裂纹识别

631、初始化粒子群的Λ个粒子，此时每个粒子每个维度寻优的初始速度和位置由matlab随机设定,取为中的最优值；由于没有比较对象，此时各粒子的初始位置即为每个粒子每个维度的

设定粒子有记忆功能，且粒子群体数Λ保持不变。

632、在步骤631的第一轮寻优结束后，第二轮C₁＝C₂＝2、r₁，r₂为matlab在[0,1]中重新随机取值，并将 代入公式(13)、(14)中得到新的速度然后Λ个粒子根据新的速度运动到新的位置根据公式(12)计算第i个粒子目标函数的大小，并在中找到最小值

此时再次比较631中第一轮的与632中第二轮的目标函数的大小，以较小者作为更新的

其它粒子按第i个粒子的方式重复Λ-1次，得到更新的个体极值，并更新得第二轮的全局极值

633、将C₁、C₂，r₁、r₂，和步骤632所得的 代入公式(13)、(14)按步骤632进行第三轮寻优并更新与的值；

634、依次循环步骤632和633至第K^*轮，当与的值不变，寻优终止，以所得的与对应的裂纹首尾端位置值作为全局最优解输出；

或者循环次数K^*＝100后，比较此时所得的与与前一轮迭代结果目标函数的大小，所得目标函数较小者对应的裂纹首尾端位置值作为全局最优解输出。

635、全局最优解的验证

将于公式(11)中对应的预测位移 带入式(15)和(16)以判断近似模型的准确性：

其中，y_j为近似模型预测位移的第j个元素；即j＝1,2,…2N，y_j依次表示y_j为步骤1实测位移的第j个元素，即公式(1)中的 为y_j中2N个元素的平均值；根据维度设定SC的最大值SC_max和EISE的最小值EISE_min，本方法粒子空间为四维，设SC_max＝0.9998；EISE_min＝0.0007。

当代入公式(15)所得结果小于SC值、且代入公式(16)所得结果大于EISE的值时，表示近似模型准确，所得的结果可靠，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息。

否则将的值带入公式(9)中求解，得到新的更新近似模型，重新初始化粒子并重复步骤631至634；直到所得的对应的位移代入公式(15)，(16)后，满足小于SC值和大于EISE值，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息。

与现有技术相比，本发明的一种矩形金属薄板直线穿透裂纹的识别方法，具有下列显著优势：1、用水平集法在仿真板上模拟裂纹位置，并引入模拟裂纹不连续面的加强函数，用扩展有限元法修改不同等参单元的加强函数以计算仿真板的不同等参单元的位移模型，建立板类结构仿真近似模型，避免了网格重划分和网格细化带来的影响；2、近似模型(RBF)代替原有的结构参数与动力响应关系，降低了仿真问题中的计算量；3、拉丁超立方采样利用分层，克服了采样样本容量大，效率低等问题,节省成本并准确的反映源模型的特性；4、通过粒子群算法反演寻优求解出与实测位移值相差最小的裂纹首尾位置值，从而准确定量预测出裂纹的相关参数，减少计算成本，耗时短结果精确。

附图说明

图1为本矩形金属薄板直线穿透裂纹的识别方法实施例整体流程图；

图2为本矩形金属薄板直线穿透裂纹的识别方法实施例步骤1实测板的几何尺寸、传感器安装位置及在上边施加拉伸力的示意图；

图3为本矩形金属薄板直线穿透裂纹的识别方法实施例步骤1实测板的几何尺寸、传感器安装位置及在上边施加剪切力的示意图；

图4为本矩形金属薄板直线穿透裂纹的识别方法实施例步骤2和3仿真板的直角坐标系、网络划分以及参照实测板A裂纹位置的三种不同类型等参单元分布示意图；

图5为本矩形金属薄板直线穿透裂纹的识别方法实施例步骤3参照实测板A裂纹位置水平集法划分仿真板的区域；

图6为本矩形金属薄板直线穿透裂纹的识别方法实施例步骤6在近似模型中粒子群算法寻优识别实测板上裂纹位置的流程图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例和附图，对本发明进一步详细说明。

本矩形金属薄板直线穿透裂纹的识别方法实施例，流程图如图1所示，主要包括如下步骤：

步骤1、实测裂纹金属板受力变形后的位移信息

在本实施例中，实测板为2024-T4铝合金板，板的长、宽都为1m、板厚为0.003m且厚度均匀。如图3所示，板的下边固定；板的弹性模量为E＝2E11Pa，泊松比ν＝0.3。

本例金属板中有一条直线穿透裂纹，位置，大小信息未知；本例目标即找到该实测板的裂纹位置；

在实测板的两条侧边上分别安装2个位移传感器，如图2和3黑色方块为位移传感器安装位置，每个侧边上的2个传感器分别位于距上边1/4和3/4边长处，两侧边的传感器相同且位置对称。

如图3所示，实测板的上边施加水平向切向力1.0MPa(对板的上边施加垂直向的拉伸力，与本例情况相似，如图2所示)，实测板在1MPa力下的变形位移，记录各传感器所测得的位移r^*为各传感器的编号，该位移用各点的直角坐标位置表示：

本例通过实验测得， 单位为米。

本例目前未有与实测板相同的仿真板近似模型存在，进入步骤2；

步骤2、建立仿真裂纹板直角坐标系，并生成位置值

步骤21根据实测板的实际尺寸建立仿真板的直角坐标系并划分网格

本例所选仿真板为与实测板相同的矩形板，长、宽都为1m，直角坐标系原点设置在仿真板的中心，二坐标轴平行于矩形两边。本例仿真板四顶点坐标依次为X_tip1(-0.5，0.5)，Y_tip1(-0.5，-0.5)，X_tip2(0.5，0.5)，Y_tip2(0.5，-0.5)。

将横坐标分成nx＝41份，纵坐标分成ny＝31份，如图4所示。

将(X_tip1，Y_tip1，X_tip2，Y_tip2，nx，ny)带入matlab中生成仿真板的网格，本例仿真板中的网格单元为矩形的等参单元，本例得到41×31个等参单元。每个等参单元的4个顶点为结点，每个结点在所述仿真板的直角坐标系中取值，即u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)

步骤22用拉丁超立方采样法生成近似模型位置值

本例采用拉丁超立方获取位置值。

拉丁超立方采样如下：

221、本例采样区间为[-0.5，0.5]，本例位置值的组数M＝20；

222、在X轴上取两个坐标值X_n,，X_n+1；同样在Y轴上取两个坐标值Y_n，Y_n+1组合成位置值本例n取值为0～19。每个包含两个坐标，这两个坐标就是裂纹的首尾坐标，即：

位置值表示裂纹的一组首尾坐标；

223、步骤222.重复20次，在matlab中筛选数据，保证所选取的20组位置值没有重复。

224、在步骤223所取的20组位置值中随机取一组值作为裂纹的首坐标(X₁,Y₁)和尾坐标(X₂,Y₂)，代入仿真板直角坐标系生成一条直线；该直线模拟一条裂纹。

步骤3、使用水平集法确定仿真板的等参单元类型

步骤31水平集法划分仿真板区域

水平集法包括裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)。如图5所示，含有一直线裂纹的仿真板上，裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)相互垂直，将Ψ(x)＝0和φ^k(x)＝0的3条直线延长到仿真板的边界将仿真板划分为4个不同区域：

Ψ(x)>0,φ^k(x)>0；

Ψ(x)<0,φ^k(x)>0；

Ψ(x)>0，φ^k(x)<0；

Ψ(x)<0,φ^k(x)<0；

Ψ(x)>0、φ^k(x)>0和Ψ(x)<0、φ^k(x)>0均分别为2个子区域，对称地位于裂纹首尾两端外侧；

步骤32水平集法判断等参单元类型

通过计算每个等参单元的4结点u_S到裂纹面Ψ(x)的距离Г和4个结点u_S到波前φ^k(x)的距离γ，等参单元的类型判断如下：

321、计算所得Г的最大值Г_max与Г的最小值Г_min之积，当所得结果小于0，且φ^k(x)大于0，此时该等参单元为裂纹贯穿单元；

322、计算所得γ的最大值与计算所得γ的最小值之积，当所得结果小于0，且φ^k(x)大于0，该等参单元为裂尖单元；

323、不符合321和322条件的是常规等参单元；

步骤4、使用扩展有限元法修改等参单元的加强函数并计算仿真板的总***移模式f

为求解等参单元内任意一点的位移，用一个插值形函数对等参单元内任一点近似插值，模拟表示其位移；用插值形函数求位移u(x)的公式如下：

其中，u_S为等参单元四个结点坐标值u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)；N_S(x)为所选择的插值形函数；N^m为等参单元的结点总集；

步骤41加强裂纹贯穿单元的阶跃函数及裂纹贯穿单元的位移模式u(x1)

裂纹贯穿单元的特征是两侧位移不连续，选取阶跃函数H(X)进行加强：

式中，X^*是考察点X在裂纹面上的投影点；n^*为裂纹面的单位外法向；

加强后的裂纹贯穿单元的位移模式u(x1)如下：

N^disc为裂纹贯穿单元的结点集合，j^{^}表示贯穿单元的第j^{^}个结点；N_j^(x)为裂纹贯穿单元的插值形函数；H(X)为阶跃函数；a_j^为结点加强变量，本法a_j^＝2。

步骤42加强裂尖单元的裂尖改进函数及裂尖单元的位移模式u(x2)

裂尖单元通过裂尖改进函数进行加强，裂尖改进函数F_k(X)表示如下：

式中，r和θ为裂纹的两个端点，即首尾坐标(X_n,Y_n,)、(X_n+1,Y_n+1)的局部极坐标，r即裂尖加强区域圆的半径，2r大于单个网格的长度、小于裂纹的长度；θ为30°～45°；

加强后的裂尖单元位移模式u(x2)如下：

式中，N^asy为裂尖加强单元k_o的集合，k表示裂尖单元k_o的第k个结点，N_ko(x)为裂尖单元对应的插值形函数，F_k(X)为裂尖单元的改进函数，b_k ^α对应结点加强变量，本例为一条裂纹设定b_k ^α＝8。

步骤43常规等参单元的位移模式u(x3)

常规等参单元的位移模式如下：

式中，N^c为常规等参单元结点；N_S(x)为常规等参单元对应的插值形函数；

步骤44仿真板的总***移模式

综合上述仿真板所包含的三种等参单元的位移模式，仿真板的总***移模式为：

f＝u(x1)+u(x2)+u(x3)即

将20组位置值分别代入上式(8)中，计算输出20组位置值的仿真板总***移f(f₁，f₂…f₂₀)；

在本实施例中，本步骤建立的仿真板如图4所示，裂纹为0.5m的位于中心的直线穿透裂纹，其中深黑色为裂纹贯穿单元，圆形区域内为裂尖单元，其他为常规等参单元。

步骤5、径向基近似模型的构建

在插值系数λ计算好之后，对于新的裂纹位置信息，只需要给出首尾坐标就可以直接用公式(9)求出位移模式。

径向基构建的裂纹近似模型用矩阵简化表示，即线性方程组：

A_ijλ＝f_ij (9)

其中， 表示裂纹首尾端的位置坐标，表示和裂纹首尾端的位置坐标有关的向量，是matlab为方便存储赋予的值，与维度相同。f_ij为该裂纹对应的仿真板总***移。i*,j*∈(0,20)

由步骤223选取的20组位置值，即M条直线裂纹的首尾端位置 作为模拟裂纹的位置坐标将由步骤44计算所得的各组位置值对应的仿真板总***移f_ij，代入公式(9)中，通过matlab计算得到一组插值系数λ＝[λ₁,λ₂…λ₂₀]。

求得λ后，对于新的模拟裂纹，即由其首尾坐标由公式(9)直接求得该模拟裂纹的位移模式f_ij。

初步完成径向基近似模型的构建。

步骤6、使用粒子群算法在近似模型中识别出实测板上裂纹的位置

本步骤识别流程如图6所示。

步骤61、建立目标函数

仿真板模拟的任意单条裂纹用如下向量表示：

其中，(X_n,Y_n)是裂纹首端坐标；(X_n+1,Y_n+1)是裂纹尾端坐标。在步骤5建立的近似模型内，每一个向量都对应该模拟直线的仿真位移计算值，用表示：

即对于不同的向量所得到的不同。

目标函数设为：

在二范数中，是步骤1实测板侧边的位移传感器实际测量结果，为仿真板上对步骤1位移传感器所处位置对应点的位移仿真计算结果；当目标函数最小，即表示仿真板上的此条模拟裂纹对应的位移与实测板上实际测量结果最接近；

步骤62、设置粒子群算法的参数

粒子群算法首先初始化Λ个粒子组成的一个随机粒子群，然后个体极值和全局极值通过更新迭代找到目标函数最优解。每轮迭代粒子的速度与位置的更新，由公式(13)、(14)计算：

i∈(1，Λ)，d＝4。本法Z_id＝(Z_i1,Z_i2…Z_id-1,Z_id)，表示粒子i在d维的位置，表示粒子i迭代了K^*代；V_id＝(V_i1,V_i2…V_id-1,V_id)表示粒子i在第d维的速度，个体极值是第i个粒子每个维度迭代K^*次经历的最优位置，全局极值表示此时其他全部粒子每个维度迭代K^*次经历的最优位置，ω是惯性权重，取值为从0.9到0.4的线性下降，其大小决定了对当前粒子速度所占的比重；K^*是迭代次数；c1和c2是学习因子，代表粒子具备自我归纳和向群体中优秀粒子学习的能力；r₁和r₂是均衡分布在[0,1]范围内的matlab随机选择的值。

本例粒子群算法(PSO)的参数设置如下：

621、初始粒子数Λ设为25；

622、粒子的寻优范围(-0.5，0.5)；

623、粒子单次运动最远距离设置为0.1；

624、加速系数C₁＝C₂＝2；

625、随机数r₁、r₂matlab于[0,1]区间分别随机取值，本例r₁≠r₂；

626、本例迭代次数K^*设置最高为100，即1≤K^*≤100；

627、惯性权重ω，取值为从0.9到0.4的线性下降。

步骤63、粒子群算法结合近似模型的裂纹识别

631、初始化粒子群的25个粒子，这25个粒子随机分布于仿真板的25个位置；此时每个粒子每个维度寻优的初始速度和位置由matlab随机设定,取为中的最优值；此时各粒子的初始位置即为每个粒子每个维度的

设定粒子有记忆功能，且粒子群体数25保持不变。

632、在步骤631的第一轮寻优结束后，第二轮将C₁＝2，C₂＝2、r₁，r₂为matlab在[0,1]中重新随机取值及 代入公式(13)、(14)中得到新的速度然后Λ个粒子根据新的速度运动到新的位置根据公式(12)计算第i个粒子目标函数的大小，并在中找到最小值

此时再次比较631中第一轮的与632中第二轮的目标函数的大小，以较小者作为更新的

其它粒子按第i个粒子的方式重复24次，得到更新的个体极值，并更新得第二轮的全局极值

633、将C₁、C₂，r₁、r₂，和步骤632所得的 代入公式(13)、(14)按步骤632进行第三轮寻优并更新与的值；

634、依次循环步骤632和633至第43轮，与的值不变，寻优终止，以所得的与对应的裂纹首尾端位置值作为全局最优解输出；

635、全局最优解的验证

将代入公式(11)中得到对应的预测位移 再带入式(15)和(16)以判断近似模型的准确性：

其中，y_j为近似模型预测位移的第j个元素，j＝1,2,…8；y_j依次表示y_j依次表示实测位移 为的平均值；本例设定SC_max＝0.9998和EISE_min＝0.0007；

当代入公式(15)所得结果为0.9981小于SC_max＝0.9998、且代入公式(16)所得结果为0.1017大于EISE_min＝0.0007，本例近似模型准确，所得的结果可靠，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息。如以下表1中实侧板A所示。

此时检查是否还有其它待识别裂纹位置的薄板，若有返回步骤1；本例还有七块相同规格的实测板，故返回步骤1在侧边相同位置安装传感器，检测受力后的各点位移数据，之后直接进入步骤6，用近似模型识别裂纹位置。

如下文表2中实测板H的数据，其裂纹实际位置为(-0.3，0.4)(-0.3，0.2)。如图2所示，其上边受垂直拉伸力，传感器安装情况与实测板A相同，各传感器实测位移结果为

实测板H在步骤633迭代次数K＝100后，得到最优解 由该最优解计算的对应点位移结果为代入(15)(16)验证其可靠性，

此时y_j依次表示(0.0048，0.0213，0.0049，0.021，-0.0062，-0.0147，0.0062，0.0146)

y_j依次表示(-0.0015，-0.0015，0，0，0.001，0.003，0.001，0.003)

是y_j的平均值，此时

此时，SC＝0.999996＞0.9998，EISE＝0.0002006875＜0.0007；

验证说明最优解不可靠，即把和代入(9)中，以21个位置值重新计算，得到新的更新近似模型，重新初始化粒子并重复步骤631至634；最后所得的(-0.3015，0.4003)(-0.2917，0.2027)对应的位移代入公式(15)，(16)后，所得SC＝0.9976和EISE＝0.00109满足小于SC_max＝0.9998和大于EISE_min＝0.0007，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息。

本实施例的识别结果分析

用相同规格的8块实测板对本实施例的近似模型进行验证，8块实测板各含一条直线裂纹，8条裂纹位置不同。表1为实测板上边受水平向剪切力情况下本法识别结果的对照表，表2为实测板上边受垂直向拉伸力情况下本法识别结果与裂纹实际值对照表，两种受力情况下，本实施例所得的识别结果都与实测板上裂纹的首尾端位置坐标十分接近，近90％的位置识别结果误差在2％之内。可见本发明方法准确地描述裂纹位置与结构的位移之间的关系，并由此准确识别了薄板直线穿透裂纹的位置。

表1实测板上边受水平向剪切力情况下本法识别结果与裂纹实际值对照表

表2实测板上边受垂直向拉伸力情况下本法识别结果与裂纹实际值对照表

上述实施例，仅为对本发明的目的、技术方案和有益效果进一步详细说明的具体个例，本发明并非限定于此。凡在本发明的公开的范围之内所做的任何修改、等同替换、改进等，均包含在本发明的保护范围之内。

Claims (7)

1.一种矩形金属薄板直线穿透裂纹的识别方法，主要包括如下步骤：

步骤1、实测裂纹金属板受力变形后的位移信息

实测板为厚度均匀的金属薄板，其有一条直线穿透裂纹；实测板的长为2L，宽为2W，L/W为1/1至2/1；板厚小于L或W的1％；

检查是否有与实测板相同的仿真板的近似模型存在，若无，进入步骤2；若已有，直接进入步骤6；

步骤2、建立仿真板的直角坐标系，并生成位置值；

步骤3、使用水平集法确定仿真板的等参单元类型；

步骤4、使用扩展有限元法修改等参单元的加强函数并计算仿真板的总***移模式，按步骤2获取的位置值一一计算各位置值对应的仿真板的总***移；

步骤5、构建径向基近似模型；

步骤6、使用粒子群算法在近似模型中识别出实测板上裂纹的位置。

2.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤1，具体执行步骤如下：

所述实测板的下边固定，在实测板的两条侧边上安装N个位移传感器，N为4、6或者8；两侧边安装的传感器型号规格和数量相同，且位置对称；在板的上边施加水平切向力或者垂直拉伸力，记录各传感器所测得的位移r^*为各传感器的编号，该位移用此点的直角坐标位置表示，

即

3.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤2具体执行步骤如下：

步骤21根据实测板的实际尺寸建立仿真板直角坐标系并划分网格；

以与步骤1的实测板长、宽和厚相同的矩形板为仿真板，直角坐标系原点设置在仿真板的中心，二坐标轴平行于矩形两边；仿真板矩形4顶点的坐标依次为X_tip1(-L，W)，Y_tip1(-L，-W)，X_tip2(L，W)，Y_tip2(L，-W)；

将横坐标分成nx份，纵坐标分成ny份，nx，ny的取值范围是30～50；

将仿真板4顶点的坐标和nx，ny带入matlab中生成仿真板的网格，仿真板得到nx×ny个网格单元，即矩形等参单元；每个矩形等参单元的4个顶点称为结点，每个结点在所述仿真板的直角坐标系中的坐标为u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)；

步骤22用拉丁超立方采样法生成20～50对近似模型位置值；

生成近似模型首先确定位置值，每组位置值就是一个裂纹的首尾坐标；采用拉丁超立方采样获取位置值，

拉丁超立方采样如下：

221、采样区间为[-L，L]，位置值的组数M为20～50；

222、在X轴上取两个坐标值X_n,，X_n+1；同样在Y轴上取两个坐标值Y_n，Y_n+1，所述n取值为0～(M-1)，(X_n,Y_n,X_n+1,Y_n+1)组合成位置值每个包含两个坐标；

223、重复步骤222M次，在matlab中筛选数据，以保所选取的M组位置值没有重复；

224、在步骤223所取的M组位置值中随机取一组值作为首坐标(X₁,Y₁)和尾坐标(X₂,Y₂)，代入仿真板直角坐标系生成一条直线，该直线模拟一条裂纹。

4.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤3，具体执行步骤如下：

步骤31水平集法划分仿真板区域

水平集法包括裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)；在含有一直线裂纹的仿真板上，裂纹面水平集函数Ψ(x)和波前水平集函数φ^k(x)相互垂直，Ψ(x)＝0和φ^k(x)＝0的3条直线延长到仿真板的边界将仿真板划分为4个不同区域：

Ψ(x)>0,φ^k(x)>0；

Ψ(x)<0,φ^k(x)>0；

Ψ(x)>0，φ^k(x)<0；

Ψ(x)<0,φ^k(x)<0；

Ψ(x)>0、φ^k(x)>0和Ψ(x)<0、φ^k(x)>0均分别为2个子区域，对称地位于裂纹首尾两端外侧；

步骤32水平集法判断裂纹单元；

计算仿真板上各个等参单元的4结点u_S到裂纹面Ψ(x)的距离Г和4个结点u_S到波前φ^k(x)的距离γ，等参单元的类型判断如下：

321、计算所得Г的最大值Г_max与Г的最小值Г_min之积，当所得结果小于0，且φ^k(x)大于0，此时该等参单元为裂纹贯穿单元；

322、计算所得γ的最大值与计算所得γ的最小值之积，当所得结果小于0，且φ^k(x)大于0，该等参单元为裂尖单元；

323、不符合321和322条件的是常规等参单元。

5.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤4具体执行步骤如下：

为求解等参单元内任意一点的位移，用一个插值形函数对等参单元内任一点近似插值，模拟表示其位移；用插值形函数求位移u(x)的公式如下：

其中，u_S为等参单元四个结点坐标值u_s1(ξ¹,η¹)、u_s2(ξ²,η²)、u_s3(ξ³,η³)、u_s4(ξ⁴,η⁴)；N_S(x)为所选择的插值形函数；N^m为等参单元的结点总集；

步骤41加强裂纹贯穿单元的阶跃函数及裂纹贯穿单元的位移模式u(x1)

裂纹贯穿单元的特征是两侧位移不连续，选取阶跃函数H(X)进行加强：

式中，X^*是考察点X在裂纹面上的投影点；n^*为裂纹面的单位外法向；

加强后的裂纹贯穿单元的位移模式u(x1)如下：

式中，N^disc为裂纹贯穿单元的结点集合，j^表示贯穿单元的第j^个结点；N_j^(x)为裂纹贯穿单元的插值形函数；H(X)为阶跃函数；a_j^为结点加强变量，本法a_j^＝2；

步骤42加强裂尖单元的裂尖改进函数及裂尖单元的位移模式u(x2)

裂尖单元通过裂尖改进函数进行加强，裂尖改进函数F_k(X)表示如下：

式中，r和θ为裂纹的两个端点，即首尾坐标(X_n,Y_n,)、(X_n+1,Y_n+1)的局部极坐标，r即裂尖加强区域圆的半径，2r大于单个网格的长度、小于裂纹的长度；θ为30°～45°；

加强后的裂尖单元位移模式u(x2)如下：

式中，N^asy为裂尖加强单元k_o的集合，k表示裂尖单元k_o的第k个结点，为裂尖单元对应的插值形函数，F_k(X)为裂尖单元的改进函数，b_k ^α对应结点加强变量，当为一条裂纹时设定b_k ^α＝8；

步骤43常规等参单元的位移模式u(x3)

常规等参单元的位移模式如下：

式中，N^c为常规单元结点；N_S(x)为常规单元对应的插值形函数；

步骤44仿真板的总***移模式

综合上述仿真板所包含的三种等参单元的位移模式，仿真板的总***移模式为：

f＝u(x1)+u(x2)+u(x3)即

将M组位置值分别代入上式(8)中，计算输出M组位置值的仿真板总***移f(f₁，f₂…f_M)。

6.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤5具体执行步骤如下：

径向基构建的裂纹近似模型用矩阵简化表示，即线性方程组：

A_ijλ＝f_ij (9)

其中， 表示裂纹首尾端的位置坐标，{X_j*}表示和裂纹首尾端的位置坐标有关的向量，是matlab赋予的值，与维度相同；f_ij为该裂纹对应的仿真板总***移，i*,j*∈(0,M)；

由步骤223选取的M组位置值，即M条直线裂纹的首尾端位置 作为裂纹的位置坐标将由步骤44计算所得的各组位置值对应的仿真板总***移f_ij，代入公式(9)中，通过matlab计算得到一组插值系数λ＝[λ₁,λ₂,…λ_M]，径向基近似模型初步建立完成。

7.根据权利要求1所述的矩形金属薄板直线穿透裂纹的识别方法，其特征在于所述步骤6具体执行步骤如下：

步骤61、建立目标函数

仿真板模拟的任意单条裂纹，用如下向量表示：

其中，(X_n,Y_n)是裂纹首端坐标，(X_n+1,Y_n+1)是裂纹尾端坐标；在步骤5建立的近似模型内，每一个向量都对应该模拟直线裂纹的仿真位移计算值，用表示：

对于不同的向量所得到的不同；

建立目标函数将步骤1传感器测量的位移值与仿真板上模拟裂纹对应点上的仿真位移值进行比较，寻找二者之间误差最小，即目标函数值最小的模拟裂纹，即识别裂纹位置；所述目标函数设为：

在二范数中，是步骤1实测板侧边的位移传感器实际测量结果，为仿真板上对步骤1位移传感器所处位置对应点位移的仿真计算结果；当目标函数Π最小，即表示仿真板上的此条模拟裂纹对应的位移与实测板上实际测量结果最接近；

步骤62、设置粒子群算法的参数

粒子群算法首先初始化Λ个粒子组成的一个随机粒子群，然后个体极值和全局极值通过更新迭代找到目标函数最优解；粒子的速度决定下次运动的方向和距离；粒子速度值由以下公式计算:

i∈(1，Λ)，d＝4；

本法Z_id＝(Z_i1,Z_i2…Z_id-1,Z_id)，表示粒子i在d维的位置，表示粒子i迭代了K^*代；V_id＝(V_i1,V_i2…V_id-1,V_id)表示粒子i在第d维的速度，个体极值是第i个粒子每个维度迭代K^*次经历的最优位置，全局极值表示此时其他全部粒子每个维度迭代K^*次经历的最优位置，ω是惯性权重，取值为从0.9到0.4的线性下降；K^*是迭代次数；c1和c2是学习因子，代表粒子具备自我归纳和向群体中优秀粒子学习的能力；r₁和r₂是均衡分布在[0,1]范围内的matlab随机选择的值；

粒子群算法的参数设置如下：

621、初始粒子数Λ设为25、50或100；

622、粒子的寻优范围即仿真板的长(-L，L)或宽(-W，W)；

623、粒子单次运动最远距离设置为(0.1L～0.2L)或(0.1W～0.2W)，即最大寻优范围的(0.05～0.1)倍；

624、加速系数C₁＝C₂＝2；

625、随机数r₁和r₂为matlab于[0,1]区间分别随机取值，r₁≠r₂；

626、本法最高迭代次数K^*设置为100，即1≤K^*≤100；

627、惯性权重ω，取值为从0.9到0.4的线性下降；

步骤63、粒子群算法结合近似模型的裂纹识别

631、初始化粒子群的Λ个粒子，此时每个粒子每个维度寻优的初始速度和位置由matlab随机设定,取为中的最优值；由于没有比较对象，此时各粒子的初始位置即为每个粒子每个维度的

设定粒子有记忆功能，且粒子群体数Λ保持不变；

632、在步骤631的第一轮寻优结束后，第二轮C₁＝C₂＝2、r₁，r₂为matlab在[0,1]中重新随机取值，并将代入公式(13)、(14)中得到新的速度然后Λ个粒子根据新的速度运动到新的位置根据公式(12)计算第i个粒子目标函数的大小，并在中找到最小值

此时再次比较631中第一轮的与632中第二轮的目标函数的大小，将小者作为更新的

其它粒子按第i个粒子的方式重复Λ-1次，得到更新的个体极值，并更新得第二轮的全局极值

633、将C₁、C₂，r₁、r₂，和步骤632所得的代入公式(13)、(14)按步骤632进行第三轮寻优并更新与的值；

634、依次循环步骤632和633至第K^*轮，当与的值不变，寻优终止，以所得的与对应的裂纹首尾端位置值作为全局最优解输出；

或者循环次数K^*＝100后，比较此时所得的与与前一轮迭代结果目标函数的大小，所得目标函数较小者对应的裂纹首尾端位置值作为全局最优解输出；

635、全局最优解的验证

将(λ₂，ζ₂)于公式(11)中对应的预测位移 带入式(15)和(16)以判断近似模型的准确性：

其中，y_j为近似模型预测位移的第j个元素；其中j＝1,2,…2N，y_j依次表示y_j为步骤1实测位移的第2N个元素，即公式(1)中的 为y_j中2N个元素的平均值；根据维度提前设定SC的最大值SC_max和EISE的最小值EISE_min，本法SC_max＝0.9998；EISE_min＝0.0007；

当代入公式(15)所得结果小于SC的值、且代入公式(16)所得结果大于EISE的值时，表示近似模型准确，所得的结果可靠，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息；

否则将的值带入公式(9)中求解，得到新的更新近似模型，重新初始化粒子并重复步骤631至634；直到所得的对应的位移代入公式(15)，(16)后，满足小于SC值和大于EISE值，输出的值即为本方法识别的实测板上裂纹的首尾端坐标位置信息。