This preview shows page 1. Sign up to view the full content.
Unformatted text preview: and consistency
implying convergence. However, in contrast, the local and global discretization errors
converge at the same rates. 8.2 The Box Scheme for FirstOrder Systems
We'll extend nite di erence methods to nonlinear rstorder BVPs of the form y = f (x y) a<x<b (8.2.1a) g (y(a)) = g (y(b)) = 0: (8.2.1b) 0 L R As in Chapter 7, y and f are mvectors, g is an lvector, and g is an rvector with
l + r = m.
We will discretize (8.2.1) by a variant of the trapezoidal rule called the \box scheme"
that was proposed by Keller 6] and developed by Lentini and Pereyra 7]. Towards this
end, introduce a partition of (a b) into N subintervals as shown in Figure 8.2.1. Let
L 3 R y yi1
yi y(x) hi
x
a = x0 x1 x i1 xi x N= b Figure 8.2.1: Mesh nomenclature for the box scheme. h =x ;x
i i (8.2.2) 1 i; and approximate (8.2.1a) using the trapezoidal rule
(y) = y ; y 1 ; f (x y ) + f (x 1 y 1) = 0
h
2
where
2
3
y0
6y 7
y = 6 ..1 7 :
6
4.7
5
i i; i i i; i; i i i = 1 2 ::: N (8.2.3a) (8.2.3b) y N The boundary conditions (8.2.1b) complete the speci cation of the discrete problem as g (y0 ) = g (y ) = 0:
L R (8.2.3c) N The BVP has been approximated by a nonlinear algebraic system (8.2.3) of dimension
mN for the N vector of unknowns y , i = 0 1 : : : N .
The midpoint rule can be used in place of the trapezoidal rule to obtain
(8.2.4)
(y) = y ; y 1 ; f (x 1 2 y + y 1 ) = 0:
h
2
The midpoint and trapezoidal rules have similar, but not identical, accuracy and stability
properties.
i i i i; i; = i i 4 i; Focusing on the trapezoidal rule, we'll linearize (8.2.3) by Newton's method. For
(8.2.3a) we have
@ (y( ) ) y( ) = ; (y( ) )
i = 1 2 ::: N
(8.2.5a)
@y
and for (8.2.3b) we have
i i (
@ g (y0 )) y( ) = ;g (y( ) )
0
0
@ y0 (8.2.5b) @ g (y( )) y( ) = ;g (y( ))
@y (8.2.5c) y( +1) = y( ) + y( ) : (8.2.5d) L L R N R N N N with Using (8.2.3a), the two nonzero contributions to the Jacobian in (8.2.5a) are @ (y( ) ) := L( ) = ; I ; fy (x 1 y( )1)
@y 1
h
2
i; i (8.2.6a) i; i i; i and @ (y( ) ) := R( ) = I ; fy (x y( ))
(8.2.6b)
@y
h
2
where I is the m m identity matrix.
For consistency and simplicity, we'll also write the Jacobians in (8.2.5b,c) as
i i i i i i L( ) :=
0 (
@ gL(y0 ) )
@ y0 (8.2.6c) (
@ gR (yN ) )
@ yN
Collecting (8.2.5) and (8.2.6), we nd the linear Newton system as
)
R(N+1 := 2
6
6
6
6
6
6
4 3
2
7
76
76
... ...
76
7
( ) R( ) 7 4
L
5 L(0 )
L(1 ) R(1 ) N )
R( +1
N 2
(
3
g (y0 ) )
6 (y( ) )
7
7 = ;6 1 .
6
..
6
... 7
6
5
4 (y( ) )
()
y
g (y ( ) )
(
y0 )
(
y1 ) N N N 5 L R N (8.2.6d)
3
7
7
7:
7
7
5 (8.2.7) We'll discuss an O(N ) algorithm for solving this block bidiagonal system in Section 8.5.
For the moment, we'll consider an application.
Example 8.2.1 ( 1], Section 5.1). Consider the secondorder linear BVP y ; 2y = 2 cos2 x + 2 2 cos 2 x 0<x<1 00 y(0) = y(1) = 0: We'll write this as a rstorder system by letting y1 = y and y2 = y= to get y1
y2 0 y1 +
y2 = 00 cos2 02
x + 2 cos 2 x : The exact solution of this problem is
" ; (1;x) + ; x
#
; cos2 x :
;
1+
y(x) = ; (1;x) ; x + sin 2 x
1+ ;
e e e e ;e e This problem is similar to one that was addressed in Section 7.1. We recall that when
1 there are boundary layers where the solution varies rapidly near both x = 0 and
x = 1. The solution varies sinusoidally in the interior (Figure 8.2.2). He...
View
Full
Document
This document was uploaded on 03/16/2014 for the course CSCI 6820 at Rensselaer Polytechnic Institute.
 Spring '14
 JosephE.Flaherty

Click to edit the document details