Homec4science
Diffusion lammps (lammm-master-deprecated) Browse

lammps/src/bond.cpp
856f78e50a84lammm-master-deprecated

Last Change
Raw File

bond.cpp
View Options

/* ----------------------------------------------------------------------
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 "string.h"
#include "bond.h"
#include "atom.h"
#include "comm.h"
#include "force.h"
#include "suffix.h"
#include "atom_masks.h"
#include "memory.h"
#include "error.h"
using namespace LAMMPS_NS;
/* -----------------------------------------------------------------------
set bond contribution to Vdwl energy to 0.0
a particular bond style can override this
------------------------------------------------------------------------- */
Bond::Bond(LAMMPS *lmp) : Pointers(lmp)
{
energy = 0.0;
writedata = 1;
allocated = 0;
suffix_flag = Suffix::NONE;
maxeatom = maxvatom = 0;
eatom = NULL;
vatom = NULL;
datamask = ALL_MASK;
datamask_ext = ALL_MASK;
}
/* ---------------------------------------------------------------------- */
Bond::~Bond()
{
memory->destroy(eatom);
memory->destroy(vatom);
}
/* ----------------------------------------------------------------------
check if all coeffs are set
------------------------------------------------------------------------- */
void Bond::init()
{
if (!allocated) error->all(FLERR,"Bond coeffs are not set");
for (int i = 1; i <= atom->nbondtypes; i++)
if (setflag[i] == 0) error->all(FLERR,"All bond coeffs are not set");
init_style();
}
/* ----------------------------------------------------------------------
setup for energy, virial computation
see integrate::ev_set() for values of eflag (0-3) and vflag (0-6)
------------------------------------------------------------------------- */
void Bond::ev_setup(int eflag, int vflag)
{
int i,n;
evflag = 1;
eflag_either = eflag;
eflag_global = eflag % 2;
eflag_atom = eflag / 2;
vflag_either = vflag;
vflag_global = vflag % 4;
vflag_atom = vflag / 4;
// reallocate per-atom arrays if necessary
if (eflag_atom && atom->nmax > maxeatom) {
maxeatom = atom->nmax;
memory->destroy(eatom);
memory->create(eatom,comm->nthreads*maxeatom,"bond:eatom");
}
if (vflag_atom && atom->nmax > maxvatom) {
maxvatom = atom->nmax;
memory->destroy(vatom);
memory->create(vatom,comm->nthreads*maxvatom,6,"bond:vatom");
}
// zero accumulators
if (eflag_global) energy = 0.0;
if (vflag_global) for (i = 0; i < 6; i++) virial[i] = 0.0;
if (eflag_atom) {
n = atom->nlocal;
if (force->newton_bond) n += atom->nghost;
for (i = 0; i < n; i++) eatom[i] = 0.0;
}
if (vflag_atom) {
n = atom->nlocal;
if (force->newton_bond) n += atom->nghost;
for (i = 0; i < n; i++) {
vatom[i][0] = 0.0;
vatom[i][1] = 0.0;
vatom[i][2] = 0.0;
vatom[i][3] = 0.0;
vatom[i][4] = 0.0;
vatom[i][5] = 0.0;
}
}
}
/* ----------------------------------------------------------------------
tally energy and virial into global and per-atom accumulators
------------------------------------------------------------------------- */
void Bond::ev_tally(int i, int j, int nlocal, int newton_bond,
double ebond, double fbond,
double delx, double dely, double delz)
{
double ebondhalf,v[6];
if (eflag_either) {
if (eflag_global) {
if (newton_bond) energy += ebond;
else {
ebondhalf = 0.5*ebond;
if (i < nlocal) energy += ebondhalf;
if (j < nlocal) energy += ebondhalf;
}
}
if (eflag_atom) {
ebondhalf = 0.5*ebond;
if (newton_bond || i < nlocal) eatom[i] += ebondhalf;
if (newton_bond || j < nlocal) eatom[j] += ebondhalf;
}
}
if (vflag_either) {
v[0] = delx*delx*fbond;
v[1] = dely*dely*fbond;
v[2] = delz*delz*fbond;
v[3] = delx*dely*fbond;
v[4] = delx*delz*fbond;
v[5] = dely*delz*fbond;
if (vflag_global) {
if (newton_bond) {
virial[0] += v[0];
virial[1] += v[1];
virial[2] += v[2];
virial[3] += v[3];
virial[4] += v[4];
virial[5] += v[5];
} else {
if (i < nlocal) {
virial[0] += 0.5*v[0];
virial[1] += 0.5*v[1];
virial[2] += 0.5*v[2];
virial[3] += 0.5*v[3];
virial[4] += 0.5*v[4];
virial[5] += 0.5*v[5];
}
if (j < nlocal) {
virial[0] += 0.5*v[0];
virial[1] += 0.5*v[1];
virial[2] += 0.5*v[2];
virial[3] += 0.5*v[3];
virial[4] += 0.5*v[4];
virial[5] += 0.5*v[5];
}
}
}
if (vflag_atom) {
if (newton_bond || i < nlocal) {
vatom[i][0] += 0.5*v[0];
vatom[i][1] += 0.5*v[1];
vatom[i][2] += 0.5*v[2];
vatom[i][3] += 0.5*v[3];
vatom[i][4] += 0.5*v[4];
vatom[i][5] += 0.5*v[5];
}
if (newton_bond || j < nlocal) {
vatom[j][0] += 0.5*v[0];
vatom[j][1] += 0.5*v[1];
vatom[j][2] += 0.5*v[2];
vatom[j][3] += 0.5*v[3];
vatom[j][4] += 0.5*v[4];
vatom[j][5] += 0.5*v[5];
}
}
}
}
/* ---------------------------------------------------------------------- */
double Bond::memory_usage()
{
double bytes = comm->nthreads*maxeatom * sizeof(double);
bytes += comm->nthreads*maxvatom*6 * sizeof(double);
return bytes;
}
c4science · Help