Home › Forum › SOFA › Using SOFA › [SOLVED] Measuring contact pressure
Tagged: 64_bits, contact, Monitor, Plugin_SoftRobots, pressure, SOFA_2006, Windows_10
 This topic has 14 replies, 4 voices, and was last updated 2 years, 10 months ago by btt.

AuthorPosts

18 January 2021 at 06:56 #18242twxuBlocked
Hi,
Has anyone had any experience with measuring contact pressure in SOFA? i.e. when object A contacts object B, measure the pressure experienced by object B at the contact area.
I have been using the Monitor object so far to measure forces on individual nodes (which has been working very well), but is there some way to convert the force measurements to pressure measurements?
For example, could the pressure on a single triangular mesh element be calculated from the force measured at the 3 triangle corners? As far as I know, the Monitor and ExtraMonitor objects do not have an option to directly measure pressure.
Any guidance is appreciated, thanks!
1 February 2021 at 11:20 #18460HugoKeymasterHi @twxu
Which contact responsive method are you using:
– penalty method
– or Lagrangemultiplier based resolution ?Depending on the method, I would point you out two different solution.
(If you are not sure about the method, could you just share your scene file?)Best
Hugo
1 February 2021 at 16:02 #18466twxuBlockedHi @Hugo,
Thanks for the response.
I am not 100% sure, but I believe I am using the Lagrangemultiplier based method.
My scene and object files are here: https://github.com/thomaswxu/SOFA
I hope it is not too annoying to use; you will have to change some file paths/names for the objects to be loaded correctly.
(i.e. change “fea_variant” variable in scene file, then change corresponding file paths in “filenames.py”)Please feel free to contact me if it doesn’t work or if you have any questions!
5 February 2021 at 12:38 #18536HugoKeymasterHey
From your scene singleFinger_wObject_XML.scn, it appears that you indeed use Lagrange multipliers: use of the FreeMotionAnimationLoop + GenericConstraintSolver.
Recovering contact forces can be done in the GenericConstraintSolver:
– activate the “computeConstaintForces” boolean
– this fills the nodal data field “constraintForces” that you can use to compute your pressureI hope this helps.
BestHugo
5 February 2021 at 15:44 #18542twxuBlockedHi @Hugo,
Thanks for your suggestion! I have activated the “computeConstraintForces” boolean, and can see the “constraintForces” data under GenericConstraintSolver in the SOFA GUI.
But, I’m not sure how to interpret this data. It appears to be a single row vector; I see in the SOFA documentation that it represents intensities of constraint forces, but I could not find further information.
Could you provide more details on how to calculate contact pressure from the “constraintForces” data? I’m not sure how to proceed with the pressure calculation currently.
26 February 2021 at 19:02 #18766HugoKeymasterHi @twxu
“constraintForces” is a vector which size equals the number of constraints (maybe contacts in your case). If you have 3 contacts you would have :
“Contact1_ForceInX Contact1_ForceInY Contact1_ForceInZ Contact2_ForceInX Contact2_ForceInY Contact2_ForceInZ Contact3_ForceInX Contact3_ForceInY Contact3_ForceInZ”
Best,
Hugo
26 February 2021 at 19:02 #18767HugoKeymasterTo compute the pressure, use directly the P = F/S relationship.
Best,Hugo
26 February 2021 at 19:25 #18769twxuBlockedHi @Hugo,
Thank you for the clarification; that makes sense.
How would you recommend computing the area term (S)? Is there a way in SOFA to calculate the area of individual mesh elements?
Or, would it be better to sum all the forces and divide by the entire area?17 March 2021 at 16:17 #18863HugoKeymasterHi @twxu
Good question.
To get such geometrical information, your C++ code must define a pointer a SetGeometryAlgorithms. For instance, you can find such computation on triangles in the MeshMatrixMass.h code :sofa::component::topology::TriangleSetGeometryAlgorithms<GeometricalTypes>* triangleGeo;
and its usage in MeshMatrixMass.inl :
triangleGeo>computeRestTriangleArea(triangleAdded[i])
I hope this helps.
BestHugo
23 March 2021 at 06:21 #18961nickljBlockedHi @Hugo
Sorry for continuing to ask questions on this “SOLVED” thread as I think my question is highly relevant to this one.
According to your reply in this post, “constraintForces” is a vector which size equals the number of constraints, and should have the form “Contact1_ForceInX Contact1_ForceInY Contact1_ForceInZ…”.
I’ve checked the source code of “GenericConstraintSolver.cpp” and “ConstraintSolverImpl.h”. It seems this “constraintForces” is nothing but read out the lambda value for Lagrange multipliers, as: “current_cp>getF()” which is actually “current_cp>f.ptr()”.
As a result, I think your reply in this thread seems different from that you posted before, as in:
which lambda also needs to be multiplied by the constraint values in the collision model to retrieve the contact force, which seems more reasonable to me.However, I’ve tried both, and in any way, I still cannot successfully extract the contact force correctly. The values seem to be more than an order lower than that I simulated in Abaqus, or measured in the experiment.
By the way, I used FreeMotionAnimationLoop and GenericConstraintSolver in my simulation.
Much appreciate if you may help to explain a bit more on this. I’ve been troubled by this problem for months…
30 March 2021 at 08:38 #19020nickljBlockedI think I may have a better understanding on the program now. lambda is the force in local contact space, and the first component is always the normal force. The “constraint” values in the “mechanicalObject” of the collision model, is actaully a transformation from local space to global coordinate.
Besides, the force obtained in constraint problem is actually dt*lambda. So in order to obtain the force in terms of unit Newton, we need to divided by dt.
Please kindly help to correct if I have any misunderstanding. Thanks.
27 May 2021 at 16:23 #19587bttBlockedHi @hugo, @nicklj,
Sorry but my comment is somehow disappeared eventhough I got reported from email that comment has been successfully summitted. I am interested in this topic as well. However, I have been struggling in either extraction of the vector of lambda or computing contact pressure as guided in thisCould anyone give me a little more detail regarding this procedure? or any other clarification would be appreciated.
Furthermore, I have the following questions:
1. let say constraint number 0 1 2 impact on node 55 then the total resulting force, socalled Fs[55], generated from collision (in x y z) on node 55 will be the summation of those calculated from constraint 0, 1, 2 right?
2. Secondly, the term F in the equation p=F/S is the summation of all contact forces Fs[i], where i are the impacted nodes, generated in the contact area, is it right?
Thank you for any help
best regards31 May 2021 at 08:23 #19595nickljBlocked@btt
For question 1. Just search the “constraint” in the collision MechanicalObject, and find all the node 55 related rows, for each row, multiples the lambda value to get the force vector in global coordinate, then add all these values together to get the total reaction force related to node 55.31 May 2021 at 10:14 #19596bttBlockedThank you @nicklj for the explanation, I did exactly what you described above however calculated values for reaction forces seem unreliable (too big). For example, the constraint from the collision model for node 8 is:
9 1 8 0.969057 0.104965 0.223409 10 1 8 0.0651265 0.981735 0.17876 11 1 8 0.238092 0.158679 0.958193 ... 15 1 8 0.971997 0.0886831 0.217617 16 1 8 0.0778631 0.995286 0.0578188 17 1 8 0.221718 0.0392554 0.97432
Lambda values for row number 9 to 11 are:
[11293731.263678096], [1448629.6310082236], [4136659.638746046]
whereas, those for row number 15 to 17 are zeros. Then the force for node 8 is
F = (constraint[9]*lamda[9]+constraint[10]*lamda[10]+constraint[10]*lamda[10])=[1.1835 0.1951 0.1700] *10^7
The unit system i used is: mm, tonne, N
Even when I multiply this value with time step dt = 0.01 as indicated in sofa documentation, the final values are still unbelievably large. Do you have any comments on this issue? Or any other help will also be appreciated.
best31 May 2021 at 14:36 #19585bttBlockedHi @hugo, @nicklj,
I am interested in this topic as well, I also want to extract the force generated by collision as pressure on contact area. I followed all related post in forum especially this one. From Sofa documentation, I understand the principle of Lagrangemultiplier based method is to create forces on two objects potentially in contact to create constraint (in this case is a nonholonomic constraint)between them. However, I have been struggling in either extraction of the vector of lambda or computing contact pressure as guided in this. Could anyone give me a little more detail regarding this procedure? or any other clarification would be appreciated.
Furthermore, I have the following questions:
1. let say constraint number 0 1 2 impact on node 550 1 55 0.172748 0.984895 0.0118137 1 1 55 0.0161293 0.0148211 0.99976 2 1 55 0.984834 0.172516 0.018446
then the total resulting force, so called Fs[55], generated from collision (in x y z) on node 55 will be the summation of those calculated from constraint 0, 1, 2 right?
2. Secondly, the term F in the equation p=F/S is the summation of all contact forces Fs[i], where i are the impacted nodes, generated in contact area, is it right?
Thank you for any help
best regards 
AuthorPosts
 You must be logged in to reply to this topic.