Dear All,
I am simulating the thermal conductance along the z axis for bilayer graphene sandwiched by Al2O3 system. This is my input setting.
variable x equal 47.59
variable y equal 41.214148966101
variable z equal 136.636
variable t equal 300.0
#Setup parameters
units real
atom_style atomic
read_data atomic_structure_atomic
mass 1 16.0
mass 2 27.0
mass 3 12.0
mass 4 12.0
velocity all create $t 87287
pair_style hybrid/overlay lj/cut 12.0 morse 6.287 tersoff shift -0.0462
pair_coeff 1 1 lj/cut 0.194633996 3.541
pair_coeff 1 2 lj/cut 0.31351524159 4.02
pair_coeff 1 3 lj/cut 0.0 1
pair_coeff 1 4 lj/cut 0.0 1
pair_coeff 2 2 lj/cut 0.50500834301 4.499
pair_coeff 2 3 morse 10.81786819 1.738 2.246
pair_coeff 2 4 morse 10.81786819 1.738 2.246
pair_coeff * * tersoff ./Copt.tersoff NULL NULL C C
neighbor 0.3 bin
neigh_modify delay 0 every 1
#layers for heat flux
region hot block INF INF INF INF 28.69356 41.810616
region cold block INF INF INF INF 96.32838 97.01156
compute Thot all temp/region hot
compute Tcold all temp/region cold
#1st equilibirum run
fix 1 all nvt temp $t $t 100.0
thermo 1
run 3000000
velocity all scale $t
unfix 1
#2nd equilibrium run
fix 1 all nve
fix hot all heat 1 23.0609 region hot
fix cold all heat 1 -23.0609 region cold
thermo_style custom step temp c_Thot c_Tcold
thermo_modify colname c_Thot Temp_hot colname c_Tcold Temp_cold
thermo 10
run 10000
#thermal conductivity calculation
compute ke all ke/atom
variable temp atom c_ke/23.0609/1.5/8.617333262145e-5
compute layers all chunk/atom bin/1d z lower 1.0 units reduced
fix 2 all ave/chunk 10 100 1000 layers v_temp file profile.heat
variable tdiff equal f_2[20][3]-f_2[36][3]
fix ave all ave/time 1 1 1000 v_tdiff ave running start 13000
variable kappa equal 23.0609/(lx*ly)*lz/f_ave
thermo_style custom step temp c_Thot c_Tcold v_tdiff f_ave
thermo_modify colname c_Thot Temp_hot colname c_Tcold Temp_cold &
colname v_tdiff dTemp_step colname f_ave dTemp
run 3000000
print "Running average thermal conductivity: $(v_kappa:%.4f)"
This is the error message I received after I launched the calculation.
LAMMPS (2 Aug 2023 - Update 1)
variable x equal 47.59
variable y equal 41.214148966101
variable z equal 136.636
variable t equal 300.0
#Setup parameters
units real
atom_style atomic
read_data atomic_structure_atomic
Reading data file ...
triclinic box = (0 0 0) to (47.59 41.214149 136.636) with tilt (-23.795 0 0)
WARNING: Triclinic box skew is large. LAMMPS will run inefficiently. (../domain.cpp:220)
2 by 2 by 8 MPI processor grid
reading atoms ...
19600 atoms
read_data CPU = 0.131 seconds
mass 1 16.0
mass 2 27.0
mass 3 12.0
mass 4 12.0
velocity all create $t 87287
velocity all create 300 87287
#pair_style hybrid/overlay lj/cut 10.0 morse 6.287 tersoff shift -0.0462 kolmogorov/crespi/z 20.0
pair_style hybrid/overlay lj/cut 12.0 morse 6.287 tersoff shift -0.0462
pair_coeff 1 1 lj/cut 0.194633996 3.541
pair_coeff 1 2 lj/cut 0.31351524159 4.02
pair_coeff 1 3 lj/cut 0.0 1
pair_coeff 1 4 lj/cut 0.0 1
pair_coeff 2 2 lj/cut 0.50500834301 4.499
pair_coeff 2 3 morse 10.81786819 1.738 2.246
pair_coeff 2 4 morse 10.81786819 1.738 2.246
pair_coeff * * tersoff ./Copt.tersoff NULL NULL C C
neighbor 0.3 bin
neigh_modify delay 0 every 1
#layers for heat flux
region hot block INF INF INF INF 28.69356 41.810616
region cold block INF INF INF INF 96.32838 97.01156
compute Thot all temp/region hot
compute Tcold all temp/region cold
#1st equilibirum run
fix 1 all nvt temp $t $t 100.0
fix 1 all nvt temp 300 $t 100.0
fix 1 all nvt temp 300 300 100.0
thermo 1
run 3000000
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 12.3
ghost atom cutoff = 12.3
binsize = 6.15, bins = 12 7 23
6 neighbor lists, perpetual/occasional/extra = 6 0 0
(1) pair lj/cut, perpetual, skip from (4)
attributes: half, newton on
pair build: skip
stencil: none
bin: none
(2) pair morse, perpetual, skip from (5)
attributes: half, newton on, cut 6.587
pair build: skip
stencil: none
bin: none
(3) pair tersoff, perpetual, skip from (6)
attributes: full, newton on, cut 2.4
pair build: skip
stencil: none
bin: none
(4) neighbor class addition, perpetual
attributes: half, newton on
pair build: half/bin/newton/tri
stencil: half/bin/3d/tri
bin: standard
(5) neighbor class addition, perpetual, trim from (4)
attributes: half, newton on, cut 6.587
pair build: trim
stencil: none
bin: none
(6) neighbor class addition, perpetual
attributes: full, newton on, cut 2.4
pair build: full/bin/atomonly
stencil: full/bin/3d
bin: standard
Per MPI rank memory allocation (min/avg/max) = 6.496 | 8.707 | 10.81 Mbytes
Step Temp E_pair E_mol TotEng Press
0 300 4.3067738e+08 0 4.306949e+08 4.4405907e+08
1 9.6397316e-73 2.5266396e+09 0 2.5266396e+09 2.5914761e+09
ERROR: Lost atoms: original 19600 current 19595 (../thermo.cpp:488)
Last command: run 3000000
I checked the unit for the potential and energy settings. It is said that the metal units would be converted into real unit by LAMMPS for the tersoff potential automatically. I did not changed the tersoff potential and only change the unit for energy from eV to Kcal/mol in the lj and mose potential.
I think that the force was not that large to blow up some atoms out of the unit cell. I do not understand why I still receive this error message.
Would anyone please give me some suggestions/comments on the parameter settings in my input file?
Thank you in advance.