CN114216582B - 一种高精度快速温度场重建方法、系统、设备及终端 - Google Patents

Abstract

本发明属于温度检测技术领域，公开了一种高精度快速温度场重建方法、系统、设备及终端，该方法包括：将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。本发明利用了温度场被测区域的尺寸，飞度时间等先验条件；通过改进型FISTA算法，快速求出粗网格点的温度，解决了温度场重建迭代过程中的不稳定现象；然后利用核极限学习机，成功地精确预测出温度场，与传统算法相比，该算法提高了平均相对误差与均方根误差。

Description

一种高精度快速温度场重建方法、系统、设备及终端

技术领域

本发明属于温度检测技术领域，尤其涉及一种高精度快速温度场重建方法、系统、设备及终端。

背景技术

目前，温度作为基本的物理量之一，通过对温度变化进行监测有助于获取多种有用的信息。因此利用声学法温度场重建技术及时准确的获取温度分布信息，在确保工业炉系统高效稳定运行、提高粮食仓储质量、大气参数远距离监测和海底热泉温度分布监测等领域具有广阔的应用前景。

随着科技的发展，我国对于电力需求日益增加，使得发电制造业飞速发展，目前我国己经赶超美国成为世界上电力装机容量第一大强国。火力发电虽然满足了国家对电力的需求，但也对环境造成了极大的污染，因而要提高火力发电厂的燃烧效率以及环保能力。煤炭燃烧时，燃烧炉内部的温度分布会影响煤炭燃烧的充分性因此研究燃烧炉内温度分布对于全世界都具有非常重要的意义。

当前常用的温度检测方法大致可以分为两类：接触式测温和非接触式测温。

接触式测温技术需要通过接触传感器和被测物体以获得传感器位置所在位置的温度。常见的传感器有：热电偶温度传感器，利用介质受热膨胀的原理制造的水银温度计，压力式温度计和双金属温度计等。接触式测温技术简单直观，同时也存在很多缺陷。由于传感器与待测对象的接触不可避免地会影响待测对象的温度分布，会给温度的测量带来一些测量误差。如果待测对象为高温或具有腐蚀特性，则会导致大大缩短传感器的使用寿命，降低温度的测量精度。接触式测温只有当传感器达到热平衡时才能准确测量温度，多以测量中容易受到待测对象自身性质和传热方式的影响，不能保证温度的实时测量。

非接触式测温可以在不与待测物体直接接触的情况下测得对象温度，是一种收集待测区域中热源或燃烧流场中丰富的物理信息，再应用这些信息间接地获取温度分布的测量技术，具有响应速度快、可以维持待测物体的温度场和流场的原貌的优点，当待测区域处于环境恶劣，温度分布情况复杂，测量区域大的情况时，非接触测温更为实用。声学法属于非接触式测量的一种。工业炉温度场监测中所应用的声学高温计采用的就是声学法温度场检测技术。声学法测温原理是根据介质温度和声波在介质中的声速存在一定的单值函数关系，通过合理布置待测区域内的声波传播路径，测量该路径上的声速，最后反演出待测区域的温度分布的方法。声学法测温的主要优点有测量范围广，响应灵敏，测量空间大，可以实现连续实时的温度监测。声学测温由于其具有非侵入式测量的大部分特点，测温范围广，响应速度快，被认为是有广阔发展前景的可视化温度测量技术。

目前超声波测温在锅炉炉膛温度场监测、海底热液口温度测量、仓储粮食温度场测量等领域应用广泛。超声波测温主要依据超声波在介质传播过程中声速与温度的相关关系，通过测量声速变化而推导出被测对象的温度信息。除了具有非接触式测温的特点之外，超声波测温技术还具有测量温度范围广、环境适应性强、实时连续等优点，且能实现空间温度分布测量。

基于声学法的温度场重建算法按照正问题模型的不同可以可分为两类：一类是把被测区域划分为若干个像素，并假设每个像素内的温度是均匀的，然后依据每个像素内声线穿过的长度建立正问题模型，这类正问题模型要求被测区域划分的像素数要小于有效声波路径数；另一类是把每个被测区域划分的像素用一个径向基函数逼近，对待测区域剖分像素数和有效声波路径数的大小关系没有限制。

求解第一类正问题模型的常见方法有代数重建算法和Landweber迭代法。

代数重建法首先给定重建温度场一个初始值，然后计算基于初值的飞渡时间与测量的飞渡时间之间的误差，再将误差用于修正初始值。通过不断地迭代，直到达到要求的精度。代数重建法容易受测量噪声影响，且重建速度与重建结果受初值影响较大。

Landweber迭代法具有良好的稳定性和抗噪能力，但存在数值解过度光滑以及未能利用被重建对象的先验信息等问题。

求解第二类问题的常见方法有截断奇异值分解法。

截断奇异值分解的重建速度快，其首先将超声波路径构成的系数矩阵进行奇异值分解，然后舍弃奇异值逼近零的项，降低小奇异值对飞渡时间测量误差的放大作用。但是由于截断奇异值分解法很容易受到奇异值的影响，导致重建效果不稳定。。因此，亟需设计一种新的高精度快速温度场重建方法、系统。

通过上述分析，现有技术存在的问题及缺陷为：

(1)由于传感器与待测对象的接触不可避免地会影响待测对象的温度分布，会给温度的测量带来一些测量误差；如果待测对象为高温或具有腐蚀特性，则会导致大大缩短传感器的使用寿命，降低温度的测量精度。

(2)接触式测温只有当传感器达到热平衡时才能准确测量温度，多以测量中容易受到待测对象自身性质和传热方式的影响，不能保证温度的实时测量。

(3)代数重建法容易受测量噪声影响，且重建速度与重建结果受初值影响较大；Landweber迭代法数值解过度光滑且未能利用被重建对象的先验信息。

(4)由于截断奇异值分解法很容易受到奇异值的影响，导致重建效果不稳定。。

解决以上问题及缺陷的难度为：

1.根据温度热平衡可知，当传感器与被测物体接触时，会发生热传递，温度高的一方会发生降温，温度低的一方会发生升温。最终达到温度近似相等。因此，传感器与待测对象的接触肯定会影响最终的测量精度，带来测量误差。并且如果待测对象具有腐蚀特性，则会大大缩短传感器的使用寿命。并且温度场状况会实时改变，使用传感器不能保证实时性测量。

2.由于代数式重建法本身固有的缺陷，迭代过程不稳定，并且迭代中存在反复的大矩阵乘积运算，导致迭代时间较长。因此，重建速度与重建结果很难得到改善。

3.由于系数矩阵的奇异值往往呈连续下降状态，很难选择合适的奇异值进行截断。

解决以上问题及缺陷的意义为：

1.保证温度场实时性的测量。温度场的测量在现实生活中大部分运用于电站锅炉当中。对电站锅炉内温度场的实时测量，有利于让我们实时了解炉内的燃烧情况，以及控制炉内的温度，保证整个锅炉系统的正常运作。

2.提高温度场重建速度与重建精度。通过本发明提出的方法可以有效地提高温度场重建迭代过程的稳定性。并且由于引入了加速公式，大大地缩短了迭代时间，使温度场重建速度与重建精度得到了显著的提高。

发明内容

针对现有超声波重建温度场算法存在的问题，本发明提供了一种高精度快速温度场重建方法、系统、设备及终端。

本发明是这样实现的，一种高精度快速温度场重建方法，所述高精度快速温度场重建方法包括：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

进一步，所述高精度快速温度场重建方法包括以下步骤：

步骤一，根据搭建的温度场重建硬件系统待测区域的几何大小、待测温度场的复杂程度，在每条边上按照特定的方式安装一定数量的超声波换能器并计算出超声波传播路径数量；

步骤二，将待测区域划分为一定数量的网格，且划分的网格数量应该大于或者等于超声波传播路径数量，权衡算法效率与重建精度；

步骤三，根据超声波传播路径数量M与划分的网格数量N，计算每一条超声波传播路径在每一个被划分的网格上的长度值，第i条超声波传播路径在第j个网格上的长度为A_i,j，形成系数矩阵A；

步骤四，通过改变超声波换能器的收发状态，测量每一条超声波传播路径的飞渡时间TOF，形成一个飞渡时间向量y；

步骤五，建立目标函数，并将改进型FISTA算法应用于求解公式(33)的x；

步骤六，将经过改进型FISTA算法得到的解x*，经过所述超声波速和气体介质温度的关系表达式，得到待测区域内每一个划分网格的中心点温度向量；

步骤七，对经过改进型FISTA求解的温度向量与对应坐标作为训练样本，带入构建好的核极限学习机进行模型训练；

步骤八，待训练完后将整个被测区域的坐标带入已经训练好的极限学习机预测，得到整个被测区域的连续温度场分布；

步骤九，通过MATLAB绘图工具，将求出的温度场数据绘图，得到等温线图与三维图。

进一步，所述步骤三中的系数矩阵A的表达式为：

所述步骤四中的渡时间向量y的表达式为：

所述步骤五中的目标函数为：

式中，矩阵A为距离矩阵；向量x为超声波传播速度的倒数向量；向量y为飞渡时间向量。

进一步，所述步骤五中的改进型FISTA算法的构建包括：

初始化x⁽⁰⁾＝0,y⁽¹⁾＝x⁽⁰⁾，t⁽¹⁾＝1，由第n次迭代计算结果x⁽ⁿ⁾到第n+1迭代计算结果x⁽ⁿ⁺¹⁾的步骤为：

(1)x⁽ⁿ⁾＝P₊(y⁽ⁿ⁾-W_tA^T(Ax⁽ⁿ⁾-y)/L)

(3)if(F(xⁿ)＞F(x^n-1))

β＝0

else

(4)y⁽ⁿ⁺¹⁾＝x⁽ⁿ⁾+β(x⁽ⁿ⁾-x^(n-1))。

所述步骤五中的将改进型FISTA算法应用于求解公式(33)的x包括：

(1)将预先确定的W矩阵，A矩阵以及y矩阵带入公式(33)；

(2)将初值x⁽⁰⁾设置为一个N*1的向量带入公式进行迭代求解，迭代次数通过调试经验所确定；

(3)通过改进型FISTA算法的迭代求解，得到最优的x*，x*即为超声波在每一个划分网格内的传播速度倒数的向量。

进一步，所述步骤六中的超声波速和气体介质温度的关系表达式为：

式中，c声速，单位m/s；γ气体介质的绝热指数，气体定压热容与定容热容的比值，无量纲量，对空气来说为1.4；R是理想气体常数，8.314^J/(mol·k)；T为气体介质的热力学温度，单位k；M为气体介质的摩尔质量，数值上等于分子量，单位kg/mol，空气为28.8*10^-3kg/mol。

所述待测区域内每一个划分网格中心点温度向量为：

进一步，所述步骤七中的核极限学习机的构建包括：

ELM的隐含层输出矩阵H改写为：

H＝[h(x₁),...,h(x_N)]；

式中，h(x_i)为x_i的非线性映射，N为训练集的个数，x_i为通过改进型FISTA得到的划分的粗网格中心点温度。若该映射的形式未知，则构造如下：

式中，K是核函数。KELM的输出是：

式中，x^*是属于测试集的样本。KELM的输出权值矩阵为：

β^*＝H^T(HH^T+λI^-1)T；

式中，λ为正则化参数，T是训练集的输出向量。

采用径向基函数，RBF函数有以下形式：

式中，δ²为核参数，决定高维数据空间的复杂性；x_i,i∈{1,2,...,N}和x^*分别是来自训练集和测试集的样本；

表示x^*的第m维。

在KELM中，δ²和λ是影响检测性能的两个关键参数；采用粒子群优化算法PSO，得到δ²和λ的最优参数。

将通过改进型FISTA算法迭代求出的划分的7*7个网格的中心点坐标与温度值作为训练样本，带入核极限学习机进行训练；待核极限学习机训练完成后，将整个被测区域的点的坐标带入提前训练好的极限学习机，求出整个温度场的温度分布。

本发明的另一目的在于提供一种应用所述的高精度快速温度场重建方法的高精度快速温度场重建系统，所述高精度快速温度场重建系统包括：

测量区域划定模块，用于将测量区域划定为一定数量的离散网格；

目标函数构建模块，用于构建基于离散网格点温度误差目标函数；

目标函数求解模块，用于基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；

温度场分布获取模块，用于将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

本发明的另一目的在于提供一种计算机设备，所述计算机设备包括存储器和处理器，所述存储器存储有计算机程序，所述计算机程序被所述处理器执行时，使得所述处理器执行如下步骤：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

本发明的另一目的在于提供一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，使得所述处理器执行如下步骤：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

本发明的另一目的在于提供一种信息数据处理终端，所述信息数据处理终端用于实现所述的高精度快速温度场重建系统。

结合上述的所有技术方案，本发明所具备的优点及积极效果为：针对现有的超声波重建温度场算法，本发明受磁共振成像与机器学习领域的启发，提出了一种新型的高精度温度场重建方法。该算法利用了温度场被测区域的尺寸，飞度时间等先验条件；通过改进型FISTA算法，快速求出粗网格点的温度，解决了温度场重建迭代过程中的不稳定现象；然后利用核极限学习机，成功地精确预测出温度场，与传统算法相比，该算法提高了平均相对误差与均方根误差。

与现有技术相比，本发明的技术方案具有以下有益技术效果：

1.保证温度场的完整性。利用最小二乘法重建温度场时，由于其算法本身的性质，只能重建出最外层网格中心点连线构成的一部份区域，无法重建出整个温度场，本发明提出的算法解决了这个问题，重建出的图像覆盖整个被测区域，保证了温度场的完整性。

2.温度场重建稳定性更好。传统的迭代算法，在增加迭代步长时，会导致温度场重建迭代曲线出现波动不稳定，本发明提出的算法通过对加速度系数重置，解决了这一问题，使温度场重建稳定性更好。

3.温度场重建精度更高。本算法提出利用核极限学习机(KELM)替代传统的插值算法，核极限学习机(KELM)是极限学习机(ELM)的改进版本，KELM继承了ELM训练速度快的优势。并且它不容易陷入局部最优，具有更强的非线性拟合能力。通过使用KELM进行被测区域的预测，提高了温度场重建精度。

4.重建时间更短。本发明提出的温度场重建算法中，使用了加速因子，使得迭代步数显著减少，与传统的迭代算法相比，重建时间得到了明显的降低。并且重建时间无限接近于非迭代算法，如最小二乘法等一步求解算法。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对本发明实施例中所需要使用的附图做简单的介绍，显而易见地，下面所描述的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下还可以根据这些附图获得其他的附图。

图1是本发明实施例提供的高精度快速温度场重建方法流程图。

图2是本发明实施例提供的高精度快速温度场重建方法原理图。

图3是本发明实施例提供的高精度快速温度场重建系统结构框图；

图中：1、测量区域划定模块；2、目标函数构建模块；3、目标函数求解模块；4、温度场分布获取模块。

图4是本发明实施例提供的待测区域以及声波飞行路径示意图。

图5是本发明实施例提供的运动示意图。

图6是本发明实施例提供的声波传播时间测量装置示意图。

图7是本发明实施例提供的KELM结构图。

图8是本发明实施例提供的高精度快速温度场重建算法重建图。

图8(a)是本发明实施例提供的二维重建图。

图8(b)是本发明实施例提供的三维重建图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

针对现有技术存在的问题，本发明提供了一种高精度快速温度场重建方法、系统、设备及终端，下面结合附图对本发明作详细的描述。

如图1所示，本发明实施例提供的高精度快速温度场重建方法包括以下步骤：

S101，将测量区域划定为一定数量的离散网格；

S102，构建基于离散网格点温度误差目标函数；

S103，基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；

S104，将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

本发明实施例提供的高精度快速温度场重建方法原理图如图2所示。

如图3所示，本发明实施例提供的高精度快速温度场重建系统包括：

测量区域划定模块1，用于将测量区域划定为一定数量的离散网格；

目标函数构建模块2，用于构建基于离散网格点温度误差目标函数；

目标函数求解模块3，用于基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；

温度场分布获取模块4，用于将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

下面结合具体实施例对本发明的技术方案作进一步描述。

1.温度场重建系统与待测区域设置

整个温度场重建系统主要由收发型声波换能器、数据采集卡、信号调理电路和通道选择电路构成。

声学法温度场重建的一般步骤可简要归纳如下：

1)考虑实际环境情况和温度分辨率的需求，选择合适的声波换能器位置以及选择有效的声学路径；

2)根据已定的声波换能器布局，选择有效的声学路径，根据待测区域的几何参数和重建算法的要求，将待测区域划分为若干粗网格；

3)通过数据采集卡的模拟输入，模拟输出功能来控制数据采集卡的输入与输出，进而驱动声波换能器发射、接收超声波信号，将采集的电压信号在MATLAB上进行信号分析、处理，由互相关时延估计算法得出每一条超声波传播路径对应的超声波飞渡时间TOF；

4)将通过互相关算法得出的每一条超声波传播路径的飞度时间TOF连同待测区域相关的参数信息一起代入温度场重建算法，重建出待测区域的温度场信息，绘制并显示温度场可视化图像，包括温度场俯视图以及三维图。

待测区域以及声波飞行路径示意图如图4所示，在矩形待测区域的每条边上以四等分的方式，一共布置16个超声波换能器在图上以1～16来表示。这种超声波换能器既可以发送超声波，也可以接受超声波，为同体式器件。两两超声波换能器之间的连线表示的是超声波传播路径，该传播路径不可能存在于同一边的超声波换能器之间。一共形成为96条传播路径。由于使用的是同体式超声波换能器，超声波换能器可以在不同时刻用作超声波发送器或者超声波接收器。超声波换能器可以用程序或者是硬件控制。图4中，当把超声波换能器1当作发射器发送信号时，超声波换能器2～16作为接收器依次接受检测声波信号；当把超声波换能器2当作发射器发送信号时，超声波换能器1，超声波换能器3～16被当作接收器依次接收检测信号。按这种方式，在每一次检测的周期里面让每一个超声波换能器依次开启和关闭收发超声波信号。在PC机上对采集到的数据进行软件分析，得出每一条超声波传播路径上超声波的渡越时间。然后将得出的渡越时间数据，待测区域的形状参数和待测区域的划分方式的数据作为基础，用温度场构建算法，得到待测区域内的温度场。

2.声学法测温基本原理

声波在气体介质中传播时，传播速度随着气体温度的变化而变化，因通过测量超声在气体介质中的传播速度，之后再加以计算就可以得到气体介质的温度，这就是超声测温的基本原理。下面将从热力学理论和波动理论对声速与温度的关系进行推理。

2.1波动方程

根假设超声沿X轴方向传播，在原点X＝0处，平面内各点的辰动状态如式(1)所示：

ζ＝ζ_m0cosωt (1)

式中ζ-t时刻的振幅；

ζ_m0-初始时刻的振幅；

声速为C，则时间t后，声波向前传播了距离x＝ct。在距离X处的振动可以写成下式：

分别对式(2)中位移x、时间t同时求二次偏导，得平面波沿X轴方向传播的波动方程，如下式所示：

2.2运动方程

假设气体介质的密度为ρ，则这段气体介质的质量可以表示为ρdx。由波动学理论可知，当声波(平面波)沿X轴方向传播时，则这段介质的右侧将受到右气体介质的作用力

，左侧将受到左方气体介质施加的作用力T，在特殊情况下作用力r即为声压p，运动示意图如图5所示。

如果这段气体介质的振动位移为ζ振动速度为μ，则其运动方程可用下式表示：

即

任何平面波的运动方程都可以用式(5)表示。

2.3速度方程

超声波属于纵波，由声学原理可知，被测空间的气体介质可以视为无限大，此时，声压力T即为声压P。长度为dx的气体介质其体积V在声压作用下会产生变化为dV；当面积不变时，体积V的相对变化表示为dV/V，实际上等于气体介质厚度的变化dζ/ζ。在声压P作用下，若气体介质体积为V，则压力增加dP，介质体积会减少dV，故根据体积弹性膜量的定义可以得到式(6)：

即

由于P＝T，则式(5)可以化为：

公式(8)即为气体介质中纵波的平面方程。

由波动方程(3)和纵波平面方程(8)可得：

对于理想气体，超声的传播过程可以看成是快速绝热过程，即值恒定，可以推出下式：

PV＝B＝γP (10)

结合理想气体状态方程

可以得到超声波速和气体介质温度的关系，如式(11)所示：

式中c--声速，单位m/s；

γ--气体介质的绝热指数，气体定压热容与定容热容的比值，无量纲量，对空气来说为1.4；R--是理想气体常数，8.314J/(mol·k)；

T--为气体介质的热力学温度，单位为k

M--为气体介质的摩尔质量(数值上等于分子量)，单位kg/mol，空气为28.8*10^{- 3}kg/mol。

若气体介质给定，则

为常数(空气通常取20.05)，因此，超声波在气体介质中的速度与气体介质的温度可以看作是单值函数关系。

2.4声波飞度时间测量原理

根据上节速度方程可知，可由声波在气体介质中的传播速度求得气体介质温度但声波的传播速度很难直接测得，因此通过声波传播路径的长度和在该路径上的飞行时间，这两个容易测得的物理量间接的求出声波的传播速度。气体介质的温度可由式(12)求得，单一路径测温装置如图6所示。

式中，d为图6所示的换能器1与换能器2间的距离，单位为m；Δt为两传声器间的声波飞行时间，单位为s；Z为空气中的气体介质常数；T为该路径上的平均温度，单位为K。

在实验中换能器的位置固定，两换能器间的距离d为定值。同一环境下的气体介质常数Z值恒定。结合式(12)可知若想得到准确的温度值T，需要先测得两换能器间的声波飞行时间Δt。在解决声波飞度时间的测量时，互相关函数法因其算法简单、计算量小、易于实现等优点被广泛应用。相关分析可以在时域上比较两信号的相似程度，利用基本互相关算法对声波信号进行互相关运算。声波换能器接收到的两信号分别如式(13)、(14)所示：

x₁(n)＝k₁s(n)+n₁(n) (13)

x₂(n)＝k₂s(n-D)+n₂(n) (14)

对两采集到的信号x₁(n)和x₂(n)做互相关可得：

实际选取的声源信号与噪声信号几乎不存在相关性，

趋近于零。当两信号噪声相关性较小时可以忽略不计，

值趋近于零。则对式(15)整理可得：

R_c(τ)＝R_ss(τ-D) (16)

由自相关函数的性质可知，当τ＝D时R_ss(τ-D)取得最大值。此时最大值点所对应的时间，就是飞度时间值。

2.5温度场重建原理

声学温度场重建系统需要在被测区域周围安装声波换能器，换能器同时具有发射和接收功能。在一个测量周期内使每个换能器依次发射声波，其余所有换能器接收该声波信号，从而形成多条有效声波传播路径。系统将采集到的声波信号经过时延算法，得到各条声波传播路径上传播时间，然后利用传播时间数据通过重建算法，将被测区域内的声速分布重建出来。最后通过式(11)中温度和声速的函数关系计算出对应的被测区域温度场。图4给出了二维测量区域的声波换能器的布置示意图。考虑到测量区域的空间维度、几何形状不同，需要合理的布置声波传感器。已知声波从不同发射源到达各接收点的传播时间，设第k条传播路径上的传播时间为t_k，则有：

式中f(x,y,z)为声速的倒数。

假设被测温度场区域均匀划分成个网格子块，系统共生成M条声波传播路径，则：

其中w_ki表示第k条路径在第i个网格子块的长度，反映了第i个子块对第k条声波路径的贡献。将其表示成线性方程组的形式：

该方程组的矩阵形式为：

Ax＝t (20)

其中，A是M*N的矩阵，M为声波传播路径总数，N为网格总数，x＝(f₁...f_N)^T表示为声速倒数向量，t＝(t₁...t_M)^T表示为声波传播时间向量。

上述即为采用离散区域法得出的温度场重建的数学模型。然后通过对方程(20)的数值求解可获得被测区域各个离散网格内的声波传播速度，再利用(11)中温度与声速之间的关系，得到每个离散网格对应的温度值，最后采用插值等数学方法获得测量区域最终的温度分布。

3.改进型FISTA算法原理

声学温度场重建实质是通过少量投影数据实现病态反问题的求解，未知变量数目常多于投影数据；从数学角度来说，所求方程的解存在不适定性。求解式(14)的经典方法是最小二乘法，通过最小化

来求得x。但若A是病态矩阵，此时条件数较大，通常需引入正则化求解，其中l1正则化方法被广泛地运用解决此问题：

其中||||₁表示为l1范数，表示为x元素绝对值之和。l1范数可以用来诱导最优解的稀疏性。针对式(15)问题，目前已经有多种算法求解，如l1-magic算法，分裂Bregman方法，迭代收缩阈值(ISTA)算法，其中ISTA方法更为简单，已成为一种极为流行的求解方法。ISTA是梯度下降法的延伸，每次迭代利用当前点的信息进行迭代估计，更新迭代点。它的迭代形式为：

x⁽ⁿ⁺¹⁾＝T_α(x⁽ⁿ⁾-2t⁽ⁿ⁾A^T(Ax⁽ⁿ⁾-y)) (22)

其中t⁽ⁿ⁾＞0为步长，一般要求

，T_α为收缩算子。

T_α(x_i)＝(|x_i|-α)₊sgn(x_i) (23)

||为绝对值，()₊表示为|x_i|-α和0中的最大值，sgn为符号函数，表示若q_i大于0，则sgn返回1；等于0，返回0；小于0，则返回-1。在声源识别中，令

为保证x能够收敛到一个极小值。

然而迭代收缩阈值方法(ISTA)的主要缺点是收敛速度比较慢，最近提出一种快速迭代收缩阈值方法(FISTA)，它在保留ISTA算法计算简单的特性的同时，还使得算法的全局收敛速度得到了很大的改善。快速迭代收缩阈值算法初始化x⁽⁰⁾＝0,y⁽¹⁾＝x⁽⁰⁾,t⁽¹⁾＝1。由第n次迭代计算结果x⁽ⁿ⁾到第n+1迭代计算结果x⁽ⁿ⁺¹⁾的具体步骤为：

其中，P₊表示为在非负象限的欧几里得投影。

FISTA相比于ISTA不同之处在于：FISTA下一个迭代点x⁽ⁿ⁺¹⁾不仅仅是依赖于前一次的迭代点x⁽ⁿ⁾，而是依赖于前两次迭代点{x⁽ⁿ⁾,x^(n-1)}的线性组合。FISTA和NNLS的主要区别在于计算步长以及中间变量y的引入。FISTA迭代的步长1/L为固定值，在进行主体循环前即可求出。梯度计算公式为：

FISTA算法与ISTA算法基本相同，但其收敛速度从O(1/T)提高到O(1/T²)。虽然FISTA比ISTA具有更好的线性收敛速度，但是该算法由于差函数及其梯度的计算中存在大矩阵乘积，计算仍然较耗时，并且在运行FISTA时，一个常见的观察是在迭代接近最优解的目标值的轨迹上出现波纹或凸起。为了提高FISTA的收敛速度，以及增加迭代过程的稳定性与鲁棒性。本发明提出了一种改进型FISTA算法，FISTA的线性收敛速率与加速度参数

密切相关，为了抑制迭代过程的纹波不稳定现象，采用参数重置法。即当式(18)在第n+1次迭代值F(xⁿ⁺¹)大于前一次的迭代值F(xⁿ)时，将加速度参数β置0，当第n+1次迭代值F(xⁿ⁺¹)小于等于前一次的迭代值F(xⁿ)时，将加速度参数β置为

，通过重置加速度β，可以很好的解决快速迭代收缩阈值算法过程中存在的波动不稳定现象，使之鲁棒性更强。但是运用FISTA算法进行迭代计算时，差函数及其梯度的计算中存在大矩阵乘积，计算仍然较较耗时，为了解决这一问题，在FISTA反卷积方法的最小化问题的梯度函数中引入正定加权矩阵，可以进一步提高收敛速度和声源识别性能与FISTA算法一样，改进型FISTA算法在每次迭代中仅需要一个矩阵向量乘积运算。

加权矩阵W通常是预先确定的，可以通过尝试法有效地得到。因此，算法在原有的基础上没有额外的计算负担。综上所述，与FISTA算法相比，本发明提出的改进型FISTA算法具有更快的收敛速度，更强的鲁棒性和更好的图像恢复能力。改进型FISTA算法步骤为：

初始化x⁽⁰⁾＝0,y⁽¹⁾＝x⁽⁰⁾，t⁽¹⁾＝1，由第n次迭代计算结果x⁽ⁿ⁾到第n+1迭代计算结果x⁽ⁿ⁺¹⁾的具体步骤为：

4.核极限学习机(KELM)原理

极限学习机(ELM)是一种基于单隐含层前馈神经网络(Single HiddenLayerFeedforward Network，SLFN)的机器学习方法。与传统的SLFN方法不同，ELM方法随机赋值输入层与隐含层之间的权值和隐含层的阈值，以最小化训练误差为目标确定输出权值。ELM方法具有非线性映射能力强、学习速度快、逼近和分类能力强、计算复杂度低及泛化能力好。与那些在训练过程中需要反复优化参数的算法(如SVM、ANN、HMMs、DNN)相比，训练效率高，训练时间大大减少。但是，由于ELM中输入层与隐含层之间的权值和隐含层的阈值是随机初始化的，故ELM的输出通常是不稳定的。因此本发明采用了ELM的一种改进形式KELM。KELM利用核函数将非线性数据映射到高维空间，便于非线性数据的融合。KELM继承了ELM训练速度快的优势。它不容易陷入局部最优，具有较强的非线性拟合能力。KELM的原理如图7所示。

ELM的隐含层输出矩阵H可以改写为：

H＝[h(x₁),...,h(x_N)] (27)

式中h(x_i)为x_i(通过改进型FISTA得到的划分的粗网格中心点温度)的非线性映射，N为训练集的个数。如果该映射的形式未知，则可以构造如下：

其中K是核函数。KELM的输出是：

其中x^*是属于测试集的样本。KELM的输出权值矩阵为：

β^*＝H^T(HH^T+λI^-1)T (30)

其中λ为正则化参数，可以提高β^*计算结果的稳定性。T是训练集的输出向量。

核函数的选择也很重要。由于径向基函数(RBF)具有较强的非线性逼近能力，因此本发明采用了径向基函数。RBF函数有以下形式：

其中δ²为核参数，它决定了高维数据空间的复杂性。x_i,i∈{1,2,...,N}和x^*分别是来自训练集和测试集的样本。KELM的结构如图7所示，

表示x^*的第m维。

在KELM中，δ²和λ是影响检测性能的两个关键参数。本发明采用粒子群优化算法(PSO)，得到了δ²和λ的最优参数。由于参数优化不是本研究的重点，此处省略。

将通过改进型FISTA算法迭代求出的划分的7*7个网格的中心点坐标与温度值作为训练样本，带入核极限学习机进行训练，待核极限学习机训练完成后，将整个被测区域的点的坐标带入提前训练好的极限学习机，求出整个温度场的温度分布。

5.温度场重建算法设计

本发明提出的高精度快速温度场重建算法的重建步骤如表1所示。

表1高精度快速温度场重建算法的重建步骤

高精度快速温度场重建算法示意图如图2所示。

6.温度场重建算法实施方式

结合温度场重建系统与重建算法，温度场重建算法实施方式可简要归纳如下：

1)根据搭建的温度场重建硬件系统待测区域的几何大小、待测温度场的复杂程度，在其每条边上按照特定的方式安装一定数量的超声波换能器并且计算出超声波传播路径数量。

2)将待测区域划分为一定数量的网格，为了权衡算法效率与重建精度，划分的网格数量应该大于或者等于超声波传播路径数量。

3)根据超声波传播路径数量(M)与划分的网格数量(N)，计算每一条超声波传播路径在每一个被划分的网格上的长度值(第i条超声波传播路径在第j个网格上的长度为A_i,j)，形成系数矩阵A。

4)通过改变超声波换能器的收发状态，测量每一条超声波传播路径的飞渡时间(TOF)，形成一个飞渡时间向量y。

5)建立基于目标函数：

上式中矩阵A为距离矩阵；向量x为超声波传播速度的倒数向量；向量y为飞渡时间向量；

6)将改进型FISTA算法应用于求解公式(33)的x：将预先确定的W矩阵，A矩阵以及y矩阵带入公式(33)。将初值x⁽⁰⁾设置为一个N*1的向量带入公式进行迭代求解，迭代次数通过调试经验所确定。最终通过改进型FISTA算法的迭代求解，得到最优的x*，x*即为超声波在每一个划分网格内的传播速度倒数的向量。

7)将经过改进型FISTA算法得到的解x*，经过公式(11)，得到待测区域内每一个划分网格的中心点温度向量：

8)对经过改进型FISTA求解的温度向量与对应坐标作为训练样本，带入构建好的核极限学习机进行模型训练。

9)待训练完后将整个被测区域的坐标带入已经训练好的极限学习机预测，得到整个被测区域的连续温度场分布。

10)通过MATLAB绘图工具，将求出的温度场数据绘图，得到等温线图与三维图。

本发明实施例提供的高精度快速温度场重建算法重建图如图8所示；

其中，图8(a)为二维重建图；图8(b)为三维重建图。

7.与传统算法温度场重建误差比较

为了便于直观理解，本发明利用数值实验途径(通过对4种常用的温度场模型进行仿真)比较典型算法与本发明提出的算法温度分布重建精度。其中，所有的计算过程均在MATLAB软件上运行。为了评价温度分布重建效果，利用平均相对误差(AverageRelativeError，ARE)和均方根误差(RootMeanSquareError，RMSE)作为重建精度的评价指标，其计算公式为：

式中n--测量区域总共划分的网格总数；

T(i)--温度分布模型中网格i的温度值；

--算法重建后网格i的温度值；

T_mean--温度分布模型的平均温度值。

各算法温度场重建误差对照表如表2～表4所示。

表2平均相对误差对照表

表3均方根误差对照表

表4温度场重建时间对照表(单位：s)

在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用全部或部分地以计算机程序产品的形式实现，所述计算机程序产品包括一个或多个计算机指令。在计算机上加载或执行所述计算机程序指令时，全部或部分地产生按照本发明实施例所述的流程或功能。所述计算机可以是通用计算机、专用计算机、计算机网络、或者其他可编程装置。所述计算机指令可以存储在计算机可读存储介质中，或者从一个计算机可读存储介质向另一个计算机可读存储介质传输，例如，所述计算机指令可以从一个网站站点、计算机、服务器或数据中心通过有线(例如同轴电缆、光纤、数字用户线(DSL)或无线(例如红外、无线、微波等)方式向另一个网站站点、计算机、服务器或数据中心进行传输)。所述计算机可读取存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。所述可用介质可以是磁性介质(例如软盘、硬盘、磁带)、光介质(例如DVD)、或者半导体介质(例如固态硬盘Solid StateDisk(SSD))等。

以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，都应涵盖在本发明的保护范围之内。

Claims (8)

1.一种高精度快速温度场重建方法，其特征在于，所述高精度快速温度场重建方法包括：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；

将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布；

所述高精度快速温度场重建方法具体包括以下步骤：

步骤一，根据搭建的温度场重建硬件系统待测区域的几何大小、待测温度场的复杂程度，在每条边上按照特定的方式安装一定数量的超声波换能器并计算出超声波传播路径数量；

步骤二，将待测区域划分为一定数量的网格，且划分的网格数量应该大于或者等于超声波传播路径数量，权衡算法效率与重建精度；

步骤三，根据超声波传播路径数量M与划分的网格数量N，计算每一条超声波传播路径在每一个被划分的网格上的长度值，第i条超声波传播路径在第j个网格上的长度为A_i,j，形成系数矩阵A；

步骤四，通过改变超声波换能器的收发状态，测量每一条超声波传播路径的飞渡时间TOF，形成一个飞渡时间向量y；

步骤五，建立目标函数，并将改进型FISTA算法应用于求解公式(33)的x；

步骤六，将经过改进型FISTA算法得到的解x*，经过超声波速和气体介质温度的关系表达式，得到待测区域内每一个划分网格的中心点温度向量；

步骤七，对经过改进型FISTA求解的温度向量与对应坐标作为训练样本，带入构建好的核极限学习机进行模型训练；

步骤八，待训练完后将整个被测区域的坐标带入已经训练好的极限学习机预测，得到整个被测区域的连续温度场分布；

步骤九，通过MATLAB绘图工具，将求出的温度场数据绘图，得到等温线图与三维图；

所述步骤七中的核极限学习机的构建包括：

ELM的隐含层输出矩阵H改写为：

H＝[h(x₁),...,h(x_N)]；

式中，h(x_i)为x_i的非线性映射，N为训练集的个数，x_i为通过改进型FISTA得到的划分的粗网格中心点温度；若该映射的形式未知，则构造如下：

式中，K是核函数；KELM的输出是：

式中，x^*是属于测试集的样本；KELM的输出权值矩阵为：

β^*＝H^T(HH^T+λI^-1)T；

式中，λ为正则化参数，T是训练集的输出向量；

采用径向基函数，RBF函数有以下形式：

式中，δ²为核参数，决定高维数据空间的复杂性；x_i,i∈{1,2,…,N}和x^*分别是来自训练集和测试集的样本；x^*m表示x^*的第m维；

在KELM中，δ²和λ是影响检测性能的两个关键参数；采用粒子群优化算法PSO，得到δ²和λ的最优参数；

将通过改进型FISTA算法迭代求出的划分的7*7个网格的中心点坐标与温度值作为训练样本，带入核极限学习机进行训练；待核极限学习机训练完成后，将整个被测区域的点的坐标带入提前训练好的极限学习机，求出整个温度场的温度分布。

2.如权利要求1所述高精度快速温度场重建方法，其特征在于，所述步骤三中的系数矩阵A的表达式为：

所述步骤四中的渡时间向量y的表达式为：

所述步骤五中的目标函数为：

式中，矩阵A为距离矩阵；向量x为超声波传播速度的倒数向量；向量y为飞渡时间向量。

3.如权利要求1所述高精度快速温度场重建方法，其特征在于，所述步骤五中的改进型FISTA算法的构建包括：

初始化x⁽⁰⁾＝0,y⁽¹⁾＝x⁽⁰⁾，t⁽¹⁾＝1，由第n次迭代计算结果x⁽ⁿ⁾到第n+1迭代计算结果x^{(n +1)}的步骤为：

(1)x⁽ⁿ⁾＝P₊(y⁽ⁿ⁾-W_tA^T(Ax⁽ⁿ⁾-y)/L)

(2)

(3)if(F(xⁿ)＞F(x^n-1))

β＝0

else

(4)y⁽ⁿ⁺¹⁾＝x⁽ⁿ⁾+β(x⁽ⁿ⁾-x^(n-1))；

所述步骤五中的将改进型FISTA算法应用于求解公式(33)的x包括：

(1)将预先确定的W矩阵，A矩阵以及y矩阵带入公式(33)；

(2)将初值x⁽⁰⁾设置为一个N*1的向量带入公式进行迭代求解，迭代次数通过调试经验所确定；

(3)通过改进型FISTA算法的迭代求解，得到最优的x*，x*即为超声波在每一个划分网格内的传播速度倒数的向量。

4.如权利要求1所述高精度快速温度场重建方法，其特征在于，所述步骤六中的超声波速和气体介质温度的关系表达式为：

式中，c声速，单位m/s；γ气体介质的绝热指数，气体定压热容与定容热容的比值，无量纲量，对空气来说为1.4；R是理想气体常数，8.314^J/(mol·k)；T为气体介质的热力学温度，单位k；M为气体介质的摩尔质量，数值上等于分子量，单位kg/mol，空气为28.8*10^-3kg/mol；

所述待测区域内每一个划分网格中心点温度向量为：

5.一种应用如权利要求1～4任意一项所述高精度快速温度场重建方法的高精度快速温度场重建系统，其特征在于，所述高精度快速温度场重建系统包括：

测量区域划定模块，用于将测量区域划定为一定数量的离散网格；

目标函数构建模块，用于构建基于离散网格点温度误差目标函数；

目标函数求解模块，用于基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；

温度场分布获取模块，用于将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

6.一种计算机设备，其特征在于，所述计算机设备包括存储器和处理器，所述存储器存储有计算机程序，所述计算机程序被所述处理器执行时，使得所述处理器执行如下步骤：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

7.一种计算机可读存储介质，存储有计算机程序，所述计算机程序被处理器执行时，使得所述处理器执行如下步骤：

将测量区域划定为一定数量的离散网格；构建基于离散网格点温度误差目标函数；基于测量的声波飞渡时间数据，采用改进型FISTA算法求解该目标函数，得到网格下的测量区域温度分布；将测量区域网格下得到的温度分布利用核极限学习机进行训练预测，得到整个测量区域的温度场分布。

8.一种信息数据处理终端，其特征在于，所述信息数据处理终端用于实现如权利要求5所述高精度快速温度场重建系统。

Images