Hi,
I am a beginner at Lammps. Recently I calculated the thermal conductivity of polyethylene terephthalate using Fourier law. I used “fix langevin tally yes” to output energy accumulation but the output data was super huge (e.g. higher than 20000 kcal/mol within the simulated time 5 ns). I checked my .in file (shown as below) and am not sure if it was wrong in some settings. Hope someone could give some suggestions.
Thanks in advance.
Initialization
units real
atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style opls
improper_style harmonic
pair_style lj/cut/coul/long 11.0 11.0
pair_modify mix geometric
special_bonds lj/coul 0.0 0.0 0.5
kspace_style pppm 0.0001
dimension 3
newton on
boundary p p p
atom_style full
neighbor 2.0 bin
neigh_modify every 1 delay 0 check yes
Data reading
read_data system.data
include “system.in.settings”
#####################################################################
Setting->system definition
velocity all create 500.0 1232
timestep 1
energy minimization
min_style cg
minimize 1e-8 1e-8 1000000 1000000
annealing
thermo 50
thermo_style custom step temp press pxx pyy pzz vol pe ke etotal
fix 1 all nvt temp 500.0 1000.0 100
dump 1 all custom 5000 heating.xyz id type x y z
run 5000
unfix 1
undump 1
reset_timestep 0
thermo 1000
thermo_style custom step temp press pxx pyy pzz vol pe ke etotal
fix 2 all nvt temp 1000.0 1000.0 100
dump 2 all custom 10000 1000k.data id mol type x y z vx vy vz
run 3000000
run 5000000
unfix 2
undump 2
reset_timestep 0
thermo 140
thermo_style custom step temp press pxx pyy pzz vol pe ke etotal
fix 3 all nvt temp 1000.0 300.0 100
dump 3 all custom 14000 cooling.xyz id mol type x y z
run 14000
unfix 3
undump 3
reset_timestep 0
equilibrium simulation
thermo 1000
thermo_style custom step temp press pxx pyy pzz vol pe ke etotal
fix 4 all npt temp 300.0 300.0 100 iso 1.0 1.0 1000
dump 4 all custom 100000 equil.data id mol type x y z vx vy vz
run 5000000
run 10000000
unfix 4
undump 4
reset_timestep 0
production simulation
thermo 1000
thermo_style custom step temp press pxx pyy pzz vol pe ke etotal
fix 5 all nvt temp 300.0 300.0 100.0
dump 5 all custom 50000 prod.data id mol type x y z vx vy vz
restart 5000000 restartequil.restart
run 5000000
run 5000000
unfix 5
undump 5
reset_timestep 0
#########################################################################
thermal conductivity calculation
compute ke1 all ke/atom
variable kb equal 1.3806504e-23 # [J/K]
variable kb equal 0.001987 # [kcal/mol/K]
variable temp1 atom c_ke1/1.5/${kb}
bins setup
variable X1 equal xlo
variable X2 equal xhi/2
variable Nlay equal 100
variable Dscale equal 1/{Nlay}
variable Len equal {X2}-{X1}
variable Dz equal {Len}/${Nlay}
heat source and sink setup coordination
variable BL1 equal {X1}
variable BL2 equal {X1}+{Dz}
variable BR1 equal {X2}
variable BR2 equal {X2}+{Dz}
region setup
region Hot block {BL1} {BL2} INF INF INF INF units box
region Cold block {BR1} {BR2} INF INF INF INF units box
group Hot region Hot
group Cold region Cold
fix 6 Hot langevin 320 320 100 14565 tally yes
fix 7 Cold langevin 280 280 100 16576 tally yes
fix 6 all momentum 1 linear 1 1 1
fix 7 all momentum 1 linear 1 1 1
fix 8 all nve
variable Time equal step
variable EL equal f_6
variable ER equal f_7
compute 18 all chunk/atom bin/1d x lower ${Dscale} units reduced
fix 9 all ave/chunk 1 200000 200000 18 v_temp1 file temp_relax.dat
dump 6 all xyz 100000 nve.xyz
thermo 20000
thermo_style custom step temp ke pe etotal press vol
run 200000
unfix 9
undump 6
reset_timestep 0
fix 10 all ave/chunk 1 10000000 10000000 18 v_temp1 file temp_equ.dat
fix 10 all ave/chunk 1 5000000 5000000 18 v_temp1 file temp_equ.dat
fix E_out all print 2000 “{Time} {EL} ${ER}” file Ener_equ.dat title “Time E1 E2” screen no
dump 7 all xyz 1000000 equ.xyz
thermo 100000
thermo_style custom step temp ke pe etotal press vol
run 5000000