单晶硅力学性能及尺寸效应的离散元模拟

1、单晶硅力学性能及尺寸效应的离散元模拟姜胜强谭援强李才等单晶硅力学性能及尺寸效应的离散元模拟姜胜强谭援强李才杨冬民湘潭大学,湘潭,411105摘要:通过建立紧密连接且随机分布在指定区域内的圆盘形颗粒离散单元的力学模型来研究单晶硅的力学性能,并采用不同尺寸的离散元模型研究力学性能尺寸效应。结果表明:泊松比、单轴抗压强度、弹性模量以及断裂韧性的尺寸效应不明显,弯曲强度则随着模型尺寸的增大而减小。建立了连续体经典断裂力学有限宽板条单边直裂纹的离散元力学模型,改变初始裂纹长度进行模拟得到的结果与单边切口梁模拟结果及文献报道的试验值一致。关键词:单晶硅;力学特性;模拟;离散元法中图分类号:TG58文章编号

2、:1004132X(2010)05058906StudyonMechanicsPropertiesandSizeEffectofMonocrystalleSiliconUsingDiscreteElementMethodJiangShengqiangTanYuanqiangLiYangXiangtanUniversity,Abstract:Themechanicsmodelsformedbyandrandomlydistrib2utedinaspecifiedregionwereusedtotofmonocrystalle-silicon,andsomedifferentsizesofoutt

3、ostudythesizeeffectofmechan2icsproperties.ThethatthereisnoobvioussizeeffectforPoissonsratio,compressivestYmodulusandfracturetoughness.However,thevalueofflexuralstrengthdecreaseshincreasingofthemodelsize.Thediscreteelementmodelofclassicalcon2tinuumfracturemechanicswasusedtosimulatethefractureprocesso

4、fsingleedgeverticalcrackplatewithfinite-widthundertensilestressininfinitedistance.Thesimulationresultswhichobtainedbychangingtheinitiatecracklengthofmodelarejustthesamewiththesimulationresultsofsingleedgenotchedbending(SENB)andtheexperimentalvaluesreportedinotherliteratures.Keywords:monocrystallesil

5、icon;mechanicsproperty;simulation;discreteelementmethod0引言随着MEMS技术的兴起和发展,硅材料MEMS器件已有大量应用。MEMS微器件及其间隔之间的尺寸为微米级甚至亚微米级,材料呈现显著的尺寸效应12。尺寸效应是脆性材料和准脆性材料的基本特性。Ba ant等3提出了很多关于脆性材料强度尺寸效应的理论,并指出离散元在研究脆性材料力学性能时存在尺寸效应。Han等4通过单晶硅的微桥梁实验得出单晶硅晶硅的断裂过程,却不能像离散元法那样模拟裂纹形成与扩展的过程7。如果构成离散元模型的颗粒太少,则模拟结果将存在明显的离散性。由于计算机能力的限制,模

6、型颗粒的数目不可能无限制扩大,因此选取一个合适的颗粒数目进行离散元模拟显得十分重要。本文用离散单元法研究单晶硅的力学性能及其尺寸效应,建立经典断裂力学有限宽板条单边直裂纹的拉伸断裂离散元模型,模拟不同初始裂纹长度下的断裂韧性。弹性模量没有明显的尺寸效应,而弯曲强度则存在明显的尺寸效应的理论。刘凯欣等5用颗粒离散元模拟了脆性材料的动态破坏过程。离散单元法在模拟脆性材料力学性能时有自身的优点,可以模拟脆性材料的大变形与微观尺度上材料断裂过程中裂纹的形成和扩展6。虽然分子动力学也能在微观尺度上模拟单收稿日期:20090429基金项目:国家自然科学基金资助项目(50875224);教育部新世纪优秀人才

7、支持计划资助项目(NCET06-0708)1模拟过程1.1离散元模型在建立单晶硅离散元BPM力学模型时,必须先确定模型中单元的尺寸。Potyondy等8提出了一个将离散元模型单元属性与材料断裂韧性联系起来的方程:RKIC=t(1)为因数,=1,=0166;式中,、t为模型中材料的拉伸强度;R为模型中单元的平均半径,模型中最大半径与最589中国机械工程第21卷第5期2010年3月上半月小半径之比取115。表2单晶硅(100面)二维离散元模型与实测主要力学性能比较力学特性文献1011的结果弹性模量E(GPa)130泊松比0.278抗压强度950c(MPa)弯曲强度3001000b(MPa)断裂韧性

8、KIC0.601.65(1/2)目前还没有一种方法或者理论能将材料宏观的力学性能与微观的参数相匹配起来,因此要通过一系列的数值模拟试验来进行校准,如单轴压缩试验、三点弯曲试验、单边切口梁试验等,最终获得匹配单晶硅力学特性的离散元颗粒单元的微观参数。m,所本文模型中的颗粒最小半径为0107用的模拟工具是Itasca公司的二维离散元软件PFC2D,利用内置AugmentedFishtank模块获得BPM力学模型9,如图1所示。离散元模拟模型尺寸如表1所示。离散元模型结果129.80.273954.37540.822(a)试验应力-应变曲线(a)单轴压缩试验(b)模型失效图图2单轴压缩离散元模拟试验

9、(b)三点弯曲试验(a)加载曲线(c)单边切口梁试验图1离散元力学模型表1离散元模拟模型尺寸模型类型单轴压缩试验三点弯曲试验单边切口梁试验xy(模型尺寸)颗粒数目118041062310623m×m)12(1326m×m)51(39.07.8m×m)51(39.07.8(b)模型失效图图3单边切口离散元模拟试验1.3离散元力学模型尺寸效应计算1.2力学特性的计算结果改变生成离散元力学模型颗粒的随机数9,每种模型分别选取10个随机数进行模拟,最后将每种模型模拟结果的平均值作为离散元模拟结果。通过离散元模拟,我们获得了单晶硅的一些主要力学性能参数,如表2所示。单晶硅三

10、种数值模拟试验力学特性的二维离散元模拟结果如图2图4所示。590在每种模型的最大尺寸和最小尺寸之间保持模型的长高比不变,即单轴压缩试验中xy=12,三点弯曲和单边切口梁试验的xy=51,再分别选取13个逐渐变化的尺寸模型(表3、表4)进行模拟。单轴压缩模型试验中尺寸较m×m)的模型选取20个随机数进小(小于612行模拟,尺寸较大的模型则选取10个随机数进行模拟。尺寸效应的模拟结果如图5所示,其结果与Han等4通过微桥梁实验得出的结论非常一单晶硅力学性能及尺寸效应的离散元模拟姜胜强谭援强李才等m(4203个颗粒)时,其值随着模型尺寸于2410的增大而快速增大,尺寸效应比较明显,当大于该

11、长度时,继续保持增长趋势,但增长变得极其缓慢。表3单轴压缩试验的模型尺寸和颗粒数目模型尺寸m×m12m×m24m×m36m×m48m×m510m×m612mm714颗粒数目692796281117174625143422(a)加载曲线模型尺寸m×m816m×m918m×m1020m×m1122m×m1224m×m1326颗粒数目44705657698484511005711804表4三点弯曲以及单边切口梁(b)模型断裂图试验的模型尺寸和颗粒数目模型尺寸m×m30.6

12、m×m61.2m×91.8120mm×m214.2颗粒数目15225122633080图4三点弯曲离散元模拟试验致,也与Koyama等12模拟岩石尺度效应的结果相似;每一模型尺寸模拟结果平均值如图6所示。从图5可以看出,抗压强度、弹性模量、泊松比在m(628个颗粒)模型长度l小于6的随机性,向于收敛;同样,当模)时,模拟结果型的长度小于15随机性很大,随着模型尺寸的增大,模拟结果逐渐收敛。从图6可以看出,随着模型尺寸的增大,抗压强度、弹性模量在模型尺寸较小时随模型尺寸m(4470的增大逐渐增大,当模型长度l大于8个颗粒)时,模拟结果平均值基本保持不变;泊松比基本不

13、随模型尺寸的增大而发生改变,而保持一稳定值;弯曲强度则随着模型尺寸的增大而减小,尺寸效应比较明显;断裂韧性在模型长度l小模型尺寸m×m244.8m275.4×m6.0m6.6mm367.2m×m397.8颗粒数目4023509162867606905210623为了分析模拟结果的随机性,本文对抗压强度、泊松比、弹性模量进行了频率统计并分析了不同模型长度时的概率密度分布函数,如图7和图8所示。从图7可以看出,模拟得出的宏观力学性能参数(抗压强度、泊松比、弹性模量)近似符合正态分布;从图8可以看出,抗压强度、泊松比、弹性模量的模拟值随着模型尺寸的增大而快速收m(565

14、7个颗粒)时,收敛敛,当模型长度大于9的速度变得缓慢,其对应的值基本上稳定的处于一个小范围内。(a)抗压强度(xy=12)(b)弹性模量(xy=12)(c)泊松比(xy=12)(d)弯曲强度(xy=51)(e)断裂韧性(xy=51)图5尺寸效应模拟结果591中国机械工程第21卷第5期2010年3月上半月(a)抗压强度(xy=12)(b)弹性模量(xy=12)(c)泊松比(xy=12)(d)弯曲强度(xy=51)(=51)图6(a)抗压强度(b)泊松比(c)弹性模量图7宏观力学性能参数出现次数的直方统计图(xy=12)(a)抗压强度(b)泊松比mmmmm1.l=12.l=33.l=64.l=95

15、.l=12(c)弹性模量图8不同模型长度时的概率密度分布函数(xy=12)通过对图8所示的概率密度分布函数进行统计分析,得出了不同模型长度时材料力学性能参数的数学期望以及标准差,如表5所示。从表5中可以得出同样的结论:随着模型尺寸的增大,抗压强度和弹性模量逐渐增大,泊松比基本保持不变;当m(69个颗粒)时,小尺寸离散性大模型长度为1的原因使得抗压强度的标准差非常大,达到20914MPa,但随着模型长度的增大,各力学参数对应的标准差随之减小,最后在一个稳定的范围内。592表5不同模型长度的材料力学性能参数的数学期望及标准差c(MPa)模型抗压强度长度l数学标准差(m)期望泊松比数学期望0.270

16、0.2780.2840.2810.280弹性模量E(GPa)数学期望121.5127.4128.4129.2129.8标准差0.0680.0270.0050.0080.009标准差12.704.462.531.941.72136912866.7917.3973.9972.0954.3209.474.836.730.722.8单晶硅力学性能及尺寸效应的离散元模拟姜胜强谭援强李才等2应力场强度因子的离散元模型拉伸断裂模拟如图9所示。线弹性断裂力学中,均匀受载的各向同性弹性连续体有限宽板条单边直裂纹的受力分析如图9a所示,其裂纹尖端的应力场强度因子为KI=cF(c/b)F(c/b)=1.12-0.2

17、31c/b+10.55(c/b)2m)处Y方向拉应力,如图行中记录模型(0,157d所示,将最大应力代入式(2)获得裂纹扩展时的应力场强度因子(材料的断裂韧性)。选取10个随机数进行模拟,模拟结果平均值为KIC=01825MPa/m1/2。另外将模型初始裂纹长度改为m,同样也选取10个随机数进行模拟,结果为2KIC=01852MPa/m1/2(2)-(3);与单边切口梁离散元模拟21.72(c/b)3+30.39(c/b)4式中,为无穷远处应力;c为裂纹初始长度;F(c/b)为有限宽板修正系数;b为有限宽板条的宽度。得到的结果(KIC=01822MPa/m1/2)相比,两者之间的误差大约是4%

18、,且都在文献报道的实验值(KIC=0161165MPa/m1/2)11之内,从而说明用离散元法建立的模型是合理可行的,其模拟有限宽板条单边直裂纹拉伸断裂的断裂过程及模拟结果与分子动力学模拟基本一致:材料以脆性方式断裂,新的裂纹一旦出现,就以很快的速度扩展,材料在瞬间发生断裂,方向,KIC7。当无穷远处拉应力逐渐增大时,初始裂纹尖端的应力场强度因子也增大,当裂纹尖端应力场强度因子达到材料的断裂韧性时,裂纹开始扩展13。利用离散元建立有限宽板条单边直裂纹的离散元力学模型,如图9b所示,模型的尺寸为m×m,初始裂纹长度为1m,颗粒数为321011175,模型的上面两层颗粒有Y度01005m

19、/s,),当只改变生成模型的随机数时,计算结果之间存在的离散性使小尺寸模型的离散性较大,随着模型尺寸的不断增大,同一模型模拟结果之间的离散性逐渐减小,最后趋于稳定。(2)随着模型尺寸的增大,抗压强度和弹性模量在模型尺寸较小时随模型尺寸的增大逐渐增m(4470个颗粒)时,模拟大,当模型长度大于8结果平均值基本保持不变;泊松比基本不随着模型尺寸的增大而发生改变,保持一稳定值;弯曲强度则随着模型尺寸的增大而减小,尺寸效应比较m(4203个明显;断裂韧性在模型长度小于2410颗粒)时,其值随着模型尺寸的增大而快速增大,尺寸效应比较明显,当大于该长度时,继续保持增长趋势,但增长变得极其缓慢。(3)数值模

20、拟得出的宏观力学性能参数(抗压强度、泊松比、弹性模量)近似符合正态分布,并且(a)有限宽板条单边直裂纹力学模型(b)单轴拉伸断裂的离散元模型(c)单轴拉伸试验断裂图模拟值随着模型尺寸的增大而快速收敛,当模型m(5657个颗粒)时,收敛变得缓慢,长度大于9其对应的值基本上稳定的处于一个小范围内。(4)随着模型颗粒数目的增大,材料的宏观力学性能参数逐渐趋向稳定;颗粒数目越小,模拟结果的离散性越大,颗粒数目越大(至少6000),其结果越可靠。(5)建立了经典断裂力学有限宽板条单边直裂纹应力场强度因子的离散元模型,通过改变离散元模型初始裂纹的长度来获得不同初始裂纹长度下材料在拉伸情况下的断裂韧性,模拟

21、结果与593(d)单轴拉伸试验应力-位移曲线图9拉伸断裂模拟的示意图中国机械工程第21卷第5期2010年3月上半月单边切口梁试验模型获得结果的误差在4%以内,且都在文献报道的实验范围之内。参考文献:1XueZ,SaifMTA,HuangY.TheStrainFradientEffectinMicroelectromechanicalSystems(MEMS)J.JournalofMicroelectromechanicalSystems,2002,11(1):2735.2SerryFM,WalliserD,MaclayGJ.TheAnharmonicCasimirCscillator(ACO)

22、theCasimirEffectinaModelMicroelectromechanicalSystemJ.JournalofMicroelectromechanicalSystems,1995,4(4):193205.3Ba antZP,ChenEr-ping.ScalingofStructuralFail2ureJ.AppliedMechanicsReviews,1997,50(10):593发表论文2篇。谭援强,男,1966年生。湘潭大学机械工程学院教授、博士研究生导师。李才,男,1982年生。湘潭大学机械工程学院硕士研究生。杨冬民,男,1984年生。湘潭大学机械工程学院助教。(上接第5

23、88页)该方法通过单元的删除,将实体单元改变为承载能力很弱且弹性模量和密度匹配的孔洞单元,防止了局部模态。双向渐进结构法对初始设计域敏度再分配技术的应用,使各单元对结构性能的贡献实现了平滑过渡,有效地抑制了棋盘格式。ANSYS与MATLAB的混合编程提高了优化效率。算例表明,方法正确和有效,具有很好的工程应用性。参考文献:1郭中泽.627.4HanGuangping,LiuKai.MechanicalPropertiesandSizeEffectsofSingleCrystalSiliconJ.ChineseJour2nalofMechanicalEngineering,2006,19(2):

24、290293.5刘凯欣,郑文刚,高凌天.脆性材料动态破坏过程的方法J.,25(8):929931.2,2006.M,StevenGP,XieYM.EvolutionaryStruc2turalOptimizationUsinganAdditiveAlgorithmJ.Fi2niteElementsinAnalysisandDesign,2000,34:291308.4XieYM,StevenGP.EvolutionaryStructuralOptimiza2tionforDynamicProblemsJ.Computers&Struc2tures,1999,58(6):10671073

25、.5XieYM,StevenGP.ASimpleEvolutionaryProcedureforStructuresUsinganEvolutionaryMethodJ.Com2puters&Structures,1993,49(5):885896.6谢忆民,杨晓英,StevenGP,等.渐进结构优化法的基数值模拟J.计算力学学报,2003,20(2)6TanYuanqiang,YangPolycrystalline23DiscreteElementJournalofMachineTools&Manufacture,2008,48(9):975982.7FujiiT,Akini

26、waY.MolecularDynamicsAnalysisforFractureBehaviorofSingleCrystalSiliconThinFilmwithMicronotchJ.ModelingandSimulationinMa2terialsScienceandEngineering,2006,14:7383.8PotyondyDO,CundallPA.ABonded-particleModelforRockJ.InternationalJournalofRockMechanics&MiningScience,2004,41(8):13291364.9PFC2DUsersManualsM.Lakewood,Colorado,USA:ItascaConsultingGroupInc.,2002.10刘加富.线锯切割单晶硅的应力场及损伤层研究D.济南:山东大学,2006.11SundararajanS,BhushanB.FractureDevelopmentofAFM-basedTechniquestoMeasureMechanicalPropertiesofNan