The Material Derivative: How a Fluid Parcel Accelerates

Post 9 closed the first series with the Fourier Neural Operator, and the running examples all the way from Post 2 were the heat equation and the wave equation. A large share of the research on neural PDE solvers, however, is demonstrated on a different family: the Navier–Stokes equations, which describe how fluids such as air and water move. This post starts a new series that derives them from scratch, the same way Post 3 derived the wave equation from \(F = ma\).

This first chapter handles only the \(a\): what it means for a fluid to accelerate, and how to compute that acceleration. This leads to two ways of describing a flow, called the Lagrangian and the Eulerian viewpoints, and ends at a quantity called the material derivative. None of these terms needs to mean anything yet; each one is built from a physical picture first, and the symbols follow.

In Post 3, the move that made the wave equation derivable was to zoom into one tiny segment of string, \([x, x+\Delta x]\), and apply \(F = ma\) to it. The same move works for a fluid, but it first needs an answer to what "a tiny piece" of fluid is.

A fluid is made of molecules, and there are far too many to track: one cubic meter of air holds about \(2.5\times10^{25}\) of them. Instead, we work with a fluid parcel: a blob of fluid that is much smaller than any feature of the flow we care about, yet still contains an enormous number of molecules. A parcel of 1 mm³ is tiny next to a 4 m car in a wind tunnel, but it still holds about \(2.5\times10^{16}\) molecules. Averaged over that many molecules, a quantity like the parcel's velocity is perfectly well defined, even though every individual molecule zips around randomly. Treating the fluid as a smooth substance with a definite velocity at every point is called the continuum assumption.

The string also needed a field to describe its state: the displacement \(u(x,t)\) from its rest position. A fluid has no rest position. Air flowing past a car's bumper never swings back to where it came from; it moves on for good, so displacement from rest is meaningless. The natural question becomes: how fast, and in which direction, is the fluid moving here, now? The answer is the velocity field:

Where \(\mathbf{x} = (x, y, z)\) denotes a position in space (m), \(t\) denotes time (s), and \(u_x\), \(u_y\) and \(u_z\) denote the fluid's velocity (m/s) along the \(x\), \(y\) and \(z\) directions at that position and time. Bold letters mark vectors.

The velocity field says where the fluid is going right now, but not how that will change. Predicting the change requires the acceleration, which means \(F = ma\), which in turn needs a force and a mass. For the string, those were the tension and the linear density. A fluid has direct counterparts:

Pressure \(p\) is the push that the surrounding fluid exerts on a parcel, measured as force per area. Only pressure differences matter: a parcel pushed equally hard from all sides feels no net force. Sucking on a straw lowers the pressure in your mouth, and the higher atmospheric pressure on the drink's surface pushes the drink up the straw.

Density \(\rho\) is mass per volume, so it measures inertia. Air has about 1.2 kg/m³ and water about 1000 kg/m³, so the same pressure difference accelerates air roughly 800 times more than water.

The forces are the topic of later chapters. This chapter is about the \(a\).

Before computing an acceleration, we have to decide what exactly we are describing, and there are two choices.

The Lagrangian viewpoint follows the matter. Think of a GPS tracker glued to a buoy drifting down a river: it reports where that one buoy is at every moment. For a fluid parcel, this is a path \(\mathbf{X}(t)\), and since the parcel moves with the fluid, its velocity is whatever the velocity field says at its current location:

Where \(\mathbf{X}(t)\) denotes the position of one specific parcel at time \(t\), \(d\mathbf{X}/dt\) denotes the parcel's velocity, and \(\mathbf{u}(\mathbf{X}(t), t)\) denotes the velocity field looked up at the parcel's current position. This has the structure of the Neural ODE from Post 1, \(dz/dt = f(z(t), t)\): the parcel's position plays the role of the state \(z\), and the velocity field plays the role of the rate function \(f\). Given a starting position and the velocity field, Euler steps trace out the parcel's path.

The Eulerian viewpoint watches fixed points instead. Think of a speed sensor bolted to a bridge pillar: it stays put and reports the velocity of whatever water happens to pass it. That is \(\mathbf{u}(\mathbf{x}, t)\) at a fixed \(\mathbf{x}\), exactly like the grid values \(u_i^n\) of Posts 2–6, where grid point \(x_i\) never moves and only its value changes. The Navier–Stokes equations are usually written in Eulerian form.

This creates a problem. \(F = ma\) is a law about one fixed lump of mass: the force acting on a lump of matter determines the acceleration of that same lump. So the \(a\) in \(F = ma\) is the Lagrangian acceleration of a parcel. A grid point, however, is not a lump of mass. Fluid constantly enters and leaves it: the water at the bridge pillar now is further downstream a second later, replaced by new water. The rate at which the velocity changes at a fixed point, \(\partial \mathbf{u}/\partial t\), is therefore not the acceleration of any single parcel.

The wave equation never ran into this. The pieces of string only oscillate up and down around their rest position, so the segment at \(x\) stays at \(x\) forever. "The segment at \(x\)" and "the point \(x\)" were the same matter, so \(\partial^2 u/\partial t^2\) at a fixed \(x\) was that segment's acceleration, and the Lagrangian law sat on the Eulerian grid for free. Fluid parcels leave for good, so we need a bridge: a way to compute a parcel's acceleration from the Eulerian velocity field.

The bridge starts from the question we want answered: as time passes, how does the velocity of this one moving piece of matter change? The recipe follows directly. Take the parcel's path \(\mathbf{X}(t)\), look up the velocity field there, \(\mathbf{u}(\mathbf{X}(t), t)\), and differentiate with respect to \(t\). By Equation (2), \(\mathbf{u}(\mathbf{X}(t), t)\) is the parcel's velocity, so its time derivative is the parcel's acceleration. It gets its own name and symbol:

Where \(D\mathbf{u}/Dt\) denotes the material derivative of the velocity, \(:=\) denotes "is defined as", and \(\mathbf{X}(t)\) is the path of the parcel we follow. The recipe works for any parcel path; nothing about a particular flow went into it.

Why \(\mathbf{x}\) becomes \(\mathbf{X}(t)\). Holding the position fixed while \(t\) advances means watching one spot as different parcels pass through it, like the bridge sensor. That answers a different question, and its answer is \(\partial \mathbf{u}/\partial t\). We asked about the same matter over time, and that matter moves, so the position where we look up the field has to move with it.

Differentiating with respect to \(t\) does not freeze the position here. A partial derivative holds the other variables fixed, but that only makes sense for variables that are independent of \(t\). Once we declare \(\mathbf{x} = \mathbf{X}(t)\), the position changes whenever \(t\) does, so it cannot be held fixed.

Three symbols, three meanings. \(\partial \mathbf{u}/\partial t\) is the rate of change of the field at a fixed point. \(d\mathbf{X}/dt\) is an ordinary derivative of a parcel's own property, its position, which depends on \(t\) alone. \(D\mathbf{u}/Dt\) follows a parcel through a field. It gets its own symbol because it silently contains the substitution \(\mathbf{x} \to \mathbf{X}(t)\); the capital \(D\) is the reminder that a field is being read along a moving path.

An Eulerian input gives a Lagrangian answer. The input is the Eulerian field \(\mathbf{u}(\mathbf{x}, t)\), and the result describes one specific parcel and its acceleration. The same field therefore answers two different questions: "what is the velocity here?", by reading it at a point, and "what is the acceleration of this parcel?", by differentiating it along the parcel's path.

A parcel's velocity can change for two reasons, and \(D\mathbf{u}/Dt\) is simply their sum. To keep things simple, start in one dimension: fluid flows only along \(x\), so the velocity is a single number \(u(x,t)\) (m/s), and the parcel's position is \(X(t)\).

1. The velocity at a fixed spot changes over time. For example, heavy rain upstream makes a river flow faster everywhere, so even a parcel that stayed in place would speed up. This is \(\partial u/\partial t\), the same \(\partial/\partial t\) as in the heat and wave posts.

2. The parcel moves to a spot where the velocity is different. For example, a river that narrows downstream flows faster there, so a drifting parcel speeds up even if nothing changes over time. This effect depends on two things: how much the velocity differs per meter, the spatial slope \(\partial u/\partial x\), and how fast the parcel moves, its drift speed \(dX/dt\). The velocity gained per second is their product. If the velocity rises by 2 m/s per meter and the parcel drifts at 3 m/s, it covers 3 m per second, each meter adds 2 m/s, so it gains 6 m/s every second. The units match an acceleration: \(\frac{\text{m/s}}{\text{m}} \cdot \frac{\text{m}}{\text{s}} = \frac{\text{m}}{\text{s}^2}\). A product is the right combination because the effect must vanish when either factor is zero: a parcel that does not move gains nothing from the slope, and a parcel moving through uniform flow gains nothing from moving.

To compute the material derivative, we differentiate Equation (3), which in one dimension reads \(\frac{d}{dt}\big[u(X(t), t)\big]\). The two reasons from the previous section are already built into this expression, because time reaches \(u\) through two paths.

The direct path. \(t\) sits in the time slot of \(u\) and changes it directly, like the rain speeding up the whole river. This contribution is \(\partial u/\partial t\), the change of the velocity over time at a fixed position.

The path through the position. \(t\) also moves the parcel, and the parcel's position \(X(t)\) feeds into the position slot of \(u\). When the output of one function is the input of another, the chain rule gives the rate of change:

Where \(A\) depends on \(B\) and \(B\) depends on \(t\), \(dA/dB\) denotes how much \(A\) changes per unit of \(B\), and \(dB/dt\) denotes how much \(B\) changes per unit of time.

Here, \(A\) is the velocity and \(B\) is the parcel's position. The first factor is the spatial slope \(\partial u/\partial x\) at the parcel's position: how much the velocity changes per meter. The slope alone does not say how much the parcel's velocity changes; that depends on how far the parcel travels along it, which is the second factor, the drift speed \(dX/dt\) in meters per second. Slope times distance travelled per second is the velocity gained per second. It is a product because the change must vanish if either factor is zero: a parcel that stays put gains nothing from the slope, and a parcel moving through uniform flow gains nothing from moving.

Adding both paths gives

Where \(\partial u/\partial t\) denotes how fast the velocity changes over time at the parcel's position, \(\partial u/\partial x\) denotes how much the velocity changes per meter at the parcel's position, and \(dX/dt\) denotes the parcel's drift speed.

This simplifies further. The parcel's drift speed \(dX/dt\) is just the fluid velocity at the parcel's position, \(u(X(t), t)\), as stated in Equation (2). Replacing \(dX/dt\) with \(u(X(t), t)\), written as just \(u\) for short, gives the material derivative in one dimension:

The velocity \(u\) now appears twice: once as the velocity whose change we measure, inside \(\partial u/\partial x\), and once as the drift speed that carries the parcel along the slope.

So far the parcel moved only along \(x\). In three dimensions, the first reason stays exactly as it is: \(\partial u/\partial t\) is the change at a fixed spot, and a fixed spot has no direction. The second reason does change, because it depends on where the parcel goes. The parcel now drifts in some combination of the \(x\), \(y\) and \(z\) directions, and the flow's velocity can change along each of them. A river, for example, flows faster downstream, slower toward the banks and slower toward the bottom. The change a parcel sees by moving therefore has to account for all three directions at once.

Take \(x\) pointing east (downstream), \(y\) north (across the river) and \(z\) up. The velocity now has three components, so start with one of them, the eastward velocity \(u_x\). How \(u_x\) changes along each direction is collected in one vector, the gradient:

Where \(\nabla\) denotes the gradient, and each entry is the spatial slope of \(u_x\) along one direction. For example, \(\partial u_x/\partial y\) is how much the eastward velocity changes per meter moved northward. The subscript says which velocity component is examined; the three entries are the three directions along which it can change.

The drift is a vector too: the parcel's velocity \(\mathbf{u} = (u_x, u_y, u_z)\), which by Equation (2) says how many meters per second the parcel moves east, north and up. In one dimension, the change by moving was slope times drift speed. In three dimensions, it is slope times drift speed for each direction, summed. That is exactly the dot product of the drift vector with the gradient:

Where \(\cdot\) denotes the dot product: multiply matching entries and add them up. Each of the three terms is the one-dimensional rule applied to one direction. For example, let the eastward velocity grow by 2 m/s per meter northward and stay the same in the other directions, so \(\nabla u_x = (0, 2, 0)\) (m/s)/m. A parcel with \(\mathbf{u} = (3, 0, 0)\) m/s gets \(3 \cdot 0 + 0 \cdot 2 + 0 \cdot 0 = 0\): it only moves east, and \(u_x\) does not change in that direction. A parcel with \(\mathbf{u} = (3, 1, 0)\) m/s gets \(3 \cdot 0 + 1 \cdot 2 + 0 \cdot 0 = 2\) m/s²: only its northward motion counts, because north is the only direction in which \(u_x\) changes.

Adding the first reason gives the rate at which the parcel's eastward velocity changes:

Where \(Du_x/Dt\) denotes the eastward component of the parcel's acceleration.

The components \(u_y\) and \(u_z\) work the same way: the same drift vector \(\mathbf{u}\) is dotted with \(\nabla u_y\) and \(\nabla u_z\). Because "\(\mathbf{u} \cdot \nabla\)" appears in front of every component, it is written once as an operator:

Where \(\partial/\partial x\), \(\partial/\partial y\) and \(\partial/\partial z\) are derivatives waiting for a field to act on. Applied to any velocity component, \((\mathbf{u} \cdot \nabla)\) gives that component's change by moving. Applying it to all three components at once gives the full material derivative:

Where \(\partial \mathbf{u}/\partial t\) denotes the local term, the change of all three velocity components at a fixed spot, and \((\mathbf{u} \cdot \nabla)\mathbf{u}\) denotes the convective term, the vector \(\big((\mathbf{u} \cdot \nabla)u_x,\; (\mathbf{u} \cdot \nabla)u_y,\; (\mathbf{u} \cdot \nabla)u_z\big)\), the change each component sees because the parcel moves. In one dimension, \(u_y = u_z = 0\) and nothing depends on \(y\) or \(z\), so the convective term reduces to \(u \, \partial u/\partial x\) and Equation (11) becomes Equation (6).

This chapter covered the \(a\) in \(F = ma\). The next chapter covers conservation of mass, which leads to the condition \(\nabla \cdot \mathbf{u} = 0\) that it builds up from scratch, and the one after that derives the pressure and viscous forces on a parcel. Assembling the pieces gives the Navier–Stokes equations in the form "density × material derivative = forces per unit volume".