Counting particle-particle collisions

Hi Felix,

yes - see compute contact/atom and compute coord/atom for details

Cheers
Christoph

I want to count the number of particle-particle collisions over all timesteps.
As far as I can see compute col contact/atom allows to track the number of collisions per timestep. However the output c_col jumps back to zero again as soon as the particles loose contact.
Is there an easy way to access the overall number of collisions per particle over the whole calculation time?

The overall number is not stored. You can modify the compute_contact_atom to sum every contact for each particle in a fix/property/atom in a way like this:
in constructor:

updFix = NULL;
updFix = static_cast(modify->find_fix_id("fppacc"));
if (updFix == NULL) error->all(FLERR,"did not find fppacc");

note: you need to add #include "fix_property_atom.h" in cpp file and to add class FixPropertyAtom* updFix; in compute_contact_atom.h to access the FixPropertyAtom Class
in ComputeContactAtom::compute_peratom() method access the new updFix like this:

...
if (rsq <= radsumsq) {
contact[i] += 1.0;
contact[j] += 1.0;
updFix->vector_atom[i] += 1; //sum contacts
updFix->vector_atom[j] += 1; //sum contacts
}

Now recompile liggghts !!

On scriptlevel you need an additional fix fppacc all property/atom fppacc scalar yes yes yes 0 in front of compute cc all contact/atom because we did not created it on c++ level (which is easiely possible).
To dump the new counter to paraview add f_fppacc to your dump custom command.
dump dmp all custom 1000 post/dump*.liggghts id type x y z vx vy vz fx fy fz omegax omegay omegaz radius c_cc f_fppacc #adds current #of contacts and sum of all contacts to dumpfile


good luck,
Christian.

ups, the forum software has eaten the "greater than" and "smaller than" symbols in the code tag block
2nd line should be:
updFix = static_cast SMALLERTHAN FixPropertyAtom*GREATERTHAN(modify->find_fix_id("fppacc"));

SMALLERTHAN <
GRATERTHAN >

Dear all,

I would like to count the overall number of collisions between different types of particles.
For example, I have particles type 1 and particle type 2, I want to know how many times
during my simulation a 1 particle collides with a 2 particle. Is there a way that I can do this?

Best,
Charlie

Hey Charlie,

just add an if (atom->type[i]!=atom->type[j]) in front of if (rsq <= radsumsq) {
to only count contact of different atom-types.

Best,
Christian.

you need to recompile liggghts src after every code change

How did you installed liggghts? Apt get or make Ubuntu?
2nd one: CD where you cloned liggghts CD Liggghts public / src Type make clean-all and than make Ubuntu or make Fedora
Than check if you symlinked from usr bin to the correct lmp_whatever executable

Sent from a mobile so some path should be uppercase..

I guess you are still starting the package liggghts version via global /usr/bin and not your self compiled one
You can see a build Date of LIGGGHTS in the first line of log
To identify your spezial Version modify liggghts_branch_version.txt
Remove all instances you get from whereis liggghts and create a new symlink from your spezial Version in /usr/bin

I have one more question as it pertains to counting the collisions between different species of particles, and I apologize for
being a bother and for my general lack of knowledge as it pertains to c++, liggghts, and ubuntu. I would like my output to give me
discrete values for the collisions between species. So if I have 3 particle types, I would like to know how many collisions between
particle type 1 and particle type 2, how many between particle type 1 and particle type 3, and how many between particle type 2
and particle type 3. Is there an easy way to implement this in the compute_contact_atom.cpp file so that when I dump to a csv file
I get a few columns that count the collisions in this manner rather than a single column with a sum of all inter-species collisions?

-Charles

easy way: duplicate ComputeContactAtom 3-time, rename to compute_contact_atom12.h and compute_contact_atom12.cpp (also 13,23)
in cpp and h file change all class-names (ComputeContactAtom::) to ComputeContactAtom12 ComputeContactAtom13 ComputeContactAtom23

change in every header file line ComputeStyle(contact/atom,ComputeContactAtom) to ComputeStyle(contact/atom12,ComputeContactAtom12)
change in cpp file line if (atom->type[i]!=atom->type[j]) to if (((atom->type[i]==0) && (atom->type[j]==1)) or ((atom->type[i]==1) && (atom->type[j]==0))) { //1-2 or 2-1, internal type starts at 0 !
recompile,
use compute contact/atom12,contact/atom13,contact/atom23 in your inputskript.

better, but more complicated: introduce two new parameters to original compute contact/atom to get wanted types from inputskript (sry, its to hard to explain how one would do this bc. you need to change the code at some other places in the cpp file)

good luck

Have a look in compute_rigid.cpp line 58-60, there the input skript value of "property" is read in, you have to adapt this logic to store your pair-combinations in two new private variables (let's call them type0_ and type1_)
Than write if (((atom->type[i]==type0_) && (atom->type[j]==type1_)) or ((atom->type[i]==type1_) && (atom->type[j]==type0_))) as condition.
Than you can use a loop in your inputskript:

variable type1 loop 1 1000
label type1loop
variable type2 loop 1 1000
label type2loop
compute ${type1}_${type2} contact/atom ${type1} ${type2}
next type2
jump SELF type2loop
next type1
jump SELF type1loop


The changes do work! But I a confused what each of the output stands for ? Is c_cc the collision count and f_fppacc is the collision frequency

Thank you!

Chai

Hi Chai,

no c_cc is compute collision count, fppacc is fix-property-particle-all-collision count, the prefix c_ and f_ is used in variable.cpp to check if the input is a fix or compute. The compute holds the collision for every atom at one timestep, while the fix stores the sum over all timesteps for all particles.

best,
Christian.

Hello

Is there a similar way to count the particle-wall collision also?

Thanks in advance

Regards
Naveen R

No. But have a look in mesh_module_stress.cpp in method add_particle_contribution(), it is easy to add an counter there with reference to previously introduced updFix = static_cast(modify->find_fix_id("fppacc"));
If you want to dump the sum of particle contacts of a mesh you would need a new mesh->prop which holds the counter. (have a look how wear is stored and dumped for mesh elements).

I did the mentioned changes in compute_contact_atom.cpp file. I did the following changes in compute_contact_atom.cpp

1. included the: #include "fix_property_atom.h"
2. defined the class: class FixPropertyAtom* updFix;
3. Added the three lines as
if(narg > 3)
{
if (narg < 5) error->all(FLERR,"Illegal compute contact/atom command");
if(strcmp("skin",arg[3])) error->all(FLERR,"Illegal compute contact/atom command, expecting keyword 'skin'");
skin = atof(arg[4]);
}

peratom_flag = 1;
size_peratom_cols = 0;
comm_reverse = 1;

nmax = 0;
contact = NULL;

updFix = NULL;
updFix = static_cast (modify->find_fix_id("fppacc"));
if (updFix == NULL) error->all(FLERR,"did not find fppacc");

4. delx = xtmp - x[j][0];
dely = ytmp - x[j][1];
delz = ztmp - x[j][2];
rsq = delx*delx + dely*dely + delz*delz;
radsum = radi + radius[j] + skin;
radsumsq = radsum*radsum;
if (rsq <= radsumsq) {
contact[i] += 1.0;
contact[j] += 1.0;
updFix->vector_atom[i] += 1; //sum of contacts
updFix->vector_atom[j] += 1; //sum of contacts
}

And following changes in file:
1. added:
fix fppacc all property/atom fppacc scalar yes yes yes 0
compute cc all contact/atom
2. dumping as:
dump dmp all custom 500 ../DEM/post1/dumpres*.liggghts id type type x y z c_2 ix iy iz vx vy vz fx fy fz omegax omegay omegaz radius f_fppacc c_cc

But I am getting different results everytime as I change the dumping interval time from 500 to 200 and so on. I am getting higher and unpractical values as I am decreasing the dumping interval.

Please suggest where did I make the mistake.

Thanks in advance.

Regards
Naveen R