simtk.org

Dear all,

I am currently reading this article, and I am quite puzzled with the amount of solutions there are to the kinematic problem. I am wondering which one of these OpenSim uses.

- It is definitely a numerical approach, rather than analytical
- There is no learning, so it is not a data-driven approach
- I don't think there is any heuristics involved
- Gauss-Newton optimization converges faster, only if they do. OpenSim is robust enough that I would think it uses a Jacobian approach instead.

However, nowhere in the OpenSim papers nor in the documentation, can I see any mention of the Jacobian. I think I gathered that it is a non-linear least-square problem, solved with global optimization with the Levenberg-Marquard algorithm. How does that imply that this uses a Jacobian method?

Subsidiary question: are forward and direct kinematics really synonym? The 6DoF way of solving joint angles (which is analytical, rather than numerical) seems to be sometimes called direct kinematics.
If so, then I suppose it is wrong to say that inverse kinematics solve joint angles from point coordinates, while forward/direct kinematics solves point coordinates from joint angles?

Thanks

Okay, I did a little further reading, and my question boils down to:
- The math formulations make me thing that global optimization differs from (local, Jacobian) Levenberg-Marquardt optimization. Is it the case?
- If so, how is it that absolutely no robotics/computer animation paper talk about global optimization, and conversely, no biomechanics paper addresses Levenber-Marquardt optimization? Is the latter faster, but also just not accurate enough?

EDIT:
There is actually a Levenverg-Marquardt function in opensim-core/OpenSim/common/lmdif.cpp, but it is only used in the SIMMtoOpenSim functions (which are quite outdated I believe).
And in any case, I don't understand OpenSim/Simulation/InverseKinematicsSolver.cpp, as I don't see any calculation nor optimization of the objective function here... Am I looking in the wrong place?

David,

These are good questions and it's good to see that you want to get to the bottom of this.

In robotics, inverse kinematics is only concerned with the position/orientation of the end effector. If the robot has more than 6 DOF, this is an underdetermined problem, with an infinite number of solutions. This is mentioned in the paper you cited. In animation, they solve the same problem when they need to have a hand or foot in a specified location.

In OpenSim, we want to find the skeleton pose that puts all markers where the motion capture system saw them. Unless the number of markers is very small, this is an overdetermined problem, with no exact solution. Therefore, we look for a least-squares solution. The least squares optimization problem is described here: https://simtk-confluence.stanford.edu:8 ... tics+Works

That page mentions a "general quadratic programming solver" but this is not a quadratic programming problem. Each function x_i(q) is a nonlinear function (it's the forward kinematics of the skeleton, which involves cos, sin, in complex combinations). Sequential quadratic programming (SQP) could be used with such a nonlinear model.

However, SQP won't be the fastest method because it does not take advantage of the special structure of the optimization objective: it is a sum of squares of functions of (q). For such nonlinear least squares problems, the Levenberg-Marquardt problem is the recommended method. Gauss-Newton is simpler and faster but less robust. With the LM or GN methods, the Hessian matrix of the optimization objective is approximated very efficiently from the Jacobian matrix of the functions in the summation [1,3].

Maybe someone from the OpenSim team can clarify which algorithm is actually used.

There is also some confusion about terminology. The early papers on marker-based IK used the term "global optimization" [2]. Here, "global" refers to optimizing the pose of multiple body segments at the same time, rather than one body segment at a time, as was done traditionally with three markers on a segment. This is an unfortunate term, because "global optimization" already has a very specific meaning in the field of optimization! The title also mentions "constraints". This refers to having joints (kinematic constraints) in the model. It's not a constrained optimization when a forward kinematic model is used. "Inverse kinematics" is a better term, but then it gets confused with the traditional end-effector IK problems. More accurate would be "least-squares inverse kinematics".

You would never want to use global optimization for IK, it is too slow, and you may end up in a local optimum. This is somewhat paradoxical, usually we use global optimization to avoid local optima. Musculoskeletal IK has many local optima, where the skeleton is twisted, but those are easy to avoid by using a neutral pose as initial guess for a local (gradient-hessian-based) optimization in the first frame of data. Each subsequent frame then uses the previous frame's solution as an initial guess. Often, you need only one or two iterations of the solver.

Apparently, the term "direct kinematics" is sometimes used to indicate the method where the (6-DOF) pose of a single segment is estimated from 3 or more markers. That's not a great term either, because "direct kinematics" is already a synonym for forward kinematics.

Personal/historical note: I worked on a proprietary marker-based IK code for Motion Analysis Corp. This was not published or patented, but commercially available in 1997 under the name "Mocap Solver". It had a lot of impact in the animation industry and I completely underestimated its potential for biomechanics research. In 1999, Vicon had a similar capability [2]. In the same year, it was independently developed by Zoran Popovic at Carnegie Mellon University who described it in one paragraph in his PhD dissertation. SIMM (a precursor of OpenSim) used Mocap Solver and made it compatible with musculoskeletal models.

--

Ton van den Bogert

[1] Numerical Recipes in Fortran, 2nd edition, 1992, chapter 15.
[2] Lu TW, O'Connor JJ (1999) Bone position estimation from skin marker co-ordinates using global optimisation with joint constraints. J Biomech 32(2):129-134.
[3] van den Bogert et al. (2013) A real-time system for mechanical analysis of human movement. Medical and Biological Engineering and Computing.

Wow, I could not hope for a clearer and more complete answer, thank you so much!
You also address a lot of the questions I was asking myself about ambiguous terms, which I did not post here in order not to make it overly long.

I am still wondering, the page "How inverse kinematics works" that you shared refers to this book for more information about their optimization procedure. Boyd, 2004, Convex Optimization. And there, I see no mention of Levenberg-Marquardt, Jacobian (nor to SQP, to be fair). Is it just a reference that is not appropriate?

You definitely think OpenSim uses a LM solver? What about Quasi-Newton approaches? Is OpenSim too fast to rely on these? Visual 3d apparently has an implementation of this.

Visual3D has 3 implementations of the solution to this optimization problem; in order of computational efficiency Levenberg-Marquardt, Quasi-Newton, and Simulated Annealing. Of particular relevance to this project is our implementation of the Quasi-Newton algorithm because this is a more stable method that we would prefer to use than the Levenberg-Marquardt algorithm, but is unfortunately slow to converge. (The Simulated Annealing algorithm is currently so slow that it is rarely used)
The Quasi-Newton algorithm is based on Newton's method for optimizing twice differentiable functions (Andrew, 2008).

Also, I really need to get better at C++, I can't believe I can't find my answer in the code...

It is extremely convoluted, but I think I slowly understand it a little better.

OpenSim/Simulation/InverseKinematicsSolver.h calls OpenSim/Simulation/AssemblySolver.h

This calls for simbody/internal/Assembler.h (on a completely different repository). According to the doc, it can be used in two ways:
- Basic assembly with assemble(), where the whole model is adjusted, and
- Repeated assembly with track(), where only generalized coordinates are adjusted. This is the one used for inverse kinematics.

Track() calls optimize(), which is defined in SimTKmath/include/simmath/Optimizer.h. However, there are several available algorithms there, and I still can't understand which one is used. The pdf doc gives a little more info, but I still don't fully understand it.

In any case, I would like to know if it is closer to a Levenberg-Marquardt algorithm, to a SQP one, or gradient, or Gauss-Newton, or anything else. And I'm still stuck

I also tried to find the answers in the C++ code and got about as far as you did. Someone must know...

The Boyd reference seems incorrect to me. This inverse kinematics problem is not a convex optimization.

If computation speed is not important, it does not matter which algorithm is used. They will all find the same answer, if the initial guess is good enough.

When I teach inverse kinematics, I keep it simple and have the students use fminsearch() from the Matlab optimization toolbox. It is fast enough for simple problems.

Ton

After having read more about inverse kinematics in biomechanics, and about optimization procedures in general, I would tend to think that the method used by OpenSim is a Quasi-Newton one. A Quasi-Newton method is a Newton method that approximates the Hessian matrix, which is computationally expensive to invert. More precisely, this would use the limited-memory BFGS approximation method, with bound constraints (option #3 here).
It is slower than Levenberg-Marquardt, but not too bad, joint limits are straightforward to implement, and it does not suffer as much from singularity issues.
An SQP method would also work, but I believe it would be slower.

To be confirmed...

Hi Pagnon David,

Sorry I might be speaking out of turn but I found your question very interesting and thought the following post might help : viewtopicPhpbb.php?f=91&t=3854&p=22173&start=0&view=

It is quite old, so I am not sure if it still valid today.

Thanks! I am not really looking into this anymore, but this is definitely a relevant post!