Atkinson K. An Introduction to Numerical Analysis
Подождите немного. Документ загружается.
_j
SYSTEMS
OF
NONLINEAR EQUATIONS 105
Let
Gn
denote the matrix in (2.10.8). Then we can rewrite this equation as
(2.10.9)
It
is convenient to introduce the Jacobian matrix for the functions g
1
and g
2
:
G(x) =
(2.10.10)
In
(2.10.9), if xn
is
close to a, then
Gn
will
be
close to
G(a).
This
will
make
the
size
or
norm
of
G(
a)
crucial in analyzing the convergence in (2.10.9). The matrix
G(
a)
plays the role of
g'(
a)
in the theory of Section 2.5. To measure the size
of
the errors a - xn and
of
the matrices
Gn
and G(
a),
we use the vector and matrix
norms
of
(1.1.16) and (1.1.19) in Chapter 1.
Theorem 2.9
Let D
be
a closed, bounded, and convex set in the plane. (We say
D
is
convex if for any two points
in
D,-the line segment joining
them
is
also in D.) Assume
that
the components of g(x) are
continuously differentiable at all points
of
D,
and
further assume
I.
g(D)
c D,
2.
A =
MaxiJG(x)IJoo
< 1
xED
Then
(a) x = g(x) has a unique solution a E D.
(2.10.11)
(2.10.12)
(b)
For
any initial point x
0
E D, the iteration (2.10.5) will
converge in D to a.
(2.10.13)
with
£n
~
0 as n
~
oo.
Proof
(a)
The
existence of a fixed point a can
be
shown by proving that the
sequence of iterates
{xn}
from (2.10.5) are convergent in D. We leave
that to a problem,
and
instead
just
show the uniqueness
of
a.
Suppose a
and
13
are both fixed points
of
g(x) in D. Then
a -13 =
g(a)-
g(l3)
(2.10.14)
106 ROOTFINDING FOR NONLINEAR EQUATIONS
Apply
the
mean value theorem to component
i,
obtaining
i = 1, 2 (2.10.15)
with
and
E(i)
ED,
on
theJine-segment
joining_(x__.and
~-
Since
JIG(x)lloo
~A.
< 1, we have-from-the-definition -of-the-norm that
XED,
i = 1, 2
Combining
this with (2.10.15),
Lsd
a)
-
!!:(~)I
~
"All«
-
~II.,.,
llg(a)-
g(f3)
lloo
~"All«-
~lloo
Combined
with (2.10.14), this yields
(2.10.16)
which is possible only if
« =
~.
showing the uniqueness
of
a.
in D.
(b)
Condition
(2.10.11) will ensure that all
xn
E D if x
0
E D. Next
subtract
xn+l
= g(xn) from
(l
= g(«), obtaining
The
result (2.10.16) applies to any two points in D. Applying this,
(2.10.17)
Inductively.
(2.10.lg)
Since A < 1, this shows X n
-+
(l
as n
-+
00.
(c)
From
(2.10.9)
and
using (1.1.21),
(2.10.19)
As n
-+
oo, the points
i;~'l
used in evaluating
Gn
will all tend to «, since
they
are
on
the line segment joining
xn
and
«.Then
IIGnJioo
-+
JIG(«)IIoo
as n
-+
oo. Result (2.1 0.13) follows from (2.1 0.19) by letting
(n =
IIGnlloo
-JIG(«)IIoo· •
SYSTEMS OF NONLINEAR EQUATIONS 107
The
preceding theorem is the generalization
to
two variables of Theorem 2.6
for functions
of
one variable. The following generalizes Theorem 2.7.
Corollary 2.10 Let a be a
fixed
point of g(x),
and
assume components of g(x)
are continuously differentiable in some neighborhood about
a.
Further assume
IIG(
a)lloo
< l
(2.10.20)
Then for
Xo
chosen sufficiently close to
a,
the iteration
xn+l
=
g(xn) will converge
to
a,
and the results of Theorem 2.9 will be
valid on some closed, bounded, convex region about a.
•
We leave the proof of this
as
a problem. Based
on
results in Chapter
7,
the
linear convergence of x
n to a
will
still be true if all eigenvalues of G
(a)
are less
·than one in magnitude, which can be shown to
be
a weaker assumption than
(2.10.20).
Example Continue the earlier example (2.10.7).
It
is straightforward to compute
and therefore
G
(-
) = [ .038920
a .008529
.000401]
-.006613
IIG(a)lloo
= 0393
Thus the condition (2.10.20) of the theorem
is
satisfied. From (2.10.13), it will be
approximately true that
for all sufficiently large
n.
Suppose that A
is
a constant nonsingular matrix of order 2 X
2.
We can then
reformulate
(2.10.1)
as
X = X + Af(x) = g(x)
(2.10.21)
I
The example (2.10.7) illustrates this procedure. To see the requirements on A,
we
produce the Jacobian matrix. Easily,
G(x) = I + AF(x)
where
F(x)
is
the Jacobian matrix of /
1
and /
2
,
F(x)
=
(2.10.22)
We
want to choose A so that (2.10.20) is satisfied.
And
for rapid convergence,
we
108 ROOTFINDING FOR NONLINEAR EQUATIONS
want IIG(
«)11
00
=
0,
or
_The matrix in (2.10.7)
was
chosen in this way using
This suggests using a continual updating of
A, say
A=
-F(xn)-
1
•
The resulting
method is
xn+l
=
Xn-
F(xn)
-lf(xn)
n~O
(2.10.23)
We consider this method in the next section.
2.11
Newton's Method for Nonlinear Systems
As with Newton's method for a single equation, there
is
more than one way of
viewing
and
deriving the Newton method for solving a system of nonlinear
equations. We begin with an analytic derivation, and then
we
give
a geometric
perspective.
Apply Taylor's theorem
for
functions of two variables to each of the equations
/;(x
1
,
x
2
)
=
0,
expanding/;(«) about x
0
:
fori=
1,
2
with
~(i)
on
the
line segment joining x
0
and
«.
If
we
drop the second-order
terms, we
obtain
the approximation
In
matrix form,
(2.11.3)
with
F(x
0
)
the Jacobian matrix of
J,
given in (2.10.22).
Solving for «,
NEWTON'S METHOD
FOR
NONLINEAR SYSTEMS 109
The
approximation x
1
should be an improvement
on
x
0
,
provided x
0
is
chosen
sufficiently close to
«. This leads to the iteration method first obtained at the end
of
the last section,
n::::::O
(2.11.4)
This
is
Newton's method for solving the nonlinear system f(x) =
0.
In actual practice,
we
do not invert F(xn), particularly for systems of more
than two equations. Instead
we
solve a linear system for a correction term to xn:
F(xJan+I
=
-f(xJ
(2.11.5)
This
is
more efficient in computation time, requiring only about one-third
as
many operations as inverting F(xn). See Sections
8.1
and
8.2 for a discussion of
the numerical solution of linear systems of equations.
There is a geometrical derivation for Newton's method, in analogy with the
tangent line approximation used with single nonlinear equations in Section
2.2.
The
graph in space of the equation
is a plane that
is
tangent to the graph of z =
/;(x
1
,
x
2
)
at
the point x
0
,
for
i = 1,
2.
If
x
0
is
near«,
then these tangent planes should be good approximations
to the associated surfaces of
z =
f;(x
1
,
x
2
),
for x =
(x
1
,
x
2
)
near «. Then the
intersection
of
the zero curves
of
the tangent planes z =
P;(x
1
,
x
2
)
should be a
good approximation to the corresponding intersection a
of
the zero curves of the
original surfaces
z = /;(x
1
,
x
2
).
This results in the statement (2.11.2). The inter-
section
of
the zero curves of z =
P;(x
1
,
x
2
),
i = 1, 2, is the point x
1
•
Example Consider the system
!
1
=
4x?
+
xi
- 4 = 0
There are only two roots, one near
(1,
0) and its reflection through the origin near
(
-1,
0).
Using (2.11.4) with x
0
= (1, 0),
we
obta:in the results given in Table 2.20.
Table 2.20
Example of Newton's method
n
Xr.n
x2.n
/r(Xn) /2(xn)
0
1.0 0.0 0.0
1.59E-
1
1
1.0
-.1029207154
1.06E-
2
4.55E-
3
2
.9986087598
-.1055307239
1.46E-
5
6.63E-
7
3
.9986069441
-
.1055304923
1.32E-
11
1.87E-
12
J
110
ROOTFINDING
FOR NONLINEAR EQUATIONS
Convergence
analysis
For
the convergence analysis
of
Newton's
method
(2.11.4),
regard
it
as
a fixed-point
iteration
method
with
g(x) =
x-
F(x) -
1
f(x)
(2.11.6)
Also
assume
Determinant
F(
a)
'4=
0
which is
the
analogue
of
assuming a is a simple
root
when
dealing with a single
equation,-as
in
Theorem
2.1.
It
can
then
be
shown
that
the
Jacobian
G(x)
of
(2.11.6) is
zero
at
x = a (see Problem 53); consequently,
the
condition
(2.10.20) is
easily satisfied.
Corollary
2.10
then
implies
that
x"
converges to
a,
provided x
0
is
chosen
sufficiently
close
to
a.
In
addition,
it
can
be
shown
that
the
iteration is quadratic.
Specifically.
the
formulas (2.11.1)
and
(2.11.4) can
be
combined
to
obtain
n
~
0 (2.11.7)
for
some
constant
B > 0.
Variations
of
Newton's
method
Newton's
method
has
both
advantages
and
disadvantages
when
compared
with
other
methods for solving systems
of
nonlin-
ear
equations.
Among
its advantages, it
is
very simple in form
and
there
is
great
flexibility
in
using
it
on
a large variety
of
problems.
If
we
do
not
want
to
bother
supplying
partial
derivatives
to
be
evaluated by a
computer
program, we can use
a difference
approximation.
For
example, we
commonly
use
(2.11.8)
with
some
very small
number
~:.
For
a
detaikd
discussion
of
the choice
of
~:,
see
Dennis
and
Schnabel
(1983, pp. 94--99).
The
first
disadvantage
of
Newton's
method
is
tliat there are
other
methods
which
are
(1) iess expensive to use,
and/or
(2) easier to use for
some
special
classes
of
problems.
For
a system of m nonlinear
equations
in m unknowns, each
iterate
for
Newton's
method requires m
2
+ m function evaluations in general. In
addition,
Newton's
method
requires the solution
of
a system
of
m linear
equations
for
each
iterate,
at
a cost of
about
fm
3
arithmetic
operations
per
linear
system.
There
are
other
methods
that
are as fast
or
almost as fast in their
mathematical
speed
of
convergence,
but
that
require fewer function evaluations
and
arithmetic
operations
per
iteration. These
are
often
referred to as Newton-like,
quasi-Newton,
and
modified Newton methods.
For
a general
presentation
of
many
of
these
methods,
see
Dennis
and
Schnabel (1983).
A
simple
modification
of
Newton's
method
is
to
fix
the
Jacobian
matrix for
several
steps,
say
k:
j=0,1,
...
,k-1
(2.11.9)
UNCONSTRAINED OPTIMIZATION
111
for r = 0, 1,
2,
....
This means the linear system in
j = 0, 1,
...
,
k-
1, (2.11.10)
can be solved much more efficiently than in the original Newton method (2.11.5).
The linear system, of order
m,
requires about fm
3
arithmetic operations for its
solution in the first case, when
j =
0.
But each subsequent case, j =
1,
...
, k -
1,
will require only 2m
2
aritl:imetic operations for its solution. See Section
8.1
for
more complete details. The speed of convergence of
(2.11.9) will be slower than
the original method
(2.11.4),
but
the actual computation time of the modified
method will often be much less. For a more detailed examination of this
question, see Petra and Ptak
(1984,
p.
119).
A second problem with Newton's method,
and
with many other methods,
is
that often x
0
must be reasonably close to «
in
order to obtain convergence.
There are modifications of Newton's method to force convergence for poor
choices of
x
0
•
For
example, define
(2-11.11)
and choose s > 0
to
minimize
m
llf(xn + sdn)
~~~
= L [.tj(xn + sdn)]
2
(2.11.12)
j=l
The
choices=
1 in (2.11.11) yields Newton's method,
but
it may not be the best
choice. In some cases,
s may need
to
be much smaller than 1,
at
least initially,
in
order to ensure convergence_ For a more detailed discussion, see Dennis and
Schnabel
(1983,
chap~
6)_
For an analysis of some current programs for solving nonlinear systems, see
Hiebert (1982). He also discusses the difficulties in producing such software.
2.12 Unconstrained Optimization
Optimization refers to finding the maximum
or
mm1mum of a continuous
function
f(x
1
,
...
,
xm)·
This is an extremely important problem, lying at the
heart of modern industrial engineering, management science, and other areas.
This section discusses some methods and perspectives for calculating the mini-
mum or maximum of a function
f(x
1
,.
_.,
xm)-
No
formal algorithms are given,
since this would require too extensive a development.
Vector notation
is
used
in
much of the presentation, to give results for a
general
numb~r
m of variables.
We
consider only the unconstrained optimization
problem,
in
which there are no limitations on
(x
1
,
•..
, xm).
For
simplicity only,
we
also assume
f(x
1
,
•••
,
xm)
is
defined for all
(x
1
,
•••
,
xm).
Because the behavior of a function
f(x)
can be quite varied, the problem must
be further limited. A point
«
is
called a strict local minimum
off
if
f(x)
>
/(
«)
112 ROOTFINDING FOR NONLINEAR EQUATIONS
for all x close to
a,
x * a.
We
limit ourselves to finding a strict local minimum of
f(x).
Generally an initial
guess
x
0
of a will
be
known, and
f(x)
will be assumed
to be twice continuously differentiable with respect to its variables x
1
,
•••
, xm.
Reformulation as a nonlinear system With the assumption of differentiability, a
necessary condition for
a
to
be a strict local minimum
is
that
aj(a)
--=0
a
xi
Thus the nonlinear system
aj(x)
--=0
ax;
i=1,
...
,m
(2.12.1)
i=l,
...
,m
(2.12.2)
can
be
solved, and each calculated solution can be checked
as
to whether it
is
a
local maximum, minimum, or neither. For notation, introduce the gradient vector
"V/(x)
=
af
Using this vector, the system (2.12.2)
is
written more compactly
as
"V/(x)=O
(2.12.3)
To
solve (2.12.3), Newton's method (2.11.4) can be used,
as
well
as
other
rootfinding methods for nonlinear systems. Using Newton's method leads
to
(2.12.4)
with H(x) the Hessian matrix of
J,
l~i,j:$;m
If
a is a strict local minimum of
J,
then Taylor's theorem (1.1.12) can be used to
show that
H(
«)
is
a nonsingular matrix; then H(x) will be nonsingular for x
close to
u.
For
convergence, the analysis of Newton's method in the preceding
section can
be
used to prove quadratic convergence of
xn
to a provided x
0
is
chosen sufficiently close to a.
The
main drawbacks with the iteration (2.12.4) are the same
as
those given in
the last section for Newton's method for solving nonlinear systems. There are
UNCONSTRAINED OPTIMIZATION 113
other, more efficient optimization methods that seek to approximate a
by
using
only
l(x)
and
\1
l(x). These methods may require more iterations, but generally
their total computing time will
be
much
less
than with Newton's method. In
addition, these methods seek to obtain convergence for a larger
·set of initial
values
x
0
•
Descent methods Suppose
we
are trying to minimize a function l(x). Most
methods for doing so are based
on
the following general two-step iteration
process.
STEP
Dl:
At
x,,
pick a direction
d,
such that
l(x)
will decrease
as
x moves
away from
x,
in the direction
d,.
STEP
D2:
Let x,+
1
=
x,
+
sd,,
with s chosen to minimize
q>(s) =
l(x,
+
sd,),
s
~
0
(2.12.5)
Usually s
is
chosen
as
the
small~st
positive relative minimum of
q>(
s
).
Such methods are called
descent=methodS'.-~With~each-iteration,
Descent methods are guaranteed to converge under more general conditions
than for Newton's method
(2.12.4). Consider
the
level surface
C = {xll(x) =
l(x
0
)}
and consider only the connected portion of it, say C', that contains x
0
•
Then if
C'
is
bounded and contains « in its interior, descent methods
will
converge under
very general conditions. This
is
illustrated for the two-variable case in Figure
2.8.
Several level curves
l(x
1
,
x
2
)
= c are shown for a set of values c approaching
I (a.). The vectors
d,
are directions in which I (
x)
is
decreasing.
There are a number of
ways.
for choosing the directions
d,,
and the best
known are as follows.
l.
The method
of
steepest descent. Here
d,
=
-\1
l(x,).
It
is
the direction in
which
l(x) decreases most rapidly when moving away from x,. It
is
a good
strategy near
x,,
but it ·usually turns out
to
be a poor strategy for rapid
convergence
to
«.
2.
QUilSi-Newton methods. These methods can be viewed as approximations
of Newton's method
(2.12.4). They
use
easily computable approximations of
H(x,)
or
H(x,)-1, and they are also descent methods. The best known
examples
are
the
Davidon-Fletcher-Powell
method
and
the
Broyden
methods.
3.
The conjugate gradient method. This
uses
a generalization of the idea of an
ortho&onal basis for a vector space
to
generate the directions
d,,
with the
114 ROOTFINDING FOR NONLINEAR EQUATIONS
Figure 2.8 Illustration
of
steepest
descent method.
directions related in an optimal way to the function
f(x)
being minimized.
In
Chapter 8, the conjugate gradient method
is
used for solving systems of
linear equations.
There are many other approaches to minimizing a function,
but
they are too
numerous to include here.
As
general references to the preceding ideas, see
Dennis
and
Schnabel (1983), Fletcher (1980), Gill et al. (1981), and Luenberger
(1984).
An
important and very different approach to minimizing a function
is
the
simplex method given in Neider and Mead
(1965), with a discussion given in Gill
et
al. (1981, p. 94) and Woods (1985, chap.
2).
This method uses only function
values (no derivative values), and it seems to
be
especially suitable for noisy
functions.
An
important project to develop programs for solving optimization problems
and nonlinear systems
is
under way at the Argonne National L::,boratory. The
program package
is
called MINPACK, and version 1
is
available [see More et al.
(1980)
and
More et al. (1984)].
It
contains routines for nonlinear systems and
nonlinear least squares problems. Future versions are intended to include pro-
grams for both unconstrained and constrained optimization problems.
Discussion of the Literature
There
is a large literature on methods for calculating the roots of a single
equation. See the books by Householder
(1970), Ostrowski (1973), and Traub
(1964) for a more extensive development than has been given here. Newton's
method
is
one
of
the most widely used methods, and its development
is
due to
many
people.
For
an historical account of contributions to it by Newton,
Raphson, and Cauchy, see Goldstine
(1977, pp. 64, 278).
For
computer programs, most people still use and individually program a
method that is especially suitable for their particular application. However, one
should strongly consider using one of the general-purpose programs that have
been developed in recent years and that are available in the commercial software
libraries. They are usually accurate, efficient,
and
easy to use. Among such
general-purpose programs, the ones based on Brent
(1973) and Dekker (1969) has