这是我的in文件，算出来的值太小了，10的负11次方，我也不知道哪里错了，就是按照example里面的例子来的，希望各位了老师帮我看一下：
#        Green-Kubo 方法计算热导率

echo                                                screen


units                                                real
#能量=Kcal/mole
#时间=fs
#长度=A


#定义参数
variable                                T  equal        86
variable                                A  equal        5.72
variable                                dt equal        0.01


dimension                                3
boundary                                p p p
newton                                        on


atom_style      atomic
neighbor                                 2.0 bin
neigh_modify                 every 2 delay 10 check yes


#初始化坐标、势函数、速度等数据（ar box）
lattice                                        fcc        ${A}
region                                        box block        0 10 0 10 0 10                #10*10*10*4=4000 原子
create_box                        1        box
create_atoms                1        box
mass                                                1        39.95


pair_style                         lj/cut 10                #LJ势的截断半径是10A
pair_coeff                        * * 3.9913e-21        3.405
velocity all create ${T} 12345          #根据温度初始化速度


#在NVT系综下平衡体系
fix                                                        NVT all nvt temp ${T} ${T} 0.5  #Nose-Hoover热浴
timestep                                 ${dt}                                                    #设置积分步长
thermo_style                custom        step        temp        press                       
thermo                                        1000                                        #1000步输出一次

#以上跑若干步，目的是使系统达到热力学平衡态
run                                                        20000


#在输出阶段用NVE系综                                                                       
unfix                                                NVT


reset_timestep  0                                                                                        #重新从0开始步长


fix                                                        NVE all nve


#在输出阶段计算热流
compute                                        KE all        ke/atom                                #计算单原子动能
compute                                        PE all        pe/atom                                #计算单原子势能
compute                                        V all          stress/atom NULL virial                #计算单原子位能
compute                                        J all          heat/flux KE PE V                        #计算单原子热流


#在输出阶段计算热流自相关函数        <J(0)J(t)>
variable                                Ns equal        10                                                #取样间隔
variable                                Nc equal        200                                                #关联数据量（最大关联时间为 Ns*Nc*Dt）
variable                                Np equal        ${Ns}*${Nc}*10                        #输出步数


variable        Jx equal c_J[1]/vol                                                            #三个方向的热流分量
variable        Jy equal c_J[2]/vol
variable        Jz equal c_J[3]/vol


fix                JJ all ave/correlate        ${Ns} ${Nc} ${Np} c_J[1] c_J[2] c_J[3] type auto file profile-ar.txt        #将每个原子热流密度的各个分量取样输出


variable                                scale equal        ${Ns}*dt/$T/$T/vol
variable                                 k11 equal trap(f_JJ[3])*${scale}
variable                                 k22 equal trap(f_JJ[4])*${scale}
variable                                 k33 equal trap(f_JJ[5])*${scale}




thermo                ${Np}
thermo_style    custom step temp v_Jx v_Jy v_Jz v_k11 v_k22 v_k33

#输出阶段
dump                                        2 all custom 100 dump2.lammpstrj id type  x y z vx vy vz fx fy fz
run                                                100000


variable                         kappa equal (v_k11+v_k22+v_k33)/3.0
print "thermal conductivity: ${kappa}"


算出来的热导率太小了，到底问题在哪里啊？是我的那几个参数的单位换算出错了吗？


PCA%_414P3N9%$_NVNRNU5T.png 返回小木虫查看更多