A Quick and Dirty Introduction to the Geometry of Surfaces

There are many ways to think about the geometry of a surface (using charts, for instance) but here’s a picture that is similar in spirit to the way we work with surfaces in the discrete setting. Imagine that our surface is a topological disk. Its geometry can be described via a map \(f: M \rightarrow \mathbb{R}^3\) from a region \(M\) in the Euclidean plane \(\mathbb{R}^2\) to a subset \(f(M)\) of \(\mathbb{R}^3\):

The differential of such a map, denoted by \(df\), tells us how to map a vector \(X\) in the plane to the corresponding vector \(df(X)\) on the surface. Loosely speaking, imagine that \(M\) is a rubber sheet and \(X\) is a little black line segment drawn on \(M\). As we stretch and deform \(M\) into \(f(M)\), the segment \(X\) also gets stretched and deformed into a different segment, which we call \(df(X)\). Later on we can talk about how to explicitly express \(df(X)\) in coordinates and so on, but it’s important to realize that fundamentally there’s nothing deeper to know about the differential than the picture you see here — the differential simply tells you how to stretch out or “push forward” vectors as you go from one space to another. For example, the length of a tangent vector \(X\) pushed forward by \(f\) can be expressed as

\[ \sqrt{ df(X) \cdot df(X) }, \]

where \(\cdot\) is the standard inner product (a.k.a. dot product or scalar product) on \(\mathbb{R}^3\). More generally, we can measure the inner product between any two tangent vectors \(df(X)\) and \(df(Y)\):

\[ g(X,Y) = df(X) \cdot df(Y). \]

The map \(g\) is called the metric of the surface, or to be more pedantic, the metric induced by \(f\).

So far we’ve been talking about tangent vectors, i.e., vectors that travel “along” the surface. We’re also interested in vectors that are orthogonal to the surface. In particular, we say that a vector \(u \in \mathbb{R}^3\) is normal to the surface at a point \(p\) if

\[ df(X) \cdot u = 0 \]

for all tangent vectors \(X\) at \(p\). For convenience, we often single out a particular normal vector \(\hat{N}\) called the unit normal, which has length one. Of course, at any given point there are two distinct unit normal vectors: \(+\hat{N}\) and \(-\hat{N}\). Which one should we use? If we can pick a consistent direction for \(\hat{N}\) then we say that \(M\) is orientable. For instance, the torus is orientable, but the Möbius band is not:

For orientable surfaces, we can actually think of \(\hat{N}\) as a continuous map \(\hat{N}: M \rightarrow S^2\) (called the Gauss map) which associates each point with its unit normal, viewed as a point on the unit sphere \(S^2 \subset \mathbb{R}^3\). In fact, why not think of \(\hat{N}\) as simply a different geometry for \(M\)? Now the map \(d\hat{N}\) (called the shape operator) tells us about the change in normal direction as we move from one point to the other. For instance, we could look at the change in normal along a particular tangent direction \(X\) by evaluating

\[ \kappa_n(X) = \frac{df(X) \cdot d\hat{N}(X)}{|df(X)|^2}. \]

(The factor \(|df(X)|^2\) in the denominator simply accounts for any “stretching out” of \(X\) that occurs as we go from the plane to the surface.) The quantity \(\kappa_n\) is called the normal curvature — we’ll have a lot more to say about curvature in the future.

Overall we end up with this picture, which captures the most fundamental ideas about the geometry of surfaces:

Conformal Parameterization

When we talked about curves, we introduced the idea of an isometric (a.k.a. arc-length or unit speed) parameterization. The idea there was to make certain expressions simpler by assuming that no “stretching” occurs as we go from the domain to \(\mathbb{R}^3\). Another way of saying this is that

\[ |df(X)| = |X|, \]

i.e., the norm of any vector \(X\) is preserved.

For surfaces, an isometric parameterization does not always exist (not even locally!). Sometimes you simply have to stretch things out. For instance, you may know that it’s impossible to flatten the surface of the Earth onto the plane without distortion — that’s why we end up with all sorts of different funky projections of the globe.

However, there is a setup that (like arc-length parameterization for curves) makes life a lot easier when dealing with certain expressions, namely conformal parameterization. Put quite simply, a map \(f\) is conformal if it preserves the angle between any two vectors. More specifically, a conformal map \(f: \mathbb{R}^2 \supset M \rightarrow \mathbb{R}^3\) satisfies

\[ df(X) \cdot df(Y) = a \langle X, Y \rangle \]

for all tangent vectors \(X, Y\), where \(a > 0\) is a positive scale factor and \(\langle \cdot, \cdot \rangle\) is the usual inner product on \(\mathbb{R}^2\). Vectors can still get stretched out, but the surface never gets sheared — for instance, orthogonal vectors always stay orthogonal:

A key fact about conformal maps is that they always exist, as guaranteed by the uniformization theorem. In a nutshell, the uniformization theorem says that any disk can be conformally mapped to the plane. So if we consider any point \(p\) on our surface \(f(M)\), we know that we can always find a conformal parameterization in some small, disk-like neighborhood around \(p\). As with unit-speed curves, it is often enough to simply know that a conformal parameterization exists — we do not have to construct the map explicitly. And, as with arc-length parameterization, we have to keep track of the least possible amount of information about how the domain gets stretched out: just a single number at each point (as opposed to, say, an entire Jacobian matrix).

The Fundamental Polygon

We’ve received a couple questions on homework problem 2.2 regarding whether a constructive proof is necessary, i.e., do you need to actually come up with examples of triangulations that achieve the desired bounds? The answer is no: you simply need to show (algebraically) that these bounds hold.

But if you were looking for a constructive proof, a nice tool to know about is something called the fundamental polygon. The basic idea is that a torus of any genus can be cut into a single disk, which can be visualized as a polygon in the plane. Visualizing a triangulation (or other data) on the fundamental polygon is often much simpler than visualizing it on the embedded torus — consider this example, for instance:

 

 

Here the arrows mark identifications of edges — for instance, the top and bottom edges get “glued” together along the positive \(x\)-axis. (Try convincing yourself that these gluings really do produce a torus — note that at some point you’ll need to make a \(180^\circ\) “twist.”) One way to visualize the identifications is to imagine that the fundamental polygon tiles the plane:

The torus is interesting because it actually admits two different funamental polygons: the hexagon and the square (corresponding to two different tilings of the Euclidean plane). So we could also visualize the torus on the square, leading to an even simpler triangulation:

 

 

(By the way, is this really what you’d call a triangulation? Each region certainly has three sides, but each “triangle”‘ has only one vertex! What about the first example? There all triangles have three vertices, but they share the same three vertices. So combinatorially these triangles are not distinct — another way of saying this is that they do not form a valid simplicial complex. What’s the smallest simplicial decomposition you can come up with for the torus?)

 

In general the fundamental polygon for a torus of genus \(g\) has \(4g\) sides with identifications \(a_1 b_1 a_1^{-1} b_1^{-1} \cdots a_n b_n a^{-1}_n b_n\), where two edges with the same letter get identified and inverse just means that the edge direction is reversed. For instance, the fundamental polygon for the double torus looks like this:

 

 

(Note that all the polygon vertices ultimately get identified with a single point on the surface.) From here it becomes easy to start playing around with tessellations — for instance, here’s how you decompose a surface of any genus into quadrilaterals, using only two irregular vertices (can you do better?):

 

 

Tiling the Euclidean plane with the fundamental polygon is impossible for a surface of genus \(g \geq 2\), since the interior angles of the fundamental polygon don’t evenly divide \(2\pi\) (proof!). Fortunately, we can still tile the hyperbolic plane, i.e., the plane with constant negative curvature. For instance, here’s a tiling of the hyperbolic plane by octagons:

 

 

From here there are all sorts of fascinating things to say about covering spaces, uniformization, and especially the fundamental group of a surface — if you’re interested I highly recommend you take a look at Allen Hatcher’s (free!) book on algebraic topology.

A Quick and Dirty Introduction to the Geometry of Curves

The picture we looked at for surfaces is actually a nice way of thinking about manifolds of any dimension. For instance, we can think of a one-dimensional curve as a map \(\gamma: I \rightarrow \mathbb{R}^3\) from a subset \(I = [0,T] \subset \mathbb{R}\) of the real line to \(\mathbb{R}^3\). Again the differential \(d\gamma\) tells us how tangent vectors get stretched out by \(\gamma\), and again the induced length of a tangent vector \(X\) is given by

\[ |d\gamma(X)| = \sqrt{d\gamma(X) \cdot d\gamma(X)}. \]

Working with curves is often easier if \(\gamma\) preserves length, i.e., if for every tangent vector \(X\) we have

\[ |d\gamma(X)| = |X|. \]

There are various names for such a parameterization (“unit speed”, “arc-length”, “isometric”) but the point is simply that the curve doesn’t get stretched out when we go from \(\mathbb{R}\) to \(\mathbb{R}^3\) — think of \(\gamma\) as a rubber band that is completely relaxed. This unit-speed view is also the right one for the discrete setting where we often have no notion of a base domain \(I\) — the curve is “born” in \(\mathbb{R}^3\) and all we can do is assume that it sits there in a relaxed state.

The Frenet Frame

Suppose we have a unit-speed curve \(\gamma\) and a unit vector \(\hat{X} \in \mathbb{R}\) oriented along the positive direction. Then

\[ T = d\gamma(\hat{X}) \]

is a unit vector in \(\mathbb{R}^3\) tangent to the curve. Carrying this idea one step further, we can look at the change in tangent direction as we move along \(\gamma\). Since \(T\) may change at any rate (or not at all!) we split up the change into two pieces:

\[ \kappa N = dT(\hat{X}). \]

The scalar \(\kappa \in \mathbb{R}\) is called the curvature and expresses the magnitude of change, while the vector \(N \in \mathbb{R}^3\) is called the principal normal and expresses the direction of change.

One thing to notice is that \(T\) and \(N\) are always orthogonal. Why? Simply because if \(T\) changed along the direction of \(T\) itself then it would no longer have unit length! Therefore one often defines a third vector \(B = T \times N\) called the binormal, which completes the orthonormal Frenet frame \((T,N,B)\) of \(\gamma\).

Finally, how does \(B\) change along \(\gamma\)? Here our story ends since the change in the binormal can be expressed in terms of the normal: \(dB(\hat{X}) = -\tau N\). The quantity \(\tau \in \mathbb{R}\) is called the torsion, and describes how much the Frenet frame “twists” around \(T\) as we travel along \(\gamma\).

Note that for a curve with torsion the normal \(N\) also twists around the tangent — in particular \(dN = -\kappa T + \tau B\). Altogether we get the Frenet-Serret formulas, which describe the motion of the Frenet frame:

\[
\begin{array}{rcc}
dT &=& \kappa N \\
dN &=& -\kappa T + \tau B \\
dB &=& -\tau N
\end{array}
\]

Since \(T\), \(N\), and \(B\) have unit length at each point, we can visualize them as curves on the unit sphere \(S^2 \subset \mathbb{R}^3\), just as we did with the Gauss map. In this context the sphere is often given a special name: the (tangent-, principal normal-, or binormal-) indicatrix. Later on we’ll see that this “indicatrix” picture has some beautiful geometric consequences.

Visualizing Curvature

Imagine that \(\gamma\) sits in the plane. At any given point there is a circle \(S\) called the osculating circle that best approximates \(\gamma\), meaning that it has the same tangent direction \(T\) and curvature vector \(\kappa N\). How do we know such a circle exists? Well, we can certainly find a point with tangent \(T\) on a circle of any radius. We can then adjust the size of \(S\) until it has the desired curvature. Alternatively, one can construct \(S\) by considering a circle passing through three consecutive points on \(\gamma\), which approaches the osculating circle \(S\) as the points move closer together. Since these points yield the same approximations of first- and second- derivatives along both \(\gamma\) and \(S\) (using finite differences), the tangents and curvatures of the two curves will agree as the points converge.

In any case, we can use the geometry of \(S\) to express the curvature of \(\gamma\). If \(S\) has radius \(r\) then it takes time \(2\pi r\) to go all the way around the circle at unit speed, and during this time the tangent turns around by an angle \(2\pi\). Of course, since \(T\) has unit length the instantaneous change in \(T\) is described exclusively by the change in angle. So we end up with \(\kappa = |\kappa N| = |dT(\hat{X})| = 2\pi/2\pi r = 1/r\): the curvature of a circle is simply the reciprocal of its radius. This fact should make some intuitive sense: if we watch a circle grow bigger and bigger from a fixed viewpoint, it eventually looks just like a straight line.

The radius and center of the osculating circle are often referred to as the radius of curvature and center of curvature, respectively. We can tell this same story for any curve in \(\mathbb{R}^3\) by considering the osculating plane \(T \times N\), since this plane contains both the tangent and the curvature vector.

For curves it makes little difference whether we consider the tilt of the tangent vector to find the curvature or the tilt of the (principal) normal vector since the two are related via a rotation by \(\pi/2\) in the osculating plane. For surfaces, however, it will make more sense to consider the change in normal vector, since we typically don’t have a distinguished tangent vector at any given point.

October 11, 2011 | Comments Closed

A Quick and Dirty Introduction to the Curvature of Surfaces

Let’s take a more in-depth look at the curvature of surfaces. The word “curvature” really corresponds to our everyday understanding of what it means for something to be curved: eggshells, donuts, and spiraling telephone cables have a lot of curvature; floors, ceilings, and cardboard boxes do not. But what about something like a beer bottle? Along one direction the bottle quickly curves around in a circle; along another direction it’s completely flat and travels along a straight line:

This way of looking at curvature — in terms of curves traveling along the surface — is often how we treat curvature in general. In particular, let \(X\) be a unit tangent direction at some distinguished point on the surface, and consider a plane containing both \(df(X)\) and the corresponding normal \(N\). This plane intersects the surface in a curve, and the curvature \(\kappa_n\) of this curve is called the normal curvature in the direction \(X\):

Remember the Frenet-Serret formulas? They tell us that the change in the normal along a curve is given by \(dN = -\kappa T + \tau B\). We can therefore get the normal curvature along \(X\) by extracting the tangential part of dN, and normalizing by any “stretching out” that occurs as we go from the domain \(M\) into \(\mathbb{R}^3\):

\[ \kappa_n(X) = \frac{-df(X) \cdot dN(X)}{|df(x)|^2}. \]

Note that normal curvature is signed, meaning the surface can bend toward the normal or away from it.

Principal, Mean, and Gaussian Curvature

At any given point we can ask: along which directions does the surface bend the most? The unit vectors \(X_1\) and \(X_2\) along which we find the maximum and minimum normal curvatures \(\kappa_1\) and \(\kappa_2\) are called the principal directions; the curvatures \(\kappa_i\) are called the principal curvatures. For instance, the beer bottle above might have principal curvatures \(\kappa_1 = 1\), \(\kappa_2 = 0\) at the marked point.

One important fact about the principal directions is that they are eigenvectors of the shape operator, in the sense that

\[ dN(X_i) = \kappa_i df(X_i). \]

(Moreover, the principal directions are orthogonal: \(df(X_1) \cdot df(X_2) = 0\) — see the appendix below for a proof.) The principal directions therefore tell us everything there is to know about the curvature at a point, since we can express any tangent vector \(Y\) in the basis \(\{X_1,X_2\}\) and easily compute the corresponding normal curvature. As we’ll discover, however, getting your hands directly on the principal curvatures is often quite difficult — especially in the discrete setting.

On the other hand, two closely related quantities — called the mean curvature and the Gaussian curvature will show up over and over again, and have some particularly nice interpretations in the discrete world. The mean curvature \(H\) is the arithmetic mean of principal curvatures:

\[ H = \frac{\kappa_1 + \kappa_2}{2}, \]

and the Gaussian curvature is the (square of the) geometric mean:

\[ K = \kappa_1 \kappa_2. \]

What should you visualize if someone tells you a surface has nonzero Gaussian or mean curvature? Perhaps the most elementary interpretation is that Gaussian curvature is like a logical “and” (is there curvature along both directions?) whereas mean curvature is more like a logical “or” (is there curvature along at least one direction?) Of course, you have to be a little careful here since you can also get zero mean curvature when \(\kappa_1 = -\kappa_2\).

It also helps to see pictures of surfaces with zero mean and Gaussian curvature. Zero-curvature surfaces are so thoroughly studied in mathematics that they even have special names. Surfaces with zero Gaussian curvature are called developable surfaces because they can be “developed” or flattened out into the plane without any stretching or tearing. For instance, any piece of a cylinder is developable since one of the principal curvatures is zero:

Surfaces with zero mean curvature are called minimal surfaces because (as we’ll see later) they minimize surface area (with respect to certain constraints). Minimal surfaces tend to be saddle-like since principal curvatures have equal magnitude but opposite sign:

The saddle is also a good example of a surface with negative Gaussian curvature. What does a surface with positive Gaussian curvature look like? The hemisphere is one example:

Note that in this case \(\kappa_1 = \kappa_2\) and so principal directions are not uniquely defined — maximum (and minimum) curvature is achieved along any direction \(X\). Any such point on a surface is called an umbilic point.

There are plenty of cute theorems and relationships involving curvature, but those are the basic facts: the curvature of a surface is completely characterized by the principal curvatures, which are the maximum and minimum normal curvatures. The Gaussian and mean curvature are simply averages of the two principal curvatures, but (as we’ll see) are often easier to get your hands on in practice.

October 12, 2011 | Comments Closed

A Quick and Dirty Introduction to Exterior Calculus — Part I: Vectors and 1-Forms

Many important concepts in differential geometry can be nicely expressed in the language of exterior calculus. Initially these concepts will look exactly like objects you know and love from vector calculus, and you may question the value of giving them funky new names. For instance, scalar fields are no longer called scalar fields, but are now called 0-forms! In many ways vector and exterior calculus are indeed “dual” to each-other, but it is precisely this duality that makes the language so expressive. In the long run we’ll see that exterior calculus also makes it easy to generalize certain ideas from vector calculus — the primary example being Stokes’ theorem. Actually, we already started using this language in our introduction to the geometry of surfaces, but here’s the full story.

Once upon a time there was a vector named \(v\):

What information does \(v\) encode? One way to inspect a vector is to determine its extent or length along a given direction. For instance, we can pick some arbitrary direction \(\alpha\) and record the length of the shadow cast by \(v\) along \(\alpha\):

The result is simply a number, which we can denote \(\alpha(v)\). The notation here is meant to emphasize the idea that \(\alpha\) is a function: in particular, it’s a linear function that eats a vector and produces a scalar. Any such function is called a 1-form (a.k.a. a covector or a cotangent).

Of course, it’s clear from the picture that the space of all 1-forms looks a lot like the space of all vectors: we just had to pick some direction to measure along. But often there is good reason to distinguish between vectors and 1-forms — the distinction is not unlike the one made between row vectors and column vectors in linear algebra. For instance, even though rows and column both represent “vectors,” we only permit ourselves to multiply rows with columns:

\[ \left[ \begin{array}{ccc} \alpha_1 & \cdots & \alpha_n \end{array} \right] \left[ \begin{array}{c} v_1 \\ \vdots \\ v_n \end{array} \right]. \]

If we wanted to multiply, say, two columns, we would first have to take the transpose of one of them to convert it into a row:

\[ v^T v = \left[ \begin{array}{ccc} v_1 & \cdots & v_n \end{array} \right] \left[ \begin{array}{c} v_1 \\ \vdots \\ v_n \end{array} \right]. \]

Same deal with vectors and 1-forms, except that now we have two different operations: sharp (\(\sharp\)), which converts a 1-form into a vector, and flat (\(\flat\)) which converts a vector into a 1-form. For instance, it’s perfectly valid to write \(v^\flat(v)\) or \(\alpha(\alpha^\sharp)\), since in either case we’re feeding a vector to a 1-form. The operations \(\sharp\) and \(\flat\) are called the musical isomorphisms.

All this fuss over 1-forms versus vectors (or even row versus column vectors) may seem like much ado about nothing. And indeed, in a flat space like the plane, the difference between the two is pretty superficial. In curved spaces, however, there’s an important distinction between vectors and 1-forms. For one thing, physical quantities “without mass” (like velocity and acceleration) are represented by vectors; quantities “with mass” (like momentum and force) are represented by 1-forms.

More generally we want to make sure that we’re taking “measurements” in the right space. For instance, suppose we want to measure the length of a vector \(v\) along the direction of another vector \(u\). It’s important to remember that tangent vectors get stretched out by the map \(f: \mathbb{R}^2 \supset M \rightarrow \mathbb{R}^3\) that takes us from the plane to some surface in \(\mathbb{R}^3\). Therefore, the operations \(\sharp\) and \(\flat\) should satisfy relationships like

\[ u^\flat(v) = g(u,v) \]

where \(g\) is the metric induced by \(f\). This way we’re really measuring how things behave in the “stretched out” space rather than the initial domain \(M\).

Coordinates

Until now we’ve intentionally avoided the use of coordinates — in other words, we’ve tried to express geometric relationships without reference to any particular coordinate system \(x_1, \ldots, x_n\). Why avoid coordinates? Several reasons are often cited (people will mumble something about “invariance”), but the real reason is quite simply that coordinate-free expressions tend to be shorter, sweeter, and easier to extract meaning from. This approach is also particularly valuable in geometry processing, because many coordinate-free expressions translate naturally to basic operations on meshes.

Yet coordinates are still quite valuable in a number of situations. Sometimes there’s a special coordinate basis that greatly simplifies analysis — recall our discussion of principal curvature directions, for instance. At other times there’s simply no obvious way to prove something without coordinates. For now we’re going to grind out a few basic facts about exterior calculus in coordinates; at the end of the day we’ll keep whatever nice coordinate-free expressions we find and politely forget that coordinates ever existed!

Let’s setup our coordinate system. For reasons that will become clear later, we’re going to use the symbols \(\frac{\partial}{\partial x^1}, \ldots, \frac{\partial}{\partial x^n}\) to represent an orthonormal basis for vectors in \(\mathbb{R}^n\), and use \(dx^i, \ldots, dx^n\) to denote the corresponding 1-form basis. In other words, any vector \(v\) can be written as a linear combination

\[ v = v^1 \frac{\partial}{\partial x^1} + \cdots + v^n \frac{\partial}{\partial x^n}, \]

and any 1-form can be written as a linear combination

\[ \alpha = \alpha_1 dx^1 + \cdots + \alpha_n dx^n. \]

To keep yourself sane at this point, you should completely ignore the fact that the symbols \(\frac{\partial}{\partial x^i}\) and \(dx^i\) look like derivatives — they’re simply collections of unit-length orthogonal bases, as depicted above. The two bases \(dx^i\) and \(\frac{\partial}{\partial x^i}\) are often referred to as dual bases, meaning they satisfy the relationship

\[ dx^i\left(\frac{\partial}{\partial x^j}\right) = \delta^i_j = \begin{cases} 1, & i = j \\ 0, & \mbox{otherwise.} \end{cases} \]

This relationship captures precisely the behavior we’re looking for: a vector \(\frac{\partial}{\partial x^i}\) “casts a shadow” on the 1-form \(dx^j\) only if the two bases point in the same direction. Using this relationship, we can work out that

\[ \alpha(v) = \sum_i \alpha_i dx^i\left( \sum_j v^j \frac{\partial}{\partial x^j} \right) = \sum_i \alpha_i v_i \]

i.e., the pairing of a vector and a 1-form looks just like the standard Euclidean inner product.

Notation

It’s worth saying a few words about notation. First, vectors and vector fields tend to be represented by letters from the end of the Roman alphabet (\(u\), \(v\), \(w\) or \(X\), \(Y\), \(Z\), repectively), whereas 1-forms are given lowercase letters from the beginning of the Greek alphabet (\(\alpha\), \(\beta\), \(\gamma\), etc.). Although one often makes a linguistic distinction between a “vector” (meaning a single arrow) and a “vector field” (meaning an arrow glued to every point of a space), there’s an unfortunate precedent to use the term “1-form” to refer to both ideas — sadly, nobody ever says “1-form field!” Scalar fields or 0-forms are often given letters from the middle of the Roman alphabet (\(f\), \(g\), \(h\)) or maybe lowercase Greek letters from somewhere in the middle (\(\phi\), \(\psi\), etc.).

You may also notice that we’ve been very particular about the placement of indices: coefficients \(v^i\) of vectors have indices up, coefficients \(\alpha_i\) of 1-forms have indices down. Similarly, vector bases \(\frac{\partial}{\partial x^i}\) have indices down (they’re in the denominator), and 1-form bases \(dx^i\) have indices up. The reason for being so neurotic is to take advantage of Einstein summation notation: any time a pair of variables is indexed by the same letter \(i\) in both the “up” and “down” position, we interpret this as a sum over all possible values of \(i\):

\[ \alpha_i v^i = \sum_i \alpha_i v^i. \]

The placement of indices also provides a cute mnemonic for the musical isomorphisms \(\sharp\) and \(\flat\). In musical notation \(\sharp\) indicates a half-step increase in pitch, corresponding to an upward movement on the staff. For instance, both notes below correspond to a “C” with the same pitch:

Therefore, to go from a 1-form to a vector we raise the indices. For instance, in a flat space we don’t have to worry about the metric and so a 1-form

\[ \alpha = \alpha_1 dx^1 + \cdots + \alpha_n dx^n \]

becomes a vector

\[ \alpha^\sharp = \alpha^1 \frac{\partial}{\partial x^1} + \cdots + \alpha^n \frac{\partial}{\partial x^n}. \]

Similarly, \(\flat\) indicates a decrease in pitch and a downward motion on the staff:

and so \(\flat\) lowers the indices of a vector to give us a 1-form — e.g.,

\[ v = v^1 \frac{\partial}{\partial x^1} + \cdots + v^n \frac{\partial}{\partial x^n}. \]

becomes

\[ v^\flat = v_1 dx^1 + \cdots + v_n dx^n. \]

October 21, 2011 | Comments Closed

A Quick and Dirty Introduction to Exterior Calculus — Part II: Differential Forms and the Wedge Product

In our last set of notes we measured the length of a vector by projecting it onto different coordinate axes; this measurement process effectively defined what we call a 1-form. But what happens if we have a collection of vectors? For instance, consider a pair of vectors \(u, v\) sitting in \(\mathbb{R}^3\):

We can think of these vectors as defining a parallelogram, and much like we did with a single vector we can measure this parallelogram by measuring the size of the “shadow” it casts on some plane:

For instance, suppose we represent this plane via a pair of unit orthogonal 1-forms \(\alpha\) and \(\beta\). Then the projected vectors have components

\[ \begin{array}{rcl}
u^\prime &=& (\alpha(u),\beta(u)), \\
v^\prime &=& (\alpha(v),\beta(v)),
\end{array} \]

hence the (signed) projected area is given by the cross product

\[ u^\prime \times v^\prime = \alpha(u)\beta(v) - \alpha(v)\beta(u). \]

Since we want to measure a lot of projected volumes in the future, we’ll give this operation the special name “\(\alpha \wedge \beta\)”:

\[ \alpha \wedge \beta(u,v) := \alpha(u)\beta(v) - \alpha(v)\beta(u). \]

As you may have already guessed, \(\alpha \wedge \beta\) is what we call a 2-form. Ultimately we’ll interpret the symbol \(\wedge\) (pronounced “wedge”) as a binary operation on differential forms called the wedge product. Algebraic properties of the wedge product follow directly from the way signed volumes behave. For instance, notice that if we reverse the order of our axes \(\alpha, \beta\) the sign of the area changes. In other words, the wedge product is antisymmetric:

\[ \alpha \wedge \beta = -\beta \wedge \alpha. \]

An important consequence of antisymmetry is that the wedge of any 1-form with itself is zero:

\[ \begin{array}{c}
\alpha \wedge \alpha = -\alpha \wedge \alpha \\
\Rightarrow \alpha \wedge \alpha = 0.
\end{array} \]

But don’t let this statement become a purely algebraic fact! Geometrically, why should the wedge of two 1-forms be zero? Quite simply because it represents projection onto a plane of zero area! (I.e., the plane spanned by \(\alpha\) and \(\alpha\).)

Next, consider the projection onto two different planes spanned by \(\alpha, \beta\) and \(\alpha, \gamma\). The sum of the projected areas can be written as

\[
\begin{array}{rcl}
\alpha \wedge \beta(u,v) + \alpha \wedge \gamma(u,v)
&=& \alpha(u)\beta(v) - \alpha(v)\beta(u) + \alpha(u)\gamma(v) - \alpha(v)\gamma(u) \\
&=& \alpha(u)(\beta(v) + \gamma(v)) - \alpha(v)(\beta(u) + \gamma(u)) \\
&=:& (\alpha \wedge (\beta + \gamma))(u,v),
\end{array}
\]

or in other words \(\wedge\) distributes over \(+\):

\[ \alpha \wedge (\beta + \gamma) = \alpha \wedge \beta + \alpha \wedge \gamma. \]

Finally, consider three vectors \(u, v, w\) that span a volume in \(\mathbb{R}^3\):

We’d like to consider the projection of this volume onto the volume spanned by three 1-forms \(\alpha\), \(\beta\), and \(\gamma\), but the projection of one volume onto another is a bit difficult to visualize! For now you can just cheat and imagine that \(\alpha = dx^1\), \(\beta = dx^2\), and \(\gamma = dx^3\) so that the mental picture for the projected volume looks just like the volume depicted above. One way to write the projected volume is as the determinant of the projected vectors \(u^\prime\), \(v^\prime\), and \(w^\prime\):

\[ \alpha \wedge \beta \wedge \gamma( u, v, w ) := \mathrm{det}\left(\left[ \begin{array}{ccc} u^\prime & v^\prime & w^\prime \end{array} \right]\right) = \mathrm{det}\left( \left[ \begin{array}{ccc} \alpha(u) & \alpha(v) & \alpha(w) \\ \beta(u) & \beta(v) & \beta(w) \\ \gamma(u) & \gamma(v) & \gamma(w) \end{array} \right] \right). \]

(Did you notice that the determinant of the upper-left 2×2 submatrix also gives us the wedge product of two 1-forms?) Alternatively, we could express the volume as the area of one of the faces times the length of the remaining edge:

Thinking about things this way, we might come up with an alternative definition of the wedge product in terms of the triple product:

\[
\begin{array}{rcl}
\alpha \wedge \beta \wedge \gamma( u, v, w ) &=& (u^\prime \times v^\prime) \cdot w^\prime \\
&=& (v^\prime \times w^\prime) \cdot u^\prime \\
&=& (w^\prime \times u^\prime) \cdot v^\prime \\
\end{array}
\]

The important thing to notice here is that order is not important — we always get the same volume, regardless of which face we pick (though we still have to be a bit careful about sign). A more algebraic way of saying this is that the wedge product is associative:

\[ (\alpha \wedge \beta) \wedge \gamma = \alpha \wedge (\beta \wedge \gamma). \]

In summary, the wedge product of \(k\) 1-forms gives us a \(k\)-form, which measures the projected volume of a collection of \(k\) vectors. As a result, the wedge product has the following properties for any \(k\)-form \(\alpha\), \(l\)-form \(\beta\), and \(m\)-form \(\gamma\):

  • Antisymmetry: \(\alpha \wedge \beta = (-1)^{kl}\beta \wedge \alpha\)
  • Associativity: \(\alpha \wedge (\beta \wedge \gamma) = (\alpha \wedge \beta) \wedge \gamma\)
and in the case where \(l=m\) we have
  • Distributivity: \(\alpha \wedge (\beta + \gamma) = \alpha \wedge \beta + \alpha \wedge \gamma\)

A separate fact is that a \(k\)-form is antisymmetric in its arguments — in other words, swapping the relative order of two “input” vectors changes only the sign of the volume. For instance, if \(\alpha\) is a 2-form then \(\alpha(u,v) = -\alpha(v,u)\). In general, an even number of swaps will preserve the sign; an odd number of swaps will negate it. (One way to convince yourself is to consider what happens to the determinant of a matrix when you exchange two of its columns.) Finally, you’ll often hear people say that \(k\)-forms are “multilinear” — all this means is that if you keep all but one of the vectors fixed, then a \(k\)-form looks like a linear map. Geometrically this makes sense: \(k\)-forms are built up from \(k\) linear measurements of length (essentially just \(k\) different dot products).

Vector-Valued Forms

Up to this point we’ve considered only real-valued \(k\)-forms — for instance, \(\alpha(u)\) represents the length of the vector \(u\) along the direction \(\alpha\), which can be expressed as a single real number. In general, however, a \(k\)-form can “spit out” all kinds of different values. For instance, we might want to deal with quantities that are described by complex numbers (\(\mathbb{C}\)) or vectors in some larger vector space (e.g., \(\mathbb{R}^n\)).

A good example of a vector-valued \(k\)-form is our map \(f: M \rightarrow \mathbb{R}^3\) which represents the geometry of a surface. In the language of exterior calculus, \(f\) is an \(\mathbb{R}^3\)-valued 0-form: at each point \(p\) of \(M\), it takes zero vectors as input and produces a point \(f(p)\) in \(\mathbb{R}^3\) as output. Similarly, the differential \(df\) is an \(\mathbb{R}^3\)-valued 1-form: it takes one vector (some direction \(u\) in the plane) and maps it to a value \(df(u)\) in \(\mathbb{R}^3\) (representing the “stretched out” version of \(u\)).

More generally, if \(E\) is a vector space then an \(E\)-valued \(k\)-form takes \(k\) vectors to a single value in \(E\). However, we have to be a bit careful here. For instance, think about our definition of a 2-form:

\[ \alpha \wedge \beta(u,v) := \alpha(u)\beta(v) - \alpha(v)\beta(u). \]

If \(\alpha\) and \(\beta\) are both \(E\)-valued 1-forms, then \(\alpha(u)\) and \(\beta(v)\) are both vectors in \(E\). But how do you multiply two vectors? In general there may be no good answer: not every vector space comes with a natural notion of multiplication.

However, there are plenty of spaces that do come with a well-defined product — for instance, the product of two complex numbers \(a + bi\) and \(c+di\) is given by \((ac-bd)+(ad+bc)i\), so we have no trouble explicitly evaluating the expression above. In other cases we simply have to say which product we want to use — in \(\mathbb{R}^3\) for instance we could use the cross product \(\times\), in which case an \(\mathbb{R}^3\)-valued 2-form looks like this:

\[ \alpha \wedge \beta(u,v) := \alpha(u) \times \beta(v) - \alpha(v) \times \beta(u). \]

A Quick and Dirty Introduction to Exterior Calculus — Part III: Hodge Duality

Previously we saw that a \(k\)-form measures the (signed) projected volume of a \(k\)-dimensional parallelpiped. For instance, a 2-form measures the area of a parallelogram projected onto some plane, as depicted above. But here’s a nice observation: a plane in \(\mathbb{R}^3\) can be described either by a pair of basis directions \((\alpha,\beta)\), or by a normal direction \(\gamma\). So rather than measuring projected area, we could instead measure how well the normal of a parallelogram \((u,v)\) lines up with the normal of our plane. In other words, we could look for a 1-form \(\gamma\) such that

\[ \gamma(u \times v) = \alpha \wedge \beta(u,v). \]

This observation captures the idea behind Hodge duality: a \(k\)-dimensional volume in an \(n\)-dimensional space can be specified either by \(k\) directions or by a complementary set of \((n-k)\) directions. There should therefore be some kind of natural correspondence between \(k\)-forms and \((n-k)\)-forms.

The Hodge Star

Let’s investigate this idea further by constructing an explicit basis for the space of 0-forms, 1-forms, 2-forms, etc. — to keep things manageable we’ll work with \(\mathbb{R}^3\) and its standard coordinate system \((x^1, x^2, x^3)\). 0-forms are easy: any 0-form can be thought of as some function times the constant 0-from, which we’ll denote “\(1\).” We’ve already seen the 1-form basis \(dx^1, dx^2, dx^3\), which looks like the standard orthonormal basis of a vector space:

What about 2-forms? Well, consider that any 2-form can be expressed as the wedge of two 1-forms:

\[ \alpha \wedge \beta = (\alpha_i dx^i) \wedge (\beta_j dx^j) = \alpha_i \beta_j dx^i \wedge dx^j. \]

In other words, any 2-form looks like some linear combination of the basis 2-forms \(dx^i \wedge dx^j\). How many of these bases are there? Initially it looks like there are a bunch of possibilities:

\[
\begin{array}{ccc}
dx^1 \wedge dx^1 & dx^1 \wedge dx^2 & dx^1 \wedge dx^3 \\
dx^2 \wedge dx^1 & dx^2 \wedge dx^2 & dx^2 \wedge dx^3 \\
dx^3 \wedge dx^1 & dx^3 \wedge dx^2 & dx^3 \wedge dx^3 \\
\end{array}
\]

But of course, not all of these guys are distinct: remember that the wedge product is antisymmetric (\(\alpha \wedge \beta = -\beta \wedge \alpha\)), which has the important consequence \(\alpha \wedge \alpha = 0\). So really our table looks more like this:

\[
\begin{array}{ccc}
0 & dx^1 \wedge dx^2 & -dx^3 \wedge dx^1 \\
-dx^1 \wedge dx^2 & 0 & dx^2 \wedge dx^3 \\
dx^3 \wedge dx^1 & -dx^2 \wedge dx^3 & 0 \\
\end{array}
\]

and we’re left with only three distinct bases: \(dx^2 \wedge dx^3\), \(dx^3 \wedge dx^1\), and \(dx^1 \wedge dx^2\). Geometrically all we’ve said is that there are three linearly-independent “planes” in \(\mathbb{R}^3\):

How about 3-form bases? We certainly have at least one:

\[ dx^1 \wedge dx^2 \wedge dx^3. \]

Are there any others? Again the antisymmetry of \(\wedge\) comes into play: many potential bases are just permutations of this first one:

\[ dx^2 \wedge dx^3 \wedge dx^1 = -dx^2 \wedge dx^1 \wedge dx^3 = dx^1 \wedge dx^2 \wedge dx^3, \]

and the rest vanish due to the appearance of repeated 1-forms:

\[ dx^2 \wedge dx^1 \wedge dx^2 = -dx^2 \wedge dx^2 \wedge dx^1 = 0 \wedge dx^1 = 0. \]

In general there is only one basis \(n\)-form \(dx^1 \wedge \cdots \wedge dx^n\), which measures the usual Euclidean volume of a parallelpiped:

Finally, what about 4-forms on \(\mathbb{R}^3\)? At this point it’s probably pretty easy to see that there are none, since we’d need to pick four distinct 1-form bases from a collection of only three. Geometrically: there are no four-dimensional volumes contained in \(\mathbb{R}^3\)! (Or volumes of any greater dimension, for that matter.) The complete list of \(k\)-form bases on \(\mathbb{R}^3\) is then

  • 0-form bases: 1
  • 1-form bases: \(dx^1\), \(dx^2\), \(dx^3\)
  • 2-form bases: \(dx^2 \wedge dx^3\), \(dx^3 \wedge dx^1\), \(dx^1 \wedge dx^2\)
  • 3-form bases: \(dx^1 \wedge dx^2 \wedge dx^3\),

which means the number of bases is 1, 3, 3, 1. In fact you may see a more general pattern here: the number of \(k\)-form bases on an \(n\)-dimensional space is given by the binomial coefficient

\[ \left( \begin{array}{c} n \\ k \end{array} \right) = \frac{n!}{k!(n-k)!} \]

(i.e., “\(n\) choose \(k\)”), since we want to pick \(k\) distinct 1-form bases and don’t care about the order. An important identity here is

\[ \left( \begin{array}{c} n \\ k \end{array} \right) = \left( \begin{array}{c} n \\ n-k \end{array} \right), \]

which, as anticipated, means that we have a one-to-one relationship between \(k\)-forms and \((n-k)\)-forms. In particular, we can identify any \(k\)-form with its complement. For example, on \(\mathbb{R}^3\) we have

\[
\begin{array}{rcl}
\star\ 1 &=& dx^1 \wedge dx^2 \wedge dx^3 \\
\star\ dx^1 &=& dx^2 \wedge dx^3 \\
\star\ dx^2 &=& dx^3 \wedge dx^1 \\
\star\ dx^3 &=& dx^1 \wedge dx^2 \\
\star\ dx^1 \wedge dx^2 &=& dx^3 \\
\star\ dx^2 \wedge dx^3 &=& dx^1 \\
\star\ dx^3 \wedge dx^1 &=& dx^2 \\
\star\ dx^1 \wedge dx^2 \wedge dx^2 &=& 1
\end{array}
\]

The map \(\star\) (pronounced “star”) is called the Hodge star and captures this idea that planes can be identified with their normals and so forth. More generally, on any flat space we have

\[ \star\ dx^{i_1} \wedge dx^{i_2} \wedge \cdots \wedge dx^{i_k} = dx^{i_{k+1}} \wedge dx^{i_{k+2}} \wedge \cdots \wedge dx^{i_n}, \]

where \((i_1, i_2, \ldots, i_n)\) is any even permutation of \((1, 2, \ldots, n)\).

The Volume Form

So far we’ve been talking about measuring volumes in flat spaces like \(\mathbb{R}^n\). But how do we take measurements in a curved space? Let’s think about our usual example of a surface \(f: \mathbb{R}^2 \supset M \rightarrow \mathbb{R}^3\). If we specify a region on our surface via a pair of unit orthogonal vectors \(u, v \in \mathbb{R}^2\), it’s clear that we don’t want the area \(dx^1 \wedge dx^2(u,v)=1\) since that just gives us the area in the plane. Instead, we want to know what a unit area looks like after it’s been “stretched-out” by the map \(f\). In particular, we said that the length of a vector \(df(u)\) can be expressed in terms of the metric \(g\):

\[ |df(u)| = \sqrt{df(u) \cdot df(u)} = \sqrt{g(u,u)}. \]

So the area we’re really interested in is the product of the lengths \(|df(u)||df(v)| = \sqrt{g(u,u)g(v,v)}\). When \(u\) and \(v\) are orthonormal the quantity \(\det(g) := g(u,u)g(v,v)-2g(u,v)\) is called the determinant of the metric, and can be used to define a 2-form \(\sqrt{\det(g)} dx^1 \wedge dx^2\) that measures any area on our surface. More generally, the \(n\)-form

\[ \omega := \sqrt{\det(g)} dx^1 \wedge \cdots \wedge dx^n \]

is called the volume form, and will play a key role when we talk about integration.

On curved spaces, we’d also like the Hodge star to capture the fact that volumes have been stretched out. For instance, it makes a certain amount of sense to identify the constant function \(1\) with the volume form \(\omega\), since \(\omega\) really represents unit volume on the curved space:

\[ \star 1 = \omega \]

The Inner Product on \(k\)-Forms

More generally we’ll ask that any \(n\)-form constructed from a pair of \(k\)-forms \(\alpha\) and \(\beta\) satisfies

\[ \alpha \wedge \star \beta = \langle\langle \alpha, \beta \rangle\rangle \omega, \]

where \(\langle\langle \alpha, \beta \rangle\rangle = \sum_i \alpha_i \beta_i\) is the inner product on \(k\)-forms. In fact, some authors use this relationship as the definition of the wedge product — in other words, they’ll start with something like, “the wedge product is the unique binary operation on \(k\)-forms such that \(\alpha \wedge \star \beta = \langle\langle \alpha, \beta \rangle\rangle \omega\),” and from there derive all the properties we’ve established above. This treatment is a bit abstract, and makes it far too easy to forget that the wedge product has an extraordinarily concrete geometric meaning. (It’s certainly not the way Hermann Grassmann thought about it when he invented exterior algebra!). In practice, however, this identity is quite useful. For instance, if \(u\) and \(v\) are vectors in \(\mathbb{R}^3\), then we can write

\[ u \cdot v = \star\left(u^\flat \wedge \star v^\flat\right), \]

i.e., on a flat space we can express the usual Euclidean inner product via the wedge product. Is it clear geometrically that this identity is true? Think about what it says: the Hodge star turns \(v\) into a plane with \(v\) as a normal. We then build a volume by extruding this plane along the direction \(u\). If \(u\) and \(v\) are nearly parallel the volume will be fairly large; if they’re nearly orthogonal the volume will be quite shallow. (But to be sure we really got it right, you should try verifying this identity in coordinates!) Similarly, we can express the Euclidean cross product as just

\[ u \times v = \star(u^\flat \wedge v^\flat)^\sharp, \]

i.e., we can create a plane with normal \(u \times v\) by wedging together the two basis vectors \(u\) and \(v\). (Again, working this fact out in coordinates may help soothe your paranoia.)

A Quick and Dirty Introduction to Exterior Calculus — Part IV: Differential Operators

Originally we set out to develop exterior calculus. The objects we’ve looked at so far — \(k\)-forms, the wedge product \(\wedge\) and the Hodge star \(\star\) — actually describe a more general structure called an exterior algebra. To turn our algebra into a calculus, we also need to know how quantities change, as well as how to measure quantities. In other words, we need some tools for differentiation and integration. Let’s start with differentiation.

In our discussion of surfaces we briefly looked at the differential \(df\) of a surface \(f: M \rightarrow \mathbb{R}^3\), which tells us something about the way tangent vectors get “stretched out” as we move from the domain \(M\) to a curved surface sitting in \(\mathbb{R}^3\). More generally \(d\) is called the exterior derivative and is responsible for building up many of the differential operators in exterior calculus. The basic idea is that \(d\) tells us how quickly a \(k\)-form changes along every possible direction. But how exactly is it defined? So far we’ve seen only a high-level geometric description.

Div, Grad, and Curl

Before jumping into the exterior derivative, it’s worth reviewing what the basic vector derivatives \(\mathrm{div}\), \(\mathrm{grad}\), and \(\mathrm{curl}\) do, and more importantly, what they look like. The key player here is the operator \(\nabla\) (pronounced “nabla”) which can be expressed in coordinates as the vector of all partial derivatives:

\[ \nabla := \left( \frac{\partial}{\partial x^1}, \ldots, \frac{\partial}{\partial x^n} \right). \]

For instance, applying \(\nabla\) to a scalar function \(\phi: \mathbb{R}^n \rightarrow \mathbb{R}\) yields the gradient

\[ \nabla\phi = \left( \frac{\partial f}{\partial x^1}, \ldots, \frac{\partial f}{\partial x^n} \right), \]

which can be visualized as the direction of steepest ascent on some terrain:

We can also apply \(\nabla\) to a vector field \(X\) in two different ways. The dot product gives us the divergence

\[ \nabla \cdot X = \frac{\partial X^1}{\partial x^1} + \cdots + \frac{\partial X^n}{\partial x^n} \]

which measures how quickly the vector field is “spreading out”, and on \(\mathbb{R}^3\) the cross product gives us the curl

\[ \nabla \times X = \left( \frac{\partial X_3}{\partial x^2} - \frac{\partial X_2}{\partial x^3}, \frac{\partial X_1}{\partial x^3} - \frac{\partial X_3}{\partial x^1}, \frac{\partial X_2}{\partial x^1} - \frac{\partial X_1}{\partial x^2} \right), \]

which indicates how much a vector field is “spinning around.” For instance, here’s a pair of vector fields with a lot of divergence and a lot of curl, respectively:

(Note that in this case one field is just a 90-degree rotation of the other!) On a typical day it’s a lot more useful to think of \(\mathrm{div}\), \(\mathrm{grad}\) and \(\mathrm{curl}\) in terms of these kinds of pictures rather than the ugly expressions above.

Think Differential

Not surprisingly, we can express similar notions using exterior calculus. However, these notions will be a bit easier to generalize (for instance, what does “curl” mean for a vector field in \(\mathbb{R}^4\), where no cross product is defined?). Let’s first take a look at the exterior derivative of 0-forms (i.e., functions), which is often just called the differential. To keep things simple, we’ll start with real-valued functions \(\phi: \mathbb{R}^n \rightarrow \mathbb{R}\). In coordinates, the differential is defined as

\[ d\phi := \frac{\partial \phi}{\partial x^1} dx^1 + \cdots + \frac{\partial \phi}{\partial x^n} dx^n. \]

It’s important to note that the terms \(\frac{\partial \phi}{\partial x^i}\) actually correspond to partial derivatives of our function \(\phi\), whereas the terms \(dx^i\) simply denote an orthonormal basis for \(\mathbb{R}^n\). In other words, you can think of \(d\phi\) as just a list of all the partial derivatives of \(\phi\). Of course, this object looks a lot like the gradient \(\nabla \phi\) we saw just a moment ago. And indeed the two are closely related, except for the fact that \(\nabla \phi\) is a vector field and \(d\phi\) is a 1-form. More precisely,

\[ \nabla \phi = (d\phi)^\sharp. \]

Directional Derivatives

Another way to investigate the behavior of the exterior derivative is to see what happens when we stick a vector \(u\) into the 1-form \(df\). In coordinates we get something that looks like a dot product between the gradient of \(f\) and the vector \(u\):

\[ df(u) = \frac{\partial f}{\partial x^1} u^1 + \cdots + \frac{\partial f}{\partial x^n} u^n. \]

For instance, in \(\mathbb{R}^2\) we could stick in the unit vector \(u = (1,0)\) to get the partial derivative \(\frac{\partial f}{\partial x^1}\) along the first coordinate axis:

(Compare this picture to the picture of the gradient we saw above.) In general, \(df(u)\) represents the directional derivative of \(f\) along the direction \(u\). In other words, it tells us how quickly \(f\) changes if we take a short walk in the direction \(u\). Returning again to vector calculus notation, we have

\[ df(u) = u \cdot \nabla f. \]

Properties of the Exterior Derivative

How do derivatives of arbitrary \(k\)-forms behave? For one thing, we expect \(d\) to be linear — after all, a derivative is just the limit of a difference, and differences are certainly linear! What about the derivative of a wedge of two forms? Harkening back to good old-fashioned calculus, here’s a picture that explains the typical product rule \(\frac{\partial}{\partial x}(f(x)g(x)) = f’(x)g(x) + f(x)g’(x) \):

The dark region represents the value of \(fg\) at \(x\); the light blue region represents the change in this value as we move \(x\) some small distance \(h\). As \(h\) gets smaller and smaller, the contribution of the upper-right quadrant becomes negligible and we can write the derivative as the change in \(f\) times \(g\) plus the change in \(g\) times \(f\). (Can you make this argument more rigorous?) Since a \(k\)-form also measures a (signed) volume, this intuition also carries over to the exterior derivative of a wedge product. In particular, if \(\alpha\) is a \(k\)-form then \(d\) obeys the rule

\[ d(\alpha \wedge \beta) = d\alpha \wedge \beta + (-1)^{k}\alpha \wedge d\beta. \]

which says that the rate of change of the overall volume can be expressed in terms of changes in the constituent volumes, exactly as in the picture above.

Exterior Derivative of 1-Forms

To be a little more concrete, let’s see what happens when we differentiate a 1-form on \(\mathbb{R}^3\). Working things out in coordinates turns out to be a total mess, but in the end you may be pleasantly surprised with the simplicity of the outcome! (Later on we’ll see that these ideas can also be expressed quite nicely without coordinates using Stokes’ theorem, which paves the way to differentiation in the discrete setting.) Applying the linearity of \(d\), we have

\[
\begin{array}{rcl}
d\alpha
&=& d(\alpha_1 dx^1 + \alpha_2 dx^2 + \alpha_3 dx^3) \\
&=& d(\alpha_1 dx^1) + d(\alpha_2 dx^2) + d(\alpha_3 dx^3).
\end{array}
\]

Each term \(\alpha_j dx^j\) can really be thought of a wedge product \(\alpha_j \wedge dx^j\) between a 0-form \(\alpha_j\) and the corresponding basis 1-form \(dx^j\). Applying the exterior derivative to one of these terms we get

\[ d(\alpha_j \wedge dx^j) = (d\alpha_j) \wedge dx^j + \alpha_j \wedge \underbrace{(ddx^j)}_{=0} = \frac{\partial \alpha_j}{\partial x^i} dx^i \wedge dx^j. \]

To keep things short we used the Einstein summation convention here, but let’s really write out all the terms:

\[
\begin{array}{rcccccc}
d\alpha &=& \frac{\partial \alpha_1}{\partial x^1} dx^1 \wedge dx^1 &+& \frac{\partial \alpha_1}{\partial x^2} dx^2 \wedge dx^1 &+& \frac{\partial \alpha_1}{\partial x^3} dx^3 \wedge dx^1 \\
&& \frac{\partial \alpha_2}{\partial x^1} dx^1 \wedge dx^2 &+& \frac{\partial \alpha_2}{\partial x^2} dx^2 \wedge dx^2 &+& \frac{\partial \alpha_2}{\partial x^3} dx^3 \wedge dx^2 \\
&& \frac{\partial \alpha_3}{\partial x^1} dx^1 \wedge dx^3 &+& \frac{\partial \alpha_3}{\partial x^2} dx^2 \wedge dx^3 &+& \frac{\partial \alpha_3}{\partial x^3} dx^3 \wedge dx^3. \\
\end{array}
\]

Using the fact that \(\alpha \wedge \beta = -\beta \wedge \alpha\), we get a much simpler expression

\[
\begin{array}{rcl}
d\alpha &=& ( \frac{\partial \alpha_3}{\partial x^2} - \frac{\partial \alpha_2}{\partial x^3} ) dx^2 \wedge dx^3 \\
&& ( \frac{\partial \alpha_1}{\partial x^3} - \frac{\partial \alpha_3}{\partial x^1} ) dx^3 \wedge dx^1 \\
&& ( \frac{\partial \alpha_2}{\partial x^1} - \frac{\partial \alpha_1}{\partial x^2} ) dx^1 \wedge dx^2. \\
\end{array}
\]

Does this expression look familiar? If you look again at our review of vector derivatives, you’ll recognize that \(d\alpha\) basically looks like the curl of \(\alpha^\sharp\), except that it’s expressed as a 2-form instead of a vector field. Also remember (from our discussion of Hodge duality) that a 2-form and a 1-form are not so different here — geometrically they both specify some direction in \(\mathbb{R}^3\). Therefore, we can express the curl of any vector field \(X\) as

\[ \nabla \times X = \left( \star d X^\flat \right)^\sharp. \]

It’s worth stepping through the sequence of operations here to check that everything makes sense: \(\flat\) converts the vector field \(X\) into a 1-form \(X^\flat\); \(d\) computes something that looks like the curl, but expressed as a 2-form \(dX^\flat\); \(\star\) turns this 2-form into a 1-form \(\star d X^\flat\); and finally \(\sharp\) converts this 1-form back into the vector field \(\left( \star d X^\flat \right)^\sharp\). The take-home message here, though, is that the exterior derivative of a 1-form looks like the curl of a vector field.

So far we know how to express the gradient and the curl using \(d\). What about our other favorite vector derivative, the divergence? Instead of grinding through another tedious derivation, let’s make a simple geometric observation: in \(\mathbb{R}^2\) at least, we can determine the divergence of a vector field by rotating it by 90 degrees and computing its curl (consider the example we saw earlier). Moreover, in \(\mathbb{R}^2\) the Hodge star \(\star\) represents a rotation by 90 degrees, since it identifies any line with the direction orthogonal to that line:

Therefore, we might suspect that divergence can be computed by first applying the Hodge star, then applying the exterior derivative:

\[ \nabla \cdot X = \star d \star X^\flat. \]

The leftmost Hodge star accounts for the fact that \(d \star X^\flat\) is an \(n\)-form instead of a 0-form — in vector calculus divergence is viewed as a scalar quantity. Does this definition really work? Let’s give it a try in coordinates on \(\mathbb{R}^3\). First, we have

\[
\begin{array}{rcl} \star X^\flat &=& \star( X_1 dx^1 + X_2 dx^2 + X_3 dx^3 ) \\
&=& X_1 dx^2 \wedge dx^3 + X_2 dx^3 \wedge dx^1 + X_3 dx^1 \wedge dx^2.
\end{array}
\]

Differentiating we get

\[
\begin{array}{rcl} d \star X^\flat &=& \frac{\partial X_1}{\partial x^1} dx^1 \wedge dx^2 \wedge dx^3 + \\
&& \frac{\partial X_2}{\partial x^2} dx^2 \wedge dx^3 \wedge dx^1 + \\
&& \frac{\partial X_3}{\partial x^3} dx^3 \wedge dx^1 \wedge dx^2, \\
\end{array}
\]

but of course we can rearrange these wedge products to simply

\[ d \star X^\flat = \left( \frac{\partial X_1}{\partial x^1} + \frac{\partial X_2}{\partial x^2} + \frac{\partial X_3}{\partial x^3} \right) dx^1 \wedge dx^2 \wedge dx^3. \]

A final application of the Hodge star gives us the divergence

\[ \star d \star X^\flat = \frac{\partial X_1}{\partial x^1} + \frac{\partial X_2}{\partial x^2} + \frac{\partial X_3}{\partial x^3} \]

as desired.

In summary, for any scalar field \(\phi\) and vector field \(X\) we have

\[
\begin{array}{rcl}
\nabla \phi &=& (d\phi)^\sharp \\
\nabla \times X &=& \left( \star d X^\flat \right)^\sharp \\
\nabla \cdot X &=& \star d \star X^\flat
\end{array}
\]

One cute thing to notice here is that (in \(\mathbb{R}^3\)) \(\mathrm{grad}\), \(\mathrm{curl}\), and \(\mathrm{div}\) are more or less just \(d\) applied to a \(0-\), \(1-\) and \(2-\) form, respectively.

The Laplacian

Another key differential operator from vector calculus is the scalar Laplacian which (confusingly!) is often denoted by \(\Delta\) or \(\nabla^2\), and is defined as

\[ \Delta := \nabla \cdot \nabla, \]

i.e., the divergence of the gradient. Although the Laplacian may seem like just yet another in a long list of derivatives, it deserves your utmost respect: the Laplacian is central to the most fundamental of physical laws (any diffusion process and all forms of wave propagation, including the Schrödinger equation); its eigenvalues capture almost everything there is to know about a given piece of geometry (can you hear the shape of a drum?). Heavy tomes and entire lives have been devoted to the Laplacian, and in the discrete setting we’ll see that this one simple operator can be applied to a diverse array of tasks (surface parameterization, surface smoothing, vector field design and decomposition, distance computation, fluid simulation… you name it, we got it!).

Fortunately, now that we know how to write \(\mathrm{div}\), \(\mathrm{grad}\) and \(\mathrm{curl}\) using exterior calculus, expressing the scalar Laplacian is straightforward: \(\Delta = \star d \star d\). More generally, the \(k\)-form Laplacian is given by

\[ \Delta := \star d \star d + d \star d \star. \]

The name “Laplace-Beltrami” is used merely to indicate that the domain may have some amount of curvature (encapsulated by the Hodge star). Some people like to define the operator \(\delta := \star d \star\), called the codifferential, and write the Laplacian as \(\Delta = \delta d + d \delta\).

One question you might ask is: why is the Laplacian for 0-forms different from the general \(k\)-form Laplacian? Actually, it’s not — consider what happens when we apply the term \(d \star d \star\) to a 0-form \(\phi\): \(\star \phi\) is an \(n\)-form, and so \(d \star \phi\) must be an \((n+1)\)-form. But there are no \((n+1)\)-forms on an \(n\)-dimensional space! So this term is often omitted when writing the scalar Laplacian.

November 1, 2011 | Comments Closed

A Quick and Dirty Introduction to Exterior Calculus — Part V: Integration and Stokes’ Theorem

In the last set of notes we talked about how to differentiate \(k\)-forms using the exterior derivative \(d\). We’d also like some way to integrate forms. Actually, there’s surprisingly little to say about integration given the setup we already have. Suppose we want to compute the total area \(A_\Omega\) of a region \(\Omega\) in the plane:

If you remember back to calculus class, the basic idea was to break up the domain into a bunch of little pieces that are easy to measure (like squares) and add up their areas:

\[ A_\Omega \approx \sum_i A_i. \]

As these squares get smaller and smaller we get a better and better approximation, ultimately achieving the true area

\[ A_\Omega = \int_\Omega dA. \]

Alternatively, we could write the individual areas using differential forms — in particular, \(A_i = dx^1 \wedge dx^2(u,v)\). Therefore, the area element \(dA\) is really nothing more than the standard volume form \(dx^1 \wedge dx^2\) on \(\mathbb{R}^2\). (Not too surprising, since the whole point of \(k\)-forms was to measure volume!)

To make things more interesting, let’s say that the contribution of each little square is weighted by some scalar function \(\phi\). In this case we get the quantity

\[ \int_\Omega \phi\ dA = \int_\Omega \phi\ dx^1 \wedge dx^2. \]

Again the integrand \(\phi\ dx^1 \wedge dx^2\) can be thought of as a 2-form. In other words, you’ve been working with differential forms your whole life, even if you didn’t realize it! More generally, integrands on an \(n\)-dimensional space are always \(n\)-forms, since we need to “plug in” \(n\) orthogonal vectors representing the local volume. For now, however, looking at surfaces (i.e., 2-manifolds) will give us all the intuition we need.

Integration on Surfaces

If you think back to our discussion of the Hodge star, you’ll remember the volume form

\[ \omega = \sqrt{\mathrm{det}(g)} dx^1 \wedge dx^2, \]

which measures the area of little parallelograms on our surface. The factor \(\sqrt{\mathrm{det}(g)}\) reminds us that we can’t simply measure the volume in the domain \(M\) — we also have to take into account any “stretching” induced by the map \(f: M \rightarrow \mathbb{R}^2\). Of course, when we integrate a function on a surface, we should also take this stretching into account. For instance, to integrate a function \(\phi: M \rightarrow \mathbb{R}\), we would write

\[ \int_\Omega \phi \omega = \int_\Omega \phi \sqrt{\mathrm{det}(g)}\ dx^1 \wedge dx^2. \]

In the case of a conformal parameterization things become even simpler — since \(\sqrt{\mathrm{det}(g)} = a\) we have just

\[ \int_\Omega \phi a\ dx^1 \wedge dx^2, \]

where \(a: M \rightarrow \mathbb{R}\) is the scaling factor. In other words, we scale the value of \(\phi\) up or down depending on the amount by which the surface locally “inflates” or “deflates.” In fact, this whole story gives a nice geometric interpretation to good old-fashioned integrals: you can imagine that \(\int_\Omega \phi\ dA\) represents the area of some suitably deformed version of the initially planar region \(\Omega\).

Stokes’ Theorem

The main reason for studying integration on manifolds is to take advantage of the world’s most powerful tool: Stokes’ theorem. Without further ado, Stokes’ theorem says that

\[ \int_\Omega d\alpha = \int_{\partial\Omega} \alpha, \]

where \(\alpha\) is any \(n-1\)-form on an \(n\)-dimensional domain \(\Omega\). In other words, integrating a differential form over the boundary of a manifold is the same as integrating its derivative over the entire domain.

If this trick sounds familiar to you, it’s probably because you’ve seen it time and again in different contexts and under different names: the divergence theorem, Green’s theorem, the fundamental theorem of calculus, Cauchy’s integral formula, etc. Picking apart these special cases will really help us understand the more general meaning of Stokes’ theorem.

Divergence Theorem

Let’s start with the divergence theorem from vector calculus, which says that

\[ \int_\Omega \nabla \cdot X dA = \int_{\partial\Omega} N \cdot X d\ell, \]

where \(X\) is a vector field on \(\Omega\) and \(N\) represents the (outward-pointing) unit normals on the boundary of \(\Omega\). A better name for this theorem might have been the “what goes in must come out theorem”, because if you think about \(X\) as the flow of water throughout the domain \(\Omega\) then it’s clear that the amount of water being pumped into \(\Omega\) (via pipes in the ground) must be the same as the amount flowing out of its boundary at any moment in time:

Let’s try writing this theorem using exterior calculus. First, remember that we can write the divergence of \(X\) as \(\nabla \cdot X = \star d \star X^\flat\). It’s a bit harder to see how to write the right-hand side of the divergence theorem, but think about what integration does here: it takes tangents to the boundary and sticks them into a 1-form. For instance, \(\int_\Omega X^\flat\) adds up the tangential components of \(X\). To get the normal component we could rotate \(X^\flat\) by a quarter turn, which conveniently enough is achieved by hitting it with the Hodge star. Overall we get

\[ \int_\Omega d \star X^\flat = \int_{\partial\Omega} \star X^\flat, \]

which, as promised, is just a special case of Stokes’ theorem. Alternatively, we can use Stokes’ theorem to provide a more geometric interpretation of the divergence operator itself: when integrated over any region \(\Omega\) — no matter how small — the divergence operator gives the total flux through the region boundary. In the discrete case we’ll see that this boundary flux interpretation is the only notion of divergence — in other words, there’s no concept of divergence at a single point.

By the way, why does \(d \star X^\flat\) appear on the left-hand side instead of \(\star d \star X^\flat\)? The reason is that \(\star d \star X^\flat\) is a 0-form, so we have to hit it with another Hodge star to turn it into an object that measures areas (i.e., a 2-form). Applying this transformation is no different from appending \(dA\) to \(\nabla \cdot X\) — we’re specifying how volume should be measured on our domain.

Fundamental Theorem of Calculus

The fundamental theorem of calculus is in fact so fundamental that you may not even remember what it is. It basically says that for a real-valued function \(\phi: \mathbb{R} \rightarrow \mathbb{R}\) on the real line

\[ \int_a^b \frac{\partial \phi}{\partial x} dx = \phi(b) - \phi(a). \]

In other words, the total change over an interval \([a,b]\) is (as you might expect) how much you end up with minus how much you started with. But soft, behold! All we’ve done is written Stokes’ theorem once again:

\[ \int_{[a,b]} d\phi = \int_{\partial[a,b]} \phi, \]

since the boundary of the interval \([a,b]\) consists only of the two endpoints \(a\) and \(b\).

Hopefully these two examples give you a good feel for what Stokes’ theorem says. In the end, it reads almost like a Zen kōan: what happens on the outside is purely a function of the change within. (Perhaps it is Stokes’ that deserves the name, “fundamental theorem of calculus!”)

November 5, 2011 | Comments Closed

Discrete Exterior Calculus

So far we’ve been exploring exterior calculus purely in the smooth setting. Unfortunately this theory was developed by some old-timers who did not know anything about computers, hence it cannot be used directly by machines that store only a finite amount of information. For instance, if we have a smooth vector field or a smooth 1-form we can’t possibly store the direction of every little “arrow” at each point — there are far too many of them! Instead, we need to keep track of a discrete (or really, finite) number of pieces of information that capture the essential behavior of the objects we’re working with; we call this scheme discrete exterior calculus (or DEC for short). The big secret about DEC is that it’s literally nothing more than the good-old fashioned (continuous) exterior calculus we’ve been learning about, except that we integrate differential forms over elements of our mesh.

Discrete Differential Forms

One way to encode a 1-form might be to store a finite collection of “arrows” associated with some subset of points. Instead, we’re going to do something a bit different: we’re going to integrate our 1-form over each edge of a mesh, and store the resulting numbers (remember that the integral of an \(n\)-form always spits out a single number) on the corresponding edges. In other words, if \(\alpha\) is a 1-form and \(e\) is an edge, then we’ll associate the number

\[ \hat{\alpha}_e := \int_e \alpha \]

with \(e\), where the use of the hat (\(\ \hat{}\ \)) is supposed to suggest a discrete quantity (not to be confused with a unit-length vector).

Does this procedure seem a bit abstract to you? It shouldn’t! Think about what this integral represents: it tells us how strongly the 1-form \(\alpha\) “flows along” the edge \(e\) on average. More specifically, remember how integration of a 1-form works: at each point along the edge we take the vector tangent to the edge, stick it into the 1-form \(\alpha\), and sum up the resulting values — each value tells us something about how well \(\alpha\) “lines up” with the direction of the edge. For instance, we could approximate the integral via the sum

\[ \int_e \alpha \approx |e|\left(\frac{1}{N} \sum_{i=1}^N \alpha_{p_i}(\hat{e})\right), \]

where \(|e|\) denotes the length of the edge, \(\{p_i\}\) is a sequence of points along the edge, and \(\hat{e} := e/|e|\) is a unit vector tangent to the edge:

Of course, this quantity tells us absolutely nothing about the strength of the “flow” orthogonal to the edge: it could be zero, it could be enormous! We don’t really know, because we didn’t take any measurements along the orthogonal direction. However, the hope is that some of this information will still be captured by nearby edges (which are most likely not parallel to \(e\)).

More generally, a \(k\)-form that has been integrated over each \(k\)-dimensional cell (edges in 1D, faces in 2D, etc.) is called a discrete differential \(k\)-form. (If you ever find the distinction confusing, you might find it helpful to substitute the word “integrated” for the word “discrete.”) In practice, however, not every discrete differential form has to originate from a continuous one — for instance, a bunch of arbitrary values assigned to each edge of a mesh is a perfectly good discrete 1-form.

Orientation

One thing you may have noticed in all of our illustrations so far is that each edge is marked with a little arrow. Why? Well, one thing to remember is that direction matters when you integrate. For instance, the fundamental theorem of calculus (and common sense) tells us that the total change as you go from \(a\) to \(b\) is the opposite of the total change as you go from \(b\) to \(a\):

\[ \int_a^b \frac{\partial\phi}{\partial x} dx = \phi(b)-\phi(a) = -(\phi(a)-\phi(b)) = -\int_b^a \frac{\partial\phi}{\partial x} dx. \]

Said in a much less fancy way: the elevation gain as you go from Pasadena to Altadena is 151 meters, so of the elevation “gain” in the other direction must be -151 meters! Just keeping track of the number 151 does you little good — you have to say what that quantity represents.

Therefore, when we store a discrete differential form it’s not enough to just store a number: we also have to specify a canonical orientation for each element of our mesh, corresponding to the orientation we used during integration. For an edge we’ve already seen that we can think about orientation as a little arrow pointing from one vertex to another — we could also just think of an edge as an ordered pair \((i,j)\), meaning that we always integrate from \(i\) to \(j\).

More generally, suppose that each element of our mesh is an oriented \(k\)-simplex \(\sigma\), i.e., a collection of \(k+1\) vertices \(p_i \in \mathbb{R}^n\) given in some fixed order \((p_1, \ldots, p_{k+1})\). The geometry associated with \(\sigma\) is the convex combination of these points:

\[ \left\{ \sum_{i=1}^{k+1} \lambda_i p_i \left| \sum_{i=1}^{k+1} \lambda_i = 1 \right. \right\} \subset \mathbb{R}^n \]

(Convince yourself that a 0-simplex is a vertex, a 1-simplex is an edge, a 2-simplex is a triangle, and a 3-simplex is a tetrahedron.)

Two oriented \(k\)-simplices have the same orientation if and only if the vertices of one are an even permutation of the vertices of another. For instance, the triangles \((p_1, p_2, p_3)\) and \((p_2, p_3, p_1)\) have the same orientation; \((p_1, p_2, p_3)\) and \((p_2, p_1, p_3)\) have opposite orientation.

If a simplex \(\sigma_1\) is a (not necessarily proper) subset of another simplex \(\sigma_2\), then we say that \(\sigma_1\) is a face of \(\sigma_2\). For instance, every vertex, edge, and triangle of a tetrahedron \(\sigma\) is a face of \(\sigma\); as is \(\sigma\) itself! Moreover, the orientation of a simplex agrees with the orientation of one of its faces as long as we see an even permutation on the shared vertices. For instance, the orientations of the edge \((p_2,p_1)\) and the triangle \((p_1,p_3,p_2)\) agree. Geometrically all we’re saying is that the two “point” in the same direction (as depicted above). To keep yourself sane while working with meshes, the most important thing is to pick and orientation and stick with it!

So in general, how do we integrate a \(k\)-form over an oriented \(k\)-simplex? Remember that a \(k\)-form is going to “eat” \(k\) vectors at each point and spit out a number — a good canonical choice is to take the ordered collection of edge vectors \((p_2 – p_1, \ldots, p_{k+1}-p_1)\) and orthogonalize them (using, say the Gram-Schmidt algorithm) to get vectors \((u_1, \ldots, u_n)\). This way the sign of the integrand changes whenever the orientation changes. Numerically, we can then approximate the integral via a sum

\[ \int_\sigma \alpha \approx \frac{|\sigma|}{N} \sum_{i=1}^N \alpha_{p_i}(u_1, \ldots, u_n) \]

where \(\{p_i\}\) is a (usually carefully-chosen) collection of sample points. (Can you see why the orientation of \(\sigma\) affects the sign of the integrand?) Sounds like a lot of work, but in practice one rarely constructs discrete differential forms via integration: more often, discrete forms are constructed via input data that is already discrete (e.g., vertex positions in a triangle mesh).

By the way, what’s a discrete 0-form? Give up? Well, it must be a 0-form (i.e., a function) that’s been integrated over every 0-simplex (i.e., vertex) of a mesh:

\[ \hat{\phi}_i = \int_{v_i} \phi \]

By convention, the integral of a function over a zero-dimensional set is simply the value of the function at that point: \(\hat{\phi}_i = \phi(v_i)\). In other words, in the case of 0-forms there is no difference between storing point samples and storing integrated quantities: the two coincide.

It’s also important to remember that differential forms don’t have to be real-valued. For instance, we can think of a map \(f: M \rightarrow \mathbb{R}^3\) that encodes the geometry of a surface as an \(\mathbb{R}^3\)-valued 0-form; its differential \(df\) is then an \(\mathbb{R}^3\)-valued 1-form, etc. Likewise, when we say that a discrete differential form is a number stored on every mesh element, the word “number” is used in a fairly loose sense: a number could be a real value, a vector, a complex number, a quaternion, etc. For instance, the collection of \((x,y,z)\) vertex coordinates of a mesh can be viewed as an \(\mathbb{R}^3\)-valued discrete 0-form (namely, one that discretizes the map \(f\)). The only requirement, of course, is that we store the same type of number on each mesh element.

The Discrete Exterior Derivative

One of the main advantages of working with integrated (i.e., “discrete”) differential forms instead of point samples is that we can easily take advantage of Stokes’ theorem. Remember that Stokes’ theorem says

\[ \int_\Omega d\alpha = \int_{\partial\Omega} \alpha, \]

for any \(k\)-form \(\alpha\) and \(k+1\)-dimensional domain \(\Omega\). In other words, we can integrate the derivative of a differential form as long as we know its integral along the boundary. But that’s exactly the kind of information encoded by a discrete differential form! For instance, if \(\hat{\alpha}\) is a discrete 1-form stored on the three edges of a triangle \(\sigma\), then we have

\[ \int_\sigma d\alpha = \int_{\partial\sigma} \alpha = \sum_{i=1}^3 \int_{e_i} \alpha = \sum_{i=1}^3 \hat{\alpha}_i. \]

In other words, we can exactly evaluate the integral on the left by just adding up three numbers. Pretty cool! In fact, the thing on the left is also a discrete differential form: it’s the 2-form \(d\alpha\) integrated over the only triangle in our mesh. So for convenience, we’ll call this guy “\(\hat{d}\hat{\alpha}\)”, and we’ll call the operation \(\hat{d}\) the discrete exterior derivative. (In the future we will drop the hats from our notation when the meaning is clear from context.) In other words, the discrete exterior derivative takes a \(k\)-form that has already been integrated over each \(k\)-simplex and applies Stokes’ theorem to get the integral of the derivative over each \(k+1\)-simplex.

In practice (i.e., in code) you can see how this operation might be implemented by simply taking local sums over the appropriate mesh elements. However, in the example above we made life particularly easy on ourselves by giving each edge an orientation that agrees with the orientation of the triangle. Unfortunately assigning a consistent orientation to every simplex is not always possible, and in general we need to be more careful about sign when adding up our piecewise integrals. For instance, in the example below we’d have

\[ (\hat{d}\hat{\alpha})_1 = \hat{\alpha}_1 + \hat{\alpha}_2 + \hat{\alpha}_3 \]

and

\[ (\hat{d}\hat{\alpha})_2 = \hat{\alpha}_4 + \hat{\alpha}_5 - \hat{\alpha}_2. \]

Discrete Hodge Star

As hinted at above, a discrete \(k\)-form captures the behavior of a continuous \(k\)-form along \(k\) directions, but not along the remaining \(n-k\) directions — for instance, a discrete 1-form in 2D captures the flow along edges but not in the orthogonal direction. If you paid attention to our discussion of Hodge duality, this story starts to sound familiar! To capture Hodge duality in the discrete setting, we’ll need to define a dual mesh. In general, the dual of an \(n\)-dimensional simplicial mesh identifies every \(k\)-simplex in the primal (i.e., original) mesh with a unique \((n-k)\)-cell in the dual mesh. In a two-dimensional simplicial mesh, for instance, primal vertices are identified with dual faces, primal edges are identified with dual edges, and primal faces are identified with dual vertices. Note, however, that the dual cells are not always simplices! (See above.)

So how do we talk about Hodge duality in discrete exterior calculus? Quite simply, the discrete Hodge dual of a (discrete) \(k\)-form on the primal mesh is an \((n-k)\)-form on the dual mesh. Similarly, the Hodge dual of an \(k\)-form on the dual mesh is a \(k\)-form on the primal mesh. Discrete forms on the primal mesh are called primal forms and discrete forms on the dual mesh are called dual forms. Given a discrete form \(\hat{\alpha}\) (whether primal or dual), its Hodge dual is typically written as \(\hat{\star} \hat{\alpha}\).

Unlike continuous forms, discrete primal and dual forms live in different places (so for instance, discrete primal \(k\)-forms and dual \(k\)-forms cannot be added to each other). In fact, primal and dual forms often have different physical interpretations. For instance, a primal 1-form might represent the total circulation along edges of the primal mesh, whereas in the same context a dual 1-form might represent the total flux through the corresponding dual edges (see illustration above).

Of course, these two quantities (flux and circulation) are closely related, and naturally leads into one definition for a discrete Hodge star called the diagonal Hodge star. Consider a primal \(k\)-form \(\alpha\). If \(\hat{\alpha}_i\) is the value of \(\hat{\alpha}\) on the \(k\)-simplex \(\sigma_i\), then the diagonal Hodge star is defined by

\[ \hat{\star} \hat{\alpha}_i = \frac{|\sigma_i^\star|}{|\sigma_i|} \hat{\alpha}_i \]

for all \(i\), where \(|\sigma|\) indicates the (unsigned) volume of \(\sigma\) (which by convention equals one for a vertex!) and \(|\sigma^\star|\) is the volume of the corresponding dual cell. In other words, to compute the dual form we simply multiply the scalar value stored on each cell by the ratio of corresponding dual and primal volumes.

If we remember that a discrete form can be thought of as a continuous form integrated over each cell, this definition for the Hodge star makes perfect sense: the primal and dual quantities should have the same density, but we need to account for the fact that they are integrated over cells of different volume. We therefore normalize by a ratio of volumes when mapping between primal and dual. This particular Hodge star is called diagonal since the \(i\)th element of the dual differential form depends only on the \(i\)th element of the primal differential form. It’s not hard to see, then, that Hodge star taking dual forms to primal forms (the dual Hodge star) is the inverse of the one that takes primal to dual (the primal Hodge star).

That’s All, Folks!

Hey, wait a minute, what about our other operations, like the wedge product (\(\wedge\))? These operations can certainly be defined in the discrete setting, but we won’t go into detail here — the basic recipe is to integrate, integrate, integrate. Actually, even in continuous exterior calculus we omitted a couple operations like the Lie derivative (\(\mathcal{L}_X\)) and the interior product (\(i_\alpha\)). Coming up with a complete discrete calculus where the whole cast of characters \(d\), \(\wedge\), \(\star\), \(\mathcal{L}_X\), \(i_\alpha\), etc., plays well together is an active and ongoing area of research, which may be of interest to aspiring young researchers like you (yes, you)!

November 9, 2011 | Comments Closed