深海热液口温度场原位在线声学检测方法.pdf

深海热液口温度场原位在线声学检测方法
    【技术领域】

    本发明属于深海装备技术领域，涉及一种深海热液口温度场原位在线声学检测方法。

    背景技术

    对现代热液活动的温度场分布原位监测是热液活动研究中的重要内容之一，也是研究深海热液喷口及其周围深部生物圈、探究海底热液成矿作用机理的基础研究工作。温度场分布又是热液系统输入海洋热通量计算的重要参数。但由于热液活动区的高温高压极端环境限制，现有的测定热液口温度的有效手段和技术非常有限。迄今为止已经开发了可以测量热液口点温度的温度计。现有的手段主要用温度传感器阵列进行测量，但传统的接触式温度测量装置具有很大的局限性：首先，接触式测量会干扰热液口的温度；其次，接触式温度测量只能测量热液口某些点的温度，无法得到热液区域二维和三维的温度场分布数据。综上所述，目前深海热液口测温的主要手段无法对热液流体的温度场准确测量。因此海洋科学家无法对脊‑洋反应区内整个温洋脊剖面的热能通量进行准确估计，只能在现有资料基础上进行一个恰当的估计。目前，国内尚未有一种能够实现热液喷口温度场原位监测，给出定量温度信息的方法面市。

    【发明内容】

    本发明的目的在于克服传统的海底热液喷口温度测量的缺点，提供了一种运行稳定可靠、测量准确快速的深海热液口温度场原位在线声学检测方法。

    本发明方法的步骤是：

    步骤(1)、将三十二个水声换能器分为两层设置在支架上，层高为0.5m，每一层有十六个水声换能器，十六个水声换能器组成围合成正方形，正方形的每边上设置四个水声换能器，每条边上水声换能器的间距为0.3m；

    步骤(2)、用扫频范围为18kHz～23kHz，声压级为69dB的扫频信号源发射卡驱动第一层任意一个水声换能器发射声波信号，同层其余15个水声换能器都作为接收水声换能器接收声信号，分别将这15个接收水声换能器接收到的声波信号与发射水声换能器发射的声波信号进行互相关函数分析

     $R\left(n\right)=\underset{k=0}{\overset{k=N}{\mathrm{\Σ}}}{f}_1^{*}\left(k\right)\·f_2(k+n)$

    其中R(n)为互相关函数，N为采样点数，可根据所需的测量精度设定，f₁(k)是发送水声换能器发送的声源信号，f₂(k)是接收水声换能器采集到的离散声波信号，互相关函数取得最大值的延迟时间即为声波的传播时间。

    再通过切换电路把按照逆时针方向确定的紧临上一个水声换能器的下一个水声换能器作为发射声波信号的水声换能器，同层其余十五个水声换能器作为接收声波信号的水声换能器，再次利用互相关函数计算声波的传播时间，以此类推，直到完成第二层中的最后一个水声换能器作为发射声波信号的水声换能器时的互相关函数分析。这样切换电路共切换三十一次形成一个完整的周期；

    步骤(3)、计算声波传播路径的平均温度值t。

    c(D，S，t)＝c(0，S，t)+(16.23+0.253t)D

                +(0.213‑0.1t)D²+[0.016

                +0.0002(S‑35)](S‑35)tD

    式中：

    c(0，S，t)＝1449.05+45.7t‑5.21t²

                +0.23t³‑(1.333‑0.126t

                +0.009t²)(S‑35)

    其中c(D，S，t)为声波在海水介质中的速度，D为被测平面的深度，单位为Km；t＝T/10，T为海水温度，单位为℃，S为含盐量的千分数。在固定的热液口盐度和被测海平面的深度均可视为常数，因此声速变化主要与温度有关；

    步骤(4)、利用数据层析还原算法重建海底热液口温度场。

    建立同层平面中十六个水声换能器之间的九十六条声波传播路径的声学测量模型。将被测正方形平面划分为8×8共六十四个网格区域。假定每块区域上的温度分布相同，网格区域标记为m，Δs为声波通过每一块区域的长度，则可建立如下关系。

     $\mathrm{Tn}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{a}_i{\mathrm{\Δs}}_{\mathrm{ni}}$

    其中，m＝64，n＝96，Tn是第n条声波传播路径上的声波飞渡时间。a_i是第i个网格区域上声速的倒数，Δs_ni是声波在第n条传播路径通过第i个网格区域的长度，a_i与Δs_ni都是常数。a_i是各路径上温度的特征函数，是要计算的值，写成矩阵形式，得到如下的方程组：

    AX＝t

    其中

     $A=\left[\begin{array}{cccc}{\mathrm{\ΔS}}_{11}& {\mathrm{\ΔS}}_{12}& \·\·\·& {\mathrm{\ΔS}}_{1m}\\ {\mathrm{\ΔS}}_{21}& {\mathrm{\ΔS}}_{22}& \·\·\·& {\mathrm{\ΔS}}_{2m}\\ \·& \·& & \·\\ \·& \·& & \·\\ \·& \·& & \·\\ {\mathrm{\ΔS}}_{n1}& {\mathrm{\ΔS}}_{n2}& \·\·\·& {\mathrm{\ΔS}}_{\mathrm{nm}}\end{array}\right]$

     $X=\left[\begin{array}{c}{a}_1\\ a_2\\ \·\\ \·\\ \·\\ a_m\end{array}\right],$$t=\left[\begin{array}{c}{T}_1\\ T_2\\ \·\\ \·\\ \·\\ T_n\end{array}\right]$

    矩阵A是已知的，矩阵t中的元素为步骤(2)中得到的声波传播时间，对于此线性方程组，由总体最小二乘方法可得到变量X，进而求出声波在每一块网格区域上传播速度的倒数，利用步骤(3)中速度与温度的关系式，可重建出海底热液口截面温度场。

    所述的总体最小二乘方法为：

        (A+ΔA)X＝t+Δt

    其中，ΔA表示声波在网格中传播路径长度的噪声矩阵，Δt表示声波在网格中传播路径上飞渡时间测量数据的噪声向量。则(A+ΔA)X＝t+Δt可以表述为：

     $[\left(\begin{array}{cc}A& t\end{array}\right)+\left(\begin{array}{cc}\mathrm{\ΔA}& \mathrm{\Δt}\end{array}\right)]\left(\begin{array}{c}X\\ -1\end{array}\right)=0$

    或者

    (W+ΔW)M＝0

    其中，W＝(At)，ΔW＝(ΔA Δt)，M＝(X‑1)^T，W为增广输入数据矩阵，ΔW是叠加于增广数据矩阵的噪声，M是增广权向量。这样，(A+ΔA)X＝t+Δt总体最小二乘解可以简单表示为求解一个解向量M_TLS，使得min||ΔA Δt||_F或min||ΔW||_F为最小。对增广输入数据矩阵W进行奇异值分解，得到：

     $W=\mathrm{U\Σ}{V}^H=\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\σ}}_j{u}_j{v}_j^H$

    其中，σ_j是第j个奇异值，∑＝diag(σ₁，…，σ_n+1)，则总体最小二乘解为：

     $\left\{\begin{array}{c}{M}_{\mathrm{TLS}}=S^{-\left(p\right)}(i+\mathrm{1,1})/S^{-\left(p\right)}\left(\mathrm{1,1}\right),i=1,\·\·\·,p\\ S^{\left(p\right)}=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}\underset{i=1}{\overset{n+1-p}{\mathrm{\Σ}}}{\mathrm{\σ}}_j^2{v}_j^i{\left(v_j^i\right)}^H\end{array}\right.$

    得出解向量M_TLS后再计算变量X，求出声波在每一块网格区域上传播速度的倒数，利用步骤(3)中速度与温度的关系式，可精确地重建出海底热液口截面温度场。

    本发明采用一种非接触式的声学热液温度场原位测量，可以克服传统接触式测量的缺点，在高温、腐蚀、多悬浮颗粒的恶劣环境下连续实时测量，提高热液口热通量测量的准确度。为科学家研究洋壳散热量、建立大洋环流模型、研究热液对临近海域海洋物理特性影响提供帮助。

    【附图说明】

    图1为同层平面的声学路径和拓扑网络图。

    【具体实施方式】

    以下结合附图对本发明的原理、结构做进一步的说明。

    本发明方法所采用的测温装置由水声换能器、机架、电缆和数据采集系统组成。安装在机架上的水声换能器通过电缆连接到驱动与前置放大腔、驱动与前置放大腔、信号采集腔和电池腔组成的数据采集系统。

    如图1所示，是同层平面的声学路径和拓扑网络图。在平面的每一边放置四个水声换能器，四个位置处于每条边五等分处。这样一个平面中共有十六个个水声换能器。

    本发明方法的具体步骤如下：

    步骤(1)、在支架上共安装两层RHS‑30水声换能器，每层RHS‑30水声换能器都处于同一水平面，每层装有十六个水声换能器，层高0.5m，每层的支架围合成正方形，正方形的每边上设置四个标准水声换能器，每条边上水声换能器的间距为0.3m，按照逆时针的顺序从T1换能器到T16换能器均匀布置。

    步骤(2).用扫频范围为18kHz～23kHz，声压级为69dB的扫频信号源发射卡驱动第T1个水声换能器发射声波信号，同层其余15个水声换能器(T2至T16)都作为接收水声换能器接收声信号，分别将T2至T16接收水声换能器接收到的声波信号与发射水声换能器T1发射的声波信号利用公式(1)进行互相关函数分析

     $R\left(n\right)=\underset{k=0}{\overset{k=N}{\mathrm{\Σ}}}{f}_1^{*}\left(k\right)\·f_2(k+n)---\left(1\right)$

    其中R(n)为互相关函数，N为采样点数，可根据所需的测量精度设定，f₁(k)是发送水声换能器发送的声源信号，f₂(k)是接收水声换能器采集到的离散声波信号，互相关函数取得最大值的延迟时间即为声波的传播时间(TOF)。

    再通过切换电路把第2个水声换能器T2作为发射声波信号的水声换能器，其余同层15个作为接收声波信号的水声换能器，再次利用公式(1)计算声波的传播时间，以此类推，直到完成第二层中的第16个水声换能器作为发射声波信号的水声换能器时的互相关函数分析。这样切换电路共切换31次才是一个完整的周期。每个测量周期的时间约为1分钟。

    步骤(3)、计算声波传播路径的平均温度值t。已知声波在海水中的传播速度与温度、盐度和深度有关，关系式如下：

    c(D，S，t)＝c(0，S，t)+(16.23+0.253t)D

                +(0.213‑0.1t)D²+[0.016             (2)

                +0.0002(S‑35)](S‑35)tD

    其中：

    c(0，S，t)＝1449.05+45.7t‑5.21t²

                +0.23t³‑(1.333‑0.126t

                +0.009t²)(S‑35)

    c(D，S，t)是声波在某种海水介质中的速度，单位为m/s；D为被测平面的深度，单位为Km；t＝T/10，液口盐度和被测海平面的深度均可视为常数，声速变化主要与温度有关。利用式(2)中温度与声速的关系得出这些有效传播路径的平均温度值。

    步骤(4)、利用数据层析还原算法重建海底热液口温度场。建立如图

    (2)所示的同层16个水声换能器，96条声波传播路径的声学测量模型。水声换能器对称分布在正方形控制面上，被测空间划分为8×8共64个网格区域。假定每块区域上的温度分布相同，网格区域标记为m，Δs为声波通过每一块区域的长度，则可建立如下关系。

     $\mathrm{Tn}=\underset{i=1}{\overset{m}{\mathrm{\Σ}}}{a}_i{\mathrm{\Δs}}_{\mathrm{ni}}---\left(3\right)$

    其中，m＝64，n＝96，Tn是第n条声波传播路径上的声波飞渡时间，a_i是第i个网格区域上声速的倒数，Δs_ni是声波在第n条传播路径通过第i个网格区域的长度，它们都是常数，取决于测温系统装置和各水声换能器的安装位置。a_i是各路径上温度的特征函数，是要计算的值，写成矩阵形式，得到如下的方程组：

    AX＝t                                            (4)

    其中，

     $A=\left[\begin{array}{cccc}{\mathrm{\ΔS}}_{11}& {\mathrm{\ΔS}}_{12}& \·\·\·& {\mathrm{\ΔS}}_{1m}\\ {\mathrm{\ΔS}}_{21}& {\mathrm{\ΔS}}_{22}& \·\·\·& {\mathrm{\ΔS}}_{2m}\\ \·& \·& & \·\\ \·& \·& & \·\\ \·& \·& & \·\\ {\mathrm{\ΔS}}_{n1}& {\mathrm{\ΔS}}_{n2}& \·\·\·& {\mathrm{\ΔS}}_{\mathrm{nm}}\end{array}\right]$

     $X=\left[\begin{array}{c}{a}_1\\ a_2\\ \·\\ \·\\ \·\\ a_m\end{array}\right],$$t=\left[\begin{array}{c}{T}_1\\ T_2\\ \·\\ \·\\ \·\\ T_n\end{array}\right]$

    方程(4)中，矩阵A是已知的，矩阵t中的元素为步骤(2)中利用公式(1)得到的声波传播时间，对于此线性方程组，由方法一即最小二乘法求解可得到变量X，求出声波在每一块网格区域上传播速度的倒数，利用式(2)中速度与温度的关系式，可重建出海底热液口截面温度场。

    最小二乘作为一种数据拟合方法，在求解矩阵方程AX＝t时，仅假设方程右端n条路径飞渡时间矩阵t中的测量数据存在误差。而在实际测量系统中，声波在网格中的传播路径长度矩阵A和飞渡时间矩阵t中同时存在误差，因此最小二乘方法具有局限性。更合理的方法是同时考虑A、t二者的扰动。本实施例中优先采用Golub G H，Loan C F V等人提出的总体最小二乘方法求解矩阵反问题，使得来自A、t的噪声扰动影响最小。总体最小二乘的求解方法如下，考虑矩阵方程：

    (A+ΔA)X＝t+Δt                                     (5)

    其中，ΔA表示声波在网格中传播路径长度的噪声矩阵，Δt表示声波在网格中传播路径上飞渡时间测量数据的噪声向量。式(5)可以表述为：

     $[\left(\begin{array}{cc}A& t\end{array}\right)+\left(\begin{array}{cc}\mathrm{\ΔA}& \mathrm{\Δt}\end{array}\right)]\left(\begin{array}{c}X\\ -1\end{array}\right)=0---\left(6\right)$

    或者

    (W+ΔW)M＝0                                      (7)

    其中，W＝(A t)，ΔW＝(ΔA Δt)，M＝(X‑1)^T，W为增广输入数据矩阵，ΔW是叠加于增广数据矩阵的噪声，M是增广权向量。这样，矩阵方程(5)的总体最小二乘解可以简单表示为求解一个解向量M_TLS，使得min||ΔA Δt||_F或min||ΔW||_F为最小。对增广输入数据矩阵W进行奇异值分解，得到：

     $W=\mathrm{U\Σ}{V}^H=\underset{j=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\σ}}_j{u}_j{v}_j^H---\left(8\right)$

    其中，σ_j是第j个奇异值，∑＝diag(σ₁，…，σ_n+1)，则总体最小二乘解为：

     $\left\{\begin{array}{c}{M}_{\mathrm{TLS}}=S^{-\left(p\right)}(i+\mathrm{1,1})/S^{-\left(p\right)}\left(\mathrm{1,1}\right),i=1,\·\·\·,p\\ S^{\left(p\right)}=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}\underset{i=1}{\overset{n+1-p}{\mathrm{\Σ}}}{\mathrm{\σ}}_j^2{v}_j^i{\left(v_j^i\right)}^H\end{array}\right.---\left(9\right)$

    得出解向量M_TLS后再计算变量X，求出声波在每一块网格区域上传播速度的倒数，利用式(2)中速度与温度的关系式，可精确地重建出海底热液口截面温度场。