- Email: [email protected]
Lattice Boltzmann simulation of gas flow over micro-scale airfoils
Computers & Fluids 38 (2009) 1675–1681
Contents lists available at ScienceDirect
Computers & Fluids j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m p fl u i d
Lattice Boltzmann simulation of gas flow over micro-scale airfoils Xiao-Dong Niu a,*, Shiaki Hyodo b, Kazuhiko Suga c, Hiroshi Yamaguchi a a b c
Energy Conversion Research Center, Department of Mechanical Engineering, Doshisha University, Kyotanabe, Kyoto 610-0321, Japan Computational Physics Lab., Toyota Central R & D Labs., Inc., Nagakute, Aichi 480-1192, Japan Department of Mechanical Engineering, Osaka Prefecture University, Sakai 599-8531, Japan
a r t i c l e
i n f o
Article history: Received 8 April 2008 Received in revised form 30 January 2009 Accepted 3 February 2009 Available online 27 February 2009
a b s t r a c t In order to understand aerodynamic issues related to design and performance of micro air aircraft, in this paper isothermal low Reynolds number flows over micro-scale airfoils are simulated using the kinetic lattice Boltzmann method [Niu et al., Phys. Rev. E 76 (2007) 036711 ]. The sample of the micro-scale airfoil investigated is a flat plate with a 5% thickness ratio. Investigation shows that low Reynolds number flows over the micro-scale airfoils are viscous and compressible, and that rarefied effects in these kinds of flows are dominant. It is also found that the lift coefficients of the micro-scale plate airfoils are always smaller than the drag coefficients of them at Reynolds numbers less than 100, and this observation is consistent with the previous studies. Ó 2009 Elsevier Ltd. All rights reserved.
1. Introduction Micro-scale aircraft have an increasing interest in special military and civil missions [1]. In recent years, many such aircraft are under development, including micro-sized unmanned aerial vehicles (UAVs) and micro air vehicles (lAVs) [2]. Due to small physical dimensions (i.e., 15 cm wingspan), the micro-scale aircraft have quite different aerodynamic performances from the manned macroscopic aircraft. For an example, steady state aerodynamics of airplanes works poorly for insect flight because, when insect wings are tested in a wind tunnel, the measured forces are substantially smaller than those required for active flight [3]. The flows over the micro-scale aircraft are gas rarefied because an important rarefied parameter, the Knudsen number ðKn ¼ k=L, where k is the molecular mean free path and L the characteristic length of the flow system), usually has a large value. It is also suggested that these kinds of flows are dominated by viscous and compressible effects, and usually at low Reynolds numbers (Re ¼ q1 U 1 L=l < 1000, where q1 ; U 1 ; l are density, velocity and viscosity of free-stream flows, respectively). Currently lots of researches and developments of the micro-scale aircraft have been conducted around the world [2–5]. However, foundamental aerodynamic studies are still in their infancy stage. There are very limited numerical and experimental studies for aerodynamics of the micro-scale aircraft [6–11]. Experimental studies of the rarefied gas flows over the micro-scale aircraft are made inherently difficult by small scales and high manufacture costs of aircraft, hence there is a strong trend to employ numerical tools to study aerodynamics of the micro-scale aircraft. * Corresponding author. E-mail address: [email protected] (X.-D. Niu). 0045-7930/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.compfluid.2009.02.003
It is well known that processes in the rarefied gas flows are described by the Boltzmann equation (BE) based on the kinetic theory [12,13]. Unfortunately, numerical solutions of the BE, either via simplified methods such as hard sphere molecule assumption [14] or via the direct simulation Monte Carlo (DSMC) approach [15,16], are difficult and a little time expensive, and numerical schemes based on the Chapman–Enskog methods [12,13,17] and the Grad’s moment approaches [18,19] have been limited to stable and boundary problems. As an alternative approach of the BE, the lattice Boltzmann method (LBM) [20,21] has been attempted to simulate the rarefied gas flows in the past few years [22–36]. Actually, the LBM is a simplified solver of the BE on a discrete lattice [37,38], and thus it losses some abilities of the BE in simulating gas fluids with compressible and rarefied effects [18,19]. Few endeavors on the LBM [39] have been made to overcome these limitations but successes are gianed by paying more complixity to the LBM scheme. Since the hydrodynamic moments of the LBE at various orders can be precisely and explicitly determined at a given order of truncations of the Hermite polynomials [37,38,40,41], adding more high-order terms to the equilibrium function of the LBM has become a more straightforward ways to considering the compressible and rarefied effects within 3 years [40–44]. Our earlire kinetic LBM model just employed this truction way to modeling the higherorder hydrodynamic moments and has been shown to be capable of capturing the Knudsen boundary layers [43]. In this paper, the goal is to test capabilities of the kinetic LBM [43] modeling the external rarefied gas flows over the micro-scale aircraft. In the LBM model [43], the distribution functions are expanded in terms of the 3rd order Hermite polynomials and the non-equilibrium distributions are regularized in the defined Hermite orthogonal space [28,32,44] to ensure the corresponding higher-order hydrodynamic moments correctly captured. In this
1676
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
paper, isothermal rarefied gas flows over a two-dimensional (2D) micro-scale flat plate with 5% thickness ratio for the Reynolds numbers in the range of 1 < Re < 200 and the free-stream Mach number (Ma ¼ U 1 =cs , where cs is the sound speed) at far field equal to 0.2 are studied, and numerical results obtained by the present LBM simulation are compared to the previous experimental, analytical and numerical results [9–11,45–47]. The rest of the paper is organized as follows. Section 2 briefly describes the kinetic lattice Boltzmann model we recently developed [43] for modeling the rarefied gas flows. Section 3 presents numerical simulations of the chosen flows using the introduced LBM model. Comparisons with the previous numerical, experimental and theoretical studies are also included in this section. Conclusions are given in Section 4.
2. Kinetic lattice Boltzmann model In the kinetic lattice Boltzmann model [42], flow systems are described by the kinetic LBE, in which a single particle distribution function fa at a set of particle discrete velocities, fna ; a ¼ 0; 1; 2; ; dg, is used as the state variable. The velocity space discretization is shown to be equivalent to project the distribution function onto a sub-space spanned by the leading N Hermite orthogonal basis, provided that na is the abscissa of a sufficiently accurate Gauss– Hermite quadrature with weight xa [37,38,40]. By using the BGK approximation [37,48] and regularization procedure [28,32,44] for the particle collision, the kinetic LBE can be written as
1 0 fa ðx; tÞ; fa ðx þ na dt; t þ dtÞ ¼ fað0Þ ðx; tÞ þ 1
s
ð1Þ
with the regularization non-equilibrium function
fa0 ¼ xa
N X n¼0
d X 1 HðnÞ ðna Þ ðfa fað0Þ ÞHðnÞ ðna Þ; n!c2n s a¼0
(
ð0Þ
fa
" # 2 na u 1 na u u2 ðx; tÞ ¼ xa q 1 þ þ RT 2 RT RT ) " # 2 1 na u u 2 na u : þ 3 6 RT RT RT
Correspondingly, the regularized non-equilibrium distribution function in the 3rd order gives
" fa0 ¼ xa
d X 1 ð2Þ H ðna Þ ðfa fað0Þ Þnai naj 4 2cs a¼0
# d X 1 ð3Þ ð0Þ þ 6 H ðna Þ ðfa fa Þnai naj nak ; 6cs a¼0
ð3Þ
where Hð2Þ ðna Þ ¼ nai naj c2s dij ; Hð3Þ ðna Þ ¼ nai naj nak c2s ðnai djk þ naj dik þnak dij Þ, and dij is the Kronecker delta function. The set of discrete velocities na can be obtained from the roots of Hermite polynomials through the optimal Gauss–Hermite quadratures. Except for a few special cases, the optimal lattice velocities usually do not coincide with the normal Cartesian coordinates. However, one can construct the quadratures on predefined abscissae by guaranteeing isotropy for tensors up to the higherorder and obtain the corresponding discrete velocities as long as a sufficient number of abscissae used, even though the number of abscissae is greater than that in the optimal Gauss quadratures. By using the above idea, a D3Q39 discrete velocity model can be obtained for the 3D Cartesian coordinates [32] and its 2D projection generates the D2Q21 model (Table 1). It is easily verified both P xa nai nan . models satisfying the six-order isotropy of tensors By using the above models, the implementation of the present kinetic LBM is eased. In the discrete velocity space, the fluid density q, velocity u are calculated as
q¼
d X
fa ; qu ¼
a¼0
where x; t and dt are the spatial coordinate vector, time and time interval, respectively; s is the relaxation parameter, HðnÞ is the nth Hermite polynomial [40]. Different from the earlier LBM models with the standard BGK collision for modeling micro-scale gas flows [22–27,29], the present model with regularized collision generates more physics. Mathematically, the regularization procedure serves as a filter and ensures the non-equilibrium distributions remain inside the defined Hermite orthogonal space by filtering out the higher-order non-equilibrium moments not supported by the defined Hermite basis. Besides improvement of stability and isotropy of the LBE (1) [28,32,44], the Knudsen layer, a region that the non-Newtonian stress/strain-rate relationship exists, can be given correctly to some extent [32–34,41–44]. This is because the nonlinear fluid behaviors in the Knudsen layer are mostly influenced by the higher-order non-equilibrium moments, and the regularization procedure ensures consistent descriptions of the higher-order hydrodynamic moments in the LBM frame. ð0Þ The equilibrium distribution function fa in Eq. (1) is the truncated Hermite expansion of the Maxwellian equilibrium [12,13] at a set of discrete velocities na . The truncation level of the Hermite ð0Þ expansion of fa determines accuracy of the LBE to approximate the BE. As revealed by the Chapman–Enskog method [12,13,16], by retaining up to the fourth order terms in the Hermite expansion, the Burnett level accuracy pertaining to the fluid momentum evolution for isothermal systems can be satisfied. When the Mach number is small, a third order Hermite expansion is said to be enough to model the momentum equations at Burnett level. Therefore, a 3rd order equilibrium distribution function is employed in present model and it is given as [32]
ð2Þ
d X
fa na :
ð4Þ
a¼0
The relaxation parameter s in the kinetic LBE (1) is determined in terms of the Knudsen number and can be given as [43]
s¼
rffiffiffiffi 2
p
Kncs
L þ 0:5; dt
ð5Þ
pffiffiffi where Kn ¼ p2 Ma . Re In the rarefied gas flows, when solid wall boundaries are presented, the slip boundary condition should be imposed. The kinetic theory states that the slip boundary condition in solving the BE can be described by the Maxwell’s diffuse boundary condition [12,13]. By using the Gauss–Hermite quadrature [32,37,38,40], this boundary condition straightforwardly yields the following diffuse-scattering boundary condition for the present LBM [23–26,43]
P 0 0 jðn uw Þ njfa0 ðx; tÞ fað0Þ fa ðx; tÞ ¼ P a 0 a ;w ðx; tÞ; ð0Þ a0 jðna uw Þ njfa0 ðx; tÞ
ð6Þ
whereðn0a uw Þ n < 0; ðna uw Þ n > 0, with n being the unit wall vector normal to the boundary surface and the subscription w meaning the wall. Previous numerical and theoretical analysis [33,35,36,44] has shown that the above equation on the higherorder ðP 2Þ is capable of capturing the Knudsen boundary layers.
3. Numerical simulations of rarefied gas flows over a microscale airfoil In this section, the isothermal rarefied gas flows over a microscale airfoil at low Reynolds numbers (Fig. 1) are simulated by using the above kinetic LBM. The airfoil is a flat plate having a
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681 Table 1 D2Q21 and D3Q39 discrete velocity models and their corresponding weighting pffiffiffiffiffiffi pffiffiffiffiffiffiffiffi functions. (For these two models, the sound speed cs ¼ RT ¼ 2=3). Model
na
xa
a
D2Q21
(0, 0) (±1, 0),(0, ±1) (±1, ±1) (±2, 0),(0, ±2) (±2, ±2) (±3, 0),(0, ± 3)
91/324 1/12 2/27 7/360 1/432 1/1620
0 1–4 5–8 9–12 13–16 17–20
D3Q39
(0, 0, 0) (±1, 0, 0),(0, ±1, 0),(0.0, ±1) (±1, ±1, ±1) (±2, 0, 0),(0, ±2, 0),(0, 0, ±2) (±2, ±2, ±2) (±3, 0, 0),(0, ±3, 0),(0, 0, ±3)
1/12 1/12 1/27 2/135 1/432 1/1620
0 1–6 7–14 15–20 21–32 33–38
1677
examine its accuracy by comparing with other theoretical, experimental and numerical works [9–11,45–47]. Fig. 2 shows the global drag coefficient C d of the micro-scale plate airfoil of zero-attack-angle flows at low Reynolds numbers from the present simulations and several other techniques, which include theoretical predictions based on the free molecular theory [45] and the Oseen–Stokes equations [46], and experimental studies [47]. The free molecular theory predicts that C d ¼ 1:35=Ma when Re < 0.1 and the Oseen–Stokes solution gives the flat plate drag coefficient C d in the slip flow (1 < Re <10) as
Cd ¼
8p 4Kn½lnð2=KnÞ þ c þ 1 ; 1 Reflnð16=ReÞ c þ 1g p½lnð16=ReÞ c þ 1
ð8Þ
The present LBM model has been numerically validated in our previous work [43] by studying shear-driven and micro-channel gas flows and comparing results with the DSMC simulations. In this work of studying flows over the micro-scale airfoils, we further
where c ¼ 0:57722. As shown in Fig. 2, the kinetic LBM gives the similar results of the slip Oseen–Stokes solution for the Reynolds number in the range of 1 < Re < 10. When Re < 0.2, the present results approach the free molecular prediction of C d ¼ 6:75, and agree well with experimental data of Schaaf and Sherman when Re > 5. In order to verify local accuracy of the kinetic LBM in the present study, we compare slip velocities, velocity profiles and pressure coefficients on the surface of the investigated micro-scale airfoils with numerical simulations of Sun and Boyd [10]. Sun and Boyd used a hybrid continuum-particle approach coupling the information preservation method [11] in the flow regions near the airfoils with the Navier–Stokes solver in the other regions. The information preservation method has been shown to be a valid particle method for modeling rarefied gas flows [9–11]. Fig. 3 shows the comparisons of slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio at the zero angle of attack and two different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10]. The comparison of the pressure coefficient distributions, C p ¼ ðp p1 Þ= ðq1 U 21 =2Þ, on the lower side of the micro-scale flat plate airfoil at the angle of attack of 20° and different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10] is displayed in Fig. 4. As shown in Figs. 3 and 4, good agreement is observed between the present simulations and those of Sun and Boyd [10], especially in the Reynolds numbers larger than 10. The differences of the present results and those of Sun and Boyd in the leading and trailing edges can be attributed to the following three reasons. First, the Burnett level high-order moments significant in rarefied gas flows are considered in the whole computational domain in the present LBM yet are only simulated in the region around the leading edge of the airfoil in simulations of Sun and Boyd [10]. Second, the diffuse-scattering boundary condition used in the present LBM is more powerful in capturing the Knudsen
Fig. 1. Sketch of computational domain of flows over the micro-scale flat plate airfoil of 5% thickness ratio.
Fig. 2. Drag coefficient of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at low Reynolds numbers from several techniques.
length L and a thickness ratio of 5%. The flow at the far field has a Mach number of 0.2. The computational domain is a 7L 4L rectangular domain with the center of the plate set in ð2L; 2LÞ. The stress integration technique [49] is employed to calculate forces acting on the plate airfoil. The diffuse-scattering boundary condition Eq. (6) is imposed on the plate surface. At the outside boundaries, a far field free-stream boundary condition is imposed and given by
fa joutside boundaries ¼ fað0Þ joutside boundaries ðq1 ; U 1 Þ:
ð7Þ
It should be noted that the open flow simulations are very sensitive to the far filed boundary condition [21,50], which should be provided in a way that does not affect the bulk flow strongly. In present simulations, the same strategy of Sun and Boyd [10] taking a domain long enough in the direction of flow (up to 7L) is adopted. All results presented below have been normalized by the characteristic length L, the free-stream density q1 and velocity U 1 . The grid independence study have been carried out first and showed that the simulation results are consistent on the uniform grids 700 400 and 1350 600. For briefness of paper we here omit the related discussing and all following results are based on the 700 400 uniform grids. 3.1. Numerical validations
1678
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
Fig. 3. Comparison of slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio at the zero angle of attack and two different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10].
boundary layer [33,35,36,44] than the first-order slip boundary condition ðU s ¼ k@u=@nÞ used by Sun and Boyd [10]. Third, Sun and Boyd considered thermal effects in their simulations while we only focus on the isothermal flows in the present work. The above three effects are coupled and the thermal effects are significant to the differences, thus in Figs. 3 and 4 it is difficult to claim accuracy of the present LBM to the approach of Sun and Boyd [10]. However, since the results obtained by the present LBM are quite close to those of Sun and Boyd, we can generally accept validation of the present LBM. 3.2. Zero-attack-angle flows over the micro-scale flat plate airfoil— viscous, compressible and rarefied effects The Reynolds number denotes the ratio of the inertia force to the viscous force. Therefore, the viscous effects are important in
Fig. 6. Pressure coefficients on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
small Reynolds number flows. The viscous effects can be observed in Figs. 3 and 4 and more demonstrated in Fig. 5, which plots velocity contours of juj ¼ 0:8 around the micro-scale plate airfoil at three Reynolds numbers 1.357, 13.57 and 135.7. Clearly seen from Fig. 5, when Re is small, the flat plate slows the flow because of the gas viscosity, and this effect increases as Re decreases. Generally, the flow is assumed to be incompressible when the Mach number is below 0.3. However, when the Reynolds number is small, a flow over the micro-sclale aircraft may become compressible even when the Mach number is small. Fig. 6 shows the pressure coefficient distribution, Cp, on the lower side of the micro-scale plate airfoil. It is found that the pressure coefficients for the cases of Re = 13.57 and 135.7 are very close, and deviate largely from that of Re = 1.357, which implies that the compressible effects are significant when the Reynolds number is below 10. Rarefied effects can also be observed from Fig. 6 and more clearly from Fig. 7. Fig. 7 shows the slip velocities on the upper side of the micro-scale plate airfoil for flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7, which correspond the Knudsen numbers of Kn = 0.28, 0.028 and 0.0028, respectively. As shown in Fig. 7, slip velocities exist for all three cases and the large Kn the large slip presents on the plate. Although the Knudsen numbers is small when the Reynolds number is larger than 10, there is a considerable amount of velocity slip near the leading and trailing edges of the micro-sclale plate airfoil. These observations indicate
Fig. 4. Comparison of pressure coefficient distributions on the lower side of the micro-scale flat plate airfoil of 5% thickness ratio at the angle of attack of 20° and three different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10].
Fig. 5. Velocity contours of juj ¼ 0:8 around the micro-scale plate airfoil for three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
Fig. 7. Slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
1679
Fig. 8. Pressure contours and velocity vectors of flows over a micro-scale flat plate airfoil of 5% thickness ratio at different Reynolds numbers when the attack-angles h ¼ 50 .
that, the rarefied effects are more dominant in the regions near the leading and trailing edges when Re > 10, and expand to larger domains when Re < 10. On the other hand, some disturbances of the pressure coefficients and slip velocities are observed near the leading and trailing edges of the micro-scale plate airfoil at Re = .357 and 13.57 in Figs. 6 and 7, and we attribute these disturbances to more demonstrations of the viscous, compressible and rarefied effects in the regions near the leading and trailing edges of the micro-scale plate airfoil. 3.3. Flows over a micro-scale flat plate airfoil at various attack-angles
Fig. 9. Pressure coefficient distributions on the lower and upper sides of the microscale plate for the different flow Reynolds numbers at attack-angles of h ¼ 50 .
Fig. 8 show the pressure contours and the velocity vectors of the flows over the micro-scale flat plate airfoil at different attack-angles of h ¼ 50 when Re = 1.357, 13.57 and 135.7, respectively. As shown in this figure, there is no flow separation at Re = 1.357 even the attack-angle is very high. However, at Re P 10, the flow has shown to separate in the region near the upper leading edge because of strong local adverse pressure gradient. By comparing flow pattern at different Reynolds numbers in Fig. 8, the flow patterns at Re = 1.357 and 13.57 are shown to demonstrate different flow behaviors from those of Re = 135.7. The flow separation is delayed and weaken at Re = 13.57, and there is even no flow separation at the attack-angle of 50° when Re = 1.357. The pressure gradient near the leading edge decreases when the flow Reynolds number
Fig. 10. Lift and drag coefficients of the flows over the micro-scale flat plate airfoil with 5% thickness ratio at different attack angles and Reynolds numbers.
1680
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
decreases. This founding can be observed clearly from Fig. 9, which shows the pressure coefficient distributions on the lower and upper sides of the plate for the different flow Reynolds numbers at attack-angles of h ¼ 50 . Shown in Fig. 9, the pressure increases where the flow faces the plate and drops where the flow leaves the plate. From Figs. 8 and 9, one can anticipate that the aerodynamic performance of the micro-scale airfoil at low Reynolds number is different from that of high Reynolds numbers. 3.4. Aerodynamic performances of a micro-scale flat plate airfoil at low Reynolds numbers The aerodynamic performances of the airfoils are usually measured by lift and drag coefficients of C l and C d . Fig. 10 shows the lift and drag coefficients of the flows over the micro-scale flat plate airfoil with 5% thickness ratio at the different attack-angles and Reynolds numbers. It is found that the drag coefficient C d decreases when the Reynolds number increases, and the lift coefficient changes little when the Reynolds number varies. Increasing of the attack angle also increases the drag coefficient. However, the lift coefficient shows a different trend via the attack-angles at the Reynolds number below and above 100. When Re < 100, the lift coefficient shows a similar monotonically increase trend to the drag coefficient when the attack angle increases. When Re > 100, the lift coefficient drops when h > 40 after first increasing at small attack-angles. Moreover, the lift coefficients are always less than the drag coefficients when the Reynolds numbers of Re = 1.357 and 13.57, which suggests that the aerodynamic performance of the micro-scale flat plate with 5% thickness ratio is poor when the flow Reynolds number is small.
4. Conclusions We have shown that the recently develope LBM [43] can quantitatively predict the fluid behaviors and aerodynamics for flows over a micro-scale flat plate airfoil with 5% thickness ratio at low Reynolds numbers below 200. The low Reynolds number flows over the micro-scale flat plate airfoils are very different from the high Reynolds number ones due to the dominant effects of viscosity, compressibility and rarefaction. The viscous effects slows flows and the flows are compressible even the Mach number is small. The rarefied effects on the flows are more dominant in the regions near the leading and trailing edges at small reynolds numbers. The drag and lift coefficients of the micro-scale flat plate also display different trend for the high and low Reynolds numbers. The lift coefficients are always less than the drag coefficients when the Re < 100, suggesting that the aerodynamic performance of the micro-scale flat plate with 5% thickness ratio is poor when the flow Reynolds number is small. Acknowledgements This work were supported by Core Research for Evolutional Science and Technology (CREST) of Japan Science and Technology (JST) Agency (No.: 228205R) and the Japan Society for the Promotion of Science through a Grant-in-Aid for Scientific Research (B) (No. 18360050).
References [1] Mueller TJ, Fixed and flapping wing aerodynamics for micro air vehicle applications. In: Mueller T, editor. AIAA Progress in Astronautics and Aeronautics, vol. 195, 2001. p. 307–39.
[2] McMichael JM, Francis MS, Micro air vehicles—toward a new dimension in flight, DARPA, Available from: http://www.darpa.mil/tto/MAV/mav_auvsi. html, 1997. [3] Dickinson MH, Lehmann FO, Sane SP. Wing rotation and the aerodynamic basis of insect flight. Science 1999;284:1954–60. [4] Johnstone R, Arntz N. CONDOR—high altitude long endurance (HALE) automatically piloted vehicle (APV). AIAA 1990:1990–3279. [5] Wang J. Dissecting insect flight. Ann Rev Fluid Mech 2005;37:183–210. [6] Null W, Shkarayev S. Effect of camber on the aerodynamics of adaptive wing micro air vehicles. In: AIAA 2004-2694, the second AIAA flow control conference. Portland, Oregon, June 29–July 1; 2004. [7] Jones KD, Bradshaw CJ, Papadopoulos J, Platzer MF. Improved performance and control of flapping-wing propelled micro air vehicles. In: AIAA 2004-0399, the AIAA 42nd aerospace sciences meeting and exhibit. Reno, Nevada, January 5–8; 2004. [8] Sunada S, Sakaguchi A, Kawachi K. Airfoil section characteristics at a low Reynolds numbers. Trans ASME J Fluids Eng 1997;119:129–35. [9] Sun Q, Boyd ID. Flat-plate aerodynamics at very low Reynolds number. J Fluid Mech 2004;502:199–206. [10] Sun Q, Boyd ID. Simulation of gas flow over micro-scale airfoils using a hybrid continuum particle approach. AIAA 2003:2003–4014. [11] Sun Q, Boyd ID, Candler GV. Numerical simulation of gas flow over micro-scale airfoils. AIAA 2001:2001–3071. [12] Cercignani C. Theory and application of the Boltzmann equation. Edinburg: Sottish Academic Press; 1975. [13] Chapman S, Cowling TG. The mathematical theory of non-uniform gases. Cambridge University Press; 1970. [14] Ohwada T. Heat flow and temperature and density distributions in a rarefied gas between parallel plates with different temperatures: finite difference analysis of the nonlinear Boltzmann equation for hard sphere molecules. Phys Fluids 1996;8:2153–60. [15] Bird G. Molecular gas dynamics and the direct simulation of gas flows. Oxford: Clarendon Press; 1994. [16] Wang M, Li Z. Simulations for gas flows in microgeometries using the direct simulation Monte Carlo method. Int J Heat Fluid Flow 2004;25(6):975–85. [17] Wang M, Li Z. An Enskog based Monte Carlo method for high Knudsen number non-ideal gas flows. Comp & Fluids 2007;36(8):1291–7. [18] Müller I, Ruggeri T. Rational extended thermodynamics. New York: Springer; 1998. [19] Struchtrup H, Torrilhon M. Regularization of Grad’s 13 moment equations: derivation and linear analysis. Phys Fluids 2003;15(9):2668–80. [20] Wolf-Gladrow DA. Lattice-gas cellular automata and lattice Boltzmann models. Berlin: Springer; 2000. [21] Succi S. The lattice Boltzmann equation: for fluid dynamics and beyond (series numerical mathematics and scientific computation). Oxford New York: Oxford University Press; 2001. [22] Nie X, Doolen G, Chen S. Lattice Boltzmann simulations of fluid flows in MEMS. J Stat Phys 2002;107:279–89. [23] Ansumali S, Karlin IV. Kinetic boundary conditions in the lattice Boltzmann method. Phys Rev E 2002;66:026311. [24] Lim CY, Shu C, Niu XD, Chew YT. Application of lattice Boltzmann method to simulate microchannel flows. Phys Fluids 2002;107:2299–308. [25] Niu XD, Shu C, Chew YT. Lattice Boltzmann BGK model for simulation of micro flows. Eur Phys Lett 2004;67:600–6. [26] Shu C, Niu XD, Chew YT. A lattice Boltzmann kinetic model for microflow and heat transfer. J Stat Phys 2005;121:239–55. [27] Niu XD, Shu C, Chew YT. Numerical simulation of isothermal micro flows by lattice Boltzmann method and theoretical analysis of the diffusive scattering boundary condition. Int J Mod Phys C 2005;16:1927–41. [28] Zhang YH, Gu XJ, Barber RW, Emerson DR. Capturing Knudsen layer phenomena using a lattice Boltzmann model. Phys Rev E 2006;74:046704. [29] Chen H, Zhang R, Staroselsky I, Jhon M. Recovery of full rotational invariance in lattice Boltzmann formulations for high Knudsen number flows. Phys A 2006;362:125–31. [30] Zhou Y, Zhang R, Staroselsky I, Chen H, Kim WT, Jhon M. Simulation of micro and nano-scale flows via the lattice Boltzmann formulations for high Knudsen number flows. Phys A 2006;362:68–77. [31] Toshi F, Succi S. Lattice Boltzmann method at finite Knudsen numbers. Eur Phys Lett 2005;69:549–55. [32] Shan X, Yuam XF, Chen H. Kinetic theory representation of hydrodynamics: a way beyond the Navier–Stokes equation. J Fluid Mech 2006;550:413–41. [33] Ansumali S, Karlin IV, Arcidiacono S, Abbas A, Prasianakis N. Hydrodynamics beyond Navier–Stokes: exact solution to the lattice Boltzmann hierarchy. Phys Rev Lett 2007;98:124502. [34] Yudistiawan WP, Ansumali S, Karlin IV. Hydrodynamics beyond Navier– Stokes: the slip flow model. Phys Rev E 2008;78:016705. [35] Kim SH, Pitsch H, Boyd LD. Slip velocity and Knudsen layer in the lattice Boltzmann method for microscale flows. Phys Rev E 2008;77:026704. [36] Kim SH, Pitsch H, Boyd LD. Accuracy of higher-order lattice Boltzmann methods for microscale flows with finite Knudsen numbers. J Comp Phys 2008;227:8655–71. [37] He X, Luo L-S. Theory of lattice Boltzmann method: from the Boltzmann equation to the lattice Boltzmann equation. Phys Rev E 1997;56(6):6811–7. [38] Shan X, He X. Discretization of the velocity space in the solution of the Boltzmann equation. Phys Rev Lett 1998;80:65–8.
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681 [39] Alexander FJ, Chen H, Chen S, Doolen GD. Lattice Boltzmann model for compressible fluids. Phys Rev A 1992;46:1967–70. [40] Zhang R, Shan X, Chen H. Efficient kinetic method for fluid simulation beyond the Navier–Stokes equation. Phys Rev E 2006;74:046703. [41] Chikatamarla S, Ansumali S, Karlin I. Entropic lattice Boltzmann models for hydrodynamics in three dimensions. Phys Rev Lett 2006;97(1):010201. [42] Ansumali S, Karlin IV, Oettinger HC. Minimal entropic kinetic models for hydrodynamics. Eur Phys Lett 2003;63:798–804. [43] Niu XD, Hyodo S, Munekata T, Suga K. Kinetic lattice Boltzmann method for micro gas flows: issues on boundary condition, relaxation time and regularization. Phys Rev E 2007;76:036711. [44] Latt J, Chopard B. Lattice Boltzmann method with regularized pre-collision distribution functions. Math Comp Simulation 2006;72:165–8.
1681
[45] Gombosi TI. Gaskinetics theory. Cambridge: Cambridge University Press; 1994. [46] Tamada K, Miura H. Slip flow past a tangential flat plate at low Reynolds numbers. J Fluid Mech 1978;85:731–42. [47] Schaaf SA, Sherman FS. Skin friction in slip flow. J Aeronaut Sci 1954;21:85–90. [48] Bhatnagar PL, Gross EP, Krook M. A model for collision processes in gases I: small amplitude processes in charged and neutral one-component systems. Phys Rev 1954;94:511–25. [49] Inamuro T, Maeba K, Ogino F. Flow between parallel walls containing the lines of neutrally buoyant circular cylinders. Int J Multiphase Flow 2000;26: 1981–2004. [50] Davis RW, Moore EF. A numerical study of vortex shedding from rectangles. J Fluid Mech 1982;116:475.
Contents lists available at ScienceDirect
Computers & Fluids j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m p fl u i d
Lattice Boltzmann simulation of gas flow over micro-scale airfoils Xiao-Dong Niu a,*, Shiaki Hyodo b, Kazuhiko Suga c, Hiroshi Yamaguchi a a b c
Energy Conversion Research Center, Department of Mechanical Engineering, Doshisha University, Kyotanabe, Kyoto 610-0321, Japan Computational Physics Lab., Toyota Central R & D Labs., Inc., Nagakute, Aichi 480-1192, Japan Department of Mechanical Engineering, Osaka Prefecture University, Sakai 599-8531, Japan
a r t i c l e
i n f o
Article history: Received 8 April 2008 Received in revised form 30 January 2009 Accepted 3 February 2009 Available online 27 February 2009
a b s t r a c t In order to understand aerodynamic issues related to design and performance of micro air aircraft, in this paper isothermal low Reynolds number flows over micro-scale airfoils are simulated using the kinetic lattice Boltzmann method [Niu et al., Phys. Rev. E 76 (2007) 036711 ]. The sample of the micro-scale airfoil investigated is a flat plate with a 5% thickness ratio. Investigation shows that low Reynolds number flows over the micro-scale airfoils are viscous and compressible, and that rarefied effects in these kinds of flows are dominant. It is also found that the lift coefficients of the micro-scale plate airfoils are always smaller than the drag coefficients of them at Reynolds numbers less than 100, and this observation is consistent with the previous studies. Ó 2009 Elsevier Ltd. All rights reserved.
1. Introduction Micro-scale aircraft have an increasing interest in special military and civil missions [1]. In recent years, many such aircraft are under development, including micro-sized unmanned aerial vehicles (UAVs) and micro air vehicles (lAVs) [2]. Due to small physical dimensions (i.e., 15 cm wingspan), the micro-scale aircraft have quite different aerodynamic performances from the manned macroscopic aircraft. For an example, steady state aerodynamics of airplanes works poorly for insect flight because, when insect wings are tested in a wind tunnel, the measured forces are substantially smaller than those required for active flight [3]. The flows over the micro-scale aircraft are gas rarefied because an important rarefied parameter, the Knudsen number ðKn ¼ k=L, where k is the molecular mean free path and L the characteristic length of the flow system), usually has a large value. It is also suggested that these kinds of flows are dominated by viscous and compressible effects, and usually at low Reynolds numbers (Re ¼ q1 U 1 L=l < 1000, where q1 ; U 1 ; l are density, velocity and viscosity of free-stream flows, respectively). Currently lots of researches and developments of the micro-scale aircraft have been conducted around the world [2–5]. However, foundamental aerodynamic studies are still in their infancy stage. There are very limited numerical and experimental studies for aerodynamics of the micro-scale aircraft [6–11]. Experimental studies of the rarefied gas flows over the micro-scale aircraft are made inherently difficult by small scales and high manufacture costs of aircraft, hence there is a strong trend to employ numerical tools to study aerodynamics of the micro-scale aircraft. * Corresponding author. E-mail address: [email protected] (X.-D. Niu). 0045-7930/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.compfluid.2009.02.003
It is well known that processes in the rarefied gas flows are described by the Boltzmann equation (BE) based on the kinetic theory [12,13]. Unfortunately, numerical solutions of the BE, either via simplified methods such as hard sphere molecule assumption [14] or via the direct simulation Monte Carlo (DSMC) approach [15,16], are difficult and a little time expensive, and numerical schemes based on the Chapman–Enskog methods [12,13,17] and the Grad’s moment approaches [18,19] have been limited to stable and boundary problems. As an alternative approach of the BE, the lattice Boltzmann method (LBM) [20,21] has been attempted to simulate the rarefied gas flows in the past few years [22–36]. Actually, the LBM is a simplified solver of the BE on a discrete lattice [37,38], and thus it losses some abilities of the BE in simulating gas fluids with compressible and rarefied effects [18,19]. Few endeavors on the LBM [39] have been made to overcome these limitations but successes are gianed by paying more complixity to the LBM scheme. Since the hydrodynamic moments of the LBE at various orders can be precisely and explicitly determined at a given order of truncations of the Hermite polynomials [37,38,40,41], adding more high-order terms to the equilibrium function of the LBM has become a more straightforward ways to considering the compressible and rarefied effects within 3 years [40–44]. Our earlire kinetic LBM model just employed this truction way to modeling the higherorder hydrodynamic moments and has been shown to be capable of capturing the Knudsen boundary layers [43]. In this paper, the goal is to test capabilities of the kinetic LBM [43] modeling the external rarefied gas flows over the micro-scale aircraft. In the LBM model [43], the distribution functions are expanded in terms of the 3rd order Hermite polynomials and the non-equilibrium distributions are regularized in the defined Hermite orthogonal space [28,32,44] to ensure the corresponding higher-order hydrodynamic moments correctly captured. In this
1676
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
paper, isothermal rarefied gas flows over a two-dimensional (2D) micro-scale flat plate with 5% thickness ratio for the Reynolds numbers in the range of 1 < Re < 200 and the free-stream Mach number (Ma ¼ U 1 =cs , where cs is the sound speed) at far field equal to 0.2 are studied, and numerical results obtained by the present LBM simulation are compared to the previous experimental, analytical and numerical results [9–11,45–47]. The rest of the paper is organized as follows. Section 2 briefly describes the kinetic lattice Boltzmann model we recently developed [43] for modeling the rarefied gas flows. Section 3 presents numerical simulations of the chosen flows using the introduced LBM model. Comparisons with the previous numerical, experimental and theoretical studies are also included in this section. Conclusions are given in Section 4.
2. Kinetic lattice Boltzmann model In the kinetic lattice Boltzmann model [42], flow systems are described by the kinetic LBE, in which a single particle distribution function fa at a set of particle discrete velocities, fna ; a ¼ 0; 1; 2; ; dg, is used as the state variable. The velocity space discretization is shown to be equivalent to project the distribution function onto a sub-space spanned by the leading N Hermite orthogonal basis, provided that na is the abscissa of a sufficiently accurate Gauss– Hermite quadrature with weight xa [37,38,40]. By using the BGK approximation [37,48] and regularization procedure [28,32,44] for the particle collision, the kinetic LBE can be written as
1 0 fa ðx; tÞ; fa ðx þ na dt; t þ dtÞ ¼ fað0Þ ðx; tÞ þ 1
s
ð1Þ
with the regularization non-equilibrium function
fa0 ¼ xa
N X n¼0
d X 1 HðnÞ ðna Þ ðfa fað0Þ ÞHðnÞ ðna Þ; n!c2n s a¼0
(
ð0Þ
fa
" # 2 na u 1 na u u2 ðx; tÞ ¼ xa q 1 þ þ RT 2 RT RT ) " # 2 1 na u u 2 na u : þ 3 6 RT RT RT
Correspondingly, the regularized non-equilibrium distribution function in the 3rd order gives
" fa0 ¼ xa
d X 1 ð2Þ H ðna Þ ðfa fað0Þ Þnai naj 4 2cs a¼0
# d X 1 ð3Þ ð0Þ þ 6 H ðna Þ ðfa fa Þnai naj nak ; 6cs a¼0
ð3Þ
where Hð2Þ ðna Þ ¼ nai naj c2s dij ; Hð3Þ ðna Þ ¼ nai naj nak c2s ðnai djk þ naj dik þnak dij Þ, and dij is the Kronecker delta function. The set of discrete velocities na can be obtained from the roots of Hermite polynomials through the optimal Gauss–Hermite quadratures. Except for a few special cases, the optimal lattice velocities usually do not coincide with the normal Cartesian coordinates. However, one can construct the quadratures on predefined abscissae by guaranteeing isotropy for tensors up to the higherorder and obtain the corresponding discrete velocities as long as a sufficient number of abscissae used, even though the number of abscissae is greater than that in the optimal Gauss quadratures. By using the above idea, a D3Q39 discrete velocity model can be obtained for the 3D Cartesian coordinates [32] and its 2D projection generates the D2Q21 model (Table 1). It is easily verified both P xa nai nan . models satisfying the six-order isotropy of tensors By using the above models, the implementation of the present kinetic LBM is eased. In the discrete velocity space, the fluid density q, velocity u are calculated as
q¼
d X
fa ; qu ¼
a¼0
where x; t and dt are the spatial coordinate vector, time and time interval, respectively; s is the relaxation parameter, HðnÞ is the nth Hermite polynomial [40]. Different from the earlier LBM models with the standard BGK collision for modeling micro-scale gas flows [22–27,29], the present model with regularized collision generates more physics. Mathematically, the regularization procedure serves as a filter and ensures the non-equilibrium distributions remain inside the defined Hermite orthogonal space by filtering out the higher-order non-equilibrium moments not supported by the defined Hermite basis. Besides improvement of stability and isotropy of the LBE (1) [28,32,44], the Knudsen layer, a region that the non-Newtonian stress/strain-rate relationship exists, can be given correctly to some extent [32–34,41–44]. This is because the nonlinear fluid behaviors in the Knudsen layer are mostly influenced by the higher-order non-equilibrium moments, and the regularization procedure ensures consistent descriptions of the higher-order hydrodynamic moments in the LBM frame. ð0Þ The equilibrium distribution function fa in Eq. (1) is the truncated Hermite expansion of the Maxwellian equilibrium [12,13] at a set of discrete velocities na . The truncation level of the Hermite ð0Þ expansion of fa determines accuracy of the LBE to approximate the BE. As revealed by the Chapman–Enskog method [12,13,16], by retaining up to the fourth order terms in the Hermite expansion, the Burnett level accuracy pertaining to the fluid momentum evolution for isothermal systems can be satisfied. When the Mach number is small, a third order Hermite expansion is said to be enough to model the momentum equations at Burnett level. Therefore, a 3rd order equilibrium distribution function is employed in present model and it is given as [32]
ð2Þ
d X
fa na :
ð4Þ
a¼0
The relaxation parameter s in the kinetic LBE (1) is determined in terms of the Knudsen number and can be given as [43]
s¼
rffiffiffiffi 2
p
Kncs
L þ 0:5; dt
ð5Þ
pffiffiffi where Kn ¼ p2 Ma . Re In the rarefied gas flows, when solid wall boundaries are presented, the slip boundary condition should be imposed. The kinetic theory states that the slip boundary condition in solving the BE can be described by the Maxwell’s diffuse boundary condition [12,13]. By using the Gauss–Hermite quadrature [32,37,38,40], this boundary condition straightforwardly yields the following diffuse-scattering boundary condition for the present LBM [23–26,43]
P 0 0 jðn uw Þ njfa0 ðx; tÞ fað0Þ fa ðx; tÞ ¼ P a 0 a ;w ðx; tÞ; ð0Þ a0 jðna uw Þ njfa0 ðx; tÞ
ð6Þ
whereðn0a uw Þ n < 0; ðna uw Þ n > 0, with n being the unit wall vector normal to the boundary surface and the subscription w meaning the wall. Previous numerical and theoretical analysis [33,35,36,44] has shown that the above equation on the higherorder ðP 2Þ is capable of capturing the Knudsen boundary layers.
3. Numerical simulations of rarefied gas flows over a microscale airfoil In this section, the isothermal rarefied gas flows over a microscale airfoil at low Reynolds numbers (Fig. 1) are simulated by using the above kinetic LBM. The airfoil is a flat plate having a
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681 Table 1 D2Q21 and D3Q39 discrete velocity models and their corresponding weighting pffiffiffiffiffiffi pffiffiffiffiffiffiffiffi functions. (For these two models, the sound speed cs ¼ RT ¼ 2=3). Model
na
xa
a
D2Q21
(0, 0) (±1, 0),(0, ±1) (±1, ±1) (±2, 0),(0, ±2) (±2, ±2) (±3, 0),(0, ± 3)
91/324 1/12 2/27 7/360 1/432 1/1620
0 1–4 5–8 9–12 13–16 17–20
D3Q39
(0, 0, 0) (±1, 0, 0),(0, ±1, 0),(0.0, ±1) (±1, ±1, ±1) (±2, 0, 0),(0, ±2, 0),(0, 0, ±2) (±2, ±2, ±2) (±3, 0, 0),(0, ±3, 0),(0, 0, ±3)
1/12 1/12 1/27 2/135 1/432 1/1620
0 1–6 7–14 15–20 21–32 33–38
1677
examine its accuracy by comparing with other theoretical, experimental and numerical works [9–11,45–47]. Fig. 2 shows the global drag coefficient C d of the micro-scale plate airfoil of zero-attack-angle flows at low Reynolds numbers from the present simulations and several other techniques, which include theoretical predictions based on the free molecular theory [45] and the Oseen–Stokes equations [46], and experimental studies [47]. The free molecular theory predicts that C d ¼ 1:35=Ma when Re < 0.1 and the Oseen–Stokes solution gives the flat plate drag coefficient C d in the slip flow (1 < Re <10) as
Cd ¼
8p 4Kn½lnð2=KnÞ þ c þ 1 ; 1 Reflnð16=ReÞ c þ 1g p½lnð16=ReÞ c þ 1
ð8Þ
The present LBM model has been numerically validated in our previous work [43] by studying shear-driven and micro-channel gas flows and comparing results with the DSMC simulations. In this work of studying flows over the micro-scale airfoils, we further
where c ¼ 0:57722. As shown in Fig. 2, the kinetic LBM gives the similar results of the slip Oseen–Stokes solution for the Reynolds number in the range of 1 < Re < 10. When Re < 0.2, the present results approach the free molecular prediction of C d ¼ 6:75, and agree well with experimental data of Schaaf and Sherman when Re > 5. In order to verify local accuracy of the kinetic LBM in the present study, we compare slip velocities, velocity profiles and pressure coefficients on the surface of the investigated micro-scale airfoils with numerical simulations of Sun and Boyd [10]. Sun and Boyd used a hybrid continuum-particle approach coupling the information preservation method [11] in the flow regions near the airfoils with the Navier–Stokes solver in the other regions. The information preservation method has been shown to be a valid particle method for modeling rarefied gas flows [9–11]. Fig. 3 shows the comparisons of slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio at the zero angle of attack and two different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10]. The comparison of the pressure coefficient distributions, C p ¼ ðp p1 Þ= ðq1 U 21 =2Þ, on the lower side of the micro-scale flat plate airfoil at the angle of attack of 20° and different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10] is displayed in Fig. 4. As shown in Figs. 3 and 4, good agreement is observed between the present simulations and those of Sun and Boyd [10], especially in the Reynolds numbers larger than 10. The differences of the present results and those of Sun and Boyd in the leading and trailing edges can be attributed to the following three reasons. First, the Burnett level high-order moments significant in rarefied gas flows are considered in the whole computational domain in the present LBM yet are only simulated in the region around the leading edge of the airfoil in simulations of Sun and Boyd [10]. Second, the diffuse-scattering boundary condition used in the present LBM is more powerful in capturing the Knudsen
Fig. 1. Sketch of computational domain of flows over the micro-scale flat plate airfoil of 5% thickness ratio.
Fig. 2. Drag coefficient of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at low Reynolds numbers from several techniques.
length L and a thickness ratio of 5%. The flow at the far field has a Mach number of 0.2. The computational domain is a 7L 4L rectangular domain with the center of the plate set in ð2L; 2LÞ. The stress integration technique [49] is employed to calculate forces acting on the plate airfoil. The diffuse-scattering boundary condition Eq. (6) is imposed on the plate surface. At the outside boundaries, a far field free-stream boundary condition is imposed and given by
fa joutside boundaries ¼ fað0Þ joutside boundaries ðq1 ; U 1 Þ:
ð7Þ
It should be noted that the open flow simulations are very sensitive to the far filed boundary condition [21,50], which should be provided in a way that does not affect the bulk flow strongly. In present simulations, the same strategy of Sun and Boyd [10] taking a domain long enough in the direction of flow (up to 7L) is adopted. All results presented below have been normalized by the characteristic length L, the free-stream density q1 and velocity U 1 . The grid independence study have been carried out first and showed that the simulation results are consistent on the uniform grids 700 400 and 1350 600. For briefness of paper we here omit the related discussing and all following results are based on the 700 400 uniform grids. 3.1. Numerical validations
1678
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
Fig. 3. Comparison of slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio at the zero angle of attack and two different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10].
boundary layer [33,35,36,44] than the first-order slip boundary condition ðU s ¼ k@u=@nÞ used by Sun and Boyd [10]. Third, Sun and Boyd considered thermal effects in their simulations while we only focus on the isothermal flows in the present work. The above three effects are coupled and the thermal effects are significant to the differences, thus in Figs. 3 and 4 it is difficult to claim accuracy of the present LBM to the approach of Sun and Boyd [10]. However, since the results obtained by the present LBM are quite close to those of Sun and Boyd, we can generally accept validation of the present LBM. 3.2. Zero-attack-angle flows over the micro-scale flat plate airfoil— viscous, compressible and rarefied effects The Reynolds number denotes the ratio of the inertia force to the viscous force. Therefore, the viscous effects are important in
Fig. 6. Pressure coefficients on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
small Reynolds number flows. The viscous effects can be observed in Figs. 3 and 4 and more demonstrated in Fig. 5, which plots velocity contours of juj ¼ 0:8 around the micro-scale plate airfoil at three Reynolds numbers 1.357, 13.57 and 135.7. Clearly seen from Fig. 5, when Re is small, the flat plate slows the flow because of the gas viscosity, and this effect increases as Re decreases. Generally, the flow is assumed to be incompressible when the Mach number is below 0.3. However, when the Reynolds number is small, a flow over the micro-sclale aircraft may become compressible even when the Mach number is small. Fig. 6 shows the pressure coefficient distribution, Cp, on the lower side of the micro-scale plate airfoil. It is found that the pressure coefficients for the cases of Re = 13.57 and 135.7 are very close, and deviate largely from that of Re = 1.357, which implies that the compressible effects are significant when the Reynolds number is below 10. Rarefied effects can also be observed from Fig. 6 and more clearly from Fig. 7. Fig. 7 shows the slip velocities on the upper side of the micro-scale plate airfoil for flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7, which correspond the Knudsen numbers of Kn = 0.28, 0.028 and 0.0028, respectively. As shown in Fig. 7, slip velocities exist for all three cases and the large Kn the large slip presents on the plate. Although the Knudsen numbers is small when the Reynolds number is larger than 10, there is a considerable amount of velocity slip near the leading and trailing edges of the micro-sclale plate airfoil. These observations indicate
Fig. 4. Comparison of pressure coefficient distributions on the lower side of the micro-scale flat plate airfoil of 5% thickness ratio at the angle of attack of 20° and three different Reynolds numbers obtained by the present LBM and approach of Sun and Boyd [10].
Fig. 5. Velocity contours of juj ¼ 0:8 around the micro-scale plate airfoil for three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
Fig. 7. Slip velocities on the upper side of the micro-scale flat plate airfoil of 5% thickness ratio for the zero-attack-angle flows at three Reynolds numbers of Re = 1.357, 13.57 and 135.7.
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
1679
Fig. 8. Pressure contours and velocity vectors of flows over a micro-scale flat plate airfoil of 5% thickness ratio at different Reynolds numbers when the attack-angles h ¼ 50 .
that, the rarefied effects are more dominant in the regions near the leading and trailing edges when Re > 10, and expand to larger domains when Re < 10. On the other hand, some disturbances of the pressure coefficients and slip velocities are observed near the leading and trailing edges of the micro-scale plate airfoil at Re = .357 and 13.57 in Figs. 6 and 7, and we attribute these disturbances to more demonstrations of the viscous, compressible and rarefied effects in the regions near the leading and trailing edges of the micro-scale plate airfoil. 3.3. Flows over a micro-scale flat plate airfoil at various attack-angles
Fig. 9. Pressure coefficient distributions on the lower and upper sides of the microscale plate for the different flow Reynolds numbers at attack-angles of h ¼ 50 .
Fig. 8 show the pressure contours and the velocity vectors of the flows over the micro-scale flat plate airfoil at different attack-angles of h ¼ 50 when Re = 1.357, 13.57 and 135.7, respectively. As shown in this figure, there is no flow separation at Re = 1.357 even the attack-angle is very high. However, at Re P 10, the flow has shown to separate in the region near the upper leading edge because of strong local adverse pressure gradient. By comparing flow pattern at different Reynolds numbers in Fig. 8, the flow patterns at Re = 1.357 and 13.57 are shown to demonstrate different flow behaviors from those of Re = 135.7. The flow separation is delayed and weaken at Re = 13.57, and there is even no flow separation at the attack-angle of 50° when Re = 1.357. The pressure gradient near the leading edge decreases when the flow Reynolds number
Fig. 10. Lift and drag coefficients of the flows over the micro-scale flat plate airfoil with 5% thickness ratio at different attack angles and Reynolds numbers.
1680
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681
decreases. This founding can be observed clearly from Fig. 9, which shows the pressure coefficient distributions on the lower and upper sides of the plate for the different flow Reynolds numbers at attack-angles of h ¼ 50 . Shown in Fig. 9, the pressure increases where the flow faces the plate and drops where the flow leaves the plate. From Figs. 8 and 9, one can anticipate that the aerodynamic performance of the micro-scale airfoil at low Reynolds number is different from that of high Reynolds numbers. 3.4. Aerodynamic performances of a micro-scale flat plate airfoil at low Reynolds numbers The aerodynamic performances of the airfoils are usually measured by lift and drag coefficients of C l and C d . Fig. 10 shows the lift and drag coefficients of the flows over the micro-scale flat plate airfoil with 5% thickness ratio at the different attack-angles and Reynolds numbers. It is found that the drag coefficient C d decreases when the Reynolds number increases, and the lift coefficient changes little when the Reynolds number varies. Increasing of the attack angle also increases the drag coefficient. However, the lift coefficient shows a different trend via the attack-angles at the Reynolds number below and above 100. When Re < 100, the lift coefficient shows a similar monotonically increase trend to the drag coefficient when the attack angle increases. When Re > 100, the lift coefficient drops when h > 40 after first increasing at small attack-angles. Moreover, the lift coefficients are always less than the drag coefficients when the Reynolds numbers of Re = 1.357 and 13.57, which suggests that the aerodynamic performance of the micro-scale flat plate with 5% thickness ratio is poor when the flow Reynolds number is small.
4. Conclusions We have shown that the recently develope LBM [43] can quantitatively predict the fluid behaviors and aerodynamics for flows over a micro-scale flat plate airfoil with 5% thickness ratio at low Reynolds numbers below 200. The low Reynolds number flows over the micro-scale flat plate airfoils are very different from the high Reynolds number ones due to the dominant effects of viscosity, compressibility and rarefaction. The viscous effects slows flows and the flows are compressible even the Mach number is small. The rarefied effects on the flows are more dominant in the regions near the leading and trailing edges at small reynolds numbers. The drag and lift coefficients of the micro-scale flat plate also display different trend for the high and low Reynolds numbers. The lift coefficients are always less than the drag coefficients when the Re < 100, suggesting that the aerodynamic performance of the micro-scale flat plate with 5% thickness ratio is poor when the flow Reynolds number is small. Acknowledgements This work were supported by Core Research for Evolutional Science and Technology (CREST) of Japan Science and Technology (JST) Agency (No.: 228205R) and the Japan Society for the Promotion of Science through a Grant-in-Aid for Scientific Research (B) (No. 18360050).
References [1] Mueller TJ, Fixed and flapping wing aerodynamics for micro air vehicle applications. In: Mueller T, editor. AIAA Progress in Astronautics and Aeronautics, vol. 195, 2001. p. 307–39.
[2] McMichael JM, Francis MS, Micro air vehicles—toward a new dimension in flight, DARPA, Available from: http://www.darpa.mil/tto/MAV/mav_auvsi. html, 1997. [3] Dickinson MH, Lehmann FO, Sane SP. Wing rotation and the aerodynamic basis of insect flight. Science 1999;284:1954–60. [4] Johnstone R, Arntz N. CONDOR—high altitude long endurance (HALE) automatically piloted vehicle (APV). AIAA 1990:1990–3279. [5] Wang J. Dissecting insect flight. Ann Rev Fluid Mech 2005;37:183–210. [6] Null W, Shkarayev S. Effect of camber on the aerodynamics of adaptive wing micro air vehicles. In: AIAA 2004-2694, the second AIAA flow control conference. Portland, Oregon, June 29–July 1; 2004. [7] Jones KD, Bradshaw CJ, Papadopoulos J, Platzer MF. Improved performance and control of flapping-wing propelled micro air vehicles. In: AIAA 2004-0399, the AIAA 42nd aerospace sciences meeting and exhibit. Reno, Nevada, January 5–8; 2004. [8] Sunada S, Sakaguchi A, Kawachi K. Airfoil section characteristics at a low Reynolds numbers. Trans ASME J Fluids Eng 1997;119:129–35. [9] Sun Q, Boyd ID. Flat-plate aerodynamics at very low Reynolds number. J Fluid Mech 2004;502:199–206. [10] Sun Q, Boyd ID. Simulation of gas flow over micro-scale airfoils using a hybrid continuum particle approach. AIAA 2003:2003–4014. [11] Sun Q, Boyd ID, Candler GV. Numerical simulation of gas flow over micro-scale airfoils. AIAA 2001:2001–3071. [12] Cercignani C. Theory and application of the Boltzmann equation. Edinburg: Sottish Academic Press; 1975. [13] Chapman S, Cowling TG. The mathematical theory of non-uniform gases. Cambridge University Press; 1970. [14] Ohwada T. Heat flow and temperature and density distributions in a rarefied gas between parallel plates with different temperatures: finite difference analysis of the nonlinear Boltzmann equation for hard sphere molecules. Phys Fluids 1996;8:2153–60. [15] Bird G. Molecular gas dynamics and the direct simulation of gas flows. Oxford: Clarendon Press; 1994. [16] Wang M, Li Z. Simulations for gas flows in microgeometries using the direct simulation Monte Carlo method. Int J Heat Fluid Flow 2004;25(6):975–85. [17] Wang M, Li Z. An Enskog based Monte Carlo method for high Knudsen number non-ideal gas flows. Comp & Fluids 2007;36(8):1291–7. [18] Müller I, Ruggeri T. Rational extended thermodynamics. New York: Springer; 1998. [19] Struchtrup H, Torrilhon M. Regularization of Grad’s 13 moment equations: derivation and linear analysis. Phys Fluids 2003;15(9):2668–80. [20] Wolf-Gladrow DA. Lattice-gas cellular automata and lattice Boltzmann models. Berlin: Springer; 2000. [21] Succi S. The lattice Boltzmann equation: for fluid dynamics and beyond (series numerical mathematics and scientific computation). Oxford New York: Oxford University Press; 2001. [22] Nie X, Doolen G, Chen S. Lattice Boltzmann simulations of fluid flows in MEMS. J Stat Phys 2002;107:279–89. [23] Ansumali S, Karlin IV. Kinetic boundary conditions in the lattice Boltzmann method. Phys Rev E 2002;66:026311. [24] Lim CY, Shu C, Niu XD, Chew YT. Application of lattice Boltzmann method to simulate microchannel flows. Phys Fluids 2002;107:2299–308. [25] Niu XD, Shu C, Chew YT. Lattice Boltzmann BGK model for simulation of micro flows. Eur Phys Lett 2004;67:600–6. [26] Shu C, Niu XD, Chew YT. A lattice Boltzmann kinetic model for microflow and heat transfer. J Stat Phys 2005;121:239–55. [27] Niu XD, Shu C, Chew YT. Numerical simulation of isothermal micro flows by lattice Boltzmann method and theoretical analysis of the diffusive scattering boundary condition. Int J Mod Phys C 2005;16:1927–41. [28] Zhang YH, Gu XJ, Barber RW, Emerson DR. Capturing Knudsen layer phenomena using a lattice Boltzmann model. Phys Rev E 2006;74:046704. [29] Chen H, Zhang R, Staroselsky I, Jhon M. Recovery of full rotational invariance in lattice Boltzmann formulations for high Knudsen number flows. Phys A 2006;362:125–31. [30] Zhou Y, Zhang R, Staroselsky I, Chen H, Kim WT, Jhon M. Simulation of micro and nano-scale flows via the lattice Boltzmann formulations for high Knudsen number flows. Phys A 2006;362:68–77. [31] Toshi F, Succi S. Lattice Boltzmann method at finite Knudsen numbers. Eur Phys Lett 2005;69:549–55. [32] Shan X, Yuam XF, Chen H. Kinetic theory representation of hydrodynamics: a way beyond the Navier–Stokes equation. J Fluid Mech 2006;550:413–41. [33] Ansumali S, Karlin IV, Arcidiacono S, Abbas A, Prasianakis N. Hydrodynamics beyond Navier–Stokes: exact solution to the lattice Boltzmann hierarchy. Phys Rev Lett 2007;98:124502. [34] Yudistiawan WP, Ansumali S, Karlin IV. Hydrodynamics beyond Navier– Stokes: the slip flow model. Phys Rev E 2008;78:016705. [35] Kim SH, Pitsch H, Boyd LD. Slip velocity and Knudsen layer in the lattice Boltzmann method for microscale flows. Phys Rev E 2008;77:026704. [36] Kim SH, Pitsch H, Boyd LD. Accuracy of higher-order lattice Boltzmann methods for microscale flows with finite Knudsen numbers. J Comp Phys 2008;227:8655–71. [37] He X, Luo L-S. Theory of lattice Boltzmann method: from the Boltzmann equation to the lattice Boltzmann equation. Phys Rev E 1997;56(6):6811–7. [38] Shan X, He X. Discretization of the velocity space in the solution of the Boltzmann equation. Phys Rev Lett 1998;80:65–8.
X.-D. Niu et al. / Computers & Fluids 38 (2009) 1675–1681 [39] Alexander FJ, Chen H, Chen S, Doolen GD. Lattice Boltzmann model for compressible fluids. Phys Rev A 1992;46:1967–70. [40] Zhang R, Shan X, Chen H. Efficient kinetic method for fluid simulation beyond the Navier–Stokes equation. Phys Rev E 2006;74:046703. [41] Chikatamarla S, Ansumali S, Karlin I. Entropic lattice Boltzmann models for hydrodynamics in three dimensions. Phys Rev Lett 2006;97(1):010201. [42] Ansumali S, Karlin IV, Oettinger HC. Minimal entropic kinetic models for hydrodynamics. Eur Phys Lett 2003;63:798–804. [43] Niu XD, Hyodo S, Munekata T, Suga K. Kinetic lattice Boltzmann method for micro gas flows: issues on boundary condition, relaxation time and regularization. Phys Rev E 2007;76:036711. [44] Latt J, Chopard B. Lattice Boltzmann method with regularized pre-collision distribution functions. Math Comp Simulation 2006;72:165–8.
1681
[45] Gombosi TI. Gaskinetics theory. Cambridge: Cambridge University Press; 1994. [46] Tamada K, Miura H. Slip flow past a tangential flat plate at low Reynolds numbers. J Fluid Mech 1978;85:731–42. [47] Schaaf SA, Sherman FS. Skin friction in slip flow. J Aeronaut Sci 1954;21:85–90. [48] Bhatnagar PL, Gross EP, Krook M. A model for collision processes in gases I: small amplitude processes in charged and neutral one-component systems. Phys Rev 1954;94:511–25. [49] Inamuro T, Maeba K, Ogino F. Flow between parallel walls containing the lines of neutrally buoyant circular cylinders. Int J Multiphase Flow 2000;26: 1981–2004. [50] Davis RW, Moore EF. A numerical study of vortex shedding from rectangles. J Fluid Mech 1982;116:475.