Direct Numerical Simulation of Flow around a Circular Cylinder Controlled Using Plasma Actuators

Flow around a circular cylinder controlled using plasma actuators is investigated by means of direct numerical simulation (DNS). The Reynolds number based on the freestream velocity and the cylinder diameter is set at Re D = 1000. The plasma actuators are placed at ±90 from the front stagnation point. Two types of forcing, that is, two-dimensional forcing and three-dimensional forcing, are examined and the effects of the forcing amplitude and the arrangement of plasma actuators are studied.The simulation results suggest that the two-dimensional forcing is primarily effective in drag reduction. When the forcing amplitude is higher, the mean drag and the lift fluctuations are suppressedmore significantly. In contrast, the three-dimensional forcing is found to be quite effective in reduction of the lift fluctuations too. This is mainly due to a desynchronization of vortex shedding. Although the drag reduction rate of the three-dimensional forcing is slightly lower than that of the two-dimensional forcing, considering the power required for the forcing, the three-dimensional forcing is about twice more efficient.


Introduction
Flow around a bluff body causes a significant amount of drag.In addition, the vortices shed from the bluff body cause flow-induced noise; the resulting lift fluctuations may damage instruments.Therefore, significant efforts have been made to control the flow around a bluff body (see, e.g., Choi et al. [1]).
Flow control methods can be classified into passive and active ones.As the passive control methods, a splitter plate [2], dimples [3], tabs [4], and wavy surfaces [5] have been used, to name a few.For the active control methods, different types of actuators have been developed, as reviewed, for example, by Cattafesta and Sheplak [6].As an example of active control of flow around a bluff body, Kim and Choi [7] performed a numerical simulation of a flow around a circular cylinder at Re  =  ∞ /] = 100 − 3900 (where ,  ∞ , and ] denote the cylinder diameter, the freestream velocity, and the kinematic viscosity, resp.)controlled by using blowing and suction and attained more than 20% drag reduction.
Recently, plasma actuators have attracted considerable attentions as novel devices for active control.The widely used single-dielectric barrier discharge (SDBD) plasma actuator is composed of two electrodes and a dielectric layer between them, as shown in Figure 1.A high frequency alternatingcurrent (AC) voltage produces plasma between these electrodes.The charged particles are accelerated through the electric field and subsequently collide with neutral particles to generate a body force and induce a wall-jet-like flow.The plasma actuator has advantages in its light weight and no moving mechanical part; therefore, it is expected to be used in a wide range of applications.For more details on the physics and the applications of SDBD plasma actuators, readers are referred to the excellent review papers by Roth and Dai [8] and Corke et al. [9].
Grundmann et al. [10] experimentally and numerically investigated the effect of plasma actuator applied to a boundary layer separation control and verified its effectiveness.Cho and Shyy [11] considered an airfoil with a plasma actuator and investigated the effect of both open-loop and closed-loop controls.In these previous studies, the so-called two-dimensional forcing was used.In the two-dimensional forcing, the body force by the plasma actuator is directed in the streamwise direction and its magnitude is uniform in the spanwise direction, so that the body force delays the flow separation by directly adding a mean momentum into the boundary layer.
In contrast, the so-called three-dimensional forcing has also been studied.In the three-dimensional forcing, the magnitude of body force varies in the spanwise direction or the body force vectors are directed in the spanwise direction, so that three-dimensional flow structure is induced to destroy the two-dimensional vortices.Kozlov and Thomas [12] experimentally investigated a flow around a circular cylinder at a high Reynolds number (Re  = 85, 000) controlled using plasma actuators.They considered not only two-dimensional forcing but also three-dimensional forcing.In their three-dimensional forcing, the plasma actuators are placed in the circumferential direction and the flow is induced in the spanwise direction.Bhattacharya and Gregory [13] experimentally examined three-dimensional forcing at Re  = 4700.In their three-dimensional forcing, the plasma actuators inducing streamwise flows are placed on a circular cylinder periodically in the spanwise direction with the wavelength of  = 4.In both studies, a better drag reduction effect than that by two-dimensional forcing was obtained.
Despite the studies introduced above, there have been fewer studies on the flow around a circular cylinder controlled using the three-dimensional forcing of plasma actuators.Detailed effects of the forcing amplitude of the twodimensional and three-dimensional forcings on the flow modifications are still to be investigated.Therefore, in the present study, we investigate the flow around a circular cylinder controlled using plasma actuators by means of direct numerical simulation (DNS).We consider both the two-dimensional forcing and the three-dimensional forcing and investigate the effect of forcing amplitudes on the flow modifications.

Direct Numerical Simulation.
We consider a circular cylinder in an incompressible uniform flow.The Reynolds number is Re  = 1000.This Reynolds number is chosen so that a direct numerical simulation (DNS) can be performed, since no turbulence models have been validated for a flow controlled using plasma actuators.The governing equations are the continuity and the Navier-Stokes equations; that is, where ⃗  and  are the velocity and the pressure, respectively, and ⃗  denotes the body force induced by the plasma actuators described later.All the quantities are made dimensionless by using the fluid density, , the freestream velocity,  ∞ , and the cylinder diameter, .
The present simulations are performed by using the DNS code of Naito and Fukagata [14].The governing equations are spatially discretized in the cylindrical coordinates by using the energy conservative second-order accurate finite difference method [16] and temporally integrated by using the lowstorage third-order Runge-Kutta/Crank-Nicolson (RK3/CN) scheme [17,18] with an SMAC-like pressure correction [19].The pressure Poisson equation is solved by using the fast Fourier transform (FFT) in the azimuthal () and the spanwise () directions and the tridiagonal matrix algorithm (TDMA) in the radial () direction.
The computational grid is shown in Figure 2. The size of computational domain (  ,   ), the number of computational cells (  ,   , and   ), and the size of a computational cell (Δ min , Δ) are summarized in Table 1.The computational grid is nonuniform in the radial () direction and uniform in  and  directions.The size of the computational domain is 70 to minimize the influence from the boundaries.
A uniform velocity  ∞ = 1 is imposed at the inlet boundary (|| ≤ (3/4), where  denotes the angle from the front stagnation point) and the convective velocity condition, that is, is used at the outlet boundary (|| > (3/4)).The periodic boundary condition is applied in  direction.On the surface of the circular cylinder, the no-slip condition is applied and the pressure gradient is set to be zero.Verification and validation for the base flow (i.e., without plasma actuators) have been reported in Naito and Fukagata [14].The computational cell has been confirmed to be sufficiently small to obtain the statistical data; for instance, change in the mean drag coefficient was less than 1% when twice finer Δ was used.Good agreement with the data available in literature has been obtained on the mean drag coefficient, root mean square (RMS) of the lift coefficient, the Strouhal number, and the base pressure coefficient.For more details, the readers are referred to Naito and Fukagata [14].

Body Force by Plasma
Actuator.The electric field formed around the electrodes is a consequence of very complex ion transport phenomena.Therefore, if one is interested in temporal variations of the body force field, the detailed ion transport should be solved.In contrast, if one is interested in the effect of body force on a relatively low-speed fluid flow, such as the present study, it is more efficient to use a timeaveraged model because the frequency of the imposed AC voltage signal is much higher than that of flow.For instance, an air flow at  ∞ = 1.5 m/s on a  = 1 × 10 −2 m cylinder results in Re  = 1000 and the vortex shedding frequency of 30 Hz, while the frequency of imposed AC voltage signal is of the order of 10 kHz.Yamamoto and Fukagata [20], Table 1: The size of computational domain (  ,   ), the number of computational cells (  ,   , and   ), and the size of a computational cell.who performed detailed numerical simulations of the ion transports in a plasma actuator, demonstrated that a timeaveraged body force field can reasonably be used when the typical frequency of flow is less than 30 times of that of the voltage signal.
As a time-averaged model, we use Shyy's model (Shyy et al. [15]), as shown in Figure 3.This model assumes a linear decay of the electric field, (, ), with its maximum at the point of the shortest distance between the electrodes, which can be expressed by where  0 is the maximum electric field strength with  1 and  2 being positive constants.Here, the values proposed by Shyy et al. [15], that is,  0 = 226.27kV/cm,  1 = 6.54 × 10 2 kV/cm 2 , and  2 = 13.08 × 10 2 kV/cm 2 , are used.The electric field vector is assumed to be unidirectional and its components   and   are modeled as The time-averaged body force, ⃗ , is computed as where   = 1.60 × 10 −19 C is the elementary charge,  (the frequency of the applied voltage signal), Δ (the discharge time), and   (the charge number density) are the parameters depending on the operating condition, and  is the collision efficiency assumed to be unity.By representing these model constants altogether by a single parameter, , the magnitude of the body force can be expressed as The purpose of the present study is to investigate the effect of forcing amplitude on the flow modification.Therefore, we set the parameter  as The forcing amplitude is varied by  3 from 0.0 to 1.2, as shown in Table 2.The resultant induced velocity will be shown in Section 3 as a function of  3 .We examine four types of arrangement as shown in Figure 4.The uncontrolled flow (i.e., base flow) is denoted as Case NC; Case 1 is a two-dimensional forcing; Cases 2 and 3 are three-dimensional forcings.In all cases, the plasma actuators are placed at ±90 ∘ from the front stagnation point.

Results and Discussion
3.1.Two-Dimensional Forcing. Figure 5 shows the maximum mean streamwise velocity around a plasma actuator,  max , as a function of the forcing amplitude,  3 .The velocity is found to increase monotonically as the amplitude parameter  3 increases.In the figure, the results from four different simulations, Cases 1a-1d, corresponding to four different time instants when the plasma actuators are turned on, are presented, because it is found from the results that the resultant values of   and    (shown below) are slightly dependent on the turn-on time.
The dependency on the turn-on time is likely to be due to the low frequency variation of lift fluctuations.Figure 6 shows the time traces of the drag coefficient,   , and the lift coefficient,   , defined, respectively, as where   and   denote the drag and the lift forces.In addition to the high frequency fluctuations corresponding to the vortex shedding ( ≃ 5), very low frequency fluctuations ( ∼ 100) are observed.The low frequency fluctuations can be explained by instantaneous vortical structures (identified by using the second invariant of deformation rate tensor) at two different time instants, as exemplified in Figure 8.At  = 50 (when the amplitude of lift fluctuations is small), fine streamwise vortices are populated in the wake and spanwise vortices are hardly seen in the far wake.On the other hand, at  = 200 (when the amplitude of lift fluctuations is large), spanwise vortices are pronounced in the far wake.Namely, the amplitude of lift fluctuations corresponds to the degree of two-dimensionality of the vortical structure.
Due to the low frequency fluctuations, the integration time should be set carefully in the computation of statistics.Figure 7 shows the mean drag,   , and the root mean square (RMS) of lift fluctuations,    , as functions of the integration time, .From this, we judge that  = 500 is sufficient as the integration time for accumulation of statistics.
Figures 9(a) and 9(b) show   and    as functions of the forcing amplitude,  3 .The mean drag,   , is found to monotonically decrease as the forcing amplitude increases; the maximum drag reduction rate is about 45% at  3 = 1.2.On the other hand,    exhibits complicated dependency on  3 ; it is at a comparable level to that of Case NC up to  3 = 0.6 and it significantly decreases for  3 > 0.6.The relationship between   and    is mapped in Figure 9(c).It is clearly observed that up to  3 = 0.6 the drag is reduced, while the lift fluctuations are not much suppressed; for  3 > 0.6, both the drag and the lift fluctuations are suppressed.From this, it is conjectured that there is a transition of flow characteristics around  3 = 0.6.
The transition of flow characteristics is illustrated by the mean streamwise velocity near the surface (Figure 10) and instantaneous vortical structures (Figure 11).On the cylinder surface, separation is found to be delayed at all amplitudes.Moreover, separation looks more synchronized in the spanwise direction at higher amplitudes.In the wake, the vortical structure is considerably different between the cases of weaker forcing ( 3 = 0.4) and stronger forcing ( 3 = 0.8 and 1.2).The vortical structure at  3 = 0.4 is not so much different from that of the uncontrolled flow.At  3 = 0.8 and 1.2, while two-dimensionality of separation is enhanced, the shed vortices seem to be too weak to generate strong streamwise vortices in the far wake.From the observations above, two different mechanisms are considered to play in the modification of flow: one is the separation delay, which reduces the drag and the amplitude of shed vortices; the other is the spanwise synchronization of the vortex shedding, which raises the lift fluctuations.The change in    suggests that these two effects compensate each other at  3 < 0.6 and the former effect exceeds at  3 > 0.6.
Figure 12 shows the distributions of pressure coefficient on the cylinder surface,   , as functions of the angle from the stagnation point, , defined as where   and  ∞ denote the mean wall pressure and the freestream pressure, respectively.The effect of plasma actuator is similar at all amplitudes, although there is a difference in its magnitude.In 30 ∘ <  < 90 ∘ , the pressure decreases due to the acceleration by the plasma actuator; at  = 90 ∘ , the pressure rapidly recovers by gaining momentum from the plasma actuator; around  = 100 ∘ , the pressure gradient reverses with flow separation; and finally on the rear surface the pressure recovers again with a larger value than that in Case NC.From the definition of drag, the greater recovery on the rear surface is most responsible for the drag reduction.
On the other hand, the change in    can be explained by the RMS pressure coefficient,    , defined as as shown in Figure 13.The dependency of    on  3 is complicated.At  3 = 0.1 and 0.2,    decreases and the peak location moves backward.From  3 = 0.2 to  3 = 0.4, the maximum value is almost unchanged, but the peak location still moves backward.The peak value increases again at  3 = 0.6 to the same level as that of Case NC; then    significantly decreases for higher  3 with the peak location nearly unchanged.This behavior is in accordance with the complicated dependency of    on  3 observed in Figure 9 and it can be explained by the two different effects discussed above: the suppression of vortex shedding amplitude and the two-dimensionalization of separation.In particular, the considerable reduction of    at  3 > 0.6 can be attributed to the suppression of    on the entire surface.

Three-Dimensional Forcing.
In the cases of three-dimensional forcing, the results are found to be independent of the initial flow conditions unlike the cases of two-dimensional forcing.
Figure 14 shows the maximum mean streamwise velocity around the upper plasma actuator in Case 2. Positions 1 and 2 denote the locations where the plasma actuator is present and absent, respectively.As  3 increases,  max at Position 1 is found to increase nearly monotonically; however, its magnitude is smaller than that in Case 1.In contrast, the velocity at Position 2 is less than that in Case NC.The difference in the mean streamwise velocity fields at Positions 1 and 2 is illustrated in Figure 15.While at Position 1 the separation is delayed similarly to the case of two-dimensional forcing (Figure 10(b)), and the velocity above the cylinder is lower than that in Case 1.The velocity at Position 1 is likely to be affected by the low velocity at Position 2. Figure 16 shows instantaneous three-dimensional vortical structures.Vortices are observed to have completely been vanished in the near wake.The far wake is turbulent; however, no clear two-dimensional structure can be observed unlike Case NC and Case 1.
Figure 17 shows   and    as functions of  3 .The dependency of   on  3 in Cases 2 and 3 is similar to each other.Unlike the case of two-dimensional forcing, there is an optimum amplitude ( 3 = 0.4) at which   is minimized; the maximum drag reduction rate is about 34%.An excellent feature of the three-dimensional forcing appears in the lift fluctuation,    ; it is drastically suppressed in a wide range of  3 due to three-dimensionalization of vortex shedding.The relationship between   and    is mapped in Figure 18.In Cases 2 and 3, both the mean drag   and the lift fluctuations    are decreased from those in Case NC.In particular, reduction in    is much larger than that in Case 1 even at low forcing amplitudes.In both cases, the maximum reduction Mathematical Problems in Engineering  rates are obtained in 0.2 ≤  3 ≤ 0.6, that is, 34% for drag and 96% for lift fluctuations.The significant reduction in the mean drag and the lift fluctuations can be explained by the distribution of the pressure coefficient,   , and the RMS pressure coefficient,    , as shown in Figure 19.The shape of the pressure distribution on Position 1 is similar to that in Case 1, while that on Position 2 is similar to that in Case NC.On both Positions 1 and 2, however, the pressure in the rear half is higher than that in Case 1.This is considered due to the enhanced mixing of momentum by the three-dimensional forcing.The pressure fluctuations are much less than those in Case NC and Case 1 due to the three-dimensionalization (i.e., desynchronization) of vortex shedding; thus, the lift fluctuations are significantly decreased.

Conclusions
We have performed direct numerical simulation of a flow around a circular cylinder at the Reynolds number of Re  = 1000 controlled using plasma actuators.A two-dimensional forcing and two types of three-dimensional forcing have been examined with different values of forcing amplitude.
In the two-dimensional forcing, higher forcing amplitude of the plasma actuators results in larger reduction of both the mean drag coefficient (  ) and the RMS lift coefficient (   ).In the best case,   and    are reduced by 44% and 66%, respectively.The mechanism for these reductions is a direct feeding of mean momentum into the boundary layer, which delays the separation; thus, the pressure recovers on the rear surface and the vortex shedding is weakened.
On the other hand, the three-dimensional forcing is found to work well too.Regardless of the forcing amplitude, greater effects are obtained in the suppression of lift fluctuations compared to the case of two-dimensional forcing.Moreover, the drag is reduced as much as that in the twodimensional forcing case, even at lower forcing amplitudes.In the best case,   and    are reduced by 34% and 96%, respectively.The mechanism is found to be different from that for the two-dimensional forcing.The reduction of the mean drag and the lift fluctuations is due to three-dimensionalizaton of vortex shedding, which enhances mixing of momentum and desynchronizes the phase of vortex shedding.
In the present study, the total length of plasma actuators in the three-dimensional forcing is half of that in the two-dimensional forcing.Therefore, considering the power required for the present three-dimensional forcing, which is also nearly half of that for the two-dimensional forcing, the three-dimensional forcing is considered nearly twice more efficient than the two-dimensional forcing.Of course, validity of the present findings is limited to the case of Re  = 1000.Although a similar effect may be expected at moderately high Reynolds numbers, the amounts of drag reduction and lift suppression should be significantly dependent on the Reynolds number.Moreover, the situation may drastically change around and above the critical Reynolds number (Re  ∼ 10 5 ).Toward practical applications, a comprehensive study on the Reynolds number dependency should be conducted as a future work.

Table 2 :
Forcing amplitude of plasma actuators.