多本构模型参数同步更新的在线数值模拟方法研究

2023-11-01刘天姿黄伟宁西占叶恒博王国元

世界地震工程 2023年4期

刘天姿,黄伟,宁西占,2,叶恒博,王国元

(1. 华侨大学土木工程学院,福建厦门 361021; 2. 中国地震局工程力学研究所地震工程与工程振动实验室,黑龙江哈尔滨 150080; 3. 芜湖市交通投资有限公司,安徽芜湖 241011)

0 引言

混合试验源于拟动力试验,是一种将物理试验和数值计算相结合来获取动力荷载作用下结构反应的先进试验方法[1]。在该方法中:待模拟结构中的关键构件或子结构在实验室中进行物理加载,称为物理子结构;而结构剩余部分以及全部的惯性力、阻尼力则在计算机中采用程序或有限元软件计算,称为数值子结构。因此,混合试验可以实现大比例尺的试验,且对加载设备无特殊需求,在过去的数十年间得到了密切关注,并取得了许多重要进展[2-5]。

为实现复杂结构的混合试验,物理子结构往往采用简化处理的边界条件,或放弃某些自由度的模拟,这样势必会降低混合试验的模拟精度,甚至导致错误的结构破坏模式。为解决物理子结构边界条件无法完全模拟的混合试验问题,WU等[6]提出了在线数值模拟方法。该方法的思想是:结构不再划分数值子结构和物理子结构,其动力反应由整体数值模型计算得到,边界条件问题自然得以避免;不完整物理子结构的实测反力不再参与运动方程的求解,而是用于数值模型本构参数估计,并以此更新整体数值模型的本构参数值,数值建模的可靠性得到了有效改善。

在线数值模拟方法的关键问题之一是本构模型参数估计。2005年YANG等[7]以实测试验数据使用神经网络在线识别物理子结构的滞回特性,并用经过物理子结构试验数据训练好的神经网络在线预测数值子结构的恢复力,以修正数值子结构数学模型与试验模型不同的影响,这是模型更新混合试验的首次应用。随后,地震工程领域学者在本构模型参数估计方面开展了大量有益探索。王涛等[8]以传统BP神经网络为基础,通过增加反馈层,提出一种在线自适应神经网络算法,并以Bouc-Wen模型为例进行了验证,研究表明:所提方法具有较高的效率和精度;王燕华等[9]将遗忘因子引入LMBP神经网络,发展了一种在线模型更新方法,并以两自由度非线性结构为例探讨了新方法的自适应性、稳定性和抗噪声能力;马天宇等[10]以最小二乘法为基础,对双折线模型进行恢复力模型更新,并开展了数值验证;CHUANG等[11]开发了一种防屈曲支撑的双屈服面宏观模型,并提出基于梯度方法的模型更新算法,完成了钢框架-防屈曲支撑结构的数值模拟;郭玉荣等[12]以Ibarra-Medina-Krawinkler模型作为钢筋混凝土的本构模型,采用UKF作为估计方法,验证了模型更新方法的有效性和优越性;WU等[13]基于UKF算法提出了截面模型更新的混合试验方法,并以钢框架结构为例验证了所提方法的有效性;MEI等[14]将UKF算法用于识别混凝土材料的本构模型参数,完成了高墩桥的模型更新混合试验,有效改善了数值建模的可靠性;ZHONG等[15]提出基于全局敏感性的模型更新混合试验方法,并以钢筋混凝土结构为例对该方法进行了系统探讨。

在已有基于模型更新的混合试验研究中,往往认为单一材料本构模型参数具有不确定性并对其进行估计更新,如钢筋或混凝土材料。但对钢筋混凝土结构,由材性试验得到的钢筋和混凝土的本构模型参数在用于结构计算时难以精确模拟其力学性能。因此,有必要对钢筋和混凝土本构模型参数同步估计以提高数值子结构的模拟精度。基于此,本文提出多本构模型参数同步更新的在线数值模拟方法,并以Levenberg-Marquardt(LM)算法作为参数估计方法,采用物理子结构的恢复力实现对钢筋和混凝土材料本构模型参数的同步估计。本文余下部分将阐述在线数值模拟方法的基本原理和基于LM算法的有限元材料本构模型参数估计方法,以一钢筋混凝土框架结构为例,探讨LM算法中阻尼因子的影响,并开展虚拟混合试验以验证所提方法的可行性。

1 在线数值模拟方法基本原理

在线数值模拟方法是一种新型混合试验方法,其既可以有效解决不完整边界物理子结构的混合试验问题,又可以改善数值建模的可靠性。与传统混合试验不同,物理子结构不再为结构反应计算提供恢复力,而是用于数值模型的本构参数估计。在线数值模拟方法包含四个模块,即:1)协调器模块,用于求解结构的动力反应。2)结构整体数值模型模块,用于确定物理子结构的边界位移并提供结构动力反应计算的整体结构静恢复力。3)物理子结构模块。4)本构模型参数估计模块,用于估计数值模型的本构参数,包含优化算法以及与物理子结构具有相同边界条件的数值模型(本文称之为“等代物理子结构”)。在线数值模拟方法的原理图见图1。

图1 在线数值模拟方法原理Fig. 1 Principle of online numerical simulation method

在线数值模拟方法的基本流程可简述如下:1)在协调器中求解结构运动方程,获取结构动力自由度方向上的位移d。2)以位移d驱动结构整体数值模型完成非线性静力分析,得到物理子结构位移dE,并将其发送给试验设备进行物理加载,获取相应的恢复力RE。3)参数估计模块以物理子结构实现的位移及其对应的恢复力,采用合适的参数估计方法估计等代物理子结构的相关本构参数,并将该参数更新结构整体数值模型中的相关参数。4)再次以位移d驱动结构整体数值模型完成非线性静力分析,提取与动力自由度相应的自由度方向上的恢复力R。5)将整体结构静恢复力反馈至协调器进行下一步运动方程的求解。6)重复1)至5)直至试验结束。

需要说明的是:在线数值模拟方法建立在对物理子结构具有一定认知的前提下,且该方法不限制本构模型的类型和参数个数,具有较强的适用性和扩展性。

2 基于LM算法的参数估计方法

近年来,有限元方法逐渐被用于混合试验以提高模拟精度,其中采用材料应力-应变关系反应结构或构件力学行为的纤维模型兼顾计算效率和精度,既能描述结构的整体行为,又能刻画结构的局部特性,因此在混合试验中得到了广泛的应用。在混合试验中:能够直接有效观测的物理量包括应变、物理子结构的力和位移,但应变易受温度和电阻影响。为此,该研究将以物理子结构的恢复力为观测量,来估计材料的微观本构模型参数。特别地,本文将以单作动器加载的物理子结构为例,对钢筋和混凝土这两种材料的本构模型参数进行估计。

2.1 本构模型参数的优化问题

参数估计就是通过不断优化等代物理子结构的本构模型参数,使由等代物理子结构计算得到的恢复力与物理子结构的实测力之差最小,可将其归结为结构参数优化问题。从有限元角度来看:等代物理子结构的恢复力F是节点力,可以表述成节点位移和本构模型参数的函数为:

F=f(x,u)

(1)

式中:f(x,u)是与历史变量相关的函数,x是本构模型参数,u是节点位移。显然,在基于有限元的本构模型参数估计问题中,恢复力是本构模型参数x的非线性函数。

假定物理子结构的实测恢复力为R,则参数估计的目标函数可表达为:

(2)

式中:k代表当前计算步。显然,以式(2)为目标函数的参数估计问题可归结为非线性最小二乘问题,需通过迭代求解。

2.2 LM算法

在试验初始阶段,由于观测量不足等因素影响,常导致目标函数的Hessian矩阵(JTJ,J表示目标函数的雅克比矩阵)为奇异阵,此时式(2)无解或存在多组解,无法实现对多本构模型参数的同步估计。为避免该问题,Levenberg和Marquardt等提出了一种阻尼高斯-牛顿方法[16],在第k个计算步的迭代格式为:

xj+1=xj+(JTJ+μjI)-1JT(f(xj,u)-R)

(3)

式中:j为迭代步,u、f和R为前k个计算步中位移、恢复力和实测力的向量表示,I表示单位矩阵,μ为大于零的阻尼因子。由式(3)可知:在式(2)的Hessian矩阵上增加正的单位矩阵,可保证式(2)解的唯一性,因此采用LM算法可实现对多本构模型参数的同步估计。当阻尼因子趋于无穷较大时,式(3)退化为最速下降法;当阻尼因子为零时,式(3)退化为高斯牛顿法。本文中,阻尼因子的更新方式为:

(4)

(5)

对式(3)和式(5)中Jacobian矩阵,采用向前差分方法近似代替微分进行求解,计算式为:

(6)

式中:φ是前向差分步长因子,是一个较小的正数。

3 数值模拟

本节将以钢筋混凝土结构为例,验证基于LM算法的多本构模型参数同步更新在线数值模拟方法的可行性和有效性。研究中:物理子结构采用有限元模拟,将其采用的本构模型参数称为真实值,并将所开展的混合试验称为虚拟混合试验。

3.1 计算模型

选取一两层两跨钢筋混凝土平面框架为研究对象,其底层层高为3.6 m,二层层高3.0 m,柱为500 mm×500 mm的矩形截面,梁宽300 mm,高550 mm,并在每侧考虑720 mm的楼板翼缘以考虑楼板对梁的约束作用;梁柱中受力筋为直径20 mm的HRB400级钢筋。结构每层质量为13 218 kg,周期为0.37 s。建模时,梁柱均采用纤维截面离散,用OpenSees中基于力的非线性梁柱单元模拟,每个单元取5个积分点,钢筋和混凝土分别用Steel02和Concrete01模型模拟[17]。结构的计算模型如图2所示。研究中:取底层边柱反弯点以下部分为物理子结构,以避免对子结构界面处弯矩的模拟[18],并考虑几何相似比和弹性模量相似比分别为1∶2和1∶1。根据《建筑抗震试验规程》(JGJ/T 101—2015)[19]的规定可知:通过恢复力估计得到的本构模型参数可直接用于原结构数值模型的计算。

图2 计算模型及物理子结构Fig. 2 Computational model and physical substructure

3.2 LM算法中阻尼因子的影响

为更好地指导虚拟混合试验的开展,首先研究了LM算法中阻尼因子的影响,并在2.2节算法基础上考虑了两种常阻尼因子,即μ=0.1和μ=0.001。该算例以物理子结构为对象开展,并采用了指定位移输入,见图3;物理子结构本构模型的真实值和等代物理子结构的初始值见表1。在计算Jacobian矩阵时,前向差分步长因子φ取为1×10-9,以使由差分计算得到的Jacobian矩阵更接近微分情况。

表1 本构模型参数及其估计值Table 1 Constitutive model parameters and their estimated values

图3 加载位移时程Fig. 3 Time history of the loading displacement

为方便对比,选用不同阻尼因子时钢筋和混凝土参数的估计结果也在表1中给出。表1显示:尽管参数估计初值与真实值具有明显的差异,三种工况下参数估计的终值与真实值的差异明显减小,表明基于LM算法的参数估计方法具有较高的估计精度。这是因为:在混合试验中,参数估计的初始值往往位于真实值的邻域内,采用LM算法可以保证估计值收敛至真实值。从表中还可以发现:除约束区混凝土的极限应力与真实值有较大差异外,其余参数几乎与真实值完全一致,其原因是:在当前最大位移幅值下,物理子结构依然具有较高的承载能力,刚度尚未有明显退化,导致Concrete01材料下降段所起作用有限,从而导致约束区混凝土极限应力未能收敛到真实值。

3.3 虚拟混合试验

为说明本文方法的优越性和有效性,开展了多本构模型同步估计的在线数值模拟和传统混合试验。地震动选为El Centro波,时间间隔通过重采样设为0.005 s,地震动峰值加速度调幅为0.31 g,持时10 s。研究中:钢筋和混凝土材料的本构模型参数初始值和真实值见表2,并以OpenSees完成的整体结构动力时程分析作为参考解以评估两种混合试验的模拟效果。值得说明的是:由于混合试验中通常会对位移和力采取滤波措施,测量的位移和力往往很光滑,因此在虚拟混合试验时未考虑噪声的影响。

表2 虚拟混合试验中材料模型参数及其估计值Table 2 Constitutive model parameters and their estimated values

在线数值模拟和传统混合试验的位移时程对比见图4,同时将定量分析误差在表3中给出。图4和表3中显示:传统混合试验与参考位移有显著的差异,而在线数值模拟方法与参考解吻合良好,显示了在线数值模拟方法的优越性。传统混合试验的误差来源有两方面:一是数值计算模型采用了材料本构模型的初始值,与真实值之间存在明显差异,进而带来模拟误差;二是在动力分析时柱的反弯点位置和轴力是发生变化的,亦即物理子结构在边界处需同时考虑弯矩和轴力边界,而试验时则忽略了这两个边界,因此造成模拟误差。对在线数值模拟方法,结构反应由整体数值模型计算得到,不存在边界条件的问题,同时本构模型参数通过物理子结构的实测反力在线估计得以改善,因此与参考解吻合良好。

表3 混合试验误差指标Table 3 Error indexes for hybrid simulation

图4 结构响应对比Fig. 4 Comparison of structural response

图5给出了底层中柱弯矩和轴力的时程对比图。从图中可以看出:无论是轴力还是弯矩,在线数值模拟方法与参考解吻合良好;传统混合试验得到的轴力幅值与参考解相差不大,但弯矩幅值有明显差异。由此可以看出:混合试验中简化的边界条件将严重改变结构的受力状态,进而影响结构性能评估的准确性。

图5 轴力与弯矩对比Fig. 5 Comparison of axial force and bending moment

图6给出了在线数值模拟方法非约束区混凝土峰值应力和钢筋屈服强度的估计结果,其余参数的估计值在表2中给出。从表中可以看出:混合试验中参数估计值与真实值几乎一致,从而表明LM算法具有较高的精度。从图6中可以看出:峰值应力在前1 s内迅速调整并逐渐趋向于真实值,随后保持不变;而屈服强度在前1.5 s内几乎没有任何变化,随后迅速调整至真实值附近,并保持该值。这是因为:在地震作用初始阶段,结构的位移幅值较小,混凝土的应变尚未达到峰值应变,但由于混凝土的应力应变关系呈现非线性,即使此时恢复力实测值与计算值的差异较小,也会使得峰值应力迅速调整,以实现恢复力差值最小;随着位移的增加,混凝土峰值应力的作用被充分激发,而此时由本构参数估计值计算得到的恢复力与实测恢复力较为接近,因此在1 s后估计值保持不变,且几乎与真实值重合。而对钢筋来讲,当位移较小时其性能主要受弹性模量影响,屈服强度的作用未被激发,因此在初始阶段屈服强度几乎没有变化;随着钢筋应变的增加,屈服强度的作用逐渐显现,导致其估计值迅速向真实值调整。当钢筋应变超过其屈服应变时,钢筋屈服强度的作用被充分激发,进而使得其估计值几乎不发生变化。由此可以看出:LM算法具有较高的模拟精度和较快的收敛速度,适合基于模型更新的混合试验。