• Author / Uploaded
• Jerson Manuel VC

Citation preview

UNIVERSIDAD NACIONAL DE INGENIER´IA FACULTAD DE INGENIER´IA CIVIL

´ GU´IA DE LABORATORIO DE METODOS ´ NUMERICOS CON MATLAB Presentado por:

LEONARDO FLORES GONZALES / CRISTINA NAVARRO FLORES ´ LIMA-PERU 2018

Cap´ıtulo 1 Laboratorio N◦1: Introducci´ on El objetivo del laboratorio, es de comprender y aplicar los diversos m´etodos num´ericos vistos en la teor´ıa del curso, para la resoluci´on de problemas, haciendo uso del programa Matlab. MATLAB puede considerarse como un lenguaje de programaci´on de c´alculo num´erico orientado a matrices. Por tanto, ser´a m´as eficiente si se dise˜ nan los algoritmos en t´erminos de matrices y vectores. Presenta las siguientes caracter´ısticas notables: La programaci´on es mucho m´as sencilla. Cuenta con una biblioteca matem´atica amplia. Tiene abundante herramientas gr´aficas, inclu´ıda funciones de interfaz con el usuario. Capacidad de vincularse con otros lenguajes de programaci´on o softwares. Una caractr´ıstica extraordinaria en Matlab cualquier variable puede contener n´ umeros de cualquier tipo sin una declaraci´on especial durante la programaci´oncon lo cual ´esta u ´ltima se hace m´as r´apida y productiva.

1.1.

Ventanas en Matlab

Matlab presenta cinco ventanas importantes que podr´ıa variar seg´ un la versi´on. Ventana de comandos, es la ventana principal que se muestra al iniciarse Matlab, que tambi´en puede llamarse ventana de ejecuci´on, y es donde se ejecuta una orden en cuanto le demos un h Enteri. Ventana del Editor, en esta ventana se escribe el c´odigo del programa a ejecutarse, que puede ser s´olo programa o funci´on creada por el usuario. 2

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil Ventana del Historial, es aqui donde se muestra las ´ordenes realizadas en la ventana de comando, y los guarda por fechas. Ventana de directorio, en esta ventana se muestran los archivos o programas alamacenados en la direcci´on de trabajo. Ventana del espacio de trabajo (workspace), contiene las variables que se crea y la almacena en la memoria s´olo durante la sesi´on de Matlab.

1.2.

Algunas funciones propias de Matlab

Se muestran las siguientes funciones predeterminadas por Matlab:

sin

seno

cos

coseno

tan

tangente

asin

arco seno

acos

arco coseno

atan

arco tangente

log

logaritmo natural

log10 logaritmo en base 10 exp

exponencial

sqrt

ra´ız cuadrada

funci´ on de Matlab

significado

abs

valor absoluto

min

devuelve el valor m´ınimo de un vector

max

devuelve el valor m´aximo de un vector

sum

devuelve la suma de los elementos de un vector

norm

norma euclidea de un vector

length

longitud de un vector

size

devuelve el orden de una matriz

det

devuelve el determinante de una matriz cuadrada

Laboratorio de M´etodos Num´ericos

3

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.3.

Vectores y Matrices en Matlab

Un vector fila en matlab es: v1×n = [v1 v2 . . . vn ], separados por un espacio o tambi´en puede usarse una coma (,). Un vector columna en matlab es: vn×1 = [v1 ; v2 ; . . . ; vn ], separados por un punto y coma (;). Luego una matriz es: Mf c = [m11 m12 . . . m1c ; m21 m22 . . . m2c ; mf 1 mf 2 . . . mf c ]

Ejemplo 1.3.1 Leer los siguientes vectores y matriz en la ventana de comandos de matlab   a) Vector fila v = 2 4 6 −4 8



2



   4       b) Vector columna w =   6     −4    8

Laboratorio de M´etodos Num´ericos

4

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil



1

2

3

4

5



   c) La matriz M =  6 7 8 9 10   11 12 13 14 15

Recordando que se usar´a el espacio en blaco o coma para separar los elementos por filas y el punto y coma para los saltos de l´ınea.

Laboratorio de M´etodos Num´ericos

5

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.3.1.

Algunas matrices predefinidas por Matlab

i) Matriz de unos.

(a)

(b)

(c)

Figura 1.1: Usos de la funci´on ones( )

ii) Matriz de ceros.

Laboratorio de M´etodos Num´ericos

6

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

(a)

(b)

(c)

Figura 1.2: Usos de la funci´on zeros( )

iii) Matriz identidad y variantes

(a)

(b)

Figura 1.3: Usos de la funci´on eye( )

Laboratorio de M´etodos Num´ericos

7

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.3.2.

Operaciones en una matriz o vector

Asumiendo que se tiene la siguiente matriz   1 2 3 4 5    M = 6 7 8 9 10   11 12 13 14 15 se puede hacer las siguientes operaciones.

Se puede adicionar un vector columna al inicio o al final de la matriz M . Se puede adicionar un vector fila al inicio o al final de la matriz M . Se puede extraer submatrices o vectores de la matriz M .

Es decir, se puede hacer las siguientes acciones en Matlab La expresi´on A = M (: , 4) se lee: en A copia de M todas las filas de la columna 4 La expresi´on A = M (2 , :) se lee: en A copia de M la fila 2 todas las columnas La expresi´on A = M (1 : 2, 3 : 5) se lee: en A copia de M desde la fila 1 hasta la fila 2, desde la columna 3 hasta la columna 5

Veamos:

Laboratorio de M´etodos Num´ericos

8

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Figura 1.4: Adicionando una columna C1 m´as a la matriz M

Laboratorio de M´etodos Num´ericos

9

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Figura 1.5: Adicionando dos filas m´as f1 y f2, al final de la matriz M

Figura 1.6: Extrae desde la fila 2 hasta la fila 3, desde la columna 2 hasta la columna 4 de la matriz A

Laboratorio de M´etodos Num´ericos

10

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Figura 1.7: Extrae una columna o fila de la matriz A

Laboratorio de M´etodos Num´ericos

11

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.4.

Expresiones matem´ aticas

S´ımbolo

significado

+/−

Suma o resta de n´ umeros o matrices

∗

Multiplicaci´on usual de matrices

/

Divisi´on

∧

potencia

\

Soluciona un sistema lineal de ecuaciones

·∗

Multiplica dos vectores o matrices de la misma dimensi´on elemento a elemento

·/

Divide dos vectores o matrices de la misma dimensi´on elemento a elemento

·∧

Para elevar a la potencia cada elemento de una vector o matriz

Ejemplo 1.4.1 Sea el n´ umero x = 2, se desea evaluar la siguiente expresi´on y = x sen(2x) + ln(x2 ) − e3x Luego en Matlab se tiene:

Figura 1.8: Para un n´ umero x = 2

Ejemplo 1.4.2 Para la misma expresi´on matem´atica y = x sen(2x) + ln(x2 ) − e3x Se quiere evaluar para una lista de n´ umeros almacenados en un vector x=



1 2 3 4 5 6 7 8 9 10



Luego en Matlab se tiene:

Laboratorio de M´etodos Num´ericos

12

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Figura 1.9: Para un vector x de 10 elementos

1.4.1.

Inversa de una matriz

La inversa de una matriz se ejecuta con la orden inv Ejemplo 1.4.3 Dada la matriz 

1 2   9 5  A=  3 1  3 4

3 4



 2 1    1 2   1 2

En matlab ser´ıa de forma directa

Laboratorio de M´etodos Num´ericos

13

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.4.2.

Matriz transpuesta de una matriz

La transpuesta de una matriz se ejecuta con el s´ımbolo ’ (ap´ostrofe) Ejemplo 1.4.4 Dada la matriz 

1 2 3 4



   9 5 2 1    A=   3 1 1 2    3 4 1 2 En Matlab ser´ıa de forma directa

Laboratorio de M´etodos Num´ericos

14

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.5.

Gr´ afica de funciones en dos dimensiones

Para graficar funciones en 2D se presentar´a dos formas, sinembargo Matlab tiene mas funciones para gr´aficos en 2D asi como para 3D.

1.5.1.

Funci´ on plot()

Esta funci´on requiere dos vectores X e Y de la misma cantidad de elementos, la primera alamcena las abscisas y la segunda las ordenadas de un punto a graficar. Es decir si la longitud de los vectores es n se tiene n puntos a graficar, plot grafica uniendo los puntos con un segmento, es por ello para que se observe una gr´afica continua y suave se requiere mayor cantidad de puntos muy cercanos. 

Sintaxis :





plot(X,Y)



Laboratorio de M´etodos Num´ericos

15

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil Ejemplo 1.5.1

Ejemplo 1.5.2

Laboratorio de M´etodos Num´ericos

16

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Adem´as que esta funci´on permite agregar detalles al gr´afico, como color de la l´ınea, forma de las l´ıneas, t´ıtulo al gr´afico, ejes, cuadros de texto interior, etc. Puede verse las diversas opciones escribiendo en la ventana de comando: help plot la opci´on hold on permite agregar un gr´afico sobre otro gr´afico. Ejemplo 1.5.3

Figura 1.10: Presenta m´as detalles en la gr´afica, adem´as se usa hold on

1.5.2.

Funci´ on fplot()

Esta funci´on es directa, y puede combinarse con la funci´on plot usando hold on 

Sintaxis :





fplot(cad,intervalo)



Donde cad: Es la expresi´on matem´atica en cadena a graficar. intervalo: Es el intervalo donde se grafica.

Laboratorio de M´etodos Num´ericos

17

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil Ejemplo 1.5.4

Observaci´ on 1.5.1 Matlab tiene una ayuda directa, si en la ventana de comando se digita help Nombre-Funcion devuelve en datalle el uso de la funci´on Nombre-Funcion a buscar.

Laboratorio de M´etodos Num´ericos

18

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.6.

Programaci´ on en Matlab

Para escribir el c´odigo de un programa vamos a la ventana del Editor, en el men´ u archivo(File), nuevo(new),archivo en blanco(Blank M-File)

Se tiene dos formas de programas: S´olo programa o´ s´olo funci´on, la diferencia est´a en que si es s´olo programa no requiere datos de ingreso ya es interno y si es funci´on requiere los datos de ingreso y es m´as flexible. Construiremos funciones. 

Sintaxis :





function [salida1, ... ,salida2]=Nombre-funcion(dato1, ...,dato2)



Y se guarda esta funci´on con el mismo nombre Nombre-funcion

Laboratorio de M´etodos Num´ericos

19

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil Ejemplo 1.6.1

Figura 1.11: Debe guardarse con el mismo nombre (pruebita). Se ejecuta en la ventana de comandos Ejecuci´on:

Figura 1.12: Para ejecutar,la ventana de comandos debe estar en el directorio donde est´a en archivo pruebita

Laboratorio de M´etodos Num´ericos

20

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

1.6.1.

Codificaci´ on

Para rutinas mas complejas se requiere usar condiciones y bucles. Matlab presenta las siguientes: if-else / si-sino Sintaxis: if hcondicioni Sentencia(1); .. . Sentencia(k); else Sentencia(1); .. . Sentencia(n); end while / mientras Sintaxis: while hcondicioni Sentencia(1); .. . Sentencia(k); end for / Desde Sintaxis:

Laboratorio de M´etodos Num´ericos

21

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

for I = ValInical : variacion : ValFinal Sentencia(1); .. . Sentencia(k); end

Ejemplo 1.6.2

Laboratorio de M´etodos Num´ericos

22

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Figura 1.13: Ejecutando la funci´on en la ventana de comandos

Funciones importantes para el mantenimiento de los programas y/o funciones

clc

Limpia la ventana de comandos

clear

Limpia las variables de la memoria

hold on

Agrega gr´aficas sobre otra ya existente

hold of

Limpia los gr´aficos

Laboratorio de M´etodos Num´ericos

23

Cap´ıtulo 2 Laboratorio N◦2 2.1.

Objetivos

Los objetivos que se pretende lograr en los dos primeros laboratorios son: Conocer y usar los comandos de Matlab. Programar en Matlab, algoritmos para una ecuaci´on no lineal: Bisecci´on, Newton y Secante

2.2.

M´ etodo de la Bisecci´ on

M´ etodo de la Bisecci´ on Siendo I = [a, b] el intervalo incial donde existe raiz. yi =

ai + bi , donde : Ii = [ai , bi ]; i = 1, 2, 3, . . . 2

Algoritmo Datos: a, b, f, tol

Salida: Raiz, niter

1. niter ← 0 2. A ← a, B ← b 3. Mientras (|A − B| > tol) 3.1 y ←

A+B 2

24

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

3.2 Si (f (A) ∗ f (y) < 0) B ← y sino A ← y 3.3 niter ← niter + 1 4. Raiz ← y

2.3.

M´ etodo de Newton Raphson

Siendo y0 una aproximaci´on inicial yi = yi−1 −

f (yi−1 ) , i = 1, 2, 3, . . . df (yi−1 )

Algoritmo Datos: y0 , f, df, tol 1. y1 ← y0 −

Salida: Raiz, niter

f (y0 ) df (y0 )

2. niter ← 1 3. Mientras (|y0 − y1 | > tol) 3.1 y0 ← y1 3.2 y1 ← y0 −

f (y0 ) df (y0 )

3.3 niter ← niter + 1 4. Raiz ← y1

2.4.

M´ etodo de Newton Raphson

M´ etodo de Newton Raphson Siendo y0 una aproximaci´on inicial yi = yi−1 −

f (yi−1 ) , i = 1, 2, 3, . . . df (yi−1 )

Laboratorio de M´etodos Num´ericos

25

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Algoritmo Datos: y0 , f, df, tol 1. y1 ← y0 −

Salida: Raiz, niter

f (y0 ) df (y0 )

2. niter ← 1 3. Mientras (|y0 − y1 | > tol) 3.1 y0 ← y1 3.2 y1 ← y0 −

f (y0 ) df (y0 )

3.3 niter ← niter + 1 4. Raiz ← y1

2.5.

M´ etodo de la Secante

M´ etodo de la Secante Siendo y0 , y1 dos aproximaciones iniciales yi = yi−1 −

f (yi−1 )(yi−1 − yi−2 ) , i = 2, 3, 4, . . . f (yi−1 ) − f (yi−2 )

Algoritmo Datos: y0 , y1 , f, tol

Salida: Raiz, niter

1. niter ← 0 2. Mientras (|y0 − y1 | > tol) 3.1 y2 ← y1 −

f (y1 )(y1 −y0 ) f (y1 )−f (y0 )

3.2 y0 ← y1 , y1 ← y2 3.3 niter ← niter + 1 3. Raiz ← y2

Laboratorio de M´etodos Num´ericos

26

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

2.6.

Aplicaci´ on

Problema 3

Por un canal trapezoidal fluye agua a una tasa de Q = 20 ms . La profundidad cr´ıtica y 2 ´ para dicho canal satisface la ecuaci´on 0 = 1 − Q 3 B donde g = 9,81 m2 , Ac : Area de la gAc

s

2

secci´on transversal (m ) y B : Ancho del canal en la superficie (m). Para este caso, el ancho y a´rea de la secci´on transversal se relacionan con la profundidad y, por medio de B = 3 + y, Ac = 3y + y 2 /2. Resuelva para la profundidad cr´ıtica con el uso de los m´etodos: 1. Bisecci´on con una tolerancia de 10−6 . 2. Newton con una tolerancia de 10−6 . 3. Secante con una tolerancia de 10−6 .

2.7.

Planteo del problema no lineal

Planteo del problema no lineal

f (y) = (9,81)[3y +

y2 3 ] − 400(y + 3) = 0 2

Si y > 0, usando bolzano tenemos y ∈ I = [1, 2]

Tambi´en se expresa su derivada como: df (y) = 3(9,81)[3y +

y2 2 ] (3 + y) − 400 2

Laboratorio de M´etodos Num´ericos

27

Cap´ıtulo 3 Laboratorio N◦3 (Evaluado) 3.0.1.

M´ etodo est´ atico, del an´ alisis s´ısmico de una edificaci´ on.

El an´alisis s´ısmico de una edificaci´on puede modelarse mediante cargas est´aticas equivalentes Pi aplicadas a cada nivel de la edificaci´on, como se muestra en la figura: Donde mi

es la masa equivalente del nivel y ki es la rigidez por entrepiso. El fuerza de corte en la base n P del edificio se calcula como V = ZC( mi )g, donde g = 9,81m/s2 es la gravedad, Z = 0,4 i=1

es el factor de zona y el factor C se calcula como: C = 1, si x > 2,5 C=

2,5 , x

si 1 ≤ x ≤ 2,5

C = 2,5, si x < 1 Donde x es el periodo de vibraci´on de la estructura. La carga equivalente Pi se calcula como Pi =

mi hi V n P mj hj

, donde hi es la altura del nivel i

j=1

medido desde la base. La fuerza de corte en la base del entrepiso i se calcula como Fi =

n P j=i

28

Pj .

UNIVERSIDAD NACIONAL DE INGENIER´IA Facultad de Ingenier´ıa Civil

Los desplazamientos laterales de los niveles se calcula como ui =

i P

Fj /kj y con estos

j=1

resultados, se puede calcular el periodo de vibraci’on x de la estructura, usando la ecuaci´on: v uP u n m u2 u i i u x = g(x), donde g(x) = 2π u i=1 (3.1) n tP Pi u i i=1

Dada la siguiente situaci´on, se tiene una edificacion de n = 5 pisos con los siguientes datos: Masas: m = [10 9 8 8 6] Alturas h = [3 5,5 8 10,5 13] Rigideces k = [10000 8000 8000 8000 8000] Se pide evaluar en x = 1,5 para los items (1, 2, 3, 4 y 5) lo siguiente: 1. Una funci´on que calcule la fuerza de corte en base V, en funci´on del periodo de vibraci´on x. 2. Una funci´on que calcule el vector de cargas equivalentes por niveles P = [Pi ], en funci´on del periodo de vibraci´on x. 3. Una funci´on que calcule el vector de cortes de entrepiso F = [Fi ], en funci’on del periodo de vibraci´on x. 4. Una funci´on que calcule el vector de los desplazamientos laterales u = [ui ], en funci´on del periodo de vibraci´on x. 5. Una funci´on que calcule el valor de g de la relaci´on (3.1) en funci´on del periodo de vibraci´on x, en base a los pasos anteriores. 6. Plantee la ecuaci´on no lineal f (x) = x − g(x) = 0, que tiene por soluci´on el periodo de vibraci´on fundamental de la estructura. Use el m´etodo de la secante para aproximar el periodo de vibraci´on (T), con un error de 10−5 usando valores iniciales de x0 = 0,5 y x1 = 1,5 7. Evalue el vector u en x = T y plotee la gr´afica u vs. h.

Laboratorio de M´etodos Num´ericos

29