M-archivos: Herramientas didácticas para procesamiento de información y visualización, Apuntes de Matemáticas

¡Descarga M-archivos: Herramientas didácticas para procesamiento de información y visualización y más Apuntes en PDF de Matemáticas solo en Docsity!

Invitación al Análisis Numérico

Carlos Enrique Mejía Salazar

Universidad Nacional de Colombia, Medellín

Escuela de Matemáticas

Julio 2002

Prefacio

En estas notas ofrecemos una variedad de temas del análisis numérico, que pueden servir de motivación y de guía a estudiantes de nuestra Universidad y a otros interesados. Se pueden considerar como unas notas de clase, pero, para la selección de temas y enfoques, no se siguió un programa de curso dado. Los más de 100 ejercicios están repartidos en pequeñas listas en medio de la exposición, ejercicios suplementarios y exámenes de en- trenamiento. Hay cerca de 30 ejemplos resueltos en detalle, que casi siempre aparecen con su rutina MATLAB acompañante. Dichas rutinas se presentan únicamente con …nes didácticos y no sustituyen software de calidad. Pretendemos publicar estas notas en la WEB y mantenerlas en evolu- ción. Por eso los aportes de usuarios son muy importantes. Comenta- rios, correcciones y sugerencias son bienvenidos en la dirección electró- nica cemejia@perseus:unalmed:edu:co: Expresamos nuestro reconocimiento a la Universidad Nacional de Colombia, por habernos concedido el año sabático durante el cual es- cribimos estas notas. Parte de este documento se utilizará como guía en un cursillo de la XIII Escuela Latinoamericana de Matemáticas, que se reunirá en Cartagena, Colombia de julio 29 a agosto 3 de 2002.

Medellín, julio 15 de 2002

8 Contenido

10 1. Presentación y repaso

numérico es un pre-requisito para la computación cientí…ca. En estas notas, hacemos de cuenta que el problema matemático a resolver ya ha sido de…nido. No nos ocupamos del problema físico o de ingeniería directamente. Aquí nos concentramos en algoritmos y uso de software. Además, es en medio de ejemplos concretos que presentamos las no- ciones de estabilidad y error, indispensables para evaluar la calidad de un método numérico. Para compensar por las obligatorias omisiones, ofrecemos bibliografía que permita profundizar en los temas a los in- teresados.

Frecuentemente, los métodos numéricos, muy en especial los que sir- ven para resolver ecuaciones diferenciales, conducen a problemas de álgebra lineal en los que las matrices tienen alguna estructura y la mayoría de sus elementos son nulos. Para estos problemas, el álgebra lineal numérica ofrece métodos especiales, en los que se está traba- jando intensamente desde hace varios años. Basta ver, por ejemplo, las Memorias de Copper Mountain Conference del año 2000 [35], que con- tienen 6 artículos sobre precondicionamiento, 6 sobre cálculo de valores propios y 3 sobre métodos multigrid. Esta es una conferencia organi- zada anualmente por University of Colorado en cooperación con SIAM (Society for Industrial and Applied Mathematics), que trata cada dos años sobre métodos multigrid y en el año intermedio sobre métodos iterativos en general. Los temas del álgebra lineal numérica son de gran importancia den- tro del análisis numérico. Dentro de estas notas veremos problemas que al discretizarlos se convierten en problemas de álgebra lineal que aquí resolvemos utilizando a MATLAB. Pero este tema es tan importante, que amerita consideración especial. Esperamos en un futuro cercano complementar estas notas con unas acerca de álgebra lineal numérica. Lo único que haremos, por ahora, es que al …nal de estas notas, comen- tamos brevemente sobre software especializado que invitamos a conocer a los interesados.

1.2 Bosquejo del contenido

La organización de estas notas es la siguiente: Enseguida, en este mismo capítulo, presentamos una sección de repaso que hasta ahora solo con-

1.3. Repaso: Normas en espacios vectoriales 11

tiene material sobre normas en espacios vectoriales. Vislumbramos que la sección de repaso crecerá un poco posteriormente. En el próximo capítulo introducimos el software MATLAB, que es el que escogimos para ejemplos y ejercicios en estas notas. En el capítulo siguiente se presenta la solución de ecuaciones no lineales por varios métodos, in- cluyendo el método de Newton. En los dos capítulos que siguen nos ocupamos de interpolación y de integración numérica y en los dos si- guientes nos concentramos en la solución numérica de ecuaciones dife- renciales. Este tema lo iniciamos con métodos de un paso para proble- mas de valor inicial para ecuaciones diferenciales ordinarias. Enseguida continuamos con métodos numéricos para problemas con valores en la frontera para ecuaciones diferenciales ordinarias. Más tarde, hace- mos una revisión de métodos numéricos para la solución de ecuaciones diferenciales parciales de tipo parabólico. En la parte …nal, hacemos re- ferencia a software de álgebra lineal numérica, con énfasis en métodos iterativos no estacionarios para la solución de sistemas de ecuaciones lineales. Al …nal de los capítulos que tratan sobre métodos numéricos, in- cluimos una lista de ejercicios suplementarios y un examen de entre- namiento. La mayoría de los problemas en estas secciones …nales están basados en problemas y ejemplos de las referencias Cheney y Kincaid (1980) [8], Golub y Ortega (1993) [12] y Johnson y Riess (1982) [19].

1.3 Repaso: Normas en espacios vectoriales

Existen nociones topológicas y geométricas que son indispensables para estudiar análisis numérico. Posiblemente la más importante sea la idea de norma, que logra generalizar a espacios abstractos los conceptos escalares de valor absoluto y módulo.

De…nición 1 Sea V un espacio vectorial no vacío, real o complejo. Una norma en V es una función N : V! R que cumple las siguientes condiciones: a. N (x) ¸ 0 para todo x 2 V y N(x) = 0 si y solo si x = 0: b. N (®x) = j®j N (x) para todo x 2 V y todo escalar ®: c. N (x + y) · N (x) + N (y) para todo x; y 2 V:

1.3. Repaso: Normas en espacios vectoriales 13

Estas …guras fueron generadas con el siguiente M-archivo. El nombre M-archivo es debido a que los programas en lenguaje MATLAB se denominan en inglés M-…les. Los nombres de archivos de programas en este lenguaje siempre tienen la letra m como extensión y son archivos de texto plano que se pueden abrir en procesadores de texto sin formato como Bloc de Notas en WINDOWS. MATLAB proporciona su propio editor de archivos en las versiones más recientes.

% normas.m Herramienta didactica, Carlos E. Mejia, 2002 % bolas unitarias con normas 1 y 2 t = 0:pi/20:2*pi; subplot(121) plot(sin(t),cos(t),’r’) axis([-1 1 -1 1]);axis square; title(’Bola Unitaria en Norma 2’); x=-1:.1:1;y=1-abs(x); subplot(122) plot(x,y,’r’,x,-y,’r’) axis([-1 1 -1 1]);axis square title(’Bola Unitaria en Norma 1’); print -deps2 .nfignnormas.eps

1.3.2 Normas de matrices

Para el espacio vectorial Rn£n^ de matrices cuadradas n £ n con ele- mentos reales, también se requiere el concepto de norma. La notación que se usa es la misma, pero, además de las propiedades a., b. y c. de la de…nición 1, se pide que se cumplan las siguientes dos condiciones:

d. kABk · kAk kBk ; para toda A; B 2 Rn£n e. Toda norma matricial debe ser compatible con alguna norma vec- torial, es decir, kAxk · kAk kxk ; para toda A 2 Rn£n; x 2 Rn^ : Las condiciones para una norma en el espacio matricial Cn£n^ son análogas a las de arriba y no las precisaremos aquí. Un concepto necesario para la de…nición de una de las normas ma- triciales más utilizadas es el de radio espectral.

De…nición 3 Sea B 2 Rn£n^ con L = f¸ 1 ; ¸ 2 ; :::; ¸r g su espectro, es decir, el conjunto de sus valores propios distintos. El radio espectral de

14 1. Presentación y repaso

B es ½(B) = max ¸j 2 L

j¸j j :

Ejemplo 4 La matriz

A =

tiene radio espectral ½ (A) = 6; pues sus valores propios son ¡ 3 y ¡ 6 :

Las normas matriciales más utilizadas son las siguientes: Sea A = [aij ] 2 Rn£n^ (de…niciones análogas para A 2 Cn£n:)

1. Norma in…nito: kAk 1 = max i=1; 2 ;:::;n

Pn j=

jaij j : (Compatible con norma

in…nito vectorial).

2. Norma 1 : kAk 1 = max j =1; 2 ;:::;n

Pn i=

jaij j : (Compatible con norma 1 vec-

torial).

3. Norma 2 : kAk 2 =

AT^ A

. (Compatible con norma 2 vecto- rial).

4. Norma de Frobenius: kAkF =

Pn i=

P^ n j =

jaij j^2

: (Compatible con

norma 2 vectorial).

Ejemplo 5 La matriz

A =

tiene kAk 1 = 18; kAk 1 = 18; kAk 2 = ½ (A) = 13: 4833 y kAkF =

p 184 :

En este ejemplo la norma 2 de A es igual a su radio espectral por ser A una matriz simétrica. Esto no es cierto en general, como puede verse con matrices sencillas como

A =

para la cual kAk 2 = 1 pero ½ (A) = 0:

16 2. Primeros pasos en MATLAB

2.1 La línea de comandos

La línea de la pantalla donde MATLAB espera comandos, empieza por el símbolo >>. A la línea se le llama línea de comandos. Para …nalizar trabajo con MATLAB, se puede escoger la opción Exit MATLAB del menú File o entrar la palabra quit en la línea de comandos. Si lo que desea es abortar los cálculos que está haciendo, presione CTRL-C. La tecla de la ‡echa " sirve para recordar las expresiones que se han escrito antes en la línea de comandos. Hay cuatro sistemas de ayuda en línea: help, lookfor, helpwin y helpdesk. El último es el sistema más completo de todos, está en formato html, es apto para cualquier explorador de internet y posee enlaces para gran cantidad de archivos con formato PDF. El primero, help, lo explicamos por medio de un ejemplo: en la línea de comandos escriba help linspace. La respuesta del sistema es mostrar las primeras líneas precedidas por % del M-archivo linspace.m que listamos a con- tinuación:

function y = linspace(d1, d2, n) %LINSPACE Linearly spaced vector. % LINSPACE(x1, x2) generates a row vector of 100 linearly % equally spaced points between x1 and x2. % % LINSPACE(x1, x2, N) generates N points between x1 and x2. % % See also LOGSPACE, :. % Copyright (c) 1984-98 by The MathWorks, Inc. % $Revision: 5.6 $ $Date: 1997/11/21 23:29:09 $ if nargin == 2 n = 100; end if n~= y = d1:(d2-d1)/(n-1):d2; else y = d2; end Si el M-archivo por el que se pide ayuda (help) fue creado por un usuario y tiene líneas de comentarios en la primera parte, también las

2.2. Los M-archivos 17

mostrará. Así que la recomendación es incluir líneas de comentario en la parte de arriba de todo M-archivo que se cree. Es una costumbre que le ahorra tiempo y esfuerzo hasta a uno mismo. Para conocer las variables que hay en un momento dado en el espa- cio de trabajo, hay dos comandos: who y whos. Ensáyelos y note la diferencia. Por cierto, también es bueno notar que en MATLAB todas las variables representan matrices. Los escalares son matrices 1 £ 1 : De hecho la palabra MATLAB es una sigla que signi…ca Matrix Laboratory. Finalmente, antes de considerar M-archivos, advertimos que el núcleo básico de MATLAB es para cálculo numérico y no simbólico. Por tanto, si a=5 y b=7+a, entonces b=12 y sigue siendo 12 aunque el valor de a cambie. MATLAB administra el cálculo simbólico a través de una de sus cajas de herramientas, que debe comprarse por aparte.

2.2 Los M-archivos

La línea de comandos deja de ser práctica cuando se desea hacer un pro- ceso con varias órdenes o comandos. Se requiere entonces crear archivos con los comandos MATLAB que se requieren para cada tarea. Se les llama M-archivos (M-…les en inglés.) Son de texto plano (ASCII) y MATLAB proporciona un editor con una agradable codi…cación de colores que ayuda mucho en el proceso de creación y depuración de dichos archivos. Los M-archivos que requieren argumentos de entrada y entregan valores de salida, se llaman funciones. Los que no tienen ni argumentos de entrada ni valores de salida se llaman guiones ( script en inglés.) El M-archivo normas.m del ejemplo 2, que consideramos en la sección 1.3, es un guión. Los M-archivos están en directorios. Para que MATLAB los encuen- tre, se debe cumplir al menos una de las siguientes condiciones: i. El directorio en el que están los archivos pertenece al search path de MATLAB. ii. El directorio en el que están los archivos es el directorio actual de MATLAB. Para saber cúal es su directorio actual, use el comando pwd. Para hacer que el directorio c:nmet sea su directorio actual, escriba la orden cd c:nmet. El search path se puede ver y actualizar por medio de la opción Set

2.2. Los M-archivos 19

function x=fib(n); % Herramienta didactica, Carlos E. Mejia 2002 % sucesion de Fibonacci, n< 1477 % uso: x=fib(n) x=ones(1,n); for j=3:n x(j)=x(j-1)+x(j-2); end Al ejecutar x=…b(10) se obtiene x = 1 1 2 3 5 8 13 21 34 55 Miremos ahora el generador de árboles: % tree1.m, Herramienta didactica, Carlos E. Mejia, 2002 % generacion aleatoria de numeros para construir un arbol clf hold on x=[0;0];y=[0;rand];line(x,y) % a lo más 10 niveles de ramas n=rand* for i=1:n, [m,longi]=size(x); angu1=(pi/2)rand(1,longi); angu2=(pi/2)rand(1,longi)+angu1; x=[x(2,:),x(2,:);x(2,:)+(2/i)randcos(angu1),... x(2,:)+(3/i)randcos(angu2)]; y=[y(2,:),y(2,:);y(2,:)+(2/i)randsin(angu1),... y(2,:)+(3/i)randsin(angu2)]; line(x,y); end title(’Arbol con 10 niveles de ramas’) print -deps2 .nfigntree1.eps Uno de los árboles que pueden obtenerse con esta rutina, aparece en la siguiente …gura.

20 2. Primeros pasos en MATLAB

-1.5^0 -1 -0.5 0 0.5 1 1.

1

2

3

4

5 Arbol con 10 niveles de ramas

Ahora es un buen momento para reposar y asimilar los comandos MATLAB que se han utilizado. También es buen momento para estu- diar algunos otros por medio de la ayuda en línea y experimentación directa. Los siguientes ejercicios pueden servir de guía para estas tareas.

2.3 Ejercicios

1. Use el sistema de ayudas de MATLAB para encontrar información sobre las siguientes constantes especiales: realmin, realmax, inf, NaN, pi, eps.
   2. Repase los M-archivos vistos hasta ahora y obtenga más informa- ción sobre los comandos utilizados en ellos, por medio de los sistemas de ayuda de MATLAB. En particular, consulte sobre for, if, end, tres puntos (...), hold y line.
   3. Intentemos calcular de forma ingenua lim x! 0 r (x) ;donde

r (x) = x¡^3

log (1 ¡ x) + x exp

³ (^) x 2

El resultado exacto es ¡

; se puede obtener por Regla de L’Hopital.