Dear all LAMMPS users,
I am new user of lammps and I am trying to find the effect of tensile strain on graphene thermal conductivity. The program runs well without applying the tensile strain.
When I applied tensile strain, I am getting error as “Divide by 0 in variable formula” I find that in the denominator the value is getting zero. In the The input script is given below. In finding the Variable Z the denominator f_5 is getting zero after 50000 run i.e after 2nd equilibration run. The v_tdiff(without the application of strain) is attached here.
I want to know how the f_5 value is getting zero in case of tensile strain which is related to v_tdiff and it calculated as f_2[96][3]-f_2[1][3]. As I am new user of lammps
f_2[96][3]-f_2[1][3] how this getting calculated which is getting stored in profile.mp which is attached here(without application of strain).
Please refer the input script below.
variable T equal 300
variable V equal vol
variable dt equal 0.001
#variable p equal 1000 # correlation length
#variable s equal 5 # sample interval
#variable d equal $p*$s # dump interval
convert from LAMMPS real units to SI
variable kB equal 1.3806504e-23 # [J/K] Boltzmann
variable eV2J equal 1.602e-19
variable A2m equal 1.0e-10
variable ps2s equal 1.0e-12
variable convert equal {eV2J}*{eV2J}/{ps2s}/{A2m}
variable alat equal 2.46
setup problem
units metal
boundary p p p
atom_style atomic
read_data Read_data_1st
mass 1 12.0
pair_style tersoff
pair_coeff * * /home/kpark39/programs/lammps-1Feb14/potentials/SiC_1994.tersoff C
neighbor 2.0 bin
neigh_modify delay 5
timestep ${dt}
variable tmp equal lx
variable lo equal ${tmp}
variable strain equal (lx-v_lo)/v_lo
variable p1 equal “-pxx/10000”
variable p2 equal “-pyy/10000”
2nd equilibration run
compute ke all ke/atom
units of Kb = J/K, units of ke = eV or J, so temp becomes K
variable temp atom c_ke*{eV2J}/(1.5*{kB})
fix 1 all nve
compute cc1 all chunk/atom bin/1d x lower 5
fix 2 all ave/chunk 10 100 1000 cc1 v_temp file profile.mp #units box
fix 3 all thermal/conductivity 10 x 190 swap 2
fix 4 all temp/berendsen 300 300 100
variable tdiff equal f_2[96][3]-f_2[1][3]
thermo_style custom step temp epair etotal f_3 v_tdiff
thermo 1000
dump 2 all cfg 1000 dump.config.*.cfg mass type xs ys zs vx vy vz x y z v_temp
run 50000
unfix 4
thermal conductivity calculation
fix 5 all ave/time 1 1 1000 v_tdiff ave running file ave.dat
variable Z equal f_3*{eV2J}*1000/2/(110*3.4)/{A2m}/{dt}/100000/{ps2s}/f_5
fix 7 all deform 1 x erate 0.0001 units box
fix def all print 1000 “{strain} {p2} ${p1}” file Graphene_tensile_data_300K.txt
thermo_style custom step temp epair etotal f_3 v_tdiff f_5 v_Z
thermo 1000
dump 7 all custom 5000 Graphene_thermal_deform*.lammpstrj id type x y z
run 100000
ave.dat (1.57 KB)
profile.mp (833 KB)