Detailed computation of unsteady droplet dynamics

Twentieth Symposium (International) on Combustion/The Combustion Institute, 1984/pp. 1743-1749

D E T A I L E D C O M P U T A T I O N OF U N S T E A D Y D R O P L E T D Y N A M I C S H. A. DWYER AND B. R. SANDERS Sandia National Laboratories Liverrrugre, California 94550

In this study, a detailed time-dependent calculation of the lifetime of a fuel droplet is presented. The calculation includes fluid flow, heat transfer, and mass transfer in both the gas and liquid phases as the droplet radius changes with time. For a starting gas phase Reynolds number of 100, the transport processes are initially dominated by convection and then become conduction and diffusion controlled during the latter part of the droplet lifetime. A comparison between the calculated local drag and heat transfer and the estimates which would be obtained using available correlations shows that there is a very large deviation and, that the correlations will lead to incorrect results.

I.

Introduction

The details of the heat transfer, vaporization, and fluid dynamics associated with fuel droplet dynamics are complex and can cover a wide range of regimes. In this paper, a detailed calculation of the history of a fuel droplet at an initial Reynolds number of 100 is carded out. The calculation is quite general and fully capable of being applied to a wide range of droplet compositions as well as conditions in the gas phase, such as pressure. The calculation is an extension of previous constant Reynolds number results, 1 and its purpose is to provide detailed information on the transfer processes which occur in practical situations. This detailed information will eventually be used to evaluate the success of simpler models, such as those discussed by Sirignano, 2 Law, 3 and Chigier. a The method employed here is the numerical solution of the two-dimensional, axisymmetric equations for momentum, energy, and mass transfer. The numerical approach is very general and allows the computational mesh to move with the changing droplet diameter. It also allows for a self-consistent drag calculation. These procedures use generalized time-dependent coordinates, a linearized coupled form of the Clausius-Clapeyron equation to allow simultaneous interaction of the surface temperature and mass fraction, and two-phase surface control volumes to calculate the conditions at the gas-liquid interface. The following sections of the paper will describe the history of a liquid fuel droplet with an initial Reynolds number of 100 in the gas flow. The properties of the droplet are similar to those of a heavy hydrocarbon, and the background pressure has been chosen as 5 atmospheres. The time development of

both the liquid and gas phases has been determined and the drag and mass transfer rates used to calculate the new droplet velocity, radius, and Reynolds number. It will be shown that the processes are very complex and that both the drag and heating rates vary greatly from typical correlation values. Also, the heat transfer varies from one dominated by convection to one controlled by conduction. Because of the large variation in fuel droplet characteristics, it is diflqcult to make generalizations to other droplet conditions without further calculations. However, it seems clear that most high Reynolds number droplet flows will involve complex flow processes, and that detailed time-dependent calculations will be necessary. A recent paper by Renksizbulut and Yuen 5 also presents numerical calculations for the flow around a single droplet vaporizing in a convective stream. However, direct comparisons cannot be made with their results as they solved for completely steady state solution, and the mechanisms and amount of mass transfer were much different. H. Methods of Analysis 1. Basic Equations

The starting point for this study is the transport equations in terms of cylindrical coordinates utilizing stream function and vorticity formulations. The equations in a compact vector form are:

aQ

--+---t at ar where

1743

--

az

+ --

ar

az

+ I~I

SPRAY COMBUSTION

1744

0 foro IprCpT]

0 f prvto ~ =

\ prvYf /

\ prYy/

=

p ,STn i s = Ors

1)

orvC T!

0 f pruo~ ~ IpruCpTI

r Or

0 r O-rrIxto fi =

rk c3T/c3r

\ pruYf /

prDfOY_--~ f dr~

,

c3~ r c~z 1

--OJ

Ixo~ + pvo~

--Ixr

~z

fi =

-koT/Oz

r

0

prDf c3YJ

where pvs is 'the surface mass flux of fuel per unit area of surface and n is the coordinate normal to the interface. For the vorticity, energy and fuel density equation a fully implicit numerical method was employed, while successive over-relaxation was used to calculate the stream function. Due to the lack of a time scale in the stream function equation and a well known sensitivity to numerical convergence, an overall iteration procedure was applied to the entire system. This iteration procedure is primarily employed in order that the stream function and vorticity be coupled and updated, and is not particularly needed for the energy and species equations. The exponential nature of the ClasiusClapeyron equation and the associated numerical difficulties were treated by the use of Newton's method. This linearization method along with the use of a block tridiagonal solution of the energy and species equations resulted in increased computational efficiency.

c3z

2. Physical Parameters of the Problem The following notation has been employed: ~ - stream function, oJ--vorticity, T--temperature, Y~. fuel mass fraction in the gas, Ix--viscosity, ,~ thermal conductivity, p--local gas density, and Cp--specific heat. Also, the dependent variables t, r, and z are the time, radial coordinate, and axial position, respectively. The formulation of the equations can be made much more useful by transformation to generalized non-orthogonal coordinates, which allow for arbitrarily shaped bodies as well as the change in droplet radius with time. The details of the transformed equations can be found in Dwyer et al. (1983).6 At the interface between the liquid and gas phases of the droplet, a special treatment of the dependent variables must be employed. The stream function and surface fuel density are determined by the mass transfer and Clausius-Clapeyron relationship, respectively. For the vorticity and energy equations a special spherical element is formed, which is part liquid and part gas. It is essential that all terms in the transport equations be retained, since at various times in the unsteady analysis they can play a key role. Also, in the present work, the mass diffusive transport equation has only been solved in the gas, since the liquid phase is homogeneous. In addition, the energy equation at the interface must be modified to include the influence of latent heat of vaporization. The amount of vaporization is determined from fuel mass conservation with diffusion, which yields:

We have chosen a set of parameters which are of practical interest. The present study will be extended in future applications to the parameter space of interest to overall droplet dynamics. The conditions chosen are as follows: GAS PRESSURE GAS TEMPERATURE INITIAL LIQUID TEMPERATURE AVERAGE VISCOSITY RATIO AVERAGE DENSITY RATIO LATENT HEAT OF VAPORIZATION OF LIQUID BOILING POINT OF LIQUID

P = 5 atm T = 1000 K TI = 400 K Ixe I~g

--=25

Pe = 300 pg L - - = 106 K Cp TB = 573 K --

= 6.675 OF LIQUID PRANDTL NUMBER OF GAS PRANDTL NUMBER OF LIQUID INITIAL GAS REYNOLDS NUMBER

4542 T Prg = 1.0 Pre = 10 UD Reg -

-

Pg

100

COMPUTATION OF UNSTEADY DROPLET DYNAMICS HI. Results and Discussion

1. Local Properties The major purpose of the present paper is to investigate the history of a vaporizing fuel droplet. To accomplish this goal, a time-dependent calculation with both a time varying radius and Reynolds number has been carried out. The calculation consists of a spherical droplet initially without internal circulation injected into a relatively hot gas stream and given an initial velocity so that the starting Reynolds number is 100. Both the gas and liquid properties have been held constant in order to better comprehend the basic processes; however, future investigations will involve variable properties. In Fig. 1 the initial coordinate system is shown for the inner part of the calculation, and there is noticeable clustering of the grid near the particle surface in both the liquid and gas phases. The flow proceeds from left to right, and the instantaneous radius Ro is given in the figure captions. The results are given in terms of a nondimensional gas time scale defined by ~g = tvg/R z, where t is time, R is initial droplet radius and vg is the gas kinematic viscosity. As the droplet radius changes in time, the entire grid rescales and the time-dependent motion of the coordinate system is included in the calculation. The isotherm distribution shown in Fig. 2 is for a time early in the droplet mass lifetime, but there has been a considerable change in the droplet Reynolds number such that it is almost one half of its initial value. The temperature distribution is typical of the high Reynolds number patterns foUnd in Ref. 1. (Notice that the isotherm spacing, which is indicated in the caption, is different away from and near the particle. Three gas isotherms with a spacing between isotherms of ATg = 100 K are shown, while the liquid and near surface gas isotherm spacing is ATe = 10 K. The influence of the decay of the Reynolds number on the internal circulation can be clearly seen by observing the isotherm patterns shown in Figs.

1745

3 through 5. At the later stages in the droplet lifetime, the conduction process in the droplet is on an equal basis with convection, and the droplet flow cannot be classified as a high Peclet number one. The history of the 530 K isotherm can be seen clearly in Figs. 3 and 4. In Fig. 3, the isotherm is just beginning to be convected away from the back of the droplet. The results of the calculations show that the convective process inside the droplet is in the transition region from a flow strongly influenced by convection to one weakly influenced by convection. Although the droplet isotherm patterns are strongly influenced by the Reynolds number decay, the droplet streamlines remain almost identical as shown in Figs. 6 and 7. This is due to the spherical shape of the droplet and also to the fact that the droplet internal recirculation region due to gas phase wake separation is very small for the parameters under consideration. 1 The droplet streamlines do not give a good indication of the decay of the internal circulation since the maximum surface velocity has decayed from 2.42 percent to .66 percent of the free-stream velocity. A good indication of the influence of surface mass transfer on the gas flow can be gained by observing the t~ = 0 streamline in Figs. 6 and 7. For a flow with zero surface mass transfer, this streamline would lie on the liquid surface, but as can be seen from the figures, this is far from being the case. The later time results show that the mass transfer has caused a major perturbation of the gas flow, and that a simple correlation for surface drag and heat transfer may not be valid. It is becoming clear that droplet flows at these intermediate Reynolds numbers are quite complex and may not lend themselves to simple modeling concepts.

2. Global Properties The results presented in Figure 2 through 7 were spatial distributions of flow properties at a given

c3

c4-

)0.

2q

2.0

O.

-2.0

-tO

0.0 X

LO

2.0

FIG. 1. Inner coordinate system utilized for droplet calculation ra = 0.

kS

tO

-'0.5

O.O X

0.5

tO

t.S

FIG. 2. Early time droplet isotherm distribution (Re = 54.8, Mass = 91.8 percent, rg = 11.92, Ro = .972, T,,,, = 470 K, ATe = 10 K, and ATg = 100 K).

1746

SPRAY COMBUSTION

~.o.

-2&l

tO

0.0

O~

I.~

L*

1.0

X

to

to

~

0.0

0~

tO

to

~.0

X

FiG. 3. Isotherm distribution (Re = 33.76, Mass = 66.2 percent, r e = 27.9, Ro = .872, Train = 530 K, ATe = 10 K, and ATe = 100 K).

FIG. 6. Streamline contours (Re = 54.8, Mass = 91.8 percent, r e = 11.92, Ro = .972).

9_

-2.0

-tO

-s

~

0.0

0~

tO

1,5

~11

x

FIG. 4. Isotherm distribution (Re = 28.5, Mass = 55.5 percent, "re = 33.4, Ro = .823, Tmi n = 5 3 0 K , ATe = 10 K, and ATe = 100 K). 9.

_~. - t . t O - t . ~ -too ~

~

-0.15

~ X

0.25 0.m

~

tOO

~

tm

FIG. 5. Isotherm distribution (Re = 21.7, Mass = 39.7 percent, r e = 41.4, Ro = .735, Tmi n = 540 K, ATe = l0 K, and ATg = 100 K).

point in the droplet lifetime. In this section global parameters of drag, mass transfer and heat transfer are given. Detailed calculations were performed which gave the spatial development in time of the pressure distribution and liquid surface friction, the distribution of mass transfer over the droplet surface, and the local heat balance between liquid and gas phases. These detailed calculations will be presented as integral quantities over the droplet sur-

X

FIG. 7. Streamline contours (Re = 21.7, Mass = 39.7 percent, r e = 41.4, Ro = .735). face, as there is not space to describe these detailed processes in this paper. Two global results are presented in Fig. 8. These are the percent of initial mass, and the instantaneous Reynolds number of the droplet based on the velocity, diameter and constant free stream gas conditions. The global mass history exhibits two distinct regimes, a relatively slow vaporization rate during the heat-up phase of the droplet followed by a nearly constant mass transfer rate during the remainder of the calculation. The Reynolds number shows a dramatic reduction in the early portion of the droplet lifetime, changing by a factor of five while the radius has been reduced by only twenty percent. The dependence of drag on the square of the instantaneous velocity is responsible for most of the early Reynolds number reduction and dictates the need of an unsteady analysis throughout most of the lifetime. The drag coefficient variation with vaporization, along with the local steady-state drag coefficient without vaporization7 is sliown in Fig. 9 . It is immediately clear that the drag coefficient with vaporization is much lower than without vaporization. The actual drag force on the sphere was calculated by applying a control volume analysis in an accelerating coordinate system with moving boundaries.

COMPUTATION OF UNSTEADY DROPLET DYNAMICS

80r}

u')c~_ I.-, Z

S I

c;o zos _..1 0 Z""

~-

0.0

20.0

I0.0

' 30.0

'' 40.0

T-GflS TIME SCALE

FIG. 8. Reynolds number and percent of initial droplet mass as a function of time. The most important and difficult parts of the force balance are the surface pressure and shear force terms. To evaluate these terms, both the second and third derivatives of the stream function were calculated. As mentioned in the introduction section, a re-

1747

cent study by Renksizbulut and Yuen 5 presented results for the steady state vaporization of a droplet without internal circulation and with much smaller values of mass transfer. Our results for Nusselt number and drag coefficient are qualitatively consistent with their results, but our drag coefficients are considerably less in magnitude. Considering the much different natures of the physical problems a detailed explanation of the differences cannot be made at this time. However, it should be mentioned that our results are consistent with the low Reynolds number calculations of Sadhal and Ayyaswamy, s which are relevant to the lower Reynolds numbers towards the end of our calculations. T h e heat transfer coefficient or Nusselt number is not as strongly influenced by vaporization as the drag coefficient, and this is probably due to the more nonlinear nature of the momentum equation compared to the energy equation. The variation between the calculated Nusselt number and that computed for a solid sphere 7 is shown in Fig. 10. A major portion of the variation is caused by using a constant temperature difference in the definition of Nusselt number and not employing a mean temperature for the droplet. It is not obvious what mean temperature to use, since the surface temperature is usually much higher than the interior of the droplet. At the final stages of the calculation, r 40, approximately 80 percent of the surface heat transfer is being used for vaporization, and the droplet is essentially isothermal, as opposed to 20 percent of the heat transfer at r ~ 10. It will be

O

o

(23

WITHOUT VAPORI ZATI ON

Z L]cD

WITHOUT 7- CO Z

L L~ LJ (D (DCO r..D ~3Z

Js WITH

:-

VAPORIZATIO ~N

zD o. fxl

WITH VAPORIZATION

CO

CZ)

O

CO I

0.0

10.0

I

20.0

I

I

30.0

40.0

T-GRS TIME SCRLE FIG. 9. Drag coefficient variation.

I

0.0

10.0

I

20.0

I

30.0

T-GFIS TIHE SCRLE FIG. 10. Nusselt number variation.

I

40.0

1748

SPRAY COMBUSTION

an objective of further work to develop a rationale for a heat transfer correlation during droplet vaporization.

ough understanding of the fuel droplet vaporization problem.

REFERENCES

Conclusions A detailed calculation of the time-dependent lifetime Of a fuel droplet has been carried out. The primary goal of the research was to understand the physical processes and to estimate the various time scales which occur in the problem as well as to develop the numerical techniques necessary to complete the calculation. The major conclusions are: 1. Both the Reynolds and Peclet numbers decay substantially during the lifetime of a fuel droplet. This decay makes the use of simple models difficult, since the process transforms from one which is convection dominated to one controlled by conduction and diffusion. 2. The calculated values of heat transfer and drag are substantially lower than the values which would be predicted from steady state correlations. However, the mass transfer rates are dominated by the rapid increase in surface concentration with temperature. 3. The present study has been carried out for fuel droplet parameters which are typical but does not cover the complete range of values. It is very probable that the present results will have to be extended to these other regions to provide thor-

1. DWYER, H. A. AND SANDERS, B. R.: "'Comparative Study of Droplet Heating and Vaporization at High Reynolds and Peclet Numbers," International Colloquium on Dynamics of Explosions and Reactive Systems, Poitiers, France, July 1983. To be published by the AIAA. 2. SIRIGNAr~O,W. A.: Progress Energy Comb. Sci., 9, 291-322 (1983). 3. LAw, C. K.: Progress in Energy and Combustion Science, 8, 169-199 (1982). 4. CmGIER, N.: Combustion and Flame, 51(2), 127139 (1983). 5. RENKSIZBULUT, M. AND YUEN, M. C.: J. of Heat

Transfer, 105, 389-397 (1983). 6. DWYER, H. A., KEE, R. J., BARR, P. K., AND SANDERS, B. R.: J. Fluids Engineering, 105, 8388 (1983). 7. BIRD, R. B., STEWART, W. E., aND LIGHTFOOT,

L. N.: Transport Phenomena, John Wiley and Sons, Inc., 1960. 8. SaDUaL, S. S., aND ArraswaMr, P. S.: J. Fluid Mech., 133, pp. 65-81 (1983). 9. LECLMR, B. P., HAMIELEE,A. E., PRUPPACHER, H. R., AND HALL, W. D.: J. of Atmospheric Sciences, 29, 728-740 (1974).

COMMENTS F. Boysan, Univ. of Sheffield, England. What kind of boundary conditions did you use at the outer gas boundary? How does the distance between the outer boundary and the droplet surface affect your calculations? Authors" Reply. The boundary condition at the outer gas boundary was the condition of uniform flow.

The outer boundary for numerical calculations should be at least 7 radii away from the droplet surface. In the present calculations 15 radii was utilized.

K. Kuo, Penn State Univ., USA. Very interesting paper. I have two questions. First, your results show that as the Reynolds number decreases, there is a conversion from an initially convective internal flow situation to a more symmetric isotherm pattern

controlled by conductive mechanism. This seems to be consistent with the one-dimensional approach just presented by Dr. Prakash at this symposium in his consideration of droplet heating process. I wonder whether you and Dr. Prakash both considered the same mechanism to reach your results or is this just a coincidence? My second question is related to the droplet shape. If the droplet shape is allowed to be nonspherical, the internal recirculation may be maintained at a high level for a longer period of time. Then, the isotherms will not quickly return to a symmetric condition. Would you please comment on the effect of droplet shape on the internal recirculation?

Authors" Reply. If the droplet lifetime is short then convection will be important. In the present calculation it was important. Droplet shape variations were not investigated in the present research.

COMPUTATION OF UNSTEADY DROPLET DYNAMICS

D. Harrington, General Motors, USA. Are you accounting for the variation in vapor pressure with droplet diameter due to surface tension? This would appear to be an important factor for the 10 to 200 p,m drop sizes associated with many fuel injectors. Although it was stated that unsteady drag coefficients were utilized, are they in fact not just steadystate drag coefficients modified by mass transfer? For example, if the droplet is accelerating (or decelerating) to the free steam velocity, it will, in general, have a different drag coefficient than when

1749

it is not accelerating, even for a rigid sphere. For a deformable droplet, the acceleration effect can be even more pronounced. The separation point for an accelerating sphere at velocity V is different from that of a non-accelerating sphere at velocity V.

Authors" Reply. The variation in vapor pressure was not accounted for in the calculation. The drag coeffcient was calculated from the flow, and full accounting was made for deceleration of the coordinate system. The deceleration term was rather small in the present work.