Optimal Drone Design

Authors: Olivier Forster, Mathieu Salz, Sylvain Pagnier

Abstract

Small, electrically powered drones are used often in recreation, structural inspections, and photography but their designs are far from optimal. This design project sought the minimization of such a drone’s weight for thermal, structural, and aerodynamic demands. The key systems were thus the battery storage, the drone structure, and the means of generating lift. Each component had a respective, individual optimization goal as follows: minimizing heat transfer resistance, maximizing the structural safety factor, and maximizing the generated lift. These goals were each achieved in individual optimizations generating three different designs characteristic of the intuitive design behind them. However, a multi-disciplinary optimization methodology was required here as all three sub-systems share competing interests. The conducted multi-disciplinary optimization yielded a design characterized by thin wall beams while satisfying the needs of each subsystem.

Design Problem statement

The overall system designed will be a small electrically powered drone with in-housing battery and vertically powered thrust. The drone ought to be capable of vertical lift and take-off, as well as satisfactory maneuverability and flight duration while in the air. Small drones are widely used in recreation, structural inspections, and photography. Thus, the overall system goal will be to minimize drone weight as it is a suitable metric to judge drone design success across all these fields as a smaller weight facilitates drone usage.

Trade-offs will likely occur between all three aforementioned sub-systems; the larger the battery, the more heat dissipation and thrust required (as weight increases). Simi- larly, the structure needs be strong enough to resist the relevant forces, but this increase in weight will require an increase in thrust generation. In addition, the heat conduction and structural integrity will both be determined by the shape or size of the drone body and arms. Therefore, they will strongly influence each other, and a balance will need to be found. It is unlikely that improving each design by itself will achieve an overall optimal goal, as the competition between systems is significant. To clarify these design interac- tions, Figure 1.2 has been created to summarize the interactions between subsystems for conceptual understanding. The interactions between the subsystems will be elucidated and explained in detail in Chapter 2.

In Figure 1.2 one may see that the goals are conflicting and will thus necessitate trade offs. As heat transfer resistance goes down, the structural strengths are impacted which also necessarily must handle the thrust of the system. Additionally, each system is typi-

Figure 1.1: An example electric drone [1].

Figure 1.2: Subsystem interaction diagram.

cally strengthened or facilitated by the adding of material. However, not only is the over- all objective to minimize weight, but also additional weight worsens the performance of the structural and thrust subsystems. This balance will define the performance of the system level performance.

Design Optimization Problem Formulation

As mentioned in Chapter 1, the system level problem dealt with the minimization of a drone’s weight while accounting for three competing subsystems: battery storage, thrust, and structure. This minimization problem may be formulated as:

$$ \min_{X} \text{mass}(X) \tag{2.1}$$

To better clarify the formulation of this design optimization problem, it is necessary to elucidate all relevant design parameters. These design parameters have been combined in Table 2.1 and represent all variables used in the attempted multi-disciplinary optimiza- tion (MDO) code. Additionally, it should be noted that the battery was chosen to generate a constant 10 Watts of heat energy with the dimensions of w_battery = 1 in, l_battery = 1 in, and hbattery = 2 in.

To better outline the formulation of the attempted MDO problem, the variables rele- vant to each subsystem is shown in a more detailed Figure 4.2. From this figure, one may see that the battery subsystem merely contains shared variables (in fact, all of its variables are necessarily shared) that the other subsystems use. The MDO coupling instead comes from the interaction between the thrust and structural systems. This is because for the structural characteristics of the system to be determined, the thrust force of the propellers

Table 2.1: Summary of the design variables used in the attempted MDO.

Variable	Subsystem	Lower bound	Upper bound	Physical de- scription
t_frame	battery	0.003175 (m)	0.0254 (m)	Drone frame thickness sur- rounding the battery
A_convection	battery	0 (m)	(wbattery ∗ hbattery + lbattery ∗ wbattery + lbattery ∗ hbattery ) (m)	Battery area open to direct convection
V_flight,avg	battery	0 (m/s)	11.176 m/s	Speed at which the drone is fly- ing
ρ_frame	battery, struc- ture	9971 (kg/m3)	139972 (kg/m3)	Drone frame density
w_arm	structure	0.01 (m)	wbattery/2 (m)	Width of the structural arm
l_arm	structure	0.5658 (m)	10 (m)	Length of the structural arm
Cross sectional diameter (D)	thrust	.001 (m2)	0.07 (m2)	Cross sectional diameter for the thrust subsys- tem
radius	thrust	0.08 (m)	0.8 (m)	Radius for the thrust subsys- tem
Helix angle	thrust	0 (deg)	75 (deg)	Helix angle for drone blades
Drone mass	structure, thrust	0 (kg)	100 (kg)	Mass of the over- all drone
Thrust force	structure, thrust	0	100000 (N)	Thrust force provided by the blades

is needed, therefore, however, for the thrust to be calculated the mass of the drone is also required.

The attempted MDO used a non hierarchical formulation to handle the relevant shared and coupled variables, with the interactions of the variables being taken from 4.2. This code was subsequently implemented into MATLAB using an altered NoHiMDO code.

Figure 2.1: A more detailed subsystem interaction diagram used to define the local, shared, and coupled variables relevant for the attempted MDO.

Additionally, it should be noted that all relevant intermediary functions are omitted from this chapter for brevity and clarity. For additional detail on the intermediary variables and functions required to solve this problem one may read Chapter 3 where the expressions are explained in detail.

Finally, comments on the results of this MDO may be found in section 4.5, though it ought to be noted that the problem was hindered by its high dimension.

Subsystem optimisation Problem

3.1 Overview

This chapter describes the relevant subsystem optimization problems that were used to form the aforementioned system level design optimization problem. These relevant subsystems are battery storage, thrust, and structural design. Each subsystem will de- scribe their respective design optimization problem formulation, their utilized analysis models, and any additional analysis that was performed to contextualize the optimized solution.

3.2 Battery Storage

This section deals with the formulation of the battery storage subsystem optimization problem. The following sub-sections will first detail the relevant subsystem design opti- mization problem and then explain the analytical theory behind it. Finally, there will be a brief discussion on whether this subsystem optimization problem was well-posed

3.2.1 Battery storage optimal design problem

The battery storage sub-system aims to design the housing for a typical on-board drone battery. Here it was chosen to use a lithium-ion battery due to the available in- formation on such a battery. The associated heat generation power has been taken to be approximately 10 Watts, based on available literature [4]. Additionally, this battery was ascribed typical dimensions of a 1 ∗ 1 ∗ 2 inch³ volume.

To best dissipate this heat, this system sought to minimize the heat transfer resistance around the battery. Thus, it was necessary to calculate the conductive and convective heat transfer resistances (radiation is ignored as the temperatures are quite low). These resistances inherently depend upon material and geometrical properties, thickness, and topology. The relevant variables can be summarized as:

$$
\text{System overall thermal resistance} = R_{\text{th,tot}} \; \left( \frac{K}{W} \right)
\tag{3.1}
$$

$$
\text{Battery surface area exposed to convective cooling} = A_{\text{convective cooling}} \; (\text{m}^2)
\tag{3.2}
$$

$$
\text{Battery surface area exposed to conductive cooling} = A_{\text{conductive cooling}} \; (\text{m}^2)
\tag{3.3}
$$

$$
\text{Frame material density} = \rho_{\text{frame}} \; \left( \frac{\text{kg}}{\text{m}^3} \right)
\tag{3.4}
$$

$$
\text{Frame material conductivity} = k_{\text{frame}} \; \left( \frac{W}{\text{m·K}} \right)
\tag{3.5}
$$

$$
\text{Battery width} = w_{\text{battery}} \; (\text{m})
\tag{3.6}
$$

$$
\text{Battery length} = l_{\text{battery}} \; (\text{m})
\tag{3.7}
$$

$$
\text{Battery height} = h_{\text{battery}} \; (\text{m})
\tag{3.8}
$$

$$
\text{Thickness of frame around battery} = t_{\text{frame,battery}} \; (\text{m})
\tag{3.9}
$$

$$
\text{Average drone flight speed} = v_{\text{flight,avg}} \; \left( \frac{\text{m}}{\text{s}} \right)
\tag{3.10}
$$

$$
\text{Convective heat transfer coefficient} = h_c \; \left( \frac{W}{\text{m}^2\text{·K}} \right)
\tag{3.11}
$$

$$
\text{Prandtl number of air} = Pr_{\text{air}} = .71
\tag{3.12}
$$

$$
\text{Air kinematic viscosity} = v_{\text{air}} = 18e{-6} \left( \frac{\text{m}^2}{\text{s}} \right)
\tag{3.13}
$$

$$
\text{Flow length scale} = L \; (\text{m})
\tag{3.14}
$$

$$
\text{Reynolds number of the flow around the battery} = Re_{\text{flow,battery}}
\tag{3.15}
$$

$$
\text{Air thermal conductivity} = k_{\text{air}} = .025 \left( \frac{W}{\text{m·K}} \right)
\tag{3.16}
$$

$$
\text{Nusselt number of the flow around the battery} = Nu_{\text{battery}}
\tag{3.17}
$$

$$
\text{Area open to convection after heat conduction normalized by the area of heat conduction} = c(\cdot)
\tag{3.18}
$$

Using these variables, one may express the optimization problem as:

$$ \min_{\mathbf{x}} \; R_{\text{th,tot}}(\mathbf{X}) \left(
\left(
\frac{t_{\text{frame,battery}}}{k_{\text{frame}} \, A_{\text{conductive cooling}}}
+
\frac{1}{h_c \, c \, * \, (A_{\text{conductive cooling}})}
\right)^{-1}
+
h_c * A_{\text{convective cooling}}
\right)^{-1}
\left( \frac{K}{W} \right)
$$

$$t_{\text{frame,battery}} < -0.125 \; \text{in}
$$

$$
t_{\text{frame,battery}} < 1 \; \text{in}
$$

$$
A_{\text{convective cooling}}
<
w_{\text{battery}} \, h_{\text{battery}}
+
l_{\text{battery}} \, h_{\text{battery}}
+
w_{\text{battery}} \, l_{\text{battery}}
+
h_{\text{battery}} \, m^2
$$

$$A_{\text{convective cooling}} < 0 \; \text{m}^2
$$

$$v_{\text{flight,avg}} < 0 \; \text{mph}
$$

$$
v_{\text{flight,avg}} < 25 \; \text{mph}
$$

$$\rho_{\text{frame}} < 9971 \; \frac{\text{kg}}{\text{m}^3}
$$

$$
\rho_{\text{frame}} < 13997 \; \frac{\text{kg}}{\text{m}^3}
$$

$$
\text{(3.19)}
$$

Reasoning for this objective function and corresponding constraints is given in sub- section 3.2.2.

3.2.2 Utilized battery storage analysis models

The relevant objective function is the minimization of R_th,tot as it was desired to reduce heat buildup within the battery as this may cause electrical failures or losses of efficiency. This overall thermal resistance was simplified to represent a model concerning only par- allel conductive and convective heat transfer (as radiation ought to be negligible at these low, near environment temperatures). Additionally, the conductive route will terminate with heat convecting just as in the parallel convective branch. However, it is important to note that this parallel heat convection branch segment will be scaled by the parameter c, where c represents the ratio of the area open to convection after heat conduction has occurred compared to the area open to conduction (A_conduction). This c value will be set to

.5 for the bulk of the simulations, but varied in a subsequent analysis to analyze the effect of this parameter. One should also note, for model simplification, that the environment temperature will be taken as a constant 25 degrees Celsius.

The constraints also stemmed from a physical nature. It was necessary to restrict the dimensions of t_frame,_battery to be within .125 and 1 inch as the former is a satisfactory lower bound gauge thickness for common laser cut parts whilst the latter simply exceeds the length scales of the drone. A_convective_cooling was also restricted to be positive to satisfy real- ism, while its upper bound was limited by the surface area required to mount unto this battery. V_flight,avg was taken to represent the drone speed and thus must be positive, but was limited by an upper bound of 25 mph representing the maximum speed at which a novice pilot is able to control their drone. It should also be noted that turbulence and streamline deflection was ignored; essentially, the heat transfer over the battery was sim- plified to the heat transfer of a plate plate with a length equal to the battery length. This is a necessary simplification as modelling turbulence analytically is quite arduous and would necessitate more complex simulations such as computational fluid dynamic sim- ulations. The bounds for the frame material density are slightly more interesting as they represent the bounds not necessarily of density, but an expression of the bounds for ther- mal conductivity that have simply been subbed in. Thermal conductivity was restricted lie between 26 (monel) and 416 (silver) ( ^W ) and related through Figure 3.3. The relevant intermediary, physical relationships may be expressed as:

Figure 3.1: Relationship between thermal conductivity and material density used in the intermediary calculations of the optimization problem [2].

$$A_{\text{conductive cooling}} = \left( 2(w_{\text{battery}} h_{\text{battery}} l_{\text{battery}} h_{\text{battery}} w_{\text{battery}} l_{\text{battery}} h_{\text{battery}}) \right) – A_{\text{convective cooling}}
\; (\text{m}^2)
\tag{3.20}
$$

$$k_{\text{frame}} = \left( 1.5 + 0.003 \cdot (\rho_{\text{frame}} – 1805) \right)
\left( \frac{W}{\text{m·K}} \right)
\tag{3.21}
$$

$$h_c = \left(
0.664 \; Re_{\text{flow,battery}}^{0.5}
\; Pr^{1/3}
\; \frac{k_{\text{fluid}}}{L}
\right)
\left( \frac{W}{\text{m}^2\text{·K}} \right)
\tag{3.22}
$$

$$Re_L = \frac{v_{\text{flight,avg}} \; l_{\text{battery}}}{v_{\text{air}}}
\tag{3.23}
$$

Comments on the nature of the optimization problem

To understand the nature of the optimization problem one may first rely on their intu- ition. The analysis models, as described in section 3.2.2, behind this specific optimization problem were simplified extensively to guarantee a logical, feasible solution. This simpli- fication facilitates not only the future construction of more complicated thermal analysis optimization models but also easy contextualization of the optimal results. The optimal battery storage solution, when analyzed in isolation, is explained in detail in section ?? but will be briefly explained here.

The optimal battery storage solution (presented in full at section 4.2), for the base case of c = .5, had the maximum allowed area open to convection as well as the highest thermal conductivity frame material. This is quite an intuitive result as the thermal resistance due to convection and conduction in series is greater than thermal convection only (even with a doubling of the heat convection area due to the heat conduction). Additionally, the max desired flight speed was desired as the greater wind velocity results in more effective thermal convection. Thus, the optimal solution results can be said to match the intuitive optimized result, supporting the well-posed nature of the problem.

More quantitatively, one may analyze the constraint activity to understand if the prob- lem was well bounded, again for the base case. MATLAB outputs a feasibility ranking of 0 for the optimal solution where the satisfied inequality constraints count as 0. Addition- ally, one may analyze the constraint activity as in table 3.1 where the constraints are given in the same order as in the design problem statement and are listed in SI units.

Table 3.1: Constraint activity for the optimal solution to the isolated battery storage opti- mization problem.

1	2	3	4	5	6	7	8
-.0032	.0032	-.0032	.0032	-11.1760	11.1760	-1.0345e5	1.0345e5

From table 3.1 one may see that either the lower or upper bound is practically satis- fied. Though, it should be recalled that the 7th and 8th constraints indirectly represent the constraints on the thermal conductivity but are written as expressions of frame density. As these constraints are active and the problem is feasible, as judged by MATLAB and visual inspection of the constraint activities, one may say that this problem was both well- bounded and well-posed. Again, as this optimization also matched the analytical intu- itive optimized values, it was not necessary to perform more rigorous and cost-intensive investigations of the nature of the optimization problem.

Drone Structure

Drone Structure optimal design problem

This section relates to the structural optimization of the drone. This subsystem aims to ensure that the drone will be as structurally sound as possible by minimizing the stresses that the drone will experience. The structure of the drone is assumed to be a rectangular body which houses the battery with four cantilevered arms that support the propellers at an angle of 45 degress as shown in Fig. 3.2. The stresses in question are the bending stress experienced by the arms due to the thrust force generated by the propellers and the torsional stress experienced by the battery housing.

Figure 3.2: Drone Structure

The relevant quantities are shown below:

arm width = w_a(m) (3.24)

arm length = L_a(m) (3.25)

arm thickness = h_a(m) (3.26)

battery width = w_b(m) (3.27)

battery length = l_b(m) (3.28)

battery width = w_b(m) (3.29) thickness of frame around battery = t_frame(m) (3.30) maximum blade radius = r(m) (3.31)

thrust force = F_t(m) (3.32)

tensile yield strength of PLA = σ_y(m) (3.33)

shear yield strength of PLA = τ_y(m) (3.34)

Young’s Modulus of PLA = E(m) (3.35)

Maximum Deflection Angle = θ_max(^◦) (3.36)

The objective function and constraints are defined below:

$$
\min_{\mathbf{x}} \; \text{Combined Stress}(\mathbf{X}) =
\left( \frac{6 F_t L_a}{w_a h_a^2 \sigma_y} \right)
\left( \frac{F_t L_a}{(w_b + t_{frame})(h_b t_{frame}) t_{frame} \sigma_y} \right)
$$

$$
1.5 – \left( \sigma_y \cdot \frac{w_a h_a^2}{6 F_t L_a} \right) < 0
$$

$$
1.5 – \left( \tau_y \cdot \frac{(w_b + t_{frame})(h_b t_{frame}) t_{frame}}{F_t L_a} \right) < 0
$$

$$
\frac{6 F_t L_a^2}{E w_a h_a^3} – 0.175 < 0
$$

$$
10 – \frac{L_a}{h_a} < 0
$$

$$
10 – \frac{L_a}{w_a} < 0
$$

$$
10 – \frac{h_b}{t_f} < 0
$$

$$
10 – \frac{w_b}{t_f} < 0
$$

$$
(3.38)
$$

The objective function is the product of both stresses. Because it is such a large number, it was normalized by dividing each stress by σ_y which will not affect teh result However, a caveat associated with the objective function being the product of both kinds of stresses is that we risk having one of the types of stresses being lower than desired, therefore a minimum safety factor of 1.5 is imposed on each individual constraint. In addition, the bending angle of the arm is constrained to be less than 0.175 radians. The remaining constraints are theory validity constraints, ensuring that the beam is indeed slender and that the battery housing is indeed thin walled.

Design variables	Fixed Parameters	Constants
w_a, L_a, h_a, t_f	r, F_t	w_b, l_b, h_b, r, τ_y, σ_y, E, θ_maz

Table 3.2: Variables, parameters and constants of the structural subsystem

The box constraints on the design variables are the following, remebering that w_b = 0.0508 and l_b, h_b = 0.0254.

$$
0.01 – w_a < 0
\tag{3.39}
$$

$$
w_a – \frac{w_b}{2} < 0
\tag{3.40}
$$

$$
\frac{\sqrt{2}}{2} \, r – L_a < 0
\tag{3.41}
$$

$$
0.01 – h_a < 0
\tag{3.42}
$$

$$
h_a – h_b < 0
\tag{3.43}
$$

$$
0.01 – t_f < 0
\tag{3.44}
$$

$$
t_f – w_b < 0
\tag{3.45}
$$

$$
\text{(3.46)}
$$

The box constraints for each variable are purely packaging constraints to ensure that the dimensions of the arms are not larger than the battery housing and that the blades of the propellers do not interfere.

Comments on the nature of the structural optimization problem

Since the objective function is monotonically increasing with the length of the arm (L_a) and monotonically decreasing with the the height (h_a) and with of the arm (w_a) as well as the thickness of the frame. L_a is bounded from below while the other three variables are bounded from above so the problem should be well bounded. The box constraints are necessary to ensure the problem is well bounded. Furthermore, the addition of the theory validity constraints could run the risk of conflicting with the box constraints but this can be avoided by ensuring that the lower bounds for the dimensions of the arms and the frame thickness are at least ten times smaller than the upper bound for L_a. The minimum safety factor constraints are also arbitrary and could conflict with the box constraints i.e. a sufficient safety factor would not be possible with the specified box constraints.

It must be noted that for the purpose of the MDO the structural problem was refor- mulated to return the mass of the drone for a given thrust. The dimensions of the arm (w_a and L_a) are solved for by setting the bending stress safety factor to 1.5 and the deflection of the arm to 0.175 radians. The width and height of the arm were set to be equal.

mulated to return the mass of the drone for a given thrust. The dimensions of the arm (w_a and L_a) are solved for by setting the bending stress safety factor to 1.5 and the deflection of the arm to 0.175 radians. The width and height of the arm were set to be equal.

3.1 Thrust from propellers

Thrust Optimal Design Problem

The thrust subsystem aims to maximise the efficiency of the propellers. It is assumed that there is a motor at the end of each arm with a set of propellers. Simple blade element theory, which decomposes the blade into differential elements, was used to model the efficiency of the propeller [5] along with the help of a vector diagram shown in Figure

3.3. The variables which will be optimized for are the radius of the propeller (from the hub), the diameter of its cross-section, the amount of blades per set of propellers, and the helix angle (ϕ) as shown in figure 3.3. As thrust is related to mass due to the structure and the battery housing, and velocity is optimized for in the battery storage subsystem, these two quantities will have to be fixed in the individual thrust subsystem optimization using a parametric study. The relevant quantities and their equations or values are displayed below.

Figure 3.3: Vector Diagram detailing reactions in the blade [3].

$$
\text{radius} = r
\tag{3.47}
$$

$$
\text{Number of Blades} = B
\tag{3.48}
$$

$$
\text{Diameter of Propeller Cross-Section} = D
\tag{3.49}
$$

$$
\text{Helix Angle} = \phi
\tag{3.50}
$$

$$
\text{Velocity} = V
\tag{3.51}
$$

$$
\text{Mass of Drone} = m_{\text{drone}}
\tag{3.52}
$$

$$
\text{Density of Air} = \rho_{\text{air}} = 1.2442 \; (\text{kg/m}^3)
\tag{3.53}
$$

$$
\text{Density of Propeller} = \rho_{\text{propeller}} = 1260 \; (\text{kg/m}^3)
\tag{3.54}
$$

$$
\gamma = 10^\circ
\tag{3.55}
$$

$$
\text{yield stress} = \sigma_y = 50 \; \text{MPa}
\tag{3.56}
$$

$$
\text{Thrust} = F_t =
\int_{0}^{r}
\frac{0.5 \, \rho_{\text{air}} \, V^2 \, B \, x \, b \, C_L \, \sin(\gamma + \phi)}
{\cos(\gamma) \, \sin^2(\phi)}
\, dr
\; (\text{N})
\tag{3.57}
$$

$$
\text{Torque} = Q =
\int_{0}^{r}
\frac{0.5 \, \rho_{\text{air}} \, V^2 \, B \, x \, r \, b \, C_L \, \cos(\gamma + \phi)}
{\cos(\gamma) \, \sin^2(\phi)}
\, dr
\; (\text{N·m})
\tag{3.58}
$$

$$
\text{angular velocity} = \omega = \frac{V}{\tan(\phi)\, r}
\; (\text{rad/s})
\tag{3.59}
$$

$$
\text{moment of inertia} = I = \frac{\pi \, D^4}{64}
\; (\text{kg·m}^2)
\tag{3.60}
$$

$$
\text{chord } b = (0.3 – 0.7 r + 4.1 r^2 – 9.1 r^3 + 8.9 r^4 – 2.9 r^5)D
\; (\text{m})
\tag{3.61}
$$

$$
\text{Cross Sectional Area} = A_s = \frac{\pi D^2}{4}
\; (\text{m}^2)
\tag{3.62}
$$

$$
\text{Lift} = L = m_{\text{drone}} \cdot 9.81
\; (\text{N})
\tag{3.63}
$$

$$
\text{coefficient of lift} = C_L =
\frac{2 \cdot L}{\rho_{\text{air}} \, V^2 \, \pi \, D \, r}
\tag{3.64}
$$

$$
\text{Wing Tip Velocity} = \omega \cdot r
\tag{3.65}
$$

$$
\text{centrifugal stress} = \sigma_c = \rho_{\text{propeller}} \, V^2
\; (\text{Pa})
\tag{3.66}
$$

$$
\text{bending stress} = \sigma_b = \frac{Q \, D / 2}{I}
\; (\text{Pa})
\tag{3.67}
$$

$$
\text{efficiency} = \eta = \frac{F_t \, V}{\omega \, Q}
\tag{3.68}
$$

The objective function and constraints are shown below.

$$
\min_{\mathbf{x}} \; -\eta(\mathbf{X}) = -\frac{F_t \, V}{\omega \, Q}
$$

$$
r \le 1.45 \; (\text{m})
$$

$$
-r \le -0.08 \; (\text{m})
$$

$$
\text{Wing Tip Velocity:} \quad \omega \, r \le 270 \; (\text{m/s})
$$

$$
B \le 8
$$

$$
-B \le -3
$$

$$
\sigma_b \le 2 \, \sigma_{\text{yield}} = 2 \cdot 50 \times 10^6 \; (\text{Pa})
\tag{3.69}
$$

$$
\sigma_c \le 2 \, \sigma_{\text{yield}} = 2 \cdot 50 \times 10^6 \; (\text{Pa})
$$

$$
D \le 0.07 \; (\text{m})
$$

$$
-D \le -0.001 \; (\text{m})
$$

$$
-2 \, \sqrt{2} \, F_t \le -L \; (\text{N})
$$

$$
\phi \le 75^\circ
$$

$$
-\phi \le 0^\circ
$$

$$
\text{Minimum Thrust: } -4 \, \frac{\sqrt{2}}{2} \, F_t \le -L
$$

Design variables	Fixed Parameters	Constants
r, D, B, ϕ	mdrone, V	ρ_propeller, ρair, γ, σ_y

Table 3.3: Variables, parameters and constants of the thrust subsystem

As this optimization problem seeks to maximise efficiency, a negative sign must be added to turn it into a minimization problem. For each design variable, there are box constraints as the geometrical properties of the propeller (r and D) must not be too big or too small for both manufacturing and application considerations. In addition, the radius must not be too large or else the chord will become a negative value, numerically reduc- ing the thrust and potentially returning a negative thrust. Similar reasoning is used for the box constraints on the number of blades, and the helix angle must be boxed so that its summation with gamma is not bigger than 90 degrees. Other constraints have to do

with making sure that the wing tip velocity is below 80% of the speed of sound, as val- ues higher than that can drastically reduce efficiency and cause shock waves increasing entropy and increasing centrifugal stresses on the blade which must be kept below twice the yield stress of the propeller). The bending stress on the propeller due to the torque it produces must also be below the same threshold. Finally, the thrust produced by all the motors (at an assume angle of 45 degrees relative to the normal to the ground) must be superior to lift needed to get the drone off the ground.

Comments on the nature of the thrust optimization problem

The efficiency is the ratio between the product of the thrust and velocity and the prod- uct of the torque and angular velocity. The equations for thrust and torque share many of the same terms, however the thrust equation has sin(γ + ϕ) while torque has cos(γ + ϕ) which means that the higher the helix angle is, past 45 degrees, the smaller the torque will be and the higher the thrust. In addition, the angular velocity is inversely propor- tional to tan(ϕ), and thus inversely proportional to ϕ. This means that the efficiency is monotonically increasing with the helix angle. Both the thrust and torque have the blade in its integral formulation, which suggests that the amount of blades is arbitrary (and can be removed from the optimization problem). This is convenient as the number of blades is a discrete parameter which causes extra difficulty in solving the problem. Finally, the torque has an extra r in its integral, meaning that the larger it is, the larger the torque is, reducing the efficiency. However, at the same time, the thrust increases as r increases because of the summation of the chord in the integral, which itself is monotonically in- creasing with the diameter of the cross section D. Thus, an increasing diameter increases the efficiency, but also increases the bending stress on the propeller. In addition, the wing tip velocity is increased with an increasing radius but has an upper bound. In essence, the term r is the most interesting variable in the optimization, as it is hard to predict in what direction r will be optimized in. Based off this analysis, the problem should be well bounded, which will be confirmed once the optimization is performed.

Numerical Optimization Results

Overview

This chapter summarizes the relevant numerical optimization results of the various system and subsystem level problems. The results will provide sufficient data to contex- tualize the results not only physically but also mathematically in terms of optimization. The numerical results of the various subsystem optimizations will first be presented fol- lowed by the system level results

Battery storage optimization results

The results for the subsystem level battery storage optimization are presented here. Information on the optimization algorithms, initial guesses, and other characteristic data will be presented. As this optimization contained many local optima, it was decided to run the optimization algorithm thousands of times with random initial guesses of the design parameters between their bounds. The choice results are provided in table 4.1 representing some of the possible local minima. Additionally, it should be noted that all optimization calculations were performed with the MATLAB fmincon code with the SQP algorithm, due to its robust nature and ease of use. The termination criteria (optimality tolerance, function tolerance, and step tolerance) of this fmincon algorithm were set to

the default values of 1e-6, 1e-6, and 1e-10, respectively. Finally, the optimal values across a range of 2000 initial guesses will be presented for a range of c values.

Table 4.1: Summary of the optimization results for the battery storage subsys- tem with c = .5. Initial guess and optimal parameter results are given as X = [tframe,battery, Aconvective area, Vflight,avg, ρframe].

Initial Guess (SI units)	Number of iterations	Optimal function value (thermal resis- tance) ^K W	Local minimum (vari- ables in SI)
[0.0227 0.0010 8.7702 75947.8035]	3	2.8698	[0.0032 0.0032 11.1760 75947.8035]
[0.0071 0.0011 6.7882 36809.9669]	4	2.8703	[0.0032 0.0032 11.1760 36809.9669]
[0.0195 0.0021 7.5887 98963.3374]	3	2.8697	[0.0032 0.0032 11.1760 98963.3374]
[0.0200 0.0025 3.9140 102633.8075]	5	2.8696	[0.0032 0.0032 11.1760 102633.8075]
[0.0262 0.0001 0.4442 148342.8249]	2	2.8695	[0.0032 0.0032 11.1760 139972.0000]

The optimal local minimum given in table 4.1 directly matches the intuitive expecta- tions behind the model thereby lending credence to the result. The optimal design param- eters include the smallest frame thickness, largest area open to convection, highest wind speed, and highest density (and thus highest thermal conductivity) possible. Again, this highest thermal conductivity value was set equal to the thermal conductivity of silver, which is quite the interesting result.

Figure 4.1: The effect of the c parametric parameter on the thermal resistance objective function.

The effect of the c parameter was quite interesting and may be understood best by viewing Table 4.2. Additionally, this was plotted in Figure 4.1 to better elucidate the re- sults. One may clearly see that as c increases the objective function decreases significantly, with nearly a 20% difference over the analyzed range. This is quite logical, as c essentially represents the degree to which the convective heat transfer resistance is decreased for the conduction branch. The interesting bifurcation in Table 4.2 occurs at c = 1. This is in fact a fascinating result and illuminates much about the physics. At this point, the conductive branch becomes more optimal and thus there was no need to take any area for purely convective area due to the low conductive resistance. This value at 1 is likely due to the symmetry in the problem and is a truly interesting physical bifurcation result.

Table 4.2: Summary of the optimization results for the battery storage subsystem with for a varying c value (benchmark at c = .5 presented in Table 4.1). Optimal parameter results are given as X = [tframe,battery, Aconvective area, Vflight,avg, ρframe].

c (-)	Optimal function value (ther- mal resistance) ^K W	Local minimum (variables in SI)
.8	2.3915	[0.0032 139972]	0.0032	11.1760
.9	2.2658	[0.0032 139972]	0.0032	11.1760
1	2.1526	[0.0032 139972]	0.0032	11.1760
1.025	2.1007	[0.0032 139972]	0.0000	11.1760
1.1	1.9575	[0.0032 139972]	0.0000	11.1760

Structural Optimization Results

The optimization was done on python using the scipy.optimize,minimize function. The algorithm chosen was ’SLSQP’ due to its ability to handle non-linear constraints. The solver was initiated at 10 randomly generated guesses within the design space. Three of the results are tabulated in Table 4.13

Initial Guess	Function value	Minimizer	Number of iterations
[0.01429, 9.61027, 0.04832, 0.01703]	0.00915	[0.0254 0.56568 0.0254 0.00254]	11
[0.02311, 6.93506, 0.24570, 0.02360]	0.01358	[0.01739 0.56568 0.0254 0.00250]	7
[0.02167, 6.37706, 0.02642, 0.02523]	0.00915	[0.0254 0.56568 0.0254 0.002539]	8

Table 4.3: Structural Subsystem Optimization results

These results are obtained with the thrust being set to 20N and a propeller radius of 0.8 m. For the first and second rows of Table 4.3, the active constraints are the lower bounds for the width (= 0.0254) of the arm and the height of the arm (= 0.0254). The Length of the arm is at its lower bound (= 0.05658) and the thickness of the frame is bounded by the first thin wall (10 − wb/t < 0) constraint. The second result converged to a local minimum where the length of the arm and the height of the arm lie on their upper bounds. To further investigate the activity of this constraint it is relaxed. Subsequently, the second thin walled constraint becomes active (10 − hb/t < 0) since h_b = w_b ∗ 2 and a new optimizer is obtained as tabulated below. The results concur with our intuition since we would expect all of the dimensions to be maximized in order to minimize stress.

Relaxed con- straint	Initial Guess	New Function Value	New Optimizer	Number of itera- tions
Thin Wall 1	[0.01247, 7.86517, 0.40609, 0.01613]	0.00400	[0.0254 0.56568543 0.0254 0.00508]	11
Thin Wall 1, Thin Wall 2	[0.01180, 7.42088, 0.04689, 0.02393]	0.00035	[0.0254 0.56568543 0.0254 0.0254 ]	4

Thrust Optimization Results

Similarly to the structural optimization, the optimization was done using Sequential Least Squares Programming (SLSQP) which is a gradient based-method. This method was chosen in order to handle non-linear constraints, which when using the scipy.optimize.minimzize function eliminates many options apart from COBYLA and trust-constr. The maximum iterations allowed were 10,000 and the default tolerance of 1e⁻6 was used. However, only SLSQP was able to converge to a result which lead to its use. In table 4.5 are the results for three guesses each with the same fixed value for mass and velocity.

Initial Guess	Function value	Minimizer (D,r,B,ϕ)	Number of iterations
[0.0345, 0.360, 2.5, 37.5]	-75.5103	[ 0.0345, 0.8, 8, 75]	8
[0.8, 0.07, 8, 75]	-75.5103	[0.07, 0.8, 8, 75]	2
[0.1, 0.01, 3, 20]	-75.5103	[0.0666073, 0.8,4, 75]	14

We see that different guesses give different values for D and B, however, this has no effect on the final efficiency of the propellers. In addition, the number of evaluations depends heavily on the initial guess, with the second requiring only two iterations (due to the fact that the guess is already very close to the final optimizer. We can see that for the second and third guess, D, r, and ϕ are all maximised, however for the first guess the value of D is closer to the center of its bounds. As the mass of the drone, and its velocity are variables solved by the other subsystems, a parametric study for both was conducted using the second guess from Table 4.5.

Table 4.6: Parametric Study of Mass

Table 4.7: Parametric Study of Velocity

Somewhat surprisingly, as can be seen in Table 4.6 and Table 4.7, varying either mass or velocity does not affect the final convergence of the design variables. However, we do see in Table 4.6 that changing the mass of the system does change the amount of thrust which makes sense. The parametric study was conducted within the bounds of the expected masses and velocities, which may be the reason why changing these fixed parameters does not affect the optimizer, as none of them are extreme enough to cause the nonlinear constraints (bending stress, centrifugal stress, wingtip velocity) to become active.

Table 4.8: Active Constraints

Table 4.8 presents which constraints are active, which is done by suppressing each constraint and reruning the problem. We find that the active constraints are the upper bound of the radius, the upper bound of the number of blades, and the upper bound for the helix angle.

Table 4.9: Thrust Subsystem Relaxation

Constraint Removed	Function value	Minimizer (D,r,B,ϕ)
Radius Upper Bound	-42676341.78845887	-2.90314460e-01, 1.41348474e+06, 4.62196811e+00, 7.55467829e+01
Blade Upper Bound	-75.51026412070432	7.00000000e-02, 8.00000000e-01, 8.00000024e+00, 7.50000000e+01
Helix Upper Bound	-1391693.9481687378	6.99902519e-02, 7.99353207e-01, 8.00002623e+00, 7.99995868e+01

The results of the relaxation are presented in table 4.9, however, when removing the first constraint, the optimizer converges to an infeasible solution (negative Diameter), even though the constraint should still be active (and the optimizer does not say that it failed) meaning that the first row can be ignored. The second constraint’s relaxation does not do much to actually improve the function value nor does it optimize to different points. However, the third row, the relaxation of the helix upper boundary, is the most interesting, as it actually does cause the value of the helix angle to increase almost up to 80 degrees. This causes the sin term in the calculation of the thrust to be maximised, and the cos term in the torque calculation to be minimized. However, this value cannot actually reach 80 degrees as or else the torque (the denominator in the efficiency function) would go to 0, causing the function value to explode, as we can see here with its very high value. As well, to be noted, efficiency is being calculated, however, values are above -1 when physically that would be impossible. This is most likely due to simplifications in the elaboration of the model.

Finally, trying to solve the optimization problem without the blade as a design variable but instead as a fixed parameter with a value of 8 was attempted. This simplified the problem as the number of blades can only take on discrete values, while the SLSQP solver returns continuous values. In order to remedy this, the problem would be solved, the number of blades would be rounded then, then every condition would be checked to see

if they were satisfied. If this was not the case, then the number of blades was increased until all the conditions were satisfied or the maximum number of blades was used. By removing B as a design variable, this loop could be entirely avoided. When this was attempted, both COBYLA and SLSQP were able to converge, providing slightly different values to those presented previously.

Table 4.10: COBYLA Parametric Study of Mass

Table 4.11: SLSQP Parametric Study of Velocity

Table 4.10 presents interesting results as we can see that the diameter increases as the mass increases, the helix angle is not maxed out but converges to 69, and the efficiency is only -24. This shows that COBYLA does not converge to the best efficiency (-75) showing that SLSQP was a good choice for an optimizer. In addition, Table 4.11 shows that the SLSQP optimizer converged not to a maximum diameter but instead to a value of 0.05.

Table 4.12: COBYLA Active Constraints

If we look at table 4.12, we see that when using COBYLA the bending stress constraint is activated, which was not the case previously.

The final optimizer chosen (when using the number of blades as a design variable) are presented in table 4.13.

Table 4.13: Thrust Subsytem Final Optimizer

D	r	B	ϕ
0.07	0.8	8	75

System level optimization results

Our whole system comprises three subsystems: thermal, structural, and thrust. The thermal subsystem seeks to reduce the heat transfer around the battery, the structural subsystem seeks to minimize the stresses that the drone will experience, and the thrust subsystem seeks to maximize the efficiency of the propellers. These subsystems are not independent however. The structural subsystem utilizes the density of the drone and the thickness of the battery frame from the thermal subsystem as box constraints, while the thrust subsystem needs the flight speed (as well from the thermal subsystem) in its calculations. However, the thermal subsystem does not need any design variables from the other two subsystems to conduct its optimization, meaning that t_frame, ρ_frame, and V are shared variables. In addition, the thrust system needs to calculate thrust, which is dependent on the mass of the drone, which itself is a function of w_arm and l_arm design variables of the structure subsystem. However, the structure subsystem calculation of stress involves thrust, which is dependent on the design variables of the thrust subsystem (D,r,ϕ). Thus, these variables are coupled, and the problem cannot be solved all-in-one due to this. This system level optimization problem is now a Multi-Disciplinary Opti-

mization problem and must be solved using different methods than the subsystems.

In order to accomplish this optimization, the structural subsystem and the thrust subsys- tem were simplified. The thrust and structural subsystems were simplified by removing the number of blades and height of the arm respectively as design variables, and instead treated them as fixed parameters with a value of 8 (as was done in the previous section) and the height of the arm equal to the width of the arm (square beam). In addition, only the box constraints were retained. Finally, the biggest change was that instead of each subsystem having its own objective function to minimize as was done when treated in- dividually, now, only the structural subsystem has an objective function which seeks to minimize mass. All the subsystems also have penalty functions which make sure that shared and coupled variables have the same values.

This setup was achieved within a Non Hierarchical Analytical Target Cascading (NHATC) architecture and implemented using the NoHiMDO package on Matlab. The NHATC ar- chitecture was suited to this particular Multi-Disciplinary Optimization problem because of the presence of coupling and the non-hierarchical nature of the problem as depicted in Figure 2.1. The problem definition is presented in table 4.14, showing whether variables are coupling variables, the link each variable has to other variables, and each variables lower and upper bounds.

The results of the MDO are presented in table 4.15.

The first thing to be noted is that since there are many changes involving the sub- systems between their independent handling and their coupled handling in the MDO, it is harder to compare the values of the design variables when solved in both instances. However, because the previous objective functions were all thrown out in favor of just optimizing for mass, it provides a great opportunity to observe the trade-offs between maximising efficiency, reducing heat transfer, reducing stresses and minimizing weight. For instance, the density of the frame of the drone is minimized in the MDO but maxi- mized when solved independently. In addition, the radius of the propellers is brought to its lower as opposed to its upper bound, and the helix angle is made very small, while the diameter does not change much. These changes contribute to reducing the thrust the mo- tors produce, allowing a smaller mass. Finally, somewhat counter-intuitively, theOur whole system comprises three subsystems: thermal, structural, and thrust. The thermal subsystem seeks to reduce the heat transfer around the battery, the structural subsystem seeks to minimize the stresses that the drone will experience, and the thrust subsystem seeks to maximize the efficiency of the propellers. These subsystems are not independent however. The structural subsystem utilizes the density of the drone and the thickness of the battery frame from the thermal subsystem as box constraints, while the thrust subsystem needs the flight speed (as well from the thermal subsystem) in its calculations. However, the thermal subsystem does not need any design variables from the other two subsystems to conduct its optimization, meaning that t_frame, ρ_frame, and V are shared variables. In addition, the thrust system needs to calculate thrust, which is dependent on the mass of the drone, which itself is a function of w_arm and l_arm design variables of the structure subsystem. However, the structure subsystem calculation of stress involves thrust, which is dependent on the design variables of the thrust subsystem (D,r,ϕ). Thus, these variables are coupled, and the problem cannot be solved all-in-one due to this. This system level optimization problem is now a Multi-Disciplinary Opti-

mization problem and must be solved using different methods than the subsystems.

The results of the MDO are presented in table 4.15.

The final objective value for the mass is 0.0302 kgs, which was achieved in 150 itera- tions, and a final maximum inconsistency of the order of 10⁻⁷. The plots of the optimiza- tion can be seen in figure 4.2.

Table 4.14: MDO Problem Definition

Table 4.15: MDO Problem Definition

Figure 4.2: Plot of Optimization of MDO

Conclusions

In this project, we conducted three individual subsystem optimizations and one sys- tem level optimization. For the individual subsystem problems, box constraints along with non-linear constraints were used to bound the design space for the design variables, and fixed parameters utilized in the optimization were varied in parametric studies to evaluate their effect on the optimizers found. The optimization algorithms used were SQP for the thermal subsystem, and SLSQP for the structural and thrust subsystems. These optimizations were able to converge to minimizers very quickly (almost always under 10 iterations).

In order to solve the system level optimization, several modifications had to be applied to each subsystem. First of all, the structural and thrust subsystems had to be simplified by removing design variables (height of the arm and number of blades). Secondly, several constraints used in the individual subsystem optimizations were relaxed. Finally, instead of each subsystem having an objective function, only the structural subsystem provided an objective function which optimized for mass. This Multi-Disciplinary Optimization problem was solved by setting up the subsystems in a Non-Hierarchical Analytical Target Cascading architecture, allowing the problem to be solved using the NoHiMDO matlab package.

The change from individual objective functions to a singular system level objective of re- ducing the mass of the system provided an opportunity to witness the tradeoffs between

the subsystem level objectives and system level objective as the final design variables in both cases differed quite heavily.

Bibliography

Mashable. Capture your summer adventures from the sky with a mini dual camera drone. 2022.
Ahmet Tandiroglu. Temperature-dependent thermal conductivity of high strength lightweight raw perlite aggregate concrete. International Journal of Thermophysics, 31:1195–1211, 06 2010.
Fred Ernest Weick. Aircraft propeller design. McGraw-Hill Book Company, Incorpo- rated, 1930.
Chunjing Lin, Fang Wang, Bin Fan, Shan Ren, Yuemeng Zhang, Liqiong Han, Shiqiang Liu, and Sichuan Xu. Comparative study on the heat generation behav- ior of lithium-ion batteries with different cathode materials using accelerating rate calorimetry. Energy Procedia, 142:3369–3374, 2017. Proceedings of the 9th International Conference on Applied Energy.
Chiranjivi Dahal, Hari Bahadur Dura, and Laxman Poudel. Design and analysis of propeller for high-altitude search and rescue unmanned aerial vehicle. International Journal of Aerospace Engineering, 2021:1–13, 2021.