Long-time simulation of catalytic MTO reaction in a fluidized bed reactor with a coarse-grained discrete particle method -EMMS-DPM

Journal Pre-proofs Long-time simulation of catalytic MTO reaction in a fluidized bed reactor with a coarse-grained discrete particle method -EMMS-DPM Xingchi Liu, Ji Xu, Wei Ge, Bona Lu, Wei Wang PII: DOI: Reference:

S1385-8947(20)30126-1 https://doi.org/10.1016/j.cej.2020.124135 CEJ 124135

To appear in:

Chemical Engineering Journal

Received Date: Revised Date: Accepted Date:

2 September 2019 10 January 2020 13 January 2020

Please cite this article as: X. Liu, J. Xu, W. Ge, B. Lu, W. Wang, Long-time simulation of catalytic MTO reaction in a fluidized bed reactor with a coarse-grained discrete particle method -EMMS-DPM, Chemical Engineering Journal (2020), doi: https://doi.org/10.1016/j.cej.2020.124135

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2020 Elsevier B.V. All rights reserved.

Long-time simulation of catalytic MTO reaction in a fluidized bed reactor with a coarse-grained discrete particle method EMMS-DPM Xingchi Liu1,2, Ji Xu1*, Wei Ge1,2,3*, Bona Lu1,4, Wei Wang1,2 1State

Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China 2School

of Chemical Engineering,

University of the Chinese Academy of Sciences, Beijing 100049, China 3Collaborative

Innovation Center of Chemical Science and Engineering (Tianjin), Tianjin 300072, China

4Dalian

National Laboratory for Clean Energy, Dalian 116023, China

*Corresponding authors, E-mail: [email protected] (J. Xu);[email protected] (W. Ge)

ABSTRACT Long-time behaviors widely exist in gas-solid reactors, e.g. the coke deposition and catalyst deactivation in the methanol to olefins (MTO) process. Due to the limitation of the computational cost, previous simulations of MTO reactors were mainly conducted with two-fluid model at short-time scales and assumed that the catalyst are perfectly mixed. Efficient simulation methods are, therefore, highly desirable, even for evaluating the assumptions through very long-time investigation. In this study, the MTO reaction kinetics is incorporated into a coarse-grained discrete particle method, namely EMMS-DPM, to simulate a pilot-scale reactor for 8 hours of physical time. The simulation is validated by comparing the obtained gaseous products and coke contents 1

with the experimental data. The perfect mixing assumption is validated by observing the residence time distribution of the catalyst for the whole reactor. The detailed coke content could provide more accurate data for accelerating steady state simulations. Furthermore, the relation between the coke content and age of the catalyst is directly obtained for the first time, which is valuable to effective use of the catalyst before discharging. The coke content is positively correlated to the age, but is also influenced by local conditions. For example, near the distributor, the catalyst particles are not well mixed and the maximum reaction rate fluctuates significantly even at large time scales. Thus, long-time particle-scale simulation would provide more dynamic details to characterize the reaction in the critical regions of the reactor, thus being very helpful to optimize operating conditions and the reactor designs.

Key words Long-time simulation; Discrete particle method; MTO; Coke distribution

1.

Introduction

Gas-solid fluidization reactors have various industrial applications including fluid catalytic cracking, coal/biomass gasification and pyrolysis, methanol to olefins, catalytic polymerization, ore calcination etc. [1] Many chemical processes need a long period of time to achieve being adequately completed [2-5]. For example, in the methanol-to-olefins (MTO) reactor, hours of physical time is required for the catalyst particles to achieve appropriate coke contents so that higher selectivity to light olefins 2

could be obtained [6-11]. However, the residence time of catalyst particles shows a much wider distribution due to complex gas-solid flows and reactions in large reactors [6]. Therefore, long-time particle-scale simulation of these reactors is essential to improve our understandings of the dynamic behaviors and variation of reactions in these long chemical processes and further helps the optimization of the process and equipment design. Industrial or large-scale gas-solid reactors are usually simulated with the two-fluid model (TFM). However, the traditional TFM is rather computationally expensive since it needs fine grid and small time step to resolve the mesoscale structures in fluidized beds [12]. Therefore, the simulations using the TFM were conducted for a relatively short time-scale, typically seconds of physical time [13, 14], which is not sufficient for investigating a long-time chemical process, e.g. coke deposition and catalyst deactivation [2, 15]. The sub-grid drag models, such as the filtered drag model [16-18] and the EMMS model [19-21], could overcome this limitation. Since these models take mesoscale structures into account to calculate the drag force, the results are reasonable even using coarse grids and large time steps. To further improve the efficiency, some researchers focus on the steady state with the pre-assumed model and/or parameters. Zhao et al. [22] applied a simplified reaction model to simulate the reaction in an MTO demo-scale reactor, but the influence of coke deposition is not taken into account. Lu et al. [23-25] and Li et al. [26] considered the influence of the coke content distribution and obtained better distribution of the gaseous products. However, the coke content

3

distribution is given under the assumption that the catalyst is perfectly mixed [23-26], which might not be the real condition of large reactors. Moreover, due to lack of the detailed coke content distribution, the coke on catalyst particles was only represented by several continuous solid species in TFM [25, 26], which does not represent the real coke content distribution of the catalyst particles, and hence limits accuracy of reactions process. The particle-scale information could be improved with the Lagrangian methods [27-31], e.g. the discrete particle methods (DPM) [32] or the combined computation fluid dynamics and discrete element method (CFD-DEM) [33]. However, huge computational cost due to the huge number of particles is a major issue for this type of methods. Thus, both the spatial and temporal scales in previous simulations are usually much smaller than the requirements of the industrial apparatus [34, 35]. For example, Yuan et al. [34] observed the perfect mixing of catalyst particles in a small fluidized bed under different gas inlets, but the average residence time of particles is only between 5 s and 30 s. Zhuang et al. [35] studied the hydrodynamics and reaction in a two-dimensional reactor containing only 40,000 particles and the physical time is only 25 seconds. To improve the capability of Lagrangian methods, newer models of lower computational cost and much more efficient computing methods are required [36]. Coarse-grain models aim to obtain reasonable results by reducing the number of modeling particles or fluid cells [37-41], so that the computational cost will be greatly

4

reduced. According to the method of calculating contact force between the coarsegrained particles, or parcels, the coarse-grain models could be classified into two kinds. One type of methods are represented by the multiphase particle-in-cell (MP-PIC) method [38], which considers the interaction between parcels through solid pressure. The other type of methods are represented by similar-particle-assembly model (SPA) [39], which treats the collision of parcels as those in the DEM [42]. These methods have already been validated against experiment [43-47] and applied in the simulations of some typical industrial-scale apparatuses [48-51]. Different from previous coarse-grain models, the EMMS-DPM [40, 41] coarse-grains the real particles based on the mesostructure in gas-solid flows which is determined by the energy-minimization multiscale (EMMS) models [52-54], so that larger coarse-graining ratio and longer time-step could be used, which will greatly reduce the computational cost. Furthermore, implementing EMMS-DPM with the heterogeneous super-computing methodologies [55-57], the large-scale high-performance computing for the gas-solid reactors is achieved with high efficiency. Thus, hours of physical time could be achieved for simulating gas-solid flows, providing possibilities to comprehensively study the longtime behaviors in the gas-solid reactor. In this study, the MTO process, in which the long-time behaviors of the catalyst particles are important for their selectivity, is chosen as an example. In Section 2, the models for gas-solid flows and the method of integrating the seven lumped MTO reaction kinetics [58] into EMMS-DPM are introduced. The configurations of the

5

reactor and the settings of the simulation are described. Then, an eight-hour physical time simulation is conducted. In Section 3, the simulation results are validated by comparison with the corresponding experiments. The residence time distribution, the particle age, the coke distribution, and their relationship are analyzed. To the best of our knowledge, it is the first time to directly correlate the particle age to the coke distribution in MTO process. The dynamic reaction in local regions is discussed, providing more information to optimize the reactor. Finally, the conclusions and the prospects are drawn. 2.

Methods and simulation settings

2.1 Models of the solid particle In EMMS-DPM, the coarse-grained particle (CGP) is tracked individually like the real particle, with its motion governed by the Newton’s second law [42]:

mCGP

dvCGP  FG  Fgp  Fc  Fdrag dt

(1)

The diameter of a CGP (dCGP) and the internal voidage of a CGP (εCGP) are two key factors in EMMS-DPM [40, 41]. dCGP should meet the condition given in equation 2, namely, dCGP should be smaller than the minimum diameter of the cluster (dcl,min) to capture the mesoscale structure and be much larger than the diameter of the real particle (dp) to greatly reduce the computational cost.

d p  d CGP  (d cl,min ) EMMS

(2)

εCGP should be chosen according to equation 3,

 CGP  max( mf ,  CGP,min ) , 6

(3)

where εmf is the voidage at the minimum fluidization, εCGP,min is the minimum voidage of a CGP. Therefore, the number of real particles inside a CGP (Np,CGP) is

N p,CGP  (

d CGP 3 ) (1   CGP )  k 3 (1   CGP ) , dp

(4)

where k is the coarse-graining ratio, namely the ratio of dCGP to dp. Thus the mass of a CGP could be presented as,

mCGP  N p,CGP mp  k 3 (1   CGP )mp .

(5)

FG and Fgp are the gravity and pressure gradient force on the CGP, respectively. Fdrag is the drag between the gas phase and CGP, which is the sum of the drag acting on all the real particles in the CGP:

Fdrag 

  3 , (u g  v CGP )(1   CGP ) d CGP 1 g 6

(6)

where ug and vCGP are the velocity of the fluid phase and CGP, respectively. β is the coefficient of the drag force. Different from the conventional drag laws [59, 60], the heterogeneous structure of solid particles at the sub-grid scale is considered in the EMMS drag model [21], which is proved to be more accurate to describe the hydrodynamics of the gas–solid flows [61-64]. Thus, the EMMS drag model is employed in this study, in which β is,



3  g (1   g )  g | u g  us | CD0 g 2.7 H D , dp 4

(7)

where HD is heterogeneity index and CD0 is the standard drag coefficient of a single particle [65]:

7

CD0

 24 0.687  Re (1  0.15 Rep )  p 0.44  Re p 

Rep  1000

(8)

Rep  1000

 g  g | u g  us | d p g

(9)

It should be noted that Rep is the Reynolds number of the real particle in the CGP. HD is adopted from the previous work [23, 24]: 1.6864  0.1443  1  ( / 0.4064) 20.4805 g   H D  (0.0999  0.9586 g )1.8008   0.08597 (1649.4591  1648.5453 g ) 

0.4   g  0.4924 0.4924   g  0.9505

(10)

0.9505   g  1

Fc is the contact force between the CGPs. In EMMS-DPM, the real particles in a CGP might be in a loose cluster, and thus, in order to model the close-packing state, the contacting diameter (dhc) is the diameter of the close packing “hard core” to calculate the contact force. The dhc and the restitution (eCGP) of the CGPs is rescaled to the structure and coarse-graining ratio of the CGPs [41]:

dhc  1 CGP  dCGP

(11)

eCGP  1  (ep2  1)k (1   CGP )1/3

(12)

1/3

Please refer to our previous work [41] for the details of the model. Additionally, the temperature is this pilot-scale reactor is assumed homogeneous according to the experiments [23]. Thus, the heat transfer between the coarse-grained particles is not taken into account. Similar to the work of Lu et al. [23], the temperature is set as 738K. 8

2.2 Models of the gas phase Gas phase is governed by the volume-averaged Navier-Stokes equations which consist of the conservation equations of mass and momentum:

( g g ) t  ( g  g u g ) t

   ( g g ug )  0

(13)

   ( g  g u g u g )    g  p     g   g  g g  (u s  u g )  (14)

where εg and ug are the volume fraction and velocity of the gas phase, respectively; τg is the viscous stress tensor and β is the drag coefficient. No energy equation is involved in this study since the temperature is assumed homogeneous in the pilot-scale reactor and the heat-transfer is not taken into account.

2.3 Reaction model The seven lumped reaction kinetics model for DMTO [58] is employed, in which seven products are directly produced from the methanol (MeOH). The formation rate of each products Ri (i=CH4, C2H4, C3H6, C3H8, C4H8, C5H10 and coke) is,

Ri  ni k i i C MeOH M i ,

(15)

where ni is the stoichiometric number, ki is the reaction rate constant of lump i (L/gcat/s), CMeOH is the methanol concentration (mol/L), Mi is the molar weight (g/mol) and φ is the deactivation function:

i 

A exp( i wc ) , 1  B exp( D  ( wc  E )) 9

(16)

where wc is the coke content (g/100gcat). A, B, D and E are the constants, which are 1, 9, 2 and 7.8, respectively. The αi is the deactivation constant. The detailed value of ki and αi for each species could be found in the literature of Lu et al. [23]. The Ri is applied in species transport equations to calculate the content of each gas species:

( g gYi ) t

  ( g g ugYi )    g Ji   g Ri ,

(17)

where Yi is the mass fraction of species i and Ji is the mass diffusion coefficient. The solid volume fraction is quite high in the current bubbling bed, at this operating condition, the turbulence has small effect on the catalytic reaction [66]. Thus, the species transport equations do not consider the turbulent effect. Coke is the only solid product, which will be deposited on the catalyst particles. Thus, it is calculated in the particle part with the formula written as:

dm coke  Rcoke m CGP dt

(18)

where Rcoke is the coke formation rate, which is determined by the equation (15). In order to obtain more detailed reacting process, the reaction is realized at the particle scale, namely the coke content on each CGP, rather than the average coke content in the fluid cells, is applied to calculate φcoke in the equation (16). In current MTO reaction model [58], the reacting rate is related to the mass distribution of catalyst particles, and thus determined by the mass distribution of CGPs, because the mass is conserved between the real particles and the CGPs in the EMMS-DPM.

10

2.4 Simulation settings As depicted in Fig. 1, the pilot-scale MTO reactor consists of three parts, namely, the lifting tube, the bubbling reactor and the sedimentation bed. The lifting tube is used for discharging the deactivated catalyst, and the sedimentation bed is designed for gassolid separation, where the reaction rarely happens in these two parts. The bubbling reactor is the main region for the reaction, where the catalyst particles have complex movements and cause complex reaction. Thus, to reduce the computational cost, only the bubbling reactor and a small part of the lifting tube are selected for the simulation.

Fig. 1. Configurations of the pilot-scale MTO reactor [23] (a); mesh and operating conditions for the simulated part (b).

11

Consistent with the experiments, 9 kilograms of fresh catalyst particles (no coke on them) are placed at the bottom of the bubbling reactor at the initial state. During the reacting process, the fresh catalyst particles are fed into the reactor at the rate of 8.8 kg/h from the top of the bubbling reactor and the coke deposited catalyst particles are discharged at the bottom of lifting tube at the same rate, so that the weight of the catalyst in the reactor is kept constant. The main reactants (Ug1 in the Table 1) consist of the methanol and steam (mass ratio 2:1), and are injected into the reactor from the distributor at the rate of 27 kg/h. The lifting stream (Ug2 in the Table 1), which is aimed to loosen the discharging particles, enters from the lifting tube at the rate of 0.12 kg/h. The effect of the coarse-graining parameters of the EMMS-DPM on the cold-flow has already been validated in our previous study [40, 41]. In order to validate the applicability of the EMMS-DPM to the reacting process, the influence of the coarsegraining parameters on the products of the reaction is studied first. According to the model proposed by Hu et al. [67], for the operating conditions of the pilot-scale MTO reactor, the (dcl,min)EMMS is 350dp., the εmf and the εCGP,min are both 0.4 under current operating conditions, thus the max(εmf, εCGP,min) is 0.4. The εCGP,min is constant for the whole simulation, because the it is determined by the operating condition which is not changed in the whole simulation. Thus, according to equations (2) and (3), three combinations of the dCGP and εCGP is designed, namely dCGP=25dp with εCGP=0.5, dCGP=30dp with εCGP=0.42, and dCGP=30dp with εCGP=0.5. The operating conditions, the properties of the particle and the simulating parameters are set according to the

12

experiments, which are given in Table 1, the thermal properties of the gas species and the catalyst could be referred to the work of Lu et al. [23]. Each case is simulated for 60 seconds of physical time for analysis.

Fig. 2. Evolution of the ratio of ethylene to propylene (C2H4/C3H6) at the outlet (a), and the average coke deposition on the particles (b) under different coarse-graining parameters. The evolution of the main gaseous products and the coke deposition on the catalyst particles are given in Fig. 2. It shows that the results related to the reaction are almost the same for the three cases, indicating that if the coarse-graining conditions for the cold gas-solid flow are met, the coarse-graining parameters have minor effect on the reaction. The reason is, as described in Section 2.3, the reacting rate is related to the mass of the catalyst, namely the distribution of the catalyst. Thus, in the reaction kinetics used in this work, proper prediction of the cold-flow will give proper reaction behaviors. Furthermore, to reduce the computational cost in the following work, the dCGP and εCGP are chosen as 30dp and 0.42, respectively. 13

Table 1. Simulation parameters of the MTO reactor Item

Variables

Value

Diameter of the real particle, dp (μm)

97

Catalyst inventory in reactor (kg)

9

Input catalyst rate (kg/h)

8.8

Density, ρp (kg/m3)

1500

Young's modulus, Y (Pa)

5×106

Thermal conductivity (W/(m K))

0.0454

Heat capacity (J/(kg K))

1220

Particle

Mass fraction of the coke on the input catalyst 0 (g/100gcat) Time step (s)

2.5×10-5

Reactor gas velocity, Ug1 (kg/h)

27

Lifting gas velocity, Ug2 (kg/h)

0.12

Gauge pressure at the exit (MPa)

0.024

Operating temperature (K)

738

Density, ρg (kg/m3)

0.4288

Dynamic viscosity, μ (Pa s)

2.43×10-5

Time step (s)

5×10-4

Gas

14

The simulation is conducted on the upgraded Mole-8.5 supercomputer [55], which features a heterogeneous architecture that couples the computing of central processing units (CPUs) and graphics processing units (GPUs). The software of EMMS-DPM is developed according to this architecture and contains three parts: the gas phase solver, which is developed based on the open source software OpenFOAM [68] and executed in CPUs; the particle simulator, which is based on our in-house developed codes DEMms [55] and executed in GPUs; the gas phase and the solid particles are coupled by the interaction model. In this study, the PISO method for the pressure velocity coupling is employed to solve the gas phase. The preconditioned conjugate gradient (PCG) solver and a diagonal incomplete-Cholesky preconditioner are employed as the solver. Moreover, for higher accuracy and numerical stability, a bi-conjugate gradient stabilized (PBiCGStab) [69, 70] is employed to solve the discretized velocity equations. The gas phase is solved on CPUs (Intel Xeon E5-2680) and the solid particles are computed on GPUs (NVIDIA Tesla K80). The overall simulation runs about 420 physical seconds per day using only 6 K80s in one node, and totally 8 hours of physical time is achieved for simulating the MTO reaction process.

3.

Results and discussion

In this section, the methods are validated by comparing the results at the steadystate between the simulation and experiment. The particle behaviors and reaction in

15

both the global reactor and key local regions are analyzed. Besides, due to the advantages of the long-time particle-scale simulation on explicitly tracking the particles, some assumptions in previous studies [25, 26] are analyzed. 3.1 Products at the steady state The ethylene and propylene are desired products in MTO process. The ratio of ethylene to propylene is important due to the market demand [71]. As shown in Fig. 3a, the ratio of ethylene to propylene, namely C2H4/C3H6, increases rapidly at the initial stage, typically before 3 hour, and then gradually increases to a plateau after 5 hour, reaching to the steady value about 1.196, which is consistent with the experiment (C2H4/C3H6=1.245) [23]. As shown in Fig. 3b, similar evolution is also observed for the coke content on the catalyst particles. This is because according to the study of Bos et al. [2], the ratio of ethylene to propylene primarily relies on the coke, namely the higher coke will result in more ethylene, and vice versa. Between 5 and 8 hour, the coke deposition reaches to the steady state. The average coke will decrease when the fresh catalyst particles are fed and the deactivated particles are discharged. At the steady state, the newly generated coke on the particles and the removed coke due to discharging of particles reach the dynamic equilibrium state, thus the average coke content fluctuates slightly along time. The average coke content on all the catalyst particles is 5.533 g/100gcat between 7 and 8 hour, which is close to the experimental data (6.0 g/100gcat). The evolutions of the C2H4/C3H6 and the coke content of the particles indicate that the MTO reaction has reached to a steady state after 5 hour, and the data between 7 and 8

16

hour will be used for the analysis of the reaction at the steady state.

Fig. 3. Evolution of the ratio of ethylene to propylene (C2H4/C3H6) at the outlet (a) and the average coke (b). The mass fractions of all the gaseous products at the outlet are shown in Fig. 4, which are in good agreement with the experimental data [23]. Detailed results are given in Table 2, which shows good consistence with the experiment [23], especially the key factors, such as the methanol conversion (ηMeOH) , the selectivity of the light olefins (S(C2H4+C3H6)) and the ratio of ethylene to propylene (C2H4/C3H6). The discrepancies of the pentene mass fraction are relatively large, might due to the lower accuracy of the kinetic model in predicting this species [11, 72].

17

Fig. 4. Normalized mass fraction of the gaseous products at the steady state (except the H2O and methanol). Table 2. The properties of the MTO reaction process obtained in the simulation and experiment [23] Simulation

Experiment.

Error (%)

YCH4

0.0144

0.01643

12.08

YC2H4

0.4378

0.4303

1.74

YC3H6

0.3662

0.3457

5.92

YC3H8

0.0293

0.03788

22.73

YC4

0.1080

0.0913

18.29

YC5

0.0443

0.07836

43.46

wc (g/100gcat)

5.533

~6

~7.78

ηMeOH (%)

99.99

99.9444

0.05

S(C2H4+ C3H6)

80.40

77.603

3.6

18

C2H4/C3H6

1.196

1.2449

3.96

The time-averaged mass fraction distributions of the ethylene, propylene and methanol at the steady state are given in Fig. 5. The axial mass fractions of the ethylene and propylene increase along the height below 0.3 m, while the axial mass fraction of the methanol decreases in such region is almost zero above 0.3 m, indicating that the reaction mainly occurs near the distributor. Besides, more than 90% particles locate in the region below 0.3 m at the steady state, and thus the region below 0.3 m is referred to as the main reaction zone. Similar spatial distribution of light olefins was also reported by Lu et al. [23] and Zhao et al. [22]. In the radial direction, especially near the distributor (height=0.06 m), the mass fractions in the centers are slightly lower than those near the walls. This is because the reactor has two gas inlets, as described in Section 2.4. One is the mixture of the steam and methanol, which is injected into the reactor from the distributor. The other is the lifting gas of the pure steam injected through the lifting tube. Thus, the concentration of methanol near the wall is slightly higher than that in the center, leading to more intense reaction and hence more products near the wall. However, no differences of the mass fractions are observed in the top region, e.g. height=0.30 m, indicating that the products will rapidly mix in the radial direction of the reactor. It indicates that the products at the outlet is homogenous in this bubbling fluidized reactor, and similar phenomenon was also reported by Lu et al. [23].

19

Fig. 5. Time-averaged mass fraction at the steady state: (a) the ethylene, propylene and methanol in the axial direction; (b) the ethylene and propylene in the radial direction at different heights. In general, the main properties of the MTO process are well predicted by coupling the reaction kinetics with the EMMS-DPM, which validates this particle-scale simulation method. Thus, with long-time simulation, this method would provide deeper understandings of the MTO reactor at the particle-scale, e.g. the residence time distribution of particles, the coke contents on particles and the reaction behaviors in local regions.

3.2 Residence time and coke content of the catalyst particles in the whole reactor Residence time of catalyst particles, which reflects the overall movement of particles in a reactor, is important for the reaction [73]. Actually, the coke depends on the particle age, which is defined by the duration of the particles from the last entrance into the reactor or some region of the rector [74]. However, it is difficult to directly 20

obtain the particle age with experiments. Thus, in the previous study [25], the coke distribution of catalyst particles at the steady state is correlated to the residence time distribution (RTD) of particles. The RTD is derived from the continuous stirred tank reactor (CSTR) model, which assumes that the catalyst particles are perfectly mixed in the reactor. However, this assumption needs further validation, because the mixing might be influenced by inner structures of the bubbling fluidized bed [75]. In the theory of bubbling fluidization, the particles of identical diameter will feature perfect mixing [76] and the characteristic time of perfect mixing is:

 perfect 

Is Gs

(19)

where Is is the total catalyst inventory in the reactor and Gs is the discharging rate of the catalyst particles. According to the inventory (9 kg) and discharging rate of the particles (8.8 kg/h) in current reactor, the characteristic time of perfect mixing (τperfect) is 3681 s.

Fig. 6. Residence time distribution of the catalyst particles. The RTD of the catalyst particles obtained from the simulation of the eight-hour 21

MTO process is given in Fig. 6. If the particles are perfectly mixed, the RTD of the catalyst particles could be expressed by the probability distribution function [77]:

p (t ) 

1

r

exp(

t

r

)

(20)

where t is the residence time, τr is the average residence time if the perfect mixing is achieved. The results in Fig. 6a could be well fitted to Eq. (20), and the characteristic time (τr) is 3579 s. This value is consistent with the characteristic time of perfect mixing (τperfect) and thus indicates that the catalyst particles are perfectly mixed. Therefore, the assumption of CSTR model in the work of Lu et al. [23, 24] is validated in this work.

Fig. 7. Axial distribution of the particle age (a) and the coke content (b) between 7 and 8 hour. To gain deeper insight of the mixing of catalyst particles in the main reaction zone of the reactor, the age and coke of the catalyst particles along the axial direction is given in Fig. 7. Apparently, their axial distributions have similar trends, indicating that they have strong relationship with each other. Both the age and the coke are identical in the 22

main reaction zone, namely below 0.3 m, and are nearly the same as the average one of all the particles in the reactor (indicated by the red dash-line in Fig. 7ab), since this region contains 90% particles. The average particle age (taverage) is 3619 s, which is consistent with the average residence time of perfect mixing (τperfect=3681 s). The average coke of the particles in the main reaction zone is 5.533 g/100gcat, which is consistent with the experiment [23]. Thus, in the main reaction zone of the reactor, the perfect mixing assumption is further validated according to the average coke and particle age. In the upper region, the age and coke decrease with the height due to the dominance of newly fed particles of shorter age and lighter coke, indicating that the newly fed particles are not well mixed with the original particles immediately. However, this phenomenon has little influence on the final production, because, as shown in Fig. 5, the reaction seldom happens in the upper region.

Fig. 8. The relationship between the coke and the particle age at the steady state.

23

Although the CSTR assumption is validated, it is more important to investigate the relationship between the coke contents and particle age. Because it will offer detailed information on obtaining appropriate coke content distribution of the catalyst particles by adjusting their residence time distribution [6]. As shown in Fig. 8, the coke is generally positively correlated to the catalyst particle age in the reactor. This is because, in current reaction kinetics model [58], the coke will not be consumed in the reaction. If the particles locate in the reactor, the coke will keep generating at a certain rate. When the particle age is shorter than 1 hour, smaller coke on the particle leads to faster generation of coke accordingly, which is consistent with the theory [58]. The coke on the particles of 1 hour is about 6.7 g/100gcat, which is slightly lower than the best operating condition of the coke (7.8 g/100gcat) [58]. However, according to the residence time distribution, 62.4% discharged particles have the age less than 1.0 hour, indicating that the catalyst is discharged too early. Thus, in current operating conditions, catalyst is not effectively used, which brings in extra burden for the regeneration of the catalyst. It should be considered in further research and development of the reactors to balance the production and energy consumption. When the particle age is longer than 4 hours, the coke increases quite slowly and then reaches to a platform for two reasons. First, according to RTD of particles (Fig. 6), only 2% particles have the age longer than 4 hours, resulting in small probability of newly generated coke. Second, even if these particles participate in the reaction, the reaction rate is very slow because higher coke will lead to much slower reaction rate. However, when these particles participate in

24

reaction, they will lower the conversion of the methanol, leading to less desired products. Furthermore, as shown in Fig. 8, the coke on the particles of the same age disperses in a certain range. This detailed coke content distribution is hard to be considered in TFM simulations [23, 24]. From this aspect of view, the EMMS-DPM simulation coupled with MTO reaction kinetics should be more accurate in predicting the MTO reactions because the coke content distribution could be considered at the particle-scale. From another aspect of view, the coke content distribution of catalyst particles obtained with the EMMS-DPM at the steady state could be provided for TFM simulations, so that the accuracy of predicting the MTO reaction would be improved and the time to solution would be greatly reduced. 3.3 Residence time and methanol concentration in local regions The analysis of the whole reactor offers the overall insights of mixing and reaction in the MTO reactor. Furthermore, to fully evaluate the movement of catalyst particles and the detailed reaction behavior, analysis for the key regions of the reactor is conducted. The bottom part of the reactor (height<0.42 m) is divided into sub-regions in the axial direction with 0.06 m interval, and the corresponding RTDs of the particles are analyzed to quantify the mixing behaviors. As shown in Fig. 9, most RTDs (0.06 m0.42 m) are unimodal. Their peaks are about 0.1s and the range is narrow, typically between 0.05 s and 0.30 s. However, the RTD of the lowest region (0 m-0.06 m) is bimodal and has broader range. The first and second peaks of the RTD is 0.1 s and 0.40

25

s, respectively, and about 24.9% particles have the residence time longer than 0.40 s, indicating poorer mixing of the particles in this region.

Fig. 9. Residence time distribution of particles in the sub-regions between 7 and 8 hour. To figure out the reason for the poor-mixing in the lowest sub-region, it needs to locate the region where the particles stay longer. The instantaneous age in sub-regions is the duration of a particle in the corresponding sub-region. It could be obtained from the difference between the time that the particle entering the sub-region and the sampling time. The distribution of particle age in sub-regions aims to figure out where the particles stay longer. Thus the radial distribution of the average particle age in the lowest sub-region between 7 and 8 hour is analyzed. As shown in Fig. 10, the age in the lowest region first decreases from the center of the region, reaching the minimum age of 0.31 s, and then increases to 0.396 s when the particles locate near the wall. It indicates that the particles stay longer in the center and periphery of the reactor, but 26

shorter in the middle. Because in the center region, the bubbling bed connects the lifting tube, which has much slower gas velocity compared with the distributor, leading to poorer mixing. The poorer mixing near the wall could be explicitly explained by the movement of particles shown in Fig. 11. The age distribution of particles at the planes of x = 0 m, y = 0 m, and z = 0.01 m are given in Fig. 11Error! Reference source not found.a. It shows that the particles of longer age (highlighted by the red circle) typically locate at the corner due to smaller velocity of particles (Fig. 11Error! Reference source not found.bc), which is consistent with the findings in Fig. 10. Near the distributor, the inlet gas would push the particles moving inclined in the radial direction. Combined with the particles falling down form the upper of reactor near the walls, the velocity of particles near the wall is relatively small, leading to different mixing of particles in the center and periphery of the reactor.

Fig. 10. Radial distribution of the particle age in the sub-region below the height of 0.06 m at the steady state. 27

Fig. 11. Instantaneous age distribution of the particles (at the planes of x = 0 m, y = 0 m, and z = 0.01 m) (a); the velocity distribution of particles near the distributor at the plane of x = 0 m (b) and y = 0 m (c). The poorer mixing near the distributor may lower the efficiency of the catalyst. If the particles of high coke contents stay in these regions too long, the reaction rate will decrease accordingly. Since the poorer mixing of particles near the distributor is caused by the slower gas velocity in the lifting tube and wall effects, increasing gas velocity or changing the location of the lifting tube might be possible methods to optimize the mixing of particles [78].

28

Fig. 12. Maximum and average mass fraction of the MeOH in different time intervals at the height of 0.01m (a) and 0.11m (b). Besides the movement of catalyst particles, the distribution of the concentration of the MeOH (YMeOH) is critical to the MTO reactions. The maximum and the average YMeOH at the steady state at the height of 0.01 m is 0.64 and 0.17, respectively, which indicates that the reaction would be non-uniform near the distributor. It is also found that the concentration of methanol has large fluctuations at different time scales, e.g. 30 s, 300 s and 1200 s. YMeOH near the distributor, e.g. the height of 0.01 m and 0.11 m, is shown in Fig. 12. At both sampling heights, the maximum and average YMeOH obtained from short time intervals have stronger fluctuation than those obtained from long time intervals, indicating that the reaction is dynamic near the distributor. Thus, current structure of the distributor might need further optimization for more homogeneous of the methanol distribution, aiming to obtain more stable operation of the reactor. Additionally, this phenomenon could not be captured in short-time simulations [23] and thus long-time simulation might be necessary to optimize the reactor structure. 29

Fig. 13. (a) Maximum and average reaction rate of the MeOH at different time intervals at the height of 0.01m and (b) the distribution of the MeOH reaction rate at the steady state. To further evaluate the reaction near the distributor, the reaction rate of methanol (RMeOH) at the height of 0.01 m is analyzed. As shown in Fig. 13a, the maximum and the average RMeOH at 0.01 m both increase at the initial stage, which is the reason of the increase of the YMeOH shown in Fig. 12a. Similar to YMeOH, the RMeOH has large fluctuations at different time scales, indicating non-uniform and dynamic MTO reaction. To quantify the probability of different reaction rate at the steady state, the distribution of the RMeOH at 0.01m is given in Fig. 13b. The reaction rate around the time-averaged one occupies 24.5% of the total time, while the largest reaction rate could triple its timeaveraged one. Although it occupies only 0.21% of the total time, it might bring about more heat and rise the local temperature, leading to faster generation of the coke and faster deactivation of the catalyst correspondingly [71]. Besides, the reaction rate less

30

than 10% of the average rate occupies 19.6% of the total time, indicating that the methanol could not be converted efficiently in this region. Furthermore, the cross section at the height of 0.01m (the region connecting the lifting tube and the distributor of the main reactor) is divided into the inner part and outer part. Apparently, the reaction mainly happens in the outer region, where the catalyst immediately contacts with the methanol from the distributor. The reaction rate in the inner region is much smaller, typically less than 10% of the average one in this region, because the methanol concentration is much lower due to the inlet of lifting steam in the inner part. Lower methanol concentration and smaller reaction rate are also the reasons for the less mass fraction of the required products in the inner region (Fig. 5b). The non-uniform reaction might cause difficulties to control the temperature or final products in the reactor. Changing the connection between the bubbling bed and the lifting tube may optimize the distribution of the methanol [79, 80]. Near the distributor of this pilot-scale MTO reactor, the particles are not well mixed and the reaction rate is non-uniform both spatially and temporally. Thus, more consideration should be taken into account for further optimization of this region in the scale-up design of the reactors. Because the mixing might be even poorer due to the complex geometry of the distributor, the complex operating condition and hence bigger difference in the behavior of catalyst particles [6, 22, 81]. For example, near the distributor of a commercial MTO reactor [24], the catalyst particles have even wider velocity distribution, whose maximum and minimum velocities are larger than 5 m/s

31

and almost 0 m/s, respectively. Thus, to optimize the geometry of the distributor or the operating conditions, the long-time particle-scale simulation with EMMS-DPM would be more helpful due to its accuracy and efficiency in modeling the catalyst particles and the reaction.

4.

Conclusions and prospects

In this work, the long-time reaction process is conducted at the particle-scale for a pilot-scale reactor by coupling EMMS-DPM with a seven lumped kinetic reaction model. The main products are consistent with the experiments, e.g. the ratio of the ethylene to propylene, the conversion of the methanol and the average coke content. For the whole reactor, perfect mixing of catalyst particles is observed and the CSTR assumption [25, 26] is validated. The relationship between the coke content and the particle age is directly obtained for the first time. For the current MTO reactor, more than 60% of the catalyst particles are discharged before they reach the required coke content for best selectivity to the light olefins. Besides, about 2% of the catalyst particles have the age longer than 4 hours and the coke reaches 8 g/100gcat, leading to slower reaction rates. The coke on the particles of the same age disperses in a certain range due to the non-uniform reactant distribution, which could not be found under the assumption of uniform reactant distribution in previous studies [25, 26]. Moreover, different behaviors of catalyst particles are observed near the distributor due to the influence of the distributor structure and the two gas inlets

32

operating mode, resulting in poorer mixing and non-uniform reaction rate. For example, in 0.21% of the operation time, the reaction rate of the methanol could triple the timeaveraged one, which will bring in more heat in the local region and might lead faster deactivation of the catalyst. This phenomenon could not be observed in the short-time simulations due to lack of samplings for such small possibility. Thus, the necessity of long-time particle-scale simulation is demonstrated. Furthermore the simulation methods developed in this work, namely implementing the EMMS-DPM with heterogeneous supercomputer and integrating reaction kinetics into gas-solid flows, and even energy equation to consider the temperature evolution in the reaction for larger reactors in the future, would be a powerful tool for studying the gas-solid reactor concerned with long-time particle-scale behaviors, such as catalyst deactivation in Fischer–Tropsch synthesis [82] and particle attrition in fluid catalytic cracking systems [83]. More details of the particle movement and the reaction for the whole reactor or the key local regions of the reactors will be obtained, which will be helpful for the equipment design or optimizing the operating conditions.

Acknowledgments This work is financially supported by the Transformational Technologies for Clean Energy and Demonstration, Strategic Priority Research Program of the Chinese Academy of Sciences (CAS) (No. XDA21030700), the Key Research Program of

33

Frontier Science, Chinese Academy of Sciences (CAS) (No. QYZDJ-SSW-JSC029), the National Key Research and Development Program of China (2017YFB0202203), the National Natural Science Foundation of China (No. 91834303), the Chinese Academy of Sciences (CAS) (No. XXH13506-301), and the Science Challenge Project (No. TZ2016001).

Nomenclature Subscripts g

gas phase

s

solid phase

p

particle

cl

cluster

CGP

coarse-grained particle

EMMS energy-minimization multi-scale MeOH methanol Variables dp

diameter of the real particle (m)

dCGP

diameter of CGP (m)

dhc

diameter of CGP hard core (m)

dcl,min

minimum cluster diameter predicted by EMMS model (m)

CD0

standard drag coefficient for a CGP

CMeOH

methanol concentration (mol/L)

e

restitution coefficient of particle

eCGP

restitution coefficient of the CGPs 34

FG

gravity force (N)

Fgp

pressure gradient force (N)

Fc

contact force (N)

Fdrag

drag force (N)

Gs

catalyst discharging rate (kg/h)

HD

heterogeneity index

Is

catalyst inventory (kg)

Ji

mass diffusion coefficient (m2/s)

k

coarse-graining ratio

ki

reaction rate constant of different products (L/gcat/s)

mCGP

mass of CGP (kg)

mp

mass of the real particle (kg)

mcoke

mass of the coke (kg)

Mi

molar weight of different products (g/mol)

ni

stoichiometric number of different products

Np

total number of real particles

Np,CGP

number of real particles inside a CGP

Ri

reaction rate of species i

Rep

Reynolds number of the real particle

S

selectivity of the light olefins

ug

velocity of gas (m/s)

us

velocity of solid (m/s)

vCGP

velocity of a CGP (m/s)

wc

coke content (g/100gcat)

Yi

mass fraction of species i 35

Greek letters

i

deactivation constant of species i



interphase friction coefficient (kg/m3 s)

CGP

voidage inside a CGP

CGP,min minimum voidage of a CGP mf

voidage at minimal fluidization

g

voidage

g

gas viscosity (Pa s)

ρg

gas density (kg/m3)

τg

viscous stress tensor

τperfect

characteristic time of perfect mixing (s)

φ

deactivation function

ηMeOH

methanol conversion

Reference

[1] W. Ge, L. Guo, X. Liu, F. Meng, J. Xu, W.L. Huang, J. Li, Mesoscience-based virtual process engineering, Computers & Chemical Engineering 126 (2019) 68-82. [2] A.N.R. Bos, P.J.J. Tromp, H.N. Akse, Conversion of methanol to lower olefins. kinetic modeling, reactor simulation, and selection, Industrial & Engineering Chemistry Research 34 (1995) 134-135. [3] M. Khorashadizadeh, H. Atashi, Modeling the kinetics of cobalt Fischer-Tropsch catalyst deactivation trends through an innovative modified Weibull distribution, 36

Physical Chemistry Chemical Physics Pccp 19 (2017). [4] J. Welt, W. Lee, F.J. Krambeck, Catalyst attrition and deactivation in fluid catalytic cracking system, Chemical Engineering Science 32 (1977) 1211-1218. [5] J. Reppenhagen, J. Werther, Catalyst attrition in cyclones, Powder Technology 113 (2000) 55-69. [6] P. Tian, Y. Wei, M. Ye, Z. Liu, Methanol to olefins (MTO): from fundamentals to commercialization, ACS Catalysis 5 (2015) 1922-1938. [7] C.D. Chang, A.J. Silvestri, The conversion of methanol and other O-compounds to hydrocarbons over zeolite catalysts, Journal of Catalysis 47 (1977) 249-259. [8] J. Tan, Z. Liu, X. Bao, X. Liu, X. Han, H.E. Changqing, R. Zhai, Crystallization and Si incorporation mechanisms of SAPO-34, Microporous & Mesoporous Materials 53 (2002) 97-108. [9] J. Li, Y. Qi, Z. Liu, G. Liu, D. Zhang, Co-reaction of ethene and methylation agents over SAPO-34 and ZSM-22, Catalysis Letters 121 (2008) 303-310. [10] J. Li, Y. Wei, G. Liu, Y. Qi, P. Tian, B. Li, Y. He, Z. Liu, Comparative study of MTO conversion over SAPO-34, H-ZSM-5 and H-ZSM-22: Correlating catalytic performance and reaction mechanism to zeolite topology, Catalysis today 171 (2011) 221-228. [11] M. Gao, H. Li, M. Yang, J. Zhou, X. Yuan, P. Tian, M. Ye, Z. Liu, A modeling study on reaction and diffusion in MTO process over SAPO-34 zeolites, Chemical Engineering Journal (2018).

37

[12] S. Sundaresan, A. Ozel, J. Kolehmainen, Toward Constitutive Models for Momentum, Species, and Energy Transport in Gas–Particle Flows, Annual Review of Chemical and Biomolecular Engineering 9 (2018) 61-81. [13] D.C. Sau, K.C. Biswal, Computational fluid dynamics and experimental study of the hydrodynamics of a gas–solid tapered fluidized bed, Applied Mathematical Modelling 35 (2011) 2265-2278. [14] F. Taghipour, N. Ellis, C. Wong, Experimental and computational study of gas– solid fluidized bed hydrodynamics, Chemical Engineering Science 60 (2005) 68576867. [15] D. Wallenstein, T. Roberie, T. Bruhin, Review on the deactivation of FCC catalysts by cyclic propylene steaming, Catalysis Today 127 (2007) 54-69. [16] X. Gao, T. Li, A. Sarkar, L. Lu, W.A. Rogers, Development and Validation of an Enhanced Filtered Drag Model for Simulating Gas-Solid Fluidization of Geldart A Particles in All Flow Regimes, Chemical Engineering Science S0009250918301726. [17] A. Sarkar, X. Sun, S. Sundaresan, Verification of sub-grid filtered drag models for gas-particle fluidized beds with immersed cylinder arrays, Chemical Engineering Science 114 (2014) 144-154. [18] Y. Igci, A.T. Andrews Iv, S. Sundaresan, S. Pannala, T. O'Brien, Filtered two-fluid models for fluidized gas-particle suspensions, AIChE Journal 54 (2008) 1431-1448. [19] J. Wang, Q. Zhou, K. Hong, W. Wang, J. Li, An EMMS-based multi-fluid model (EFM) for heterogeneous gas–solid riser flows: Part II. An alternative formulation from

38

dominant mechanisms, Chemical Engineering Science 75 (2012) 349-358. [20] W. Wang, J. Li, Simulation of gas-solid two-phase flow by a multi-scale CFD approach - Extension of the EMMS model to the sub-grid level, Chemical Engineering Science 62 (2007) 208-231. [21] N. Yang, W. Wang, W. Ge, L. Wang, J. Li, Simulation of heterogeneous structure in a circulating fluidized-bed riser by combining the two-fluid model with the EMMS approach, Industrial & Engineering Chemistry Research 43 (2004) 5548-5561. [22] Y. Zhao, H. Li, M. Ye, Z. Liu, 3D numerical simulation of a large scale MTO fluidized bed reactor, Industrial & Engineering Chemistry Research 52 (2013) 1135411364. [23] B. Lu, H. Luo, H. Li, W. Wang, M. Ye, Z. Liu, J. Li, Speeding up CFD simulation of fluidized bed reactor for MTO by coupling CRE model, Chemical Engineering Science 143 (2016) 341-350. [24] B. Lu, J. Zhang, H. Luo, W. Wang, H. Li, M. Ye, Z. Liu, J. Li, Numerical simulation of scale-up effects of methanol-to-olefins fluidized bed reactors, Chemical Engineering Science 171 (2017) 244-255. [25] J. Zhang, B. Lu, F. Chen, H. Li, M. Ye, W. Wei, Simulation of a large methanolto-olefins fluidized bed reactor with consideration of coke distribution, Chemical Engineering Science 189 (2018) 212-220. [26] H. Li, X. Yuan, M. Gao, M. Ye, Z. Liu, Study of Catalyst Coke Distribution Based on Population Balance Theory: Application to Methanol to Olefins Process, AIChE

39

Journal 65 (2019) 1149-1161. [27] W. Ge, L. Wang, J. Xu, F. Chen, G. Zhou, L. Lu, Q. Chang, J. Li, Discrete simulation of granular and particle-fluid flows: from fundamental study to engineering application, Reviews in Chemical Engineering 33 (2017) 551-623. [28] R. Bing, W. Zhong, C. Yu, C. Xi, B. Jin, Z. Yuan, L. Yong, CFD-DEM simulation of spouting of corn-shaped particles, Particuology 10 (2012) 562-572. [29] C.L. Wu, K. Nandakumar, A.S. Berrouk, H. Kruggel-Emden, Enforcing mass conservation in DPM-CFD models of dense particulate flows, Chemical Engineering Journal 174 (2011) 475-481. [30] A.V. Patil, Y.M. Lau, J.A.M. Kuipers, Modeling 3-D bubble heat transfer in gassolid fluidized beds using CFD-DEM, Industrial & Engineering Chemistry Research 54 (2015) 121-131. [31] F. Li, F. Song, S. Benyahia, W. Wei, J. Li, MP-PIC simulation of CFB riser with EMMS-based drag model, Chemical Engineering Science 82 (2012) 104-113. [32] B.H. Xu, A.B. Yu, Numerical simulation of the gas-solid flow in a fluidized bed by combining discrete particle method with computational fluid dynamics, Chemical Engineering Science 52 (1997) 2785-2809. [33] Y. Tsuji, T. Kawaguchi, T. Tanaka, Discrete particle simulation of twodimensional fluidized bed, Powder Technology 77 (1993) 79-87. [34] X. Yuan, H. Li, M. Ye, Z. Liu, Study of the coke distribution in MTO fluidized bed reactor with MP-PIC approach, Canadian Journal of Chemical Engineering (2018).

40

[35] Y.-Q. Zhuang, X.-M. Chen, Z.-H. Luo, J. Xiao, CFD–DEM modeling of gas–solid flow and catalytic MTO reaction in a fluidized bed reactor, Computers & Chemical Engineering 60 (2014) 1-16. [36] W. Ge, Q. Chang, C. Li, J. Wang, Multiscale structures in particle–fluid systems: Characterization, modeling, and simulation, Chemical Engineering Science 198 (2019) 198-223. [37] M. Sakai, S. Koshizuka, Large-scale discrete element modeling in pneumatic conveying, Chemical Engineering Science 64 (2009) 533-539. [38] M.J. Andrews, P.J. O'Rourke, The multiphase particle-in-cell (MP-PIC) method for dense particulate flows, International Journal of Multiphase Flow 22 (1996) 379402. [39] K. Kuwagi, The similar particle assembly (SPA) model, an approach for largescale discrete element (DEM) simulation, Fluidization (2004) 243-250. [40] L. Lu, J. Xu, W. Ge, G. Gao, Y. Jiang, M. Zhao, X. Liu, J. Li, Computer virtual experiment on fluidized beds using a coarse-grained discrete particle method—EMMSDPM, Chemical Engineering Science 155 (2016) 314-337. [41] L. Lu, J. Xu, W. Ge, Y. Yue, X. Liu, J. Li, EMMS-based discrete particle method (EMMS–DPM) for simulation of gas–solid flows, Chemical Engineering Science 120 (2014) 67-87. [42] P.A. Cundall, O.D.L. Strack, A discrete numerical model for granular assemblies, geotechnique 29 (1979) 47-65.

41

[43] M.A. Mokhtar, K. Kuwagi, T. Takami, H. Hirano, M. Horio, Validation of the similar particle assembly (SPA) model for the fluidization of Geldart's group A and D particles, AIChE Journal 58 (2012) 87-98. [44] M. Sakai, M. Abe, Y. Shigeto, S. Mizutani, H. Takahashi, A. Viré, J.R. Percival, J. Xiang, C.C. Pain, Verification and validation of a coarse grain model of the DEM in a bubbling fluidized bed, Chemical Engineering Journal 244 (2014) 33-43. [45] Y. Liang, Y. Zhang, T. Li, C. Lu, A critical validation study on CPFD model in simulating gas–solid bubbling fluidized beds, Powder technology 263 (2014) 121-134. [46] A. Lanza, M.A. Islam, H. de Lasa, CPFD modeling and experimental validation of gas–solid flow in a down flow reactor, Computers & Chemical Engineering 90 (2016) 79-93. [47] X. Gao, T. Li, A. Sarkar, L. Lu, W.A. Rogers, Development and validation of an enhanced filtered drag model for simulating gas-solid fluidization of Geldart A particles in all flow regimes, Chemical Engineering Science 184 (2018) 33-51. [48] L. Lu, S. Benyahia, Advances in Coarse Discrete Particle Methods With Industrial Applications, Advances in Chemical Engineering, Elsevier2018, pp. 53-151. [49] A.S. Berrouk, C. Pornsilph, S.S. Bale, Y. Du, K. Nandakumar, Simulation of a large-scale FCC riser using a combination of MP-PIC and four-lump oil-cracking kinetic models, Energy & Fuels 31 (2017) 4758-4770. [50] K. Chu, J. Chen, A. Yu, Applicability of a coarse-grained CFD–DEM model on dense medium cyclone, Minerals Engineering 90 (2016) 43-54.

42

[51] C. Hu, K. Luo, S. Wang, L. Sun, J. Fan, Influences of operating parameters on the fluidized bed coal gasification process: A coarse-grained CFD-DEM study, Chemical Engineering Science 195 (2019) 693-706. [52] J. Li, M. Kwauk, Particle-fluid two-phase flow: the energy-minimization multiscale method, Metallurgical Industry Press1994. [53] J. Li, W. Huang, From Multiscale to Mesoscience: Addressing Mesoscales in Mesoregimes of Different Levels, Annual Review of Chemical and Biomolecular Engineering 9 (2018) 41-60. [54] J. Li, Mesoscale spatiotemporal structures: opportunities from challenges, National Science Review 4 (2017) 787-787. [55] J. Xu, H. Qi, X. Fang, L. Lu, W. Ge, X. Wang, M. Xu, F. Chen, X. He, J. Li, Quasireal-time simulation of rotating drum using discrete element method with parallel GPU computing, Particuology 9 (2011) 446-450. [56] J. Xu, X.W. Wang, X.F. He, Y. Ren, W. Ge, J.H. Li, Application of the Mole-8.5 supercomputer:Probing the whole influenza virion at the atomic level, Chinese Science Bulletin 56 (2011) 2114-2118. [57] F. Chen, W. Ge, L. Guo, X. He, B. Li, J. Li, X. Li, X. Wang, X. Yuan, Multi-scale HPC system for multi-scale discrete simulation—Development and application of a supercomputer with 1 Petaflops peak performance in single precision, Particuology 7 (2009) 332-335. [58] L. Ying, X. Yuan, M. Ye, Y. Cheng, X. Li, Z. Liu, A seven lumped kinetic model

43

for industrial catalyst in DMTO process, Chemical Engineering Research and Design 100 (2015) 179-191. [59] S. Ergun, Fluid flow through packed columns, Chemical Engineering Progress 48 (1952) 89-94. [60] D. Gidaspow, R. Bezburuah, J. Ding, Hydrodynamics of circulating fluidized beds: Kinetic theory approach, United States, 1991. [61] X. Liu, Y. Jiang, C. Liu, W. Wang, J. Li, Hydrodynamic modeling of gas–solid bubbling fluidization based on energy-minimization multiscale (EMMS) theory, Industrial & Engineering Chemistry Research 53 (2014) 2800-2810. [62] A. Nikolopoulos, A. Stroh, M. Zeneli, F. Alobaid, N. Nikolopoulos, J. Ströhle, S. Karellas, B. Epple, P. Grammelis, Numerical investigation and comparison of coarse grain CFD – DEM and TFM in the case of a 1MWth fluidized bed carbonator simulation, Chemical Engineering Science 163 (2017) 189-205. [63] A. Stroh, A. Daikeler, M. Nikku, J. May, F. Alobaid, M. von Bohnstein, J. Ströhle, B. Epple, Coarse grain 3D CFD-DEM simulation and validation with capacitance probe measurements in a circulating fluidized bed, Chemical Engineering Science 196 (2019) 37-53. [64] J. Wang, Y. Liu, EMMS-based Eulerian simulation on the hydrodynamics of a bubbling fluidized bed with FCC particles, Powder Technology 197 (2010) 241-246. [65] P.N. Rowe, Drag forces in a hydraulic model of a fluidised bed-Part II, Transactions of the Institution of Chemical Engineers 39 (1961) 175-180.

44

[66] J. Bruchmüller, B.G.M. van Wachem, S. Gu, K.H. Luo, R.C. Brown, Modeling the thermochemical degradation of biomass inside a fast pyrolysis fluidized bed reactor, AIChE Journal 58 (2012) 3030-3042. [67] S. Hu, X. Liu, N. Zhang, J. Li, W. Ge, W. Wang, Quantifying cluster dynamics to improve EMMS drag law and radial heterogeneity description in coupling with gassolid two-fluid method, Chemical Engineering Journal 307 (2017) 326-338. [68] T.O. Foundation, OpenFOAM v5 User Guide, (2017). [69] R. Barrett, M.W. Berry, T.F. Chan, J. Demmel, J. Donato, J. Dongarra, V. Eijkhout, R. Pozo, C. Romine, H. Van der Vorst, Templates for the solution of linear systems: building blocks for iterative methods, Siam1994. [70] H. van der Vorst, Bi-CGSTAB: A fast and smoothly converging variant of Bi-CG for the solution of nonsymmetric linear systems, SIAM Journal on Scientific and Statistical Computing 13 (1992) 631-644. [71] X. Wu, M.G. Abraha, R.G. Anthony, Methanol conversion on SAPO-34: reaction condition for fixed-bed reactor, Applied Catalysis A General 260 (2004) 63-69. [72] X. Yuan, H. Li, M. Ye, Z. Liu, Kinetic modeling of methanol to olefins process over SAPO-34 catalyst based on the dual-cycle reaction mechanism, AIChE Journal 65 (2019) 662-674. [73] Y.-Q. Zhuang, X. Gao, Y.-p. Zhu, Z.-h. Luo, CFD modeling of methanol to olefins process in a fixed-bed reactor, Powder Technology 221 (2012) 419-430. [74] N.E. Monsen, J.E. Cloern, L.V. Lucas, S.G. Monismith, A comment on the use of

45

flushing time, residence time, and age as transport time scales, Limnology and Oceanography 47 (2002) 1545-1553. [75] G. Zhao, X. Shi, Y. Wu, M. Wang, M. Zhang, J. Gao, X. Lan, 3D CFD simulation of gas-solids hydrodynamics and bubbles behaviors in empty and packed bubbling fluidized beds, Powder Technology 351 (2019) 1-15. [76] D. Kunii, O. Levenspiel, Fluidization engineering, Elsevier2013. [77] J.J. Zacca, J.A. Debling, W.H. Ray, Reactor residence-time distribution effects on the multistage polymerization of olefins—II. Polymer properties: bimodal polypropylene and linear low-density polyethylene, Chemical Engineering Science 52 (1997) 1941-1967. [78] K. Papadikis, S. Gu, A.V. Bridgwater, Geometrical Optimization of a Fast Pyrolysis Bubbling Fluidized Bed Reactor Using Comutational Fluid Dynamics, Energy Fuels 24 (2010) 5634-5651. [79] D. Christensen, J. Nijenhuis, J.R. van Ommen, M.O. Coppens, Influence of Distributed Secondary Gas Injection on the Performance of a Bubbling Fluidized-Bed Reactor, Industrial & Engineering Chemistry Research 47 (2008) 3601-3618. [80] T. Li, K. Pougatch, M. Salcudean, D. Grecov, Numerical simulation of single and multiple gas jets in bubbling fluidized beds, Chemical Engineering Science 64 (2009) 4884-4898. [81] L.T. Zhu, L. Xie, J. Xiao, Z.H. Luo, Filtered model for the cold-model gas–solid flow in a large-scale MTO fluidized bed reactor, Chemical Engineering Science 143

46

(2016) 369-383. [82] A.M. Saib, D.J. Moodley, I.M. Ciobîcă, M.M. Hauman, B.H. Sigwebela, C.J. Weststrate, J.W. Niemantsverdriet, J. van de Loosdrecht, Fundamental understanding of deactivation and regeneration of cobalt Fischer–Tropsch synthesis catalysts, Catalysis Today 154 (2010) 271-282. [83] D. Wu, Z. Gu, Y. Li, Attrition of catalyst particles in a laboratory-scale fluidizedbed reactor, Chemical Engineering Science 135 (2015) 431-440.

Declaration of interests

☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:

47

A pilot-scale MTO reactor is simulated for 8 hours of physical time using EMMSDPM. The perfect mixing of the catalyst particles in the whole reactor is observed. The correlation between the coke content distribution and age of the catalyst is directly obtained for the first time. In the local regions of the reactor, the poor mixing of catalyst particles is observed and the reaction fluctuates at large time-scales. Thus, the longtime particle-scale simulation would provide detailed information to characterize the behaviors in the critical regions of the reactor, and would be helpful to optimize the operating conditions and the reactor design.

48

1. Flow field and reactions in a pilot-scale MTO reactor is well predicted by EMMS-DPM. 2. Eight-hour particle-scale simulation reveals the perfect mixing of the catalyst globally, but not in some locations. 3. Dependence of coke deposition on catalyst age is directly obtained for the first time. 4. Reaction rate fluctuates significantly even at large time scales.

Table 3. Simulation parameters of the MTO reactor Item

Variables

Value

Diameter of the real particle, dp (μm)

97

Catalyst inventory in the reactor (kg)

9

Input catalyst rate (kg/h)

8.8

Density, ρp (kg/m3)

1500

Young's modulus, Y (Pa)

5×106

Thermal conductivity (W/(m K))

0.0454

Heat capacity (J/(kg K))

1220

Mass fraction of the coke on the input

0

Particle

49

catalyst (g/100gcat) Time step (s)

2.5×10-5

Reactor gas velocity, Ug1 (kg/h)

27

Lifting gas velocity, Ug2 (kg/h)

0.12

Gauge pressure at the exit (MPa)

0.024

Operating temperature (K)

738

Density, ρg (kg/m3)

0.4288

Gas

2.43×10Dynamic viscosity, μ (Pa s) 5

5×10-4

Time step (s)

Table 2. The properties of the MTO reaction process obtained in the simulation and experiment [23] Simulation

Experiment.

Error (%)

YCH4

0.0144

0.01643

12.08

YC2H4

0.4378

0.4303

1.74

YC3H6

0.3662

0.3457

5.92

YC3H8

0.0293

0.03788

22.73

YC4

0.1080

0.0913

18.29

50

YC5

0.0443

0.07836

43.46

wc (g/100gcat)

5.533

~6

~7.78

ηMeOH (%)

99.99

99.9444

0.05

S(C2H4+ C3H6)

80.40

77.603

3.6

C2H4/C3H6

1.196

1.2449

3.96

51