add a compute bond that allows accessing per sub-style bond energies.

doc/compute_bond.txt

@@ -0,0 +1,55 @@
+"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
+
+compute bond command :h3
+
+[Syntax:]
+
+compute ID group-ID bond bstyle :pre
+
+ID, group-ID are documented in "compute"_compute.html command
+bond = style name of this compute command
+bstyle = style name of a bond style that calculates additional values :ul
+
+[Examples:]
+
+compute 1 all bond harmonic
+compute 2 all bond morse :pre
+
+[Description:]
+
+Define a computation that extracts additional values calculated by a
+bond style, sums them across processors, and makes them accessible for
+output or further processing by other commands.  The group specified
+for this command is ignored.
+
+The specified {bstyle} must be a bond style used in your simulation
+either by itself or as a sub-style in a "bond_style hybrid"_bond_hybrid.html
+command.
+
+This compute is useful when using "bond_style hybrid"_bond_hybrid.html
+if you want to know the portion of the total energy contributed by one
+sub-style.
+
+[Output info:]
+
+This compute calculates a global scalar which is the contribution of
+the named bond style to the bond energy {ebond}.
+
+The scalar value calculated by this compute is "extensive" and will be
+in energy "units"_units.html.
+
+[Restrictions:] none
+
+[Related commands:]
+
+"compute pe"_compute_pe.html, "compute pair"_compute_pair.html
+
+[Default:]
+
+The default for {evalue} is {ebond}.

doc/compute_pair.txt

@@ -76,7 +76,7 @@ see the doc page for the pair style for details.
 
 [Related commands:]
 
-"compute pe"_compute_pe.html
+"compute pe"_compute_pe.html, "compute bond"_compute_bond.html
 
 [Default:]
 

src/compute_bond.cpp

@@ -0,0 +1,89 @@
+/* ----------------------------------------------------------------------
+   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.
+------------------------------------------------------------------------- */
+
+#include <mpi.h>
+#include <string.h>
+#include "compute_bond.h"
+#include "update.h"
+#include "force.h"
+#include "bond.h"
+#include "error.h"
+
+using namespace LAMMPS_NS;
+
+/* ---------------------------------------------------------------------- */
+
+ComputeBond::ComputeBond(LAMMPS *lmp, int narg, char **arg) :
+  Compute(lmp, narg, arg)
+{
+  if (narg < 4 || narg > 5) error->all(FLERR,"Illegal compute bond command");
+  if (igroup) error->all(FLERR,"Compute bond must use group all");
+
+  scalar_flag = 1;
+  extscalar = 1;
+  peflag = 1;
+  timeflag = 1;
+
+  int n = strlen(arg[3]) + 1;
+  if (lmp->suffix) n += strlen(lmp->suffix) + 1;
+  bstyle = new char[n];
+  strcpy(bstyle,arg[3]);
+
+  // check if bond style with and without suffix exists
+
+  bond = force->bond_match(bstyle);
+  if (!bond && lmp->suffix) {
+    strcat(bstyle,"/");
+    strcat(bstyle,lmp->suffix);
+    bond = force->bond_match(bstyle);
+  }
+  if (!bond)
+    error->all(FLERR,"Unrecognized bond style in compute bond command");
+
+  vector = NULL;
+}
+
+/* ---------------------------------------------------------------------- */
+
+ComputeBond::~ComputeBond()
+{
+  delete [] bstyle;
+}
+
+/* ---------------------------------------------------------------------- */
+
+void ComputeBond::init()
+{
+  // recheck for bond style in case it has been deleted
+
+  bond = force->bond_match(bstyle);
+
+  if (!bond)
+    error->all(FLERR,"Unrecognized bond style in compute bond command");
+}
+
+/* ---------------------------------------------------------------------- */
+
+double ComputeBond::compute_scalar()
+{
+  invoked_scalar = update->ntimestep;
+  if (update->eflag_global != invoked_scalar)
+    error->all(FLERR,"Energy was not tallied on needed timestep");
+
+  double eng;
+  eng = bond->energy;
+
+  MPI_Allreduce(&eng,&scalar,1,MPI_DOUBLE,MPI_SUM,world);
+  return scalar;
+}
+

src/compute_bond.h

@@ -0,0 +1,67 @@
+/* -*- 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 COMPUTE_CLASS
+
+ComputeStyle(bond,ComputeBond)
+
+#else
+
+#ifndef LMP_COMPUTE_BOND_H
+#define LMP_COMPUTE_BOND_H
+
+#include "compute.h"
+
+namespace LAMMPS_NS {
+
+class ComputeBond : public Compute {
+ public:
+  ComputeBond(class LAMMPS *, int, char **);
+  ~ComputeBond();
+  void init();
+  double compute_scalar();
+
+ private:
+  int nbond;
+  char *bstyle;
+  class Bond *bond;
+};
+
+}
+
+#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: Compute bond must use group all
+
+Bond styles accumlate energy on all atoms.
+
+E: Unrecognized bond style in compute bond command
+
+Self-explanatory.
+
+E: Energy was not tallied on needed timestep
+
+You are using a thermo keyword that requires potentials to
+have tallied energy, but they didn't on this timestep.  See the
+variable doc page for ideas on how to make this work.
+
+*/