CN114386296B - 一种混响水池中三维声场的数值计算方法 - Google Patents

Abstract

本发明公开一种混响水池中三维声场的数值计算方法。步骤1：构建混响水池模型；步骤2：在混响水池中设有一个声源S₀，在池壁表面放置一个矢量水听器R₀，基于步骤1的混响水池模型，直达声线为轴，对混响水池中的声场进行平面划分，得到声场的截面；步骤3：基于步骤2的声场截面，利用虚源法绘制声场平面的声线传播情况，构建出声线在声场平面的虚源图像及其反射声线模型；步骤4：基于步骤3的虚源图像及其反射声线模型构建出声源位于混响水池中心位置时的三维声场。本发明针对现有技术中对对水池内的声场进行数值计算时，用时长，效率低的问题。

Description

一种混响水池中三维声场的数值计算方法

技术领域

本发明涉及声场的构建技术领域，具体涉及一种混响水池中三维声场的数值计算方法。

背景技术

在水声工程领域，常使用水池对声纳或换能器性能进行测试与评估。由于建设周期和成本限制，多数水池的池壁上都未覆盖吸声橡胶，而是直接使用瓷砖或防水涂料。这种池壁会对声信号进行多次反射后，水池中的声场较为复杂。对水池中的声场进行计算时，不光要考虑声线在水面与池底的反射，还需要考虑声线在三维空间场中与四周池壁的反射。

现有的方法多使用有限元分析软件可以对水池内的声场进行数值计算，其计算过程较为复杂，计算时间较长。单次仿真仅能对单个频率下若干个点上的声压幅值进行计算，计算效率不高。

发明内容

本发明提供一种混响水池中三维声场的数值计算方法，针对现有技术中对对水池内的声场进行数值计算时，用时长，效率低的问题。

本发明通过以下技术方案实现：

一种混响水池中三维声场的数值计算方法，所述数值计算方法包括以下步骤：

步骤1：构建混响水池模型，所述混响水池为一个底面为正方形，四周侧面为矩形的长方体且混响水池壁内表面光滑的方体混响水池；

步骤2：在混响水池中设有一个声源S₀，在池壁表面放置一个矢量水听器R₀，基于步骤1的混响水池模型，直达声线为轴，对混响水池中的声场进行平面划分，得到声场的截面；

步骤3：基于步骤2的声场截面，利用虚源法绘制声场平面的声线传播情况，构建出声线在声场平面的虚源图像及其反射声线模型；

步骤4：基于步骤3的虚源图像及其反射声线模型构建出声源位于混响水池中心位置时的三维声场。

进一步的，所述步骤2具体包括以下两种情况：

第一种：混响水池正中心有一个声源S₀，在其正对的某一个池壁表面中心放置一个矢量水听器R₀，以直达声线为轴，对混响水池中的声场进行平面划分，得到声场的水平截面与垂直截面；

第二种：混响水池中有一个声源S₀，放置一个矢量水听器R₀，即声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分。

进一步的，所述步骤1具体为，声源在其中发射时，声波经过各个侧壁的反射与叠加，形成混响声场，混响水池的长、宽、高分别为l_x、l_y、l_z，混响水池的四面侧壁设为平面Ⅰ1、平面Ⅱ2、平面Ⅲ3和平面Ⅳ4，上表面为平面Ⅴ5，池底为平面Ⅵ6。

进一步的，基于水平截面与垂直截面，设真实声源位置为S₀，其关于反射界面的镜像虚源为S_0n；矢量水听器为R₀，各个虚源到达矢量水听器的距离为R_0n，各条声线的声程长，到达时延为τ_0n，

若假设水池声速为c，直达声时延为τ₀，声源S₀与四个侧壁的距离为k₁，与池底及水面的距离为k₂，声线经水池侧壁的反射次数为m₀₁，声线经池底与水面的反射次数为m₀₂，由勾股定理计算出各条声线的虚源与接收点的距离为：

R₀S₀＝c·τ₀ (1)

进一步的，所述声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分，得到声场的截面l₁与截面t₁，其中直线y₁为平面l₁与平面t₁的相交线，S₁与R₁均位于相交线上，

矢量水听器的真实位置为R₁，真实声源位置为S₁；R₁S₁为直达声线，直达时延为τ₁，直达声线R₁S₁与x轴的夹角为θ；真实声源关于反射界面的镜像虚源为S_1n，矢量水听器与各个镜像虚源的距离为R_1n，即各条声线的声程，时延设为τ_1n；

假设声源S₁与边界2的距离为h₁，与边界4的距离为h₂，与边界3的距离为m₁，矢量水听器R₁与边界1的间距为m₂；截面l₁中水池边界的长宽长度分别为a＝h₁+h₂与b＝m₁+m₂+R₁S₁·cosθ；另外，根据声线在六个平面的反射次数的不同，分别将声线与平面Ⅰ1、平面Ⅱ2、平面Ⅲ3、平面Ⅳ4、平面Ⅴ5和平面Ⅵ6的反射次数用参数n₁、n₂、n₃、n₄、n₅以及n₆来表示；根据每条声线与不同边界面的反射次数次数的不同，由勾股定理得到不同的声线的声程为：

(9)直达声线R₁S₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝0，R₁S₁＝c·τ₁＝L₁

(10)R₁S₁₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝1

R₁₁＝(L₁ cosθ)²+(L₁ sinθ+2h₂)²＝(L₁ cosθ)²+(L₁ sinθ+n₄·2h₂)²

(11)R₁S₁₂，即n₁＝0，n₂＝1，n₃＝0，n₄＝0

R₁₂＝(L₁ cosθ)²+(2h₁-L₁ sinθ)²＝(L₁ cosθ)²+(n₂·2h₁-L₁ sinθ)²

(12)R₁S₁₃，即n₁＝0，n₂＝0，n₃＝1，n₄＝0

R₁₃＝(L₁ cosθ+2m₁)²+(L₁ sinθ)²＝(L₁ cosθ+n₃·2m₁)²+(L₁ sinθ)²

(13)R₁S₁₄，即n₁＝1，n₂＝0，n₃＝0，n₄＝0

R₁₄＝(L₁ cosθ+2m₂)²+(L₁ sinθ)²＝(L₁ cosθ+n₁·2m₂)²+(L₁ sinθ)²

(14)R₁S₁₅，即n₁＝0，n₂＝1，n₃＝0，n₄＝1

R₁₅＝(L₁ cosθ)²+(L₁ sinθ+2h₁+2h₂)²＝(L₁ cosθ)²+(L₁ sinθ+n₂·2h₁+n₄·2h₂)²

(15)R₁S₁₆，即n₁＝0，n₂＝1，n₃＝1，n₄＝0

R₁₆＝(L₁ cosθ+2m₁)²+(2h₁-L₁ sinθ)²＝(L₁ cosθ+n₃·2m₁)²+(n₂·2h₁-L₁ sinθ)²

(16)R₁S₁₇，即n₁＝0，n₂＝1，n₃＝1，n₄＝1

R₁₇＝(L₁cosθ+2m₁)²+(L₁sinθ+2h₁+2h₂)²

＝(L₁cosθ+n₃·2m₁)²+(L₁ sinθ+n₂·2h₁+n₄·2h₂)²

(17)R₁S₁₈，即n₁＝1，n₂＝0，n₃＝1，n₄＝1

R₁₈＝(L₁cosθ+2m₁+2m₂)²+(L₁ sinθ+2h₂)²

＝(L₁cosθ+n₃·2m₁+n₁·2m₂)²+(L₁ sinθ+n₄·2h₂)²

(18)R₁S₁₉，即n₁＝1，n₂＝1，n₃＝1，n₄＝2

进一步的，结合虚源法所列出的公式(1)-(8)，推测在平面l₁内，混响水池中虚源与接收点距离的表达式P_n如下：

P_n＝(2n₃m₁+2n₁m₂+L₁ cosθ)²+(2n₂h₁+2n₄h₂±L₁ sinθ)²

其中L₁ sinθ项的系数为正或为负由最终反射界面决定。

进一步的，根据平面l₁的声场情况，可得出平面t₁中虚源与接收点距离的表达式Q_n如下：

Q_n＝(2n₃m₃+2n₁m₄+L₁ cosα)²+(2n₅h₃+2n₆h₄±L₁ sinα)²。

进一步的，所述步骤4具体为，将得到的声线长度代入到虚源辐射声压的计算公式里并结果进行叠加，便可得到这平面t₁与平面l₁内的所有镜像虚源在接收点的总声压；

对接收点的每条声线与直达声线轴所形成的空间角进行积分的方式，来达到对整个空间的虚源的辐射声压进行叠加的目的。

进一步的，所述对接收点的每条声线与直达声线轴所形成的空间角进行积分具体包括以下步骤：

步骤4.1：设R₂为接收水听器，S₂为声源，两者位于水平面l₂与垂直平面t₂的相交直线上；S’₂与S”₂为声源关于边界的镜像虚源；

步骤4.2：BC与直线R₂S₂垂直，R₂A与直线S’₂S₂垂直，R₂S₂与直线S”₂S₂垂直，DB垂直平面l₂；R₂S’₂为平面l₂上声线的声程，R₂S”₂为平面t₂上反射声线的声程；在空间上存在多个声场平面，而每个平面上，所有经过多次声线反射所产生的镜像虚源均位于平面r₂上；

步骤4.3：基于步骤4.1与步骤4.2的定义可知，若令∠AR₂S₂＝θ₀，由于R₂A与直线S’₂S₂垂直且BC与直线R₂S₂垂直，因此∠CBS₂＝θ₀，同样的，设声源R₂与平面Ⅱ2和平面Ⅳ4的距离分别为h₁与h₂，与平面Ⅴ5和平面Ⅵ6的距离分别为h₃与h₄，R₂S₂＝d₀，R₂S’₂＝d₁，R₂S”₂＝d₂；则由几何关系，得到如下等式：

上式看作虚源与接收点距离通式的变形；

步骤4.4：若平面l₂以直线R₂S₂为轴进行旋转，则B点的运动轨迹为沿着BD方向的直线，BC即为DC在平面l₂上的映射，BS₂即为DS₂在平面l₂上的映射，而∠DBC则为旋转后的平面l’₂与原平面l₂的夹角，若令∠DBC＝β，∠CDS₂＝θ，DS₂＝h，则根据几何关系有：

步骤4.5：对一次反射声线而言，基于公式(20)和(21)等式在任何声场平面上均成立；因此若以h代替声源与边界2的距离h₁，并将以上两式代入其中，得到虚源与接收点距离关于角度θ的等式如下：

步骤4.6：在考虑反射系数的情况下，单个虚源在接收点的辐射声压为：

Ω(φ)＝kdsinφsinε Ω_z(φ)＝kdcosφcosε

其中J₀(·)为零阶第一类贝塞尔函数，V(θ)为液体，弹性界面的反射系数计算公式；k为波束矢量，k＝ω/c，ω为角频率，c为介质声速；Z₁、Z₂为反射界面两侧的声阻抗，

为声线入射角，/>

为折射角；φ为声线传播矢量与空间直角坐标系的z轴的夹角，ε为声线传播矢量在xOy面的投影与x轴的夹角，d为虚源到达接收点的距离；

步骤4.7：基于步骤4.6的单个虚源在接收点的辐射声压为只是一个声场平面的一个虚源在接收点的辐射声压，若要得到整个声场的虚源在接收点的声压总和，将所有由通式算得的声线长度，代入到虚源辐射声压的计算公式，并对声场平面绕直达声线旋转的旋转角β进行积分，即

其中p_n(d_n)为单个虚源辐射声压，θ₁与θ₂为积分角θ的积分范围。

步骤4.8：将步骤4.6的公式中cosθ与sinθ做如下替换：

而由平面l₁的旋转，则确定角度β的积分范围分别为

和/>

对三个角度范围分别进行积分，便得到整个混响场中的虚源在接收点R处的声压；/>

算出空间中每个点的声场声压，从而构建出三维声场。

本发明的有益效果是：

本发明计算简单，运算量小，可对宽带声场进行快速计算。

本发明为混响水池的数值构建提供了一个较为简便的理论模型，一定程度上解决了混响水池的声场声压较难计算的问题。

本发明将二维的虚源法拓展到三维的空间中，以空间几何的方式来直观的分析声线在混响水池的反射情况，简化了混响水池的边界问题。

本发明避免了建立理想水下扩散场等较为复杂的问题，仅以几何方式来分析整个声场的声传播。

附图说明

图1是本发明所构建的混响水池模型。

图2是本发明对声源位于水池正中间的情况进行三维声场平面划分示意图。

图3是本发明中混响声场平面l₀的虚源图像及其反射声线。

图4是本发明中混响声场平面t₀的虚源图像及其反射声线。

图5是本发明中声源位于特殊位置时的水池声场的平面划分。

图6是本发明中混响声场平面l₁的虚源图像及其反射声线。

图7是本发明中混响声场平面t₁的虚源图像及其反射声线。

图8是本发明中混响水池的三维空间声线与声场平面示意图。

图9是本发明所构建混响水池中声源与矢量水听器相对位置示意图。

图10是本发明中水池的接收点声压随频率变化曲线。

具体实施方式

下面将结合本发明实施例中的附图对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

一种混响水池中三维声场的数值计算方法，所述数值计算方法包括以下步骤：

步骤1：构建混响水池模型，所述混响水池为一个底面为正方形，四周侧面为矩形的长方体且混响水池壁内表面光滑的方体混响水池；

步骤2：在混响水池中设有一个声源S₀，在池壁表面放置一个矢量水听器R₀，基于步骤1的混响水池模型，直达声线为轴，对混响水池中的声场进行平面划分，得到声场的截面；

步骤3：基于步骤2的声场截面，利用虚源法绘制声场平面的声线传播情况，构建出声线在声场平面的虚源图像及其反射声线模型；

步骤4：基于步骤3的虚源图像及其反射声线模型构建出声源位于混响水池中心位置时的三维声场。

进一步的，所述步骤2具体包括以下两种情况：

第一种：混响水池正中心有一个声源S₀，在其正对的某一个池壁表面中心放置一个矢量水听器R₀，以直达声线为轴，对混响水池中的声场进行平面划分，得到声场的水平截面与垂直截面；

第二种：混响水池中有一个声源S₀，放置一个矢量水听器R₀，即声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分。

对于声源与接收水听器均位于一般位置时，不同的声线，其声线长度的计算不仅与声源和各个壁面的间距有关，还与声线在不同界面的反射次数有关。下面便针对这种情况，在声源位于水池中心时所得结果的基础之上，利用虚源法进行三维声场构建方法的推导。

进一步的，所述步骤1具体为，声源在其中发射时，声波经过各个侧壁的反射与叠加，形成混响声场，混响水池的长、宽、高分别为l_x、l_y、l_z，混响水池的四面侧壁设为平面Ⅰ1、平面Ⅱ2、平面Ⅲ3和平面Ⅳ4，上表面为平面Ⅴ5，池底为平面Ⅵ6；

一般情况下水池内壁以及底部都有瓷砖，但在这里所构建的水池模型，其内表面边界设为水——混凝土边界，其相对声导纳为α＝ρ₀c₀/Z(ρ₀，c₀，Z分别为声波在水中传播的速度及水池壁的声阻抗)，即在后续的计算过程中需要考虑声线在经过池壁反射后所产生的能量损失；而水池的上表面为自由边界，以绝对软边界来近似。

进一步的，基于水平截面与垂直截面，设真实声源位置为S₀，其关于反射界面的镜像虚源为S_0n；矢量水听器为R₀，各个虚源到达矢量水听器的距离为R_0n，各条声线的声程长，到达时延为τ_0n，

若假设水池声速为c，直达声时延为τ₀，声源S₀与四个侧壁的距离为k₁，与池底及水面的距离为k₂，声线经水池侧壁的反射次数为m₀₁，声线经池底与水面的反射次数为m₀₂，由勾股定理计算出各条声线的虚源与接收点的距离为：

R₀S₀＝c·τ₀ (1)

进一步的，所述声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分，得到声场的截面l₁与截面t₁，其中直线y₁为平面l₁与平面t₁的相交线，S₁与R₁均位于相交线上，

矢量水听器的真实位置为R₁，真实声源位置为S₁；R₁S₁为直达声线，直达时延为τ₁，在平面l₁，也就是图6上的直达声线R₁S₁与x轴的夹角为θ；图中绘制了9条声线的反射情况，与之前类似，真实声源关于反射界面的镜像虚源为S_1n，矢量水听器与各个镜像虚源的距离为R_1n，即各条声线的声程，时延设为τ_1n；

假设声源S₁与边界2(是指平面Ⅱ2上的一条直线，图6中的平面l₁是声场的截面，这样平面l₁中所示的边界2其实是平面l₁与平面Ⅱ2，这两个平面相交的一条直线，而这个距离其实是声源S₁到这条相交直线的距离)的距离为h₁，与边界4(是指平面l₁与平面Ⅳ4相交而形成的相交线)的距离为h₂，与边界3(是指平面l₁与平面Ⅲ3相交而形成的相交线)的距离为m₁，矢量水听器R₁与边界1(是指平面l₁与平面Ⅰ1相交而形成的相交线)的间距为m₂；截面l₁中水池边界的长宽长度分别为a＝h₁+h₂与b＝m₁+m₂+R₁S₁·cosθ；另外，根据声线在六个平面的反射次数的不同，分别将声线与平面Ⅰ1、平面Ⅱ2、平面Ⅲ3、平面Ⅳ4、平面Ⅴ5和平面Ⅵ6的反射次数用参数n₁、n₂、n₃、n₄、n₅以及n₆来表示；根据每条声线与不同边界面的反射次数次数的不同，由勾股定理得到图6中的不同的声线的声程为：

(9)直达声线R₁S₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝0，R₁S₁＝c·τ₁＝L₁

(10)R₁S₁₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝1

R₁₁＝(L₁ cosθ)²+(L₁ sinθ+2h₂)²＝(L₁ cosθ)²+(L₁ sinθ+n₄·2h₂)²

(11)R₁S₁₂，即n₁＝0，n₂＝1，n₃＝0，n₄＝0

R₁₂＝(L₁ cosθ)²+(2h₁-L₁ sinθ)²＝(L₁ cosθ)²+(n₂·2h₁-L₁ sinθ)²

(12)R₁S₁₃，即n₁＝0，n₂＝0，n₃＝1，n₄＝0

R₁₃＝(L₁ cosθ+2m₁)²+(L₁ sinθ)²＝(L₁ cosθ+n₃·2m₁)²+(L₁ sinθ)²

(13)R₁S₁₄，即n₁＝1，n₂＝0，n₃＝0，n₄＝0

R₁₄＝(L₁ cosθ+2m₂)²+(L₁ sinθ)²＝(L₁ cosθ+n₁·2m₂)²+(L₁ sinθ)²

(14)R₁S₁₅，即n₁＝0，n₂＝1，n₃＝0，n₄＝1

R₁₅＝(L₁ cosθ)²+(L₁ sinθ+2h₁+2h₂)²＝(L₁ cosθ)²+(L₁ sinθ+n₂·2h₁+n₄·2h₂)²

(15)R₁S₁₆，即n₁＝0，n₂＝1，n₃＝1，n₄＝0

R₁₆＝(L₁ cosθ+2m₁)²+(2h₁-L₁ sinθ)²＝(L₁ cosθ+n₃·2m₁)²+(n₂·2h₁-L₁ sinθ)²

(16)R₁S₁₇，即n₁＝0，n₂＝1，n₃＝1，n₄＝1

R₁₇＝(L₁cosθ+2m₁)²+(L₁sinθ+2h₁+2h₂)²

＝(L₁cosθ+n₃·2m₁)²+(L₁ sinθ+n₂·2h₁+n₄·2h₂)²

(17)R₁S₁₈，即n₁＝1，n₂＝0，n₃＝1，n₄＝1

R₁₈＝(L₁cosθ+2m₁+2m₂)²+(L₁ sinθ+2h₂)²

＝(L₁cosθ+n₃·2m₁+n₁·2m₂)²+(L₁ sinθ+n₄·2h₂)²

(18)R₁S₁₉，即n₁＝1，n₂＝1，n₃＝1，n₄＝2

进一步的，结合虚源法所列出的公式(1)-(8)，推测在平面l₁内，混响水池中虚源与接收点距离的表达式P_n如下：

P_n＝(2n₃m₁+2n₁m₂+L₁ cosθ)²+(2n₂h₁+2n₄h₂±L₁ sinθ)²

其中L₁ sinθ项的系数为正或为负由最终反射界面决定。

进一步的，根据平面l₁的声场情况，可得出平面t₁中虚源与接收点距离的表达式Q_n如下：

Q_n＝(2n₃m₃+2n₁m₄+L₁ cosα)²+(2n₅h₃+2n₆h₄±L₁ sinα)²。

为了验证以上表达式是否正确，现利用数学归纳法对其进行推导证明，推导过程中，引入参数k来代替参数n。

(1)首先，假设当n＝1时(此处为了便于推导，以单独的参数n来代替整体的多个反射参数)，也即n₁＝1，n₂＝1，n₃＝1，n₄＝1时，令以下等式成立

P₁＝(2m₁+2m₂+L₁ cosθ)²+(2h₁+2h₂±L₁ sinθ)²

(2)若令n＝k时，P_k-1及P_k也成立，即

P_k-1＝[2(k₃-1)m₁+2(k₁-1)m₂+L₁ cosθ]²+[2(k₂-1)h₁+2(k₄-1)h₂±L₁ sinθ]²

(3)则当n＝k+1时，P_k+1的表达式如下所示

最后，得证，通式P_n成立。

进一步的，所述步骤4具体为，将得到的声线长度代入到虚源辐射声压的计算公式里并结果进行叠加，便可得到这平面t₁与平面l₁内的所有镜像虚源在接收点的总声压；

但是对于整个声场空间来说，仅仅是两个平面的虚源的辐射声压进行叠加是远远不够的，但又不可能将每个声场平面以及其声线反射情况都绘制出来。

对接收点的每条声线与直达声线轴所形成的空间角进行积分的方式，来达到对整个空间的虚源的辐射声压进行叠加的目的。

不过对于不同的声线，其与直达声线的空间角则较难求出，也不易表示，因此需要利用空间几何对其进行代换或者寻找一个新的积分量，以简化积分计算。

进一步的，所述对接收点的每条声线与直达声线轴所形成的空间角进行积分具体包括以下步骤：

在单个平面中，声源到达接收点的各条声线所形成的镜像虚源，分布在与反射边界垂直且互相平行的轴线上，其中的每个虚源与接收点的距离均可由P_n得出。但是如果放在整个空间平面上来看的话，这些虚源应位于互相平行的空间平面上，并且反射边界相同的声线(不需要声程也完全一样，仅需要反射边界相同即可)，其虚源也位于同样一平面。每条空间中的声线均可看作位于某一声场平面之上，其声程仍可利用该平面的通式来计算；而不同平面的声程计算通式可利用平面间夹角，建立几何关系来互相表示。由此，可认为，每个声场平面可看作是由另一个声场平面绕直达声线轴旋转了角度β得到的。

步骤4.1：建立三维声场示意图如图8所示，设R₂为接收水听器，S₂为声源，两者位于水平面l₂与垂直平面t₂的相交直线上；S’₂与S”₂为声源关于边界的镜像虚源；

步骤4.2：B点，C点(图8中，B点是直线S₂S₂’与平面Ⅱ2的交点。而直线S₂S₂’是水平面l₂与网格平面r₂的交线。而C点则是从B点向直线R₂S₂做垂线的垂足)BC与直线R₂S₂垂直，R₂A与直线S’₂S₂垂直，R₂S₂与直线S”₂S₂垂直，B点，D点(D点是B点的运动轨迹上任何一个点的统称；当平面l₂以直线R₂S₂为轴进行旋转时，直线S₂S₂’的位置也会发生相应的改变，但始终是平面l₂与网格平面r₂的交线，网格平面r₂的位置是不变的；直线S₂S₂’的位置变化也就造成B点的变化，而BD就是整个B点的变化轨迹，D可以是这条轨迹上的任何一个点。而轨迹BD是平面Ⅱ2上的一条垂直于水平面l₁的，或者更准确的说，是一条垂直于底部平面Ⅵ6的垂线，垂足为B)DB垂直平面l₂；R₂S’₂为平面l₂上声线的声程，R₂S”₂为平面t₂上反射声线的声程；在空间上存在多个声场平面，而每个平面上，所有经过多次声线反射(这里所说的多次反射仅指与截面声场中的上下边界有关的反射声线，即与平面2、5、4、6有关的反射声线)所产生的镜像虚源均位于平面r₂上；

步骤4.3：基于步骤4.1与步骤4.2的定义可知，若令∠AR₂S₂＝θ₀，由于R₂A与直线S’₂S₂垂直且BC与直线R₂S₂垂直，因此∠CBS₂＝θ₀，同样的，设声源R₂与平面Ⅱ2和平面Ⅳ4的距离分别为h₁与h₂，与平面Ⅴ5和平面Ⅵ6的距离分别为h₃与h₄，R₂S₂＝d₀，R₂S’₂＝d₁，R₂S”₂＝d₂；则由几何关系，得到如下等式：

上式看作虚源与接收点距离通式的变形；

步骤4.4：若平面l₂以直线R₂S₂为轴进行旋转，则B点的运动轨迹为沿着BD方向的直线，BC即为DC在平面l₂上的映射，BS₂即为DS₂在平面l₂上的映射，而∠DBC则为旋转后的平面l’₂与原平面l₂的夹角，若令∠DBC＝β，∠CDS₂＝θ，DS₂＝h，则根据几何关系有：

步骤4.5：对一次反射声线而言，基于公式(20)和(21)等式在任何声场平面上均成立；因此若以h代替声源与边界2的距离h₁，并将以上两式代入其中，得到虚源与接收点距离关于角度θ的等式如下：

步骤4.6：由于本发明中所涉及的声源为点声源，其声传播在空间上近似为球面波，因此若要求得单个虚源在接收点处的辐射声场，就考虑球面波在水平界面上的反射情况，即在考虑反射系数的情况下，单个虚源在接收点的辐射声压为：

Ω(φ)＝kdsinφsinε Ω_z(φ)＝kdcosφcosε

其中J₀(·)为零阶第一类贝塞尔函数，V(θ)为液体，弹性界面的反射系数计算公式；k为波束矢量，k＝ω/c，ω为角频率，c为介质声速；Z₁、Z₂为反射界面两侧的声阻抗，

为声线入射角，/>

为折射角；φ为声线传播矢量与空间直角坐标系的z轴的夹角，ε为声线传播矢量在xOy面的投影与x轴的夹角，d为虚源到达接收点的距离；

步骤4.7：基于步骤4.6的单个虚源在接收点的辐射声压为只是一个声场平面的一个虚源在接收点的辐射声压，若要得到整个声场的虚源在接收点的声压总和，将所有由通式算得的声线长度，代入到虚源辐射声压的计算公式，并对声场平面绕直达声线旋转的旋转角β进行积分，即

其中p_n(d_n)为单个虚源辐射声压，θ₁与θ₂为积分角θ的积分范围。

步骤4.8：对于本发明来说，θ的变化范围不易得出，而β的变化则更容易判断，因此根据前边已推导出的几何关系，可以将步骤4.6的公式中cosθ与sinθ做如下替换：

而由平面l₁的旋转，则确定角度β的积分范围分别为

和/>

对三个角度范围分别进行积分，便得到整个混响场中的虚源在接收点R处的声压；

算出空间中每个点的声场声压，从而构建出三维声场。

下面将本发明方法的计算结果与有限元分析软件Comsol的仿真结果来进行对比。计算水池材质为混凝土，考虑反射边界为水——混凝土的情况，选用混凝土参数：密度为ρ₀＝2.6×103kg/m3，介质声速为v＝4000m/s。水池池底的外表面长宽均为11m，水池高为6m，池壁厚度为0.5m，水池内表面池底长宽均为10m，水域高度为5m。现以水池内表面的其中的一个底角为坐标原点，建立空间直角坐标系Oxyz，并在水池深为3m，与距离侧壁分别为2m与8m的位置放置一个功率为1W的声源S，即声源坐标为S(8，2，2),同时分别在该水池中的四个不同的位置分别放置一个水听器，其位置坐标分别为Rs₁(2，2，2)，Rs₂(2，8，2)，Rs₃(8，8，2)，Rs₄(2，8，4)。最终所建水池的模型模型示意图如下图9所示(为了较为清楚的看到图中各个水听器与声源的相对位置，因此将图中前面的两个侧壁已隐藏)。

现分别利用三维空间虚源法以及有限元分析软件Comsol，对四个接收水听器所处位置的声压进行计算。本发明在进行虚源辐射声压计算时，仅考虑声线反射次数小于等于4次的情况，忽略其他反射声线，因此其计算结果与Comsol软件算得的结果如图10所示。

图10中，实线为Comsol软件的声压计算结果，虚线为本发明方法所计算结果。从图中可以看出，本发明方法计算得到的辐射声压相较于利用Comsol直接算得的声场声压吻合较好。

Claims (8)

1.一种混响水池中三维声场的数值计算方法，其特征在于，所述数值计算方法包括以下步骤：

步骤1：构建混响水池模型，所述混响水池为一个底面为正方形，四周侧面为矩形的长方体且混响水池壁内表面光滑的方体混响水池；

步骤2：在混响水池中设有一个声源S₀，在池壁表面放置一个矢量水听器R₀，基于步骤1的混响水池模型，直达声线为轴，对混响水池中的声场进行平面划分，得到声场的截面；

步骤3：基于步骤2的声场截面，利用虚源法绘制声场平面的声线传播情况，构建出声线在声场平面的虚源图像及其反射声线模型；

步骤4：基于步骤3的虚源图像及其反射声线模型构建出声源位于混响水池中心位置时的三维声场。

2.根据权利要求1所述一种混响水池中三维声场的数值计算方法，其特征在于，所述步骤1具体为，声源在其中发射时，声波经过各个侧壁的反射与叠加，形成混响声场，混响水池的长、宽、高分别为l_x、l_y、l_z，混响水池的四面侧壁设为平面Ⅰ(1)、平面Ⅱ(2)、平面Ⅲ(3)和平面Ⅳ(4)，上表面为平面Ⅴ(5)，池底为平面Ⅵ(6)。

3.根据权利要求1所述一种混响水池中三维声场的数值计算方法，其特征在于，所述步骤2具体包括以下两种情况：

第一种：混响水池正中心有一个声源S₀，在其正对的某一个池壁表面中心放置一个矢量水听器R₀，以直达声线为轴，对混响水池中的声场进行平面划分，得到声场的水平截面与垂直截面；

第二种：混响水池中有一个声源S₀，放置一个矢量水听器R₀，即声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分。

4.根据权利要求3所述一种混响水池中三维声场的数值计算方法，其特征在于，基于水平截面与垂直截面，设真实声源位置为S₀，其关于反射界面的镜像虚源为S_0n；矢量水听器为R₀，各个虚源到达矢量水听器的距离为R_0n，各条声线的声程长，到达时延为τ_0n，

若假设水池声速为c，直达声时延为τ₀，声源S₀与四个侧壁的距离为k₁，与池底及水面的距离为k₂，声线经水池侧壁的反射次数为m₀₁，声线经池底与水面的反射次数为m₀₂，由勾股定理计算出各条声线的虚源与接收点的距离为：

R₀S₀＝c·τ₀ (1)

5.根据权利要求2或3所述一种混响水池中三维声场的数值计算方法，其特征在于，所述声源S₀与矢量水听器R₀距离各个池壁以及水面的间距都不相同，声源S₀与矢量水听器R₀也未处于同一水平面，以直达声线为轴，对声场进行平面划分，得到声场的截面l₁与截面t₁，其中直线y₁为平面l₁与平面t₁的相交线，S₁与R₁均位于相交线上，

矢量水听器的真实位置为R₁，真实声源位置为S₁；R₁S₁为直达声线，直达时延为τ₁，直达声线R₁S₁与x轴的夹角为θ；真实声源关于反射界面的镜像虚源为S_1n，矢量水听器与各个镜像虚源的距离为R_1n，即各条声线的声程，时延设为τ_1n；

假设声源S₁与边界2的距离为h₁，与边界4的距离为h₂，与边界3的距离为m₁，矢量水听器R₁与边界1的间距为m₂；截面l₁中水池边界的长宽长度分别为a＝h₁+h₂与b＝m₁+m₂+R₁S₁·cosθ；另外，根据声线在六个平面的反射次数的不同，分别将声线与平面Ⅰ(1)、平面Ⅱ(2)、平面Ⅲ(3)、平面Ⅳ(4)、平面Ⅴ(5)和平面Ⅵ(6)的反射次数用参数n₁、n₂、n₃、n₄、n₅以及n₆来表示；根据每条声线与不同边界面的反射次数次数的不同，由勾股定理得到不同的声线的声程为：

直达声线R₁S₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝0，R₁S₁＝c·τ₁＝L₁ (9)

R₁S₁₁，即n₁＝0，n₂＝0，n₃＝0，n₄＝1

R₁₁＝(L₁cosθ)²+(L₁sinθ+2h₂)²＝(L₁cosθ)²+(L₁sinθ+n₄·2h₂)² (10)

R₁S₁₂，即n₁＝0，n₂＝1，n₃＝0，n₄＝0

R₁₂＝(L₁cosθ)²+(2h₁-L₁sinθ)²＝(L₁cosθ)²+(n₂·2h₁-L₁sinθ)² (11)

R₁S₁₃，即n₁＝0，n₂＝0，n₃＝1，n₄＝0

R₁₃＝(L₁cosθ+2m₁)²+(L₁sinθ)²＝(L₁cosθ+n₃·2m₁)²+(L₁sinθ)² (12)

R₁S₁₄，即n₁＝1，n₂＝0，n₃＝0，n₄＝0

R₁₄＝(L₁cosθ+2m₂)²+(L₁sinθ)²＝(L₁cosθ+n₁·2m₂)²+(L₁sinθ)² (13)

R₁S₁₅，即n₁＝0，n₂＝1，n₃＝0，n₄＝1

R₁₅＝(L₁cosθ)²+(L₁sinθ+2h₁+2h₂)²＝(L₁cosθ)²+(L₁sinθ+n₂·2h₁+n₄·2h₂)² (14)

R₁S₁₆，即n₁＝0，n₂＝1，n₃＝1，n₄＝0

R₁₆＝(L₁cosθ+2m₁)²+(2h₁-L₁sinθ)²＝(L₁cosθ+n₃·2m₁)²+(n₂·2h₁-L₁sinθ)² (15)

R₁S₁₇，即n₁＝0，n₂＝1，n₃＝1，n₄＝1

R₁S₁₈，即n₁＝1，n₂＝0，n₃＝1，n₄＝1

R₁S₁₉，即n₁＝1，n₂＝1，n₃＝1，n₄＝2

6.根据权利要求5所述一种混响水池中三维声场的数值计算方法，其特征在于，结合虚源法所列出的公式(1)-(8)，推测在平面l₁内，混响水池中虚源与接收点距离的表达式P_n如下：

P_n＝(2n₃m₁+2n₁m₂+L₁cosθ)²+(2n₂h₁+2n₄h₂±L₁sinθ)²

其中L₁sinθ项的系数为正或为负由最终反射界面决定。

7.根据权利要求5所述一种混响水池中三维声场的数值计算方法，其特征在于，根据平面l₁的声场情况，可得出平面t₁中虚源与接收点距离的表达式Q_n如下：

Q_n＝(2n₃m₃+2n₁m₄+L₁cosα)²+(2n₅h₃+2n₆h₄±L₁sinα)²

对接收点的每条声线与直达声线轴所形成的空间角进行积分的方式，来达到对整个空间的虚源的辐射声压进行叠加的目的。

8.根据权利要求7所述一种混响水池中三维声场的数值计算方法，其特征在于，所述对接收点的每条声线与直达声线轴所形成的空间角进行积分具体包括以下步骤：

步骤4.1：设R₂为接收水听器，S₂为声源，两者位于水平面l₂与垂直平面t₂的相交直线上；S’₂与S”₂为声源关于边界的镜像虚源；

步骤4.2：BC与直线R₂S₂垂直，R₂A与直线S’₂S₂垂直，R₂S₂与直线S”₂S₂垂直，DB垂直平面l₂；R₂S’₂为平面l₂上声线的声程，R₂S”₂为平面t₂上反射声线的声程；在空间上存在多个声场平面，而每个平面上，所有经过多次声线反射所产生的镜像虚源均位于平面r₂上；

步骤4.3：基于步骤4.1与步骤4.2的定义可知，若令∠AR₂S₂＝θ₀，由于R₂A与直线S’₂S₂垂直且BC与直线R₂S₂垂直，因此∠CBS₂＝θ₀，同样的，设声源R₂与平面Ⅱ(2)和平面Ⅳ(4)的距离分别为h₁与h₂，与平面Ⅴ(5)和平面Ⅵ(6)的距离分别为h₃与h₄，R₂S₂＝d₀，R₂S’₂＝d₁，R₂S”₂＝d₂；则由几何关系，得到如下等式：

上式看作虚源与接收点距离通式的变形；

步骤4.4：若平面l₂以直线R₂S₂为轴进行旋转，则B点的运动轨迹为沿着BD方向的直线，BC即为DC在平面l₂上的映射，BS₂即为DS₂在平面l₂上的映射，而∠DBC则为旋转后的平面l’₂与原平面l₂的夹角，若令∠DBC＝β，∠CDS₂＝θ，DS₂＝h，则根据几何关系有：

步骤4.5：对一次反射声线而言，基于公式(20)和(21)等式在任何声场平面上均成立；因此若以h代替声源与边界2的距离h₁，并将等式(20)与(21)代入等式(19)，得到虚源与接收点距离关于角度θ的等式如下：

步骤4.6：在考虑反射系数的情况下，单个虚源在接收点的辐射声压为：

Ω(φ)＝kdsinφsinε Ω_z(φ)＝kdcosφcosε

其中J₀(·)为零阶第一类贝塞尔函数，V(θ)为液体，弹性界面的反射系数计算公式；k为波束矢量，k＝ω/c，ω为角频率，c为介质声速；Z₁、Z₂为反射界面两侧的声阻抗，

为声线入射角，

为折射角；φ为声线传播矢量与空间直角坐标系的z轴的夹角，ε为声线传播矢量在xOy面的投影与x轴的夹角，d为虚源到达接收点的距离；

步骤4.7：基于步骤4.6的单个虚源在接收点的辐射声压为只是一个声场平面的一个虚源在接收点的辐射声压，若要得到整个声场的虚源在接收点的声压总和，将所有由通式算得的声线长度，代入到虚源辐射声压的计算公式，并对声场平面绕直达声线旋转的旋转角β进行积分，即

其中p_n(d_n)为单个虚源辐射声压，θ₁与θ₂为积分角θ的积分范围；

步骤4.8：将步骤4.6的公式中cosθ与sinθ做如下替换：

而由平面l₁的旋转，则确定角度β的积分范围分别为

和

对三个角度范围分别进行积分，便得到整个混响场中的虚源在接收点R处的声压；

算出空间中每个点的声场声压，从而构建出三维声场。

Images