Jonke, A. A., Levenson, M., Levitz, N. M., Steindler, M. J., Vogel, R. C., Nucleonics, 25 (5), 58-63 (1967a). Jonke, A. A,, Levenson, M., Steunenberg, R. K., Vogel, R. C., Symposium Proceedings: Plutonium as a Reactor Fuel, International Atomic Energy Agency, Vienna, Austria, 196713, pp 537. Levenspiel, O., “Chemical Reaction Engineering,” Wiley, New York, N. Y., 1962, p 348.. Levitz, N. M., Vogel, G. J., Carls, E. L., USAEC Report ANL-7225, 102, Argonne, Ill., 1966. Mukherjee, S. P., Doraiswamy, L. K., Brit. Chem. Eng., 12 (l),70 (1967). Ramaswami, D., Anastasia, L. J., Levitz, N. M., Mecham, W. J., Jonke, A. A., USAEC Report ANL-7339, Argonne, Ill., 1968. Rand, M. H., At. Energy Reu., 4, Special Issue No. 1, 7 (1966). Rand, M. H., Kubaschewski, O., “The Thermochemical Properties of Uranium Compounds,” Wiley, New York, N. Y., 1963.
Stein, L., J . Phys. Chem., 66, 288 (1962). Steindler, M. J., Steidl, D. V., T!SAEC Report ANL5759, Argonne, Ill., 1957. Steindler, M. J., Steidl, D . V., Steunenberg, R. K., Nucl. Sci. Eng., 6, 333 (1959). Steindler, M. J., USAEC Report ANL-6753, Argonne, Ill., 1963. Trevorrow, L. E., Shinn, W. A., Steunenberg, R. K., J . Phys. Chem., 65, 398 (1961). Vandenbussche, G., Commissariat a 1’Energie Atomique, Report No. CEA-R-2859, Paris, France, 1966. Wen, C. Y., Ind. Eng. Chem., 60 (9), 34 (1968).
RECEIVED for review June 6, 1969 ACCEPTED November 30, 1970. Work performed under the auspices of the U. S. Atomic Energy Commission.
Modeling of the Air Flow Pattern in a Countercurrent Spray-Drying Tower John R. Paris’, Phillip N. ROSS,Jr., Soli P. Dastur, and Robert 1. Morris’ Miami Valley Laboratories, The Procter & Gamble Co., Cincinnati, Ohio 45239 A model for the air flow pattern in a synthetic detergent spray-drying unit i s developed by by of by
deconvolution of residence-time distribution data. Experimental data are obtained injecting a pulse of helium a t the unit air inlet and monitoring the concentration tracer in the effluent gas by on-line mass spectroscopy. Deconvolution i s performed several analog and digital techniques which are compared.
T h i s paper discusses the development and interpretation of a model for the air flow in a spray-drying tower used for the commercial production of detergents. The parameters of the model apply to actual operating conditions. There are two fundamental approaches to the mathematical modeling of the flow of a fluid in a vessel. The first is the now-classical transport phenomena method which consists of writing mass, momentum, and energy balances over some volume which, depending on the degree of details in which the flow pattern can be described, may vary from an elementary volume to the whole vessel. The equations in their differential or integral form are then solved to predict the flow patterns of the fluid. These models are very desirable when they can be obtained, though, unfortunately, they are restricted to relatively simple systems. When the nature of the flow is complex, particularly if its pattern varies from region to region inside the vessel or if the relationship between the transport parameters and the geometrical coordinates cannot be expressed continuously, one must turn to the second approach, that is, stochastic modeling. I n these instances the flowing material is considered as a collection of countable entities (particles of solid, small elements of liquid or gas) which
’ To whom correspondence should be addressed. ‘Present address, Kraftco Corp.. Glenview, Ill. 60025
are treated statistically. Their flow is expressed by various age probability functions, such as the age distribution of the elements of fluid in the vessel at a given time, or the age distribution of the elements of fluid leaving the vessel at a given time, and known as the residencetime distribution function ( R T D F or the f-function). More complete definitions of these functions and others have been given by Levenspiel (1962) as well as by Himmelblau and Bischoff (1968). The age spread of elements of fluid leaving the vessel is indicative of the gross flow patterns in the vessel (short-circuiting, channeling, back-mixing, etc.). Danckwerts (1953) developed an analysis of the R T D F curve by introducing the concept of longitudinal dispersion where the spread in residence time is related to an effective diffusion process counter to the bulk flow of fluid. The equivalence between the longitudinel dispersion model and a series of equally sized. perfectly mixed tanks was later established (Levenspiel and Smith, 1957; Van der Laan, 1958), adding flexibility to the mathematical treatment. Dispersion models have been widely used in connection with tubular reactors, packed beds, and packed columns (Levenspiel and Bischoff, 1963), but they have not been successful in representing complex flow patterns where micromixing and macromixing occur simultaneously. Attempts have been made to describe such systems by empirical models, generally referred to as “mixed” models Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2 , 1971
157
AirOut
,(:);iter , I
Air In
Principles of the RTDF Analysis
SprayDirection
I
-I. -Air
Flow direction
in, t h e Plenum
inlection Equipment Granules
Figure 1. Schematic of spray-drying tower with tracer injection and sampling points
because they are formed by a network of ideal flow components arranged to reproduce the R T D F of the real unit. By ideal flow component a system is meant which has a specifically characterized flow pattern, usually taken as one of the following: plug flow (no mixing); “perfectly mixed” flow (complete mixing on a macroscopic scale); or dispersed plug flow (plug flow plus axial diffusion). The unsteady material balance for each of the ideal flow components is a linear differential equation and the mixed model results in a system of linear differential equations. The determination of a mixed model from the R T D F alone is not unique. Several models having different numbers and arrangements of ideal flow components may produce the same residence-time distribution. Hopefully, one of the models may be sufficiently close t o the physical system that it can be used for simulation purposes. Its validity must be checked by independent physical measurements such as velocity and temperature distributions. This paper presents the development of such a model for the air flow in a spray-drying tower. Comprehensive studies of mixed models, of their construction and utilization, have been the object of a number of publications (Van Deemter, 1961; Wen and Chung, 1967; Levelspiel and Bischoff. 1963; Adler et al., 1963). The spray-drying tower used in this work (Figure 1) was described previously by Chaloud et al. (1957); it is also similar in design and dimensions to the unit studied by Place et al. (19591, its overall height being 80 ft and the diameter of‘ the cylindrical section, 20 ft. Velocity traverses were determined by Chaloud et al. (1957), but in the absence of spray, a condition which has been found since to dramatically alter the pattern of flow of the air. Place et al. (1959) conducted tracer experiments with several injection points within the tower and one sampling point a t the outlet of the tower. These experiments were also performed in the absence of spray. Furthermore, the R T D F of the sampling and analyzing equipment was neglected, a simplification which, in light of the present work, must be questioned. The results from these two earlier experiments were very useful, however, as guidelines for the selection of a feasible model among several theoretically possible ones, and for the physical interpretation of the selected model. 158
Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
The residence-time distribution of a fluid in a vessel can be determined by injecting a pulse, 6 ( t ) , of tracer into the inlet stream and measuring the tracer concentration in the effluent stream as a function of time. By definition, the normalized tracer concentration curve is the residence-time distribution function of the fluid in the vessel, it is denoted f r ( t ) . I n practice, the equipment used to sample and analyze the effluent (for convenience, this equipment will be referred t o hereafter as “detector”) is characterized by an R T D F of its own, and what is actually measured is the overall R T D F for the vessel and the detector combined. If the response of this total system-Le., vessel plus detector-is linear, Duhamel’s theorem states that its response is independent of the serial order of the elements comprising it. The linearity of a real system can be experimentally checked since the normalized response curve must be independent of the quantity of tracer injected as a pulse. Figure 2 illustrates the problem solution for a linear system. Diagram a represents the resolved problem, that is, the R T D F of the vessel-Le., spray-drying tower in the case a t hand-fr(t) is the response to an input 6 ( t ) . Diagram b shows that in the actual experiment f ~ , ) ( t )
Figure 2. Relationships between the RTDF’s of tower, detector, and overall system
Figure 3. (detector)
Helium sampling
1. 2.
and
analyzing
equipment
Sampling probe, %J-in. id stoinless steel tubing Injection point for determination of fr, 3. Counter-flow stotion, to prevent plugging of the sampling probe between runs 30. Three-way valve: 1, during runs; 2, between runs Drying station: cold traps; liquid NJ 4. 5 . Vapor superheating station, heating tope coiled around %-in. copper tubing. Beyond this point, no water may condense 6. Filter, pore diameter = 3 microns 7. High vacuum sampling probe (C.E.C. adjustoble probe 24-042) 8. Leak detector (C.E.C.-24-120B)
9.
Pumping station
is also the tower response to the R T D F of the detector ) itself. Since this function can also be measured, f ~ ( t can be obtained by application of the deconvolution theorem t o f n and f r u
A convenient and established way of solving the deconvolution problem is t o assume a flow model with a sufficient number of differential equations and parameters to describe the system response. With the measured f D ( t ) as the input function, the parameters of the model are varied until the solution to the model equations matches the measured overall R T D F . The final set of parameters and differential equations can then be used to produce the true system response, f T ( t ) ,to a pulse input. As pointed out by Moser and Cupit (1966), the model used for the deconvolution may be arbitrary; that is, it does not have to represent the physical flow through the system. For the purposes of this work we will use a model believed t o represent the physical flow, thus solving two problems simultaneously: the deconvolution of the tower response, and the evaluation of the parameters in a flow model. Experimental
A diagram of the spray-drying tower showing the positioning of the tracer injection and detection equipment is given in Figure 1. The input signal was generated by the discharge of a pressure vessel filled with helium. The discharge was realized through a nozzle designed t o spread the helium across one of the openings from the plenum chamber into the tower a t an angle with the louvres so as to attenuate the jet effect and to generate good mixing between the helium and the air entering the tower a t this opening. The discharge system was designed to inject a quantity of helium sufficient for adequate detection in the output signal within as short a time as possible while delivering a jet linear velocity of the same order of magnitude as the velocity of the air a t the same point. The pulse was released by the action of an air-operated valve actuated by a solenoid. The discharge time was approximately 0.35 sec, well within the range where the assumption of a true impulse could be accepted since the mean residence of the overall system (tower plus detector) was about 60 sec. The concentration of helium in the air leaving the tower was continuously monitored by a helium analyzer (C.E.C.24-120B-Helium Leak Detector) capable of measuring helium t o the ppm level. The air was sampled in two steps (Figure 3). A very small fraction of the gas leaving the tower (-0.6 scfm) was first sampled through a 1 2 in. id probe located at the center of the outlet duct. Part of the primary gas was then diverted toward the leak detector through a high vacuum probe, that is, the precision needle valve furnished with the helium analyzer and normally used during leak detection. The design of this valve eliminated the dead space on the high pressure side, present in conventional high vacuum throttling valves, which could not be tolerated for this experiment. The air leaving the tower contained detergent dust and water vapor which condenses as the gas is cooled. The deposition and accumulation of the liquid and solid material a t various points of the equipment would be undesirable. Furthermore, the most essential parts of the system
x 0
\-
TOWER PLUS DETECTOR
10 20 30 40 50 60 70 80 90 100 110 120130140 TIME (SEC 1
Figure 4. Actual response curves,
fD(t)
and
fTD(t)
(high vacuum probe, helium sensing unit) could not function properly in the presence of dust or condensing vapors. Purification of the gas was achieved by means of a series of baffled cold traps followed by a micropore filter. This set-up allowed for a constant flow rate and pressure in the detector over a sufficient period of time for the recording of several experiments. However, it resulted in a relatively large mean residence time for the detector ( - 40 sec) compared to the spray-drying tower. The helium analyzer itself had a mean residence time of only 3.3 sec. No absorption and later release of the helium by any part of the equipment was observed during experiments. The residence-time distribution of the detector was recorded under actual sampling conditions during the experiment by injecting a pulse of helium (actual duration -0.15 sec) at the base of the sampling probe while the spray-drying unit and sampling pumps were operating at their normal steady state. The output of the helium analyzer was recorded on a Brush Mark 240 recorder, together with injection time. and pressure in the helium analyzer, which was to remain strictly constant. Helium discharges a t the inlet to the tower and a t the inlet to the sampler were triggered by a common mechanism which allowed for recording of a series of alternate residence-time distributions for the overall system (tower plus detector) and for the detector alone. The ultimate criterion used t o evaluate the validity of the experimental data was a very high degree of reproducibility. This was achieved by comparing the curves, after having reduced them to common peak ordinate. I t was possible to record a series of several such curves following a start-up period of the sampling equipment, when an equilibrium was achieved and before erratic behavior was observed because of excessive loading of the cold traps. This procedure was also a check that the spray-drying tower had been operated a t a fixed steady state during the whole experiment. This precaution was necessary because the operation of the tower was not in the control of the experimenter and operating parameters could vary without the experimenter being aware of it. Since the experiments were conducted during actual plant production periods, the process parameters could not be varied at will. It was possible to establish. however, that the experiment was reproducible when operating conditions were maintained constant. The non-normalized, experimental R T D F ' s on which the published model parameters are based are given in Figure 4. T o compute accurately a mean residence time. the tails of the curves were approximated by an exponential decay. A semilog plot of the curves showed that they contained at least two vanishing exponentials and Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
159
Figure 5. Model of the air flow pattern in the spray-drying tower
tional analog elements. Analog elements were used to reproduce the exponential decay while the VDFG was used to simulate the rest of the curve by a 10-straightsegments approximation. The simulated and experimental curves are shown in Figure 6. Since a 10-sec lag was observed between the two R T D F curves, a corresponding plug flow was incorporated into the model, but it was divided between an initial plug flow involving the whole stream of air and a by-pass. The analog simulation of plug flow elements is discussed in detail in the Appendix. This modeling of the bulk of the flow by two stirred tanks in parallel with a plug flow by-pass is thought to reflect the existence of a rapidly ascending central stream surrounded by an annular zone of intense turbulence (Chaloud et al., 1957). The third perfectly stirred tank added to produce a better match of the tail of the overall response seems to represent a region of relative stagnation. The system of equations describing the flow of tracer through the system is:
-.-.experimental
- model
10 20 30
50
75
160
125
Time (sec 1 Figure 6. Experimental and simulated RTDF curves
gave a range for the region which could successfully be approximated by an exponential. The exponential coefficient and the intercept were determined by a linear curve-fitting technique performed on a digital computer. For each curve, the curve fitting was repeated with a successively increased number of data points in the tail. This method was very satisfactory for determining the segment of the curve over which a constant and accurate fit could be obtained. With too few points, noise would overshadow the time decay, and erratic data were observed when increasing the number of points. When too many points were used, a nonrandom deviation between experimental data and computed data was observed. For the two curves used in the simulation, the exponential approximations are: detector
yl = 1.175 x lO-‘exp - 0.266 ( t - 75.2)
tower + detectory? = 1.175 x lO-‘exp
= 5.5 sec
72
= 22.5 sec
r 3 = 4.5 sec
74
= 19.9 sec
T~
- 0.069 (t - 112.4)
Mean residence times can then be computed; for the tower it is: t = 22.3 & 1.0 sec. The error on 7 was estimated from the values obtained with the various pairs of curves taken t o check the reproducibility. This value corresponds to a mean flow rate of 1075 f t ’ per sec. Analog Simulation
The model shown in Figure 5 was developed on a TR-48 EA1 analog computer by the process described below. The input signal was generated in part by a variable diode function generator (VDFG) and in part by conven160
The subscripts of the various functions correspond to the numbered points on the model (Figure 5 ) , and fe is the general response to some input, f,,. The residence-time distribution function for the model corresponds to the solution of Equations 1 to 7 for f e ( t ) given a set of values for the parameters T ) , T ~ ,T ? , T ~ , a , a. The solution was obtained on the analog computer and f e ( t ) was compared to the experimental R T D F , f ~ ~ ( t ) , which was predrawn on the chart of an X-Y plotter where the output was recorded. The parameter values were adjusted until the best fit by eye of the model R T D F to the experimental R T D F was obtained. Two of the six parameters can be constrained a priori by the following conditions imposed by the physical system: the sum of the two plug flows must equal 10 sec; the average residence time of the model must be equal to 23.3 sec. The four parameters, T ~ T, ~ a, , and a, were selected as independent parameters. A best fit was obtained with the following set of parameters:
Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
= 0.93
/3 = 0.07
u = 0.13
b = 0.87
01
However, the model response was particularly sensitive to changes in some of the parameters making the evaluation of the lack of fit extremely difficult and tedious. In addition, there appeared to be a number of different sets of parameter values which gave “good” fits. I t was concluded that a quantitative evaluation of the error was necessary and that digital curve-fitting methods should be employed.
Use of Statistical Methods with Digital Simulation
Use of a digital solution of the model equations permits a more quantitative measure of the lack of fit between the model R T D F and the digitized, experimental RTDF. Three different error criteria were used to evaluate the lack of fit between the curves. Let the deviation between the curves a t time t, be dt
= fm(tt)
- fe(ti)
then the error criteria selected were the “total sum of squares”
” ti
df
= i = 1
the “total relative error”
(9) and the “total absolute error” n
€
3
~ 1 = 1
Since feasible upper and lower bounds could be imposed on the four parameters, a central, composite, statistical design was constructed by which the allowable region of interest in parameter space could be investigated in detail. Initially, 25 sets of parameter values which covered the allowable range were chosen using a central, composite design. The tower model Equations 1 to 7 were solved numerically on an IBM 1800 using the simulation language CSMP (continuous systems modeling program) for each of these 25 sets of parameters and the values of c l , e L , and were computed. These results were used to fit the three “responses”-cl, € 2 , and cj-using full, secondorder, polynomial models in the four independent parameters, CY, a, T ~ ,7 3 . The least squares method of estimation was employed. From these prediction equations, contours for each of the error criteria were print-plotted using the on-line printer of the 1800 computer system. An example of these contours is shown in Figure 7. Visual examination of the contour plots indicated that there was a region in parameter space that would minimize the error for .O
-
9 8
+
I
6
:
I
I
t”
I
;
I I
I
designotion
C5 C14 c25 Analog solution
Error criteria
x 10’
x lo2
cy
a
TI
TI
0.168 0.168 0.251
0.631 0.631 0.398
7.94 10.0 7.94
10.0 14.0 18.0
0.994 0.290 0.455
1.695 1.446 1.724
1.3iO 0.963 1.280
0.93
0.13
4.5
17.9
0.543
4.128
1.459
f,
f?
f t
each of the three criteria; the parameter set that produced the best visual fit on the analog computer was in this region. Thirty sets of parameters were chosen from this region and values of CY, a, T % , and r4 were computed; three sets of parameters gave a better fit for all the error criteria than the analog solution (Table I ) . For the best parameter values (C14 in Table I ) r1 is 10 (therefore, T~ = O ) , indicating that the initial plug flow in the model was not needed. The results also indicated that additional searches around the set C14 were required to find the true optimum values. Although this technique was better than a visual method using the TR-48 analog computer, the whole procedure was quite time consuming. I t was decided that an automatic-directed search method should be applied t o find the true optimum values. Method of Howland-Vaillencourt with Digital Simulation
The parameters of the flow model were then determined on the digital computer by use of the generalized Newton’s method described by Howland and Vaillencourt (1961), also referred to as the method of differential corrections. In the sense of a problem in parameter estimation, f e ( t ) may be regarded as a function of the parameter values. Let any set of parameter values be denoted by the vector 6 and let the state variables of the model be represented by f so that
f= fe
I
L - - l
The state Equations 1 to 7 are then in vector form
i. = @ ( t ,f, 6)
I
OL = 0 . 2 5
Independent model parameters
I
\
I
I
I
Point
I
POlNTNO.CZ5
I I I
I I
7
p
/
I
Table I. Best Solutions for Points Selected from Contour Plots
I II I
r4=18.0sec I I
3
(11)
where @ is the tensor of arithmetic operators of the right-hand sides of Equations 1 to 7. Any solution of the state Equation 11, numerical or analytical, may be represented by:
If the model is physically correct, then corresponding to the observed tracer distribution, f r D ( t ) ,there exists some set of parameters 6* such that f < ( t ) is close in some sense t o f ~ , ) ( t ) .The truncated Taylor series expansion of fc, in the neighborhood of 6* is
2
4
1
.1
.2
.3
.4
.5 .6
.7 .8 .91.0 i = 1
a-
f ~ ~1 , (J t = 1,.. . . N, are available, then Equation 13 relates the error in each point,
If N experimental points, Figure
7. Examples of
contours of error surface
Ind. Eng.
Chem. Process Des. Develop., Vol.
10, No. 2, 1971
161
t,, of the curve fit to a correction in the parameter set, given by
G * D= E
(14)
where
G=
The elements of the matrix, G, are determined by differentiating the state Equation 11 implicitly to form the “equations of variation.” Defining the compound vectors
VI VZ
U1
uz U(t) =
V(t) = U7
V7
where
ui
vi =
=
the “equations of variation” may then be written in general form as
u=V
(16)
Equation 16 represents a set of first-order differential equations which may be solved simultaneously with the state e quat ions
f=d
6”-’
The iteration was continued until either D was arbitrarily small or until the error criterion chosen was satisfied. The relaxation factor, A, was introduced to reduce excessive oscillations in the parameter values. If the initial guess, 6 ’ , was not “close” to the optimal set, 6‘, the vector D would produce wild oscillations in 6 on the first few iterations, thereby slowing down the convergence. Therefore, a small value of A ( ~ 0 . 1 0 was ) used in the first few iterations and the value increased as the solution converged. The error criterion used in this method was the total relative error defined in Equation 9. T o be sure that the search led to a global optimum rather than a local optimum, several starting values covering a wide range in parameter space were used. As shown in Table 11, this method yielded a unique global optimum very close to one of the three solutions obtained in the previous search, but different from the original TR-48 solution. Furthermore, when any of these three solutions was used as initial value, the global optimum was still generated. Interpretation of the Flow Model
The optimum values for the parameters imply several important points with regard to the physical interpretation of the model. The initial plug flow section was insignificant. The first stirred tank was initially conceived as a relatively stagnant air pocket in the cone of the tower with a fairly large residence time. The parameters indicated that 83% of the air flow entered this well-mixed area with a mean residence time of only 13 sec. Since the R T D F has a linear property, the two major portions of the model may be considered in any order. Therefore, the first stirred tank is believed to represent a highly turbulent region a t the top of the tower where the spray is atomized. Even though this zone is close to the outlet of the tower, only slight channeling occurs (17%). The pure plug flow section, through which 60% of the air flows immediately on entrance to the tower, may represent the forced vortex swirl found by Chaloud et al. (1957) in the tower without spray. The two tanks in-series may represent a large-scale, free vortex flow in the region around the walls of the tower. These hypothesized flow patterns are pictured in Figure 8. The velocity profiles in the unloaded tower-Le., no spray-from the previous study by Chaloud et al. (1957) suggest several significant effects of the spray on the air flow pattern: the turbulent region around the spray nozzles is absent in the unloaded tower; the forced vortex air core is much smaller in the unloaded tower; and, there
(17)
to produce the response function f e ( t ) and the gradients
Table II. Reduction of Parameters to the Final Values by N e w t o n Iterations
a f e ( t ) / a 6 , ,i = 1, . . 4.
An iterative procedure was established for searching through the parameter space based on Equation 13. The procedure was as follows: An initial guess was made for the parameter set, 6 ; the state Equation 16 and the variational Equation 17 were solved for the f e ( t i ) and af,(t,)/aij,, i = 1,. , . 4 , j = 1,. . . N ; the matrix, G, and the vectors, D and E, were computed; the direction and magnitude by which each parameter value was to be adjusted was computed by solving Equation 13. 162
D = (GTE).(GTG)-’ = 6 “ +XD,O < X < 1
Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
Point
urged for initial value
c5 C14 C25 Analog solution
No. of iterations
a
U
T?
74
29.70 28.70 28.70
1.212 1.189 1.189
29.70
1.212
21 3 53
0.169 0.167 0.167
0.610 0.607 0.607
9.37 10.00 10.00
71
0.169
0.610
9.37
is a significant amount of “backflow” in the central portion of the loaded tower that was not observed in the unloaded tower.
a. Generation of fl
FunctionGeneration
*fl
Conclusions
A model of the air flow in the spray-drying tower has been developed. This model corresponds to actual production conditions, in presence of the spray. I t can be physically interpreted in terms of flow patterns which are compatible with velocity profiles obtained independently. This model will be used for further studies of the process and of the plant. How flexible and reliable a tool it will be in these studies will, of course, be the ultimate test of its validity. Analog simulation on the TR-48 was adequate to construct an arrangement of ideal flow components providing a reasonable interpretation of the actual air flow as well as a satisfactory fit to the experimental data. The limitations of the equipment, however, did not permit a definitive estimation of the parameters. The statistical analysis with digital simulation was an excellent method for locating the regions of a possible global minimum in the error, but further searching for this minimum required an excessive number of solutions of the differential equations. By far the fastest and most acceptable method was Howland and Vaillencourt’s modification of the generalized Newton’s search method with digital simulation. This method converged very rapidly to a local minimum for an appropriate initial guess. The only difficulty with the method arose when estimates were very different from the optimal values, in which case the method tended to oscillate wildly. This problem was overcome by introducing the relaxation parameter, ,I, and, by an appropriate definition of a sequence { A, },the method converged to a local optimum nearly monotonically. Appendix. Simulation of the Plug Flows in the Model
To generate f l ( t ) and f 4 ( t ) in Equations 2 and 4, a plug flow (transport delay) must be simulated on the analog or digital computer. The function f l ( t ) is identical t o f ! , ( t ) delayed by a time T ~ . The only form of f,,(t) which is of practical
b. Generation of f4
(t-r1-r3)
Function Generation for fin(t1
Figure 9. Simulation of the plug flows in the model
interest is the R T D F of the detector f ~ ( t )obtained experimentally as a continuous curve. Therefore, on the digital computer, this curve can be used to generate f,, ( t ) in tabular form using any of the well-known digitization techniques. For the analog computer, a number of strategic points can be selected from the curve and used to set up a variable diode function generator, resulting in the continuous generation of f L n ( t )u p to a point beyond which the curve is approximated by an exponential decay generated by conventional analog elements. In block notation, the generation of f l ( t ) is shown in Figure 9, case a. The same technique can be extended further to generate f 4 ( t ) using the same function generator. Referring to Equations 3 and 4, to generate f 4 ( t ) we must generate f 3 ( t - 74, f d t - 74, f l ( t - 74, and ft,(t - T~ - 74. The resulting diagram is shown in Figure 9 as case b. Acknowledgment
The authors are indebted to Edward F. Leonard, Columbia University, for his consultation throughout the course of this work, and to Richard K. Sakulich who programmed the calculations for the Howland-Vaillencourt method. They are also indebted to many of their colleagues a t Procter & Gamble Co. whose expertise on the spraydrying process and generous advice made this work possible. Nomenclature
Spray Direction
1 enum hamber I
-0-
i
Air Flowdirection in the Plenum
[Grzies Figure 8. Sketch of the air flow pattern in the tower
f T = tower residence-time distribution function f D = detector residence-time distribution function f r n = combined tower-detector system residencetime distribution function f t = intermediate response of the flow model to any input function, f,,, a t points correspondingly labeled (Figure 5) f? = response of the model to input f,, ii = average residence time for each element o the flow model, sec a, b , CY, d = flow splits in the model t , , t 2 , c 3 = error criteria = relaxation parameter Literature Cited
Adler, R. J., Long, W. M., Rooze, J., World Petrol. Congr Proc., 6th (1963). Chaloud, J. H., Baker, J. S., Chem. Eng. Progr., 53, 5936 (1957). Danckwerts, P. V., Chem. Eng. Sci., 2, 1 (1953). Himmeblau, D. M., Bischoff, K. B., “Process Analysis and Simulation,” Wiley, New York, N . Y . , 1968, pp 61-6. Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
163
Howland, J. L., Vaillencourt, R., J . SOC.Znd. Appl. Math., 9 ( 2 ) , 165-8 (1961). Levenspiel, O., “Chemical Reaction Engineering,” Wiley, New York, N . Y., pp 244-8. Levenspiel, O., Bischoff, K. B., Aduan. Chem. Eng., 4, 95-198, (1963). Levenspiel, O., Smith, W. K., Chem. Eng. Sci., 6, 227 (1957). Moser, J. H., Cupit, C. R., Chem. Eng. Progr., 62 (6), 60-5 (1966).
Place, G., Ridgway, K., Danckwerts, P. V., Trans. Inst. Chem. Eng., 37, 268-76 (1959). Van Deemter, J. J., Chem. Eng. Sci., 13 (3), 143-54 (1961). Van der Laan, E. Th., ibid., 7, 187 (1958). Wen, C. Y., Chung, S. F., Can. J . Chem. Eng., 43 (3) 101 (1966). RECEIVED for review August 7, 1969 ACCEPTED December 7, 1970
Analysis of Nonisothermal Moving Bed for Noncatalytic Solid-Gas Reactions Masaru lshida and C. Y. Wen’ Department o f Chemical Engineering, West Virginia Uniuersity, Morgantown, W . Va. 26506 The design and performance of a moving bed for solid-gas reactions under nonisothermal conditions are discussed based on unreacted-core shrinking model. Numerical solutions based on unsteady-state analysis are compared with graphical solutions based on pseudosteady-state analysis. The pseudosteady-state analysis i s satisfactory if no abrupt temperature change occurs in the bed. For an exothermic reaction, thermal instability may occur. In such cases, the initial temperature of solid will considerably affect the reactor performance. Moreover, the transition of rate controlling steps occurs more readily in a moving bed reactor than in a single particle.
M o v i n g beds, rotary kilns, and fluidized beds have been widely used in industries for solid-gas reactions. In spite of this, the design and performance of these reactors are not thoroughly based on fundamental principles of solidgas reactions. A nonisothermal analysis of both irreversible and reversible reactions between gas and a single particle was presented both from theoretical and experimental points of view in previous papers (Ishida and Shirai, 1969, 1970; Ishida and Wen, 1968; Wen, 1968; Wen and Wang, 1970). The important role played by heat of reaction in solidgas systems, was discussed. If the basic information obtained from a single particle can be applied directly to the design of the multiparticle systems, the optimal conditions of an industrial-size reactor can be selected. Such attempts have been made by several investigators. Meissner and Schora (1960) analyzed a multistage reactor for reduction of iron ore by comparing it with a distillation column. Yagi et al. (1968) analyzed a moving bed reactor in which isothermal solid-gas reaction took place. Horio et al. (1969) extended this study to nonisothermal cases assuming that the temperature of the solid is uniform and the reaction proceeds under the chemical reaction controlling regime. Jackman and Aris (1968) discussed the optimal control of pyrolysis reactors in which a coke layer builds u p decreasing the heat transfer coefficient a t the wall. Since these analyses are made either under
’ To whom correspondence should be addressed. 164
Ind. Eng. Chem. Process Des. Develop., Vol. 10, No. 2, 1971
isothermal condition or in the absence of diffusional effect within the solids, they can only be applied to restricted cases. This paper will present nonisothermal steady-state analysis of a moving bed utilizing the knowledge gained in the single particle studies. I n these studies, the heat of reaction and simultaneous mass and heat transfer within the particle are taken into consideration. Thermal and transitional instabilities of exothermic reactions are also discussed. The particle follows the unreacted-core shrinking model, while solids move in the bed according to plug flow. Equation Derivation
Consider the following reversible solid-gas reaction
S,
+ aA
p1
Z SI,
+ bB,
(1)
In an unreacted-core shrinking particle, the reactibn rate per unit reaction surface area, m 4 , is assumed to be uniquely determined by a set of intensive properties a t the reaction surface as follows: mA
= mA(TC, C A c , C B ~C, I ~ )
(2)
Furthermore, the particle size is assumed to be constant during the reaction. Figure 1 shows a schematic diagram of a moving bed. The gas at temperature ( T o ) ois fed from z = 0 a t a molar flow rate, G, and leaves the reactor a t t = 2 a t a temperature, ( T o ) l .I n cocurrent flow, the solids, S,
