基于空间谱单元的声学振动模拟预测方法、存储介质及设备与流程

1.本发明属于声学振动分析技术领域，尤其是涉及一种基于空间谱单元的声学振动模拟预测方法、存储介质及设备。


背景技术：

2.声学与振动的数值模拟和预测广泛应用于各大关系国计民生的工业领域中，比如建筑、桥梁、船舶、火车、航空航天、医疗器械等。目前该领域使用最广泛的数值模拟方法还是有限元法，对于某些高频情况使用统计能量法，对于某些大空间的辐射和传播使用边界元法。这些方法各有优势和缺陷。
3.有限元法理论简单，模型成熟，应用最为广泛，但是有限元要求每个波长内需有一定数量的结点，以保证结果的精度，这就使得有限元矩阵在计算大型模型高频情况时，会因为模型太大而需要极其强大的运算资源。即使是低频情况，如果模型较大，比如船舶工业上的整船，有限元就有些力不从心了。统计能量法作为替代方法可以在一定程度上弥补有限元方法在高频运算上的缺陷。但是，统计能量法主要基于声振传播理论，在理论上的简化很大，结果很难由实验验证，更多的是作为一个在统计层面上的指标。
4.另外，有限元不擅长处理大空间内的声辐射和传播，因为有限元需要将整个空间划分网格。对于这些工况，边界元方法可以使用较少的单元数达到与有限元类似的精度，从而在一定程度上替代有限元，但是由于边界元本身理论的原因，其系统矩阵往往不是对称方阵，格林函数在某些情况下不收敛或者耗时极长。另外，对于一些工况，比如折角和拐点，边界元不能很好的处理。但是边界元有它不可替代的优势，那就是只需要考虑边界，不管空间大小。这是有限元不能及的。所以，边界元作为弥补有限元在预测大空间内的声辐射和传播缺陷的替代方法更为合适。
5.再有，有限元和边界元单元内的值都是由结点值插值得到的，这就使得任意点值在一定程度上受到了插值函数的影响，从而引入误差。谱单元法也是解声学偏微分方程的一种方法，它不需要像有限元这样划分很细的单元，以保证每个波长内的单元数量，因此单元数量可以大大减少。谱单元法的未知量是谱的系数，因此在这些系数确定后，单元内任一点的值都是确定的，而非差值。这都是谱元强大的优势。谱元的缺陷是理论比较复杂，针对每个单元类型及其相关耦合开发出解析解费时费力。如果由数值解代替解析解，则运算量成指数倍增长。正是因为理论推导的复杂性，之前的谱单元都是正规的矩形、直角三角形等。虽然有人利用与有限元类似的形函数方法对单元进行改造，但是这样就失去了谱元不受波长内单元数的限制的优势。如果将谱元划分到与有限元单元相同的大小，由于其系统矩阵的复杂性，它就不如有限元的效率高了。也是因此，谱元在处理细微的结构细节时不如有限元。同时，适用于谱单元的相关声振辐射引发的附加质量和辐射效率的模型之前也没有被系统开发过，这些都限制了谱元在声学和振动预测领域的应用。


技术实现要素：

6.本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种能够突破谱单元形状限制、计算效率高、并能够直接在系统方程中考虑声辐射对结构振动影响，不需要后期处理或迭代的基于空间谱单元的声学振动模拟预测方法、存储介质及设备。
7.本发明的目的可以通过以下技术方案来实现：
8.一种基于空间谱单元的声学振动模拟预测方法，包括以下步骤：
9.对待预测对象进行单元划分，将大型结构划分为空间谱单元，其余细微结构划分为有限元，所述大型结构包括梁、板、壳或管道，基于单元划分结果对所述待预测对象进行建模；
10.根据待预测对象的实际情况设置边界条件、激励、响应点及各单元间的耦合，加入已构建的模型中；
11.计算当前模型中各谱单元的系数矩阵，并由系数矩阵组装得到系统振动方程中的系统矩阵；
12.根据实际需求对所述系统振动方程进行求解，获得求解结果，所述求解结果包括固有频率和系统响应；
13.基于所述求解结果对构建的模型进行可视化。
14.进一步地，设置各单元间的耦合时，以弹性构件作为耦合件，所述耦合包括点耦合、面耦合或线耦合。
15.进一步地，所述谱单元的系数矩阵的计算具体为：
16.基于利用傅里叶级数获取振动解；
17.基于所述振动解获取动能、势能和外力功；
18.利用rayleigh-ritz方法和hamilton方程，考虑系统稳定状态，推导获得系数矩阵。
19.进一步地，若单元的材料和特性随位置变化，则使用数值积分方法计算系数矩阵。
20.进一步地，若存在单元间的耦合，则耦合刚度的影响加入到势能项中。
21.进一步地，该方法还包括：
22.构建辐射模型，该辐射模型将每个单元划分为多个区域，每个区域看作是一个声源，根据该单元在该区域振动幅度的大小，推导出该声源对模型中某一点处声压的贡献，将所有声源的贡献进行叠加，则得到该点处的总声压；
23.在计算各单元的系数矩阵时考虑该单元上各点的所述总声压。
24.进一步地，在精度要求高、计算资源充裕的情况下可采用边界元作为上述辐射模型的替代和补充。
25.进一步地，所述声源为点声源或偶极子声源。
26.进一步地，在计算各单元的系数矩阵时考虑该单元上各点的所述总声压具体为：
27.将所述总声压视为一个系统附加的外力，获取该外力所产生的外力功，基于该外力功的叠加计算系数矩阵。
28.本发明还提供一种计算机可读存储介质，包括供电子设备的一个或多个处理器执行的一个或多个程序，所述一个或多个程序包括用于执行如上所述基于谱单元的声学振动模拟预测方法的指令。
29.本发明还提供一种电子设备，包括一个或多个处理器、存储器和被存储在存储器中的一个或多个程序，所述一个或多个程序包括用于执行如上所述基于谱单元的声学振动模拟预测方法的指令。
30.与现有技术相比，本发明具有以下有益效果：
31.1、本发明给出了梁、任意形状板(含加筋板，可有空洞)、任意形状壳(含加筋壳，可有空洞)、任意截面形状管道(含有流体填充管道)在内的多种谱单元以及相关耦合单元的解析解。将以往结构离散再耦合的过程融入到了获得系统矩阵的积分过程中，从而突破了之前谱单元对于单元形状的限制，减小了由不必要的单元耦合带来的误差和计算量，大大提高了谱单元的适用性和计算效率。用谱元替代有限元，在预测精度相同的情况下，矩阵自由度可以成倍减小，从而提高运算效率，大大提高模型频率上限。同时给出了与解析方法相对应的数值方法，以处理材料和单元特性随位置变化的问题。
32.2、本发明给出了可适用于谱单元的声振辐射模型以及由辐射阻抗所引发的附加质量和辐射效率问题的解决方法，并可直接在系统方程中考虑声辐射对结构振动影响，不需要后期处理或迭代，使得利用谱单元解决船舶工业中极为重要的附连水问题和火车工业中需要的轨道辐射效率问等问题成为可能，从而扩展了谱元法在相关重要领域中的应用。
33.3、本发明将谱元与有限元和边界元进行耦合，利用有限元处理结构小细节，利用边界元作为辐射模型的补充处理大空间内的声辐射和传播，使各个方法取长补短，从而能够在全频域模拟和预测声学振动以及相关的辐射问题。
附图说明
34.图1为本发明的流程示意图；
35.图2为实施例中钢弹簧浮置板轨道平面与剖面图；
36.图3为实施例中响应点垂向位移幅值曲线。
具体实施方式
37.下面结合附图和具体实施例对本发明进行详细说明。本实施例以本发明技术方案为前提进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。
38.如图1所示，本发明提供一种基于空间谱单元的声学振动模拟预测方法，包括以下步骤：
39.1)对待预测对象进行单元划分，将大型结构划分为谱单元，其余细微结构划分为有限元，所述大型结构包括梁、任意形状板(含加筋板，可有空洞)、任意形状壳(含加筋壳，可有空洞)、任意截面形状管道(含有流体填充管道)，基于单元划分结果对所述待预测对象进行建模；
40.2)判断是否需要辐射模型，若是，则根据问题对于声辐射和声场预测的要求，结合精度和速度方面的考量，选择辐射模型，然后执行下一步，若否，则直接执行下一步；
41.3)根据待预测对象的实际情况设置边界条件、激励、响应点及各单元间的耦合，加入已构建的模型中；
42.4)计算当前模型中各单元的系数矩阵，并由各单元的系数矩阵组装得到系统振动
方程中的系统矩阵；
43.5)根据问题需求求解系统振动方程，得到固有频率、系统响应以及其它辐射参数；
44.6)基于所述求解结果对构建的模型进行可视化。
45.上述方法的具体创新点如下：
46.1)给出了包括梁、任意形状板(含加筋板，可有空洞)、任意形状壳(含加筋壳，可有空洞)、任意截面形状管道(含有流体填充管道)在内的多种谱单元以及相关耦合单元的解析解，并给出了与解析方法相对应的数值方法，以处理材料和单元特性随位置变化的问题。
47.解析解由理论推导出来，为可以根据所给的参数直接计算出所需要的方程结果的公式。本方法中，各单元的解析解的推导思路为：给出单元振动解的傅里叶级数形式，由该解的形式推导出单元动能、势能、外力功的表达式；再将这些带入hamilton方程，利用rayleigh-ritz方法求得单元振动方程的矩阵形式。这个过程中涉及很多积分才能得到振动方程中的矩阵，需要对每一个积分求得解析解，才能得到最终系统矩阵的解析解形式。以往结构离散再耦合的过程也在这里被融入到了获得系统矩阵的积分过程中，从而突破了之前谱单元对于单元形状的限制，减小了由不必要的单元耦合带来的误差和计算量，大大提高了谱单元的适用性和计算效率。
48.此处以一维梁单元只考虑垂向挠度为例给出解析解推导的方法，其他类型谱单元以及其他自由度方向可遵从相同思路获得，不一一赘述。
49.利用傅里叶级数的概念给出一维梁挠度方向振动的解为：
[0050][0051]
其中w为挠度，m表示x方向傅里叶级数阶数，m比表示傅里叶级数展开所取阶数的最大值，a表示傅里叶级数的系数，p表示平滑辅助函数，l表示梁的长度。这里p的选取可以任意，但是需要保证解在边界上的高阶连续。
[0052]
利用解(1)可以给出单元的动能，势能，和外力功
[0053][0054][0055][0056]
其中，v为势能，t为动能，w为外力功，d为抗弯刚度，t为时间变量，ρ为密度，s为梁截面积，f为所考虑挠度方向的外力。
[0057]
利用rayleigh-ritz方法和hamilton方程，可知在系统稳定的状态下，需满足
[0058]
δ(t
–
v+w)＝0
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0059]
利用(1)-(5)，将(1)-(4)全部代入(5)中，可以整理推导出该单元的系统振动方程的矩阵形式为：
[0060]
(k-ω2m)a＝f
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0061]
其中，k、m、f分别为刚度矩阵、质量矩阵和力矩阵，ω是圆频率，a是系统傅里叶级数展开系数和平滑函数系数组成的向量，也是问题要求解的未知数。
[0062]
对于每一个梁单元，都可以由此得到该单元的系数矩阵。同理可以得到其他类型
单元的系数矩阵。这些单元包括但不限于，欧拉梁、铁木申科梁、任意形状板(含加筋板，可有空洞)、任意形状壳(含加筋壳，可有空洞)、任意截面形状中空管道、任意截面形状流体填充管道。这里如果该单元的材料和单元特性随位置变化，比如变截面梁和非均匀密度板，则公式(2)到(5)中的积分可使用数值积分方法以便考虑每一个积分点的具体参数特性。如果使用数值积分，则积分区间将被划分为多个小区间，认为每个小区间内参数值是一个常数，叠加所有小区间的贡献，可以得到最后的积分结果。
[0063]
各单元之间的可以通过类似弹簧的构件进行耦合，弹簧的刚度决定了耦合的强度。根据单元类型的不同，相应耦合可以是点、线、面的形式。耦合刚度的影响需要加入到方程中势能项中，以影响系统的总刚度矩阵。
[0064]
在所有单元的系数矩阵和耦合矩阵都得到后，利用每个单元的系数矩阵，可以组装得到整个系统振动方程的矩阵形式。总系统方程与(6)形式相同，只是其中的矩阵，已经考虑到了所有元素的影响。
[0065]
若不需要耦合有限元和边界元，则该方程可以直接求解。如果有有限元和边界元，则系统矩阵还需加入相关耦合项后再求解。
[0066]
2)给出基于谱单元的声振辐射模型
[0067]
如果声辐射对于系统振动的影响需要被考虑，则根据问题对于声辐射和声场预测的要求，结合精度和速度方面的考量，选择辐射模型，可选择的模型包括边界元或者基于谱单元的快速声振辐射模型。由于边界元会大大降低模型的计算效率，基于谱单元的辐射模型是系统默认模型，而边界元可作为该模型的替代和补充在精度要求高、计算资源充裕的情况下采用。无论是使用哪一种，这里声辐射对结构振动影响是可以直接在系统方程中考虑的，不需要后期处理或迭代。
[0068]
本发明设计了可以适用于谱单元的快速声振辐射模型。该模型的具体思路是将每个单元看作是一个声源，根据该单元上振动幅度的大小，推导出该声源对模型中某一点处声压的贡献。根据具体情况不同，该声源可以看作是无限大障板上的点声源，比如对板或壳类型的谱单元(7)，或者偶极子点声源，比如对梁或管类型的谱单元(8)。
[0069][0070][0071]
其中p为某点处声压，ρ0为媒质密度，w为声源处位移，r为某点到声源的距离，ω是圆频率，k是波长，ds是声源单元面积，l是偶极子声源间距，θ是偶极子声源与声压点连线的夹角。将所有声源的作用进行叠加，可得到该点处的总声压。
[0072]
对于有限元单元，这个过程是很简单的，因为结点位移就是有限元方程的未知变量。但是对于谱单元，由于其未知变量是位移在空间域谱的系数，该方法不能直接适用，必须要根据每类单元的定义分别推导。另外，因为谱单元往往几何尺寸比较大，为保证精度，需要将单元划分成多个辐射区域。如果还以前面提到的梁单元为例，具体步骤如下：
[0073]
1、将该单元划分为多个小的段，每段视为一个声源点，该声源点的位置在该端的中心，该声源点的声能是该段辐射能量的总和。
[0074]
2、对于每一个这样的声源，利用式(8)可以求得该声源在空间任意点处的声压贡献。
[0075]
3、针对结构上的某一点，将该系统中所有这样的声源对于该点的声压贡献加在一起，则可以得到在整个系统在振动的状态下，因为系统振动对于该点所产生的压力。
[0076]
4、将该压力视为一个系统附加的外力，因这个外力所产生的外力功，也可以由式(4)表达。将该外力功与其他外力功一起带入式(5)，可以导出考虑了声振辐射所产生影响(附加质量和辐射阻尼)后的单元振动方程。
[0077]
5、请注意，这里式(4)、(7)、(8)的未知量都是w，而w是由式(1)表达的，因此，最后的未知量是傅里叶级数展开系数和平滑函数系数，与该单元振动方程(6)中的未知量相同。如果该系统中有多个单元，则每两个单元之间以及该单元对自身的辐射影响都需要考虑。
[0078]
3)给出谱单元与有限元和边界元耦合的模型
[0079]
如果有限元和边界元需要与谱单元进行耦合，则可以将系统方程和相应方程进一步表示为。
[0080][0081]
其中k、m、a、f分别代表系统刚度矩阵、质量矩阵、变量向量、力矩阵。下标s、f、b分别代表谱元、有限元和边界元。这里k
couple
代表由耦合产生的耦合刚度项。谱单元与其他单元间的耦合是利用类似弹簧的构件来实现的，具体矩阵形式可以按照之前所述的求解谱单元的步骤求得。至于有限元和边界元系统矩阵的推导这里不赘述。
[0082]
在获得系统振动方程后，通过对振动方程的求解，能够得到固有频率、系统响应以及其它辐射参数，具体地：采用特征值法可以求出系统的固有频率；对系统方程进行求解，可求得系统傅里叶级数展开系数和平滑函数系数，将与某单元对应的傅里叶级数展开系数和平滑函数系数带回到该单元振动解的傅里叶级数形式中，比如式(1)，可以得到该单元任意点的响应；利用求得的单元的振动幅值，可以由公式(7)或(8)预测该单元在空间任意点上对声压的贡献，辐射阻抗和效率之类的参数都是这之后通过公式利用某点的系统总声压和结构的振动幅值信息求得。
[0083]
实施例
[0084]
本实施例将上述方法应用于悬置板火车轨道振动问题中，具体模拟预测过程包括：
[0085]
1)将轨道进行单元划分。钢轨选用t60轨，扣件采用dtvi2型，浮置板具体几何参数如图2，材料参数如表1。这里整块浮置板可划为一个矩形板谱单元，两跟铁轨可以划分为两个铁木申科梁谱单元。
[0086]
2)在相应点给出浮置板的弹簧支撑边界条件，保持铁轨边界自由。在左轨1.5米处设置激励点，给与1千牛的垂向激励。在左轨2.7米处设置响应点。
[0087]
表1 t60轨和dtvi2扣件参数
[0088][0089]
3)在梁单元和板单元之间利用点耦合模拟轨道与浮置板的扣件连接。
[0090]
4)组装系统矩阵并求解。可得到该结构前10阶的固有频率如表2，响应点0-1500赫垂向位移幅值如图3。
[0091]
表2浮置板前10阶固有频率
[0092]
阶数频率(赫兹)14.2283e-0524.6056e-0536.484847.134657.3455638.0354741.8748891.1985991.623010103.0473
[0093]
以上详细描述了本发明的较佳具体实施例。应当理解，本领域的普通技术人员无需创造性劳动就可以根据本发明的构思作出诸多修改和变化。因此，凡本技术领域中技术人员依本发明的构思在现有技术的基础上通过逻辑分析、推理或者有限的实验可以得到的技术方案，皆应在由权利要求书所确定的保护范围内。