Biomedical Engineering Reference
In-Depth Information
This contrasts with the level set method, where each time step requires an
additional pass over the mesh to evolve the level set function in time.
The implementation of the fast marching method also uses the numerical
flux functions discussed in Section 4.2.2; however, in this case, only one-sided
differences such as Godunov's method may be used. For example, suppose the
values of
φ
i
−
1
,
j
,
φ
i
,
j
+
1
are already determined, and we wish to compute
φ
ij
. Then
Eq. 4.21 is discretized using one-sided differences to obtain
F
ij
((
D
−
x
φ
ij
)
2
+
(
D
+
y
φ
ij
)
2
)
=
1
.
(4.22)
This equation can be rewritten as a quadratic in terms of the unknown
φ
ij
:
1
x
2
+
y
2
−
2
φ
i
−
1
,
j
i
−
1
,
j
x
2
+
i
,
j
+
1
y
2
−
φ
φ
1
F
2
=
0
.
(4.23)
In most cases, solving Eq. 4.23 will produce two solutions, one which is less than
the values of
φ
i
−
1
,
j
,
φ
i
,
j
+
1
, and one which is greater. The larger of the two values
is always chosen because of the causality assumption made by this method;
values that are unknown are always greater than the known values.
Occasionally, Eq. 4.23 will not have any real roots. In that case, each of the
coordinate directions is considered separately. For example, if we consider the
x
-direction, we assume that
∂φ/∂
y
=
0, and then discretize Eq. 4.21 to get
1
x
2
+
φ
i
,
j
+
1
i
,
j
φ
φ
i
,
j
+
y
2
F
ij
D
−
x
φ
ij
=
1
.
(4.24)
This equation is linear and is easily solved for
φ
ij
. Similarly, the
y
-direction is
considered, and the smaller of the two solutions is taken as the new estimate
for
φ
ij
.
The key to solving Eq. 4.21 in one pass is to traverse the mesh in the proper
order. The grid points must be evaluated in the order of increasing
t
. This is
accomplished by using a sorted heap which always keeps track of which grid
point is to be evaluated next. To begin, the set of grid points is divided into
three disjoint sets, the
accepted
points
A
, the
tentative
points
T
, and the
distant
points
D
. The accepted points in
A
are the points
x
ij
for which the computed
value of
φ
ij
is already determined. The tentative points in
T
are the points
x
ij
for which a tentative value for
φ
ij
is computed. The remainder of the points
are in the set
D
. One by one, points in
T
are taken, in order of increasing
value of
φ
ij
, from the set
T
into
A
. Each time, points
φ
ij
in
D
which become
adjacent to points in the set
A
are moved into the set
T
and a tentative value
