【摘要】本文从简单实用的角度论述了空间杯系结构的几何非线性分析理论。文中分析了非线性有限元方法的求解过程，特别强调决定几何非线性收敛结果的关键问题，即由节点位移增量计算单元的内力增量。通过引入随转坐标系，论述了平面和空间梁单元小应变变形时单元内力增量的计算问题。用本文方法可以分析大跨度
桥梁结构的六位移大旋转问题。并且用实桥算例进行了验证。
关键词 大跨度桥梁 几何非线性 实用分析 非线性有限元 小应变理论 江阴 长江大桥
一.引言．
现代大跨度桥梁等工程结构的柔性特征已十分明显，对于这些结构考虑几何非线性的影响己必不可少。并且，计算机能力的大大提高也使得分析大型复杂结构的非线性问题成为可行。80年代国外对几何非线性问题的发展已相当完善［1，2］，国内在这方面也做了不少的工作［4－6］
在工程结构几何非线性分析中，按照参考构形的不同可分为TL（Total Lagranrian)法和UL（Updated Lagrangian)法[1]。后来，引入随转坐标系后又分别得出 CR（Co-rotational)-TL法和CR－LU法[2,3]，在工程中UL（或CR-UL）法应用较多。以前的文献大都对结构的几何刚度矩阵进行了复杂而详细的推导。从文中的分析可以发现，结构几何刚度矩阵的精确与否并不实质性地影响迭代收敛的最终结果，求解几何非线性问题的关键在于如何由节点位移增量准确地计算出单元的内力增量，而这一点以前文献都没有提到过。因此，本文的重点放在论述单元内力增量的计算上。
工程上很早就开始使用拖动坐标系来求解大跨度桥梁结构的大挠度问题，本文则把它应用到单元内力增量的计算中。从实质上说，这里的拖动坐标系与上面提到的随转坐标系没有区别。因此，在理论方法上，目前文中的方法可以归类到CR－UL法。但由于本文重点不在于详细介绍这种方法的理论体系，所以论述中均不再使用该名词。本文的目的主要是通过简化复杂的几何非线性分析方法，推广该方法在实际工程中的应用。
二、非线性商限元求解过程
对于工程结构的非线性问题，用有限元方法求解时的非线性平衡方程可写成以下的一般形式：
Fs(δ）-P0（δ）=0 （l）
其中，为节点的位移向量；Fs(δ)为结构的等效节点抗力向量，它随节点位移及单元内力而变化；PO（δ）为外荷载作用的等效节点荷载向量，为方便起见，这里暂时假定它不随节点位移而变化。
由于式（l）中的等效节点抗力一般无法用节点位移显式表示，故不可能直接对非线性平衡方程进行求解。但实际结构的整体切向刚度容易得到，所以通常应用Newton-Raphson迭代方法求解该问题。结构的整体切向刚度矩阵KT可表示如下
dPO＝ KTdδ （2）
式中，KT＝ KE十KG，其中KE为结构的整体弹性刚度矩阵，KG为几何刚度矩阵。
用混合Newton－Raphson迭代方法求解结构非线性问题的基本过程如下：
(1）将等效节点荷载PO分成n步，ΔP0＝PO／n，计算并组集结构的整体切向刚度矩阵，进入加载步循环；
（2）求解节点位移增量；
（3）计算各单元内力增量，修正单元内力；
（4）更新节点坐标，计算节点不平衡力R；
（5）判断节点不平衡力R是否小于允许值，如满足条件，则进入下一个加载步；如不满足条件，重新计算结构的整体切向刚度矩阵，用R代替ΔP0，回到第2步；
(6）全部加载步完成之后，结束。
从上述求解过程中可见，最为关键的一步是第3步，即由节点位移增量计算单元的内力增量。也可以说是由这一步决定了最终的收敛结果，以下将对此着重论述。其实结构的整体切向刚度矩阵对结果并无实质性的影响，修正的NetwRaphson方法正是利用这一点来节省迭代计算的时间。
以前的文献对空间梁单元几何刚度矩阵的推导方面论述较多，都建立在一些假定的基础上，这里就不详细说明。考虑到结构的整体切向刚度矩阵精确与否并不改变最终结果，仅影响迭代收敛的速度，并且不是越精确的整体切向刚度矩阵迭代收敛越快。所以，本文建议计算时的单元几何刚度矩阵来用以下最普遍的形式：
式中，N为单元的轴向拉力。
三、小应变时单元内力增百计算
在一般情况下，工程结构的几何非线性都属于小应变大位移（大平移、大转动）问题。对于这类问题，单元内力增量的计算比较简单。平面梁单元是空间梁单元发展的基础，故这里先分析平面梁单元的情况。
图1表示了平面梁单元在整体坐标系（OXY）下从t到t十Δt时刻的变形情况。定义随转坐标系（oxy）的原点固定在单元的一端（i端），x轴始终保持沿i→j的直线方向。可见，在随转坐标系中平面梁单元的自由度减少为三个（uxθiθj），并有如下的关系式
式中，L（t）和L（t十Δt）分别为 t和t十Δt时刻梁单元的直线长度。
从随转坐标系中的三个自由度可以看出，它们反映的是单元的真实变形情况，与单元所经历的刚性位移无关。在用有限元方法求解非线性问题时，只要将单元尺寸划分得适当小，整体坐标系下的小应变大位移问题在单元随转坐标系中就转化为小应变小位移问题，这一点可从非线性连续介质力学给出证明。这样，随转坐标系下的受力变形情况就可近似地接线性处理，单元内力增量的计算也就与线性情况一样，这里不再赘述。同时也正说明了工程中常用拖动坐标法计算平面结构大变形问题的正确性。
四、算例分析
结合以上论述，编制了相应的非线性有限元计算程序。为验证本文方法和有限元程序，下面首先分析了45度弯梁空间弯扭大位移问题。大跨度悬索桥在施工阶段的几何非线性比较明显，因此，必须准确地考虑，否则计算结果可能不正确。作为实桥算例，对江阴长江大桥在20％拼装率施工阶段的几何非线性问题进行了分析，并与Ansys程序的计算结果相比较。
1．45度弯梁空间弯扭六位移分析
本例是ADINA中的45度弯梁大位移分析考题，梁的形状和截面尺寸见图5。该梁位于X一y平面内，梁根固定，在自由端沿Z方向受一个集中荷载的作用，梁因此发生空间弯扭大变形。分析时将梁划分为8个单元，每步加载量为10.0，有关参数见图5中所示。分别用ADINA，AnsyS和本文程序计算了60个加载步，各计算结果均基本上一致。梁自由端无量纲位置坐标在初始时刻，加载30步与加载60步时的比较列于表1，可见三者相互较吻合。为了进行对比，都没有考虑剪切影响。为简洁起见，这里不指定专门的量纲单位。
2．江阴长江大桥非线性分析
江阴长江大桥是我国目前建成的最大跨度悬索桥，见图6所示。主跨跨度为1385m。主梁为宽36.9m，高3.0m的扁平状闭口钢箱梁。主缆相距32.5米，吊杆间距为1.6m，矢跨比为 l／10.5。桥塔为门式框架结构，南北桥塔高分别为187m和184m。桥面波置为R＝27710m的竖曲线。
根据设计资料，建立了江阴长江大桥的计算模型。在成桥状态下，单根主缆的水平内力约为23878t，单根吊杆的内力约为144t，考虑到悬索桥在施工时主缆与塔顶有相对位移，计算模型中主缆与塔顶在顺桥向可自由移动，而其他方向均耦联。
悬索桥施工过程中分段安装主梁，小拼装率时各主梁段之间相互饺接。由于悬索桥在成桥状态的位置和内力一般为已知，施工状态均从成桥状态通过拆除梁段的方法确定。图7绘出了江阴长江大桥在成桥状态拆除两端梁段后，但未发生变形之前20％拼装率的初始状态。由于该初始状态的节点位置和单元内力均为己知，用以上的非线性有限元程序可得出20％拼装率变形后的施工平衡状态。见图8所示。从图中可见，在变形后的施工平衡状态下，跨中梁段随主缆发生了较大的烧曲，主跨两端的主缆形状比成桥状态时变化较明显。跨中竖向向下的位移为2.595m，塔顶处的主缆发生向外0.669m的位移。为验证该程序的计算结果，对上述同样的工况用AnsyS程序进行了分析。计算结果中跨中向下的位移为2.581m，塔顶处主缆向外的位移0.664m。从本文方法的变形和内力结果与Ansys程序计算结果的比较发现，两者均较吻合，这就验证了本文方法和非线性有限元程序的可靠性和有效性。
五、结语
以上从简单实用的角度论述了空间杯系结构的几何非线性分析理论。通过对有限元求解几何非线性问题过程的分析，特别强调了用选代方法求解杯系结构几何非线性问题中的关键问题，即由节点的位移增量计算单元内力增量的重要性。在引人随转坐标系之后，论述了小应变问题中单元内力增量的计算。从论述中可知，随转坐标系下的受力变形情况可近似地接线性处理，单元内力增量的计算也与线性情况一样，同时也说明了工程中常用拖动坐标法计算大跨度桥梁结构大变形问题的正确性。
对空间杆系结构用数值算例对本文方法进行了验证。为保证分析结果的正确性，用多个程序进行相互校核。
对江阴长江大桥在20％拼装率施工阶段的几何非线性问题进行了分析，分析结果与
AnsyS程序的计算结果吻合。从分析中可见，在小拼装率施工阶段，悬索桥跨中梁段随主缆
发生了较大的挠曲，主跨两端的主缆形状比成桥状态时变化比较明显。
参考文献
[1]K．J．Bathe＆ S．Bofourchi，Large displacement analysis of three－dimensional beam structures，Int．J．Num．Meth．Eng，V14（1979（， 961－986
［2］K．M．Hsiao，H．J．Horng ＆Y．R．Chen，A corotational procedure that handle large rotations of Spstial beam structures Comput.Struct．27（1987），769－781
［3］G． Narayanan， C．S．Krihnamoorthy，An investigation of geometric non－linear formulation for 3－D beam element， Int． J．Non-linear Mechanics，1990，25（6）， 643-662
［4］陈政清，曾庆元，颜全胜．空间杆系结构大挠度问题内力分析 UL列式法．土木工程学报，1992，25（2），34-43
[5]黄文，李明瑞，黄文彬.杆系结构的几何非线性分析一Ⅱ三维问题．计算结构力学及其应用，1995，12（2），133-141
［6］潘家英，程庆国．大跨度悬索桥有限位移分析．土木工程学报，1994，27（1），l-10
桥梁结构的六位移大旋转问题。并且用实桥算例进行了验证。
关键词 大跨度桥梁 几何非线性 实用分析 非线性有限元 小应变理论 江阴 长江大桥
一.引言．
现代大跨度桥梁等工程结构的柔性特征已十分明显，对于这些结构考虑几何非线性的影响己必不可少。并且，计算机能力的大大提高也使得分析大型复杂结构的非线性问题成为可行。80年代国外对几何非线性问题的发展已相当完善［1，2］，国内在这方面也做了不少的工作［4－6］
在工程结构几何非线性分析中，按照参考构形的不同可分为TL（Total Lagranrian)法和UL（Updated Lagrangian)法[1]。后来，引入随转坐标系后又分别得出 CR（Co-rotational)-TL法和CR－LU法[2,3]，在工程中UL（或CR-UL）法应用较多。以前的文献大都对结构的几何刚度矩阵进行了复杂而详细的推导。从文中的分析可以发现，结构几何刚度矩阵的精确与否并不实质性地影响迭代收敛的最终结果，求解几何非线性问题的关键在于如何由节点位移增量准确地计算出单元的内力增量，而这一点以前文献都没有提到过。因此，本文的重点放在论述单元内力增量的计算上。
工程上很早就开始使用拖动坐标系来求解大跨度桥梁结构的大挠度问题，本文则把它应用到单元内力增量的计算中。从实质上说，这里的拖动坐标系与上面提到的随转坐标系没有区别。因此，在理论方法上，目前文中的方法可以归类到CR－UL法。但由于本文重点不在于详细介绍这种方法的理论体系，所以论述中均不再使用该名词。本文的目的主要是通过简化复杂的几何非线性分析方法，推广该方法在实际工程中的应用。
二、非线性商限元求解过程
对于工程结构的非线性问题，用有限元方法求解时的非线性平衡方程可写成以下的一般形式：
Fs(δ）-P0（δ）=0 （l）
其中，为节点的位移向量；Fs(δ)为结构的等效节点抗力向量，它随节点位移及单元内力而变化；PO（δ）为外荷载作用的等效节点荷载向量，为方便起见，这里暂时假定它不随节点位移而变化。
由于式（l）中的等效节点抗力一般无法用节点位移显式表示，故不可能直接对非线性平衡方程进行求解。但实际结构的整体切向刚度容易得到，所以通常应用Newton-Raphson迭代方法求解该问题。结构的整体切向刚度矩阵KT可表示如下
dPO＝ KTdδ （2）
式中，KT＝ KE十KG，其中KE为结构的整体弹性刚度矩阵，KG为几何刚度矩阵。
用混合Newton－Raphson迭代方法求解结构非线性问题的基本过程如下：
(1）将等效节点荷载PO分成n步，ΔP0＝PO／n，计算并组集结构的整体切向刚度矩阵，进入加载步循环；
（2）求解节点位移增量；
（3）计算各单元内力增量，修正单元内力；
（4）更新节点坐标，计算节点不平衡力R；
（5）判断节点不平衡力R是否小于允许值，如满足条件，则进入下一个加载步；如不满足条件，重新计算结构的整体切向刚度矩阵，用R代替ΔP0，回到第2步；
(6）全部加载步完成之后，结束。
从上述求解过程中可见，最为关键的一步是第3步，即由节点位移增量计算单元的内力增量。也可以说是由这一步决定了最终的收敛结果，以下将对此着重论述。其实结构的整体切向刚度矩阵对结果并无实质性的影响，修正的NetwRaphson方法正是利用这一点来节省迭代计算的时间。
以前的文献对空间梁单元几何刚度矩阵的推导方面论述较多，都建立在一些假定的基础上，这里就不详细说明。考虑到结构的整体切向刚度矩阵精确与否并不改变最终结果，仅影响迭代收敛的速度，并且不是越精确的整体切向刚度矩阵迭代收敛越快。所以，本文建议计算时的单元几何刚度矩阵来用以下最普遍的形式：
式中，N为单元的轴向拉力。
三、小应变时单元内力增百计算
在一般情况下，工程结构的几何非线性都属于小应变大位移（大平移、大转动）问题。对于这类问题，单元内力增量的计算比较简单。平面梁单元是空间梁单元发展的基础，故这里先分析平面梁单元的情况。
图1表示了平面梁单元在整体坐标系（OXY）下从t到t十Δt时刻的变形情况。定义随转坐标系（oxy）的原点固定在单元的一端（i端），x轴始终保持沿i→j的直线方向。可见，在随转坐标系中平面梁单元的自由度减少为三个（uxθiθj），并有如下的关系式
式中，L（t）和L（t十Δt）分别为 t和t十Δt时刻梁单元的直线长度。
从随转坐标系中的三个自由度可以看出，它们反映的是单元的真实变形情况，与单元所经历的刚性位移无关。在用有限元方法求解非线性问题时，只要将单元尺寸划分得适当小，整体坐标系下的小应变大位移问题在单元随转坐标系中就转化为小应变小位移问题，这一点可从非线性连续介质力学给出证明。这样，随转坐标系下的受力变形情况就可近似地接线性处理，单元内力增量的计算也就与线性情况一样，这里不再赘述。同时也正说明了工程中常用拖动坐标法计算平面结构大变形问题的正确性。
四、算例分析
结合以上论述，编制了相应的非线性有限元计算程序。为验证本文方法和有限元程序，下面首先分析了45度弯梁空间弯扭大位移问题。大跨度悬索桥在施工阶段的几何非线性比较明显，因此，必须准确地考虑，否则计算结果可能不正确。作为实桥算例，对江阴长江大桥在20％拼装率施工阶段的几何非线性问题进行了分析，并与Ansys程序的计算结果相比较。
1．45度弯梁空间弯扭六位移分析
本例是ADINA中的45度弯梁大位移分析考题，梁的形状和截面尺寸见图5。该梁位于X一y平面内，梁根固定，在自由端沿Z方向受一个集中荷载的作用，梁因此发生空间弯扭大变形。分析时将梁划分为8个单元，每步加载量为10.0，有关参数见图5中所示。分别用ADINA，AnsyS和本文程序计算了60个加载步，各计算结果均基本上一致。梁自由端无量纲位置坐标在初始时刻，加载30步与加载60步时的比较列于表1，可见三者相互较吻合。为了进行对比，都没有考虑剪切影响。为简洁起见，这里不指定专门的量纲单位。
2．江阴长江大桥非线性分析
江阴长江大桥是我国目前建成的最大跨度悬索桥，见图6所示。主跨跨度为1385m。主梁为宽36.9m，高3.0m的扁平状闭口钢箱梁。主缆相距32.5米，吊杆间距为1.6m，矢跨比为 l／10.5。桥塔为门式框架结构，南北桥塔高分别为187m和184m。桥面波置为R＝27710m的竖曲线。
根据设计资料，建立了江阴长江大桥的计算模型。在成桥状态下，单根主缆的水平内力约为23878t，单根吊杆的内力约为144t，考虑到悬索桥在施工时主缆与塔顶有相对位移，计算模型中主缆与塔顶在顺桥向可自由移动，而其他方向均耦联。
悬索桥施工过程中分段安装主梁，小拼装率时各主梁段之间相互饺接。由于悬索桥在成桥状态的位置和内力一般为已知，施工状态均从成桥状态通过拆除梁段的方法确定。图7绘出了江阴长江大桥在成桥状态拆除两端梁段后，但未发生变形之前20％拼装率的初始状态。由于该初始状态的节点位置和单元内力均为己知，用以上的非线性有限元程序可得出20％拼装率变形后的施工平衡状态。见图8所示。从图中可见，在变形后的施工平衡状态下，跨中梁段随主缆发生了较大的烧曲，主跨两端的主缆形状比成桥状态时变化较明显。跨中竖向向下的位移为2.595m，塔顶处的主缆发生向外0.669m的位移。为验证该程序的计算结果，对上述同样的工况用AnsyS程序进行了分析。计算结果中跨中向下的位移为2.581m，塔顶处主缆向外的位移0.664m。从本文方法的变形和内力结果与Ansys程序计算结果的比较发现，两者均较吻合，这就验证了本文方法和非线性有限元程序的可靠性和有效性。
五、结语
以上从简单实用的角度论述了空间杯系结构的几何非线性分析理论。通过对有限元求解几何非线性问题过程的分析，特别强调了用选代方法求解杯系结构几何非线性问题中的关键问题，即由节点的位移增量计算单元内力增量的重要性。在引人随转坐标系之后，论述了小应变问题中单元内力增量的计算。从论述中可知，随转坐标系下的受力变形情况可近似地接线性处理，单元内力增量的计算也与线性情况一样，同时也说明了工程中常用拖动坐标法计算大跨度桥梁结构大变形问题的正确性。
对空间杆系结构用数值算例对本文方法进行了验证。为保证分析结果的正确性，用多个程序进行相互校核。
对江阴长江大桥在20％拼装率施工阶段的几何非线性问题进行了分析，分析结果与
AnsyS程序的计算结果吻合。从分析中可见，在小拼装率施工阶段，悬索桥跨中梁段随主缆
发生了较大的挠曲，主跨两端的主缆形状比成桥状态时变化比较明显。
参考文献
[1]K．J．Bathe＆ S．Bofourchi，Large displacement analysis of three－dimensional beam structures，Int．J．Num．Meth．Eng，V14（1979（， 961－986
［2］K．M．Hsiao，H．J．Horng ＆Y．R．Chen，A corotational procedure that handle large rotations of Spstial beam structures Comput.Struct．27（1987），769－781
［3］G． Narayanan， C．S．Krihnamoorthy，An investigation of geometric non－linear formulation for 3－D beam element， Int． J．Non-linear Mechanics，1990，25（6）， 643-662
［4］陈政清，曾庆元，颜全胜．空间杆系结构大挠度问题内力分析 UL列式法．土木工程学报，1992，25（2），34-43
[5]黄文，李明瑞，黄文彬.杆系结构的几何非线性分析一Ⅱ三维问题．计算结构力学及其应用，1995，12（2），133-141
［6］潘家英，程庆国．大跨度悬索桥有限位移分析．土木工程学报，1994，27（1），l-10
