Numerical approximation of principal curvature

curvaturediscrete geometrynumerical linear algebra

I have a surface given by z-values on an xy-grid (a 2D-array of values).
To calculate surface tension, I need to calculate mean curvature in every point.
To calculate mean curvature, I need to calculate principal curvature.
Since principal curvature lies in a plane containing the surface-normal-vector, this is where it gets tricky.

My current approach is very brute force:

  1. take a point and it's surrounding 8 points
  2. calculate gradient at central point
  3. rotate all points so that gradient is now in the xy-plane (moving all points off the grid!)
  4. calculate intersection of triangles defined by these 8 points and a sample of planes normal to the xy-plane
  5. find maximum

This is of course very cumbersome and slow. I'm furthermore not sure how many normal plains I should sample and weather there are special planes that give better estimates (e.g. planes containing the rotated points).
Is there a more elegant way to do this?

This may be a duplicate: Numerical computation of surface curvature
However, that question didn't receive an answer.

Best Answer

Finally found the time to do the legwork and calculate the curvature by hand.

Curvature defined as $\nabla \vec{e}_n$ can be calculated with a suitably defined $\vec n$. Having a grid of Z-values, a vector normal to the surface can be found by using the gradient: $$\vec n = \binom{1}{\vec \nabla Z } $$ With $\vec \nabla z$ supplying $x$ and $y$ values of $\vec n$ and $z=1$. Thereby: $$\vec e _n = {\vec{n}\over{|\vec n|}} = {\vec{n}\over \sqrt{1+(\vec \nabla Z)^2}}$$

Since the z-component doesn't depend on z, the divergence is nearly the same, the only difference being the added one in the square root. Therefore, the curvature will be overestimated slightly, if the method portrayed in this linked paper (Chapter 4) is used.

Related Question