基于HCA的改進(jìn)收斂準(zhǔn)則優(yōu)化算法研究

向曉波 田啟華，2 杜義賢，2，3

（1.三峽大學(xué) 機(jī)械與材料學(xué)院，湖北 宜昌 443002；2.三峽大學(xué) 水電機(jī)械設(shè)備設(shè)計與維護(hù)湖北省重點(diǎn)實(shí)驗(yàn)室，湖北 宜昌 443002；3.三峽區(qū)域能源裝備三峽大學(xué)協(xié)同創(chuàng)新中心，湖北 宜昌 443002）

隨著結(jié)構(gòu)優(yōu)化技術(shù)研究的深入發(fā)展，結(jié)構(gòu)拓?fù)鋬?yōu)化已成為結(jié)構(gòu)優(yōu)化領(lǐng)域的研究熱點(diǎn)之一［1］.與尺寸優(yōu)化、形狀優(yōu)化相比，結(jié)構(gòu)拓?fù)鋬?yōu)化具有更大的設(shè)計自由度和更廣闊的設(shè)計空間以滿足苛刻的設(shè)計要求［2］.目前，常用的連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化方法大多基于梯度優(yōu)化算法，在優(yōu)化模型的求解中需要反復(fù)進(jìn)行敏度分析，不利于求解復(fù)雜的結(jié)構(gòu)優(yōu)化問題.Tovar等［3］將元胞自動機(jī)理論引入結(jié)構(gòu)拓?fù)鋬?yōu)化方法中，提出了混合元胞自動機(jī)（Hybrid Cellular Automata，HCA）算法.周珍珍等［4］研究混合元胞自動機(jī)方法處理結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計的問題.王安麟［5］提出了基于進(jìn)化建立元胞自動機(jī)局部直接規(guī)則的方法，有效解決了局部間接規(guī)則存在的局限性問題.

在求解數(shù)學(xué)中的最優(yōu)化問題及工程中結(jié)構(gòu)拓?fù)鋬?yōu)化問題中，能夠直接解析求解得到結(jié)果的非常有限，一般都需要通過迭代法進(jìn)行求解［6］.而在算法的迭代過程中，收斂準(zhǔn)則的選用將直接影響優(yōu)化設(shè)計的精度及結(jié)果.目前，采用收斂準(zhǔn)則求解問題的思想在人工智能和自動化控制等領(lǐng)域都獲得較為廣泛的應(yīng)用［7］.因此，研究收斂準(zhǔn)則具有十分重要的理論意義和應(yīng)用價值.

本文在總結(jié)現(xiàn)有收斂準(zhǔn)則算法的基礎(chǔ)上，結(jié)合范數(shù)矩陣的理論，提出了一種基于HCA的改進(jìn)收斂準(zhǔn)則優(yōu)化算法，來抑制拓?fù)鋬?yōu)化過程中可能出現(xiàn)的灰度單元問題.并通過經(jīng)典二維拓?fù)鋬?yōu)化數(shù)值算例，探討所提出的收斂準(zhǔn)則在解決拓?fù)鋬?yōu)化中灰度單元問題上的可靠性和有效性.

1 連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化建模

1.1 密度插值函數(shù)模型

變密度法中常用的插值模型是固體各向同性懲罰微結(jié)構(gòu)模型（solid Isotropic microstructure with penalization，SIMP）［8］，SIMP材料插值方法主要通過引入懲罰因子對單元的中間密度值進(jìn)行懲罰，本文基于各向同性材料，建立非線性的冪函數(shù)來表示材料的彈性模型與單元相對密度之間的顯式對應(yīng)關(guān)系

式中，E0和ρ0分別是實(shí)體材料的彈性模量和密度，ρe和Ee是單元材料的密度和相對密度，Ee和Emin是固體和孔洞材料的彈性模量，ΔE0＝E0-Emin，Emin＝E0／1000，p是懲罰因子.

1.2 基于HCA的拓?fù)鋬?yōu)化數(shù)學(xué)模型

HCA算法是將CA理論和有限元方法結(jié)合起來的一種無梯度優(yōu)化算法.在CA模型迭代演化時，將單元的平均應(yīng)變能密度值與目標(biāo)值之間的差值作為反饋信號，優(yōu)化的目的是使反饋信號值最小化.本文采用變密度法中SIMP材料插值模型，建立基于HCA的結(jié)構(gòu)拓?fù)鋬?yōu)化模型為

2 改進(jìn)的算法收斂準(zhǔn)則

收斂準(zhǔn)則就是判斷在一定的精度要求下，數(shù)值求解過程中所得到的迭代點(diǎn)是否正是需要的.不合適的收斂準(zhǔn)則可能會導(dǎo)致最終的迭代結(jié)果無法收斂或者得到不好的優(yōu)化結(jié)果.

理想的迭代收斂準(zhǔn)則是將結(jié)構(gòu)的質(zhì)量變化作為收斂準(zhǔn)則，質(zhì)量的表達(dá)式為

式中，ρe（t）和ρ0分別是變質(zhì)量密度和固體材料密度，M（t）是第t次迭代時結(jié)構(gòu)的質(zhì)量.

目標(biāo)函數(shù)的收斂準(zhǔn)則取決于設(shè)計變量更新時所采用局部控制規(guī)則的類型，Tovar在文獻(xiàn)［3］中將迭代過程中結(jié)構(gòu)的質(zhì)量變化達(dá)到穩(wěn)定狀態(tài)作為優(yōu)化模型的收斂準(zhǔn)則，為避免目標(biāo)函數(shù)出現(xiàn)過早收斂，將連續(xù)兩次迭代時質(zhì)量絕對變化的平均值小于給定值作為收斂條件，具體的數(shù)學(xué)表達(dá)式為

式中，ΔM（t）＝M（t）-M（t-1）是第t次迭代時結(jié)構(gòu)質(zhì)量的變化量，ε是給定容差.

范數(shù)（norm）是數(shù)學(xué)中的一種基本概念，在泛函分析中，范數(shù)是一種定義在賦范線性空間中函數(shù)，滿足相應(yīng)條件后的函數(shù)都可以被稱為范數(shù)［9-10］.最常用的范數(shù)就是P范數(shù).若x＝［x1，x2，x3，…，xn］T，那么

本文引入P范數(shù)理論的概念，通過判斷在這些方向上充分接近t＋1迭代點(diǎn)處，達(dá)到判別最優(yōu)點(diǎn)的目的，基于這種思想，提出了一種新的算法收斂準(zhǔn)則.當(dāng)相鄰兩次設(shè)計變量的迭代結(jié)果足夠小時，收斂準(zhǔn)則的數(shù)學(xué)表達(dá)式為

式中，M（t＋1）是第t＋1次迭代時結(jié)構(gòu)質(zhì)量.

本文提出的收斂準(zhǔn)則引入P范數(shù)后，使前后兩個迭代點(diǎn)充分絕對接近，具有更多下一時刻t＋1迭代點(diǎn)周圍的信息，因?yàn)樵诘^程中能夠搜索到極小值點(diǎn)，使收斂準(zhǔn)則有更高的可靠性.

3 基于HCA的結(jié)構(gòu)拓?fù)鋬?yōu)化建模與求解

3.1 HCA拓?fù)鋬?yōu)化方法

HCA算法是一種無梯度的優(yōu)化方法，將元胞自動機(jī)方法與有限元分析結(jié)合進(jìn)行結(jié)構(gòu)優(yōu)化，最初用于模擬骨重塑過程，根據(jù)骨骼功能自適應(yīng)理論進(jìn)行骨骼重新塑造，使骨骼具有均勻應(yīng)變能密度.

HCA模型的基本組成有元胞柵格、每個元胞的狀態(tài)集合和狀態(tài)更新規(guī)則集合.每個元胞的狀態(tài)由當(dāng)前時刻t狀態(tài)及其鄰域元胞的狀態(tài)決定下一個時刻t＋1該元胞的狀態(tài).元胞第j種的狀態(tài)更新規(guī)則可以表示為

式中，αe（t），αe＋Δ1（t），αe＋ΔN（t）表示元胞e的鄰居元胞.

CA本質(zhì)上是時空離散化的數(shù)學(xué)模型，由有限個元胞單元構(gòu)成，每個元胞單元具有特定的狀態(tài)值.在某一時刻（或某次迭代），每個元胞的狀態(tài)值僅與與其相鄰的幾個元胞的狀態(tài)值相關(guān).狀態(tài)更新規(guī)則需要收集每個元胞領(lǐng)域元胞的狀態(tài)信息.本文采用包含8個鄰居元胞的Moore鄰居類型，如圖1（a）所示.為了定義位于設(shè)計域邊界上的元胞更新規(guī)則，必須先定義不同的邊界條件.本文選用固定邊界條件，即設(shè)計域邊界外的狀態(tài)為零，如圖1（b）所示.

圖1 元胞自動機(jī)的領(lǐng)域和邊界條件

拓?fù)鋬?yōu)化迭代過程中，每個元胞的設(shè)計變量可定義為單元密度xe，每個元胞的狀態(tài)場變量可以定義為應(yīng)力、應(yīng)變、應(yīng)變能、互應(yīng)變能或者它們的函數(shù)Se.這時，元胞的狀態(tài)可以表示為

式中，αe（t）是元胞的狀態(tài)，xe（t）是元胞的設(shè)計變量，Se（t）是元胞的狀態(tài)場變量.

式中，Ue（t）為元胞應(yīng)變能密度值；Ui為鄰居元胞的應(yīng)變能密度值；N為元胞鄰居數(shù)量.

3.2 局部控制規(guī)則

設(shè)計域中的材料分布由局部控制規(guī)則確定，它尋求狀態(tài)場變量平均值與狀態(tài)場變量設(shè)定值差值減小到最小.局部控制規(guī)則有兩位置、線性、積分和微分控制.較復(fù)雜的控制策略是比例-積分-微分（Proportional-Integral-Derivative，PID）反饋控制，該控制策略已廣泛應(yīng)用于控制應(yīng)用問題中，數(shù)學(xué)表達(dá)式為

式中，cP，cI，cD分別是比例，積分和微分的權(quán)系數(shù).本文采用局部比例控制.

3.3 HCA算法的求解流程

本文采用HCA算法求解優(yōu)化模型，求解流程如圖2所示.在迭代計算時，CA狀態(tài)的平衡條件-＝0與滿應(yīng)力設(shè)計準(zhǔn)則一致，使非空元胞處于滿應(yīng)力狀態(tài).

Step1：定義設(shè)計域、材料屬性、載荷等參數(shù)，初始化元胞值x（0），離散設(shè)計域?yàn)橛邢迋€CA單元；

Step2：采用有限元方法進(jìn)行結(jié)構(gòu)分析，根據(jù)式（11）分別計算元胞單元在第t次迭代時應(yīng)變能Ue和平均應(yīng)變能

Step3：采用基于反饋控制的CA局部控制規(guī)則更新設(shè)計域內(nèi)材料的分布；

Step4：根據(jù)改進(jìn)的收斂準(zhǔn)則和給定容差ε，判斷迭代是否收斂.若迭代滿足收斂準(zhǔn)則，則得到設(shè)計域內(nèi)材料的最優(yōu)分布，迭代終止；否則，返回到Step2.

圖2 混合元胞自動機(jī)算法求解流程

4 經(jīng)典數(shù)值算例及結(jié)果分析

圖3為懸臂梁的設(shè)計域，結(jié)構(gòu)初始設(shè)計域離散為80×50個單元.材料彈性模量E＝210GPa，泊松比μ＝0.3；懸臂梁的左端為固支約束，豎直向下的集中載荷F＝1kN作用在右邊界中心位置.單元相對密度的初始值取為1.0（實(shí)體材料）；利用有限元法求解初始結(jié)構(gòu)的總應(yīng)變能U0＝19.2Nmm，單元的局部應(yīng)變能密度目標(biāo)值＝U0／（80×50）＝4.8×10-3Nmm／mm3.

圖3 懸臂梁的設(shè)計域

1）范數(shù)值對優(yōu)化結(jié)果的影響

針對懸臂梁結(jié)構(gòu)拓?fù)鋬?yōu)化模型，采用本文提出的改進(jìn)收斂準(zhǔn)則，在相同的初始條件和約束下，分別選用不同的范數(shù)值，通過HCA算法對優(yōu)化模型進(jìn)行求解，得到圖4所示的優(yōu)化結(jié)果.

圖4 不同范數(shù)下拓?fù)鋬?yōu)化結(jié)果

從圖4可知，隨著范數(shù)值的不斷增大，本文提出的收斂準(zhǔn)則得到的拓?fù)浣Y(jié)構(gòu)的應(yīng)變能和灰度單元會逐漸降低.當(dāng)范數(shù)值小于20時，拓?fù)鋱D形會出現(xiàn)嚴(yán)重的灰度單元現(xiàn)象.當(dāng)范數(shù)值在60和100之間時，應(yīng)變能和灰度單元數(shù)量開始出現(xiàn)略微下降.因此，當(dāng)采用本文提出的收斂準(zhǔn)則時，選用合理的范數(shù)值能得到更低的應(yīng)變能值且能更好地抑制灰度單元的出現(xiàn).

2）灰度單元問題

針對本文提出的收斂準(zhǔn)則，采用3種不同的范數(shù)值，分別為20，80，100.對模型進(jìn)行優(yōu)化求解，與采用文獻(xiàn)［3］中的收斂準(zhǔn)則得到的拓?fù)鋬?yōu)化結(jié)果進(jìn)行對比，見表1.

表1 兩種收斂準(zhǔn)則下拓?fù)鋬?yōu)化結(jié)果對比

由表1可知，當(dāng)范數(shù)值為20時，得到的優(yōu)化結(jié)果中出現(xiàn)大量的灰度單元.當(dāng)范數(shù)值為80或100時，得到的拓?fù)鋱D形跟經(jīng)典拓?fù)鋱D形相似.當(dāng)范數(shù)值為100時，得到應(yīng)變能值20.983和灰度單元數(shù)目116，比采用文獻(xiàn)［3］中收斂準(zhǔn)則得到的應(yīng)變能值21.006和灰度單元數(shù)目130，分別減少了0.11%和10.8%.因此，采用本文提出的收斂準(zhǔn)則得到的應(yīng)變能更小，并能更好地抑制灰度單元的出現(xiàn).

5 結(jié) 論

本文對CA和拓?fù)鋬?yōu)化的結(jié)合技術(shù)進(jìn)行了探討，借助范數(shù)矩陣的理論，提出了一種改進(jìn)的拓?fù)鋬?yōu)化收斂準(zhǔn)則優(yōu)化算法.通過經(jīng)典二維拓?fù)鋬?yōu)化算例，分析了不同的范數(shù)值對拓?fù)鋬?yōu)化結(jié)果的影響，并與采用文獻(xiàn)［3］中收斂準(zhǔn)則得到的拓?fù)鋬?yōu)化結(jié)果進(jìn)行對比，證明了本文提出的收斂準(zhǔn)則能夠有效地抑制灰度單元的出現(xiàn)，取得更好的應(yīng)變能目標(biāo)值，并最終能夠得到輪廓較清晰的拓?fù)鋬?yōu)化圖形.

［1］左孔天.連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化理論與應(yīng)用研究［D］.武漢：華中科技大學(xué)，2004.

［2］杜義賢，王 偉，李 然，等.結(jié)構(gòu)拓?fù)鋬?yōu)化中敏度過濾技術(shù)研究［J］.工程設(shè)計學(xué)報，2012，19（1）：20-24.

［3］Tovar A，Patel N M，Niebur G L，et al.Topology Optimization Using a Hybrid Cellular Automaton Method with Local Control Rules［J］.Journal of Mechanical Design，2006，128（6）：1205-1216.

［4］周珍珍.基于混合元胞自動機(jī)方法的結(jié)構(gòu)拓?fù)鋬?yōu)化研究［D］.武漢：華中科技大學(xué)，2009.

［5］王安麟，姜 濤.基于進(jìn)化元胞自動機(jī)的結(jié)構(gòu)拓?fù)鋬?yōu)化［J］.機(jī)械工程學(xué)報，2005，41（2）：1-5.

［6］孫云龍.關(guān)于最優(yōu)化問題的算法收斂準(zhǔn)則［J］.四川師范大學(xué)學(xué)報：自然版，1999，22（3）：323-325.

［7］唐煥文，秦學(xué)志.實(shí)用最優(yōu)化方法［M］.大連：大連理工大學(xué)出版社，2000.

［8］Bendse M P，Sigmund O.Topology Optimization：Theory，Methods，and Applications［M］.Berlin Heidelberg，New York：Springer，2003.

［9］趙瑞珍，秦 周，胡紹海.一種基于特征值分解的測量矩陣優(yōu)化方法［J］.信號處理，2012，28（5）：653-658.

［10］徐 仲，張凱院.矩陣論簡明教程［M］.北京：科學(xué)出版社，2005.