Development of a Flexible Framework Multi-Design Optimization Scheme for a Hand Launched Fuel Cell-Powered UAV

: This paper presents different methods for the design of a hand-launchable, ﬁxed wing, fuel cell-powered unmanned aerial vehicle (UAV) to maximize ﬂight endurance during steady level ﬂight missions. The proposed design methods include the development of physical models for different propulsion system components. The performance characteristics of the aircraft are modeled through empirical contributing analyses in which each analysis corresponds to an aircraft subsystem. The contributing analyses are collected to form a design structure matrix which is included into a multi-disciplinary analysis to solve for the design variables over a deﬁned design space. The optimal solution is found using a comprehensive optimization tool developed for long endurance ﬂight missions. Optimization results showed a signiﬁcant improvement in UAV ﬂight endurance that reached up to 475 min with take-off ratio equals to 59 min/kg. Wind tunnel and bench-top tests and HiL simulation tests are performed to validate the results obtained from the optimization tools. Validated optimization results showed an increase of the overall UAV ﬂight endurance by 19.4% compared to classical approaches in design methods.


Introduction
The UAV and fuel cell industries are seeking to enhance the capability and performance of miniature unmanned aerial vehicles (MUAVs) powered by fuel cells as propulsion systems [1][2][3]. MUAVs are primarily used for surveillance, target engagement and tracking missions [4]. Such missions require the MUAV to possess long endurance performance characteristics capable of maximizing the UAV flight time [5]. Hence, there will be more fuel and a larger UAV propulsion system size in terms of engines, fuel storage and other components. Similarly, electric propulsion systems that utilize batteries as the main power source also require more batteries to increase the endurance of the UAV, and hence, increase UAV weight. This leads to an increase in UAV size and thus more drag force. At a certain point, the UAV will reach a threshold in which adding more fuel or batteries will increase the total weight and size, and in turn, will decrease the endurance. Moreover, such a solution is not effective, as it will increase both the initial and the running costs. Second, there is a need for new power sources that lower the dependence on fossil fuel and achieves more sustainability. Third, an environmentally friendly power source such as solar cells or fuel cells will be considered as a better alternative to internal combustion engines [6,7].
Hand launched unmanned aerial vehicles are considered to be a good choice for small scale UAVs designed for long endurance flight missions as they have many advantages over other types of MUAVs [8][9][10]. First, hand launched UAVs do not require runways or mechanical launching aids. Second, they require less training and less setup time. Third, they can have longer endurance because they require less power for take-off during flight missions. Figure 1 shows a sample of the most popular MUAVs in which endurance versus method of the HiL simulation tests. The new developed MDO can be successfully implemented by other researches and be used in MATLAB optimization toolbox easily.
To meet this paper objectives, the paper is organized as follows: Section 2 describes the technical approach. Section 3 includes the UAV conceptual design which defines the UAV baseline configuration. Section 4 includes the preliminary UAV design which includes the development of the multi-disciplinary design analysis (optimization of inner loop modeling). Section 5 includes the detail design which includes the development of the multi-disciplinary design optimization (optimization of outer loop modeling). Section 6 includes the results of the analysis and the validation of results through wind-tunnel tests and hardware in the loop simulation tests.

Technical Approach
The design process is divided into three main parts; UAV conceptual design level, UAV preliminary design level and then a detail design level. In the UAV conceptual design level, the design requirements and constraints are specified. After that, the design space exploration process defines a number of alternatives through a matrix of alternatives that meets the design requirements. Then, the base line aircraft configuration is analyzed through a qualitative and quantitative study for the matrix of alternatives. Based on that, design variables and design matrices are defined.
Aircraft design consists of subsystems which are described as multi-disciplinary complex systems by nature such as fuel cell and propeller subsystems. This requires simplifications and decomposition schemes to provide a well-structured method for analyzing the MDA which are done in the preliminary design stage. Figure 2 illustrates the process of obtaining a logical MDA for aircraft design. First, different design variables and Deployment of PEMFC hand launched UAVs emerges some challenges that can lead to better performance characteristics. There is a contradiction between the takeoff mass requirements and long endurance flight missions which results in high take-off/endurance ratio such as in [11,12]. Increasing the endurance means increasing the energy storage. Increasing the energy storage results in the increase of the takeoff mass. However, the takeoff mass is limited to the ability of the user. There is a need for a well-developed flexible optimization scheme to give the optimal solution among the tradeoff between flight endurance and takeoff mass.
Most MUAVs shown in Figure 1 use batteries as the main power supply [13,14]. PEM fuel cells have proposed a very good solution due their high energy densities that can reach up to five times compared to lithium ion batteries [15,16]. Talking about PEM fuel cell includes the study and analysis of fuel cell subsystems as well. PEMFC is a very sensitive and costly device that should be analyzed correctly. PEMFC models still need further validation and analysis when integrated in MUAV propulsion systems because they tend to show difficulties with temperature change and transient response [17][18][19][20].
A high-level conceptual design study of small-scale long-endurance aircraft was performed in [20]. It was found that fuel cell model that depends on polarization curve resulted in accurate fuel cell voltage and current predictions. Lindhal and his team introduced a mathematical derivation, simulation and validation of a simple propulsion system for UAV at steady state conditions [13]. The propulsion system components used were a solid oxide fuel cell (SOFC) as the power source, motor controller, brushless motor and propeller. The results validation showed a good agreement between simulation and wind tunnel test results. A team from Colorado State University succeeded in designing, building and flying a fuel cell powered UAV [11]. The goal of the study was to achieve 24 h endurance. The team was able to achieve up to 28 h of flight but the aircraft had a high take-off/endurance ratio. Wagner et al. in [14] have developed an in-depth mathematical formulation for a hand launched UAV to achieve long endurance. The optimization technique used to solve the DSM was a genetic algorithm. Results showed that a notable improvement of the UAV endurance from 2 h to almost 6 h. However, the main power source was a lithium polymer battery and not a fuel cell. Moreover, a basic fuel cell powered UAV was carried by a team from Sydney University [17]. It included a preliminary design with a demonstrator to validate the UAV performance. The work included a good aerodynamic modeling. However, fuel cell system model was poor because the reference presented only the rate of hydrogen consumption of the PEMFC without further investigation. Furrutter et al. designed a small PEM fuel cell powered UAV. The fuel cell was used in the steady flight condition and it had a power supply of 100 W. Mathematical formulation and derivation was presented in [18] with the achievement of total endurance of 150 s and the designed UAV showed a high drag force. Lapena Rey et al. in [19] succeeded to achieve a 4 h theoretical endurance for 200 W PEM fuel cell powered UAV. The team used a hydrogen generator kit as a fuel source. The MUAV succeeded to fly for 130 min and the hydrogen generator faced thermal activity that limited the UAV flight endurance. Another fuel cell powered UAV was introduced in [12] using a 550 W PEM fuel cell integrated with Li-Po battery. A qualitative and quantitative analysis for hydrogen storage tank was implemented. Thermal management that depended on a fin-and-tube radiator with an external fan was designed instead of having traditional radiators that would have a high surface area. The UAV succeeded to fly for 23 h and 19 min with the achievement of specific energy of 1170 Wh/kg. This paper presents the design, development and validation of a flexible optimization scheme that can be used in the design of any PEMFC powered UAV of a small scale with validated design results. The major contributions of this research concentrate on the integration of PEMFC and its components into the optimization scheme. This in term includes hydrogen subsystem and oxygen subsystem. In addition, the integration of the wing/tail design/selection into the optimization scheme are considered in the analysis. The paper also developed a flexible selection scheme that can be utilized by other interested researches in the field. This scheme composed of hundreds of commercial off-the-shelf components that includes propellers, motors, hydrogen tanks and PEMFCs in which the selected components in the optimization will produce the optimum flight endurance. The results are validated through both wind tunnel tests and the development of a novelty method of the HiL simulation tests. The new developed MDO can be successfully implemented by other researches and be used in MATLAB optimization toolbox easily.
To meet this paper objectives, the paper is organized as follows: Section 2 describes the technical approach. Section 3 includes the UAV conceptual design which defines the UAV baseline configuration. Section 4 includes the preliminary UAV design which includes the development of the multi-disciplinary design analysis (optimization of inner loop modeling). Section 5 includes the detail design which includes the development of the multi-disciplinary design optimization (optimization of outer loop modeling). Section 6 includes the results of the analysis and the validation of results through wind-tunnel tests and hardware in the loop simulation tests.

Technical Approach
The design process is divided into three main parts; UAV conceptual design level, UAV preliminary design level and then a detail design level. In the UAV conceptual design level, the design requirements and constraints are specified. After that, the design space exploration process defines a number of alternatives through a matrix of alternatives that meets the design requirements. Then, the base line aircraft configuration is analyzed through a qualitative and quantitative study for the matrix of alternatives. Based on that, design variables and design matrices are defined.
Aircraft design consists of subsystems which are described as multi-disciplinary complex systems by nature such as fuel cell and propeller subsystems. This requires simplifications and decomposition schemes to provide a well-structured method for analyzing the MDA which are done in the preliminary design stage. Figure 2 illustrates the process of obtaining a logical MDA for aircraft design. First, different design variables and objectives should be defined. These design variables are modeled in their contribution analyses. A design structure matrix is needed for managing the flow of information among different contributing analyses. DSM should be fitted into MDA framework solution scheme that guarantees that all design objectives are met. MDA flexible framework includes compatibility equations and nonlinear solvers. So, the preliminary design level allows for in depth analysis for the baseline UAV defined in the conceptual design. In addition, all powertrain components are modeled in this stage. components and design variables to the given PEM fuel cell and fuselage including the total aircraft mass. In this stage, all components' models are defined but the components themselves are not selected yet because there are many commercial off-the-shelf components that are visualized in the optimization space. As a result, an optimization scheme is developed in this part to select the components of the powertrain and wing geometry and obtain all their design variables that will result in the optimum aircraft endurance. The optimization scheme best described as "mixed discrete continuous" because the solution is structured based on DSM rather than a close form mathematical model for the objective function. Genetic algorithm stands to be one among the best solution methodologies for such optimization scheme. In addition, results comparison and validation through conducting wind tunnel test and hardware in the simulation (HiL) tests are carried on in this section. The UAV conceptual design task is to define a baseline aircraft configuration through defining a design space and defining an initial UAV structural design.

Design Space Definition
The aircraft is a hand-launched UAV which is considered to be a small scale aircraft with maximum take-off weight equals to 85 N. This weight is considered to guarantee hand-launchability. In addition, design requirements are necessary to satisfy the UAV to fly at low altitude and low cruising speed of 200 m and 20 m/s, respectively. This flight speed ensures that the Reynolds number will stay in the laminar region.
A quantitative study of the market was achieved to select the fuel cells that can be integrated with the UAV. Three fuel cells of interest were found; 200 W, 400 W Aeropack Horizon fuel cells and 500 W Ion Tiger fuel cell. Since the design is for long endurance flight mission, it was decided to choose the 200 W fuel cell as the total aircraft weight will not exceed 85 N. Hence, the design of the UAV should have a high aerodynamic efficiency, low drag, low pitching moment and ease of manufacturing. Based on these design requirements, a matrix of alternatives is built to down select the conceptual aircraft design attributes as shown in Table 1. Based on quantitative and qualitative calculations and statistics as well as previous aircraft projects, the final conceptual design was selected from hundreds of possible combinations from the matrix of alternatives.
The selected vehicle configuration is conventional which is used widely for small scale UAVs due its proper gripping surface. Aircraft mass will vary during the optimization process based on the selection process of different propulsion system components The detail design part main objective is to optimize the different propulsion system components and design variables to the given PEM fuel cell and fuselage including the total aircraft mass. In this stage, all components' models are defined but the components themselves are not selected yet because there are many commercial off-the-shelf components that are visualized in the optimization space. As a result, an optimization scheme is developed in this part to select the components of the powertrain and wing geometry and obtain all their design variables that will result in the optimum aircraft endurance. The optimization scheme best described as "mixed discrete continuous" because the solution is structured based on DSM rather than a close form mathematical model for the objective function. Genetic algorithm stands to be one among the best solution methodologies for such optimization scheme. In addition, results comparison and validation through conducting wind tunnel test and hardware in the simulation (HiL) tests are carried on in this section. The UAV conceptual design task is to define a baseline aircraft configuration through defining a design space and defining an initial UAV structural design.

Design Space Definition
The aircraft is a hand-launched UAV which is considered to be a small scale aircraft with maximum take-off weight equals to 85 N. This weight is considered to guarantee hand-launchability. In addition, design requirements are necessary to satisfy the UAV to fly at low altitude and low cruising speed of 200 m and 20 m/s, respectively. This flight speed ensures that the Reynolds number will stay in the laminar region.
A quantitative study of the market was achieved to select the fuel cells that can be integrated with the UAV. Three fuel cells of interest were found; 200 W, 400 W Aeropack Horizon fuel cells and 500 W Ion Tiger fuel cell. Since the design is for long endurance flight mission, it was decided to choose the 200 W fuel cell as the total aircraft weight will not exceed 85 N. Hence, the design of the UAV should have a high aerodynamic efficiency, low drag, low pitching moment and ease of manufacturing. Based on these design requirements, a matrix of alternatives is built to down select the conceptual aircraft design attributes as shown in Table 1. Based on quantitative and qualitative calculations and statistics as well as previous aircraft projects, the final conceptual design was selected from hundreds of possible combinations from the matrix of alternatives. The selected vehicle configuration is conventional which is used widely for small scale UAVs due its proper gripping surface. Aircraft mass will vary during the optimization process based on the selection process of different propulsion system components and wing airfoil. As a result, canard vehicle configuration was not selected due sensitivity to weight changes and weight distribution. While flying wing selection contradicts with the hand launchability design requirement.
For better stability and control of the aircraft, tapered planform with conventional tail was selected. For the wing type, a high wing was selected based on historical data of hand launched fixed wing aircrafts and taking into consideration ease of hand-launchability. In addition, the wing is placed in a way that the aircraft center of gravity is below the wing neutral point, as this will cause an opposing moment due any pitch disturbance.
In addition, a high aspect ratio equals initially to 10 was selected because of high aerodynamic efficiency design requirement criteria. Fuselage shape was cylindrical as it's easier to manufacture taking into consideration propulsion system components' sizes. The landing gear selection resulted in a simple retractable landing gear. Finally, a tractor propeller position was selected for lower drag and pitching moment.
Chemical hydride technology was not feasible in the selection of hydrogen storage because it has a high cost and it needs a continuous maintenance and change of its special catalyst. Carbon fiber showed better weight specifications than aluminum composites tanks. As a result, it was decided to use carbon fiber material in the hydrogen storage tank design process.
Aircraft structural material is composed mainly of carbon fiber due its high strength and stress specifications with lightweight specifications compared to other composites. The wing tubes are made of carbon fiber sandwiched structure with balsa wood due its softness specifications.

Initial Aircraft Structural Design Sizing
The total lift and drag coefficients are obtained by defining the cruise speed and maximum weight. Then, fuselage is designed taking into consideration different components' sizes.
A parametric study was conducted to set initial values for the coefficient of lift, wing surface area and UAV total weight assuming an altitude of 1000 m as shown in Figure 3 using integrated Lift equation. The parametric study resulted in air cruise speed equals to 17.8 m/s, which is close to the design requirement. The coefficient of lift is set equal to 0.4, wing surface area equals to 0.9 m 2 and UAV weight equals to 75 N. Even though lower cruise speed is possible, it will lead to higher coefficient of lift, which is considered less practical and harder to achieve. While using the integrated Drag equation and knowing that maximum power available is 200 W, another parametric study was conducted to specify the maximum coefficient of drag. The value resulted is CD max = 0.067.  Total drag coefficient is composed of induced drag and parasitic drag covered by the following equations: where is total coefficient of lift, represents Oswald efficiency factor ranges from 0.85 to 0.95 depending on how much elliptic the shape is [21]. Parasitic drag is drag generated from other non-lifting parts such as fuselage. Hence, total drag coefficient is: Based on the selected values of the aspect ratio and UAV weight and wing surface area, the wing geometry design specifications are shown in Table 2. For vertical and horizontal wings, sizing is based on coefficient of vertical and horizontal wings covered by the following equations and results summarized in Table 3. Finally, cylindrical streamlined fuselage shape is designed with specifications shown in Table 4.  Total drag coefficient is composed of induced drag and parasitic drag covered by the following equations: where C L is total coefficient of lift, e represents Oswald efficiency factor ranges from 0.85 to 0.95 depending on how much elliptic the shape is [21]. Parasitic drag is drag generated from other non-lifting parts such as fuselage. Hence, total drag coefficient is: Based on the selected values of the aspect ratio and UAV weight and wing surface area, the wing geometry design specifications are shown in Table 2. For vertical and horizontal wings, sizing is based on coefficient of vertical and horizontal wings covered by the following equations and results summarized in Table 3. Finally, cylindrical streamlined fuselage shape is designed with specifications shown in Table 4.

UAV Preliminary Design
The preliminary design task is to develop a multi-disciplinary analysis scheme that allows better analyzing and designing of UAV powertrain and aerodynamics by decomposing the design into smaller subsystems called contributing analyses. Components' models are developed under steady state conditions because dynamic response occurs for a small period of time that will have a negligible effect on total endurance.

Powertrain Contributing Analyses
The powertrain system is decomposed into six CA's: hydrogen tank CA, PEM fuel cell and its subsystems CA, propeller CA, motor CA, and UAV performance analysis.

Hydrogen Tank Contributing Analysis
Hydrogen is stored in a composite fully wrapped pressure vessel, carbon fiber cylinder. The hydrogen tank shape and configuration depends on several criteria such as fuselage and type of tank [22]. The final shape is a cylindrical geometry with hemispherical end caps. The cylinder thickness of the overwrap should be enough to resist both the axial stress and the hoop stress that are caused by the pressure load. The total cylinder thickness can be obtained from the following relation [23,24]: and the total tank mass is obtained as: where x f s is linear load sharing factor of safety to yield, f mount is tank mounting/bosses/tubing mass fraction, m liner is hydrogen tank liner mass, m comp is hydrogen tank composite overwrap mass, m H 2 is the hydrogen mass, m reg regulator mass, r is the cylinder radius. The volume (V) required by this amount of hydrogen stored at a pressure (P) is calculated using the Redlich-Kwong relation [25].

PEMFC Contributing Analysis
The fuel cell model consists of three main parts; fuel cell output voltage, fuel cell hydrogen consumption and fuel cell oxygen consumption. PEMFC output voltage is represented by its polarization curve obtained experimentally by connecting the PEMFC to the load directly and operated until it reaches its full output power [26] as shown in Figure 4. Hydrogen consumption can be calculated from Faraday's law [27]: where U H2 is the hydrogen utilization factor. It is defined as the ratio between the fuel flow that reacts and the input fuel flow. The typical value used is ranged from 0.8 to 0.9 [28]. Also, the fuel, oxidant and water generated output should be taken into consideration. Oxygen consumption can be calculated from Faraday's law as well: .
where U O2 is the oxidant utilization factor which is defined as the ratio between the oxidant flow that reacts and the input air flow rate. The typical value used is 0.21.
where UH2 is the hydrogen utilization factor. It is defined as the ratio between that reacts and the input fuel flow. The typical value used is ranged from 0 Also, the fuel, oxidant and water generated output should be taken into co Oxygen consumption can be calculated from Faraday's law as well: where UO2 is the oxidant utilization factor which is defined as the ratio betw dant flow that reacts and the input air flow rate. The typical value used is 0.2

Air Supply Subsystem
Air will be supplied for both cooling and oxygen reactant supply. 2 axia used to supply air from atmosphere to the fuel cell system. The reason behin ble fans instead of a single fan is mainly because of sizing issues. In addition, allow for better distribution of air through the fuel cell which results in bett factor and better cooling. The power needed is considered as a loss or paras as follows: where ̇ is air flow rate, −1 , is specific heat, J −1 −1 .

Hydrogen Supply Subsystem
For hydrogen supply subsystem, the dead-end mode technique is used ing the hydrogen supply to the stack as shown in Figure 5. This technique i as the simplest mode of connecting the hydrogen supply to the stack especial of a simple small low temperature PEM stack [28]. In dead-end mode, one o sure regulators are connected between the pressure supply and the stack to l drogen supplied pressure to the required operating point.

Air Supply Subsystem
Air will be supplied for both cooling and oxygen reactant supply. 2 axial DC fans are used to supply air from atmosphere to the fuel cell system. The reason behind using double fans instead of a single fan is mainly because of sizing issues. In addition, two fans will allow for better distribution of air through the fuel cell which results in better utilization factor and better cooling. The power needed is considered as a loss or parasitic load [27] as follows: where . m airin is air flow rate, gs −1 , c p is specific heat, J g −1 K −1 .

Hydrogen Supply Subsystem
For hydrogen supply subsystem, the dead-end mode technique is used for connecting the hydrogen supply to the stack as shown in Figure 5. This technique is considered as the simplest mode of connecting the hydrogen supply to the stack especially in the case of a simple small low temperature PEM stack [28]. In dead-end mode, one or more pressure regulators are connected between the pressure supply and the stack to lower the hydrogen supplied pressure to the required operating point.

Electric Motor Contributing Analysis
A brushless direct current electric motor (BLDC) will be used because it has higher efficiency than brushed DC electric motors. BLDC motor model treats the motor subsystem as steady state with assumption that the controller has the ability to state the position of the rotor [29]. Hence, an equivalent motor circuit with conventional lumped parameter is used where motor voltage and torque are obtained as follows [30,31]: where is the BLDC resistance per winding appearing with a factor of two due two windings for each commutating sequence, is the motor speed constant, is motor drag coefficient and is the torque of electric origin. Equations (10) and (11) can be converted to "no load" terms of current I0 and voltageV0 data usually supplied by manufacturers by calculating a drag coefficient.
Taking into consideration the speed controller losses, motor output voltage and current are covered by the following equations [32]: where represents the speed controller duty cycle, is the speed controller output current and is the speed controller output voltage. is the speed controller efficiency and is the speed controller internal resistance.

Propeller Contributing Analysis
Propeller model is based on the non-dimensional coefficients defined by the vortex theory [33,34]. To calculate for the propeller-fuselage interference, a corrected factor is used that is derived based on experimental data done by Lowry in [35,36] as follows (Figure 6): where is the fuselage diameter and dp is the propeller diameter.

Electric Motor Contributing Analysis
A brushless direct current electric motor (BLDC) will be used because it has higher efficiency than brushed DC electric motors. BLDC motor model treats the motor subsystem as steady state with assumption that the controller has the ability to state the position of the rotor [29]. Hence, an equivalent motor circuit with conventional lumped parameter is used where motor voltage and torque are obtained as follows [30,31]: where R m is the BLDC resistance per winding appearing with a factor of two due two windings for each commutating sequence, K v is the motor speed constant, β is motor drag coefficient and τ e is the torque of electric origin. Equations (10) and (11) can be converted to "no load" terms of current I 0 and voltageV 0 data usually supplied by manufacturers by calculating a drag coefficient. Taking into consideration the speed controller losses, motor output voltage and current are covered by the following equations [32]: where D t represents the speed controller duty cycle, I c is the speed controller output current and V c is the speed controller output voltage. η controller is the speed controller efficiency and R controller is the speed controller internal resistance.

Propeller Contributing Analysis
Propeller model is based on the non-dimensional coefficients defined by the vortex theory [33,34]. To calculate for the propeller-fuselage interference, a corrected factor is used that is derived based on experimental data done by Lowry in [35,36] as follows ( Figure 6): where d f is the fuselage diameter and dp is the propeller diameter.

Performance Contributing Analysis
Performance contributing analysis calculates the propeller efficiency, motor efficiency, fuel cell efficiency and endurance of the aircraft as follows:

UAV Simulation Model
The aerodynamic simulation for the proposed aircraft was done using XFLR5 for a range of angle of attacks at a fixed Reynold number of 400,000 as shown in Figure 7

Performance Contributing Analysis
Performance contributing analysis calculates the propeller efficiency, motor efficiency, fuel cell efficiency and endurance of the aircraft as follows:

UAV Simulation Model
The aerodynamic simulation for the proposed aircraft was done using XFLR5 for a range of angle of attacks at a fixed Reynold number of 400,000 as shown in Figure 7. The resulted relations of the coefficients of lift and drag at an angle of for three low speed airfoils; NACA23010, SD7062, NACA 23012 are as follows:  NACA23012 : NACA23010 : Figure 8 represents a parametric analysis of the effect of changing the angle of attack on UAV performance characteristics.

UAV Weight Tabulation CA
The aircraft weight breakdown consists of two main parts. The first part includes the UAV fixed masses which are the airframe that includes the fuselage weight, wing and tail weight, PEMFC weight and landing gear weight. While the second part includes the components that needs to be specified which are the hydrogen tank mass, hydrogen content mass, electric motor mass and propeller mass. Those components correspond to 44.7% of the total UAV weight as shown in Figure 9.

UAV Weight Tabulation CA
The aircraft weight breakdown consists of two main parts. The first part includes the UAV fixed masses which are the airframe that includes the fuselage weight, wing and tail weight, PEMFC weight and landing gear weight. While the second part includes the components that needs to be specified which are the hydrogen tank mass, hydrogen content mass, electric motor mass and propeller mass. Those components correspond to 44.7% of the total UAV weight as shown in Figure 9.

Aerodynamic Contributing Analysis
The meeting point between the propulsion system and aerodynamic system is power. The power available from the propulsion system must equal the power required by aerodynamic in steady flight condition as follows: Available thrust can then be obtained as follows: where C T is the coefficient of thrust, ρ is air density, d p is propeller diameter and ω is rotational speed of the propeller.
To maximize the range, thrust should be minimized. As a result, C L C D should be maximized, while to maximize endurance, power should be minimized. As a result, it is necessary to maximize

Aerodynamic Contributing Analysis
The meeting point between the propulsion system and aerodynamic system is power. The power available from the propulsion system must equal the power required by aerodynamic in steady flight condition as follows: Available thrust can then be obtained as follows: where is the coefficient of thrust, is air density, is propeller diameter and is rotational speed of the propeller.
To maximize the range, thrust should be minimized. As a result, should be maximized, while to maximize endurance, power should be minimized. As a result, it is necessary to maximize 3 2 [37].

Design Structure Matrix
The flow of inputs and outputs between the different CAs requires a manageable system which is described as the design structure matrix. CAs are organized diagonally from upper left to lower right. Flow of data is shown as the lines connecting the CAs in which CAs are executed serially from left to right. The interaction between different CAs requires a flexible framework to execute all of CAs and guarantees that all feedback variables are in harmony with other outputs as shown in Figure 10.

MDA with Simultaneous Analysis
The technique used here is a MATLAB function "fsolve". It evaluates the convergance or divergance percentage between x i and x i+1 by minimizaing the merit function using the the trust region dogleg method for solving nonlinear systems.
The resulted MDA is shown in Figure 11. The solver variables are the fuel cell current and propeller speed. These variables are set to some initial values. As the analysis begins. The new resulted variable values are fed back through the loop to the MDA again.

MDA with Simultaneous Analysis
The technique used here is a MATLAB function "fsolve". It evaluat gance or divergance percentage between xi and xi+1 by minimizaing the me ing the the trust region dogleg method for solving nonlinear systems.
The resulted MDA is shown in Figure 11. The solver variables are the f and propeller speed. These variables are set to some initial values. As the a The new resulted variable values are fed back through the loop to the MDA For a predefined wing airfoil and propulsion system components suc tank, fuel cell, motor and propeller, aircraft total mass and hydrogen cont tained. Coefficient of lift and drag can be obtained based on the calculate mass and the predefined angle of attack. The required cruise speed and t are then obtained from the aerodynamic box. With the guessed value of ro the advance ratio could be obtained and propeller corrected coefficients ca Using the required thrust obtained and coefficient of thrust, the new valu speed can be obtained and fed back. The propeller speed and Cp obtained get the propeller torque. The motor current and voltage can then be obta fuel cell current initial value, the fuel cell voltage can be obtained using th curve equation. Neglecting the motor speed controller losses will lead to th ered from the fuel cell be the same as the motor power. For a predefined wing airfoil and propulsion system components such as hydrogen tank, fuel cell, motor and propeller, aircraft total mass and hydrogen content can be obtained. Coefficient of lift and drag can be obtained based on the calculated total aircraft mass and the predefined angle of attack. The required cruise speed and thrust required are then obtained from the aerodynamic box. With the guessed value of rotational speed, the advance ratio could be obtained and propeller corrected coefficients can be obtained. Using the required thrust obtained and coefficient of thrust, the new value of propeller speed can be obtained and fed back. The propeller speed and Cp obtained will be used to get the propeller torque. The motor current and voltage can then be obtained. With the fuel cell current initial value, the fuel cell voltage can be obtained using the polarization curve equation. Neglecting the motor speed controller losses will lead to the power delivered from the fuel cell be the same as the motor power. Dividing the calculated motor power over the new fuel cell voltage will result in the new fuel cell current. Finally, efficiencies of different components are calculated and the aircraft endurance is obtained.
The new variables of fuel cell current and propeller speed are fed to the compatibility equations before passing them to the fsolve solver. The solver keeps iterating the variables until the compatibility equations are values are less than the presetted tolerance value. The compatibility equations defined are: Those equations are set equal to zero. If the compatibility equation solution is less than the preset tolerance value, the nonlinear system of equations is implicitly solved. This flexible MDA can solve all the UAV design variables for a given UAV geometry and propulsion system. MDA development analysis is important for the development of the multi-disciplinary design optimization (MDO) scheme as it represents the main structure for the MDO and its evaluation.

UAV Detailed Design (Optimization)
Genetic algorithms are commonly used technique to search for optimal operational points especially in cases of complex problems where there is no fixed mathematical formulation frame to solve the optimization problem [38][39][40]. GA optimization technique depends on evolutionary programming to search for optimal results [41]. GA flowchart is shown in Figure 12.

UAV Detailed Design (Optimization)
Genetic algorithms are commonly used technique to search for optimal operational points especially in cases of complex problems where there is no fixed mathematical formulation frame to solve the optimization problem [38][39][40]. GA optimization technique depends on evolutionary programming to search for optimal results [41]. GA flowchart is shown in Figure 12.

Genetic Algorithm Problem Formulation
The objective function of interest is the aircraft endurance. The optimization scheme main task is to maximize the objective function given in Equation (29): The negative sign in Equation (29) is introduced to maximize the endurance. Design variables and constraints are illustrated in Figure 13. The optimization scheme involves discrete and continuous design variables. These design variables are the airfoil index, motor index, propeller index, and hydrogen tank index. Design variables associated with the component indexes are illustrated in Figure 13. Since commercially off-the-shelf (COTS) ready components will be used, thousands of possible combinations of wing geometries, hydrogen tanks, motors and propellers could be visualized. For this study, 3 selected airfoils, 42 available BLDC motors' models supplied from Hacker Motors and AXI Motors are selected for the analysis and 57 propellers ranges from 12-27 inch diameter manufactured by APC are selected and the PEMFC current ranged from 0 to 13 A. For hydrogen tank, the selection process will be done after getting the hydrogen tank design variables solution but the maximum allowable diameter is 180 mm and length of 400 mm to guarantee its fit in the fuselage.
Design constraints are defined to guarantee the feasibility of the proposed solutions. These constraints are the aircraft mass, fuel cell current, propeller tip Mach number, motor and propeller efficiency, and angle of attack as presented in Figure 13.

Genetic Algorithm Problem Formulation
The objective function of interest is the aircraft endurance. The optimization scheme main task is to maximize the objective function given in Equation (29): The negative sign in Equation (29) is introduced to maximize the endurance. Design variables and constraints are illustrated in Figure 13. The optimization scheme involves discrete and continuous design variables. These design variables are the airfoil index, motor index, propeller index, and hydrogen tank index. Design variables associated with the component indexes are illustrated in Figure 13. Since commercially off-the-shelf (COTS) ready components will be used, thousands of possible combinations of wing geometries, hydrogen tanks, motors and propellers could be visualized. For this study, 3 selected airfoils, 42 available BLDC motors' models supplied from Hacker Motors and AXI Motors are selected for the analysis and 57 propellers ranges from 12-27 inch diameter manufactured by APC are selected and the PEMFC current ranged from 0 to 13 A. For hydrogen tank, the selection process will be done after getting the hydrogen tank design variables solution but the maximum allowable diameter is 180 mm and length of 400 mm to guarantee its fit in the fuselage.
The developed genetic algorithm code will be used is the GA toolbox and direct search box provided by MATLAB. The evaluation of the fitness function used is the same that is developed in the design structure matrix as indicated in Figure 13. After defining the optimization scheme, objective functions, all design variables and constraints; the multi-disciplinary design optimization is obtained as presented in Figure 13. This optimization scheme stands as the baseline flexible framework that can vary depending on the implementation.

Optimization Results
The problem formulation of the multi-disciplinary design optimization analysis is fully summarized in Tables 5 and 6. This analysis was conducted assuming steady level flight performance using genetic algorithm as the optimization routine. The evolution history of the optimization scenarios are also presented in Figures 14 and 15 in which the top graph represents the penalty function value while the bottom one represents the stopping criteria. Figure 14 shows the obtained results for the maximum endurance while Figure  15 shows the results for the maximum range. Design constraints are defined to guarantee the feasibility of the proposed solutions. These constraints are the aircraft mass, fuel cell current, propeller tip Mach number, motor and propeller efficiency, and angle of attack as presented in Figure 13.
The developed genetic algorithm code will be used is the GA toolbox and direct search box provided by MATLAB. The evaluation of the fitness function used is the same that is developed in the design structure matrix as indicated in Figure 13. After defining the optimization scheme, objective functions, all design variables and constraints; the multi-disciplinary design optimization is obtained as presented in Figure 13. This optimization scheme stands as the baseline flexible framework that can vary depending on the implementation.

Optimization Results
The problem formulation of the multi-disciplinary design optimization analysis is fully summarized in Tables 5 and 6. This analysis was conducted assuming steady level flight performance using genetic algorithm as the optimization routine. The evolution history of the optimization scenarios are also presented in Figures 14 and 15 in which the top graph represents the penalty function value while the bottom one represents the stopping criteria. Figure 14 shows the obtained results for the maximum endurance while Figure 15 shows the results for the maximum range. airfoil. While in order to maximize range, should be maximized which corresponds to an angle of attack equals to 7° for the same UAV. UAV with NACA23010 gave better performance characteristics in terms of endurance and range than other airfoils. The hydrogen tank diameter is 144 mm and length of 322 mm. The calculated weight is 1.65 kg. The nearest commercially available carbon fiber tank has specifications of 130 mm diameter, 300 mm length and 1.55 kg weight. The total UAV weight constraint limited the hydrogen tank size because the bigger the hydrogen tank, the more compressed hydrogen will be.
Hydrogen consumption is always a function of the generated fuel cell current. Thus, the lower the fuel cell current, the lower the hydrogen consumption and a higher endurance is achieved. The resulted fuel cell current value is the lowest possible operational value in both cases. Motor index chosen for both cases is the same due its lower Kv compared to other motors indices. A motor with low Kv produces lower rotational speed than a high Kv motor at the same voltage. In order to compensate for that, a propeller with bigger diameter relatively is required. Hence, the resulted propeller indices for both scenarios are relatively high especially in the first scenario where higher thrust is required.    Other performance parameters are listed in Table 8. Power produced in the second scenario is higher than the first one. This is because the rotational speed in the second scenario is higher. As a result, the hydrogen flow rate in the second scenario is higher as well. Finally, the coefficients of lift and drag for the first scenario are higher than the second scenario. This is due to the fact that angle of attack that produces the longest flight endurance is higher than the angle of attack that produces the maximum flight range.  The penalty function graph indicated in blue dots represents the mean value of the penalty function while the black dots represent the best value at each generation. The penalty function evaluation consists of two main parameters namely constraints and fitness functions. For each parameter, the penalty function picks a design entity from the created population mentioned above at each generation for checking its consistency with the design parameter. In the case of consistency, the fitness function value is set to be equal to the value of the penalty function. Otherwise, the value of the penalty function is also set to be equal to the maximum fitness value among the worst feasible members of the population. For the subsequent generation, the values chosen are those with the minimum penalty function values. As the optimization routine process is carried on, each new generation is converged to the best penalty function value. This is because of the reproduction process of the new generation that takes into consideration the relation of the best penalty function and new parents selection process.
Bottom graph included the stopping criteria in which if any of the termination criteria are met, it will cause the optimization process to be stopped. The stall time criteria evaluates the improvement of the objective function over a specified time interval. The time interval is set equal to the time criteria, which is set to be equal to infinity. The stall generations criteria evaluate the weighted average change in fitness value for a stall generation that is equal to almost 50. The specified value for the change is sit equal to a tolerance value of 1 × 10 −6 . Finally, the number of generations is set to be 100. If the optimization process reaches the maximum number of generations, the process is terminated. As can be seen from Figures 14 and 15, the optimization process for maximum endurance was terminated after 59 generations and 2403 fitness functions evaluations. As for the range optimization process, a total number of 2841 fitness function values were evaluated and 70 generations were created. For both scenarios, the stall (G) generation criteria was satisfied as the change in the stall generation was too small, less than 1 × 10 −6 , and the penalty function value converged to the minimum objective function value.
Selected indices for both scenarios resulted from the optimization process are listed in Table 7. The selected airfoil that resulted in both maximum endurance and maximum range is NACA23010. As stated previously, to maximize endurance, attack equals to 7 • for the same UAV. UAV with NACA23010 gave better performance characteristics in terms of endurance and range than other airfoils. The hydrogen tank diameter is 144 mm and length of 322 mm. The calculated weight is 1.65 kg. The nearest commercially available carbon fiber tank has specifications of 130 mm diameter, 300 mm length and 1.55 kg weight. The total UAV weight constraint limited the hydrogen tank size because the bigger the hydrogen tank, the more compressed hydrogen will be.
Hydrogen consumption is always a function of the generated fuel cell current. Thus, the lower the fuel cell current, the lower the hydrogen consumption and a higher endurance is achieved. The resulted fuel cell current value is the lowest possible operational value in both cases. Motor index chosen for both cases is the same due its lower Kv compared to other motors indices. A motor with low Kv produces lower rotational speed than a high Kv motor at the same voltage. In order to compensate for that, a propeller with bigger diameter relatively is required. Hence, the resulted propeller indices for both scenarios are relatively high especially in the first scenario where higher thrust is required.
Other performance parameters are listed in Table 8. Power produced in the second scenario is higher than the first one. This is because the rotational speed in the second scenario is higher. As a result, the hydrogen flow rate in the second scenario is higher as well. Finally, the coefficients of lift and drag for the first scenario are higher than the second scenario. This is due to the fact that angle of attack that produces the longest flight endurance is higher than the angle of attack that produces the maximum flight range.

Wind Tunnel Tests
To validate the optimization tool results, a wind tunnel experimental setup was developed in an open circuit low speed wind tunnel that has a cross section of 90 by 72 cm and a maximum airspeed of 40 m/s as shown in Figure 16. The test stand is composed of a load cell that measures thrust force produced by the propeller. Power source used is the Aeropack fuel cell that can supply up to 200 W. Hydrogen tank chosen was ultra-weight carbon fiber 1.25 L that includes high purity compressed hydrogen with a two step down pressure regulators.
a load cell that measures thrust force produced by the propeller. Power source used is the Aeropack fuel cell that can supply up to 200 W. Hydrogen tank chosen was ultra-weight carbon fiber 1.25 L that includes high purity compressed hydrogen with a two step down pressure regulators.
Selected propeller was 16 by 11 inches and electric motor Hacker C50 15 XL electric motor with kv equals to 173 RPM/volt. The chosen propeller and motor were available in the lab and has the closest performance characteristics to those resulted from the optimization scheme.
Calibrated instruments and their accuracies used in the experimental setup are described in Table 9. Uncertainty analysis was conducted for all measured and calculated variables using the Kline-McClintok method explained in [35] as follows: where R is the result, is the uncertainty result, are the dependent variables and are the associated uncertainty of those variables. Experiments were conducted in a lab environment and under conditions where the airspeed equals 16 m/s in which the used fuel cell input current and voltage values are obtained from optimization scheme results.  Selected propeller was 16 by 11 inches and electric motor Hacker C50 15 XL electric motor with kv equals to 173 RPM/volt. The chosen propeller and motor were available in the lab and has the closest performance characteristics to those resulted from the optimization scheme.
Calibrated instruments and their accuracies used in the experimental setup are described in Table 9. Uncertainty analysis was conducted for all measured and calculated variables using the Kline-McClintok method explained in [35] as follows: where R is the result, w R is the uncertainty result, xn are the dependent variables and w n are the associated uncertainty of those variables. Experiments were conducted in a lab environment and under conditions where the airspeed equals 16 m/s in which the used fuel cell input current and voltage values are obtained from optimization scheme results. Table 10 compares results obtained from optimization techniques and those obtained from experiments. All obtained results are acceptable with error range less 10%. The least error is in the thrust required while the maximum is in the propeller advance ratio. This due the fact that cruise velocity contributes to the advance ratio with a factor of 2. That is the double. Propeller model that uses a slow factor proofed to be a very good choice when designing a small UAV. In addition, fuel cell model that depends on the polarization curve had also showed a good agreement with the experimental results.  Fuel cell electrical DC output power was measured to be 96.31 W. The fuel cell used 4.87 W as the feedback power for the FC microcontrollers and small DC motor. The electric input motor was 91.44 W. The motor produced a net mechanical power as a rotational speed equals to 81.28 W. This represents a 10.15 W losses due speed controller and motor internal resistance. Propeller output power as a propulsive power was measured to be 74.58 W with 6.71 W losses. System total utilization of power is 77.43% taking into consideration the fuel cell as input power and propeller propulsive power as output power. The power utilization is considered to be high with compared to other propulsion systems.

HiL Simulation Test
It was decided to perform a complete HiL simulation test to validate the flight endurance obtained from optimization before building the demonstration aircraft. The PEMFC is considered as a very expensive component that should be treated carefully. Performing a HiL simulation test would allow to better study the propulsion system behavior especially the fuel cell behavior during the entire flight without the need to perform a real flight test. The HiL test is provided with transitional hard stops and fail-safes to avoid system failure or components damage.
HiL simulation consists of three main components; hardware simulation, software simulation and interface components. The design of HiL simulation requires a tradeoff between the hardware components and software components.
The proposed design is shown in Figure 17. The hardware simulation includes the 200 W Horizon Aeropack PEMFC system (Horizon, Singapore), compressed hydrogen supplied from ultra-weight carbon fiber tank with 3.0 L high purity hydrogen content, a Hacker BLDC motor model number A60 24S (Hacker Motor GmbH, Ergolding, Germany) with a voltage constant of 200 RPM/V and 16.5 × 14 inch propeller supplied by APC (APC, Woodland, CA, USA). The software simulation includes the aircraft aerodynamics, UAV control unit and the propeller model. Interface components regulate and monitor the hardware components and feed the inputs to the software simulation and delivers throttle command to the electric motor. 458.7 min with 475.3 min obtained from optimization results. The associated error percentage is 3.62% which considered very acceptable. Even though a total hydrogen content of 24.2 mole was included in the hydrogen tank, only 20.88 mole was utilized. As the hydrogen content was decreasing due hydrogen supply to the fuel cell, the hydrogen pressure inside the vessel was decreasing as well. When the hydrogen pressure decreased to less than 0.4 bar which is the required pressure by the fuel cell, the fuel cell suffered difficulties in operation which led to stop the HiL simulation test. The UAV cruise velocity average is 13.25 m/s which is also close to the cruise velocity resulted from optimization.  Figure 17. Schematic of the HiL simulation. The PEMFC stack supplies the required electrical power to rotate the electric motor shaft via a DC electrical bus connection between the PEMFC and BLDC motor. The shaft rotates causing the propeller to prop with a specific rotational speed controlled by the speed controller which in turns receives a throttle command from the UAV control unit via signal generator. Generated thrust values are passed to the aircraft model to calculate the actual drag and lift forces and actual airspeed velocity with the angle of attack being defined previously. In addition, the propeller model calculates the coefficient of thrust and the corresponding thrust required which is used to obtain the required velocity. Both the actual velocity and required velocity are fed into the UAV control unit. Based on the difference between the actual and required velocities, the velocity hold controller generates a throttle command to the speed controller of the electric motor.
HiL simulation test results are shown in Figure 18. The measured endurance was 458.7 min with 475.3 min obtained from optimization results. The associated error percentage is 3.62% which considered very acceptable. Even though a total hydrogen content of 24.2 mole was included in the hydrogen tank, only 20.88 mole was utilized. As the hydrogen content was decreasing due hydrogen supply to the fuel cell, the hydrogen pressure inside the vessel was decreasing as well. When the hydrogen pressure decreased to less than 0.4 bar which is the required pressure by the fuel cell, the fuel cell suffered difficulties in operation which led to stop the HiL simulation test. The UAV cruise velocity average is 13.25 m/s which is also close to the cruise velocity resulted from optimization.

Discussion and Comparison of Results
The UAV propulsion system design using the classical approach depends on selecting the propeller-motor combination based on their combined efficiency as shown in Figure 19. The classical approach design resulted in selecting APC propeller of 20 by 13 inch Figure 18. HiL simulation test results.

Discussion and Comparison of Results
The UAV propulsion system design using the classical approach depends on selecting the propeller-motor combination based on their combined efficiency as shown in Figure 19. The classical approach design resulted in selecting APC propeller of 20 by 13 inch for the Hacker A60 24S motor, the one resulted from optimization. After specifying the operational point, the fuel cell current and voltage is obtained. The calculated endurance for this set of propeller-motor combination is 383 min. This represents a total of 19.4% difference from the endurance resulted from the developed optimization routine.

Discussion and Comparison of Results
The UAV propulsion system design using the classical approach depends on selecting the propeller-motor combination based on their combined efficiency as shown in Figure 19. The classical approach design resulted in selecting APC propeller of 20 by 13 inch for the Hacker A60 24S motor, the one resulted from optimization. After specifying the operational point, the fuel cell current and voltage is obtained. The calculated endurance for this set of propeller-motor combination is 383 min. This represents a total of 19.4% difference from the endurance resulted from the developed optimization routine. Going back to Figure 1, the PUMA fuel cell aircraft achieved the maximum endurance for small scale fuel cell powered UAV. The take-off mass ratio to endurance was 20 min/kg. The maximum obtained result from the optimization scheme was 475 min. This represents a ratio of 59 min/kg which is almost three times PUMA's ratio.

Conclusions
In this study, a flexible optimization tool was successfully developed for the design of a hand launched fixed wing fuel cell-powered UAV. The design was decomposed into three levels. In the conceptual design level, the aircraft baseline configuration was defined. While in the preliminary design level, the aircraft design was decomposed into subsystems defined as contributing analyses. For the hydrogen storage vessels, it was found that Going back to Figure 1, the PUMA fuel cell aircraft achieved the maximum endurance for small scale fuel cell powered UAV. The take-off mass ratio to endurance was 20 min/kg. The maximum obtained result from the optimization scheme was 475 min. This represents a ratio of 59 min/kg which is almost three times PUMA's ratio.

Conclusions
In this study, a flexible optimization tool was successfully developed for the design of a hand launched fixed wing fuel cell-powered UAV. The design was decomposed into three levels. In the conceptual design level, the aircraft baseline configuration was defined. While in the preliminary design level, the aircraft design was decomposed into subsystems defined as contributing analyses. For the hydrogen storage vessels, it was found that storing hydrogen in a composite fully wrapped vessels such as carbon-fiber increased the hydrogen content because it has better performance characteristics than other types of tanks in terms of pressure, volume and weight. In addition, a fuel cell model that used polarization curve fit showed acceptable results in terms of obtaining the required voltage and current. A propeller model that includes the effect of propeller-fuselage interference showed consistent results. The inclusion of components' masses in the design process resulted in the increase of the endurance-takeoff mass ratio.
The contributing analyses were combined in a multidisciplinary analysis through a design structure matrix. The developed MDA is considered as a flexible framework scheme for the design of small-scale fuel cell-powered UAVs and can be modified based on the design requirements. The MDA is solved iteratively by setting the non-linear compatibility equations equal to zero. The solver used in the analysis was "fsolve" function in Matlab that uses the trust region dogleg method for solving nonlinear systems.
In the detailed design process, an optimization tool for the design of the UAV powertrain was developed and analyzed. The optimization technique used was genetic algorithm in which hydrogen tank, motor and propeller combinations were varied until an optimal solution was defined. The genetic algorithm is considered to be a very good choice to such design processes because it has flexibility in defining many near optimal designs. Results showed that the lower the fuel cell current, the better aircraft endurance is achieved. In addition, a motor-propeller combination in which bigger diameter propeller and lower motor voltage constant proofed to improve the UAV endurance.
The design tool results were verified using a wind tunnel experimental setup and HiL simulation tests. The setup was developed in which it simulated the real flying conditions in terms of air speed velocity, temperature and ambient pressure. The result error margin was less than 10% showing the optimization results and experimental results were consistent. Further work will include real time flight tests for the aircraft to validate the aerodynamic models.
The measured endurance was 458.7 min as compared with 475.3 min obtained from the optimization results. The associated error percentage is about 3.62% which is considered to be very acceptable for this type of work.

Acknowledgments:
We would like to thank Mohamed Zulkefle for his valuable contribution to the work especially his assistance in conducting the wind tunnel tests.

Conflicts of Interest:
The authors declare no conflict of interest. Nomenclature a 0 zero lift slope for the airfoil, rad −1 b w wing span, m C D l total coefficient of drag C L total coefficient of lift C T total coefficient of thrust C HT horizontal tail coefficient of volume C VT vertical tail coefficient of volume c p specific heat J g −