非线性索单元（共4篇）
非线性索单元 篇1
0 引言
缆索单元在工程领域有着广泛的应用,如桥梁拉索、起重拉索、运输索道等。缆索的运动方程为强非线性的偏微分方程,直接求解需要很大的计算量。很多文献用有限元法来求解绳索在外力作用下偏离平衡位置时的变形应力及动态响应[1,2,3,4]。应用有限元法时,必须有足够的约束条件,且不能发生刚体运动。绳索基本是柔性体,抗弯曲能力很小,在约束条件下绳索受外力时,沿绳索方向发生弹性变形,而沿法线方向运动却不受绳索的弹性限制,易出现大变形情况,必须考虑几何非线性。文献[3,4]从多体理论出发,提出计算绳索动力学的有限段模型,这种方法将绳索离散成为一系列铰接刚性绳段组成的多体系统,即用一系列具有不同几何物理参数的刚性绳段近似无限自由度的绳索,应用多体理论求解[5]。有限段方法的有效性在于它可以近似模拟绳索的轮廓形状,保持原有系统的质量分布特性,并可以采用多体系统中的刚体接触模型,求解起重和索道过程中绳索与驱动轮之间的接触关系。对于弹性索,李晓平等[6]在绳段间引入弹簧阻尼,并将多体运动力学和弹性力学相关理论结合起来,提出了采用铰接弹性段组成的多体模型,对弹性的处理效果类似于考虑了几何非线性的杆单元有限元法,可以处理绳索大范围的运动响应。对于预张紧桥梁拉索,文献[7]提出了非线性有限元模型,并在ABAQUS软件中进行了仿真。
客运索道中的拉索同样存在张力,由于工作过程中索道除了弹性变形,还需要做刚体运动,并与驱动轮进行接触,属于多体动力学范畴。ADAMS是多体动力学典型的仿真软件,但它不支持柔索单元。因此本文根据预张力非线性索单元模型,推导出绳索之间的连接力,将绳索离散成直径相等的圆柱小段,段之间用力连接。为了保证位移相容性,增加连接点扭矩和阻尼的联系,防止突变角位移的突变。将等效模型在ADAMS中进行了仿真,并通过计算得出了预张力钢丝绳的频率、振幅与张力的关系,验证了等效多体动力学模型的有效性。
1 非线性索单元模型
设索单元是单向受力构件,随着应变的非线性增大,索力也呈非线性增大。在三维索单元计算中,坐标x、y、z和位移u、v、w的变量表达式为[7]
式中,i、j均为节点编号。
应变公式为
式中,L为索单元的长度。
索的张力为
F=ε SE+F0 (3)
式中,S为截面面积;E为弹性模量;F0为初始张力。
在总体坐标下,单元刚度阵为
其中,单元刚度矩阵的子矩阵k3×3分别由线性和非线性矩阵项组成:
索单元的节点质量为
索单元的质量矩阵为
结构运动方程为
式中,F为作用在结构上的外力;u为结构位移;M为总体刚度矩阵;K为总质量矩阵。
在不断变化的索道中求解该运动方程,得到节点的位移值。
该单元模型已经被应用到有预紧力的荆州长江大桥斜拉桥的计算分析中,并证明了结果是可靠的[7]。在索道中,也存在有预紧力的钢丝绳,钢丝绳与绷轮、驱动轮、从动轮之间是通过摩擦进行传动的,属于多体动力学分析范畴。ADAMS是出色的多体动力学仿真软件,但它没有非线性索单元模型,因此需要根据非线性索单元模型进行等效多体动力学建模。
2 等效多体动力学模型
根据单元位移法,可以求出x、y、z三个方向的等效刚度。在ADAMS中,弹性连接是小位移假设,不适合于索单元等大变形问题。但ADAMS提供通用力连接,其中类轴衬选项可以模拟弹性连接,只需给出连接力的方程即可。索道离散刚体等效动力学模型如图1所示。
根据非线性索单元的单元刚度矩阵,列出等效连接力如下:
离散后的钢丝绳每一小段有12个自由度,总自由度数为12n(n为离散单元数)。与有限元法不同,单元的力连接并不减少自由度方程数,只是建立了相连节点的受力平衡关系,并不能保证位移(特别是角位移)的相容性。因此,需要增加连接点转矩T和扭转刚度KT与阻尼CT的联系:
式中,θji为单元的j节点与另一个单元i节点的夹角;
防止突变角位移的突变。
通过引入弯曲刚度,建立的旋转自由度的受力联系,在线性(逐步加载)受载和运动范围内,可以保证位移的相容性。
3 非线性索道张力下的振动分析
为了验证预张紧钢丝绳等效多体动力学模型的有效性,建立了长为20m、半径为8mm,且两端铰接的一段钢丝绳模型。每50mm离散成一段刚体,共有400个单元。所用材料为钢,钢丝绳拉伸时的弹性模量为110GPa。对不同张力下的运动情况进行分析,线性加载时间为1s,一共分析了5s内的受力情况(图2中示出了0～2.5s的情况)和运动情况(图3中示出了0～2.5s的情况),振动频率由2～5s内稳态响应曲线经过傅里叶变换后的频谱分析得到(图3b)。张力取不同的值,得到的分析结果列于表1。
3.1 频率与钢丝绳张力的关系
从表1数据中,绘出频率与张力的关系,并通过曲线拟合,得出拟合曲线,如图4所示。可以看到,半径为8mm、长为20m钢丝绳的振动频率与张力呈指数函数关系f=0.6201F0.5072,与理论解y=ax0.5(a为系数)相符。
3.2 振幅与钢丝绳张力的关系
同样,绘出振幅与张力的关系,并通过曲线拟合,得出拟合曲线如图5所示。可以看到,半径为8mm,长为20m钢丝绳的振幅与张力呈指数函数关系A=3513.1F-0.5697,符合实际情况。
将创建预张力索道等效动力学模型方法编写成程序,并生成ADAMS软件可执行的命令文件,再运用到本文所述的非线性索道单元等效动力学模型应用与试验索道多体动力学分析中,得到的动力学模型如图6所示。
4 结语
本文建立了非线性索单元的多体动力学等效模型,为索道整体的多体动力学仿真提供了依据。在ADAMS中建立了20m预张力钢丝绳模型,通过多体动力学仿真,得到了不同张力下的钢丝绳振动频率和振幅,并分析了振动频率、振幅与张力的关系。本文结论为索道系统的多刚体动力学仿真提供了依据,也为ADAMS对带传动、链传动等非线性单元的仿真提供了参考。
摘要：根据预应力非线性索单元本构方程,推导了钢丝绳多刚体离散等效节点连接力,将索道离散成刚体小段,在刚体之间施加等效载荷,并引入了扭转刚度和阻尼来保证位移相容性。利用该模型编制了ADAMS程序代码,通过预紧力索道多体动力学进行仿真,揭示了钢丝绳振动频率与张力的关系,为索道多体动力学仿真提供了依据。
关键词：张紧索道,非线性索单元,多体动力学模型,位移相容性
参考文献
[1]丁虎.轴向运动梁横向非线性振动建模、分析和仿真[D].上海:上海大学,2008.
[2]王连华.斜拉索的非线性动力学分析[D].长沙:湖南大学,2001.
[3]Kammen J W,Huston R L.Modeling of VariableLength Towed and Tethered Cable System[J].Journal of Guidance,Control,and Dynamics,1999,22(4):602-608.
[4]Huston R L,Kammen J W.Validation of FiniteSegment Cable Models[J].Computers and Struc-tures,1982,15(6):653-660.
[5]休斯顿,刘又午.多体系统动力学[M].天津:天津大学出版社,1991.
[6]李晓平,王树新,何曼丽.海洋缆索的动力学仿真研究[J].海洋技术,2005,24(1):52-57.
[7]庄茁,由小川,廖剑晖,等.基于ABAQUS的有限元分析和应用[M].北京:清华大学出版社,2009.
非线性索单元 篇2
索穹顶结构属于柔性结构，结构的刚度主要是由预应力提供，结构在荷载作用下呈较强的几何非线性特性。在索穹顶结构施加预应力的过程中，预应力的张拉是很难控制的，要想使实际工程中所施加到结构上的预应力值与设计预应力值接近，就需要对所施加的预应力值进行反复的调试，这就是索穹顶结构不可避免的预应力分布差异现象的产生。这种预应力分布差异的存在是否会影响结构的正常使用，针对这个问题的出现，并结合索穹顶结构在荷载作用下呈现的几何非线性特性，本文编制了相应的计算程序来解决此问题。
2 索穹顶结构静力分析方法
索穹顶结构静力分析方法主要有线性和非线性两种方法。由于索穹顶在荷载作用下呈较强的几何非线性特性，通常采用非线性方法进行求解。对于非线性的问题通常使用增量分析的方法，将荷载分成若干个荷载步进行求解，这样将一个复杂的非线性问题分成若干个线性荷载步。
增量法是指荷载是以增量的形式逐级加至结构上，在每级荷载增量作用的过程中，假定结构的刚度是不变的。节点位移和杆端力是由每级荷载增量区间内区间起点的结构刚度计算出的，然后利用所求的节点位移和杆端力求出区间终点变形后的结构刚度，此刚度作为下一个荷载增量起点的刚度。在第j级荷载增量作用下平衡方程可表示为：
式中，[K]j表示第j级荷载增量区间的起点刚度；{Δε}j+1表示第j级荷载增量引起的节点位移增量；{ΔP}j表示第j级荷载增量；{ε}j、{ε}j+1表示第j级荷载始末的节点位移。
在增量法求解过程中，每次计算结束，都应检查计算结果是否收敛到误差范围内，如果收敛误差太松，最终得到的解是不精确的；如果收敛误差太紧，会造成精度计算量的浪费。适当的收敛准则对于增量求解是至关重要的。一般非线性有限元分析中判断结果是否达到给定的精度，通常采用位移、不平衡力或是两者同时控制的收敛准则。
1）位移收敛准则是指当某次增量迭代中位移增量与当前总位移相比足够小时，认为已经收敛：
式中，‖Δεj+1‖表示位移增量向量范数；βd表示位移收敛容差。
2）不平衡力收敛准则是指某次增量迭代后的不平衡力向量与当前节点外荷载相比，二者相差足够小时，认为已经收敛。
式中，‖ΔPj+1‖表示不平衡力向量范数；‖P‖表示外荷载向量范数；βp表示不平衡力收敛容差。
本文的程序在范数选取上采用的是不平衡力收敛准则。根据Levy式索穹顶结构失效形式[1]为：上斜索内力值减小，索出现松弛而退出工作，丧失对压杆的约束能力这一特点，即第一圈上斜索在荷载作用下内力值小于零。应用本文编制的静力计算程序对实际施加到结构上的预应力值进行判别。若上斜索（共8根）出现负值，即索出现松弛而退出工作，索穹顶结构在正常使用荷载作用下失效，实际施加到结构中的预应力值需要进行调试。反之，若上斜索（共8根）全是正值，则此组施加到结构中的预应力值不需要进行调试即可满足正常使用要求。
本文的实际施加到结构中的预应力值（模拟施工预应力值）采用文献[7]的设计方法进行设计。
3 非线性静力计算程序的实现
3.1 程序设计框图（见图1)
虽然利用几何非线性静力分析方法在理论上是正确的，但是在编程求解过程中，仍然遇到了很大的困难，也就是遇到了求解精度不足的问题。在排除了程序中的缺陷后，这个问题仍然存在，作者先后试验了调整程序中变量精度，仍无法解决该问题。在文献[3]中，也指出用高斯消元法求解平衡方程的精度很差，且计算结果随着变量预设精度变化很大，但文献中并没有提出适用于编程的算法。通过查阅资料与反复试验，利用数学软件MATLAB编写了解平衡方程的程序CALC以及自动累加程序减小计算误差，最终解决了精度问题，实现了几何非线性静力计算。
3.2 程序编写
通过查阅资料与反复试验，利用数学软件MATLAB编写了解平衡方程的程序CALC才最终解决了精度问题。采用了MATLAB程序CALC.m,程序如下：
function calc
load A.dat;
load B.dat;
C=AB;
fid1=fopen('D:matlabworktest.dat','wt');
fprintf(fid1,'n%fn',C(:,1));
fclose(fid1)。
通过执行此程序，用命令C=A/B求解方程，接着程序生成一个数据文件test.dat，并将解向量C中的数据写入文件test.dat，再利用Visual C++中编制的自动累加程序将计算出的内力和坐标增量自动累加，本文采用Visual C++与MATLAB混合编程的方式，直接调用MATLAB的库函数，实现C++程序对MATLAB程序的自动调用，最终完成了多次调用MATLAB的求解过程。
3.3 程序实现
3.3.1 前处理程序块
前处理程序主要包括结构几何尺寸、材料信息、约束条件、各索杆实际施工预应力值、荷载信息等。主要录入的数据有：总圈数（由内向外排列），各圈压杆数、个圈半径（单位：cm）、压杆与斜索的夹角、各压杆弹性模量（单位：N/cm2）、各杆件截面积（单位：cm2）、各索杆预应力值（单位：N）、荷载信息（单位：N）、等，数据录入时都是由内向外的顺序，半径、压杆高度和压杆与斜索的夹角录入都是先内圈后外圈。程序将进行杆件编码、节点编码、计算节点坐标等工作，并最终生成可输出的前处理文件。
3.3.2 计算并输出结果程序块
主要实现数据的输出，例如：节点信息文件、杆件预应力信息文件、材料信息文件等。
4 算例分析
所选结构为2圈8压杆的Levy式索穹顶，外径为400cm，内径为200cm，结构的尺寸见图2，各索、杆截面尺寸及弹性模量见表1，节点荷载值为P=-10kN，模拟施工预应力值见表2，用所编制的几何非线性静力计算程序对该结构模拟施工预应力进行求解计算结果见表3。（由于数据量大仅列出上斜索的相关值）。
根据索穹顶结构在荷载作用下呈现的几何非线性特性，笔者分别利用高斯消元方法以及调用MATLAB库函数的方法分别编制了计算程序。对得出的数据进行分析，最终找出能够体现几何非线性特性的编程方法，并对计算结果和编程方法进行分析。
5 结论
1）根据表3中的计算数据，对于索穹顶结构呈现的几何非线性特性高斯消元法是不适合的，需要采用调用MATLAB库函数的方法进行求解计算。原因在于该非线性方程组行列式值接近于零，该矩阵近似于奇异，属于病态方程组，才会出现求解精度偏大的情况。
2）由计算结果可以看出此组模拟施工预应力值在荷载作用下内力值小于零，是需要对其进行调试，将调试结果多次代入本文所编制的程序中进行计算，直到结果中8根上斜索内力值全部大于零，表明此组施加到结构中的预应力值（模拟施工预应力值）不需要进行调试即可满足正常使用要求。
摘要：针对索穹顶结构在施工过程中不可避免的会产生预应力分布差异的现象,提出合理的方法即几何非线性方法对索穹顶结构进行静力分析,并编制相应的计算程序,判断索穹顶结构在正常使用状态下结构是否有效。此程序为施工预应力是否需要调试提供判别的依据。
关键词：索穹顶结构,施工预应力,设计预应力,几何非线性
参考文献
【1】曹喜.张拉整体索穹顶结构的设计理论与实验研究[D].天津大学博士学位论文,1997(1):69-71.
【2】张建华.索穹顶结构施工成形理论及试验研究[D].北京工业大学博士学位论文,2008(3):13-14.
【3】詹伟东.葵花形索穹顶结构的理论分析与试验研究[D].浙江大学博士学位论文,2004(5):25-26.
【4】曹喜,张娜,王利辉.索穹顶结构预应力设计及结构找形分析的研究[D].内蒙古工业大学硕士学位论文,2008(11):22-24.
【5】曹喜,王利辉,张娜.索穹顶结构施工预应力水平评估的研究[D].内蒙古工业大学硕士学位论文,2008(11):18-25.
【6】史勇.张拉整体索穹顶结构分析程序TCDS的研究与开发[D].内蒙古工业大学硕士学位论文,2006(11):46-48.
非线性索单元 篇3
目前,国内胡建华、王连华等学者通过对索结构单元进行柔性迭代分析,依据内力改变与位移增量的关系,建立了一种具有较高精度的索结构几何非线性分析方法———悬链线索单元法[2]; 西南交通大学的冉志红,李乔等人在考虑垂度和斜度的基础上,利用奇异摄动解法,推导出了斜拉索振动频率和振型的解析式[3]。国外,Kiisa Martti、Idnurm Juhan和Idnurm Siim三位学者在集中载荷的作用下,对弹性索进行了离散分析[4]; 而Vairo Giuseppe教授根据应变与弦长变化的关系,针对斜拉索的非线性响应提出了一个精致的闭式二阶模型[5]。但是,这些研究都是针对悬索桥或斜拉索桥的单一索结构,并且也没有将影响索结构几何非线性的因素进行完全考虑。悬索管桥的索结构作为一个体系,不但各索之间相互耦合作用,多种几何非线性行为也相互耦合,具有独特的结构特性和受力特点,笔者将在悬索管桥正常运营的现状态下,充分考虑几何非线性因素耦合作用的基础上,研究其索结构体系的应力分布规律,为以后悬索管桥的优化设计、可靠性分析以及风险评价提供参考依据。
1 索结构体系
一座悬索式管桥的完整索结构体系一般应该由1 对主索、1 对风索、1 对下稳定索、若干对吊索及系索组成。其中最重要的是主索,主索两端与塔架相连,整体几何形状在竖直平面内呈抛物线型; 它在整个悬索管桥的结构中起着至关重要的作用,是最主要的受力结构,基本上承载着整个管桥的所有内在载荷和外在载荷; 一旦主索发生失效破坏,将对悬索管桥产生致命的影响,会导致整个管桥的倒塌、断裂。风索的两端固定在水泥制成的锚固墩上,在水平面内呈抛物线型,在竖直平面内通常伴随着桥桁架的形状,可能呈拱形分布,也可能水平分布。由于悬索管桥的具有大柔性、大跨度的特点,桥身结构在横向极其不稳定,受到横向载荷作用时,容易出现大幅度的漂移,导致管桥毁坏,发生事故,风索的作用就是加强管桥在横向的约束,是管桥承受横向载荷的主要结构。下稳定索的两端也固定锚固墩上,索体与桥桁架底部相连,因此沿着桥桁架分布,在横向和纵向都起着稳定整个管桥的作用。吊索是连接主索和桥桁架的结构,是除主索以外的主要受力索。系索则连接着风索和桥桁架,正常情况下与吊索的对数相同,和风索一起对管桥进行横向约束。
2 理论分析
2. 1 几何非线性行为
在结构分析中,非线性问题包括材料非线性、几何非线性和状态非线性三类; 其实质是求解过程中结构的刚度矩阵不断变化。几何非线性是指结构经受大变形,几何形状的变化引起结构响应的非线性。悬索管桥的索结构属于大跨度、低阻尼、超柔性体系,变形特征为大位移、小应变,是典型的几何非线性结构[6]。虽然索结构的刚度与位移呈明显的非线性关系,但应力、应变仍符合虎克定律,即材料仍是线性的。索结构体系的几何非线性行为主要表现为:
2. 1. 1 结构大位移
这是影响索结构几何非线性最主要的因数。在外载荷作用下,主索、风索、下稳定索及系索会发生下挠,吊索将会出现伸长和倾斜,这些位移通常情况下比较大,不可避免地将对索结构的内力产生影响。由于力和变形是非线性的,因此在对索结构进行分析时,要依据变形后的几何位置来建立力的平衡方程[7]。
2. 1. 2 初始应力的影响
为了使柔性的索结构保持一定的几何形状,在悬索管桥成桥时,会对各个索施加一定的初始应力,这必然会使索结构具有初始刚度。在受到后续的载荷作用时,初始应力产生的初始刚度对后续的变形存在一定的抵抗力,而后续载荷建立新的平衡时必须打破原有平衡,这个过程中必然存在着能量消耗,也就是说初始应力可以提高结构刚度。
2. 1. 3 自重垂度的作用
进行有限元分析时,索单元通常选取直杆单元。但实际上索结构在自重的作用下将会产生一定的垂度,受到轴向力作用时,其弦长变化由材料的弹性伸长和结构的垂度变化两部分组成,其中后者使轴向力与弦长形成几何非线性关系[8]。
2. 2 耦合刚度矩阵
对悬索管桥索结构进行分析时,充分考虑以上三种非线性行为的耦合作用,采用有限位移理论[9],其耦合的刚度矩阵表达式如下:
式( 1) 中Ke为结构弹性刚度矩阵,KL为大位移刚度矩阵,R为恒活载内力,Kσ为结构几何刚度矩阵,δ为节点位移向量。
选取两节点组成的直杆单元来模拟索单元。取一个杆单元,长度为l ,横截面积为A ,弹性模量为E ,节点位移的列阵形式为qe= {ui,vi,uj,vj},由于索结构在实际状态下除了自重,只受轴向拉力,所以杆单元也只承受轴向拉力的作用,其线性应变根据经典的弹性理论得:
但是索结构是非线性的,横向位移v会引起单元的附加伸长,因此杆单元的应变可表达为:
式( 3) 中为应变的非线性项。不考虑杆单元的弯曲变形,那么位移便是x的线性函数:
式( 4) 中I = I2为二阶的单位矩阵,Ni、Nj为形函数,,由公式( 3) 和公式(4) 可以得到:
式( 5) 中T为体载荷向量,。因为杆在变形后,转角相对较小,故。
根据胡克定律,杆的拉伸应力为:
考虑到索结构还受自重垂度的非线性影响,所以上式的弹性模量E应依据Ernst公式[10]进行换算,采用等效弹性模量Eeq:
式( 7) 中W为单位长度索的重力,l' 为索弦长在水平面的投影长度,F为索内的张力。故杆的应力变为:
将各个参数值代入式( 8) 可得,结构弹性刚度矩阵为:
几何刚度矩阵为:
因为F是索的初始轴向力,即[K]σ也称为初始应力刚度矩阵,体现了成桥时初始轴向力对横向位移的作用。
大位移刚度矩阵为:
故最终索结构体系的耦合刚度矩阵为:
综上所述,所构建的刚度矩阵综合考虑了索结构的大位移、初始应力以及自重垂度三种几何非线性行为的影响,理论上具有足够的精度,可以用来进行下面的有限元分析。
3 有限元模型
3. 1 模型的建立
选取索结构体系比较完整的黄河悬索管桥为建模对象,管桥横跨黄河两岸,一端位于陕西省府谷县,另一端位于山西省保德县。桥跨长度270 m,塔架高44 m,桥面为桁架结构,两边高、中间低,呈微拱形,高度差为6 m,输气管道为 Φ660 × 14. 3 mm,最大输气压力6. 4 MPa,钢丝索包括1 对主索、1 对风索、1 对下稳定索、53 对吊索和53 对系索,每对吊、系索间相隔5 m。根据管桥各部分结构的受力特点和材料特性[11]选择单元类型,具体的单元材料特性如表1 所示。
结合工程实际情况,进行结构模型的简化及优化和边界条件的约束,依据施工图纸以及现场测量结果,运用ANSYS软件建立管桥的三维空间有限元仿真模型[12],为了方便研究,管桥桥跨中点为坐标原点,从府谷端到桥跨中点的X轴坐标为负值,从桥跨中点到保德端的X轴坐标为正值,并沿着管道内部天然气介质输送的方向逐渐增大,同时将每对吊索和系索进行编号,从小到大依次标号1 到53。图为1 全桥模型,图为2 索结构体系模型。
3. 2 索结构找形
以悬索管桥的成桥状态为模型初始平衡状态[13]进行索结构找形,因为此时管桥各部分结构的几何位置是已知的,所以模型在受到成桥恒载作用时应保持与管桥成桥状态相同的几何位置,这可以通过不断调整索结构的初始应变来近似地实现,假设主索的初始应变为 ε主0,风索的初始应变为 ε风0,下稳定索的初始应变为 ε下0,吊索的初始应变为ε吊0,系索的初始应变为 ε系0,则管桥模型成桥状态最大位移值f0max是 ε主0、ε风0、ε下0、ε吊0、ε系0的函数,记为f0max= f( ε主0,ε风0,ε下0,ε吊0,ε系0) ,那么当f0max等于零或接近于零时,所得到的 ε主0、ε风0、ε下0、ε吊0、ε系0值便是各个索结构在成桥状态下的初始应变值,所得模型便是后续分析的理想模型。考虑到主索、风索、下稳定索、吊索、系索的耦合作用,找形计算[14]通过循环程序实现,找形结果见表2。
4 模拟分析与验证
在管桥正常运营的工况下对索结构体系进行几何非线性有限元分析。该管桥已运行15 年,各部分结构都存在着一定程度的腐蚀损耗,因此必须依据现场实测的结构现尺寸数据对管桥模型进行修正[15]。正常运营工况是指只考虑输气管道内输送的天然气介质重力载荷和管桥结构自重,不考虑其他外载荷。通过加载计算,得到了主索、风索、下稳定索、吊索、系索五种索结构的应力分布规律,并将模拟结果与工程检测[16]数据进行对比验证。其中主索、风索、下稳定索的实测数据来自于该悬索管桥的安全监测预警系统[17],每索9 个测点。吊索和系索由于数量较多,首先通过SL - 30T型钢丝绳张力测量仪现场逐一测得每对索的轴向张力,然后换算成轴向应力。另外,由于主索、风索、下稳定索在整个索结构体系中起着主要作用,而且各只有1 对,所以模拟分析它们整条索不同位置的应力变化情况,图为3 为主索轴向应力曲线,图4 为风索轴向应力曲线,图5 为下稳定索轴向应力曲线; 对于吊索和系索,因为数量众多,长短不一,且每对索所处的位置也不同,不考虑每对索具体的应力变化情况,而是把它们作为一个整体来研究,图6 为吊索体系应力分布图,图7 为系索体系应力分布图。
从以上各图可以看出,仿真模拟结果与现场实测数据吻合良好。其中主索的轴向应力呈“U”型分布,在管桥中点处最小,距离管桥中点越远,其值越大,最大值为290 MPa左右; 风索的轴向应力值也呈现中间低、两边高的走势,但变化较为平缓,最大值约80 MPa,最小值约70 MPa; 而下稳定索的轴向应力基本保持不变,为一个恒定值,仿真结果为94. 4MPa,现场实测结果稍大,稳定在100 MPa左右; 对于吊索体系,管桥两端的几对应力值明显高于其它吊索,且越接近桥端,其值越大,其中第1 对和第53对的轴向应力值达到130 MPa,剩余吊索的应力变化比较平缓,保持在55 MPa左右; 系索体系,每对索的轴向应力越接近管桥中心越高,但峰值并不是出现在中间的第27 对系索上,而是出现在中间附近两边的系索上。
5 结论
( 1) 悬索管桥各个索的应力值都相对较小,远小于其标准破坏强度1 570 MPa,说明该索结构体系确实是小应变、大位移的几何非线性结构。
( 2) 由于管桥结构的对称性,每条主索、风索、下稳定索的轴向应力沿桥中心呈对称分布,吊索体系、系索体系的各对索的轴向应力以各自的第27 对索为中心,也呈对称分布。
( 3) 系索的最大应力没有出现在第27 对索上,而是在其两侧附近的系索上,可能是受管桥桥面桁架微拱形结构的影响。
非线性索单元 篇4
1 斜拉索的模拟
目前大多数有限元分析软件都没有精确的索模型来模拟斜拉索的性质，所以在选择索模型时，需要注意各种模型的优缺点以及适用范围。目前常用的斜拉索的模拟方法有等效弹性模量法、多段杆单元法、曲线索单元法等。
1.1 等效弹性模量法（切线模量法）
在有限元分析中，若将斜拉桥的拉索直接简化成直杆单元，就无法考虑垂度影响。为了考虑垂度效应，工程界广泛采用了由Ernst提出的等效弹性模量的方法。Ernst修正的实质是将拉索的几何非线性问题转化为材料非线性问题。拉索弹性模量的Ernst修正公式为：
式中：Eeq为等效弹性模量;E为弹性模量;γ为容重;l为拉索水平投影长度;σ为拉索拉应力。
当索内应力水平较高、拉索水平投影长度较小、应力幅值不大时，等效弹性模量法具有较高精度。但是对于低应力或是大位移状态下的长索，切线模量法的误差很大。切线模量法由于无法求出不平衡力，所以只能采用类似简单增量法求解，不可避免地产生漂移误差。切线模量法的计算精度与计算所采用的荷载步个数密切相关。尤其对于较长的斜拉索，增加荷载步个数，可以提高切线模量法的计算精度。在基于ANSYS平台进行斜拉桥分析时，需要不断获取拉索应力状态，并反馈给ANSYS，实现这一过程相当繁琐。
1.2 多段杆单元法[1]
ANSYS可以考虑结构的几何非线性的特性，采用多段直杆法，即使用多个杆单元(如LINK8等)来模拟拉索。因为杆单元不约束转角自由度，不能承担弯曲载荷，拉索在自重作用下就会下垂，接近其工作状态，即通过多个杆单元节点来定义拉索中间点的运动，从而模拟索的非线性行为。一般来说，所取杆单元个数够多的话，计算精度能满足要求。
为了满足计算精度，多段杆单元法引入了许多自由度，而且是几何非线性的，这大大增加了计算机存储量和计算时间，要求输入的数据量大，节点坐标、初始应变输入的准确性要求很高，有时难以收敛。
1.3 曲线索单元[2]
除了以上2种方法，还有多节点曲线索单元、悬链线索单元等方法，但这些方法使用较为复杂。
1.4 σ-ε法
ANSYS具有强大的非线性分析功能，其中就包括非线性弹性材料。在定义实际材料性质时，ANSYS提供了各种可以自行定义的项目，即自己定义σ-ε关系。
基于σ-ε的索模拟方法，根据弹性索的解析式导出索元的刚度矩阵，且能适应各种极端情形，因此是一个精度高、通用性强的方法。它的主要优点是：计入垂度影响，能模拟各种极端情况;分析模型的自由度个数少;通过仅仅给定索力(初始应变)和端点坐标即可方便地确定出结构的施工状态和使用状态下的位移和内力，保持了结构的连续性和计算精度。
与Ernst修正弹性模量法相比，割线模量法同样具有很高的精度，公式简单，而且割线模量法对应的刚度为割线刚度，不会产生计算漂移等误差。因此建议在建立斜拉索σ-ε关系时采用割线模量法，见公式(2)和(3)，对应初应变按照公式(4)计算。
根据笔者粗略统计，一般情况下，拉索第1次张拉和第2次张拉的应力水平大约为130～350 MPa和350～500 MPa。图1为某根斜拉索在初应力水平下σ-ε曲线，斜拉索水平投影长度l=200 m，高度c=100 m，弹性模量E=2.0×105 MPa，斜拉索面积A=0.011 m2，重度γ=88.3 kN/m3，初应力σ0分别为130 MPa、350 MPa。由图1可见：(1)初应力不同，相同索力时名义应变相差较大;(2)不同初应力水平下，σ-ε曲线上段平行;(3)初应力水平越高，σ-ε曲线与材料本身的弹性模量E线越接近。
2 斜拉索张拉的实现
在斜拉桥的施工过程中，斜拉索张拉改变的是拉索的长度，其实质是改变无应力索长。在ANSYS中，可采取改变拉索温度的方法来等效这一过程。但从精确、简便要求来衡量，“标准”方法更可取，即斜拉索张拉分2步。第1步：斜拉索张拉力的施加;第2步：去掉第1步施加的力，激活给定初始应变的斜拉索单元。
3 算例与结果分析
为验证σ-ε法的有效性，对某斜拉索进行了计算，斜拉索水平投影长度分别为l=200 m、l=400 m，对应高度分别为c=100 m、c=200 m，弹性模量E=2.0×105MPa，斜拉索面积A=0.011 m2，重度γ=88.3kN/m3，初应力σ0分别为130 MPa、350 MPa。约束塔端水平、竖向位移和梁端竖向位移，在梁端施加一竖向力，分别用多段杆单元法、σ-ε法计算得出主梁梁端竖向位移。不同应力水平下计算结果见表1、表2。
从表1可见，本文的σ-ε方法计算得到的拉索的力学性能和多段直杆法所得的结果一致，误差极小。不论是在高应力状态下还是在低应力状态下，也不论水平投影长度大小和位移变化大小，σ-ε法计算结果精度均很高。另外，多段直杆法计算时应考虑大变形效应和应力刚化，否则位移误差较大。
4 结语
本文基于在ANSYS平台上模拟斜拉索，提出了自定义σ-ε关系模拟拉索的垂度效应，并编制了与ANSYS相匹配的功能模块，通过算例验证了本文计算方法简便、准确。
摘要：为了在ANSYS平台上模拟斜拉索,比较了几种常用的方法,提出了自定义σ-ε关系模拟拉索的垂度效应,编制了与ANSYS相匹配的功能模块,并通过算例验证了这些方法的有效性。
关键词：桥梁工程,斜拉索单元,垂度效应,σ-ε法,ANSYS
参考文献
[1]周强,杨文兵,杨新华.斜拉桥索力调整在ANSYS中的实现[J].华中科技大学学报(城市科学版),2005,22(增刊):27-28,38.
[2]施溪溪,李鸿晶.斜拉桥拉索单元模型及其计算模拟[J].钢结构,2005,20(5):53-56,89.