**BVP of ODE**

^{1}

**Boundary-Value Problems for Ordinary** **Differential Equations**

**NTNU**

**Tsung-Min Hwang**
**December 20, 2003**

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{2}

1 – Mathematical Theories . . . 4

2 – Finite Difference Method For Linear Problems . . . 15

2.1 – The Finite Difference Formulation . . . 15

2.2 – Convergence Analysis . . . 19

3 – Shooting Methods . . . 22

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{3}

The two-point boundary-value problems (BVP) considered in this chapter involve a

second-order differential equation together with boundary condition in the following form:

###

###

###

### y

^{00}

### = f (x, y, y

^{0}

### )

### y(a) = α, y(b) = β

(1)

The numerical procedures for finding approximate solutions to the initial-value problems can not be adapted to the solution of this problem since the specification of conditions involve two different points,

### x = a

^{and}

### x = b

. New techniques are introduced in this chapter for handling problems (1) in which the conditions imposed are of a boundary-value rather than an initial-value type.**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{3}

The two-point boundary-value problems (BVP) considered in this chapter involve a

second-order differential equation together with boundary condition in the following form:

###

###

###

### y

^{00}

### = f (x, y, y

^{0}

### )

### y(a) = α, y(b) = β

(1)

The numerical procedures for finding approximate solutions to the initial-value problems can not be adapted to the solution of this problem since the specification of conditions involve two different points,

### x = a

^{and}

### x = b

^{.}

New techniques are introduced in this chapter for handling problems (1) in which the conditions imposed are of a boundary-value rather than an initial-value type.

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{3}

The two-point boundary-value problems (BVP) considered in this chapter involve a

second-order differential equation together with boundary condition in the following form:

###

###

###

### y

^{00}

### = f (x, y, y

^{0}

### )

### y(a) = α, y(b) = β

(1)

The numerical procedures for finding approximate solutions to the initial-value problems can not be adapted to the solution of this problem since the specification of conditions involve two different points,

### x = a

^{and}

### x = b

. New techniques are introduced in this chapter for handling problems (1) in which the conditions imposed are of a boundary-value rather than an initial-value type.**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{4}

**1 – Mathematical Theories**

Before considering numerical methods, a few mathematical theories about the two-point boundary-value problem (1), such as the existence and uniqueness of solution, shall be discussed in this section.

**Theorem 1 Suppose that**

### f

*in (1) is continuous on the set*

### D = {(x, y, y

^{0}

### )|a ≤ x ≤ b, −∞ < y < ∞, −∞ < y

^{0}

### < ∞} ,

*and that*

### ∂f

### ∂y

^{and}### ∂f

### ∂y

^{0}

*are also continuous on*

### D

^{. If}*1.*

### ∂f

### ∂y (x, y, y

^{0}

### ) > 0

^{for all}### (x, y, y

^{0}

### ) ∈ D

^{, and}*2. a constant*

### M

*exists, with*

### ∂f

### ∂y

^{0}

### (x, y, y

^{0}

### )

### ≤ M

^{,}### ∀ (x, y, y

^{0}

### ) ∈ D

^{,}*then (1) has a unique solution.*

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{5}

When the function

### f (x, y, y

^{0}

### )

has the special form### f (x, y, y

^{0}

### ) = p(x)y

^{0}

### + q(x)y + r(x),

the differential equation become a so-called linear problem. The previous theorem can be simplified for this case.

**Corollary 1 If the linear two-point boundary-value problem**

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x) y(a) = α, y(b) = β

*satisfies*

*1.*

### p(x), q(x)

^{, and}### r(x)

*are continuous on*

### [a, b]

^{, and}*2.*

### q(x) > 0

^{on}### [a, b]

^{,}*then the problem has a unique solution.*

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{6}

Many theories and application models consider the boundary-value problem in a special form as follows.

###

###

###

### y

^{00}

### = f (x, y)

### y(0) = 0, y(1) = 0

We will show that this simple form can be derived from the original problem by simple

techniques such as changes of variables and linear transformation. To do this, we begin by changing the variable. Suppose that the original problem is

###

###

###

### y

^{00}

### = f (x, y)

### y(a) = α, y(b) = β

^{(2)}

where

### y = y(x)

^{. Now let}

### λ = b − a

and define a new variable### t = x − a

### b − a = 1

### λ (x − a).

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{7}

That is,

### x = a + λt

. Notice that### t = 0

corresponds to### x = a

^{, and}

### t = 1

corresponds to### x = b

. Then we define### z(t) = y(a + λt) = y(x)

with

### λ = b − a

. This gives### z

^{0}

### (t) = d

### dt z(t) = d

### dt y(a + λt) = d

### dx y(x) d

### dt (a + λt)

### = λy

^{0}

### (x)

and, analogously,

### z

^{00}

### (t) = d

### dt z

^{0}

### (t) = λ

^{2}

### y

^{00}

### (x) = λ

^{2}

### f (x, y(x)) = λ

^{2}

### f (a + λt, z(t)).

Likewise the boundary conditions are changed to

### z(0) = y(a) = α

^{and}

### z(1) = y(b) = β.

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{8}

With all these together, the problem (2) is transformed into

###

###

###

### z

^{00}

### (t) = λ

^{2}

### f (a + λt, z(t))

### z(0) = α, z(1) = β

^{(3)}

Thus, if

### y(x)

is a solution for (2), then### z(t) = y(a + λt)

is a solution for the boundary-value problem (3). Conversely, if### z(t)

is a solution for (3), then### y(x) = z(

_{λ}

^{1}

### (x − a))

is a solution for (2).**Example 1 Simplify the boundary conditions of the following equation by use of changing***variables.*

###

###

###

### y

^{00}

### = sin(xy) + y

^{2}

### y(1) = 3, y(4) = 7

*Solution: In this problem*

### a = 1, b = 4

^{, hence}

### λ = 3

. Now define the new variable### t =

^{1}

_{3}

### (x − 1)

^{, hence}

### x = 1 + 3t

^{, and let}

### z(t) = y(x) = y(1 + 3t)

^{. Then}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{9}

### λ

^{2}

### f (a + λt, z) = 9 sin(1 + 3t)z + z

^{2}

### ,

and the original equation is reduced to

###

###

###

### z

^{00}

### (t) = 9 sin((1 + 3t)z) + 9z

^{2}

### z(0) = 3, z(1) = 7

To reduce a two-point boundary-value problem

###

###

###

### z

^{00}

### (t) = g(t, z)

### z(0) = α, z(1) = β

to a homogeneous system, let

### u(t) = z(t) − [α + (β − α)t]

then

### u

^{00}

### (t) = z

^{00}

### (t)

^{, and}

### u(0) = z(0) − α = 0

^{and}

### u(1) = z(1) − β = 0

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{10}

Moreover,

### g(t, z) = g(t, u + α + (β − α)t) ≡ h(t, u).

The system is now transformed into

###

###

###

### u

^{00}

### (t) = h(t, u)

### u(0) = 0, u(1) = 0

**Example 2 Reduce the system**

###

###

###

### z

^{00}

### = [5z − 10t + 35 + sin(3z − 6t + 21)]e

^{t}

### z(0) = −7, z(1) = −5

*to a homogeneous problem by linear transformation technique.*

*Solution: Let*

### u(t) = z(t) − [−7 + (−5 + 7)t] = z(t) − 2t + 7.

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{11}

Then

### z(t) = u(t) + 2t − 7

^{, and}

### u

^{00}

### = z

^{00}

### = [5z − 10t + 35 + sin(3z − 6t + 21)]e

^{t}

### = [5(u + 2t − 7) − 10t + 35 + sin(3(u + 2t − 7) − 6t + 21)]e

^{t}

### = [5u + sin(3u)]e

^{t}

The system is transformed to

###

###

###

### u

^{00}

### (t) = [5u + sin(3u)]e

^{t}

### u(0) = u(1) = 0

**Example 3 Reduce the following two-point boundary-value problem**

###

###

###

### y

^{00}

### = y

^{2}

### + 3 − x

^{2}

### + xy y(3) = 7, y(5) = 9

*to a homogeneous system.*

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{12}

*Solution: In the original system,*

### a = 3, b = 5, α = 7, β = 9

^{. Let}

### λ = b − a = 2

^{and}

define a new variable

### t = 1

### 2 (x − 3) =⇒ x = 2t + 3.

Let the function

### z(t) = y(x) = y(2t + 3)

^{. Then}

### z

^{00}

### (t) = λ

^{2}

### y

^{00}

### (2t + 3) = λ

^{2}

### f (2t + 3, u)

### = 4[z

^{2}

### + 3 − (2t + 3)

^{2}

### + (2t + 3)z]

### = 4[z

^{2}

### + 3z + 2tz − 4t

^{2}

### − 12t − 6]

The original problem is first transformed into

###

###

###

### z

^{00}

### (t) = 4[z

^{2}

### + 3z + 2tz − 4t

^{2}

### − 12t − 6]

### z(0) = 7, z(1) = 9

Next let

### u(t) = z(t) − [7 + 2t],

or equivalently,### z(t) = u(t) + 2t + 7.

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{13}

Then

### u

^{00}

### (t) = 4[(z + 2t + 7)

^{2}

### + 3(u + 2t + 7) + 2t(u + 2t + 7) − 4t

^{2}

### − 12t − 6]

### = 4[u

^{2}

### + 6tu + 17u + 4t

^{2}

### + 36t + 64].

The original problem is transformed into the following homogeneous system

###

###

###

### u

^{00}

### (t) = 4[u

^{2}

### + 6tu + 17u + 4t

^{2}

### + 36t + 64]

### u(0) = u(1) = 0

**Theorem 2 The boundary-value problem**

###

###

###

### y

^{00}

### = f (x, y)

### y(0) = 0, y(1) = 0

*has a unique solution if*

### ∂f

### ∂y

*is continuous, non-negative, and bounded in the strip*

### 0 ≤ x ≤ 1

^{and}### −∞ < y < ∞

^{.}**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{14}

**Theorem 3 If**

### f

*is a continuous function of*

### (s, t)

*in the domain*

### 0 ≤ s ≤ 1

^{and}### −∞ < t < ∞

^{such that}### |f (s, t

_{1}

### ) − f (s, t

_{2}

### )| ≤ K|t

_{1}

### − t

_{2}

### |, (K < 8).

*Then the boundary-value problem*

###

###

###

### y

^{00}

### = f (x, y)

### y(0) = 0, y(1) = 0

*has a unique solution in*

### C[0, 1]

^{.}**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{15}

**2 – Finite Difference Method For Linear Problems**

We consider finite difference method for solving the linear two-point boundary-value problem of the form

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x)

### y(a) = α, y(b) = β.

^{(4)}

Methods involving finite differences for solving boundary-value problems replace each of the derivatives in the differential equation by an appropriate difference-quotient approximation.

**2.1 – The Finite Difference Formulation**

First, partition the interval

### [a, b]

^{into}

### n

equally-spaced subintervals by points### a = x

_{0}

### < x

_{1}

### < . . . < x

n### < x

n### = b

. Each mesh point### x

i can be computed by### x

i### = a + i ∗ h, i = 0, 1, . . . , n,

^{with}

### h = b − a n

where

### h

is called the mesh size.**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{15}

**2 – Finite Difference Method For Linear Problems**

We consider finite difference method for solving the linear two-point boundary-value problem of the form

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x)

### y(a) = α, y(b) = β.

^{(4)}

Methods involving finite differences for solving boundary-value problems replace each of the derivatives in the differential equation by an appropriate difference-quotient approximation.

**2.1 – The Finite Difference Formulation**

First, partition the interval

### [a, b]

^{into}

### n

equally-spaced subintervals by points### a = x

_{0}

### < x

_{1}

### < . . . < x

n### < x

n### = b

. Each mesh point### x

i can be computed by### x

i### = a + i ∗ h, i = 0, 1, . . . , n,

^{with}

### h = b − a n

where

### h

is called the mesh size.**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{15}

**2 – Finite Difference Method For Linear Problems**

We consider finite difference method for solving the linear two-point boundary-value problem of the form

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x)

### y(a) = α, y(b) = β.

^{(4)}

Methods involving finite differences for solving boundary-value problems replace each of the derivatives in the differential equation by an appropriate difference-quotient approximation.

**2.1 – The Finite Difference Formulation**

First, partition the interval

### [a, b]

^{into}

### n

equally-spaced subintervals by points### a = x

_{0}

### < x

_{1}

### < . . . < x

n### < x

n### = b

. Each mesh point### x

i can be computed by### x

i### = a + i ∗ h, i = 0, 1, . . . , n,

^{with}

### h = b − a n

where

### h

is called the mesh size.**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{15}

**2 – Finite Difference Method For Linear Problems**

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x)

### y(a) = α, y(b) = β.

^{(4)}

**2.1 – The Finite Difference Formulation**

First, partition the interval

### [a, b]

^{into}

### n

equally-spaced subintervals by points### a = x

_{0}

### < x

_{1}

### < . . . < x

n### < x

n### = b

^{.}

Each mesh point

### x

i can be computed by### x

i### = a + i ∗ h, i = 0, 1, . . . , n,

^{with}

### h = b − a n

where

### h

is called the mesh size.**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{15}

**2 – Finite Difference Method For Linear Problems**

###

###

###

### y

^{00}

### = p(x)y

^{0}

### + q(x)y + r(x)

### y(a) = α, y(b) = β.

^{(4)}

**2.1 – The Finite Difference Formulation**

First, partition the interval

### [a, b]

^{into}

### n

equally-spaced subintervals by points### a = x

_{0}

### < x

_{1}

### < . . . < x

n### < x

n### = b

. Each mesh point### x

i can be computed by### x

i### = a + i ∗ h, i = 0, 1, . . . , n,

^{with}

### h = b − a n

where

### h

is called the mesh size.**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{16}

At the interior mesh points,

### x

i, for### i = 1, 2, . . . , n − 1

, the differential equation to be approximated satisfies### y

^{00}

### (x

i### ) = p(x

i### )y

^{0}

### (x

i### ) + q(x

i### )y(x

i### ) + r(x

i### ).

^{(5)}

The central finite difference formulae

### y

^{0}

### (x

i### ) = y(x

i+1### ) − y(x

i−1### )

### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### ),

^{(6)}

for some

### η

i in the interval### (x

i−1### , x

i+1### )

^{, and}

### y

^{00}

### (x

i### ) = y(x

i+1### ) − 2y(x

i### ) + y(x

i−1### )

### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ),

^{(7)}

for some

### ξ

i in the interval### (x

i−1### , x

i+1### )

, can be derived from Taylor’s theorem by expanding### y

^{about}

### x

i.Let

### u

i denote the approximate value of### y

i### = y(x

i### )

^{. If}

### y ∈ C

^{4}

### [a, b]

, then a finite difference method with truncation error of order### O(h

^{2}

### )

can be obtained by using the approximations**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{16}

At the interior mesh points,

### x

i, for### i = 1, 2, . . . , n − 1

, the differential equation to be approximated satisfies### y

^{00}

### (x

i### ) = p(x

i### )y

^{0}

### (x

i### ) + q(x

i### )y(x

i### ) + r(x

i### ).

^{(5)}

The central finite difference formulae

### y

^{0}

### (x

i### ) = y(x

i+1### ) − y(x

i−1### )

### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### ),

^{(6)}

for some

### η

i in the interval### (x

i−1### , x

i+1### )

^{,}

and

### y

^{00}

### (x

i### ) = y(x

i+1### ) − 2y(x

i### ) + y(x

i−1### )

### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ),

^{(7)}

for some

### ξ

i in the interval### (x

i−1### , x

i+1### )

, can be derived from Taylor’s theorem by expanding### y

^{about}

### x

i.Let

### u

i denote the approximate value of### y

i### = y(x

i### )

^{. If}

### y ∈ C

^{4}

### [a, b]

, then a finite difference method with truncation error of order### O(h

^{2}

### )

can be obtained by using the approximations**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{16}

At the interior mesh points,

### x

i, for### i = 1, 2, . . . , n − 1

, the differential equation to be approximated satisfies### y

^{00}

### (x

i### ) = p(x

i### )y

^{0}

### (x

i### ) + q(x

i### )y(x

i### ) + r(x

i### ).

^{(5)}

The central finite difference formulae

### y

^{0}

### (x

i### ) = y(x

i+1### ) − y(x

i−1### )

### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### ),

^{(6)}

for some

### η

i in the interval### (x

i−1### , x

i+1### )

^{, and}

### y

^{00}

### (x

i### ) = y(x

i+1### ) − 2y(x

i### ) + y(x

i−1### )

### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ),

^{(7)}

for some

### ξ

i in the interval### (x

i−1### , x

i+1### )

, can be derived from Taylor’s theorem by expanding### y

^{about}

### x

i.Let

### u

i denote the approximate value of### y

i### = y(x

i### )

^{. If}

### y ∈ C

^{4}

### [a, b]

, then a finite difference method with truncation error of order### O(h

^{2}

### )

can be obtained by using the approximations**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{16}

At the interior mesh points,

### x

i, for### i = 1, 2, . . . , n − 1

, the differential equation to be approximated satisfies### y

^{00}

### (x

i### ) = p(x

i### )y

^{0}

### (x

i### ) + q(x

i### )y(x

i### ) + r(x

i### ).

^{(5)}

The central finite difference formulae

### y

^{0}

### (x

i### ) = y(x

i+1### ) − y(x

i−1### )

### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### ),

^{(6)}

for some

### η

i in the interval### (x

i−1### , x

i+1### )

^{, and}

### y

^{00}

### (x

i### ) = y(x

i+1### ) − 2y(x

i### ) + y(x

i−1### )

### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ),

^{(7)}

for some

### ξ

i in the interval### (x

i−1### , x

i+1### )

, can be derived from Taylor’s theorem by expanding### y

^{about}

### x

i.Let

### u

i denote the approximate value of### y

i### = y(x

i### )

^{.}

If

### y ∈ C

^{4}

### [a, b]

, then a finite difference method with truncation error of order### O(h

^{2}

### )

can be obtained by using the approximations**Department of Mathematics – NTNU**

**Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{16}

At the interior mesh points,

### x

i, for### i = 1, 2, . . . , n − 1

, the differential equation to be approximated satisfies### y

^{00}

### (x

i### ) = p(x

i### )y

^{0}

### (x

i### ) + q(x

i### )y(x

i### ) + r(x

i### ).

^{(5)}

The central finite difference formulae

### y

^{0}

### (x

i### ) = y(x

i+1### ) − y(x

i−1### )

### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### ),

^{(6)}

for some

### η

i in the interval### (x

i−1### , x

i+1### )

^{, and}

### y

^{00}

### (x

i### ) = y(x

i+1### ) − 2y(x

i### ) + y(x

i−1### )

### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ),

^{(7)}

for some

### ξ

i in the interval### (x

i−1### , x

i+1### )

, can be derived from Taylor’s theorem by expanding### y

^{about}

### x

i.Let

### u

i denote the approximate value of### y

i### = y(x

i### )

^{. If}

### y ∈ C

^{4}

### [a, b]

, then a finite difference method with truncation error of order### O(h

^{2}

### )

can be obtained by using the approximations**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively.Furthermore, let

### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

. Equation (8) can be written in theform

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{. In (8),}

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

linear equations to be solved. The resulting system of linear equations can be expressed in the matrix form

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively. Furthermore, let### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

. Equation (8) can be written in theform

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{. In (8),}

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

linear equations to be solved. The resulting system of linear equations can be expressed in the matrix form

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively. Furthermore, let### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

^{.}

Equation (8) can be written in the

form

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{. In (8),}

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

linear equations to be solved. The resulting system of linear equations can be expressed in the matrix form

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively. Furthermore, let### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

. Equation (8) can be written in theform

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{.}

In (8),

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively. Furthermore, let### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

. Equation (8) can be written in theform

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{. In (8),}

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

linear equations to be solved.

The resulting system of linear equations can be expressed in the matrix form

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{17}

### y

^{0}

### (x

i### ) ≈ u

i+1### − u

i−1### 2h

^{and}

### y

^{00}

### (x

i### ) ≈ u

i+1### − 2u

i### + u

i−1### h

^{2}

for

### y

^{0}

### (x

i### )

^{and}

### y

^{00}

### (x

i### )

, respectively. Furthermore, let### p

i### = p(x

i### ), q

i### = q(x

i### ), r

i### = r(x

i### ).

The discrete version of equation (4) is then

### u

i+1### − 2u

i### + u

i−1### h

^{2}

### = p

i### u

i+1### − u

i−1### 2h + q

i### u

i### + r

i### , i = 1, 2, . . . , n − 1,

^{(8)}

together with boundary conditions

### u

_{0}

### = α

^{and}

### u

n### = β

. Equation (8) can be written in theform

### 1 + h 2 p

i### u

i−1### − 2 + h

^{2}

### q

i### u

i### +

### 1 − h 2 p

i### u

i+1### = h

^{2}

### r

i### ,

^{(9)}

for

### i = 1, 2, . . . , n − 1

^{. In (8),}

### u

_{1}

### , u

_{2}

### , . . . , u

n−1 are the unknown, and there are### n − 1

### Au = f,

^{(10)}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{18}

where A =

−2 − h^{2}q_{1} 1 −

h
2 p_{1}
1 + ^{h}

2p_{2} −2 − h^{2}q_{2} 1 −

h
2p_{2}

. .. . .. . ..

. .. . .. . ..

1 + ^{h}

2p^{n}_{−2} −2 − h^{2}q^{n}_{−2} 1 −

h

2p^{n}_{−2}
1 + ^{h}

2p^{n}_{−1} −2 − h^{2}q^{n}_{−1}

,

### u =

###

###

###

###

###

###

###

###

###

###

### u

_{1}

### u

_{2}

.. .

### u

n−2### u

n−1###

###

###

###

###

###

###

###

###

###

### ,

^{and}

### f =

###

###

###

###

###

###

###

###

###

###

### h

^{2}

### r

_{1}

### − 1 +

^{h}

_{2}

### p

_{1}

### α h

^{2}

### r

_{2}

.. .

### h

^{2}

### r

n−2### h

^{2}

### r

n−1### − 1 −

^{h}

_{2}

### p

n−1### β

###

###

###

###

###

###

###

###

###

###

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{19}

Since the matrix

### A

is tridiagonal, this system can be solved by a special Gaussian elimination in### O(n

^{2}

### )

^{flops.}

**Theorem 4 Suppose that**

### p(x), q(x)

^{, and}### r(x)

*in (4) are continuous on*

### [a, b]

^{, and}### q(x) > 0

^{on}### [a, b]

*. Then (10) has a unique solution provided that*

### h < 2/L

^{, where}### L = max

_{a≤x≤b}

### |p(x)|

^{.}**2.2 – Convergence Analysis**

We shall analyze that when

### h

converges to zero, the solution### u

i of the discrete problem (8) converges to the solution### y

i of the original continuous problem (5).### y

i satisfies the following system of equations### y

i+1### − 2y

i### + y

i−1### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ) = p

i### y

i+1### − y

i−1### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### )

### + q

i### y

i### + r

i### ,

(11) for

### i = 1, 2, . . . , n − 1

. Subtract (8) from (11) and let### e

i### = y

i### − u

i, the result is### e

_{i+1}

### − 2e

i### + e

_{i−1}

### h

^{2}

### = p

i### e

_{i+1}

### − e

_{i−1}

### 2h + q

i### e

i### + h

^{2}

### g

i### , i = 1, 2, . . . , n − 1,

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{19}

Since the matrix

### A

is tridiagonal, this system can be solved by a special Gaussian elimination in### O(n

^{2}

### )

^{flops.}

**Theorem 4 Suppose that**

### p(x), q(x)

^{, and}### r(x)

*in (4) are continuous on*

### [a, b]

^{, and}### q(x) > 0

^{on}### [a, b]

*. Then (10) has a unique solution provided that*

### h < 2/L

^{, where}### L = max

_{a≤x≤b}

### |p(x)|

^{.}**2.2 – Convergence Analysis**

We shall analyze that when

### h

converges to zero, the solution### u

i of the discrete problem (8) converges to the solution### y

i of the original continuous problem (5).### y

i satisfies the following system of equations### y

i+1### − 2y

i### + y

i−1### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ) = p

i### y

i+1### − y

i−1### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### )

### + q

i### y

i### + r

i### ,

(11) for

### i = 1, 2, . . . , n − 1

. Subtract (8) from (11) and let### e

i### = y

i### − u

i, the result is### e

_{i+1}

### − 2e

i### + e

_{i−1}

### h

^{2}

### = p

i### e

_{i+1}

### − e

_{i−1}

### 2h + q

i### e

i### + h

^{2}

### g

i### , i = 1, 2, . . . , n − 1,

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{19}

Since the matrix

### A

is tridiagonal, this system can be solved by a special Gaussian elimination in### O(n

^{2}

### )

^{flops.}

**Theorem 4 Suppose that**

### p(x), q(x)

^{, and}### r(x)

*in (4) are continuous on*

### [a, b]

^{, and}### q(x) > 0

^{on}### [a, b]

*. Then (10) has a unique solution provided that*

### h < 2/L

^{, where}### L = max

_{a≤x≤b}

### |p(x)|

^{.}**2.2 – Convergence Analysis**

We shall analyze that when

### h

converges to zero, the solution### u

i of the discrete problem (8) converges to the solution### y

i of the original continuous problem (5).### y

i satisfies the following system of equations### y

i+1### − 2y

i### + y

i−1### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ) = p

i### y

i+1### − y

i−1### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### )

### + q

i### y

i### + r

i### ,

(11) for

### i = 1, 2, . . . , n − 1

. Subtract (8) from (11) and let### e

i### = y

i### − u

i, the result is### e

_{i+1}

### − 2e

i### + e

_{i−1}

### h

^{2}

### = p

i### e

_{i+1}

### − e

_{i−1}

### 2h + q

i### e

i### + h

^{2}

### g

i### , i = 1, 2, . . . , n − 1,

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{19}

Since the matrix

### A

is tridiagonal, this system can be solved by a special Gaussian elimination in### O(n

^{2}

### )

^{flops.}

**Theorem 4 Suppose that**

### p(x), q(x)

^{, and}### r(x)

*in (4) are continuous on*

### [a, b]

^{, and}### q(x) > 0

^{on}### [a, b]

*. Then (10) has a unique solution provided that*

### h < 2/L

^{, where}### L = max

_{a≤x≤b}

### |p(x)|

^{.}**2.2 – Convergence Analysis**

We shall analyze that when

### h

converges to zero, the solution### u

i of the discrete problem (8) converges to the solution### y

i of the original continuous problem (5).### y

i satisfies the following system of equations### y

i+1### − 2y

i### + y

i−1### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ) = p

i### y

i+1### − y

i−1### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### )

### + q

i### y

i### + r

i### ,

(11) for

### i = 1, 2, . . . , n − 1

^{.}

Subtract (8) from (11) and let

### e

i### = y

i### − u

i, the result is### e

_{i+1}

### − 2e

i### + e

_{i−1}

### h

^{2}

### = p

i### e

_{i+1}

### − e

_{i−1}

### 2h + q

i### e

i### + h

^{2}

### g

i### , i = 1, 2, . . . , n − 1,

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{19}

Since the matrix

### A

is tridiagonal, this system can be solved by a special Gaussian elimination in### O(n

^{2}

### )

^{flops.}

**Theorem 4 Suppose that**

### p(x), q(x)

^{, and}### r(x)

*in (4) are continuous on*

### [a, b]

^{, and}### q(x) > 0

^{on}### [a, b]

*. Then (10) has a unique solution provided that*

### h < 2/L

^{, where}### L = max

_{a≤x≤b}

### |p(x)|

^{.}**2.2 – Convergence Analysis**

We shall analyze that when

### h

converges to zero, the solution### u

i of the discrete problem (8) converges to the solution### y

i of the original continuous problem (5).### y

i satisfies the following system of equations### y

i+1### − 2y

i### + y

i−1### h

^{2}

### − h

^{2}

### 12 y

^{(4)}

### (ξ

i### ) = p

i### y

i+1### − y

i−1### 2h − h

^{2}

### 6 y

^{(3)}

### (η

i### )

### + q

i### y

i### + r

i### ,

(11) for

### i = 1, 2, . . . , n − 1

. Subtract (8) from (11) and let### e

i### = y

i### − u

i, the result is### e

_{i+1}

### − 2e

i### + e

_{i−1}

### h

^{2}

### = p

i### e

_{i+1}

### − e

_{i−1}

### 2h + q

i### e

i### + h

^{2}

### g

i### , i = 1, 2, . . . , n − 1,

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{20}

where

### g

i### = 1

### 12 y

^{(4)}

### (ξ

i### ) − 1

### 6 p

i### y

^{(3)}

### (η

i### ).

After collecting terms and multiplying by

### h

^{2}

^{, we have}

### 1 + h 2 p

i### e

i−1### − 2 + h

^{2}

### q

i### e

i### +

### 1 − h 2 p

i### e

i+1### = h

^{4}

### g

i### , i = 1, 2, . . . , n − 1.

Let

### e = [e

_{1}

### , e

_{2}

### , . . . , e

_{n−1}

### ]

^{T}

^{and}

### |e

k### | = kek

_{∞}

^{. Then}

### 2 + h

^{2}

### q

k### e

k### =

### 1 + h 2 p

k### e

k−1### +

### 1 − h 2 p

k### e

k+1### − h

^{4}

### g

k### ,

and, hence

### 2 + h

^{2}

### q

k### |e

k### | ≤

### 1 + h 2 p

k### |e

k−1### | +

### 1 − h 2 p

k### |e

k+1### | + h

^{4}

### |g

k### |

### ≤

### 1 + h 2 p

k### kek

_{∞}

### +

### 1 − h 2 p

k### kek

_{∞}

### + h

^{4}

### kgk

_{∞}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{20}

where

### g

i### = 1

### 12 y

^{(4)}

### (ξ

i### ) − 1

### 6 p

i### y

^{(3)}

### (η

i### ).

After collecting terms and multiplying by

### h

^{2}

^{, we have}

### 1 + h 2 p

i### e

i−1### − 2 + h

^{2}

### q

i### e

i### +

### 1 − h 2 p

i### e

i+1### = h

^{4}

### g

i### , i = 1, 2, . . . , n − 1.

Let

### e = [e

_{1}

### , e

_{2}

### , . . . , e

_{n−1}

### ]

^{T}

^{and}

### |e

k### | = kek

_{∞}

^{. Then}

### 2 + h

^{2}

### q

k### e

k### =

### 1 + h 2 p

k### e

k−1### +

### 1 − h 2 p

k### e

k+1### − h

^{4}

### g

k### ,

and, hence

### 2 + h

^{2}

### q

k### |e

k### | ≤

### 1 + h 2 p

k### |e

k−1### | +

### 1 − h 2 p

k### |e

k+1### | + h

^{4}

### |g

k### |

### ≤

### 1 + h 2 p

k### kek

_{∞}

### +

### 1 − h 2 p

k### kek

_{∞}

### + h

^{4}

### kgk

_{∞}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**

**BVP of ODE**

^{20}

where

### g

i### = 1

### 12 y

^{(4)}

### (ξ

i### ) − 1

### 6 p

i### y

^{(3)}

### (η

i### ).

After collecting terms and multiplying by

### h

^{2}

^{, we have}

### 1 + h 2 p

i### e

i−1### − 2 + h

^{2}

### q

i### e

i### +

### 1 − h 2 p

i### e

i+1### = h

^{4}

### g

i### , i = 1, 2, . . . , n − 1.

Let

### e = [e

_{1}

### , e

_{2}

### , . . . , e

_{n−1}

### ]

^{T}

^{and}

### |e

k### | = kek

_{∞}

^{.}

Then

### 2 + h

^{2}

### q

k### e

k### =

### 1 + h 2 p

k### e

k−1### +

### 1 − h 2 p

k### e

k+1### − h

^{4}

### g

k### ,

and, hence

### 2 + h

^{2}

### q

k### |e

k### | ≤

### 1 + h 2 p

k### |e

k−1### | +

### 1 − h 2 p

k### |e

k+1### | + h

^{4}

### |g

k### |

### ≤

### 1 + h 2 p

k### kek

_{∞}

### +

### 1 − h 2 p

k### kek

_{∞}

### + h

^{4}

### kgk

_{∞}

**Department of Mathematics – NTNU** **Tsung-Min Hwang December 20, 2003**