DESARROLLO EN MATLAB SIMULINK DEL MODELO MATEMATICO PARA SIMULADOR DE UNA CALDERA Fredy Antonio Estupiñán1, Jaime Tenjo García2 1,2 Universidad ECCI – Maestría Maestría en Ingeniería Cr 19 No. 49-20 1
[email protected] [email protected] [email protected] 2 jaimete
RESUMEN El proyecto consiste en realizar el desarrollo teórico de la generación de potencia de una caldera, la cual estará ubicada en Bucaramanga. En una primera etapa del diseño se desarrolla y analiza la reacción química de todos los elementos que intervienen en el subsistema de combustión, luego en una segunda etapa se desarrolla y analiza los elementos elementos que intervienen en en el subsistema de hogar, los cálculos de calor generados y perdidos durante el proceso y el balance energético del mismo, por último se elabora las ecuaciones de los bancos de convección, los sobrecalentadores y los cálculos c álculos de eficiencia y gasto de vapor de la caldera. Durante todo el proyecto, se elaboran varias ecuaciones donde intervinieron las diferentes entalpias de cada proceso; se generan los cálculos para cada variable generada que intervienen en los diferentes subsistemas. Con estos datos se elaboraron bloques donde se consignaban los valores calculados de cada uno de estos elementos, se anexan al sistema de combustión previamente desarrollado, para que finalmente se realizara la simulación del funcionamiento real de los subsistemas empleando el software Matlab-Simulink. El objetivo principal es conocer el comportamiento del sistema cuando es perturbado con diferentes valores de entrada como humedad del combustible, temperatura ambiente, temperatura del agua de alimento, presión de trabajo, humedad del aire, temperatura de sobrecalentamiento del vapor, delta de temperatura del sobrecalentador, caudal de combustible
PALABRAS CLAVE: caldera, entalpia, hogar, balance, simulación. 1. INTRODUCCIÓN El propósito de este proyecto es lograr realizar la simulación del subsistema de hogar de una caldera, según las siguientes condiciones: ubicación de la caldera en Bucaramanga, y porcentajes del peso de cada uno de los elementos los cuales se describirán posteriormente. Para lograr esta simulación se debe desarrollar previamente las ecuaciones de balance de energía, las cuales deben indicar los calores generados y perdidos dentro del subsistema hogar. Los datos generados dentro del banco de convección y la eficiencia lograda con los parámetros establecidos por el autor. Se seleccionaron un conjunto de datos de entrada que corresponden a un estado estable inicial de referencia: la combustión se utiliza chips de madera a una razón de 3500 kg/h de composición y 25% de humedad, 30% de aire en exceso con una humedad del 30%, presión de trabajo de 15 bares. Las temperaturas para el aire ambiente de 22°C, agua de alimentación 105°C, de sobre-calentamiento 220°C y delta de T en el pre-calentador de aire de 130°C.
2. OBJETIVO GENERAL Señalar el desarrollo matemático que se requiere para el subsistema de combustión de una caldera, el subsistema de hogar, los bancos de convección y todas las ecuaciones necesarias para obtener la eficiencia de la caldera, para finalmente finalmente demostrar mediante mediante los resultados resultados obtenidos obtenidos en la la simulación hecha con el software Matlab-Simulink, la coherencia de los resultados obtenidos en todas las ecuaciones del sistema con los obtenidos en la simulación.
3. OBJETIVOS ESPECÍFICOS - Aplicar y comprender conceptos y principios de entalpia, entropía y balance de materia. - Identificar la ecuación química del subsistema de combustión. - Identificar las ecuaciones de balance térmico de la caldera. - Realizar la simulación de los diferentes subsistemas con el software Matlab-Simulink.
4. DESARROLLO Para el desarrollo de este proyecto se tuvo en cuenta las siguientes condiciones iniciales: -
Ubicación de la caldera en la ciudad de Bucaramanga. Temperatura ambiente 22° Humedad del aire 83% El combustible a utilizar es Madera. H = 6.2% N = 0.13% C = 44.59% O = 48.78% S = 0% Ash = 0% Área de Hogar empleada para cálculos 51m 2 Caudal de combustible 3500 KgFuel/Hora Emisividad 0.85 0. 85 Porcentaje de la humedad 25%, Debido a que la humedad influye negativamente en el poder calorífico [1], se empleó la figura 1 para determinar el poder calorífico de la leña:
Fig. 1: Poder calorífico superior (Fernández, R. A, 2009).
Según figura 1, el poder calorífico de la leña con un 25% de humedad [2] es de 14068 KJ/Kg. Para la simulación en Matlab-Simulink se establece programar la función de determinar el poder calorífico con base del porcentaje de humedad ingresado. Al tener una función lineal en este cálculo, se emplea la función lineal y=mx + b, donde según la figura 1 se obtiene m=-44 y b=15168, de acuerdo a datos interpolados de la función.
Cálculo de Qps: Para este cálculo teniendo como base el porcentaje de humedad del combustible y aplicando la función lineal descrita (y = -44x+15168), se obtiene el bloque de la figura 2.
Fig. 2: Cálculo poder calorífico del combustible (Autor Matlab-Simulink R2016a) Para obtener el porcentaje del exceso de aire con relación a los gases que resultan de la combustión de la madera, se empleó la figura 3:
Fig. 3: Exceso de aire en la combustión (Thomasset, C.W Pag 128)
El proyecto se desarrolla considerando una combustión incompleta y de acuerdo a la Figura 3, con un exceso de aire del 30% (zona amarilla), se obtienen un 14,74% CO2 y un 4,647% O2. Cabe resaltar que según la literatura para este proceso se sugiere un exceso de aire óptimo entre 30% y 50% de exceso. [3]
4.1.
COMBUSTIÓN
La ecuación que se usará para la reacción química es la siguiente: C + H2 + O2 + N2 + S + (O2 + 3.76N2)
CO2 + H2O + N2 + O2 + SO2 + CO
Se calcula el flujo molar: C = 12.01 gr/mol H = 1 gr/mol O = 16 gr/mol N = 14 gr/mol S = 32.06 gr/mol Dónde
Datos de entrada: Se determina dejar variables estos datos de entrada como se observa en la figura 4, al emplear un software de simulación y empleando las ecuaciones dadas por el autor [2], se implementa el programa con la función de realizar los cálculos para diferentes valores de entrada y no limitarse a una condición estable.
Fig. 4: Variables datos de entrada (Autor Matlab-Simulink R2016a)
Además de los datos del combustible, se dan como datos de entrada el porcentaje de exceso de aire, la humedad del combustible y el porcentaje de carbono no quemado. 3.71C + 3.1H2 + 1.52O2 + 0.0046N2 + 0S + α(O2 + 3.76N2) βCO2 + γH2O + μ N2 + νO2 + δSO2 + εCO
ENTRA
SALE
1) 3.71 = β + ε ε=1 2.71 = β
2) 3.1 = γ 3) 1.52 + α = β + γ/2 + ν + δ + ε/2 α = 2.71 + 3.1/2 + 1 + 0 + 1/2 - 1.52 α = 4.226
4)
0.0046 + α x 3.76 = μ 0.0046 + 4.22 * 3.76 = μ 15.89 = μ
5) ν= 1
Cálculo de flujo molar : Para este cálculo se emplean constantes para especificar el peso atómico de cada elemento, un bloque sumador para tener en cuenta el factor del porcentaje de carbono no quemado y bloques divisores para determinar finalmente el cálculo requerido según se observa en la figura 5.
Fig. 5: Cálculo del flujo molar (Autor Matlab-Simulink R2016a)
Cálculo de alpha: Este es uno de los cálculos más importantes en el proceso del balance de la ecuación, debido a que esta variable es clave para el balance del oxígeno y el nitrógeno. En la figura 6 se muestra los bloques empleados para el desarrollo de la ecuación 3.
Fig. 6: Cálculo de alpha (Autor Matlab-Simulink R2016a) Se verifica el exceso de aire en el proceso, se utiliza la ecuación dada por el autor: [2]
) ( % Para obtener el gasto de aire utilizado en la reacción se calcula de la siguiente forma:
α ) ( ) ( Para obtener el gasto de humos húmedos calculamos: GHH = 1 + Ga - % Ash - %C GHH = 1 + 583.18 – 0 – 44.59 GHH = 539.54
Para cálculos posteriores se consideran a los humos sin la humedad del combustible, por lo que definimos el gasto de humos semi-húmedos con la siguiente ecuación: GH9H = GHH - % (H2O) f GH9H = 539.54 - 25 GH9H = 514.54
Cálculo de gastos de salida: Como se observa en la figura 7 se emplean bloques sumadores, bloques de producto, constantes y valores obtenidos de ecuaciones anteriores.
Fig. 7: Cálculo de gastos de salida (Autor Matlab-Simulink R2016a)
Cálculo de salida de CO y O 2: Como se observa en la figura 8 se tienen unas entrada calculadas en procesos anteriores, se emplean 2 bloques sumadores, 2 bloques divisores y dos generadores de ganancia para dar las respuestas de salidas en porcentajes
Fig. 8: Cálculo de porcentaje salida de CO y O 2 (Autor Matlab-Simulink R2016a)
Resultados balance ecuación entrada: De acuerdo a la figura 9 se tienen los valores para el balance de la ecuación de entrada, estos valores cambian automáticamente si se cambia alguna variable
Fig. 9: Resultados balance productos de entrada C+H 2+O2+N2+S2 (Autor Matlab-Simulink R2016a)
Resultados balance ecuación salida: De acuerdo a la figura 10 se tienen los valores para el balance de la ecuación de salida, estos valores cambian automáticamente si se cambia alguna variable
Fig. 10: Resultados balance productos de salida CO 2+H2O+N2 O2+SO2+CO (Autor Matlab-Simulink R2016a)
Salida Qps, Ga y GH9H: De acuerdo a la figura 11 se tienen los resultados de salida para los parámetros de poder calorífico del combustible con 25% de humedad, gastos de aire y gastos de humos semi húmedos que se utilizaran en el desarrollo del subsistema de hogar.
Fig. 11: Resultados Qps, Ga y GH9H (Autor Matlab-Simulink R2016a) La ecuación a balancear en este subsistema hogar es la siguiente, donde cada uno de los componentes está dado en KJ/KgFuel
A continuación se desarrollará cada uno de los componentes de la ecuación y se verificará el balance energético del subsistema hogar.
4.2 CALCULO DE CALOR DE HOGAR (Qh). Para el desarrollo de la ecuación de calor de hogar o intercambio radiativo entre pared y hogar calculado se emplea la siguiente ecuación.
̇ Dónde SB corresponde a la constante de Stefan-Boltzmann 20.42x10 -8 KJ/h.m2.k 4 ε corresponde a la emisividad = 0. 85
Ah corresponde al área de hogar establecida en 51m 2 Wf es el caudal de combustible = 3500 KgFuel/Hora. Th temperatura de hogar. Tp temperatura de pared. La temperatura de pared es comúnmente aproximada a la temperatura de saturación mas 150K, esta temperatura de saturación es obtenida de tablas termodinámicas e interpolando el valor en condiciones de presión y temperatura para la ciudad de Bucaramanga se obtiene un valor de 471.42K, por lo que el resultado de la temperatura de pared es 621.42K. Según la ecuación del cálculo de calor de hogar, la temperatura de hogar debe ser mayor a la temperatura de pared, para evitar que de un valor negativo en el resultado, para garantizar lo anterior se maneja una segunda ecuación donde Tp = (Th + Tsat) / 2, para lo cual se debe manejar un bloque de saturación dinámica en simulación que cumpla con tal condición como se muestra en la figura 12.
Fig. 12: Bloque cálculo calor de hogar (Autor Matlab-Simulink R2016a)
Para el cálculo de la temperatura de hogar se debe iterar este valor con la finalidad de que el resultado de la sumatoria presente en el subsistema se iguale a cero, el autor [2] especifica que los valores de emisividad, área del hogar y porcentaje de pérdidas pueden ser variados para corresponder a las necesidades del modelo. Se emplea un bloque algebraic constraint encargado de realizar esta función y realizar un ajuste al balance energético del sistema como se muestra en la figura 13.
Fig. 13: Ajuste Th con bloque algebraic constraint (Autor Matlab-Simulink R2016a) El bloque de restricción algebraica limita la señal de entrada f (z) a cero y genera un estado algebraico z. El bloque emite el valor que produce un cero en la entrada. La salida debe afectar la entrada a través de una trayectoria de realimentación directa, es decir, la trayectoria de realimentación contiene solo bloques con avance directo. El resultado al reemplazar los valores en la ecuación es de 43714 KJ7Kg
4.3 CÁLCULO CALOR DE HUMOS SEMIHÚMEDOS (QH9H) Las pérdidas de calor por los humos semihúmedos de combustión se obtienen por la diferencia del calor de los productos de combustión y el calor total obtenido en aire frio que ingresa a la caldera.
Donde GH9H corresponde al gasto de húmedos semihúmedos obtenidos en la primera etapa del proceso y que corresponden a 5.15% hH9H es la entalpia de humos semihúmedos a temperatura K, valor que se obtiene de tablas termodinámicas y corresponde a -17.31 KJ/Kg To corresponde al valor de temperatura ambiente, en este caso es de 22° correspondientes a la ciudad de Bucaramanga, equivalentes a 295.15K.
El resultado de este cálculo es de -40794 KJ/Kg, el símbolo negativo es debido a un proceso exotérmico de calor en el sistema. El bloque empleado para este cálculo se muestra en la figura 14.
Fig. 14: Cálculo de calor de humos semihúmedos (Autor Matlab-Simulink R2016a)
4.4 CÁLCULO DE CALOR DE AIRE SECO DE HUMOS DE COMBUSTIÓN (Qair)
La ecuación para este cálculo es la siguiente:
Donde Ga es el gasto de aire calculado previamente en la etapa 1 y dado en 5.83% Wa es la humedad del aire en este caso 83% correspondientes a la ciudad de Bucaramanga. Cpa es la capacidad calorífica promedio del aire entre Ta y To = 1.017 KJ/Kg°C Ta corresponde a la temperatura del aire precalentado (125°C) o To si no hay precalentamiento. Qair al reemplazar estos valores en la ecuación igual a 104 KJ/Kg. El bloque empleado para este cálculo se muestra en la figura 15
Fig. 15: Cálculo de calor de aire seco (Autor Matlab-Simulink R2016a)
4.5 CÁLCULO DE CALOR DEL COMBUSTIBLE SECO (Qfuel)
Fig. 16: Cálculo de calor de combustible seco (Autor Matlab-Simulink R2016a)
Para determinar este valor se emplea la siguiente ecuación:
Donde (H2O)f corresponde a la humedad del combustible dada en 25% cpf es la capacidad calorífica promedio del combustible entre T f y To correspondiente según tablas termodinámicas a 1.38KJ/Kg°C. Al reemplazar estos valores en la ecuación realizando el bloque como se muestra en la figura 16, se obtiene como resultado para Qfuel = 107KJ/Kg
4.6 CÁLCULO DE CALOR DE CAMBIO DE FASE DEL AGUA (Qhfg) Según el autor se obtiene con la siguiente ecuación que es simulada según la figura 17.
Fig. 17: Cálculo de calor de cambio de fase de agua (Autor Matlab-Simulink R2016a) Donde H corresponde a la cantidad de hidrógeno en el combustible según el último análisis equivalente a 6.2%, se utiliza en bloque gain para especificar este valor igual a 0.062. hfg corresponde a la entalpia de cambio de estado del agua a temperatura T, este valor se obtiene de las tablas termodinámicas y equivalea Tch = 2852.8 KJ/Kg. Al reemplazar todos los valores en la ecuación dada se obtiene Qhfg = 1592KJ/Kg
4.6 CÁLCULO DE CALOR DE AGUA DE COMBUSTIBLE (Q(H2O)f ) La ecuación empleada en este cálculo es la siguiente:
Acá se requiere cpv y cpw correspondientes a la capacidad calorífica promedio del agua entre Th y To y Tf y To respectivamente. De las tablas termodinámicas se obtiene cpv = 2.63KJ/Kg°C y cpw = 1.96KJ/Kg°C.
Fig. 18 Cálculo de calor de combustible (Autor Matlab-Simulink R2016a) El valor de entalpia hfg se importa desde el anterior bloque, así como los valores correspondientes a la humedad del combustible. En la figura18 se muestra el bloque creado para realizar el cálculo correspondiente dando como resultado Q(H2O)f =968KJ/Kg.
4.6 CÁLCULO DE CALOR DE HUMEDAD DEL AIRE (Qwa) La ecuación empleada en este cálculo es la siguiente:
Donde Wa corresponde al porcentaje de humedad del aire correspondiente a 83%. En la figura 19 se muestra el bloque correspondiente a este cálculo.
Al reemplazar los valores obtenidos da como resultado Qwa = 5827 KJ/Kg , nótese que los valores de cpv, Ga y Th – To fueron tomados de otros bloques en el modelo.
Fig. 19 Cálculo de calor de la humedad del aire (Autor Matlab-Simulink R2016a)
4.7 CÁLCULO DE CALOR DE CO Y CARBONO NO QUEMADO (QCO+C) La ecuación empleada en este cálculo es la siguiente:
Dónde: 23958 es una constante que corresponde al cociente entre el poder calorífico del CO expresado en KJ/KMolCO y el peso molar de C expresado en KJ/KgFuel. CO: Es la cantidad de CO resultante de la combustión, según el cálculo de la primera etapa corresponde a 1 KMol/KgFuel. CO2: Es la cantidad de CO 2 resultante de la combustión, según el cálculo de la primera etapa corresponde a 2.71 KMol/KgFuel. C: Corresponde al porcentaje de carbono en combustible, según el último análisis equivale a 44.59%. C’: Corresponde al porcentaje de carbono sin quemar, según el último análisis equivale a 0%.
En la figura 20 se muestra el bloque correspondiente a este cálculo.
Al reemplazar los valores obtenidos da como resultado QCO+C = 2834 KJ/Kg
Fig. 20 Cálculo de calor de CO y carbono no quemado (Autor Matlab-Simulink R2016a)
4.8 CÁLCULO DE CALOR PERDIDO (%Perd.Qps) Es el porcentaje del Qps (poder calorífico superior del combustible) perdido a través de las paredes del ambiente mediante radiación y convección. Según el autor se aproxima comúnmente a un 2%.
Fig. 21 Cálculo de porcentaje de pérdidas de Qps (Autor Matlab-Simulink R2016a) Teniendo un Qps = 14068KJ/Kg y de acuerdo a la figura 21 este valor es 281KJ/Kg
Datos de entrada: Se determina dejar variables estos datos de entrada como se observa en la figura 22, al emplear un software de simulación y empleando las ecuaciones dadas por el autor [2], se implementa el programa con la función de realizar los cálculos para diferentes valores de entrada y no limitarse a una condición estable.
Fig. 22 Bloque subsistema hogar (Autor Matlab-Simulink R2016a)
Por cálculo teórico = 43714 + (-40794) + 104 + 107 + 1592 + 968 + 5827 + 2834 + (-281) Qps Teorico = 14071 KJ/Kg Qps Simulación = 14065 KJ/Kg Qps Real = 14068 KJ/Kg
4.9 BANCOS DE CONVECCIÓN : Los bancos de convección fueron resueltos mediante un método iterativo igualando NTU-eficiencia Dimensiones de diseño del banco de convección. Las propiedades de los componentes de los humos fueron determinadas con correlación en función de la T de entrada y a la T de salida del sistema para el primer y segundo paso por tubos.
4-10 EFICIENCIA DE LA CALDERA: Se obtiene por la relación entre el calor realmente producido (calor útil) en la caldera y la capacidad de producir calor del combustible empleado. Para el cálculo de esta eficiencia se emplea al método directo y la ecuación se muestra a continuación y la simulación se muestra en la figura 23
Dónde h2 corresponde a la entalpia del vapor saliente, vapor saturado o sobrecalentado a la presión de trabajo expresada en (KJ/Kg) para este valor se tomó de tablas termodinámicas 2792.15 KJ/Kg. Ahora, h1 corresponde a la entalpia del agua entrante a T determinado por datos técnicos o precalentada expresada en (KJ/Kg) para este valor se tomó de tablas termodinámicas 440.13 KJ/Kg. En cuanto al gasto de vapor las variables que más lo afectan en forma directa son la temperatura de alimentación del agua, presión de trabajo y el incremento de temperatura del aire e n el precalentador, la primera con mayores efectos favorables. Las variables que afectan disminuyendo su valor son humedad de la madera y la temperatura de sobrecalentamiento del vapor (mv). Se especifica las libras de vapor producido con un valor de 15 toneladas por hora, se debe realizar la conversión a kg/s obteniendo como resultado 4.16kg/s = mv Se tiene mc que corresponde al consumo de combustible expresado en kilogramos por hora, se debe realizar una conversión sencilla para pasarlo a kg/s, para el ejercicio se emplea inicialmente los valores dados por el autor 3500kg/h correspondientes a 0.97kg/s, pero estos valores para efecto de simulación se dejan variables.
Fig. 23 Bloque cálculo de eficiencia (Autor Matlab-Simulink R2016a)
Fig. 24 Resultado eficiencia de la caldera (Autor Matlab-Simulink R2016a) Según la figura 24 se puede mostrar el resultado para el cálculo de valor de eficiencia de la caldera, según la simulación y los valores ingresados que corresponden a los especificados por el autor y a las condiciones de la ciudad de Bucaramanga se tiene como resultado una eficiencia del 71.65%
Fig. 25 Resultado eficiencia de la caldera (Autor Matlab-Simulink R2016a) Por último en la figura 25 se muestra el desarrollo de los bloques en Matlab-Simulink de las ecuaciones de cada uno de los subsistemas que se tuvieron en cuenta a la hora de realizar la simulación.
5. CONCLUSIONES La temperatura ambiente es poco sensible al proceso en general en cuanto la temperatura de sobrecalentamiento del vapor, es una variable del proceso que no genera un cambio de fase por lo que sólo afecta la calidad del vapor que se va a producir y por ende no genera cambio en las demás unidades del sistema. De igual forma la temperatura de alimentación del agua, ya que ésta no afecta el calor que producen los humos. Se consideraba que el exceso de aire sería la variable de mayor impacto en el funcionamiento de la caldera por lo que afecta directamente la combustión y de allí toda la producción de vapor pero no fue así, modificó algunos comportamientos del sistema pero no fue la más importante. La mayor eficiencia de la caldera se obtiene a una menor humedad del combustible, debido a que se modifica el poder calorífico superior del mismo, afectando de manera directa el resultado obtenido. El ejercicio de realizar el balance de energías es un tema complejo en cuanto a ajustes, debido a que intervienen diversos factores que alteran el resultado y que no se pueden controlar con facilidad, aunque el autor sugiere iterar con el valor de la temperatura hogar hay otros datos que pueden incidir directamente con las variaciones presentadas tales como el área de hogar, la emisividad o las pérdidas generadas por las paredes del ambiente. Los resultados obtenidos con el desarrollo matemático, difieren con los resultados obtenidos en la simulación por factores de decimales, lo que indica que tanto el desarrollo matemático como la simulación están correctos, y que no siempre la teoría concuerda con los datos obtenidos de ensayos reales. El software Matlab con su herramienta simulink, es una herramienta confiable y eficaz a la hora de realizar una simulación en la realización de este tipo de proyectos.
6. BIBLIOGRAFÍA [1] Y. A. Cengel and M.A Boles, Termodinamica.: McGraw Hill, 2006. [2] Fernando Alem Rodriguez, "Desarrollo de un simulador de calderas en MATLAB-SIMULINK," Memoria de trabajos de difusión cienttífica y técnica , no. 7, 2009. [3] Carlos Walter Thomaset, Pequeño manual del foguista , Décima ed. Montevideo, Uruguay, 2011.
