MATLAB实现管道漏磁曲线拟合方法

1.拟合数据的提取

利用ANSYS有限元分析软件，选取半椭圆形缺陷模型。它可以近似地表征出腐蚀缺陷的形态特征，与矩形缺陷模型相比其更接近实际腐蚀缺陷的形状。在建立二维仿真实体模型时，设置不同的缺陷几何参数，半椭圆形缺陷模型可以用半椭圆的长轴和短轴的值来表述缺陷的几何参数，其中，半椭圆的长轴为缺陷的长，半椭圆的短半轴为缺陷的深，如图9-2所示。

图9-2 缺陷的几何参数

图9-2中半椭圆形缺陷的几何参数为缺陷的长度l和缺陷的深度d。由于在ANSYS有限元分析中建立实体模型时，管道壁模型的厚度设置为10mm，因此缺陷的几何参数时是有一定的取值范围的，半椭圆形缺陷的深度的取值范围为1～9mm，长度的取值范围为1～10mm，取值间隔均为1mm，这样就有不同深度和长度的缺陷模型共90个。分别对这90个缺陷模型进行管道漏磁内检测二维有限元仿真分析，共得到90组缺陷漏磁场磁通密度径向分量B_x的数据，这就为下一步缺陷漏磁信号径向分量B_x特征量的提取提供了基础数据。

然后，利用MATLAB软件编写程序，将这90组数据导入MATLAB，以待对其进行数据的处理和分析。编写算法程序，求出这90组缺陷漏磁信号径向分量B_x数据中的特征量B_x峰峰值y和B_x峰峰间距z，分别存入数组y（k）和z（k），k=1，2，…，90。再将半椭圆形缺陷的几何参数即缺陷的深度d和长度l，与这90组漏磁信号径向分量B_x特征数据相对应，也分别存入数组d（k）和l（k）。缺陷的几何参数数据d（k）和l（k），与漏磁信号径向分量B_x特征数据y（k）和z（k），共同组成待拟合数据，截取其中的13组，见表9-1。

表9-1 拟合数据

2.拟合方程的建立

在进行管道漏磁内检测二维有限元仿真数据的最小二乘曲线拟合之前，要根据实际要求建立拟合模型，即建立关于缺陷深度d和长度l的二元函数方程。由于拟合出来的方程要为缺陷量化分析提供数学模型，因此要在符合精度要求的前提下，尽量减少计算量。根据管道漏磁内检测缺陷量化的要求，以及二元函数最小二乘曲线拟合原理，分别建立漏磁信号径向分量B_x特征量的二元一次拟合方程和二元二次拟合方程。

（1）建立一次拟合方程

径向分量B_x峰峰值P₁为

P₁=α₀d+α₁l+α₂ （9-11）

径向分量B_x峰峰间距Q₁为

Q₁=β0d+β₁l+β₂ （9-12）

（2）建立二次拟合方程

径向分量B_x峰峰值P₂为

P₂=a₀d²+a₁d+a₂dl+a₃l+a₄l²+a₅ （9-13）

径向分量B_x峰峰间距Q₂为

Q₂=b₀d²+b₁d+b₂dl+b₃l+b₄l²+b₅ （9-14）

式中，α₀～α₂、β₀～β₂、a₀～a₅、b₀～b₅为待求的拟合系数。

3.拟合方程的求解

根据二元函数最小二乘法曲线拟合的基本原理，对管道漏磁内检测二维有限元仿真得到的径向分量B_x特征数据与缺陷几何参数数据进行曲线拟合，利用MATLAB软件编写程序，实现拟合系数的求解。将求得的拟合系数代入所建立的拟合函数方程，得到径向分量B_x特征的拟合函数方程。因为自变量为两个，即缺陷的深度d和长度l，所以，径向分量B_x峰峰值和径向分量B_x峰峰间距的差值经过最小二乘曲线拟合得到的拟合方程的图形均应为三维图形。

用MATLAB实现二元函数最小二乘曲线拟合的具体步骤如下：

（1）输入待拟合数据 缺陷的深度d（k）、长度l（k）、径向分量B_x峰峰值y（k）、径向分量B_x峰峰间距z（k），均为1行90列的矢量，将其输入MATLAB软件以待拟合。

（2）编写算法子程序 将待拟合数据分别代入方程组（9-7）、式（9-10）中，形成三元一次方程组和六元一次方程组，利用MATLAB软件中的solve函数求解这两个方程组，经过计算得到拟合方程的系数。

（3）绘制拟合曲面 在主程序中调用（2）中的算法子程序，将求解出的各系数代入到相对应的拟合方程中，绘制出径向分量B_x特征方程的三维拟合曲面。

求得的径向分量B_x特征的二元一次拟合方程组为

式中，求解出的各系数均取四位有效数字。

在MATLAB软件中，根据方程组（9-15）绘制出径向分量B_x特征的一次函数方程组的三维拟合曲面，如图9-3所示。(www.xing528.com)

图9-3 漏磁信号径向分量B_x特征的一次函数三维拟合曲面

a）径向分量B_x峰峰值的拟合曲面 b）径向分量B_x峰峰间距的拟合曲面

其中，图9-3a所示为漏磁信号径向分量B_x峰峰值的拟合曲面，图9-3b所示为漏磁信号径向分量B_x峰峰间距的拟合曲面。图中的散点为原始的待拟合数据点。因为漏磁信号径向分量B_x特征的二元一次拟合方程组（9-15）中的自变量d和l均为一次，所以拟合的结果均为平面。

根据图9-3漏磁信号径向分量B_x特征的二元一次三维拟合曲面可以看出，二元一次函数拟合具有一定的局限性，拟合结果不够精确，根据一次拟合方程组所绘制出的平面的拟合效果较差，很多原始数据点都离拟合的曲面较远，只能大致地描绘出原始数据点的分布情况。

根据MATLAB实现二元函数最小二乘曲线拟合的具体步骤，求得的漏磁信号径向分量B_x特征的二元二次拟合方程组为

式中，求解出的各系数均取四位有效数字；d为缺陷的深度；l为缺陷的长度；P₂为二次拟合方程的漏磁信号径向分量B_x峰峰值；Q₂为二次拟合方程的峰峰间距的差值。

在MATLAB软件中，根据方程组（9-16）绘制出漏磁信号径向分量B_x特征的二次函数方程组的三维拟合曲面，如图9-4所示。

图9-4 漏磁信号径向分量B_x特征的二元二次三维拟合曲面

a）径向分量B_x峰峰值的拟合曲面 b）径向分量B_x峰峰间距的拟合曲面

其中，图9-4a所示为漏磁信号径向分量B_x峰峰值的拟合曲面，图9-4b所示为漏磁信号径向分量B_x峰峰间距的拟合曲面。图中的散点为原始的待拟合数据点。因为漏磁信号径向分量B_x特征（峰峰值和峰峰间距）的二元二次拟合方程组（9-16）中的自变量d和l均为二次，所以拟合结果为二次曲面。

根据图9-4漏磁信号径向分量B_x特征的二元二次三维拟合曲面可以看出，二元二次函数拟合的局限性较小，原始数据点都比较贴合二次拟合的曲面，尤其是漏磁信号径向分量B_x峰峰值的三维拟合曲面，少许离开曲面的散点与拟合曲面也都十分接近，更能准确地描绘出原始数据点的分布情况。因此，根据二元二次拟合方程组（9-16）所绘制出的漏磁信号径向分量B_x特征曲面的拟合效果较好。

通过对图9-3所示二元一次三维拟合曲面和图9-4所示二元二次三维拟合曲面进行对比，可以得知，二元二次函数曲线拟合比二元一次函数曲线拟合的拟合效果更好，能更准确地描绘出原始拟合数据点的分布情况。这与最小二乘法曲线拟合原理中的曲线拟合函数方程的次数越高、拟合效果越好的特性是一致的。

4.拟合误差的分析

针对管道漏磁内检测二维有限元仿真的缺陷处漏磁信号径向分量B_x特征的二元一次函数拟合结果和二元二次函数拟合结果，进行拟合误差的分析，从误差数值的角度进一步对拟合结果进行比较和分析。

采用误差分析中的SSE和RMSE来对拟合结果进行拟合误差的分析。

SSE（和方差）是拟合结果数据和原始拟合数据对应点误差的平方和，其计算公式为

式中，y_i和分别为拟合结果数据和原始拟合数据。

SSE越接近于0，说明所建立的模型拟合效果越好。

RMSE（均方根）是拟合结果数据和原始拟合数据对应点误差平方和均值的平方根，也称为拟合标准差，其计算公式为

根据管道漏磁内检测半椭圆形缺陷二维有限元仿真分析所得到的漏磁信号径向分量B_x特征值数据，可以知道，漏磁信号径向分量B_x峰峰值和峰峰间距值的原始拟合数据点均为90个。因此，在进行误差分析的计算时，式（9-17）和式（9-18）中的n=90。

根据式（9-17）和式（9-18），对管道漏磁内检测半椭圆形缺陷二维有限元仿真分析所得到的漏磁信号径向分量B_x特征拟合进行误差分析，分析结果见表9-2。

表9-2 拟合结果数据和原始拟合数据的误差

由误差分析结果可以明显地看出：二元二次函数拟合结果的误差比二元一次函数拟合结果的误差要小，SSE值和RMSE值都更接近于0。二元二次函数曲线拟合模型的拟合效果比较好，拟合模型可靠，能够得到很好的拟合结果，可以满足二维缺陷量化分析的要求。