Graphics Reference
In-Depth Information
Add body forces to get
u
B
=
u
A
+Δ
tg
.
•
Solve for pressure to get a divergence-free
u
n
+1
= project(
u
B
,
Δ
t
)
•
satisfying the inviscid boundary conditions.
Note that there is the possibility that the final velocity field
u
n
+1
will no
longer exactly satisfy the viscous boundary conditions, just the inviscid
ones—however, it will be off by at most
O
(Δ
t
), which we can live with.
Now that we have a full time step in place, let's determine how to
actually apply viscous forces. Once again staggered grids are our friend.
Let's take a look at the contribution to the horizontal component of velocity,
given the viscous stress tensor
τ
:
∂τ
11
∂x
.
+
∂τ
12
∂y
+
∂τ
13
∂z
=
u
n
+
Δ
t
ρ
u
V
Since
u
is located in the grid at, say, (
i
+1
/
2
,j,k
), it's natural to ask
for
τ
11
to be at grid cell centers (
i, j, k
), for
τ
12
to be at the edge-center
(
i
+1
/
2
,j
+1
/
2
,k
), and for
τ
13
at (
i
+1
/
2
,j,k
+1
/
2). This gives an elegant
discretization:
τ
11
τ
11
i,j,k
i
+1
,j,k
−
u
i
+1
/
2
,j,k
=
u
i
+1
/
2
,j,k
+
Δ
t
ρ
Δ
x
.
τ
12
τ
12
i
+1
/
2
,j−
1
/
2
,k
τ
13
τ
13
i
+1
/
2
,j,k−
1
/
2
i
+1
/
2
,j
+1
/
2
,k
−
i
+1
/
2
,j,k
+1
/
2
−
+
+
Δ
x
Δ
x
Similarly for the other components of velocity:
τ
12
τ
12
i−
1
/
2
,j
+1
/
2
,k
i
+1
/
2
,j
+1
/
2
,k
−
v
i,j
+1
/
2
,k
=
v
i,j
+1
/
2
,k
+
Δ
t
ρ
Δ
x
,
τ
23
τ
23
i,j
+1
/
2
,k−
1
/
2
+
τ
22
i,j
+1
,k
−
τ
22
i,j,k
i,j
+1
/
2
,k
+1
/
2
−
+
Δ
x
Δ
x
τ
13
τ
13
i−
1
/
2
,j,k
+1
/
2
i
+1
/
2
,j,k
+1
/
2
−
w
i,j,k
+1
/
2
=
w
i,j,k
+1
/
2
+
Δ
t
ρ
Δ
x
+
τ
23
i,j
+1
/
2
,k
+1
/
2
−
τ
23
i,j−
1
/
2
,k
+1
/
2
τ
33
τ
33
i,j,k
i,j,k
+1
−
+
.
Δ
x
Δ
x
Note that these formulas make use of the symmetry of
τ
, e.g.,
τ
12
=
τ
21
.In
2D, they simplify the obvious way. But how do we determine the values of
τ
on the staggered grid, and how do we simplify when viscosity is constant?