改进的区间参数结构频响函数迭代解法

范芷若， 姜　东， 董萼良， 费庆国

(1.东南大学 土木工程学院 工程力学系, 南京　210096;2. 江苏省工程力学分析重点实验室, 南京　210096; 3. 南京林业大学 机械电子工程学院, 南京　210096)

改进的区间参数结构频响函数迭代解法

范芷若1, 2， 姜东3， 董萼良1, 2， 费庆国1, 2

(1.东南大学 土木工程学院 工程力学系, 南京210096;2. 江苏省工程力学分析重点实验室, 南京210096; 3. 南京林业大学 机械电子工程学院, 南京210096)

针对区间不确定性系统频响函数求解提出了一种改进的迭代方法。采用区间分析方法描述不确定性结构参数，将频响函数求解转化为区间线性方程组形式，使用定点迭代方法，构造迭代格式；为了解决迭代方法中出现的收敛性问题，结合子区间概念提出改进的迭代方法，计算结构频响函数的区间。以弹簧-质量系统为仿真算例，研究了子区间数目对计算结果的影响。最后，将方法应用于文献中的工程结构，对比研究表明：该方法适用于全频段频响函数分析、具有较高的精度和效率。

不确定性；区间分析；区间线性方程；子区间；频响函数

结构参数不确定性是工程中不可避免的问题。参数不确定性产生的原因很多[1]，结构物理量存在离散性，如杨氏模量、泊松比、密度等；在设计初始阶段因掌握资料不足产生的不可预测性，如边界条件和连接方式难以准确参数化。不确定性的存在对结构动态特性和结构动响应的影响不容忽视。

针对不确定性问题的概率方法研究较为广泛[2]，这类方法在掌握不确定参数的概率密度函数或均值、方差等统计量后开展研究，其中最有效的方法是蒙特卡洛模拟法[3]，但这种方法要求大量确定性计算，在实际工程结构的分析应用中存在难度。此外，工程中难以得到不确定参数的概率分布或统计量，概率方法并不适用于描述某些特殊的不确定性。为了克服概率方法的缺点，出现了一些非概率的方法，其中区间分析与其他方法在理论上有密切联系。Moore等[4]最早提出了区间分析的概念，用一个带上下边界的区间来定义参数的不确定性。区间有限元方法[5]是在考虑所有参数的不确定性后，分析结构响应的取值范围。这类方法因其所需信息少、计算初始条件简洁等特点，在不确定性结构响应分析中具有广阔的应用前景。Ben-Haim等[6]提出了不确定性结构力学分析的凸模型方法。Chen等[7]采用区间量描述系统矩阵的不确定性，并提出区间摄动法，以减少区间参数的相关性在区间运算中导致的区间扩张问题。邱志平等[8]将区间摄动法应用于不确定性结构振动问题特征值求解中。在摄动法基础上，李金平等[9]提出基于Neumann展开与摄动理论相结合的线性方程组的解法，提高了摄动法的准确性，但未解决摄动法在不确定性水平较高时的精度不足问题。苏静波等[10]将子区间方法与摄动法联合应用，提高不确定性水平较高时的区间运算精度，并通过计算得出精度较高时子区间的数目。Muhanna等[11]提出一种结构静力响应的区间线性方程解法，该算法是对线性方程形式的一种特定迭代解法。

频响函数求解在结构动力学响应计算中至关重要。针对于带区间参数结构频响函数的求解，Moens等[12]提出基于模态叠加法原理的模态矩形法，用于计算带区间参数结构频响函数，其结果精度高，但是不适用于全频段分析，且算法较复杂；Dessombz等[13]提出使用区间算法计算结构频响函数区间。

本文提出一种改进的不确定性结构动态特性频响函数的求解方法。将区间线性方程方法应用于带区间参数结构频响函数包络值的求解，结合子区间概念提出改进解法以提高计算精度，采用弹簧-质量系统为仿真算例验证方法的有效性，并将方法应用于工程结构分析。

1　区间频响函数

1.1区间参数

区间参数是实数域R内的一个封闭集合，对于给定的区间参数xI∈ IR，可以将它表示为

(1)

(2)

令

(3)

得到区间参数的另一种常见表达形式

xI=x0(1+αI)

(4)

xI○yI=[min(xi○yi),max(xi○yi)],

○∈{+,-,×,/}

(5)

1.2不确定线性系统频响函数

在不确定系统中，由于参数具有不确定性，将不确定性采用区间参数描述，则结构刚度矩阵KI∈IRn×n、质量矩阵MI∈IRn×n及阻尼矩阵CI∈IRn×n，IRn×n表示n×n维区间矩阵，代入多自由度线性系统的动力学方程

(-ω2MI+iωCI+KI)HI=I

(6)

式中：ω为圆频率，HI为区间频响函数，I为单位阵。

由此可将不确定线性系统频响函数求解转化为求解形如Ax=b的区间线性方程的解集

Σ(A,b)={x∈Rn|∃A∈A,∃b∈b:Ax=b}

(7)

其中

(8)

2　定点迭代解法

在区间线性方程的求解方面目前已有一些较成熟的方法[11]，例如区间高斯消元、区间Gauss-Seidel迭代、Krawczyk迭代和定点迭代等。本文采用定点迭代方法构造迭代格式求解区间线性方程。

2.1迭代格式

对方程进行如下处理，

R=inv(mid(A))=A0-1

(9)

解的估值为

x0=R×mid(b)

(10)

定义迭代初始条件

g=R×(b-A×x0)=-αRA1x0

(11)

G=[I]-R×A=-αRA1

(12)

迭代格式如下

(13)

设ε为迭代系数，对迭代格式做处理

(14)

将迭代格式进行展开，在进行了n次迭代后，解的格式可以表示为

(15)

(16)

由式(15)，(16)可看出，迭代最终结果中包含区间向量的乘方运算。区间运算规则决定了区间运算后会出现误差累积，因此在迭代次数过多时，该方法会出现过度估计的结果。需要尽量减少迭代次数，以提高算法精度。这就要求迭代矩阵G具有强收敛性。

2.2改进的迭代解法

上述解决区间线性方程的迭代算法，是基于定点迭代方法的推广，因此迭代矩阵G必须满足收敛性条件ρ(|G|)<1，才能得到算法的保守解。

由式(12)可知，不确定参数区间宽度较小时，迭代矩阵G收敛性较好，线性方程迭代解法可以获得较高精度，而对于参数区间宽度过大的系统，其迭代形式的收敛性就会降低，甚至出现不收敛的情况，迭代结果与真实结果相比误差相应增大。

为解决这一问题，结合子区间概念对迭代格式进行改进，提出一种改进的迭代解法，可以有效的提高迭代精度，解决迭代收敛问题。

i=1,2,…,L

(17)

则区间参数矩阵A相应分为L个

(18)

对于每个Ai，区间宽度由α变为αi，A0的值变为A=A0+mid(αi)A1，A1不变。相应对于每个子区间线性方程，其迭代格式变为

Ri=inv(mid(Ai))

(19)

x0=Ri×mid(b)

(20)

gi=Ri×(b-Ai×x0)=-αiRiA1x0

(21)

Gi=[I]-Ri×Ai=-αiRiA1

(22)

(23)

利用结合子区间参数矩阵修改后的迭代格式进行求解，待求解的区间线性方程改写为

Aixi=bx=∪xii=1,2,…,L

(24)

对于带多个区间参数的系统，可将迭代中各参数做相应推广(带N个区间参数)

(25)

3　算例研究

3.1弹簧质量系统

考虑如图1所示的带有阻尼的6自由度弹簧质量系统，结构阻尼系数η=0.06。激励点与激励方向如图1所示，确定性参数如表1所示。

表1　6自由度系统确定性参数

图1　6自由度弹簧质量系统Fig.1 Spring mass system of 6 DOFs

结构系统的部分刚度k(N/m),质量m(kg)具有误差或有界不确定性。设这些区间参数分别为

k1=[6 400, 9 600], m2=[0.009 6, 0.014 4]

对该系统计算质点1振动频响函数。图2、图3分别给出了频率从0～250 Hz时，区间迭代方法及其改进算法计算得到的结构频响函数幅值、实部和虚部的上下界，以及蒙特卡洛方法计算得到的频响函数上下界，其中改进算法中将每个区间参数分别划分为16个子区间，蒙特卡洛抽样样本为20 000，蒙特卡洛的收敛解可被视为准确解[15]。

对比图2(a)、图3(a)系统频响函数幅值图，在频响函数非峰值处，三种计算方法得到频响函数幅值区间结果相近；在频响函数峰值处，与蒙特卡洛方法对比，改进的迭代方法精度明显高于区间线性方程法。

若分别考虑系统频响函数的实部与虚部区间，如图2(b)、2(c)、3(b)、3(c)，在对应频响函数峰值的频率点处，区间线性方程法与蒙特卡洛方法计算结果差距极大，而改进的迭代方法结果为蒙特卡洛方法吻合度较高。

图2　区间线性方程法与蒙特卡洛方法结果对比Fig.2 Results comparison between Liner interval equation method and Monte-carlo method

图3　改进的迭代方法(16子区间)与蒙特卡洛方法对结果比Fig.3 Results comparison between modified iteration method and Monte-carlo method

使用子区间摄动法[10]对质点1振动频响函数区间进行求解作为参照。该将频响函数峰值点处四种方法数值区间提取得表2。

在f=20 Hz频率点处，随着子区间划分数目变化，改进的迭代方法结果与蒙特卡洛方法对比的误差如图4所示。

表2　四种方法在特定频率点处幅值区间

由图4可知，在子区间划分数目超过8以后，改进的区间线性方程方法计算所得结果与蒙特卡洛方法计算结果误差降低到5%以内，且随着子区间数目增多，计算结果收敛。

表3　计算时间比较

图4　20 Hz处误差随子区间划分数目变化趋势Fig.4 Change of deviation against the number of sub-intervals at 20 Hz

计算所需时间对比如表3，由时间对比可得改进的区间线性方程迭代解法运算时间远少于蒙特卡洛方法。

综上算例结果可得，区间线性方程解法的改进解法可以得到足够的精度，计算时间相较于蒙特卡洛方法大幅减少，可以有效地提高计算效率。

3.2工程结构

将区间线性方程方法应用于文献[12]中的工程结构，进行频响函数区间包络值求解。图5为相应的确定性结构，为一个27自由度的二维梁模型。模型的各参数的属性在表3中列出。其中节点3的质量以及单元6的截面积为区间不确定参数，取值区间为

m3=[144, 216]kg,A6=[3.0, 5.0]×10-4m2

表4　单元参数

表5　节点参数

对该模型节点3竖向自由度的频响函数进行求解。图6、图7分别将区间线性方程法及其改进方法求出的频响函数包络值与蒙特卡洛方法结果进行对比。由于此为无阻尼系统，因此对频响函数最终结果进行处理，去掉结果在(0, 10-11)以及(10-5,5)数值段的数据，保留其他数据进行对比分析。

图6、图7进一步证明了改进的区间线性方程方法的精度。与文献[12]方法相比，改进的迭代方法在保证精度的情况下，不会出现文献[12]中结构部分频段频响函数不能求解的情况，因此在区间参数结构频响函数求解上具有更大的优越性。

图5　二维梁模型Fig.5 Two dimension beam model

图6　区间线性方程法与蒙特卡洛方法结果对比Fig.6 Results comparison between Liner interval equation method and Monte-Carlo method

图7　改进的迭代方法(16子区间)与蒙特卡洛方法对结果比Fig.7 Results comparison between modified iteration method and Monte-Carlo method

4　结　论

本文对区间不确定工程结构的动力特性频响函数的求解，给出了一种改进的迭代方法。将求解频响函数的问题转化为区间线性方程的形式，并按照区间的运算规则，运用定点迭代的方法，得出结构频响函数上、下界的求解方法。运用子区间的思想，对迭代格式进行改进，保证迭代的精度和收敛性。通过对算例计算结果的对比和分析表明，文中的改进方法能求出较精确的解、计算效率高且能得到全频段的响应区间。

[1] 苏静波, 邵国建. 基于区间分析的工程结构不确定性研究现状与展望[J]. 力学进展, 2005, 35(3):338-344.

SU Jing-bo, SHAO Guo-jian. Current research and prospects on interval analysis in engineering structure uncertainty analysis [J]. Advances in Mechanics, 2005, 35(3): 338-344.

[2] 邱志平, 王靖. 不确定参数结构特征值问题的概率统计方法和区间分析方法的比较[J]. 航空学报, 2007, 28(3):590-592.

QIU Zhi-ping, WANG Jing. Comparison of probabilistic Statistical method and interval analysis method for Eigenvalue problem of structure system with uncertain parameters [J]. Acta Aeronautica Et Astronautica Sinica-Series,2007, 28(3):590-592.

[3] Metropolis N, Ulam S. The monte carlo method [J]. Journal of the American statistical association,1949,44(247):335-341.

[4] Moore R E, Kearfott R B, Cloud M J. Introduction to interval analysis [M]. Society for Industrial and Applied Mathematics, Philadelphia, 2009.

[5] Moens D, Hanss M. Non-probabilistic finite element analysis for parametric uncertainty treatment in applied mechanics: Recent advances [J]. Finite Elements in Analysis and Design, 2011, 47(1): 4-16.

[6] Ben-Haim Y, Elishakoff I. Convex models of uncertainty in applied mechanics [M]. Elsevier, 2013.

[7] Chen S H, Qiu Z P, Song D. A new method for computing the upper and lower bounds on frequencies of structures with interval parameters [J]. Mechanics Research Communications, 1994, 21(6): 583-592.

[8] 邱志平, 陈塑寰, 刘中生. 区间参数结构振动问题的矩阵摄动法[J]. 应用数学和力学, 1994(6): 519-527.

QIU Zhi-ping, CHEN Su-huan, LIU Zhong-sheng. Matrix perturbation method of structure vibration problems with interval parameters [J]. Applied Mathematics and Mechanics, 1994(6): 519-527.

[9] 李金平, 陈建军, 朱增青,等. 结构区间有限元方程组的一种解法[J]. 工程力学, 2010, 27(4): 79-83.

LI Jin-ping, CHEN Jian-jun, ZHU Zeng-qing,et al. A method for solving the structural interval finite element equations [J]. Engineering Mechanics, 2010, 27(4): 79-83.

[10] 苏静波, 邵国建, 刘宁. 基于单元的子区间摄动有限元方法研究[J]. 计算力学学报, 2007, 24(4): 524-528.

SU Jing-bo, SHAO Guo-jian, LIU Ning. Static sub-interval perturbed finite element method based on the elements [J]. Chinese Journal of Computational Mechanics, 2007, 24(4): 524-528.

[11] Muhanna R L, Zhang H, Mullen R L. Interval finite elements as a basis for generalized models of uncertainty in engineering mechanics [J]. Reliable Computing, 2007, 13(2): 173-194.

[12] Moens D, Vandepitte D. An interval finite element approach for the calculation of envelope frequency response functions [J]. International Journal for Numerical Methods in Engineering, 2004, 61(14): 2480-2507.

[13] Dessombz O, Thouverez F, Lané J P, et al. Analysis of mechanical systems using interval computations applied to finite element methods [J]. Journal of Sound and Vibration, 2001, 239(5): 949-968.

[14] Rump S M. On the solution of interval linear systems [J]. Computing, 1992, 47(3/4): 337-353.

[15] 梁震涛, 陈建军, 王小兵. 不确定性结构区间分析的改进MonteCarlo方法[J]. 系统仿真学报,2007,19(6):1220-1223.

LIANG Zhen-tao, CHEN Jian-jun, WANG Xiao-bing. Improved Monte Carlo method for interval analysis of uncertain structures [J]. Journal of System Simulation, 2007, 19(6): 1220-1223.

A modified iteration method to solve frequency response function of structures with interval parameters

FAN Zhi-ruo1,2, JIANG Dong3, DONG E-liang1,2, FEI Qing-guo1,2

(1. Department of Engineering Mechanics, Southeast University, Nanjing 210096, China; 2. Jiangsu Key Laboratory of Engineering Mechanics, Nanjing 210096, China;3. College of Mechanical and Electronic Engineering, Nanjing Forestry University, Nanjing 210096, China)

A modified iteration method for solving frequency response function of a structure with uncertain parameters was presented. Interval parameters were adopted to represent structural uncertain parameters. The equation to solve a structural frequency response function was converted to a group of interval linear equations. The iteration matrix was obtained by using the fixed point iteration theorem. Due to a convergence problem in iteration, sub-intervals were introduced and a modified iteration method was proposed to get the upper and lower bounds of the frequency response function. The accuracy of the solution was improved and the convergence problem was solved after modification. A spring-mass system was adopted in simulation to explore the influence of sub-interval’s quantity. At last, a comparative study was conducted using a structure in reference. The results showed that the proposed method is suitable for a full frequency band frequency response analysis with better accuracy and computation efficiency.

uncertainty; interval analysis; interval linear equation; sub-interval; frequency response function

10.13465/j.cnki.jvs.2016.13.004

国家自然科学基金(10902024);教育部新世纪优秀人才支持计划(NCET-11-0086)

2015-04-27修改稿收到日期：2015-06-26

范芷若 女，硕士生，1991年生

董萼良 男，副教授，1966年生

TB122

A