I add “compute tc command” into lammps, the usage for it is like this:
compute tc command
compute ID group-ID tc c_temp-ID c_flux-ID v_factor_ac v_factor_tc direction portion M tc_ntimestep_start tc_ntimestep_freq
- ID, group-ID are documented in compute command
- tc = style name of this compute command
- temp-ID = ID of a compute that calculates temperature
- flux-ID = ID of a compute that calculates heat flux vector
- factor_ac = ID of a variable whose value is used to transform the autocorrelation unit defined in “units command” to the units desired
- factor_tc = ID of a variable whose value is used to transform the thermal conductivity unit defined in “units command” to the units desired
- direction = x or y or z or iso, direction for thermal conductivity computation
- portion = first(total heat flux: jc + jv) or second(just conductive portion: jc only) for heat flux contribution
- M = M value in discrete Green Kubo thermal conductivity expression
- tc_ntimestep_start = starting timestep for “ac.dat” and “tc.dat” output(timestep units, not time units)
- tc_ntimestep_freq = output frequency for “ac.dat” and “tc.dat” (timestep units, not time units)
The command can output two files “ac.dat”(autocorrelation) and “tc.dat”(thermal conductivity), every “tc_ntimestep_freq” timesteps, and start at “tc_ntimestep” timestep.
And the value returned from this command is thermal conductivity value at each timestep, so you can obtain a file like “tc_time.dat” using fix ave/time command, then you can get a curve, whose vertical value is thermal conductivity, and horizontal value is timestep.
“v_factor_ac” and “v_factor_tc” can be used to transform dimensionless units defined by “units command” in lammps to SI units or other units you desire.
Citing: discrete GK expression is cited from PK Schelling’s article “Comparison of atomic-level simulation methods for computing thermal conductivity” (2002), as attached.
Note: Because thermal conductivity computation is accumulated as timestep increases, so you can cite the compute tc command only once.
compute myKE all ke/atom
compute myPE all pe/atom
compute myStress all stress/atom virial
variable factor_ac equal 1.0
variable factor_tc equal 1.3806504e-23*sqrt(1.67e-21/6.633e-26)/3.405e-10^2
compute jflux all heat/flux myKE myPE myStress
compute tc all tc c_thermo_temp c_jflux v_factor_ac v_factor_tc iso first 10000 900000 100000
fix tc_out all ave/time 1 1 1 c_tc file tc_time.dat
compute_tc.cpp (14.9 KB)
compute_tc.h (1.23 KB)
Schellingetal-2002.pdf (147 KB)