Nine components, not six

The second interactive companion to our paper:

Part one ended with a recipe in one dimension: average a conservation law over a grid cell, and what the coarse simulation leaves unresolved is a single flux that can be written down exactly. A stress tensor in three dimensions looks like the same thing with more indices.

It is not quite. The textbook grants the unresolved stress two properties: it is symmetric, and it is local. Neither survives the trip to a grid. Both are lost for the same reason, and the reason is the shape of a box.

A box has more than one kind of face

In one dimension a cell has two end points, and the filter-swap traded a derivative for the difference between them. In three dimensions a cell has six faces, in three pairs, and the same trade is the divergence theorem. Write $\sigma_{ij}$ for the stress, the flux of $i$-momentum in direction $j$. Then, for each $i$,

On the left, the bar is the average over a cell of the force on the fluid. On the right, $\langle \cdot \rangle_j$ is the average over a face normal to direction $j$, and $\partial_j^h$ is the difference between the two opposite faces, divided by $h$. Nothing is approximated. The momentum in a cell changes by what passes through its faces.

Now look at which component sits on which face. Each term has its own: $\sigma_{i1}$ is averaged over the faces normal to $x_1$, and $\sigma_{i2}$ over those normal to $x_2$. A grid cell never sees “the stress”. It sees each column of it through a different window.

Same stress, two windows

The stress of a fluid is symmetric. $\sigma_{12}$ and $\sigma_{21}$ are the same number at every point.

Then the box averages them. $\sigma_{12}$ enters the equation for $u_1$ through a face normal to $x_2$. $\sigma_{21}$ enters the equation for $u_2$ through a face normal to $x_1$. On a staggered grid these two faces cross at right angles, at a corner of the cell. Here is a patch of that one field, taken from a simulation, with the two faces drawn on it.

The two numbers are averages of the same field around the same point, and they disagree. They agree when the cell shrinks to a single value, which is the textbook limit, and wherever the field happens to look the same along both faces. Switch the map to see where they differ: wherever the stress has structure that is finer than the cell and longer in one direction than in the other. In a turbulent flow, that is most places.

So the unresolved stress on a grid has $\tau_{12} \neq \tau_{21}$. In three dimensions it has nine independent components, not six.

The pressure will not stay in its lane

An incompressible flow comes with a constraint, $\partial_j u_j = 0$, and a pressure $p$ that enforces it. That pressure is not set by the fluid around a point. It is set by the whole flow at once. Change the stress in one place, and the pressure adjusts everywhere, so that the velocity stays divergence-free.

Most of the response stays close to its source. But it fades slowly, and there is no distance beyond which it is gone.

In the textbook derivation this is harmless: filter the equations, and the filtered pressure becomes the pressure of the filtered flow. On a grid it is not harmless, for the reason this post keeps returning to. A coarse simulation computes its own pressure, to keep its own field divergence-free on its own grid. That pressure can absorb the part of the true one that looks the same on every face of a cell. But the box gives each diagonal component its own face, $\langle p \rangle_1$ for the first and $\langle p \rangle_2$ for the second, and the difference between those is nobody’s pressure. It has to go into the stress:

This is the whole result. The first term is the true stress, pressure included, with each component averaged over the face that the box assigns to it. The second is the stress that your scheme computes from the coarse field $\bar{u}$, exactly as in part one. And since $p$ depends on the flow everywhere, so does $\tau$. The exact unresolved stress is not local.

(The same theorem explains why the coarse field needs a pressure of its own. Applied to $\partial_j u_j = 0$, it says that the face averages of the velocity are divergence-free on the grid. Momentum conservation, above, evolves cell averages, and those are not. So $\bar{u}$ is the cell average, projected back onto fields that the coarse grid considers divergence-free.)

The ladder, live

Does any of this matter, or are these corrections to corrections? The test is the one from part one, moved up a dimension. Below, a two-dimensional turbulent flow is simulated directly, and six coarse simulations run alongside it on a grid five times coarser in each direction. All six use the same scheme. They differ only in what they are fed from the DNS.

Read the rows from the top.

• Nothing and Textbook stress end close together. The classical subgrid stress, computed here from the full DNS and so as good as the classical theory can make it, barely helps.
• … plus your scheme is the largest single step. This is part one again: the numerical stress belongs inside the residual.
• Exact, made symmetric is the exact stress with its two shear components replaced by their mean. This is as close as a symmetric closure can get to the exact stress, and the coarse run it drives is still off by several percent.
• Exact, pressure left out keeps the asymmetry but leaves the pressure to the coarse grid alone.
• Exact stays on the filtered DNS to round-off.

In two dimensions the stress has four components instead of nine, and a symmetric one has three instead of six. The argument does not change. The paper runs this experiment in three dimensions, at a resolution no browser will manage, and finds the same ladder: 26% with nothing, 21% with the textbook stress, 6.4% with the scheme inside, 2.9% for the exact stress made symmetric, and machine precision for the exact one. The pressure row is not in the paper. It is new here, and so far only in two dimensions.

What this asks of a closure

Nine outputs. A closure that models the stress itself, for a finite-volume code, should predict $\tau_{12}$ and $\tau_{21}$ separately. Building symmetry into the model, as the tensor-basis closures do by construction, means the best it can represent is the stress of the fourth row. (This is symmetry of the tensor, $\tau_{ij} = \tau_{ji}$. It is a different thing from the rotation and reflection symmetries of another post.)

Possibly a longer reach. Through the pressure, the exact stress depends on the whole flow. How much of that a model needs is open. The paper measures what symmetry costs, not what locality costs, and in the run above leaving the pressure out costs less than making the stress symmetric.

The same caveat as before. This $\tau$ belongs to one grid and one scheme. A staggered grid puts the two shear components at the same corner, which is what made the comparison above so clean. Another layout will hand out its faces differently, and its exact stress has to be derived for it.

One more box

An interval gave a flux that knows your scheme. A cube gave a stress with nine components and a long reach. The last box is a brick in space-time: one cell wide, one time step tall. Forward Euler turns out to be the same identity along the time axis, and it adds exactly one more term. That is part three.