Due to my previous stupidity, I limited the scope of the project. The previous project needed to ensure that the number of atoms on the left and right ends of the illusion interface was equal. But in principle it is not required.
Notes:
When use this method to calculation some material with very low thermal conductivity, maybe one can get some very strange results (the spectral heat flux will turn to negative or positive). I don't know why, but one can increase the temperature difference between the heat source and sink to address this problem. I would appreciate it if someone could email me the mechanism behind it.
If you use the code, the following citations are highly recommended.
- K. Sääskilahti, J. Oksanen, J. Tulkki, and S. Volz, Phys. Rev. B 90, 134312 (2014)
- K. Sääskilahti, J. Oksanen, S. Volz, and J. Tulkki, Phys. Rev. B 91, 115426 (2015)
- Xu K, Deng S, Liang T, et al, Efficient mechanical modulation of the phonon thermal conductivity of Mo6S6 nanowires. Nanoscale, 2022
- Yao Y, Ren G, Yu Y, et al., Thermal conduction mechanism of ferroelastic Zr‐Y‐Yb‐Ta‐Nb‐O high‐entropy oxides with glass‐like thermal conductivity. Journal of the American Ceramic Society, 2022
- relax_thermal.in (input file for lammps, and use to get the atoms velocity)
- forces.in (input file for force_calculate.py, in order to call the python interface of lammps)
- SHC_generate.py (Main program, useful for calculating spectral decomposition thermal conductance)
1. According to one specific needs, modify the potential format in the files (relax_thermal.in and forces.in)
- Modify in relax_thermal.in and forces.in file
# ************ Modify in relax_thermal.in and forces.in file ************
# ********************* Potential function setting *****************
pair_style tersoff
pair_coeff * * BNC.tersoff B C N
# ************ Modify in relax_thermal.in and forces.in file ************
# ************************* For Y direction **************************************
variable y_max equal ly # Depends on the direction of heat transport
variable P equal ${y_max}/2-100
variable P1 equal ${y_max}/2+100
variable tmp1 equal ${P1}-${P}
variable tmp equal ${tmp1}/40
variable L1 equal ${P}+3*${tmp}
variable R1 equal ${P1}-3*${tmp}
- Modify in relax_thermal.in file
variable dmid equal 6 ## Set to 6 (A) here, one can modify it
variable middle equal ${y_max}/2
variable mid_left equal ${middle}-${dmid}
variable mid_right equal ${middle}+${dmid}
- Modify in relax_thermal.in file (For
dn
parameter)
# ************************** Collect Velocities for the calculation of force constants *************************
variable dn equal 15
dump vels interface custom ${dn} vels.dat id type vx vy vz
dump_modify vels format line "%d %d %0.8g %0.8g %0.8g"
dump_modify vels sort id
postprocessor = SHCPostProc(compact_velocities_file,
fileprefix, # forces.in
dt_md = md_timestep, # timestep
dn = dn, # Dump interval
scaleFactor = unit_scaling_factor,
LAMMPSDumpFile = atomic_velocities_file, # velocity file
widthWin = frequency_window_width,
LAMMPSInFile = in_file,
in_plane = False,
out_of_plane = False,
reCalcVels = True,
reCalcFC = True)
g++ compactify_vels.cpp -o compactify_vels
-
lmp_mpi < relax_thermal.in
-
python SHC_generate.py
If one want to run lammps and python in parallel, one can use following command (need to set
if_MPI = True
in SHC_generate.py file):mpirun -np 10 python SHC_generate.py # for generate the force constant in parallel
Maybe one need to modify the input file according to the error report
- tips: check the
generate_force.log
file is a good option
If you are lucky, maybe you don’t need to modify anything
- tips: check the