forked from lijiext/lammps
178 lines
4.5 KiB
Plaintext
178 lines
4.5 KiB
Plaintext
# This script reproduces stress trajectories from Fig. 1 in
|
|
# Ravelo, Holian, Germann, and Lomdahl, PRB 70 014103 (2004)
|
|
#
|
|
# Three thermostatting scenarios are visited: undamped (nodrag),
|
|
# damped (drag) and Nose-Hoover chain (nhchains).
|
|
#
|
|
# The axial and shear stress trajectories are printed to the
|
|
# file "stress_vs_t.dat". For the damped case, the original figure
|
|
# seems to be a plot of 2*tau, rather than tau.
|
|
#
|
|
# The script also demonstrates how to
|
|
# orient a crystal along <110>,
|
|
# and how to use the lj/cubic pair style.
|
|
|
|
units lj
|
|
boundary p p p
|
|
|
|
atom_style atomic
|
|
|
|
# Set up FCC lattice with z axis along <110>
|
|
|
|
lattice fcc 1.4142136 &
|
|
orient x 0 0 1 &
|
|
orient y 1 -1 0 &
|
|
orient z 1 1 0
|
|
|
|
region mycell block 0.0 5.0 0.0 5.0 0.0 5.0 units lattice
|
|
create_box 1 mycell
|
|
mass * 1.0
|
|
create_atoms 1 box
|
|
|
|
# Using units of Rmin, so sigma = 2^-1/6 = 0.8908987
|
|
|
|
pair_style lj/cubic
|
|
pair_coeff * * 1.0 0.8908987
|
|
|
|
# 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 0.01 87287 mom yes rot yes dist gaussian
|
|
write_restart restart.equil
|
|
|
|
# Start Run #1
|
|
|
|
log log.nodrag
|
|
|
|
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
|
|
|
|
fix myhug all nphug temp 1.0 1.0 10.0 z 40.0 40.0 70.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
|
|
|
|
# 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 stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (no drag)'
|
|
|
|
#dump id all atom 500 dump.hugoniostat
|
|
|
|
#dump 1 all image 500 image.*.jpg type type &
|
|
# axes yes 0.8 0.02 view 60 -30
|
|
#dump_modify 1 pad 5
|
|
|
|
run 10000
|
|
|
|
# Start Run #2
|
|
|
|
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 1.0 1.0 1.0 z 40.0 40.0 70.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
|
|
|
|
# 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 stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (with drag)'
|
|
|
|
run 10000
|
|
|
|
# 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 1.0 1.0 1.0 z 40.0 40.0 70.0
|
|
|
|
# Specify reference state from paper, times 1000 atoms
|
|
|
|
fix_modify myhug e0 -6334.0 p0 0.0 v0 680.73519
|
|
|
|
# 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 stress all print 10 "${time} ${pzz} ${tau} " screen no append stress_vs_t.dat title '#time pzz tau (Nose-Hoover chain)'
|
|
|
|
run 10000
|
|
|