Asymmetric Hessians

I encountered a really counter-intuitive result this past week regarding the Hessians of a simple 2D Harmonic Oscillator:

U(r(x_0, y_0, x_1, y_1), b) = (r-b)^2/2 =(\sqrt{(x_0-x_1)^2+(y_0-y_1)^2}-b)^2/2

It is clear that the energies for r+\delta and r-\delta are the same, and the forces are identical in magnitude except opposite in direction.

However, while the internal Hessian \dfrac{\partial^2 U}{\partial r^2} = 1 and therefore positive semi-definite (PSD) every where w.r.t. r, this is not the case for the 4×4 Cartesian Hessian \dfrac{\partial^2 U}{\partial x_i \partial x_j}. In fact, the Cartesian Hessian is PSD only if r \geq b, that is, it has only positive eigenvalues or zero. We can actually analytically show that the eigenvalues of the Cartesian Hessian are 0, 0, 2 and \lambda where:

\lambda= 2(1-\dfrac{b}{r^2})

So if r < b then the eigenvalue is negative (no longer PSD), and if r>b then the eigenvalue is positive (PSD).

To see this for yourself in sympy, the code is:

import sympy as sp
from sympy.functions import Abs, sqrt
from IPython.display import display, Math
x0, x1, y0, y1, b = sp.symbols('x0 x1 y0 y1 b', real=True)
U = ((sqrt((x0-x1)**2+(y0-y1)**2)-b)**2)/2
sp.init_printing()
args = (x0, y0, x1, y1)
H = sp.Matrix(list(sp.derive_by_array(sp.derive_by_array(U, args), args)))
H = H.reshape(4,4)
H = sp.simplify(H)
H.eigenvals()

Leave a comment

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: