Von Neumann’s Q, often referred to as ‘artificial viscosity’ requires knowing difference in the volume of zones between cycles. The volume is directly calculated from the vertex positions. Between sweeps in the current code the positions have been updated by the latest velocities. The volume change is thus not available without extra storage of 2mn floats. Many nit-picking little changes! This is the sort of change that makes the promise of languages such as Sisal so great. I think that I shall write a version that displays the logic of the equations first.