共軛梯度法在膜結構分析中的應用

高振宇

1 概述

力密度法應用于膜結構分析時建立的方程組如下:

其中，C為結構的拓撲矩陣;Q為力密度矩陣;Px，Py和Pz均為節點荷載，具體推導詳見文獻[1][2]。隨著結構網格劃分的增加，總結點數變大，力密度法方程組階數升高，膜結構分析的最終問題就歸結為大型線性方程組的求解上。C矩陣本身是含有大量零元素的稀疏矩陣，并且本文通過推導證明，該方程組具有對稱正定的特點，所以，本文提出可利用其方程組的結構特性進行算法優化，以提高計算程序的效率。

一般來說，解線性方程組 Ax=b有兩種數值分析方法:1)直接法。但直接法涉及到系數矩陣的分解，當 A為大型稀疏矩陣時，因為所要求的分解因子是稠密的，直接法的效率不高甚至難以實現。2)迭代法。這類方法是生成近似解序列{x(k)}，矩陣 A只在矩陣向量乘法時才用到，而稀疏矩陣的乘法已相當成熟。數值方法的一個基本原則是:求解任一問題都應利用其方程組的結構特性，下面將介紹一種迭代解法——共軛梯度法，其優點正好體現于求解稀疏、對稱及正定的方程組的效率上。

2 共軛梯度法在力密度方程組求解中的應用

2.1 力密度方程組的對稱正定性

根據式(1)，設任意n為向量p，且力密度矩陣 Q為對角矩陣，其中對角元素 q(i，i)＞0，(i=1，n)，則:

所以方程組系數矩陣(CTQC)為對稱正定矩陣。

2.2 共軛梯度法

表1 共軛梯度法與Gauss消去直接法找形結果對比表

表2 共軛梯度法與Gauss消去直接法計算耗時對比表

3 算例分析

根據以上算法，本文編制了相應的程序，并進行了算例分析。

本算例為馬鞍形膜結構，采用矩形網格劃分，劃分間距為0.5 m，膜面預應力為1.0 kN/m，膜網內部力密度值為1.0 kN/m，邊索力密度值為14.0 kN/m。結點總數為552，其中4個為固定結點，坐標分別為 99 000 001(12，16，4)，99 000 002(3，10，0)，99 000 003(12，1，4)，99 000 004(20，9，0)。單元總數 1 000，其中“T”單元總數為128。

本文編制的程序計算迭代166次，滿足終止準則‖r‖2＜1×e-6。圖1，圖2分別給出結構網格劃分、幾何零狀態時與找形后的平衡形狀的對比，表1，表2分別給出了共軛梯度法和Gauss消去直接法的結點坐標計算結果及其計算耗時的比較。算例分析表明，本文算法精度高、收斂快，是求解力密度方程組的可靠、高效的數值分析方法，為力密度法應用于大型膜結構分析奠定了堅實的基礎。

[1] H.J.Schek.The force density methodfor form finding and computation of general network[J].Computer Methods in Applied Mechanics and Engineer，1974(3):115-134.

[2] 王 勇，魏德敏，高振宇.張拉膜結構力密度法混合找形分析[J].計算力學學報，2005(10):629-632.

[3] 鄭 超，張建海.預處理共軛梯度法在巖土工程有限元中的應用[J].巖石力學與工程學報，2007(7):2820-2825.

[4] G.H.戈盧布，C.F.范洛恩.矩陣計算[M].袁亞湘，譯.北京:科學出版社，2001.

[5] 李慶揚，關 治，白峰杉.數值計算原理[M].第2版.北京:清華大學出版社，2001.