# ViscoEl contact normal stiffness kn in the source file is not consistent with that in the documentation

Hello,

I am trying to use the visco-elastic model to do the simulation. I checked the kn calculation in the source file and that in the documentation, then I found that there is a factor 2 difference.

In the documentation, the stiffness is kn=2*kn1*kn2/(kn1+kn2), as shown in the following link:

However, in the source file, the Kn is calculated by kn= kn1*kn2/(kn1+kn2), as shown in source file line 279 and 315-323:
279: phys->kn = contactParameterCalculation(kn1, kn2)

315: Real contactParameterCalculation(const Real& l1, const Real& l2)
316: {
317: // If one of paramaters > 0. we DO NOT return 0
318: Real a = (l1 ? 1 / l1 : 0) + (l2 ? 1 / l2 : 0);
319: if (a)
320: return 1 / a;
321: else
322: return 0;
323: }

Why they are different?

I know in the DEM configuration, Fn=kn1*dl1=kn2*dl2=kn*(dl1+dl2), and then we get kn=kn1*kn2/(kn1+kn2) by derivation , but this has a requirement that (dl1+dl2) is the length corresponds to kn. However, the penetration depth in YADE is calculated by dl=R1+R2-L where L is the length from sphere1 center to sphere2 center, in this case, dl=1/2*(dl1+dl2). Then F=kn*(dl1+dl2)=kn*(2*dl)=2*kn*dl=2*kn1*kn2/(kn1+kn2)*dl, and the new kn = 2*kn1*kn2/(kn1+kn2), that is exactly the kn expression in documentation. However, in the source file, this factor 2 is missing. Am I thinking correctly?

Kind Rrgards,
Yuxuan Wen

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Jérôme Duriez
Solved:
Last query:
 Revision history for this message Launchpad Janitor (janitor) said on 2021-01-03: #1

This question was expired because it remained in the 'Open' state without activity for the last 15 days.

 Revision history for this message Jérôme Duriez (jduriez) said on 2021-01-05: #2

Hi,

If your question was "there is a mismatch between documentation and source code, which should I trust?" the answer is "source code".
Behavior of YADE executable is obviously the one of its source code (the one of the used version). Documentation just attempts to retranscript this behavior in a non-C++, more friendly, language. Sometimes the attempt fails.. (if that's the cas, it may be fixed in a near future)

If your question is "how to justify the kn constant expression, with or without a 2", this may have been discussed few times on the mailing list in the last 10 years (..) and I think there is simply is no answer.
As far as I'm concerned, I can imagine no bullet-proof theoretical development for such a linear approximation (why Fn = kn dl1 would hold ? This is not what Hertz theory says). And there are then several (size-independent, at the macro scale) expressions that are all equally as good, with or without multiplying factors 2, or pi, or whatever number you could think of.

 Revision history for this message Yuxuan Wen (wenyuxuan) said on 2021-06-28: #3

Thank you Jerome!

So do you mean that a constant factor such as 2, or pi, doesn't really affect the results of simulation? Because the fundamental equation in DEM Fn=kn dl1 is also not strictly proved and is just a approximation.

Kind Regards,
Yuxuan

 Revision history for this message Jérôme Duriez (jduriez) said on 2021-06-29: #4

Even though having a different kn will affect the simulation, I mean that you can always get the same kn (and the same simulation results), as soon as you change your kn-related user-parameter (that will combine with 2 or pi to give kn).

And there is no physical grounds to discriminate between one user-parameter value (the one combined with 2) or another (the one combined with pi) => no reason to prefer the expression with 2 or the expression with pi.

 Revision history for this message Yuxuan Wen (wenyuxuan) said on 2021-06-29: #5

Thanks Jérôme Duriez, that solved my question.