Topology optimization for 3D metamaterials with negative Poisson's ratio
HU Tiannan1, GUO Honghu2 (1. DepartmentofechanicalEngneeringandience,Kyoto UnversityKyoto6553o,Jan;2.Facultyofied Engineering, Waseda University, Kyoto l69-8555, Japan)
Abstract: Negative Poisson's ratio metamaterials have significant application potential in aerospace, biomedical engineering, and flexible electronics due to their anomalous mechanical behavior of expanding laterally under tension and contracting laterally under compression. However, existing research predominantly focuses on 2D or 3D isotropic designs, which struggle to meet practical engineering demands for direction-dependent material properties. This paper proposed a density-based topology optimization method for designing 3D orthotropic negative Poisson's ratio metamaterials. By constructing a novel multi-objective optimization function combined with homogenization theory, negative Poisson's ratio characteristics in three orthogonal directions for unit cell structures was achieved. First, based on the SIMP (solid isotropic material with penalization) material interpolation model, geometric constraints were introduced to ensure unit cell symmetry while eectively reducing computational scale. Second, a new objective function and optimization model were established through penalty functions. Finally, the equivalent mechanical properties were calculated using finite element homogenization under periodic boundary conditions, with design variables updated through sensitivity analysis. Numerical examples demonstrate that the optimized unit cells exhibit negative Poisson's ratio behavior in all three principal directions. This study provides a theoretical foundation for the controllable fabrication of 3D orthotropic auxetic metamaterials and expands the design and application scope of mechanical metamaterials.
Keywords: topology optimization; negative Poisson's ratio; metamaterial; homogenization method; density method
負泊松比材料作為一種力學超材料,因其在單軸載荷下呈現(xiàn)反常的橫向變形行為(拉伸時橫向膨脹,壓縮時橫向收縮),近年來在航空航天、生物醫(yī)學、防護裝備及智能傳感器等領(lǐng)域展現(xiàn)出廣闊的應用前景[1]。與傳統(tǒng)材料相比,負泊松比材料在抗剪切性、抗斷裂性和能量吸收率等方面更有優(yōu)勢[2]。然而,現(xiàn)有研究多集中于二維結(jié)構(gòu)或基于經(jīng)驗設(shè)計的周期性胞元,三維負泊松比超材料的設(shè)計仍面臨幾何自由度受限、力學性能可調(diào)性不足及制造成本高昂等挑戰(zhàn)[3]。針對這些問題,提出結(jié)合拓撲優(yōu)化方法與均勻化理論,系統(tǒng)研究三維負泊松比超材料的優(yōu)化設(shè)計方法,旨在突破傳統(tǒng)設(shè)計框架,實現(xiàn)高性能、可定制化的負泊松比超材料結(jié)構(gòu)設(shè)計。
拓撲優(yōu)化作為一種基于數(shù)學模型的材料分布優(yōu)化方法,能夠通過目標函數(shù)與約束條件的靈活配置,高效探索材料與結(jié)構(gòu)的最優(yōu)構(gòu)型[4]。均勻化法通過表征微觀結(jié)構(gòu)與宏觀等效性能的關(guān)系,為超材料的多尺度設(shè)計提供了理論框架。Sigmund 等[5]系統(tǒng)總結(jié)了拓撲優(yōu)化在多物理場超材料設(shè)計中的應用,指出其在實現(xiàn)負泊松比、負熱膨脹等多功能耦合方面的潛力。Anaya-Jaimes等提出了基于雙向漸進結(jié)構(gòu)優(yōu)化的無約束設(shè)計框架,結(jié)合鄰近材料相間插值、靈敏度穩(wěn)定化策略及動態(tài)添加比率調(diào)控,實現(xiàn)了二維正交各向異性超材料的熱膨脹系數(shù)精準定制。Murai等運用拓撲優(yōu)化技術(shù)實現(xiàn)了可捕捉局部共振效應,并精準調(diào)控多頻電磁波行為的超材料設(shè)計。Paulino等構(gòu)建了針對功能梯度材料的多尺度架構(gòu)設(shè)計框架,并在一維、二維梯度及極端性能(近零剪切模量、負泊松比)材料中驗證了方法的普適性與有效性。Ha等[9]提出了一種基于手性負泊松比晶格(內(nèi)凹六邊形或手性三角形結(jié)構(gòu))的材料設(shè)計方法,通過采用兩種具有不同熱膨脹系數(shù)的材料,構(gòu)建雙材料梁單元,在宏觀尺度上實現(xiàn)了對材料整體熱膨脹系數(shù)的主動調(diào)控。Ai等[1基于參數(shù)化水平集方法和緊支徑向基函數(shù)無網(wǎng)格方法,提出了一種二維微結(jié)構(gòu)機械超材料拓撲優(yōu)化設(shè)計機制,實現(xiàn)了負泊松比拉脹超材料設(shè)計,并可擴展至設(shè)計具有負熱膨脹系數(shù)或頻率帶隙的超材料。Qin 等[1]以最大化初始機械阻抗水平為目標,提出了一種用于提升減振性能的超材料拓撲優(yōu)化設(shè)計方法。Ye等[12]設(shè)計了一種自適應漸進變硬機械超材料,克服了傳統(tǒng)均勻負泊松比材料常存在的剛度突變或應力局限性,為開發(fā)具有定制化非線性力學響應的超材料提供了設(shè)計理論與方法。Jebellat等[13]研究提出了一種由智能材料構(gòu)成的超結(jié)構(gòu)設(shè)計,該結(jié)構(gòu)不僅具有形狀記憶效應,還具備負泊松比性能,可應用于新型傳感器、執(zhí)行器及生物醫(yī)學領(lǐng)域。Agrawal等[14]基于密度法,提出了一種魯棒拓撲優(yōu)化(robust topologyoptimization,RTO)設(shè)計方法,以應對材料不確定性對負泊松比超材料性能的影響。
雖然運用拓撲優(yōu)化實現(xiàn)負泊松比結(jié)構(gòu)設(shè)計已十分成熟,三維負泊松比超材料的優(yōu)化設(shè)計研究卻相對較少,限制了超材料的工程應用。Evans等[15]提出了三維正交內(nèi)凹蜂窩負泊松比結(jié)構(gòu)。基于拉一扭耦合效應,Duan 等[提出了一種新型三維負泊松比材料設(shè)計方法,通過將多層二維負泊松比材料與三維手性鏈結(jié)構(gòu)耦合連接,構(gòu)建出兼具負泊松比和正泊松比性能的三維材料體系。然而,這些三維超材料設(shè)計都并非基于拓撲優(yōu)化設(shè)計,難以通過合理的材料分布實現(xiàn)最優(yōu)的機械性能。Clausen等[17]運用幾何非線性拓撲優(yōu)化技術(shù),首次實現(xiàn)了泊松比在 -0.8~0.8 范圍內(nèi)可編程調(diào)控的三維超材料設(shè)計,并通過參數(shù)化超級橢圓結(jié)構(gòu)優(yōu)化與多噴嘴連續(xù)打印策略,解決了傳統(tǒng)拓撲優(yōu)化超材料結(jié)構(gòu)難以規(guī)?;圃斓碾y題。Zhang等[18]利用GPU并行計算優(yōu)化設(shè)計了高分辨率的三維負泊松比微結(jié)構(gòu)。Vogiatzis等[9]采用調(diào)和水平集方法實現(xiàn)了多材料負泊松比超材料的優(yōu)化設(shè)計,顯著提升了結(jié)構(gòu)的剛度與能量吸收性能。雖然基于拓撲優(yōu)化的三維負泊松比超材料設(shè)計已逐漸獲得關(guān)注,如何便捷有效地實現(xiàn)三維正交異性負泊松比性能仍然需要進一步研究。
針對上述挑戰(zhàn),本文提出一種基于拓撲優(yōu)化與均勻化理論的三維負泊松比超材料設(shè)計方法。首先,通過密度法構(gòu)建參數(shù)化設(shè)計空間,結(jié)合均勻化理論計算等效力學性能;其次,引入調(diào)節(jié)因子調(diào)控結(jié)構(gòu)的負泊松比行為,實現(xiàn)拉壓載荷下力學響應的對稱性優(yōu)化;最后,通過光固化(stereolithography,SLA)技術(shù)驗證設(shè)計方案的可行性。本文的創(chuàng)新點包括:a.提出了一種兼顧幾何自由度與計算效率的三維拓撲優(yōu)化框架;b.構(gòu)建了新的三維負泊松比拓撲優(yōu)化目標函數(shù);c.通過增材制造驗證了優(yōu)化結(jié)構(gòu)的可制造性。研究結(jié)果將為高性能負泊松比超材料的定制化設(shè)計提供理論支撐與技術(shù)路徑。
1 三維負泊松比超材料拓撲優(yōu)化設(shè)計方法
1.1 均勻化法
利用均勻化法計算超材料單胞的等效彈性矩陣,通過施加周期性邊界條件對微結(jié)構(gòu)單胞進行有限元分析,結(jié)構(gòu)平衡方程為
KU=f
式中: K 為單胞的整體剛度矩陣; U 為單胞的位移場; f 為由均勻應變場引起的外力,可以表示為
式中: B 為應變-位移矩陣; Y 為單胞; Ye 為單胞中第 e 個單元的體積; De 為單元 e 的彈性矩陣;m 為單元總個數(shù)。針對三維問題,應變 ε=[ε1,ε2,ε3 ε4, ε5 , ε6] , ε1=[1,0,0,0,0,0]T , ε2=[0,1,0,0,0,0]T ,ε3=[0,0,1,0,0,0]T , ε4=[0,0,0,1,0,0]T , ε5=[0,0,0,0,1,0]T ,ε6=[0,0,0,0,0,1]T 。
得到位移場后,利用均勻化法,超材料的等效彈性矩陣 DH 可以表示為
式中: |Y| 為單胞的總體積; ue0 為通過施加單位應變計算得到的節(jié)點位移場; Ue 為通過施加全局單位應變計算得到的位移場[20]; ke 為單元剛度矩陣。
keue0=fe
矩陣 DH 中的彈性張量 Drs 可以表示為
式中: r 對應應力分量, r=1,2,3,4,5,6 s 對應應變分量, s=1,2,3,4,5,6
DH=
式中: c 為彈性系數(shù)。對于施加對稱性約束后的三維正交各向異性材料而言, C1122H=C2211H ,C1133H=C3311H , C2233H=C3322H°
通過對彈性矩陣中相關(guān)項的計算,可得到3個方向的泊松比的計算公式為
1.2 優(yōu)化數(shù)學模型
實現(xiàn)正交各向異性負泊松比超材料設(shè)計的關(guān)鍵是彈性矩陣中 為負數(shù),理論上,通過最小化 C1122H+C1133H+C2233H 可實現(xiàn)負泊松比設(shè)計。因此,負泊松比超材料設(shè)計的核心轉(zhuǎn)變?yōu)槿绾瓮ㄟ^迭代,實現(xiàn) C1122 ! C1133 和 C2233 的最小化設(shè)計。為了改善迭代過程的收斂性,基于文獻[14]和[21],引人松弛項 wβl(C1111H+C2222H+C3333H) ,構(gòu)建新的目標函數(shù)。最終,通過引入松弛函數(shù)構(gòu)造的新目標函數(shù)為
式中: β 為常數(shù),取值范圍為(0,1);指數(shù) l 為迭代次數(shù); w 為調(diào)節(jié)參數(shù)。
為了調(diào)節(jié)3個正交方向上的泊松比數(shù)值,引入額外的3個調(diào)節(jié)參數(shù) w1,w2 及 w3 ,將式(8)修改為
將固體材料體積設(shè)定為上限值,其作為不等式約束,具體上限值可根據(jù)實際應用需求進行靈活調(diào)整。因此,最終的優(yōu)化數(shù)學模型為
式中: u 是實體材料的體積; uf 是體積分數(shù); u0 是微結(jié)構(gòu)的總體積; x 是設(shè)計變量; xmin 是為了避免數(shù)值不穩(wěn)定現(xiàn)象設(shè)置的1個小的正值。體積分數(shù)約束定義了微觀結(jié)構(gòu)中材料的最大體積?;赟IMP法,單元剛度矩陣 ke 可以表示為
ke=xep(ks-k0)+k0
式中: ks 為實體單元的剛度矩陣; k0 為空材料單元對應的剛度矩陣; p 為懲罰因子。
此外,為了保證最終設(shè)計結(jié)果的正交各向異性,在 xy 、yz、 zx 面對單胞結(jié)構(gòu)施加幾何對稱約束,如圖1所示。由于結(jié)構(gòu)是對稱的,因此,拓撲優(yōu)化對應的獨立設(shè)計變量為圖1藍色部分所示的單胞微結(jié)構(gòu)1/8部分。對稱約束的施加,不僅保證了設(shè)計結(jié)果的正交性,也減小了設(shè)計變量,提高了計算效率。
1.3 并行計算
三維結(jié)構(gòu)拓撲優(yōu)化相較于二維結(jié)構(gòu),其計算規(guī)模呈指數(shù)增加。為了實現(xiàn)三維負泊松比超材料的高分辨率設(shè)計,采用開源并行框架對上述問題進行求解。使用開源軟件FreeFEM ++ 對均勻化問題進行離散,采用PETSc并行求解線性方程組。并行求解算法邏輯如下:首先,使用圖分割包,如METIS,對網(wǎng)格進行分解,并指派給CPU各核心;隨后,由FreeFEM執(zhí)行偏微分方程離散,并將離散后的線性代數(shù)方程組傳遞給PETSc進行并行求解;最后,將分別求解得到的數(shù)據(jù)進行回傳處理,得到最終分析結(jié)果。本文所有數(shù)值算例均在DELL桌面級工作站Precision3660上完成,其主要配置如下:中央處理器為Inteli913900K,基頻為 3.0GHz ;內(nèi)存為128GB;顯卡為NVIDIAGeForceRTX4080。
1.4 優(yōu)化算法及流程
優(yōu)化算法的流程如圖2所示。其中,有限元離散和線性方程求解通過并行計算處理,兩者均基于MPI(message passing interface)實現(xiàn)。
優(yōu)化算法的描述如下:
a.建立設(shè)計域和參數(shù)初始化
b.對設(shè)計域進行離散,并網(wǎng)格進行分解。
c.有限元法求解控制方程。
d.計算目標函數(shù),進行收斂判斷。如果目標
函數(shù)收斂,則優(yōu)化結(jié)束,否則跳轉(zhuǎn)至步驟e進行迭代直至收斂。收斂條件:優(yōu)化迭代的目標函數(shù)連續(xù)兩次的差值小于收斂容差,或者迭代次數(shù)到達設(shè)定的最大值 Nc
e.計算靈敏度。
f.基于移動漸近線算法(methodofmovingasymptotes,MMA)更新設(shè)計變量,并跳轉(zhuǎn)至步驟c。
2 數(shù)值算例
2.1 算例設(shè)置與優(yōu)化過程
為了驗證所提設(shè)計方法的有效性,選取邊長為1的三維單胞微結(jié)構(gòu)作為設(shè)計域,采用結(jié)構(gòu)化的四面體網(wǎng)格進行離散,單邊網(wǎng)格數(shù)為30,網(wǎng)格總數(shù)為 32.4×104 。為了實現(xiàn)結(jié)果的對稱性,在單胞構(gòu)型上施加三向?qū)ΨQ面 (xy,yz,zx) 約束并構(gòu)建網(wǎng)格,如圖3所示。優(yōu)化數(shù)學模型中的參數(shù)設(shè)定為:松弛因子 β=0.8 ,調(diào)節(jié)參數(shù) w=1000 , w1= 3, w2=1 , w3=4 ,體積分數(shù) uf=0.5 。材料屬性為:彈性模量 E0=0.91 ,泊松比 μ=0.3 。
迭代歷程及優(yōu)化結(jié)果如圖4所示。在初始階段,迭代步 llt;50 時,設(shè)計域內(nèi)材料呈均勻分布,3 個正交方向的泊松比 μxy,μyz 和 μzx 相等,為材料本身的泊松比0.3。隨著迭代推進,迭代步 l 為50~200 時,目標函數(shù)通過懲罰項動態(tài)調(diào)整彈性矩陣項 C1111H+C2222H+C3333H ,驅(qū)動泊松比向負值演變。
如圖5所示,至迭代步 l=200 時, μzx 率先下降到 -0.5 ,表明 zx 方向上率先形成穩(wěn)定的拉脹結(jié)構(gòu)。最終,迭代步 lgt;300 時,目標函數(shù)趨于穩(wěn)定,3個方向上的泊松比分別達到 μxy=-0.558 ,μyz=-0.161,μzx=-0.801,實現(xiàn)三維正交各向異性負泊松比超材料結(jié)構(gòu)逆向設(shè)計??梢钥吹剑罱K獲得的3個方向的負泊松比數(shù)值間的比值關(guān)系基本滿足設(shè)計要求 w1=3 , w2=1 , w3=4 。
2.2 優(yōu)化構(gòu)型與力學機理
圖6為優(yōu)化后的單胞結(jié)構(gòu)及其與經(jīng)典二維手性超材料的對比。圖 6(a)~(b) 所示的優(yōu)化結(jié)果呈現(xiàn)出獨特的空間螺旋拓撲結(jié)構(gòu): z 軸方向通過交叉肋條形成類似二維手性結(jié)構(gòu)的扭轉(zhuǎn)耦合效應,而x、 y 軸方向通過弧形空腔實現(xiàn)橫向膨脹抑制。這種三維耦合機制使得單軸拉伸時,交叉肋條的彈性變形引發(fā)多向收縮響應。值得注意的是,優(yōu)化結(jié)構(gòu)的材料分布與二維手性超材料具有幾何相似性,如圖6(c所示,驗證了目標函數(shù)對負泊松性能的捕捉能力。
由于優(yōu)化是基于密度法進行的,從圖6(b)可以看到,優(yōu)化后的單胞有較為明顯的鋸齒形邊界面,這一現(xiàn)象可以通過引入額外的技術(shù)來消除,如貼體網(wǎng)格技術(shù)。然而,考慮到單胞的幾何尺寸,圖6(b)所示的鋸齒形邊界面在實際制造中并不會變得明顯,且對結(jié)構(gòu)的性能影響較小,因此,在本文中不作額外的處理
進一步將優(yōu)化單胞進行 3×3 周期擴展,如圖7所示,宏觀等效泊松比呈現(xiàn)顯著方向依賴性:沿主對角線方向( μzx=-0.801 )的拉脹效應最強,而yz方向 μyz=-0.161 )受對稱約束影響較弱。這種各向異性可通過調(diào)整體積分數(shù)和幾何約束進行定制,為多軸載荷場景提供設(shè)計靈活性。
3 超材料增材制造
為了驗證三維負泊松比超材料的可制造性,使用3D技術(shù)對設(shè)計結(jié)果進行增材制造。首先,將文中的優(yōu)化結(jié)果導出為STL文件;然后,將STL文件導入切片軟件以生成G代碼;最后,將G代碼導入3D打印機打印結(jié)構(gòu)。3D打印機的打印最高分辨率為 16μm ,打印精度為: 20~85μm (尺寸小于 50mm 零件); 200μm (全尺寸模型, 260mm 零件),成型材料為光敏樹脂。3D打印的結(jié)果如圖8所示。從圖中可以看出,打印的結(jié)構(gòu)具有很高的分辨率,與優(yōu)化結(jié)果基本一致,驗證了優(yōu)化結(jié)果的可制造性。
4結(jié)論
本文提出了一種基于密度法的三維正交各向異性負泊松比超材料并行拓撲優(yōu)化框架,通過理論推導、數(shù)值仿真與增材制造驗證,得出以下結(jié)論:a.通過融合幾何對稱約束與動態(tài)懲罰函數(shù),成功實現(xiàn)了三維單胞在3個正交方向的負泊松比各向異性調(diào)控 'μxy=-0.558 , μyz=-0.161 , μzx=-0.801 。優(yōu)化結(jié)構(gòu)呈現(xiàn)空間螺旋拓撲,其力學行為與二維手性超材料具有機理一致性,驗證了目標函數(shù)的物理合理性。
b.采用PETSc并行求解器與MPI任務分配策略,有效縮短了計算時間,相較傳統(tǒng)串行算法效率提升8倍,為高分辨率設(shè)計提供技術(shù)支撐。
c.增材制造樣件證實了優(yōu)化結(jié)構(gòu)的可制造性,其多向拉脹特性可應用于定制化緩沖材料、柔性傳感器及航空航天輕量化結(jié)構(gòu)。通過調(diào)節(jié)體積分數(shù)與對稱約束,可進一步擴展至負熱膨脹、能量吸收等多功能耦合設(shè)計。
未來工作將聚焦于非線性大變形超材料拓撲優(yōu)化算法開發(fā),并探索金屬基超材料的SLA成形工藝,以提升結(jié)構(gòu)的承載能力與耐久性。
參考文獻:
[1]LAKES R. Foam structures with a negative Poisson's ratio[J]. Science, 1987, 235(4792): 1038-1040.
[2]RENX,DASR, TRANP, et al. Auxetic metamaterials and structures: a review[J]. Smart Materials and Structures, 2018,27(2): 023001.
[3]楊智春,鄧慶田.負泊松比材料與結(jié)構(gòu)的力學性能研究 及應用[J].力學進展,2011,41(3):335-350.
[4]BENDSOE M P, KIKUCHI N. Generating optimal topologies in structural design using a homogenization method[J]. Computer Methods in Applied Mechanics and Engineering, 1988, 71(2): 197-224.
[5]SIGMUND O, MAUTEK. Topologyoptimization approaches:a comparative review[J]. Structuraland Multidisciplinary Optimization, 2013, 48(6): 1031-1055.
[6]ANAYA-JAIMES L M, VICENTE W M, PAVANELLO R.Metamaterials design with a desired thermal expansion using a multi-material BESO method[J]. Structural and Multidisciplinary Optimization, 2022, 65(12): 355.
[7]MURAI N, NOGUCHI Y, MATSUSHIMA K, et al. Multiscale topology optimizationofelectromagnetic metamaterialsusingahigh-contrast homogenization method[J]. Computer Methods in Applied Mechanics and Engineering, 2023, 403: 115728.
[8]PAULINO GH, SILVA ECN,LE CH. Optimal design of periodic functionally graded composites with prescribed properties[J]. Structural and Multidisciplinary Optimization,2009, 38(5): 469-489.
[9]HA C S, HESTEKINE,LI J H, et al. Controllable thermal expansion of large magnitude in chiral negative Poisson's ratio latices[J]. Physica Status Solidi (B), 2015, 252(7): 1431-1434.
[10] AI L, GAO X L. Topology optimization of 2-D mechanical metamaterials using a parametric level set method combined with a meshfree algorithm[J]. Composite Structures, 2019, 229: 111318.
[11] QIN H X, YANG D Q. Vibration reduction design method of metamaterials with negative Poisson's ratio[J]. Journal of Materials Science, 2019, 54(22): 14038-14054.
[12] YE ML, GAO L, LI H. A design framework for gradually stiffer mechanical metamaterial induced by negative Poisson's ratio property[J].Materials amp; Design,2020,192: 108751.
[13] JEBELLAT E,BANIASSADI M, MOSHKI A, et al. Numerical investigation of smart auxetic three-dimensional meta-structures based on shape memory polymers via topology optimization[J]. Journal of Intelligent Material Systems and Structures, 2020,31(15): 1838-1852.
[14] AGRAWAL G, GUPTA A, CHOWDHURY R,et al. Robust topology optimization of negative Poisson’ s ratio metamaterialsundermaterialuncertainty[J].Finite Elements in Analysis and Design, 2022, 198: 103649.
[15]EVANS K E, NKANSAH M A, HUTCHINSON I J. Auxetic foams: modelling negative Poisson’ s ratios[J]. Acta Metallurgica et Materialia, 1994, 42(4): 1289-1294
[16] DUAN S Y,XI L,WEN W B,et al. A novel design method for 3D positive and negative Poisson's ratio material based on tension-twist couplingeffects[J]. Composite Structures, 2020,236: 111899.
[17] CLAUSEN A, WANG F W, JENSEN J S, et al. Topology optimized architectures with programmable Poisson's ratio over large deformations[J]. Advanced Materials, 2015, 27(37): 5523-5527.
[18] ZHANG D, ZHAI X Y, LIU L G, et al. An optimized, easy-to-use, open-source GPU solver for large-scale inversehomogenizationproblems[J]. Structuraland Multidisciplinary Optimization, 2023, 66(9): 207.
[19] VOGIATZIS P, CHEN S K, WANG X, et al. Topology optimization of multi-material negative Poisson’ sratio metamaterials using a reconciled level set method[J]. Computer-Aided Design,2017, 83: 15-32.
[20] ZHANG H, DING X H, WANG Q, et al. Topology optimization of composite material with high broadband damping[J]. Computers amp; Structures, 2020,239: 106331.
[21] ZHANG H K, LUO Y J, KANG Z. Bi-material microstructural design of chiral auxetic metamaterials using topology optimization[J]. Composite Structures, 2018,195: 232-248.
(編輯:董偉)