|
LAMMPS讲解49-纳米流体Green-Kubo法计算热导率和粘度Green-Kubo方法的计算热导率的原理是利用系统热流的时间自相关函数进行计算，其公式为
由于体系是各向同性的所以xyz三个方向上的热流自相关函数计算出来的热导率应该是相同的，所以可以计算出来三个方向上的热导率系数然后取平均。在固液共存相中这里没有考虑焓的影响，所以计算结果有问题。因此内容可选择性借鉴。在LAMMPS中利用compute heat/flux来计算系统的热流。这个compute会输出6个全局矢量。该矢量的前三个分量分别是全局热流，除以系统体积就为热流了。以x方向为例，Jx equal c_flux[1]/vol。把这个代入上面公式中的
粘度的计算与热导类似，使用的是应力的自相关。其公式为
建模分为两个部分第一部分可以借鉴上一节的过程： ######基本模拟规则设置######### units metal dimension 3 boundary p p p neighbor 3.0 bin neigh_modify delay 0 every 1 chech yes atom_style atomic pair_style hybrid sw eam lj/cut 10 #######定义几何模型########## #使用LAMMPS内部命令建立模型需要三个命令联用：region命令，create_box命令和create_atoms命令 region box block -90 90 -90 90 -90 90 units box create_box 2 box region sphere_cu 0.0 0.0 0.0 45 units box side in lattice fcc 3.615 create_atoms 1 sphere_cu region wat subtract 2 box sphere_cu lattic fcc 5.10873 create_atoms 2 wat ########定义物理模型############## pair_style * * sw NULL WT pair_style 1 1 eam Cu_u6.eam pair_style 1 2 lj/cut 0.006 3.0 #水和铜的参数，可以通过一些实验参数进行调节 write_data pre.data run 1 通过控制水区域和铜区域的比例就可以调控纳米颗粒的体积分数。写出data文件后重新读入data文件开始模拟。in文件如下： ######基本模拟规则设置######### units metal dimension 3 boundary p p p neighbor 3.0 bin neigh_modify delay 0 every 1 chech yes atom_style atomic pair_style hybrid sw eam lj/cut 10 #######定义几何模型########## read_data pre.data replicate 5 5 5 #利用复制命令建立更大的模型，具体模型多大需要经验。具体可以跑一下算例 ########定义物理模型############## pair_style * * sw NULL WT pair_style 1 1 eam Cu_u6.eam pair_style 1 2 lj/cut 0.006 3.0 #水和铜的参数，可以通过一些实验参数进行调节 #######定义分组信息######### group cu type 1 group wat type 2 #######弛豫部分####### timestep 0.002 thermo 1000 thermo_style custom step temp etotal press thermo_modify flush yes velocity all create 300 3265456 rot yes dist gaussian loop local dump 1 all custom 10000 dump.lammpstrj id type x y z fix fxnvt move nvt temp 300 300 0.2 run 1000000 write_restart restart.equ ########粘度和热导率计算之Green-Kubo法######### variable p equal 20000 # correlation length variable s equal 5 # sample interval variable d equal $p*$s # dump interval variable t equal 300 variable kB equal 1.3806504e-23 variable A2m equal 1.0e-10 variable ps2s equal 1.0e-12 variable bar2pa equal 100000 variable eV2J 1.60218E-19 variable convert_eta equal ${bar2pa}*${bar2pa}*${A2m}*${A2m}*${A2m}*${ps2s} variable convert_kappa equal ${eV2J}*${eV2J}/${ps2s}/${A2m} reset_timestep 0 compute myKE all ke/atom compute myPE all pe/atom compute myStress all stress/atom NULL virial compute flux all heat/flux myKE myPE myStress variable Jx equal c_flux[1]/vol variable Jy equal c_flux[2]/vol variable Jz equal c_flux[3]/vol unfix fxnvt fix 1 all nve fix JJ all ave/correlate $s $p $d c_flux[1] c_flux[2] c_flux[3] type auto file profile.heatflux ave running variable scale equal $s*dt/$t/$t/${kB}/vol variable k11 equal trap(f_JJ[3])*${scale}*${convert_kappa} variable k22 equal trap(f_JJ[4])*${scale}*${convert_kappa} variable k33 equal trap(f_JJ[5])*${scale}*${convert_kappa}
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 profile.gk ave running variable scale equal $s*dt/$t/${kB}*vol variable diagfac equal 2-2/3 variable vxy equal trap(f_SS[3])*${scale}*${convert_eta} variable vxz equal trap(f_SS[5])*${scale}*${convert_eta} variable vyz equal trap(f_SS[7])*${scale}*${convert_eta}
thermo $d thermo_style custom step temp press pe ke v_k11 v_k22 v_k33 v_vxy v_vxz v_vyz thermo_modify flush yes run 1000000 计算结果分析
感谢鲍路瑶老师的分享，内容来自于鲍老师分享出来的资料 如有需要添加微信：lmp_zhushou 进入微信群，帮助他人，共建社区 获取完整版lammps讲义可以加微信lmp_zhushou或加入QQ群994359511 |
