SOLUTION of NONLINEAR SYSTEMS
  Background:
       • A nonlinear system with n equations, n variables
                       f (x ,x ,...,x ) = 0
                        1  1  2      n
                       f (x ,x ,...,x ) = 0
                        2  1  2      n
                                         .
                                         .
                                         .
                       f (x ,x ,...,x ) = 0,
                        n  1  2      n
         or f(x) = 0. Solution(s)?
         NumericalmethodsusegeneralizationsofNewtonmethod.
       • Examples
         a)
                       2    2
                     x +2x −x −2x = 0
                       1    2    2     3
                           2    2
                         x −8x +10x = 0
                           1    2      3
                               2
                              x −7x x = 0.
                               1     2 3
         b) three intersecting radius-1 spheres:
                   (u−1)2+(v−1)2+w2 = 1
                   (u−1)2+v2+(w−1)2 = 1
                   u2 +(v −1)2 +(w−1)2 = 1.
                              1
                SOLUTION of NONLINEAR SYSTEMS
    Multivariable Newton Method
            • Use Taylor approximation of f near approximate x :
                                                                                                         0
                                                                                                      2
                    f(x) = f(x ) + Df(x )(x −x )+O(||x−x || ),
                                       0                 0              0                         0
                with Jacobian matrix
                                                   ∂f1 ∂f1 ... ∂f1 
                                                   ∂x1 ∂x2                  ∂xn 
                                                   ∂f2 ∂f2 ... ∂f2 
                                        Df =∂x1 ∂x2                         ∂xn 
                                                         .       .             .
                                                    .           .             .    
                                                         .       .    . . .    .
                                                       ∂fn ∂fn ... ∂fn
                                                       ∂x1 ∂x2               ∂xn
            • Linear approximation produces an iterative method:
                for each iteration, solve
                                  f(x ) + Df(x )(x                       −x )=0,
                                         k                 k      k+1          k
                to get multivariate Newton iteration
                                                                            −1
                                    x         =x −(Df(x )) f(x ).
                                       k+1           k                 k              k
                Note: (Df)−1 should not be explicitly computed.
                   Typical algorithm solves Df(x )s = −f(x ),
                                                                         k    k                k
                       then uses the update x                         =x +s .
                                                               k+1           k        k
            • Convergence Theory: if r is a solution,
                Df(r) is nonsingular and x0 is close enough to r,
                   then {x } converges quadratically to r.
                                 k
                This means ||x                   −r|| = O(||x −r||2).
                                          k+1                            k
                                                         2
          SOLUTION of NONLINEAR SYSTEMS
       • Multivariate Newton Examples
            2     2
        a) x +2x −x −2x =0,
            1     2    2     3
            2    2
          x −8x +10x =0,
            1    2      3
            2
          x −7x x =0, so
            1    2 3
        Df =
        Matlab
        f=@(x)[x(1)^2+2*x(2)^2-x(2)-2*x(3); ...
             x(1)^2-8*x(2)^2+10*x(3); x(1)^2-7*x(2)*x(3)];
        Df=@(x)[2*x(1) 4*x(2)-1 -2; 2*x(1) -16*x(2) 10;...
             2*x(1) -7*x(3) -7*x(2)];
        x = [1 1 1]’; k = 6;
        for i = 1:k, x = x-Df(x)\f(x); disp([i x’]), end
        disp(norm(f(x)));
        1        0.72791        0.73488        0.33023
        2        0.61379        0.53493        0.16056
        3        0.55262        0.43269        0.11125
        4        0.53125         0.4026        0.10077
        5        0.52917        0.40002             0.1
        6        0.52915            0.4             0.1
           2.4858e-09
        x = [.1 .1 .1]’; k = 6;
        for i = 1:k, x = x-Df(x)\f(x); disp([i x’]), end
        1         -36.35           -9.5           -0.8
        2        -18.125        -4.7032       -0.34947
        3        -9.0129        -2.3113       -0.12306
        ...
        Solution at 0 not found because of singular Df.
                              3
         SOLUTION of NONLINEAR SYSTEMS
                 2         2    2
        b) (x −1) +(x −1) +x =1,
            1         2         3
                 2   2         2
          (x −1) +x +(x −1) =1,
            1        2    3
           2         2         2
          x +(x −1) +(x −1) =1.
           1    2         3
        Matlab
        f=@(x)[(x(1)-1)^2+(x(2)-1)^2+x(3)^2-1; ...
               (x(1)-1)^2+x(2)^2+(x(3)-1)^2-1; ...
               x(1)^2+(x(2)-1)^2+(x(3)-1)^2-1];
        Df=@(x)[ 2*x(1)-2 2*x(2)-2 2*x(3); ...
                  2*x(1)-2 2*x(2) 2*x(3)-2; ...
                  2*x(1) 2*x(2)-2 2*x(3)-2 ];
        x = [0 0 0]’; k = 4;
        for i = 1:k, x = x-Df(x)\f(x); disp([i x’]), end
        disp(norm(f(x)));
        1          0.25          0.25          0.25
        2         0.325         0.325         0.325
        3      0.33323       0.33323       0.33323
        4      0.33333       0.33333       0.33333
           5.3649e-08
        x = [3 2 -1]’; k = 6;
        for i = 1:k, x = x-Df(x)\f(x); disp([i x’]), end
        disp(norm(f(x)));
        1          3.25          3.25          3.25
        2        1.9798        1.9798        1.9798
        3        1.3656        1.3656        1.3656
        4        1.0956        1.0956        1.0956
        5        1.0107        1.0107        1.0107
        6        1.0002        1.0002        1.0002
           0.00057171
                            4