應(yīng)該為結(jié)構(gòu)拓撲優(yōu)化的設(shè)計變量正名和處理方法順言

彭細榮， 隋允康

(1.湖南城市學(xué)院土木工程學(xué)院， 湖南 益陽 413000；2.北京工業(yè)大學(xué)工程數(shù)值模擬中心， 北京 100124)

?

應(yīng)該為結(jié)構(gòu)拓撲優(yōu)化的設(shè)計變量正名和處理方法順言

彭細榮1， 隋允康2

(1.湖南城市學(xué)院土木工程學(xué)院， 湖南 益陽 413000；2.北京工業(yè)大學(xué)工程數(shù)值模擬中心， 北京 100124)

為了澄清拓撲化概念及促進發(fā)展，從基本概念入手，闡述了為拓撲變量正名、為有關(guān)映射順言的必要性和迫切性. 連續(xù)體結(jié)構(gòu)拓撲優(yōu)化方面的研究，以往學(xué)者們不注意厘清拓撲變量的概念，沒有探討其獨立地位，而是依附于結(jié)構(gòu)優(yōu)化低層次，采用懲罰函數(shù)使變量趨于0/1，有研究者于2004年意識到應(yīng)當(dāng)引進Heaviside函數(shù)取代懲罰函數(shù). 其實，筆者在1996年提出的獨立連續(xù)映射(independent continuous mapping，ICM)方法中已經(jīng)徹底解決了這些問題：定義獨立層次的拓撲變量；采用階躍函數(shù)及其逆函數(shù)的近似逼近. 本文闡述了ICM方法的概念，并以頻率約束拓撲優(yōu)化為例詳述了建模及求解. 算例表明ICM方法因概念清晰顯示出明顯優(yōu)勢. 《論語·子路篇》曰：“名不正則言不順，言不順則事不成.”正名方可順言，是時候該為結(jié)構(gòu)拓撲優(yōu)化的設(shè)計變量正名和處理方法順言了.

連續(xù)體結(jié)構(gòu)拓撲優(yōu)化；拓撲設(shè)計變量的獨立化；近似逼近函數(shù)；懲罰函數(shù)；階躍函數(shù)；Heaviside函數(shù)；獨立連續(xù)映射(independent continuous mapping，ICM)方法；變密度法

連續(xù)體結(jié)構(gòu)拓撲優(yōu)化是指通過子域的“有”或“無”的拓撲布局，在滿足一定約束條件下，求得使某目標最優(yōu)的結(jié)構(gòu)構(gòu)型[1]. 約束及目標可取為結(jié)構(gòu)強度、剛度等性能指標或材料用量等經(jīng)濟指標.

目前使用的連續(xù)體結(jié)構(gòu)拓撲優(yōu)化大部分是構(gòu)建在基結(jié)構(gòu)(ground structure)之上的[2]，其拓撲優(yōu)化模型在本質(zhì)上屬于大型離散規(guī)劃問題. 對離散模型的求解目前沒有高效算法，智能算法，如模擬退火算法(simulated annealing，SA)[3]、遺傳算法(genetic algorithm，GA)[4]、粒子群優(yōu)化算法(partical swarm optimization algorithm，PSOA)[5]等，其普遍特點是收斂速度都較慢. 啟發(fā)式的進化結(jié)構(gòu)優(yōu)化(evolutionary structural optimization，ESO )[6]法按某種準則，對單元進行刪除或增加，其收斂速度也較慢，且最優(yōu)構(gòu)型受準則參數(shù)影響大等嚴重困擾[2,7].

為了克服離散變量優(yōu)化問題的求解困難，出現(xiàn)了用連續(xù)變量“代替”或“逼近”離散變量的做法. 連續(xù)變量“代替”離散變量的大部分做法是將拓撲層次優(yōu)化問題降低為低層次尺寸或形狀優(yōu)化問題，用連續(xù)的材料性能變量、截面變量或形狀變量取代離散的拓撲變量. 主要可分為兩大類方法：一類方法是依據(jù)設(shè)計變量在設(shè)計區(qū)域內(nèi)的分布場來確定，典型方法有均勻化方法[8]、變密度法[9]、變厚度法[10]等；另一類方法是在設(shè)計區(qū)域內(nèi)通過孔洞的形成、合并及孔洞與邊界的形狀改變來確定，典型方法有水平集法[11]、相場法[12]及拓撲導(dǎo)數(shù)法[13-14]等.

由于不同方法間的參考與借鑒，不同方法間呈現(xiàn)交叉、融合的趨勢. 如變密度法的密度過濾所采用的投影技術(shù)類似于水平集方法中的水平集函數(shù)[15]，可看作一種參數(shù)化的水平集方法[16]. 使用人工假想材料描述邊界的水平集方法，“空”域內(nèi)單元指定一個很小的人工密度，“實”域內(nèi)單元指定人工密度為1，邊界附近單元人工密度依據(jù)Heaviside 投影確定取適當(dāng)?shù)闹虚g密度值，其實質(zhì)是在用水平集函數(shù)控制人工密度場分布，而非真實的有限元網(wǎng)格邊界[17]，可看作一種變密度法[16]. 基于敏度分析及過濾技術(shù)改進的雙向進化結(jié)構(gòu)優(yōu)化(bi-directional evolutionary structural optimization，BESO)方法，已不完全是初始的基于生物進化啟發(fā)而提出的完全離散的優(yōu)化算法[6,18]，可看作一種采取離散更新策略[16]的固體各向同性材料懲罰(solid isotropic material with penalization，SIMP)[9]法.

從連續(xù)體結(jié)構(gòu)拓撲優(yōu)化方法的發(fā)展來看，無論是基于投影技術(shù)的變密度法，還是基于人工材料描述邊界的水平集方法，材料方式描述拓撲或形狀方式描述拓撲已經(jīng)開始融合，基于設(shè)計變量分布場來確定結(jié)構(gòu)拓撲布局的概念已經(jīng)形成. 然而，不幸的是，由于歷史發(fā)展的原因，這個設(shè)計變量分布場卻生硬地依附在人工材料密度或水平集函數(shù)所描述的邊界形狀等一些低層次材料或形狀優(yōu)化變量上，使得在問題描述及建模上出現(xiàn)了各種所謂的雜交方法. 既然各種方法開始融合，是時候應(yīng)努力尋求各種方法歸一的可能，即各類初始從不視角提出的方法，是否有可能在某種理論框架下進行統(tǒng)一.

突破口應(yīng)回歸到連續(xù)體結(jié)構(gòu)拓撲優(yōu)化問題的本源上來. 一是選取什么量作為設(shè)計變量？二是如何利用基于連續(xù)變量的優(yōu)化算法來解決大規(guī)劃離散優(yōu)化的困難？

本質(zhì)上，上述概念和處理方法相關(guān)的問題，早在1996年由隋允康提出的拓撲優(yōu)化獨立連續(xù)映射(independent contimous mapping，ICM)方法已經(jīng)解決：一是ICM方法定義了獨立層次的拓撲變量；二是通過階躍函數(shù)及其逆函數(shù)的近似逼近，在連續(xù)拓撲變量的基礎(chǔ)上建立及求解拓撲優(yōu)化模型[19].

連續(xù)體結(jié)構(gòu)拓撲優(yōu)化發(fā)展之初，對于這2個問題，大多方法將設(shè)計變量依附于低層次的變量上，如板的厚度、微結(jié)構(gòu)尺寸、人工材料密度、水平集函數(shù)描述的隱式邊界形狀等，通過這樣處理，以達到將一個大規(guī)劃離散優(yōu)化問題轉(zhuǎn)化為連續(xù)變量優(yōu)化問題求解的目的. Sigmund等[16]提到基于密度的各類拓撲優(yōu)化方法中，設(shè)計變量場為純數(shù)學(xué)意義上的概念. 以人工密度來表征子域的“有”或“無”，采用懲罰函數(shù)的概念來實現(xiàn)大部分設(shè)計變量取0/1，其實，這種處理方式是算不上純數(shù)學(xué)的拓撲優(yōu)化變量. Heasivide投影法[15,20]用于變密度法中過濾處理以得到無棋盤格及網(wǎng)格依賴現(xiàn)象、拓撲清晰(即過渡單元很少)的最優(yōu)拓撲，其用連續(xù)函數(shù)近似Heasivide函數(shù)(其為階躍函數(shù)的別名)類似于ICM法中磨光函數(shù)逼近于階躍函數(shù)，但其拓撲優(yōu)化建模仍是變密度法，設(shè)計變量仍是人工密度. 到目前為止，綜合國內(nèi)外相關(guān)研究文獻來看，除ICM方法外，還沒有其他方法提出過獨立層次的拓撲變量及相應(yīng)的建模處理方法.

近年來，一些學(xué)者敏銳地意識到還原拓撲變量獨立地位的重要意義，ICM方法得到他們的關(guān)注，并且做出了可喜的研究工作，如龍凱等[21-22]基于ICM方法提出了基于節(jié)點獨立變量的方法和混合插值建模方法，后來定義物質(zhì)點拓撲變量為[0 , 1]上用于表征物質(zhì)點及其領(lǐng)域存在與否的實數(shù)[23]，以拓撲變量來表征一個“物質(zhì)點”的“有”或“無”，從概念上來講還是ICM方法所提出的獨立連續(xù)拓撲變量. 鄧果等[24]基于ICM方法，研究了位移線性近似式和應(yīng)力約束轉(zhuǎn)換表達式以及有效結(jié)構(gòu)信息到結(jié)構(gòu)最大設(shè)計域信息的映射轉(zhuǎn)換方法.

連續(xù)體結(jié)構(gòu)拓撲優(yōu)化發(fā)展到今天，可謂方興未艾，卻在最基本概念上如此良莠紛繁而未能趨同. 然而，在建筑結(jié)構(gòu)形態(tài)創(chuàng)構(gòu)領(lǐng)域，研究者卻自然應(yīng)用了拓撲變量的概念來描述拓撲優(yōu)化問題[25]，作為結(jié)構(gòu)拓撲優(yōu)化的一個延展應(yīng)用的領(lǐng)域，卻在基本概念的使用上更為自覺及合理，這不得不引起結(jié)構(gòu)拓撲優(yōu)化領(lǐng)域研究者的反思乃至自?。菏菚r候該為結(jié)構(gòu)拓撲優(yōu)化設(shè)計變量正名了，同時也該為相應(yīng)處理方法順言了.

1 獨立層次拓撲設(shè)計變量定義及拓撲優(yōu)化模型

求解連續(xù)體結(jié)構(gòu)拓撲優(yōu)化問題，旨在確定設(shè)計區(qū)域內(nèi)各子域的“有”或“無”，應(yīng)該怎樣定義拓撲設(shè)計變量呢？

ICM方法解決這一問題時，首先考察了低層次的截面及形狀優(yōu)化設(shè)計變量的取法：在截面優(yōu)化中，需要設(shè)計的是桿件或板的截面，其設(shè)計變量通常取桿件的截面尺寸或面積、板的厚度等；在形狀優(yōu)化中，需要設(shè)計的是結(jié)構(gòu)的形狀，其設(shè)計變量通常取骨架類結(jié)構(gòu)的連接節(jié)點坐標或連續(xù)體結(jié)構(gòu)的邊界或孔洞的有限元網(wǎng)格節(jié)點坐標或形狀描述函數(shù)參數(shù)，亦即取“形狀”作為設(shè)計變量. 到了拓撲優(yōu)化層次，設(shè)計變量不再是具體的幾何量或物理量，而是各子域的“有”或“無”，很自然，各子域“有”或“無”的拓撲狀態(tài)就可以分別用0和1兩個數(shù)來表征，即：

子域——一個幾何點鄰域內(nèi)足夠小的區(qū)域.

子域的拓撲變量是離散量0或1，分別代表此子域為“無”與“有”.

不難看出，上述定義的拓撲設(shè)計變量是一個純粹的數(shù)學(xué)量. 對連續(xù)體結(jié)構(gòu)拓撲優(yōu)化問題，作為子域“有”或“無”的表征，也可以稱為子域拓撲分布場的平均值. 可稱為離散拓撲設(shè)計變量或拓撲分布場的離散平均值. 這里強調(diào)離散，是為了區(qū)分ICM方法定義的連續(xù)拓撲變量.

可以看到，由于定義了表征子域“有”或“無”的拓撲狀態(tài)量，結(jié)構(gòu)拓撲的描述就變得非常清晰及自然. 選取其他物理量作為拓撲優(yōu)化模型的設(shè)計變量，如變厚度法中的板厚度、均勻化方法中的微結(jié)構(gòu)尺寸、變密度法中的人工密度等均沒有直接反映拓撲優(yōu)化的本質(zhì)，應(yīng)是拓撲優(yōu)化發(fā)展初期，為解決大規(guī)模離散優(yōu)化模型的“組合爆炸”困難，將問題轉(zhuǎn)化為基于連續(xù)變量的規(guī)劃模型而采取的不得已的一些物理轉(zhuǎn)化方式.

ICM方法則不依賴物理轉(zhuǎn)化，而是把離散的0/1變量延拓到[0,1]，0與1之間的中間實數(shù)值類似于模糊集合論中的隸屬度，表示對于0或1的靠近程度. 這就是連續(xù)的拓撲變量.

另一方面，發(fā)展到現(xiàn)在，變密度法由于其廣泛的研究及應(yīng)用，“密度”一詞一直沿用，但其所表達的意義已由初始的人工密度的起始概念悄悄發(fā)生了變化，其意義已非常接近于ICM方法中的連續(xù)拓撲變量. 對比一下變密度法中用“密度”這個詞匯來表征子域“有”或“無”的拓撲狀態(tài)，其值取值在[0,1]，利用此人工密度變量，建立單元體積及單元彈性模量的函數(shù)關(guān)系. 很明顯，變密度法在此處使用密度這一詞匯是極易與真正的材料密度概念產(chǎn)生混亂的. 為什么要用“密度”來說一個子域的“有”或“無”，讓人感覺非常奇怪，難道真正的材料密度必須用人造的、人工的、假的或偽的“密度”表達嗎？

更重要的是，伴隨著純數(shù)學(xué)的拓撲設(shè)計變量由1變成0，單元所具備的所有的幾何量、物理量到從有變成無，不僅僅是密度從有變成無，因此，單獨關(guān)注密度，并不是公正的. 可見，最合理的做法是采用純數(shù)學(xué)的拓撲設(shè)計變量，ICM方法稱之為獨立連續(xù)的拓撲變量.

選定了設(shè)計變量為離散拓撲變量，連續(xù)體結(jié)構(gòu)拓撲優(yōu)化的提法可表達為

(1)

式中：ti為子域i對應(yīng)的離散拓撲設(shè)計變量，取值0或1，t為由tj組成的設(shè)計變量向量；c(t)為目標函數(shù)；Aj(t)及Bk(t)為等式及不等式約束條件；J及K分別為等式及不等式約束數(shù)目；N為設(shè)計變量數(shù)目.

這是個離散數(shù)學(xué)規(guī)劃問題，在基于有限單元法的結(jié)構(gòu)拓撲優(yōu)化中，通常取一個單元對應(yīng)于一個子域，即對應(yīng)于一個拓撲設(shè)計變量，由于有限元單元數(shù)量通常成千上萬，因此，此離散數(shù)學(xué)規(guī)劃的設(shè)計變量數(shù)通常很大，而求解大規(guī)模離散數(shù)學(xué)規(guī)劃問題遭遇極為棘手的“組合爆炸”困難. 如引言中所述，為了克服離散變量優(yōu)化問題的求解困難，通常用連續(xù)變量“代替”或“逼近”離散變量的做法，均勻化方法、變厚度法、變密度法及水平集法等采用“代替”的做法，唯有ICM法從函數(shù)逼近的視角采用“逼近”的做法.

2 離散及連續(xù)拓撲變量轉(zhuǎn)化——階躍函數(shù)及其逆函數(shù)近似逼近

2.1 離散拓撲變量與階躍函數(shù)

考察變厚度法、均勻化方法或變密度法等中的設(shè)計變量如板厚度、微結(jié)構(gòu)尺寸或人工材料密度等，這些在結(jié)構(gòu)拓撲優(yōu)化過程中變化的物理量，無論其如何小，對應(yīng)的拓撲變量為1，而當(dāng)其等于0時，拓撲變量突變?yōu)?；物理量是連續(xù)變化的，但拓撲變量卻是在零點跳躍的. ICM法引入階躍函數(shù)描述上述關(guān)系，函數(shù)圖形如圖1所示.

階躍函數(shù)為

(2)

2.2 連續(xù)拓撲變量與磨光函數(shù)

(3)

磨光函數(shù)逼近階躍函數(shù)，使拓撲變量由取值0/1的離散變量擴展為到取值在[0, 1]的連續(xù)變量. 此函數(shù)逼近過程稱為磨光近似映射. 連續(xù)拓撲變量取(0,1)值時反映了對應(yīng)子域“有”與“無”的程度.

函數(shù)圖形如圖2所示.

2.3 基于連續(xù)拓撲變量的拓撲化建模與求解

取連續(xù)拓撲變量作為設(shè)計變量，連續(xù)體結(jié)構(gòu)拓撲優(yōu)化的提法可表達為

(4)

式中ti為子域i對應(yīng)的連續(xù)拓撲設(shè)計變量，取值[0, 1]. 這是個連續(xù)變量數(shù)學(xué)規(guī)劃問題，可應(yīng)用基于連續(xù)函數(shù)的高效求解算法求解，如對偶序列二次規(guī)劃法、移動漸近線法(method of moving asymptotes,MMA)等.

2.4 離散連續(xù)拓撲變量間的映射與求解

求解基于連續(xù)拓撲設(shè)計變量建立的結(jié)構(gòu)拓撲優(yōu)化模型(如式(4)所示)，得到最優(yōu)的設(shè)計變量分布場t*，稱之為連續(xù)拓撲分布場. 最后取連續(xù)值的拓撲變量需要回歸到取0/1離散值的拓撲變量，在階躍函數(shù)的逆函數(shù)(為引用方便，以下均稱之為跨欄函數(shù))中，拓撲變量只取離散值0/1，如圖3所示，因此，要使得拓撲變量取中間值，需要一個連續(xù)可導(dǎo)函數(shù)(稱為過濾函數(shù)，如圖4所示)逼近跨欄函數(shù)

(5)

過濾函數(shù)實現(xiàn)將連續(xù)拓撲變量離散化，同時在建模時，起到了對相應(yīng)的子域物理量進行識別的作用：

(6)

過濾函數(shù)為磨光函數(shù)的逆函數(shù)，跨欄函數(shù)為階躍函數(shù)的逆函數(shù)，過濾函數(shù)近似逼近跨欄函數(shù)，而磨光函數(shù)近似逼近階躍函數(shù)，4種函數(shù)的關(guān)系如圖5所示.

過濾函數(shù)或磨光函數(shù)可以為多樣函數(shù)形式，到目前為止，ICM法研究過的過濾函數(shù)大致有3類：冪函數(shù)、修正的Sigmoid函數(shù)和指數(shù)函數(shù)[7].

綜上所述，ICM法定義了獨立拓撲變量，應(yīng)用函數(shù)逼近理論，對階躍函數(shù)及其逆函數(shù)進行連續(xù)可導(dǎo)化的函數(shù)逼近近似，從而實現(xiàn)用基于連續(xù)變量的優(yōu)化模型求解原本離散的大規(guī)劃優(yōu)化問題.

(7)

所示Heaviside函數(shù)表示這種關(guān)系，然后，用一個連續(xù)光滑的函數(shù)近似逼近此函數(shù)，如

(8)

所示.

β取不同值時的近似程度如圖6所示，可以看到，此與ICM法中磨光函數(shù)對階躍函數(shù)的近似是類似的(見圖2)，然而其應(yīng)用卻是完全不同的. 兩者之間的共同點是均用了階躍函數(shù)來表征單元的“有”或“無”，均用了連續(xù)光滑函數(shù)對階躍函數(shù)逼近. 但ICM法是用階躍函數(shù)的近似逼近來將離散拓撲變量連續(xù)化，將不可導(dǎo)函數(shù)變換為連續(xù)光滑的可導(dǎo)函數(shù)，從而用連續(xù)變量規(guī)劃模型來求解原離散規(guī)劃模型，而Heaviside投影法是應(yīng)用在拓撲邊界的尖銳化處理上，迫使過渡單元的人工密度盡量取值為0或1.

從上述對Heaviside投影法及ICM法的對比分析可以看到，發(fā)展到現(xiàn)在，人們逐漸意識到階躍函數(shù)與單元“有”或“元”表達間的關(guān)系及應(yīng)用光滑連續(xù)函數(shù)進行建模處理的好處，但還處于朦朧階段，在此發(fā)展階段，深入闡述ICM法的獨立層次拓撲變量及其基于階躍函數(shù)及其逆函數(shù)近似的建模處理方法，對加深對拓撲優(yōu)化問題本質(zhì)的認識及促進其發(fā)展是有非常重要的意義.

3 ICM法合理化建模及求解

選取連續(xù)拓撲變量作為設(shè)計變量，可以建立如式(4)形式的連續(xù)體結(jié)構(gòu)拓撲優(yōu)化模型. 這里還有2個問題：

問題1 式(4)所示模型中目標函數(shù)c(x)，等式及不等式約束條件Aj(x)及Bk(x)的選取問題.

問題2 如何將目標函數(shù)及約束條件表達為連續(xù)拓撲設(shè)計變量的近似顯式函數(shù)的問題.

對問題1，將結(jié)構(gòu)優(yōu)化問題涉及的物理量劃分為兩大類指標：經(jīng)濟指標，如結(jié)構(gòu)總重量(或總體積)、總造價等；性能指標，如結(jié)構(gòu)應(yīng)力、位移、頻率等. 則結(jié)構(gòu)優(yōu)化模型可分為兩大類：1) 在滿足性能指標下結(jié)構(gòu)經(jīng)濟指標最小化問題；2) 在結(jié)構(gòu)經(jīng)濟指標限制下的結(jié)構(gòu)性能最優(yōu)化問題.

連續(xù)體結(jié)構(gòu)拓撲優(yōu)化中研究得最多的模型是以結(jié)構(gòu)總體積(即對應(yīng)于結(jié)構(gòu)總重量)約束下結(jié)構(gòu)柔順度極小化問題. Sigmund[16]指出：應(yīng)實際工程問題的需要，可以加入其他性能約束條件，但體積約束即使不重要，也作為一個基本約束條件加入在優(yōu)化模型中以提高求解過程收斂的穩(wěn)定性.

考察低層次的截面及形狀優(yōu)化問題，可以發(fā)現(xiàn)通常建立的優(yōu)化模型是第1類以經(jīng)濟指標最小化為目標的模型，而發(fā)展到拓撲優(yōu)化層次，模型的主流卻發(fā)生了變化，轉(zhuǎn)變成了第2類模型. 拓撲優(yōu)化以體積約束作為一個基本約束，建立第2類模型， 存在如下不足：

1) 體積約束值的確定是沒有依據(jù)的、先驗的，指定不同的體積約束值會得到不同的最優(yōu)拓撲構(gòu)型.

2) 取性能為目標，在很多優(yōu)化問題中不得不建立多目標優(yōu)化模型，而多目標間的權(quán)重系數(shù)也是先驗的.

3) 與低層次的截面及形狀優(yōu)化的模型不一致，不便于建立截面、形狀及拓撲集成優(yōu)化模型.

隋允康等[26-27]研究表明：第1種模型比第2種模型更加合理. 在ICM法中，均是建立第1類優(yōu)化模型.

對問題2，在過濾函數(shù)逼近跨欄函數(shù)的過程中，建立了連續(xù)拓撲設(shè)計變量與各物理量的顯式函數(shù)關(guān)系，如式(6)所示，在此基礎(chǔ)上，通過敏度分析求得各物理量對拓撲設(shè)計變量的偏導(dǎo)數(shù)，通過一階Taylor展式，可以建立目標函數(shù)或約束與連續(xù)拓撲設(shè)計變量間的近似顯式函數(shù).

以重量極小受結(jié)構(gòu)頻率約束的拓撲優(yōu)化問題為例，其優(yōu)化模型為

(9)

上述3個過濾函數(shù)取為冪函數(shù)形式，分別為

fw(ti)=tαw,fk(ti)=tαk,fm(ti)=tαm

(10)

為避免當(dāng)拓撲變量ti→0時，質(zhì)量剛度比趨于無窮大，取αw=a,αk=αm=b.

各單元重量、單元剛度矩陣及單元質(zhì)量矩陣與拓撲設(shè)計變量間的關(guān)系分別用對應(yīng)的過濾函數(shù)進行識別：

(11)

式中上標k指第k次迭代的值.

記

令

則所有類型頻率約束可統(tǒng)一表述為

由此得到頻率約束下、重量最輕為目標的連續(xù)體結(jié)構(gòu)拓撲優(yōu)化模型為

(12)

式(12)可分離變量的非線性規(guī)劃問題，可采用對偶系列二次規(guī)劃法求解.

迭代求解的收斂準則為

ΔW=|(W(k+1)-W(k))/W(k+1)|≤ε

(13)

式中：W(k)及W(k+1)為前輪與本輪迭代的結(jié)構(gòu)總重量；ε為收斂精度，本文取ε=0.001.

其他類型性能約束的拓撲優(yōu)化問題建模及求解的過程類似[8].

4 數(shù)值算例

算例1 多點位移約束拓撲優(yōu)化

如圖7所示，基本結(jié)構(gòu)為16 mm×10 mm的平面，初始厚度為1 mm，材料彈性模量為200 GPa，泊松比為0.3. 劃分為64×40個矩形單元，載荷F=1 000 N. 左、右下腳點采用固定約束.A點處豎直向下位移約束為0.25 mm，B、C點處豎直向下位移約束為0.2 mm.

結(jié)構(gòu)初始總體積為160 mm3，A、B及C三點豎向的位移分別為-0.100 6、-0.940 2及-0.940 2 mm. 經(jīng)39次迭代后得到的拓撲圖形，如圖8所示，目標及約束的歷史迭代曲線如圖9及10所示. 最優(yōu)點總體積35.852 1 mm3，A、B及C三點最優(yōu)點位移分別為-0.249 7、-0.199 9及-0.199 9 mm.

算例2 多個頻率約束拓撲優(yōu)化

如圖11所示，60 mm×20 mm×2 mm的矩形板，2個角點固定，在下邊界布置3個集中質(zhì)量塊，1/2跨處的集中質(zhì)量塊為8 g，1/4及3/4跨處的集中質(zhì)量塊為4 g，集中質(zhì)量塊只有Y方向的慣性，彈性模量E=1×106MPa，泊松比為0.3， 材料密度ρ=1 mg/mm3. 用矩形單元，劃分為120×40=4 800個單元.

為比較動態(tài)加入頻率約束的方式解決模態(tài)交換的效果，計算2種情況：

1) 只考慮指定的頻率約束.

2) 除考慮指定的頻率約束外，動態(tài)加入防止模態(tài)交換的頻率約束.

計算情況1)不收斂，得到的部分迭代過程目標及約束的歷史曲線如圖12、13所示，計算情況2)經(jīng)過51次迭代后收斂，得到的最優(yōu)拓撲圖形如圖14，目標及約束的迭代歷史曲線如圖15、16所示，最優(yōu)點處體積為512.12 mm3， 前三階頻率為15 007.78、15 770.05、30 058.07 Hz.

從迭代歷史曲線中可以看出：在沒有加入防止模態(tài)交換的頻率約束的情況下，在迭代進行到一定程度后，第一階頻率值與第二階頻率值非常接近，從而發(fā)生模態(tài)交換，使得敏度計算得到錯誤的信息，從而引起迭代過程中出現(xiàn)振蕩，最終不收斂；加入防止模態(tài)交換的頻率約束后，迭代過程非常平穩(wěn). 由此可見，動態(tài)加入頻率約束的方式確實可以有效地防止迭代過程中模態(tài)交換現(xiàn)象的發(fā)生，從而達到穩(wěn)定收斂過程的目的.

算例3 多點頻率響應(yīng)位移幅值約束拓撲優(yōu)化

如圖17所示，基結(jié)構(gòu)尺寸為520 mm×260 mm×6 mm, 材料彈性模量為68.890 GPa，泊松比為0.3，密度為7.8×103kg/m3. 集中載荷F=7 000×sin (2πft)N分別作用于下邊界A、B及C三點，激勵頻率為500 Hz. 底部2個角點固定. 位移約束A點許用Y方向位移為0.5 mm，B及C點許用Y方向位移為0.45 mm. 計算無阻尼及結(jié)構(gòu)阻尼為0.3兩種情況.

在無阻尼情況下，41次迭代后計算得到的最優(yōu)拓撲圖形如圖18所示，有阻尼情況下，43次迭代后計算得到的最優(yōu)拓撲圖形如圖19所示. 初始基結(jié)構(gòu)的位移及最優(yōu)結(jié)構(gòu)的位移比較如表1所示，位移約束迭代歷史迭代曲線及目標迭代歷史曲線如圖20～22所示. 從最優(yōu)拓撲圖形中可以看出：在相同約束條件下，無阻尼結(jié)構(gòu)及有阻尼結(jié)構(gòu)得到的最優(yōu)拓撲圖形不同.

結(jié)構(gòu)A點B點C點基結(jié)構(gòu)(無阻尼)0.488260.444540.44454最優(yōu)結(jié)構(gòu)(無阻尼)0.499860.449150.44915基結(jié)構(gòu)(有阻尼)0.418480.381090.38108最優(yōu)結(jié)構(gòu)(有阻尼)0.499010.449570.44957

5 結(jié)論

1) 以往連續(xù)體結(jié)構(gòu)拓撲優(yōu)化沒有能夠厘清拓撲變量基本概念，而是用依附于低層次結(jié)構(gòu)優(yōu)化上，而ICM方法定義拓撲變量表征子域“有”或“無”的拓撲狀態(tài)，自然而清晰，不再需要將拓撲優(yōu)化問題降格為材料優(yōu)化、尺寸優(yōu)化或形狀優(yōu)化問題. 可惜，連續(xù)體結(jié)構(gòu)拓撲優(yōu)化研究領(lǐng)域內(nèi)還沒有對此引起足夠重視，例如變密度方法中“密度”一詞極易與真實的材料密度概念混淆.

2) 不同于采用懲罰函數(shù)進行“有”或“無”拓撲狀態(tài)的變密度方法處理，ICM方法通過對階躍函數(shù)及其逆函數(shù)的逼近，定義了磨光函數(shù)及過濾函數(shù)，識別出各物理量與連續(xù)拓撲變量間的函數(shù)關(guān)系，解決了目標函數(shù)與約束條件按拓撲設(shè)計變量近似顯式化的建模問題. 而Heaviside投影法盡管有了用階躍函數(shù)表征子域的“有”或“無”的意識，也應(yīng)用到了用連續(xù)光滑函數(shù)逼近階躍函數(shù)，但其應(yīng)用卻不是用于建立連續(xù)變量規(guī)劃模型以解決大規(guī)模離散規(guī)劃模型的求解困難，而是用在人工密度的過濾處理方面，僅起到對存在過渡單元的拓撲邊界進行銳化的作用.

3) ICM方法以經(jīng)濟指標為目標函數(shù)、性能指標為約束函數(shù)建立的拓撲優(yōu)化模型可以避免通常的以體積為約束的拓撲優(yōu)化模型的不足，與低層次的截面及形狀優(yōu)化在目標取法上保持了一致，更有利于各層次優(yōu)化間的集成.

[1] ROZVANY G I N. A critical review of established methods of structural topology optimization [J]. Structural Multidisciplinary Optimization, 2009, 37(3): 217-237.

[2] BENDSOE M P, SIGMUND O. Topology optimization: theory, methods and applications [M]. Berlin: Springer, 2003.

[3] RAO S S. Simulated annealing and parallel processing: an implementation for constrained global design optimization [J]. Engineering Optimization, 2000, 32(5): 659-685.

[4] CHEN S Y, RAJAN S D. Robust genetic algorithm for structural optimization [J]. Structural Engineering and Mechanics, 2000, 10(4): 313-336.

[5] FOURIE P C, GROENWOLD A A. Particle swarms in topology optimization [C]∥Proceedingsof Fouthworld Congress of Structural and Multidisciplinary Optimization. Dalian: WCSMO_4, 2001.

[6] HUANG X, XIE Y M. A further review of ESO type methods for topology optimization [J]. Structural and Multidisciplinary Optimization, 2010, 41: 671-683.

[7] 隋允康, 葉紅玲. 連續(xù)體結(jié)構(gòu)拓撲優(yōu)化的ICM方法[M]. 北京: 科學(xué)出版社, 2013: 50-80.

[8] BENDSOE M P, KIKUCHI N. Generating optimal topologies in structural design using a homogenization method [J]. Computer Methods in Applied Mechanics and Engineering, 1988, 69: 197-224.

[9] MLEJNEK H P. Some aspects of the genesis of structures [J]. Structural Optimization, 1992, 5: 64-69.

[10] 程耿東, 張東旭. 受應(yīng)力約束的平面彈性體的拓撲優(yōu)化[J]. 大連理工大學(xué)學(xué)報, 1995, 35(1): 1-9. CHENG G D, ZHANG D X. Topology optimization of plane elastic continuum with stress constraints [J]. Journal of Dalian University of Technology, 1995, 35(1): 1-9. (in Chinese)

[11] OSHER S, SETHIAN J. Fronts propagating with curvature dependent speed-algorithms based on Hamilton-Jacobi formulations [J]. J Comput Phys, 1988, 79(1): 12-49.

[12] BOURDIN B, CHAMBOLLE A. Design-dependent loads in topology optimization [J]. ESAIM: Control, Optimisation and Calculus of Variations, 2003, 9(8): 19-48.

[13] ESCHENAUER H A, KOBELEV V V, SCHUMACHER A. Bubble method for topology and shape optimization of structures [J]. Structural Optimization, 1994, 8(1): 42-51.

[14] NORATO J, BENDSOE O, HABER R, et al. A topological derivative method for topology optimization [J]. Structural and Multidisciplinary Optimization, 2007, 33(4/5): 375-386.

[15] GUEST J, PREVOST J, BELYTSCHKO T. Achieving minimum length scale in topology optimization using nodal design variables and projection functions [J]. Int J Numer Methods Eng, 2004, 61(2): 238-254.

[16] SIGMUND O, MAUTE K. Topology optimization approaches — a comparative review [J]. Structural Multidisciplinary Optimization, 2013, 48(6): 1031-1055.

[17] PINGEN G, WAIDMANN M, EVGRAFOV A, et al. A parametric level-set approach for topology optimization of flow domains [J]. Structural and Multidisciplinary Optimization, 2010, 41: 117-131.

[18] HUANG Y, XIE Y M. Evolutionary topology optimization of geometrically and materially nonlinear structures under prescribed design load [J]. Struct Eng Mech, 2010, 34(5): 581-595.

[19] 隋允康. 建?！ぷ儞Q·優(yōu)化——結(jié)構(gòu)綜合方法新進展[M]. 大連：大連理工大學(xué)出版社, 1996: 90-100.

[20] KAWAMOTO A, MATSUMORI T, YAMASAKI S, et al. Heaviside projection based topology optimization by a PDE-filtered scalar function [J]. Structural and Multidisciplinary Optimization, 2011, 44: 19-24.

[21] 龍凱, 左正興, 閆清東. 靜動態(tài)多目標下的連續(xù)體結(jié)構(gòu)拓撲優(yōu)化[J]. 宇航學(xué)報, 2008, 29(2): 456-460. LONG K, ZUO Z X, YAN Q D. Static and dynamic multi-objective topological optimization of continuum structure [J]. Journal of Astronautics, 2008, 29(2): 456-460. ( in Chinese)

[22] 龍凱, 左正興. 基于不同過濾函數(shù)的ICM拓撲優(yōu)化混合建模方法[J]. 中國機械工程, 2007, 18(19): 2303- 2306. LONG K, ZUO Z X. ICM method for topological Optimization of continuum structure based on different filter functions [J]. Chinese Journal of Mechanical Engineering, 2007, 18(19): 2303-2306. (in Chinese)

[23] 龍凱, 左正興. 物質(zhì)點拓撲變量法在柔性機構(gòu)設(shè)計中的應(yīng)用[J]. 計算機輔助設(shè)計與圖形學(xué)學(xué)報, 2010, 22(1): 158-164. LONG K, ZUO Z X. Topology optimization for compliant mechanism using material point topological variable[J]. Journal of Computer-Aided Design & Computer Graphics, 2010, 22(1): 158-164. (in Chinese)

[24] 鄧果, 榮見華, 邢曉娟. 一種考慮位移和應(yīng)力約束的結(jié)構(gòu)拓撲優(yōu)化方法[J]. 長沙理工大學(xué)學(xué)報, 2008, 5(3): 26-33. DENG G, RONG J H, XING X J. A structural topological optimization method considering displacement and stress constraints [J]. Journal of Changsha University of Science and Technology, 2008, 5(3): 26-33. (in Chinese)

[25] 吳志海. 基于遺傳算法的網(wǎng)殼結(jié)構(gòu)形態(tài)創(chuàng)構(gòu)方法研究[D]. 哈爾濱: 哈爾濱工業(yè)大學(xué), 2013. WU Z H. Research on form finding method of mesh shell structures based ongenetic algorithm [D]. Harbin: Harbin Institute of Technology, 2013. (in Chinese)

[26] YI G L, SUI Y K. Different effects of economic and structural performance indicators on model construction of structural topology optimization[J]. Acta Mechanica Sinia, 2015, 9: 1-12.

[27] 彭細榮, 隋允康. 對連續(xù)體結(jié)構(gòu)拓撲優(yōu)化合理模型的再探討[J]. 固體力學(xué)學(xué)報, 2016, 37(1): 1-11. PENG X R, SUI Y K. A further discussion on rational topology optimization models for continuum structures [J]. Chinese Journal of Solid Mechanics, 2016, 37(1): 1-11. (in Chinese)

(責(zé)任編輯 楊開英)

Name Correction for Design Variable of Structural Topology Optimization and Presentation of Its Corresponding Treatment Method

PENG Xirong1, SUI Yunkang2

(1.School of Civil Engineering, Hunan City University, Yiyang 413000, Hunan, China;2.Numerical Simulation Center for Engineering, Beijing University of Technology, Beijing 100124, China)

To clarify the concepts and promote developments of the structural topology optimization by analyzing the basic concepts, the necessity and urgency of name correction of the design variable of structural topology optimization and presenting its corresponding mapping method are discussed. For the topology optimization of continuum structures, the concept of topology variable was not clarified. The independent level of the topology variable was not discussed and the topology variable was attached to the variables of the structural optimization in lower levels. A penalty function method was adopted to assure the design variables of lower levels to approach to 0/1. In 2004, the Heaviside function was introduced to replace the penalty function by some researchers. Actually, the above problems and their treatments were solved completely in the independent continuous mapping (ICM) method, which was put forward in 1996. The one was the definition of the topology variable with independent level. The other was the approximations of the step function and its inverse function. In this paper, the basic concepts in the ICM method were first introduced. The process of modeling and solution was described by taking the topology optimization problem with frequency constraints as an example. Some examples show that the ICM method has obvious advantages in the aspects of modeling and solution because of its clear concepts. The Zi lu article in the Analects said that a speech will not be in right order if a name is not correct, and nothing will be successful if a speech is not in right order. It is time to correct the name of design variable of structural topology optimization and present its corresponding treatment method.

topology optimization of continuum structures; Independent topology variable; approximation functions; penalty function; step function; heaviside function; independent continuous mapping (ICM) method; variable density method

2016- 07- 24

國家自然科學(xué)基金資助項目(11672103)；湖南省自然科學(xué)基金資助項目(2016JJ6016)

彭細榮(1972—)， 男， 副教授， 主要從事結(jié)構(gòu)優(yōu)化方面的研究， E-mail：pxr568@163.com

隋允康(1943—)，男，教授，博士生導(dǎo)師，主要從事結(jié)構(gòu)優(yōu)化方面的研究，E-mail：ysui@bjut.edu.cn

O 343.1

A

0254-0037(2016)12-1787-11

10.11936/bjutxb2016070072