Abstract
A reduced order modeling method based on algorithm fusion and multifidelity framework for nonlinear unsteady aerodynamics is proposed to obtain a low-cost and high-precision unsteady aerodynamic model. This method integrates the traditional algorithm, intelligent algorithm, and multifidelity data fusion algorithm. In this method, the traditional algorithm is based on separated flow theory, the intelligent algorithm refers to the nonlinear autoregressive (NARX) method, and the multifidelity data fusion algorithm uses different fidelity data for aerodynamic modeling, which can shorten the time cost of data acquisition. In the process of modeling, firstly, a multifidelity model with NARX description provides a general intelligent algorithm framework for unsteady aerodynamics. Then, based on the separated flow theory, the correction equation from low-fidelity model to high-fidelity result is constructed, and the cuckoo algorithm based on chaos optimization is used to identify the parameters. In order to verify the effectiveness of the method, an unsteady aerodynamic model of NACA0012 airfoil is established. Three kinds of data with low, medium, and high fidelity are used for modeling. The low-fidelity and medium-fidelity data is obtained from the CFD-Euler solver and CFD-RANS solver, respectively, while the high-fidelity data comes from the experimental results. Then, the model is established, and its prediction of unsteady aerodynamic coefficients is in good agreement with the CFD results and the experimental data. After that, the model is applied to a two-dimensional aeroelastic system, and the bifurcation and limit cycle response analysis are compared with the experimental results, which further shows that the model can accurately capture the main flow characteristics in the flow range of low speed and high angle of attack. In addition, the convergence of the model is studied; the accuracy and generalization ability as well as applicability scope of the model are compared with other aerodynamic models and finally discussed.
1. Introduction
Most of the problems in aircraft aerodynamic design are closely related to unsteady aerodynamic forces. Aerodynamic modeling is a major foundation for flight overall parameter design, flight trajectory design, and flight maneuverability and stability analysis. The accuracy of the aerodynamic model directly affects the performance of the control system and the reliability of flight simulation [1]. There are many factors that can affect the aerodynamic force such as flight altitude, airspeed, movement form, and rudder interference, which leads to a strong nonlinear aerodynamic system. However, the traditional aerodynamic modeling method is based on the assumption of quasisteady and linearization theory, which cannot meet the requirements for advanced civil and military aircraft design. So, in recent years, unsteady and nonlinear aerodynamic modeling methods have been paid more and more attention. There are two major problems in the unsteady aerodynamic modeling, one is the calculation of dynamic derivative, and the other is the aerodynamic modeling at a high angle of attack. As the former directly determines flight quality analysis and flight control system design, the latter is one of the important ways to evaluate the performance of aircraft at a high angle of attack. Therefore, it is of great value to study the flight dynamic derivative identification methods as well as nonlinear and unsteady aerodynamic modeling at high angles of attack.
The unsteady aerodynamic modeling method appeared at the beginning of the 20th century, and a variety of modeling methods had been widely developed since the 1980s. These methods can be divided into two categories: traditional mathematical expression method and intelligent learning methods. The traditional mathematical expression methods are a kind of method which obtain mathematical relationships between aerodynamic force, flight states, control inputs, and other parameters according to the physical characteristics and statistical laws of aerodynamic force. The traditional aerodynamic model is based on a dynamic stability derivative at the beginning, which is used to predict the change of flight load after the flight parameters change. These kinds of models regard each dynamic stability derivative as a constant value independent of time and other variables, so these models are linear aerodynamic models and are only applicable to the linear flight states with a small angle of attack [2]. In order to expand the application range of these models, the nonlinear algebraic model [3], Fourier series model [4–6], integral equation model [7], state space model [8, 9], etc., were developed. However, generally speaking, the traditional mathematical expression method of aerodynamic modeling is more suitable for aerodynamic modeling with limited flight state range changes, such as the linear/steady aerodynamic modeling under certain altitude, speed, angle of attack, and other conditions. The relationship between aerodynamic force and flight height, velocity, angle of attack, sideslip angle, angular velocity, rudder angle, and other parameters is highly nonlinear for the aircraft with large changes in airspace, airspeed, flight attitude, and control surface deflection range. It is difficult for traditional mathematical modeling methods to give an accurate description of all the relationships with acceptable errors.
With the rapid development of computer technology and the integration of interdisciplinary ideas, many other mathematical methods have also been introduced into the field of unsteady aerodynamic modeling, such as the autoregressive model [10–12], Volterra series model [13–15], fuzzy logic method [16], neural network method [17–21], and support vector machine method [22–24], collectively referred to as the intelligent algorithm. The unsteady aerodynamic model based on these algorithms does not pay attention to the physical principles; that is to say, the relationship between the aerodynamic force and the input state quantity is regarded as a “black box.” The aerodynamic model is established by a large number of high-quality sample data, so this kind of model is also called the data-driven model.
Although these methods can give accurate prediction within the range of sample flight states, the ability of extrapolation is usually very poor. At present, the research of aerodynamic modeling based on the intelligent algorithm is mostly based on the results of the wind tunnel dynamic test or CFD calculation. However, the model based on wind tunnel and flight test is difficult to reflect the aerodynamic characteristics of all the flight states or motion forms we are interested in, and the high-fidelity CFD calculation may take a lot of time to obtain enough data, so it is difficult to directly use these aerodynamic models in flight simulation. In addition, when the parameters in the aerodynamic system increase, the number of sample flight states and the calculation requirement for modeling also increase, which makes many parameters cannot be identified. Compared with the traditional methods, these intelligent algorithm models are based on “sample data” and do not care about the flow characteristics, so they do not have a physical background and cannot provide more information for the subsequent application. Therefore, it is necessary to improve the intelligent algorithm to make it suitable for practical application.
To sum up, the aerodynamic model based on the traditional mathematical method needs full understanding of the flow characteristics. The stability of the model can meet the requirements, but the modeling process is very complex. Although the aerodynamic model based on the intelligent algorithm does not rely on physical equations, the generalization ability of the intelligent model is often unsatisfactory. Both of the two methods have their own shortcomings, which makes them difficult to be applied to engineering practice.
To solve this problem, this paper proposes a reduced order model (ROM) for unsteady aerodynamics, which combines the intelligent algorithm with a traditional mathematical model and can consider the aerodynamic nonlinearity at a high angle of attack. Specifically, the intelligent algorithm based on NARX is used to build the structure of nonlinear unsteady aerodynamic ROM. Then, the correction functions based on the aerodynamic physical characteristic assumption of the modified Leishman-Beddoes dynamic stall (LB for short) model [25] are added into the intelligent model, and the parameters of the correction functions are identified by aerodynamic data and intelligent algorithm. Finally, the nonlinear aerodynamic ROM is obtained.
For most of the existing data-driven aerodynamic ROM, the dataset comes from the same data source: numerical simulation, experiment, or flight testing. However, aerodynamic data often comes from different sources in practice, and the cost and accuracy are also different. In order to handle the balance between model accuracy and data generation cost, data-driven modeling driven by data fusion is receiving more and more attention [26]. Here, data fusion refers to the process of combining data and information from multiple sources in order to refine, estimate, and get a better understanding of the data [27], and data-driven modeling based on data fusion mainly refers to the multifidelity modeling method based on aerodynamic data with different fidelity. For aerodynamic data, generally, high-fidelity aerodynamic data are obtained through the flight test, wind tunnel test, or direct numerical simulation, while low-fidelity data usually comes from numerical simulation and simplified, such as rough discretization, relaxation convergence tolerance, low precision, and omission of physics [28]. For data with different fidelity, the cost of data acquisition becomes higher with the increase of data fidelity, and due to the limitation of cost, the amount of available data decreases with the increase of fidelity. In addition, it is difficult for high-fidelity data to include all interesting information, which makes the aerodynamic modeling based on these data expensive and unsatisfactory. For example, most experimental data come from harmonic motion rather than random motion. That is because random motion is vulnerable to external noise interference, and the data quality of harmonic motion is better. However, the dynamic information contained in harmonic motion is limited. Based on the above problem, the multifidelity method has been developed and widely used, which provides a general framework for fusing data with different fidelity in the modeling process. The basic idea of the commonly used multifidelity framework is to use low-fidelity data or models (physically driven or data-driven) to provide the dynamic trend, while high-fidelity data is used to correct the model to better reproduce the high-fidelity results. In conclusion, the multifidelity model can obtain high-fidelity results with less cost.
Multifidelity methods have been used in many fields, such as mechanism analysis, optimization design [29–31], statistical inference, and uncertainty quantification [32]. In the field of fluid dynamics, some scholars have applied the multifidelity method to flapping wing dynamics analysis [33, 34], aerodynamic optimization [35–38], flight simulation [28, 39, 40], hypersonic aerodynamic load prediction [41], low-fidelity turbulence model correction [42, 43], and uncertainty quantification of fluid dynamics system [44–47]. It is worth noting that most of these studies focus on the evaluation of steady-state aerodynamics to achieve rapid simulation and optimization, while there are few studies on unsteady multifidelity aerodynamic modeling. Because the generation of unsteady aerodynamic data is much more expensive than steady aerodynamic data, it is particularly important to explore the effectiveness of the unsteady aerodynamic multifidelity method. Ghoreyshi et al. [48] firstly showed that using both Euler and RANS data for modeling can improve the prediction results of aerodynamic ROM, and the results are better than that of the model based on a single data source. High-fidelity, low-fidelity, and multifidelity ROM are built by the NARX framework. The multifidelity ROM is constructed by adding low-fidelity time delay prediction as additional inputs of high-fidelity ROM. After that, Kou and Zhang [49] took the cokriging concept as the underlying theory and systematically introduced the problem of unsteady aerodynamic data fusion modeling. They built the multifidelity aerodynamic ROM with improved stability of NACA0012 airfoil under transonic conditions by using the NARX framework. The correction term of the model is constructed by the multicore neural network. Euler and unsteady RANS solvers are used to provide low-fidelity and high-fidelity data, respectively. Results show that only three groups of high-fidelity training cases of harmonic motion are needed to establish a high-fidelity model and obtain high-fidelity results.
In this paper, in order to improve the modeling efficiency of the ROM, the multifidelity data fusion modeling method is also introduced, using a small amount of high-fidelity data and a large amount of low-fidelity data to build the high-fidelity model. Using data fusion modeling is aimed at improving the efficiency of data acquisition as the cost of high-fidelity data is usually pretty high, which can also improve the modeling efficiency at the same time. In this paper, three kinds of data with different fidelity are used for fusion modeling. The low-fidelity data is obtained by the CFD solver based on Euler. It should be noted that these low-fidelity data can also reflect the real situation under certain restrictions. That is to say, although it is low-fidelity data, it still meets the accuracy requirements. The medium-fidelity data is obtained by the CFD solver based on RANS. The high-fidelity data comes from the existing experimental data. For the modeling process, firstly, the low-fidelity aerodynamic ROM is established by using low-fidelity data. Then, the nonlinear correction function is constructed based on the flow separation principle [25] and is added to the low-fidelity ROM. After that, the parameters of the nonlinear correction function are preliminarily identified by using the medium-fidelity data and the intelligent parameter identification method. Finally, the high-fidelity data from experiments is used to optimize the correction parameters, and the final high-fidelity nonlinear aerodynamic ROM is obtained. In order to verify the proposed high-fidelity nonlinear aerodynamic ROM, the lift and moment coefficients of NACA0012 airfoil in pitching and plunging motion are predicted and compared with the experimental data and CFD-RANS results. Then, the high-fidelity ROM is applied to the two-dimensional aeroelastic system. After that, the convergence of the model is analyzed, and its performance is compared with the traditional modified LB model and the intelligent aerodynamic model based on the RBF neural network.
2. Multifidelity Framework for Nonlinear Unsteady Aerodynamic ROM Materials and Methods
In this part, the multifidelity framework for nonlinear unsteady aerodynamic ROM based on traditional and intelligent algorithm fusion is introduced. It can be divided into the following parts. First of all, it is necessary to select a suitable multifidelity modeling method to develop the relationship between the low-fidelity model and the high-fidelity model. Secondly, the low-fidelity model is established. Thirdly, correction functions are constructed for the difference between the high-fidelity model and the low-fidelity model, and the parameters are identified. Finally, the nonlinear unsteady multifidelity aerodynamic model suitable for a certain speed range and high angle of attack can be obtained.
2.1. Multifidelity Modeling Method
The multifidelity model, that is, the models with different fidelity and computational efficiency, is combined in a certain way, so that the new model can obtain the accuracy of the high-fidelity model with less computational cost. The low-fidelity model can be divided into three categories, including the simplified model, data fitting model, and projection-based model [32]. Here, we use the data fitting model which is also a kind of data-driven model. To combine the low-fidelity model and the high-fidelity model, the method based on correction is adopted; that is, the low-fidelity model is modified by the data generated by the high-fidelity model. Currently, the main correction methods can be briefly divided into three categories [50]: (1)Additive and multiplicative corrections:(2)Comprehensive corrections:(3)Space mapping (input correction):where is the high-fidelity model and is the low-fidelity model. represents the multiplicative corrector while represents the additive corrector. is the input vector.
In this paper, comprehensive corrections are used. However, there are some differences with the standard form. We divide the input factors into two categories according to their effects on the nonlinear response of the system: . We think represents the input factor which can cause strong nonlinearity of the system response, is also a part of input factors, but the nonlinear response of the system caused by is very small. According to the above assumption, we divide the high-fidelity model and the low-fidelity model into two parts, respectively: where and are the total response of the system calculated by the high-fidelity model and the low-fidelity model, respectively. and are the response of the system caused by , and and are the response of the system caused by . Superscripts 1 and 2 are identifiers without mathematical meaning which represent the response part caused by and , respectively.
Because the previous assumption that the nonlinear response of the system caused by is very small, we think that the main difference between the high-fidelity model and the low-fidelity model is the difference between and , so only this part needs to be modified. And there is little difference between and , so we regard . Then, the relationship between the high-fidelity model and the low-fidelity model can be written as the following equation: where the multiplicative corrector and the additive corrector are both time-varying functions, which can be identified from high-fidelity data.
2.2. Reduced Order Modeling Methodology
2.2.1. Nonlinear Autoregressive with Exogenous Input (NARX) Description
Considering the unsteady time delay effect of aerodynamic force, the nonlinear autoregressive model (NARX) modeling method with exogenous input for typical nonlinear dynamic systems is used [51, 52]. For MIMO discrete-time systems, the NARX model can be written as follows: where is the time scale, is the system output at time , and is the system input at time . NARX is aimed at establishing the nonlinear relationship between the outputs and the inputs at the current and several previous time steps.
In the current study, the nonlinear aerodynamic ROM for two-dimensional airfoil with pitching and plunging freedom at low speed and high angle of attack is modeled; the two-dimensional aeroelastic system is shown in Figure 1. In this case, . and represent the lift coefficient and pitching moment coefficient at time , respectively. , represents the pitching angle of airfoils, while represents the plunging displacement. Let , can be also written as .
2.2.2. Model Structure
It is considered that the nonlinearity of the aerodynamic model mainly comes from the nonlinear relationship between the input and output of the static part, while the dynamic part can be represented by a linear system. Therefore, the input variables in formula (7) are divided into two groups, and the input variables related to the static part are separated:
It should be noted that when in is sufficiently larger, the time delayed in the first- and second-order derivatives in may not be needed. However, in the modeling process, with the increase of the value of , the model structure will become more complicated, and the parameters that need to be identified will increase too, which will influence the modeling efficiency and model stability. Therefore, generally, of the model is not necessarily a large value. At this time, the first-order and second-order derivatives are still needed as supplementary parameters to show the dynamic characteristics of the system.
According to the multifidelity modeling method, low-fidelity and high-fidelity models can be written as follows: where and are the multiplicative and additive correction functions of the high-fidelity model, respectively, and is the set of correction parameters to be identified. and are low-fidelity linear models, and is high-fidelity nonlinear models. In the current study, we can further expand them as the following functions: where and are the parameters to be identified corresponding to different input variables. and represent the order of the model, reflecting the unsteady effect of the system. Formula (9) can be written as
It should be noted that in order to improve the accuracy and efficiency of identification, the parameters of the low-fidelity model and the high-fidelity model are identified, respectively. The parameters of the low-fidelity model are identified firstly. The parameter set of the low-fidelity ROM can be written as follows:
And the output of the low-fidelity ROM can be written as where and represent the lift coefficient generated by input parameters and , respectively, while and represent the pitching moment coefficient generated by input parameters and , respectively.
2.2.3. Nonlinear Correction Model Based on Flow Separation
For the multifidelity model, some scholars use a neural network to model the nonlinear correction term [49]. However, data-driven models generally have the problems of easy overfitting, poor generalization ability, and poor robustness, and a large number of high-precision calculation data are required to ensure the accuracy and generalization ability of the model, which reduces the modeling efficiency, so it is difficult to be applied to engineering practice. The traditional mathematical expression modeling method can avoid the above problems, because it is based on the physical characteristics of the flow field. Therefore, according to the physical characteristics of the flow field, the nonlinear aerodynamic correction terms and at low speed and high angle of attack are modeled.
For the case of a high angle of attack, its nonlinearity is mainly caused by the separation of airflow. The influence of flow separation on airfoil aerodynamics can be divided into two aspects. The first is that the flow separation can directly change the aerodynamics, and the second is that flow separation causes the generation of vortex which leads to additional aerodynamics. and correspond to the two aspects, respectively. In addition, when the pitching motion with a high angle of attack is carried out at a low Mach number, the separation and reattachment between the flow and wing surface will be delayed, which will bring additional aerodynamic force called “overshoot” [25, 53]. Here, the influence of “overshoot” is represented by . To sum up, for the high-fidelity ROM model at a higher Mach number (generally and ), it can be written as formula (11). For the high-fidelity ROM model at a lower Mach number (generally ), it can be written as
, , and can be further expanded according to the output:
(1) Modeling for . Some scholars have found that the aerodynamic correction term directly caused by flow separation can be written as a function of the separation point [25]. The form of the function is shown in formula (16). Here, we directly use the form of this function, but the parameters in the function are not fixed and need to be identified according to the actual situation: where are the parameters to be identified related to the airfoil section. is the separation point at time , which can be written as the following iterative form: where , . is the velocity of the incoming flow, is the time step, and is the chord length of the airfoil. is the time constant of the dynamic separation point and is one of the time parameters to be identified. can be written in the following forms: where are the parameters to be identified according to different airfoils, and corresponds to the critical separation angle of attack of airfoils. are the parameters to be identified related to the separation point. is the effective angle of attack at time .
(2) Modeling for . The correction term of the additional lift and moment generated by the separation vortex can also be written as a function of the separation point, which is also affected by the time parameter . The additional lift generated by the separated vortex can be calculated by the following equation: where
is the time parameter of separation vortex failure. is the time constant of the vortex moving on the whole airfoil.
The correction term for the pitching moment coefficient generated by vortex can be given by the following formula: where is the parameter to be identified.
(3) Modeling for . For the special “overshoot” lift at a low Mach number, the lift correction factor of this part is also the function of the separation point, which can be written as where is the parameter to be identified related to the airfoil, is the initial angle at which the airflow readheres to the airfoil when the angle of attack decreases, and are the shape functions of the lift coefficient generated by the “overshoot.” The expression is as follows [54]:
The additional moment coefficient on the wing surface caused by the overshoot effect can be calculated by the following formula: where and are the parameters to be identified.
To sum up, the nonlinear correction parameter set can be written as , and the parameters to be identified can be divided into three parts: airfoil parameter , separation point equation parameter , and separation-related time coefficient . Let
Then, the parameter set to be identified for high-fidelity ROM is .
In conclusion, the construction of the multi-fidelity model can be divided into two steps. Firstly, identify by low-fidelity data and then obtain the parameter set of the correction equation based on high-fidelity data. In addition, it should be also noted that when identifying and , airfoil chord length and flow velocity are used for dimensionless processing. Therefore, chord length and flow velocity are also used as additional input parameters of ROM, which broadens the application scope of the model. The parameter identification algorithm can be used to identify the above parameters. The next section will introduce the parameter identification algorithm in detail.
2.3. Parameter Identification Methodology
In this paper, a parameter optimization algorithm based on cuckoo search (CS) is used for parameter identification. CS, proposed by Yang and Deb in 2009 [55], is a typical swarm intelligence optimization algorithm with iterative search characteristics. The CS algorithm originated from egg parasitism of cuckoo species, and it simplified the biological behavior into three idealized hypotheses: (I) each cuckoo lays one egg at a time and places it in randomly selected host nests, (II) high-quality eggs and the best nest will be passed on to the next generation, and (III) the number of available host nests is fixed, and the probability that host birds may find foreign eggs is .
Each egg is represented as a solution, and the last hypothesis indicates that some host nests will be replaced by new random nests. For the maximization problem, the quality or fitness of the solution can be simply proportional to the value of the objective function. Other forms of fitness can be defined in a way similar to fitness function in a genetic algorithm. In terms of the search algorithm, the CS algorithm does not search data through a simple isotropic random walk but updates data through the combination of local random walk and global random walk. Controlled by switch parameters , local random walk can be written as follows: where and are two solutions generated by random walk, is the heaviside function, is the random number with uniform distribution, and is the step size. is the step size scaling factor. The global random walk adopts the Levy flight mode, which can be written as follows: where can be written as
However, the above algorithm uses a random number to initialize the position of the bird’s nest, which is blindness. The position of the random solution may be far away from the correct solution. Meanwhile, the data search and update mechanism of Levy flight reduces the convergence speed in the iterative process, and it may not be able to find the solution with the required accuracy in a certain period of time. In order to solve this problem, chaos optimization is introduced. The randomness and ergodicity of chaos can avoid the problem of stopping searching when the local minimum is found, which can overcome the shortcomings of traditional optimization algorithms. In this paper, the logistic map, which is commonly used in chaos optimization, is introduced into the solution updating process of the CS algorithm: where is the chaos control variable. Before chaos optimization, the variables to be optimized need to be mapped to .
The cuckoo search algorithm with chaos optimization can be written in the following steps: (1)Define the number of initial nests and the initial solution (2)Obtain a new stochastic solution by Levy flight and calculate its fitness (3)Compare with the fitness of the other solution in the nest. If , the solution in the nest is replaced by the new solution ; otherwise, the solution in the nest is retained. Find the best solution of this iteration as (4)Project to and obtain the chaotic sequence by Logistic mapping. Then inversely map the chaotic sequence to the original solution space to obtain a new solution and calculate the fitness of the new solution. Repeat the above process several times and find the optimal solution in as . By comparing with , the final best solution of this iteration is obtained(5)Discard a certain number of nests by probability and replace them with new nests. Retain the optimal solution after the replacement(6)Judge whether the optimal solution meets the accuracy requirements. If it meets the requirements, stop the iteration process and output . Otherwise, return to step 2
The whole process can also be illustrated by the flow chart shown in Figure 2.
It should be noted that for parameters with physical meaning, the optimization range of these parameters will be given according to the conclusions of some literatures such as Ref. [25, 54] to meet the physical meaning.
3. Model Training and Validation
The NACA0012 airfoil is taken as an example for model training and validation. Firstly, different fidelity data sources are introduced; then, the model training process is represented, and finally, the model is verified in aerodynamic prediction and aeroelastic analysis.
3.1. High-Fidelity, Medium-Fidelity, and Low-Fidelity Data
Theoretically, establishing a multifidelity model only requires two kinds of datasets with low and high fidelity and two steps. For the first step, the parameters of the low-fidelity model are identified through the low-fidelity data. For the second step, the parameters of the high-fidelity model are obtained according to the high-fidelity data. However, for the multifidelity aerodynamic model in this paper, the high-fidelity data obtained by experiments is limited, and the airfoil motion corresponding to the experimental data is usually harmonic motion, which makes the dynamic information of the system in high-fidelity data insufficient. In addition, there are many parameters to be identified in the high-fidelity model; directly using experimental data to identify the nonlinear parameters may cause overfitting. The generalization ability of the model will be very poor. So, it is not suitable to use limited experimental data alone for nonlinear modeling. Therefore, it is necessary to obtain more high-fidelity aerodynamic data with more dynamic information of the system for modeling. The CFD-RANS model is a good choice to obtain more aerodynamic data with relatively high fidelity. On the one hand, the CFD results based on the RANS solver can reflect the nonlinear aerodynamic characteristics of a high angle of attack to a certain extent. On the other hand, the airfoil motion form can be given randomly according to the actual requirements, which is not limited by the experimental conditions. As a result, CFD data based on the RANS model is chosen as the third kind of data with relatively high fidelity (here defined as medium fidelity) for multifidelity aerodynamic model modeling, which can be used for the preliminary identification of nonlinear high-fidelity model parameters. For the multifidelity aerodynamic model established in this paper, the modeling process is divided into three stages. Firstly, the low-fidelity model is identified through the low-fidelity data; then, the high-fidelity model parameters are identified for the first time through the medium-fidelity data, and finally, based on the results of the first identification, the high-fidelity model parameters are modified for the second time through the high-fidelity data from experiments.
As just described, in this study, three kinds of data with different fidelity are used in modeling, which are called high-fidelity data, medium-fidelity data, and low-fidelity data, respectively. Here, low-fidelity data is obtained by CFD-Euler solvers. Medium-fidelity data is obtained by the CFD-RANS solver. High-fidelity data is obtained by experiments. It should be noted that the final model will be modified by high-fidelity experimental data, so the CFD-RANS model only needs to reflect the trend of aerodynamic force at a high angle of attack. The result of the RANS solver does not need to be in perfect agreement with the experimental data, which further saves the time to adjust the parameters and mesh quality of the RANS solver, as well as improve the modeling efficiency. Therefore, the aerodynamic force calculated by RANS is taken as the relatively high-fidelity data and here is defined as medium-fidelity data.
In this study, the ALE (Arbitrary Lagrange Euler) finite volume method is used as the governing equation. In this way, the motion of the material structure boundary can be effectively tracked in the calculation process, and the grid position can be appropriately adjusted according to the defined parameters in the calculation process without serious distortion. The double time step implicit time discretization method is used for time discretization, which can improve computational efficiency. For the turbulence model, SST is selected, which can be applied to the calculation of flow separation and reverse pressure gradient. The scheme of the computational grid deformation is based on diffusion, and the aerodynamic grid for CFD-Euler and CFD-RANS is shown in Figure 3. For the mesh for the RANS solver, the grid height of the first boundary layer is 10-3 mm, and the growth rate is 1.2.
(a) Mesh for the Euler solver
(b) Mesh for the RANS solver
3.2. Model Training
The training data of the unsteady aerodynamic model are usually obtained by two methods. One is to design a training signal with specified excitation (such as random signal, 3211 signal, and chrip signal), which needs to contain the dynamic information of the main motion state in the actual flight situation. The other is to select some representative states in the sampling space, where each state corresponds to a specific motion type, such as harmonic motion with constant frequency and amplitude [49]. Both of the two methods can be easily applied to numerical simulation, but the first method cannot be applied to the experiment due to the noise. For the nonlinear unsteady aerodynamic ROM established in this paper, on the one hand, it needs a large amount of data with a wide range of characteristics to meet the generalization ability requirement; on the other hand, it needs experimental data for secondary correction, so both of the two methods are used to obtain the training data.
As described in Section 3.1, the model training process can be divided into three stages. The first stage is the low-fidelity linear ROM training. The Euler solver is used to obtain the aerodynamic training data. It should be noted that for the Euler solver which ignores the viscosity, the motion form corresponding to the training signal should be limited in a certain range. Otherwise, due to the influence of viscous force, the obtained aerodynamic coefficients cannot meet the actual situation and cannot be used in aerodynamic modeling. The RANS solver does not have this problem because it considers the influence of viscosity. However, whether the training data is obtained by the Euler solver or RANS solver, it should meet the following requirement to achieve a better training result; that is, the training signal needs to contain a wider range of angle of attack (initial angle of attack is also included) and reduced frequency (the reduced frequency is defined as , where is the angular frequency of pitching and plunging motion). Therefore, several groups of chirp signals with variable amplitude as well as initial value are selected as training signals, one of which is shown in Figure 4. The frequency range of the signal is 0-5 Hz, and the amplitude range is set according to the actual situation; for the Euler solver, the amplitude of the effective angle of attack (defined as ) is limited to . Each group of training signals also corresponds to multiple working conditions with different incoming velocities, which can further widen the reduced frequency range. The aerodynamic training data corresponding to the signal in Figure 4 is shown in Figure 5.
(a) Lift coefficient under different inlet velocities
(b) Moment coefficient under different inlet velocities
The second stage is the first time training of the high-fidelity nonlinear ROM, which is also the initial identification of nonlinear correction function parameters. In this stage, the RANS solver is used to obtain the training data. The training signals are the same as described above which are also the chirp signals with variable amplitude, but the effective angle of attack changes from –30° to 30°. CFD data based on the RANS model is used for the preliminary identification of nonlinear correction parameters, and one or more groups of nonlinear correction parameters with stronger generalization ability can be obtained. A group of training data is shown in Figure 6. It can be seen that under the same inflow velocity, the nonlinearity of aerodynamic force is stronger at the motion with high amplitude and low frequency. In addition, the aerodynamic nonlinearity becomes more significant with the increase of the velocity, which may be the result of the decrease of the reduced frequency.
(a) Lift coefficient under different inlet velocities
(b) Moment coefficient under different inlet velocities
The third stage is the second time training of the high-fidelity nonlinear ROM, which is also a second time of optimization for nonlinear correction function. Based on the nonlinear parameters identified in the second stage, the experimental data in McAlister et al.’s [56], Sheng et al.’s [25], and Leishman’s [57] paper is used as the training data, and the cuckoo search algorithm with chaos optimization is used for the optimization. Finally, the high-fidelity nonlinear ROM with high accuracy and strong generalization ability is obtained.
For low-fidelity and medium-fidelity training datasets, six velocity points are selected as the training velocity, which are 10 m/s, 20 m/s, 40 m/s, 60 m/s, 80 m/s, and 100 m/s, respectively. There are one Euler training dataset and two RANS training datasets at each training speed. The training signal of each dataset is the chirp signal with variable amplitude, and the maximum amplitude of the signal is 5°, 20°, and 30°, respectively. The unit time step is 0.001 s, and the number of time points of each training set is 5000. The range of the Reynolds number of the training datasets is . For high-fidelity experimental training datasets, three groups of experimental data are used for modeling, and the parameters contained in the experimental data are listed in Table 1, corresponding to Case 3, Case 5, and Case 6. The reason for selecting Case 3, Case 5, and Case 6 as the training dataset is that the high-fidelity model needs to reflect the strong nonlinear characteristics of aerodynamic force. Therefore, three sets of data with stronger nonlinearity are selected as the training datasets. The dataset of Case 4 with relatively weak nonlinearity is used as the test dataset to test the generalization ability of the model.
Through the above training progress, the parameters of the low-fidelity model and high-fidelity model can be identified. For the low-fidelity model, the delay orders and must be determined first. Here, is determined by a trial-and-error method. After that, the parameters in Equations (12) and (26) can be identified by the cuckoo search algorithm with chaos optimization. In total, there are 39 parameters needed to be identified by multifidelity aerodynamic data. After the identification of the parameters, the low-fidelity and high-fidelity models are established.
3.3. Model Validation
3.3.1. Validation in Aerodynamic Prediction
By estimating the aerodynamic loads of NACA0012 airfoil under different pitching and plunging motions at low speed, the effectiveness of the unsteady nonlinear multifidelity modeling method based on data fusion is verified.
Firstly, several simple harmonic pitching motions are selected as test cases. The information and the working parameters of the test cases are listed in Table 1. Among them, case3-6 have experimental data and case1-2 are the test cases without experimental data. Use RMSE (root mean squared error) to measure the performance of different aerodynamic models. The calculation formula of RMSE is as follows: where is the time scale and means the number of time points of a group of test data. and represent the aerodynamic coefficients obtained by high-fidelity models/data and low-fidelity models/data, respectively.
Table 2 lists the RMSE of test case1-2 compared with the CFD-RANS result. Table 3 lists the RMSE of test case3-6, which is obtained by comparing with the experimental data. Figures 7 and 8 show the comparison of aerodynamic coefficients predicted by different aerodynamic models.
(a) Lift coefficient versus pitching angle of Case 1
(b) Moment coefficient versus pitching angle of Case 1
(c) Lift coefficient versus pitching angle of Case 2
(d) Moment coefficient versus pitching angle of Case 2
(a) Lift coefficient versus pitching angle of Case 3
(b) Moment coefficient versus pitching angle of Case 3
(c) Lift coefficient versus pitching angle of Case 4
(d) Moment coefficient versus pitching angle of Case 4
(e) Lift coefficient versus pitching angle of Case 5
(f) Moment coefficient versus pitching angle of Case 5
(g) Lift coefficient versus pitching angle of Case 6
(h) Moment coefficient versus pitching angle of Case 6
It can be seen from Figure 7 that for linear test Case 1 and Case 2 which pitching angle changes from [-10°, 10°] and reduced frequency is 0.46 and 0.23, respectively, the shape of the hysteresis loop is almost elliptical, and all models can effectively capture the unsteady characteristics of CL (lift coefficient) and CM (moment coefficient), which verifies the low-fidelity ROM and the high-fidelity ROM in the linear range. In addition, it can be seen from Table 2 and Figure 7 that the results calculated by the low-fidelity model and CFD-Euler model are almost the same. And Figures 7(b) and 7(d) show that the curve of the low-fidelity model is closer to the result of CFD-Euler, while the curve of the high-fidelity model corrected by CFD-RANS data is closer to the result of CFD-RANS, which also proves the effectiveness of the nonlinear correction model to a certain extent.
For Cases 3-6, the aerodynamic nonlinearity is very strong, and the low-fidelity model cannot meet the accuracy requirements. CFD-RANS can reflect the trend of nonlinear aerodynamics to a certain extent, but they are not accurate enough to fit the details. The proposed high-fidelity ROM improves this defect, as its parameters are optimized based on experimental data. Results show that the prediction results of the proposed high-fidelity ROM have the best consistency with the experimental data, which can be seen from the fitting degree of the hysteresis loop in Figure 8 as well as the RMSE data in Table 1. In addition, we also found from Table 1 and Figure 8 that with the increase of the initial angle of attack, the aerodynamic nonlinear characteristics will become more significant, and the consistency between CFD results and experimental data will become worse. This may imply that when the initial angle of attack or total motion amplitude increases to a certain value, the computational accuracy of CFD will no longer meet the requirements. The proposed high-fidelity ROM based on multifidelity data can avoid the above problem to a certain degree, because it can modify itself based on the experimental data, which makes it more extensible for different conditions.
In order to prove the generalization ability of the model, more random test signals are needed. We choose a group of fwgn signals with random parameters as the test signal. The motion of the test case is a combination of pitching motion and plunging motion. The test signal is shown in Figure 9. Figure 10 shows the comparison among the output results from different models and CFD solvers. It can be seen from Figure 10 that the established high-fidelity ROM well reflects the nonlinear aerodynamic characteristics, which proves its good generalization and extrapolation ability.
(a) FWGN signal for pitching motion
(b) FWGN signal for plunging motion
(a) Lift coefficient versus time
(b) Moment coefficient versus time
As shown in the test case, the proposed aerodynamic modeling method only needs a group of low-fidelity and medium-fidelity CFD data and a few high-fidelity data with a specified harmonic motion to establish high-fidelity aerodynamic ROM with satisfactory generalization ability.
3.3.2. Verification in Two-Dimensional Aeroelastic Response
In order to further verify the validity of the model, the nonlinear ROM is applied to two-dimensional aeroelastic analysis, and the results are compared with the experimental results in Ref. [58].
The two-dimensional aeroelastic system is already shown in Figure 1, which has two degrees of freedom: pitching motion and plunging motion. The important experiment parameters in Ref. [58] are listed in Table 3.
The dynamic equation of the system can be written as follows: where is the pitch angle and is the amplitude of plunge motion. and are spring stiffness of plunging and pitching freedom, respectively. and represent structural damping. is the mass per unit length. represents the mass static moment, and is the moment of inertia at quarter chord length. and are external aerodynamic forces and aerodynamic moments, which can be calculated by the proposed nonlinear aerodynamic ROM with the following equation: where and are aerodynamic coefficient and aerodynamic moment coefficient, respectively, which are calculated by nonlinear aerodynamic ROM. is chord length and is span length. is the distance between elastic wheelbase and quarter chord line. represents air density. represents flow velocity.
The Runge-Kutta method is used to solve Equation (32) in the time domain, and the bifurcation and limit cycle analysis can be carried out, as shown in Figures 11 and 12. The -axis parameter in Figure 12 can be expressed by pitching angle and plunging velocity .
(a) m/s
(b) m/s
It can be seen from Figures 11 and 12 that the results of two-dimensional aeroelastic bifurcation and limit cycle analysis based on nonlinear aerodynamic ROM proposed in this paper are in good agreement with the experimental results, which further verifies the effectiveness of the aerodynamic reduced order model established in this paper.
4. Model Analysis and Comparison
In the previous section, the proposed high-fidelity aerodynamic ROM based on the multifidelity framework is verified by the accurate prediction of the nonlinear unsteady aerodynamic coefficients as well as the simulation of two-dimensional limit cycle response induced by nonlinear aerodynamics. However, the performance of the model still needs to be investigated. This section studies the performance of the model from three aspects, respectively. Firstly, the convergence of model error with the velocity interval of training samples is studied. Then, the calculation accuracy of the proposed model is compared with the modified dynamic stall model (also called the modified LB model proposed in Ref. [25]) and nonlinear aerodynamic model based on RBF neural network. After that, the application scope and generalization ability of the proposed ROM model is explored and compared with other aerodynamic models.
4.1. Convergence Analysis
As described in Section 3.2, the modeling process of the nonlinear ROM can be divided into three stages. The first stage is to identify the low-fidelity model, and the second and the third stages are to identify the high-fidelity correction equation. However, due to the limited training data in the third stage, only the convergence of the first two stages is studied. The study of the convergence is divided into two parts also, according to different training stages.
The study is conducted as follows. Firstly, we need to define the velocity range of the training model. In this study, the velocity range is 10-60 m/s. Then, the training cases should be grouped according to different velocity intervals. Here, four training groups with different velocity intervals are set up; the corresponding speed intervals are 50 m/s, 25 m/s, 10 m/s, and 5 m/s. The velocity included in each training group is shown in Table 4. After that, the training data of each group are calculated, and the different ROM models are obtained after training. In order to eliminate the influence of signal forms, the training signals for each velocity point are the same, as shown in Figure 4. Finally, a set of test cases are selected, and the aerodynamic coefficients are predicted by different ROM models. The signal form of the test cases is shown in Figure 13, and the velocities of the test cases are listed in Table 5.
(a) Test signal for low-fidelity ROM
(b) Test signal for high-fidelity ROM
As can be seen from Figure 14, for the low-fidelity ROM, the RMSE of different test cases all decreases with the decreasing of velocity interval. However, it should be noted that reducing the velocity interval means more training data need to be obtained, which will increase the time cost for modeling. So, it is necessary to balance the accuracy and efficiency. For the case of this paper, the change of the RMSE is very limited when the interval of training velocity points is less than 10 m/s. Therefore, it is the best to select training group 3 as the training case for low-fidelity ROM.
(a) RMSE of lift coefficient
(b) RMSE of moment coefficient
For the high-fidelity ROM, as shown in Figure 15, RMSE of different test cases change differently with the decreasing of training velocity interval. For test Cases 3~5 where inlet velocity is less than 30 m/s, RMSE increases first and then reduces, while for other test cases with higher inlet velocity, the opposite is true. This is because the influence of gas viscosity on aerodynamic force is different from the change in flow velocity. When the flow velocity is pretty low, the viscous effect of the gas will be stronger and the influence on the aerodynamic force will be greater. Therefore, for the same pitching motion at 10 m/s and 60 m/s, there are some differences in aerodynamic characteristics. If the training data contains more cases at a higher speed, the aerodynamic characteristics reflected by the training data are more consistent with the aerodynamic situation at a higher speed. Therefore, for test Cases 1 and 2, the result of training group 2 is the best, as it contains only one case whose inlet velocity is less than 30 m/s, while the other two cases are both at a higher speed. After reducing the velocity interval and increasing lower incoming flow speed as training points, the model is adjusted to fit the aerodynamic characteristics at a lower velocity, so that the fitting degree at 30-60 m/s is reduced. However, compared with the improvement of model generalization ability, the loss of accuracy is acceptable. In addition, it can be seen from Figure 15(a) that for the lift coefficient, RMSE increases slightly when the velocity interval decreases to 5 m/s. This may be due to the overfitting of the model, which should be avoided as far as possible. Therefore, for the high-fidelity ROM, the best training velocity interval is 10 m/s, corresponding to training group 3.
(a) RMSE of lift coefficient
(b) RMSE of moment coefficient
4.2. Calculation Accuracy Analysis and Comparison
In this part, the calculation accuracy of the proposed high-fidelity ROM is compared with that of the modified LB model and aerodynamic ROM based on the RBF neural network. The training data of the aerodynamic ROM based on the RBF neural network is exactly the same as that of the high-fidelity model proposed in this paper. The time delayed orders of the RBF model are also the same with the high-fidelity ROM. The comparison cases are the cases with experimental data listed in Table 1. The comparison of three aerodynamic ROM is shown in Figure 16 and Table 6.
(a) Lift coefficient versus pitching angle of Case 3
(b) Moment coefficient versus pitching angle of Case 3
(c) Lift coefficient versus pitching angle of Case 4
(d) Moment coefficient versus pitching angle of Case 4
(e) Lift coefficient versus pitching angle of Case 5
(f) Moment coefficient versus pitching angle of Case 5
(g) Lift coefficient versus pitching angle of Case 6
(h) Moment coefficient versus pitching angle of Case 6
It can be seen from Figure 16 and Table 6 that the aerodynamic ROM based on the RBF neural network has the best fitting degree for the test case. This is because the RBF neural network has the ability to fit any nonlinear function, and the test cases are also in the training set, so the results of the RBF model are the best. However, the generalization ability of the RBF model needs to be further verified because of the limited training data, which will be explained in the next subsection.
The fitting degree of the proposed algorithm of test cases is basically equivalent to that of the traditional modified LB algorithm, which can be seen from Figure 16 and Table 6. However, this traditional modified LB model is based on a large amount of experimental data, and parameter determination is based on subjective experience, which costs a lot of time and brings a lot of uncertainty for modeling. Moreover, after changing the experimental parameters and airfoil conditions, it is necessary to redesign the experiment and obtain the data, so the modeling process of the LB model will cost a pretty long time. The aerodynamic ROM proposed in this paper only needs a group of CFD data with low-fidelity or medium-fidelity and few experimental data for modeling, which improves the modeling efficiency.
4.3. Generalization Ability and Application Scope Analysis and Comparison
In this section, the generalization ability and the application scope of the proposed high-fidelity aerodynamic ROM are explored and compared with the RBF model. Firstly, the training data are obtained in a certain speed range (10-60 m/s), and the unsteady aerodynamic ROM is established by the RBF neural network algorithm and the algorithm proposed in this paper, respectively. The signal of training data is the same as that in Figure 4. Five velocities (10 m/s, 20 m/s, 30 m/s, 40 m/s, 50 m/s, and 60 m/s) are selected as training points. Then, the two kinds of ROM are applied to the test cases beyond the training speed range, and the prediction results are compared with the CFD-RANS calculation results. There are two main reasons for taking CFD-RANS calculation results as the reference. Firstly, the experimental data at a higher speed is insufficient. CFD-RANS data is a good supplement as the relatively high-fidelity data, and it has been proved in the article that the CFD-RANS model can reflect the trend of nonlinear aerodynamic force to a certain extent. And the generalization ability can be also seen as the ability to follow the response of the aerodynamic trend under different conditions. Therefore, we think it is reasonable to take the CFD-RANS data as the reference, because CFD-RANS data can reflect the trend of aerodynamic force. Another reason for taking CFD data as a reference is that CFD data is the easiest to obtain. In order to more clearly reflect the variation of aerodynamic nonlinear characteristics with different incoming flow velocities, the chirp signal with variable amplitude is also used as the test case.
Figures 17 and 18 show the comparison of the aerodynamic coefficients predicted by the CFD-RANS solver and the two kinds of ROM under different inlet velocities of 100-250 m/s. It can be seen that the aerodynamic coefficients predicted by high-fidelity ROM proposed in this paper are in good agreement with the CFD-RANS results, which shows that the model has a strong generalization ability in the speed range of 100-250 m/s. However, the fitting degree of the RBF model to CFD data decreases rapidly with the increase of airspeed, which indicates that the generalization ability of the RBF model is really poor, and the RBF model needs more data in modeling. Therefore, the modeling efficiency of the proposed ROM is higher than that of the aerodynamic model based on the RBF neural network.
(a) 100 m/s
(b) 150 m/s
(c) 200 m/s
(d) 250 m/s
(a) 100 m/s
(b) 150 m/s
(c) 200 m/s
(d) 250 m/s
It is worth noting in Figures 17 and 18 that with the increase of the incoming flow velocity, the fitting degree of ROM proposed in this paper will become lower, which means that the ROM proposed in this paper may no longer be applicable when the incoming flow velocity is higher than 250 m/s. This is because when the velocity of the incoming flow continues to grow, shock waves will be generated locally on the wing surface, which will affect the aerodynamic characteristics. The proposed aerodynamic models cannot reflect the effect of shock wave, which means it is no longer applicable. Moreover, when the inflow velocity is less than 10 m/s, the adaptability of the ROM model also becomes worse due to the enhanced viscous effect of the gas. In general, the prediction accuracy of the lift coefficient is higher than that of the moment coefficient, because the nonlinearity of the moment coefficient is stronger. Based on the above analysis, we can draw the conclusion that the ROM based on algorithm fusion and multifidelity framework established in this paper is suitable for the airspeed condition at 0.03~0.8 Ma.
5. Conclusions
A multifidelity modeling framework for low-speed nonlinear unsteady aerodynamic ROM is presented in this paper, which combines traditional and intelligent aerodynamic modeling methods. The intelligent algorithm based on NARX is used to build the nonlinear aerodynamic model framework, and the traditional LB separation equation based on the flow characteristics is used to build the modified equation. The framework integrates data with different fidelity and extends the multifidelity method of nonlinear system identification to simulate the unsteady effects of nonlinear aerodynamics. The low-fidelity data and the medium-fidelity data are obtained by the CFD-Euler solver and CFD-RANS solver, respectively, and the high-fidelity data is obtained from the experimental results. The main conclusions drawn from the results can be summarized as follows: (1)The proposed high-fidelity ROM performs well in aerodynamic coefficient prediction and two-dimensional nonlinear aeroelastic analysis of NACA0012 airfoil at low speed and high angle of attack, which indicates that the proposed aerodynamic modeling method only needs a group of low-fidelity or medium-fidelity CFD data and few high-fidelity data with a specified harmonic motion to establish the high-fidelity ROM with satisfactory generalization ability. The results are in good agreement with the experimental data(2)Compared with the traditional modified LB aerodynamic model, the proposed high-fidelity ROM is more flexible and more adaptable(3)Compared with the aerodynamic model based on the RBF neural network, the proposed high-fidelity ROM has stronger generalization ability, higher stability, and modeling efficiency(4)Data fusion can further reduce the modeling cost and improve modeling efficiency. Thus, the efficiency of calculation and simulation is greatly improved compared with high-precision CFD and experiments(5)Future work includes improving the modeling and identification algorithm to better represent the correction term and involving it in three-dimensional configurations
Data Availability
The numerical simulation result data used to support the findings of this study are available from the corresponding author upon request.
Disclosure
This paper is a part of a research project supported by the Department of Aeroelasticity, School of Aeronautic Science and Engineering, Beihang University.
Conflicts of Interest
The authors declare that there is no conflict of interest regarding the publication of this paper.