Simulation of active flow control actuator using CFD with application to rotor blade vibration reduction

Paper 69

SIMULATION OF ACTIVE FLOW CONTROL ACTUATOR USING CFD WITH APPLICATION TO

ROTOR BLADE VIBRATION REDUCTION

Ryan P. Patterson*_{and Peretz P. Friedmann}†

rppat@umich.edu peretzf@umich.edu

Department of Aerospace Engineering, University of Michigan, Ann Arbor, MI, USA Yuehan Tan*_{and Ari Glezer}‡

Woodruff School of Mechanical Engineering, Georgia Institute of Technology, Atlanta, GA, USA

Abstract

A computational fluid dynamics (CFD) model is developed to determine the unsteady aerodynamic effects of active flow control implemented by combustion-powered actuation (COMPACT) on a two-dimensional airfoil. Previous work with COMPACT actuators employed pulsed-jet actuation at moderate to high angles of attack to control flow separation and dynamic stall. However, the focus of the present study is to imple-ment actuation in the moderate to low angle of attack range suitable for helicopter rotor vibration control. At low angles of attack, the effect of actuation on the sectional aerodynamic forces of the airfoil diminishes. Therefore, modifications to the COMPACT actuator configuration are considered to enhance the actuation authority for the relevant operating range. Modifications include placing a ramp on the airfoil upstream of the actuator and relocating the actuator near the airfoil trailing edge. These are considered using both CFD simulations and wind tunnel experiments, which are compared to validate the CFD model. Finally, a surrogate-based reduced-order modeling technique is described to address the high computational cost of the CFD simulations. The reduced-order model is used to accurately reproduce full-order CFD results for the unsteady changes in lift, moment, and drag due pulsed actuation on a static airfoil. Furthermore, the results show that the reduced-order model represents a feasible method for representing the unsteady aerodynamic effects of COMPACT actuation. This will be refined in future work and implemented in a com-prehensive aeroelastic code for helicopter rotor vibration reduction.

NOMENCLATURE

c

Airfoil chord length

C

l

, C

m

, C

d 2D lift, moment, and drag coefficients

C

L

, C

M

, C

D Experimentally measured lift,

mo-ment, and drag coefficients

C

q Volume flow rate ratio jets/freestream

∆C

x ,qs Quasi-steady increment in

C

x due to

repeated actuation;

_x

_{≡ l , m,}

or

_d

*_{Ph.D. Candidate} †

François-Xavier Bagnoud Professor of Aerospace Engineer-ing

‡

Professor, Woodruff School of Mechanical Engineering

Copyright Statement

The authors confirm that they, and/or their company or or-ganization, hold copyright on all of the original material included in this paper. The authors also confirm that they have obtained permission, from the copyright holder of any third party material included in this paper, to publish it as part of their paper. The authors confirm that they give per-mission, or have obtained permission from the copyright holder of this paper, for the publication and distribution of this paper as part of the ERF proceedings or as individual offprints from the proceedings and for inclusion in a freely accessible web-based repository.

f

act Actuation frequency, Hz

F

+ Non-dimensional reduced actuation frequency,

f

act

c /U

∞

k

Turbulent kinetic energy

m, n

Number of previous time steps re-quired to account for time-history ef-fects in SBRF

M

∞ Freestream Mach number

P

act Actuator combustion chamber pres-sure

P

r Pressure ratio of actuator combustion chamber,

P

act

/P

∞

P

r,peak Peak pressure ratio of actuator com-bustion chamber

P

∞ Freestream static pressure

Re

c Reynolds number based on chord length,

U

∞

c /µ

∞

t

Dimensional time

t

0 Initial time of actuator pulse in CFD simulations

T

act,pulse Actuator pulse duration in CFD simu-lations

u

Representative actuator input signal for SBRF model of COMPACT actua-tion

u

Vector of external inputs

U

∞ Freestream flow speed

x

State vector for SBRF model

ˆ

x

Estimated state vector for SBRF model

y

“Exact” solution for a given output

variable; e.g.,

C

l,

C

m, or

C

d

ˆ

y

Approximate solution for a given out-put variable

α

Angle of attack

 (ˆ

y , y )

Normalized mean absolute error of

ˆ

y (t)

relative to

y (t)

µ

∞ Molecular dynamic viscosity

µ

t Turbulent dynamic viscosity

Φ

Non-linear mapping function

ˆ

Φ

Surrogate mapping function

ˆ

Φ

qs Quasi-steady surrogate mapping function for SBRF model

ω

Specific rate of dissipation of turbu-lent kinetic energy

1. INTRODUCTION, BACKGROUND, AND OBJEC-TIVES

The principal thrust of research in active flow con-trol for rotorcraft applications has been to concon-trol dynamic stall, which is a major source of vibrations in a helicopter rotor during high-speed flight. Dy-namic stall is characterized by unsteady flow sep-aration and reattachment that result from the pe-riodic variations in effective angle of attack of the blade sections. For high forward flight speeds, the inboard sections of the retreating blade encounter large angles of attack that exceed the static stall an-gle of the blade. Near the maximum anan-gle of attack, a dynamic stall vortex is shed from the blade sur-face. This induces a large negative pitching moment on the blade that contributes to the overall vibra-tory loads on the rotor.

Active flow control has been implemented in wind tunnel experiments to alleviate the effects of dynamic stall on pitching airfoils. A typical ap-proach is to inject momentum from the airfoil sur-face when the airfoil is at moderate to high angles of attack.1,2,3,4 Previous studies demonstrated that actuation during the upstroke of the dynamic stall pitching cycle can suppress the formation of the dy-namic stall vortex and effectively delay lift-stall to higher angles of attack. Additionally, actuation dur-ing the downstroke can accelerate the flow reat-tachment process after the dynamic stall vortex has been shed. Note that the benefits of active flow con-trol were achieved in these studies because of the interaction withseparated flows.

Computational fluid dynamics (CFD) studies have also characterized the influence of active flow con-trol on dynamic stall and flow separation. Un-steady Reynolds-averaged Navier-Stokes (URANS) simulations have been implemented to produce engineering-level accuracy suitable for design sen-sitivity studies5,6 and optimization7 of the actu-ators. In Ref. 8, the accuracy of URANS simula-tions was compared with higher-fidelity LNS (limited numerical scales) simulations. The study demon-strated that the lower-fidelity URANS simulations were sufficient to reproduce experimental results for low frequencies of actuation (i.e., for

F

+

< 3.3

), whereas the LNS simulations were more accurate up to

F

+

_{= 10}

_{. In the literature, the actuation} fre-quency is commonly expressed in non-dimensional form:

F

+

= f

act

c /U

∞. Here the actuation fre-quency

f

actis scaled with the convective time scale,

T

conv

= c /U

∞, which is the time for disturbances to travel over the chord length of the airfoil. Large eddy simulations (LES) have also been implemented to produce highly accurate unsteady results for ac-tive flow control on a statically stalled airfoil.9 How-ever, the application of LES to active flow control modeling is limited, due its high computational cost. Although dynamic stall is an important consider-ation in rotorcraft aeroelasticity, it is not the only source of vibrations in a helicopter rotor. Active control studies for rotor vibration reduction (re-viewed in Ref. 10) have applied a systematic ap-proach to reducing vibrations, rather than target-ing only dynamic stall. As a result, the control ef-fort in these studies was applied primarily in tached flows, at low to moderate blade angles of at-tack. Among these studies, the AVINOR aeroelastic simulation code11 has been used for active vibra-tion control implemented by active flaps and micro-flaps.12 The AVINOR code is a comprehensive sim-ulation code for helicopter rotor vibration reduc-tion. It uses a reduced-order model, based on the rational function approximation,13,14 to model the unsteady aerodynamic loads during vibration con-trol. The computational efficiency and accuracy of the reduced-order model are critical for performing closed-loop active vibration control studies.

The overall objective of the current study is to develop an accurate representation of active flow control implemented by combustion-powered actu-ation (COMPACT)4. The specific objectives of this pa-per are:

1) Develop a CFD simulation methodology for modeling the unsteady effects of COMPACT actua-tion on a two-dimensional (2D) airfoil operating in the moderate to low angle of attack range that is relevant for helicopter rotor vibration control.

wind tunnel experiments for actuation on a static airfoil.

3) Develop a surrogate-based reduced-order model for representing the full-order CFD solutions, capturing changes in cross-sectional aerodynamic coefficients (lift, moment, and drag) due to actua-tion on the 2D airfoil.

These objectives are important steps towards evaluating active flow control for helicopter rotor vibration reduction when using a comprehensive aeroelastic simulation code.

2. CFD SIMULATION OF COMPACT ACTUATION The COMPACT actuator produces a high-velocity jet that is generated by a combustion reaction inside a cubic-centimeter-scale chamber. This causes flow to accelerate when exiting through a narrow slot. An overview of the experimental development and characterization of the actuators was presented in Ref. 15. In the current study CFD simulations are used to represent a VR-12 airfoil with COMPACT ac-tuators embedded in the airfoil surface. This con-figuration simulates a wind tunnel experiment at Georgia Tech that also represents a configuration used in Ref. 3.

2.1. Description of Flow Solver

Simulations of the COMPACT actuators are per-formed using two-dimensional (2D) RANS CFD sim-ulations. The code used in this study is CFD++,16 which solves the compressible unsteady RANS equations using a finite volume formulation. The CFD++ code is a general-purpose code that is ca-pable of modeling a variety of mesh types, in-cluding structured, unstructured, hybrid, patched, and overset grids through a unified grid method-ology. Spatial discretization is based on a second-order total variation diminishing scheme. For time-dependent simulations, a first-order implicit algo-rithm is implemented, with second-order accuracy achieved by dual time-stepping. Additionally, multi-grid convergence acceleration is utilized in the sim-ulations. Several different turbulence models are available in CFD++. The Menter

k

-

ω

shear stress transport (SST) model was selected for the current study, assuming a fully turbulent boundary layer.

2.2. Actuator Boundary Conditions

The pulsed jet is modeled in the CFD simulations by specifying a time-varying inlet boundary condi-tion for the total pressure within the actuator, as illustrated in Fig. 1. The inlet boundary represents

Combustion chamber (not modeled)

Pact(t)

Figure 1: Actuator inlet boundary condition for CFD simulation of COMPACT pulse actuation

a simplification of the internal flow, since the com-bustion chamber is not modeled. This approach to modeling the jet is similar to that used in a previous CFD study with COMPACT actuation.6

The time varying pressure within the actuator is represented by an analytical function that repro-duces the rise, peak, and drop in combustion cham-ber pressure observed in experiments,

(1)

P

act

(t) = P

r

(t)P

∞

,

where

P

∞is the freestream static pressure and the pressure ratio is,

(2)

P

r

(t) = 1+

P

r,peak

− 1

2

×



cos



2π

t

− t

0

T

act,pulse

− π



+ 1



.

The cosine function in Eq. (2) is used for

P

r

(t)

for

0

≤ t − t

0

≤ T

act,pulse, where

T

act,pulse

= 7

× 10

−4 seconds. Otherwise, the actuator inlet boundary is treated as a wall boundary condition.

Within the CFD++ simulation environment, this approach is implemented with a built-in boundary type. The “Time Varying Total Pressure and Tem-perature” boundary is a reservoir inlet boundary condition for which the stagnation pressure, tem-perature, and turbulence quantities are specified as functions of time.17 For the current simulations, only the time-varying pressure is specified, while the values of temperature (

T

) and turbulence vari-ables (

k

and

ω

) are held constant at the inlet bound-ary. The quantities

T

,

k

, and

ω

at the inlet bound-ary are prescribed to be equal to the values that are specified at the far-field boundary of the do-main. These correspond to a temperature of 300 K, a turbulence intensity of 1%, and a turbulent-to-molecular viscosity ratio (

µ

t

/µ

∞) of 50.

Figure 2: Computational grid for baseline clean VR-12 airfoil

For repeated pulse actuation, the actuator strength is frequency dependent, such that

_P

_r,peak decreases with increasing actuation frequency,

f

act. This dependence was characterized in Ref. 18, and is taken into account in our CFD simulations. For example, the peak pressure ratio of a single, iso-lated pulse is 3.3. Repeated actuation at 50 Hz (pulses/second) yields a peak pressure ratio of 2.3. In Ref. 18, data for

_P

_r,peak versus

_f

_actare available in 5-Hz intervals (up to 85 Hz) and are linearly inter-polated in the current study.

2.3. Grids for Modeling the Airfoil/Actuator Combination

An overset, structured grid approach is employed in the vicinity of the actuator in order to capture the pulsed jet with sufficient resolution. First, a base-line grid for the clean VR-12 airfoil was generated. The computational domain is a C-grid, as shown in Fig. 2, with the far-field boundary extending to 50 chord-lengths from the airfoil. A distance of 30–50 chord lengths is typically used in order to minimize the effects of numerical reflection from the far-field boundary.

Next, an overset grid for the COMPACT actua-tor was generated, consisting of the actuaactua-tor slot and the domain that immediately surrounds the ac-tuator. This includes a small region upstream (5% chord) and a larger region downstream (45% chord) from the actuator. The height of the overset grid ex-tends approximately 12% chord normal to the air-foil surface and overlaps a portion of the baseline clean grid, as shown in Fig. 3. For the configura-tion shown, the COMPACT actuator is located at 10% chord from the airfoil leading edge.

A grid refinement study was performed with the 2D grids in order to determine the mesh density

Figure 3: Overset grid for COMPACT actuator lo-cated at 10%-chord on VR-12 airfoil

required to accurately capture the unsteady pulse-actuation. Three levels of mesh refinement for the combined baseline and overset grids were tested – coarse, nominal, and fine. The coarse mesh consists of a baseline airfoil mesh of approximately 85,000 quadrilateral cells, combined with an overset actu-ator mesh of approximately 25,000 cells. The nom-inal mesh uses the same 85,000-cell baseline mesh and a finer overset mesh of approximately 40,000 cells. Finally, the fine mesh consists of a baseline mesh of 160,000 cells, combined with the 40,000-cell overset mesh.

The unsteady lift, moment, and drag coefficients (

C

l,

C

m,

C

d) were computed on the coarse, nomi-nal, and fine grids for a single actuator pulse. For this comparison, the airfoil was fixed at a static

₈

◦ angle of attack, with freestream flow conditions cor-responding to

M

∞

= 0.06

and

Re

c

= 535, 000

. The unsteady coefficients were compared over an interval of

5T

conv after the pulse (

T

conv

= c /U

∞). The relative error was computed using the follow-ing metric: (3)

 (ˆ

y , y ) = 100

×

1 N

P

N n=1

|ˆ

y (t

n

)

− y (t

n

)|

y

max

− y

min

.

This metric is the mean absolute error of

y (t)

ˆ

rel-ative to

y (t)

, where

y (t)

represents the “exact” response and

_{y (t)}

_ˆ

represents an approximate re-sponse. The quantity

y

is used to represent a single output of interest, and is replaced by

C

l,

C

m, and

C

d, respectively, and

N

is the number of time sam-ples from the CFD responses that are used to com-pute the relative error. In the grid refinement study, the responses calculated on the fine grid represent the “exact” solutions and the responses calculated on the coarse and nominal grids represent approx-imate solutions.

The unsteady errors of the coarse and nominal grids relative to the fine grid were calculated using Eq. (3) and are summarized in Table 1. The response of the nominal grid approximates the response of the fine grid within a reasonable margin of error: 1.9% in

C

l, 0.6% in

C

m, and 0.5% in

C

d. Therefore, the nominal grid is used for the remainder of the CFD computations in this study.

2.4. Modification of Actuation for Low Angles of Attack

Active flow control has been shown in the litera-ture to be effective for separated flows considered at moderate to high angles of attack. However, the principal region of interest for helicopter rotor vi-bration control is in attached flows at lower angles of attack. Reference 19 presented a clear example of the diminished authority of active flow control at low angles of attack (

α < 10

◦), when compared to higher angles of attack (up to

30

◦). Motivated by the need to increase the effectiveness of COMPACT actuators at low angles of attack, several modifica-tions of the basic configuration are considered.

The first modification to COMPACT consisted of placing a ramp upstream of the actuator. By placing a ramp with an aft-facing step on the airfoil surface, a localized region of “trapped vorticity” surrounds the actuator. This localized region of trapped vortic-ity is subsequently manipulated by COMPACT actu-ation to effect changes in aerodynamic forces. How-ever, it is noted that this actuation approach would only provide a unidirectional change (net reduction) in the lift force relative to the uncontrolled base-line. This modification was implemented in both the wind tunnel experiments and in the CFD simula-tions. The addition of the ramp was accounted for in the CFD simulations by generating a new over-set mesh for the actuator/ramp combination. This will be referred to as the COMPACT+Ramp configu-ration and is shown in Fig. 4 for a ramp height of 1.3% chord.

The second type of modification to COMPACT consisted of relocating the actuator. The overset grid method for generating the CFD meshes made the investigation of different locations straightfor-ward in simulations, because only the overset ac-tuator mesh needed to be regenerated for each new configuration. This flexibility of the simulations is not present in actual wind tunnel experiments. Therefore, the results examining alternate locations of the COMPACT actuator are obtained exclusively from simulations.

Figure 4: Overset grid for COMPACT actuator lo-cated at 10%-chord on VR-12 airfoil with added 1.3%-chord ramp

3. RESULTS

The CFD simulation results are presented in this section for COMPACT actuation on a static airfoil. First, the effects of actuation at low and high an-gles of attack are compared in order to illustrate the fundamental differences between the flow regimes. Next, the effect of adding a ramp to the airfoil to enhance the actuation at low angles is shown. Fi-nally, alternate locations for the COMPACT actua-tor are investigated to identify locations where the actuators are most effective. The CFD simulations are used throughout this section and are compared with the wind tunnel experiments where available.

3.1. Actuator Operation at Low and High Angles of Attack

The performance of the COMPACT actuator is de-pendent on the ambient flow conditions, and is ef-fective in controlling separated flow. When the base flow is fully attached the actuation can be utilized indirectly by regulating the presence of trapped vor-ticity concentration that can be engendered using a passive element. This is illustrated by simulat-ing pulse actuation at several angles of attack us-ing 2D simulations. The cases considered are

4

◦,

8

◦,

12

◦, and

16

◦, with freestream flow conditions corresponding to the experimental test conditions (

_M

_∞

_{= 0.06}

,

_Re

_c

_{= 535, 000}

). Before actuation,

α = 4

◦ and

8

◦ have fully attached flows under steady state conditions. At

α = 12

◦partial flow sep-aration is encountered near the airfoil trailing edge and at

α = 16

◦there is substantial flow separation. This behavior is identified from contour plots of the Mach number, provided in Fig. 5 for each case.

Table 1: Relative error of unsteady aerodynamic coefficients in 2D CFD simulations of a single pulse Grid # of cells baseline # of cells overset

 ˆ

C

l

, C

l

 ˆ

C

m

, C

m

 ˆ

C

d

, C

d

coarse 86,304 24,689 3.6% 1.6% 1.4% nominal 86,304 40,399 1.9% 0.6% 0.5% fine 158,404 40,399 – – – (a)α = 4◦ (b)α = 8◦ (c)_{α = 12}◦ (d)_{α = 16}◦

Figure 5: Mach number contours for VR-12 airfoil at

α = 4

◦

, 8

◦

, 12

◦

,

and

₁₆

◦, before COMPACT actua-tion;

M

∞

= 0.06

,

Re

c

= 535, 000

actuation was applied for each angle of attack. The actuation consisted of seven consecutive pulses ap-plied at a reduced frequency

F

+ of 1.67. The raw CFD data for the lift, moment, and drag coefficients due to actuation contain high-frequency oscillations after each pulse. Therefore, in order to compare the quasi-steady trends in the unsteady data, the phase average of each response was computed over inter-vals of

1/F

+. Figure 6 illustrates how the computed phase average compares with raw CFD data for the lift coefficient. Note that in the figure the raw CFD data are clipped to be within the prescribed axis bounds. In this case, the range of

∆C

l

(t)

is -0.8 to 0.5. Also shown at the bottom of the figure is the actuation signal. Here the magnitude of the actua-tion signal is arbitrary and is plotted only to show the timing of each pulse relative to the response.

The phase average highlights the quasi-steady nature of the lift response. As each pulse is applied, the lift coefficient incrementally increases towards a converged value. After the final pulse, the lift co-efficient gradually returns to the initial steady-state

Actuation Signal

Figure 6: Phase average of lift coefficient data from 2D CFD simulation of COMPACT actuation at

F

+

=

1.67

,

α = 8

◦

value that existed before the actuation. For this par-ticular flow condition, the post-actuation relaxation takes place over approximately

_8T

_conv.

The phase averages of the responses for

α =

4

◦,

8

◦,

12

◦, and

16

◦ are plotted in Fig. 7 in terms of the increments in lift, moment, and drag coeffi-cients. The most straightforward comparison is ob-tained by comparing the responses for the lift co-efficient. The peak value of lift due to actuation in-creases greatly from

_{α = 12}

◦ to

_{α = 16}

◦. Ad-ditionally, the post-actuation relaxation period for

α = 16

◦ takes place over

30T

conv– approximately three times longer than any of the lower angles. Fur-thermore, the overshoot phenomenon in the lift re-sponse for

α = 16

◦ is not observed at the lower angles of attack.

The greater receptivity of the flow to actuation at

α = 16

◦is due to the natural instability of the sep-arated flow domain. In previous experiments with COMPACT actuators,20the transient effects of actu-ation in a fully separated flow were characterized via particle image velocimetry (PIV). In the study, PIV frames were used to identify that the immedi-ate effect of a pulse from the actuators is to sever the separated vorticity layer. Consequently, a large-scale vortex is shed which causes the progressive reattachment of flow to the upper surface of the air-foil. The result is a dramatic increase in the

circula-Actuation Signal

(a) Lift coefficient

Actuation Signal

(b) Moment coefficient

Actuation Signal

(c) Drag coefficient

Figure 7: Phase-averaged increments in lift, moment, and drag coefficients due to COMPACT actuation for

F

+

= 1.67

at

_{α = 4}

◦,

₈

◦,

₁₂

◦, and

₁₆

◦

tion of the airfoil. After the vortex passes the airfoil trailing edge, the attached vorticity layer gradually lifts off the airfoil surface, thus reestablishing the baseline separated flow domain. For low to moder-ate angles of attack, the lack of a substantial sepa-rated flow domain weakens the influence of the ac-tuation.

3.2. Influence of COMPACT+Ramp for Enhanced Actuation

A ramp was added to the airfoil/actuator configu-ration to alleviate the diminished performance at low angles of attack. The rear portion of the ramp is an aft-facing step, which induces local flow sep-aration near the COMPACT actuator and generates a region of “trapped vorticity”. Although in the ab-sence of COMPACT actuation the local separation bubble formed downstream of the passive ramp can lead to a decrease in lift relative to the base-line airfoil, this reduction can be subsequently aug-mented by COMPACT actuation. The actuation can increase the lift performance such that it is nearly at the level of the base airfoil. The magnitude of the lift may be optimized by the streamwise posi-tion of the actuators, although in the present ex-periments it is a bit lower than the lift of the base airfoil. Contour plots of the Mach number from CFD simulations, provided in Fig. 8, show the influence of a 1.3%-chord ramp height on the flow field. Both the original COMPACT configuration and the COM-PACT+Ramp configuration are depicted prior to ac-tuation.

The CFD simulations are compared with experi-mental results for unsteady lift, moment, and drag coefficients. It is important to note that the wind tunnel model used in the experiments is only equipped with the COMPACT actuators and ramps over the center 21% of the model span. This consti-tutes theactive section of the model, as the

remain-(a) COMPACT (b) COMPACT+Ramp Figure 8: Mach number contours from 2D CFD sim-ulation of VR-12 airfoil with COMPACT and COM-PACT+Ramp actuator configurations, before actua-tion;

_{α = 8}

◦,

_M

_∞

_{= 0.06}

,

_Re

_c

_{= 535, 000}

ing outer sections are clean airfoil sections. Because the CFD simulations are representative of a fully-actuated, infinite span, a direct comparison of the calculated 2D coefficients from simulations with the measured coefficients from experiments is not a one-to-one comparison. Therefore, an approximate 2D comparison is used to scale the experimental re-sults in an appropriate manner.

The 2D comparison consists of scaling the incre-ment in the measured coefficients by 1/0.21, as if COMPACT actuation were applied across the entire span of the wind tunnel model. To illustrate the scal-ing, denote the unsteady lift coefficient measured in experiments as

C

_LEXP

(t)

, the baseline, steady-state value of the lift as

_C

EXP

L,0 , and the unsteady increment as the difference between the two, i.e.,

∆C

_LEXP

(t) = C

_LEXP

(t)

− C

EXP

L,0 . The measured in-crement in lift is scaled to yield the approximate 2D increment in lift on the active section of the wind tunnel model:

(4)

∆C

_lEXP

(t) =

∆C

EXP L

(t)

0.21

.

the increment in lift obtained in the CFD simu-lations:

_∆C

CFD l

(t) = C

CFD l

(t)

− C

CFD l ,0 . The ex-perimental data for moment and drag coefficients (

C

_MEXP and

C

_DEXP) are also scaled for comparison with CFD results in a similar manner. The incre-ments in aerodynamic coefficients are all calculated relative to the original baseline COMPACT configu-ration.

To demonstrate the approximate 2D scaling, a 2D CFD simulation with the COMPACT+Ramp configu-ration is compared with experimental results. Thir-teen consecutive pulses were applied at a reduced actuation frequency

F

+ of 0.634. The flow condi-tions correspond to a fixed angle of attack of

8

◦with

M

∞

= 0.06

and

Re

c

= 535, 000

. The unscaled lift, moment, and drag coefficients are plotted in Figs. 9a-c and the scaled, relative coefficients are shown in Figs. 9d-f.

It is evident from Fig. 9 that much better agree-ment is obtained between the CFD simulations and the experiments when the 2D scaling is applied to the experimental data. This is particularly evident in the lift coefficient, as both quasi-steady and tran-sient effects are accurately represented. The con-verged value of the lift due to the repeated actu-ation is denoted as the quasi-steady lift increment,

∆C

l ,qs. From the CFD simulations

∆C

l ,qsCFD

=

−0.163

and from the scaled experimental data

_∆C

EXP

l ,qs

=

−0.165

. This demonstrates remarkable agreement for the combined effects of the ramp and the ac-tuation. Similarly, the quasi-steady increments in moment and drag coefficient were calculated and are summarized in Table 2. For comparison, results are also included for the configuration without the ramp.

The aerodynamic effectiveness of the ramp on the airfoil is limited by several practical considera-tions. First, when the COMPACT actuators are not in use, the lift coefficient of the airfoil is 0.51 lower than the unmodified airfoil, while drag coefficient is increased by 0.084. The control authority (force excursion) of the actuator is greatly enhanced with the addition of the ramp, as the change in lift coeffi-cient in the absence and presence of COMPACT ac-tuation is 0.35, compared to 0.05 of the base airfoil. However, the combined effect of COMPACT and the ramp still produces a net decrease in the lift coeffi-cient of 0.16 and net increase in drag coefficoeffi-cient of 0.045, which is undesirable in practical applications at low angles of attack.

3.3. Alternate Actuator Locations

Although the presence of the passive ramp leads to significant improvement in the control authority of COMPACT, the performance penalties associated

with the presence of the ramp exceed levels that are deemed practically acceptable. Therefore, consider-able optimization must be undertaken to minimize the loss in the aerodynamic performance in the ab-sence of actuation.

The location of the actuator in the simulations and experiments has been fixed at 10% chord from the leading edge on the upper surface of the airfoil. This resembles the configuration used in previous experimental3,4and computational6studies where the objective was to control dynamic stall. However, it may not be the best location for active flow con-trol in the low to moderate angle of attack range. The flexibility of the overset grid method for model-ing the COMPACT actuator in CFD simulations facil-itates testing of alternate actuator locations. There-fore, several new locations for the COMPACT actua-tor are explored and compared in this section.

The new COMPACT configurations include an ac-tuator at 75% chord on the upper surface, 90% chord on the upper surface, and 90% chord on the lower surface of the airfoil. These configurations are located closer to the trailing edge of the airfoil, in contrast to the original configuration where actu-ation was located near the leading edge. The ex-ploration of locations near the trailing edge is mo-tivated by manipulating the unsteady aerodynamic loads by affecting the Kutta condition. For each con-figuration, seven consecutive pulses were applied at

F

+

= 1.67

. As before, the flow conditions corre-spond to a fixed angle of attack of

8

◦ with

M

∞

=

0.06

and

Re

c

= 535, 000

. The phase-averaged re-sponses for the increments in lift, moment, and drag coefficients are plotted in Fig. 10 for the three new configurations as well as the original configu-ration located at 10%-chord from the leading edge. The lift responses of the 10% and 75% configu-rations are nearly identical. However, actuation at 75% chord induces larger increments in the pitching moment and drag. The two 90% configurations af-fect changes in opposite directions for both lift and moment. From these results, there is no definitive improvement over the original actuator location. However, the responses of the two actuators at 90% chord demonstrate that bi-directional control of lift and moment may be achievable by utilizing actua-tion on both the upper and lower surfaces of the airfoil. This could potentially lead to an implemen-tation that resembles the active trailing-edge flap for control of vibrations on a rotor blade.10To sim-ulate positive and negative active flap deflections, the COMPACT actuators could be alternately pulsed on the upper and lower surfaces of the blade.

Preliminary experiments with steady blowing jets at the trailing edge support the bi-directional lift

Actuation Signal (a) Actuation Signal (b) Actuation Signal (c) Actuation Signal (d) Actuation Signal (e) Actuation Signal (f)

Figure 9: CFD and experimental data for multiple COMPACT pulses at

F

+

= 0.634

for COMPACT+Ramp configuration; absolute coefficients (a)-(c) and scaled relative coefficients (d)-(f) for unsteady lift, moment, and drag.

Table 2: Quasi-steady coefficients due to repeated COMPACT actuation at

_F

+

_{= 0.634}

,

_{α = 8}

◦ Configuration Method

∆C

l ,qs

∆C

m,qs

∆C

d ,qs

COMPACT 2D CFD

+0.0495

−0.0017

−0.0006

EXP

+0.0274

−0.0050

+0.0009

COMPACT+Ramp 2D CFD

−0.1632

+0.0092

+0.0452

EXP

_−0.1655

_−0.0017

_+0.0578

and moment control predicted in CFD simulations. The trailing-edge blowing jets were placed at 87.5%-chord locations on both the upper and lower sur-faces of the wind tunnel model. In contrast to the CFD simulations and previous experiments, the jets were directed upstream relative to the oncoming flow (

₆₅

◦relative to chord), as shown in Fig. 11.

Jet flow rates of up to 600 liters/minute (

C

q

=

0.00052

) were tested for angles of attack ranging from

₋₂

◦ to

₁₂

◦. The lift, moment, and drag co-efficients due to the steady blowing are shown in Fig. 12 for both the upper and lower surface jets. The control authority of the jets to affect lift and moment scale linearly with the jet flow rate and is nearly unaffected by the angle of attack. The strength of the jet, and thereby the changes in forces, can be varied continuously for application to

vibration control.

In future experiments, COMPACT actuators will replace the steady blowing jets at the trailing edge of the airfoil. The steady jets currently act as a substitute for the COMPACT actuators because the steady jets could be integrated with the wind tunnel model in a more convenient manner. However, it is difficult to predict the effects that COMPACT actu-ation will have on the aerodynamic loads by con-sidering only the performance of the steady jets. Nevertheless, the experiments and the CFD simu-lations suggest that the bi-directional response of the lift and moment can be anticipated with confi-dence. Additional CFD simulations will also be per-formed in order to investigate the effect of directing the COMPACT actuator jets upstream relative to the oncoming flow.

Actuation Signal

(a) Lift coefficient

Actuation Signal

(b) Moment coefficient

Actuation Signal

(c) Drag coefficient

Figure 10: Phase-averaged coefficients due to COMPACT actuation at 10% chord, 75% chord, and 90% chord on airfoil upper surface and 90% chord on airfoil lower surface

Figure 11: Drawing of steady blowing jets located at 87.5% chord on upper and lower surfaces of VR-12 airfoil model in wind tunnel experiments

4. DEVELOPMENT OF A SURROGATE-BASED REDUCED-ORDER MODEL

The computational expense of the CFD approach for modeling COMPACT actuation is not suitable for direct integration with a comprehensive rotor-craft code. In order to study rotor vibration reduc-tion, an efficient realization of the relationship be-tween the actuator input and the induced cross-sectional force and moment coefficients on the air-foil is needed. Therefore, a surrogate approach is used to generate a reduced-order model (ROM) of the unsteady response due to actuation. Surrogates provide an efficient global approximation by inter-polation of a limited number of samples from ex-pensive, full-order results. A surrogate model can then be evaluated many times in optimization or control studies at a fraction of the cost when com-pared to evaluating the full-order model directly. Thus surrogates can reproduce the results of high-fidelity CFD at a fraction of the cost.

A surrogate-based recurrence framework (SBRF) developed in Ref. 21 is adapted in the current study to model CFD solutions of COMPACT actuation. The SBRF was initially developed to model the unsteady aerodynamic loads induced on an airfoil by

com-bined pitching and plunging motions.21 The SBRF is suitable for predicting unsteady aerodynamic re-sponses, because time history effects associated with unsteadiness are incorporated directly into the model structure. This approach was implemented in Ref. 22 to predict dynamic stall characteristics of an airfoil undergoing realistic rotor-blade motions, including the effects of time-varying Mach number. The previous studies demonstrated that the SBRF approach is suitable for predicting unsteady aero-dynamic loads across a wide range of operating conditions, from attached to separated flows, with significantly less computational expense than full-order CFD simulations.

In the present study, the SBRF is constructed based on the time-dependent responses for the air-foil lift, moment, and drag coefficients calculated us-ing the 2D CFD simulation methodology described in Sec. 2. The design of the SBRF closely follows the approach described in Ref. 21. Therefore, to avoid repeating all the details available in Ref. 21, only the essential modifications made in order to model COMPACT actuation are presented here.

4.1. Surrogate-Based Recurrence Framework The SBRF approximates the response of a non-linear dynamical system in which the output at any time

t

is not only a function of the instantaneous input at

t

, but is also a function of the input time history. The mapping from external inputs,

u(t)

, to a single output of interest,

_{y (t)}

, is expressed in the following form:

(5)

y (t) = Φ (u(t), u(t

− ∆t), . . . , u(t − m∆t),

y (t

− ∆t), . . . , y (t − n∆t)) ,

where

Φ

is an unknown function that, in general, may be computationally expensive to evaluate.

m

(a) Lift coefficient (b) Moment coefficient (c) Drag coefficient

Figure 12: Time-averaged lift, moment, and drag coefficients due to steady trailing-edge blowing jets on upper (

₄

) and lower (

₅

) surface of airfoil model in wind tunnel experiments

and

_n

are integers representing the number of pre-vious inputs and outputs required to account for time-history effects.

The SBRF constructs a surrogate mapping func-tion

_Φ

ˆ

that emulates the response of the full-order system. The SBRF approximation is written as (6)

y (t) = ˆ

ˆ

Φ (ˆ

x(t)) ,

where the state vector

ˆ

x(t)

contains the short-term memory of the external inputs and outputs of the surrogate model:

(7)

ˆ

x(t) = [ u(t) u(t

− ∆t) . . . u(t − m∆t)

ˆ

y (t

− ∆t) . . . ˆ

y (t

− n∆t) ] .

The termrecurrence refers to the property that past outputs from the surrogate model,

y

ˆ

, are fed back as inputs to

Φ

ˆ

to improve the prediction of the re-sponse at the current time

t

.

Note that responses corresponding to

_n

previ-ous time steps are required as inputs to the sur-rogate mapping function. The SBRF is initialized by constructing an additional surrogate model that is based on only the instantaneous value of the input signal:

(8)

_{y (t}

_ˆ

_j

_{) = ˆ}

_Φ

_qs

_u(t

_j

)
,

for

_{j = 1, . . . , n ,}

where

Φ

ˆ

qs denotes that the surrogate can be con-sidered a steady approximation. The quasi-steady surrogate model corresponds to an SBRF model with

m

=

n

=

0

. Following the initialization, the response is predicted using an SBRF model with se-lected values of

m

and

n

, utilizing the previous val-ues generated by the quasi-steady approximation.

The surrogate model for

_Φ

ˆ

is constructed based on kriging interpolation of sample data. In kriging, the unknown function is assumed to be of the form (9)

Φ(x) = f (x) + Z(x) ,

where

_{f (x)}

is an assumed function and

_Z(x)

is a stochastic (random) process that is assumed to be Gaussian. The regression function

f (x)

is a global approximation of

Φ(x)

, while

Z(x)

accounts for lo-cal deviations from

_{f (x)}

that ensure that the kriging model interpolates the sample points exactly. In this study,

f (x)

is assumed to be a second-order poly-nomial and

_Z(x)

is based on a linear spatial corre-lation function provided in Ref. 21. The surrogates were constructed using an available Matlab kriging toolbox,23which determines the optimal fitting pa-rameters for the kriging predictor

_Φ(x)

ˆ

based on maximum likelihood estimation.24

In this study, the outputs of interest from the CFD simulations are the unsteady lift, moment, and drag coefficients of the airfoil. Therefore, each of the out-puts (

C

l,

C

m, or

C

d) corresponds to a mapping func-tion (

Φ

l,

Φ

m, or

Φ

d) for which a surrogate mapping function (

Φ

ˆ

l,

Φ

ˆ

m, or

Φ

ˆ

d) was generated.

4.2. Time-Dependent Sample Data From CFD Simulations

The sample data required for surrogate construc-tion are generated from time-dependent CFD sim-ulations of COMPACT actuation on a static airfoil. This represents the simplest case for the demon-stration of the SBRF, because the only external input is the actuator input signal,

u(t)

. The SBRF originally used in Ref. 21 was for an airfoil with pitching and plunging motions, which required the use of sev-eral external input variables, including the pitch an-gle

_θ(t)

, plunge velocity

_h(t)

˙

, and their time deriva-tives, such that,

(10)

_{u(t) =}

 θ(t) ˙θ(t) ¨

_θ(t)

_{h(t) ¨}

˙

_h(t)

 .

The approach was also extended to represent time-varying Mach number,

_M(t)

, as an external input. However, for the static airfoil considered in the

present study, the vector of external inputs to the SBRF model simplifies to a scalar input, i.e.,

_{u(t) =}

u(t)

. For the current application the actuator input signal is described by a sinusoidal function,

(11)

u(t) =

1

2

(1

− cos (2πf

act

t)) ,

defined by a single parameter,

f

act. This representa-tive input signal

u(t)

is used to correlate the CFD response due to pulses applied at frequency

_f

_act with a reference signal of the same frequency. Al-though actuation is applied in the CFD simulations as discrete, short-duration pulses – i.e., Eq. (2) – it was found that a sinusoidal function such as Eq. (11) results in better performance of the SBRF model. That is, when a pulse input signal with pulse width

T

act,pulse was used to construct the SBRF, the pre-dicted response typically diverged after the first pulse. However, when Eq. (11) is used as a repsentative input signal, the predicted response re-produces the full-order CFD results over the entire simulation period.

As mentioned in Sec. 3.1, the unsteady lift, mo-ment, and drag coefficient data from the CFD sim-ulations contain high-frequency oscillations asso-ciated with each pulse. In order to transform the data from the CFD simulations into a more tractable form for constructing the SBRF models, the CFD data were filtered and sampled with larger time steps. A moving-average filter was applied to the CFD data to smooth the high-frequency compo-nents of the responses. Figure 13 illustrates how fil-tering the CFD data removes the oscillations while preserving the unsteady trends. Next, the filtered data were sampled at time-steps of

∆t = 0.05T

conv to reduce the overall number of data points. This re-sults in two orders of magnitude fewer data points for each training case, since the time step sizes re-quired for accurate CFD results are

O(T

conv

×10

−4

)

.

4.3. Demonstration of the SBRF for COMPACT Actuation on a Static Airfoil

The SBRF approach is used to reproduce full-order CFD results for unsteady changes in lift, moment, and drag coefficients due to COMPACT actuation. The demonstration of the approach is limited to actuation on a static airfoil in this study. The CFD simulations required to construct the SBRF models were performed using the original actuator located at 10% chord on the upper surface. For the case considered, the angle of attack of the airfoil was set at

α = 8

◦, with the same freestream flow con-ditions that have been used throughout this study (

M

∞

= 0.06

,

Re

c

= 535, 000

).

Figure 13: Filtered CFD data for lift coefficient due to COMPACT actuation for reduced-order model con-struction

Several training cases were used to construct the SBRF models, spanning the range of actuation fre-quencies from

f

act

= 10

Hz to

50

Hz. For each training case, actuation was applied from

_{t = 0}

to

6T

conv, a time period sufficient for the unsteady lift coefficient to reach a converged quasi-steady-state value. After the final pulse, a relaxation period of

6T

convwas included in the simulations to allow the lift, moment, and drag coefficients to return to their initial values.

Two different SBRF models were constructed for each the lift, moment, and drag: one based on data from 11 training cases and the other based on 21 training cases. For both models, the training cases were equally spaced across the frequency range of 10-50 Hz. The accuracy of each of the models was compared for four test cases that were not included in either of the training data sets. The test cases were

f

act

= 15, 25, 35,

and

45

Hz, corresponding to

F

+

= 0.75, 1.25, 1.75,

and

2.25

.

The values of

_m

and

_n

for the SBRF models were selected after trying various combinations of the parameters, with values up to 5. The values that produced the best overall accuracy for the test cases were

_{(m, n) = (1, 4)}

for

_C

_l,

_{(2, 4)}

for

_C

_m, and

(2, 4)

for

C

d. Figure 14 shows an example of the unsteady responses predicted by the SBRF models for the lift, moment, and drag due to actuation at

f

act

= 15

Hz (

F

+

= 0.75

). These were obtained us-ing the 21-trainus-ing-case SBRF model.

The errors of the SBRF models relative to the full-order CFD responses were calculated using Eq. (3) for each test case examined. A comparison of the errors for both the 11-training-case and the 21-training-case SBRF models is presented in Ta-ble 3. There is a clear improvement in the error of the lift coefficient when comparing the 11- and 21-training-case SBRF models, as the average unsteady

(a) Lift coefficient (b) Moment coefficient (c) Drag coefficient

Figure 14: Predicted unsteady lift, moment, and drag coefficients due to COMPACT actuation at

F

+

= 0.75

from SBRF models and 2D CFD simulation

error decreases from 12.5% to 6.8%. However, there is only marginal improvement in the errors of the moment and drag coefficients when using the 21-training-case model. Furthermore, the SBRF mod-els reproduced the full-order CFD results within a reasonable margin of error, considering the limited number of training data sets that the models were constructed from. The SBRF models required a frac-tion of a second on a single computer processor to compute the predictions for each test case, while each of the CFD results required 50-60 hours of computer time for calculation on four parallel pro-cessors.

The demonstration of the SBRF for actuation on a static airfoil represents a necessary step towards implementing active flow control for rotor-blade vi-bration reduction. The SBRF was shown in the previ-ous studies to accurately predict unsteady lift, mo-ment, and drag for realistic rotor blade motions at reduced computational cost. Therefore, the model-ing approach is clearly appropriate for application to helicopter rotor aeroelastic computations. How-ever, further development of the SBRF is needed in order to model the combined effects of COMPACT actuation and unsteady rotor blade motions. The ultimate goal is to implement the reduced-order model into a comprehensive code for helicopter ro-tor vibration reduction.

5. CONCLUSIONS

The principal conclusions of this study are:

1) The unsteady aerodynamic effects of COM-PACT actuation on the lift, moment, and drag co-efficients of a 2D airfoil are accurately represented using the computational methods described in this paper. The 2D CFD simulations were correlated with wind tunnel experiments for pulse actuation on a static airfoil, and the agreement with the

experi-mental data was acceptable.

2) The 2D CFD method presented is suitable for studying actuator modifications needed to enhance control authority at low angles of attack. The modi-fications consisted of ramp placement on the airfoil as well as actuator relocation.

3) No definite improvement was found from the modifications. However, actuation at the airfoil trail-ing edge indicated the potential for bi-directional control of lift and moment when actuators were used on both the upper and lower surfaces.

4) The surrogate-based reduced order model re-produces full-order CFD results for COMPACT actu-ation with a reasonable degree of accuracy and at a fraction of the computational cost.

ACKNOWLEDGEMENT

The authors thank Professor M. Smith (Georgia Tech) for her advice in formulating the mesh struc-ture used for modeling the COMPACT actuators in the CFD simulations. This study is supported under the NRTC Vertical Lift Rotorcraft Center of Excellence (VLRCOE) from the U.S. Army Aviation and Missile Research, Development and Engineer-ing Center (AMRDEC) under Technology Investment Agreement W911W6-17-2-0002, entitled Georgia Tech Vertical Lift Research Center of Excellence (GT-VLRCOE) with Dr. Mahendra Bhagwat as the Pro-gram Manager.

REFERENCES

[1] Gardner, A. D., Richter, K., Mai, H., and Neuhaus, D., “Experimental investigation of high-pressure pulsed blowing for dynamic stall control,”CEAS Aeronautical Journal, Vol. 5, No. 2, Feb. 2014, pp. 185–198.

Table 3: Unsteady error of SBRF models constructed from 11 and 21 training data sets for 4 test cases SBRF model

_{ ˆ}

_C

_l

_{, C}

_l

_{ ˆ}

_C

_m

_{, C}

_m

_{ ˆ}

_C

_d

_{, C}

_d

11 training cases

4.0

− 17.1%

(avg

12.5%

)

3.8

− 5.1%

(avg

4.4%

)

4.1

− 6.3%

(avg

5.0%

) 21 training cases

_4.0

_{− 9.6%}

(avg

_6.8%

)

_2.8

_{− 5.8%}

(avg

_4.5%

)

_3.7

_{− 5.6%}

(avg

_4.7%

)

[2] Bons, J., Frankhouser, M., and Gregory, J., “Syn-chronized Flow Control of Dynamic Stall under Coupled Pitch and Freestream Oscillations,” AHS International 73rd Annual Forum, American Helicopter Society International, Fort Worth, TX, May 2017.

[3] Tan, Y., Crittenden, T. M., and Glezer, A., “Aero-dynamic Control of a Dynamically Pitching VR-12 Airfoil using Discrete Pulsed Actuation,” 54th AIAA Aerospace Sciences Meeting, American Institute of Aeronautics and Astronautics, San Diego, CA, Jan. 2016, AIAA Paper 2016-0321. [4] Matalanis, C., Bowles, P., Lorber, P., Crittenden,

T., Glezer, A., Schaeffler, N., Min, B.-Y., Jee, S., Kuczek, A., and Wake, B., “High-Speed Experi-ments on Combustion-Powered Actuation for Dynamic Stall Suppression,” _{AHS 72nd Annual} Forum, American Helicopter Society Interna-tional, West Palm Beach, FL, May 2016.

[5] Duraisamy, K. and Baeder, J., “Active Flow Con-trol Concepts for Rotor Airfoils Using Synthetic Jets,”1st Flow Control Conference, American In-stitute of Aeronautics and Astronautics, St. Louis, MO, June 2002, AIAA Paper 2002-2835. [6] Jee, S., Bowles, P., Matalanis, C., Min, B.-Y.,

Wake, B., Crittenden, T., and Glezer, A., “Com-putations of Combustion-Powered Actuation for Dynamic Stall Suppression,”AHS 72nd An-nual Forum, American Helicopter Society Inter-national, West Palm Beach, FL, May 2016. [7] Duvigneau, R., Hay, A., and Visonneau, M.,

“Study on the Optimal Location of a Synthetic Jet for Stall Control,”_{3rd AIAA Flow Control} Con-ference, American Institute of Aeronautics and Astronautics, San Francisco, CA, June 2006, AIAA Paper 2006-3679.

[8] Parekh, D., Palaniswamy, S., and Goldberg, U., “Numerical Simulation of Separation Control via Synthetic Jets,”1st Flow Control Conference, American Institute of Aeronautics and Astro-nautics, St. Louis, MO, June 2002, AIAA Paper 2002-3167.

[9] You, D. and Moin, P., “Active Control of Flow Separation over an Airfoil using Synthetic Jets,” Journal of Fluids and Structures, Vol. 24, No. 8, Nov. 2008, pp. 1349–1357.

[10] Friedmann, P. P., “On-Blade Control of Rotor Vi-bration, Noise, and Performance: Just Around

the Corner? The 33rd Alexander Nikolsky Hon-orary Lecture,” _{Journal of the American} Heli-copter Society, Vol. 59, No. 4, Oct. 2014, pp. 1– 37.

[11] Glaz, B., Friedmann, P., Liu, L., Kumar, D., and Cesnik, C., “The AVINOR Aeroelastic Simula-tion Code and Its ApplicaSimula-tion to Reduced Vi-bration Composite Rotor Blade Design,”50th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference, Ameri-can Institute of Aeronautics and Astronautics, Palm Springs, CA, May 2009, AIAA Paper 2009-2601.

[12] Padthe, A. K. and Friedmann, P. P., “Simul-taneous Blade-Vortex Interaction Noise and Vibration Reduction in Rotorcraft Using Mi-croflaps, Including the Effect of Actuator Sat-uration,”Journal of the American Helicopter So-ciety, Vol. 60, No. 4, Oct. 2015.

[13] Myrtle, T. F. and Friedmann, P. P., “Applica-tion of a New Compressible Time Domain Aerodynamic Model to Vibration Reduction in Helicopters Using an Actively Controlled Flap,”_{Journal of the American Helicopter Society,} Vol. 46, No. 1, Jan. 2001, pp. 32–43.

[14] Liu, L., Padthe, A. K., Friedmann, P. P., Quon, E., and Smith, M. J., “Unsteady Aerodynam-ics of an Airfoil/Flap Combination on a Heli-copter Rotor Using Computational Fluid Dy-namics and Approximate Methods,”_{Journal of} the American Helicopter Society, Vol. 56, No. 3, July 2011.

[15] Crittenden, T., Woo, G. T., and Glezer, A., “Combustion Powered Actuators for Separa-tion Control,”6th AIAA Flow Control Conference, American Institute of Aeronautics and Astro-nautics, New Orleans, LA, June 2012, AIAA Pa-per 2012-3135.

[16] Chakravarthy, S., Peroomian, O., Goldberg, U., and Palaniswamy, S., “The CFD++ Computa-tional Fluid Dynamics Software Suite,”_{AIAA and} SAE, 1998 World Aviation Conference, American Institute of Aeronautics and Astronautics, Ana-heim, CA, Sept. 1998, AIAA Paper 1998-5564. [17] Metacomp Technologies, Inc., Agoura Hills, CA,

USA,CFD++: Version 17.1, User Manual, 2017. [18] Crittenden, T., Warta, B., and Glezer, A.,

for Flow Control,”_{3rd AIAA Flow Control} Confer-ence, American Institute of Aeronautics and As-tronautics, San Francisco, CA, June 2006, AIAA Paper 2006-2864.

[19] An, X., Williams, D. R., Eldredge, J., and Colo-nius, T., “Modeling Dynamic Lift Response to Actuation,” 54th AIAA Aerospace Sciences Meet-ing, American Institute of Aeronautics and As-tronautics, San Diego, CA, Jan. 2016, AIAA Pa-per 2016-0058.

[20] Brzozowski, D. P., Woo, G. T. K., Culp, J. R., and Glezer, A., “Transient Separation Control Us-ing Pulse-Combustion Actuation,”AIAA Journal, Vol. 48, No. 11, Nov. 2010, pp. 2482–2490. [21] Glaz, B., Liu, L., and Friedmann, P. P.,

“Reduced-Order Nonlinear Unsteady Aerodynamic Mod-eling Using a Surrogate-Based Recurrence Framework,”_{AIAA Journal, Vol. 48, No. 10, Oct.} 2010, pp. 2418–2429.

[22] Glaz, B., Liu, L., Friedmann, P. P., Bain, J., and Sankar, L. N., “A Surrogate-Based Approach to Reduced-Order Dynamic Stall Modeling,” Jour-nal of the American Helicopter Society, Vol. 57, No. 2, April 2012.

[23] Lophaven, S. N., Nielsen, H. B., and Sønder-gaard, J., “DACE - A Matlab Kriging Toolbox, Ver-sion 2.0,” Tech. rep., Informatics and Mathe-matical Modelling, Technical University of Den-mark, DTU, Lyngby, DenDen-mark, 2002, IMM-TR-2002-12.

[24] Sacks, J., Welch, W., Mitchell, T., and Wynn, H., “Design and Analysis of Computer Experi-ments,”_{Statistical Science, Vol. 4, No. 4, 1989,} pp. 409–435.