Biomedical Engineering Reference
In-Depth Information
for
k
=
i
,
i
+
1,
=
j
,
j
+
1. This gives 16 equations for the 16 unknown coef-
ficients
a
m
,
n
. Solving for
a
m
,
n
makes
p
(
x
,
y
) a bicubic interpolating function
of
φ
(
x
,
y
) on the rectangle bounded by the corners
x
i
,
j
,
x
i
+
1
,
j
,
x
i
,
j
+
1
, and
x
i
+
1
,
j
+
1
.
Since
φ
is only known on the mesh points, the values for the derivatives of
φ
must be approximated. We use second-order finite difference approximations
for the derivatives of
φ
:
1
2
x
(
φ
(
x
m
+
1
,
n
)
−
φ
(
x
m
−
1
,
n
))
∂φ
∂
x
(
x
m
,
n
)
≈
∂φ
∂
y
(
x
m
,
n
)
≈
1
2
y
(
φ
(
x
m
,
n
+
1
)
−
φ
(
x
m
,
n
−
1
))
2
1
4
x
y
(
φ
(
x
m
+
1
,
n
+
1
)
−
φ
(
x
m
−
1
,
n
+
1
)
−
φ
(
x
m
+
1
,
n
−
1
)
+
φ
(
x
m
−
1
,
n
−
1
))
φ
∂
x
∂
y
(
x
m
,
n
)
≈
∂
for
m
=
i
,
i
+
1 and
n
=
j
,
j
+
1. Thus, construction of the interpolant
p
requires
all the points shown in Fig. 4.7. Higher order local approximations can be made
using higher order finite difference approximations and using a larger set of grid
points around the box where the interpolant is used.
Now, given the interpolating function
p
(
x
,
y
) in the domain [
x
i
,
x
i
+
1
]
×
[
y
j
,
y
j
+
1
], and given a point (
x
0
,
y
0
) in that domain, we compute the distance
between (
x
0
,
y
0
) and the zero level curve of
p
(
x
,
y
). The point (
x
1
,
y
1
)onthe
zero level curve closest to (
x
0
,
y
0
) must satisfy two conditions:
p
(
x
1
,
y
1
)
=
0
,
(4.27)
∇
p
(
x
1
,
y
1
)
×
((
x
0
,
y
0
)
−
(
x
1
,
y
1
))
=
0
.
(4.28)
Equation 4.27 is a requirement that (
x
1
,
y
1
) must be on the interface. Equa-
tion 4.28 is a requirement that the interface normal, given by
∇
p
(
x
1
,
y
1
),
must be aligned with the line through the points (
x
0
,
y
0
) and (
x
1
,
y
1
). Equa-
tions 4.27 and 4.28 are solved simultaneously using Newton's method. Typ-
ically, less than five iterations are necessary in order to achieve sufficient
accuracy.
Given the front speed
F
(
x
1
,
y
1
) and the initial distance to the front,
d
=
(
x
1
,
y
1
)
−
(
x
0
,
y
0
)
, the initial value for a point adjacent to the initial front
for the general fast marching method solving Eq. 4.21 is
d
/
F
.
