new 17Nov16 patch and stable

compute pressure doc change
Merge pull request #278 from akohlmey/pair-agni
2016-11-17 16:51:35 -07:00 · 2016-11-17 16:11:30 -07:00 · 2016-11-17 09:04:31 -07:00 · 2016-11-17 09:02:27 -07:00 · 2016-11-17 09:01:46 -07:00 · 2016-11-17 09:01:37 -07:00
179 changed files with 7689 additions and 321 deletions
--- a/doc/src/Eqs/pair_agni.jpg
+++ b/doc/src/Eqs/pair_agni.jpg
--- a/doc/src/Eqs/pair_tersoff_mod_c.jpg
+++ b/doc/src/Eqs/pair_tersoff_mod_c.jpg
--- a/doc/src/Eqs/pair_tersoff_mod_c.tex
+++ b/doc/src/Eqs/pair_tersoff_mod_c.tex
@ -0,0 +1,10 @@
+\documentclass[12pt]{article}
+\pagestyle{empty}
+
+\begin{document}
+
+\begin{eqnarray*}
+  V_{ij}  & =  & f_C(r_{ij}) \left[ f_R(r_{ij}) + b_{ij} f_A(r_{ij}) + c_0 \right]
+\end{eqnarray*}
+
+\end{document}
--- a/doc/src/Eqs/pressure.jpg
+++ b/doc/src/Eqs/pressure.jpg
--- a/doc/src/Eqs/pressure.tex
+++ b/doc/src/Eqs/pressure.tex
@ -3,7 +3,7 @@
 \begin{document}

 $$
-   P = \frac{N k_B T}{V} + \frac{\sum_{i}^{N} r_i \bullet f_i}{dV}
+   P = \frac{N k_B T}{V} + \frac{\sum_{i}^{N'} r_i \bullet f_i}{dV}
 $$

 \end{document}
--- a/doc/src/Eqs/pressure_tensor.jpg
+++ b/doc/src/Eqs/pressure_tensor.jpg
--- a/doc/src/Eqs/pressure_tensor.tex
+++ b/doc/src/Eqs/pressure_tensor.tex
@ -4,7 +4,7 @@

 $$
   P_{IJ} = \frac{\sum_{k}^{N} m_k v_{k_I} v_{k_J}}{V} + 
-   \frac{\sum_{k}^{N} r_{k_I} f_{k_J}}{V}
+   \frac{\sum_{k}^{N'} r_{k_I} f_{k_J}}{V}
 $$

 \end{document}
--- a/doc/src/Manual.txt
+++ b/doc/src/Manual.txt
@ -1,7 +1,7 @@
 <!-- HTML_ONLY -->
 <HEAD>
 <TITLE>LAMMPS Users Manual</TITLE>
-<META NAME="docnumber" CONTENT="27 Oct 2016 version">
+<META NAME="docnumber" CONTENT="17 Nov 2016 version">
 <META NAME="author" CONTENT="http://lammps.sandia.gov - Sandia National Laboratories">
 <META NAME="copyright" CONTENT="Copyright (2003) Sandia Corporation.  This software and manual is distributed under the GNU General Public License.">
 </HEAD>
@ -21,7 +21,7 @@
 <H1></H1>

 LAMMPS Documentation :c,h3
-27 Oct 2016 version :c,h4
+17 Nov 2016 version :c,h4

 Version info: :h4

--- a/doc/src/Section_commands.txt
+++ b/doc/src/Section_commands.txt
@ -886,6 +886,8 @@ KOKKOS, o = USER-OMP, t = OPT.
 "body"_pair_body.html,
 "bop"_pair_bop.html,
 "born (go)"_pair_born.html,
+"born/coul/dsf"_pair_born.html,
+"born/coul/dsf/cs"_pair_born.html,
 "born/coul/long (go)"_pair_born.html,
 "born/coul/long/cs"_pair_born.html,
 "born/coul/msm (o)"_pair_born.html,
@ -979,11 +981,12 @@ KOKKOS, o = USER-OMP, t = OPT.
 "table (gko)"_pair_table.html,
 "tersoff (gkio)"_pair_tersoff.html,
 "tersoff/mod (gko)"_pair_tersoff_mod.html,
+"tersoff/mod/c (o)"_pair_tersoff_mod.html,
 "tersoff/zbl (gko)"_pair_tersoff_zbl.html,
 "tip4p/cut (o)"_pair_coul.html,
 "tip4p/long (o)"_pair_coul.html,
 "tri/lj"_pair_tri_lj.html,
-"vashishta (o)"_pair_vashishta.html,
+"vashishta (ko)"_pair_vashishta.html,
 "vashishta/table (o)"_pair_vashishta.html,
 "yukawa (go)"_pair_yukawa.html,
 "yukawa/colloid (go)"_pair_yukawa_colloid.html,
@ -993,6 +996,7 @@ These are additional pair styles in USER packages, which can be used
 if "LAMMPS is built with the appropriate
 package"_Section_start.html#start_3.

+"agni (o)"_pair_agni.html,
 "awpmd/cut"_pair_awpmd.html,
 "buck/mdf"_pair_mdf.html,
 "coul/cut/soft (o)"_pair_lj_soft.html,
--- a/doc/src/Section_intro.txt
+++ b/doc/src/Section_intro.txt
@ -366,11 +366,11 @@ complementary modeling tasks.
 "DL_POLY"_dlpoly
 "Tinker"_tinker :ul

-:link(charmm,http://www.scripps.edu/brooks)
-:link(amber,http://amber.scripps.edu)
+:link(charmm,http://www.charmm.org)
+:link(amber,http://ambermd.org)
 :link(namd,http://www.ks.uiuc.edu/Research/namd/)
 :link(nwchem,http://www.emsl.pnl.gov/docs/nwchem/nwchem.html)
-:link(dlpoly,http://www.cse.clrc.ac.uk/msi/software/DL_POLY)
+:link(dlpoly,http://www.ccp5.ac.uk/DL_POLY_CLASSIC)
 :link(tinker,http://dasher.wustl.edu/tinker)

 CHARMM, AMBER, NAMD, NWCHEM, and Tinker are designed primarily for
--- a/doc/src/compute_pressure.txt
+++ b/doc/src/compute_pressure.txt
@ -37,12 +37,18 @@ The pressure is computed by the formula

 where N is the number of atoms in the system (see discussion of DOF
 below), Kb is the Boltzmann constant, T is the temperature, d is the
-dimensionality of the system (2 or 3 for 2d/3d), V is the system
-volume (or area in 2d), and the second term is the virial, computed
-within LAMMPS for all pairwise as well as 2-body, 3-body, and 4-body,
-and long-range interactions.  "Fixes"_fix.html that impose constraints
-(e.g. the "fix shake"_fix_shake.html command) also contribute to the
-virial term.
+dimensionality of the system (2 or 3 for 2d/3d), and V is the system
+volume (or area in 2d).  The second term is the virial, equal to
+-dU/dV, computed for all pairwise as well as 2-body, 3-body, 4-body,
+manybody, and long-range interactions, where r_i and f_i are the
+position and force vector of atom i, and the black dot indicates a dot
+product.  When periodic boundary conditions are used, N' necessarily
+includes periodic image (ghost) atoms outside the central box, and the
+position and force vectors of ghost atoms are thus included in the
+summation.  When periodic boundary conditions are not used, N' = N =
+the number of atoms in the system.  "Fixes"_fix.html that impose
+constraints (e.g. the "fix shake"_fix_shake.html command) also
+contribute to the virial term.

 A symmetric pressure tensor, stored as a 6-element vector, is also
 calculated by this compute.  The 6 components of the vector are
@ -62,8 +68,9 @@ compute temperature or ke and/or the virial.  The {virial} keyword
 means include all terms except the kinetic energy {ke}.

 Details of how LAMMPS computes the virial efficiently for the entire
-system, including the effects of periodic boundary conditions is
-discussed in "(Thompson)"_#Thompson.
+system, including for manybody potentials and accounting for the
+effects of periodic boundary conditions are discussed in
+"(Thompson)"_#Thompson.

 The temperature and kinetic energy tensor is not calculated by this
 compute, but rather by the temperature compute specified with the
--- a/doc/src/info.txt
+++ b/doc/src/info.txt
@ -14,7 +14,7 @@ info args :pre

 args = one or more of the following keywords: {out}, {all}, {system}, {communication}, {computes}, {dumps}, {fixes}, {groups}, {regions}, {variables}, {styles}, {time}, or {configuration}
     {out} values = {screen}, {log}, {append} filename, {overwrite} filename
-     {styles} values = {all}, {angle}, {atom}, {bond}, {compute}, {command}, {dump}, {dihedral}, {fix}, {improper}, {integrate}, {kspace}, {minimize}, {region} :ul
+     {styles} values = {all}, {angle}, {atom}, {bond}, {compute}, {command}, {dump}, {dihedral}, {fix}, {improper}, {integrate}, {kspace}, {minimize}, {pair}, {region} :ul

 [Examples:]

@ -70,8 +70,9 @@ The {variables} category prints a list of all currently defined
 variables, their names, styles, definition and last computed value, if
 available.

-The {styles} category prints the list of styles available in LAMMPS. It
-supports one of the following options to control what is printed out:
+The {styles} category prints the list of styles available in the
+current LAMMPS binary. It supports one of the following options
+to control which category of styles is printed out:

 all
 angle
@ -86,6 +87,7 @@ improper
 integrate
 kspace
 minimize
+pair
 region :ul

 The {time} category prints the accumulated CPU and wall time for the
--- a/doc/src/lammps.book
+++ b/doc/src/lammps.book
@ -59,6 +59,7 @@ dump_h5md.html
 dump_image.html
 dump_modify.html
 dump_molfile.html
+dump_nc.html
 echo.html
 fix.html
 fix_modify.html
@ -152,6 +153,7 @@ fix_colvars.html
 fix_controller.html
 fix_deform.html
 fix_deposit.html
+fix_dpd_energy.html
 fix_drag.html
 fix_drude.html
 fix_drude_transform.html
@ -272,6 +274,7 @@ fix_viscosity.html
 fix_viscous.html
 fix_wall.html
 fix_wall_gran.html
+fix_wall_gran_region.html
 fix_wall_piston.html
 fix_wall_reflect.html
 fix_wall_region.html
@ -390,6 +393,7 @@ compute_voronoi_atom.html
 compute_xrd.html

 pair_adp.html
+pair_agni.html
 pair_airebo.html
 pair_awpmd.html
 pair_beck.html
@ -622,3 +626,4 @@ USER/atc/man_unfix_flux.html
 USER/atc/man_unfix_nodes.html
 USER/atc/man_write_atom_weights.html
 USER/atc/man_write_restart.html
+
--- a/doc/src/pair_agni.txt
+++ b/doc/src/pair_agni.txt
@ -0,0 +1,128 @@
+"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
+
+:link(lws,http://lammps.sandia.gov)
+:link(ld,Manual.html)
+:link(lc,Section_commands.html#comm)
+
+:line
+
+pair_style agni command :h3
+pair_style agni/omp command :h3
+
+[Syntax:]
+
+pair_style agni :pre
+
+[Examples:]
+pair_style      agni
+pair_coeff      * * Al.agni Al
+
+[Description:]
+
+Style {agni} style computes the manybody vectorial force components for
+an atom as
+
+:c,image(Eqs/pair_agni.jpg)
+
+{u} labels the individual components, i.e. x, y or z, and {V} is the
+corresponding atomic fingerprint. {d} is the Euclidean distance between
+any two atomic fingerprints. A total of N_t reference atomic
+environments are considered to construct the force field file. {alpha_t}
+and {l} are the weight coefficients and length scale parameter of the
+non-linear regression model.
+
+The method implements the recently proposed machine learning access to
+atomic forces as discussed extensively in the following publications -
+"(Botu1)"_#Botu2015adaptive and "(Botu2)"_#Botu2015learning. The premise
+of the method is to map the atomic enviornment numerically into a
+fingerprint, and use machine learning methods to create a mapping to the
+vectorial atomic forces.
+
+Only a single pair_coeff command is used with the {agni} style which
+specifies an AGNI potential file containing the parameters of the
+force field for the needed elements. These are mapped to LAMMPS atom 
+types by specifying N additional arguments after the filename in the 
+pair_coeff command, where N is the number of LAMMPS atom types:
+
+filename
+N element names = mapping of AGNI elements to atom types :ul
+
+See the "pair_coeff"_pair_coeff.html doc page for alternate ways
+to specify the path for the force field file.
+
+An AGNI force field is fully specified by the filename which contains the
+parameters of the force field, i.e., the reference training environments
+used to construct the machine learning force field. Example force field 
+and input files are provided in the examples/USER/misc/agni directory. 
+
+:line
+
+Styles with {omp} suffix is 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 5"_Section_accelerate.html
+of the manual.  The accelerated style takes the same arguments and
+should produce the same results, except for round-off and precision
+issues.
+
+The accelerated style is part of the USER-OMP.  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 style 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 5"_Section_accelerate.html of the manual for
+more instructions on how to use the accelerated styles effectively.
+
+:line
+
+[Mixing, shift, table, tail correction, restart, rRESPA info]:
+
+This pair style does not support the "pair_modify"_pair_modify.html
+shift, table, and tail options.
+
+This pair style does not write its information to "binary restart
+files"_restart.html, since it is stored in potential files.  Thus, you
+need to re-specify the pair_style and pair_coeff commands in an input
+script that reads a restart file.
+
+This pair style can only be used via the {pair} keyword of the
+"run_style respa"_run_style.html command.  It does not support the
+{inner}, {middle}, {outer} keywords.
+
+:line
+
+[Restrictions:]
+
+Currently, only elemental systems are implemented. Also, the method only
+provides access to the forces and not energies or stresses. However, one
+can access the energy via thermodynamic integration of the forces as
+discussed in "(Botu3)"_#Botu2016construct.  This pair style is part
+of the USER-MISC package. It is only enabled if LAMMPS was built with
+that package. See the "Making LAMMPS"_Section_start.html#start_3 section
+for more info.
+
+The AGNI force field files provided with LAMMPS (see the
+potentials directory) are parameterized for metal "units"_units.html.
+You can use the AGNI potential with any LAMMPS units, but you would need
+to create your own AGNI potential file with coefficients listed in the
+appropriate units if your simulation doesn't use "metal" units.
+
+[Related commands:]
+
+"pair_coeff"_pair_coeff.html
+
+[Default:] none
+
+:line
+
+:link(Botu2015adaptive)
+[(Botu1)] V. Botu and R. Ramprasad, Int. J. Quant. Chem., 115(16), 1074 (2015).
+
+:link(Botu2015learning)
+[(Botu2)] V. Botu and R. Ramprasad, Phys. Rev. B, 92(9), 094306 (2015).
+
+:link(Botu2016construct)
+[(Botu3)] V. Botu, R. Batra, J. Chapman and R. Ramprasad, https://arxiv.org/abs/1610.02098 (2016).
--- a/doc/src/pair_born.txt
+++ b/doc/src/pair_born.txt
@ -19,6 +19,8 @@ pair_style born/coul/msm/omp command :h3
 pair_style born/coul/wolf command :h3
 pair_style born/coul/wolf/gpu command :h3
 pair_style born/coul/wolf/omp command :h3
+pair_style born/coul/dsf command :h3
+pair_style born/coul/dsf/cs command :h3

 [Syntax:]

@ -37,7 +39,11 @@ args = list of arguments for a particular style :ul
  {born/coul/wolf} args = alpha cutoff (cutoff2)
    alpha = damping parameter (inverse distance units)
    cutoff = global cutoff for non-Coulombic (and Coulombic if only 1 arg) (distance units)
-    cutoff2 = global cutoff for Coulombic (optional) (distance units) :pre
+    cutoff2 = global cutoff for Coulombic (optional) (distance units)
+  {born/coul/dsf} or {born/coul/dsf/cs} args = alpha cutoff (cutoff2)
+    alpha = damping parameter (inverse distance units)
+    cutoff = global cutoff for non-Coulombic (and Coulombic if only 1 arg) (distance units)
+    cutoff2 = global cutoff for Coulombic (distance units) :pre

 [Examples:]

@ -62,6 +68,10 @@ pair_style born/coul/wolf 0.25 10.0 9.0
 pair_coeff * * 6.08 0.317 2.340 24.18 11.51
 pair_coeff 1 1 6.08 0.317 2.340 24.18 11.51 :pre

+pair_style born/coul/dsf 0.1 10.0 12.0
+pair_coeff * *   0.0 1.00 0.00 0.00 0.00
+pair_coeff 1 1 480.0 0.25 0.00 1.05 0.50 :pre
+
 [Description:]

 The {born} style computes the Born-Mayer-Huggins or Tosi/Fumi
@ -90,10 +100,14 @@ term.
 The {born/coul/wolf} style adds a Coulombic term as described for the
 Wolf potential in the "coul/wolf"_pair_coul.html pair style.

+The {born/coul/dsf} style computes the Coulomb contribution with the
+damped shifted force model as in the "coul/dsf"_pair_coul.html style.
+
 Style {born/coul/long/cs} is identical to {born/coul/long} except that
 a term is added for the "core/shell model"_Section_howto.html#howto_25
 to allow charges on core and shell particles to be separated by r =
-0.0.
+0.0. The same correction is introduced for {born/coul/dsf/cs} style
+which is identical to {born/coul/dsf}.

 Note that these potentials are related to the "Buckingham
 potential"_pair_buck.html.
@ -116,9 +130,10 @@ The second coefficient, rho, must be greater than zero.
 The last coefficient is optional.  If not specified, the global A,C,D
 cutoff specified in the pair_style command is used.

-For {born/coul/long} and {born/coul/wolf} no Coulombic cutoff can be
-specified for an individual I,J type pair.  All type pairs use the
-same global Coulombic cutoff specified in the pair_style command.
+For {born/coul/long}, {born/coul/wolf} and {born/coul/dsf} no
+Coulombic cutoff can be specified for an individual I,J type pair.
+All type pairs use the same global Coulombic cutoff specified in the
+pair_style command.

 :line

--- a/doc/src/pair_cs.txt
+++ b/doc/src/pair_cs.txt
@ -8,19 +8,24 @@

 pair_style born/coul/long/cs command :h3
 pair_style buck/coul/long/cs command :h3
+pair_style born/coul/dsf/cs command :h3

 [Syntax:]

 pair_style style args :pre

-style = {born/coul/long/cs} or {buck/coul/long/cs}
+style = {born/coul/long/cs} or {buck/coul/long/cs} or {born/coul/dsf/cs}
 args = list of arguments for a particular style :ul
  {born/coul/long/cs} args = cutoff (cutoff2)
    cutoff = global cutoff for non-Coulombic (and Coulombic if only 1 arg) (distance units)
    cutoff2 = global cutoff for Coulombic (optional) (distance units)
  {buck/coul/long/cs} args = cutoff (cutoff2)
    cutoff = global cutoff for Buckingham (and Coulombic if only 1 arg) (distance units)
-    cutoff2 = global cutoff for Coulombic (optional) (distance units) :pre
+    cutoff2 = global cutoff for Coulombic (optional) (distance units)
+  {born/coul/dsf/cs} args = alpha cutoff (cutoff2)
+    alpha = damping parameter (inverse distance units)
+    cutoff = global cutoff for non-Coulombic (and Coulombic if only 1 arg) (distance units)
+    cutoff2 = global cutoff for Coulombic (distance units) :pre

 [Examples:]

@ -32,6 +37,10 @@ pair_style buck/coul/long/cs 10.0 8.0
 pair_coeff * * 100.0 1.5 200.0
 pair_coeff 1 1 100.0 1.5 200.0 9.0 :pre

+pair_style born/coul/dsf/cs 0.1 10.0 12.0
+pair_coeff * *   0.0 1.00 0.00 0.00 0.00
+pair_coeff 1 1 480.0 0.25 0.00 1.05 0.50 :pre
+
 [Description:]

 These pair styles are designed to be used with the adiabatic
@ -39,7 +48,7 @@ core/shell model of "(Mitchell and Finchham)"_#MitchellFinchham.  See
 "Section 6.25"_Section_howto.html#howto_25 of the manual for an
 overview of the model as implemented in LAMMPS.

-These pair styles are identical to the "pair_style
+The styles with a {coul/long} term are identical to the "pair_style
 born/coul/long"_pair_born.html and "pair_style
 buck/coul/long"_pair_buck.html styles, except they correctly treat the
 special case where the distance between two charged core and shell
@ -63,6 +72,14 @@ where C is an energy-conversion constant, Qi and Qj are the charges on
 the core and shell, epsilon is the dielectric constant and r_min is the
 minimal distance.

+The pair style {born/coul/dsf/cs} is identical to the
+"pair_style born/coul/dsf"_pair_born.html style, which uses the
+the damped shifted force model as in "coul/dsf"_pair_coul.html
+to compute the Coulomb contribution. This approach does not require
+a long-range solver, thus the only correction is the addition of a
+minimal distance to avoid the possible r = 0.0 case for a
+core/shell pair.
+
 [Restrictions:]

 These pair styles are part of the CORESHELL package.  They are only
--- a/doc/src/pair_snap.txt
+++ b/doc/src/pair_snap.txt
@ -15,7 +15,7 @@ pair_style snap :pre
 [Examples:]

 pair_style snap
-pair_coeff * * snap InP.snapcoeff In P InP.snapparam In In P P :pre
+pair_coeff * * InP.snapcoeff In P InP.snapparam In In P P :pre

 [Description:]

@ -27,9 +27,9 @@ it uses bispectrum components
 to characterize the local neighborhood of each atom
 in a very general way. The mathematical definition of the
 bispectrum calculation used by SNAP is identical
-to that used of "compute sna/atom"_compute_sna_atom.html.
+to that used by "compute sna/atom"_compute_sna_atom.html.
 In SNAP, the total energy is decomposed into a sum over
-atom energies. The energy of atom {i} is
+atom energies. The energy of atom {i } is
 expressed as a weighted sum over bispectrum components.

 :c,image(Eqs/pair_snap.jpg)
@ -183,8 +183,7 @@ LAMMPS"_Section_start.html#start_3 section for more info.
 :line

 :link(Thompson2014)
-[(Thompson)] Thompson, Swiler, Trott, Foiles, Tucker, under review, preprint
-available at "arXiv:1409.3880"_http://arxiv.org/abs/1409.3880
+[(Thompson)] Thompson, Swiler, Trott, Foiles, Tucker, J Comp Phys, 285, 316 (2015).

 :link(Bartok2010)
 [(Bartok2010)] Bartok, Payne, Risi, Csanyi, Phys Rev Lett, 104, 136403 (2010).
--- a/doc/src/pair_tersoff_mod.txt
+++ b/doc/src/pair_tersoff_mod.txt
@ -7,32 +7,43 @@
 :line

 pair_style tersoff/mod command :h3
+pair_style tersoff/mod/c command :h3
 pair_style tersoff/mod/gpu command :h3
 pair_style tersoff/mod/kk command :h3
 pair_style tersoff/mod/omp command :h3
+pair_style tersoff/mod/c/omp command :h3

 [Syntax:]

 pair_style tersoff/mod :pre

+pair_style tersoff/mod/c :pre
+
 [Examples:]

 pair_style tersoff/mod
 pair_coeff * * Si.tersoff.mod Si Si :pre

+pair_style tersoff/mod/c
+pair_coeff * * Si.tersoff.modc Si Si :pre
+
 [Description:]

-The {tersoff/mod} style computes a bond-order type interatomic
-potential "(Kumagai)"_#Kumagai based on a 3-body Tersoff potential
-"(Tersoff_1)"_#Tersoff_1, "(Tersoff_2)"_#Tersoff_2 with modified
-cutoff function and angular-dependent term, giving the energy E of a
-system of atoms as
+The {tersoff/mod} and {tersoff/mod/c} styles computes a bond-order type
+interatomic potential "(Kumagai)"_#Kumagai based on a 3-body Tersoff
+potential "(Tersoff_1)"_#Tersoff_1, "(Tersoff_2)"_#Tersoff_2 with
+modified cutoff function and angular-dependent term, giving the energy
+E of a system of atoms as

 :c,image(Eqs/pair_tersoff_mod.jpg)

 where f_R is a two-body term and f_A includes three-body interactions.
 The summations in the formula are over all neighbors J and K of atom I
 within a cutoff distance = R + D.
+The {tersoff/mod/c} style differs from {tersoff/mod} only in the
+formulation of the V_ij term, where it contains an additional c0 term.
+
+:c,image(Eqs/pair_tersoff_mod_c.jpg)

 The modified cutoff function f_C proposed by "(Murty)"_#Murty and
 having a continuous second-order differential is employed. The
@ -69,10 +80,11 @@ are placeholders for atom types that will be used with other
 potentials.

 Tersoff/MOD file in the {potentials} directory of the LAMMPS
-distribution have a ".tersoff.mod" suffix.  Lines that are not blank
-or comments (starting with #) define parameters for a triplet of
-elements.  The parameters in a single entry correspond to coefficients
-in the formula above:
+distribution have a ".tersoff.mod" suffix. Potential files for the
+{tersoff/mod/c} style have the suffix ".tersoff.modc". Lines that are
+not blank or comments (starting with #) define parameters for a triplet
+of elements.  The parameters in a single entry correspond to
+coefficients in the formulae above:

 element 1 (the center atom in a 3-body interaction)
 element 2 (the atom bonded to the center atom)
@ -93,13 +105,15 @@ c1
 c2
 c3
 c4
-c5 :ul
+c5
+c0 (energy units, tersoff/mod/c only):ul

 The n, eta, lambda2, B, lambda1, and A parameters are only used for
 two-body interactions.  The beta, alpha, c1, c2, c3, c4, c5, h
 parameters are only used for three-body interactions. The R and D
-parameters are used for both two-body and three-body interactions. The
-non-annotated parameters are unitless.
+parameters are used for both two-body and three-body interactions.
+The c0 term applies to {tersoff/mod/c} only. The non-annotated
+parameters are unitless.

 The Tersoff/MOD potential file must contain entries for all the elements
 listed in the pair_coeff command.  It can also contain entries for
--- a/doc/src/pair_vashishta.txt
+++ b/doc/src/pair_vashishta.txt
@ -8,6 +8,7 @@

 pair_style vashishta command :h3
 pair_style vashishta/omp command :h3
+pair_style vashishta/kk command :h3
 pair_style vashishta/table command :h3
 pair_style vashishta/table/omp command :h3

--- a/doc/src/pairs.txt
+++ b/doc/src/pairs.txt
@ -6,6 +6,7 @@ Pair Styles :h1
   :maxdepth: 1

   pair_adp
+   pair_agni
   pair_airebo
   pair_awpmd
   pair_beck
--- a/examples/USER/misc/addtorque/in.addtorque
+++ b/examples/USER/misc/addtorque/in.addtorque
@ -0,0 +1,64 @@
+
+units           lj
+lattice         fcc 0.8442
+boundary        f f f
+
+region          ball   sphere   0.0 0.0 0.0 5.0
+region          box    block  -10 10 -10 10 -10 10
+region          half   block  -10  0 -10 10 -10 10
+
+# add molecule ids so we can use chunk/atom
+fix             0   all property/atom mol ghost yes
+
+create_box      2  box
+create_atoms    1  region ball
+
+pair_style      lj/cut 4.0
+pair_coeff      * *  1.0 1.0
+mass            *    1.0
+set             group all mol 1
+
+# label half the sphere with a different type for better visualization
+set             region half   type 2
+
+# use a dynamic group (may require a patch to fix addtorque with older versions of LAMMPS)
+group           ball   dynamic all   region ball
+
+neigh_modify    delay 2  every 1  check yes
+
+minimize        0.0  0.0    1000  10000
+reset_timestep  0
+
+velocity        all create 1.2 12351235
+
+fix             1   all  nve
+fix             2   all  wall/reflect  xlo EDGE xhi EDGE ylo EDGE yhi EDGE zlo EDGE zhi EDGE
+
+compute 1 all chunk/atom molecule
+compute 2 ball omega/chunk 1
+compute 3 ball inertia/chunk 1
+
+# compute rotational kinetic energy: 1/2 * I * omega**2
+variable r_ke_x equal c_2[1][1]*c_2[1][1]*c_3[1][1]*0.5
+variable r_ke_y equal c_2[1][2]*c_2[1][2]*c_3[1][2]*0.5
+variable r_ke_z equal c_2[1][3]*c_2[1][3]*c_3[1][3]*0.5
+
+# output moments of inertia for x,y,z and angular velocity as well as rotational kinertic energy
+
+thermo_style custom step ke pe c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z
+thermo 500
+
+# dump 1 all atom 100 dump.lammpstrj
+
+# dump 2 all movie 10 ball.mp4 type mass
+
+# equilibration w/o torque added
+run 1000 post no
+
+# start spinning the ball. rotation around z should increase and Erot_z should grow
+fix 4 ball addtorque 0.0 0.0 200.0
+run 10000 upto  post no
+
+# continue without adding more torque. rotation continues at fixed speed
+unfix 4
+run 5000
--- a/examples/USER/misc/addtorque/log.9Nov16.addtorque.g++.1
+++ b/examples/USER/misc/addtorque/log.9Nov16.addtorque.g++.1
@ -0,0 +1,202 @@
+LAMMPS (9 Nov 2016)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:90)
+  using 1 OpenMP thread(s) per MPI task
+
+units           lj
+lattice         fcc 0.8442
+Lattice spacing in x,y,z = 1.6796 1.6796 1.6796
+boundary        f f f
+
+region          ball   sphere   0.0 0.0 0.0 5.0
+region          box    block  -10 10 -10 10 -10 10
+region          half   block  -10  0 -10 10 -10 10
+
+# add molecule ids so we can use chunk/atom
+fix             0   all property/atom mol ghost yes
+
+create_box      2  box
+Created orthogonal box = (-16.796 -16.796 -16.796) to (16.796 16.796 16.796)
+  1 by 1 by 1 MPI processor grid
+create_atoms    1  region ball
+Created 2123 atoms
+
+pair_style      lj/cut 4.0
+pair_coeff      * *  1.0 1.0
+mass            *    1.0
+set             group all mol 1
+  2123 settings made for mol
+
+# label half the sphere with a different type for better visualization
+set             region half   type 2
+  1142 settings made for type
+
+# use a dynamic group (may require a patch to fix addtorque with older versions of LAMMPS)
+group           ball   dynamic all   region ball
+dynamic group ball defined
+
+neigh_modify    delay 2  every 1  check yes
+
+minimize        0.0  0.0    1000  10000
+WARNING: Resetting reneighboring criteria during minimization (../min.cpp:168)
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 0 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 4.3
+  ghost atom cutoff = 4.3
+  binsize = 2.15 -> bins = 16 16 16
+Memory usage per processor = 3.68322 Mbytes
+Step Temp E_pair E_mol TotEng Press 
+       0            0   -6.2285099            0   -6.2285099  -0.38871568 
+    1000            0   -7.3616908            0   -7.3616908 -9.1828951e-16 
+Loop time of 12.4181 on 1 procs for 1000 steps with 2123 atoms
+
+100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+Minimization stats:
+  Stopping criterion = max iterations
+  Energy initial, next-to-last, final = 
+        -6.22850993032     -7.36169083402     -7.36169083402
+  Force two-norm initial, final = 197.762 3.3539e-12
+  Force max component initial, final = 7.88704 1.52475e-13
+  Final line search alpha, max atom move = 1 1.52475e-13
+  Iterations, force evaluations = 1000 1994
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 12.282     | 12.282     | 12.282     |   0.0 | 98.91
+Neigh   | 0.06315    | 0.06315    | 0.06315    |   0.0 |  0.51
+Comm    | 0.000443   | 0.000443   | 0.000443   |   0.0 |  0.00
+Output  | 2.7e-05    | 2.7e-05    | 2.7e-05    |   0.0 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 0.07231    |            |       |  0.58
+
+Nlocal:    2123 ave 2123 max 2123 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    270799 ave 270799 max 270799 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 270799
+Ave neighs/atom = 127.555
+Neighbor list builds = 11
+Dangerous builds = 0
+reset_timestep  0
+
+velocity        all create 1.2 12351235
+
+fix             1   all  nve
+fix             2   all  wall/reflect  xlo EDGE xhi EDGE ylo EDGE yhi EDGE zlo EDGE zhi EDGE
+
+compute 1 all chunk/atom molecule
+compute 2 ball omega/chunk 1
+compute 3 ball inertia/chunk 1
+
+# compute rotational kinetic energy: 1/2 * I * omega**2
+variable r_ke_x equal c_2[1][1]*c_2[1][1]*c_3[1][1]*0.5
+variable r_ke_y equal c_2[1][2]*c_2[1][2]*c_3[1][2]*0.5
+variable r_ke_z equal c_2[1][3]*c_2[1][3]*c_3[1][3]*0.5
+
+# output moments of inertia for x,y,z and angular velocity as well as rotational kinertic energy
+
+thermo_style custom step ke pe c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z
+thermo 500
+
+# dump 1 all atom 100 dump.lammpstrj
+
+# dump 2 all movie 10 ball.mp4 type mass
+
+# equilibration w/o torque added
+run 1000 post no
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 4.3
+  ghost atom cutoff = 4.3
+  binsize = 2.15 -> bins = 16 16 16
+Memory usage per processor = 2.68359 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+       0    1.7991521   -7.3616908    51125.265    51125.265    51125.265 0.0034331372 -0.0045852283 0.0091015032   0.30129221   0.53743693     2.117541 
+     500    0.8882476     -6.43927    52585.112    52205.936    52112.109 0.0033058719 -0.0043753253 0.0089502237    0.2873458   0.49970143    2.0872595 
+    1000   0.79407357    -6.337372    53754.145    54093.977    53787.624 0.0033214912 -0.0042125031 0.0090171557   0.29651603   0.47995375    2.1867116 
+Loop time of 4.26119 on 1 procs for 1000 steps with 2123 atoms
+
+100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+# start spinning the ball. rotation around z should increase and Erot_z should grow
+fix 4 ball addtorque 0.0 0.0 200.0
+run 10000 upto  post no
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Memory usage per processor = 2.68359 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+    1000   0.79407357    -6.337372    53754.145    54093.977    53787.624 0.0033214912 -0.0042125031 0.0090171557   0.29651603   0.47995375    2.1867116 
+    1500   0.68685306   -6.2226287    55026.889    54809.958    55224.858 0.0026096779 -0.0039390202   0.01797948   0.18737807   0.42521238    8.9260406 
+    2000   0.65370325   -6.1832475    54914.897    55655.542    55405.781 0.0029310978 -0.0040761978  0.025816421   0.23589612   0.46236922    18.463634 
+    2500   0.69337585   -6.2170462     54604.66    54800.001    54488.865 0.0028821313 -0.0045216915  0.035781895   0.22679174   0.56021203    34.882245 
+    3000   0.76778067   -6.2850756    53423.198    53620.349    53692.133  0.004088872 -0.004451787  0.044703139   0.44658786   0.53133496     53.64839 
+    3500   0.79707238   -6.3044974    53055.225    53071.129    52927.537 0.0036751739 -0.0037584362  0.054889715    0.3583059    0.3748372    79.732181 
+    4000   0.80746429   -6.3010044    53519.853    53643.284    54058.105  0.003813517 -0.0041637733  0.062983015   0.38916725   0.46500703    107.22047 
+    4500   0.81206394   -6.2884719    53371.354    53883.202    53854.559   0.00385001 -0.003643712  0.072544638   0.39555051   0.35769386    141.71085 
+    5000   0.81648687   -6.2737414    53776.764    54233.367    53762.425 0.0024484228 -0.0043310965  0.080673643   0.16118978   0.50866551    174.94929 
+    5500   0.81888245   -6.2572796     53908.22    53502.342    54717.506 0.0037110524  -0.00327586  0.088836946   0.37120958   0.28707375    215.91536 
+    6000   0.86533749   -6.2804248    53687.533    53571.135    53536.171 0.0025223465 -0.0047753224  0.099646475   0.17078626   0.61081016    265.79156 
+    6500   0.88029206   -6.2719195     53344.67    53291.068    53298.665  0.003937416 -0.0033910578   0.10778737   0.41350774   0.30640427    309.61504 
+    7000    0.9066019   -6.2714707    53928.646    53524.142    54003.175 0.0028500736 -0.0039730774   0.11855069   0.21902903    0.4224485     379.4875 
+    7500   0.94601421   -6.2819912    53534.525    53547.598    53851.344 0.0028610722 -0.0049440438   0.12716472   0.21910969    0.6544472    435.41142 
+    8000    0.9562253   -6.2626222    53486.577    53033.175    53858.803 0.0025501008 -0.0048075887   0.13526164   0.17391198   0.61287549    492.69254 
+    8500   0.99679401   -6.2691139    53566.766    53721.523    53831.283 0.0034137155   -0.0039265   0.14392854    0.3121189   0.41412316    557.56894 
+    9000    1.0371056   -6.2761647    53415.879    53605.078     53985.02 0.0029321914 -0.0046362889    0.1518846    0.2296281   0.57612526     622.6884 
+    9500    1.0598491     -6.26216    53493.003    53049.859    53699.774 0.0032127038 -0.0050624912   0.16002437   0.27606311   0.67980256    687.56651 
+   10000    1.1014855   -6.2654945     53418.49    53251.964     53867.56 0.0037330249 -0.0056278017   0.17103133   0.37220602   0.84330216    787.85924 
+Loop time of 43.025 on 1 procs for 9000 steps with 2123 atoms
+
+100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+# continue without adding more torque. rotation continues at fixed speed
+unfix 4
+run 5000
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Memory usage per processor = 2.68359 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+   10000    1.1014855   -6.2654945     53418.49    53251.964     53867.56 0.0037330249 -0.0056278017   0.17103133   0.37220602   0.84330216    787.85924 
+   10500      1.07259   -6.2363698    52825.233    52896.327    55753.551 0.0035928335 -0.0050843805   0.16344484   0.34094601   0.68370948    744.70621 
+   11000    1.0644214   -6.2282099    52016.795    51950.497    54922.101 0.0047316668 -0.0050149935   0.16196531   0.58229343   0.65328165    720.37919 
+   11500    1.0887339   -6.2525622    52602.789    52903.594     54461.53 0.0044295697 -0.0046710153   0.16698036   0.51606197   0.57713546    759.26022 
+   12000    1.0707466    -6.234719    52785.654    52997.192    54943.066 0.0057389353 -0.0030340721   0.16553451   0.86925773    0.2439353    752.76594 
+   12500    1.0758302   -6.2397291    52375.734    52783.309    55011.986 0.0047029783 -0.0023526884   0.16493895   0.57922337   0.14608158    748.29657 
+   13000    1.0716957   -6.2354221    52711.183    52788.224    55234.737 0.0034033406  -0.00206283   0.16427569   0.30526962   0.11231401    745.29615 
+   13500     1.083667   -6.2475953    52698.902    52203.324    55102.881 0.0032012254 -0.0021366488   0.16381832   0.27002507   0.11916109    739.38261 
+   14000     1.085106   -6.2490655    52767.613    52353.974    55225.438 0.0025647931 -0.0022235227    0.1636534   0.17355699   0.12942041    739.53587 
+   14500    1.0838261   -6.2477856    52292.343    51995.567    54735.836  0.001794954 -0.0029396951   0.16409339  0.084239299   0.22466783    736.92607 
+   15000    1.0827419   -6.2468971    51917.584    51388.833    54481.681 0.0017979486 -0.0025793756   0.16196568  0.083914884   0.17094953    714.60575 
+Loop time of 22.7848 on 1 procs for 5000 steps with 2123 atoms
+
+Performance: 94800.138 tau/day, 219.445 timesteps/s
+100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 18.891     | 18.891     | 18.891     |   0.0 | 82.91
+Neigh   | 3.5735     | 3.5735     | 3.5735     |   0.0 | 15.68
+Comm    | 0.005778   | 0.005778   | 0.005778   |   0.0 |  0.03
+Output  | 0.001862   | 0.001862   | 0.001862   |   0.0 |  0.01
+Modify  | 0.27476    | 0.27476    | 0.27476    |   0.0 |  1.21
+Other   |            | 0.0379     |            |       |  0.17
+
+Nlocal:    2123 ave 2123 max 2123 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    242472 ave 242472 max 242472 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 242472
+Ave neighs/atom = 114.212
+Neighbor list builds = 560
+Dangerous builds = 0
+Total wall time: 0:01:22
--- a/examples/USER/misc/addtorque/log.9Nov16.addtorque.g++.4
+++ b/examples/USER/misc/addtorque/log.9Nov16.addtorque.g++.4
@ -0,0 +1,202 @@
+LAMMPS (9 Nov 2016)
+OMP_NUM_THREADS environment is not set. Defaulting to 1 thread. (../comm.cpp:90)
+  using 1 OpenMP thread(s) per MPI task
+
+units           lj
+lattice         fcc 0.8442
+Lattice spacing in x,y,z = 1.6796 1.6796 1.6796
+boundary        f f f
+
+region          ball   sphere   0.0 0.0 0.0 5.0
+region          box    block  -10 10 -10 10 -10 10
+region          half   block  -10  0 -10 10 -10 10
+
+# add molecule ids so we can use chunk/atom
+fix             0   all property/atom mol ghost yes
+
+create_box      2  box
+Created orthogonal box = (-16.796 -16.796 -16.796) to (16.796 16.796 16.796)
+  1 by 2 by 2 MPI processor grid
+create_atoms    1  region ball
+Created 2123 atoms
+
+pair_style      lj/cut 4.0
+pair_coeff      * *  1.0 1.0
+mass            *    1.0
+set             group all mol 1
+  2123 settings made for mol
+
+# label half the sphere with a different type for better visualization
+set             region half   type 2
+  1142 settings made for type
+
+# use a dynamic group (may require a patch to fix addtorque with older versions of LAMMPS)
+group           ball   dynamic all   region ball
+dynamic group ball defined
+
+neigh_modify    delay 2  every 1  check yes
+
+minimize        0.0  0.0    1000  10000
+WARNING: Resetting reneighboring criteria during minimization (../min.cpp:168)
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 0 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 4.3
+  ghost atom cutoff = 4.3
+  binsize = 2.15 -> bins = 16 16 16
+Memory usage per processor = 3.77014 Mbytes
+Step Temp E_pair E_mol TotEng Press 
+       0            0   -6.2285099            0   -6.2285099  -0.38871568 
+    1000            0   -7.3616908            0   -7.3616908 -9.7399049e-16 
+Loop time of 8.29581 on 4 procs for 1000 steps with 2123 atoms
+
+89.4% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+Minimization stats:
+  Stopping criterion = max iterations
+  Energy initial, next-to-last, final = 
+        -6.22850993032     -7.36169083402     -7.36169083402
+  Force two-norm initial, final = 197.762 3.40861e-12
+  Force max component initial, final = 7.88704 1.60379e-13
+  Final line search alpha, max atom move = 1 1.60379e-13
+  Iterations, force evaluations = 1000 1994
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 5.8462     | 6.8198     | 8.0872     |  34.0 | 82.21
+Neigh   | 0.021211   | 0.028936   | 0.035891   |   3.8 |  0.35
+Comm    | 0.10672    | 1.3842     | 2.3694     |  76.2 | 16.69
+Output  | 2.8e-05    | 3.35e-05   | 5e-05      |   0.2 |  0.00
+Modify  | 0          | 0          | 0          |   0.0 |  0.00
+Other   |            | 0.0629     |            |       |  0.76
+
+Nlocal:    530.75 ave 543 max 514 min
+Histogram: 1 0 0 0 0 0 2 0 0 1
+Nghost:    1058.25 ave 1075 max 1046 min
+Histogram: 1 0 0 2 0 0 0 0 0 1
+Neighs:    67699.8 ave 82013 max 55681 min
+Histogram: 1 1 0 0 0 0 1 0 0 1
+
+Total # of neighbors = 270799
+Ave neighs/atom = 127.555
+Neighbor list builds = 11
+Dangerous builds = 0
+reset_timestep  0
+
+velocity        all create 1.2 12351235
+
+fix             1   all  nve
+fix             2   all  wall/reflect  xlo EDGE xhi EDGE ylo EDGE yhi EDGE zlo EDGE zhi EDGE
+
+compute 1 all chunk/atom molecule
+compute 2 ball omega/chunk 1
+compute 3 ball inertia/chunk 1
+
+# compute rotational kinetic energy: 1/2 * I * omega**2
+variable r_ke_x equal c_2[1][1]*c_2[1][1]*c_3[1][1]*0.5
+variable r_ke_y equal c_2[1][2]*c_2[1][2]*c_3[1][2]*0.5
+variable r_ke_z equal c_2[1][3]*c_2[1][3]*c_3[1][3]*0.5
+
+# output moments of inertia for x,y,z and angular velocity as well as rotational kinertic energy
+
+thermo_style custom step ke pe c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z
+thermo 500
+
+# dump 1 all atom 100 dump.lammpstrj
+
+# dump 2 all movie 10 ball.mp4 type mass
+
+# equilibration w/o torque added
+run 1000 post no
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 4.3
+  ghost atom cutoff = 4.3
+  binsize = 2.15 -> bins = 16 16 16
+Memory usage per processor = 2.77052 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+       0    1.7991521   -7.3616908    51125.265    51125.265    51125.265 0.0040198435 -0.0024237256 0.0011306245    0.4130702   0.15016629  0.032677011 
+     500   0.90968729   -6.4603004    52559.075    52593.598    52370.876 0.0039073429 -0.0023505431 0.0011032256   0.40121833   0.14529121  0.031870471 
+    1000   0.80269935   -6.3456699    53966.859     54192.03    54298.931 0.0040875347 -0.0023098179 0.00077010445   0.45083753   0.14456425  0.016101285 
+Loop time of 2.3842 on 4 procs for 1000 steps with 2123 atoms
+
+89.0% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+# start spinning the ball. rotation around z should increase and Erot_z should grow
+fix 4 ball addtorque 0.0 0.0 200.0
+run 10000 upto  post no
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Memory usage per processor = 2.77052 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+    1000   0.80269935   -6.3456699    53966.859     54192.03    54298.931 0.0040875347 -0.0023098179 0.00077010445   0.45083753   0.14456425  0.016101285 
+    1500   0.68922642   -6.2266953    54966.109    54681.607    54947.256 0.0038798951 -0.0016325797 0.0099837591   0.41371853  0.072871893    2.7384461 
+    2000   0.64105167   -6.1739656    54782.995    54570.486    54910.316 0.0034645832 -0.0027243304  0.017763588   0.32878936   0.20251043    8.6633395 
+    2500   0.69929409   -6.2287557     54307.47    53952.925    54538.409 0.0035199565 -0.0022538001  0.028279733   0.33643732   0.13703004    21.808363 
+    3000   0.77013056   -6.2945597    53916.743     53801.81    53955.863 0.0039732989 -0.0024543292  0.037182149   0.42559463   0.16204384    37.297319 
+    3500   0.80807105   -6.3246605    53451.163    53387.178    53474.789 0.0043137676 -0.0020556348  0.047270147   0.49732542   0.11279735    59.743821 
+    4000   0.81007199   -6.3142362    53334.317    53243.089    53685.963 0.0038548696 -0.0031009535  0.055811043    0.3962745   0.25599044    83.612467 
+    4500    0.7850954   -6.2735818    53738.002    53682.367    53639.471 0.0033046861 -0.0018472801  0.065975851   0.29343502  0.091594032    116.74129 
+    5000   0.77992741   -6.2508277    53864.644    53804.867    53877.025 0.0038258186 -0.0028703189  0.073848203   0.39420539    0.2216419    146.91071 
+    5500   0.79428302   -6.2467907    54353.329    53987.578    54234.062 0.0034715133 -0.0030161617  0.082746129   0.32751699   0.24556875    185.66819 
+    6000   0.82211943   -6.2549436    54273.545    53195.299    54061.645 0.0030929934 -0.0031282112  0.090458895   0.25960687   0.26027676    221.18813 
+    6500   0.87630771   -6.2870988    54042.229    53505.982    54492.239 0.0026486452 -0.0024783378   0.10046947   0.18956181    0.1643211    275.02546 
+    7000   0.88614639   -6.2727488    53701.993    52682.206    53948.142 0.0035384498 -0.0035905797   0.11030427   0.33619131   0.33959641    328.19439 
+    7500   0.92102182   -6.2787281    53410.068    52577.853    53132.511 0.0033084315 -0.0031776605   0.11973875    0.2923058   0.26545312     380.8902 
+    8000   0.94010525   -6.2697758    53732.562    53384.271    53795.933 0.0046460009 -0.0032755993   0.12803006    0.5799174   0.28639462    440.90328 
+    8500   0.95359399    -6.252319    53444.305    53558.444    53789.691 0.0037919474 -0.0035729209   0.13616177   0.38423423   0.34185722    498.63117 
+    9000   0.98631627     -6.25197     53115.01    54017.327    53955.578 0.0045084495 -0.0034368377   0.14488919   0.53981096   0.31902236    566.34143 
+    9500    1.0328238   -6.2632597    53287.675    53682.978    53769.692 0.0044595175 -0.0025931203   0.15416363    0.5298739    0.1804895    638.95665 
+   10000    1.0741755   -6.2667785    53103.657    53319.569    53601.523 0.0059237675 -0.0019561182   0.16445664   0.93173079   0.10201096    724.85302 
+Loop time of 21.8398 on 4 procs for 9000 steps with 2123 atoms
+
+90.3% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+# continue without adding more torque. rotation continues at fixed speed
+unfix 4
+run 5000
+WARNING: One or more dynamic groups may not be updated at correct point in timestep (../fix_group.cpp:153)
+Memory usage per processor = 2.77052 Mbytes
+Step KinEng PotEng c_3[1][1] c_3[1][2] c_3[1][3] c_2[1][1] c_2[1][2] c_2[1][3] v_r_ke_x v_r_ke_y v_r_ke_z 
+   10000    1.0741755   -6.2667785    53103.657    53319.569    53601.523 0.0059237675 -0.0019561182   0.16445664   0.93173079   0.10201096    724.85302 
+   10500    1.0588571   -6.2509381    52374.303    52131.544    55020.367 0.0048843769 -0.0018716797   0.15729531   0.62475047  0.091313217    680.65188 
+   11000    1.0554911   -6.2471863    52803.345    52177.891    55200.756 0.0042073234 -0.0024283269   0.15742315   0.46735107   0.15384055    683.99392 
+   11500    1.0559499   -6.2469955    54031.351    52919.728     54882.35 0.0046703444 -0.0016225764   0.15994353   0.58926908   0.06966232    701.99832 
+   12000    1.0311805   -6.2224471     52812.51     52444.13    55356.101 0.0044986993 -0.0019225732   0.15672327   0.53441759  0.096924293    679.83334 
+   12500    1.0423882   -6.2339087    52000.562    52043.013    55003.272 0.0038688875 -0.0022935053    0.1548654   0.38917977   0.13687746    659.57977 
+   13000    1.0548857   -6.2465445    52196.499    52869.024    54622.553 0.0036650563 -0.0025542156    0.1590498   0.35056832   0.17245921    690.88895 
+   13500    1.0443009   -6.2360149    51921.746    53124.078    54750.325 0.0052756473 -0.0011658069   0.15689119   0.72255483  0.036100621    673.83538 
+   14000    1.0505583    -6.241923    51861.696    52153.234    54321.531 0.0037119391 -0.00045576703   0.15738082   0.35728798 0.0054167284    672.73745 
+   14500     1.040343   -6.2316147    52035.588    51680.479    54443.305 0.0026177168 -0.0014795729   0.15428968    0.1782854  0.056567797    648.01978 
+   15000    1.0404962   -6.2322338    52376.795    51327.453    54677.693 0.0025711787 -0.0021695312   0.15403509   0.17313044   0.12079571    648.66363 
+Loop time of 11.9695 on 4 procs for 5000 steps with 2123 atoms
+
+Performance: 180458.440 tau/day, 417.728 timesteps/s
+89.7% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 7.0313     | 8.4017     | 9.6932     |  35.5 | 70.19
+Neigh   | 0.97886    | 1.3633     | 1.7406     |  26.8 | 11.39
+Comm    | 0.34002    | 2.0227     | 3.7911     |  94.9 | 16.90
+Output  | 0.001198   | 0.0015522  | 0.002578   |   1.5 |  0.01
+Modify  | 0.12841    | 0.131      | 0.13336    |   0.5 |  1.09
+Other   |            | 0.04924    |            |       |  0.41
+
+Nlocal:    530.75 ave 534 max 529 min
+Histogram: 1 0 2 0 0 0 0 0 0 1
+Nghost:    1038 ave 1048 max 1029 min
+Histogram: 1 0 1 0 0 0 1 0 0 1
+Neighs:    60482.5 ave 72547 max 49124 min
+Histogram: 1 0 1 0 0 0 1 0 0 1
+
+Total # of neighbors = 241930
+Ave neighs/atom = 113.957
+Neighbor list builds = 555
+Dangerous builds = 0
+Total wall time: 0:00:44
--- a/examples/USER/misc/agni/Al_prb.agni
+++ b/examples/USER/misc/agni/Al_prb.agni
--- a/examples/USER/misc/agni/adatom.data
+++ b/examples/USER/misc/agni/adatom.data
@ -0,0 +1,197 @@
+Data File created from VASP POSCAR
+
+181 atoms
+1 atom types
+
+0 17.121440767 xlo xhi
+0 14.8276026536 ylo yhi
+0 39.3197318979 zlo zhi
+
+Masses
+
+1 26.982
+
+Atoms
+
+1 1 0.0 1.64751140595 15.0
+2 1 2.85357346116 1.64751140595 15.0
+3 1 5.70714692232 1.64751140595 15.0
+4 1 8.56072038348 1.64751140595 15.0
+5 1 11.4142938446 1.64751140595 15.0
+6 1 14.2678673058 1.64751140595 15.0
+7 1 1.42678673058 4.11877851488 15.0
+8 1 4.28036019174 4.11877851488 15.0
+9 1 7.1339336529 4.11877851488 15.0
+10 1 9.98750711406 4.11877851488 15.0
+11 1 12.8410805752 4.11877851488 15.0
+12 1 15.6946540364 4.11877851488 15.0
+13 1 0.0 6.59004562381 15.0
+14 1 2.85357346116 6.59004562381 15.0
+15 1 5.70714692232 6.59004562381 15.0
+16 1 8.56072038348 6.59004562381 15.0
+17 1 11.4142938446 6.59004562381 15.0
+18 1 14.2678673058 6.59004562381 15.0
+19 1 1.42678673058 9.06131273274 15.0
+20 1 4.28036019174 9.06131273274 15.0
+21 1 7.1339336529 9.06131273274 15.0
+22 1 9.98750711406 9.06131273274 15.0
+23 1 12.8410805752 9.06131273274 15.0
+24 1 15.6946540364 9.06131273274 15.0
+25 1 0.0 11.5325798417 15.0
+26 1 2.85357346116 11.5325798417 15.0
+27 1 5.70714692232 11.5325798417 15.0
+28 1 8.56072038348 11.5325798417 15.0
+29 1 11.4142938446 11.5325798417 15.0
+30 1 14.2678673058 11.5325798417 15.0
+31 1 1.42678673058 14.0038469506 15.0
+32 1 4.28036019174 14.0038469506 15.0
+33 1 7.1339336529 14.0038469506 15.0
+34 1 9.98750711406 14.0038469506 15.0
+35 1 12.8410805752 14.0038469506 15.0
+36 1 15.6946540364 14.0038469506 15.0
+37 1 0.0 0.0 17.3299329745
+38 1 2.85357346116 0.0 17.3299329745
+39 1 5.70714692232 0.0 17.3299329745
+40 1 8.56072038348 0.0 17.3299329745
+41 1 11.4142938446 0.0 17.3299329745
+42 1 14.2678673058 0.0 17.3299329745
+43 1 1.42678673058 2.47126710893 17.3299329745
+44 1 4.28036019174 2.47126710893 17.3299329745
+45 1 7.1339336529 2.47126710893 17.3299329745
+46 1 9.98750711406 2.47126710893 17.3299329745
+47 1 12.8410805752 2.47126710893 17.3299329745
+48 1 15.6946540364 2.47126710893 17.3299329745
+49 1 0.0 4.94253421786 17.3299329745
+50 1 2.85357346116 4.94253421786 17.3299329745
+51 1 5.70714692232 4.94253421786 17.3299329745
+52 1 8.56072038348 4.94253421786 17.3299329745
+53 1 11.4142938446 4.94253421786 17.3299329745
+54 1 14.2678673058 4.94253421786 17.3299329745
+55 1 1.42678673058 7.41380132679 17.3299329745
+56 1 4.28036019174 7.41380132679 17.3299329745
+57 1 7.1339336529 7.41380132679 17.3299329745
+58 1 9.98750711406 7.41380132679 17.3299329745
+59 1 12.8410805752 7.41380132679 17.3299329745
+60 1 15.6946540364 7.41380132679 17.3299329745
+61 1 0.0 9.88506843572 17.3299329745
+62 1 2.85357346116 9.88506843572 17.3299329745
+63 1 5.70714692232 9.88506843572 17.3299329745
+64 1 8.56072038348 9.88506843572 17.3299329745
+65 1 11.4142938446 9.88506843572 17.3299329745
+66 1 14.2678673058 9.88506843572 17.3299329745
+67 1 1.42678673058 12.3563355446 17.3299329745
+68 1 4.28036019174 12.3563355446 17.3299329745
+69 1 7.1339336529 12.3563355446 17.3299329745
+70 1 9.98750711406 12.3563355446 17.3299329745
+71 1 12.8410805752 12.3563355446 17.3299329745
+72 1 15.6946540364 12.3563355446 17.3299329745
+73 1 1.42678673058 0.823755702976 19.6598659489
+74 1 4.28036019174 0.823755702976 19.6598659489
+75 1 7.1339336529 0.823755702976 19.6598659489
+76 1 9.98750711406 0.823755702976 19.6598659489
+77 1 12.8410805752 0.823755702976 19.6598659489
+78 1 15.6946540364 0.823755702976 19.6598659489
+79 1 0.0 3.29502281191 19.6598659489
+80 1 2.85357346116 3.29502281191 19.6598659489
+81 1 5.70714692232 3.29502281191 19.6598659489
+82 1 8.56072038348 3.29502281191 19.6598659489
+83 1 11.4142938446 3.29502281191 19.6598659489
+84 1 14.2678673058 3.29502281191 19.6598659489
+85 1 1.42678673058 5.76628992084 19.6598659489
+86 1 4.28036019174 5.76628992084 19.6598659489
+87 1 7.1339336529 5.76628992084 19.6598659489
+88 1 9.98750711406 5.76628992084 19.6598659489
+89 1 12.8410805752 5.76628992084 19.6598659489
+90 1 15.6946540364 5.76628992084 19.6598659489
+91 1 0.0 8.23755702976 19.6598659489
+92 1 2.85357346116 8.23755702976 19.6598659489
+93 1 5.70714692232 8.23755702976 19.6598659489
+94 1 8.56072038348 8.23755702976 19.6598659489
+95 1 11.4142938446 8.23755702976 19.6598659489
+96 1 14.2678673058 8.23755702976 19.6598659489
+97 1 1.42678673058 10.7088241387 19.6598659489
+98 1 4.28036019174 10.7088241387 19.6598659489
+99 1 7.1339336529 10.7088241387 19.6598659489
+100 1 9.98750711406 10.7088241387 19.6598659489
+101 1 12.8410805752 10.7088241387 19.6598659489
+102 1 15.6946540364 10.7088241387 19.6598659489
+103 1 0.0 13.1800912476 19.6598659489
+104 1 2.85357346116 13.1800912476 19.6598659489
+105 1 5.70714692232 13.1800912476 19.6598659489
+106 1 8.56072038348 13.1800912476 19.6598659489
+107 1 11.4142938446 13.1800912476 19.6598659489
+108 1 14.2678673058 13.1800912476 19.6598659489
+109 1 0.0 1.64751140595 21.9897989234
+110 1 2.85357346116 1.64751140595 21.9897989234
+111 1 5.70714692232 1.64751140595 21.9897989234
+112 1 8.56072038348 1.64751140595 21.9897989234
+113 1 11.4142938446 1.64751140595 21.9897989234
+114 1 14.2678673058 1.64751140595 21.9897989234
+115 1 1.42678673058 4.11877851488 21.9897989234
+116 1 4.28036019174 4.11877851488 21.9897989234
+117 1 7.1339336529 4.11877851488 21.9897989234
+118 1 9.98750711406 4.11877851488 21.9897989234
+119 1 12.8410805752 4.11877851488 21.9897989234
+120 1 15.6946540364 4.11877851488 21.9897989234
+121 1 0.0 6.59004562381 21.9897989234
+122 1 2.85357346116 6.59004562381 21.9897989234
+123 1 5.70714692232 6.59004562381 21.9897989234
+124 1 8.56072038348 6.59004562381 21.9897989234
+125 1 11.4142938446 6.59004562381 21.9897989234
+126 1 14.2678673058 6.59004562381 21.9897989234
+127 1 1.42678673058 9.06131273274 21.9897989234
+128 1 4.28036019174 9.06131273274 21.9897989234
+129 1 7.1339336529 9.06131273274 21.9897989234
+130 1 9.98750711406 9.06131273274 21.9897989234
+131 1 12.8410805752 9.06131273274 21.9897989234
+132 1 15.6946540364 9.06131273274 21.9897989234
+133 1 0.0 11.5325798417 21.9897989234
+134 1 2.85357346116 11.5325798417 21.9897989234
+135 1 5.70714692232 11.5325798417 21.9897989234
+136 1 8.56072038348 11.5325798417 21.9897989234
+137 1 11.4142938446 11.5325798417 21.9897989234
+138 1 14.2678673058 11.5325798417 21.9897989234
+139 1 1.42678673058 14.0038469506 21.9897989234
+140 1 4.28036019174 14.0038469506 21.9897989234
+141 1 7.1339336529 14.0038469506 21.9897989234
+142 1 9.98750711406 14.0038469506 21.9897989234
+143 1 12.8410805752 14.0038469506 21.9897989234
+144 1 15.6946540364 14.0038469506 21.9897989234
+145 1 0.0 0.0 24.3197318979
+146 1 2.85357346116 0.0 24.3197318979
+147 1 5.70714692232 0.0 24.3197318979
+148 1 8.56072038348 0.0 24.3197318979
+149 1 11.4142938446 0.0 24.3197318979
+150 1 14.2678673058 0.0 24.3197318979
+151 1 1.42678673058 2.47126710893 24.3197318979
+152 1 4.28036019174 2.47126710893 24.3197318979
+153 1 7.1339336529 2.47126710893 24.3197318979
+154 1 9.98750711406 2.47126710893 24.3197318979
+155 1 12.8410805752 2.47126710893 24.3197318979
+156 1 15.6946540364 2.47126710893 24.3197318979
+157 1 0.0 4.94253421786 24.3197318979
+158 1 2.85357346116 4.94253421786 24.3197318979
+159 1 5.70714692232 4.94253421786 24.3197318979
+160 1 8.56072038348 4.94253421786 24.3197318979
+161 1 11.4142938446 4.94253421786 24.3197318979
+162 1 14.2678673058 4.94253421786 24.3197318979
+163 1 1.42678673058 7.41380132679 24.3197318979
+164 1 4.28036019174 7.41380132679 24.3197318979
+165 1 7.1339336529 7.41380132679 24.3197318979
+166 1 9.98750711406 7.41380132679 24.3197318979
+167 1 12.8410805752 7.41380132679 24.3197318979
+168 1 15.6946540364 7.41380132679 24.3197318979
+169 1 0.0 9.88506843572 24.3197318979
+170 1 2.85357346116 9.88506843572 24.3197318979
+171 1 5.70714692232 9.88506843572 24.3197318979
+172 1 8.56072038348 9.88506843572 24.3197318979
+173 1 11.4142938446 9.88506843572 24.3197318979
+174 1 14.2678673058 9.88506843572 24.3197318979
+175 1 1.42678673058 12.3563355446 24.3197318979
+176 1 4.28036019174 12.3563355446 24.3197318979
+177 1 7.1339336529 12.3563355446 24.3197318979
+178 1 9.98750711406 12.3563355446 24.3197318979
+179 1 12.8410805752 12.3563355446 24.3197318979
+180 1 15.6946540364 12.3563355446 24.3197318979
+181 1 7.1339336529 4.11877851488 26.7197318979
+
--- a/examples/USER/misc/agni/in.adatom
+++ b/examples/USER/misc/agni/in.adatom
@ -0,0 +1,23 @@
+
+processors * * 1
+units metal
+boundary p p f
+read_data adatom.data
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+
+neighbor 0.3 bin 
+neigh_modify delay 2 check yes 
+
+timestep 0.0005
+velocity all create 500 12345
+fix 1 all nvt temp 250 250 0.2 
+fix 5 all momentum 1 linear 1 1 1
+
+thermo 100
+thermo_style custom step ke temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
--- a/examples/USER/misc/agni/in.vacancy
+++ b/examples/USER/misc/agni/in.vacancy
@ -0,0 +1,23 @@
+
+units metal
+boundary p p p
+read_data vacancy.data
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+neighbor 0.3 bin 
+neigh_modify delay 2 check yes 
+
+
+timestep 0.0005
+velocity all create 1000 12345
+fix 1 all nvt temp 900 900 200
+fix 5 all momentum 1 linear 1 1 1
+
+
+thermo 100
+thermo_style custom step ke etotal temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
--- a/examples/USER/misc/agni/log.9Nov16.adatom.g++.1
+++ b/examples/USER/misc/agni/log.9Nov16.adatom.g++.1
@ -0,0 +1,82 @@
+LAMMPS (9 Nov 2016)
+  using 1 OpenMP thread(s) per MPI task
+
+processors * * 1
+units metal
+boundary p p f
+read_data adatom.data
+  orthogonal box = (0 0 0) to (17.1214 14.8276 39.3197)
+  1 by 1 by 1 MPI processor grid
+  reading atoms ...
+  181 atoms
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+Reading potential file Al_prb.agni with DATE: 2016-11-11
+
+neighbor 0.3 bin
+neigh_modify delay 2 check yes
+
+timestep 0.0005
+velocity all create 500 12345
+fix 1 all nvt temp 250 250 0.2
+fix 5 all momentum 1 linear 1 1 1
+
+thermo 100
+thermo_style custom step ke temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 8.3
+  ghost atom cutoff = 8.3
+  binsize = 4.15 -> bins = 5 4 10
+Memory usage per processor = 2.37049 Mbytes
+Step KinEng Temp 
+       0    11.633413          500 
+     100     4.049974    174.06646 
+     200    5.8983472    253.50889 
+     300    5.3667309    230.66021 
+     400    4.9343935     212.0785 
+     500    5.4054496    232.32432 
+     600    6.1779127    265.52452 
+     700    6.3749266     273.9921 
+     800    6.0701481    260.89283 
+     900    6.4582394    277.57286 
+    1000    6.4047444    275.27366 
+Loop time of 20.8273 on 1 procs for 1000 steps with 181 atoms
+
+Performance: 2.074 ns/day, 11.571 hours/ns, 48.014 timesteps/s
+100.0% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 20.79      | 20.79      | 20.79      |   0.0 | 99.82
+Neigh   | 0.022742   | 0.022742   | 0.022742   |   0.0 |  0.11
+Comm    | 0.0040836  | 0.0040836  | 0.0040836  |   0.0 |  0.02
+Output  | 0.00011086 | 0.00011086 | 0.00011086 |   0.0 |  0.00
+Modify  | 0.0089345  | 0.0089345  | 0.0089345  |   0.0 |  0.04
+Other   |            | 0.001819   |            |       |  0.01
+
+Nlocal:    181 ave 181 max 181 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    562 ave 562 max 562 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  18810 ave 18810 max 18810 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 18810
+Ave neighs/atom = 103.923
+Neighbor list builds = 33
+Dangerous builds = 0
+
+Please see the log.cite file for references relevant to this simulation
+
+Total wall time: 0:00:20
--- a/examples/USER/misc/agni/log.9Nov16.adatom.g++.4
+++ b/examples/USER/misc/agni/log.9Nov16.adatom.g++.4
@ -0,0 +1,82 @@
+LAMMPS (9 Nov 2016)
+  using 1 OpenMP thread(s) per MPI task
+
+processors * * 1
+units metal
+boundary p p f
+read_data adatom.data
+  orthogonal box = (0 0 0) to (17.1214 14.8276 39.3197)
+  2 by 2 by 1 MPI processor grid
+  reading atoms ...
+  181 atoms
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+Reading potential file Al_prb.agni with DATE: 2016-11-11
+
+neighbor 0.3 bin
+neigh_modify delay 2 check yes
+
+timestep 0.0005
+velocity all create 500 12345
+fix 1 all nvt temp 250 250 0.2
+fix 5 all momentum 1 linear 1 1 1
+
+thermo 100
+thermo_style custom step ke temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 8.3
+  ghost atom cutoff = 8.3
+  binsize = 4.15 -> bins = 5 4 10
+Memory usage per processor = 2.48695 Mbytes
+Step KinEng Temp 
+       0    11.633413          500 
+     100     4.049974    174.06646 
+     200    5.8983472    253.50889 
+     300    5.3667309    230.66021 
+     400    4.9343935     212.0785 
+     500    5.4054496    232.32432 
+     600    6.1779127    265.52451 
+     700    6.3749266     273.9921 
+     800    6.0701481    260.89283 
+     900    6.4582394    277.57286 
+    1000    6.4047444    275.27366 
+Loop time of 5.96868 on 4 procs for 1000 steps with 181 atoms
+
+Performance: 7.238 ns/day, 3.316 hours/ns, 167.541 timesteps/s
+99.9% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 5.4176     | 5.4602     | 5.5505     |   2.3 | 91.48
+Neigh   | 0.0056074  | 0.0060464  | 0.0062635  |   0.3 |  0.10
+Comm    | 0.39544    | 0.48696    | 0.53111    |   7.9 |  8.16
+Output  | 0.0001545  | 0.00015736 | 0.0001595  |   0.0 |  0.00
+Modify  | 0.010492   | 0.011565   | 0.012588   |   0.9 |  0.19
+Other   |            | 0.003794   |            |       |  0.06
+
+Nlocal:    45.25 ave 47 max 42 min
+Histogram: 1 0 0 0 0 0 0 0 2 1
+Nghost:    374.75 ave 380 max 373 min
+Histogram: 3 0 0 0 0 0 0 0 0 1
+Neighs:    0 ave 0 max 0 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+FullNghs:  4702.5 ave 4868 max 4389 min
+Histogram: 1 0 0 0 0 0 0 1 1 1
+
+Total # of neighbors = 18810
+Ave neighs/atom = 103.923
+Neighbor list builds = 33
+Dangerous builds = 0
+
+Please see the log.cite file for references relevant to this simulation
+
+Total wall time: 0:00:05
--- a/examples/USER/misc/agni/log.9Nov16.vacancy.g++.1
+++ b/examples/USER/misc/agni/log.9Nov16.vacancy.g++.1
@ -0,0 +1,82 @@
+LAMMPS (9 Nov 2016)
+  using 1 OpenMP thread(s) per MPI task
+
+units metal
+boundary p p p
+read_data vac.data
+  orthogonal box = (0 0 0) to (8.07113 8.07113 8.07113)
+  1 by 1 by 1 MPI processor grid
+  reading atoms ...
+  31 atoms
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+Reading potential file Al_prb.agni with DATE: 2016-11-11
+neighbor 0.3 bin
+neigh_modify delay 2 check yes
+
+
+timestep 0.0005
+velocity all create 1000 12345
+fix 1 all nvt temp 900 900 200
+fix 5 all momentum 1 linear 1 1 1
+
+
+thermo 100
+thermo_style custom step ke etotal temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 8.3
+  ghost atom cutoff = 8.3
+  binsize = 4.15 -> bins = 2 2 2
+Memory usage per processor = 2.40842 Mbytes
+Step KinEng TotEng Temp 
+       0    3.8778043    3.8778043         1000 
+     100    2.8126642    2.8126642    725.32391 
+     200    3.7110413    3.7110413     956.9955 
+     300    3.2361084    3.2361084    834.52081 
+     400    3.4625769    3.4625769    892.92201 
+     500    3.4563307    3.4563307    891.31126 
+     600    2.8486344    2.8486344    734.59982 
+     700    3.1183057    3.1183057    804.14208 
+     800    2.9164818    2.9164818    752.09618 
+     900     3.464416     3.464416    893.39629 
+    1000    3.5954546    3.5954546    927.18825 
+Loop time of 3.86777 on 1 procs for 1000 steps with 31 atoms
+
+Performance: 11.169 ns/day, 2.149 hours/ns, 258.547 timesteps/s
+99.9% CPU use with 1 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 3.8463     | 3.8463     | 3.8463     |   0.0 | 99.44
+Neigh   | 0.011294   | 0.011294   | 0.011294   |   0.0 |  0.29
+Comm    | 0.0057271  | 0.0057271  | 0.0057271  |   0.0 |  0.15
+Output  | 0.00014257 | 0.00014257 | 0.00014257 |   0.0 |  0.00
+Modify  | 0.0029459  | 0.0029459  | 0.0029459  |   0.0 |  0.08
+Other   |            | 0.001395   |            |       |  0.04
+
+Nlocal:    31 ave 31 max 31 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    878 ave 878 max 878 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    0 ave 0 max 0 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+FullNghs:  4334 ave 4334 max 4334 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 4334
+Ave neighs/atom = 139.806
+Neighbor list builds = 51
+Dangerous builds = 0
+
+Please see the log.cite file for references relevant to this simulation
+
+Total wall time: 0:00:03
--- a/examples/USER/misc/agni/log.9Nov16.vacancy.g++.4
+++ b/examples/USER/misc/agni/log.9Nov16.vacancy.g++.4
@ -0,0 +1,82 @@
+LAMMPS (9 Nov 2016)
+  using 1 OpenMP thread(s) per MPI task
+
+units metal
+boundary p p p
+read_data vac.data
+  orthogonal box = (0 0 0) to (8.07113 8.07113 8.07113)
+  1 by 2 by 2 MPI processor grid
+  reading atoms ...
+  31 atoms
+
+pair_style agni
+pair_coeff * * Al_prb.agni Al
+Reading potential file Al_prb.agni with DATE: 2016-11-11
+neighbor 0.3 bin
+neigh_modify delay 2 check yes
+
+
+timestep 0.0005
+velocity all create 1000 12345
+fix 1 all nvt temp 900 900 200
+fix 5 all momentum 1 linear 1 1 1
+
+
+thermo 100
+thermo_style custom step ke etotal temp
+
+# dump MyDump all custom 250 dump.atoms id type x y z vx vy vz fx fy fz
+
+run 1000
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 2 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 8.3
+  ghost atom cutoff = 8.3
+  binsize = 4.15 -> bins = 2 2 2
+Memory usage per processor = 2.3974 Mbytes
+Step KinEng TotEng Temp 
+       0    3.8778044    3.8778044         1000 
+     100    2.8126642    2.8126642    725.32391 
+     200    3.7110413    3.7110413    956.99549 
+     300    3.2361084    3.2361084    834.52081 
+     400    3.4625769    3.4625769    892.92201 
+     500    3.4563307    3.4563307    891.31126 
+     600    2.8486343    2.8486343    734.59981 
+     700    3.1183057    3.1183057    804.14208 
+     800    2.9164819    2.9164819    752.09618 
+     900    3.4644161    3.4644161    893.39631 
+    1000    3.5954546    3.5954546    927.18824 
+Loop time of 1.11007 on 4 procs for 1000 steps with 31 atoms
+
+Performance: 38.916 ns/day, 0.617 hours/ns, 900.843 timesteps/s
+99.8% CPU use with 4 MPI tasks x 1 OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 0.906      | 0.99469    | 1.0291     |   5.1 | 89.61
+Neigh   | 0.0026186  | 0.0027665  | 0.0028622  |   0.2 |  0.25
+Comm    | 0.066125   | 0.10079    | 0.1896     |  16.2 |  9.08
+Output  | 0.00012875 | 0.00014615 | 0.00018787 |   0.2 |  0.01
+Modify  | 0.0080338  | 0.0083079  | 0.00861    |   0.2 |  0.75
+Other   |            | 0.003372   |            |       |  0.30
+
+Nlocal:    7.75 ave 8 max 7 min
+Histogram: 1 0 0 0 0 0 0 0 0 3
+Nghost:    623.5 ave 630 max 616 min
+Histogram: 1 0 0 0 0 2 0 0 0 1
+Neighs:    0 ave 0 max 0 min
+Histogram: 4 0 0 0 0 0 0 0 0 0
+FullNghs:  1083.5 ave 1131 max 988 min
+Histogram: 1 0 0 0 0 0 0 1 1 1
+
+Total # of neighbors = 4334
+Ave neighs/atom = 139.806
+Neighbor list builds = 51
+Dangerous builds = 0
+
+Please see the log.cite file for references relevant to this simulation
+
+Total wall time: 0:00:01
--- a/examples/USER/misc/agni/vacancy.data
+++ b/examples/USER/misc/agni/vacancy.data
@ -0,0 +1,46 @@
+Data File created from VASP POSCAR
+
+31 atoms
+1 atom types
+
+0 8.071125 xlo xhi
+0 8.071125 ylo yhi
+0 8.071125 zlo zhi
+
+Masses
+
+1 26.9815
+
+Atoms
+
+1 1 8.05711986217 3.20498589607 7.09652861184
+2 1 8.05262255028 3.62006786258 3.16719841667
+3 1 2.08891866821 1.38430927213 3.14852514324
+4 1 4.25446836692 3.27689661974 3.35678388118
+5 1 7.92524269451 7.20500664579 3.03232792051
+6 1 6.04056771113 7.24499020906 1.11223380379
+7 1 2.32585852889 5.29910389395 7.31500292009
+8 1 2.09613190567 1.27658214906 7.44277603054
+9 1 3.96852985867 7.2805082905 3.37568009522
+10 1 0.0773420461671 1.29964047903 5.27451616984
+11 1 7.96501442334 1.24471347504 1.17853896176
+12 1 2.13035246804 5.36148411996 3.3817805118
+13 1 2.06211525033 7.25482811482 1.52039033766
+14 1 3.99735704234 7.4099829467 7.05753768668
+15 1 3.84113228596 5.1855444403 1.41642147402
+16 1 0.231862769544 5.38528175164 5.51171817022
+17 1 0.12718452785 5.35814065671 1.11669573581
+18 1 8.05303937039 7.38861123542 7.41398359808
+19 1 1.88506066609 3.17578974033 1.20929473631
+20 1 4.33739926831 1.37976783613 5.28141762358
+21 1 2.23200994743 3.12419127088 5.36881641316
+22 1 6.22871004896 1.34968648416 7.24032447626
+23 1 6.08380394159 1.16222146146 3.30535465675
+24 1 6.16629028099 5.22806528503 3.7675179617
+25 1 4.30194966153 1.14526017671 1.45054175732
+26 1 6.24221620153 5.05377575942 7.17573714759
+27 1 3.92820642281 2.9627641757 7.71515743722
+28 1 4.33686872315 4.73096617728 5.57649231331
+29 1 6.05033104136 3.51389714904 1.34127903322
+30 1 6.27311587476 7.19257797516 5.46814369382
+31 1 1.81274009101 7.47392095028 5.35484578074
--- a/examples/coreshell/in.coreshell.dsf
+++ b/examples/coreshell/in.coreshell.dsf
@ -0,0 +1,71 @@
+# Testsystem for core-shell model compared to Mitchel and Finchham
+# Hendrik Heenen, June 2014
+
+# ------------------------ INITIALIZATION ----------------------------
+
+units 		metal
+dimension	3
+boundary	p	p	p
+atom_style	full
+
+# ----------------------- ATOM DEFINITION ----------------------------
+
+fix csinfo all property/atom i_CSID
+read_data data.coreshell fix csinfo NULL CS-Info
+
+group cores type 1 2
+group shells type 3 4
+
+neighbor 2.0 bin
+comm_modify vel yes
+
+# ------------------------ FORCE FIELDS ------------------------------
+
+pair_style   born/coul/dsf/cs 0.1 20.0 20.0    # A, rho, sigma=0, C, D 
+pair_coeff   * *      0.0 1.000   0.00  0.00   0.00
+pair_coeff   3 3    487.0 0.23768 0.00  1.05   0.50 #Na-Na
+pair_coeff   3 4 145134.0 0.23768 0.00  6.99   8.70 #Na-Cl
+pair_coeff   4 4 405774.0 0.23768 0.00 72.40 145.40 #Cl-Cl
+
+bond_style harmonic
+bond_coeff 1 63.014 0.0
+bond_coeff 2 25.724 0.0
+
+# ------------------------ Equilibration Run -------------------------------
+
+reset_timestep 0
+
+thermo 50
+thermo_style custom step etotal pe ke temp press &
+             epair evdwl ecoul elong ebond fnorm fmax vol
+
+compute CSequ all temp/cs cores shells 
+
+# output via chunk method
+
+#compute prop all property/atom i_CSID
+#compute cs_chunk all chunk/atom c_prop
+#compute cstherm all temp/chunk cs_chunk temp internal com yes cdof 3.0
+#fix ave_chunk all ave/time 100 1 100 c_cstherm file chunk.dump mode vector 
+
+thermo_modify temp CSequ
+
+# velocity bias option
+
+velocity all create 1427 134 dist gaussian mom yes rot no bias yes temp CSequ
+velocity all scale 1427 temp CSequ
+
+fix thermoberendsen all temp/berendsen 1427 1427 0.4
+fix nve all nve
+fix_modify thermoberendsen temp CSequ
+
+# 2 fmsec timestep
+
+timestep 0.002
+run 500
+
+unfix thermoberendsen
+
+# ------------------------ Dynamic Run -------------------------------
+
+run 1000
--- a/examples/coreshell/log.9Nov16.coreshell.dsf.1
+++ b/examples/coreshell/log.9Nov16.coreshell.dsf.1
@ -0,0 +1,185 @@
+LAMMPS (27 Oct 2016)
+# Testsystem for core-shell model compared to Mitchel and Finchham
+# Hendrik Heenen, June 2014
+
+# ------------------------ INITIALIZATION ----------------------------
+
+units 		metal
+dimension	3
+boundary	p	p	p
+atom_style	full
+
+# ----------------------- ATOM DEFINITION ----------------------------
+
+fix csinfo all property/atom i_CSID
+read_data data.coreshell fix csinfo NULL CS-Info
+  orthogonal box = (0 0 0) to (24.096 24.096 24.096)
+  1 by 1 by 1 MPI processor grid
+  reading atoms ...
+  432 atoms
+  scanning bonds ...
+  1 = max bonds/atom
+  reading bonds ...
+  216 bonds
+  1 = max # of 1-2 neighbors
+  0 = max # of 1-3 neighbors
+  0 = max # of 1-4 neighbors
+  1 = max # of special neighbors
+
+group cores type 1 2
+216 atoms in group cores
+group shells type 3 4
+216 atoms in group shells
+
+neighbor 2.0 bin
+comm_modify vel yes
+
+# ------------------------ FORCE FIELDS ------------------------------
+
+pair_style   born/coul/dsf/cs 0.1 20.0 20.0    # A, rho, sigma=0, C, D
+pair_coeff   * *      0.0 1.000   0.00  0.00   0.00
+pair_coeff   3 3    487.0 0.23768 0.00  1.05   0.50 #Na-Na
+pair_coeff   3 4 145134.0 0.23768 0.00  6.99   8.70 #Na-Cl
+pair_coeff   4 4 405774.0 0.23768 0.00 72.40 145.40 #Cl-Cl
+
+bond_style harmonic
+bond_coeff 1 63.014 0.0
+bond_coeff 2 25.724 0.0
+
+# ------------------------ Equilibration Run -------------------------------
+
+reset_timestep 0
+
+thermo 50
+thermo_style custom step etotal pe ke temp press              epair evdwl ecoul elong ebond fnorm fmax vol
+
+compute CSequ all temp/cs cores shells
+
+# output via chunk method
+
+#compute prop all property/atom i_CSID
+#compute cs_chunk all chunk/atom c_prop
+#compute cstherm all temp/chunk cs_chunk temp internal com yes cdof 3.0
+#fix ave_chunk all ave/time 100 1 100 c_cstherm file chunk.dump mode vector
+
+thermo_modify temp CSequ
+
+# velocity bias option
+
+velocity all create 1427 134 dist gaussian mom yes rot no bias yes temp CSequ
+Neighbor list info ...
+  1 neighbor list requests
+  update every 1 steps, delay 10 steps, check yes
+  max neighbors/atom: 2000, page size: 100000
+  master list distance cutoff = 22
+  ghost atom cutoff = 22
+  binsize = 11 -> bins = 3 3 3
+velocity all scale 1427 temp CSequ
+
+fix thermoberendsen all temp/berendsen 1427 1427 0.4
+fix nve all nve
+fix_modify thermoberendsen temp CSequ
+
+# 2 fmsec timestep
+
+timestep 0.002
+run 500
+Memory usage per processor = 7.04355 Mbytes
+Step TotEng PotEng KinEng Temp Press E_pair E_vdwl E_coul E_long E_bond Fnorm Fmax Volume 
+       0   -635.80596   -675.46362    39.657659         1427   -21302.622   -675.46362    1.6320365   -677.09565            0            0 1.3517686e-14 2.942091e-15      13990.5 
+      50    -633.9898   -666.02679     32.03699    1152.7858   -4578.5681   -668.50431    37.800204   -706.30452            0    2.4775226    14.568073    4.3012389      13990.5 
+     100   -631.89604   -661.96148    30.065442    1081.8436   -3536.6738   -664.61798     39.18583   -703.80381            0    2.6564973    14.677968    3.9051029      13990.5 
+     150   -630.08723   -662.95879    32.871559     1182.816   -109.19506   -665.76772    46.247821   -712.01554            0    2.8089226    15.270039    2.9328953      13990.5 
+     200   -628.55895   -663.97376    35.414806    1274.3296     -1748.35   -666.58439    41.738552   -708.32294            0    2.6106349    14.148282    3.1047826      13990.5 
+     250   -627.28761   -661.92274    34.635123    1246.2743   -1280.4899     -664.917    43.045475   -707.96247            0    2.9942594    14.248617    2.4694705      13990.5 
+     300    -626.6163   -663.65651    37.040209    1332.8164   -1887.9043   -666.35215     40.84964   -707.20179            0    2.6956373    13.142643    1.9263242      13990.5 
+     350   -625.76781   -664.66441    38.896607    1399.6151    -1839.482   -667.47659    40.999206   -708.47579            0    2.8121749    13.601238    1.9262698      13990.5 
+     400   -625.02586   -661.46042    36.434568    1311.0236    -868.2031   -664.40231     43.21398   -707.61629            0    2.9418875    14.945389    2.7493413      13990.5 
+     450    -624.3278   -660.50844    36.180639    1301.8865   -2203.3944   -663.49896    40.008669   -703.50763            0    2.9905179    14.158866    1.7299899      13990.5 
+     500   -623.56254   -661.33839    37.775849    1359.2869   -810.50736   -664.11652    42.993999   -707.11052            0    2.7781274     13.68709    2.9115277      13990.5 
+Loop time of 10.7162 on 1 procs for 500 steps with 432 atoms
+
+Performance: 8.063 ns/day, 2.977 hours/ns, 46.658 timesteps/s
+99.9% CPU use with 1 MPI tasks x no OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 10.478     | 10.478     | 10.478     |   0.0 | 97.78
+Bond    | 0.0029511  | 0.0029511  | 0.0029511  |   0.0 |  0.03
+Neigh   | 0.14159    | 0.14159    | 0.14159    |   0.0 |  1.32
+Comm    | 0.074382   | 0.074382   | 0.074382   |   0.0 |  0.69
+Output  | 0.00054097 | 0.00054097 | 0.00054097 |   0.0 |  0.01
+Modify  | 0.010588   | 0.010588   | 0.010588   |   0.0 |  0.10
+Other   |            | 0.007748   |            |       |  0.07
+
+Nlocal:    432 ave 432 max 432 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    9280 ave 9280 max 9280 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    297636 ave 297636 max 297636 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 297636
+Ave neighs/atom = 688.972
+Ave special neighs/atom = 1
+Neighbor list builds = 20
+Dangerous builds = 0
+
+unfix thermoberendsen
+
+# ------------------------ Dynamic Run -------------------------------
+
+run 1000
+Memory usage per processor = 7.04355 Mbytes
+Step TotEng PotEng KinEng Temp Press E_pair E_vdwl E_coul E_long E_bond Fnorm Fmax Volume 
+     500   -623.56254   -661.33839    37.775849    1359.2869   -810.50736   -664.11652    42.993999   -707.11052            0    2.7781274     13.68709    2.9115277      13990.5 
+     550    -623.5004   -660.74472    37.244326    1340.1611   -1413.4326   -663.99669    41.875014    -705.8717            0    3.2519651    15.097948     2.278405      13990.5 
+     600   -623.46963   -659.61729    36.147655    1300.6997   -521.50578   -662.54994    43.956071   -706.50601            0    2.9326492     14.99649    2.6334959      13990.5 
+     650   -623.49291   -661.50698    38.014069    1367.8588   -1230.0925   -664.21074    42.027844   -706.23859            0    2.7037578    13.982308    1.6247207      13990.5 
+     700    -623.4913   -660.11564     36.62434    1317.8522   -727.89052   -663.24921    43.413397   -706.66261            0    3.1335699    15.009937    2.0563966      13990.5 
+     750   -623.50292   -657.95982      34.4569    1239.8613    636.46644   -661.16971    46.539267   -707.70898            0    3.2098934     15.25993    2.1864622      13990.5 
+     800    -623.5176   -659.92032    36.402711    1309.8773   -912.75799   -662.84989    42.668309    -705.5182            0    2.9295708    13.577516    2.0006099      13990.5 
+     850   -623.44098   -660.92727    37.486295    1348.8679   -550.40358   -664.08308    43.667245   -707.75033            0    3.1558098    14.836208     2.279198      13990.5 
+     900   -623.46361   -661.21737    37.753765    1358.4923    1267.8647   -664.52195     47.67284   -712.19479            0    3.3045765    15.058502     1.886141      13990.5 
+     950   -623.50114   -660.58464    37.083492    1334.3739    1754.7359   -663.48186     48.70363   -712.18549            0     2.897226    15.519042    2.2654928      13990.5 
+    1000   -623.50161   -660.02915    36.527539     1314.369    228.76104   -663.31152    45.374099   -708.68562            0    3.2823685    14.783709    2.4201134      13990.5 
+    1050   -623.45985   -660.57417    37.114321    1335.4832    -1490.604   -663.75391    41.258878   -705.01279            0    3.1797391    14.250262    2.3153255      13990.5 
+    1100   -623.51051   -661.20338    37.692871    1356.3011    1791.7899   -664.01042    48.626451   -712.63687            0     2.807039    15.559872     3.184101      13990.5 
+    1150   -623.51067   -663.19545    39.684776    1427.9758    1023.0584   -666.07723      46.5628   -712.64003            0    2.8817804    13.895322    2.3950292      13990.5 
+    1200   -623.49625    -659.6715    36.175253    1301.6927    1600.2805   -662.62259    48.522365   -711.14495            0    2.9510854    15.567834    2.1677651      13990.5 
+    1250   -623.48282   -660.56735    37.084533    1334.4113   -871.67341   -663.86673    42.560699   -706.42743            0    3.2993759    14.569539    2.0093709      13990.5 
+    1300   -623.47744   -663.63125    40.153811     1444.853    1343.7147   -666.39564    47.104842   -713.50048            0    2.7643857    14.186019    1.4599359      13990.5 
+    1350   -623.49121   -661.42731    37.936096    1365.0531    589.73669   -664.46099    45.947687   -710.40867            0    3.0336821    14.801223    2.7486556      13990.5 
+    1400   -623.50803   -660.03912     36.53109    1314.4968    362.97431   -663.24909    45.772904   -709.02199            0    3.2099708    14.566488    1.9170714      13990.5 
+    1450   -623.51243   -659.65548    36.143052     1300.534    2853.0755    -663.0534    51.355353   -714.40875            0    3.3979157    15.890282    2.5251359      13990.5 
+    1500   -623.51621   -661.87741    38.361201    1380.3496    740.04973   -665.00896    46.208742    -711.2177            0    3.1315492    15.168927    2.4710846      13990.5 
+Loop time of 22.2766 on 1 procs for 1000 steps with 432 atoms
+
+Performance: 7.757 ns/day, 3.094 hours/ns, 44.890 timesteps/s
+99.9% CPU use with 1 MPI tasks x no OpenMP threads
+
+MPI task timing breakdown:
+Section |  min time  |  avg time  |  max time  |%varavg| %total
+---------------------------------------------------------------
+Pair    | 21.8       | 21.8       | 21.8       |   0.0 | 97.86
+Bond    | 0.005852   | 0.005852   | 0.005852   |   0.0 |  0.03
+Neigh   | 0.30423    | 0.30423    | 0.30423    |   0.0 |  1.37
+Comm    | 0.14388    | 0.14388    | 0.14388    |   0.0 |  0.65
+Output  | 0.0010855  | 0.0010855  | 0.0010855  |   0.0 |  0.00
+Modify  | 0.0064189  | 0.0064189  | 0.0064189  |   0.0 |  0.03
+Other   |            | 0.01527    |            |       |  0.07
+
+Nlocal:    432 ave 432 max 432 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Nghost:    9318 ave 9318 max 9318 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+Neighs:    297131 ave 297131 max 297131 min
+Histogram: 1 0 0 0 0 0 0 0 0 0
+
+Total # of neighbors = 297131
+Ave neighs/atom = 687.803
+Ave special neighs/atom = 1
+Neighbor list builds = 44
+Dangerous builds = 0
+Total wall time: 0:00:33
--- a/examples/threebody/in.threebody
+++ b/examples/threebody/in.threebody
@ -114,3 +114,35 @@ neighbor        1.0 bin
 neigh_modify    every 1 delay 10 check yes
 run             100

+# Test Tersoff/Mod model for Si
+
+clear
+read_restart	restart.equil
+
+pair_style      tersoff/mod
+pair_coeff 	* * Si.tersoff.mod 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
+
+# Test Tersoff/Mod/C model for Si
+
+clear
+read_restart	restart.equil
+
+pair_style      tersoff/mod/c
+pair_coeff 	* * Si.tersoff.modc 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
+
--- a/potentials/Al_prb.agni
+++ b/potentials/Al_prb.agni
--- a/potentials/Si.tersoff.modc
+++ b/potentials/Si.tersoff.modc
@ -0,0 +1,10 @@
+# DATE: 2016-11-09 CONTRIBUTOR: Ganga P Purja Pun (George Mason University, Fairfax) CITATION: Unknown
+#
+# Format:
+# element1 element2 element3
+# beta alpha h eta
+# beta_ters lam2 B R D lam1 A
+# n c1 c2 c3 c4 c5 C
+Si Si Si 3.00000000 1.80536502 -0.38136087 2.16152496
+1 1.39343356 117.78072440 2.87478837 0.33090566 3.18011795 3198.51383127
+1.98633876 0.20123243 614230.04310619 996439.09714140 3.33560562 25.20963770 -0.00592042
--- a/src/.gitignore
+++ b/src/.gitignore
@ -539,6 +539,8 @@

 /pair_adp.cpp
 /pair_adp.h
+/pair_agni.cpp
+/pair_agni.h
 /pair_airebo.cpp
 /pair_airebo.h
 /pair_airebo_morse.cpp
@ -976,6 +978,8 @@
 /fix_ttm_mod.h
 /pair_born_coul_long_cs.cpp
 /pair_born_coul_long_cs.h
+/pair_born_coul_dsf_cs.cpp
+/pair_born_coul_dsf_cs.h
 /pair_buck_coul_long_cs.cpp
 /pair_buck_coul_long_cs.h
 /pair_coul_long_cs.cpp
--- a/src/ASPHERE/compute_temp_asphere.cpp
+++ b/src/ASPHERE/compute_temp_asphere.cpp
@ -38,7 +38,8 @@ enum{ROTATE,ALL};
 /* ---------------------------------------------------------------------- */

 ComputeTempAsphere::ComputeTempAsphere(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg),
+  id_bias(NULL), tbias(NULL), avec(NULL)
 {
  if (narg < 3) error->all(FLERR,"Illegal compute temp/asphere command");

--- a/src/BODY/body_nparticle.cpp
+++ b/src/BODY/body_nparticle.cpp
@ -28,7 +28,7 @@ enum{SPHERE,LINE,TRI};           // also in DumpImage
 /* ---------------------------------------------------------------------- */

 BodyNparticle::BodyNparticle(LAMMPS *lmp, int narg, char **arg) :
-  Body(lmp, narg, arg)
+  Body(lmp, narg, arg), imflag(NULL), imdata(NULL)
 {
  if (narg != 3) error->all(FLERR,"Invalid body nparticle command");

--- a/src/BODY/compute_body_local.cpp
+++ b/src/BODY/compute_body_local.cpp
@ -33,7 +33,7 @@ enum{ID,TYPE,INDEX};
 /* ---------------------------------------------------------------------- */

 ComputeBodyLocal::ComputeBodyLocal(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg), which(NULL), index(NULL), avec(NULL), bptr(NULL)
 {
  if (narg < 4) error->all(FLERR,"Illegal compute body/local command");

--- a/src/BODY/compute_temp_body.cpp
+++ b/src/BODY/compute_temp_body.cpp
@ -37,7 +37,7 @@ enum{ROTATE,ALL};
 /* ---------------------------------------------------------------------- */

 ComputeTempBody::ComputeTempBody(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg), id_bias(NULL), tbias(NULL), avec(NULL)
 {
  if (narg < 3) error->all(FLERR,"Illegal compute temp/body command");

--- a/src/CORESHELL/Install.sh
+++ b/src/CORESHELL/Install.sh
@ -31,8 +31,10 @@ action () {
 action compute_temp_cs.cpp
 action compute_temp_cs.h
 action pair_born_coul_long_cs.cpp pair_born_coul_long.cpp
+action pair_born_coul_dsf_cs.cpp pair_born_coul_dsf.cpp
 action pair_buck_coul_long_cs.cpp pair_buck_coul_long.cpp
 action pair_born_coul_long_cs.h pair_born_coul_long.h
+action pair_born_coul_dsf_cs.h pair_born_coul_dsf.h
 action pair_buck_coul_long_cs.h pair_buck_coul_long.h
 action pair_coul_long_cs.cpp pair_coul_long.cpp
 action pair_coul_long_cs.h pair_coul_long.h
--- a/src/CORESHELL/compute_temp_cs.cpp
+++ b/src/CORESHELL/compute_temp_cs.cpp
@ -39,7 +39,7 @@ using namespace LAMMPS_NS;
 /* ---------------------------------------------------------------------- */

 ComputeTempCS::ComputeTempCS(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg), vint(NULL), id_fix(NULL), fix(NULL)
 {
  if (narg != 5) error->all(FLERR,"Illegal compute temp/cs command");

--- a/src/CORESHELL/pair_born_coul_dsf_cs.cpp
+++ b/src/CORESHELL/pair_born_coul_dsf_cs.cpp
@ -0,0 +1,162 @@
+/* ----------------------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+/* ----------------------------------------------------------------------
+   Contributing author: Ariel Lozano (arielzn@gmail.com)
+   References: Fennell and Gezelter, JCP 124, 234104 (2006)
+------------------------------------------------------------------------- */
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include "pair_born_coul_dsf_cs.h"
+#include "atom.h"
+#include "comm.h"
+#include "force.h"
+#include "neighbor.h"
+#include "neigh_list.h"
+#include "math_const.h"
+#include "memory.h"
+#include "error.h"
+#include "math_special.h"
+
+using namespace LAMMPS_NS;
+using namespace MathConst;
+
+#define EPSILON 1.0e-20
+
+
+/* ---------------------------------------------------------------------- */
+
+PairBornCoulDSFCS::PairBornCoulDSFCS(LAMMPS *lmp) : PairBornCoulDSF(lmp)
+{
+  writedata = 1;
+  single_enable = 0;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void PairBornCoulDSFCS::compute(int eflag, int vflag)
+{
+  int i,j,ii,jj,inum,jnum,itype,jtype;
+  double qtmp,xtmp,ytmp,ztmp,delx,dely,delz,evdwl,ecoul,fpair;
+  double r,rsq,r2inv,r6inv,forcecoul,forceborn,factor_coul,factor_lj;
+  double prefactor,erfcc,erfcd,arg;
+  double rexp;
+  int *ilist,*jlist,*numneigh,**firstneigh;
+
+  evdwl = ecoul = 0.0;
+  if (eflag || vflag) ev_setup(eflag,vflag);
+  else evflag = vflag_fdotr = 0;
+
+  double **x = atom->x;
+  double **f = atom->f;
+  double *q = atom->q;
+  int *type = atom->type;
+  int nlocal = atom->nlocal;
+  double *special_lj = force->special_lj;
+  double *special_coul = force->special_coul;
+  int newton_pair = force->newton_pair;
+  double qqrd2e = force->qqrd2e;
+
+  inum = list->inum;
+  ilist = list->ilist;
+  numneigh = list->numneigh;
+  firstneigh = list->firstneigh;
+
+  // loop over neighbors of my atoms
+
+  for (ii = 0; ii < inum; ii++) {
+    i = ilist[ii];
+    qtmp = q[i];
+    xtmp = x[i][0];
+    ytmp = x[i][1];
+    ztmp = x[i][2];
+    itype = type[i];
+    jlist = firstneigh[i];
+    jnum = numneigh[i];
+
+  // self coulombic energy
+    if (eflag) {
+      double e_self = -(e_shift/2.0 + alpha/MY_PIS) * qtmp*qtmp*qqrd2e;
+      ev_tally(i,i,nlocal,0,0.0,e_self,0.0,0.0,0.0,0.0);
+    }
+
+    for (jj = 0; jj < jnum; jj++) {
+      j = jlist[jj];
+      factor_lj = special_lj[sbmask(j)];
+      factor_coul = special_coul[sbmask(j)];
+      j &= NEIGHMASK;
+
+      delx = xtmp - x[j][0];
+      dely = ytmp - x[j][1];
+      delz = ztmp - x[j][2];
+      rsq = delx*delx + dely*dely + delz*delz;
+      jtype = type[j];
+
+      if (rsq < cutsq[itype][jtype]) {
+        rsq += EPSILON; // Add Epsilon for case: r = 0; Interaction must be removed by special bond;
+        r2inv = 1.0/rsq;
+
+        if (rsq < cut_coulsq) {
+          r = sqrt(rsq);
+          prefactor = qqrd2e*qtmp*q[j] / r;
+	  arg = alpha * r ;
+          erfcd = MathSpecial::expmsq(arg);
+          erfcc = MathSpecial::my_erfcx(arg) * erfcd;
+          forcecoul = prefactor * (erfcc/r + 2.0*alpha/MY_PIS * erfcd +
+                                   r*f_shift) * r;
+          if (factor_coul < 1.0) forcecoul -= (1.0-factor_coul)*prefactor;
+        } else forcecoul = 0.0;
+
+        if (rsq < cut_ljsq[itype][jtype]) {
+          r6inv = r2inv*r2inv*r2inv;
+          r = sqrt(rsq);
+          rexp = exp((sigma[itype][jtype]-r)*rhoinv[itype][jtype]);
+          forceborn = born1[itype][jtype]*r*rexp - born2[itype][jtype]*r6inv
+            + born3[itype][jtype]*r2inv*r6inv;
+        } else forceborn = 0.0;
+
+        fpair = (forcecoul + factor_lj*forceborn) * r2inv;
+
+        f[i][0] += delx*fpair;
+        f[i][1] += dely*fpair;
+        f[i][2] += delz*fpair;
+        if (newton_pair || j < nlocal) {
+          f[j][0] -= delx*fpair;
+          f[j][1] -= dely*fpair;
+          f[j][2] -= delz*fpair;
+        }
+
+        if (eflag) {
+          if (rsq < cut_coulsq) {
+            ecoul = prefactor * (erfcc - r*e_shift - rsq*f_shift);
+            if (factor_coul < 1.0) ecoul -= (1.0-factor_coul)*prefactor;
+          } else ecoul = 0.0;
+          if (rsq < cut_ljsq[itype][jtype]) {
+            evdwl = a[itype][jtype]*rexp - c[itype][jtype]*r6inv +
+              d[itype][jtype]*r6inv*r2inv - offset[itype][jtype];
+            evdwl *= factor_lj;
+          } else evdwl = 0.0;
+        }
+
+        if (evflag) ev_tally(i,j,nlocal,newton_pair,
+                             evdwl,ecoul,fpair,delx,dely,delz);
+      }
+    }
+  }
+
+  if (vflag_fdotr) virial_fdotr_compute();
+}
+
--- a/src/CORESHELL/pair_born_coul_dsf_cs.h
+++ b/src/CORESHELL/pair_born_coul_dsf_cs.h
@ -0,0 +1,59 @@
+/* -*- c++ -*- ----------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+#ifdef PAIR_CLASS
+
+PairStyle(born/coul/dsf/cs,PairBornCoulDSFCS)
+
+#else
+
+#ifndef LMP_PAIR_BORN_COUL_DSF_CS_H
+#define LMP_PAIR_BORN_COUL_DSF_CS_H
+
+#include "pair_born_coul_dsf.h"
+
+namespace LAMMPS_NS {
+
+class PairBornCoulDSFCS : public PairBornCoulDSF {
+ public:
+  PairBornCoulDSFCS(class LAMMPS *);
+  virtual void compute(int, int);
+};
+
+}
+
+#endif
+#endif
+
+/* ERROR/WARNING messages:
+
+E: Illegal ... command
+
+Self-explanatory.  Check the input script syntax and compare to the
+documentation for the command.  You can use -echo screen as a
+command-line option when running LAMMPS to see the offending line.
+
+E: Incorrect args for pair coefficients
+
+Self-explanatory.  Check the input script or data file.
+
+E: Pair style born/coul/dsf requires atom attribute q
+
+The atom style defined does not have this attribute.
+
+E: All pair coeffs are not set
+
+All pair coefficients must be set in the data file or by the
+pair_coeff command before running a simulation.
+
+*/
--- a/src/Depend.sh
+++ b/src/Depend.sh
@ -53,6 +53,7 @@ fi

 if (test $1 = "CLASS2") then
  depend GPU
+  depend KOKKOS
  depend USER-OMP
 fi

--- a/src/GRANULAR/fix_pour.cpp
+++ b/src/GRANULAR/fix_pour.cpp
@ -47,7 +47,10 @@ enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED};    // several files
 /* ---------------------------------------------------------------------- */

 FixPour::FixPour(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), radius_poly(NULL), frac_poly(NULL), 
+  idrigid(NULL), idshake(NULL), onemols(NULL), molfrac(NULL), coords(NULL), 
+  imageflags(NULL), fixrigid(NULL), fixshake(NULL), recvcounts(NULL), 
+  displs(NULL), random(NULL), random2(NULL)
 {
  if (narg < 6) error->all(FLERR,"Illegal fix pour command");

--- a/src/GRANULAR/fix_wall_gran.cpp
+++ b/src/GRANULAR/fix_wall_gran.cpp
@ -46,7 +46,7 @@ enum{NONE,CONSTANT,EQUAL};
 /* ---------------------------------------------------------------------- */

 FixWallGran::FixWallGran(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), idregion(NULL), shearone(NULL), fix_rigid(NULL), mass_rigid(NULL)
 {
  if (narg < 4) error->all(FLERR,"Illegal fix wall/gran command");

--- a/src/GRANULAR/fix_wall_gran_region.cpp
+++ b/src/GRANULAR/fix_wall_gran_region.cpp
@ -44,7 +44,8 @@ enum{HOOKE,HOOKE_HISTORY,HERTZ_HISTORY,BONDED_HISTORY};
 /* ---------------------------------------------------------------------- */

 FixWallGranRegion::FixWallGranRegion(LAMMPS *lmp, int narg, char **arg) :
-  FixWallGran(lmp, narg, arg)
+  FixWallGran(lmp, narg, arg), region(NULL), region_style(NULL), ncontact(NULL), 
+  walls(NULL), shearmany(NULL), c2r(NULL)
 {
  restart_global = 1;
  motion_resetflag = 0;
--- a/src/KOKKOS/Install.sh
+++ b/src/KOKKOS/Install.sh
@ -59,6 +59,8 @@ action bond_harmonic_kokkos.cpp bond_harmonic.cpp
 action bond_harmonic_kokkos.h bond_harmonic.h
 action comm_kokkos.cpp
 action comm_kokkos.h
+action comm_tiled_kokkos.cpp
+action comm_tiled_kokkos.h
 action compute_temp_kokkos.cpp
 action compute_temp_kokkos.h
 action dihedral_charmm_kokkos.cpp dihedral_charmm.cpp
@ -167,6 +169,8 @@ action pair_reax_c_kokkos.cpp pair_reax_c.cpp
 action pair_reax_c_kokkos.h pair_reax_c.h
 action pair_sw_kokkos.cpp pair_sw.cpp
 action pair_sw_kokkos.h pair_sw.h
+action pair_vashishta_kokkos.cpp pair_vashishta.cpp
+action pair_vashishta_kokkos.h pair_vashishta.h
 action pair_table_kokkos.cpp
 action pair_table_kokkos.h
 action pair_tersoff_kokkos.cpp pair_tersoff.cpp
--- a/src/KOKKOS/atom_vec_kokkos.h
+++ b/src/KOKKOS/atom_vec_kokkos.h
@ -83,8 +83,13 @@ class AtomVecKokkos : public AtomVec {
                   std::is_same<typename ViewType::execution_space,LMPDeviceType>::value,
                   Kokkos::CudaHostPinnedSpace,typename ViewType::memory_space>::type,
                 Kokkos::MemoryTraits<Kokkos::Unmanaged> > mirror_type;
-    if(buffer_size < src.capacity())
+    if (buffer_size == 0) {
+       buffer = Kokkos::kokkos_malloc<Kokkos::CudaHostPinnedSpace>(src.capacity());
+       buffer_size = src.capacity();
+    } else if (buffer_size < src.capacity()) {
       buffer = Kokkos::kokkos_realloc<Kokkos::CudaHostPinnedSpace>(buffer,src.capacity());
+       buffer_size = src.capacity();
+    }
    return mirror_type( buffer ,
                             src.dimension_0() ,
                             src.dimension_1() ,
@ -104,8 +109,13 @@ class AtomVecKokkos : public AtomVec {
                   std::is_same<typename ViewType::execution_space,LMPDeviceType>::value,
                   Kokkos::CudaHostPinnedSpace,typename ViewType::memory_space>::type,
                 Kokkos::MemoryTraits<Kokkos::Unmanaged> > mirror_type;
-    if(buffer_size < src.capacity())
+    if (buffer_size == 0) {
+       buffer = Kokkos::kokkos_malloc<Kokkos::CudaHostPinnedSpace>(src.capacity()*sizeof(typename ViewType::value_type));
+       buffer_size = src.capacity();
+    } else if (buffer_size < src.capacity()) {
       buffer = Kokkos::kokkos_realloc<Kokkos::CudaHostPinnedSpace>(buffer,src.capacity()*sizeof(typename ViewType::value_type));
+       buffer_size = src.capacity();
+    }
    mirror_type tmp_view( (typename ViewType::value_type*)buffer ,
                             src.dimension_0() ,
                             src.dimension_1() ,
--- a/src/KOKKOS/fix_reaxc_bonds_kokkos.cpp
+++ b/src/KOKKOS/fix_reaxc_bonds_kokkos.cpp
@ -103,6 +103,7 @@ void FixReaxCBondsKokkos::Output_ReaxC_Bonds(bigint ntimestep, FILE *fp)
    ((PairReaxCKokkos<LMPDeviceType>*) reaxc)->PackBondBuffer(k_buf,nbuf_local);
  else
    ((PairReaxCKokkos<LMPHostType>*) reaxc)->PackBondBuffer(k_buf,nbuf_local);
+  buf[0] = nlocal;

  // Receive information from buffer for output
  RecvBuffer(buf, nbuf, nbuf_local, nlocal_tot, numbonds_max);
@ -110,6 +111,8 @@ void FixReaxCBondsKokkos::Output_ReaxC_Bonds(bigint ntimestep, FILE *fp)
  memory->destroy_kokkos(k_buf,buf);
 }

+/* ---------------------------------------------------------------------- */
+
 double FixReaxCBondsKokkos::memory_usage()
 {
  double bytes;
--- a/src/KOKKOS/fix_reaxc_bonds_kokkos.h
+++ b/src/KOKKOS/fix_reaxc_bonds_kokkos.h
@ -32,7 +32,7 @@ class FixReaxCBondsKokkos : public FixReaxCBonds {
  void init();

 private:
-  double nbuf;
+  int nbuf;
  void Output_ReaxC_Bonds(bigint, FILE *);
  double memory_usage();
 };
--- a/src/KOKKOS/neighbor_kokkos.cpp
+++ b/src/KOKKOS/neighbor_kokkos.cpp
@ -21,6 +21,11 @@
 #include "atom_masks.h"
 #include "error.h"
 #include "kokkos.h"
+#include "force.h"
+#include "bond.h"
+#include "angle.h"
+#include "dihedral.h"
+#include "improper.h"

 using namespace LAMMPS_NS;

@ -601,7 +606,19 @@ void NeighborKokkos::build_topology_kokkos() {
    k_anglelist.modify<LMPDeviceType>();
    k_dihedrallist.modify<LMPDeviceType>();
    k_improperlist.modify<LMPDeviceType>();
-  } else {
+
+    // Transfer topology neighbor lists to Host for non-Kokkos styles
+ 
+    if (force->bond && force->bond->execution_space == Host)
+      k_bondlist.sync<LMPHostType>();
+    if (force->angle && force->angle->execution_space == Host)
+      k_anglelist.sync<LMPHostType>();
+    if (force->dihedral && force->dihedral->execution_space == Host)
+      k_dihedrallist.sync<LMPHostType>();
+    if (force->improper && force->improper->execution_space == Host)
+      k_improperlist.sync<LMPHostType>();
+
+   } else {
    neighbond_host.build_topology_kk();

    k_bondlist = neighbond_host.k_bondlist;
--- a/src/KOKKOS/pair_tersoff_kokkos.cpp
+++ b/src/KOKKOS/pair_tersoff_kokkos.cpp
@ -126,26 +126,26 @@ void PairTersoffKokkos<DeviceType>::setup_params()
  for (i = 1; i <= n; i++)
    for (j = 1; j <= n; j++)
      for (k = 1; k <= n; k++) {
-	m = elem2param[i-1][j-1][k-1];
-	k_params.h_view(i,j,k).powerm = params[m].powerm;
-	k_params.h_view(i,j,k).gamma = params[m].gamma;
-	k_params.h_view(i,j,k).lam3 = params[m].lam3;
-	k_params.h_view(i,j,k).c = params[m].c;
-	k_params.h_view(i,j,k).d = params[m].d;
-	k_params.h_view(i,j,k).h = params[m].h;
-	k_params.h_view(i,j,k).powern = params[m].powern;
-	k_params.h_view(i,j,k).beta = params[m].beta;
-	k_params.h_view(i,j,k).lam2 = params[m].lam2;
-	k_params.h_view(i,j,k).bigb = params[m].bigb;
-	k_params.h_view(i,j,k).bigr = params[m].bigr;
-	k_params.h_view(i,j,k).bigd = params[m].bigd;
-	k_params.h_view(i,j,k).lam1 = params[m].lam1;
-	k_params.h_view(i,j,k).biga = params[m].biga;
-	k_params.h_view(i,j,k).cutsq = params[m].cutsq;
-	k_params.h_view(i,j,k).c1 = params[m].c1;
-	k_params.h_view(i,j,k).c2 = params[m].c2;
-	k_params.h_view(i,j,k).c3 = params[m].c3;
-	k_params.h_view(i,j,k).c4 = params[m].c4;
+        m = elem2param[map[i]][map[j]][map[k]];
+        k_params.h_view(i,j,k).powerm = params[m].powerm;
+        k_params.h_view(i,j,k).gamma = params[m].gamma;
+        k_params.h_view(i,j,k).lam3 = params[m].lam3;
+        k_params.h_view(i,j,k).c = params[m].c;
+        k_params.h_view(i,j,k).d = params[m].d;
+        k_params.h_view(i,j,k).h = params[m].h;
+        k_params.h_view(i,j,k).powern = params[m].powern;
+        k_params.h_view(i,j,k).beta = params[m].beta;
+        k_params.h_view(i,j,k).lam2 = params[m].lam2;
+        k_params.h_view(i,j,k).bigb = params[m].bigb;
+        k_params.h_view(i,j,k).bigr = params[m].bigr;
+        k_params.h_view(i,j,k).bigd = params[m].bigd;
+        k_params.h_view(i,j,k).lam1 = params[m].lam1;
+        k_params.h_view(i,j,k).biga = params[m].biga;
+        k_params.h_view(i,j,k).cutsq = params[m].cutsq;
+        k_params.h_view(i,j,k).c1 = params[m].c1;
+        k_params.h_view(i,j,k).c2 = params[m].c2;
+        k_params.h_view(i,j,k).c3 = params[m].c3;
+        k_params.h_view(i,j,k).c4 = params[m].c4;
      }

  k_params.template modify<LMPHostType>();
--- a/src/KOKKOS/pair_tersoff_mod_kokkos.cpp
+++ b/src/KOKKOS/pair_tersoff_mod_kokkos.cpp
@ -125,27 +125,27 @@ void PairTersoffMODKokkos<DeviceType>::setup_params()
  for (i = 1; i <= n; i++)
    for (j = 1; j <= n; j++)
      for (k = 1; k <= n; k++) {
-	m = elem2param[i-1][j-1][k-1];
-	k_params.h_view(i,j,k).powerm = params[m].powerm;
-	k_params.h_view(i,j,k).lam3 = params[m].lam3;
-	k_params.h_view(i,j,k).h = params[m].h;
-	k_params.h_view(i,j,k).powern = params[m].powern;
-	k_params.h_view(i,j,k).beta = params[m].beta;
-	k_params.h_view(i,j,k).lam2 = params[m].lam2;
-	k_params.h_view(i,j,k).bigb = params[m].bigb;
-	k_params.h_view(i,j,k).bigr = params[m].bigr;
-	k_params.h_view(i,j,k).bigd = params[m].bigd;
-	k_params.h_view(i,j,k).lam1 = params[m].lam1;
-	k_params.h_view(i,j,k).biga = params[m].biga;
-	k_params.h_view(i,j,k).cutsq = params[m].cutsq;
-	k_params.h_view(i,j,k).c1 = params[m].c1;
-	k_params.h_view(i,j,k).c2 = params[m].c2;
-	k_params.h_view(i,j,k).c3 = params[m].c3;
-	k_params.h_view(i,j,k).c4 = params[m].c4;
-	k_params.h_view(i,j,k).c5 = params[m].c5;
-	k_params.h_view(i,j,k).ca1 = params[m].ca1;
-	k_params.h_view(i,j,k).ca4 = params[m].ca4;
-	k_params.h_view(i,j,k).powern_del = params[m].powern_del;
+        m = elem2param[map[i]][map[j]][map[k]];
+        k_params.h_view(i,j,k).powerm = params[m].powerm;
+        k_params.h_view(i,j,k).lam3 = params[m].lam3;
+        k_params.h_view(i,j,k).h = params[m].h;
+        k_params.h_view(i,j,k).powern = params[m].powern;
+        k_params.h_view(i,j,k).beta = params[m].beta;
+        k_params.h_view(i,j,k).lam2 = params[m].lam2;
+        k_params.h_view(i,j,k).bigb = params[m].bigb;
+        k_params.h_view(i,j,k).bigr = params[m].bigr;
+        k_params.h_view(i,j,k).bigd = params[m].bigd;
+        k_params.h_view(i,j,k).lam1 = params[m].lam1;
+        k_params.h_view(i,j,k).biga = params[m].biga;
+        k_params.h_view(i,j,k).cutsq = params[m].cutsq;
+        k_params.h_view(i,j,k).c1 = params[m].c1;
+        k_params.h_view(i,j,k).c2 = params[m].c2;
+        k_params.h_view(i,j,k).c3 = params[m].c3;
+        k_params.h_view(i,j,k).c4 = params[m].c4;
+        k_params.h_view(i,j,k).c5 = params[m].c5;
+        k_params.h_view(i,j,k).ca1 = params[m].ca1;
+        k_params.h_view(i,j,k).ca4 = params[m].ca4;
+        k_params.h_view(i,j,k).powern_del = params[m].powern_del;
      }

  k_params.template modify<LMPHostType>();
--- a/src/KOKKOS/pair_tersoff_zbl_kokkos.cpp
+++ b/src/KOKKOS/pair_tersoff_zbl_kokkos.cpp
@ -136,30 +136,30 @@ void PairTersoffZBLKokkos<DeviceType>::setup_params()
  for (i = 1; i <= n; i++)
    for (j = 1; j <= n; j++)
      for (k = 1; k <= n; k++) {
-	m = elem2param[i-1][j-1][k-1];
-	k_params.h_view(i,j,k).powerm = params[m].powerm;
-	k_params.h_view(i,j,k).gamma = params[m].gamma;
-	k_params.h_view(i,j,k).lam3 = params[m].lam3;
-	k_params.h_view(i,j,k).c = params[m].c;
-	k_params.h_view(i,j,k).d = params[m].d;
-	k_params.h_view(i,j,k).h = params[m].h;
-	k_params.h_view(i,j,k).powern = params[m].powern;
-	k_params.h_view(i,j,k).beta = params[m].beta;
-	k_params.h_view(i,j,k).lam2 = params[m].lam2;
-	k_params.h_view(i,j,k).bigb = params[m].bigb;
-	k_params.h_view(i,j,k).bigr = params[m].bigr;
-	k_params.h_view(i,j,k).bigd = params[m].bigd;
-	k_params.h_view(i,j,k).lam1 = params[m].lam1;
-	k_params.h_view(i,j,k).biga = params[m].biga;
-	k_params.h_view(i,j,k).cutsq = params[m].cutsq;
-	k_params.h_view(i,j,k).c1 = params[m].c1;
-	k_params.h_view(i,j,k).c2 = params[m].c2;
-	k_params.h_view(i,j,k).c3 = params[m].c3;
-	k_params.h_view(i,j,k).c4 = params[m].c4;
-	k_params.h_view(i,j,k).Z_i = params[m].Z_i;
-	k_params.h_view(i,j,k).Z_j = params[m].Z_j;
-	k_params.h_view(i,j,k).ZBLcut = params[m].ZBLcut;
-	k_params.h_view(i,j,k).ZBLexpscale = params[m].ZBLexpscale;
+        m = elem2param[map[i]][map[j]][map[k]];
+        k_params.h_view(i,j,k).powerm = params[m].powerm;
+        k_params.h_view(i,j,k).gamma = params[m].gamma;
+        k_params.h_view(i,j,k).lam3 = params[m].lam3;
+        k_params.h_view(i,j,k).c = params[m].c;
+        k_params.h_view(i,j,k).d = params[m].d;
+        k_params.h_view(i,j,k).h = params[m].h;
+        k_params.h_view(i,j,k).powern = params[m].powern;
+        k_params.h_view(i,j,k).beta = params[m].beta;
+        k_params.h_view(i,j,k).lam2 = params[m].lam2;
+        k_params.h_view(i,j,k).bigb = params[m].bigb;
+        k_params.h_view(i,j,k).bigr = params[m].bigr;
+        k_params.h_view(i,j,k).bigd = params[m].bigd;
+        k_params.h_view(i,j,k).lam1 = params[m].lam1;
+        k_params.h_view(i,j,k).biga = params[m].biga;
+        k_params.h_view(i,j,k).cutsq = params[m].cutsq;
+        k_params.h_view(i,j,k).c1 = params[m].c1;
+        k_params.h_view(i,j,k).c2 = params[m].c2;
+        k_params.h_view(i,j,k).c3 = params[m].c3;
+        k_params.h_view(i,j,k).c4 = params[m].c4;
+        k_params.h_view(i,j,k).Z_i = params[m].Z_i;
+        k_params.h_view(i,j,k).Z_j = params[m].Z_j;
+        k_params.h_view(i,j,k).ZBLcut = params[m].ZBLcut;
+        k_params.h_view(i,j,k).ZBLexpscale = params[m].ZBLexpscale;
      }

  k_params.template modify<LMPHostType>();
--- a/src/KOKKOS/pair_vashishta_kokkos.cpp
+++ b/src/KOKKOS/pair_vashishta_kokkos.cpp
@ -0,0 +1,961 @@
+/* ----------------------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+/* ----------------------------------------------------------------------
+   Contributing author: Stan Moore (SNL)
+                        Anders Hafreager (UiO), andershaf@gmail.com
+------------------------------------------------------------------------- */
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include "pair_vashishta_kokkos.h"
+#include "kokkos.h"
+#include "pair_kokkos.h"
+#include "atom_kokkos.h"
+#include "neighbor.h"
+#include "neigh_request.h"
+#include "force.h"
+#include "comm.h"
+#include "memory.h"
+#include "neighbor.h"
+#include "neigh_list_kokkos.h"
+#include "memory.h"
+#include "error.h"
+#include "atom_masks.h"
+#include "math_const.h"
+
+using namespace LAMMPS_NS;
+using namespace MathConst;
+
+#define MAXLINE 1024
+#define DELTA 4
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+PairVashishtaKokkos<DeviceType>::PairVashishtaKokkos(LAMMPS *lmp) : PairVashishta(lmp)
+{
+  respa_enable = 0;
+
+  atomKK = (AtomKokkos *) atom;
+  execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
+  datamask_read = X_MASK | F_MASK | TAG_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
+  datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
+}
+
+/* ----------------------------------------------------------------------
+   check if allocated, since class can be destructed when incomplete
+------------------------------------------------------------------------- */
+
+template<class DeviceType>
+PairVashishtaKokkos<DeviceType>::~PairVashishtaKokkos()
+{
+  if (!copymode) {
+    memory->destroy_kokkos(k_eatom,eatom);
+    memory->destroy_kokkos(k_vatom,vatom);
+    eatom = NULL;
+    vatom = NULL;
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+void PairVashishtaKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
+{
+  eflag = eflag_in;
+  vflag = vflag_in;
+
+  if (neighflag == FULL) no_virial_fdotr_compute = 1;
+
+  if (eflag || vflag) ev_setup(eflag,vflag);
+  else evflag = vflag_fdotr = 0;
+
+  // reallocate per-atom arrays if necessary
+
+  if (eflag_atom) {
+    memory->destroy_kokkos(k_eatom,eatom);
+    memory->create_kokkos(k_eatom,eatom,maxeatom,"pair:eatom");
+    d_eatom = k_eatom.d_view;
+  }
+  if (vflag_atom) {
+    memory->destroy_kokkos(k_vatom,vatom);
+    memory->create_kokkos(k_vatom,vatom,maxvatom,6,"pair:vatom");
+    d_vatom = k_vatom.d_view;
+  }
+
+  atomKK->sync(execution_space,datamask_read);
+  if (eflag || vflag) atomKK->modified(execution_space,datamask_modify);
+  else atomKK->modified(execution_space,F_MASK);
+
+  x = atomKK->k_x.view<DeviceType>();
+  f = atomKK->k_f.view<DeviceType>();
+  tag = atomKK->k_tag.view<DeviceType>();
+  type = atomKK->k_type.view<DeviceType>();
+  nlocal = atom->nlocal;
+  newton_pair = force->newton_pair;
+  nall = atom->nlocal + atom->nghost;
+
+  inum = list->inum;
+  const int ignum = inum + list->gnum;
+  NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
+  d_ilist = k_list->d_ilist;
+  d_numneigh = k_list->d_numneigh;
+  d_neighbors = k_list->d_neighbors;
+
+  k_list->clean_copy();
+  copymode = 1;
+
+  EV_FLOAT ev;
+  EV_FLOAT ev_all;
+
+  int max_neighs = d_neighbors.dimension_1();
+
+  if ((d_neighbors_short_2body.dimension_1() != max_neighs) ||
+     (d_neighbors_short_2body.dimension_0() != ignum)) {
+    d_neighbors_short_2body = Kokkos::View<int**,DeviceType>("Vashishta::neighbors_short_2body",ignum,max_neighs);
+  }
+  if (d_numneigh_short_2body.dimension_0()!=ignum) {
+    d_numneigh_short_2body = Kokkos::View<int*,DeviceType>("Vashishta::numneighs_short_2body",ignum);
+  }
+  
+  if ((d_neighbors_short_3body.dimension_1() != max_neighs) ||
+     (d_neighbors_short_3body.dimension_0() != ignum)) {
+    d_neighbors_short_3body = Kokkos::View<int**,DeviceType>("Vashishta::neighbors_short_3body",ignum,max_neighs);
+  }
+
+  if (d_numneigh_short_3body.dimension_0()!=ignum) {
+    d_numneigh_short_3body = Kokkos::View<int*,DeviceType>("Vashishta::numneighs_short_3body",ignum);
+  }
+
+  Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagPairVashishtaComputeShortNeigh>(0,neighflag==FULL?ignum:inum), *this);
+
+
+
+  // loop over neighbor list of my atoms
+
+  if (neighflag == HALF) {
+    if (evflag)
+      Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeHalf<HALF,1> >(0,inum),*this,ev);
+    else
+      Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeHalf<HALF,0> >(0,inum),*this);
+    ev_all += ev;
+  } else if (neighflag == HALFTHREAD) {
+    if (evflag)
+      Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeHalf<HALFTHREAD,1> >(0,inum),*this,ev);
+    else
+      Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeHalf<HALFTHREAD,0> >(0,inum),*this);
+    ev_all += ev;
+  } else if (neighflag == FULL) {
+    if (evflag)
+      Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeFullA<FULL,1> >(0,inum),*this,ev);
+    else
+      Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeFullA<FULL,0> >(0,inum),*this);
+    ev_all += ev;
+
+    if (evflag)
+      Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeFullB<FULL,1> >(0,ignum),*this,ev);
+    else
+      Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairVashishtaComputeFullB<FULL,0> >(0,ignum),*this);
+    ev_all += ev;
+  }
+
+  if (eflag_global) eng_vdwl += ev_all.evdwl;
+  if (vflag_global) {
+    virial[0] += ev_all.v[0];
+    virial[1] += ev_all.v[1];
+    virial[2] += ev_all.v[2];
+    virial[3] += ev_all.v[3];
+    virial[4] += ev_all.v[4];
+    virial[5] += ev_all.v[5];
+  }
+
+  if (vflag_fdotr) pair_virial_fdotr_compute(this);
+
+  if (eflag_atom) {
+    k_eatom.template modify<DeviceType>();
+    k_eatom.template sync<LMPHostType>();
+  }
+
+  if (vflag_atom) {
+    k_vatom.template modify<DeviceType>();
+    k_vatom.template sync<LMPHostType>();
+  }
+
+  copymode = 0;
+}
+
+template<class DeviceType>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeShortNeigh, const int& ii) const {
+    const int i = d_ilist[ii];
+    const int itype = d_map[type[i]];
+    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_2body = 0;
+    int inside_3body = 0;
+    for (int jj = 0; jj < jnum; jj++) {
+      int j = d_neighbors(i,jj);
+      j &= NEIGHMASK;
+      const int jtype = d_map[type[j]];
+      const int ijparam = d_elem2param(itype,jtype,jtype);
+      
+      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 < d_params[ijparam].cutsq) {
+        d_neighbors_short_2body(i,inside_2body) = j;
+        inside_2body++;
+      }
+
+      if (rsq < d_params[ijparam].cutsq2) {
+        d_neighbors_short_3body(i,inside_3body) = j;
+        inside_3body++;
+      }
+    }
+    d_numneigh_short_2body(i) = inside_2body;
+    d_numneigh_short_3body(i) = inside_3body;
+}
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeHalf<NEIGHFLAG,EVFLAG>, const int &ii, EV_FLOAT& ev) const {
+
+  // The f array is atomic
+
+  Kokkos::View<F_FLOAT*[3], typename DAT::t_f_array::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > a_f = f;
+
+  F_FLOAT delr1[3],delr2[3],fj[3],fk[3];
+  F_FLOAT evdwl = 0.0;
+  F_FLOAT fpair = 0.0;
+
+  const int i = d_ilist[ii];
+  const tagint itag = tag[i];
+  const int itype = d_map[type[i]];
+  const X_FLOAT xtmp = x(i,0);
+  const X_FLOAT ytmp = x(i,1);
+  const X_FLOAT ztmp = x(i,2);
+
+  // two-body interactions, skip half of them
+
+  const int jnum = d_numneigh_short_2body[i];
+
+  F_FLOAT fxtmpi = 0.0;
+  F_FLOAT fytmpi = 0.0;
+  F_FLOAT fztmpi = 0.0;
+
+  for (int jj = 0; jj < jnum; jj++) {
+    int j = d_neighbors_short_2body(i,jj);
+    j &= NEIGHMASK;
+    const tagint jtag = tag[j];
+
+    if (itag > jtag) {
+      if ((itag+jtag) % 2 == 0) continue;
+    } else if (itag < jtag) {
+      if ((itag+jtag) % 2 == 1) continue;
+    } else {
+      if (x(j,2) < ztmp) continue;
+      if (x(j,2) == ztmp && x(j,1) < ytmp) continue;
+      if (x(j,2) == ztmp && x(j,1) == ytmp && x(j,0) < xtmp) continue;
+    }
+
+    const int jtype = d_map[type[j]];
+
+    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;
+
+    const int ijparam = d_elem2param(itype,jtype,jtype);
+    
+    twobody(d_params[ijparam],rsq,fpair,eflag,evdwl);
+
+    fxtmpi += delx*fpair;
+    fytmpi += dely*fpair;
+    fztmpi += delz*fpair;
+    a_f(j,0) -= delx*fpair;
+    a_f(j,1) -= dely*fpair;
+    a_f(j,2) -= delz*fpair;
+
+    if (EVFLAG) {
+      if (eflag) ev.evdwl += evdwl;
+      if (vflag_either || eflag_atom) this->template ev_tally<NEIGHFLAG>(ev,i,j,evdwl,fpair,delx,dely,delz);
+    }
+  }
+
+  const int jnumm1 = d_numneigh_short_3body[i];
+
+  for (int jj = 0; jj < jnumm1-1; jj++) {
+    int j = d_neighbors_short_3body(i,jj);
+    j &= NEIGHMASK;
+    const int jtype = d_map[type[j]];
+    const int ijparam = d_elem2param(itype,jtype,jtype);
+    delr1[0] = x(j,0) - xtmp;
+    delr1[1] = x(j,1) - ytmp;
+    delr1[2] = x(j,2) - ztmp;
+    const F_FLOAT rsq1 = delr1[0]*delr1[0] + delr1[1]*delr1[1] + delr1[2]*delr1[2];
+    
+    F_FLOAT fxtmpj = 0.0;
+    F_FLOAT fytmpj = 0.0;
+    F_FLOAT fztmpj = 0.0;
+
+    for (int kk = jj+1; kk < jnumm1; kk++) {
+      int k = d_neighbors_short_3body(i,kk);
+      k &= NEIGHMASK;
+      const int ktype = d_map[type[k]];
+      const int ikparam = d_elem2param(itype,ktype,ktype);
+      const int ijkparam = d_elem2param(itype,jtype,ktype);
+
+      delr2[0] = x(k,0) - xtmp;
+      delr2[1] = x(k,1) - ytmp;
+      delr2[2] = x(k,2) - ztmp;
+      const F_FLOAT rsq2 = delr2[0]*delr2[0] + delr2[1]*delr2[1] + delr2[2]*delr2[2];
+
+      threebody(d_params[ijparam],d_params[ikparam],d_params[ijkparam],
+                rsq1,rsq2,delr1,delr2,fj,fk,eflag,evdwl);
+
+      fxtmpi -= fj[0] + fk[0];
+      fytmpi -= fj[1] + fk[1];
+      fztmpi -= fj[2] + fk[2];
+      fxtmpj += fj[0];
+      fytmpj += fj[1];
+      fztmpj += fj[2];
+      a_f(k,0) += fk[0];
+      a_f(k,1) += fk[1];
+      a_f(k,2) += fk[2];
+
+      if (EVFLAG) {
+        if (eflag) ev.evdwl += evdwl;
+        if (vflag_either || eflag_atom) this->template ev_tally3<NEIGHFLAG>(ev,i,j,k,evdwl,0.0,fj,fk,delr1,delr2);
+      }
+    }
+
+    a_f(j,0) += fxtmpj;
+    a_f(j,1) += fytmpj;
+    a_f(j,2) += fztmpj;
+  }
+
+  a_f(i,0) += fxtmpi;
+  a_f(i,1) += fytmpi;
+  a_f(i,2) += fztmpi;
+}
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeHalf<NEIGHFLAG,EVFLAG>, const int &ii) const {
+  EV_FLOAT ev;
+  this->template operator()<NEIGHFLAG,EVFLAG>(TagPairVashishtaComputeHalf<NEIGHFLAG,EVFLAG>(), ii, ev);
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeFullA<NEIGHFLAG,EVFLAG>, const int &ii, EV_FLOAT& ev) const {
+
+  F_FLOAT delr1[3],delr2[3],fj[3],fk[3];
+  F_FLOAT evdwl = 0.0;
+  F_FLOAT fpair = 0.0;
+
+  const int i = d_ilist[ii];
+
+  const tagint itag = tag[i];
+  const int itype = d_map[type[i]];
+  const X_FLOAT xtmp = x(i,0);
+  const X_FLOAT ytmp = x(i,1);
+  const X_FLOAT ztmp = x(i,2);
+
+  // two-body interactions
+
+  const int jnum = d_numneigh_short_2body[i];
+
+  F_FLOAT fxtmpi = 0.0;
+  F_FLOAT fytmpi = 0.0;
+  F_FLOAT fztmpi = 0.0;
+
+  for (int jj = 0; jj < jnum; jj++) {
+    int j = d_neighbors_short_2body(i,jj);
+    j &= NEIGHMASK;
+    const tagint jtag = tag[j];
+
+    const int jtype = d_map[type[j]];
+
+    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;
+
+    const int ijparam = d_elem2param(itype,jtype,jtype);
+
+    twobody(d_params[ijparam],rsq,fpair,eflag,evdwl);
+
+    fxtmpi += delx*fpair;
+    fytmpi += dely*fpair;
+    fztmpi += delz*fpair;
+
+    if (EVFLAG) {
+      if (eflag) ev.evdwl += 0.5*evdwl;
+      if (vflag_either || eflag_atom) this->template ev_tally<NEIGHFLAG>(ev,i,j,evdwl,fpair,delx,dely,delz);
+    }
+  }
+
+  const int jnumm1 = d_numneigh_short_3body[i];
+
+  for (int jj = 0; jj < jnumm1-1; jj++) {
+    int j = d_neighbors_short_3body(i,jj);
+    j &= NEIGHMASK;
+    const int jtype = d_map[type[j]];
+    const int ijparam = d_elem2param(itype,jtype,jtype);
+    delr1[0] = x(j,0) - xtmp;
+    delr1[1] = x(j,1) - ytmp;
+    delr1[2] = x(j,2) - ztmp;
+    const F_FLOAT rsq1 = delr1[0]*delr1[0] + delr1[1]*delr1[1] + delr1[2]*delr1[2];
+
+    for (int kk = jj+1; kk < jnumm1; kk++) {
+      int k = d_neighbors_short_3body(i,kk);
+      k &= NEIGHMASK;
+      const int ktype = d_map[type[k]];
+      const int ikparam = d_elem2param(itype,ktype,ktype);
+      const int ijkparam = d_elem2param(itype,jtype,ktype);
+
+      delr2[0] = x(k,0) - xtmp;
+      delr2[1] = x(k,1) - ytmp;
+      delr2[2] = x(k,2) - ztmp;
+      const F_FLOAT rsq2 = delr2[0]*delr2[0] + delr2[1]*delr2[1] + delr2[2]*delr2[2];
+
+      threebody(d_params[ijparam],d_params[ikparam],d_params[ijkparam],
+                rsq1,rsq2,delr1,delr2,fj,fk,eflag,evdwl);
+
+      fxtmpi -= fj[0] + fk[0];
+      fytmpi -= fj[1] + fk[1];
+      fztmpi -= fj[2] + fk[2];
+
+      if (EVFLAG) {
+        if (eflag) ev.evdwl += evdwl;
+        if (vflag_either || eflag_atom) this->template ev_tally3<NEIGHFLAG>(ev,i,j,k,evdwl,0.0,fj,fk,delr1,delr2);
+      }
+    }
+  }
+
+  f(i,0) += fxtmpi;
+  f(i,1) += fytmpi;
+  f(i,2) += fztmpi;
+}
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeFullA<NEIGHFLAG,EVFLAG>, const int &ii) const {
+  EV_FLOAT ev;
+  this->template operator()<NEIGHFLAG,EVFLAG>(TagPairVashishtaComputeFullA<NEIGHFLAG,EVFLAG>(), ii, ev);
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeFullB<NEIGHFLAG,EVFLAG>, const int &ii, EV_FLOAT& ev) const {
+
+  F_FLOAT delr1[3],delr2[3],fj[3],fk[3];
+  F_FLOAT evdwl = 0.0;
+  F_FLOAT fpair = 0.0;
+
+  const int i = d_ilist[ii];
+
+  const int itype = d_map[type[i]];
+  const X_FLOAT xtmpi = x(i,0);
+  const X_FLOAT ytmpi = x(i,1);
+  const X_FLOAT ztmpi = x(i,2);
+
+  const int jnum = d_numneigh_short_3body[i];
+
+  F_FLOAT fxtmpi = 0.0;
+  F_FLOAT fytmpi = 0.0;
+  F_FLOAT fztmpi = 0.0;
+
+  for (int jj = 0; jj < jnum; jj++) {
+    int j = d_neighbors_short_3body(i,jj);
+    j &= NEIGHMASK;
+    if (j >= nlocal) continue;
+    const int jtype = d_map[type[j]];
+    const int jiparam = d_elem2param(jtype,itype,itype);
+    const X_FLOAT xtmpj = x(j,0);
+    const X_FLOAT ytmpj = x(j,1);
+    const X_FLOAT ztmpj = x(j,2);
+
+    delr1[0] = xtmpi - xtmpj;
+    delr1[1] = ytmpi - ytmpj;
+    delr1[2] = ztmpi - ztmpj;
+    const F_FLOAT rsq1 = delr1[0]*delr1[0] + delr1[1]*delr1[1] + delr1[2]*delr1[2];
+
+    const int j_jnum = d_numneigh_short_3body[j];
+
+    for (int kk = 0; kk < j_jnum; kk++) {
+      int k = d_neighbors_short_3body(j,kk);
+      k &= NEIGHMASK;
+      if (k == i) continue;
+      const int ktype = d_map[type[k]];
+      const int jkparam = d_elem2param(jtype,ktype,ktype);
+      const int jikparam = d_elem2param(jtype,itype,ktype);
+
+      delr2[0] = x(k,0) - xtmpj;
+      delr2[1] = x(k,1) - ytmpj;
+      delr2[2] = x(k,2) - ztmpj;
+      const F_FLOAT rsq2 = delr2[0]*delr2[0] + delr2[1]*delr2[1] + delr2[2]*delr2[2];
+
+      if (vflag_atom)
+        threebody(d_params[jiparam],d_params[jkparam],d_params[jikparam],
+                  rsq1,rsq2,delr1,delr2,fj,fk,eflag,evdwl);
+      else
+        threebodyj(d_params[jiparam],d_params[jkparam],d_params[jikparam],
+                  rsq1,rsq2,delr1,delr2,fj);
+
+      fxtmpi += fj[0];
+      fytmpi += fj[1];
+      fztmpi += fj[2];
+
+      if (EVFLAG)
+        if (vflag_atom || eflag_atom) ev_tally3_atom(ev,i,evdwl,0.0,fj,fk,delr1,delr2);
+    }
+  }
+
+  f(i,0) += fxtmpi;
+  f(i,1) += fytmpi;
+  f(i,2) += fztmpi;
+}
+
+template<class DeviceType>
+template<int NEIGHFLAG, int EVFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::operator()(TagPairVashishtaComputeFullB<NEIGHFLAG,EVFLAG>, const int &ii) const {
+  EV_FLOAT ev;
+  this->template operator()<NEIGHFLAG,EVFLAG>(TagPairVashishtaComputeFullB<NEIGHFLAG,EVFLAG>(), ii, ev);
+}
+
+/* ----------------------------------------------------------------------
+   set coeffs for one or more type pairs
+------------------------------------------------------------------------- */
+
+template<class DeviceType>
+void PairVashishtaKokkos<DeviceType>::coeff(int narg, char **arg)
+{
+  PairVashishta::coeff(narg,arg);
+
+  // sync map
+
+  int n = atom->ntypes;
+
+  DAT::tdual_int_1d k_map = DAT::tdual_int_1d("pair:map",n+1);
+  HAT::t_int_1d h_map = k_map.h_view;
+
+  for (int i = 1; i <= n; i++)
+    h_map[i] = map[i];
+
+  k_map.template modify<LMPHostType>();
+  k_map.template sync<DeviceType>();
+
+  d_map = k_map.d_view;
+}
+
+/* ----------------------------------------------------------------------
+   init specific to this pair style
+------------------------------------------------------------------------- */
+
+template<class DeviceType>
+void PairVashishtaKokkos<DeviceType>::init_style()
+{
+  PairVashishta::init_style();
+
+  // irequest = neigh request made by parent class
+
+  neighflag = lmp->kokkos->neighflag;
+  int irequest = neighbor->nrequest - 1;
+
+  neighbor->requests[irequest]->
+    kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value &&
+    !Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
+  neighbor->requests[irequest]->
+    kokkos_device = Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
+
+  // always request a full neighbor list
+
+  if (neighflag == FULL || neighflag == HALF || neighflag == HALFTHREAD) {
+    neighbor->requests[irequest]->full = 1;
+    neighbor->requests[irequest]->half = 0;
+    neighbor->requests[irequest]->full_cluster = 0;
+    if (neighflag == FULL)
+      neighbor->requests[irequest]->ghost = 1;
+    else
+      neighbor->requests[irequest]->ghost = 0;
+  } else {
+    error->all(FLERR,"Cannot use chosen neighbor list style with pair vashishta/kk");
+  }
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+void PairVashishtaKokkos<DeviceType>::setup_params()
+{
+  PairVashishta::setup_params();
+
+  // sync elem2param and params
+
+  tdual_int_3d k_elem2param = tdual_int_3d("pair:elem2param",nelements,nelements,nelements);
+  t_host_int_3d h_elem2param = k_elem2param.h_view;
+
+  tdual_param_1d k_params = tdual_param_1d("pair:params",nparams);
+  t_host_param_1d h_params = k_params.h_view;
+
+  for (int i = 0; i < nelements; i++)
+    for (int j = 0; j < nelements; j++)
+      for (int k = 0; k < nelements; k++)
+        h_elem2param(i,j,k) = elem2param[i][j][k];
+
+  for (int m = 0; m < nparams; m++)
+    h_params[m] = params[m];
+
+  k_elem2param.template modify<LMPHostType>();
+  k_elem2param.template sync<DeviceType>();
+  k_params.template modify<LMPHostType>();
+  k_params.template sync<DeviceType>();
+
+  d_elem2param = k_elem2param.d_view;
+  d_params = k_params.d_view;
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::twobody(const Param& param, const F_FLOAT& rsq, F_FLOAT& fforce,
+                     const int& eflag, F_FLOAT& eng) const
+{
+  F_FLOAT r,rinvsq,r4inv,r6inv,reta,lam1r,lam4r,vc2,vc3;
+  r = sqrt(rsq);
+  rinvsq = 1.0/rsq;
+  r4inv = rinvsq*rinvsq;
+  r6inv = rinvsq*r4inv;
+  reta = pow(r,-param.eta);
+  lam1r = r*param.lam1inv;
+  lam4r = r*param.lam4inv;
+  vc2 = param.zizj * exp(-lam1r)/r;
+  vc3 = param.mbigd * r4inv*exp(-lam4r);
+
+  fforce = (param.dvrc*r
+      - (4.0*vc3 + lam4r*vc3+param.big6w*r6inv
+         - param.heta*reta - vc2 - lam1r*vc2)
+      ) * rinvsq;
+
+  if (eflag) eng = param.bigh*reta + vc2 - vc3 - param.bigw*r6inv - r*param.dvrc + param.c0;
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::threebody(const Param& paramij, const Param& paramik, const Param& paramijk,
+                       const F_FLOAT& rsq1, const F_FLOAT& rsq2,
+                       F_FLOAT *delr1, F_FLOAT *delr2,
+                       F_FLOAT *fj, F_FLOAT *fk, const int& eflag, F_FLOAT& eng) const
+{
+  F_FLOAT r1,rinvsq1,rainv1,gsrainv1,gsrainvsq1,expgsrainv1;
+  F_FLOAT r2,rinvsq2,rainv2,gsrainv2,gsrainvsq2,expgsrainv2;
+  F_FLOAT rinv12,cs,delcs,delcssq,facexp,facrad,frad1,frad2,pcsinv,pcsinvsq,pcs;
+  F_FLOAT facang,facang12,csfacang,csfac1,csfac2;
+
+  r1 = sqrt(rsq1);
+  rinvsq1 = 1.0/rsq1;
+  rainv1 = 1.0/(r1 - paramij.r0);
+  gsrainv1 = paramij.gamma * rainv1;
+  gsrainvsq1 = gsrainv1*rainv1/r1;
+  expgsrainv1 = exp(gsrainv1);
+
+  r2 = sqrt(rsq2);
+  rinvsq2 = 1.0/rsq2;
+  rainv2 = 1.0/(r2 - paramik.r0);
+  gsrainv2 = paramik.gamma * rainv2;
+  gsrainvsq2 = gsrainv2*rainv2/r2;
+  expgsrainv2 = exp(gsrainv2);
+
+  rinv12 = 1.0/(r1*r2);
+  cs = (delr1[0]*delr2[0] + delr1[1]*delr2[1] + delr1[2]*delr2[2]) * rinv12;
+  delcs = cs - paramijk.costheta;
+  delcssq = delcs*delcs;
+  pcsinv = paramijk.bigc*delcssq + 1.0;
+  pcsinvsq = pcsinv*pcsinv;
+  pcs = delcssq/pcsinv;
+
+  facexp = expgsrainv1*expgsrainv2;
+
+  facrad = paramijk.bigb * facexp * pcs;
+  frad1 = facrad*gsrainvsq1;
+  frad2 = facrad*gsrainvsq2;
+  facang = paramijk.big2b * facexp * delcs/pcsinvsq;
+  facang12 = rinv12*facang;
+  csfacang = cs*facang;
+  csfac1 = rinvsq1*csfacang;
+
+  fj[0] = delr1[0]*(frad1+csfac1)-delr2[0]*facang12;
+  fj[1] = delr1[1]*(frad1+csfac1)-delr2[1]*facang12;
+  fj[2] = delr1[2]*(frad1+csfac1)-delr2[2]*facang12;
+
+  csfac2 = rinvsq2*csfacang;
+
+  fk[0] = delr2[0]*(frad2+csfac2)-delr1[0]*facang12;
+  fk[1] = delr2[1]*(frad2+csfac2)-delr1[1]*facang12;
+  fk[2] = delr2[2]*(frad2+csfac2)-delr1[2]*facang12;
+
+  if (eflag) eng = facrad;
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::threebodyj(const Param& paramij, const Param& paramik, const Param& paramijk,
+                       const F_FLOAT& rsq1, const F_FLOAT& rsq2, F_FLOAT *delr1, F_FLOAT *delr2, F_FLOAT *fj) const
+{
+  F_FLOAT r1,rinvsq1,rainv1,gsrainv1,gsrainvsq1,expgsrainv1;
+  F_FLOAT r2,rinvsq2,rainv2,gsrainv2,gsrainvsq2,expgsrainv2;
+  F_FLOAT rinv12,cs,delcs,delcssq,facexp,facrad,frad1,frad2,pcsinv,pcsinvsq,pcs;
+  F_FLOAT facang,facang12,csfacang,csfac1,csfac2;
+
+  r1 = sqrt(rsq1);
+  rinvsq1 = 1.0/rsq1;
+  rainv1 = 1.0/(r1 - paramij.r0);
+  gsrainv1 = paramij.gamma * rainv1;
+  gsrainvsq1 = gsrainv1*rainv1/r1;
+  expgsrainv1 = exp(gsrainv1);
+
+  r2 = sqrt(rsq2);
+  rinvsq2 = 1.0/rsq2;
+  rainv2 = 1.0/(r2 - paramik.r0);
+  gsrainv2 = paramik.gamma * rainv2;
+  gsrainvsq2 = gsrainv2*rainv2/r2;
+  expgsrainv2 = exp(gsrainv2);
+
+  rinv12 = 1.0/(r1*r2);
+  cs = (delr1[0]*delr2[0] + delr1[1]*delr2[1] + delr1[2]*delr2[2]) * rinv12;
+  delcs = cs - paramijk.costheta;
+  delcssq = delcs*delcs;
+  pcsinv = paramijk.bigc*delcssq + 1.0;
+  pcsinvsq = pcsinv*pcsinv;
+  pcs = delcssq/pcsinv;
+
+  facexp = expgsrainv1*expgsrainv2;
+
+  facrad = paramijk.bigb * facexp * pcs;
+  frad1 = facrad*gsrainvsq1;
+  frad2 = facrad*gsrainvsq2;
+  facang = paramijk.big2b * facexp * delcs/pcsinvsq;
+  facang12 = rinv12*facang;
+  csfacang = cs*facang;
+  csfac1 = rinvsq1*csfacang;
+
+  fj[0] = delr1[0]*(frad1+csfac1)-delr2[0]*facang12;
+  fj[1] = delr1[1]*(frad1+csfac1)-delr2[1]*facang12;
+  fj[2] = delr1[2]*(frad1+csfac1)-delr2[2]*facang12;
+}
+
+/* ---------------------------------------------------------------------- */
+
+template<class DeviceType>
+template<int NEIGHFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::ev_tally(EV_FLOAT &ev, const int &i, const int &j,
+      const F_FLOAT &epair, const F_FLOAT &fpair, const F_FLOAT &delx,
+                const F_FLOAT &dely, const F_FLOAT &delz) const
+{
+  const int VFLAG = vflag_either;
+
+  // The eatom and vatom arrays are atomic for half/thread neighbor list
+
+  Kokkos::View<E_FLOAT*, typename DAT::t_efloat_1d::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_eatom = k_eatom.view<DeviceType>();
+  Kokkos::View<F_FLOAT*[6], typename DAT::t_virial_array::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_vatom = k_vatom.view<DeviceType>();
+
+
+  if (eflag_atom) {
+    const E_FLOAT epairhalf = 0.5 * epair;
+    v_eatom[i] += epairhalf;
+    if (NEIGHFLAG != FULL)
+      v_eatom[j] += epairhalf;
+  }
+
+  if (VFLAG) {
+    const E_FLOAT v0 = delx*delx*fpair;
+    const E_FLOAT v1 = dely*dely*fpair;
+    const E_FLOAT v2 = delz*delz*fpair;
+    const E_FLOAT v3 = delx*dely*fpair;
+    const E_FLOAT v4 = delx*delz*fpair;
+    const E_FLOAT v5 = dely*delz*fpair;
+
+    if (vflag_global) {
+      if (NEIGHFLAG != FULL) {
+        ev.v[0] += v0;
+        ev.v[1] += v1;
+        ev.v[2] += v2;
+        ev.v[3] += v3;
+        ev.v[4] += v4;
+        ev.v[5] += v5;
+      } else {
+        ev.v[0] += 0.5*v0;
+        ev.v[1] += 0.5*v1;
+        ev.v[2] += 0.5*v2;
+        ev.v[3] += 0.5*v3;
+        ev.v[4] += 0.5*v4;
+        ev.v[5] += 0.5*v5;
+      }
+    }
+
+    if (vflag_atom) {
+      v_vatom(i,0) += 0.5*v0;
+      v_vatom(i,1) += 0.5*v1;
+      v_vatom(i,2) += 0.5*v2;
+      v_vatom(i,3) += 0.5*v3;
+      v_vatom(i,4) += 0.5*v4;
+      v_vatom(i,5) += 0.5*v5;
+
+      if (NEIGHFLAG != FULL) {
+        v_vatom(j,0) += 0.5*v0;
+        v_vatom(j,1) += 0.5*v1;
+        v_vatom(j,2) += 0.5*v2;
+        v_vatom(j,3) += 0.5*v3;
+        v_vatom(j,4) += 0.5*v4;
+        v_vatom(j,5) += 0.5*v5;
+      }
+    }
+  }
+}
+
+/* ----------------------------------------------------------------------
+   tally eng_vdwl and virial into global and per-atom accumulators
+   called by SW and hbond potentials, newton_pair is always on
+   virial = riFi + rjFj + rkFk = (rj-ri) Fj + (rk-ri) Fk = drji*fj + drki*fk
+ ------------------------------------------------------------------------- */
+
+template<class DeviceType>
+template<int NEIGHFLAG>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::ev_tally3(EV_FLOAT &ev, const int &i, const int &j, int &k,
+          const F_FLOAT &evdwl, const F_FLOAT &ecoul,
+                     F_FLOAT *fj, F_FLOAT *fk, F_FLOAT *drji, F_FLOAT *drki) const
+{
+  F_FLOAT epairthird,v[6];
+
+  const int VFLAG = vflag_either;
+
+// The eatom and vatom arrays are atomic for half/thread neighbor list
+
+  Kokkos::View<E_FLOAT*, typename DAT::t_efloat_1d::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_eatom = k_eatom.view<DeviceType>();
+  Kokkos::View<F_FLOAT*[6], typename DAT::t_virial_array::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_vatom = k_vatom.view<DeviceType>();
+
+  if (eflag_atom) {
+    epairthird = THIRD * (evdwl + ecoul);
+    v_eatom[i] += epairthird;
+    if (NEIGHFLAG != FULL) {
+      v_eatom[j] += epairthird;
+      v_eatom[k] += epairthird;
+    }
+  }
+
+  if (VFLAG) {
+    v[0] = drji[0]*fj[0] + drki[0]*fk[0];
+    v[1] = drji[1]*fj[1] + drki[1]*fk[1];
+    v[2] = drji[2]*fj[2] + drki[2]*fk[2];
+    v[3] = drji[0]*fj[1] + drki[0]*fk[1];
+    v[4] = drji[0]*fj[2] + drki[0]*fk[2];
+    v[5] = drji[1]*fj[2] + drki[1]*fk[2];
+
+    if (vflag_global) {
+      ev.v[0] += v[0];
+      ev.v[1] += v[1];
+      ev.v[2] += v[2];
+      ev.v[3] += v[3];
+      ev.v[4] += v[4];
+      ev.v[5] += v[5];
+    }
+
+    if (vflag_atom) {
+      v_vatom(i,0) += THIRD*v[0]; v_vatom(i,1) += THIRD*v[1];
+      v_vatom(i,2) += THIRD*v[2]; v_vatom(i,3) += THIRD*v[3];
+      v_vatom(i,4) += THIRD*v[4]; v_vatom(i,5) += THIRD*v[5];
+
+      if (NEIGHFLAG != FULL) {
+        v_vatom(j,0) += THIRD*v[0]; v_vatom(j,1) += THIRD*v[1];
+        v_vatom(j,2) += THIRD*v[2]; v_vatom(j,3) += THIRD*v[3];
+        v_vatom(j,4) += THIRD*v[4]; v_vatom(j,5) += THIRD*v[5];
+
+        v_vatom(k,0) += THIRD*v[0]; v_vatom(k,1) += THIRD*v[1];
+        v_vatom(k,2) += THIRD*v[2]; v_vatom(k,3) += THIRD*v[3];
+        v_vatom(k,4) += THIRD*v[4]; v_vatom(k,5) += THIRD*v[5];
+      }
+    }
+  }
+}
+
+/* ----------------------------------------------------------------------
+   tally eng_vdwl and virial into global and per-atom accumulators
+   called by SW and hbond potentials, newton_pair is always on
+   virial = riFi + rjFj + rkFk = (rj-ri) Fj + (rk-ri) Fk = drji*fj + drki*fk
+ ------------------------------------------------------------------------- */
+
+template<class DeviceType>
+KOKKOS_INLINE_FUNCTION
+void PairVashishtaKokkos<DeviceType>::ev_tally3_atom(EV_FLOAT &ev, const int &i,
+          const F_FLOAT &evdwl, const F_FLOAT &ecoul,
+                     F_FLOAT *fj, F_FLOAT *fk, F_FLOAT *drji, F_FLOAT *drki) const
+{
+  F_FLOAT epairthird,v[6];
+
+  const int VFLAG = vflag_either;
+
+  if (eflag_atom) {
+    epairthird = THIRD * (evdwl + ecoul);
+    d_eatom[i] += epairthird;
+  }
+
+  if (VFLAG) {
+    v[0] = drji[0]*fj[0] + drki[0]*fk[0];
+    v[1] = drji[1]*fj[1] + drki[1]*fk[1];
+    v[2] = drji[2]*fj[2] + drki[2]*fk[2];
+    v[3] = drji[0]*fj[1] + drki[0]*fk[1];
+    v[4] = drji[0]*fj[2] + drki[0]*fk[2];
+    v[5] = drji[1]*fj[2] + drki[1]*fk[2];
+
+    if (vflag_atom) {
+      d_vatom(i,0) += THIRD*v[0]; d_vatom(i,1) += THIRD*v[1];
+      d_vatom(i,2) += THIRD*v[2]; d_vatom(i,3) += THIRD*v[3];
+      d_vatom(i,4) += THIRD*v[4]; d_vatom(i,5) += THIRD*v[5];
+    }
+  }
+}
+
+namespace LAMMPS_NS {
+template class PairVashishtaKokkos<LMPDeviceType>;
+#ifdef KOKKOS_HAVE_CUDA
+template class PairVashishtaKokkos<LMPHostType>;
+#endif
+}
+
--- a/src/KOKKOS/pair_vashishta_kokkos.h
+++ b/src/KOKKOS/pair_vashishta_kokkos.h
@ -0,0 +1,162 @@
+/* ----------------------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+#ifdef PAIR_CLASS
+
+PairStyle(vashishta/kk,PairVashishtaKokkos<LMPDeviceType>)
+PairStyle(vashishta/kk/device,PairVashishtaKokkos<LMPDeviceType>)
+PairStyle(vashishta/kk/host,PairVashishtaKokkos<LMPHostType>)
+
+#else
+
+#ifndef LMP_PAIR_VASHISHTA_KOKKOS_H
+#define LMP_PAIR_VASHISHTA_KOKKOS_H
+
+#include "pair_vashishta.h"
+#include "pair_kokkos.h"
+
+template<int NEIGHFLAG, int EVFLAG>
+struct TagPairVashishtaComputeHalf{};
+
+template<int NEIGHFLAG, int EVFLAG>
+struct TagPairVashishtaComputeFullA{};
+
+template<int NEIGHFLAG, int EVFLAG>
+struct TagPairVashishtaComputeFullB{};
+
+struct TagPairVashishtaComputeShortNeigh{};
+
+namespace LAMMPS_NS {
+
+template<class DeviceType>
+class PairVashishtaKokkos : public PairVashishta {
+ public:
+  enum {EnabledNeighFlags=FULL};
+  enum {COUL_FLAG=0};
+  typedef DeviceType device_type;
+  typedef ArrayTypes<DeviceType> AT;
+  typedef EV_FLOAT value_type;
+
+  PairVashishtaKokkos(class LAMMPS *);
+  virtual ~PairVashishtaKokkos();
+  virtual void compute(int, int);
+  virtual void coeff(int, char **);
+  virtual void init_style();
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeHalf<NEIGHFLAG,EVFLAG>, const int&, EV_FLOAT&) const;
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeHalf<NEIGHFLAG,EVFLAG>, const int&) const;
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeFullA<NEIGHFLAG,EVFLAG>, const int&, EV_FLOAT&) const;
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeFullA<NEIGHFLAG,EVFLAG>, const int&) const;
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeFullB<NEIGHFLAG,EVFLAG>, const int&, EV_FLOAT&) const;
+
+  template<int NEIGHFLAG, int EVFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeFullB<NEIGHFLAG,EVFLAG>, const int&) const;
+
+  KOKKOS_INLINE_FUNCTION
+  void operator()(TagPairVashishtaComputeShortNeigh, const int&) const;
+
+  template<int NEIGHFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void ev_tally(EV_FLOAT &ev, const int &i, const int &j,
+      const F_FLOAT &epair, const F_FLOAT &fpair, const F_FLOAT &delx,
+                  const F_FLOAT &dely, const F_FLOAT &delz) const;
+
+  template<int NEIGHFLAG>
+  KOKKOS_INLINE_FUNCTION
+  void ev_tally3(EV_FLOAT &ev, const int &i, const int &j, int &k,
+            const F_FLOAT &evdwl, const F_FLOAT &ecoul,
+                       F_FLOAT *fj, F_FLOAT *fk, F_FLOAT *drji, F_FLOAT *drki) const;
+
+  KOKKOS_INLINE_FUNCTION
+  void ev_tally3_atom(EV_FLOAT &ev, const int &i,
+            const F_FLOAT &evdwl, const F_FLOAT &ecoul,
+                       F_FLOAT *fj, F_FLOAT *fk, F_FLOAT *drji, F_FLOAT *drki) const;
+
+ protected:
+  typedef Kokkos::DualView<int***,DeviceType> tdual_int_3d;
+  typedef typename tdual_int_3d::t_dev_const_randomread t_int_3d_randomread;
+  typedef typename tdual_int_3d::t_host t_host_int_3d;
+
+  t_int_3d_randomread d_elem2param;
+  DAT::t_int_1d_randomread d_map;
+
+  typedef Kokkos::DualView<Param*,DeviceType> tdual_param_1d;
+  typedef typename tdual_param_1d::t_dev t_param_1d;
+  typedef typename tdual_param_1d::t_host t_host_param_1d;
+
+  t_param_1d d_params;
+
+  virtual void setup_params();
+  void twobody(const Param&, const F_FLOAT&, F_FLOAT&, const int&, F_FLOAT&) const;
+  void threebody(const Param&, const Param&, const Param&, const F_FLOAT&, const F_FLOAT&, F_FLOAT *, F_FLOAT *,
+                 F_FLOAT *, F_FLOAT *, const int&, F_FLOAT&) const;
+  void threebodyj(const Param&, const Param&, const Param&, const F_FLOAT&, const F_FLOAT&, F_FLOAT *, F_FLOAT *,
+                 F_FLOAT *) const;
+
+  typename ArrayTypes<DeviceType>::t_x_array_randomread x;
+  typename ArrayTypes<DeviceType>::t_f_array f;
+  typename ArrayTypes<DeviceType>::t_tagint_1d tag;
+  typename ArrayTypes<DeviceType>::t_int_1d_randomread type;
+
+  DAT::tdual_efloat_1d k_eatom;
+  DAT::tdual_virial_array k_vatom;
+  DAT::t_efloat_1d d_eatom;
+  DAT::t_virial_array d_vatom;
+
+  DAT::t_int_1d_randomread d_type2frho;
+  DAT::t_int_2d_randomread d_type2rhor;
+  DAT::t_int_2d_randomread d_type2z2r;
+
+  typename ArrayTypes<DeviceType>::t_neighbors_2d d_neighbors;
+  typename ArrayTypes<DeviceType>::t_int_1d_randomread d_ilist;
+  typename ArrayTypes<DeviceType>::t_int_1d_randomread d_numneigh;
+  //NeighListKokkos<DeviceType> k_list;
+
+  int neighflag,newton_pair;
+  int nlocal,nall,eflag,vflag;
+
+  int inum;
+  Kokkos::View<int**,DeviceType> d_neighbors_short_2body;
+  Kokkos::View<int*,DeviceType> d_numneigh_short_2body;
+  Kokkos::View<int**,DeviceType> d_neighbors_short_3body;
+  Kokkos::View<int*,DeviceType> d_numneigh_short_3body;
+  friend void pair_virial_fdotr_compute<PairVashishtaKokkos>(PairVashishtaKokkos*);
+};
+
+}
+
+#endif
+#endif
+
+/* ERROR/WARNING messages:
+
+E: Cannot use chosen neighbor list style with pair vashishta/kk
+
+Self-explanatory.
+
+*/
--- a/src/KOKKOS/region_block_kokkos.h
+++ b/src/KOKKOS/region_block_kokkos.h
@ -33,10 +33,10 @@ template<class DeviceType>
 class RegBlockKokkos : public RegBlock {
  friend class FixPour;

+ public:
  typedef DeviceType device_type;
  typedef ArrayTypes<DeviceType> AT;

- public:
  RegBlockKokkos(class LAMMPS *, int, char **);
  ~RegBlockKokkos();
  void match_all_kokkos(int, DAT::t_int_1d);
--- a/src/KSPACE/ewald.cpp
+++ b/src/KSPACE/ewald.cpp
@ -40,14 +40,20 @@ using namespace MathConst;

 /* ---------------------------------------------------------------------- */

-Ewald::Ewald(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
+Ewald::Ewald(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg),
+  kxvecs(NULL), kyvecs(NULL), kzvecs(NULL), ug(NULL), eg(NULL), vg(NULL),
+  ek(NULL), sfacrl(NULL), sfacim(NULL), sfacrl_all(NULL), sfacim_all(NULL),
+  cs(NULL), sn(NULL), sfacrl_A(NULL), sfacim_A(NULL), sfacrl_A_all(NULL),
+  sfacim_A_all(NULL), sfacrl_B(NULL), sfacim_B(NULL), sfacrl_B_all(NULL),
+  sfacim_B_all(NULL)
 {
+  group_allocate_flag = 0;
+  kmax_created = 0;
  if (narg != 1) error->all(FLERR,"Illegal kspace_style ewald command");

  ewaldflag = 1;
  group_group_enable = 1;
-  group_allocate_flag = 0;
-
+  
  accuracy_relative = fabs(force->numeric(FLERR,arg[0]));

  kmax = 0;
--- a/src/KSPACE/ewald_disp.cpp
+++ b/src/KSPACE/ewald_disp.cpp
@ -45,7 +45,9 @@ enum{GEOMETRIC,ARITHMETIC,SIXTHPOWER};   // same as in pair.h

 /* ---------------------------------------------------------------------- */

-EwaldDisp::EwaldDisp(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
+EwaldDisp::EwaldDisp(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg),
+  kenergy(NULL), kvirial(NULL), energy_self_peratom(NULL), virial_self_peratom(NULL),
+  ekr_local(NULL), hvec(NULL), kvec(NULL), B(NULL), cek_local(NULL), cek_global(NULL)
 {
  if (narg!=1) error->all(FLERR,"Illegal kspace_style ewald/n command");

--- a/src/KSPACE/msm.cpp
+++ b/src/KSPACE/msm.cpp
@ -44,7 +44,19 @@ enum{REVERSE_RHO,REVERSE_AD,REVERSE_AD_PERATOM};
 enum{FORWARD_RHO,FORWARD_AD,FORWARD_AD_PERATOM};
 /* ---------------------------------------------------------------------- */

-MSM::MSM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
+MSM::MSM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg),
+  factors(NULL), delxinv(NULL), delyinv(NULL), delzinv(NULL), nx_msm(NULL),
+  ny_msm(NULL), nz_msm(NULL), nxlo_in(NULL), nylo_in(NULL), nzlo_in(NULL),
+  nxhi_in(NULL), nyhi_in(NULL), nzhi_in(NULL), nxlo_out(NULL), nylo_out(NULL),
+  nzlo_out(NULL), nxhi_out(NULL), nyhi_out(NULL), nzhi_out(NULL), ngrid(NULL),
+  active_flag(NULL), alpha(NULL), betax(NULL), betay(NULL), betaz(NULL), peratom_allocate_flag(0),
+  levels(0), world_levels(NULL), qgrid(NULL), egrid(NULL), v0grid(NULL), v1grid(NULL),
+  v2grid(NULL), v3grid(NULL), v4grid(NULL), v5grid(NULL), g_direct(NULL),
+  v0_direct(NULL), v1_direct(NULL), v2_direct(NULL), v3_direct(NULL), v4_direct(NULL),
+  v5_direct(NULL), g_direct_top(NULL), v0_direct_top(NULL), v1_direct_top(NULL),
+  v2_direct_top(NULL), v3_direct_top(NULL), v4_direct_top(NULL), v5_direct_top(NULL),
+  phi1d(NULL), dphi1d(NULL), procneigh_levels(NULL), cg(NULL), cg_peratom(NULL),
+  cg_all(NULL), cg_peratom_all(NULL), part2grid(NULL), boxlo(NULL)
 {
  if (narg < 1) error->all(FLERR,"Illegal kspace_style msm command");

@ -96,8 +108,6 @@ MSM::MSM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
  egrid = NULL;
  v0grid = v1grid = v2grid = v3grid = v4grid = v5grid = NULL;

-  levels = 0;
-
  peratom_allocate_flag = 0;
  scalar_pressure_flag = 1;
  warn_nonneutral = 0;
--- a/src/KSPACE/msm_cg.cpp
+++ b/src/KSPACE/msm_cg.cpp
@ -42,7 +42,8 @@ enum{FORWARD_RHO,FORWARD_AD,FORWARD_AD_PERATOM};

 /* ---------------------------------------------------------------------- */

-MSMCG::MSMCG(LAMMPS *lmp, int narg, char **arg) : MSM(lmp, narg, arg)
+MSMCG::MSMCG(LAMMPS *lmp, int narg, char **arg) : MSM(lmp, narg, arg),
+  is_charged(NULL)
 {
  if ((narg < 1) || (narg > 2))
    error->all(FLERR,"Illegal kspace_style msm/cg command");
@ -53,7 +54,6 @@ MSMCG::MSMCG(LAMMPS *lmp, int narg, char **arg) : MSM(lmp, narg, arg)
  else smallq = SMALLQ;

  num_charged = -1;
-  is_charged = NULL;
 }

 /* ----------------------------------------------------------------------
--- a/src/KSPACE/pppm.cpp
+++ b/src/KSPACE/pppm.cpp
@ -64,8 +64,20 @@ enum{FORWARD_IK,FORWARD_AD,FORWARD_IK_PERATOM,FORWARD_AD_PERATOM};

 /* ---------------------------------------------------------------------- */

-PPPM::PPPM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
+PPPM::PPPM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg), 
+  factors(NULL), density_brick(NULL), vdx_brick(NULL), vdy_brick(NULL), vdz_brick(NULL),
+  u_brick(NULL), v0_brick(NULL), v1_brick(NULL), v2_brick(NULL), v3_brick(NULL),
+  v4_brick(NULL), v5_brick(NULL), greensfn(NULL), vg(NULL), fkx(NULL), fky(NULL),
+  fkz(NULL), density_fft(NULL), work1(NULL), work2(NULL), gf_b(NULL), rho1d(NULL),
+  rho_coeff(NULL), drho1d(NULL), drho_coeff(NULL), sf_precoeff1(NULL), sf_precoeff2(NULL),
+  sf_precoeff3(NULL), sf_precoeff4(NULL), sf_precoeff5(NULL), sf_precoeff6(NULL),
+  acons(NULL), density_A_brick(NULL), density_B_brick(NULL), density_A_fft(NULL),
+  density_B_fft(NULL), fft1(NULL), fft2(NULL), remap(NULL), cg(NULL), cg_peratom(NULL),
+  part2grid(NULL), boxlo(NULL)
 {
+  peratom_allocate_flag = 0;
+  group_allocate_flag = 0;
+
  if (narg < 1) error->all(FLERR,"Illegal kspace_style pppm command");

  pppmflag = 1;
@ -108,9 +120,6 @@ PPPM::PPPM(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
  nmax = 0;
  part2grid = NULL;

-  peratom_allocate_flag = 0;
-  group_allocate_flag = 0;
-
  // define acons coefficients for estimation of kspace errors
  // see JCP 109, pg 7698 for derivation of coefficients
  // higher order coefficients may be computed if needed
--- a/src/KSPACE/pppm_cg.cpp
+++ b/src/KSPACE/pppm_cg.cpp
@ -48,7 +48,8 @@ enum{FORWARD_IK,FORWARD_AD,FORWARD_IK_PERATOM,FORWARD_AD_PERATOM};

 /* ---------------------------------------------------------------------- */

-PPPMCG::PPPMCG(LAMMPS *lmp, int narg, char **arg) : PPPM(lmp, narg, arg)
+PPPMCG::PPPMCG(LAMMPS *lmp, int narg, char **arg) : PPPM(lmp, narg, arg), 
+  is_charged(NULL)
 {
  if ((narg < 1) || (narg > 2))
    error->all(FLERR,"Illegal kspace_style pppm/cg command");
@ -59,7 +60,6 @@ PPPMCG::PPPMCG(LAMMPS *lmp, int narg, char **arg) : PPPM(lmp, narg, arg)
  else smallq = SMALLQ;

  num_charged = -1;
-  is_charged = NULL;
  group_group_enable = 1;
 }

--- a/src/KSPACE/pppm_disp.cpp
+++ b/src/KSPACE/pppm_disp.cpp
@ -64,7 +64,48 @@ enum{FORWARD_IK, FORWARD_AD, FORWARD_IK_PERATOM, FORWARD_AD_PERATOM,

 /* ---------------------------------------------------------------------- */

-PPPMDisp::PPPMDisp(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg)
+PPPMDisp::PPPMDisp(LAMMPS *lmp, int narg, char **arg) : KSpace(lmp, narg, arg),
+  factors(NULL), csumi(NULL), cii(NULL), B(NULL), density_brick(NULL), vdx_brick(NULL), 
+  vdy_brick(NULL), vdz_brick(NULL), density_fft(NULL), u_brick(NULL), v0_brick(NULL), 
+  v1_brick(NULL), v2_brick(NULL), v3_brick(NULL), v4_brick(NULL), v5_brick(NULL), 
+  density_brick_g(NULL), vdx_brick_g(NULL), vdy_brick_g(NULL), vdz_brick_g(NULL), 
+  density_fft_g(NULL), u_brick_g(NULL), v0_brick_g(NULL), v1_brick_g(NULL), v2_brick_g(NULL), 
+  v3_brick_g(NULL), v4_brick_g(NULL), v5_brick_g(NULL), density_brick_a0(NULL), 
+  vdx_brick_a0(NULL), vdy_brick_a0(NULL), vdz_brick_a0(NULL), density_fft_a0(NULL), 
+  u_brick_a0(NULL), v0_brick_a0(NULL), v1_brick_a0(NULL), v2_brick_a0(NULL), 
+  v3_brick_a0(NULL), v4_brick_a0(NULL), v5_brick_a0(NULL), density_brick_a1(NULL), 
+  vdx_brick_a1(NULL), vdy_brick_a1(NULL), vdz_brick_a1(NULL), density_fft_a1(NULL), 
+  u_brick_a1(NULL), v0_brick_a1(NULL), v1_brick_a1(NULL), v2_brick_a1(NULL), 
+  v3_brick_a1(NULL), v4_brick_a1(NULL), v5_brick_a1(NULL), density_brick_a2(NULL), 
+  vdx_brick_a2(NULL), vdy_brick_a2(NULL), vdz_brick_a2(NULL), density_fft_a2(NULL), 
+  u_brick_a2(NULL), v0_brick_a2(NULL), v1_brick_a2(NULL), v2_brick_a2(NULL), 
+  v3_brick_a2(NULL), v4_brick_a2(NULL), v5_brick_a2(NULL), density_brick_a3(NULL), 
+  vdx_brick_a3(NULL), vdy_brick_a3(NULL), vdz_brick_a3(NULL), density_fft_a3(NULL), 
+  u_brick_a3(NULL), v0_brick_a3(NULL), v1_brick_a3(NULL), v2_brick_a3(NULL), 
+  v3_brick_a3(NULL), v4_brick_a3(NULL), v5_brick_a3(NULL), density_brick_a4(NULL), 
+  vdx_brick_a4(NULL), vdy_brick_a4(NULL), vdz_brick_a4(NULL), density_fft_a4(NULL), 
+  u_brick_a4(NULL), v0_brick_a4(NULL), v1_brick_a4(NULL), v2_brick_a4(NULL), 
+  v3_brick_a4(NULL), v4_brick_a4(NULL), v5_brick_a4(NULL), density_brick_a5(NULL), 
+  vdx_brick_a5(NULL), vdy_brick_a5(NULL), vdz_brick_a5(NULL), density_fft_a5(NULL), 
+  u_brick_a5(NULL), v0_brick_a5(NULL), v1_brick_a5(NULL), v2_brick_a5(NULL), 
+  v3_brick_a5(NULL), v4_brick_a5(NULL), v5_brick_a5(NULL), density_brick_a6(NULL), 
+  vdx_brick_a6(NULL), vdy_brick_a6(NULL), vdz_brick_a6(NULL), density_fft_a6(NULL), 
+  u_brick_a6(NULL), v0_brick_a6(NULL), v1_brick_a6(NULL), v2_brick_a6(NULL), 
+  v3_brick_a6(NULL), v4_brick_a6(NULL), v5_brick_a6(NULL), density_brick_none(NULL), 
+  vdx_brick_none(NULL), vdy_brick_none(NULL), vdz_brick_none(NULL), 
+  density_fft_none(NULL), u_brick_none(NULL), v0_brick_none(NULL), v1_brick_none(NULL), 
+  v2_brick_none(NULL), v3_brick_none(NULL), v4_brick_none(NULL), v5_brick_none(NULL), 
+  greensfn(NULL), vg(NULL), vg2(NULL), greensfn_6(NULL), vg_6(NULL), vg2_6(NULL), 
+  fkx(NULL), fky(NULL), fkz(NULL), fkx2(NULL), fky2(NULL), fkz2(NULL), fkx_6(NULL), 
+  fky_6(NULL), fkz_6(NULL), fkx2_6(NULL), fky2_6(NULL), fkz2_6(NULL), gf_b(NULL), 
+  gf_b_6(NULL), sf_precoeff1(NULL), sf_precoeff2(NULL), sf_precoeff3(NULL), 
+  sf_precoeff4(NULL), sf_precoeff5(NULL), sf_precoeff6(NULL), sf_precoeff1_6(NULL), 
+  sf_precoeff2_6(NULL), sf_precoeff3_6(NULL), sf_precoeff4_6(NULL), sf_precoeff5_6(NULL), 
+  sf_precoeff6_6(NULL), rho1d(NULL), rho_coeff(NULL), drho1d(NULL), drho_coeff(NULL), 
+  rho1d_6(NULL), rho_coeff_6(NULL), drho1d_6(NULL), drho_coeff_6(NULL), work1(NULL),
+   work2(NULL), work1_6(NULL), work2_6(NULL), fft1(NULL), fft2(NULL), fft1_6(NULL), 
+   fft2_6(NULL), remap(NULL), remap_6(NULL), cg(NULL), cg_peratom(NULL), cg_6(NULL), 
+   cg_peratom_6(NULL), part2grid(NULL), part2grid_6(NULL), boxlo(NULL)
 {
  if (narg < 1) error->all(FLERR,"Illegal kspace_style pppm/disp command");

--- a/src/KSPACE/pppm_stagger.cpp
+++ b/src/KSPACE/pppm_stagger.cpp
@ -52,7 +52,8 @@ enum{FORWARD_IK,FORWARD_AD,FORWARD_IK_PERATOM,FORWARD_AD_PERATOM};
 /* ---------------------------------------------------------------------- */

 PPPMStagger::PPPMStagger(LAMMPS *lmp, int narg, char **arg) :
-  PPPM(lmp, narg, arg)
+  PPPM(lmp, narg, arg), 
+  gf_b2(NULL)
 {
  if (narg < 1) error->all(FLERR,"Illegal kspace_style pppm/stagger command");
  stagger_flag = 1;
--- a/src/MAKE/OPTIONS/Makefile.icc_mpich
+++ b/src/MAKE/OPTIONS/Makefile.icc_mpich
@ -7,7 +7,7 @@ SHELL = /bin/sh
 # specify flags and libraries needed for your compiler

 CC =		mpicxx -cxx=icc
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -restrict
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

--- a/src/MAKE/OPTIONS/Makefile.icc_mpich_link
+++ b/src/MAKE/OPTIONS/Makefile.icc_mpich_link
@ -7,7 +7,7 @@ SHELL = /bin/sh
 # specify flags and libraries needed for your compiler

 CC =		icc
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -restrict
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

--- a/src/MAKE/OPTIONS/Makefile.icc_openmpi
+++ b/src/MAKE/OPTIONS/Makefile.icc_openmpi
@ -8,7 +8,7 @@ SHELL = /bin/sh

 export OMPI_CXX = icc
 CC =		mpicxx
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -restrict
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

--- a/src/MAKE/OPTIONS/Makefile.icc_openmpi_link
+++ b/src/MAKE/OPTIONS/Makefile.icc_openmpi_link
@ -7,7 +7,7 @@ SHELL = /bin/sh
 # specify flags and libraries needed for your compiler

 CC =		icc
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -restrict
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

--- a/src/MAKE/OPTIONS/Makefile.icc_serial
+++ b/src/MAKE/OPTIONS/Makefile.icc_serial
@ -7,7 +7,7 @@ SHELL = /bin/sh
 # specify flags and libraries needed for your compiler

 CC =		icc
-CCFLAGS =	-g -O3
+CCFLAGS =	-g -O3 -restrict
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

--- a/src/MAKE/OPTIONS/Makefile.kokkos_phi
+++ b/src/MAKE/OPTIONS/Makefile.kokkos_phi
@ -6,12 +6,12 @@ SHELL = /bin/sh
 # compiler/linker settings
 # specify flags and libraries needed for your compiler

-CC =		mpiicpc
+CC =		mpicxx
 CCFLAGS =	-g -O3
 SHFLAGS =	-fPIC
 DEPFLAGS =	-M

-LINK =		mpiicpc
+LINK =		mpicxx
 LINKFLAGS =	-g -O3
 LIB = 
 SIZE =		size
@ -20,7 +20,7 @@ ARCHIVE =	ar
 ARFLAGS =	-rc
 SHLIBFLAGS =	-shared
 KOKKOS_DEVICES = OpenMP
-KOKKOS_ARCH = KNC
+KOKKOS_ARCH = KNL

 # ---------------------------------------------------------------------
 # LAMMPS-specific settings, all OPTIONAL
--- a/src/MANYBODY/fix_qeq_comb.cpp
+++ b/src/MANYBODY/fix_qeq_comb.cpp
@ -39,7 +39,7 @@ using namespace FixConst;
 /* ---------------------------------------------------------------------- */

 FixQEQComb::FixQEQComb(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg),
-fp(NULL), qf(NULL), q1(NULL), q2(NULL)
+  fp(NULL), comb(NULL), comb3(NULL), qf(NULL), q1(NULL), q2(NULL)
 {
  if (narg < 5) error->all(FLERR,"Illegal fix qeq/comb command");

--- a/src/MANYBODY/pair_tersoff.h
+++ b/src/MANYBODY/pair_tersoff.h
@ -49,6 +49,7 @@ class PairTersoff : public Pair {
    double ZBLcut,ZBLexpscale;
    double c5,ca1,ca4;           // added for TersoffMOD
    double powern_del;
+    double c0;                   // added for TersoffMODC
  };

  Param *params;                // parameter set for an I-J-K interaction
--- a/src/MANYBODY/pair_tersoff_mod.h
+++ b/src/MANYBODY/pair_tersoff_mod.h
@ -30,7 +30,7 @@ class PairTersoffMOD : public PairTersoff {
  ~PairTersoffMOD() {}

 protected:
-  void read_file(char *);
+  virtual void read_file(char *);
  virtual void setup_params();
  double zeta(Param *, double, double, double *, double *);

--- a/src/MANYBODY/pair_tersoff_mod_c.cpp
+++ b/src/MANYBODY/pair_tersoff_mod_c.cpp
@ -0,0 +1,196 @@
+/* ----------------------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+/* ----------------------------------------------------------------------
+   Contributing author: Ganga P Purja Pun (George Mason University, Fairfax)
+------------------------------------------------------------------------- */
+
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include "pair_tersoff_mod_c.h"
+#include "atom.h"
+#include "neighbor.h"
+#include "neigh_list.h"
+#include "neigh_request.h"
+#include "force.h"
+#include "comm.h"
+#include "memory.h"
+#include "error.h"
+
+#include "math_const.h"
+
+using namespace LAMMPS_NS;
+using namespace MathConst;
+
+#define MAXLINE 1024
+#define DELTA 4
+
+/* ---------------------------------------------------------------------- */
+
+void PairTersoffMODC::read_file(char *file)
+{
+  int params_per_line = 21;
+  char **words = new char*[params_per_line+1];
+
+  memory->sfree(params);
+  params = NULL;
+  nparams = maxparam = 0;
+
+  // open file on proc 0
+
+  FILE *fp;
+  if (comm->me == 0) {
+    fp = force->open_potential(file);
+    if (fp == NULL) {
+      char str[128];
+      sprintf(str,"Cannot open Tersoff potential file %s",file);
+      error->one(FLERR,str);
+    }
+  }
+
+  // read each line out of file, skipping blank lines or leading '#'
+  // store line of params if all 3 element tags are in element list
+
+  int n,nwords,ielement,jelement,kelement;
+  char line[MAXLINE],*ptr;
+  int eof = 0;
+
+  while (1) {
+    if (comm->me == 0) {
+      ptr = fgets(line,MAXLINE,fp);
+      if (ptr == NULL) {
+	    eof = 1;
+	    fclose(fp);
+      } else n = strlen(line) + 1;
+    }
+    MPI_Bcast(&eof,1,MPI_INT,0,world);
+    if (eof) break;
+    MPI_Bcast(&n,1,MPI_INT,0,world);
+    MPI_Bcast(line,n,MPI_CHAR,0,world);
+
+    // strip comment, skip line if blank
+
+    if ((ptr = strchr(line,'#'))) *ptr = '\0';
+    nwords = atom->count_words(line);
+    if (nwords == 0) continue;
+
+    // concatenate additional lines until have params_per_line words
+
+    while (nwords < params_per_line) {
+      n = strlen(line);
+      if (comm->me == 0) {
+        ptr = fgets(&line[n],MAXLINE-n,fp);
+        if (ptr == NULL) {
+	      eof = 1;
+	      fclose(fp);
+        } else n = strlen(line) + 1;
+      }
+      MPI_Bcast(&eof,1,MPI_INT,0,world);
+      if (eof) break;
+      MPI_Bcast(&n,1,MPI_INT,0,world);
+      MPI_Bcast(line,n,MPI_CHAR,0,world);
+      if ((ptr = strchr(line,'#'))) *ptr = '\0';
+      nwords = atom->count_words(line);
+    }
+
+    if (nwords != params_per_line)
+      error->all(FLERR,"Incorrect format in Tersoff potential file");
+
+    // words = ptrs to all words in line
+
+    nwords = 0;
+    words[nwords++] = strtok(line," \t\n\r\f");
+    while ((words[nwords++] = strtok(NULL," \t\n\r\f"))) continue;
+
+    // ielement,jelement,kelement = 1st args
+    // if all 3 args are in element list, then parse this line
+    // else skip to next line
+
+    for (ielement = 0; ielement < nelements; ielement++)
+      if (strcmp(words[0],elements[ielement]) == 0) break;
+    if (ielement == nelements) continue;
+    for (jelement = 0; jelement < nelements; jelement++)
+      if (strcmp(words[1],elements[jelement]) == 0) break;
+    if (jelement == nelements) continue;
+    for (kelement = 0; kelement < nelements; kelement++)
+      if (strcmp(words[2],elements[kelement]) == 0) break;
+    if (kelement == nelements) continue;
+
+    // load up parameter settings and error check their values
+
+    if (nparams == maxparam) {
+      maxparam += DELTA;
+      params = (Param *) memory->srealloc(params,maxparam*sizeof(Param),
+					  "pair:params");
+    }
+
+    params[nparams].ielement = ielement;
+    params[nparams].jelement = jelement;
+    params[nparams].kelement = kelement;
+    params[nparams].powerm = atof(words[3]);
+    params[nparams].lam3 = atof(words[4]);
+    params[nparams].h = atof(words[5]);
+    params[nparams].powern = atof(words[6]);
+    params[nparams].beta = atof(words[7]);
+    params[nparams].lam2 = atof(words[8]);
+    params[nparams].bigb = atof(words[9]);
+    params[nparams].bigr = atof(words[10]);
+    params[nparams].bigd = atof(words[11]);
+    params[nparams].lam1 = atof(words[12]);
+    params[nparams].biga = atof(words[13]);
+    params[nparams].powern_del = atof(words[14]);
+    params[nparams].c1 = atof(words[15]);
+    params[nparams].c2 = atof(words[16]);
+    params[nparams].c3 = atof(words[17]);
+    params[nparams].c4 = atof(words[18]);
+    params[nparams].c5 = atof(words[19]);
+    params[nparams].c0 = atof(words[20]);
+
+    // currently only allow m exponent of 1
+
+    params[nparams].powermint = int(params[nparams].powerm);
+
+    if (
+	params[nparams].lam3 < 0.0 || params[nparams].powern < 0.0 ||
+	params[nparams].beta < 0.0 || params[nparams].lam2 < 0.0 ||
+	params[nparams].bigb < 0.0 || params[nparams].bigr < 0.0 ||
+	params[nparams].bigd < 0.0 ||
+                               params[nparams].bigd > params[nparams].bigr ||
+	params[nparams].lam3 < 0.0 || params[nparams].biga < 0.0 ||
+	params[nparams].powerm - params[nparams].powermint != 0.0 ||
+    (params[nparams].powermint != 3 && params[nparams].powermint != 1))
+      error->all(FLERR,"Illegal Tersoff parameter");
+
+    nparams++;
+  }
+
+  delete [] words;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void PairTersoffMODC::repulsive(Param *param, double rsq, double &fforce,
+                                int eflag, double &eng)
+{
+  double r,tmp_fc,tmp_fc_d,tmp_exp;
+
+  r = sqrt(rsq);
+  tmp_fc = ters_fc(r,param);
+  tmp_fc_d = ters_fc_d(r,param);
+  tmp_exp = exp(-param->lam1 * r);
+  fforce = -param->biga * tmp_exp * (tmp_fc_d - tmp_fc * param->lam1) / r - param->c0 * tmp_fc_d / r;
+  if (eflag) eng = tmp_fc * param->biga * tmp_exp + param->c0 * tmp_fc;
+}
+
--- a/src/MANYBODY/pair_tersoff_mod_c.h
+++ b/src/MANYBODY/pair_tersoff_mod_c.h
@ -0,0 +1,66 @@
+/* -*- c++ -*- ----------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   http://lammps.sandia.gov, Sandia National Laboratories
+   Steve Plimpton, sjplimp@sandia.gov
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+#ifdef PAIR_CLASS
+
+PairStyle(tersoff/mod/c,PairTersoffMODC)
+
+#else
+
+#ifndef LMP_PAIR_TERSOFF_MOD_C_H
+#define LMP_PAIR_TERSOFF_MOD_C_H
+
+#include "pair_tersoff_mod.h"
+
+namespace LAMMPS_NS {
+
+class PairTersoffMODC : public PairTersoffMOD {
+ public:
+  PairTersoffMODC(class LAMMPS *lmp) : PairTersoffMOD(lmp) {};
+  ~PairTersoffMODC() {}
+
+ protected:
+  void read_file(char *);
+  void repulsive(Param *, double, double &, int, double &);
+};
+
+}
+
+#endif
+#endif
+
+/* ERROR/WARNING messages:
+
+E: Cannot open Tersoff potential file %s
+
+The specified potential file cannot be opened.  Check that the path
+and name are correct.
+
+E: Incorrect format in Tersoff potential file
+
+Incorrect number of words per line in the potential file.
+
+E: Illegal Tersoff parameter
+
+One or more of the coefficients defined in the potential file is
+invalid.
+
+E: Potential file has duplicate entry
+
+The potential file has more than one entry for the same element.
+
+E: Potential file is missing an entry
+
+The potential file does not have a needed entry.
+
+*/
--- a/src/MANYBODY/pair_vashishta.cpp
+++ b/src/MANYBODY/pair_vashishta.cpp
@ -64,6 +64,8 @@ PairVashishta::PairVashishta(LAMMPS *lmp) : Pair(lmp)

 PairVashishta::~PairVashishta()
 {
+  if (copymode) return;
+  
  if (elements)
    for (int i = 0; i < nelements; i++) delete [] elements[i];
  delete [] elements;
--- a/src/MANYBODY/pair_vashishta.h
+++ b/src/MANYBODY/pair_vashishta.h
@ -34,7 +34,6 @@ class PairVashishta : public Pair {
  double init_one(int, int);
  void init_style();

- protected:
  struct Param {
    double bigb,gamma,r0,bigc,costheta;
    double bigh,eta,zi,zj;
@ -45,7 +44,7 @@ class PairVashishta : public Pair {
    double lam1rc,lam4rc,vrcc2,vrcc3,vrc,dvrc,c0;
    int ielement,jelement,kelement;
  };
-
+ protected:
  double cutmax;                // max cutoff for all elements
  int nelements;                // # of unique elements
  char **elements;              // names of unique elements
--- a/src/MC/fix_atom_swap.cpp
+++ b/src/MC/fix_atom_swap.cpp
@ -57,7 +57,7 @@ using namespace MathConst;

 FixAtomSwap::FixAtomSwap(LAMMPS *lmp, int narg, char **arg) :
  Fix(lmp, narg, arg),
-  idregion(NULL), type_list(NULL), mu(NULL), qtype(NULL), sqrt_mass_ratio(NULL), random_equal(NULL), random_unequal(NULL)
+  idregion(NULL), type_list(NULL), mu(NULL), qtype(NULL), sqrt_mass_ratio(NULL), local_swap_iatom_list(NULL), local_swap_jatom_list(NULL), local_swap_atom_list(NULL), random_equal(NULL), random_unequal(NULL), c_pe(NULL)
 {
  if (narg < 10) error->all(FLERR,"Illegal fix atom/swap command");

--- a/src/MC/fix_bond_break.cpp
+++ b/src/MC/fix_bond_break.cpp
@ -37,7 +37,7 @@ using namespace FixConst;

 FixBondBreak::FixBondBreak(LAMMPS *lmp, int narg, char **arg) :
  Fix(lmp, narg, arg),
-  partner(NULL), finalpartner(NULL), distsq(NULL), broken(NULL), copy(NULL), random(NULL)
+  partner(NULL), finalpartner(NULL), distsq(NULL), probability(NULL), broken(NULL), copy(NULL), random(NULL)
 {
  if (narg < 6) error->all(FLERR,"Illegal fix bond/break command");

--- a/src/MC/fix_bond_create.cpp
+++ b/src/MC/fix_bond_create.cpp
@ -40,8 +40,8 @@ using namespace FixConst;

 FixBondCreate::FixBondCreate(LAMMPS *lmp, int narg, char **arg) :
  Fix(lmp, narg, arg),
-  bondcount(NULL), partner(NULL), finalpartner(NULL), distsq(NULL), 
-  created(NULL), copy(NULL), random(NULL)
+  bondcount(NULL), partner(NULL), finalpartner(NULL), distsq(NULL),
+  probability(NULL), created(NULL), copy(NULL), random(NULL), list(NULL)
 {
  if (narg < 8) error->all(FLERR,"Illegal fix bond/create command");

--- a/src/MC/fix_bond_swap.cpp
+++ b/src/MC/fix_bond_swap.cpp
@ -53,7 +53,8 @@ static const char cite_fix_bond_swap[] =

 FixBondSwap::FixBondSwap(LAMMPS *lmp, int narg, char **arg) :
  Fix(lmp, narg, arg),
-  tflag(0), alist(NULL), id_temp(NULL)
+  tflag(0), alist(NULL), id_temp(NULL), type(NULL), x(NULL), list(NULL),
+  temperature(NULL), random(NULL)
 {
  if (lmp->citeme) lmp->citeme->add(cite_fix_bond_swap);

--- a/src/MC/fix_tfmc.cpp
+++ b/src/MC/fix_tfmc.cpp
@ -38,7 +38,7 @@ using namespace FixConst;

 FixTFMC::FixTFMC(LAMMPS *lmp, int narg, char **arg) :
  Fix(lmp, narg, arg),
-  xd(NULL), rotflag(0)
+  xd(NULL), rotflag(0), random_num(NULL)
 {
  if (narg < 6) error->all(FLERR,"Illegal fix tfmc command");

--- a/src/MISC/compute_ti.cpp
+++ b/src/MISC/compute_ti.cpp
@ -36,7 +36,8 @@ enum{PAIR,TAIL,KSPACE};
 /* ---------------------------------------------------------------------- */

 ComputeTI::ComputeTI(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg), nterms(0), which(NULL), ivar1(NULL), ivar2(NULL),
+  ilo(NULL), ihi(NULL), var1(NULL), var2(NULL), pptr(NULL), pstyle(NULL)
 {
  if (narg < 4) error->all(FLERR,"Illegal compute ti command");

--- a/src/MISC/dump_xtc.cpp
+++ b/src/MISC/dump_xtc.cpp
@ -53,7 +53,8 @@ int xdr3dfcoord(XDR *, float *, int *, float *);

 /* ---------------------------------------------------------------------- */

-DumpXTC::DumpXTC(LAMMPS *lmp, int narg, char **arg) : Dump(lmp, narg, arg)
+DumpXTC::DumpXTC(LAMMPS *lmp, int narg, char **arg) : Dump(lmp, narg, arg), 
+  coords(NULL)
 {
  if (narg != 5) error->all(FLERR,"Illegal dump xtc command");
  if (binary || compressed || multifile || multiproc)
--- a/src/MISC/fix_deposit.cpp
+++ b/src/MISC/fix_deposit.cpp
@ -45,7 +45,9 @@ enum{DIST_UNIFORM,DIST_GAUSSIAN};
 /* ---------------------------------------------------------------------- */

 FixDeposit::FixDeposit(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), idregion(NULL), idrigid(NULL),
+  idshake(NULL), onemols(NULL), molfrac(NULL), coords(NULL), imageflags(NULL),
+  fixrigid(NULL), fixshake(NULL), random(NULL)
 {
  if (narg < 7) error->all(FLERR,"Illegal fix deposit command");

--- a/src/MISC/fix_efield.cpp
+++ b/src/MISC/fix_efield.cpp
@ -41,7 +41,8 @@ enum{NONE,CONSTANT,EQUAL,ATOM};
 /* ---------------------------------------------------------------------- */

 FixEfield::FixEfield(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), xstr(NULL), ystr(NULL), zstr(NULL),
+  estr(NULL), idregion(NULL), efield(NULL)
 {
  if (narg < 6) error->all(FLERR,"Illegal fix efield command");

--- a/src/MISC/fix_evaporate.cpp
+++ b/src/MISC/fix_evaporate.cpp
@ -35,7 +35,7 @@ using namespace FixConst;
 /* ---------------------------------------------------------------------- */

 FixEvaporate::FixEvaporate(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), idregion(NULL), list(NULL), mark(NULL), random(NULL)
 {
  if (narg < 7) error->all(FLERR,"Illegal fix evaporate command");

--- a/src/MISC/fix_gld.cpp
+++ b/src/MISC/fix_gld.cpp
@ -43,7 +43,8 @@ using namespace FixConst;
 ------------------------------------------------------------------------- */

 FixGLD::FixGLD(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), 
+  step_respa(NULL), prony_c(NULL), prony_tau(NULL), s_gld(NULL), random(NULL)
 {
  int narg_min = 8;
  // Check to make sure we have the minimal number of inputs
--- a/src/MISC/fix_orient_bcc.cpp
+++ b/src/MISC/fix_orient_bcc.cpp
@ -57,7 +57,8 @@ static const char cite_fix_orient_bcc[] =
 /* ---------------------------------------------------------------------- */

 FixOrientBCC::FixOrientBCC(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), 
+  xifilename(NULL), chifilename(NULL), order(NULL), nbr(NULL), sort(NULL), list(NULL)
 {
  if (lmp->citeme) lmp->citeme->add(cite_fix_orient_bcc);

--- a/src/MISC/fix_orient_fcc.cpp
+++ b/src/MISC/fix_orient_fcc.cpp
@ -54,7 +54,9 @@ static const char cite_fix_orient_fcc[] =
 /* ---------------------------------------------------------------------- */

 FixOrientFCC::FixOrientFCC(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg),
+  xifilename(NULL), chifilename(NULL), order(NULL), nbr(NULL), 
+  sort(NULL), list(NULL)
 {
  if (lmp->citeme) lmp->citeme->add(cite_fix_orient_fcc);

--- a/src/MISC/fix_thermal_conductivity.cpp
+++ b/src/MISC/fix_thermal_conductivity.cpp
@ -36,7 +36,8 @@ using namespace FixConst;

 FixThermalConductivity::FixThermalConductivity(LAMMPS *lmp,
                                               int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg),
+  index_lo(NULL), index_hi(NULL), ke_lo(NULL), ke_hi(NULL)
 {
  if (narg < 6) error->all(FLERR,"Illegal fix thermal/conductivity command");

--- a/src/MISC/fix_ttm.cpp
+++ b/src/MISC/fix_ttm.cpp
@ -40,7 +40,12 @@ using namespace FixConst;
 /* ---------------------------------------------------------------------- */

 FixTTM::FixTTM(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg), 
+  random(NULL), fp(NULL), fpr(NULL), nsum(NULL), nsum_all(NULL), 
+  T_initial_set(NULL), gfactor1(NULL), gfactor2(NULL), ratio(NULL), 
+  flangevin(NULL), T_electron(NULL), T_electron_old(NULL), sum_vsq(NULL), 
+  sum_mass_vsq(NULL), sum_vsq_all(NULL), sum_mass_vsq_all(NULL), 
+  net_energy_transfer(NULL), net_energy_transfer_all(NULL)
 {
  if (narg < 15) error->all(FLERR,"Illegal fix ttm command");

--- a/src/MISC/fix_viscosity.cpp
+++ b/src/MISC/fix_viscosity.cpp
@ -37,7 +37,8 @@ using namespace FixConst;
 /* ---------------------------------------------------------------------- */

 FixViscosity::FixViscosity(LAMMPS *lmp, int narg, char **arg) :
-  Fix(lmp, narg, arg)
+  Fix(lmp, narg, arg),
+  pos_index(NULL), neg_index(NULL), pos_delta(NULL), neg_delta(NULL)
 {
  if (narg < 7) error->all(FLERR,"Illegal fix viscosity command");

--- a/src/MOLECULE/fix_cmap.cpp
+++ b/src/MOLECULE/fix_cmap.cpp
@ -62,7 +62,10 @@ using namespace MathConst;

 /* ---------------------------------------------------------------------- */

-FixCMAP::FixCMAP(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
+FixCMAP::FixCMAP(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg),
+  crosstermlist(NULL), num_crossterm(NULL), crossterm_type(NULL), crossterm_atom1(NULL),
+  crossterm_atom2(NULL), crossterm_atom3(NULL), crossterm_atom4(NULL), crossterm_atom5(NULL),
+  g_axis(NULL), cmapgrid(NULL), d1cmapgrid(NULL), d2cmapgrid(NULL), d12cmapgrid(NULL)
 {
  if (narg != 4) error->all(FLERR,"Illegal fix cmap command");

--- a/src/PERI/compute_damage_atom.cpp
+++ b/src/PERI/compute_damage_atom.cpp
@ -32,7 +32,7 @@ using namespace LAMMPS_NS;
 /* ---------------------------------------------------------------------- */

 ComputeDamageAtom::ComputeDamageAtom(LAMMPS *lmp, int narg, char **arg) :
-  Compute(lmp, narg, arg)
+  Compute(lmp, narg, arg), damage(NULL)
 {
  if (narg != 3) error->all(FLERR,"Illegal compute damage/atom command");

@ -40,7 +40,6 @@ ComputeDamageAtom::ComputeDamageAtom(LAMMPS *lmp, int narg, char **arg) :
  size_peratom_cols = 0;

  nmax = 0;
-  damage = NULL;
 }

 /* ---------------------------------------------------------------------- */
--- a/Show More
+++ b/Show More
Author	SHA1	Message	Date
Steve Plimpton	4b51719e67	new 17Nov16 patch and stable	2016-11-17 16:51:35 -07:00
Steve Plimpton	25d7be5f3d	compute pressure doc change	2016-11-17 16:11:30 -07:00
sjplimp	74c0e4dd5c	Merge pull request #278 from akohlmey/pair-agni Implementation of the AGNI manybody potential	2016-11-17 09:04:31 -07:00
sjplimp	073e8a0524	Merge pull request #276 from akohlmey/doc-updates Small bugfixes and updates	2016-11-17 09:02:27 -07:00
sjplimp	5320bbf585	Merge pull request #275 from andeplane/IP_VORONOI Initializing pointers in VORONOI	2016-11-17 09:01:46 -07:00
sjplimp	4448819824	Merge pull request #274 from andeplane/IP_POEMS Initialize pointers in POEMS	2016-11-17 09:01:37 -07:00
sjplimp	300ac30332	Merge pull request #273 from akohlmey/auto-memalign turn on -DLAMMPS_MEMALIGN=64 automatically when USER-INTEL is installed	2016-11-17 09:01:27 -07:00
sjplimp	2535e44991	Merge pull request #271 from akohlmey/tersoff-modc pair style tersoff/mod/c	2016-11-17 09:00:52 -07:00
sjplimp	cdae794383	Merge pull request #242 from andeplane/vashishta_kokkos Added KOKKOS vashishta	2016-11-17 09:00:03 -07:00
Steve Plimpton	8756a1017d	Kokkos updates by Stan	2016-11-17 08:58:22 -07:00
Axel Kohlmeyer	22d93fe8fb	add restrict to CCFLAGS for makefiles intended for intel compilers	2016-11-16 14:31:10 -05:00
Axel Kohlmeyer	57124b9b25	update documentation metadata files for recent changes	2016-11-12 09:18:21 -05:00
Axel Kohlmeyer	03b3834fe3	add documentation for pair style agni	2016-11-12 09:07:42 -05:00
Axel Kohlmeyer	d0124eac95	optimized data access and using approximate exponential for USER-OMP version	2016-11-12 08:36:27 -05:00
Axel Kohlmeyer	5685131fe2	add USER-OMP version of pair style agni	2016-11-11 19:08:01 -05:00
Axel Kohlmeyer	22fc92f9d8	use special ev_tally() function suitable for this kind of force compute	2016-11-11 18:32:55 -05:00
Axel Kohlmeyer	b9770766a8	add adatom and vacancy examples for AGNI pair style	2016-11-11 18:19:29 -05:00
Axel Kohlmeyer	9cc0c8badd	error exit when requested element is not in potential file	2016-11-11 18:10:57 -05:00
Axel Kohlmeyer	6e1492a86c	add potential file for pair style AGNI	2016-11-11 18:10:26 -05:00
Axel Kohlmeyer	9b0987d8c4	first complete implementation of AGNI pair style	2016-11-11 17:32:47 -05:00
Axel Kohlmeyer	e453adaf81	implemented parser for 1 element potential files	2016-11-11 15:53:37 -05:00
Axel Kohlmeyer	8e0fd88697	add example demonstrating the use of fix addtorque	2016-11-11 10:52:50 -05:00
Axel Kohlmeyer	fdcabd7d1d	fix addtorque is compatible with dynamic groups	2016-11-11 09:01:18 -05:00
Axel Kohlmeyer	c5c8c50e97	initialize 'nper'	2016-11-11 07:47:42 -05:00
Axel Kohlmeyer	72b0841b28	Merge branch 'doc-updates' of github.com:akohlmey/lammps into doc-updates	2016-11-11 07:44:44 -05:00
Axel Kohlmeyer	801111a7ab	dummy framework implementation for AGNI pair style	2016-11-10 15:00:36 -05:00
Axel Kohlmeyer	bfc478c320	simpler variant of the segfault workaround, that does not offend Clang c++.	2016-11-10 14:12:02 -05:00
Axel Kohlmeyer	2b75ee761d	avoid segmentation fault, when creating a LAMMPS instance from the library interface. arg[] may be NULL.	2016-11-10 13:28:32 -05:00
Axel Kohlmeyer	352e177fcd	Merge branch 'master' into small-fixes-and-updates	2016-11-10 13:26:00 -05:00
Anders Hafreager	c20ee34c7b	Initializing pointers in VORONOI	2016-11-10 09:30:07 +01:00
Anders Hafreager	95a7f7160e	Initialize pointers in POEMS	2016-11-10 09:26:52 +01:00
Axel Kohlmeyer	1f38e1a771	Merge branch 'master' into doc-updates	2016-11-09 18:18:05 -05:00
Steve Plimpton	9806da69f3	Stan bug fixes for fix reaxc/bonds/kk	2016-11-09 15:47:50 -07:00
Axel Kohlmeyer	fec87c070d	simplify compiling USER-INTEL package, by defaulting to -DLAMMPS_MEMALIGN=64 in case it is not set	2016-11-09 16:53:39 -05:00
Axel Kohlmeyer	3d3a99c082	added missing potential for tersoff/mod/c	2016-11-09 16:50:34 -05:00
Axel Kohlmeyer	3e36ec3754	remove unused class member	2016-11-09 16:17:46 -05:00
Axel Kohlmeyer	9ed5c4f0fa	Merge branch 'master' into kokkos-vashishta	2016-11-09 15:15:10 -05:00
Axel Kohlmeyer	c55fd502e0	correct typo in formula	2016-11-09 15:04:24 -05:00
Axel Kohlmeyer	71ee2ecaa1	integrate pair style tersoff/mod/c contributed by Ganga P Purja Pun (GMU) This includes docs, added testing and inclusion of USER-OMP support.	2016-11-09 14:52:39 -05:00
sjplimp	bfea3dce7d	Merge pull request #268 from arielzn/born_dsf pair styles born/coul/dsf and born/coul/dsf/cs added	2016-11-09 11:57:09 -07:00
sjplimp	eef862ee1c	Merge pull request #267 from akohlmey/pager-help use pager for help message, if connected to stdout	2016-11-09 11:53:43 -07:00
sjplimp	0cc2fbf1d6	Merge pull request #266 from andeplane/IP_USER_OMP Initializing pointers in USER-OMP	2016-11-09 11:52:05 -07:00
sjplimp	ae00666994	Merge pull request #265 from andeplane/IP_DIFF_DPD Initializing pointers in USER-DIFFRACTION and USER-DPD	2016-11-09 11:51:53 -07:00
sjplimp	51b3b5fb35	Merge pull request #264 from andeplane/IP_SNAP_SRD Initialize pointers in SNAP and SRD	2016-11-09 11:51:44 -07:00
sjplimp	176f2c3aa1	Merge pull request #263 from andeplane/IP_RIGID_SHOCK Initialize pointers in RIGID and SHOCK	2016-11-09 11:51:36 -07:00
sjplimp	3f71bfb185	Merge pull request #262 from andeplane/IP_PERI_QEQ_REPLICA Initializing pointers in PERI, QEQ and REPLICA	2016-11-09 11:51:25 -07:00
sjplimp	cf3ab51679	Merge pull request #261 from andeplane/IP_MISC_MOLECULE Initialized pointers in MISC and MOLECULE	2016-11-09 11:51:09 -07:00
sjplimp	59922f894b	Merge pull request #260 from andeplane/IP_MANYBODY_MC Initialize pointers in MANYBODY and MC	2016-11-09 11:51:01 -07:00
sjplimp	5e2b9d8bf3	Merge pull request #259 from andeplane/IP_KSPACE Initialize pointers in KSPACE	2016-11-09 11:50:50 -07:00
sjplimp	2d132cad6b	Merge pull request #258 from andeplane/IP_GRANULAR Initialize pointers in GRANULAR	2016-11-09 11:50:41 -07:00
sjplimp	ef6801f8bf	Merge pull request #257 from andeplane/IP_CORESHELL Initialize pointers in CORESHELL	2016-11-09 11:50:31 -07:00
sjplimp	c81a723642	Merge pull request #256 from andeplane/IP_BODY Initialized pointers in BODY	2016-11-09 11:50:23 -07:00
sjplimp	f9eb2a99ce	Merge pull request #255 from andeplane/IP_ASPHERE Initialize pointers in ASPHERE	2016-11-09 11:50:11 -07:00
sjplimp	16a02ef27d	Merge pull request #254 from andeplane/IP_root Initialized pointers in src folder	2016-11-09 11:47:21 -07:00
Steve Plimpton	2c801320c2	fixed links in Section_intro.txt	2016-11-09 11:46:10 -07:00
Anders Hafreager	d20b32092e	Building correct shortlists and removed rsq test in force loops	2016-11-08 18:57:27 +01:00
arielzn	9de1a2a08f	added input using born/coul/dsf/cs to examples/coreshell	2016-11-08 18:27:44 +01:00
arielzn	cdb5d47e9f	add FLERR argument to force->bounds() in born/coul/dsf	2016-11-08 18:24:54 +01:00
Steve Plimpton	a23b287a7a	sync with SVN for creation of stable release	2016-11-08 09:05:50 -07:00
Steve Plimpton	31204aab6a	sync with SVN	2016-11-08 08:57:51 -07:00
arielzn	25e7d074cf	documentation added for born/coul/dsf and born/coul/dsf/cs styles	2016-11-08 16:51:54 +01:00
arielzn	667f4dfe28	pair style born/coul/dsf added with its coreshell version	2016-11-08 11:32:38 +01:00
Axel Kohlmeyer	21694ca3a8	improve help and it through a pager, when screen == stdout	2016-11-07 17:10:12 -05:00
Axel Kohlmeyer	9b910d5511	make name of the actual executable (i.e. arg[0]) accessible	2016-11-07 17:07:40 -05:00
Anders Hafreager	054ab6bff3	Initializing pointers in USER-OMP	2016-11-07 21:07:10 +01:00
Anders Hafreager	616420cda8	Initializing pointers in USER-DIFFRACTION and USER-DPD	2016-11-07 20:51:36 +01:00
Anders Hafreager	fb3ac9afba	Initialize pointers in SNAP and SRD	2016-11-07 20:30:38 +01:00
Anders Hafreager	7cd7cda2d4	Initialize pointers in RIGID and SHOCK	2016-11-07 20:22:04 +01:00
Anders Hafreager	db0524278a	Initializing pointers in PERI, QEQ and REPLICA	2016-11-07 19:58:27 +01:00
Anders Hafreager	1ff75eaba2	Initialized pointers in MISC and MOLECULE	2016-11-07 17:15:48 +01:00
Anders Hafreager	30dede867a	Initialize pointers in MANYBODY and MC	2016-11-07 17:02:32 +01:00
Anders Hafreager	a5c6104d64	Initialize pointers in KSPACE	2016-11-07 16:54:59 +01:00
Anders Hafreager	c5869bdee2	Initialize pointers in GRANULAR	2016-11-07 16:33:50 +01:00
Anders Hafreager	e7a2c6b5d1	Initialize pointers in CORESHELL	2016-11-07 16:26:16 +01:00
Anders Hafreager	06959a9c59	Initialized pointers in BODY	2016-11-07 16:21:08 +01:00
Anders Hafreager	cd65d44d95	Initialize pointers in ASPHERE	2016-11-07 16:15:39 +01:00
Anders Hafreager	45f2e86dd6	NULLed ptrs in files	2016-11-07 16:07:37 +01:00
Anders Hafreager	f8226e8ae5	NULL ptrs in dump_custom and dump_image	2016-11-07 15:56:47 +01:00
Anders Hafreager	b221b15d24	NULLing ptrs in comm_brick and dump*.cpp	2016-11-07 15:50:18 +01:00
Axel Kohlmeyer	3a3d96b877	info styles also prints out pair styles	2016-11-04 18:18:40 -04:00
Anders Hafreager	f333d659c2	Using short neighborlists in neigh full	2016-10-29 22:54:43 +02:00
Anders Hafreager	51e2313fac	Using short neighbor lists	2016-10-29 22:35:29 +02:00
Anders Hafreager	e37d2b5c94	Calculating short neighbor lists	2016-10-29 22:20:37 +02:00
Axel Kohlmeyer	3870780894	Merge branch 'master' into kokkos-vashishta	2016-10-28 10:43:05 -04:00
stamoor	21619f6a2f	Recommitting reverted change git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@15794 f3b2605a-c512-4ea7-a41b-209d697bcdaa (cherry picked from commit `c0b98f5299`)	2016-10-25 15:25:39 -04:00
Anders Hafreager	039bda9b61	Added updated vashishta for KOKKOS support (cherry picked from commit 96089a42547f625e70aa2ac3933d248d2731b731)	2016-10-25 15:07:10 -04:00
Anders Hafreager	6929603eef	Added KOKKOS vashishta (cherry picked from commit 5edc474bf0be574ddba96d00bb63894edf400ddb)	2016-10-25 15:07:10 -04:00