Práctica 1. Carolina López Sánchez..doc

Laboratorio de Computación Científica 14

Curso 2013-

PRÁCTICA 1. Introducción a Matlab Alumno: Carolina López Sánchez MATLAB diferencia entre las mayúsculas y minúsculas. No es lo mismo escribir EJEMPLO que Ejemplo o ejemplo. Todas tus respuestas deben ser escritas en color rojo para diferenciarlas de los enunciados. Cuando termines esta práctica la debes enviar al email de tu profesor. Las secciones de la 1 a la 4 incluida se realizarán como máximo en 6 sesiones de laboratorio. SECCIÓN 1. Aprendiendo operadores en Matlab

comandos,

tipos

de

variables,

funciones

y

Las siguientes órdenes pueden serte útiles para el desarrollo de esta práctica:  help: Te ofrece ayuda sobre una función, incluyendo algún ejemplo sobre su uso y otras funciones relacionadas. Por ejemplo, teclea >>help clc 

lookfor: Lista todas las funciones en cuya ayuda aparece una palabra clave. Por ejemplo, si quisiéramos saber qué funciones tiene Matlab para manipular directorios podríamos teclear >>lookfor directory y a continuación teclear >>help pwd, >>help dir, etc. para saber cómo funciona cualquier función que nos interese.

En la siguiente tabla se indican algunas órdenes útiles en MATLAB. ORDEN SIGNIFICADO DE LA ORDEN HELP  Ayuda sobre órdenes y funciones internas de Matlab.  Ayuda sobre nuestras propias funciones .m WHAT Da una lista de funciones .m en el directorio especificado (Por ejemplo: >> what a:\ejemplos) CD  Cambia al directorio padre o directorio anterior (>> cd ..)  Cambia al directorio que se especifique (>> cd a:\ejemplos) DIR Da una lista del contenido del directorio en el que me encuentro PWD ¿Qué hace este comando? Varias órdenes: Dir p*: te muestra todos los archivos que hay y empiezan por p en el directorio; Dir pr*: aquellos que empiezan por pr; P tab: para ver las funciones internas además de las tuyas que comienzan por p; Dir: muestra los archivos del directorio; Edit: editor

1

1. Haciendo uso del comando HELP puedes saber las funciones de los comandos diary, who, clear, hold, format y edit. ¿Qué diferencia existe entre los comando edit y type?.

2

Las siguientes órdenes son interesantes para poder conocer los operadores definidos en Matlab. Algunos de ellos los vas a utilizar frecuentemente, sería útil que grabaras la información que tienes como respuesta. (a)help \ >> help \ Operadores y carácteres especiales Operadores aritméticos plus - Plus + uplus - Unary plus + minus - Minus uminus - Unary minus mtimes - Matrix multiply * times - Array multiply .* mpower - Matrix power ^ power - Array power .^ mldivide - barra invertida o dividir matriz de la izquierda \ mrdivide - (Raya vertical o dividir la matriz derecha ) Slash or right matrix divide ldivide - Left array divide .\ rdivide - Right array divide ./ kron - Kronecker tensor product kron Operadores relaciones eq - Igual ne - No igual lt - menor que gt - mayor que le - menor que o igual ge - mayor que o igual

/

== ~= < > <= >=

Los operadores lógicos relop - Short-circuit logical AND && relop - Short-circuit logical OR || and - Element-wise logical AND & or - Element-wise logical OR | not - Logical NOT ~ xor - Logical EXCLUSIVE OR any - True if any element of vector is nonzero all - True if all elements of vector are nonzero Special characters. colon - Colon : paren - Parentheses and subscripting () paren - Brackets [] paren - Braces and subscripting {} punct - Function handle creation @ punct - Decimal point . punct - Structure field access . punct - Parent directory .. punct - Continuation ... punct - Separator , punct - Semicolon ; punct - Comment % punct - Invoke operating system command ! punct - Assignment = punct - Quote ' transpose - Transpose .' ctranspose - Complex conjugate transpose ' horzcat - Horizontal concatenation [,] vertcat - Vertical concatenation [;] subsasgn - Subscripted assignment ( ),{ },. subsref - Subscripted reference ( ),{ },. subsindex - Subscript index Bitwise operators. bitand - Bit-wise AND. bitcmp - Complement bits. bitor - Bit-wise OR. bitmax - Maximum floating point integer. bitxor - Bit-wise XOR. bitset - Set bit. bitget - Get bit. bitshift - Bit-wise shift. Establecer los operadores. union - establecer la unión (Set union)

3

unique - conjunto único (Set unique) intersect - establecer intersección (Set intersection) setdiff - establecer diferencia (Set difference) setxor - establecer o-exclusiva (Set exclusive-or) ismember - cierto para miembro del conjunto (True for set member) See also arith, relop, slash, function_handle.

(b)help arith. >> help arith Los operadores aritméticos.(Arithmetic operators) + Plus. X + Y añade matrices X e Y. X e Y deben tener el mismo dimensiones que el que no un escalar (una matriz de 1 por 1). Un escalar se puede añadir a cualquier cosa. + Plus. X + Y adds matrices X and Y. X and Y must have the same dimensions unless one is a scalar (a 1-by-1 matrix). A scalar can be added to anything. - Menos. X - Y resta matriz X de Y. X e Y debe tener el mismo dimensiones que el que no un escalar. A escalar se puede restar de nada. - Minus. X - Y subtracts matrix X from Y. X and Y must have the same dimensions unless one is a scalar. A scalar can be subtracted from anything. * La multiplicación de matrices. X * Y es el producto de la matriz de X e Y. Un escalar (una matriz de 1 por 1) puede multiplicar cualquier cosa. De lo contrario, el número de columnas de X debe ser igual al número de filas de Y. * Matrix multiplication. X*Y is the matrix product of X and Y. Any scalar (a 1-by-1 matrix) may multiply anything. Otherwise, the number of columns of X must equal the number of rows of Y. . * Matriz de multiplicación X. * Y denota multiplicación elemento a elemento. X e Y debe tener las mismas dimensiones a menos que uno sea un escalar. A escalar se puede multiplicar por todo. * Array multiplication X.*Y denotes element-by-element multiplication. X and Y must have the same dimensions unless one is a scalar. A scalar can be multiplied into anything. ^ Poder Matrix. Z = X ^ y X es elevado a y, si y es un escalar y X es cuadrada. Si y es un número entero mayor que uno, la potencia se calcula mediante la repetida multiplicación. Para otros valores de y el cálculo implica valores y vectores propios. Z = x ^ Y es x a la potencia Y, si Y es una matriz cuadrada y x es un escalar, se calcula usando valores propios y vectores propios. Z = X ^ Y, donde X e Y son matrices, es un error. ^ Matrix power. Z = X^y is X to the y power if y is a scalar and X is square. If y is an integer greater than one, the power is computed by repeated multiplication. For other values of y the calculation involves eigenvalues and eigenvectors. Z = x^Y is x to the Y power, if Y is a square matrix and x is a scalar, computed using eigenvalues and eigenvectors. Z = X^Y, where both X and Y are matrices, is an error. . ^ Potencia matriz. Z = X ^ Y denota la operación elemento por elemento. X e Y deben tener las mismas dimensiones a menos que uno sea un escalar. Un escalar puede operar cualquier cosa. . ^ Array power. Z = X.^Y denotes element-by-element powers. X and Y must have the same dimensions unless one is a scalar. A scalar can operate into anything.

(c) help slash. >> help slash Matrix division. \ Barra invertida o división por la izquierda . A \ B es la división de la matriz de A entre B , que es aproximadamente el mismo que inv ( A) * B , excepto que se calcula de una manera diferente .

4

Si A es una matriz N por N y B es un vector columna con N componentes , o una matriz con varias de estas columnas , Entonces X = A \ B es la solución a la ecuación A * X = B calculada por Eliminación de Gauss. Un mensaje de advertencia aparecerá si A es gravemente reducido o casi singular. A \ OJO (SIZE ( A) ) produce el inverso de A. Si A es una matriz M - por - N con M < o > N y B es una columna vector con componentes M o una matriz con varias de estas columnas , entonces X = A \ B es la solución en los mínimos cuadrados sentido para la sub o sistema sobredeterminado de ecuaciones A * X = B. Los rango eficaz , K , de A se determina a partir de la descomposición QR con pivotante . Una solución X se calcula que tiene a lo sumo K componentes distintos de cero por columna . Si K < N esto generalmente no lo hará ser la misma solución que PINV ( A) * B. A \ OJO (SIZE ( A) ) produce una inversa generalizada de A. \ Backslash or left division. A\B is the matrix division of A into B, which is roughly the same as INV(A)*B , except it is computed in a different way. If A is an N-by-N matrix and B is a column vector with N components, or a matrix with several such columns, then X = A\B is the solution to the equation A*X = B computed by Gaussian elimination. A warning message is printed if A is badly scaled or nearly singular. A\EYE(SIZE(A)) produces the inverse of A. If A is an M-by-N matrix with M < or > N and B is a column vector with M components, or a matrix with several such columns, then X = A\B is the solution in the least squares sense to the under- or overdetermined system of equations A*X = B. The effective rank, K, of A is determined from the QR decomposition with pivoting. A solution X is computed which has at most K nonzero components per column. If K < N this will usually not be the same solution as PINV(A)*B. A\EYE(SIZE(A)) produces a generalized inverse of A. / La raya vertical o división derecha. B / A es la matriz de la división de A entre B , que es aproximadamente igual que B * INV ( A) , excepto que se calcula de una manera diferente . Más precisamente, B / A = (A ' \ B') ' . Ver \ / Slash or right division. B/A is the matrix division of A into B, which is roughly the same as B*INV(A) , except it is computed in a different way. More precisely, B/A = (A'\B')'. See \. . / Matriz división derecha. B. / A significa división elemento a elemento . A y B deben tener las mismas dimensiones a menos que uno sea un escalar . Un escalar se puede dividir con cualquier cosa . ./ Array right division. B./A denotes element-by-element division. A and B must have the same dimensions unless one is a scalar. A scalar can be divided with anything. . \ Matriz división izquierda. A. \ B. denota división elemento a elemento . A y B debe tener las mismas dimensiones a menos que uno seaun escalar . un escalar se puede dividir con cualquier cosa . .\ Array left division. A.\B. denotes element-by-element division. A and B must have the same dimensions unless one is a scalar. A scalar can be divided with anything.

(d)help relop. >> help relop Relational operators. < > Operadores relacionales. Los seis operadores relacionales son <, < = ,>, > =, == y ~ = . A Relational operators. The six relational operators are <, <=, >, >=, ==, and ~=. A < B does element by element comparisons between A and B and returns a matrix of the same size with elements set to logical 1 (TRUE) where the relation is true and elements set to logical 0 (FALSE) where it is not. A and B must have the same dimensions

5

(or one can be a scalar).

Ejemplo: >> A=[1,2;0,0];B=[1,0;9,0]; >> A
0 0

& Element-wise Logical AND. A & B es una matriz cuyos elementos son 1 lógico (verdad), donde tanto A como B tienen elementos distintos de cero , y 0 lógico (falso) , donde uno de ellos tiene un elemento cero . A y B deben tener las mismas dimensiones (o una ser un escalar ) . A & B is a matrix whose elements are logical 1 (TRUE) where both A and B have non-zero elements, and logical 0 (FALSE) where either has a zero element. A and B must have the same dimensions (or one can be a scalar).

Ejemplo: >> A=[1,2;0,0];B=[1,0;9,0]; >> A&B ans = 1 0

0 0

&& Short-Circuit Logical AND. A && B es un valor escalar que es el AND lógico del escalar A y B. Esta es una operación de " cortocircuito " en que MATLAB evalúa B solamente si el resultado no está totalmente determinado por A. Por ejemplo , si A es igual a 0, entonces toda la expresión se evalúa como 0 lógico (FALSO) , independientemente del valor de B. Bajo estas circunstancias , no hay necesidad de evaluar B porque el resultado ya es conocido A && B is a scalar value that is the logical AND of scalar A and B. This is a "short-circuit" operation in that MATLAB evaluates B only if the result is not fully determined by A. For example, if A equals 0, then the entire expression evaluates to logical 0 (FALSE), regardless of the value of B. Under these circumstances, there is no need to evaluate B because the result is already known. | Element -wise O lógico . A | B es una matriz cuyos elementos son 1 lógico (Verdad) cuando sea A o B tiene un elemento distinto de cero , y 0 lógico ( falso ) , donde ambos tienen cero elementos . A y B deben tener las mismas dimensiones (o una ser un escalar ) . | Element-wise Logical OR. A | B is a matrix whose elements are logical 1 (TRUE) where either A or B has a non-zero element, and logical 0 (FALSE) where both have zero elements. A and B must have the same dimensions (or one can be a scalar).

Ejemplo: >> A=[1,2;0,0];B=[1,0;9,0]; >> A|B ans = 1 1

1 0

| | Corto - Circuito lógico OR. A | | B es un valor escalar que es el OR lógico de escalar A y B. Esta es una operación de " cortocircuito " en que MATLAB evalúa B solamente si el resultado no está totalmente determinado por A. Por ejemplo , si A es igual 1, entonces toda la expresión se evalúa como 1 lógico (verdad) , independientemente del valor de B. Bajo estas circunstancias , no hay necesidad para evaluar B porque el resultado ya es conocido . || Short-Circuit Logical OR. A || B is a scalar value that is the logical OR of scalar A and B. This is a "short-circuit" operation in that MATLAB evaluates B only if the result is not fully determined by A. For example, if A equals 1, then the entire expression evaluates to logical 1 (TRUE), regardless of the value of B. Under these circumstances, there is no need

6

to evaluate B because the result is already known. ~ Complemento lógico ( NOT) . ~ A es una matriz cuyos elementos son 1 lógico (verdad) en la que A tiene cero elementos , y lógica 0 (FALSO ) en la que A tiene elementos distintos de cero . ~ Logical complement (NOT). ~A is a matrix whose elements are logical 1 (TRUE) where A has zero elements, and logical 0 (FALSE) where A has non-zero elements.

Ejemplo: >> A=[1,2;0,0]; >> ~A ans = 0 1

0 1

xor OR exclusivo . XOR ( A, B ) es 1 lógico ( VERDADERO ) donde A o B , pero no ambos , es distinto de cero . Ver XOR . xor Exclusive OR. xor(A,B) is logical 1 (TRUE) where either A or B, but not both, is non-zero. See XOR.

Ejemplo: >> A=[1,2;0,0];B=[1,0;9,0]; >> xor(A,B) ans = 0 1

1 0

2. Haz lo siguiente usando las líneas de comando de Matlab (indicadas por >>) Operación Ejemplo Indica los resultados de las operaciones Suma + >>format short Resultados de a, b y suma: >>a=7 >> format short >>b=2; %variable escalar >> a=7;b=2;suma=a+b; (puedes combertirlo en un >> suma comentario, y lo que hay antes, suma = a la derecha, no se ejecuta) >>suma=a+b; %El ; hace que 9 no se vea el resultado pero si se ejecuta la orden o sentencia. >>suma Resta >>resta=a-b; Resultado resta: >>resta >> resta=a-b; >> resta resta = Multiplicaci ón *

>>multiplicar=a*b

5 Resultado multiplicar: >> multiplicar=a*b multiplicar =

División /

>>div1=a/b

14 Resultado div1: >> div1=a/b div1 =

7

División \

3.5000 >>div2=a\b Resultado div2: (divide b entre a) >>div3=b\a %¿Es igual que >> div2=a\b div1? El resultado es el mismo porque \ significa que divides lo div2 = que hay en la derecha entre lo que hay en la izquierda, en este 0.2857 caso a entre b, al igual que en Resultado div3: div1 >> div3=b\a div3 =

Potencia ^ >>a=3;potencia=a^2 (primera en ejecución)

3.5000 % es lo mismo ya que b\a=a/b Resultado potencia: >> a=3;potencia=a^2 potencia =

>>a=9; >>a^3

9 Resultado: >> a=9; >> a^3 ans = 729 ¿Cómo se llama la variable a la que se le asigna el resultado? ans porque si no le asignamos algún nombre de una variable nosotros, Matlab le asigna ans al resultado de la operación.

>>ans >>c=7;d=8;ans=1 >>(c*d)/(ans+7) >>(c*d)/ans+7

¿Cuáles son los valores que va tomando la variable por defecto ans? 1. Ans=729 2. Ans=1 3. Ans=7 4. Ans=15 (a) Ejecuta en la línea de comandos el comando who. ¿Qué aparece en pantalla?. Las variables que hay activas en ese momento. (b) Ejecuta: >> clear b. Pregunta por el valor de b. ¿Qué sale en pantalla?. Clear b, elimina la asignación que le habiamos hecho a b, por lo tanto aparece que no está definida. (c) Ejecuta el comando who. ¿Cuál es la diferencia que encuentras con la ejecución anterior? Aparecen todas las variables activas en ese momento, la diferencia aparece cuando nos fijamos en b, en este caso no aparece ya que en el apartado b la hemos eliminado.

8

3. Busca el significado de las siguientes funciones internas de Matlab, pon un ejemplo y su resultado. Averigua cuál se puede escribir en minúsculas, mayúsculas o indistintamente. Normalmente en el segundo párrafo se especifica la sintaxis de la función. En general, en cualquier lenguaje de alto nivel las funciones trigonométricas deben tener el argumento o variable de entrada en radianes. Matlab también dispone de funciones trigonométricas cuyo argumento se defina en grados. Búscalas. Función Ejemplo ABS >> abs(-120) ABS (X) es el valor >> abs(sqrt(-1)) absoluto de los elementos de X. Cuando X es complejo, ABS (X) es el módulo complejo (magnitud) de los elementos de X.

Resultado ejemplo ans =

SQRT Sqrt (x) es la raíz cuadrada de los elementos de X. Complejo resultados se producen si X no es positivo.

>> sqrt(25)

ans =

>> sqrt(-20)

5 ans =

>> sqrt(-1)

0 + 4.4721i ans =

>> rand(2)

0 + 1.0000i ans =

RAND RAND distribuye uniformemente números pseudoaleatorios. R = RAND (N) devuelve una matriz N por N que contiene valores pseudoaleatorios elaborado a partir de una distribución uniforme en el intervalo unidad. RAND (M, N) o RAND ([m, n]) devuelve una matriz M-por-N. RAND sin argumentos devuelve un escalar. RAND (SIZE (A)) devuelve una matriz del mismo

120 ans = 1

0.8147 0.9058 >> rand(1,2)

0.1270 0.9134

ans = 0.6324

0.0975

>> rand ans = >> A=[1 2;3 4];rand(size(A))

0.2785 ans = 0.5469 0.9575

0.9649 0.1576

9

tamaño como A. SIN Sin (x) es el seno de los elementos de X.

>> sin(pi) >> sin(pi/2)

SIND >> SIND(180) SIND (X) es el seno de los elementos de X, >> SIND(30) expresada en grados. Al ponerlo en mayúsculas Matlab te da un mensaje de aviso, pero te devuelve el resultado, en este caso. TAN TAN tangente del argumento en radianes. Tan (x) es la tangente de los elementos de X. ASIN ASIN (X) es el arcoseno de los elementos de X. Si el argumento es mayor de 1 te devuelve el número complejo.

>> tan(pi)

ans = 1.2246e-016 ans = 1 ans = 0 ans = 0.5000

ans = -1.2246e-016

>> asin(0)

ans =

>> asin(1)

0 ans =

>> asin(sqrt(3)/2)

1.5708 ans = 1.0472

SINH >> sinh(pi) SINH (X) es el seno hiperbólico de los elementos de X. EXP >> exp(1) EXP (X) es la exponencial de los elementos de X, e elevado a X. Por complejo z = X + i * Y, EXP (Z) = exp (X) * (COS (Y) + i * sen (y)). LOG >> log(2) Log (x) es el logaritmo natural

ans = 11.5487 ans = 2.7183

ans = 0.6931

10

de los elementos de X. Se producen resultados complejos si X no es positivo. LOG10 Log10 (X) es el logaritmo en base 10 de los elementos de X. Resultados complejos se producen si X no es positivo. Se puede utilizar en mayúsculas REM REM (x, y) es x -. N * Y donde n = f (x / y) si y ~ = 0. El resto que resulta al dividir x entre y ROUND Round(X) redondea los elementos de X hacia el entero más cercano EPS EPS(X) separación de los números de punto flotante. X puede ser de doble precisión o precisión simple. Para todoX, EPS (X) es igual a EPS (ABS (X)). PI Es el número pi

>> log(-2)

ans = 0.6931 + 3.1416i

>> log10(20)

ans =

>> log10(10)

1.3010 ans = 1

>> rem(7,-2)

ans =

>> rem(2,2)

1 ans = 0

>> round(9/3.5)

ans =

>> round(1.2456)

3 ans =

>> eps(3.5)

1 ans =

>> eps(3)

4.4409e-016 ans = 4.4409e-016

PI

3.1415926535897....

4. Introducción a un script. Un script es un programa que construye el usuario de Matlab. Este programa se construye con el lenguaje de alto nivel de Matlab. Todas las asignaciones que se realicen en el script o/y en la ventana de comandos se comparten. Todo script tiene que tener un nombre y una extensión fija .m. Esta extensión indica que el script o fichero está escrito en el lenguaje de Matlab. El script se puede crear con un editor (p.e. edit) o con la pestaña File del menú superior del entorno de matlab. Cuando se haya terminado la edición grabarlo (save as) con nombre y extensión .m (la extensión normalmente se pone por defecto, es decir no hay que ponerla). Tras la grabación el programa estará en el directorio donde se ha grabado (en nuestro caso en d:\alumnos\work). Al ejecutarlo en la ventana de comandos solo se escribe el nombre del script. Fíjate en los detalles. Los apóstrofos están en la tecla del cierre de interrogación.

11

Ejemplo de script % Script ejemplo_script.m disp(´Cuánto vale A_script y B_vc en ejemplo_script´) A_script=37.5 disp([B_vc]) disp(´Abandono el script´)

Ejecución de un script en la ventana de comandos >>edit ejemplo_script.m >>B_vc=5555; >> ejemplo_script %ejecución de un script >>disp(´Cuánto vale a_script en ventana de comandos´) >>A_script

IMPORTANTE: Matlab diferencia entre las letras mayúsculas y minúsculas como ya se indicó. Luego, en el ejemplo anterior si se pregunta por A_script dará el valor asignado 37.5, mientras que si se pregunta por a_script nos dirá que esta variable no existe. IMPORTANTE: Como en todo lenguaje de alto nivel hay palabras reservadas, p.e. disp., que no se pueden usar como nombres de variable para realizar asignaciones. IMPORTANTE: El separador de la parte entera y de la parte decimal en los números reales es un punto no una coma. IMPORTANTE: El directorio, por defecto, en donde se graba un script lo puedes saber ejecutando >>pwd o >>cd. Si quieres ir al directorio padre teclea >>cd … Si quieres volver al directorio en donde trabajabas teclea >>cd nombredeldirectorio. Si quieres trabajar con tu pen usando el comando cd indica a continuación el nombre con el que se identifica tu pen. 5. Introducción básica a funciones. Una función al igual que un script la tiene que crear un programador o usuario de Matlab. En Matlab hay dos tipos de funciones, función interna, como p.e. la función ya conocida sin, y función externa o creada por el usuario de Matlab. Centrándonos en esta última, una función tiene que tener un nombre y la extensión .m, se crea con un editor y se graba en un directorio. CUIDADO: La función se debe grabar con el mismo nombre que está indicado en la primera línea (la línea donde está la palabra reservada function). Se ejecuta en la ventana de comandos. Otro ejemplo de palabra reservada es function. Vas a realizar muchas funciones a lo largo del curso; a continuación sólo se muestra un ejemplo que debes hacer. Presta atención. Ejemplo de función Ejecución de una función en la ventana de comandos function [c_fun,d_fun]=ejemplo_fun(a,b) %ejemplo_fun.m >>a=10; b=5; >> [C_fun,D_fun]=ejemplo_fun(a,b) disp('Cuánto vale a,b, c_fun y c_fun dentro >> disp(´Los valores numéricos de las de ejemplo_fun') variables c_fun y d_fun definidas disp(a) dentro de la función se han asignado a disp(b) las variables C_fun y D_fun en la ventana de comandos tras ejecutar la c_fun=a+b; función ejemplo_fun en dicha ventana d_fun=a*b; ´) >> disp([C_fun D_fun]) disp([c_fun d_fun]) % Fíjate en este disp. >>% las variables a y b son las disp('Abandono la función') denominadas variables o parámetros de entrada a la función. >>% las variables C_fun y D_fun son las denominadas variables o parámetros de salida, obtenidos tras la ejecución de la función. >> disp(´Sin embargo, en la ventana de comandos no conoce las variable

12

c_fun y d_fun, ya que únicamente están definidas dentro de la función. La siguiente ejecución así lo dirá´) >> disp([c_fun d_fun]) 6. Escribe un script que calcule la superficie y volumen de un cilindro de radio 3 m y altura 5 m. La asignación a variables escalares hazlas fuera del script. A este script llámale sup_vol_script.m. Escribe el script y las órdenes en la ventana de comandos. Script Líneas de ejecución del script en la %sup_vol_script ventana de comandos de Matlab s=pi*r*r*2+2*pi*r*h; >> r=3;h=5;sup_vol_script v=pi*r*r*h; >> s,v s= 150.7964 v= 141.3717 7. Escribe una función que tenga como variables de entrada el radio y la altura, y como variables de salida el volumen y la superficie. Llama a esta función sup_vol_fun.m. Escribe la función y las órdenes. Escribe la función Líneas de ejecución de la función en function [s,v]=sup_vol_fun(r,h) la ventana de comandos de Matlab %sup_vol_fun >> r=3;h=5; s=pi*r*r*2+2*pi*r*h; >>[sup,vol]=sup_vol_fun(r,h) v=pi*r*r*h; %disp([r,h,s,v])

sup = 150.7964 vol = 141.3717

8. Introducción a vectores. Como en todos los lenguajes de alto nivel las variables numéricas se clasifican básicamente en variables escalares (un número), variables vectoriales (vector) o variables matriciales. En ejemplos anteriores se ha visto la asignación numérica escalar, ahora vamos a ver una de las formas para asignar valores/elementos a un vector. Busca el significado de las órdenes LINSPACE usando el comando help, y ejecuta >> x=linspace(0,2*pi,3) %Linspace (X1, X2,n) crea un vector fila de n puntos entre X1 y X2, si n<2 te devuelve X2. Nota. Esta línea de comandos permite asignar valores numéricos a la variable llamada x. A diferencia de las asignaciones anteriores la variable x es un vector. Sin embargo, hay que aclarar que realmente Matlab trabaja con variables con estructura matricial,

13

en donde p.e. x será una matriz de 1 fila y 30 columnas. (a)Escribe el resultado de la asignación anterior. >> x=linspace(0,2*pi,3) x= 0 3.1416 6.2832 (b)Sea la siguiente orden >> y=sin(x) Escribe los valores numéricos de la variable y >> y=sin(x) y= 1.0e-015 * 0 0.1225 -0.2449 (c) Busca información sobre el comando o la orden plot. Ejecuta >>plot(x,y,x,y,´r*´). Como ves se ha generado una figura con los elementos de dos vectores. Plot(x,y) dibuja el vector x frente al vector y. >> plot(x,y,x,y,'r*') -16

1.5

x 10

1 0.5 0 -0.5 -1 -1.5 -2 -2.5

0

1

2

3

4

5

6

7

(d)Para determinar el número de elementos de un vector usa la función interna length. Indica la orden que has ejecutado para saber el número de elementos de los vectores x e y. >> x,y x= 0

3.1416

6.2832

y= 1.0e-015 * 0

0.1225 -0.2449

>> length(x),length(y) ans = 3

14

ans = 3 (e)La siguiente ejecución es incorrecta: >> z=x^2. Modifícala para que tengamos el cuadrado de cada uno de los elementos de vector x. Recuerda el ejercicio 1. Se le debe añadir un punto delante de ^ para que haga la operación elemento a elemento. >> z=x.^2 z= 0

9.8696 39.4784

(f) [opcional]Otra forma de asignar valores numéricos al vector con nombre x es por ejemplo: x=[10:-1:1] Ejecuta la siguiente función interna que permite construir una base de datos, s=struct('a', x, 'b', 'Nombre', 'c', int16([x;x])). Visualiza su contenido en el workspace y responde a las siguientes preguntas: ¿Que tipos de datos son x y s?. Cambia el contenido del campo b de s. >> x=[10:-1:1] x= 10 9 8 7 6 5 4 3 2 1 >> s=struct('a', x, 'b', 'Nombre', 'c', int16([x;x])) s= a: [10 9 8 7 6 5 4 3 2 1] b: 'Nombre' c: [2x10 int16] %s es una variable en la que se guardan las variables ‘a’ cuyo valor es la variable que ya teníamos guardada x, un vector fila, a ‘b’ se le asigna una palabra, y a ‘c’ una matriz formada por dos vectores x. Se crea una matriz de estructura con los campos y los valores especificados >> s=struct('a', x, 'b', x, 'c', int16([x;x])) s= a: [10 9 8 7 6 5 4 3 2 1] b: [10 9 8 7 6 5 4 3 2 1] c: [2x10 int16] >> s=struct('a',int16([y;y]),'b', y, 'c', int16([x;x])) s= a: [2x1 int16] b: 9 c: [2x10 int16] 9. Indica los resultados tras ejecutar las siguientes órdenes. Si hubiera errores corrígelos. Fíjate en las órdenes que han dado lugar a matrices o variables matriciales. >> X=[ 1 2 9 4 5]

X=

15

ans = 1 2 9 4 >> X=[1,2,9,4,5] X= 1

2

9

4

5 1 2 5

>>X=20:-2:1 x=20:-2:1

>> X=[ 1;2;9;4;5] X=

x=

1 2 9 4 5

2

9

>>X(3) %le término 3 ans =

4 está

5 pidiendo

el

9

18 2

16

14

12

8

7 6

8 7

Columns 17 through 32 8 9 10 11 12 3 5 6 7 8 9 10 12 13

4 11

2 >>longitud=length(X) Longitud= 10 >>tamano =size(X) tamano = 1

10

>>[fil,col]=size(X) fil =

>>X(0) Matlab te dice que necesita un número entero y lógico, el cero para él no lo es. Probamos con X(1), ya que no existe el término cero del vector. >> X(1) ans =

1 col = 10

>>a=[1:10,2:12,3:13] a=

1 >>X([1,3,5]) %debe tener un ‘=’ X=([1,3,5]) X= 1

20 4

>>X(end) ans =

>>X' ans = 1

10

6

3

5

>>X(1:4) da error dado que el vector X solo tiene 3 elementos. X(1:4) significa que te de los elementos del vector del 1 al 4 con un salto (por defecto) de 1. Para que nos de una respuesta ponemos X(1:2): >> X(1:2)

Columns 1 through 16 1 9

2 10

3 2

4 3

5 4

6 5

>> a=[1 2 7 568 3 2 6] a= 1

2

7 16

5 3

6 2

8 6

d= 1 4 6

>> b=[1 2 3] b= 1

2

3

2 6 2 2

3 0 1

>>e=[1 2 3;4 3 0;6 8 1] e=

>> c=[a;b] c= 1 5 3 1

2 3 8

1 4 6 7 8 6 3

2 3 8

3 0 1

>>f=e’, g=b’; f= 1 2 3

>>d=[1,2,3;4,3,0;6,8,1]

4 3 0

6 8 1

SECCIÓN 2. Variables 1. Crear un vector x de 4 componentes equi-espaciado entre los valores 6 y 7. >> x=linspace(6,7,4) x= 2.

6.0000 6.3333 6.6667 7.0000 Sumar 1 al tercer elemento del vector x. >> x(3)=x(3)+1 x=

6.0000 6.3333 7.6667 7.0000 3. Crear un vector y de las mismas dimensiones que x con los primeros números impares. y=linspace(1,7,4) y= 1

3

5

7

4. Crea un vector v que contenga la tabla de multiplicar del 3, cuyo primer elemento sea cero y el último sea 30. >> v=(0:3:30) v= 0 5.

3

6

9

12

15

18

21

24

27

30

Crear un vector y con la diferencia entre sucesivos elementos de v. y=linspace(v(1)-v(2),v(end-1)-v(end),length(v)-1) y=

17

-3

-3

-3

-3

-3

-3

-3

-3

-3

-3

6. Sumar 5 a los elementos del vector v con índice par 5.+[v(1),v(3),v(5),v(7),v(9),v(11)] ans = 5

11

17

23

29

35

7. Crear los vectores a=[0 6 8 3] y b=[-1 7 8*2 4]. (a) Crear una matriz A de tamaño 2x4 cuyas filas sean los vectores a y b y (b) Crear una matriz 4x2 cuyas columnas sean a y b. >> A=[a;b] A= 0 -1

6 7

8 16

3 4

>> B=A' B= 0 6 8 3

-1 7 16 4

8. Escribe un vector z que vaya del 0 al 1 con un espaciado de 0.1. Calcular un vector con los valores de la función interna sin aplicada a los elementos de z. Hazlo también para la función interna exp. >> z=(0:0.1:1) z= Columns 1 through 7 0

0.1000

0.2000

0.3000

0.4000

0.5000

0.6000

0.4794

0.5646

Columns 8 through 11 0.7000

0.8000

0.9000

1.0000

>> g=sin(z) g= Columns 1 through 7 0

0.0998

0.1987

0.2955

0.3894

Columns 8 through 11 0.6442

0.7174

0.7833

0.8415

>> h=exp(z)

18

h= Columns 1 through 7 1.0000

1.1052

1.2214

1.3499

1.4918

1.6487

1.8221

Columns 8 through 11 2.0138

2.2255

2.4596

2.7183

9. Crear un vector x con las fracciones 1, 1/2 , 1/3, 1/4, … 1/10. v=1./[1:1:10] v= Columns 1 through 8 1.0000

0.5000

0.3333

0.2500

0.2000

0.1667

0.1429

0.1250

Columns 9 through 10 0.1111 0.1000 function n=[1:1:10] V=1./n 10. Generar los primeros 10 números de la serie: 0, 1/2, 2/3, 3/4, … w=[0:1:9]./[1:1:10] w= Columns 1 through 8 0

0.5000

0.6667

0.7500

0.8000

0.8333

0.8571

0.8750

Columns 9 through 10 0.8889 11.

0.9000

Crear la serie (1) n con n= 0, …, 10

n=(0:1:10),(-1).^n n= 0

1

2

3

4

5

6

7

8

9

10

-1

1

-1

1

-1

1

-1

1

-1

1

ans = 1

19

( 1) n desde n = 0, …, 10 . Usa la 2n  1 función interna sum para sumar todos los elementos del vector. >> n=(0:1:10);w=((-1).^n)./(2.*n+1),s=sum(w) 12.

Crear un vector con los números de la serie:

w= Columns 1 through 7 1.0000 -0.3333

0.2000 -0.1429

0.1111 -0.0909

0.0769

Columns 8 through 11 -0.0667

0.0588 -0.0526

0.0476

s= 0.8081 Veamos cómo calcular un valor aproximado de  . Suma los primeros 100  (1) n 1 1 1 1 1        ...  . Repetir el paso términos de la fórmula de Leibniz:  1 3 5 7 9 4 n  0 2n  1 anterior con 10000 términos. ¿Cuál de los dos resultados está más próximo al valor de  ?. >> n=(0:1:100);w=((-1).^n)./(2.*n+1);s=sum(w) 13.

s= 0.7879 >> n=(0:1:10000);w=((-1).^n)./(2.*n+1);s=sum(w) s= 0.7854 14. Crear una matriz B cuya primera fila corresponda al vector a, cuya segunda fila sea un vector de ceros, la tercera un vector de unos y la cuarta el contenido de b. Ver ones y zeros. >>length(a),length(b) ans = 4 ans = 4 >> B=[a;zeros(1,4);ones(1,4);b] B= 0

6

8

3

20

0 1 -1

0 1 7

0 1 16

0 1 4

15. Definir un vector k cuyos elementos sean iguales a los dos primeros elementos de la diagonal de la matriz B. Utiliza la función interna diag. k=[diag(B(1)),diag(B(2))] k= 0

0

16. Crea una submatriz C de tamaño 3x3 que sea el contenido de las tres primeras filas y columnas de la matriz B. >> C=[B(1:3,1:3)] C= 0 0 1 17.

6 0 1

8 0 1

Crear dos matrices A y B: (a) Resultados de las dos operaciones suma y resta entre matrices: C = A+B D=A-B. >> A=[1 2;4 -1];B=[4 -2;-6 3]; >> C=A+B,D=A-B

y

C= 5 -2

0 2

D= -3 4 10 -4 (b) Calcula los determinantes de ambas matrices y sus inversas. >>det(A),det(B) ans = -9 ans = 0 >>inv(A) ans = 0.1111 0.2222 0.4444 -0.1111 inv(B)

21

Warning: Matrix is singular to working precision. ans = Inf Inf Inf Inf %Nos advierte de que B es una matriz singular, es decir su determinante es cero y por lo tanto no tiene inversa. (c) Resultados de: E=A*B y F=B*A. >> E=A*B,F=B*A E= -8 4 22 -11 F= -4 10 6 -15 (d) Empleando la operación elemento a elemento obtén la matriz G que sale de la siguiente expresión. G=A-B*A2/3 Si observas algo extraño averigua el por qué. >> G=A-B.*A.^(2/3) G= -3.0000 5.1748 19.1191 0.5000 - 2.5981i %En esta matriz podemos observar que en la columna 2 y 2 fila el número es complejo esto se debe a que A(2,2)=-1 y (-1)^(2/3) es un número complejo. >> (A(2,2))^(2/3) ans = -0.5000 + 0.8660i (e) Tiene sentido realizar las siguientes operaciones A/B y A\B? Razona tu respuesta >> A/B Warning: Matrix is singular to working precision. ans = Inf Inf Inf Inf >> A\B ans = -0.8889 0.4444 2.4444 -1.2222 No tiene ningún sentido la primera operación, porque significa

que

22

multiplicamos la matriz A por la inversa de B (B^(-1)) y anteriormente hemos dicho que las matrices cuyos determinantes son cero no tienen inversa, por lo que esta operación no tiene sentido, mientras que A\B si se puede realizar porque aquí multiplicamos B por la inversa de A, ambas matrices cuadradas, y esta última tiene un determinante distinto de cero. 18. Resolución de sistemas de ecuaciones. Sea el siguiente sistema de ecuaciones de la forma A·x=b, resolverlo mediante los operadores matriciales estudiados en clase

>> A=[1 3 6;8 9 1;3 10 9]; >> b=[1;3;4]; >>x=b/A x= 0.4176 19.

0.0110

0.1648

Veamos ahora el indexado lógico de matrices y vectores: Sea el vector x=[-1 0 2 4 -2 3 1 4 7 0 -3 -1]; . Si queremos asignar a cero los elementos de x mayores de cero hacemos en la línea de comandos x(x>0) = 0. Recuerda el ejercicio 1 de la sección 1. (a) Poner a -5 los elementos de x menores de 0 >> x(x<0)=-5 x= -5

0 2 4 -5 3 1 4 7 0 -5 -5 Poner a 3 los elementos de x mayores de 1 y menores de 5: >> x((x>1)&(x<5))=3

(b) x=

-1

0

3

3

-2

3

1

3

7

0

-3

-1

(c) Crear un vector y con los elementos de x menores o iguales a 3: >> y=x((x<=3)) y= -1 (d)

0

3

3

-2

3

1

3

0

-3

-1

Sumar 2 a los elementos mayores de 5 >> x(x>5)+2

ans = 9 (e) Crear un nuevo vector y que tenga 1's en las posiciones de los elementos de x que sean mayores que la media (usar la función interna mean) y 0's en las de los menores que la media. NOTA: primero haz >>y=x %para evitar modificaciones del vector x y, por tanto, también de su media. >> x(xmean(x))=1

23

x= 0

0

5

5

0

5

5

5

7

0

0

0

x= 0

0

1

1

0

1

1

1

1

0

0

0

(f) Asignar a cero los elementos pares (usar función interna rem) >> x(rem(x,2)==0)=0 x= -1 0 0 0 0 3 1 0 7 0 -3 -1 (g) Cambiar el signo de los valores de x verificando: 2<=x(i) < 5: >> x=(-1)*(x(2<=x<5)) x= 1 0 0 0 0 -3 -1 SECCIÓN 3. Funciones y scripts.

0

-7

0

3

1

1. Crear la función suma_vectores que tome como parámetros de entrada los vectores a y b, devolviendo un vector c que realice su suma. Hazlo también con un script, llámalo suma_vectores_s.. ¿Que diferencia aprecias entre suma_vectores y suma_vectores_s? ¿Las variables tienen el mismo ámbito? ¿La forma de invocarlos es la misma?. Función: function [c]=suma_vectores(a,b) c=a+b;

>> a=[1 2 3];b=[2 3 1]; >> c=suma_vectores(a,b) c= 3 5 Scrip:

4

%suma_vectores_s c=a+b

>> a=[1 2 3];b=[2 3 1]; >> suma_vectores_s c= 3

5

4

La diferencia entre la función y el scrip, la diferencia se encuentra en la sintaxis, mientras que la función necesita variables de entrada y salida, el scrip es una secuencia de comandos. Además la función no está conectada con Matlab, esto quiere decir que si le preguntas a Matlab por una variable de la función te va a decir que está indefinida, o si a la función le pides una variable que has dado en la línea de comandos y no en la variable de entrada la función no la reconoce, sin embargo un scrip reconoce las variables del entorno de Matlab.

24

Las variables en los scrips son de ámbito global, mientras que las de las funciones son de ámbito local. No se invocan de la misma manera, para invocar un scrip asignamos las variables y llamamos al scrip (>>suma_vectores_s), pero en la función debes escribir el nombre del fichero y las variables de entrada.

2.

Escribir una función que reciba como argumento un valor x y un entero n y

devuelva el valor de la función matemática:

f (x ) 

1 1 x

n

. Construye una figura que

muestre cómo varía f(x) si x=3 y n varía entre -40 y 40. Pon etiquetas a los ejes usando los comando xlabel e ylabel. %abs(x)es el valor absoluto de x elevado a n %recibe como argumento un valor x y un entero n, y devuelve el valor de fx Edit: function [fx]=argumento(x,n) fx=1./(1+((abs(x)).^n))

Línea de comandos: >> n=[-40:40]; >> x=3; >> fx=argumento(x,n); fx = Columns 1 through 9 1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

1.0000

0.9999

0.9998

0.9995

0.9986

0.9959

0.7500

0.5000

0.2500

0.1000

0.0357

0.0122

0.0002

0.0001

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

Columns 10 through 18 1.0000

1.0000

1.0000

Columns 19 through 27 1.0000

1.0000

1.0000

Columns 28 through 36 1.0000

1.0000

1.0000

Columns 37 through 45 0.9878

0.9643

0.9000

Columns 46 through 54 0.0041

0.0014

0.0005

Columns 55 through 63 0.0000

0.0000

0.0000

Columns 64 through 72

25

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

Columns 73 through 81 0.0000

0.0000

0.0000

>> plot(fx,n) >> xlabel('f(x) en x=3') >> ylabel('n entre -40 y 40') 40

30

20

n entre -40 y 40

10

0

-10

-20

-30

-40 -1

-0.5

0

0.5 f(x) en x=3

1

1.5

2

3. Crear una función llamada alcance.m que determine el alcance máximo horizontal (xmax), de un objeto lanzado con una velocidad inicial vo en el plano xy. Considera como variables de entrada el ángulo θ y la velocidad de origen vo del lanzamiento del objeto, ambas serán variables escalares. La variable de salida será el alcance máximo.

Edit: function [xmax]=alcance(z,vo) %z es el ángulo %vo es la velocidad de origen del lanzamiento g=9.8; %alcance.m xmax=(vo^2*sin(2*z))/g; disp([xmax,z,vo]) %disp('determina el alcance máximo horizontal, de un objeto lanzado con una %velocidad inicial vo en xy') %disp('variable entrada z, ángulo en radianes, velocidad origen vo') >>z=0.7854; vo=40; >> [xmax]=alcance(z,vo) 163.2653 0.7854 40.0000

26

xmax = 163.2653 Modificamos la función alcance: >> vo=linspace(40,90,3); >> alcance(z,vo) ??? Error using ==> mpower Matrix must be square. Error in ==> alcance at 5 xmax=(vo^2*sin(2*z))/g; %Lo que hacemos para que nos deje de dar error es ponerle delante de los signos de multiplicación, potencia y división %un punto para que la operación se haga elemento a elemento >> alcance(z,vo) 163.2653 431.1224 826.5306

0.7854 40.0000 65.0000 90.0000

ans = 163.2653 431.1224 826.5306 Longitud: >> length(vo) ans = 3 >> length(z) ans = 1

4. Haz lo mismo que en el ejercicio 3 con un script, llámale alcance_s.m Edit: %alcance_s g=9.8; %alcance.m xmax=(vo.^2.*sin(2.*z))./g; disp('xmax');disp([xmax]) %disp('determina el alcance máximo horizontal, de un objeto lanzado con una %velocidad inicial vo en xy') %disp('variable entrada z, ángulo en radianes, velocidad origen vo')

Línea de comandos: >> z=0.7854; vo=40; >> alcance_s xmax 163.2653 431.1224 826.5306 5. Crea una función tiempo.m que determine el tiempo máximo del lanzamiento. Variables de entrada: ángulo y velocidad inicial. Variable de salida: tiempo máximo.

Edit: function [tm]=tiempo(z,vo)

27

g=9.8; tm=(2.*vo.*sin(z))./g; %z es el ángulo en radianes %vo es la velocidad inicial %tm es la variable de salida (tiempo máximo), de entrada el ángulo en radianes (a) y la velocidad inicial en m/s(vo) >> z=pi/4;vo=40; >> [tm]=tiempo(z,vo) tm = 5.7723

6. Crear una función llamada alcance_tiempo.m que llame a las funciones alcance.m y tiempo.m. Las variables de salida serán el alcance y el tiempo máximo. function xtmax(z,vo) %disp('hemos creado una función, llamando a otras dos funciones, posible porque se encuentran en el mismo directorio'); [xmax]=alcance(z,vo) [tm]=tiempo(z,vo) xmax2=vo*cos(z)*tm %xmax2 fue añadido en clase, un ejercicio añadido. Línea de comandos: >>z=pi/4;vo=40;

>> xtmax(z,vo) 163.2653 0.7854 40.0000 xmax = 163.2653 tm = 5.7723 xmax2 = 163.2653 7. Crear una función llamada alcance_vect.m con variable de salida una matriz que contenga todos los ángulos, velocidades iniciales, tiempos en alcanzar de nuevo el suelo y los alcances máximos de todas las trayectorias generadas. Las variables de entrada serán el ángulo y la velocidad.(a) Aplícalo a un ángulo de 45º y una velocidad de entre 40 y 60 m/s con un salto de 10 m/s. (b) Aplícalo a ángulos entre 10º y 90º con un salto de 5º. (c) Aplícalo a una velocidad entre 20 y 40 m/s con un salto de 10 m/s y ángulos entre 10º y 30º con un salto de 10º. (d) Aplícalo a los ángulos de 10º a 90º saltando de 5º en 5º y una velocidad inicial del lanzamiento de 40 m/s a 60 m/s con un salto de 10 m/s. Como verás no es posible. Este último apartado lo tendrás que hacer con bucles en el ejercicio 11 de la Sección 4, a este programa llámale alcance_bucle.m. a) Edit: function [A]=alcance_vect(z,vo) disp('xmax,tm,vo,z')

28

g=9.8; [xmax]=alcance(z,vo); [tm]=tiempo(z,vo); A=[xmax,tm,vo,z]; %llamamos a otras funciones Línea de comandos >> z=45; >> [A]=alcance_vect(z,vo) xmax,tm,vo,z 145.9586 228.0604 328.4069 45.0000 40.0000 50.0000 60.0000 A= Columns 1 through 9 145.9586 228.0604 328.4069

6.9462

8.6827 10.4192 40.0000 50.0000 60.0000

Column 10 45.0000 b) Le falta la velocidad de inicio, ya que si usamos la del apartado anterior estaríamos realizando el apartado (d), por lo que supondremos que la velocidad de inicio es 40m/s. Edit: function [A]=alcance_vect(z,vo) %disp('xmax,tm,vo,z') g=9.8; [tm]=tiempo(z,vo); [xmax]=alcance(z,vo); A=[xmax tm vo z]; Línea de comandos: >> [A]=alcance_vect(z,vo) Columns 1 through 9 149.0523 -161.3113 121.6511 -42.8367 -49.7650 126.3495 -162.2675 145.9586 -82.6719 Columns 10 through 18 -7.2233 94.7937 -151.8540 160.0391 -116.7145 35.8245 56.5958 -130.8004 10.0000 Columns 19 through 27 15.0000 20.0000 25.0000 30.0000 35.0000 40.0000 45.0000 50.0000 55.0000 Columns 28 through 35 60.0000 65.0000 70.0000 75.0000 80.0000 85.0000 90.0000 40.0000 A= Columns 1 through 9 149.0523 -161.3113 121.6511 -42.8367 -49.7650 126.3495 -162.2675 145.9586 -82.6719 Columns 10 through 18

29

-7.2233 94.7937 -151.8540 160.0391 -116.7145 35.8245 56.5958 -130.8004 -4.4410 Columns 19 through 27 5.3085

7.4526 -1.0804 -8.0656 -3.4954

6.0826

6.9462 -2.1418 -8.1613

Columns 28 through 36 -2.4882

6.7496

6.3175 -3.1656 -8.1134 -1.4374

7.2979 40.0000 10.0000

Columns 37 through 45 15.0000 20.0000 25.0000 30.0000 35.0000 40.0000 45.0000 50.0000 55.0000 Columns 46 through 52 60.0000 65.0000 70.0000 75.0000 80.0000 85.0000 90.0000

c) Para este apartado tuve que buscar un comando que crease una matriz por bloques de dimensiones, “repmat(A,[m,n])” con copias de la matriz A. Edit: function[A]=alcance_vect1(vo,z) g=9.8; r=z.*pi./180; vo1=repmat(vo,1,length(r)); z1=repmat(r,1,length(vo)); xmax=((vo1.^2).*(sin(2.*z1)))./g; tm=(2.*vo1.*sin(z1))./g; A=[xmax; vo1; z1]; Línea de comandos: >> [A]=alcance_vect1(z,vo) A= Columns 1 through 9 6.5591 19.8832 40.1962 40.9941 79.5329 123.1010 104.9449 178.9491 251.2265 10.0000 15.0000 20.0000 25.0000 30.0000 35.0000 40.0000 45.0000 50.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 Columns 10 through 18 198.4115 318.1318 424.5727 321.3938 497.0809 643.1398 473.8919 715.7965 10.0491 55.0000 60.0000 65.0000 70.0000 75.0000 80.0000 85.0000 90.0000 10.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 Columns 19 through 27 14.7579 35.3480 62.8066 59.0315 108.2532 160.7849 132.8209 220.9248 303.9840 15.0000 20.0000 25.0000 30.0000 35.0000 40.0000 45.0000 50.0000 55.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 Columns 28 through 36 236.1261 373.3630 492.4039 368.9470 565.5676 726.0445 531.2836 8.8370 22.6104 60.0000 65.0000 70.0000 75.0000 80.0000 85.0000 90.0000 10.0000 15.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 Columns 37 through 45

30

26.2362 55.2312 90.4415 80.3485 141.3919 203.4934 163.9764 267.3191 361.7661 20.0000 25.0000 30.0000 35.0000 40.0000 45.0000 50.0000 55.0000 60.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981 Columns 46 through 51 277.1202 433.0127 565.2596 419.7797 638.4728 813.9738 65.0000 70.0000 75.0000 80.0000 85.0000 90.0000 0.3491 0.5236 0.6981 0.3491 0.5236 0.6981

SECCIÓN 4. Bucles y condicionales Se tiene que usar bucles o/y condicionales en todos los programas de esta sección. 1. Escribe una función potencia.m que reciba como argumentos de entrada un número x y un entero n y devuelva el valor de la n-ésima potencia de x (x^n). El segundo argumento de entrada puede ser opcional. Si se omite debe tomarse como n=2 (devolviendo el cuadrado de x). Ver la función interna nargin. Usa la estructura if. EDIT: function [nesima]=potencia(x,n) if nargin==1 n=2 end nesima=(x^n)

Líneas de commando: >> potencia(3) n= 2 nesima = 9 ans = 9 2. Crea una función que determine la suma de los inversos de los impares menores de 1000. Usa las funciones tic y toc para determinar el tiempo que se tarda. Usa la estructura for. 1 1 1 1 Serie  1     ...  3 5 7 999 function [serie]=inversos(N) serie=0; tic for i=1:2:N; serie=serie+(1/i); end toc >> inversos(1000) Elapsed time is 0.000005 seconds. ans = 4.0891

3. Los números de Fibonacci, F(n) son una serie de números enteros, inicializados con F(1)=1 y F(2)=2. Los demás términos se tienen como la suma de los dos

31

anteriores: F(n) = F(n-1) + F(n-2). Escribe un script que genere los primeros 100 números de la serie de Fibonacci, almacenándolos en un vector F. Usa la estructura while. %script=fibonacci_s f(1)=1; f(2)=2; n=3; N=100; %clear f while (n> fibonacci_s N= 100 f= 1.0e+020 * Columns 1 through 9 0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0000

0.0001

0.0001

0.0001

0.0002

0.0004

0.0006

Columns 10 through 18 0.0000

0.0000

0.0000

Columns 19 through 27 0.0000

0.0000

0.0000

Columns 28 through 36 0.0000

0.0000

0.0000

Columns 37 through 45 0.0000

0.0000

0.0000

Columns 46 through 54 0.0000

0.0000

0.0000

Columns 55 through 63 0.0000

0.0000

0.0000

Columns 64 through 72 0.0000

0.0000

0.0000

Columns 73 through 81 0.0000

0.0000

0.0000

Columns 82 through 90

32

0.0010

0.0016

0.0026

0.0042

0.0068

0.0110

0.0178

0.0288

0.0466

0.3194

0.5168

0.8362

1.3530

2.1892

3.5422

Columns 91 through 99 0.0754

0.1220

0.1974

4. Escribir las líneas de código para que dado un vector x, genere otro vector con el orden de los elementos invertido. Edit: %vector_s for i=0:1:length(x); if i~=length(x); n=x(length(x)-i); %n=x(end)-(x-i); y=n; disp(y) end end Línea de comandos: >> x=[1 2 -4 6 7] x= 1

2

-4

6

7

>> vector_s 7 6 -4 2 1

5. Implementar una función que tome como entrada un vector x y que devuelva la posición y el contenido del primer elemento negativo. Sintaxis de llamada: [elem, pos] = buscaelem(x). (a) Añadir una condición de parada si no hubiese ningún elemento negativo, que muestre por pantalla el mensaje: No hay ningún elemento negativo Edit: function [elem,pos]=bucaelem(x) %x es un vector for i=1:1:length(x) if x(i)>=0 disp('no hay ningun elemento negativo') else elem=x(i) pos=i break end end Línea de comandos: >> x=[1 -3 5 -3];buscaelem(x) no hay ningun elemento negativo elem =

33

-3 pos = 2 ans = -3

6. Implementar una función que dado un vector devuelva el menor valor del vector y la posición donde se encuentra sin usar la función min o max. Por ejemplo v=[ 1 2 3 4 5 6 -1 -3 -2] devolvería el valor de -3 y la posición 8. function maximoyminimo(v) for i=1:length(v); if v(i)<=v; disp(v(i)) disp(i) end if v(i)>=v; disp(v(i)) disp(i) end end

>> maximoyminimo(v) 6 6 -3 8 7. Crear una función que dados los valores escalares de n, a 1, vector de n componentes tal que a n = an-1 - 3an-2. function [x]=vector(n,a1,a2) x(1)=a1; x(2)=a2; N=10; n=3; length(x)=n; %clear f while n> vector(3,1,4) 1 4 1 -11 -14 19 61

a2 devuelva un

4 -179

8. Crear un script que obtenga el valor medio y desviación media del vector x descrito como:

Compara los resultados con las funciones internas mean y std. Compara ahora el resultado de std con tu variable desv si la divides entre (n-1) en vez de n.

34

%mediaydesviacion_s n=length(x); suma1=0; for i=1:1:n suma1=suma1+x(i); end media=suma1/n desviacion = sqrt(sum((x-media).^2)/n) >> x=[1 2 3 9 8 7]; >> mediaydesviacion_s >> n,media,desviacion n= 6 media = 5 desviacion = 3.1091 %Lo comparamos >> mean(x),std(x) ans = 5 ans = 3.4059 (no da exactamente la desviación igual, ya que Matlab no usa la misma fórmula que nosotros hemos utilizado)

9. Implementar una función que tome como entrada una matriz A y que devuelva la suma de todos los elementos de las 2 diagonales de dicha matriz. No utilizar la función diag(). Edit: function [z]=matriz(A) [nf,nc]=size(A); n=0; p=0; %solo matrices cuadradas for i=1:1:nf n=n+A(i,(nc-i+1)); for j=1:1:nc if i==j p=p+A(i,j); end end end z=n+p Líneas de comando: >> A=[1 2 3;4 5 6;7 8 9]; diagon(A) z= 30

10. Implementar una función que dadas dos matrices A y B definidas como (a ij)mxn y (bij)nxp, implemente la multiplicación de matrices C=A·B mediante bucles, definida como (cij)mxp donde

.

function [C]=matriz2(A,B) [nf1,nc1]=size(A); [nf2,nc2]=size(B); if nc1==nf2 for i=1:1:nf1 for j=1:1:nc2 C(i,j)=sum(A(i,:).*B(:,j)');

35

end

end

end C=A*B >> A=[2 3;4 6]; >> B=[1 4;2 7]; >> matriz2(A,B) 8 29 16 58

11. Crear una función llamada alcance_bucle.m que determine el alcance máximo horizontal y el tiempo máximo de un objeto lanzado con una velocidad inicial vo y un ángulo θ según lo indicado en el apartado d del ejercicio 7 de la sección anterior. Edit: function[xmax,tm]=alcance_bucle(vo,z1) g=9.8; r=z1.*pi./180; for i=1:length(vo) for j=1:length(r) xmax(i,j)=((vo(i)^2)*(sin(2*r(j))))/g; tm(i,j)=2*(vo(i))*(sin((r(j))))/g; end end Línea de commandos:

>> z1=[10:5:90];vo=[40:10:60]; >> [xmax,tm]=alcance_bucle(vo,z1) xmax = Columns 1 through 9

55.8400 81.6327 104.9449 125.0685 141.3919 153.4192 160.7849 163.2653 160.7849 87.2500 127.5510 163.9764 195.4195 220.9248 239.7175 251.2265 255.1020 251.2265 125.6401 183.6735 236.1261 281.4041 318.1318 345.1932 361.7661 367.3469 361.7661 Columns 10 through 17 153.4192 141.3919 125.0685 104.9449 81.6327 55.8400 28.3507 0.0000 239.7175 220.9248 195.4195 163.9764 127.5510 87.2500 44.2980 0.0000 345.1932 318.1318 281.4041 236.1261 183.6735 125.6401 63.7891 0.0000 tm = Columns 1 through 9 1.4175 1.7719 2.1263

2.1128 2.6410 3.1692

2.7920 3.4900 4.1880

3.4499 4.3124 5.1749

4.0816 5.1020 6.1224

4.6823 5.8528 7.0234

5.2472 6.5591 7.8709

5.7723 7.2154 8.6585

7.6710 9.5887

7.8851 8.0392 8.1322 8.1633 9.8564 10.0491 10.1653 10.2041

6.2534 7.8168 9.3801

Columns 10 through 17 6.6870 8.3587

7.0696 8.8370

7.3984 9.2480

36

10.0304 10.6044 11.0976 11.5064 11.8277 12.0589 12.1983 12.2449 12.

Haz el ejercicio 17.d de la Sección 2.

Empleando la operación elemento a elemento obtén la matriz expresión. G=A-B*A2/3

G que sale de la siguiente

%matriz3 [nf,nc]=size(A); %deben ser cuadradas y de las mismas dimensiones, porque si no va a dar error if size(A)==size(B) for i=1:nf for j=1:nc G(i,j)=(A(i,j))-(B(i,j))*(A(i,j))^(2/3); end end end disp(G) línea de comandos: >> matriz3 -3.0000 5.1748 19.1191 0.5000 - 2.5981i >> G=A-B.*A.^(2/3) G= -3.0000 19.1191

5.1748 0.5000 - 2.5981i

SECCIÓN 5. Depurador (esta sección es opcional) 1. Utiliza la funcionalidad del depurador de Matlab estudiando el comportamiento del código de la función calcula_phi.m. Descarga del campus virtual calcula_phi.m y dibuja_valores.m. Ejecuta la calcula_phi(0.01) a. Añade un punto de ruptura (breakpoint) en la línea 13. b. Ejecuta una introducción pasa a paso (Step). ¿Cómo cambia el valor de phi y por qué? c. ¿Como va evolucionando el valor de E y por qué? d. Coloca un breakpoint en la línea 16: “dibuja_valores(iter,phi,'b*')” i. ¿Qué diferencia observas entre Step, Step In y Step Out? ii. ¿Y con Continue? SECCIÓN 6. Representación gráfica 1. Escribe la función dibujo_parabola.m que muestre paso a paso en una figura la trayectoria del tiro parabólico para vo 40 m/s y θ=45º y para vo igual a 60 m/s y θ=60º. En cada instante se mostrará por un circulo rojo la posición (x,y). Llama a alguna de las funciones anteriores. Usa la función plot y hold on

(a) Añade el nombre de cada eje con xlabel, ylabel y el título de la figura con title (b) Emplea los comandos de la ventana que contiene el gráfico para cambiar el grosor y el color de las líneas de tres trayectorias. (no hacer, poner para que

37

la pelota aparezca poco a poco, dicho en clase) Editor: function [x,y]=dibujo_parabola(vo,z)%O en grados n=30; %n=al número de puntos g=9.8; tmax=(2*vo*sind(z))/g; for t=linspace(0,tmax,n); x=vo*cos(z)*t; y=vo*sin(z)*t-1/2*g*t^2; plot(x,y,'r o') hold on pause(tmax/n) end end Línea de comandos: >> dibujo_parabola(40,45),dibujo_parabola(60,60) ans = 163.2653 ans = 318.1318 Tiro parabólico

350 300

altura

250 200 150 100 50 0

2.

0

50

100

150

Longitud

200

250

300

350

Usa únicamente sola figura para representar en diferentes panales (usa subplot)

las siguientes funciones. En el panel 1, dibujar la función en el intervalo [-2 2] empleando 100 datos equiespaciados y en los paneles 3 y 5, respectivamente, las funciones y . Cada función debe ir en un color y en un estilo de línea diferente. Añade sólo en la figura del panel 5 el comando grid on. En el panel 6 se representará la última función en escala logarítmica en uno de los ejes. En el panel 2 se representará la primera función empleando stem. En el panel 4 se representará la segunda función empleando hist. ¿Tiene sentido lo que observas? Editor: function [y]=figura(x) y=exp(-x.^2); subplot(3,2,1)

38

plot(x,y,'y*') y=exp(-x.^2) subplot(3,2,2) stem(x,y,'k') y=exp(-(x./2).^2) subplot(3,2,3) plot(x,y,'g.') y=exp(-(x./2).^2) subplot(3,2,4) hist(y,x) y=exp(-(2.*x).^2) subplot(3,2,5) plot(x,y,'bp') grid on y=exp(-(2.*x).^2) subplot(3,2,6) semilogx(x,y,'cd') Línea de comandos: >>x=linspace(-2,2,100);figure(x) 1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0 -2

-1.5

-1

-0.5

0

0.5

1

1.5

2

0 -2

1

20

0.8

15

0.6

10

0.4

5

0.2 -2

-1.5

-1

-0.5

0

0.5

1

1.5

2

0 -3

1

1

0.8

0.8

0.6

0.6

0.4

0.4

0.2

0.2

0 -2

-1.5

-1

-0.5

0

0.5

1

1.5

2

0 -2 10

-1.5

-2

-1

-0.5

-1

-1

10

0

0

0.5

1

1

1.5

2

0

10

2

3

1

10

%la última barra del histograma no tiene mucho sentido dado que no le estamos pidiendo eso

3. Crear una función dibujatriangulo que tenga como argumento de entrada los 3 vértices. Editor: function figuratriangulo(a,b,c) x=[a(1),b(1),c(1),a(1)]; y=[a(2),b(2),c(2),a(2)]; plot(x,y) area(x,y,'facecolor','c')

39

2 1.8 1.6 1.4 1.2 1 0.8 0.6 0.4 0.2 0

1

1.2

1.4

1.6

1.8

2

2.2

2.4

2.6

2.8

3

4. Crear una función dibujaMUC que dibuje la trayectoria circular de una partícula conociendo su velocidad angular ω y el radio de giro r. Recuerda el significado del periodo de un movimiento circular. Las ecuaciones cartesianas son:

Editor: function dibujomcu(r,w) t1=(2*pi)/w; for t=linspace(0,t1,100) x=r*cos(w*t); y=r*sin(w*t); plot(x,y,'ko') hold on pause(t/100) end Línea de comandos: >> dibujomcu(5,62.8) Figure: 10 8 6 4 2 0 -2 -4 -6 -8 -10 -10

-8

-6

-4

-2

0

2

4

6

8

10

(a) Añadir a la trayectoria en cada punto el vector velocidad con el comando quiver, mediante las ecuaciones de la velocidad:

Editor: function dibujomcu(r,w) n=15; t1=(2*pi)/w; for t=linspace(0,t1,n) x=r.*cos(w.*t); y=r*sin(w*t); vx=-r.*w.*sin(w.*t);

40

vy=r.*w.*cos(w.*t); plot(x,y,'o') hold on quiver(x,y,vx,vy) end Línea de comandos: >> dibujomcu(5,5) Figure: 25 20 15 10 5 0 -5 -10 -15 -20 -25 -25

-20

-15

-10

-5

0

5

10

15

20

25

5. Crear una función dibujatrayectoria3D que tenga como parámetros de entrada el radio r y la velocidad angular ω sabiendo que la partícula se mueve en el eje z con una velocidad constate vz.. function dibujotrayectoria3d(r,w) %n puntos %N vueltas %v es la velocidad en z n=100; N=4; v=8; tm=(2*pi*N)/w; t=linspace(0,tm,n); z=v.*t; y=r.*sin(w.*t); x=r.*cos(w.*t); plot3(x,y,z,'g') hold on

12 10 8 6 4 2 0 20 10 0 -10 -20

-20

-10

0

10

20

6. Crear una función dibujatipoparabolico que represente el tiro parabólico en 3D sabiendo que la función tiene como entrada los ángulos φ y θ de acuerdo a las

41

siguiente expresiones:

A la hora de representar la función, emplea el comando plot3 y pause para ir viendo como se dibuja la trayectoria parabólica. Calcula el alcance máximo en cada eje previamente y fija de antemano los valores de los ejes empleando los comandos xlim ,ylim y zlim respectivamente. function dibujotipoparabolico(j,i) %i y j son los ángulos de entrada g=9.8; n=100; vo=60; tm=(2*vo*sind(j))/g; for t=linspace(0,tm,n) x=vo*cosd(j)*cosd(i)*t; y=vo*cosd(j)*sind(i)*t; z=vo*sind(j)*t-(1/2)*g*t^2; plot3(x,y,z,’ro’) pause(tm/n) hold on end xm=((vo.^2).*sind(2.*j).*cosd(i))./g; ym=((vo.^2).*sind(2.*j).*sind(i))./g; zm=((vo.*sind(j)).^2)./(2.*g); xlim([min(min(xm),0) max(xm)]); ylim([min(min(ym),0) max(ym)]); zlim([min(min(zm),0) max(zm)]); grid on Líneas de comando: >> dibujotipoparabolico(10,20)

5

4

3

2

1

0 40 20 0

0

20

40

60

80

100

7. Usando el Comando meshgrid, crea una retícula cuadrada en el intervalo x=[-1 1] y=[-2 2] emplea para ello un paso de malla de valor 0.1. Crea una matriz de ceros del tamaño de las matrices que definen la retícula. Representa, empleando el comando mesh, la matriz de ceros creada sobre la retícula. Representa gráficamente la superficie

. Dibuja las curvas de nivel usando el comando contour

42

Línea de comandos: >> y=[-2:0.1:2];x=[-1:0.1:1]; >> B=meshgrid(x,y); >> A=zeros(size(B)); >> mesh(A)

1 0.5 0 -0.5 -1 60 40 20 0

0

15

10

5

20

25

>> [A,B]=meshgrid(x,y); >> Z=exp(-((A.^2)+(B.^2))); >> mesh(A,B,Z)

1 0.8 0.6 0.4 0.2 0 2 1 0 -1 -2

>>contour(A,B,Z)

-0.5

-1

0.5

0

1

2 1.5 1 0.5 0 -0.5 -1 -1.5 -2 -1

-0.8

-0.6

-0.4

-0.2

0

0.2

0.4

0.6

0.8

1

CONTOUR+MESH >> y=[-2:0.1:2];x=[-1:0.1:1]; >>B=meshgrid(x,y); >>A=zeros(size(B)); >> [A,B]=meshgrid(x,y);

43

>>Z=exp(-((A.^2)+(B.^2))); >> figure >> mesh(A) >> meshc(Z)

1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2

25 20

0.1 0 45

15 40

35

10 30

25

20

15

5 10

5

0

0

SURF+CONTOUR >> y=[-2:0.1:2];x=[-1:0.1:1]; >> B=meshgrid(x,y); >> A=zeros(size(B)); >> mesh(A) >> [A,B]=meshgrid(x,y); >> Z=exp(-((A.^2)+(B.^2))); >> surfc(Z)

44

1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2

20

0.1

15

0

10 40

30

5 20

10

0

0

SECCIÓN 7. Entrada y salida 1. Usando el comando fprintf escribe en la pantalla el número pi con 10 decimales. Mira previamente las opciones de formato de fprintf. >> fprintf('%11.10f',pi) 3.1415926536

2. Usando fprintf escribe el número e con 2 decimales en una línea y con 8 decimales en la siguiente. >> fprintf('%3.2f\n%9.8f\n',exp(1),exp(1)) 2.72 2.71828183

3. Sea la variable x = 174, usa fprintf para ver en una sola línea en pantalla el valor de x como:  un entero (%d), >> fprintf('%d\n',x) 174



un entero reservando 4 columnas (%4d)



un entero, reservando 4 columnas y rellenando con 0's (%04d)

>> fprintf('%4d\n',x) 174 >> fprintf('%04d\n',x) 0174



un entero con signo (%+d)



un número real con 2 decimales (%.2f)

>> fprintf('%+d\n',x) +174 >> fprintf('%.2f\n',x) 174.00



un número real en notación científica (%e ) >> fprintf('%e\n',x)

45

1.740000e+002

4. Usa el comando input para asignar el valor 222 a una variable. Pregunta a continuación por el valor de dicha variable en la línea de comando de matlab. %input le asocia a la variable 'a' aquel número que le das cuando te %devuelve lo que le has puesto en (' ') >> a=input('dame un número:'); dame un número:222 >> a a= 222

5. Crear un programa que dado un número introducido por el usuario, devuelva por pantalla la tabla de multiplicar del número introducido. Edit: %tablademultiplicar a=input('¿Cuál es el número cuya tabla de multiplicar desea saber?'); for x=1:1:10 fprintf('%3i %1s %1i %1s %1i\n' ,a,'x',x,'=',(a*x)) end

Línea de comandos: >> tablademultiplicar ¿Cuál es el número cuya tabla de multiplicar desea saber?4 4 4 4 4 4 4 4 4 4 4

x x x x x x x x x x

1=4 2=8 3 = 12 4 = 16 5 = 20 6 = 24 7 = 28 8 = 32 9 = 36 10 = 40

6. Crear 5 vectores aleatorios denominados v1, v2, v3, v4 y v5 de tamaños diferentes y guardar su contenido en un fichero usando el comando save. Recuperar la información usando load verificando que la información guardada y leída es la misma. (a) Los vectores v deben tener la misma dimensión, ejecuta >> save ale1.txt v* – ascii, donde el nombre del fichero es ale1.txt en este fichero se almacenará la información contenida en todas las variables que empiezan por v y su formato será legible por eso se usa la opción -ascii. Usa type para visualizar el contenido de este fichero. Recupera la información usando load, es decir ejecuta >> load ale1.txt y asigna, por ejemplo, a una variable vectorial a el contenido de la segunda fila de las columnas de las 2 a la 4 ejecutando >>a=ale1(2,[2:4]). Borra el fichero usando >>delete ale1.txt >> v1=rand(1,3);v2=rand(1,3);v3=rand(1,3);v4=rand(1,3);v5=rand(1,3); >> save aleatorios1.txt v* -ascii >> type aleatorios1.txt 1.4188634e-001 7.9220733e-001 3.5711679e-002 6.7873515e-001 3.9222702e-001

4.2176128e-001 9.5949243e-001 8.4912931e-001 7.5774013e-001 6.5547789e-001

9.1573553e-001 6.5574070e-001 9.3399325e-001 7.4313247e-001 1.7118669e-001

>> a=aleatorios1(1,:)

46

??? Undefined variable aleatorios1. >> who Your variables are: v1 v2 v3 v4 v5 >> load aleatorios1.txt >> who Your variables are: aleatorios1 v2 v1 v3

v5

v4

>> a=aleatorios1(1,:) a= 0.1419 0.4218 0.9157 >> delete aleatorios1.txt

(b) Los vectores v deben tener diferente longitud. Al usar load verás que no es posible utilizarlo. Ejecuta las siguientes órdenes >>save ale2 v* %no lo graba en formato legible ya que no usas la opción –ascii, el fichero sin extensión se llamará ale2 y v* indica, como sabes,que se grabaran en este fichero las variables que comienzan con v, >>clear v* % para borrar todas las variables que empiezan por v del workspace, >> who %para ver si existen variables que comienzan con v, >> dir ale* % para ver los ficheros que comienzan con ale, ¿qué ves?, >>load ale %recupera información junto con asignación, >>v1 %verás que recupera la información asignada a v1. >> v1=rand(1,2);v2=rand(1,3);v3=rand(1,4);v4=rand(1,5);v5=rand(1,6); >> save aleatorios2.txt v* -ascii >> load aleatorios2.txt ??? Error using ==> load Number of columns on line 1 of ASCII file C:\Documents and Settings\CLS\Mis documentos\MATLAB\aleatorios2.txt must be the same as previous lines. %Nos da error debido a que los vectores v* no tienen las mismas dimensiones y no es posible la matriz. >> >> >> >>

delete aleatorios2.txt save aleatorios2 v* clear c* who

Your variables are: a

aleatorios1

>> dir aleatorios* aleatorios2.mat aleatorios1.txt >> load aleatorios2 >> who Your variables are: a v1 aleatorios1 v2

v3

v4

v5

47

%Si no lo guardamos como “.txt” y sin “-ASCII” se guarda como .mat., y es un fichero sin extensión, esto nos permite cargar el documento con “load” %sin tener que hacer previamente “import data”.

(c) Volvemos a generar cinco vectores de longitud diferente y grabamos su información legible en un fichero ale3.txt con save. Si usas load no podrás recuperar la información. Vete a la pestaña de File y activa la opción Import Data hasta terminar. Mira si existe alguna variable con el nombre ale3 en el workspace y visualiza su contenido, ¿qué ves?, ¿podrías utilizar la información que contiene?. >> v1=rand(1,4),v2=rand(1,2),v3=rand(1,3),v4=rand(1,5),v5=rand(1,3) v1 = 0.9575

0.9649

v2 = 0.9572

0.485

v3 = 0.8003

0.1419

0.4218

v4 = 0.9157

0.7922

0.9595

v5 = 0.8491

0.9340

0.6787

0.1576

0.9706

0.6557

0.0357

>> save aleatorios.txt v* -ascii >> type aleatorios.txt 9.5750684e-001 9.5716695e-001 8.0028047e-001 9.1573553e-001 8.4912931e-001

9.6488854e-001 4.8537565e-001 1.4188634e-001 7.9220733e-001 9.3399325e-001

1.5761308e-001 9.7059278e-001 4.2176128e-001 9.5949243e-001 6.5574070e-001 3.5711679e-002 6.7873515e-001

>> a=aleatorios(1,:) ??? Undefined variable aleatorios. >> load aleatorios.txt ??? Error using ==> load Number of columns on line 1 of ASCII file C:\Documents and Settings\CLS\Mis documentos\MATLAB\aleatorios.txt must be the same as previous lines. >> who Your variables are: v1 v2

v3 v4

v5

%Nos da error debido a que los vectores v* no tienen las mismas dimensiones y no es posible la matriz. %Cuando esto ocurre en Matlab le damos a ‘file/import data’ seleccionamos el documento, en nuestro %caso “aleatorios.txt”, y nos aparece una ventana donde se encuentran nuestros vectores todos con la %misma dimensión, en los huecos aparece “NAN”, y le damos a “next” y “finish”. %Después de realizar esto ya podemos llevar a cabo, dado que se convierte en una variable: >> a=aleatorios(1,:) a= 0.9575 0.9649 >> aleatorios aleatorios =

0.1576

0.9706

48

0.9575 0.9572 0.8003 0.9157 0.0357 0.8491

0.9649 0.4854 0.1419 0.7922 NaN 0.9340

0.1576 NaN 0.4218 0.9595 NaN 0.6787

0.9706 NaN NaN 0.6557 NaN NaN

>> who Your variables are: a v1 aleatorios v2

v3

v4

v5

7. Ayuda en el uso de fopen, fprintf y fclose en ficheros de texto. fid=fopen(‘nombrefichero’,’modo’) fid es el identificador del fichero que se utilizará en el resto de las operaciones. Es un número natural que asigna matlab al nombre del fichero que le hayamos dado en nombrefichero. Es decir, fid es un nombre de variable a la que se le va a asignar el anterior valor numérico, de tal forma que internamente matlab sabe que el número asignado a fid se refiere al nombrefichero. modo. son las opciones identificativas del fichero r – se refiere a un fichero existente del que se va a leer la información que contiene. w– fichero que queremos crear. Si nos equivocamos y existiera se borraría toda su información. rt para leer el fichero en modo texto, es decir, legible wt para escribir un fichero en modo textofclose(fid) Cierra el fichero identificado con fid Ejemplo: fprintf(fid, 'El valor es: Escribe en el fichero identificado con fid %8.2f\n', y) una cadena de caracteres, el valor de la variable y según el formato indicado y finalmente salta de línea. La orden fprintf sólo se usara para crear ficheros con modo w. Ejemplo: fscanf sirve para ver la información de [A,contador]=fscanf(fid,’formato’) tipo texto, es decir legible, que contiene el fichero identificado con fid. Esta información se almacena en un vector A y en contador se almacena el número total de datos que contiene el fichero. Ejemplo: La opción [14 inf] estructura la variable A [A,contador]=fscanf(fid,’formato’,[14 como una matriz de 14 filas hasta el final inf]) de los datos

49

8. Realizar la misma operación del ejercicio 6 creando un fichero de texto plano (pone extensión .dat) haciendo uso de la ayuda del punto 7 y visualízalo con un editor de textos. >> v1=rand(1,3);v2=rand(1,3);v3=rand(1,3);v4=rand(1,3);v5=rand(1,3); >> alea3=fopen('aleatorios3.dat','w') alea3 = 3 >> fprintf(alea3,'%8.7f %8.7f %8.7f %8.7f %8.7f \n',v1,v2,v3,v4,v5) ans = 153 >> fclose(alea3) ans = 0

9. Leer el fichero de precipitaciones anual precip.txt que contiene los campos (Año, Enero, Febrero, .... Noviembre, Diciembre, Total) para ello usa el comando type. ¿observas algún error en los datos?. Indicalos. A continuación realiza las siguiente operaciones: (a) Extraer las precipitaciones en cada uno de los meses con la ayuda indicada en el punto 6. Si se usara la función interna load se haría de la siguiente forma: >> load(‘precipt.txt’). Si además se quisiera extraer la información de la precipitación total anual en una variable llamada agno y las precipitaciones en cada uno de los meses en una variable llamada mes se haría: >>agno=precip(:,1); mes=precip(:,2:14). Hazlo ahora usando los comandos de la ayuda del punto 7 de esta sección >> load precip.txt >> agno=precip(:,1); >> mes=precip(:,2:14);

(b) Dibujar en una gráfica de la evolución de la precipitación anual del mes de enero en los diferentes años. >> plot(agno,mes(:,1))

50

200

180

160

p re c ip it a c ió n e n e ro

140

120

100

80

60

40

20

0 1980

1985

1990

1995 años

2000

2005

2010

(c) Dibujar en una única gráfica (usa el comando figure antes de hacerlo) la evolución de la precipitación anual de todos los meses y del total de los diferente años. >> >> >> >>

figure hold on plot(agno,mes(:,14),'k') plot(agno,mes) Precipitaciones 1983-2010

8000

7000

6000

precipitación

5000

4000

3000

2000

1000

0 1980

1985

1990

1995 años

2000

2005

2010

%La función de color negro es la de las precipitaciones totales a lo largo de los diferentes años, en la gráfica también podemos observar esos errores que anteriormente hemos visto a simple vista en el documento.

(d) Chequear si la columna de precipitación total corresponde a la suma total de las precipitaciones de los meses. Si no fuera así construye otro fichero en el que se escriba el dato de precipitación correcto para un mes concreto cuando la discrepancia entre la suma mensual difiere mucho de la precipitación total, para ello considera que la precipitación total es correcta, y si la precipitación total no difiere mucho de la suma de las precipitaciones mensuales para un año sustituye el valor de esta suma en el dato de precipitación total. Hazlo en un script precip_bien.m. En este script debes incorporar la visualización por pantalla del año donde se ha detectado un error, de la suma total mensual y de la precipitación total. Si la suma mensual difiere mucho de la precipitación total que índique además el mes y la precipitación de ese mes que es

51

errónea. Llama a este fichero de datos con los errores corregidos precip_bien.dat. EDIT %precip_bien1 [nf,nc]=size(mes); for i=1:1:nf A=mes(i,1:nc-1); z=sum(mes(i,1:nc-1)); total=mes(i,nc); for j=1:1:nc-1 if z~=total & A(j)>total %consideramos que total está bien y el error se encuentra %en A(j) disp(['ERROR EN UN MES']) disp(['año:']) disp([agno(i)]) disp(['mes:']) disp([j]) disp(['precipitación del mes erronea:']) disp([A(j)]) disp(['precipitación total anual erronea:']) disp([z]) disp(['precipitación total anual correcta:']) disp([total]) t=z-total; mes(i,j)=A(j)-t; z=total; end end if z~=total %el error se encuentra en el total disp(['ERROR EN LA SUMA']) disp(['año:']) disp([agno(i)]) disp(['precipitación total anual correcta:']) disp([z]) disp(['precipitación total anual erronea:']) disp([total]) mes(i,nc)=z; end end precip_bien=[agno,mes]; %save precip_bien.dat precip_bien –ascii Línea de comandos: >> load precip.txt >> mes=precip(:,2:14); >> agno=precip(:,1); >> precip_bien1 ERROR EN UN MES año: 1985 mes: 10 precipitación del mes erronea: 7205 precipitación total anual erronea:

52

7643 precipitación total anual correcta: 510.5000 ERROR EN LA SUMA año: 1994 precipitación total anual correcta: 223.5000 precipitación total anual erronea: 220.5000 ERROR EN LA SUMA año: 1997 precipitación total anual correcta: 318.5000 precipitación total anual erronea: 336.5000 ERROR EN UN MES año: 1998 mes: 10 precipitación del mes erronea: 765.5000 precipitación total anual erronea: 1.0485e+003 precipitación total anual correcta: 358.5000 ERROR EN LA SUMA año: 2002 precipitación total anual correcta: 724 precipitación total anual erronea: 845.5000 ERROR EN LA SUMA año: 2007 precipitación total anual correcta: 576 precipitación total anual erronea: 617

53

>> save precip_bien.dat precip_bien –ascii >> precip_bien1

(e) Usando las órdenes de fopen, etcétera indicadas en la ayuda del punto7, crea otro fichero con el nombre precip_modificado.dat en el que se añada la fila media que corresponde a la media aritmética de la precipitación mensual y como última columna la media aritmética mensual de todos esos años. Crea un script para hacerlo e incorpora en este script la representación gráfica, en dos figuras, de los valores medios indicados. EDIT: %media_precip [nf,nc]=size(mes); for i=1:1:nf media=mean(mes(i,1:nc-1)); m(i,1)=media'; end for j=1:1:nc-1 media1=mean(mes(1:nf,j)); m2(1,j)=media1; end precip_modificado=([agno,mes,m]); subplot(2,1,1) plot(agno,m) t=[1:12]; subplot(2,1,2) plot(t,m2)

media aritmética de las precipitaciones por mes

media de las precipitaciones por año

>> media_precip 80 70 60 50 40 30 20 10 1980

1985

1990

1995 años 1983-2010

2000

2005

2010

2

4

6 meses Enero(1)-Diciembre(12)

8

10

12

70 60 50 40 30 20 10 0

0

>> p1=fopen('precip_modificado.dat','w'); >> fprintf(p1,'%6.5f %6.5f \n',precip_modificado,m2); >> fclose(p1);

54

SECCIÓN 8. Formato IEEE754 y errores 1. Un algoritmo es inestable cuando los errores de redondeo se acumulan degradando la exactitud del resultado final, a pesar de que en aritmética exacta el algoritmo sea correcto. La sucesión siguiente es un ejemplo de algoritmo numérico inestable,

Se puede comprobar que el término general de esta sucesión es igual a:

Si suponemos que x'n son los valores exactos de la sucesión se puede estimar el error de redondeo de xn. Escribir un programa en MATLAB (.m) que resuelva esta sucesión y visualice en pantalla, desde n igual a 1 hasta 20 iteraciones, xn, x'n, el error absoluto y el error relativo

. Escribe el programa y los resultados que

obtienes. %En la teoría en la última diapositiva de aritmética aparece la solución a este problema pero en vez de 20 términos aparecen 23, al hacerlo nosotros nos damos cuenta de que nos aparecen los 19 primeros términos iguales a los de la diapositiva y el 23 de la diapositiva como nuestro término 20 pero el 20, 21 y 22 de la diapositiva (en negro) en nuestro programa no aparecen, esto es debido a que en la diapositiva se repiten en esos términos los términos 17,18 y 19 (en verde).

Edit:

55

%algoritmo clear %Pongo el clear ya que si le pides que te haga un vector de 22 y luego uno %de 20 con este mismo scrip devuelve al pedirle el de 20 el de 22, sin %embargo al poner el clear limpia lo anterior y te devuelve el vector de %los elementos que le pides x(1)=1; x(2)=1/3; n=1; N=20; while n<=N if n<(N-1) x(n+2)=((13/3)*x(n+1))-((4/3)*x(n)); end x1(n)=(1/3)^(n-1); errorabsoluto(n)=abs(x(n)-x1(n)); errorrelativo(n)=(errorabsoluto(n))/(abs(x1(n))); n=n+1; end disp('x=');disp(x) disp('x1=');disp(x1) disp('error absoluto:');disp(errorabsoluto) disp('error relativo:');disp(errorrelativo) Línea de comandos >> algoritmo x= Columns 1 through 5 1.000000000000000 0.012345679012343

0.333333333333333

0.111111111111111

0.037037037037036

0.001371742112432

0.000457247370625

0.000152415789465

0.000005644977344

0.000001881468722

0.000000626394672

-0.000000029940803

-0.000000204941979

-0.000000848160840

0.333333333333333

0.111111111111111

0.037037037037037

0.001371742112483

0.000457247370828

0.000152415790276

0.000005645029269

0.000001881676423

0.000000627225474

Columns 6 through 10 0.004115226337436 0.000050805260180 Columns 11 through 15 0.000016935074827 0.000000205751947 Columns 16 through 20 0.000000056398875 -0.000003402107668 x1= Columns 1 through 5 1.000000000000000 0.012345679012346 Columns 6 through 10 0.004115226337449 0.000050805263425 Columns 11 through 15 0.000016935087808 0.000000209075158

56

Columns 16 through 20 0.000000069691719 0.000000000860392

0.000000023230573

0.000000007743524

0.000000002581175

error absoluto: 1.0e-005 * Columns 1 through 5 0

0 0.000000000016653 0.000000000077022 0.000000000316240

Columns 6 through 10 0.000000001267389 0.000000324532319

0.000000005070727

0.000000020283233

0.000000081133066

0.000005192517181

0.000020770068726

0.000083080274904

0.005317137593855

0.021268550375422

0.085074201501686

Columns 11 through 15 0.000001298129294 0.000332321099616 Columns 16 through 20 0.001329284398464 0.340296806006744 error relativo: 1.0e+003 * Columns 1 through 5 0

0 0.000000000000000 0.000000000000000 0.000000000000000

Columns 6 through 10 0.000000000000003 0.000000000063878

0.000000000000037

0.000000000000444

0.000000000005323

0.000000009198388

0.000000110380661

0.000001324567931

0.002288853385213

0.027466240622556

0.329594887470678

Columns 11 through 15 0.000000000766532 0.000015894815175 Columns 16 through 20 0.000190737782101 3.955138649648134

2. Matlab/Octave utilizan el formato representación de números reales:

a.

IEEE754 en

doble

precisión

para

la

Teniendo en cuenta que la mantisa es de 52 bits, corroborar que la precisión del

57

computador es eps=2-52 >> eps ans = 2.2204e-016 >> 2^(-52) ans = b.

2.2204e-016 Realiza la siguiente operación a=1+2-52, b=1-a. ¿Qué valor se obtiene?¿Por qué? >> a=1+(2^(-52)) a= 1.0000 >> b=1-a b=

-2.2204e-016 Este error se debe al error de redondeo, ya que 2^(-55) es un número muy pequeño. El error es del orden del eps del ordenador. La operación la tiene en cuenta pero no representa el formato del resultado. c. Realiza la siguiente operación a=1+2-53, b=1-a. ¿Qué valor se obtiene?¿Por qué? ¿Qué tipo de error se produce y por qué? >> a=1+2^(-53) a= 1 >> b=1-a b=

d.

0 Matlab considera (2^(-53)) como si fuese cero debido a su reducido tamaño, ya que está por debajo de la precisión. Este error se debe al error de redondeo, se produce cuando el número tiene una representación binaria infinita. No se puede representar debido a su largo excesivo. Ejecuta el siguiente script. ¿Son correctos los resultados? Razona tu respuesta >> ejercicio %ejercicio k=1; while ((1.0+2^(-k))>1.0) k=k+1; end fprintf(1,'El número de bits de la mantisa es %d\n', k-1); fprintf(1,'El eps calculado es %e\n', 2^-(k-1));

58

El número de bits de la mantisa es 52 El eps calculado es 2.220446e-016 >> eps ans = 2.220446049250313e-016 e. ¿Qué tipo de errores se producen cuando se teclea las siguiente cantidades 10^308 y 10^309? >> 10^308 ans = 1.0000e+308 >> 10^309 ans = Inf

%Error de desbordamiento, el número es más alto que la capacidad de Matlab, alcanza esta capacidad. El primer número representado es el número más alto que entra dentro de la capacidad de Matlab. Números en binario 3. Crea la función dec2bin_entero que ante la entrada de un número natural devuelva un vector con el contenido de los números en binario. Usa las funciones internas mod y floor. Ejemplo: x=input('¿número decimal?'); k=1; while x/2 > 0 %ya que cuando sea cero ese resto es el primer número en binario, no hay más n(k)=mod(x,2); x=floor(x/2); k=k+1; end n(k)=x; %el vector que nos da es el numero binario pero al contrario por lo que %debemos crear un vector opuesto. p=n((end-1):-1:1)

>> decimal ¿número decimal?19 p= 1 4.

0 0 1 1 Crea la función bin2dec_entero que realice el funcionamiento inverso, dado un

>> dec2bin_entero(19) ans = 1 0 0 1 1 vector en binario devuelva su valor en decimal. 5.

59

>> bin2dec_entero([1 0 0 1 1]) ans = >> dec2bin_decimal(0.125) 19 ans = 0 0 1 6. Crea la función dec2bin_decimal que ante la entrada de un numero fraccionario (sin parte entera, es decir, decimal positivo y menor que cero) y devuelva un vector con el contenido de los números en binario. Ejemplo: 7. Crea la función dec_to_ieee754 que implemente la conversión al formato IEEE754 en simple precisión y normalización 1.m. La función tendrá la siguiente cabecera [ieee, signo, mantisa, ex]=dec_ieee754(x) donde: - ieee es el vector con el número codificado en dicho formato - signo corresponde al bit de signo en dicho formato - mantisa corresponde al vector mantisa en dicho formato - exponente corresponde al campo exponente en dicho formato

>> [i,s,m,e] = dec_to_ieee754(-0.125) i=1 0 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 s= 1 m= 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 e= 0 1 1 1 1 1 0 0 Esta función debe llamar a las funciones que previamente has creado: dec2b_decimal, dec2bin_entero.m 8. Crear la función análoga ieee754_to_dec que dado un número en formato IEEE754 en simple precisión devuelva el número convertido en decimal. Llama a tu finción bin2dec_entero.m

60