Document
AIAA 2003-0653
Wind Tunnel Database Development
using Modern Experiment Design and
Multivariate Orthogonal Functions
Eugene A. Morelli
NASA Langley Research Center
Hampton, VA
Richard DeLoach
NASA Langley Research Center
Hampton, VA
41 st AIAA Aerospace Sciences Meeting and Exhibit
January 6-9, 2003 / Reno, NV
For permission to copy or to republish, contact the American Institute of Aeronautics and Astronautics, 1801 Alexander Bell Drive, Suite 500, Reston, VA, 20191-4344 AIAA-2003-0653 WIND TUNNEL DATABASE DEVELOPMENT USING MODERN EXPERIMENT DESIGN AND MULTIVARIATE ORTHOGONAL FUNCTIONS Eugene A. Morelli* and Richard DeLoacht NASA Langley Resem'ch Center Hampton, Virginia US,_ 23681- 2199 Nomenclature Abstract a parameter vector A wind tunnel experiment for characterizing the aerodynamic and propulsion forces and moments lift, drag, and side force coefficients Q,CD, q.
acting on a research model airplane is described. The rolling moment coefficient Q model airplane, called the Free-flying Airplane for pitching moment coefficient Sub-scale Experimental Research (FASER), is a C,.
modified off-the-shelf radio-controlled model yawing moment coefficient C.
airplane, with 7 ft wingspan, a tractor propeller covariance matrix Coy driven by an electric motor, and aerobatic capability.
J cost function FASER was tested in the NASA Langley 12-foot MDOE Modem Design Of Experiments Low-Speed Wind Tunnel. using a combination of traditional sweeps and modem experiment design. number of model terms n Power level was included as an independent variable N total number of data points in the wind tunnel test, to allow characterization of One Factor At a Time OFAT power effects on aerodynamic forces and moments.
PSE predicted squared error A modeling technique that employs multivariate power level, percent pwr orthogonal functions was used to develop accurate analytic models for the aerodynamic and propulsion modeling function vector P force and moment coefficient dependencies from the thrust force, Ibf T wind tunnel data. Efficient methods for generating x independent variable vector orthogonal modeling functions, expanding the measured output vector Y orthogonal modeling functions in terms of ordinary Oc angle of attack, deg polynomial functions, and analytical orthogonal sideslip angle, deg blocking were developed and discussed. The resulting models comprise a set of smooth, aileron deflection, deg differentiable functions for the non-dimensional elevator deflection, deg
8e
aerodynamic force and moment coefficients in terms flap deflection, deg of ordinary polynomials in the independent variables, _f suitable for nonlinear aircraft simulation.
rudder deflection, deg
8r
cr 2 variance * Research Engineer, Senior Member AIAA ordinary polynomial function vector _"Senior Research Scientist Copyright © 2003 by the American Institute of Aeronautics and superscripts Astronautics, Inc. No copyright is asserted in the United States T transpose under Title 17, US. Code The U.S. Government has a royalty- free license to exercise all rights under the copyright claimed estimate herein for Governmental purposes. All other rights are reserved by -1 matrix inverse the copyright owner.
normalized American Institute of A erona_tics and Astronautics subscripts FASER was designed so that the flight vehicle max maximum could be installed in the wind tunnel, see Figure 1.
min minimum This avoids any Reynolds number or scaling effects, and ensures that there are no physical differences o nominal between the wind tunnel model and the flight vehicle.
In contrast, full scale flight tests and drop model Introduction tests are expensive and are sometimes separated by months (and even years) for a particular research Modem aeronautical research involves expanded activity. The number of these tests is always tightly flight envelopes, as a result of the desire for improved constrained by budget. There is also a substantial fighter maneuverability for tactical advantages, and the difference in the cost of overhead if the aircraft is to be desire to improve flight safety. The expanded flight kept in flyable condition. Since FASER is inexpensive envelopes involve nonlinear aerodynamics which must and unmanned, risks can be taken in research and be modeled accurately.
development that could never be tolerated in a piloted Since nonlinear aerodynamics are much more flight test or even in a drop model test. Advances in complex than linear aerodynamics, more sophisticated instrumentation now make it possible to instrument a experimentation is required to accurately characterize sub-scale model aircraft with research-quality, miniaturized flight test instrumentation at a reasonable the functional dependencies. Nonlinear aerodynamics cost.
violate linear modeling assumptions such as superposition, quasi-steady flow, and no This paper describes the experiment design, data interdependence of the effects of states and controls. In analysis, and mathematical modeling involved in addition, aircraft designs have evolved with increasing developing a wind tunnel database for the FASER numbers of control effectors. Traditional wind tunnel aircraft. Accuracy of this database is critical for the testing methods would set each control effector to development of a high-fidelity nonlinear simulation to different fixed levels while sweeping through angle of be used for control law design, flight envelope attack and sideslip angle, for example. With such an expansion, flight experiment design, and pilot training.
approach, the number of data points required increases A preliminary version of the nonlinear simulation for exponentially with the number of control effectors, if FASER has already been developed, using U.S. Air information on control surface interaction effects is Force DATCOM to generate an aerodynamic model, desired. These considerations highlight the need to with experimentally-determined values for the mass develop more efficient wind tunnel testing and and inertia characteristics of FASER 1. The work modeling techniques to accurately characterize described in this paper will upgrade the aerodynamic nonlinear aerodynamics, with possible interaction model with analytic models derived from wind tunnel effects among a large number of independent variables.
data, add an engine thrust model, and include propulsion effects on the aerodynamics. Since FASER At NASA Langley, the Free-flying Airplane for was intended to be a research vehicle from the outset, Sub-scale Experimental Research (FASER) is being developed to study problems such as that described the approach used for the experiment design and data above. FASER is a modified off-the-shelf radio control analysis for the wind tunnel testing was not traditional.
model called the Ultra-Stick, manufactured by Hangar This paper explains how the wind tunnel testing was 9, Ltd., see Figure 1. FASER has a conventional high- done, and examines the results. The paper also wing and tail configuration with 7 ft wingspan, a describes a method for generating orthogonal modeling foldable tractor propeller driven by an electric motor, functions based on the independent variable data, along and acrobatic capability. subsequent expansion of the orthogonal modeling functions in terms of ordinary multivariate The purpose of FASER is to provide an polynomials. This method is slightly different from inexpensive aircraft for developing and demonstrating that described in Refs. [2] and [3], and represents an advanced experiment design, data analysis and evolutionary improvement of the technique. In Ref.
modeling techniques, and control law design methods.
[2], the concept of response surface modeling using As long as the goal is technology demonstration or multivariate orthogonal functions was successfully basic research, a sub-scale model that is not applied to inference subspaces for limited ranges of dynamically scaled for a specific full-scale aircraft is a angle of attack and Mach number with fixed sideslip completely acceptable test vehicle for these purposes.
angle and control surface deflections. This paper American Institute of Aeronautics and Astronautics extends the multivariate orthogonal function modeling nt;n-dimensional coefficients, get an idea of the concept to identify aerodynamic models for a large response levels, collect basic static stability and trim infbrmation, and define the boundaries of the flight envelope, with more independent variables.
independent variable subspaces to be used for further experimentation. OFAT sweeps were used because the Experiment Design data acquisition system in the NASA Langley 12-foot For this wind tunnel test, the fundamental Low-Speed Wind Tunnel is set up to collect OFAT data objective was to find a mathematical description for the in an automated fashion, making the sweeps very dependence of non-dimensional aerodynamic force and efficient in terms of collecting data points in minimum moment coefficients on independent variables that are time. However, the initial OFAT sweeps are really varied during the experiment. Each mathematical only intended for qualitative use, namely to define the description or model can be thought of geometrically as b_,tmdaries for subspaces that will be the focus of a hyper-surface, also called a response surface. Critical detailed experimentation and modeling in procedure issues for successfully identifying an adequate response four. One advantage of operating in this manner is that surface model from experimental data include the any issues related to instrumentation, data collection, or experiment design (or, how the independent variable experimental procedures can be worked out during the values are set when measuring the output responses), OFAT sweeps without impact on the experiment.
noise level on the measured outputs, identification of a because the data from the OFAT sweeps is being used mathematical model structure that can capture the for qualitative purposes only. This approach also functional dependence of the output variables on the provides a good rough overview of the landscape' that independent variables, accurate estimates of unknown characterizes the dependence of force and moment parameters in the identified model structure, and the coefficients on the independent variables.
ability of the identified model to predict outputs for The independent variables for the FASER wind data that was not used to identify the response surface model. ttmnel test were angle of attack a, sideslip angle ft, p_)wer level pwr, elevator deflection 8 e , aileron The experiment design used for FASER wind deflection _a, rudder deflection _r, and flap tunnel testing was a hybrid design broken down into a series of procedures. The procedures are listed in deflection _y-. The response variables were Table 1. Randomization 4-6 was used throughout the n_n-dimensional aerodynamic coefficients for lift, drag, testing, to separate independent variable effects from a:ld side forces (Cz,CD, andCy), and rolling, time-dependent systematic errors.
pitching, and yawing moments (CI.C m, and Cn).
The first procedure consisted of randomized Each data point produced measured values for all engine power sweeps with the wind tunnel air off, to independent and response variables.
determine the static thrust from the electric motor and the propeller. All of these runs were made with the The experiment was designed assuming an), of the model at zero angle of attack and zero sideslip angle, so independent variables could influence any of the that the thrust measurement was obtained from the response variables. It was assumed (as an initial guess longitudinal force measured by the balance mounted in oaly) that the dependencies could be modeled the model, Figure 2 shows the static thrust plotted as a accurately with polynomial terms in the independent function of pulses per second from a Hall effect variables of order 3 or less within each independent transducer on the electric motor, which is proportional variable subspace. In addition, it was assumed that to the propeller RPM. The model shown in Figure 2 is 14mgitudinal controls ( 8 e , _f. and pwr) do not interact the result of a simple least squares fit of thrust to pulse with lateral/directional controls (t_ a and 8 r ). In count, using a quadratic model structure. This model practical terms, this meant that longitudinal and structure was identified automatically from the data, lateral/directional controls were not varied using the orthogonal function modeling technique simultaneously to allow their mutual interaction effects described later.
to be quantified. As a result, the subspaces were called In the second and third procedures, the approach "longitudinal" if the longitudinal controls were moved, was to use One Factor At a Time (OFAT) sweeps, and "lateral/directional" if the lateral/directional wherein one independent variable is changed with all controls were moved.
others held constant, to characterize the general Based on experience with similar airplanes, it was topology of the response surfaces for the l,nown that the dependence of non-dimensional American Institute of A eronautics and Astronautics aerodynamic force and moment coefficients on control independent variable values are found by mapping the surface deflections could be modeled with low order independent variable values in engineering units for polynomials for the entire range of control surface each subspace onto the interval [-1,1]. The deflections. With that assumption, it was not necessary normalization of each independent variable was to vary the control surface deflections to search for implemented by inference subspace boundaries along the dimensions of the independent variable space corresponding to control X -- Xmm ) (1) surface deflections. Inference subspace boundaries "r =-1+ 2 (xma x -Xmin ) were therefore sought only for _ fl, and pwr. These independent variables were varied using OFAT sweeps where _ was the normalized value of the independent to identify the inference subspace boundaries.
variable, and the independent variable range in Figure 3 shows an OFAT sweep on angle of engineering units was [xmi,,Xm_ ]. The inverse attack. The vertical lines mark the selected subspace transformation was boundaries in angle of attack, which are intended to mark the boundaries of regions where the character of (J+l) the response surfaces change. There is a trade-off in (2) X = Xmin +T( Xmax --Xmi,1 ) selecting the subspace boundaries, in that more subspaces mean more individual experimentation regions in procedure four, while fewer subspaces All modeling for the inference subspaces was generally require more resources in the data collection done using normalized values of the independent and modeling for each subspace. Figures 4 and 5 show variables. The fmal models used for prediction were OFAT sweeps on sideslip angle and power level, with written in terms of engineering units.
the selected subspace boundaries marked by vertical Second-order central composite design 5,6 in four lines. All control surface deflections were zero for the independent variables (for the power-off subspaces) or sweeps shown in Figures 3-5.
five independent variables (for the power-on The full independent variable ranges tested are subspaces) was used in each subspace, augmented with listed in Table 2. Tables 3 and 4, and corresponding a 3ra order D-optimal design 5,6 in the appropriate Figures 6, 7, and 8, show the definitions of the number of independent variables. A two-dimensional inference subspaces in terms of boundary values of projection of this constellation of data points in angle of attack, sideslip angle, and power level. The normalized independent variable space is shown in symmet D, of the vehicle was used to justify omitting Figure 9. The central composite design points occupy testing in most of the subspaces with high negative the comers, the centers of each face, and the center of sideslip angles, see Figures 6, 7, and 8. All the the normalized subspace, while the D-optimal points subspaces together comprised the full inference space, generally fill in between. Although some of the defined by the full range of the independent variables in D-optimal points are coincident with the central Table 2. All control surfaces were tested over their full composite design points, the number of times that this physical deflection ranges for each subspace.
happens is not represented accurately in Figure 9, because of the projection onto two dimensions. This In the fourth and final procedure, Modem Design experiment design enabled identification of up to 3 ra Of Experiment (MDOE) techniques 4-6 were applied to order functional dependencies and interaction effects.
each defined subspace in order to obtain the most Provisions were made to augment the designed accurate and complete characterization of the functional experiment if the data indicated a lack of fit that dependencies, and also to make sure all interaction required modeling functions with higher than 3_aorder.
effects were adequately modeled. Refs. [2], [4]-[6] describe and demonstrate in detail the advantages of the MDOE approach compared to traditional OFAT for Instrumentation and Data Collection detailed modeling of the functional dependencies, both FASER was used as the wind tunnel model and in terms of the modeling accuracy and in the economy tested at a nominal flight speed, thus avoiding scaling, of experimentation resources required to arrive at an Reynolds number, or geometric dissimilarity issues for acceptable result.
comparisons with future flight test data. All control Within each subspace, independent variables were surfaces were instrumented with potentiometers. Air set according to normalized values. Normalized data vanes and airspeed pinwheels were installed on American Institute of Aeronautics and Astronautics Modeling, booms attached at each wing tip and extending 1 chord length in front of the wing. The air data sensors, which Typically, once the experimental data are will be used for flight testing, were calibrated as part of collected, polynomials in the independent variables are the wind tunnel experiment, since aerodynamic used to model the functional dependence of the output incidence angles and airspeed were carefully controlled variables on the independent variables, and the model and measured in the wind tunnel. The wind tunnel parameters are estimated from the measured data using balance was installed near the e.g. of the airplane in the least squares linear regression 5,6. The question of space normally occupied by the accelerometer and rate which polynomial terms should be included in the gyro package during flight test operations.
model for a given set of data, called model structure Control surface deflections and power level were determination, gets more difficult with increases in the automatically set to the values required by the number of independent variables, increased ranges for experiment design via a serial port interface from a the independent variables, or increased complexity of laptop computer in the control room to a the underlying functional dependency.
servomechanism controller onboard on the airplane.
Various hypothesis testing techniques 6,7 can be The same onboard equipment will be used to command used to identify an adequate model structure, but these the control surfaces and power level during flight methods are iterative and require the involvement of an operations. Angle of attack and sideslip angle were set experienced analyst. Neural networks using radial from the control room using servomechanisms driving basis functions with subspace partitioning, or back the movable C-strut in the test section. The angle of propagation with layered and interconnected nonlinear attack and sideslip angles were set automatically during activation functions, have also been applied to the OFAT sweeps, but had to be set manually for the response surface modeling problem 8. For this type of MDOE data points, using joystick controllers and a approach, there is a loss of physical insight and a measurement feedback to the control room. Dynamic danger of over-fitting the data, because the model pressure in the tunnel was regulated to 2 psf by an structures contain many parameters, typically with no automatic closed loop control on the wind tunnel fan mechanism for limiting the size of the model other than motor speed.
the judgment of the analyst.
The experimental set-up was designed to In this work, a nonlinear multivariate orthogonal accommodate the MDOE approach, which typically modeling technique 2,3 was used to model response requires changes in more than one independent variable surfaces from wind tunnel data. The technique for successive data points. Since the control surface generates nonlinear orthogonal modeling functions and power level settings were automated using the from the independent variable data, and uses those laptop computer, each data point required only manual modeling functions with a predicted squared error setting of angle of attack and sideslip angle using the metric to determine appropriate model structure. The joysticks in the control room. Each data point was orthogonal functions are generated in a manner that taken as the average of a ten-second dwell using a allows them to be decomposed without ambiguity into sampling rate of 100 Hz.
an expansion of ordinary multivariate polynomials.
The data for power-on subspaces was collected by This allows the identified orthogonal function model to interleaving power-on points with power-off points bc converted to a multivariate ordinary polynomial from other subspaces, in order to keep the engine expansion in the independent variables, which provides temperature at low levels for extended testing periods.
physical insight into the identified functional This was necessary to avoid damage to the electric dependencies.
motor. A regulated DC power supply was used to The next section briefly describes the multivariate power the electric motor, so that batteries would not orthogonal function modeling approach. In the Results have to be swapped in and out. The power-on section, the multivariate orthogonal function modeling subspaces were limited to relatively low angle of attack method is applied to identify response surface models and sideslip angle (see Figure 7), because of excessive fi,r non-dimensional aerodynamic force and moment vibration of the wind tunnel rig and model for powered coefficients for inference subspaces, based on data runs at high angles of attack and/or high sideslip from the FASER wind tunnel test.
angles.
The wind tunnel experiment described above was conducted over 4 weeks in May 2002, in the NASA Langley 12-foot Low-Speed Wind Tunnel.
American Institute of Aeronautics and Astronautics Multivariate Orthogonal Functions Assume an N-dimensional vector of response variable values, y = [Yl ,Y2 ..... YN ]r, modeled in terms where E is the expectation operator, and the error variance o -2 can be estimated from the residuals, of a linear combination of n modeling functions p j, j = 1,2 ..... n. Each pj is an N-dimensional vector v = y - Pti (9) which in general depends on the independent variables.
Then.
Y = al Pl +a2 P2 + .,. + an Pn +£ (3) _2_ 1 Pd)]- vWv (N-n) [(y_pd:)T (y_ (N-n) (I0) The a s , j = 1, 2 ..... n are constant model parameters to and n is the number of elements in parameter vector a.
be determined, and _ denotes the modeling error vector.
Parameter standard errors are computed as the square Eq. (3) represents the usual mathematical model used to root of the diagonal elements of the Coy(d) matrix fit a response surface to measured data from an experiment. We put aside for the moment the from Eq. (8), using 6 -2 from Eq. (10).
important questions of determining how candidate Estimated model output is modeling functions Ps should be computed from the independent variables, as well as which candidate )=P_i (11) modeling functions should be included in Eq. (3), which implicitly determines n. Now define an Nxn For response surface modeling, the modeling matrix P, functions (columns of P) are often chosen as polynomials in the measured independent variables.
P = [t_, P2 ..... Pn ] (4) This approach corresponds to using the terms of a Taylor series expansion to approximate the functional and let a =[al,a 2 ..... an] T. Eq. (3) can be written as a dependence of the output response variable on the independent variables.
standard least squares regression problem, If the modeling functions are instead multivariate y =Pa+e (5) orthogonal functions generated from the measured independent variable data, advantages accrue in the where y is a vector of measured dependent variable model structure determination for response surface values, P is a matrix whose columns contain modeling modeling. After the model structure is determined functions of the measured independent variables, and a using the multivariate orthogonal modeling functions, is a vector of unknown parameters. The variable each retained modeling function can be decomposed represents a vector of errors that are to be minimized in into an expansion of ordinary polynomials in the independent variables. Combining like terms from this a least squares sense. The goal is to determine a that final step puts the final model in the form of a Taylor minimizes the least squares cost function series expansion. It is this latter form of the model that provides the physical insight, particularly in the case of d = (y -Pa) T (y -Pa)=cTe (6) modeling non-dimensional aerodynamic force and moment coefficients. This is the reason that aircraft The parameter vector estimate that minimizes this cost dynamics and control analyses are nearly always fimction is computed from 3,5-7 conducted with the assumption of this form for the dependence of the non-dimensional aerodynamic force and moment coefficients on independent variables such it =IP T P1-1 pT y (7) as angle of attack and sideslip angle.
Ref. [3] describes a procedure for using the The estimated parameter covariance matrix is independent variable data to generate multivariate American Institute of Aeronautics and Astronautics orthogonal modeling functions pj, which have the P = __ G -1 (17) following important property: The columns of G -1 contain the coefficients for p,lpj = 0 i_j, i,j =1, 2 ..... n (12) expansion of each column of P (i.e., each multivariate orthogonal function) in terms of the ordinar 3' It is also possible to generate multivariate polynomial functions contained in the columns of --.
orthogonal functions by first generating ordinary Eq. (17) can be used to express each multivariate multivariate polynomials in the independent variables, orthogonal function in terms of ordinary multivariate then orthogonalizing these functions using polynomials.
Gram-Schmidt orthogonalization. The process begins by choosing one of the ordinary multivariate The orthogonal functions are generated in a polynomial functions as the first orthogonal function: manner that allows them to be decomposed without _lbiguity into an expansion of ordinary multivariate p_=_ (13) polynomials. The orthogonalization process can be repeated using arbitrary ordinal, multivariate where 41 is the ordinary multivariate polynomial polynomials to generate orthogonal functions of arbitrary order in the independent variables, subject function chosen to be the first orthogonal function.
only to limitations related to the information contained Then to make each subsequent ordinary multivariate in the data. For the FASER wind tunnel data response polynomial function orthogonai to the preceding surface modeling, the orthogonal modeling functions orthogonal function(s), define the f,h orthogonal were generated in the manner described above.
function pj as: If an additional independent variable is introduced j-I to represent blocking in the experiment, the orthogonal PJ=_J--EYkJPk j=2 ..... n (14) k=l function modeling can be used to orthogonalize the block effects with respect to the other independent w_iable effects. A blocking variable is typically used where 4i is the jth ordinary multivariate polynomial to indicate some change in the experimental conditions vector, and the Ykj are scalars determined from that cannot be controlled by the experimenter. The blocking variable to the first power can simply be made one of the ordinary polynomial vectors el. The p_'¢j k = 1, 2 ..... j -1 (15) orthogonalization procedure in Eqs. (13)-(15) makes Yk; - pT[ Pk j = 2 ..... n the blocking variable orthogonal to the other orthogonal modeling functions, which are generated where n is the total number of ordinary multivariate from ordinary multivariate polynomial functions. This polynomials used as raw material for generating the approach allows arbitrary blocking of the data points in multivariate orthogonal functions. Eq. (15) results the experiment using analytical means alone.
from multiplying both sides of Eq. (13) by p[, Normally, experiment designers try to arrange the invoking the mutual orthogonality of Pk, k = 1..... j, normalized settings of the independent variables so that the blocking variable is orthogonal to the other and solving for Ykj. It can be seen from independent variables and their polynomial Eqs. (13)-(15) that each orthogonal function can be combinations. This separates the block effect from the expressed in terms of ordinary polynomial functions. If oiher model terms and allows identification of a block the pj vectors and the Cj vectors are arranged as effect independent of the other terms in the model.
columns of matrices P and ---, respectively, and the However, with the analytical orthogonal blocking described above, the blocking variable is made y,j are elements in the k th row and f,h column of an orthogonal to the other modeling functions analytically, upper triangular matrix G with ones on the diagonal, so that arranging the experiment so that independent then variable settings and their polynomial combinations are orthogonal to the blocking variable is not necessary.
--=PG (16) All the models identified in this work included an awmlytic blocking variable that accounted for drift in which leads to American Institute of A eronaatics and Astronautics experimental conditions between the time when the where o-2,a., is the maximum variance of elements in central composite design points were collected and the time when the D-optimal design points were collected. the error vector tr, assuming the correct model structure, and n is the number of model terms. The With the modeling functions orthogonal, using PSE in Eq. (22) depends on the mean squared fit error Eqs. (4) and (12) in Eq. (7), the j_ element of the J/N, and a term proportional to the number of terms estimated parameter vector _i is given by in the model, n. The latter term prevents over-fitting the data with too many model terms, which is (18)
+=(,;.)/i,; ,'.)
detrimental to model prediction accuracy 9. The factor of 2 in the model over-fit penalty term accounts for the Combining Eqs. (4), (6), and ( 11 )-(12), fact that the PSE is being used when the model structure is not correct, i.e., during the model structure determination stage. Ref. [9] contains further justifying j= yT),__ .7 (19) statistical arguments and analysis for the form of PSE )=1 in Eqs. (21)-(22). Note that while the mean squared fit error J/N must decrease with the addition of each or, using Eq. (18), orthogonal modeling function to the model (by Eq. (19) or (20)), the over-fit penalty term o-2m, u n/N increases
)
(20)
T , " T =/,-,.., p, > p+ p+
with each added model term (n increases). Introducing j=l the orthogonal modeling functions into the model in order of most effective to least effective in reducing the Eq. (20) shows that when the modeling functions mean squared fit error (quantified by h_(pTs pj) for are orthogonal, the reduction in the estimated cost the fh orthogonal modeling function) means that the resulting from including the term aj pj in the model PSE metric will always have a single global minimum depends only the dependent variable data y and the value. Figure 10 depicts this graphically, using actual added orthogonal modeling function pj. This modeling results from Ref. [2]. Ref. [9] contains decouples the least squares modeling problem, and details on the statistical properties of the PSE metric, makes it possible to evaluate each orthogonal modeling including justification for its use in modeling problems.
function in terms of its ability to reduce the least For wind tunnel testing, repeated runs at the same squares model fit to the data, regardless of which other orthogonai modeling functions are present in the test conditions are often available. If o-_ is the output model. When the modeling functions p; are instead variance estimated from measurements of the output for polynomials in the independent variables (or any other repeated runs at the same test conditions, then o-2,a_ non-orthogonal function set), the least squares problem can be estimated as is not decoupled, and iterative analysis is required to find the subset of modeling functions for an adequate "_ 2 o-,7,,ax '-- 25 o-o (23) model structure.
The orthogonal modeling functions to be included If the output errors were Gaussian, Eq. (23) would in the model are chosen to minimize predicted squared correspond to conservatively placing the maximum error PSE, defined by 9 output variance at 25 times the estimated value (corresponding to a 5o- 0 maximum deviation).
However, the estimate of o-o 2 may not be very good,
ese=(y-ea) (y- ea) +2O- L
(21) N N because of relatively few repeated runs available or errors in the independent variable settings or drift errors or when duplicating test conditions for the repeat runs.
PSE = + 2 o-_,a_-- (22) The 5o- o value was found to give the most accurate N N models in model identification algorithm testing done using simulated data. In Ref. [2], the model structure American Institute of Aeronautics and Astronautics determined using PSE was found to be virtually the Results same for o',;,cL_ in the range: Figures I1, 12, and 13 show results for the response surface model fits to measured lift, drag, and < _ 2 pilching moment coefficient data obtained by applying 9 o',2 _ cr,_,a x -< 100 tyo (24) the multivariate orthogonal function nonlinear mc_deling technique to experimental data from This happens because the plot of mean square fit longitudinal subspace 16 for the FASER wind tunnel error versus added orthogonal modeling functions is test. The crosses shown in the top plot of the figures typically very fiat in the region of minimum PSE, see arc measured data, and the circles are values computed Figure 10.
fr_m the identified response surface models. Values Using orthogonal functions to model the response fo:" cr_ were found from twelve repeated runs at the variable made it possible to evaluate the merit of normalized center point of the independent variable including each modeling function individually as part ranges, using the method described above. Model of the model, using the predicted squared error PSE.
structure determination and parameter estimation were Since the goal is to select a model structure with done automatically using the orthogonal function minimum PSE, and the PSE always has a single global modeling technique. The orthogonal function modeling minimum for orthogonal modeling functions, the model software allows manual override by the analyst in the structure determination was a well-defined and model structure determination stage, but this option was straightforward process that could be (and was) not used for the results presented here. All data automated.
analysis and modeling was done on a 1.2 GHz personal After the orthogonal modeling functions that computer running MATLAB ® 6.5. The orthogonal minimized PSE were selected, each retained orthogonai modeling technique described above was implemented function was expanded into an ordinary polynomial as an m-file function.
expression, and common terms in the ordinary The model residuals shown in the middle plot of polynomials were combined using double precision Figures I 1. 12, and 13 show a random character with arithmetic to arrive finally at a multivariate model using magnitudes on the order of the noise level estimated only ordinary polynomials in the independent variables.
from repeated runs. The dashed lines in the middle and Ordinary polynomial terms that contributed less than lower plots in the figures represent +_or 0 , the square 0.1 percent of the final model root-mean-square root of the estimated output variance, computed from magnitude were dropped.
repeated runs at the subspace center point. The lower Orthogonal modeling functions are useful in plots in the figures show that the prediction residual determining the model structures for the response magnitudes are also within the estimated noise levels variables using the PSE metric, by virtue of the for the output responses. Results shown for this properties of orthogonal modeling functions and the subspace were typical.
resultant decoupling of the associated least squares Table 5 contains the identified model structures problem. The subsequent decomposition of the retained orthogonal functions is done to express the for C_, C z, and C,,, for this inference subspace.
results in physically meaningful terms and to allow Based on the data, the identification algorithm analytic differentiation for partial derivatives of the determined that no sideslip angle dependence was response variable with respect to the independent nceded for the longitudinal coefficients models in variables.
subspace 16, so the final models did not include sideslip angle.
The results shown in Figures 11, 12. and 13 and the associated analysis of residuals indicated that the filnctional dependencies were accurately captured by the response surface models identified using the experiment design and modeling techniques described.
Similar results were generated for the other longitudinal and lateral/directional subspaces listed in Tables 3 and 4 and depicted in Figures 6.7. and 8.
American Institute of A erona,,_tics and Astronautics References Figure 14 shows a screen capture of a 3D plotting tool developed to inspect the response surface fits to 1. Monzon, B.R. _'Development of a Nonlinear measured data. The symbols represent the measured Simulation for a Research Model Airplane," data points, and the smooth surface is the identified MSE Thesis, George Washington University response surface. The viewpoint for the 3D plot can be JIAFS, NASA Langley Research Center, rotated and zoomed in or out, under the control of the September 2001.
analyst, and the two independent variables for the 3D 2. Morelli, E.A. and DeLoach, R. "Response plotting can be selected from the complete set of Surface Modeling using Multivariate Orthogonal independent variables. This tool gives a good overview Functions," AIAA paper 2001-0168, 39 th of the response surface topology, using 3D slices AIAA Aerospace Sciences Meeting and Exhibit, through the inference space.
Reno, Nevada, January" 2001.
Concluding Remarks 3. Morelli, E.A., "'Global Nonlinear Aerodynamic Modeling using Multivariate Orthogonal This paper describes and demonstrates an efficient Functions," Journal of Aircraft, Vol. 32, No. 2, and effective approach to developing a wind tunnel March-April 1995, pp. 270-77.
database for a research model airplane. The use of OFAT sweeps for subspace definition proved useful in 4. DeLoach, R. "'Applications of Modem def'ming subspace boundaries and as trial runs to Experiment Design to Wind Tunnel Testing at identify and resolve problems related to experimental NASA Langley Research Center," AIAA procedure, data collection, and instrumentation. The 2000-2639, 36 'h AIAA Aerospace Sciences use of multivariate orthogonal modeling functions Meeting and Exhibit, Reno, NV, January 1998.
simplified the model structure determination problem 5. Box, G.E.P., Hunter, W.G., and Hunter, J.S., and allowed arbitrary orthogonal blocking. The quality Statistics for Experim enters - An Introduction to of the modeling and prediction results from the Design, Data Analysis, and Model Building, experiment design and modeling approach described John Wiley & Sons, Inc., New York, NY, 1978, and used were excellent. However, because the wind Chapter 14.
tunnel was not set up for MDOE experimentation, the overall efficiency of the testing was not what it could 6. Montgomery', D.C., Design and Analysis of have been.
Experiments, 5th Ed., John Wiley & Sons, New York, NY, 2001.
Significant efficiency improvement would accrue if the data could be transferred to a computer for 7. Klein, V., Batterson, J.G., and Murphy, P.C., analysis in real time, and if the wind tunnel rig could be "Determination of Airplane Model Structure set automatically to arbitrary angles of attack and from Flight Data by Using Modified Stepwise sideslip angles. In this case, the modeling software Regression," NASA TP-1916, 1981.
could automatically determine the model structure and 8. Lo, C.F., Zhao, J.L., and DeLoach, R., parameter values in real time as each data point is "Application of Neural Networks to Wind taken, then compute model fit metrics and make 3D Tunnel Data Response Surface Methods," plots. The augmentation of the MDOE experiment AIAA 2000-2639, 21 st AIAA Aerodynamic designs for higher order model fits, or the re-definition Measurement Technology and Ground Testing of subspace boundaries could be automated as well, Conference, Denver, CO, June 2000.
resulting in a testing operation where the human analyst would be needed only to provide high-level oversight 9. Barron. A.R., "Predicted Squared Error : A of the entire operation.
Criterion for Automatic Model Selection," Self-Organizing Methods in Modeling, Fariow, S.J., Ed., Marcel Dekker, Inc., New York, NY, 1984, pp. 87-104.
10. Using MATLAB, Version 6.1, The MathWorks, Inc., Natick, MA, 2001.
American Institute of Aeronautics and Astronautics Table 4 Lateral/Directional Inference Subspaces Table 1 FASER Wind Tunnel Experiment Procedures Inference Angle of Power Sideslip Procedure Description Subspace Attack Level Angle (deg) 1 Power sweep with wind off for (deg) (percent) static thrust modeling min max min max mm max 2 Angle of attack sweeps at min, max, and 7 -7.5 10 0 0 -10 10 zero elevator, for static stability and trim 8 -7.5 10 0 0 10 30 3 Angle of attack, sideslip angle, and power sweeps for inference subspace 26 -7.5 10 0 0 -30 -10 identification 19 10 20 0 0 -10 10 4 2 nd order central composite design and 20 20 40 0 0 -10 10 3rd order D-optimal design for inference subspace modeling 21 10 20 0 0 10 30 22 20 40 0 0 10 30 11 40 80 0 0 -10 10 10 40 80 0 0 10 30 Table 2 Independent Variable Ranges Variable a fl pwr 8 e 8 a 8 r 8f Min -7.5 -30 0 -25 -25 -30 0 Table 5 Longitudinal Subspace Identified Model Max 80 30 100 25 25 30 30 Model Structure Inference Units deg deg % deg deg deg deg Subspace C D = a 0 + ala 2 + a20t8 e + a38 f C L = b 0 +blae +b2_Sf +b3otcSe +b4 a'2 Table 3 Longitudinal Inference Subspaces + bsa6 _ Inference Angle of Power Sideslip Subspace Attack Level Angle (deg) C m = c O+ c18e + c2a'2 + c3_ f + c4£_¢_ e (deg) (percent) + c5a82 + c6ot min max min max min max 13 -7.5 10 0 0 -10 10 23 -7.5 10 0 100 -10 10 6 -7.5 10 0 0 10 30 15 10 20 0 0 -10 10 17 10 20 0 0 10 30 24 10 20 0 100 -10 10 16 20 40 0 0 -10 10 18 20 40 0 0 10 30 5 40 80 0 0 -10 10 4 40 80 0 0 10 30 American Institute of A eronautics and Astronautics 8 [ ,_ta t -1 o 200 4oo-- 60o 80o _doo _2oo Pul_ Count (pul_s I wc) Fi2ure 2 Static Engine Thrust, Wind Tunnel Air Off i CD 1_ o -20 0 20 40 60 80 100 .;L_ J. _ f j __ _ i -20 0 20 40 60 80 1 O0 i 0_ Cm I i -2 -20 0 20 40 60 80 100 Angle of Attack (dig) Figure 3 Angle of Attack Sweep for Subspace Boundary Definition, fl=O, pwr=0 FiLi.q_re 1 Free-flying Airplane for Sub-scale Experimental Research (FASER) in the NASA Langley 12 ft Low-Speed Wind Tunnel American Institute of A eronautics and Astronautics 05 ; r 3O SS6 SS17 SSI8 SS4 l0 -30 -20 - 10 0 10 20 30 SSl3 SSI5 SSI6 SS5 0.2, ' i -10 -30 -O -7.5 10 20 40 80 -30 -20 -10 0 10 20 30 a (deg) _ Longitudinal Subspaces, Power Off 43.02_'_1_414'_ " _ J -30 -20 - 10 0 10 20 30 3O Sideslip Angle (deg) Figure 4 Sideslip Angle Sweep for Subspace Boundary Definition, a=80, pwr=0 $523 SS24 -10 1 7 • -30 -7.5 10 20 40 80 CD O. _ - ; a (deg) Figure 7 Longitudinal Subspaces, Power On 0 _ 40 _ _ 1_ 1.5 -- , ..... _ --- CL o._ SS8 SS21 SS22 SSI0 0 20 40 60 80 100 SSII SS7 SS19 SS20 -10 SS26 Cm_.6J_y_"4_'_ _i _ ...... .......
-30 -7.5 l0 20 40 80 0 20 40 60 80 100 Power Level (percent) (deg) Fi2ure 5 Power Sweep for Subspace Boundary Lateral/Directional Subspaces, Power Off Definition, a=40, fl=0 American Institute of Aerona.tics and Astronautics _ o _ ....
1.5, . Data i i / i
i
X Central Composite Design Points i , _0± MOdel j i O D-Optimal Points j 1 + - --4_ ,® - CD i o5 i o o i !
0 J _ • 0 10 20 30 40 50 60 70 !
0 04- - _ M_el Residual i ' ' I-- NotseLevet , ' i
-0.5i O .... C _ + • ...+.."s" oO',° :'° t I
I i nn,+| 0 4_ " ; , °_ ° .0. ° l , 0 , = , 0 , -004! _ ; , 0 10 20 30 40 50 60 70
-_ -o._ _ o o.s _t
0.04 i . [ • PredictionResiduai - - Noise Le_t ' i i FiEure 9 Experiment Design Projection onto Two
L , , +--
0.02 +-+ _ +.- '+ "+ "_ ....
Dimensions = 0 , ° , I i = O+ 0 - ° ---_ - T - 4 I ' ' -o.o2- ; , + , * : I
I+ "--'+'+'"" .... +'""-- "+---'-----""
-0.04 L + . j _ _ • + 0 2 4 6 8 10 12 Run Number Fi2ure 11 Modeling and Prediction Results for Drag Coefficient, Longitudinal Subspace 16 x 10 -s Mean Square Error O_Nt Penalty , ' ! 1.4 , - - I PSE i 1.2L _ : : ' +O Mode ll j CL i_ I I i Minimum PSE l ® _, _, ! I / 0.8 _ ' _ - _ 0 10 20 30 40 50 60 70 0.1 _ O Model Residual- / ', --- Noise Le,_el I 1 2 3 4 5 6 7 8 9 0.05- . -4_+ , -- " .....
: = -&+ 0=--- . _ .
Orthogonal Function Number _,d_' 0 e", .6° re' .d_ '4r .. ,0 v_i Fiszure 10 Predicted Squared Error (PSE) Components i 0 10 20 30 40 50 60 70 0.05r-- 4 0 Prediction Residual _ + • -- Noise Le_=l O / o o o_ • ; ° ...... + - : + ....
i i -0.05 t _, 0 2 4 6 8 10 12 Run Number Fii_ure 12 Modeling and Prediction Results for Lift Coefficient, Longitudinal Subspace 16 American Institute o_A eronautics and Astronautics Data 0 * i x i i ;0 mod_ -0.6! == _'_ _ " _-'_ _ ' J'_ i o 1o 20 so 4o 5o 6o m 0.04 r • # Model Residual I - + Noise Level • O.O2
f
' • _' • ,m_-- --':" --"-%-'*---- _- -- <+-a,. ....
_+ _ _ L 2J0 L _ L o lo 30 40 so 60 70 0.04 ,'.... _+ i_redictio_P, esiclual 0.021- i . _ -.-- Noise Level , , • _ • _ , 4, - o- • * * + + i -0.02 • - J • -0D4 _ J _ - 0 2 4 6 8 10 12 Run Number Fiszure 13 Modeling and Prediction Results for Pitching Moment Coefficient. Longitudinal Subspace 16 ,2iID_U S31ZE_m_I/_I • Ui : =.Lgi ,_l F_xed In depllnder.! Variables =, F . "[ , zJ._l .._d ,_I.L + ,_:r_ ,,j ,_] .,r Lj.2 -d -04]'0"3J 'X JJ__J+.. ,_I ,,I Lj-: 'q .o.7"°a[ '_ _ I+5 _-++_. / >...,,, ,<I" 20 30 _,_.._ 1£1 Figure 14 3D Plotting Tool using Pitching Moment Cofficient for Longitudinal Subspace 16 American Institute of Aerona.,'ztics and Astronautics