Hi all,
I am simulating a group of particle bonding together ( like a string). I divide a string (or fibre) into a number of balls. By this way, I need to define the bonding peroperties (kn; ks) among these balls so that they can be created and also be broken. I got a nice suggestion from Christian (Richti) in terms of granular bond package. Could anyone share with me more information about this package? how can I get and install it? Does it work on public version of LIGGGHTS? DEM simulation would be strongly limited without bonding.
Thanks,
richti83 | Mon, 12/01/2014 - 08:21
HowTo
You need these files:
https://github.com/Polyun/LIGGGHTS-PUBLIC/tree/master/bondspackage
add them to src folder.
I think you need to hack some error->{warning|all|error}("message") to error->{warning|all|error}(FLERR,"message") statements to compile it (bc. this code is for liggghts2 and not maintained anymore.
Be carefull in bond_gran.cpp line 286 the force-dissapation factor is set to 0.995, this will lead to loose forces over time, I would suggest to set it to 1 or develop a velocity-dependant damping function there ..
As far as I remember the code was not stable in multicore usage when atoms leaves the simulation domain, so make sure you have a bounding box of walls.
Good Luck,
Christian.
b.saprativ | Mon, 12/01/2014 - 08:31
DEM simulation would be strongly limited without bonding.
Great question..Even I would like to know more about bonded particles and effect of forces on those assembly of particles. Is it available with vanilla liggghts and how to use the fix command
richti83 | Mon, 12/01/2014 - 12:41
Literatur
A good starting point is Potyondya & Cundall "A bonded-particle model for rock" http://www.sciencedirect.com/science/article/pii/S1365160904002874, the implementation in the source files I linked above seems to be equal to this paper and so it is equal to the implementation in PFC.
Another Bond Model is the V-Model, unfurtunately I could not find any source files: http://tm.spbstu.ru/Implementation_of_V_model_in_LIGGGHTS
When you have successfull compiled the bond_gran package into LIGGGHTS you can use this example as starting point:
https://github.com/Polyun/LIGGGHTS-PUBLIC/blob/master/bondspackage/bonds...
Line 28 defines the kn/ks stiffness an normal/shear strength of one bond,
Line 65 establish bonds between atoms within the bond-cuttoff-distance (0.0051m)
Hope this helps.
When you are not familar with C++ you can use the old Liggghts version:
git clone https://github.com/Polyun/LIGGGHTS-PUBLIC.git LIGGGHTS-BONDS
cd LIGGGHTS-BONDS/src/
make fedora -j
I prefer to create a symlink in /usr/bin/ to start liggghts-with-bonds from everywhere:
cd /usr/bin/
sudo ln -s /home/YOURUSERNAME/LIGGGHTS-BONDS/src/lmp_fedora lmp4bonds
Then
cd $HOME
cd LIGGGHTS-BONDS/bondspackage/bondsfalloncone/
lmp4bonds -in in.template
Good Luck
Christian
apple | Sat, 04/22/2017 - 11:27
Could you provide me with the bond_gran package?
Thanks very much! I have read your reply. I think the link"https://github.com/Polyun/LIGGGHTS-PUBLIC/blob/master/bondspackage/bonds..." is a great example about bond. But I don't have the successfull compiled the bond_gran package into LIGGGHTS. Could you give me? Thank you very much. I'm looking forward to your reply.
Best regards
appple
apple | Sat, 04/22/2017 - 12:00
My problem
I have installed the LIGGGTHTs 3.3 and bonded several particles into a ellipsoid. Now I want this ellipsoid move in particle system(unbonded). But I failed,the feedback is :
ERROR: Fix bond/create can only be used together with dedicated 'granular' bond styles (../fix_bond_create_gran.cpp:225)
Could you help me ? Thanks very much.
NTT1508 | Wed, 12/03/2014 - 03:50
Compiling of bondpackage fail
Hi Christian,
Thank very much. I am happy to see the file grand_bond with ks,kn and rb. I havenot spent much time on this file but I think mechanism should be the same for displacement-force relationship bonding that used in PFC. In fact I havenot been able to run your example "ballfalloncone" successfully due to compiling error. Note I am using LIGGGHTS 3.0.6
First, according to your guide, I downloaded bondpakage and copied them to LIGGGHTS-PUBLIC/src. I then compiled again by typing make fedora. It came out with error ( unfortunately I havenot known how to crack error as you mentioned: >{warning|all|error}("message") to error->{warning|all|error}(FLERR,"message").). I am not quite sure how I can show you log file for that compiling.
Second, I then downloaded old version of LIGGGHTS (clone github) as you pointed out in the link. Here I found LIGGGHTS-BOND with bondpackage available. Compiling by make fedora -j didnot inform any problem however I could not run example because it states " lmp4bond: command not found". I think compiling was fail or in conflict with my previous LIGGGHTS (without bondpackage). When I check current version of LIGGGHTS by simply typing liggghts, it shows my current version is 3.0.6. That means LIGGGHTS without bondpackage is still active. Do I need to remove LIGGGHTS without bond before compiling LIGGGHTS-BOND.
Furthermore, I prefer to adding bondpackage to my current LIGGGHTS 3.0.6 which is in relation to CFDEM and OpenFOAM. Could you be patient to tell me again how to compile the bondpackage to my current LIGGGHTS 3.0.6? I am not really aware of why I have to use old version of LIGGGHTS for bondpackage. I am Geo-engieer and not good at C++ unfortunately.
And I am really surprise to hear that bondpackage is not provided and updated anymore in recent version of LIGGGHTS. I believe bonding among particles is one of the key features of modelling that is normally required in a wide range of topics.
Thank Christian for your kind and generous contribution.
Regards,
Nathan,
richti83 | Wed, 12/03/2014 - 09:13
Hi Nathan,
Hi Nathan,
The code in the bondspackage is not campatible to the LIGGGHTS3 Framework, for example the error class has changed. So you need to read the error messages of the compiling and change the lines in the source files the compiler tells you. Without C++ knowledge you should ask a software engineer to help you.
When the old version compiles fine you should get a lmp_fedora binary in LIGGGHTS-BONDS/src/ directory, I predict you missed the 2nd step and did not created a symbolic link in your /usr/bin directory, so the bash can not find lmp4bonds executable ...
sudo ln -s /home/YOURUSERNAME/LIGGGHTS-BONDS/src/lmp_fedora lmp4bonds
As I said in my last post.
I'm not related to DCS GmbH or CFDEM Project, so I can not tell you why the bondspackage is not supported officially.
About your question how to add bonds to your CFDEM coupling: when you have a working LIGGGHTS3-standalone version just copy the bond files to LIGGGHTS/LIGGGHTS-PUBLIC/src/ and run cfdemCompLiggghts again to recompile the library used by cfdem-coupling. I'm not familiar with CFD-stuff but I think the solver wont give correct results for bonded conglomerates of particles (wrong voidfraction, c_d value ....)
richti83 | Wed, 12/03/2014 - 11:37
Marry Christmas
I spent this morning to hack the original code to work with LIGGGHTS 3.0.6.
See this video for an eye-candy:
http://youtu.be/o3KRwzNKBmg
The code is now here for public use:
https://github.com/richti83/LIGGGHTS-WITH-BONDS
(note: I'm NOT the author of the bond package, I just copied the files from polyun - fork to current tree and adapted the code by adding missing methods and changing some LOC to fit to the new api of L3)
And the sample inputskript is here:
https://github.com/richti83/LIGGGHTS-WITH-BONDS/blob/master/examples/bon...
(note: I tried to use values from a PFC tutorial to compare the result later, so the units are a little bit funny for European)
Marry Christmas Time,
Christian.
NTT1508 | Wed, 12/03/2014 - 22:48
Fantastic specialist of LIGGGHTS and C++
Wow Christian, It is amazing work. I feel you're coding like having fun ;). I wish to have your skills and knowledge of C++ beside.
Anyway, I have not tried how it work but it sounds bright. What I need to do now is to replace my currently installed LIGGGHTS with your updated one then compile again with OpenFOAM and CFDEM. Hope to share with you some result soon. And I am excited to know you have intention of comparing LIGGGHTS and PFC. It would be highly expected game for many.
Thank very much and Marry Christmas,
Regards,
Nathan,
NTT1508 | Thu, 12/04/2014 - 09:19
some minor errors in example
Hi Christian,
I think I am successful to compile LIGGGHTS 3.0.6 with bond and it works. I tried to run your examples (ballfalloncone and beam.bonds). It run with some syntax errors such as mis-placing comments # and also warning about beta version of bonds. Finally I got stuck at one errors below that I checked manual but could not determine exactly where it is in script. Some values are out of defined range.
bond_coeff 1 0.2032 271447137401.575 65147312976.3781 1 1e10 1e10
ERROR: Numeric index is out of bounds (../force.cpp:755)
Could you double check again your script? That error appeared in both examples (ballfalloncone and beam). Basically I think it works now, not quite sure if there is any limitation of application.
Thanks,
Regards,
richti83 | Thu, 12/04/2014 - 13:23
ups
Argh, sorry I did forget to commit the last code changes to github, it is now fixed.
Please git pull again from richti83 repository and recompile.
Correct example is
atom_style hybrid granular bond/gran n_bondtypes 1 bonds_per_atom 6
....
create_box 2 reg
The ballfalloncone is not mine, but should be easy to adapt to new syntax.
NTT1508 | Fri, 12/05/2014 - 01:13
Now it works
Hi Christian,
I have updated LIGGGHTS with bond and now it works greatly. Thank very much for your interest. I will build my personal case in which I will import coordinates of a quite big number of particles from a file and assign them bonding. Hope to share with you and get your advices soon.
Cheer,
Nathan,
JoshuaP | Mon, 01/12/2015 - 16:12
hi
is it possible to take just the necessary files for the bond function? Because I extended my code already and dont wanna lose this. If so, which files are necessary for the bond_gran?
Thanks
Joshua
JoshuaP | Mon, 01/12/2015 - 16:13
btw
is it possible to break the bonds? Is there a limit value which can be set?
Thanks
richti83 | Mon, 01/12/2015 - 18:09
Patch
Hi,
you need at least these 8 files from src/ folder:
https://github.com/richti83/LIGGGHTS-WITH-BONDS/commit/fb3287e629331633c...
I would pull from the fork to a new folder and copy the files by hand my local master to have the newest files.
Yes bonds can break (after commit https://github.com/richti83/LIGGGHTS-WITH-BONDS/commit/7530ce80595e352be...),
have a look in https://github.com/richti83/LIGGGHTS-WITH-BONDS/blob/master/src/bond_gra... Line 305 - 339
Best,
Christian
NTT1508 | Tue, 01/13/2015 - 04:59
Break_Stresses and Radius
Hi Chistian,
I have looked at beam.csv to obtain development of tip (lastone) atom's displacement over time. After 10000 steps stable without pull_out force (fx=0) activated, the displacement of such atom is linearly increased over time due to pull_out force fx. However, if I change break_stresses ( 0.001 to 1.e+10), I cannt see the breaking point in displacement-time graph. The same trend can be shown if I increase pull_out force (from (10 to 1000)*F). What I can see is that all atoms except freezed one are moved equally along x axis. The higher pull-out force is applied, the farther displacement of atoms is.
I understand the theory of this kind of bonding that when stresses of bond disk reach break_stress values, the bonding will be broken. By that time, we should observe free movement of tip atom (that should be seen on displacement-time graph or Paraview).
Breaking condition by r_break (breaking type 0) is more understandable. When I try with a suitable value of r_break, the bonding is broken at a certain time and atoms move freely. I am wrong somewhere, am I?
Many thanks, Christian. You have inspired us much.
Regards,
Nathan,
NTT1508 | Tue, 01/13/2015 - 05:51
Visualize bond stress
Moreover, Is that possible to use your bond_reader plugin to visualize the development of bond stress over time? For example, in stead of placing a force at the tip of beam, assign the lastone (tip) of beam a velocity ( constant strain). In that case, it is neccessary to see how stress developed. Can dump/local and your liggghts_forcechain_reader work out this issue, Christian?
Thanks,
Regards,
Nathan,
richti83 | Tue, 01/13/2015 - 10:23
Break Stress
Hi Nathan,
about your first question: What you see there is the influence of the force dissipate factor in bond_gran line 286. Set it to 1.0 and the tip ball will not move. (Otherwise every timestep the bond force is damped by 0.5%). We need to develop a velocity-depending damping function. But as there is no funding for this, it is hard to do in less freetime.
Please see TODO in bond_gran.cpp:
/*NP
large TODO list for granular bonds: (could be a diploma thesis?)
+ need a better dissipative formulation than the hardcoded
'dissipate' value which produces plastic deformation
need some vel-dependant damping
+ need to carefully debug and validate this bond style
valiation against fix rigid
+ check whether run this bond style w/ or w/o gran pair style active,
(neigh_modify command)
+ need to store bond radii per particle, not per type
+ parallel debugging and testing not yet done
+ need evtally implemetation
*/
Here is a graph with and without damping:
http://www.richtisoft.de/transfer/vrgl_dissipate.png
About your second question: when you increase break_stress to 1e10 the bondstress will not be lager than the break condition. I tested with 1e4 and 1000 kip and git the following result:
http://www.richtisoft.de/transfer/bond_breaks_at_t0p1.png
The bond stress depends on force AND stiffness, so you realy need to calibrate BOTH parameter - I spent half a year to find good values for a specefic material, so this is a kind of magic and very time consuming.
At the moment you can not directly access bond stress. But you can try to dump bondforce with compute property/local bforceX bforceY bforceZ and dump it and than use a calculator filter to devide it by your bond-area to get normal stress.
Best,
Christian.
NTT1508 | Wed, 01/14/2015 - 12:09
Your verification is important
Hi Christian,
Well I thought I had tried to calibrate both parameters force and stiffness simultaneously. It is interesting to see dissipative consideration here. I carried out several tensile tests for natural fibres that could be a good source to verify.
Thank Christian,
Best regards,
Nathan,
NTT1508 | Fri, 01/16/2015 - 12:16
LIGGGHTS 3.07 with Bonds
HI Christian,
Could you please compile bond_gran for LIGGGHTS3.07 just released (like you did with ver3.06)? I have tried version 3.07 with mesh/surface/stress/servo and that seems to kill the previous bug completely. It would be great for my model with bond and mesh/surface/stress/servo in hand.
Or you can briefly tell me how to do it, I will try.
Thank you much,
Have a good weekend,
Ragards,
Nathan,
richti83 | Fri, 01/16/2015 - 13:25
pull repo
Hi Nathan,
well, I merged with 3.0.7 allready 9 days ago. This fixes this problem. You should only need to git pull from richti83 fork and recompile src.
Here is the commit:
https://github.com/richti83/LIGGGHTS-WITH-BONDS/commit/626b365369441d7a8...
(see changes in src/multi_node_mesh_I.h at the bottom of the page, this fixes the servo bug)
In Addition I fixed the sign error in forcechain output mentioned here:
http://www.cfdem.com/forums/sign-error-fixwallgran-relating-cwl-addwall1 with my last comit to https://github.com/richti83/LIGGGHTS-WITH-BONDS
(you asked how about visualisize particle wall contacts here, see my in_v.txt in the other thread how it's done. My forcechainreader can read particle-wall dumpfiles too)
Best,
Christian.
NTT1508 | Sun, 01/18/2015 - 04:24
Thanks
Hi Christian,
That's great. Thank for your active and close update.
Have a nice week,
Regards,
NTT1508 | Wed, 01/21/2015 - 04:16
Cannot be used with write_restart
Hi Christian,
I found that bond_grand currently cannot be used with write_restart command. fix_bond_propagate_gran.cpp states that Restart functionality not yet tested for granular bonds ...". Is that possible to enable wrtite_restart with bond_grand on current version of LIGGGHTS?
Thank you,
regards,
Nathan,
richti83 | Wed, 01/21/2015 - 13:13
not sure
The question is, what should bond/propagate save to restart ?
Try uncommend line 214 and 243 and recompile.
Mybe Chris can tell us more about the reason for inserting these error messages.
NTT1508 | Sun, 01/25/2015 - 00:16
Recompiling fail due to uncommending
Hi Christian,
Thanks, I tried to uncommend and re-compile but re-compiling fail due to this uncommending. I will contact with Chris to clarify about this.
Best regards,
richti83 | Mon, 01/26/2015 - 11:32
I did it
Hi Nathan,
not sure what you have commented, but I tried by my self and found three other stumbling blocks. However I got it to work, see this commit:
https://github.com/richti83/LIGGGHTS-WITH-BONDS/commit/e2270cd3bd43028a6...
An Example is there too.
Be carefull to define bond_style gran BEFORE you use read_restart - took me one hour to track this problem down.
Best,
Christian.
NTT1508 | Thu, 01/29/2015 - 00:47
How awsome you have done
Hi Christian,
Sorry for my late response. I have tested your new commit. It works fine. I checked what you had done with fix_bond_propagate_gran.cpp/h. what I can see is uncommending //error->warning(FLERR,"Restart functionality not yet tested for granular bonds..."). I think you had done more than that. Could I have learned from your action, Christian?
I did git pull on your github to my current directory, it updated 8 new files including fix_bond_propagate_gran.cpp but It then fail to compile. However when I did git clone your github and it compiles and run successfully. Possibly I did something wrong with git-pull. To be honest, I am a new guy in LINUX based system and Open-source softwares but it's really lucky for me to have your great advices constantly especially regarding LIGGGHTS. Without you and this forum, I would be in the dump.
Thank very much,
Best regards,
Nathan,
richti83 | Thu, 01/29/2015 - 13:57
what I've done:
Could I have learned from your action, Christian?
Offcourse, just disabling the error message does not solve the problem. There was a change in the place where the bond parameters are set.
Let's have a look at the changes in AtomVecBondGran Class:
first off all, read restart calls ->settings() but settings() expects 4 args. My work-a-round is to do nothing when no params are given. This will lead to a crash when one forget to give the bondparameters in atom style command and it is not a restart case, but that's why this is a beta version :-)
Second problem was that the functions read_restart_settings and write_restart_settings has been missing in src/atom_vec_bond_gran.cpp , I reverse engineered the pair_style table which uses parameters in restart case and adapted this to the atom_vec_bond_gran.cpp
Third of all a variable name changed in AtomVecBondGran::size_restart() which leads to a segmentation fault because the pointer pointed to not existing data.
Last but not least I found out that the bond gran commend has to come before restart command to set correct bond history size, allocate memory and so on. That's why I added an error->all message in AtomVecBondGran::grow where the memory is allocated for bond history.
Hope it's now clear,
Christian.
richti83 | Thu, 01/29/2015 - 14:03
make clean-all
I did git pull on your github to my current directory, it updated 8 new files including fix_bond_propagate_gran.cpp but It then fail to compile.
you should make clean-all everytime you updated the code, before make fedora -j otherwise the old object files persists in the source tree and the linker complains about not found blabla. (the -j NUMBER_OF_CORES specifies the number of processor to use to compile, this significant reduces build time)
NTT1508 | Fri, 01/30/2015 - 13:38
Thank your kindness
Hi Christian,
Thank for your kind and detailed explanation. It definitely helps me understand the logic of coding alot.
Have a good weekend,
Best regards,
Nathan,
NTT1508 | Tue, 02/03/2015 - 05:25
Constant Strain Rate
Hi Christian,
There is one point that's quite strange in tension test using gran_bond that I am trying (I think I am wrong somewhere). Suppose I have a string of atoms. One end if fixed while another is pulled by a constant velocity (velocity or fix/move) (tension by constant strain rate). Due to velocity (displacement of atoms), bond should be broken at some points. However problem is that velocity not come into effect with bonding. No atom moving when using velocity or fix/move an one end. Do I understand correctly.
Thank much, Christian,
Regards,
NTT1508 | Thu, 02/05/2015 - 05:36
My Mistake
Hi Christian,
Sorry for disturbance, I made a mistake. I got it worked now.
Regards,
richti83 | Fri, 02/06/2015 - 18:23
mind sharing your experience
mind sharing your experience what causes the problem ? i.e what you expected and where you made the mistake ?
NTT1508 | Sun, 02/08/2015 - 13:41
Dump Bond Force
Hi Christian,
I have spent recent 4 days to calibrate coefficients of bond/gran for my material but have not succeeded yet. There are several unclear points that would like hear your advices. Note that my intention is to simulate a tension test with constant strain rate. To make atom string (around 10 atoms) tensile, I apply a constant velocity to the atom at one end while another end of string is fixed. I try with both commands: velocity and fix/move because they donot give the same results in my case.
(1) I create atoms by several methods and see how they are different with respect to bond/gran and tension force. First I create string of atoms by atom_create in lattice region, second by atom_create as single one in specific coordinates and finally create atom by read_data. And it is quite interesting that the velocity recorded is not constant.
- For lattice region: velocity recorded at the tip of string (which is applied velocity directly) is always far different from input value in velocity/set or fix/move commands. For example, input velocity of 1 m/s results in recorded value at 10^(-3) m/s. I cannt explain why they are much different. (velocity is recorded by fix/print).
- Creating atom by read_data and atom_create/single give the same results which are divided into two cases. In the first case (1), one end of string is free (not fixed), all atoms in string move in the same velocity which is the same as input value (this makes sense to me). In the second case (2), the top of string is fixed while the bottom atom is moved. Now the velocity is recorded at the bottom atom is different from input value (say 0.7 m/s compared to input value of 1m/s). If I reduce time step, the recorded velocity comes closer to input one. Note that all time steps selected satisfy Rayleigh and Hert condition. Could you understand this?
(2) For bonded atoms without freezing one side, If I use velocity to make tip atom moved, the actual recorded velocity on tip particle is always equal to (input velocity)/(number of particles bonded). On the other hand, fix/move generates actual velocity on atom as prescribed (input) value. I cannt see any explanation of velocity/set in manual for this respect.
(3) Using compute property/local for bonding force doesnot give expected outcome. Though relative distance between particles in a bond is increasing step by step due to constant movement (fix/move) bond force recorded by dump/local is constant for all steps. Bond force can be calculated based on distance, bond-disk area and stiffness by spread sheet. But this calculated value is not matched to one output from compute/property/local.
Sorry for my long words but this is complicated investigation. Your patience and understanding are appreciated. Could you try the tension test under constant velocity (strain rate) on your beam.bond example then extract the bond force until the bond broken? I hope what I got is due to my scripting problem, not due to the bugs of codes. Otherwise I would be in trouble.
Thank very much, Christian for your interest and active contribution.
Best regards,
Nathan,
richti83 | Thu, 02/12/2015 - 19:58
I played arround with this
I played arround with this problem last weekend, but can not give you a clear answer.
First I checked with fixed velocity:
variable vx_tip equal 1 #pull with 1m/svelocity tipball set ${vx_tip} NULL NULL #expected: vel is once set and than the integrator overules the value [confirmed]
#final velocity = initial velocity - F_bond/m*dt (equilibrium of F)
#(gives recvel=0.508*vx_tip)
As expected the velocity is once set and than the integrator overcomes.
2nd I checked with equal style variable velocity setting:
velocity tipball set v_vx_tip NULL NULL #expected: vel is set every ts [confirmed] (gives recvel=vx_tip)
In this case the command velocity seems to overrule the integrator (which makes sense when the command is executed after all other fixes)
same is with
set atom ${lastone} vx 1After that I decided to not manipulate the velocity because in the two-step velocity-verlet algorithm used in ligggghts (have a look in wiki !!) it's a bad idea to overule the integrator-start-condition.
Instead I made a script-based force-controller:
#force 'controller'
variable mball equal mass[${lastone}]
variable acc equal (${vx_tip}-vx[${lastone}])/dt #recalculate acc every step
variable Fpull equal ${mball}*v_acc #calculate force to accelerate tipball
fix tipF tipball addforce v_Fpull 0 0 #apply acc-force every step
This gives a reasonable velocity of 1m/s as expected.
Now it's your turn to investigate what happens to the bonded chain. In my case the first bond breaks and the tip ball moves with constant velocity to right border of domain or it does not break and my beam is enlarged until hell freezes over (unbreakable bond because I use break-stress condition with very very high sigma_break). Maybe there is still a bug in the bond-algorithm, it seems that bonds break, but are never deleted.
Here is a nice animation of the force-controller-strain-bond-test:
http://youtu.be/2BgvgeNcoVU
skript:
http://www.richtisoft.de/transfer/beam_setVelTip.bonds
Note: fix print is a liar about particle velocity in some cases when it comes before the velocity command it dumps 0 or e^-200, you need to have a look in ParaView because the dumped velocities seems right.
I'm not totaly sure about the workflow but something like this seems to happen:
command
computes
fixes
dumps
command
....
NTT1508 | Sun, 02/15/2015 - 08:20
Extract the bond force and delete the bond
Hi Christian,
Many thanks for your passion and hard work, you would be a great professor in the future. Your force controller works like a charm to manage speed of atom that's clearly better than fix/move and velocity/set. Using fix/move/linear can generate expected velocity for atom but it is very sensitive to the size of time step as far as I tried.
More importantly, breaking condition of bond is a doubtable issue. In previous discussion, when tip force is increased until a certain value, bond looks like broken on Paraview. To make sure whether bond completely broken or not, I unfix the tip force and continue to run, the tip ball move back to its neighbour but it behaves like a simple contact, not a bonding. However, look at the bond*.bond of properties/local dump files, the bond between tip ball and its neighbour is still existing. How can we trust these files?
Besides, I carry out a simple test for 3 atoms bonded subjected to constant velocity (using your force controller). The tip ball moving constantly under a certain velocity pulls the middle ball move together due to bond force. The farther distance between two atoms, the velocity of neighbour becomes higher due to the bigger bond force generated. However, bond seem to never be broken even I reduce the break stress to very small value (e.g., 1e+1). If you have a look at dump file by compute/property/local/bforce X/Y/Z, you can see bond force is always constant though time varying and distance between atoms increasing. My question is how to extract the bond force conrrectly? I cannt see how the bond force developed. If I calculate bond force based on distance; kn and bond-disk area, the bond force is really big and far different from what recorded on dump files.
See the link download for script of my simple test of 3 atoms below.
http://www.datafilehost.com/d/46fd4407
I have seen the bond_gran.cpp and I believe all formulas used to calculate the bond force are correct (for this type of bond). There is only one thing that I am confused is at breaking condition point. What happened next after conditional lines?
bool nstress = sigman_break[type] < (+nforce_mag/A + 2.*ttorque_mag/J*rb[type]);
bool tstress = tau_break[type] < (+tforce_mag/A + ntorque_mag/J*rb[type]);
Which lines after this to deactivate the bond? Can we delete the bond if the conditions above satisfied?
Thank very much Christian. Really exciting whenever see any of your responses in this forum.
Best regards,
naren | Mon, 03/09/2015 - 07:08
Hi Christian
I have a couette geomety (of parallel plates) system, in which bottom plate is fixed and top plate is moving with a fixed velocity. This system is filled with granular particles. I have an agglomerate (particles with bonds) lump in between the free flowing granular particles. I want to see that, how the bonds between the particles in the agglomerate breaks with the shear in the system. Can i use this LIGGGHTS-WITH-BONDS for my system? Any suggestions and advice will be much helpful for me.
Thanks,
Naendra
richti83 | Mon, 03/09/2015 - 08:53
Yes you can
but be aware that this is a beta version and is not maintained by the liggghts core developer. I do this in my freetime.
There are some restrictions you should know about:
/*NP
large TODO list for granular bonds:
+ need a better dissipative formulation than the hardcoded
'dissipate' value which produces plastic deformation
need some vel-dependant damping (1)
+ need to carefully debug and validate this bond style
valiation against fix rigid
+ check whether run this bond style w/ or w/o gran pair style active,
(neigh_modify command)
+ need to store bond radii per particle, not per type (2)
+ parallel debugging and testing not yet done
+ need evtally implemetation
*/
(1) this is the main limitation, with dissipative "damping" the bond forces get lost over time, but without (set factor to 1) the behaviour is like glass
(2) this could be easily implemented like in PFC3D by introducing a bond radius scaling factor lambda and than say
rbmin=rb[type]*MIN(radius[i1],radius[i2]); //lamda * min(rA,rB) see P.Cundall, "A bondet particle model for rock"A = M_PI * rbmin* rbmin;
J = A * 0.5 * rbmin * rbmin;
(note I misapply rb[type] to route lamda from inputdeck to code for a quick&dirty test of this idea).
As there is no funding for me to test and debug I will not improve the code itself, but give hints how it could be done.
Please note: the bond code has been initial released by CKloss in 2012 and has been extended by PFodor later [2]. I only made it compatibel to LIGGGHTS 3.x code base, so it is not "my" intellectual property.
[1] http://www.cfdem.com/comment/12090#comment-12090
[2] https://github.com/Polyun/LIGGGHTS-PUBLIC/commit/a73e80ecb971fb633716453...
naren | Tue, 03/10/2015 - 11:12
Thanks
Thank you Christian. I am trying to understand the code. Do you mind sharing ballfalloncone example with me? Just to get familiar with the script.
Thanks,
Narendra
richti83 | Tue, 03/10/2015 - 12:40
ballfalloncone
https://github.com/Polyun/LIGGGHTS-PUBLIC/tree/master/bondspackage/bonds...
it is allready out there. Compare syntax with this:
https://github.com/richti83/LIGGGHTS-WITH-BONDS/tree/master/examples/bon...
(U need to swap n_bondtypes and bonds_per_atom from create box to atom_style line and change pair style to L3 syntax)
naren | Wed, 03/11/2015 - 06:34
Thank you
Thanks for your time Christian. I have one more small doubt.
With the atom_style
hybrid granular bond/gran n_bondtypes 1 bonds_per_atom 6
Can we define the bonded particle in between the free particles?
Thanks,
Narendra
richti83 | Wed, 03/11/2015 - 15:07
yes of course
You can do this with 2 different atom types, say atom type 1 should not be bonded but atom type 2 should be bonded with bond-type 1:
fix bondcr all bond/create/gran 1 2 2 0.0051 1 6 #every 1, type 2 to type 2, cutoff 0.0051m, bondtype 1, maximum 6 new bonds per timestep
naren | Thu, 03/12/2015 - 06:02
Thanks
Thank you Christian.
NTT1508 | Mon, 05/18/2015 - 11:40
Output information of bonding
Hi Christian,
Are we able to dump the bond force of particle or by any means of exporting the bond force to an external file? It is impossible to verify the bond force based on the total force dumped by the normal DUMP commands. Currently I can check the bond force by a bit editing the bond_gran.cpp but donot know how to write specific information of bonding ( e.x. the bond force for a certain particle with its ID).
Thank you,
Regards,
richti83 | Mon, 05/18/2015 - 13:26
yes
there is a un-documented AddOn in my github fork for compute_property_local:
else if (strcmp(arg[iarg],"bforceX") == 0)
pack_choice[i] = &ComputePropertyLocal::pack_bforceX;
else if (strcmp(arg[iarg],"bforceY") == 0)
pack_choice[i] = &ComputePropertyLocal::pack_bforceY;
else if (strcmp(arg[iarg],"bforceZ") == 0)
this means, make a compute bforce all property/local bforceX bforceY bforceZ
and a dump local c_bforce[1]..[3]
should do the trick.
Be carefull, this compute is only known to work in single-processor execution bc. in mpi mode it is possible that the bonded atoms are in different domains and than it will crash bc. it doesnt know the bondforce history array of the other proc.
Cheers,
Christian.
NTT1508 | Tue, 05/19/2015 - 10:04
Deviation of Velocity
Hello Christian,
I have checked bond information (force and ID) that is dumped by your compute_property_local. It is matched exactly to what I got from bond_gran.cpp. A bit trouble is that the dump/local commend exports information in a complicated format that might require some data filtering before using (e.x. just graphing bond force over time)
Additionally, I found a significant difference between velocity used in bond_gran.cpp and velocity of particles that obtained by dumping or fix/print. this is not understandable for me.
Many thanks,
Regards,
richti83 | Tue, 05/19/2015 - 12:45
Bond Reader
You can use a modified version of the bond reader:
https://github.com/richti83/ParaView_Reader_for_LIGGGHTS/blob/master/rea...
and than PlotSelectionOverTimer in PV.
Just add 3 more fields after line 215 to catch the force-output and store it in a new cellarray.
I'm on holliday this week, when you can wait until next Thuesday I can make the changes for you.
About fix print: this fix seems to be invoked after the 2nd half verlet step so it is after velocity integration. But the force models (pair,bond style) using the first half step value as far as I know.
Best,
Christian.
NTT1508 | Mon, 05/25/2015 - 01:21
Compile Paraview 4.3.1-source
Hi Christian,
I have used Matlab to filter the output. Additionally I tried with your Plug-In on Paraview4.3.1-source. However I could not carried out compiling Paraview successfully. This happened last time with Paraview 4.3.0, so that I had to use your pre-compiled plug-in.
Following is the error come when I tried to make in Paraview/build
[ 36%] Built target vtkFiltersGeometry
Built target vtkIOPLY
[ 36%] Built target vtkFiltersStatistics
Built target vtkParallelCore
[ 36%] Built target vtkImagingGeneral
[ 36%] [ 36%] [ 36%] Built target vtkFiltersVerdictCxxTests
[ 36%] Built target vtkPVCommonCxxTests
Built target vtkCommonExecutionModelCxxTests
Built target vtkIOMovieCxxTests
Linking CXX static library ../../lib/libvtkIOImageCS-pv4.3.a
[ 36%] Built target vtkIOImageCS
Linking CXX static library ../../lib/libvtkFiltersCoreCS-pv4.3.a
[ 36%] Built target vtkFiltersCoreCS
Linking CXX shared library ../../../lib/libvtkIOXML-pv4.3.so
[ 36%] Built target vtkIOXML
Linking CXX shared library ../../../lib/libvtkFiltersGeneral-pv4.3.so
[ 36%] Built target vtkFiltersGeneral
make: *** [all] Error 2
Thank you much,
Regards,
richti83 | Tue, 05/26/2015 - 14:52
Update
I just pushed the enhanced precompiled bond reader to github.
https://github.com/richti83/ParaView_Reader_for_LIGGGHTS/tree/master/pre...
You need to use a dump local like this:
compute b1 all property/local batom1x batom1y batom1z batom2x batom2y batom2z batom1 batom2 btype bforceX bforceY bforceZ
compute b2 all bond/local mol
dump bnd all local 1000 post/bonds*.bond c_b1[1] c_b1[2] c_b1[3] c_b1[4] c_b1[5] c_b1[6] c_b1[7] c_b1[8] c_b1[9] c_b2 c_b1[10] c_b1[11] c_b1[12] #x1 y1 z1 x2 y2 z2 id1 id2 typeofbond molID fbx fby fbz
About your compile error: I can not see what's going wrong, maybe better ask the PV mailing list and attach a full compile log.
Sometimes it helps to start make again with make -j1 (to only use one proc).
Was the ccmake step without errors ?
Cheers,
Christian.
NTT1508 | Wed, 05/27/2015 - 11:09
Nice Bond Visualization
Hi Christian,
Well, I have checked bond force output by compute property/local and dump/local. They are all right. And plug-in bondreader to visualize bond on Paraview works well. But I think it cannt represent the magnitude of bond force
Only strange point is a significant conflict of velocity. Velocity calculated by (Position variation)/timestep is much different from the velocities obtained by Fix/print or Dump. Difference of 1 half verlet step should not result in much variation. But Bond_force is quite close to my analytical calculation.
Steps on ccmake including configure and generate look Ok. No error reported. I will try to understand this.
Pages