ADVANCED HOMOTOPY STOCHASTIC FINITE ELEMENT METHOD FOR STRUCTURAL ELASTIC STABILITY ANALYSIS
-
摘要: 结构参数的不确定性将对结构稳定性产生不可忽视的影响，实现随机结构屈曲荷载与屈曲模态的高效高精度求解和统计分析，对结构设计与安全评估有重要意义。该文基于随机残差最小化法改进现有的同伦随机有限元法，并利用新方法高效高精度地求解了大变异随机参数结构的屈曲特征值和屈曲模态。将随机参数结构的屈曲特征值和屈曲模态以同伦级数的形式表达，并给出了同伦级数中任意阶系数的显式递推表达式；在此基础上，定义了关于弹性屈曲方程近似解的随机残余误差，通过使该随机残余误差最小化，得到了优化的随机屈曲特征值和屈曲模态的同伦级数展开表达式。该文提出的改进的同伦随机有限元法能够实现同伦级数展开的自动寻优，有效避免了现有同伦随机有限元法(HSFEM)计算精度易受样本选点影响的缺点。当随机参数变异性较大时，随着级数展开阶数的增加，该文方法计算结果除了能保持良好的收敛性外，相比于HSFEM具有更好的稳定性，而摄动随机有限元法则可能出现发散现象；与蒙特卡洛模拟法相比，新方法具有很高的求解效率。通过强非线性函数算例以及变截面轴心受压杆和框架结构的弹性稳定性分析说明了该文方法的有效性。Abstract: The uncertainty in structural parameters has notable impact on stability analysis of structures. It is of great significance for structural design and safety evaluation to obtain random buckling loads and buckling modes with efficiency and high-precision. The homotopy stochastic finite element method is improved to effectively solve the structural elastic stability problem involving large fluctuation of random parameters, and the statistical properties of the buckling eigenvalue and buckling mode are obtained. The buckling eigenvalue and buckling mode of the structure involving random parameters are expressed using the homotopy series, and the arbitrary order coefficients of the homotopy series are given as explicit recursive relationship formulas. Further, the stochastic residual error with respect to the buckling governing equation is defined, and the optimal form of the homotopy series is determined by minimizing the stochastic residual error. The proposed advanced homotopy stochastic finite element method can automatically realize the optimization process, which overcomes the drawbacks of the existing homotopy stochastic finite element method, named HSFEM, that the computational accuracy is affected by the selected samples and relies on empirical knowledge. In addition, for a structure involving large fluctuation of random parameters, the proposed method has better stability than HSFEM when higher-order terms in homotopy series are employed, while the results from the traditional perturbation method based on Taylor series may diverge. And the proposed method has an excellent computational efficiency compared with the Monte Carlo simulation method. The validity of the proposed method is verified through a strong-nonlinearity function, the stability analysis of a variable cross-section column subjected to axial force and the stability analysis of a 7-story frame structure.
-
索膜结构受力高效、造型优美，集建筑学、结构力学、材料科学、计算机技术等于一体，长期以来得到设计人员与建设单位的关注与青睐[1]。索膜结构设计理论与建造技术的飞速发展应用相比[2-3]，在役索膜结构的检测、检查、评估、鉴定技术相对落后，难以为此类工程结构的安全运营提供长期有效的技术保障。虽然部分国内外专家学者在这些方面进行了针对性的研究，但迄今为止，在役索膜结构的损伤识别与性能评估尚有诸多问题需要解决。
在役索膜结构的安全性评估，无法直接参照国内现行的结构可靠性鉴定标准。根据《民用建筑可靠性鉴定标准》(GB 50292−2015)[4]、《工业建筑可靠性鉴定标准》(GB 50144−2019)[5]中的规定，钢构件的承载力、变形、缺陷、锈蚀均对安全性等级评定均有影响；《高耸与复杂钢结构检测与鉴定标准》(GB 51008−2016)[6]规定，大跨度及空间钢结构的安全性鉴定应按结构整体性和结构承载安全性两个项目分别评定等级，上述标准并无针对索膜结构的专门规定。从结构特点上讲，索膜结构属于柔性结构，施工阶段的初始张拉应力大小与索力分布决定了结构的刚度和形状，在结构服役阶段，索、膜的损伤与张拉应力的损失同样会显著影响结构位形。因此，服役阶段的索膜结构不能直接使用设计模型进行结构性能分析与安全评价，需要全面、准确地测量结构的实际位形以及索力、损伤等指标，采用基于实测数据的模型修正方法进行结构性能分析。然而，工程实践与文献资料中未见成熟可靠的在役索膜结构模型修正技术与方法。
本文提出基于空间位形参数的在役索膜结构有限元模型修正方法以实现结构的性能评估。采用现场三维扫描获取结构的真实空间位形信息，以模型中拉索的初始预应力作为优化参数、以结构关键节点坐标的匹配作为模型修正目标函数，设计ABAQUS与MATLAB联合仿真优化程序，通过全局搜索与局部优化实现精细有限元模型的迭代修正。修正后模型与实际结构受力状况更为接近，索力的实测与对比结果表明，有限元模型的索力由修正前最大偏差10%~30%降低至10%以下。进一步，利用修正后有限元模型对在役结构进行状态评估并得出可靠结论，为实际工程中同类结构的安全性检测鉴定、损伤识别与性能评估以及相关研究提供参考。
1 研究现状
1.1 空间位形参数获取
在役索膜结构的位形与结构的内部应力状态具有直接联系，因此变形检测是对索膜结构进行评估的重要项目，目前获取索膜结构三维信息方法与其他空间结构并无很大差别，以下几种方法均可用于获取索膜结构的三维信息：
1) 常规大地测量法。利用常规的大地测量仪器测量方向、角度、边长、高差等技术来测定变形或三维信息，该方法灵活度、精度较高，但效率较低。此外，由于构件本身具有尺寸参数，该方法选取的目标点一般是构件表面的特征点，对于大型空间结构来说，有时通过选取的特征点不能准确地表达构件的空间位置信息。
2) 三维激光扫描。三维激光扫描能大范围、高精度、高分辨率地，以非接触的方式快速获取目标表面每个采样点的三维坐标数据。目前三维激光扫描被广泛地用于工程中的变形测量[7]、三维建模[8]等领域。该方法具有较高的采样率及较高的精度，足够用于工程测量分析。
3) 数字图像法多视角几何三维重建。采集目标图像信息，可使用基于数字图像的多视角几何三维重建方法，获取结构三维点云模型与关键节点空间坐标。现场操作相对简单、设备要求低，在测绘工作与三维空间建模[9]中应用广泛。当前各种算法的实现，逐渐提高了多视角几何重建法的精度[10]。
对于体量较大的索膜结构，采用数字图像法多视角几何三维重建的后期处理工作量大、计算需求高，因此，本文采用三维激光扫描方法直接获取索膜结构的空间信息。
1.2 有限元模型修正方法
有限元模型修正已成为在结构工程各个领域中改进数值模型的一种广泛使用的方法。如今，修正、优化后的有限元模型，已成为评估整个结构使用寿命期间结构的承载能力、安全性和可维修性的必不可少的工具。对在役索膜结构进行安全性评估的一个关键点在于建立符合结构当前服役性态的精确有限元模型，以便开展准确的承载力计算。
目前，对在役索膜结构进行损伤识别或安全性分析，部分研究以及工程实践仍是基于原始设计数据建立的结构有限元模型[11]。由于索膜结构本身的特点，实际建成结构的空间位形可能与设计模型存在差异，尤其是对于服役时间较长的在役索膜结构，采用设计模型可能会出现很大误差。对于类似索膜结构的空间预应力结构有限元模型的修正方法，既有的研究较少，而且难以对修正后模型进行评估[12]。因此，需要建立适合索膜结构的有限元模型修正方法。
有限元模型修正一般基于结构的静力或动力试验数据，通过调整模型的参数(构件的刚度、材料密度等)使得模拟结果与试验结果基本吻合，同时使得有限元模型能够更精确地反映结构的实际特性。修正方法多是采用迭代的方法对模型中的相关参数进行调整，优化的一般流程如图1所示。迭代方法是将模型修正问题转化为优化问题，其目标函数一般是根据数值模拟结果与试验结果的差异来定义的，如：
min (1) 式中：
{Y_0} 为试验测量值；{Y_1}(X) 为有限元模型计算值；X 为优化参数组成的向量。目前，对于有限元模型的修正方法的研究很普遍，也有较多的优秀成果，从修正方式到内部算法，基本上形成了较为完善的系统，但对于具体的结构形式，需要更加细致地调整修正方法。从修正方式来看，一般基于静力学或动力学数据进行优化，追求优化效率与准确度。翁顺、朱宏平[13]提出了基于子结构有限元模型修正方法，可大幅提高大型复杂结构的有限元模型的优化效率。Wendy 等[14]基于实测的固有频率对斜拉桥有限元模型进行修正，提高了有限元模型的分析效果。Hendrik等[15]引入了多响应目标函数，根据实测数据优化有限元模型，显著提高了有限元模型的质量。沈雁彬等[12]提出一种基于索力敏感度分析的有限元模型修正方法，利用有限元软件计算得到索力敏感度矩阵信息，将其代入结构零状态和初始态循环迭代的过程，最终得到能够准确反映结构实际受力状态的数值模型。从优化算法来看，非线性优化算法的使用更为普遍[16-17]，同时，新算法的尝试与改进也在进行。Javier等[18]提出基于两种优化算法(和声搜索和有效集算法)及人工神经网络的有限元模型修正方法，能够显著降低模拟计算时间。
对于索膜结构，结构位形与结构内力联系紧密，因此对索膜结构有限元模型的优化，也可以认为是通过对相关参数的优化调整，使得结构位形与结构内力的关系与实测结果相吻合。基于这一点，可以通过试验得到的结构位形信息，对模型进行修正，并以结构内力作为参照对优化结果进行评估；反之亦可用实测结构内力进行模型修正。
2 有限元模型修正方法与实现
本文提出的索膜结构有限元模型修正方法主要是基于实测的结构位形指标，选择结构中拉索的预应力作为优化参数，因为造成索膜结构位形发生变化的主要是荷载、施工时拉索预应力大小以及结构运营时的预应力损失。其中，荷载造成的影响可以通过在有限元分析过程中根据试验时的实际情况进行合理的假定来减小或消除，如根据试验现场的风力情况、气温等环境状况，对有限元分析时的风荷载等条件进行相应的设置；将拉索预应力作为优化参数则可认为，施工时施加预应力及之后的预应力损失同时考虑，优化结果可视为等效的预应力值。
2.1 有限元模型修正方法
如图2所示，索膜结构由施工至运营可以简化为3个阶段：初始状态为施加预应力前结构的状态，仅存在约束条件，结构内力为0；施加预应力后，结构由于受牵拉位形发生变化，结构内力改变；在外荷载(重力荷载、风荷载等)作用下，结构位形与结构内力再次改变。图中，T为各索力组成的向量，X为结构特征点的位形指标组成的向量。为了清楚表明结构位形变化，以虚线为等高线作为参照。
基于对三阶段过程的分析，可以建立相应的有限元模型修正方法，如下所述。
首先选取结构中各拉索的预应力T作为优化参数，在确定结构施加的外荷载之后，便可求解得到结构特征点的位形指标，因此X可看做是关于T的函数，基于结构特征点的位形指标X构建形如式(1)的目标函数：
{\text{min }}F({\boldsymbol{T}}) = \sum {\left| {{X_A}({\boldsymbol{T}})(i) - {X_M}(i)} \right|} (2) 式中：
{X_M}(i) 为试验测量值的第i个分量；{X_A}({\boldsymbol{T}})(i) 为在预应力为T时的分析计算值的第i个分量。在确定优化参数与目标函数值后，需要针对优化问题的特点选取合适的优化算法。首先，索膜结构一般规模较大，拉索数量以及可选取的位形特征点较多；另外，结构内力与结构位形之间具有复杂的关系。可以确定该问题为多元非线性优化问题，对于类似的工程优化问题的求解方法有多种，如牛顿-高斯迭代法[19]、序列二次规划(SQP)[20]等算法。本文提出的优化方法基于SQP算法，该算法基本原理是将复杂的非线性问题转化为较简单的二次规划问题，是目前公认的求解约束非线性优化问题最有效的方法之一。针对以上优化问题，其优化原理如下所述。
对第k次迭代设优化参数，即预应力为
{{\boldsymbol{T}}_k} = ({t_{k1}},{t_{k2}}, \cdot \cdot \cdot ,{t_{kn}}) 。将目标函数在{{\boldsymbol{T}}_k} 处简化为二次函数：F({\boldsymbol{T}}) = \frac{1}{2}{[{\boldsymbol{T}} - {{\boldsymbol{T}}_k}]^{\rm T}}{\boldsymbol{H}}({{\boldsymbol{T}}_k})[{\boldsymbol{T}} - {{\boldsymbol{T}}_k}] + {\boldsymbol{G}}({{\boldsymbol{T}}_k})[{\boldsymbol{T}} - {{\boldsymbol{T}}_k}] (3) 式中，
{\boldsymbol{G}}({{\boldsymbol{T}}_k}) 、{\boldsymbol{H}}({{\boldsymbol{T}}_k}) 分别为目标函数在{{\boldsymbol{T}}_k} 处的梯度向量与海森矩阵。{\boldsymbol{G}}({{\boldsymbol{T}}_k}) = {\left( {\frac{{\partial F}}{{\partial {t_1}}},\frac{{\partial F}}{{\partial {t_2}}}, \cdot \cdot \cdot ,\frac{{\partial F}}{{\partial {t_n}}}} \right)^{\rm T}} (4) {{\boldsymbol{H}}}({{\boldsymbol{T}}_k}) = \left[ {\begin{array}{*{20}{c}} {\dfrac{{{\partial ^2}F}}{{\partial {t_1}\partial {t_1}}}}&{\dfrac{{{\partial ^2}F}}{{\partial {t_1}\partial {t_2}}}}&{...}&{\dfrac{{{\partial ^2}F}}{{\partial {t_1}\partial {t_n}}}} \\ {\dfrac{{{\partial ^2}F}}{{\partial {t_2}\partial {t_1}}}}&{\dfrac{{{\partial ^2}F}}{{\partial {t_2}\partial {t_2}}}}&{...}&{\dfrac{{{\partial ^2}F}}{{\partial {t_2}\partial {t_n}}}} \\ {...}&{...}&{...}&{...} \\ {\dfrac{{{\partial ^2}F}}{{\partial {t_n}\partial {t_1}}}}&{\dfrac{{{\partial ^2}F}}{{\partial {t_n}\partial {t_2}}}}&{...}&{\dfrac{{{\partial ^2}F}}{{\partial {t_n}\partial {t_n}}}} \end{array}} \right] (5) \begin{split} & {\boldsymbol{G}}({{\boldsymbol{T}}}_{k})\approx \left(\frac{F({{\boldsymbol{T}}}_{k}+\Delta {{\boldsymbol{T}}}_{1})-F({{\boldsymbol{T}}}_{k})}{\Delta {{\boldsymbol{T}}}_{1}(1)},\right.\\&\left.\frac{F({{\boldsymbol{T}}}_{k}+\Delta {{\boldsymbol{T}}}_{2})-F({{\boldsymbol{T}}}_{k})}{\Delta {{\boldsymbol{T}}}_{2}(2)}, \cdot \cdot \cdot ,\frac{F({{\boldsymbol{T}}}_{k}+\Delta {{\boldsymbol{T}}}_{n})-F({{\boldsymbol{T}}}_{k})}{\Delta {{\boldsymbol{T}}}_{n}(n)}\right)^{\rm T} \end{split} (6) 其中：
\Delta {{\boldsymbol{T}}_1} = (\Delta {t_1},0,\cdots,0);\;\Delta {{\boldsymbol{T}}_2} = (0,\Delta {t_2},\cdots,0);\; \cdots; \Delta {{\boldsymbol{T}}_n} = (0,0,\cdots,\Delta {t_n}) 。而本优化问题并没有显式的优化函数用以求得梯度向量与海森矩阵，因此需要用到有限差分法(式(6))，进行近似求解。
问题转化为求解式(3)最小值的二次优化问题，且将
{\boldsymbol{T}} - {{\boldsymbol{T}}_k} 作为优化参数，得到优化结果\Delta {\boldsymbol{T}} = {\boldsymbol{T}} - {{\boldsymbol{T}}_k} 后，将\Delta {\boldsymbol{T}} 作为下一步迭代方向，选取合适的步长参数α，计算第k+1次的优化参数：{{\boldsymbol{T}}_{k + 1}} = {{\boldsymbol{T}}_k} + \alpha \Delta {\boldsymbol{T}} (7) 最终，当满足给定的精度要求时，输出当前迭代得到的优化参数作为最优解。
以上方法是进行优化时的迭代过程，但有时并不能直接通过迭代过程得到最优解，往往需要使用其他算法进行预处理。对于优化参数较多的有限元模型，由于可能存在较多的局部最优值，并不能直接使用上述原理进行迭代，否则可能出现不收敛或求解结果为局部最优值的情况。对于类似的工程优化问题，一般需要预先通过全局搜索[21](包括粒子群算法、遗传算法等)的方法得到若干起点或缩小优化参数范围，之后进行局部最优化。
2.2 有限元模型修正程序实现
基于上述原理与方法，本文作者实现了ABAQUS与MATLAB联合仿真优化程序，基于索膜结构的空间位形指标对结构有限元模型进行参数修正，程序结构如图3所示。该程序主要通过MATLAB脚本实现，使用MATLAB中的优化工具箱(optimization toolbox)进行全局搜索与优化迭代过程。每次迭代的优化参数由MATLAB写入ABAQUS的输入文件，目标函数值通过MATLAB调用ABAQUS程序进行数值分析并读取输出文件进而计算得到。
3 试验验证
为了验证上述有限元模型修正方法的可行性，对北京首都国际机场南线收费大棚索膜结构进行了现场试验。该收费大棚于2008年4月16日竣工并投入使用，建筑物平面呈曲壳型，索膜结构水平投影面积3275 m2，桅杆纵向间距为166.14 m，边环索锚座横向间距为53.4 m。马路两边立有4根钢桅杆，前、后为边环索束，每根桅杆顶部交叉下拉7根吊索，拉住边环索，膜顶中间设置1根谷索，谷索与边环索之间由联系索联结，每根桅杆顶部有2根下拉索与外侧的锚座节点相连保持平衡。结构所处的地理位置与现场照片如图4、图5所示。试验当天天气状况良好，故进行有限元分析时可忽略风荷载的影响，试验内容包括使用全站仪测量结构位形、使用激光扫描仪获取结构三维点云模型、对部分拉索进行索力测量及其他结构状况检测项目。
![]()
图 4 首都国际机场南线收费大棚索膜结构地理位置Figure 4. Location of the toll station of south line of Beijing Capital International Airport3.1 量测设备
试验通过全站仪(KON-XSY-211)以及三维扫描仪(PENTAX S-3180)获取结构的位形信息。基于本结构的跨度，全站仪测距误差小于5 mm，结构特征点位置信息的测量误差在毫米级别。三维扫描仪测距、精度等参数如表1所示，现场测量照片如图6所示。
结构索力采用振动法进行测量，使用无线智能加速度计(HCF400-A0)与数据采集仪(HGL400-A1 4G)获取拉索振动信息，进而得到索力。振动法，是基于拉索索力和自振频率之间的关系，通过测量自振频率间接得到索力，可以适当考虑钢索的松弛和弯曲[22]，误差一般小于3%。
表 1 三维扫描仪参数Table 1. 3D scanner parameters项目 参数 型号 PENTAX S-3180 测量距离 187.3 m 最小测距 0.3 m 分辨率 0.1 mm 数据获取速率 ≤1 016 027点/秒 50 m线性误差 ≤1 mm 3.2 试验测量结果
通过全站仪测量的目标点主要是结构环索与谷索上的特征点，如图7所示，由三维激光扫描获取的点云模型如图8所示。提取点云模型中的谷索与环索的坐标信息：首先对密集的点云模型进行5%~10%的采样，可在一定程度上去除噪声，然后针对目标点，对一定球型范围内的点进行记录，求所有记录点的重心位置作为目标点坐标，可大幅减小少量噪声点的影响。其中谷索坐标信息如表2所示，与全站仪测量结果统一坐标系，得到的对比结果如图9所示。通过对比可以发现，北环索与谷索两种测量结果吻合良好，之后将主要采用三维扫描得到的结构三维指标对有限元模型进行优化。
表 2 谷索坐标信息Table 2. Coordinate information of valley cables测点 纵向 X/m 横向 Y/m 高度 H/m Pt1 −9.801 −15.809 3.971 Pt2 −13.238 −16.088 4.959 Pt3 −23.950 −16.031 7.361 Pt4 −34.293 −15.953 9.565 Pt5 −43.739 −16.180 11.454 Pt6 −53.253 −15.950 12.264 Pt7 −61.468 −15.946 13.011 Pt8 −68.883 −15.906 13.432 Pt9 −75.777 −15.937 13.607 Pt10 −82.497 −16.000 13.341 Pt11 −89.819 −16.083 12.947 Pt12 −98.199 −16.248 12.280 Pt13 −107.559 −16.015 11.171 Pt14 −117.057 −15.926 9.474 Pt15 −127.271 −15.992 7.230 Pt16 −138.130 −16.008 4.688 Pt17 −142.312 −16.282 4.085 索力检测，无论是在膜结构的建设过程中还是在其日常维护检测中，都具有举足轻重的地位。索力是否处在合理的范围内，将直接影响结构的整体受力状态和线形的平顺程度。试验使用振动法对结构索力进行测量，结果见表3。
表 3 索力测量结果Table 3. Cable force measurement results拉索编号 索力/kN 截面积/m2 应力/MPa LS21 702.5 0.003 91 179.668 LS22 639.0 0.003 91 163.427 LS31 5538.5 0.014 59 379.609 LS32 6161.0 0.014 59 422.276 LS41 5076.5 0.014 59 347.944 LS42 2917.0 0.014 59 199.931 LS51 411.5 0.003 91 105.243 LS52 686.0 0.003 91 175.448 3.3 有限元模型修正
3.3.1 建立有限元模型
考虑到通用有限元软件ABAQUS强大的非线性分析功能，有限元建模使用了ABAQUS软件，模型如图10所示。整个模型由梁单元、桁架单元、膜单元构成：两侧桅杆采用了梁单元进行建模，每个杆件划分为1个单元；所有的拉索采用T3D2的桁架单元进行模拟，每根拉索简化划分为1个单元；膜采用了壳分类中的膜单元进行建模，与拉索的桁架单元节点协调，每个三角形膜结构使用1个三角形单元。桅杆与接地的拉索均与地面刚接。实际结构建成时间较长，无法准确获取结构的老化损伤信息，按照设计时的材料属性进行定义。
通过温差法对拉索施加预应力，因此在之后的优化过程中，选取温度作为优化参数，本质上还是将拉索预应力作为优化参数，结合现场试验时无风，不计风荷载，仅考虑重力荷载。
3.3.2 目标函数
由于特征点数目较多，选取部分具有代表性的点构建目标函数，同时未使用到的特征点可以用来对修正后的有限元模型进行评估。分别选取两条环索与一条谷索上各3个特征点，如图11所示，通过其位置信息建立目标函数，该目标函数反映了有限元计算结果与实测结果的吻合程度：
\min {\text{ }}F = \sum {\left| {\frac{{{X_A}(i) - {X_M}(i)}}{{{X_M}(i)-10}}} \right|} (8) 式中：
{X_M} 为试验测量得到的高度坐标；{X_A} 为有限元计算得到的高度坐标，i = 1,2,\cdots,9 。3.3.3 修正流程
1) 启动优化算法。输入试验测量结果并设定优化参数初始值与可行域，设置优化算法相关参数，如收敛限值、有限差分最小步长、最大迭代次数等。
2) ABAQUS运行。通过MATLAB脚本将优化参数写入输入文件(INP文件)并调用ABAQUS执行输入文件，输出结果文件(ODB文件)。
3) 求解目标函数值。MATLAB运行python脚本读取二进制结果文件，得到相应计算值，与实测结果进行计算得到目标函数值。
4) 优化处理。利用目标函数值进行相关运算，计算修正梯度向量、海森矩阵，得到下一个优化参数继续步骤2)循环计算，直至满足收敛限值或迭代次数限值条件。
首次局部优化完成后，程序执行全局搜索算法，基于优化参数可行域，生成一系列均匀分布的优化参数，并对其进行计算筛选出部分优质的(目标函数值相对较小)参数作为迭代起点进行局部优化，最终得到全局最优解。在这一过程中，可以通过观察全局搜索过程结果，人工选取迭代起点并缩小优化参数可行域，提高收敛速度。
3.3.4 全局搜索
由于该索膜结构较为复杂、优化参数较多，可能存在多个局部最优解，因此初始值的选取对优化结果影响重大。可通过全局搜索得到多个优质的迭代起点，同时可根据全局搜索结果手动缩小优化参数的可行域，提高迭代过程的收敛速度。全局搜索过程中，目标函数的变化如图12所示，根据变化情况可选取部分结果较好的点进行下一步迭代计算。
3.3.5 局部优化
基于以上全局搜索的结果，缩小优化参数的可行域，选取优质点作为迭代起点，经过300次左右的函数计算(主要是用于求解梯度向量、修正海森矩阵)，完成27次迭代，最终达到预先设定的迭代精度，得到最优解，过程如图13所示。可以发现，经过全局搜索后得到的迭代起点，在迭代过程中得到进一步优化，使得目标函数值大幅下降。目标函数的下降可以说明，前述结构中的特征点的位形信息与试验测量结果更加吻合。具体情况可以通过图14的对比看出，修正后的变形与实测变形相差很小。在施加温差应力模拟张拉后，再对结构施加外部荷载，拉索上的特征点经历变形增大、变形减小两个过程，实测以及修正后有限元的结果中特征点1变形较小，因此在差值较小的情况下，计算偏差百分比数值较大，未在偏差对比图中绘出。其他特征点的修正后偏差均远小于修正前偏差，有限元模型优化工作已经基本完成。
3.3.6 模型修正结果的验证
对修正后的有限元模型进行验证，判断模型与实际结构的吻合程度，Friswell和Mottershead[23]建议将实测结果一部分用于模型修正，另一部分用于模型评估。根据试验测量得到的索力数据，对修正后的模型进行独立验证，同时也是对本文提出的有限元模型修正方法的可行性验证。试验时共测得了8根拉索的索力值，拉索位置如图15所示，通过实测索力值计算拉索应力，与修正前后有限元模型相应构件进行对比，结果如图16所示。与模型修正前相比，修正后的应力计算值与实测结果偏差普遍大幅减小，其中LS21、LS22、LS31、LS42、LS51效果明显，有限元模型的索力由修正前最大偏差10%~30%降低至10%以下。
另外，未参与模型修正的特征点也可以用来对模型进行评估，选取部分特征点对修正后的模型进行验证，如图17所示。可以发现，与修正前相比，修正后的模型变形与实测值更加接近。从应力、应变两个方面来看，修正后模型与实际结构更加接近。因此本文提出的有限元模型修正方法可显著优化有限元模型，使其反映结构真实的服役性态，从而进一步对结构进行性能评价。
3.4 结构安全性评估
基于修正后的有限元模型，对结构进行承载力分析，在最不利的荷载布置下，得到结构的应力分布如图18所示，膜面变形云图如图19所示，由于膜结构应力较小且会遮挡拉索的应力显示，因此应力分布图对膜单元进行了消隐处理。图中显示，应力最大的拉索处于大约1/3跨连接边环索与谷索的位置，为1035 MPa，小于抗拉强度1770 MPa，因此该结构各拉索构件的安全性满足要求。
![]()
图 19 膜面变形云图(竖向挠度) /mFigure 19. Membrane surface deformation cloud map (vertical displacement)4 结论
本文提出一种基于空间位形的在役索膜结构模型修正方法，通过现场试验验证了该方法的可行性，采用修正后的模型实现结构的安全评估，主要结论如下：
(1) 以结构位形构建目标函数，以拉索预应力为修正参数，通过全局搜索与局部优化可实现基于空间位形的在役柔性索膜结构的有限元模型修正。
(2) 利用现场试验进行方法验证，经过迭代修正，有限元模型的空间位形信息与实测数据吻合良好，更准确地反映了结构真实的服役性态。利用实测索力值进行独立验证，结果表明方法有效可靠且精度较高。
(3) 进行索膜结构有限元模型修正，由于优化参数较多且结构非线性特征明显，可能会求解得到局部最优值，可以通过全局搜索筛选迭代起点、缩减优化参数可行域，进而提高优化效率与准确度。
注：该文在第30届结构工程学术会议(2021广州)应邀作特邀报告
-
![]()
图 1 随机残余误差关于辅助参数h值变化曲线
Figure 1. Stochastic residual error with respect to the auxiliary parameter h
![]()
图 3 变截面柱屈曲特征值前4阶统计矩(工况1)
Figure 3. First four order statistical moments of the buckling eigenvalue for Case 1
![]()
图 4 变截面柱屈曲特征值前4阶统计矩(工况2)
Figure 4. First four order statistical moments of the buckling eigenvalue in Case 2
![]()
图 5 辅助参数h值和趋近函数随δ变化关系图(工况2)
Figure 5. Values of the auxiliary parameter h and the approaching function in the homotopy series expansion for Case 2
![]()
图 7 一阶屈曲模态前4阶统计矩
Figure 7. First four order statistical moments of the 1st-order buckling mode
表 1 同伦级数中需要计算系数的个数
Table 1 Number of terms in homotopy series
随机变量个数 2阶展开 3阶展开 4阶展开 一维 二维 完整 一维 二维 完整 一维 二维 完整 1 3 3 3 4 4 4 5 5 5 5 11 21 21 16 46 56 21 81 126 10 21 66 66 31 166 286 41 311 1001 15 31 136 136 46 361 816 61 691 3876 20 41 231 231 61 631 1771 81 1221 10 626 
下载: 导出CSV
表 2 不同方法计算Y的前4阶统计矩
Table 2 First four order statistical moments of Y calculated by different methods
统计矩 蒙特卡洛 泰勒级数 PC法 HSFEM方法 本文方法 1σ 2σ 3σ 均值 0.43 0.42 0.43 0.43 0.43 0.43 0.43 均方差 0.07 0.06 0.07 0.06 0.07 0.07 0.07 偏度 1.63 −0.82 1.88 −0.03 2.57 1.34 1.56 峰度 5.58 71.98 6.68 40.42 12.35 5.36 5.12
下载: 导出CSV
表 3 框架结构屈曲特征值的前4阶统计矩( {\delta _{{E_{\rm{c}}}}} =0.3)
Table 3 First four order statistical moments of the buckling eigenvalue ( {\delta _{{E_{\rm{c}}}}} =0.3)
统计矩 蒙特卡洛 2阶展开 3阶展开 4阶展开 摄动法 HSFEM 本文方法-1 本文方法-2 摄动法 HSFEM 本文方法-1 本文方法-2 摄动法 HSFEM 本文方法-1 本文方法-2 均值 16.22 16.22 16.21 16.22 16.22 16.22 16.21 16.22 16.22 16.22 16.21 16.22 16.22 均方差 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 0.05 偏度 1.47 1.34 1.41 1.49 1.49 1.58 1.51 1.49 1.49 1.44 1.46 1.47 1.47 峰度 7.31 6.16 6.62 7.22 7.23 8.90 7.92 7.48 7.48 6.85 7.06 7.26 7.26
下载: 导出CSV
表 4 框架结构屈曲特征值的前4阶统计矩( {\delta _{{E_{\rm{c}}}}} =0.6)
Table 4 First four order statistical moments of the buckling eigenvalue ( {\delta _{{E_{\rm{c}}}}} =0.6)
统计矩 蒙特卡洛 2阶展开 3阶展开 4阶展开 摄动法 HSFEM 本文方法-1 本文方法-2 摄动法 HSFEM 本文方法-1 本文方法-2 摄动法 HSFEM 本文方法-1 本文方法-2 均值 16.34 16.33 16.32 16.32 16.32 16.34 16.33 16.33 16.33 16.34 16.33 16.33 16.33 均方差 0.09 0.09 0.09 0.09 0.09 0.10 0.09 0.09 0.09 0.09 0.09 0.09 0.09 偏度 1.28 1.00 1.11 1.31 1.31 2.67 1.91 1.34 1.34 −0.55 1.04 1.29 1.29 峰度 6.02 4.68 4.94 5.95 5.97 53.29 21.56 6.38 6.39 101.88 8.51 6.04 6.05
下载: 导出CSV
表 5 各方法计算时间
Table 5 Computation time for different methods
/s 蒙特卡洛 2阶展开 3阶展开 4阶展开 摄动法 本文方法-1 本文方法-2 摄动法 本文方法-1 本文方法-2 摄动法 本文方法-1 本文方法-2 710730 71 64 65 106 83 87 182 94 131
下载: 导出CSV
表 6 框架结构屈曲特征值的前4阶统计矩( \delta_{E_{\rm{c}}} =0.6, δb=0.3)
Table 6 First four order statistical moments of the buckling eigenvalue ( \delta_{E_{\rm{c}}} =0.6, δb=0.3)
统计矩 蒙特
卡洛2阶展开 3阶展开 4阶展开 本文方
法-1本文方
法-2本文方
法-1本文方
法-2本文方
法-1本文方
法-2均值 15.97 15.98 15.95 15.98 15.95 15.98 15.96 均方差 0.69 0.71 0.70 0.73 0.69 0.73 0.69 偏度 −0.27 −0.53 −0.24 −0.55 −0.26 −0.55 −0.26 峰度 2.72 3.26 2.74 3.30 2.71 3.31 2.71
下载: 导出CSV
-
[1] GODOY L A, ELISHAKOFF I. The experimental contribution of petrus van musschenbroek to the discovery of a buckling formula in the early 18th century [J]. International Journal of Structural Stability and Dynamics, 2020, 20(5): 2050063. doi: 10.1142/S0219455420500637
[2] 范浩, 王新, 董卫国, 等. 基于屈曲模态的弹性压弯构件平面内二阶弯矩计算[J]. 工程力学, 2022, 39(2): 37 − 50. doi: 10.6052/j.issn.1000-4750.2020.12.0900 FAN Hao, WANG Xin, DONG Weiguo, et al. Calculation of second-order in-plane bending moment of elastic compression-bending members based on buckling modes [J]. Engineering Mechanics, 2022, 39(2): 37 − 50. (in Chinese) doi: 10.6052/j.issn.1000-4750.2020.12.0900
[3] 钮鹏, 李旭, 李世荣, 等. 弹性地基上复合材料夹层梁的热过屈曲[J]. 工程力学, 2017, 34(增刊): 26 − 30. doi: 10.6052/j.issn.1000-4750.2016.04.S023 NIU Peng, LI Xu, LI Shirong, et al. The thermal buckling of composite sandwich beams on elastic foundations [J]. Engineering Mechanics, 2017, 34(Suppl): 26 − 30. (in Chinese) doi: 10.6052/j.issn.1000-4750.2016.04.S023
[4] AARON P G, CARLOS V P, MAURICIO A H. Assessing the influence of side friction and intermediate end-boundary conditions on the critical buckling load in piles using the differential transform method [J]. Engineering Structures, 2021, 239: 112269. doi: 10.1016/j.engstruct.2021.112269
[5] ZHENG L, DENG T, LIU Q J. Buckling of piles in layered soils by transfer matrix method [J]. International Journal of Structural Stability and Dynamics, 2021, 21(8): 2150109. doi: 10.1142/S0219455421501091
[6] 赵桂平, 赵锺斗. 基于改进弧长法的层压复合壳后屈曲反应分析[J]. 工程力学, 2003, 20(1): 106 − 111. doi: 10.3969/j.issn.1000-4750.2003.01.021 ZHAO Guiping, ZHAO Chongdou. Post-buckling analysis of laminated composite shell based on an improved arc-length method [J]. Engineering Mechanics, 2003, 20(1): 106 − 111. (in Chinese) doi: 10.3969/j.issn.1000-4750.2003.01.021
[7] ANSARI Q M, TRINH L C, ZUCCO G, et al. Effect of elastic support on the linear buckling response of quasi-isotropic cylindrical shells under axial compression [J]. Engineering Structures, 2021, 244: 112796. doi: 10.1016/j.engstruct.2021.112796
[8] ZAFARABADI M M, AGHDAM M M. Semi-analytical solutions for buckling and free vibration of composite anisogrid lattice cylindrical panels [J]. Composite Structures, 2021, 275: 114422. doi: 10.1016/j.compstruct.2021.114422
[9] 李振眠, 余杨, 余建星, 等. 基于向量有限元的深水管道屈曲行为分析[J]. 工程力学, 2021, 38(4): 247 − 256. doi: 10.6052/j.issn.1000-4750.2020.06.0357 LI Zhenmian, YU Yang, YU Jianxing, et al. Buckling analysis of deep water pipelines by vector form intrinsic finite element method [J]. Engineering Mechanics, 2021, 38(4): 247 − 256. (in Chinese) doi: 10.6052/j.issn.1000-4750.2020.06.0357
[10] 余杨, 李振眠, 余建星, 等. 穿越平移断层海底埋地管道屈曲失效分析[J]. 工程力学, 2022, 39(9): 242 − 256. doi: 10.6052/j.issn.1000-4750.2021.05.0391 YU Yang, LI Zhenmian, YU Jianxing, et al. Buckling failure analysis of subsea buried pipeline crossing strike-slip fault [J]. Engineering Mechanics, 2022, 39(9): 242 − 256. (in Chinese) doi: 10.6052/j.issn.1000-4750.2021.05.0391
[11] ELISHAKOFF I. Probabilistic methods in the theory of structures: Strength of materials, random vibrations, and random buckling [M]. Singapore: World Scientific, 2017.
[12] ROMES A B, LUIZ F F, THIAGO P S, et al. Stochastic eigenfrequency and buckling analyses of plates subjected to random temperature distributions [J]. Mechanical Systems and Signal Processing, 2021, 147: 107088. doi: 10.1016/j.ymssp.2020.107088
[13] SHINOZUKA M, ASTILL C J. Random Eigenvalue Problems in Structural Analysis [J]. AIAA Journal, 1972, 10(4): 456 − 462. doi: 10.2514/3.50119
[14] SONG D, ELLINGWOOD B R, COX J V. Solution methods and initialization techniques in SFE analysis of structural stability [J]. Probabilistic Engineering Mechanics, 2005, 20(2): 179 − 187. doi: 10.1016/j.probengmech.2005.05.002
[15] CHEN L, LIU S W, ZHANG Z J, et al. Efficient algorithm for elastic buckling of corroded I-section steel members with Monte Carlo simulation [J]. Thin-Walled Structures, 2022, 175: 109216. doi: 10.1016/j.tws.2022.109216
[16] JEONG G D. Critical buckling load statistics of an uncertain column [C]// Probabilistic Mechanics and Structural and Geotechnical Reliability. Proceedings of the Sixth Specialty Conference. Denver, Colorado, ASCE, 1992: 563 − 566.
[17] ZHANG J, ELLINGWOOD B. Effects of uncertain material properties on structural stability [J]. Journal of Structural Engineering-ASCE, 1995, 121(4): 705 − 716. doi: 10.1061/(ASCE)0733-9445(1995)121:4(705)
[18] RAMU S A, GANESAN R. Response and stability of a stochastic beam-column using stochastic FEM [J]. Computers & Structures, 1995, 54(2): 207 − 221.
[19] GRAHAM L L, SIRAGY E F. Stochastic finite-element analysis for elastic buckling of stiffened panels [J]. Journal of Structural Engineering-ASCE, 2001, 127(1): 91 − 97.
[20] ALTUS E, TOTRY E M. Buckling of stochastically heterogeneous beams, using a functional perturbation method [J]. International Journal of Solids & Structures, 2003, 40(23): 6547 − 6565.
[21] GUPTA A, ARUN C O. Stochastic meshfree method for elastic buckling analysis of columns [J]. Computers & Structures, 2018, 194: 32 − 47.
[22] KAMIŃSKI M, ŚWITA P. Generalized stochastic finite element method in elastic stability problems [J]. Computers & Structures, 2011, 89(11): 1241 − 1252.
[23] WU D, GAO W, TANGARAMVONG S, et al. Robust stability analysis of structures with uncertain parameters using mathematical programming approach [J]. International Journal for Numerical Methods in Engineering, 2015, 100(10): 720 − 745.
[24] WU D, GAO W, SONG C, et al. Probabilistic interval stability assessment for structures with mixed uncertainty [J]. Structural Safety, 2016, 58: 105 − 118. doi: 10.1016/j.strusafe.2015.09.003
[25] LI K, WU D, GAO W. Spectral stochastic isogeometric analysis for linear stability analysis of plate [J]. Computer Methods in Applied Mechanics and Engineering, 2019, 352: 1 − 31.
[26] ALIBRANDI U, IMPOLLONIA N, RICCIARDI G. Probabilistic eigenvalue buckling analysis solved through the ratio of polynomial response surface [J]. Computer Methods in Applied Mechanics and Engineering, 2010, 199(9): 450 − 464.
[27] 徐亚洲, 白国良. 考虑混凝土材料变异性的超大型冷却塔随机屈曲承载力分析[J]. 工程力学, 2012, 29(8): 208 − 212. doi: 10.6052/j.issn.1000-4750.2010.11.0853 XU Yazhou, BAI Guoliang. Analysis of random buckling bearing capacity of super-large cooling towers considering the variation of the property of concrete [J]. Engineering Mechanics, 2012, 29(8): 208 − 212. (in Chinese) doi: 10.6052/j.issn.1000-4750.2010.11.0853
[28] SU C, FAN X, MA H, et al. Green’s Function Method for Stability Analysis of Stochastic Structures [J]. Journal of Engineering Mechanics-ASCE, 2015, 141(3): 04014121-1 − 04014121-11. doi: 10.1061/(ASCE)EM.1943-7889.0000842
[29] ZHANG H, HUANG B. A new homotopy-based approach for structural stochastic analysis [J]. Probabilistic Engineering Mechanics, 2019, 55: 42 − 53. doi: 10.1016/j.probengmech.2018.09.001
[30] 张衡, 王鑫, 陈辉, 等. 基于同伦分析方法的随机结构静力响应求解[J]. 工程力学, 2019, 36(11): 27 − 33, 61. doi: 10.6052/j.issn.1000-4750.2018.11.0619 ZHANG Heng, WANG Xin, CHEN Hui, et al. Solution for static response of structure with random parameters based on homotopy analysis method [J]. Engineering Mechanics, 2019, 36(11): 27 − 33, 61. (in Chinese) doi: 10.6052/j.issn.1000-4750.2018.11.0619
[31] HUANG B, ZHANG H, PHOON K K. Homotopy approach for random eigenvalue problem [J]. International Journal for Numerical Methods in Engineering, 2018, 113(3): 450 − 478. doi: 10.1002/nme.5622
[32] 廖世俊. 超越摄动: 同伦分析方法基本思想及其应用[J]. 力学进展, 2008, 38(1): 1 − 34. doi: 10.3321/j.issn:1000-0992.2008.01.001 LIAO Shijun. Beyond perturbation: The basic concepts of the homotopy analysis method and its applications [J]. Advances in Mechanics, 2008, 38(1): 1 − 34. (in Chinese) doi: 10.3321/j.issn:1000-0992.2008.01.001
[33] LIAO S J. An optimal homotopy-analysis approach for strongly nonlinear differential equations [J]. Communications in Nonlinear Science and Numerical Simulation, 2010, 15: 2003 − 2016. doi: 10.1016/j.cnsns.2009.09.002
[34] 黄斌, 贺志赟, 张衡. 随机桁架结构几何非线性问题的混合摄动-伽辽金法求解[J]. 力学学报, 2019, 51(5): 1424 − 1436. doi: 10.6052/0459-1879-19-099 HUANG Bin, HE Zhiyun, ZHANG Heng. Hybrid perturbation-Galerkin method for geometrical nonlinear analysis of truss structures with random parameters [J]. Chinese Journal of Theoretical and Applied Mechanics, 2019, 51(5): 1424 − 1436. (in Chinese) doi: 10.6052/0459-1879-19-099
[35] 廖世俊, 刘曾. 同伦分析方法进展综述[J]. 力学进展, 2019, 49: 201902. doi: 10.6052/1000-0992-18-005 LIAO Shijun, LIU Zeng. A brief review of the homotopy analysis method [J]. Advances in Mechanics, 2019, 49: 201902. (in Chinese) doi: 10.6052/1000-0992-18-005
[36] NI P H, XIA Y, LI J, et al. Multi-scale stochastic dynamic response analysis of offshore risers with lognormal uncertainties [J]. Ocean Engineering, 2019, 189: 106333. doi: 10.1016/j.oceaneng.2019.106333
-
期刊类型引用(13)
1. 林楷奇，郑俊浩，陆新征. 数字孪生技术在土木工程中的应用：综述与展望. 哈尔滨工业大学学报. 2024(01): 1-16 . 
百度学术
2. 张帆，孙楚津，覃思中，赵馨怡. 基于BIM和深度学习点云分割的施工检查方法模拟研究. 工程力学. 2024(02): 194-201 . 本站查看
3. 秦志浩，李志强，李永辉. 城镇老旧小区建筑结构检测及安全性评价分析. 科技创新与生产力. 2024(03): 57-60 . 百度学术
4. 宋华霖，姜绍飞. 数字孪生在空间结构生命周期管理中的应用与挑战. 福州大学学报(自然科学版). 2024(03): 369-378 . 百度学术
5. 卢耀光. 大空间张拉索膜结构分析与优化. 新疆有色金属. 2024(05): 104-105 . 百度学术
6. 刘慕广，乔磊，王雷，余先锋，张春生，谢壮宁，张丽. 桅杆结构有限元模型修正与参数敏感性研究. 振动工程学报. 2024(07): 1107-1114 . 百度学术
7. 王超，孙小刚，张李军，罗晓群. 寒冷地区单层网格支承膜结构天线罩工程检测研究. 建筑结构. 2024(21): 110-116 . 百度学术
8. 王磊，黄杜康，马亚飞，黄可. 基于元学习的结构损伤定位与量化. 中国公路学报. 2024(11): 38-51 . 百度学术
9. 梁龙. 桥隧串接洞门施工结构稳定性分析研究. 铁道建筑技术. 2023(04): 151-154 . 百度学术
10. 申跃奎，李洋，朱新雨，刘晨迪，李星. 不同索网加劲的气承式膜结构受力特性研究. 计算力学学报. 2023(06): 1029-1035 . 百度学术
11. 霍林生，丁宁，李庆伟，李伯遥. 基于空间位形的FAST索网有限元模型修正方法研究. 建筑结构. 2023(S2): 557-561 . 百度学术
12. 刘宇飞，齐玉，李保罗，冯楚乔，丁一凡，聂鑫，樊健生. 工程结构运维多场景的敏感服役指标智能识别. 建筑结构学报. 2022(10): 1-15 . 百度学术
13. 曹辉，许涛，申金刚，王浩，穆明. 钢骨架索膜结构储煤棚工程设计研究. 煤炭工程. 2022(12): 24-27 . 百度学术
其他类型引用(6)
计量
- 文章访问数: 352
- HTML全文浏览量: 69
- PDF下载量: 63
- 被引次数: 19
