【LAMMPS学习】八、基础知识(3.7)计算粘度

8. 基础知识

此部分描述了如何使用 LAMMPS 为用户和开发人员执行各种任务。术语表页面还列出了 MD 术语,以及相应 LAMMPS 手册页的链接。 LAMMPS 源代码分发的 examples 目录中包含的示例输入脚本以及示例脚本页面上突出显示的示例输入脚本还展示了如何设置和运行各种模拟。

8.1.通用基础知识

8.2. 设置入门

8.3. 分析入门

8.3.1. LAMMPS 的输出(thermo, dumps, computes, fixes, variables

8.3.2.使用chunks计算系统属性 

8.3.3.使用分布式网格 

8.3.4.计算温度

8.3.5.计算弹性常数

8.3.6.计算热导率

8.3.7.计算粘度

可以使用 LAMMPS 中的各种选项以至少 6 种方式测量流体的剪切粘度 eta。请参阅examples/VISCOSITY 目录,了解实现此处讨论的简单 Lennard-Jones 流体模型的 5 种方法和 SPC/E 水模型的 1 种方法的脚本。另请参阅计算热导率( calculating thermal conductivity)页面,了解热导率的类似讨论。

Eta 是流体在与速度或动量流方向垂直的方向上传递动量的倾向的量度。或者,它是流体必须被剪切的阻力。它由下式给出

J = -eta grad(Vstream)

其中 J 是动量通量,单位为每面积每时间的动量。 grad(Vstream) 是沿另一个方向移动的流体速度的空间梯度,该方向垂直于动量流过的区域。因此,粘度具有压力-时间单位。


第一种方法是通过 fix deform命令剪切模拟框来执行非平衡 MD (NEMD) 模拟,并使用 fix nvt/sllod 命令通过 SLLOD 运动方程对流体进行恒温。


或者,作为第二种方法,可以使用一个或多个移动壁来剪切它们之间的流体,同样使用某种恒温器,仅修改速度的热(非剪切)分量以防止流体升温。

在这两种情况下,通过此过程在流体中设置的速度分布可以通过 fix ave/chunk 命令进行监控,该命令确定上面等式中的 grad(Vstream)。例如。流体运动 Vx 分量在 y 方向的导数或 grad(Vstream) = dVx/dy。还可以监视由 compute pressure 命令计算得出的压力或应力张量的 Pxy 非对角分量,它是上式中的 J 项。有关 NEMD 模拟的详细信息,请参阅Howto nemd ​​​​​​​ 页面。


第三种方法是使用 fix viscosity命令执行反向非平衡MD模拟,该命令实现了Muller-Plathe的rNEMD算法。一维的动量在不同维度的模拟盒的两个不同层中的原子之间交换。这会产生速度梯度,可以使用  fix ave/chunk ​​​​​​​命令进行监控。该修复会计算其执行的累积动量传递。有关详细信息,请参阅 fix viscosity命令。


第四种方法基于 Green-Kubo (GK) 公式,该公式将应力/压力张量自相关的整体平均值与 eta 联系起来。这可以在完全平衡的模拟中完成,这与前面两种非平衡方法相反,在非平衡方法中动量连续流过模拟箱。

Here is an example input script that calculates the viscosity of liquid Ar via the GK formalism:
下面是一个示例输入脚本,通过 GK 形式计算液体 Ar 的粘度:

# Sample LAMMPS input script for viscosity of liquid Ar

units       real
variable    T equal 200.0       # run temperature
variable    Tinit equal 250.0   # equilibration temperature
variable    V equal vol
variable    dt equal 4.0
variable    p equal 400     # correlation length
variable    s equal 5       # sample interval
variable    d equal $p*$s   # dump interval

# convert from LAMMPS real units to SI

variable    kB equal 1.3806504e-23    # [J/K] Boltzmann
variable    atm2Pa equal 101325.0
variable    A2m equal 1.0e-10
variable    fs2s equal 1.0e-15
variable    convert equal ${atm2Pa}*${atm2Pa}*${fs2s}*${A2m}*${A2m}*${A2m}

# setup problem

dimension    3
boundary     p p p
lattice      fcc 5.376 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region       box block 0 4 0 4 0 4
create_box   1 box
create_atoms 1 box
mass         1 39.948
pair_style   lj/cut 13.0
pair_coeff   * * 0.2381 3.405
timestep     ${dt}
thermo       $d

# equilibration and thermalization

velocity     all create ${Tinit} 102486 mom yes rot yes dist gaussian
fix          NVT all nvt temp ${Tinit} ${Tinit} 10 drag 0.2
run          8000

# viscosity calculation, switch to NVE if desired

velocity     all create $T 102486 mom yes rot yes dist gaussian
fix          NVT all nvt temp $T $T 10 drag 0.2
#unfix       NVT
#fix         NVE all nve

reset_timestep 0
variable     pxy equal pxy
variable     pxz equal pxz
variable     pyz equal pyz
fix          SS all ave/correlate $s $p $d &
             v_pxy v_pxz v_pyz type auto file S0St.dat ave running
variable     scale equal ${convert}/(${kB}*$T)*$V*$s*${dt}
variable     v11 equal trap(f_SS[3])*${scale}
variable     v22 equal trap(f_SS[4])*${scale}
variable     v33 equal trap(f_SS[5])*${scale}
thermo_style custom step temp press v_pxy v_pxz v_pyz v_v11 v_v22 v_v33
run          100000
variable     v equal (v_v11+v_v22+v_v33)/3.0
variable     ndens equal count(all)/vol
print        "average viscosity: $v [Pa.s] @ $T K, ${ndens} atoms/A^3"

第五种方法与上述 Green-Kubo 方法相关,但使用爱因斯坦公式,类似于自扩散率的爱因斯坦均方位移公式。时间积分动量通量起到笛卡尔坐标的作用,其均方位移在足够长的时间内随时间线性增加。


第六种是周期摄动法,这也是一种非平衡MD方法。然而,它不是测量响应于施加的速度梯度的动量通量,而是测量响应于施加的应力的速度分布。通过 fix accelerate/cos ​​​​​​​ 命令向系统添加余弦形状的周期性加速度,并使用compute viscosity/cos ​​​​​​​命令监视生成的速度分布并在恒温之前消除速度偏差。

Ref:

 

  1. (Daivis and Todd) Daivis and Todd, Nonequilibrium Molecular Dynamics (book), Cambridge University Press, Nonequilibrium Molecular Dynamics, (2017).
  2. (Hess) Hess, B. The Journal of Chemical Physics 2002, 116 (1), 209-217.

370f1e7a41bb4f9f8259b9ae65b1b457.png

 

相关推荐

最近更新

  1. TCP协议是安全的吗?

    2024-04-23 20:14:03       18 阅读
  2. 阿里云服务器执行yum,一直下载docker-ce-stable失败

    2024-04-23 20:14:03       19 阅读
  3. 【Python教程】压缩PDF文件大小

    2024-04-23 20:14:03       18 阅读
  4. 通过文章id递归查询所有评论(xml)

    2024-04-23 20:14:03       20 阅读

热门阅读

  1. Redis雪崩

    2024-04-23 20:14:03       18 阅读
  2. python多线程详解

    2024-04-23 20:14:03       16 阅读
  3. Ubuntu搭建RP2040开发环境-1

    2024-04-23 20:14:03       18 阅读
  4. Springboot2.7解决静态资源302问题

    2024-04-23 20:14:03       14 阅读
  5. LeetCode 42. 接雨水 - PHP

    2024-04-23 20:14:03       12 阅读
  6. 2023年图灵奖揭晓

    2024-04-23 20:14:03       13 阅读