懸索橋空間纜索實用找形方法

馬召宇， 許福友， 檀永剛

(大連理工大學 建設工程學部，大連 116024)

1 引 言

隨著懸索橋跨度的增加，其橫向剛度和動力穩定性不斷降低，因此空間主纜得到廣泛采用，以提高懸索橋橫向剛度和穩定性，如韓國永宗橋和杭州江東大橋等。另外，一些較窄的人行懸索橋，其主梁剛度和抗風穩定性相對較差，需要加設空間風纜，如河北平山紅崖谷玻璃吊橋和日本九重夢大橋等。相比于平面纜索，空間纜索的線形計算應考慮主纜和吊桿的耦合效應，計算過程更為復雜。目前主要有分段懸鏈線法[1,2]SCM(Segmental Catenary Method)和有限元法兩類方法。

羅喜恒等[3]以空間的主纜索段為研究對象，將SCM擴展至空間纜索的線形計算。李傳習等[4]推導了懸索橋空間主纜線形變化剛度矩陣的遞推公式，提高了SCM的計算效率。文獻[5-8]基于空間分段懸鏈線理論對SCM加以修正，改善了其計算性能。SCM迭代參數少，求解速度快，是目前計算精度最高，應用最普遍的方法。但該方法對迭代初值敏感，迭代過程容易發散，且需要推導復雜的懸鏈線方程，編程求解非線性方程組，不便于一般工程師和科研人員使用。有限元法[9]更具有通用性，在方便應用、收斂性等方面相對SCM更有優勢，也可以實現空間纜索的線形計算。文獻[10,11]將懸鏈線理論和有限元法相結合計算纜索線形，計算過程過于復雜，耗時較多。Sun等[12,13]將計算拱橋合理拱軸線的思想用于懸索橋空間主纜的線形計算，提出了坐標迭代法。該方法收斂穩定，但仍需要編程實現，不便在通用有限元軟件中應用。Xiao 等[14,15]將吊桿力等效為集中荷載，建立主纜的有限元模型，并通過ANSYS進行簡單的代數運算和流程控制，實現了空間纜索的線形計算。該方法不需要復雜的理論推導和編程，但仍存在以下不足，(1) 當主纜跨度較大或空間效應明顯時，由拋物線法指定迭代初值時，主纜節點坐標可能與真實值偏差較大，容易導致第一步迭代計算發散，為確保收斂需多次人為調整迭代初值； (2) 當主纜坐標更新量較大時，由于耦合效應吊桿力會明顯變化，容易導致下一步迭代計算發散，為確保收斂，迭代時需采取增加臨時約束或減小坐標更新量等人為干預措施，不便于操作； (3) 忽略了吊桿的垂度效應。

針對以上問題，基于ANSYS提出一種懸索橋空間纜索的實用找形方法。本文方法建立空間纜索純索體系模型，以纜索節點坐標和單元內力為未知參數，采用小彈性模量法確定其初值，通過內循環參數更新外循環坐標修正的嵌套循環迭代計算空間纜索線形。通過2個算例驗證本文方法用于橫橋向和縱橋向斜吊桿空間纜索線形計算的精度和收斂的穩定性，并通過1個工程實例驗證了本文方法用于工程實踐的有效性。

2 空間纜索設計參數

懸索橋的風纜與空間主纜的受力特點和求解方法是一致的，但通常情況下，風纜橫向矢跨比大于空間主纜，故風纜的空間特性更明顯，線形計算難度更大，因此本文以懸索橋風纜為工程背景展開論述。如圖1所示定義空間坐標系，x軸為順橋向，y軸為橫橋向，z軸為豎向。風纜有y和z兩個方向的垂度，y方向垂度通常事先設定，z方向垂度需要計算；對于吊桿，在確定A1～An - 1坐標前，其線形和內力無法確定，B1～Bn - 1吊點處吊桿力y方向分量可以事先設定，z方向分量需要計算。因此，計算空間纜索線形時，設計參數如下。

圖1 風纜分析模型

(1) 風纜和吊桿的彈性模量、截面面積和材料密度等參數。

(2) 風纜A0和An兩錨固點x，y和z坐標，A1～An - 1吊點的x坐標，控制點Ai的y坐標ycon。

(3) 吊桿在主梁上吊點B1～Bn -1的x，y和z坐標。

(4)B1～Bn -1吊點處吊桿力y方向分量設計值F1 y～F(n - 1)y。

3 空間纜索線形計算

采用有限元計算纜索線形屬于典型的大變形小應變的幾何非線性問題，通常采用迭代法計算，而線形計算過程中，有限元建模方法、迭代初值的確定和迭代算法均十分關鍵，任何一步操作不合理都可能造成計算發散。

3.1 基本假定

為了簡化操作流程，避免復雜的理論推導，結合空間纜索的特點，本文基本假定如下。

(1) 采用桿單元(Link10)模擬風纜和吊桿，單元呈直線線形，且只能受拉而不能彎壓，單元材料符合胡克定律。

(2) 將沿單元均勻分布的自重荷載轉化為等效的集中荷載作用于單元節點上。

(3) 每根吊桿及風纜相鄰吊點間的索段均采用一個或多個桿單元模擬。

3.2 建模方法

有限元模型是基于ANSYS計算空間纜索線形的基礎，建模方法的優劣將直接關系到迭代過程的穩定性和結果的準確性。文獻[14，15]采用僅包含主纜的模型，每步迭代時，根據吊桿力豎向分量和主纜節點坐標更新吊桿力，并施加于主纜節點上。采用該方法建模，當主纜坐標更新量較大時，由于耦合效應，吊桿力會發生明顯變化，容易導致迭代發散。另外，將吊桿力等效為集中荷載忽略了吊桿的垂度效應，而引入懸鏈線理論又會增加算法的難度，不便于工程應用。針對以上問題，本文建模方法如下。

(1) 如圖1所示，建立風纜和吊桿單元且均采用Link10單元模擬。

(2) 固結風纜A0和An兩點，約束B1～Bn -1吊點x和z方向自由度，釋放y方向自由度。

(3) 將吊桿力y方向分量的設計值F1 y～F(n -1)y作為集中荷載施加于吊點B1～Bn -1上。

(4) 索夾重量采用質量單元施加于風纜節點A1～An -1上。

基于以上方法建立的模型迭代計算空間纜索線形具有以下明顯優勢。

(1) 模型建立了吊桿單元，每一步迭代時直接將上步吊桿力計算值更新為下次迭代初值，避免了吊桿力更新變化導致的迭代發散。

(2) 模型自身可以考慮吊桿的垂度效應以及風纜和吊桿的耦合效應。

(3) 模型除了約束A0和An兩點外，B1～Bn -1吊點的x和z方向均施加了約束，可以提高計算的穩定性。

3.3 迭代初值的確定

基于以上建模方法，以風纜和吊桿單元內力和節點坐標為未知參數，迭代計算空間纜索線形。迭代計算時，首先應賦予未知參數合理的迭代初值，否則很可能會造成第一步迭代計算的發散。文獻[6，7]采用拋物線法確定迭代初值，對于空間效應明顯的風纜，該方法確定的初值往往與真實值相差較大，容易導致發散，往往需要多次人為調整，不便于實際應用。為解決以上問題，采用小彈性模量法計算確定未知參數的迭代初值。所謂小彈性模量法，即在有限元計算時將索單元的彈性模量取小于真實值的2～3個數量級，使索單元獲得較強的變形能力和穩定的收斂性。基本步驟如下。

(1) 如圖2所示，以A0An直線為風纜初始線形建立初始有限元模型。

圖2 初始有限元模型

(2) 賦予單元的截面積、線密度和小彈性模量(縮小103倍)。

(3) 賦予風纜和吊桿各單元相同且較大的初應變。

(4) 非線性求解，獲得風纜和吊桿單元內力以及風纜節點y和z坐標。

(5) 計算控制點Ai的y坐標與設計值ycon的差值Δ1，若Δ1>Δcon1 -1，返回步驟(4)；若|Δ1|<Δcon1 -1，則輸出風纜和吊桿單元內力和風纜節點y和z坐標。

上述步驟確定的未知參數初值只需確保第一步迭代計算的收斂，對參數取值的精度要求不高，因此步驟(5)收斂容許值Δcon1 -1可以在0.1 m～0.5 m 范圍內取值。工程人員一般試算2～3次即可確定步驟(3)合理的初應變值，進而獲得滿足收斂要求的迭代初值。相比于拋物線法，本文方法具有以下優勢。

(1) 迭代初值由有限元計算確定，滿足小彈性模量下力的平衡條件，相比于由拋物線法確定的初值更接近于真實值，故收斂性更好，不需要多次人為調整。

(2) 初始模型風纜和吊桿的線形均指定為直線，操作簡便快捷，便于掌握。

3.4 迭代算法

基于有限元迭代計算空間纜索線形時，迭代方法的優劣將直接關系到計算的斂散性和求解速度的快慢。本文以風纜和吊桿各單元內力及節點坐標為未知參數，采用嵌套循環展開迭代計算，基本步驟如下。

(1) 根據3.3節輸出結果更新有限元模型，恢復索單元的真實彈性模量。

(2) 非線性求解，獲得風纜和吊桿各單元內力及各節點y和z坐標。

(3) 如圖3所示，根據風纜控制點Ai的y坐標設計要求值ycon，對上次計算獲得的風纜y坐標a進行修正，設c為過A0，(xi，Δ1)，An三點的拋物線，將a減去c即為修正后的風纜y坐標b。

圖3 風纜節點y坐標修正簡圖

(4) 更新各單元內力，以及風纜和吊桿各節點y和z坐標。

(5) 非線性求解，獲得風纜和吊桿各單元內力及各節點y和z坐標。

(6) 獲得B1～Bn -1節點y方向位移Δ2，若|Δ2|>Δcon 2，返回步驟(4)；若|Δ2|<Δcon 2，則進入步驟(7)。

(7) 計算控制點Ai的y坐標與設計值ycon的差值Δ1，若|Δ1|>Δcon1 - 2，返回步驟(3)；若|Δ1|<Δcon1 -2，則結束迭代，輸出結果。

步驟(6，7)的收斂容許值Δcon 2和Δcon1 -2一般取10-6m即可滿足工程精度要求，以上迭代算法在ANSYS中基于APDL語言實現。

本文迭代算法的內循環通過風纜和吊桿各單元內力及各節點y和z坐標的更新，確保B1～Bn -1節點y方向位移快速收斂于Δcon 2，內循環結束時|Δ1|可能不滿足收斂條件，即風纜控制點Ai的y坐標與設計值ycon相差較大，故進入外循環進行風纜y坐標修正，由于坐標修正打破了上次內循環迭代結束時纜索的平衡狀態，因此需要再次進入內循環迭代更新，由此內外循環依次交替進行，|Δ1|和|Δ2|即可快速滿足收斂條件，進而獲得滿足已知設計條件的空間纜索線形。

4 算 例

采用上述方法對3個算例進行分析。其中算例1對橫向斜吊桿空間纜索進行分析，并通過與SCM的計算結果進行對比，驗證了本文方法的精度；基于算例1，進一步開展了24種工況的線形計算，并通過與SCM的對比，驗證了本文方法的穩定性。算例2對縱橋向斜吊桿空間纜索進行了分析，并通過有限元軟件Midas驗證了本文方法計算縱橋向斜吊桿空間纜索線形的準確性。算例3對一座人行懸索橋風纜找形分析，驗證了本文方法在工程實踐中的有效性。

4.1 橫向斜吊桿空間纜索

如圖1所示，風纜兩錨固點坐標分別為A0(-100，0，0) m和A15(80，25，20) m，吊桿吊點B1(-65，100，60) m和B14(65，100，60) m，控制點A8的y坐標為60 m，吊桿沿x方向間距為 10 m，吊桿在B1～B14點沿y方向的張拉力設計值為 45 kN。風纜和吊桿的截面積、彈性模量和線密度分別為62.80 cm2，158 GPa，528.78 N/m和 5.22 cm2，165 GPa，45.21 N/m，計算該空間纜索系統的線形和內力。

風纜和吊桿均采用Link10單元模擬，風纜沿x方向每2.5 m劃分1個單元，共72個單元，每根吊桿劃分5個單元，整個纜索系統共142個單元。本算例外循環共迭代10次，內循環累計迭代 92次，完成迭代。表1給出了本文方法確定的風纜節點初始坐標及最終坐標。為了便于對比驗證，根據SCM編制Matlab程序分析算例1，SCM的計算結果以及本文方法與SCM的計算差值一并給出。

由表1可知，兩種方法計算得到的風纜各節點y坐標幾乎完全相等，z坐標最大相差0.11 mm，可見兩種方法計算結果相差很小，說明本文方法有效可行,且精度與SCM一致。

表1 風纜節點坐標Tab.1 Coordinates of wind cable nodes

為了驗證本文方法的計算效率和穩定性，以算例1為基礎，將風纜控制點A8的y坐標、吊桿橫橋向張拉力、A15的y坐標和A15的z坐標分別設為 6種不同的數值，并在修改其中一個參數時，其余參數取值不變，共計得到24種工況，列入表2。采用本文方法對以上24種工況下的空間纜索線形進行計算，表2給出了內循環迭代次數，括號內第1個數值為采用本文方法的循環迭代次數，第2個數值為SCM的迭代次數。為了便于對比驗證，SCM迭代次數也一并給出。

由表2可知，采用2種方法計算以上24種工況均全部收斂。SCM迭代次數較多，且迭代過程中為了保證收斂，需要多次調整迭代參數初值和修改松弛因子等人為參與。本文方法確定的初值不需人為調整即可保證第一步迭代的收斂，迭代過程中不需人為參與，具有良好的收斂性。可見本文方法人為參與程度低，計算效率高，具有很強的穩定性和可靠性，可以減少計算時間。

4.2 縱橋向斜吊桿空間纜索

如圖4所示，目前一些懸索橋的風纜吊桿采用縱向傾斜的布置形式，為驗證本文方法對縱向斜吊桿空間纜索線形計算的適用性，增加以下算例。風纜兩錨固點A0和A9坐標分別為(-100，0，0) m和(80，25，20) m，吊點B1和B7坐標分別為(-65，100，60) m和(55，100，60) m，控制點A5的y坐標為60 m，吊點A1和A8的x坐標分別為-75 m 和65 m，吊點沿x方向間距為20 m，吊桿在B1～B7點沿y方向張拉力設計值為45 kN。風纜和吊桿的截面積、彈性模量和線密度分別為 62.80 cm2，158 GPa，528.78 N/m 和 5.22 cm2，165 GPa，45.21 N/m，計算該空間纜索系統的線形和內力。

圖4 縱橋斜吊桿風纜分析模型

對于縱橋向斜吊桿的空間纜索，為使B1～B7吊點處x方向吊桿力相互平衡，采用本文方法開展迭代計算時，有限元模型應釋放B1～B7吊點x方向的約束，本算例外循環共迭代10次，內循環累計迭代120次，完成迭代。表3給出了采用本文方法確定的風纜節點初始坐標及最終坐標。

表3 風纜節點坐標

根據本文方法計算得到的風纜和吊桿各單元無應力長度及各節點坐標，建立該空間纜索純索體系的Midas模型，進行靜力分析。由Midas模型靜力分析得到的各節點位移均小于10-6m，且B1～B7吊點橫橋向反力為45 kN，說明本文方法得到的計算結果是準確的，進而表明采用本文方法計算縱橋向斜吊桿空間纜索線形的可行性。

4.3 工程實例

一座單跨地錨式人行懸索橋，跨度為445 m。主纜理論矢高為45.0 m，矢跨比為1/9.89(圖5)，為提高橋梁抗風性能，橋面下方增設抗風纜。風纜具有以下設計要求，(1) 風纜和吊桿的截面積、彈性模量和線密度分別為 62.80 cm2，158 GPa，528.78 N/m 和5.22 cm2，165 GPa，45.21 N/m； (2) 上游側風纜橫橋向矢高為34.2 m，矢跨比為1/12.28，兩錨固點順橋向豎向和橫橋向坐標分別為(-210，-25，40) m和(210，-25，40) m； (3) 風纜共設23根順橋向間距為13.32 m的吊桿，吊桿在主梁上吊點的豎向和橫橋向坐標為(-5.896，1.800) m，且橫橋向張拉力為45 kN。依據以上設計要求，通過本文方法，風纜線形的計算列入表4。

圖5 設風纜懸索橋平面及立面布置

表4 風纜節點坐標

工程師在設計懸索橋風纜時，通常需要多次試算風纜吊桿橫橋向張拉力，每次試算時風纜線形的計算都有較大挑戰，尤其對于圖6所示吊桿縱橋向傾斜設置的風纜更為突出。通過本文方法計算空間纜索線形，工程人員只需給出風纜設計參數，操作簡便、快捷且穩定，為設置空間纜索的懸索橋的設計、施工和監控提供參考和借鑒。

圖6 縱橋向風纜斜吊桿懸索橋平面及立面布置

5 結 論

(1) 基于通用有限元軟件ANSYS提出了一種懸索橋空間纜索線形計算的實用方法，并給出了建立模型、確定迭代初值和迭代計算的詳細步驟。

(2) 本文方法計算精度與SCM一致，但比SCM迭代次數少、穩定性好、計算效率更高，還適用于縱向斜吊桿空間纜索的線形計算。

(3) 本文方法操作簡便快捷，不需要復雜的理論推導，人為參與程度低，具有較強的穩定性和可靠性，可供設計、施工和科研部門參考使用。