王 涛, 吴 斌, 孟丽岩, 许国山, 张 健, 尹晓黎
(1.黑龙江科技大学 建筑工程学院,哈尔滨150022;2.哈尔滨工业大学 土木工程学院,哈尔滨150090)
确定结构恢复力模型及其模型参数是进行结构非线性地震反应分析、健康监测及结构混合实验的关键技术。迄今为止,研究者们已提出一系列的恢复力模型[1]来描述力与变形之间的关系。从形式上可分为折线形或曲线形,从应用的对象可分为层间层次、构件层次、截面层次及材料层次的恢复力模型。其中,Bouc 等[2]建立了经典的Bouc -Wen模型。之后,Baber 等[3]进一步发展了Bouc-Wen 模型,能够考虑强度退化、刚度退化及滑移捏缩效应。
在过去的几十年里,EKF 被广泛地应用在大量土木工程应用中的结构损伤识别和非线性结构参数识别等领域[4]。由于EKF 需要在每一个递推步中求解相应的雅克比矩阵,当非线性系统不可导时,则很难求解雅克比矩阵。为了解决EKF 存在的问题,Julier 等[5]提出了隐性卡尔曼滤波器(Unscented Kalman filter,UKF)。与EKF 不同,UKF 采用确定性采样策略去逼近状态量的统计特征值,即隐性变换(Unscented transformation,UT),然后利用卡尔曼滤波的预测-更新递推算法框架进行状态估计。UKF 同时具有很强的非线性状态估计能力和较高的计算效率。
2008年,Wu 等[6]应用UKF 在线识别考虑退化、滑移的非线性滞回模型参数,研究表明,UKF 算法具有较高的识别精度和计算效率,能够实现实时状态估计和参数识别。Xie 等[7]提出迭代UKF 算法,用来在线识别非线性结构参数,该算法提高了UKF 对观测噪声的鲁棒性。2013年,王涛等[8-9]在UKF 算法基础上提出了约束UKF(Constrained UKF,CUKF)算法,并将此算法应用于结构混合实验中。
采用UKF 或CUKF 进行参数识别时,需要事先给出初始估计误差协方差、过程噪声协方差和观测噪声协方差等滤波器的初始参数,这些初始参数的设置会直接影响参数识别结果。对于Bouc -Wen模型这样的非线性系统,采用解析的方式来分析滤波器初始参数对模型参数识别影响规律将变得异常困难。为此,笔者提出Bouc-Wen 模型在线参数识别方法,采用数值模拟方式分析约束UKF 滤波器的初始参数对Bouc -Wen 模型参数识别精度和收敛速度的影响规律。
UKF 的基本思想就是将UT 变换和卡尔曼滤波框架相结合,进行非线性状态估计。在卡尔曼滤波框架下,预测步中通过确定采样的方法生成2n1+1个Sigma 点,将这些Sigma 点通过状态方程,然后加权统计变换的Sigma 点,可以得到先验状态估计和协方差;在更新步,将状态Sigma 点通过观测方程,然后加权统计变换后的观测Sigma 点,得到先验的观测先验估计和协方差及状态和观测的互协方差,进而可得到卡尔曼增益矩阵,然后通过卡尔曼更新方程,得到状态后验估计和协方差矩阵。
UKF 方法具体算法[10],将一个非线性动力系统表达为离散状态方程:
观测方程为
式(1)、(2)中,xk为系统状态向量为状态向量维数;yk+1为系统观测向量为观测向量维数;uk为系统输入;vk为过程噪声,假定为高斯白噪声,其均值为零,协方差矩阵为Qk;wk+1为观测噪声,假定为高斯白噪声,其均值为零,协方差矩阵为Rk+1。
状态估计的目的就是根据当前第k 步及第k 步之前的观测信息,及系统的状态方程和观测方程来计算第k+1 步系统状态估计值。计算第k 步2n1+1 个确定性Sigma 点Xk,i和其相应的权重:
上面得到的第k +1 步系统状态估计均值和协方差,并没有考虑到第k +1 步系统观测yk+1的影响。故可由观测方程式(2)计算观测向量Sigma点
应用卡尔曼滤波更新方程,计算第k+1 步系统状态估计均值和协方差Pk+1:
式(15)、(16)中,Kk+1为卡尔曼滤波增益矩阵,
可以看出,标准的UKF 算法并未考虑有界限约束,或其他数学约束条件的状态估计问题。若状态估计值违反了约束条件,将会影响在线模型更新准确性,甚至会影响非线性模型本身的稳定性。为了实现约束状态估计,采用基于Sigma 点调整的约束UKF 算法[9],该算法在标准UKF 算法基础上进行以下两个方面的改进:第一,在预测步中,将违反界限约束的Sigma 点沿着协方差的方向移至界限约束边界上,同时对在约束界限内的Sigma 点作相应移动,以保持重新调整后的一组Sigma 点仍具有对称性;第二,在校正步中,采用状态更新方程产生变换后的Sigma 点,当更新的状态估计违反约束时,将违反约束的Sigma 点垂直投影到约束边界上,通过加权统计的方式计算后验的状态估计均值和协方差。
经典Bouc-Wen 滞回模型的数学表达式[2]:
式(19)、(20)中,r 为恢复力;d、v、z 分别为位移、速度、滞变位移;k、α、β、γ、n 分别为控制滞回环形状和尺寸的五个模型参数。当n =∞时,Bouc -Wen模型骨架曲线近似为双折线模型,此时,k 和α 可以被理解为双折线模型的初始刚度和第二刚度系数。然而,Bouc -Wen 模型参数的物理意义并不明确,这使得确定这些模型参数变得十分困难。我们通过在线参数识别来得到Bouc-Wen 模型的参数,即k、β、γ、n、α 和时变的滞变位移z。
采用约束UKF 方法在线识别模型参数,需要定义系统的状态向量x:
系统观测为恢复力,即
位移d 和速度v 为系统的已知输入,则系统的连续状态方程为
系统的观测方程
设Bouc-Wen 模型参数真实值分别为k =135 kN/mm,β=0.2,γ =0.2,n =1,α =0.02。设真实系统恢复力观测噪声wk+1为离散的高斯白噪声序列,噪声标准差为真实系统输出恢复力标准差的2%,其大小相当于真实恢复力峰值的0.84%。位移输入为曾完成的防屈曲支撑混合实验的位移数据,位移输入时程如图1a 所示;速度输入采用位移向前差分方法计算,得到速度输入时程如图1b 所示。
图1 Bouc-Wen 模型位移和速度输入Fig.1 Displacement and velocity inputs of Bouc-Wen model
模型参数界限约束为:k≥0,β≥-γ,n≥1,0≤α≤1。过程噪声协方差设为Q =10-6I6;观测噪声协方差设为R=7.839 6 kN2。
状态初始估计值为
初始状态估计误差协方差为
建立了系统状态空间模型和滤波器初始参数之后,就可以采用递推形式的约束UKF 算法来在线识别Bouc-Wen 模型参数。
为了分析初始协方差对参数识别的影响,设观测噪声协方差P0=kPP,分别对比kP为0.1、1.0、10.0 三种情况下模型参数识别结果,此时其余滤波器参数保持不变,即不同P0下参数识别值时程曲线对比如图2 所示。参数识别终值及其相对误差对比如表1 所示。
图2 P0 对Bouc-Wen 模型参数识别值影响Fig.2 Estimated parameters of Bouc-Wen model with different P0
表1 不同P0 下参数识别终值及其相对误差Table 1 Final identification values and relative errors with different P0
图2 的结果表明:在无模型误差情况下,随着kP的增加参数识别值收敛速度明显加快,识别终值的相对误差绝对值变小。当kP=10.0 时,除了参数α 以外其余参数识别终值的相对误差均小于2.5%。同时也可以看出,kP从0.1 变化到10.0,约束UKF 均可以正常工作,可见识别算法对P0的选择具有较好的鲁棒性。因此,当我们对参数初值信心不是很大时,选择较大的P0可以提高参数识别收敛速度和精度。
设观测噪声协方差Rk= kRR,分别对比kR为0.1、1.0、10.0 三种情况下模型参数识别结果,此时其余滤波器参数保持不变,即Q。参数识别终值及其相对误差对比见表2。不同Rk下参数识别值时程曲线对比如图3 所示。
表2 不同Rk 下参数识别终值及其相对误差Table 2 Final identification values and relative errors with different Rk
图3 的结果表明:在无模型误差情况下,当kR=1.0时,即滤波器的观测噪声协方差与真实系统观测噪声协方差相等时,参数识别值收敛速度更快,识别终值的相对误差绝对值变小。当kR=1.0 时,除了参数α 以外其余参数识别终值的相对误差均小于7.5%。kR越小,说明我们更相信观测,这样前期参数识别受观测噪声的影响更大,识别值的波动性会加大;相反,kR越大,我们会更相信初值选择,识别值的收敛速度会变慢。同时也可以看出,kR从0.1 变化到10.0,约束UKF 均可以正常工作,可见识别算法对Rk的选择具有较好的鲁棒性。因此,最好能够事先统计真实的系统观测噪声协方差,将此作为滤波器的观测噪声协方差,这样,可以提高参数识别收敛速度和精度。
图3 Rk 对Bouc-Wen 模型参数识别值影响Fig.3 Estimated parameters of Bouc-Wen model with different Rk
设过程噪声协方差Qk=kQQ,分别对比kQ为0、10.0、100.0 三种情况下模型参数识别结果,此时其余滤波器参数保持不变,即参数识别终值及其相对误差对比见表3。不同Qk下参数识别值时程曲线对比如图4 所示。
表3 不同Qk 下参数识别终值及其相对误差Table 3 Final identification values and relative errors with different Qk
图4 的结果表明:在无模型误差情况下,当kQ=0 时,参数识别值收敛速度更快,识别终值的相对误差绝对值变小。当kQ=0 时,除了参数α 以外其余参数识别终值的相对误差均小于5.5%。kQ越大,参数识别值的波动性会加大,收敛速度会变慢。同时也可以看出,kQ从0 变化到10.0,约束UKF 均可以正常工作,而且对识别结果影响并不显著,可见识别算法对Qk的选择具有较好的鲁棒性。因此,当真实系统无模型不确定性问题时,可以不用考虑过程噪声,这样可以提高参数识别收敛速度和精度。对于实际情况,系统模型一般不可能完全精确描述真实系统,可以通过设置较小的过程噪声协方差来减小模型误差影响。
图4 Qk 对Bouc-Wen 模型参数识别值影响Fig.4 Estimated parameters of Bouc-Wen model with different Qk
图5 对Bouc-Wen 模型参数识别值影响Fig.5 Estimated parameters of Bouc-Wen model with different
表4 不同 下参数识别终值及其相对误差Table 4 Final identification values and relative errors with different
表4 不同 下参数识别终值及其相对误差Table 4 Final identification values and relative errors with different
注:相对误差(%)=100% × |识别值-真实值|/真实值。
参数 真实值 初始值kx =0.8识别值 误差/%kx=1.0识别值 误差/%kx =1.2识别值 误差/%k 135.00 115.0 132.471 1.87 133.080 1.41 134.134 0.64 β 0.20 0.5 0.200 0 0.200 0 0.201 0.50 γ 0.20 0.5 0.178 11.00 0.185 -7.50 0.196 2.00 n 1.00 2.0 1.083 8.30 1.069 6.90 1.051 5.10 α 0.02 0.1 0.035 75.00 0.038 90.00 0.043 115.00
图5 的结果表明:在无模型误差情况下,kx从0.8 变化到1.2,约束UKF 均可以较好地识别模型参数,除了参数α 以外其余参数识别终值的相对误差均不大于8.3%,可见识别算法对^x0的选择具有较好的鲁棒性。整体上,初始参数值与真实值越接近,识别值收敛会越快,识别终值的相对误差绝对值也相应变小。因此,当对参数初值^x0的信心更有把握时,可以提高参数识别收敛速度和精度。
(1)与标准的UKF 算法相比,约束UKF 法在不增加计算负荷的基础上,可以解决考虑有界限约束的非线性状态估计问题。可以减小参数识别值的前期波动幅度,提高识别值的收敛速度,有助于提高在线模型更新准确性。文中给出基于约束UKF 的Bouc-Wen 模型参数识别方法。
(2)数值模拟研究表明,在无模型误差的情况下,滤波器对于初始参数具有较好的鲁棒性,初始参数可以具有较宽的取值范围。
(3)给出在实际应用中滤波器参数选择建议,通过适当的增大初始状态估计协方差,减小过程噪声,采用真实系统观测噪声协方差,以及减小初始参数值与真实值的偏差,可以有效提高参数识别收敛速度和精度。
(4)该研究可以为解决土木结构非线性模型在线参数识别和模型更新混合实验提供参考。
[1]郭子雄,杨 勇.恢复力模型研究现状及存在问题[J].世界地震工程,2004,20(4):47 -51.
[2]WEN Y K.Method for random vibration of hysteretic systems[J].Journal of the Engineering Mechanics Division ASCE,1976,102(2):249 -263.
[3]BABER T T,WEN Y K.Random vibration hysteretic,degrading systems[J].Journal of the Engineering Mechanics Division ASCE,1981,107(6):1069 -1087.
[4]YANG J N,LIN S,HUANGL H,et al.An adaptive extended Kalman filter for structural damage identification[J].Structural Control and Health Monitoring,2006,13(4):849 -867.
[5]JULIER S J,UHLMANN J K.Unscented filtering and nonlinear estimation[J].Proceendings of the IEEE,2004,92(3):401-422.
[6]WU M,SMYTH A.Real-time parameter estimation for degrading and pinching hysteretic models[J].International Journal of Non-Linear Mechanics,2008,43(9):822 -833.
[7]XIE Z B,FENG J C.Real-time nonlinear structural system identification via iterated unscented Kalman filter[J].Mechanical Systems and Signal Processing,2012,28:309 -322.
[8]王 涛,吴 斌.基于UKF 模型更新的混合试验方法[J].振动与冲击,2013,32(5):138 -143.
[9]王 涛,吴 斌.基于约束UKF 模型更新的混合试验方法[J].地震工程与工程振动,2013,33(5):100 -109.
[10]WAN E,VAN DER MERWE R.Kalman filtering and neural networks[M].Wiley:[s.n.],2001:221 -280.