Unformatted text preview: ...
... ... ...
(6.3.26)
6
74 . 5 6
7
..
6
7
6
7
2
4
4
bN ;1 aN ;1
cN ;1 5 y
h rN ;1 5
N
2
2 aN ; 2hCcN
h rN ; 2hcN B
which may be solved by the tridiagonal algorithm of Figure 6.3.2.
Now let us return to the original nonlinear problem (6.3.12). Most iterative schemes
for solving nonlinear algebraic equations can be used to determine the solution, but we'll
illustrate the use of Newton's method, which is the most popular. To begin, let us write
the nite di erence system (6.3.12a) in the form
;
Fi(y) = yi;1 ; 2yi + yi+1 ; h2 f (xi yi yi+1 2h yi;1 ) = 0
i = 1 2 ::: N ; 1
(6.3.27)
subject to the Dirichlet boundary conditions (6.3.12b) and the de nition of the vector of
unknowns y given by (6.3.18c).
Newton's iteration involves solving Fy (y( ) )(y( +1) ; y( ) ) = ;F(y( ) )
where 2
F1(y)
6 F2 (y)
F(y) = 6 ..
6
4
.
FN ;1(y) 3
7
7
7
5 2
6
6
Fy (y) = 6
6
4 (6.3.28a) @F1
@y1
@F2
@y1 @F1
@y2
@F2
@y2 @ F1
@yN ;1
@ F2
@yN ;1 @FN ;1
@y1 @FN ;1
@y2 @ FN ;1
@yN ;1 Di erentiating (6.3.27), the Jacobian Fy (y) is
8
>
>
>
@Fi = <
@yj >
>
>
: = 0 1 ::: ... ... ... 3
7
7
7:
7
5 (6.3.28b) if j = i ; 1
if j = i
:
if j = i + 1
otherwise (6.3.29) y( ) ; y( )
@f
b(i ) = 1 + h @y0 (xi yi( ) i+1 2h i;1 )
2 (6.3.30a) @f
1 + h @y0
2
;2 ; h2 @f
@y
@f
1 ; h @y0
2
0 Letting 23 )
@f (x y( ) yi(+1 ; yi(;)1 )
ai = ;2 ; h @y i i
2h (6.3.30b) y( ) ; y( )
@f
c(i ) = 1 ; h @y0 (xi yi( ) i+1 2h i;1 ):
2 (6.3.30c) () gives 2 2 () ()
ac
6 b(1 ) a1 ) c( )
(
62
2
2
6
()
... ... ...
6
Fy (y ) = 6
6
)
)
)
b(N ;2 a(N ;2 c(N ;2
4
()
() 3
7
7
7
7:
7
7
5 (6.3.30d) bN ;1 aN ;1
Each Newton iteration requires the solution of a tridiagonal system. The Jacobian of
this system need not be reevaluated and factored after each Newton step thus, only the
forward and backward substitution steps of the tridiagonal algorithm shown in Figure
6.3.2 need be performed at each iterative step. The derivatives @[email protected] and @[email protected] can
be approximated by nite di erences.
Convergence of Newton's method is typically quadratic except at a bifurcation point
where it is often linear. The use of nite di erence approximations in the Jacobian also
slows the convergence rate.
Example 6.3.3. Consider the elastica problem described in Section 6.1 and repeated
here using the notation of this Section as
y00 + P sin y = 0
Hence, and y(0) = y(1=2) = 0: f (x y y0) = ;P sin y
@f = ;P cos y
@f = 0
@y
@y0
b(i ) = c(i ) = 1 a(i ) = ;2 + h2 P cos yi( ): Using the convergence criteria that jjF(y( ))jj1 = 1 max;1 jFi(y( ) )j 10;9
iN
24 and setting P = 40, we found the number of Newton iterations and solution at x = 0:25
to be as recorded in Table 6.3.2. The number of Newton iterations decreases as the mesh
becomes ner. This is a result of the solution appearing to be smoother. The convergence
rate seems to be nearly quadratic. h Ki
yi
0.1 5 5 0.41240948
0.05 4 10 0.34795042
0.025 3 20 0.32985549
Table 6.3.2: Number of iterations K to reach convergence and the approximate solution
at x = 0:25 for Example 6.3.3.
The development and description of nite di erence equations may be simpli ed by
introducing a set of nite di erence operators as shown in Table 6.3.3.
The next several examples illustrate some applications of these nite di erence operators.
Example 6.3.4. The centered di erence formula (6.3.7) can be expressed in terms of
the central di erence and averaging operators and as
yi = (yi+1=2 ; yi;1=2) = yi+1 ; yi;1 :
h
h
2h
Example 6.3.5. An operator raised to a positive integer power is iterated, e.g.,
2 yi = (yi+1=2 ; yi;1=2) = yi+1 ; 2yi + yi;1: Thus, the centered second di erence approximation (6.3.8) of the second derivative can
be written as
2
yi00 = hyi :
2
Example 6.3.6. Expanding y (xi+1) in a Taylor's series about xi yields
2
3
y(xi+1) = y(xi) + hy0(xi ) + h y00(xi) + h y000 (xi) + : : : :
2
6 Using the derivative operator D de ned in Table 6.3.3 y(xi+1) = 1 + hD + h D2 + : : : ]y(xi):
2
2 25 Operator
Forward Di erence Symbol Backward Di erence r De nition
yi := yi+1 ; yi ryi := yi ; yi;1 Central Di erence yi := yi+1=2 ; yi;1=2 Average yi := (yi+1=2 + yi;1=2 )=2 Shift E Eyi := yi+1
Dyi := yi0 Derivative
D
Table 6.3.3: De nition of nite di erence operators.
This suggests the shorthand operator notation Ey(xi) = y(xi+1) = ehD y(xi)
where E is the shift operator (Table 6.3.3). We, thus, infer the identity between the shift,
exponential, and derivative operators E = ehD : (6.3.31) Additional relationships can be obtained by noting that yi = (E ; 1)yi, which implies
that = E ; 1 or E = 1 + . Using this with (6.3.31) gives hD = ln E = ln(1 + ) = ; 1
2 2 +1
3 3 ; ::: (6.3.32a) where the series expansion of ln(1+ x), jxj < 1, has been used. A similar relation in terms
of the backward di erence operator can be constructed by noting that r = 1 ; E ;1 thus,
1
hD = ln E = ; ln(1 ; r) = r + 2 r2 + 1 r3 + : : :
3 (6.3.32b) These identities can be used to derive highorder nite di erence...
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