The rootfinding problem becomes much more difficult when multiple variables and equations are involved.
Particular problems are often posed using scalar variables and equations.
While the equations of Example 4.5.1 are easy to solve by hand, in practice even establishing the existence and uniqueness of solutions for any particular system is typically quite difficult.
To extend rootfinding methods to systems, we will keep to the basic philosophy of constructing easily managed models of the exact function. As usual, the starting point is a linear model. We base it on the multidimensional Taylor series,
Because of the Jacobian’s role in (4.5.3), we may write J(x) as f′(x). Like any derivative, it is a function of the independent variable x.
The terms f(x)+J(x)h in (4.5.3) represent the linear part of f near x. If f is actually linear, i.e., f(x)=Ax−b, then the Jacobian matrix is the constant matrix A and the higher-order terms in (4.5.3) disappear.
With a method in hand for constructing a linear model for the vector system f(x), we can generalize Newton’s method. Specifically, at a root estimate xk, we set h=x−xk in (4.5.3) and get
Note that J−1f now plays the role that f/f′ had in the scalar case; in fact, the two are the same in one dimension. In computational practice, however, we don’t compute matrix inverses.
An extension of our series analysis of the scalar Newton’s method shows that the vector version is also quadratically convergent in any vector norm, under suitable circumstances and when the iteration converges at all.
An implementation of Newton’s method for systems is given in Function 4.5.2. Other than computing the Newton step using backslash and taking vector magnitudes with norm, Function 4.5.2 is virtually identical to the scalar version Function 4.3.2 presented earlier.
✍ Suppose that f(x)=Ax−b for a constant n×n matrix A and constant n×1 vector b. Show that Newton’s method converges to the exact root in one iteration.
Two curves in the (u,v) plane are defined implicitly by the equations ulogu+vlogv=−0.3 and u4+v2=1.
(a) ✍ Write the intersection of these curves in the form f(x)=0 for two-dimensional f and x.
(b) ✍ Find the Jacobian matrix of f.
(c) ⌨ Use Function 4.5.2 to find an intersection point starting from u=1, v=0.1.
(d) ⌨ Use Function 4.5.2 to find an intersection point starting from u=0.1, v=1.
Two elliptical orbits (x1(t),y1(t)) and (x2(t),y2(t)) are described by the equations
(b) ✍ Write out a 2×2 nonlinear system of equations that describes an intersection of these orbits. (Note: An intersection is not the same as a collision—they don’t have to occupy the same point at the same time.)
(c) ✍ Write out the Jacobian matrix of this nonlinear system.
(d) ⌨ Use Function 4.5.2 to find all of the unique intersections.
⌨ Suppose one wants to find the points on the ellipsoid x2/25+y2/16+z2/9=1 that are closest to and farthest from the point (5,4,3). The method of Lagrange multipliers implies that any such point satisfies
(a) Write out this system in the form f(u)=0. (Note that the system has four variables to go with the four equations.)
(b) Write out the Jacobian matrix of this system.
(c) Use Function 4.5.2 with different initial guesses to find the two roots of this system. Which is the closest point to (5,4,3), and which is the farthest?
⌨ Any three noncollinear points in the plane determine a unique circle. Suppose the points are given as (xi,yi) for i=1,2,3. We can define the circle in terms of its center (a,b) and radius r. Then