CN101604019B - 一种海洋环境与声场不确实性的表征和传递的快速计算方法 - Google Patents

Abstract

本发明涉及声学领域，具体涉及一种海洋环境与声场不确实性的表征和传递的快速计算方法，解决的问题是不确实环境参量和不确实声压场的随机谱表征，以及不确实性如何从环境参量传递到声压场。本发明对不确实环境参量和不确实声压场分别建立以已知量为中心的领域空间，并对不确实量作概率密度函数描述，根据概率密度函数作不确实量的多项式混沌展开随机谱表征，进一步将多项式混沌随机谱表征与确定波方程相结合，导出嵌入不确实性的随机波方程，获得不确实性系数的偏微分方程组，最后得到由不确实性系数加权的多项式混沌基函数线性叠加构成的不确实声压场。本发明有益的效果：得到的不确实声压场与MonteCarlo得到的结果完全一致，其计算速度提高了10倍以上。

Description

一种海洋环境与声场不确实性的表征和传递的快速计算方法

所属技术领域

本发明涉及声学领域，具体涉及一种海洋环境与声场不确实性的表征和传递的快速计算方法。

背景技术

由于海洋环境参量例如海深、声速剖面、海底地声特性等不确实，即使发射确定信号，海洋中接收点声压场也表现为随机起伏场：在空间和时间上都具有不确实性和变异性。而先进的声纳信号处理都基于传播模型对海洋中传播的信号作预测，再与接收的信号作拷贝相关处理。由于处理中使用确定的传播模型，即输入一组特定的环境参量，得到确定的预测声压场。这样的传播模型不能反映空时变异和不确实的实际声场，往往导致声纳信号处理的性能不稳定。因而，寻求海洋环境和海洋声场的不确实性表征和传递规律是研究减弱、消除和进一步利用不确实性对声纳处理性能影响的基础。

数学上，不确实性定义为误差的概率密度函数(probability density function PDF)。误差是指实际场与估计/测量场的差。在基于模型的预测问题中，初始条件误差、模型及边界条件的误差都会影响预测的精度，引起预测声场的不确实性，是初始误差和模型集成过程引入的误差的综合影响。在海洋环境中，对于近、中距离，通常我们认为声传播模型的精度已足够高，对预测声场不确实性的贡献可以忽略，因而认为引起预测声场不确实性的主要因素是作为初始条件和边界条件的海洋环境知识的不确实性。环境知识不确实性主要来自于水体，内波、潮、锋、涡、流及其动力特性引起的声速梯度变化。边界不确实性主要来自于界面，表面波，海底粗糙度、地形变化、海底构成等。本发明中，不确实性者，指不完备的环境知识及所导致的声场的随机性，它们的定量度量用一种正交随机多项式展开(PCE)所形成的随机谱作表征。

本发明主要针对假定声场预测模型是精确的，由于海洋环境知识不完备所引起的不确实性导致声场不确实性的问题。不确实声场建模将海洋环境及在其中传播的声场两部分合在一起考虑，将不确实性作为一个必要的组成嵌入，并提供它对声场影响的概率度量。波方程是单变量(声压场)和单参量(声速场)方程，环境部分不确实，通过声速场耦合到声场部分，建立一组耦合差分方程，既描述声场的传播，又描述所关联的不确实性的传播。

发明内容

本发明所要解决的技术问题之一是环境参量不确实性和声场不确实性的表征方法。用于驱动声传播模型的环境参量包括水体声速梯度剖面，海深，沉积层声速梯度、沉积层密度、沉积层衰减、沉积层厚度，基底声速、基底密度、基底衰减，用一个向量表示。不确实环境参量由一个确知的向量及其邻域构成，邻域的大小决定不确实性程度。同样一个基阵采集的随机声场也由一个确知的声场向量与邻域构成。不确实性的领域由多项式混沌展开作随机谱表征。

本发明所要解决的技术问题之二是环境参量的不确实性如何传递到声场不确实性。声传播模型，例如射线模型、抛物方程模型、简正波模型等，能够解决由确定的环境参量如何产生确定的声场问题，而不能直接解决不确实性由环境参量到声场的传递，因而也无法直接给出用于不确实信号检测所需的不确实驾驶向量空间。将不确实性嵌入到传播模型产生随机的传播模型既描述声场的传播也实现不确实性从环境参量到声场的传递。

为解决以上技术问题，本发明是提出以下技术方案实现的：

本发明所述的这种海洋环境与声场不确实性的表征与传递的实现方法，步骤如下：

1.海洋环境不确实性和声场不确实性的表征：

(1)海洋环境不确实性表征

用向量b表征环境参量(由水中声速梯度，海深，沉积层声速梯度，沉积层密度，沉积层衰减，沉积层厚度，基底声速，基底密度，基底衰减构成)。用b₀表示我们利用测量或反演方法获得的已知环境参量，作为不确实环境参量的均值，不确实环境参量采用参量不确实性域

来表征，它表示了环境参量落在b₀附近的一个领域之中，用概率密度函数p(b)表征不确实环境参量的分布。

表征了环境参量不确实性的范围。下面对不确实环境参量作多项式混沌展开随机谱表征。

步骤一：将不确实性环境参量空间

表示成随机变量b；

步骤二：将不确实环境参量b分解为确实性部分和随机扰动部分：

b(r，z)＝b₀+δb(r，z)                         (1)

其中b₀为确定性部分，δb为随机扰动部分，r，z分别表示距离和深度。

步骤三：根据环境参量b的概率概率密度分布p(b)，对随机扰动部分作多项式混沌展开：

$\mathrm{\δb}(r,z;\mathrm{\θ})=\underset{s=0}{\overset{\∞}{\mathrm{\Σ}}}{\mathrm{\α}}_s(r,z){\mathrm{\Λ}}_s\left(\mathrm{\ξ}\left(\mathrm{\θ}\right)\right),---\left(2\right)$

式中Λ_s(ξ(θ))是对应环境参量b概率密度分布的随机多项式混沌基函数，θ为随机参量空间。例如，b是高斯分布，Λ_s(ξ(θ))为Gauss-Hermite多项式；b是均匀分布，Λ_s(ξ(θ))为Legendre多项式；b是γ分布，Λ_s(ξ(θ))为Laguerre多项式，b是β分布，Λ_s(ξ(θ))为Jacobi多项式，b是泊松分布，Λ_s(ξ(θ))为Charlier多项式，等等。α_s(r，z)是随机扰动环境参量投影到随机多项式混沌基函数得到的确定性系数。

步骤四：对无穷多阶多项式混沌展开作有限截取

对不同分布的环境参量利用对应的多项式混沌基函数展开，展开的阶数是有限的。对于高斯分布分布的b，只要一阶就够了，其中α₀(r，z)等于0，α₁(r，z)＝σ，式中σ是随机扰动环境参量的标准离差。具有类似的特性，其它分布的环境参量也可以利用投影方法得到有限阶(阶数为S)的确定性系数。

步骤五：多项式混沌项数确定

随机环境参量b由多个元素构成，b＝[b₁b₂...b_n]，根据最小均方误差准则，n维随机变量b所必需的多项式项数为(nS)！/n！S！。

下面的过程我们针对n＝1时的情况进行讨论。其中的环境参量选取水层的声速梯度，令

b(r，z，θ)＝c(r，z，θ)＝c₀(z)+δc(r，z，θ)           (3)

(2)声场不确实性表征

①确实的环境参量b₀驱动声传播模型产生确实的声场p₀；

②不确实环境参量空间

对应不确实声场空间

③用概率密度函数p(p)表征不确实声场的分布；

④表征了声场不确实性的范围。

同样，对不确实声场也作多项式混沌展开随机谱表征。

步骤一：不确实声场作多项式混沌展开：

$p(r,z;\mathrm{\θ})=\underset{q=0}{\overset{\∞}{\mathrm{\Σ}}}{\mathrm{\γ}}_q(r,z){\mathrm{\Λ}}_q\left(\mathrm{\ξ}\left(\mathrm{\θ}\right)\right)---\left(4\right)$

步骤二：对无穷多阶多项式混沌展开作有限截取

根据最小均方误差准则，对无穷多阶多项式混沌展开作有限阶Q的截取：

$p(r,z;\mathrm{\θ})\widetilde{\‾}\underset{q=0}{\overset{Q}{\mathrm{\Σ}}}{\mathrm{\γ}}_q(r,z){\mathrm{\Λ}}_q\left(\mathrm{\ξ}\left(\mathrm{\θ}\right)\right)---\left(5\right)$

式中p(r，z；θ)为预测的随机声压场。

2.环境参量不确实性到声场不确实性的传递方法：

将环境参量与预测声场的多项式随机谱表征与确实性波方程相结合，得到嵌入不确实性的随机波方程，实现环境参量不确实性到声场不确实性的传递。具体步骤如下：

(1)确定性波方程到随机波方程

步骤一：基于简正波传播模型的确定性波方程

${\▿}^{2}p+{\mathrm{\κ}}^{2}p=0---\left(6\right)$

式中

表示拉普拉氏算子，κ是波数，κ＝ω/c(z)，ω为角频率。

步骤二：随机扰动声速作二项式展开

对随机扰动声速作二项式展开，并保留一阶项，则不确实声速的波数可近似为，

${\widetilde{\mathrm{\κ}}}^2=\frac{{\mathrm{\ω}}^2}{c_0^2{(1+\mathrm{\δc})}^2}\≈{\mathrm{\κ}}^2(1-\frac{2\mathrm{\δc}(r,z;\mathrm{\θ}}{c_0\left(z\right)})---\left(7\right)$

步骤三：确定性波方程推广到随机波方程

将近似的波数代入到确定性波方程，将确定性波方程推广到随机波方程：

$\frac{{\∂}^{2}p(r,z;\mathrm{\θ})}{{\∂r}^2}+\frac{{\∂}^{2}p(r,z;\mathrm{\θ})}{{\∂z}^2}+{\mathrm{\κ}}^{2}p-2{\mathrm{\κ}}^2\frac{\mathrm{\δc}(r,z;\mathrm{\θ})}{c_0(r,z)}p=0---\left(8\right)$

将式(3)和式(4)代入(8)式，并对等式两边同时乘以Λ₁(ξ(θ))，然后通过系总平均将各项投影到基空间，可得

$\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_q(r,z)}{{\&PartialD;r}^2}<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>+\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_q(r,z)}{{\&PartialD;z}^2}<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>+{\mathrm{\&kappa;}}^2\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>---\left(9\right)$

$-\frac{2{\mathrm{\&kappa;}}^2}{c_0\left(z\right)<{\mathrm{\&Lambda;}}_1^2>}\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\underset{s=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s(r,z){\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_s{\mathrm{\&Lambda;}}_1>=0$

将混沌基函数的正交性应用于上述方程，得到一组声场不确实性系数的耦合差分方程，

$\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;\mathrm{\&gamma;}}_1^2(r,z)}{{\&PartialD;r}^2}+\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;z}^2}+{\mathrm{\&kappa;}}^2\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_1(r,z)-\frac{{2\mathrm{\&kappa;}}^2}{c_0\left(z\right)<{\mathrm{\&Lambda;}}_1^2>}\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\underset{s=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s(r,z){\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_s{\mathrm{\&Lambda;}}_1>=0---\left(10\right)$

以高斯分布随机声速为例进行下面公式的推导。δc(r，z，θ)＝α₀He₀+α₁He₁+…＝σξ，ξ是均值为0，方差为1的高斯随机变量。将δc代入(10)式给出第1个不确实系数的耦合方程：

$\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;z}^2}+{\mathrm{\&kappa;}}^2{\mathrm{\&gamma;}}_1(r,z)-\frac{2{\mathrm{\&kappa;}}^2\mathrm{\&sigma;}}{c_0<{\mathrm{He}}_1^2>}\underset{q=0}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z)<{\mathrm{He}}_q\mathrm{\&xi;}{\mathrm{He}}_1>=0---\left(11\right)$

利用Hermite基函数的回归特性和正交性，计算耦合项的系总平均，得到系数γ₁的偏微分回归方程：

$\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{\&PartialD;r^2}+\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{\&PartialD;z^2}+D{\mathrm{\&gamma;}}_1(r,z)-E\sqrt{1}{\mathrm{\&gamma;}}_{l-1}(r,z)-E\sqrt{l+1}{\mathrm{\&gamma;}}_{l+1}(r,z)=0---\left(12\right)$

式中D≡κ²，E≡2κ²σ/c₀。用向量-矩阵形式，上述偏微分方程可化简为

$\frac{{\&PartialD;}^2\mathrm{\&gamma;}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2\mathrm{\&gamma;}}{{\&PartialD;z}^2}+\mathrm{A\&gamma;}=0---\left(13\right)$

式中，A为环境参量不确实性到声场不确实性的耦合方程。A是对称的三角阵。

(2)求解随机波方程得到不确实声场

步骤一：耦合矩阵A特征值分解

$A=\underset{q=1}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_q{g}_q---\left(14\right)$

式中G＝[{g₁}，{g₂}，…]为耦合矩阵A的特征向量，GG^-1＝I，(13)式可以表示为

$\frac{{\&PartialD;}^2{G}^{-1}\mathrm{\&gamma;}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{G}^{-1}\mathrm{\&gamma;}}{{\&PartialD;z}^2}+G^{-1}{\mathrm{AGG}}^{-1}\mathrm{\&gamma;}=0---\left(15\right)$

令 $\widehat{\mathrm{\&gamma;}}=G^{-1}\mathrm{\&gamma;}$ ，上式可以简化为

$\frac{{\&PartialD;}^2\widehat{\mathrm{\&gamma;}}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2\widehat{\mathrm{\&gamma;}}}{{\&PartialD;z}^2}+\mathrm{\&Omega;}\widehat{\mathrm{\&gamma;}}=0---\left(16\right)$

式中，Ω＝G^-1AG为对角阵，对角阵元素为特征值λ₁，λ₂，…。利用特征值和特征向量分解技术，得到第1个确定系数的波方程

$\frac{{\&PartialD;}^2{\widehat{\mathrm{\&gamma;}}}_1}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{\widehat{\mathrm{\&gamma;}}}_1}{{\&PartialD;z}^2}+{\mathrm{\&lambda;}}_1{\widehat{\mathrm{\&gamma;}}}_1=0---\left(17\right)$

步骤二：确定系数求解

等式(17)与标准波方程(6)有相同的形式。因而求解过程与标准的波方程求解过程是一样的。

的解可以表示为

${\widehat{\mathrm{\&gamma;}}}_u(r,z)=\frac{i}{\sqrt{8\mathrm{\&pi;r}}}{e}^{-i\frac{\mathrm{\&pi;}}{4}}{a}_{u0}\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&psi;}}_m\left(z_s\right){\mathrm{\&psi;}}_m\left(z\right)}{\sqrt{k_{\mathrm{rm}}}}{e}^{i{k}_{\mathrm{rm}}r}---\left(18\right)$

式中k_rm、ψ_m(z)分别是水平波数和模深度函数，由简正波传播模型求解得到。a_u0表示矩阵A的特征向量的逆的第一列中的第u个元素。在声速梯度不确实情况下，模深度函数没有受到不确实声速的影响，但声速梯度不确实导致了水平波数的不确实， $k_{\mathrm{rm}}=\sqrt{{\mathrm{\&lambda;}}_1-k_{\mathrm{zm}}^2}$ ，其中k_zm是深度波数。根据γ₁与

的关系式，系数γ₁是

的线性组合

${\mathrm{\&gamma;}}_1=\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{\mathrm{lu}}{\widehat{\mathrm{\&gamma;}}}_u(r,z)---\left(19\right)$

步骤三：Q个多项式混沌基函数的确定性系数加权累加得到不确实性声场

最后，在远场近似下的不确实性声场可表示为：

$p(r,z;\mathrm{\&theta;})=\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\left[\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{\mathrm{qu}}{\widehat{\mathrm{\&gamma;}}}_u(r,z)\right]{\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)$

$=\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\left[\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{\mathrm{qu}}\left(\frac{i}{\sqrt{8\mathrm{\&pi;r}}}{e}^{-i\frac{\mathrm{\&pi;}}{4}}{a}_{u0}\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&psi;}}_m\left(z_s\right){\mathrm{\&psi;}}_m\left(z\right)}{\sqrt[4]{{\mathrm{\&lambda;}}_u-k_{\mathrm{zm}}^2}}{e}^{i\sqrt{{\mathrm{\&lambda;}}_u-k_{\mathrm{zm}}^2}r}\right)\right]{\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)---\left(20\right)$

若SSP是确实的，则上式中只有γ₀是不等于0的，收敛到声速确实时的声场简正模表达式。

步骤四：不确实性声场的均值和相关函数

由多项式混沌展开的不确实声场是随机声场，其均值可表示为

$<p(r,z;\mathrm{\&theta;})>={\widehat{\mathrm{\&gamma;}}}_0(r,z)---\left(21\right)$

相关函数为

$R_{\mathrm{pp}}(r_1,z_1;r_2,z_2)=\underset{q=1}{\overset{Q}{\mathrm{\&Sigma;}}}{\widehat{\mathrm{\&gamma;}}}_q*(r_1,z_1){\widehat{\mathrm{\&gamma;}}}_q(r_2,z_2)---\left(22\right)$

本发明能带来以下有益效果：

(1)对海洋环境不确实性和声场不确实性进行了定量化表征，分别表示为海洋参量不确实邻域空间

和声场不确实性邻域空间

，由概率密度函数p(b)和p(p)定量化描述，并将不确实性部分表征为多项式混沌随机谱，为减弱、消除和利用不确实性的信号处理所需的驾驶向量邻域空间

提供输入。

(2)不确实性海洋环境参量和不确实性的声场的多项式混沌随机谱表示与简正波声传播模型相结合，得到嵌入不确实性的随机波方程，提供了不确实性从海洋环境参量到声场的耦合和传递。通过这种传递方法比通常采用的Monte Carlo方法可以提高计算效率，节省计算时间。利用Monte Carlo方法获得(21)和(22)的结果，需要进行至少3000次传播模型的运行，而利用多项式混沌展开方法，只需要运行Q次传播模型，而Q＜＜3000。

附图说明

图1海洋环境与声场不确实性的表征和传递结构框图

图2(a)确实声场与高斯分布不确实声速N(0，2.5²)导致的100Hz不确实声场的MC结果与PCE结果比较：确实声场

图2(b)确实声场与高斯分布不确实声速N(0，2.5²)导致的100Hz不确实声场的MC结果与PCE结果比较：不确实声场MC仿真

图2(c)确实声场与高斯分布不确实声速N(0，2.5²)导致的100Hz不确实声场的MC结果与PCE结果比较：不确实声场PCE仿真

图3(a)高斯分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：均值

图3(b)高斯分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：方差

图3(c)高斯分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：归一化自相关

图4(a)确实声场与高斯分布不确实声速N(0，2.5²)导致的200Hz不确实声场的MC结果与PCE结果比较：确实声场

图4(b)确实声场与高斯分布不确实声速N(0，2.5²)导致的200Hz不确实声场的MC结果与PCE结果比较：不确实声场MC仿真

图4(c)确实声场与高斯分布不确实声速N(0，2.5²)导致的200Hz不确实声场的MC结果与PCE结果比较：不确实声场PCE仿真

图5(a)高斯分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：均值

图5(b)高斯分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：方差

图5(c)高斯分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：归一化自相关

图6(a)确实声场与高斯分布不确实声速N(0，5²)导致的100Hz不确实声场的MC结果与PCE结果比较：确实声场

图6(b)确实声场与高斯分布不确实声速N(0，5²)导致的100Hz不确实声场的MC结果与PCE结果比较：不确实声场MC仿真

图6(c)确实声场与高斯分布不确实声速N(0，5²)导致的100Hz不确实声场的MC结果与PCE结果比较：不确实声场PCE仿真

图7(a)高斯分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：均值

图7(b)高斯分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：方差

图7(c)高斯分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：归一化自相关

图8(a)均匀分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：均值

图8(b)均匀分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：方差

图8(c)均匀分布不确实声速c～N(0，2.5²)100Hz声场的均值、方差和归一化自相关：归一化自相关

图9(a)均匀分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：均值

图9(b)均匀分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：方差

图9(c)均匀分布不确实声速c～N(0，2.5²)200Hz声场的均值、方差和归一化自相关：归一化自相关

图10(a)均匀分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：均值

图10(b)均匀分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：方差

图10(c)均匀分布不确实声速c～N(0，5²)100Hz声场的均值、方差和归一化自相关：归一化自相关

具体实施方式

下面结合具体实施例子和附图对本发明做进一步的描述：

图1是本发明的海洋环境与声场不确实性的表征和传递结构框图。从图中可以看出本发明的实施过程：对接收阵接收的阵数据，首先进行放大、滤波，然后进行A/D转换，将接收的模拟信号进行数字采样。对采样后的离散随机声压场与预测随机声压场作拷贝相关处理，处理后的结果送给显示，给出距离和深度两维的能量分布。最后，根据处理结果进行有无目标的判决和目标距离/深度的估计。其中预测随机声压场由本发明提供：不确实海洋环境参量经多项式混沌展开后与确定性波方程相结合，得到将不确实性由海洋环境参量耦合到声场的随机波方程，对随机波方程的耦合矩阵作特征值分解，并作适当变换，得到一组确定性系数方程，确定性的解与多项式混沌线性组合得到不确实声场，不确实声场作为不确实驾驶向量的输入。

实施例子采用的浅海环境参量为Pekeris波导，深度200米，垂直布阵40元，阵元间距5米，声源位于距离5km，深度50m。共实现3000次Monte Carlo仿真。

实施例子一：不确实声速剖面服从高斯分布

当不确实声速剖面服从高斯分布时，式(13)中的矩阵A为对称三角阵，其中对角线给出确定性部分的信息，而非对角线给出了不确实环境参量与不确实声场相耦合的信息。例如取不确实声压场的五阶多项式混沌展开，则A为6×6的矩阵，

$A=\left[\begin{array}{cccccc}D& -E& 0& 0& 0& 0\\ -E& D& -\sqrt{2}E& 0& 0& 0\\ 0& -\sqrt{2}E& D& -\sqrt{3}E& 0& 0\\ 0& 0& -\sqrt{3}E& D& -2E& 0\\ 0& 0& 0& -2E& D& -\sqrt{5}E\\ 0& 0& 0& 0& -\sqrt{5}E& D\end{array}\right]---\left(23\right)$

图2给出了f＝100Hz时确实声场与不确实性声速服从高斯分布的垂直阵采集的随机声场的实部与虚部，与Monte Carlo仿真结果的比较。(a)声速为1500m/s时垂直阵采集的声场(实部与虚部)，(b)不确实声速服从高斯分布N(0，2.5²)时Monte Carlo仿真结果，(c)不确实声速服从高斯分布N(0，2.5²)时100阶PCE展开的仿真结果。

图3给出了图2中接收水听器在50m深度处不确实声场的均值、方差和归一化自相关。

图4给出了f＝200Hz时确实声场与不确实性声速服从高斯分布N(0，2.5²)时的垂直阵采集的随机声场的实部与虚部，与Monte Carlo仿真结果的比较，Q＝300。

图5给出了图4中接收水听器在50m深度处不确实声场的均值、方差和归一化自相关。

图6给出了f＝100Hz时确实声场与不确实性声速服从高斯分布N(0，5²)时垂直阵采集的随机声场的实部与虚部，与Monte Carlo仿真结果的比较，Q＝200。

图7给出了图2中接收水听器在50m深度处不确实声场的均值、方差和归一化自相关。

表1给出了f＝100Hz当不确实声速服从不同方差高斯分布时所需的最小多项式混沌阶数。

实施例子二：不确实声速剖面服从均匀分布

式(13)给出的随机波方程不仅适合高斯分布，也适合其它分布，如均匀分布，所不同的只是A的取值。不确实声速均匀分布时，取不确实声压场的五阶多项式混沌展开，A仍为6×6的矩阵，

$A=\left[\begin{array}{cccccc}D& -\sqrt{1/3}E& 0& 0& 0& 0\\ -\sqrt{1/3}E& D& -\sqrt{4/5}E& 0& 0& 0\\ 0& -\sqrt{4/5}E& D& -\sqrt{9/35}E& 0& 0\\ 0& 0& -\sqrt{9/35E}& D& -\sqrt{16/63}E& 0\\ 0& 0& 0& -\sqrt{16/63}E& D& -\sqrt{25/99}E\\ 0& 0& 0& 0& -\sqrt{25/99}E& D\end{array}\right]---\left(24\right)$

图8给出了f＝100Hz、不确实性声速服从均匀分布U(0，2.5²)、接收水听器在50m深度处时不确实声场的均值、方差和归一化自相关，Q＝25。

图9给出了f＝200Hz、不确实性声速服从均匀分布U(0，2.5²)、接收水听器在50m深度处时不确实声场的均值、方差和归一化自相关，Q＝50。

图10给出了f＝100Hz、不确实性声速服从均匀分布U(0，5²)、接收水听器在50m深度处时不确实声场的均值、方差和归一化自相关，Q＝50。

表2给出了f＝100Hz不确实声场服从不同方差均匀分布时所需的最小多项式混沌阶数。

表3给出了当均匀分布不确实声场标准离差为5m/s时，不同频率所需的最小多项式混沌展开阶数。

从以上处理结果可以看出，无论从不确实声场的实部和虚部，或是不确实声场的均值、方差和归一化自相关，多项式混沌展开的不确实声场与Monte Carlo仿真结果都是一致的。而多项式混沌展开方法所需的计算量要较Monte Carlo方法小得多，在上述示例中，不确实声速取值大到标准离差5m/s时，在频率f＝200Hz时，所需的多项式混沌阶数为300，计算量是3000次Monte Carlo方法的1/10。如果环境参量由1增大到n，Monte Carlo方法的计算量是3000×2ⁿ，而多项式混沌方法则是(n+Q)！/n！Q！，计算的有效性更明显。同时本发明提供的方法还提供了不确实声速是如何传递到不确实声场的内在机理。

表1不确实声速服从高斯分布时不确实声场所需的最小多项式混沌阶数

标准离差σ(m/s)	1	2	3	4	5	6	7	8	9	10
											最小多项式阶数	15	35	67	110	160	220	290	380	470	570

表2不确实声速服从均匀分布时不确实声场所需的最小多项式混沌阶数

标准离差σ(m/s)	1	2	3	4	5	6	7	8	9	10
											最小多项式阶数	12	20	30	40	50	60	70	80	90	100

表3不确实声速服从均匀分布时不同频率声场所需的最小多项式混沌阶数

频率(Hz)	50	100	150	200	250	300
							最小多项式阶数	25	50	75	100	125	150

以上对本发明的描述不具有限制性，如果本领域的普通技术人员受其启示，在不脱离本发明权利要求的保护的情况，作出本发明的其它结构变形和实施方式，均属于本发明的保护范围。

Claims (3)

1.一种海洋环境与声场不确实性的表征和传递的快速计算方法，其特征在于：

一、海洋环境不确实性和声场不确实性的表征：

(1)海洋环境不确实性表征

用向量b表征环境参量，用b₀表示我们利用测量或反演方法获得的已知环境参量，作为不确实环境参量的均值，不确实环境参量采用参量不确实性域

来表征，它表示了环境参量落在b₀附近的一个领域之中，用概率密度函数p(b)表征不确实环境参量的分布，

表征了环境参量不确实性的范围，下面对不确实环境参量作多项式混沌展开随机谱表征；

步骤一：将不确实性环境参量空间

表示成随机变量b；

步骤二：将不确实环境参量b分解为确实性部分和随机扰动部分：

b(r，z)＝b₀+δb(r，z)    (1)

其中b₀为确定性部分，δb为随机扰动部分，r，z分别表示距离和深度；

步骤三：根据环境参量b的概率概率密度分布p(b)，对随机扰动部分作多项式混沌展开：

$\mathrm{\&delta;b}(r,z;\mathrm{\&theta;})=\underset{s=0}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s(r,z){\mathrm{\&Lambda;}}_s\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right),---\left(2\right)$

式中Λ_s(ξ(θ))是对应环境参量b概率密度分布的随机多项式混沌基函数，θ为随机参量空间α_s(r，z)是随机扰动环境参量投影到随机多项式混沌基函数得到的确定性系数；

步骤四：对无穷多阶多项式混沌展开作有限截取

对不同分布的环境参量利用对应的多项式混沌基函数展开，展开的阶数是有限的；

步骤五：多项式混沌项数确定

随机环境参量b由多个元素构成，b＝[b₁b₂...b_n]，根据最小均方误差准则，n维随机变量b所必需的多项式项数为(n+S)！/n！S！；其中S为多项式混沌的阶数；

环境参量b选取为水层的声速梯度，则

b(r，z，θ)＝c(r，z，θ)＝c₀(z)+δc(r，z，θ)    (3)

(2)声场不确实性表征

①确实的环境参量b₀驱动声传播模型产生确实的声场p₀；

②不确实环境参量空间

对应不确实声场空间

③用概率密度函数p(p)表征不确实声场的分布；

④

表征了声场不确实性的范围；

同样，对不确实声场也作多项式混沌展开随机谱表征；

步骤一：不确实声场作多项式混沌展开：

$p(r,z;\mathrm{\&theta;})=\underset{q=0}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z){\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)---\left(4\right)$

步骤二：对无穷多阶多项式混沌展开作有限截取

根据最小均方误差准则，对无穷多阶多项式混沌展开作有限阶Q的截取：

$p(r,z;\mathrm{\&theta;})\&ap;\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z){\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)---\left(5\right)$

式中p(r，z；θ)为预测的随机声压场；

二、环境参量不确实性到声场不确实性的传递方法：

将环境参量与预测声场的多项式随机谱表征与确实性波方程相结合，得到嵌入不确实性的随机波方程，实现环境参量不确实性到声场不确实性的传递，具体步骤如下：

(1)确定性波方程到随机波方程

步骤一：基于简正波传播模型的确定性波方程

${\&dtri;}^{2}p+{\mathrm{\&kappa;}}^{2}p=0---\left(6\right)$

式中表示拉普拉氏算子，κ是波数，κ＝ω/c(z)，ω为角频率；

步骤二：随机扰动声速作二项式展开

对随机扰动声速作二项式展开，并保留一阶项，则不确实声速的波数可近似为，

${\widetilde{\mathrm{\&kappa;}}}^2=\frac{{\mathrm{\&omega;}}^2}{c_0^2{(1+\mathrm{\&delta;c})}^2}\&ap;{\mathrm{\&kappa;}}^2(1-\frac{2\mathrm{\&delta;c}(r,z;\mathrm{\&theta;})}{c_0\left(z\right)})---\left(7\right)$

步骤三：确定性波方程推广到随机波方程

将近似的波数代入到确定性波方程，将确定性波方程推广到随机波方程：

$\frac{{\&PartialD;}^{2}p(r,z;\mathrm{\&theta;})}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^{2}p(r,z;\mathrm{\&theta;})}{{\&PartialD;z}^2}+{\mathrm{\&kappa;}}^{2}p-{2\mathrm{\&kappa;}}^2\frac{\mathrm{\&delta;c}(r,z;\mathrm{\&theta;})}{c_0(r,z)}p=0---\left(8\right)$

将式(3)和式(4)代入(8)式，并对等式两边同时乘以Λ₁(ξ(θ))，然后通过系总平均将各项投影到基空间，

可得

$\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_q(r,z)}{{\&PartialD;r}^2}<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>+\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_q(r,z)}{{\&PartialD;z}^2}<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>+{\mathrm{\&kappa;}}^2\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_1>$

(9)

$-\frac{2{\mathrm{\&kappa;}}^2}{c_0\left(z\right)<{\mathrm{\&Lambda;}}_1^2>}\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\underset{s=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s(r,z){\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_s{\mathrm{\&Lambda;}}_1>=0$

将混沌基函数的正交性应用于上述方程，得到一组声场不确实性系数的耦合差分方程，

$\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{\&PartialD;{\mathrm{\&gamma;}}_1^2(r,z)}{{\&PartialD;r}^2}+\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;z}^2}+{\mathrm{\&kappa;}}^2\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_1(r,z)-\frac{{2\mathrm{\&kappa;}}^2}{c_0\left(z\right)<{\mathrm{\&Lambda;}}_1^2>}\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\underset{s=0}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_s(r,z){\mathrm{\&gamma;}}_q(r,z)<{\mathrm{\&Lambda;}}_q{\mathrm{\&Lambda;}}_s{\mathrm{\&Lambda;}}_1>=0---\left(10\right)$

以高斯分布随机声速进行下面公式的推导，δc(r，z，θ)＝α₀He₀+α₁He₁+…＝σξ，ξ是均值为0，方差为1的高斯随机变量；将δc代入(10)式给出第1个不确实系数的耦合方程：

$\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;z}^2}+{\mathrm{\&kappa;}}^2{\mathrm{\&gamma;}}_1(r,z)-\frac{{2\mathrm{\&kappa;}}^2\mathrm{\&sigma;}}{c_0<{\mathrm{He}}_1^2>}\underset{q=0}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{\mathrm{\&gamma;}}_q(r,z)<H{e}_q\mathrm{\&xi;}{\mathrm{He}}_1>=0---\left(11\right)$

利用Hermite基函数的回归特性和正交性，计算耦合项的系总平均，得到系数γ₁的偏微分回归方程：

$\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{\mathrm{\&gamma;}}_1(r,z)}{{\&PartialD;z}^2}+{\mathrm{D\&gamma;}}_1(r,z)-E\sqrt{l}{\mathrm{\&gamma;}}_{l-1}(r,z)-E\sqrt{l+1}{\mathrm{\&gamma;}}_{l+1}(r,z)=0---\left(12\right)$

式中D≡κ²，E≡2κ²σ/c₀；用向量-矩阵形式，上述偏微分方程可化简为

$\frac{{\&PartialD;}^2\mathrm{\&gamma;}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2\mathrm{\&gamma;}}{{\&PartialD;z}^2}+\mathrm{A\&gamma;}=0---\left(13\right)$

式中，A为环境参量不确实性到声场不确实性的耦合方程，A是对称的三角阵；

(2)求解随机波方程得到不确实声场

步骤一：耦合矩阵A特征值分解

$A=\underset{q=1}{\overset{Q}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_q{g}_q{g}_q^T---\left(14\right)$

式中，上标T表示转置，G＝[{g₁}，{g₂}，…]为耦合矩阵A的特征向量矩阵，GG^-1＝I，(13)式可以表示为

$\frac{{\&PartialD;}^2{G}^{-1}\mathrm{\&gamma;}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{G}^{-1}\mathrm{\&gamma;}}{{\&PartialD;z}^2}+G^{-1}{\mathrm{AGG}}^{-1}\mathrm{\&gamma;}=0---\left(15\right)$

令

上式可以简化为

$\frac{{\&PartialD;}^2\widehat{\mathrm{\&gamma;}}}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2\widehat{\mathrm{\&gamma;}}}{{\&PartialD;z}^2}+\mathrm{\&Omega;}\widehat{\mathrm{\&gamma;}}=0---\left(16\right)$

式中，Ω＝G^-1AG为对角阵，对角阵元素为特征值λ₁，λ₂，…，利用特征值和特征向量分解技术，得到第1个确定系数的波方程

$\frac{{\&PartialD;}^2{\widehat{\mathrm{\&gamma;}}}_1}{{\&PartialD;r}^2}+\frac{{\&PartialD;}^2{\widehat{\mathrm{\&gamma;}}}_1}{{\&PartialD;z}^2}+{\mathrm{\&lambda;}}_1{\widehat{\mathrm{\&gamma;}}}_1=0,---\left(17\right)$

步骤二：确定系数求解

等式(17)与标准波方程(6)有相同的形式，因而求解过程与标准的波方程求解过程是一样的；

的解表示为

${\widehat{\mathrm{\&gamma;}}}_u(r,z)=\frac{i}{\sqrt{8\mathrm{\&pi;r}}}{e}^{-i\frac{\mathrm{\&pi;}}{4}}{a}_{u0}\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&psi;}}_m\left(z_s\right){\mathrm{\&psi;}}_m\left(z\right)}{\sqrt{k_{\mathrm{rm}}}}{e}^{{\mathrm{ik}}_{\mathrm{rm}}r}---\left(18\right)$

式中k_rm、ψ_m(z)分别是水平波数和模深度函数，由简正波传播模型求解得到，a_u0表示矩阵A的特征向量的逆的第一列中的第u个元素，在声速梯度不确实情况下，模深度函数没有受到不确实声速的影响，但声速梯度不确实导致了水平波数的不确实，

其中k_zm是深度波数，根据γ₁与

的关系式，系数γ₁是的线性组合

${\mathrm{\&gamma;}}_1=\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{1u}{\widehat{\mathrm{\&gamma;}}}_u(r,z)---\left(19\right)$

步骤三：Q个多项式混沌基函数的确定性系数加权累加得到不确实性声场

最后，在远场近似下的不确实性声场可表示为：

$p(r,z;\mathrm{\&theta;})=\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\left[\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{\mathrm{qu}}{\widehat{\mathrm{\&gamma;}}}_u(r,z)\right]{\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)$

(20)

$=\underset{q=0}{\overset{Q}{\mathrm{\&Sigma;}}}\left[\underset{u=0}{\overset{Q}{\mathrm{\&Sigma;}}}{g}_{\mathrm{qu}}\left(\frac{i}{\sqrt{8\mathrm{\&pi;r}}}{e}^{-i\frac{\mathrm{\&pi;}}{4}}{a}_{u0}\underset{m=1}{\overset{M}{\mathrm{\&Sigma;}}}\frac{{\mathrm{\&psi;}}_m\left(z_s\right){\mathrm{\&psi;}}_m\left(z\right)}{\sqrt[4]{{\mathrm{\&lambda;}}_u-k_{\mathrm{zm}}^2}}{e}^{i\sqrt{{\mathrm{\&lambda;}}_u-k_{\mathrm{zm}}^2}r}\right)\right]{\mathrm{\&Lambda;}}_q\left(\mathrm{\&xi;}\left(\mathrm{\&theta;}\right)\right)$

若SSP是确实的，则上式中只有γ₀是不等于0的，此时式(20)收敛到声速确知时的声场简正模表达式；

步骤四：不确实性声场的均值和相关函数

由多项式混沌展开的不确实声场是随机声场，其均值可表示为

$<p(r,z;\mathrm{\&theta;})>={\widehat{\mathrm{\&gamma;}}}_0(r,z)---\left(21\right)$

相关函数为

$R_{\mathrm{pp}}(r_1,z_1;r_2,z_2)=\underset{q=1}{\overset{Q}{\mathrm{\&Sigma;}}}{\widehat{\mathrm{\&gamma;}}}_q*(r_1,z_1){\widehat{\mathrm{\&gamma;}}}_q(r_2,z_2)---\left(22\right).$

2.根据权利要求1所述的海洋环境与声场不确实性的表征和传递的快速计算方法，其特征在于：b是高斯分布，对应的Λ_s(ξ(θ))为Gauss-Hermite多项式；b是均匀分布，对应的Λ_s(ξ(θ))为Legendre多项式；b是γ分布，对应的Λ_s(ξ(θ))为Laguerre多项式，b是β分布，对应的Λ_s(ξ(θ))为Jacobi多项式，b是泊松分布，对应的Λ_s(ξ(θ))为Charlier多项式。

3.根据权利要求2所述的海洋环境与声场不确实性的表征和传递的快速计算方法，其特征在于：对于高斯分布的b，不确实环境参量的展开只要一阶就够，其中α₀(r，z)等于0，α₁(r，z)＝σ，式中σ是随机扰动环境参量的标准离差；具有类似的特性，其它分布的环境参量也可以利用投影方法得到有限S阶的确定性系数。

Images