Compute Hessian of MbP dynamics function
If I have a function f(u, q, v, v̇, λ) = M * v̇ + C(q, v) - Bu - J'λ, I wish to compute the Hessian of this function f w.r.t u, q, v, v̇, λ. Some of the Hessian terms are easy to compute, but to compute ∂² f / ∂ q², it requires computing the Hessian of ∂² M / ∂ q² and ∂² C / ∂ q², which I don't know how to do yet.
I guess the easiest solution might be to instantiate MBP with nested AutoDiffScalar.
@sherm1
Since f is a vector function, ∂f / ∂q is a matrix. That makes ∂²f / ∂q² a hypermatrix. Is that what you're looking for?
I should have been more clear, I mean ∂²f(i) / ∂q², namely the second order derivative of the i'th entry of f w.r.t q.
We don't have a way to calculate just a single generalized force f(i) so the nested AutoDiff calculation would give { ∂²f(i)/∂q² | i ∈ 0..n-1 }. (That's a way to represent the 3D hypermatrix.) I'm not seeing an easier way to get that though.
You don't have u in the functional dependence of f. Should it be there?
You don't have u in the functional dependence of f. Should it be there?
Yes u should be there also.
We don't have a way to calculate just a single generalized force f(i) so the nested AutoDiff calculation would give { ∂²f(i)/∂q² | i ∈ 0..n-1 }. (That's a way to represent the 3D hypermatrix.) I'm not seeing an easier way to get that though.
Makes sense. I think if we use nested autodiff, we can still return f as an Eigen::Vector<NestedAutoDiff>, where the i'th entry of this f contains the value of f, the gradient and the Hessian.
Yes, each entry would contain its value, gradient, and Hessian. Added u in the issue description.
BTW @edrumwri wrote a document on AutoDiffing that has a Hessian section.
Assigning to @rpoyner-tri for now as vaguely related to #10991 -- this is an application that would really stretch AutoDiff performance!
Given where our thinking about autodiff is going, I'm handing this back to @hongkai-dai for further consideration.
This is something @shrutigarg914 I would find very useful -- there are various cases where we might want to compute the kinematic Hessian, which this would encompass.
(Of course we can do it symbolically or analytically for now.)