Dear all,

I am implementing a granular simulation, where the rattlers are identified and removed to test the force balance. After utilizing “compute rattlers/atom”, I label the non-rattler particles and group them together for later computing, as shown below

``````compute         rat all rattlers/atom radius 4 10000
run             0

variable        pp atom "c_rat == 0"
group           nonrattler variable pp

compute         f nonrattler pair/local force
compute         af nonrattler property/atom fx fy fz
``````

I have tested these commands with a sample containing 100 particles (with 99 being non-rattlers). I have confirmed that group nonrattler does contain 99 particles. However, using “compute af nonrattler …” and “aforce = lmp.numpy.extract_compute(‘af’, 1, 2)”, I will obtain an array with 100 elements, instead of 99 elements. Specifically, one element corresponding to the rattler particle is exactly zero vector.

Now I want to calculate the mean total force of non-rattler particles and mean force of neighbor pairs (should be non-rattler particles), and the termination for new packing generation is set to be

``````npair = lmp.extract_compute('f', 2, 0)
nonRat = lmp.get_natoms() - int(lmp.extract_compute('rat',0, 0))
``````

For new packing generation, some thing weird happened: when system is nearly jammed, both pair force and total force of particles within nonrattlers are zero. For example,

``````compute         af all property/atom fx fy fz
compute         f all pair/local force
then
average pair force: 9.600610516428262e-07 average total force: 1.500483357014579e-07

compute         af nonrattler property/atom fx fy fz
compute         f all pair/local force
then
average pair force: 1.0062173713140847e-06 average total force: 0.0

compute         af nonrattler property/atom fx fy fz
compute         f nonrattler pair/local force
then
npair ==0, average total force: 0.0
``````

Can you please specify what exactly you mean by “weird?” What behavior do you expect from LAMMPS and how does it differ from what LAMMPS returns? If you are below jamming, are all particles just rattlers?

It is also generally hard to give any detailed advice without information on your specific system and input script. Some advice for posts can be found here. For instance, I don’t know what pair style you are using. There is a caveat for granular pair styles on the compute pair/local doc page which may or may not be relevant.

Can you please specify what exactly you mean by “weird?” What behavior do you expect from LAMMPS and how does it differ from what LAMMPS returns? If you are below jamming, are all particles just rattlers?

All particles will be rattlers if system is below jamming. For an already jammed system (a sample), when I calculate within the nonrattler group containing 99 particles, I would expect the “compute atom” command give me an array with len(array) == 99, instead of 100 (the total particle numbers).

Also, I tried to generate new packing samples, during which I can monitor the average pair force and average total force of all particles. As shown in my previous post, when the system is nearly jammed, “compute nonrattler” gives zero average total force and no pair force, but “compute all” yields reasonable results. So I think it is weird.

It is also generally hard to give any detailed advice without information on your specific system and input script. Some advice for posts can be found here. For instance, I don’t know what pair style you are using. There is a caveat for granular pair styles on the compute pair/local doc page which may or may not be relevant.

The simulation concerns frictionless granular particles interacting with hooke law, i.e.,

``````pair_style      granular
pair_coeff      * * hooke 1 \${gamma} damping mass_velocity limit_damping tangential linear_history 1 1 0
``````

I have attached necessary files (one packing sample, one analysis script, and one packing-generation script, respectively). Thank you for your prompt response.
test.in (2.8 KB)
test.data (8.7 KB)
nph_mono.in (4.1 KB)

The length of per-atom output is based on the total number of atoms. Calculations are skipped for atoms outside of the group, or zero, as you noted.

Unfortunately, I still do not know what you mean by weird. How exactly does this observed behavior differ from what you expect? If all particles are rattlers, does that not explain why there are no forces on non-rattlers? Have you confirmed that there should be non-rattlers in your system? Does ‘nearly jammed’ mean below jamming?

Thank you for sharing some files, but the data file does not seem to be the system you are concerned with. It only has one rattler and non-zero forces between non-rattlers.