劉海峰,李正光

(1. 西安交通大學 數(shù)學與統(tǒng)計學院,西安 710049; 2. 吉林大學 數(shù)學學院,長春 130012)

結(jié)構靜力重分析預處理共軛梯度法的一種有效改進

劉海峰1,李正光2

(1. 西安交通大學 數(shù)學與統(tǒng)計學院,西安 710049; 2. 吉林大學 數(shù)學學院,長春 130012)

考慮預處理共軛梯度法在結(jié)構靜力重分析中的改進. 利用增量剛度矩陣相對結(jié)構修改后剛度矩陣稀疏的特點,修改處理共軛梯度法的實施過程,以達到降低計算量、 節(jié)省計算時間的目的. 數(shù)值算例驗證了該方法的有效性.

預處理共軛梯度法； 靜力重分析； 稀疏性； Cholesky分解

在飛行器設計、 橋梁隧道和能源與動力等領域,有許多大型復雜的工程結(jié)構需要進行設計. 結(jié)構設計通常需要進行多次修改,每次修改都需要對結(jié)構進行分析,計算量較大. 為改善這種情況,人們提出了結(jié)構靜力重分析問題. 結(jié)構靜力重分析的目的是盡可能多地利用初始分析的相關信息,設計高效算法計算修改后結(jié)構的響應,以降低計算量,節(jié)省計算時間,加快結(jié)構設計進程[1].

目前,重分析方法主要分為兩類： 逼近重分析法和精確重分析法. 逼近重分析法又分為迭代逼近法、 組合逼近法、 多點逼近法和單點逼近法[2]. 其中最有效的逼近重分析法是迭代逼近法中的預處理共軛梯度法,該算法用給定的預處理矩陣作用結(jié)構修改后的平衡方程,使系數(shù)矩陣的條件數(shù)大幅度減小,從而達到加快收斂的目的. 此外,預處理共軛梯度法還可以處理結(jié)構節(jié)點數(shù)目發(fā)生變化的重分析問題[3]. 文獻[4]給出了預處理共軛梯度法在重分析問題中的最新研究進展.

精確重分析法又稱為直接重分析法,在不考慮計算機舍入誤差的情況下,可計算出結(jié)構修改后位移向量的精確解. 這類重分析方法一般都建立在Sherman-Morrison-Woodbury(SMW)低秩校正公式基礎上. 文獻[5]將直接重分析法應用到裂縫增長計算中,取得了較好的效果.

本文考慮預處理共軛梯度法在結(jié)構靜力重分析問題中的實施問題,利用增量剛度矩陣比結(jié)構修改后的剛度矩陣更稀疏的特點,修改了預處理共軛梯度法的執(zhí)行過程,從而降低了計算量,節(jié)省了計算時間.

1 問題的數(shù)學描述

設初始結(jié)構自由度數(shù)目為m,給定其剛度矩陣K0∈m×m和荷載向量R0∈m,則位移向量x0可使用如下平衡方程：

計算,其中K0是一個對稱正定矩陣. 根據(jù)初始分析,它的Cholesky分解

假設結(jié)構修改后,其剛度矩陣變?yōu)?/p>

其中ΔK稱為增量剛度矩陣. 則相應的平衡方程變?yōu)?/p>

其中:R∈m表示修改后結(jié)構的荷載; 剛度矩陣K∈m×m也是一個對稱正定矩陣. 靜力重分析的目的是盡量多地利用結(jié)構修改前的初始信息計算修改后結(jié)構的位移向量x,以降低計算量.

2 靜力重分析中的預處理共軛梯度法

在逼近重分析法中,預處理共軛梯度法是最有效的方法. 在結(jié)構修改程度較小的情況下,選取初始結(jié)構的剛度矩陣做預處理矩陣,使用少數(shù)幾步迭代即可得到精度較高的近似解.

算法1[6-8]

假設初始猜測x=0

k=0

r=R

δ0=rTr

t=‖R‖2

SolveK0z=rforz

ρk+1=rTz

k=k+1

ifk=1

p=z

else

p=z+βp

end

w=Kp

x=x+αp

r=r-αw

δk=rTr

end

其中ε是事先給定的相對誤差精度. 算法1的計算量主要體現(xiàn)在兩個步驟： 解方程K0z=r和計算矩陣向量乘積w=Kp. 由于K0的Cholesky分解已知,使用前回代法和后回代法即可對K0z=r進行求解. 下面研究w=Kp的計算問題.

在第k步迭代中,分別記w,p,z,r為wk,pk,zk,rk,則由式(3)在第k步迭代中,有

因此,w=Kp可寫成另一種形式：

w=r+ΔKz+βw.

由于在靜力重分析問題中,增量剛度矩陣ΔK的非零元素個數(shù)通常遠少于矩陣K的非零元素個數(shù),所以ΔKz的計算量遠小于w=Kp的計算量. 因此可將算法1的執(zhí)行過程進行修改,得到如下算法.

算法2

假設初始猜測x=0

k=0

r=R

δ0=rTr

t=‖R‖2

SolveK0z=rforz

ρk+1=rTz

k=k+1

ifk=1

p=z

else

p=z+βp

end

w=r+ΔKz+βw

x=x+αp

r=r-αw

δk=rTr

end

算法1稱為標準的預處理共軛梯度法,算法2稱為改進的預處理共軛梯度法.

3 數(shù)值算例

下面用數(shù)值算例驗證算法2的有效性. 電腦配置為Pentium 4,3.2 GHz CPU,2 Gb RAM. 使用軟件為MATLAB 7.2.0.232(R2006a).

考慮如圖1所示的橋體結(jié)構,橋長660 m,寬24 m,每個橋墩高24 m,6個橋塔均高出橋面36 m. 該橋可離散成一個具有4 420個節(jié)點和4 626個單元的有限單元模型. 除18個約束節(jié)點外,每個節(jié)點有6個自由度,因此整個結(jié)構的自由度為26 412. 結(jié)構單元有3種類型： 索單元(180個),梁單元(486個),板單元(3 960個). 索和梁單元使用材料的彈性模量均為E1=2×1011Pa,梁單元使用的材料Poisson比為υ1=0.3； 板單元使用材料的彈性模量為E2=3×1010Pa,Poisson比為υ2=0.2. 索單元的橫截面積為0.031 416 m2,梁單元的橫截面為1 m×1 m,板厚0.2 m. 橋面每個節(jié)點受垂直向下的力,大小為P=1×104N. 為加固此橋體結(jié)構,將索單元的橫截面積增加到0.045 239 m2,梁單元的橫截面變?yōu)?.2 m×1.2 m. 分別使用算法1和算法2對修改后的結(jié)構進行計算,相對誤差精度ε=10-6. 由算法1和算法2計算橋結(jié)構的計算時間分別為1.095 313 s和0.814 062 5 s. 本文所提出方法的計算時間約為標準預處理共軛梯度法的75%.

圖1 橋結(jié)構Fig.1 Structure of the bridge

綜上所述,本文對靜力重分析中預處理共軛梯度法的執(zhí)行問題進行了研究,利用靜力重分析問題中增量剛度矩陣ΔK比結(jié)構修改后剛度矩陣K稀疏的特點,將預處理共軛梯度法的執(zhí)行過程進行修改,節(jié)約了計算時間,保留了預處理共軛梯度法的全部優(yōu)點. 同時,由于本文提出的算法不需要存儲K,因此所需的存儲空間也小于標準預處理共軛梯度法.

[1]Kassim A M A,Topping B H V. Static Reanalysis: A Review [J]. Journal of Structural Engineering,1987,113(5): 1029-1045.

[2]Kirsch U. Reanalysis of Structures [M]. Dordrecht: Springer,2008: 121-158.

[3]Wu B S,Lim C W,Li Z G. A Finite Element Algorithm for Reanalysis of Structures with Added Degrees of Freedom [J]. Finite Elements in Analysis and Design,2004,40(13/14): 1791-1801.

[4]Voormeeren S,Rixen D. Updating Component Reduction Bases of Static and Vibration Modes Using Preconditioned Iterative Techniques [J]. Computer Methods in Applied Mechanics and Engineering,2013,253(1): 39-59.

[5]Pais M J,Yeralan S N,Davis T A,et al. An Exact Reanalysis Algorithm Using Incremental Cholesky Factorization and Its Application to Crack Growth Modeling [J]. International Journal for Numerical Methods in Engineering,2012,91(12): 1358-1364.

[6]Golub G H,Loan C F,van. Matrix Computations [M]. 3rd ed. Baltimore: The Johns Hopkins University Press,1996: 532-544.

[7]李正光,吳柏生. 關于結(jié)構剛度修改的一個新算法 [J]. 吉林大學學報： 理學版,2003,41(1): 36-39. (LI Zhengguang,WU Baisheng. A New Algorithm for Modifications of Structural Stiffness [J]. Journal of Jilin University: Science Edition,2003,41(1): 36-39.)

[8]Kirsch U,Kocvara M,Zowe J. Accurate Reanalysis of Structures by a Preconditioned Conjugate Gradient Method [J]. International Journal for Numerical Methods in Engineering,2002,55(2): 233-251.

(責任編輯： 趙立芹)

EfficientImprovementofPreconditionedConjugateGradientMethodforStructuralStaticReanalysis

LIU Haifeng1,LI Zhengguang2
(1.SchoolofMathematicsandStatistics,Xi’anJiaotongUniversity,Xi’an710049,China;
2.CollegeofMathematics,JilinUniversity,Changchun130012,China)

This paper deals with the implementation of preconditioned conjugate gradient method (PCG) for structural static reanalysis. The implementation of PCG is revised by exploring sparsity of incremental stiffness matrices so that the computational cost and time can be reduced. Numerical example is given to show the effectiveness of the method.

PCG method; static reanalysis; sparsity; Cholesky factorization

2013-12-27.

劉海峰(1982—),男,漢族,博士,講師,從事計算數(shù)學的研究,E-mail： lhflkcc@126.com. 通信作者： 李正光(1975—),男,漢族,博士,教授,從事計算力學的研究,E-mail： lizg@jlu.edu.cn.

國家自然科學基金(批準號： 51005096).

O342

A

1671-5489(2014)05-0971-04