Pei-qing Liu, Yue Zhao, Qiu-lin Qu, Tian-xiang Hu
Lu Shijia Laboratory, Beihang University, Beijing 100191, China
Abstract: For correct identification of vortices, this paper first analyzes the properties of the rigid vortex core and its induced flow field given by the Rankine vortex model, and it is concluded that the concentrated vortex structure should consist of the vortex core and the induced flow field (the potential flow region with a weak shear layer). Then the vortex structure is analyzed by using the Oseen vortex model. Compared with the Rankine vortex, the Oseen vortex is a concentrated vortex with a deformed vortex core. The vortex structure consists of the vortex core region, the transition region and the shear layer region (or the potential flow region). The transition region reflects the properties of the resultant vorticity of the same magnitude and the resultant deformation rate of the shear layer, and the transition region also determines the boundary of the vortex core. Finally, the evolution of leading-edge vortices of the double-delta wing is numerically simulated. And with different vortex identification methods, the shape and the properties of the leading-edge vortices identified by each method are analyzed and compared. It is found that in the vorticity concentration region, the vortices obtained by using ω,2λ,Ω criteria and Q criteria are basically identical when appropriate threshold values are adopted. However, in the region where the vorticity is dispersed, due to the influence of the flow viscous effect and the adverse pressure gradient, the results obtained by different vortex identification methods can be quite different, as well as the related physical properties, which need to be further studied.
Key words: Vortex identification methods, Rankine vortex model, Oseen vortex model, leading-edge vortex, double-delta wing
The vortex is a common flow phenomenon in nature, closely related to various fluid mechanics fields, including turbulence[1-3], atmospheric phenomena[4-5], ocean physics[6-8], lift augmentation[9-10],drag reduction[11-12], flow control[13-15], and aeroacoustics[16-17]. For example, the trailing-edge vortices and the circulation of the flow around a wing have a dominant impact on the aerodynamic loading of the aircraft. In turbulence studies, it is expected that the mechanism of the turbulence generation and dissipation can be revealed through the study of vortices, and the relationship between the turbulence coherent structure and the flow parameters can be established[1]. The British scientist Küchemann, a student of the German scientist Prandtl, once said(1965) that the vortices are the tendons of the fluid motion[18]. Another student of Prandtl, Shijia Lu, a Chinese scientist, said that the essence of the fluid is the vortices, because the fluid cannot stand "rubbing",and once being “rubbed” the vortices come out. At present, it is generally accepted that vortices consist of a group of fluid micro-clusters rotating around a common center[19]. The spatial structure of the vortices is usually of three kinds: the spiral vortex, the disc vortex and the columnar vortex (Fig. 1)[20].
By definition, the vorticity is two times the rotational speed of the fluid micro-cluster in the flow field around its own axis during the motion[21]. This is a purely kinematic physical quantity, which represents the rotation speed of the fluid micro-cluster. The mathematical expression of the vorticity is
The vortex flux refers to the sum of the vorticity passing through an arbitrary surface (Fig. 2). For three-dimensional vortices, the vorticity flux is expressed as
Fig. 1 (Color online) Spatial structure of vortex
Fig. 2 Definition of vortex flux
The velocity circulation refers to the line integral of the velocity along an arbitrary closed curve in the flow field, which represents the accumulated vorticity in a region. According to the Stokes formula, the relationship between the accumulated vorticity and the circulation of the flow field is
Currently, most of the researches of the vortex identification focus on describing the vortical feature based on some subjective conjecture without an in-depth understanding of the physical properties of the vortex. In order to see the shape and the properties of a coherent vortex, we need to categorize the types of the vortex. First, the vortex with a “rigid core” and the induced flow field properties is represented by the Rankine vortex model. It is known that a coherent vortex consists of the vortex core and the induced flow field (the shear layer). Secondly, the Oseen vortex model describes a deformed vortex, with the vortex core region, the transition region and the shear layer region. Based on above two models of vortices,the shape and the characteristics of (numerically simulated) leading-edge vortices for the double-delta wing are analyzed, in which different identification methods are employed to describe the leading-edge vortices and to make related evaluations. The aim of the present work is to help future researches in choosing an appropriate method to identify the vortex in the flow field.
British scientist W. J. M. Rankine (1820-1872)firstly studied the rigid vortex core and its induced flow field[22]. It is assumed that the vortex core is a rigid cylinder with constant vorticity and the velocity field outside the vortex core satisfies the induced velocity field generated based on the Biot-Savart law,and a combined vortex model is proposed, which consists of a rigid vortex core and its induced shear layer, so it is called the Rankine vortex model. This combined vortex model satisfies the exact solution of the Navier-Stokes equations. The vorticity in the vortex core is uniformly distributed, and the induced velocity field outside the vortex core is in the hyperbolic function distribution, therefore with only velocity gradient but with no vorticity. The basic flow field can be described as:
In the vortex core (the radius of the vortex core isR), the circumferential velocity is
The vorticity distribution in the vortex core is uniform,and its size is
The accumulated vorticity (the circulation) within the core region is
Outside the vortex core, the induced tangential velocity is
According to the definition of the vorticity
Substituting Eq. (7) into Eq. (8), the accumulated vorticity is obtained to be equal to zero, which means that the vorticity outside the vortex core is equal to zero, and the induced velocity field is a potential flow field. In a polar coordinate system, the two times of the shear deformation rate of the fluid filament is
Substituting Eq. (7) into Eq. (9), it can be seen that outside the vortex core the shear deformation ratezγis
It is shown that the shear deformation rate is nonzero for the velocity field induced outside the vortex core. Ifp∞is the farfield pressure, the difference between the farfield pressure and the pressure at the center of the vortex core can be calculated as
The structure of the Rankine vortex and its pressure distribution along the radius are shown in Fig.3(a). The velocity distribution of the Rankine vortex model is shown in Fig. 3(a). The circular region in the middle is the rigid vortex core region. The circumferential velocity of the fluid micro-clusters in this region increases linearly with the distance from the center. Outside the circular region is the shear layer region, in which the circumferential velocity of the fluid micro-clusters is inversely proportional to the distance from the vortex center. The pressure distribution of the Rankine vortex model is shown in Fig. 3(b).The pressure distribution is funnel-shaped, and the pressure in the center of the vortex is the lowest.Figure 3(c) shows an example of the Rankine vortex schematically.
According to the Rankine vortex model and the conservation law of the Helmholtz (1821-1894)vorticity motion[23], the structure of tornado vortices formed by the three-dimensional rising (or descending)airflow can be constructed as shown in Fig. 4. For the concentrated vortices, the following qualitative conclusions can be reached by analyzing the Rankine vortex:
Fig. 3 (Color online) Rankine vortex model
Fig. 4 (Color online) Tornado structure caused by rising (hot)or descending (cold) airflow (I-Deceleration and pressure decrease region, II-Acceleration and pressure increase region, III-Static pressure distribution)
(1) There is a group of fluid micro-clusters rotating around a common axis, which is also the basic property of the coherent vortex.
(2) The concentrated vortices are composed of the vortex core region and its induced region.
(3) Inside the vortex core region, the fluid micro-clusters rotate around the common axis. In the induced region outside the vortex core, the shear deformation rate is dominant, the vorticity tends to be zero, and the flow field can be seen as a potential flow area.
(4) Inside the vortex core region, there is a pressure drop due to the rotational motion of the flow filaments.
Along with the Rankine combined vortex model,in 1912, the Swedish physicist Oseen described the viscous flow field with consideration of the vortex diffusion behavior. By solving the Navier-Stokes equations, a vortex model was developed to characterize the vortex diffusive phenomenon of viscous flow, namely, the Oseen vortex model[22]. The vorticity transport equation is obtained by taking the curl of the Navier-Stokes equations. For an axisymmetric two-dimensional flow, the vorticity diffusion equation of the vorticity in thezdirection can be expressed as
Fig. 5 (Color online) Radial distribution of dimensionless characteristic quantities
According to the Oseen vortex model, the basic characteristics of the deformed vortex core and its induced flow field can be characterized as follows:
(1) Comparing with the rigid vortex (the Rankine vortex) of ideal fluids, the vortex in viscous fluids is deformed with radius proportional to the square root of the time and is diffused in accordance with the law of the viscous diffusion.
(2) For the viscous fluid, the coherent vortical structure consists of the vortex core region, the transition region and the shear layer region.
(3) In the vortex core region, the vorticity is predominant, and the shear deformation rate is small.Outside the vortex core, there is the vortex core induced region, including the transition region and the pure shear region (the potential flow region).
(4) The shear layer is located in the further outside of the vortex core, and the shear deformation rate is predominant, and the vorticity tends to zero.
(5) The vorticity in the transition region, the boundary of the vortex core, is of the same magnitude as the shear deformation rate.
(6) At the vortex axis, the shear deformation rate is zero, the vorticity is the largest, and the fluid micro-clusters are in pure rotational motion. Based on this feature, Prof. Chaoqun Liu proposed the Liutex vortex identification method[24-25](the method of pure rotation of the vortex axis).
(7) Vortex core is the tendon of the fluid movement and the soul in controlling the fluid movement.
(8) The generation of vortices depends on“rubbing”, and the dissipation of vortices also depends on “rubbing”. The “rubbing” can not only produce vortices, but also dissipate them.
In an incompressible flow, the angular rotation rate tensor is a physical quantity that represents the degree of rotation of the fluid micro-clusters, and the strain rate tensor is a physical quantity that represents the strain degree of the fluid micro-clusters.According to the Helmholtz's law of the velocity field decomposition[23], the velocity gradient tensor is defined as
It can be uniquely partitioned into symmetric and antisymmetric parts
Ωis the angular rotation rate tensor of the fluid micro-clusters, which is defined as
Ωis the norm of the angular rotation rate tensorΩ,indicating the degree of rotation of the fluid microcluster
Sis the strain rate tensor of the fluid micro-cluster,which is defined as
Sis the norm of the strain rate tensorS, indicating the strain degree of the fluid micro-cluster,
From the Oseen vortex model, it can be seen that the core region is dominated byΩ, withStaking the secondary part, and the norm of the strain rate tensor is zero at the vortex axis, where the fluid micro-clusters are in the pure rotational motion. In the outer shear layer induced by the vortex core,Sis dominant andΩtends to zero. In the boundary region of the vortex core (the transition region), one sees the transition of the fluid micro-clusters from the rotation to the shear, andΩandSare in the same order of magnitude, as shown in Fig. 6.
Fig. 6 (Color online) Vortex core boundary region
In the fluid mechanics, the identification method of vortices is one of the hot topics. For the leadingedge vortices in the case of double-delta wings, the rationality of various methods of identifying the equivalent surface vortices is discussed in terms of the shape and the properties of the leading-edge vortices.The DDES[26-27]turbulence model is used to simulate the flow field of the static double-delta wing, and the results are verified by the water tunnel force measurement experiment. Then, the leading-edge vortices for the double-delta wing atα=20° are extracted by theωcriterion, theQcriterion, the 2λcriterion and theΩmethod, respectively. The shape and the physical properties of the leading-edge vortices extracted by different vortex identification methods are analyzed, including the leading-edge vorticity intensity, the vorticity distribution, the radius of the vortex core and the circumferential velocity distribution.
The experiment was conducted in the low Reynolds number recirculation water channel at Lu Shijia Laboratory of Beihang University, as shown in Fig. 7(a). The model used in the experiment is shown in Fig. 7(b). The maximum chord length of the double delta wing model is 188 mm, the leading-edge sweep angle of the strake wing is 75°, and the leading-edge sweep angle of the main wing is 50°, the thickness of the model is 4mm, the leading-edge and the trailingedge of the model are chamfered at 45° on both sides,the flow velocity at infinity is 0.15 m/s, and the Reynolds number of the incoming flow with the chord length as the reference length is 1.8×104.
The model used in the numerical simulation is shown in Fig. 8(a). In order to reduce the number of grids and save the computing resources, a half model is used in the numerical simulation of the flow field of the double delta wing in this study. The flow velocity at the infinity in the calculation is the same as that in the water channel experiment, which is 0.15 m/s. In Fig. 8(b), from the vorticity distribution in different sections of the double-delta wing atα=20°, one sees leading-edge vortices on both the strake wing and the main wing.
Fig. 7(a) (Color online) low Reynolds number recirculation water channel
Fig. 7(b) Experiment model
Fig. 8(a) Double-delta wing half model
Fig. 8(b) (Color online) Cross section vorticity distribution(α = 20°)
Fig. 8(c) (Color online) Lift coefficient verification
The comparison between the experiment lift results and the numerical simulation results is shown in Fig. 8(c). In the figure, “CFD-c” represents the calculation result when the number of grids is 3×106,and “CFD-f” represents the calculation result when the number of grids is 7×106. It can be seen from Fig.8(c) that the results of experiment and calculation are in good agreement, indicating that the numerical simulation enjoys a high accuracy.
This method is based on the conditions of vorticity dominance in the vortex core area, so the vortex core boundary is determined by using the condition that the resultant vorticity is larger than a threshold[21]. The choice of the vorticity thresholds is a problem, and with different thresholds, different radii of the vortex core will be obtained. This criterion is relatively simple and does not take into account the characteristics of the shear deformation rate. This method is invalid for concentrated vortices when the maximum vorticity is not at the center of the vorticity axis. As shown in Fig. 9, differentωthreshold iso-surfaces have a great influence on the shape of vortices. If the maximum vorticity is located in the boundary layer, it should be eliminated in the process of the vortex identification. For the iso-surface ofω= 30, the vorticity distributions on the cross sectionsx/c=0.3, 0.6 and 0.9 are shown in Fig. 10.On all cross-sections of the leading-edge vortex, the different vorticity iso-lines are quite different.
Fig. 9 (Color online) Vortex structure of different vorticity threshold iso-surfaces (α = 20°)
Fig. 10 (Color online) Vortex shape at ω = 30 and crosssection vorticity distribution (α = 20°)
The method is mainly based on the idea that the resultant vorticity and the resultant shear deformation rate are of the same magnitude at the boundary of the vortex core[28-30]. For incompressible fluid motion, theQvalue is defined as
Qcriterion requiresQ>threshold (generally greater than zero), with the pressure funnel in the vortex core region. The vortices obtained with differentQvalues are shown in Fig. 11 (at the angle of attack of 20°). It can be seen that differentQvalues have different effects on the vortices in different regions.Vortex lines are mainly distributed in the iso-surface,and the results of the identification for the concentrated vortex are reasonable. In the case of the angle of attack of 20°, the cross-section vorticity distribution withinQ=100 iso-line and theQ-value distributions on the cross-sectionsx/c=0.3, 0.6 and 0.9 are shown in Fig. 12. Most of vorticity is distributed within theQ=100 iso-line. In the front of the wing, the leading-edge vorticity is concentrated,and the vortex structures obtained by taking differentQvalues are not different. In the rear of the wing,due to the influence of the viscous and adverse pressure gradient, the leading-edge vorticity diffuses, and the vortex shapes obtained by taking differentQvalues are quite different. In this case, the selection ofQthreshold has a great influence on the vortex shapes.
Fig. 11 (Color online) Vortex structure of different Q threshold iso-surface (α = 20°)
Fig. 12 (Color online) Vortex shape at Q=100 , cross- section vorticity distribution withinQ=100 iso-line at(a)x/ c=0.3,(b) x/ c=0.6,(c) x/ c=0.9 and crosssection Q-value distribution at(d) x/ c=0.3,(e)x/ c=0.6 and (f) x/ c=0.9 (α = 20°)
Letλ1,λ2,λ3be three eigenvalues of the matrixΩ2+S2, andλ1≥λ2≥λ3. For incompressible flow, theQcriterion can be expressed as
According to Eq. (22), theQcriterion is not exactly the same as that of the2λcriterion. The specific difference is shown in Table 1. In some cases, they give agreed results and in other cases they do not,depending on whether a point belongs to a vortex core.So the difference of the vortex identification results between theQcriterion and the2λcriterion does not all come from the threshold value, and the difference of the method will also cause a difference of results. In addition, in terms of the physical essence,theQcriterion holds that there are vortices in the region where the rotation is greater than the strain, but 2λcriterion is based on the notion that a pressure minimum should exist across a vortex core. Therefore,there are differences between the two vortex identification methods.
Fig. 13 (Color online) Vortex structure of different λ2 threshold iso-surface (α = 20°)
Table 1 Comparison between Q criterion and λ2 criterion
Fig. 14 (Color online) Vortex shape at λ2 =-100, cross-section vorticity distribution within λ2 =-100 iso-line at (a) x/ c=0.3,(b) x/ c=0.6 and (c) x/ c=0.9 and cross-section λ2-value distribution at (d) x/ c=0.3, (e) x/ c=0.6 and (f) x/ c=0.9 (α= 20°)
TheΩcriterion indicates that in the boundary region of the vortex core, the resultant vorticity and the resultant shear deformation rate are of the same magnitude[32-35]. For the incompressible fluid motion,theΩvalue is defined as
In the boundary region of the vortex core, Liu suggested that Ω = 0.52[36]. For ease of comparison,Fig. 15 shows the iso-surface ofΩ= 0.52,Q=100 andλ2= -1 00 atα=20°, showing that the vortices are basically identical in structure. Similarly,in the front of the wing, the leading-edge vorticity is concentrated, and the vorticity structures obtained by the three methods are basically the same. In the rear of the wing, the vorticity diffuses along the leading-edge due to the influence of the viscous and adverse pressure gradient, and the vortex shape obtained by different vortex identification methods is different.
Liutex method is a new method of the vortex identification proposed by Prof. Chaoqun Liu in recent years[25]. This method mainly identifies the axis of the vortex by finding the pure rolling region. The definitions are as follows
Fig. 15 (Color online) Comparison of vortex shapes obtained by different vortex identification methods (α = 20°)
whereQis the transformation matrix that makes the following formula hold
The Rotex vector obtained by the Liutex method is the axis of the rotation at the center of the vortex core. The core and the shear layer of the vortex revolve around the axis in the direction of the Rotex vector.With this method, one can only identify the revolving center line of the vortex, but not the radius of the core and the region of the vortex. It should be used in practice together with other identification methods.The verification example used in this paper is the identification of the concentrated vortices on the leading edge of double delta wings. For the concentrated vortices on the leading edge, the location and the shape of the vortex axis are not the key points of the identification. So this method is not suitable for the identification of leading-edge vortices.
The main characteristics and the vorticity distribution along the leading edge of the delta wing are shown in Fig. 16. The vorticity layer gradually coils from the wing head to the trailing edge, and the vortex core gradually enlarges[9]. According to Fig. 16,the vorticity distribution curves at different sections are shown in Fig. 17. The maximum vorticity in the vortex core region decreases gradually along the path,and the radius of the vortex core increases gradually.Figure 18 shows the circumferential velocity distribution around the vortex core at the cross-sectionx/c= 0.4, in agreement with the distribution law of the Oseen vortex.
Using the data of Fig. 17, as shown in Fig. 19,taking the central point of the vortex core as the center,the vorticity integral of the circular region with different radii is obtained at the cross sectionx/c=0.4, and the curve of the variation of the velocity circulation with the vortex radius is obtained.It is obvious that the radius of the vortex determined by the three methods does not reach the maximum circulation value. We all know that neither the vorticity nor the circulation can be used to define the vortex, but for the leading edge vorticity of double delta wings, the vorticity integral in the circular region with different radius represents the concentration and the intensity of vorticity in this region. The circular region with the concentrated vorticity can be regarded as the region where the leading-edge vortex exists.Before the vortex breakdown, the leading edge vortex of the double-delta wing can be regarded as a rigid rotating concentrated vortex, and the boundary of the vortex should be near the place of the maximum value of the circulation. As shown in Fig. 19, the maximum circulation value of the three vortex identification methods is not reached, which may be due to the selection of the threshold value of the vortex identification method.
Fig. 16(a) Winding process of leading-edge vortex
Fig. 16(b) vorticity distribution of cross-section[9]
Fig. 17 (Color online) Vorticity distribution of different crosssections
Fig. 18 (Color online) The circumferential velocity distribution around the vortex core at the cross section
Fig. 19 (Color online) Integral of vorticity of different radii at the cross section of x/ c=0.4
Fig. 20 (Color online) Integral of vortex intensity along different sections of three different vortex identification methods
For the strake wing vortex, whenx/c< 0.6, the vorticity in the free shear layer separated from the leading-edge continuously enters the vortex core, and the vorticity intensity increases gradually as the vortex develops downstream. Whenx/c> 0.6, the edge vortex core moves over the main wing with no free shear layer involved, so the vorticity does not increase,as shown in Fig. 20. It can be seen that the vorticity is concentrated in the front of the wing, and the vortex shape and the distribution of the vorticity intensity along the wing are basically the same as obtained by different vorticity identification methods. However, in the rear of the wing, the vorticity diffusion increases,and the vorticity intensity obtained by different methods of vortex identification is quite different. The axial velocity distribution of the concentrated vortices is shown in Fig. 21. Whenx/cis less than 0.6, the axial velocity distribution of the edge vortices is of jet type, which is less affected by the inverse pressure gradient. Whenx/cis greater than 0.6, the axial velocity distribution of the edge vortices gradually changes from the jet type to the wake type, which is greatly affected by the inverse pressure gradient, as shown in Fig. 22.
Fig. 21 (Color online) Distribution of vortex axis velocity at different chord positions (front wing)
With the numerical simulation results of leadingedge vortices of double-delta wing, the shape and the characteristics of the vortices are analyzed by different methods of vortex identification. It is shown that the results obtained by usingω,Q,2λ, andΩcriteria are similar and basically identical for the vorticity concentrated vortices as long as appropriate threshold values in the methods are chosen. For the vorticity dispersed vortices, the results are highly dispersed. The following conclusions can be drawn:
(1) Physically, theωcriterion highlights the vorticity in the vortex core region, theQcriterion and theΩmethod require that the resultant vorticity and the resultant shear deformation rate have the same order of magnitude at the boundary of the vortex core;and the2λcriterion requires that there exists a minimum pressure region in the center of the vortex core. Various methods of vortex identification take different thresholds, and the shape and the properties of vortices are different.
Fig. 22 (Color online) Distribution of vortex axis velocity at different chord positions (rear wing)
(2) In the front of the wing, the leading-edge vortices of the double-delta wing identified by theQcriterion, the2λcriterion and theΩmethod are basically identical in the vorticity concentration region.
(3) In the rear of the delta wing. the leading-edge vortices of the double-delta wing identified by theQcriterion, the2λcriterion and theΩmethod are dispersed in the vorticity dispersion region. The shape and the physical properties of the vortices are highly dispersed due to the influence of the viscous and adverse pressure gradients. The accuracy of each vortex identification method cannot be rightly judged.
(4) The thresholds of theQcriterion and the 2λcriterion are difficult to choose, while the threshold of theΩmethod can be given (0.52).
水動(dòng)力學(xué)研究與進(jìn)展 B輯2020年5期