2012II+Manual-Uso Del Matlab
2012II+Manual-Uso Del Matlab
2012II+Manual-Uso Del Matlab
Editor:
EDITADO EN HUANCAYO, DICIEMBRE, 2008. Hecho el Depsito Legal en la Biblioteca Nacional del Per N 2008- 15609
INDICE
Introduccin Fundamentos bsicos Arreglos Operaciones con arrays Arreglos bi dimensionales Construccin de arrays Ejercicios propuestos de arrays Matrices Definicin de matrices desde el teclado Conceptos bsicos sobre matrices Operaciones con matrices Suma y Resta de matrices Multiplicacin de matrices Productos de un nmero real por una matriz Matriz inversa Manipulacin de matrices comandos relacionados con ellos Funciones que actan sobre matrices Ejercicios propuestos sobre matrices Grficos con Matlab Grficos simples Tipos y colores de lneas Ejemplos relacionados con Grficos Ley de Boyle Ejemplo Calculo del volumen con la Ecuacin de van Der Waals Ejercicios propuestos con grficos Grficos 3D Plot3 Grficos de malla: mesh(X,Y,Z) meshz(X,Y,Z) meshc(X,Y,Z) Grficos de contorno Polgonos en tres dimensiones Fundamentos tericos de programacin Procesamiento de datos Algoritmos y programas Partes de un programa Diagramas de flujo Aplicaciones de los diagramas de flujo Elementos bsicos de un programa Ejercicios Comando If..ElseEnd Operadores necesarios algunas variedades del comando IfEnd Cdigo de programa para determinar si un nmero es par o impar Cdigo de programa para redondear un nmero de hasta 4 cifras Uso de if mltiple Uso de comando Switch ...End Ejercicios propuesto sobre comando if Comando For...End - Estructura Ejercicios Cdigo de programa produccin de amoniaco: Conversin vs Temperatura a diferentes presiones Comando WhileEnd - estructura y relaciones
7 9 19 21 24 25 20 30 31 33 34 35 37 42 43 47 47 50 51 52 60 66 75 77 80 81 83 84 85 86 88 89 82 85 101 105 107 108 110 112 115 117 118 121 126 133
Cdigo de programa para calcular la presin de rocio Ejercicios propuestos Ecuaciones lineales Sistema de ecuaciones lineales parte terica Aplicaciones con programas mtodo de biseccin Mtodo de la Regla Falsa Mtodo de Newton Raspn Solucin grfica Ejercicios propuestos Ajuste de curvas Mtodo mnimos cuadrados digo de programa grficos C Ajuste de curvas Diferenciacin numrica Polinomios Funciones orientadas al clculo con polinomios Ajuste polinomial Diferenciacin Uso del Comando diff Estadstica Descriptiva : media, mediana Desviacin estndar ejercicios Integracin numrica Mtodo de Simpson y mtodo abreviado
136 140 145 146 147 148 149 150 151 152 154 157 158 159 163 164 165 167 169
Bibliografa
173
Prlogo
El hombre que ha dejado de aprender, no merece deambular libremente en estos das tan peligrosos M.M. Coady
El presente aporte es como parte del proceso de aprendizaje en las aulas universitarias, dirigido a los estudiantes de ingeniera de Qumica y otras facultades y/o escuelas profesionales.
La presentacin del texto se encuentra plasmada en la parte introductoria del presente texto.
Dedicatoria
A la memoria de Francisco Or Alcntara JJOC
Introduccin
El presente texto constituye un aporte para los estudiantes de la Facultad de Ingeniera Qumica, la misma que contiene problemas de los exmenes del curso de LENGUAJES DE PROGRAMACION de los aos 2006 y 2007, as mismo se encuentran las prcticas de laboratorio correspondientes a las ocho primeras semanas establecidas en el sylabus correspondiente.
El texto preparado con fines didcticos, escrito usando cdigo del programa MATLAB pueden ser trasladados a otros programas con pequeas modificaciones, que se deja al estudiante para estimularlos en su desarrollo y que exploren otros medios de solucin para un problema.
Se considera que el presente texto ser tomado como una versin preliminar y su perfeccionamiento se har en base a las recomendaciones de los usuarios que sern recibidos en el siguiente correo electrnico: erojasza@hotmail.com, quedando sumamente agradecido por sus aportes.
LOS AUTORES
Agradecimiento
Nuestros sinceros agradecimientos a los docentes de la facultad de Ingeniera de Sistemas por las constantes recomendaciones para la mejora del presente trabajo. As mismo, al encargado del Centro de Cmputo de la mencionada Facultad por su gentileza en el apoyo con la tecnologa de informacin correspondiente.
Los Autores
OBJETIVO
Usar el MatLab en modo interactivo con las funciones establecidas en el Lenguaje de programacin citado.
10
y observaremos la siguiente pantalla, donde podemos escribir algunos cdigo para la ejecucin en modo interactivo. Se denomina modo interactivo, debido a que estamos interactuando con el computador haciendo uso del lenguaje de programacin MATLAB. En la pantalla adjunta se puede observar el promt o seala de ingreso de datos, funciones o palabras reservadas por el MATLAB.
11
Tambin; podemos carga el HELP o ayuda del matlab y ubicarnos las principales funciones que pose el MATLAB para hacer las operaciones que desea el usuario, estas mismas funciones se usan para escribir cdigos de programa. O presionando la tecla F1 y se observa la figura siguiente figura.
12
Entre ellas los operadores aritmticos, operadores relacionales y operadores lgicos, tambin el significado de algunos caracteres especiales usados por el MATLAB.
Una de las primeras funciones que se muestra es abs que sirve para tratar cualquier magnitud o arreglo. Tambin se muestran las funciones acos y acosd con las cuales se calculan la inversa de la funcin coseno y los resultados se muestran en radianes y grados respectivamente.
Cada funcin tiene su parte explicativa y uso para lo cual haremos clik sobre una funcin determinada, por ejemplo; para ver la sintaxis de la funcin sqrt, el cual nos permite sacar la raz cuadrada de cualquier nmero, hacemos click sobre la funcin y observamos la sintaxis y tiene la siguiente forma:
13
Podemos usar dicha funcin en la VENTANA DE COMANDO del MATLAB y tendremos el siguiente resultado:
Ahora; usaremos la funcin trigonomtrica seno (sin en matlab) para ver los resultados (primero se observar la sintaxis de la funcin, tambin se observa el uso de la funcin PLOT que nos permite graficar dicha funcin en un determinado rango)
14
Tambin podemos hacer uso de operadores aritmticos para realizar operaciones ms complejas
OPERADORES + * / ^ ()
Significado Suma o adicin Resta o sustraccin Producto o Multiplicacin Divisin o cociente Potenciacin Orden de operacin
15
OPERADORES < <= > >= == Menor que Menor que o igual a Mayor que Mayor que o igual a Igual a No igual a
~=
Tambin podemos hacer uso de operadores lgicos para tomar decisiones; tales como:
A&B AB ~A
El smbolo &, y ~ son operadores lgicos AND, OR y NOT, ellos trabajan elemento por elemento en los arreglos. Los operadores lgicos retornan un arreglo lgico con elementos (1) para verdadero y elementos (0) para falso.
POR EJEMPLO: Dado A=3, B=2 y C=5; vamos a calcular las siguientes expresiones: y A^B+C
A + B + C, A * B C
SOLUCION: Debemos tener en cuenta la asignacin de valores a las variables y proceder al clculo de las expresiones aritmticas OJO: Matlab es sensible a la declaracin de variables y asignacin de valores. Si asignamos el valor de 3 para la variable A, 2 para la variable B y 5 para la variable C, son esas variables las que tiene en cuenta Matlab. Vamos a usar la variable indistintamente la variable a en lugar de la variable A, matlab muestra un mensaje de error Podemos observar los resultados obtenidos:
16
2. Resolver los siguientes ejercicios de manera interactiva con los comandos del Matlab
F( x, y) =
G ( x, y ) =
Con x = 1, 3, 5, 7
y = 7, 5, 3, 1
17
Solucin
x=1:2:7; y=7:-2:1;
f=(sin(x).*cos(x))./(log10(x).*log(x))+(sqrt(exp(x)+exp(y))).*(log10(y).*log(y))./(sin(y).*cos(y)) g=sqrt((x+y)./(x.*y)+log10(10^5^exp(1)))+(sqrt(sin(x).*sin(y))./sqrt(cos(y).*cos(x)))*exp(exp(1))
(5 + 7
* 3 + 4 * 7/3 )6 + 4567
6/8
333 3 158
(95 * 7 3 + 4 * 7 / 3 ) 6
45678 + log 12
9 / 23
3456781 158
12
log ( 1234 ) 6
+ 3457
6 /8
3345453 158
(4565
26/8
32345 3 111158
18
61/2
3. Defina las variable x=3.6, y=14.2 y z=456 y calcule las siguientes expresiones
( 3 * x + 7 * y + 3 * z + 4 * x * y * z * 7 * y/3 )6 * x * y
(5 * x * y * z + 7 * 3 * Log(x * y * z) + 4 * 7/3)6
4. Sabiendo que p=543 y r=234 y definida las variables: a, b, c y d como: a=3*p, b=4*a-b; c=3*b-c y d=3,5*(a*c*r/b), calcular las siguientes expresiones:
a * b ( a + d) 2 a + c a *b
6 + 4567 3 * a
33 * b3 158 * z
a * b ( a + d) 2 Log(a) * Ln(a * b) + c a *b
6 + ( a b c) 3 * a
c * z * b3 a * b * 158 * z
5. Acudimos a su imaginacin para que pueda realizar 30 ejercicios del nivel de su grado de estudios, la misma que ser evaluada rigurosamente. (Los ejercicios deben incluir operaciones con funciones algebraicas, trigonomtricas: seno, coseno, tangente, etc. Tambin deben contener operadores relacionales, lgicos)
19
ARREGLOS
OBJETIVO
Aplicar el concepto de arreglo y sus funciones establecidas en el MatLab para manipular variables y aplicarlos en la solucin de problemas especficos
20
ARREGLOS:
Variables de arreglo unidimensional: Estas variables tienen forma de fila o columna y estn relacionadas con los vectores y las matrices. En MATLAB:
arreglo de fila es lo mismo que vector de fila y arreglo de columna es lo mismo que vector de columna.
La variable x puede definirse como vector de fila especificando sus elementos; por ejemplo: x = [0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9] ;
La variable x tiene 10 elementos, si desea imprimir o visualizar un elemento cualquiera, digite: x con el subndice correspondiente. Por ejemplo, si digita x(8) como un comando se mostrar el siguiente resultado: Ans = 0.7 Una forma equivalente de definir la misma x es: for i = 1 : 10 x(i) = (i - 1) * 0.1 end ; Otra forma de escribir un arreglo de fila, asignado a una variable, con un incremento o decremento fijo es: x = 6 : -0.4 : 0 Que produce el siguiente resultado: x = 6.0 5.6 5.2 4.8 4.4 4.0 3.6 3.2 2.8 2.4 2.0 1.6 1.2 0.8 0.4 0.0 o si digitamos la siguiente expresin se tiene lo siguiente: >>
[x]'
ans = 6.0000 5.6000 5.2000 4.8000 4.4000 4.0000 3.6000 3.2000 2.8000 2.4000 2.0000 1.6000 1.2000 0.8000 0.4000 0.0000
21
La definicin de un arreglo de columna es similar a la de un arreglo de fila, excepto que los elementos se separan mediante signos de punto y coma; por ejemplo, z = [0.0; 0.1; 0.2; 0.3; 0.4; 05; 06 ] Cuyo resultados es: z= 0 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000
Una forma alternativa para definir lo anterior es agregar un apstrofo a un arreglo de fila:
z=[0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 06 ]' ; Cuyo resultado es: z= 0 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000
c(9) = 15
Se supondr c(i) = 0 para i = 1 hasta 7. Por tanto, si teclea c como comando obtendr c= 0 0 0 0 0 0 0 0 15
OPERACIONES CON ARRAYS Los operadores de arrays en MATLAB son los siguientes:
+ adicin o suma sustraccin o resta \ divisin-izquierda / divisin-derecha .* producto elemento a elemento ./ y .\ divisin elemento a elemento .^ elevar a una potencia elemento a elemento
Cuando y y x tienen la misma longitud y la misma forma (fila y columna), los vectores y y x se pueden sumar, restar, multiplicar y dividir empleando los operadores aritmticos de arreglos
z=x+y z=x-y z = x .* y z = x ./ y
Que equivalen a:
Listado 1.12 For i=1 : 6 ; z(i) For i=1 : 6 ; z(i) For i=1 : 6 ; z(i) For i=1 : 6 ; z(i)
+ * /
22
Se puede comprobar digitando la siguiente sentencia o haciendo correr el siguiente programa: x=[1 2 3 4 5 6] y=[2 4 6 8 10 12] z = x ./ y for i=1 : 6 ; z(i) = x(i) ./ y(i) end
El operador potenciacin en arreglos se puede ilustrar con: G = z. ^ 1.2; Donde z es un vector de longitud 6, se coloca un punto antes del operador (^) y g se convierte en un vector de la misma longitud. El enunciado anterior equivale a % prueba de una funcin exponencial clear all; clear memory; clear command history; clc for i=1:6 g(i) = i .^ 1.2 end
Produce el siguiente resultado: g= 1 g = 1.0000 g = 1.0000 g = 1.0000 g = 1.0000 g = 1.0000 2.2974 2.2974 2.2974 2.2974 2.2974 3.7372 3.7372 3.7372 3.7372 5.2780 5.2780 5.2780 6.8986 6.8986 8.5858
g = 1.0000
2.2974
3.7372
23
As mismo:
El tamao de un arreglo puede incrementarse anexndole un elemento o un vector (o vectores). Por ejemplo, suponga: >> x = [ 4 , 8 ] x= 4
El tamao de un arreglo puede incrementarse anexndole un elemento o un vector (o vectores). Por ejemplo, suponga: >> m = [ 24 , 81 ] m= 24 81
El comando que sigue anexa 15 a m y hace que su longitud sea 3 >> m = [ m , 7 ] lo que devuelve m= 24 81 15
81
15
Podemos anexar un vector o varios vectores a un vector columna. Suponga que y es un vector columna, >> y = [ 2 ; 3; 6 ] y = 2 3 6 entonces; >> y = [ y ; 9] produce y= 2 3 6 9
Adjunto; 9 se aade al final del vector columna. Observe que se utiliza un signo de punto y coma para anexar a un vector de columna. De forma similar, [ -4. y ] produce y= -4 2 3 6 9
Un procedimiento inverso consiste en extraer una parte de un vector. Con la y anterior, >> w = y (2:3) Define a w que equivale al tercer y cuarto elemento de y, a saber: w = 2 3
24
Tambin
Para >> z = [9; 2; 3; 5; 8; 10; 12 ] >> z = [9, 2, 3, 5, 8, 10, 12 ] >> size( z ) devolver ans = 4 1
size.
Por ejemplo,
donde la primera cifra es el nmero de filas y la segunda es el nmero de columnas. Esta respuesta nos dice que y es un arreglo de 3 por 1, es decir, un vector de columna de longitud 3.
un arreglo bidimensional es lo mismo que una matriz en MATLAB, Por ejemplo, un arreglo de 3 por 3 se puede definir mediante
m = [0.1,
0.2,
0.3;
0.4,
0.5,
0.6;
0.7,
0.8,
0.9 ] ;
Observe que los elementos de una fila terminan con un signo de punto y coma. Desde luego, todas las filas deben tener el mismo nmero de elementos; si no es as, la definicin no ser aceptada. Si digitamos m obtenemos m= 0.1000 0.4000 0.7000 0.2000 0.5000 0.8000 0.3000 0.6000 0.9000
Los arreglos bidimensionales se pueden sumar, restar, multiplicar y dividir con los operadores aritmticos de arreglos: c=a+b c=ab c = a .*b c = a ./ b
25
Argumentos de arreglo: La mayor parte de las funciones de MATLAB puede aceptar vectores y matrices como argumentos. Por ejemplo, si x= 1 9 entonces: sin(x) producir: ans = 0.8415 04121 0.9093 0.9894 0.1411 0.6570 2 9 3 7
El cual es una matriz del mismo tamao que x. El clculo realizado equivale a: % Ejercicio 1. clear all; clear memory; clear command history; clc x=[1 2 3; 9 9 7] sin(x) % Tambien for i = 1 : 2 for j = 1 : 3 x(i,j) = sin(x(i,j)) end end
CONSTRUCCION DE ARRAYS % CONSTRUCCION DE ARRAYS X = (0 : 1 : 12)' CC = [X X+X X.^2 X.^3 sqrt(X) ] Rpta: X X= 0 1 2 3 4 5 6 7 8 9 10 11 12 CC = 1.0e+003* 0 0.0010 0.0020 0.0030 0.0040 0.0050 0.0060 0.0070 0.0080 0.0090 0.0100 0.0110 0.0120 0 0.0020 0.0040 0.0060 0.0080 0.0100 0.0120 0.0140 0.0160 0.0180 0.0200 0.0220 0.0240 0 0.0010 0.0040 0.0090 0.0160 0.0250 0.0360 0.0490 0.0640 0.0810 0.1000 0.1210 0.1440 0 0.0010 0.0080 0.0270 0.0640 0.1250 0.2160 0.3430 0.5120 0.7290 1.0000 1.3310 1.7280 0 0.0010 0.0014 0.0017 0.0020 0.0022 0.0024 0.0026 0.0028 0.0030 0.0032 0.0033 0.0035 X+X X^2 X^3 Sqrt(X)
26
X = mean(S)
X = max(S)
X = min(S)
[a,b] = max(S)
[c,b] = min(S)
X = sum(S)
X = cumsum(S)
X = prod(S)
X = cumprod(S)
X = sort(S)
[a,b] = sort(S)
X = sortrows(S) X = rand(1,S)
X = randperm(S)
27
FUNCIONES DE POSICION DE ELEMENTOS, FILAS Y COLUMNAS FUNCION disp display isequal : A(:,j) A(i: ) A(:,:) A(j,k) A(:) DESCRIPCION Muestra un arreglo Muestra un arreglo Verdadero si los arreglos son idnticos Operaciones con arreglos Es la j th columna de A Es la i TH fila de A Es el equivalente a un arreglo de dos dimensiones Es A(j), A(j+1), ..A(k) A(:) muestra los elementos de A, mostrados en una columna simple
Ejercicios de Aplicacin - Usos de la funcin (sort) En arreglos o vectores >> x=[ 5 2 1 4 8 9 3 6 5 7 3] x= 5 2 1 4 8 9 >> x1=sort(x) x1 = 1 2 3 En matrices >> x=magic(3) x= 8 1 6 3 5 7 4 9 2 >> x1=sort(x) x1 = 3 1 4 5 8 9
>> x1=sort(x,'ascend') x1 = 1 2 3 3 4
2 6 7
>> x2=sort(x,'descend') x2 = 9 8 7 6 5
>> [a,b]=sort(x,'ascend') a= 1 2 3 3 4 5 b= 3 2 7 11 4 1
9 >> [a,b]=sort(x,'ascend') a= 3 1 2 b= 4 5 6 8 9 7 2 3 1 1 2 3 3 1 2
10
>> [a,b]=sort(x,'descend') a= 9 8 7 6 5 5 b= 6 5 10 8 9 1
11
>> [a,b]=sort(x,'descend') a= b= 8 9 7 1 3 4 5 6 3 2 3 1 2 2 1
2 1 3
28
Ejercicios de Aplicacin - Usos de las funciones length(x), max(x), min(x), mean(x), sum(x), prod(x), cumsum(x) y cumprod(x) En arreglos o vectores >> x=[3 5 8 5 8 2 6] x= 3 5 8 5 >> x1=length(x) x1 = 7 >> x2=mean(x) x2 = 5.28571428571429 >> x3=max(x) x3 = 8 >> x4=min(x) x4 = 2 >> [a,b]=max(x) a= 8 b= 3 >> [a,b]=min(x) a= 2 b= 6 >> x=magic(3) x= 8 1 6 3 5 7 4 9 2 >> c1=mean(x) c1 = 5 5 5 >> c2=max(x) c2 = 8 9 7 >> c3=min(x) c3 = 3 1 2 >> [a,b]=max(x) a= 8 9 7 >> [a,b]=min(x) a= 3 1 2 >> c4=cumprod(c) c4 = 1 3 4 6 >> c=magic(2) c= 1 3 4 2 >> c1=sum(c) c1 = 5 5 En matrices >> v = [4 6 5 8 9 3 4] v= 4 6 5 8 >> c11=sum(v) c11 = 39 >> c12=prod(v) c12 = 103680 >> c13=cumsum(v) c13 = 4 10 15
23
32
35
39
>> c2=prod(c) c2 = 4 6
>> c3=cumsum(c) c3 = 1 3 5 5
b= 1
b= 2
29
2.
3.
4.
5.
7.
8.
9.
10. Teniendo en consideracin que los arreglos sirven para guardar variables de memoria, defina dos arreglos llamados presin (P) y volumen (V), pero el volumen puede generarse como resultado de la siguiente ecuacin PV = nRT (Considere: n=1 mol, R=0.09205 y T=300) 11. En general hacer 20 ejercicio usando arrays y operadores aritmticos, relacionales y lgicos pero relacionado con la carrera de ingeniera qumica
30
MATRICES
OBJETIVO
Conocer las funciones establecidas para el tratamiento de matrices y su relacin con el concepto de arreglos
31
MATRICES: MATLAB es un programa especializado en clculo matricial. DEFINICIN DE MATRICES DESDE TECLADO Las matrices y vectores son variables que tienen nombres. Una matriz de orden mn es un conjunto rectangular de elementos aij dispuestos en m lneas horizontales (filas) y n lneas verticales (columnas) de la forma:
Las matrices se definen por filas; los elementos de una misma fila estn separados por blancos o comas, mientras que las filas estn separadas por pulsaciones intro o por caracteres punto y coma (;).
A partir de este momento la matriz A est disponible para hacer cualquier tipo de operacin con ella. Por ejemplo, una operacin con A es matriz hallar su traspuesta. En MATLAB el apstrofo (') es el smbolo de trasposicin matricial. Para calcular digitamos siguiente: >> A' ans = 1 4 2 5 3 6 7 8 9 A' lo
>> A = [ 1 2 3; 4 5 6; 7 8 9 ]
Como el resultado de la operacin no ha sido asignado a ninguna otra matriz, MATLAB utiliza un nombre de variable por defecto (ans, de answer), que contiene el resultado de la ltima operacin.
Ahora ya estn definidas las matrices A y B, y es posible seguir operando con ellas. Por ejemplo, se puede hacer el producto B * A (deber resultar una matriz simtrica): >> B * A ans = 66 78 90
del la
2 5 8
3 6 9
La variable ans puede ser utilizada como operando en la siguiente expresin que se introduzca.
78 93 108
90 108 126
Tambin podra haberse asignado el resultado a otra matriz llamada B: >> B = A' B= 1 2 3
4 5 6
7 8 9
32
En MATLAB se accede a los elementos de un vector poniendo el ndice entre parntesis ( por ejemplo: x( 3 ) x( i ) ). Los elementos de las matrices se acceden poniendo los dos ndices entre parntesis, separados por una coma (por ejemplo A(1,2) A(i,j)). Las matrices se almacenan por columnas., y teniendo en cuenta esto puede accederse a cualquier elemento de una matriz con un slo subndice. Por ejemplo, si A es una matriz (3x3) se obtiene el mismo valor escribiendo A(1,2) que escribiendo A(4).
Ahora se va a calcular la inversa de A y el resultado se asignar a B. Para ello basta hacer uso de la funcin inv( ): B = inv(A) B = 0.1803 0.2213 -0.1885 0.13110 0.0246 0.0902 -0.0984 0.1066 0.0574
Para comprobar que este resultado es correcto basta multiplicar A por B; >> B * A ans = 1.0000 0.0000 0.0000 0.0000 1.0000 0.0000 0.0000 0.0000 1.0000
De forma anloga a las matrices, es posible definir un vector fila x en la forma siguiente (si los tres nmeros estn separados por blancos o comas, el resultado ser un vector fila): >> x = [ 10 20 30 ]
Por el contrario, si los nmeros estn separados por intros o puntos y coma (;) se vector obtendr un columna: >> y = [ 11; 12; 13 ] % vector columna y=
MATLAB tiene en cuenta la diferencia entre vectores fila y vectores columna. Por ejemplo, si se intenta sumar los vectores x e y se obtendr el siguiente mensaje de error:
Estas dificultades desaparecen si se suma x con el vector traspuesto de y: >> x+y' ans =
>> x+y 11 12 13
% vector fila x= 10
21 32 43
20
30
>> x + z ans = 44 65 52
45
22
Nota: MATLAB considera comentarios todo lo que va desde el carcter tanto por ciento (%).
33
ALGUNOS CONCEPTOS BASICOS: Matriz cuadrada Matriz Nula Matriz Identidad : : : Una matriz de M = N Todo los elementos de la matriz nula son ceros Es una matriz cuadrada donde los elementos diagonales son la unidad y
el resto de los elementos son cero. La orden en matlab para definir una matriz identidad es eye Matriz transpuesta : La transupuesta de la Matriz A = [ a i,j ] es At = [ a
j,i
Matriz de Permutacin :
OPERACIONES CON MATRICES MATLAB operar con matrices por medio de operadores y por medio de funciones. Se han visto ya los operadores suma (+), producto (*), traspuesta (') e invertir inv( ). MATLAB son los siguientes: Los operadores matriciales de
+ * ' ^
\ / .* ./ y .\ .^
divisin-izquierda divisin-derecha producto elemento a elemento divisin elemento a elemento elevar a una potencia elemento a elemento
Estos operadores se aplican tambin a las variables o valores escalares, aunque con algunas diferencias. Nota: no se puede por sumar matrices que no sean del mismo tamao. Si los operadores no se usan de modo correcto se obtiene un mensaje de error. Los operadores anteriores se pueden aplicar tambin de modo mixto, es decir con un operando escalar y otro matricial. En este caso la operacin con el escalar se aplica a cada uno de los elementos de la matriz. Considrese el siguiente ejemplo:
>> A = [ 1 2; 3 4 ] A= 1 2 3 4
>> A * 2 ans = 2 6
4 8
>> A - 4 ans = -3 -1
-2 0
Los operadores de divisin requieren una cierta explicacin adicional. Considrese el siguiente sistema de ecuaciones lineales, Ax = b (1)
Donde x y b son vectores columna, y A una matriz cuadrada invertible. La resolucin de este sistema de ecuaciones se puede escribir en las 2 formas siguientes: x = inv(A) * b x = A\ b (2a) (2b)
34
SUMA Y RESTA DE MARICES: Podemos sumar una matriz a otra o restarla de otra si ambas tienen el mismo tamao (el mismo nmero de columnas y de filas)
EJEMPLOS
>> A = [ 1 2 4; 3 1 2; 4 1 3 ] A= 1 3 4
>> B = [ 7 3 1; 2 3 5; 8 1 6 ] B=
>> A + B ans =
2 1 1
4 2 3
7 2 8
3 3 1
1 5 6
8 5 12 >> x + y ans =
>> x = [ 1; 4; 2 ] x= 1 4 2
>> y = [ 3; 9; 4 ] y= 3 9 4
4 13 6
N-M ans = 11 8 44 44
35
MULTIPLICACION Si el nmero de columnas de A y el nmero de filas de B son idnticos, las matrices pueden multiplicarse como: C = AB, el nmero de filas de C es igual al de A y el nmero de columnas de C es igual al de B, en otras palabras, si A es una matriz de p por q y C es una matriz de q por r, entonces C es una matriz de p por r , mostramos algunos ejemplos tpicos: Suponga que B y C son matrices:
Para recordar el mtodo para multiplicar matrices, existe una disposicin ms llamativa: se sube la matriz de la derecha, y se escribe el producto debajo de ella:
36
EJEMPLOS:
>> M = [ 8 2 9; 3 9 4; 6 9 4 ] M= 8 3 6 2 9 9 9 4 4
>> N = [ 7 3 1 ; 2 3 5 ; 8 1 6 ] N= 7 2 8 3 3 1 1 5 6
>> y = [ 7; 9; 3 ] y= 7 9 3
Ojo: Para poder multiplicar dos matrices, la primera debe tener el mismo nmero de columnas que filas de la segunda
>> A = [ 1 2; 4 3; 0 2 ] A= 1 4 0
>> A = [ 2 1 7 ] A= 2
>> A = [ 8 1 3 ; 1 5 2 ] A= 8 1
>> A = [ 1 2; 4 3; 0 2 ] A= 1 4 0
2 3 2
1 5
3 2
2 3 2
>> B = [ 5 ; 1 ] B= 5 1 >> A * B
>> B = [ 1 2 ; 4 3 ; 0 2 ] B= 1 4 0 2 3 2
>> B = [ 1 2 ; 4 3 ; 0 2 ] B= 1 4 0 2 3 2
>> B = [ 8 1 3 ; 1 5 2 ] B= 8 1 1 5 3 2
ans = 7 23 2
7 18 4
37
Para multiplicar un escalar por una matriz se multiplica el escalar por todos los elementos de la matriz, obtenindose otra matriz del mismo orden.
EJEMPLOS:
>> A = [ 5 ] A= 5
>> M = [ 10 ] M= 10
>> x = [ 1 3; 4 2; 6 2 ] x= 1 4 6 3 2 2
>> B = [ 7 3 1 ; 2 3 5 ; 8 1 6 ] B= 7 2 8 3 3 1 1 5 6
>> N = [ 8 4 7 ; 8 3 8 ; 2 0 7 ] N= 8 8 2
4 3 0
7 8 7
>> A * B = B * A ans = 35 10 40
>> A * x = x * A ans = 5 15 20 10 30 10
15 15 5
5 25 30
40 70 30 80 0 70
38
DIVISION DE MATRICES: La divisin de matrices se define como el producto del numerador multiplicado por la matriz inversa del denominador. Es decir, sean las matrices A y B tal que A/B = AB-1 .. Si una matriz est dividida entre un escalar, todos los trminos de la matriz quedarn divididos por ese escalar. As:
EJEMPLOS
>> A = [ 1 2 4 ; 3 1 2 ; 4 1 3 ] A= 1 3 4
>> B = [ 7 3 1 ; 2 3 5 ; 8 1 6 ] B=
>> x = [ 1; 4; 2 ] x=
2 1 1
4 2 3
7 2 8
3 3 1
1 5 6
1 4 2
2 3 5
8 1 6
Nota: La transpuesta de un producto de matrices es el producto de las trasnpuestas de las matrices en orden invertido. Por ejemplo:
3 1 2
4 1 3
3 1 2
4 1 3
11
14
39
CONCEPTO DE MATRIZ INVERSA Una de las aplicaciones del mtodo de Gauss-Jordan, es el clculo de matrices inversas. Recordamos primero la definicin de matriz inversa. Definicin. Sea A una matriz de que:
nxn
nxn
tal
Se escribe B = A para denotar la matriz inversa. Cuando la matriz inversa existe, es nica, pro no siempre existe la matriz inversa. Un resultado de algebra lineal prueba que la matriz inversa
A 1
Es decir, en una matriz comenzamos por escribir la matriz A, y a su derecha agregamos la matriz identidad
In
de Gauss-Jordan para hacer los ceros y unos y obtener del lado izquierdo la matriz identidad lado derecho lo que obtendremos ser la matriz inversa de A.
I n . Del
Ejemplo 1. Usar el mtodo de Gauss-Jordan para calcular la matriz inversa de la siguiente matriz:
I2 :
a11 = 4
1 4
40
a 22 = 0.25 .
11 0.25
Finalmente, hacemos los 1s en la diagonal principal. Para ello, multiplicamos el rengln 1 por rengln 2 por
1 0.25
1 4
y el
Ejemplo 2.
Solucin.
a11 = 2
debajo de este elemento. Para ello multiplicamos el rengln 1 por tambin, multiplicamos el mismo rengln 1 por
0.3125 2
y lo sumamos al rengln 2;
41
Para elegir el segundo elemento pivote, debemos escoger el elemento mayor (con valor absoluto) entre a 22 = 0.2 y a 32 = 1.25 , el cual obviamente es ste ltimo. Por lo tanto, debemos intercambiar el rengln 2 y el rengln 3. Tenemos entonces:
Procedemos a hacer ceros arriba y abajo de nuestro segundo elemento pivote; para ello, multiplicamos el rengln 2 por rengln 2 por
0 .2 1.25
1.4 25
a33 = 0.4
3.125 0.4
Finalmente, hacemos los 1s en la diagonal principal. Para ello multiplicamos el rengln 1, 2 y 3 por
1 2
1 1.25
1 0.4
42
eje(n), as:
EJEMPLOS
eye(3) ans = 1 0 0 0 1 0 0 0 1
J = [ 1 3 5; 5 6 7 ; 1 1 1 ] J= 1 5 1 3 6 1 5 7 1
Z=inv(J) Z = 1.0e+015 * 0.5004 -1.0008 4.5036 -1.0008 2.0016 -9.0072 0.5004 -1.0008 4.5036
M= eye(4) M= 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1
J = [ 1 2 3 4; 5 6 7 8 ; 1 1 1 1; 2 2 2 2] J= 1 5 1 2 2 6 1 2 3 7 1 2 4 8 1 2
S=M * J S= 1 5 1 2 2 6 1 2 3 7 1 2 4 8 1 2
>> A = [ 3 6 4; 3 8 5; 1 4 2 ] A= 3 3 1
6 8 4
4 5 2
>> B = inv(A) B= 2.0000 -2.0000 1.0000 0.5000 -1.0000 1.5000 -2.0000 3.0000 -3.0000
0 0 1.0000
0 1.0000 0
0 0 1.0000
43
>> x = zeros(5) x= 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0
0 0
0 0
0 0
Tambin:
>> d = ones(3) d= 1 1 1 1 1 1 1 1 1
>> r = ones(3,3) r= 1 1 1 1 1 1 1 1 1
MANIPULACION DE MATRICES
MATLAB proporciona funciones que manipulan matrices (mas a detalle se pueden encontrar en el Help del Matlab):
Devuelve el nmero de filas y de columnas de la matriz A. Si la matriz es cuadrada basta recoger el primer valor de retorno calcula el nmero de elementos de un vector x forma una matriz diagonal A cuyos elementos diagonales son los elementos de un vector ya existente x (puedes usarlo para generar una matriz diagonal con los componentes del vector existente)
A=diag(x)
x=[ 1 4 6 7 ]
x= 1 4 6 7
44
0 4 0 0
0 0 6 0
0 0 0 7
0 0 0 8 0
0 0 0 0 4
Tambin se puede usar el comando diag para buscar las diagonales de una matriz >> A = magic(5) A= 17 24 1 8 15 23 5 7 14 16 4 6 13 20 22 10 12 19 21 3 11 18 25 2 9
45
x=diag(A)
2 5 8
3 6 9
diag(diag(A))
2 11 7 14
3 13 10 8 6 12 15 1
crea una matriz diagonal de sub-matrices a partir de las matrices que se le pasan como argumentos
triu(A)
Forma una matriz triangular superior a partir de una matriz A (no es necesario que sea cuadrada).
tril(A)
rot90(A,k) flipud(A)
Gira k*90 grados la matriz rectangular A en sentido anti-horario. k es un entero que puede ser negativo. Si se omite, se supone k=1 halla la matriz simtrica de A respecto de un eje horizontal >> A = [ 1 3 5; 5 6 7 ] A= 1 3 5 5 6 7
46
fliplr(A)
halla la matriz simtrica de A respecto de un eje vertical Cambia el tamao de la matriz A devolviendo una matriz de tamao m n cuyas columnas se obtienen a partir de un vector formado por las columnas de A puestas una a continuacin de otra. Si la matriz A tiene menos de m n elementos se produce un error. >> a=[10 22 31; 45 53 62; 78 86 97; 14 26 34]
a= 10 reshape(A,m,n) 45 78 14
22 53 86 26
31 62 97 34
>> k=reshape( a , 2 , 6) k= 10 45
78 14
22 53
86 26
31 62
97 34
Generamos una matriz aleatoria con: >> g = rand(2,4) g= 0.9501 0.2311 0.6068 0.4860 0.8913 0.7621 0.4565 0.0185
0 0 0
0 0 0
se puede generar la matriz nula y la matriz unidad del mismo tamao de f con la siguientes ordenes:
>> h = rand(3,3) h= 0.8214 0.4447 0.6154 0.7919 0.9218 0.7382 0.1763 0.4057 0.9355
1 1 1
1 1 1
>> h = rand(5,2) h= 0.9501 0.7621 0.2311 0.4565 0.6068 0.0185 0.4860 0.8214 0.8913 0.4447
47
FUNCIONES QUE ACTAN SOBRE MATRICES Las siguientes funciones exigen que el/los argumento/s sean matrices. En este grupo aparecen algunas de las funciones ms tiles y potentes de MATLAB.
Calcula la traspuesta (conjugada) de la matriz A Calcula la traspuesta (sin conjugar) de la matriz A Devuelve un vector v con los coeficientes del polinomio caracterstico de la matriz cuadrada A Devuelve la traza t (suma de los elementos de la diagonal) de una matriz cuadrada A Devuelve el nmero de filas m y de columnas n de una matriz rectangular A Devuelve el tamao de una matriz cuadrada A Devuelve el nmero de filas de A Devuelve el nmero de columnas de A
2.
Dada la siguiente matriz A, escribir la orden usando funciones conocidas para el tratamiento de matrices y tener la parte resaltada con rojo en el arreglo BB y la parte resaltada con verde en el arreglo CC: Nota: no copie los elementos de la matriz
64 9 17 40 32 41 49 8
2 55 47 26 34 23 15 58
3 54 46 27 35 22 14 59
61 12 20 37 29 44 52 5
60 13 21 36 28 45 53 4
6 51 43 30 38 19 11 62
7 50 42 31 39 18 10 63
57 16 24 33 25 48 56 1
48
3.
Sin usar la forma comn de definir una matriz; es decir: a=[2 1 1; 1 2 1; 1 1 2], solo usando comando conocidos para el tratamiento de matrices definir la siguiente matriz:
2 1 1 1 2 1 1 1 2
4. Dada la siguiente matriz A, escribir la orden usando funciones conocidas para el tratamiento de matrices y transformar dicha matriz en la matriz inferior. As mismo; la parte resaltada con rojo y lila debe transferirse a una matriz CC : Nota: no copie los elementos de la matriz
64 9 17 40 32 41 49 8
2 55 47 26 34 23 15 58
3 54 46 27 35 22 14 59
61 12 20 37 29 44 52 5
60 13 21 36 28 45 53 4
6 51 43 30 38 19 11 62
7 50 42 31 39 18 10 63
57 16 24 33 25 48 56 1
64 9 17 40 32 41 49 8
2 55 47 26 34 23 15 58
3 54 0 0 0 0 14 59
61 12 0 1 1 0 52 5
60 13 0 1 1 0 53 4
6 51 0 0 0 0 11 62
7 50 42 31 39 18 10 63
57 16 24 33 25 48 56 1
0 0 0 0
5.
0 1 1 0
0 1 1 0
0 0 0 0
Dada la siguiente matriz A, escribir la orden usando funciones conocidas para el tratamiento de matrices para EXTRAER los elementos con color a arreglos AA y BB: Nota: no copie los elementos de la matriz
49
6.
Dada la siguiente matriz A, escribir la orden usando funciones conocidas para el tratamiento de matrices para EXTRAER los elementos de color rojo a arreglos AA y BB, luego calcular la suma y producto de dichos elementos: Nota: no copie los elementos de la matriz
35 3 31 8 30 4
1 32 9 28 5 36
6 7 2 33 34 29
26 21 22 17 12 13
19 23 27 10 14 18
24 25 20 15 16 11
7.
Dada la siguiente matriz A, escribir las ordenes usando funciones conocidas para el tratamiento de matrices para tener como resultado la matriz B: Nota: no copie los elementos de la matriz
17 23 4 10 11
24 5 6 12 18
1 7 13 19 25
8 14 20 21 2
15 16 22 3 9
9 3 22 16 15
2 1 20 14 8
25 19 13 7 1
18 12 6 5 4
11 10 4 3 7
8.
Hacer un programa en matlab para establecer campos para ingresar datos de una matriz ( El usuarios debe preguntar Cuntas filas y cuntas columnas debe tener la matriz ), de acuerdo a ello se debe abrir una tabla para ingresar datos y operar los siguientes clculos (suma, resta, multiplicacin y divisin)
9.
5 5 5 5 5 5
5 5 5 5 5 5
5 5 5 5 5 5
4 4 4 4 4 4
4 4 4 4 4 4
4 4 4 4 4 4
10. Dibuje un rbol de navidad usando comandos del MATLAB (relacionado matrices u otros)
50
51
ESQUEMATIZACION SIMPLE Suponga que desea graficar un conjunto de puntos de datos, (xi, yi), i=1, 2, ...n. Es necesario preparar x y y en forma de arreglos con igual nmero de puntos, es decir, convertirlos en arreglos de fila o de columna de la misma longitud. Los datos se grafican con plot. Por ejemplo, y = sen(x)exp(-0.4x), -1<= i >= 10, se grafica con el listado 2.1
EJEMPLOS
Listado 2.1 x = -1 : 0.08 : 10; y = sin (x) .* exp(-0.3.* x); plot(x,y) xlabel( 'x' ); ylabel('y') obteniendo el siguiente resultados
Otra forma es utilizar vectores de columna en los argumentos de plot, como se muestra en el siguiente guin: Listado 2.2 x = ( -1 : 0.05 :10)' ; y = sin(x).*exp(-0.3.* x); plot(x,y) xlabel( 'x' ) ; ylabel( 'y' )
52
USO DE MARCAS: Los datos pueden graficarse slo con marcas sin estar conectados por lneas. Se dispone de cinco tipos de marcas o letras: Tipo de marca Punto Ms Estrella Simbolo . + * Tipo de marca Circulo Marca x Simbolo o x
Si desea graficar con un solo tipo de marca, coloque el smbolo de la marca como una cadena despus de las coordenadas en los argumentos de plot. La grfica del listado 2.4 se muestra en la figura 2.3
Listado 2.3 x = (0 : 0.4 : 10 ) ; y = sin(x) . * exp(-0.4 * x) ; plot ( x , y , + ) xlabel( 'x' ) ; ylabel( 'y' ) obteniendo el resultado adjunto
Tambien
53
Smbolo --
Smbolo : -.
El tipo de lnea por omisin es el continuo. Si deseas graficar con un tipo de lnea en particular, especifique la marca de lnea despus de las coordenadas; por ejemplo:
Smbolo r y m c
Smbolo g b w k
Utilice el smbolo del color igual que los tipos de lnea en el argumento de plot; por ejemplo: Plot( x , y , g ) Tambin es posible combinar marcas y colores: Plot( x , y , + g ) grafica los datos con marcas + de color verde
54
%Tambin se pueden utilizar %vectores de columna en los %argumentos de plot, como se %muestra a continuacin
%PROG23 p=0: 0.05:8*pi; z = (cos(p) + i*sin(2*p)).*exp(0.05*p) + 0.01*p; plot(real(z), imag(z)) xlabel('Re(z)') ; ylabel('Im(z)') title('PROG23') Nota: Re(z) : Parte real de un nmero complejo Imag(z) : Parte imaginaria de un nmero complejo
55
%PROG24 x=(0:0.4:10)'; y=sin(x).*exp(-0.4*x); plot(x,y,'*') xlabel('x');ylabel('y') title('PROG24: grafico con un solo tipo de marca')
%PROG25 x=(0:0.2:10)'; y=sin(x).*exp(-0.4*x); plot(x,y) grid on xlabel('x');ylabel('y') title('PROG25: Figura con reticula')
56
%PROG27 t=0.1:0.1:3; x=exp(t); y=exp(t.*sinh(t)); loglog(x,y) grid xlabel('x'); ylabel('y'); title('PROG27: Grafica semilogaritmica...')
57
%PROG32 x=0:0.05:5; y=sin(x); z=cos(x); plot(x,y ,'. -r' , x,z ,'*b') ; title('PROG32: MULTIPLES CURVAS')
58
%PROG33 clear; clf; hold off x=0:0.05:5; y=sin(x); plot(x,y,'*b'); hold on z=cos(x); plot(x,z,'--r') xlabel('x'); ylabel('y(-), z(--)'); hold off
%PROG34 M=[0:0.01:1]' k=1.4; p0_entre_p = (1 + (k1)/2+M.^2).^(k/(k-1)); plot(M,p0_entre_p) xlabel('M,numero de Mach') ylabel('p0/p') title('Relacion de presion, p(estancamiento/p(estatica)')
%PROG34-A clear; clf; hold off x=0:0.05:5; y=tan(x); plot(x,y,'*b'); hold on z=cot(x); plot(x,z,'--r'); title('Grfico Tangente - Cotangente') xlabel('x'); ylabel('y(-), z(--)'); hold off
59
%PROG35 %RELACION DE PRESION vs NUMERO DE MACH clear; clf; hold off; M=[0:0.01:0.9]'; k=1.4; p0_entre_p = (1 + (k-1)/2+M.^2) .^ (k/(k-1)); hold on axis('square'); %Hace que la grafica sea cuadrada plot(M,p0_entre_p) xlabel('M,numero de Mach') ylabel('p0/p') title('Relacion de presion, p(estancamiento/p(estatica)') text(0.35, 5.55, 'compresible') Mb = [0:0.01:0.9]'; p0_entre_pb = 1 + k/2*Mb.^2; plot(Mb, p0_entre_pb,'*') text(0.5, 0.8, 'Incompresible')
60
EJEMPLOS ADICIONALES
PROGRAMA PARA LA LEY DE BOYLE %Ley de Bolye-Mariotte PV = K %PV = nR.T n=1 y R=0.08205 t=25+273=298K clc hold off clear all; clear memory; clear command history; clc for T=298:100:698; K = 1.*0.08205.*T ; P=1:0.01:10; V = K./P ; switch T case 298 xox='-r'; case 398 xox='-k'; case 498 xox='-b'; case 598 xox='-m'; case 698 xox='-.R'; end hold on; grid on; title(' Perfiles de la Ley de Boyle - Mariotte ') ylabel(' Presin - atm '); xlabel(' Volmen - Lts '); plot(V,P,xox); end h = legend('T=298K','T=398K','T=498K','T=598K','T=698K',1); axis([0 60 0 10]) % Para el caso del grfico (P vs 1/V) figure; hold on; grid on for T=298:100:698; K = 1.*0.08205.*T ; P=1:0.01:10; V = K./P ; switch T case 298 xox='-r'; case 398 xox='-k'; case 498 xox='-b'; case 598 xox='-m'; case 698 xox='-.R'; end title(' Perfiles de la Ley de Boyle - Mariotte ') ylabel(' Presin - atm '); xlabel(' Volmen (1/V) Lts -1 '); plot((1./V),P,xox) end h = legend('T=298K','T=398K','T=498K','T=598K','T=698K',4); axis([0 0.45 0 10])
61
% Para el caso del grfico (P.V vs 1/V ) figure; hold on; grid on for T=298:100:698; K = 1.*0.08205.*T ; P=1:0.01:10; V = K./P ; PV=P.*V ; switch T case 298 xox='-r'; case 398 xox='-k'; case 498 xox='-b'; case 598 xox='-m'; case 698 xox='-.R'; end title(' Perfiles de la Ley de Boyle - Mariotte ') ylabel(' Presin*Volumen - atm*L '); xlabel(' Volmen (1/V) Lts -1 '); plot((1./V),PV,xox) end h = legend('T=298K','T=398K','T=498K','T=598K','T=698K',1);
RESULTADOS
62
63
PROGRAMA PARA LA LEY DE CHARLES %Ley de Charles V/T = K Proceso a Presin constante %PV = nR.T n=1 y R=0.08205 P=1 atm clc hold off clear all; clear memory; clear command history; clc for P=1:5:21; K = 1.*0.08205./P ; V=0.05:0.05:10; T = V./K ; switch P case 1 xox='-r'; case 6 xox='-k'; case 11 xox='-b'; case 16 xox='-m'; case 21 xox='-.R'; end hold on; grid on; title(' Perfiles de la Ley de Charles ') ylabel('Volmen - Lts '); xlabel(' Temperatura K '); plot(T,V,xox); end h = legend('P=1 atm','P=6 atm','P=11 atm','P=16 atm','P=21 atm',4); axis([0 3000 0 10]) % Para el caso del grfico (P vs 1/T) figure; hold on; grid on for P=1:5:21; K = 1.*0.08205./P ; V=0.05:0.05:10; T = V./K ; switch P case 1 xox='-r'; case 6 xox='-k'; case 11 xox='-b'; case 16 xox='-m'; case 21 xox='-.R'; end hold on; grid on; title(' Perfiles de la Ley de Charles ') ylabel('Volmen - Lts '); xlabel(' Temperatura (1/T) K-1 '); plot((1./T),V,xox); end h = legend('P=1 atm','P=6 atm','P=11 atm','P=16 atm','P=21 atm',1); axis([0 0.035 0 10])
64
% Para el caso del grfico (P.V vs 1/V ) figure; hold on; grid on for P=1:5:21; K = 1.*0.08205./P ; V=0.01:0.01:10; T = V./K ; switch P case 1 xox='-r'; case 6 xox='-k'; case 11 xox='-b'; case 16 xox='-m'; case 21 xox='-.R'; end hold on; grid on; title(' Perfiles de la Ley de Charles ') ylabel('Volmen - Lts '); xlabel(' Temperatura (t) C '); plot((T-273),V,xox) end h = legend('P=1 atm','P=6 atm','P=11 atm','P=16 atm','P=21 atm',4); axis([-500 2500 0 10])
RESULTADOS
65
66
PROGRAMA PARA EL CLCULO DEL VOLUMN CON LA ECUACION DE VAN DER WAALS % Clculo del volumen con Ecuacin de Van Der Waals Metano clear all; clear memory; clear command history; clc T=273; % K P=100; % atm Tc=190.6 % K Pc=45.4 % atm R=0.08205 b=R*Tc/(8*Pc) a=27*R*R*Tc*Tc / (64*Pc) Vm=0.1:0.1:100; PP = ((R.*T./(Vm - b)) - (a./(Vm.*Vm))); plot(Vm,PP); grid on nombre_f='( (R.*T./(x - b)) - (a./(x.*x)) )-P'; i=1; delta=0.005; x0=0.05; e=1; while e>3E-12 & i<=18 x=x0; fx0=eval(nombre_f); x=x0-delta; df1=eval(nombre_f); x=x0+delta; df2=eval(nombre_f); dfx0=(df2-df1)./(2.*delta); r=x0-(fx0./dfx0); e=abs((r-x0)./r); x0=r; i=i+1; end fprintf('la raiz es:%10.9f\n',x0); axis([0 4 0 120])
67
PROGRAMA %TRATAMIENTO TERCIARIO DE AGUAS RESIDUALES CON CARBON ACTIVADO clc;clear memory;clear command history;clear all format short g p1=[0.76 1.52 2.28]; t1=[740 1789 2780]; y1=polyfit(p1,t1,1); m1=y1(1); %pendiente ; b1=y1(2); %ordenada al origen p11=[0.0:0.01:2.40]; %plot(p11,p11*m1+b1); hold on %---------------p2=[0.76 1.52 3.04]; t2=[180 560 1330]; y2=polyfit(p2,t2,1);m2=y2(1); b2=y2(2); p22=[0.0:0.01:3.10]; %plot(p22,p22*m2+b2); %-------------p3=[1.52 3.04 4.56]; t3=[170 500 830]; y3=polyfit(p3,t3,1);m3=y3(1); b3=y3(2); p33=[0.0:0.01:5.00]; %plot(p33,p33*m3+b3); %-----------------p4=[1.52 4.56 7.60]; t4=[60 390 730]; y4=polyfit(p4,t4,1); m4=y4(1); b4=y4(2); p44=[0.0:0.01:8.00]; %plot(p44,p44*m4+b4); grid on plot(p11,p11*m1+b1,p22,p22*m2+b2,p33,p33*m3+b3,p44,p44*m4+b4); text(2.2,2560,'Experimento 1','BackgroundColor',[.7 .9 .7]);text(3,1060,'Experimento 2','BackgroundColor',[.7 .9 .7]); text(4.7,700,'Experimento 3','BackgroundColor',[.7 .9 .7]), text(5.7,400,'Experimento 4','BackgroundColor',[.7 .9 .7]); xlabel('Profundidad D (m)'); ylabel('Tiempo t (h)'); title('REPRESENTACION DE TIEMPO-t vs PROFUNDIDAD DE LECHO-D','BackgroundColor',[.7 .9 .7]) hold off %--------------------%Caudal en m3/hr y diametro del lecho de prueba piloto en mts qa1=0.00248; qa2=0.0049; qa3=0.0099; qa4=0.0198; diamlecho=0.0254; %Caudal en L/(min*m2) - el factor de conversion de conversion de m3/hr a L/(min*m2) es 21.2207 qa11=(4*1000/(60*pi))*qa1/(diamlecho^2); qa22=(4*1000/(60*pi))*qa2/(diamlecho^2); qa33=(4*1000/(60*pi))*qa3/(diamlecho^2); qa44=(4*1000/(60*pi))*qa4/(diamlecho^2); %Calculo de la velocidad en mts/hr v1=0.06*qa11; v2=0.06*qa22; v3=0.06*qa33; v4=0.06*qa44; %Calculo de la concentracion inicial de contaminante en mg/L - el factor de %conversion es 0.001 para que este dado en Kg/m3 - debe ingresar con input Coo y Cee Coo=20; Cee=1; Co=0.001*Coo; Ce=0.001*Cee; %Calculo de No en Kg/m3 No1=Co*v1*m1; No2=Co*v2*m2; No3=Co*v3*m3;
68
No4=Co*v4*m4; %Calculo de la constante de velocidad - m3/(Kg de carbon * hr) %Se trabaja con Coo=20 y Cee=1 Ambos deben ingresar con imput k1=-log((Co/Ce)-1)/(b1*Co); k2=-log((Co/Ce)-1)/(b2*Co); k3=-log((Co/Ce)-1)/(b3*Co); k4=-log((Co/Ce)-1)/(b4*Co); %Calculo de la profundidad crtica (mts) Do1=log((Co/Ce)-1)*v1 / (k1*No1); Do2=log((Co/Ce)-1)*v2 / (k2*No2); Do3=log((Co/Ce)-1)*v3 / (k3*No3); Do4=log((Co/Ce)-1)*v4 / (k4*No4); %Grafico de los datos encontrados figure Q=[qa11 qa22 qa33 qa44]; KK=[k1 k2 k3 k4]; DD=[Do1 Do2 Do3 Do4]; NN=[No1 No2 No3 No4]; subplot(2,2,1);plot(Q,KK,'*-');xlabel('Caudal-L/(min*m2)');ylabel('K - m3/(Kg*hr)');grid on subplot(2,2,2);plot(Q,DD,'O-');xlabel('Caudal-L/(min*m2)');ylabel('Do - m)');grid on subplot(2,2,3);plot(Q,NN,'*-');xlabel('Caudal-L/(min*m2)');ylabel('No - Kg/m3)');grid on %---%Ajuste de curva para encontrar la ecuacion de cada una de ellas Q1Q=[81.5:50:681.5]; k2 =polyfit(Q,KK,2); % ES IGUAL A [p,S] = polyfit(Q,KK,2) ko2=k2(1).*Q1Q.*Q1Q + k2(2).*Q1Q + k2(3); cof11=k2(1); cof22=k2(2); cof33=k2(3); %k3=polyfit(Q,KK,3) %Existe mucha dispersion, por lo que se desecha %ko3=k3(1).*Q1Q.*Q1Q.*Q1Q + k3(2).*Q1Q.*Q1Q + k3(3).*Q1Q + k3(4) n2 =polyfit(Q,NN,2); % Es igual a [p,s]=polyfit(Q,NN,2) no2=n2(1).*Q1Q.*Q1Q + n2(2).*Q1Q + n2(3); cof1=n2(1); cof2=n2(2); cof3=n2(3); %n3 =polyfit(Q,NN,3); %Existe mucha dispersion, por lo que se desecha %no3=n3(1).*Q1Q.*Q1Q.*Q1Q + n3(2).*Q1Q.*Q1Q + n3(3).*Q1Q + n3(4) %no33=(no2+no3)/2 %plot(Q,NN,'*-',Q1Q,no2,'o-',Q1Q,no3,'*-') d2 =polyfit(Q,DD,2); do2=d2(1).*Q1Q.*Q1Q + d2(2).*Q1Q + d2(3); cof111=d2(1); cof222=d2(2); cof333=d2(3); %---------INCREMENTADO figure subplot(2,2,1);plot(Q1Q,ko2,'*-');xlabel('Caudal-L/(min*m2)');ylabel('K - m3/(Kg*hr)');grid on subplot(2,2,2);plot(Q1Q,do2,'O-');xlabel('Caudal-L/(min*m2)');ylabel('Do - m)');grid on subplot(2,2,3);plot(Q1Q,no2,'*-');xlabel('Caudal-L/(min*m2)');ylabel('No - Kg/m3)');grid on %-------FIN DE INCREMENTO %DATOS DE PROBLEMA REAL %Caudal en m3/hr y diametro del lecho de prueba piloto en mts - profundidad - prof en mts qaa1=4; diamlecho1=0.76; prof=1.82; %Caudal en L/(min*m2) - el factor de conversion de conversion de m3/hr a L/(min*m2) es 21.2207 qaa11=(4*1000/(60*pi))*qaa1/(diamlecho1^2); area=(pi/4)*diamlecho1^2; %Calculo de la velocidad en mts/hr ve1=0.06*qaa11; %Modulo para halla K ko22 en m3/(kg*hr) y No no22 en kg/m3 no22=cof1*qaa11*qaa11 + cof2*qaa11 + cof3; ko22=cof11*qaa11*qaa11 + cof22*qaa11 + cof33;
69
do22=cof111*qaa11*qaa11 + cof222*qaa11 + cof333; %Calculo del tiempo de servicio t=((no22*prof)/(ve1*Co)) - (log((Co/Ce)-1)/(ko22*Co)); %horas/ciclo %Parte2:Calculo del numero de cambios de carbon por ao y volumen anual de carbon NumCamCar=365*24 / t ; VolAnualCarbon = prof*area*NumCamCar ; %Parte3:Estimacion de eliminacion del soluto (Kg/ao) %Kg soluto separado por ciclo=kg soluto en afluente-Kg soluto en efluente VolAR=qaa1*t; %Kg de soluto en el afluente por ciclo SolutAfluen=VolAR * Co; %Calculo del soluto residual - Kg de soluto por ciclo SolutEfluen=VolAR * Ce; SolutSeparado = SolutAfluen - SolutEfluen; SolutSepanual = SolutSeparado*NumCamCar; %Parte4: Rendimiento de la adsorcion CapTotAdsor=no22*VolAnualCarbon; rendimiento = (SolutSepanual/CapTotAdsor)*100; %Rendimiento basandose en Do para un caudal de 148 L/(min*mt2) Rendimiento = ((prof - do22)/prof)*100; %Los resultados disp(' '); disp('------------------------------------------------------------------') disp('RESULTADOS '); disp('------------------------------------------------------------------') disp(' APLICACION DE LOS DATOS DE LABORATORIO ') disp(' AL DISEO DE UNA PLANTA A ESCALA REAL ') disp('------------------------------------------------------------------') disp('PARTE1: CALCULO DEL TIEMPO DE SERVICIO') disp('------------------------------------------------------------------') fprintf('Capacidad de Adsorcin (No) = %g (Kg/m3)\n' ,no22) fprintf('Constante de Velocidad (K) = %g (m3/(kg*hr)\n',ko22) fprintf('Tiempo de servicio (t) = %g (hr/ciclo)\n' ,t) disp(' ') disp('-------------------------------------------------------------------') disp('PARTE2: CAMBIO DE CARBON POR AO Y VOLUMEN ANUAL DE CARBON') disp('-------------------------------------------------------------------') fprintf('# de Cambios de Carbn / ao (NumCamCar) = %g (Ciclo/ao)\n',NumCamCar) fprintf('Volmen anual de carbn (VolAnualCarbon) = %g (m3)\n',VolAnualCarbon) disp(' ') disp('-------------------------------------------------------------------') disp('PARTE3: ESTIMACION DE LA ELIMINACION DE SOLUTO (kg/ao)') disp('kg de soluto separado=Kg soluto en afluente - Kg soluto en efluente') disp('-------------------------------------------------------------------') fprintf('Volmen de Agua Residual (VolAR) = %g (m3/ciclo)\n',VolAR) fprintf('Eliminacin de soluto (SolutSeparado) = %g (kg/ciclo)\n',SolutSeparado) fprintf('Soluto separado anual (SolutSepanual) = %g (kg/ciclo)\n',SolutSepanual) disp(' ') disp('-------------------------------------------------------------------') disp('PARTE4: RENDIMIENTO DE LA ADSORCION') disp('-------------------------------------------------------------------') disp('Basado en No: Capacidad de adsorcin') fprintf('Capacidad Total de Adsorcin (CapTotAdsor) = %g (m3/ciclo)\n',CapTotAdsor)
70
fprintf('Rendimiento (rendimiento) = %g (kg/ciclo)\n',rendimiento) disp('Basado en Do: Profundidad Crtica') fprintf('Rendimiento (Rendimiento) = %g (m3/ciclo)\n',Rendimiento) disp(' ') disp('-------------------------------------------------------------------') disp(' ') disp(' ')
RESULTADOS
71
72
PROGRAMA %CAPACIDAD DE ADSORCION DE CARBON GRANULADO clc;clear memory;clear command history;clear all dc=[0 0.4 0.8 1.6 3.2 9.75]; dcA=[0 0.08 0.16 0.32 0.64 1.95]; dcB=[0 0.08 0.16 0.32 0.64 3.00]; imreA=[400 248 164 92 34 10]; imreB=[400 280 200 120 60 10]; %subplot(1,2,1) plot(dcA,imreA,'o-',dcB,imreB,'*-') title('Isoterma de adsorcion');xlabel('Dosage de carbn - gr/200mL'); ylabel('Impurezas remanentes - mg/L') text(0.8,110,'O - Carbn A');text(0.8,85,' * - Carbn B');grid on %subplot(1,2,2) figure imreA=[400 248 164 92 34 10]; imreB=[400 280 200 120 60 10]; XMA=[0 380 295 193 114 40]; XMB=[0 300 250 175 106 26]; LnimreA=log(imreA); LnimreB=log(imreB); LnXMA=log(XMA); LnXMB=log(XMB); %figure loglog(imreA,XMA,'o-',imreB,XMB,'-*') title('Isoterma de adsorcion'); xlabel('Impurezas remanentes - mg/L'); ylabel('X/M Impurezas adosrbidas/peso de carbon usado - mgr/gr'); text(100,50,'O - Carbn A');text(100,40,'* - Carbn B'); grid on %Determinacion de carga de carbon desde los perfiles de las isotermas a un %nivel de impureza final deseada - Mtodo de Interpolacin SeA=20 ; SeB=20; Carga_carA = interp1(imreA,XMA,SeA); Carga_carB = interp1(imreB,XMB,SeB); %Cculo de la impureza eliminada - Se conoce So=400 mg/L y Se=20 Mg/L SoA=400; SoB=400; % mg/L Impu_elimA = SoA - SeA; Impu_elimB = SoB - SeB; %Promedio de uso de carbon Prom_uso_carbonA = Impu_elimA / Carga_carA; Prom_uso_carbonB = Impu_elimB / Carga_carB; clc; %Los resultados disp(' '); disp('------------------------------------------------------------------') disp('RESULTADOS '); disp('------------------------------------------------------------------') disp('CAPACIDAD DE ADSORCION DE CARBON GRANULADO ') disp('------------------------------------------------------------------') disp(' ') disp('Clculo de carga de carbn desde las isotermas:')
73
disp('-----------------------------------------------') fprintf('Carbn A (X/M) = %g (mg/g)\n' ,Carga_carA) %g es el formato mas corto fprintf('Carbn B (X/M) = %g (mg/g)\n' ,Carga_carB) %g es el formato mas corto disp(' ') disp('Promedio de uso de Carbn :') disp('---------------------------') fprintf('Uso de Carbn A = %g (mg/L)\n',Prom_uso_carbonA) fprintf('Uso de Carbn B = %g (mg/L)\n',Prom_uso_carbonB) disp(' ') disp('----------------------------------------------------------------------------') disp('De los resultados calculados se concluye que se requiere de aproximadamente') fprintf('%g \n',Prom_uso_carbonB - Prom_uso_carbonA ) disp('mas de carbn B para levar a cabo el mismo tratamiento. El precio unitario ') disp('de cada tipo de carbn puede ser usado para determinarel mejor caso de ') disp('costo-efectividad para una aplicacin de purificacin ') disp('----------------------------------------------------------------------------')
RESULTADOS
74
75
1.
2. Mostrar las ordenes para obtener el siguiente grfico. Entre -3 y 3 hay 50 puntos
3. Dibujar en una sola ventana con dos subventanas las funciones y = x2 -3*x 2 y z = x3 2*x + 1 con x variando entre -5 y 5 a intervalos de 0.02. La primera en lneas verdes continuas, la segunda en lneas azules discontinuas. Marcar con un '+' rojo el punto (2,1) en ambas grficas. Escribir las instrucciones
5.
Se han usado diverso comandos no tratados en la clase correspondientes a grficos con Matlab, es necesario que Ud. Con ayuda el HELP del MATLAB, mencione el significado
76
7. El mecanismo de un pistn conectado mediante una varilla a una manivela es un problema clsico utilizado frecuentemente en aplicaciones de ingeniera. En el mecanismo que se muestra en la figura adjunta, la manivela tiene una velocidad constante de rotacin de 500 rpm. Calcular y representar grficamente la posicin, la velocidad y aceleracin del pistn para una de las revoluciones de la manivela. Represente tres grficos distintos. Consideras theta=0 en el instante t=0
8. Un barco A viaja hacia el sur a una velocidad de 8 millas por hora, mientras que otro barco B viaja hacia el este a una velocidad de 16 millas por hora. A la 7 am los barcos estan a las distancias que se indica en la figura adjunta. Represente la distancia entre los barcos en funcin del tiempo para las siguientes cuatro horas. El eje horizontal debe mostrar el tiempo actual del da, comenzando por las 7 am, mientras que el eje vertical debe mostrar la distancia. Ponga titulo y etiquetas a los ejes. Si la visibilidad es de 8 millas, estime la hora a partir de la cual las personas de un barco pueden ver a las del otro. NOTA: los problemas 6 y 7 fueron tomados del texto MATLAB una introduccin con ejemplos prcticos - Amos Gilat
9. Grafique diez funciones relacionados con los cursos de Ingeniera Qumica, Fsica u otro, deben tener cierto grado de dificultad con conceptos de simulacin
77
GRAFICOS 3D
OBJETIVO
Generar grficos tri-dimensionales de diversas funciones establecidas con el comando plot3 del MatLab
78
GRAFICOS 3D
La orden plot que utilizbamos en 2 dimensiones podemos tambin usarla en 3 dimensiones. El formato ser el mismo que en 2D, pero en vez de representar los datos en parejas, lo haremos en tripletes. El formato generalizado de la orden plot en 3D es: plot3(x,y,z,s) "x", "y", "z", sern los tres vectores, y "s" ser una cadena de caracteres adicional para indicar color y tipo de lnea. Si queremos representar varias curvas en un solo grfico, la orden ser de la forma:
plot3(x1,y1,z1,s1,x2,y2,z2,s2,...xn,yn,zn,sn)
Por ejemplo, si queremos representar una hlice en tres dimensiones: % Programa para representar una hlice en tres dimensiones clc; clear memory; clear command history; clear all; t=0:0.1:30; plot3(sin(t),cos(t),2*t) title('Hlice') xlabel('Seno(t)'); ylabel('Coseno(t)'),zlabel('2t')
Hlice
Para fijar manualmente los lmites de los ejes, usaremos la misma orden que en 2D, pero teniendo en cuenta al eje z: >>axis([xmin xmax ymin ymax zmin zmax] ) Del mismo modo, podemos introducir texto en la grfica, usando coordenadas en 3D: >>text(x,y,z,'texto')
79
Por ejemplo, si en la grfica anterior queremos limitar los ejes al doble de su valor actual, e introducir la etiqueta "Grfica en 3D" en las coordenadas x=0.25, y=0.25, z=20: % Programa para representar una hlice en tres dimensiones clc; clear memory; clear command history; clear all; t=0:0.1:30; plot3(sin(t),cos(t),2*t) title('Hlice') xlabel('Seno (t)'); ylabel('Coseno (t)'),zlabel ('2 t') axis([-2 2 -2 2 0 120]); text(0.25,0.25,100,'GRAFICA EN 3D') grid
Hlice
Grficos de Malla Un grfico tridimensional de malla viene definido por una funcin z=f(x,y), de tal forma que los puntos de la superficie se representan sobre una rejilla, resultado de levantar los valores de z dados por f(x,y) sobre los correspondientes puntos (x,y) del plano. Los puntos adyacentes se unen mediante lneas rectas, dndolo al grfico un aspecto de malla. Para representar un grfico de malla, se utiliza el comando mesh y sus variantes. Estos grficos son muy tiles a la hora de visualizar grandes matrices o representar funciones de dos variables. El primer paso para representar una funcin de dos variables z=f(x,y), es utilizar el comando meshgrid, que bsicamente define la matriz de puntos (X,Y) sobre los cuales se evala la funcin de dos variables para hacer su representacin grfica. Su sintaxis es la siguiente: [X,Y]=meshgrid(x,y) Transforma el campo de definicin dado de las variables "x" e "y", de la funcin a representarz=f(x,y), en argumentos matriciales utilizables por el comando mesh para obtener el grfico de malla.
80
Disponemos de tres tipos de grficos de malla: mesh(X,Y,Z) Representa el grfico de malla de la funcin z=f(x,y). Podemos aadir un argumento opcional C que indique el color de la rejilla.
meshz(X,Y,Z) Representa el grfico de malla de la funcin z=f(x,y) con una especie de teln o cortina en la parte inferior
meshc(X,Y,Z) Representa el grfico de malla de la funcin z=f(x,y) junto con el grfico de contorno correspondiente (curvas de nivel proyectadas sobre el plano x,y)
El programa que introduciremos en MATLAB ser: % Programa para representar grfico tipo malla clc; clear memory; clear command history; clear all; [X,Y]=meshgrid(-7.5:0.5:7.5); Z=sin(sqrt(X.^2+Y^2))./(sqrt(X.^2+Y.^2)+eps); % Le sumaremos esp al denominador para no dividir por cereo % (eps es la cantidad ms pequea que maneja MATLAB, mayor que cero) subplot(2,2,1) mesh(X,Y,Z) title('Comando mesh') subplot(2,2,2) meshz(X,Y,Z) title('Comando meshz') subplot(2,2,3) meshc(X,Y,Z) title('Comando meshc')
81
Grficos de Contorno Las grficas de contorno o sistema de planos acotados muestran lneas de elevacin o altura constante. As, dibujando diferentes curvas de nivel correspondientes a alturas constantes, se puede describir un mapa de lneas de nivel de la superficie. Son las grficas usadas en los mapas topogrficos. Los grficos de contorno pueden representarse en dos o en tres dimensiones En estos grficos se observa, por tanto, la variacin de z=f(x,y) con respecto a "x" y a "y". Cuando el espacio entre las curvas de nivel es grande, significa que la variacin de la variable "z" es lenta, mientras que un espacio pequeo indica un cambio rpido de "z". Los comandos que utiliza MATLAB para la representacin de grficos de contorno son los siguientes (todos ellos aceptan el argumento adicional sobre el color de las lneas): contour(Z) Dibuja el grfico de contorno para la matriz Z. El nmero de lneas de contorno a utilizar se elige automticamente contour(Z,n) Dibuja el grfico de contorno para la matriz Z usando n lneas de contorno
contour(x,y,Z,n) Dibuja el grfico de contorno para la matriz Z usando en los ejes X e Y el escalado definido por los vectores "x" e "y", con n lneas de contorno contour3(Z), contour3(Z,n), contour3(x,y,Z,n) dimensiones Dibujan los grficos de contorno en tres
82
pcolor(X,Y,Z) Dibuja un grfico de contorno para la matriz (X,Y,Z) utilizando una representacin basada en densidades de colores. Suele denominarse grfico de densidad. Como ejemplo, vamos a representar las curvas de nivel que responden a la ecuacin: z=Sen(x)Sen(y) con -2 < x,y < 2
En la primera grfica vemos su grfico de contorno en 2D, con 20 lneas. En la segunda, su grfico de contorno tridimensional, con 40 lneas. En la tercera dibujaremos su grfico de densidad: % Programa para representar grfico tipo contorno clc; clear memory; clear command history; clear all; [X,Y]=meshgrid(-2:0.1:2); Z=sin(X).*sin(Y); subplot(2,2,1) contour(Z,20,'b') title('GRAFICA DE CONTORNO 2D') subplot(2,2,2) contour3(Z,40,'b') title('GRAFICA DE CONTORNO 3D') subplot(2,2,3) pcolor(X,Y,Z) title('GRAFICA DE DENSIDAD')
83
Polgonos en Tres Dimensiones MATLAB permite dibujar polgonos en tres dimensiones. Para ello utiliza los siguientes comandos: fill3(X,Y,Z,C) Dibuja el polgono compacto cuyos vrtices son las tripletas de componentes (Xi,Yi,Zi) de los vectores columna X, Y y Z. C es un vector de la misma dimensin de X, Y y Z, que contiene los colores Ci de cada punto(Xi,Yi,Zi). Si C es solo un carcter, se pintarn todos los puntos del polgono del color correspondiente al carcter.
Un ejemplo de este tipo de grficos tridimensionales sera: % Programa para REPRESENTAR polgono en tres dimensiones clc; clear memory; clear command history; clear all; x=cos(0:0.01:8*pi); y=sin(0:0.01:8*pi); z=(0:0.01:8*pi); fill3(x,y,z,'r')
84
85
2.
3.
PROCESAMIENTO DE DATOS
En la actualidad es muy comn solucionar problemas con ayuda del computador, sabemos que existen muchos datos que son posibles tratarlos para poderlos tener como informacin que pueda ayudarnos para tomar decisiones importantes en la vida diaria
DATOS
INFORMACION
PROBLEMA
DISEO ALGORITMO
DEL
PROGRAMA DE COMPUTADORA
86
ALGORITMOS Y PROGRAMAS
La resolucin de problemas exige al menos los siguientes pasos: 1. 2. 3. 4. Definicin o Anlisis del problema Diseo del Algoritmo Transformacin del algoritmo en un programa Ejecucin y validacin del programa
ESPECIFICACIONES DE ENTRADA
ESPECIFICACIONES DE SALIDA
Es una frmula o mtodo para la resolucin de problemas, podemos ayudarnos de diagramas de flujo CARACTERISTICA DE LOS ALGORITMOS:
Debe ser preciso e indicar el orden de ejecucin de cada paso Debe ser definido. Si se hace dos veces se debe obtener el mismo resultado Debe ser finito. Debe terminar en algn momento
Por ejemplo:
CALCULAR
SUPERFICIE (CIRCULO) LONGITUD (CIRCUNFERENCIAS)
CALCULO DE SUPERFICIE
S = PI * R * R
CALCULO DE LA CIRCUNFERENCIA S
SALIDA DE RESULTADOS R L
L = 2 * PI * R
ENTRADA PROCESO SALIDA
87
DATOS
PROCESO
INFORMACION
BASURITA
PROCESO
BASURITA
PROCESO DE PROGRAMACION
D O C U M E N T A C I O N
M A N T E N I M I E N T O
DEPURACION Y VERIFICACION
88
PARTES DE UN PROGRAMA
INSTRUCCION 1 INSTRUCCION 2 INSTRUCCION 3 INSTRUCCION 4 INSTRUCCION 5 INSTRUCCION 6 ................... .................. .................. INSTRUCCION 7 INSTRUCCION 8
INSTRUCCIONES DE ENTRADA INSTRUCCIONES CONDICIONALES INSTRUCCIONES DE PROCESO INSTRUCCIONES REPETITIVAS INSTRUCCIONES DE SALIDA
NOTA: EN ALGUNOS CASOS EXISTEN PROGRAMAS CON INSTRUCCIONES DE UNA SOLA PASADA
INSTRUCCIONES CONDICIONALES
ES CIERTO ?
NO
SI
Instruccione 1 Instruccione 2 Instruccione 3 Instruccione 4
O
Instruccione 1 Instruccione 2 Instruccione 3 Instruccione 4
INSTRUCCIONES REPETITIVAS
FIN DE LA CONDICION
INSTRUCCION n-2 INSTRUCCION n-1 INSTRUCCION n
FIN
89
DIAGRAMAS DE FLUJO
Los diagramas de flujo representan la forma de especificar y documentar los detalles algortmicos de un producto de programacin. Tambin podemos decir que el diagrama de flujo es la fotografa de la solucin de un problema. Estos diagramas utilizan cajas rectangulares para especificar las acciones, cajas en forma de rombos para las proposiciones de decisin, arcos dirigidos para las interconexiones entre las diversas cajas, as como una variedad de formas especiales para denotar las entradas, las salidas, los almacenamientos,
entre datos y en funcin del resultado de la misma determina cual de las alternativas se tomar. Normalmente tiene dos salidas
diagrama de flujo, el primero sirve para conectar dos funciones dentro de la misma pgina, mientras que el segundo conecta diagramas de pginas diferentes
sirve para indicar que los resultados de la operacin se observar impreso en un papel
FUNCION SALIDA.-
sirve para indicar que los resultados de la operacin se observar en el monitor del computador
FUNCION SALIDA.-
90
DIAGRAMAS DE FLUJO En un inicio, la programacin estructurada fue desarrollada en sus principios por Edsgar W. Dijkstra en sus Notes on Structured Programming y se basa en el denominado Teorema de la Estructura desarrollado en 1966 por Bmh y Jacopini, que se ratific con los trabajos de Charlan D. Mills. En la actualidad existen diversas definiciones de estos diagramas, pero todas ellas giran alrededor del teorema de estructura que, como ya hemos dicho, se debe a Bmh y Jacopini que inician todo esto con esta tcnica de programacin a travs de mdulos o bloques. Para un buen entendimiento del teorema mencionado, se realiza una definicin previa de algunos de los conceptos que trata el teorema:
1.
Diagrama Propio. Es aquel que posee un solo punto de entrada y uno de salida.
C1
S B
A
N
2.
INICIO
Programa Propio. Es aquel programa que cumple las siguientes condiciones: Posee un solo inicio y un solo fin. Todo elemento del programa es accesible, es decir, existe al menos un camino desde el inicio al fin que pasa a travs de l. No posee bucles infinitos.
C1
A
N
FIN
91
CARACTERISTICAS DE LOS DIAGRAMAS DE FLUJO En los distintos departamentos de informtica existentes no siempre se dispone de los mismos programadores con respecto al tiempo que se pretende que dure una aplicacin, por lo cual es de suma importancia que un programa realizado por una persona sea fcil de modificas u mantener por otra. En este sentido, la diagramacin estructurada ofrece muchas ventajas para logras estos objetivos. Con esto podemos decir que: Un diagrama estructurado es: Fcil de leer y comprender. Fcil de codificar en una amplia gama de lenguajes y en diferentes sistemas. Fcil de mantener, Eficiente, aprovechando al mximo los recursos de la computadora. Es posible tener los programas de una aplicacin por mdulos.
Los diagramas de flujo estructurados difieren de los diagramas tradicionales en que los primeros tienen restriccin en cuanto a las formas de uso; con las actuales se obtiene que la grfica obtenida sea un equivalente grfico de la descripcin del problema por medio del seudo cdigo estructurado; un ejemplo de las formas comunes y de los equivalentes en seudo cdigo son:
A S
F
P V S
V F P S1
S2
Los diagramas estructurados poseen una entrada nica y una salida nica; as estas formas pueden ser anidadas dentro de otras formas hasta el nivel deseado de anidamiento.
92
APLICACIONES DE LOS DIAGRAMAS DE FLUJO. Los diagramas de flujo estructurados son actualmente caracterizados como una herramienta de la programacin estructurada. Gracias a esta herramienta podemos interpretar cada accin de un programa y representarlo grficamente.
2.
Estructura Alternativa. Es una estructura con una sola entrada y una sola salida en la cual se realiza una accin de entre varias, segn una condicin o se realiza una accin segn el cumplimiento o no de una determinada condicin. Esta condicin puede ser simple o compuesta. Las estructuras alternativas pueden ser:
De dos salidas, en la que una de ellas puede ser la accin nula. representado por el comando IF .then..else end
En el Matlab esta
S
COND
COND
A N
Alternativa Simple
Alternativa doble
93
De tres o ms salidas, que tambin se llama mltiple. En el Matlab esta representado por el comando switch
Sintaxis COND. switch switch_expr case case_expr statement,...,statement case case_expr statement,...,statement ... otherwise statement,...,statement end
Alternativa Mltiple
3.
Estructura Repetitiva. Es una estructura con una entrada y una salida en la cual se repite una accin un nmero determinado o indeterminado de veces, dependiendo en este caso del cumplimiento de una condicin. Las estructuras repetitivas pueden ser:
94
4.
Estructura MIENTRAS (WHILE). En esta estructura se repite una accin mientras se cumpla la condicin que controla el bucle. La caracterstica principal de esta estructura es la de que la condicin es evaluada siempre antes de cada repeticin. El nmero de repeticiones oscila entre 0 e infinito, dependiendo de la evaluacin de la condicin, cuyos argumentos en los casos de repeticin, al menos una vez, debern modificarse dentro del bucle, pues de no ser as el nmero de repeticiones ser infinito y nos encontraremos en un bucle sin salida.
Y = Y + 1 Y = Y + 1
95
Es una variable cuyo valor se incrementa o decrementa en una cantidad constante en cada iteracin. Tiene la funcin de contar los sucesos o acciones del bucle. Una forma de controlar un bucle es mediante un contador N = N+1
CONCEPTO DE ACUMULADOR
S = S + N
CONCEPTO DE BUCLE
Es un segmento de programa o conjunto de sentencias u ordenes de un programa que se repiten un nmero determinado de veces mientras se cumple una determinada condicin
BUCLES INDEPENDIENTES
BUCLES ANIDADOS
96
EJERCICIO 1: Disear el diagrama de flujo que permita escribir los nmeros del 1 al 1000
INICIO
I=1
SI
Es I <=1000 NO Escribir I
I = I + 1
FIN
INICIO
I=1 S=0
SI
Es I <=1000 S = S + I I = I + 1
NO
Escribir S
FIN
97
EJERCICIO 3: Escribir el diagrama de flujo para calcular el factorial de un nmero, ingresado por teclado sea mayor que 1 considerar que el nmero
INICIO
Leer N
SI
NO
P = P * N N = N - 1
NO
Es N = 1
SI
ESCRIBIR P
1 FIN
98
Ejercicio 1-1 Preparar un algoritmo que represente la receta que se muestra a continuacin: Ingredientes: 1. 2. 3. 4. 5. 6. 7. 8. 9. 1 taza de almendras picadas. de libra de chocolate en bloque para hornear de libra de malvaviscos cortados a la mitad 3 tazas de azcar taza de leche evaporada taza de miel de maz 1 cucharadita de vainilla libra de mantequilla cucharadita de sal
Virtase la leche y aada la miel de maz, el azcar, el chocolate y la sal en un recipiente de 1 litro, y cocnese en flama alta, mezclando constantemente hasta que hierva la mezcla. Redzcase, a flama moderada y continese hirviendo y revolviendo hasta que una gota de la miel forme una pelota suave en un vaso de agua fra. Qutese de la flama y djesele enfriar durante 10 minutos. Mzclese la mantequilla y la vainilla hasta que estn completamente incorporadas. Agrguense las almendras. Distribyanse las mitades de malvavisco en el fondo de una charola de horneado, de 30 centmetros por lado, engrasada. Virtase la miel sobre los malvaviscos. Djesele enfriar durante 10 minutos. Crtese en cuadros y srvase. TOMADO DE: http://www.uady.mx/~contadur/fundamentos_de_programacion/ejercicios1.htm
Ejercicio 2.1 Un maestro asign a sus estudiantes el problema de construir un diagrama de flujo en la forma siguiente. La entrada consista en las longitudes y anchuras de varios rectngulos. El objeto era producir una lista con lneas numeradas consecutivamente que mostrasen la longitud L, el ancho W y el rea A, nicamente de aquellos rectngulos cuyo permetro fuese mayor que 12. Los estudiantes entregaron los diagramas de flujo que se encuentran en los diagramas adjunto del 1 al 8, como soluciones al problema. a. Aplicar la prueba a cada uno de los diagramas (debers entregar la evidencia de cada prueba) y dgase cules soluciones son correctas y cules incorrectas (para nuestro objeto una solucin correcta es la que produce el resultado deseado. Sin embargo, puede no ser la ms eficiente). En cuanto a las soluciones incorrectas, cules son los errores? Construir un diagrama de flujo que utilice las caractersticas ms eficientes de cada uno (el programa ms eficiente es el que requiere el menor nmero de clculos).
b. c.
1 2 3 4
99
100
Tomado para clases de: http://www.uady.mx/~contadur/fundamentos_de_programacion/ejemplos/Ejemplo_Comparar%20Dia gramas1.jpg NOTA: Ambos smbolos representan entrada de datos
101
Ejercicios propuestos
1. 2. 3. 4. 5. 6. 7. 8.
Hacer un diagrama para calcular el rea de un triangulo. Hacer un diagrama para convertir de grados centgrados a grados Fahrenheit. Hacer un diagrama para imprimir la suma de los nmeros del 1 al 100. Hacer un diagrama que te pida un nmero y te diga si es par, es non y/o es primo. Hacer un diagrama para imprimir la sucesin de Fibonacci. Hacer un diagrama que pida 10 nmeros y muestre el promedio. Hacer un diagrama que pida 3 nmeros y diga cual es el mayor. Hacer un diagrama que pida la edad y despliegue si es menor de edad (<18), mayor (>=18) o si pertenece a la 3 edad.(>=60)
9.
Hacer un diagrama que te pida un nmero y te diga si es par, es non y/o es primo.
10. Hacer un diagrama para calcular el factorial de un nmero. 11. Hacer un diagrama que calcule e imprima N nmeros primos. 12. Hacer un diagrama que solicite 4 calificaciones y diga si est reprobado o no, segn las reglas de tu escuela. 13. Hacer un diagrama que pida un nmero N y despliegue todas las combinaciones de dos nmeros que sumados den N. 14. Hacer un diagrama que despliegue la tabla de multiplicar de un nmero X. 15. Hacer un diagrama que calcule la probabilidad de que dos dados lanzados sumen 7. 16. Hacer un diagrama que pida 100 nmeros y diga cual es la mediana. 17. Hacer un diagrama que solicite los datos de una matriz de 4x4 y la muestre invertida. 18. Hacer un diagrama que pida 3 nmeros y calcule el comn denominador. 19. Hacer un diagrama que llene una matriz de 3x3 y despliegue los valores de la diagonal principal. 20. Hacer un diagrama que pida 2 matrices y despliegue el producto cruz de las mismas.
102
TOMADO DE:
http://www.uady.mx/~contadur/fundamentos_de_programacion/ejemplos/Ejemplo_Ciclos1.jpg
103
Cinta magntica
Disco magntico
Conector de pagina
Lneas de flujo
Anotacin
104
106
COMANDO IF.END
Sirve para realizar clculos con consultas o toma de decisiones, ello de acuerdo con la respuesta de la consulta que puede ser VERDADERA o FALSA. La consulta puede ser con respecto a una, dos o tres variables
if
EXPRESION
Else
end
107
OPERADORES NECESARIOS:
Less than less than or equal to greater than greater than or equal to equal to not equal to Element-wise Element-wise (<), (<=), (>), (>=), (==), (~=) AND (&) OR (|)
A S
V
P
F
P S1 S2
S
V
If p = 20 Hacer S end
108
PROBLEMA 1 Escribir el cdigo de programa y el diagrama de flujo para determinar si un nmero ingresado por el teclado es par o impar SOLUCION %Hacer un programa que permita ingresar un nmero cualquiera por teclado y %determine si dicho nmero es par o impar, tener en consideracin que no deben ingresar 0 ni nmeros negativos clc;clear memory;clear command history; ni=input('ingrese numero inicial : '); par=0; impar=0; %mensaje=' '; if ni<=0 mensaje='numero ingresado no es mayor que cero'; else if rem(ni,2)~=0 mensaje='IMPAR'; else mensaje='PAR'; end end disp(' ');disp(' '); [mensaje]
109
INICIO
INGRESA VALOR DE N
DEFINIR VARIABLES
n=VAL(TxtN)
SI
Es N <= 0
NO
SI IF REM(N, 2) = 0
NO
TxtRpta = Es Par
TxtRpta = Es Impar
110
PROBLEMA 2
Escribir un programa que permita ingresar cuatro dgitos: A, B, C y D de un entero positivo N . Se desea redondear N a la centena ms prxima y visualizar el resultado. Por ejemplo: Si A es 2, B es 3, C es 6 y D es 2, entonces N ser 2362 y el resultado redondeado ser 2400. Si N 2342 el resultado ser 2300 y si N es 2962 el resultado ser 3000.
SOLUCION
% PROGRAMA DE REDONDEO clear all; clear memory; clear command history; clc a = input('Ingrese el primer digito '); b = input('Ingrese el segundo digito '); c = input('Ingrese el tercer digito '); d = input('Ingrese el cuarto digito '); if a > 0 x=1000*a + 100*b + 10*c + d k=10*c + d; if k>49 y = x + (100-k); else y = x - k; end end [(y)]
111
INICIO
SI X = 1000 * a + 100 * b + 10 * c + d k = 10 * c + d
SI IF k>49 Y = X + ( 100 k )
NO
Y = X - k
NUMERO REDONDEADO
FIN
112
PROBLEMA 3 Supongamos que Ud. Solicita bolas de billar de diferentes colores pero que tendra que ayudarse de un programa para verificar el color, tal que: 1=rojo, 2=azl, 3=verde, 4=amarillo y 5=negro Escribir el cdigo de programa:
SOLUCION Se digitar el cdigo de programa y tendr las siguientes ordenes % Uso de if multiple % 1=rojo, 2= azul, 3=verde, 4=amarillo 5=negro clear all; clear memory; clear command history; clc; color = input('Ingrese el numero que represente un color....'); if color == 1 color1='rojo' else if color == 2 color1='azul' else if color == 3 color1='verde' else if color == 4 color1 = 'amarillo' else if color == 5 color1 = 'negro' else color1 = 'no dado' end end end end end [color1]
113
INICIO
INGRESE BOLA
N
COND .
S S
COND .
N
COND .
ROJO
AZUL
S
VERDE
COND .
S
AMARILLO
N
COND .
S
NEGRO
NO DADO
COLOR1
FIN
114
OTRA FORMA DE CODIGO: % Uso de if multiple % 1=rojo, 2= azul, 3=verde, 4=amarillo 5=negro clear all; clear memory; clear command history; clc; color = input('Ingrese el numero que represente un color....'); switch color case 1 color1='rojo' case 2 color1='azul' case 3 color1='verde' case 4 color1='amarillo' case 5 color1='negro' otherwise color1='no dado' end [color1] DIAGRAMA DE FLUJO:
INICIO
COND.
FIN
115
OTRA FORMA
% Uso de if multiple % 1=rojo, 2= azul, 3=verde, 4=amarillo 5=negro clear all; clear memory; clear command history; clc; color = input('Ingrese el numero que represente un color....'); if color == 1 color1='rojo' end if color == 2 color1='azul' end if color == 3 color1='verde' end if color == 4 color1 = 'amarillo' end if color == 5 color1 = 'negro' else color1 = 'no dado' end [color1]
116
INICIO
INGRESE BOLA
COND .
S
ROJO
COND .
S
AZUL
COND .
S
VERDE
COND .
S
AMARILLO
N
COND .
S
NEGRO
NO DADO
COLOR1
FIN
117
1.
Elaborar el diagrama de flujo y el programa en MATLAB para que reciba dos 2 nmeros enteros determine: El mayor de los dos nmeros y lo imprima. El menor de los dos nmeros y lo imprima.
2.
Elaborar el diagrama de flujo y el programa en MATLAB para que reciba 3 nmeros enteros determine: El mayor de los tres nmeros y muestre en pantalla. El menor de los tres nmeros y muestre en pantalla. El intermedio de los tres nmeros y muestre en pantalla.
3.
Elaborar el diagrama de flujo y el programa en MATLAB para que reciba 4 nmeros enteros (teniendo en consideracin que los cuatro pueden ser iguales, tres pueden ser iguales y dos pueden ser iguales) y determine: El mayor de los cuatro nmeros y muestre en pantalla. El menor de los cuatro nmeros y muestre en pantalla. Los intermedios en orden ascendente o descendente y muestre en pantalla.
4.
Elaborar el diagrama de flujo y el programa en MATLAB que determine si tres nmeros enteros (a, b y c) mayores que 0 representan los lados de un tringulo. El programa deber de imprimir SI si los lados forman un tringulo y en caso contrario, se deber imprimir NO. Nota: tener en consideracin las expresiones siguientes: (a-c)<b<(a+c), (a-b)<c<(a+b) y (b-c)<a<(b+c)
5.
Elaborar el diagrama de flujo y el programa en MATLAB para que dado los 3 lados de un tringulo (a, b yc), imprima R si el tringulo es rectngulo, I si el tringulo es issceles y E si el tringulo es escaleno. Nota: tener en consideracin las expresiones siguientes: (a-c)<b<(a+c), (a-b)<c<(a+b) y (b-c)<a<(b+c)
6.
Elaborar el diagrama de flujo y el programa en MATLAB para que reciba 4 notas de las evaluaciones de un curso e imprima APROBADO si el promedio de notas es mayor o igual a 10.5, en caso contrario se imprimir DESAPROBADO, tener en cuenta lo siguiente: El Promedio de notas es igual a la suma de las 4 notas dividido entre 4. El Promedio de notas es igual a la suma de las 4 notas menos la ms baja, dividido entre 3.
Derechos reservados
118
COMANDO FOREND
OBJETIVO
Conocer y aplicar el concepto de estructura repetitiva y multi- repetitiva (Comando for-end) en la solucin de problemas especficos.
119
COMANDO FOR
Sirve para realizar clculos repetitivos, repitiendo las sentencias que se encuentran dentro del bucle un nmero especfico de veces. Segn las necesidades de solucin del problema puede haber un solo bucle o puede haber varios bucles for-end
for
VARIABLE = 1: N: M
end
120
Estructura Repetitiva. Es una estructura con una entrada y una salida en la cual se repite una accin un nmero determinado o indeterminado de veces, dependiendo en este caso del cumplimiento de una condicin. Las estructuras repetitivas pueden ser: Estructura para (FOR) simple y anidado
El caso dos es muy usado para construir la tabla de multiplicacin del 1 al 12 en un solo paso. El caso tres podra darse en un proceso de simulacin simple.
121
PROBLEMA:
Hacer un programa y el diagrama de flujo que permita calcular la suma de la serie: ingresando el valor de n S = (1/2)+(2/2*2)+(3/2*2*2)+...+(n/2*2*2*2*.*n)
SOLUCION % Uso del comando FOR clear all; clear memory; clear command history; clc; n = input('ingrese el valor de "n" : '); su=0; if n > 0 for I = 1:n ; su = su + I / (2 ^ I); end disp(' ') disp(' ') [su] else disp(' ') disp(' ') disp('Ingresar valores de "n" mayores de cero') end
INICIO
NO IMPRIMIR
Es N>0
SI
ERROR
FOR I = 1 to N
SU = SU + I / (2I)
IMPRIMIR
SU
DIAGRAMA DE FLUJO
FIN
122
PROBLEMAS 2 Hacer el diagrama de flujo y escribir el cdigo de programa para hacer la tabla del 1 al 12, SOLUCION % Uso del comando FOR % Tabla de multiplicar clear all; clear memory; clear command history; clc; for n=1:12; disp(' ') fprintf(' TABLA DEL = %g \n ',n) ; disp(' ') for m=1:12 ; p=m*n; disp([m,n,p]); end end
INICIO
FOR n = 1 to 12
FOR m = 1 to 12
p = m * n
IMPRIMIR
m, n, p
FIN
123
PROBLEMA
% solubilidad del oxigeno disuelto en funcion % de la temperatura y salinidad en KCl - mg/L clc clear clear memory clear history xlabel('Temperatura - C'), ylabel('Concentracin de oxigeno disuelto, mg/L') title('solubilidad del oxigeno disuelto en funcin de la temperatura y salinidad en KCl - mg/L') for cl=0:5:25; T=273.15:0.5:323.15; C=-139.344+(1.5757e+5)./T-(6.6423e+7)./(T.^2)+(1.2438e+10)./(T.^3)-(8.6219e+11)./(T.^4)-cl.* (3.1929e-2 -(19.428./T) + 3867.3./(T.^2)); C=exp(C); hold on plot(T-273.15,C) end grid
124
PROBLEMA: Un nmero perfecto es un entero positivo que es igual a la suma de todos los enteros positivos (excluido el mismo) que son divisores del nmero. El primer nmero perfecto es 6, ya que los divisores son 1, 2, 3 y 1 + 2 + 3 = 6. Esquematice el Diagrama de flujo y escribir un programa que encuentre los nmeros perfectos en el rango de 2 al nmero n ingresado por teclado
% calcula la cantidad de nmeros perfectos clear memory; clear command history; clear all; clc; XX=input('Ingrese un numero cualquiera : '); if XX>0 np = 0; % Calcula la cantidad de numeros perfectos entre 1 y XX for II = 1:XX; suma = 1; z=round(II/2); for I = 2:z; %Realiza un bucle para saber cuantos divisores if rem(II,I)== 0 %Tiene el nmero que se esta probando suma = suma + II / I; %Acumulacin de todo los divisores end end if II == suma %Verifica si la suma de los divisores es Idntico al nmero probado a ser perfecto np = np + 1; %Cuenta la cantidad de nmeros perfectos existentes en el rango disp([II]) end end disp(['la cantidad de numeros perfectos es = ',num2str(np)]) fprintf('La cantidad de numeros perfectos es = %g \n' ,np)%np representa la cantidad de nmeros perfectos en un rango else disp(' ') disp(' ') disp(['debe ingresar cantidad mayor que cero']) en
125
INICIO XX
SI Es XX<= 0 NO
II = 1 to XX
SUMA = 1
I = 2 to round(II/I)
SI
Es II mod I = 0 NO
SUMA = SUMA + II / I
Es SUMA = II
SI
NP = NP + 1
FIN
126
PROBLEMA DE AMONIACO %PROGRAMA PARA EL PROCESO DE AMONIACO clear all; clear memory; clear command history; clc; ce1 = 2; ce2 = 3; ce3 = 1; T=298; To=298; P=1; for P=1:1:5; T1=300:20:700; a1 = 2.73.*10; b1 = 2.38E-2; c1 = 1.71E-5; d1 = -1.19E-8; a2 = 2.71.*10; b2 = 9.27E-3; c2 = -1.38E-5; d2 = 7.65E-9; a3 = 3.12.*10; b3 = -1.36E-2; c3 = 2.68E-5; d3 = -1.17E-8; H=-45700; syms T dCp = (ce1*a1-ce2*a2-ce3*a3)+(ce1*b1-ce2*b2-ce3*b3)*T+(ce1*c1-ce2*c2ce3*c3)*T^2+(ce1*d1-ce2*d2-ce3*d3)*T^3; I=H-((ce1*a1-ce2*a2-ce3*a3)*To+(ce1*b1-ce2*b2-ce3*b3)/2*To^2+(ce1*c1-ce2*c2ce3*c3)/3*To^3+(ce1*d1-ce2*d2-ce3*d3)/4*To^4); DH=I+(ce1*a1-ce2*a2-ce3*a3)*T+(ce1*b1-ce2*b2-ce3*b3)/2*T^2+(ce1*c1-ce2*c2ce3*c3)/3*T^3+(ce1*d1-ce2*d2-ce3*d3)/4*T^4; %DH1=I+(ce1.*a1-ce2.*a2-ce3.*a3).*T1+(ce1.*b1-ce2.*b2-ce3.*b3)./2.*T1.^2+(ce1.*c1-ce2.*c2ce3.*c3)./3.*T1.^3+(ce1.*d1-ce2.*d2-ce3.*d3)./4.*T1.^4; R = 8.314; % En joules Go = -16160; % En joules To = 298; ko = exp(-Go/(R*To)); Ka=exp(33.3398+(3642.8./T1)-6.97.*log(T1)+0.0020086.*T1+0.0000009776.*T1.^20.00000000035041.*T1.^3); Tc1 = 405.7; Tc2 = 33.19; Tc3 = 126.2; Pc1 = 111.32; Pc2 = 12.96; Pc3 = 33.55; Pr1 = P./Pc1; Pr2 = P./Pc2; Pr3 = P./Pc3; Tr1 = T1./Tc1; Tr2 = T1./Tc2; Tr3 = T1./Tc3; w1=0.253; w2=-0.216; w3=0.038; Bo1=0.083-0.422./Tr1.^1.6; Bo2=0.083-0.422./Tr2.^1.6; Bo3=0.083-0.422./Tr3.^1.6; B11=0.139-0.172./Tr1.^4.2; B12=0.139-0.172./Tr2.^4.2; B13=0.139-0.172./Tr3.^4.2; Fi1=2.7182.^((Pr1./Tr1).*(Bo1+w1.*B11)); Fi2=2.7182.^((Pr2./Tr2).*(Bo2+w2.*B12)); Fi3=2.7182.^((Pr3./Tr3).*(Bo3+w3.*B13)); Fit=Fi1.^ce1./((Fi2.^ce2).*(Fi3.^ce3)); nombre_f='(4.*Fit.*(x.^2).*(4-2.*x).^2./((1-x).*(P.^2).*(3-3.*x).^3))- Ka'; i=1; delta=0.05; x0=0.88; e=1; while e>3E-12 & i<=18 x=x0; fx0=eval(nombre_f); x=x0-delta; df1=eval(nombre_f); x=x0+delta; df2=eval(nombre_f); dfx0=(df2-df1)./(2.*delta); r=x0-(fx0./dfx0); e=abs((r-x0)./r); x0=r; i=i+1;
127
end %fprintf('la raiz es:%10.9f\n',x0); hold on switch P case 1 xox='-r'; case 2 xox='-k'; case 3 xox='-b'; case 4 xox='-m'; case 5 xox='-R'; end plot(T1,x0,xox); h = legend('P=1 atm','P=2 atm','P=3 atm','P=4 atm','P=5 atm',1); xlabel('Temperatura - K'); ylabel('Extensin de la reaccin'); title('EXTENSION VS TEMPERTAURA PARA EL PROCESO AMONIACO') axis([250 700 0 1]) grid on P,[T1' x0' Fit'] end %fprintf('la raiz es:%10.9f\n',x0)
RESULTADOS
128
PRESION= 1 atm P= 1 300.0000 320.0000 340.0000 360.0000 380.0000 400.0000 420.0000 440.0000 460.0000 480.0000 500.0000 520.0000 540.0000 560.0000 580.0000 600.0000 620.0000 640.0000 660.0000 680.0000 700.0000 0.8267 0.7710 0.7079 0.6393 0.5677 0.4960 0.4271 0.3633 0.3062 0.2565 0.2142 0.1787 0.1493 0.1251 0.1053 0.0891 0.0758 0.0648 0.0559 0.0485 0.0424
PRESION = 2 atm P= 2 300.0000 0.8771 320.0000 0.8365 340.0000 0.7895 360.0000 0.7368 380.0000 0.6796 400.0000 0.6192 420.0000 0.5575 440.0000 0.4964 460.0000 0.4375 480.0000 0.3824 500.0000 0.3320 520.0000 0.2869 540.0000 0.2473 560.0000 0.2128 580.0000 0.1833 600.0000 0.1581 620.0000 0.1366 640.0000 0.1185 660.0000 0.1031 680.0000 0.0901 700.0000 0.0791
PRESION = 3 atm P=3 300.0000 320.0000 340.0000 360.0000 380.0000 400.0000 420.0000 440.0000 460.0000 480.0000 500.0000 520.0000 540.0000 560.0000 580.0000 600.0000 620.0000 640.0000 660.0000 680.0000 700.0000 0.8998 0.8664 0.8274 0.7831 0.7343 0.6818 0.6268 0.5707 0.5149 0.4607 0.4092 0.3615 0.3179 0.2787 0.2440 0.2135 0.1870 0.1640 0.1441 0.1270 0.1122
PRESION= 4 atm P=4 300.0000 320.0000 340.0000 360.0000 380.0000 400.0000 420.0000 440.0000 460.0000 480.0000 500.0000 520.0000 540.0000 560.0000 580.0000 600.0000 620.0000 640.0000 660.0000 680.0000 700.0000 0.9135 0.8845 0.8504 0.8116 0.7683 0.7214 0.6715 0.6198 0.5674 0.5154 0.4649 0.4168 0.3719 0.3306 0.2931 0.2594 0.2295 0.2032 0.1800 0.1598 0.1421
PROBLEMA: Disear el diagrama de flujo y cdigo de programa: para determinar la cantidad de nmero primos existentes en un rango (entre 1 y el nmero ingresado por teclado)
129
INICIO
NO
Es NN>0
SI
XX = 0
n = 1 To NN
X = 0
Es X = 0
SI
XX = XX + 1
FIN
130
% Determinar la cantidad de nmeros primos % Entre 1 y un nmero ingresado por el teclado % Numero primo es aquel que es divisible por la unidad y el mismo numero clc; clear memory; clear command history; clear all; NN=input('Ingrese un numero cualquiera : '); if NN > 0 XX = 0; for n = 1:NN % Bucle para buscar nmeros primos en el rango de 1 a nn X = 0; for p = 2:(n-1) % Busca otro divisor aparte del 1 y el mismo numero if rem(n,p)== 0 % Determinacin del numero de divisores X = X + 1; % Si x>0 el numero cuenta con mas de dos divisores end % por lo tanto no es primo end if X == 0 % Si el numero de divisores es cero, entonces es primo XX = XX + 1; % Cuenta los nmeros primos disp([n]) end end fprintf('La cantidad de nmeros primos es = %g \n' ,XX) %np representa la cantidad de nmeros perfectos en un rango else disp(' ') disp(' ') disp(['debe ingresar cantidad mayor que cero']) end
131
sen ( x ) = x
x3 x5 x7 x9 + + .... 3! 5! 7! 9! x2 x4 x6 x8 + + .... 2! 4! 6! 8!
cos ( x ) = 1
2.
Elaborar el diagrama de flujo y el programa correspondiente para que pueda leer un nmero entero entre 11 y 9999 y como resultado deber mostrar los dgitos que forman dicho nmero, indicando si el dgito es par o impar. Guiarse por el ejemplo mostrado. Si N=8638 es el nmero ingresado 8 par 6 par 3 impar 8 par
3.
Elaborar el diagrama de flujo y el programa correspondiente para que lea un nmero entre 1 y 9999999 y una posicin determinada. El objetivo es determinar el dgito que est ocupando la posicin indicada. Por ejemplo: 856384 y 5; sean el nmero y posicin ingresada, la respuesta indicar: El dgito que ocupa la posicin (5) es el (8)
8 1
5 2
6 3
3 4
8 5
4 6
NOTA: El programa deber tener la consistencia entre el nmero ingresado y las posiciones para evitar errores, es decir; no se puede ingresar 4567 y digitar la posicin 7, porque no existe.
4.
Elaborar el diagrama de flujo y el programa correspondiente para que lea un nmero entre 1 y 9999999, una posicin determinada y un dgito entre 0 y 9. El objetivo es: una vez identificado el dgito que ocupa la posicin determinada deber cambiarse dicho dgito por el dgito ingresado por el teclado y escribir el nuevo nmero generado. Por ejemplo: 947390, 3 y 8; sean el nmero, la posicin de un digito de dicho nmero y el nuevo dgito, la respuesta indicar: El dgito que ocupa la posicin (3) es el (7) y deber cambiarse por (8) y el nuevo nmero es: 948390
9 1 9
4 2 4
7 3 8
3 4 3
9 5 9
0 6 0
132
NOTA: El programa deber tener la consistencia entre el nmero ingresado y las posiciones para evitar errores, es decir; no se puede ingresar 45678 y digitar la posicin 8, porque no existe.
5.
Elaborar el diagrama de flujo y el programa para mostrar lo k elementos de la serie de Fibonacci, as como la suma de dicho elementos. Sabemos que la serie de Fibonacci se establece teniendo en consideracin que los dos primeros elementos de la serie son el cero (0) y uno(1). El resto se calcula como la suma de los dos elementos que los anteceden. Ello puede observarse en la tabla adjunta. 0 1er ter 1 2do ter 1 1=0+1 2 2=1+1 3 3=2+1 5 5=2+3 8 8=3+5 13 13=5+8 21 21=8+13 34 34=13+21
6.
Elaborar un diagrama de flujo y el programa correspondiente para determinar si un nmero entero entre cero (0) y 9999999 es un nmero Amstrong. En caso que sea deber mostrar en pantalla el resultado de la siguiente forma: El nmero (N) es un nmero Amstrong. Y en el caso que no lo sea deber mostrar: El nmero (N) no es un nmero Amstrong NOTA: Un nmero (N) es Amstrong si la suma de los cubos de sus dgitos es igual a (N)
7.
Elaborar el diagrama de flujo y el programa correspondiente para saber cuntos nmeros Amstrong existen en el rango establecido como datos de entrada al iniciar el programa
8.
Elaborar el diagrama de flujo y el programa correspondiente para saber si un nmero ingresado por teclado de la PC es un nmero Amstrong. Si el nmero ingresado es Amstrong se mostrar el siguiente resultado: El nmero (N) es un nmero Amstrong y el prximo nmero Amstrong esta despus de (x) nmeros Ejemplo: N = 371
371
407
36 RESPUESTA: El nmero 371 es Amstrong y el siguiente se encuentra despus de 36 nmeros En el caso: que el nmero ingresado no sea un nmero Amstrong, encontrar el siguiente nmero Amstrong ms cercano, indicando cuanto se tiene que restar o sumar para hallar el siguiente nmero Amstrong. Ejemplo: N = 398
371
398
407
27 9 RESPUESTA: Al nmero ingresado hay que sumar 9 para encontrar el siguiente Nmero Amstrong ms cercano
133
COMANDO WHILE
OBJETIVO
Conocer y aplicar el concepto de estructura repetitiva y multi- repetitiva incluyendo expresiones de control (Comando While -end) en la solucin de problemas especficos.
134
COMANDO WHILE
El bucle while es similar al bucle for en el que permite repetir la ejecucin de las sentencias Matlab, la sintaxis del bucle while tiene la siguiente forma:
while expression % MATLAB command 1 % MATLAB command 2 % Mas comandos para ejecutar repetidamente hasta que la expresin sea no verdadera end
Donde la expresin es una expresin lgica y que puede ser falsa o verdadera, por ejemplo: consideramos el siguiente bucle while n=1 while n < 3 n = n+1 end Este cdigo crea la siguiente salida: n= n= n= 1 2 3
Tambin
count = 1; num = 2; while count < 16 cnum= [count num] disp(cnum) num = num*2; count = count + 1; if count >= 10 break end end % expression is modified % forms a vector % pretty printing
135
NO
Condicin Verdadera ?
SI
FIN
Menor que Mayor que Menor o igual que Mayor o igual que Igual No igual
& ~
Y O No
136
PROBLEMA PARA CALCULAR LA PRESION DE ROCIO Encontrar la presin de roco y la composicin en la fase lquida para una mezcla vapor de 55% de acetona (1) y 45% de hexano (2) a la temperatura de 20C. Datos: Pc (bar) 47,0 30,1 Vc cm3/mol 209 370
Zc 0.232 0,264
0.304 0,299
A 7,11714 6,91058
B 1210,595 1189,64
C 229,664 226,280
B , T: C y Psat : mmHg. ( T + C)
137
ALGORITMO
Leer temperatura (T1) y fraccin molar gaseosa (Y1) Constantes Hacer i = 1.0
sat
Gammai = 1.0
Calcular P
P = 1 / ( Yi * i
/ Gammai* Pi
sat
) )
Xi = ( Yi * i * P ) / (Gammai * Pi
sat
Evaluar i
sat
No
Es Gammai < 0
Calcular Presin:
Si
P = 1 / ( Yi * i
/ Gammai* Pi
sat
No
Es P < 0
Si
Pi Xi
138
% Programa para el problema - (1) Acetona (2) Hexano clear all; clear command history; clc; y1 = 0.55; y2 = 0.45;
T = 20 + 273;
% Clculo de la Psat Con Constantes de ANTOINE A1 = 7.11714; B1 = 1210.595; C1 = 229.664; A2 = 6.91058; B2 = 1189.64; C2 = 226.280; Psat1 = exp(A1 - B1 / (T + C1)); Psat2 = exp(A2 - B2 / (T + C2)); % Calculo de la presion fi1=1; fi2=1; gama1=1.0; gama2=1.0; P=1/( y1*fi1/(gama1*Psat1) + y2*fi2/(gama2*Psat2) ) % Calculo de X1 y X2 x1=y1*fi1*P/(gama1*Psat1);
x2=y2*fi2*P/(gama2*Psat2);
% Calculo de gama1 y gama2 MARGULES A12 = 1.7448; A21 = 1.8012; gama1 = exp(x2 ^ 2 * (A12 + 2 * x1 * (A21 - A12))) gama2 = exp(x1 ^ 2 * (A21 + 2 * x2 * (A12 - A21))) % Calcular P PP=1/( y1*fi1/(gama1*Psat1) + y2*fi2/(gama2*Psat2));
gama11=0; gama22=0;
while (PP-P)>0.00001
% Evaluar fi1 y fi2 Tc1 = 508.1; Tc2 = 507.5; Pc1 = 47. * 750.0617; Pc2 = 30.1 * 750.0617; w1 = 0.304; w2 = 0.299; Tr1 = T / Tc1; Tr2 = T / Tc2; Pr1 = PP / Pc1; Pr2 = PP / Pc2; B01 = 0.083 - 0.422 / Tr1 ^ (1.6); B02 = 0.083 - 0.422 / Tr2 ^ (1.6); B11 = 0.139 - 0.172 / Tr1 ^ (4.2); B12 = 0.139 - 0.172 / Tr2 ^ (4.2); fi1 = exp(Pr1 / Tr1 * (B01 + w1 * B11)); fi2 = exp(Pr2 / Tr2 * (B02 + w2 * B12));
while gama11-gama1>0.00000001
%Nuevo calculo de x1 y x2 xx1 = y1*fi1*PP/(gama1*Psat1); xx2 = y2*fi2*PP/(gama2*Psat2) xtot = xx1 + xx2; xx1 = xx1/(xx1+xx2); xx2=xx2/(xx1+xx2); %forzando a 1.0 xtot1 = xx1 + xx2; x1 = xx1; x2 = xx2 %Recalculo de gama1 y gama2 gama1=gama11; gama2=gama22; A12 = 1.7448; A21 = 1.8012; gama11 = exp(x2 ^ 2 * (A12 + 2 * x1 * (A21 - A12))) gama22 = exp(x1 ^ 2 * (A21 + 2 * x2 * (A12 - A21)))
end
P=PP; PP = 1/( y1*fi1/(gama1*Psat1) + y2*fi2/(gama2*Psat2))
end
[ x1 x2 P ]'
RESULTADOS= COMPOSICION DE ACETONA COMPOSICION DE HEXANO PRESION DE ROCIO : : : 0.5048 0.4952 176.6582
139
140
Donde f1 es el factor de horas-hombre laboradas en un ao para un trabajador y f2 es el factor del salario anual del trabajador vs. el salario anual de todos los trabajadores de la empresa. F1 = hF / (hF x cF + hlT1 x clT1 + hlT2 x clT2) F2 = sF / (sF x cF + clT1 x slT1 + clT2 x slT2) Siendo cada una de las variables que intervienen en las frmulas anteriores: hF: Horas por da que trabaja un funcionario sF: Salario que recibe el funcionario cF: Cantidad de funcionarios que existen en la empresa hIT1: Horas por da que trabaja un ingeniero tipo 1 cIT1: Cantidad de ingenieros tipo 1 que existen en la empresa sIT1: Salario que recibe el ingeniero tipo 1 hIT2: Horas por da que trabaja un ingeniero tipo 2 cIT2: Cantidad de ingenieros tipo 2 sIT2: Salario que recibe el ingeniero tipo 2
2.
Implementar un programa en VBA, que nos diga el da de la semana que corresponda a una fecha ingresada en forma de dia (xx), mes (yy), ao (yyyy) da mes ao da de la semana
(Respuesta)
Un mtodo para determinar qu da de la semana correspondiente a un da, mes y ao determinado es el siguiente: A. Tome los dos ltimos dgitos del ao y divdalo entre 4, quedndose con la parte entera como resultado. B. Aada el da ingresado. C. Luego, aada uno de los siguientes valores, de acuerdo al mes ingresado: - Abril, julio: 0
141
- Enero, octubre: 1 - Mayo: 2 - Agosto: 3 - Febrero, marzo, noviembre: 4 - Junio: 5 - Septiembre, diciembre: 6 D. Reste 1 si el mes ingresado es enero o febrero y el ao ingresado es bisiesto. E. De acuerdo al ao ingresado (asuma que el ao ingresado se encuentra en alguno de los rangos anteriores), aada: - Para los 1700 - 1799 4 - Para los 1800 - 1899 2 - Para los 1900 - 1999 0 - Para los 2000 - 2999 6 F. Aada al resultado los ltimos dos dgitos del ao. G. Por ltimo, divida entre 7 y tome el resto. El resultado ser un nmero entre 0 y 6, que corresponder a un da de la semana: 0 sbado, 1 domingo, 2 lunes, 3 martes, 4 mircoles, 5 jueves, 6 viernes. NOTA: Un ao es bisiesto si ocurre que: Es divisible por 4 y no por 100 Es divisible por 400.
3.
Se tiene una diana (bull) para el juego de dardos, formada por 5 zonas delimitadas por circunferencias concntricas de radio 2, como se observa en la figura. Cada punto Px, est conformado por (x, y) o datos de ingreso.
Ejemplos: Observe el grfico y la Tabla mostrados anteriormente. El punto P1 est en la zona 1, circunferencia de radio 2, entonces el puntaje ser 50 puntos. El punto P2 est en la zona 3, limitada por la circunferencia de radio 4 y 6, entonces el puntaje ser 30 puntos. El punto P3 est en el lmite de la zona 4 y 5 (se considera la menor zona, zona 4), entonces el puntaje ser 20 puntos. El punto P4 est fuera de la diana, entonces el puntaje ser 0 (cero) puntos. Para lo cual se pide desarrollar un Diagrama de flujo y el Programa correspondiente para que
142
nos devuelva el mensaje correspondiente, segn los casos siguientes: Si x e y son nmeros positivos, nos devuelva las coordenadas ingresadas pertenecen al primer cuadrante. Si x e y son nmeros negativos, nos devuelva las coordenadas ingresadas pertenecen al tercer cuadrante. Si x es un nmero positivo e y es un nmero negativo, nos devuelva las coordenadas ingresadas pertenecen al cuarto cuadrante. Si y es un nmero positivo y x es un nmero negativo, nos devuelva las coordenadas ingresadas pertenecen al segundo cuadrante. Y el puntaje alcanzado en cada tiro (ingreso de dato)
4. Hacer un programa para saber si un nmero es feliz o infeliz. A partir de ello puede hacerse
variaciones, por ejemplo: cuantos nmeros felices hay en un rango establecido por el usuario o Cuntos nmeros felices o infelices hay entre 0 y 1000000. NOTA: Un nmero natural es Nmero feliz: cuando cumple lo siguiente: si sumamos los cuadrados de sus dgitos y seguimos el proceso con los resultados obtenidos el resultado es 1. Por ejemplo, el nmero 203 es un nmero feliz ya que 22+02+32=13; 12+32=10; 12+02=1. En caso contrario es nmero infeliz.
7. Hacer un programa para saber cuntos nmeros afortunados existen entre 0 y 100000 0 y
1000000 NOTA: los nmeros que sobreviven despus de una serie de operaciones se llama Nmeros afortunados: por ejemplo: tomemos la secuencia de todos los naturales a partir del 1: 1, 2, 3, 4, 5, Tachemos los que aparecen en las posiciones pares. Queda: 1, 3, 5, 7, 9, 11, 13, Como el segundo nmero que ha quedado es el 3 tachemos todos los que aparecen en las posiciones mltiplo de 3. Queda: 1, 3, 7, 9, 13, Como el siguiente nmero que qued es el 7 tachamos ahora todos los que aparecen en las posiciones mltiplos de 7. As sucesivamente. Los nmeros que sobreviven se denominan nmeros afortunados. Tambin; un nmero afortunado, que fue el nombre que le dio Reo Franklin Fortune, es un nmero primo que puede resultar de la expresin: Q = (Pn) q, donde (Pn), es el producto de los primeros (n) primos y (q) es el nmero primo ms pequeo, pero mayor que (Pn) + 1.
143
Segn la conjetura de Fortune el nmero Q siempre ser primo, pero no todos los nmeros primos pueden resultar de esta frmula. Slo los primos que pueden tomar el valor de Q, se les llama, Nmeros afortunados. Si n = 7, los primeros siete primos son (2, 3, 5, 7, 11, 13 y 17), cuyo producto es 510510, (P7) = 510510; el primo menor, pero ms grande que 510511 es 510529, q = 510529, luego Q = 510510 510529 = 19, entonces 19 es un nmero afortunado
9. Hacer el Diagrama de flujo y el programa correspondiente para calcular la suma de los cubos de
los dgitos de un nmero natural ingresado por el teclado, el nmero ingresado puede ser de 3, 4, 5 y 6 dgitos, pero no de 1 ni de 2.
10. Hacer 10 problemas relacionado con la Ingeniera Qumica haciendo uso de comando WHILE
144
ECUACIONES LINEALES
OBJETIVO
Aplicacin del MatLab en la solucin de ecuaciones lineales incluyendo mtodos conocidos para tal fin.
145
SISTEMA DE ECUACIONES LINEALES Consideramos un conjunto de ecuaciones de n incgnitas dado por: a1,1 x1 + a1,2 x2 + a1,3 x3 + + a1,n xn = y1 a2,1 x1 + a2,2 x2 + a2,3 x3 + + a2,n xn = y2 . am,1 x1 + am,2 x2 + am,3 x3 + + am,n xn = ym Donde: ai,j son coeficientes conocidos, x i son incgnitas yi trminos conocidos (trminos no homogneos) Las ecuaciones lineales anteriores se pueden expresar en forma compacta como: A x = y
A =
x =
a1,1 a1,n : a a x1 : xn y1 : yn
Podr solucionarse mediante: x=A\y x = inv(A) * y z=y/A inv(A) * y A ^ (-1) * y Y * inv(A) (El resultado esta en forma de vector columna) (El resultado esta en forma de vector columna) (El resultado esta en forma de vector fila ) (El resultado esta en forma de vector columna) (El resultado esta en forma de vector columna) (El resultado esta en forma de vector fila)
y =
= 14 = 10 = 9
>> A = [ 1 2 3; 2 1 2; 2 1 1] A= 1 2 3 2 1 2 2 1 1
>> y = [14 ; 10 ; 9 ] y= 14 10 9
146
PARTE TEORICA
Un sistema de ecuaciones lineales es un conjunto de ecuaciones lineales que podemos escribir de forma tradicional as:
un sistema as expresado tiene "m" ecuaciones y "n" incgnitas, donde aij son nmeros reales, llamados coeficientes del sistema, los valores bm son nmeros reales, llamados trminos son las variables del sistema, independientes del sistema, las incgnitas xj y la solucin del sistema es un conjunto ordenado de nmeros reales (s1, s2, ..., sn) tales que al sustituir las incgnitas x1, x2, ... , xn por los valores s1, s2, ..., sn se verifican a la vez las "m" ecuaciones del sistema. Este mismo sistema de ecuaciones lineales en notacin matricial tiene esta forma:
Donde:
Llamamos matriz del sistema a la matriz de dimensin mn formada por los coeficientes del sistema, y la designamos por A. Designamos por X a la matriz columna formada por las incgnitas. Denotamos por B a la matriz columna formada por los trminos independientes.
As mismo; llamamos matriz ampliada de dimensin m(n+1) a la matriz que se obtiene al aadir a la matriz del sistema (= matriz de coeficientes) la columna de los trminos independientes, y la denotamos por A*, es decir
147
APLICACIONES CON PROGRAMAS Programas Aplicados a Los Sistemas de Ecuaciones Lineales %mtodo de la biseccin para los sistemas de ecuaciones fprintf(\n); nombre_f=input( ingrese funcin asociada f(x)= ,s); a=input( ingrese limite inferior: ); b=input( ingrese limite superior: ); fprintf(\n); fprintf( it a b aprox error\n); i=1;e=1; r=0; while e>=3E-6 & i<=10 va=r; r=(a+b)/2; x=a;fa=eval(nombre_f); x=b;fb=eval(nombre_f); x=r;fr=eval(nombre_f); fprintf( %3.0f %10.6f %10.6f %10.6f,I,a,b,r); if fa*fr<=0 b=r; e=abs((r-va)/r); fprintf( %10.6f\n,e); else a=r; e=abs((r-va)/r); fprintf( %10.6f\n,e); end i=i+1; end fprintf(\n); fprintf(la raz es: %10.9f\n, r);
SALIDA: ingrese funcin asociada f(x)= x^3+3*x ingrese lmite inferior: 4 ingrese lmite superior: 8 it 1 2 3 4 5 6 7 8 9 10 a 4.000000 6.000000 7.000000 7.500000 7.750000 7.875000 7.937500 7.968750 7.984375 7.992188 b 8.000000 8.000000 8.000000 8.000000 8.000000 8.000000 8.000000 8.000000 8.000000 8.000000 aprox 6.000000 7.000000 7.500000 7.750000 7.875000 7.937500 7.968750 7.984375 7.992188 7.996094 error 1.000000 0.142857 0.066667 0.032258 0.015873 0.007874 0.003922 0.001957 0.000978 0.000489
148
EJERCICIO %mtodo de la regla falsa fprintf('\n'); nombre_f=input(' ingrese funcin asociada f(x)= ','s'); a=input(' ingrese lmite inferior: '); b=input(' ingrese lmite superior: '); fprintf('\n'); fprintf(' it a b aprox error\n'); i=1;e=1; r=0; while e>=3E-6 & i<=10 va=r; r=(a+b)/2; x=a;fa=eval(nombre_f); x=b;fb=eval(nombre_f); x=r;fr=eval(nombre_f); fprintf(' %3.0f %10.6f
%10.6f %10.6f',i,a,b,r);
if fa*fr<=0 b=r; e=abs((r-va)/r); fprintf(' %10.6f\n',e); else a=r; e=abs((r-va)/r); fprintf(' %10.6f\n',e); end i=i+1; end fprintf('\n'); fprintf('la raz es: %10.9f\n', r);
SALIDA: ingrese funcin asociada f(x)= sqrt(x)+2*x ingrese lmite inferior: 2 ingrese lmite superior: 6 it 1 2 3 4 5 6 7 8 9 10 a 2.000000 4.000000 5.000000 5.500000 5.750000 5.875000 5.937500 5.968750 5.984375 5.992188 b 6.000000 6.000000 6.000000 6.000000 6.000000 6.000000 6.000000 6.000000 6.000000 6.000000 aprox 4.000000 5.000000 5.500000 5.750000 5.875000 5.937500 5.968750 5.984375 5.992188 5.996094 error 1.000000 0.200000 0.090909 0.043478 0.021277 0.010526 0.005236 0.002611 0.001304 0.000651
149
EJERCICIO % mtodo de Newton Raphson clc; nombre_f = input(' ingrese funcin asociada f(x) = ','s'); x0 = input(' ingrese valor inicial : '); fprintf('\n'); fprintf(' it aprox g(x) error\n'); i=1; e=1; delta=0.001;
while e>=3E-12 & i<=18 x=x0; fx0=eval(nombre_f); x=x0-delta; df1=eval(nombre_f); x=x0+delta; df2=eval(nombre_f); dfx0=(df2-df1)/(2*delta); r=x0-(fx0/dfx0); e=abs((r-x0)/r); fprintf(' %3.0f %10.6f %10.6f %10.6f\n',i,x0,r,e); x0=r; i=i+1; end fprintf('la raz es: %10.9f\n', x0);
SALIDA: Ingrese funcin asociada f(x) = sqrt(x)+2*x^3 Ingrese valor inicial : 12 it aprox 1 12.000000 2 7.996659 3 5.324967 4 3.538690 5 2.338336 6 1.520431 7 0.941642 8 0.489185 9 0.055119 10 -0.054335 11 0.054324 12 -0.053579 13 0.053535 14 -0.052828 15 0.052757 16 -0.052087 17 0.051994 18 -0.051359 g(x) 7.996659 5.324967 3.538690 2.338336 1.520431 0.941642 0.489185 0.055119 -0.054335 0.054324 -0.053579 0.053535 -0.052828 0.052757 -0.052087 0.051994 -0.051359 0.051247 error 0.500627 0.501729 0.504785 0.513337 0.537942 0.614659 0.924921 7.875092 2.014432 2.000056 2.013915 2.000502 2.013405 2.000890 2.012908 2.001225 2.012426 2.001513
150
Tambin podr solucionarse el problema de manera grfica, as % Grfico funcin clear all; clear memory; clear command history; clc; x=10.9:0.1:11.2 y=0.9.*x.^2+8.*x-200 plot(x,y) grid on
SOLUCION
3 2 1 0 -1 -2 -3 -4 -5 -6 10.9
10.95
11
11.05
11.1
11.15
11.2
11.25
151
152
153
%Y variable dependiente: Y = b1 + m1*X b=(Sy*Sx2-Sx*Sxy)/(N*Sx2-Sx^2) m=(N*Sxy-Sx*Sy)/(N*Sx2-Sx^2) xx=-3:1:14; plot(xx,b+m*xx,x,y,'r-') text(5,2.8,' \leftarrow Y=b+m*X ','EdgeColor','red','LineWidth',2); text(9.5,3.5,['b = ',num2str(b)]) text(9.5,2.4,['m = ',num2str(m)]) xlabel('X'); ylabel('Y') axis([-6 15 -10 15])
154
yy=-5:1:10; plot(yy,b1+m1*yy,'b',y,x,'k*-') text(-3.5,-7,' \leftarrow X=b1+m1*Y ','EdgeColor','black','LineWidth',2); text(1.5,-6.4,['b = ',num2str(b1)]) text(1.5,-7.4,['m = ',num2str(m1)]) xlabel('X'); ylabel('Y') axis([-6 14 -10 15])
CUYOS RESULTADOS SE DAN EN LA FIGURA ADJUNTA: SE PIDE CONFIRMAR LOS RESULTADOS, explicando paso por paso lo que dice cada sentencia del programa
155
156
EJERCICIOS PROPUESTOS: 1. Es necesario revisar el Help del MatLab para verificar el significado de cada uno de los comandos usados en tcnica de ajuste de curvas usando el mtodo de mnimos cuadrados.
157
158
POLINOMIOS - Funciones para clculos con polinomios En MATLAB un polinomio se puede definir mediante un vector de coeficientes. Por ejemplo, el polinomio:
X4 8X2 + 6X 10 = 0
se puede representar mediante el vector [1, 0, -8, 6, -10]. MATLAB puede realizar diversas operaciones sobre l, como por ejemplo evaluarlo para un determinado valor de x (funcin polyval()) y calcular las races (funcin roots())
>> pol = [ 1 0 -8 6 10 ] pol = 1 0 -8 6 >> roots(pol) ans = -3.2800 2.6748 0.3026 + 1.0238i 0.3026 - 1.0238i
-10
Para calcular el producto de polinomios MATLAB utiliza una funcin llamada conv() (de producto de convolucin). En el ejemplo se va a ver cmo se multiplica un polinomio de segundo grado por otro de tercer grado, asi mismo; Para dividir polinomios existe otra funcin llamada deconv(). Por ejemplo; sean los polinomios siguientes: Pol1 = X2 2X + 4 y pol2 = X3 + 3X 4
>> pol1 = [ 1 -2 4 ] pol1 = 1 -2 4 >> pol2 = [ 1 0 3 4 ] pol2 = 1 0 3 -4 >> pol3 = conv(pol1,pol2) pol3 = 1 -2 7 -10
El polinomio ser:
159
polinomio caracterstico de la matriz A races del polinomio pol evaluacin del polinomio pol para el valor de x. Si x es un vector, pol se evala para cada elemento de x evaluacin del polinomio pol de la matriz A producto de convolucin de dos polinomios p1 y p2 divisin del polinomio p por el polinomio q. En c se devuelve el cociente y en r el resto de la divisin descompone el cociente entre p1 y p2 en suma de fracciones simples calcula la derivada de un polinomio calcula la derivada de producto de polinomios calcula los coeficientes de un polinomio p(x) de grado n que se ajusta a los datos p(x(i)) ~= y(i), en el sentido de mnimo error cuadrtico medio. calcula el valor interpolado para la abscisa x a partir de un conjunto de puntos dado por los vectores xp e yp. como el anterior, pero permitiendo especificar tambin el mtodo de interpolacin. La cadena de caracteres m admite los valores 'nearest', 'linear', 'spline', 'pchip', 'cubic' y 'v5cubic'.
polyfit(x,y,n)
interp1(xp,yp,x)
interp1(xp,yp,x,'m')
AJUSTE POLINOMIAL:
Si se dispone de n+1 puntos (xi,yi), se pueden encontrar un polinomio nico de orden n . Estos coeficientes ser determinados por la funcin polyfit . Sean los puntos (x,y): (1,2), (2,5), (3,10), (4,17), (5,26), la aplicacin estar dado por:
EDU>> x = [ 1 2 3 4 5 ] ; EDU>> y = [ 2 5 10 17 26 ] ; EDU>> a = polyfit( x , y , length(x) 1 ) a= 0.0000 0.0000 1.0000 0.0000 1.0000
y = x2 + 1
160
El concepto nos sirve para ajustar curvas cuando tenemos un conjunto de puntos tales como:
Xi = [590 1100 1750 2400 3200 4100 5800 6700 7800 8900 11000]; Vsz = [3.306 2.71 2.13 1.65 1.16 0.82 0.38 0.30 0.22 0.13 0.08];
Graficando siguiente:
obtendremos
lo
Para hacer un ajuste manual: hacemos click en Tools y clic en Basic Fitting y aparecer el men que aparece en la siguiente pgina
161
El estudiante podr probar ajustes de 5th, 6th y 7th grado polinomial y escribir sus conclusiones Tambin podemos escribir cdigo de programa para ajustar curvas, un caso particular es el mtodo de mnimos cuadrados para un ajuste lineal. A continuacin se muestra un ejemplo tpico: %Ejemplo - Ajuste Lineal %Y variable dependiente: Y = b + m*X clc;clear memory; clear command history x=[1 3 4 6 8 9 11 14] y=[1 2 3.8 4 5 7 8 9] plot(x,y,'o-r') xlabel('X'); ylabel('Y') axis([-6 15 -10 15]) grid on figure plot(y,x,'o-b') xlabel('X'); ylabel('Y') axis([-6 14 -10 15]) grid on figure Sx=sum(x); Sy=sum(y); Sx2=sum(x.^2); Sy2=sum(y.^2); Sxy=sum(x.*y); N=length(x) %Y variable dependiente: Y = b1 + m1*X b=(Sy*Sx2-Sx*Sxy)/(N*Sx2-Sx^2) m=(N*Sxy-Sx*Sy)/(N*Sx2-Sx^2) xx=-3:1:14; plot(xx,b+m*xx,x,y,'r-') text(5,2.8,' \leftarrow Y=b+m*X ','EdgeColor','red','LineWidth',2); text(9.5,3.5,['b = ',num2str(b)]) text(9.5,2.4,['m = ',num2str(m)])
162
xlabel('X'); ylabel('Y') axis([-6 15 -10 15]) figure %hold on %X variable dependiente: X = b1 + m1*Y b1=(Sx*Sy2-Sy*Sxy)/(N*Sy2-Sy^2) m1=(N*Sxy-Sx*Sy)/(N*Sy2-Sy^2) yy=-5:1:10; plot(yy,b1+m1*yy,'b',y,x,'k*-') text(-3.5,-7,' \leftarrow X=b1+m1*Y ','EdgeColor','black','LineWidth',2); text(1.5,-6.4,['b = ',num2str(b1)]) text(1.5,-7.4,['m = ',num2str(m1)]) xlabel('X'); ylabel('Y') axis([-6 14 -10 15])
% falta calcular el error % verificar en la bibliografa y adicionar....TAREA % Tarea para el estudiante: completar el programa con la parte de clculo de error en el ajuste
163
DIFERENCIACION
La derivada de un polinomio es
>> q = polyder(p) q= 4 6 6 4
Tambin
( 3X2 + 6X + 9 )*( x2 + 2X )
36
42
18 kk=polyder(c) kk = 12 36 42 18
164
EJERCICIOS Tambin; Uso del comando diff: syms x d=3*x^2+6*x+9 diff(d) ans = 6*x+6
OTRO:
165
ESTADISTICA DESCRIPTIVA MEDIA M=mean(A) retorna el valor medio de los elementos de arreglos de diferentes dimensiones Si A es un vector, mean(A) retorna el valor medio de A Si A es una matriz, mean(A) trata las columnas de A como vectores, retornando un vector fila con los valores de la media, por ejemplo:
>> a=[2 4 6 8 3 6 8] a= 2 4 6 8 3 6 8
>> A = [1 2 3; 3 3 6; 4 6 8; 4 7 7] A= 1 3 4 4 2 3 6 7 3 6 8 7
MEDIANA M=median(A) retorna el valor de la mediana de los elementos de un arreglo de cualquier dimensin. Si A es un vector, median(A) retorna el valor de la mediana de A. Si A es una matriz, median(A) trata a las columnas de A como un vector y retorna un vector fila con los valores de la mediana, por ejemplo:
>> A = [1 2 4 4; 3 4 6 6; 5 6 8 8; 5 6 8 8] A= 1 3 5 5 2 4 6 6 4 6 8 8 4 6 8 8
166
DESVIACION ESTNDAR Sintxis S=std(X) S=std(X,flag) S=std(X,flan,dim) Donde: n es el nmero de elementos en la muestra. Las dos formas de las ecuaciones difieren en n-1 y n en el divisor. s = std(X), donde X es un vector, retorna la desviacin estndar usando la ecuacin (1) de 2 arriba. Si X es una muestra al azar de datos desde una distribucin normal, s es la mejor respuesta estimada de su varianza. Si X es una matriz, std(X) retorna un vector fila conteniendo la desviacin estndar de los elementos de cada columna de X. Si X es un arreglo multidimensional, std(X) es la desviacin estndar de los th elementos a lo largo de la primera dimensin de X. S=std(X , flag) para flag=0, es la misma como std(X). Para flag=1, std(X,1) retorna la desviacin estndar usando la ecuacin (2), produciendo el segundo momento cerca a la mitad Para srd(X,flag,dim) calcula la desviacin estndar a lo largo de la dimensin de X especificado por el escalar dim. Por ejemplo: >> X=[1 5 9; 7 15 22] X= 1 5 7 15 9 22 % usando la ecuacin (1) de arriba 7.0711 9.1924
7.0711
5.0000
167
>> s0=std(X,1,3)
s0 =
dx ( x 2 + 3 x) = log x dt senx
dy ( x 3 + y 3 ) = sen( x) cos( y ) dx 2y dy ( y 3 + y 2 ) 3 y = + dz (1 y ) ln y
dy x + y 2 log( x + y ) =( ) + sen( x) cos( y ) dx x y e x y
x2 +1 f ( x) = x +1 g ( x) = x 3 * sen( x) x2 + 2
168
INTEGRACION NUMERICA
OBJETIVO
Usar la funcin int del MatLab y hacer uso del concepto de programacin para usar algunos mtodos de integracin numrica.
169
METODOS DE INTEGRACION NUMERICA % Integracion Simpson 3/8 forma abreviada function simpson38 clc; clear all; clear memory; clear command history; f=input(' Ingrese la funcin a integrar f(x) = ','s'); a=input(' Ingrese limite inferior : '); b=input(' Ingrese limite superior : '); n=input(' Ingrese numero de trapecios a considerar en la integracin : %f n=2*n; xmin=a; xmax=b; h=(b-a)/n; x=a:h:b; fx=eval(f);y=abs(fx); suma1=y(1)+y(n+1); suma2=3*sum(y(2:3:n-1)); suma3=3*sum(y(3:3:n)); suma4=2*sum(y(4:3:n-2)); suma= suma1+suma2+suma3+suma4; integral=(3/8)*h*suma; fprintf('El rea es : %10.9f\n ',integral); % grfico hold on xp= xmin:0.1:xmax; x=xp; yp=eval(f); x3=0:0.01:0.95; x=x3; yp1=eval(f); plot(x3,yp1,'r',xp,yp,'bo') x=a:0.1:b; y=eval(f); bar(x,y,'y'); axis([0 1 0 30]) % Metodo Abreviado syms x yy=int('3.2824e+5*x^10-1.5407e+6*x^9+3.1342e+6*x^8-3.6235e+6*x^7+2.6261e+6*x^61.2424e+6*x^5+3.8818e+5*x^4-79512*x^3+10465*x^2-855.51*x+44.084','x',0.25,0.75) %yy=3.2824e+5.*x.^10 - 1.5407e+6.*x.^9 + 3.1342e+6.*x.^8 - 3.6235e+6.*x.^7 + 2.6261e+6.*x.^6 1.2424e+6.*x.^5 + 3.8818e+5.*x.^4 - 79512.*x.^3 + 10465.*x.^2 - 855.51.*x + 44.084 y=(3.2824e+5.*x.^10 - 1.5407e+6.*x.^9 + 3.1342e+6.*x.^8 - 3.6235e+6.*x.^7 + 2.6261e+6.*x.^6 1.2424e+6.*x.^5 + 3.8818e+5.*x.^4 - 79512.*x.^3 + 10465.*x.^2 - 855.51.*x + 44.084);
');
170
RESPUESTA: Introducir la funcin sgte: y=(3.2824e+5.*x.^10 - 1.5407e+6.*x.^9 + 3.1342e+6.*x.^8 3.6235e+6.*x.^7 + 2.6261e+6.*x.^6 -1.2424e+6.*x.^5 + 3.8818e+5.*x.^4 - 79512.*x.^3 + 10465.*x.^2 - 855.51.*x + 44.084);
Se puede observar los resultados: mtodo Simpson y el mtodo abreviado El rea es 5.361569895 El rea es yy = 5.3937307293604290674603174603175 Mtodo Simpson Mtodo abreviado
OJO.Explique la diferencia
30
25
20
15
10
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
171
Tambin; Programas Aplicados a la integracin de Funciones % Este programa sirve para integrar por integrales numricas % por el mtodo de %SIMPSON una funcin dada % y para graficar dicha funcin. clc disp(' INTEGRACION NUMERICA CON EL METODO DE SIMPSON') f=input('ingrese la funcin a integrar f(x)= ','s');%FUNCION a=input('ingrese limite inferior: '); b=input('ingrese limite superior: '); n=input('ingrese numero de trapecios: '); n=2*n;xmin=a-1;xmax=b+1; h=(b-a)/n; x=a:h:b; fx=eval(f);y=abs(fx); suma1=y(1)+y(n+1); suma2=4*sum(y(2:2:n)); suma3=2*sum(y(3:2:n-1)); suma=suma1+suma2+suma3; integral=(h/3)*suma; fprintf('el rea es :%10.9f\n',integral); %grafica xp=xmin:0.2:xmax; x=xp; yp=eval(f); plot(xp,yp,'g'); hold on x=a:0.05:b; y=eval(f); bar(x,y,'r');
INTEGRACION NUMERICA CON EL METODO DE SIMPSON ingrese la funcin a integrar f(x)= sqrt(x.^2+x.*3+2) ingrese limite inferior: 2 ingrese limite superior: 9 ingrese nmero de trapecios: 10 12 el rea es :48.862387973
10
10
172
El primer problema puede solucionarse con el siguiente comando abreviado: syms x yy=int('3.2824e+5*x^10-1.5407e+6*x^9+3.1342e+6*x^8-3.6235e+6*x^7+2.6261e+6*x^61.2424e+6*x^5+3.8818e+5*x^4-79512*x^3+10465*x^2-855.51*x+44.084','x',0.25,0.75)
ESQUEMA DE LA SENTENCIA U ORDEN AL COMPUTADOR USANDO MATLAB syms x yy=int( funcin,x,limite inferior, limite mximo)
x2 +1 f ( x) = x +1 x 3 * sen( x) x2 + 2
g ( x) =
f ( x) =
Mtodo del trapecio
x2 +1 x +1
g ( x) =
x 3 * sen( x) x2 + 2
3.
Escribir 10 funciones para su integracin, con los dos mtodos (programa y orden abreviada)
173
BIBLIOGRAFIA
1. 2. 3. 4.
Vasquez Paragulla, Julio - Diseo de Programacin Edit. San Marcos 2000 Perez, Cesar - MATLAB Aplicaciones en las Ciencias y la Ingeniera Edit. Prentice Hall 2002 Nakamura, Shoichiro - MATLAB Anlisis Numrico - Edit. Prentice Hall 1997 Morales Marchena, Hern - Mtodos numricos y Visualizacin Grfica - Edit. Megabit 2005 Chapra, Steven - Mtodos Numricos para Ingenieros - Edit. Mc. Graw Hill 2004 Lecca, Eduardo - Mi Primer Matlab - Ed. Faffo-Lecca Editores - 2001
5. 6.
174