units           metal
boundary        p p p

atom_style	atomic
atom_modify sort 0  0
# Set up FCC lattice with z axis along <110>

lattice                fcc 4.04 &
                       orient x 1 0 0 &
                       orient y 0 1 0 &
                       orient z 0 0 1

region         mycell block 0.0 20.0 0.0 20.0 0.0 20.0 units lattice
create_box     1 mycell
mass           1 26.98
create_atoms   1 box

# Using units of Rmin, so sigma = 2^-1/6 = 0.8908987

pair_style eam/alloy
pair_coeff * * Al.eam.alloy Al 

#pair_style adp
#pair_coeff * * AL-cu.adp Al Cu

# Relax box dimensions

fix             3 all box/relax aniso 0.0 vmax 1.0e-4 nreset 100

thermo          100
thermo_style	custom step temp pe etotal pxx pyy pzz lx ly lz 

min_modify      line quadratic
minimize	0.0 1.0e-6 10000 100000

# Define initial velocity

velocity        all create 300.0 87287 mom yes rot yes dist gaussian
write_restart restart.equil

# Start Run #1


log log.drag

clear
read_restart restart.equil

neighbor        0.2 bin
neigh_modify    every 1 delay 0 check yes
timestep        0.001
reset_timestep 0

# Pzz = 40.0, drag/damping term on

fix myhug all nphug temp 300.0 300.0 200.0 z 200000.0 200000.0 1000.0 drag 200.0 tchain 1 pchain 0 

# Specify reference state from paper, times 1000 atoms

#fix_modify myhug e0 -6334.0 p0 0.0 v0 680.73519
#fix_modify myhug e0 -129600.0 p0 0.0 v0 288542.70
#fix_modify myhug e0 -2073600000.0 p0 0.0 v0 288542.70

# Add fix energy to ouput etotal
 
#fix_modify myhug energy yes 

# Define output

variable dele  equal f_myhug[1] # energy delta [temperature]
variable us    equal f_myhug[2] # shock velocity [distance/time]
variable up    equal f_myhug[3] # particle velocity [distance/time]
variable pzz   equal pzz  # axial stress
variable tau equal 0.5*(pzz-0.5*(pxx+pyy)) # shear stress
variable time equal dt*step

thermo          1000
thermo_style custom step temp ke epair etotal pzz  lz v_us v_up
fix pr2 all print 10 "${time} ${us} ${up} " screen no append us-up-t-withdrag.dat title '#time us up (with drag)'

fix stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (with drag)'

run		1000000


log log.nodrag

clear
read_restart restart.equil

neighbor        2.0 bin
neigh_modify    every 1 delay 0 check yes

timestep        0.001
reset_timestep 0

# Pzz = 40.0, drag/damping term off

fix myhug all nphug temp 300.0 300.0 100.0 z 2000000.0 2000000.0 1000.0 drag 0.0 tchain 1 pchain 0 

# Specify reference state from paper, times 1000 atoms

#fix_modify myhug e0 -6334.0 p0 0.0 v0 680.73519
fix_modify myhug e0 -129600.0 p0 0.0 v0 288542.70

# Add fix energy to ouput etotal
 
fix_modify myhug energy yes 

# Define output

variable dele  equal f_myhug[1] # energy delta [temperature]
variable us    equal f_myhug[2] # shock velocity [distance/time]
variable up    equal f_myhug[3] # particle velocity [distance/time]
variable pzz   equal pzz  # axial stress
variable tau equal 0.5*(pzz-0.5*(pxx+pyy)) # shear stress
variable time equal dt*step

thermo          1000
#thermo_style custom step temp ke epair etotal pzz v_tau lz f_myhug v_dele v_us v_up
thermo_style custom step temp ke epair etotal pzz  lz v_up

fix pr1 all print 10 "${time} ${us} ${up} " screen no append us-up-t-nodrag.dat title '#time us up (no drag)'


fix stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (no drag)'

#dump		id all atom 5000 dump.hugoniostat

#dump		2 all image 500 image.*.jpg type type &
#		axes yes 0.8 0.02 view 60 -30
#dump_modify	2 pad 5

#dump		3 all movie 500 movie.mpg type type &
#		axes yes 0.8 0.02 view 60 -30
#dump_modify	3 pad 5

run		100000

# Start Run #2

# Start Run #3

log log.nhchains

clear
read_restart restart.equil

neighbor        0.2 bin
neigh_modify    every 1 delay 0 check yes
timestep        0.001
reset_timestep 0

# Pzz = 40.0, drag/damping term off, Nose-Hoover chains

fix myhug all nphug temp 300.0 300.0 10.0 z 2000000.0 2000000.0 70.0

# Specify reference state from paper, times 1000 atoms

#fix_modify myhug e0 -6334.0 p0 0.0 v0 680.73519
fix_modify myhug e0 -129600.0 p0 0.0 v0 288542.70

# Add fix energy to ouput etotal
 
fix_modify myhug energy yes 

# Define output

variable dele  equal f_myhug[1] # energy delta [temperature]
variable us    equal f_myhug[2] # shock velocity [distance/time]
variable up    equal f_myhug[3] # particle velocity [distance/time]
variable pzz   equal pzz  # axial stress
variable tau equal 0.5*(pzz-0.5*(pxx+pyy)) # shear stress
variable time equal dt*step

thermo          1000
thermo_style custom step temp ke epair etotal pzz v_tau lz f_myhug v_dele v_us v_up
fix pr3 all print 10 "${time} ${us} ${up} " screen no append us-up-t-Nosehoover.dat title '#time us up (Nose-Hoover chain)'
fix stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (Nose-Hoover chain)'

run		100000