PDF de programación - GUIA MATLAB - SOLUCION DE ECUACIONES NO LINEALES Y SISTEMAS LINEALES

SOLUCION DE ECUACIONES NO LINEALES Y SISTEMAS

GUIA MATLAB

LINEALES

En este taller usaremos el programa MATLAB con el fin de resolver ecuaciones no lineales

y sistemas de ecuaciones lineales, de manera rápida y fácil. Se usarán tanto las herramientas

propias de MATLAB, como rutinas creadas por el usuario que nos llevarán paso a paso a la

solución de problemas.

PRIMERA PARTE: SOLUCIÓN DE ECUACIONES NO LINE-

ALES DE LA FORMA F (X) = 0

Comencemos revizando unos fáciles ejemplos que nos mostrarán cuales son los pasos para re-

solver una ecuación no lineal empleando MATLAB, por medio de rutinas creadas por el usuario

y que con anterioridad deben ser creadas por usted en la forma como se explicó en la inducción

(y en la guía que aparece en la página del curso).

• Ejemplo: Determine valores aproximados de las soluciones positivas de la ecuación

1
2

exp³ x

3´ − sin (x) = 0

Solución: Lo primero que debemos hacer es graficar la función f (x) = 1

2 exp¡ x

3¢−sin (x)

para identificar las raices de la ecuación, para ello empleamos la instruccion fplot que

permite graficar funciones simbólicas. La sintaxis es:

fplot( ‘funcion’,[xmin xmax]) = grafica la función para los valores de x en el intervalo

[xmin,xmax].

fplot( ‘funcion’,[xmin xmax ymin ymax]) = grafica la función para los valores de x

en el intervalo [xmin,xmax] y las imágenes comprendidas en el intervalo [ymin,ymax].

1

Ambas instrucciones permiten agregar color adicionando ’color’ antes de cerrar el parén-

tesis.

Para distinguir las raíces adicionamos la instrucción gird on que activa una

cuadrícula a la gráfica de la función.

Ejecutamos las siguientes instrucciones

fplot(‘0.5*exp(x/3)-sin(x)’,[-10 10 -1 1]), grid on)

fplot(‘0.5*exp(x/3)-sin(x)’,[0 2 -0.3 0.3]), grid on)

con la primera identificamos el intervalo donde se encuentran las raices positivas y con el

segundo obtenemos la figura:

Para aproximar las raices podemos emplear el método de Bisección en los intervalos [0.5, 1]

y [1.5, 2], clamamente la función es continua por ser la suma de funciones continuas

(exponencial y trigonometrica). También podemos emplear el método de Newton (puesto

que la función, la primera y la segunda derivada son continuas en todos los reales), pero

debemos tomar el valor inicial muy cerca de la raíz para no caer en una zona de divergencia

de la sucesión (es decir, que una aproximación sea el valor donde la función alcanza el

mínimo).

Para empezar, apliquemos el método de bisección con una tolerancia de 10−6 (es decir,

2

queremos que cuando la distancia entre dos aproximaciones sea menor a 10−6 paremos),
para ello utilizamos la función "bisección" en MATLAB, que debe ser creada por el usuario

de la misma manera como se indicó en la inducción (Una muestra sencilla de un programa

típico se incluye al final del la guía). Digitando

biseccion(’0.5*exp(x/3)-sin(x)’,0.5,1,0.000001)

en la ventana Commad Window y oprimiento la tecla ENTER, obtenemos los siguientes

resultados

3

Y así mismo para la otra aproximación

Para emplear el método de Newton tomaremos la misma tolerancia para la distancia entre

las aproximaciones y emplearemos una función del método de Newton en MATLAB, que

debe ser creada por el usuario (una muestra está al final la guía). Si digitamos en la

ventana Command Window:

newton(’0.5*exp(x/3)-sin(x)’,0,0.000001)

4

obtenemos los siguientes resultados para la primer raíz

y para la segunda raíz

• Ejemplo: El factor de fricción f para los fluidos turbulentos en una tubería está dado

por

1
√f

= 1.14 − 2 log10µ e

D

+

9.35

Re√f¶

llamada correlación de Colebrook, donde Re es el número de Reynolds, e es la aspereza

de la superficie de la tubería y D es el diámetro de la tubería. Resolver la ecuación para

f utilizando el método de punto fijo para los siguientes casos:

1. (a) D = 0.1m, e = 0.0025m, Re = 3 × 104
(b) D = 0.1m, e = 0.0005m, Re = 5 × 106

Solución: Si queremos resolver el problema empleando el método de punto fijo debemos

5

llevar la ecuación a la forma

para ello llamemos x =

1
√f

, y así

x = g (x) ,

luego la función:

+

9.35
Re

D

x = 1.14 − 2 log10µ e
g (x) = 1.14 − 2 log10µ e

D

x¶
x¶ .

+

9.35
Re

puede ser una posible función de iteración de punto fijo para f.

Como queremos encontrar un punto fijo graficamos la función g y la recta y = x para

tomar un valor inicial. Empecemos con el caso (a) donde

g (x) = 1.14 − 2 log10µ 0.0025
3 × 104 x¶
g (x) = 1.14 − 2 log10¡0.025 + 3.1167 × 10−4x¢

9.35

0.1

+

Para graficar utilizamos las instrucciones

fplot(’1.14-2*log10*(0.025+3.1167*10^(-4))’,[-2 2 -2 2]),grid on

hold on,fplot(’1*x’,[-2 2 -2 2])

donde hold on permite graficar varias funciones en un mismo sistema coordenado.

6

Obtenemos:

De la gráfica anterior podemos ver que la función g cumple las condiciones del Teorema

Fundamental de Punto Fijo (¿dónde y por qué?). Para emplear el método tomaremos la

misma tolerancia para la distancia entre las aproximaciones y emplearemos una función

de punto fijo para MATLAB creada por el usuario (un ejemplo de ésta se muestra al final

de la guía). Si digitamos:

puntofijo(’1.14-2*log10*(0.025+3.1167*10^(-4))’,1,0.000001)

obtenemos

7

luego el valor para f lo obtenemos de

1
√f

= 1.0850701544

f = (1.0850701544)−2

f = 0.84934544855

En forma similar se resuelve para el caso (b).

De manera más fácil podemos utilizar las rutinas internas de MATLAB, que encuentran

raíces de ecuaciones no lineales. Algunas de estas son:

• fzero(fun,x0): Encuentra una raíz de la función f (x) = f un, que debe ser definida
antes o allí mismo con la instrucción inline. Esta instrucción busca un cero de fun cerca

del punto x0 especificado por el usuario, y se basa en los métodos de la secante, bisección

e interpolación cuadrática inversa, por lo que la función ingresada debe ser continua.

Por ejemplo, para el ejemplo 1 podemos digitar:

>> x=fzero(inline(’0.5*exp(x/3)-sin(x)’),0)

con la cual obtenemos: x = 6.7721e-001. Que es el mismo resultado que se obtuvo con

las rutinas de bisección y Newton que se usaron en el ejemplo 1.

• roots(p): Encuentra todas las raíces de un polinomio p, tanto reales como complejas.
Para usarla es necesario tener en cuenta que en MATLAB, un polinomio se representa

por medio de un vector de coeficientes.

Por ejemplo, para calcular todas las raíces del polinomio P (x) = x5 + 3x3 − 2x + 1,
debemos digitar la instrucción:

>> r=roots([1 0 3 0 -2 1])

(donde el vector [1 0 3 0 -2 1] representa al

polinomio P).

con la cual obtenemos:

r =

8

-3.3865e-002 +1.8892e+000i

-3.3865e-002 -1.8892e+000i

-9.0261e-001

4.8517e-001 +2.7374e-001i

4.8517e-001 -2.7374e-001i

que como vemos coincide con lo que esperábamos: Las raíces complejas siempre vienen

por pares conjugados y hay al menos una raíz real.

SEGUNDA PARTE: SOLUCION DE SISTEMAS DE ECUA-

CIONES NO LINEALES

Queremos resolver un sistema de la forma F (X) = 0, donde F es un campo vectorial. Para ello

emplearemos el método de Newton para Sistemas no lineales; para empezar consideraremos el

caso de un sistema de dos ecuaciones no lineales.

Ejemplo: Dado el sistema no lineal

x =

y =

8x − 4x2 + y2 + 1
2x − x2 + 4y − y2 + 3

8

4

encuentre una aproximación a la solución del sistema.

Solución: Para resolver el sistema primero lo llevamos a la forma F (X) = 0, para ello

expresamos cada ecuación igualada a cero

4x2 − y2 − 1 = 0
−2x + x2 + y2 − 3 = 0

9

definimos

así

f1 (x, y) = 4x2 − y2 − 1
f2 (x, y) = −2x + x2 + y2 − 3

F

: R2 → R2

X = (x, y) → F (X) =

f1 (x, y)

f2 (x, y)  .

queremos resolver F (X) = 0 por el método de Newton para sistemas donde la ecuación de

iteración está dada por:

X (n+1) = X (n) + H (n)

con

H(n) = −JF³X (n)´ F³X (n)´

Ahora realizamos una gráfica para identificar un valor inicial, empleando la instrucción de

Matlab para graficar funciones implicitas contour. Antes de esta instrucción es necesario

crear un arreglo bidimensional para la gráfica con la instrucción meshgrid. Si ejecutamos las

instrucciones

»xa=-3:0.1:3; ya=-3:0.1:3; [x,y]=meshgrid(xa,ya);

»f1=4*x.^2-y.^2-1; f2=-2*x+x.^2+y.^2-3;

»contour(x,y,f1,[0,0],’k’); hold on; grid on; contour(x,y,f1,[0,0],’k’);

10

obtenemos

en la primera línea indicamos la región inicial sobre la cual queremos ver las graficas de las

funciones y generamos el arreglo bidimensional sobre el cual aparecerá la gráfica, en la segunda

línea definimos las funciones que queremos graficar (las operaciones deben ir acompañadas por

punto ya que son operaciones entre vectores) y en la tercera línea aparece la instrucción de

contour para graficar, hold on para que ambas esten en el mismo plano y grid on para poder

tomar los valores iniciales.

De la gráfica vemos que existen 4 soluciones, busquemos una de ellas, tomando por valor
inicial x0 = −0.8 y y0 = −1, es decir, X (0) = [−0.8,−1]T . Aplicamos el método de Newton
no Lineal con una tolerancia de 10−6 (es decir, queremos que cuando la distancia entre dos
aproximaciones sea menor a 10−6 paremos), para ello utilizamos la función de newtonnl en
matlab (esta al final del la guía) digitando

newtonnl(4*x^2-y^2-1,-2*x+x^2+y^2-3,[-0.8,-1]’,0.000000001,10)

11

obtenemos los siguientes resultados

salió porque cumplio la condición de la tolerancia.

TERCERA PARTE: SOLUCION DE SISTEMAS DE ECUA-

CIONES LINEALES

En esta parte veremos como usar MATLAB para resolver sistemas de ecuaciones lineales. Al

igual que en la sección anterior usaremos tanto las rutinas internas de MATLAB como rutinas

propias creadas por el usuario, como las que se encuentran disponibles en la guía básica de

MATLAB que se encuentra en la página del curso. Pero además usaremos una combionación

de los dos procedimientos.

Comencemos con las fáciles y sencillas rutinas internas de MATLAB: Si queremos resolver

un sistema de ecuaciones lineales, lo primero que debemos hacer es escribirlo en la forma

matricial Ax = b. Por ejemplo consideremos el sistema de ecuaciones:

2x − y − z = 2
2x + 3y − z = 6
x + 2y + 3z