I am trying to get equation of state of silicon in sw potential for T =0.05 and density =0.445 but I am unable to get right pressure at these density could anyone tell me what’s wrong in my script.

```
#variable
variable sigma equal 2.0951
variable epsilon equal 2.1683
variable tempScale equal 2.5173e4
variable densityScale equal 5.0571
variable T equal 0.05
variable t equal $(v_T*v_tempScale)
variable dt equal 0.005
variable rho equal 0.489
variable a equal (8.0*$(v_sigma)*$(v_sigma)*$(v_sigma)/${rho})^(1.0/3.0)
#main script
units metal
atom_style atomic
boundary p p p
lattice custom $a &
a1 1.0 0.0 0.0 &
a2 0.0 1.0 0.0 &
a3 0.0 0.0 1.0 &
basis 0.0 0.0 0.0 &
basis 0.0 0.5 0.5 &
basis 0.5 0.0 0.5 &
basis 0.5 0.5 0.0 &
basis 0.25 0.25 0.25 &
basis 0.25 0.75 0.75 &
basis 0.75 0.25 0.75 &
basis 0.75 0.75 0.25
region myreg block 0 4 &
0 4 &
0 4
create_box 8 myreg
create_atoms 1 region myreg &
basis 1 1 &
basis 2 2 &
basis 3 3 &
basis 4 4 &
basis 5 5 &
basis 6 6 &
basis 7 7 &
basis 8 8
mass * 28.0855
velocity all create $t 42342 loop geom
pair_style sw
pair_coeff * * Si.sw Si Si Si Si Si Si Si Si
neighbor 2.0 bin
neigh_modify every 1 delay 0 check no
timestep ${dt}
thermo 10
thermo_style custom step density ke pe etotal vol temp press
fix 1 all nvt temp $t $t 0.1
#output variable
compute myTemp all temp
compute Press all pressure myTemp
variable curT equal temp
variable curStep equal step
variable curV equal vol
variable curH equal enthalpy
variable curKE equal ke
variable curPE equal pe
variable curE equal etotal
variable curP equal (c_Press[1]+c_Press[2]+c_Press[3])/3.0
variable curD equal density
fix 2 all print 1 "${curStep} ${curD} ${curKE} ${curPE} ${curE} ${curV} ${curH} ${curT} ${curP}" screen no file Si_rho${rho}_T0.05_metal.dat title "#Step Density KE PE TotalE Vol Enthalpy Temp Press"
dump 1 all custom 10 Si_rho${rho}_T0.05_.dump id type xu yu zu x y z vx vy vz
log Si_rho${rho}_T0.05_.log
run 100000
```

Thanks