¡Descarga caracterizar un reactor y más Ejercicios en PDF de Matemáticas Aplicadas solo en Docsity!
Caracterización de la estática y la estabilidad de un reactor de tipo
tanque agitado con camisa de enfriamiento
Resumen
En el presente trabajo se estudió como es el comportamiento de un reactor, como caso específico
un reactor tipo tanque agitado con camisa de enfriamiento, se analizarán cuáles son los puntos de
equilibrio del sistema, así como su estabilidad. Se utilizaron principios de conservación, así como la
caracterización de la estática, también se emplearon algunos métodos numéricos y gráficos que
nos ayuden a comprender más lo que sucede en el sistema, además se determinó la estabilidad
del sistema mediante el criterio de Lyapunov. Al finalizar se logró obtener tres puntos de
equilibrio del sistema (x 1
, x
2
y x
3
), con un diagrama de Van Heerden se determinó la estabilidad de
los estados estacionarios del sistema, esto al graficar las curvas de calor generado y removido.
Posteriormente, mediante el criterio de Lyapunop, se logró determinar analíticamente la
estabilidad de los estados estacionarios y como resultado se obtuvo que el sistema cuenta con dos
estados estacionarios estables (x 1
y x
3
) y un punto inestable (x
2
). Este ultimo resultado
complementa al análisis que se realizo con el diagrama de Van Heerden.
Estados estacionarios de un reactor
Se muestra un reactor tipo tanque agitado con camisa de enfriamiento (figura 1) en el cual se
producirán ciertas toneladas de producto anualmente. Se necesita conocer algunas condiciones de
operación (concentración y temperatura principalmente) para que el reactor opere de manera
adecuada, estos son estados estacionarios o puntos de equilibrio.
Figura 1. Reactor químico tipo tanque agitado con camisa de enfriamiento
Con respecto a la figura 1 se realizaron los siguientes balances dinámicos de materia y energía
correspondientes a la masa total, masa reactiva, calor del reactor y calor en la camisa.
Balance de materia total
d
d t
a
[ Vρ ] =q
e
ρ
e
−qρ
(1)
Balance de materia en el tanque (reactante)
d
d t
a
[ CV ] =q
e
C
e
−qC−VR(C ,T )
(2)
Balance de calor en el tanque
d
d t
a
[ VρCpT ] =q
e
ρ
e
Cp
e
T
e
−q ρCpT −UAr
T −T
e
−Δ H
r
VR (C , T )
(3)
Balance de calor en la camisa de enfriamiento
d
d t
a
[
V
c
ρ
c
Cp
c
T
]
=q
c
ρ
c
CpT −q
c
ρ
c
CpT
c
+UAr (T −T
c
(4)
Se presenta la siguiente nomenclatura (tabla 1) utilizada en las ecuaciones (1), (2), (3) y (4).
Tabla 1._ Nomenclatura
Símbolo Significado Símbolo Significado
t
a
Tiempo actual Cp Calor especifico
V Volumen Cp
e
Calor especifico de entrada
Ρ Densidad Cp
c
Calor especifico en la camisa
ρ
e
Densidad de entrada T Temperatura
q
e
Flujo de entrada
T
e
Temperatura de entrada
q Flujo de salida
T
c
Temperatura de camisa
C Concentración Δ
H
r
Delta H de reacción
C
e
Concentración de entrada U Coef. Global de
transferencia de calor
Ar Área transversal
R
Función Arrhenius
De acuerdo a las características del reactor se presenta una lista de suposiciones:
El flujo en la entrada y salida del reactor es constante.
El flujo en la entrada y salida de la camisa de enfriamiento es constante.
Se considera un mezclado perfecto.
Consideramos un líquido incompresible, por lo tanto, la densidad permanece constante.
El calor específico es constante en el tanque.
En la tabla 2 se enlistan las variables y parámetros presentes en el reactor.
Tabla 4. números adimensionales identificados en los balances de materia y energía del reactor
Numero
Adimensional
Variables
involucradas
Interpretación Física
Stanton
UA
VρCp
Tiempo de residencia que alude a un tiempo característico de
intercambio de calor por la pared del reactor
Damköhler
k
o
V
q
r
Es el cociente de dos tiempos característicos, uno relacionado
con la convección de masa y el otro a la rapidez de reacción.
ϵ
E
a
R
g
T
r
Tiempo característico de la cinética de Reacción
Los balances dinámicos se pueden escribir en forma vectorial de tal manera que:
f
c , τ
[
f
c
(c , τ)
f
τ
(c , τ )
]
[
θ ( c
e
−c )−cα ( τ)
θ
τ
e
−τ
−υ
τ−τ
e
]
(7.a1)
Dónde:
f
c
( c , τ )=
dc
dt
f
τ
( c , τ )=
dτ
dt
, c (0)=
C
0
y τ (0)=
τ
0
Considerando la estática del reactor, con f ( c , τ ) = 0 , y considerando que:
θ=υ=c
e
, y que
τ
e
=τ
c
, también
a= 25 , ϵ= 50 (8.a1)
Esto es que la derivada de los balances dinámicos (ecuación 5.a1 y 6.a1) es igual a cero,
entonces se puede despejar la concentración adimensional (por conveniencia) y obtener
dos graficas correspondientes a las ecuaciones (2.d1) y (3.d1). Para ello se utilizaron los
datos de (8.a1). Entonces:
c
1
θ c
e
θ+e
a−ϵ /τ
(9.a1)
c
2
υ
τ−τ
e
−θ
τ
e
−τ
e
a−ϵ / τ
(10.a1)
El único dato faltante es τ por lo que se optó por evaluarlo, en un intervalo de 7/4 a 11/4,
esto ya que el punto mínimo es
τ
e
, por lo que el punto máximo es el mínimo más el
incremento adiabático (en este caso es 1) por lo que
τ
max
Con base en los cálculos realizados se obtuvo la gráfica de temperatura contra
concentración (grafica 1).
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
2
c
c
concentración (c)
Temperatura (τ)
Grafica 1 ._ temperatura (τ) contra concentración (
c
1
y c
2
)
La grafica 1 se puede localizar los puntos de intersección, señalados con un respectivo marcador.
Esto nos indica que el reactor puede operar de tres formas distintas, por lo que los puntos
solución son:
x
1
[
]
,
x
2
[
]
y
x
3
[
]
Dichos puntos representan que el reactor puede operar de tres maneras distintas, esto es que
presenta multiplicidad.
Método de Newton-Raphson
En el método de Newton-Raphson para una ecuación se utiliza la derivada de una función para
calcular las intersecciones con el eje de la variable independiente, que en otras palabras son las
raíces de la función. El método está basado en la expansión de la serie de Taylor de primer orden.
f
x
i + 1
=f
x
i
x
i + 1
−x
i
f ' (x
i
(11.a1)
Dónde:
x
i
=es el valor inicial de la raiz ( propuesto)
x
i+ 1
=es el valor en el cuálla recta tangente intersecta al eje x
Parte del método de consiste en igualar la ecuación (11.a1) a cero, con lo cual los puntos donde la
función se convierta en cero serán los puntos solución, de este modo la ecuación (11.a1) se
manipulo de la siguiente forma:
x
1
=π r
2
decir que se ha encontrado la solución, esto es cuando se calcula el error del punto iterado por
medio de la ecuación 12.a
|
x
i+ 1
− x
i
|
|
x
i
|
< 1 x 10
− 4
(12.a2)
De esta manera se calcula el error del punto iterado, si el valor obtenido es menor que el lado
derecho de la igualdad se encontró una solución del sistema, de lo contrario se continúa iterando
con el nuevo punto resultante hasta satisfacer la ecuación 12.a2.
Regresando al problema del reactor, debido a la complejidad que presentan las ecuaciones
resultantes de los balances de materia y energía, se utilizó este método para determinar los
estados estacionarios de una manera mas exacta, aprovechando que de manera gráfica se
obtuvieron los valores aproximados de las soluciones y que este método funciona de mejor
manera cuando se comienza a iterar cerca de los puntos solución.
Se encontraron 3 puntos solución, tal y como se había visto en el método gráfico, sin embargo, con
una aproximación más exacta de la solución. A continuación, se muestran los puntos obtenidos
por medio del método de Newton Raphson modificado.
x
1
[
]
,
x
2
[
]
y
x
3
[
]
Estáti ca y dinámica (estabilidad) del reactor
Los reactores químicos son unidades que presentan frecuentemente grandes dificultades de
control, en particular si las reacciones son rápidas y exotérmicas, pudiendo ser peligrosas y
potencialmente explosivas. Que el fenómeno sea controlable depende de la velocidad de
producción y de eliminación de calor del sistema. Si la velocidad de producción es superior a la de
eliminación el reactor se irá calentando, pudiendo llegar a fundir y/o producirse una explosión.
Consecuentemente se procederá a estudiar como es el calor producido sobre el eliminado en cada
uno de los posibles puntos estacionarios para determinar su estabilidad (diagrama de Van
Heerden ).
Si utilizamos la ecuación (9.a1) y sustituimos en (6.1) igualada a cero obtenemos que
0 =θ
τ
e
−υ
τ−τ
e
θ c
e
α (τ )
θ+α (τ)
(13.a1)
Haciendo un reacomodo de la ecuación (13.a1) obtenemos:
[
( θ+υ ) τ−(θ τ
e
+υ τ
c
]
θ c
e
α (τ )
θ+α (τ )
(13.a2)
En la ecuación (13.a2) la parte encerrada en corchetes cuadrados se puede ver como una recta
(curva de remoción de calor ƞ(τ)), donde el término dependiente de τ es la pendiente y el segundo
término es la ordenada al origen. Si graficamos esta recta con los valores de τ con los valores
previamente seleccionados contra la función no lineal que esta fuera del corchete de la ecuación
13.a2 (curva de generación de calor Ω(τ)) se obtiene la gráfica 2.
1.75 1.8 1.85 1.9 1.95 2 2.05 2.1 2.15 2.2 2.25 2.
0
1
Ω(τ)
ƞ(τ)
ƞ(τ), Ω(τ)
1.75 1.8 1.85 1.9 1.95 2 2.05 2.1 2.15 2.2 2.25 2.
0
1
c(τ)
τ
c
Grafica 2. Grafica de τ vs. Ω (τ), ƞ(τ) (diagrama de Van Heerden), comparada con la conversión de la
concentración dependiente de τ
Por medio de la estática se encontraron tres puntos donde el reactor puede operar, esto
al igualar las dos ecuaciones (para dejar en términos de una variable) y graficar la parte
lineal y la no lineal contra la variable de la temperatura en el diagrama de Van Heerden
(gráfica 2), se puede observar que existe multiplicidad en el reactor.
Al modificar algunas condiciones de operación se notó que el punto de equilibrio de en
medio no tiende a regresar a su estado original, como lo hacen el punto caliente (x 1
) y frio
(x 2
), lo cual nos indica que el punto de equilibrio es inestable. Por otro lado, se encontró
que al moverse por las vecindades de x 1,
al aumentar τ, la curva de remoción de calor es
mayor a la de generación, por lo que este punto tiende a regresar a x 1
. Ahora al disminuir
De esta manera, podemos pensar que el sistema (1.c) se encuentra próximo al sistema
lineal,
(
x
y ˙
)
= A
(
x−x
0
y− y
0
)
(4.c)
Cuando
(x , y ) está cerca de
(x
0
, y
0
, y por consiguiente, es natural esperar que el
comportamiento de las trayectorias del sistema (1.c) cerca del punto crítico
(x
0
, y
0
, sea
similar al de las trayectorias del sistema linealizado (4.c). En general, el punto de
equilibrio (x 0
, y
0
) del sistema autónomo (1.c) hereda la estabilidad, y en algunos casos la
naturaleza, del punto de equilibrio
( x
0
, y
0
para el sistema lineal (4.c).
Teorema 4.1 Teorema de Linealización de Lyapunov y Poincaré.
- El punto crítico( x
0
, y
0
) del sistema (1.c) es asintóticamente estable, si y sólo si,
todos los valores propios de la matriz A poseen parte real negativa
- El punto crítico
( x
0
, y
0
del sistema (1.c) es inestable, si y sólo si, la matriz A del
sistema posee un valor propio con parte real positiva
Más aún, si los valores propios de A son distintos entre sí y distintos de cero se puede
decir lo siguiente:
Si λ1 < λ2 < 0, entonces (x
0
, y
0
)es un nodo asintóticamente estable.
Si λ1 > λ2 > 0, entonces
(x
0
, y
0
es un nodo inestable.
Si λ1 < 0 < λ2, entonces (x
0
, y
0
) es un punto de silla.
Si λ1 no es real y Re (λ1) < 0, entonces
(x
0
, y
0
es un foco asintóticamente estable.
Si λ1 no es real y Re (λ1) > 0, entonces (x
0
, y
0
) es un foco inestable
Estabilidad de los estados estacionarios del reactor de
tanque agitado con camisa de enfriamiento basado en
el criterio de Lyapunov
Basándonos en la representación de Taylor, y en los criterios de la estabilidad de estados
estacionario de Lyapunov para sistemas dinámicos no lineales, se procede a analizar la
estabilidad del reactor con camisa de enfriamiento controlado. Se presentan los balances
de materia y energía del reactor controlado.
c ˙=θ
c
e
−c
−c α ( τ )=f
c
(3.c1)
τ =θ
τ
e
−τ
−υ(τ−τ
e
)+cα ( τ )=φ
τ
(3.c2)
Donde
α ( τ )=e
a−
ϵ
τ
De acuerdo con la representación de Taylor (2.c1) en las localidades del estado
estacionario (c , τ) y al considerar el truncamiento lineal obtenemos:
c ˙=
(c−c¿)=f
c
( c , τ )+
(
∂ f
c
∂ c
)
(c ,τ )
( c−c ) +
(
∂ f
c
∂ τ
)
(c , τ )
( τ −τ ) +… … …+ 0 ¿
(3.c3)
τ =
(τ−τ ¿)=φ
τ
( c , τ )+
(
∂ φ
τ
∂ c
)
(c , τ )
( c−c) +
(
∂ φ
τ
∂ τ
)
(c , τ)
( τ−τ ) +… … …+ 0 ¿
(3.c4)
El primer término de las ecuaciones (3.c3) y (3.c4) es equivalente a cero, de tal modo que
al ordenar la matriz Jacobiana (matriz de derivadas) se obtiene que:
A=
(
(
∂ f
c
∂ c
)
( c , τ
)
(
∂ f
c
∂ τ
)
( c , τ
)
(
∂ φ
τ
∂ c
)
(c , τ)
(
∂ φ
τ
∂ τ
)
(c , τ)
)
(4.c1)
Dónde:
(
∂ f
c
∂ c
)
=−θ−e
a−
ϵ
τ
… … … … … .. a
11
(
∂ f
c
∂ τ
)
−c ϵe
a−
ϵ
τ
τ
2
… … … … … … a
12
(
∂ φ
τ
∂ c
)
=e
a−
ϵ
τ
… … … … … … … … a
21
(
∂ φ
τ
∂ τ
)
=−(θ +ν ) +
c ϵe
a−
ϵ
τ
τ
2
… … … … … … a
22
De la matriz A se obtendrán los valores propios para poder analizar la estabilidad por
medio del criterio de Lyapunov, esto es:
( A−λI )=
(
a
11
−λ a
12
a
21
a
22
−λ
)
(4.c2)
Al aplicar el determinante a (4.c2) y desarrollar términos se llega a obtener el polinomio
característico de donde se obtendrán los valores propios para conocer si es estable en el
punto o no de acuerdo con el criterio de Lyapunov.
P
λ
=λ
2
[
a
11
22
]
λ+(a
11
a
22
−a
12
a
21
) (4.c3)
S
n− 3
c
1
c
2
c
3
c
4
c
n
Donde los coeficientes
b
1
, b
2
, b
3
, … … , c
1
, c
2
, c
3
, … … , d
1
, d
2
, etc., se evalúan de la
siguiente manera:
El criterio de estabilidad de Routh-Hurwitz plantea que el número de raíces de la ecuación
con partes reales positivas es igual al número de cambios de signo de los coeficientes de la
primera columna del arreglo. La condición necesaria y suficiente para que todas las raíces
de la ecuación se encuentren en el semiplano izquierdo del plano s es que todos los
coeficientes de la ecuación sean positivos y que todos los términos de la primera columna
del arreglo tengan signo positivo.
Anexo I
Suponiendo que V es constante y que el T c
es un dato conocido, además tomando en cuenta las
suposiciones antes mencionadas, las ecuaciones (2) y (3) se pueden expresar de la siguiente
manera (al sacar las constantes de la derivada):
dC
d t
a
q C
e
−qC−VR(C ,T )
V
(2.a1)
dT
d t
a
qρCp T
e
−qρCpT −UAr
T −T
e
− Δ H
r
VR(C , T )
VρCp
(3.a1)
Eliminando términos comunes y factorizando algunos términos obtenemos:
dC
d t
a
q
V
C
e
−C
−R(C ,T )
(2.a2)
dT
d t
a
q
V
T
e
−T
UAr
T −T
e
VρCp
−∆ H
r
R(C , T )
ρCp
(3.a2)
Donde en ambas ecuaciones (2.a2) y (3.a2) el término R=(k
0
e
− E
a
/ R
g
T
)C , que es la función
Arrhenius, aquí la
E
a
=energía de activación ,
R
g
=cte. de los gases , la función se sustituirá en
las mismas ecuaciones obteniendo:
dC
d t
a
q
V
C
e
−C
−(k
0
e
−E
a
/ R
g
T
)C
(2.a3)
dT
d t
a
q
V
T
e
−T
UAr
T −T
e
VρCp
−∆ H
r
(k
0
e
−E a
/ R g
T
)C
ρCp
(3.a3)
Las ecuaciones (2.a3) y (3.a3) se adimensionaran para pasar de una ecuación compleja a una más
sencilla, utilizando las siguientes variables adimensionales:
t=
t
a
t
r
,
c=
C
C
r
,
c
e
C
e
C
r
,
τ =
T
T
r
,
τ
e
T
e
T
r
, θ=¿ , α ( τ )=¿
Donde las variables adimensionales θ y α (τ) se definirán más adelante. De tal modo que
despejando las variables originales de las variables adimensionales y sustituyendo en las
ecuaciones (2.a3) y (3.a3) obtenemos:
dc C
r
dt t
r
q
V
c
e
C
r
−c C
r
−( k
0
e
−E
a
/ R
g
τ T
r
)c C
r
(2.b1)
dτ T
r
dt t
r
q
V
τ
e
T
r
−τ T
r
UAr
τT
r
−τ
e
T
r
VρCp
−∆ H
r
(k
0
e
−E
a
/ R
g
τ T
r
)c C
r
ρCp
(3.b1)
En las ecuaciones (2.b1) y (3.b1) se considera a
C
r
,T
r
, t
r
como constantes por lo cual se
pueden despejar de las ecuaciones, por lo tanto:
dc
dt
q
V
c
e
C
r
−c C
r
−(k
0
e
− E a
/ R g
τ T r
) c C
r
)t
r
C
r
(2.b2)
dτ
dt
q
V
τ
e
T
r
−τ T
r
t
r
T
r
UAr ( τT
r
−τ
e
T
r
) t
r
VρCp T
r
(−∆ H
r
) t
r
(k
0
e
− E
a
/ R
g
τ T
r
) c C
r
ρCp T
r
(3.b2)
Ahora se pueden eliminar algunos términos en las ecuaciones (2.b2) y (3.b2) obteniendo
las siguientes expresiones:
dc
dt
q t
r
V
( c
e
−c )−(k
0
e
− E a
/ R g
τ T r
) t
r
c
(2.b3)
dτ
dt
q t
r
V
τ
e
−τ
UAr t
r
τ −τ
e
VρCp
−∆ H
r
t
r
(k
0
e
−E
a
/ R
g
τ T
r
) c C
r
ρCp T
r
(3.b3)
Ahora bien en la ecuación (3.c2) se agruparan los términos
UAr
ρCp q
r
como upsilon (υ), el cual
tendrá que ver como un tiempo de transporte, y
−∆ H
r
C
r
ρCp T
r
como beta (β), el cual hace
referencia al incremento adiabático de temperatura (en este caso su valor es 1), por lo tanto:
dτ
dt
=θ
τ
e
−τ
−υ
τ−τ
e
ln Da−ϵ / τ
) c (3.c3)
Si el logaritmo natural del número de Damköhler (ln (Da)) lo representamos con una
constante (a), entonces las ecuaciones 2.c2 y 3.c3 se expresan como:
dc
dt
=θ
c
e
−c
−(e
a−ϵ / τ
)c (2.d1)
dτ
dt
=θ
τ
e
−τ
−υ
τ−τ
e
+(e
a−ϵ / τ
) c (3.d1)
Por comodidad la expresión e
a−ϵ /τ
se denotara como alpha (α (τ)), esta porque hace alusión
hacia a la función Arrhenius, la cual tienen que ver con la cinética de reacción respecto a la
temperatura, entonces:
dc
dt
=θ
c
e
−c
−c α (τ ) (5.a1)
dτ
dt
=θ
τ
e
−τ
−υ
τ−τ
e
+cα ( τ) (6.a1)
De esta manera las ecuaciones (2.a3) y (3.a3) pasaron a ser adimensionadas en las
ecuaciones (5.a1) y (6.a1)
