added examples/threebody, fix reaxc/speceies/kk

This commit is contained in:
Steve Plimpton
2016-11-04 10:56:04 -06:00
parent d9891abdf4
commit f48b71f46b
28 changed files with 1180 additions and 337 deletions

View File

@ -69,8 +69,8 @@ velocity for each atom. Note that if there is only one atom in the
bin, its thermal velocity will thus be 0.0. bin, its thermal velocity will thus be 0.0.
After the spatially-averaged velocity field has been subtracted from After the spatially-averaged velocity field has been subtracted from
each atom, the temperature is calculated by the formula KE = (dim/2 N each atom, the temperature is calculated by the formula KE = (dim*N
- dim*Nx*Ny*Nz) k T, where KE = total kinetic energy of the group of - dim*Nx*Ny*Nz) k T/2, where KE = total kinetic energy of the group of
atoms (sum of 1/2 m v^2), dim = 2 or 3 = dimensionality of the atoms (sum of 1/2 m v^2), dim = 2 or 3 = dimensionality of the
simulation, N = number of atoms in the group, k = Boltzmann constant, simulation, N = number of atoms in the group, k = Boltzmann constant,
and T = temperature. The dim*Nx*Ny*Nz term are degrees of freedom and T = temperature. The dim*Nx*Ny*Nz term are degrees of freedom

View File

@ -8,6 +8,7 @@
fix reax/bonds command :h3 fix reax/bonds command :h3
fix reax/c/bonds command :h3 fix reax/c/bonds command :h3
fix reax/c/bonds/kk command :h3
[Syntax:] [Syntax:]
@ -47,6 +48,31 @@ commands"_Section_howto.html#howto_15. No parameter of this fix can
be used with the {start/stop} keywords of the "run"_run.html command. be used with the {start/stop} keywords of the "run"_run.html command.
This fix is not invoked during "energy minimization"_minimize.html. This fix is not invoked during "energy minimization"_minimize.html.
:line
Styles with a {gpu}, {intel}, {kk}, {omp}, or {opt} suffix are
functionally the same as the corresponding style without the suffix.
They have been optimized to run faster, depending on your available
hardware, as discussed in "Section_accelerate"_Section_accelerate.html
of the manual. The accelerated styles take the same arguments and
should produce the same results, except for round-off and precision
issues.
These accelerated styles are part of the GPU, USER-INTEL, KOKKOS,
USER-OMP and OPT packages, respectively. They are only enabled if
LAMMPS was built with those packages. See the "Making
LAMMPS"_Section_start.html#start_3 section for more info.
You can specify the accelerated styles explicitly in your input script
by including their suffix, or you can use the "-suffix command-line
switch"_Section_start.html#start_7 when you invoke LAMMPS, or you can
use the "suffix"_suffix.html command in your input script.
See "Section_accelerate"_Section_accelerate.html of the manual for
more instructions on how to use the accelerated styles effectively.
:line
[Restrictions:] [Restrictions:]
The fix reax/bonds command requires that the "pair_style The fix reax/bonds command requires that the "pair_style

View File

@ -7,6 +7,7 @@
:line :line
fix reax/c/species command :h3 fix reax/c/species command :h3
fix reax/c/species/kk command :h3
[Syntax:] [Syntax:]
@ -129,6 +130,31 @@ No parameter of this fix can be used with the {start/stop} keywords of
the "run"_run.html command. This fix is not invoked during "energy the "run"_run.html command. This fix is not invoked during "energy
minimization"_minimize.html. minimization"_minimize.html.
:line
Styles with a {gpu}, {intel}, {kk}, {omp}, or {opt} suffix are
functionally the same as the corresponding style without the suffix.
They have been optimized to run faster, depending on your available
hardware, as discussed in "Section_accelerate"_Section_accelerate.html
of the manual. The accelerated styles take the same arguments and
should produce the same results, except for round-off and precision
issues.
These accelerated styles are part of the GPU, USER-INTEL, KOKKOS,
USER-OMP and OPT packages, respectively. They are only enabled if
LAMMPS was built with those packages. See the "Making
LAMMPS"_Section_start.html#start_3 section for more info.
You can specify the accelerated styles explicitly in your input script
by including their suffix, or you can use the "-suffix command-line
switch"_Section_start.html#start_7 when you invoke LAMMPS, or you can
use the "suffix"_suffix.html command in your input script.
See "Section_accelerate"_Section_accelerate.html of the manual for
more instructions on how to use the accelerated styles effectively.
:line
[Restrictions:] [Restrictions:]
The fix species currently only works with The fix species currently only works with

View File

@ -94,6 +94,13 @@ Of course this is also possible by not using any suffix commands, and
explictly appending or not appending the suffix to the relevant explictly appending or not appending the suffix to the relevant
commands in your input script. commands in your input script.
NOTE: The default "run_style"_run_style.html verlet is invoked prior to
reading the input script and is therefore not affected by a suffix command
in the input script. The KOKKOS package requires "run_style verlet/kk",
so when using the KOKKOS package it is necessary to either use the command
line "-sf kk" command or add an explicit "run_style verlet" command to the
input script.
[Restrictions:] none [Restrictions:] none
[Related commands:] [Related commands:]

View File

@ -0,0 +1,74 @@
# DATE: 2013-03-21 CONTRIBUTOR: Cem Sevik CITATION: Kinaci, Haskins, Sevik and Cagin, Phys Rev B, 86, 115410 (2012)
# Tersoff parameters for B, C, and BN-C hybrid based graphene like nano structures
# multiple entries can be added to this file, LAMMPS reads the ones it needs
# these entries are in LAMMPS "metal" units:
# A,B = eV; lambda1,lambda2,lambda3 = 1/Angstroms; R,D = Angstroms
# other quantities are unitless
# Cem Sevik (csevik at anadolu.edu.tr) takes full blame for this
# file. It specifies B-N, B-C, and N-C interaction parameters
# generated and published by the reseacrh group of Prof. Tahir Cagin.
# 1. Physical Review B 84, 085409 2011
# Characterization of thermal transport in low-dimensional boron nitride nanostructures,
#
# 2. Physical Review B 86, 075403 2012
# Influence of disorder on thermal transport properties of boron nitride nanostructures
#
# 3. Physical Review B 86, 075403 2012, Please see for further information about B-C and N-C parameters
# Thermal conductivity of BN-C nanostructures
#
# The file also specifies C-C, interaction parameters
# generated and published by the reseacrh group of Dr. D. A. Broido
# Physical Review B 81, 205441 2010
# Optimized Tersoff and Brenner empirical potential parameters for
# lattice dynamics and phonon thermal transport in carbon nanotubes and graphene
# Users in referring the full parameters can cite the full parameter paper (3) as:
# A. Kinaci, J. B. Haskins, C. Sevik, T. Cagin, Physical Review B 86, 115410 (2012)
# Thermal conductivity of BN-C nanostructures
#
# format of a single entry (one or more lines):
# element 1, element 2, element 3,
# m, gamma, lambda3, c, d, costheta0, n, beta, lambda2, B, R, D, lambda1, A
N B B 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
N B N 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
N B C 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
B N B 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
B N N 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
B N C 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.199 340.00 1.95 0.05 3.568 1380.0
N N B 3.0 1.0 0.0 17.7959 5.9484 0.00000 0.6184432 0.019251 2.6272721 138.77866 2.0 0.1 2.8293093 128.86866
N N N 3.0 1.0 0.0 17.7959 5.9484 0.00000 0.6184432 0.019251 2.6272721 138.77866 2.0 0.1 2.8293093 128.86866
N N C 3.0 1.0 0.0 17.7959 5.9484 0.00000 0.6184432 0.019251 2.6272721 138.77866 2.0 0.1 2.8293093 128.86866
B B B 3.0 1.0 0.0 0.52629 0.001587 0.5 3.9929061 1.6e-6 2.0774982 43.132016 2.0 0.1 2.2372578 40.0520156
B B N 3.0 1.0 0.0 0.52629 0.001587 0.5 3.9929061 1.6e-6 2.0774982 43.132016 2.0 0.1 2.2372578 40.0520156
B B C 3.0 1.0 0.0 0.52629 0.001587 0.5 3.9929061 1.6e-6 2.0774982 43.132016 2.0 0.1 2.2372578 40.0520156
C C C 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2119 430.00 1.95 0.15 3.4879 1393.6
C C B 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2119 430.00 1.95 0.15 3.4879 1393.6
C C N 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2119 430.00 1.95 0.15 3.4879 1393.6
C B B 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
C B N 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
C B C 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
C N B 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78
C N N 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78
C N C 3.0 1.0 0.0 3.8049e4 4.3484 -0.93000 0.72751 1.5724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78
B C C 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
B C B 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
B C N 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 339.068910 1.95 0.10 3.5279 1386.78
N C C 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78
N C B 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78
N C N 3.0 1.0 0.0 25000 4.3484 -0.89000 0.72751 1.25724e-7 2.2054 387.575152 1.95 0.10 3.5279 1386.78

View File

@ -0,0 +1,233 @@
### DATE: 2013-08-09 CONTRIBUTOR: X. W. Zhou, xzhou@sandia.gov, CITATION: Zhou, Ward, Martin, van Swol, Cruz-Campa, and D. Zubia, Phys. Rev. B, 88, 085309 (2013).
#
# Note that the way the parameters can be entered is not unique.
# As one way, we assume that eps_ijk is equal to eps_ik and
# lambda_ijk is equal to sqrt(lambda_ij*eps_ij*lambda_ik*eps_ik)/eps_ik,
# and all other parameters in the ijk line are for ik.
#
# The twobody ik pair parameters are entered on the i*k lines, where *
# can be any species. This is consistent with the LAMMPS requirement
# that twobody ik parameters be defined on the ikk line. Entries on all
# the other i*k lines are ignored by LAMMPS
#
# These entries are in LAMMPS "metal" units: epsilon = eV;
# sigma = Angstroms; other quantities are unitless
#
# cutoff distance = 4.632
# eps sigma a lambda gamma cos(theta) A B p q tol
Cd Cd Cd 1.182358e+00 2.663951e+00 1.527956e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Cd Te 1.385284e+00 2.352141e+00 1.810919e+00 3.002537e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Cd Zn 6.908179e-01 2.238699e+00 1.812616e+00 4.251831e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Cd Se 1.352371e+00 2.045165e+00 1.953387e+00 3.038855e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Cd Hg 4.881231e-01 2.432694e+00 1.677987e+00 5.058167e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Cd S 1.300376e+00 1.804151e+00 2.124568e+00 3.099013e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te Cd 1.182358e+00 2.663951e+00 1.527956e+00 3.517858e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te Te 1.385284e+00 2.352141e+00 1.810919e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te Zn 6.908179e-01 2.238699e+00 1.812616e+00 4.602259e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te Se 1.352371e+00 2.045165e+00 1.953387e+00 3.289311e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te Hg 4.881231e-01 2.432694e+00 1.677987e+00 5.475051e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Te S 1.300376e+00 1.804151e+00 2.124568e+00 3.354428e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn Cd 1.182358e+00 2.663951e+00 1.527956e+00 2.484224e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn Te 1.385284e+00 2.352141e+00 1.810919e+00 2.295069e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn Zn 6.908179e-01 2.238699e+00 1.812616e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn Se 1.352371e+00 2.045165e+00 1.953387e+00 2.322829e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn Hg 4.881231e-01 2.432694e+00 1.677987e+00 3.866344e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Zn S 1.300376e+00 1.804151e+00 2.124568e+00 2.368813e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se Cd 1.182358e+00 2.663951e+00 1.527956e+00 3.475816e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se Te 1.385284e+00 2.352141e+00 1.810919e+00 3.211159e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se Zn 6.908179e-01 2.238699e+00 1.812616e+00 4.547256e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se Se 1.352371e+00 2.045165e+00 1.953387e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se Hg 4.881231e-01 2.432694e+00 1.677987e+00 5.409618e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Se S 1.300376e+00 1.804151e+00 2.124568e+00 3.314338e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg Cd 1.182358e+00 2.663951e+00 1.527956e+00 2.088207e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg Te 1.385284e+00 2.352141e+00 1.810919e+00 1.929206e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg Zn 6.908179e-01 2.238699e+00 1.812616e+00 2.731909e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg Se 1.352371e+00 2.045165e+00 1.953387e+00 1.952541e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg Hg 4.881231e-01 2.432694e+00 1.677987e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd Hg S 1.300376e+00 1.804151e+00 2.124568e+00 1.991194e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd S Cd 1.182358e+00 2.663951e+00 1.527956e+00 3.408343e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.674460e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd S Te 1.385284e+00 2.352141e+00 1.810919e+00 3.148823e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd S Zn 6.908179e-01 2.238699e+00 1.812616e+00 4.458985e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd S Se 1.352371e+00 2.045165e+00 1.953387e+00 3.186911e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Cd S Hg 4.881231e-01 2.432694e+00 1.677987e+00 5.304605e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Cd S S 1.300376e+00 1.804151e+00 2.124568e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd Te 1.849775e+00 2.905254e+00 1.594353e+00 2.812506e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd Zn 1.546239e+00 2.056363e+00 1.907922e+00 3.076200e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd Se 1.295053e+00 2.231716e+00 1.809645e+00 3.361313e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd Hg 1.204715e+00 2.135591e+00 1.892491e+00 3.485063e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Cd S 1.450015e+00 2.297301e+00 1.726905e+00 3.176630e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Te Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.755548e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Te Te 1.849775e+00 2.905254e+00 1.594353e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Te Zn 1.546239e+00 2.056363e+00 1.907922e+00 3.554713e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Te Se 1.295053e+00 2.231716e+00 1.809645e+00 3.884177e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Te Hg 1.204715e+00 2.135591e+00 1.892491e+00 4.027176e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Te S 1.450015e+00 2.297301e+00 1.726905e+00 3.670765e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.433620e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn Te 1.849775e+00 2.905254e+00 1.594353e+00 2.971408e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn Zn 1.546239e+00 2.056363e+00 1.907922e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn Se 1.295053e+00 2.231716e+00 1.809645e+00 3.551222e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn Hg 1.204715e+00 2.135591e+00 1.892491e+00 3.681964e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Zn S 1.450015e+00 2.297301e+00 1.726905e+00 3.356105e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Se Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.142373e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Se Te 1.849775e+00 2.905254e+00 1.594353e+00 2.719366e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Se Zn 1.546239e+00 2.056363e+00 1.907922e+00 2.974328e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Se Se 1.295053e+00 2.231716e+00 1.809645e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Se Hg 1.204715e+00 2.135591e+00 1.892491e+00 3.369652e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Se S 1.450015e+00 2.297301e+00 1.726905e+00 3.071433e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.030791e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg Te 1.849775e+00 2.905254e+00 1.594353e+00 2.622805e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg Zn 1.546239e+00 2.056363e+00 1.907922e+00 2.868714e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg Se 1.295053e+00 2.231716e+00 1.809645e+00 3.134597e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg Hg 1.204715e+00 2.135591e+00 1.892491e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te Hg S 1.450015e+00 2.297301e+00 1.726905e+00 2.962370e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te S Cd 1.385284e+00 2.352141e+00 1.810919e+00 3.325065e+01 1.200000e+00 -3.333333e-01 7.049600e+00 8.861252e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te S Te 1.849775e+00 2.905254e+00 1.594353e+00 2.877465e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.307283e-01 4.000000e+00 0.000000e+00 0.000000e+00
Te S Zn 1.546239e+00 2.056363e+00 1.907922e+00 3.147250e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te S Se 1.295053e+00 2.231716e+00 1.809645e+00 3.438949e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te S Hg 1.204715e+00 2.135591e+00 1.892491e+00 3.565557e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Te S S 1.450015e+00 2.297301e+00 1.726905e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd Cd 6.908179e-01 2.238699e+00 1.812616e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd Te 1.546239e+00 2.056363e+00 1.907922e+00 2.172335e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd Zn 1.392961e+00 2.367650e+00 1.525521e+00 2.288736e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd Se 1.691181e+00 2.028827e+00 1.836907e+00 2.077161e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd Hg 4.951616e-01 2.239186e+00 1.761363e+00 3.838766e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Cd S 2.208390e+00 2.323783e+00 1.589241e+00 1.817721e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te Cd 6.908179e-01 2.238699e+00 1.812616e+00 4.862279e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te Te 1.546239e+00 2.056363e+00 1.907922e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te Zn 1.392961e+00 2.367650e+00 1.525521e+00 3.424146e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te Se 1.691181e+00 2.028827e+00 1.836907e+00 3.107611e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te Hg 4.951616e-01 2.239186e+00 1.761363e+00 5.743124e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Te S 2.208390e+00 2.323783e+00 1.589241e+00 2.719467e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn Cd 6.908179e-01 2.238699e+00 1.812616e+00 4.614993e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn Te 1.546239e+00 2.056363e+00 1.907922e+00 3.084711e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn Zn 1.392961e+00 2.367650e+00 1.525521e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn Se 1.691181e+00 2.028827e+00 1.836907e+00 2.949563e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn Hg 4.951616e-01 2.239186e+00 1.761363e+00 5.451040e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Zn S 2.208390e+00 2.323783e+00 1.589241e+00 2.581160e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se Cd 6.908179e-01 2.238699e+00 1.812616e+00 5.085067e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se Te 1.546239e+00 2.056363e+00 1.907922e+00 3.398914e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se Zn 1.392961e+00 2.367650e+00 1.525521e+00 3.581039e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se Se 1.691181e+00 2.028827e+00 1.836907e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se Hg 4.951616e-01 2.239186e+00 1.761363e+00 6.006272e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Se S 2.208390e+00 2.323783e+00 1.589241e+00 2.844072e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg Cd 6.908179e-01 2.238699e+00 1.812616e+00 2.751535e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg Te 1.546239e+00 2.056363e+00 1.907922e+00 1.839156e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg Zn 1.392961e+00 2.367650e+00 1.525521e+00 1.937704e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg Se 1.691181e+00 2.028827e+00 1.836907e+00 1.758578e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg Hg 4.951616e-01 2.239186e+00 1.761363e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn Hg S 2.208390e+00 2.323783e+00 1.589241e+00 1.538930e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn S Cd 6.908179e-01 2.238699e+00 1.812616e+00 5.810847e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.010632e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn S Te 1.546239e+00 2.056363e+00 1.907922e+00 3.884033e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.255846e+00 4.000000e+00 0.000000e+00 0.000000e+00
Zn S Zn 1.392961e+00 2.367650e+00 1.525521e+00 4.092153e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.676279e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn S Se 1.691181e+00 2.028827e+00 1.836907e+00 3.713865e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn S Hg 4.951616e-01 2.239186e+00 1.761363e+00 6.863534e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Zn S S 2.208390e+00 2.323783e+00 1.589241e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd Cd 1.352371e+00 2.045165e+00 1.953387e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd Te 1.295053e+00 2.231716e+00 1.809645e+00 3.321142e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd Zn 1.691181e+00 2.028827e+00 1.836907e+00 2.906271e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd Se 2.400781e+00 2.789002e+00 1.544925e+00 2.439242e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd Hg 1.299758e+00 2.113406e+00 1.831821e+00 3.315126e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Cd S 1.307592e+00 2.229392e+00 1.747782e+00 3.305180e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Te Cd 1.352371e+00 2.045165e+00 1.953387e+00 3.180382e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Te Te 1.295053e+00 2.231716e+00 1.809645e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Te Zn 1.691181e+00 2.028827e+00 1.836907e+00 2.844016e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Te Se 2.400781e+00 2.789002e+00 1.544925e+00 2.386992e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Te Hg 1.299758e+00 2.113406e+00 1.831821e+00 3.244113e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Te S 1.307592e+00 2.229392e+00 1.747782e+00 3.234380e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn Cd 1.352371e+00 2.045165e+00 1.953387e+00 3.634382e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn Te 1.295053e+00 2.231716e+00 1.809645e+00 3.713938e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn Zn 1.691181e+00 2.028827e+00 1.836907e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn Se 2.400781e+00 2.789002e+00 1.544925e+00 2.727735e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn Hg 1.299758e+00 2.113406e+00 1.831821e+00 3.707211e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Zn S 1.307592e+00 2.229392e+00 1.747782e+00 3.696088e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Se Cd 1.352371e+00 2.045165e+00 1.953387e+00 4.330238e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Se Te 1.295053e+00 2.231716e+00 1.809645e+00 4.425026e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Se Zn 1.691181e+00 2.028827e+00 1.836907e+00 3.872260e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Se Se 2.400781e+00 2.789002e+00 1.544925e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Se Hg 1.299758e+00 2.113406e+00 1.831821e+00 4.417011e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Se S 1.307592e+00 2.229392e+00 1.747782e+00 4.403758e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg Cd 1.352371e+00 2.045165e+00 1.953387e+00 3.186153e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg Te 1.295053e+00 2.231716e+00 1.809645e+00 3.255898e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg Zn 1.691181e+00 2.028827e+00 1.836907e+00 2.849177e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg Se 2.400781e+00 2.789002e+00 1.544925e+00 2.391323e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg Hg 1.299758e+00 2.113406e+00 1.831821e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se Hg S 1.307592e+00 2.229392e+00 1.747782e+00 3.240249e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se S Cd 1.352371e+00 2.045165e+00 1.953387e+00 3.195742e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.116149e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se S Te 1.295053e+00 2.231716e+00 1.809645e+00 3.265696e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.005396e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se S Zn 1.691181e+00 2.028827e+00 1.836907e+00 2.857751e+01 1.200000e+00 -3.333333e-01 7.049600e+00 9.510930e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se S Se 2.400781e+00 2.789002e+00 1.544925e+00 2.398520e+01 1.200000e+00 -3.333333e-01 7.917000e+00 7.672131e-01 4.000000e+00 0.000000e+00 0.000000e+00
Se S Hg 1.299758e+00 2.113406e+00 1.831821e+00 3.259780e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Se S S 1.307592e+00 2.229392e+00 1.747782e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd Cd 4.881231e-01 2.432694e+00 1.677987e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd Te 1.204715e+00 2.135591e+00 1.892491e+00 2.068740e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd Zn 4.951616e-01 2.239186e+00 1.761363e+00 3.226819e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd Se 1.299758e+00 2.113406e+00 1.831821e+00 1.991668e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd Hg 1.272807e+00 2.699097e+00 1.498503e+00 2.012643e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Cd S 1.531211e+00 2.025045e+00 1.833708e+00 1.834976e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te Cd 4.881231e-01 2.432694e+00 1.677987e+00 5.105765e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te Te 1.204715e+00 2.135591e+00 1.892491e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te Zn 4.951616e-01 2.239186e+00 1.761363e+00 5.069347e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te Se 1.299758e+00 2.113406e+00 1.831821e+00 3.128919e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te Hg 1.272807e+00 2.699097e+00 1.498503e+00 3.161872e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Te S 1.531211e+00 2.025045e+00 1.833708e+00 2.882756e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn Cd 4.881231e-01 2.432694e+00 1.677987e+00 3.273348e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn Te 1.204715e+00 2.135591e+00 1.892491e+00 2.083602e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn Zn 4.951616e-01 2.239186e+00 1.761363e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn Se 1.299758e+00 2.113406e+00 1.831821e+00 2.005976e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn Hg 1.272807e+00 2.699097e+00 1.498503e+00 2.027102e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Zn S 1.531211e+00 2.025045e+00 1.833708e+00 1.848159e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se Cd 4.881231e-01 2.432694e+00 1.677987e+00 5.303345e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se Te 1.204715e+00 2.135591e+00 1.892491e+00 3.375766e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se Zn 4.951616e-01 2.239186e+00 1.761363e+00 5.265518e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se Se 1.299758e+00 2.113406e+00 1.831821e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se Hg 1.272807e+00 2.699097e+00 1.498503e+00 3.284228e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Se S 1.531211e+00 2.025045e+00 1.833708e+00 2.994311e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg Cd 4.881231e-01 2.432694e+00 1.677987e+00 5.248074e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg Te 1.204715e+00 2.135591e+00 1.892491e+00 3.340584e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg Zn 4.951616e-01 2.239186e+00 1.761363e+00 5.210641e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg Se 1.299758e+00 2.113406e+00 1.831821e+00 3.216129e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg Hg 1.272807e+00 2.699097e+00 1.498503e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg Hg S 1.531211e+00 2.025045e+00 1.833708e+00 2.963105e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg S Cd 4.881231e-01 2.432694e+00 1.677987e+00 5.756205e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.250999e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg S Te 1.204715e+00 2.135591e+00 1.892491e+00 3.664028e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.445180e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg S Zn 4.951616e-01 2.239186e+00 1.761363e+00 5.715148e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.461167e-01 4.000000e+00 0.000000e+00 0.000000e+00
Hg S Se 1.299758e+00 2.113406e+00 1.831821e+00 3.527522e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.150200e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg S Hg 1.272807e+00 2.699097e+00 1.498503e+00 3.564673e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.211532e+00 4.000000e+00 0.000000e+00 0.000000e+00
Hg S S 1.531211e+00 2.025045e+00 1.833708e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Cd Cd 1.300376e+00 1.804151e+00 2.124568e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Cd Te 1.450015e+00 2.297301e+00 1.726905e+00 3.077737e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Cd Zn 2.208390e+00 2.323783e+00 1.589241e+00 2.493905e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Cd Se 1.307592e+00 2.229392e+00 1.747782e+00 3.241019e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Cd Hg 1.531211e+00 2.025045e+00 1.833708e+00 2.995023e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Cd S 2.434871e+00 2.423171e+00 1.711097e+00 2.375088e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Te Cd 1.300376e+00 1.804151e+00 2.124568e+00 3.431904e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Te Te 1.450015e+00 2.297301e+00 1.726905e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Te Zn 2.208390e+00 2.323783e+00 1.589241e+00 2.633490e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Te Se 1.307592e+00 2.229392e+00 1.747782e+00 3.422421e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Te Hg 1.531211e+00 2.025045e+00 1.833708e+00 3.162656e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Te S 2.434871e+00 2.423171e+00 1.711097e+00 2.508023e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Zn Cd 1.300376e+00 1.804151e+00 2.124568e+00 4.235326e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Zn Te 1.450015e+00 2.297301e+00 1.726905e+00 4.010837e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Zn Zn 2.208390e+00 2.323783e+00 1.589241e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Zn Se 1.307592e+00 2.229392e+00 1.747782e+00 4.223622e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Zn Hg 1.531211e+00 2.025045e+00 1.833708e+00 3.903046e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Zn S 2.434871e+00 2.423171e+00 1.711097e+00 3.095161e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Se Cd 1.300376e+00 1.804151e+00 2.124568e+00 3.259006e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Se Te 1.450015e+00 2.297301e+00 1.726905e+00 3.086266e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Se Zn 2.208390e+00 2.323783e+00 1.589241e+00 2.500815e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Se Se 1.307592e+00 2.229392e+00 1.747782e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Se Hg 1.531211e+00 2.025045e+00 1.833708e+00 3.003322e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Se S 2.434871e+00 2.423171e+00 1.711097e+00 2.381670e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Hg Cd 1.300376e+00 1.804151e+00 2.124568e+00 3.526684e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Hg Te 1.450015e+00 2.297301e+00 1.726905e+00 3.339756e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Hg Zn 2.208390e+00 2.323783e+00 1.589241e+00 2.706220e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Hg Se 1.307592e+00 2.229392e+00 1.747782e+00 3.516939e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S Hg Hg 1.531211e+00 2.025045e+00 1.833708e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S Hg S 2.434871e+00 2.423171e+00 1.711097e+00 2.577288e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00
S S Cd 1.300376e+00 1.804151e+00 2.124568e+00 4.447203e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.540087e+00 4.000000e+00 0.000000e+00 0.000000e+00
S S Te 1.450015e+00 2.297301e+00 1.726905e+00 4.211484e+01 1.200000e+00 -3.333333e-01 7.049600e+00 7.794685e-01 4.000000e+00 0.000000e+00 0.000000e+00
S S Zn 2.208390e+00 2.323783e+00 1.589241e+00 3.412585e+01 1.200000e+00 -3.333333e-01 7.049600e+00 4.643181e-01 4.000000e+00 0.000000e+00 0.000000e+00
S S Se 1.307592e+00 2.229392e+00 1.747782e+00 4.434914e+01 1.200000e+00 -3.333333e-01 7.049600e+00 6.932325e-01 4.000000e+00 0.000000e+00 0.000000e+00
S S Hg 1.531211e+00 2.025045e+00 1.833708e+00 4.098300e+01 1.200000e+00 -3.333333e-01 7.049600e+00 1.184541e+00 4.000000e+00 0.000000e+00 0.000000e+00
S S S 2.434871e+00 2.423171e+00 1.711097e+00 3.250000e+01 1.200000e+00 -3.333333e-01 7.917000e+00 1.049688e+00 4.000000e+00 0.000000e+00 0.000000e+00

View File

@ -0,0 +1,38 @@
# DATE: 2015-10-14 CONTRIBUTOR: Aidan Thompson, athomps@sandia.gov CITATION: Branicio, Rino, Gan and Tsuzuki, J. Phys Condensed Matter 21 (2009) 095002
#
# Vashishta potential file for InP, Branicio, Rino, Gan and Tsuzuki,
# J. Phys Condensed Matter 21 (2009) 095002
#
# These entries are in LAMMPS "metal" units:
# H = eV*Angstroms^eta; Zi, Zj = |e| (e = electronic charge);
# lambda1, lambda4, rc, r0, gamma = Angstroms;
# D = eV*Angstroms^4; W = eV*Angstroms^6; B = eV;
# other quantities are unitless
# element1 element2 element3
# H eta Zi Zj lambda1 D lambda4
# W rc B gamma r0 C cos(theta)
In In In 273.584 7 -1.21 -1.21 4.5 0.0 2.75
0.0 6.0 0.0 0.0 0.0 0.0 0.0
P P P 1813.06 7 1.21 1.21 4.5 52.7067 2.75
0.0 6.0 0.0 0.0 0.0 0.0 -0.333333333333
In P P 4847.09 9 1.21 -1.21 4.5 26.3533 2.75
270.105 6.0 4.34967 1.0 3.55 7.0 -0.333333333333
P In In 4847.09 9 1.21 -1.21 4.5 26.3533 2.75
270.105 6.0 4.34967 1.0 3.55 7.0 -0.333333333333
In In P 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0
In P In 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0
P In P 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0
P P In 0.0 0.0 0.0 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0

18
examples/threebody/Si.sw Normal file
View File

@ -0,0 +1,18 @@
# DATE: 2007-06-11 CONTRIBUTOR: Aidan Thompson, athomps@sandia.gov CITATION: Stillinger and Weber, Phys Rev B, 31, 5262, (1985)
# Stillinger-Weber parameters for various elements and mixtures
# multiple entries can be added to this file, LAMMPS reads the ones it needs
# these entries are in LAMMPS "metal" units:
# epsilon = eV; sigma = Angstroms
# other quantities are unitless
# format of a single entry (one or more lines):
# element 1, element 2, element 3,
# epsilon, sigma, a, lambda, gamma, costheta0, A, B, p, q, tol
# Here are the original parameters in metal units, for Silicon from:
#
# Stillinger and Weber, Phys. Rev. B, v. 31, p. 5262, (1985)
#
Si Si Si 2.1683 2.0951 1.80 21.0 1.20 -0.333333333333
7.049556277 0.6022245584 4.0 0.0 0.0

View File

@ -0,0 +1,116 @@
# Simple regression tests for threebody potentials
# NOTE: These are not intended to represent real materials
units metal
atom_style atomic
atom_modify map array
boundary p p p
atom_modify sort 0 0.0
# temperature
variable t equal 1800.0
# cubic diamond unit cell
variable a equal 5.431
lattice custom $a &
a1 1.0 0.0 0.0 &
a2 0.0 1.0 0.0 &
a3 0.0 0.0 1.0 &
basis 0.0 0.0 0.0 &
basis 0.0 0.5 0.5 &
basis 0.5 0.0 0.5 &
basis 0.5 0.5 0.0 &
basis 0.25 0.25 0.25 &
basis 0.25 0.75 0.75 &
basis 0.75 0.25 0.75 &
basis 0.75 0.75 0.25
region myreg block 0 4 &
0 4 &
0 4
create_box 8 myreg
create_atoms 1 region myreg &
basis 1 1 &
basis 2 2 &
basis 3 3 &
basis 4 4 &
basis 5 5 &
basis 6 6 &
basis 7 7 &
basis 8 8
mass * 28.06
velocity all create $t 5287287 mom yes rot yes dist gaussian
# Equilibrate using Stillinger-Weber model for silicon
pair_style sw
pair_coeff * * Si.sw Si Si Si Si Si Si Si Si
thermo 10
fix 1 all nvt temp $t $t 0.1
fix_modify 1 energy yes
timestep 1.0e-3
neighbor 1.0 bin
neigh_modify every 1 delay 10 check yes
run 100
write_restart restart.equil
# Test Stillinger-Weber model for Cd/Te/Zn/Se/Hg/S
clear
read_restart restart.equil
pair_style sw
pair_coeff * * CdTeZnSeHgS0.sw Cd Zn Hg Cd Te S Se Te
thermo 10
fix 1 all nvt temp $t $t 0.1
fix_modify 1 energy yes
timestep 1.0e-3
neighbor 1.0 bin
neigh_modify every 1 delay 10 check yes
run 100
# Test Vashishta model for In/P
clear
read_restart restart.equil
pair_style vashishta
pair_coeff * * InP.vashishta In In In In P P P P
thermo 10
fix 1 all nvt temp $t $t 0.1
fix_modify 1 energy yes
timestep 1.0e-3
neighbor 1.0 bin
neigh_modify every 1 delay 10 check yes
run 100
# Test Tersoff model for B/N/C
clear
read_restart restart.equil
variable fac equal 0.6
change_box all x scale ${fac} y scale ${fac} z scale ${fac} remap
pair_style tersoff
pair_coeff * * BNC.tersoff N N N C B B C B
thermo 10
fix 1 all nvt temp $t $t 0.1
fix_modify 1 energy yes
timestep 1.0e-3
neighbor 1.0 bin
neigh_modify every 1 delay 10 check yes
run 100

View File

@ -18,9 +18,11 @@
# sqrt(lambda_ij*epsilon_ij*lambda_ik*epsilon_ik)/lambda_ik, and the # sqrt(lambda_ij*epsilon_ij*lambda_ik*epsilon_ik)/lambda_ik, and the
# results are directly entered in this table. Obviously, this # results are directly entered in this table. Obviously, this
# conversion does not change the two-body parameters epsilon_ijj. # conversion does not change the two-body parameters epsilon_ijj.
# All other ik pair parameters are entered on the i*k line, where *
# can be any species. This is consistent with the requirement of # The twobody ik pair parameters are entered on the i*k lines, where *
# the ik parameter being on the ikk line. # can be any species. This is consistent with the LAMMPS requirement
# that twobody ik parameters be defined on the ikk line. Entries on all
# the other i*k lines are ignored by LAMMPS
# These entries are in LAMMPS "metal" units: epsilon = eV; # These entries are in LAMMPS "metal" units: epsilon = eV;
# sigma = Angstroms; other quantities are unitless # sigma = Angstroms; other quantities are unitless

View File

@ -1,9 +1,17 @@
### DATE: 2013-08-09 CONTRIBUTOR: X. W. Zhou, xzhou@sandia.gov, CITATION: Zhou, Ward, Martin, van Swol, Cruz-Campa, and D. Zubia, Phys. Rev. B, 88, 085309 (2013). ### DATE: 2013-08-09 CONTRIBUTOR: X. W. Zhou, xzhou@sandia.gov, CITATION: Zhou, Ward, Martin, van Swol, Cruz-Campa, and D. Zubia, Phys. Rev. B, 88, 085309 (2013).
# #
# Note that the way the parameters can be entered is not unique. As one way, we assume that eps_ijk is equal to eps_ik and lambda_ijk is equal to # Note that the way the parameters can be entered is not unique.
# sqrt(lambda_ij*eps_ij*lambda_ik*eps_ik)/eps_ik, and all other parameters in the ijk line are for ik. # As one way, we assume that eps_ijk is equal to eps_ik and
# lambda_ijk is equal to sqrt(lambda_ij*eps_ij*lambda_ik*eps_ik)/eps_ik,
# and all other parameters in the ijk line are for ik.
# #
# These entries are in LAMMPS "metal" units: epsilon = eV; sigma = Angstroms; other quantities are unitless; # The twobody ik pair parameters are entered on the i*k lines, where *
# can be any species. This is consistent with the LAMMPS requirement
# that twobody ik parameters be defined on the ikk line. Entries on all
# the other i*k lines are ignored by LAMMPS
#
# These entries are in LAMMPS "metal" units: epsilon = eV;
# sigma = Angstroms; other quantities are unitless
# #
# cutoff distance = 4.632 # cutoff distance = 4.632
# eps sigma a lambda gamma cos(theta) A B p q tol # eps sigma a lambda gamma cos(theta) A B p q tol

View File

@ -37,7 +37,7 @@
#include "math_const.h" #include "math_const.h"
#include "memory.h" #include "memory.h"
#include "error.h" #include "error.h"
#include "pair_kokkos.h" #include "pair_reax_c_kokkos.h"
using namespace LAMMPS_NS; using namespace LAMMPS_NS;
using namespace FixConst; using namespace FixConst;
@ -62,6 +62,11 @@ FixQEqReaxKokkos<DeviceType>::FixQEqReaxKokkos(LAMMPS *lmp, int narg, char **arg
nmax = nmax = m_cap = 0; nmax = nmax = m_cap = 0;
allocated_flag = 0; allocated_flag = 0;
reaxc = (PairReaxC *) force->pair_match("reax/c/kk",1);
use_pair_list = 0;
if (reaxc->execution_space == this->execution_space)
use_pair_list = 1;
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
@ -80,11 +85,35 @@ void FixQEqReaxKokkos<DeviceType>::init()
atomKK->k_q.modify<LMPHostType>(); atomKK->k_q.modify<LMPHostType>();
atomKK->k_q.sync<LMPDeviceType>(); atomKK->k_q.sync<LMPDeviceType>();
FixQEqReax::init(); //FixQEqReax::init();
{
if (!atom->q_flag) error->all(FLERR,"Fix qeq/reax requires atom attribute q");
ngroup = group->count(igroup);
if (ngroup == 0) error->all(FLERR,"Fix qeq/reax group has no atoms");
// need a half neighbor list w/ Newton off and ghost neighbors
// built whenever re-neighboring occurs
if (!use_pair_list) {
int irequest = neighbor->request(this,instance_me);
neighbor->requests[irequest]->pair = 0;
neighbor->requests[irequest]->fix = 1;
neighbor->requests[irequest]->newton = 2;
neighbor->requests[irequest]->ghost = 1;
}
init_shielding();
init_taper();
if (strstr(update->integrate_style,"respa"))
nlevels_respa = ((Respa *) update->integrate)->nlevels;
}
neighflag = lmp->kokkos->neighflag; neighflag = lmp->kokkos->neighflag;
int irequest = neighbor->nrequest - 1;
if (!use_pair_list) {
int irequest = neighbor->nrequest - 1;
neighbor->requests[irequest]-> neighbor->requests[irequest]->
kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value && kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value &&
!Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value; !Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
@ -104,6 +133,7 @@ void FixQEqReaxKokkos<DeviceType>::init()
neighbor->requests[irequest]->full_cluster = 0; neighbor->requests[irequest]->full_cluster = 0;
neighbor->requests[irequest]->ghost = 1; neighbor->requests[irequest]->ghost = 1;
} }
}
int ntypes = atom->ntypes; int ntypes = atom->ntypes;
k_params = Kokkos::DualView<params_qeq*,Kokkos::LayoutRight,DeviceType> k_params = Kokkos::DualView<params_qeq*,Kokkos::LayoutRight,DeviceType>
@ -213,6 +243,8 @@ void FixQEqReaxKokkos<DeviceType>::pre_force(int vflag)
k_shield.template sync<DeviceType>(); k_shield.template sync<DeviceType>();
k_tap.template sync<DeviceType>(); k_tap.template sync<DeviceType>();
if (use_pair_list)
list = reaxc->list;
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list); NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
d_numneigh = k_list->d_numneigh; d_numneigh = k_list->d_numneigh;
d_neighbors = k_list->d_neighbors; d_neighbors = k_list->d_neighbors;

View File

@ -146,7 +146,7 @@ class FixQEqReaxKokkos : public FixQEqReax {
double memory_usage(); double memory_usage();
protected: protected:
int inum; int inum,use_pair_list;
int allocated_flag; int allocated_flag;
typedef Kokkos::DualView<int***,DeviceType> tdual_int_1d; typedef Kokkos::DualView<int***,DeviceType> tdual_int_1d;

View File

@ -71,12 +71,6 @@ void FixReaxCSpeciesKokkos::init()
"pair_style reax/c/kk"); "pair_style reax/c/kk");
FixReaxCSpecies::init(); FixReaxCSpecies::init();
int irequest = neighbor->request(this,instance_me);
neighbor->requests[irequest]->pair = 0;
neighbor->requests[irequest]->fix = 1;
neighbor->requests[irequest]->newton = 2;
neighbor->requests[irequest]->ghost = 1;
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
@ -86,12 +80,20 @@ void FixReaxCSpeciesKokkos::FindMolecule()
int i,j,ii,jj,inum,itype,jtype,loop,looptot; int i,j,ii,jj,inum,itype,jtype,loop,looptot;
int change,done,anychange; int change,done,anychange;
int *mask = atom->mask; int *mask = atom->mask;
int *ilist;
double bo_tmp,bo_cut; double bo_tmp,bo_cut;
double **spec_atom = f_SPECBOND->array_atom; double **spec_atom = f_SPECBOND->array_atom;
inum = list->inum; inum = reaxc->list->inum;
ilist = list->ilist; typename ArrayTypes<LMPHostType>::t_int_1d ilist;
if (reaxc->execution_space == Host) {
NeighListKokkos<LMPHostType>* k_list = static_cast<NeighListKokkos<LMPHostType>*>(reaxc->list);
k_list->k_ilist.sync<LMPHostType>();
ilist = k_list->k_ilist.h_view;
} else {
NeighListKokkos<LMPDeviceType>* k_list = static_cast<NeighListKokkos<LMPDeviceType>*>(reaxc->list);
k_list->k_ilist.sync<LMPHostType>();
ilist = k_list->k_ilist.h_view;
}
for (ii = 0; ii < inum; ii++) { for (ii = 0; ii < inum; ii++) {
i = ilist[ii]; i = ilist[ii];

View File

@ -41,8 +41,6 @@ class FixWallReflectKokkos : public FixWallReflect {
void operator()(TagFixWallReflectPostIntegrate, const int&) const; void operator()(TagFixWallReflectPostIntegrate, const int&) const;
protected: protected:
class AtomKokkos *atomKK;
typename AT::t_x_array x; typename AT::t_x_array x;
typename AT::t_v_array v; typename AT::t_v_array v;
typename AT::t_int_1d_randomread mask; typename AT::t_int_1d_randomread mask;

View File

@ -164,6 +164,7 @@ if (GHOST) {
list->gnum = 0; list->gnum = 0;
} }
list->k_ilist.template modify<DeviceType>();
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */

View File

@ -49,8 +49,9 @@ void NeighListKokkos<Device>::grow(int nmax)
if (nmax <= maxatoms) return; if (nmax <= maxatoms) return;
maxatoms = nmax; maxatoms = nmax;
d_ilist = k_ilist =
typename ArrayTypes<Device>::t_int_1d("neighlist:ilist",maxatoms); DAT::tdual_int_1d("neighlist:ilist",maxatoms);
d_ilist = k_ilist.view<Device>();
d_numneigh = d_numneigh =
typename ArrayTypes<Device>::t_int_1d("neighlist:numneigh",maxatoms); typename ArrayTypes<Device>::t_int_1d("neighlist:numneigh",maxatoms);
d_neighbors = d_neighbors =

View File

@ -71,7 +71,8 @@ public:
void clean_copy(); void clean_copy();
void grow(int nmax); void grow(int nmax);
typename ArrayTypes<Device>::t_neighbors_2d d_neighbors; typename ArrayTypes<Device>::t_neighbors_2d d_neighbors;
typename ArrayTypes<Device>::t_int_1d d_ilist; // local indices of I atoms typename DAT::tdual_int_1d k_ilist; // local indices of I atoms
typename ArrayTypes<Device>::t_int_1d d_ilist;
typename ArrayTypes<Device>::t_int_1d d_numneigh; // # of J neighs for each I typename ArrayTypes<Device>::t_int_1d d_numneigh; // # of J neighs for each I
typename ArrayTypes<Device>::t_int_1d d_stencil; // # of J neighs for each I typename ArrayTypes<Device>::t_int_1d d_stencil; // # of J neighs for each I
typename ArrayTypes<LMPHostType>::t_int_1d h_stencil; // # of J neighs per I typename ArrayTypes<LMPHostType>::t_int_1d h_stencil; // # of J neighs per I

View File

@ -426,8 +426,6 @@ class PairReaxCKokkos : public PairReaxC {
typename AT::t_ffloat_2d_dl d_sum_ovun; typename AT::t_ffloat_2d_dl d_sum_ovun;
typename AT::t_ffloat_2d_dl d_dBOpx, d_dBOpy, d_dBOpz; typename AT::t_ffloat_2d_dl d_dBOpx, d_dBOpy, d_dBOpz;
class AtomKokkos *atomKK;
int neighflag,newton_pair, maxnumneigh, maxhb, maxbo; int neighflag,newton_pair, maxnumneigh, maxhb, maxbo;
int nlocal,nall,eflag,vflag; int nlocal,nall,eflag,vflag;
F_FLOAT cut_nbsq, cut_hbsq, cut_bosq, bo_cut, thb_cut, thb_cutsq; F_FLOAT cut_nbsq, cut_hbsq, cut_bosq, bo_cut, thb_cut, thb_cutsq;

View File

@ -121,6 +121,18 @@ void PairSWKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
EV_FLOAT ev; EV_FLOAT ev;
EV_FLOAT ev_all; EV_FLOAT ev_all;
// build short neighbor list
int max_neighs = d_neighbors.dimension_1();
if ((d_neighbors_short.dimension_1() != max_neighs) ||
(d_neighbors_short.dimension_0() != ignum)) {
d_neighbors_short = Kokkos::View<int**,DeviceType>("SW::neighbors_short",ignum,max_neighs);
}
if (d_numneigh_short.dimension_0()!=ignum)
d_numneigh_short = Kokkos::View<int*,DeviceType>("SW::numneighs_short",ignum);
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagPairSWComputeShortNeigh>(0,neighflag==FULL?ignum:inum), *this);
// loop over neighbor list of my atoms // loop over neighbor list of my atoms
if (neighflag == HALF) { if (neighflag == HALF) {
@ -174,6 +186,38 @@ void PairSWKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
copymode = 0; copymode = 0;
} }
/* ---------------------------------------------------------------------- */
template<class DeviceType>
KOKKOS_INLINE_FUNCTION
void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeShortNeigh, const int& ii) const {
const int i = d_ilist[ii];
const X_FLOAT xtmp = x(i,0);
const X_FLOAT ytmp = x(i,1);
const X_FLOAT ztmp = x(i,2);
const int jnum = d_numneigh[i];
int inside = 0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj);
j &= NEIGHMASK;
const X_FLOAT delx = xtmp - x(j,0);
const X_FLOAT dely = ytmp - x(j,1);
const X_FLOAT delz = ztmp - x(j,2);
const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cutmax*cutmax) {
d_neighbors_short(i,inside) = j;
inside++;
}
}
d_numneigh_short(i) = inside;
}
/* ---------------------------------------------------------------------- */
template<class DeviceType> template<class DeviceType>
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
@ -196,14 +240,14 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeHalf<NEIGHFLAG,EVFLAG>
// two-body interactions, skip half of them // two-body interactions, skip half of them
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT fxtmpi = 0.0; F_FLOAT fxtmpi = 0.0;
F_FLOAT fytmpi = 0.0; F_FLOAT fytmpi = 0.0;
F_FLOAT fztmpi = 0.0; F_FLOAT fztmpi = 0.0;
for (int jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const tagint jtag = tag[j]; const tagint jtag = tag[j];
@ -245,7 +289,7 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeHalf<NEIGHFLAG,EVFLAG>
const int jnumm1 = jnum - 1; const int jnumm1 = jnum - 1;
for (int jj = 0; jj < jnumm1; jj++) { for (int jj = 0; jj < jnumm1; jj++) {
int j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const int jtype = d_map[type[j]]; const int jtype = d_map[type[j]];
const int ijparam = d_elem2param(itype,jtype,jtype); const int ijparam = d_elem2param(itype,jtype,jtype);
@ -260,7 +304,7 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeHalf<NEIGHFLAG,EVFLAG>
F_FLOAT fztmpj = 0.0; F_FLOAT fztmpj = 0.0;
for (int kk = jj+1; kk < jnum; kk++) { for (int kk = jj+1; kk < jnum; kk++) {
int k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
const int ktype = d_map[type[k]]; const int ktype = d_map[type[k]];
const int ikparam = d_elem2param(itype,ktype,ktype); const int ikparam = d_elem2param(itype,ktype,ktype);
@ -331,14 +375,14 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeFullA<NEIGHFLAG,EVFLAG
// two-body interactions // two-body interactions
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT fxtmpi = 0.0; F_FLOAT fxtmpi = 0.0;
F_FLOAT fytmpi = 0.0; F_FLOAT fytmpi = 0.0;
F_FLOAT fztmpi = 0.0; F_FLOAT fztmpi = 0.0;
for (int jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const tagint jtag = tag[j]; const tagint jtag = tag[j];
@ -368,7 +412,7 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeFullA<NEIGHFLAG,EVFLAG
const int jnumm1 = jnum - 1; const int jnumm1 = jnum - 1;
for (int jj = 0; jj < jnumm1; jj++) { for (int jj = 0; jj < jnumm1; jj++) {
int j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const int jtype = d_map[type[j]]; const int jtype = d_map[type[j]];
const int ijparam = d_elem2param(itype,jtype,jtype); const int ijparam = d_elem2param(itype,jtype,jtype);
@ -380,7 +424,7 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeFullA<NEIGHFLAG,EVFLAG
if (rsq1 >= d_params[ijparam].cutsq) continue; if (rsq1 >= d_params[ijparam].cutsq) continue;
for (int kk = jj+1; kk < jnum; kk++) { for (int kk = jj+1; kk < jnum; kk++) {
int k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
const int ktype = d_map[type[k]]; const int ktype = d_map[type[k]];
const int ikparam = d_elem2param(itype,ktype,ktype); const int ikparam = d_elem2param(itype,ktype,ktype);
@ -438,14 +482,14 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeFullB<NEIGHFLAG,EVFLAG
const X_FLOAT ytmpi = x(i,1); const X_FLOAT ytmpi = x(i,1);
const X_FLOAT ztmpi = x(i,2); const X_FLOAT ztmpi = x(i,2);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT fxtmpi = 0.0; F_FLOAT fxtmpi = 0.0;
F_FLOAT fytmpi = 0.0; F_FLOAT fytmpi = 0.0;
F_FLOAT fztmpi = 0.0; F_FLOAT fztmpi = 0.0;
for (int jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
if (j >= nlocal) continue; if (j >= nlocal) continue;
const int jtype = d_map[type[j]]; const int jtype = d_map[type[j]];
@ -461,10 +505,10 @@ void PairSWKokkos<DeviceType>::operator()(TagPairSWComputeFullB<NEIGHFLAG,EVFLAG
if (rsq1 >= d_params[jiparam].cutsq) continue; if (rsq1 >= d_params[jiparam].cutsq) continue;
const int j_jnum = d_numneigh[j]; const int j_jnum = d_numneigh_short[j];
for (int kk = 0; kk < j_jnum; kk++) { for (int kk = 0; kk < j_jnum; kk++) {
int k = d_neighbors(j,kk); int k = d_neighbors_short(j,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
if (k == i) continue; if (k == i) continue;
const int ktype = d_map[type[k]]; const int ktype = d_map[type[k]];

View File

@ -34,6 +34,8 @@ struct TagPairSWComputeFullA{};
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
struct TagPairSWComputeFullB{}; struct TagPairSWComputeFullB{};
struct TagPairSWComputeShortNeigh{};
namespace LAMMPS_NS { namespace LAMMPS_NS {
template<class DeviceType> template<class DeviceType>
@ -75,6 +77,9 @@ class PairSWKokkos : public PairSW {
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
void operator()(TagPairSWComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const; void operator()(TagPairSWComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const;
KOKKOS_INLINE_FUNCTION
void operator()(TagPairSWComputeShortNeigh, const int&) const;
template<int NEIGHFLAG> template<int NEIGHFLAG>
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
void ev_tally(EV_FLOAT &ev, const int &i, const int &j, void ev_tally(EV_FLOAT &ev, const int &i, const int &j,
@ -136,6 +141,9 @@ class PairSWKokkos : public PairSW {
int nlocal,nall,eflag,vflag; int nlocal,nall,eflag,vflag;
int inum; int inum;
Kokkos::View<int**,DeviceType> d_neighbors_short;
Kokkos::View<int*,DeviceType> d_numneigh_short;
friend void pair_virial_fdotr_compute<PairSWKokkos>(PairSWKokkos*); friend void pair_virial_fdotr_compute<PairSWKokkos>(PairSWKokkos*);
}; };

View File

@ -12,7 +12,7 @@
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
Contributing author: Ray Shan (SNL) Contributing author: Ray Shan (SNL) and Christian Trott (SNL)
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
#include <math.h> #include <math.h>
@ -191,7 +191,7 @@ void PairTersoffKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
nall = atom->nlocal + atom->nghost; nall = atom->nlocal + atom->nghost;
newton_pair = force->newton_pair; newton_pair = force->newton_pair;
const int inum = list->inum; inum = list->inum;
const int ignum = inum + list->gnum; const int ignum = inum + list->gnum;
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list); NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
d_numneigh = k_list->d_numneigh; d_numneigh = k_list->d_numneigh;
@ -204,6 +204,18 @@ void PairTersoffKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
EV_FLOAT ev; EV_FLOAT ev;
EV_FLOAT ev_all; EV_FLOAT ev_all;
// build short neighbor list
int max_neighs = d_neighbors.dimension_1();
if ((d_neighbors_short.dimension_1() != max_neighs) ||
(d_neighbors_short.dimension_0() != ignum)) {
d_neighbors_short = Kokkos::View<int**,DeviceType>("Tersoff::neighbors_short",ignum,max_neighs);
}
if (d_numneigh_short.dimension_0()!=ignum)
d_numneigh_short = Kokkos::View<int*,DeviceType>("Tersoff::numneighs_short",ignum);
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagPairTersoffComputeShortNeigh>(0,neighflag==FULL?ignum:inum), *this);
if (neighflag == HALF) { if (neighflag == HALF) {
if (evflag) if (evflag)
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffComputeHalf<HALF,1> >(0,inum),*this,ev); Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffComputeHalf<HALF,1> >(0,inum),*this,ev);
@ -257,6 +269,35 @@ void PairTersoffKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
template<class DeviceType>
KOKKOS_INLINE_FUNCTION
void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeShortNeigh, const int& ii) const {
const int i = d_ilist[ii];
const X_FLOAT xtmp = x(i,0);
const X_FLOAT ytmp = x(i,1);
const X_FLOAT ztmp = x(i,2);
const int jnum = d_numneigh[i];
int inside = 0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj);
j &= NEIGHMASK;
const X_FLOAT delx = xtmp - x(j,0);
const X_FLOAT dely = ytmp - x(j,1);
const X_FLOAT delz = ztmp - x(j,2);
const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cutmax*cutmax) {
d_neighbors_short(i,inside) = j;
inside++;
}
}
d_numneigh_short(i) = inside;
}
/* ---------------------------------------------------------------------- */
template<class DeviceType> template<class DeviceType>
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
@ -273,21 +314,22 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
const int itype = type(i); const int itype = type(i);
const int itag = tag(i); const int itag = tag(i);
int j,k,jj,kk,jtag,jtype,ktype;
F_FLOAT rsq1, cutsq1, rsq2, cutsq2, rij, rik, bo_ij;
F_FLOAT fi[3], fj[3], fk[3]; F_FLOAT fi[3], fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
//const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i); //const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
// repulsive // repulsive
for (jj = 0; jj < jnum; jj++) { F_FLOAT f_x = 0.0;
j = d_neighbors(i,jj); F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
jtag = tag(j); const int jtag = tag(j);
if (itag > jtag) { if (itag > jtag) {
if ((itag+jtag) % 2 == 0) continue; if ((itag+jtag) % 2 == 0) continue;
@ -315,9 +357,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
(tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r; (tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r;
const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp; const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp;
a_f(i,0) += delx*frep; f_x += delx*frep;
a_f(i,1) += dely*frep; f_y += dely*frep;
a_f(i,2) += delz*frep; f_z += delz*frep;
a_f(j,0) -= delx*frep; a_f(j,0) -= delx*frep;
a_f(j,1) -= dely*frep; a_f(j,1) -= dely*frep;
a_f(j,2) -= delz*frep; a_f(j,2) -= delz*frep;
@ -330,35 +372,35 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
delx1 = xtmp - x(j,0); const F_FLOAT delx1 = xtmp - x(j,0);
dely1 = ytmp - x(j,1); const F_FLOAT dely1 = ytmp - x(j,1);
delz1 = ztmp - x(j,2); const F_FLOAT delz1 = ztmp - x(j,2);
rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1; const F_FLOAT rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1;
cutsq1 = paramskk(itype,jtype,jtype).cutsq; const F_FLOAT cutsq1 = paramskk(itype,jtype,jtype).cutsq;
bo_ij = 0.0; F_FLOAT bo_ij = 0.0;
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); const F_FLOAT rij = sqrt(rsq1);
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2); bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2);
} }
@ -369,16 +411,16 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij); const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij);
const F_FLOAT fatt = -0.5*bij * dfa / rij; const F_FLOAT fatt = -0.5*bij * dfa / rij;
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa;
a_f(i,0) += delx1*fatt; f_x += delx1*fatt;
a_f(i,1) += dely1*fatt; f_y += dely1*fatt;
a_f(i,2) += delz1*fatt; f_z += delz1*fatt;
a_f(j,0) -= delx1*fatt; F_FLOAT fj_x = -delx1*fatt;
a_f(j,1) -= dely1*fatt; F_FLOAT fj_y = -dely1*fatt;
a_f(j,2) -= delz1*fatt; F_FLOAT fj_z = -delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
const F_FLOAT eng = 0.5*bij * fa;
if (eflag) ev.evdwl += eng; if (eflag) ev.evdwl += eng;
if (vflag_either || eflag_atom) if (vflag_either || eflag_atom)
this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1); this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1);
@ -386,29 +428,29 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
a_f(i,0) += fi[0]; f_x += fi[0];
a_f(i,1) += fi[1]; f_y += fi[1];
a_f(i,2) += fi[2]; f_z += fi[2];
a_f(j,0) += fj[0]; fj_x += fj[0];
a_f(j,1) += fj[1]; fj_y += fj[1];
a_f(j,2) += fj[2]; fj_z += fj[2];
a_f(k,0) += fk[0]; a_f(k,0) += fk[0];
a_f(k,1) += fk[1]; a_f(k,1) += fk[1];
a_f(k,2) += fk[2]; a_f(k,2) += fk[2];
@ -420,7 +462,13 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeHalf<NEIGHFL
if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik); if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik);
} }
} }
a_f(j,0) += fj_x;
a_f(j,1) += fj_y;
a_f(j,2) += fj_z;
} }
a_f(i,0) += f_x;
a_f(i,1) += f_y;
a_f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -450,12 +498,15 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2; X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
//const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i); //const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
// repulsive // repulsive
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const int jtype = type(j); const int jtype = type(j);
@ -475,9 +526,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
(tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r; (tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r;
const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp; const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp;
f(i,0) += delx*frep; f_x += delx*frep;
f(i,1) += dely*frep; f_y += dely*frep;
f(i,2) += delz*frep; f_z += delz*frep;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -490,7 +541,7 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); jtype = type(j);
@ -506,7 +557,7 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -530,9 +581,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) += delx1*fatt; f_x += delx1*fatt;
f(i,1) += dely1*fatt; f_y += dely1*fatt;
f(i,2) += delz1*fatt; f_z += delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) ev.evdwl += 0.5*eng; if (eflag) ev.evdwl += 0.5*eng;
@ -544,7 +595,7 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -559,9 +610,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
f(i,0) += fi[0]; f_x += fi[0];
f(i,1) += fi[1]; f_y += fi[1];
f(i,2) += fi[2]; f_z += fi[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrij[3], delrik[3]; F_FLOAT delrij[3], delrik[3];
@ -571,6 +622,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullA<NEIGHF
} }
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -599,12 +653,16 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
F_FLOAT fj[3], fk[3]; F_FLOAT fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2; X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
if (j >= nlocal) continue; if (j >= nlocal) continue;
jtype = type(j); jtype = type(j);
@ -619,10 +677,10 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); rij = sqrt(rsq1);
j_jnum = d_numneigh[j]; j_jnum = d_numneigh_short[j];
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -648,9 +706,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) -= delx1*fatt; f_x -= delx1*fatt;
f(i,1) -= dely1*fatt; f_y -= dely1*fatt;
f(i,2) -= delz1*fatt; f_z -= delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -662,7 +720,7 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -677,9 +735,9 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
rik = sqrt(rsq2); rik = sqrt(rsq2);
ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fj,fk); rik,delx2,dely2,delz2,fj,fk);
f(i,0) += fj[0]; f_x += fj[0];
f(i,1) += fj[1]; f_y += fj[1];
f(i,2) += fj[2]; f_z += fj[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrji[3], delrjk[3]; F_FLOAT delrji[3], delrjk[3];
@ -692,11 +750,14 @@ void PairTersoffKokkos<DeviceType>::operator()(TagPairTersoffComputeFullB<NEIGHF
const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij); const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij);
ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2, ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2,
rij,delx1,dely1,delz1,fk); rij,delx1,dely1,delz1,fk);
f(i,0) += fk[0]; f_x += fk[0];
f(i,1) += fk[1]; f_y += fk[1];
f(i,2) += fk[2]; f_z += fk[2];
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -749,8 +810,9 @@ double PairTersoffKokkos<DeviceType>::bondorder(const int &i, const int &j, cons
const F_FLOAT costheta = (dx1*dx2 + dy1*dy2 + dz1*dz2)/(rij*rik); const F_FLOAT costheta = (dx1*dx2 + dy1*dy2 + dz1*dz2)/(rij*rik);
if (int(paramskk(i,j,k).powerm) == 3) arg = pow(paramskk(i,j,k).lam3 * (rij-rik),3.0); const F_FLOAT param = paramskk(i,j,k).lam3 * (rij-rik);
else arg = paramskk(i,j,k).lam3 * (rij-rik); if (int(paramskk(i,j,k).powerm) == 3) arg = param*param*param;//pow(paramskk(i,j,k).lam3 * (rij-rik),3.0);
else arg = param;
if (arg > 69.0776) ex_delr = 1.e30; if (arg > 69.0776) ex_delr = 1.e30;
else if (arg < -69.0776) ex_delr = 0.0; else if (arg < -69.0776) ex_delr = 0.0;
@ -780,7 +842,6 @@ KOKKOS_INLINE_FUNCTION
double PairTersoffKokkos<DeviceType>:: double PairTersoffKokkos<DeviceType>::
ters_dgijk(const int &i, const int &j, const int &k, const F_FLOAT &cos) const ters_dgijk(const int &i, const int &j, const int &k, const F_FLOAT &cos) const
{ {
const F_FLOAT ters_c = paramskk(i,j,k).c * paramskk(i,j,k).c; const F_FLOAT ters_c = paramskk(i,j,k).c * paramskk(i,j,k).c;
const F_FLOAT ters_d = paramskk(i,j,k).d * paramskk(i,j,k).d; const F_FLOAT ters_d = paramskk(i,j,k).d * paramskk(i,j,k).d;
const F_FLOAT hcth = paramskk(i,j,k).h - cos; const F_FLOAT hcth = paramskk(i,j,k).h - cos;
@ -838,9 +899,9 @@ double PairTersoffKokkos<DeviceType>::ters_dbij(const int &i, const int &j,
const int &k, const F_FLOAT &bo) const const int &k, const F_FLOAT &bo) const
{ {
const F_FLOAT tmp = paramskk(i,j,k).beta * bo; const F_FLOAT tmp = paramskk(i,j,k).beta * bo;
if (tmp > paramskk(i,j,k).c1) return paramskk(i,j,k).beta * -0.5*pow(tmp,-1.5); if (tmp > paramskk(i,j,k).c1) return paramskk(i,j,k).beta * -0.5/sqrt(tmp*tmp);//*pow(tmp,-1.5);
if (tmp > paramskk(i,j,k).c2) if (tmp > paramskk(i,j,k).c2)
return paramskk(i,j,k).beta * (-0.5*pow(tmp,-1.5) * return paramskk(i,j,k).beta * (-0.5/sqrt(tmp*tmp) * //*pow(tmp,-1.5) *
(1.0 - 0.5*(1.0 + 1.0/(2.0*paramskk(i,j,k).powern)) * (1.0 - 0.5*(1.0 + 1.0/(2.0*paramskk(i,j,k).powern)) *
pow(tmp,-paramskk(i,j,k).powern))); pow(tmp,-paramskk(i,j,k).powern)));
if (tmp < paramskk(i,j,k).c4) return 0.0; if (tmp < paramskk(i,j,k).c4) return 0.0;
@ -883,15 +944,16 @@ void PairTersoffKokkos<DeviceType>::ters_dthb(
fc = ters_fc_k(i,j,k,rik); fc = ters_fc_k(i,j,k,rik);
dfc = ters_dfc(i,j,k,rik); dfc = ters_dfc(i,j,k,rik);
if (int(paramskk(i,j,k).powerm) == 3) tmp = pow(paramskk(i,j,k).lam3 * (rij-rik),3.0); const F_FLOAT param = paramskk(i,j,k).lam3 * (rij-rik);
else tmp = paramskk(i,j,k).lam3 * (rij-rik); if (int(paramskk(i,j,k).powerm) == 3) tmp = param*param*param;//pow(paramskk(i,j,k).lam3 * (rij-rik),3.0);
else tmp = param;
if (tmp > 69.0776) ex_delr = 1.e30; if (tmp > 69.0776) ex_delr = 1.e30;
else if (tmp < -69.0776) ex_delr = 0.0; else if (tmp < -69.0776) ex_delr = 0.0;
else ex_delr = exp(tmp); else ex_delr = exp(tmp);
if (int(paramskk(i,j,k).powerm) == 3) if (int(paramskk(i,j,k).powerm) == 3)
dex_delr = 3.0*pow(paramskk(i,j,k).lam3,3.0) * pow(rij-rik,2.0)*ex_delr; dex_delr = 3.0*param*param*paramskk(i,j,k).lam3*ex_delr;//pow(rij-rik,2.0)*ex_delr;
else dex_delr = paramskk(i,j,k).lam3 * ex_delr; else dex_delr = paramskk(i,j,k).lam3 * ex_delr;
cos = vec3_dot(rij_hat,rik_hat); cos = vec3_dot(rij_hat,rik_hat);
@ -951,15 +1013,16 @@ void PairTersoffKokkos<DeviceType>::ters_dthbj(
fc = ters_fc_k(i,j,k,rik); fc = ters_fc_k(i,j,k,rik);
dfc = ters_dfc(i,j,k,rik); dfc = ters_dfc(i,j,k,rik);
if (int(paramskk(i,j,k).powerm) == 3) tmp = pow(paramskk(i,j,k).lam3 * (rij-rik),3.0); const F_FLOAT param = paramskk(i,j,k).lam3 * (rij-rik);
else tmp = paramskk(i,j,k).lam3 * (rij-rik); if (int(paramskk(i,j,k).powerm) == 3) tmp = param*param*param;//pow(paramskk(i,j,k).lam3 * (rij-rik),3.0);
else tmp = param;//paramskk(i,j,k).lam3 * (rij-rik);
if (tmp > 69.0776) ex_delr = 1.e30; if (tmp > 69.0776) ex_delr = 1.e30;
else if (tmp < -69.0776) ex_delr = 0.0; else if (tmp < -69.0776) ex_delr = 0.0;
else ex_delr = exp(tmp); else ex_delr = exp(tmp);
if (int(paramskk(i,j,k).powerm) == 3) if (int(paramskk(i,j,k).powerm) == 3)
dex_delr = 3.0*pow(paramskk(i,j,k).lam3,3.0) * pow(rij-rik,2.0)*ex_delr; dex_delr = 3.0*param*param*paramskk(i,j,k).lam3*ex_delr;//pow(paramskk(i,j,k).lam3,3.0) * pow(rij-rik,2.0)*ex_delr;
else dex_delr = paramskk(i,j,k).lam3 * ex_delr; else dex_delr = paramskk(i,j,k).lam3 * ex_delr;
cos = vec3_dot(rij_hat,rik_hat); cos = vec3_dot(rij_hat,rik_hat);
@ -1012,15 +1075,16 @@ void PairTersoffKokkos<DeviceType>::ters_dthbk(
fc = ters_fc_k(i,j,k,rik); fc = ters_fc_k(i,j,k,rik);
dfc = ters_dfc(i,j,k,rik); dfc = ters_dfc(i,j,k,rik);
if (int(paramskk(i,j,k).powerm) == 3) tmp = pow(paramskk(i,j,k).lam3 * (rij-rik),3.0); const F_FLOAT param = paramskk(i,j,k).lam3 * (rij-rik);
else tmp = paramskk(i,j,k).lam3 * (rij-rik); if (int(paramskk(i,j,k).powerm) == 3) tmp = param*param*param;//pow(paramskk(i,j,k).lam3 * (rij-rik),3.0);
else tmp = param;//paramskk(i,j,k).lam3 * (rij-rik);
if (tmp > 69.0776) ex_delr = 1.e30; if (tmp > 69.0776) ex_delr = 1.e30;
else if (tmp < -69.0776) ex_delr = 0.0; else if (tmp < -69.0776) ex_delr = 0.0;
else ex_delr = exp(tmp); else ex_delr = exp(tmp);
if (int(paramskk(i,j,k).powerm) == 3) if (int(paramskk(i,j,k).powerm) == 3)
dex_delr = 3.0*pow(paramskk(i,j,k).lam3,3.0) * pow(rij-rik,2.0)*ex_delr; dex_delr = 3.0*param*param*paramskk(i,j,k).lam3*ex_delr;//pow(paramskk(i,j,k).lam3,3.0) * pow(rij-rik,2.0)*ex_delr;
else dex_delr = paramskk(i,j,k).lam3 * ex_delr; else dex_delr = paramskk(i,j,k).lam3 * ex_delr;
cos = vec3_dot(rij_hat,rik_hat); cos = vec3_dot(rij_hat,rik_hat);

View File

@ -39,6 +39,8 @@ struct TagPairTersoffComputeFullA{};
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
struct TagPairTersoffComputeFullB{}; struct TagPairTersoffComputeFullB{};
struct TagPairTersoffComputeShortNeigh{};
template<class DeviceType> template<class DeviceType>
class PairTersoffKokkos : public PairTersoff { class PairTersoffKokkos : public PairTersoff {
public: public:
@ -77,6 +79,9 @@ class PairTersoffKokkos : public PairTersoff {
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const; void operator()(TagPairTersoffComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const;
KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffComputeShortNeigh, const int&) const;
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const; double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const;
@ -184,6 +189,7 @@ class PairTersoffKokkos : public PairTersoff {
// hardwired to space for 15 atom types // hardwired to space for 15 atom types
//params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1]; //params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1];
int inum;
typename AT::t_x_array_randomread x; typename AT::t_x_array_randomread x;
typename AT::t_f_array f; typename AT::t_f_array f;
typename AT::t_int_1d_randomread type; typename AT::t_int_1d_randomread type;
@ -203,10 +209,12 @@ class PairTersoffKokkos : public PairTersoff {
typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh; typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh;
//NeighListKokkos<DeviceType> k_list; //NeighListKokkos<DeviceType> k_list;
class AtomKokkos *atomKK;
int neighflag,newton_pair; int neighflag,newton_pair;
int nlocal,nall,eflag,vflag; int nlocal,nall,eflag,vflag;
Kokkos::View<int**,DeviceType> d_neighbors_short;
Kokkos::View<int*,DeviceType> d_numneigh_short;
friend void pair_virial_fdotr_compute<PairTersoffKokkos>(PairTersoffKokkos*); friend void pair_virial_fdotr_compute<PairTersoffKokkos>(PairTersoffKokkos*);
}; };

View File

@ -191,7 +191,7 @@ void PairTersoffMODKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
nall = atom->nlocal + atom->nghost; nall = atom->nlocal + atom->nghost;
newton_pair = force->newton_pair; newton_pair = force->newton_pair;
const int inum = list->inum; inum = list->inum;
const int ignum = inum + list->gnum; const int ignum = inum + list->gnum;
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list); NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
d_numneigh = k_list->d_numneigh; d_numneigh = k_list->d_numneigh;
@ -204,6 +204,18 @@ void PairTersoffMODKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
EV_FLOAT ev; EV_FLOAT ev;
EV_FLOAT ev_all; EV_FLOAT ev_all;
// build short neighbor list
int max_neighs = d_neighbors.dimension_1();
if ((d_neighbors_short.dimension_1() != max_neighs) ||
(d_neighbors_short.dimension_0() != ignum)) {
d_neighbors_short = Kokkos::View<int**,DeviceType>("Tersoff::neighbors_short",ignum,max_neighs);
}
if (d_numneigh_short.dimension_0()!=ignum)
d_numneigh_short = Kokkos::View<int*,DeviceType>("Tersoff::numneighs_short",ignum);
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagPairTersoffMODComputeShortNeigh>(0,neighflag==FULL?ignum:inum), *this);
if (neighflag == HALF) { if (neighflag == HALF) {
if (evflag) if (evflag)
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffMODComputeHalf<HALF,1> >(0,inum),*this,ev); Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffMODComputeHalf<HALF,1> >(0,inum),*this,ev);
@ -257,6 +269,35 @@ void PairTersoffMODKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
template<class DeviceType>
KOKKOS_INLINE_FUNCTION
void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeShortNeigh, const int& ii) const {
const int i = d_ilist[ii];
const X_FLOAT xtmp = x(i,0);
const X_FLOAT ytmp = x(i,1);
const X_FLOAT ztmp = x(i,2);
const int jnum = d_numneigh[i];
int inside = 0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj);
j &= NEIGHMASK;
const X_FLOAT delx = xtmp - x(j,0);
const X_FLOAT dely = ytmp - x(j,1);
const X_FLOAT delz = ztmp - x(j,2);
const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cutmax*cutmax) {
d_neighbors_short(i,inside) = j;
inside++;
}
}
d_numneigh_short(i) = inside;
}
/* ---------------------------------------------------------------------- */
template<class DeviceType> template<class DeviceType>
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
@ -273,21 +314,22 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
const int itype = type(i); const int itype = type(i);
const int itag = tag(i); const int itag = tag(i);
int j,k,jj,kk,jtag,jtype,ktype;
F_FLOAT rsq1, cutsq1, rsq2, cutsq2, rij, rik, bo_ij;
F_FLOAT fi[3], fj[3], fk[3]; F_FLOAT fi[3], fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
//const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i); //const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
// repulsive // repulsive
for (jj = 0; jj < jnum; jj++) { F_FLOAT f_x = 0.0;
j = d_neighbors(i,jj); F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
jtag = tag(j); const int jtag = tag(j);
if (itag > jtag) { if (itag > jtag) {
if ((itag+jtag) % 2 == 0) continue; if ((itag+jtag) % 2 == 0) continue;
@ -315,9 +357,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
(tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r; (tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r;
const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp; const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp;
a_f(i,0) += delx*frep; f_x += delx*frep;
a_f(i,1) += dely*frep; f_y += dely*frep;
a_f(i,2) += delz*frep; f_z += delz*frep;
a_f(j,0) -= delx*frep; a_f(j,0) -= delx*frep;
a_f(j,1) -= dely*frep; a_f(j,1) -= dely*frep;
a_f(j,2) -= delz*frep; a_f(j,2) -= delz*frep;
@ -330,35 +372,35 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
delx1 = xtmp - x(j,0); const F_FLOAT delx1 = xtmp - x(j,0);
dely1 = ytmp - x(j,1); const F_FLOAT dely1 = ytmp - x(j,1);
delz1 = ztmp - x(j,2); const F_FLOAT delz1 = ztmp - x(j,2);
rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1; const F_FLOAT rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1;
cutsq1 = paramskk(itype,jtype,jtype).cutsq; const F_FLOAT cutsq1 = paramskk(itype,jtype,jtype).cutsq;
bo_ij = 0.0; F_FLOAT bo_ij = 0.0;
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); const F_FLOAT rij = sqrt(rsq1);
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2); bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2);
} }
@ -369,16 +411,16 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij); const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij);
const F_FLOAT fatt = -0.5*bij * dfa / rij; const F_FLOAT fatt = -0.5*bij * dfa / rij;
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa;
a_f(i,0) += delx1*fatt; f_x += delx1*fatt;
a_f(i,1) += dely1*fatt; f_y += dely1*fatt;
a_f(i,2) += delz1*fatt; f_z += delz1*fatt;
a_f(j,0) -= delx1*fatt; F_FLOAT fj_x = -delx1*fatt;
a_f(j,1) -= dely1*fatt; F_FLOAT fj_y = -dely1*fatt;
a_f(j,2) -= delz1*fatt; F_FLOAT fj_z = -delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
const F_FLOAT eng = 0.5*bij * fa;
if (eflag) ev.evdwl += eng; if (eflag) ev.evdwl += eng;
if (vflag_either || eflag_atom) if (vflag_either || eflag_atom)
this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1); this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1);
@ -386,29 +428,29 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
a_f(i,0) += fi[0]; f_x += fi[0];
a_f(i,1) += fi[1]; f_y += fi[1];
a_f(i,2) += fi[2]; f_z += fi[2];
a_f(j,0) += fj[0]; fj_x += fj[0];
a_f(j,1) += fj[1]; fj_y += fj[1];
a_f(j,2) += fj[2]; fj_z += fj[2];
a_f(k,0) += fk[0]; a_f(k,0) += fk[0];
a_f(k,1) += fk[1]; a_f(k,1) += fk[1];
a_f(k,2) += fk[2]; a_f(k,2) += fk[2];
@ -420,7 +462,13 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeHalf<N
if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik); if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik);
} }
} }
a_f(j,0) += fj_x;
a_f(j,1) += fj_y;
a_f(j,2) += fj_z;
} }
a_f(i,0) += f_x;
a_f(i,1) += f_y;
a_f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -450,12 +498,15 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2; X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
//const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i); //const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
// repulsive // repulsive
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const int jtype = type(j); const int jtype = type(j);
@ -475,9 +526,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
(tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r; (tmp_fcd - tmp_fce*paramskk(itype,jtype,jtype).lam1) / r;
const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp; const F_FLOAT eng = tmp_fce * paramskk(itype,jtype,jtype).biga * tmp_exp;
f(i,0) += delx*frep; f_x += delx*frep;
f(i,1) += dely*frep; f_y += dely*frep;
f(i,2) += delz*frep; f_z += delz*frep;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -490,7 +541,7 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); jtype = type(j);
@ -506,7 +557,7 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -530,9 +581,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) += delx1*fatt; f_x += delx1*fatt;
f(i,1) += dely1*fatt; f_y += dely1*fatt;
f(i,2) += delz1*fatt; f_z += delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) ev.evdwl += 0.5*eng; if (eflag) ev.evdwl += 0.5*eng;
@ -544,7 +595,7 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -559,9 +610,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
f(i,0) += fi[0]; f_x += fi[0];
f(i,1) += fi[1]; f_y += fi[1];
f(i,2) += fi[2]; f_z += fi[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrij[3], delrik[3]; F_FLOAT delrij[3], delrik[3];
@ -571,6 +622,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullA<
} }
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -599,12 +653,16 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
F_FLOAT fj[3], fk[3]; F_FLOAT fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2; X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
if (j >= nlocal) continue; if (j >= nlocal) continue;
jtype = type(j); jtype = type(j);
@ -619,10 +677,10 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); rij = sqrt(rsq1);
j_jnum = d_numneigh[j]; j_jnum = d_numneigh_short[j];
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -648,9 +706,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) -= delx1*fatt; f_x -= delx1*fatt;
f(i,1) -= dely1*fatt; f_y -= dely1*fatt;
f(i,2) -= delz1*fatt; f_z -= delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -662,7 +720,7 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -677,9 +735,9 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
rik = sqrt(rsq2); rik = sqrt(rsq2);
ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fj,fk); rik,delx2,dely2,delz2,fj,fk);
f(i,0) += fj[0]; f_x += fj[0];
f(i,1) += fj[1]; f_y += fj[1];
f(i,2) += fj[2]; f_z += fj[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrji[3], delrjk[3]; F_FLOAT delrji[3], delrjk[3];
@ -692,11 +750,14 @@ void PairTersoffMODKokkos<DeviceType>::operator()(TagPairTersoffMODComputeFullB<
const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij); const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij);
ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2, ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2,
rij,delx1,dely1,delz1,fk); rij,delx1,dely1,delz1,fk);
f(i,0) += fk[0]; f_x += fk[0];
f(i,1) += fk[1]; f_y += fk[1];
f(i,2) += fk[2]; f_z += fk[2];
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>

View File

@ -39,6 +39,8 @@ struct TagPairTersoffMODComputeFullA{};
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
struct TagPairTersoffMODComputeFullB{}; struct TagPairTersoffMODComputeFullB{};
struct TagPairTersoffMODComputeShortNeigh{};
template<class DeviceType> template<class DeviceType>
class PairTersoffMODKokkos : public PairTersoffMOD { class PairTersoffMODKokkos : public PairTersoffMOD {
public: public:
@ -77,6 +79,9 @@ class PairTersoffMODKokkos : public PairTersoffMOD {
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffMODComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const; void operator()(TagPairTersoffMODComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const;
KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffMODComputeShortNeigh, const int&) const;
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const; double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const;
@ -184,6 +189,7 @@ class PairTersoffMODKokkos : public PairTersoffMOD {
// hardwired to space for 15 atom types // hardwired to space for 15 atom types
//params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1]; //params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1];
int inum;
typename AT::t_x_array_randomread x; typename AT::t_x_array_randomread x;
typename AT::t_f_array f; typename AT::t_f_array f;
typename AT::t_int_1d_randomread type; typename AT::t_int_1d_randomread type;
@ -203,10 +209,12 @@ class PairTersoffMODKokkos : public PairTersoffMOD {
typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh; typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh;
//NeighListKokkos<DeviceType> k_list; //NeighListKokkos<DeviceType> k_list;
class AtomKokkos *atomKK;
int neighflag,newton_pair; int neighflag,newton_pair;
int nlocal,nall,eflag,vflag; int nlocal,nall,eflag,vflag;
Kokkos::View<int**,DeviceType> d_neighbors_short;
Kokkos::View<int*,DeviceType> d_numneigh_short;
friend void pair_virial_fdotr_compute<PairTersoffMODKokkos>(PairTersoffMODKokkos*); friend void pair_virial_fdotr_compute<PairTersoffMODKokkos>(PairTersoffMODKokkos*);
}; };

View File

@ -205,7 +205,7 @@ void PairTersoffZBLKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
nall = atom->nlocal + atom->nghost; nall = atom->nlocal + atom->nghost;
newton_pair = force->newton_pair; newton_pair = force->newton_pair;
const int inum = list->inum; inum = list->inum;
const int ignum = inum + list->gnum; const int ignum = inum + list->gnum;
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list); NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
d_numneigh = k_list->d_numneigh; d_numneigh = k_list->d_numneigh;
@ -218,6 +218,18 @@ void PairTersoffZBLKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
EV_FLOAT ev; EV_FLOAT ev;
EV_FLOAT ev_all; EV_FLOAT ev_all;
// build short neighbor list
int max_neighs = d_neighbors.dimension_1();
if ((d_neighbors_short.dimension_1() != max_neighs) ||
(d_neighbors_short.dimension_0() != ignum)) {
d_neighbors_short = Kokkos::View<int**,DeviceType>("Tersoff::neighbors_short",ignum,max_neighs);
}
if (d_numneigh_short.dimension_0()!=ignum)
d_numneigh_short = Kokkos::View<int*,DeviceType>("Tersoff::numneighs_short",ignum);
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagPairTersoffZBLComputeShortNeigh>(0,neighflag==FULL?ignum:inum), *this);
if (neighflag == HALF) { if (neighflag == HALF) {
if (evflag) if (evflag)
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffZBLComputeHalf<HALF,1> >(0,inum),*this,ev); Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairTersoffZBLComputeHalf<HALF,1> >(0,inum),*this,ev);
@ -271,6 +283,35 @@ void PairTersoffZBLKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
template<class DeviceType>
KOKKOS_INLINE_FUNCTION
void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeShortNeigh, const int& ii) const {
const int i = d_ilist[ii];
const X_FLOAT xtmp = x(i,0);
const X_FLOAT ytmp = x(i,1);
const X_FLOAT ztmp = x(i,2);
const int jnum = d_numneigh[i];
int inside = 0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors(i,jj);
j &= NEIGHMASK;
const X_FLOAT delx = xtmp - x(j,0);
const X_FLOAT dely = ytmp - x(j,1);
const X_FLOAT delz = ztmp - x(j,2);
const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
if (rsq < cutmax*cutmax) {
d_neighbors_short(i,inside) = j;
inside++;
}
}
d_numneigh_short(i) = inside;
}
/* ---------------------------------------------------------------------- */
template<class DeviceType> template<class DeviceType>
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
@ -287,21 +328,22 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
const int itype = type(i); const int itype = type(i);
const int itag = tag(i); const int itag = tag(i);
int j,k,jj,kk,jtag,jtype,ktype;
F_FLOAT rsq1, cutsq1, rsq2, cutsq2, rij, rik, bo_ij;
F_FLOAT fi[3], fj[3], fk[3]; F_FLOAT fi[3], fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
//const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i); //const AtomNeighborsConst d_neighbors_i = k_list.get_neighbors_const(i);
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
// repulsive // repulsive
for (jj = 0; jj < jnum; jj++) { F_FLOAT f_x = 0.0;
j = d_neighbors(i,jj); F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (int jj = 0; jj < jnum; jj++) {
int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
jtag = tag(j); const int jtag = tag(j);
if (itag > jtag) { if (itag > jtag) {
if ((itag+jtag) % 2 == 0) continue; if ((itag+jtag) % 2 == 0) continue;
@ -359,9 +401,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
eng = (1.0 - fermi_k(itype,jtype,jtype,r)) * eng_z + eng = (1.0 - fermi_k(itype,jtype,jtype,r)) * eng_z +
fermi_k(itype,jtype,jtype,r) * eng_t; fermi_k(itype,jtype,jtype,r) * eng_t;
a_f(i,0) += delx*frep; f_x += delx*frep;
a_f(i,1) += dely*frep; f_y += dely*frep;
a_f(i,2) += delz*frep; f_z += delz*frep;
a_f(j,0) -= delx*frep; a_f(j,0) -= delx*frep;
a_f(j,1) -= dely*frep; a_f(j,1) -= dely*frep;
a_f(j,2) -= delz*frep; a_f(j,2) -= delz*frep;
@ -374,35 +416,35 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (int jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); int j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); const int jtype = type(j);
delx1 = xtmp - x(j,0); const F_FLOAT delx1 = xtmp - x(j,0);
dely1 = ytmp - x(j,1); const F_FLOAT dely1 = ytmp - x(j,1);
delz1 = ztmp - x(j,2); const F_FLOAT delz1 = ztmp - x(j,2);
rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1; const F_FLOAT rsq1 = delx1*delx1 + dely1*dely1 + delz1*delz1;
cutsq1 = paramskk(itype,jtype,jtype).cutsq; const F_FLOAT cutsq1 = paramskk(itype,jtype,jtype).cutsq;
bo_ij = 0.0; F_FLOAT bo_ij = 0.0;
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); const F_FLOAT rij = sqrt(rsq1);
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2); bo_ij += bondorder(itype,jtype,ktype,rij,delx1,dely1,delz1,rik,delx2,dely2,delz2);
} }
@ -413,16 +455,16 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij); const F_FLOAT bij = ters_bij_k(itype,jtype,jtype,bo_ij);
const F_FLOAT fatt = -0.5*bij * dfa / rij; const F_FLOAT fatt = -0.5*bij * dfa / rij;
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa;
a_f(i,0) += delx1*fatt; f_x += delx1*fatt;
a_f(i,1) += dely1*fatt; f_y += dely1*fatt;
a_f(i,2) += delz1*fatt; f_z += delz1*fatt;
a_f(j,0) -= delx1*fatt; F_FLOAT fj_x = -delx1*fatt;
a_f(j,1) -= dely1*fatt; F_FLOAT fj_y = -dely1*fatt;
a_f(j,2) -= delz1*fatt; F_FLOAT fj_z = -delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
const F_FLOAT eng = 0.5*bij * fa;
if (eflag) ev.evdwl += eng; if (eflag) ev.evdwl += eng;
if (vflag_either || eflag_atom) if (vflag_either || eflag_atom)
this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1); this->template ev_tally<NEIGHFLAG>(ev,i,j,eng,fatt,delx1,dely1,delz1);
@ -430,29 +472,29 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < jnum; kk++) { for (int kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); int k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); const int ktype = type(k);
delx2 = xtmp - x(k,0); const F_FLOAT delx2 = xtmp - x(k,0);
dely2 = ytmp - x(k,1); const F_FLOAT dely2 = ytmp - x(k,1);
delz2 = ztmp - x(k,2); const F_FLOAT delz2 = ztmp - x(k,2);
rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2; const F_FLOAT rsq2 = delx2*delx2 + dely2*dely2 + delz2*delz2;
cutsq2 = paramskk(itype,jtype,ktype).cutsq; const F_FLOAT cutsq2 = paramskk(itype,jtype,ktype).cutsq;
if (rsq2 > cutsq2) continue; if (rsq2 > cutsq2) continue;
rik = sqrt(rsq2); const F_FLOAT rik = sqrt(rsq2);
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
a_f(i,0) += fi[0]; f_x += fi[0];
a_f(i,1) += fi[1]; f_y += fi[1];
a_f(i,2) += fi[2]; f_z += fi[2];
a_f(j,0) += fj[0]; fj_x += fj[0];
a_f(j,1) += fj[1]; fj_y += fj[1];
a_f(j,2) += fj[2]; fj_z += fj[2];
a_f(k,0) += fk[0]; a_f(k,0) += fk[0];
a_f(k,1) += fk[1]; a_f(k,1) += fk[1];
a_f(k,2) += fk[2]; a_f(k,2) += fk[2];
@ -464,7 +506,13 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeHalf<N
if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik); if (vflag_either) this->template v_tally3<NEIGHFLAG>(ev,i,j,k,fj,fk,delrij,delrik);
} }
} }
a_f(j,0) += fj_x;
a_f(j,1) += fj_y;
a_f(j,2) += fj_z;
} }
a_f(i,0) += f_x;
a_f(i,1) += f_y;
a_f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -498,8 +546,11 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
// repulsive // repulsive
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
const int jtype = type(j); const int jtype = type(j);
@ -549,9 +600,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
eng = (1.0 - fermi_k(itype,jtype,jtype,r)) * eng_z + eng = (1.0 - fermi_k(itype,jtype,jtype,r)) * eng_z +
fermi_k(itype,jtype,jtype,r) * eng_t; fermi_k(itype,jtype,jtype,r) * eng_t;
f(i,0) += delx*frep; f_x += delx*frep;
f(i,1) += dely*frep; f_y += dely*frep;
f(i,2) += delz*frep; f_z += delz*frep;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -564,7 +615,7 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
jtype = type(j); jtype = type(j);
@ -580,7 +631,7 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -604,9 +655,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(itype,jtype,jtype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) += delx1*fatt; f_x += delx1*fatt;
f(i,1) += dely1*fatt; f_y += dely1*fatt;
f(i,2) += delz1*fatt; f_z += delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) ev.evdwl += 0.5*eng; if (eflag) ev.evdwl += 0.5*eng;
@ -618,7 +669,7 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
for (kk = 0; kk < jnum; kk++) { for (kk = 0; kk < jnum; kk++) {
if (jj == kk) continue; if (jj == kk) continue;
k = d_neighbors(i,kk); k = d_neighbors_short(i,kk);
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -633,9 +684,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthb(itype,jtype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fi,fj,fk); rik,delx2,dely2,delz2,fi,fj,fk);
f(i,0) += fi[0]; f_x += fi[0];
f(i,1) += fi[1]; f_y += fi[1];
f(i,2) += fi[2]; f_z += fi[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrij[3], delrik[3]; F_FLOAT delrij[3], delrik[3];
@ -645,6 +696,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullA<
} }
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>
@ -673,12 +727,16 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
F_FLOAT fj[3], fk[3]; F_FLOAT fj[3], fk[3];
X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2; X_FLOAT delx1, dely1, delz1, delx2, dely2, delz2;
const int jnum = d_numneigh[i]; const int jnum = d_numneigh_short[i];
F_FLOAT f_x = 0.0;
F_FLOAT f_y = 0.0;
F_FLOAT f_z = 0.0;
// attractive: bond order // attractive: bond order
for (jj = 0; jj < jnum; jj++) { for (jj = 0; jj < jnum; jj++) {
j = d_neighbors(i,jj); j = d_neighbors_short(i,jj);
j &= NEIGHMASK; j &= NEIGHMASK;
if (j >= nlocal) continue; if (j >= nlocal) continue;
jtype = type(j); jtype = type(j);
@ -693,10 +751,10 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
if (rsq1 > cutsq1) continue; if (rsq1 > cutsq1) continue;
rij = sqrt(rsq1); rij = sqrt(rsq1);
j_jnum = d_numneigh[j]; j_jnum = d_numneigh_short[j];
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -722,9 +780,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij); const F_FLOAT prefactor = 0.5*fa * ters_dbij(jtype,itype,itype,bo_ij);
const F_FLOAT eng = 0.5*bij * fa; const F_FLOAT eng = 0.5*bij * fa;
f(i,0) -= delx1*fatt; f_x -= delx1*fatt;
f(i,1) -= dely1*fatt; f_y -= dely1*fatt;
f(i,2) -= delz1*fatt; f_z -= delz1*fatt;
if (EVFLAG) { if (EVFLAG) {
if (eflag) if (eflag)
@ -736,7 +794,7 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
// attractive: three-body force // attractive: three-body force
for (kk = 0; kk < j_jnum; kk++) { for (kk = 0; kk < j_jnum; kk++) {
k = d_neighbors(j,kk); k = d_neighbors_short(j,kk);
if (k == i) continue; if (k == i) continue;
k &= NEIGHMASK; k &= NEIGHMASK;
ktype = type(k); ktype = type(k);
@ -751,9 +809,9 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
rik = sqrt(rsq2); rik = sqrt(rsq2);
ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1, ters_dthbj(jtype,itype,ktype,prefactor,rij,delx1,dely1,delz1,
rik,delx2,dely2,delz2,fj,fk); rik,delx2,dely2,delz2,fj,fk);
f(i,0) += fj[0]; f_x += fj[0];
f(i,1) += fj[1]; f_y += fj[1];
f(i,2) += fj[2]; f_z += fj[2];
if (vflag_atom) { if (vflag_atom) {
F_FLOAT delrji[3], delrjk[3]; F_FLOAT delrji[3], delrjk[3];
@ -766,11 +824,14 @@ void PairTersoffZBLKokkos<DeviceType>::operator()(TagPairTersoffZBLComputeFullB<
const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij); const F_FLOAT prefactor_jk = 0.5*fa_jk * ters_dbij(jtype,ktype,itype,bo_ij);
ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2, ters_dthbk(jtype,ktype,itype,prefactor_jk,rik,delx2,dely2,delz2,
rij,delx1,dely1,delz1,fk); rij,delx1,dely1,delz1,fk);
f(i,0) += fk[0]; f_x += fk[0];
f(i,1) += fk[1]; f_y += fk[1];
f(i,2) += fk[2]; f_z += fk[2];
} }
} }
f(i,0) += f_x;
f(i,1) += f_y;
f(i,2) += f_z;
} }
template<class DeviceType> template<class DeviceType>

View File

@ -39,6 +39,8 @@ struct TagPairTersoffZBLComputeFullA{};
template<int NEIGHFLAG, int EVFLAG> template<int NEIGHFLAG, int EVFLAG>
struct TagPairTersoffZBLComputeFullB{}; struct TagPairTersoffZBLComputeFullB{};
struct TagPairTersoffZBLComputeShortNeigh{};
template<class DeviceType> template<class DeviceType>
class PairTersoffZBLKokkos : public PairTersoffZBL { class PairTersoffZBLKokkos : public PairTersoffZBL {
public: public:
@ -77,6 +79,8 @@ class PairTersoffZBLKokkos : public PairTersoffZBL {
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffZBLComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const; void operator()(TagPairTersoffZBLComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const;
KOKKOS_INLINE_FUNCTION
void operator()(TagPairTersoffZBLComputeShortNeigh, const int&) const;
KOKKOS_INLINE_FUNCTION KOKKOS_INLINE_FUNCTION
double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const; double ters_fc_k(const int &i, const int &j, const int &k, const F_FLOAT &r) const;
@ -190,6 +194,7 @@ class PairTersoffZBLKokkos : public PairTersoffZBL {
// hardwired to space for 15 atom types // hardwired to space for 15 atom types
//params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1]; //params_ters m_params[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1];
int inum;
typename AT::t_x_array_randomread x; typename AT::t_x_array_randomread x;
typename AT::t_f_array f; typename AT::t_f_array f;
typename AT::t_int_1d_randomread type; typename AT::t_int_1d_randomread type;
@ -209,10 +214,12 @@ class PairTersoffZBLKokkos : public PairTersoffZBL {
typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh; typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh;
//NeighListKokkos<DeviceType> k_list; //NeighListKokkos<DeviceType> k_list;
class AtomKokkos *atomKK;
int neighflag,newton_pair; int neighflag,newton_pair;
int nlocal,nall,eflag,vflag; int nlocal,nall,eflag,vflag;
Kokkos::View<int**,DeviceType> d_neighbors_short;
Kokkos::View<int*,DeviceType> d_numneigh_short;
// ZBL // ZBL
F_FLOAT global_a_0; // Bohr radius for Coulomb repulsion F_FLOAT global_a_0; // Bohr radius for Coulomb repulsion
F_FLOAT global_epsilon_0; // permittivity of vacuum for Coulomb repulsion F_FLOAT global_epsilon_0; // permittivity of vacuum for Coulomb repulsion

View File

@ -176,7 +176,8 @@ void PairMorse::settings(int narg, char **arg)
void PairMorse::coeff(int narg, char **arg) void PairMorse::coeff(int narg, char **arg)
{ {
if (narg < 5 || narg > 6) error->all(FLERR,"Incorrect args for pair coefficients"); if (narg < 5 || narg > 6)
error->all(FLERR,"Incorrect args for pair coefficients");
if (!allocated) allocate(); if (!allocated) allocate();
int ilo,ihi,jlo,jhi; int ilo,ihi,jlo,jhi;