橢圓鞍形單層正交索網結構解析計算方法與形態優化

王博厚，孫 遜，方立新

(1.東南大學土木工程學院，江蘇，南京 210096；2.東南大學建筑設計研究院有限公司，江蘇，南京 210096；3.東南大學建筑學院，江蘇，南京 210096)

索網結構因其在充分利用材料強度實現結構功能的同時，能保持整體結構的輕巧簡潔，能適應建筑功能及外觀需求等優點，受到工程結構界的廣泛關注。單層正交索網，由相互正交、曲率相反的兩組鋼索連接而形成，形成負高斯曲率的曲面。索網懸掛在強大的邊緣構件(外環梁)上。單層正交索網可按其平面投影形狀分類，如平面投影為菱形、橢圓形、圓形的單層正交索網結構。

1952 年建成的道頓競技館被公認為是第一個具有現代意義的大跨度索網屋蓋[1]。隨后在國內外的大型場館、公共建筑中得到廣泛應用，如倫敦奧運自行車館[2]、蘇州奧體中心游泳館[3]、國家速滑館[4]等。對于索網結構的研究，主要包括找形與荷載分析[5]、解析計算[6-7]、靜力與動力性能研究[8]、形態優化[9-11]、斷索分析[12]、施工模擬等方面。

用于索網結構初始形態確定和荷載狀態分析的主要方法有：非線性有限元法[13]、動力松弛法[14]和力密度法[15]。而非線性有限元法是一種更為精確的數值計算方法[16]，且隨著有限元計算軟件的發展，該方法在找形及荷載分析中應用最廣。

沈世釗等[6]、金問魯[7]將鞍形索網連續化為薄膜，分別利用平衡與變形協調關系、能量原理建立以豎向位移為未知量的非線性方程，其所得解析解精度極大程度依賴于所假定的位移函數。

肖康麗等[8]對索網結構靜力性能進行參數分析，研究矢跨比、邊界形狀、索初始預張力、索直徑等參數對變形、內力的影響，其結論可用于指導結構優化設計。目前在索膜結構設計中，大多根據建筑造型要求，在所允許的范圍內調整幾何參數以確定最終形態，在此過程中雖然會進行一定的方案比選，但該過程大多依賴于設計者的經驗，很難保證最終形態在給定約束條件下最優[17]。

對索網結構形態的優化，關鍵問題在于選取合適的優化目標。形態優化研究根據優化目標可分為單目標優化、多目標優化，目前學者所提出的優化目標可概括為材料用量最少、豎向剛度最大兩大類[9-11]，對于材料用料的優化主要針對索和下部支撐結構。已有形態優化方法需借助有限元、優化算法自編程序實現，計算過程繁雜且不易收斂，不能快速判斷結構形態是否合理。

索網屋蓋結構中，外環梁材料用量大、受力復雜，其截面極大地影響著結構的經濟性和安全性，在優化設計中未受到廣泛關注。外環梁大多處壓彎狀態，重力荷載方向的彎矩由支承條件決定，而鞍面的形態對外環梁的面內彎矩有極大影響。外環梁越接近軸壓構件，其材料利用越充分，則可相應地減小外環截面尺寸和變形[18]。

針對以上問題，本文首先利用平衡方程、變形協調方程，以豎向位移為未知量推導了正交單層索網連續化計算方法，提出一個新的位移函數，采用伽遼金法求解。其次以外環梁平面內彎矩最小為優化目標，以中央承重索垂度為優化變量，以索力、豎向位移和矢跨比為約束條件。利用外環梁上任一點彎矩計算公式，求解方程得到結構最優形態。最后對一個簡單算例分別采用有限元法、解析法計算，驗證位移、索力、外環彎矩計算公式和形態優化的準確性。在結構初步設計階段，利用本方法僅需求解一元方程，便可快速得到結構的最優形態，避免了有限元編程、反復試算的過程。

1 索網連續化計算方法

目前，單層正交索網的計算理論主要有離散化理論和連續化理論[6]。連續化理論通過假設索網在受力變形后的位移函數并以之為未知量，進而求出索網面內的索力分布，避免了離散化理論中大量方程、未知量造成的求解困難。

本文利用連續化理論，通過平衡微分方程、變形協調方程建立以豎向位移為未知量的方程，以變分原理求解，推導正交索網的解析公式。

1.1 基本假定

索網計算時采用如下基本假設：① 索是理想柔性的，只能承受拉力；② 索材料符合胡克定律；③ 兩個方向的索之間永遠保持接觸，能相互傳遞豎向力的作用；④ 只考慮小垂度問題，且僅考慮豎向均布荷載作用；⑤x方向索為承重索，y方向索為穩定索，相同方向索的初始預拉力、截面面積、間距均相同。

1.2 索網的三種形態

根據索受力狀態和外荷載的不同，索網結構分為3 種狀態：零狀態、初始平衡態、荷載態。假設結構在初始態下的曲面形狀為z0(x,y)，索自重折算為豎向均布荷載q0(x,y)，兩個方向索折算到單位寬度內的薄膜預張力為Hx0、Hy0，折算到單位寬度內的拉伸剛度為EAx、EAy。

假設對結構施加的外荷載為Δq(x,y)，豎向荷載變為q=q0+Δq。索網在外荷載作用下產生三個方向的位移分量u、v、w，但在小垂度情況下水平位移分量u、v與豎向位移分量w相比是高階微量，可以忽略其對曲面形狀的影響，變形后索網形狀變為z=z0+w。索水平張力獲得增量ΔHx、ΔHy，達到Hx=Hx0+ΔHx、Hy=Hy0+ΔHy。

索網結構在初始態時的各項參數均為已知，在荷載態時兩正交方向的水平拉力Hx(或ΔHx)、Hy(或ΔHy)及位形z(或w)是未知的，求解三個未知量需要三個方程，即x、y方向的變形協調方程和z方向的平衡微分方程。索在受到外荷載作用從初始態轉換為荷載態的過程中，索中應變增量很小，但位移較大。這時的平衡條件應建立在變形后的位形上，應變表達式也應包含位移的二次項，即在建立方程時應考慮幾何非線性的影響。

1.3 索網平衡方程

在荷載態，從連續化的索網中任取微分單元dxdy，受力情況如圖1 所示。

圖1 索網微分單元受力簡圖Fig.1 Stress diagram of cable-net differential element

考慮x、y方向的平衡條件，可得Hx=Hx(y)、Hy=Hy(x)，即Hx在x方向不變，Hy在y方向不變。再考慮z方向平衡條件并將以上結論代入，可得：

式(1)即為荷載態時索網的平衡微分方程。同理可得初始態時索網的平衡微分方程：

若將q=q0+Δq、z=z0+w、Hx=Hx0+ΔHx、Hy=Hy0+ΔHy代入式(1)，并減去式(2)，可得索網平衡方程的增量形式：

1.4 索網變形協調方程

以x方向任意位置的承重索MN為例，假設MN從初始態轉換為荷載態時，兩端產生的位移為(uM、vM、wM)、(uN、vN、wN)，索溫度變化為Δtx。

1) 由幾何關系推導索長變化

考察在初始態長度為dsx0的微元AB，在荷載態移動到A′B′，其長度伸長至dsx，如圖2 所示。

圖2 由幾何關系計算索長簡圖Fig.2 Schematics of calculating cable length from geometric relations

由幾何關系可知，微元AB的伸長量為：

在小垂度問題中，略去水平位移u、v的高階微量，將根號利用級數展開并保留級數前兩項，并將z=z0+w代入，可得：

承重索MN的總伸長量為：

式中：M、N為承重索的兩個邊界點；uN、uM表示這兩點在x方向的位移。

2) 由物理關系推導索長變化

由物理關系可知，索的伸長由拉伸變形和溫度變化引起，即：

式中：lx為索MN的水平投影長度；ΔHx為索MN索力增量的水平分量；α 為材料線膨脹系數；Δtx為索MN溫度變化。在小垂度問題中，認為(dz0/dx)2與1 相比是微量。則式(7)可簡化為：

文獻[6]表明：當f/l≤0.1 時，這種近似引起的索力誤差在5%左右。令式(6)與式(8)相等，可得x方向承重索MN的變形協調方程，即：

同理，y方向穩定索PQ的變形協調方程為：

1.5 索網近似解

從理論上講，當給定索網初始狀態、邊界條件、荷載增量時，即已知z0、Hx0、Hy0、q0、Δq，由變形協調方程和平衡方程，可解出ΔHx、ΔHy和w三個未知量。但由于三個未知量是函數式，三個方程是積分或微分方程，直接求解并不容易。

實際應用中可求助于基于變分原理的各種近似解法，建立簡便的實用計算公式。常用的近似解法有配點法、力矩法、最小二乘法和伽遼金法等，因伽遼金法求解精度相對較高，采用伽遼金法對結構進行近似求解。具體步驟如下：

1) 根據索網的邊界條件和荷載分布情況，假設一個位移試函數w(x,y)，其中包含有限個待定常數，一般形式如下：

式中：fi(x,y) (i=0, 1, …,n)為選定的已知函數，f0在邊界處滿足給定的邊界條件；ci(i=1, 2, …,n)為n個待定參數。

2) 將假定的位移函數w代入變形協調方程，可得索水平張力增量ΔHx、ΔHy的函數，可寫為：

3) 將w、z0、ΔHx和ΔHy代入式(3)，得包含各待定參數ci的函數S。

4) 按平衡條件要求，函數S應恒等于0，但由于位移函數是事先假定的，不可能精確。此方法要求函數S與所選位移函數中的各已知函數fi正交，即：

可得n個方程，近似求解出n個待定參數ci。

2 橢圓鞍形正交索網解析解

平面投影為橢圓的鞍形索網由馬鞍面和橢圓柱面相切而成。馬鞍面與橢圓柱面相交而成的空間曲線即為結構的外環，其在xoy面上的投影即橢圓柱面在xoy面上的投影。索網結構的平面及側面投影、幾何參數、坐標體系如圖3 所示，側視圖中實線代表索初始態位置，虛線代表索荷載態位置。初始態下馬鞍面、橢圓柱面方程如下：

圖3 索網結構平面及側視圖Fig.3 Plan and side view of cable-net structure

式中：f1為中央承重索的跨中垂度；f2為中央穩定索的跨中拱度，設H為結構最高點、最低點之間的高差，則H=f1+f2；a、b分別為橢圓的兩個半軸長度。可見，f1、f2、a、b這4 個參數將決定索網結構的形態。

由于索網結構跨度大、自重輕，在初始態忽略索網自重，令q0=0，根據初始態平衡方程，可得：

即兩個方向索的初始預張力與曲面在相應方向的曲率成反比。若已知一個方向的預張力，則可確定滿足平衡方程的另一方向的預張力。

2.1 位移函數的選取

位移函數w需滿足以下條件：① 索網承受對稱荷載時，位移函數應同樣具有對稱的性質；② 索網的邊緣構件支承于下部結構之上，本身剛度很大，可認為索網邊界上的豎向位移為0。位移函數的選取將極大地影響位移解的精度，此前文獻均采用如下形式的位移函數，區別在于所選取的待定參數和已知函數的個數不同。

式中：wmn為需求解的待定參數。對于承受豎向均布荷載的情況，沈世釗等[6]取1 個待定參數w00，其位移函數精度差，不能準確描述索網的變形，但僅需求解一元方程。金問魯[7]取3 個待定參數w00、w20、w02，其位移函數可準確地描述索網變形，與有限元解十分接近；但需聯立求解方程組，很難寫出各待定參數的表達式。

因此，本文提出一個新的位移函數，該函數既滿足對稱性和邊界條件，又在僅取一個待定參數的條件下，有較高的精度，以求解索網豎向位移及內力 。位移函數形式如下：

式中：w0為待定參數，表示索網中點的豎向位移。

2.2 單參數三次方程解

采用伽遼金法求解位移函數中的待定參數，對馬鞍面方程式(15)、位移函數式(19)代入變形協調方程式(9)與式(10)，沿x、y方向積分，可得索力增量：

將馬鞍面方程式(15)、豎向位移函數式(19)及式(20)、式(21)代入平衡方程式(3)，經整理后得：

式中：A0、A1、A2、均為x、y的函數。

由式(13)，伽遼金變分方程為：

將函數S代入伽遼金變分方程，利用橢圓積分運算并整理后，可得用于確定待定參數w0的一元三次方程，即：

一元三次方程可用卡丹求根公式求解，但表達式過于復雜，實際中采用迭代法求解。

2.3 單參數一次方程解

索網在外荷載作用下的中點豎向位移與跨度的比值w0/a、w0/b通常是微量。在不考慮支座位移和溫度變化的條件下，為使一元三次方程便于求解與表達，忽略式(24)中w0/a及w0/b的二次項、三次項，將其簡化為一次方程：

得待定參數w0的線性表達式：

與三次方程相比，此表達式大大簡化，可快速計算出索網中點在荷載態下的豎向位移，并在此基礎上計算各節點豎向位移及各索索力增量。索網中任一節點的位移可由式(19)計算，任一根索的索力增量計算公式便簡化為：

3 索網形態優化

外環梁大多處壓彎狀態，重力荷載方向的彎矩由支承條件決定，而鞍面的形態對外環梁的面內彎矩有極大影響，外環梁越接近軸壓構件，其材料利用越充分，可相應地減小外環截面尺寸和變形。因此在結構平面投影尺寸a和b、高差H、承重索預拉力Hx0不變的條件下，以外環梁平面內彎矩最小為優化目標，以跨中承重索垂度f1為優化變量，以索力、位移及矢跨比為約束條件，采用解析方法求解結構最優形態。

3.1 外環荷載

為簡化計算做出以下假定：① 實際工程中馬鞍面高差與跨度相比較小，將空間曲梁簡化為平面構件；② 只考慮外環所受平面內彎矩；③ 外環橫截面處處相等；④ 由于外環剛度較大，計算時不考慮幾何非線性的影響；⑤ 將索對外環的集中力作用轉換為分布荷載。

馬少華[19]將索對外環的作用轉換為雙向均布荷載，根據平衡條件推導出外環任一點軸力、彎矩的表達式。但其未考慮索力的真實分布情況，所得彎矩表達式精度較差。以外側受拉為正，其外環上任一點彎矩M表達式為：

本文考慮索力的真實分布情況，由于對稱性，僅取1/4 外環進行分析，外環受力如圖4 所示。

圖4 1/4 外環受力簡圖Fig.4 Load condition of 1/4 outer-ring

承重索、穩定索初始索力Hx0、Hy0呈均勻分布，索力增量ΔHx、ΔHy呈曲線分布，索力增量按式(27)、式(28)計算。如何選取形態優化時采用的外荷載Δq，關系到形態優化的優劣，目前主要有兩種思路：① 取結構全壽命周期內發生概率最大的荷載組合作為外荷載；② 取各組合中最大、最小荷載增量的平均值作為外荷載[20]。

3.2 外環彎矩及其分布

考慮y方向平衡條件，可得：

式中：NA表示A點軸力。再對B點取距并將NA代入，可得B點彎矩值：

以外側受拉為正，結合式(29)和式(31)，假設外環上任一點彎矩M表達式為：

式(32)中的彎矩分布符合邊界條件，即θ=0°時，M=MA；θ=90°時，M=MB。MA為A點彎矩，可由等截面閉合環的彎矩沿周長積分為零的條件確定，即解得：

將MA代回式(32)，便得到外環上任一點的彎矩M表達式，式中 γ 的含義同式(29)。

1/4 外環彎矩分布曲線如圖5 所示，彎矩在A、B兩點達到最大，即結構最高點與最低點。

圖5 1/4 外環彎矩圖Fig.5 Bending moment of 1/4 outer-ring

3.3 索網最優形態求解

若式(34)中括號內結果等于0，外環上任一點彎矩便為0。以f1為變量，固定a、b、H、Hx0保持不變，改變索網空間形態。將式(17)、式(26)代入，可得關于f1的一元四次方程式(35)。求解得到f1后，便得到外環平面內彎矩處處為0 時的索網結構形態。

3.4 約束條件

索網結構設計時還需遵循一定的約束條件，也就是根據安全和設計要求對變量的限制，主要以索力、豎向位移、矢跨比為約束條件。

1) 索力。在任意荷載組合下，拉索不能松弛，索力應始終為正，且最大索力應滿足材料強度限值，以索水平力近似代替索力。索力約束條件為：

式中：αc為考慮材料安全儲備的安全系數；[σ]為材料容許應力。

2) 豎向位移。結構的豎向位移應滿足《索結構技術規程》[21]的要求，即最大撓度不宜大于跨度的1/250。位移約束條件為：

式中，[uz]為節點豎向位移限值。

3) 矢跨比。文獻[21]給出了承重索垂度、穩定索拱度與其跨度比值的適宜取值區間，還應根據建筑功能、造型具體考慮其限值。矢跨比約束條件為：

式中，[f/l]max、[f/l]min分別為矢跨比上限、下限值，應分別考慮承重索與穩定索。

綜上，用解析法對橢圓單層正交索網結構進行優化設計的完整步驟為：① 給定初始參數，包含橢圓長半軸長度a、短半軸長度b、高差H、承重索初始索力Hx0，確定索力、位移和矢跨比的限值；② 求解式(35)，計算使外環上任一點彎矩為0 時的承重索垂度f1，得到結構最優形態；③ 驗算所得形態在各種荷載組合下是否滿足約束條件；④ 若滿足，則該形態為結構的最優形態；若不滿足，可調整初始索力Hx0或索截面面積Ax、Ay，返回步驟② 再次求解。

4 有限元驗證

本節將前文所推導的索網位移和索力計算公式、外環彎矩計算公式、形態優化結果，與有限元分析結果進行對比，以驗證解析計算方法的準確性。為便于比較，取與文獻[7]P212 算例相同的幾何參數、荷載條件，具體參數如表1 所示。

表1 橢圓單層正交鞍形索平網結構參數Table 1 Parameters of elliptic saddle single-layer orthogonal cable-net structure

4.1 有限元找形、荷載分析方法

采用ANSYS APDL 參數化設計語言編寫找形及荷載分析程序。由于索中預拉力遠大于自重引起的張力，索可用只受拉的桿單元模擬，采用Link10 單元模擬懸索，初始預拉力通過初應變施加；采用Surf154 單元模擬膜面施加均布荷載；采用Beam189 單元模擬外環梁。程序分為以下3 部分：

1) 找形分析

對索網結構進行荷載分析之前，首先應通過找形分析，建立初始態時索網結構模型。找形分析將初始預張力分布、邊界條件作為已知，尋求對應的平衡形態。本文采用一種改進的綜合找形法[22]進行初始形態分析。該方法在支座位移法找形后，約束邊界點，采用節點平衡法迭代求解，此時形態變化較大，控制位移達到要求精度；再恢復真實彈模、初應變，再度采用節點平衡法迭代求解，此時索力變化較大，控制索力達到要求精度。共3 次找形過程，具體分析流程見圖6(a)。

圖6 有限元分析流程圖Fig.6 Finite element analysis process

對案例進行分析，找形后各節點z坐標值與理論雙曲拋物面最大相差3.5 mm，索力值與預設初始預拉力最大相差0.47 kN。該方法找形結果準確，且找形過程共經歷7 次計算，收斂快、計算效率高。

2) 索網荷載分析

在找形后的模型上，借助表面效應單元施加豎向均布荷載，并將其轉換為節點荷載施加于索網自由節點上，具體分析流程見圖6(b)。

3) 外環荷載分析

在荷載態模型基礎上，記錄邊界節點上的節點力，刪除索網并建立外環，將所記錄的節點力施加給外環梁，具體流程見圖6(c)。

4.2 索網位移、索力對比

將有限元所得位移、索力結果與上文及文獻[6 - 7]的解析計算結果進行對比，近似求解方法均采用伽遼金法。在所有節點中，中點附近位移最大，其值決定了結構是否能滿足撓度限值，因此僅對比中點豎向位移結果，如表2 所示。在所有拉索中，中央承重索、穩定索的索力變化量最大，其索力值決定了所有索是否松弛與是否滿足強度要求，因此僅對比中央承重索、穩定索索力增量，如表3 所示。

表2 不同求解方式中點豎向位移解對比Table 2 Comparison of vertical displacement of midpoints in different solutions

表3 不同求解方式中央索索力增量對比Table 3 Comparison of central cable force increment in different solutions

應用本文位移函數的單參數一次方程解所得位移、索力結果，與本文單參數三次解相比，結果相差很小，說明將三次方程簡化為一次方程是可行的。與有限元相比，其誤差在可接受范圍內；與文獻[6]單參數解相比，位移解結果精度大大提高；與文獻[7]三參數解相比，雖然精度稍有差距，但所需運算量大大減小且表達式簡便，便于后續求解結構最優形態。因此可使用本文所提出的計算方法快速估算索網結構在各種荷載作用下的位移、內力。

4.3 外環彎矩對比

1) 不同位置彎矩值對比

取表1 算例，對1/4 外環計算結果進行對比，每15°取一測點，分別采用有限元法、式(29)和式(34)計算，外環彎矩值隨角度θ 的變化如圖7所示。相比于式(29)的計算方法，本文的彎矩計算公式所得結果與有限元解更接近，在結構最高點(θ=0°)、最低點(θ=90°)誤差分別為1.1%、3.4%。

圖7 不同角度θ 時外環彎矩值對比Fig.7 Comparison of outer-ring moment at different angles θ

2) 不同形態彎矩值對比

固定H=5 m，使f1在2.8 m～4.0 m 區間內變化，f2、Hy0隨之改變，其余參數保持不變，分別采用有限元法、式(34)計算外環A、B點彎矩值。彎矩值隨f1的變化如圖8 所示，解析解與有限元結果基本吻合。A點彎矩隨著f1的增加，從內側受拉轉變為外側受拉，B點則相反；A、B點彎矩均接近0 時，外環越接近軸壓狀態。

圖8 不同垂度f1 時外環彎矩值對比Fig.8 Comparison of outer-ring moment at different sags f1

綜上，通過與有限元結果對比，驗證了本文所提出的外環彎矩計算公式的準確性。

4.4 最優形態對比

采用表1 案例，固定H=5 m，改變f1、f2以求解結構最優形態。通過式(35)解得f1=3.52 m，則f2=1.48 m，即為解析方法求得的索網最優形態。

實際情況中并不存在一個使外環梁彎矩處處為0 的形態，在有限元模型中以外環彎矩絕對值的最大值最小為目標，以f1為變量，尋求最優形態。求得f1=3.53 m 時，為有限元法得到的結構最優形態。

分別取表1 案例、解析法與有限元法所得最優形態的幾何參數，采用有限元法進行荷載分析，外環梁上最大彎矩、軸力值對比如表4 所示。

表4 不同形態下最大彎矩、軸力值對比Table 4 Comparison of maximum bending moment and axial force under different shapes

與表1 案例相比，解析法與有限元法所得最優形態下，外環梁上最大彎矩均大幅度減小并接近于0，最大軸力小幅度增加。外環中應力大大減小，可減小外環截面面積，使得材料強度充分利用。解析法與有限元法形態優化結果僅相差0.01 m，且外環上彎矩、軸力最大值相差很小，解析法所得形態優化結果具有相當高的精度。

5 結論

本文推導了橢圓鞍形正交單層索網連續化計算方法和外環梁上任一點彎矩計算公式。隨后，提出一種快速求解此類結構最優形態的方法，該方法以外環面內彎矩最小為優化目標，以中央承重索垂度為優化變量，以索力、豎向位移和矢跨比為約束條件，利用外環彎矩計算公式，得到外環彎矩最小時的結構形態。本文得到的主要結論如下：

(1) 對橢圓鞍形單層正交索網結構的連續化計算方法中所假設的位移函數進行了改進，利用變形協調方程、平衡微分方程結合伽遼金法推導出結構荷載態時的位移、內力。改進后的位移函數既具有較高的精度又使得表達式簡便，并采用有限元驗證了計算公式的準確性，可利用本文公式快速計算出給定結構的位移、內力。

(2) 利用所得索力，對外環施加更符合真實情況的荷載作用，改進了外環上任一點彎矩的計算公式。通過與有限元對比分析，本文外環彎矩計算公式具有較高的精度。

(3) 提出一種快速求解索網結構最優形態的方法，該方法以外環梁面內彎矩最小為優化目標，以中央承重索垂度為優化變量，以索力、豎向位移和矢跨比為約束條件。利用外環彎矩計算公式，得到以優化變量為未知數的一元方程，求解方程便快速得到結構的最優形態。解析方法所得結構最優形態與有限元結果基本一致，優化結果具有相當高的精度。在結構初步設計階段，利用本文解析計算方法，可快速、準確的獲得結構最優形態，避免了有限元法編程、反復試算的繁雜過程。