I’m looking for a lammps command which can do similar thing like what I described below:
The compute centro/atom command can be used to calculate the centro-symmetry parameter for each atom in the simulation box. Now let’s say that I want to calculate the average centro-symmetry parameter for each atom, which is defined by averaging the centro-symmetry parameters of the central atom and its 10 neighbors.
Is there a similar command in lammps? At the moment I’m trying to extent lammps to do such a thing. If there is a similar command, then I can look at its source and header files and change them to fulfill my requirement.
Thanks and look for your response.
As far as coding it the easiest place to start is probably from compute_centro_atom. The “centro” array pointer is the atom vector for the compute and you can probably reassign to its contents the average you desire at the end of the original calculation. I’m not sure if the ten neighbors you’re referring to is the variable “nhalf” limiting the loop.
Thank you very much for your quick response. In fact I’m aksing my friend who is good at coding but knows nothing about lammps to help me extent lammps. We first defined a new compute command which has a similar function like the compue centro/atom. This goal has been achieved. Let’s say the new compute command calculates a scalar value for each atom, and we call this per-atom scalar as “feature”. Then we change our code to calculate a new value for each atom. The new value (let’s call it “average_feature”) is the avarage of 11 “feature” values which come from the central atom and its 10 neighbors. However, this time problem occurred. The “average_feature” for each atom was totally wrong. My friend guesses that this is due to parallel computing, namely, sometimes a “average_feature” is calculated before the 11 "feature"s have been calculated. So we think maybe we can first define a new compute command to calculate the “feature”, and then define another compute command to calcualte the “average_feature”.
So our questions are: (1) Is the problem due to parallel computing? We realized that the parallel function is not defined in the source code of compute centro/atom command. It’s defined outside. (2) Is it possible to compute the “average_feature” through just one compute command?
Thanks again and best wishes,
There are a number of issues at hand here.
Parallel computing, or rather more accurately, domain decomposition is one of them. You need to understand the concept of “local” and “ghost” atoms in LAMMPS. computes generally provide computed properties only for the local atoms. but for your averaging over neighbors, you may need to average over atoms that are “ghost” atoms, i.e. copies of atoms from neighboring subdomains. for those you usually cannot access the computed data, since the per-atom property arrays are including only local atoms (i.e. those belonging to a sub-domain). to make that data accessible, you would have to increase the storage (to use atom->nmax in size, which is the maximum of atom->nlocal + atom->nghost across all processors) and then do what is called a “forward communication”, i.e. an exchange of data where properties from local atoms are sent to the corresponding storage of their “ghosts” in the neighboring subdomains.
once that is done, you can loop over the neighbors like before. however, if you want a fixed limit of 10 (closest) neighbors, you would have to loop over the neighbor list and create a sublist that sorts entries by distance and then pick the N closest. generally, it is easier to just pick a cutoff radius and then average over all atoms within that distance.
and to more specifically answer your second question. yes, you can do the averaging as a second step within the same compute, but it would actually be much more general to implement a generic fix ave/sphere that would do this averaging within a given spherical radius (or number of neighbors, if your prefer) and have this implemented as a generic feature, so it could be used with any compute, that computes a per-atom property.
Thank you very much for your responses. They really help. We’ll check these these further and might contact you again later this week.
Since it might take a lot of efforts to tackle with the domain decomposition problem. We are considering another way to do this.
(1) Step1: At the beginning of the simulation, we use “atom_style atomic” and perform the simulation, getting a dump file which contains the ID, x, y, z, ix, iy, iz of each atom.
(2) Step2: We use the “rerun command” to calculate the “feature” of each atom using the new “compute feature/atom” command which we add to lammps by modifying the original "compute centro/atom” command. We output the ID, x, y, z, ix, iy, iz, feature of each atom.
(3) Step3: We first change the dump file obtained from Step 2 by replacing all the “feature” words (which should have the form of c_ID) with “q”. Then we use “atom_style charge” and define a simulatio box with the same number of atoms as the dump file. lastly, we use the “rerun command” again. Similar to “compute centro/atom”, if we can pass the charge of atom into a new compute command which we define, we can compute the average_feature.
Is this possible? If it is, how to pass the change to a compute command? Looking forward to your response.
I second Axel’s last suggestion, which is to make this a new command, separate
from compute centro/atom. It could be a compute which takes any other
per-atom compute (e.g. compute centro/atom) as an input and averages the
values for atoms nearby (within a cutoff? or N nearest?), to calculate a
spatial average for each atom. The output of the compute could
also be used as input to fix ave/atom to perform further time-averaging of
the quantity. The new compute could do the needed ghost communication
that Axel also mentioned.
Note that compute chunk/atom and fix ave/chunk can do something similar
which is to compute an average over the atoms in a small bin (i.e. geometric
region), e.g. for any per-atom quantity (e.g. centro-symmetry).
Not exactly what you’re asking, but it is close enough to what you want.