Story Transcript
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
35
2 La transformada de Laplace En esta sección vamos a mencionar algunos aspectos relativos a la utilización de la transformada de Laplace que pueden ser tratados con ayuda de Matlab. El primer detalle a tratar es la resolución de un sistema algebraico de ecuaciones lineales. Este problema se presenta frecuentemente en análisis de circuitos, especialmente en régimen sinusoidal permanente. No es estrictamente un problema de utilización de transformadas de Laplace, pero puede estar relacionado indirectamente con él. El planteamiento del problema consiste, básicamente, en reducirlo a uno del tipo para el que Matlaab fue desarrollado inicialmente (manipulación de matrices y vectores). El segundo apartado de la sección está dedicado a presentar distintas operaciones con polinomios, de las cuales una de las más importantes es la determinación de sus raíces, que es un problema básico en el cálculo de transformadas inversas de Laplace. Finalmente, en el tercer apartado se hace referencia a distintas posibilidades de representar el módulo y la fase de una función resultante de particularizar la transformada de Laplace de la función de transferencia de un circuito para el caso del régimen sinusoidal permanente.
2.1 Resolución de un sistema de ecuaciones La forma en la que generalmente se presenta un sistema de ecuaciones en un circuito en régimen sinusoidal permanente es [V] = [Z] x [I] donde [V] es un vector columna de n elementos cada uno de los cuales es la suma algebraica de las fuentes independientes de tensión de las n mallas; [I] es un vector columna de n elementos cada uno de los cuales es la corriente de una malla, y [Z] es una matriz cuadrada n x n cuyos elementos son las impedancias de las distintas mallas (los elementos del tipo Z ij con i = j) y las impedancias compartidas entre mallas diferentes (los elementos del tipo Z ij con i ≠ j). Ocasionalmente, [I] puede incluir tensiones introducidas por fuentes dependientes. En todo caso, [I] debe reunir todas las incógnitas del sistema, las cuales no pueden aparecer ni en [V], ni en [Z]. Algo similar ocurriría si quisiéramos resolver el circuito utilizando la técnica de análisis por nudos en lugar de la de análisis por mallas.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
36
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
La resolución de este sistema (es decir, la obtención de los elementos de [I]) mediante Matlab requiere los siguientes pasos. En primer lugar hay que definir los vectores y las matrices que aparecen en el problema. En el caso general, una matriz se define mediante la instrucción M =[M11 M12 ... M1n; M21 M22 ... M2n; ...; Mk1 Mk2 ... Mkn]
% Definición de una matriz
% Mij: elemento correspondiente a la fila i y a la columna j % k: número de filas de la matriz k x n % n: número de columnas de la matriz k x n
Como puede observarse, los elementos de una misma fila van separados por espacios en blanco; el final de una fila se separa del comienzo de la siguiente mediante un punto y coma. Con la misma instrucción podemos definir vectores fila y vectores columna. F =[F11 F12 ... F1n] C =[C11; C21; ...; Ck1]
% Definición de un vector fila de n elementos % Definición de un vector columna de k elementos
Los elementos de la matriz [Z] serán complejos en general. Un número complejo en Matlab se define mediante el comando real+imag*i
% Definición de un número complejo
% real: parte real del número complejo % imag: parte imaginaria del número complejo
Obsérvese que lo que es una j en teoría de análisis de circuitos se convierte aquí en una i (también vale la j). Obsérvese también que no se han dejado blancos ni antes, ni después del signo más; esto tiene por objeto evitar confusiones en la definición de vectores y matrices que incluyen números complejos (recuérdese que, como se indicó antes, los elementos de una fila van separados por blancos). Un ejemplo de la utilización de Matlab para obtener las corrientes de malla en un circuito en régimen sinusoidal permanente es el siguiente. Supóngase que el vector de tensiones es [V (V)] = [1 - j2; -5; 15 - j2] y que la matriz de impedancias es [Z (Ω)] = ([3, -2 + j2, 0; -1, 4 + j2, 1 - j3; 0, 1 - j2, 8 + j5]. El programa que obtiene las corrientes es
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
%%%%% RESOLUCIÓN DE UN SISTEMA DE ECUACIONES clear all;
37
%%%%%
% Elimina variables utilizadas en otras rutinas
V = [1-2*i; -5; 15-2*i] Z = [3 -2+2*i 0; -1 4+2*i 1-3*i; 0 1-2*i 8+5*i] I = Z\V clear all;
% Elimina las variables utilizadas en esta rutina
Lo que se obtiene en la ventana de comandos de Matlab es V = 1.0000- 2.0000i -5.0000 15.0000- 2.0000i Z = 3.0000 -1.0000 0
-2.0000+ 2.0000i 4.0000+ 2.0000i 1.0000- 2.0000i
0 1.0000- 3.0000i 8.0000+ 5.0000i
I = 1.0000 0.0000+ 1.0000i 1.0000- 1.0000i
Es decir, las corrientes buscadas son [I (A)] = (1; j; 1 - j). Obsérvese que los resultados anteriores aparecen en la ventana de comandos porque las distintas instrucciones no han sido terminadas con puntos y comas. De haber incluido esta terminación en una o más sentencias, el programa habría conservado la información internamente, pero no la habría reflejado en la ventana de comandos. Obsérvese que, contra lo que podría parecer lógico ([I] = [V] / [Z]), la instrucción de división tiene en este caso la forma [I] = [Z] \ [V]. En realidad, esta forma es equivalente a realizar la operación [I] = inv[Z] x [V], que ya se corresponde con lo que parecería razonable. Esta sentencia viene impuesta por la forma de operar de Matlab y se denomina división a la izquierda, como contraposición a la habitual división a la derecha. De hecho, podríamos haber definido el sistema de ecuaciones como [V] = [I] x [Z] En ese caso el programa para obtener las corrientes sería
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
38
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
%%%%% RESOLUCIÓN DE UN SISTEMA DE ECUACIONES clear all;
%%%%%
% Elimina variables utilizadas en otras rutinas
V = [1-2*i -5 15-2*i] Z = [3 -2+2*i 0; -1 4+2*i 1-3*i; 0 1-2*i 8+5*i] I = V/Z.' clear all;
% Elimina las variables utilizadas en esta rutina
y los resultados en la ventana de comandos de Matlab serían V = 1.0000- 2.0000i
-5.0000
15.0000- 2.0000i
Z = 3.0000 -1.0000 0
-2.0000+ 2.0000i 4.0000+ 2.0000i 1.0000- 2.0000i
0 1.0000- 3.0000i 8.0000+ 5.0000i
0.0000+ 1.0000i
1.0000- 1.0000i
I = 1.0000
Obsérvese que ahora el vector [V] se ha definido como fila y no como columna, lo cual hace que el vector [I] también sea presentado en la misma forma. Además, en lugar de utilizar la matriz [Z], se utiliza su transpuesta, definida mediante la instrucción. Z.’
% Obtiene la matriz transpuesta de la matriz Z
Téngase presente que la obtención de la transpuesta de una matriz de números complejos requiere un procedimiento matemático de cierta complicación.
2.2 Transformada inversa de Laplace Raíces de un polinomio En el cálculo de transformadas inversas de Laplace un problema frecuente es la determinación de las raíces de un polinomio en s. En general, un polinomio de este tipo es de la forma
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
39
P(s) = aksk + ak-1sk-1 + ... + a0 en la que los coeficientes de los términos exponenciales son números reales. Sin embargo, esta forma no es la más adecuada para utilizar el procedimiento convencional de obtener la transformada inversa de Laplace. Es preferible expresar el polinomio como el producto de una serie de términos en los que aparecen las raíces de aquél. Es decir, P(s) = aksk + ak-1sk-1 + ... + a0 = (s + s1)(s + s2) ... (s + sj)q ... (s + sm) A cualquiera de los términos -sk se le denomina raíz del polinomio. sj (y cualquier otra que figure en un término exponencial) se le denomina raíz múltiple de orden q (el exponente del término exponencial). El orden del polinomio es k (es decir, el exponente de la mayor potencia de s que figure en el polinomio). k es igual al número total de raíces del polinomio; en el cómputo de este número cada raíz múltiple cuenta tantas veces como su propio orden. Las raíces de un polinomio pueden ser reales o complejas; si un número complejo es raíz de un polinomio representativo de un circuito, su complejo conjugado también lo es. En el problema de determinar las raíces de un polinomio pueden presentarse varios casos muy sencillos, como son los que se mencionan seguidamente. Polinomio P(s)
Raíces
s+r
-r
(s + r)n
-r (raíz múltiple de orden n)
as2 + bs + c (polinomio cuadrático)
- b ± b 2 - 4ac r 1, 2 = 2a
Al margen de estos casos, no existe ningún procedimiento analítico sencillo para obtener las raíces de polinomios de orden superior al segundo. Es necesario recurrir a procedimientos numéricos para realizar dicha tarea. En Matlab esta operación se realiza con la instrucción roots (pol)
% Obtiene las raíces del polinomio pol
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
40
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
Obsérvese que la aplicación de esta instrucción requiere definir previamente el polinomio del que se quieren obtener las raíces. El polinomio se define como un vector fila cuyos elementos son los coeficientes ordenados de los términos exponenciales, empezando por el correspondiente a la potencia mayor. Si faltan algunos de los términos exponenciales, esta circunstancia es tenida en cuenta asignando el valor nulo al coeficiente correspondiente. Recuérdese que el término independiente del polinomio, si lo hay, corresponde a una potencia cuyo exponente es nulo; por tanto, dicho término debe ser incluido en la definición del polinomio. Así, la instrucción % Definición de un polinomio
pol = [a(k) a(k-1) ... a(0)]
corresponde al polinomio P(s) = aksk + ak-1sk-1 + ... + a0. Si el polinomio no va a ser utilizado con posterioridad a la obtención de sus raíces, pueden combinarse las dos instrucciones que acabamos de mencionar en una sola y utilizar la sentencia % Obtiene las raíces de un polinomio
roots ([a(k) a(k-1) ... a(0)])
A título de ejemplo, considérese el polinomio P(s) = s4 + 6s3 + 14s2 + 16s + 8. La utilización de la rutina %%%%% RAÍCES DE UN POLINOMIO clear all;
%%%%% % Elimina variables utilizadas en otras rutinas
roots ([1 6 14 16 8]) clear all;
% Elimina las variables utilizadas en esta rutina
proporciona los siguientes resultados en la ventana de comandos de Matlab: -2.0000+ -2.0000-1.0000+ -1.0000-
0.0000i 0.0000i 1.0000i 1.0000i
Es decir, el polinomio tiene una raíz real doble (-2 + j0) y dos raíces complejas (-1 + j, -1 - j), de las cuales una es la conjugada de la otra.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
41
Valor de un polinomio Otra operación de cierto interés que puede realizarse con Matlab es la determinación del valor particular que toma un polinomio cuando la variable independiente tiene un valor dado. Ello puede realizarse utilizando la instrucción polyval (pol, k)
% Obtiene el valor del polinomio pol cuando la variable vale k
A título de ejemplo, la rutina siguiente permite hallar el valor que toma el polinomio P(s) definido anteriormente cuando s es igual a -2. %%%%% VALOR DE UN POLINOMIO clear all;
%%%%% % Elimina variables utilizadas en otras rutinas
P = [1 6 14 16 8]; polyval (P, -2) clear all;
% Elimina las variables utilizadas en esta rutina
El resultado es ans = 0
Es un resultado lógico, por cuanto hemos pretendido calcular el valor del polinomio cuando la variable es igual a una de sus raíces y, por definición, éstas son los valores que anulan el polinomio. Producto de polinomios Otra operación de interés es el producto de dos polinomios. Se realiza utilizando la instrucción conv (pol1, pol2)
% Efectúa el producto de los polinomios pol1 y pol2
Así, por ejemplo, si multiplicamos los polinomios P1(s) = s2 + 2s + 2 y P2(s) = s2 + 4s + 4 mediante la rutina
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
42
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
%%%%% PRODUCTO DE POLINOMIOS clear all;
%%%%% % Elimina variables utilizadas en otras rutinas
P1 = [1 2 2]; P2 = [1 4 4]; conv (P1, P2) clear all;
% Elimina las variables utilizadas en esta rutina
obtenemos el siguiente resultado: ans = 1
6
14
16
8
Es decir, el polinomio resultante del producto está definido por sus coeficientes, al igual que cualquier otro polinomio. Obsérvese que el resultado obtenido en este ejemplo coincide con el polinomio P(s) definido anteriormente. Obsérvese también que la instrucción que se utiliza para el producto (conv) significa convolución. Sin embargo, aunque esta instrucción tiene un papel relevante en la obtención de la convolución de señales discretas, no guarda una relación inmediata con la convolución de señales continuas, tal y como se indicó en la sección precedente de este manual. División de polinomios Se trata de la operación inversa de la que se acaba de mencionar. Se realiza utilizando la instrucción % Divide el polinomio pol1 entre el polinomio pol2
deconv (pol1, pol2)
Así, por ejemplo, si dividimos el polinomio P1(s) = s3 + 3s2 + 2s + 2 entre el polinomio P2(s) = s2 + 4s + 4 mediante la rutina %%%%% DIVISIÓN DE POLINOMIOS clear all;
%%%%% % Elimina variables utilizadas en otras rutinas
P1 = [1 3 2 2]; P2 = [1 4 4]; [cociente, resto] = deconv (P1, P2) clear all;
% Elimina las variables utilizadas en esta rutina
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
43
obtenemos los siguientes valores para el cociente y el resto de la operación: cociente = 1
-1
resto = 0
0
2
6
Ambas variables definen sendos polinomios. El segundo es el numerador de una fracción cuyo denominador coincide con el del cociente original. En otras palabras, el resultado obtenido nos indica que P 1(s) = s - 1 + 2s + 6 P 2(s) s 2 + 4s + 4 Esta instrucción tiene un interés especial cuando una transformada de Laplace está definida por una fracción impropia (es decir, el orden del numerador es mayor o igual que el del denominador). Recuérdese que la aplicación del procedimiento convencional para obtener la transformada inversa requiere efectuar el cociente de forma que el resto sea una fracción propia (el orden del numerador es menor que el del denominador). Descomposición de un cociente en suma de fracciones simples Podemos utilizar todo lo que acabamos de exponer para obtener transformadas inversas de Laplace. Después del cociente de polinomios, quedaría por obtener las raíces del denominador del resto (roots) y formular tantas fracciones como raíces tenga dicho denominador (recuérdese que una raíz múltiple cuenta tantas veces como su orden). Sin embargo, hay un procedimiento directo para descomponer un cociente en una suma de fracciones simples cuyos denominadores corresponden a las raíces del denominador original. Ese procedimiento se ejecuta en Matlab mediante la instrucción [num, raíces, cociente] = residue (pol1, pol2)
% Expande en suma de fracciones simples % el cociente entre los polinomios pol1 y pol2
% num: vector columna en la que aparecen los numeradores de las distintas fracciones % raíces: vector columna de las raíces del denominador ordenadas igual que los numeradores % cociente: cociente de los dos polinomios (sólo si la fracción original es impropia)
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
44
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
No es preciso invocar explícitamente estos tres parámetros en la rutina dedicada a la descomposición de un cociente; si no se pone un punto y coma tras la instrucción residue, tales parámetros aparecen directamente en la ventana de comandos. En cualquier caso, son guardados internamente por el programa y puede hacerse uso de ellos en operaciones posteriores. Al ejecutar esta instrucción una raíz múltiple no aparece señalada explícitamente como tal. La indicación de semejante característica la proporciona el hecho de que la raíz múltiple aparece repetida tantas veces como su orden. A cada una de ellas le corresponde el numerador que ocupa la misma posición que la raíz (o su repetición) en ambos vectores. Recuérdese que, al formular las fracciones simples originadas por una raíz múltiple, el denominador de la primera de éstas aparece elevado a la unidad; el segundo, al cuadrado, y así sucesivamente hasta alcanzar el orden de multiplicidad de la rutina. A título de ejemplo, considérese la rutina %%%%% EXPANSIÓN DE UN COCIENTE EN FRACCIONES SIMPLES clear all;
%%%%%
% Elimina variables utilizadas en otras rutinas
P1 = [1 3 2 2]; P2 = [1 4 4]; [num, raíces, cociente] = residue (P1, P2) clear all;
% Elimina las variables utilizadas en esta rutina
El resultado que aparece en la ventana de comandos de Matlab es num = 2 2 raíces = -2 -2 cociente = 1
-1
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
45
Es decir, se ha realizado la siguiente operación P 1(s) 2 =s-1+ 2 + P 2(s) s + 2 (s + 2) 2 siendo P1(s) y P2(s) los polinomios definidos anteriormente. La operación inversa de la que acabamos de describir permite recomponer el cociente de dos polinomios a partir de los valores del cociente, las raíces y los numeradores. Así, por ejemplo, la rutina %%%%% RECOMPOSICIÓN DE UN COCIENTE clear all;
%%%%%
% Elimina variables utilizadas en otras rutinas
num = [2; 2]; raíces = [-2; -2]; cociente = [1; -1]; [P1, P2] = residue (num, raíces, cociente) clear all;
% Elimina las variables utilizadas en esta rutina
devuelve el resultado P1 = 1
3
2
1
4
4
2
P2 =
que, como puede comprobarse, corresponde a los polinomios a los que nos estamos refiriendo. Evidentemente, todo lo que acabamos de exponer es aplicable igualmente en el caso de que la expansión conduzca a fracciones en las que aparecen raíces complejas. Así, por ejemplo, si se trata del cociente P 1(s) 0.5s = P 2(s) s 2 + s + 0.5 la ejecución de la correspondiente instrucción residue conduce a
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
46
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
num = 0.2500+ 0.2500i 0.2500- 0.2500i raíces = -0.5000+ 0.5000i -0.5000- 0.5000i
En términos matemáticos, este resultado puede expresarse como P 1(s) 0.5s = = 0.25 + j0.25 + 0.25 - j0.25 = s +Ks + K* P 2(s) s 2 + s + 0.5 s + 0.5 - j0.5 s + 0.5 + j0.5 1 s + s *1 con lo que la transformada inversa de Laplace de este cociente es directamente un único término en el que aparecen el módulo y la fase de K. Éstos pueden obtenerse por cálculo manual, o utilizando instrucciones de Matlab a las que se hace referencia más adelante.
2.3 Respuesta en régimen sinusoidal permanente La respuesta de un circuito en régimen sinusoidal permanente puede ser determinada a partir del conocimiento de su función de transferencia expresada en función de la variable s; para ello basta sustituir dicha variable por jω. Esto supone una diferencia importante con relación a lo expuesto en el apartado anterior. En él, s era una variable genérica, mientras que ahora pasa a ser una variable continua. Dado que es imposible reproducir exactamente en Matlab una variable continua, la transformación implica el paso de la variable genérica s a la variable discreta jω. En otras palabras, en lo que sigue s va a comportarse como un vector que almacena los distintos valores que puede tomar ω. Por otro lado, una vez que se ha efectuado la transformación indicada, la función de transferencia H(s) pasa a ser una función variable con la frecuencia (H(jω)). La observación de la variación de esta función proporciona información muy útil acerca del comportamiento del circuito; por ejemplo, permite determinar su tipo de respuesta cuando el circuito es tratado como un filtro. Dado que la nueva función H(jω) es en general compleja, dicha variación queda determinada a partir de las variaciones de su módulo y su fase con ω. En lo que sigue abordaremos precisamente estos aspectos: el tratamiento de s como vector, la definición del rango de variación de este vector, y la obtención de módulos y fases de funciones complejas.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
47
En definitiva, se trata de un problema cuyo planteamiento puede resumirse en las siguientes equivalencias matemáticas H(s) =
P 1(s) P (jω) P 1(jω) ∠θ1(ω) ⇒ H(jω) = 1 = = H(jω) ∠θ(ω) P 2(s) P 2(jω) P 2(jω) ∠θ2(ω) H(jω) =
P 1(jω) P 2(jω)
θ(ω) = θ 1(ω) - θ 2(ω)
En otras palabras, hay distintas posibilidades de afrontarlo. A continuación consideraremos algunas de ellas. Rango de frecuencias En los problemas abordados en este apartado lo que se persigue es representar la variación de una magnitud frente a la frecuencia angular (ω). Por consiguiente, lo primero que hay que hacer en cualquiera de dichos problemas es establecer el rango de valores que debe tomar ω. Este aspecto es, en principio, formalmente idéntico al de establecer una base de tiempos que hemos considerado en la sección anterior. Por tanto, podemos definir el rango de frecuencias mediante cualquiera de los dos procedimientos que se mencionaron en dicha sección, procedimientos que habían sido denominados lineales; basta sustituir la variable t por la variable ω para definir el rango de frecuencias. Ahora bien, al abordar problemas de respuesta en frecuencia de circuitos funcionando en régimen sinusoidal permanente el rango de variación de la frecuencia suele cubrir varias décadas, con lo que establecer como lineal tal variación no suele ser lo más conveniente. Una instrucción más adecuada es instrucción logspace (d1, d2, puntos)
% Define un rango de valores de variación logarítmica
% d1: parámetro que define el inicio del rango (10d1) % d2: parámetro que define el final del rango (10d2) % puntos: número de puntos totales entre el inicio y el final del rango (la separación entre dos puntos consecutivos es igual en coordenadas logarítmicas)
Así, por ejemplo, la instrucción w = logspace (1, 3, 1000) genera un vector de 1000 puntos comprendidos entre 10 y 103.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
48
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
Representación de números complejos Hay cuatro instrucciones básicas para caracterizar un número complejo. En el caso particular que estamos tratando, dicho número tendrá la forma genérica z(ω) = a(ω) + jb(ω). Las instrucciones a las que nos referimos son: % Obtiene la parte real de z % Obtiene la parte imaginaria de z % Obtiene el valor absoluto de z % Obtiene la fase de z (en rad)
real (z) imag (z) abs (z) angle (z)
A título de ejemplo, consideremos la siguiente rutina. %%%%% NÚMEROS COMPLEJOS %%%%% clear all;
% Elimina variables utilizadas en otras rutinas
w = logspace (-1, 3, 1000);
% Rango de frecuencias
a = 1; b = (w - 1./w); z = a + i*b;
% Número complejo
subplot (2, 2, 1); % Parte real semilogx (w, real(z), 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Parte real', 'FontName', 'Times', 'Fontsize', 14); % Ord axis ([10^(-1) 10^3 (1/3)*min(real(z)) (3/2)*max(real(z))]); % Área de dibujo grid on; % Malla subplot (2, 2, 2); % Parte imaginaria semilogx (w, imag(z), 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Parte imaginaria', 'FontName', 'Times', 'Fontsize', 14); % Ord axis ([10^(-1) 10^3 (1/3)*min(imag(z)) (3/2)*max(imag(z))]); % Área de dibujo grid on; % Malla subplot (2, 2, 3); % Módulo semilogx (w, abs(z), 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Módulo', 'FontName', 'Times', 'Fontsize', 14); % Ord axis ([10^(-1) 10^3 (1/3)*min(abs(z)) (3/2)*max(abs(z))]); % Área de dibujo grid on; % Malla subplot (2, 2, 4); % Fase semilogx (w, angle(z), 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Fase (rad)', 'FontName', 'Times', 'Fontsize', 14); % Ord axis ([10^(-1) 10^3 -3*pi/4 3*pi/4]); % Área de dibujo grid on; % Malla clear all;
% Elimina las variables utilizadas en esta rutina
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
49
1.5
1 b=ω-ω
Parte imaginaria
a=1
Parte real
1400
1
0.5
z = a + jb
1200 1000 800 600 400 200
0
10
0
2
10
Frecuencia angular (rad/s)
Datos:
0
10
2
10
Frecuencia angular (rad/s)
1400
2
Módulo
a
Fase (rad)
1200 1000 800 600 400 200
1 0 -1 -2
0
10
2
10
Frecuencia angular (rad/s)
0
10
2
10
Frecuencia angular (rad/s)
Como puede observarse, la rutina representa la variación con ω de la parte real, la parte imaginaria, el módulo y la fase del número complejo que aparece definido a la izquierda de la figura. Algunas particularidades de esta rutina se mencionan seguidamente. En primer lugar, cabe destacar que la habitual instrucción plot ha sido remplazada. Como se ha visto en ejemplos anteriores, dicha instrucción permitía obtener una representación lineal de ambos ejes coordenados. Teniendo en cuenta que la sensibilidad de un circuito a la variación de la frecuencia suele ser exponencial, si utilizamos dicha instrucción obtendremos una representación muy pobre, con la consiguiente pérdida de información. Se sugiere al lector que observe este efecto analizando el comportamiento de la rutina cuando se utiliza la instrucción plot. A cambio, en este caso (y en prácticamente todos los que implican variaciones con la frecuencia) se ha utilizado la instrucción semilogx
% El eje de abscisas se representa en coordenadas semilogarítmicas
Esta instrucción sencillamente sustituye a la de plot, con lo que puede ser dotada de los mismos atributos que aquélla. Hay otras instrucciones similares, que permiten representar semilogarítmicamente el eje de ordenadas o bien utilizar coordenadas logarítmicas en ambos ejes, pero carecen de interés desde la perspectiva de este manual.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
50
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
En segundo lugar, obsérvese que la fase del número complejo se calcula y se representa en radianes. Ésa es la representación de Matlab por defecto, pero existen otras posibilidades, que serán expuestas más adelante. También es importante destacar que, en la definición del número complejo, se utiliza el operador ./, que indica que el numerador debe ser dividido por todos y cada uno de los elementos del vector que constituye el denominador. De no utilizarse este operador, se produciría una indefinición (con la consiguiente paralización del programa), ya que Matlab ignoraría qué valor concreto del vector denominador debe ser utilizado en el cociente. Finalmente, en la rutina se utilizan de forma accesoria las instrucciones % Obtiene el valor máximo del vector x % Obtiene el valor mínimo del vector x
max (x) min (x)
Números complejos expresados como fracciones Aunque útil para ilustrar algunos aspectos de la representación de números complejos, el ejemplo que acabamos de considerar es raramente representativo de lo que realmente interesa en el análisis de la respuesta en frecuencia de circuitos funcionando en régimen sinusoidal permanente. Lo más habitual es que la función de transferencia esté caracterizada por el cociente de dos polinomios de s; al sustituir esta variable por jω queda un número complejo que varía con la frecuencia. Es a este caso al que nos referiremos en este epígrafe. La siguiente rutina permite obtener el módulo y la fase de una función de transferencia representada por una fracción propia (el orden del numerador es menor que el del denominador), Hs, en la que el numerador y el denominador son dos polinomios de s. Como puede observarse, tanto el módulo y la fase se obtienen aplicando directamente las instrucciones abs y angle a la fracción en lugar de hacerlo al numerador y al denominador por separado. %%%%% NÚMEROS COMPLEJOS COMO FRACCIONES clear all;
%%%%%
% Elimina variables utilizadas en otras rutinas
inicial = -2; final = - inicial; puntos = 1000; % Rango de frecuencias w = logspace (inicial, final, puntos); s = i*w; Hs = 0.5*s./(s.^2 + 2*s + 2); moduloh = abs(Hs); fase = (180/pi)*angle(Hs);
% Número complejo como fracción
subplot (2, 1, 1); % Módulo semilogx (w, moduloh, 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Módulo', 'FontName', 'Times', 'Fontsize', 14); % Ord
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
51
axis ([w(1) w(puntos) (1/3)*min(moduloh) (3/2)*max(moduloh)]); % Área de dibujo grid on; % Malla title ('Respuesta en frecuencia', 'FontName', 'Times', 'Fontsize', 24); % Título subplot (2, 1, 2); % Fase semilogx (w, fase, 'b', 'LineWidth', 2); % Curva en azul de grosor 2 xlabel ('Frecuencia angular (rad/s)', 'FontName', 'Times', 'Fontsize', 14); % Abs ylabel ('Fase (º)', 'FontName', 'Times', 'Fontsize', 14); % Ord axis ([w(1) w(puntos) -180 180]); % Área de dibujo grid on; % Malla modmax = max(moduloh) m = find (moduloh == max(moduloh)); wcentral = w(m) fasecentral = fase(m)
% Frecuencia central
positivo = find (fase >= 0); n = length (positivo); wfasenula = w(n) modfasnul = moduloh(n) fasenula = fase(n)
% Frecuencia de fase nula
clear all;
% Elimina las variables utilizadas en esta rutina
Respuesta en frecuencia 0.35
H(s) =
0.5s s 2 + 2s + 2
Módulo
0.3 0.25 0.2 0.15 0.1 0.05 -2
10
-1
10
0
10
1
10
2
10
Frecuencia angular (rad/s)
Datos: 150
inicial Fase (º)
100 50 0 -50 -100 -150 -2
10
-1
10
0
10
1
10
2
10
Frecuencia angular (rad/s)
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
52
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
Obsérvese, en la definición de la fracción, la utilización de los operadores ./ y .ˆ, ya que las operaciones aludidas son efectuadas entre vectores. Nótese también que la fase de la función de transferencia ha sido expresada en º. La rutina se completa con un conjunto de instrucciones finales destinadas a proporcionar ciertos datos que no pueden obtenerse (al menos con precisión) de la observación de las curvas. El primer grupo de tales instrucciones está destinado a obtener la frecuencia central (es decir, aquélla para la que se obtiene el máximo del módulo de la función de transferencia), y los valores del módulo y la fase correspondientes a dicha frecuencia. El segundo grupo tiene por objeto determinar la frecuencia para la que se anula la fase de la función de transferencia, y los consiguientes valores de módulo y fase. Los resultados obtenidos se indican más abajo. Es conveniente que el lector repare en la forma en la que se ha establecido el cálculo de la segunda frecuencia. El máximo del módulo es un valor concreto del vector que contiene los valores de aquél; por tanto, es posible hallar directamente la posición que ocupa tal valor en el vector. En cambio, la condición de que la fase sea nula es algo impuesto externamente. En otras palabras, el vector que almacena los valores de la fase puede no contener el valor 0 exacto. Por consiguiente, es preciso determinar de forma aproximada la posición dentro del vector que ocupa el valor de la fase más próximo a 0. modmax = 0.2500 wcentral = 1.4130 fasecentral = 0.0681 wfasenula = 1.4130 modfasnul = 0.2500 fasenula = 0.0681
Obsérvese que en este caso las dos frecuencias son idénticas en apariencia. Ello puede ser así realmente, o bien el resultado de las aproximaciones y redondeos internos del programa en el caso particular de la función de transferencia considerada. Sin embargo, debe recordarse que, en general, tales frecuencias pueden diferir significativamente una de otra. Consideraciones sobre la representación de la fase En ocasiones la representación de la fase presenta una discontinuidad (por ejemplo, pasando de -180 º a 180 º) a una frecuencia dada. Esta circunstancia puede deberse a que el
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
53
circuito se comporta de forma que realmente se produce dicha discontinuidad, o bien puede ser originada por la manipulación matemática interna de Matlab. A este respecto, recuérdese que, a fin de cuentas, los ángulos α y 180 º + α tienen la misma tangente y que es la tangente la función que suele utilizarse en la determinación de ángulos en ordenadores y calculadoras. Para tener la seguridad de que el comportamiento representado es el que corresponde realmente al circuito puede utilizarse la instrucción unwrap (fase)
% Elimina discontinuidades debidas al cálculo en la representación de la fase
La instrucción unwrap utiliza y devuelve ángulos expresados en radianes. Esta circunstancia debe ser tenida en cuenta cuando se pretende obtener representaciones en grados sexagesimales. A título de ejemplo de la utilización de esta instrucción considérese el siguiente fragmento de rutina. s = i*w; Hs = s./(s.^4 + 2*s.^3 + 4*s.^2 + 2*s + 1); fas = angle(Hs); fase1 = (180/pi)*unwrap(fas); fase2 = (180/pi)*angle(Hs); semilogx (w, fase1, 'b', 'LineWidth', 2); % Sin discontinuidades de cálculo (azul) hold on; semilogx (w, fase2, 'r', 'LineWidth', 2); % Con posibles discontinuidades (rojo)
Con ayuda de este fragmento pueden obtenerse curvas como las representadas en la figura siguiente. 200
150
100
50
0
Fase (º)
s H(s) = 4 3 s + 2s + 4s 2 + 2s + 1
-50
-100
-150
-200
-250
-300 -1 10
0
10
1
10
Frecuencia angular (rad/s)
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
54
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
Por otro lado, la figura siguiente muestra un caso en el que la discontinuidad es propia del circuito y, por tanto, no puede ser eliminada con ayuda de la instrucción unwrap. 100
80
60
40
+ 0.5 s 2 + 0.5s + 1
20
Fase (º)
H(s) =
0.5s 2
0
-20
-40
-60
-80
-100 -1 10
0
10
1
10
Frecuencia angular (rad/s)
Esta segunda figura ha sido obtenida con ayuda del mismo fragmento de rutina considerado anteriormente, aunque cambiando la instrucción que define la función de transferencia. Determinación simultánea de módulo y fase En los epígrafes anteriores hemos presentado rutinas que tratan separadamente la representación de módulos y fases de funciones de transferencia de circuitos funcionando en régimen sinusoidal permanente. Sin embargo, existe una instrucción que permite efectuar dicha representación simultáneamente. Se trata de la instrucción freqs (pol1, pol2, w)
% Obtiene la variación con w de la fracción en la que pol1 y pol2 % son los polinomios del numerador y el denominador, respectivamente
La rutina y la figura siguientes recogen un ejemplo de la utilización de esta instrucción. Como puede observarse, su uso implica una reducción del número de instrucciones necesarias (antes eran precisas abs y angle), así como una disminución del riesgo de equivocarse en la utilización de operadores precedidos de punto. Incluso, como se hace en la rutina siguiente, no es necesario definir explícitamente los polinomios que constituyen el numerador y el denominador de la fracción a representar.
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
55
A cambio, la utilización de esta instrucción comporta algunos inconvenientes. Por ejemplo, no es posible editar las figuras (la edición la hace automáticamente Matlab). Además, el eje de ordenadas de la representación del módulo es de tipo logarítmico, lo cual da lugar a una forma de presentación muy poco habitual en este tipo de problemas. En este sentido, puede observarse que la función considerada en este ejemplo ya ha sido representada en uno anterior y comparar ambos tipos de representación. %%%%% MÓDULO Y FASE DE UNA FRACCIÓN
%%%%%
clear all;
% Elimina variables utilizadas en otras rutinas
w = logspace (-2, 3, 1000); s = i*w; freqs ([0.5 0], [1 2 2], w);
% Rango de frecuencias % Módulo y fase del cociente de dos polinomios
clear all;
% Elimina las variables utilizadas en esta rutina 0
10
-1
Magnitude
10
-2
10
-3
10
-4
10
0.5s H(s) = s 2 + 2s + 2
-2
10
-1
10
0
1
0
1
10 10 Frequency (radians)
2
10
3
10
Phase (degrees)
100
50
0
-50
-100 -2 10
-1
10
10 10 Frequency (radians)
2
10
3
10
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO
56
Prácticas de circuitos como sistemas lineales. Ejercicios sencillos con Matlab
2.4 Ejercicios propuestos Ejercicio 6
- V1 +
El circuito de la figura, en cuya representación se ha utilizado notación fasorial, funciona en régimen sinusoidal permanente a una frecuencia angular dada. Se desea obtener los valores de los fasores de las corrientes de malla y de la tensión en la fuente dependiente.
R1
M VG
C1 I1
aV1 L1
C2 I2
R3 C3 L2
L3
I3
VG = 0.5 - j1.5 V, a = 2, ω = 1 krad/s R1 = 1 Ω, R3 = 1 Ω L 1 = 1 mH, L2 = 2 mH, L3 = 2.5 mH, M = 1 mH C1 = 1 mF, C2 = 0.5 mF, C3 = 0.4 mF
Ejercicio 7 Obtener la transformada de Laplace (expresándola como cociente de dos polinomios) de la función f(t) = 3δ(t) - 2te-5t + 4e-2tcos(8t) Ejercicio 8 Obtener la transformada inversa de Laplace de la función 6 5 4 3 2 L(s) = 2s + 2.3s - 6.39s - 6.755s - 1.9925s - 0.0275s + 0.005 s 6 + 3.1s 5 + 2.55s 4 + 0.725s 3 + 0.05s 2
Ejercicio 9 Representar la variación con la frecuencia del módulo y la fase de la función de transferencia de un circuito funcionando en régimen sinusoidal permanente. H(s) =
2×10 6s 2 + 10 18 (s + 10 6) (s 2 + 2×10 6s + 10 12)
Dpto. Teoría de la Señal y Comunicaciones. Escuela Técnica Superior de Ing. Telecomunicación. UNIVERSIDAD DE VIGO