摘 要: 采用非平衡态分子动力学( NEMD) 仿真了近场热记录系统的热响应过程. 对探针施加不同的作用力时, 随着作用力的增加, 基体的热响应明显加快, 系统达到热平衡所需时间相应减少, 可有效提高热机械数据存储速率. 采用数值差分方法对一维探针- 基体的传热模型进行数值求解, 数值解和分子动力学模拟结果总体趋势基本吻合.
关键词: 热响应; 非平衡态; 分子动力学; 数值差分
Gerd 等人于1985 年研制出了原子力显微镜Atomic Fo rce Micro scope ( AFM) [ 1] , AFM 主要应用于微观领域的成像技术和加工, 在此基础上, BinNIng等提出基于AFM 的热机械数据存储技术, 用于数据的存储、读写以及擦除, 针尖直径在50 纳米时, 数据存储密度可望达到100GB/ in2 . 热机械数据存储系统有可能成为令人满意的数据存储技术[ 2-5] .热机械写过程是在聚合物上施加一个局部力, 并加热探针针尖来软化聚合物层, 从而在基底上写下一个数据位, 如图1[ 5] 所示. 在探针加热过程中, 局部高温的探针针尖与基体构成了典型的一维热传导过程, 当针尖直径在纳米量级时, 所构成的微尺度热传导模型对了解声子输运、界面热阻机理有重要理论研究价值.
Peterson R. B. [ 6] 从德拜晶体模型中获得了基本声子频率散射以及采用很简单的声子碰撞模拟,模拟了连续有限体上的一维温度曲线. Ho- Ki Lyeo等人[ 7] 用超高真空扫描式热电显微镜研究半导体纳米结构的局部热电能, 揭示出当应用一个p-n 结时,在距结2 nm 内热电能图突然改变方向. Maruyama和Kimur a[ 8] 对固液界面热阻的计算表明, 固液界面热阻大约相当于5~ 20 nm 厚度的液体的热阻, 并随着液体对固体润湿能力的增加而迅速减少.Lukes 等[ 9] 采用NEMD 方法对薄膜的一半热传导特性进行了仿真, 仿真结果表明薄膜的热传导性能随薄膜厚度的增加而增加, 并最终可超过体态( bulk) 实验值30%左右. Li Shi、Arun Majumdar 等人[ 10-11] 通过实验方法, 证实当热源尺寸减小时, 固体和液体之间的传导变得更为重要, 液体薄膜的存在可以提高热传导. 然而, 由于针尖直径极其微小, 它和聚合物薄膜之间的热传导是很弱的, 影响了它的存储速度. 可以通过提高两者的作用时间或对针尖施加作用力来提高到聚合物的热传导, 但是如果针尖作用力太大或者作用时间过长就会导致聚合物上的数据点的尺寸过大, 影响存储密度. 所以如何准确的控制针尖作用力及作用时间是热机械写技术的一个重要问题.
为了进一步研究这个问题, 本文采用分子动力学Mo lecular dy namics( MD) 方法, 模拟针尖通过液体薄膜向基底传热的一维过程, 得出系统中各个位置的温度响应. 在此基础上, 对针尖施加作用力, 比较在不同的力作用下的温度响应, 讨论针尖施加的作用力的大小以及作用时间对传热过程的影响. 在纳米尺度范围内, 作为实验方法的一个有效补充手段, MD 方法已被广泛应用于热传导的研究工作中.MD 方法主要分为基于Green-Kubo 关系的平衡态分子动力学( EMD) 和基于Fourier 定律的非平衡态分子动力学( NEMD) 方法[ 12] . 本文采用的是NEMD 方法. 最后用数值差分方法[ 13 ] 求针尖传热的温度响应的数值解, 和分子模拟所得的解作比较.
1 分子动力学模型
本文选取Ar 作为针尖和基底的材料, 采用Lennard-Jones ( LJ) 两体作用势描述原子间作用力[ 1 4] , 即
NEMD 的分子动力学模型如图2 所示, 粒子按照面心立方( FCC) 晶格排列, 选用恒温墙[ 15] 的方法来施加温度梯度. 系统初期平衡完成后, 在保持热( 冷) 域内粒子净动量不变的基础上, 通过调整热( 冷) 域内粒子的动量使热( 冷) 域在模拟过程中恒定在给定的温度TH ( TC ) [ 16] 上( T H , TC 为热域、冷域的温度, TH > T C) . 绝热壁粒子在模拟过程中静止在各自平衡位置, 并通过LJ 势同其他粒子相互作用.绝热壁的引入降低了热( 冷) 域表面粒子受力的不均衡程度, 有效地防止了在模拟过程中表面粒子的蒸发. 由于原子间作用力是短程的( r Frc ) , 因此, 选用L 4 为2 个晶胞( unit cell , UC) 长度来模拟无穷厚度的绝热壁. 为了避免绝热壁的引入对模拟结果造成任何影响, 必须将绝热壁的边界效应限制在热( 冷)域中, 即热( 冷) 域x 向厚度L 3 必须保证绝热壁粒子同针尖和基底粒子不能直接发生相互作用, 为此选用L3 = 2 UC. 另外, 选取基体和针尖的x 向厚度L 2= 4 U C, 薄膜的x 向厚度粒子L 1= 1 U C. 对基体的y , z 两个方向上施加周期性边界条件, 对针尖不施加周期性边界条件.
采用Velocity Verlet[ 17] 算法计算粒子运动方程的数值积分, 积分步长为0. 1 fs, 总模拟时间为0. 2 ns( 200 万步) . 由于Ar 的熔点为83. 7 K, 为了得到稳定的固体结构, 取热域和冷域温度分别为60 K 和20 K,模拟的平均温度为40 K.
2 一维热传导模型
直角坐标系下导热微分方程为:
系统的初始温度为40 K, 基体左端的温度取20K, 针尖右端的温度取60 K. 由于本文的目的是探讨针尖传热的普遍特性, 而不是求解某种具体的针尖基体材料的热传导, 所以并没有选取物体实际的热扩散率, 而只是选一个比例, 这里令A1 = 10A2 = 5A3 ,x 1 : x 2 : x 3= 4: 1: 4
3 分子模拟结果与讨论
采用基于经典Bolt zmann 分布的能量均分定理, 可以得到垂直于X 轴的粒子平面的温度, 即
其中, mA 为粒子质量, N p 为该平面中粒子数量, v i 为第个粒子的速度, kB 为Bo lt zman 常数. 公式成立的条件为: ①每个X 平面可以达到区域热平衡; ②模拟温度高于材料的Deby e 温度. 所以这里选取的模拟的横截面为5 U C x 5 U C, 而Ar 的Debye温度为92 K, 高于我们模拟的平均温度( 40 K) .X 向共有17 个晶胞, 因为是FCC 结构, 所以可以将有粒子的位置在垂直于X 轴的方向分成34 个平面, 每隔1000 步记录下各个平面的温度.
图4 显示了一维热传导的温度随时间的响应图, 分别选取了2. 5 万步( 2. 5 ps) , 5 万步( 5 ps) , 20万步( 20 ps) , 75 万步( 75 ps) , 150 万步( 1. 5 ns) 时的温度分布图, 曲线反映了从传热开始到传热逐渐稳定过程的温度变化. 记录基体、针尖、薄膜的温度,共18 个平面, 其中最前面的8 个平面为基体, 中间的两个平面为薄膜, 最右侧的8 个平面为针尖, 曲线上每个点代表一个平面的温度. 图4( a) 中没有对针尖施加作用力, 在2. 5 万步和5 万步处薄膜附近位置温度大约为系统的初始温度40 K, 只有冷热域附近的粒子才开始传热, 到75 万步处系统基本趋于平衡. 对于基体和针尖两段曲线跟Peterson R. B. 得到的连续有限体上的一维的温度曲线[ 6] 相吻合. 在系统达到平衡以后, 由于界面热阻的影响, 在薄膜处有较大的温度梯度, 跟Maruyama 等人[ 8] 所得到的固液界面热阻较大的结论相一致.
图4( b) 中, 对针尖施加了0. 1 个UC 的作用力, 即针尖插入基体0. 1 个UC, 将图4( b) 和( a) 作比较, 可以明显地看出在2. 5 万步和5 万步的两条曲线中, ( b) 图中曲线各个点的温度比( a) 图中的要高, 到系统趋于平衡后, 两个图的曲线差别很小. 图4( c) 为对针尖作用0. 2 个UC 力, 前两条曲线的温度要比a、b 图中的要高. 通过对图4 中3 个图的比较, 可以得出随着对针尖的作用力的增加, 针尖对基体的传热越快, 系统达到稳定所需的时间也相应的减少. 所以对针尖施加适当的作用力, 可以提高数据存储的速率, 且不影响其存储密度.
4 一维热传导数值解
用数值差分方法求方程的数值解, 求得解如图5, 和Peterso n R. B. 得到的连续有限体上的一维的温度曲线[ 6] 相吻合. 五条曲线的时间分别为t= 100,t= 500, t= 1000, t= 5000, t= 20000. 从结果可以看出, 数值解的结果跟用分子动力学模拟出来的结果的变化趋势是很接近的. 由于在数值解中没有考虑界面热阻的作用, 当系统稳定时, 温度成线性分布,与模拟结果有所差异.
5 结 论
采用NEMD 算法模拟了探针针尖通过聚合物薄膜向基底传热的一维过程, 得到了系统各个位置的温度响应, 当针尖没有受到外界作用力时, 基体热响应相对缓慢, 在薄膜位置界面热阻较大, 温度梯度较大. 随着作用力的增加, 针尖对基体的热传导明显增强, 有利于针尖对基底的热机械写数据, 系统温度平衡所需时间也相应减少, 增大数据存储速率. 用数值差分方法对一维的针尖传热所得到的数值解和分子动力学模拟结果的趋势基本上吻合.
参考文献:
[ 1] Binning G, Quat e C. F, Gerber Ch. Physical Review Lett ers[ J] . 1982. 56: 930-933.
[ 2] Mamin H. J , Rugar D. Applied Phy sics Let ters [ J] . 1992. 61:1003- 1005.
[ 3] Binnin g G, Despont M, Drech sler U, et al . Applied Ph ysicsLet ters[ J ] . 1999. 74: 1329-1331.
[ 4 ] Lutwy che M. I, Des pont M, Drech sler. U , et al. AppliedPhysi cs Let t ers[ J] . 2000. 77: 3299-3301.
[ 5] Wil liam P. King, T homas W. Kenny, Kenn eth E. Goods on ,et al. Applied Physi cs Let t ers [ J] . 2001. 78: 1300-1302.
[ 6] Peterson R B. Journal of Heat Transfer [ J] . 1994. 116: 815- 822.
[ 7] H o-Ki Lyeo, Khajetoorians, Li S hi, et al. Science[ J] . 2004.303: 816- 818.
[ 8] Maruyama S , Kimura T . T herm al Science& En gineering [ J ] .1999, 7: 63- 68.
[ 9] Luk es J R, Li D Y, Liang X G, et al. Journal of Heat Transfer[ J] . 2000, 122: 536-543
[ 10] Li Shi. Arunava M ajum dar. Journal of H eat Transf er [ J ] .2002. 124: 329-337.
[ 11] Ohmyoung Kw on, Li Sh i. Aru n Majumdar. J ournal of H eatT ran sf er[ J] . 2003. 125: 156- 163.
[ 12] 冯晓利, 李志信, 过增元. 工程热物理学报[ J] . 2001. 22: 195- 198.
[ 13] 孙志忠, 袁慰平, 闻震初. 数值分析[ M] ( 第二版) , 南京: 东南大学出版社, 2002. 318-346.
[ 14] All en M P, T ildesley D J. Compu ter simul at ion of liquids[ M] . Oxford: Claren don Press, 1987. 9- 50.
[ 15] Chou F C, Luk es J R, Liang X G, et al. Annual Review ofH eat T ransf er[ J ] . 1999, 10: 141-176.
[ 16] M ou nt ain R D, MacDonald R A. Phys ical Review B [ J ] .1983, 28: 3022-3025.
[ 17] Rapaport D C. T he art of m ol ecu lar dynam ics sim ulati on[ M] . Cambridge: Universit y Press, 1995: 50-118.
