1 Introduction

Nanofluids are engineered colloidal suspensions of nanoparticles in water,ethylene glycol, oil,or other base fluids. Nanofluids have attracted an increasing attention because of extensive applications in energy saving,heat transfer,microelectronics,transportation,pharmaceutical processes and so on. Nanofluids have some special properties of flow and heat transfer,e.g.,they have higher ability of thermal conductivity and convective heat transfer^[1],result in significant increase in viscosity than the base fluid^[2],and exhibit viscosity decrease when the shear rate is increased^[3].

However,the special properties of flow and heat transfer properties for nanofluids depend on the flow state. For example,the heat transfer coefficient is very large^[4],and the proportional increase in pressure drop is very low in turbulent flows^[5]. Thus,the transition of nanofluids from laminar to turbulent is an important research topic. Nanofluids usually experience instabilities when they are driven away from the thermal or mechanical equilibrium. There are some studies related to the thermal instability in nanofluids^{[6, 7, 8, 9, 10, 11]},while the hydrodynamic instability produced by mechanical non-equilibrium in fluids is quite normal. There are few of researches on the hydrodynamic instability produced by mechanical non-equilibrium in nanofluids. Therefore,we explore the hydrodynamic instability of nanofluids in jet flow in the present study.

We choose jet flow because of its numerous industrial applications. For many industrial processes,the distribution of the nanoparticle is an important factor which controls the efficiency and the instability of the processes. Several studies have been devoted to the hydrodynamic instability of particulate flow with micron-sized particles and gas-liquid flow in a jet. Sykes and Lyell^[12] showed that the particles enhance the stability of the flow. Parthasarathy^[13] indicated that based on the temporal stability analysis,the particles leads to decrease in the wave amplification and increase in the wave velocity,respectively. Wave amplification and velocity decrease with the increase of the particle mass loading. The particles decrease the wave amplification rate based on the spatial stability analysis. Lin and Zhou^[14] indicated that the particles have an effect on the flow instability significantly in moving particulate jet flow with dense particles. DeSpirito andWang^[15] found that the particles can make the flow unstable when the Stokes number is small,while make the flow stable when the Stokes number changes from an intermediate value to a large one. When the particle mass loading and wavenumber remain unchanged,the Stokes number on the order of 1 corresponds to a maximum stability, which also depends on the wavenumber. Chan et al.^[16] showed that the disturbance with the axisymmetric mode is usually more stable than that with the first azimuthal mode. The particles enhance the flow stability,and the flow will become more stable with the increase of the particle mass loading. Xie et al.^[17] showed that the critical Stokes number which affects the flow stability is about 2.

Although there are some works on the hydrodynamic instability of particulate flow with micron-sized particles in a jet as mentioned above,the literature related to the hydrodynamic instability of particulate flow with nano-sized particles in jet flow has not been found. Therefore, the aim of this study is to explore the instability properties of nanofuilds in a jet,and to assess the effects of the shape factors of the velocity profile,Reynolds number,particle loading, Knudsen number,and Stokes number on the hydrodynamic instability of nanofluids. 2 Mathematical equations 2.1 Primitive equation

The particulate flow with nano-sized particles is shown in Fig. 1,where the cylindrical coordinate system is used. The flow is assumed to be impressible,and the effect of particle on the base fluid can be neglected because the particle mass loading is low. The particle diameter is so small that the Cunningham slip correction can be used to modify the Stokes drag expression^{[13, 15]}. The basic equations for the axisymmetric jet are as follows:

where p is the base fluid pressure,u and v are the base fluid velocity and particle velocity, respectively,μ is the base fluid viscosity,and ρ_f and ρ_f are the fluid and nanoparticle density, respectively. The momentum coupling is described by the term 3πμNd (u−v)/(ρ_fC_c)^[18] on the right-hand side of (2),where N is the particle number density,d is the nanoparticle diameter, and C_c=1+1.591Kn is the Cunningham slip correction for the Stokes drag with Kn being the Knudsen number (where Kn=2λ/d,and λ is the mean free path length of the gas molecules). The particle volume fraction α is Nπd³/6.

Using the base fluid velocity U₀ at the centerline,the radius r₀ of the jet core,and the mean particle volume fraction α₀ as the reference variables,we can non-dimensionalize (1)-(4) as follows:

where Re=U₀r₀/ν is the Reynolds number,Z is the mean particle mass loading (=α₀ρ_p/ρ_f ), and the Stokes number is given by St=(ρ_pd²/(18νρ_f ))/(r₀/U₀). 2.2 Stability equations

The linear stability theory is still an important approach to study the instability of jet flow. Based on the linear stability theory,the basic variables consist of the mean component and the perturbation component (with ) as follows:

where P is the steady mean pressure of the flow,U and V are the steady mean velocity of the base fluid and particle,respectively,and α₀ is equal to 1. In a steady state flow,the nanoparticles move with the same velocity as the fluid because the sedimentation velocity can be neglected for nanoparticles. Thus,we have U=V=(0,0,U_z). Substituting (9) into (5)-(8), then subtracting the mean flow terms from the equations,and finally neglecting the higher order terms of the perturbation,we have the linear stability equations as follows:

We can give the perturbations of a traveling wave as follows:

where n is the azimuthal mode,u_r(r),u_θ(r),u_z(r),v_r(r),v_θ(r),v_z(r),and p(r) are the amplitudes, β is the axial wavenumber,c is the wave amplification factor,and u'_r is taken to be proportional to iu_r(r) because the phase of u'_r differs by π/2 from those of u'_θ,u'_z,and v'_r. There are two kinds of instability modes,i.e.,the spatial instability and temporal instability. Here, we consider the temporal instability,therefore β is a real quantity,and c = c_r +ic_i is a complex number. Here,c_r is the phase speed,and c_i is the amplification factor,c_i > 0 means that the perturbation will grow with time,and the flow is unstable. c_i < 0 corresponds to a decaying perturbation,and hence the flow is stable. c_i = 0 means a neutral perturbation.

Substituting (18) into (10)-(17) and rearranging the equations,we obtain

where D and D_* are defined as follows:

Equations (19)-(25) form a closed eigenvalue problem and can be simulated numerically with the corresponding boundary conditions. 2.3 Boundary conditions for equation

The boundary conditions of particle perturbation are the same as that of the base fluid phase. Therefore,the boundary conditions of perturbation for different azimuthal modes are as follows.

3 Computation 3.1 Computational method

Equations (19)-(25) form a system of linear homogeneous equations shaped as Mφ=cLφ, where M and L are the constant coefficient matrices,and φ=(u_r,u_θ,v_r,v_θ,v_z)^T is a vector matrix to be solved. The flow instability study can be converted to a generalized eigenvalue problem of ordinary differential equations because there should be a nonzero solution for the homogeneous linear equation system. To solve the equations numerically,the continuous equations and boundary conditions need to be discretized,and then a central difference scheme is used to solve the discretized equations. To solve the eigenvalue problem,we use the QZ algorithm proposed firstly by Moler and Stewart^[19] and accepted widely as a valid serial algorithm to solve the eigenvalue problem of asymmetric matrices. 3.2 Velocity profile shape factor

In the hydrodynamic stability analysis of flow,a velocity distribution for steady state flow should be given in advance. As far as jet flow is concerned,a normalized velocity profile is^[13]

where B=R/θ (R and θ are the middle line and the momentum boundary layer thickness of the jet shear layer,respectively) is the shape factor of the velocity profile,and represents the velocity profiles at different axial locations. The smaller B is,the farther the axial location is from the jet exit. Figure 2 shows the velocity profiles for different shape factors. 3.3 Definition of outer boundary

The outer boundary in the radial direction r is located at infinity theoretically,but it cannot be reached in the actual numerical simulation. Therefore,the following assumption is made, i.e.,the outer boundary is defined at r = R. The wave amplifications of a base fluid jet are shown in Fig. 3,from which it can be seen that the value of c_i at R =6 is almost the same as that at R=10,which illustrates that the outer boundary can be defined at R=6.

4 Results and discussion 4.1 Verification of formulations and schemes

In order to verify the numerical formulations and schemes,we have simulated numerically the top-hat jet flow of base fluid and compared the numerical results with the theoretical results given by Batchelor and Gill^[20] as shown in Fig. 4,where the relationship of wave amplification and perturbation wavenumber for different azimuthal modes is

with where I_n and K_n are the first and second kinds of modified Bessel functions,respectively. From Fig. 4,we can see that the numerical results agree with the theoretical ones excellently. 4.2 Effect of shape factor of velocity profile

The relation of wave amplification and wavenumber for different shape factors of the velocity profile and azimuthal modes is shown in Fig. 5. It can be seen that the wave amplification increases with increasing B,which means the flow is more unstable when B is larger. This may be attributed to,as shown in Fig. 2,a larger shear force because of a sharper velocity changes from the centerline of the jet to the periphery,which makes the flow more unstable. We can see that the wave amplification increases first and then decreases with the increase of wavenumber. The range of perturbation wavenumber making the flow unstable increases with increasing B. There exists a maximum of wave amplification at a certain wavenumber,which implies that the flow is most unstable at this wavenumber. The wavenumber corresponding to the maximum of wave amplification becomes large with the increase of B. In other words,the more rapidly the velocity changes from the centerline of the jet to the periphery,the larger the wavenumber making the flow more unstable is. From Figs. 5 and 6,we can see that the curves of c_i-β are different for n=0 and n=1. The maximum of wave amplification,for a fixed value of B,is larger for n=1 than that for n=0,which means that the perturbation with the first azimuthal mode makes the flow unstable more easily than that with axisymmetric azimuthal mode. This conclusion was also derived in the single-phase jet flow^[21]. For n=0,the region of wavenumber corresponding to the maximum of wave amplification is in the interval (0.5,1.5) when the value of B varies from 2.5 to 10. While for n=1,the corresponding interval is (0.25,0.75),which indicates that the wavenumbers corresponding to the maximum of wave amplification are more concentrated for the perturbation with axisymmetric azimuthal mode.

4.3 Effect of Reynolds number

Figure 7 shows the relation of wave amplification and wavenumber for different Reynolds numbers Re and azimuthal modes n. The wave amplification increases and gradually attains the inviscid solution with increasing Re,i.e.,the flow is more unstable with larger Re,which is in accord with the mechanism of hydrodynamic instability. The region of wavenumber making the flow unstable increases with increasing Re. The larger Re is,the smaller the wavenumber making the flow unstable is,which is different from the case with different shape factors of the velocity profile as shown in Fig. 5. Comparing Figs. 7(a) and 7(b),we can see that the curves of c_i-β are different for n=0 and n=1. The maximum of wave amplification,for a fixed Re,is larger for n=1 than that for n=0. For n=0,the region of wavenumber corresponding to the maximum of wave amplification is in the interval (0.5,1.0) when Re varies from 200 to 2 000. However,for n=1,the corresponding interval is (0,0.5).

The relation of wave amplification and Re for different azimuthal modes of perturbation is shown in Fig. 8. The wave amplification increases first drastically and then slowly with the increase of Re,and the variation of wave amplification is sensitive to Re when Re is less than 400. It is obvious that the perturbation with the first azimuthal mode makes the flow unstable more easily than that with the axisymmetric azimuthal mode when Re is the same.

4.4 Effect of particle mass loading

The relation of wave amplification and wavenumber for different values of particle loading Z and azimuthal mode n is shown in Fig. 9. The wave amplification decreases with the increase of Z,which means that particles suppress the flow instability and the flow becomes more stable with increasing Z. For the given values of parameters as shown in Fig. 9 (i.e.,Re=1000,B=5, St=1,and Kn=1),the wave amplifications reduce to zero when Z increases to 1.25,which means flow instability will never happen. Thus,there exists a critical particle loading beyond which the flow is stable. The range of wavenumber making the flow unstable decreases with increasing Z. The larger Z is,the smaller the wavenumber making the flow unstable is. From Figs. 9 and 10,it can be seen that the maximum of wave amplification,for a fixed Z,is larger for n=1 than that for n=0. Therefore,the perturbation with the first azimuthal mode makes the flow unstable more easily than that with axisymmetric azimuthal mode when Z is the same. For n=0,the region of wavenumber corresponding to the maximum of wave amplification is in the interval (0.5,1.0) when Z varies from 0 to 0.3,while the corresponding interval is (0,0.5) for n=1.

4.5 Effect of Stokes number

The relation of wave amplification and Stokes number St for different azimuthal modes is shown in Fig. 11,from which we can see that the relation of c_i with St depends on the value of St. The wave amplification decreases first with St when St changes from 0.01 to 1,and then increases with St when St changes from 1 to 100. There exists a turning point St_c (St_c=1 for Re=1 000,B=5,Z=0.3,Kn=1,and β=1),at which the value of wave amplification is the smallest,i.e.,the flow is relatively stable. The perturbation with the first azimuthal mode makes the flow unstable more easily than that with the axisymmetric azimuthal mode when other parameters remain unchanged.

4.6 Effect of Knudsen number

The Knudsen number can be used to determine whether the Stokes drag should be modified. For a large Knudsen number,i.e.,the particle sizes are comparable with the mean free path length of the gas molecules,the Cunningham slip correction is usually used to modify the Stokes drag. The relation of wave amplification and wavenumber for different Knudsen numbers Kn and azimuthal modes n is shown in Fig. 12. It can be seen that the variation of wave amplification with Kn depends on the value of Kn. The wave amplification increases with Kn when Kn changes from 0.5 to 1,while decreases with Kn when Kn changes from 1 to 2. There exists a turning point Kn_c (Kn_c=1 for Re=1 000,B=5,Z=0.3,and St=1),at which the value of wave amplification is the largest,i.e.,the flow is most unstable. Because Kn is related to the mean free path length of the gas molecules and the particle diameter,while the former usually remains unchanged,it implies that there exists a particle size with which the flow is most unstable. Comparing Figs. 12(a) and 12(b),one can find that the perturbation with the first azimuthal mode makes the flow unstable more easily than that with axisymmetric azimuthal mode when Kn is the same. The regions of perturbation wavenumber making the flow unstable are almost the same for n=0 and n=1. For n=0,the region of wavenumber corresponding to the maximum of wave amplification is in the interval (0.5,1.0) when Kn varies from 0 to 2, while the corresponding interval is (0.2,0.5) for n=1.

Figure 13 shows the wave amplification as a function of Kn for different azimuthal modes of perturbation. The wave amplification decreases first and then increases with the increase of Kn,and the minimum of wave amplification appears around Kn=1.

5 Conclusions

An analysis for the hydrodynamic instability of nanofluids in a jet flow is performed. The instability equations are derived and solved numerically. Some present results are compared with the available theoretical results. The effects of shape factors of the velocity profile,Reynolds number,particle loading,Knudsen number,and Stokes number on the flow instability are analyzed. The main conclusions can be summarized as follows.

(Ⅰ) The presence of nanoparticles suppresses the flow instability. There exists a critical particle loading beyond which the flow is stable. As the shape factor of the velocity profile and the Reynolds number are increased,the perturbation wave amplification and the region of wavenumber making the flow unstable increase so that the flow becomes more unstable. However,as the particle mass loading is increased,the wave amplification and the region of wavenumber making the flow unstable decrease so that the flow becomes more stable. The wavenumber corresponding to the maximum of wave amplification becomes large with the increase of the shape factor of the velocity profile and with the decrease of the particle mass loading and the Reynolds number.

(Ⅱ) The variations of wave amplification with Stokes number and Knudsen number depend on the values of the Stokes number and the Knudsen number,respectively. There exists a turning point of the Knudsen number Kn_c and the Stokes number St_c. The flow is relatively stable and most unstable when St_c=1 and Kn_c=1,respectively,when other parameters remain unchanged.

(Ⅲ) The perturbation with the first azimuthal mode makes the flow unstable more easily than that with the axisymmetric azimuthal mode. The wavenumbers corresponding to the maximum of wave amplification are more concentrated for the perturbation with the axisymmetric azimuthal mode.