馬召宇, 許福友, 檀永剛
(大連理工大學(xué) 建設(shè)工程學(xué)部,大連 116024)
隨著懸索橋跨度的增加,其橫向剛度和動力穩(wěn)定性不斷降低,因此空間主纜得到廣泛采用,以提高懸索橋橫向剛度和穩(wěn)定性,如韓國永宗橋和杭州江東大橋等。另外,一些較窄的人行懸索橋,其主梁剛度和抗風(fēng)穩(wěn)定性相對較差,需要加設(shè)空間風(fēng)纜,如河北平山紅崖谷玻璃吊橋和日本九重夢大橋等。相比于平面纜索,空間纜索的線形計算應(yīng)考慮主纜和吊桿的耦合效應(yīng),計算過程更為復(fù)雜。目前主要有分段懸鏈線法[1,2]SCM(Segmental Catenary Method)和有限元法兩類方法。
羅喜恒等[3]以空間的主纜索段為研究對象,將SCM擴展至空間纜索的線形計算。李傳習(xí)等[4]推導(dǎo)了懸索橋空間主纜線形變化剛度矩陣的遞推公式,提高了SCM的計算效率。文獻[5-8]基于空間分段懸鏈線理論對SCM加以修正,改善了其計算性能。SCM迭代參數(shù)少,求解速度快,是目前計算精度最高,應(yīng)用最普遍的方法。但該方法對迭代初值敏感,迭代過程容易發(fā)散,且需要推導(dǎo)復(fù)雜的懸鏈線方程,編程求解非線性方程組,不便于一般工程師和科研人員使用。有限元法[9]更具有通用性,在方便應(yīng)用、收斂性等方面相對SCM更有優(yōu)勢,也可以實現(xiàn)空間纜索的線形計算。文獻[10,11]將懸鏈線理論和有限元法相結(jié)合計算纜索線形,計算過程過于復(fù)雜,耗時較多。Sun等[12,13]將計算拱橋合理拱軸線的思想用于懸索橋空間主纜的線形計算,提出了坐標迭代法。該方法收斂穩(wěn)定,但仍需要編程實現(xiàn),不便在通用有限元軟件中應(yīng)用。Xiao 等[14,15]將吊桿力等效為集中荷載,建立主纜的有限元模型,并通過ANSYS進行簡單的代數(shù)運算和流程控制,實現(xiàn)了空間纜索的線形計算。該方法不需要復(fù)雜的理論推導(dǎo)和編程,但仍存在以下不足,(1) 當(dāng)主纜跨度較大或空間效應(yīng)明顯時,由拋物線法指定迭代初值時,主纜節(jié)點坐標可能與真實值偏差較大,容易導(dǎo)致第一步迭代計算發(fā)散,為確保收斂需多次人為調(diào)整迭代初值; (2) 當(dāng)主纜坐標更新量較大時,由于耦合效應(yīng)吊桿力會明顯變化,容易導(dǎo)致下一步迭代計算發(fā)散,為確保收斂,迭代時需采取增加臨時約束或減小坐標更新量等人為干預(yù)措施,不便于操作; (3) 忽略了吊桿的垂度效應(yīng)。
針對以上問題,基于ANSYS提出一種懸索橋空間纜索的實用找形方法。本文方法建立空間纜索純索體系模型,以纜索節(jié)點坐標和單元內(nèi)力為未知參數(shù),采用小彈性模量法確定其初值,通過內(nèi)循環(huán)參數(shù)更新外循環(huán)坐標修正的嵌套循環(huán)迭代計算空間纜索線形。通過2個算例驗證本文方法用于橫橋向和縱橋向斜吊桿空間纜索線形計算的精度和收斂的穩(wěn)定性,并通過1個工程實例驗證了本文方法用于工程實踐的有效性。
懸索橋的風(fēng)纜與空間主纜的受力特點和求解方法是一致的,但通常情況下,風(fēng)纜橫向矢跨比大于空間主纜,故風(fēng)纜的空間特性更明顯,線形計算難度更大,因此本文以懸索橋風(fēng)纜為工程背景展開論述。如圖1所示定義空間坐標系,x軸為順橋向,y軸為橫橋向,z軸為豎向。風(fēng)纜有y和z兩個方向的垂度,y方向垂度通常事先設(shè)定,z方向垂度需要計算;對于吊桿,在確定A1~An - 1坐標前,其線形和內(nèi)力無法確定,B1~Bn - 1吊點處吊桿力y方向分量可以事先設(shè)定,z方向分量需要計算。因此,計算空間纜索線形時,設(shè)計參數(shù)如下。
圖1 風(fēng)纜分析模型
(1) 風(fēng)纜和吊桿的彈性模量、截面面積和材料密度等參數(shù)。
(2) 風(fēng)纜A0和An兩錨固點x,y和z坐標,A1~An - 1吊點的x坐標,控制點Ai的y坐標ycon。
(3) 吊桿在主梁上吊點B1~Bn -1的x,y和z坐標。
(4)B1~Bn -1吊點處吊桿力y方向分量設(shè)計值F1 y~F(n - 1)y。
采用有限元計算纜索線形屬于典型的大變形小應(yīng)變的幾何非線性問題,通常采用迭代法計算,而線形計算過程中,有限元建模方法、迭代初值的確定和迭代算法均十分關(guān)鍵,任何一步操作不合理都可能造成計算發(fā)散。
為了簡化操作流程,避免復(fù)雜的理論推導(dǎo),結(jié)合空間纜索的特點,本文基本假定如下。
(1) 采用桿單元(Link10)模擬風(fēng)纜和吊桿,單元呈直線線形,且只能受拉而不能彎壓,單元材料符合胡克定律。
(2) 將沿單元均勻分布的自重荷載轉(zhuǎn)化為等效的集中荷載作用于單元節(jié)點上。
(3) 每根吊桿及風(fēng)纜相鄰吊點間的索段均采用一個或多個桿單元模擬。
有限元模型是基于ANSYS計算空間纜索線形的基礎(chǔ),建模方法的優(yōu)劣將直接關(guān)系到迭代過程的穩(wěn)定性和結(jié)果的準確性。文獻[14,15]采用僅包含主纜的模型,每步迭代時,根據(jù)吊桿力豎向分量和主纜節(jié)點坐標更新吊桿力,并施加于主纜節(jié)點上。采用該方法建模,當(dāng)主纜坐標更新量較大時,由于耦合效應(yīng),吊桿力會發(fā)生明顯變化,容易導(dǎo)致迭代發(fā)散。另外,將吊桿力等效為集中荷載忽略了吊桿的垂度效應(yīng),而引入懸鏈線理論又會增加算法的難度,不便于工程應(yīng)用。針對以上問題,本文建模方法如下。
(1) 如圖1所示,建立風(fēng)纜和吊桿單元且均采用Link10單元模擬。
(2) 固結(jié)風(fēng)纜A0和An兩點,約束B1~Bn -1吊點x和z方向自由度,釋放y方向自由度。
(3) 將吊桿力y方向分量的設(shè)計值F1 y~F(n -1)y作為集中荷載施加于吊點B1~Bn -1上。
(4) 索夾重量采用質(zhì)量單元施加于風(fēng)纜節(jié)點A1~An -1上。
基于以上方法建立的模型迭代計算空間纜索線形具有以下明顯優(yōu)勢。
(1) 模型建立了吊桿單元,每一步迭代時直接將上步吊桿力計算值更新為下次迭代初值,避免了吊桿力更新變化導(dǎo)致的迭代發(fā)散。
(2) 模型自身可以考慮吊桿的垂度效應(yīng)以及風(fēng)纜和吊桿的耦合效應(yīng)。
(3) 模型除了約束A0和An兩點外,B1~Bn -1吊點的x和z方向均施加了約束,可以提高計算的穩(wěn)定性。
基于以上建模方法,以風(fēng)纜和吊桿單元內(nèi)力和節(jié)點坐標為未知參數(shù),迭代計算空間纜索線形。迭代計算時,首先應(yīng)賦予未知參數(shù)合理的迭代初值,否則很可能會造成第一步迭代計算的發(fā)散。文獻[6,7]采用拋物線法確定迭代初值,對于空間效應(yīng)明顯的風(fēng)纜,該方法確定的初值往往與真實值相差較大,容易導(dǎo)致發(fā)散,往往需要多次人為調(diào)整,不便于實際應(yīng)用。為解決以上問題,采用小彈性模量法計算確定未知參數(shù)的迭代初值。所謂小彈性模量法,即在有限元計算時將索單元的彈性模量取小于真實值的2~3個數(shù)量級,使索單元獲得較強的變形能力和穩(wěn)定的收斂性?;静襟E如下。
(1) 如圖2所示,以A0An直線為風(fēng)纜初始線形建立初始有限元模型。
圖2 初始有限元模型
(2) 賦予單元的截面積、線密度和小彈性模量(縮小103倍)。
(3) 賦予風(fēng)纜和吊桿各單元相同且較大的初應(yīng)變。
(4) 非線性求解,獲得風(fēng)纜和吊桿單元內(nèi)力以及風(fēng)纜節(jié)點y和z坐標。
(5) 計算控制點Ai的y坐標與設(shè)計值ycon的差值Δ1,若Δ1>Δcon1 -1,返回步驟(4);若|Δ1|<Δcon1 -1,則輸出風(fēng)纜和吊桿單元內(nèi)力和風(fēng)纜節(jié)點y和z坐標。
上述步驟確定的未知參數(shù)初值只需確保第一步迭代計算的收斂,對參數(shù)取值的精度要求不高,因此步驟(5)收斂容許值Δcon1 -1可以在0.1 m~0.5 m 范圍內(nèi)取值。工程人員一般試算2~3次即可確定步驟(3)合理的初應(yīng)變值,進而獲得滿足收斂要求的迭代初值。相比于拋物線法,本文方法具有以下優(yōu)勢。
(1) 迭代初值由有限元計算確定,滿足小彈性模量下力的平衡條件,相比于由拋物線法確定的初值更接近于真實值,故收斂性更好,不需要多次人為調(diào)整。
(2) 初始模型風(fēng)纜和吊桿的線形均指定為直線,操作簡便快捷,便于掌握。
基于有限元迭代計算空間纜索線形時,迭代方法的優(yōu)劣將直接關(guān)系到計算的斂散性和求解速度的快慢。本文以風(fēng)纜和吊桿各單元內(nèi)力及節(jié)點坐標為未知參數(shù),采用嵌套循環(huán)展開迭代計算,基本步驟如下。
(1) 根據(jù)3.3節(jié)輸出結(jié)果更新有限元模型,恢復(fù)索單元的真實彈性模量。
(2) 非線性求解,獲得風(fēng)纜和吊桿各單元內(nèi)力及各節(jié)點y和z坐標。
(3) 如圖3所示,根據(jù)風(fēng)纜控制點Ai的y坐標設(shè)計要求值ycon,對上次計算獲得的風(fēng)纜y坐標a進行修正,設(shè)c為過A0,(xi,Δ1),An三點的拋物線,將a減去c即為修正后的風(fēng)纜y坐標b。
圖3 風(fēng)纜節(jié)點y坐標修正簡圖
(4) 更新各單元內(nèi)力,以及風(fēng)纜和吊桿各節(jié)點y和z坐標。
(5) 非線性求解,獲得風(fēng)纜和吊桿各單元內(nèi)力及各節(jié)點y和z坐標。
(6) 獲得B1~Bn -1節(jié)點y方向位移Δ2,若|Δ2|>Δcon 2,返回步驟(4);若|Δ2|<Δcon 2,則進入步驟(7)。
(7) 計算控制點Ai的y坐標與設(shè)計值ycon的差值Δ1,若|Δ1|>Δcon1 - 2,返回步驟(3);若|Δ1|<Δcon1 -2,則結(jié)束迭代,輸出結(jié)果。
步驟(6,7)的收斂容許值Δcon 2和Δcon1 -2一般取10-6m即可滿足工程精度要求,以上迭代算法在ANSYS中基于APDL語言實現(xiàn)。
本文迭代算法的內(nèi)循環(huán)通過風(fēng)纜和吊桿各單元內(nèi)力及各節(jié)點y和z坐標的更新,確保B1~Bn -1節(jié)點y方向位移快速收斂于Δcon 2,內(nèi)循環(huán)結(jié)束時|Δ1|可能不滿足收斂條件,即風(fēng)纜控制點Ai的y坐標與設(shè)計值ycon相差較大,故進入外循環(huán)進行風(fēng)纜y坐標修正,由于坐標修正打破了上次內(nèi)循環(huán)迭代結(jié)束時纜索的平衡狀態(tài),因此需要再次進入內(nèi)循環(huán)迭代更新,由此內(nèi)外循環(huán)依次交替進行,|Δ1|和|Δ2|即可快速滿足收斂條件,進而獲得滿足已知設(shè)計條件的空間纜索線形。
采用上述方法對3個算例進行分析。其中算例1對橫向斜吊桿空間纜索進行分析,并通過與SCM的計算結(jié)果進行對比,驗證了本文方法的精度;基于算例1,進一步開展了24種工況的線形計算,并通過與SCM的對比,驗證了本文方法的穩(wěn)定性。算例2對縱橋向斜吊桿空間纜索進行了分析,并通過有限元軟件Midas驗證了本文方法計算縱橋向斜吊桿空間纜索線形的準確性。算例3對一座人行懸索橋風(fēng)纜找形分析,驗證了本文方法在工程實踐中的有效性。
如圖1所示,風(fēng)纜兩錨固點坐標分別為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方向的張拉力設(shè)計值為 45 kN。風(fēng)纜和吊桿的截面積、彈性模量和線密度分別為62.80 cm2,158 GPa,528.78 N/m和 5.22 cm2,165 GPa,45.21 N/m,計算該空間纜索系統(tǒng)的線形和內(nèi)力。
風(fēng)纜和吊桿均采用Link10單元模擬,風(fēng)纜沿x方向每2.5 m劃分1個單元,共72個單元,每根吊桿劃分5個單元,整個纜索系統(tǒng)共142個單元。本算例外循環(huán)共迭代10次,內(nèi)循環(huán)累計迭代 92次,完成迭代。表1給出了本文方法確定的風(fēng)纜節(jié)點初始坐標及最終坐標。為了便于對比驗證,根據(jù)SCM編制Matlab程序分析算例1,SCM的計算結(jié)果以及本文方法與SCM的計算差值一并給出。
由表1可知,兩種方法計算得到的風(fēng)纜各節(jié)點y坐標幾乎完全相等,z坐標最大相差0.11 mm,可見兩種方法計算結(jié)果相差很小,說明本文方法有效可行,且精度與SCM一致。
表1 風(fēng)纜節(jié)點坐標Tab.1 Coordinates of wind cable nodes
為了驗證本文方法的計算效率和穩(wěn)定性,以算例1為基礎(chǔ),將風(fēng)纜控制點A8的y坐標、吊桿橫橋向張拉力、A15的y坐標和A15的z坐標分別設(shè)為 6種不同的數(shù)值,并在修改其中一個參數(shù)時,其余參數(shù)取值不變,共計得到24種工況,列入表2。采用本文方法對以上24種工況下的空間纜索線形進行計算,表2給出了內(nèi)循環(huán)迭代次數(shù),括號內(nèi)第1個數(shù)值為采用本文方法的循環(huán)迭代次數(shù),第2個數(shù)值為SCM的迭代次數(shù)。為了便于對比驗證,SCM迭代次數(shù)也一并給出。
由表2可知,采用2種方法計算以上24種工況均全部收斂。SCM迭代次數(shù)較多,且迭代過程中為了保證收斂,需要多次調(diào)整迭代參數(shù)初值和修改松弛因子等人為參與。本文方法確定的初值不需人為調(diào)整即可保證第一步迭代的收斂,迭代過程中不需人為參與,具有良好的收斂性??梢姳疚姆椒ㄈ藶閰⑴c程度低,計算效率高,具有很強的穩(wěn)定性和可靠性,可以減少計算時間。
如圖4所示,目前一些懸索橋的風(fēng)纜吊桿采用縱向傾斜的布置形式,為驗證本文方法對縱向斜吊桿空間纜索線形計算的適用性,增加以下算例。風(fēng)纜兩錨固點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方向張拉力設(shè)計值為45 kN。風(fēng)纜和吊桿的截面積、彈性模量和線密度分別為 62.80 cm2,158 GPa,528.78 N/m 和 5.22 cm2,165 GPa,45.21 N/m,計算該空間纜索系統(tǒng)的線形和內(nèi)力。
圖4 縱橋斜吊桿風(fēng)纜分析模型
對于縱橋向斜吊桿的空間纜索,為使B1~B7吊點處x方向吊桿力相互平衡,采用本文方法開展迭代計算時,有限元模型應(yīng)釋放B1~B7吊點x方向的約束,本算例外循環(huán)共迭代10次,內(nèi)循環(huán)累計迭代120次,完成迭代。表3給出了采用本文方法確定的風(fēng)纜節(jié)點初始坐標及最終坐標。
表3 風(fēng)纜節(jié)點坐標
根據(jù)本文方法計算得到的風(fēng)纜和吊桿各單元無應(yīng)力長度及各節(jié)點坐標,建立該空間纜索純索體系的Midas模型,進行靜力分析。由Midas模型靜力分析得到的各節(jié)點位移均小于10-6m,且B1~B7吊點橫橋向反力為45 kN,說明本文方法得到的計算結(jié)果是準確的,進而表明采用本文方法計算縱橋向斜吊桿空間纜索線形的可行性。
一座單跨地錨式人行懸索橋,跨度為445 m。主纜理論矢高為45.0 m,矢跨比為1/9.89(圖5),為提高橋梁抗風(fēng)性能,橋面下方增設(shè)抗風(fēng)纜。風(fēng)纜具有以下設(shè)計要求,(1) 風(fēng)纜和吊桿的截面積、彈性模量和線密度分別為 62.80 cm2,158 GPa,528.78 N/m 和5.22 cm2,165 GPa,45.21 N/m; (2) 上游側(cè)風(fēng)纜橫橋向矢高為34.2 m,矢跨比為1/12.28,兩錨固點順橋向豎向和橫橋向坐標分別為(-210,-25,40) m和(210,-25,40) m; (3) 風(fēng)纜共設(shè)23根順橋向間距為13.32 m的吊桿,吊桿在主梁上吊點的豎向和橫橋向坐標為(-5.896,1.800) m,且橫橋向張拉力為45 kN。依據(jù)以上設(shè)計要求,通過本文方法,風(fēng)纜線形的計算列入表4。
圖5 設(shè)風(fēng)纜懸索橋平面及立面布置
表4 風(fēng)纜節(jié)點坐標
工程師在設(shè)計懸索橋風(fēng)纜時,通常需要多次試算風(fēng)纜吊桿橫橋向張拉力,每次試算時風(fēng)纜線形的計算都有較大挑戰(zhàn),尤其對于圖6所示吊桿縱橋向傾斜設(shè)置的風(fēng)纜更為突出。通過本文方法計算空間纜索線形,工程人員只需給出風(fēng)纜設(shè)計參數(shù),操作簡便、快捷且穩(wěn)定,為設(shè)置空間纜索的懸索橋的設(shè)計、施工和監(jiān)控提供參考和借鑒。
圖6 縱橋向風(fēng)纜斜吊桿懸索橋平面及立面布置
(1) 基于通用有限元軟件ANSYS提出了一種懸索橋空間纜索線形計算的實用方法,并給出了建立模型、確定迭代初值和迭代計算的詳細步驟。
(2) 本文方法計算精度與SCM一致,但比SCM迭代次數(shù)少、穩(wěn)定性好、計算效率更高,還適用于縱向斜吊桿空間纜索的線形計算。
(3) 本文方法操作簡便快捷,不需要復(fù)雜的理論推導(dǎo),人為參與程度低,具有較強的穩(wěn)定性和可靠性,可供設(shè)計、施工和科研部門參考使用。