D01:10.13374j.isml00103x2006.卫.040 第28卷第12期 北京科技大学学报 Vol.28 No.12 2006年12月 Journal of University of Science and Technology Beijing Dec.2006 纳米尺度孔隙内气体导热系数的分子动力学模拟 刘育松张欣欣于帆 北京科技大学机械工程学院北京100083 摘要应用LennardJone作用势，在300K和01MPa条件下，对边长20nm的立方体孔隙内氨 气的导热系数进行了平衡分子动力学模拟.结果得出分子分速度和速率的分布与统计力学得到的 Maxwel速度和速率分布曲线基本一致，并且分子的平均自由程受到孔隙壁的严格限制.通过 Green-Kubo关系式计算得出了孔隙内氮气的导热系数.并与文献中的结果进行了比较模拟结果 接近于实验值.仅为同样条件下自由空间氮气的导热系数的V3左右. 关键词氨气：导热系数：纳米尺度：孔隙：分子动力学模拟 分类号TK1240414 近年来，具有纳米或亚微米尺度结构的介孔 Lennard-Jone作用势，通过平衡MD方法对纳米 材料由于其极高的绝热性能成为传热领域的一个 孔隙内的氮气的运动行为进行模拟和分析，通过 研究热点).这些介孔材料中的传热一般有三种 Green-Kubo关系式计算了纳米孔隙内氮气的导 方式，即固相导热、气相导热和辐射传热.在常规 热系数. 介孔材料中，气相导热要远高于固相导热和辐射 1 模拟方法 传热，占材料总传热量的50%以上，甚至达到 80%.由于具有纳米或亚微米尺度结构的介孔 对于所研究的氮气，其分子间作用势 材料与常规介孔材料相比具有更高的孔隙率和更 Lennard-Jone(LJ)为： 小的孔隙尺寸，因此研究该尺度范围内的气相导 u(r）= 4 r (1) 热对于理解该类材料的传热机理极为重要. 目前，微/纳米尺度下的气相导热系数的直接 其中，G为L山作用势的尺寸参数.e为LJ作用势 测量还难以进行.随着计算机技术的发展，分子 的能量参数，r是分子i和分子j之间的距离.对 动力学模拟(molecular dy namics以下简称MD) 于氮气，o=0331nm,e/kg=37.3K,kB为玻尔 方法，在计算物理、化学、材料和化工等领域己经 兹曼常数网 得到广泛的应用，成为比较成熟的模拟实验研究 假设氮气位于一个开放性的立方体孔隙内， 手段.在热科学领域中，采用MD方法对物质 孔隙壁边长L取为20nm,孔隙的六个壁面均由 的传热性质的研究正逐步开展起来.Fernandez 一层以fc点阵形式排列的氮气分子组成.如图I 等采用平衡MD方法计算了Ar,Kr和CH4以 所示.所模拟的整个系统共包含864个氮气分 及它们的混合气体的导热系数：Takashi等采用 子，其中孔隙内分子数为256个，六个壁面内的分 ● 非平衡MD方法对液态O2、CO、CS2、Cb和Br2等 的导热系数进行了模拟计算；Volz等9应用平衡 MD方法计算了硅纳米线的导热系数：Konstantin 等通过平衡MD方法模拟计算了固态氩的导 热系数.这些模拟都取得了很好的结果，展现了 。 MD方法的潜力和良好的应用前景.本文采用 ● 收稿日期：2005-09-08修回日期：200605-12 基金项目：国家自然科学基金资助项目(N0.50276003) 作者简介：刘有松(1974一人男.博士研究生：张欣欣(1957一)， 图1模拟的立方体好隙结构示意图 男，教授，博士生导师 Fig.I Diagrammatic sketch of pore structure in simulation
纳米尺度孔隙内气体导热系数的分子动力学模拟 刘育松 张欣欣 于 帆 北京科技大学机械工程学院, 北京 100083 摘 要 应用 Lennard-Jone 作用势 , 在 300 K 和 0.1 MPa 条件下, 对边长 20 nm 的立方体孔隙内氮 气的导热系数进行了平衡分子动力学模拟.结果得出分子分速度和速率的分布与统计力学得到的 Maxwell 速度和速率分布曲线基本一致, 并且分子的平均自由程受到孔隙壁的严格限制.通过 Green-Kubo 关系式计算得出了孔隙内氮气的导热系数, 并与文献中的结果进行了比较, 模拟结果 接近于实验值, 仅为同样条件下自由空间氮气的导热系数的 1/ 3 左右. 关键词 氮气;导热系数;纳米尺度;孔隙;分子动力学模拟 分类号 TK 124;O 414 收稿日期:2005 09 08 修回日期:2006 05 12 基金项目:国家自然科学基金资助项目(No .50276003) 作者简介:刘育松(1974—), 男, 博士研究生;张欣欣(1957—), 男, 教授, 博士生导师 近年来 ,具有纳米或亚微米尺度结构的介孔 材料由于其极高的绝热性能成为传热领域的一个 研究热点[ 1] .这些介孔材料中的传热一般有三种 方式 ,即固相导热 、气相导热和辐射传热.在常规 介孔材料中 ,气相导热要远高于固相导热和辐射 传热, 占材料总传热量的 50 %以上, 甚至达到 80 %[ 2] .由于具有纳米或亚微米尺度结构的介孔 材料与常规介孔材料相比具有更高的孔隙率和更 小的孔隙尺寸, 因此研究该尺度范围内的气相导 热对于理解该类材料的传热机理极为重要. 目前 ,微/纳米尺度下的气相导热系数的直接 测量还难以进行 .随着计算机技术的发展, 分子 动力学模拟(molecular dy namics, 以下简称 M D) 方法 ,在计算物理 、化学 、材料和化工等领域已经 得到广泛的应用, 成为比较成熟的模拟实验研究 手段[ 3] .在热科学领域中, 采用 MD 方法对物质 的传热性质的研究正逐步开展起来 .Fernández 等[ 4] 采用平衡 M D 方法计算了 Ar 、Kr 和 CH4 以 及它们的混合气体的导热系数 ;Takashi 等[ 5] 采用 非平衡M D 方法对液态 O2 、CO 、CS2 、Cl2 和Br2 等 的导热系数进行了模拟计算 ;Volz 等[ 6] 应用平衡 M D 方法计算了硅纳米线的导热系数;Konstantin 等 [ 7] 通过平衡 MD 方法模拟计算了固态氩的导 热系数.这些模拟都取得了很好的结果 ,展现了 M D 方法的潜力和良好的应用前景 .本文采用 Lennard-Jone 作用势 ,通过平衡 MD 方法对纳米 孔隙内的氮气的运动行为进行模拟和分析 ,通过 Green-Kubo 关系式计算了纳米孔隙内氮气的导 热系数. 1 模拟方法 对于 所 研 究 的 氮 气, 其 分 子 间 作 用 势 Lennard-Jone(LJ)为: u(rij)=4ε σ rij 12 - σ rij 6 (1) 其中, σ为 LJ 作用势的尺寸参数, ε为 LJ 作用势 的能量参数, rij是分子 i 和分子j 之间的距离 .对 于氮气, σ=0.331 nm , ε/ kB =37.3 K , kB 为玻尔 兹曼常数[ 8] . 图 1 模拟的立方体孔隙结构示意图 Fig.1 Diagrammatic sketch of pore structure in simulation 假设氮气位于一个开放性的立方体孔隙内, 孔隙壁边长 L 取为 20 nm ,孔隙的六个壁面均由 一层以 fcc 点阵形式排列的氮气分子组成,如图 1 所示.所模拟的整个系统共包含 864 个氮气分 子,其中孔隙内分子数为 256 个,六个壁面内的分 第 28 卷 第 12 期 2006 年 12 月 北 京 科 技 大 学 学 报 Journal of University of Science and Technology Beijing Vol.28 No.12 Dec.2006 DOI :10.13374/j .issn1001 -053x.2006.12.040
Vol.28 No.12 刘育松等：纳米尺度孔隙内气体导热系数的分子动力学摸拟 。1183。 子数为608个.在分子动力学模拟的过程中，壁 0.1MPa条件下进行的，采用Vedet算法的速度 面内的氮气分子均保持固定不动并且和孔隙内 形式来求解分子的牛顿运动方程，L山作用势的截 的氮气分子相互之间具有作用势.由于所研究的 断半径取为10o.模拟首先在NVT系统中进行， 对象是孔隙内的氮气分子，所以在X,Y和Z三 即模拟系统中的分子数(N)、体积(V)和温度 个方向上采用周期性边界条件， (T)要保持不变用以确定各部分能量及热流. 采用Green-Kubo关系式来计算孔隙内氮气 各个分子初速度设定为零，然后通过一调温系数 的导热系数，即： 来自动调节温度，时间步长△t由△tNe/m/6= 1 k=3WBTJ。<J0J》d (2) 00001计算.当NVT系统达到平衡，即在一定 的模拟时间内，温度T保持恒定，则启动NVE系 其中，V是模拟系统的体积T为系统的温度，符 统下的平衡MD模拟.在NVE系统下，模拟系统 号()表示的是在整个模拟时间内系统的平均值， 中的分子数(N)、体积(V)和能量(E)都保持不 J(t)是热流 变.图2为NVT系统达到平衡后，各个分子在不 =空+片是，5》 (3) 同的分速度和速率区间分子数的分布.'m为分 子的平均速率.从图中可以看出，通过平衡MD 其中，：是分子i的速度矢量，F是分子i由于 方法模拟得到的分子分速度和速率的分布与统计 L山作用势的作用受到邻近分子j的作用力.e 力学得到的Maxwell速度和速率分布曲线基本一 可以通过下式给定： 致.从模拟结果计算得出，此时分子的平均速率 (4) 为460ms1.通过对分子在不同时间的位置和 其中，m是氮气分子质量. 速度的模拟计算可知，分子的平均自由程受立方 体的边长严格限制，远小于自由空间分子的平均 2 模拟结果与讨论 自由程. 图3为NVE系统的能量变化曲线.可以看 平衡MD模拟计算是在温度300K和压力 35[( 30 30l (b) ■模拟结果 25 ■模拟结果 -Maxwe分布 -Maxwell分布 20 ■ 15 15 ◆◆ ■ 10 10 ◆ 5 ■ 要■ -0.6 -0.2020.61.0 86 020.40.60.81.01.21.4 (瑞) 图2氨气分子的分速度和速率的分布.()氨气分子沿x方向速度分布：(b)氮气分子的速率分布 Fig.2 Distributions of the velocity component and speed of nitrogen molecues:(a velocity component distribution in x direction:(b) speed distribution 出，在经过1.4ps后，系统各部分能量达到平衡. Zeng等9依据理想气体模型和动力学理论 J(t)以及时间相关函数的变化曲线分别如图4 推导出了多孔介质中气相导热系数： 和5所示.在此基础上，采用式(2)计算得出孔隙 kg= 内氮气的导热系数为6.8×103W·(mK)1. Q225r-1.25)0.461mg (P/knT)(8kBT/m)2 MD模拟结果仅为同样条件下自由空间的氮气导 0.255.pr (P/kg T d 热系数(2.4X102W(mK)-1)的1/3左右，其 (5) 原因主要有：一是由于孔隙壁的存在，使得分子平 其中，S,是固相的比表面积Ppr为密度，中为气 均自由程减小限制了分子输运能量的能力：二是 孔率，P为压力，mg为分子质量，dg为分子直径 由于分子之间以及分子与孔隙壁分子之间的相互 =cp/cm,cp和cr分别为气体比定压热容和比 作用势，使分子在孔隙壁上的漫反射得到加强. 定容热容.为了比较MD模拟结果、式(5)的计算
子数为 608 个 .在分子动力学模拟的过程中, 壁 面内的氮气分子均保持固定不动, 并且和孔隙内 的氮气分子相互之间具有作用势 .由于所研究的 对象是孔隙内的氮气分子, 所以在 X , Y 和Z 三 个方向上采用周期性边界条件 . 采用Green-Kubo 关系式来计算孔隙内氮气 的导热系数[ 8] ,即 : k = 1 3 Vk B T 2∫ ∞ 0 dt (2) 其中 , V 是模拟系统的体积, T 为系统的温度 ,符 号〈〉表示的是在整个模拟时间内系统的平均值, J(t)是热流 . J(t)= ∑i viεi +1 2 i ,∑ j , i ≠j rij(Fi jvi) (3) 其中, vi 是分子 i 的速度矢量, Fij 是分子 i 由于 LJ 作用势的作用受到邻近分子 j 的作用力 .εi 可以通过下式给定: εi = 1 2 m vi 2 + 1 2 ∑i ≠j u(rij) (4) 其中 , m 是氮气分子质量. 2 模拟结果与讨论 平衡 MD 模拟计算是在温度 300 K 和压力 0.1 MPa 条件下进行的 , 采用 Verlet 算法的速度 形式来求解分子的牛顿运动方程, LJ 作用势的截 断半径取为 10σ.模拟首先在 NV T 系统中进行, 即模拟系统中的分子数(N )、体积(V)和温度 (T)要保持不变, 用以确定各部分能量及热流. 各个分子初速度设定为零 ,然后通过一调温系数 来自动调节温度,时间步长 Δt 由 Δt ε/ m/ σ= 0.000 1计算.当 NV T 系统达到平衡 , 即在一定 的模拟时间内 ,温度 T 保持恒定 ,则启动 NVE 系 统下的平衡M D 模拟 .在 NVE 系统下 ,模拟系统 中的分子数(N)、体积(V)和能量(E)都保持不 变.图 2 为 NVT 系统达到平衡后 ,各个分子在不 同的分速度和速率区间分子数的分布 .v m 为分 子的平均速率.从图中可以看出 , 通过平衡 M D 方法模拟得到的分子分速度和速率的分布与统计 力学得到的 Maxwell 速度和速率分布曲线基本一 致.从模拟结果计算得出, 此时分子的平均速率 为460 m·s -1 .通过对分子在不同时间的位置和 速度的模拟计算可知 , 分子的平均自由程受立方 体的边长严格限制 , 远小于自由空间分子的平均 自由程. 图3为NVE系统的能量变化曲线 .可以看 图 2 氮气分子的分速度和速率的分布.(a)氮气分子沿 x 方向速度分布;(b)氮气分子的速率分布 Fig.2 Distributions of the velocity component and speed of nitrogen molecul es:(a)velocity component distribution in x direction;(b) speed distribution 出,在经过 1.4 ps 后 ,系统各部分能量达到平衡. J (t)以及时间相关函数的变化曲线分别如图 4 和5 所示.在此基础上, 采用式(2)计算得出孔隙 内氮气的导热系数为 6.8 ×10 -3 W·(m·K)-1 . M D 模拟结果仅为同样条件下自由空间的氮气导 热系数(2.4 ×10 -2 W·(m·K)-1)的 1/3 左右 ,其 原因主要有:一是由于孔隙壁的存在,使得分子平 均自由程减小, 限制了分子输运能量的能力;二是 由于分子之间以及分子与孔隙壁分子之间的相互 作用势,使分子在孔隙壁上的漫反射得到加强 . Zeng 等[ 9] 依据理想气体模型和动力学理论 推导出了多孔介质中气相导热系数 : k g = (2.25 r -1.25)0.461 mg(P/ k B T)(8 k B T/πm g) 1/2 0.25S s ρpo r -1 + 2(P/ k B T)πd 2 g c V (5) 其中 , S s 是固相的比表面积, ρpor为密度, 为气 孔率 , P 为压力, mg 为分子质量, dg 为分子直径, r =cp/ cV , cp 和cV 分别为气体比定压热容和比 定容热容 .为了比较 MD 模拟结果 、式(5)的计算 Vol.28 No.12 刘育松等:纳米尺度孔隙内气体导热系数的分子动力学模拟 · 1183 ·
。1184。 北京科技大学学报 2006年第12期 12.068 (a) 12.066 12.064 12.062 2 12.060 12.058 12.056 0 12.05%62002040.60810121416 0 0.51012.02.5 时间ps 时间ps 0.008 (b) 图5J八)时间相关函数变化曲线 0.006 Fig 5 Curve of the time correlation function of J(t) 0.004 结果以及实验测量结果，以一种硅气凝胶材料为 0.002 例，其主要参数为：S,=797m2g1,P=110 0 kgm3,=94%,平均孔隙尺寸20m.在温度 -0.002 为300K和压力0.1MPa下，采用热线法间接测 -0.04%20.02040.60.81012141.6 量得到气相导热系数为8X103W·(mK)1 时间s (孔隙内气体为空气)0.在同样条件下，通过式 12.06292 (c) (5)求解得到该硅气凝胶材料气相导热系数为 12.06290 9X103W(mK)1.平衡MD模拟结果为 12.06288 68X103W·(mK)一1，和实验测量结果比较接 12.06286 近，但小于式(5)的计算结果.这是由于Zemg采 酒 12.06284 用的是理想气体模型，忽略了分子之间的相互作 12.06282 用势 12.06280 -0.200.20.40.60.81.012141.6 3结论 时间ps 采用了Lennard--Jone作用势，通过平衡MD 图3NVE系统的能量变化曲线.(a)动能变化曲线；(b)势 方法在温度300K和压力0.1MPa条件下对边长 能变化曲线：（©）总能量变化曲线 Fig.3 Relationship curves of energy and time in the NVE:(a) 为20nm的立方体孔隙内的氮气的运动行为进行 kinetic energy curve (b)potential energy curve;(c)total ener- 了模拟，得出了分子分速度和速率的分布，与统计 gy curve 力学得到的Maxwell速度和速率分布曲线基本一 致，由此计算出氮气分子平均速率为460ms1. 12 通过Gen-Kubo关系式计算得出了立方体孔隙 1.0 内氮气导热系数为68×103W·(mK)-1,仅为 同样条件下自由空间氮气的导热系数的1/3 Iu-OIKI 0.6 左右. 0.2 参考文献 1]Cahill D G.Fond W K.Goodson K E.et al.Nanoscale ther 00.20.40.60.81.012141.6 时间/ps mal transport.J Appl Phys 2003.93(2):793 I24 徐烈。方荣生。马庆芳.绝热技术。北京：国防工业出版社， 图4J()变化曲线 1990 Fig.4 Curve of J(t) [3 Heermann D W.Computer simulation methods in theoretical physics.London:Springer-Verlag-Bedin Heidelberg.1986
图 3 NVE 系统的能量变化曲线.(a)动能变化曲线;(b)势 能变化曲线;(c)总能量变化曲线 Fig.3 Relationship curves of energy and time in the NVE:(a) kinetic energy curve;(b)potential energy curve ;(c)total energy curve 图 4 J(t)变化曲线 Fig.4 Curve of J(t) 图 5 J(t)时间相关函数变化曲线 Fig.5 Curve of the time correlation function of J(t) 结果以及实验测量结果 ,以一种硅气凝胶材料为 例,其主要参数为:Ss =797 m 2 ·g -1 , ρpor =110 kg·m -3 , =94 %, 平均孔隙尺寸 20 nm .在温度 为 300 K 和压力 0.1 M Pa 下 ,采用热线法间接测 量得到气相导热系数为 8 ×10 -3 W ·(m·K) -1 (孔隙内气体为空气)[ 10] .在同样条件下, 通过式 (5)求解得到该硅气凝胶材料气相导热系数为 9 ×10 -3 W·(m·K) -1 .平衡 MD 模拟 结 果为 6.8 ×10 -3 W·(m·K)-1 ,和实验测量结果比较接 近,但小于式(5)的计算结果 .这是由于 Zeng 采 用的是理想气体模型 , 忽略了分子之间的相互作 用势 . 3 结论 采用了 Lennard-Jone 作用势 , 通过平衡 M D 方法在温度 300 K 和压力 0.1 MPa 条件下对边长 为 20 nm 的立方体孔隙内的氮气的运动行为进行 了模拟,得出了分子分速度和速率的分布,与统计 力学得到的 Maxwell 速度和速率分布曲线基本一 致,由此计算出氮气分子平均速率为460 m·s -1 . 通过 Green-Kubo 关系式计算得出了立方体孔隙 内氮气导热系数为 6.8 ×10 -3W·(m·K) -1 ,仅为 同样条件下自由空间氮气的导热系数的 1/3 左右 . 参 考 文 献 [ 1] Cahill D G , Ford W K , Goodson K E , et al.Nanoscale thermal transport .J Appl Phys, 2003 , 93(2):793 [ 2] 徐烈, 方荣生, 马庆芳.绝热技术.北京:国防工业出版社, 1990 [ 3] Heermann D W .Comput er simulation methods in theoretical physics.London:Sp ringer-Verlag-Berlin Heidelberg , 1986 · 1184 · 北 京 科 技 大 学 学 报 2006 年第 12 期
Vol.28 No.12 刘育松等：纳米尺度孔隙内气体导热系数的分子动力学摸拟 。1185。 [4]Femndez G A.Vmbec J.Hasse H.A molecular simulation gon from mokcular dynamics simulations.J Chem Phys. study of shear and bulk viscosity and thermal conductivity of 2004,120(8):3765 simple real fhids.Flud Phase Equilib 2004 221:157 [8 Allen M P.Tildesey D J.Computer simulation of liquids. [5]Takashi T.Kenjm K.Molecular dynamics study for the ther Oxford:Clarendon Press 1987 mal conductivity of diatomic liqid.Superlattices Microstruct. [9 Zeng S Q.Hunt AJ.Greif R.Mean free path and apparent 200435:217 themal conductivity of a gas in a porous medium.ASME J [6]Volz S.Chen G.Mokcular dynamics simulation of themal Heat Transfer.1995.117:758 conductivity of silicon nanow ires.Appl Phys Lett,1999,75: 10]Zeng S Q Hunt A J,Greif R.Geometric structure and ther 2651 mal conductivity of porous medium silica aerogd.ASME J [7]Tretiakow K V.Scandolo S.Thermal conductivity of solid ar Heat Transfer.1995,117:1055 Molecular dynamics simulation of thermal conductivity of a gas in nanoscale pores LIU Yusong,ZHANG Xinxin,YU Fan Mechanical Engineering School University of Science and Technology Beijing.Beijing 100083 China ABSTRACT Using Lennard-Jones pair potential the thermal conductivity of nitrogen in a cubic pore of 20 nm in length was calculated by equilibrium molecular dynamics (M D)simulation at 300K and 0.IM Pa. The distributions of velocity component and speed of molecules calculated by the simulation and the Maxwell distribution equation fit very well.The mean free path of the molecules was strictly confined by the wall size of the nanoscale pore.The thermal conductivity of nitrogen in the nanoscale pore was deter- mined by the GreemKubo method and separately compared with the experimental and theo retic data in the literatures.The simulation result is close to the experimental one and only about one third of the value in the free space under the same conductions. KEY WORDS nitrogen;thermal conductivity;nanoscale;pore;molecular dynamics simulation
[ 4] Ferná ndez G A , Vrabec J, Hasse H .A molecular simulation study of shear and bulk viscosity and thermal conductivity of simple real fluids.Fluid Phase Equilib, 2004 , 221:157 [ 5] Takashi T , Kenjiro K .Molecular dynamics study f or the thermal conductivity of diat omic liquid.Superlattices Microstruct , 2004 , 35:217 [ 6] Volz S , Chen G .Molecular dynamics simulation of thermal conductivit y of silicon nanow ires.Appl Phys Lett , 1999 , 75: 2651 [ 7] Tretiakov K V , S candolo S.Thermal conductivity of solid argon from molecular dynamics simulations.J Chem Phys, 2004 , 120(8):3765 [ 8] Allen M P, Tildesley D J .Computer simulation of liquids. Oxford:Clarendon Press, 1987 [ 9] Zeng S Q , Hunt A J, Greif R.Mean free path and apparent thermal conductivity of a gas in a porous medium .ASME J Heat Transfer, 1995 , 117:758 [ 10] Zeng S Q, Hunt A J, Greif R.Geometric structure and thermal conductivit y of porous medium silica aerogel.ASME J Heat Transfer, 1995 , 117:1055 Molecular dynamics simulation of thermal conductivity of a gas in nanoscale pores LIU Y usong , ZHANG X inxin , Y U Fan Mechanical Engineering School, University of Science and Technology Beijing , Beijing 100083 , China ABSTRACT Using Lennard-Jones pair potential, the thermal conductivity of nitrogen in a cubic po re of 20 nm in leng th w as calculated by equilibrium molecular dynamics (M D)simulatio n at 300K and 0.1M Pa. The distributions of velocity component and speed of molecules calculated by the simulation and the M axwell distributio n equatio n fit very well.The mean free path of the molecules w as strictly confined by the w all size of the nanoscale pore.The thermal conductivity of nitrogen in the nanoscale pore w as determined by the Green-Kubo method and separately compared with the experimental and theo retic data in the literatures .The simulation result is close to the experimental one and only about one third of the value in the free space under the same conductions . KEY WORDS nitrogen ;thermal conductivity ;nanoscale ;pore ;molecular dynamics simulation Vol.28 No.12 刘育松等:纳米尺度孔隙内气体导热系数的分子动力学模拟 · 1185 ·
