Gas Axial Dispersion Behavior and Gas Residence ... - ACS Publications

Sep 7, 2017 - behavior as well as the gas residence time by CFD simulation and theoretical ... the porous model32 in the Fluent 6.3.26 software.33 Tab...
0 downloads 12 Views 5MB Size
Article pubs.acs.org/IECR

Gas Axial Dispersion Behavior and Gas Residence Time in the Radial Flow Bed Ruojin Wang, Yiping Fan,* and Chunxi Lu* College of Chemical Engineering, State Key Laboratory of Heavy Oil Processing, China University of Petroleum (Beijing), Beijing 102249, China S Supporting Information *

ABSTRACT: The gas flow pattern in the radial flow bed is simulated by porous and transport models combined with a 6-lump kinetic reaction model (C4 hydrocarbon catalytic pyrolysis), which are verified by the published experimental results. Aside from the simulation method, the gas axial dispersion behavior is well explained by introducing the “common pressure pool” (CPP) and the “constant pressure drop lines” (CPDL) methods. The most severe gas axial dispersion behavior occurs at the turning point of (r − r1)/(r2 − r1) = 0.62. However, the gas flow pattern can still be roughly seen as a plug flow in whole bed for the small influences of the gas axial dispersion behavior. The appropriate gas residence time can then be determined by the plug flow model and the simulation method to achieve different production targets. reaction process simulation (the flow pattern is usually assumed to be plug flow24). However, although noticed by some researchers,21,25 the gas axial dispersion behavior (or the gas end effect) has not yet been quantified and fully investigated. It closely relates to the above undesirable phenomena and flow pattern. The gas axial dispersion behavior reflects the magnitude of gas flow rate entering the ends of the bed and appears in the radial flow bed inevitably due to the existence of the feed and the discharge influence zones.6,9 It is different from gas diffusion, which needs to be studied in the microscopic view and ignored under conditions of high Peclet number.6 In the feed influence zone, for instance, the so-called “solid seal” should be high enough in the case of gas short circuit.19,21 In this sense, the gas flow pattern is not exactly the plug flow as part of the gas phase flows toward the ends of the bed. As a result, the radial pressure drop reduces gradually at the ends of the bed,9,17 which indicates the appearance of the gas axial dispersion behavior. There are two major adverse effects of the gas axial dispersion behavior. First, the gas residence time varies with the axial position. In detail, compared with the middle section of the bed (the gas−solid contact zone), the gas tends to stay in the bed for a longer time at the ends of the bed (the feed and the discharge influence zones).6,13 In most reaction/filtration processes, the gas residence time should stay within a quite

1. INTRODUCTION The radial flow bed has some advantages, such as high gas flow rate and low pressure drop. It is used in many chemical engineering processes including gas filtration,1 naphtha reforming,2 and catalytic cracking.3,4 The radial flow bed has four different gas flow modes: inward Z, outward Z, inward π, and outward π. The gas phase flows through the bed centripetally in the inward mode and centrifugally in the outward mode. The directions of the gas flow rate in the distribution and collection channels are the same in the Z mode and opposite in the π mode.5 In order to improve the bed operating flexibility, the outward gas flow mode is then preferred when the pinning is more severe than the cavity under the same conditions.6 The gas flow mode is determined to be Z or π according to the axial uniformity of the gas flow rate.5,7−9 In this paper, the nonuniformity of the gas flow rate is anticipated to be controlled by introducing several gas inlets and outlets.10 The gas flow pattern in all mentioned modes can be further calculated based on that in the outward mode. In detail, the gas flow pattern in the inward mode can be derived from that in the outward mode. As for Z and π modes, the influences of the locations of the gas inlets and outlets should also be considered. In recent years, the flow pattern in the radial flow bed is investigated by many researchers, including some undesirable phenomena (e.g., cavity,6,11−14 pinning,6,14−16 air lock,6,14 gas nonuniformity distribution5,7−9), the simulation methods (e.g., single phase,7,9,17,18 Eulerian−Eulerian,13,19,20 Eulerian−Lagrangian21 methods), the structural optimization (e.g., trapezoidal/ conical bed,6,14 baffles/Johnson-net/louvers1,22,23), and the © XXXX American Chemical Society

Received: Revised: Accepted: Published: A

July 10, 2017 August 28, 2017 September 7, 2017 September 7, 2017 DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research narrow range; otherwise the anticipated products are difficult to obtain.26−28 The gas axial dispersion behavior may produce undesirable effects on the reaction/filtration efficiency of the bed. In addition, the gas axial dispersion behavior may cause the cavity to grow. The bubble fluidized bed can be taken as a similar example. In the bubble fluidized bed, the bubble forms at the bottom of the bed, and then goes forward to break up in the surface of the bed.29 The size of the bubble increases when the gas flow rate (or the radial pressure drop) rises. When the gas flow rate is high enough, the bubble fluidized bed even becomes a slug fluidized bed.30 Corresponding to the bubbles in the bubble fluidized bed, the cavity in a radial flow bed also experiences the recycle processes of formation, growth, and breakup.13,19 Moreover, the size of the cavity grows when the gas flow rate in the feed influence zone increases. On the other hand, the air lock phenomenon occurs when the gas flow rate in the feed influence zone is sufficiently high.14 Notably, at the ends of the bed, the gas flow rate in the feed influence zone rises only when the gas axial dispersion behavior aggravates it. Thus, the gas axial dispersion behavior relates closely to the cavity and the air lock. The nonuniformity of the gas flow rate always exists due to the gas axial dispersion behavior phenomenon. However, the influence of the gas axial dispersion behavior is limited to the small areas of the feed and the discharge influence zones where the gas axial dispersion behavior phenomenon appears. The gas flow can approximately be considered as plug flow in the bed. As mentioned above, the gas residence time is an important factor in the reaction.26−28 Nowadays, with the increasing demands for low-carbon olefins, it becomes popular to investigate the catalytic pyrolysis reaction of light hydrocarbons to low-carbon olefins. The intermediates are required in the series reaction such as the catalytic pyrolysis. Meanwhile, the chemical reaction is greatly affected by the temperature.3,4 The appropriate gas residence time in the catalytic pyrolysis reaction is then studied in this paper as well as the appropriate temperature. The gas axial dispersion behavior in a moving bed has close correlation to the uniformity of the gas residence time and the cavity/air lock. Besides, the gas residence time determines the reaction or filtration efficiency of the bed. The calculation of the appropriate gas residence time requires an assumption of plug flow, which requires investigation of the influence of the gas axial dispersion behavior on the uniformity of the gas flow rate. Thus, this paper puts emphasis on the gas axial dispersion behavior as well as the gas residence time by CFD simulation and theoretical analysis.

Figure 1. Schematic diagram of the radial flow bed: (a) geometry; (b) grid.

The gas phase is air at atmospheric pressure and room temperature. The gas flow rate is 183−641 m3/h. The solid phase is the supporter 3861 for the catalytic reforming reaction, whose diameter and density are 1.65 mm and 912 kg/m3.6 The gas phase passes through the bed horizontally, while the solid phase moves vertically.

3. SIMULATION 3.1. Mesh Generation. As seen in Figure 1b, the radial flow bed is modeled by the GAMBIT 2.4.6 software. The original rounds are replaced by the rectangles in both of the gas inlets and outlets for simplification. To improve the mesh quality, different blocks and mesh generation methods are used in the bed. In detail, the bed is made of hexagonal mesh (the catalyst bed and the collection channel), pentahedral mesh (the distribution channel), and tetrahedral mesh (others). 3.2. Simulation without Chemical Reaction. Considering that the magnitude order of the solid velocity is mm/s, the gas flow pattern in the moving bed differs little from that in the fixed bed.12,31 The gas flow characteristics are then simulated by the porous model32 in the Fluent 6.3.26 software.33 Table S1 shows the main governing equations. The boundary conditions are determined as follows: the mass inlet boundary condition is conducted in the gas inlets; the pressure outlet boundary condition is conducted in the gas outlets; the no-slip condition is conducted in the wall. 3.3. Simulation with Chemical Reaction. Some settings are changed when the catalytic pyrolysis reaction of light hydrocarbons to low-carbon olefins3,4 is considered. However, its reaction kinetic parameters (Table 1) can be directly used in the simulation models as they vary little with the size of the radial flow bed structure.30 It is because the gas flow pattern is nearly plug flow and the solid holdup remains almost unchanged in the radial flow bed. Moreover, the catalyst deactivation and coking deposition can be neglected in this catalytic pyrolysis reaction process.4 The reaction process involves a number of components and a series of complicated subreactions. It is very hard to simulate

2. RADIAL FLOW BED A semicylindrical radial flow bed6 is chosen to be studied in this paper. As shown in Figure 1, the bed is divided into three parts: the distribution channel, the catalyst bed, and the collection channel. These three parts are separated by the Johnson net, which prevents only the solid phase from passing through it. Several gas inlets and outlets are applied in the semicylindrical radial flow bed to obtain a more uniform gas flow rate distribution in the catalyst bed. As to the catalyst bed, the inside and outside radius are 54 mm and 246 mm. Its height is 1300 mm, which includes the height of the feed influence zone (200 mm), the gas−solid contact zone (1000 mm), and the discharge influence zone (100 mm).6 B

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Table 1. Reaction Kinetic Parametersa

a

reaction

k1

k2

k3

k4

k5

k6

k7

k8

k9

apparent activation energy, 106 J·kmol−1 pre-exponential factor, cm3·g−1·h−1

−28.5 16.25

0.2 1.82 × 104

43.5 2.47 × 106

59.5 5.55 × 107

−24.7 12.03

56.2 3.52 × 106

−27.6 207.47

−21.4 146.18

61.6 3.84 × 108

Reprinted in part with permission from ref 4. Copyright 2010 American Chemical Society.

the real reaction process. Thus, the lump kinetic model, one typical simplified model, is employed to simplify the calculation. The model involves the simplifications of the components, the reaction types, and the reaction mechanism. The process mainly talks about the cracking reaction of butene to propylene and ethylene. Some other reactions are also discussed involving butane, liquid + coke, and hydrogen + light alkanes. Nevertheless, many reactions are neglected such as the cracking reaction and some reversible reactions. Coke deposition is also excluded for the extremely low coke yield.4 As illustrated in Figure 2, only six components and nine important reactions

Table 2. Mass Fraction of Each Component in the Gas Inletsa feedstock butene propylene ethylene butane liquid + coke hydrogen + light alkane nitrogen

content of original,b content of whole gas mixture, wt % wt % 81.35 0.23 0.02 17.66 0.51 0.23

32.8764 0.0930 0.0081 7.1370 0.2061 0.0930 59.5865

a

Dilution ratio is 2.97. bReprinted in part with permission from ref 4. Copyright 2010 American Chemical Society.

The pressure in the near-wall area of the bed is then investigated in the simulation. Figure 3a indicates that the pressure drop distribution assumes almost the same profile under different mesh sizes. Considering the small size of the gas collection channel, the mesh size of 10 mm is then used to improve the mesh quality. 4.2. Turbulent Models. The turbulent pulsation should be taken into account when the turbulent flow occurs. However, the turbulent model requires large computational capacity due to its complexity. The Reynolds number is used to distinguish between laminar and turbulent flows. In the catalyst bed, the calculated Reynolds ρvd number is small enough (Re = μ , its value is less than

Figure 2. Reaction network of 6-lump kinetic model Reprinted with permission from ref 4. Copyright 2010 American Chemical Society.

remain to be considered. The reaction process can be simplified as eq 1. Furthermore, the reaction order is assumed to be of the first order, which is indicated in the right side of the mass balance equation (eq 2).

approximately 282) to use the laminar model. However, in the distribution and the collection channels, the Reynolds number is fairly high, and then the turbulent model is needed. Considering that high Reynolds number appears only in small areas, the laminar equation is anticipated to cause negligible calculation error. In order to choose a proper model, five common models are all studied including the laminar model, the standard κ−ε turbulent model, the RNG κ−ε turbulent model, the realizable κ−ε turbulent model, and the Reynolds stress model (RSM). 4.2.1. Equations of Turbulent Models. For the laminar flow, the gas flow pattern is modeled by the basic momentum conservation equation mentioned above. Every parameter has its own mean value. However, in the case of the turbulent flow, every parameter has mean and the fluctuating values (eq 3). The momentum conservation equation is then revised into eq 4.

ki

Ii → vijIj ⎛ ∂ρCi ⎞ ⎛ ∂C ⎞ ⎜ ⎟ + Gv ⎜ i ⎟ = ⎝ ∂t ⎠x ⎝ ∂x ⎠

(1)

∑ (±vjikjCjρ)

ρs ε

(2)

When the 6-lump kinetic model is included, some modifications have to be made in the simulation, such as the gas material and the operating condition. As shown in Table 2 and Tables S2 and S3, the mass fraction of each component is determined as well as the physical parameters of each component and the gas mixture parameters.33 Table S4 shows the governing equations in the condition of chemical reaction. In particular, there are six mass conservation equations for each component. The energy conservation equation has to be considered. The momentum conservation equation remains to be one, but the gas mixture parameters are introduced into the equation.

(3) v = v̅ + v′ ∂(ρg v) + ∇·(ρg vv) = −∇p + ∇·τ + ρg g + ∇·( − ρg v′v′) ∂t

4. MODEL VERIFICATION 4.1. Mesh Independent. Three different mesh sizes (5 mm, 10 mm, and 15 mm) were tested according to similar simulation cases.7,18,20,27 The corresponding total numbers of the mesh cells are 1 782 336, 372 738, and 94 677, respectively. The maximum equisize skew is approximately 0.84, 0.79, and 0.77, while the maximum aspect ratio is 3.8, 3.3, and 4.0. The measurement points are placed in the near-wall area of the bed in the experiment.6

(4)

Compared with the laminar flow, the Reynolds stresses are considered in the momentum conservation equation in the turbulent flow. In the κ−ε turbulent models, the Reynolds stresses are calculated via the Boussinesq hypothesis (eq 5). This hypothesis simplifies the calculation by the isotropic assumption. C

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 3. Radial pressure drop distribution under (a) three different mesh sizes and (b) turbulent models (near-wall).

The Reynolds stress model, one anisotropy model, considers the effects of every possible factor on the Reynolds stresses. However, it makes the calculation more complicated.33 ⎛ ∂v ⎛ ∂vj ⎞ ∂v ⎞ ⎟⎟ − 2 ⎜ρκ + μt k ⎟δij −ρg vi′vj′ = μt ⎜⎜ i + ∂xi ⎠ ∂xk ⎠ 3⎝ ⎝ ∂xj

(5)

4.2.2. Parameter Settings. In the gas inlets and outlets, both the hydraulic diameter (eq 6) and the turbulence intensity (eqs 7 and 8) are confirmed. L=

4A′ P

(6)

It = 0.16Re−0.125

(7)

ρvL μ

(8)

Figure 4. Radial pressure drop distribution under different solid holdup (the experimental radial pressure drop comes from the calculation of the pressure in x/l = 0.078, 0.781) Adapted with permission from ref 6. Copyright 2014 China university of petroleum (Beijing).

4.2.3. Different Turbulent Models. As shown in Figure 3b, the pressure drop distribution maintains almost the same configuration under different turbulent models. Particularly, the gas flow pattern in the catalyst bed, rather than the distribution and collection channels, is the focus in this paper. The gas will distribute uniformly and collect in the distribution and collection channels due to multiple gas inlets and outlets under different turbulent models. In this sense, the laminar model can be selected for maintaining accuracy and saving computational resources. 4.3. Solid Holdup and Simulation Validation without Chemical Reaction. According to the Ergun equation, the gas velocity has great effects on the radial pressure drop as well as the solid holdup. As indicated in Figure 4, the radial pressure drop remains almost unvaried in the gas−solid contact zone in the near-wall of the bed. It also grows with the rising of the solid holdup. When the solid holdup equals 0.58, the pressure drop distribution computed by the simulation method agrees well with the experimental data.6 4.4. Simulation Validation with Chemical Reaction. The gas flow rate is set as 549 m3/h, and the temperature is 500 °C. The mole dilution ratio of nitrogen to C4 hydrocarbons is 2.95 in the case RA,4 while it is 0 in the case RB. As illustrated in Table 3, the simulation results agree well with the experimental data4 including the conversion, the product yields, and the mass fraction of the main components. It validates the accuracy of simulation models. Moreover, the mole dilution ratio can be set as zero since it is almost independent from the conversion and the product yields.

5. RESULTS AND DISCUSSION 5.1. Different Gas Flow Rate. Figure 5 indicates that the radial pressure drop in the middle section of the bed is higher than that at the ends of the bed. It is also the criterion for dividing the catalyst bed into three parts including the feed influence, the gas−solid contact, and the discharge influence zones.6,13 In addition, the radial pressure drop increases with increase of the gas flow rate, which agrees well with Ergun’s equation.27 Moreover, the radial pressure drop seems to be zero in two axial positions. The lower one is placed at y/H = 0 in the discharge influence zone; while the upper one is located at y/H = 0.93 in the feed influence zone, near the bottom of the feed tube. That is to say, the gas is restricted to flow between the bottom of the feed tube and the bottom of the bed in the axial direction. Moreover, the radial pressure drop remains nearly invariant with the axial position in the gas−solid contact zone. However, it reduces gradually at the ends of the bed. Based on the assumption of plug flow, the pressure drop of the semicylindrical radial flow bed is computed by eqs 9 and 10.6,32 As seen in Table 4, the calculated radial pressure drop of the whole bed (Ergun equation) is slightly lower than the simulation and the experimental6 results. It is because the gas flow rate is assumed to be uniform in the axial direction as mentioned above, which is inconsistent with the facts. However, the nonuniformity caused by the gas axial dispersion behavior can be neglected for the computation of the radial pressure drop of the whole bed due to the small areas of the feed and the discharge influence zones. Thus, the radial pressure drop

Re =

D

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research Table 3. Conversion and Product Yields in the Simulation and the Experiment b

experiment (%) simulation of case RA (%) deviation between RA and expt (%) simulation of case RB (%) deviation between RB and expt (%) a

butenea

propylene

ethylene

butane

liquid + coke

hydrogen + light alkane

76.31 75.72 0.59 76.04 0.27

28.45 24.99 3.46 24.62 3.83

6.51 8.53 2.02 8.82 2.31

1.78 4.16 2.38 4.38 2.6

16.77 19.08 2.31 19.08 2.31

8.57 4.17 4.4 4.42 4.15

Note: is the conversion, others is the product yields. bReprinted in part with permission from ref 4. Copyright 2010 American Chemical Society.

changes in the radial direction. It is then discussed at the position of the turning point. The gas axial dispersion behavior in the radial flow bed needs to be investigated as well as the turning point. Before further discussion, the upward (in the feed influence zone) or the downward (in the discharge influence zone) is specified as the positive direction in this paper. The outward is also assigned as the positive direction. 5.2.1. Simulation Results. Seven virtual surfaces are selected in the bed to calculate its area-weighted average values, which include the whole bed, the feed influence, the gas−solid contact, and the discharge influence zones. Each surface has its own area-weighted average values in different zones. As illustrated in Figure 7a, the average gas radial velocity decreases from the center of the bed toward the outside of the bed. It grows only in the feed/discharge influence zone near the inside of the bed. The average gas radial velocity in the gas−solid contact zone is higher than those in the feed and the discharge influence zones. The average gas radial velocity computed by eq 11 agrees well with the simulation data in the whole bed. However, this equation cannot be used to compute the average gas radial velocity in other zones of the bed. In addition, the ratio of the average gas radial velocity in the feed/discharge influence zone to that in the gas−solid contact zone gets its peak in the turning point, which locates at approximately (r − r1)/(r2 − r1) = 0.6. The direction of the gas velocity changes from upward to downward at this position in the feed influence zone.

Figure 5. Radial pressure drop distribution under different gas flow rate.

of the whole bed has the same trend as the gas flow rate between the experiment, the simulation, and the theoretical (Ergun) results.6 The plug flow can roughly represent the gas flow pattern in the moving bed. C2ρg ⎛ Q ⎞2 dp μ Q ⎜ ⎟ = + 2 ⎝ πrH ⎠ α πrH dr

Δp =

(9)

C2ρg ⎛ Q ⎞2 ⎛ 1 r μ Q 1⎞ ⎜ ⎟ ⎜ ln 2 + − ⎟ 2 ⎝ πH ⎠ ⎝ r1 r2 ⎠ α πH r1

(10)

5.2. The Influence of the Gas Axial Dispersion Behavior in the Radial Direction. Figure 6 shows the gas flow streamlines in the feed influence zone. The pressure at the ends of the bed remains almost unchanged in the radial direction due to its small gas velocity and gas flow resistance. In addition, part of gas flows toward the feed influence zone. The gas streamlines can be simplified in the below theoretical analysis. As to the rectangular bed, it is obvious that the turning point at which the direction of the gas axial flow rate or the average gas axial velocity changes (e.g., from upward to downward) is in the middle of the bed. However, the gas velocity changes with the radial position in the semicylindrical radial flow bed. The turning point is not in the middle of the semicylindrical radial flow bed. The most severe gas axial dispersion behavior occurs at the turning point. The gas flow rate along the axial direction

vr =

Q πrH

(11)

The gas phase flows in both radial and axial directions in the feed and discharge influence zones. The gas axial velocity is then also used to explain the gas axial dispersion behavior. Figure 7b indicates that the average gas axial velocity almost equals zero in both the whole bed and the gas−solid contact zone. The average gas axial velocity presents x-axial symmetry between the feed and the discharge influence zones. It signals that the gas axial dispersion behaviors in these two zones are almost the same. The average gas axial velocity in the wall becomes zero as the result of the no-slip boundary condition in the simulation. The position of the turning point is likely to be

Table 4. Radial Pressure Drop of Whole Bed under Different Gas Flow Rate gas flow rate, Q (m3/h) radial pressure drop, Δp (kPa)

experimenta simulationb Ergun eq CPDL

183 0.14 0.13 0.11 0.12

275 0.22 0.22 0.18 0.20

336 0.33 0.32 0.23 0.29

458 0.46 0.43 0.35 0.40

549 0.57 0.56 0.46 0.52

641 0.70 0.70 0.58 0.65

a

Reprinted in part with permission from ref 6. Copyright 2014 China university of petroleum (Beijing). bNote: comes from the near-wall area in y/H = 0.5. E

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 6. Pressure contour and gas streamlines in the feed influence zone in (a) CPP and (b) CPDL methods.

Figure 7. Average gas velocity and gas radial velocity ratio: (a) radial velocity; (b) axial velocity.

at (r − r1)/(r2 − r1) = 0.62, which is exactly corresponding to the position determined by the gas radial velocity. More importantly, the average gas axial velocity near the inside of the bed has positive correlation to the gas force, which prevents the particles from moving downward to form the cavity in the feed influence zone. The average gas axial velocity is then preferred to be small near the inside of the bed. 5.2.2. Results of Theoretical Analysis. As mentioned above, the pressure keeps almost unchanged in the radial direction at the top and the bottom of the bed. As the gas axial dispersion behaviors are nearly the same in the feed and the discharge influence zones, the gas axial dispersion behavior is only discussed in the feed influence zone. The “common pressure pool” method (CPP) is put forward to calculate the average gas axial velocity. The detailed calculation process is shown in Figure 8. In particular, some viewpoints need to be clarified during this theoretical calculation process of CPP. First, the net gas flow rate entering the feed influence zone equals zero according to the mass conservation equation of the gas phase in the feed influence zone. Second, as illustrated in Figure 6a, the impetus for the formation of the gas axial velocity is assumed to be the pressure difference between the top and the bottom of the feed influence zone. The pressure can be calculated by eq 12 at the top of the feed influence zone. Meanwhile, the gas axial velocity is computed by eq 13. Third, the pressure at the top of the feed influence zone is invariant

with the radial position. In this way, the average gas axial velocity can be computed as well as the pressure at the top of feed influence zone. Δp|r =

C2ρg ⎛ Q ⎞2 ⎛ 1 r μ Q 1⎞ ⎜ ⎟ ⎜ ln 2 + − ⎟ 2 ⎝ πH ⎠ ⎝ r r2 ⎠ α πH r

±(Δp|r − Δp|gas ) HL

=

C2ρg ⎛ Q ′ ⎞2 μ Q′ ⎜ ⎟ + 2 ⎝ πrdr ⎠ α π r dr

(12)

(13)

The CPP method can be used to determine many parameters such as the pressure at the top of the feed influence zone, the average gas axial flow rate, the average gas axial velocity, and the position of the turning point. Notably, the air lock closely relates to the pressure at the top of the feed influence zone. When the pressure is high enough, the particles in the solid feed tube is hardly to move downward, and the air lock appears. The cavity aggravates when the particles are pushed to move upward, which is caused by the increased average gas axial flow rate near the cavity areas. The average gas axial velocity tends to zero at the turning point. Thus, the gas axial dispersion behavior is most severe at the turning point. The pressure is irrelevant to the height of the feed/discharge influence zone (HL) when the height is high enough (great than 100 mm) as well as the position of the turning point. F

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

rises when the gas flow rate increases. The position of the turning point is not affected by the gas flow rate. 5.3. The Influence of the Gas Axial Dispersion Behavior in the Axial Direction. The gas radial flow rate needs to be further analyzed to explain the gas axial dispersion behavior clearly. The average gas velocity in the turning point of (r − r1)/(r2 − r1) = 0.62 is explained. 5.3.1. Simulation Results. As the gas velocity changes little in the circumferential direction (Figure S1), the results in a full cylindrical bed are close to those in the semicylindrical bed as the influence of the no-slip boundary condition is limited. The gas velocity is small in the wall of the bed due to the no-slip boundary condition. Therefore, only the gas velocity in the near-wall is discussed. Based on above discussion, the gas flow rate remains unvaried in the axial direction in the gas−solid contact zone and decreases in other two zones, which is one significant mark of the appearance of the gas axial dispersion behavior. In addition, the gas flow rate in the near-wall is almost equivalent to that in the center of the bed regardless of the small discrepancy caused by the no-slip boundary condition of the wall. 5.3.2. Theoretical Analysis Results. The “constant pressure drop lines” method (CPDL) is proposed to calculate the average gas radial velocity in the axial direction. Figure 10 illustrates the theoretical calculation process of CPDL. Notably, the CPDL method requires some assumptions to simplify the calculation. First, the simplified streamlines of the gas phase are shown in Figure 6b. The value of α′, which ranges from 0° to 90°, is proportional to the axial distance between the streamline and the Johnson-net. Second, the gas velocity remains invariant with the radial position along each streamline. Finally, the pressure drop of each streamline is equal to each other. Under these assumptions, the radial pressure drop of the whole bed can be solved by eqs 14 and 15. The average gas radial velocity is then obtained according to the radial pressure drop of the whole bed.

Figure 8. Flowchart of the computation of the average gas axial velocity in the feed influence zone.

5.2.3. Comparison between Simulation and Theoretical Results. The theoretical results are validated with the simulation data in Figure 9. There exist small differences on the average gas axial velocity caused by the assumption of the formation of the gas axial velocity in the CPP method. However, the average gas axial velocity has the same trend with the radial position between the simulation and the theoretical results. Moreover, it grows with increase of the gas flow rate. The simulation results of the pressure at the top of the influence zone and the position of the turning point agree well with the theoretical data by the CPP method. The pressure at the top of the influence zone

C2ρg ⎛ Q ⎞2 dp μ Q ⎜ ⎟ = + dr /cos α′ 2 ⎝ πrH ⎠ α πrH

(14)

⎡μ Q C2ρg ⎛ Q ⎞2 ⎛ 1 1 ⎞⎤ r ⎜ ⎟ ⎜ ln + Δp = ⎢ − ⎟⎥ /cos α′ ⎢⎣ α πH 2 ⎝ πH ⎠ ⎝ r1 r ⎠⎥⎦ r1 (15)

Figure 9. Comparison of the gas axial dispersion between the simulation and the CPP theoretical results: (a) average gas axial velocity; (b) pressure at the top of the feed influence zone. G

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

obtained by the CPDL method is a little higher than the simulation result at the ends of the bed while it is lower than that in the gas−solid contact zone. It is because the gas radial velocity is assumed to be unchanged in the gas−solid contact zone based on the CPDL method. Meanwhile, the gas flow rate is uniform in the gas−solid contact zone, nevertheless, it is relatively small at the ends of the bed. Thus, the gas flow pattern can be roughly seen by plug flow in the whole bed sequentially. In addition, as illustrated in Table 4 the radial pressure drop of the whole bed deduced by the CPDL and the simulation methods agrees well with the experimental results.6 5.4. The Influence of Chemical Reaction on the Gas Axial Dispersion Behavior. The gas volumetric flow rate grows with increasing reaction time in the catalytic pyrolysis reaction process of C4 hydrocarbons. The following discussion focuses on the influence of the chemical reaction on the gas axial dispersion behavior. 5.4.1. Radial Direction. Figure 12 illustrates the variation of the average gas radial velocity in the radial direction. The average gas radial velocity grows when the chemical reaction is introduced. The cavity phenomenon is little affected since the gas velocity remains unchanged in the inside of the bed. Meanwhile, the variation of the gas velocity is not very large especially at the ends of the bed. Thus, the CPP method can also be used to compute the trend of the average gas axial velocity for simplification. The position of the turning point changes from 0.62 to 0.56. It is because the chemical reaction is carried out more completely at the ends of the bed than that in the middle section of the bed due to its small gas flow rate. Thus, the gas flow rate has greater growth at the ends of the bed. The gas has the tendency to flow into the middle section of the bed in the feed and the discharge influence zones, which results in the decrease of the location of the turning point. 5.4.2. Axial Direction. As the gas flow rate has not changed much in the catalytic pyrolysis reaction process, the CPDL method can also be used to compute the average gas radial velocity for simplification. 5.5. Gas Residence Time. The gas residence time is an important factor in the chemical reactions. In the catalytic pyrolysis reaction process of C4 hydrocarbons, the influence of the gas residence time on the mass fraction of each component needs to be further investigated. As the gas axial dispersion behavior only appears in local areas, its influence can be neglected for the calculation of the general parameters of the whole bed. The plug flow model is used to calculate the mass fraction of each component and the gas residence time. The appropriate gas residence time is computed by both the CFD simulation and the theoretical analysis (plug flow and 6-lump kinetic model) in this paper. 5.5.1. Mass Fraction Change of Each Component with the Gas Residence Time. In the catalytic pyrolysis reaction process of C4 hydrocarbons, most of the butene is converted into propylene and ethylene. The mass fraction of liquid + coke increases linearly in the radial direction as does that of the hydrogen + light alkanes. Moreover, the mass fraction of the propylene reaches its maximum in the middle section of the bed. The mass fraction contour of each component is given in Figure S2. It means that the appropriate gas residence time for maximizing propylene is in a narrow range. 5.5.2. Theoretical Calculation Results. The mathematical eq 16 for the plug flow is used to calculate the mass fraction change of each component.

Figure 10. Flowchart of the computation of the average gas radial velocity.

In the conical bed, the bed width varies with the axial position. Its radial pressure drop is then difficult to calculate. The radial pressure drop distribution can also be calculated based on the CPDL method. 5.3.3. Comparison between Simulation and Theoretical Results. As indicated in Figure 11, the gas radial velocity has the

Figure 11. Comparison of the gas radial velocity at the turning point in the axial direction between the simulation and theoretical results (near-wall).

same trend with the gas flow rate between the simulation and the theoretical results (CPDL). The gas radial velocity rises with increase of the gas flow rate. The gas radial velocity H

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research

Figure 12. Comparison of the average gas velocity between the nonreaction and the reaction processes: (a) average gas radial velocity; (b) average gas axial velocity.

⎛ ∂ρg Ci ⎞ ⎛ ∂C ⎞ ⎜⎜ ⎟⎟ + Gv ⎜ i ⎟ = ⎝ ∂x ⎠ ⎝ ∂t ⎠x

∑ (±vjkjCjρg )

ρp ε

(16)

When the bed reaches a stable state, the mass fraction change of each component is computed by eqs 17−19). ρp dCi ,r = ∑ (±vjkjCj ,r) dr πrHugε (17) dYi ,r dr

=

ug, r =

ρp πrHug, r ε

∑ (±kjYj ,r)

(18)

Q 0M̅ 0 πrHM r

Figure 13. Comparison of the mass fraction of each component between the simulation and the theoretical (plug flow and 6-lump kinetic models) results (near-wall).

(19)

As mentioned above, the bed is also divided into a number of concentric cylindrical elements with a thickness dr (Figure 6). In addition, according to the original mass fraction and the mass fraction variation of each component, the mass fraction of each component in any radial position is obtained by eq 20. Meanwhile, the gas residence time is obtained by eqs 21 and 22. Yi ,r + dr = Yi ,r +

tr + dr = tr + r2

t=

∫r1

dYi ,r dr

dr

dr ug,r

rdr = ug,r

The gas residence time is greatly influenced by both the bed structure and the gas flow rate. In order to obtain different gas residence times in simulation, the gas flow rate rather than the bed structure is changed for simplification. The appropriate gas residence time should be determined for obtaining different production targets (Figure S3), which is discussed as follows. (1) Maximizing propylene. The mass fraction of the propylene reaches its maximum in a narrow range of gas residence time. As for maximizing the propylene yield, the appropriate range of the gas residence time is computed to be about 0.089−0.151 s (1032−1799 m3/h) by theoretical method (plug flow and 6-lump kinetic models) and 0.089−0.177 s (875−1799 m3/h) by simulation method. In particular, the mass fraction obtained by the simulation method is almost the same as the theoretical calculation result. (2) Maximizing ethylene. For maximizing the ethylene yield, the appropriate range of the gas residence time is about 0.467−1.209 s (143−232 m3/h) by the theoretical method. (3) Maximizing whole profits. The optimized function is needed for maximizing the whole profits. For example, the optimized function becomes ∑(Yi − Yi0)pri when the profit of each component is given as pri. The appropriate gas residence time can then be roughly determined by the theoretical calculation (plug flow and 6-lump kinetic models) method and verified by the simulation method.

(20)

(21) r2

∫r1

πrHrdr πH 2 = (r2 − r12) Q 2Q

(22)

5.5.3. Comparison between Simulation and Theoretical Calculation Results. As indicated in Figure 13, the theoretical calculation method, the combination of the plug flow and the 6-lump kinetic models, is proved to be accurate in modeling the mass fraction of each component under different radial position (or the gas residence time). 5.5.4. Determination of Appropriate Gas Residence Time. Figure 13 also illustrates the mass fraction of each component under different gas residence times. The mass fraction of most components, for example, butane, propylene, ethylene, and liquid + coke, reaches its peak in a given gas residence time. Moreover, propylene reaches the maximum value of appropriately 28% at 500 °C, which has good accordance with the experimental result.3 I

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research



5.5.5. Effects of the Temperature on the Appropriate Gas Residence Time. In the catalytic pyrolysis reaction process of C4 hydrocarbons, the temperature greatly affects the reaction process as well as the appropriate gas residence time. The mass fraction of each component is then investigated under different temperature. Take the propylene for example. The appropriate gas residence time reaches its peak at about 650 °C for maximizing propylene as well as propylene. High temperature accelerates the transition of propylene to hydrogen + light alkanes. Low temperature is undesirable for the conversion of butene to propylene. Besides, the difference of the maximum mass fraction of the propylene is less than 1% when the temperature ranges from 590 to 720 °C computed by the theoretical calculation (plug flow and 6-lump kinetic models) method and from 550 to 650 °C in the experiment.3 The appropriate gas residence time also remains unvaried in this range for maximizing propylene. More details can be seen in Figure S4. The propylene selectivity reaches the highest 56% at about 600 °C obtained by the theoretical calculation (plug flow and 6-lump kinetic models) method and 58.87% at about 570 °C in the experiment.3 These also validate the accuracy of the 6-lump kinetic model.

Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.iecr.7b02817. Tables with the governing equations and the gas mixture parameters and figures with the results excluded from the full paper including the gas velocity magnitude contour and flow streamlines, mass fraction contour and mass fraction change of components (PDF)



AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. ORCID

Ruojin Wang: 0000-0001-9750-9783 Chunxi Lu: 0000-0002-9803-9119 Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS This work was supported by NKBRDP (the National Key Basic Research Development Program) of China (973, 2012CB215000).

6. CONCLUSIONS The porous model is applied to simulate the gas flow pattern in the radial flow bed combined with the transport and the 6-lump kinetic models. The simulation models are verified by the experimental data considering the radial pressure drop, the conversion and the production yields. According to above discussion, the laminar model is enough for simulating the gas flow pattern in the catalyst bed here. The gas axial dispersion behavior is first quantified by the average gas velocity, which is determined by the CFD simulation and theoretical calculation (CPP and CPDL) methods. When the gas flow rate increases, the gas axial dispersion behavior aggravates for increasing the average gas velocity. However, the turning point remains almost unchanged with the gas flow rate; the location of (r − r1)/(r2 − r1) = 0.62. In addition, the effect of the reaction on the gas axial dispersion behavior can be neglected if the gas flow rate has not changed much in this reaction process. Besides, the theoretical calculation results can also be used for other purposes. For instance, the pressure at the top of the feed influence zone clarifies the air lock phenomenon; the velocity near the upstream of the bed affects the cavity; the gas flow rate distribution in the axial direction denotes the severity of the gas axial dispersion behavior. The gas flow pattern can also be computed in other radial flow beds for structural optimization purposes, for example, the trapezoidal and the conical bed. The appropriate gas residence time should be determined for different production targets such as maximizing a single component and the whole profits. The appropriate gas residence time can then be computed by the theoretical analysis (plug flow and 6-lump kinetic models) and verified by the CFD simulation. Besides, the appropriate gas residence time is greatly affected by the temperature as well as the mass fraction of each component. The intermediates in the chemical reaction reach their maximum at appropriate temperature. The gas flow pattern needs further calculation based on the ideas proposed when the gas flow rate changes considerably in the reaction processes.

NOMENCLATURE A pre-exponential factor, s−1 A′ cross sectional area, m2 C2 inertial resistance factor, m−1 Ci molecular weight of lump i in one gram of the gas, mol/g Cp species heat capacity, J/(kg·k) D dispersion coefficient, m2/s d particle diameter, mm e restitution coefficient between particles, unitless E activation energy for reaction, J/kmol g gravitational acceleration, m2/s Gv mass velocity in the cross section of the reactor, kg/(m2·s) h enthalpy, J/kg H height of the catalyst bed, mm H′ solid-seal height, mm HF height of the feed influence zone, mm HJ height of the Johnson net, mm HL height of the feed/discharge influence zone, mm Ii species i, unitless It turbulence intensity, unitless identity matrix, unitless I J dispersion flux, kg/(m2·s) ki, kj reaction rate constant of species, s−1 kt thermal conductivity, W/(m·k) km mass diffusivity, m2/s Keff effective thermal conductivity, W/(m·k) L hydraulic diameter, m Le Lewis number, unitless M molecular weight, kg/kmol P wetted perimeter, m p pressure, kPa pr profit of component, unitless Q gas flow rate, m3/s Q′ assumed gas flow rate, m3/s r radius, mm r̅ gas velocity ratio, unitless J

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research R Ri Re S T u u̅gz v v̅ v′ vij y Y Si t

(10) Ahn, Y. K.; Chen, H. C.; Fan, L. T.; Wan, C. G. A modified moving bed grain dryer. Ind. Eng. Chem. Process Des. Dev. 1964, 3, 96− 100. (11) Song, X. Q.; Jin, Y.; Yu, Z. Q. Influence of outward radial gas flow on particle movement in an annular moving bed. Powder Technol. 1994, 79, 247−256. (12) Chen, Y. H.; Zhu, X. D.; Wu, Y. Q.; Dong, Y. H.; Zhu, Z. B. Phenomenon of cavity and its size in rectangular moving-bed. J. Chem. Ind. Eng. 2006, 57, 731−737 (in Chinese). (13) Wang, R. J.; Xu, J.; Long, W. Y.; Fan, Y. P.; Lu, C. X. Experiment and simulation of dense gas-solid two-phase flow in a rectangular cross-flow moving bed. Chin. J. Process. Eng. 2015, 361−367 (in Chinese). (14) Long, W. Y.; Xu, J.; Fan, Y. P.; Lu, C. X. Pinning and cavity in two types of cross-flow moving beds. Powder Technol. 2015, 269, 66− 74. (15) Pilcher, K. A.; Bridgwater, J. Pinning in a rectangular moving bed reactor with gas cross-flow. Chem. Eng. Sci. 1990, 45, 2535−2542. (16) Doyle, F. J.; Jackson, R.; Ginestra, J. C. The phenomenon of pinning in an annular moving bed reactor with crossflow of gas. Chem. Eng. Sci. 1986, 41, 1485−1495. (17) Takahashi, H.; Ohno, O.; Takeuchi, T. Analysis of gas flow in a solid-gas cross-flow moving bed by the finite element method. Kagaku Kogaku Ronbunshu 1990, 16, 116−122. (18) Gao, X.; Zhu, Y. P.; Luo, Z. H. CFD modeling of gas flow in porous medium and catalytic coupling reaction from carbon monoxide to diethyl oxalate in fixed-bed reactors. Chem. Eng. Sci. 2011, 66, 6028−6038. (19) He, Y. J.; Xiao, F. Z.; Luo, Z. H. Numerical modeling of the cavity phenomenon and its elimination way in rectangular radial moving bed reactor. Powder Technol. 2015, 274, 28−36. (20) Xiao, F. Z.; Luo, Z. H. A Two-Phase CFD Modeling Approach to Investigate the Flow Characteristics in Radial Flow Moving-Bed Reactors. Int. J. Chem. React. Eng. 2014, 12, 497−512. (21) Lu, C. X.; Long, W. Y.; Lin, Q. Z. Numerical simulation of gassolid flow characteristics in radial moving bed. J. Chem. Eng. Chinese Univ. 2015, 29, 49−57 (in Chinese). (22) Hsiau, S. S.; Smid, J.; Chyou, Y. P.; Liu, T. C.; Huang, T. C.; Hsu, C. J. Impact of flow-corrective insert on flow patterns in twodimensional moving bed. Chem. Eng. Process. 2013, 73, 7−15. (23) Chen, Y. H.; Zhu, X. D.; Wu, Y. Q.; Zhu, Z. B. Investigation of the effect of a dividing wall in a moving bed. Chem. Eng. Technol. 2007, 30, 1028−1035. (24) Rahimpour, M. R.; Iranshahi, D.; Bahmanpour, A. M. Dynamic optimization of a multi-stage spherical, radial flow reactor for the naphtha reforming process in the presence of catalyst deactivation using differential evolution (DE) method. Int. J. Hydrogen Energy 2010, 35, 7498−7511. (25) Mu, Z. Z.; Wang, J. F.; Wang, T. F.; Jin, Y. Optimum design of radial flow moving-bed reactors based on a mathematical hydrodynamic model. Chem. Eng. Process. 2003, 42, 409−417. (26) Hsu, C. J.; Hsiau, S. S. A study of filtration performance in a cross-flow moving granular bed filter: the influence of gas flow uniformity. Powder Technol. 2015, 274, 20−27. (27) Xiao, F. Z.; Chen, H. Y.; Luo, Z. H. Numerical evaluation and improvement efficiency of radial flow moving-bed reactors for catalytic pyrolysis of light hydrocarbons to low carbon olefins. Can. J. Chem. Eng. 2015, 93, 1033−1043. (28) Sadrameli, S. M. Thermal/catalytic cracking of liquid hydrocarbons for the production of olefins: A state-of-the-art review II: Catalytic cracking review. Fuel 2016, 173, 285−297. (29) Sun, X. D.; Kim, S.; Ishii, M.; Beus, S. G. Modeling of bubble coalescence and disintegration in confined upward two-phase flow. Nucl. Eng. Des. 2004, 230, 3−26. (30) Hoomans, B. P. B.; Kuipers, J. A. M.; Briels, W. J.; Van Swaaij, W. P. M. Discrete particle simulation of bubble and slug formation in a two-dimensional gas-fluidised bed: A hard-sphere approach. Chem. Eng. Sci. 1996, 51, 99−118.

universal gas constant, J/(kmol·K) net rate of production of species i by chemical reaction, kg/(m3·s) Reynolds number, unitless momentum source term, kg/(m3·s) temperature, K superficial velocity, m/s Average gas axial velocity, m/s instant velocity, m/s mean velocity, m/s fluctuating velocity, m/s Stoich coefficient for the reaction of species i to species j axial position, mm mass fraction, % momentum source term, kg/(m2·s2) time step, s

Greek Symbols

α permeability resistance factor, m2 ρg, ρp gas and particle density, kg/m3 μ dynamic viscosity, kg/(m·s) μt turbulent viscosity, kg/(m·s) shear stress, N/m2 τ κ turbulence kinetic energy, m2/s2 ε bed porosity, unitless Subscripts

i, j opt r 0 1, 2

lump species, i and j are in the same reaction process optional value radial position original value inside and outside of the catalyst bed

Abbreviations

CFD computational fluid dynamics CPP common pressure pool CPDL constant pressure drop lines



REFERENCES

(1) Hsiau, S. S.; Smid, J.; Tsai, S. A.; Tzeng, C. C.; Yu, Y. J. Flow of filter granules in moving granular beds with louvers and sublouvers. Chem. Eng. Process. 2008, 47, 2084−2097. (2) Iranshahi, D.; Karimi, M.; Amiri, S.; Jafari, M.; Rafiei, R.; Rahimpour, M. R. Modeling of naphtha reforming unit applying detailed description of kinetic in continuous catalytic regeneration process. Chem. Eng. Res. Des. 2014, 92, 1704−1727. (3) Li, L.; Gao, J. S.; Xu, C. M.; Meng, X. H. Reaction behaviors and mechanisms of catalytic pyrolysis of C4 hydrocarbons. Chem. Eng. J. 2006, 116, 155−161. (4) Meng, X. H.; Xu, C. M.; Li, L.; Gao, J. S. Kinetic study of catalytic pyrolysis of C4 hydrocarbons on a modified ZSM-5 zeolite catalyst. Energy Fuels 2010, 24, 6233−6238. (5) Song, X. Q.; Wang, Z. W.; Jin, Y.; Gong, M. S. Investigations on hydrodynamics of radial flow moving bed reactors. Chem. Eng. Technol. 1993, 16, 383−388. (6) Long, W. Y. Hydrodynamic characteristics of cross-flow moving bed. Ph.D. Thesis, College of Chemical Engineering, China University of Petroleum, Beijing, 2014 (in Chinese). (7) Kareeri, A. A.; Zughbi, H. D.; Al-Ali, H. H. Simulation of flow distribution in radial flow reactors. Ind. Eng. Chem. Res. 2006, 45, 2862−2874. (8) Zhang, C. F.; Zhu, Z. B.; Xu, M. S.; Zhu, B. C. Design for uniform distribution of stream in radial reactor. J. Chem. In. 1979, 1, 67−90. (9) Wang, J. F.; Jing, S.; Wang, T. F.; Jin, Y.; Ma, X. Q.; Gao, L. P. Mathematical modeling and flow field characteristics of radial flow moving-bed reactors. J. Chem. Eng. Chinese Univ. 1999, 13, 435−441 (in Chinese). K

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Article

Industrial & Engineering Chemistry Research (31) Chen, Y. H.; Zhu, X. D.; Wu, Y. Q.; Zhu, Z. B. The investigation of pressure drop in a moving bed. Chin. J. Process. Eng. 2006, 697−702 (in Chinese). (32) Ergun, S. Fluid flow through packed columns. Chem. Eng. Prog. 1952, 48, 89−94. (33) Fluent Inc. Fluent 6.3 User’s Guide 2006.

L

DOI: 10.1021/acs.iecr.7b02817 Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX