ANALYTICAL AND NUMERICAL SOLUTIONS FOR THE FLOW IN MICROTUBE WITH THREE-DIMENSIONAL CORRUGATED SURFACE, PART 1: STEADY FLOW*

WANG Hao-li, YANG Meng

College of Metrological Technology and Engineering, China Jiliang University, Hangzhou 310018, China, E-mail: whl@cjlu.edu.cn

WANG Yuan

Department of Fluid Engineering, School of Energy and Power Engineering, Xi’an Jiaotong University, Xi’an 710049, China

ANALYTICAL AND NUMERICAL SOLUTIONS FOR THE FLOW IN MICROTUBE WITH THREE-DIMENSIONAL CORRUGATED SURFACE, PART 1: STEADY FLOW*

WANG Hao-li, YANG Meng

College of Metrological Technology and Engineering, China Jiliang University, Hangzhou 310018, China, E-mail: whl@cjlu.edu.cn

WANG Yuan

Department of Fluid Engineering, School of Energy and Power Engineering, Xi’an Jiaotong University, Xi’an 710049, China

(Received October 20, 2009, Revised July 2, 2010)

The influences of three-dimensional corrugated wall on the fully-developed steady no-slip flows in microtube are studied by analytical and numerical methods in this article. Detailed analytical solutions for the space-averaged equations and the numerical method for the solutions of the disturbance equations are given. An iterative arithm of coupled equations with respect to space-averaged velocities and disturbance velocities is suggested. The study shows that a three-dimensional subsidiary stress layer exists in the near-wall region. The relative roughness, wavenumber and Reynolds number are three important parameters influencing the subsidiary stresses and the space-averaged pressure drop. The space-averaged pressure drop subject to an invariable flow rate mainly depends on the position of datum surface. When the datum surface is taken at the balance position of wall function, the value of pressure drop is determined by the dynamic characteristics of the flow.

microtube, three-dimensional corrugated surface, space-averaged method, subsidiary stress

1. Introduction

Flows at microscales have been widely concerned subjects for their abroad application in biology, chemistry and industry since early 1980’s[1]. Compared with the flows at macroscales, the flows at microscales (the scales between 1 μm to 1 mm) are more complicated[2,3]. In order to research the effects of scales on microflows, the Knudsen number, Kn, has been introduced to estimate gas flow states at microscales, in which no-slip, slip, transition and free molecular flow regimes are divided. The gas flows will drop into no-slip regime at the scales ofcollocation method. The coupled equations between the space-averaged and the disturbance flows are solved iteratively.

Fig.1 Schematic of microtube with three-dimensional corrugated wall

2. Mathematic formulation

Consider an infinitely long microtube whose inner surface have distributions of the three-dimensional corrugated elements (see Fig.1). The wall rough function, F( x ,θ) , is assumed as a two-dimensional periodic function. Choose their radius, R, as the characteristic length and the bulk averaged velocity, U = Q /(π R2) as the characteristic velocity (Q, the bulk flow rates), then the time-independent dimensionless Navier-Stokes equations of laminar flow in microtube and the boundary conditions are given in Eqs.(1) and (2)

where u=(u, v, w)is the velocity vector in cylindrical coordinates, and three components are theaxial, radial and azimuthal velocities with respect to the coordinates x, r and θ, respectively, p the pressure, p = p( x, r,θ) , F( x,θ) the function of rough wall, F( x,θ) = O (1), ε the amplitude of rough function ε?1, Re the Reynolds number, Δ the Laplace operator in cylindrical coordinates. The non-physical wall surface at r=1 are termed as the datum surface in this study.

The wall rough function, F( x,θ), is represented in terms of a constant plus an oscillatory function, given as

where C0denotes the distance between the datum surface and the balance position of rough function, f( x,θ) the oscillatory function.

The velocity and pressure of flow in the rough-wall microtube can be decomposed into space-averaged components plus their disturbance components as follows:

where U0=(U0,V0)and Φ are the space-averaged velocity and pressure drop calculated in Eq.(5), u′, p′ are the disturbance velocity and pressure, respectively, p0( r) the space-averaged pressure on the reference section of x=0. It is known that p0=const for the Hagen-Poiseuille flow.

As ε?1, the boundary conditions are expanded into the Taylor series at datum wall surface r=1 as follows:

Using the standard methods, the velocity and pressure are decomposed into perturbation series as

Substituting Eqs.(3) and (4) into Eq.(7), we have

We can see that the second term on the left hand side of Eq.(9) is a space-averaged quantity termed as slip velocity and the third as well as the forth terms are disturbance velocities. The fifth term is a product of two disturbance quantities which still can be decomposed into a space-averaged quantity plus a disturbance one. According to the Fourier method, disturbance quantities are decomposed as

where α, β are the wavenumbers in the axial and azimuthal directions respectively and c.c. the complex conjugate. Thus we can gain the following relationship:

where u~mn= [u~mn,v~mn,w~mn]is the Fourier coefficient of disturbance velocity. We can see that the first term on the right hand side of Eq.(11) is also the space-averaged velocity to be termed as additional slip velocity, and the second term is additional disturbance velocity. The additional slip velocity is denoted as uwhere.

All the terms in the Navier-Stokes Eqs.(1) are spatially averaged according to Eq.(6), and we obtain the space-averaged equations

The space-averaged boundary condition at r=1 can be written as

Subtracting the space-averaged equations from the time-independent Navier-Stokes Eqs.(1) and eliminating the terms of higher thanεyields the disturbance equations and their boundary condition

Equations (12)-(16) compose a set of coupled equations.

3. Solution methods

3.1 Solution for space-averaged velocities

By expanding Eq.(12), the space-averaged equations are reduced to a set of ODEs in terms of space-averaged axial, radial velocities and pressure as follows:

where U0, V0are the space-averaged axial and radial velocities, respectively.

The space-averaged radial velocity V0( r) has the zero solution from the continuity equation and the boundary conditions. Thus Eqs.(17a) and (17b) could be reduced to

Substituting Eq.(21) into (22) gives thespace-averaged pressure drop Φ

The radial distribution of space-averaged pressure p0( r) is approached by integrating Eq.(20) and the solution is

3.2 Solutions for disturbance equations

Expanding Eqs.(15)-(16), we can obtain the explicit forms of disturbance equations

The explicit forms of boundary conditions are deduced as follows:

The solutions for the disturbance equations can be approached numerically by the infinite difference method or spectral method. Due to its higher accuracy, the spectral collocation method is employed in this study. The Fourier coefficients of disturbance velocity and pressure are decomposed in spectral space as

where J is the number of collocation points. By using the transform of r =(1 -ζ) /2, the variable, r, in physical space is transformed into that of the Chebyshev space, ζ.

In this study, three-dimensional rough functions are given by the cosine functions as

where C0is the space of the datum surface deviating from equilibrium position of rough wall, A the amplitude of rough wall function. As the relative wall roughness defined here is the difference between the peak and trough of rough wall, the constant, A, is taken as the value of 1/2 in this study for the sake of being equal to the small parameter, ε.

From the boundary conditions given in Eq.(16), it can be seen that there are not the sub-harmonic components for this kind of rough walls, so according to Eqs.(10) and (27) the disturbance velocity and pressure are written as

Substituting Eq.(30) into Eqs.(15) and (16) and separating the Fourier modes lead to the discretization equations and their boundary conditions in spectralspace. The discretization equations at collocation points have total 4J+4 unknown quantities. The explicit forms of boundary condition at r=0 are dependent on the values of β as follows:

The equation of the radial velocity is expanded as

which is considered as the boundary condition ofat r=1 in this study.

3.3 Solution for the coupled equations

The coupled equations include two ODEs with respect to the axial and swirling velocities and four PDEs. Solutions for the coupled equations are found by calculating u′ and p′ using a typical parabolic velocity profile, i.e., U?=2(1 - r2), and the initial

0iterative solutions u ′?, p′?, τ?and Φ?can be obtained respectively. Substituting τ?and Φ?into the space-averaged equations yields the correctional solutions of space-averaged velocities U0??. In turn, the next iterative solutions, u′??, p′??, τ??and Φ??are approached. This process is repeated until all the quantities simultaneously attain some high degrees of accuracy. In the process of calculating the space-averaged velocities through Eqs.(21), (23) and (26), the numerical integration method is employed. The convergent criteria for coupled equations is given as where Θ denotes one of the physical quantities, k the iterative index, ω the convergent accuracy which is taken as 10?6in this study.

Calculations show that the convergent speeds are satisfactory at low and moderate Reynolds numbers (Re＜ 103) and small relative wall roughness (ε＜10%). In general, Re could be fairly large for smaller relative wall roughness and vise versa. The convergent speeds are proved quite high when the relative wall roughness taken as lower than 5%and the Reynolds number lower than 500 no matter what values of wavenumbers and C0are taken. However, it is difficult to converge for a larger Reynolds number as relative roughness increases and vice versa. In addition, the convergent speeds are tightly related to the parameter C0if the relative wall roughness and Reynolds number are taken comparatively large (equal to or greater than 7.5% for relative wall roughness and 500 for the Reynolds number, say). For a non-negative value of C0, the convergent speed is comparatively high, and the computational results can reach its accuracy no more than twenty steps. For a negative value of C0, however, the convergence becomes very slow. For examples, if C0= -1 /2 and the relative wall roughness and Reynolds number are taken as 7.5% and 1 000, respectively, the iterative steps are more than 100. This indicates that an increase of velocity based on an invariable flow rate will lead to an increase of non-linear effect.

4. Results and discussion

4.1 Subsidiary stress layers and space-averaged pressures

The distributions of subsidiary stresses in the axial, radial and azimuthal directions are illustrated in Figs.2-4, respectively. Because the influences of rough wall on flows mainly exist in the regions near wall, subsidiary stress layers exist near the regions of wall as a matter of fact. In subsidiary stress layers, the fluid gain its subsidiary drive or drag forces in the three coordinate directions, which makes the flow patterns greatly different from those flows in smooth or two-dimensional rough tubes. At the first glance, we can see that there are different configurations of stress layers in different components. As is illustrated in Fig.2, the axial stresses go through two steps in the process far away from the datum surface: all the axial stress curves slope down from zero value to extreme values and then slope up to zero. Hence, the values of all axial stresses drop into a negative region, which indicates drag forces impose on the axial flows. As is illustrated in Fig.3, the radial stresses go through three steps in the process far away from datum surface: all the radial stress curves slopes up from zero to positiveextreme values and then slope down into the negative region. After arriving at negative extreme values, the radial stresses begin to slope up and soon disappear. Figure 4 presents the distributions of azimuthal stresses, and we can see that the azimuthal stress curves are different from radial stress curves in the process far from datum surface.

Fig.2 Influences of the wall roughness on subsidiary stresses in axial direction

The amplitudes and thicknesses of subsidiary stress layers with different wall rough parameters and Reynolds numbesr can be analyzed from these figures. As are illustrated in Figs.2(a), 3(a) and 4(a), the amplitudes of subsidiary stresses increase with the augments of wall rough amplitudes, but the thicknesses of subsidiary stress layeres are invariable fundamentally. This indicates that the amplitudes of subsidiary stresses are sensitive to the wall rough amplitudes, while the thicknesses of stress layers are not. It can be seen from Figs. 2(b), 3(b) and 4(b) that the amplitudes of subsidiary stresses increase with the augments of the wall wavenumbers, but the subsidiary stress layers become thin gradually. On one hand, the increases of the wall wavenumbers lead to the augments of shear rates in the near-wall region due to their high disturbance frequencies, and on the other hand the average kinetic energy of main flow increases with increase of wavenumber. The former makes the subsidiary stresses increase in the near-wall regions and the latter makes wall disturbances tend to be weakened. The influences of the Reynolds numbers on the subsidiary stresses are illustrated in Fig.4. It can be seen that the subsidiary stress layers become thin gradually with increasing Reynolds numbers. This can be explained as that the increase of the Reynolds numbers leads to the increase of the average kinetic energy of the main flow, so that the effects of rough wall on the flow field are weakened.

Fig.3 Influences of the wall roughness on subsidiary stresses in radial direction

It is known that the cross-section pressure distribution is a constant for the Hagen-Poiseuille flow. For the laminar flows in rough microtubes, however, the space-averaged pressure distributions (pressure distributions in short) are the functions of radius. The pressure distribution functions for different wall and flow parameters are presented in Fig.5, where the y-coordinate denotes the pressure distributions. We can see that the pressuredistributions slope down from a positive region to a negative region at first, and then slope up to zero (the cross section is the section of x=0). The pressures arrive at their maximum values on the datum surface.

Fig.4 Influences of wall roughness on subsidiary stresses in azimuthal direction

4.2 Space-averaged pressure drop

Despite the fundamental simplicities of laminar flows in straight microchannels, experimental studies of microscale flows have often failed to reveal the expected relationship between the pressure drop (or friction factors) and Reynolds number. Further, flow discrepancies are neither consistently higher nor lower than macroscopic predictions. Sharp[19]presented some typical experimental results to illustrate the relationship between the pressure drop, Φ?, defined as Eq.(34) and the Reynolds number. The experimental data indicate that the pressure drop falls into a wide range of 0.5-2.5 under the conditions of the Reynolds numbers of 10?3to 103.

Fig.5 Influences of wall roughness on pressure distribution at cross section

In order to present a reasonable explanation for the diversity of measurement results of pressure drop, the Space-Averaged Pressure Drops (SAPDs) are calculated for different datum surfaces based on a constant flow rate. We can see from the Eq. (23) that the SAPD, Φ, composes two parts. One is only related to the geometrical configuration of the wall rough function, written as

And another is the dynamic part,DΦ, written as

where the quantities of O(ε3) in Eq.(23) have been ignored. We have

Fig.6 Influences of wall roughness on the space-averaged pressure drop

Obviously, if C0=0, the SAPD is only related with the dynamic correlated part,DΦ.

The variations of SAPD with the increase of the Reynolds number for different parameters of wall roughness are illustrated in Fig.6, where the x-coordinate is the logarithm of the Reynolds number and y-coordinate denotes the SAPD. At first, we study

The case of C0=0 is another concerned problem in this study, because the SAPD in this case is merely dependent on the dynamic correlation part. C0=0 means the datum surface to be taken at the balance position of the wall rough function. We can see from Fig.6(a) that the values of SAPD are approximately equal to the theoretical solution of the Hagen-Poiseuille flow for the small relative roughness (ε≤5%) and small wavenumbers (α = β= 5). Even for large wavenumbers, large variations of SAPD are still not found for small relative roughness, which can be seen in Fig.6(c), where the wavenumber is taken as 30. Hence, we can see that the dynamic correlation part of SAPD is insensitive to the wall parameters for the small relative roughness. However, for larger relative roughnesses (ε=7.5% and 10%, say), there are some differences. We can see from Fig.6(b) that the lines of SAPD have deviated from that of the Hagen-Poiseuille flow for small wavenumbers (α = β =5). From Fig.6(d) we can see that the values of SAPD for the relative roughness of 10% have great augments from theoretical solution with the increase of wavenumber, which can be attributed to the following two aspects. One is that the drag forces of flowing around the rough elements have some increases due to non-linear effect intensified under the condition of larger relative roughness, so the SAPD increases even under the small wall wavenumbers. Another reason is that greater wavenumbers lead to higher disturbance frequency, so that the effect of theviscous dissipation resulting from a higher wall shear action increases. This brings a conjecture for us that the wall roughness will give rise to an important effect on the viscous heating for microtube flows because a virtual surface includes numerous high wavenumber components according to Fourier’s theory.

5. Conclusion

The influences three-dimensional rough wall on the laminar microtube flows have been analyzed in this article. All physical quantities are decomposed into space-averaged flows plus disturbances owning to the presence of three-dimensional rough wall. The space-averaged equations and the disturbance equations are solved by the analytical methods and the spectral collocation method, respectively. A set of coupled equations are approached in an iterative arithmatic. Analytical and numerical results indicate that flows in three-dimensional rough wall microtubes have the following three main characteristics: three-dimensional subsidiary stress layers exist near wall; the laminar pressure drops in microtubes may be present three possibilities, i.e., higher than, equal to and even less than the solution of Hagen-Poisueille flow. These characteristics are influenced by the wall rough parameters in terms of wall relative roughness, wavenumbers as well as the Reynolds number.

[1] WU P. Y., LITTLE W. A. Measurement of friction factors for the flow of gases in very fine channels used for microminiature Joule-Thomson refrigerators[J]. Cryogenics, 1983, 23(5) : 273-277.

[2] GUO Qiang, CHENG Rui and SILBER-LI Zhan-hua. Influence of capillarity on nano-liter flowrate measuremet with displacemet method[J]. Journal of Hydrodynamics, Ser. B, 2007, 19(5): 594-600.

[3] LIU Zhong-chun, YUE Xiang-an and HOU Ji-rui et al. Experimental study of microscale flow for micro molecule liquid and polymer solution[J]. Journal of Hydrodynamics, Ser. B, 2005,17(3): 352-357.

[4] GAD-EL-HAK M. The MEMS handbook[M]. Boca Raton, USA: CRC Press, 2002.

[5] MALA G. M., LI D. Q. Flow characteristics of water in microtubes[J]. Int. J. Heat Fluid Flow, 1999, 20(2): 142-148.

[6] PAPAUTSKY I., BRAZILE J. and AMEEL T. et al. Laminar fluid behavior in microchannels using micropolar fluid theory[J]. Sensors and Actuators A, 1999, 73(1): 101-108.

[7] CELATA G. P., CUMO M. and GULIELMI M. et al. Experimental investigation of hydraulic and single phase heat transfer in 0.130 mm capillary tube[J]. Nanoscale and Microscale Thermophysical Engineering, 2002, 6(2): 85-97.

[8] LI Z. X., DU D. X. and GUO Z. Y. Experimental study on flow characteristics of liquid in circular microtubes[C]. Proceedings of the International conference on Heat Transfer and Transport Phenomena in Microscale. Ban, Canada, 2000, 162-167.

[9] STANLEY R. S. Two-phase flow in microchannels[C]. Ph. D. Thesis, Ruston, L.A., USA: Louisiana Technology University, 1997.

[10] CHEN Zhou, QIAN Jia-zhong and LUO Shao-he et al. Experimental study of friction factor for groundwater flow in a single rough fracture[J]. Journal of Hydrodynamics, 2009, 21(6): 820-825.

[11] KLEINSTREUER C., KOO L. Computational analysis of wall roughness effects for liquid flow in micro-conduits[J]. J. Fluids Eng., 2004, 126(1): 1-9.

[12] KOO L., KLEINSTREUER C. Liquid flow in micro-channels: Experimental observations and computational analyses of microfluidics effects[J]. J. Micromech. Microeng., 2003, 13(5): 568-579.

[13] RAWOOL A. S., MITRA S. K. and KANDLIKAR S. G. Numerical simulation of flow through microchannels with designed roughness[J]. Microfluids and Nanofluids, 2006, 2(3): 215-221.

[14] TUCK E. O., KOUZOUBOV A. A laminar roughness boundary condition[J]. J. Fluid Mech., 1995, 300: 59-70.

[15] WANG Hao-li, WANG Yuan and ZHANG Jia-zhong. Influence of ribbon structure rough wall on the microscale Poiseuille flow[J]. J. Fluids Eng., 2005, 127(6): 1140-1145

[16] WANG Hao-li, WANG Yuan. Flow in microchannels with rough walls: Flow pattern and pressure drop[J]. J. Micromech. Microeng., 2007, 17(3): 586-596

[17] WANG Hao-li, WANG Yuan. Influence of threedimensional wall roughness on the laminar flow in microtube[J]. Int. J. Heat Fluid Flow, 2007, 28(2): 220-228.

[18] GAMRAT G., FAVRE-MARINET M. and PERSON S. L. An experimental study and modelling of roughness effects on laminar flow in microchannels[J]. J. Fluid Mech., 2008, 594: 399-423.

[19] SHARP K. V. Experimental investigation of liquid and particle-laden flows in microtubes[C]. Ph. D. Thesis, Urbana-Champaign, USA: University of Illinois, 2001.

10.1016/S1001-6058(09)60099-8

* Project supported by the National Natural Science Foundation of China (Grant No. 10702066), Natural Science Foundation of Zhejiang Province (Grant No. Y7080398).

Biography: WANG Hao-li (1972-), Male, Ph. D., Associate Professor