Trabajo realizado con el apoyo del Programa UNAM-DGAPA-PAPIME PE101019

• Autor:
• Miguel Angel Pérez León
• Rev: jue jun 19 17:08:27 CDT 2020

Introducción

Múltiples problemas de las ciencias exactas implican la resolución de sistemas de N ecuaciones no lineales con N incógnitas.

Un método es linealizar y resolver, repetidamente. Esta es la misma estrategia que usa el método de Newton para resolver una sola ecuación no lineal.

Así que es natural pensar en una forma de extender esta forma de resolver ecuaciones no lineales, de manera que sea posible resolver un sistema de ecuaciones no lineales.

En el caso general, un sistema de N ecuaciones no lineales con N incógnitas $x_{i}$ se puede presentar en la forma:

$$Sistema = \begin{cases} f_{1}(x_{1},x_{2},\ldots,x_{n}) & =0\\ f_{2}(x_{1},x_{2},\ldots,x_{n}) & =0 \\ \vdots \\ f_{n}(x_{1},x_{2},\ldots,x_{n}) & =0 \end{cases}$$

En caso de que el sistema fuera de ecuaciones lineales, podemos representar este sistema de N ecuaciones con N incognitas mediante su forma matricial. $$A\vec{x}=\vec{0}$$

Con $A\in M_{n\times n}$ y $\vec{x}, \vec{b} \in \mathbb{R}^{n}$, es decir: $$\left(\begin{array}{ccccccc} a_{11} & a_{12} & \cdots & \cdots & \cdots & \cdots & a_{1n}\\ a_{21} & a_{22} & \cdots & \cdots & \cdots & \cdots & a_{2n}\\ \vdots & \ddots & \ddots & \ddots & \ddots & \ddots & \vdots\\ \vdots & \ddots & \ddots & \ddots & \ddots & \ddots & \vdots\\ \vdots & \ddots & \ddots & \ddots & \ddots & \ddots & \vdots\\ a_{n1} & \cdots & \cdots & \cdots & \cdots & \cdots & a_{nn} \end{array}\right)\left(\begin{array}{c} x_{1}\\ x_{2}\\ \vdots\\ \vdots\\ \vdots\\ x_{n} \end{array}\right)=\left(\begin{array}{c} 0\\ 0\\ \vdots\\ \vdots\\ \vdots\\ 0 \end{array}\right)$$

Usando notación vectorial, podemos reescribir el sistema en una forma más elegante: $$ F(\vec{X})=\vec{0}$$

Definiendo vectores columna como: $$F = \left[f_{1},f_{2},\ldots,f_{N}\right]^{T}$$$$ \vec{X} = \left[x_{1,}x_{2},\ldots,x_{N}\right]$$$$ \vec{0} = \left[x_{1,}x_{2},\ldots,x_{N}\right]^{T}$$

Si recordamos la expresión del método de Newton para ecuaciones no lineales es. $$ x_{n+1}=x_{n}-\frac{f\left(x_{n}\right)}{f^{'}\left(x_{n}\right)}$$

No es muy difícil extender la expresión anterior a sistemas de ecuaciones no lineales, lo que nos da como resultado la forma iterativa para resolver sistemas de ecuaciones no lineales mediante el metodo de Newton. $$\vec{X}^{\left(k+1\right)}=\vec{X}^{\left(k\right)}-\left[F'\left(\vec{X}^{(k)}\right)\right]^{-1}F\left(\vec{X}^{(k)}\right)$$

Desarrollo

El desarrollo para llegar al método de Newton para sistemas de ecuaciones lineales se mostrará a continuación y para ello tomaremos las siguientes 3 ecuaciones no lineales $$Sistema =\begin{cases} f_{1}(x_{1},x_{2},x_{3}) & =0\\ f_{2}(x_{1},x_{2},x_{3}) & =0 \quad \tag{1}\\ f_{3}(x_{1},x_{2},x_{3}) & =0 \end{cases}$$

Si se aplica el desarrollo del polinomio de Taylor centrado en $\vec{0}$ para 3 variables $i=1,2,3$, se tiene una solución aproximada para cada $f_{i}$. $$P_{i}(x_{1},x_{2},x_{3}) \thickapprox f_{i}(0,0,0)+\frac{\delta f_{i}}{\delta x_{1}}x_{1}+\frac{\delta f_{i}}{\delta x_{2}}x_{2}+\frac{\delta f_{i}}{\delta x_{3}}x_{3}$$

Supongamos que el vector $\vec{X}^{(0)}=(x_{1}^{(0)},x_{2}^{(0)},x_{3}^{(0)})$ es una solución aproximada de (1). Sea $\vec{H}= \left(h_{1},h_{2},h_{3}\right)$ una corrección de la solución inicial, de modo que: $$\vec{X}^{(1)}=\vec{X}^{(0)}+\vec{H}=(x_{1}^{(0)}+h_{1},x_{2}^{(0)}+h_{2},x_{3}^{(0)}+h_{3})$$

La aproximación $\vec{X}^{(1)}$ es una mejor aproximación y además: $$\vec{H}=\vec{X}^{(1)}-\vec{X}^{(0)}=(x_{1}^{(1)}-x_{1}^{(0)},x_{2}^{(1)}-x_{2}^{(0)},x_{3}^{(1)}-x_{3}^{(0)})$$

Descartando los términos de orden superior en el desarrollo de Taylor, una mejor aproximación esta dada por: $$f_{i}(x_{1}^{(1)},x_{2}^{(1)},x_{3}^{(1)})=f_{i}(x_{1}^{(0)},x_{2}^{(0)},x_{3}^{(0)})+\frac{\delta f_{i}}{\delta x_{1}}(x_{1}^{(1)}-x_{1}^{(0)})+\frac{\delta f_{i}}{\delta x_{2}}(x_{2}^{(1)}-x_{2}^{(0)})+\frac{\delta f_{i}}{\delta x_{3}}(x_{3}^{(1)}-x_{3}^{(0)})$$

Reescribiendo la ecuación anterior de manera más compacta. $$f_{i}(x_{1}+h_{1},x_{2}+h_{2},x_{3}+h_{1})=f_{i}(x_{1},x_{2},x_{3})+\frac{\delta f_{i}}{\delta x_{1}}h_{1}+\frac{\delta f_{i}}{\delta x_{2}}h_{2}+\frac{\delta f_{i}}{\delta x_{3}}h_{3}\tag{2} $$

Ya que en $(1)$ el sistema esta igualado a $\vec{0}$, entonces: $$f_{i}(x_{1}+h_{1},x_{2}+h_{2},x_{3}+h_{1})\thickapprox 0$$

En notación vectorial, se tiene: $$\vec{0}\thickapprox F\left(\vec{X}^{(0)}+\vec{H}\right)\thickapprox F\left(\vec{X}^{(0)}\right)+F'\left(\vec{X}^{(0)}\right)\vec{H} \tag{3}$$

Donde la matriz Jacobiana ó Jacobiano, está definida por $$F'\left(\vec{X}^{(0)}\right)=\left[\begin{array}{ccc} \frac{\delta f_{1}}{\delta x_{1}} & \frac{\delta f_{1}}{\delta x_{2}} & \frac{\delta f_{1}}{\delta x_{3}}\\ \frac{\delta f_{2}}{\delta x_{1}} & \frac{\delta f_{2}}{\delta x_{2}} & \frac{\delta f_{2}}{\delta x_{3}}\\ \frac{\delta f_{3}}{\delta x_{1}} & \frac{\delta f_{3}}{\delta x_{2}} & \frac{\delta f_{3}}{\delta x_{3}} \end{array}\right] $$

En la primera iteración todas las derivadas parciales se evalúan en $\vec{X}^{(0)}$, es decir: $$ \frac{\delta f_{i}}{x_{j}}=\frac{\delta f_{i}\left(\vec{X}^{(0)}\right)}{x_{j}}$$

También suponemos que la matriz Jacobiana es no singular, por lo que su inversa existe. Por lo tanto resolviendo para H en (3), se tiene: $$\vec{0}\thickapprox F\left(\vec{X}^{(0)}\right)+F'\left(\vec{X}^{(0)}\right)\vec{H}$$$$ \Rightarrow F'\left(\vec{X}^{(0)}\right)\vec{H} \thickapprox -F\left(\vec{X}^{(0)}\right) $$$$ \Rightarrow \left[F'\left(\vec{X}^{(0)}\right)\right]^{-1}F'\left(\vec{X}^{(0)}\right)\vec{H} \thickapprox-\left[F'\left(\vec{X}^{(0)}\right)\right]^{-1}F\left(\vec{X}^{(0)}\right) $$$$\Rightarrow \vec{H}\thickapprox-\left[F'\left(\vec{X}^{(0)}\right)\right]^{-1}F\left(\vec{X}^{(0)}\right) $$

De tal manera que $\vec{H}$ nos dice cuánto se tiene que modificar la solución anterior para encontrar una nueva mejor solución y el término: $$\left[F'\left(\vec{X}^{(0)}\right)\right]^{-1}$$

Es la inversa de la matriz Jacobiaba. En general el método de Newton se ve de la siguiente manera: $$ \vec{X}^{(k+1)}=\vec{X}^{(k)}-\left[F'\left(\vec{X}^{(k)}\right)\right]^{-1}F\left(\vec{X}^{(k)}\right) $$

En la práctica el método de Newton no implica invertir la matriz Jacobiana, solo resolver los sistemas lineales jacobianos. $$ \left[F'\left(\vec{X}^{(k)}\right)\right]\vec{H}^{(k)}=-F\left(\vec{X}^{(k)}\right) \tag{4}$$

Podemos asegurar que $(4)$ representa un sistema lineal, ya que la incógnita es $\vec{H}$ y los coeficientes de la matriz jacobiana no dependen de $\vec{H}$.

Esto significa que para resolver un sistema de ecuaciones no lineales necesitamos resolver un sistema de ecuaciones lineal, es decir, que para resolver un problema que parece complejo, se resuelve un problema equivalente pero menos complejo.

Ejemplo

Emplea el método de Newton para resolver el siguiente sistema de ecuaciones con una tolerancia de $10^{-8}$. $$2x_{1}-x_{2}-e^{-x_{1}}=0, \quad -x_{1}+2x_{2}-e^{-x_{2}}=0$$

Solución

Este sistema de ecuaciones se puede expresar en notación vectorial de la siguiente manera: $$F\left(\vec{X}\right)=\begin{cases} 2x_{1}-x_{2}-e^{-x_{1}} & =0\\ -x_{1}+2x_{2}-e^{-x_{2}} & =0 \end{cases}$$

De tal manera que la matriz jacobina asociada a este sistema está definida por: $$F'(\vec{X})=\begin{cases} 2+e^{-x_{1}} & -1\\ -1 & 2+e^{-x_{2}} \end{cases}$$

Para calcular la inversa de $F'(\vec{X})$ es necesario calcular el determinante, el cual está dado por: $$det(F'(\vec{X}))=\left(2+e^{-x_{1}}\right)\left(2+e^{-x_{2}}\right)-\left(-1\times-1\right)=4+2e^{-x_{2}}+2e^{-x_{1}}+e^{-x_{1}-x_{2}}-1$$ Y la inversa de $F'(\vec{X})$ se define como:

$$'(\vec{X})^{-1}=\frac{1}{det(F'(\vec{X}))}\left(\begin{array}{cc} 2+e^{-x_{2}} & 1\\ 1 & 2+e^{-x_{1}} \end{array}\right)\\=\left(\begin{array}{cc} \frac{2+e^{-x_{2}}}{4+2e^{-x_{2}}+2e^{-x_{1}}+e^{-x_{1}-x_{2}}-1} & \frac{1}{4+2e^{-x_{2}}+2e^{-x_{1}}+e^{-x_{1}-x_{2}}-1}\\ \frac{1}{4+2e^{-x_{2}}+2e^{-x_{1}}+e^{-x_{1}-x_{2}}-1} & \frac{2+e^{-x_{1}}}{4+2e^{-x_{2}}+2e^{-x_{1}}+e^{-x_{1}-x_{2}}-1} \end{array}\right)$$

Ahora solo nos hace falta una primera aproximación $$ \vec{X}^{(0)}=(0,0) $$

No forzosamente tiene que ser el vector $\vec{0}$, pero este método no tiene restricciones al respecto.

De tal manera que ahora podemos aplicar la definición del método de Newton para sistemas de ecuaciones y resolver para $$\vec{X}^{(1)}=\vec{X}^{(0)}-\left[F'\left(\vec{X}^{(0)}\right)\right]^{-1}F\left(\vec{X}^{(0)}\right)\tag{5}$$

Al evaluar$$\vec{X}^{(0)}=(0,0)$$

en (5), se tiene $$\vec{X}^{(1)}=\left(0,0\right)-\left[F'\left(0,0\right)\right]^{-1}F\left(0,0\right)$$$$=\left(0,0\right)-\left(\begin{array}{cc} \frac{3}{8} & \frac{1}{8}\\ \frac{1}{8} & \frac{3}{8} \end{array}\right) \left(\begin{array}{c} -1\\ -1 \end{array}\right) $$$$ =\left(0,0\right)-\left(\begin{array}{c} \left(\frac{3}{8}\cdot-1\right)+\left(\frac{1}{8}\cdot-1\right)\\ \left(\frac{1}{8}\cdot-1\right)+\left(\frac{3}{8}\cdot-1\right) \end{array}\right) $$$$=\left(0,0\right)-\left(\begin{array}{c} -\frac{1}{2}\\ -\frac{1}{2} \end{array}\right) $$$$\vec{X}^{(1)}=\left(\begin{array}{c} \frac{1}{2}\\ \frac{1}{2} \end{array}\right)$$

Después de $k$ iteraciones con el método de Newton podemos ver que la solución converge al valor $$\vec{X}^{(k)}=\left(\begin{array}{c} 5.6714329040978384\times10^{-1}\\ 5.6714329040978384\times10^{-1} \end{array}\right)$$

Lo que significa que si evaluamos $F\left(\vec{X}^{(k)}\right)$ obtendremos un vector muy cercano al vector $\vec{0}$ $$F\left(\vec{X}^{(k)}\right)=\left(\begin{array}{c} 0.000000001\\ 0.000000001 \end{array}\right)$$

Refeencias

• Butt, R. (2009). Introduction to Numerical Analysis Using MATLAB®. Jones & Bartlett Learning.
• Cheney, W., & Kincaid, D. (2010). Métodos numéricos y computación (6a. Ed.). Cengage Learning Editores S.A. de C.V.
• Burden, R. L., Faires, J. D., & C, S. M. (1985). Análisis numérico. Grupo Editorial Iberoamérica.
• Skiba, Y. N. (2001). Introducción a los métodos numéricos. UNAM, Dirección General de Publicaciones y Fomento Editorial.