Hooke/Herz force and the problem of negative (attractive) forces

Dear LAMMPS Developers,

I am simulating frictionless granular spheres with

pair_style gran/hooke 1 0 1 0 0 0
pair_coeff * *

Then I dump distance and forces between pairs as follows:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force

dump dmp2 all local 10000 pairs/dump.pair.* index c_1[1] c_1[2] c_2[1] c_2[2]

In this model particles interact only if they overlap.

I have two questions:

  1. I get a lot of zero forces for particles that do not overlap. These pairs appear in the dump file because they are in one another’s neighbouring list.

e.g.:
1 74 135 1.30228 0
2 32 79 1.23184 0
3 32 148 1.28095 0
4 79 62 1.24551 0
5 79 155 1.1986 0
6 62 88 1.31448 0

Is there any flag in LAMMPS to tell it not to dump particles that do not overlap? This will save a lot of space when one uses millions of particles.

Re: your Q about discarding output when the particles are not

touching. See the “cutoff radius” option on the compute pair local

command. It should do what you want.

Re: negative forces - I have not seen this before. I’m CCing Dan

Bolintineanu who can comment.

Steve

If the particle velocity during unloading is large, the damping can cause a negative (attractive) force. This is not a LAMMPS implementation issue as much as an issue with the classic spring/dashpot model. We’ve usually found that for typical stiffness values, dynamics and damping set to near critical values, it’s a very minor effect. One solution (as implemented in LIGGGHTS) is to zero out the force if it becomes negative, which is just a single line to implement in pair/gran/hooke and hertz. Also, if you’re finding a lot of negative force values, it may be that your time step is too large, causing large overlaps->high force->high unloading velocities?

Dan

Dan S. Bolintineanu, Ph. D.

Fluid and Reactive Processes Department, 1516

Sandia National Laboratories

Ph: 505 284 5460

FYI, the suggested modification to the current granular models in LAMMPS so that there are no attractive forces due to the damping are now included in pull request #510 and should become part of the the next LAMMPS patch. https://github.com/lammps/lammps/pull/510

axel.

Sorry to throw a wrench in the machine of progress.

I personally do not see the possibility of having an attractive force in the linear spring dash-pot model as a bad thing - even if unphysical on a microscopic level. (note that even the linear spring dashpot model is unphysical, a point that T. Poeschel referenced by Habib has also been quite adamant about) As Dan says the linear model does not often substantially affect the kinematics or stress in most granular systems. Moreover, you lose the extremely nice feature of a constant coefficient of restitution. This constant CoR has allowed many people, including myself, to develop continuum theories and study the affect of combining models (collision time-scales, energy loss, aggregation rate, etc.) with granular models - such as van der waals cohesion and electrostatics. Studies on the constant coefficient of restitution systems then allow us to look for similar scalings in the hertz/kono-kuwabara models - also the zeroed out force alternatve.

I suggest an alternative - adding a flag in the pair style to decide whether or not the collision force (combined dashpot and spring forces) should be zeroed if they become negative.

I do agree with Eric

to give a possibility to users to opt in/out the attractive force in granular spring models.

http://newton.kias.re.kr/~habib/

Dear Steve,

Thanks for your reply,

Dear Steve,

Thanks for your reply,

Re: your Q about discarding output when the particles are not
touching. See the "cutoff radius" option on the compute pair local
command. It should do what you want.

I did try this by:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force cutoff radius

However I get the following error:

ERROR: Invalid keyword in compute pair/local command
(../compute_pair_local.cpp:64)

My lammps is: LAMMPS (17 Nov 2015)

​your version is too old. the "cutoff" option exists only since the
22July2016 version.
http://lammps.sandia.gov/bug2016.html

please always refer to the docs matching your specific version and not the
online docs (those always correspond to the very latest development
version).

axel.​

Thanks Alex for your helpful response!

I updated my LAMMPS with:

sudo apt-get install lammps-daily 

and now run my simulation with:

lammps-daily  -in  in.shearLE.2D.lammps -echo screen

with this compute:

compute 1 all property/local patom1 patom2

compute	2 all pair/local  dist  force cutoff radius

I get now the new following error:

Per MPI rank memory allocation (min/avg/max) = 11.14 | 11.14 | 11.14 Mbytes
Step Temp E_pair E_mol TotEng Press Volume 
       1 5.8491046e+17            0            0 0.00050068525 0.00047900632          100 
ERROR on proc 0: Dump local count is not consistent across input fields (../dump_local.cpp:310)

Thanks Alex for your helpful response!

I updated my LAMMPS with:

sudo apt-get install lammps-daily

and now run my simulation with:

lammps-daily -in in.shearLE.2D.lammps -echo screen

with this compute:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force cutoff radius

I get now the new following error:

Per MPI rank memory allocation (min/avg/max) = 11.14 | 11.14 | 11.14 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
       1 5.8491046e+17 0 0 0.00050068525 0.00047900632 100
ERROR on proc 0: Dump local count is not consistent across input fields (../dump_local.cpp:310)
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------

However when I remove the "cutoff radius" flag the program runs quite well:

i.e.:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force

Is this again a version problem?

​no, this is the expected behavior of dump local.​

axel

Thanks Alex for your helpful response!

I updated my LAMMPS with:

sudo apt-get install lammps-daily

and now run my simulation with:

lammps-daily -in in.shearLE.2D.lammps -echo screen

with this compute:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force cutoff radius

I get now the new following error:

Per MPI rank memory allocation (min/avg/max) = 11.14 | 11.14 | 11.14 Mbytes
Step Temp E_pair E_mol TotEng Press Volume
       1 5.8491046e+17 0 0 0.00050068525 0.00047900632 100
ERROR on proc 0: Dump local count is not consistent across input fields (../dump_local.cpp:310)
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.
--------------------------------------------------------------------------

However when I remove the "cutoff radius" flag the program runs quite well:

i.e.:

compute 1 all property/local patom1 patom2

compute 2 all pair/local dist force

Is this again a version problem?

​no, this is the expected behavior of dump local.​

axel

How can I solve this problem? I searched in the mailing list, could not
find any solution.

​don't use the cutoff option or use two different dump local commands to
output the different size local vectors into different files.

axel.​

Dear Axel,

Thanks for your kind help.

I use cutoff radius and now print into two separate files like:

dump dmpl2 all local ${dumpSteps} pairs/dump.pairList.* index c_1[1] c_1[2]

dump dmpl3 all local ${dumpSteps} pairs/dump.pairForce.* index c_2[1] c_2[2]

and it runs properly without termination and error.

Thanks again for your reply.

Yours sincerely,

Dear all,

in our last communication Axel wrote me that

Dear all,

in our last communication Axel wrote me that

FYI, the suggested modification to the current granular models in LAMMPS
so that there are no attractive forces due to the damping are now included
in pull request #510 and should become part of the the next LAMMPS patch.
https://github.com/lammps/lammps/pull/510

axel.

the change has been (so far) rejected pending further feedback from an
expert on granular media.

How can I get now this modification into my lammps?

Is there an automatic way to do that?

just issue the following commands:

git clone https://github.com/lammps/lammps.git
cd lammps/
git checkout -b granular-no-attractive-force master
git cherry-pick 066123007cd9277d4c03669707cf34ab6c05bc06

this will create a checkout of the current LAMMPS git repository, then
create a custom branch (so you don't mess up the master branch) called
granular-no-attractive-force.
into that branch the specific commit with the reverted change is
included selectively with git cherry-pick

now you have a correspondingly customize LAMMPS version and can
install the desired packages (none are installed by default) and
compile with the required flags.

axel.

p.s.: isn't git fantastic? nothing ever gets lost.