Microscale Simulation of Catalyst Deactivation
|Introduction||Solve the PDE Model|
|Mass Transport Model||Appendix|
|Initial and Boundary Conditions|
A catalyst is a substance that increases the rate of a chemical reaction without being consumed itself. Although it is not consumed the catalyst can be deactivated by secondary processes. Two types of catalyst deactivation are common. First, the loss of a catalyst's reactivity, the ability to increase the rate of a reaction. Second, the loss of selectivity, the ability to direct a reaction to yield a particular product over time . To avoid insufficient catalyst activity and selectivity the deactivation of catalysts should be understood and monitored. This study is to simulate a gas-solid catalyst reaction within a porous structure and model the deactivation process over time.
The type of catalyst deactivation we look at in this study is due to coke accumulation . In this case coke accumulates over catalyst sites, , and reduces the catalyst's reactivity. In this model a gas-phase reactant, , is passed over a porous structure of a catalyst particle, and is converted into a gas-phase product, . This product will then react with the catalyst as the secondary process, and form coke, , on the catalyst sites. This coke accumulation eventually deactivates catalyst sites and reduces the rate of the reaction possible at the sites.
The catalyst particle is roughly assumed to be spherical. As such the reactant's concentration will decrease with radius . Because of symmetry it is then sufficient to use a 2D cutout of the catalyst particle for the simulation.
In this model we assume that there is no external fluid flow involved, which means the reactant is transported by diffusion only; there is no convective mass transfer.
To measure the level of the catalyst deactivation, we can compute the amount of coke that accumulated on the catalyst sites . The net increase in the coke concentration is computed by integrating over the porous structure of the catalyst particle.
The symbols and corresponding units used here are summarized in the Nomenclature section.
Please refer to the information provided in the "Mass Transport" theory manual for a more general theoretical background for mass transport/chemical reaction analysis.
The non-conservative mass transport equation (4) is used to solve for the concentration distribution in a mass transport model:
Here and are the reaction rate constants of this two-steps process. First we will simulate the concentration evolution of reactant , catalyst sites and the resulting product . Once we have computed the concentration profile of , and , we can proceed to solve for the accumulation of coke based on the second step of the reaction in (5).
During the first step the concentration of the reactant denoted as will reduce while the product concentration increase. It is important to note that the concentration of catalyst sites will only get consumed in the second reaction, resulting in a raise of the coke concentration .
The reaction terms could very well also me modeled with a "MassReactionRate" model parameter in place of the "MassSource" term. Then a dependent variable needs to be factored out of the reaction rates and the sign of the equations negated. The choice is purely a matter of taste.
In this model we assume that the catalyst reaction occurs homogeneously from the particle surface into the core of the particle. As such the variation of the species concentration in the azimuthal direction and the polar direction can be ignored. Then a two dimensional model is sufficient to represent the 3D catalyst particle. Here we make use of the lower-left quarter of the catalyst particle as the modeling domain .
To build a realistic geometry that represents the catalyst microstructure, we will use a microscopy image of the catalyst . Note that the image needs to be binarized in order to distinguish the gas phase region (shown in white) and the solid phase region (shown in black).
The binarized image has a finely resolved catalyst boundary. When we discretize this image into a mesh a large number of elements is expected. Since a large number of elements will result in a longer run time and a larger memory usage one may consider to defeatured the domain. Defeaturing is done by making use of a less accurate catalyst boundary while still preserving the important characteristics of the original structure. The details about defeaturing the domain and its construction is presented in the appendix section: Construction of De-featured Modeling Domain. For the sake of accuracy we make use of the original region.
To generate the full element mesh, we have to set up element markers to differentiate the gas phase region and the solid phase region. These markers can then be used to specify, for example, the species diffusivity in different regions of the simulation domain.
More information on markers and their generation in meshes can be found in the Element Mesh Generation: Markers documentation.
In this model the diffusivity for the reactant and the product are different between the gas and solid domains. The diffusivity of the catalyst sites is set to zero to represent its stationary feature.
A simple way to model the discontinuous initial concentration is by using the If statement with element markers:
However, for the default second order mesh a discontinuous initial condition may result in overshoot/undershoot values for on the gas-solid interface due to the quadratic interpolation. This issue is explained in a separate tutorial: Finite Element Method Usage Tips.
To improve the accuracy we can use a linear interpolation instead. To do so we create a first order mesh and extract the original (second order) interpolation result on first order nodes. Then we compute the linear interpolation based on these values.
There are three types of the boundary conditions involved in this mass transport model: the concentration boundary condition, the surrounding flux boundary condition, and the symmetric boundary condition.
A default symmetric boundary condition is implicitly applied on the axis of symmetry.
A concentration boundary condition is used to model the supply of the reactant on the open surface. To avoid a discontinuous jump of reactant at the boundary a helper function that smoothly ramps up the reactant supplied over time is created.
During the catalyst reaction the product is free to diffuse out of the simulation domain through the open surface. This is modeled by a surrounding flux boundary condition.
As mentioned in the previous section, due to the large number of elements used we expect a longer run time and a larger memory usage to find a solution than in typical application examples shown. The number of the degrees of freedom is proportional to the number of elements and gives a measure of the effort needed to solve a PDE.
During the reaction the reactant gradually diffuses into the porous structure of the catalyst particle, and is converted into the product . This product then reacted within the porous structure and forms the coke on the catalyst sites . The coke accumulation eventually deactivated catalyst sites by decreasing the concentration of the catalyst sites.
The coke gradually accumulates throughout the porous structure of the catalyst particle over time. To measure the level of the coke accumulation, the increase in the mean coke concentration: is computed. This is done by integrating over the solid phase domain :
Due the large number of elements used, it took a significant amount of time/memory to solve this PDE model. To reduce the number of elements and the time/memory usage, we can defeatured the domain instead, that is, we model a less accurate catalyst boundary but still preserving the important characteristics of the original structure.
To build a defeatured domain we will smooth out the catalyst boundary in the binarized image with DeleteSmallComponents.
Compared with the result from the original domain: here and here, it is seen that the catalyst reaction took a slightly longer time to propagate to the inner region through the porous structure. This is because some minor channels within particle has been smoothen out and treated as a solid phase region with much smaller diffusivity, which makes it more difficult for the reactant to pass through.
3. Ciesielski, P. N., Robichaud, D., Donohoe, B., Nimlos M., Microscale Simulation of Catalyst Deactivation during Gas-Phase Upgrading of Biomass Pyrolysis Vapors, Biosciences Center, National Renewable Energy Laboratory (2015).