本发明航天器力学环境预示、力学环境条件设计及声振响应预示等研究领域，特别是一种航天器中低频声振响应预示方法。
背景技术：
在航天任务过程中，航天器是否能够经受住发射段残酷的力学环境是决定任务成败的关键因素之一，比较典型的包括级间分离载荷、气动载荷、火工冲击载荷等。上述载荷主要通过两种途径传递到航天器结构或者关键部件上：一是通过星箭界面或者安装界面传递，这种力学环境一般通过正弦、随机振动和冲击试验进行考核；二是通过整流罩内的噪声环境直接作用在航天器结构与设备上，一般通过噪声试验进行考核。航天器经历的噪声环境十分残酷，一方面中高频段的噪声环境会造成机电设备、光学仪器、电子电路等的失效和精度下降，另一方面，在实际工程中，中低频噪声环境也会偶尔出现太阳翼、天线和馈源等轻质薄璧结构的破坏案例。
目前航天工程中主要针对不同的频段采用不同的分析方法进行声振力学环境预示：在低频段采用结构有限元-声学边界元方法，此方法中结构采用有限元方法建模，航天器的内外声场采用直接边界元方法建模，然而直接边界在建立内外部声场模型时需要建立封闭的边界，但是对于航天器结构而言，由于其构型复杂，太阳翼、天线等结构建立封闭声腔的难度很大，因此直接边界元法的应用受限制较多；另外噪声载荷模型也致关重要，目前比较流行的是采用扩散声场载荷建模方法(daf)，该模型认为噪声空间为完全混响空间，然后基于直接场与扩散声场的互易关系，将声压转化为结构网格节点上的力谱，然后直接施加在结构上求解，该模型主要特点是建模方法简单，能够描述声场载荷的空间相关性，但是其空间相关性仅在施加的区域内存在，不同区域之间的声载荷无相关性，无法刻画混响声空间中由于结构存在导致的声场不均匀性；另外该模型主要施加在结构上，无法考虑结构的声辐射、散射特征，需要结合半无限声场系统模拟结构的声辐射等特征。
技术实现要素：
本发明解决的技术问题是：克服现有技术的不足，提供了一种航天器中低频声振响应预示方法，针对航天器中低频声振响应预示的问题给出了有限元-间接边界元的混合建模方法，并基于平面波叠加理论，提出能够直接施加在声学间接边界元模型上的混响声场载荷模型，并给出了确定平面波数目的方法，可为后续的工程应用提供借鉴。
本发明的技术解决方案是：一种航天器中低频声振响应预示方法，包括如下步骤：
(1)建立包含平面波载荷的声学间接边界元模型；
(2)建立航天器结构有限元模型；
(3)建立基于平面波叠加建立混响声场载荷模型；
(4)建立航天器中低频声振响应混合模型；
(5)根据航天器中低频声振响应混合模型进行中低频声振响应预示。
所述的建立包含平面波载荷的声学间接边界元模型为：
其中，ρa为声场空气密度，ω为角频率，u为位移向量，μ和δμ分别为声场边界上双层势函数向量及双层势函数变分向量，cμu为双层势函数与位移间的耦合矩阵，dμμ为双层势函数的自相关矩阵，finc为由平面波载荷产生的载荷，<>为行向量，{}为列向量。
所述的建立航天器结构有限元模型为：
<δu>([k]{u}-ω2[m]{u}+[cuμ]{μ}-{f})＝0
其中，u和δu分别表示离散点上对应的航天器结构的位移值和位移变分值，k为航天器结构的刚度矩阵，m为航天器结构的质量矩阵，cuμ为航天器结构位移与表面双层势函数的耦合矩阵，f为航天器结构的外载荷矩阵。
所述的基于平面波叠加建立混响声场载荷模型的建立方法为：
(1)将整个声空间建立在球坐标系中，将球面按照经度和纬度进行切分，在纬度上分为m层，则每层之间间隔角度为
(2)假设经度剖分数目n和纬度剖分数目m相等，则
p00＝8.259；
p10＝-0.01022；
p01＝-1.456；
p20＝1.695e-05；
p11＝0.009679；
p02＝0.3041；
p30＝-1.024e-08；
p21＝3.783e-07；
p12＝-4.932e-05；
p03＝-0.01996
其中f为圆频率，c为声场的声速，r0为航天器结构的最大包括尺寸。
所述的建立航天器中低频声振响应混合模型为
所述的根据航天器中低频声振响应混合模型进行中低频声振响应预示的方法为：
其中，p(x)为声场内非边界上空间内位置x处的声压，
一种计算机可读存储介质，所述计算机可读存储介质存储有计算机程序，，所述计算机程序被处理器执行时实现如权利要求1-权利要求6任意所述方法的步骤。
本发明与现有技术相比的优点在于：
本发明针对现有技术的不足，提出结构有限元-声学间接边界元的中低频声振响应预示方法，便于对具有复杂边界的航天器结构建立声学模型；另外针对声学间接边界元，提出了基于平面波叠加的混响声场载荷建模方法，并给出了确定平面波数目的方法，便于工程应用。
附图说明
图1为流程框图；
图2为结构域的声固耦合示意图；
图3为平面位置目及其有效辐射面积示意图；
图4为平面波数目及其拟合结果
具体实施方式
本发明针对航天器中低频声振响应预示这一问题，提出了一种航天器结构有限元-声学间接边界元的混合建模方法，在此基础上提出结合基于平面波叠加的方法模拟混响声场载荷，并给出了确定平面波数目的方法，可为后续的航天工程应用提供借鉴，进而完成航天器中低频声响响应预示。
本发明的具体实施步骤如图1所示，主要包括如下步骤：
(1)建立包含平面波载荷的声学间接边界元模型
由于大多数航天器结构并不具备封闭区域，采用直接边界元建立模型非常困难，这里采用间接积分格式建立声场模型。对于航天器结构的声固耦合问题，声场的边界为纽曼边界，所以采用双层势，可通过双极子源对声场中位置x处的声压p进行重构：
其中p(x)为x点的声压，
如果考虑入射平面波激励，在非边界上有：
其中pinc(x)为平面在位置x处的声压。
在边界
其中上标“+”代表的是边界上法向的正向边界，上标“-”代表的是边界上法向的负向边界，nx为边界x处的法向。当结构与声场耦合时，耦合处满足法向的位移协调条件：
其中ρa为声场的密度，ω为圆频率，un(x)为x处的法向位移。
对式(3)在nx方向上求微分，并带入(4)式有：
式(5)的右边积分项为超奇异积分，无法直接积分，因此对式(5)两边乘以变分δμ(x)并在整个边界上以变量x进行积分有：
其中dsx为位置x处的面积微元。
上式右边的积分为正则积分，经有限元离散后进行积分可得到如下矩阵形式的表达式：
其中<>表示行向量，{}表示列向量，dμμ为双层势函数的自相关矩阵，可通过直接在一系列离散单元上计算式(7)的积分得到，μ和δμ分别表示离散点上对应的双层势函数的值及其变分值。
其余积分可写为：
其中cμu为双层势函数μ与位移u的耦合矩阵，finc为平面波引起的等效载荷向量。综合式(7)、(8)和(9)有：
(2)建立航天器结构有限元模型
考虑如图2所示的系统，如图3所示为平面位置目及其有效辐射面积示意图，其中u为结构的位移场，
采用迦辽金弱积分格式可建立结构的动力学变分格式：
其中“:”为张量双点积，“·”为点积。有限元离散并积分后有：
其中δu和u分别表示离散点上对应的位移变分值和位移值，k为结构的刚度矩阵，m为质量矩阵，cup为位移与声压的耦合矩阵，f为结构的外载荷矩阵。
因此结构的动力学响应可写为：
<δu>([k]{u}-ω2[m]{u}+[cup]{p}-{f})＝0(12)
上式中考虑了结构位移与声压的耦合关系，可实现结构与直接边界元的耦合，但是在间接边界元模型中，自变量为双层势函数μ，因此必建立结构周围声压p和双层势函数μ的关系。为了得到该关系，假设结构为薄壁结构，而且两边均与声场接触，在这种情况下，双层势函数μ可以描述为声压跨过
假设边界
采用单元离散后进行积分有：
其中cuμ为位移u与双层势函数μ的耦合矩阵。
则结构的动力学方程可写为：
<δu>([k]{u}-ω2[m]{u}+[cuμ]{μ}-{f})＝0(15)
(2)基于平面波叠加建立混响声场载荷模型
理想混响声场在各个方传的声能传播概率相同，空间内声能密度处处相等。理论分析表明，在各个平面波相位随机的情况下，平面波的声能密度叠加式满足线性叠加原理的。
假设将整个声空间建立在球坐标系中，将球面按照经度和纬度进行切分，在纬度上分为m层，则每层之间间隔角度为
rij＝[sin(θi)cos(ψj),sin(θi)sin(ψj),cos(θi)](16)
平面波将球体分为m×n份，平面波pij作用在球面上的面积假设为sij，在纬度方向，面积元的球面宽度为δθ。因此平面波pij在球面上的作用面积近似为：
δsij≈sin(θi)δψδθ(17)
为了保证混响场在球面上的能量均匀分布，有：
其中c为常数。
对于具有均匀分布随机相位的平面波，满足线性叠加关系，即：
其中p为声场的混响载荷声压。
这样就可以基于(19)式确定每个平面波的幅值：
(3)依据结构尺寸和分析频率确定平面波数目
当平面波数目趋向于+∞时，平面波叠加得到的混响声场载荷的空间相关性与理论空间相关性完全一致：
其中k为声场的波数，r为空间中任意两点的距离,γref为该两点位置之间的声场载荷的空间相关系数。在中低频，声场载荷的空间相关性对响应分析结果的影响非常显著，对于实际工程而言，给出如何确定平面波的数目具有很重要的工程应用价值。
平面波的数目与分析频率的上限和航天器结构的特征尺寸(最大包络尺寸)密切相关，以分析频率上限和结构最大包络尺寸(r0)为参数，以空间理论相关性为目标函数，以经度与纬度的平面波数(二者相等)为自变量，确定使得目标函数满足给定误差时需要的平面波数目，由于容差函数是一个具有多个零点的衰减振荡函数，在零点处无法定义相对误差。考虑到空间相关性最大值为1，并随着分析频率与空间距离的增加逼近于0，因此取绝对误差，这里仅给出绝对误差为0.01的结果：
|γ(k,r,m,n)-γref(k,r)|≤0.01(22)
其中γ(k,r,m,n)为声场中平面波经度、纬度划分数目分别为m和n时，距离为r的空间中两点在波数为k时的空间相关性。图4为容差为0.01时平面波数目与最大包络尺寸和分析频率上限的关系，为了便于工程应用，这里假设经度的剖分与纬度的剖分数目相同，采用二次曲面进行拟合，平面波数目与频率和最大包络尺寸满足如下关系：
其中f为圆频率，c为声场的声速，r0为航天器结构的最大包括尺寸。
(4)依据(1)、(2)和(4)步建立混合模型
在第(1)、(2)和(4)步建立好间接边界模型、结构模型和混响载荷模型后，组装形成混合声振响应分析模型。对于间接边界元，式(9)对于任意的δμ成立，则有：
ρ0ω2[cμu(ω)]{u}-[dμμ(ω)]{μ}＝{f'inc}(25)
而结构的响应对于任意的δu满足式(15)，因此有：
[k]{u}-ω2[m]{u}+[cuμ]{μ}＝{f}(26)
两者耦合可写为：
(5)求解得到结构响应与声场响应
基于式(27)可求解得到结构的位移响应{u}以及间接边界元的双层势{μ}，然后依据式(2)和式(3)可求解得到声场内任意一点的声压值p
本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。