含螺栓連接結(jié)構的非線性靜力學拓撲優(yōu)化

張曉頔， 郝志明， 李一堃， 劉遠東

(1.中國工程物理研究院總體工程研究所， 四川 綿陽 621900； 2.成都大學高等研究院， 成都 610106)

連接結(jié)構在航空航天、軍事裝備和核電等領域中均有廣泛應用[1-2]，大型裝備如航天器、軍用車輛、反應堆壓力容器等是由不同的零部件通過螺栓預緊連接、鉚釘連接、楔環(huán)連接等多種連接方式連接而成。在結(jié)構系統(tǒng)的拓撲優(yōu)化設計中，認為連接結(jié)構的接觸表面不發(fā)生相對運動，將其簡化為線性模型分析，計算效率更高。實際上連接結(jié)構的剛度和阻尼具有不連續(xù)性、遲滯性等非線性特性，分布于結(jié)構系統(tǒng)不同位置的剛度不同，對整體結(jié)構的力學特性的影響也不同[3]。由于連接接觸界面在外載荷作用下發(fā)生微觀滑移現(xiàn)象甚至宏觀滑移現(xiàn)象，引起能量耗散并產(chǎn)生干摩擦阻尼，占結(jié)構整體阻尼的90%；同時連接接觸界面的滑移還會引起結(jié)構局部剛度非線性變化。因此，在對復雜工程結(jié)構開展拓撲優(yōu)化研究時，考慮其連接部分的非線性具有重要的意義。

結(jié)構拓撲優(yōu)化研究經(jīng)歷了從離散體結(jié)構到連續(xù)體結(jié)構的發(fā)展過程。從早期由Michell的桁架理論和Dorn等人提出的基結(jié)構法建立起離散結(jié)構的拓撲優(yōu)化方法，發(fā)展到由Bends?e和Kikuchi提出的結(jié)構拓撲優(yōu)化設計均勻化方法，標志著連續(xù)體結(jié)構拓撲優(yōu)化設計的誕生。Xie[4]又提出了一種連續(xù)體拓撲優(yōu)化方法——漸進結(jié)構優(yōu)化方法(Evolutionary Structural Optimization，ESO)，與均勻化方法相比，它不產(chǎn)生0-1之間的中間密度，拓撲形式清晰，也避免了多變量的數(shù)學規(guī)劃求解，這種方法已成功應用于大量工程優(yōu)化問題的求解。目前關于線彈性材料和小變形的結(jié)構拓撲優(yōu)化研究較為充分，但考慮結(jié)構非線性行為的研究仍相對較少[5]。Ryu等[6]率先進行了非線性響應的設計靈敏度分析，開辟了非線性結(jié)構優(yōu)化研究的先例。Buhl等[7]研究了幾何非線性結(jié)構的拓撲優(yōu)化，指出非線性和線性情況下優(yōu)化結(jié)果的區(qū)別。Yuge等[8]將材料非線性引入結(jié)構拓撲優(yōu)化，改進彈塑性均勻化方法，提高了計算效率。Jung等[9]在優(yōu)化中同時考慮了幾何非線性和材料非線性，進行了數(shù)值驗證。考慮界面非線性的結(jié)構優(yōu)化研究文獻很少，僅在尺寸優(yōu)化和形狀優(yōu)化方面有少量工作[10-11]，尚未見界面非線性問題的拓撲優(yōu)化研究報道。

開展含螺栓連接結(jié)構的非線性拓撲優(yōu)化設計，需模擬連接界面的接觸行為。傳統(tǒng)摩擦模型[12]如庫倫摩擦模型在有限元分析中主要采用接觸算法，計算較為耗時且不易收斂，難以應用于拓撲優(yōu)化分析。因此，尋求一種更加合理的模型描述連接界面的接觸行為十分必要。連接接觸本構模型[13]是一類能描述摩擦接觸非線性遲滯行為的模型，如Iwan模型、Bouc-Wen模型、Valanis模型和剪切層模型等。其中Iwan模型[14-15]由一系列理想彈塑性元件組成，最初用以描述金屬材料的彈塑性力學行為。Segalman[16]最早采用并聯(lián)-串聯(lián)形式的Iwan模型描述連接接觸非線性力學行為，提出了一個四參數(shù)Iwan模型，能夠準確反映微觀滑移能量耗散的冪次關系，但不能描述發(fā)生宏觀滑移后的殘余剛度現(xiàn)象。在四參數(shù)Iwan模型的基礎上，Li等[17-18]提出了一個含截斷冪律分布和雙脈沖的非均勻密度函數(shù)的六參數(shù)Iwan模型，可以同時描述微觀滑移階段能量耗散冪次關系和宏觀滑移階段殘余剛度現(xiàn)象。此外，該模型引入有限元分析，不需要采用接觸算法，計算易于收斂[19-20]。

本文旨在建立一種含界面強非線性的靜力學拓撲優(yōu)化方法。采用六參數(shù)Iwan模型描述含螺栓連接的非線性力學行為，ESO方法結(jié)合結(jié)構非線性分析，用ANSYS軟件的參數(shù)化設計語言(APDL語言)編寫命令流文件，形成一套優(yōu)化算法，應用于含螺栓連接圓筒結(jié)構的拓撲優(yōu)化設計。

1 螺栓連接的非線性力學模型

Iwan模型由若干Jenkins單元(彈簧-滑塊)構成。圖1所示為并聯(lián)-串聯(lián)型Iwan模型，其中ki為彈簧單元的剛度，fi為滑塊單元的屈服力，F(xiàn)和u為外力和位移。令每個Jenkins單元的剛度均為k，位移載荷為x，單元屈服力f以密度函數(shù)r(f) 表示，可以得到Iwan模型表達式如下：

(1)

將式(1)對x求二階導數(shù)后整理得

(2)

再令φ=f/k,ρ(φ)=k2r(f)，則式(1)可改寫為

(3)

(4)

式中ρ(φ)為密度函數(shù)。Li等[17-18]提出一個含截斷冪律分布和雙脈沖的非均勻密度函數(shù)，如圖2所示，其表達式為

(5)

將式(5)代入式(3)中積分可以得到六參數(shù)Iwan模型的力-位移曲線關系，其示意圖見圖3。O-A段力-位移關系呈線性；A-B段為微觀滑移階段，力-位移關系呈非線性；從B點開始模型進入宏觀滑移階段，力-位移關系再次呈線性。

將六參數(shù)Iwan模型應用于數(shù)值計算時，將其離散為n+2個并聯(lián)的Jenkins單元，如圖4所示[20]。前n個Jenkins單元J1～Jn描述微觀滑移階段連接結(jié)構剛度的變化；第n+1個Jenkins單元Jn+1描述宏觀滑移時刻剛度的變化量，其剛度為K2；第n+2個Jenkins單元Jn+2為描述宏觀滑移后接觸界面殘余剛度現(xiàn)象的彈簧單元，其剛度為K∞。

2 考慮接觸非線性的漸進結(jié)構優(yōu)化

一個理想的結(jié)構，其每一部分的應力應該接近于相同的安全水平[21]。而受力分析表明結(jié)構中的應力分布并不均勻，可以認為應力水平較高的區(qū)域是結(jié)構破壞的主要區(qū)域，而在應力水平較低的區(qū)域材料處于沒有充分利用的狀態(tài)。以單元的等效應力作為刪除判據(jù)，將低應力的材料逐漸刪除，更新結(jié)構設計，使得優(yōu)化后的結(jié)構應力水平更加均勻，材料構成趨于合理，同時也提高了材料的利用率。對于各向同性材料可采用von Mises應力作為單元刪除判據(jù)。

(6)

式中Sij為應力偏張量，且i=j。

(7)

則認為該單元處于低應力狀態(tài)，是無效的或低效的單元，對結(jié)構整體性能貢獻最小，應從結(jié)構中刪除。

保持當前刪除率不變，重復刪除單元，直至再無單元滿足式(7)，即優(yōu)化結(jié)果對于當前刪除率RR已達到穩(wěn)定狀態(tài)。此時引入一個進化率ER，將刪除率修改為

RR=RR+ER

(8)

以新的刪除率重復刪除和更新刪除率過程，直到結(jié)構滿足給定的約束條件為止。

含螺栓連接結(jié)構的拓撲優(yōu)化的約束條件應包括兩個方面：1)結(jié)構的強度要求，即結(jié)構的應力不超過許用應力；2)結(jié)構的連接剛度要求，即一般不允許螺栓連接界面產(chǎn)生宏觀滑移，造成螺栓受剪切載荷作用。以結(jié)構的質(zhì)量最小為目標函數(shù)，結(jié)構的最大等效應力滿足要求、連接不產(chǎn)生宏觀滑移為約束條件，以設計區(qū)域材料單元的存在狀態(tài)為設計變量，建立如下拓撲優(yōu)化數(shù)學模型：

(9)

基于應力水平的含螺栓連接結(jié)構的非線性漸進結(jié)構優(yōu)化方法的步驟如下：

1) 建立含螺栓連接的組合結(jié)構非線性有限元模型。利用六參數(shù)Iwan模型描述螺栓連接界面的非線性特性；定義載荷及邊界約束條件；

2) 確定優(yōu)化準則。本文采用應力準則，如式(7)；

3) 定義初始刪除率RR和進化率ER；明確結(jié)構的許用應力和螺栓連接處相對位移(滑移量)允許值；

4) 采用增量法分步施加載荷，對當前離散結(jié)構進行非線性靜力學求解；依據(jù)當前刪除率RR，將滿足式(7)的低應力單元組成集合Ωσ，刪除集合中的單元；

5) 保持當前刪除率不變，重復步驟4)，直至再無單元滿足優(yōu)化準則，即優(yōu)化結(jié)果對前刪除率RR已達到穩(wěn)定狀態(tài)。此時引入進化率ER，按式(8)對刪除率進行修改；

6) 以新的刪除率重復步驟4)～5)，直到單元最大等效應力達到許用應力或螺栓連接處的相對位移達到允許值，優(yōu)化完成。

3 算例

3.1 模型的建立

圖5(a)所示為一含螺栓連接的圓筒組合結(jié)構，由圓筒和底座兩個部分組成。底座上半部分為一內(nèi)徑120 mm、外徑370 mm、厚度10 mm的圓盤，其上表面與圓筒連為一體，下半部分為直徑370 mm、厚度10 mm的圓板。三對單螺栓連接組件以120°間隔均勻排布在直徑為300 mm的圓周上，將底座上下部分連接起來。螺栓預緊力為7 117 N(圖5中未畫出連接螺栓)；圓筒高為335 mm，內(nèi)徑為230mm，外徑為234 mm，厚為2 mm。結(jié)構材料為鋼材，彈性模量=200 GPa，泊松比ν=0.3，密度ρ=7 850 kg/m3，許用應力σ*=105 MPa。

有限元建模采用SOLID185單元，以六面體網(wǎng)格將結(jié)構劃分為31 412個單元，如圖5(b)所示。

文獻[23]對預緊力為7 117 N時的螺栓連接實驗結(jié)果進行了參數(shù)辨識，具體值見表1。將表1的參數(shù)代入六參數(shù)Iwan模型，獲得如圖3所示的單個螺栓連接的力-相對位移關系(骨線方程)。圖3B點為宏觀滑移的起始點，對應的力為4 430 N，相對位移φ2=1.02×10-5m。若取不發(fā)生宏觀滑移的安全系數(shù)[n]=1.2，則允許的最大相對滑移量Smax≤0.85×10-5m。

表1 預緊力7 117 N的參數(shù)辨識結(jié)果

3.2 拓撲優(yōu)化

在ANSYS有限元分析中，每對螺栓連接組件之間均以六參數(shù)Iwan模型描述，將每對連接離散為13個COMBIN40單元和1個COMBIN14單元(圖4)，各單元的單元剛度和屈服力的選取參見文獻[23]，APDL命令流如下。

/prep7

!定義單元類型!

et,1,40

et,2,14

!定義單元的剛度和屈服力!

R,1,1.31E+07,,,,2,,

R,2,1.31E+07,,,,8,,

R,3,1.31E+07,,,,15,,

R,4,1.31E+07,,,,24,,

R,5,1.31E+07,,,,35,,

R,6,1.31e+07,,,,46,,

R,7,1.31e+07,,,,57,,

R,8,1.31e+07,,,,70,,

R,9,1.31e+07,,,,83,,

R,10,1.31e+07,,,,97,,

R,11,1.31e+07,,,,111,,

R,12,1.31e+07,,,,126,,

R,13,3.48e8,,,,3549.6,,

R,14,2.065e7,,,,,,

!將1-13號單元設置為類型1!

*do,rin,1,13

type,1,

real,rin

e,m,n

*enddo

!將14號單元設置為類型2!

type,2,

real,14,

e,m,n

(說明：命令流中m和n分別為每對螺栓連接的兩個接觸面中心節(jié)點的編號。)

約束底座下端面，對筒體上端面施加7.5 MPa的均勻壓力，并選定初始刪除率RR=0.2%，進化率ER=0.1%。利用編寫的命令流輸入文件，完成非線性求解、優(yōu)化循環(huán)及后處理等過程[22]。

優(yōu)化后的圓筒結(jié)構如圖6所示，最大等效應力為105 MPa(對應刪除率為60%)。最終優(yōu)化構型中肋條的上部沿圓周均勻分布，承受加載在筒體上端面的均布載荷；下部與螺栓連接結(jié)構的排布相對應，呈120°分別集中，、向下方的連接部分傳力。很明顯，優(yōu)化后得到拓撲構型很好地反映了結(jié)構的受力和傳力特征。

表2中列出了隨刪除率的增大，筒體的單元最大von Mises應力的變化。從表2可以看出，筒體的最大等效應力隨刪除率增加而逐漸增大，當刪除率達60%時，筒體最大等效應力達到許用應力，此時連接單元的相對滑移量Smax=8.12×10-6m，小于允許的最大相對滑移量8.5×10-6m，即連接處沒有產(chǎn)生宏觀滑移。

表2 筒體應力隨刪除率的變化

4 結(jié)論

本文采用六參數(shù)Iwan模型描述螺栓連接的非線性力學行為，再基于漸進結(jié)構優(yōu)化方法和非線性有限元分析，建立含螺栓連接圓筒結(jié)構的非線性靜力學拓撲優(yōu)化方法。利用ANSYS軟件中的APDL語言編寫命令流輸入文件，實現(xiàn)優(yōu)化流程。通過含螺栓連接圓筒組合結(jié)構的優(yōu)化算例，驗證了含界面強非線性的靜力學拓撲優(yōu)化方法可以應用于含螺栓連接復雜工程結(jié)構的靜力學拓撲優(yōu)化設計。