ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' into hybrid-svector-fix

Browse Source
This commit is contained in:
Axel Kohlmeyer
2019-09-16 22:28:48 -04:00
parent e8f94447f7 357ca4605c
commit f81c1aae26
15 changed files with 1820 additions and 5 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
doc/src/Commands_bond.txt
View File

@ -108,7 +108,7 @@ OPT.
"class2 (ko)"_dihedral_class2.html,
"cosine/shift/exp (o)"_dihedral_cosine_shift_exp.html,
"fourier (io)"_dihedral_fourier.html,
"harmonic (io)"_dihedral_harmonic.html,
"harmonic (iko)"_dihedral_harmonic.html,
"helix (o)"_dihedral_helix.html,
"multi/harmonic (o)"_dihedral_multi_harmonic.html,
"nharmonic (o)"_dihedral_nharmonic.html,

1
doc/src/dihedral_harmonic.txt
View File

@ -8,6 +8,7 @@
dihedral_style harmonic command :h3
dihedral_style harmonic/intel command :h3
dihedral_style harmonic/kk command :h3
dihedral_style harmonic/omp command :h3
[Syntax:]

1
examples/README
View File

@ -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

176
examples/steinhardt/log.13Sept18.bcc.g++.1 Normal file
View File

@ -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

176
examples/steinhardt/log.13Sept18.bcc.g++.4 Normal file
View File

@ -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

172
examples/steinhardt/log.13Sept18.fcc.g++.1 Normal file
View File

@ -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

172
examples/steinhardt/log.13Sept18.fcc.g++.4 Normal file
View File

@ -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

235
examples/steinhardt/log.13Sept18.icos.g++.1 Normal file
View File

@ -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

235
examples/steinhardt/log.13Sept18.icos.g++.4 Normal file
View File

@ -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

2
src/KOKKOS/Install.sh
View File

@ -91,6 +91,8 @@ action dihedral_charmm_kokkos.cpp dihedral_charmm.cpp
action dihedral_charmm_kokkos.h dihedral_charmm.h
action dihedral_class2_kokkos.cpp dihedral_class2.cpp
action dihedral_class2_kokkos.h dihedral_class2.h
action dihedral_harmonic_kokkos.cpp dihedral_harmonic.cpp
action dihedral_harmonic_kokkos.h dihedral_harmonic.h
action dihedral_opls_kokkos.cpp dihedral_opls.cpp
action dihedral_opls_kokkos.h dihedral_opls.h
action domain_kokkos.cpp

537
src/KOKKOS/dihedral_harmonic_kokkos.cpp Normal file
View File

@ -0,0 +1,537 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing author: Stan Moore (SNL)
------------------------------------------------------------------------- */
#include "dihedral_harmonic_kokkos.h"
#include <cmath>
#include <cstdlib>
#include "atom_kokkos.h"
#include "comm.h"
#include "neighbor_kokkos.h"
#include "domain.h"
#include "force.h"
#include "update.h"
#include "memory_kokkos.h"
#include "error.h"
#include "atom_masks.h"
using namespace LAMMPS_NS;
#define TOLERANCE 0.05
#define SMALL 0.001
#define SMALLER 0.00001
/* ---------------------------------------------------------------------- */
template<class DeviceType>
DihedralHarmonicKokkos<DeviceType>::DihedralHarmonicKokkos(LAMMPS *lmp) : DihedralHarmonic(lmp)
{
atomKK = (AtomKokkos *) atom;
neighborKK = (NeighborKokkos *) neighbor;
execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
datamask_read = X_MASK | F_MASK | Q_MASK | ENERGY_MASK | VIRIAL_MASK;
datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
k_warning_flag = DAT::tdual_int_scalar("Dihedral:warning_flag");
d_warning_flag = k_warning_flag.view<DeviceType>();
h_warning_flag = k_warning_flag.h_view;
}
/* ---------------------------------------------------------------------- */
template<class DeviceType>
DihedralHarmonicKokkos<DeviceType>::~DihedralHarmonicKokkos()
{
if (!copymode) {
memoryKK->destroy_kokkos(k_eatom,eatom);
memoryKK->destroy_kokkos(k_vatom,vatom);
}
}
/* ---------------------------------------------------------------------- */
template<class DeviceType>
void DihedralHarmonicKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
{
eflag = eflag_in;
vflag = vflag_in;
ev_init(eflag,vflag,0);
// reallocate per-atom arrays if necessary
if (eflag_atom) {
memoryKK->destroy_kokkos(k_eatom,eatom);
memoryKK->create_kokkos(k_eatom,eatom,maxeatom,"dihedral:eatom");
d_eatom = k_eatom.view<DeviceType>();
}
if (vflag_atom) {
memoryKK->destroy_kokkos(k_vatom,vatom);
memoryKK->create_kokkos(k_vatom,vatom,maxvatom,6,"dihedral:vatom");
d_vatom = k_vatom.view<DeviceType>();
}
k_k.template sync<DeviceType>();
k_cos_shift.template sync<DeviceType>();
k_sin_shift.template sync<DeviceType>();
k_sign.template sync<DeviceType>();
k_multiplicity.template sync<DeviceType>();
x = atomKK->k_x.view<DeviceType>();
f = atomKK->k_f.view<DeviceType>();
neighborKK->k_dihedrallist.template sync<DeviceType>();
dihedrallist = neighborKK->k_dihedrallist.view<DeviceType>();
int ndihedrallist = neighborKK->ndihedrallist;
nlocal = atom->nlocal;
newton_bond = force->newton_bond;
h_warning_flag() = 0;
k_warning_flag.template modify<LMPHostType>();
k_warning_flag.template sync<DeviceType>();
copymode = 1;
// loop over neighbors of my atoms
EV_FLOAT ev;
if (evflag) {
if (newton_bond) {
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagDihedralHarmonicCompute<1,1> >(0,ndihedrallist),*this,ev);
} else {
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagDihedralHarmonicCompute<0,1> >(0,ndihedrallist),*this,ev);
}
} else {
if (newton_bond) {
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagDihedralHarmonicCompute<1,0> >(0,ndihedrallist),*this);
} else {
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagDihedralHarmonicCompute<0,0> >(0,ndihedrallist),*this);
}
}
// error check
k_warning_flag.template modify<DeviceType>();
k_warning_flag.template sync<LMPHostType>();
if (h_warning_flag())
error->warning(FLERR,"Dihedral problem",0);
if (eflag_global) energy += ev.evdwl;
if (vflag_global) {
virial[0] += ev.v[0];
virial[1] += ev.v[1];
virial[2] += ev.v[2];
virial[3] += ev.v[3];
virial[4] += ev.v[4];
virial[5] += ev.v[5];
}
if (eflag_atom) {
k_eatom.template modify<DeviceType>();
k_eatom.template sync<LMPHostType>();
}
if (vflag_atom) {
k_vatom.template modify<DeviceType>();
k_vatom.template sync<LMPHostType>();
}
copymode = 0;
}
template<class DeviceType>
template<int NEWTON_BOND, int EVFLAG>
KOKKOS_INLINE_FUNCTION
void DihedralHarmonicKokkos<DeviceType>::operator()(TagDihedralHarmonicCompute<NEWTON_BOND,EVFLAG>, const int &n, EV_FLOAT& ev) const {
// The f array is atomic
Kokkos::View<F_FLOAT*[3], typename DAT::t_f_array::array_layout,DeviceType,Kokkos::MemoryTraits<Kokkos::Atomic|Kokkos::Unmanaged> > a_f = f;
const int i1 = dihedrallist(n,0);
const int i2 = dihedrallist(n,1);
const int i3 = dihedrallist(n,2);
const int i4 = dihedrallist(n,3);
const int type = dihedrallist(n,4);
// 1st bond
const F_FLOAT vb1x = x(i1,0) - x(i2,0);
const F_FLOAT vb1y = x(i1,1) - x(i2,1);
const F_FLOAT vb1z = x(i1,2) - x(i2,2);
// 2nd bond
const F_FLOAT vb2x = x(i3,0) - x(i2,0);
const F_FLOAT vb2y = x(i3,1) - x(i2,1);
const F_FLOAT vb2z = x(i3,2) - x(i2,2);
const F_FLOAT vb2xm = -vb2x;
const F_FLOAT vb2ym = -vb2y;
const F_FLOAT vb2zm = -vb2z;
// 3rd bond
const F_FLOAT vb3x = x(i4,0) - x(i3,0);
const F_FLOAT vb3y = x(i4,1) - x(i3,1);
const F_FLOAT vb3z = x(i4,2) - x(i3,2);
// c,s calculation
const F_FLOAT ax = vb1y*vb2zm - vb1z*vb2ym;
const F_FLOAT ay = vb1z*vb2xm - vb1x*vb2zm;
const F_FLOAT az = vb1x*vb2ym - vb1y*vb2xm;
const F_FLOAT bx = vb3y*vb2zm - vb3z*vb2ym;
const F_FLOAT by = vb3z*vb2xm - vb3x*vb2zm;
const F_FLOAT bz = vb3x*vb2ym - vb3y*vb2xm;
const F_FLOAT rasq = ax*ax + ay*ay + az*az;
const F_FLOAT rbsq = bx*bx + by*by + bz*bz;
const F_FLOAT rgsq = vb2xm*vb2xm + vb2ym*vb2ym + vb2zm*vb2zm;
const F_FLOAT rg = sqrt(rgsq);
F_FLOAT rginv,ra2inv,rb2inv;
rginv = ra2inv = rb2inv = 0.0;
if (rg > 0) rginv = 1.0/rg;
if (rasq > 0) ra2inv = 1.0/rasq;
if (rbsq > 0) rb2inv = 1.0/rbsq;
const F_FLOAT rabinv = sqrt(ra2inv*rb2inv);
F_FLOAT c = (ax*bx + ay*by + az*bz)*rabinv;
const F_FLOAT s = rg*rabinv*(ax*vb3x + ay*vb3y + az*vb3z);
// error check
if ((c > 1.0 + TOLERANCE || c < (-1.0 - TOLERANCE)) && !d_warning_flag())
Kokkos::atomic_fetch_add(&d_warning_flag(),1);
if (c > 1.0) c = 1.0;
if (c < -1.0) c = -1.0;
const int m = d_multiplicity[type];
F_FLOAT p = 1.0;
F_FLOAT ddf1,df1;
ddf1 = df1 = 0.0;
for (int i = 0; i < m; i++) {
ddf1 = p*c - df1*s;
df1 = p*s + df1*c;
p = ddf1;
}
p = p*d_cos_shift[type] + df1*d_sin_shift[type];
df1 = df1*d_cos_shift[type] - ddf1*d_sin_shift[type];
df1 *= -m;
p += 1.0;
if (m == 0) {
p = 1.0 + d_cos_shift[type];
df1 = 0.0;
}
E_FLOAT edihedral = 0.0;
if (eflag) edihedral = d_k[type] * p;
const F_FLOAT fg = vb1x*vb2xm + vb1y*vb2ym + vb1z*vb2zm;
const F_FLOAT hg = vb3x*vb2xm + vb3y*vb2ym + vb3z*vb2zm;
const F_FLOAT fga = fg*ra2inv*rginv;
const F_FLOAT hgb = hg*rb2inv*rginv;
const F_FLOAT gaa = -ra2inv*rg;
const F_FLOAT gbb = rb2inv*rg;
const F_FLOAT dtfx = gaa*ax;
const F_FLOAT dtfy = gaa*ay;
const F_FLOAT dtfz = gaa*az;
const F_FLOAT dtgx = fga*ax - hgb*bx;
const F_FLOAT dtgy = fga*ay - hgb*by;
const F_FLOAT dtgz = fga*az - hgb*bz;
const F_FLOAT dthx = gbb*bx;
const F_FLOAT dthy = gbb*by;
const F_FLOAT dthz = gbb*bz;
const F_FLOAT df = -d_k[type] * df1;
const F_FLOAT sx2 = df*dtgx;;
const F_FLOAT sy2 = df*dtgy;;
const F_FLOAT sz2 = df*dtgz;;
F_FLOAT f1[3],f2[3],f3[3],f4[3];
f1[0] = df*dtfx;
f1[1] = df*dtfy;
f1[2] = df*dtfz;
f2[0] = sx2 - f1[0];
f2[1] = sy2 - f1[1];
f2[2] = sz2 - f1[2];
f4[0] = df*dthx;
f4[1] = df*dthy;
f4[2] = df*dthz;
f3[0] = -sx2 - f4[0];
f3[1] = -sy2 - f4[1];
f3[2] = -sz2 - f4[2];
// apply force to each of 4 atoms
if (NEWTON_BOND || i1 < nlocal) {
a_f(i1,0) += f1[0];
a_f(i1,1) += f1[1];
a_f(i1,2) += f1[2];
}
if (NEWTON_BOND || i2 < nlocal) {
a_f(i2,0) += f2[0];
a_f(i2,1) += f2[1];
a_f(i2,2) += f2[2];
}
if (NEWTON_BOND || i3 < nlocal) {
a_f(i3,0) += f3[0];
a_f(i3,1) += f3[1];
a_f(i3,2) += f3[2];
}
if (NEWTON_BOND || i4 < nlocal) {
a_f(i4,0) += f4[0];
a_f(i4,1) += f4[1];
a_f(i4,2) += f4[2];
}
if (EVFLAG)
ev_tally(ev,i1,i2,i3,i4,edihedral,f1,f3,f4,
vb1x,vb1y,vb1z,vb2x,vb2y,vb2z,vb3x,vb3y,vb3z);
}
template<class DeviceType>
template<int NEWTON_BOND, int EVFLAG>
KOKKOS_INLINE_FUNCTION
void DihedralHarmonicKokkos<DeviceType>::operator()(TagDihedralHarmonicCompute<NEWTON_BOND,EVFLAG>, const int &n) const {
EV_FLOAT ev;
this->template operator()<NEWTON_BOND,EVFLAG>(TagDihedralHarmonicCompute<NEWTON_BOND,EVFLAG>(), n, ev);
}
/* ---------------------------------------------------------------------- */
template<class DeviceType>
void DihedralHarmonicKokkos<DeviceType>::allocate()
{
DihedralHarmonic::allocate();
int n = atom->ndihedraltypes;
k_k = DAT::tdual_ffloat_1d("DihedralHarmonic::k",n+1);
k_cos_shift = DAT::tdual_ffloat_1d("DihedralHarmonic::cos_shift",n+1);
k_sin_shift = DAT::tdual_ffloat_1d("DihedralHarmonic::sin_shift",n+1);
k_sign = DAT::tdual_int_1d("DihedralHarmonic::sign",n+1);
k_multiplicity = DAT::tdual_int_1d("DihedralHarmonic::multiplicity",n+1);
d_k = k_k.template view<DeviceType>();
d_cos_shift = k_cos_shift.template view<DeviceType>();
d_sin_shift = k_sin_shift.template view<DeviceType>();
d_sign = k_sign.template view<DeviceType>();
d_multiplicity = k_multiplicity.template view<DeviceType>();
}
/* ----------------------------------------------------------------------
set coeffs for one type
------------------------------------------------------------------------- */
template<class DeviceType>
void DihedralHarmonicKokkos<DeviceType>::coeff(int narg, char **arg)
{
DihedralHarmonic::coeff(narg, arg);
int n = atom->ndihedraltypes;
for (int i = 1; i <= n; i++) {
k_k.h_view[i] = k[i];
k_cos_shift.h_view[i] = cos_shift[i];
k_sin_shift.h_view[i] = sin_shift[i];
k_sign.h_view[i] = sign[i];
k_multiplicity.h_view[i] = multiplicity[i];
}
k_k.template modify<LMPHostType>();
k_cos_shift.template modify<LMPHostType>();
k_sin_shift.template modify<LMPHostType>();
k_sign.template modify<LMPHostType>();
k_multiplicity.template modify<LMPHostType>();
}
/* ----------------------------------------------------------------------
proc 0 reads coeffs from restart file, bcasts them
------------------------------------------------------------------------- */
template<class DeviceType>
void DihedralHarmonicKokkos<DeviceType>::read_restart(FILE *fp)
{
DihedralHarmonic::read_restart(fp);
int n = atom->ndihedraltypes;
for (int i = 1; i <= n; i++) {
k_k.h_view[i] = k[i];
k_cos_shift.h_view[i] = cos_shift[i];
k_sin_shift.h_view[i] = sin_shift[i];
k_sign.h_view[i] = sign[i];
k_multiplicity.h_view[i] = multiplicity[i];
}
k_k.template modify<LMPHostType>();
k_cos_shift.template modify<LMPHostType>();
k_sin_shift.template modify<LMPHostType>();
k_sign.template modify<LMPHostType>();
k_multiplicity.template modify<LMPHostType>();
}
/* ----------------------------------------------------------------------
tally energy and virial into global and per-atom accumulators
virial = r1F1 + r2F2 + r3F3 + r4F4 = (r1-r2) F1 + (r3-r2) F3 + (r4-r2) F4
= (r1-r2) F1 + (r3-r2) F3 + (r4-r3 + r3-r2) F4
= vb1*f1 + vb2*f3 + (vb3+vb2)*f4
------------------------------------------------------------------------- */
template<class DeviceType>
//template<int NEWTON_BOND>
KOKKOS_INLINE_FUNCTION
void DihedralHarmonicKokkos<DeviceType>::ev_tally(EV_FLOAT &ev, const int i1, const int i2, const int i3, const int i4,
F_FLOAT &edihedral, F_FLOAT *f1, F_FLOAT *f3, F_FLOAT *f4,
const F_FLOAT &vb1x, const F_FLOAT &vb1y, const F_FLOAT &vb1z,
const F_FLOAT &vb2x, const F_FLOAT &vb2y, const F_FLOAT &vb2z,
const F_FLOAT &vb3x, const F_FLOAT &vb3y, const F_FLOAT &vb3z) const
{
E_FLOAT edihedralquarter;
F_FLOAT v[6];
// The eatom and vatom arrays are atomic
Kokkos::View<E_FLOAT*, typename DAT::t_efloat_1d::array_layout,DeviceType,Kokkos::MemoryTraits<Kokkos::Atomic|Kokkos::Unmanaged> > v_eatom = k_eatom.view<DeviceType>();
Kokkos::View<F_FLOAT*[6], typename DAT::t_virial_array::array_layout,DeviceType,Kokkos::MemoryTraits<Kokkos::Atomic|Kokkos::Unmanaged> > v_vatom = k_vatom.view<DeviceType>();
if (eflag_either) {
if (eflag_global) {
if (newton_bond) ev.evdwl += edihedral;
else {
edihedralquarter = 0.25*edihedral;
if (i1 < nlocal) ev.evdwl += edihedralquarter;
if (i2 < nlocal) ev.evdwl += edihedralquarter;
if (i3 < nlocal) ev.evdwl += edihedralquarter;
if (i4 < nlocal) ev.evdwl += edihedralquarter;
}
}
if (eflag_atom) {
edihedralquarter = 0.25*edihedral;
if (newton_bond || i1 < nlocal) v_eatom[i1] += edihedralquarter;
if (newton_bond || i2 < nlocal) v_eatom[i2] += edihedralquarter;
if (newton_bond || i3 < nlocal) v_eatom[i3] += edihedralquarter;
if (newton_bond || i4 < nlocal) v_eatom[i4] += edihedralquarter;
}
}
if (vflag_either) {
v[0] = vb1x*f1[0] + vb2x*f3[0] + (vb3x+vb2x)*f4[0];
v[1] = vb1y*f1[1] + vb2y*f3[1] + (vb3y+vb2y)*f4[1];
v[2] = vb1z*f1[2] + vb2z*f3[2] + (vb3z+vb2z)*f4[2];
v[3] = vb1x*f1[1] + vb2x*f3[1] + (vb3x+vb2x)*f4[1];
v[4] = vb1x*f1[2] + vb2x*f3[2] + (vb3x+vb2x)*f4[2];
v[5] = vb1y*f1[2] + vb2y*f3[2] + (vb3y+vb2y)*f4[2];
if (vflag_global) {
if (newton_bond) {
ev.v[0] += v[0];
ev.v[1] += v[1];
ev.v[2] += v[2];
ev.v[3] += v[3];
ev.v[4] += v[4];
ev.v[5] += v[5];
} else {
if (i1 < nlocal) {
ev.v[0] += 0.25*v[0];
ev.v[1] += 0.25*v[1];
ev.v[2] += 0.25*v[2];
ev.v[3] += 0.25*v[3];
ev.v[4] += 0.25*v[4];
ev.v[5] += 0.25*v[5];
}
if (i2 < nlocal) {
ev.v[0] += 0.25*v[0];
ev.v[1] += 0.25*v[1];
ev.v[2] += 0.25*v[2];
ev.v[3] += 0.25*v[3];
ev.v[4] += 0.25*v[4];
ev.v[5] += 0.25*v[5];
}
if (i3 < nlocal) {
ev.v[0] += 0.25*v[0];
ev.v[1] += 0.25*v[1];
ev.v[2] += 0.25*v[2];
ev.v[3] += 0.25*v[3];
ev.v[4] += 0.25*v[4];
ev.v[5] += 0.25*v[5];
}
if (i4 < nlocal) {
ev.v[0] += 0.25*v[0];
ev.v[1] += 0.25*v[1];
ev.v[2] += 0.25*v[2];
ev.v[3] += 0.25*v[3];
ev.v[4] += 0.25*v[4];
ev.v[5] += 0.25*v[5];
}
}
}
if (vflag_atom) {
if (newton_bond || i1 < nlocal) {
v_vatom(i1,0) += 0.25*v[0];
v_vatom(i1,1) += 0.25*v[1];
v_vatom(i1,2) += 0.25*v[2];
v_vatom(i1,3) += 0.25*v[3];
v_vatom(i1,4) += 0.25*v[4];
v_vatom(i1,5) += 0.25*v[5];
}
if (newton_bond || i2 < nlocal) {
v_vatom(i2,0) += 0.25*v[0];
v_vatom(i2,1) += 0.25*v[1];
v_vatom(i2,2) += 0.25*v[2];
v_vatom(i2,3) += 0.25*v[3];
v_vatom(i2,4) += 0.25*v[4];
v_vatom(i2,5) += 0.25*v[5];
}
if (newton_bond || i3 < nlocal) {
v_vatom(i3,0) += 0.25*v[0];
v_vatom(i3,1) += 0.25*v[1];
v_vatom(i3,2) += 0.25*v[2];
v_vatom(i3,3) += 0.25*v[3];
v_vatom(i3,4) += 0.25*v[4];
v_vatom(i3,5) += 0.25*v[5];
}
if (newton_bond || i4 < nlocal) {
v_vatom(i4,0) += 0.25*v[0];
v_vatom(i4,1) += 0.25*v[1];
v_vatom(i4,2) += 0.25*v[2];
v_vatom(i4,3) += 0.25*v[3];
v_vatom(i4,4) += 0.25*v[4];
v_vatom(i4,5) += 0.25*v[5];
}
}
}
}
/* ---------------------------------------------------------------------- */
namespace LAMMPS_NS {
template class DihedralHarmonicKokkos<LMPDeviceType>;
#ifdef KOKKOS_ENABLE_CUDA
template class DihedralHarmonicKokkos<LMPHostType>;
#endif
}

104
src/KOKKOS/dihedral_harmonic_kokkos.h Normal file
View File

@ -0,0 +1,104 @@
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
#ifdef DIHEDRAL_CLASS
DihedralStyle(harmonic/kk,DihedralHarmonicKokkos<LMPDeviceType>)
DihedralStyle(harmonic/kk/device,DihedralHarmonicKokkos<LMPDeviceType>)
DihedralStyle(harmonic/kk/host,DihedralHarmonicKokkos<LMPHostType>)
#else
#ifndef LMP_DIHEDRAL_HARMONIC_KOKKOS_H
#define LMP_DIHEDRAL_HARMONIC_KOKKOS_H
#include "dihedral_harmonic.h"
#include "kokkos_type.h"
namespace LAMMPS_NS {
template<int NEWTON_BOND, int EVFLAG>
struct TagDihedralHarmonicCompute{};
template<class DeviceType>
class DihedralHarmonicKokkos : public DihedralHarmonic {
public:
typedef DeviceType device_type;
typedef EV_FLOAT value_type;
typedef ArrayTypes<DeviceType> AT;
DihedralHarmonicKokkos(class LAMMPS *);
virtual ~DihedralHarmonicKokkos();
void compute(int, int);
void coeff(int, char **);
void read_restart(FILE *);
template<int NEWTON_BOND, int EVFLAG>
KOKKOS_INLINE_FUNCTION
void operator()(TagDihedralHarmonicCompute<NEWTON_BOND,EVFLAG>, const int&, EV_FLOAT&) const;
template<int NEWTON_BOND, int EVFLAG>
KOKKOS_INLINE_FUNCTION
void operator()(TagDihedralHarmonicCompute<NEWTON_BOND,EVFLAG>, const int&) const;
//template<int NEWTON_BOND>
KOKKOS_INLINE_FUNCTION
void ev_tally(EV_FLOAT &ev, const int i1, const int i2, const int i3, const int i4,
F_FLOAT &edihedral, F_FLOAT *f1, F_FLOAT *f3, F_FLOAT *f4,
const F_FLOAT &vb1x, const F_FLOAT &vb1y, const F_FLOAT &vb1z,
const F_FLOAT &vb2x, const F_FLOAT &vb2y, const F_FLOAT &vb2z,
const F_FLOAT &vb3x, const F_FLOAT &vb3y, const F_FLOAT &vb3z) const;
protected:
class NeighborKokkos *neighborKK;
typename AT::t_x_array_randomread x;
typename AT::t_f_array f;
typename AT::t_int_2d dihedrallist;
DAT::tdual_efloat_1d k_eatom;
DAT::tdual_virial_array k_vatom;
typename ArrayTypes<DeviceType>::t_efloat_1d d_eatom;
typename ArrayTypes<DeviceType>::t_virial_array d_vatom;
int nlocal,newton_bond;
int eflag,vflag;
DAT::tdual_int_scalar k_warning_flag;
typename AT::t_int_scalar d_warning_flag;
HAT::t_int_scalar h_warning_flag;
DAT::tdual_ffloat_1d k_k;
DAT::tdual_ffloat_1d k_cos_shift;
DAT::tdual_ffloat_1d k_sin_shift;
DAT::tdual_int_1d k_sign;
DAT::tdual_int_1d k_multiplicity;
typename AT::t_ffloat_1d d_k;
typename AT::t_ffloat_1d d_cos_shift;
typename AT::t_ffloat_1d d_sin_shift;
typename AT::t_int_1d d_sign;
typename AT::t_int_1d d_multiplicity;
void allocate();
};
}
#endif
#endif
/* ERROR/WARNING messages:
*/

2
src/MOLECULE/dihedral_harmonic.cpp
View File

@ -42,7 +42,7 @@ DihedralHarmonic::DihedralHarmonic(LAMMPS *lmp) : Dihedral(lmp)
DihedralHarmonic::~DihedralHarmonic()
{
if (allocated) {
if (allocated && !copymode) {
memory->destroy(setflag);
memory->destroy(k);
memory->destroy(sign);

6
src/MOLECULE/dihedral_harmonic.h
View File

@ -29,16 +29,16 @@ class DihedralHarmonic : public Dihedral {
DihedralHarmonic(class LAMMPS *);
virtual ~DihedralHarmonic();
virtual void compute(int, int);
void coeff(int, char **);
virtual void coeff(int, char **);
void write_restart(FILE *);
void read_restart(FILE *);
virtual void read_restart(FILE *);
void write_data(FILE *);
protected:
double *k,*cos_shift,*sin_shift;
int *sign,*multiplicity;
void allocate();
virtual void allocate();
};
}

4
src/compute_pe.cpp
View File

@ -26,6 +26,7 @@
#include "modify.h"
#include "domain.h"
#include "error.h"
#include "atom_masks.h"
using namespace LAMMPS_NS;
@ -65,6 +66,9 @@ ComputePE::ComputePE(LAMMPS *lmp, int narg, char **arg) :
iarg++;
}
}
datamask_read = EMPTY_MASK;
datamask_modify = EMPTY_MASK;
}
/* ---------------------------------------------------------------------- */