Effect of NACA0012 Airfoil Pitching Oscillation on Flow Past a Cylinder

: The ﬂow past a cylinder is a classical problem in ﬂow physics. In a certain range of Reynolds number, there will be Karman vortex street phenomenon in the wake of a cylinder, which will greatly increase the pressure drag of the cylinder. By controlling the vortex shedding phenomenon, drag reduction of the cylinder could be effectively realized. In this paper, a NACA0012 airfoil with pitching oscillation is placed downstream of the cylinder. Based on the tight coupling method, kinematics equations and Navier–Stokes equations in the arbitrary Lagrangian–Eulerian form are solved. Firstly, the effect of airfoil oscillation period and the distance between airfoil leading edge and cylinder center ( x / D ) are studied respectively, especially considering the aspects of vortex shedding and drag reduction effect. Besides, the vortex interaction in the ﬂow ﬁeld around the airfoil and cylinder is analyzed in detail. It is found that the NACA0012 airfoil with pitching oscillation can change the period of vortex shedding. Moreover, it can also increase the drag reduction rate to as high as 50.5%, which presents a certain application prospect in the engineering drag reduction ﬁeld, e.g., for launch vehicles, masts, submarine etc.


Introduction
Flow past a cylinder is a common phenomenon of flow physics in nature and engineering practices, e.g., in launch vehicles [1,2], ship masts [3,4], submarine pipelines [5,6], and fluidics at the micro and nano scale [7,8]. Within the range of Reynolds number common in engineering practices, an alternately shedding Karman vortex street can be generated in the wake of a cylinder. The alternately shedding vortex creates a low-pressure area at the trailing edge of the cylinder, leading to a sharp increase in the pressure difference between the front and back walls of the cylinder. Pressure drag contributes most sources of resistance to the cylinder. On this basis, drag reduction can be effectively achieved in the cylinder via controlling the vortex shedding phenomenon.
Recently, achieving drag reduction through flow control has aroused extensive attention. Choi et al. [9] divided the control of flow past a cylinder into boundary layer control and wake adjustment. As mentioned in the literature [10], flow past a cylinder controlled using electromagnetic force can delay the separation of the turbulent boundary layer, thus decreasing resistance and suppressing lift pulsation. Jeon [11] and Kimura [12] made separation bubbles appear on the cylinder surface by means of blow-suction disturbance and adding cross grooves, so that drag reduction can be achieved with the separated-andreattaching flow of separation bubbles. Boundary layer control presents requirements regarding the object shape and Reynolds number. By comparison, it is simpler and more efficient to achieve drag reduction through controlling the cylindrical wake.
By changing the flow state of the wake, the pressure distribution on the cylinder surface in the subsonic flow can be adjusted, thereby markedly lowering front-back pressure drag of the object. Moreover, adding a circumferential flow plate is also an effective flow control method. As Figure 1a shows, many scholars [13][14][15][16][17][18][19] added a spoiler in the wake of a cylinder, and found that the drag reduction effect and the vortex shedding frequency are strongly associated with the spoiler length. Anderson et al. [17] discovered that the drag coefficient can be reduced by roughly 30% at most when the length of the circumferential flow plate is equivalent to the cylinder diameter according to numerical simulation. Besides, as shown in Figure 1b, Roshko et al. [20] carried out wind tunnel tests of cases existing a gap between the cylinder and spoiler. They concluded that the drag coefficient and vortex shedding frequency could achieve minimum values when the distance (between the spoiler leading edge and the cylinder trailing edge) is 2.3 times the cylinder diameter. Likewise, Hwang [21] found the same law during numerical simulation. Ozono [22] and Bearman et al. [23] also achieved the purpose of preventing vortex shedding by conducting a wind tunnel test on asymmetrically arranged cylinders and spoilers. In addition, Rathakrishnan et al. [24] revealed that h/4 (h is the height of the rectangular cylinder) was the optimum distance between the cylinder and the plate, and the best drag reduction effect could be achieved when the spoiler was placed horizontally based on their experimental study of the rectangular cylinder.
Energies 2021, 14, x FOR PEER REVIEW 2 of 17 By changing the flow state of the wake, the pressure distribution on the cylinder surface in the subsonic flow can be adjusted, thereby markedly lowering front-back pressure drag of the object. Moreover, adding a circumferential flow plate is also an effective flow control method. As Figure 1a shows, many scholars [13][14][15][16][17][18][19] added a spoiler in the wake of a cylinder, and found that the drag reduction effect and the vortex shedding frequency are strongly associated with the spoiler length. Anderson et al. [17] discovered that the drag coefficient can be reduced by roughly 30% at most when the length of the circumferential flow plate is equivalent to the cylinder diameter according to numerical simulation. Besides, as shown in Figure 1b, Roshko et al. [20] carried out wind tunnel tests of cases existing a gap between the cylinder and spoiler. They concluded that the drag coefficient and vortex shedding frequency could achieve minimum values when the distance (between the spoiler leading edge and the cylinder trailing edge) is 2.3 times the cylinder diameter. Likewise, Hwang [21] found the same law during numerical simulation. Ozono [22] and Bearman et al. [23] also achieved the purpose of preventing vortex shedding by conducting a wind tunnel test on asymmetrically arranged cylinders and spoilers. In addition, Rathakrishnan et al. [24] revealed that h/4 (h is the height of the rectangular cylinder) was the optimum distance between the cylinder and the plate, and the best drag reduction effect could be achieved when the spoiler was placed horizontally based on their experimental study of the rectangular cylinder. In addition to the spoiler, methods of placing small-sized cylinders and symmetric airfoils similar to the spoiler in the cylindrical wake [25,26] have also aroused the attention of scholars in recent years. Strykowski et al. [27] achieved the control of vortex shedding in the wake of the cylinder through placing a secondary small cylinder behind the cylinder, and disclosed that variations in the local flow field in the wake might affect the global characteristics. According to the literature [28], the addition of airfoil might affect the aerodynamic force and vortex street shape. Zhang et al. [29] found that not only were characteristics related the cylinder affected, but the lift-drag characteristics of the static airfoil in the wake of the cylinder also varied significantly. Furthermore, the maximum lift coefficient of the airfoil can be increased by 34% by adjusting the longitudinal distance between the airfoil leading edge and the center line of the wake.
The dynamic airfoil has a more significant impact on the wake of a cylinder in the cylindrical-airfoil tandem configuration apart from the static airfoil. Gopalkrishnan et al. [30] detected three varied flow modes in the airfoil wake from the experimental study on the phase difference between the airfoil motion and the vortex street. After the numerical simulation, Streitlien et al. [31] discovered that the cylinder-airfoil tandem configuration presented the least drag when the flapping phase was consistent with the vortex phase in the case of the airfoil carrying out parallel and pitching oscillations. Chao et al. [32] performed numerical simulation on a swing airfoil placed on the center line of two parallel In addition to the spoiler, methods of placing small-sized cylinders and symmetric airfoils similar to the spoiler in the cylindrical wake [25,26] have also aroused the attention of scholars in recent years. Strykowski et al. [27] achieved the control of vortex shedding in the wake of the cylinder through placing a secondary small cylinder behind the cylinder, and disclosed that variations in the local flow field in the wake might affect the global characteristics. According to the literature [28], the addition of airfoil might affect the aerodynamic force and vortex street shape. Zhang et al. [29] found that not only were characteristics related the cylinder affected, but the lift-drag characteristics of the static airfoil in the wake of the cylinder also varied significantly. Furthermore, the maximum lift coefficient of the airfoil can be increased by 34% by adjusting the longitudinal distance between the airfoil leading edge and the center line of the wake.
The dynamic airfoil has a more significant impact on the wake of a cylinder in the cylindrical-airfoil tandem configuration apart from the static airfoil. Gopalkrishnan et al. [30] detected three varied flow modes in the airfoil wake from the experimental study on the phase difference between the airfoil motion and the vortex street. After the numerical simulation, Streitlien et al. [31] discovered that the cylinder-airfoil tandem configuration presented the least drag when the flapping phase was consistent with the vortex phase in the case of the airfoil carrying out parallel and pitching oscillations. Chao et al. [32] performed numerical simulation on a swing airfoil placed on the center line of two parallel semi-cylinders and discovered that the propulsion efficiency was affected by lateral gap of semi-cylinders, the lateral distance between the semi-cylinder and the airfoil, as well as the swing frequency of the airfoil. Xiao et al. [33] added an airfoil-like foil behind the cylinder and made the foil swing in accordance with the "subcarangiform" motion model, and through numerical simulation, they found that the swing frequency was a vital factor affecting the lift-drag performance of the cylinder.
The abovementioned studies mainly focus on the aerodynamic characteristics and propulsion efficiency of the airfoil. However, studies on the aerodynamic performance of the cylinder and the drag reduction effect of combined configuration remain relatively scarce, and there is no detailed analysis on the vortex interaction in the flow field. Therefore, the supplement of the research on cylinder-airfoil wake is of great significance, and this paper adopts numerical simulations to conduct a detailed analysis about that. A NACA0012 airfoil, one of classical symmetrical wing models of aircraft (the airfoil profile is defined in Equation (1), x ∈ [0, 1]), with pitching oscillation of small attack angle, is placed in the wake of a cylinder. With the dynamic airfoil, a novel drag reduction method is introduced, which is conducive to the engineering drag reduction. Then variations in the vortex shedding frequency and the drag coefficient of tandem configuration are studied through changing the oscillation frequency of the airfoil as well as the distance between airfoil leading edge and cylinder center (x/D). On this basis, the interaction between the cylinder and the vortex around the airfoil is also analyzed. The mechanisms found in this paper may be helpful for the future engineering applications, e.g., in launch vehicles, ship masts, submarine pipelines, etc.
This paper is organized as follows. To begin with, physical models and numerical methods are introduced in Section 2 together with the independent analysis of mesh and verification of numerical methods. Then, the influences of the airfoil oscillation period on the vortex shedding and drag of the cylinder are analyzed in Section 3. Moreover, here, we analyze the situation whereby the distance between airfoil leading edge and cylinder center is different and investigate the vortex interaction in the flow field. Finally, conclusions are made in Section 4.

Physical Model
This paper studied the two-dimensional cylinder-airfoil tandem configuration, and the computational domain is shown in Figure 2. Specifically, the cylinder diameter D is 100mm, and the airfoil chord length is equal to the cylinder diameter (l/D = 1.0). Besides, NACA0012 airfoil chord along the x-axis and the distance between airfoil leading edge and cylinder center (x/D) can be adjusted. Furthermore, the radius of the computational domain is set as 20D, so as to remove the influence of boundary conditions of subsonic velocity. The NACA0012 airfoil located downstream oscillates along the 1/4 chord point in pitch, and its attack angle is changed as in Equation (2): where α is the average attack angle and α m is the oscillation amplitude.

Numerical Methods
The NACA0012 airfoil with forced pitching oscillation is related to the motion of mesh boundary. For accurately describing the flow of moving boundary, a coordinate system is established based on the arbitrary Lagrangian-Eulerian (ALE) method. In the ALE coordinate system, the governing equation in the integral form can be expressed as Equation (3): where U represents the conservation vector; F and F ν denote the non-viscous flux and viscous flux, respectively; νg symbolizes the velocity of moving boundary; dS and dV stand for the surface area and volume of the control volume; and n points to the direction of exterior normal of the control volume surface. The calculation in this paper is based on the finite volume method, which is more flexible and fast-solved than other discrete methods, and this method is applied by many popular CFD software programs, such as FLUENT, CFX, etc. Besides, the one-equation SA turbulence model [34], suitable for simulating the turbulence on unsteady dynamic motion characteristics at low speed, is adopted for flow simulation, and the computation quantity of the model is lower compared to other turbulent models. Besides, the non-viscous and viscous fluxes are differentially split in the Roe format [35]. In addition, the dualtime-step and LU-SGS implicit marching method [36] are used for the time term by setting the number of sub-iteration steps as 80.
The Navier-Stokes (NS) equations in the ALE coordinate system interact with the dynamics/ kinematics equations, controlling the motion of the airfoil by tight coupling. Regarding the classical problem of flow past a cylinder, flow in the case of 300 ≤ Re < 2 × 10 5 is normally considered as in the sub-critical zone. In this study, the Reynolds number is set as 5000. Under this condition, the laminar flow in the upstream boundary layer of the cylinder has converted to turbulent flow, together with the wake transforming into turbulent vortex street, which can promote the separation of the airfoil boundary layer for better flow observation, and this can be the inflow condition for the downstream airfoil. Specific parameters of the flow field are presented in Table 1.

Numerical Methods
The NACA0012 airfoil with forced pitching oscillation is related to the motion of mesh boundary. For accurately describing the flow of moving boundary, a coordinate system is established based on the arbitrary Lagrangian-Eulerian (ALE) method. In the ALE coordinate system, the governing equation in the integral form can be expressed as Equation (3): where U represents the conservation vector; F and F ν denote the non-viscous flux and viscous flux, respectively; ν g symbolizes the velocity of moving boundary; dS and dV stand for the surface area and volume of the control volume; and n points to the direction of exterior normal of the control volume surface. The calculation in this paper is based on the finite volume method, which is more flexible and fast-solved than other discrete methods, and this method is applied by many popular CFD software programs, such as FLUENT, CFX, etc. Besides, the one-equation SA turbulence model [34], suitable for simulating the turbulence on unsteady dynamic motion characteristics at low speed, is adopted for flow simulation, and the computation quantity of the model is lower compared to other turbulent models. Besides, the non-viscous and viscous fluxes are differentially split in the Roe format [35]. In addition, the dual-time-step and LU-SGS implicit marching method [36] are used for the time term by setting the number of sub-iteration steps as 80.
The Navier-Stokes (NS) equations in the ALE coordinate system interact with the dynamics/ kinematics equations, controlling the motion of the airfoil by tight coupling. The tightly coupled marching-iteration strategy is presented in Figure 3, where CFD (computational fluid dynamics) is the solution of NS equations and RBD (rigid body dynamics) is the solution of flight dynamic and kinematic equations. These methods have been integrated into our unsteady RANS solver. For more details, please refer to the work of Zhang and Wang [37] and Chang et al. [38]. The tightly coupled marching-iteration strategy is presented in Figure 3, where CFD (computational fluid dynamics) is the solution of NS equations and RBD (rigid body dynamics) is the solution of flight dynamic and kinematic equations. These methods have been integrated into our unsteady RANS solver. For more details, please refer to the work of Zhang and Wang [37] and Chang [38] et al.

Mesh Generation
The dynamic overlapping grid technique [39] is utilized to handle the moving boundary problem during computation. Essentially, since there might be overlapping and nesting areas between mesh blocks generated separately by varied objects, necessary pretreatment should be performed on mesh blocks before computation. Then, a transitive relationship should be established on the boundary of the overlapping area, so that the boundary information can be updated through exchanging data interpolation during computation. On this basis, the solution of the total flow field can be obtained when time marching has been conducted.
The grids are composed of a cylinder mesh, a mesh near the airfoil, and a background mesh. Among them, the background mesh covers the whole computational domain. Furthermore, the mesh can be divided into structured and unstructured hybrid meshes in accordance with Figure 4a,b. Specifically, quadrilateral mesh and triangular mesh are adopted for the body-fitted mesh near the cylinder as well as airfoil walls, respectively. To seize the viscous boundary layer more accurately, the body-fitted quadrilateral mesh is equipped with more than 20 layers. Moreover, the mesh near the object is properly encrypted in the background mesh, so that the flow field around the object can be displayed more subtly. Besides, the mesh is processed into blocks before computation for conducting parallel calculation and improving the computational efficiency. More precisely speaking, cylinder, airfoil, and background meshes are partitioned into 16 blocks, respectively, as shown in Figure 4c.

Mesh Generation
The dynamic overlapping grid technique [39] is utilized to handle the moving boundary problem during computation. Essentially, since there might be overlapping and nesting areas between mesh blocks generated separately by varied objects, necessary pre-treatment should be performed on mesh blocks before computation. Then, a transitive relationship should be established on the boundary of the overlapping area, so that the boundary information can be updated through exchanging data interpolation during computation. On this basis, the solution of the total flow field can be obtained when time marching has been conducted.
The grids are composed of a cylinder mesh, a mesh near the airfoil, and a background mesh. Among them, the background mesh covers the whole computational domain. Furthermore, the mesh can be divided into structured and unstructured hybrid meshes in accordance with Figure 4a,b. Specifically, quadrilateral mesh and triangular mesh are adopted for the body-fitted mesh near the cylinder as well as airfoil walls, respectively. To seize the viscous boundary layer more accurately, the body-fitted quadrilateral mesh is equipped with more than 20 layers. Moreover, the mesh near the object is properly encrypted in the background mesh, so that the flow field around the object can be displayed more subtly. Besides, the mesh is processed into blocks before computation for conducting parallel calculation and improving the computational efficiency. More precisely speaking, cylinder, airfoil, and background meshes are partitioned into 16 blocks, respectively, as shown in Figure 4c.

Verification of Numerical Method
When carrying out the verification of mesh independence, three sets of meshes are adopted for simulating the shape of a single cylinder. The first layer height of mesh and mesh volume of the corresponding wall are displayed in Table 2. With the mesh being encrypted, the vortex shedding period in the wake of a cylinder is unchanged, presenting certain convergence. To improve the computational efficiency, the medium-density mesh is applied in this paper on the premise of safeguarding the accuracy of results.
To seize the viscous boundary layer more accurately, the body-fitted quadrilateral mesh is equipped with more than 20 layers. Moreover, the mesh near the object is properly encrypted in the background mesh, so that the flow field around the object can be displayed more subtly. Besides, the mesh is processed into blocks before computation for conducting parallel calculation and improving the computational efficiency. More precisely speaking, cylinder, airfoil, and background meshes are partitioned into 16 blocks, respectively, as shown in Figure 4c.    [40], are adopted for verifying the dynamic computational performance of the numerical method in this paper.
Firstly, the steady flow field at the average attack angle is calculated for the initial field of unsteady simulation. The lift coefficient varying with the attack angle is shown in Figure 5, incorporating the experimental results of Landon [40] and numerical results of other researchers [41]. The computational result of this paper agrees well with the numerical simulation results in the literature, and is basically identical with the distribution of the experimental results. This proves that the numerical method used in this paper can accurately simulate the pitching oscillation of the NACA0012 airfoil at a small attack angle.

Verification of Numerical Method
When carrying out the verification of mesh independence, three sets of meshes are adopted for simulating the shape of a single cylinder. The first layer height of mesh and mesh volume of the corresponding wall are displayed in Table 2. With the mesh being encrypted, the vortex shedding period in the wake of a cylinder is unchanged, presenting certain convergence. To improve the computational efficiency, the medium-density mesh is applied in this paper on the premise of safeguarding the accuracy of results.  [40], are adopted for verifying the dynamic computational performance of the numerical method in this paper.
Firstly, the steady flow field at the average attack angle is calculated for the initial field of unsteady simulation. The lift coefficient varying with the attack angle is shown in Figure 5, incorporating the experimental results of Landon [40] and numerical results of other researchers [41]. The computational result of this paper agrees well with the numerical simulation results in the literature, and is basically identical with the distribution of the experimental results. This proves that the numerical method used in this paper can accurately simulate the pitching oscillation of the NACA0012 airfoil at a small attack angle.

Effect of T/T 0
When the Reynolds number is between 500 and 2 × 10 5 , the Strouhal number (St) of flow past a cylinder will basically remain unchanged at 0.2. The vortex shedding period T 0 can be calculated by Equation (4), and the airfoil oscillation period is T.
where D represents the cylinder diameter and U ∞ notes the velocity of free inflow. This section discusses the influence of T/T 0 (the ratio of the airfoil oscillation period to the vortex shedding period of a single cylinder) on the tandem configuration. Computational results and analysis are presented in the following two subsections.

Influence on Vortex Shedding Frequency
Aerodynamic characteristics of the airfoil and cylinder are changed periodically under the combined effect of the airfoil oscillation and the cylindrical wake. Concretely, the aerodynamic coefficients' frequencies of the cylinder and airfoil vary with the increase of T/T 0 , as shown in Figure 6a. When T/T 0 is less than roughly 0.2, the frequency oscillation modes of cylinder and airfoil also change, and the smaller the T/T 0 , the greater the difference. Besides, when T/T 0 is greater than the critical value, the frequencies of both are coupled, which is called "aerodynamic coupling" in this paper. When the Reynolds number is between 500 and 2 × 10 5 , the Strouhal number (St) of flow past a cylinder will basically remain unchanged at 0.2. The vortex shedding period T0 can be calculated by Equation (4), and the airfoil oscillation period is T.
where D represents the cylinder diameter and U∞ notes the velocity of free inflow. This section discusses the influence of T/T0 (the ratio of the airfoil oscillation period to the vortex shedding period of a single cylinder) on the tandem configuration. Computational results and analysis are presented in the following two subsections.

Influence on Vortex Shedding Frequency
Aerodynamic characteristics of the airfoil and cylinder are changed periodically under the combined effect of the airfoil oscillation and the cylindrical wake. Concretely, the aerodynamic coefficients' frequencies of the cylinder and airfoil vary with the increase of T/T0, as shown in Figure 6a. When T/T0 is less than roughly 0.2, the frequency oscillation modes of cylinder and airfoil also change, and the smaller the T/T0, the greater the difference. Besides, when T/T0 is greater than the critical value, the frequencies of both are coupled, which is called "aerodynamic coupling" in this paper. The coupling frequency of lift coefficient varying with various airfoils and x/D is shown in Figure 6b, comparing with the variation frequency of the static cylinder lift coefficient. As T/T0 increases, the coupling frequency will decrease. The coupling frequency and static results have little difference under high-frequency oscillation of the airfoil, whereas the former approaches the airfoil oscillation frequency at low-frequency oscillation. It is evident that spectral characteristics in the wake of a cylinder have a greater influence on the airfoil low-frequency oscillation.
To investigate critical factors of aerodynamic coupling, fast Fourier analysis (FFT) is performed on the instantaneous lift coefficients of the cylinder and airfoil. Furthermore, the amplitudes corresponding to the fundamental frequencies of both are calculated, the ratio of which is presented in Figure 7. The results show that the ratio gradually decreases with the increase of T/T0, and a turning point appears when T/T0 = 0.3, that is, the variation The coupling frequency of lift coefficient varying with various airfoils and x/D is shown in Figure 6b, comparing with the variation frequency of the static cylinder lift coefficient. As T/T 0 increases, the coupling frequency will decrease. The coupling frequency and static results have little difference under high-frequency oscillation of the airfoil, whereas the former approaches the airfoil oscillation frequency at low-frequency oscillation. It is evident that spectral characteristics in the wake of a cylinder have a greater influence on the airfoil low-frequency oscillation.
To investigate critical factors of aerodynamic coupling, fast Fourier analysis (FFT) is performed on the instantaneous lift coefficients of the cylinder and airfoil. Furthermore, the amplitudes corresponding to the fundamental frequencies of both are calculated, the ratio of which is presented in Figure 7. The results show that the ratio gradually decreases with the increase of T/T 0 , and a turning point appears when T/T 0 = 0.3, that is, the variation trend of the ratio slackens in the presence of aerodynamic coupling. By inference, the amplitude ratios of the cylinder and airfoil aerodynamic coefficients corresponding to the fundamental frequency represent a critical factor affecting aerodynamic coupling. When the ratio is less than a critical value (it is about 1.5 when x/D = 1.5; note that the critical value is varied along with the variation of x/D), aerodynamic coupling will occur in the corresponding flow field. trend of the ratio slackens in the presence of aerodynamic coupling. By inference, the amplitude ratios of the cylinder and airfoil aerodynamic coefficients corresponding to the fundamental frequency represent a critical factor affecting aerodynamic coupling. When the ratio is less than a critical value (it is about 1.5 when x/D = 1.5; note that the critical value is varied along with the variation of x/D), aerodynamic coupling will occur in the corresponding flow field.

Influence on Drag Reduction Characteristics
When T/T0 > 0.2, aerodynamic coupling is observed in the flow field. At this time, airfoil oscillation has a significant impact on the wake of a cylinder. Subsequently, the drag reduction effect of tandem configuration (T/T0 ∈ [0.4, 10]) is analyzed, and the mean drag coefficient value of a single cylinder is computed (1.6087). Figure 8 indicates that the relationship between the drag reduction rate of tandem configuration and T/T0 is related to x/D (the distance between the airfoil leading edge and the cylinder center). Moreover, two completely different tendencies before x/D = 1.5 and after x/D = 1.75 are presented in the curve of drag reduction. To be specific, when x/D ≤ 1.5, the drag reduction rate varies slightly, which can reach above 45% within the studied range. When x/D ≥ 1.75, the drag reduction rate is relatively small, reaching 30% at most. As T/T0 increases, the drag reduction rate increases at first, and then decreases with the maximum value taken around T/T0 = 2. The reasons are shown in Section 3.2.
The above variations are explained below in combination with flow fields when x/D = 1.0 and x/D = 2.5, respectively. As can be observed from Figure 9, there is a relatively short distance between the cylinder and the airfoil when x/D = 1.0. Thus, the airfoil is surrounded by the wake vortex of the cylinder and no Karman vortex street is formed between both. Moreover, with the strong wake vortex of the cylinder, the airfoil oscillation has a relatively small effect on the flow field. Therefore, the drag coefficient varies little when T/T0 changes.
Pressure contour and pressure coefficient of the cylinder wall at x/D = 2.5 are displayed in Figure 10 (φ = 0 represents the cylindrical wall when x = −0.5, and φ increases clockwise). A low-pressure area can be observed at the trailing edge of the cylinder (partially y < 0). When T/T0 ≤ 2.0, the increase in the oscillation period will suppress the development of the low-pressure area. The pressure difference between the first half (φ ∈ [0°, 90°]∪[270°, 360°]) and second half of the cylinder (φ ∈ [90°, 270°]) gets lowered gradually. On the contrary, when T/T0 ≥ 2.0, the increase in the oscillation period will promote

Influence on Drag Reduction Characteristics
When T/T 0 > 0.2, aerodynamic coupling is observed in the flow field. At this time, airfoil oscillation has a significant impact on the wake of a cylinder. Subsequently, the drag reduction effect of tandem configuration (T/T 0 ∈ [0.4, 10]) is analyzed, and the mean drag coefficient value of a single cylinder is computed (1.6087). Figure 8 indicates that the relationship between the drag reduction rate of tandem configuration and T/T 0 is related to x/D (the distance between the airfoil leading edge and the cylinder center). Moreover, two completely different tendencies before x/D = 1.5 and after x/D = 1.75 are presented in the curve of drag reduction. To be specific, when x/D ≤ 1.5, the drag reduction rate varies slightly, which can reach above 45% within the studied range. When x/D ≥ 1.75, the drag reduction rate is relatively small, reaching 30% at most. As T/T 0 increases, the drag reduction rate increases at first, and then decreases with the maximum value taken around T/T 0 = 2. The reasons are shown in Section 3.2.
The above variations are explained below in combination with flow fields when x/D = 1.0 and x/D = 2.5, respectively. As can be observed from Figure 9, there is a relatively short distance between the cylinder and the airfoil when x/D = 1.0. Thus, the airfoil is surrounded by the wake vortex of the cylinder and no Karman vortex street is formed between both. Moreover, with the strong wake vortex of the cylinder, the airfoil oscillation has a relatively small effect on the flow field. Therefore, the drag coefficient varies little when T/T 0 changes. the development of the low-pressure area. In addition, together with the rising pressure difference between the front and back surfaces of the cylinder, the drag reduction rate will decline.  the development of the low-pressure area. In addition, together with the rising pressure difference between the front and back surfaces of the cylinder, the drag reduction rate will decline. Pressure contour and pressure coefficient of the cylinder wall at x/D = 2.5 are displayed in Figure 10 (ϕ = 0 represents the cylindrical wall when x = −0.5, and ϕ increases clockwise). A low-pressure area can be observed at the trailing edge of the cylinder (partially y < 0). When T/T 0 ≤ 2.0, the increase in the oscillation period will suppress the development of the low-pressure area. The pressure difference between the first half (ϕ ∈ [0 • , 90 • ]∪[270 • , 360 • ]) and second half of the cylinder (ϕ ∈ [90 • , 270 • ]) gets lowered gradually. On the contrary, when T/T 0 ≥ 2.0, the increase in the oscillation period will promote the development of the low-pressure area. In addition, together with the rising pressure difference between the front and back surfaces of the cylinder, the drag reduction rate will decline.

Effect of x/D
Through analyzing the influence of T/T 0 , it is apparent that x/D (the distance between the airfoil leading edge and cylinder center) has a significant influence on the flow characteristics.

Effect of x/D
Through analyzing the influence of T/T0, it is apparent that x/D (the distance between the airfoil leading edge and cylinder center) has a significant influence on the flow characteristics.

Influence on the Flow Mechanism
According to computational results listed in Section 3.1.2, drag characteristics of tandem configuration before x/D = 1.5 and after x/D = 1.75 are completely different. Hence, evolutions of the flow field when x/D = 1.0 and x/D = 1.75 are selected for following discussion.
The vorticity contour when x/D = 1.0 is displayed in Figure 11. At this moment, there is no complete Karman vortex street between the cylinder and the airfoil. Specifically, t0-t5 represent different evolution stages of the flow field respectively.
The wake vortex of the cylinder in the flow field can be called the "primary vortex" for its dominating effect. Taking the local flow field of y < 0 as an example, when the primary vortex meets the NACA0012 airfoil, since the vortex length is far greater than the airfoil chord length, the latter is wrapped by the primary vortex, as shown in Figure 11b. When the primary vortex marches to the middle of the airfoil, a low-pressure area will form, causing the induced vortices generated in the trailing and leading edges of the airfoil successively, as shown in Figure 11c. In the meantime, the airfoil, like a vortex separator, can divide the primary vortex into two parts. As can be observed from Figure 11d, a part of vortex is fused with the induced vortex generated by the upper airfoil in the trailing edge, while the other part is fused with the induced vortex developed at the lower airfoil in the leading edge. After that, the new primary vortex flows past the upper airfoil, as shown in Figure 11e. The direction of the primary vortex at the upper airfoil is opposite to the residual vortex at the lower airfoil, and the primary vortex at the upper airfoil expands and pushes the vortex at the lower airfoil to move along the opposite direction of inflow and collide with the cylinder surface, as shown in Figure 11f. Next, the flow separation at the y < 0 section is aggravated. Hence, the vortex that is pushed to the surrounding wall is accumulated, forming a new primary vortex to repeat the above process.
When x/D = 1.75, the flow mechanism will be changed significantly. For example, a complete Karman vortex street is formed between the cylinder and the airfoil. In Figure   Figure 10. Pressure contour and pressure coefficient of the cylinder wall at the moment of drag coefficient amplitude (x/D = 2.5).

Influence on the Flow Mechanism
According to computational results listed in Section 3. The vorticity contour when x/D = 1.0 is displayed in Figure 11. At this moment, there is no complete Karman vortex street between the cylinder and the airfoil. Specifically, t 0 -t 5 represent different evolution stages of the flow field respectively. at the lower airfoil is divided into two vortices. Then, the first half is fused with the induced vortex, and the second half is fused with the small vortex shedding from the upper surface. Meanwhile, another induced vortex arises at the tail region of airfoil, as shown in Figure 12c. Afterwards, the trailing edge induced vortices gradually develop and trigger the vortex-breakdown at the lower side. By this time, vorticity of the induced vortex after fusion at the leading edge gets increased, as shown in Figure 12d. Influenced by this, the near-wall vortex rolls up near the airfoil. As shown in Figure 12e, the near-wall vortex spreads to the opposite side of the airfoil and fused with that on the opposite side, while the vortex in the lower surface of the airfoil spreads downstream. In Figure 12f, the induced vortex shedding from the trailing edge is fused with that in the upper airfoil, and the vortex breakdown at the lower airfoil is fused again after shedding and will advance downstream. The wake vortex of the cylinder in the flow field can be called the "primary vortex" for its dominating effect. Taking the local flow field of y < 0 as an example, when the primary vortex meets the NACA0012 airfoil, since the vortex length is far greater than the airfoil chord length, the latter is wrapped by the primary vortex, as shown in Figure 11b. When the primary vortex marches to the middle of the airfoil, a low-pressure area will form, causing the induced vortices generated in the trailing and leading edges of the airfoil successively, as shown in Figure 11c. In the meantime, the airfoil, like a vortex separator, can divide the primary vortex into two parts. As can be observed from Figure 11d, a part of vortex is fused with the induced vortex generated by the upper airfoil in the trailing edge, while the other part is fused with the induced vortex developed at the lower airfoil in the leading edge. After that, the new primary vortex flows past the upper airfoil, as shown in Figure 11e. The direction of the primary vortex at the upper airfoil is opposite to the residual vortex at the lower airfoil, and the primary vortex at the upper airfoil expands and pushes the vortex at the lower airfoil to move along the opposite direction of inflow and collide with the cylinder surface, as shown in Figure 11f. Next, the flow separation at the y < 0 section is aggravated. Hence, the vortex that is pushed to the surrounding wall is accumulated, forming a new primary vortex to repeat the above process.
When x/D = 1.75, the flow mechanism will be changed significantly. For example, a complete Karman vortex street is formed between the cylinder and the airfoil. In Figure 12, red and blue contours represent the counterclockwise vortex and clockwise vortex respectively. Firstly, the vortex street breaks down resulting from the effect of the airfoil. After that, most of the vortex street continues to march along the lower side of the airfoil, while a small part is offset by the vortex in the opposite direction of the upper airfoil. When the vortex is in the middle of the airfoil, induced vortices will be generated in the leading edge of the airfoil, as shown in Figure 12b. Under the effect of airfoil, the vortex at the lower airfoil is divided into two vortices. Then, the first half is fused with the induced vortex, and the second half is fused with the small vortex shedding from the upper surface. Meanwhile, another induced vortex arises at the tail region of airfoil, as shown in Figure 12c. Afterwards, the trailing edge induced vortices gradually develop and trigger the vortex-breakdown at the lower side. By this time, vorticity of the induced vortex after fusion at the leading edge gets increased, as shown in Figure 12d. Influenced by this, the near-wall vortex rolls up near the airfoil. As shown in Figure 12e, the near-wall vortex spreads to the opposite side of the airfoil and fused with that on the opposite side, while the vortex in the lower surface of the airfoil spreads downstream. In Figure 12f, the induced vortex shedding from the trailing edge is fused with that in the upper airfoil, and the vortex breakdown at the lower airfoil is fused again after shedding and will advance downstream. at the lower airfoil is divided into two vortices. Then, the first half is fused with the induced vortex, and the second half is fused with the small vortex shedding from the upper surface. Meanwhile, another induced vortex arises at the tail region of airfoil, as shown in Figure 12c. Afterwards, the trailing edge induced vortices gradually develop and trigger the vortex-breakdown at the lower side. By this time, vorticity of the induced vortex after fusion at the leading edge gets increased, as shown in Figure 12d. Influenced by this, the near-wall vortex rolls up near the airfoil. As shown in Figure 12e, the near-wall vortex spreads to the opposite side of the airfoil and fused with that on the opposite side, while the vortex in the lower surface of the airfoil spreads downstream. In Figure 12f, the induced vortex shedding from the trailing edge is fused with that in the upper airfoil, and the vortex breakdown at the lower airfoil is fused again after shedding and will advance downstream.

Influence on the Vortex Shedding Frequency
Since the airfoil oscillation upon aerodynamic coupling has a greater influence on the wake of the cylinder, three groups of varied T/T0 values after the critical period are selected for the study in this section, that is, T/T0 = 0.5, 1.0, 2.0. The variation of the coupling frequency with x/D is shown in Figure 13. It can be seen that the three groups of varied

Influence on the Vortex Shedding Frequency
Since the airfoil oscillation upon aerodynamic coupling has a greater influence on the wake of the cylinder, three groups of varied T/T 0 values after the critical period are selected for the study in this section, that is, T/T 0 = 0.5, 1.0, 2.0. The variation of the coupling frequency with x/D is shown in Figure 13. It can be seen that the three groups of varied T/T 0 present the same variation trend. When x/D ≤ 1.6, the coupling frequency will gradually decrease as x/D increases, with the minimum coupling frequency at roughly x/D = 1.6. Besides, when x/D ≥ 1.625, the coupling frequency will be on the rise, and when x/D > 3.5, the coupling frequency will become stabilized.
Energies 2021, 14, x FOR PEER REVIEW 12 of 17 1.6. Besides, when x/D ≥ 1.625, the coupling frequency will be on the rise, and when x/D > 3.5, the coupling frequency will become stabilized. When x/D ≤ 1.6, the wake vortex of the cylinder is fully developed with the rising x/D. In this case, longer time for the vortex interaction with the airfoil is required as the length of vortex increases, so as to form a relatively stable vortex to shed from the airfoil trailing edge, as shown in Figure 14a-c. Consequently, the distance between the initial shedding position of the vortex and the airfoil trailing edge increases gradually, together with the decline in coupling frequency. When x/D ≥ 1.625, a stable Karman vortex street forms between the cylinder and the airfoil. Moreover, as x/D increases, the vortex street is dissipated before flowing past the airfoil, leading to a decrease in the strength of the vortex street, along with shorter interaction time with the airfoil. As depicted in Figure 14d-f, the distance between adjacent vortices in the airfoil wake reduces with the growth of coupling frequency.
(a) (b) (c)  1.6. Besides, when x/D ≥ 1.625, the coupling frequency will be on the rise, and when x/D > 3.5, the coupling frequency will become stabilized. When x/D ≤ 1.6, the wake vortex of the cylinder is fully developed with the rising x/D. In this case, longer time for the vortex interaction with the airfoil is required as the length of vortex increases, so as to form a relatively stable vortex to shed from the airfoil trailing edge, as shown in Figure 14a-c. Consequently, the distance between the initial shedding position of the vortex and the airfoil trailing edge increases gradually, together with the decline in coupling frequency. When x/D ≥ 1.625, a stable Karman vortex street forms between the cylinder and the airfoil. Moreover, as x/D increases, the vortex street is dissipated before flowing past the airfoil, leading to a decrease in the strength of the vortex street, along with shorter interaction time with the airfoil. As depicted in Figure 14d-f, the distance between adjacent vortices in the airfoil wake reduces with the growth of coupling frequency. When x/D ≤ 1.6, the wake vortex of the cylinder is fully developed with the rising x/D. In this case, longer time for the vortex interaction with the airfoil is required as the length of vortex increases, so as to form a relatively stable vortex to shed from the airfoil trailing edge, as shown in Figure 14a-c. Consequently, the distance between the initial shedding position of the vortex and the airfoil trailing edge increases gradually, together with the decline in coupling frequency. When x/D ≥ 1.625, a stable Karman vortex street forms between the cylinder and the airfoil. Moreover, as x/D increases, the vortex street is dissipated before flowing past the airfoil, leading to a decrease in the strength of the vortex street, along with shorter interaction time with the airfoil. As depicted in Figure 14d-f, the distance between adjacent vortices in the airfoil wake reduces with the growth of coupling frequency.

Influence on Drag Reduction Characteristics
The relationship between the drag reduction rate and x/D of the cylinder-airfoil tandem configuration under three different T/T 0 conditions is presented in Figure 15. The results show that the drag characteristics varying with x/D is basically unrelated to T/T 0 . When x/D ≤ 1.6, the drag reduction rate is on the rise with the increase of x/D, reaching more than 49% when x/D = 1.6. After x/D ≥ 1.625, the drag reduction rate falls abruptly and then becomes stabilized. Around x/D = 1.6, a large variation is witnessed in the drag reduction rate of tandem configuration, which corresponds to the conclusion made in Section 3.1.2.

Influence on Drag Reduction Characteristics
The relationship between the drag reduction rate and x/D of the cylinder-airfoil tandem configuration under three different T/T0 conditions is presented in Figure 15. The results show that the drag characteristics varying with x/D is basically unrelated to T/T0. When x/D ≤ 1.6, the drag reduction rate is on the rise with the increase of x/D, reaching more than 49% when x/D = 1.6. After x/D ≥ 1.625, the drag reduction rate falls abruptly and then becomes stabilized. Around x/D = 1.6, a large variation is witnessed in the drag reduction rate of tandem configuration, which corresponds to the conclusion made in Section 3.1.2.  Figure 16 depicts the pressure contour of flow fields when x/D = 1.0, 1.25, and 1.5. Apparently, the distance between the primary vortex core and the cylindrical wall on the same side increases with the increase of x/D, together with a low pressure intensity at the center of vortex core. Consequently, when x/D ≤ 1.6, the pressure drag will decrease gradually with the rising drag reduction rate. Evolutionary conditions of the flow field are analyzed in Section 3.2.1. The flow mechanism in the case of x/D ≤ 1.6 is similar to that of x/D = 1.0. Figure 16 depicts the pressure contour of flow fields when x/D = 1.0, 1.25 and 1.5. Apparently, the distance between the primary vortex core and the cylindrical wall on the same side increases with the increase of x/D, together with a low pressure intensity at the center of vortex core. Consequently, when x/D ≤ 1.6, the pressure drag will decrease gradually with the rising drag reduction rate. Comparing flow fields of x/D ≥ 1.625 and x/D ≤ 1.6, the main difference will be witnessed between the cylinder and the airfoil. When x/D ≥ 1.625, the phenomenon of Karman vortex street shedding is found. Since the vortex street has a great effect on the lowpressure area behind the cylinder, there will be an abrupt increase in the time-averaged drag coefficient of tandem configuration near x/D = 1.6.
The pressure contour of flow fields in different x/Ds when T/T0 = 0.5 is displayed in Figure 17. As x/D increases, the strength of the Karman vortex street also increases. It can be seen that the time-average drag coefficient increases slightly after x/D = 1.625. When x/D increases further, the Karman vortex street will expand downstream and become stabilized. At the same time, its resulting low-pressure area has a weakening effect on the cylinder. Although the space between the cylinder and the airfoil is large enough to accommodate the third low-pressure area at x/D = 4.0, the drag reduction effect is less affected for its large gap with the trailing edge of the cylinder. Therefore, the drag reduction rate of the tandem configuration remains unchanged after a small rise.

Conclusions
A NACA0012 airfoil of pitching oscillation is placed downstream of the cylinder. Then, the flow state in the wake of the cylinder changes by changing the pitch oscillation of the airfoil, so as to adjust the distribution of the low-pressure area downstream of the cylinder. The influence of T/T0 (the ratio of the airfoil oscillation period to the Strouhal period of the cylinder) and x/D (the ratio of the distance between airfoil leading edge and cylinder center to the cylinder diameter) on the shedding frequency in the wake of a cylinder and the drag reduction effect of tandem configuration were studied, particularly on the flow mechanism of the flow field. Conclusions are made as follows: 1. When T/T0 is greater than a critical value, the aerodynamic coefficients of the cylinder and the airfoil will vary in accordance with the same frequency, i.e., the phenomenon Comparing flow fields of x/D ≥ 1.625 and x/D ≤ 1.6, the main difference will be witnessed between the cylinder and the airfoil. When x/D ≥ 1.625, the phenomenon of Karman vortex street shedding is found. Since the vortex street has a great effect on the low-pressure area behind the cylinder, there will be an abrupt increase in the time-averaged drag coefficient of tandem configuration near x/D = 1.6.
The pressure contour of flow fields in different x/Ds when T/T 0 = 0.5 is displayed in Figure 17. As x/D increases, the strength of the Karman vortex street also increases. It can be seen that the time-average drag coefficient increases slightly after x/D = 1.625. When x/D increases further, the Karman vortex street will expand downstream and become stabilized. At the same time, its resulting low-pressure area has a weakening effect on the cylinder. Although the space between the cylinder and the airfoil is large enough to accommodate the third low-pressure area at x/D = 4.0, the drag reduction effect is less affected for its large gap with the trailing edge of the cylinder. Therefore, the drag reduction rate of the tandem configuration remains unchanged after a small rise. Comparing flow fields of x/D ≥ 1.625 and x/D ≤ 1.6, the main difference will be witnessed between the cylinder and the airfoil. When x/D ≥ 1.625, the phenomenon of Karman vortex street shedding is found. Since the vortex street has a great effect on the lowpressure area behind the cylinder, there will be an abrupt increase in the time-averaged drag coefficient of tandem configuration near x/D = 1.6.
The pressure contour of flow fields in different x/Ds when T/T0 = 0.5 is displayed in Figure 17. As x/D increases, the strength of the Karman vortex street also increases. It can be seen that the time-average drag coefficient increases slightly after x/D = 1.625. When x/D increases further, the Karman vortex street will expand downstream and become stabilized. At the same time, its resulting low-pressure area has a weakening effect on the cylinder. Although the space between the cylinder and the airfoil is large enough to accommodate the third low-pressure area at x/D = 4.0, the drag reduction effect is less affected for its large gap with the trailing edge of the cylinder. Therefore, the drag reduction rate of the tandem configuration remains unchanged after a small rise.

Conclusions
A NACA0012 airfoil of pitching oscillation is placed downstream of the cylinder. Then, the flow state in the wake of the cylinder changes by changing the pitch oscillation of the airfoil, so as to adjust the distribution of the low-pressure area downstream of the cylinder. The influence of T/T0 (the ratio of the airfoil oscillation period to the Strouhal period of the cylinder) and x/D (the ratio of the distance between airfoil leading edge and cylinder center to the cylinder diameter) on the shedding frequency in the wake of a cylinder and the drag reduction effect of tandem configuration were studied, particularly on the flow mechanism of the flow field. Conclusions are made as follows: 1. When T/T0 is greater than a critical value, the aerodynamic coefficients of the cylinder and the airfoil will vary in accordance with the same frequency, i.e., the phenomenon

Conclusions
A NACA0012 airfoil of pitching oscillation is placed downstream of the cylinder. Then, the flow state in the wake of the cylinder changes by changing the pitch oscillation of the airfoil, so as to adjust the distribution of the low-pressure area downstream of the cylinder. The influence of T/T 0 (the ratio of the airfoil oscillation period to the Strouhal period of the cylinder) and x/D (the ratio of the distance between airfoil leading edge and cylinder center to the cylinder diameter) on the shedding frequency in the wake of a cylinder and the drag reduction effect of tandem configuration were studied, particularly on the flow mechanism of the flow field. Conclusions are made as follows: 1.
When T/T 0 is greater than a critical value, the aerodynamic coefficients of the cylinder and the airfoil will vary in accordance with the same frequency, i.e., the phenomenon of "aerodynamic coupling" will occur. According to our simulation results, the ratio of