% Activate the following line by filling in the right side. If for example the name of the root file is Main.tex, write
% "...root = Main.tex" if the chapter file is in the same directory, and "...root = ../Main.tex" if the chapter is in a subdirectory.
 
%!TEX root = ../Metodos_Numericos.tex  

\chapter[Sistemas de Ecuaciones]{Solución de sistemas de ecuaciones lineales y no-lineales}

\section{Solución de Sistemas lineales}

Se llama sistema de ecuaciones lineales a un conjunto de ecuaciones de la forma:

donde $x_1$, ..., $x_n$ son las incógnitas,  $b_1$, ..., $b_n$ se denominan términos independientes y los números $a_{i,j}$ se llaman coeficientes de las incógnitas, formando una matriz que denominaremos $A$, matriz de coeficientes. 
\begin{equation*}
\left[
\begin{array}{cccc}
a_{1,1} & a_{1,2} & \cdots  & a_{1,n}\\
a_{2,1} & a_{2,2} & \cdots  & a_{2,n} \\
\vdots & \vdots  & \cdots & \vdots\\	
a_{n,1} & a_{n,2} &  \cdots & a_{n,n}
\end{array}
\right]
\left[
\begin{array}{c}
x_1 \\
x_2 \\	
\vdots \\
x_n
\end{array}
\right]
= 	
\left[
\begin{array}{c}
b_1 \\
b_2 \\
\vdots\\
b_n
\end{array}
\right]
\end{equation*}

En forma compacta 

\begin{equation*}
A x = b
\end{equation*}

En las siguientes subsecciones se analizarán algunos de los métodos para resolver el sistema de ecuaciones.

\subsection{Método iterativo de Jacobi}

Consideremos el siguiente sistema de ecuaciones.

\begin{equation*}
\left[
\begin{array}{ccc}
a_{1,1} & a_{1,2} & a_{1,3}\\
a_{2,1} & a_{2,2} & a_{2,3} \\	
a_{3,1} & a_{3,2} & a_{3,3}
\end{array}
\right]
\left[
\begin{array}{c}
x_1 \\
x_2 \\	
x_3
\end{array}
\right]
= 	
\left[
\begin{array}{c}
b_1 \\
b_2 \\
b_3
\end{array}
\right]
\end{equation*}

Vamos a representar cada una de las variables en términos de ellas mismas.

\begin{eqnarray*} 
x_1 & = & \frac{ b_1 - a_{1,2} x_2 - a_{1,3} x_3}{a_{1,1}} \\
x_2 & = & \frac{ b_2 - a_{2,1} x_1 - a_{2,3} x_3}{a_{2,2}} \\
x_3 & = & \frac{ b_3 - a_{3,1} x_1 - a_{3,2} x_2}{a_{3,3}} \\
\end{eqnarray*}

Lo cual nos sugiere el siguiente esquema iterativo de solución.

\begin{eqnarray*} 
x_1^{(t+1)} & = & \frac{ b_1 - a_{1,2} x_2^{(t)} - a_{1,3} x_3^{(t)}}{a_{1,1}} \\
x_2^{(t+1)} & = & \frac{ b_2 - a_{2,1} x_1^{(t)} - a_{2,3} x_3^{(t)}}{a_{2,2}} \\
x_3^{(t+1)} & = & \frac{ b_3 - a_{3,1} x_1^{(t)} - a_{3,2} x_2^{(t)}}{a_{3,3}} \\
\end{eqnarray*}
 
En general podemos escribir como

\begin{equation*} 
x_k^{(t+1)} = \frac{b_k – \sum_{l=1, l \neq k}^N a_{k,l}x_l^{(t)}}{a_{k,k}} 
\end{equation*}

La implementación en Matlab es

\begin{verbatim}
function y = Jacobi(A, x, b)

N = length(x);
y = zeros(N,1);

for iter=1:100000
    for k = 1:N
        suma =0;
        for l= 1:N
            if k ~= l
                suma = suma + A(k,l)*x(l);
            end;
        end;
        y(k) = (b(k) - suma)/A(k,k);
    end;
    if sqrt(norm(x-y)) < 1e-6 
        break;
    else
       x = y;
    end;
end;\end{verbatim}

\subsection{Algoritmo iterativo de Gauss-Seidel}

El cambio que debemos hacer respecto al de Jacobi, es que las variables nuevas son utilizadas una vez que se realiza el cálculo de ellas así, para un sistema de tres ecuaciones tendremos:

\begin{eqnarray*} 
x_1^{(t+1)} & = & \frac{ b_1 - a_{1,2} x_2^{(t)} - a_{1,3} x_3^{(t)}}{a_{1,1}} \\
x_2^{(t+1)} & = & \frac{ b_2 - a_{2,1} x_1^{(t+1)} - a_{2,3} x_3^{(t)}}{a_{2,2}} \\
x_3^{(t+1)} & = & \frac{ b_3 - a_{3,1} x_1^{(t+1)} - a_{3,2} x_2^{(t+1)}}{a_{3,3}} \\
\end{eqnarray*}

La implementación en Matlab es

\begin{verbatim}
function y = Gauss_Seidel(A, x, b)

N = length(x);
y = zeros(N,1);

for iter=1:100000
    for k = 1:N
        suma =0;
        for l= 1:N
            if k ~= l
                suma = suma + A(k,l)*x(l);
            end;
        end;
        x(k) = (b(k) - suma)/A(k,k);
    end;
    if sqrt(norm(x-y)) < 1e-6 
        break;
    else
        y=x;
    end;
end;
\end{verbatim}

\subsubsection{Ejemplo}

Resolver el siguiente sistema de ecuaciones utilizando el método de Jacobi y comparar con el método de Gauss-Seidel.
 
\begin{equation*}
\left[
\begin{array}{ccc}
 4  &  -1 &  1  \\
 4  &  -8 &  1 \\
-2  &   1 &  5
\end{array}
\right]
\left[
\begin{array}{c}
x_1 \\
x_2 \\	
x_3
\end{array}
\right]
= 	
\left[
\begin{array}{c}
  7 \\
-21 \\
 15
\end{array}
\right]
\end{equation*}

Las primeras 20 iteraciones del algoritmo de Jacobi son :

\begin{center}
\begin{tabular}{|r|r|r|r|}
\hline
$k$ & $x_1$ & $x_2$ & $x_3$ \\
\hline
1 & 1.0000 & 2.0000 & 2.0000 \\
2 & 1.7500 & 3.3750 & 3.0000 \\
3 & 1.8438 & 3.8750 & 3.0250\\
4 & 1.9625 & 3.9250 & 2.9625\\
5 & 1.9906 & 3.9766 & 3.0000\\
6 & 1.9941 & 3.9953 & 3.0009\\
7 & 1.9986 & 3.9972 & 2.9986\\
8 & 1.9996 & 3.9991 & 3.0000\\
9 & 1.9998 & 3.9998 & 3.0000\\
10 & 1.9999 & 3.9999 & 2.9999 \\
11 & 2.0000 & 4.0000 & 3.0000 \\
12 & 2.0000 & 4.0000 & 3.0000 \\
\hline
\end{tabular}
\end{center}

Para correr hacer 

\begin{verbatim}
>> Jacobi([4 -1 1; 4 -8 1; -2 1 5], [1,2,2]', [7,-21, 15]')

ans =

    2.0000
    4.0000
    3.0000
\end{verbatim}

La solución utilizando Gauss-Seidel es :

\begin{center}
\begin{tabular}{|r|r|r|r|}
\hline
$k$ & $x_1$ & $x_2$ & $x_3$ \\
\hline
1 & 1.0000 & 2.0000 & 2.0000 \\
2 & 1.7500 & 3.7500 & 2.9500 \\
3 & 1.9500 & 3.9688 & 2.9863 \\
4 & 1.9956 & 3.9961 & 2.9990 \\
5 & 1.9993 & 3.9995 & 2.9998 \\
6 & 1.9999 & 3.9999 & 3.0000 \\
7 & 2.0000 & 4.0000 & 3.0000 \\
8 & 2.0000 & 4.0000 & 3.0000 \\
\hline
\end{tabular}
\end{center}
 
Note que Gauss-Seidel requiere de 7 iteraciones mientras Jacobi de 11, para convergir. Para correr hacer

\begin{verbatim}
>> Gauss_Seidel([4 -1 1; 4 -8 1; -2 1 5], [1,2,2]', [7,-21, 15]')

ans =

    2.0000
    4.0000
    3.0000
\end{verbatim}

\subsection{Ejemplo matrices dispersas} 
Resolver el sistema de ecuaciones utilizando matrices dispersas y los métodos de Jacobi y Gauss-Seidel

\[
A = 
\left[
\begin{array}{ccccc}
    10 &   -1  &  -1  &  -1  &   1 \\
     1  &  4   &  0   &  0   &  0\\
     1  &  0   &  4   &  0   &  0\\
     1  &   0  &   0  &   4  &   0\\
     1  &   0  &   0  &   0   &  3
\end{array}
\right]
\left[
\begin{array}{c}
x_1\\
x_2\\
x_3\\
x_4\\
x_5
\end{array}
\right]
= 
\left[
\begin{array}{c}
1\\
2\\
3\\
4\\
5
\end{array}
\right]
\]

\begin{verbatim}
A = sparse(5,5);
A(1,1) = 10;
A(1,2) = -1;
A(1,3) = -1;
A(1,4) = -1;
A(1,5) = 1;
A(2,1) = 1;
A(3,1) = 1;
A(4,1) = 1;
A(5,1) = 1;
A(2,2) = 4;
A(3,3) = 4;
A(4,4) = 4;
A(5,5) = 3;
disp(A);

size(A)

b = [1,2,3,4,5]';

disp(b);
Jacobi(A,[0,0,0,0,0]', b)
Gauss_Seidel(A,[0,0,0,0,0]', b)
\end{verbatim}

La solución para ambos métodos es:

\begin{verbatim}
ans =

    0.1520
    0.4620
    0.7120
    0.9620
    1.6160


ans =

    0.1520
    0.4620
    0.7120
    0.9620
    1.6160
 \end{verbatim}

\subsection[Eliminación Gaussiana]{Solución de sistemas lineales de ecuaciones por el método de Eliminación Gaussiana}

Consideremos que tenemos un sistema lineal Ax=b, donde la matriz A no tiene las condiciones de ser triangular superior.
$$
\left[%
\begin{array}{ccc}
a_{1,1} & a_{1,2} & a_{1,3} \\
a_{2,1} & a_{2,2} & a_{2,3} \\
a_{3,1} & a_{3,2} & a_{3,3} \\
\end{array}%
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c}
  b_1 \\
  b_2 \\
  b_3 \\
\end{array}
\right]
$$

Comenzaremos por despejar  $x_1$ de la ecuación (I)

$$
x_1 = \frac{b_1 - a_{1,2}x_2 - a_{1,3}x_3}{a_{1,1}}
$$

y sustituimos en las ecuaciones (II)

$$
a_{2,1} \left( \frac{b_1 - a_{1,2}x_2 - a_{1,3}x_3}{a_{1,1}} \right) + a_{2,2} x_2
+ a_{2,3}x_3 = b_2
$$

agrupando términos semejantes tenemos:

$$
\left(a_{2,2}-\frac{a_{2,1} *a_{1,2}}{a_{1,1}} \right) x_2 + \left (a_{2,3}-a_{2,1}
\frac{a_{1,3}}{a_{1,1}} \right) x_3 = \left (b_2 - a_{2,1} \frac{b_1}{a_{1,1}} \right)
$$

Ahora sustituimos en la ecuación (III)

$$
\left(a_{3,1}\frac{b_1 - a_{1,2}x_2 - a_{1,3}x_3}{a_{1,1}} \right ) + a_{3,2}x_2
+ a_{3,3}x_3 = b_3
$$
$$ \left(a_{3,2}-
\frac{a_{3,1} a_{1,2}}{a_{1,1}} \right)x_2+ \left(a_{3,3}-
\frac{a_{3,1}a_{1,3}}{a_{1,1}}\right) x_3 = \left(b_3-a_{3,1} \frac{b_1}{a_{1,1}}\right)
$$

Lo cual nos da un nuevo sistema simplificado dada por

$$
\left[
\begin{array}{ccc}
a_{1,1} & a_{1,2} & a_{1,3} \\
      0 & a'_{2,2} & a'_{2,3} \\
      0 & a'_{3,2} & a'_{3,3} \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] = \left[
\begin{array}{c}
  b_1 \\
  b'_2 \\
  b'_3 \\
\end{array}
\right ]
$$

donde los valores de $a'_{i,j}$ los calculamos como:
\begin{eqnarray*}
a'_{i,j} & = & a_{i,j}-a_{i,k} \frac{a_{k,j}}{a_{k,k}}\\
b'_i & = & b_i - a_{i,k} \frac{b_k}{a_{k,k}}
\end{eqnarray*}

 Si repetimos el procedimiento, ahora para $x_1$, tendremos un sistema dado como:

$$
\left[
\begin{array}{ccc}
a_{1,1} & a_{1,2} & a_{1,3} \\
      0 & a'_{2,2} & a'_{2,3} \\
      0 & 0 & a''_{3,3} \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c}
  b_1 \\
  b'_2 \\
  b''_3 \\
\end{array}
\right]
$$

Solución de un sistema Triangular superior

\begin{equation*}
\left[
\begin{array}{cccc}
a_{1,1} & a_{1,2} & a_{1,3} & a_{1,4}\\
         0 & a_{2,2} & a_{2,3}  & a_{2,4}\\
         0 &          0 & a_{3,3}  & a_{3,4} \\
         0 &          0 &           0  & a_{4,4}
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right] =
\left[
\begin{array}{c}
  b_1 \\
  b_2 \\
  b_3 \\
  b_4\\
\end{array}
\right]
\end{equation*}

Podemos resolverlo empezando con la ecuación 4 cuya solución es la más simple, para luego solucionar 3 y así sucesivamente. Este esquema de solución queda

\begin{eqnarray*}
x_4 & = & \frac{b_4}{a_{4,4}}\\
x_3 & = & \frac{b_3 - a_{3,4} x_4}{a_{3,3}}\\
x_2 & = & \frac{b_2 - a_{2,3} x_3 - a_{2, 4} x_4}{a_{2,2}}\\
x_1 & = & \frac{b_1 - a_{1,2} x_2 - a_{1,3} x_3 - a_{1,4} x_4}{a_{1,1}}
\end{eqnarray*}

En general podemos calcular la solución haciendo

\begin{equation*}
x_k = \frac{b_k - \sum_{i=k+1}^N a_{k,i} x_i}{a_{k,k}}
\end{equation*}


La implementación en Matlab queda

\begin{verbatim}
function x = Eliminacion_Gaussiana(A, b)

N = length(b);

x = zeros(N,1);
    
for k =1:1:N-1
    for n = k+1:1:N
        b(n) = b(n) - A(n,k)*b(k)/A(k,k);
        for m=N:-1:k
            A(n,m) = A(n,m) - A(n,k)*A(k,m)/A(k,k);
        end;
    end;
end;

for k=N:-1:1
    suma = 0;
    for m=k+1:1:N
        suma = suma +A(k,m)*x(m);
    end;
    x(k) = (b(k)-suma)/A(k,k);
end;
\end{verbatim}

\noindent 
que corresponde a un sistema triangular superior que podemos solucionar utilizando sustitución hacia atrás.

\subsubsection{Ejemplo 1}

Dado el sistema lineal de ecuaciones, calcular el sistema Triangular Superior utilizando el algoritmo de Eliminación Gaussiana

\begin{equation*}
\left[
\begin{array}{r r r r}
    10  &  -1 &    2  &   1 \\
    -1   & 15 &   -3  &    1\\
     2   & -3  &    6  &  -3\\
     1   &  1  &   -3  &   7\\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right]
=
\left[
\begin{array}{r }
  1 \\
  2 \\
  2 \\
  1\\
\end{array}
\right]
\end{equation*}

Primer paso $k=1$

\begin{equation*}
\left(
\begin{array}{cccc}
 10 & -1 & 2 & 1 \\
 \\
 0 & \frac{149}{10} & -\frac{14}{5} & \frac{11}{10} \\
 \\
 0 & -\frac{14}{5} & \frac{28}{5} & -\frac{16}{5} \\
 \\
 0 & \frac{11}{10} & -\frac{16}{5} & \frac{69}{10}
\end{array}
\right)\left(
\begin{array}{c}
 x_1 \\
 \\
 x_2 \\
 \\
 x_3 \\
 \\
 x_4
\end{array}
\right)\text{= } \left(
\begin{array}{c}
 1 \\
 \\
 \frac{21}{10} \\
 \\
 \frac{9}{5} \\
 \\
 \frac{9}{10}
\end{array}
\right)\end{equation*}

\begin{equation*}
\left[
\begin{array}{r r r r}
   10.0000   & -1.0000  &  2.0000  &  1.0000 \\
         0   & 14.9000 &  -2.8000  &  1.1000\\
         0   & -2.8000  &  5.6000  &  -3.2000\\
         0   & 1.1000 &  -3.2000   &  6.9000\\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right]
=
\left[
\begin{array}{r }
    1.0000 \\
    2.1000\\
    1.8000\\
    0.9000\\
\end{array}
\right]
\end{equation*}

Segundo paso $k=2$

\begin{equation*}
\left(
\begin{array}{cccc}
 10 & -1 & 2 & 1 \\
 \\
 0 & \frac{149}{10} & -\frac{14}{5} & \frac{11}{10} \\
 \\
 0 & 0 & \frac{756}{149} & -\frac{446}{149} \\
 \\
 0 & 0 & -\frac{446}{149} & \frac{1016}{149}
\end{array}
\right)\left(
\begin{array}{c}
 x_1 \\
 \\
 x_2 \\
 \\
 x_3 \\
 \\
 x_4
\end{array}
\right)\text{= } \left(
\begin{array}{c}
 1 \\
 \\
 \frac{21}{10} \\
 \\
 \frac{327}{149} \\
 \\
 \frac{111}{149}
\end{array}
\right)
\end{equation*}


\begin{equation*}
\left[
\begin{array}{r r r r}
   10.0000  &  -1.0000  &  2.0000  &  1.0000\\
         0  & 14.9000 &  -2.8000   & 1.1000\\
         0    &     0 &   5.0738 &  -2.9933\\
         0    &     0  &  -2.9933 &   6.8188\\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right]
=
\left[
\begin{array}{r }
    1.0000\\
    2.1000\\
    2.1946\\
    0.7450\\
\end{array}
\right]
\end{equation*}

Tercer paso $k=3$

\begin{equation*}
\left(
\begin{array}{cccc}
 10 & -1 & 2 & 1 \\
 \\
 0 & \frac{149}{10} & -\frac{14}{5} & \frac{11}{10} \\
 \\
 0 & 0 & \frac{756}{149} & -\frac{446}{149} \\
 \\
 0 & 0 & 0 & \frac{955}{189}
\end{array}
\right)\left(
\begin{array}{c}
 x_1 \\
 \\
 x_2 \\
 \\
 x_3 \\
 \\
 x_4
\end{array}
\right)\text{= } \left(
\begin{array}{c}
 1 \\
 \\
 \frac{21}{10} \\
 \\
 \frac{327}{149} \\
 \\
 \frac{257}{126}
\end{array}
\right)
\end{equation*}

\begin{equation*}
\left[
\begin{array}{r r r r}
   10.0000  & -1.0000 &   2.0000  &  1.0000\\
         0   &14.9000  & -2.8000 &  1.1000\\
         0   &      0  &  5.0738  & -2.9933\\
         0   &     0   &      0  &  5.0529\\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right]
=
\left[
\begin{array}{r }
    1.0000\\
    2.1000\\
    2.1946\\
    2.0397\\
\end{array}
\right]
\end{equation*}

Sustitución hacia atras

Paso 1 $k=4$

\begin{eqnarray*}
x_4 & = & \frac{b_4}{a_{4,4}}\\
x_4 & = & \frac{2.0397}{5.0529} = 0.4037
\end{eqnarray*}

Paso 2 $k=3$
\begin{eqnarray*}
x_3 & = & \frac{b_3 - a_{3,4} x_4}{a_{3,3}}\\
x_3 & = & \frac{2.1946 - (-2.9933) (0.4037)}{5.0738} = 0.6707
\end{eqnarray*}

Paso 3 $k=2$

\begin{eqnarray*}
x_2 & = & \frac{b_2 - a_{2,3} x_3 - a_{2, 4} x_4}{a_{2,2}}\\
x_2 & = & \frac{2.1000 - (-2.8000 )(0.6707) - (1.1000) (0.4037)}{14.9000 } = 0.2372
\end{eqnarray*}

Paso 4 $k=1$
\begin{eqnarray*}
x_1 & = & \frac{b_1 - a_{1,2} x_2 - a_{1,3} x_3 - a_{1,4} x_4}{a_{1,1}}\\
x_1 & = & \frac{1.0000 - (-1.0000) (0.2372) - (2.0000) (0.6707) - (1.0000)(0.4037)}{10.0000} = - 0.0508
\end{eqnarray*}

La solución en Matlab es 

\begin{verbatim}
>> Eliminacion_Gaussiana([10 -1 2 1; -1 15 -3 1; 2 -3 6 -3; 1 1 -3 7], [1,2,2,1]')

   10.0000   -1.0000    2.0000    1.0000
         0   14.9000   -2.8000    1.1000
         0         0    5.0738   -2.9933
         0         0         0    5.0529

    1.0000    2.1000    2.1946    2.0397

ans =

   -0.0508
    0.2372
    0.6707
    0.4037
\end{verbatim}

\subsubsection{Ejemplo 2}

Calcular el sistema triangular superior utilizando eliminación Gaussiana.

$$
5x + 2y + 1z = 3
$$
$$
2x + 3y - 3z = -10
$$
$$
1x - 3y + 2z = 4 \\
$$

Primer paso

$$
5x +      2y +      1z = 3
$$
$$
0  + (11/5)y-(17/5)z =-(56/5)
$$
$$
0  - (17/5)y +  (9/5)z = (17/5)
$$

Segundo paso

$$
5x +      2y +      1z   = 3
$$
$$
0x + (11/5)y-( 17/5)z =-(56/5)
$$
$$
0x -      0y -(38/11)z =-(153/11)
$$

La solución en Matlab es

\begin{verbatim}
>> Eliminacion_Gaussiana([5 2 1; 2 3 -3; 1 -3 2], [3, -10, 4])
    5.0000    2.0000    1.0000
         0    2.2000   -3.4000
         0         0   -3.4545

    3.0000
  -11.2000
  -13.9091


ans =

   -0.6579
    1.1316
    4.0263
\end{verbatim}

\subsubsection{Ejemplo 3}

Resolver el sistema de ecuaciones

$$
\left[
\begin{array}{rrr}
3 &  -1 &  -1\\
-1 &  1 &  0 \\
-1 &  0 &  1\\
\end{array}
\right]
\left[
\begin{array}{r}
  x \\
  y \\
  z \\
\end{array}
\right] = \left[
\begin{array}{r}
  0 \\
  1 \\
  1 \\
\end{array}
\right ]
$$

Aplicando Eliminación Gaussiana nos queda.

$$
\left[
\begin{array}{rrr}
3 &  -1 &  -1 \\
0 &  2/3 &  -1/3\\
0 &  0 &  1/2 \\
\end{array}
\right]
\left[
\begin{array}{r}
  x \\
  y \\
  z \\
\end{array}
\right] = \left[
\begin{array}{r}
  0 \\
  1 \\
  3/2 \\
\end{array}
\right ]
$$

La solución en Matlab es

\begin{verbatim}

>> Eliminacion_Gaussiana([3 -1 -1; -1 1 0; -1 0 1], [0 1 1])
    3.0000   -1.0000   -1.0000
         0    0.6667   -0.3333
         0         0    0.5000

         0
    1.0000
    1.5000


ans =

    2.0000
    3.0000
    3.0000
\end{verbatim}

\subsubsection{Problema con la Eliminación Gaussiana}

Un ejemplo de sistema donde es necesario hacer un cambio de renglón por renglón para que tenga solución es el siguiente sistema. 

$$
\left[
\begin{array}{rrr}
1 & 2 & 6 \\
4 & 8 & -1 \\
-2 & 3 &  5\\
\end{array}
\right]
$$

Aplicando el primer paso de la eliminación gaussiana tenemos:

$$
\left[
\begin{array}{rrr}
1 & 2 & 6 \\
0 & 0 & -25 \\
0 & 7 & 17 \\
\end{array}
\right]
$$

\noindent
note, que aparece un cero en el elemento 22, lo cual nos da un sistema sin solución. Permutando los renglones 2 y 3 el sistema tiene solución.

$$
\left[
\begin{array}{rrr}
1 &  2 & 6 \\
-2 &  3 &  5 \\
4 &  8 & -1 \\
\end{array}
\right]
$$

La solución en Matlab tenemos
\begin{verbatim}
>> Eliminacion_Gaussiana([1 2 6; 4 8 -1; -2 3 5], [1 1 1])
     1     2     6
     0     0   -25
     0   NaN   Inf

     1
    -3
   Inf


ans =

   NaN
   NaN
   NaN
\end{verbatim}

Haciendo el cambio de renglones tenemos

\begin{verbatim}
>> Eliminacion_Gaussiana([1 2 6;  -2 3 5; 4 8 -1], [1 1 1])
     1     2     6
     0     7    17
     0     0   -25

     1
     3
    -3


ans =

    0.0057
    0.1371
    0.1200
\end{verbatim}

\subsection{Gauss-Jordan}

El método de Gauss-Jordan es una variación de la eliminación de Gauss. La principal diferencia consiste em que cuando una incógnita  se elimina en el método de Gauss-Jordan, esta es eliminada de todas las ecuaciones en lugar de hacerlo en la subsecuentes. Además todos los renglones se normalizan al dividirlos entre su elemento privote. De esta forma, el paso de eliminación genera una matriz identidad en vez de una matriz triangular. En consecuencia no es necesario usar la sustitución hacia atrás para obtener la solución. 

Comenzaremos planteando un sistema de 3 ecuaciones con 3 incógnitas dado por 

\begin{equation*}
\left[
\begin{array}{c|c|c}
a_{1,1} & a_{1,2} & a_{1,3} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} \\
\hline
a_{3,1} & a_{3,2} & a_{3,3} \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c}
  b_1 \\
  b_2 \\
  b_3 \\
\end{array}
\right]
\end{equation*}

Este sistema de ecuaciones lo podemos escribir como:

\begin{equation*}
Ax = I b
\end{equation*}

donde $I$ es la matriz identidad

\begin{equation*}
\left[
\begin{array}{ccc}
a_{1,1} & a_{1,2} & a_{1,3} \\
a_{2,1} & a_{2,2} & a_{2,3} \\
a_{3,1} & a_{3,2} & a_{3,3} \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{ccc}
1 & 0 & 0 \\
0 & 1 & 0 \\
0 & 0 & 1 \\
\end{array}
\right]
\left[
\begin{array}{c}
  b_1 \\
  b_2 \\
  b_3 \\
\end{array}
\right]
\end{equation*}

Lo cual da lugar a la matriz aumentada que hemos utilizado en el procedimiento de Gauss-Jordan

\begin{equation}
\left[
\begin{array}{ccc|ccc|c}
a_{1,1} & a_{1,2} & a_{1,3} & 1 & 0 & 0 & b_1 \\
a_{2,1} & a_{2,2} & a_{2,3} & 0 & 1 & 0 & b_2\\
a_{3,1} & a_{3,2} & a_{3,3} & 0 & 0 & 1 & b_3\\
\end{array}
\right]  \label{eq:Gauss_Jordan}
\end{equation}

\textbf{Primer paso}

Vamos a despejar la variable $x_1$ del sistema y la sustituimos en las otras dos ecuaciones

\begin{equation} 
\frac{b_1 - a_{1,2} x_2 - a_{1,3} x_3} {a_{1,1}} = x_1   \label{eq:eli_par_1}
\end{equation}

Sustituyendo en la ecuación 2

\begin{equation*}
a_{2,1}\left( \frac{b_1 - a_{1,2} x_2 - a_{1,3} x_3} {a_{1,1}} \right)+ a_{2,2} x_2 + a_{2,3} x_3 = b_2
\end{equation*}
\begin{equation}
0 x_1 + \left(\frac{a_{2,1}}{a_{1,1}}\right) b_1 + \left( a_{2,2} - \frac{a_{2,1}a_{1,2} }{a_{1,1}} \right) x_2 + \left( a_{2,3} -\frac{a_{2,1} a_{1,3}} {a_{1,1}} \right) x_3 = b_2  \label{eq:eli_par_2}
\end{equation} 

Sustituyendo en la ecuación 3

\begin{equation*}
a_{3,1}\left( \frac{b_1 - a_{1,2} x_2 - a_{1,3} x_3} {a_{1,1}} \right)+ a_{3,2} x_2 + a_{3,3} x_3 = b_3
\end{equation*}

\begin{equation}
0 x_1 + \left(\frac{a_{3,1}}{a_{1,1}}\right) b_1 + \left( a_{3,2} - \frac{a_{3,1}a_{1,2} }{a_{1,1}} \right) x_2 + \left( a_{3,3} - \frac{a_{3,1} a_{1,3}} {a_{1,1}} \right) x_3 = b_3  \label{eq:eli_par_3}
\end{equation}

En forma matricial podemos escribir las ecuaciones \ref{eq:eli_par_1}, \ref{eq:eli_par_2} y \ref{eq:eli_par_3}

\begin{equation*}
\left[
\begin{array}{c|c|c}
1 & \frac{a_{1,2}}{a_{1,1}} & \frac{a_{1,3}}{a_{1,1}} \\
\hline
0 & \left( a_{2,2} - \frac{a_{2,1}a_{1,2} }{a_{1,1}} \right) &\left( a_{2,3} - \frac{a_{2,1} a_{1,3}} {a_{1,1}} \right) \\
\hline
0 & \left( a_{3,2} - \frac{a_{3,1}a_{1,2} }{a_{1,1}} \right) &\left( a_{3,3} - \frac{a_{3,1} a_{1,3}} {a_{1,1}} \right)  \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c|c|c}
\frac{1}{a_{1,1}} & 0 & 0 \\
\hline
-\frac{a_{2,1}}{a_{1,1}} & 1 & 0 \\
\hline
-\frac{a_{3,1}}{a_{1,1}} & 0 & 1 \\
\end{array}
\right]
\left[
\begin{array}{c}
  b_1 \\
  b_2 \\
  b_3 \\
\end{array}
\right]
\end{equation*}

La matriz aumentada dada por (\ref{eq:Gauss_Jordan} ) tenemos:

\begin{equation*}
\left[
\begin{array}{c|c|c||c|c|c||c}
1 & \frac{a_{1,2}}{a_{1,1}}  & \frac{a_{1,3}}{a_{1,1}} & \frac{1}{a_{1,1}}  & 0 & 0 & \frac{1}{a_{1,1}}  b_1\\
\hline
0 & \left( a_{2,2} - \frac{a_{2,1}a_{1,2} }{a_{1,1}} \right) & \left( a_{2,3} - \frac{a_{2,1} a_{1,3}} {a_{1,1}} \right) & -\frac{a_{2,1}}{a_{1,1}}  & 1 & 0 & b_2 - \frac{a_{2,1}}{a_{1,1}} b_1\\
\hline
0 & \left( a_{3,2} - \frac{a_{3,1}a_{1,2} }{a_{1,1}} \right) &\left( a_{3,3} - \frac{a_{3,1} a_{1,3}} {a_{1,1}} \right)  & -\frac{a_{3,1}}{a_{1,1}} & 0 & 1 & b_3 - \frac{a_{3,1}}{a_{1,1}}  b_1\\
\end{array}
\right]  
\end{equation*}

\textbf{Segundo Paso}

Para un sistema equivalente 

\begin{equation}
\left[
\begin{array}{c|c|c}
1 & a'_{1,2} & a'_{1,3} \\
\hline
0 & a'_{2,2} & a'_{2,3}  \\
\hline
0 & a'_{3,2} & a'_{3,3}   \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c|c|c}
c'_{1,1} & 0 & 0 \\
\hline
c'_{2,1} & 1 & 0 \\
\hline
c'_{3,1} & 0 & 1 \\
\end{array}
\right]
\left[
\begin{array}{c}
  b'_1 \\
  b'_2 \\
  b'_3 \\
\end{array}
\right]
\end{equation}

Despejamos $x_2$ de la ecuación 2 

\begin{equation*}
x_2 = \frac{c'_{2,1} b'_1 + b'_2 - a'_{2,3} x_3}{a'_{2,2}}
\end{equation*}

Sustituyendo en la ecuación 1

\begin{equation*}
x1 + a'_{1,2} \left( \frac{c'_{2,1} b'_1 + b'_2 - a'_{2,3} x_3}{a'_{2,2}} \right)+ a'_{1,3} = c'_{1,1}b'_1
\end{equation*}


\begin{equation*}
x1 + 0 x_2 + \left (a'_{1,3}  - \frac{a'_{1,2}a'_{2,3}}{a'_{2,2}} x_3\right )= \left( c'_{1,1}b'_1 -  \frac{a'_{1,2}c'_{2,1}}{a'_{2,2}} \right) - \frac{a'_{2,1}} {a'_{2,2}} b'_2 
\end{equation*}

Sustituyendo en la ecuación 3

\begin{equation*}
0 x_1 + 0 x_2 + \left( a'_{3,3} - \frac{a'_{3,2} a'_{2,3}}{a'_{2,2}} \right) x_3 = \left ( c'_{3,1} - \frac{a'_{3,2} c'_{2,1}}{a'_{2,2}} \right) b'_1- \frac{a'_{3,2}}{a'_{2,2}} b'2 + b'3
\end{equation*}

En forma matricial

\begin{equation*}
\left[
\begin{array}{c|c|c}
1 & 0 & a'_{1,3} - \frac{a'_{1,2}a'_{2,3}}{a'_{2,2}}\\
\hline
0 & 1 & \frac{a'_{2,3}}{a'_{2,2}}  \\
\hline
0 & 0 & a'_{3,3} - \frac{a'_{3,2} a'_{2,3}}{a'_{2,2}}   \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c|c|c}
c'_{1,1} - \frac{a'_{1,2} c'_{2,1}}{a'_{2,2}} & \frac{a'_{1,2}}{a'_{2,2}} & 0 \\
\hline
\frac{c'_{2,1}}{a'_{2,2}} & \frac{1}{a'_{2,2}} & 0 \\
\hline
c'_{3,1} - \frac{a'_{3,2} c'_{2,1}}{a'_{2,2}} & -\frac{a'_{3,2}}{a'_{2,2}} & 1 \\
\end{array}
\right]
\left[
\begin{array}{c}
  b'_1 \\
  b'_2 \\
  b'_3 \\
\end{array}
\right]
\end{equation*}

\textbf{Tercer Paso}

Dado el sistema equivalente 

\begin{equation*}
\left[
\begin{array}{c|c|c}
1 & 0 & a''_{1,3}\\
\hline
0 & 1 &  a''_{2,3} \\
\hline
0 & 0 & a''_{3,3}   \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c|c|c}
c''_{1,1}  & c''_{1,2} & 0 \\
\hline
c''_{2,1} & c''_{2,2} & 0 \\
\hline
c''_{3,1}  & c''_{3,2} & 1 \\
\end{array}
\right]
\left[
\begin{array}{c}
  b''_1 \\
  b''_2 \\
  b''_3 \\
\end{array}
\right]
\end{equation*}

Despejamos $x_3$ de la ecuación 3 y sustituimos en la ecuaciones 1 y 2

\begin{equation*}
\left[
\begin{array}{c|c|c}
1 & 0 & 0\\
\hline
0 & 1 &  0 \\
\hline
0 & 0 & 1   \\
\end{array}
\right]
\left[
\begin{array}{c}
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right] =
\left[
\begin{array}{c|c|c}
c''_{1,1}  & c''_{1,2} - \frac{a''_{1,3} c''_{3, 2}}{a''_{3,3}}& -\frac{a''_{1,3}}{a''_{3,3}} \\
\hline
c''_{2,1} & c''_{2,2} - \frac{a''_{2,3}}{a''_{3,3}}& -\frac{a''_{2,3}}{a''_{3,3}} \\
\hline
\frac{c''_{3,1}}{a''_{3,3}}  & \frac{c''_{3,2}}{a''_{3,3}} & \frac{1}{a''_{3,3}} \\
\end{array}
\right]
\left[
\begin{array}{c}
  b''_1 \\
  b''_2 \\
  b''_3 \\
\end{array}
\right]
\end{equation*}
Como resultado tenemos que la matrix $c''$ es la inversa de nuestro sistema y $b''$ la solución del sistema de ecuaciones

\textbf{Resumen}

Dado el sistema 

\begin{equation*}
\left[
\begin{array}{ccc|c|c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

\textbf{Paso 1}

Dividimos la primer ecuación entre $a_{1,1}$

\begin{equation*}
\begin{array}{c}
\frac{1}{a_{1,1}}\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 1 por $a_{2,1}$ y la restamos a la ecuación 2

\begin{equation*}
\begin{array}{c}
a_{2,1}\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
1 & \frac{a_{1,2}}{a_{1,1}} & \frac{a_{1,3}}{a_{1,1}}  & \cdots &  \frac{a_{1,M}}{a_{1,1}} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots  &  a_{2,M}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots  &  a_{3,M}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
1 & \frac{a_{1,2}}{a_{1,1}} & \frac{a_{1,3}}{a_{1,1}} & \cdots & \frac{a_{1,M}}{a_{1,1}} \\
\hline
0 & a_{2,2} - a_{2,1}\left( \frac{a_{1,2}}{a_{1,1}} \right )& a_{2,3} - a_{2,1}\left(\frac{a_{1,3}}{a_{1,1}} \right) & \cdots &  a_{2,M} - a_{2,1} \left( \frac{a_{1,M}}{a_{1,1}} \right )\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots  & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 1 por $a_{3,1}$ y la restamos a la ecuación 3

\begin{equation*}
\begin{array}{c}
a_{3,1}\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
1 & \frac{a_{1,2}}{a_{1,1}} & \frac{a_{1,3}}{a_{1,1}} & \cdots  & \frac{a_{1,M}}{a_{1,1}} \\
\hline
0 & a_{2,2} - a_{2,1}\left( \frac{a_{1,2}}{a_{1,1}} \right )& a_{2,3} - a_{2,1}\left(\frac{a_{1,3}}{a_{1,1}} \right) & \cdots &  a_{2,M} - a_{2,1} \left( \frac{a_{1,M}}{a_{1,1}} \right )\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}


\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
1 & \frac{a_{1,2}}{a_{1,1}} & \frac{a_{1,3}}{a_{1,1}} & \cdots & \frac{a_{1,M}}{a_{1,1}} \\
\hline
0 & a_{2,2} - a_{2,1}\left( \frac{a_{1,2}}{a_{1,1}} \right )& a_{2,3} - a_{2,1}\left(\frac{a_{1,3}}{a_{1,1}} \right) & \cdots &  a_{2,M} - a_{2,1} \left( \frac{a_{1,M}}{a_{1,1}} \right )\\
\hline
0   & a_{3,2} - a_{3,1} \left(\frac{a_{1,2}}{a_{1,1}}\right )  & a_{3,3} -a_{3,1} \left(\frac{a_{1,3}}{a_{1,1}} \right) & \cdots & a_{3,M}-  a_{3,1} \left(\frac{a_{1,M}}{a_{1,1}} \right )\\
\end{array}
\right]  
\end{equation*}

\textbf{Paso 2}

Dividimos la segunda ecuación entre $a_{2,2}$

\begin{equation*}
\begin{array}{c}
\\
\frac{1}{a_{2,2}}\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
\frac{a_{2,1}}{a_{2,2}}  & 1 & \frac{a_{2,3}}{a_{2,2}} & \cdots & \frac{a_{2,M}}{a_{2,2}}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 2 por $a_{1,2}$ y la restamos a la primer ecuación

\begin{equation*}
\begin{array}{c}
\\
a_{1,2}\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
\frac{a_{2,1}}{a_{2,2}}  & 1 & \frac{a_{2,3}}{a_{2,2}} & \cdots & \frac{a_{2,M}}{a_{2,2}}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,2} \left( \frac{a_{2,1}}{a_{2,2}} \right )& 0 & a_{1,3}- a_{1,2} \left(\frac{a_{2,3}}{a_{2,2}} \right) & \cdots & a_{1,M} - a_{1,2} \left(\frac{a_{2,M}}{a_{2,2}} \right )\\
\hline
\frac{a_{2,1}}{a_{2,2}}  & 1 & \frac{a_{2,3}}{a_{2,2}} & \cdots & \frac{a_{2,M}}{a_{2,2}}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 2 por $a_{3,2}$ y se la restamos a la ecuación 3

\begin{equation*}
\begin{array}{c}
\\
a_{3,2}\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,2} \left( \frac{a_{2,1}}{a_{2,2}} \right )& 0 & a_{1,3}- a_{1,2} \left(\frac{a_{2,3}}{a_{2,2}} \right) & \cdots & a_{1,M} - a_{1,2} \left(\frac{a_{2,M}}{a_{2,2}} \right )\\
\hline
\frac{a_{2,1}}{a_{2,2}}  & 1 & \frac{a_{2,3}}{a_{2,2}} & \cdots & \frac{a_{2,M}}{a_{2,2}}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}  \\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,2} \left( \frac{a_{2,1}}{a_{2,2}} \right )& 0 & a_{1,3}- a_{1,2} \left(\frac{a_{2,3}}{a_{2,2}} \right) & \cdots & a_{1,M} - a_{1,2} \left(\frac{a_{2,M}}{a_{2,2}} \right )\\
\hline
\frac{a_{2,1}}{a_{2,2}}  & 1 & \frac{a_{2,3}}{a_{2,2}} & \cdots & \frac{a_{2,M}}{a_{2,2}}\\
\hline
a_{3,1} - a_{3,2} \left( \frac{a_{2,1}}{a_{2,2}} \right )& 0 & a_{3,3}- a_{3,2} \left(\frac{a_{2,3}}{a_{2,2}} \right) & \cdots & a_{3,1} - a_{3,2} \left(\frac{a_{2,M}}{a_{2,2}} \right )\\
\end{array}
\right]  
\end{equation*}

\textbf{Paso 3}

Dividimos la tercera ecuación entre $a_{3,3}$

\begin{equation*}
\begin{array}{c}
\\
\\
\frac{1}{a_{3,3}}\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
a_{3,1} & a_{3,2} & a_{3,3} & \cdots & a_{3,M}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
\frac{a_{3,1}}{a_{3,3}} & \frac{a_{3,2}}{a_{3,3}} & 1 & \cdots & \frac{a_{3,M}}{a_{3,3}}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 3 por $a_{1,3}$ y la restamos a la primer ecuación

\begin{equation*}
\begin{array}{c}
\\
\\
a_{1,3}\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} & a_{1,2} & a_{1,3} & \cdots & a_{1,M} \\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
\frac{a_{3,1}}{a_{3,3}} & \frac{a_{3,2}}{a_{3,3}} & 1 & \cdots & \frac{a_{3,M}}{a_{3,3}}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,3} \left(\frac{a_{3,1}}{a_{3,3}}\right) & a_{1,2} - a_{1,3}\left(\frac{a_{3,2}}{a_{3,3}}\right)& 0 & \cdots & a_{1,M} - a_{1,3}\left(\frac{a_{3,M}}{a_{3,3}}\right)\\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
\frac{a_{3,1}}{a_{3,3}} & \frac{a_{3,2}}{a_{3,3}} & 1 & \cdots & \frac{a_{3,M}}{a_{3,3}}\\
\end{array}
\right]  
\end{equation*}

Multiplicamos la ecuación 3 por $a_{2,3}$ y la restamos a la segunda ecuación

\begin{equation*}
\begin{array}{c}
\\
\\
a_{2,3}\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,3} \left(\frac{a_{3,1}}{a_{3,3}}\right) & a_{1,2} - a_{1,3}\left(\frac{a_{3,2}}{a_{3,3}}\right)& 0 & \cdots & a_{1,M} - a_{1,3}\left(\frac{a_{3,M}}{a_{3,3}}\right)\\
\hline
a_{2,1} & a_{2,2} & a_{2,3} & \cdots & a_{2,M}\\
\hline
\frac{a_{3,1}}{a_{3,3}} & \frac{a_{3,2}}{a_{3,3}} & 1 & \cdots & \frac{a_{3,M}}{a_{3,3}}\\
\end{array}
\right]  
\end{equation*}

\begin{equation*}
\begin{array}{c}
\\
\\
\\
\end{array}
\left[
\begin{array}{c|c|c||c||c}
a_{1,1} - a_{1,3} \left(\frac{a_{3,1}}{a_{3,3}}\right) & a_{1,2} - a_{1,3}\left(\frac{a_{3,2}}{a_{3,3}}\right)& 0 & \cdots & a_{1,M} - a_{1,3}\left(\frac{a_{3,M}}{a_{3,3}}\right)\\
\hline
a_{2,1} - a_{2,3} \left( \frac{a_{3,1}}{a_{3,3}} \right) & a_{2,2}- a_{2,3} \left(\frac{a_{3,2}}{a_{3,3}} \right) & 0 & \cdots & a_{2,M} - a_{2,3} \left(\frac{a_{3,M}}{a_{3,3}}\right)\\
\hline
\frac{a_{3,1}}{a_{3,3}} & \frac{a_{3,2}}{a_{3,3}} & 1 & \cdots & \frac{a_{3,M}}{a_{3,3}}\\
\end{array}
\right]  
\end{equation*}

Note la similitud del método con el método de sumas y restas. No olvidar que este método es equivalente al calculo parcial de la matriz inversa. El método se ilustra mejor con un ejemplo.

\subsubsection{Ejemplo 1}

Utilice la técnica de Gauss-Jordan para resolver el siguiente sistema de ecuaciones.

\begin{eqnarray*}
30 x_1 -  x_2 - 2 x_3 & = & 78\\
 x_1 + 70 x_2  - 3 x_3 & = & -193\\
3 x_1 - 2 x_2 + 100 x_3 & = & 714
\end{eqnarray*}

Primero exprese los coeficientes y el lado derecho como una matriz aumentada

\[
\left(
\begin{array}{rrr|rrr|r}
 30 & -1 & -2 & 1 & 0 & 0 & 78 \\
 1 & 70 & -3 & 0 & 1 & 0 & -193 \\
 3 & -2 & 100 & 0 & 0 & 1 & 714
\end{array}
\right)
\]

\textbf{Primer iteración}

Luego se normaliza el primer renglón, al dividirlo entre el elemento pivote, 30, se obtiene

\[
\left[
\begin{array}{rrr|rrr|r}
 1 & -\frac{1}{30} & -\frac{1}{15} & \frac{1}{30} & 0 & 0 & \frac{13}{5} \\
 &  & & & & &  \\
 1 & 70 & -3 & 0 & 1 & 0 & -193 \\
 &&&&&&\\
 3 & -2 & 100 & 0 & 0 & 1 & 714
\end{array}
\right]
\]

El término $x_1$ se puede eliminar del segundo renglón restando el primer renglón multiplicado por $1$ del segundo renglón. En forma similar restando el primer renglón multiplicando por $3$ eliminará el término $x_1$ del tercer renglón

\[
\left[
\begin{array}{rrr|rrr|r}
 1 & -\frac{1}{30} & -\frac{1}{15} & \frac{1}{30} & 0 & 0 & \frac{13}{5} \\
 &&&&&&\\
 0 & \frac{2101}{30} & -\frac{44}{15} & -\frac{1}{30} & 1 & 0 & -\frac{978}{5} \\
 &&&&&&\\
 0 & -\frac{19}{10} & \frac{501}{5} & -\frac{1}{10} & 0 & 1 & \frac{3531}{5}
\end{array}
\right]
\]

\[
\left[
\begin{array}{rrr|rrr|r}
 1. & -0.0333 & -0.0667 & 0.0333 & 0. & 0. & 2.6 \\
 0. & 70.0333 & -2.9333 & -0.0333 & 1. & 0. & -195.6 \\
 0. & -1.9 & 100.2 & -0.1 & 0. & 1. & 706.2
\end{array}
\right]
\]

\textbf{Segunda Iteración}

En seguida, se normaliza el segundo renglón dividiéndolo entre $\frac{2101}{30}$

\[
\left[
\begin{array}{rrr|rrr|r}
 1 & -\frac{1}{30} & -\frac{1}{15} & \frac{1}{30} & 0 & 0 & \frac{13}{5} \\
 &&&&&&\\
 0 & 1 & -\frac{8}{191} & -\frac{1}{2101} & \frac{30}{2101} & 0 & -\frac{5868}{2101} \\
 &&&&&&\\
 0 & -\frac{19}{10} & \frac{501}{5} & -\frac{1}{10} & 0 & 1 & \frac{3531}{5}
\end{array}
\right]
\]

\[
\left[
\begin{array}{rrr|rrr|r}
 1. & -0.0333 & -0.0667 & 0.0333 & 0. & 0. & 2.6 \\
 0. & 1. & -0.0419 & -0.0005 & 0.0143 & 0. & -2.793 \\
 0. & -1.9 & 100.2 & -0.1 & 0. & 1. & 706.2
\end{array}
\right]
\]   

Al reducir los términos $x_2$ de las ecuaciones 1 y 2 tenemos
\[
\left[
\begin{array}{rrr|rrr|r}
 1 & 0 & -\frac{13}{191} & \frac{70}{2101} & \frac{1}{2101} & 0 & \frac{5267}{2101} \\
 &&&&&&\\
 0 & 1 & -\frac{8}{191} & -\frac{1}{2101} & \frac{30}{2101} & 0 & -\frac{5868}{2101} \\
 &&&&&&\\
 0 & 0 & \frac{19123}{191} & -\frac{212}{2101} & \frac{57}{2101} & 1 & \frac{1472577}{2101}
\end{array}
\right]
\]

\[
\left[
\begin{array}{rrr|rrr|r}
 1. & 0. & -0.0681 & 0.0333 & 0.0005 & 0. & 2.5069 \\
 0. & 1. & -0.0419 & -0.0005 & 0.0143 & 0. & -2.793 \\
 0. & 0. & 100.12 & -0.1009 & 0.0271 & 1. & 700.893
\end{array}
\right]
\]  


\textbf{Tercer iteración}

El tercer renglón se normaliza entonces al dividirlo entre $\frac{19123}{191}$

\[
\left[
\begin{array}{rrr|rrr|r}
 1 & 0 & -\frac{13}{191} & \frac{70}{2101} & \frac{1}{2101} & 0 & \frac{5267}{2101} \\
 &&&&&&\\
 0 & 1 & -\frac{8}{191} & -\frac{1}{2101} & \frac{30}{2101} & 0 & -\frac{5868}{2101} \\
 &&&&&&\\
 0 & 0 & 1 & -\frac{212}{210353} & \frac{57}{210353} & \frac{191}{19123} & \frac{1472577}{210353}
\end{array}
\right]
\]  

\[
\left[
\begin{array}{ccccccc}
 1. & 0. & -0.0681 & 0.0333 & 0.0005 & 0. & 2.5069 \\
 0. & 1. & -0.0419 & -0.0005 & 0.0143 & 0. & -2.793 \\
 0. & 0. & 1. & -0.001 & 0.0003 & 0.01 & 7.0005
\end{array}
\right]
\]

Por último, los términos $x_3$ se pueden reducir de la primera y segunda ecuación para obtener

 \[
\left[
\begin{array}{rrr|rrr|r}
 1 & 0 & 0 & \frac{538}{16181} & \frac{8}{16181} & \frac{1}{1471} & \frac{48274}{16181} \\
 &&&&&& \\
 0 & 1 & 0 & -\frac{109}{210353} & \frac{3006}{210353} & \frac{8}{19123} & -\frac{525828}{210353} \\
 &&&&&& \\
 0 & 0 & 1 & -\frac{212}{210353} & \frac{57}{210353} & \frac{191}{19123} & \frac{1472577}{210353}
\end{array}
\right]
\]  

\[
\left[
\begin{array}{ccccccc}
 1. & 0. & 0. & 0.0332 & 0.0005 & 0.0007 & 2.9834 \\
 0. & 1. & 0. & -0.0005 & 0.0143 & 0.0004 & -2.4997 \\
 0. & 0. & 1. & -0.001 & 0.0003 & 0.01 & 7.0005
\end{array}
\right]
\]

El código en Matlab es

\begin{verbatim}
function [x, Ainv] = Gauss_Jordan(A, b)
N = length(b);

B  = [A,  diag(ones(N,1)), b];

for n = 1:N;
    B(n,:) = B(n, :)./B(n,n);
    for m = 1: N
        if n ~= m
            B(m,:) = B(m,:) - B(n,:)*B(m,n);
        end;
    end;
end;
x = B(:,2*N+1);
Ainv = B(:, (N+1):(N+N));
\end{verbatim}


Como resultado tenemos en primer termino la matriz identidad resultado de las eliminaciones, en segundo lugar la matriz inversa y finalmente la solución del sistema de ecuaciones.


Para el ejemplo hacemos

\begin{verbatim}
>> [b, Ainv] = Gauss_Jordan([3, -0.1, -0.2; 0.1 7 -0.3; 0.3 -0.2 10],
 [7.85; -19.3; 71.4])

x =

    3.0000
   -2.5000
    7.0000


Ainv =

    0.3325    0.0049    0.0068
   -0.0052    0.1429    0.0042
   -0.0101    0.0027    0.0999
\end{verbatim}     

\subsubsection{Ejemplo 2}

Dado el sistema lineal de ecuaciones calcular la solución utilizando el método de Gauss-Jordan

\begin{equation*}
\left[
\begin{array}{r r r }
     1   &   1 &    2  \\
     1   &   3 &   -6  \\
     2   &  -1  &   2  \\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
\end{array}
\right]
=
\left[
\begin{array}{r }
  -1 \\
  7 \\
  0 \\
\end{array}
\right]
\end{equation*}

La matriz aumentada para este sistema es:

\begin{equation*}
\left[
\begin{array}{r r r | r }
     1   &   1 &    2  & -1\\
     1   &   3 &   -6  &  7\\
     2   &  -1  &   2  & 0 \\
\end{array}
\right] 
\end{equation*}

\textbf{Primer iteración}

\begin{equation*}
\left[
\begin{array}{r r r | r }
     1   &   1 &    2  & -1\\
     0   &   2 &   -8  &  8\\
     0   &  -3  &  -2  & 2 \\
\end{array}
\right] 
\end{equation*}

\textbf{Segunda iteración}

\begin{equation*}
\left[
\begin{array}{r r r | r }
     1   &   0 &    6  & -5\\
     0   &   1 &   -4  &  4\\
     0   &  0  &  -14  &  14 \\
\end{array}
\right] 
\end{equation*}

\textbf{Tercer iteración}

\begin{equation*}
\left[
\begin{array}{r r r | r }
     1   &   0 &  0  & 1\\
     0   &   1 &  0  &  0\\
     0   &  0  &  1  &  -1 \\
\end{array}
\right] 
\end{equation*}


\subsubsection{Ejemplo 3}

Dado el sistema lineal de ecuaciones, calcular la solución y la matriz inversa utilizando el método de Gauss-Jordan

\begin{equation*}
\left[
\begin{array}{r r r r}
    10  &  -1 &    2  &   1 \\
    -1   & 15 &   -3  &    1\\
     2   & -3  &    6  &  -3\\
     1   &  1  &   -3  &   7\\
\end{array}
\right] 
\left[
\begin{array}{r }
  x_1 \\
  x_2 \\
  x_3 \\
  x_4\\
\end{array}
\right]
=
\left[
\begin{array}{r }
  1 \\
  2 \\
  2 \\
  1\\
\end{array}
\right]
\end{equation*}

El sistema inicial aumentado es
\[
\left[
\begin{array}{rrrr|rrrr|r}
    10 &   -1   &  2  &   1  &   1  &   0  &   0   &  0   &  1\\
    -1 &   15   & -3   &  1  &   0  &   1  &   0   &  0  &   2\\
     2 &   -3   &  6   & -3   &  0   &  0   &  1   &  0  &   2\\
     1  &   1  &  -3   &  7   &  0   &  0   &  0   &  1  &   1\\
 \end{array}
\right]
\]
    
\textbf{Primer iteración}

\[
\left(
\begin{array}{rrrr|rrrr|r}
 1 & -\frac{1}{10} & \frac{1}{5} & \frac{1}{10} & \frac{1}{10} & 0 & 0 & 0 & \frac{1}{10} \\
 &&&&&&&\\
 0 & \frac{149}{10} & -\frac{14}{5} & \frac{11}{10} & \frac{1}{10} & 1 & 0 & 0 & \frac{21}{10} \\
 &&&&&&&\\
 0 & -\frac{14}{5} & \frac{28}{5} & -\frac{16}{5} & -\frac{1}{5} & 0 & 1 & 0 & \frac{9}{5} \\
 &&&&&&&\\
 0 & \frac{11}{10} & -\frac{16}{5} & \frac{69}{10} & -\frac{1}{10} & 0 & 0 & 1 & \frac{9}{10}
\end{array}
\right)
\]

\[
\left[
\begin{array}{rrrr|rrrr|r}
    1.0000 &  -0.1000  &  0.2000  &  0.1000  &  0.1000  &       0   &      0  &       0  &  0.1000\\
         0     &  14.9000  & -2.8000 &   1.1000 &   0.1000 &   1.0000    &     0  &       0  &  2.1000\\
         0     &  -2.8000  &  5.6000  & -3.2000  & -0.2000  &       0   & 1.0000   &      0   & 1.8000\\
         0     &   1.1000  & -3.2000  &  6.9000  & -0.1000   &      0    &     0  &  1.0000  &  0.9000\\
\end{array}
\right]
\]

\textbf{Segunda iteración}

\[
\left[
\begin{array}{rrrr|rrrr|r}
 1 & 0 & \frac{27}{149} & \frac{16}{149} & \frac{15}{149} & \frac{1}{149} & 0 & 0 & \frac{17}{149} \\
  &&&&&&&\\
 0 & 1 & -\frac{28}{149} & \frac{11}{149} & \frac{1}{149} & \frac{10}{149} & 0 & 0 & \frac{21}{149} \\
  &&&&&&&\\
 0 & 0 & \frac{756}{149} & -\frac{446}{149} & -\frac{27}{149} & \frac{28}{149} & 1 & 0 & \frac{327}{149} \\
  &&&&&&&\\
 0 & 0 & -\frac{446}{149} & \frac{1016}{149} & -\frac{16}{149} & -\frac{11}{149} & 0 & 1 & \frac{111}{149}
\end{array}
\right]
\]

\[
\left[
\begin{array}{rrrr|rrrr|r}
    1.0000     &    0 &   0.1812 &   0.1074  &  0.1007  &  0.0067     &    0    &     0  &  0.1141\\
         0   & 1.0000 &  -0.1879 &   0.0738  &  0.0067  &  0.0671     &    0    &     0  &  0.1409\\
         0   &      0  &  5.0738  & -2.9933 &  -0.1812  &  0.1879  &  1.0000     &    0   & 2.1946\\
         0    &     0 &  -2.9933 &  6.8188  & -0.1074 &   -0.0738   &      0  &  1.0000  &  0.7450\\
\end{array}
\right]
\]

\textbf{Tercer iteración}

\[
\left(
\begin{array}{rrrr|rrrr|r}
 1 & 0 & 0 & \frac{3}{14} & \frac{3}{28} & 0 & -\frac{1}{28} & 0 & \frac{1}{28} \\
   &&&&&&&\\
 0 & 1 & 0 & -\frac{1}{27} & 0 & \frac{2}{27} & \frac{1}{27} & 0 & \frac{2}{9} \\
   &&&&&&&\\
 0 & 0 & 1 & -\frac{223}{378} & -\frac{1}{28} & \frac{1}{27} & \frac{149}{756} & 0 & \frac{109}{252} \\
   &&&&&&&\\
 0 & 0 & 0 & \frac{955}{189} & -\frac{3}{14} & \frac{1}{27} & \frac{223}{378} & 1 & \frac{257}{126}
\end{array}
\right)
\]

\[
\left[
\begin{array}{rrrr|rrrr|r}
    1.0000    &     0     &    0  &  0.2143 &   0.1071   & 0.0000 & -0.0357   &      0   & 0.0357\\
         0   & 1.0000     &    0  & -0.0370   &      0  &  0.0741   & 0.0370      &   0  &  0.2222 \\
         0   &      0   & 1.0000  & -0.5899 &  -0.0357  &  0.0370  &  0.1971   &      0  &  0.4325\\
         0    &     0   &      0  &  5.0529  & -0.2143  &  0.0370  &  0.5899   & 1.0000  &  2.0397\\
\end{array}
\right]
\]

\textbf{Cuarta iteración}

\[
\left[
\begin{array}{rrrr|rrrr|r}
 1 & 0 & 0 & 0 & \frac{111}{955} & -\frac{3}{1910} & -\frac{58}{955} & -\frac{81}{1910} & -\frac{97}{1910} \\
    &&&&&&&\\
 0 & 1 & 0 & 0 & -\frac{3}{1910} & \frac{71}{955} & \frac{79}{1910} & \frac{7}{955} & \frac{453}{1910} \\
    &&&&&&&\\
 0 & 0 & 1 & 0 & -\frac{58}{955} & \frac{79}{1910} & \frac{254}{955} & \frac{223}{1910} & \frac{1281}{1910} \\
    &&&&&&&\\
 0 & 0 & 0 & 1 & -\frac{81}{1910} & \frac{7}{955} & \frac{223}{1910} & \frac{189}{955} & \frac{771}{1910}
\end{array}
\right]
\]

\[
\left[
\begin{array}{rrrr|rrrr|r}
    1.0000    &     0    &     0    &     0  &  0.1162 &  -0.0016 &  -0.0607 &  -0.0424 &  -0.0508\\
         0  &  1.0000    &     0    &     0  & -0.0016 &   0.0743 &   0.0414 &   0.0073 &   0.2372\\
         0   &      0  &  1.0000    &     0  & -0.0607  &  0.0414  &  0.2660  &  0.1168 &   0.6707\\
         0   &      0   &      0   & 1.0000  & -0.0424  &  0.0073  &  0.1168  &  0.1979  &  0.4037\\
\end{array}
\right]
\]

\begin{verbatim}
>> [x, Ainv] = Gauss_Jordan([10 -1 2 1; -1 15 -3 1; 2 -3 6 -3; 1 1 -3 7], 
[1,2,2,1]')

b =

   -0.0508
    0.2372
    0.6707
    0.4037


Ainv =

    0.1162   -0.0016   -0.0607   -0.0424
   -0.0016    0.0743    0.0414    0.0073
   -0.0607    0.0414    0.2660    0.1168
   -0.0424    0.0073    0.1168    0.1979
\end{verbatim}


\section{Métodos para sistemas no lineales}

Una ecuación no lineal es aquella que tiene una forma diferente a $f(x) = a_0 + a_1 x$ en cuyo caso calcular la solución consiste en resolver despejando  de $x=-a_0/a_1$. Pero el caso es que queremos resolver un sistema de ecuaciones no lineales de la forma

\begin{eqnarray*}
f_1(x_1, x_2, \cdots, x_n) & = & 0\\ 
f_2(x_1, x_2, \cdots, x_n) & = & 0\\
\vdots & = & \vdots \\
f_n(x_1, x_2, \cdots, x_n) & = & 0
\end{eqnarray*}

donde $n$ es el número de ecuaciones y se tiene el mismo números de funciones $f_i$  que variables $x_i$. En esta sección veremos dos de estos métodos que son la iteración de punto fijo y el método de Newton-Raphson, los cuales se pueden utilizar para resolver sistemas de ecuaciones no lineales


\subsection{Método Gráfico para sistemas no lineales en dos dimensiones}

En el caso de tener un sistema de ecuaciones de dos variables, es posible hacer la gráfica de las curvas de nivel. Por ejemplo en el caso de tener dos funciones como el siguiente 

\begin{eqnarray*}
f_1(x, y) & = & x^2 + x y -10 = 0\\
f_2(x, y) & = & y + 3 x y^2 - 57  =  0 \\
\end{eqnarray*}

podemos graficar la curva de nivel con la función contour de matlab. El código correspondiente es el siguiente y en la Figura \ref{fig:SNL-01} se muestra las curvas de nivel obtenidas.

\begin{verbatim}
xa = [-10:0.01:10];
ya = [-10:0.01:10];
[x,y] = meshgrid(xa, ya);

f1 = x.^2 + x.*y - 10;
f2 = y + 3*x.*y.^2-57;

contour(x,y, f1,[0,0],'k')
hold on
contour(x,y, f2,[0,0],'r')
\end{verbatim}


\subsection{Método de iteración de punto fijo para sistemas}

Considerando un sistema no lineal de ecuaciones

\begin{eqnarray*}
f_1(x_1, x_2, \cdots, x_n) & = & 0\\ 
f_2(x_1, x_2, \cdots, x_n) & = & 0\\
\vdots & = & \vdots \\
f_n(x_1, x_2, \cdots, x_n) & = & 0
\end{eqnarray*}

El método de iteración de punto fijo intentara despejar de cada una de la las ecuaciones $f_i(x_1, x_2, \cdots, x_n)$ la $i$-esima variable tal que

\begin{eqnarray*}
x_1  -  g_1(x_1, x_2, \cdots, x_n) & = & f_1(x_1, x_2, \cdots, x_n) \\ 
x_2 -  g_2(x_1, x_2, \cdots, x_n) & = & f_2(x_1, x_2, \cdots, x_n) \\
\vdots & = & \vdots \\
x_n -  g_n(x_1, x_2, \cdots, x_n) & = &f_n(x_1, x_2, \cdots, x_n) 
\end{eqnarray*}

Para resolver de manera iterativa

\begin{eqnarray*}
x_1^{(t+1)}  & = &  g_1(x_1^{(t)}, x_2^{(t)}, \cdots, x_n^{(t)})  \\ 
x_2^{(t+1)} & = &  g_2(x_1^{(t)}, x_2^{(t)}, \cdots, x_n^{(t)})  \\
\vdots & = & \vdots \\
x_n^{(t+1)} & = &  g_n(x_1^{(t)}, x_2^{(t)}, \cdots, x_n^{(t)})
\end{eqnarray*}

La implementación en Matlab es

\begin{verbatim}
function r = Punto_Fijo(g, x1) 

while 1
    x2 = g(x1)
    error = abs(norm((x2-x1)/x2));
    if(error < 0.0001) break
    else x1 = x2;
    end;
end;
r = x2;
\end{verbatim}

\subsubsection{Ejemplo 1} 

Utilice el método de iteración de punto fijo para determinar las raíces del sistema de ecuaciones dado. Considere como valores iniciales $x_1^{(0)} = 1.5$ y $x_2^{(0)}=3.5$.

\begin{eqnarray*}
x_1^2 + x_1 x_2 -10 & = & 0\\
x_2 + 3 x_1 x_2^2 - 57 & = & 0 
\end{eqnarray*}

De acuerdo con lo descrito despejamos $x_1$ y $x_2$ 

\begin{eqnarray*}
f_1(x_1, x_2) & = & x_1 -  \sqrt{10 -  x_1 x_2 }\\
f_2(x_1, x_2) & = & x_2 - \sqrt{\frac{57 - x_2}{3 x_1} }
\end{eqnarray*}

Despejando tenemos el siguiente sistema iterativo de ecuaciones

\begin{eqnarray*}
x_1^{(t+1)} & = &  \sqrt{10 -  x_1^{(t)} x_2^{(t)}  }\\
x_2^{(t+1)} & = & \sqrt{\frac{57 - x_2^{(t)}}{3 x_1^{(t+1)}} }
\end{eqnarray*}

La implementación en Matlab es

\begin{verbatim}
function x = g1(x)

x(1) = sqrt(10 - x(1)*x(2));
x(2) = sqrt((57 - x(2))/(3 * x(1)));
\end{verbatim}

El proceso iterativo será:

\textbf{Iteración 1}

Con $x^{(0)} = [1.5, 3.5]^T$ tenemos

\begin{eqnarray*}
x_1^{(1)} & = &  \sqrt{10 -  x_1^{(0)} x_2^{(0)}  } = \sqrt{10 - (1.5)(3.5)} = 2.1794\\
x_2^{(1)} & = & \sqrt{\frac{57 - x_2^{(0)}}{3 x_1^{(1)}} } = \sqrt{\frac{57 - (3.5)}{3 (2.1794)} } = 2.8605
\end{eqnarray*}

\textbf{Iteración 2}

Con $x^{(1)} = [2.1794, 2.8605]^T$ tenemos

\begin{eqnarray*}
x_1^{(2)} & = &  \sqrt{10 -  x_1^{(1)} x_2^{(1)}  } = \sqrt{10 - (2.1794)(2.8605)} = 1.9405\\
x_2^{(2)} & = & \sqrt{\frac{57 - x_2^{(0)}}{3 x_1^{(1)}} } = \sqrt{\frac{57 - (2.8605)}{3 (1.9405)} } = 3.0496
\end{eqnarray*}

El resumen el proceso iterativo se muestra en la siguiente tabla, donde podemos ver que la solución es $x^* = [2, 3]^T$ en 9 iteraciones.

\begin{center}
\begin{tabular}{|r|r|r|}
    \hline
    $k$ & $x_1^{(k)}$ & $x_2^{(k)}$ \\
    \hline
    0  &  1.5000  &  3.5000\\
    1  &  2.1794  &  2.8605\\
    2  &  1.9405  &  3.0496\\
    3  &  2.0205  &  2.9834\\
    4  &  1.9930  &  3.0057\\
    5  &  2.0024  &  2.9981\\
    6  &  1.9992  &  3.0007\\
    7  &  2.0003  &  2.9998\\
    8  &  1.9999  &  3.0001\\
    9  &   2.0000 &   3.0000\\
    \hline
\end{tabular}
\end{center}

Para realizar la ejecución dar en Matlab 

\begin{verbatim}
Punto_Fijo(@g1, [1.5; 3.5])
\end{verbatim}

\subsubsection{Ejemplo 2} 

Utilice el método de iteración de punto fijo para determinar las raíces del sistema de ecuaciones dado. Considere como valores iniciales $x_1^{(0)} = 0.0$ y $x_2^{(0)}=1.0$.

\begin{eqnarray*}
x_1^2 - 2 x_1 - x_2 + 0.5 & = & 0\\
x_1^2 + 4 x_2^2 - 4 & = & 0 
\end{eqnarray*}

De acuerdo con lo descrito despejamos $x_1$ y $x_2$ 

\begin{eqnarray*}
f_1(x_1, x_2) & = & x_1 -  \sqrt{2 x_1 + x_2 - 0.5  }\\
f_2(x_1, x_2) & = & x_2 - \sqrt{\frac{4 - x_1^2}{4} }
\end{eqnarray*}

Despejando tenemos el siguiente sistema iterativo de ecuaciones

\begin{eqnarray*}
x_1^{(t+1)} & = &  \sqrt{2 x_1^ {(t)} +  x_2^{(t)} - 0.5  }\\
x_2^{(t+1)} & = & \sqrt{\frac{4 - \left(x_1^{(t+1)} \right)^2}{4} }
\end{eqnarray*}

La implementación en Matlab es

\begin{verbatim}
function x = g2(x)

x(1) = sqrt(2*x(1) + x(2) - 0.5);
x(2) = sqrt(4 -x(1)^2)/2;  
\end{verbatim}

El proceso iterativo será:

\textbf{Iteración 1}

Con $x^{(0)} = [0,1]^T$ tenemos

\begin{eqnarray*}
x_1^{(1)} & = &  \sqrt{2 x_1^ {(0)}+  x_2^{(0)} - 0.5  } = \sqrt{2 (0) + (1) - 0.5  } = 0.7071 \\
x_2^{(1)} & = & \sqrt{\frac{4 - \left(x_1^{(1)} \right)^2}{4} }= \sqrt{\frac{4 - (0.7071)^2}{4} } = 0.9354
\end{eqnarray*}

\textbf{Iteración 2}

Con $x^{(1)} = [0.7071, 0.9354]^T$ tenemos

\begin{eqnarray*}
x_1^{(2)} & = &  \sqrt{2 x_1^ {(1)}+  x_2^{(1)} - 0.5  } = \sqrt{2 (0.7071) + (0.9354) - 0.5  } = 1.3600 \\
x_2^{(2)} & = & \sqrt{\frac{4 - \left(x_1^{(2)} \right)^2}{4} }= \sqrt{\frac{4 - (1.3600)^2}{4} } = 0.7332
\end{eqnarray*}

El resumen el proceso iterativo se muestra en la siguiente tabla, donde podemos ver que la solución es $x^* = [1.9007, 0.3112]^T$ en 8 iteraciones.

\begin{center}
\begin{tabular}{|r|r|r|}
    \hline
    $k$ & $x_1^{(k)}$ & $x_2^{(k)}$ \\
    \hline
    0  &  0.0000  & 1.0000\\
    1  &  0.7071  &  0.9354\\
    2  & 1.3600   &  0.7332\\
    3  & 1.7185  &   0.5116\\
    4  & 1.8570  &   0.3713\\
    5  & 1.8935  &  0.3220\\
    6  & 1.8997  &  0.3127\\
    7  & 1.9006  &  0.3114 \\
    8  & 1.9007  &  0.3112\\
    \hline
\end{tabular}
\end{center}

Para realizar la ejecución dar en Matlab 

\begin{verbatim}
Punto_Fijo(@g2, [0; 1])
\end{verbatim}

\subsection{Método de Newton-Raphson para sistemas}

Consideremos el sistema de ecuaciones no lineales

\begin{eqnarray*}
f_1(x_1, x_2, \cdots, x_n) & = & 0\\ 
f_2(x_1, x_2, \cdots, x_n) & = & 0\\
\vdots & = & \vdots \\
f_n(x_1, x_2, \cdots, x_n) & = & 0
\end{eqnarray*}

Utilizando la serie de Taylor podemos hacer una aproximación lineal para una función $f_i(x_1^{(t)} +\delta x_1^{(t)}, x_2^{(t)} +\delta x_2^{(t)}, \cdots, x_n^{(t)} + \delta x_n^{(t)})$ en un incremento $\delta x_i^{(t)} = x_i^{(t+1)}- x_i^{(t)}$ como:

\begin{eqnarray*}
f_1(x_1^{(t)} +\delta x_1^{(t)}, x_2^{(t)} +\delta x_2^{(t)}, \cdots, x_n^{(t)} + \delta x_n^{(t)}) \\
= f_1(x_1^{(t)}, x_2^{(t)} + \cdots, + x_n^{(t)})  + \frac{\partial f_1}{\partial x_1} \delta x_1^{(t)} + \frac{\partial f_1}{\partial x_2} \delta x_2^{(t)} + \cdots + \frac{\partial f_1}{\partial x_n} \delta x_n^{(t)}\\
f_2(x_1^{(t)} +\delta x_1^{(t)}, x_2^{(t)} +\delta x_2^{(t)}, \cdots, x_n^{(t)} + \delta x_n^{(t)}) \\
= f_2(x_1^{(t)}, x_2^{(t)} + \cdots, + x_n^{(t)})  + \frac{\partial f_2}{\partial x_1} \delta x_1^{(t)} + \frac{\partial f_2}{\partial x_2} \delta x_2^{(t)} + \cdots + \frac{\partial f_2}{\partial x_n} \delta x_n^{(t)}\\
\vdots \\
f_n(x_1^{(t)} +\delta x_1^{(t)}, x_2^{(t)} +\delta x_2^{(t)}, \cdots, x_n^{(t)} + \delta x_n^{(t)}) \\
= f_n(x_1^{(t)}, x_2^{(t)} + \cdots, + x_n^{(t)})  + \frac{\partial f_n}{\partial x_1} \delta x_1^{(t)} + \frac{\partial f_n}{\partial x_2} \delta x_2^{(t)} + \cdots + \frac{\partial f_n}{\partial x_n} \delta x_n^{(t)}\\
\end{eqnarray*}

Si escribimos el sistema en forma matricial tenemos

\begin{eqnarray*}
\left[
\begin{array}{c}
f_1(x_1^{(t)} + \delta x_1^{(t)}, x_2^{(t)}+ \delta x_2^{(t)}, \cdots,  x_n^{(t)}+ \delta x_n^{(t)} )  \\
f_2(x_1^{(t)} + \delta x_1^{(t)}, x_2^{(t)}+ \delta x_2^{(t)}, \cdots,  x_n^{(t)}+ \delta x_n^{(t)} )  \\
\vdots\\ 
f_n(x_1^{(t)} + \delta x_1^{(t)}, x_2^{(t)}+ \delta x_2^{(t)}, \cdots,  x_n^{(t)}+ \delta x_n^{(t)} )  \\
\end{array}
\right] = \\
\left[
\begin{array}{c}
f_1(x_1^{(t)}, x_2^{(t)}, \cdots,  x_n^{(t)} )  \\
f_2(x_1^{(t)}, x_2^{(t)}, \cdots,  x_n^{(t)} )  \\
\vdots\\ 
f_n(x_1^{(t)}, x_2^{(t)}, \cdots,  x_n^{(t)} )  \\
\end{array}
\right] +
\left[
\begin{array}{cccc}
\frac{\partial f_1(x)}{\partial x_1}  &  \frac{\partial f_1(x)}{\partial x_2}  & \cdots & \frac{\partial f_1(x)}{\partial x_n}   \\
\frac{\partial f_2(x)}{\partial x_1}  &  \frac{\partial f_2(x)}{\partial x_2}  & \cdots & \frac{\partial f_2(x)}{\partial x_n}   \\
\vdots  & \vdots  & \cdots & \vdots \\
\frac{\partial f_n(x)}{\partial x_1}  &  \frac{\partial f_n(x)}{\partial x_2}  & \cdots & \frac{\partial f_n(x)}{\partial x_n}   \\  
\end{array}
\right]
\left[
\begin{array}{c}
\delta x_1^{(t)} \\
\delta x_2^{(t)} \\
\vdots \\
\delta x_n^{(t)} \\
\end{array}
\right]
\end{eqnarray*}

En forma compacta

\begin{eqnarray*}
f(x^{(t)} + \delta x^{(t)}) &=& f(x^{(t)}) + J(x^{(t)}) \delta x^{(t)} \\
f(x^{(t)} + \delta x^{(t)}) &=& f(x^{(t)}) + J(x^{(t)}) (x^{(t+1)} - x^{(t)} ) \\
\end{eqnarray*}


donde $J(x)$ es la matriz de primeras derivada o Jacobiano.

Como queremos encontrar el cero de la función,  la aproximación lineal que debemos resolver es:

\begin{equation*}
0 = f(x^{(t)}) + J(x^{(t)})  (x^{(t+1)} - x^{(t)} )\\
\end{equation*} 

donde las actualizaciones las hacemos como

\begin{equation*}
x^{(t+1)}  =  x^{(t)}  -  \left[J(x^{(t)})\right]^{-1} f(x^{(t)})
\end{equation*} 

La implementación del algoritmo en Matlab es:

\begin{verbatim}
function r = Newton_Raphson(f, J, x1)

while 1
    x2 = x1 - inv(J(x1))* f(x1);
    if norm((x1-x2)./x2) < 0.0001 break;
    else x1 = x2;
    end;
end;
r = x2;
\end{verbatim}

\subsubsection{Ejemplo 1}

Resolver el siguiente sistema de ecuaciones dado por y cuya gráfica se muestra en la Figura \ref{fig:SNL-01}

\begin{figure}[!htp]
   \centering
   \centerline{\psfig{figure=Capitulo_4/figuras/SNL-01.eps,width=8.0cm}}
   \caption{Curvas con la solución del sistema del ejemplo 1}
   \label{fig:SNL-01}
\end{figure}

\begin{eqnarray*}
f_1(x_1,x_2) &=& x_1^2 + x_1 x_2 -10\\
f_2(x_1,x_2) &= & x_2 + 3 x_1 x_2^2 -57
\end{eqnarray*}

El Jacobiano es
\[
J(x) = 
\left[
\begin{array}{c|c}
2 x_1 + x_2  &  x_1   \\
\hline
3 x_2^2    &   1 + 6 x_1 x_2    
\end{array}
\right]
\]
y el arreglo de funciones

\[
f(x) = 
\left[
\begin{array}{c}
x_1^2 + x_1 x_2 -10\\
x_2 + 3 x_1 x_2^2 -57
\end{array}
\right]
\]

Considerando como valores iniciales $x^{(0)} = [1.5, 3.5]^T$ tenemos:
 
\textbf{Primer iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(1)}    \\
x_2^{(1)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
    1.5000 \\
    3.5000     
\end{array}
\right] -
\left[
\begin{array}{rr}
    6.5000  &  1.5000 \\
   36.7500 &  32.5000
\end{array}
\right]^{-1}
\left[
\begin{array}{l}
   -2.5000 \\
    1.6250   
\end{array}
\right] =
\left[
\begin{array}{r}
    2.0360 \\
    2.8439
\end{array}
\right]
\end{equation*}


\textbf{Segunda iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(2)}    \\
 x_2^{(2)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
    2.0360 \\
    2.8439     
\end{array}
\right] -
\left[
\begin{array}{rr}
    6.9159  &  2.0360 \\
   24.2629 &  35.7413\\
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
   -0.0644 \\
   -4.7562 
\end{array}
\right] =
\left[
\begin{array}{r}
    1.9987\\
    3.0023
\end{array}
\right]
\end{equation*}

\textbf{Tercer iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(3)}    \\
x_2^{(3)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
    1.9987 \\
    3.0023     
\end{array}
\right] -
\left[
\begin{array}{rr}
    6.9997  &  1.9987 \\
   27.0412  & 37.0041
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
   -0.0045 \\
    0.0496
\end{array}
\right] =
\left[
\begin{array}{r}
    2.0000 \\
    3.0000
\end{array}
\right]
\end{equation*}

\textbf{Cuarta iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(4)}    \\
 x_2^{(4)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
    2.0000\\
    3.0000    
\end{array}
\right] -
\left[
\begin{array}{rr}
    7.0000  &  2.0000 \\
   27.0000 &  37.0000\\
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
   -0.00000129\\
   -0.00002214
\end{array}
\right] =
\left[
\begin{array}{r}
    2.0000 \\
    3.0000
\end{array}
\right]
\end{equation*}

Note que solamente 4 iteraciones son necesarias para llegar a la solución $x = [2,3]^T$

La implementación en Matlab para este ejemplo son:

Para la función tenemos

\begin{verbatim}
function f = f1(x)

n = length(x);
f = zeros(n,1);

f(1) = x(1)^2 + x(1)*x(2) - 10; 
f(2) = x(2) + 3*x(1)*x(2)^2 - 57;
end
\end{verbatim}

El Jacobiano

\begin{verbatim}
function J = J1(x)
n = length(x);
J = zeros(n, n);

J(1,1) = 2*x(1) + x(2);
J(1,2) = x(1);
J(2,1) = 3*x(2)^2;
J(2,2) = 1 + 6*x(1)*x(2);
end
\end{verbatim}

y la ejecución

\begin{verbatim}
>> Newton_Raphson(@f1, @J1, [1.5,3.5]')

ans =

    2.0000
    3.0000

\end{verbatim}

\subsubsection{Ejemplo 2}

Resolver el siguiente sistema de ecuaciones dado por las ecuaciones y cuya solución gráfica se muestra en la figura \ref{fig:SNL-02}

\begin{figure}[!htp]
   \centering
   \centerline{\psfig{figure=Capitulo_4/figuras/SNL-02.eps,width=8.0cm}}
   \caption{Curvas con la solución del sistema del ejemplo 2}
   \label{fig:SNL-02}
\end{figure}


\begin{eqnarray*}
f_1(x_1,x_2) &=& x_1^2 – 2x_1 – x_2 + 0.5\\
f_2(x_1,x_2) &= & x_1^2 + 4x_2^2 – 4
\end{eqnarray*}

El Jacobiano es
\[
J(x) = 
\left[
\begin{array}{c|c}
2 x_1 -2  &  -1   \\
\hline
2 x_1    &   8x_2    
\end{array}
\right]
\]
y el arreglo de funciones

\[
f(x) = 
\left[
\begin{array}{c}
x_1^2 – 2x_1 – x_2 + 0.5\\
x_1^2 + 4x_2^2 – 4
\end{array}
\right]
\]

Considerando como valores iniciales $x^{(0)} = [0, 1]^T$ tenemos:
 
\textbf{Primer iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(1)}    \\
x_2^{(1)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
0    \\
1     
\end{array}
\right] -
\left[
\begin{array}{rr}
-2  & -1     \\
 0  &  8  
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
-0.5    \\
0.0   
\end{array}
\right] =
\left[
\begin{array}{r}
-0.25    \\
1.00  
\end{array}
\right]
\end{equation*}


\textbf{Segunda iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(2)}    \\
x_2^{(2)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
-0.2500    \\
1.0000      
\end{array}
\right] -
\left[
\begin{array}{rr}
2 .5000 & -1.0000     \\
-0.5000  &  8.0000  
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
0.0625    \\
0.0625  
\end{array}
\right] =
\left[
\begin{array}{r}
-0.2226   \\
 0.9939  
\end{array}
\right]
\end{equation*}

\textbf{Tercer iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(3)}    \\
x_2^{(3)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
   -0.2226\\
    0.9939    
\end{array}
\right] -
\left[
\begin{array}{rr}
   -2.4451  & -1.0000 \\
   -0.4451  &  7.9512\\
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
   -0.0008\\
    0.0009 
\end{array}
\right] =
\left[
\begin{array}{r}
   -0.2222 \\
    0.9938
\end{array}
\right]
\end{equation*}


\textbf{Cuarta iteración}

\begin{equation*}
\left[
\begin{array}{r}
x_1^{(4)}    \\
x_2^{(4)}     
\end{array}
\right] = 
\left[
\begin{array}{r}
   -0.2222 \\
    0.9938   
\end{array}
\right] -
\left[
\begin{array}{rr}
   -2.4444  &  -1.0000\\
   -0.4444   &  7.9505
\end{array}
\right]^{-1}
\left[
\begin{array}{r}
   -0.00002300\\
    0.00000038
\end{array}
\right] =
\left[
\begin{array}{r}
   -0.2222 \\
    0.9938
\end{array}
\right]
\end{equation*}

La implementación en Matlab para este ejemplo son:


Para la función tenemos

\begin{verbatim}
function f = f2(x) 

n = length(x);
f = zeros(n,1);

f(1) = x(1)*x(1) - 2*x(1) - x(2) + 0.5; 
f(2) = x(1)*x(1) + 4*x(2)*x(2) - 4;
end
\end{verbatim}

El Jacobiano

\begin{verbatim}
function y = J2(x)
n = length(x);
y = zeros(n, n);
y(1,1) = 2*x(1) -2;
y(1,2) = -2;
y(2,1) = 2*x(1);
y(2,2) = 8*x(2);
end
\end{verbatim}

y la ejecución

\begin{verbatim}
>> Newton_Raphson(@f2, @J2, x)

ans =

   -0.2222
    0.9938

\end{verbatim}

\subsection{Ejemplo}

Para el circuito que se muestra en la figura \ref{fig:sistema}, esta constituido por dos mallas y dos cargas. Se sea calcular la corriente que circula por cada uno de los elementos y la potencia que debe suministrar la fuente de voltaje.

\begin{figure}[!htp]
   \centering
   \centerline{\psfig{figure=Capitulo_4/figuras/sistema.eps,width=8.0cm}}
   \caption{Circuito eléctrico de dos mallas}
   \label{fig:sistema}
\end{figure}

La ecuaciones de Voltaje para cada una de las mallas es:

\begin{eqnarray*}
V - R_1 (i_1 + i_2) - L_1/i_1 = 0\\
V - R_1 (i_1 + i_2) - R_2 i_2 - L_2/i_2 =0
\end{eqnarray*}

Reorganizando las funciones tenemos

\begin{equation*}
F = \left[
\begin{array}{c}
V i_1 - R_1 (i_1 + i_2) i_1 - L_1\\
V i_2 - R_1 (i_1 + i_2) i_2- R_2 i_2^2 - L_2
\end{array}
\right] = 0
\end{equation*}

\begin{verbatim}
function f = Fun(i)
N = length(i);
f = zeros(N,1);

V = 100;
R1 = 2;
R2 = 3;
L1 = 100;
L2 = 200;

f(1) = V*i(1) - R1*(i(1)+i(2))*i(1) - L1;
f(2) = V*i(2) - R1*(i(1)+i(2))*i(2) - R2*i(2)^2 - L2;
end
\end{verbatim}

El Jacobiano para este sistema es:

\begin{equation*}
J = \left[
\begin{array}{c|c}
V  - R_1 (2 i_1 + i_2) \  & - i_1 R_1\\
\hline
-i_2 R1 & V_1  - R_1 (i_1 + 2 i_2)  - 2 R_2 i_2   \\
\end{array}
\right]
\end{equation*}

\begin{verbatim}
function J = JacSistema(i)
V = 100;
R1 = 2;
R2 = 3;
%L1 = 100;
%L2 = 200;

J(1,1) = V - R1*(2*i(1) + i(2)) ;
J(1,2) = -i(1)*R1;
J(2,1) = -i(2)*R1;
J(2,2) = V - (i(1) + 2*i(2))*R1 - 2*i(2)*R2 ;
end
\end{verbatim}

La solución considerando un valor inicial de corrientes $I^{(0)} = [0,0]^T$:

\begin{verbatim}

>> I = NewtonRaphson(@FunSistema, @JacSistema, [0,0]')
1.-      1     2

2.-     1.0720    2.3114

3.-     1.0728    2.3185

4.-     1.0728    2.3185


I =

    1.0728
    2.3185
\end{verbatim}

Potencia de la fuente $P_V =  V(i_1+i_2) = 100(1.0728 +  2.3185) = + 339.1282$\\ 
Potencia consumida por $R_1$  es $P_{R_1} = -R_1(i_1+i_2)^2 = -2(1.0728 +  2.3185)^2  = - 23.0016$\\
Potencia consumida por $R_2$  es $P_{R_1} = -R_2i_2^2 = - 3 (2.3185) ^2 = - 16.1266$\\
La potencia consumida por las cargas es $L1 + L2 = -300$.\\
La suma de la potencia es cero.
 
\subsubsection{Desventajas del Método de Newton}

Aunque el método de Newton-Raphson en general es muy eficiente, hay situaciones en que se comporta en forma deficiente. Un caso especial, raíces múltiples.




\subsubsection{Ejemplo}

Determinar la raíz de la función $f(x) = x^{10} –1$.

La solución utilizando el método de Newton-Raphson queda:

\begin{equation*}
x^{(k+1)} = x^{(k)} - \frac{f(x^{(k)})} {f'(x^{(k+1)})}
\end{equation*}

Sustituyendo valores tenemos

\begin{equation*}
x^{(k+1)} = x^{(k)} - \frac{\left(x^{(k)} \right)^{10}-1} {10 \left(x^{(k)} \right)^9}
\end{equation*}

Y la solución numérica es:

\begin{center}
\begin{tabular}{|r|r|r|}
\hline
$x^{(k)}$ & $f(x^{(k)})$ & $f'(x^{(k)})$ \\
\hline
0.5000 & -0.9990 & 0.0195 \\
51.6500 & 135114904483914000.0000 & 26159710451871000.0000 \\
46.4850 & 47111654129711500.0000 & 10134807815362300.0000 \\
41.8365 & 16426818072478500.0000 & 3926432199748670.0000 \\
37.6529 & 5727677301318310.0000 & 1521180282851980.0000 \\
33.8876 & 1997117586819850.0000 & 589336409039672.0000\\
30.4988 & 696351844868619.0000 & 228320999775654.0000\\
27.4489 & 242802875029547.0000 & 88456233382052.8000\\
24.7040 & 84660127717097.5000 & 34269757191973.2000\\
22.2336 & 29519161271064.1000 & 13276806089225.7000\\
20.0103 &10292695105054.7000 & 5143706707446.1600\\
18.0092 & 3588840873655.1100 & 1992777367871.5700\\
16.2083 & 1251351437592.9200 & 772042782329.1500\\
14.5875 & 436319267276.5290 & 299105192259.1190\\
13.1287 & 152135121499.2910 & 115879479847.7330\\
11.8159 & 53046236848.5329 & 44894084747.9692\\
10.6343 & 18496079117.2577 & 17392888266.5936\\
9.5708 & 6449184014.3077 & 6738361277.7304\\
8.6138 & 2248691421.7628 & 2610579221.6818\\
7.7524 & 784070216.9426 & 1011391879.0870\\
\hline
\end{tabular}
\end{center}
 

