About mixing of coefficients of pair_style granular

Dear all
hope you are all well recently, I hope to obtain some help for the mixing of coefficients of pair_style granular. I hope to have two different materials, but both the “hertz/material” model.

for material 2, i hope to have “rolling” part, while material 1 dose not have “rolling” part.
but this cannot really work in the lastest lammps by default mixing by lammps, therefore, i wish to give the “pair_coeff 1 2” by myself and want to give pair_coeff 1 2 without “rolling” part, but the remaining parameters just follow the equation for lammps documentations:
however, i do meet a problem for Young’s modulus of pair_coeff 1 2. Firstly, for the equations of effective young’s modulus, there are two different equations for E_eff,ij; see equation 2 and 3 in the image.

i attached the image (from lammps documentations) for these two equations. I am not pretty sure whether this is a typo for the equation 3, whether the square is missing in the equation 3.
So I assume that the equation 2 is right. So in my case, i tried to check this without “rolling” part;
so input script 1, the mixing coefficients are automatically obtained by LAMMPS; input script 2, i give
the mixing coefficients based on the equations from LAMMPS documentations. I think i should get the 100% same results.

for my material 1, E1 = 7e6, v1 = 0.5; material 2, E2 = 7e10, v2 = 0.2, so based on equation 1, the Eeff,ij = 9.33e6; then i should give the suitable Eij so that i obtain the same Eeff,ij (9.33e6), based on my calculation from equation 2, when Eij is 1.68e7, the Eeff,ij is 9.33e6, which should be ok, but it seems that i got totally different simulations results for these input script 1 and input script 2:
see input 1, the key sentences are:

pair_coeff 2 2 hertz/material 7.0e10 0.1 0.2 tangential mindlin_rescale/force NULL 1.0 1.0 
pair_coeff 1 1 hertz/material 7e6 0.5 0.5 tangential mindlin_rescale/force NULL 1.0 0.001

see input 2, the key sentences are :

pair_coeff 2 2 hertz/material 7.0e10 0.1 0.2 tangential mindlin_rescale/force NULL 1.0 1.0 
pair_coeff 1 1 hertz/material 7e6 0.5 0.5 tangential mindlin_rescale/force NULL 1.0 0.001
pair_coeff 1 2 hertz/material 1.68e7 0.316228 0.316228 tangential mindlin_rescale/force NULL 1.0 0.031623

all parameters are geometric averaging, and Eij is calculated based on equation 1,2 and 3.

so does anyone has the idea that how should i correctly give the values for mixing of coefficients, i just want to obtain the same results compared to the default results obtained from lammps

thanks for your time,
hope you all have a nice day


This looks like it’s just a typo in the documentation. There should be a square exponent in Eq. 3, which is present in the code.

This looks correct when compared to the code, so your argument for pair_coeff 1 2 should be 9.33e6 not 1.68e7. I’ll have to check with the primary author of pair granular whether there’s a disconnect in the documentation.

Should the mixed value of Poisson’s ratio be ~0.22?

Dear jtclemm

Thanks for your response, i have also checked the source code of pair_style granular, it seems that the equation 2 is the right one. it could be a typo for equation 3.

for the calculation of pair_coeff 1 2, it seems that the lammps automatically calculate the effective Eij based on equation 1 if we do not give the Eij.
so if we give the exact value of Eij, the effective Eij should be calculated based on equation 2.
it seems that the Eij is different compared to effective Eij.
but for the mixing, i have a idea for using mixing, i can use rolling sds 0 0 0. Although we have rolling command, do not really have this rolling

Hi, I have tried to get help from one expert from uk, he is very nice and help to check the source code. He mentioned that the stable version of lammps 23 Jun 2022 may have some potential bugs for the granular pair style;
I cited this for his email.
In the version you are using, I already found one clear bug. In PairGranular::init_one, there is a condition that looks like the following:

if (normal_model[i][j] == HERTZ || normal_model[i][j] == HOOKE) {
[use the geometric mixing rule for Young’s modulus]
} else {
[use equation 2 in your e-mail]

Well, the problem is that normal_model[i][j] == HOOKE (surprisingly) as the mixed i/j value hasn’t been set properly. As you would expect, normal_model[i][i] == normal_model[j][j] == HERTZ_MATERIAL. So this is not right.
For example, you can write out values of E_eff from the compute by copy–pasting the line of code below into PairGranular::compute at around line 270 in pair_granular.cpp (just after E has been retrieved in the function):

fprintf(screen,“E in compute is %1.5e between %i and %i\n”,E,itype,jtype);

When you do that, you will get output in the terminal like the following:

E in compute is 6.22222e+06 between 1 and 1

E in compute is 5.49972e+08 between 1 and 2

E in compute is 5.49972e+08 between 2 and 1

E in compute is 4.86111e+10 between 2 and 2

These 1–1 and 2–2 values of E_eff are correct; the 1–2 (or 2–1) value is not.
I therefore hope to update this in this community in case anyone has the similar problem.

Hi Tom_lammps, thanks for checking the code. This is an old bug that was fixed in the current development version of LAMMPS.

However, the new code may have a different inconsistency with mixing. I’ve patched it in another branch but need to discuss it with my coworker to confirm.

Dear jtclemm
many thanks for your help and time for that. It is really helpful.
would you mind sharing any details for the new inconsistency? I am planning to use the new version of the LAMMPS, 28 Mar 2023. your message may let me feel a little bit anxious for that. do you mean that the automatic mixing of coefficient 1 and 2 shows inconsistency. the automatic calcuated values are not right?
so as a compromise at this stage. is this possbible for me to directly give the coefficient between 1 and 2 to cover the automatic calculation from lammps so that i can continue the research.

hope you have a nice day
best regards

Hi Tom_lammps, in the new code it’s possible that unspecified i-j coefficients (for i != j) are first mixed using Ei, Ej, nui, and nuj and then subsequently mixed with themselves. Details, and a patch, are in this pull request: Misc minor patches/features in BPM/Granular packages by jtclemm · Pull Request #3807 · lammps/lammps · GitHub

Currently, the best work around is to simply specify coefficients for all i-j pairs. Alternatively, you could try cloning the branch in that pull request and confirm it works as you anticipate. Hopefully, we’ll be able to finish checking the patch and merge it to the main LAMMPS branch relatively quickly.

Thank you again for thoroughly checking the mixing of coefficients and reporting potential issues.

Hi jtclemm
Many thanks for your detailed response, it is really helpful!. i will try to look this from the pull request, or maybe its safe to wait for official release of the new version. It always take much effort to get help of engineer to install new softare in hpc.
would you mind share the idea that when the new verision could be released ?

by the way, i have identified a thing (which i am not really clear) in your previous answer. If you do not mind, i cited that:
This looks correct when compared to the code, so your argument for pair_coeff 1 2 should be 9.33e6 not 1.68e7.
if i understand this correctly. it should be a disconnect between the documentation and code. but i am not pretty sure about this.
the documenation stated that if Eij is specified, the Eeff,ij would be calculated following this equation (this image). it is clear that Eij did not equal to Eeff,ij.
so, just for instance, in my case, the Eeff,ij should be 9.33e6, i should use this value Eeff,ij to recalcuate the right Eij, so that i can put this Eij in the argument (which is 1.68e7 in the calculation). I think this is what i can see from documentation. do i make any silly mistake for understanding this documentation?

hope you have a nice day!!

You can save a whole lot of this effort by simply compiling LAMMPS for yourself. There is nothing in LAMMPS that requires having system administrator privilege.

Compiling is not that difficult once your have figured out the basics, especially when using the CMake based build procedure. You have complete control over which features and settings are included into your LAMMPS binaries. When you use a git checkout, you can also follow the development directly and can retrieve different versions easily with a suitably crafted “git pull” command. It may bit a bit of a learning curve, but it quickly pays off. There are plenty of large scale clusters and similar machines, where you are required to compile any application software on your own.

Hi akohlmey

thanks for your suggestions!. I just have sucessfully installed and compile lammps in hpc. it works. thank you very much for that. but currently i am not sure how to use a git checkout and git pull. i directly upload the latest lammps.zip to hpc, unzip and complie. it is worthwhile for me to learn how to use these command. So once the Misc minor patches/features in BPM/Granular packages is merged into develop branch. I can directly download this and complie this, i may not need wait for the release of the stable version.

thank you for your time and suggestions. it is really helpful.
hope you have a nice day

To learn about using git, you can go to https://git-scm.com/

Unfortunately I cannot say. It depends on when other people have time to review and check the code.

Yes, if you specify E_ij and nu_ij for type i-j LAMMPS will use the equation you quoted to compute an effective stiffness. So you just need to figure out which value of E_ij and nu_ij is needed to give you whatever effective stiffness your application requires for every combination of i-j types.

Only when you don’t specify pair coefficients for some i-j types will LAMMPS derive a set of mixed coefficients from other pair coefficients which may lead to a (potential) bug.

hi jtclemm

thanks your clarity. it is really helpful. thank for your time. I will check this git website frequently to say what will happen. Thanks for your contribution.

its just i have identified another weird thing for bpm package. it just happen in LAMMPS Feature Release 28 Mar 2023. I have reported this weired thing in new matsci. Since you are the developer of this bpm. would you mind share any idea for this weird thing?

so if we get distinct gap between two version of lammps. dose this suggest potential bugs or problems?

hope you have a nice day

hi akohlmey

Many thanks for your time and response, I will try to learn these useful things