Curvature and geodesics

Curvatures

Mathematical definitions

and the Mean Curvature is 
 

where S is the Shape operator (negative derivative of the unit vector field, called also Second Fundamental Tensor), and k1(p), k2(p) are the Principal Curvatures (max and min of the normal curvature of a surface).

See also John Oprea, Differential Geometry and its applications, Prentice Hall, 1997
 

More intuitively

• Normal curvature is similar to the 2D case, if you slice the surface with a plane in a given direction. It is a normal acceleration for the curve. The mean curvature is the mean of the two extremes.
• Gaussian curvature is the ratio between the area of a patch on the surface and the area of the same solid angle on the Gaussian sphere.

Geodesics

Definition

See reference book for equations
 

Maple demo

with (plots):

dp := proc(X,Y)
        X[1] * Y[1] + X[2] * Y[2] + X[3] * Y[3];
     end:

nrm := proc(X)
        sqrt( dp(X,X) );
        end:
 

xp := proc(X,Y)
        local a,b,c;
        a = X[2] * Y[3] - X[3] * Y[2];
        b = X[3] * Y[1] - X[1] * Y[3];
        c = X[1] * Y[2] - X[2] * Y[1];
        [a,b,c];
        end:

Jacf := proc(X)
        local Xu, Xv;
        Xu := [ diff(X[1],u) , diff(X[2],u), diff(X[3],u) ];
        Xv := [ diff(X[1],v) , diff(X[2],v), diff(X[3],v) ];
        simplify( [Xu,Xv] );
        end:

EFG := proc(X)
        local E,F,G,Y;
        Y := Jacf(X);
        E := dp( Y[1], Y[1] );
        F := dp( Y[1], Y[2] );
        G := dp( Y[2], Y[2] );
        simplify( [E,F,G] );
        end:
 

geoeq := proc(X)
        local M, eq1, eq2;
        M := EFG( X );
        eq1 := diff(u(t), t$2) + subs({u=u(t), v=v(t)},
                diff(M[1], u)/(2*M[1])) * diff(u(t),t)^2
                +subs({u=u(t), v=v(t)}, diff(M[1],v)/(M[1]))
                *diff(u(t),t)* diff(v(t),t)
                -subs({u=u(t), v=v(t)}, diff(M[3],u)/(2*M[1]))
                *diff(v(t),t)^2 = 0;
        eq2 :=  diff(v(t),t$2) - subs({u=u(t), v=v(t)},
                diff(M[1], v)/(2*M[3])) * diff(u(t),t)^2
                +subs({u=u(t), v=v(t)}, diff(M[3],u)/(M[3]))
                *diff(u(t),t) * diff(v(t),t)
                +subs({u=u(t), v=v(t)}, diff(M[3],v)/(2*M[3]))
                *diff(v(t),t)^2 = 0;
        eq1, eq2;
        end:
 

plotgeo := proc( X, ustart, uend, vstart, vend, u0, v0,
                Du0, Dv0, n, gr, theta, phi )
        local sys, desys, dequ, deqv, listp, j, geo, plotX;
        sys:= geoeq(X);
        desys := dsolve( {sys, u(0) =u0, v(0) = v0,
                        D(u)(0) = Du0, D(v)(0) = Dv0 },
                        { u(t), v(t)}, type = numeric,
                        output = listprocedure );
        dequ := subs(desys, u(t) ); deqv := subs(desys, v(t) );
        listp := [seq(subs({u=dequ(j/n[3]), v= deqv(j/n[3]) },X),
                j=n[1]..n[2] ) ];
        geo := spacecurve( {listp}, color=black, thickness=2 ):
        plotX := plot3d( X, u=ustart..uend, v=vstart..vend,
                                grid=[gr[1],gr[2]] ) :
        display( {geo, plotX}, style = wireframe, scaling = constrained,
                        orientation = [theta,phi ] );
        end:
 
 

Tor := [(5+cos(u)) * cos(v), (5+cos(u)) * sin(v), sin(u) ];

geoeq( Tor );
 

plotgeo(Tor, 0,2*Pi, 0, 2*Pi, Pi/2, 0, 0, 1, [0,240,15], [20,20], 0, 60 );

Huge thanks to Neel !
 

Eric Amram, Feb 1998, MIT 6.838 Advanced Topics in Computer Graphics - Prof. Seth Teller.