Abstract
In the current study, numerical simulations are performed in order to investigate effects of incidence angle and airfoil thickness on alternating vortex pattern of symmetric airfoils at Re = 1000. This alternating vortex pattern is found to be significantly varying in shape as the incidence angle increases. The results are obtained with 1° increment from 0° to 41° and then with 10° increment from 40° until 180°. The instantaneous and mean vortex patterns are investigated around 2% thick (NACA 0002) and 12% thick (NACA 0012) airfoils. The mean lift, drag, and pitching moment coefficients in addition to Strouhal number are computed and compared with the results available in literature. It is found that the wake behind the airfoil exhibits a continuous vortex shedding pattern below 8° incidence angle for NACA 0002 and below 7° incidence angle for NACA 0012 at Re = 1000. The wake structures are classified into five different modes according to their pattern obtained from instantaneous and mean vorticity fields by also taking into account the amplitude spectrum of the lift coefficient, the instantaneous and mean aerodynamic force coefficients, velocity fields, and the longitudinal and lateral vortex spacings.
Introduction
The Reynolds number for species such as hummingbirds is Re = O (104), while for smaller insects, such as fruit flies or honeybees, Re = O (102−103). 1 In low Reynolds number aerodynamics, thickness of the wings become as small as possible.
It is known that bird or insect wings have thin cambered airfoils at low Reynolds numbers. Different researches on airfoil thickness at these low Reynolds numbers reveal also the advantage of thin airfoils on aerodynamic forces compared to thick airfoils. However, there are very few research on vortex shedding pattern which leads to different flow structures as the airfoil thickness is decreased and the ones which exist are limited to low angles of attack.
Different numerical simulations are carried out by researchers around airfoils for Re = O(102−104) but mainly for low angles of attack ranges.
From Stanford University, Kunz and Kroo 2 used a two-dimensional incompressible Navier–Stokes solver INS2D which utilizes the artificial compressibility method. Analysis were done at angles of attack lower than about 10° with Reynolds number ranging from 1000 to 12,000 for different thickness and cambered NACA profiles.
Mateescu and Abdo 3 used a pseudo-time integration method using artificial compressibility to accurately solve the Navier–Stokes equations for Reynolds number ranging from 400 to 6000 at low angles of attack. Their results are very close to Kunz and Kroo 2 solutions.
Ahuja and Rowley 4 have performed a study on a flat plate at Re = 100 from 0° to 90° incidence angles. Khalid and Akhtar, 5 Liu et al., 6 Hoarau et al., 7 and Mittal and Tezduyar 8 have shown numerical simulations around NACA 0012 airfoil at Re = 1000 for different angles of attack. Their results are compared with a previous work of author in Kurtulus. 9
The vortex street formations at low Reynolds numbers are investigated in literature mostly around bluff bodies and in particularly around circular cylinders.10–13 Bluff bodies at certain Reynolds numbers exhibit alternate vortex shedding, called as von Karman vortex, consisting of alternating vortices of opposite signs. The time averaged velocity profile for Karman vortex structure is a momentum deficit profile. The bluff body creates drag and the mean velocity at the wake is lower than the free-stream velocity.
Schaefer and Eskinazi 14 defined three basic regions of vortex street behavior behind a bluff body as a formation region, where the vortex street is developed, a stable region where the vortices display a periodic laminar behavior, and unstable region where the street disappears and turbulence develops for 50<Re<125.
Durgin and Karlsson 15 investigated a von Karman vortex street subject to a deceleration which causes change in the ratio of longitudinal to lateral spacing between vortices. The deceleration of the flow is performed by inserting a large cylinder in the downstream of a vertically placed thin vortex shedding cylinder. From mean velocity and vorticity profiles obtained from hot-wire measurements, they observed three regions behind the shedding cylinder as primary vortex region, calm region (in the field where no velocity fluctuations were visible), and secondary vortex region. They denoted that the calm region obtained from mean velocity profiles consists of two adjacent layers of vorticity of opposite signs. They noted (as was also explained by von Karman) that one vortex can be deformed by the combined action of all the others. In addition, the shape of a vortex deviated from its idealized circular shape can be further changed by its own induced velocities. They also indicate that similar processes may be present in turbulent flows.
Strykowski and Sreenivasan 16 shows that vortex formation at low Reynolds numbers can be completely suppressed by suitably placing in the wake of the vortex shedding cylinder a smaller cylinder (called as control cylinder). In their experiments at Re = 80, a circular cylinder of 1/7 times diameter of the shedding cylinder was placed at the downstream of the shedding cylinder and hydrogen-bubble visualization technique in addition to numerical simulations reveals suppression of the vortex street formation at this Reynolds number. An increase in the Reynolds number by the increase of the free-stream velocity reveals once more the vortex street. It is important to understand the vortex street formations behind airfoils at low Reynolds numbers for the performance and control of micro air vehicles (MAV). This knowledge could lead the application of different control problems for MAV applications.
Unsteady vortex pattern at the wake of the airfoils is mainly investigated around the pitching/flapping airfoils at low Reynolds numbers. The characteristics of the wake patterns are not investigated as much as force measurements on oscillating wings.
Bratt 17 investigated the flow patterns in the wake of an oscillating NACA 0016 airfoil by using smoke visualization technique. He has also performed some theoretical calculations for the interpretation of the smoke traces. He describes that the wake can be treated as a system of discrete vortices and the elements of a continuous vortex sheet could be subject to normal induced velocities, which would tend to break the sheet up. The induced velocities for a uniform infinite vortex sheet would be zero.
Koochesfahani 18 investigated experimentally the vortical flow patterns in the wake of NACA 0012 airfoil pitching at small amplitudes in a water channel. The experiments are performed at Re = 12,000. It is concluded that oscillation frequency, amplitude, and waveform shape affect the structure of the wake. At low frequencies, an undulating vortex pattern is observed. At higher frequencies, an alternating vortex pattern is observed opposite to typical Karman vortex street which corresponds to a jet pattern. This vortex arrangement corresponds to the thrust generating airfoil. He also observed a case with no momentum deficit in the wake (no drag) when the alternating vortices are positioned exactly on a straight line. Similarly, Freymuth 19 and Lai and Platzer 20 experimentally show thrust producing flapping airfoil configurations depending on flapping frequency and amplitude. Huang et al. 21 investigated experimentally in a towing water tank, the formation and shedding of the vortex system on the suction surface of an impulsively started NACA 0012 wing.
The present study aims to analyze numerically the unsteady flow structure around NACA 0002 in comparison with NACA 0012 airfoil 9 at Re = 103 for angles of attack in the range of 0°–180° using a finite volume-based package program ANSYS Fluent. 22 The flow structure analysis includes the vortex formation, wake shedding frequency, and their effects on aerodynamic coefficients.
Numerical method
Governing equations and geometry
The governing equations are Navier–Stokes equations for incompressible, laminar, and two-dimensional flow
The flow around NACA airfoil is obtained at Re = 1000 steady external conditions. The thickness distribution of NACA four digit airfoils,
The angle of attack, defined to be positive in the clockwise direction, was varied from 0 to 180° with an increment step of 1 for α < 41° and 10° for 40° < α < 180°. The pivot point was at the quarter chord location (0.25 c) from the leading edge of the airfoil. It should be noted that pitching moment calculations are performed with respect to this pivot point for the current study.
Grid and time refinement studies have been carried out for NACA 0012 airfoil in a previous study of the author.
9
The same number of grid and time increments is also used for NACA 0002 airfoil using the MATLAB program written to create a journal file for grid generation in GAMBIT software; 300 nodes are used around the airfoil for the medium mesh selected after the grid refinement study compared to 150 nodes at coarse grid and 500 nodes at fine grid. The first cell spacing of the boundary layer for medium mesh has value of 0.0015 c. Total number of elements are in the order of 200,000. The mesh around the airfoil is composed of two parts, namely inner region and outer region. The inner region above the boundary layer is constructed with a semi-circle having radius of 1.5c centered at the airfoil c/4 location at the upstream and a rectangular region with width of 1.5c at the downstream of the airfoil. The inner region mesh is unstructured triangular grid. The outer domain is of C-type structured mesh with 15c radius (Figure 1). With the preprocessor program, the airfoil is rotated for the given angle of attack in the inner region by keeping the outer domain and wake region fixed with the same structured grid. The time refinement study carried out with time increments of Δt = 0.0025 s, Δt = 0.005 s, and Δt = 0.01 s in Kurtulus
9
reveals that Δt = 0.005 s in comparison with Δt = 0.0025 s give almost same results; therefore, the study is performed with time increment of Δt = 0.005 s.
Computational domain around NACA 0002 airfoil.
In the current paper, the results are simulated until t = 100 s which corresponds to a non-dimensional time of t* = t.Uinf/c = 146, where Uinf = 0.146 m/s is the free-stream velocity and c = 0.01 m is the chord length of the airfoil. The computational time interval is t*∈[0 146]. The mean contours and values are obtained by averaging the parameters in the interval of 73 ≤ t* ≤ 146 which corresponds to 50 s ≤ t ≤ 100 s. Cl, Cd, and Cm denote the lift, drag, and moment coefficients for the airfoil, respectively. The corresponding mean values are calculated in the interval of t*∈[73 146] and denoted with Cl, Cd, and Cm.
Results
Aerodynamic coefficients
The mean lift and drag coefficients are compared with the results of Kunz and Kroo
2
which are given for angle of attack values of 0° to 8°, and the results are found to be very similar as can be seen from Figure 2. The aerodynamic force coefficients are also compared with a previous study for NACA 0012 airfoil.
9
From this comparison, it is concluded that thin airfoils have very high Cl/Cd values at low Reynolds numbers of order 1000. Drag change does not change as much as lift change as angle of attack increases for low angles of attack until 30° and this results a high Cl/Cd = 4.5 value at α = 7° for NACA 0002 airfoil. The maximum Cl/Cd value reached for NACA 0012 is 2.55 at α = 11°. It is also noted that the mean moment coefficient with respect to quarter chord point is negative until α = 11° for NACA 0012 airfoil. For NACA 0002 airfoil, mean Cm is always positive with respect to quarter chord location which shows an unstable characteristic at zero angle of attack for NACA 0002 airfoil at this Reynolds number.
Mean aerodynamic coefficients at Re = 1000 for NACA 0002 airfoil and their comparison with literature data and NACA 0012 airfoil (0° ≤ α ≤ 30°).
Different from previous study,
9
NACA 0012 airfoil solutions for α > 90° with 10° increments are also added in the current study and the aerodynamic coefficients ranging between 0° and 180° are compared with respect to the thinner airfoil NACA 0002 as shown in Figure 3.
Mean aerodynamic coefficients at Re = 1000 for NACA 0002 and NACA 0012 airfoils (0° ≤ α ≤ 180°).
For low angles of attack (α ≤ 7°), there is a significant reduction in lift curve slope (Figure 4) which is also observed by Kunz and Kroo.
2
Then there is a sudden increase in Clα at α = 9° for NACA 0002 and α = 8° for NACA 0012. Amplitude spectrum of Cl coefficient reveals vortex shedding with a Strouhal number of 1.027 at α = 9° for NACA 0002. Similarly, α = 8° is the angle of attack where the vortex shedding is started to be observed for NACA 0012.
9
Strouhal number based on the amplitude spectrum of lift coefficient is shown in Figure 5. Maximum Strouhal number is reached at α = 9° for NACA 0002 airfoil with St = 1.027. It is noted that the leading edge separation is delayed in thin sections.
2
The results for Strouhal number in Figure 5 are also compared with the simulations in literature for NACA 0012 airfoil7,8 at Re = 1000 and are found to be in good agreement.
Clα for NACA 0002 and NACA 0012 at Re = 1000. Strouhal number at Re = 1000 for NACA 0002 and NACA 0012 airfoils (0° ≤ α ≤ 180°).

The lift curve slope of NACA 0002 airfoil and NACA 0012 airfoil at Re = 1000 are shown in Figure 4. It is found that the thickness plays an important role at these Reynolds numbers and lift drastically increases at small angles of attack. It is seen that the curve slope was equal to π value at α = 0° and decreases almost linearly for NACA 0012 until 7° angle of attack. Until this angle of attack, the vortex sheets are not alternating and a continuous vortex sheet mode was observed in the flowfield. Similarly, the lift curve slope for NACA 0002 is observed to decrease until 8° angle of attack. The lift curve slope for NACA 0002 is found to decrease very slightly until 3° and then decreases steeply. The increase in Clα as the thickness of the airfoil increases results also very high Cl/Cd values at low angles of attack as is shown in Figure 2.
In the previous study, 9 for NACA 0012, α = 8° was found to be the angle of attack where a low amplitude vortex shedding is observed by investigating the Strouhal number based on the amplitude of the Fast Fourier Transform of the instantaneous lift coefficient. Similarly, in the current study this angle of attack is found to be α = 9° for NACA 0002 airfoil. This shows that the thickness of the airfoil affects the vortex shedding angle of attack at these low Reynolds numbers.
The Strouhal number is defined by equation (4).
Vortex pattern modes
The vortex pattern at the downstream of NACA 0002 and NACA 0012 airfoils are investigated at different angles of attack for Re = 1000. The intensity and location of the shedding vortices influence the wake pattern. Vortices shed into the wake are also interacting with each other as they move downstream. The vortex shedding patterns are classified into five different vortex shedding modes according to the vortex pattern, amplitude spectrum of lift coefficient, and the longitudinal and lateral spacing of the vortices.
Mode 1: “continuous vortex sheet mode” (α ≤ 8° for NACA 0002 and α ≤ 7° for NACA 0012)
It is observed that angle of attach highly influences the vortex shedding phenomena and periodic variation of the flowfield. A MATLAB code is written to detect vortex cores in the flowfields. At each x station, the vorticity distribution is analyzed and maximum positive vorticity values (red line corresponding to counter-clockwise vortices) and maximum negative vorticity values (blue line corresponding to clockwise vortices) are drawn as a function of x/c as shown in Figure 6. It is observed that the vorticity distribution for continuous vortex sheet mode does not show any oscillatory behavior; however, other modes due to the formation of vortices show oscillations. The peak values of these vorticity distributions are detected to locate the core of the vortices. In the meantime, this approach is used to obtain the instantaneous and mean velocity of the vortex cores. This technique is found to be very satisfactory to detect aligned vortex core locations in the flowfield. The clockwise vortices are marked with blue circles and the counter-clockwise vortices are marked with red circles. In the meantime, longitudinal and lateral spacing of the vortices are also calculated for each individual vortex detected.
Example vorticity distribution along x/c for Mode 1 and Mode 2 cases. (a) α = 8° for NACA 0002 (Mode 1). (b) α = 9° for NACA 0002 (Mode 2).
Mode 1 is “ Vortex street for Mode 1. (a) α ≤ 8° for NACA 0002. (b) α ≤ 7° for NACA 0012.
A transition is observed between Modes 1 and 2 at α = 8° for NACA 0012 airfoil. This transition should occur between α = 8° and α = 9° for NACA 0002 airfoil. In this transition angle of attack, the wake of the instantaneous vorticity contours is in continuous vortex sheet mode close to the airfoil region and the successive counter-rotating vortices are observed at the downstream (Figure 8 top). The vortex street transforms from a continuous vortex sheet mode to an alternating rows of clockwise vortices at the top and counter-clockwise vortices at the bottom. The vortex shape is also observed to be changed from a droplet like shape, towards a triangular shape and then it stretches in horizontal direction by transforming to an elliptical shape at the downstream. Funakoshi
25
shows that lateral spacing and longitudinal spacing of Karman vortex street influence the evolution of uniform or non-uniform circular vorticity regions. Figure 8 also shows the mean vorticity contour for α = 8° of NACA 0012 airfoil. It is observed that the negative vorticity regions (blue contours) and positive vorticity regions (red contours) deviate towards up and down respectively at the downstream.
Vortex street for the transition between Mode 1 and Mode 2 at α = 8° for NACA 0012.
The transition from continuous mode to alternating mode is also visible from instantaneous aerodynamic force coefficients (Figure 10) and amplitude spectrums of the lift coefficient (Figure 9). As it was previously stated, the highest Cl/Cd ratio of 4.5 is obtained at α = 7° for NACA 0002 airfoil. The moment coefficient is seen to be positive for NACA 0002 for these low angles of attack but it is negative for NACA 0012 as can be seen from the last row of Figure 10. A peak in the amplitude spectrum is observed in Figure 9 at α = 9° for NACA 0002 and at α = 8° for NACA 0012 airfoil. These are the transition angles from the continuous vortex shedding mode (Mode 1) to alternating vortex shedding mode (Mode 2) for the vortex patterns at the downstream of the airfoil.
Amplitude spectrum of the lift coefficient at α = 1°, 7°,8° and 9° for NACA 0002 (left column) and NACA 0012 (right column). Instantaneous aerodynamic force coefficients at α = 1°, 7°,8° and 9° for NACA 0002 (left column) and NACA 0012 (right column).

The variation of the lateral and longitudinal spacing of the shed vortices and the vortex formation with the cylinder oscillation at Longitudinal and lateral spacings of vortex-street.
The formation of vortex patterns is investigated locally in terms of the vorticity magnitude and velocity magnitudes at the core of the vortices. The first two rows of Figure 12 show the mean and instantaneous vorticity field at α = 8° for NACA 0012. This angle is found to be a transition angle between the continuous vortex sheet mode and the alternating vortex shedding mode. The continuous vortex sheet is stretched and the instantaneous vortex shape is found to change from a droplet like shape to a more triangular shape and at the end to a stretched elliptic shape in the horizontal direction. The vortices grow as they travel downstream.
Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 8° for NACA 0012 (
The vorticity distribution shown in the third row of Figure 12 (as it was explained in Figure 6) results in oscillatory behavior which allows detection of the vortex core locations by the program developed.
The fourth and fifth rows in Figure 12 correspond to the instantaneous x-velocity (U/Uinf) and y-velocities (V/Uinf) at the core of the vortices, where Uinf is the free-stream velocity. On the same figures, the mean velocities at the corresponding vortex core locations are also presented by square symbols in order to visualize the mean induced velocities at the locations investigated. It is observed that the U/Uinf velocity is increasing towards downstream and an induced velocity (V/Uinf) in vertical direction is observed until x/c = 5.8 where the flow transforms from the continuous mode to the alternating mode. The induced velocity goes to zero towards downstream. At the downstream of x/c = 5.8, the mean velocity in the streamwise direction (U/Uinf) of the vortex cores converges to 92.5% of the free-stream velocity for α = 8°. The mean streamwise velocities of vortex cores decrease to 90% at α = 9° and to 89% at α = 12°. It is observed that the velocity of the vortex cores decreases until 20° incidence angle for NACA 0012 airfoil and reaches approximately 77% of the free-stream velocity at α = 20°. Bearman 13 has performed an experiment with 2D blunt trailing edge body at Re = 2.3 × 104 and Re = 4.1 × 104 based on base height of the body. He noted that in the initial part of the wake, the vortices are accelerating and the vortices were travelling between 88% and 89% of free-stream velocity in the stable region.
The sixth row in Figure 12 shows the longitudinal spacing (a/c) of the successive clockwise (blue circles) and counter-clockwise (red circles) vortices. The longitudinal vortex spacing is found by subtracting the horizontal coordinates of the two consecutive vortices having the same sign (rotating in the same direction). The x/c location in the horizontal axis of this figure is the mean x coordinate of the two successive vortices. In addition to the local longitudinal spacing of the vortices, the mean values of a/c ratios of the domain are also presented as a solid line on the same figure. The region where longitudinal spacing (a/c) is found constant is referred as stable region by Bearman 13 and Schaefer and Eskinazi. 14
For NACA 0012 at α = 8°, the longitudinal vortex spacing (a) is observed to be approximately constant and equal to 1.07 times chord length at the downstream of the airfoil after x/c > 5.3 (Figure 12), where stable region is reached. The location where the longitudinal vortex spacing becomes constant is found to be a transition location to the stable region for the shedding pattern of the vortex.
The last row in Figure 12 shows the longitudinal to lateral spacing (b/a) ratio of the vortex street. In this figure, the black solid line corresponds to b/a = 0.281 value obtained by von Karman. 13 Interestingly, it is observed that the transition from attached mode to alternating vortex shedding mode starts after x/c > 5.3 locations where b/a reaches 0.2817 as can be seen from Figure 12. The second transition from droplet like shape to triangular shape is observed at x/c > 10.6 where b/a ratio is found to be decreasing to 0.277 value and then increasing once more from this value.
At the downstream of x/c>11.7, b/a ratio is always bigger than 0.3 and reaches 0.36 at x/c = 17. Around this region, the vortex shape is observed to be triangular. Funakoshi 25 reported a value of 0.365 for triangular vortex blob observation and a value between 0.3<b/a/0.4 was reported by Tsuboi and Oshima. 27 Tsuboi and Oshima 27 and Funakoshi 25 investigated effect of the increase in longitudinal and lateral distance in the shape of vorticity region analytically. They defined merging and non-merging domains and shape of the vorticities in plane of (b, a). It should also be noted that if a vorticity region is initially circular, it will be deformed under shear flow.
The elliptical horizontally stretched shape of the vortices at the downstream is obtained when b/a ratio is increased to the values of 0.47 for clockwise (blue) vortex at far downstream and 0.478 for counter-clockwise (red) vortex at far downstream at the downstream of x/c = 17.
Both clockwise (blue) and counter-clockwise (red) vortices are observed to have almost constant longitudinal spacing in stable region. The lateral spacing increases at the downstream which also results in the deviation of the vortex street from the horizontal line towards up and down. The increase in spacing ratio (b/a) is a result of an increase in b (lateral spacing) since it is increasing with distance along the wake.
Mode 2: “alternating vortex shedding mode” (α = 9°–11°, α = 14°–20° for NACA 0002 and α = 8°–22° for NACA 0012)
Mode 2 is “
Mode 2 could also be investigated into two categories. The first category includes the angles of attack between 9° and 11° for NACA 0002 and angles between 9°and 12° for NACA 0012, where the vortex shedding pattern looks like a droplet shape close to the airfoil, then turns out to be triangle and afterwards the horizontally stretched elliptical shape. In the first category of Mode 2, the separation point moves upstream on the airfoil surface starting from the trailing edge
9
as the angle of attack increases. The trailing edge vortex (counter-clockwise) and a suction-side surface vortex (clockwise) are alternating by separating from the surface of the airfoil.
Instantaneous and Mean Vortex street of Mode 2 for NACA 0002. (9° ≤ α ≤ 11° and 14° ≤ α ≤ 20°). Instantaneous and Mean Vortex street of Mode 2 for NACA 0012. (9° ≤ α ≤ 22°). Amplitude spectrum of the lift coefficient at α = 9°, 14°, 17° and 20° for NACA 0002 (left column) and NACA 0012 (right column). Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 9° for NACA 0002. Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 14° for NACA 0002. Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 17° for NACA 0002.





The second category is where the shedding pattern close to the airfoil starts with a mushroom like vortex pair structure and returns to be an alternating vortex type at the downstream. The mushroom like vortex pair is formed of two counter-rotating vortices which are almost in circular shape initially as can be seen clearly from the instantaneous vorticity contours of α = 14° for NACA 0002 airfoil (Figure 13). The clockwise vortices (red contours) keep its circular shape more than the counter-clockwise vortices (blue contours). In addition, the alternating pattern turns out to be continuous vortex street at the downstream for the incidence angles between 15° and 19° for NACA 0002 (Figure 13) and for the incidence angels between 14° and 21° for NACA 0012 (Figure 14). It is known that co-rotating vortices have a tendency to merge, and counter-rotating vortices have a tendency to repel. Vorticity regions of Karman vortex can merge to continuous vorticity layers and circular shape can turn out to be horizontally elongated shape due to the combination of deformation and rotation in the flowfield.26–35
An important point is that, the vortices are observed to grow in size in the downstream for this specific Category II of Mode 2. For other modes, it is observed that the vortices grow before escaping from the airfoil surface or near the airfoil.
Panda and Zaman 28 have also observed mushroom like vortex structures experimentally for NACA 0012 airfoil at α = 15° for Re = 44,000.
The mean vorticity contours for this second category reveals also a secondary trace of the counter-clockwise vortices (blue contours) at the bottom row (see for example mean vorticity contours for α = 20° in Figure 13 or for α = 21° in Figure 14) due to the trace of counter-clockwise vortex which is stretching around the clockwise vortex (red contours).
It should be noted that, the mushroom type vortex pair structure is observed to be tilted upstream which results drag. For plunging airfoils, this mushroom type vortex pair is observed to be tilted downstream which resulted in a jet type wake and thrust generating condition.18–21,29 No momentum deficit in the wake (no drag) can be resulted also by plunging airfoils if the mushroom vortex pair is tilted vertically (not downwash or upwash) and when the alternating vortices are positioned exactly on a straight line. 18
Figure 15 shows the amplitude spectrum of the lift coefficient for angles of attack classified as Mode 2 in alternating vortex shedding mode. There is only one dominant peak in the amplitude spectrum which is the fundamental frequency for Mode 2 cases and the other amplitudes are decreasing for the harmonics of the fundamental frequency.
As can be investigated also from the variation of Strouhal number in Figure 5 (right column), the frequency corresponding to the maximum amplitude of the spectrum decreases continuously with the increase of the incidence angle in Mode 2. The flow is periodic with a single dominant frequency for Mode 2. With further increase of incidence angle, the period doubling occurs.
The initial growth of the disturbances is followed by the generation and growth of harmonics and a sub-harmonic of the fundamental frequency, which is indicative of nonlinear interactions between the disturbances.30,31 Subsequent merging of the roll-up vortices produces a second peak in the spectra centered at the sub-harmonic of the fundamental frequency 0.5 f0. 31
A transition occurs between the Modes 2 and 3 at α = 20° for NACA 0002 and α = 22° for NACA 0012. The longitudinal to lateral spacing (b/a) and vorticity distributions are almost some for both cases. The results for NACA 0002 at α = 20° are presented in Figure 19. It is observed that b/a is close to zero until x/c = 13–14 at α = 20° for NACA 0002 which concludes that the clockwise and counter -clockwise vortices are aligned individually along the horizontal axis for each mushroom type vortex pair. However, there is a deflection of the overall vortex pattern toward downwards.
Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 20° for NACA 0002.
Mode 3: “alternating vortex pair shedding mode” (α = 12°–13°, α = 21°–26°, α = 34°–41°for NACA 0002 and α = 23°–26°, α = 27°–28°, α = 31°–41° for NACA 0012)
Mode 3 is “ Vortex street for Mode 3, 12° ≤ α ≤ 13° for NACA 0002 (Category I). Amplitude spectrum of the lift coefficient at α = 11° (Mode 1), 12° (Mode 2), 13°(Mode 2) and 14° (Mode 1) for NACA 0002. Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 13° for NACA 0002. Vortex street for Mode 3 (Category II). (a) 21° ≤ α ≤ 26° for NACA 0002. (b) 23° ≤ α ≤ 28° for NACA 0012. Amplitude spectrum of the lift coefficient at α = 23°, 24°, 25°, and 26° for NACA 0002 (left column) and NACA 0012 (right column). Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 23° for NACA 0002. Vortex street for Mode 3 (Category III). (a) 34° ≤ α ≤ 41° for NACA 0002. (b) 31° ≤ α ≤ 41°for NACA 0012.






The transition locations between modes can also be observed from Strouhal number versus angle of attack curve in Figure 5. At α = 12°, the period doubling results in very low frequencies which shows a transition from Modes 2 to 3 for NACA 0002. There is also a low frequency in the amplitude spectrum at α = 25° for NACA 0002, a sudden decrease in Strouhal number reveals a transition from Modes 3 to 4.
Mode 3 is investigated in three sub-categories. The first category is the case of α = 12°–13° (Figure 20) for NACA 0002 airfoil for a transition from first category of Mode 2 to second category of Mode 2. From the amplitude spectrum, it is observed that the period at α = 11° (Mode 2) is quadrupled at α = 12° and α = 13° and then turns out to be a harmonic distribution at α = 14° (Mode 2) with a fundamental frequency. For this case, two vortex pair is visible in streamwise direction and alternatingly one vortex pair is in vertical ascent. This behaviour is more visible for α = 13° in Figure 20. The deformed couple seems to be responsible for this vertical ascent as shown in Figure 20. As the deformed vortex recovers its shape, the vortex structures becomes more aligned for α = 13°. This transition cannot be observed in NACA 0012 airfoil among the integer angles investigated.
It is also noted by Wu et al. 32 that the uneven spacing between the wake vortices implies that they form a row of short-spaced vortex couples rather than a Karman street. Each vortex couple, which in the co-moving frame of reference appears as a Kelvin oval, carries along a body of fluid.
Longitudinal spacings (a/c) between vortices in Mode 3 are not as uniform as in Mode 2 (sixth row in Figure 22). The spacing ratios are almost constant for vortex pairs in streamwise direction but the spacing ratio is close to zero for pairs uplifted. It should be noted that the counter-clockwise vortex (red contours at top) of the vortex pair in vertical ascent seems to rotate around the clockwise vortex as it goes towards downstream for α = 12°, 13° and 24°, 25°. Therefore, some of the b/a ratios are in negative sign since the clockwise vortex (blue vortex) is below the counter-clockwise vortex (red vortex) for the vortex couple in ascent. As can be seen from the fourth row of Figure 22, the vertical velocities of the ascending vortices (at approximately x/c = 4, x/c = 8 and x/c = 12.5) are high compared to other vortex pair in the flow field. However, the mean vertical velocities at the given locations are observed to be close to zero (square symbols in the fourth row of Figure 22).
A transition occurs between Modes 3 and 4 at α = 28° for NACA 0012. The single successive vortex at the top formed in Mode 4 (at α = 29°–30° for NACA 0012) becomes visible at this angle of attack. In comparison to α = 27°, it is seen how the mushroom like vortex pair is pushed towards the bottom row rather than being at the top row (Figure 23(b)). A similar transition phase is observed at α = 26° from Modes 3 to 4 for NACA 0002 airfoil. In the previous study 9 and also in the study of Khalid et al., 5 α = 27° was found to be chaos for NACA 0012 airfoil based on the amplitude spectrum of the lift coefficient. Similarly, α = 25° was found to be chaos for NACA 0002 airfoil at Re = 1000 (Figure 24). The decrease in the thickness affects the incidence angle where the wide range of frequencies are observed in the amplitude spectrum of the lift coefficient.
In Figure 25, the longitudinal spacing of counter-clockwise (red) vortices are found to be almost constant at α = 23° for NACA 0002 airfoil. The longitudinal spacings between clockwise (blue) vortices reveal also double wake pattern. It is observed that uplifted vortex pair has b/a ratios close to zero implying that the vortex pairs are aligned with horizontal at α = 23° for NACA 0002 airfoil.
At Category III of Mode 3, a constant Strouhal number (St = 0.233) is observed at incidence angles of 34° ≤ α ≤ 41° for NACA 0002 (see Figures 5 and 27). Similarly, Strouhal number is constant and equal to 0.247 at incidence angles of 31° ≤ α ≤ 39° for NACA 0012 and turns out to be equal to St = 0.233 for α = 40° and α = 41°. As the Strouhal number decreases, the period of the vortex shedding increases and the number of vortices at the domain investigated decreases also in the same time. For NACA 0002 at α = 37°, the spacing ratio of vortices is observed to become zero at the downstream which shows the horizontal alignment of the vortex pairs (Figure 28). The double vortex structure is observed to diminish its strength at the bottom row (see Figures 26(a) and 28 for mean vorticity contours at α = 37° for NACA 0002).
Amplitude spectrum of the lift coefficient at α = 34°, 37°, 38°, and 40° for NACA 0002 (left column) and NACA 0012 (right column). Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 37° for NACA 0002.

Mode 4: “alternating single vortex with vortex pair shedding mode” (α = 27°–33° for NACA 0002 and α = 29°–30° for NACA 0012
In Mode 4, during the shedding of the counter-rotating vortex pairs towards downstream, simultaneously small clockwise vortices form a row at the top (Figure 29). Same vortex formation as Mode 4 was visualized with NACA 0015 airfoil at Re = 5200 in plunge motion at 5° angle of attack with a plunge frequency of 6 Hz and amplitude of 10% of the chord length.
29
Vortex street for Mode 4. (a) 27° ≤ α ≤ 33° for NACA 0002. (b) 29° ≤ α ≤ 30°for NACA 0012.
There is a transition phase at α = 33° from Modes 4 to 3 for NACA 0002 (Figure 29(a)). This transition can also be visualized from the amplitude spectrum of the lift coefficient (Figure 30), where the magnitude of the second harmonic grows and the fundamental frequency decreases. Then at α = 34°, Mode 3 shows wide range of frequency as is shown in Figure 27. Mode 4 is visible in a larger incidence range (27° ≤ α ≤ 33°) for NACA 0002 airfoil compared to NACA 0012 airfoil, where just α = 29° and α = 30° exhibit this behavior in the wake pattern (Figure 29(b)). The longitudinal spacing of the vortices is decreasing from 27° to 33° for NACA 0002 airfoil. The average streamwise velocities of the vortex cores are observed to decrease to 70% of the free-stream velocity at α = 31° (Figure 31).
Amplitude spectrum of the lift coefficient at α = 27°, 29°, 31,° and 33° for NACA 0002. Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 31° for NACA 0002.

It should be noted that the streamlines of the mean velocity field around the NACA 0012 airfoil revealed an abrupt decrease in the size of the primary clockwise leading edge vortex on the upper surface of the airfoil when increasing the incidence angle from α = 30° to α = 31°. 9
Mode 5: “bluff body vortex shedding mode” (α = 50°–90° for NACA 0002 and NACA 0012)
Mode 5 is “bluff body vortex shedding mode”. The vortex shedding at large incidence angles behaves like a bluff body (Figure 32). The vortex pattern shows vortices which are amalgamated with time before it is shed to the wake. The vortices grow very slowly. The vortex separated from the trailing edge is slightly bigger than that separating from leading edge due to the shape difference between leading edge (rounded) and trailing edge (flat).
21
The amplitude spectrums reveals a very low frequency of vortex shedding with a strong fundament frequency amplitude at this mode (Figure 33). The longitudinal spacing is found to be very high for Mode 5 (Figure 34). The streamwise component of the vorticity cores is decreased to approximately 50% of the free-stream velocity at α = 70° for NACA0012.
Vortex street for Mode 5. (a) 50° ≤ α ≤ 90°for NACA 0002 (b) 50° ≤ α ≤ 90°for NACA 0012. Amplitude spectrum of the lift coefficient at α = 50°, 70°, 80°, and 90° for NACA 0002 (left column) and NACA 0012 (right column). Mean and instantaneous vortex streets, vorticity distribution, vortex core streamwise, and vertical velocities, vortex core longitudinal spacing (a/c) and b/a ratio at α = 70° for NACA 0012.


Figure 35 shows average longitudinal spacings of clockwise and counter-clockwise at different incidence angles. The error bars in the same figure vortices show also the maximum and minimum longitudinal spacings obtained from the flowfield. As a result, Figure 35 also shows the interval of longitudinal spacings obtained even at double wake structure cases. There is no discrepancy between the longitudinal spacings of clockwise and counter-clockwise vortices for Mode 2 and first category of Mode 3 (Figure 35). Second and third category of Modes 3 and 4 show a wide range of longitudinal spacings in the field due to the double wake structures formed. It is important to note that Mode 3 looks like a transition phase between different modes.
Mean longitudinal spacings (a) NACA 0002. (b) NACA 0012.
Conclusion
Effect of incidence angle and airfoil thickness on vortex pattern at the wake of symmetric airfoils are investigated in the current study. The instantaneous and mean flow fields are obtained for 2% thick (NACA 0002) and 12 % thick (NACA 0012) airfoils.
Five different vortex pattern modes are identified according to the shape of the vortices, amplitude spectrum of the lift coefficient, instantaneous and mean velocity fields and longitudinal and lateral spacings between vortex cores at the wake of the airfoil.
For each mode, the location where each individual longitudinal vortex spacing becomes constant is found to be a transition location to the stable region for the shedding pattern of the vortex. It is observed that the lateral spacing increases at the downstream which also results the deviation of the vortex street from the horizontal line towards up and down. The increase in spacing ratio (b/a) is a result of an increase in the lateral spacing due to the wake deflection.
The flow is found to be periodic with a single dominant frequency for Mode 2 and with the further increase of the angle of attack, the period doubling occurs.
It is interestingly noted that Mode 3 looks like a transition phase between the passages at different modes. Period quadrupling and chaos condition occur during Mode 3 phase of the vortex pattern.
A transition occurs from Modes 2 to 3 at α = 20° for NACA 0002 and α = 22° for NACA 0012, where the longitudinal to lateral spacing (b/a) are very close to zero which states that each pair of counter-rotating vortices are aligned horizontally.
The vortices are observed to grow in size in the downstream for Category II of Mode 2. For other modes, it is observed that the vortices grow before escaping from the airfoil surface or near the airfoil.
Previously, α = 27° was found to be chaos for NACA 0012 airfoil based on the amplitude spectrum of the lift coefficient. 9 Similarly, α = 25° was found to be chaos for NACA 0002 airfoil at Re = 1000. The decrease in the thickness, diminishes also the angle of attack of this transition location where wide range of frequencies are observed in the amplitude spectrum.
At Category III of Mode 3, a constant Strouhal number (St = 0.233) is observed at incidence angles of 34° ≤ α ≤ 41° for NACA 0002. Similarly, Strouhal number is constant and equal to 0.247 at incidence angles of 31° ≤ α ≤ 39° for NACA 0012 and turns out to be equal to St = 0.233 for α = 40° and α = 41°.
Current study gives a collection of the flow patterns obtained at the wake of symmetric airfoils at Re = 1000 and their classification according to their shape, frequency spectrums, and vortex core spacings.
Footnotes
Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work is supported TUBITAK 213M327 project. The financial support from Turkish Academy of Science-Young Scientific Award (TUBA GEBIP) is gratefully acknowledged by the author.
