9.2 Finite Difference Methods
291
The PDE is valid at all spatial points x ∈ Ω, but we may relax this condition and
demand that it is fulfilled at the internal mesh points only, x 1 , . . . , x N−1 :
∂u(x i , t)
∂t
= β
∂ 2 u(x i , t)
∂x 2
+ g(x i , t), i = 1, . . . , N − 1 .
(9.5)
Now, at any point x i we can approximate the second-order derivative by a finite
difference:
∂ 2 u(x i , t)
∂x 2
≈
u(x i+1 , t) − 2u(x i , t) + u(x i−1 , t)
Δx 2
.
(9.6)
It is common to introduce a short notation u i (t) for u(x i , t), i.e., u approximated at
some mesh point x i in space. With this new notation we can, after inserting (9.6)
in (9.5), write an approximation to the PDE at mesh point (x i , t) as
du i (t)
dt
= β
u i+1 (t) − 2u i (t) + u i−1 (t)
Δx 2
+ g i (t), i = 1, . . . , N − 1 .
(9.7)
Note that we have adopted the notation g i (t) for g(x i , t) too.
What is (9.7)? This is nothing but a system of ordinary differential equations in
N −1 unknowns u 1 (t), . . . , u N−1 (t)! In other words, with aid of the finite difference
approximation (9.6), we have reduced the single PDE to a system of ODEs, which
we know how to solve. In the literature, this strategy is called the method of lines.
We need to look into the initial and boundary conditions as well. The initial
condition u(x, 0) = I (x) translates to an initial condition for every unknown
function u i (t): u i (0) = I (x i ), i = 0, . . . , N. At the boundary x = 0 we need
an ODE in our ODE system, which must come from the boundary condition at
this point. The boundary condition reads u(0, t) = s(t). We can derive an ODE
from this equation by differentiating both sides: u
0 (t) = s (t). The ODE system
above cannot be used for u
0 since that equation involves some quantity u
−1 outside
the domain. Instead, we use the equation u
0 (t) = s (t) derived from the boundary
condition. For this particular equation we also need to make sure the initial condition
is u 0 (0) = s(0) (otherwise nothing will happen: we get u = 283 K forever).
We remark that a separate ODE for the (known) boundary condition u 0 = s(t)
is not strictly needed. We can just work with the ODE system for u 1 , . . . , u N , and
in the ODE for u 0 , replace u 0 (t) by s(t). However, these authors prefer to have an
ODE for every point value u i , i = 0, . . . , N, which requires formulating the known
boundary at x = 0 as an ODE. The reason for including the boundary values in the
ODE system is that the solution of the system is then the complete solution at all
mesh points, which is convenient, since special treatment of the boundary values is
then avoided.
The condition ∂u/∂x = 0 at x = L is a bit more complicated, but we can
approximate the spatial derivative by a centered finite difference:
∂u
∂x
i=N
≈
u N+1 − u N−1
2Δx
= 0 .
291
The PDE is valid at all spatial points x ∈ Ω, but we may relax this condition and
demand that it is fulfilled at the internal mesh points only, x 1 , . . . , x N−1 :
∂u(x i , t)
∂t
= β
∂ 2 u(x i , t)
∂x 2
+ g(x i , t), i = 1, . . . , N − 1 .
(9.5)
Now, at any point x i we can approximate the second-order derivative by a finite
difference:
∂ 2 u(x i , t)
∂x 2
≈
u(x i+1 , t) − 2u(x i , t) + u(x i−1 , t)
Δx 2
.
(9.6)
It is common to introduce a short notation u i (t) for u(x i , t), i.e., u approximated at
some mesh point x i in space. With this new notation we can, after inserting (9.6)
in (9.5), write an approximation to the PDE at mesh point (x i , t) as
du i (t)
dt
= β
u i+1 (t) − 2u i (t) + u i−1 (t)
Δx 2
+ g i (t), i = 1, . . . , N − 1 .
(9.7)
Note that we have adopted the notation g i (t) for g(x i , t) too.
What is (9.7)? This is nothing but a system of ordinary differential equations in
N −1 unknowns u 1 (t), . . . , u N−1 (t)! In other words, with aid of the finite difference
approximation (9.6), we have reduced the single PDE to a system of ODEs, which
we know how to solve. In the literature, this strategy is called the method of lines.
We need to look into the initial and boundary conditions as well. The initial
condition u(x, 0) = I (x) translates to an initial condition for every unknown
function u i (t): u i (0) = I (x i ), i = 0, . . . , N. At the boundary x = 0 we need
an ODE in our ODE system, which must come from the boundary condition at
this point. The boundary condition reads u(0, t) = s(t). We can derive an ODE
from this equation by differentiating both sides: u
0 (t) = s (t). The ODE system
above cannot be used for u
0 since that equation involves some quantity u
−1 outside
the domain. Instead, we use the equation u
0 (t) = s (t) derived from the boundary
condition. For this particular equation we also need to make sure the initial condition
is u 0 (0) = s(0) (otherwise nothing will happen: we get u = 283 K forever).
We remark that a separate ODE for the (known) boundary condition u 0 = s(t)
is not strictly needed. We can just work with the ODE system for u 1 , . . . , u N , and
in the ODE for u 0 , replace u 0 (t) by s(t). However, these authors prefer to have an
ODE for every point value u i , i = 0, . . . , N, which requires formulating the known
boundary at x = 0 as an ODE. The reason for including the boundary values in the
ODE system is that the solution of the system is then the complete solution at all
mesh points, which is convenient, since special treatment of the boundary values is
then avoided.
The condition ∂u/∂x = 0 at x = L is a bit more complicated, but we can
approximate the spatial derivative by a centered finite difference:
∂u
∂x
i=N
≈
u N+1 − u N−1
2Δx
= 0 .
