連續(xù)介質力學和微納觀理論的研究進展

0求解宏觀—引言隨著計算機的快速發(fā)展，數值計算方法（如金元法、邊界源法、有限差分法和無網格法）已成為解決各種復雜問題的有力工具。如何建立一套有效的多尺度關聯分析方法求解從宏觀到微觀問題是值得探討的課題。若沒有大尺度計算為小尺度計算的界面提供邊界條件,內部的多尺度描述就不會表述清楚;若小尺度計算所得的物理量不能通過一定的等效方法向大尺度計算傳遞,系統(tǒng)整體的物理、力學性能就不會被準確了解。因此,有必要建立尺度間的相互作用理論和多尺度級進計算方法,以實現不同尺度間的有效關聯。目前,主要的求解宏觀—微觀問題的多尺度數值計算方法有多尺度連續(xù)介質方法、連續(xù)介質—分子動力學耦合方法、準連續(xù)介質方法。本文根據所涉獵的有限領域討論了宏觀—微觀多尺度數值計算方法以及發(fā)展前景。1多尺度有限元計算方法隨著新技術、新材料的不斷發(fā)展,許多復雜問題用一般的連續(xù)介質方法,往往需要作網格密化或節(jié)點加密處理,但這樣勢必會增大計算量,使問題難以解決。于是,基于連續(xù)介質的多尺度方法得到了發(fā)展,對整體求解域作均勻化、大尺度的計算,對需要精細計算的局部區(qū)域采用局部密化或尺度關聯的小尺度計算,使上述問題的求解得到較好的實現。崔俊芝等提出了基于雙尺度漸近分析的有限元方法,并給出了詳細的計算步驟,用于解決在局部區(qū)域內間斷且跳躍性很大、區(qū)域內含有周期性洞穴或裂縫且周期很小的復合材料或周期結構的力學問題。曹禮群等研究了二維多孔復相介質穩(wěn)態(tài)溫度場的多尺度漸近數值計算方法。例如對二維多孔復相介質(圖1)問題:式中,ε為一個小周期參數,0<ε<l/L?1,取ε=1/6;l為基本周期單元(單胞)的細觀尺寸;ω為n維歐氏空間R中具有1—周期結構的無界區(qū)域。基本單胞Q∩ω如圖2所示,Q為單位方體。陰影部分表示周期孔洞,孔洞邊界上滿足第一類邊界條件。對區(qū)域Ωε和Q∩ω實行一致正方形剖分,如圖1、圖2所示,其剖分尺寸分別是h=1/36?ˉh=1/18h=1/36?hˉ=1/18,由于式(1)的精確解無法得到,可用較細網格上的有限元解作為原始問題溫度場uε(x)的近似解。設f(x)=36×106[x1x2(1-x1)(1-x2)]3,圖3是uε(x)的有限元解與多尺度有限元方法耦合解之間的絕對誤差曲面,無論是邊界附近或內部點均有較好的逼近。羅劍蘭等在以上研究基礎上,進一步解決了隨機多孔復相介質穩(wěn)態(tài)溫度場問題。Ladeveze等提出了一種可考慮摩擦效應的接觸問題的雙尺度計算模型。對于在相對于整體結構有較小尺度區(qū)域或局部區(qū)域,具有高梯度解,對不同特征尺度組成的異構材料如復合材料的問題,采用多尺度計算是一條可行的途徑。Ladeveze等提出的方法的基本思想是將結構分解成若干界面相關聯的子結構,子結構的界面可傳遞位移分布和力分布信息。在界面上的未知量(位移和力)可分解為S=SM+Sm(2)式中,SM為宏觀量;Sm為微觀量。在界面上采用力平衡準則強制邊界條件,這樣對于接觸邊界可較方便地引入摩擦力的作用。在具體計算中,采用稱為LATIN方法的宏觀—微觀計算模型組裝各子結構和界面,并對問題進行求解。在此研究中,僅考慮了線彈性情況。如圖4所示,矩形板(90mm×90mm)上受壓縮和彎曲載荷各100MPa,在板中有若干裂紋。在彎曲過程中,裂紋由于變形會發(fā)生接觸和摩擦,設摩擦系數f=1,楊氏模量E=130GPa,泊松比ν=0.2,平面應變問題。矩形板被分解為若干子結構。圖5給出了矩形板的變形。在計算中,宏觀尺度解反映了邊界力的總趨勢,而微觀尺度解則改善了裂紋尖端的計算結果。由于采用了子結構的多尺度關聯,這一方法適于并行計算。O?nn?ate提出了采用有限計算(finitecalculus,FIC)的多尺度有限元計算方法。FIC方法主要解決典型的多尺度問題,如具有局部高梯度的流體力學對流—擴散問題、固體或流體力學中的不可壓縮問題、固體的剪切帶預測或可壓縮流體的沖擊波問題等。此方法通過有限元計算對全局離散系統(tǒng)的節(jié)點提供了宏觀均勻化解,而對高梯度區(qū)域則通過遞階(hierarchical)耦合的多尺度有限元近似方法得到局部精確解。O?nn?ate給出的遞階耦合插值函數為Φ=n∑i=1ΝiΦi+m∑j=n+1Νh(huán)jahjΦ=∑i=1nNiΦi+∑j=n+1mNhjahj(3)式中,Ni為有限元形函數;Φi為節(jié)點位移;Nhjhj為遞階耦合形函數;ahjhj為遞階耦合的局部密化節(jié)點位移。此外,Juanes等提出了求解三相介質流動問題的多尺度分析方法,顯示了求解高非線性問題的潛力,并用此方法解決滲流區(qū)域的汽油過濾和碳氫化合物存儲器的水氣注射等實際問題。Li等將可變多尺度方法(variationalmultiscalemethod)和等效特征應變定律(equivalenteigenstrainprinciple)結合,提出了可變特征應變多尺度有限元方法(variationaleigenstrainmultiscalefiniteelementmethod),使多尺度計算中的粗糙尺度計算可獲得更好的精度。Masud等提出了多尺度穩(wěn)定有限元方法求解對流—擴散流體力學問題的解決方法。Gravemeier等提出了三層有限元多尺度分析方法求解了不可壓縮流體的Navier-Stokes方程?？梢?隨著當前工程技術問題的復雜度的不斷增加,考慮非線性因素的程度越來越高,多尺度分析方法的發(fā)展是必不可少的。2算例2—連續(xù)介質—分子動力學耦合方法在多尺度計算中,連續(xù)介質—分子動力學耦合方法是較常用且有效的方法。這種方法將不同空間或時間尺度對象的材料特性、物理特性、力學特性分別用適用于不同尺度的計算方法求解,通過建立一條宏觀—微觀—納觀尺度間的有效關聯,對多尺度問題進行求解。Liu等提出了一種多尺度遞階耦合的計算方法,對整體求解對象,采用連續(xù)介質的計算方法如有限元或無網格方法作均勻化計算,而在微結構與宏觀材料的交界面和微小結構內部,采用分子動力學或者能量密度函數方法計算其微/納觀性能,通過一定的關聯方法實現對問題的求解。這種方法可應用于材料力學性能的研究,為新材料的設計提供了有效的預測和指導手段。合金材料通常含有不同尺寸的雜質。其中有些雜質的引入在材料的生成過程中是不可避免的,它們可能會對材料剛度、強度、硬度等性能產生不良影響;另一些雜質的引入則是在材料生成過程中人為的加入以提高其力學性能。如何對這些材料進行合理的計算,確定引入微結構的種類、微結構尺寸的大小、微結構在材料中的排列方式等的影響,對于指導新材料的設計具有重要的現實意義。圖6所示為一個運用上述多尺度方法的簡單算例,對于含雜質(黑色部分)的材料,考察其承受剪切應力后的力學性能。圖7給出了受力過程中材料的應力—應變關系的變化,其中Ef表示不同材料之間的分離能?？梢钥闯?由于微結構的存在,首先在微結構與主材料界面間發(fā)生斷裂,使應力—應變曲線發(fā)生突變,表明材料性能在遭受破壞后迅速減弱。但通過合理的選擇所要加入雜質的類型以及調整相應的參數,將可以避免這種現象的發(fā)生并且可以提升材料的性能。Liu等進一步提出了一種多尺度并發(fā)(concurrent)耦合的關聯方法,用于實現研究對象中含有不同時間和空間尺度的材料或器件的有效聯接,或者在研究對象中需要考慮局部區(qū)域的微/納觀特性,如位錯、結合鍵分離等,而對其余區(qū)域作一般性考慮的問題。與以往提出的連續(xù)介質—分子動力學耦合方法相比,這一方法不需要在連續(xù)介質和分子動力學計算域之間建立重疊區(qū)域,或對局部區(qū)域進行加密計算,而采用多尺度映射的關聯方法實現多尺度計算。在計算中,對整個求解區(qū)域運用無網格方法,而在感興趣的區(qū)域嵌入分子結構,運用分子動力學進行計算。在此方法中,通過引入映射算子P,微/納觀尺度的位移被映射到宏觀尺度的位移形函數中,有uα=∑ΙΝΙ(xα)dΙ+dα-Ρdαuα=∑INI(xα)dI+dα?Pdα(4)式中,uα為原子位移;NI(xα)為計算原子α的節(jié)點I的有限元或無網格形函數;dI為節(jié)點位移;dα為分子動力位移。式中的Pdα是分子動力位移dα的映射,被稱為關聯尺度,用于唯一的分離宏觀和微觀尺度。通過合適的映射,可以只通過內力和邊界條件的信息傳遞實現尺度關聯,而剛度矩陣不需要進行耦合處理。這樣可以直接利用已有的有限元/無網格和分子動力學代碼對問題進行求解。圖8給出了多層碳納米管的純彎曲多尺度計算。其中,圖8a為無網格方法—分子動力學模擬耦合方法,圖8b為實驗測量的多層碳納米管的變形,圖8c為數值計算的多層碳納米管的變形。在計算中,將分子動力學粒子嵌入無網格離散系統(tǒng)。實際的多層碳納米管包含了上百萬個原子,而在此模擬中,總共使用了約30000個節(jié)點和80000個原子,大大減少了運算規(guī)模。經計算,多層碳納米管在彎曲過程中出現皺褶現象,與實驗結果相吻合。此外,Nattavut等研究了微觀非均勻化的粘土材料的多尺度方法,結合分子行為和宏觀連續(xù)介質現象,提出了分子動力學模擬和多尺度均勻化分析方法結合的計算流程。Horstemeyer等使用三種不同尺度的數值計算方法:嵌入原子方法(embeddedatommetohod)、分子動力學模擬和有限元方法,計算了金屬鎳面心立方(FCC)的簡單剪切問題。Gungor等運用分子動力學模擬和邊界元方法的遞階耦合研究了金屬薄膜的失效問題。Muralidharan等提出了有限差分時域方法/分子動力學混合方法(finitedifferencetimedomainmethod/moleculardynamicshybridmethod),研究了彈性波的傳播問題,考察了原子尺度的波動效應對材料響應的影響,在兩種方法的重疊部分采用了并發(fā)耦合。3能量密度函數準連續(xù)介質方法是一種將粒子系統(tǒng)簡化,不需要考慮大量的原子、分子數目而又能反映這些粒子信息和物理特性的數學方法。其基本思想是采用一種稱為Cauchy-Born的映射準則,將粒子系統(tǒng)視為由許多特定的圍繞表征原子的半徑為Rc的微晶結構等效成的連續(xù)介質,在微晶結構中包含了真實的原子或分子(圖9)。套用連續(xù)介質力學的本構方程,引入適當的原子模型勢能函數。通過在宏觀力學中已經被廣泛應用的有限元方法等求解問題。Tadmor等結合Cauchy-Born準則提出了晶體有限變形的思想,將有限元方法應用到晶格缺陷和復雜晶體的研究當中。這一準則假設系統(tǒng)中微晶結構滿足變形協調關系,通過變形梯度函數F將初始微晶結構和變形后的微晶結構之間的變形關系聯系起來(圖10):a=FA(5)F=?x?X(6)式中,a和A分別為變形后和初始狀態(tài)的晶格向量;x和X分別為變形后和初始狀態(tài)的材料坐標向量。因為變形梯度是材料變形前后空間坐標的函數,而變形后的空間坐標與初始狀態(tài)的空間坐標在時間上是相關的,所以變形梯度也是初始狀態(tài)空間坐標和時間的函數。原子模型通過適當的勢能函數可以獲得較好的近似。準連續(xù)介質方法的本構方程也必須引入這些勢函數近似,建立相應的應力、應變關系式。這里,可以直接利用非線性有限元或無網格方法的本構方程。已知勢函數U,有W=UA0Τ(7)式中,W為與粒子空間坐標相關的能量密度函數;A0為粒子系統(tǒng)未變形晶格的面積;T為粒子系統(tǒng)的厚度。第一Piola-Kirchhoff應力P和拉格朗日切向剛度張量C(Lagrangiantangentstiffnesstensor)分別為Ρ=?W?FΤC=?2W?FΤ?FΤ(8)由上可知,勢函數可將能量密度函數近似地引入到系統(tǒng)本構關系中,因此可以利用有限元等方法對問題求解。與傳統(tǒng)分子動力學模擬相比,準連續(xù)介質方法通過與一些數值方法如有限元方法等結合,可以簡化求解模型。例如,將三維分子動力學模型等效為二維有限元或無網格離散模型進行求解,同時可以用較少的節(jié)點信息表征大量原子、分子信息,節(jié)省了計算機資源,為求解大規(guī)模分子動力學問題提供了一條可能的途徑。人們已通過一些數值算例驗證了這一方法的可行性,并將其用于晶體斷層、晶粒邊界之間的相互作用、晶體結構斷裂、碳納米管各項力學性能等方面的模擬當中,展現了較廣闊的應用前景。Shenoy等將準連續(xù)介質方法用于求解多晶結構問題,給出了晶界結構和能量的有效處理方法。Miller等進一步用此方法求解原子尺度的晶體斷裂問題,在計算中考慮了晶體內部的斷裂和晶界之間斷裂傳遞的雙重效應。Knap等將準連續(xù)介質方法的應用擴展到三維晶體結構的計算當中并分析了它的計算精確性和收斂性。Smith等對二維和三維硅納米壓痕(Nanoindentation)作了模擬,并與實驗結果相對比,定性地再現了實驗過程,并從微觀角度反映了在壓痕過程中材料的相變和位錯(圖11),給出了應用stillingger-weber(SW)勢函數近似和非正交tight-binding(TB)哈密頓函數計算壓痕過程所得的材料位錯和相變。Dong等運用準連續(xù)介質方法和無網格方法結合分子動力學模擬研究了C60(或稱為富勒烯)和碳納米管之間的力學作用,討論了C60分子可以進入哪些尺寸的碳納米管,為制造具有較高強度的帶C60的納米管提供了參考。Arroyo等修正了納米薄膜的Cauchy-Born準則,利用準連續(xù)介質-有限元方法計算了納米管力學特性。由于碳納米管可視為石墨片卷曲而成,筆者采用tersoff-brenner勢函數近似的準連續(xù)介質—有限元方法對單層石墨片拉伸變形作了計算,得到材料等效楊氏模量E=1.030TPa、泊松比ν=0.326,與Dong等采用準連續(xù)介質—無網格方法計算得到的等效楊氏模量E=0.989TPa、泊松比ν=0.367接近。文獻給出的單壁碳納米管等效楊氏模量E=1.021TPa、泊松比ν=0.19。泊松比與前兩者的結果相差較大,Dong等認為這是勢函數選擇不合適的緣故。圖12給出了單層石墨片變形前后各原子的相對位置的變化。近來,Miller等對準連續(xù)介質方法作了綜述,介紹