Hi,
I’m using lammps to induce a temperature gradient in a Silicon Sample. My output data shows that my the gradient ends up ranging from ~100K to ~4000K instead of my targeted 1000-2000. Am I initializing it wrong, or using fix npt incorrectly.
Here is my input script:
Silicon
#Initialization
units metal
dimension 3
boundary p p p
atom_style atomic
lattice diamond 5.430710
region box block 0 20 0 5 0 5
create_box 1 box
create_atoms 1 box
mass 1 28.085
pair_style sw
pair_coeff * * Si.sw Si
Define Regions
region 1 block 0 2 INF INF INF INF
region 2 block 10 12 INF INF INF INF
group hot region 1
group cold region 2
region 3 block 2 10 INF INF INF INF
group m region 3
region 4 block 12 20 INF INF INF INF
group n region 4
group other union m n
#set up eqilibrium
velocity all create 1683 21564
fix eq all npt temp 1683 1683 .1 iso 0 0 1
#Compute
variable kb equal 8.617339e-5 #eV/K
variable mass_g equal mass(all)/6.02e23
variable vol_cm3 equal vol1.0e-24
variable mydensity equal v_mass_g/v_vol_cm3
compute k all ke/atom
compute 1 all temp
variable t equal c_1
variable temp atom c_k/(1.5v_kb)
#Outputs
dump id all custom 100 dump.Si.06 id x y z c_k
log log.Si.06
thermo_style custom step temp press v_mydensity enthalpy
thermo 50
run 6000
thermo_modify lost ignore flush yes
#set temp bc’s
unfix eq
velocity hot create 2000 13544
velocity cold create 1000 18450
fix NPT hot npt temp 2000 2000 .1 iso 0 0 1
fix NPT2 cold npt temp 1000 1000 .1 iso 0 0 1
fix NVE other nve
Correlate
fix a all ave/spatial 10 20 200 x lower 1 c_k &
v_temp file cor.dat6
thermo_style custom step temp press v_mydensity enthalpy f_NPT f_NPT2
run 50000
Thanks so much!
Kelly