Graphics Reference
In-Depth Information
This is the compatibility condition, without which the PDE for pressure
has no solution. (Note of course that if there is a free surface as well, no
condition is needed.)
More precisely, we want the discrete linear system to have a solution.
Even if the solid wall velocity does exactly integrate to zero, once discretized
on a grid there may be numerical errors that cause the linear solve to fail.
Thus we need a routine which can “correct” the right-hand side of the
linear system to make it consistent.
We can view this as a pure linear algebra problem, as outlined by Guen-
delman et al. [Guendelman et al. 05]. That is, given a singular
A
and a
vector
b
,wewanttoremovefrom
b
any components which lie outside the
range of
A
: the modified linear system
Ax
=
b
then has a solution. If
A
is symmetric, which it is in our case, this is equivalent to removing com-
ponents from
b
that lie in the null-space of
A
. Luckily, for the pressure
problem on a MAC grid, this null-space is very simple: for every set of
connected grid cells with no free surface conditions, there is a null-space
basis vector consisting of unit pressures in those voxels. By “connected,”
we mean with respect to the graph of the matrix
A
: two grid cells are con-
nected if there is a non-zero entry in
A
corresponding to the pair; finding
these connected sets is a simple graph traversal problem. Note that the sets
are necessarily disjoint, so each of the null-space vectors we find is orthog-
onal to the rest. Therefore we can eciently remove the incompatibilities
in
b
by orthogonal projection. Boiled down to basics, this simply means
adding up the entries of
b
in the set, dividing by the number of entries, and
subtracting this mean off of each entry. It is imperative here to use double
precision floating-point in calculating the sum, as otherwise the accumu-
lated rounding errors can seriously undermine the method (necessitating
repeating the process and/or subtracting the sum off a single randomly
chosen entry). Also note that this process works for both the classic vox-
elized solve and the more accurate variational method from the previous
section. From a physical standpoint, this causes us to find a solution that
strikes a balance between respecting the boundary conditions and staying
divergence-free.
