units METAL variable T equal 70 variable V equal vol variable dt equal 0.004 variable p equal 200 variable s equal 10 variable d equal $p*$s variable kB equal 1.3806504e-23 variable ev2J equal 1.602176565e-19/6.02214e23 #ev to J coeficent convertor 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 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 $T 102486 mom yes rot yes dist gaussian fix NVT all nvt temp $T $T 10 drag 0.2 run 8000 # thermal conductivity calculation, switch to NVE if desired #unfix NVT #fix NVE all nve reset_timestep 0 compute myKE all ke/atom compute myPE all pe/atom compute myStress all stress/atom 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])*${scale} thermo_style custom step temp 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"