Bonded Particle with Damping (Extending richti83's LIGGGHTS-With-Bonds)

Submitted by mschramm on Tue, 10/17/2017 - 20:56

Hello,
I'm posting what I have done to incorporate damping into a flexible fiber model.
I am using Yu Guo's paper "Granular Shear Flows of Flexible Rod-like Particles"
https://www.researchgate.net/publication/258741920_Granular_Shear_Flows_...
This release currently works for LIGGGHTS 3.7.0
The only thing that i would like to explore further is this bond model interaction with the fix_template_multiplespheres command.
As this would allow the user to insert directly into the simulation rather than read in the particle and bond information from a text file.
I will also be uploading a FORTRAN code that allows the insertion of flexible fibers into a cylindrical region.

Here is a link to richti83's original bond code
https://github.com/richti83/LIGGGHTS-WITH-BONDS

Here is my link
https://github.com/schrummy14/LIGGGHTS_Flexible_Fibers

EDIT 08/18/2020
I have dynamic insertion of bonds working for a single processor applications. There is currently an issue during information sharing across processors.
You can certainly uses this feature to first spawn in your fibers then write a restart file using a single processor. Then read in the restart file with multiple processors. Once this is working in multi-core applications, I will merge the branch to the master. Here is a teaser video.
https://www.youtube.com/watch?v=h5n0SNLrwes&feature=youtu.be

Forums:

mschramm | Thu, 09/03/2020 - 20:42

I have a flag for this that you can signal during the atom_style call.
See the 3-point-bending example to see how the flag is used.
The flag seems like it will only disable the normal contact force but it currently disables all contact forces between two bonded particles.
https://github.com/schrummy14/LIGGGHTS_Flexible_Fibers/blob/master/src/p...

embastias | Sat, 09/05/2020 - 19:12

Hi, I would like to know if this works from a restart file. I am generating a beam similar to the one in your example, but with elements of two different sizes and two types of bonds. When I create the bonds they start to break without applying any load, this happens with sigma_break and tau_break less than 15e7.

As I understand it independently of the resistance I use in the bonds, they should have 0 stress when created. To ensure this I have put the command:
set group all vx 0 vy 0 vz 0 omegax 0 omegay 0 omegaz 0

before creating the bonds.

I don't know if you can known why this might be happening.

Thanks for your help
Emmanuel Bastias

mschramm | Sun, 09/06/2020 - 19:40

Hello,
I just tried with my shear cell example and experienced no issues.
Added "disableNormalContact 1" to the atom style call.
Set the time step to 1.0e-7 (has the potential to be higher but I did not play with it).
Added 3 spheres to the multiplesphere template.
Reduced the bond skin.

Added the following calls after each run command to reduce particle movement until bonds are created.
set group all vx 0 vy 0 vz 0 omegax 0 omegay 0 omegaz 0

This is where I wish I have finished my dynamic insertion code for multiplesphere templates. For this exact example, It would work as you are doing insertion on a single core, writing a restart and finishing the simulation. Currently, the dynamic insertion code breaks in parallel usage.

embastias | Mon, 09/07/2020 - 05:10

Hi, I think I found my problem. I was watching if when I run the simulation the code is executed after the line 289
if (isBonded) continue;

write down the elements that this line of code is passed through and display them in paraview, I found that it was calculating the forces for elements that are far away, not in contact.

I don't know if this is a problem in the code or if I am missing any commands in my script to update the neighbors of each element.
I hope you can help me, since I don't know in which file to check which elements are neighbors to each other.

In my script I have put these commands
neighbor 0.001 bin
neigh_modify delay 0

and the diameters of the elements are 15.4 mm and 5mm.

thank you for your help

Emmanuel Bastias

embastias | Mon, 09/07/2020 - 06:05

it had nothing to do with the code it was my model that had elements in contact with each other without a bond between them.
thank you for your support with the above answers

mschramm | Mon, 09/07/2020 - 16:53

I should state that the flag isn't perfect as it only checks the particles that are directly bonded to a specific sphere.
So it works as intended if you are only using kissing spheres to assemble your bonded shape. But as soon as we
start to make a "smooth" particle we run into trouble as the contact equations start to be used again. Looking at
this again, I need to add another list to the atom class that holds the tag of an atom that the contact forces should
not be applied to.
Using the flag with spheres that are kissing should not have any issues.

benkla | Tue, 06/29/2021 - 15:59

Hello,
thanks a lot for sharing your code. I tried to build your code with MS VisualStudio 14 (2015) for a windows application, it stopped with error C2131. It seems MSVS doesn´t like variable-length arrays ;)

I had to change following parts:

jkr_lookup_table.cpp, Line 206: "std::vector <double> op(Degree+1),zeroi(Degree),zeror(Degree); // Coefficient vectors"

PolynominalRootFinter.cpp, Line 34: "void rpoly_ak1(std::vector<double> &op, int* Degree, std::vector<double> &zeror, std::vector<double> &zeroi){"
+ change definition in .h-File

Not fully tested yet but the executable works. Maybe that´s helpful to other users.

Best reguards,
Benedikt

mschramm | Wed, 06/30/2021 - 06:49

Hello,
You may want to try simply replacing the calls to creating op, zeroi, and zeror with
double * op = new double[Degree+1];
double * zeroi = new double[Degree];
double * zeror = new double[Degree];

This creates the variables in a dynamic way.

Don't forget the corresponding deletes at the end of the function.
delete [] op;
delete [] zeroi;
delete [] zeror;

This was tested using the jkr examples and found to have no change.
This also allows you to keep the originally written code in rpoly_ak1.

Pages