Efficient Variable-Pitch Propeller Aerodynamic Model Development for Vectored-Thrust eVTOL Aircraft

Introduction

Efficient Variable-Pitch Propeller Aerodynamic Model

Development for Vectored-Thrust eVTOL Aircraft

∗ Benjamin M. Simmons NASA Langley Research Center, Hampton, Virginia 23681 This paper describes a variable-pitch propeller testing and aerodynamic modeling method- ology for electric vertical takeoff and landing (eVTOL) aircraft using wind tunnel testing.

Propellers used for eVTOL aircraft propulsion systems experience a wide range of operating conditions resulting in significant variation of axial thrust and torque, as well as off-axis forces and moments, across a typical flight envelope. An experimental design is developed to collect informative wind tunnel data in an efficient manner using response surface methods.

System identification methods are then applied to isolated propeller wind tunnel data to develop a mathematical model of the propeller valid throughout the transition envelope of a vectored-thrust eVTOL aircraft. Multiple explanatory variable definitions are postulated and compared using modeling and prediction performance metrics. Modeling results validated against data withheld from the modeling process indicate good predictive capability and agree with theoretical expectations. All identified propeller models are provided, which allows the models to be used in future eVTOL aircraft flight dynamics simulations.

Nomenclature 𝐶 , 𝐶 , 𝐶 = propeller force coefficients 𝑇 𝑇 𝑇 𝑥 𝑦 𝑧 𝐶 , 𝐶 , 𝐶 = propeller moment coefficients 𝑄 𝑄 𝑄 𝑥 𝑦 𝑧 𝐷 = propeller diameter, ft 𝑖 = propeller incidence angle, rad 𝑝 𝐽 = advance ratio 𝐽 = normal advance ratio 𝑥 𝐽 = edgewise advance ratio 𝑧 𝑛 = propeller and motor rotational speed, revolutions/s 𝑄 , 𝑄 , 𝑄 = propeller moments, ft · lbf 𝑥 𝑦 𝑧 ¯ 𝑞 = freestream dynamic pressure, lbf/ft 𝑇 , 𝑇 , 𝑇 = propeller forces, lbf 𝑥 𝑦 𝑧 𝑉 = freestream velocity, ft/s 𝛿 = collective pitch (blade root) angle, rad 𝑐 𝜂 = motor pulse width modulation command, 𝜇 s 𝑚 𝜂 = collective pitch servo-actuator pulse width modulation command, 𝜇 s 𝑐 𝜌 = air density, slug/ft I. Introduction lectric vertical takeoff and landing (eVTOL) aircraft are currently drawing significant interest in the aerospace E industry as an enabling technology for future Urban Air Mobility (UAM) transportation missions. Mission requirements for eVTOL aircraft necessitate operations in a wide variety of flight conditions spanning hover, transition, and forward flight. Consequently, eVTOL aircraft propellers experience aerodynamic conditions that significantly differ from conventional propeller and rotor operation in a typical flight profile. Propeller aerodynamics are conventionally modeled in the axial airflow condition where data tables or functional representations of axial thrust and torque coefficients are sufficient to model the propeller aerodynamics. At high incidence angles, axial thrust and torque deviate from nominal axial airflow values and off-axis propulsive forces and moments become significant.

∗ Research Engineer, Flight Dynamics Branch, MS 308, Member AIAA.

Next: Background →

Background

Several works have employed methods for theoretically and computationally predicting propeller aerodynamics at incidence [ 1 – 10 ]; however, experimental techniques provide the most accuracy in revealing the highly complex and nonlinear behavior of high incidence angle propeller aerodynamics. Previous experimental studies of propellers operating at high incidence angles include Refs. [ 11 – 19 ]. In particular, the work described in this paper builds on experimental design and modeling techniques for fixed-pitch, eVTOL aircraft propellers described in Refs. [ 16 – 19 ].

The significant new contributions of this work include characterization of a variable-pitch eVTOL aircraft propeller, an extended experimental design methodology, a comparison of modeling results using different explanatory variable definitions, and distribution of variable-pitch propeller aerodynamics models valid through a significant portion of an eVTOL aircraft flight envelope.

The objective of this work is to identify a model of the propeller aerodynamics enabling accurate flight dynamics simulation development for a new vectored-thrust eVTOL vehicle. Vectored thrust in the context of eVTOL aircraft includes both tilt-wing and tilt-rotor configurations. The intended aircraft is similar in scale and utility to the Langley Aerodrome No. 8 (LA-8) tandem tilt-wing, eVTOL, subscale wind-tunnel/flight test (SWFT) aircraft pictured in Fig. 1 [ 20 ]. Design of experiments (DOE) and response surface methodology (RSM) techniques are leveraged to develop an efficient, statistically-rigorous wind tunnel experiment [ 21 , 22 ]. Response surface models are then developed using propeller aerodynamic model structure insight [ 16 ] and system identification algorithms [ 23 ]. The proposed propeller model identification strategies can be applied for many tilt-wing and tilt-rotor eVTOL vehicles currently under ∗ development and future eVTOL vehicles. This work is complemented by recent eVTOL aircraft modeling research leveraging isolated propeller models for flight dynamics simulation development [24, 25].

Fig. 1 LA-8 mounted in the NASA Langley 12-Foot Low-Speed Tunnel.

This paper is organized as follows. Section II provides background information on propeller aerodynamics.

Section III describes the wind tunnel experimental setup. The experimental design leveraging DOE/RSM techniques is developed in Sec. IV. Section V outlines the overall propeller aerodynamic modeling approach, which is followed by an overview of the employed model identification methods given in Sec. VI. Propeller aerodynamic modeling results are presented in Sec. VII. Overall conclusions are summarized in Sec. VIII.

II. Background Propeller aerodynamics are well-defined and well-researched for nominal operating conditions in axial flow where aerodynamic predictions can be made analytically and/or experimentally. Theoretical techniques include momentum theory, blade element methods, and vortex theories [ 26 ]. Experimental techniques typically consist of developing data tables or functional representations from wind tunnel data. For a propeller operating in airflow normal to the propeller disk, only a net thrust force and a net aerodynamic torque acting along the axis of rotation are generated [ 27 ]. The individual propeller blades can be thought of as rotating wings which each produce a lift force perpendicular to the relative flow direction and a drag force parallel to the relative flow direction [ 6 ]. The summed lift forces produced by the propeller blades is the propeller thrust 𝑇 . The summed drag forces result in a net moment about the propeller shaft 𝑥 opposite to the direction of rotation, which is the propeller aerodynamic torque 𝑄 .

𝑥 ∗ Information available online at https://evtol.news/aircraft [retrieved 29 October 2021] Propeller data are generally nondimensionalized and presented in terms of thrust coefficient 𝐶 and torque coefficient 𝑇 𝑥 𝐶 (or equivalently by power coefficient 𝐶 = 2 𝜋𝐶 where 𝑃 = 2 𝜋𝑛𝑄 ). The thrust and torque coefficients are 𝑄 𝑃 𝑄 𝑥 𝑥 𝑥 defined as: 𝑇 𝑥 𝐶 = (1) 𝑇 𝑥 2 4 𝜌𝑛 𝐷 𝑄 𝑥 𝐶 = (2) 𝑄 𝑥 2 5 𝜌𝑛 𝐷 The thrust and torque coefficients can be shown through dimensional analysis to have a functional dependence on advance ratio 𝐽 , propeller blade Reynolds number 𝑅𝑒 , and propeller tip Mach number 𝑀 for a given fixed-pitch tip propeller design [27, 28]. For propellers with variable-pitch hubs, thrust and torque also vary with blade pitch angle.

Advance ratio 𝐽 , which relates to the linear distance traveled by the propeller in one revolution, is defined as: 𝑉 𝐽 = (3) 𝑛𝐷 Advance ratio generally has the largest effect on fixed-pitch propeller aerodynamics, and accordingly, thrust and torque coefficient representations are commonly expressed as only a function of advance ratio. Representing propeller aerodynamics only as a function of advance ratio for a fixed blade pitch angle requires that airflow is parallel to the propeller axis of rotation as well as the assumptions that viscous and compressibility effects are negligible [27].

Reynolds number is a dimensionless quantity which corresponds to the ratio of inertial to viscous forces acting on a body. For large propellers, the propeller blade Reynolds number effects are minimal and can generally be neglected.

For small propellers, the Reynolds number is lower, indicating that the viscous forces become important and results in propeller performance degradation [ 29 , 30 ]. Following the definition given in Ref. [ 30 ], the propeller blade Reynolds number 𝑅𝑒 is 𝜌𝑉 𝑐 𝑝 𝑅𝑒 = (4) 𝜇 where 𝑐 is the propeller chord at 75% blade length, 𝜇 is the dynamic viscosity, and 𝑉 = 0 . 75 𝜋𝑛𝐷 is the propeller blade 𝑝 linear speed at 75% blade length. Mach number is the ratio of flow speed to the speed of sound 𝑎 , which physically represents the ratio of inertial forces to forces related to compressibility of the fluid [ 31 ]. The propeller tip Mach number, which quantifies the averaged compressibility effects, is defined as [32]: 𝜋𝑛𝐷 𝑀 = (5) tip 𝑎 For the propeller studied in this paper, the propeller blade Reynolds number and tip Mach number effects were small compared to advance ratio and blade pitch angle effects; however, the influence of propeller blade Reynolds number and/or tip Mach number, which are both proportional to rotational speed 𝑛 , was found to be beneficial to consider for model development. This will be discussed further in Sec. VII.

When the airflow relative to a propeller is not parallel to the axis of rotation, the propeller will produce auxiliary forces and moments other than the axial thrust and torque [27]. In this condition, periodic variation in propeller blade local angle of attack results in a non-uniform load distribution on the propeller disk. Thus, in a general case of arbitrary flow direction relative to the propeller disk, propeller forces and moments will also be dependent on angle between the freestream velocity and propeller axis of rotation, in addition to advance ratio, propeller blade pitch angle, propeller blade Reynolds number, tip Mach number, and the propeller design. This angle between the freestream airflow and propeller rotation axis is the propeller incidence angle, 𝑖 , shown in Fig. 2. The value of 𝑖 is zero when airflow is 𝑝 𝑝 normal to the propeller disk, opposing the direction of axial thrust. Previous work has demonstrated the benefit of using the normal and edgewise (tangential) component of advance ratio 𝑉 cos 𝑖 𝑝 𝐽 = = 𝐽 cos 𝑖 (6) 𝑥 𝑝 𝑛𝐷 𝑉 sin 𝑖 𝑝 𝐽 = = 𝐽 sin 𝑖 (7) 𝑧 𝑝 𝑛𝐷 as explanatory variables for modeling propeller aerodynamics at high incidence angles [ 16 ]. Similar advance ratio component definitions are used for rotorcraft [31, 32]; this representation of advance ratio was also used in a previous propeller modeling effort which developed lookup tables for propulsive forces and moments for a quadrotor vehicle [ 33 ].

Next: Wind Tunnel Test Setup →

Wind Tunnel Test Setup

𝑥 𝑖 𝑉 𝑝 𝑧 Fig. 2 Propeller incidence angle definition and coordinate system.

The propeller side force 𝑇 , normal force 𝑇 , pitching moment 𝑄 , and yawing moment 𝑄 can be non-dimensionalized 𝑦 𝑧 𝑦 𝑧 in a manner similar to the thrust 𝑇 and torque 𝑄 [ 6 ]. The propeller normal force coefficient 𝐶 , side force coefficient 𝑥 𝑥 𝑇 𝑧 𝐶 , pitching moment coefficient 𝐶 , and yawing moment coefficients 𝐶 , are defined as: 𝑇 𝑄 𝑄 𝑦 𝑦 𝑧 𝑇 𝑄 𝑦 𝑇 𝑦 𝑄 𝑧 𝑧 𝐶 = , 𝐶 = , 𝐶 = , 𝐶 = (8) 𝑇 𝑇 𝑄 𝑄 𝑦 𝑧 𝑦 𝑧 2 4 2 4 2 5 2 5 𝜌𝑛 𝐷 𝜌𝑛 𝐷 𝜌𝑛 𝐷 𝜌𝑛 𝐷 The propeller force and moment sign convention used in this work follows the right-handed propeller coordinate system shown in Fig. 2, where the 𝑦 -axis is pointed into the page.

III. Wind Tunnel Test Setup The test article for this study was a variable-pitch, three-bladed, 19.5-inch diameter, clockwise rotating propeller.

Assuming that all clockwise and counterclockwise rotating propellers on an aircraft have perfectly mirrored designs, only one orientation needs to be tested. The propeller data and models would be identical, except that the signs of the lateral propeller force and moment components, 𝑄 , 𝑇 , and 𝑄 , would be reversed for the counterclockwise rotating 𝑥 𝑦 𝑧 propeller. Figure 3 shows the assembled propeller and variable-pitch hub on the wind tunnel sting. The propeller was powered by a 220 KV electric motor and 100 amp electronic speed control (ESC). A variable-pitch mechanism was designed and fabricated to equip the flight vehicle with high-bandwidth variable-pitch control for aircraft stabilization.

The ability to leverage variable propeller pitch as a control effector becomes important for certain full-scale eVTOL aircraft configurations because motor/propeller rotational speed control becomes slower as propeller diameter and inertia increase. As configured for testing, the custom variable-pitch mechanism allowed setting the collective pitch angle 𝛿 𝑐 (blade root pitch angle) between − 9 . 6 and + 7 . 2 deg. A multi-exposure image showing propeller blade collective pitch movement is shown in Fig. 4.

The isolated propeller wind tunnel experiment was performed in the NASA Langley Research Center 12-Foot † Low-Speed Tunnel. The test was designed to cover a significant portion of the operational flight envelope of a future eVTOL SWFT aircraft. For operational convenience and test facility implementation consideration, test factors selected to control during wind tunnel testing were freestream velocity 𝑉 (or dynamic pressure ¯ 𝑞 ), incidence angle 𝑖 , motor pulse 𝑝 width modulation (PWM) command 𝜂 , and collective pitch servo-actuator PWM command 𝜂 . The factor settings to 𝑚 𝑐 test were determined using DOE/RSM experimental design techniques that will be discussed in Sec. IV. At each test point, the six force and moment components were measured using a strain gauge balance and a propeller rotational speed measurement was provided by the ESC. A polynomial calibration model between the collective servo-actuator PWM command 𝜂 and the collective pitch angle 𝛿 was used to determine the tested collective pitch angle.

𝑐 𝑐 Each individual wind tunnel run was executed by automatically commanding incidence angle, motor PWM command, and collective PWM command at each test point throughout the run. Incidence angle variation is depicted in the multi-exposure image displayed in Fig. 5a. The propeller test coordinate system and incidence angle measurement are shown in Fig. 5b. The wind tunnel dynamic pressure setting had to be manually changed from the control room and was slower to change compared to the other test factors, but could be adjusted during the course of a run without † Information available online at https://researchdirectorate.larc.nasa.gov/12-foot-low-speed-tunnel-12-ft-lst/ [retrieved 29 October 2021] Fig. 3 Experimental propeller mounted on the wind tunnel sting.

Fig. 4 Multi-exposure image showing propeller blade pitch movement.

shutting off tunnel airflow. In accordance with the ease of changing the settings for each test factor, the incidence angle, motor PWM command, and collective servo-actuator PWM command were designated as easy-to-change (ETC) factors; freestream velocity (dynamic pressure) was designated as a hard-to-change (HTC) factor. These considerations informed the experimental design discussed in the next section.

𝑉 𝑉 𝑦 𝑧 𝑥 𝑖 𝑝 Sting Rotation (a) Multi-exposure image showing incidence angle variation (b) Propeller test coordinate system Fig. 5 Experimental propeller mounted in the NASA Langley 12-Foot Low-Speed Tunnel (viewed from above).

Each of the test factors were varied to encompass the anticipated transition envelope for the future eVTOL vehicle.

Although the test point conditions were specified using motor PWM command and collective pitch servo-actuator PWM command, these settings were used to indirectly sweep variables more pertinent to propeller aerodynamics: propeller rotational speed and collective pitch angle. Freestream velocity and propeller rotational speed were then further reduced to calculate propeller advance ratio. Changes in flow velocity drive changes in propeller advance ratio. Changes in motor PWM command change the propeller rotational speed, which has the primary effect of changing the propeller

Next: Experimental Design →

Experimental Design

blade Reynolds number and tip Mach number, but also contributes to changes in the propeller advance ratio.

IV. Experimental Design As mentioned previously, the propeller test matrix was designed using four variable factors: freestream velocity 𝑉 (HTC), incidence angle 𝑖 (ETC), motor command 𝜂 (ETC), and collective servo-actuator command 𝜂 (ETC).

𝑝 𝑚 𝑐 For the wind tunnel used for this study, changing 𝑉 required more time and effort than to change 𝑖 , 𝜂 , and 𝜂 .

𝑝 𝑚 𝑐 Consequently, for testing efficiency, 𝑉 was treated as a HTC factor held at a constant value for several consecutive test points, rather than being varied between each test point. Four different test regions were designed and tested: 1) a 2-factor Hover Region varying 𝜂 and 𝜂 with no airflow, 𝑚 𝑐 2) a 4-factor Low Incidence Region with a higher airspeed range (high-speed transition), 3) a 4-factor High Incidence Region with a lower airspeed range (low-speed transition), and 4) a 3-factor Descent Region surrounding the vortex ring state (VRS) condition at 𝑖 = 180 deg.

𝑝 ‡ The factor ranges tested for each region are listed in Table 1, with a visual overview shown in Fig. 6. These test ranges were selected based on the expected transition flight envelope, test hardware restrictions, and insight from a previous study [ 16 , 17 ], where it was expected that a cubic response surface model would be able to adequately describe the propeller force and moment coefficient variation in each region. Note that a lower collective command setting 𝜂 𝑐 corresponds to a higher collective pitch angle 𝛿 , and vice versa.

𝑐 Table 1 Factor ranges for each propeller test region Factor Units Hover Region Low Incidence Region High Incidence Region Descent Region 𝑉 (HTC) ft/s 0 (fixed) 14.5 to 71.1 14.5 to 54.3 14.5 to 45.9 𝑖 (ETC) deg N/A (fixed at 0) 0 to 65 50 to 100 180 (fixed) 𝑝 𝜂 (ETC) 𝜇 s 1500 to 1960 1500 to 1960 1500 to 1960 1500 to 1960 𝑚 𝜂 (ETC) 𝜇 s 1000 to 2000 1000 to 1500 1200 to 1600 1200 to 1600 𝑐 Fig. 6 Freestream airspeed and incidence angle range for each propeller test region.

® Modeling test matrix design was accomplished with the aid of Design-Expert , a commercially available statistical § software package. The test factor settings in each region were determined using I-Optimal response surface designs.

I-optimal designs minimize the integrated prediction variance for a predefined model order over the range of test factors, which reduces prediction error for the identified models [ 21 , 22 ]. The primary goal for this study was to develop propeller aerodynamic models with low prediction error over the region of operability, so an I-optimal design is a good response surface design choice. The Hover Region was run as a completely randomized design, meaning that each test factor was varied at each new test point. The Low Incidence Region, High Incidence Region, and Descent Region employed split-plot designs [ 34 , 35 ], or designs with restricted randomization for HTC factor(s). Split-plot designs are used to conduct efficient experiments with HTC factors, while still executing a statistically-rigorous designed experiment. The HTC factors (whole-plot factors) are held at a constant value for several consecutive test points while ‡ Propeller incidence angle 𝑖 is undefined in hover ( 𝑉 = 0 ). For Hover Region testing, the sting orientation angle was fixed at 0 deg, as indicated 𝑝 in Table 1. The Hover Region is plotted from 0 to 180 deg in Fig. 6 to reflect that hover can be approached from all incidence angles at low speeds.

§ Information available online at https://www.statease.com/software/design-expert/ [retrieved 29 October 2021] the ETC factors (sub-plot factors) are varied between each test point. The set of consecutive test points where the settings for HTC factor(s) are held constant is referred to as a group. The different amount of randomization for HTC factors and ETC factors results in different whole-plot and sub-plot variance components, which must be considered for experiment design and data analysis [ 22 ]. All split-plot designs used for this work were ordinary least-squares equivalent [ 22 , 36 – 38 ], meaning that the ordinary least-squares parameter estimates are equal to parameter estimates obtained using generalized least-squares with sub-plot/whole-plot variance component estimates. In other words, the design performance and estimates of model parameters are decoupled from the sub-plot and whole-plot variance components.

A set of completely randomized and split-plot I-optimal experiment designs were developed with different design model complexities to investigate the adequacy of fitting a cubic model. Assessment of the prediction variance of a response surface design for a given model structure provides insight into its precision of prediction and allows comparison of different response surface designs. For a completely randomized design, the variance of the predicted response is − 1 2 𝑇 𝑇 Var [ ˆ 𝑦 ( 𝒙 )] = 𝜎 𝒙 𝑿 𝑿 𝒙 (9) 0 0 where ˆ 𝑦 ( 𝒙 ) is the predicted response evaluated at the design space location 𝒙 expanded to the form of the model 0 0 structure, 𝑿 is a matrix composed of the designed test points in the form of the model structure, and 𝜎 is the measurement error variance [ 22 ]. As can be seen by examining Eq. (9) , the prediction variance is a function of the experimental design, the model structure, the location in the design space, and the measurement facility error variance.

The scaled prediction variance (SPV) and unscaled prediction variance (UPV), which remove the dependence on 𝜎 , are used to compare experimental designs prior to experimentation. SPV is defined as − 1 𝑁 Var [ ˆ 𝑦 ( 𝒙 )] 𝑇 𝑇 SPV = = 𝑁 𝒙 𝑿 𝑿 𝒙 (10) 𝜎 where the number of test points 𝑁 penalizes a larger design size [ 22 ]. The SPV considers the prediction accuracy as well as the expense of test points when comparing designs. The UPV, defined as − 1 Var [ ˆ 𝑦 ( 𝒙 )] 𝑇 𝑇 UPV = = 𝒙 𝑿 𝑿 𝒙 (11) 𝜎 provides an assessment of the prediction precision independent from the size of the experimental design. For split-plot designs, restricted randomization for HTC factors and the consequent compound error structure results in a different prediction variance computation [39]. The prediction variance for split-plot designs is calculated as: − 1 𝑇 𝑇 − 1 Var [ ˆ 𝑦 ( 𝒙 )] = 𝒙 𝑿 𝚺 𝑿 𝒙 (12) 0 0 The SPV and UPV for split-plot designs are: − 1 𝑇 𝑇 − 1 𝑁 𝒙 𝑿 𝚺 𝑿 𝒙 𝑁 Var [ ˆ 𝑦 ( 𝒙 )] 0 0 0 SPV = = (13) 2 2 2 2 𝜎 + 𝜎 𝜎 + 𝜎 sp wp sp wp − 1 𝑇 − 1 𝑇 𝒙 𝑿 𝚺 𝑿 𝒙 Var [ ˆ 𝑦 ( 𝒙 )] 0 UPV = = (14) 2 2 2 2 𝜎 + 𝜎 𝜎 + 𝜎 sp wp sp wp 2 2 Here, 𝜎 is the sub-plot variance and 𝜎 is the whole-plot variance, with the measurement error variance of a single sp wp 2 2 data point being 𝜎 + 𝜎 . For 𝑛 total whole-plots (groups), the covariance matrix 𝚺 is wp sp wp   𝚺 . . . 0  1    . .

  . . .

𝚺 = .

. .

      0 . . . 𝚺 𝑛 wp   with the covariance matrix for each 𝑗 th whole plot given as [39]:  2 2 2  𝜎 + 𝜎 . . . 𝜎   sp wp wp   . .

  . . .

𝚺 = , 𝑗 = 1 , 2 , ..., 𝑛 𝑗 . wp . .

    2 2 2   𝜎 . . . 𝜎 + 𝜎 wp sp wp   Graphical presentation of the distribution of prediction variance throughout the design space is an effective way to assess experimental designs. Each candidate design used for this work was compared using fraction of design space (FDS) plots [ 40 ]. FDS graphs depict the prediction variance distribution over the design space in a concise manner, where the prediction variance metrics are plotted against the FDS encompassing a prediction variance less than or equal to a particular value. It is also useful to consider the FDS including a particular model precision, quantified by the confidence interval half-width 𝛿 [ 41 – 43 ]. The model precision 𝛿 normalized by the response standard deviation 𝜎 plotted against FDS provides further insight into the prediction capability of the model developed from a particular experiment design, prior to conducting the experiment. For this study, a design was deemed to be adequate for fitting a particular model if 𝛿 / 𝜎 was less than two for greater than 95% of the design space. The prediction variance threshold ∗ PV used to determine the FDS within a given model precision level is 𝛿 / 𝜎 ∗ PV = (15) 𝑡 𝛼 / 2 , 𝑁 − 𝑛 𝑝 where 𝑁 is the number of test points, 𝑛 is the number of parameters in the model, and 𝛼 is the significance level 𝑝 (chosen as 𝛼 = 0 . 05 for this study).

A. Hover Region: 2-Factor Completely Randomized Design In the Hover Region, only the motor PWM command and the collective PWM command were varied with no tunnel airflow and sting incidence angle fixed. A comparison of candidate 2-factor completely randomized I-optimal designs using a cubic, quartic, quintic, and sixth-order design model complexity is shown in Table 2 and Fig. 7. Each candidate design contains five model points beyond the minimum number of test points needed to fit a full model of the specified order, as well as one center point. Table 2 lists the number of test points, as well as the FDS with 𝛿 / 𝜎 ≤ 2 for a cubic and quartic evaluation model. Due to the small number of factors, increasing the design model complexity only marginally increases the number of test points. An adequate FDS ( FDS ≥ 0 . 95 ) for a normalized model precision 𝛿 / 𝜎 ≤ 2 is obtained with a design order one power larger than the evaluation model order. Figure 7 shows the UPV, SPV, and 𝛿 / 𝜎 threshold values against FDS for a cubic evaluation model. The UPV and 𝛿 / 𝜎 threshold curves decrease in value and become more uniform (flat) as the design order increases. The SPV curve is similar for each design order.

Table 2 Comparison of candidate Hover Region experimental designs Cubic Model Quartic Model Design Order Points FDS with 𝛿 / 𝜎 ≤ 2 FDS with 𝛿 / 𝜎 ≤ 2 Cubic Design 16 0.921 0.000 Quartic Design 21 0.998 0.865 Quintic Design 27 1.000 0.998 Sixth-Order Design 34 1.000 0.999 Fig. 7 Hover Region FDS plots for a cubic evaluation model.

Although the quartic design meets the FDS criteria for a cubic evaluation model, the quintic design was selected to use for the experiment due to the increased data density and the relatively low expense of executing the additional test points for the higher order design. Prior to running the experiment, three additional center points were added and spaced throughout the test matrix (resulting in four total center points). The center points aid in further reduction of the prediction variance within the experimental region and allow estimation of pure error [ 21 ]. Six additional test points selected using a random number generator were added throughout the design to use as validation data withheld from model identification. The commanded Hover Region motor and collective PWM factor settings are shown in Fig. 8a.

The observed values of the propeller rotational speed and collective pitch angle are shown in Fig. 8b, where it can be seen that a higher collective pitch angle setting results in a lower propeller speed. This is due to increased loading on the propeller at higher collective pitch settings and the open-loop nature of the motor commands used for the experiment.

Although propeller speed could be commanded directly and regulated using a feedback controller, using motor PWM settings as the commanded test factor is more operationally convenient from an experimental design and wind tunnel implementation perspective, with minimal impact on the modeling results [24].

(a) Motor command against collective command (b) Observed rotational speed against collective pitch angle Fig. 8 Designed and observed test variable values for the Hover Region.

B. Low and High Incidence Regions: 4-Factor Split-Plot Design In the Low Incidence Region and High Incidence Region, all four test factors were varied to characterize the propeller aerodynamics across the transition flight envelope for the intended eVTOL vehicle. Each region used the same experiment design in coded units, but translated to the factor settings in engineering units shown in Table 1. A comparison of candidate 4-factor split-plot I-optimal designs using a cubic, quartic, quintic, and sixth-order design model complexity is shown in Table 3 and Fig. 9. Each candidate design contains five model points beyond the minimum number of test points and two groups beyond the minimum number of groups needed to fit a full model of the specified order. Table 3 lists the number of test points and number of groups, as well as the FDS with 𝛿 / 𝜎 ≤ 2 for a cubic and quartic evaluation model. An adequate FDS ( FDS ≥ 0 . 95 ) for a normalized model precision 𝛿 / 𝜎 ≤ 2 is obtained with a design order two powers larger than the evaluation model order. Figure 9 shows the UPV, SPV, and 𝛿 / 𝜎 threshold values against FDS for a cubic evaluation model. The UPV and 𝛿 / 𝜎 threshold curves decrease in value and become more uniform as the design order increases. The SPV curve increases in value as the design order increases.

Table 3 Comparison of candidate Low/High Incidence Region experimental designs Cubic Model Quartic Model Design Order Points Groups FDS with 𝛿 / 𝜎 ≤ 2 FDS with 𝛿 / 𝜎 ≤ 2 Cubic Design 40 6 0.000 0.000 Quartic Design 75 7 0.922 0.000 Quintic Design 131 8 0.988 0.785 Sixth-Order Design 215 9 0.999 0.982 Fig. 9 Low and High Incidence Region FDS plots for a cubic evaluation model.

Balancing design quality to support identifying a cubic model and test efficiency, the quintic split-plot design was selected for the Low Incidence and High Incidence Region experiments. The base quintic design was augmented with three center groups containing three center points each (resulting in nine total center points) near the start, middle, and end of the test matrix to allow estimation of pure error and further reduce prediction variance. The ETC factor settings were reset between center points within the center point groups by commanding randomized changes in the ETC factors, which also functioned as validation test points. Furthermore, two test points selected using a random number generator were appended to each non-center point group to use as additional validation data withheld from model identification.

These validation data were found to provide a good estimate of prediction error while remaining a modest number of test points for the present application.

The experimental test points for the Low Incidence Region are shown in Fig. 10. Figure 10a shows the commanded freestream velocity (HTC factor) against test point number. Figure 10b shows the incidence angle against freestream velocity. The commanded motor and collective factor settings are shown in Fig. 10c. The observed values of the propeller rotational speed and collective pitch angle are shown in Fig. 10d. The coded test matrix for the high incidence angle region was identical, and thus, followed the same test matrix as shown in Figs. 10a-10c with different values of test factors in engineering units.

(a) Freestream velocity (HTC factor) against test point number (b) Incidence angle against freestream velocity (c) Motor command against collective command (d) Observed rotational speed against collective pitch angle Fig. 10 Designed and observed test variable values for the Low Incidence Region.

C. Descent Region: 3-Factor Split-Plot Design In the Descent Region, the freestream velocity, motor PWM command, and collective PWM command were varied at a fixed incidence angle of 𝑖 = 180 deg to characterize the location and average thrust losses associated with the 𝑝 VRS condition. A comparison of candidate 3-factor split-plot I-optimal designs using a cubic, quartic, quintic, and sixth-order design model complexity is shown in Table 4 and Fig. 11. Each candidate design contains five model points beyond the minimum number of test points and two groups beyond the minimum number of groups needed to fit a full model of the specified order. The candidate design comparison interpretation is similar to that described in Sec. IV.B.

Table 4 Comparison of candidate Descent Region experimental designs Cubic Model Quartic Model Design Order Points Groups FDS with 𝛿 / 𝜎 ≤ 2 FDS with 𝛿 / 𝜎 ≤ 2 Cubic Design 25 6 0.000 0.000 Quartic Design 40 7 0.837 0.000 Quintic Design 61 8 0.965 0.492 Sixth-Order Design 89 9 0.997 0.956 Fig. 11 Descent Region FDS plots for a cubic evaluation model.

Due to the possibility of increased thrust sensitivity near the VRS condition, the sixth-order design was selected to run for the experiment for the increased data density throughout the experimental region. Center points and validation points were added throughout the test matrix following the same approach and justification as discussed in Sec. IV.B The commanded and observed test variable values for the Descent region are shown in Fig. 12.

(a) Freestream velocity (HTC factor) against test point number (b) Motor command against freestream velocity (c) Motor command against collective command (d) Observed rotational speed against collective pitch angle Fig. 12 Designed and observed test variable values for the Descent Region.

D. Final Design Summary The final design qualities for each propeller test region are summarized in Table 5 and Fig. 13 following the same format as the design comparisons presented in Secs. IV.A-IV.C. The validation points are included in the total number of points but are withheld from the FDS calculations. The prediction variance distribution is relatively uniform and has a low overall value for each design. The model precision is also adequate for fitting at least a full cubic model.

Table 5 Final experimental design properties for each test region Cubic Model Quartic Model Design Order Points Groups FDS with 𝛿 / 𝜎 ≤ 2 FDS with 𝛿 / 𝜎 ≤ 2 Hover Region 36 N/A 1.000 0.999 Low/High Incidence Region 162 11 0.993 0.864 Descent Region 122 12 0.997 0.963 Fig. 13 Final design FDS plots for a cubic evaluation model in each test region.

Figure 14 shows a four-dimensional representation of the observed freestream velocity, incidence angle, propeller rotational speed, and collective pitch angle for each experimental region with non-zero freestream velocity. Figure 15 shows the 𝐽 and 𝐽 coverage for each experimental region. Table 6 shows the observed ranges of several propeller test 𝑥 𝑧 variables in each experimental region. Excluding experimental setup time, execution of all test points described in this paper took approximately 3.5 hours. This test strategy leveraging DOE/RSM experimental design techniques is significantly faster compared to a similar propeller characterization study conducted using a one-factor-at-a-time (OFAT) test matrix [ 16 , 17 ]. Additionally, the randomized test point execution makes the designs more robust to time-varying systematic instrumentation errors and extraneous variables [ 21 ]. These errors are transferred to the parameter variance as opposed to corrupting the parameter estimates.

Fig. 14 Observed propeller test variable values for the Low Incidence, High Incidence, and Descent Regions.

Fig. 15 Observed normal and edgewise advance ratio values in each experimental region.

Table 6 Observed ranges of propeller test variables Variable Hover Region Low Incidence Region High Incidence Region Descent Region 𝑛 , rpm [ 2240 , 5450 ] [ 2230 , 5640 ] [ 2220 , 5350 ] [ 2290 , 5190 ] 𝛿 , deg [ − 9 . 56 , 7 . 23 ] [ − 0 . 387 , 7 . 23 ] [ − 2 . 33 , 4 . 9 ] [ − 2 . 3 , 4 . 77 ] 𝑐 𝐽 0 [ 0 . 11 , 1 . 19 ] [ 0 . 103 , 0 . 905 ] [ 0 . 123 , 0 . 661 ] 𝐽 0 [ 0 . 0476 , 0 . 697 ] [ − 0 . 157 , 0 . 438 ] [ − 0 . 661 , − 0 . 123 ] 𝑥 𝐽 0 [ 0 , 1 . 08 ] [ 0 . 0974 , 0 . 891 ] 0 𝑧 𝑅𝑒 [ 98300 , 239000 ] [ 96300 , 246000 ] [ 95300 , 231000 ] [ 98700 , 224000 ] 𝑀 [ 0 . 17 , 0 . 415 ] [ 0 . 169 , 0 . 429 ] [ 0 . 168 , 0 . 405 ] [ 0 . 174 , 0 . 393 ] tip

Next: Propeller Aerodynamic Modeling Approach →

Propeller Aerodynamic Modeling Approach

V. Propeller Aerodynamic Modeling Approach The main goal of this effort was to develop a propeller aerodynamic model relating the propeller states to output propulsive forces and moments. The response variables were the dimensionless propeller force and moment coefficients, 𝐶 , 𝐶 , 𝐶 , 𝐶 , 𝐶 , and 𝐶 , defined previously in Eqs. (1) , (2) and (8) . Response surface models developed 𝑇 𝑇 𝑇 𝑄 𝑄 𝑄 𝑥 𝑦 𝑧 𝑥 𝑦 𝑧 from designed experiments conventionally evaluate the factors under test directly as candidates for explanatory variables.

For this work, analysis was instead performed by redefining certain explanatory variables for modeling due to the unique characteristics of propeller aerodynamics at incidence and test facility integration considerations. As mentioned in Secs. III-IV, wind tunnel testing was performed by varying motor PWM command 𝜂 and collective PWM command 𝑚 𝜂 directly; however, modeling was performed using the measured propeller rotational speed 𝑛 and collective pitch 𝑐 angle 𝛿 , which are more directly related to propeller aerodynamics. Three different explanatory variable definitions 𝑐 were considered and compared for this work: ˆ ˆ ˆ 1) Explanatory Variable Definition I (EV-I): Using 𝑉 , 𝑖 , ˆ 𝑛 , and 𝛿 as explanatory variables, which is the closest 𝑝 𝑐 explanatory variable representation to the factors under test ˆ ˆ ˆ 2) Explanatory Variable Definition II (EV-II): Using 𝐽 , 𝐽 , and 𝛿 as explanatory variables, which was 𝑥 𝑧 𝑐 hypothesized to better represent propeller aerodynamics when 𝑅𝑒 and 𝑀 effects are negligible tip ˆ ˆ ˆ 3) Explanatory Variable Definition III (EV-III): Using 𝐽 , 𝐽 , ˆ 𝑛 , and 𝛿 as explanatory variables, which was 𝑥 𝑧 𝑐 hypothesized to better represent propeller aerodynamics when 𝑅𝑒 and/or 𝑀 effects are significant.

tip The hat ( ˆ ) notation is used to indicate that the variables are centered on the median value within the modeling data. It is important to perform modeling with explanatory variables expressed in coded units, or in engineering units with the explanatory variables centered on a reference value, to maintain low correlation among candidate regressors. Note that the variables not directly varied in the Hover and Descent Regions are omitted for their respective modeling analyses.

The utility of each of these three parameterizations are compared in Sec. VII.

Since the transformations between commanded test factors and the explanatory variables used for modeling were not linear, data collinearity assessment was performed to assess the impact on model identification. Data collinearity is defined as a correlation between regressors high enough to cause corrupted model identification [ 23 ]. Data collinearity will cause difficulty in both model structure determination and parameter estimation because the effects of certain regressors on the response cannot be distinguished. Collinearity assessment is useful for confirming that a choice of modeling candidate regressors from a given experiment design or data set are sufficiently decorrelated for model identification. Correlation between two candidate regressors can be assessed using the pairwise correlation coefficient 𝑇 ¯ ¯ ( 𝝃 − 𝜉 ) ( 𝝃 − 𝜉 ) 𝑖 𝑗 𝑖 𝑗 𝑟 = √︃ (16) 𝑖 𝑗 √︁ 𝑇 𝑇 ¯ ¯ ¯ ¯ ( 𝝃 − 𝜉 ) ( 𝝃 − 𝜉 ) ( 𝝃 − 𝜉 ) ( 𝝃 − 𝜉 ) 𝑖 𝑖 𝑗 𝑗 𝑖 𝑖 𝑗 𝑗 ¯ ¯ where 𝝃 and 𝝃 are two regressor measurement histories, with means denoted 𝜉 and 𝜉 , respectively. A correlation 𝑖 𝑗 𝑖 𝑗 coefficient value of zero means the signals are uncorrelated, or orthogonal, and an absolute correlation coefficient of one indicates that the signals are completely correlated. A correlation coefficient between regressors with magnitude greater than 0.9 indicates that data collinearity problems may be encountered [23].

Figure 16 shows the pairwise correlation coefficient values between all candidate regressors sorted in ascending order within a full cubic model structure in the Low Incidence Region. The candidate regressors are assembled using the coded test factors as well as the three explanatory variable definitions centered on each respective median value. The highest correlation for each case is between the linear and cubic candidate regressors for each respective explanatory variable, which are challenging to decorrelate for any modeling problem. As would be expected, the overall lowest correlation is observed for the candidate regressors assembled from the test factors, or the designed test matrix. Low candidate regressor correlation, however, is still sufficiently maintained for a large majority of the candidate regressors for each of the explanatory variable definitions (EV-I, EV-II, and EV-III). This supports the validity of using the transformed explanatory variable definitions for model identification. Additionally, the model structure selection method used for this work, discussed next, automatically avoids adding highly correlated regressors into the model, thus, providing additional protection against data collinearity [23].

Next: Model Identification Methodology →

Model Identification Methodology

Fig. 16 Pairwise correlation between candidate regressors for a cubic model in the Low Incidence Region.

VI. Model Identification Methodology As highlighted previously, most of the data used for model development in this paper were obtained from split-plot experiments [ 34 , 35 ], or experiments with restriction to randomization for certain test factor(s). This results in a compound error structure with different whole-plot and sub-plot variance estimates. Generally, restricted maximum likelihood (REML) analysis is used to estimate the variance components and generalized least-squares is used for parameter estimation when developing models for split-plot experiments [ 44 ]. However, because the split-plot experimental designs used for this work were ordinary least-squares equivalent, ordinary least-squares regression can be used for parameter estimation and provides the same solution as generalized least-squares regression, independent of the variance component estimates [ 22 , 36 – 38 ]. Model structure determination, however, is still complicated by the compound error structure [ 44 ]. Statistical software packages generally use REML analysis to estimate variance components and significance testing leveraging methods described in Ref. [45] to refine the model structure.

For this work, an alternative model structure determination approach was sought that does not rely on significance testing or the variance component estimation. An approach meeting these requirements allows more model development flexibility for when the explanatory variables used for modeling are transformed versions of the designed test factors, which was necessary for this study for the reasons explained in Sec. V. For example, for explanatory variable definitions EV-II and EV-III, the whole-plot and sub-plot factors are combined to calculate the advance ratio components which precludes using conventional model structure development techniques used for split-plot experiments. This approach represents an example of an engineering solution used to develop an improved, aerodynamically-informed model taking priority over statistical convention. Although it is not a statistically-based choice, using explanatory variable definitions different from the test factors was necessary to accommodate the practical considerations associated with this experimental study and yield the most useful mathematical model for the intended engineering application.

Multivariate orthogonal function (MOF) modeling, described in Refs. [ 23 , 46 ], meets the aforementioned model structure identification requirements—it does not rely on significance testing or variance component estimation to determine the model structure. MOF modeling instead selects regressors to include in the model based on their independent ability to improve characterization of the response variable. MOF modeling has also been successfully used in previous modeling studies using wind tunnel testing to characterize fixed-wing aircraft [ 47 , 48 ] and propeller aerodynamics at incidence [ 16 ]. Consequently, MOF modeling was selected as the model structure determination technique for this work.

The MOF modeling approach [ 23 , 46 ] is initiated by orthogonalizing a predefined set of candidate regressors using an algorithm such as Gram-Schmidt orthogonalization or QR decomposition. Orthogonal regressors are convenient for model structure development because of the ability to independently assess the potential of the orthogonalized candidate regressors to model the response variable—this facilitates selecting only the model terms that significantly contribute to model effectiveness. For model selection, the orthogonal regressors are ranked from highest to lowest decrease in the mean squared fit error (MSFE), 𝑇 MSFE = ( 𝒛 − ˆ 𝒚 ) ( 𝒛 − ˆ 𝒚 ) (17) 𝑁 where ˆ 𝒚 is the length 𝑁 model response vector and 𝒛 is the length 𝑁 measured response vector. In other words, this ranks the regressors from highest to lowest ability to improve the model. Candidate regressors are brought into the model structure in this order. Note that this procedure prevents data collinearity because if there is high correlation between candidate regressors, after the first candidate regressor is orthogonalized, any other highly correlated candidate regressors will be close to zero after passing through the orthogonalization process, which will prevent the latter model terms from being included in the model structure [23].

Deciding which terms to include in the final model can be done using one or more statistical metrics. A common threshold for MOF modeling is to minimize the predicted squared error (PSE) [ 23 , 49 ]. The PSE is the sum of the MSFE for a model and a complexity penalty term related to the number of terms included in the model 𝑛 𝑝 PSE = MSFE + 2 𝜎 (18) max 𝑁 where 𝑛 is the number of terms in the current model structure and 𝜎 is an estimate of the upper-bound of mean 𝑝 max squared error for the model prediction of data not used to develop the model. The model complexity penalty term contains a factor of 2 because the PSE metric is being employed for model structure identification, where inadequate model forms are being evaluated [ 47 , 48 ]. In wind tunnel testing, the quantity 𝜎 can be estimated using the variance max of measured response between repeated data points or from the variance between the measured response 𝒛 and mean measured response ¯ 𝒛 . The quantity 𝜎 is estimated for this work using the measurement error variance between max repeated data points ˆ 𝜎 as [23, 47, 48]: 2 2 𝜎 = 25 ˆ 𝜎 (19) max Another statistical metric that has been used as a stopping criterion for MOF modeling is the coefficient of 2 2 determination 𝑅 [50, 51]. The 𝑅 metric, calculated as 𝑇 2 ˆ 𝒚 𝒛 − 𝑁 ¯ 𝒛 𝑅 = (20) 𝑇 2 𝒛 𝒛 − 𝑁 ¯ 𝒛 quantifies the amount of variation in the response variable about its mean value that is described by the model. Because 𝑅 always increases with the addition of more model terms, it is important that each model term added on the basis 2 2 of the 𝑅 metric significantly increases its value. A common 𝑅 increase constituting a significant increase with the addition of a new model term is Δ 𝑅 = 0 . 5% [ 23 ]. This means that the model term describes a minimum of 0 . 5% of the total variation about the mean response.

Both PSE and 𝑅 were used as a cutoff threshold for candidate model terms to include in the final model structure.

After the orthogonal regressors were ranked by their ability to reduce the MSFE, the cutoff for model term addition was chosen to be either the candidate model term that minimized the PSE or the last term to increase 𝑅 by 0 . 5% , whichever admitted more terms into the model. This selection was made because in certain circumstances, PSE minimization was found to admit too few terms into the model due to having a rough estimated value of 𝜎 .

max After determining the model terms to include in the model structure, the final parameter values were estimated using ordinary least-squares regression in ordinary regressor space as: − 1 𝑇 𝑇 ˆ 𝜽 = 𝑿 𝑿 𝑿 𝒛 (21) ˆ Here, 𝜽 is a vector of 𝑛 estimated parameters, 𝑿 is a 𝑁 × 𝑛 matrix consisting of column vectors of regressors, and 𝒛 is 𝑝 𝑝 the measured response variable. The software implementing the MOF modeling and least-squares regression algorithm ¶ was from the System IDentification Programs for AirCraft (SIDPAC) software toolbox.

After identifying the model structure and parameter estimates, model adequacy was examined using data withheld from the model development process. Regression methods minimize the summation of squared modeling residuals between modeled and measured response, so inspection of modeling fit metrics and modeling residuals alone does not provide information about the model predictive capability. Assessment of model performance using validation data not used for modeling provides a more reliable estimate of model prediction accuracy. Validation assessment can be performed by analyzing the residuals between the measured response 𝒛 and predicted response ˆ 𝒚 for the same explanatory variable inputs. Comparison of modeling and prediction residuals is useful because a significant increase in the spread of prediction residuals compared to modeling residuals is a way of diagnosing an improper model fit.

For this work, modeling and prediction residuals are compared using the normalized root-mean-square modeling error (NRMSE) metric: √︄ 𝑇 1 ( 𝒛 − ˆ 𝒚 ) ( 𝒛 − ˆ 𝒚 ) NRMSE = (22) range ( 𝒛 ) 𝑁 ¶ Information available online at https://software.nasa.gov/software/LAR-16100-1 [retrieved 29 October 2021]

Next: Modeling Results →

Modeling Results

Range normalization, as opposed to other normalization metrics such as the mean or maximum absolute value of the response variable, provides the fairest comparison between prediction error metrics for this work because axial propeller responses are generally biased above or below zero and off-axis propeller responses are centered about zero [16].

A prediction error metric defined using critical binomial analysis of validation residuals is also useful as a quantitative measure of the model adequacy. For this analysis, each validation data point is considered to either pass or fail relative to a prediction error threshold. Failed trials can indicate model inadequacy or measurement error. The binomial test provides a threshold to determine when the number of failures is statistically significant. For this metric, the prediction error level associated with the number of successful trials being equal to the critical binomial number quantifies the model prediction capability. The process of computing critical binomial analysis of residuals prediction error metric ∗ ( 𝑒 ) is shown in Ref. [ 24 ] and further explanation of critical binomial analysis of residuals and justification for using 𝑐𝑣 this metric to assess prediction error are given in Ref. [52].

VII. Modeling Results Separate propeller aerodynamic models were identified for each test region using each of the three explanatory ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ variable definitions [EV-I ( 𝑉 , 𝑖 , ˆ 𝑛 , 𝛿 ), EV-II ( 𝐽 , 𝐽 , 𝛿 ), and EV-III ( 𝐽 , 𝐽 , ˆ 𝑛 , 𝛿 )] postulated in Sec. V. To analyze 𝑝 𝑐 𝑥 𝑧 𝑐 𝑥 𝑧 𝑐 the adherence to regression modeling assumptions and model fit adequacy, sample residual plots for 𝐶 in the Low 𝑇 𝑥 Incidence Region are shown in Fig. 17 for the EV-I model and Fig. 18 for the EV-III model. Figures 17a and 18a (a) Residual history (b) Normal probability plot (c) Residuals against predicted response Fig. 17 Externally studentized residual diagnostics for the EV-I 𝑪 model in the Low Incidence Region.

𝑻 𝒙 (a) Residual history (b) Normal probability plot (c) Residuals against predicted response Fig. 18 Externally studentized residual diagnostics for the EV-III 𝑪 model in the Low Incidence Region.

𝑻 𝒙 show the externally studentized modeling residuals, 𝑡 , against data point number [ 22 , 53 ]. Figures 17b and 18b show 𝑖 normal probability plots for the externally studentized modeling residuals. Figures 17c and 18c show the externally studentized modeling residuals against the predicted response. These plots show that the residuals for each model definition are reasonably independent, normally distributed, and have constant variance, satisfying the regression modeling assumptions and suggesting that the identified model parameters are meaningful. There is a slight trend observed in Figs. 17a and 18a, but this was deemed to be small enough as to not significantly impact the modeling results—parameter estimation results are robust to time-varying instrumentation errors due to the randomization in the experimental design. Residuals with similar character were obtained for the other response variables, test regions, and explanatory variable definitions.

Figures 19-23 compare modeling and prediction error metrics for each modeling region and explanatory variable definition. Tables 7-11 at the end of this paper list the numerical values shown in the figures. These metrics are used to assess the quality of the modeling results and compare the utility of the different explanatory variable combinations.

Fig. 19 Coefficient of determination, 𝑹 , for each local model.

Fig. 20 Number of model parameters identified for each local model.

Fig. 21 Modeling data NRMSE for each local model response.

Figure 19 (Table 7) shows the coefficient of determination ( 𝑅 ) value for each local model. The thrust coefficient 𝐶 models have an 𝑅 value above 98% indicating that a large majority of the variation in the response is described by 𝑇 𝑥 the models. The 𝐶 models for EV-III generally have the highest 𝑅 values, although the difference in value is small 𝑇 𝑥 compared to the EV-I and EV-II models. For the other propeller force and moment coefficients, it is observed that EV-II most often has the lowest 𝑅 value compared to EV-I and EV-III. This suggests that Reynolds number and/or tip Mach Fig. 22 Validation data NRMSE for each local model response.

∗ Fig. 23 Binomial analysis of residuals prediction error metric, 𝒆 , for each local model.

𝒄𝒗 number effects become significant for characterizing the variation of the torque coefficient 𝐶 and the off-axis force 𝑄 𝑥 and moment coefficients. The 𝐶 and 𝐶 models for EV-I and EV-III have 𝑅 values above 97%. The 𝐶 and 𝐶 𝑇 𝑄 𝑄 𝑄 𝑧 𝑥 𝑦 𝑧 2 2 models for EV-I and EV-III have 𝑅 values above 82%, with models in the High Incidence Region having 𝑅 values above 94%. The models for 𝐶 have the lowest 𝑅 values, but this is expected since 𝐶 is the weakest response for 𝑇 𝑇 𝑦 𝑦 propellers at incidence [16].

Figure 20 (Table 8) shows the number of model terms identified in each local model. Generally, the models corresponding to EV-I have the largest number of model parameters for the stronger response variables compared to respective EV-II and EV-III models. The models for EV-II are seen to generally have the fewest number of model terms, but as explained above, they omit variables that characterize propeller blade Reynolds number and tip Mach number effects leading to a worse model fit. This suggests that propeller blade Reynolds number effects or tip Mach number effects are important to consider for the propeller used in this study.

Figures 21-22 (Tables 9-10) show the NRMSE for modeling data ( NRMSE ) and validation data ( NRMSE ), 𝑚 𝑣 respectively. Although 𝑅 and NRMSE reflect the model fit quality, model prediction metrics calculated using 𝑚 validation data not used in the modeling process are generally considered a more reliable indicator of modeling success.

The NRMSE values calculated using modeling and validation data for each respective local model are similar and low-valued indicating that a quality model has been identified. The binomial analysis of residuals prediction error metric ∗ ( 𝑒 ) value for each model is shown in Figure 23 (Table 11). In the Hover, Low Incidence, and High Incidence regions, 𝑐𝑣 the prediction error for thrust coefficient is less than 3%; in the Descent region prediction error for thrust coefficient is around 5% or less. The models for 𝐶 and 𝐶 in all respective regions, as well as 𝐶 and 𝐶 in the High 𝑇 𝑄 𝑄 𝑄 𝑧 𝑥 𝑦 𝑧 Incidence region, have prediction error values of approximately 10% or less. This indicates that high-quality propeller aerodynamic models have been developed. The 𝐶 models and the 𝐶 , 𝐶 models in the Low Incidence region have 𝑇 𝑄 𝑄 𝑦 𝑦 𝑧 prediction error values of roughly 20% or less, which are reasonable values because these models characterize relatively weak responses with a lower signal-to-noise ratio.

After investigating the three different explanatory variable definitions, one definition needed to be selected for the final model. Since the EV-II models were suspected to lack the parameterization to best describe the aerodynamics for the propeller described in this work (i.e., EV-II models cannot characterize propeller blade Reynolds number or tip Mach number effects), EV-II was not selected as the final model explanatory variable type. Comparing EV-I and ˆ ˆ ˆ EV-III in general, and especially for stronger propeller force and moment coefficients, using EV-III ( 𝐽 , 𝐽 , ˆ 𝑛 , and 𝛿 as 𝑥 𝑧 𝑐 explanatory variables) results in a better modeling fit and lower prediction error, while also having a lower number ˆ ˆ of modeling terms. Additionally, using dimensionless propeller variables 𝐽 , 𝐽 , as opposed to 𝑉 and 𝑖 extends the 𝑥 𝑧 𝑝 ˆ ˆ ˆ extrapolation capabilities of the model. For these reasons, the models developed with 𝐽 , 𝐽 , ˆ 𝑛 , and 𝛿 as explanatory 𝑥 𝑧 𝑐 ˆ variables (EV-III) were selected as the final models. The EV-III model parameter estimates 𝜃 and parameter standard ˆ errors 𝑠 ( 𝜃 ) for each propeller force and moment coefficient in each test region are shown in Tables 12-27. The model parameters are presented in the order in which they were added to the model, so terms appearing first are most significant to the model. As an example of how the polynomial models would appear for usage, the 𝐶 model in the Low Incidence 𝑇 𝑥 Region in Table 14 is expressed in polynomial equation form as: ˆ ˆ ˆ ˆ ˆ ˆ 𝐶 = 𝐶 + 𝐶 𝐽 + 𝐶 𝛿 + 𝐶 𝐽 + 𝐶 𝐽 + 𝐶 𝐽 𝛿 𝑇 𝑇 𝑇 𝑥 𝑇 𝑐 𝑇 𝑇 𝑧 𝑇 𝑥 𝑐 𝑥 𝑥𝑜 𝑥 𝑥 𝑥 𝑥 𝑥 𝑥 ˆ ˆ 2 ˆ ˆ ˆ 𝐽𝑥 𝛿𝑐 ˆ 𝐽𝑧 𝐽𝑥 𝛿𝑐 𝐽 𝑥 Recall that the explanatory variables are all centered on a reference value in the model equations (see Sec. V). The median values used to center each explanatory variable are given in Table 28. As an example for 𝐽 , the centered 𝑥 explanatory variable is ˆ 𝐽 = 𝐽 − 𝐽 𝑥 𝑥 𝑥 𝑜 where the value of 𝐽 , as well as the respective reference values for the other explanatory variables, are listed in 𝑥 𝑜 Table 28.

Figures 24-26 show the modeled response for axial thrust coefficient 𝐶 in each test region compared to modeling 𝑇 𝑥 and validation wind tunnel data. For the Hover Region (Fig. 24), the 𝐶 model and data are plotted against rotational 𝑇 𝑥 Fig. 24 Hover Region 𝑪 response surface model Fig. 25 Low and High Incidence Region 𝑪 models 𝑻 𝑻 𝒙 𝒙 compared to measured data. compared to measured data.

speed 𝑛 and collective pitch angle 𝛿 , which were the only two explanatory variables varied in the experiment. For the 𝑐 Low and High Incidence Regions (Fig. 25), the 𝐶 model is plotted against normal advance ratio 𝐽 and edgewise 𝑇 𝑥 𝑥 advance ratio 𝐽 at a collective pitch angle of 𝛿 = 4 deg and a propeller rotational speed of 𝑛 = 3500 rpm. The modeling 𝑧 𝑐 and validation test data shown in the figure have a collective pitch angle between 3 deg and 5 deg (i.e., near the nominal value used to display the modeled response). All tested propeller rotational speed values are shown due to the relatively weak influence of propeller speed on the 𝐶 response compared to the other explanatory variables. As can be seen in 𝑇 𝑥 Table 14 and Table 20, the 𝐶 model in the Low Incidence Region is not dependent on 𝑛 and the 𝐶 model in the High 𝑇 𝑇 𝑥 𝑥 Incidence Region only has one relatively weak model term dependent on 𝑛 . The 𝑛 effects, attributed to propeller blade Reynolds number and/or tip Mach effects, were more influential in other propeller force and moment coefficients. Small differences seen between the model and data are partially attributed to the collective pitch angle and propeller rotational speed settings for the displayed test data not perfectly matching the displayed modeled response surface value. For the Descent Region (Fig. 25), the 𝐶 model and data are plotted against normal advance ratio 𝐽 and collective pitch angle 𝑇 𝑥 𝑥

Next: Conclusions →

Conclusions

Fig. 26 Descent Region 𝑪 response surface model Fig. 27 Modeled 𝑪 response in the Descent Region 𝑻 𝑻 𝒙 𝒙 compared to measured data. at different collective pitch angle settings.

of 𝛿 . As can be seen in Table 26, the 𝐶 model in the Descent Region is independent of direct 𝑛 effects (propeller 𝑐 𝑇 𝑥 blade Reynolds number and/or tip Mach effects). The modeled response in each region shows overall good agreement with the modeling and validation data points indicating that the models are suitable for their intended application of flight dynamics simulation development. Note that models identified over different ranges of certain variables require blending methods to eliminate discontinuities for simulation, where it is desirable to have continuous, differentiable transition between modeling regions. Reference [ 16 ] presents an approach to smoothly blend propeller models for this purpose, in addition to other practical aspects of using the identified propeller models in a flight dynamics simulation.

Figure 27 shows the modeled 𝐶 response in the Descent Region against normal advance ratio at several different 𝑇 𝑥 collective pitch angle settings shown in the legend. The minimum 𝐶 for each collective pitch angle on the plot is the 𝑇 𝑥 approximate location of the VRS at each collective pitch angle setting. It can be seen that the VRS location occurs at a more negative 𝐽 value (greater total advance ratio 𝐽 ) as the collective pitch angle increases.

𝑥

VIII. Conclusions

A method for variable-pitch propeller aerodynamic model development suitable for use in flight dynamics simulations for vectored-thrust eVTOL aircraft has been presented. Wind tunnel data were collected using designed experiments allowing identification of accurate cubic response surface propeller aerodynamic models in four regions of flight: hover, high-incidence/low-speed transition, low-incidence/high-speed transition, and descent. Application of I-optimal completely randomized and split-plot experiment designs enabled efficient collection of high-quality wind tunnel data for model identification. The propeller model identification approach leveraged multivariate orthogonal function modeling to identity local polynomial models, which was a suitable and effective choice for the selected experimental designs and modeling problem definitions. Three different explanatory variable formulations were compared and the results indicated that the variables derived from propeller aerodynamics theory provided the best models. The final models describe the variation of the dimensionless force and moment coefficients in each test region as a function of normal advance ratio, edgewise advance ratio, propeller rotational speed, and collective pitch angle. Collectively, the models identified in the four different test regions describe propeller aerodynamics over a wide range of flight conditions seen in operational flight for eVTOL vehicles. Model validation assessment indicated that the models are high quality and sufficient for the purpose of supporting flight dynamics model development for a future eVTOL aircraft. Because eVTOL vehicle dynamics are highly dependent on propulsive effects, accurate propeller aerodynamic modeling is essential for high-fidelity simulator development. This paper demonstrated several novel propeller modeling techniques useful for modeling future eVTOL vehicles and provides progress in this new area of aerodynamic modeling research.

Acknowledgments

This research was funded by the NASA Aeronautics Research Mission Directorate (ARMD) Transformational Tools and Technologies (TTT) project. Wind tunnel testing support was provided by Gregory Howland, Ronald Busan, Wes O’Neal, George Altamirano, Rose Weinstein, and Clinton Duncan. Wes O’Neal and David Hatke added a new capability to the 12-Foot Low-Speed Tunnel software allowing dynamic pressure to be changed in an automated data collection run. The automatic propeller variable-pitch mechanism was developed by Gregory Howland. Photography support was provided by Lee Pollard. Additional team members contributing to the wind tunnel testing effort included Earl Harris, Stephen Riddick, Sue Grafton, Jacob Cook, David North, Brian Duvall, Matthew Gray, Steven Geuther, Jason Welstead, and Siena Whiteside.

Modeling and Validation Metric Tables

ˆ ˆ

𝑥 𝑐

ˆ ˆ

Descent Region EV-III: 𝐽 , ˆ 𝑛 , 𝛿 5.02 — — 3.51 — —

𝑥 𝑐

EV-III Parameter Estimate Tables

Table 12 Parameter estimates for 𝑪 in the Hover Table 14 Parameter estimates for 𝑪 in the Low Inci- 𝑻 𝑻 𝒙 𝒙 2 2

Region ( 𝑹 = 99 . 97% ) dence Region ( 𝑹 = 99 . 40% )

ˆ ˆ ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 ) Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 7.054E − 02 2.643E − 04 𝐶 + 5.102E − 02 4.433E − 04

𝑇 𝑇 𝑥𝑜 𝑥𝑜

𝐶 + 4.877E − 01 5.187E − 03 𝐶 − 2.317E − 01 1.832E − 03

𝑇 𝑇 𝑥 𝑥 ˆ ˆ 𝛿𝑐 𝐽𝑥

𝐶 + 2.844E − 01 2.273E − 02 𝐶 + 5.686E − 01 6.237E − 03

𝑇 𝑇 𝑥 𝑥 2 ˆ ˆ 𝛿𝑐 𝛿 𝑐

𝐶 − 2.301E − 01 1.242E − 02

𝑇

𝐶 − 1.696E + 00 2.918E − 01 𝑥

𝑇 𝑥 2 ˆ 3 𝐽 ˆ 𝛿 𝑥 𝑐

𝐶 + 1.927E − 02 1.414E − 03

𝑇 𝑥 ˆ 𝐽𝑧

𝐶 + 3.828E − 01 4.092E − 02

𝑇 𝑥 ˆ ˆ 𝐽𝑥 𝛿𝑐

Table 13 Parameter estimates for 𝑪 in the Hover

𝑸 𝒙

Region ( 𝑹 = 99 . 91% )

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 − 3.706E − 03 2.478E − 05

𝑄 𝑥𝑜

𝐶 − 2.757E − 02 1.761E − 04

𝑄 𝑥 ˆ 𝛿𝑐

𝐶 − 8.434E − 02 1.872E − 03

𝑄 𝑥 ˆ 𝛿 𝑐

𝐶 + 6.931E − 06 1.232E − 06

𝑄 𝑥 ˆ 𝑛 Table 15 Parameter estimates for 𝑪 in the Low Inci- Table 17 Parameter estimates for 𝑪 in the Low In- 𝑻 𝑸 𝒚 𝒙

cidence Region ( 𝑹 = 98 . 12% )

dence Region ( 𝑹 = 86 . 13% )

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 − 5.958E − 03 4.606E − 05

𝐶 − 1.894E − 03 1.572E − 04 𝑄

𝑇 𝑥𝑜 𝑦𝑜

𝐶 − 3.964E − 02 6.993E − 04

𝐶 − 8.849E − 03 6.870E − 04 𝑄

𝑥 𝑇 𝑦 ˆ 𝛿𝑐 ˆ 𝐽𝑧

𝐶 + 1.064E − 02 2.020E − 04

𝐶 + 8.063E − 03 5.899E − 04 𝑄

𝑇 𝑥 𝑦 ˆ ˆ 𝐽𝑥 𝐽𝑥

𝐶 + 4.634E − 05 2.593E − 06

𝐶 + 8.497E − 05 1.053E − 05 𝑄

𝑥 𝑇 𝑦 ˆ 𝑛 ˆ 𝑛

𝐶 + 2.048E − 02 1.218E − 03

𝑄

𝐶 − 4.461E − 02 8.483E − 03

𝑇 𝑥 𝑦 2 ˆ ˆ ˆ 𝐽 𝐽𝑧 𝛿𝑐 𝑥

𝐶 + 3.002E − 04 4.551E − 05

𝐶 + 2.554E − 02 4.590E − 03

𝑇 𝑄 𝑦 𝑥 ˆ ˆ ˆ 𝛿𝑐 ˆ 𝑛 𝐽𝑥 𝐽𝑧

𝐶 + 9.493E − 05 1.538E − 05

𝐶 − 6.945E − 02 1.630E − 02

𝑄 𝑇 𝑥 𝑦 2 ˆ ˆ ˆ 𝐽 ˆ 𝑛 𝐽𝑥 𝐽 𝑧 𝑧

𝐶 − 1.378E + 00 4.614E − 01

𝑇 𝑦 ˆ 𝛿 𝑐

𝐶 − 6.985E − 03 3.081E − 03

𝑇 𝑦 ˆ 𝛿 ˆ 𝑛 𝑐

𝐶 + 1.253E − 01 2.599E − 02

𝑇 𝑦 ˆ ˆ 𝐽 𝐽𝑧 𝑥

𝐶 − 3.445E − 04 9.993E − 05

𝑇 𝑦 ˆ 𝐽 ˆ 𝑛 𝑧

𝐶 − 1.376E − 06 4.495E − 07

𝑇 𝑦 ˆ 𝑛

𝐶 + 7.399E − 04 3.290E − 04

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝐽𝑧 ˆ 𝑛

𝐶 + 3.368E − 02 1.154E − 02

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝛿𝑐

Table 18 Parameter estimates for 𝑪 in the Low In-

𝑸 𝒚

𝐶 − 1.032E − 02 3.554E − 03

𝑇 𝑦 ˆ 𝐽 𝑥

cidence Region ( 𝑹 = 84 . 92% )

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 − 2.664E − 03 1.221E − 03

𝑄 𝑦 ˆ 𝐽𝑥

𝐶 + 9.627E − 05 3.513E − 05

𝑄 𝑦 ˆ 𝐽𝑧 ˆ 𝑛

𝐶 + 5.446E − 05 7.954E − 06

𝑄 𝑦 ˆ 𝑛

𝐶 − 7.585E − 04 1.221E − 04

𝑄 𝑦𝑜

𝐶 + 6.799E − 03 2.169E − 03

𝑄 𝑦 ˆ 𝛿𝑐

𝐶 − 1.915E − 02 3.680E − 03

𝑄 𝑦 ˆ ˆ 𝐽𝑥 𝐽𝑧

Table 16 Parameter estimates for 𝑪 in the Low Inci-

𝑻 𝒛

𝐶 + 1.586E − 02 3.095E − 03

𝑄

dence Region ( 𝑹 = 98 . 27% ) 𝑦

ˆ 𝐽 𝑥

𝐶 + 5.184E − 05 3.617E − 05

𝑄 𝑦 ˆ 𝐽𝑥 ˆ 𝑛

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 7.998E − 03 2.447E − 03

𝑄 𝑦 ˆ 𝐽 𝑧

𝐶 − 3.929E − 02 1.094E − 03

𝑇 𝑧 ˆ 𝐽𝑧

𝐶 − 7.827E − 06 1.701E − 06

𝑄 𝑦 ˆ 2 𝐽𝑧 ˆ 𝑛

𝐶 − 6.643E − 03 2.001E − 04

𝑇 𝑧𝑜

𝐶 − 1.540E − 06 3.586E − 07

𝑄 𝑦

𝐶 − 1.883E − 01 2.088E − 02 2

𝑇 𝑧 ˆ 𝑛 ˆ ˆ 𝐽𝑧 𝛿𝑐

𝐶 + 5.352E − 02 1.550E − 02

𝑄 𝑦

𝐶 + 1.137E − 02 9.038E − 04 ˆ ˆ

𝑇 𝐽 𝛿𝑐 𝑧 𝑧 ˆ 𝐽𝑥

𝐶 − 9.327E − 04 2.494E − 04

𝑄

𝐶 − 2.922E − 02 3.387E − 03 𝑦

𝑇 ˆ 2 𝑧 𝐽 ˆ 𝑛 ˆ 𝑥 𝛿𝑐

𝐶 − 1.389E − 05 7.449E − 06

𝐶 + 7.735E − 05 1.025E − 05

𝑄 𝑇 𝑦 𝑧 ˆ 2 ˆ 𝑛 𝛿𝑐 ˆ 𝑛

𝐶 − 1.551E − 02 3.114E − 03

𝑇 𝐶 − 8.027E − 06 2.750E − 06

𝑧 𝑄 2 𝑦 ˆ 2 𝐽 ˆ 𝐽𝑥 ˆ 𝑛 𝑧

𝐶 − 2.945E − 06 6.773E − 07 𝐶 − 8.403E − 04 2.783E − 04

𝑄 𝑇 𝑦 𝑧 ˆ ˆ 2 𝐽𝑥 𝐽𝑧 ˆ 𝑛 ˆ 𝑛

𝐶 − 9.038E − 02 1.668E − 02

𝐶 − 4.838E − 02 1.733E − 02

𝑇 𝑄 𝑧 𝑦 2 3 ˆ ˆ ˆ 𝐽𝑥 𝐽 𝐽 𝑧 𝑥

𝐶 + 1.919E − 01 4.802E − 02 𝐶 − 3.084E − 02 1.477E − 02

𝑇 𝑄 𝑧 𝑦 2 2 ˆ ˆ ˆ ˆ 𝐽 𝛿𝑐 𝐽𝑥 𝐽 𝑧 𝑧 Table 19 Parameter estimates for 𝑪 in the Low In- Table 21 Parameter estimates for 𝑪 in the High 𝑸 𝑻 𝒛 𝒚

cidence Region ( 𝑹 = 83 . 26% )

Incidence Region ( 𝑹 = 82 . 68% )

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 − 1.266E − 02 1.142E − 03

𝑄 𝐶 − 3.768E − 03 9.358E − 05

𝑧 𝑇 𝑦𝑜 ˆ 𝐽𝑧

𝐶 − 2.000E − 02 1.746E − 03

𝐶 − 9.130E − 06 4.602E − 06

𝑇 𝑄 𝑦 𝑧 ˆ ˆ 2 𝛿𝑐 𝐽𝑥 ˆ 𝑛

𝐶 − 3.965E − 01 6.011E − 02

𝐶 − 9.283E − 05 1.064E − 05 𝑇

𝑦 𝑄 𝑧 ˆ ˆ ˆ ˆ 𝑛 𝐽𝑥 𝐽𝑧 𝛿𝑐

𝐶 + 2.317E − 06 5.943E − 07 𝐶 − 1.692E − 03 5.718E − 04

𝑄 𝑧 𝑇 𝑦 ˆ 𝑛 ˆ 𝐽𝑧

𝐶 − 1.583E − 01 2.820E − 02

𝑄 𝐶 − 6.493E − 06 3.799E − 05

𝑧 𝑇 2 𝑦 ˆ ˆ ˆ 𝐽 𝛿𝑐 𝑧 𝐽𝑧 ˆ 𝑛

𝐶 − 8.719E − 02 1.859E − 02 𝐶 − 9.258E − 04 7.534E − 04

𝑄 𝑇 𝑧 𝑦 3 ˆ ˆ ˆ 𝐽𝑧 𝛿𝑐 ˆ 𝑛 𝐽 𝑥

𝐶 − 6.420E − 02 7.830E − 03

𝐶 + 1.611E − 02 4.730E − 03 𝑇

𝑦 𝑄 𝑧 3 2 ˆ ˆ 𝐽 𝐽 𝑧 𝑥

𝐶 + 2.734E − 02 3.522E − 03

𝐶 − 8.568E − 03 2.738E − 03

𝑇 𝑄 𝑦 𝑧 ˆ 2 ˆ 2 𝐽 𝐽 𝑧 𝑧

𝐶 + 2.527E − 01 5.708E − 02

𝐶 − 1.241E − 03 4.331E − 04 𝑇

𝑄 𝑦 𝑧 ˆ ˆ ˆ ˆ 𝐽𝑥 𝐽𝑧 ˆ 𝑛 𝐽 𝛿𝑐 𝑧

𝐶 + 1.803E − 02 3.241E − 03

𝑇 𝑦 ˆ 2 𝐽 𝑥

𝐶 − 3.029E − 02 4.665E − 03

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝐽𝑧

𝐶 + 1.321E − 01 2.068E − 02

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝐽 𝑧

𝐶 + 1.166E − 03 3.400E − 04

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝐽𝑧 ˆ 𝑛

𝐶 − 3.342E − 02 8.854E − 03

𝑇 𝑦 ˆ ˆ 𝐽𝑧 𝛿𝑐

𝐶 + 2.014E − 05 5.943E − 06

𝑇 𝑦 ˆ 𝑛

𝐶 + 1.519E − 01 6.317E − 02

𝑇 𝑦 ˆ ˆ 𝐽 𝛿𝑐 𝑥

𝐶 + 6.441E − 02 2.982E − 02

𝑇 𝑦 ˆ 𝛿 𝑐

𝐶 + 2.395E − 03 7.666E − 04

𝑇 𝑦 ˆ 𝐽𝑥

𝐶 − 6.251E − 01 2.434E − 01

𝑇 𝑦 ˆ ˆ 𝐽𝑥 𝛿 𝑐

Table 20 Parameter estimates for 𝑪 in the High

𝑻 𝒙

𝐶 − 4.123E − 03 2.078E − 03

2 𝑇 𝑦

Incidence Region ( 𝑹 = 99 . 87% ) ˆ

𝛿 ˆ 𝑛 𝑐

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 8.638E − 02 2.227E − 04

𝑇 𝑥𝑜

𝐶 − 1.667E − 01 1.936E − 03

𝑇 𝑥 ˆ 𝐽𝑥

𝐶 + 4.844E − 01 2.783E − 03

𝑇 𝑥 ˆ 𝛿𝑐

𝐶 + 6.300E − 02 1.394E − 03

𝑇 𝑥 ˆ 𝐽𝑧

𝐶 − 3.751E − 01 1.330E − 02

𝑇 𝑥 ˆ ˆ

𝐽𝑥 𝐽𝑧 Table 22 Parameter estimates for 𝑪 in the High In-

𝑻 𝒛

𝐶 + 1.544E − 01 6.887E − 03

cidence Region ( 𝑹 = 97 . 76% )

𝑇 𝑥 ˆ 𝐽 𝑧

𝐶 + 3.234E − 01 3.279E − 02

𝑇 𝑥 ˆ ˆ

𝐽𝑥 𝐽 ˆ ˆ

𝑧 Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 3.205E − 01 2.174E − 02

𝑇 𝑥

ˆ ˆ 𝐶 − 1.132E − 02 8.805E − 05

𝐽𝑥 𝛿𝑐 𝑇 𝑧𝑜

𝐶 − 2.401E − 01 1.832E − 02

𝑇

𝑥 𝐶 − 3.037E − 02 7.763E − 04

3 𝑇 ˆ 𝑧 𝐽 ˆ 𝑧 𝐽𝑧

𝐶 − 4.661E − 01 5.525E − 02

𝐶 − 3.567E − 02 1.866E − 03

𝑇 𝑥 𝑇 𝑧 ˆ ˆ 𝐽 𝛿𝑐 𝑥

𝐶 + 2.554E − 01 5.030E − 02 𝐶 − 5.633E − 02 4.439E − 03

𝑇 𝑇 𝑧 𝑥 ˆ ˆ 2 𝐽𝑥 𝐽𝑧 ˆ ˆ 𝐽 𝐽𝑧 𝑥

𝐶 − 5.915E − 02 6.191E − 03

𝐶 − 2.828E − 06 6.887E − 07 𝑇

𝑧 𝑇 𝑥 3 ˆ 2 𝐽 ˆ 𝑛 𝑧

𝐶 + 6.817E − 02 1.876E − 02

𝐶 − 4.666E − 03 8.625E − 04

𝑇 𝑇 𝑥 𝑧 ˆ ˆ ˆ 𝐽𝑧 𝛿𝑐 𝐽𝑥 Table 23 Parameter estimates for 𝑪 in the High Table 25 Parameter estimates for 𝑪 in the High 𝑸 𝑸 𝒙 𝒛 2 2

Incidence Region ( 𝑹 = 98 . 09% ) Incidence Region ( 𝑹 = 96 . 39% )

ˆ ˆ ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 ) Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 − 4.968E − 03 4.436E − 05 𝐶 − 1.915E − 02 5.020E − 04

𝑄 𝑄 𝑥𝑜 𝑧 ˆ 𝐽𝑧

𝐶 − 3.383E − 02 7.032E − 04

𝑄

𝑥 𝐶 − 2.770E − 03 9.572E − 05

ˆ 𝑄 𝛿𝑐 𝑧𝑜

𝐶 + 3.540E − 04 1.044E − 04

𝐶 − 5.641E − 02 4.267E − 02

𝑄 𝑄 𝑥 𝑧 ˆ 2 2 𝐽 ˆ 𝑛 ˆ ˆ 𝐽 𝛿𝑐 𝑧 𝑧

𝐶 + 1.605E − 02 1.064E − 03 𝐶 − 2.216E − 02 1.776E − 03

𝑄 𝑥 𝑄 𝑧 2 ˆ ˆ 𝐽 𝛿𝑐 𝑥

ˆ 3 𝑻 𝒙 𝐽 𝑧

Region ( 𝑹 = 98 . 43% )

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 8.906E − 02 6.646E − 04

𝑇 𝑥𝑜

𝐶 + 1.126E + 00 1.815E − 01

𝑇 𝑥 ˆ 𝐽 𝑥

Table 24 Parameter estimates for 𝑪 in the High

𝑸 𝒚

𝐶 + 4.380E − 01 1.109E − 02

2 𝑇 𝑥 ˆ

Incidence Region ( 𝑹 = 95 . 70% ) 𝛿𝑐

𝐶 − 1.445E − 01 6.936E − 03

𝑇 𝑥 ˆ 𝐽𝑥

ˆ ˆ

Parameter 𝜃 𝑠 ( 𝜃 )

𝐶 + 8.904E − 01 4.206E − 02

𝑇 𝑥 ˆ 𝐽 𝑥

𝐶 − 2.780E − 02 7.570E − 04

𝑄 𝑦 ˆ

𝐶 + 6.459E − 01 8.246E − 02

𝐽𝑥 𝑇 𝑥 ˆ ˆ 𝐽𝑥 𝛿𝑐

𝐶 + 2.879E − 03 9.625E − 05

𝑄 𝑦𝑜

𝐶 + 1.403E − 02 7.421E − 04

𝑄 𝑦 ˆ 𝐽𝑧

Table 27 Parameter estimates for 𝑪 in the Descent

𝐶 + 5.117E − 04 1.395E − 04

𝑄 𝑸 𝑦 𝒙 ˆ 𝐽 ˆ 𝑛 𝑧

Region ( 𝑹 = 98 . 18% )

𝐶 + 8.761E − 05 6.092E − 06

𝑄 𝑦 ˆ 𝑛

𝐶 + 2.827E − 04 3.718E − 05

𝑄

𝑦 ˆ ˆ

ˆ Parameter 𝜃 𝑠 ( 𝜃 )

𝐽𝑧 ˆ 𝑛

𝐶 + 3.307E − 02 4.561E − 03

𝑄

𝑦 𝐶 − 4.858E − 03 3.650E − 05

ˆ 2 𝑄 𝑥𝑜 𝐽 𝑥

𝐶 + 6.980E − 03 1.543E − 03 𝐶 − 3.438E − 02 5.012E − 04

𝑄 𝑄 𝑥 𝑦 ˆ ˆ 𝛿𝑐 𝛿𝑐

𝐶 + 1.424E − 02 8.473E − 04

𝐶 + 1.328E − 01 1.812E − 02 𝑄

𝑄 𝑥 𝑦 2 ˆ ˆ ˆ 𝐽 𝐽𝑥 𝐽 𝑥 𝑧

𝐶 − 1.066E − 01 1.466E − 02

𝐶 − 4.272E − 02 6.805E − 03 𝑄

𝑄 𝑥 𝑦 2 ˆ ˆ ˆ 𝛿 𝐽𝑥 𝐽𝑧 𝑐 Table 28 Median reference values used to center explanatory variables for each modeling region

Variable Hover Region Low Incidence Region High Incidence Region Descent Region

𝐽 — + 3.320E − 01 + 8.333E − 02 − 3.038E − 01

𝑥

𝐽 — + 1.693E − 01 + 3.215E − 01 —

𝑧

𝑛 , rev/s + 6.249E + 01 + 6.345E + 01 + 6.155E + 01 + 6.144E + 01

𝛿 , rad − 6.747E − 03 + 6.992E − 02 + 2.676E − 02 + 2.550E − 02

𝑐

References

[1] Ribner, H. S., “Propellers in Yaw,” NACA TR-820, 1945.

[2] Ribner, H. S., “Formulas for Propellers in Yaw and Charts of the Side-Force Derivative,” NACA TR-819, 1945.

[3] Crigler, J. L., and Gilman Jr., J., “Calculation of Aerodynamic Forces on a Propeller in Pitch or Yaw,” NACA TN-2585, Jan.

1952.

[4] De Young, J., “Propeller at High Incidence,” Journal of Aircraft , Vol. 2, No. 3, 1965, pp. 241–250. https://doi.org/10.2514/3.

43646.

[5] Smith, H. R., “Engineering Models of Aircraft Propellers at Incidence,” Ph.D. thesis, University of Glasgow, Jan. 2015.

nd [6] Phillips, W. F., Mechanics of Flight , 2 ed., John Wiley & Sons, Hoboken, NJ, 2010. Chap. 2.

[7] McCormick, B. W., Aljabri, A. S., Jumper, S. J., and Martinovic, Z. N., “The Analysis of Propellers Including Interaction Effects,” NASA CR-158111, 1979.

[8] Theys, B., Dimitriadis, G., Hendrick, P., and De Schutter, J., “Experimental and Numerical Study of Micro-Aerial-Vehicle Propeller Performance in Oblique Flow,” Journal of Aircraft , Vol. 54, No. 3, 2017, pp. 1076–1084. https://doi.org/10.2514/1.

C033618.

[9] Leng, Y., Yoo, H., Jardin, T., Bronz, M., and Moschetta, J.-M., “Aerodynamic Modeling of Propeller Forces and Moments at High Angle of Incidence,” AIAA SciTech 2019 Forum , AIAA Paper 2019-1332, Jan. 2019. https://doi.org/10.2514/6.2019-1332.

[10] Fei, X., Litherland, B. L., and German, B. J., “Development of an Unsteady Vortex Lattice Method to Model Propellers at Incidence,” AIAA Journal , Vol. 60, No. 1, 2022, pp. 176–188. https://doi.org/10.2514/1.J060133.

[11] McLemore, H. C., and Cannon, M. D., “Aerodynamic Investigation of a Four-Blade Propeller Operating through an Angle-of- ◦ ◦ Attack Range from 0 to 180 ,” NACA TN-3228, Jun. 1954.

◦ [12] Yaggy, P. F., and Rogallo, V. L., “A Wind-Tunnel Investigation of Three Propellers Through an Angle-of-Attack Range from 0 ◦ to 85 ,” NASA TN D-318, May 1960.

[13] Kuhn, R. E., and Draper, J. W., “Investigation of the Aerodynamic Characteristics of a Model Wing-Propeller Combination and ◦ of the Wing and Propeller Separately at Angles of Attack up to 90 ,” NACA TR-1263, 1956.

[14] Theys, B., Dimitriadis, G., Andrianne, T., Hendrick, P., and De Schutter, J., “Wind Tunnel Testing of a VTOL MAV Propeller in Tilted Operating Mode,” 2014 International Conference on Unmanned Aircraft Systems (ICUAS) , May 2014, pp. 1064–1072.

https://doi.org/10.1109/ICUAS.2014.6842358.

[15] Leng, Y., Jardin, T., Bronz, M., and Moschetta, J.-M., “Experimental Analysis of Propeller Forces and Moments at High Angle of Incidence,” AIAA SciTech 2019 Forum , AIAA Paper 2019-1331, Jan. 2019. https://doi.org/10.2514/6.2019-1331.

[16] Simmons, B. M., “System Identification for Propellers at High Incidence Angles,” Journal of Aircraft , Vol. 58, No. 6, 2021, pp.

1336–1350. https://doi.org/10.2514/1.C036329.

[17] Simmons, B. M., and Hatke, D. B., “Investigation of High Incidence Angle Propeller Aerodynamics for Subscale eVTOL Aircraft,” NASA TM-20210014010, May 2021.

[18] Stratton, M., and Landman, D., “Wind Tunnel Test and Empirical Modeling of Tilt-Rotor Performance for eVTOL Applications,” AIAA SciTech 2021 Forum , AIAA Paper 2021-0834, Jan. 2021. https://doi.org/10.2514/6.2021-0834.

[19] Stratton, M. C., “Empirical Modeling of Tilt-Rotor Aerodynamic Performance,” Master’s thesis, Old Dominion University, 2021. https://doi.org/10.25777/mh0a-e343.

[20] North, D. D., Howland, G., and Busan, R. C., “Design and Fabrication of the Langley Aerodrome No. 8 Distributed Electric Propulsion VTOL Testbed,” AIAA SciTech 2021 Forum , AIAA Paper 2021-1188, Jan. 2021. https://doi.org/10.2514/6.2021-1188.

th [21] Montgomery, D. C., Design And Analysis of Experiments , 8 ed., John Wiley & Sons, Inc., Hoboken, NJ, 2013.

[22] Myers, R. H., Montgomery, D. C., and Anderson-Cook, C. M., Response Surface Methodology: Process and Product th Optimization Using Designed Experiments , 4 ed., John Wiley & Sons, Hoboken, NJ, 2016.

nd [23] Morelli, E. A., and Klein, V., Aircraft System Identification: Theory and Practice , 2 ed., Sunflyte Enterprises, Williamsburg, VA, 2016.

[24] Simmons, B. M., and Murphy, P. C., “Aero-Propulsive Modeling for Tilt-Wing, Distributed Propulsion Aircraft Using Wind Tunnel Data,” Journal of Aircraft , published online 2 March 2022. https://doi.org/10.2514/1.C036351.

[25] Simmons, B. M., Buning, P. G., and Murphy, P. C., “Full-Envelope Aero-Propulsive Model Identification for Lift+Cruise Aircraft Using Computational Experiments,” AIAA AVIATION 2021 Forum , AIAA Paper 2021-3170, Aug. 2021. https: //doi.org/10.2514/6.2021-3170.

nd [26] McCormick, B. W., Aerodynamics, Aeronautics, and Flight Mechanics , 2 ed., John Wiley & Sons, New York, NY, 1995.

Chap. 6.

[27] Dommasch, D. O., Elements of Propeller and Helicopter Aerodynamics , Pitman Publishing Corporation, 1953. Chaps. 1-3.

[28] Greitzer, E. M., Spakovszky, Z. S., and Waitz, I. A., Thermodynamics and Propulsion , Lecture Notes, Massachusetts Institute of Technology, Cambridge, MA, 2006. Chap. 11, http://web.mit.edu/16.unified/www/FALL/thermodynamics/notes/node86.html, [retrieved March 2022].

[29] Brandt, J., and Selig, M., “Propeller Performance Data at Low Reynolds Numbers,” 49th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition , AIAA Paper 2011-1255, Jan. 2011. https://doi.org/10.2514/6.2011- 1255.

[30] Deters, R. W., Ananda, G. K., and Selig, M. S., “Reynolds Number Effects on the Performance of Small-Scale Propellers,” 32nd AIAA Applied Aerodynamics Conference , AIAA Paper 2014-2151, Jun. 2014. https://doi.org/10.2514/6.2014-2151.

nd [31] Leishman, J. G., Principles of Helicopter Aerodynamics , 2 ed., Cambridge University Press, Cambridge, England, U.K., 2016. Chaps. 2, 7.

[32] Johnson, W., Rotorcraft Aeromechanics , Cambridge University Press, Cambridge, England, U.K., 2013. Chaps. 2, 6.

[33] Foster, J. V., and Hartman, D., “High-Fidelity Multi-Rotor Unmanned Aircraft System (UAS) Simulation Development for Trajectory Prediction Under Off-Nominal Flight Dynamics,” 17th AIAA Aviation Technology, Integration, and Operations Conference , AIAA Paper 2017-3271, Jun. 2017. https://doi.org/10.2514/6.2017-3271.

[34] Jones, B., and Nachtsheim, C. J., “Split-Plot Designs: What, Why, and How,” Journal of Quality Technology , Vol. 41, No. 4, 2009, pp. 340–361. https://doi.org/10.1080/00224065.2009.11917790.

[35] Cortes, L. A., Simpson, J. R., and Parker, P. A., “Response surface split-plot designs: A literature review,” Quality and Reliability Engineering International , Vol. 34, No. 7, 2018, pp. 1374–1389. https://doi.org/10.1002/qre.2350.

[36] Vining, G. G., Kowalski, S. M., and Montgomery, D. C., “Response Surface Designs within a Split-Plot Structure,” Journal of Quality Technology , Vol. 37, No. 2, 2005, pp. 115–129. https://doi.org/10.1080/00224065.2005.11980310.

[37] Parker, P. A., Kowalski, S. M., and Vining, G. G., “Classes of Split-plot Response Surface Designs for Equivalent Estimation,” Quality and Reliability Engineering International , Vol. 22, No. 3, 2006, pp. 291–305. https://doi.org/10.1002/qre.771.

[38] Parker, P. A., Kowalski, S. M., and Vining, G. G., “Construction of Balanced Equivalent Estimation Second-Order Split-Plot Designs,” Technometrics , Vol. 49, No. 1, 2007, pp. 56–65. https://doi.org/10.1198/004017006000000462.

[39] Liang, L., Anderson-Cook, C. M., and Robinson, T. J., “Fraction of Design Space Plots for Split-plot Designs,” Quality and Reliability Engineering International , Vol. 22, No. 3, 2006, pp. 275–289. https://doi.org/10.1002/qre.770.

[40] Zahran, A., Anderson-Cook, C. M., and Myers, R. H., “Fraction of Design Space to Assess Prediction Capability of Response Surface Designs,” Journal of Quality Technology , Vol. 35, No. 4, 2003, pp. 377–386. https://doi.org/10.1080/00224065.2003.

11980235.

[41] Anderson, M. J., and Whitcomb, P. J., RSM Simplified: Optimizing Processes Using Response Surface Methods for Design of nd Experiments , 2 ed., CRC Press, Boca Raton, FL, 2017. https://doi.org/10.1201/9781315382326.

[42] Anderson, M. J., Adams, W. F., and Whitcomb, P. J., “How to Properly Size Response Surface Method Experiment (RSM) Designs for System Optimization,” 2016. URL https://cdn.statease.com/media/public/documents/how_to_properly_size_ experiments_aimed_at_system_optimization.pdf, [retrieved April 2022].

[43] Whitcomb, P., “FDS—A Power Tool for Designers of Optimization Experiments,” Stat-Teaser Newsletter , Stat-Ease, Inc., Sep.

2008, pp. 1–3. URL https://cdn.statease.com/media/public/documents/statteaser-0908.pdf, [retrieved April 2022].

[44] Goos, P., and Jones, B., Optimal Design of Experiments , John Wiley & Sons, West Sussex, United Kingdom, 2011. Chaps. 7, 10.

[45] Kenward, M. G., and Roger, J. H., “Small Sample Inference for Fixed Effects from Restricted Maximum Likelihood,” Biometrics , Vol. 53, No. 3, 1997, pp. 983–997. https://doi.org/10.2307/2533558.

[46] Morelli, E. A., “Global Nonlinear Aerodynamic Modeling Using Multivariate Orthogonal Functions,” Journal of Aircraft , Vol. 32, No. 2, 1995, pp. 270–277. https://doi.org/10.2514/3.46712.

[47] Morelli, E. A., and DeLoach, R., “Response Surface Modeling using Multivariate Orthogonal Functions,” 39th AIAA Aerospace Sciences Meeting and Exhibit , AIAA Paper 2001-0168, Jan. 2001. https://doi.org/10.2514/6.2001-168.

[48] Morelli, E. A., and DeLoach, R., “Wind Tunnel Database Development Using Modern Experiment Design and Multivariate Orthogonal Functions,” 41st Aerospace Sciences Meeting and Exhibit , AIAA Paper 2003-653, Jan. 2003. https://doi.org/10.

2514/6.2003-653.

[49] Barron, A. R., “Predicted Squared Error: A Criterion for Automatic Model Selection,” Self-Organizing Methods in Modeling , Farlow, S. J., Ed., Marcel Dekker, Inc., New York, NY, 1984, pp. 87–104.

[50] Morelli, E. A., “Real-Time Global Nonlinear Aerodynamic Modeling for Learn-To-Fly,” AIAA Atmospheric Flight Mechanics Conference , AIAA Paper 2016-2010, Jan. 2016. https://doi.org/10.2514/6.2016-2010.

[51] Morelli, E. A., “Practical Aspects of Real-Time Modeling for the Learn-To-Fly Concept,” 2018 Atmospheric Flight Mechanics Conference , AIAA Paper 2018-3309, Jun. 2018. https://doi.org/10.2514/6.2018-3309.

[52] DeLoach, R., “Assessment of Response Surface Models Using Independent Confirmation Point Analysis,” 48th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition , AIAA Paper 2010-741, Jan. 2010.

https://doi.org/10.2514/6.2010-741.

th [53] Montgomery, D. C., Peck, E. A., and Vining, G. G., Introduction to Linear Regression Analysis , 5 ed., John Wiley & Sons, Hoboken, New Jersey, 2012. Chap. 4.

Overview

Introduction View in PDF