【LAMMPS如何系列】计算热膨胀系数

热膨胀系数(thermal expansion coefficent)是物质因温度改变时,体积发生变化的趋势。

维基百科的定义:Thermal expansion is the tendency of matter to change in volume in response to a change in temperature.  (Link:thermal expansion)

热膨胀系数有体膨胀系数 β=ΔV/(V*ΔT) 和 线膨胀系数 α=ΔL/(L*ΔT) 。

从这个定义大致就可以获得计算热膨胀系数的思路:计算不同温度下物质的体积或长度,那么体积或长度随着温度的变化曲线的斜率就是ΔL/ΔT了。所以,在模拟的过程中,体积是变化的,那么自然就会想到使用NPT系综。在使用NPT系综之前和之后,建议使用NVT系综进行平衡一定的步数。

举例计算热膨胀系数

下面仍然以铜为例,具体介绍如何使用LAMMPS计算线膨胀系数α。

输入脚本:in.thermal.expansion

# This input script is used to calculate
# the thermal expansion of copper.
# Powered by Xianbao Duan
# Email: xianbao.d@gmail.com
# Website: http://www.52souji.net/

units           metal
boundary        p p p
atom_style      atomic

variable	i loop 10
variable        x equal 200+100*$i

lattice		fcc 3.62
region		box block 0 8 0 8 0 8
create_box	1 box
create_atoms 	1 box

pair_style      eam
pair_coeff      1 1 Cu_u3.eam

variable        N equal step
variable        pote equal pe
variable        Etotal equal etotal
variable        T equal temp
variable        Press equal press
variable        V equal vol

velocity        all create 2.5 825577 dist gaussian

timestep	0.002
thermo		1000

fix		1 all nvt temp 2.5 2.5 0.2
run		2000

unfix		1

fix             2 all npt temp 2.5 $x 0.2 iso 0 0 10
run             120000

unfix		2

fix             extra all print 100 "${N} ${T} ${V} ${pote} ${Etotal} ${Press}" append data

fix		3 all nvt temp $x $x 0.2
run		2000

unfix		3
unfix		extra

clear
next		i
jump		in.thermal.expansion

对以上脚本进行简单解释:

1. 第12行定义的变量x为温度,即计算这些温度下体系的晶格常数。

2. 第34-35行将体系在2.5K下进行NVT平衡,为NPT升温过程做准备。

3. 第39-40行将体系在NPT系综下进行升温。

4. 第44行定义输出有关变量。注意这里使用的是append关键字,而不是file,是为了让脚本在循环执行时不重写这个文件。

5. 第46-47行将升温后的体系在NVT下进行适当平衡,在这个过程中输出了有关变量到文件中。

得到的data文件即为结果数据。将其中的温度和体积复制出来列表,绘图,如下图。

thermal expansion coefficient - temperature vs lattice constant

图中曲线的斜率是6.89444e-5,即ΔL/ΔT= 6.89444e-5 A/K。又初始的晶格常数大致为3.6346A,所以可以算出线膨胀系数α=ΔL/(L*ΔT) =18.97 10-6 K-1

实验结果为17.5 10-6 K-1,两者基本吻合。

材料的热膨胀系数可以参考:热膨胀系数查询表

标签: lammps, input file

相关文章推荐

添加新评论 (无需注册,可直接评论)

已有 28 条评论

  1. 11

    data文件里面没有数据

  2. 可否加一下联系方式,诚心向大神请教lammps相关问题

  3. 我用lammps跑了一下,为什么结果跟博主的不一致,向博主请教

  4. 没有晶格常数怎么办

    1. 非晶体物质怎么算呢

      1. 我没有做过非晶物质的建模,可能需要使用一些软件或者编程来产生。

  5. guang

    博主你好,问下为什么到最后出现错误:cannot open input script in.thermal.expansion

    1. 是不是你把脚本的文件名改了? 找不到文件了,就会出现这个提示。

  6. 我套用你的模版去做其他元素的,可是结果与实际值差好多。是什么原因

  7. mjj

    请问博主,体积和晶格常数如何换算?多谢了

  8. zhongyanan

    请问可以由温度和体积计算体膨胀系数吗

    1. zhongyanan

      博主您好,我用你的程序算了一下金元素的,得到了十个温度和十个体积,初始温度是273K,得到的体积和温度的斜率是44.145,初始的体积是35319,那斜率和初始体积相乘得到的数据也太大了吧,这个体积单位是不是要换算一下呢?请博主指导一下

  9. minlei

    请问如何用得到的温度和体积,来得到下图的温度和晶格常数的关系曲线,

    1. 我爱搜集网博主

      将体积换算成晶格常数就可以啦~

  10. 冯闯

    用最新版的LAMMPS跑了下你的程序,跑一个循环后显示错误“unbalanced quotes in input file",这个什么原因?

    1. 冯闯

      少输入了一个引号,见笑了。

  11. charles

    代码非常简洁明了。想请教博主几个问题:

    图中的纵坐标晶格常数是输出data文件中的那个参数?是体积吗? 这个例子当中的热胀系数跟温度没有关系,如果有的材料热胀系数跟温度有关系,该怎么处理呢?
    期待您的回复。多谢了。
    1. 我爱搜集网博主
      可以通过体积换算,也可以通过修改脚本直接输出晶格常数 可以考虑将温度控制在一个较小的范围内
      1. charles

        多谢,我试试。

  12. chy

    您好,我现在在学习计算弹性模量,不知道用lammps计算,请教一下请您。

    1. 我爱搜集网博主

      请问你想要计算什么模量?
      如果是弹性常数,请参考本文。
      如果是体积模量请参考
      http://www.52souji.net/how-to-calculate-bulk-modulus-using-lammps/

      1. chy

        谢谢,我是想计算杨氏模量,和剪切模量,但是一直出问题。

  13. 付佩

    脚本写的简洁明了,但是我想问下,我运行之后,改变npt的步数,发现温度下对应的体积会发生变化,然后热膨胀系数也会变化,请问那么如何确定正确的npt步数呢?

  14. 高士才

    为什么没写质量呢?

    1. 我爱搜集网博主

      EAM势函数不需要设置质量的。
      the EAM potential files list atomic masses; thus you do not need to use the mass command to specify them.
      link: http://lammps.sandia.gov/doc/pair_eam.html

  15. 高士才

    脚本写的简洁好看,但是为什么没写质量呢?

  16. donghaipark

    为什么NVT和NPT的温度是2.5K?

    1. 我爱搜集网博主

      这个温度差不多是随便给的,是为了将体系从一个比较低的温度开始升温。