Utilizando el modelamiento matemático en función del tiempo, Laplace y transformada Z de

un horno. Simular los resultados en Matlab, Simulink y Control System.

0.1998
𝐺(𝑠) =
𝑠 + 0.03772
1. MATLAB
Modelamiento en función del tiempo
%****************** MODELACION DE UN HORNO *******************************
%***************** en función del tiempo *******************************
clear all; close all; clc;

%Pametros de simmulacion

T0 = 0.1; %Tiempo de muestreo

tsim = 3000; %Tiempo de simulacion
t = [0:T0:tsim]; %Vector tiempo
%Parametros de simulacion del entorno, consideramos las siguientes constantes

%Calculo de C1=m*c
m=13; %masa de la caga(Kg)
c=0.75; %Calor especifido del producto
C1=m*c; %Capacitancia de la carga

%Calculo de C2=cp*p*v
cp=0.28; %Calor especifico del aire a 172°C
p=0.748; %Densidad kg/m3
v=7.95; %Volumen interno del horno m3
C2=cp*v*p; %Capacitancia de la fuente de calor

%Calculo de C3
C3=70.84; %Calor especifico de la lamina de las paredes
internas, masa=554.38kg
rt=0.1;
V=438; %Voltaje de entrada

%Condiciones iniciales
temp(1)=25; %Temperatura ambiente

%Modelo de cambio de la potencia de entrada que define la temperatura final.

Q = [242*ones(1,2000) 242*ones(1,2000) 242*ones(1,2000)...
242*ones(1,2000) 242*ones(1,2000) 242*ones(1,2000)...
242*ones(1,2000) 242*ones(1,2000) 242*ones(1,length(t)-16000)];

%************************* MODELACION DE UN HORNO** ************************

for k=1:length(t)
% 2) Modelo no Lineal

temp(k)=(V/(C1+C2+C3))*Q(k)*(1-exp(-k/((C1*rt)*(C2+C3))/(C1+C2+C3)));
end
%******************** ANIMACION DE LOS ESTADOS ***********************
%**************************************************************************
figure ('Name','Modelado del horno en funcion del tiempo')
graf1 = plot(t(1),temp(1),'r','LineWidth',2); hold on;
for i=1:50:length(t)

delete(graf1)
xlabel ('Tiempo');
ylabel ('Temperatura (°C)');
graf1 = plot(t(1:i),temp(1:i),'r','LineWidth',2); hold on;
drawnow;
end
 Ilustración 1 Gráfica del modelamiento de un horno en función del tiempo

Modelamiento en Laplace

%% ************************************************************************
%************************* MODELACION DE UN HORNO** ************************
%**************************************************************************
for k=1:length(s)
% 2) Modelo no Lineal

temp(k)=Q(k)*(V/(C1+C2+C3))/((((C1*rt)*(C2+C3))/(C1+C2+C3))*k+1);
end
%% ************************************************************************
%******************** ANIMACION DE LOS ESTADOS ***********************
%**************************************************************************
figure ('Name','Modelado de la horno en funcion de la frecuencia')
graf1 = plot(s(1),temp(1),'r','LineWidth',2); hold on;
for i=1:50:length(s)

delete(graf1)
xlabel ('Frecuencia');
ylabel ('Temperatura (°C)');
graf1 = plot(s(1:i),temp(1:i),'r','LineWidth',2); hold on;
drawnow;
end
 Ilustración 2 Gráfica del modelamiento de un horno en el dominio de la frecuencia

Modelamiento en Transformada Z

%% ************************************************************************
%************************* MODELACION DE UN HORNO** ************************
%**************************************************************************
for k=1:length(z)
% 2) Modelo no Lineal

temp(k)=Q(k)*(V/(C1+C2+C3))/(C1*rt)*(k/(k-exp(-
(C1+C2+C3)/((C2+C3)*(C1*rt))*T0)));
end
%% ************************************************************************
%******************** ANIMACION DE LOS ESTADOS ***********************
%**************************************************************************
figure ('Name','Modelado de la horno en funcion de la frecuencia')
graf1 = plot(z(1),temp(1),'r','LineWidth',2); hold on;
for i=1:50:length(z)

delete(graf1)
xlabel ('Frecuencia');
ylabel ('Temperatura (°C)');
graf1 = plot(z(1:i),temp(1:i),'r','LineWidth',2); hold on;
drawnow;
end
 Ilustración 3 Gráfica del modelamiento de un horno en transformada Z

2. SIMULINK
Función de transferencia

Ilustración 4 Función de transferencia de un horno en Simulink

Ilustración 5 Gráfica de la función de transferencia de un horno

Transformada Z

Ilustración 6 Transformada Z en Simulink

Ilustración 7 Gráfica en Simulink

3. Control System
Para guardar en Matlab la función de transferencia del horno se realizó el siguiente
programa:
clc, clear all
num=[0.199835];
den=[1 0.037717];
fprintf('Función de Transferencia \n')
Gs=tf(num,den)
step(Gs)
 Ilustración 8 Ejecución del programa

Ilustración 9 Gráfica de la función de transferencia

Ilustración 10 Programa en control system

 1. OBJETIVO

Realizar la simulación en Matlab, Simulink y Control System correspondiente a la tarea

del "Modelo matemático en función del tiempo, Laplace, Transformada Z."

2. DESARROLLO

2.1. DESARROLLO MATEMÁTICO

Figura 1 Variables que intervienen en un motor DC

En la figura 1, donde:

R= Resistencia eléctrica de la armadura del motor.

L= Inductancia de la bobina de la armadura del motor.

V= Tensión generada por la armadura del motor.

i = Corriente de la armadura.

ω = Velocidad angular del motor.

Tm = Torque generado por el motor.

J = Momento de Inercia.
B = Coeficiente de rozamiento viscoso.

ea = Fuerza contraelectromotriz

if = Corriente de Campo.

Para el presente modelo matemático las variables de interés son:

ω = Velocidad angular del motor

V= Tensión generada por la armadura del motor.

Puesto que lo que se pretende es realizar el control de la velocidad del motor mediante la
variación en la tensión

Modelo Eléctrico

Aplicando la Segunda Ley de Kirchoff en el circuito se tiene:

La segunda ley de Kirchoff nos dice que la suma algebraica de todos los voltajes en una
malla o bucle cerrado debe ser igual a cero. Expresada matemáticamente, la segunda ley
de Kirchoff se resume de la siguiente forma.

∑ Voltajes del bucle cerrado = 0

Sabiendo previamente que la caída de voltaje en una resistencia se representa como:

V𝐑 = Ri(t) (𝟏)

Se sabe que el voltaje en un inductor se calcula como:

di(t)
VL = L (𝟐)
dt
 VR + VL + Ea (t) − V(t) = 0 (𝟑)

Reemplazando (1) y (2) en (3)

di(t)
Ri(t) + L + Ea (t) − V(t) = 0 (𝟒)
dt

di(t)
Despejando L de (3) se tiene:
dt

di(t)
L = V(t) − Ri(t) − Ea (t) (𝟓)
dt

Donde Ea (t) es una tensión generada que resulta cuando los conductores se mueven a
través del flujo de campo establecido por la corriente de campo if .

Modelo Mecánico

La ecuación de la parte mecánica proviene del cálculo del torque del motor como se
muestra a continuación:

El cálculo de del momento de inercia se lo describe como

Mm = Jα

También se conoce que la aceleración angular es la derivada de la velocidad angular ω,

por lo tanto

dω(t)
Mm = J (𝟔)
dt
Se tiene que considerar que el coeficiente de fricción del motor respecto a la carga que
ejerce un momento el cual calcula de la siguiente manera:

Mr = βω(t) (𝟕)

El cálculo del torque del motor en el cual se hace uso de las ecuaciones (6) y (7) se lo
realiza como:

dω(t)
Tm (t) = J + Bω(t) (𝟖)
dt

En donde:

Tm (t) es el torque del motor.

B es el coeficiente de fricción cuando el motor se encuentra con carga.

J es el momento de inercia total del rotor y de la carga

ω(t) es la velocidad angular del motor

dω(t)
Despejando J de la ecuación (8) se tiene:
dt

dω(t)
J = Tm (t) − Bω(t) (𝟗)
dt

Las ecuaciones que van a regir el sistema son la ec. (5) y la ec. (9) y se asumen ciertas
relaciones descritas a continuación para linealizar el modelo.
Se asume que existe una relación proporcional, K a , entre el voltaje inducido en la
armadura y la velocidad angular del eje del motor. Se supone que cuando la armadura
está girando, se induce en ella un voltaje proporcional al producto del flujo por la
velocidad angular, para un flujo constante el voltaje inducido es directamente
proporcional a la velocidad angular por lo cual se tiene la siguiente ecuación:

E a ( t) = K a ω ( t) (𝟏𝟎)

Además, también se supone que existe una relación electromecánica la cual establece que
el torque mecánico es proporcional, K m a la corriente eléctrica. Se considera que, para
una corriente de campo constante, el flujo se vuelve constante y el par desarrollado por el
motor es directamente proporcional a la corriente de la armadura por lo tanto esto resulta
en:

Tm (t) = K m i(t) (𝟏𝟏)

A continuación, se aplica la transformada de Laplace en las ecuaciones (5), (9), (10) y

(11) obteniéndose así las siguientes ecuaciones:

Lsi(s) = v(s) − Ri(s) − Ea (s) (𝟏𝟐)

Jsω(s) = Tm (s) − Bω(s) (𝟏𝟑)

Ea (s) = K a ω(s) (𝟏𝟒)

Tm (s) = K m i(s) (𝟏𝟓)

De manera seguida se reemplaza la ecuación (14) y (15) en (12)

Tm (s) Tm (s)
Ls = v(s) − R − K a ω(s)
Km Km

(R + Ls)Tm (s)
v(s) = + K a ω(s) (𝟏𝟔)
Km

De la ec. (6) se puede obtener el torque del motor:

Tm (s) = (Js + B)ω(s) (𝟏𝟕)

Reemplazando la ec. (17) en (16) se obtiene:

(R + Ls)(Js + B)ω(s) + K a K m ω(s)

v(s) =
Km

(RJs + RB + LJs2 + LBs)ω(s) + K a K m ω(s)

v(s) =
Km

Obteniéndose así la función de transferencia velocidad angular – voltaje en dominio S

ω(s) Km
= 2
v(s) LJs + (RJ + LB)s + RB + K m K a
Tabla 1: Parámetros del motor utilizados para la simulación

Parámetro Símbolo Valor

Constante de torque Km 0.01

Constante de F.E.M. Ka 0.01

Resistencia de armadura R 1Ω

Inductancia L 0.5 H

Kg. m2
Inercia del motor con carga J 0.01
s2

Coeficiente de fricción B 0.1 N ∗ m

Voltaje de Alimentación V 1v

Función de Transferencia Dominio s

ω(s) 0.01
= 2
v(s) 0.005s + 0.06s + 0.1001

Transformada al Tiempo

2
G(s) =
s 2 + 12s + 20.02

Aplicamos fracciones parciales

 2
(s − 2.002500782)(s − 9.997499218)s
A B C
= + +
s s − 2.002500782 s − 9.997499218

2 = A(s − 2.002500782)(s − 9.997499218) + B(s − 9.997499218)s + C(s

− 2.002500782)s

2 = A(s 2 + 12s + 20.02) + Bs 2 − Bs ∗ 9.997499218 + Cs2 − Cs ∗ 2.002500782

Establecemos el siguiente sistema de ecuaciones

A+B+C = 0
{12 A − 9.997499218 B − 2.002500782 C = 0
20.02 A = 2

Resolviendo el sistema, tenemos el siguiente valor para las constantes

100
A=
1001

B = 0.1749657912

C = −0.2748658911

Reemplazamos en la ecuación 1

100
1001 + 0.1749657912 + −0.2748658911
s s − 2.002500782 s − 9.997499218
Aplicamos la transformada inversa de Laplace mediante el uso de tablas y tenemos la
siguiente ecuación:

100
+ 0.1749657912e2.002500782∗t − 0.2748658911e9.997459218∗t
1001

Transformada a Z

A partir de la ecuación de tiempo y aplicando las siguientes transformaciones de tabla de

transformación:

z z
e−at = y u ( t) =
z − e−aT z−1

Obtenemos:

100 z z
∗ + 0.1749657912 ∗ 2.002500782∗T
− 0.2748658911
1001 z − 1 z−e
z
∗ 9.997459218∗T
z−e

2.2. SIMULACIONES

2.2.1. SIMULACIÓN EN MATLAB

Con la cual se tiene la siguiente función de transferencia:

w(s) 2
= 2
v(s) s + 12 s + 20.02
Para realizar la simulación en Matlab y obtener las gráficas de la función de transferencia
en el dominio del tiempo, frecuencia y Z se utilizó el siguiente programa:

clc, clear all

disp(' UNIVERSIDAD DE LAS FUERZAS ARMADAS');

disp(' ESPE ');

fprintf('\n')

disp('CONTROL DIGITAL');

disp('NRC: 2057 ');

disp('Docente: Ing. Milton Pérez');

fprintf('\n \n')

syms s t k z

num=[2];

den=[1 12 20.02];

fprintf('Función de Transferencia en dominio de Laplace \n')

Gs=tf(num,den)

Gs_1=(poly2sym(num,s)/poly2sym(den,s));

Gs_2=Gs_1*(1/s);

fprintf('Funcion de transferencia en el dominio del Tiempo \n ');

H(t)=ilaplace(Gs_1)

fprintf('Funcion de transferencia en el dominio del Tiempo por escalon unitario \n ');

H2(t)=ilaplace(Gs_2)

fprintf('Transformada Z \n');

Gz=ztrans(H2)

Gk=iztrans(Gz,k);

%Gráficos de las funciones

subplot(2,2,1),fplot(heaviside(t-0)*H2, 'g','linewidth',3),xlim([0 10]),ylim([0

0.1]),title('Dominio del tiempo por escalon');

grid on;

subplot(2,2,2),fplot(H),xlim([0 10]),ylim([0 0.15]),title('Dominio del tiempo');

grid on

subplot(2,2,3),fplot(Gk, 'r','linewidth',2),xlim([0 10]),ylim([0 0.1]),title('Transformada

Z');

grid on;

subplot(2,2,4),step(Gs),xlim([0 10]),ylim([0 0.1]),title('Dominio de Laplace');

grid on;

Una vez compilado el programa se obtienen las siguientes gráficas:

Figura 2 Gráficas resultantes de la programación en Matlab

SIMULINK

Lo primero implementar el modelo matemático en función del tiempo en simulink

empleando las constantes del motor previamente definidas utilizando los diferentes
bloques que nos proporciona la herramienta Simulink como se muestra a continuación:
Figura 3 Modelo matemático en función del tiempo Velocidad Angular/Voltaje
(Simulink)

A continuación se realizó la implementación de la función de transferencia en dominio s

para poder observar su grafica posteriormente.

Figura 4 Función de Transferencia en Dominio S Velocidad Angular/Voltaje

A continuación, se realizó la implementación de la función de transferencia en dominio

Z para poder observar su grafica posteriormente.
Figura 5 Función de Transferencia en Dominio Z Velocidad Angular/Voltaje

Obteniéndose finalmente la implementación del modelo en dominio t, s y Z

Se ha utilizado un tiempo de muestreo de 10 para todas las simulaciones en simulink.

Figura 6: Gráfica de Respuesta del Sistema Dominio del Tiempo

Figura 7: Gráfica de Respuesta del Sistema Dominio S

 Figura 8: Gráfica de Respuesta del Sistema Dominio Z

2.2.2. CONTROL SYSTEM DESIGNER

Función de transferencia o transformada de Laplace del modelo matemático:

Para la respectiva simulación de la función de transferencia de un motor primero se

realizó un programa en el editor de Matlab con el fin de que nuestra función de
transferencia se guarde en una variable.

Programa en Editor de Matlab:

clc, clear all

disp(' UNIVERSIDAD DE LAS FUERZAS ARMADAS');

disp(' ************ ESPEL ************');

fprintf('\n')

disp('CONTROL DIGITAL');

disp('Modelamiento del motor');

disp('Docente: Ing. Milton Pérez');

disp('NRC: 2057 ');

fprintf('\n \n');

num=[2];

den=[1 12 20.02];

fprintf('Función de Transferencia \n')

Gs=tf(num,den)

step(Gs)

Resultado de la Ejecución del Programa

Figura 9 Función de transferencia ingresada

Grafica de la función
Figura 10 Gráfica resultante de la función de transferencia ingresada

A continuación, se abrió la extensión Control System, seleccionando la estructura del

sistema de control y asignando a Gs como función de transferencia de trabajo, así:

Figura 11 Localización de la aplicación Control System Designer en Matlab

Una vez abierto no saldrá la siguiente ventana para empezar a trabajar

 Figura 12 Interfaz de la aplicación Control System Designer

Seleccionamos la estructura del sistema de control en Edit Architecture nos saldrá una
nueva ventana:

Figura 13 Opción para elegir el tipo de diagrama de bloque

Asignación de la variable Gs como función de transferencia del sistema, presionamos en
la flecha color verde para asignar dicha función:

Figura 14 Interfaz de selección de diagrama de bloques

Asignamos la función de transferencia y le damos clic en Import.

Figura 15 Interfaz para la selección de los datos que se alojarán en un bloque

Visualizamos la simulación

Figura 16 Gráfica de la función de transferencia obtenida en Control System

Designer

También podemos obtener más gráficas, en este caso obtendremos las gráficas de impulso
y también la de los polos y ceros.

Figura 17 Opción de selección de otras gráficas disponibles dentro de Control

System Designer
Figura 18 Gráficas adicionales de la función de transferencia ingresada

RESULTADOS Y DISCUSIÓN

2.3. MATLAB

Se puede apreciar que tanto la gráfica en dominio del tiempo multiplicada por la función
escalón unitario, la transformada Z y la gráfica en el dominio de Laplace son similares,
teniendo la misma amplitud y estabilizándose en un tiempo similar. Por otro lado, la
gráfica en el dominio del tiempo sin ser multiplicada por la función escalón unitario,
corresponde a una gráfica totalmente diferente, pues solo posee un pulso hasta un máximo
de 0,12 y después de 0 es igual a 0 en el infinito del tiempo.

2.4. SIMULINK

Se puede apreciar que las tres graficas presentan similitud tanto en sus amplitudes que
son de 0.1 para los tres casos, así como también en su tiempo de establecimiento
notándose una pequeña diferencia en la Transformada Z debido a que es una señal
discretizada esa muestra ha sido tomada en un tiempo muy pequeño por delante de las
tras graficas en dominio del tiempo y dominio s. Se puede apreciar que la reconstrucción
de la gráfica en el dominio Z se lo realizó correctamente acercándose en gran medida a la
gráfica de los otros dominios
 2.5. CONTROL SYSTEM DESIGNER

Con estás gráficas resultantes podemos observar que la función de transferencia es igual
a la halla en los otros métodos, tanto como los ploteos en Matlab así como en los scope´s
en Simulink, teniendo una amplitud 0.1 y un tiempo de establecimiento aproximadamente
en 3. Otra gráfica que podemos comparar es la del impulso que también fue obtenido en
Matlab a partir de la función en el tiempo.

3. CONCLUSIONES

Se pudo comprobar que el modelo matemático utilizado trabaja correctamente empleando

la simulación en Simulink mostrando sus gráficas de respuesta del sistema en las cuales
se pudo observar una gran similitud y sobre todo la discretización de la señal en el
dominio Z.

Control System Designeer es una extensión de Matlab que nos permite elaborar y diseñar
un sistema de control a partir de parámetro establecidos por nosotros, observa el
comportamiento de la planta y del controlador como tal, siendo así una herramienta muy
útil en el campo de la ingeniería.
Realizar la simulación en Matlab, Simulink y Control System correspondiente a la tarea
del "Modelo matemático en función del tiempo, Laplace, Transformada Z de una Planta
de Nivel"

Función de Transferencia

0.0636
𝐺(𝑠) =
𝑠 + 0.0266

1. MATLAB

Programa:

clc, clear all

disp(' UNIVERSIDAD DE LAS FUERZAS ARMADAS');
disp(' ************ ESPEL ************');
fprintf('\n')
disp('CONTROL DIGITAL');
disp('Docente: Ing. Milton Pérez');
disp('NRC: 2057 ');
fprintf('\n \n')
syms s t k z
num=[0.0636];
den=[1 0.0266];
fprintf('Función de Transferencia \n')
Gs=tf(num,den)
Gs1=(poly2sym(num,s)/poly2sym(den,s));
Os=Gs1*(1/s);
fprintf('Modelo matemático o Transformada inversa de
Laplace \n ')
M(t)=ilaplace(Os);
pretty(M);
fprintf('Transformada Z \n')
Gz=ztrans(M);
pretty(Gz);
Gk=iztrans(Gz,k)
%Gráficos
fplot(heaviside(t-0)*M, 'g','linewidth',3), hold on
fplot(Gk, 'r','linewidth',2); hold on
step(Gs)
hold on
grid on
xlim([0 200]);
xlabel('t,s,z ')
ylim([0 2.5]);
ylabel('M, Gs, Gz ')
title('Gráficos en Matlab')
legend('Modelo Matemático','Transformada
Z','Laplace','Location','south')

Gráficos:

Comentario: Se observa que las tres gráficas inciden en la misma trayectoria por lo
cual para su apreciación se ha realizado las líneas con distinto grosor.
 2. SIMULINK

Modelo Matemático

 Desarrollo del esquema de trabajo

 Grafica obtenida en el scope tras ejecutar la simulacion

Comentario: Se observa que el modelo en Simulink se mantiene al generado en Matlab

Funcion de transferencia o transformada de Laplace

 Esquema de trabajo realizado

 Resultado obtenido tras ejecutar la simulacion

Comentario: La trayectoria de la función de transferencia se mantienen a la del modelo

matemático.
Transformada Z

 Esquema de trabajo realizado

 Resultado obtenido tras ejecutar la simulacion

Comentario: Se observa que la gráfica incide en la misma trayectoria que el modelo

generado en Matlab.
 3. CONTROL SYSTEM

Función de transferencia o transformada de Laplace del modelo matemático

Para generar la simulación de la función de transferencia en control system designer en primer

lugar se realizó un pequeño script, de modo que la función de transferencia deseada, sea
guardada en la memoria de Matlab, así:

clc, clear all

disp(' UNIVERSIDAD DE LAS FUERZAS ARMADAS');
disp(' ************ ESPEL ************');
fprintf('\n')
disp('CONTROL DIGITAL');
disp('Docente: Ing. Milton Pérez');
disp('NRC: 2057 ');
fprintf('\n \n')
num=[0.0636];
den=[1 0.0266];
fprintf('Función de Transferencia \n')
Gs=tf(num,den)
step(Gs)

Resultado de la ejecución del programa:

 Command Window
  Grafica de la función

A continuación, se abrió y genero la simulación en Control System, seleccionando la estructura

del sistema de control y asignando a Gs como función de transferencia de trabajo, así:

 Selección de la estructura del sistema de control

 Asignación de la variable Gs como función de transferencia del sistema

 Al ser un sistema en lazo abierto, se le asigna un valor de cero a la variable H

  Generación de la simulación al aceptar las configuraciones con el botón OK

Comentario: Se logra visualizar que tanto la gráfica generada en el script de Matlab como la
desarrollada en Control System son exactamente iguales, con lo cual se comprueba que el
programa es perfectamente funcional.

 Comprobacion de funcionalidad y similitud de las simulaciones en Matlab, Simulink

y Control System