[lammps-users] (no subject)

Dear mailinglist,

Currently I'm implementing a new compute style, which should calculate
some averages in a local environment of an atom. More concrete: within a
sphere of radius force->pair->cutforce averaged values will be
calculated. One of them is a local stress:
(virial[i][X] + virial[i][Y] + virial[i][Z])

The sekeleton (run through neighbors etc) of the new compute is ripped of
the compute_centro_atom.cpp. The virial is accessed from the compute

The problem which occurs is the following: The serial version of the
code works smoothly, but the || code leads to artificial structures (one
can clearly see the edges of the decomposition). The stress for the
ghost atoms is too high. An example:

stress_avg_noghostsinneighborhood=60280 (real units)

From within the same new compute style local temperature or density can
be calculated without problems. Also if I visualize the pressure/atom
from compute_stress_atom.cpp itself, there are no artifacts. There seems
to be a problem in my code with accessing properties of ghost atoms
calculated by another compute style.

My question to the developers is now: Is there anything, which I have to
take into account, when I want to access the pressure tensor (or virial)
from a new compute style, which does a summation over neighbor lists? Do
the ghost atoms contain the current values, after I call

If there are more detailled questions, here are the relevant(?)


void ComputeAvgAtom::init()
==> request neighbor lists as in compute_centro_atom.cpp


double ComputeAvgAtom::ComputeStress (int i)
  if (force->pair)
     stress += force->pair->vatom[i][X] + force->pair->vatom[i][Y] +
  ==> etc for all other force types
return stress


void ComputeAvgAtom::compute_peratom ()
==> neighbor list stuff etc
for (ii = 0; ii < inum; ii++) {
    i = ilist[ii];
    ==> init some vars
    stress = ComputeStress (i);
    if (mask[i] & groupbit) {
       jlist = firstneigh[i];
       jnum = numneigh[i];
       for (jj = 0; jj < jnum; jj++) {
          j = jlist[jj];
          if (j >= nall) j %= nall;
          rsq = SQ(x[i][X]-x[j][X]) + SQ(x[i][Y]-x[j][Y]) + SQ(x[i][Z]-x[j][Z]);
          if (rsq < rcut_pot_sq) {
             stress += ComputeStress (j);
       vector_per_atom[i][0] = stress;

If there are even more question, I could pm a snapshot of my repository.

Best regards:

The virial terms are computed by potentials, not by computes. So you
need to insure the potentials compute them on the proper
timesteps you need them. This is typically done by the code
that calls the compute. See fix ave/atom for example.

Ghost atoms will not know about any
per-atom virial info.