Graphics Reference
In-Depth Information
we increment the negative divergence measured on the right-hand side with
a term involving the difference between fluid and solid velocity. This can
be implemented in code as an additional loop to modify
rhs
,asshownin
Figure 4.3.
4.3.1
Putting It In Matrix-Vector Form
We have now defined a large system of linear equations for the unknown
pressure values. We can conceptually think of it as a large coecient
matrix,
A
, times a vector consisting of all pressure unknowns,
p
,equaltoa
vector consisting of the
negative
divergences in each fluid grid cell,
b
(with
appropriate modifications at solid wall boundaries):
Ap
=
b.
(4.12)
In the implementation we have discussed so far, of course,
p
and
b
are
logically stored in a two or three-dimensional grid structure (since each
entry corresponds to a grid cell).
We needn't store
A
directly as a matrix. Notice that each row of
A
corresponds to one equation, i.e., one fluid cell. For example, if grid cell
(
i, j, k
) is fluid, then there will be a row of the matrix that we can label
with the indices (
i, j, k
). The entries in that row are the coecients of all
the pressure unknowns in that equation: almost all of these are zero except
possibly for the seven entries corresponding to
p
i,j,k
and its six neighbors,
p
i±
1
,j,k
,
p
i,j±
1
,k
,and
p
i,j,k±
1
. (In two dimensions there are at most four
neighbors of course.) We only have non-zeros (
i, j, k
) and its fluid cell
neighbors. It is of course pointless to store all the zeros: this is a
sparse
matrix.
Let's take a closer look at
A
. In the equation for (
i, j, k
), the coecients
for neighboring fluid cells are all equal to
Δ
t/
(
ρ
Δ
x
2
), and if there are
n
i,j,k
fluid- or air-cell neighbors the coecient for
p
i,j,k
is
n
i,j,k
Δ
t/
(
ρ
Δ
x
2
).
One of the nice properties of the matrix
A
is that it is symmetric. For
example,
A
(
i,j,k
)
,
(
i
+1
,j,k
)
, the coecient of
p
i
+1
,j,k
in the equation for grid
cell (
i, j, k
), has to be equal to
A
(
i
+1
,j,k
)
,
(
i,j,k
)
. Eitherit'szeroifoneof
those two cells is not fluid, or it's the same non-zero value. This symmetry
property will hold even with the more advanced discretization at the end
of the chapter. Thus we only have to store half of the non-zero entries in
A
, since the other half are just copies!
This leads us to the following structure for storing
A
. In two dimen-
sions, we will store three numbers at every grid cell: the diagonal entry
A
(
i,j
)
,
(
i,j
)
and the entries for the neighboring cells in the positive directions,
−
