{'fraction_non_alphanumeric': 0.030082466833990677, 'fraction_numerical': 0.010935819290068125, 'mean_word_length': 4.365717583686033, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Lift Coefficient Estimation for a Rapidly Pitching Airfoil

A method for the lift coefficient estimation over a rapidly pitching NACA0009 wing is proposed that contains three components. First, we establish that the Goman-Khrabrov model is in fact, a linear parameter-varying (LPV) system, therefore it is suitable for a Kalman filter without any linearization. In the second part we attempt to estimate the lift coefficient by measuring the surface pressure from four pressure sensors located on the suction side of the wing. We demonstrate that four pressure sensors alone, are not sufficient to capture the lift coefficient variation during the rapidly pitching maneuvers, and this results in non-Gaussian error. In the last part we demonstrate the non-Gaussian error from the pressure estimated lift coefficient introduces additional errors into the estimator when we employ the conventional Kalman filter design. To address this issue, we propose a new method of coupling the model and the measurement through the Kalman filter. It is shown that the proposed Kalman filter is capable of estimating the lift coefficient accurately on a NACA 0009 wing that is undergoing rapidly pitching maneuvers.


Introduction
Rapidly pitching maneuvers are commonly seen on the helicopter rotational blades, airplanes and their control surfaces (e.g. elevators). As an example, during forward flight, in order to counterbalance the relative incoming flow speed difference over the blades between the port and starboard sides, helicopters have to continuously adjust the incline angle of their rotor blades within one rotation cycle. Failing to predict and control the lift force variation on the blades during this type of rapidly pitching maneuver will result in the asymmetry of lift between the port and starboard sides, which causes an undesired rolling moment on the helicopter. In fact, the asymmetry of lift plays an important role in limiting the helicopter forward flying speed. On the other hand, in order to adjust the flight attitude during landing approaches or maintain the flight attitude when flying through wind gusts, airplanes have to move their control surfaces in the manner of fast pitching. Under these circumstances, the control surfaces could undergo rapidly pitching maneuvers even the airplanes themselves do not. Meanwhile, high-performance airplanes are capable of changing their pitch angle in a short time, which leads to the entire airplane undergoing rapidly pitching maneuvers. To achieve better flight performance for both helicopters and airplanes, model estimators that are capable of accurately predicting the lift force in real time are desired for the flight control applications.
In contrast to the quasi-steady pitching maneuvers, rapidly pitching maneuvers could trigger complicated flow responses, including the dynamic stall vortex formation and vortex shedding [3], which make the lift variation difficult to estimate in real time. High fidelity methods, such as Computational Fluid Dynamics (CFD), are too costly in terms of computation time for realtime control systems, so that lower-order models are desired for real-time controllers. The classical Wagner model [14], which is a low-order approach employing the potential flow assumption, accounts for the bound vorticity, the vorticity in the wake as well as the Kutta condition. Theodorsen [13] extended Wagner's approach to a multi-plate configuration. The drawback of these methods is that the flow separation is not modeled, so that they lose accuracy when the boundary layer separates from the plate's surface.
More recently, using experimental data a linear frequency response model developed by Kerstens, et al. [9] was used to capture the unsteady lift variation due to periodic and quasi-random time-varying freestream flow acting on an airfoil at high angles of attack. The wing's angle of attack was always higher than the static stall angle, so the flow on the suction side remains separated all the time. But this black-box model was not capable of capturing the transitional behavior between the attached and separated flow states.
The ability of a low-dimensional model to capture the transition between the attached and separated flow is essential for a useful unsteady pitching wing aerodynamic models. Hemati, et al. [7], Brunton, et al. [2], Dawson [4] and Provost, et al. [10] introduced similar linear parameter-varying models, which showed good performance for aerodynamic loads tracking for rapidly pitching wings. However, these models are purely data-driven and the models' physical insight remains to be investigated. Moreover, unstable eigenvalues might be identified from the training data set, due to some phenomena within the system's dynamic that are not modeled (eg. turbulence effect). Meanwhile, Goman & Khrabrov [5] proposed a model, which we refer to as the G-K model, utilizing the static measurement as a forcing term. This model is capable of predicting the lift force during rapid maneuvers and the model's stability is guaranteed by constraining the time constant values. Grimaud [6] and Williams, et al. [16] modified the G-K model in such a way that the lift as a function of an inner state variable and angle of attack. The lift function is directly generated from the static lift curve. This is a critical modification on the G-K model that allows one to systematically generate the lift function rather than the trial and error approach in the original G-K model. Readers interested in the evolution of the G-K type models please refer to a detailed review article given by Williams and King [17] Despite the good performance of the G-K model, this low-order model neglects some other important features in the flowfield, such as, trailing-edge vortex separation, natural vortex shedding and some other unmodeled disturbances [12]. This motivates us to employ real-time measurement to account for the unmodeled features. For example, surface pressure has been widely used for lift estimation. Comparing to the force measurement, surface pressure measurement provides more detailed information about the surrounding flow field and it is more feasible to measure during actual flight than lift. Some investigators [4] [10] have employed the Kalman filter [8] to assimilate pressure measurements into low-order models (other than the G-K model) for estimation of the real-time aerodynamic loads variation in response to different types of wing/aircraft maneuvers.
In reality, it is not feasible to measure the pressure on the entire surface of the wing/aircraft. In fact, sparse distributions of the pressure sensors are required for most of the applications. An [1] showed that with a limited number of pressure sensors, it is possible to project the state variable (pressure distribution along the entire airfoil) onto its sub-space (sparse pressure measurements), which then leads to colored (non-Gaussian) noise for the lift that is estimated by the sparse pressure measurements. Such colored noise can be even nonlinear. We refer to this colored noise as the measurement error for the remainder of this paper, in contrast to the white (Gaussian) measurement noise. The measurement error can be problematic when the Kalman filter is implemented. Unlike white noise, with little knowledge of the measurement error (nonlinear colored noise), it is difficult for a Kalman filter to reduce the measurement error. Therefore, a new way of coupling the model and the measurement is proposed in the present work that enables the Kalman filter to reduce the measurement error.
In this paper, we will show that the G-K model is in fact an LPV model with a nonlinear forcing term, so that the model is suitable for use in a conventional Kalman filter without any linearization of the G-K model. This is the original Kalman filter, which can be only used for linear systems Then we will exhibit the measurement error (colored noise) within the lift estimated by the sparsely distributed pressure sensors. Finally, a state estimator design based on the Kalman filter will be used to couple the LPV model and the pressure measurements to reduce both the modeling and the measurement errors. The rest of this paper is organized as follows. The experimental setup is described in Sec 2. The derivation of the LPV form is given in Sec 3. The measurement error of the lift coefficient due to the sparsely distributed pressure sensors is shown in Sec 4. The design of the Kalman filter is discussed in Sec 5, and the conclusions are given in Sec 7.

Experimental Setup
The experiments were conducted in the Andrew Fejer Unsteady Flow Wind Tunnel at Illinois Institute of Technology. The test section of the wind tunnel has cross-section dimensions 600mm × 600mm. A nominally two-dimensional NACA0009 wing with a wingspan b = 596mm and chord length c = 245mm was used as the test article (Fig. 1a). The freestream speed was U ∞ = 3m/s, corresponding to a convective time t convect = c U∞ ≈ 0.08s, t + = t tconvect , and chord-based Reynolds number 49,000. The reduced frequency k is defined as k = πf c U∞ where f is the frequency in Hz. The wing was mounted to an ATI Nano-17 force/moment transducer that was connected to the pitch-plunge mechanism consisting of two computer-controlled Copley servo tubes. The two servo tubes enable the pivot point for the pitching motion to be changed. For the results presented in this paper, the pivot point was at the location, x/c = 0.15. Pitch rates were restricted to 2 Hz or less (k ¡ 0.55) to avoid overstressing the force balance (FB). Forces were measured with the force balance located inside the model at 30% of the chord, which is the center of gravity of the wing. Four pressure sensors (All Sensors 1inchD2P4Vmini) are located on the upper surface of the wing along its chord line. The locations of the force balance and the pressure sensors are shown in Fig. 1b. 3 Lift response to the pitching maneuver

Nonlinear G-K model
The Goman-Khrabrov [5] (G-K) state-space model was reported to successfully predict the lift variation over a wing or the entire aircraft during arbitrary pitching maneuvers [11]. The G-K model (Eq. 1 and Eq. 2) uses two time constants (τ 1 and τ 1 ), and an internal dynamic variable x that nominally represents the degree of flow attachment over a wing. Fully attached flow corresponds to x = 1, and fully separated flow is x = 0. The quasi-steady position of the separation point is given by the function x 0 (α), which is shown in Fig. 2. Since the computational cost of G-K model is very low, the two time constants are often obtained by running through all the possible values  of the time constants to find the values that minimize the mean square error between the model and the training data sets of a dynamic pitching motion. For the current test conditions, the time lag associated with dynamic stall vortex formation and its convection over the wing is represented by τ 2 = 4.375t + , which is approximately four convective times. The relaxation time constant is τ 1 = 3.75t + . The time constants and function X 0 (α) are determined once from the training data for a specific airfoil. The Euler method of integration is then used to compute X(t) from Eq. 1 during real-time experiments. The instantaneous lift coefficient is found from Eq. 2. Eq. 2 was initially identified by Grimaud [6], and then extended to a more generalized form by Williams, et al. [18], which is shown in Eq. 3.
where C 1 is the α − C L slope, dC L dα , when the flow is fully attached, C 3 is the zero-lift angle, C 2 is the α − C L slope, dC L dα for fully separated flow and C 4 is the C L value at the smallest α when the flow is fully separated. We will refer Eq. 1 and Eq. 3 as the modified G-K model, or mG-K model for abbreviation. Eq. 3 enables a systematic method for C L (α, x) function generation. In the original G-K model, however, it is a trial and error approach.

The relation between the linear parameter-varying model and the mG-K model
By substituting Eq. 3 into Eq. 1 and discretizing it, one obtains Eq. 4, can be reorganized into the form where, An LPV system is a linearly evolving system with time varying coefficients (e.g. gain) that depends on some measurable parameters. Therefore, Eq. 5 is a linear parameter-varying (LPV) dynamic system, because the function A only depends on the input variables α(k) and α(k − 1) but not any output variables C L . The last two terms on the right-hand side are nonlinear input functions (forcing terms), which only act on the system input variables.
This model can be applied directly to the conventional Kalman filter as part of a feedback controller. A more detailed proof of the application of Kalman filter on the LPV system is given in Appendix. Since the LPV model and the mG-K models are exactly the same in the rest of this paper, we will not differentiate the LPV model and the mG-K model. Next, the mG-K model will be validated by fast (k ≥ 0.05) periodic and random pitching maneuvers.

Periodic motion
When the wing is pitching over some range in α, the lift coefficient deviates from the quasi-steady values and hysteresis loops are formed. The ability of the mG-K model to predict the lift hysteresis that occurs during periodic pitching motions is shown for four different pitching cases, k = 0.05 from 2 o to 8 o , k = 0.1 from 2.  The mG-K model is capable of tracking the changes in C L during periodic pitching motions. Even the dynamic stall in Fig. 3a and Fig. 3b is captured where the flow is attached in the quasi-steady case shown in Fig. 2.

Quasi-random motion
The ability of the mG-K model to predict the lift produced by a quasi-random pitching motion is shown in Fig. 4. The quasi-random pitching motion was constructed by superposing 10 sinusoidal signals with random initial phases relative to each other. The highest frequency sinusoidal signal is k = 0.51. The mG-K model prediction for the random pitching is shown by the solid red lines, and is compared to the experimental measurement shown as the solid blue lines. The mG-K model prediction closely tracks the experimental data, and the correlation coefficient between them is 0.956. However, some errors still exist and some small amplitude fluctuations are not modeled by the mG-K model. To improve the model for the C L estimation we will include sparsely distributed real-time pressure measurements in a Kalman filter.  Prior to building a Kalman filter that incorporates pressure measurements, a relation between the pressure measurements and the lift coefficient has to be established. To estimate the lift coefficients with the pressure measurement data, a set of weighting coefficients along with an offset are employed to formulate where C L (k) is the C L at time instant k, α(k) is the angle of attack at time instant k, p i (k) is the pressure reading on the ith pressure sensor at time instant k, w i is a weighting coefficient for the ith pressure sensor and N is the number of pressure sensors. Note that p N +1 (k) = 1 is a constant offset. Since we only have pressure sensors on the suction side of the wing, p N +1 (k) = 1 takes pressure force on the pressure side of the wing into account. To simplify Eq. 6, we let P i (k) = cos (α(k))p i (k). Then, a more compact form of C Lpressure relation can be obtained by rewriting Eq. 6 in matrix form, It is very important to point out that the pressure P has already taken the α(k) into account, thus the weight w is independent of α(k). The weights can be solved by the pseudoinverse of Eq. 7, Therefore, w can be solved offline once and for all with some experimental training data sets that include both pressure and C L measurements. This means that C L can be estimated using the measured pressure by Eq. 6 after C L is identified.
In the current research, there are four pressure sensors on the surface of the airfoil's suction side (Fig. 1b). The reason for placing all four pressure sensors on the suction side of the airfoil is that the present research focuses on positive angles of attack, and most of the complex fluid dynamic events (e.g. flow separation, leading-edge vortex shedding) happen on the suction side. A training data set consisting of 4076 data points from an airfoil that pitches from 13 o to 19 o at k = 0.13 was used to solve for w. Fig. 5 shows the lift coefficient estimated by the pressure measurements when the airfoil is pitching between 13 o < α < 19 o at a frequency k = 0.13. The C L estimation using pressure fits the data measured by the force balance very well for the training case itself. The correlation coefficient between the pressure estimated C L and directly measured C L is 0.906. The α effect is included in P as mentioned in Eq. 7, so the weight w function should not be a strong function of α). The strong correlation between the measured lift coefficient and the pressurebased estimate shown in Fig. 5 suggests that this simple approach would be viable on its own as a model for lift. However, there are two major drawbacks when using pressure to estimate C L directly. First, Fig. 6 indicates that the pressure C L estimation could produce more measurement error (colored noise) for the non-training cases. It can be seen that there exists a bias for the steady state between 0t + and 200t + , which produces a time-varying error after the pitching motion is started. The other drawback of using the pressure C L prediction directly is that it would also introduce white noise from the real-time measurement. To alleviate these additional sources of noise, the idea of coupling the pressure measurement and the mG-K model utilizing a Kalman filter will be given in the next section. Fig. 6: Comparison of C L measured by the force balance and C L reconstructed by pressure sensors for the non-training case.

Linear Quadratic Estimator design
The common approach to Kalman filter design would be to use the C L directly predicted by mG-K model in combination with the C L estimated by the pressure measurements as shown in the schematic of the Kalman filter in Fig. 7. However, as it is shown in Fig. 8, this Kalman filter architecture is not capable of accurately estimating the C L due to the measurement error that was mentioned in Sec. 4.  Recall that the mG-K model is capable of modeling the C L with no bias error due to its utilization of the static lift measurement. This model only has angle of attack as an input, and has no information about the unmodeled disturbances. On the other hand, the pressure-predicted C L has a measurement error as shown in Fig. 6. The pressure-based model can provide real-time information about the flowfield that is not modeled by the mG-K model. This motivates us to modify the Kalman filter design in a way that couples the mG-K model and pressure measurements to take the advantage of both the model prediction and the update with a measurement.
To achieve this goal, we start by representing the C L by all pressure measurements. Reorganizing Eq. 6 into a more compact form allows C L to be expressed as where i is the index of the pressure sensor number, k denotes the time instant. In this case, M = 5, since we have four pressure sensors and one offset constant. P i is the corrected pressure which is the production of the pressure and cosine of the time-varying α, so that the weighting coefficients are independent of α. If C L is replaced with C L,mG−K that is computed by the mG-K model, then each pressure sensor value can be expressed by C L,mG−K and the sum of the weighted pressures on the other pressure sensors, as shown in Eq. 10.
Writing the entire system combining the pressure and mG-K model in the matrix form, we have Note that A is the time-varying coefficient in Eq. 5, w 5 is used to track the error for the pressure estimated C L . The corresponding state vectors arê P of f is a state variable that is used to track the error of pressure predicted C L . Here ". " denotes the a posterior state estimation given the measurement. The time-evolving equation for this system than becomes Coefficients B and C are computed by Eq. 5. They can be treated as constant at each time step, and they are not part of matrix A K . Thus, the standard form of Kalman filter update algorithm can be applied. The matrix A K needs an online update for each time step. The measurement matrix is expressed as Eq. 15, since we only have access to data from four pressure sensors.
In this work, defining ω as the processing noise and v as the measurement noise, the covariance matrices Q = E(ωω T ) and R = E(vv T ) are chosen to be diagonal matrices with appropriate equal diagonal elements. The values of the diagonal elements of R matrix can be 10 3 to 10 4 times larger than Q to reduce the measurement noise. The conventional Kalman filter prediction step and measurement update are used here, which is the common approach used in numerous control applications. Readers interested in the details of the Kalman filter should refer to [8] The architecture of the Kalman filter is shown in Fig. 9, the input to the Kalman filter are mG-K predicted C L , 4 pressure measurements as well as the pressure computed from Eq. 10. The major advantage of this new Kalman filter algorithm (Fig. 9) compared to the original Kalman filter design (Fig. 7) is that the pressure is directly coupled with the mG-K model thoughĈ L . Since the mG-K model is forced by the static C L measurement, it is unlikely to have colored noise within the model. The colored noise of the pressure estimated C L can be suppressed by mG-K model, even if the colored noise is nonlinear. Fig. 9: Block digram of the Kalman filter architecture.

Validation of the Kalman filter
The proposed Kalman filter was validated using a randomly pitching NACA0009 airfoil. Two random pitching signals that were different from the initial training data were used. The first random pitching signal is the same as the one used in Sec. 3. The second random pitching signal is different from the first one and will be discussed in detail later.
The validation of the Kalman filter against the first random pitching motion is exhibited in Fig. 10. It shows that the C L estimation from the Kalman filter is tracking the experimental force balance data "true" C L very well. Both the unmodeled fluctuation and the main trend of C L are captured, and the C L noise (both colored and white) level is also reduced by the Kalman filter. The correlation coefficient between the experimental data and the Kalman filter output is 0.9637, which is higher than either the mG-K model (0.959) or pressure C L estimation (0.9277) alone.
One could argue that the good performance of the Kalman filter is only due to the accurate mG-K model. To further exhibit the ability of this Kalman filter, two types of artificial errors were added into the mG-K model. The first modeling error (case 1) is simulated through an error onα, in which theα input to the system was reduced by 80% from the actual value, and then multiplied by an extra error term of sin(t/0.05). The second modeling error (case 2) is to use incorrect time constants in the mG-K model to simulate the time response error. In this case τ 1 is increased by 55% from the original value and τ 2 is reduced to 20% of its original value. The results are shown in Fig. 11 and Fig.  12. In both cases, the Kalman filter is tracking the force balance measured "true" C L signal well. The correlation coefficients between the Kalman filter output and the experimental force balance data are 0.9489 for the former case and 0.9390 for the latter case. Therefore, it has been demonstrated that this Kalman filter design achieved good error-tolerant features for both the processing model and the pressure measurement. It is worth of pointing out that all the three figures ( Fig. 10 to Fig. 12) were used to show the validation of the Kalman filter against the first random pitching motion.   In order to further validate the Kalman filter, we also applied the Kalman filter on the second random motion. The major difference between this random pitching motion and the previous one this that the smaller pitching amplitude results in worse performance of the mG-K model. To be more specific, a series of vortex motions, that are not modeled by the mG-K model, play a more important role in the C L variation for smaller pitching amplitude. Fig. 13 exhibits the pitching motion of the second random pitching signal. The demonstration of the validation of the Kalman filter against the second random pitching motion is shown in Fig. 14. It is obvious that the Kalman filter outperform both the mG-K model and the pressure measurements C L prediction. To quantify the performance of the Kalman filter against the second random pitching, the correlation coefficient between the experimental force balance data "true" C L and the Kalman filter is 0.7070, which is higher than either the mG-K model (0.6362) or pressure C L estimation (0.6458) alone.

Conclusion
In the current work, we have shown that the mG-K model is a linear parametervarying (LPV) model with nonlinear input forcing terms. Lift coefficient estimation using pressure measurements was performed by giving a weight to each pressure sensor. The measurement error (colored noise) of the lift coefficient estimation from the sparsely distributed pressure sensors was also exhibited. Then, a Kalman filter framework that couples the mG-K model and pressure measurement was introduced for a rapidly pitching airfoil to take advantage of both the mG-K model and the pressure measurements. The Kalman filter was validated using a randomly pitching NACA0009 airfoil. Good performance of C L estimation was observed, and both the modeling error and the measurement error (colored noise) and noise (white noise) are reduced. It also has been shown that even with the significant artificial and real modeling error, this Kalman filtering approach still provides accurate estimates of the timevarying C L .

Appendix
The proof of Kalman filter is applicable on LPV systems with nonlinear input Starting from a LPV dynamic system with nonlinear input where X k is the state X ∈ R n at time instant k, f (u) is the nonlinear input function and ω is the white Gaussian processing noise. All the subscripts in the remaining of this section denote the time instants. The measurement Z ∈ R m at time instant k + 1 is Where H is the measurement matrix and v is the measurement noise. Following a similar algorithm proposed by Kalman, et. al. [8] and Welch and Bishop [15], the time update of the discrete Kalman filter can be then expressed aŝ here,X − k+1 is the a priori estimate at (k + 1)th time step, P − k+1 is the a priori estimate error covariance, e − k+1 is the a priori estimate error, E[.] denotes the expectation, andX k is the a posteriori state estimate at time step k as a linear combination of the a priori estimateX − k and a weighted difference between the actual measurement Z k and a measurement prediction HX − k . Hence, following a similar procedure proposed by Kalman,et. al. [8], the discrete Kalman filter measurement update equations are where K k+1 is the Kalman gain at time instant k + 1, R is the measurement noise covariance and P k+1 is the a posteriori estimate error covariance.
It is important to point out that at this stage, P − k+1 in Eq. 19 is the only term that contains the nonlinear input function f (u k ), since the a priori estimate error e − k+1 can be expressed as However, it is obvious that the f (u k ) terms are canceled out. Thus, the a priori estimate error covariance where Q is process noise covariance. Therefore, by replacing Eq. 19 with Eq. 25, it can be seen that the Kalman filter algorithm for the LPV dynamic system is the same as the Original Kalman filter despite the time-varying A k .