Instituto Politécnico Nacional

Escuela Superior de Ingeniería Química e Industrias

Extractivas
Departamento de Ingeniería Química Petrolera

Ingeniería de Reactores Heterogéneos

Proyecto del tercer periodo

Oxidación de O-xileno a Anhídrido Ftálico
INTEGRANTES:

- Granados Rojas Jorge Iván

- López Vázquez Bianca
- Sosa Perea Ana Gabriela
Grupo: 3PM61
Profesor: Dr. Edgar Ramírez Jiménez

07 Junio 2022

Contenido
Planteamiento del problema........................................................................................................................3
Ecuaciones.........................................................................................................................................................3
Estrategia de solución....................................................................................................................................5
Balance de materia.....................................................................................................................................6
Balance de energía.....................................................................................................................................7
Planteamiento del problema
Este proceso altamente exotérmico se lleva a cabo en reactores multitubulares
enfriados mediante sal fundida (molten salt) que se encuentra circulando
alrededor de los tubos y que intercambia calor hacia una caldera interna o
externa. La longitud de los tubos es de 3 m y su diámetro interno de 2.54 cm. Un
reactor puede contener 2500 tubos en paralelo y hasta 10000 o más en su
versión más grande. En procesos alemanes el catalizador es oxido de vanadio
V2O5 soportado en gel de sílice, y el intervalo de temperatura de operación varía
desde 335°C hasta 415°C. El diámetro de la partícula es de 3 mm, y la densidad
aparente es de 1300 kg/m3. El hidrocarburo es vaporizado y mezclado con aire
antes de entrar al rector. Su concentración se mantiene por debajo del 1 %mol,
con el fin de mantenerse debajo del límite de explosión.
La presión de operación es cercana a la atmosférica. La producción de anhídrido
ftálico en un reactor con 2500 tubos es de 1650 toneladas/año. De aquí se sabe
que este catalizador tiene un flux másico de la mezcla de gasta igual a 4684
kg/m2h. Con una densidad media del fluido de 1.293 kg/m 3, esto siguiendo una
velocidad superficial del fluido de 3600 m/h. El calor típico de reacción es de
1285409 kJ/kmol y el calor específico es de 0.992 kJ/kg K. En este problema, la
ecuación cinética para la conversión del hidrocarburo será considerada en una
primera aproximación a ser pseudo-primer orden, debido al gran exceso de
oxigeno:
0
rA=k p B p A

Donde p0B = 0.211 bar representa la presión parcial del oxígeno. Donde k está dado
por:
13636 −1 −2 −1
ln=19.837 − mol A g × atm ×ℎ
T

2
Ecuaciones
La ecuación de continuidad del hidrocarburo, en términos de presiones parciales y
de la ecuación de energía, puede escribirse, para una densidad constante:

dpA M m ptρb 0
vs + k p B p A =0 → ⟨ 1 ⟩
dz ρg
dT 4U
( T −T r )=0→ ⟨ 2 ⟩
0
v s ρgC p − ( −∆ H ) ρb k pB p A +
dz dt

Con
p A = p0 en z=0

T =T 0=T r en z=0

La presión total se considera constante e igual a 1 bar. El coeficiente global de

transferencia de calor U puede considerarse de 0.096 kJ/m 2sK. Se elige que Tr es
de 352°C.
Estas ecuaciones constituyen ecuaciones diferenciales ordinarias no lineales
acopladas de primer orden que pueden ser pueden resolverse convenientemente
mediante el procedimiento Runge-Kutta.

3
Estrategia de solución
Para la solución utilizaremos variables de referencia, las ecuaciones de balance
de materia y energía estarán en función de las siguientes variables para llevar a
cabo el modelo adimensional.

 Longitud adimensional
z
ξ=
L
 Temperatura adimensional
T
y=
T0

 Presiones parciales
Pa
x=
Pt

 Número de Arrhenius
Ea
γ=
RgT 0

 Numero de Damkohler
−γ
A 0 e ρ b M m Pt P0 L
Da=
vs ρ g

 Beta
−∆ H
β=
CP M mT 0

 Delta
4 UL
δ=
v s ρ gC p dt

Se conoce que el coeficiente cinético es función de la temperatura, también se

sabe que K0 es función T0, por lo tanto:
− Ea
Rg T 0
K 0= A 0 e

4
  Ecuación de Arrhenius:

( ) ( )
K − Ea 1 1
ln = −
K0 Rg T T 0

Se resuelve la ecuación anterior para encontrar K

( )→ ( )
− Ea 1 1 − Ea T 0 1 1

( ) ( )
− −
K Rg T T 0 K Rg T 0 T T 0
=e =e
K0 K0

( T − T )→ K =e −γ ( T − T )
1 1 T0 T0

( ) ( )
−γ T
K 0

=e
0 0

K0 K0

− γ ( − 1) (T T
T0 1 T0
)
( ) ( )
−γ −
K T K
=e → =e
0

K0 K0

−γ( ) (1 −y y )
1−y

( )
K −γ
y
=e → K=e K0
K0

−γ ( 1 −y y ) − Ea
R g T0
−γ ( 1 −y y ) −γ
K=e A0 e → K=e A0 e

Balance de materia
Sustituir en la ecuación 1 los números adimensionales propuestos:
Ec. 1

dpA M m ptρb 0
vs + k p B p A =0
dz ρg

Sustitución

dpA M m ptρb 0
vs + k p B p A =0
dz ρg
0
dx Pt M m P t ρb 0 d x Pt − M m Pt ρ B KP B P A
vs + K p B p A =0 → =
d ( ξL ) ρg d ( ξL ) ρ g vs
0 0
P t dx − M m Pt ρB KP B P A dx − M m Pt ρ B KP B P A L
= → =
Ldξ ρ g vs dξ ρ g v s Pt

Sustituyendo K

5
 0
−γ ( 1 −y y ) −γ 0
dx − M m Pt ρ B KP B P A L dx − M m P t ρB e A 0 e PB P A L
= → =
dξ ρ g v Pt dξ ρ g v s Pt

Para Damkohler:
−γ ( 1 −y y )
dx Da e PA
=
dξ Pt

Sustituyendo x:
−γ( )
1− y
dx y
=Dae x =0
dξ
Balance de energía
De ec. 2
dT 0 4U
v s ρgC p
dz
− ( −∆ H ) ρ K P B P A +
dt
( T −T r )=0

Se sabe que:
z
ξ= → z=ξ × L→ dz=dξ × L
L
T
y= →T = y ×T 0 → dT =dy × T 0
T0

dT 0 4U v s ρg C p T 0 dy 0 4U
v s ρgC p
dz
− ( −∆ H ) ρ K P B P A +
dt
( T −T r ) →
L dξ
− ( −∆ H ) ρb K P B P A +
dt
(T −T r )

vs ρ g C pT 0
Dividiendo entre :
L
v s ρg C p T 0 dy 0 4U
− ( − ∆ H ) ρb K P B P A + (T − T r) 0
L dξ dt dy − ∆ H ρb K P B P A L 4 UL
→ − + (T −T r )
v s ρg C p T 0 dξ v s ρg C p T 0 dt vs ρ g C pT 0
L
De:
PA
X A= → P A =PT × X A
PT

Sustituyendo:

6
 0
dy − ∆ H ρb K PB PT X A L 4 UL
− + (T − T r)
dξ vs ρ gC pT 0 dt v s ρ g C p T 0

Con Arrhenius se obtiene K0.

Se sabe que:

K= A 0 e e
−γ
−
( REaT − REaT )
g g 0

T
y= →T = y ×T 0
T0

Sustituyendo:

K= A 0 e e −γ
−
( REaT − REaT ) → K =A
g g 0
e −γ
e
−
( Ea
R yT
g
−
R T )
Ea
0 g 0
0

K= A 0 e e
−γ
−
( REaT − REaT ) → K =A
g g 0 −γ
e e
− ( γy − γ)
0

−γ
− ( γy −γ ) −γ
−γ ( 1y −1) −γ
−γ ( 1 −y y )
K= A 0 e e → K =A 0 e e → K = A0 e e

Sustituyendo en ecuación de balance de energía:

0 −γ
−γ (1 −y y ) 0
dy − ∆ H ρb K P PT X A L
B 4 UL dy − ∆ H ρb A 0 e e PB PT X A L 4 UL
− + (T − T r)→→ − + (T −T
dξ vs ρ gC pT 0 dt v s ρ g C p T 0 dξ vs ρ gC p T0 d t v s ρg C p T 0

Mm
Se multiplica por
Mm

( ( 1−y y )
)
−γ
−γ 0
Mm dy − ∆ H ρb A 0 e e PB PT X A L 4 UL
− + (T −T r )
Mm dξ v s ρg C p T 0 dt vs ρ gC pT 0

Desarrollar el número de Damkohler, el número β y δ:

Da:

7
 −γ
A 0 e ρ b M m Pt P0 L
Da=
vs ρ g

−γ ( 1y − 1) 0
dy −∆ H e PB PT X A L 4 UL
Da= − D a + (T −T r )
dξ v s ρg C p T 0 dt vs ρ gC p T 0

β:
−∆ H
β=
CP M mT 0

− γ ( −1)
1
dy 4 UL
− Da β e y X A + (T −T r )
dξ dt vs ρ g C p T 0

δ:
4 UL
δ=
v s ρ gC p dt

− γ ( −1)
1
dy
dξ
− Da β e
y
X A +δ
1
T0 ( )
(T −T r )

Como:

T Tr
y= → y r= →T r=T 0 × y r → T = y × T 0
T0 T0

Por lo tanto:
− γ ( −1) − γ ( −1)
1 1
dy
dξ
− Da β e
y
( )
X A +δ
1
T0
dy
(T −T r )→ − Da β e
dξ
y
X A+ δ
1
T0 ( )
(T 0 × y − T 0 × y r )

− γ ( −1) − γ ( −1)
1 1
dy
dξ
− Da β e
y
X A +δ( )
1
T0
dy
(T 0 × y −T 0 × y r )→ − D a β e
dξ
y
X A+ δ ( y − yr )

− γ ( −1)
1
dy
− D a β e y X A +δ ( y − y r )
dξ

Tabla de datos

8
 L=3m Cp= 0.992 kJ/kg K
0
Diámetro = 0.0254 m PB = 0.211 bar = 0.2082 atm

Nt = 2500 tubos Pt= 1 bar = 0.9869 atm

T = 335 °C a 415 °C U= 0.096 kJ/m2sK = 345.6 kJ/m2hK

NA= 4684 kg/m2h Tr=T0= 352 °C = 625.15 K

Ea= (13636)(8.314)= 113369.7040
ρb= 1300 kg/m3
KJ/Kmol
vs= 3600 m/h A0= exp
(19.837)
= 4.12x108 kmol /kg atm2h

∆H= 1285409 kJ/kmol R= 8.314 kJ/kmol K

Calculo de Mm:

(
M m= 0.211× 32
kg
kmol)(
+ 0.01× 106.16
kg
kmol
+ 0.779 ×28 )(
kg
kmol
=29.6256
kg
kmol )
Calculo de ρg:
kg
4684
2
NA m ℎ kg
ρ g= = =1.3011 3
vs m m
3600
ℎ
Calculo de 𝛶:
kJ
113369.7040
Ea kmol
γ= = =21.8123
RgT 0
( 8.314
kj
kmol K )
(625.15 K )

Calculo Da:
−γ
A 0 e ρ b M m Pt P0 L
Da=
vs ρ g

9
 Da=
( 4.12 x 10
kg atm ℎ )
kmol
8
2
(e −21.8123
( m )( kmol )
) 1300 kg3 29.6256 kg ( 0.9869 atm ) ( 0.2082 atm ) ( 3 m )
=0.7027
( 3600
m
ℎ )( kg
1.3011 3
m )
Beta:
−∆ H
β=
CP M mT 0

kJ
1285409
kmol
β= =69.9645
( 0.992 kgkJK )( 29.6256
kg
kmol )
(625.15 K)

Delta:
4 UL
δc=
vs ρ gC p dt

δc=
(
4 345.6
kJ
m ℎK
(3 m)
2 ) =35.1395
( 3600
m
ℎ )( kg
1.3011 3 0.992
m
kJ
kg K )(
( 0.0254 m ) )
Ya obtenido nuestro modelo adimensional se usara el método numérico Runge-
Kutta de 4to orden para resolver nuestro sistema de ecuaciones diferenciales
simultaneas no lineales por medio de un programa hecho en FORTRAN.
program Oxileno

implicit none

!Declaración de variables

real*8 k1,k2,k3,k4,c1,c2,c3,c4,Mm,Pt,Db,Dg,Vs,U,dt,Cp,Hr,L,T,Pa,T0,Tr,yR,ya,x,G,Da,D,B,yb,Pb,x0,H,xfi,y,Xa

integer::i,n

!Valores Constantes

!Valores Constantes

Mm=29.6330!Kg/kmol peso del gas

Pt=0.9869 !atmh

Db=1300!kg/m^3

10
Pb=0.2082 !atm

Dg=1.3011 !kg/m^3

Vs=3600!m/h

U=345.6!Kg/m^2*h*k

dt=0.0254!m

Cp=0.992!Kj/kg*k

Hr=1285409 !Kj/kmol

L=3!m

x0=0

xfi=3

write(*,*)'Ingresa valor para n'

read(*,*)n

Write(*,*)'Ingresa valor inicial de Ya'

read(*,*)y

Write(*,*)'Ingrese el valor incial de Presion de A'

read(*,*)Pa

Write(*,*)'Ingrese el valor inicial de T en Kelvin'

Read(*,*)T0

Write(*,*)"Ingrese el valor de T"

Read(*,*)T

Write(*,*)'Ingrese el valor de Tr'

read(*,*)Tr

!Declaracion de numero adimensionales

h=(xfi-x0)/n

Da=(Mm*Pt*Db)/(Dg*3600)*exp(19.837)*exp(-G)*L*Pb

G=13636/T0

B=Hr/(Cp*T0*Mm)

D=(4*U*L)/(dt*Vs*Dg*Cp)

yR=Tr/T0

11
Xa=Pa/Pt

y=T/T0

!Impresion de numeros adimensionales

Write (*,*) 'valor de Da',Da

Write(*,*)'el valor de D',D

Write(*,*)'el valor de B',B

Write(*,*)'el valor de G',G

!Calculo Runge-Kutta cuarto Orden

Do i=1,n

k1=-Ya*Da*exp(-G*((1-y)/y))

c1=B*Da*(exp(-G*((1-y)/y)))*ya-D*(y-yb)

k2=-(ya+(k1*h)/2)*Da*exp(-G*((1-(y+(h*c1)/2))/(y+(h*c1/2))))

c2=B*Da*(exp(-G*((1-(y+(h*c1)/2))/(y+(h*c1)/2))))*(ya+(h*k1)/2)-D*((y+(h*c1)/2)-yb)

k3=-(ya+(k2*h)/2)*Da*exp(-G*((1-(y+(h*c2)/2))/(y+(h*c2/2))))

c3=B*Da*(exp(-G*((1-(y+(h*c2)/2))/(y+(h*c2)/2))))*(ya+(h*k2)/2)-D*((y+(h*c2)/2)-yb)

k4=-(ya+(k3*h))*Da*exp(-G*((1-(y+(h*c3)))/(y+(h*c3))))

c4=B*Da*(exp(-G*((1-(y+(h*c3)))/(y+(h*c3)))))*(ya+(h*k3))-D*((y+(h*c3))-yb)

!Obtención de valores nuevos de x0,y0,z0

xa=xa+(h/6)*(k1+2*k2+2*k3+k4)

y=y+(h/6)*(c1+2*c2+2*c3+c4)

x=x+h

!Impresión de Resultados

T=T0*y

Pa=xa*Pt

write(*,*)'iteración numero',i

write(*,*)'el valor de k1 es',k1

write(*,*)'el valor de c1 es',c1

write(*,*)'el valor de k2 es',k2

12
 write(*,*)'el valor de c2 es',c2

write(*,*)'el valor de k3 es',k3

write(*,*)'el valor de c3 es',c3

write(*,*)'el valor de k4 es',k4

write(*,*)'el valor de c4 es',c4

Write(*,*)'el valor de la T adimensional es',y

write(*,*)'el valor de la presion adimensional es',Ya

write(*,*)'el valor de T',T

Write(*,*)'el valor de Pa es',T

Write(*,*)'el valor de Pa es', Pa

end do

end program Oxileno

13