張東東， 欒福強， 趙禮輝， 鄭 玲

（1. 上海理工大學 機械工程學院，上海 200093；2. 上海市新能源汽車可靠性評價公共技術平臺，上海 200093；3. 重慶大學 機械與運載工程學院，重慶 400044）

引 言

振動與噪聲的控制水平是衡量航天飛行器、大型飛機、艦艇、高速列車等現代裝備發(fā)展的重要技術指標.這些裝備日益向高速、重載、輕質等方向發(fā)展，由此帶來的振動與噪聲問題更加突出.阻尼復合結構具有結構簡單、減振降噪效果良好、可靠性高等優(yōu)勢，被廣泛應用于各種大型裝備薄壁結構的振動噪聲控制，特別是對材料用量有苛刻要求的薄壁承載結構[1].

阻尼復合結構的減振性能取決于材料分布和材料屬性.借助于拓撲優(yōu)化方法，對阻尼復合材料的分布進行優(yōu)化，可在控制材料用量的約束下達到高效減振降噪目的.Ansari 等[2]采用水平集法（level set method,LSM）對薄板結構表面約束阻尼材料的位置和形狀進行尋優(yōu)設計，并通過與試驗結果對比驗證了優(yōu)化方法的有效性.Fang 等[3]以最小化諧振響應為目標函數，采用漸進優(yōu)化方法（evolutionary structural optimization, ESO）尋找板結構表面約束阻尼材料的最優(yōu)布局.Pang 等[4]以最大化加權模態(tài)損耗因子為目標對自由阻尼板結構進行了優(yōu)化設計.在上述阻尼復合結構的宏觀尺度拓撲優(yōu)化中，一般都預先確定了阻尼材料的屬性.

然而，材料屬性對阻尼復合結構的減振降噪性能也有著重要影響，比如夾層阻尼結構主要依靠阻尼材料的橫向剪切變形耗散振動能量，而且不同方向上剪切模量的差異還會引起抑振性能的變化.事實上，針對不同的抑振目標或工作環(huán)境，往往要求不同的阻尼材料屬性，如夾層阻尼板梁結構不同階抑振目標所要求的阻尼材料最優(yōu)剪切模量都不相同[5].在工程應用中，不可能按照材料屬性需求制造出所有相應的阻尼材料.此時，通過材料微觀結構的優(yōu)化設計獲得期望的材料性能是一種非常有效實用的路徑.

一方面，可以宏觀結構的性能為目標設計材料微結構.Yang 等[6]以宏觀結構的輻射聲功率最小化為目標，對阻尼材料的微結構進行了優(yōu)化設計.Chen 等[7]基于均勻化理論和實體各向同性材料懲罰（solid isotropic material with penalization, SIMP）插值方法設計了阻尼材料的微結構拓撲構型，旨在獲得約束阻尼梁結構最大模態(tài)損耗因子.Fang 等[8]以約束阻尼結構模態(tài)損耗因子最大化為目標，基于能量法預測微結構的等效性能，設計了阻尼材料三維微結構拓撲構型.另一方面，則是直接以期望性能為目標開展阻尼材料的微結構設計.Chen 等[5]采用序列線性規(guī)劃（sequential linear programming, SLP）方法研究了具有期望性能的阻尼材料微結構設計，并將其應用到約束阻尼梁結構.Andreasen 等[9]首先分析了在給定體積約束下的阻尼復合材料等效性能上下界，研究討論了考慮制造性約束時阻尼材料微結構設計能否達到等效性能邊界的問題.Huang 等[10]則以期望的材料損耗因子或彈性模量為目標，基于均勻化理論對阻尼復合材料的微結構進行了優(yōu)化設計，將優(yōu)化后阻尼材料的儲能模量和耗能模量與Hashin 和Shtrikman[11-12]提出的理論界限（H-S 界限理論）進行了對比，結果表明以最大化剛度為目標優(yōu)化設計得到的材料的儲能模量接近H-S 理論的上限，以最大化阻尼為目標優(yōu)化設計得到的材料的儲能模量接近H-S 理論的下限.上述文獻對阻尼材料微結構優(yōu)化設計研究大都面向二維微結構問題，可以滿足主要依靠拉壓變形來耗散振動能量的自由阻尼結構的多尺度設計要求.但在工程實踐中，具備更大減振降噪潛力的夾層阻尼復合結構，主要依靠阻尼材料橫向剪切變形來耗散振動能量，那么，夾層阻尼復合結構的多尺度設計就需考慮阻尼材料的橫向剪切模量.因此，將作為夾層的阻尼材料微結構優(yōu)化設計定義為三維問題則是非常有必要的.另外，上述文獻的優(yōu)化過程都預先設定了阻尼材料的體積分數，并可獲得預期材料性能，但對應的材料用量并非最少，難以保證充分利用材料.

復合材料有效模量的界限理論，如Voigt-Reuss 界限理論[13]和H-S[11-12]界限理論等，都是基于不同材料相的體積分數和本構參數構建多相復合材料等效性能極限的估計表達式，可用來指導材料微結構的設計.基于變分法的H-S 界限理論比較準確地描述了復合材料用量和等效模量界限的關系，應用較為廣泛.為了實現最小材料用量獲得期望的材料性能，本文將阻尼材料微結構等效為由實體材料和空材料構成的兩相復合材料，基于H-S 界限理論提出一種變體積約束的阻尼材料微結構拓撲優(yōu)化設計方法.基于均勻化理論，以阻尼材料期望的剪切模量為目標函數，構建了阻尼材料三維結構拓撲優(yōu)化問題.通過逆用H-S 界限理論估計對應于期望等效剪切模量的材料體積分數限，并考慮到界限理論的精度局限性，提出了一種材料體積分數限移動準則，將原優(yōu)化問題轉化為體積約束下最大化等效模量的優(yōu)化問題.采用最優(yōu)準則法進行求解得到阻尼材料三維最優(yōu)微結構，實現最小材料用量獲得期望材料性能.通過典型數值算例驗證了該方法的可行性和有效性，著重討論了初始微構型、網格密度和彈性模量等對阻尼材料微結構的影響.

1 阻尼材料均勻化

1.1 阻尼材料復模量模型

那么，阻尼材料復數模量的張量形式可以重新表示為

1.2 阻尼材料均勻化方法

典型的夾層阻尼復合結構，如圖1（a）所示，主要依靠黏彈性阻尼材料的橫向剪切變形（圖1（b））耗散能量實現減振降噪.圖1 中的uc，uv，ub分別為約束層中性面、阻尼層中性面和基板層中性面沿x方向的位移，w為夾層阻尼復合結構的橫向位移，?x為阻尼層繞y軸的轉角，?w/?x為夾層阻尼復合結構繞x軸轉動的角度.因此，對于夾層阻尼結構的阻尼材料微結構拓撲設計，需將其考慮為三維結構（圖1（c））的優(yōu)化設計問題.本文將阻尼材料單胞定義為正交各向異性材料，其本構關系為

圖1 含阻尼材料微結構的夾層阻尼結構示意圖：(a) 夾層阻尼結構；(b)阻尼材料變形圖；(c)阻尼材料三維微結構Fig. 1 Schematic drawings of a sandwich damping structure with a damping core microstructure: (a) the sandwich damping structure;(b) the deformation of the damping core; (c) the 3D microstructure of the damping core

其中，σij和εkl為黏彈性材料的應力和應變張量.

均勻化理論是基于微觀結構和材料組分來預測宏觀結構性能的一種方法，它具有嚴格的數學理論基礎，其性能預測結果更接近真實值[15].假設材料單胞尺寸遠小于材料宏觀結構尺寸，基于能量均勻化方法，阻尼材料單胞的等效復模量可以寫為[16]

式中，vi表示虛擬位移場，屬于可允許位移空間，yj表示微觀尺度坐標系下的坐標.將式（1）代入式（5），則阻尼材料等效復模量可以寫為

因此，基于單元的交互性能量，式（7）可轉化為另一種新形式[17]：

本文中，周期性邊界條件為給三維單胞邊界施加約束方程，約束方程用周期性單胞頂點、邊和面上的節(jié)點之間的位移差來表示[18].設阻尼材料單胞的微觀單元材料是各向同性的，采用SIMP 插值模型函數來表示阻尼材料微觀單元復模量的實部，則第e個微觀單元剛度矩陣的實部為

2 變體積約束的阻尼材料微結構拓撲優(yōu)化

2.1 傳統(tǒng)拓撲優(yōu)化模型

對于以期望等效材料性能為目標的阻尼材料微結構拓撲優(yōu)化問題，通?？梢员硎鰹?/p>

2.2 變體積約束拓撲優(yōu)化模型及求解

本文提出一種變體積約束的阻尼材料微結構拓撲優(yōu)化設計方法.逆用H-S 界限理論估計對應于期望等效剪切模量的材料體積分數限，并考慮到界限理論的精度局限性，提出一種材料體積分數限移動準則，將傳統(tǒng)優(yōu)化問題轉化為體積約束下最大化等效模量的優(yōu)化問題.采用優(yōu)化準則（optimal criterion，OC）法進行求解，得到阻尼材料三維最優(yōu)微結構，實現最小材料用量獲得期望材料性能.

2.2.1 目標函數

對于夾層阻尼結構，主要依靠阻尼材料的橫向剪切變形來耗散能量.為最大化利用材料且充分挖掘阻尼材料的阻尼特性，將式（14）表達的拓撲優(yōu)化問題轉化為體積約束下最大化等效模量的優(yōu)化問題.這里的體積約束限在下一小節(jié)由H-S 界限理論公式給出，是理論上獲得期望等效模量的體積約束限.本文中最大化目標函數取為阻尼材料微結構單胞三個方向剪切模量的平均值，表示如下：

2.2.2 變體積約束

1) 初始體積約束限

如引言中所述，復合材料的界限理論描述了組分材料體積分數與等效性能極限之間的近似關系，應用較廣泛的是H-S 界限理論.阻尼材料微結構可認為是由兩相材料組成的復合材料，第一相為空材料，第二相為阻尼材料.那么，在給定阻尼材料體積下，基于H-S 界限理論公式可以預測阻尼材料微結構的等效剪切模量上限Gu為[11-12]

其中，G1，K1分別為空材料的剪切模量和體積模量，G2，K2分別為阻尼材料的剪切模量和體積模量，f1，f2分別為空材料和阻尼材料的體積分數.顯然K1

本文旨在采用最小的材料用量獲得具有期望等效性能的阻尼材料微結構.逆用H-S 界限理論，獲得期望等效性能對應的阻尼材料體積限，這個體積即為理想情況下獲得期望等效材料性能的最小體積.設Gexp為已知的期望阻尼材料剪切模量，則由方程（16）可獲得Gu=Gexp時空材料和阻尼材料的體積分數f1和f2之間的關系如式（17）所示：

設空材料和阻尼材料的體積分數之和為1，令f1= 1–f2時，則可得到阻尼材料期望等效性能對應的阻尼材料體積分數限：

由于H-S 界限理論是一個理想理論，那么由式（18）的材料體積分數計算得到的最大剪切模量與期望剪切模量可能存在一定的誤差.因此，本文進一步提出變體積約束限移動準則，保證優(yōu)化后的阻尼材料微結構的等效剪切模量盡可能達到期望的剪切模量.將式（18）得到的阻尼材料初始體積分數f2記為所提出優(yōu)化問題中的阻尼材料初始體積分數約束限f(0).

2) 體積約束限的移動準則

本文中的初始體積約束限為f(0)，在迭代過程中通過改變阻尼材料的體積分數約束限使優(yōu)化后阻尼材料單胞的等效剪切模量達到期望的剪切模量.材料體積分數約束限f(k)的迭代格式可表示為

2.2.3 變體積約束拓撲優(yōu)化模型及求解

以等效彈性模量構造目標函數是材料微結構設計中廣泛采用的一種方法，尤其是以材料等效剪切模量構造目標函數，可使優(yōu)化過程具有非常高的穩(wěn)健性.基于SIMP 插值方法，以式（15）為目標函數、材料的移動體積分數限為約束進行材料微結構優(yōu)化設計，提出變體積約束的阻尼材料微結構拓撲優(yōu)化模型如下：

式中，f(k)是由式（19）定義的在第k迭代步的阻尼材料體積分數約束限.采用直接求導法得到目標函數和約束條件對變量的靈敏度后，通過優(yōu)化準則法對該優(yōu)化問題進行求解.為了保證微結構的連續(xù)性，采用靈敏度過濾技術消除棋盤格式[19].

3 數值算例與分析

筆者基于上述所提出的方法編制了變體積約束的阻尼材料微結構拓撲優(yōu)化程序.阻尼材料三維微結構的初始構型如圖2 所示，為方便描述，將其命名為初始構型A.將其離散為20 × 20 × 20 個有限元單元，每個單元的尺寸為0.05 mm × 0.05 mm × 0.05 mm，8 個角點單元為空材料（密度變量等于0），其余單元為實體阻尼材料單元（密度變量等于1）.實體阻尼材料采用復常數模量模型，彈性模量Ev為12 MPa、Poisson 比為0.48、材料損耗因子為0.5.

圖2 阻尼材料單胞初始構型AFig. 2 Initial configuration A of the damping material unit cell

考慮阻尼材料微結構單胞的期望剪切模量Gexp分別為1 MPa，2 MPa 和3 MPa，那么由式（18）可得對應于3 個期望剪切模量的阻尼材料體積分數約束值分別為0.360，0.624 和0.830.此時為了獲得3 種期望的剪切模量值，阻尼材料微結構的拓撲優(yōu)化設計問題就轉化為理想情況下在體積約束分別為0.360，0.624 和0.830 下的剪切模量最大化問題，目標函數如式（15）所示.在優(yōu)化過程中，為了與求解最小化問題的優(yōu)化算法相適應，將目標函數修改為最小化問題.

對于上述3 個不同的期望剪切模量，優(yōu)化后的阻尼材料微結構拓撲構型、迭代歷程以及等效剪切模量如圖3 所示.優(yōu)化后的阻尼材料微結構仍呈軸對稱特征，迭代過程中阻尼材料的體積分數限呈現出移動變化的特征.對應的等效剪切模量分別為1.02 MPa，2.04 MPa 和3.04 MPa，非常接近期望剪切模量且誤差都能控制在2%以內.優(yōu)化后阻尼材料的體積分數分別為0.395，0.614 和0.830，與逆用H-S 界限理論得到的阻尼材料體積分數限相比，分別增加了9.73%，?1.61%，0%，與初始的材料體積分數限有一定的差異.這是由于H-S 界限理論本身是一種理想的近似理論，那么逆用H-S 界限理論獲得的對應于期望等效剪切模量的材料體積分數限也是一個理想的估計值，必然存在一定的精度局限.因此，本文構造了材料體積分數限移動準則，在優(yōu)化過程中通過增大或者減小材料的體積分數約束限來獲得最小材料用量下具有期望剪切模量的阻尼材料微結構構型，實現材料的充分利用.

圖3 基于初始構型A 獲得阻尼材料單胞微結構拓撲優(yōu)化結果Fig. 3 Topology optimization results of the damping material unit cell microstructure for initial configuration A

3.1 優(yōu)化結果對初始構型的依賴性分析

在阻尼材料微結構單胞尺寸和網格密度不變的情況下，改變阻尼材料單胞的初始構型，分析它們對優(yōu)化結果的影響.定義了兩種不同的初始構型B 和C，如圖4 所示.初始構型B 的中心8 個單元為空材料，其余單元為實體阻尼材料；初始構型C 的中心8 個單元為實體阻尼材料，其余單元為密度等于初始體積分數限值的阻尼材料.

圖4 阻尼材料初始構型：(a) 阻尼材料初始構型B； (b) 阻尼材料初始構型CFig. 4 The initial configurations of the damping material unit cell：(a) initial configuration B of the damping material unit cell;(b) initial configuration C of the damping material unit cell

基于初始構型B 的阻尼材料微結構最優(yōu)拓撲構型及等效剪切模量如圖5 所示.優(yōu)化后的阻尼材料微結構仍呈軸對稱特征，對應的等效剪切模量分別為1.03 MPa，2.04 MPa 和3.04 MPa，可以認為達到了期望的剪切模量.優(yōu)化后阻尼材料的體積分數分別為0.385，0.620 和0.828，與逆用H-S 界限理論得到的阻尼材料體積分數限相比，分別增加了6.95%，?0.60%，?0.24%.

圖5 基于初始構型B 獲得阻尼材料單胞微結構拓撲優(yōu)化結果Fig. 5 Topology optimization results of the damping material unit cell microstructure for initial configuration B

基于初始構型C 的阻尼材料微結構最優(yōu)拓撲構型及剪切模量如圖6 所示.優(yōu)化規(guī)整后的等效剪切模量分別為1.03 MPa，2.05 MPa 和3.05 MPa，與期望的剪切模量相比，誤差控制在3%以內，仍可認為達到了期望的剪切模量.優(yōu)化后阻尼材料的體積分數分別為0.388，0.604 和0.820，分別增加了7.78%，?3.20%，?1.2%.

圖6 基于初始構型C 獲得阻尼材料單胞微結構拓撲優(yōu)化結果Fig. 6 Topology optimization results of the damping material unit cell microstructure for initial configuration C

分析上述優(yōu)化結果可以看出：在不同的初始構型下，都可獲得具有期望剪切模量的阻尼材料微結構，且材料用量（微結構體積分數）差別不大，但優(yōu)化構型存在差異.本文采用相似度描述優(yōu)化構型的差異性，定義為：不同初始構型的優(yōu)化構型與參考構型的單元重合數量之比[20].以初始構型A 對應的優(yōu)化構型（如圖3 所示）為

參考構型，可獲得初始構型B 下3 個優(yōu)化構型的相似度分別為0.162，0.475 和0.794；初始構型C 下3 個優(yōu)化構型的相似度分別為0.830，0.917 和0.968.因此，可以認為以初始構型B 獲得的優(yōu)化構型與初始構型A 差異大；以初始構型C 獲得的優(yōu)化構型與初始構型A 相似度高.

3.2 優(yōu)化結果對網格密度的依賴性

討論了阻尼材料單胞的網格密度對優(yōu)化結果的影響.以初始構型A 為例，修改阻尼材料單胞網格密度分別為15 × 15 × 15，25 × 25 × 25 和30 × 30 × 30，對應的優(yōu)化結果如圖7～9 所示.

圖7 單胞初始構型A 優(yōu)化結果(15 × 15 × 15)Fig. 7 Topology optimization results of initial configuration A (15 × 15 × 15)

通過對比分析可以發(fā)現：①不同的網格密度下獲得的阻尼材料單胞等效剪切模量與期望剪切模量的誤差都在3%以內；②與逆用H-S 界限理論得到的阻尼材料理想體積分數相比，在期望剪切模量較小時（1 MPa），不同的網格密度下阻尼材料體積分數略微增加，其余兩種期望剪切模量下的阻尼材料體積分數則變化很小；③不同的網格密度下，相對于阻尼材料的微結構最優(yōu)構型具有相似的的材料分布特征，且隨著網格數量的增加，微結構最優(yōu)構型越光滑.

3.3 優(yōu)化結果對材料彈性模量的依賴性分析

基于初始構型A，阻尼材料微結構的網格密度為20 × 20 × 20，每個微觀單元的尺寸為0.05 mm × 0.05 mm ×0.05 mm .阻尼材料的Poisson 比和損耗因子分別為0.48 和0.5，阻尼材料的彈性模量分別選取10 MPa 和15 MPa.

阻尼材料的彈性模量為10 MPa 時，由式（18）可得對應于3 個期望剪切模量的阻尼材料體積分數約束值分別為0.417，0.712 和0.931.優(yōu)化后的阻尼材料微結構構型及等效剪切模量如圖10 所示.對應的等效剪切模量分別為1.01 MPa，2.03 MPa 和3.01 MPa，可以認為達到了期望的剪切模量；優(yōu)化后阻尼材料的體積分數分別為0.436，0.701 和0.929，與逆用H-S 界限理論得到的阻尼材料體積分數限相比，分別增加了4.36%，?1.57%，?0.22%.

圖9 單胞初始構型A 優(yōu)化結果(30 × 30 × 30)Fig. 9 Topology optimization results of initial configuration A (30 × 30 × 30)

圖10 基于初始構型A 的拓撲優(yōu)化結果(Ev=10 MPa)Fig. 10 Topology optimization results of initial configuration A (Ev=10 MPa)

阻尼材料的彈性模量為15 MPa 時，對應于3 個期望剪切模量的阻尼材料體積分數約束值分別為0.295，0.526 和0.712.優(yōu)化后的阻尼材料微結構構型及等效剪切模量如圖11 所示.對應的等效剪切模量分別為1.02 MPa，2.04 MPa 和3.03 MPa，與期望的剪切模量相比，誤差控制在2%以內，仍可認為達到了期望的剪切模量；優(yōu)化后阻尼材料的體積分數分別為0.333，0.531 和0.699，與逆用H-S 界限理論得到的阻尼材料體積分數限相比，分別增加了11.41%，?0.94%，?1.86%.

對比不同彈性模量下阻尼材料的優(yōu)化結果，可以看出：阻尼材料的彈性模量變化時，仍可獲得具有期望剪切模量的微結構；當期望的剪切模量相同時，阻尼材料彈性模量越大，優(yōu)化后的阻尼材料微結構體積分數越小.

4 結 論

本文以最小材料用量獲得期望阻尼材料特性為目標，提出了基于變體積約束的阻尼材料微結構拓撲優(yōu)化設計方法.以獲得阻尼材料的期望剪切模量為目標，逆用H-S 界限理論估計對應于期望等效剪切模量的材料體積分數限，并提出體積分數約束限的移動準則，將阻尼材料微結構拓撲優(yōu)化問題轉化為材料體積約束下的阻尼材料剪切模量最大化拓撲優(yōu)化問題.通過數值算例，得出了以下結論：

1) 該方法可以很好地實現在最小材料用量下獲得具有期望性能的阻尼材料微結構；

2) 優(yōu)化后阻尼材料微結構的拓撲構型受初始構型的影響較大，但基于不同初始構型優(yōu)化后的阻尼材料微結構體積分數差別不大；

3) 優(yōu)化后阻尼材料微結構的拓撲構型和體積分數對阻尼材料單胞網格密度的依賴性很小，且優(yōu)化后的拓撲構型隨著網格密度的增加而變得光滑；

4）優(yōu)化后阻尼材料微結構的拓撲構型對彈性模量的依賴性較小，阻尼材料彈性模量越大，優(yōu)化后的阻尼材料微結構體積分數越小.

參考文獻（ References ） ：

［1］曾 昭陽, 范紅偉, 焦映厚, 等. 基于三明治夾層約束阻尼結構的潛艇降噪[J]. 科學技術與工程, 2020, 20（22）: 8975-8982. （ZENG Zhaoyang, FANG Hongwei, JIAO Yinghou, et al. Noise reduction of submarine based on sandwich constrained layer damping structure[J].Science Technology and Engineering, 2020, 20（22）: 8975-8982.（in Chinese））

［2］A NSARI M, KHAJEPOUR A, ESMAILZADEH E. Application of level set method to optimal vibration control of plate structures[J].Journal of Sound and Vibration, 2013, 332（4）: 687-700.

［3］F ANG Z P, ZHENG L. Topology optimization for minimizing the resonant response of plates with constrained layer damping treatment[J].Shock and Vibration, 2015, 2015（2）: 1-11.

［4］P ANG J, ZHENG W G, YANG L, et al. Topology optimization of free damping treatments on plates using level set method[J].Shock and Vibration, 2020, 2020: 5084167.

［5］C HEN W J, LIU S T. Topology optimization of microstructures of viscoelastic damping materials for a prescribed shear modulus[J].Structural and Multidisciplinary Optimization, 2014, 50（2）: 287-296.

［6］Y ANG R Z, DU J B. Microstructural topology optimization with respect to sound power radiation[J].Structural and Multidisciplinary Optimization, 2013, 47（2）: 191-206.

［7］C HEN W J, LIU S T. Microstructural topology optimization of viscoelastic materials for maximum modal loss factor of macrostructures[J].Structural and Multidisciplinary Optimization, 2016, 53（1）: 1-14.

［8］F ANG Z P, LEI Y, TIAN S X, et al. Microstructural topology optimization of constrained layer damping on plates for maximum modal loss factor of macrostructures[J].Shock and Vibration, 2020, 2020: 8837610.

［9］A NDREASEN C S, ANDREASSEN E, JENSEN J S, et al. On the realization of the bulk modulus bounds for two-phase viscoelastic composites[J].Journal of the Mechanics and Physics of Solids, 2014, 63（1）: 228-241.

［10］H UANG X D, ZHOU S W, SUN G Y, et al. Topology optimization for microstructures of viscoelastic composite materials[J].Computer Methods in Applied Mechanics and Engineering, 2015, 283: 503-516.

［11］H ASHIN Z, SHTRIKMAN S. A variational approach to the theory of the elastic behaviour of multiphase materials[J].Journal of the Mechanics and Physics of Solids, 1963, 11（2）: 127-140.

［12］H ASHIN Z. Complex moduli of viscoelastic composites, Ⅰ: general theory and application to particulate composites[J].International Journal of Solids and Structures, 1970, 6（5）: 539-552.

［13］楊 大鵬, 劉新田. 復合材料有效彈性模量的上、下限的求解[J]. 鄭州大學學報(工學版), 2002, 23（2）: 106-109.（YANG Dapeng, LIU Xintian. Solution to the upper and lower bounds of the effective elastic modulus of composite materials[J].Journal of Zhengzhou University(Engineering Edition), 2002, 23（2）: 106-109.（in Chinese））

［14］L IU Q M, RUAN D, HUANG X. Topology optimization of viscoelastic materials on damping and frequency of macrostructures[J].Computer Methods in Applied Mechanics and Engineering, 2018, 337: 305-323.

［15］程 可朋, 王憲杰, 張洵安, 等. 周期性復合材料構型及結構一體化優(yōu)化[J]. 應用數學和力學, 2015, 36（7）: 725-732.（CHENG Kepeng, WANG Xianjie, ZHANG Xun’an, et al. Collaborative optimization of structures with periodic composite materials[J].Applied Mathematics and Mechanics, 2015, 36（7）: 725-732.（in Chinese））

［16］S IGMUND O. Materials with prescribed constitutive parameters: an inverse homogenization problem[J].International Journal of Solids and Structures, 1994, 31（17）: 2313-2329.

［17］X IA L, BREITKOPF P. Design of materials using topology optimization and energy-based homogenization approach in Matlab[J].Structural and Multidisciplinary Optimization, 2015, 52（6）: 1229-1241.

［18］G AO J, LI H, GAO L, et al. Topological shape optimization of 3D micro-structured materials using energy-based homogenization method[J].Advances in Engineering Software, 2018, 116: 89-102.

［19］S IGMUND O. Morphology-based black and white filters for topology optimization[J].Structural and Multidisciplinary Optimization, 2007, 33（4）: 401-424.

［20］Z HANG D D, WU Y H, LU X, et al. Topology optimization of constrained layer damping plates with frequencyand temperature-dependent viscoelastic core via parametric level set method[J].Mechanics of Advanced Materials and Structures, 2021, 29（1）: 154-170.