Saltar al contenido
Control Automático Educación

Modelo Matemático de un Sistema Térmico

En esta entrada vamos a realizar el modelo matemático de un sistema térmico el cual puede ser reproducido a sistemas clásicos como hornos, refrigeradores, calderas, etc. Y vamos a comparar el modelo matemático con un sistema real para validar el modelo.

Antes de comenzar, te hago la invitación de que mires nuestro CURSO GRATUITO DE ANÁLISIS DE SISTEMA.

Y que te suscríbas al canal de YouTube para seguir aprendiendo más sobre estos temas.

Transferencia de Calor

Existen diferentes formas de transferencia de calor: por radiación, convección y conducción.  

La transferencia de calor por radiación proviene del movimiento de partículas y se emite como fotones electromagnéticos principalmente como radiación infrarroja.  

La transferencia de calor por convección se da a través del fluido circundante, como por ejemplo el aire. La convección puede ser forzada, por ejemplo usando un ventilador o puede ocurrir de forma natural.  

La transferencia de calor por conducción es a través del contacto físico de dos materiales donde ocurre la transferencia de energía, pero, a diferencia de la convección, el material de contacto es estacionario

Modelos Térmicos

Vamos a ver el modelado térmico de un proceso, el cual vamos posteriormente modificar o controlar su temperatura.  

Existen diversos procesos térmicos que podríamos analizar: horno, duchas, refrigeradores o calentamiento de calderas.  

Sistemas Térmicos
Diferentes ejemplos de sistemas térmicos

Para este caso, vamos a tomar un componente electrónico y vamos modelar su comportamiento dinámico.

Vamos a utilizar un transistor que será usado como un calentador, en este caso aprovechamos para usar una tarjeta conocida como el Laboratorio de control de temperatura, diseñada por el  Profesor John D. Hedengren de la universidad de  Brigham Young University.

Posteriormente, en otra entrada vamos a realizar un control de temperatura con Arduino usando esta tarjeta.

Toda la información del modelo matemático que veremos a continuación junto con información adicional de la tarjeta, implementaciones de controladores en Matlab o Python pueden ser encontradas en la página de APMonitor o en el canal de YouTube de APMonitor.

Modelo Matemático de Sistema Térmico – Transistor BJT

En este proceso, vamos a necesitar conocer todo tipo de temperatura que entra y sale del sistema.  

T. Entrada: La temperatura del transistor aumenta a medida que la corriente eléctrica fluye a través del pequeño dispositivo  

T. Salida: La energía que se dispersa fuera del transistor es provocada por dos mecanismos principales: convección y radiación.

\dfrac{\Delta Q}{\Delta t}=Q_{in}-Q_{out}

La temperatura de entrada:

Q_{in}=\alpha Q_i

Q_i es el porcentaje de salida del calentador. El parámetro α es un factor que relaciona la salida del calentador (0-100%) con la potencia disipada por el transistor en vatios.

La temperatura de salida es dada por:

Q_{out}=k_T(T-T_{\infty})+\epsilon \sigma A (T^4-T_{\infty}^4)

El primer término de la ecuación se conoce como la Ley de enfriamiento de newton y el segundo como la Ley de Stefan-Boltzmann

k_T constante de perdida para el ambiente, T temperatura del transistor,
T_\infty temperatura del ambiente, ϵ  emisividad, σ es la constante de Stefan-Boltzmann, A es el área

donde:

k_T=UA

U es el coeficiente de transferencia de calor, A es el área

Reemplazando valores tenemos que:

\dfrac{\Delta Q}{\Delta t}=\alpha Q_i+UA(T_{\infty}-T)+\epsilon \sigma A (T_{\infty}^4-T^4)

Utilizando el concepto de derivada

\dfrac{d Q}{d t}=\alpha Q_i+UA(T_{\infty}-T)+\epsilon \sigma A (T_{\infty}^4-T^4)

Podemos calcular el calor liberado o absorbido utilizando el calor específico, c_p, la masa, m, del componente y el cambio en temperatura (T-T_{ref}) en la ecuación:

Q=mc_p(T-T_{ref})

reemplazando tenemos:

mc_p\dfrac{d T}{d t}-mc_p\dfrac{d T_{ref}}{d t}=\alpha Q_i+UA(T_{\infty}-T)+\epsilon \sigma A (T_{\infty}^4-T^4)

Considerando que la temperatura de referencia es constante, por lo tanto la derivada da CERO.

mc_p\dfrac{d T}{d t}=\alpha Q_i+UA(T_{\infty}-T)+\epsilon \sigma A (T_{\infty}^4-T^4)

Parámetros del Modelo

Los parámetros del modelo fueron tomados de la pagina de APMonitor, modificando el coeficiente de transferencia de calor y el factor del calentador para Mi caso específico:

QuantityValue
Temperatura inicial (T0)296.15 K (23oC)
Temperatura Ambiente (T∞)296.15 K (23oC)
Salida del Calentador (Q)0 to 1 W (0%-100%)
Factor del Calentador (α)0.014 W/(% heater)
Capacidad Calorifica (Cp)500 J/kg-K
Area de la superficie (A)1.2×10-3 m2 (12 cm2)
Masa (m)0.004 kg (4 gm)
Coeficiente de transferencia de calor(U)5 W/m2-K
Emisividad (ε)0.9
Stefan Boltzmann Constant (σ)5.67×10-8 W/m2-K4

Linealización del Modelo Térmico

Por causa del término de T^4 agregado por la ecuación de radiación de la lay de Stefan Boltzmann, el modelo fenomenológico de la temperatura del transistor es NO lineal.

Procedemos a linealizar el modelo en un punto de equilibrio usando la expansion por series de Taylor.

Comenzamos definiendo el modelo como una función f:

f(T,Q_i)=\dfrac{d T}{d t}=\dfrac{\alpha}{mc_p} Q_i+\dfrac{UA}{mc_p}(T_{\infty}-T)+\dfrac{\epsilon \sigma A}{mc_p} (T_{\infty}^4-T^4)

Punto de Equilibrio

En el punto de equilibrio la derivada es nula, debido a que ya no existe más variación de temperatura dado que el sistema se ha estacionado.

0=\alpha Q_{is}+UA(T_{\infty s}-T_s)+\epsilon \sigma A (T_{\infty s}^4-T_s^4)

De esa forma llegamos a una ecuación algebraica de orden 4:

(\epsilon \sigma A)T_s^4+(UA)T_s-(\alpha Q_{is}+UAT_{\infty s}+\epsilon \sigma A T_{\infty s}^4)=0

Debemos encontrar las raíces de la ecuación anterior (serán 4 raíces) y seleccionar únicamente la raíz positiva. Note que la ecuación anterior depende de la entrada Q_{is} que será la cantidad de potencia que vamos a aplicar al transistor para que genere calor.

La siguiente función en MATLAB puede ayudarlo a encontrar la raíz positiva de ese polinomio y entregar su resultado en grados Kelvin o en grados Celsius. Ta debe ser ingresada en Kelvin.

function [TK,TC] = T_ss(eps,sigma,A,U,Qi,Ta,alpha)
%Equilibrio
    p=[eps*sigma*A 0 0 U*A -alpha*Qi-U*A*Ta-eps*sigma*A*Ta^4];
    rp=roots(p);
    Ts=[];
    for i=1:4
        if isreal(rp(i)) && rp(i)>0
            Ts=[Ts;rp(i)];
        end
    end
  TK=Ts;
  TC=TK-273.15;

Series de Taylor

Aplicando las series de Taylor:

\dfrac{dT}{dt}=\left.\begin{matrix}\dfrac{\partial f}{\partial T}\end{matrix}\right|_{T_s,Q_{is}}(T-T_s)+\left.\begin{matrix}\dfrac{\partial f}{\partial Q_i}\end{matrix}\right|_{T_s,Q_{is}}(Q_i-Q_{is})+f(T_s,Q_{is})

Realizando las derivadas parciales con relación a la temperatura y a la entrada evaluando todo en los puntos de equilibrio encontrados anteriormente llegamos a la siguiente ecuación:

\dfrac{dT}{dt}=\left(-\dfrac{UA}{mc_p}-\dfrac{4\epsilon \sigma A}{mc_p} T_s^3\right)\Delta T+\left( \dfrac{\alpha}{mc_p}\right)\Delta Q+\left.\begin{matrix}\dfrac{dT}{dt}\end{matrix}\right|_{T_s,Q_{is}}

\dfrac{d\Delta T}{dt}=\left(-\dfrac{UA}{mc_p}-\dfrac{4\epsilon \sigma A}{mc_p} T_s^3\right)\Delta T+\left( \dfrac{\alpha}{mc_p}\right)\Delta Q

Función de Transferencia

Aplicando la transformada de Laplace al modelo lineal encontrado usando las series de Taylor:

\Delta T(s)s+\left(\dfrac{UA}{mc_p}+\dfrac{4\epsilon \sigma A}{mc_p} T_s^3\right)\Delta T(s)=\left( \dfrac{\alpha}{mc_p}\right)\Delta Q(s)

Ordenando tenemos que:

\dfrac{\Delta T(s)}{\Delta Q(s)}=\dfrac{\left( \dfrac{\alpha}{mc_p}\right)}{s+\left(\dfrac{UA}{mc_p}+\dfrac{4\epsilon \sigma A}{mc_p} T_s^3\right)}

Reemplazando los valores de la tabla y encontrando un punto de equilibrio cuando NO estoy inyectando potencia en el transistor Q_{is}=0 voy a encontrar que la temperatura en equilibrio será igual a la temperatura del ambiente T_s=23C:

G(s)=\dfrac{0.007}{s + 0.006213}

Representando el sistema de primer orden en su forma estándar:

G(s)=\dfrac{1.127}{160 s + 1}

Retardo de Tiempo

Sin embargo, para que la temperatura en el transistor comience a ser evidenciado por el sensor, puede tomar algún tiempo.  

Eso se conoce como el retardo del sistema. En nuestro caso, para dejar nuestro modelo mucho mas fiel al proceso real, podemos incrementar un retardo el cual podremos determinar experimentalmente (En mi caso encontré un retardo de 10 segundos).

Realizando experimentos, también opté por reducir un poco la ganancia de la función de transferencia.

G(s)=\dfrac{1.04}{160 s + 1}e^{-10s}

Validación del Modelo

A continuación veremos el comportamiento dinámico del sistema real vs los modelos lineales y no lineales que hemos desarrollado aquí:

Modelo Termodinámico
Validación del Modelo Termodinámico

Descargar Código en Matlab

A continuación te dejo el código en matlab del modelo matemático del sistema térmico del transistor BJT donde se hace la comparación entre el modelo lineal y el modelo no lineal.

[sociallocker id=”948″]

[/sociallocker]

Eso es todo por la entrada del dia de hoy, espero les haya gustado y hayan aprendido algo nuevo. Si te ha servido el contenido de esta entrada, de los videos y los códigos de implementación y deseas apoyar mi trabajo invitandome a un café super barato, puedes hacerlo en el siguiente link:

👉 Invitar a Sergio a un Café ☕️

Que esten muy bien, nos vemos en la siguiente entrada.