設計相關動載荷作用下的水下耐壓結構拓撲優化

詹 可，蔣垣騰，趙 敏

(上海交通大學 船舶海洋與建筑工程學院 海洋工程國家重點實驗室，上海 200240)

潛水器是進行深海勘探和海洋科學研究的重要裝備，也是一個國家科學技術水平的重要體現。耐壓殼作為潛水器的關鍵結構，其設計對潛水器的安全性和空間利用率等性能具有重要影響，在保障設備工作和人員安全方面起著重要的作用。目前學者們往往按照傳統的金屬材料耐壓結構設計方法開展水下耐壓殼結構設計，即將靜水壓力作為環境載荷作用于耐壓結構上，然后開展安全性設計。但是，水下耐壓結構也會面臨一種特別的沖擊載荷——內爆。當水下耐壓結構不能承受外部水壓而被壓潰塌陷時，流場靜水壓力轉化為流體動能，水流壓縮結構至最小限度時，會發生水錘型的沖擊，水流動能轉化為沖擊波壓力對周圍結構造成破壞，這種沖擊載荷就是內爆載荷。單體結構內爆誘發多體結構內爆的現象就是殉爆。美國于2009年研制成功的11 000米級“海神號”無人潛水器搭載了大量陶瓷耐壓結構，其結構如圖1所示[1]，該潛水器在探索位于新西蘭的世界第二深海溝克馬德克海溝時下潛至9 990 m處由于單體陶瓷耐壓結構內爆引發多體殉爆，導致整個潛水器損毀，所有設備幾乎成為碎片。

圖1 “海神號”上陶瓷耐壓結構Fig. 1 Ceramic pressure hulls on hybrid remotely operated vehicle Nereus

美國海軍[2]2004年發布的水下無人潛水器(UUV)計劃要求潛水器的設計必須嚴格考慮耐壓支承結構內爆所產生的影響。Turner[3]使用中空玻璃耐壓結構在6.996 MPa的壓力筒中進行水下內爆試驗，得到了沖擊波壓力曲線；陳鋒華和趙敏[4]對該耐壓結構進行內爆模擬，取得了與試驗值吻合的數值計算結果，并分析了流場壓力波動特性，如圖2所示，圖2 (a) 為測試玻璃球體及監測裝置，圖2 (b) 為數值模擬7 MPa下內爆后的壓力曲線。

圖2 測試玻璃球體及監測裝置和數值模擬7 MPa下內爆后的壓力曲線Fig. 2 Test stand with glass sphere, blade and three pressure sensors numerical simulation of pressure curve after implosion at 7 MPa

從圖2(b)可以看到，內爆后產生的壓力大小隨著時間在不斷的發生變化，由于存在多次壓縮與反彈過程，壓力出現多個周期性峰值，對于周圍耐壓結構而言，環境載荷不再保持為靜載。因此，在進行耐壓殼優化設計時，要考慮動態載荷即載荷隨時間變化對優化結果的影響。圖2(b)中內爆后產生的壓力存在多個周期性峰值，從載荷的類別上，可以將內爆后產生的壓力以最高峰值為界分為兩個階段，從防護的角度研究水下耐壓結構優化設計，由于單體內爆瞬間產生極大水錘型沖擊，第一階段可通過設置耐壓結構防護裝置，防止陶瓷耐壓結構單體內爆導致的多體殉爆。第二階段可將周期性壓力變化簡化成簡諧載荷的形式，開展設計相關動載荷作用下的水下耐壓結構拓撲優化理論及方法研究。圖2(b)右上角為峰值后的壓力曲線放大圖，可以看到其載荷幅值波動在2 MPa以內，其載荷變化形式類似于簡諧載荷。從橫軸的時間尺度來看，周期約為1 ms，頻率約為6 000 rad/s，即動載荷的頻率較高。在多耐壓結構共存的情況下，若某一水下耐壓結構內爆產生周期性動態載荷，將會對其他耐壓結構造成影響，僅考慮靜載的耐壓結構傳統優化設計方法在該工況下難以適用。因此，研究設計相關動載荷作用下的水下耐壓結構優化設計，是具有重要的理論意義和工程價值的。

目前，動力學相關的拓撲優化主要集中在最大化結構特征頻率設計[5-7]，以及以頻率為約束的優化問題[8]。但是，關于結構動響應比如以結構的動柔順度為目標的研究較少。Ma等[9]首次將動柔順度的概念引入到結構動力學拓撲優化中，并以動柔順度為目標，利用均勻化方法研究結構動力學問題。Olhoff和Du[10-11]提出了頻率漸變的IF方法和GIF方法研究動柔順度相關問題，并且采用質量修正的插值模型去避免虛假模態問題，其結果表明該方法簡單有效，求得的頻率值準確。Liu等[12]給出了簡諧載荷下不同列式的拓撲優化問題，用ESO算法給出了相應的以動柔順度為目標的結構動力學優化設計，得到了較為理想的結果。Silva等[13]研究了不同因素對穩態強迫振動下的單材料結構拓撲優化的影響。以上研究均是關于固定載荷的優化設計問題，然而本文物理背景下研究的問題是設計相關動載荷作用下的優化問題，關于這類問題的研究相對較少，Olhoff和Du[10]和張暉等[14]利用其開發的邊界搜索算法來研究此類問題。張暉等[14]進行了受內壓的容器和水箱的優化設計，其采用RAMP(rational approximation of material properties)插值模型來避免虛假模態，并對結構拓撲的變化給出了定性的解釋。但是，上述文獻并未將設計相關動載荷作用下的優化問題引入到水下耐壓結構設計問題中，并且對于得到的結構沒有給出更為深入的解釋。

本文研究聚焦于設計相關動載荷下的水下耐壓結構拓撲優化問題，其中設計相關載荷的邊界搜索是區別于固定載荷動力學問題的關鍵環節。目前關于加載面搜索算法的研究主要分為兩類：一類是通過對節點密度的處理來搜索加載面，Du與Olhoff[15]提出利用貝塞爾樣條曲線代表加載邊界，由節點相對密度值插值得到的等值點來作為樣條曲線控制點的等值線方法和Zhang等[16]提出的單元節點加載法均屬于此類，該類方法的優點是易于實現，而且避免了載荷靈敏度分析；另一類通過構造多物理場來跟蹤壓力載荷的變化，Chen和Kikuchi[17]、Bourdin和Chambolle[18]分別提出的通過流體流動來識別壓力載荷表面的技術。Kumar等[19]利用巖石力學中的Darcy定律來處理設計相關載荷。Wang和Qian[20]提出的基于密度梯度的方法來跟蹤加載面。多物理場的優勢在于壓力直接作用于兩相交界面上，無需進行加載面搜索，但其有限元建模相對復雜。除上述兩類搜索算法外，也有一些學者進行了水平集方法下的加載面搜索算法研究，郭旭和趙康[21]，Jiang和Zhao[22]利用水平集演化技術，提出了一種設計相關載荷作用邊界的搜索方法，通過適當的數學變換，方便的處理施加在結構上的拓撲相關載荷，避免了在密度法中繁瑣的邊界提取工作。

采用Ibhadode等[23]提出的BILE(boundary identification-load evolution)模型作為邊界搜索算法，這一邊界識別和載荷演化的模型通過設置單元偽密度閾值來搜索加載面，該閾值在每次迭代中通過特定的方式遞增以進行邊界識別。同時，引入一個參數來定義兩個邊界識別步驟之間載荷演化的迭代次數，也可以控制優化的速度。BILE模型屬于第一類加載面搜索算法，無需進行載荷的靈敏度分析。因此，該模型易于應用，而且在80～100次迭代步后得到的數值結果都具有不錯的計算速度和可行性[23]。

通過將內爆后耐壓結構所受載荷簡化為周期性載荷，研究設計相關動載荷作用下水下耐壓結構的拓撲優化理論及方法，采用SIMP(solid isotropic microstructures with penalization)方法，選擇Ibhadode等[23]提出的BILE模型作為邊界搜索算法，并結合杜建鑌[24]提出的修正的材料插值模型來避免局部模態現象，采用靜柔度權重因子的拓撲優化數學模型來解決高頻下的水下耐壓結構拓撲優化設計，給出在不同頻率下的水下耐壓結構最優構型，指導耐壓結構的概念設計，為結構設計人員提供思路。

1 問題提出

1.1 有限元分析

考慮受簡諧變化的壓力設計相關動載荷拓撲優化問題，用有限元方法離散連續體結構時，結構在隨時間變化的力作用下的動力響應平衡方程可寫成：

(1)

不考慮阻尼C=0，且在諧響應的情況下可以得到:

(K-w2M)u=F

(2)

其中，M,K分別代表了結構的總質量矩陣和剛度矩陣，利用有限元方法可得：

(3)

其中，n是設計域中單元的總數。mi和ki是單元質量和剛度矩陣，通過式(4)、(5)計算：

(4)

(5)

其中，V是單元的體積，ρ和D分別是質量密度和本構矩陣，N和B表示單元形狀函數和應變位移矩陣，mc為單元的一致質量矩陣，ml為單元的集中質量矩陣。本文采用形如式(4)中的組合形式來減少計算誤差。

1.2 拓撲優化數學模型

基于密度法的連續體結構拓撲優化問題，將設計變量松弛為0～1之間任意值的連續變量優化問題。由于這樣的松弛方式將導致最優結構中出現中間密度的材料，為了減少灰度單元，常采用懲罰因子來對單元偽密度進行懲罰。最常采用的是SIMP方法，設給定設計區域為Ω，用N個單元離散，每個有限單元的相對密度為：xi(i=1,2……N)。對于各項同性材料，由SIMP模型可以得到單元的彈性模量：

(6)

式中：E0為固體材料的彈性模量；p為懲罰系數，在文中取3；Emin是為了避免總體剛度矩陣奇異性而引入的非零常數，取0.001。

本研究給定體積約束下最小化結構的動柔順度問題，當使用標準的SIMP法進行動力學優化時，存在局部模態(也稱虛假模態)現象。局部模態現象是指在優化過程中，設計域內會出現一些低密度區域，由于這些區域的剛度與質量比遠低于其他位置，會導致最后求解的振動模態出現劇烈的振蕩。因此，在動力學拓撲優化中需要對局部模態進行處理。Lazarus和Hagiwara[25]建議采用0.1作為低密度區的門檻值，Pedersen[6]提出一種修正的剛度插值模型；杜建鑌[24]借鑒Pedersen[6]提出的插值形式的思路，提出對低密度區的質量插值模型進行修正，采用杜建鑌[24]提出的質量修正模型：

me=xeρe,ρ≥0.1

(7)

me=xeqρe,ρ<0.1

(8)

其中，me為單元的質量，xe為單元的相對密度，ρe為單元的物理密度,q是進行單元質量插值時的懲罰因子，文中q取為 6。因此，優化問題的離散數學模型為：

(9)

當頻率接近或高于一階固有頻率時，采用上述優化模型無法得到可靠的結果(3.1節會給出相關算例來說明)，需引入靜柔度約束，或目標函數中引入含靜柔度的權重因子[13]，采用第二種方式，其數學模型為：

(10)

其中，F為設計相關動載荷的載荷幅值向量，U為位移幅值向量，K、M分別為整體剛度和整體質量矩陣，w為激勵頻率，Cd為結構動柔順度，Cs為靜柔度，計算KU=F得到，V為實際材料體積，V0為設計域體積，f為給定材料體分比，n為設計域中單元的個數。

1.3 靈敏度分析

為了解決式(9)中列出的優化問題，需要目標函數對設計變量的導數。采用優化準則法對優化問題進行求解。目標函數對設計變量的求導為：

(11)

根據KdU=F，可得：

(12)

(13)

代入式(11)可得：

(14)

對于BILE模型(第2節中將給出具體的搜索方法)，單元密度被用于定義加載面，壓力載荷直接作用于加載節點，每個等效力在每次迭代中大致保持相同的大小。在兩次迭代之間，由于加載面的長度和加載節點的數量同時增加或減少，每個等效力的大小都不會有顯著變化。因此，在整個優化過程中，每個等效節點力的大小保持不變[23]。載荷對設計變量的伴隨靈敏度為0，則靈敏度為：

(15)

即：

(16)

2 BILE模型

在設計相關載荷作用下結構拓撲優化的研究難點在于如何確定設計相關載荷作用的加載面。采用Ibhadode等[23]提出的BILE模型作為加載面搜索算法，因為其收斂速度快，能夠在80～100個迭代步給出與其他算法相似的結構，且本文為應用研究，BILE模型能夠更快給出一個較為準確的結果，有效提高工作效率。

該方法的基本流程為：

1)確定每次迭代的邊界節點。根據單元密度閾值以及每個節點相鄰單元的平均密度值確定邊界節點如圖3(a)和3(b)，第一次迭代時由于各單元密度值相等，邊界節點為設計域四周節點。而迭代過程中，某一節點n的選取公式為:

(17)

其中，xn,i為節點周圍的4個單元，xsn,l為當前節點鄰接的4個節點Sn的四周單元，若當前節點四周的單元密度符合公式的選取規則，則當前節點被標記為邊界節點，某一次的選取如圖3(b)所示。

2)從生成的邊界節點中選擇加載節點[如圖3(c)]形成加載面。在BILE模型中，選取最外側邊界點作為加載節點，具體是否選取某一邊的加載節點取決于初始邊界條件和邊界條件與加載節點關系。

圖3 節點n的鄰接節點以及某一次優化過程中的邊界節點和載荷節點Fig. 3 Adjacent node of node n, boundary node in a certain optimization process, load node in a certain optimization process

3)施加等效節點力。對于一個加載節點，該點上力的角度由連接兩個最近的加載節點到垂直正方向的直線的法線逆時針來計算。

由于是通過單元偽密度閾值并人為設定一定規則來搜索壓力加載面，不必求解力的伴隨靈敏度項，這樣不僅簡化了目標函數靈敏度的求解過程，而且也加快求解的速度，提高了整體的計算效率。

3 數值算例

由于實際結構均存在阻尼，阻尼對高頻響應有很強的抑制作用，對低頻響應的影響較小[24]，這里主要關注小于或略高于結構一階固有頻率的低頻共振優化問題，研究不同設計域和約束以及頻率條件對結構拓撲的影響。

Synergy HT酶標儀（美國BioTeK公司）；Mini-PROTEAN Tetra蛋白電泳儀、Trans-Blot SD半干轉膜系統（美國Bio-rad公司）；Direct-Q超純水儀（美國Millipore公司）；5417R高速冷凍離心機（德國Eppendorf公司）；IKA T18 basic勻漿器（德國ULTRATURRAX公司）；AX70顯微照相系統（日本Olympus公司）；ImageQuant LAS 4000全自動圖像分析系統（美國GE公司）。

采用優化準則法進行設計變量的更新，利用MATLAB中的eigs函數求解結構頻率，其求解對稱矩陣采用Lanczos算法。在本文算例中，Sigmund[26]、Sigmund和Maute[27]提出的靈敏度過濾方法被用來消除計算中的棋盤格，網格依賴性等數值不穩定現象。

3.1 受內壓的拱形結構

作為承受壓力載荷的結構優化問題中的經典算例，選擇內部受壓的拱形結構優化設計作為第一個算例，用以驗證搜索算法在動態問題上的適用性。初始的設計域及約束條件等如圖4所示。動載荷幅值為1.0，實體材料的材料屬性均采用無量綱化的取值，彈性模量取值為1.0，泊松比取值為0.3，質量密度ρ=1.0×10-6，本節算例的上述參數取值相同。設計區域離散成 40×20個正方形單元，且最大允許的材料體積占初始總體積的體積分數為50%。

圖4 初始設計域及邊界條件Fig. 4 Design domain and boundary conditions

此前，已有許多學者進行了設計相關載荷下的拓撲優化研究。圖5表明本文優化結果的拓撲形式與Zhang等[16]研究以及Xia等[28]研究的結果類似。表1開展了本文優化結果的柔度與Xia等[28]研究和Zhang等[16]研究結果的比較，當楊氏模量采用1.0，網格為100×50時，本文柔度值為10.02，小于表1中Xia等[28]研究的結果；當楊氏模量采用100，網格為40×20時，本文的柔度值為0.07，略小于表1中Zhang等[28]研究的結果。本文的方法獲得了較小柔度的拓撲結構，優于參考文獻的結果。

圖5 本文結果與不同文獻結果的對比Fig. 5 Comparison of optimization results between this paper and other references

表1 設計相關載荷下拓撲優化的數值示例對比：受內壓拱形結構Tab. 1 Comparison of numerical examples of topology optimization under design-dependent loads: arch structure subjected to internal pressure

圖6給出了采用式(9)的優化模型時隨激勵頻率變化目標函數以及拓撲形式的變化。可以看到，當激勵頻率低于500 rad/s時，隨著激勵頻率的增大，雖然優化結果的動柔度略有增大，但拓撲結構形式并無明顯差異。當激勵頻率為200 rad/s時，本文優化結果與張暉等[14]的研究結果對比如圖7(a)和7(c)所示，拓撲形式幾乎一致。隨著激勵頻率增大到600 rad/s時，由于頻率過大，求解并不收斂，中間結果如圖6所示。圖7(b)為激勵頻率800 rad/s時的中間結果，雖然與圖7 (d)張暉等[14]的研究結果類似，但繼續優化后求解并不收斂，最終結果發散；張暉等[14]在其論文中已經說明其結果不滿足體積約束，非可行解。由此可見，傳統的以單一動柔度為目標的拓撲優化數學模型存在缺陷。

圖6 不同激勵頻率下的優化結果對比曲線Fig. 6 Comparison curve of optimization results under different excitation frequencies

圖7 不同激勵頻率下本文優化結果與參考文獻[14]的對比Fig. 7 Comparison of optimization results between this paper and reference [14] with different excitation frequencies

因此，引入了如式(10)的目標函數——含靜柔度權重因子的優化模型來計算高頻下的結果，其有效性已在Silva等[13]研究中得到證明。該方法中η的選取對結果影響較大。表2給出了激勵頻率為800 rad/s時，不同η下優化結果，其中nlter為迭代步數。

表2 激勵頻率為800 rad/s時不同的η對結果的影響Tab. 2 The effect of different η on the result when excitation fequency is 800 rad/s

當η大于0.6時，拱形結構中會存在一些灰度單元，且結果并不對稱；隨著η的減小，優化結果呈現“0～1”分布，且η為0.4時結構的一階固有頻率為808.2 rad/s，比靜載時的結構頻率750.5 rad/s高；因此，在利用式(10)的優化模型進行受內壓的拱形結構優化時，建議η在0.2～0.4取值。

3.2 設計相關動載荷下的水下耐壓結構平面應變模型優化設計

圖8 設計域與邊界條件Fig. 8 The design domain and boundary conditions

圖9和圖10分別為不同激勵頻率下的優化結果以及激勵頻率200 rad/s時的優化過程。

圖9 六種激勵頻率下的最優拓撲Fig. 9 Optimal topology under six excitation frequencies

圖10 設計問題的優化過程Fig. 10 Optimization history for the design problem

圖9(a)中為靜水壓(激勵頻率為0 rad/s)下的優化結果，與以往學者[30-31]所得到靜水壓下的結果基本一致，其一階固有頻率為418.1 rad/s。圖9(b)～(f)中隨著激勵頻率的增大，由于算法振蕩，導致優化結構的拓撲構型在局部邊界處灰度單元減少，但優化結構的整體構型基本保持不變。在設計域中無固定端等強約束時，圖8所示的初始設計域在低于靜水壓下結構的一階固有頻率的不同激勵頻率下，圓環均為最優的結構形式，進一步驗證了圓柱殼作為耐壓結構的合理性。當激勵頻率高于靜水壓下結構的一階固有頻率(418.1 rad/s)時，采用式(10)的模型來進行優化。激勵頻率為500 rad/s時不同η下的結果對比如表3。

表3 激勵頻率為500 rad/s時η對結果的影響Tab. 3 The effect of η on the result when excitation fequency is 500 rad/s

從表4可以發現，當η為0.80時，表3中優化結果的一階固有頻率較圓環型耐壓結構小，且并不呈現出中心對稱的構型，不可取；隨著η減小到0.40，優化結果呈現較為清晰的“0～1”分布，且表4中顯示結構的一階固有頻率普遍高于圓環型耐壓結構的一階固有頻率，但大于0.30時，結構對稱性沒有小于0.20時好，在利用式(10)的優化模型進行設計相關動載荷下的水下耐壓結構平面應變模型優化設計時，建議η在0.05～0.20內取值。

表4 激勵頻率為500 rad/s時不同η下結果的第一階固有頻率Tab. 4 Fundamental frequency of results under different η when excitation fequency is 500 rad/s

3.3 設計相關動載荷下的水下耐壓結構軸對稱模型優化設計

除了上述開展耐壓結構平面應變模型的研究外，一些學者[30-31]進行了靜載環境軸對稱下水下耐壓結構的設計研究。在前人的基礎上，開展設計相關動載荷作用下水下耐壓結構的拓撲優化研究，其設計域與邊界條件等如圖11(a)所示。考慮到圖11(a)中的模型及邊界條件的對稱性，取 1/4 模型如圖11(b)所示作為設計域，并離散成120×48個正方形單元，且最大允許的材料體積占初始總體積的體積分數為35%。

圖11 設計域與邊界條件Fig. 11 The design domain and boundary condition

圖12(a)中結構的第一階固有頻率值為1 253.8 rad/s。從圖12(b)、12(c)及圖13優化結果可以看到，隨著激勵頻率的增大，固定約束處的材料聚集明顯，整體材料分布有靠近約束的趨勢。

圖12 三種激勵頻率下的1/4軸對稱模型最優拓撲Fig. 12 Optimal topology under three excitation frequencies of 1/4 axisymmetric model

圖13 設計問題的優化過程Fig. 13 Optimization history for the design problem

定義固定端約束兩側5列單元密度和為約束處面積，表5及圖14反映出隨著激勵頻率增大，約束處面積不斷增大，這與杜建鑌[24]所做的流道設計中隨著激勵頻率變大，頂部稍薄，支撐處有所加強有相似之處。為了解釋這種現象，考慮不同頻率下最終優化構型的第一階模態剛度與模態質量比值的變化：

圖14 不同激勵頻率下固定約束處體積Fig. 14 The volume of the fixed restraint under different excitation frequencies

表5 不同頻率下約束處體積Tab. 5 Restricted volume at different frequencies

(18)

其中，wn為結構的圓頻率，ur為結構的第r階模態向量，K、M分別為結構總的剛度矩陣和質量矩陣。由于頻率不同，結構的第r階模態剛度和第r階模態質量均在發生變化，無法很好地反映材料變化對結構模態剛度和結構模態質量的影響。因此，利用振動力學模態質量正則化方法[32]，可得相應的正則模態剛度為：

(19)

通過比較正則模態剛度的變化，可以比較材料分布變化對結構一階固有頻率的影響。圖15給出了不同激勵頻率下優化結構的第一階正則模態剛度變化情況。隨著激勵頻率的增大，優化結構的正則模態剛度也在不斷增大。因此，拓撲形式上的不同，實際上反映的是優化過程中材料的增加和去除對結構一階模態剛度與一階模態質量比值的影響，若某位置材料的增加或去除能使其比值變大，則能夠提高結構的一階固有頻率。圖16優化結果的一階固有頻率隨激勵頻率的變化情況也驗證了上述觀點。因此，對于多球交接的水下耐壓結構形式，在多球交接的位置采用環向加強肋可以提高結構的動力學特性，使耐壓殼能夠適應更加復雜的水下環境，提高潛水器抗沖擊性能。

圖15 優化結果一階正則模態剛度隨激勵頻率變化Fig. 15 The first-order canonical mode stiffness varies with the excitation frequency

圖16 優化結果一階固有頻率隨激勵頻率變化Fig. 16 The first-order frequency varies with the excitation frequency

當激勵頻率高于靜水壓下結構的第一階固有頻率(1 253.8 rad/s)時，采用式(10)的數學模型來進行優化。表6給出了激勵頻率為1 300 rad/s時不同η的拓撲優化結果。

表6 激勵頻率為1 300 rad/s時不同η對結果的影響Tab. 6 The effect of different η on the result when excitation fequency is 1 300 rad/s

當η大于0.20時，優化結果出現較多的灰度單元，結果不可取；隨著η的減小，優化結果呈現較為清晰的“0～1”分布。因此，在利用式(10)的優化模型進行設計相關動載荷作用下的水下耐壓結構軸對稱模型優化設計時，建議η在0.01～0.15內取值。

為了進一步說明優化結構的應用可行性，圖17(a)、17(b)分別給出了MIT 水下潛水器設計團隊設計的一種多球交接耐壓殼體模型[33]及張建等[34]提出的多蛋交接耐壓殼三蛋形殼截面，本文得到的多球交接耐壓結構如圖17(c)所示。可見，本文中得到的多球交接耐壓結構球殼交接處與圖17(a)、17(b)非常相似，在連接處會進行局部加強，結構一階固有頻率更高，動力學性能更好，說明本文對于水下耐壓結構的探索具有一定的工程應用價值。

圖17 考慮設計相關動載荷得到的優化結果與其他文獻結果的比較Fig. 17 Comparison of optimization results between this paper under design-dependent dynamic loading and other references

4 結 語

通過將內爆后耐壓結構第二階段所受高頻率的周期性變化載荷簡化為簡諧載荷，基于Ibhadode等提出的BILE加載面搜索算法，研究了設計相關動載荷作用下的水下耐壓結構拓撲優化研究及方法。同時，采用了靜柔度權重因子的拓撲優化數學模型來解決高頻下的水下耐壓結構拓撲優化設計，數值算例驗證算法有效性。

著重研究了設計相關動載荷作用下的水下耐壓結構拓撲優化問題，探索了水下耐壓結構新形式。低頻時，圓環型耐壓結構與靜水壓下相似，但對于多球交接耐壓結構，隨著頻率增加，材料往約束處聚集以提升結構的動力學性能；接近或高于一階頻率時，圓環型耐壓結構與多球交接型耐壓結構形式均與靜水壓下的結果存在明顯差異，多球交接型耐壓結構與現有的耐壓結構概念設計結果相似。最終，探究了其在工程實際中的應用，對新型水下耐壓結構的概念設計具有積極的意義。

本文的研究只是對設計相關動載荷作用的水下耐壓結構優化設計的初步探究，在未來的工作中，將致力于該方法的三維拓展及屈曲約束下的水下耐壓結構優化設計以及其對應的三維拓撲優化問題。