Thermal Conductivity - Green-kubo


I submitted a request called “Thermal Conductivity - Green-kubo Argon Example” and tried to cancel my request with the link provided because I figured out what was the issue. Now I would like to ask why is this method not to work properly when trying to find the thermal conductivity of Aluminum. I utilized a eam potential from NIST and changed the units and conversion variables appropriately and I’m getting an output of 13.45 W/mK @ 20 C when it should be much higher. Are there parameters that need to be altered? Below is the script that I’m using.


Sample LAMMPS input script for thermal conductivity of solid Ar

#units real
units metal
variable T equal 293.15
variable V equal vol
#variable dt equal 4.0
variable dt equal 0.004
variable p equal 200 # correlation length
variable s equal 10 # 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 kCal2J equal 4186.0/6.02214e23
#variable A2m equal 1.0e-10
#variable fs2s equal 1.0e-15
#variable convert equal {kCal2J}*{kCal2J}/{fs2s}/{A2m}
variable eV2J equal 1.602176565e-19 #[J]
variable A2m equal 1.0e-10
variable ps2s equal 1.0e-12
variable convert equal {eV2J}*{eV2J}/{ps2s}/{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
lattice fcc 4.0 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
pair_style eam/alloy
pair_coeff * * Al99.eam.alloy Al
timestep ${dt}
thermo $d

equilibration and thermalization

velocity all create $T 102486 mom yes rot yes dist gaussian
fix NVT all nvt temp $T $T 0.01 drag 0.2
run 80000

thermal conductivity calculation, switch to NVE if desired

#unfix NVT
#fix NVE all nve

reset_timestep 0
dump 1 all atom 1000 dump.output
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
fix JJ all ave/correlate $s $p d & c_flux[1] c_flux[2] c_flux[3] type auto file J0Jt.dat ave running variable scale equal {convert}/${kB}/$T/$T/$V*s*{dt}
variable k11 equal trap(f_JJ[3]){scale} variable k22 equal trap(f_JJ[4])*{scale}
variable k33 equal trap(f_JJ[5])
thermo_style custom step temp density v_Jx v_Jy v_Jz v_k11 v_k22 v_k33
run 100000
variable k equal (v_k11+v_k22+v_k33)/3.0
variable ndens equal count(all)/vol
print “average conductivity: $k[W/mK] @ T K, {ndens} /A^3”

GK only calculates the phonon contribution to thermal conductivity, in which case 13.5 W/m-K, while not correct, is within a factor of ~3.

Hello Mr. Deangelis,

Thank you for the information and help, I really appreciate it. I had another question regarding utilizing the green-kubo method. What parameters should I alter or change when trying to use this method in calculating the thermal conductivity of a polymer? Are there certain parameters that I should alter when dealing with different types of molecules from the sample script that is given on lammps. Obviously I would have different pair_styles and coefficients, I’m specifically asking about the simulation parameters. Also, any suggestions on how to use this method for accurately calculating thermal conductivity would help.


If you get answers on the mail list, great.
But you’d probably be better off finding papers
that use GK on a system similar to your, and
trying to reproduce their results using the
details of the GK params they provide.