专利名称：计算硅锗超晶格材料界面热阻的方法
技术领域：
本发明涉及硅锗超晶格材料界面热阻的领域。
背景技术：
硅锗超晶格材料是ー种由硅、锗分子按照周期性结构排列而形成的复合材料，由于它具有良好的热电性质，能显著提高热电设备效率，在计算机芯片、MEMS器件、航空航天等领域有广泛应用。硅锗超晶格材料的界面热阻是微尺度传热研究領域的热点问题，由于它反应的是界面层的电子与声子的传输性质，因此直接影响材料的热传导性能，从而对微/纳米器件的设计和热优化产生影响。
目前通过微尺度传热理论分析和实验方法尚未对硅锗超晶格材料界面热阻的内部机制和变化规律建立完善的理论体系和研究方法，存在许多亟待澄清和解决的科学问题。从实验研究的角度而言难点主要是实验材料难以制备，试件大多为多晶体，无法独立讨论材料的界面热阻效应。另外实验条件、试件的材料和制作エ艺的差异，也使实验结果没有可比性。
近几年发展起来的分子动力学模拟是理论和实验有效的补充手段，它可以从微观角度求解有相互作用的各个分子的运动方程，得到每个分子空间位置和运动状态随时间的演进状况，从而统计出材料的宏观行为特性。分子动力学模拟的应用范围涉及机械工程、化学、材料科学等众多领域，但利用分子动力学计算硅锗超晶格材料界面热阻的方法，目前国内外未见报道。
发明内容
为了解决目前通过理论和实验方法不能准确分析硅锗超晶格材料传热机理的问题，提供ー种计算硅锗超晶格材料界面热阻的方法。
计算硅锗超晶格材料界面热阻的方法，具体步骤如下 步骤ー建立硅锗超晶格材料的非平衡态分子动力学导热模型； 步骤ニ 设定步骤一中导热模型中粒子的初始状态，所述初始状态包括粒子的初始位置和速度，硅、锗粒子的初始位置分布在面心立方结构的晶格点阵上，而粒子的初始速度的分布采用麦克斯韦速度分布；将系统温度调整到要求温度，重新标定硅、锗超晶格材料中粒子在要求温度下的速度； 步骤三根据步骤ニ中获得的粒子速度，采用Stillinger-Webber多体势能函数来计算硅、锗分子之间的作用势能； 步骤四根据分子动力学分析方法，用步骤三中计算得到的势能函数的梯度来计算硅锗超晶格材料中粒子间作用力； 步骤五由步骤四获得的粒子间作用力，求得硅锗超晶格材料中粒子的运动方程，采用Verlet蛙跳算法对运动方程进行积分，计算粒子在t时刻的新位置和速度； 步骤六采用各向异性非平衡态分子动力学方法来计算界面热阻，选用恒温墙方法对系统施加温度梯度，计算导热方向每ー层粒子的温度，并且对温度进行量子化修正；利用公式
来计算硅锗超晶格结构的界面热阻，其中AT表示温度梯度方向的
温差，Q为热流量，An为面YOZ的横截面积，ATx为X方向的温差Jx为X方向的热流。
本发明建立硅锗超晶格材料的非平衡态热传导模型，在分子水平上探讨界面热阻的变化规律，本发明方法简单易行，结果稳定可靠，对超晶格材料热导率研究具有指导意义。
图I是计算硅锗超晶格材料界面热阻的方法流程图；图2是硅锗超晶格材料的非平衡态分子动力学导热模型；图3是硅锗界面的温度示意图；图4是不同温度下界面热阻的曲线图；图5是不同周期长度的界面热阻曲线图。
具体实施例方式
具体实施方式
一、结合图I说明本实施方式，计算硅锗超晶格材料界面热阻的方法，具体步骤如下 步骤ー建立硅锗超晶格材料的非平衡态分子动力学导热模型； 步骤ニ 设定步骤一中导热模型中粒子的初始状态，所述初始状态包括粒子的初始位置和速度，硅、锗粒子的初始位置分布在面心立方结构的晶格点阵上，而粒子的初始速度的分布采用麦克斯韦速度分布；将系统温度调整到要求温度，重新标定硅、锗超晶格材料中粒子在要求温度下的速度； 步骤三根据步骤ニ中获得的粒子速度，采用Stillinger-Webber多体势能函数来计算硅、锗分子之间的作用势能； 步骤四根据分子动力学分析方法，用步骤三中计算得到的势能函数的梯度来计算硅锗超晶格材料中粒子间作用力； 步骤五由步骤四获得的作用力，求得硅锗超晶格材料中粒子的运动方程，采用Verlet蛙跳算法对运动方程进行积分，计算粒子在t时刻的新位置和速度； 步骤六采用各向异性非平衡态分子动力学模拟方法来仿真界面热阻，选用恒温墙方法对系统施加温度梯度，计算导热方向每ー层粒子的温度，并且对温度进行量子化修正；利用公式及来计算硅锗超晶格结构的界面热阻,其中AT表示温度梯度方向的温差，Q为热流量，An为面YOZ的横截面积，ATx为X方向的温差Jx为X方向的热流。 具体实施方式
ニ、本实施方式是对具体实施方式
一中步骤ー的进ー步说明以硅锗材料交替生长方向为X方向，X方向、Y方向与Z方向满足右手定则，在X方向布置随机恒温热墙，以建立温度梯度，高温热墙和低温热墙的材料均与各自临近的超晶格材料相同 ’导热模型的最外层设置绝热壁；在Y方向、Z方向施加周期性边界条件，并且面YOZ的截面面积大于3UCX3UC(Unit Cell,晶胞长度)，选取YOZ横截面积An为4UCX4UC。
Si的晶格參数是
的晶格參数aGe =5.657A = 0.5657nm。高温热墙和低温热墙至少设置3层晶格单元，即厚度至少为3UC，绝热壁也需要设置2层晶格単元，即厚度至少为2UC。在热流方向上(即X方向)采用固定边界条件，通过改变晶格单元数目，实现不同的周期长度和周期数，需要保证垂直热流方向导热层的总厚度小于lOOnm。YOZ横截面积An选取4UCX4UC的目的是减少计算量。
导热模型的最外层设置绝热壁的作用是减少导热层内的粒子蒸发，并且设定该区域粒子的速度为O。 具体实施方式
三、本实施方式是对具体实施方式
一中步骤ニ的进ー步说明硅、锗粒子的初始位置分布在面心立方结构的晶格点阵上，设定晶体単元的边长是2，其初始位置是 (0,0,0) ;(1，1，0) ;(1，0，1) ;(0，1，1)； (I. 5,0. 5,0. 5) ； (0. 5,1. 5,1. 5) ； (0. 5,0. 5,1. 5) ； (I. 5,1. 5,1. 5)； 粒子的初始速度的分布采用麦克斯韦速度分布，即每个粒子在X，Y，Z三个方向的分速度在区间(-0. 5,0. 5)之间均匀随机分布的速度； 根据经典的Boltzmann统计，将系统温度调整到要求温度，标定硅锗超晶格材料中粒子在要求温度下的速度，分两步进行 步骤_.一、保证系统的总动量为零 公式一
其中，N表示总的粒子数，vixold表示设定温度前粒子i在X方向的速度分量，vixnew表示满足设定温度后的粒子i在X方向的标定的速度分量； 步骤ニニ、标定硅锗超晶格材料中粒子速度，来满足给定温度
其中，Tref是參考温度，Treq是系统的要求温度，m表示粒子的总质量。 具体实施方式
四、本实施方式是对具体实施方式
一中步骤三的进ー步说明采用Stillinger-Webber多体势能函数来描述娃、锗分子之间的作用势能，N个粒子的总势能为
两体式f2、三体式f3分别用下式表达
其中rx、ry、rz表示粒子X、Y、Z方向的位置向量,r表示位置向量,r = (rx,ry,rz)；O为平衡常数；e为势阱常数；r。表示截断半径，r^.表示粒子i和粒子j之间的距离，r^-= YjiJjk表示粒子j和粒子k之间的距离，rJk = _rkJ，rik表示粒子i和粒子k之间的距离，rik = _rki, 0 iJk为rid和rik的夹角，9 Jik为rdi和rjk的夹角；9 kiJ为rki和rkj的夹角，、
和れ为势能參数，び=-:S> 2 (k，e = (esieGe)1/2。
娃错Stillinger-Webber势能參数由表I列出 表I
为了减小计算中的误差，必须对物理量进行无量纲化处理，长度单位CTsi =2.095U , (Tck =2.181 A ； ￡ Si = 3. 4723X IO^19J, e Ge = 3. 085XKT19J ;力的単位 f =
e /o ;时间单位T= [mo2/e]1/2 ;温度单位T*= e/Kb (Kb = I. 3805 X 10'为玻尔兹曼常数);速度単位V* = [e /m]1,2 ;热流单位J* = m(v*)3/ o 2 ;热阻单位R* = J/T*。 具体实施方式
五、本实施方式是对具体实施方式
一中步骤四的进ー步说明根据分子动力学仿真分析方法，用势能函数的梯度AiQ来计算カFi,由于粒子间相互作用カ随着距离增加而迅速衰减，选定ー个截断半径r。，只计算以截断半径为球体内的粒子间的作用力，而当粒子间的距离超过截断半径，则不考虑它们之间的作用。
本发明采用Cell列表法为Velet列表法建立列表，截断半径r。= 2. 5 o，用Velet列表法计算原子间作用力，这样可以提高计算速度。 具体实施方式
六、本实施方式是对具体实施方式
一中步骤五的进ー步说明积分所有粒子的Newton运动方程，获得下一个时间步长内粒子的位置。本发明采用Verlet推导的蛙跳(Leapfrog)算法对运动方程进行积分
任意时间层上的速度
其中r (t)表示粒子t时刻的位置，加一点代表微分一次，得到速度，加两个点代表微分两次，得到加速度。
蛙跳算法中离散方程的数值积分时间步长为A t = 0. 005 ^ ^ 0. Olps0 具体实施方式
七、结合图2说明本实施方式，本实施方式是对具体实施方式
一中步骤六的进ー步说明采用各向异性非平衡态分子动力学模拟方法来仿真界面热阻，选用恒温墙方法对系统施加温度梯度，任一平面的温度为 公式四
其中，Np是第p层粒子数。
当热域或冷域的粒子域与其他粒子发生作用,温度势能必会发生改变。要使接触面的温度和热域的温度保持一致，在内部形成稳定的温度梯度，需要通过公式一至公式三的速度标定法来实现。
公式四成立的前提条件之一是系统温度远高于材料的debye温度0 D，由于平衡温度低于硅的Debye温度(645K)，因此需对系统的局域稳定进行量子化修正，才能获得超晶格结构热导率的真实值，设定系统总能量是动能的2倍，通过求解积分方程得到系统的真头温度 公式五
其中，公式五右边是系统中粒子的总能量，D(co)为声子密度函数；《为声子频率；〈n>表示对应于热平衡温度T的声子平均占有数，该占有数满足Planck分布，即Bose-Einstein 统计； 界面热阻的表达式为
其中，ATx为X方向的温差；Jx为X方向的热流，A T表示温度梯度方向的温差。Q为热流量，An为面YOZ的横截面积。
热流Jx的计算公式如下
其中V表示系统体积；Vi表示第i个粒子速度；ru表示第i个粒子和第j个粒子的距离ばu表示第i个粒子和第j个粒子的作用力；ち表示第i个粒子的总能量，可用公式表示为
其中Iiii表示粒子i的质量。 具体实施方式
八、结合图3、图4和图5说明本实施方式，本实施方式是利用本发明的方法求解硅锗超晶格材料界面热阻实例，高低温热墙的温度设为Thrt = T+20K及Trald= T-20K。图3是硅锗超晶格结构X方向导热层的温度分布曲线，受界面热阻效应的影响在Si、Ge界面处存在很明显的温度跳跃，并且在温度较高的第一个界面的突变要比温度较低界面的跳跃更加明显，图4是周期长度为IOUC的硅锗超晶格结构温度与界面热阻的变化关系曲线，可以看出随着仿真温度的升高，界面热阻越来越小；图5是平均温度为400K时，周期长度与界面热阻的变化关系曲线，随着超晶格结构周期长度的増加，界面热阻逐渐减小。总的计算步长数为7X107，其中前IO6步使系统趋于平衡。
权利要求
1.计算硅锗超晶格材料界面热阻的方法，其特征在于具体步骤如下 步骤一建立硅锗超晶格材料的非平衡态分子动力学导热模型； 步骤二 设定步骤一中导热模型中粒子的初始状态，所述初始状态包括粒子的初始位置和速度，硅、锗粒子的初始位置分布在面心立方结构的晶格点阵上，而粒子的初始速度的分布采用麦克斯韦速度分布；将系统温度调整到要求温度，重新标定硅、锗超晶格材料中粒子在要求温度下的速度； 步骤三根据步骤二中获得的粒子速度，采用Stillinger-Webber多体势能函数来计算硅、锗分子之间的作用势能； 步骤四根据分子动力学分析方法，用步骤三中计算得到的势能函数的梯度来计算硅锗超晶格材料中粒子间作用力； 步骤五由步骤四获得的粒子间作用力，求得硅锗超晶格材料中粒子的运动方程，采用Verlet蛙跳算法对运动方程进行积分，计算粒子在t时刻的新位置和速度； 步骤六采用各向异性非平衡态分子动力学方法来计算界面热阻，选用恒温墙方法对系统施加温度梯度，计算导热方向每一层粒子的温度，并且对温度进行量子化修正；利用公式
来计算硅锗超晶格结构的界面热阻，其中ΔΤ表示温度梯度方向的温差，Q为热流量，An为面YOZ的横截面积，ATx为X方向的温差Jx为X方向的热流。
2.根据权利要求I所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤一具体为以硅锗材料交替生长方向为X方向，X方向、Y方向与Z方向满足右手定则，在X方向布置随机恒温热墙，以建立温度梯度，高温热墙和低温热墙的材料均与各自临近的超晶格材料相同；导热模型的最外层设置绝热壁；在Y方向、Z方向施加周期性边界条件，并且面YOZ的截面面积大于3UCX3UC(Unit Cell,晶胞长度)，选取YOZ横截面积An为4UCX4UC。
3.根据权利要求I所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤二具体为硅、锗粒子的初始位置分布在面心立方结构的晶格点阵上，设定晶体单元的边长是.2，其初始位置是(0,0,0) ;(1，1，0) ;(1，0，1) ;(0，1，1)； (I. 5,0. 5,0. 5) ；(0. 5,1. 5,1. 5) ; (O. 5，O. 5，I. 5) ； (I. 5,1. 5,1. 5)； 粒子的初始速度的分布采用麦克斯韦速度分布，即每个粒子在X，Y，Z三个方向的分速度在区间(-0. 5,0. 5)之间均匀随机分布的速度； 根据经典的Boltzmann统计，将系统温度调整到要求温度，标定硅锗超晶格材料中粒子在要求温度下的速度，分两步进行 步骤二一、保证系统的总动量为零 公式一
其中，N表示总的粒子数，vix°ld表示设定温度前粒子i在X方向的速度分量，Vixn 表示满足设定温度后的粒子i在X方向的标定的速度分量； 步骤二二、标定硅锗超晶格材料中粒子速度，来满足给定温度 公式二
其中，Tref是参考温度，Treq是系统的要求温度，m表示粒子的总质量。
4.根据权利要求I所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤三具体为采用Stillinger-Webber多体势能函数来描述娃、锗分子之间的作用势能，N个粒子的总势能为
两体式f2、三体式f3分别用下式表达
其中rx、ry、rz表示粒子X、Y、Z方向的位置向量，r表示位置向量，r = (rz, ry, rz)；σ为平衡常数；ε为势阱常数；r。表示截断半径，riJ表示粒子i和粒子j之间的距离，ru=Tji, rJk表示粒子j和粒子k之间的距离，rJk = -rkJ, rik表示粒子i和粒子k之间的距离，rik = -rki, Θ iJk为rij和rik的夹角，Θ Jik为rji和rjk的夹角；Θ kiJ为rki和rkj的夹角，Y 和 Η 为势能参数，2 °e，ε = ( ε Si ε Ge)1/2。
5.根据权利要求I所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤四具体为根据分子动力学仿真分析方法，用势能函数的梯度来计算力Fi,由于粒子间相互作用力随着距离增加而迅速衰减，选定一个截断半径r。，只计算以截断半径为球体内的粒子间的作用力，而当粒子间的距离超过截断半径，则不考虑它们之间的作用。
6.根据权利要求I所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤五具体为积分所有粒子的Newton运动方程，获得下一个时间步长内粒子的位置，本发明采用Verlet推导的娃跳(Leapfrog)算法对运动方程进行积分
任意时间层上的速度
其中r(t)表不粒子t时刻的位置,加一点代表微分一次,得到速度,加两个点代表微分两次，得到加速度。
7.根据权利要求I或3所述的计算硅锗超晶格材料界面热阻的方法，其特征在于步骤六具体为采用各向异性非平衡态分子动力学模拟方法来仿真界面热阻，选用恒温墙方法对系统施加温度梯度，任一平面的温度为
其中，Np是第P层粒子数； 需要通过公式一至公式三的速度标定来实现接触面的温度和热域的温度保持一致，在内部形成稳定的温度梯度； 对系统的局域稳定进行量子化修正，获得超晶格结构热导率的真实值，系统总能量是动能的2倍，通过求解积分方程得到系统的真实温度 公式BNkBTMD = ][ Ζ)(ω)( (ω,Γ)) /ω 其中，公式五右边是系统中粒子的总能量，D(co)为声子密度函数；ω为声子频率；<η>表示对应于热平衡温度T的声子平均占有数，该占有数满足Planck分布，即Bose-Einstein统计； 界面热阻的表达式为
其中，ATxS X方向的温差JxSX方向的热流，Λ T表示温度梯度方向的温差；Q为热流量，K为面YOZ的横截面积； 热流Jx的计算公式如下
其中V表示系统体积；Vi表示第i个粒子速度；ru表示第i个粒子和第j个粒子的距离表示第i个粒子和第j个粒子的作用力Wi表示第i个粒子的总能量，可用公式表示为
其中Hii表示粒子i的质量。
全文摘要
计算硅锗超晶格材料界面热阻的方法，涉及计算硅锗超晶格材料界面热阻，解决了目前通过理论和实验方法不能准确分析硅锗超晶格材料传热机理的问题，具体步骤如下A建立硅锗超晶格材料的非平衡态分子动力学导热模型；B设定导热模型粒子的初始状态，将系统温度调整到要求温度，标定硅锗超晶格材料中粒子在要求温度下的速度；C确定粒子间的作用势能；D计算硅锗超晶格材料中粒子间的作用力；E根据牛顿第二定律，求得硅锗超晶格材料中粒子的运动方程，积分运动方程，求得硅锗超晶格材料中粒子的运动速度；F运用非平衡态分子动力学方法求解硅锗超晶格结构的界面热阻。本发明适用于硅锗超晶格材料热导率研究领域，评价不同因素对界面热阻的影响。
文档编号C09K5/14GK101760183SQ20091021746
公开日2010年6月30日 申请日期2009年12月30日 优先权日2009年12月30日
发明者孙兆伟, 张兴丽 申请人:哈尔滨工业大学