Contact Modeling in OpenSim

[Milad Ghafary]

Hi
I am studying and trying to understand how contact forces are modeled in the Simbody, especially Hunt and Crossley force. The model that I am studying is a sphere bouncing on a plane. The stiffness coefficient between the sphere and plane is 10000000. At the first step of understanding, I have neglected the dissipation and the effect of friction in this model. I have used the Forward Dynamic tool in OpenSim to simulate the contact between the sphere and the plane to obtain the contact force between them.
To validate the contact forces, I have subtracted the radius of the sphere from the displacement of the center of mass of the sphere in the Y direction. In instances in which this subtraction is negative, it means that the sphere is in contact with the plane and the contact force exists. I have calculated this contact force with the relation which is given on the SimTK website. k is the stiffness coefficient which is 10000000 in this model. x is obtained from the absolute magnitude of the subtraction that is mentioned above. And exponent n is 1.5 here. The c coefficient is assumed zero because of the neglection of the dissipation. When I calculate the contact force with this relation, the calculated force for each instance is 9.21525727104773 times the contact force that is calculated in the Simbody engine.
Could you please help me to understand what is missing here?
I need to validate this relation.
I have attached the contact forces calculated in Simbody and calculated from the mentioned relation in an Excel file.
Thanks
Milad

[Nicos Haralabidis]

Hello Mihad,

Is this the code you were referring too? https://github.com/simbody/simbody/blob ... ontact.cpp

This article also explains how the forces are calculated using the Hunt-Crossley model:
https://www.sciencedirect.com/science/a ... 3811000241

I hope that helps in some way!

Cheers.

Nicos

[Michael Sherman]

Hi, Milad.

Given a compliant sphere of radius R and elastic modulus E in contact with a rigid halfspace, Hertz contact force is

F = 4/3 sqrt(R) E x^(3/2)

where x is the penetration depth. With no dissipation or friction, that should give the total contact force. So the relevant questions are:

• What is the radius R?
• How does OpenSim map the "stiffness" to the elastic modulus E?

From your data, it looks like your sphere has radius R=.053 so 4/3 sqrt(R) is .307. That accounts for a factor of 3, but not the factor of 9 you are seeing. However, if instead I treat .053 as the sphere diameter so the radius is R=.0265 then 4/3 sqrt(R) is .21705 accounting for a factor of 4.60721 which is exactly 1/2 of the 9.21 factor you see. That suggests that we have E=k/2 and R=.0265 which would account for the numbers you are seeing. My speculation would be that the given OpenSim "stiffness" k is being treated as the elastic modulus of both the sphere and the halfspace, so that the composite elastic modulus is 1/2 of that.

That's my speculation but I don't actually know how OpenSim is mapping your model into a Simbody calculation; likely someone else on this forum can give you a more definitive answer.

Regards,
Sherm