From 790d7d9fae8c6cc6e9bb4ad54953550a769df608 Mon Sep 17 00:00:00 2001
From: Aidan Thompson <athomps@sandia.gov>
Date: Fri, 13 Sep 2019 15:59:25 -0600
Subject: [PATCH] Added log files and updated README

---
 examples/README                             |   1 +
 examples/steinhardt/log.13Sept18.bcc.g++.1  | 176 +++++++++++++++
 examples/steinhardt/log.13Sept18.bcc.g++.4  | 176 +++++++++++++++
 examples/steinhardt/log.13Sept18.fcc.g++.1  | 172 ++++++++++++++
 examples/steinhardt/log.13Sept18.fcc.g++.4  | 172 ++++++++++++++
 examples/steinhardt/log.13Sept18.icos.g++.1 | 235 ++++++++++++++++++++
 examples/steinhardt/log.13Sept18.icos.g++.4 | 235 ++++++++++++++++++++
 7 files changed, 1167 insertions(+)
 create mode 100644 examples/steinhardt/log.13Sept18.bcc.g++.1
 create mode 100644 examples/steinhardt/log.13Sept18.bcc.g++.4
 create mode 100644 examples/steinhardt/log.13Sept18.fcc.g++.1
 create mode 100644 examples/steinhardt/log.13Sept18.fcc.g++.4
 create mode 100644 examples/steinhardt/log.13Sept18.icos.g++.1
 create mode 100644 examples/steinhardt/log.13Sept18.icos.g++.4

diff --git a/examples/README b/examples/README
index dbfdb3363b..46148aea8b 100644
--- a/examples/README
+++ b/examples/README
@@ -105,6 +105,7 @@ shear:    sideways shear applied to 2d solid, with and without a void
 snap:     use of SNAP potential for Ta
 srd:      stochastic rotation dynamics (SRD) particles as solvent
 snap:     NVE dynamics for BCC tantalum crystal using SNAP potential
+steinhardt: Steinhardt-Nelson Q_l and W_l parameters usng orientorder/atom
 streitz:  Streitz-Mintmire potential for Al2O3
 tad:      temperature-accelerated dynamics of vacancy diffusion in bulk Si
 threebody: regression test input for a variety of manybody potentials
diff --git a/examples/steinhardt/log.13Sept18.bcc.g++.1 b/examples/steinhardt/log.13Sept18.bcc.g++.1
new file mode 100644
index 0000000000..b0a82fca90
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.bcc.g++.1
@@ -0,0 +1,176 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for BCC
+
+variable  	rcut equal 3.0
+
+boundary	p p p
+
+atom_style	atomic
+neighbor	0.3 bin
+neigh_modify	delay 5
+
+# create geometry
+
+lattice		bcc 1.0
+Lattice spacing in x,y,z = 1.25992 1.25992 1.25992
+region		box block 0 3 0 3 0 3
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.77976 3.77976 3.77976)
+  1 by 1 by 1 MPI processor grid
+create_atoms	1 box
+Created 54 atoms
+  create_atoms CPU = 0.000282049 secs
+
+mass		1 1.0
+
+# LJ potentials
+
+pair_style	lj/cut ${rcut}
+pair_style	lj/cut 3
+pair_coeff	* * 1.0 1.0 ${rcut}
+pair_coeff	* * 1.0 1.0 3
+
+# 14 neighbors, perfect crystal
+
+compute 	qlwlhat all orientorder/atom degrees 6 2 4 6 8 10 12 nnn 14 wl/hat yes
+compute 	avql all reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5] c_qlwlhat[6]
+compute 	avwlhat all reduce ave c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10] c_qlwlhat[11] c_qlwlhat[12]
+
+thermo_style    custom step temp epair etotal c_avql[*] c_avwlhat[*]
+
+run		0
+WARNING: No fixes defined, atoms won't move (../verlet.cpp:52)
+Neighbor list info ...
+  update every 1 steps, delay 5 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 3.3
+  ghost atom cutoff = 3.3
+  binsize = 1.65, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 3.746 | 3.746 | 3.746 Mbytes
+Step Temp E_pair TotEng c_avql[1] c_avql[2] c_avql[3] c_avql[4] c_avql[5] c_avql[6] c_avwlhat[1] c_avwlhat[2] c_avwlhat[3] c_avwlhat[4] c_avwlhat[5] c_avwlhat[6] 
+       0            0   -7.8104466   -7.8104466 8.6570408e-17  0.036369648   0.51068823   0.42932247   0.19519122   0.40479919            0   0.15931737  0.013160601  0.058454791 -0.090130212 -0.049573639 
+Loop time of 9.53674e-07 on 1 procs for 0 steps with 54 atoms
+
+209.7% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 9.537e-07  |            |       |100.00
+
+Nlocal:    54 ave 54 max 54 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    1187 ave 1187 max 1187 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    4536 ave 4536 max 4536 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  9072 ave 9072 max 9072 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 9072
+Ave neighs/atom = 168
+Neighbor list builds = 0
+Dangerous builds = 0
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q2ref equal 0.0
+variable 	q4ref equal 0.036
+variable 	q6ref equal 0.511
+variable 	q8ref equal 0.429
+variable 	q10ref equal 0.195
+variable 	q12ref equal 0.405
+
+variable 	q2 equal c_avql[1]
+variable 	q4 equal c_avql[2]
+variable 	q6 equal c_avql[3]
+variable 	q8 equal c_avql[4]
+variable 	q10 equal c_avql[5]
+variable 	q12 equal c_avql[6]
+
+print "q2     = $(v_q2:%10.6f) delta = $(v_q2-v_q2ref:%10.4f)"
+q2     =   0.000000 delta =     0.0000
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.036370 delta =     0.0004
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.510688 delta =    -0.0003
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.429322 delta =     0.0003
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.195191 delta =     0.0002
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.404799 delta =    -0.0002
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w4hatref equal 0.159317
+variable 	w6hatref equal 0.013161
+variable 	w8hatref equal -0.058455
+variable 	w10hatref equal -0.090130
+
+variable 	w4hat equal c_avwlhat[2]
+variable 	w6hat equal c_avwlhat[3]
+variable 	w8hat equal c_avwlhat[4]
+variable 	w10hat equal c_avwlhat[5]
+
+print "w4hat  = $(v_w4hat:%10.6f) delta = $(v_w4hat-v_w4hatref:%10.6f)"
+w4hat  =   0.159317 delta =   0.000000
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =   0.013161 delta =  -0.000000
+print "w8hat  = $(v_w8hat:%10.6f) delta = $(v_w8hat-v_w8hatref:%10.6f)"
+w8hat  =   0.058455 delta =   0.116910
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.090130 delta =  -0.000000
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+Total wall time: 0:00:00
diff --git a/examples/steinhardt/log.13Sept18.bcc.g++.4 b/examples/steinhardt/log.13Sept18.bcc.g++.4
new file mode 100644
index 0000000000..ee1f0c69c1
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.bcc.g++.4
@@ -0,0 +1,176 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for BCC
+
+variable  	rcut equal 3.0
+
+boundary	p p p
+
+atom_style	atomic
+neighbor	0.3 bin
+neigh_modify	delay 5
+
+# create geometry
+
+lattice		bcc 1.0
+Lattice spacing in x,y,z = 1.25992 1.25992 1.25992
+region		box block 0 3 0 3 0 3
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.77976 3.77976 3.77976)
+  1 by 2 by 2 MPI processor grid
+create_atoms	1 box
+Created 54 atoms
+  create_atoms CPU = 0.000533104 secs
+
+mass		1 1.0
+
+# LJ potentials
+
+pair_style	lj/cut ${rcut}
+pair_style	lj/cut 3
+pair_coeff	* * 1.0 1.0 ${rcut}
+pair_coeff	* * 1.0 1.0 3
+
+# 14 neighbors, perfect crystal
+
+compute 	qlwlhat all orientorder/atom degrees 6 2 4 6 8 10 12 nnn 14 wl/hat yes
+compute 	avql all reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5] c_qlwlhat[6]
+compute 	avwlhat all reduce ave c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10] c_qlwlhat[11] c_qlwlhat[12]
+
+thermo_style    custom step temp epair etotal c_avql[*] c_avwlhat[*]
+
+run		0
+WARNING: No fixes defined, atoms won't move (../verlet.cpp:52)
+Neighbor list info ...
+  update every 1 steps, delay 5 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 3.3
+  ghost atom cutoff = 3.3
+  binsize = 1.65, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 3.722 | 3.722 | 3.722 Mbytes
+Step Temp E_pair TotEng c_avql[1] c_avql[2] c_avql[3] c_avql[4] c_avql[5] c_avql[6] c_avwlhat[1] c_avwlhat[2] c_avwlhat[3] c_avwlhat[4] c_avwlhat[5] c_avwlhat[6] 
+       0            0   -7.8104466   -7.8104466 8.6382997e-17  0.036369648   0.51068823   0.42932247   0.19519122   0.40479919            0   0.15931737  0.013160601  0.058454791 -0.090130212 -0.049573639 
+Loop time of 3.99351e-06 on 4 procs for 0 steps with 54 atoms
+
+87.6% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 3.994e-06  |            |       |100.00
+
+Nlocal:    13.5 ave 15 max 12 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Nghost:    819.5 ave 821 max 818 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+Neighs:    1134 ave 1260 max 1008 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+FullNghs:  2268 ave 2520 max 2016 min
+Histogram: 2 0 0 0 0 0 0 0 0 2
+
+Total # of neighbors = 9072
+Ave neighs/atom = 168
+Neighbor list builds = 0
+Dangerous builds = 0
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q2ref equal 0.0
+variable 	q4ref equal 0.036
+variable 	q6ref equal 0.511
+variable 	q8ref equal 0.429
+variable 	q10ref equal 0.195
+variable 	q12ref equal 0.405
+
+variable 	q2 equal c_avql[1]
+variable 	q4 equal c_avql[2]
+variable 	q6 equal c_avql[3]
+variable 	q8 equal c_avql[4]
+variable 	q10 equal c_avql[5]
+variable 	q12 equal c_avql[6]
+
+print "q2     = $(v_q2:%10.6f) delta = $(v_q2-v_q2ref:%10.4f)"
+q2     =   0.000000 delta =     0.0000
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.036370 delta =     0.0004
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.510688 delta =    -0.0003
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.429322 delta =     0.0003
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.195191 delta =     0.0002
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.404799 delta =    -0.0002
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w4hatref equal 0.159317
+variable 	w6hatref equal 0.013161
+variable 	w8hatref equal -0.058455
+variable 	w10hatref equal -0.090130
+
+variable 	w4hat equal c_avwlhat[2]
+variable 	w6hat equal c_avwlhat[3]
+variable 	w8hat equal c_avwlhat[4]
+variable 	w10hat equal c_avwlhat[5]
+
+print "w4hat  = $(v_w4hat:%10.6f) delta = $(v_w4hat-v_w4hatref:%10.6f)"
+w4hat  =   0.159317 delta =   0.000000
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =   0.013161 delta =  -0.000000
+print "w8hat  = $(v_w8hat:%10.6f) delta = $(v_w8hat-v_w8hatref:%10.6f)"
+w8hat  =   0.058455 delta =   0.116910
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.090130 delta =  -0.000000
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+Total wall time: 0:00:00
diff --git a/examples/steinhardt/log.13Sept18.fcc.g++.1 b/examples/steinhardt/log.13Sept18.fcc.g++.1
new file mode 100644
index 0000000000..a34a33d7b1
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.fcc.g++.1
@@ -0,0 +1,172 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for FCC
+
+variable  	rcut equal 3.0
+
+boundary	p p p
+
+atom_style	atomic
+neighbor	0.3 bin
+neigh_modify	delay 5
+
+# create geometry
+
+lattice		fcc 1.0
+Lattice spacing in x,y,z = 1.5874 1.5874 1.5874
+region		box block 0 3 0 3 0 3
+create_box	1 box
+Created orthogonal box = (0 0 0) to (4.7622 4.7622 4.7622)
+  1 by 1 by 1 MPI processor grid
+create_atoms	1 box
+Created 108 atoms
+  create_atoms CPU = 0.000289917 secs
+
+mass		1 1.0
+
+# LJ potentials
+
+pair_style	lj/cut ${rcut}
+pair_style	lj/cut 3
+pair_coeff	* * 1.0 1.0 ${rcut}
+pair_coeff	* * 1.0 1.0 3
+
+# 12 neighbors, perfect crystal
+
+compute 	qlwlhat all orientorder/atom wl/hat yes
+compute 	avql all reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5]
+compute 	avwlhat all reduce ave c_qlwlhat[6] c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10]
+
+thermo_style    custom step temp epair etotal c_avql[*] c_avwlhat[*]
+
+run		0
+WARNING: No fixes defined, atoms won't move (../verlet.cpp:52)
+Neighbor list info ...
+  update every 1 steps, delay 5 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 3.3
+  ghost atom cutoff = 3.3
+  binsize = 1.65, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 3.767 | 3.767 | 3.767 Mbytes
+Step Temp E_pair TotEng c_avql[1] c_avql[2] c_avql[3] c_avql[4] c_avql[5] c_avwlhat[1] c_avwlhat[2] c_avwlhat[3] c_avwlhat[4] c_avwlhat[5] 
+       0            0   -8.1295091   -8.1295091   0.19094065   0.57452426   0.40391456  0.012857043   0.60008302  -0.15931737 -0.013160601  0.058454791 -0.090130212  0.087390889 
+Loop time of 9.53674e-07 on 1 procs for 0 steps with 108 atoms
+
+209.7% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 9.537e-07  |            |       |100.00
+
+Nlocal:    108 ave 108 max 108 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    1580 ave 1580 max 1580 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    7560 ave 7560 max 7560 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  15120 ave 15120 max 15120 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 15120
+Ave neighs/atom = 140
+Neighbor list builds = 0
+Dangerous builds = 0
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q4ref equal 0.190
+variable 	q6ref equal 0.575
+variable 	q8ref equal 0.404
+variable 	q10ref equal 0.013
+variable 	q12ref equal 0.600
+
+variable 	q4 equal c_avql[1]
+variable 	q6 equal c_avql[2]
+variable 	q8 equal c_avql[3]
+variable 	q10 equal c_avql[4]
+variable 	q12 equal c_avql[5]
+
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.190941 delta =     0.0009
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.574524 delta =    -0.0005
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.403915 delta =    -0.0001
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.012857 delta =    -0.0001
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.600083 delta =     0.0001
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w4hatref equal -0.159316
+variable 	w6hatref equal -0.013161
+variable 	w8hatref equal 0.058454
+variable 	w10hatref equal -0.090128
+
+variable 	w4hat equal c_avwlhat[1]
+variable 	w6hat equal c_avwlhat[2]
+variable 	w8hat equal c_avwlhat[3]
+variable 	w10hat equal c_avwlhat[4]
+
+print "w4hat  = $(v_w4hat:%10.6f) delta = $(v_w4hat-v_w4hatref:%10.6f)"
+w4hat  =  -0.159317 delta =  -0.000001
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =  -0.013161 delta =   0.000000
+print "w8hat  = $(v_w8hat:%10.6f) delta = $(v_w8hat-v_w8hatref:%10.6f)"
+w8hat  =   0.058455 delta =   0.000001
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.090130 delta =  -0.000002
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+Total wall time: 0:00:00
diff --git a/examples/steinhardt/log.13Sept18.fcc.g++.4 b/examples/steinhardt/log.13Sept18.fcc.g++.4
new file mode 100644
index 0000000000..8573488be6
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.fcc.g++.4
@@ -0,0 +1,172 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for FCC
+
+variable  	rcut equal 3.0
+
+boundary	p p p
+
+atom_style	atomic
+neighbor	0.3 bin
+neigh_modify	delay 5
+
+# create geometry
+
+lattice		fcc 1.0
+Lattice spacing in x,y,z = 1.5874 1.5874 1.5874
+region		box block 0 3 0 3 0 3
+create_box	1 box
+Created orthogonal box = (0 0 0) to (4.7622 4.7622 4.7622)
+  1 by 2 by 2 MPI processor grid
+create_atoms	1 box
+Created 108 atoms
+  create_atoms CPU = 0.000549078 secs
+
+mass		1 1.0
+
+# LJ potentials
+
+pair_style	lj/cut ${rcut}
+pair_style	lj/cut 3
+pair_coeff	* * 1.0 1.0 ${rcut}
+pair_coeff	* * 1.0 1.0 3
+
+# 12 neighbors, perfect crystal
+
+compute 	qlwlhat all orientorder/atom wl/hat yes
+compute 	avql all reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5]
+compute 	avwlhat all reduce ave c_qlwlhat[6] c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10]
+
+thermo_style    custom step temp epair etotal c_avql[*] c_avwlhat[*]
+
+run		0
+WARNING: No fixes defined, atoms won't move (../verlet.cpp:52)
+Neighbor list info ...
+  update every 1 steps, delay 5 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 3.3
+  ghost atom cutoff = 3.3
+  binsize = 1.65, bins = 3 3 3
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 3.737 | 3.737 | 3.737 Mbytes
+Step Temp E_pair TotEng c_avql[1] c_avql[2] c_avql[3] c_avql[4] c_avql[5] c_avwlhat[1] c_avwlhat[2] c_avwlhat[3] c_avwlhat[4] c_avwlhat[5] 
+       0            0   -8.1295091   -8.1295091   0.19094065   0.57452426   0.40391456  0.012857043   0.60008302  -0.15931737 -0.013160601  0.058454791 -0.090130212  0.087390889 
+Loop time of 3.09944e-06 on 4 procs for 0 steps with 108 atoms
+
+88.7% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0          | 0          | 0          |   0.0 |  0.00
+Neigh   | 0          | 0          | 0          |   0.0 |  0.00
+Comm    | 0          | 0          | 0          |   0.0 |  0.00
+Output  | 0          | 0          | 0          |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 3.099e-06  |            |       |100.00
+
+Nlocal:    27 ave 27 max 27 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+Nghost:    1053 ave 1053 max 1053 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+Neighs:    1890 ave 1890 max 1890 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+FullNghs:  3780 ave 3780 max 3780 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 15120
+Ave neighs/atom = 140
+Neighbor list builds = 0
+Dangerous builds = 0
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q4ref equal 0.190
+variable 	q6ref equal 0.575
+variable 	q8ref equal 0.404
+variable 	q10ref equal 0.013
+variable 	q12ref equal 0.600
+
+variable 	q4 equal c_avql[1]
+variable 	q6 equal c_avql[2]
+variable 	q8 equal c_avql[3]
+variable 	q10 equal c_avql[4]
+variable 	q12 equal c_avql[5]
+
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.190941 delta =     0.0009
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.574524 delta =    -0.0005
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.403915 delta =    -0.0001
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.012857 delta =    -0.0001
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.600083 delta =     0.0001
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w4hatref equal -0.159316
+variable 	w6hatref equal -0.013161
+variable 	w8hatref equal 0.058454
+variable 	w10hatref equal -0.090128
+
+variable 	w4hat equal c_avwlhat[1]
+variable 	w6hat equal c_avwlhat[2]
+variable 	w8hat equal c_avwlhat[3]
+variable 	w10hat equal c_avwlhat[4]
+
+print "w4hat  = $(v_w4hat:%10.6f) delta = $(v_w4hat-v_w4hatref:%10.6f)"
+w4hat  =  -0.159317 delta =  -0.000001
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =  -0.013161 delta =   0.000000
+print "w8hat  = $(v_w8hat:%10.6f) delta = $(v_w8hat-v_w8hatref:%10.6f)"
+w8hat  =   0.058455 delta =   0.000001
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.090130 delta =  -0.000002
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+Total wall time: 0:00:00
diff --git a/examples/steinhardt/log.13Sept18.icos.g++.1 b/examples/steinhardt/log.13Sept18.icos.g++.1
new file mode 100644
index 0000000000..f42fd346db
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.icos.g++.1
@@ -0,0 +1,235 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for icosahedral cluster
+# W_6_hat is sensitive to icosohedral order
+
+variable        rcut equal 1.2 # a bit bigger than LJ Rmin
+variable        rcutred equal 0.75 # a bit bigger than 1/sqrt(2)
+
+# create a perfect fcc crystallite
+
+atom_style	atomic
+boundary	s s s
+lattice		fcc 1.0 # neighbors at LJ Rmin
+Lattice spacing in x,y,z = 1.5874 1.5874 1.5874
+region		box block 0 2 0 2 0 2
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.1748 3.1748 3.1748)
+  1 by 1 by 1 MPI processor grid
+create_atoms	1 box
+Created 63 atoms
+  create_atoms CPU = 0.000341177 secs
+mass		1 1.0
+
+region 		centralatom sphere 1 1 1 0.0 side in
+group 		centralatom region centralatom
+1 atoms in group centralatom
+
+region 		mysphere sphere 1 1 1 ${rcutred} side out
+region 		mysphere sphere 1 1 1 0.75 side out
+delete_atoms 	region mysphere
+Deleted 50 atoms, new total = 13
+
+# LJ potential
+
+pair_style	lj/cut 100.0
+pair_coeff	* * 1.0 1.0 100.0
+
+# define output for central atom
+
+compute 	qlwlhat all orientorder/atom wl/hat yes cutoff ${rcut} nnn NULL
+compute 	qlwlhat all orientorder/atom wl/hat yes cutoff 1.2 nnn NULL
+compute 	avql centralatom reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5]
+compute 	avwlhat centralatom reduce ave c_qlwlhat[6] c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10]
+variable	q6 equal c_avql[2]
+variable	w6hat equal c_avwlhat[2]
+
+compute 	mype all pe/atom
+compute 	centralatompe centralatom reduce ave c_mype
+
+# gently equilibrate the crystallite
+
+velocity     	all create 0.001 482748
+fix	     	1 all nve
+neighbor	0.3 bin
+neigh_modify	every 1 check no delay 0
+timestep	0.003
+thermo_style    custom step temp epair etotal c_centralatompe v_q6 v_w6hat
+thermo		10
+
+run		10
+Neighbor list info ...
+  update every 1 steps, delay 0 steps, check no
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 100.3
+  ghost atom cutoff = 100.3
+  binsize = 50.15, bins = 1 1 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 60.87 | 60.87 | 60.87 Mbytes
+Step Temp E_pair TotEng c_centralatompe v_q6 v_w6hat 
+       0        0.001    -3.134107   -3.1327224           -6   0.57452426 -0.013160601 
+      10 0.0021974351   -3.1357656    -3.132723   -5.9995795   0.57450739 -0.013160482 
+Loop time of 0.210631 on 1 procs for 10 steps with 13 atoms
+
+Performance: 12305.887 tau/day, 47.476 timesteps/s
+99.5% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 1.6928e-05 | 1.6928e-05 | 1.6928e-05 |   0.0 |  0.01
+Neigh   | 0.13335    | 0.13335    | 0.13335    |   0.0 | 63.31
+Comm    | 1.4782e-05 | 1.4782e-05 | 1.4782e-05 |   0.0 |  0.01
+Output  | 0.033121   | 0.033121   | 0.033121   |   0.0 | 15.72
+Modify  | 2.8849e-05 | 2.8849e-05 | 2.8849e-05 |   0.0 |  0.01
+Other   |            | 0.0441     |            |       | 20.94
+
+Nlocal:    13 ave 13 max 13 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    78 ave 78 max 78 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  156 ave 156 max 156 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 156
+Ave neighs/atom = 12
+Neighbor list builds = 10
+Dangerous builds not checked
+
+# quench to icosehedral cluster
+
+minimize      	1.0e-10 1.0e-6 100 1000
+WARNING: Using 'neigh_modify every 1 delay 0 check yes' setting during minimization (../min.cpp:168)
+Per MPI rank memory allocation (min/avg/max) = 113 | 113 | 113 Mbytes
+Step Temp E_pair TotEng c_centralatompe v_q6 v_w6hat 
+      10 0.0021974351   -3.1357656    -3.132723   -5.9995795   0.57450739 -0.013160482 
+      20 0.0021974351   -3.1449631   -3.1419205   -5.9766502   0.57452794  -0.01317299 
+      30 0.0021974351    -3.377441   -3.3743984   -5.6930377   0.65479437  -0.16922776 
+      40 0.0021974351   -3.4096335   -3.4065909   -5.6325443   0.66328926  -0.16975382 
+      49 0.0021974351    -3.409754   -3.4067114   -5.6323333   0.66332496  -0.16975389 
+Loop time of 0.0821278 on 1 procs for 39 steps with 13 atoms
+
+99.9% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+Minimization stats:
+  Stopping criterion = energy tolerance
+  Energy initial, next-to-last, final = 
+        -3.13576562743     -3.40975395481     -3.40975395503
+  Force two-norm initial, final = 6.44841 0.000945077
+  Force max component initial, final = 1.43234 0.000348946
+  Final line search alpha, max atom move = 1 0.000348946
+  Iterations, force evaluations = 39 94
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0.00020051 | 0.00020051 | 0.00020051 |   0.0 |  0.24
+Neigh   | 0.024762   | 0.024762   | 0.024762   |   0.0 | 30.15
+Comm    | 1.2398e-05 | 1.2398e-05 | 1.2398e-05 |   0.0 |  0.02
+Output  | 0.048263   | 0.048263   | 0.048263   |   0.0 | 58.77
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 0.00889    |            |       | 10.82
+
+Nlocal:    13 ave 13 max 13 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    78 ave 78 max 78 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  156 ave 156 max 156 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 156
+Ave neighs/atom = 12
+Neighbor list builds = 3
+Dangerous builds not checked
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q4ref equal 0.0
+variable 	q6ref equal 0.663
+variable 	q8ref equal 0.0
+variable 	q10ref equal 0.363
+variable 	q12ref equal 0.585
+
+variable 	q4 equal c_avql[1]
+variable 	q6 equal c_avql[2]
+variable 	q8 equal c_avql[3]
+variable 	q10 equal c_avql[4]
+variable 	q12 equal c_avql[5]
+
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.000002 delta =     0.0000
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.663325 delta =     0.0003
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.000003 delta =     0.0000
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.362951 delta =    -0.0000
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.585423 delta =     0.0004
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w6hatref equal -0.169754
+variable 	w10hatref equal -0.093967
+
+variable 	w4hat equal c_avwlhat[1]
+variable 	w6hat equal c_avwlhat[2]
+variable 	w8hat equal c_avwlhat[3]
+variable 	w10hat equal c_avwlhat[4]
+variable 	w12hat equal c_avwlhat[5]
+
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =  -0.169754 delta =   0.000000
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.093968 delta =  -0.000001
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+
+Total wall time: 0:00:00
diff --git a/examples/steinhardt/log.13Sept18.icos.g++.4 b/examples/steinhardt/log.13Sept18.icos.g++.4
new file mode 100644
index 0000000000..d64e179c8f
--- /dev/null
+++ b/examples/steinhardt/log.13Sept18.icos.g++.4
@@ -0,0 +1,235 @@
+LAMMPS (7 Aug 2019)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:93)
+  using 1 OpenMP thread(s) per MPI task
+# Steinhardt-Nelson bond orientational order parameters for icosahedral cluster
+# W_6_hat is sensitive to icosohedral order
+
+variable        rcut equal 1.2 # a bit bigger than LJ Rmin
+variable        rcutred equal 0.75 # a bit bigger than 1/sqrt(2)
+
+# create a perfect fcc crystallite
+
+atom_style	atomic
+boundary	s s s
+lattice		fcc 1.0 # neighbors at LJ Rmin
+Lattice spacing in x,y,z = 1.5874 1.5874 1.5874
+region		box block 0 2 0 2 0 2
+create_box	1 box
+Created orthogonal box = (0 0 0) to (3.1748 3.1748 3.1748)
+  1 by 2 by 2 MPI processor grid
+create_atoms	1 box
+Created 63 atoms
+  create_atoms CPU = 0.000592947 secs
+mass		1 1.0
+
+region 		centralatom sphere 1 1 1 0.0 side in
+group 		centralatom region centralatom
+1 atoms in group centralatom
+
+region 		mysphere sphere 1 1 1 ${rcutred} side out
+region 		mysphere sphere 1 1 1 0.75 side out
+delete_atoms 	region mysphere
+Deleted 50 atoms, new total = 13
+
+# LJ potential
+
+pair_style	lj/cut 100.0
+pair_coeff	* * 1.0 1.0 100.0
+
+# define output for central atom
+
+compute 	qlwlhat all orientorder/atom wl/hat yes cutoff ${rcut} nnn NULL
+compute 	qlwlhat all orientorder/atom wl/hat yes cutoff 1.2 nnn NULL
+compute 	avql centralatom reduce ave c_qlwlhat[1] c_qlwlhat[2] c_qlwlhat[3] c_qlwlhat[4] c_qlwlhat[5]
+compute 	avwlhat centralatom reduce ave c_qlwlhat[6] c_qlwlhat[7] c_qlwlhat[8] c_qlwlhat[9] c_qlwlhat[10]
+variable	q6 equal c_avql[2]
+variable	w6hat equal c_avwlhat[2]
+
+compute 	mype all pe/atom
+compute 	centralatompe centralatom reduce ave c_mype
+
+# gently equilibrate the crystallite
+
+velocity     	all create 0.001 482748
+fix	     	1 all nve
+neighbor	0.3 bin
+neigh_modify	every 1 check no delay 0
+timestep	0.003
+thermo_style    custom step temp epair etotal c_centralatompe v_q6 v_w6hat
+thermo		10
+
+run		10
+Neighbor list info ...
+  update every 1 steps, delay 0 steps, check no
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 100.3
+  ghost atom cutoff = 100.3
+  binsize = 50.15, bins = 1 1 1
+  2 neighbor lists, perpetual/occasional/extra = 1 1 0
+  (1) pair lj/cut, perpetual
+      attributes: half, newton on
+      pair build: half/bin/atomonly/newton
+      stencil: half/bin/3d/newton
+      bin: standard
+  (2) compute orientorder/atom, occasional
+      attributes: full, newton on
+      pair build: full/bin/atomonly
+      stencil: full/bin/3d
+      bin: standard
+Per MPI rank memory allocation (min/avg/max) = 60.76 | 60.76 | 60.76 Mbytes
+Step Temp E_pair TotEng c_centralatompe v_q6 v_w6hat 
+       0        0.001    -3.134107   -3.1327224           -6   0.57452426 -0.013160601 
+      10 0.0021821015   -3.1357444    -3.132723   -5.9994885   0.57450756 -0.013160431 
+Loop time of 0.15387 on 4 procs for 10 steps with 13 atoms
+
+Performance: 16845.410 tau/day, 64.990 timesteps/s
+99.2% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 9.7752e-06 | 1.7285e-05 | 2.9087e-05 |   0.0 |  0.01
+Neigh   | 0.055218   | 0.065482   | 0.073829   |   2.9 | 42.56
+Comm    | 0.0061202  | 0.018279   | 0.028283   |   6.3 | 11.88
+Output  | 0.020282   | 0.020287   | 0.020298   |   0.0 | 13.18
+Modify  | 2.7895e-05 | 2.8968e-05 | 2.9802e-05 |   0.0 |  0.02
+Other   |            | 0.04978    |            |       | 32.35
+
+Nlocal:    3.25 ave 4 max 2 min
+Histogram: 1 0 0 0 0 1 0 0 0 2
+Nghost:    9.75 ave 11 max 9 min
+Histogram: 2 0 0 0 0 1 0 0 0 1
+Neighs:    19.5 ave 35 max 7 min
+Histogram: 2 0 0 0 0 0 0 1 0 1
+FullNghs:  39 ave 48 max 24 min
+Histogram: 1 0 0 0 0 1 0 0 0 2
+
+Total # of neighbors = 156
+Ave neighs/atom = 12
+Neighbor list builds = 10
+Dangerous builds not checked
+
+# quench to icosehedral cluster
+
+minimize      	1.0e-10 1.0e-6 100 1000
+WARNING: Using 'neigh_modify every 1 delay 0 check yes' setting during minimization (../min.cpp:168)
+Per MPI rank memory allocation (min/avg/max) = 112.9 | 112.9 | 112.9 Mbytes
+Step Temp E_pair TotEng c_centralatompe v_q6 v_w6hat 
+      10 0.0021821015   -3.1357444    -3.132723   -5.9994885   0.57450756 -0.013160431 
+      20 0.0021821015   -3.1449613   -3.1419399   -5.9764731   0.57452454  -0.01316152 
+      30 0.0021821015   -3.3366586   -3.3336372   -5.7717004   0.63679987  -0.16411081 
+      40 0.0021821015   -3.4097262   -3.4067048   -5.6321229   0.66331761  -0.16975374 
+      49 0.0021821015    -3.409754   -3.4067326   -5.6323347   0.66332496  -0.16975389 
+Loop time of 0.0932837 on 4 procs for 39 steps with 13 atoms
+
+97.2% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+Minimization stats:
+  Stopping criterion = energy tolerance
+  Energy initial, next-to-last, final = 
+        -3.13574438249     -3.40975395519     -3.40975395529
+  Force two-norm initial, final = 6.46606 0.000429142
+  Force max component initial, final = 1.45519 0.000196446
+  Final line search alpha, max atom move = 1 0.000196446
+  Iterations, force evaluations = 39 98
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 5.4121e-05 | 9.0897e-05 | 0.00014281 |   0.0 |  0.10
+Neigh   | 0.019662   | 0.02379    | 0.03176    |   3.0 | 25.50
+Comm    | 0.0072601  | 0.014304   | 0.019575   |   3.8 | 15.33
+Output  | 0.033646   | 0.033654   | 0.03368    |   0.0 | 36.08
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 0.02144    |            |       | 22.99
+
+Nlocal:    3.25 ave 6 max 1 min
+Histogram: 1 0 0 0 2 0 0 0 0 1
+Nghost:    9.75 ave 12 max 7 min
+Histogram: 1 0 0 0 0 0 2 0 0 1
+Neighs:    19.5 ave 49 max 1 min
+Histogram: 2 0 0 0 1 0 0 0 0 1
+FullNghs:  39 ave 72 max 12 min
+Histogram: 1 0 0 0 2 0 0 0 0 1
+
+Total # of neighbors = 156
+Ave neighs/atom = 12
+Neighbor list builds = 5
+Dangerous builds not checked
+
+# check Q_l values
+
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+print "Comparison with reference values of Q_l "
+Comparison with reference values of Q_l 
+print "   [Table I in W. Mickel, S. C. Kapfer,"
+   [Table I in W. Mickel, S. C. Kapfer,
+print "   G. E. Schroeder-Turkand, K. Mecke, "
+   G. E. Schroeder-Turkand, K. Mecke, 
+print "   J. Chem. Phys. 138, 044501 (2013).]"
+   J. Chem. Phys. 138, 044501 (2013).]
+print " "
+ 
+
+variable 	q4ref equal 0.0
+variable 	q6ref equal 0.663
+variable 	q8ref equal 0.0
+variable 	q10ref equal 0.363
+variable 	q12ref equal 0.585
+
+variable 	q4 equal c_avql[1]
+variable 	q6 equal c_avql[2]
+variable 	q8 equal c_avql[3]
+variable 	q10 equal c_avql[4]
+variable 	q12 equal c_avql[5]
+
+print "q4     = $(v_q4:%10.6f) delta = $(v_q4-v_q4ref:%10.4f)"
+q4     =   0.000001 delta =     0.0000
+print "q6     = $(v_q6:%10.6f) delta = $(v_q6-v_q6ref:%10.4f)"
+q6     =   0.663325 delta =     0.0003
+print "q8     = $(v_q8:%10.6f) delta = $(v_q8-v_q8ref:%10.4f)"
+q8     =   0.000002 delta =     0.0000
+print "q10    = $(v_q10:%10.6f) delta = $(v_q10-v_q10ref:%10.4f)"
+q10    =   0.362951 delta =    -0.0000
+print "q12    = $(v_q12:%10.6f) delta = $(v_q12-v_q12ref:%10.4f)"
+q12    =   0.585423 delta =     0.0004
+
+# check W_l_hat values
+
+print " "
+ 
+print "Comparison with reference values of W_l_hat"
+Comparison with reference values of W_l_hat
+print "   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, "
+   [Table I in P. Steinhardt, D. Nelson, and M. Ronchetti, 
+print "   Phys. Rev. B 28, 784 (1983).]"
+   Phys. Rev. B 28, 784 (1983).]
+print " "
+ 
+
+variable 	w6hatref equal -0.169754
+variable 	w10hatref equal -0.093967
+
+variable 	w4hat equal c_avwlhat[1]
+variable 	w6hat equal c_avwlhat[2]
+variable 	w8hat equal c_avwlhat[3]
+variable 	w10hat equal c_avwlhat[4]
+variable 	w12hat equal c_avwlhat[5]
+
+print "w6hat  = $(v_w6hat:%10.6f) delta = $(v_w6hat-v_w6hatref:%10.6f)"
+w6hat  =  -0.169754 delta =   0.000000
+print "w10hat = $(v_w10hat:%10.6f) delta = $(v_w10hat-v_w10hatref:%10.6f)"
+w10hat =  -0.093968 delta =  -0.000001
+print " "
+ 
+print "*******************************************************************"
+*******************************************************************
+print " "
+ 
+
+Total wall time: 0:00:00