LIGGGHTS Public v3.8.0 with JKR contact model

Bias's picture
Submitted by Bias on Thu, 08/22/2019 - 14:43

Hi all,
I decided to share the code for adhesive (JKR) contact models I have been working on for the previous months. The code also includes negative overlap after contact.
The implementation is based on the article "Discrete-element modeling of particulate aerosol flow" by J.S. Marshall 2009 in Journal of Computational Physics 2009. For the normal model the contact force is calculated using equation 40 while for the tangential model equation 43 is used.

The files can be found here: https://github.com/eidevag/LIGGGHTS-PUBLIC-JKR

For running the normal and tangential jkr model the following files are needed:
PolynomialRootFinder.cpp
PolynomialRootFinder.h
jkr_lookup_table.cpp
jkr_lookup_table.h
normal_model_jkr.h
tangential_model_jkr_tan.h

The code is using a lookup table that is calculated before the simulation and then the contact force is interpolated. The number of points is hardcoded in the file "jkr_lookup_table.h" with an arbitrary but for most cases sufficient value of 3e5, this can be changed by the user if wanted. The resolution of lookup table is a value set in the input file. If an overlap will exceed the maximum overlap value in the lookup table the code will give the warning "Simulation outside LUT table, increase the resolution for speed! delta_delta_c = " and then start solving the fourth-order polynomial for each iteration which significantly slows down the simulation. Usually, this happens for an ill-conditioned case with an extreme overlap at start or a very large timestep.

To run the normal model the following properties needs to be set (See example in examples/LIGGGHTS/jkr_cases/)
youngsModulus peratomtype
poissonsRation peratomtype
coefficientRestitution peratomtypepair
workOfAdhesion peratomtypepair
resolutionJKR scalar

The coefficientRestitution is a heritage from the hertz model which the normal model is coded from. The workOfAdhesion is in J/m2 (for SI) and for similar materials is equal to twice the surface energy. The scalar resolutionJKR specifies the spacing for the lookup table. Lower values gives higher resolution but lower range of displacement values. Default recommended value is 1e-4.

To run the tangential model
coefficientFriction peratomtypepair
is also needed.

Word of caution: It should be noted that the JKR model has a zero contact force Fn at the equilibrium point with displacement a0. Therefore the normal model should be combined with caution with other tangential and rolling models.
For example, tangential model history and rolling model epsd both has a resistance to sliding/rolling that is proportional to sidata.Fn. This will cause a particle to free slide/roll at the equilibrium point which is most likely not what is intended! For rolling the cdt model is recommended instead since it calculates Fn=kn*delta_n instead of using sidata.Fn

Old related forum post:
https://www.cfdem.com/forums/how-add-new-contact-model

I hope someone can get some of this code and feel free to ask me if anything is uncertain!
Special thanks to Mohammadreza Tamadondar for help in the coding process.

Best regards
Tobias Eidevåg

Bias's picture

Bias | Sat, 08/24/2019 - 13:42

Thank you for pointing that out.

Yes I think you're right about the limitForce. It's a heritage from the Hertz model I based the code on.
I never used it and it should most likely not be used for this model.

Best regards
Tobias

Bias's picture

Bias | Thu, 02/06/2020 - 15:42

Hi,
Just wanted to update this thread that I also added a rolling model that models the rolling resistance due to adhesion.
https://github.com/eidevag/LIGGGHTS-PUBLIC-JKR

Here is also our latest article describing the contact models:
Eidevåg, T., Abrahamsson, P., Eng, M., & Rasmuson, A. (2019). Modeling of dry snow adhesion during normal impact with surfaces. Powder Technology, 361, 1081–1092. https://doi.org/10.1016/j.powtec.2019.10.085

Bias's picture

Bias | Mon, 02/10/2020 - 16:55

Hi mschramm,
No I don't mind just great if the code comes to use so feel free to add it to your fibers repo.

Syedabulhassan | Sat, 10/09/2021 - 19:54

Can you please post the syntax to run the tangential model
I am using pair_style gran model jkr tangential history but it is not working

Shakir | Thu, 02/02/2023 - 11:43

I'm working on a project regarding simulation of nanoindentation using JKR model, I wanted to simulate the adhesion of sample material on the indenter over the Nano indenter. Can I use this model in commercial software like Ansys or abaqcus?
Thank you.