From 4dec943a5e0ecd9edc06fb04b34c08a36a90cdce Mon Sep 17 00:00:00 2001 From: sjplimp Date: Tue, 22 Mar 2016 13:44:36 +0000 Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14777 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- src/USER-TALLY/compute_pe_mol_tally.cpp | 129 ++++++++++++++++++++++++ src/USER-TALLY/compute_pe_mol_tally.h | 59 +++++++++++ 2 files changed, 188 insertions(+) create mode 100644 src/USER-TALLY/compute_pe_mol_tally.cpp create mode 100644 src/USER-TALLY/compute_pe_mol_tally.h diff --git a/src/USER-TALLY/compute_pe_mol_tally.cpp b/src/USER-TALLY/compute_pe_mol_tally.cpp new file mode 100644 index 0000000000..09ee04d57a --- /dev/null +++ b/src/USER-TALLY/compute_pe_mol_tally.cpp @@ -0,0 +1,129 @@ +/* ---------------------------------------------------------------------- + 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 +#include "compute_pe_mol_tally.h" +#include "atom.h" +#include "group.h" +#include "pair.h" +#include "update.h" +#include "memory.h" +#include "error.h" +#include "force.h" + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +ComputePEMolTally::ComputePEMolTally(LAMMPS *lmp, int narg, char **arg) : + Compute(lmp, narg, arg) +{ + if (narg < 4) error->all(FLERR,"Illegal compute pe/mol/tally command"); + + igroup2 = group->find(arg[3]); + if (igroup2 == -1) + error->all(FLERR,"Could not find compute pe/mol/tally second group ID"); + groupbit2 = group->bitmask[igroup2]; + + vector_flag = 1; + size_vector = 4; + timeflag = 1; + + extvector = 1; + peflag = 1; // we need Pair::ev_tally() to be run + + did_compute = invoked_vector = -1; + vector = new double[size_vector]; +} + +/* ---------------------------------------------------------------------- */ + +ComputePEMolTally::~ComputePEMolTally() +{ + if (force && force->pair) force->pair->del_tally_callback(this); + delete[] vector; +} + +/* ---------------------------------------------------------------------- */ + +void ComputePEMolTally::init() +{ + if (force->pair == NULL) + error->all(FLERR,"Trying to use compute pe/mol/tally with no pair style"); + else + force->pair->add_tally_callback(this); + + if (atom->molecule_flag == 0) + error->all(FLERR,"Compute pe/mol/tally requires molecule IDs."); + + if (force->pair->single_enable == 0 || force->pair->manybody_flag) + error->all(FLERR,"Compute pe/mol/tally used with incompatible pair style."); + + if ((comm->me == 0) && (force->bond || force->angle || force->dihedral + || force->improper || force->kspace)) + error->warning(FLERR,"Compute pe/mol/tally only called from pair style"); + + did_compute = -1; +} + + +/* ---------------------------------------------------------------------- */ +void ComputePEMolTally::pair_tally_callback(int i, int j, int nlocal, int newton, + double evdwl, double ecoul, double, + double, double, double) +{ + const int * const mask = atom->mask; + const tagint * const molid = atom->molecule; + + // do setup work that needs to be done only once per timestep + + if (did_compute != update->ntimestep) { + did_compute = update->ntimestep; + + etotal[0] = etotal[1] = etotal[2] = etotal[3] = 0.0; + } + + if ( ((mask[i] & groupbit) && (mask[j] & groupbit2)) + || ((mask[i] & groupbit2) && (mask[j] & groupbit)) ){ + + evdwl *= 0.5; ecoul *= 0.5; + if (newton || i < nlocal) { + if (molid[i] == molid[j]) { + etotal[0] += evdwl; etotal[1] += ecoul; + } else { + etotal[2] += evdwl; etotal[3] += ecoul; + } + } + if (newton || j < nlocal) { + if (molid[i] == molid[j]) { + etotal[0] += evdwl; etotal[1] += ecoul; + } else { + etotal[2] += evdwl; etotal[3] += ecoul; + } + } + } +} + +/* ---------------------------------------------------------------------- */ + +void ComputePEMolTally::compute_vector() +{ + invoked_vector = update->ntimestep; + if ((did_compute != invoked_vector) || (update->eflag_global != invoked_vector)) + error->all(FLERR,"Energy was not tallied on needed timestep"); + + // sum accumulated energies across procs + + MPI_Allreduce(etotal,vector,size_vector,MPI_DOUBLE,MPI_SUM,world); +} + diff --git a/src/USER-TALLY/compute_pe_mol_tally.h b/src/USER-TALLY/compute_pe_mol_tally.h new file mode 100644 index 0000000000..b2c5ffab7f --- /dev/null +++ b/src/USER-TALLY/compute_pe_mol_tally.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 COMPUTE_CLASS + +ComputeStyle(pe/mol/tally,ComputePEMolTally) + +#else + +#ifndef LMP_COMPUTE_PE_MOL_TALLY_H +#define LMP_COMPUTE_PE_MOL_TALLY_H + +#include "compute.h" + +namespace LAMMPS_NS { + +class ComputePEMolTally : public Compute { + + public: + ComputePEMolTally(class LAMMPS *, int, char **); + virtual ~ComputePEMolTally(); + + void init(); + void compute_vector(); + + void pair_tally_callback(int, int, int, int, + double, double, double, + double, double, double); + + private: + bigint did_compute; + int igroup2,groupbit2; + double etotal[4]; +}; + +} + +#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. + +*/