螺栓法兰连接结构有限元模型参数确定方法

蒋国庆，马 斌，陈万华

(中国空气动力研究与发展中心， 四川 绵阳 621000)

螺栓法兰连接结构广泛应用于各种工业装备中，起着重要的连接作用。该结构的特点是构造简单、维护简便、成本低廉、可操作性好等，但由于其固有的结构不连续性和几何突变性，在外力作用下结构可能出现应力分布不均匀甚至局部应力过大的现象。为对螺栓法兰连接结构进行给定载荷作用下的强度校核和刚度分析，需要得到高精度的有限元模型。

Gray等[1]利用壳单元模拟螺栓连接板件和螺栓，McCarthy等[2]利用三维实体单元对螺栓以及板件进行了精细网格划分，这些模型的准确性较高，有时甚至可以代替试验研究。Pedersen等[3]从弹性能量的角度来确定螺栓连接的连接刚度，Kim等[4]建立了实体螺栓模型、耦接螺栓模型、无螺栓模型、蜘蛛型螺栓模型四种有限元模型并定义了这些模型的连接刚度。王建民等[5]针对工程上建立的螺栓对接结构线性模型无法反映时变非线性特征这一不足，建立了考虑接触和预紧力等因素的非线性解析模型，其分析效率得到显著提高。杨敏[6]通过实验研究建立了剪切式螺栓连接件的唯象等效模型并将其推广至典型螺栓连接结构的简化建模。张文元等[7]通过对螺栓连接进行弹性、滑移和屈服等阶段的受力机理研究，得到了用于代替螺栓的连接件的本构关系，并成功运用于钢结构高强度螺栓连接的简化。江金锋等[8]应用Global/Local技术建立了螺栓连接结构的有限元模型并进行了带EWK断裂子程序的非线性分析，分析效率得到了明显提升。张红艳等[9]建立了二维螺栓连接桁架的完全固连模型、螺栓固连模型和非线性摩擦接触模型，结果表明螺栓固连模型比完全固连模型更接近非线性接触模型。张琪昌等[10]利用改进后的正交模型正交模态法对某类螺栓连接结构的刚性有限元模型进行了修正，修正后模型的计算精度和可靠度显著提高。陈宏威等[11]建立了组合结构的梁-双弹簧等效计算模型，分析了不同材料垫片、螺栓预紧力、几何尺寸等对轴向刚度非线性特性的影响。上述文献几乎涵盖了螺栓连接研究的方方面面，但具体到每篇文献而言，都或多或少存在一定不足，例如有些单纯只开展了数值仿真研究而没有进行试验验证、有些考虑的因素不够全面。本文从有限元网格参数和接触参数等出发开展螺栓法兰连接结构有限元模型参数确定方法研究，综合采用数值仿真手段和试验手段来获得高精度的有限元模型。

1 螺栓法兰连接结构简介

图1为典型螺栓法兰连接结构及其主要几何尺寸(单位均为mm)，由图可知该结构主要由上、下部段和若干螺栓螺母组件组成，其中部段材料为铝合金LY12，螺栓螺母材料为45号钢。螺栓法兰连接结构所用材料的主要力学性能参数如表1所示。

图1 典型螺栓法兰连接结构Fig.1 Typical bolted flange joint

表1 材料主要力学性能参数

图1中螺栓法兰连接结构几何模型与实际结构高度接近，因而直接对其进行网格划分而得到的有限元模型的计算精度主要取决于模型的网格参数和接触参数。下面的研究主要围绕这些参数的确定而展开。

2 网格参数确定

螺栓法兰连接结构主要由螺栓和部段组成，为使模型的网格数目最少，需要分别确定各个部位的网格参数。网格参数的确定步骤如下：

1)根据螺栓法兰连接结构不同部件的结构形式设计合理的计算模型。

2)以单元边长为变量设置不同的种子数，分析不同单元边长下螺栓法兰连接结构不同部件有限元模型计算结果的差异性。当相邻边长所对应的模型的计算结果相差接近1%时，可以认为该部件的网格参数已经达到合理水平。

2.1 螺栓网格参数确定

为确定螺栓网格参数，设计如图2所示计算模型，该模型的边界载荷条件设置为底端固支、顶端施加垂向的拉力(2 kN)或者压力(4 kN)。图3为螺栓单元边长对螺栓最大轴向位移的影响规律。

图2 螺栓网格参数确定计算模型Fig.2 Model for getting mesh parameter of the bolt

(a) 拉力(a) Tensile force

(b) 压力(b) Pressure force图3 单元边长对螺栓最大轴向位移的影响Fig.3 Element length′s influence on the bolt′s max axial displacement

由图3可知，不管外载荷是拉力还是压力，随着单元边长的减小，螺栓有限元模型的计算精度都越来越高。单元边长为0.6 mm和0.5 mm时，所对应模型的计算结果最大相差1.2%，而后者单元规模是前者单元规模的1.92倍。可见，0.6 mm是一个比较理想的单元边长。

2.2 部段网格参数确定

确定部段的网格参数时，可直接选用图1中的结构作为计算模型。对该模型进行网格划分前，考虑到部段的竖直部位在部段中占有更大的体积百分比，划分网格时为避免沙漏现象，需要在壁厚方向至少划分三层单元，这样部段的单元规模会变得异常庞大。因此可以在法兰面对部段进行分割，分割后的竖直部分采用壳单元建模，而水平部分依然采用体单元建模，如图4所示。计算模型的边界载荷条件为底端固支、顶端施加50 kN的拉力或者100 kN的压力，提取的响应为顶端最大轴向位移。

(a) 竖直部位(a) Vertical part

(b) 水平部位(b) Horizontal part图4 部段网格参数确定计算模型Fig.4 Model for getting mesh parameter of the section

图5为单元边长对部段竖直部位最大轴向位移的影响规律。由图可知，不管外载荷是拉力还是压力，竖直部位的单元边长对最大轴向位移的影响都较小。单元边长为20 mm和单元边长为5 mm所对应模型的计算结果最大相差1.2%，而后者的单元规模是前者单元规模的14.5倍。因此，对于部段的竖直部位而言，单元边长20 mm是一个比较合理的网格参数。

(a) 拉力(a) Tensile force

(b) 压力(b) Pressure force图5 单元边长对竖直部位最大轴向位移的影响Fig.5 Element length′s influence on the vertical section′s max axial displacement

图6为单元边长对部段水平部位最大轴向位移的影响规律。由图可知，不管外载荷是压力还是拉力，水平部位的单元边长对最大轴向位移的影响都较大。这一方面是因为单元边长会影响水平部位与竖直部位的壳体耦合面积，另一方面是因为水平部位存在圆孔等特殊结构，当单元边长较大时，在圆孔附近会出现较多的畸形单元。单元边长为4 mm和3 mm时所对应模型的计算结果相差在1.0%以内。同时，单元边长为4 mm时，水平部位的厚度方向共有3层单元，这是弯矩载荷作用下避免出现沙漏现象的最少单元层数，因而厚度方向已无降低单元规模的可能。因此，对于部段的水平部位而言，单元边长4 mm是一个比较理想的网格参数。

(a) 拉力(a) Tensile force

(b) 压力(b) Pressure force图6 单元边长对水平部位最大轴向位移的影响Fig.6 Element length′s influence on the horizontal section′s max axial displacement

3 接触参数确定

接触参数是螺栓法兰连接结构有限元模型中的关键参数，不仅决定着计算结果的准确性，还影响着模型的收敛性。接触参数主要包括接触刚度和摩擦系数，其中接触刚度一般采用系统默认的硬接触模型计算求得，摩擦系数一般根据工程经验选取，但是有时难以获得合理的工程经验数值。铝合金相对偏软，铝合金与铝合金之间的摩擦系数与多种因素有关，如载荷、加工精度等。本文中的摩擦主要来源于铝合金与铝合金之间，不宜直接选用某一特定情况下获得的工程经验数值。下面采用数值仿真与试验相结合的方法来确定摩擦系数。

图7 螺栓法兰连接结构静力学试验Fig.7 Static test for the bolted flange joint

图7为螺栓法兰连接结构静力学试验示意图，主要由万能试验机、电子应变引伸计、结构多点变形测量仪和结构多点数据采集仪等组成。其中电子应变引伸计布置如图8所示，图中圆点表示该位置有螺栓，而方框则为电子应变引伸计。静力学试验包括拉伸试验和压缩试验，其中拉伸试验的最大拉力为50 kN，压缩试验的最大压力为100 kN。为减小误差，分别进行了4次拉伸试验和4次压缩试验。对各次试验的数据进行相应的均化处理后统计如图9所示，图中螺栓处为①、③、⑤、⑦，其余为相邻螺栓中间处。

图8 电子应变引伸计布置示意Fig.8 Distribution diagram of electronic strain extensometer

(a) 拉伸试验(a) Tensile test

(b) 压缩试验(b) Pressure test图9 电子应变引伸计测量数据Fig.9 Data from electronic strain extensometer

为确定模型中的摩擦系数，首先要分析摩擦系数对计算结果的影响规律，如图10所示。图10中测量部位与图7中电子应变引伸计两刀刃之间的部位相对应。

(a) 50 kN拉力(a) Tensile force is equal to 50 kN

(b) 100 kN压力(b) Pressure force is equal to 100 kN图10 摩擦系数对测量部位轴向位移的影响Fig.10 Friction coefficient′s influence on axial displacement of the tested part

由图10可知，随着摩擦系数μ的增加，连接结构的抗拉刚度随之增加，但增加幅度越来越小，而抗压刚度基本保持不变。采用最小二乘法拟合可知，当摩擦系数μ取0.16时，模型计算精度较高，如图11和图12所示。由图11和图12可知，50 kN拉力作用下，测量部位在螺栓处和相邻螺栓中间处模型计算与试验测量的轴向位移误差率分别为6.41%、-5.65%；100 kN压力作用下，测量部位在螺栓处和相邻螺栓中间处模型计算与试验测量的轴向位移误差率分别为-1.88%、1.59%。

(a) 螺栓处(a) The place where bolts are located

(b) 相邻螺栓中间处(b) The place which locate at the middleposition of two adjacent bolts图11 μ为0.16时拉力作用下模型计算结果与试验测量结果对比曲线Fig.11 Comparing curves between computing data and testing data under tensile force when μ=0.16

(a) 螺栓处(a) The place where bolts are located

(b) 相邻螺栓中间处(b) The place which locate at the middle position of two adjacent bolts图12 μ为0.16时压力作用下模型计算结果与试验测量结果对比曲线Fig.12 Comparing curves between computing data and testing data under compression force when μ=0.16

4 结论

综合采用数值仿真手段和试验手段研究了螺栓法兰连接结构有限元模型中网格参数和接触参数的确定方法，得到的主要结论如下：

1)确定结构的网格参数时可以遵循从整体到局部的思想，在保证计算精度的前提下最大限度地降低单元规模，进而可为模型的后续使用提供便利；

2)选取结构的摩擦系数时，对于铝合金而言，不宜直接选用某一特定情况下获得的工程经验数值，可以采用数值仿真与试验相结合的方法予以确定；

3)网格参数和接触参数确定后的螺栓法兰连接结构有限元模型不仅计算效率得到较大幅度提高，且计算精度能够保持在合理水平。