Alexander Vervust steam cracking units Coupled furnace-reactor CFD simulations of industrial

143 Pages·2015·7.06 MB·English

Coupled furnace-reactor CFD simulations of industrial
steam cracking units
Alexander Vervust
Supervisor: Prof. dr. ir. Kevin Van Geem
Counsellor: Pieter Reyniers
Master's dissertation submitted in order to obtain the academic degree of
Master of Science in Chemical Engineering
Department of Chemical Engineering and Technical Chemistry
Chairman: Prof. dr. ir. Guy Marin
Faculty of Engineering and Architecture
Academic year 2014-2015

Coupled furnace-reactor CFD simulations of industrial
steam cracking units
Alexander Vervust
Supervisor: Prof. dr. ir. Kevin Van Geem
Counsellor: Pieter Reyniers
Master's dissertation submitted in order to obtain the academic degree of
Master of Science in Chemical Engineering
Department of Chemical Engineering and Technical Chemistry
Chairman: Prof. dr. ir. Guy Marin
Faculty of Engineering and Architecture
Academic year 2014-2015
Acknowledgements
Several people deserve a word of gratitude for their support during what has been an incredible journey
into the world of steam cracking and computational fluid dynamics.
I would like to thank my promotor prof. dr. ir. Kevin M. Van Geem for giving me the opportunity to
conduct this master thesis under his patronage and for his counsel and guidance throughout the past few
months. I am most grateful for his support in my endeavor to continue my research on steam cracking.
My gratitude also goes out to prof. dr. ir. Guy B. Marin for allowing me to commence my thesis as well
as continue my career at the Laboratory for Chemical Technology.
Even though he was not able to coach me all the way until the end, I want to express my thanks to Carl
for passing on his knowledge to me. He has made the goal of this master thesis crystal clear, which
allowed me to continue my work during his absence.
Special thanks go out to Pieter for taking on the role of coach during the second half of my thesis. Thanks
to his instructions I have been able to sharpen both my research and reporting skills. His help in
correcting my report has been invaluable. I would also like to thank David for sharing his insights on
steam cracker reactor modeling and Yu for answering my questions on radiation modeling.
For all the fun conversations and activities in and outside the LCT, I give thanks to my fellow student and
in particular Jessica, Jühl and colleague Jens who I shared an office with. Together we have endured the
many hardships caused by the close proximity of our office to the pilot and the lack of windows to
provide us with fresh air. They always succeeded at lightening the mood, especially Jühl who frequently
graced us with a wonderful song or a memorable one-liner. Jenoff also deserves a mention in this
acknowledgement for providing a meal when we were working late. But more importantly for his
friendship and support over the past five years.
I would not be writing this acknowledgement if it were not for my mother. She has always urged me to
aim as high as possible and I am happy to have followed her advise. Her support of my choice to become
a chemical engineer has helped me come this far.
Last but most certainly not least I would like to thank Eva. She has always been there for me when I
needed a morale boost or when I needed to vent my frustrations. I appreciate her for enduring the long
hours inherent to being a chemical engineer and for always greeting me with a smile when I come home.
Alexander Vervust
22 mei 2015
Deze pagina is niet beschikbaar omdat ze persoonsgegevens bevat.
Universiteitsbibliotheek Gent, 2021.
This page is not available because it contains personal information.
Ghent University, Library, 2021.
Coupled furnace-reactor CFD simulations of
industrial steam cracking units

Alexander Vervust

Master's dissertation submitted in order to obtain the academic degree of
Master of Science in Chemical Engineering
Academic year 2014-2015

Promotor: prof. dr. ir. Kevin M. Van Geem
Coach: Ir. Pieter A. Reyniers

GHENT UNIVERSITY
Faculty of Engineering and Architecture
Department of Chemical Engineering and Technical Chemistry
Laboratory for Chemical Technology
Chairman: Prof. dr. ir. Guy B. Marin

Abstract
Steam cracking of fossil feedstocks is the predominant process for the production of light olefins, such as
ethene, propene and 1,3-butadiene. As it is the single most energy-consuming process in chemical
industry, improving the thermal efficiency of the radiant section of a steam cracking unit can lead to
substantial economic and environmental benefits. An important tool in evaluating existing and new steam
cracker radiant section and reactor tube designs is computational fluid dynamics. Three-dimensional
simulations of both the radiant section and the reactors have been performed using the open-source CFD
code OpenFOAM®. Coupling these two simulations provides a powerful tool for the evaluation of
existing and new steam cracking technologies. A dedicated solver, called edmSimpleFoam, is developed
to model gas-fired furnaces and has been validated by means of the Sandia Flame D. The simulations
show that using the eddy dissipation concept alongside the discrete ordinates model and exponential wide
band model provides the highest accuracy. However as the eddy dissipation concept is sensitive to initial
conditions with respect to stability, eddy dissipation model simulations serve as initial guess for the
simulations with the eddy dissipation concept. The edmSimpleFoam solver is used to perform a full-scale
three-dimensional simulation of a Kellogg Millisecond propane cracking furnace coupled with
one-dimensional COILSIM1D reactor simulations. An in house developed solver is adopted for the
three-dimensional reactor simulations of a bare reactor and two swirl flow tube reactors. The mild swirl
flow tube reactor showed superior performance to the other reactors as the maximum external tube metal
temperature is 71 K lower at the cost of a pressure drop ratio of 1.36 compared to the bare reactor.
Keywords: Steam cracking, computational fluid dynamics, Kellogg Millisecond furnace, Sandia flame D,
swirl flow tube reactor, radiative heat transfer, OpenFOAM®
Coupled furnace-reactor CFD simulations of
industrial steam cracking units
Alexander Vervust
Promotor: prof. dr. ir. Kevin M. Van Geem
Coach: ir. Pieter A. Reyniers
Abstract: Steam cracking of fossil feedstocks is the Plehiers et al. 6 were amongst the first to simulate the radiant
predominant process for the production of light olefins, such as section in three-dimensions, followed by Nozawa et al. 7.
ethene, propene and 1,3-butadiene. As it is the single most Schietekat et al. 8, 9 performed three-dimensional simulations
energy-consuming process in chemical industry, improving the on enhanced reactor tube geometries, while Zhang et al. 10,
thermal efficiency of the radiant section of a steam cracking unit Hu et al. 11 and Stefanidis et al. 12 performed
can lead to substantial economic and environmental benefits. An
three-dimensional simulations of the radiant section coupled
important tool in evaluating existing and new steam cracker
with one-dimensional reactor simulations.
radiant section and reactor tube designs is computational fluid
The aim of this work is to perform three-dimensional
dynamics. Three-dimensional simulations of both the radiant
section and the reactors have been performed using the simulations of the radiant section of the steam cracking
open-source CFD code OpenFOAM®. Coupling these two furnace as well as three-dimensional simulations of the steam
simulations provides a powerful tool for the evaluation of cracking reactors using the open source computational fluid
existing and new steam cracking technologies. A dedicated dynamics code OpenFOAM®. Simulations of the Sandia
solver, called edmSimpleFoam, is developed to model gas-fired flame D serve as validation for the edmSimpleFoam solver.
furnaces and has been validated by means of the Sandia Flame
Full-scale three-dimensional simulations of a Kellogg
D. The simulations show that using the eddy dissipation concept
Millisecond furnace coupled with one-dimensional reactor
alongside the discrete ordinates model and exponential wide
simulations are performed as well as three-dimensional steam
band model provides the highest accuracy. However as the eddy
cracker reactor simulations of a bare reactor and two swirl
dissipation concept is sensitive to initial conditions with respect
to stability, eddy dissipation model simulations serve as initial flow tube reactors.
guess for the simulations with the eddy dissipation concept. The
edmSimpleFoam solver is used to perform a full-scale II. LITERATURE STUDY: FURNACE MODELING
three-dimensional simulation of a Kellogg Millisecond propane
Developing accurate predictive models for non-premixed
cracking furnace coupled with one-dimensional COILSIM1D
turbulent combustion proves to be a complex assignment.
reactor simulations. An in house developed solver is adopted for
the three-dimensional reactor simulations of a bare reactor and Two major challenges are accurately modeling the chemical
two swirl flow tube reactors. The mild swirl flow tube reactor source terms and the energy source term related to radiative
showed superior performance to the other reactors as the transfer. The latter involves both the solution of the radiative
maximum external tube metal temperature is 71 K lower at the transfer equation as well as calculating the radiative properties
cost of a pressure drop ratio of 1.36 compared to the bare of the participating medium. In what follows, several models
reactor.
available in literature for each of these source terms are
Keywords: Steam cracking, computational fluid dynamics, compared.
Kellogg Millisecond furnace, Sandia flame D, swirl flow tube The turbulence-chemistry interaction of the chemical source
reactor, radiative heat transfer, OpenFOAM®
term can be modeled using the eddy dissipation model, which
focusses on the fast and slow chemistry limit relative to
I. INTRODUCTION mixing, or the Eddy dissipation concept, based on the energy
Steam cracking of fossil feedstocks is the predominant cascade model, which is applicable in the full range of
process for the production of light olefins, such as ethene, chemistry time scales. The probability density function
propene and 1,3-butadiene. 1, 2 With a global annual method requires the solution of transport equations for the
production of 167 million tons in 2014, ethene is the most probability density function of the reactive scalars, shifting the
important steam cracker product with major applications in closure problem from the chemical source term to the closure
the polymer industry, e.g. for the production of polyethene models needed for the probability density function transport
and other derived polymers. 3, 4 equation. By calculating the average values of the reactive
Steam cracking is the single most energy-consuming scalars conditioned on a well-chosen conditioning variable,
process in chemical industry. 5 Improving the thermal the conditional moment closure method attempts to reduce the
efficiency of the radiant section can not only lead to an error on the chemical source term by effectively accounting
increased profit but it can also reduce the emissions of for turbulent fluctuations. The flamelet model decouples the
environmentally harmful gases such as CO and NO . An flame structure from the turbulent flow, allowing the use of a
2 x
important tool in evaluating existing and new steam cracker flamelet library to calculate averaged values of the reactive
radiant section and reactor tube designs is computational fluid scalars and reaction rate.
dynamics. Models used to solve the radiative transfer equation are the
Significant effort has already been dedicated to the P-1, Rosseland and discrete ordinates model, which transform
simulation of the furnace section of steam cracking units. the radiative transfer equation into a set of partial differential
equations to facilitate the computation of the spectral
intensity, while the discrete transfer radiation model solves the  u~
 k 0 (1)
radiative transfer equation to obtain a recurrence relation for
t x
use in between boundary conditions and along a set of k
predefined representative rays.
~ ~ ~
Two major groups of gas-phase radiative property models  u  u u  p   uu
l  k l   kl  F  k l (2)
exist, namely the gray and the non-gray models. The first t x x x k x
assume radiative properties to be wavelength independent k k k k
while the second account for changes in radiative properties
amloondgel tihs ea scpoemctmraol nrlyan ugsee.d Tghrea y wmeoigdhetle. dT hseu mno no-fg rgarya ym goadseelss Y~i u~kY~i  Jki R~ ukYi (3)
t x x i x
can again be subdivided in two groups, namely narrow band k k k
models and wide band models, depending on the spectral
rEelssoalsusteiro n moof detlh e amndo detlh. e Nasrtarotiwst icbaal ndm omdoedl,e ls,a ree .g.m othree h~u~kh~ Jkh ukhi  q (4)
t x x x rad
fundamental models, with an inherently high computational k k k
cost, while wide band models, e.g. the box model and the In these equations 𝜏 are the viscous stresses, modeled with
𝑘𝑙
exponential wide band model, are more suited for engineering the molecular viscosity model, and 𝐹 are the external forces
𝑘
applications as they provide accurate results with a moderate acting on the system, e.g. the gravitational force. The
computational cost. Reynolds stresses −〈𝜌𝑢′′𝑢′′〉 are modeled by a turbulence
𝑘 𝑙
Simulations on steam cracking furnaces performed by model. The terms 𝐽𝑖 are the molecular fluxes of species,
Stefanidis et al. 12, Hu et al. 11, Habibi et al. 13 and Zhang et al. 𝑘
which are calculated using Fick’s law for species diffusion,
10, show that a combination of the eddy dissipation model or and the molecular fluxes of enthalpy 𝐽ℎ are modeled with
eddy dissipation concept to model the turbulence chemistry 𝑘
Fourier’s law of heat conduction. The chemical source term 𝑅̃
interaction, the discrete ordinates model to solve the radiative 𝑖
is computed using a combustion model, which comprises a
transfer equation and the exponential wide band model to
reaction network and a model to handle the turbulence-
calculate the radiative properties of the absorbing and emitting
chemistry interaction. −〈𝑞̇ 〉 is the volumetric rate of
gas-phase species, is the most suited for the simulation of 𝑟𝑎𝑑
heating due to radiation, which is calculated by solving the
industrial scale steam cracking units.
radiative transfer equation. The terms 〈𝜌𝑢′′𝑌′′〉 and 〈𝜌𝑢′′ℎ′′〉
𝑘 𝑖 𝑘 𝑖
are approximated using the gradient diffusion assumption.
III. REACTOR AND FURNACE MODELS

2) Reactor model
A. One-dimensional reactor model
The reactor model uses a reaction network specifically
One-dimensional reactor simulations are performed using developed for propane cracking by reducing the full single-
COILSIM1D, a fundamental simulation model for calculating event microkinetic CRACKSIM model to its relevant core.
product yields of steam cracker reactors which has been This resulted in a reaction network of 203 reactions between
extensively validated with pilot plant experiments and 26 species, of which 13 radical species are modeled with the
industrial data. The two main parts of COILSIM1D are the pseudo steady state approximation. 9
reactor model and the single event reaction network. 14 The k-omega SST model is used to model turbulence in the
The reactor model comprises a momentum equation, species reactor model to accurately predict the near-wall flow in the
equation and an energy equation, to describe the wall-resolved simulations.
non-isothermal, non-adiabatic and non-isobaric character of
steam cracking. The equations are solved using a dedicated 3) Furnace model
solver developed by Dente et al. 15which adequately accounts A dedicated solver, called edmSimpleFoam, has been
for the inherent stiffness of the Jacobian of the set of developed to model industrial scale gas-fired furnaces.
discretized equations. Turbulence is modeled with the standard k-epsilon model.
The single event reaction network is based on the radical The implemented kinetic model for the combustion of
chemistry of the cracking process. It consists of a methane, developed by Westbrook and Dryer 16, consists of
monomolecular μ network and a mono- and bimolecular β the following two reaction steps.
network. The β network considers all reactions possible for
species with 5 carbon atoms or less. It contains 1324 CH + 1.5 O  CO + 2 H O
4 2 2
reversible elementary reactions between 51 molecules and 43
β radicals. The μ network adds another 676 molecules to the
network, resulting in a total of 770 components considered by CO + 0.5 O2  CO2
COILSIM1D. 14
The eddy dissipation model with finite rate approximation is
adopted to account for the effect of turbulent mixing on the
B. Three-dimensional models
reaction rates. In this model the reaction rate is determined by
1) Governing equations the minimum of the chemical reaction rate determined from
The Favre-averaged conservation equations of mass, the rate expressions and the rate of mixing, given by equation
momentum, species and energy are given in Einstein notation (5) 17, in which 𝐴 is an empirical constant with value 4.
by equations (1) to (4) respectively.
R  Aminy ,yox (5)
fu k  fu r 
 fu 
The discrete ordinates is used to solve the radiative transfer
equation. In this model the solid angles range of 4𝜋 is
subdivided into 𝑛 directions for which the radiative transfer
equation is solve. Consequently the radiative transfer equation
is approximated by a set of 𝑛 differential equations (6).
sˆ I (r,sˆ )( )I (r,sˆ ) I
i  i  s  i  b
 n (6)
 s w I (r,sˆ)(sˆ sˆ)
4 j  j j
j1

The radiative properties of the absorbing and emitting gas-
Figure 1. Temperature along the central axis of the burner predicted
phase species CO2 and H2O are calculated using the by the edmSimpleFoam solver and experimental data for the Sandia
implementation of the exponential wide band model flame D. gray / one-step. gray / two-step.
developed by Zhang et al. 10. This implementation accounts non-gray / one-step. non-gray / two-step,
for five absorption bands, two of which are overlapping and experimental data.
thus considered as one, and five transparent spectral windows
Species mass fraction plots in axial and radial direction
which are present in-between these absorption bands.
show an overprediction of the reaction rate resulting in an
overprediction of the temperature of up to 734 K close to the
IV. EDMSIMPLEFOAM VALIDATION: SANDIA FLAME D burner (x/d < 45), where the temperature is determined by the
rate of heat release due to the combustion reactions. A good
A. Geometry and operating conditions performance of the non-gray model is observed further away
The Sandia flame D is an often used validation case for from the burner (x/d > 45), where the temperature is mainly
radiation models. It is an unconfined CH /air flame, stabilized determined by the convective transfer of heat from the flame
4
on a piloted burner. The burner comprises a main jet and an and the absorbed and emitted radiation.
annular shaped pilot. The dimensions and operating The gray model predicts total absorption coefficients
conditions of the Sandia flame D are adopted from Barlow et between 0.4 m-1 and 3.5 m-1, while the more accurate
al. 18. non-gray model calculates total absorption coefficients in the
range 0-0.74 m-1. The gray model predicts the highest total
B. Solution procedure absorption coefficient in the hottest part of the flame, while
the total absorption coefficient field calculated by the
The calculations have been performed using the P1 model
together with the WSGGM, used by Lilleberg et al. 19, on the non-gray model shows a maximum in the zone around the
flame. Since the amount of radiation emitted by a gas is also
one hand (gray model) and using the discrete ordinates model
determined by the value of the total absorption coefficient, the
(DOM) together with the exponential wide band model
difference in calculated total absorption coefficient explains
(EWBM) on the other hand (non-gray model). Both
the 207 K difference in temperature prediction between the
simulations were carried out with the one-step reaction
mechanism used by Lilleberg et al. 19 on the one hand and the non-gray model and the gray model in the zone close to the
burner (x/d < 45), as seen in Figure 1.
two-step mechanism on the other hand, resulting in a total of
Comparison of the simulation results with those performed
four Sandia flame D simulations.
by Lilleberg et al 19 with the eddy dissipation concept shows
that the combination of a detailed reaction mechanism with
C. Simulation results and discussion
the eddy dissipation concept to model the
Figure 1 shows the temperature along the central axis of the
turbulence-chemistry interaction, the discrete ordinates model
burner predicted by the edmSimpleFoam solver as well as
to solve the radiative transfer equation and the exponential
experimental data. The edmSimpleFoam solver with the gray
wide band model to compute the radiative properties of the
model predicts the correct peak temperature value along the
gas would give the best simulation results. However as the
central axis of the burner, but not at the right location. The
eddy dissipation concept is sensitive to initial conditions with
edmSimpleFoam solver with the non-gray model overpredicts respect to stability 19, simulations with the eddy dissipation
the peak temperature by 221 K along the central axis of the
model should serve as initial guess for the simulations with
burner, however the maximum is located at the right location
the eddy dissipation concept. As neither the eddy dissipation
and the temperature is accurately predicted beyond the
model nor the eddy dissipation concept where implemented in
maximum. The same observations can be made for the OpenFOAM® 2.2 prior to this work, priority was given to the
temperature in the radial direction.
implementation of the eddy dissipation model.
V. KELLOGG MILLISECOND FURNACE: temperature profiles of each reactor, which serve as input for
3D FURNACE - 1D REACTOR COUPLED SIMULATIONS the furnace simulations. This procedure is repeated until the
difference between the previous and the updated maximum
A. Geometry and operating conditions tube metal temperature for every reactor is less than 1 K.
Due to the symmetry in the radiant section of the Kellogg
Millisecond furnace, it is sufficient to simulate a part of the
furnace, with a length of 0.328 m, width of 2.180 m and a
height of 10.556 m, comprising four reactor tubes with half a
burner at each side. Figure 2 (left) shows a schematic top view
of the simulated furnace floor. Located at the center of the
furnace are the reactor tubes. Their length matches the height
of the furnace and they pass through the furnace only once.
0.04365 kg/s of a propane/steam mixture, comprising 75.4
wt% propane and 24.6 wt% steam, is fed to each reactor at a
temperature of 903 K.
For reasons of computational cost, the burner design is
simplified: methane is supplied as fuel through a rectangular
shaped opening at the center of the burner with a flow rate of
0.0185 kg/s and air is supplied through the rectangular
opening surrounding the fuel inlet with a flow rate of
0.4308 kg/s. A general overview of the simulated furnace is
given in Figure 2 (right).

Figure 3. Flow chart of the solution algorithm for the coupled 3D
furnace – 1D reactor simulations of the Kellogg Millisecond furnace.
C. Simulation results and discussion
The results and discussion will be included in the final
version of the thesis.
VI. THREE-DIMENSIONAL STEAM CRACKER REACTOR
SIMULATIONS
A. Geometry and operating conditions
A bare reactor and two swirl flow tube reactors are
simulated. The dimensions of the reactors have been adapted
from the works of Schietekat et al. 8, 9. All reactor have a
length of 10.556 m, which matches the height of the Kellogg
Millisecond furnace. The swirl flow tube reactors have a
helical shape. The mild swirl flow tube reactor (SFT-M) has

an amplitude of 0.007 m and a pitch of 0.309 m and the high
Figure 2. Schematic view floor (left) and general overview (right) of swirl flow tube reactor (SFT-H) has an amplitude of 0.009 m
the simulated Kellogg millisecond furnace. and a pitch of 0.235 m. The reactor feedstock is propane with
a steam dilution of 0.326 kg/kg, the coil inlet temperature is
B. Solution procedure
903 K and the coil outlet pressure is 203 000 Pa, as is the case
The Kellogg Millisecond furnace was simulated in three for the reactors in the Kellogg Millisecond furnace. The heat
dimensions with the edmSimpleFoam solver. Figure 3 shows flux profile is determined based on one-dimensional
the solution algorithm used for the coupled 3D furnace – 1D CHEMKIN® simulations.
reactor simulations. The boundary condition for the reactors is
the least straightforward as this should be coupled with an B. Solution procedure
explicit reactor simulation in order to reach a realistic steady The solver developed by Schietekat et al. 9 is adopted for
state situation in the furnace. In the first furnace iteration, an
the simulation of the steam cracking reactors. In the bare
azimuthally averaged tube metal temperature profile as
reactor turbulence is model with the standard k-epsilon model
function of the axial position in the reactor is obtained based
and wall functions are used to describe the near-wall behavior
on 1D simulations in CHEMKIN. The boundary condition for
of the flow, while the k-omega SST model is used to model
the four reactor tubes are identical. Once the first iteration of
turbulence in the swirl flow tube simulations, which are
the furnace simulation has converged, axial heat flux profiles
wall-resolved. The external heat flux is corrected for the
for each reactor tube are obtained by averaging the calculated
difference in external surface area between the reactors, such
heat flux in azimuthal direction. These heat flux profiles are
that the total heat flux is the same for all simulations.
used as input for the one-dimensional reactor simulations,
with dimensions and operating conditions analogous to those
of the simulated steam cracking reactors. The reactor
simulations in turn give an update for the external tube metal