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.
# ************ 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}
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}
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
generate_force.log
file is a good optionIf you are lucky, maybe you don’t need to modify anything
2. Accumulated conductance