On the interplay of body-force distributions and flow speed for dielectric-barrier discharge plasma actuators

The dielectric-barrier discharge plasma actuator is a well-established device commonly operated in boundary-layer airflows for active flow control. In the present experimental investigation, their ability to cause momentum transfer to the surrounding fluid is analyzed by means of spatio-temporal body-force distributions in both quiescent air and external airflow conditions. The work is motivated by the limitation to quiescent-air operating conditions of frequent previous efforts. Available analytical velocity-information-based force derivation approaches are contrasted to investigate the actuator performance under conditions of their area of application. Results of body force in quiescent air, in agreement with literature, confirm the major taken assumption for Navier–Stokes-based body-force formulations—a negligible pressure gradient. However, the previous circumstance turns out as an invalid


Introduction
Ever since political agenda considers climate change a relevant topic, different technological domains push for environmental-friendly solutions, where aircraft industries is likewise in charge to reduce fuel consumption and greenhousegas emissions. Various active and passive means to control and accordingly optimize the flow scenario at hand have been developed and evaluated for flow control applications [1].
In the past two decades, also dielectric barrier discharge (DBD) plasma actuators (PAs) have considerably matured to competing devices for active flow control (AFC) [2]. In order to comprehensively evaluate the performance of PAs for AFC, numerous characterization studies have been performed, as continuously summarized in review articles [3][4][5][6][7]. Similarly, a vast amount of PA-force modeling efforts with varying degree of accuracy and (accordingly) complexity have been undertaken, which range from straight forward phenomenological to advanced first principle models [8][9][10][11][12].
However, such studies still leave open questions on the details of the underlying physical working principle of PAs for AFC. A good summary of the respective strategies to determine the momentum transfer from the charged particles to neutral gas has been provided by Jayaraman and Shyy [13]. More recent model development studies largely build upon these earlier concepts and address various modification strategies of both phenomenological [14][15][16] and first-principle [17][18][19][20][21] models.
In addition, complementary efforts to determine the spatiotemporal body-force distributions of PAs on the grounds of experimentally obtained velocity fields have been undertaken to contribute towards PA-force-field estimates, which mostly applied particle image velocimetry (PIV) to achieve the required velocity information. These approaches do not rely on physical descriptions of the gas discharge, but only apply Navier-Stokes or vorticity equations (NSE, VE) on the measured wall jet in proximity of the PA to solve for the unknown force fields [22][23][24][25][26][27].
Furthermore, Maden et al [28] converted PIV-based force fields into an empirical model, in order to allow a more straight forward implementation of the experimentally derived PA force in CFD. Likewise, further (semi-)empirical models have been proposed to predict the force distributions under absence of spatially resolved velocity information, but by means of various discharge-specific quantities [29][30][31][32].
The direct comparison of different phenomenological and PIV-based models demonstrated reasonable accuracy at sufficient downstream distance from the discharge zone. In contrast, considerable differences from the reference flow field became salient in immediate vicinity of the discharge zone. Depending on the target application, such remaining model uncertainties might be acceptable. However, given the spatiotemporal body-force distribution to be of utmost importance for the flow-control success-such as, e.g. topology pinning [33] or virtual wall oscillations [34,35]-this limited spatial model accuracy might imply an oversimplified representation of the momentum transfer.
Moreover, the vast majority of the above-reported research into PA characteristics as yet addresses quiescent-air studies, which imply the assumption that the PA manipulates the airflow and not vice versa, where research into the influence of wind has received very little attention so far.
First investigations of Baughn et al [36] into the interplay of PA force and small free-stream velocities U ∞ ⩽ 6.8 m s −1 (M < 0.02) revealed no significant airflow influence on the PA performance. The little number of publications (e.g. [37][38][39][40][41]) on this influence above U ∞ = 6.8 m s −1 , however, provide strong evidence that the above quiescent-air assumption in fact might turn out as an oversimplification.
Early studies by Pavon et al [37] uncovered changes of the DBD towards a more filamentary character with significantly smaller current peaks, when the plasma discharge encounters transonic and low supersonic airflow conditions. In addition to these morphology changes, Kriegseis et al [38] identified significant reductions of plasma extent ∆x and underlying actuator power consumption P A , when the DBD PA was installed in co-flow configuration at low and moderate Mach numbers (M ⩽ 0.78). The observed performance drop was in the range of 4% at M = 0.1 and 30% at M = 0.4. The influence of pressure on the exerted body force [42] is particularly important to be considered in scenarios of plasma-airflow interaction, giving rise to a competition between the implications of pressure and airflow on the plasma discharges [40]. Furthermore, the interplay of plasma discharge and external airflow has been studied by Fan et al [43] and Wang et al [44], where the authors found a reduction in current-peak intensity and a simultaneous increase of current-pulse density for increasing airflow speed.
To advance beyond only discharge-specific studies of the PA-performance drop, Pereira et al [39] implemented PAs in both co-and counter-flow configuration for M ∞ ⩽ 0.2 and measured the resulting actuator thrust with a load cell in addition to P A . Even though no significant electrical performance drop was determined within the uncertainty margin, the load cell measurements uncover an important interplay between the momentum of PA and the boundary layer. Interestingly, for co-flow forcing the actuator thrust remained constant for M ∞ < 0.1 and was found to increase for higher external flow speeds. Counter-flow forcing instead was found to result in a constant actuator thrust throughout the range of investigated flow speeds M ∞ ⩽ 0.2.
The observed increase of net thrust has been attributed to the changing wall shear [39], where the increasing wall shear accordingly implies a diminishing contribution of self-induced drag [45] to the overall momentum transfer into the boundary layer. It is worth to mention that this insight indicates an improved fluid-mechanic effectiveness and/or efficiency [46][47][48] of PAs under airflow operation-despite the observed reductions of electric performance. A deeper insight into this interplay, however, requires knowledge of the body-force distributions during airflow encounter, which remains yet to be investigated.
The first pioneering contribution towards such information has been published by Dörr and Kloker [41], who applied the so-called Maden model [28] in direct numerical simulations (DNSs) of both a quiescent-air wall jet and a M ∞ = 0.1 boundary layer. Analysis of the resulting velocity fields with the above-mentioned NSE and VE approaches [22,23] and subsequent comparison with the implemented force fields revealed close agreement between modelled and extracted force fields for the quiescent-air case as well as a negligible pressure gradient in the DNS data. Additionally, the VE approach extracted almost identical force fields from the M ∞ = 0.1-case. In contrast, the NSE-based estimations changed significantly, which resulted in unphysical distributions of body force and corresponding pressure gradient.
Despite the fact that this DNS study is an inherently loosely-coupled approach, the outcomes particularly indicate that the NSE approach might oversimplify the problem for U ∞ > 0, as a result of the increasing pressure gradient. The main motivation for the present work, therefore, is an immediate continuation of Dörr and Kloker's work [41] on the grounds of experimental data of a DBD-manipulated boundary layer.
The present investigation, accordingly, comprises several interrelated objectives. First and most importantly, the quiescent-air assumptions made to determine the actuator body force from both NSE and VE approaches [22,23] require further verification in such a way that the body-force distributions are directly extracted from flow fields, where the DBD encounters a range of flow speeds U ∞ . These results will allow to uncover possible changes of the force distributions under airflow influence.
Additionally, it remains yet to fully clarify how the airflow onset manipulates the electrical and/or fluid-mechanical performance of the PA. As such, electrical information will be recorded simultaneously to the PIV measurements. Further analysis of the individual terms of NSE and VE approach is foreseen to provide valuable information regarding advanced efforts to derive more sophisticated empirical PA models. In order to provide a more robust reference quantity, the integrated force distributions are furthermore compared to the (phase-resolved) integral approaches [36,42,49].
Finally, as an additional technical objective the present study attempts to provide an outline of requirements and limitations of the standard PIV method to derive the underlying velocity fields for the subsequent body-force estimation efforts, where dynamic range and motion blur are in direct competition to e.g. the mandatory high spatial resolution.

Experimental setup and procedure
The time-dependent momentum transfer of the PA to the surrounding fluid in quiescent air and under presence of an external airflow is investigated using planar high-speed PIV. The experiments are performed in a blower-type wind tunnel at ISTM, KIT, which features optical access through a transparent test-section module. The dimensions of the test section in length, width and height are 700 × 320 × 220 mm 3 , respectively. A flat plate with elliptic leading edge has been horizontally mounted in the center of the test section, to provide support for a co-flow oriented PA of L = 150 mm length in spanwise z direction (wall-jet along mean-flow direction x). The dimensions of the flat-plate model in length, width and height are 420 × 320 × 20 mm 3 , respectively.
In continuation of earlier studies [25,27] the electrodes were made of copper and five layers of Kapton (0.4 mm thickness) were chosen as dielectric material. The width of exposed and encapsulated electrodes was set to 2.5 and 10 mm, respectively, without an inter-electrode gap. A high-voltage generator (Minipuls 2, GBS Elektronik GmbH) is used to supply the PA with V pp = 12 kV peak-to-peak voltage at an alternatecurrent (AC) frequency f ac = 10 kHz (sine wave). The wind tunnel was either switched off (U ∞ = 0) or has been operated at free-stream velocities of U ∞ = 2.5, 5, 7.5, 10, 15, 20 and 30 m s −1 , which corresponds to Re ∞ = 1.9 · 10 6 1 m −1 for the fastest flow speed of U ∞ = 30 m s −1 . The PA installed on the flat plate and all utilized measurement equipment is sketched in figure 1.

PIV
The PIV high-speed system is comprised of a Nd:YLF dualcavity laser (Quantronix Darwin Duo) and a Photron Fastcam SA4 camera. As indicated in figure 1, the laser illuminated the field of view (FOV) from the tunnel outlet into the test section. Also, a picture of the arrangement of experimental components is added in figure 2 for clarity. The camera has been placed outside the test section and was equipped with a Nikon AF micro Nikkor 200 mm f/4D IF-ID lens (f # = 8) and an additional extension ring of 108 mm total length, in order to simultaneously obtain a suitable object distance and a sufficiently high spatial resolution of 96 px mm −1 . The resulting FOV spans 8.5 × 2 mm 2 (x × y) in the x-y plane and is located at mid-span (z = 0) on both the plate and the PA, i.e. the x axis is along the mean-flow direction, y is the wall-normal coordinate and z is along the span. The coordinate-system origin x = y = z = 0 is located at the downstream edge of the exposed electrode on the flat-plate surface 57 mm downstream of the elliptical leading edge. . The PA is mounted on a flat plate with elliptical leading edge; phase-resolved high-speed planar PIV is applied to record a field of view (FOV) in proximity of the DBD; electrical equipment is indicated, where either the probe capacitor (Cp) or the shunt resistor (Rs) is connected to the electric circuit. The frame of reference (x = y = z = 0) sits at the exposed electrode's trailing edge in the spanwise center of the PA. The flow was seeded with d p = 0.9 µm di-ethylhexylsebacat (DEHS) tracer particles (density ρ p = 912 kg m −3 , dynamic viscosity η, particle-response time τ p = ρ p d 2 p /18η = 2.25 µs), resulting in a Stokes number of Stk ∞ = τ p U ∞ /l c ⩽ 6.75 × 10 −2 for PIV measurements in the considered range of free-stream velocities of U ∞ ⩽ 30 m s −1 and a characteristic length scale of l c = 1 mm (cp. [27,50]). For PIV measurements in quiescent air only the test-section chamber was seeded, whereas for PIV measurements with external airflow, the seeding was supplied to the inlet nozzle of the wind tunnel. Phase-resolved velocity information over the plasma-discharge cycle has been acquired in doubleframe mode with 24 equidistant phases per discharge cycle (phase resolution N ϕ = 24), which resulted in a phase-tophase spacing and time of ∆ϕ = π/12 and of ∆t ϕ = 4.17 µs, respectively.
For each PIV run 5400 image pairs were recorded with f s = 6.570 kHz camera frame rate and 4 µs inter-frame time or laser pulse distance, which corresponds to 3285 double-frames per second and 225 velocity fields per phase. Furthermore, all runs were repeated three times to achieve 675 fields per phase in total and thus ensure both repeatability and statistical significance of the phase-resolved data. Note that also the base-flow velocity field was recorded for each U ∞ for reference.

Data processing and measurement uncertainty
The raw images were pre-processed with Matlab, where subtraction and surface masking was applied. Note that also the first two pixel rows above the surface were masked to eliminate reflections at the actuator electrodes. Subsequently, the velocity fields were cross-correlated with PIVview2C software (version 3.8.0, PIVTEC GmbH) in a multigrid/ multipass approach with an initial interrogation area (IA) size of 96 × 96 px 2 . The final IA size was 32 × 16 px 2 (x × y) with overlaps of 75% (x) and 50% (y). Outliers (<1%) were eliminated based on the normalized median test [51]. Furthermore, statistical significance was verified for characteristic locations (fastest velocity •, shear layer ♢, wall jet □ and suction △) in the wall jet according to Kriegseis et al [25], see figure 3. Note that identical locations have been chosen for quiescent air and U ∞ = 30 m s −1 for consistency purposes.
The relative standard deviation of the time-and phaseaveraged velocities is shown in figures 3(b) and (c), respectively. Figure 3(b) indicates that at least the order of N = 10 2 samples is mandatory to achieve a statistically meaningful basis, where also the footprint of the AC character can be identified from the high fluctuation level at the suction region of the wall jet (cp. also [25]). Note that all other cases with 0 < U ∞ < 30 m s −1 also converged within this range. Furthermore, figure 3(c) reveals that the phase-resolved fields converge for N > 200 snap shots, which essentially corresponds to the same number of passed discharge cycles. At two phase positions, the data for U ∞ = 30 m s −1 shows a significant increase of the standard deviation. However, the figure inset indicates reasonable convergence level for N > 200, as well.
The accuracy of the PIV-velocity fields was additionally evaluated with an uncertainty-quantification strategy, reported by Sciacchitano and Wieneke [52]. This approach determines the uncertainty of the time-averaged velocity data from the standard deviation of the velocity field, being computed with fluctuation components and measurement errors (systematic error sources can not be quantified). Individual analysis of the uncertainty associated with each single PIV case with plasma actuation resulted in a maximum uncertainty of <2% for quiescent air and <4% for U ∞ = 30 m s −1 . The uncertainty of all other cases is established within this range (not shown for brevity).
Despite this low net uncertainty, it is important to mention motion blur of the particle images as an additional uncertainty. Motion blur leads to non-Gaussian image conditions [53], which imply a diminished signal-to-noise ratio [54] of the inherently broadened Airy disks and correlation peaks [55] such that additional issues arise for the corresponding displacement estimation [56]. This issue becomes particularly sensible for particle acceleration during exposure, which in turn further skews the smeared particle-image intensity [57]. For the given setup (0.12 µs laser pulse width, 96 px mm −1 spatial resolution) and the investigated velocity range (U ∞ ⩽ 30 m s −1 ) the motion blur of the particle image remains below 0.35 px for all experiments and is, therefore,   expected to mostly affect the results for high velocities in combination with high local accelerations.

Base-flow boundary layer
The velocity fields of the base flow without plasma actuation were evaluated by means of boundary-layer thickness δ 99 and shape factor H 12 (see e.g. [58]) in proximity of the PA to identify the local boundary-layer condition. The extracted boundary-layer properties at x = 0 are summarized in table 1. The shape factor H 12 is above 2.59 (laminar Blasius flow) for all free-stream velocities, which is attributed to a slightly decelerating flow and correspondingly small adverse pressure gradient along the streamwise direction. The freestream velocity was limited to U ∞ = 30 m s −1 , to ensure laminar flow in the FOV. To provide additional information on the effective flow speed as encountered by the discharge and the exerted body force, the characteristic mean velocity estimate U 0 = ⟨u(x = 0, y ⩽ 1.5)⟩ y across the FOV-height has also been added to the table together with an accordingly diminished Stokes number Stk 0 < Stk ∞ . Based on the U 0 range in table 1, the before-indicated upper bound of 0.35 px motion blur accordingly accommodates values between 0.01 and 0.24 px in the near-wall region (y ⩽ 1.5 mm). The time-averaged flow field u(x, y) and a corresponding laminar velocity profile at x = 0 are plotted together in figure 4(a). The inflection point of the profile above the wall indicates the adverse pressure gradient ∂p/∂x > 0 in the boundary layer, which has been quantified from the nonactuated reference velocity fields u(x, y) by means of the steady two-dimensional (2D) NSE without any body-force term f x (cp. section 4, equation (1a)). The resulting gradient field is shown in figure 4(b), which ranges around the same orders of magnitude (10 0 to 10 1 kN m −3 ) as typical body forces from earlier studies [25,26]. Consequently, the distributions of streamwise pressure gradient ∂p/∂x are determined for all flow speeds considered and included for the respective force derivations, as will be discussed in section 4.1.

Electrical characteristics
Direct recordings of the individual current peaks provide sensible information on the periodic discharge onset after each dark period [59,60]. Earlier studies have, however, shown that the Lissajous-figure-based approach provides a more robust means to quantify the actuator power consumption P A [61,62]. As such, the combination of both measurements provides sufficient information to study the interplay of discharge characteristics, outer airflow velocity and resulting body-force distributions [26,38,63].
As indicated in figure 1, a charge-probe capacitor with C p = 22 nF was used for the Lissajous-figure analysis, which was temporarily replaced by a shunt resistor (R s = 2.2 Ω) to identify the regimes of discharge collapse of the plasma from the current histories I(t) [63]. The choice for C p = 22 nF is based on the study of Kriegseis et al [61], where the authors used the identical DBD PA as for the current work. An Agilent oscilloscope (InfiniiVision 2000 X-Series) with 9 bit vertical resolution and adjusted acquisition rates of 5 MHz and 50 MHz was used for capacitor and resistor measurements, respectively. It is to be noted that the acquired current signal is not used for quantitative analysis in the present study, but rather serves a phenomenological orientation within the evaluated discharge cycle. Figure 5 shows the quiescent-air Lissajous figure with superimposed current signal as a red solid line to the cyclogram. The current peaks in the positive-going half cycle of the discharge (lower branch) reveal large current peaks, which is typical for the filamentary-dominated discharge regime [60]. The negative-going half cycle (upper branch) indicates similar trends at weaker magnitude, represented by a more diffuse and uniform discharge, as typical for surface discharges [64] (cp. also Benard and Moreau [65] for intensified charge-coupled device images of both discharge regimes). Furthermore, the boxes in the cyclogram indicate the dark periods of no discharge, which have been demonstrated to transfer no momentum to the air [27]. The comparison of the quiescent-air and U ∞ = 30 m s −1 current signals gives rise to a constant phase position of the discharge onset across the investigated U ∞ range.
Interestingly, the measured power consumption as well appears to be quasi independent from the airflow impact, where the observed ±2% fluctuations of the relative actuator performance Π PA = P A /P A,0 [38] around the quiescent-air result of P A,0 /L = 89.4 W m −1 are within the uncertainty margin of the utilized equipment; see figure 6. This result resembles the findings of Pereira et al [39] in co-flow configuration. For comparison purposes to earlier studies by Kriegseis et al [38], results for identical actuators and operating conditions are added to the diagram (▲). The mild performance drop of Ψ PA = 1 − Π PA < 3% for the range up to M ∞ = 0.1 suggests that the seeming contradiction falls back to the limited temporal and vertical resolution of the oscilloscope as used by Pereira et al [39] and in the current study. As such, the corner in the cyclograms (where the discharge ceases [59]) are not sufficiently well resolved to determine the enclosed area within an uncertainty margin below 1%.
However, the near-proportional relationship of actuator power and resulting production of force and thrust under quiescent air conditions [25,[65][66][67][68] provides evidence that small changes of power consumption of a few percent due to external  airflow similarly only lead to mild changes of the resulting body-force magnitude.

Body-force estimation
Body-force estimations on the grounds of experimentally obtained planar velocity information either build upon the 2D NSE where ρ is the fluid density, µ is the dynamic viscosity of the fluid, the added source terms f x , f y refer to the body-force components along the indicated space coordinate and ω = ω z = ∂v/∂x − ∂u/∂y is the (out of plane) vorticity of the 2D flow field. Since either approach has more unknowns (on the left side) than equations, Wilke [22] estimated that the pressure-gradients ∂p/∂x i are orders of magnitude weaker as compared to the respective body-force components f i and accordingly proposed to exclude the pressure term from (1). Likewise, Albrecht et al [23] argued that the the curl of the force is predominated by ∂f x /∂y and proposed to exclude the counter-part ∂f y /∂x. Both simplifications have been verified [25,26] and have, therefore, been proven to be valid assumptions in a quiescent-air environment. Earlier studies [26,41] furthermore demonstrated that the horizontal component f x of the force dominates the problem by orders of magnitude. Consequently, only equation (1a) is considered for the unsteady NSE-approach, which reduces to under absence of the pressure term. The cumulative integration of equation (2) along y from ∞ down to the wall [23] and only consideration of the f x contribution leads to for the unsteady VE-approach. As a technical note, it is important to mention that a moving average filter (7 × 1) has been applied to the raw velocity data to minimize the typical scratchy pattern [26], which results from the third derivatives of the velocity required in equation (4). The contribution due to local acceleration (acc.), convection (conv.) and diffusion (diff.) are indicated in either approach. The acceleration term has been approximated for all phases ϕ via central differences according to with m = 1, …, N ϕ (N ϕ : phase resolution of discharge cycle) and where ∆ϕ corresponds to the considered phase-to-phase spacing [69] and ϑ represents either u or ω [27]. To additionally test the airflow influence on the quasi-steady force fields, the acceleration term (acc.) is furthermore eliminated from equations (3) and (4) for the accordingly time-averaged force estimates.
Furthermore, the force fields are integrated over the control volume (CV) in order to achieve a measure for the possible drop of the force magnitude where F = F x has already been implicitly applied to the equation, since the body force is predominated by the horizontal force component and thus F x is considered the only meaningful contributor to the momentum-transfer magnitude [22,25,41]. The such derived integral body-force values are also compared to more robust phase-resolved momentumbalance-based force estimations, where the approaches of Baughn et al [36] and Debien et al [49] are combined to determine the body force from the boundary-layer flow.
For quiescent-air environment, Debien et al [49] identified the acceleration term inside the CV as the major contributor for an accurate estimation of the phase-resolved body force F from the momentum balance ρV reduced method, Kotsonis et al. [24] ∂u j ∂t dV + ρS u j u i n i dS Versailles et al. [42] = −S pn j dS +S τ ij n i dS Versailles et al. [42] + F j , as recently summarized by Hehner et al [69] with reference to the respective contributors [24,42]. As for actuated boundarylayer flows, Baughn et al [36] determined the time-averaged body force from the momentum-balance differences between actuated (PA) and reference (0) cases in streamwise direction (j = x). Further distinction of acceleration force F ACC , shear force F S and thrust F T as contributors leads to the total force Note that consideration of a sufficiently large 2D CV allows to neglect possible local pressure changes (cp. Kotsonis et al [24]) in equation (8) and F ACC,0 is zero as per definition for a steady reference flow and is accordingly excluded as well. Further approximation of the wall shear stress τ w ≈ µ∆u/∆y| 1.IA at the first IA above the wall [25,42] leads to for the body-force magnitude, which will be referred to as the CV approach for the remainder of the present work. Recall from above that the acceleration term cancels out for the timeaveraged considerations.

Surface-drag determination
To provide a more generalized measure for the shearcontribution analysis, a normalized version of the integrand of the last term in equation (9) according to is proposed as a valuable measure to evaluate the changing impact of self-induced drag [39,45] across the considered velocity range. Note that the normalization scales the wall shear stresses to the range 0 ⩽ S ⩽ 1 as per definition. The limits of S = 0 and S = 1 correspond to cases where no and all shear at the surface results from the effect of self-induced drag.

Performance quantification
A quantitative measure for the fluid-mechanic power P FM is determined as which has been proposed by Giepman and Kotsonis [47] for PA-power estimations. Note that both phase-resolved and time-averaged integral forces F are determined from equation (6), whereas only the time-averaged fluid-mechanic power P FM is extracted via equation (11). Finally, further evaluation of the determined PA performance with a more universal figure of merit is possible, if the ratio of either quantities, i.e. fluid-mechanic power P FM and integral body force F, together with the (electrical) actuator power P A is considered to provide measures for the fluid-mechanic effectiveness and efficiency respectively [48].

Results and discussion
The body-force estimation for the time-averaged velocity fields based on the previously introduced methodology (see section 4.1) is discussed in section 5.1. The role of the unknown pressure gradient for the NSE and CV approaches will be assessed in section 5.2. A more detailed elaboration of spatial phase-resolved body-force distributions and magnitudes follows in section 5.3. The high phase resolution N ϕ = 24 of the discharge cycle for the present experiment further raises the question for the effects of N ϕ on the previous body-force analysis, which is addressed in section 5.4. Finally, the performance of the PA under external airflow impact will be evaluated in section 5.5.

Time-averaged body force
First, the computed body-force fields for quiescent-air conditions are shown in figure 7(a) to provide a comparison to earlier PIV-based force-estimation efforts [25][26][27]. In agreement with these earlier investigations, the similarity of the forcedistribution patterns f x (x, y) as calculated with the NSE and VE approaches confirms the validity of the above-recapped respective assumptions included in equations (3) and (4) under absence of external airflow impact. Furthermore, integration of either force distribution inside the indicated CV (---) according to equation (6) and additional comparison with the CV approach (9) reveal comparable force magnitudes to earlier studies with similar operating conditions [25,27]. The influence of external airflow impact for increasing freestream velocity U ∞ on the body-force distribution, however, reveals dissimilar outcomes for NSE and VE approaches. For NSE (see figures 7(b) and (c), top row), the body force appears to vanish as the magnitude of the main part for quiescent air reduces from 12 kN m −3 (U ∞ = 0) to < 4 kN m −3 (U ∞ = 20 m s −1 ). Furthermore, domains of f x < 0 start to grow both in space and magnitude. In contrast, the body-force fields computed with the VE approach for U ∞ > 0 (see figures 7(b) and (c), bottom row) undergo a deformation, where a f x > 0 contribution grows on top of the observed bulk-force domain (y < 1.5 mm) in quiescent air. Locally, a f x < 0 region arises more downstream, when U ∞ increases from 0 to 20 m s −1 . The occurrence of the stressed f x > 0 part is unexpected, since the body force is exerted in proximity of the wall. While the influence of pressure gradients on the VE approach (see equation (4)) was eliminated, the findings raise the question, if negligence of ∂f y /∂x in equation (2) plays an adverse role for the body-force computation f x via equation (4). In that context, Kriegseis et al [25] have demonstrated that within the bulkforce domain, i.e. f x ⩾ 0.1f x,max [25], the expression holds true. However, above that domain, corresponding to the f x > 0 contribution in figure 7 (VE), the ratio in (14) is about 1 and, therefore, might impact the VE approach. Further indications of the observed body-force deformation are also given by the DNS study of Dörr and Kloker [41], where the authors highlight a significant f x > 0 body-force contribution in the same location as for the present results (figures 7(b) and (c), VE). Unveiling the cause for these findings, however, requires more investigations, as the mandatory assumption for the VE approach (see equation (2) and corresponding discussion) cannot be verified for U ∞ > 0.
To further study the global effects of external airflow impact on the body force, the integral magnitudes F are compared for

NSE
, VE (•) and CV approaches in figure 8. Among them, the consistent body-force magnitude for F(U ∞ = 0) confirms validity of the made assumptions. Saliently, the NSE and CV results reveal a significant drop of F, which even persists across the sign flip to F < 0 values. These outcomes underlie no physical mechanism involved, which might lead to such reversal of force direction. Moreover, this result retroactively provides an experimental validation of the DNSbased conclusion by Dörr and Kloker [41], i.e. that exclusion of the pressure gradient ∂p/∂x from the NSE approach (3) is an oversimplification of the problem unless applied to quiescent-air conditions. As such, the assumption of negligence of the pressure gradient under plasma actuation in presence of an external airflow will be revisited in section 5.2.
The VE (•) results, in contrast, reveal a mild F reduction once the discharge encounters external airflow. Interestingly, however, the mildly negative slope ceases beyond free stream velocities U ∞ = 10 m s −1 (correspondingly U 0 ⩾ 5.59 m s −1 ), which requires further elaboration of the data to examine the hypothetical changes of self-induced drag [39,45]. Therefore, the surface distribution of the normalized shear ratio S(x, U ∞ ) is shown in figure 9 along the wall boundary of the CV (cp. figure 7) for all investigated velocities U ∞ .
Obviously, for velocities U ∞ > 0 the shear ratio S drops immediately upstream of the actuator, which is due to the absence of discharge at x < 0. However, the impact of the momentum transfer on the effect of self-induced drag also decreases considerably for increasing free-stream velocities in the discharge region (x > 0).
Since the boundary-layer profile is steeper for faster flow speeds in vicinity of the discharge, the direct modification of the velocity gradients in immediate wall proximity becomes weaker and the footprint of the imparted momentum likewise shrinks for higher free stream velocities. This effect of shearratio reduction, therefore, counteracts any observable airflow impact on the electrical performance measures [38] and, in turn, explains the near-constant net-force in the flow [39].
Above U ∞ = 10 m s −1 , the force results from the VE approach indicate a clear trend towards the quiescent-air F value, while a sudden increase of F for U ∞ = 30 m s −1 (VE) is observed, which is unexpected as P A was constant (cp. figure 6). Recall, that the uncertainty margin of the PIV data increases with U ∞ due to increasing near-wall velocity gradients in the IAs. Moreover, the increasing impact of motion blur for the fastest investigated flow speeds influences the result. Since the VE approach furthermore implies third spatial derivatives of the recorded velocity information, the increasing uncertainty for faster flow speeds accordingly propagates heavily into the VE results. Therefore, the data for U ∞ = 30 m s −1 are only considered as trends rather than precise numbers throughout the discussion. The significance of F ≈ constant, in turn, is anticipated to hold for further increase of U ∞ to 30 m s −1 . Accordingly, the -symbol was added to the graph of figure 8 to indicate this hypothetical yet more realistic trend for the force, which is chosen identical to F(U ∞ = 20 m s −1 ).

Remarks on the unknown pressure gradient
The force-estimation approaches NSE (3) and likewise CV (9) rely on the necessary assumption that the pressure gradient ∂p/∂x is negligible [22,24] when computing the spatial distributions f x and integral magnitudes F of the body force, respectively. Accordingly, the body force exerted by the PA in quiescent air mainly drives momentum transfer, which results in a 2D velocity field (cp. figure 2(a)). The VE approach is instead independent of pressure variations in the flow field. As a consequence, the comparison of F in figure 8 in quiescent air (U ∞ = 0) asserts the unknown pressure gradient to be negligible (cp. also Kriegseis et al [25] and Dörr and Kloker [41]), showing excellent agreement between NSE, VE and CV approaches.
Superposition of an external flow has yet introduced significant and rather unexpected differences for F between the force-estimation approaches, as illustrated in figure 8. Considering the simultaneously monitored electrical discharge characteristics (see figure 5) and PA power consumption P A (see figure 6) to remain approximately constant, the exerted body-force magnitude F must be conjectured to undergo the same trend. However, F drops significantly for the NSE and CV approaches and even accommodates nonphysical values of F < 0 for U ∞ > 15 m s −1 . As further clarified in figure 7, the NSE-estimated body force f x (x, y) deforms substantial and even turns partly negative. These experimental findings are in strong agreement with the DNS study of Dörr and Kloker [41], where the authors have implemented a spatial body-force distribution f x (x, y), to compute the resulting velocity and pressure fields in a laminar M ∞ = 0.1 boundary-layer flow. The subsequent reproduction of the body-force field f x (x, y) via both NSE and VE, likewise revealed a large domain of f x < 0 for the NSE approach. Retrieving the pressure field from the DNS, Dörr and Kloker [41] attributed the modified solution of applying the NSE approach in a boundary-layer flow (U ∞ > 0) to the locally induced pressure gradient, driven by the exerted body force.
Within the scope of the present investigation the effect of fluid-flow heating on the applied body-force approaches (see equations (3) and (4)) and on the resulting spatial distributions was verified as insignificant [70].
The outcomes of the current work clearly manifest the role of the unknown pressure gradient as too significant to be neglected in the NSE approach. Therefore, either of the timeaveraged body-force fields shown in figures 7(b) and (c) for NSE represent a combination of subtracting f x and ∂p/∂x, where the true body-force field is deformed and becomes indistinguishable from the existing pressure gradient. For the CV approach, likewise, increasing U ∞ increasingly disrupts uniformity of pressure on the borders of the actuator-enclosing CV. Following the recommendations of Dörr and Kloker [41] for NSE, both NSE and CV approaches are to be disregarded in external airflow conditions as the inherent pressure assumption is violated. Consequently, in such conditions only the VE approach, in agreement with the findings of Dörr and Kloker [41], is rendered a meaningful method for spatial and integral body-force estimation. It is, however, to be noted that further verification of the f x > 0 contribution in figure 8 (VE) (see also section 5.1) above the bulk-force domain will be required in the scope of future studies.

Phase-resolved body force
The time-averaged analysis of the PA body force is expanded to a phase-wise comparison, where first the integral phaseresolved magnitude F(ϕ) and involved contributions separated by single terms in equations (3), (4) and (8) will be discussed in section 5.3.1. Subsequently, the spatial distributions of the body force along the discharge cycle will be presented in section 5.3.2 for the VE approach only, as delivered results by the other methods (NSE and CV) are related to nonphysical force values.

Integral magnitude.
The phase-resolved integral magnitudes F(ϕ) of the exerted body force are contrasted in figures 10(a)-(c) for the NSE, VE and CV approaches, respectively. The total force magnitude F(ϕ) (first row in figure 10) was computed through integration of f x (x, y, ϕ) according to equation (6) for NSE and VE, or via direct application of equation (9) for CV. The phase-wise progression of F(ϕ) reflects the periodic and unsteady discharge character, as previously shown by Benard et al [26], Kuhnhenn et al [27] and Debien et al [49].
The dark period [59], related to F = 0, after the positivegoing half cycle (0.5π < ϕ < π) is well captured by the applied force-estimation approaches. The second dark period (1.5π < ϕ < 2π), however, gives rise to a small negative force F < 0 for NSE and CV, which remains unclear. In agreement for NSE, VE and CV, the results of F(ϕ) during dark periods are deemed independent from U ∞ . Therefore, the pressure gradient neglected for NSE and CV plays an insignificant role during dark periods, considering F = 0. Apart from U ∞ = 30 m s −1 , the independence of U ∞ holds further true when the filamentary discharge is active (−0.2π < ϕ < 0.6π or 1.8π < ϕ < 2.6π ). In contrast, the integral force magnitude drops significantly with increasing U ∞ for NSE and CV, once the diffuse discharge commences (0.8π < ϕ < 1.6π ). In summary, the identified behavior in figure 8, rendering NSE and CV illposed for U ∞ > 0, mainly arises from the negative-going half cycle, which is representative of a strong F > 0 contribution to the time-averaged body force.
To further identify the impact of local acceleration, convective acceleration and viscous diffusion on the force magnitude, specific contributions are shown in the second, third and fourth rows of figure 10, respectively. Independent from U ∞ , the temporal average of local acceleration remains below O(10 −15 ) mN m −1 , which demonstrates that a time-averaged flow-field analysis is sufficient to gain quasi-steady net forces.
For NSE and CV, the results of F(ϕ) (first row in figure 10) are equivalently influenced by both local and convective acceleration terms (note the order-of-magnitude difference between the given terms). However, as outlined in section 5.2, the first row in figure 10 represents an indistinguishable combination of the spatial integration of f x and ∂p/∂x. In contrast, for VE the dismantled integral force magnitudes F into individual terms (see equation (4)) delivers equivalent results for local and convective acceleration terms for the investigated U ∞ range; note that U ∞ = 30 m s −1 was considered an outlier previously (see section 5.1). The magnitude of viscous diffusion, for clarity, as shown in the bottom row in figure 10 is O(10 −2 ) lower than F and remains insignificant.
While NSE and CV approaches have to be avoided for estimation of the exerted body force by the PA in presence of an external airflow, spatio-temporal body-force distributions are discussed in the subsequent part alongside the VE approach to provide further spatio-temporal insights.

Spatial distribution.
The temporal evolution of the spatial body-force distribution based on the VE approach is shown in figures 11 and 12 by means of phase-resolved force fields f x (x, y, ϕ) for the positive-and negative-going discharge half cycles, respectively. The first column in each figure represents a phase position ϕ during the dark period (as indicated by ) prior to the correspondingly presented discharge half cycle. Note that U ∞ for the computed force fields increases from top (0) to bottom (30 m s −1 ), as depicted on the right of the figures. Furthermore, in contrast to previous works [25,26], the FOV was cropped to −1 ⩽ x ⩽ 5 mm and 0 ⩽ y ⩽ 1.5 mm, so as to give a detailed picture of the relevant domain of the exerted body force.
The force fields in the first column of figure 11 (ϕ = 22π/12) clearly refer to the dark period, where the integral magnitude was F = 0 (cp. figure 10(b)). The subsequent ignition of the filamentary discharge (second column, ϕ ⩾ 1π/12) triggers a significant f x < 0 component (cp. also [26,27]). Moving along ϕ the force magnitudes become weaker, approaching zero for the fourth column (ϕ = 5π/12). The U ∞ -dependent development of the spatial distribution of f x Figure 11. Phase-resolved body-force distribution fx(x, y, ϕ) based on the VE approach for increasing free-stream velocity U∞. Selected phase-wise positions refer to dark period (ϕ = 22π/12, ) and subsequent positive-going half cycle (ϕ = 1π/12, 3π/12 and 5π/12), i.e. filamentary discharge (negative-going half cycle see figure 12). U∞ is indicated on the right of the figure. during the filamentary discharge obviously experiences a shape change for increasing U ∞ , where the most salient yet peculiar modification is an intensifying tilt of the curved force pattern in the second (ϕ = 1π/12) and third (ϕ = 3π/12) column towards the surface. Furthermore, for increasing U ∞ this part of the body-force contribution considerably weakens in magnitude.
The force fields in the first column of figure 12 (ϕ = 10π/12) show-likewise as for figure 11 (first column, ϕ = 22π/12) -the dark period prior to the diffuse discharge regime or negative-going cycle, rendering F = 0 (cp. figure 10(b)). The ignition of the discharge triggers a very strong f x > 0 component (second column, ϕ = 13π/12), much larger and more pronounced than for the filamentary discharge regime (cp. also [26,27]). For the phase-wise development of f x , the shape of the body force in fourth column (17π/12) resembles the one of the second column in the filamentary discharge regime (figure 11, ϕ = 1π/12) for opposite signs of f x . The previously discussed tilt of the curved body-force pattern towards the surface likewise occurs in the diffuse discharge regime for increasing U ∞ . The diminishing effect on the force magnitudes during the filamentary discharge regime (figure 11) is obviously much weaker during the diffuse discharge. In terms of differences of the force-magnitude variation between either discharge cycle, the phase-wise shape change of the force with increasing U ∞ , consequently implies the resultant modified time-averaged body-force fields for different U ∞ (see figure 8, VE).

Influence of phase resolution on the estimated integral force magnitude
The strong temporal variation of the exerted body force as presented in figures 10(b), 11 and 12 implies the accuracy Figure 12. Phase-resolved body-force distribution fx(x, y, ϕ) based on the VE approach for increasing free-stream velocity U∞. Selected phase-wise positions refer to dark period (ϕ = 10π/12, ) and subsequent negative-going half cycle (ϕ = 13π/12, 15π/12 and 17π/12), i.e. diffuse discharge (positive-going half cycle see figure 11). U∞ is indicated on the right of the figure. of the phase-resolved body-force magnitude F(ϕ) to be very sensitive to the selected phase resolution N ϕ of a discharge cycle of the PA. The current phase resolution of N ϕ = 24 allows to further analyze the effect of reduced phase resolutions on F(ϕ) (cp. also Benard et al [26] and Kuhnhenn et al [27]). Figure 13(a) provides the comparison of F(ϕ) in quiescent air (U ∞ = 0) for phase resolutions of N ϕ = 24 for the VE approach. The phase-to-phase spacing ∆ϕ, accordingly, is inversely proportional to N ϕ . The force magnitudes F(ϕ) in figure 13(a) were estimated for ϕ 1 = 0 (see figure 5 for orientation). Obviously, both peaks F min and F max reduce in magnitude for decreasing N ϕ . As further illustrated by figure 13(b), reducing N ϕ from 24 to 4, yields six options of choosing ϕ 1 , which in turn leads to a horizontal shift of the peak magnitude F max . Such a shift consequently introduces an uncertainty on the estimation of the phase-wise peak-force location ϕ(F max ). The errorbars in figure 13(c) indicate this increasing uncertainty for decreasing N ϕ .
The impact of N ϕ on the time-averaged body-force magnitude F is shown in figure 14. The resulting deviation compared to N ϕ = 24 remains below 2% for N ϕ ⩾ 6 and appears to be dominated by the above-mentioned measurement uncertainty for increasing velocities. Lower phase resolution of N ϕ = 4, however, reveals considerably increasing errors for increasing U ∞ , culminating in ≈7% error for U ∞ = 30 m s −1 . While the phase-wise convective acceleration shown in figure 10(b) remains largely constant for U ∞ < 20 m s −1 , a variation occurs for larger U ∞ .
It is furthermore worth to mention that the phase-resolved force magnitude F(ϕ) does not immediately mimic a harmonic oscillator due to the different characters of either half cycle. That is, a sufficient sampling frequency for either part of the signal is N ϕ f ac /2, which leads to an effective Nyquist limit for  N ϕ = 4. As such, neither half cycle F(ϕ) < 0 and F(ϕ) > 0 of the phase-resolved force magnitude F(ϕ) can be reconstructed accurately for phase resolutions N ϕ < 4, which explains the sudden increase of the error for the time-averaged force magnitude F(N ϕ = 4).
In summary, an accurate representation of the phaseresolved body force requires sufficiently large N ϕ , while complexity of the experiment can be reduced if pure time-averaged body forces are considered.

Actuator performance
Both fluid-mechanic effectiveness η * FM (see equation (12)) and efficiency η FM (see equation (13)) of the considered PA configuration are determined for the investigated U ∞ range. In quiescent air, the effectiveness of a PA yields an optimum peak-to-peak voltage V pp for a given operating frequency f ac [48]. Based on the power consumption P A and the integral magnitudes F shown in figures 6 and 8, respectively, the effectiveness η * FM is presented in figure 15 (•). These results clarify that the PA remains as effective as in quiescent air under external airflow impact. The deviation of the integral force magnitude F(U ∞ = 30 m/s), which is included in the numerator of η * FM (cp. equation (12)), terminated the trend of F(U ∞ ) ≈ constant (cp. figure 8). However, since this outcome is rather attributed to the increasing effect of measurement uncertainty (see section 5.1) the hypothesized trend is also added to the diagram similar to figure 8 for clarity.
Since the fluid-mechanic power P FM (see equation (11)) is determined from the supplied body force multiplied with the velocity of the surrounding fluid of the (local) momentum transfer, the observed quasi-constant force magnitude directly leads to increased fluid-mechanic power for increasing flow speeds. Given that electrical actuator power and supplied force scale linearly, the such increased P FM , therefore, immediately increases the fluid-mechanic efficiency η FM (see equation (13)) of the PA. That is, more power is pumped into the airflow at the same expense of consumed power P A . Typically, the fluid-mechanic PA efficiency η FM is O(10 −1 )% [47,48,50], where the local flow velocity was entirely generated by the PA momentum transfer. The increase of U ∞ to up to 20 m s −1 , consequently, involves an increase of η FM by a factor 3, as shown in figure 15 (•). Recall however, that the •-symbol for U ∞ = 30 m s −1 is considered an outlier (see discussion of figure 8), representing an overestimation of the increase of η FM . As above, this data point is supplemented with the hypothesized value as determined from a quasiconstant body force. Therefore, exerting a constant body force improves η FM , as one of the individually considered chain links in the power-flow diagram of flow-control efficiency [48], and is in favor of maintaining the fluid-mechanic efficacy of the PA. However, the global implications on AFC efficiency are not as easily determined, considering that the flow-control effect decreases for PA AFC at higher airflow speeds, while e.g. power savings scale with the external flow velocity. Even though further increase of U ∞ is conjectured to further increase η FM , it has to be mentioned, however, that P A is expected to drop significantly for higher flow speeds [38,40], which requires to re-adjust/increase the input power of the overall system. The present findings, as yet, cannot be extrapolated to higher free-stream velocities. Once the external flow velocity is comparable to the drift velocity 6 of the ions, an impact on the performance and, thus, a change on the PA forcing mechanism might occur [39].
The phase resolution N ϕ of gathered velocity information along the plasma discharge cycle is of utmost importance for the accuracy of computed body forces. According to the tested scenarios, N ϕ ⩾ 24 is recommended. Lowering N ϕ has shown to adversely affect the estimation of integral force magnitudes and to give rise to an uncertainty of the phase position related to minimum/maximum forces. The worst case of N ϕ < 4, has to be avoided, since then the discharge phenomenon of an entire half cycle is smeared into a single ∂u/∂t estimate, rendering any gained body-force results meaningless.
If for a study only time-averaged or phase-resolved integral force magnitudes are requested, the acquisition of body-forcerelevant data in quiescent air is sufficient. In contrast, either supported by FCT through IDMEC, under LAETA, Project UIDB/50022/2020. Gonçalo Coutinho acknowledges the PhD scholarship 2021.04780.B D attributed by Fundação para a Ciênciae Tecnologia (FCT).