Graphics Reference
In-Depth Information
scale = dt / (rho*dx*dx);
loop over i,j,k:
if cell(i,j,k)==FLUID and cell(i+1,j,k)==FLUID:
Adiag(i,j,k) += scale;
Adiag(i+1,j,k) += scale;
Aplusi(i,j,k) = -scale;
else if cell(i,j,k)==FLUID and cell(i+1,j,k)==EMPTY:
Adiag(i,j,k) += scale;
if cell(i,j,k)==FLUID and cell(i,j+1,k)==FLUID:
Adiag(i,j,k) += scale;
Adiag(i,j+1,k) += scale;
Aplusj(i,j,k) = -scale;
else if cell(i,j,k)==FLUID and cell(i,j+1,k)==EMPTY:
Adiag(i,j,k) += scale;
if cell(i,j,k)==FLUID and cell(i,j,k+1)==FLUID:
Adiag(i,j,k) += scale;
Adiag(i,j,k+1) += scale;
Aplusk(i,j,k) = -scale;
else if cell(i,j,k)==FLUID and cell(i,j,k+1)==EMPTY:
Adiag(i,j,k) += scale;
Figure 4.4.
Setting up the matrix entries for the pressure equations.
A ( i,j ) , ( i +1 ,j ) and A ( i,j ) , ( i,j +1) . We could call these entries Adiag(i,j) ,
Aplusi(i,j) ,and Aplusj(i,j) in our code. In three dimensions, we
would similarly have Adiag(i,j,k) , Aplusi(i,j,k) , Aplusj(i,j,k) ,and
Aplusk(i,j,k) . When we need to refer to an entry like A ( i,j ) , ( i− 1 ,j ) we
use symmetry and instead refer to A ( i− 1 ,j ) , ( i,j ) = Aplusi(i-1,j) .SeeFig-
ure 4.4 for pseudocode to set up the matrix in this structure.
4.3.2 The Conjugate Gradient Algorithm
The matrix A is a very well-known type of matrix, sometimes referred to
as the five- or seven-point Laplacian matrix, in two or three dimensions,
respectively. It has been exhaustively studied, serves as the subject of
countless numerical linear algebra papers, and is the prototypical first ex-
ample of a sparse matrix in just about any setting. More effort has been
put into solving linear systems with this type of matrix than probably
all other sparse matrices put together! We won't look very far into this
Search WWH ::




Custom Search