Analysis of Thermal Effects in a Cavitating Inducer Using Rayleigh Equation

A simple model based on the resolution of Rayleigh equation is used to analyze thermal effects in cavitation. Two different assumptions are considered for the modeling of heat transfer toward the liquid/vapor interface. One is based upon a convective type approach using a convection heat transfer coefficient or the equivalent Nusselt number. The other one is based upon the resolution of the heat diffusion equation in the liquid surrounding the bubble. This conductive-type approach requires one to specify the eddy thermal diffusivity or the equivalent Peclet number. Both models are applied to a cavitating inducer. The basic pressure distribution on the blades is determined from a potential flow computation in a two-dimensional cascade of flat plates. The sheet cavity, which develops from the leading edge, is approximated by the envelope of a hemispherical bubble traveling on the suction side of the blade. Cavity shape and temperature distribution predicted by both models are compared. The evolutions of cavity length with the cavitation number for cold water (without thermal effects) and for Refrigerant 114 at two different temperatures is compared to experimental data. Such a simple model is easy to apply and appears to be quite pertinent for the analysis of thermal effects in a cavitating inducer.

Introduction

In a number of fluids, the development of cavitation goes with significant thermal effects. The so-called thermodynamic effect has been extensively studied since the early works of Stahl et al. ͓1͔ and Stepanoff ͓2͔. The phenomenon is characterized, in particular, by a temperature drop inside two-phase regions. It is due to the latent heat of vaporization, which is taken to the liquid surrounding the cavitating zone whose temperature is then lowered.

The effect is negligible in water at room temperature but may be significant at a higher temperature or in other liquids as, e.g., cryogenic fluids used for the propulsion of space rockets. Several parameters have been introduced to estimate the importance of thermal effects in cavitating flows. Stepanoff ͓2͔ defined a characteristic temperature drop by

⌬T * = v L ᐉ c pᐉ ͑1͒
This is the temperature drop to be applied to a unit volume of liquid to supply the latent heat required for the vaporization of a unit volume of vapor. Experiments show that this parameter is a relevant order of magnitude of the temperature drop in cavitating flows. For water at room temperature, it is quite small and equal to 0.01 K, whereas it exceeds 1 K for liquid hydrogen at 22 K. This big difference demonstrates the major influence of the liquid properties on the phenomenon. The well-known B factor of Stepanoff ͓2͔ is the nondimensional parameter:

B = ⌬T ⌬T * ͑2͒
where ⌬T is the actual temperature drop in the cavitating flow.

Using the slope dp v / dT of the vapor pressure curve, a characteristic pressure difference can be built on the basis of the previous characteristic temperature difference ⌬T * :

⌬p v * = dp v dT ⌬T * ͑3͒
This is the drop in vapor pressure associated with the drop in temperature ⌬T*. Since ⌬T * is a relevant order of magnitude of the temperature drop, ⌬p v * can be considered as a relevant order of magnitude of the pressure drop inside the cavities. This is true as long as the thermodynamic equilibrium condition is achieved, which is generally the case. Nonequilibrium is expected only at the very final stage of collapse when the kinetics of phase change can become a limiting factor because of the very high velocity of the interface.

Brennen ͓3,4͔ suggests the use of the thermodynamic parameter

⌺ = ⌬p v * ᐉ ͱ ␣ ᐉ ͑4͒
to measure the importance of thermal effects in cavitating flows. The ⌺ parameter was derived from the analysis of the dynamics of a single bubble with the assumption that heat transfer to the interface is provided by conduction through the surrounding liquid.

The thickness of the thermal boundary layer then plays an important role and accounts for the presence of the square root of the liquid thermal diffusivity ␣ ᐉ . The ⌺ parameter has the units m/s 3/2 . Considering the case of a sheet cavity and still assuming heat transfer to the cavity interface by conduction through the neighboring liquid, Kato ͓5͔ introduced a similar parameter defined by

␣ = ͱ ᐉ v ⌺ ͑5͒
The importance of thermal effects increases with ␣ and ⌺. From a systematic analysis of Hord's experiments ͓6͔, Kato ͓5͔ showed a good correlation of the nondimensional pressure drop 2⌬p v / ᐉ V 2 inside the cavity with the nondimensional parameter ␣ ͱ c / V 3 , where c and V are a reference length scale and a reference velocity scale of the flow.

It is of major concern to observe that, in such models based on a well-defined interface between the liquid and either a cavitation bubble or a sheet cavity, heat transfer by conduction characterized For interface models, the volume of liquid surrounding vapor is almost infinite and heat transfer is supposed to be limited by the growth rate of the thermal boundary layer. Conversely, two-phase models assume an intimate mixing of liquid and vapor and the volume of liquid surrounding a given bubble is naturally limited by the surrounding bubbles. Hence, thermal effects are essentially regulated by the limited volume of liquid able to supply the heat for vaporization and not by heat diffusion as for interface models. Thermal effects are then basically controlled by the evolution of void fraction without any influence of the liquid thermal diffusivity.

The difference between both types of models can easily be demonstrated by comparing basic expressions of the B factor. In the case of a bubble whose radius grows from 0 to R during time t, the B factor is given by ͑see, e.g., Ref.

͓21͔͒ B Ϸ R ͱ ␣ ᐉ t

͑6͒

On the other hand, for a two-phase mixture of void fraction ␣ v , the heat conservation equation leads to ͑see Ref. ͓22͔͒:

B Ϸ ␣ v 1-␣ v ͑7͒
In the case of Eq. ͑6͒, thermal effects depend upon the liquid thermal diffusivity, which is obviously not the case for Eq. ͑7͒.In addition, Eq. ͑7͒ becomes singular when the void fraction approaches unity, which shows that the two-phase model is valid for rather dilute mixtures but fails when the cavitating region becomes pure vapor. If so, it is necessary to switch to interface-type models. Hence, it is of primary importance, in practice, to recognize the type of cavities in order to be able to identify the appropriate physical mechanism, which actually limits the supply of heat necessary for vaporization. The present paper is devoted to an analysis of thermal effects on the basis of the Rayleigh equation. It belongs to the class of interface models for which heat transfer through the thermal boundary layer is supposed to be the limiting physical phenomenon. The present analysis is very similar to that conducted by Watanabe et al. ͓7͔ except that it considers bubble cavitation instead of sheet cavitation. It is also applied to the case of a cavitating inducer of a space rocket turbopump investigated experimentally by Franc et al. ͓22͔.

The analysis comprises two successive steps. First, the pressure distribution in a two-dimensional ͑2D͒ noncavitating cascade is computed using a boundary element method. Second, the Rayleigh equation is solved using the previous pressure distribution in order to compute the evolution of a cavitation bubble on the blade suction side.

The original aspect of the present contribution essentially lies in the thermal term added to the Rayleigh equation in order to account for thermal aspects. Two different forms based either on a convective-or a conductive-type approach are considered. Both approaches are evaluated by comparison with the experimental results ͓22͔. The main objective is to evaluate whether this kind of simple model is able to account for the trends observed experimentally as the increase of the B factor with the cavity length.

In the present approach, the sheet cavity observed experimentally is compared to the computed envelope of a hemispherical bubble, which grows and collapses on the suction side of the blades. This kind of simple model obviously suffers from limitations. One is the lack of any feedback of cavitation on the pressure distribution since the bubble evolution is computed using the noncavitating pressure distribution. The 2D nature of the computation of the pressure distribution is also a limitation since inducers may be subjected to backflow, which may significantly change angles of attack. Nevertheless, the present work supplies a theoretical support for the physical interpretation and understanding of experimental results.

Pressure Distribution Computation

The inducer is represented by a 2D cascade of flat plates, as shown in Fig. 1. This 2D configuration is obtained from a cylindrical section of the inducer at radius r. Enlargements of the leading edge and trailing edge are shown in Figs. 2 and3. Calculation of the bubble evolution was started at the back end of the bevel, from point B where the pressure is minimum ͑see also Fig. 4͒. This location corresponds to the observed point of detachment of cavities. As for the trailing edge, its shape was chosen in order to allow a smooth detachment of the flow in a direction parallel to the blades.

The relative thickness of the blade is about 1% and the solidity c / H is 2. The incidence angle ␣ a and the stagger angle ␥ of the cascade are such that the blade bevel is parallel to the j direction of the grid. The inlet velocity V u i is supposed perpendicular to the grid ͑no prerotation͒. At upstream infinity, the reference relative velocity is

W u = ͱ V u 2 + 2 r 2 ͑8͒
where is the angular speed of rotation. Any radial flow normal to the cylinder surface is ignored and the flow is then supposed to be 2D in the ͑i , j͒ plane. In addition, it is assumed incompressible, irrotational, and inviscid. Under these assumptions, the pressure coefficient is given by

C p =1- W 2 -2 r 2 W u 2 -2 r 2 ͑9͒
The pressure distribution on the blade is computed using the boundary element method. A double distribution of sources and vortices ͑ * , ␥ * ͒ is distributed on the blade surface with the space periodicity H in the j direction to take into account the infinite series of blades. This distribution is chosen so that the relative velocity at upstream infinity is in the I direction. The boundary conditions consist in a zero normal velocity on the blade ͑slip condition͒ and the continuity of pressure at the trailing edge. Conservation of flowrate through the inducer requires that the component of the relative velocity W • i is the same at upstream and downstream infinity. The flow at downstream infinity is deviated with an angle given by tan =

͑sin ␥͒ Ͷ blade ␥ * dx HW u + ͑cos ␥͒ Ͷ blade ␥ * dx

͑10͒

The double distribution used here introduces more unknowns than equations to be solved. Additional conditions are imposed in order to minimize the discretization errors. This allows us to compute, with an improved accuracy, the pressure coefficient on blades of small relative thickness. The computational method has been validated by comparison with the results obtained by Watanabe et al. ͓7͔. Both computations agree fairly well, as shown in Fig. 4. The differences, which are observed particularly at the leading edge near the bevel, are due to the nonlinear nature of the present computation contrary to Ref. ͓7͔, which assumes blades of zero thickness.

Basic Equations

Nondimensional Form of the Rayleigh Equation.

Once the pressure distribution has been computed, Rayleigh equation

R d 2 R dt 2 + 3 2 ͩ dR dt ͪ 2 = p v ͑T c ͒ -p͑t͒ ᐉ ͑11͒
is solved in order to determine the evolution of the radius R of a bubble traveling on the suction side of the blades. p͑t͒ is the time dependent pressure to which the bubble is subjected as it moves along the blades. The effects of viscosity and gas content, which generally become negligible as soon as the original microbubble becomes a macroscopic cavitation bubble, are disregarded and the present work is focused on the thermodynamic effect. To account for thermal effects, it is necessary to consider in Rayleigh equation the vapor pressure corresponding to the actual temperature T c inside the bubble, which is different from the liquid temperature T ϱ far from the bubble. Using the pressure coefficient previously calculated

C p = p -p u ͑1/2͒ ᐉ W u 2 ͑12͒
and the usual cavitation number defined on the basis of the vapor pressure at the liquid temperature at infinity

p v ͑T ϱ ͒ v = p u -p v ͑T ϱ ͒ ͑1/2͒ ᐉ W u 2 ͑13͒ Rayleigh equation ͑11͒ becomes R d 2 R dt 2 + 3 2 ͩ dR dt ͪ 2 = p v ͑T c ͒ -p v ͑T ϱ ͒ ᐉ -W u 2 C p + v 2 ͑14͒
The first term on the right hand side accounts for thermal effects. d dt

= W u ͱ 1-C p d dx ͑15͒
In addition, lengths are made nondimensional using the foil chord length c. Then, Eq. ͑14͒ takes the following nondimensional form:

͓1-C p ͔ͫR ¯d2 R dx ¯2 + 3 2 ͩ dR d x ¯ͪ2 ͬ-1 2 dC p dx ¯R ¯dR d x ¯+ 1 ᐉ W u 2 dp v dT ͑T ϱ -T c ͒ =- C p + v 2 , ͑16͒
where bars denote nondimensional variables.

Convective Approach.

The problem is now to estimate the temperature inside the bubble T c . Two different approaches are considered. In the present section, heat transfer at the bubble wall is supposed to be of convective nature. In other words, the heat flux at the interface is supposed to be proportional to the temperature difference T ϱ -T c :

= h͑T ϱ -T c ͒͑ 17͒
where h is the convection heat transfer coefficient. The validity of the convective approach is open to discussion especially when considering that the bubble is traveling at the same velocity as the liquid. From this viewpoint, the conductive approach considered in Sec. 3.3 might appear more appropriate. Nevertheless, the statement that the heat flux at the bubble wall is proportional to the temperature difference appears physically reasonable. In the absence of any validated correlation applicable to cavitating flows, the heat transfer coefficient is considered here as a tuning parameter.

Under this assumption, the heat balance for the bubble is writen as

d dt ͫ 4 3 R 3 vͬ L =4R 2 h͑T ϱ -T c ͒͑ 18͒
This equation allows us to compute the temperature inside the bubble:

T ϱ -T c = v L h dR dt ͑19͒
Equation ͑19͒ shows that the thermal effect ⌬T = T ϱ -T c is proportional to dR/ dt. In the growing phase, the bubble temperature is then lower than the liquid temperature and is higher during collapse. Using form ͑19͒ of the thermal effect, Rayleigh equation ͑16͒ becomes

͓1-C p ͔ͫR ¯R ¯+ 3 2 R ¯2ͬ - 1 2 C ˙pR ¯R ¯+ 1 Nu ⌬p v * ᐉ ␣ ᐉ ͱ 1-C p R ¯=- C p + v 2 ͑20͒
where Nu is the Nusselt number defined by Nu = hc ᐉ ͑21͒ and = c / W u is the transit time.

According to this approach, the thermodynamic effect introduces an additional term in Rayleigh equation, which is proportional to the bubble growth rate R ¯. The magnitude of thermal effects depends upon the factor:

1 Nu ⌬p v * ᐉ ␣ ᐉ ͑22͒
It is the ratio of two characteristic times, the transit time and the following time:

T =Nu ᐉ ␣ ᐉ ⌬p v * ͑23͒
which can be considered as a thermal time since it essentially contains the information on heat transfer at the interface.

If the thermal time is much larger than the transit time, thermal effects are negligible and Eq. ͑20͒ reduces to the usual Rayleigh equation. Conversely, if the thermal time is much smaller than the transit time, thermal effects are predominant. A suitable criterion for estimating the magnitude of thermal effects is then based on the comparison of these two characteristic times. Besides the Nusselt number, the thermal time depends on ⌬p v * / ᐉ ␣ ᐉ , which is very similar to Brennen's parameter ͑4͒. Equation ͑20͒ will be solved to compute the bubble evolution with thermal effects. As already mentioned, the convection heat transfer coefficient or the corresponding nondimensional Nusselt number is a free parameter of the model. We can also imagine correlating the Nusselt number with the Reynolds ͑Re= cW u ᐉ / ᐉ ͒ and Prandtl ͑Pr= ᐉ / ␣ ᐉ ᐉ ͒ numbers as in traditional convection heat transfer theory, but this has not yet been done for cavitating flows.

Let us observe that Eq. ͑19͒ takes the following nondimensional form:

B = Pe ͱ 1-C p Nu R ¯͑24͒
where Pe is the Peclet number ͑see Eq. ͑27͒͒.

Conductive Approach.

The second approach tested here is the conductive approach. It requires to solve the heat diffusion equation in the liquid surrounding the bubble. An explicit expression of the bubble wall temperature has been proposed by Plesset and Zwick ͓23-25͔͑ see also Brennen ͓3͔͒ under the assumption that the thickness of the thermal boundary layer is much smaller than the bubble radius. The validity of this assumption is evaluated at the end of the present section. If so, bubble wall temperature T c can be calculated at each instant t by the following expression:

B = T ϱ -T c ⌬T* = 1 ͱ ␣ ᐉ ͵ u=0 u=t R 2 ͑u͒͑dR/dt͒͑u͒du ͱ͵ v=u v=t R 4 ͑v͒dv ͑25͒
This equation is the solution of the unsteady heat diffusion problem with a moving spherical boundary. It includes the enhancement coefficient discussed in Sec. 1. It is similar to some extent to Eq. ͑19͒, although more complicated. Both approaches will be compared in detail in Sec. 4.1.

Using the nondimensional procedure given in Sec. 3.1 and switching to a space dependency rather than a time dependency, Eq. ͑25͒ becomes

B = T ϱ -T c ⌬T * = ͱ Pe J͑x ¯͒͑ 26͒
where Pe is the Peclet number defined by

Pe = Re Pr = W u c ␣ ᐉ ͑27͒
and J is the integral

J͑x ¯͒ = 1 ͱ ͵ u=0 u=x ¯R ¯2͑u͒dR ¯/dx ¯͑u͒du ͱ͵ v=u v=x R ¯4͑v͒/ ͱ 1-C p dv

͑28͒

By introducing expression ͑26͒ of the interface temperature into Rayleigh equation ͑16͒, we get the final following equation:

͓1-C p ͔ͫR ¯R ¯+ 3 2 R ¯2ͬ - 1 2 C ˙pR ¯R ¯+ ⌬p v * ᐉ W u 2 ͱ Pe J =- C p + v 2 ͑29͒
When using Brennen's parameter ⌺ defined in Eq. ͑4͒, previous equation takes the equivalent form:

͓1-C p ͔ͫR ¯R ¯+ 3 2 R ¯2ͬ - 1 2 C ˙pR ¯R ¯+ ⌺ͱ c W u 3 1 ͱ J =- C p + v 2

͑30͒

Bubble evolution including thermal effects is obtained by solving Eq. ͑30͒ together with Eq. ͑28͒. The numerical procedure is presented in Sec. 3.4. As mentioned previously, the consistency of the approach depends upon the relative order of magnitude of the thermal boundary layer thickness compared to bubble size, which is typically of the order of the chord length. The order of magnitude of the boundary layer thickness is

␦ Ϸ ͱ ␣ ᐉ ͑31͒ so that ␦ c Ϸ ͱ ͱ Pe = ͱ ͱ Re Pr

͑32͒

A typical value of the Reynolds number for rocket turbopump inducers is 10 6 . For Refrigerant 114 considered later in this work, we have PrХ 5. A typical value for the amplifying coefficient is 10 4 ͑see Sec. 4͒. Hence, we have ␦ / c Ϸ 0.04 so that the hypothesis of a thin thermal boundary layer compared to the characteristic length scale of the cavities appears appropriate. The larger the , the larger the effective thermal diffusivity and the less accurate this hypothesis.

In a few cases and particularly when thermal effects become important, the conductive model presented above may lead to an unrealistic shape of the cavity characterized by oscillations. This effect is generally observed around cavity closure when the bubble starts to collapse, i.e., when the local pressure exceeds the vapor pressure. It is due to the fact that the temperature increases because of condensation as well as the corresponding vapor pressure. If thermal effects are large enough, the vapor pressure may exceed the local pressure. The bubble will then grow again so that the cavity interface oscillates.

To avoid this effect, which is not observed experimentally, it is chosen to set the temperature difference T c -T ϱ to zero as soon as it becomes positive. This is equivalent to assuming that the bubble temperature cannot exceed the liquid temperature. In other words, it is assumed that heat diffusion is considerably enhanced during the collapse phase so that the temperature becomes uniform almost instantaneously. The instability of the spherical shape together with the splitting of the vapor phase into small scale structures during collapse can justify this assumption. The entire growth phase is not affected by this effect. A similar assumption is made by Watanabe et al. ͓7͔, who neglected the heat flux due to condensation near the trailing edge of the cavity, assuming that condensation occurs far downstream of the cascade because of the shedding of vapor structures by the cavity.

For coherence and to make relevant comparisons between the conductive and the convective model, the same assumption of a uniform temperature in the condensation zone is also made for the convective model. The consequence on the temperature distribution in the cavity is clearly visible on Fig. 5͑b͒, for instance.

Numerical Aspects.

Equation ͑20͒ is a nonlinear differential equation of the second order, whereas Eq. ͑30͒ is an integrodifferential equation because of integral ͑28͒. To solve these equations, a combination of two Runge-Kutta methods of orders 4 and 5, respectively, is used. The first provides the solution, and the second an estimate of the error. This combination, due to Fehlberg and improved by Kash and Carp ͓26͔, allows us to control the error during iterations. Indeed, for a specified tolerance, the order of magnitude of the error is calculated and the integration step is modified so that the precision criterion is satisfied: the step value is increased if the error is too small and decreased if the error is too large. Thus, the integration step is variable and the computing time is considerably reduced while keeping an excellent relative precision ͑10 -6 ͒. This method is extended to the resolution of a system of two first order differential equations, using R and dR/ dx as unknown functions. In this case, the relative error is simply the maximum value of the two relative errors.

The iterative resolution procedure makes it possible to transform the integrodifferential equation ͑30͒ into an ordinary differential equation. Indeed, at the current iteration, integral J is calculated using estimates obtained at the previous iteration.

Mechanical and thermal equilibriums are assumed as initial conditions for the bubble. So, the interface velocity ͑dR/ dx͒ 0 and acceleration ͑d 2 R / dx 2 ͒ 0 are initially zero, and the temperature is equal to T ϱ and J 0 =0 ͑Eq. ͑25͒͒. To correctly initiate the iterative procedure, a nonzero value of R 0 is required. It has been checked that the subsequent behavior of the bubble is insensitive to the value of its initial radius within a very large range of variation.

Under these assumptions, strictly speaking, Eq. ͑20͒ or ͑30͒ imposes that C p0 + v = 0 at the initial time. The cavity starting point being fixed, this last condition is generally not satisfied, except for the very special case v =-C p0 . Nevertheless, for all values of v , the preceding initial conditions are preserved. Hence, the interface velocity presents a discontinuity at the initial time. The numerical results show that this discontinuity does not affect the resolution of the differential equation. diffusivity to the molecular one for the conductive model have been chosen so that the cavity length is the same. Figure 5͑b͒ shows that temperature distributions are quite different for both cases.

Results and Discussion

For the convective model, temperature exhibits a sudden drop at the leading edge of the cavity. The temperature depression is maximum there and the temperature increases regularly along the cavity. As shown by Eq. ͑19͒, the convective model assumes that the temperature depression is proportional to the bubble growth rate. Since the latter is maximum at the cavity leading edge where pressure is minimum, the temperature depression is also maximum at this location.

For the conductive model, the temperature depression is zero at the cavity leading edge. The temperature inside the cavity first decreases and then increases so that a minimum is observed. This smooth behavior is the consequence of the solution of the heat diffusion equation in the liquid, which does not allow any jump in temperature.

Because of the large temperature drop at the cavity leading edge, the initial bubble growth rate is smaller for the convective approach. In the conductive case, it is the same as that without thermal effects ͑see Fig. 9͒, since temperature depression is initially zero. As a consequence, the cavity is somewhat thinner for the convective model.

The maximum cavity thickness is located closer to the cavity trailing edge for the convective approach ͑see Fig. 5͑a͒͒. Its location corresponds exactly to the point where the B factor vanishes, still because of the proportionality between B and R ˙. For the conductive model, the relationship between B and R ˙is more complicated ͑see Eq. ͑25͒͒, since the cavity temperature at any time depends upon the whole temperature history of the surrounding liquid. As a consequence, the conductive model predicts a shift between the location of maximum cavity thickness and the point where the B factor vanishes, as shown in Fig. 5͑a͒.

It is difficult to know which model is more appropriate. On the one hand, the instantaneous initial jump in temperature predicted by the convective model may be considered as unrealistic from a physical viewpoint. On the other hand, experimental results ͑see Hord ͓6͔ and Fruman et al. ͓8͔͒ tend to prove that the temperature depression is maximum at the cavity leading edge and decreases regularly downstream. However, the maximum observed experimentally is usually not as pronounced as the convective model predicts and measured temperature profiles are generally rather flat just downstream of the cavity detachment. On the whole, the conductive model may be preferred although no definite conclusion can be drawn from the present work. It will be shown in next section that the differences between both models decrease when thermal effects increase and that both models tend toward the same limit for large enough thermal effects.

Effect of Nu or on Thermal Effects.

Both models include a free parameter, which controls the amplitude of thermal effects. The tuning parameter is the Nusselt number for the convective approach and the ratio for the conductive one. Thermal effects are zero when Nu or are infinite and increase with decreasing Nu or .

The similar influence of Nu and can be understood by comparing both models. The heat flux to the interface is given by

q = h • ⌬T ͑33͒
for the convective model, and by approximately

q Ϸ t • ⌬T ␦ ͑34͒
for the conductive one since ⌬T / ␦ is the order of magnitude of the temperature gradient. The identification of forms ͑33͒ and ͑34͒ of the heat flux together with the estimate ͑31͒ of the boundary layer thickness leads to the following nondimensional relation:

Nu Ϸ ͱ Pe ͑35͒
This equation suggests a link between the conductive and the convective models. To some extent, it can be expected that the influence of the Nusselt number in the convective model is similar to that of ͱ Pe for the conductive one. Figure 6͑a͒ confirms that the variations of cavity length with Nu and ͱ Pe for the convec- tive and the conductive model, respectively, are comparable.

When thermal effects are increased, i.e., when either Nu or ͱ Pe are decreased, cavity length decreases in a similar way. Moreover, both models predict the same minimum cavity length for very large thermal effects.

This minimum cavity length is not zero. However, Fig. 6͑b͒ indicates, and detailed computations confirm, that the cavity thickness as well as the cavity volume comes to zero when thermal effects are increased. The bubble growth rate R ˙together with R Fig.

also vanishes. Then, the inertial term on the left hand side of Rayleigh equation ͑11͒ tends to zero. As a consequence, the right hand side term should also vanish. This implies that, whatever the thermal model may be the vapor pressure p v ͑T c ͒ based on the actual cavity temperature T c tends toward the local pressure at any point along the cavity. This is confirmed on Fig. 7, which clearly shows that the opposite of the cavitation numberc based on the actual cavity pressure approaches everywhere the pressure coefficient C p when thermal effects are increased. In the limit case corresponding to the suppression of cavitation by thermal effects, the temperature distribution inside the cavity is determined only by the pressure distribution on the blade and the vapor pressure curve of the fluid. In particular, the maximum value of the B factor is independent of the model, as shown by Fig. 6͑c͒, and determined by the minimum pressure p min on the blade. The minimum value of the temperature inside the cavity T c min is determined by p v ͑T c min ͒ = p min . The existence of a minimum value of the cavity length when thermal effects become predominant can be qualitatively understood on the basis of the convective model. The limit shape of the cavity is given by the solution of Rayleigh equation ͑20͒ in which the inertial terms are neglected as previously mentioned. This simplified equation writes

R ¯Х - T C p + v 2 ͱ 1-C p ͑36͒ Since ͵ 0 ᐉ R ˙dx =0 ͑37͒
the minimum cavity length ᐉ is determined by the following equation:

͵ 0 ᐉ C p + v ͱ 1-C p dx =0 ͑38͒
It depends only upon the pressure distribution and the cavitation number. In practice, the situation is more complicated since the B factor is forced to zero in the condensation region, where Eq. ͑36͒ is then no longer valid. The situation is even more complex in the conductive case for which equations are much more complicated. Anyway, the previous simplified argument makes the idea of a minimum cavity length quite understandable when thermal effects are increased even though an exact equation for the quantitative prediction of this minimum cavity length could not be derived.

Comparison With Experiments.

Figure 8 presents the evolution of cavity length with the cavitation number for water and Refrigerant 114. The cavitation number plotted here is the inducer cavitation number inducer . It is deduced from the cavita- tion number of the 2D cascade defined in Eq. ͑13͒ by the following relation:

inducer Х 1 2 ͫ r R ͬ 2 ͑39͒
where r is the local radius between hub and casing where the 2D cascade is considered and R is the inducer peripheral radius used as reference length scale for the computation of the inducer cavitation number.

In the case of water, i.e., without thermal effects, the predicted curve ᐉ͑͒ is in reasonable agreement with experiments. Because of several shortcuts in the modeling mentioned in Sec. 1, together with the oversimplification of the real shape of the inducer blades by flat plates, discrepancies are unavoidable. However, they are not prejudicial to the present discussion, which is focused on thermodynamic effect.

The influence of thermal effects computed from both models on cavity length ᐉ͑͒ is shown in Fig. 8 for R114 at two different temperatures ͑20 and 40°C͒. The computation was conducted with the value = 5000 for the conductive model. As for the convective one, the Nusselt number was adjusted to give the same cavity length as the conductive model in the particular case inducer = 0 for R114 at 20°C. This value of the Nusselt number ͑2.085ϫ 10 6 ͒ was then kept constant for all computations including at 40°C.

The comparison of experimental and results shows that, with these values of and Nu, both models lead to a reduction in cavity length due to the thermal effects of the same order of magnitude as that observed experimentally. Both models can then be considered as acceptable for the modeling of thermal effects although further computations are required for an advanced validation. In particular, the choice of the Nusselt number or the parameter remains an open question. The use of a specific correlation between Nusselt number and Reynolds and Prandtl numbers as for classical convection would make the method predictive and give a basis for transposition of experimental results to other fluids and/or other operating conditions.

Figure 9 presents a comparison of cavity shapes obtained with the conductive model, in water ͑without thermal effects͒ and in Refrigerant 114 for two different temperatures. For each computation, the value was adjusted in order to keep the cavity length approximately constant. As already observed, the cavity becomes thinner as thermal effects increase. This is due to a decrease in the bubble growth rate ͑and also in the collapse rate͒. Similarly, Ahuja and Hosangadi ͓16͔ noticed differences in the cavity content computed from their two-phase model. They report that the cavity becomes more porous and frothy as thermal effects increase as it can reasonably be expected assuming a smaller growth rate due to thermal effects.

From Fig. 9, it can be concluded that there is no similarity in cavity shape when the amplitude of thermal effects changes. Strictly speaking, cavities with and without thermodynamic effect are not similar as well as cavities for the same fluid but at two different temperatures. If thermal effects are not negligible, exact scaling requires the conservation of an additional parameter, characterizing thermal effects ͑as ⌺ ͱ c / V 3 introduced in Sec. 1͒ in addition to the usual scaling law. Results can be transposed rigorously only if such an additional scaling law is satisfied. The same conclusion holds for the temperature distribution inside the cavity which cannot be transposed from water to R114 neither between two different temperatures in R114.

The computed maximum temperature drop inside the cavity is plotted as a function of cavity length on Fig. 10 and compared to the experimental results. The experimental procedure for computing the B factor consisted in determining the equivalent vapor pressure drop ͑and so the corresponding temperature drop͒, which would lead to the observed reduction in cavity length in comparison with reference tests in cold water ͓22͔. In the absence of any direct measurement of the temperature distribution inside the cavities, this experimental procedure based on visualizations is the only way to get a practical estimate of the temperature depression in the rotating machinery and make a possible comparison with computational results.

The predicted temperature depression increases with cavity length, as shown by experiments. The conductive approach predicts a rapid increase in the B factor for short cavities, whereas the experimental trend is an almost linear increase of the B factor with cavity length over the whole range of variation of cavity length investigated experimentally. As for the influence of the fluid temperature on thermal effects, the present approach appears to be in reasonable agreement with experiments, since it predicts a negligible influence of the fluid temperature on the nondimensional temperature depression B at least for small enough cavities as observed experimentally.

The conductive model is suitable for a comparison with experimental data in terms of B factor since the temperature distribution in the cavity is rather flat ͑see Fig. 5͒. The minimum temperature can then be considered as representative of a kind of mean cavity temperature as determined experimentally from visualizations of cavitation extent. This is not the case for the convective approach because of the sudden jump in temperature at the cavity leading edge ͑see Fig. 5͒. The temperature is far from being uniform and the comparison with experiments in terms of B factor is difficult. The mean cavity temperature, for example, proved to be inadequate for the comparison. Nevertheless, it has to be noticed that the experimental evolution of the B factor with cavity length shown in Fig. 10 was deduced from the ᐉ͑͒ curves presented in Fig. 8 and then contains basically the same information. Hence, in the absence of direct measurements of temperature inside the cavities, it seems better to base the confrontation of computational and experimental results on the curves ᐉ͑͒ rather than B͑ᐉ͒.

Concluding Remarks

The present work is devoted to an analysis of thermal effects in a cavitating inducer. The rotating machinery is modeled by a 2D cascade of flat plates. The potential flow is computed using a nonlinear approach based on the boundary element method. The development of cavitation is predicted by solving Rayleigh equation using the pressure distribution deduced from the fully wetted flow computation. The interface of the sheet cavity is approximated by the envelope of a hemispherical bubble traveling on the suction side of the blades. Thermal effects in the Rayleigh equation are taken into account by considering the vapor pressure at the actual bubble temperature T c , which is different from the liquid temperature T ϱ far from the bubble.

Two different approaches are evaluated to model heat transfer at the bubble wall and compute the cavity temperature T c .I h e convective approach, the heat flux is supposed to be proportional to the temperature difference ⌬T = T ϱ -T c . The amplitude of thermal effects is controlled by the magnitude of the convective heat transfer coefficient or the corresponding nondimensional Nusselt number. The second approach is based upon the resolution of the heat diffusion equation in the liquid surrounding the bubble. For this conductive-type model, thermal effects are controlled by the ratio of the eddy thermal diffusivity to the liquid molecular thermal diffusivity.

The main conclusions of the present work are the following.

͑1͒

The convective model predicts a sudden drop in temperature at the cavity leading edge followed by a gradual increase up to T ϱ . For the conductive model, the temperature decreases smoothly from T ϱ at cavity detachment, and exhibits a minimum before gradually increasing to T ϱ . Such differences in temperature distributions induce differences in cavity shapes. The convective approach predicts a thinner cavity and the location of maximum thickness is also different. ͑2͒ Thermal effects are increased by decreasing either the Nusselt number for the convective model or the heat transfer enhancement factor for the conductive one. Cavity length decreases and reaches a minimum value, which is the same for both models but which is not zero. The cavity thickness also decreases but tends toward zero. Both models tend toward the same solution when thermal effects are predominant. This solution corresponds to the suppression of cavitation by thermal effects, which is achieved by a local vapor pressure p v ͑T c ͒ equal everywhere to the local pressure on the blade. ͑3͒ The evolution of cavity length with cavitation number is compared to experimental data obtained for an inducer ͓22͔ in the case of cold water ͑i.e., without thermal effect͒ and R114, which exhibits a significant thermal delay at 20 and 40°C. On the whole, both models lead to an acceptable reduction in cavity length as compared to experiments when the Nusselt number or the eddy thermal diffusivity is kept constant.

Such a simple model based on a 2D potential flow computation coupled to the resolution of the Rayleigh equation presents clear limitations due to simplifying assumptions. However, this work shows that they are easy to run and can supply a support for the interpretation of experimental data relative to thermal delays, in particular, in inducers. They can also help in the prediction of basic trends following changes in fluid, operating conditions and geometry.

Nomenclature

  by the classical liquid thermal diffusivity ␣ ᐉ is far too small to account for the measured values of the temperature drop. For a correct prediction, Kato ͓5͔ and Watanabe et al. ͓7͔ have shown that the liquid thermal diffusivity should be multiplied by an amplifying coefficient of the order of 8000 ͓7͔ up to 10 5 ͓5͔. This large value is supposed to account for the enhancement of heat transfer due to both turbulence and an increase in the interface area probably due to the small-scale instabilities of the interface. The importance of the so-called interface roughness was confirmed by Fruman et al. ͓8,9͔ on the basis of Billet's entrainment theory ͓10-12͔.Besides interface models, two-phase models have been developed. The two-phase mixture is generally considered as a pseudofluid, which obeys the Navier-Stokes equations. A suitable constitutive equation for the mixture as the barotropic law can be specified ͑Rapposelli et d'Agostino ͓13͔, Rolland et al. ͓14͔͒.Another option consists in solving an additional transport equation for the vapor phase with a source term modeling the cavitation phenomenon ͑see, e.g., Tani and Nagashima ͓15͔ and Hosangadi and Ahuja ͓16-18͔͒. The vaporization or condensation process can be modeled either by ad hoc equations as in Ref. ͓16͔ or by seeding the flow field with microbubbles and computing their evolution using the Rayleigh-Plesset equation, as proposed originally by Kubota et al. ͓19͔͑ see, e.g., Singhal et al. ͓20͔͒. To address thermal effects, it is necessary to solve the energy equation for the mixture together with the mass and momentum conservation equations. These equations are essentially coupled through the dependence of the vapor pressure with the temperature.

Fig. 1

 1 Fig. 1 Scheme of the 2D cascade "vertical scale has been expanded…

  The bubble is supposed to travel on the blade with the local fluid velocity W u ͱ 1-C p . Time derivatives in the original Rayleigh equation are then transposed into space derivatives using the following equation:

Fig. 2

 2 Fig. 2 Close view of the blade leading edge

4. 1

 1 Comparison Between the Convective and the Conductive Model. Cavity shape and temperature distribution obtained by both approaches are compared on Fig. 5. Values of Nusselt number Nu for the convective model and of the ratio of the eddy

Fig. 5

 5 Fig. 5 Comparison between the convective "Nu= 7.31Ã 10 6 … and the conductive " =5Ã 10 4 … models for the same cavity length. "a… cavity shape. "b… Temperature distribution "R114 40°Cv =0….

Fig. 7

 7 Fig. 7 Effect of increasing thermal effects on the distribution of the cavitation number c based on the vapor pressure at the local cavity temperature p v "T c …. Comparison with the pressure coefficient distribution "R114 40°Cv =0….

Fig. 8

 8 Fig. 8 Cavity length versus cavitation number. Comparison between computation and experiments "Ref. †22 ‡… for cold water and R114 at two different temperatures "20 and 40°C…. Cavity length is made nondimensional using blade spacing H "convective model, Nu= 2.085Ã 10 6 ; conductive model, = 5000….

Fig. 9

 9 Fig. 9 Influence of thermodynamic effect on the shape of the cavity. Results for water and R114 at two different temperatures "20 and 40°C… are compared. For the three cases, the cavity length is kept constant H AE 0.5 "conductive model, = 5000….

  B ϭ B factor of Stepanoff ͑Eq. ͑2͒͒ c ϭ chord length of the blade C p ϭ pressure coefficient c pᐉ ϭ liquid heat capacity e ϭ blade thickness h ϭ convection heat transfer coefficient ᐉ ϭ cavity length L ϭ latent heat of vaporization Nu ϭ Nusselt number based on chord length ͑Eq. ͑21͒͒ p ϭ pressure p u ϭ reference upstream pressure p v ϭ vapor pressure Pe ϭ Peclet number ͑Eq. ͑27͒͒ Pr ϭ Prandtl number R ϭ bubble radius R ¯ϭ nondimensional bubble radius R / c Re ϭ Reynolds number based on chord length t ϭ time T ϭ temperature T c ϭ local temperature in the cavity T ϱ ϭ liquid temperature at infinity V ϭ flow velocity V u ϭ upstream flow velocity W ϭ relative flow velocity W u ϭ upstream relative flow velocity x ϭ curvilinear abscissa on the blade x ¯ϭ nondimensional distance x / c ␣ ϭ Kato's thermodynamic parameter ͑Eq. ͑5͒͒ ␣ a ϭ angle of attack ␣ ᐉ ϭ thermal diffusivity of the liquid ␣ t ϭ turbulent thermal diffusivity ␣ v ϭ void fraction ␦ ϭ thermal boundary layer thickness ⌬p v * ϭ characteristic vapor pressure drop ͑Eq. ͑3͒͒ ⌬T ϭ temperature drop in the cavitating flow T ϱ -T c ⌬T * ϭ characteristic temperature drop ͑Eq. ͑1͒͒ ϭ heat transfer enhancement coefficient ␣ t / ␣ ᐉ or t / ᐉ ᐉ ϭ thermal conductivity of the liquid t ϭ turbulent thermal conductivity ᐉ ϭ water molecular viscosity v ϭ vapor density ᐉ ϭ liquid density ⌺ϭBrennen's thermodynamic parameter ͑Eq. ͑4͒͒ c ϭ local cavitation number based on p v ͑T c ͒

v ϭ cavitation number ͑Eq. ͑10͒͒ ϭ transit time c / W

Effect of Nu "convective model… or "conductive model… on "a… cavity length, "b… cavity shape, and "c… maximum B factor "R114 40°Cv =0….

Acknowledgment

This research was supported by CNES/SNECMA ͑Contract No. 2005-015-I͒. The authors are particularly grateful to H. Horiguchi and Y. Tsujimoto for their help in the validation of the potential flow computation and to B. Chanut for his contribution to this work.