基于malab的fe-aa混合方法的研究

基于malab的fe-aa混合方法的研究

基于波的fe-ea混合方法在高频振动和噪声负荷的作用下，每种复杂结构通常表现出相对复杂的动态特性。究其原因主要是由于各类不确定性引起的,如复杂结构高频参数的不确定性、结构制造工艺误差、振动载荷及声场的不确定性等。根据结构特征尺寸与结构内波长关系的不同使得结构不确定性对结构响应的影响也有所不同[1~3]。结构特征尺寸小于或相当于结构内波长时,即结构的模态稀疏时,结构不确定性对响应的影响很小,甚至低频时可忽略不计;反之,结构特征尺寸大于结构内波长时,即结构的模态密集时,不确定性对响应的影响较大,在高频段结构响应往往呈现出随机的特性。传统的有限元方法(FiniteElementMethod,FEM)能很好地反映确定性对结构响应的影响而无法反映不确定性的影响,因此该方法仅适用于低频段响应分析;而传统的统计能量分析(StatisticalEnergyAnalysis,SEA)仅能反映不确定性对结构的影响,因此该方法仅适用于高频段。对于复杂结构,各部组件特征尺寸与其结构波长的关系不同,如运载火箭上的桁条蒙皮结构,桁条的结构波长与其特征尺寸相当,而蒙皮的结构波长远小于其特征尺寸,此时各部组件的不确定性对整体结构响应影响不一致,有的呈现出低频特征,而有的呈现出高频特征。单独使用FEM或SEA均不能有效地对此类复杂结构的响应进行预示。对此类介于低频FEM和高频SEA的问题统称为中频问题。FE-SEA混合方法正是针对中频问题的特点提出的一种方法。该方法可根据子结构特征尺寸与结构波长关系的不同选取不同方法进行建模,如子结构特征尺寸小于或相当于其结构内波长时,应用FEM或边界元方法进行建模,称其为确定性子结构(或确定性子系统);而子结构特征尺寸大于结构内波长时,应用SEA进行建模,并称其为随机子结构(或随机子系统)。建模完成以后,通过确定性子结构与随机子结构连接处直接场和混响场之间的互易原理将子结构重新连接,得到整体结构的响应。FE-SEA既是传统FEM与SEA结合体,同时也是对传统方法的改进,填补了传统分析方法在中频段应用的空白。而且,当频率从中频逐渐向低频和高频两极趋近时,FE-SEA逐渐退化为FEM和SEA,也可以说FEM和SEA是FE-SEA的两个特例。1999年,Langley和Bremner提出了FE-SEA混合方法的基本理论。在应用模态叠加法的同时,将振型及相应的主坐标分解为整体模态集和局部模态集,并对与两个模态集相关的子系统分别应用有限元和统计能量分析进行建模。2005年,Langley和Shorter等在上述基于模态的FE-SEA混合法的基础上,提出了基于波动理论的FE-SEA混合方法。在各类混合方法中,基于波动理论的FE-SEA混合法可以处理各类结构间载荷双向传递问题,克服其他混合方法只能考虑能量单向传递的不足,使得该方法可应用于更为复杂的结构形式。同年,法国ESI集团也以基于波的FE-SEA混合法为理论基础推出了PAM-VAOne商业软件。此后,Cotoni和Shorter分别利用数值方法和试验方法对FE-SEA方法进行了验证。目前,伴随VAOne商业软件的推广,基于波的FE-SEA混合法逐步应用到轿车随机振动及噪声、飞机舱内噪声和火箭整流罩内噪声分析及航天器声振响应分析上[10~12],展现出良好的工程应用前景。本文首先阐述了基于波动理论的FE-SEA混合方法的基本理论,建立了该方法完整的分析流程。而后,应用Matlab开发了相应的计算程序并进行了校验。最后,结合MonteCarlo仿真分析和基于有限元的能量流分析,通过实例对基于波的FE-SEA混合方法进行了仿真验证。1基于波动理论的fe-ea混合方法1.1各子系统间边界划分考虑整体系统中的一个子系统,如图1所示,其区域为K,边界为Γ。若稳态谐波激励加载在该子系统边界区域上,则其边界上的时域位移响应可表示为ub(t)=Re{ubexp(ikt)},其中ub为与广义边界自由度对应的幅值,通常为复数向量。按照边界元法的思想,子系统边界区域的响应可由格林函数积分求得。为简化计算,可将边界离散为一系列广义坐标qb,此时边界位移可由式(1)求得,即式中ub(x∈Γ)为边界x处的位移,Okb为边界上的基函数,qkb为第k个边界广义坐标。只考虑子系统间的相互作用关系,忽略外载荷对边界的作用,子系统边界的动力学方程可表示为式中qb为子系统的边界广义坐标,fb为边界广义力(只表示其他子结构给予的反作用力)。若子系统中,一部分边界上所有的物理性质已知,而另一部分边界物理性质只是部分已知或完全无法确定,如图2所示,则可将子系统的边界划分为确定性边界Γd和随机边界Γr。子系统均可通过确定性边界Γd和随机边界Γr与其他子系统相连并通过这两类边界进行能量交换,确定性边界Γd和随机边界Γr也均可承受外部载荷。此外,由图2可知,确定性边界Γd无需满足边界连续的条件,任意子系统间的连接均可应用在这两类边界上。子系统边界进行能量交换的实质就是通过边界的位移场将能量传递到其他子系统上或接受其他子系统传递给自身的能量。对应子系统边界的分类,可将边界的位移场分为两类:直接场(DirectField)与混响场(ReverberantField)。直接场仅满足确定性边界上的边界条件,并不考虑随机边界上的边界条件。通过直接场可精确表示出确定性边界上的输出位移场,即通过直接场子系统可向外辐射能量;但通过直接场表示的随机边界上的位移场一般很复杂,暂且记为O,如图3所示。为满足随机边界上的边界条件,定义第二个位移场,即混响场。混响场的作用就是将混响场中的边界条件与直接场中的边界条件线性叠加后使确定性边界条件和随机边界条件同时得到满足。因此,混响场必须满足两个条件:1)混响场中确定性边界上的位移为0;2)与直接场线性叠加后随机边界条件得到满足,例如,若随机边界为固支,则混响场中确定性边界上的位移为0,随机边界上的位移为-O。对应于边界分类及直接场和混响场的定义,可将边界广义坐标qb分离为qdb和qrb,相应地将边界广义力fb分离为fdb和fbrev,则子系统边界的动力学方式可分解为由于随机边界上部分信息未知,因此有必要对式(3)进行化简。由式(3)的方程(1)可得式中Dbdir定义为直接场的动刚度阵并用于求解确定性边界上的位移,其表达式为对于每一项均为实数的基函数,矩阵Dbdir为复数对称矩阵。Dbdir通常可由边界元法求得,但对于某些规则的连接方式,如点、线和面连接,与此类规则的连接方式相关的直接场的动刚度阵Dbdir可由理论公式直接推导出解析解,即利用边界连接处的位移协调关系及波在结构中传播的性质推导出解析解[14~18]。1.2混响场建模与约束在基于波动理论的FE-SEA混合方法中,一个复杂系统通常可分为若干子系统。据子系统的特征尺寸与其系统中波长的关系可将子系统分为两类:若其特征尺寸与其系统中的波长相当,即该子系统的刚度很大、模态稀疏,系统中的一些不确定因素,如制造公差等,对其响应不产生影响,则该子系统为确定性子系统,可用FEM建模;若其特征尺寸大于其系统中的波长,即该子系统的柔性很大、模态密集,其响应对系统中不确定因素的变化很敏感,如相同设计同一批次的两个结构,由于制造公差等因素使两个结构的响应变化很大,则此类子系统为随机子系统,可用SEA建模。同时,对应地在数值上对各系统做离散化处理,系统响应可由一系列位移广义坐标q表示,其中确定性子系统的广义坐标为qd。类比有限元频响分析,经傅里叶变换FE-SEA混合方法中系统的动力学方程为式中fext为确定性子系统所受外力的向量;frev(m)为混响场受挡力;m为与确定性子系统相连的随机子系统的个数;Dtot为系统的总动刚度阵,是确定性子系统的动刚度阵Dd与直接场动刚度阵Ddir(m)的线性叠加Dd可直接由有限元求得,而Ddir(m)可由相应的Dbdir进行坐标变换求得。在式(6)中,直接场动刚度阵Ddir(m)表示确定性子系统通过直接场向随机子系统传递的能量以及随机子系统间通过确定性子系统进行的间接能量传输,frev(m)表示随机子系统通过混响场向确定性子系统传递的能量。由于结构不确定性的存在,类比统计能量分析中求集合平均的概念,将式(6)改写为互谱的形式并求集合平均得其中符号·-H表示矩阵的共轭转置并求逆的运算。假设系统具有最大熵(MaximumEntropy),即系统的平均值为所有可能出现样本的平均,简言之系统具有最大平均信息量。此时,混响场表示所有混响场样本的集合平均,随机子系统在混响场中的受挡力frev(m)与系统不确定性因素的变化无关,因此与frev(m)相关的集合平均趋近于下面的极限值式中Tm为与混响场振幅相关的比例常数,其求解如下式所示式中Em和nm分别为第m个随机子系统在混响场中的所具有的能量和第m个随机子系统的模态密度。通过子系统间的功率平衡关系可求解各随机子系统的能量Em。若系统中包含p+q个随机子系统,其中p个随机子系统与确定性子系统相连,其余的q个子系统只与随机子系统相连而与确定性子系统不相连,则随机子系统间的功率平衡方程为式中htot,m表示第m个子系统在混响场中的能量损耗系数,即随机子系统通过混响场向确定性子系统输出能量;hmn为随机子系统m与n通过确定性子系统建立耦合关系的耦合损耗因子;Zmk为能量从第m个随机子系统直接传递到第k个随机子系统时的耦合损耗因子;Zm为第m个随机子系统的内损耗因子;Pin(m,0)表示确定性子系统通过直接场传递到第m个随机子系统上的功率;Pin(m,1)表示外界直接对第m个随机子系统的输入功率;Em和En分别为第m个随机子系统和第n个随机子系统的能量;nm和nn分别为第m个随机子系统和第n个随机子系统的模态密度。htot,m,hmn及Pin(m,0)分别为1.3随机子系统相关参数分析综上所述,对一个计算频段基于波动理论的FE-SEA混合方法的分析过程可大体分为7步(如图4所示)。(1)系统的划分首先,将系统划分为确定性子系统和随机子系统;而后,划分子系统的确定性边界和随机边界并确定边界上的连接处及连接方式;最后,选取确定性广义坐标qd。(2)边界连接处的动刚度阵的获取规则点、线和面连接处的直接场的动刚度阵Dbdir可直接得到解析解。据确定性广义坐标qd的选取,动刚度阵Dbdir需进行相关的坐标变换。(3)与确定性子系统相关的参数分析应用FEM计算与确定性子系统相关的动刚度阵Dd,据选取的确定性广义坐标qd,相应地对动刚度阵Dd进行坐标变换;然后,将动刚度阵Dd与Ddir(m)线性叠加求得系统的总动刚度阵Dtot。同时,将确定性子系统的外载荷fext变换为互功率谱Sffext并进行相应的坐标变换。(4)与随机子系统相关的参数分析与随机子系统相关的参数分析包括:计算确定性子系统对与其相连的随机子系统的输入功率;计算各类随机子系统间的耦合损耗因子;计算随机子系统在混响场中的能量损耗系数;计算所有随机子系统由于自身的内损耗而损失的功率;将直接作用在随机子系统上的外力向量变换为相应输入功率Pin(m,1)。(5)随机子系统响应的求解据步骤(4)中的所得参数建立随机子系统的功率平衡方程,求得各随机子系统能量的集合平均。(6)确定性子系统响应的求解据步骤(5)中的所得随机子系统的能量求得随机子系统通过混响场作用在确定性子系统上的力谱;然后据步骤(3)建立的动力学求解确定性广义坐标qd的集合平均值〈Sd,qq〉。(7)重复步骤(2)~(6),计算下一个频段响应的集合平均。2计算机器的开发根据上述的FE-SEA混合方法分析流程,以Matlab为平台开发了相应的计算程序。而后,通过实例将混合方法的结果与MonteCarlo仿真分析的结果相比较对基于波的FE-SEA混合方法进行仿真验证。2.1响应系数和响应曲线如图5所示的板梁组合结构,结构的具体参数如表1所示,具有16个圆孔的Z型梁分别通过两个点和3个点与平板1和2相连。结构中Z型梁的一端施加固支边界条件,在板1的板面垂直方向施加空间位置不相关的宽频随机分布式激励。然后分别应用FE-SEA混合方法和MonteCarlo分析计算出板1受激励时随机子系统的能量影响系数(EnergyInfluenceCoefficient,EIC)并进行比较。其中,计算频段取值为1Hz的等带宽,中心频率范围为1~1500Hz。两类方法的计算过程中均选用模态为广义坐标,模态截断位置为3500Hz,即选用3500Hz之前结构的所有模态进行计算。由于在3500Hz内,Z型梁仅有28阶模态,模态非常稀疏,因此将梁建为确定性子系统;当中心频率为200Hz时,如表2所示,由板1的波数推导出板1弯曲、拉伸和剪切方向的波长分别为0.222,27.165和15.198m,与板1的对角线长度相比,弯曲方向的波长约为对角线的一半,而拉伸和剪切方向远大于对角线长度。对于该结构中的薄板,当某一方向的波长小于对角线长度的2/3时,该方向的位移可见为随机子系统。因此,可将板1面外弯曲建为随机子系统,而面内的拉伸和剪切方向的位移建为确定性子系统。同理,对板2各方向的位移进行处理。综合考虑该组合结构,则可将两块板弯曲方向的位移分别建为两个随机子系统1和2,面内位移与Z型梁组合建为确定性子系统。2.2不同集中质量对模拟结构的影响FE-SEA混合方法的结果是集合平均的结果,因此应用MonteCarlo仿真分别对300个样本进行分析,然后对这300次的结果进行平均得到结构的随机平均值,并与FE-SEA的结果比较来验证FE-SEA方法。对于任意单次样本,作用在板1上的输入功率P1和两个随机子系统的能量均是通过基于有限元的能量流分析求得。在MonteCarlo仿真过程中,随着结构中随机样本数的不断增加,样本平均结果的稳定性越来越强,即当样本具有最大包络时仿真中的集合平均结果与结构的不确定性无关,集合平均结果唯一确定。由于结构中的不确定性的影响主要是由随机子系统的不确定性引起的,因此在增加样本的包络方面主要是模拟随机子系统不确定性的变化。为模拟结构样本的最大包络选择在随机子系统的随机位置上添加集中质量。集中质量主要作用是增大随机子系统低阶模态对外界扰动的敏感性从而增大整个随机子系统对外界扰动的敏感性(中高频时,随机子系统本身对外界扰动已经很敏感)。建立MonteCarlo仿真中使用的FE模型(未添加集中质量),如图6所示,整个组合结构上应用5148个Tria3壳单元,总的节点自由度数为16866。对组合结构进行模态分析并提取出3500Hz以内结构的前342阶模态,图7所示为FE模型(未添加集中质量)的第100阶模态。考虑计算成本,仿真过程中选择FEM中的模态法并对模型进行300次MonteCarlo仿真。仿真的每个样本即对图6所示FE模型,分别在两平板上添加10个集中质量,所有集中质量的总和为整个组合结构质量的25%,所在位置为FE模型中两板的节点位置。令每块平板上集中质量所在节点的节点号的概率分布为均匀分布,应用Matlab分别生成每块板上的300组节点号随机数。据随机数每次生成集中质量的节点位置并据该节点位置相应地将集中质量添加到组合结构的FE模型中,然后对样本进行模态分析并提取出3500Hz以内样本的所有模态。在此基础上,对样本FE模型,在板1每个节点上(除点连接所在节点)施加垂直于板面方向的外力自谱,其中外力自谱在300次仿真中保持不变,而后利用模态法求得样本响应。最后利用基于有限元的能量流分析推导出两个随机子系统的能量E1和E2、输入功率P1及相应的EIC。300次MonteCarlo仿真时间约为5天,计算周期较长。图8所示为仿真过程中板2上的随机子系统EIC的两个单次样本。图9和10分别显示了300次MonteCarlo仿真的结果,从图中也可看出单次仿真结果间的差异明显。图9和10中,由300次仿真的平均结果可以看出在900Hz附近,即第100阶模态附近,两子系统间的能量传输明显,板1(激励板)的EIC下降即板1的能量向外传输时,板2(接收板)的EIC显著增大即板2接受了外部传入其子系统上的能量。2.3fe-ea混合法对FE-SEA混合方法,为得到确定性子系统的FE模型,将图7中FE模型上两块板的壳单元换为膜单元,梁上的壳单元不变。为了保证板上面内位移的质量阵,只考虑板的面内质量需单独建立膜单元的质量阵,因此在FE模型中可直接在节点相应的面内位移上添加集中质量自行建立膜单元的质量阵。求解之前通过Nastran输出确定性子系统的刚度阵和质量阵,并导入Matlab程序作为输入;同时还需将MonteCarlo仿真中所有的节点外力自谱通过点连接的动刚度阵将其转化为FE-SEA混合法中随机子系统的外部输入功率。在此基础上,应用FE-SEA混合法通过Matlab程序求得随机子系统的能量和相应的EIC。FE-SEA混合方法的计算过程中使用的确定性子系统的模态为3500Hz内的共28阶模态,整个的分析过程约40s左右,两项指标均远小于MonteCarlo仿真。如图9和10所示,由FE-SEA混合法所得的集合平均值与MonteCarlo仿真的样本平均值结果基本一致,尤其是中高频段,如900Hz附近,两类结果吻合得很好。但同时,低频段两类结果间存在误差,可能的原因包括:a)FE-SEA求解过程中多次求逆引入