Radial Effects in a Tubular Reactor

Introduction

Chemical reactors are at the core of chemical engineering processes. It is in these reactors that chemical reactions are made to take place. Chemical engineers are concerned about the design of the reactors seen that the desired output product should be produced with high efficiency.

In addition to continuous-stirred tank reactors and batch reactors, tubular reactors are also commonly used in the industry. Tubular reactors are made up of a cylindrical pipe. In the reactor, a variety of reactions can occur and many approaches can be taken in modeling a tubular reactor.

One aspect to consider when modeling tubular reactors is the flow profile in the reactor. Traditionally, to keep things simple, a plug-flow was assumed. A plug-flow assumes a constant velocity through the cross section of the tubular reactor. In an adiabatic case this has the benefit that no radial contributions, for example temperature, need to be considered and a 1D model is sufficient in this case.

In this notebook we will first show a simple 1D plug flow based model and then move on to a more elaborate 2D model.

While tubular reactors are typically used for gas phase reactions, liquid phase reactors are not uncommon either. This example will consider a liquid phase reactor.

In this tubular reactor a liquid phase, exothermic and irreversible reaction will take place first in a plug flow regime and then in a laminar flow regime. The reaction taking place inside the reactor is a hydrolysis reaction between propylene oxide product , and water product to form propylene glycol, product . This process can take place at room temperature when catalyzed by sulfuric acid. This kind of reaction is exothermic which can be approximated as a first-order reaction given that the reaction takes place in an excess of water [Fogler, 2016]. First-order reactions are reactions that entirely depend on the concentration, such that by doubling the concentration the product is also doubled.

Show the chemical reaction:

This notebook illustrates how to model an exothermic reaction in a tubular reactor and see the concentration of the species and the temperature distribution in the domain.

Load the finite element package:

1 Dimensional Adiabatic Plug Flow

A schematic of a classical tubular reactor model with a plug flow profile and reactants A and B and resultant C is shown. The rector is oriented along the -axis.

8.gif

In this first example we assume a plug flow and we assume that we have an adiabatic system. Assuming an adiabatic system means that no energy is flowing into or out of the system, say, through a cooling mechanism. These assumptions allow us to neglect all radial effects. As a consequence the model will be 1D.

Multiphysics Model Setup

Since this problem considers more than one physics domain, a multiphysics model is to be constructed. As a first step we look at the contributing physics domains independently. A heat transfer model and a mass transport model are set up to solve for the temperature field and concentration field within the reactor.

Initially, the flow field velocity within the reactor will be approximated by a plug-flow velocity profile in the -direction given by the following equation that represents the average velocity:

where is the total flow rate and is the radius of the reactor.

As a first step we define 1D steady-state heat transfer and a mass transport model.

Heat Transfer Model

The heat equation is used to solve for the temperature field in a heat transfer model:

where is the mass density in , is the specific heat capacity in , is the flow velocity in , is the thermal conductivity in , and is the heat source in .

Since the reaction is an exothermic process, the heat generated by the reaction is modeled by a heat source , which depends on both the chemical reaction rate, , and the heat of reaction, :

Set up a 1D steady-state heat transfer model:

Mass Transport Model

For a steady-state mass transport model the concentration field within the reactor is described by the following equation.

where denote the diffusion coefficient in , is the flow velocity, and is the reaction rate.

The chemical reaction rate depends on both the reactant concentration and the reaction rate constant, .

Note that for consumption reactions the reaction rate, , is negative.

The rate of reaction is temperature dependent according to the Arrhenius equation:

where is a frequency factor in , the activation energy in [J/mol], and the gas constant in .

Note that the standard notation of the reaction constant is denoted by , however, to avoid the conflict with thermal conductivity in this case, we use the non-standard symbol instead.

Set up a 1D steady-state mass transport model:

Parameters

Model variables

Besides the concentration variable , which gives us the concentration of the species , we can also compute the conversion of species and the concentration of species .

The conversion ratio of species is given by:
Define the concentration of species :

Material Parameters

The feed into the reactor consists of two streams. One stream is an equivolumetric mixture of propylene oxide and methanol, and the other stream is water containing 0.1 wt% (weight percent) of sulfuric acid. The molar flow rate of propylene oxide fed to the tubular reactor is 0.1 . The water is fed at a volumetric rate 3.5 times larger than the propylene oxide-methanol feed [Fogler, 2016].

To solve any of the proposed multiphysics models, we need to specify material properties of the reactants and the mixture.

First we define the density , heat capacity and molar weight of each of the species.

Set up the material properties of Propylene Oxide using as abbreviation:
Set up the material properties of Methanol using as abbreviation:
Set up the material properties of the fluid medium water:
Set up the frequency factor , the activation energy , and the heat of reaction for Propylene Oxide:
Set up the volumetric flow ratio:
Set up the molar and volumetric flow rate at the inlet:
Set up the total flow rate:
Set up the concentrations at the inlet:
Set up the thermal conductivity and the diffusion coefficient of the mixture using as abbreviation:
Set up the gas constant :
The density of the solution is defined as the weight sum using the initial concentrations of the species:
The heat capacity of the solution is defined as the weight sum using the concentration of and as a variable:
Set up the average velocity:
For further computation the magnitudes of the concentrations are extracted:

Model parameters

Specify the multiphysics model parameters.

Set up the reaction rate and the heat source term :
Set up the parameters to use:

Domain

The domain can be represented as a line that goes through the main axis of the tubular reactor from 0 to .

Specify the parameters of the geometry:
Create a 1D mesh of a Line primitive:

Boundary Conditions

Both, the temperature and the concentration are known at the inlet boundary .

Set up the HeatTemperatureCondition boundary condition at :

Solve the PDE Model

In the following section, the fully coupled multiphysics model will be solved.

Specify the multiphysics PDE with the model parameters.

The PDE model presented here is nonlinear. To solve nonlinear PDEs the solver needs an initial seed for the solver. The default initial seed is 0 for every dependent variable. Now, since the model uses a reaction rate that has a expression in it, a 0 initial seed for the temperature is not possible. The remedy for this is easy. Here we set the initial seed values as . This assumes there is no concentration of and the temperature of the reactor is that of the inlet.

Solve the 1D coupled PDE model and monitor the total time and memory used.

Verify that there are no radial effects.

Visualize the temperature distribution:
Visualize the temperature changes at the inlet, middle section, and outlet:

Extending the model

The second tubular reactor model presented here will additionally have a cooling jacket surrounding the reactor. The cooling jacket will introduce radial effects as the temperature varies from the cooling jacket to the center of the reactor. Since the concentration is also linked to the temperature the concentration will also see a radial effect.

To consider radial effects it is necessary to make use of a 2D model using cylindrical coordinates . This can be done because the model is rotationally symmetric about the -axis, so the cylindrical coordinate variable disappears.

Since radial effects will be present due to the cooling jacket a laminar flow profile can also be considered. The radial effects are actually further increased if a laminar flow profile is present. That is because a laminar flow field has a parabolic form. Considering the cooling jacket at the wall and a low flow velocity at the wall compared to the center of the cylinder the laminar flow profile will increase the radial effects. The model shown is based on [Fogler, 2016] and is usually known as a Non-adiabatic laminar flow model.

Multiphysics Model Setup

The flow field velocity within the reactor will be approximated by a laminar-flow velocity profile in the -direction given by the following equation, taking as an assumption that velocity in the -direction is zero.

Set up the laminar velocity profile:
Verify that represents a laminar flow profile:

Heat Transfer Model

For an axisymmetric steady-state heat transfer model the temperature distribution is described by the following equation.

The HeatTransferPDEComponent function can produce the axisymmetric form of the heat transfer equation. To do so the parameter "RegionSymmetry" is set to "Axisymmetric".

Set up a 2D axisymmetric steady-state heat transfer model:

Mass Transport Model

For an axisymmetric steady-state mass transport model the concentration field within the reactor is described by the following equation.

The MassTransportPDEComponent function can produce the axisymmetric form of the mass transport equation. To do so the parameter "RegionSymmetry" is set to "Axisymmetric".

Set up a 2D axisymmetric steady state mass transport model:

Parameters

We redefine the parameters and the model variables, because now they are and dependent.

The conversion of species is now given by:
Define the concentration of species :
Set up the reaction rate and the heat source term :
Set up the parameters to use:

Domain

The chemical reaction takes place inside a tubular reactor of 1 in length and 0.2 in diameter. The reactor has the property that it is rotationally symmetric about the -axis, so the geometry can be approximated by 2D rectangle which represents a cross-section of the reactor in the plane. The sketches from below show the simulation domain.

107.gif

In order to get a good result, a finer grid needs to be set up near the reactor inlet and outlet, and the reactor wall. This will be done by using ToGradedMesh.

Create two graded 1D meshes from a Line primitive:
Create the full 2D mesh:
Visualize part of the mesh:

Note, how the mesh is finer from top to bottom and from left to right.

Boundary Conditions

Since this multiphysics simulation contains a heat transfer model and a mass transport model, the boundary conditions for each physics mode are set up respectively.

Heat Transfer Boundary Conditions

There are four types of boundary conditions involved in the heat transfer model, one for each boundary.

At the inlet a HeatTemperatureCondition models an immediate temperature rise caused by the heat produced by the mixing of the two feed streams. The temperate at the inlet is set to 312 .

Set up the HeatTemperatureCondition boundary condition at :

At the reactor wall a convective boundary condition is in place. This models the heat exchange between the reactor and the cooling jacket. The heat transfer coefficient is 1300 and the temperature of the cooling jacket is 277 .

Set up the HeatTransferValue boundary condition at :

At the outlet an outflow boundary condition is applied and at the symmetry line a symmetric boundary condition is set up. Since both the outflow boundary condition and the symmetric boundary condition are Neumann zero conditions, they are implicitly applied without further setup.

Mass Transport Boundary Conditions

There are four types of boundary conditions involved in the mass transfer model, one for each boundary.

At the inlet the concentration of the specie is fixed at .

On the wall, , the mass flux of the species is assumed to be zero, and is modeled by an impermeable boundary condition.

An outflow boundary condition and a symmetric boundary condition are applied on the flow outlet and the symmetric axis respectively. However, since they are Neumann zero conditions, they are implicitly applied without further setup.

Solve the PDE Model

In the following section, the fully coupled multiphysics model will be solved.

Define the fully coupled heat transfer-mass transport model:

Here we set the initial seed values as . This assumes there is no concentration of and the temperature of the reactor is that of the cooling jacket.

Solve the fully coupled PDE model and monitor the total time and memory used:

Visualization

To see the effects of the temperature in the process of the reaction one can visualize the cross-section of the temperature profile and of the conversion profile in the reactor.

An alternative definition of the conversion of species which can be used for the visualizations:
Visualize the temperature and conversion distribution:

To visualize the full 3D solutions from the axisymmetric model one can apply the interpolating function to visualize the data in a 3D domain. In other words, do a revolution of the data.

Use RegionPlot3D to make a plot showing the three-dimensional region in which pred is the equation of a cylinder.

Visualize the concentration of species from the axisymmetric model in 3D using RegionPlot3D:

In order to see more details of the radial effects one can extract the changes from different locations along the length of the reactor: at the inlet, at the middle section, and at the outlet.

Extract the temperature and conversion changes from the three locations and visualize them in a single plot:

Nomenclature

References

1.  Fogler, H.S. Elements of Chemical Reaction Engineering, 5th Edition, Prentice-Hall Inc., New Jersey. (2016).