MODELO MATEMÁTICO DEL REACTOR DE SÍNTESIS DE AMONIACO DE LA EFNC Ramón Eduardo García García*, Iván Rodríguez Rico**, Leidy Morales Fundora***, Raúl Larramendi Cabrerizo***
*Universidad de Matanzas, **Universidad Central de Las Villas, ***Empresa de Fertilizantes Nitrogenados de Cienfuegos (EFNC) Se presenta un modelo matemático pseudohomogéneo que describe con buena precisión el comportamiento comportamiento en estado estacionario del reactor de síntesis de amoniaco de la Empresa de Fer tiliz ti liz ant es Nit rog ena dos de Cie nfueg nf ueg os, Cuba. Cu ba. Se bri nda det all es de la imple imp lemen menta tació ció n del mismo en el resolvedor de ecuaciones TK-Solver-Plus. Se incluye restricciones del modelo, su validación y los resultados de una corrida con datos reales, lo que permite el cálculo de la actividad del catalizador. Palabras clave: modelo matemático, reactores químicos, síntesis de amoniaco, TK-Sol ver-Plus. ___ ___ ___ ___ ___ ___ ___
It is presented prese nted a pseu doh omo gen eou s mat hemati hem ati cal mode l whi ch descr ibes ib es with good prec isi on the ste ady state sta te beha vio r of an ammoni amm oni a synthesi synth esi s rea cto r ins tal led at Cienf ueg os, Cuba Cub a Nitro gen ate d Fertil Fer til izer ize r Ent erpris erp rise. e. Th ey are ar e giv en de tai ls ab out i ts imple im ple mentat men tat ion o n the e qua tio n sol ver TKT KSolver-Plus. They are included model restrictions, restrictions, validation and a running with real data to calculate the catalyst activity. Key words: mathematical model, chemical reactors, ammonia synthesis, TK-Solver-Plus.
Introducción La síntesis de amoniaco ha sido muy estudiada con vistas a mejorar el diseño de nuevas plantas y la operación de las ya existentes. Debido a la gran capacidad y velocidad de cálculo de las computadoras modernas, se han desarrollado modelos matemáticos complejos para la generalidad de los tipos de reactores diseñados con este propósi pro pósi to. Elnasha Eln ashaie ie /3/ / 3/ hace un resum r esumen en de d e los lo s casos considerados en la literatura hasta 1987, pero ni en este est e perí odo ni hasta has ta la fecha fech a se ha obtenido información acerca de la modelación matemática de un reactor de síntesis de amoniaco con las características de éste en particular, que pertenece pert enece a la firma fir ma Helder He lder-Tho -Thompso mpsoe, e, Ingl I nglateaterra. Díaz /1/ organizó un sistema de ecuaciones que conforma un modelo homogéneo que fue implementado en computadora, y llegó a la conclusión de que el reactor no se podía describir con precisi prec isión ón a parti pa rtirr de este es te tipo ti po de model m odel o, lo cual coincide con el criterio emitido por Elnashaie /3/. Giger /5/ propone un modelo heterogéneo para describir el comportamiento dinámico de los reactores adiabáticos donde ocurren reacciones en
28
fase gaseosa catalizadas por sólidos (torres em pacadas pac adas), ), que descr iben ibe n la propag pro pagaci ación ón de ondas de temperatura y concentración (en la corriente gaseosa y en el sólido) a lo largo del equipo. El objetivo del modelo que aquí se desarrolla es describir los perfi les x, T en cada cama catalítica y los perfiles de t emperatura en el intercambiador de calor para diversas condiciones de operación, con la finalidad de determinar la zona de operación permisible en cada caso. Para describir con exactitud el comportamiento del reactor en estado estacionario Gaines /4/ y Elnashaie /3/ proponen utilizar un modelo pseudo homogéneo, que se basa en los principios de modelación del reactor como sistema homogéneo y la corrección de la velocidad de reacción con un término E denominado factor de efectividad de l a cama, que se calcula a partir de una ecuación empírica obtenida experimentalmen te en función de las condiciones de operación en cada punto de la misma. Morales /8/ realizó un estudio detallado acerca de las condiciones de operación del reactor, obtuvo datos de proceso, fundamentalmente, durante la arrancada de la planta y recogió los criterios que posee el personal vinculado a la operación de este. Dzingayi /2/ obtuvo un modelo
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
que describe el comportamiento en estado estacionario de la cama I, y llegó a la conclusión de que ésta define prácticamente el éxito de la operación posterior (cama II e intercambiador), ya que en ella se genera casi el 70 % del amoniaco y se establece los niveles de temperatura corres pondientes a una operación eficiente.
Desarrollo El modelo matemático que se describe a continuación se basa, fundamentalmente, en la aplicación a este reactor en particular de los resultados obtenidos por Gaines /4/y Elnashaie /3/, y extiende el modelo presentado por Dzingayi /2/ a todo el reactor.
Descripción del reactor La descripción detallada del reactor ha sido brindada por Dzingayi /2/ procedente de los manuales de operación de la planta instalada /7/. Por la importancia que esto tiene para la comprensión del modelo se incluye aquí una breve explicación: El convertidor de síntesis de amoniaco está instalado en una planta que tiene una capacidad de 700 toneladas de amoniaco por día. La planta cuenta con una sección de preparación del gas de síntesis, formado por H 2, N2, Ar y CH 4 básicamente, que entra como corriente de reposición a un lazo cerrado donde están dispuestos el reactor y un conjunto de equipos encargados de recuperar calor, separar parte del amoniaco procedente del
Esquema del reactor de síntesis de NH 3. Empresa de fertilizantes nitrogenados de Cienfuegos
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
29
reactor y preparar los gases restantes para su recirculación. Debido a las características pro pias de la reacción de síntesis, es conveniente extraer el exceso de inertes del lazo en forma de purga.
Estructura interna y circulación de gases El reactor está formado por dos camas catalíticas aisladas térmicamente, con blanqueo intermedio y un intercambiador de tubos y coraza situado en la zona inferior, colocados dentr o de una carcaza de presión aislada térmicamente; por fuera está la pared del reactor y entre ambas existe un espacio anular. El 70 % de los gases se alimenta por la entrada principal del reactor y desciende a través del espacio anular exterior con el objetivo de mantener refrigerada la pared exterior del reactor. Se ha comprobado en la práctica que el cambio de temperatura de estos gases es despreciable; entran al intercambiador de calor por el lado de la coraza, se calientan hasta la temperatura de reacción y ascienden por un tubo central aislado térmicamente hasta la parte superior del reactor, donde está la cama I; circulan a través de ella radialmente hacia fuera, se mezclan con el 30 % de los gases frescos restantes, se enfrían y se enriquecen en reactantes, luego atraviesan la cama II radialmente hacia adentro, descienden po r u n esp aci o an ul ar i n t eri or hast a el intercambiador de calor, circulan a través de los tubos, intercambian calor con la corriente de alimentación y salen del reactor.
La temperatura de entrada de los gases al reactor debe ser de 130 °C aproximadamente y la de salida de 300 °C, para garantizar la eficiencia del resto del lazo; por ello el reactor es de tipo autotérmico, o sea, el calor generado por la reacción química se utiliza para precalentar los gases de alimentación. Se debe notar el aislamiento casi perfecto con que cuenta este reactor, fundamentalmente, en el tubo central y la canasta, así como un sistema de distribución y recolección de gases que determina el éxito de la operación. La circulación radial establecida en las camas permite operar con flujos elevados sin caídas de presión apreciables y con pastillas de pequeño diámetro, que facilitan el contacto íntimo entre el sólido y los gases. La caída de presión a través de todo el reactor, excluyendo los efectos entrada-salida, es sólo de 2,5 kg/cm 2 (manométrica) cuando el reactor opera con una presión de entr ada de 254 kg/ cm 2 (manométrica).
Modelo matemático Consideraciones 1. 2. 3. 4.
Camas adiabáticas. Flujo radial uniforme. Tubo central con aislamiento perfecto. Los gases en la coraza (ánulo exterior) no absorben calor. 5. Caída de presión despreciable en todo el equi po.
Camas I y II
Reacción química única 2 N2(g) + 3 H2(g)
2 NH3(g)
Para describir el comportamiento de las camas I (superior) y II (inferior) se utiliza un sist ema de ecuaciones diferenciales ordinarias, obtenido a partir de balances de materiales y energía /6/, adaptados al caso de circulación radial, no reportado en la literatura consultada. Los perfiles x, T se calculan para un flujo de alimentación determinado, F.
30
∆ H 298 K = -10,96 kcal/mol de NH 3 formado π r c r H/F10)
(1)
π ∆Hr r c r H)/(F Cp(T,P)
(2)
dx/dr = (2
dT/dr = -(4/3
La velocidad de reacción se calcula a partir de la ecuación de Temkin-Pyzev modificada, aporta-
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
da por Gaines/4/, corregida a partir de los resultados de Walas /10/, que incluye un término Ka2 en el denominador del segundo bloque de actividades. Esto se comprobó mediante varias corridas preliminares. a 1−a 3 2 rc
=
3
2
I Pγ a I a M − G J G J M a K H Ka a K P H M P N Q
K 2 E Ka 2 a 2
1 2 3
3 2 3 1
Gases de blanqueo Gases a la cama II
(3) Gases a la cama II
Las expresiones de cálculo de Ka, K , las 2 actividades de los componentes involucrados en la reacción y E se han tomado con excelentes resultados del artículo de Gaines /4/. µ es un parámetro ajustable en el modelo que permite regular la velocidad de producción de amoniaco durante el proceso de validación del modelo. γ representa la actividad del catalizador 1 en condiciones de diseño (cargas catalíticas nuevas) y es también un parámetro ajustable en el modelo cuando se realiza corridas con datos reales ( γ <1). El cálculo de estos perfiles requiere que se plantee un balance de materiales en cada etapa de integración para calcular los flujos puntuales y la composición de todos los componentes de la mezcla reaccionante y el cálculo puntual de Cp, ∆Hr y de las actividades.
Punto de mezcla
Intercambiador de calor Los perfiles de temperatura en los tubos y la coraza se calculan a partir del siguiente sistema de ecuaciones diferenciales ordinarias: dT i dz
dT sh
Zona de quench
dz
Se realiza un balance de materiales donde se suma los flujos por componentes de las corrientes de salida de la cama I y quench, para obtener los flujos por componentes de la corriente de entrada a la cama II y un balance de energía donde se incluye las entalpias de las corrientes de entrada a la cama II y quench para el cálculo de la temperatura de entrada de los gases a la cama II. Para el cálculo de TeI I se obtuvo la expresión siguiente: Cp (TeI I, P) TeI I - B = 0
UA Wi Cpi
UA W sh Cpsh
( Ti
( Ti
− T sh )
− T sh )
(tubos) (6)
(coraza) (7)
U se calcula como: U =
hio hi
hio
+ ho + (hio ho Rd )
(8)
z=L
dz
(4)
donde:
=
=
z=0 B
=
+ (1 − m) Froc CprocT roc F si + (1 − m) F roc
F si Cp si
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
(5)
Elemento diferencial. Intercambiador
31
Rd es un tercer parámetro ajustable en el modelo. h io se calcula a partir de la ecuación de Colburn /4/. h o se calcula también a partir de la ecuación de Colburn eliminando la relación Di/Do, y utilizando el diámetro equivalente en la coraza.
misma fuente para las condiciones de alta presión con que se opera. µ m y K m se calculan puntualmente en función de la temperatura y la composición de los gases, tanto en los tubos como en la coraza.
Las propiedades físicas µ m, K m se obtienen con las ecuaciones generales de Perry /9/, corregidas a partir de los nomogramas generalizados de esa
Condiciones de contorno del modelo. Tablas 1y2
Tabla 1 Condiciones de contorno del modelo (Camas catalíticas)
r = r int Cama I Cama II
r = r ext
T = TeI, x = 0 T = TeII, x = 0
T = T sI, x = x sI T = TsII, x = xsII
Tabla 2 Condiciones de contorno del modelo (Intercambiador)
Fondo, intercambiador, z=0 Tope, intercambiador, z=L
Solución del modelo La no linealidad considerable del modelo pro puesto impide su solución analíticamente, por ello se recurrirá a la solución numérica utilizada por Díaz /1/, Gaines /4/ y otros. En todos los trabajos anteriores se ha elaborado softwares o programas específicos para la implementación en com putadora digital del modelo por resolver. Díaz /1/ reconoce las dificultad es de este método, que han sido comprobadas por los autores en un intento inicial de elaborar un programa con este fin en lenguaje Turbo Pascal 5.0. La presencia de ecuaciones de diferentes fuentes, empíricas en muchos casos, con gran variedad dimensional, requiere la comprobación minuciosa de cada una de ellas mediante la evaluación y comparación de los resultados obtenidos con la información numérica procedente de la literatura o la industria; esto
32
Tubos Tt = Ts Tt = TsII
Coraza Tsh = Trec Tsh = TeI
permite evit ar el solapamiento de los errores cuyas fuentes pueden ser difíciles de detectar. Existe una tendencia muy reciente de incluir una etapa previa a la programación que consiste en implementar el modelo en un software resolvedor de ecuaciones, lo cual resulta muy útil, y en ocasiones tal implementación puede constituir de hecho el software específico para la explotación del modelo. En este trabajo se implementa el modelo sobre el paquete profesional TK-Solver-Plus, que es un software resolvedor de ecuaciones basado en la pseu doprogramación. Entre l os conocidos en el país éste ha resul tado el más apropiado.
Esquema de cálculo 1. Cama I. 2. Zona de quench.
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
3. Cama II. 4. Intercambiador de calor.
Orden de cálculo 1. Conocida Trec se supone un valor TeI (si se dispone de datos industriales, T eI es conocido) y se calcula los perfiles x, T en la cama I. 2. Conocida la composición de la corriente quench (igual a la de los gases de alimentación) se calcula los flujos de entrada a la cama II. Se calcula, además T eII (temperatura de entrada a la cama II) a partir de la solución de la ecuación (4) aplicando el método de aproximaciones sucesivas. 3. Se calcula los perfiles x, T en la cama II. 4. Se calcula los perfiles T t, T sh en el intercam biador de calor.
5. Si T sh(z=0) - Trec >1 °C, se debe repetir el cálculo desde el punto 1 suponiendo otro valor de T eI o modificando los parámetros ajustables del modelo (en dependencia del propósito de la corrida), de acuerdo con la diferencia de tem peraturas obtenida en el pu nto 5. Nota: Los sistemas de ecuaci ones diferenciales planteados se resuelven mediante el procedimiento RK4TWO.tk, que es una implementación del método clásico Runge-Kutta 4 to orden. Se realiza la validación del modelo utilizando los datos de diseño (tabla 3). Como T eI es conocido en el punto 5 del orden de cálculo se ajusta el parámetro α para γ =1 (considerando cargas catalíticas nuevas). Se obtuvo una diferencia modular T sh(z=0)-T rec igual a 0,16 K. En esta corrida la temperatura de salida del reactor T t(z=0) difiere de la de diseño
Tabla 3 Datos de proceso Diseño
Operación*
243,70
143,32
5 839 709,50
2 897 205,80
Ts (K)
573,00
544,00
Te (K)
403,00
418,00
TeI (K)
673,00
685,00
15,42
13,02
P (atm, abs) Frec (lbmol/h)
NH3s (%V)
Composición a la entrada (% molar) H2
63,67
69,12
N2
21,23
20,45
NH3
3,10
4,53
3,19
2,61
8,81
3,29
Ar CH4
Parámetros del modelo Cama I
Cama II
Cama I
Cama II
α
1,239
1,300
1,239
1,300
γ
1,000
1,000
0,130
0,135
0
0
0
0
Rd
* Los datos de operación corresponden a valores medios provenientes de varias arrancadas (17-18 h), donde se manifiestan las mayores inestabilidades.
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
33
x 0,18 0,15 0,12 0,09 0,06 0,03 0 0,7
1,1
1,5
1,9
2,3
2,7
r(pie)
Fig. 1 Perfil de conversión. Cama I (diseño).
en 8,03 K, que representa un error relativo de 1,4 %, la diferencia de composición de amoniaco en la corriente de salida del reactor (0,147 8) con respecto a la de diseño (0,154 2) es 0,006 4, que representa un error relativo de 4,5 %, lo cual se considera válido dentro del grado de exactitud esperada para datos de instalaciones industriales, que en este caso se brinda como rangos y no como valores precisos. La desviación se atribuye, además, a las consideraciones del modelo establecidas, al conocimiento incompleto de la ecuación de velocidad de reacción y a la variedad de fuentes que se han utilizado para la construcción del modelo. Algunas ecuaciones (por ejemplo, la del factor de efectividad del catalizador) se han obtenido para condiciones no exactamente iguales a las de este equipo. Durante estas corridas se obtuvo los valores de α para cada cama (tabla 3).
34
En las figuras 1, 2, 3 y 4 se muestra los perfiles x, T en cada cama, donde se observa un comportamiento ascendente de x, T casi constante, lo que evidencia un buen dimensionamiento de las camas, aunque en la zona cercana a la salida se observa una pequeña disminución en las velocidades de crecimiento dx/dr, dT/dr, que se justifica en la práctica por el desplazamiento del equilibrio en sentido inverso, debido al exceso de amoniaco ya presente en las corrientes de gases y al nivel de temperatura alcanzado, que es termodinámicamente perjudicial para la verificación de la reacción directa, ya que ∆ Hr<0. En la figura 5 se muestra los perfiles T t, T sh en función de la altura del intercambiador. Se observa una diferencia de temperaturas adecuada en ambos cabezales y la no existencia de un cruce de temperaturas, que sí pudiera ocurrir si se opera fuera de los rangos permisibles de operación.
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
Fig. 2 Perfil de temperatura. Cama I (diseño).
Fig. 3 Perfil de conversión. Cama II (diseño).
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
35
Fig. 4 Perfil de temperatura. Cama II (diseño).
Fig. 5 Perfil de temperatura. Intercambiador (diseño).
36
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
Esto ha sido comprobado a través de sucesivas corridas del programa. Se realizó corridas con datos reales de proceso (T eI conocida) y se ajustó γ para cada cama (tabla 3). Los perfiles obtenidos son muy similares a los de la cama I. Se nota que el intercambiador no tiene la longitud suficiente para precalentar todo el flujo de gases, lo que resulta lógico ya que durante esta operación se alimenta una parte de los gases de entrada directamente a través del tubo central para regular la temperatura de entrada a la cama I (Cold-Shot). Esta distribución de flujo no ha sido considerada en el modelo. Des pués de estas corridas en el modelo se ha determinado todos los parámetros ajustables y por tanto, el mismo queda listo para su explotación industrial. El cálculo de γ se debe realizar cada cierto período de tiempo, ya que el catalizador pierde actividad durante la oper ación, lo cual permite además, determinar en qué momento se debe cambiar las cargas catalíticas.
Conclusiones y recomendaciones En este trabajo se ha presentado un modelo matemático que describe con precisión el com portami ento en estado estacionario del reactor de síntesis de amoniaco de la EFNC. El mismo permite calcular los perfiles x, T a través del reactor y estimar la actividad del catalizador. Se recomienda analizar la posibilidad de utilización del modelo en la industria o la elaboración de un software específico ejecutable a partir de los resultados obtenidos, ya que la implementación realizada por los autores requiere para su explotación el dominio del resolvedor de ecuaciones utilizado y el uso de un procesador matemático de alta velocidad.
Nomenclatura A a B Cp
area de la sección transversal de la cama (pie 2) actividad de los componentes de la mezcla de gases (atm) variable auxiliar calor específico de la mezcla de gases (BTU/ lbmol-K)
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000
∆ Hr variación de entalpia de la reacción (BTU/ lbmol NH3) D i diámetro interior de los tubos del intercambiador (pie) D o diámetro exterior de los tubos del intercambiador (pie) E factor de efectividad (0-1) F flujo molar (lbmol/h) H altura de la cama (pie) I cama I II cama II i componente i K conductividad térmica (BTU/h pie °F) Ka constante de equilibrio de la reacción directa (atm -2 ) K 2 constante específica de velocidad de la reacción inversa (lbmol H 2/h-pie 3) L longitud del intercambiador de calor (pie) m factor de distribución de gases (0,7) máx máximo mín mínimo P presión de trabajo (atm) Rd factor de obstrucción en el intercambiador (h-pie 3/BTU) r radio, posición radial en la cama (pie) rc velocidad de la reacción química (lbmol NH3/ h-pie 3) T temperatura de los gases (K) U coeficiente total de transferencia de calor (BTU/h-pie 3) W flujo molar superficial en el intercambiador (lbmol/h-pie2) x fracción convertida de hidrógeno (lbmol de hidrógeno convertido/lbmol de hidrógeno alimentado) z profundid ad dentro de la cama (pie) a parámetro cinético γ actividad del catalizador m viscosidad ( µ P)
Subíndices e ext int io m o
entrada exterior interior lado de los tubos (referido al área exterior) mezcla lado de la coraza
37
re c s sh t 0 1 2 3 4 5
recirculación salida coraza tubos inicial hidrógeno nitrógeno amoniaco argón metano
Bibligrafía 1. Díaz, I., Simulación de un reactor de síntesis de amoniaco, Tesis de Grado, Facultad de Ingeniería Química, Biblioteca EFNC, Santa Clara, 1975. 2. Dzingayi, E.; García, R. E. ;Rodríguez, I., Análisis de estabilidad del reactor de síntesis de amoniaco de la EFNC, Trabajo de Diploma, Facultad de Ingeniería Química, Universidad de Matanzas "Camilo Cienfuegos", 1994.
38
3. Elnashaie, S. S., "Digital Simulation of an Industrial Ammonia Reactor", en Chem. Eng. Sci., vol. XXXII, núm. 11, diciembre 1987. 4. Gaines, P., "Simulation of an Ammonia Synthesis Loop", en Chem. Eng. Sci., vol. XXX, núm. 8, octubre 1987. 5. Giger, G. K., "Control of Temperature Peaks in Adiabatic Fixed-Bed Tubular Reactors", en Ind. & Eng. Chem. Fund., vol. XIX, núm. 4, noviembre 1980, 389 págs. 6. Levenspiel, O., Chemical Reaction Engineering , La Habana, Editorial Pueblo y Educación, 1972. 7. Manual, "Operaciones. Planta Amoniaco", Biblioteca EFNC, 1970. 8. Morales, F. L.; García, R. E.; Rodríguez, I.; Larramendi, R., Estudio preliminar sobre la estabilidad operacional del convertidor de síntesis de amoniaco de la E.F.N.C., Trabajo de Diploma, Facultad de Química-Farmacia, Universidad Central de Las Villas, Santa Clara, 1992. 9. Perry, R.; Green, D., Chemical Engineers’ Handbook , 6ta edic., La Habana, Edición Revolucionaria, 1984. 10. Walas, S. M., Reacti on Ki net ic s for Chemical Engineers, New York, McGraw-Hill, 1982.
TECNOLOGÍA QUÍMICA Vol. XX, No. 2, 2000