New analytical results on the study of aircraft performance with velocity dependent forces

Abstract

Most aeronautical accidents happen during takeoff and landing. The main objective when studying those phases of the flight mission is to answer a seemingly simple questions: can the airplane safely takeoff and land on the stipulated runway dimensions with the intended weight? The main objective of the present paper is to obtain new analytical answer to those questions, for fixed wing airplanes. To our present knowledge such a solution, with the degree of generalisation proposed here, is new in the literature. Regarding previous studies, first a new power unit traction equation is employed to explicitly consider the influence of air density, angular velocity and diameter of the propeller. Then a new method for calculating the maximum weight is proposed. Next, the use of breaks is modeled and analysed and an equation to calculate the static gliding wind velocity is proposed. Finally, a toolbox created to perform the calculations is described. A thorough analysis of the influence of the airplane design parameters on the behavior of the motion equations is made, with special attention to the use of brakes. Numerical results are successfully compared with experimental data from two models of a commercial airplane, the Cessna 172 Skyhawk models N and S, and four UAV prototypes. The methodology employed uses simple laws of classical mechanics allied to basic calculus and is easy to understand by first year students of physics, engineering or mathematics.

Keywords
Flight mechanics; aircraft performance; takeoff and landing distances; maximum weight

1. Introduction

The application of the principles of classical physics to mechanical systems marked the beginning of what is currently understood as modern engineering, dating back to the the scientific revolution on the 16^th and 17^th centuries. In the field of aeronautical engineering the use of classical mechanics can be traced back to the aviation pioneers, around the start of the 20^th century, when machines heavier than the air, capable of flying by its own means with fixed wings started to be studied.

As soon as the aeronautical designers understood that airplanes should have propeller systems that were independent of their lift generating systems¹ (different from bird and helicopters), the need of precise calculations of a number of parameters related to its performance appeared. After all, flying is an activity that leaves very little space for trial and error, as so many early accidents demonstrated.

From the beginning, two crucial phases of the flight were carefully considered: the takeoff and the landing. It is not to say that the other phases were not important, but most accidents happen around or during takeoff and landing. According to a statistical survey by the Boeing company, 49% of all fatal accidents happen during the final descent and landing, while 14% of all fatal accidents happen during the takeoff and the initial climb.

Thus, the main objective when studying takeoff and landing is to answer the seemingly simple questions: can the airplane takeoff and land on the stipulated runway dimensions with the intended weight? Obviously all texts in aeronautical engineering offer an answer to this question. Nevertheless, it still offers some interesting aspects to be studied.

The question first received attention from the authors in 2014, as they became involved with the design of a small unmanned aerial vehicle (UAV) to take part in that year’s SAE AeroDesign competition. The contest challenges students from all around the country to conceive, design and build cargo UAVs to perform specific missions on a tournament held every spring in São Paulo, Brazil.

This kind of project typically begins with flight mechanics calculations in general, and specifically with an analysis of the takeoff and landing phases of the flight mission. The idea is to calculate the necessary runway length to takeoff and to land safely given the total weight of the airplane, or conversely to calculate the maximum weight, given the runway length available for each phase.

In any case, calculations starts from Newton’s second law of motion, and the distances are calculated as the result of an integral equation involving these forces, which all depend on velocity, except for the weight. Traditional aeronautical literature ([¹] based on the results of [²]) simplifies such integration by using average values of the forces. This simplification works well in practice but has several known limitations [³], the most noteworthy being the fact that the approximate integration does not yield an equation of motion, but only the total takeoff and landing distances.

Other approach relies on statistics ([⁴], for example), correlating experimental data of takeoff distance and landing with weight, wing area, engine power and lift coefficient for different models of commercial aircraft. Of course, the problem is simple enough to be tackled numerically ([¹] and [⁴], for example) but again the dynamic details of the maneuver are not explicit. Only one author [⁵] presents an analytical solution for the problem, but under a rather particular condition that, more often than not, fails to model a real aircraft.

With this motivation in mind, a 2015 study [⁶] proposed a closed form analytical solution to the problem. The study was originally intended to provide students with a simple analytic model to be embedded in their optimization routines at the beginning of each year’s project for the SAE AeroDesign competition. In the years that followed, the paper achieved almost 14,000 downloads, suggesting a large interest on the part of engineering and physics undergraduates and, probably, of the other AeroDesign teams. Therefore, came the idea of revisiting this work to include as many new results as possible, to correct one unnoticed mistake and to make the results available to a greater audience, by publishing it in English.

Therefore the main objective of the present paper is to obtain new analytical solutions for the takeoff and landing phases of the fixed wing airplanes flight, comparing the predicted results to measurements takes for four UAVs and two distinct models of a commercial size airplanes. To our present knowledge, such a solution is new in the literature. The main additions to the results presented previously in [⁶] are:

• A new power unit traction equation that explicitly considers the influence of air density, angular velocity and diameter of the propeller. In the new model, the traction equation can be used for any power regime, including reverse thrust, and not only for maximum power and idle.²;

• Generalization of the equations to calculate the takeoff and the landing distances to turn them into the equations of motion;

• A new method for calculating the maximum weight;

• Analysis of the average solution due to [²] and comparison with the present results;

• An equations describing the use of breaks during landing for different kinds of braking profiles;

• An equation to calculate the static gliding wind velocity and the associated analysis of its relevance during the landing;

• Updated values of the rolling friction coefficients based on recent measurements [⁸];

• Application of the results to three more UAV prototypes;

• A much broader analysis of the influence of the airplane design parameters on the behavior of the motion equations, with special attention to the use of brakes;

• Description of a toolbox created to perform the calculations³.

2. Takeoff

In this section we present the study of the takeoff phase of the flight. Some material from [⁶] is included for the sake of clarity.

Consider an aircraft accelerating on a horizontal runaway, with velocity V in relation to the runaway, Figure 1. For a wind velocity U₀ the aircraft velocity relative to the Earth is

(1) $$V=\upsilon +U_0$$

where $\upsilon$ is the velocity of the aircraft relative to the air. Values of $U_0<0$ are called headwinds and $U_0>0$ are called tailwinds.

Figure 1
Takeoff maneuver.

Figure 2 illustrates the forces acting on an aircraft during takeoff, where T is the power unit traction (also called thrust), D is the total drag force, L is the lift force, W is the aircraft weight and R is the rolling friction between all tires of the landing gear and the runway. All forces except the weight depend on the aircraft velocity relative to the air, $\upsilon$.

Figure 2
Forces acting on an airplane during takeoff or landing.

The takeoff maneuver starts from rest, or nearly so, and the aircraft must accelerate to a velocity V_to where the lift force is larger than the aircraft weight. Thus, taking off should be simply a matter of having the necessary runway length. However, as the dissipative forces, D and R, increase with velocity while the traction decreases, this condition may not be satisfied. In symbols, the condition for takeoff is

(2) $$T\ge D+R\mathit{}\text{for}\{\begin{array}{cc}V\in [0,V_{to}]\hfill & \\ L_{max}>W\hfill & \end{array}$$

However, the power unit may not be able to accelerate the aircraft to the necessary velocity, however long the runway might be and the airplane may not be able to take off.

To understand the situation a little better, the reader is referred to Figure 3. It shows the relation between the lift and drag forces generated by the airfoil NACA43012A and the angle of attack. The relation is shown through the lift and drag coefficients, C_L and C_D respectively, to be formally defined soon. The lift coefficient is seen to increase almost linearly with $\alpha$ up to a maximum and then to drop sharply. This point is known as the stall point, where a rapid loss of lift occurs due to the separation of the boundary layer over most of the airfoil upper surface. Figure 3 also shows that the minimum value of C_D occurs at $\alpha \approx 0$.

Figure 3
Lift and drag coefficients for a conventional airfoil. Reproduced from [¹⁰].

In a real takeoff maneuver first the airplane must accelerate from rest to a velocity where the maximum lift force L_max, occurring at $\approx {15}^{\circ }$ in Figure 3, equals the weight. To reduce the takeoff distance, the ground roll must be done with minimum drag force, using $\alpha \approx 0$.⁴ At this point, however, the aircraft do not yet leave the ground, because at $\alpha \approx 0$ the lift is still insufficient. Shortly after this moment (for safety reasons) the pilot rotates the aircraft (Figure 1) pulling its nose up to increase the lift to its maximum, making the airplane finally takes off.

The whole maneuver is performed against the wind, so that the relative velocity between the air and the aircraft increases, thus increasing the lift force for a given velocity in relation to the ground. Of course this also increases the drag force, but this is not a problem as long as the power necessary is available.

During the ground roll, Newton’s second law of motion with the forces shown in Figure 2 in the (considered) inertial frame attached to the runway is

(3) $$T-D-R=m\frac{dV}{dt}$$

for $0\le V\le V_{to}$, where t is time.

By the chain rule, $dV/dt=(dV/ds)(ds/dt)=V(dV/ds)$, where s is the displacement. Integrating Eqn. (3) between V=0 and any given velocity $V\in [0,V_{to}]$ yields

(4) $$s=\frac{W}{g}{\int }_0^V\frac{VdV}{T-D-R}$$

where W is the total weight of the aircraft (payload plus aircraft weight), considered constant, as the only weight that varies during take-off is the fuel’s, but the decrease is not significant.

The following subsections present details about the three forces in Eq. (4).

2.1. Propeller traction

Many authors, as [¹] and the previous study of [⁶], use a simple quadratic model to describe the traction, i.e., $T=a_1{\upsilon }^2+a_2\upsilon +a_3$. Written in this form, we recognize that the coefficients a₁, a₂ and a₃ must depend on a number of parameters, as T also does: the propeller diameter⁵ D, the angular velocity n, and the air density $\rho$. To explicitly access this dependence, we can use the expression [¹¹]

(5) $$T=\rho n^2{D}^4{C}_T$$

where C_T is the thrust coefficient.

This result comes from dimensional analysis and is very popular in the propeller design area. Comparing the previous result to the quadratic equation, it follows that C_T must also depend quadratically on V, i.e.,

(6) $$C_T=a{\upsilon }^2+b\upsilon +C_{T0}$$

as shown for example by [¹²]. Subtituting this relation back in equation (6) yields

(7) $$T=\rho n^2{D}^4(a{\upsilon }^2+b\upsilon +C_{T0})$$

The main advantage of this equation in relation to the $T=a_1{\upsilon }^2+a_2\upsilon +a_3$ form is the fact that D, $\rho$ and specially n are represented explicitly.

Using the original quadratic form it is necessary to have values of a₁, a₂ and a₃ for the takeoff, the landing and any intermediate power regime necessary to the flight mission. Using equation (5) on the other hand, it is only necessary to know a, b and $C_{T0}$, which can be determined by just one experiment. The idle regime can be represented by n=0, full power by n=1 and any other intermediate power regime by the adequate value of $0<n<1$. Even the use of reverse thrust can be represented by a negative⁶ value of $(n^2)$.

Analytical solutions found in the literature [⁵] solve the integral on equation (4) considering a₂ = 0 (and, thus, b = 0) to simplify the calculations. This simplification is not necessary because as shown before by [⁶] the integral have a simple analytical closed form solution. Furthermore, experimental data [¹³] and our experiments confirms that $b\ne 0$ for at least some commercial size engines and for some internal combustion and electric UAVs motors (Figure 4).

Figure 4
Normalized traction data for typical UAVs and commercial size engines.

Our experiments also shows that the parabola is concave downward ($a<0$) for the UAV engines tested. For commercial aircraft piston engines, the parabola generally presents a slight concavity upwards ($a>0$), as Figure 4 also exemplifies. In all cases, the traction decreases monotonically with increasing velocity, thus presenting a negative derivative at V=0 and, consequently, $b<0$.

2.2. Rolling friction

According to the modern understanding of tribology, there are two main causes of the energy dissipation which give rise to the resistance in the rolling contact between surfaces [¹⁴]. The first is associated with the work involved in the making and the breaking of adhesion bonds formed in the points of contact. The second is the cyclic deformation of the bodies.

For cases of free rolling (no slipping), where the rolling body deforms but the surface is rigid, R can be approximately modeled by [¹⁵]

(8) $$R={\mu }_r(W-L)$$

where ${\mu }_r$ is the rolling friction coefficient and $(W-L)$ is the normal load, resulting from the difference between the total weight of the aircraft and the lift force during takeoff.

It is known that rolling friction depends on the velocity, the torque applied to the wheels and their diameter [¹⁵]. To simplify the problem, however, only the dependency on the surface materials will be considered here.

2.3. Drag and lift

The drag and lift forces are traditionally modeled as [¹⁶]

(9) $$D=\frac{1}{2}\rho {\upsilon }^2{C}_{Dtot}S$$

(10) $$L=\frac{1}{2}\rho {\upsilon }^2{C}_{L}S$$

where $C_{Dtot}$ is the total drag coefficient of the airplane, C_L is the lift coefficient of the wing, and S is the vertical projected wing area.

Equations (9) and (10) also come from dimensional analysis. The coefficients $C_{Dtot}$ and C_L depend on the non dimensional geometrical parameters of the body, collectively denoted by G here, the angle of attack, the Mach number and the Reynolds number, $\mathrm{Re}=\rho \upsilon c/\mu$, where $\mu$ is the dynamic viscosity of the air and c is the aerodynamic mean chord of the airfoil. According to dimensional analysis, thus, we have

(11) $$C_n=C_n(\alpha ,G,\mathrm{Ma},\mathrm{Re})$$

for $n=D,L$ or any other force of interest.

Therefore, for a given geometry and angle of attack in subsonic flows, D and L depend only on $\mathrm{Re}$, and consequently only on $\upsilon$, explicitly through $({\upsilon }^2)$ and implicitly, through $C_{Dtot}$ and C_L.

The dependence of C_L with $\mathrm{Re}$, and thus with $\upsilon$, is simple, as seen in Figure 3: C_L is practically independent of $\mathrm{Re}$ in the linear part of the curve. The dependence of $C_{Dtot}$ with $\upsilon$ is a bit more pronounced as seen if Figure 3, but will be ignored here to simplify the analysis, as most aeronautical texts do.

As the lift force is practically all generated by the wing, its projected area appears naturally as the reference surface, S in Eqn. (10). The total drag, however, has several components due to the different parts of the aircraft and the interaction between them, and each one should be associated, in principle, with its own reference area. Nevertheless, for simplicity, it is common practice to refer the total drag force also to S.

2.4. Equations of motion

Substitution eqs. (7)–(10) in equation(4) yields

(12) $$s=\frac{W}{g}{\int }_0^V\frac{VdV}{A_{to}{V}^2+B_{to}V+C_{to}}$$

where

(13) $$A_{to}=\rho \left[a{n}^2{D}^4+\frac{S}{2}\left({\mu }_r{C}_L-C_{Dtot}\right)\right]$$

(14) $$B_{to}=b\rho n^2{D}^4-2A_{to}{U}_0$$

(15) $$C_{to}=\left(C_{T0}-b{U}_0\right)\rho n^2{D}^4-{\mu }_{r}W+A_{to}{U}_0^2$$

are independent of V. The result of the integration is shown in Eqs. (16) and (17). These are the equations of motion for the airplane during takeoff for the cases where $4A_{to}{C}_{to}-B_{to}^2>0$ and $4A_{to}{C}_{to}-B_{to}^2<0$ respectively and constitute one of the main results of the present work.

(16) $$s=\frac{W}{2g{A}_{to}}\left[\ln \right|\frac{A_{to}{V}^2+B_{to}V+C_{to}}{C_{to}}|$$ $$-\frac{2B_{to}}{\sqrt{4A_{to}{C}_{to}-B_{to}^2}}(\arctan \frac{2A_{to}V+B_{to}}{\sqrt{4A_{to}{C}_{to}-B_{to}^2}}$$ $$-\arctan \frac{B_{to}}{\sqrt{4A_{to}{C}_{to}-B_{to}^2}})]$$

(17) $$s=\frac{W}{2g{A}_{to}}\left[\ln \right|\frac{A_{to}{V}^2+B_{to}V+C_{to}}{C_{to}}|-\frac{B_{to}}{\sqrt{B_{to}^2-4A_{to}{C}_{to}}}$$ $$\times (\ln \frac{2A_{to}V+B_{to}+\sqrt{B_{to}^2-4A_{to}{C}_{to}}}{2A_{to}V+B_{to}-\sqrt{B_{to}^2-4A_{to}{C}_{to}}}$$ $$-\ln \frac{B_{to}+\sqrt{B_{to}^2-4A_{to}{C}_{to}}}{B_{to}-\sqrt{B_{to}^2-4A_{to}{C}_{to}}})]$$

for $4A_{to}{C}_{to}-B_{to}^2>0$ and $B_{to}^2-4A_{to}{C}_{to}>0$ respectively.

To obtain the takeoff distance from the preceding equations, it is only necessary to substitute the value of the takeoff velocity in it. The velocity adopted for this is the minimum velocity at which an aircraft can fly, known as stall velocity, with a security margin. The stall velocity can be obtained by equating the airplane weight to the maximum lift force

(18) $$W=\frac{1}{2}\rho {\upsilon }_{stall}^2{C}_{Lmax}S$$

It follows that the relative stall velocity is ${\upsilon }_{stall}={(2W/\rho C_{Lmax}S)}^{1/2}$ and, thus, that

(19) $$V_{to}={\text{SF}}_{to}\sqrt{\frac{2W}{\rho C_{Lmax}S}}+U_0$$

where ${\text{SF}}_{to}$ is the safety margin for takeoff. According to accepted safety standards ([¹⁷] for example), ${\text{SF}}_{to}=1.2$ for takeoff and ${\text{SF}}_{\mathrm{\ell }}=1.3$ for landing.

Mathematically, Eqns. (16) and (17) are not defined for $4A_{to}{C}_{to}=B_{to}^2$. Physically, if $4A_{to}{C}_{to}=B_{to}^2$, then equation $A_{to}{V}^2+B_{to}V+C_{to}=0$ has two identical real roots, $-B_{to}/2A_{to}$. In this case, since $A_{to}{V}^2+B_{to}V+C_{to}=T-D-R=\sum F_x$ is the sum of the horizontal forces in the takeoff direction (by Eqns. (3) and (12)) we have $\sum F_x=0$ at some point during the ground roll. Obviously, the airplane is not able to take off if $\sum F_x=0$ except at the exact moment of rolling. This implies that $A_{to}{V}^2+B_{to}V+C_{to}=0$ cannot have roots in the problem domain, consistently with the natural logarithm arguments in Eqns. (16) and (17).

Equations (16) and (17) can be simplified if the headwind is neglected and traction is given by $T=a_1{V}^2+a_3$. In this case, we would have B=0 and Eqns. (16) or (17) would reduce to the simpler logarithmic form presented by [⁵].

2.5. Averaged motion equation

As mentioned before, the classical analytical solution for the problem under study unnecessarily considers all forces to be constant when integrating equation (12), a simplification known to work well in practical terms. We revise this result here to understand why.

Shevell [²] proposes that an average total force can represent the variable forces at the denominator of equation (12). According to the author, calculating the forces at $V=\kappa V_{to}$ yields a correct result for $\kappa =0.7$.

To investigate this claim, we substitute $V=\kappa V_{to}$ at the denominator of equation (12) to obtain the average forces, integrate in V and evaluate the result at $V=V_{to}$. The result is

(20) $$s_{to}=\frac{W{V}_{to}^2/2g}{{\kappa }^2{A}_{to}{V}_{to}^2+\kappa B_{to}{V}_{to}+C_{to}}$$

This result can be equated to Eqs. (16) or (17) at $V=V_{to}$, yielding a second degree equation in $\kappa$

(21) $${\kappa }^2(A_{to}{V}_{to}^2)+\kappa (B_{to}{V}_{to})+C_{to}-\frac{W{V}_{to}^2}{2g{s}_{to}}=0$$

which solution is

(22) $$\kappa =\frac{-B_{to}{V}_{to}\pm \sqrt{{(B{V}_{to})}^2-4A_{to}{V}_{to}^2{C}_{\kappa }}}{2A_{to}{V}_{to}^2}$$

where $C_{\kappa }=C_{to}-W{V}_{to}^2/2g{s}_{to}$.

Equation (22) suggests that $\kappa$ is probably not constant, for it depends on every input variable of the problem. The sensitivity of $\kappa$ to parameters as a, b, $C_{T0}$ and n will be evaluated in a later section.

2.6. Maximum weight and MTOW

In many situations regarding the operations of UAVs, the operators may prefer to calculate the maximum weight for the airplane to take off given the available runaway length, rather than the takeoff distance given its weight.

The maximum weight, hereafter denoted by $W_{max}$, is the weight for which the airplane reaches the takeoff velocity just at the end of the runaway. For any $W>W_{max}$ a longer runaway must be available, provided the airplane have the power to keep accelerating, i.e, that $T>D+R$. Thus, calculating $W_{max}$ is only a matter of inverting the equation of motion, i.e., of finding a relation

(23) $$W_{max}=W_{max}(s_{run},V_{to})$$

where $s=s_{run}$ is the available length of the runaway.

Unfortunately, this task is easier said than done, as $V_{to}$ and $C_{to}$ also depend on W. Therefore, we could find no analytical way to invert those equations.

An approximate method to calculate $W_{max}$ was proposed in [⁶] but our comparisons of its results with a larger set of data shows that they are not always accurate, because they assume that $V_{to}$ occur at $s_{run}$and when the sum of the horizontal forces acting on the airplane are zero. This is simply not always the case.

Assuming that Eqs. (16) and (17) are bijective and can, thus, be inverted, we may calculate $W_{max}$ numerically. The task, although analytically impossible (we believe) is not a problem for a simple numerical code. We simply give the code a tentative value of W, then calculate $V_{to}$ and $C_{to}$ (that depends on W), after that calculate $A_{to}$ and $B_{to}$ (that do not), and finally obtain $s_{to}$ from Eqs. (16) or (17). If $s_{to}=s_{run}$, within some prescribed tolerance, then $W=W_{max}$.

Although this is a analytical paper, we chose to present this method here, both because of its practical importance and to point out the mistake in the analytical method presented by [⁶].

At this point, it is worthy to explain the difference between $W_{max}$ and the maximum takeoff weight (MTOW), used in commercial and military aviation. The MTOW is defined as the maximum weight the airplane is certified to take off with, due to structural and “other limits”, and does not take the runway length into consideration. Being a safety threshold, it gives priority regard to structural limitations. The value of $W_{max}$ calculated here nevertheless assumes the airplane is capable of resisting to all loads and simply informs the maximum weight it can take off with, for a given runaway length.

3. Landing

The same forces acting on the airplane during the takeoff ground roll, act during the landing roll with the exception of the braking force. The engine traction is sometimes considered null during the landing roll, since the pilot reduces the engine to idle at touch down, but not every landing is done this way. Landings on difficult weather conditions or in aircraft carriers, for example, are made with substantial engine power reserve in case the aircraft needs to abort landing. Also, most large aircraft uses reverse thrust to reduce landing distances.

Landing an aircraft means slowing it down from touchdown velocity, $V_{\mathrm{\ell }},$ to rest or nearly so. Figure 5 illustrates the process. The approach to the runway is made at a constant descent rate, in low traction regime and keeping the velocity above stall. When crossing the runway threshold, the pilot reduces engine power and rotates to increase drag and lift, reducing the velocity and the rate of descent. The airplane is kept “gliding” a few decimeters above the runway, while it losses velocity and lift, until it touches down.

Figure 5
Landing maneuver.

Right after touchdown, the pilot keeps the aircraft nose up for a few seconds to avoid impacts on the front landing gear. With the subsequent reduction in velocity, the nose pulls down, reducing the angle of attack and, thus, the lift. Deceleration continues until rest or taxi speed is achieved. During this stage, brakes, reverse thrust or lift spoilers may be used to reduce the landing distance.

The landing, as the take-off, is also carried out against the wind, so that the necessary lift is obtained at a lower velocity relative to the runway. This gives the pilot better control and also increases the drag, reducing the landing distance.

The “rough” landings everyone who flights frequently has experienced more than once, means that, for some reason, one or more of the procedures just described have failed. For example, in cases of heavy rain pilots intentionally approach the runway with a steeper rate of descent, making the touch down more abrupt to break the water film and avoid hydroplaning.

Small aircraft land in a similar way, but usually touch the runway with ${\text{SF}}_{\mathrm{\ell }}<1.3$. They also give the pilot considerably more freedom to chose the aircraft landing configuration, as runaways tend to be more than long enough for them. For example, the Pilot’s Operating Handbook for the Cessna 172N [¹⁸] says that “normal landing approaches can be made with power-on or power-off with any flap setting desired”. Small aircraft are also not generally equipped with spoilers for reasons of space and use brakes to reduce velocity. For the Cessna 172S the recommendation is to apply the brakes “as needed” in normal landings but to “apply maximum brake pressure without sliding the tires” in short field landings.

To model the landing, we return to equation 4 including a braking force, F_B. The integral then becomes

(24) $$s=\frac{W}{g}{\int }_{V_{\mathrm{\ell }}}^V\frac{V}{T-D-R-F_B}\mathit{d}V$$

where $V_{\mathrm{\ell }}$ is tha landing velocity and $V\in [V_{\mathrm{\ell }},0]$.

3.1. Braking force

Modeling the braking force in general is difficult because, as we have just pointed out, the pilot can apply the brake as necessary. It is, however, possible to establish a working relation for the braking force, multiplying the maximum braking force by a braking intensity coefficient ${\iota }_B$, as in [¹⁹], i.e.

(25) $$F_B={\iota }_B{\mu }_e(W-L)$$

where ${\mu }_e$ is the coefficient of static friction between the tires and the runaway and $(W-L)$ is the normal load. The term $F_{Bmax}={\mu }_e(W-L)$ represents the maximum force that can be applied by the brakes on the wheels without sliding. Any value larger than that will reduce the breaking efficiency, since the coefficient of dynamic friction is smaller than the coefficient of static friction. This is avoided automatically by the ABS system, when present, or by the pilot’s intervention when not.

The braking intensity, $0\le {\iota }_B\le 1$, quantifies the pressure exerted by the pilot on the brake pedal and the limits ${\iota }_B=0$ and ${\iota }_B=1$ correspond to no braking and maximum braking, respectively. Nevertheless, the actual braking process is described by the braking profile, i.e, by the function ${\iota }_B={\iota }_B(V)$.

In an optimization study, [²⁰] presents the typical breaking profile used in commercial aviation, reproduced in our Figure 6. Braking begins just after touchdown with maximal intensity. Then, after a rapid substantial velocity reduction, the aircraft is allowed to roll freely, coming to rest by the influence of drag and rolling friction. Figure 6 shows that braking time consists of 8.5% of the total time, which corresponds to 17.9% of the total distance and to a reduction of 50% in the velocity.

Figure 6
Typical braking profile for commercial size aircraft. Reproduced from [²⁰].

The formal similarity between eqs. (8) and (25) show that the use of brakes acts mathematically as an increment to the rolling friction. However, while typical values of ${\mu }_r$ for car tires on asphalt are around 0.01, those of ${\mu }_e$ are around 1.0 [¹⁵]. This explains why the use of brakes is so efficient, and indeed essential, in most landings.⁷

3.2. Equations of motion

As in the take-off analysis, some authors [¹] approximate the integral in equation(24) using mean values for the forces, also calculated at 70% of the landing velocity, while others [¹, ⁴] solve the integral numerically. In spite of the fact that there is no practical problem associated with that, it is completely unnecessary.

Substituting equation (25) into equation (24) and modeling the forces as before yields

(26) $$s=\frac{W}{g}{\int }_{V_{\mathrm{\ell }}}^V\frac{V}{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}\mathit{d}V$$

where

(27) $$A_{\mathrm{\ell }}=\rho \left[a{n}^2{D}^4+\frac{S}{2}\left(({\mu }_r+{\iota }_B{\mu }_e)C_L-C_{Dtot}\right)\right]$$

(28) $$B_{\mathrm{\ell }}=b\rho n^2{D}^4-2A_{\mathrm{\ell }}{U}_0$$

(29) $$C_{\mathrm{\ell }}=\left(C_{T0}-b{U}_0\right)\rho n^2{D}^4-({\mu }_r+{\iota }_B{\mu }_e)W+A_{\mathrm{\ell }}{U}_0^2$$

are independent of V if ${\iota }_B$ is constant. In this case, this integral is identical to equation(12) except for the limits of integration. The result appears in Eqs. (30) and (31).

(30) $$s=\frac{W}{2g{A}_{\mathrm{\ell }}}\left[\ln \right|\frac{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}{A_{\mathrm{\ell }}{V}_{\mathrm{\ell }}^2+B_{\mathrm{\ell }}{V}_{\mathrm{\ell }}+C_{\mathrm{\ell }}}|+\frac{2B_{\mathrm{\ell }}}{\sqrt{4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}-B_{\mathrm{\ell }}^2}}$$ $$\times (\arctan \frac{2A_{\mathrm{\ell }}{V}_{\mathrm{\ell }}+B_{\mathrm{\ell }}}{\sqrt{4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}-B_{\mathrm{\ell }}^2}}$$ $$-\arctan \frac{2A_{\mathrm{\ell }}V+B_{\mathrm{\ell }}}{\sqrt{4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}-B_{\mathrm{\ell }}^2}})]$$

and

(31) $$s=\frac{W}{2g{A}_{\mathrm{\ell }}}\left[\ln \right|\frac{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}{A_{\mathrm{\ell }}{V}_{\mathrm{\ell }}^2+B_{\mathrm{\ell }}{V}_{\mathrm{\ell }}+C_{\mathrm{\ell }}}|+\frac{B_{\mathrm{\ell }}}{\sqrt{B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}}}$$ $$\times (\ln \frac{2A_{\mathrm{\ell }}{V}_{\mathrm{\ell }}+B_{\mathrm{\ell }}+\sqrt{B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}}}{2A_{\mathrm{\ell }}{V}_{\mathrm{\ell }}+B-\sqrt{B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}}}$$ $$-\ln \frac{2A_{\mathrm{\ell }}V+B_{\mathrm{\ell }}+\sqrt{B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}}}{2A_{\mathrm{\ell }}V+B_{\mathrm{\ell }}-\sqrt{B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}}})]$$

for $4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}-B_{\mathrm{\ell }}^2>0$ and $B_{\mathrm{\ell }}^2-4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}>0$ respectively.

The landing distance can be obtained substituting V=0 into those equations.

The landing velocity, analogously to Eq.(19), is

(32) $$V_{\mathrm{\ell }}={\text{SF}}_{\mathrm{\ell }}\sqrt{\frac{2W}{\rho C_{Lmax}S}}+U_0$$

The integral in equation (26) has obviously no solution for an arbitrary braking profile. Therefore, rather than solving it for a few cases where an analytical solution is possible, we present the solution for a generalisation of Figure 6, i.e, a piecewise constant distributions of ${\iota }_B$ of the form

(33) $${\iota }_B(V)=\{\begin{array}{cc}{\iota }_1,\text{for}{V}_1\le V\le V_{\mathrm{\ell }}\hfill & \\ {\iota }_2,\text{for}{V}_2\le V<V_1\hfill & \\ \mathrm{⊊}\hfill & \\ {\iota }_m,\text{for}\mathit{}0\le V<V_{m-1}\hfill & \end{array}$$

where $V_{\mathrm{\ell }}>V_1>V_2>\mathrm{⊋}{V}_{m-1}>0,$ and ${\iota }_1$, ${\iota }_2,\mathrm{⊋}$, ${\iota }_m$ are constants.

Figure 6 depicts a case where ${\iota }_B=1.0$ for $0.5V_{\mathrm{\ell }}<V\le V_{\mathrm{\ell }}$ and ${\iota }_B=0$ for $0\le V\le 0.5V_{\mathrm{\ell }}$. In the general case, the integral in equation (26) must simply be solved by pieces between $V_{\mathrm{\ell }}$ and V yielding

(34) $$\begin{array}{c}s=\frac{W}{g}[{\int }_{V_{\mathrm{\ell }}}^{V_1}\frac{VdV}{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}\hfill \\ \hfill +{\int }_{V_1}^{V_2}\frac{VdV}{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}\\ \hfill +\mathrm{⊋}+{\int }_{V_k}^V\frac{VdV}{A_{\mathrm{\ell }}{V}^2+B_{\mathrm{\ell }}V+C_{\mathrm{\ell }}}]\end{array}$$

with k any integer in the range $[1,m]$. Each integral on the preceding summation has the exact same solution presented in Eqs. (30) or (31) substituting the appropriate limits of integration, but the whole expression is a bit large to fit in this paper. The total landing distance can be obtained substituting k=m into the last integral.

The physical interpretation of the case $4A_{\mathrm{\ell }}{C}_{\mathrm{\ell }}=B_{\mathrm{\ell }}^2$ is similar to the takeoff case. The second degree polynomial here represents the total resistive force during landing and must be strictly negative. The existence of one or more velocities where $\sum F_x=0$ in the interval $[V_{\mathrm{\ell }},0]$ would imply that the aircraft stops decelerating, thus never coming to a full stop.

In spite of all similarity between Eqs. (12) and (26), our calculations show that takeoff distances are considerably shorter than landing distances for commercial size airplanes when no brakes are used. Indeed, if there is no rolling friction and no brakes are applied, the aircraft slows down indefinitely but never stops, because as the velocity decreases, drag decreases faster as it is proportional to the square of the velocity. If the rolling friction is included in the model, it makes the resistive forces reduce at a smaller rate and is able to finally take the plane to rest. Even with rolling friction, however, landing distances with no brakes are generally large because ${\mu }_r$ is typically two orders of magnitude smaller than ${\mu }_e$ for commercial airplanes and one order fot the UAVs tested [⁸].

3.3. Static gliding wind velocity

We can equate the lift force to the total aircraft weight to calculate a headwind velocity where the airplane “floats” in the air in idle for a given angle of attack and is, consequently, unable to land. This velocity, hereafter denoted static gliding wind velocity is

(35) $$U_{0glide}=\sqrt{\frac{2W}{\rho C_{L}S}}$$

where $C_L$ is the lift coefficient for that particular angle of attack.

This equation, and the term static gliding wind velocity, are seldom, if ever, found in the aeronautical literature, because for commercial airplanes a very strong headwind would be necessary to satisfy this condition. Even for gliders, such headwinds are not usual.

Of course, being unable to land a commercial aircraft is not possible, because the pilot always has the option of pulling the plane’s nose down until C_L = 0. For UAVs however, wind gusts often achieve static gliding wind velocity. Thus, the pilot have to adjust the attitude of the aircraft (nose up or down) in a extremely rapid fashion to land.

4. Traction Measurements

Values for the coefficients a, b and $C_{T0}$ in equation (7) can obtained from the literature for a small number of commercial size aircraft⁸. For the UAVs, this kind of data should be made available by the power unit manufacturer, but our experience proved them to be unreliable. Therefore, a low cost traction measurement apparatus was built.

To build the apparatus, two vertical supports were attached onto a plywood base (Figure 7), one fixed and one pivoted. Onto the fixed one, a single point load cell (HBM model PW6DC3) was attached, and onto the other, the power unit to be tested was attached. The power unit was connected to the load cell by a steel cable just long enough to keep the pivoted support vertical. To measure the relative air speed, a pitot tube (Eagle Tree Systems) was installed in a lateral support, connected to a differential pressure sensor (Freescale Semiconductor MPXV5004DP). The output signals were acquired by a micro-controller (8-bit Atmel model ATmega328) and transmitted to a laptop where they were converted to velocity and traction.

Figure 7
Traction measurement apparatus.

The whole apparatus was fixed to a larger base attached to a pickup truck (Figure 8). A highway near the campus was used to perform the experiment on a windless day, varying the velocity of the pickup from rest to approximately 80 km/h, recording measurements during non-accelerated parts of the run. Measurements agreed well with a parabola, as expected. The parabola coefficients for all UAVs tested as well as the two models of the Cessna 172 (data obtained from [¹³]) are shown in Table 1.

Figure 8
Dynamic traction measurement system.

Table 1
Takeoff and landing distances and W_max.

	Unit	Commercial aircraft		UAV
		Cessna172N	Cessna172S	2014	2017	2018	2019
$\rho$	kg/m3	1.2250	1.2250	1.1226	1.1240	1.1226	1.0564
n	rpm	2700	2700	2500	2500	2500	2500
D	m	1.905	1.905	0.305	0.356	0.4064	0.4826
$a\times {10}^8$	${\text{(s/m)}}^2$	2.60	2.60	$-861$	$-878$	$-827$	$-879$
$b\times {10}^5$	s/m	$-1.44$	$-1.44$	$-53.5$	$-22.8$	$-24.8$	$-26.4$
$C_{T0}\times {10}^3$	–	189	189	21.7	18.4	24.5	26.1
M tot	kg	1043.00	1156.60	3.13	6.12	11.80	7.24
C Dtot	–	0.0320	0.0320	0.0646	0.0760	0.0771	0.0715
C L	–	0.41	0.41	0.44	0.77	0.85	0.95
C Lmax	–	2.1	2.1	1.4	2.1	2.2	2.2
S	m2	16.07	16.07	0.34	0.42	0.50	0.65
${\mu }_r$	–	0.045	0.045	0.11	0.11	0.13	0.13
${\mu }_e$	–	1.0	1.0	0.4	0.4	0.4	0.4
Sto exp.	m	219	292	45	39	60	40
Sto Eq. (17)	m	209	301	37	48	56	27
Difference	–	4%	3%	18%	22%	$-6$%	32%
Sto Eq. (20)	m	210	301	36	48	56	28
Wmax exp.	N	10230.9	10898.1	31.4	62.4	115.7	70,9
W max	N	11220.9	11220.9	34.1	63.5	110.6	87.8
Difference	–	4%	3%	10%	2%	$-4$%	24%
$S_{\ell }$ exp.	m	175	175	50	52	57	42
$S_{\ell }$ no brakes	m	799	927	82	104	143	74
$S_{\ell }$ brakes	m	182	211	35	43	61	32
Difference	–	4%	17%	30%	18%	7%	25%

5. Results

To demonstrate the wide range of aircraft to which the proposed results applies, the takeoff distance was calculated for the Cessna 172 Skyhawk models S and N, and for four UAVs used in the AeroDesign competitions of 2014, 2017, 2018 and 2019. Table 1 summarizes the input data, the takeoff and landing distances and the maximum weight.

For the Cessnas, the input data were obtained from [¹³, ¹⁸, ²¹] and for the UAVs they came from the team’s own measurements. Unfortunately, the absence of a telemetry system tends to make the later a little imprecise. During takeoff, the pilot cannot be sure that the airplane is at takeoff velocity and consequently tends to rotate a little later, reporting overestimated values of $S_{to}$. During the landing there is no way to be sure about the landing velocity either. Therefore, the experimental values in Table 1 represent our best estimates. This problem, which does not disqualify the methodology proposed, is currently under consideration, as new experiments and telemetry systems are being planed.

Another important detail about the input values involves the friction coefficients. The extensive use of our analytical results since 2015 in a toolbox developed by the authors showed that $S_{to}$ and $S_{\mathrm{\ell }}$ are rather sensitive to the value of ${\mu }_e$ used. To tackle the problem, we developed a simple methodology [⁸] to measure ${\mu }_r$ for the UAVs tested. The values obtained are considerably larger than the available values for large pneumatic tires over concrete used for the Cessna (truck tires, actually).

Table 1 shows that the values calculated for $s_{to}$ by equation (17) and the approximate solution equation (20) agree almost exactly for $\kappa =0.7$. The results were found not to be very sensitive to $\kappa$ for the Cessnas: a variation in $\kappa$ of $\pm 20\%$ resulted in a variation in $s_{to}$ of just $\pm 6\%$. For the UAVs however, the same variation in $\kappa$ resulted in a variation in $s_{to}$ of $\pm 12\%$.

Equation (22) was used to calculate $\kappa$ exactly for the 2014 UAV and for the Cessna 172N to access its dependence on parameters as $C_{Dtot}$ and $U_0$. For $0.04\le C_{Dtot}\le 0.09$ variations in $\kappa$ were never superior to $\pm 2\%$ and for $-10\le U_0\le 10$ m/s) never superior to $\pm 5\%$. Figure 9 illustrates the point for U₀. The variation is clearly seen to be more accentuated for the UAV.

Figure 9
Factor $\kappa$ dependence on the wind velocity.

Table 1 shows that the values calculated for $S_{to}$ by equation(16) underestimate the experimental values for the Cessna 172N and 172S by 4.3% and 2.8% respectively. For the UAV our model underestimates the measured value of $S_{to}$ by 17.8%. For the other UAVs the largest difference occurred for the 2019 prototype, 32%, and is probably due to poor measurement of engine’s angular velocity.

To access the impact of measurements errors in $S_{to}$ we supposed that the pilot rotated the UAV only 0.5 seconds after $V_{to}$ was achieved. In this case, we would have $S_{to}\approx 39$ m, yielding a difference of 6% between calculated and experimental values.

Table 1 also shows that the proposed methodology to calculate $W_{max}$ gives good predictions for all airplanes tested, except again for the 2019 prototype. For the Cessnas the difference is smaller than 4% and for the UAVs smaller than 10%.

Regarding the values calculated for $S_{\mathrm{\ell }}$ by Eqns. (31) and (34), the values obtained with no use of brakes largely overestimate the experimental values in all cases. Unfortunately, as we have mentioned before, there is no way to be sure if and how the Cessna pilots applied the brakes (remember the recommendation to apply the brakes “as needed”). Therefore, just to have some reference, we calculated the landing distances using full breaks until 40% of the landing velocity, a profile very similar to that in Figure 6. We obtained $S_{\mathrm{\ell }}=182$ m for the Cessna 172N and $S_{\mathrm{\ell }}=211$ m for the Cessna 172S (4% lower and 17 % larger than the measured values), in much better agreement with measurements.

For the UAVs with the same braking profile, we obtained a difference between calculated and experimental $S_{\mathrm{\ell }}$ between 7% and 30%. However, no brakes were available in the prototypes. Then, why assume a braking profile in this case?

Due to constructions limitations, our UAV’s wheels tend to be a bit unstable, making the airplanes woble during landing. Each time they woble, a small friction force appears between the wheels and the runaway, acting as a brake. Of course, assuming a braking profile to emulate the landing instabilities is definitively not the best practice, but we choose to include this result to illustrate the importance of introducing a braking system in future projects. However we are absolutely sure that the landing distance was never larger than 60, m simply because in this case the flight would have been disqualified, which never happened (by the competition rules, UAVs should land in less than 60 m).

To further illustrate the behavior of our analytical solution, relevant parameters for the operation of UAVs were varied for the 2014 prototype. The results are presented in the following subsections.

5.1. Takeoff results

Figure 10 shows the relation between velocity and displacement obtained from equation (17) with U₀ = 0, for different payloads around the value given on Table 1. Here we see one of the advantages of the present solution in comparison with the traditional approximation given by equation (20): the former yields the motion equation whereas the later only gives the value of $S_{to}$.

Figure 10
Motion equation for takeoff with $U_0=0$. $M_{tot}$ from Table 1.

Figure 10 shows that in the first 10–15% of the takeoff run, the airplane accelerates rather quickly because the power unit traction is near to its maximum value, while the dissipative forces are near their minimum. As the plane gains velocity the relation is progressively inverted, to the point takeoff occurs. The last point of each curve represents the takeoff distance which increases as total mass increases, as expected.

Figure 11 shows the dependence of the takeoff distance on the total weight of the aircraft for various wind velocities. The lower limit in the calculations was set as its empty mass. In addition to the expected reduction of $S_{to}$ with the increasing of the headwind and the decreasing of the mass, our results show that for certain values of $(W,U_0)$ the takeoff distance tends to zero, i.e., the aircraft takes off by itself as predicted by equation (35). For the UAV considered, equation (35) yields the values of $U_{0glide}$ presented in Table 2, where $M_{empty}$ is the mass of the empty aircraft.

Figure 11
Relation between takeoff distance and total weight. Negative wind velocities correspond to headwinds.

Table 2
Static gliding headwind velocity, U_0glide (m/s).

	M empty	M tot
C Lmax	5.5	10.7
C L	9.9	19.0

According to Table 2, for our UAV to get airborne on its on, a headwind of only 5.5 m/s would be enough, provided the wind also pulled the airplane nose up, which is not at all uncommon. For large commercial airplanes, on the other hand, this is rare, as the wind velocities involved would be typical of a hurricane. Nevertheless, a headwind of 21.8 m/s (less than 80 km/h) may float an empty Cessna 172N with its nose up.

Figure 12 shows the dependence between takeoff distance, total weight and total drag coefficient. This is an important result at the optimization phase of the project, where eventual adjustments of the fuselage and cargo bay directly impact on $C_{Dtot}$ while leaving all other features of the aircraft unchanged.⁹ Figure 12 shows that the dependence with $C_{Dtot}$ is not strong, with for example a 20% reduction in $C_{Dtot}$ causing a 7% reduction in $s_{to}$ for the considered mass of 3.13 kg.

Figure 12
Relation between takeoff distance and total weight for different total drag coefficients.

5.2. Landing results

Figure 13 shows the relation between velocity and displacement, obtained from equation (31), for different payloads, with no use of brakes, for landings in idle (13a) and with 20% of power unit thrust (13b). As we can see, after touchdown the airplane loses velocity at a fairly constant rate, being decelerated mainly by drag because the normal load, and thus the rolling friction, are still small due to the lift force being almost as large as the weight. As the plane loses velocity, friction gets larger and the deceleration increases until the airplane comes to a halt.

Figure 13
Motion equation for landing with U₀ = 0 for two different power unit thrust regimes.

Figure 14 shows the relation between velocity and displacement, obtained from Eqs. (31) and (34), for landings in idle, with five different braking profiles similar with the one depicted in Figure 6 and two distinct braking intensities: ${\iota }_B=1.0$ (14a) and ${\iota }_B=0.5$ (14b), applied from touchdown to braking limit velocities V₁ varying from $V_{\mathrm{\ell }}$ (no braking) to $0.4V_{\mathrm{\ell }}$. It is evident that even a conservative braking attitude, $V_1=0.8V_{\mathrm{\ell }}$, can affect the landing distance significantly. Also worth of note is the strong variation of the deceleration, represented by the derivative of the curves, as the aircraft goes from braking to free roll. Of course the effect is less dramatic for smaller braking intensities, as seen in Figure 14b for 50% use of the brakes.

Figure 14
Motion equation for landing with different braking profiles. Braking limit velocity V₁ defined at equation (33).

Figure 15 shows the dependence of the landing distance on the total weight for five distinct wind velocities. Again, the trend is as expected: $s_{\mathrm{\ell }}$ increases with the total mass and the tail wind. As before, the headwind may reach a point where the airplane glides statically and beyond this point it is almost impossible to land the airplane, which will be flying backwards. This phenomena is familiar to most participants of the AeroDesign competitions, where high lift small weigh UAVs where seen standing still in the air and even moving backwards¹⁰.

Figure 15
Relation between landing distance and total weight. Negative velocities correspond to headwinds.

6. The APT Toolbox

As mentioned, both the design and the operation of an aircraft require multiple analysis of the aircraft’s flight mechanics, and in particular of its performance. However, such an analysis requires input data related to the aerodynamic characteristics of the aircraft, its propulsion system, its dimensions and weight, the runway length available for takeoff and landing, the thermodynamic properties of the air at the operating location, etc. In addition, it inevitably ends up generating a large amount of output data.

Despite the fact that the physics used in the analysis of an aircraft’s performance are well described in the traditional literature, with some recent contributions as the one presented in this work, the large quantity of input and output data involved favors the use of a computational tool.

With the AeroDesign competition in mind and the result previously obtained [⁶], in the year of 2014 the authors set out to obtain a computational tool to aid the team in the calculations. Their experience showed that, in general, the computational packages developed for performance analysis must meet the following characteristics:

• Be supported by simple mathematical models but allow a complete assessment of the aircraft’s performance characteristics;

• Have a graphical interface that allows the creation or loading of a database containing the information of the aircraft to be analyzed, and allowing the selection of the types of analyzes to be performed;

• Generate a clear and flexible data output, allowing for the comparison of two or more aircraft. Also, the resulting graphics must have layouts that can be easily adapted to the user’s needs;

• Be able to be implemented in a computational environment disseminated in the academic environment and be freely distributed.

Furthermore, the computational tools must be able to describe, at least, the following flight phases: take-off, climb, unaccelerated level flight, level turn, gliding and landing.

Unfortunately, a large part of aeronautical science is developed by private companies, that keep confidentiality about the tools developed by them or charge for their use. Thus, just a few computational packages developed for aircraft performance analysis are presently available and the ones available do not meet some of the requirements mentioned. Therefore, the authors realized that the only way to assure that the team would have access to an efficient and complete computation tool, within an affordable cost, would be to develop their own. Thus, the Aircraft Performance Toolbox (APT) was created. It was designed to be a complete and free-to-use tool to help students calculate the most relevant aircraft performance parameters and help them explore the numerous possibilities associated with different aircraft configurations.

A number of Matlab scripts were created to perform the calculations, and a graphical interface was created to unify them and to facilitate the end-user operation, making it clear and intuitive. The software included our main results, Eqs. (16), (17), (30), (31) and the methofology to calculate $W_{max}$, together with a number of other important formulas for variables such as climb rate, curve radius, angle for maximum reach in gliding, and many others. To check the tool, the APT was used in the analysis of the UAV designed to the SAE Brasil AeroDesign competition of the year of 2014, the same used for comparisons here, and showed good results. Since then, the code has been extensively used by the students in the following year’s competitions and the evaluations returned were always extremely positive.

6.1. Software design

The free student version of the Matlab^® platform was the chosen development environment due to its renowned ability to deal with large numbers of variables and to perform calculations with large indexed variables. In spite of the fact that Matlab^® is a commercial software, due to its huge numerical problem solving capability and simplicity of coding, it is widely spread in the academic community, being available at most universities in its student edition and higher.

To ensure that the didactic objectives of the tool were fulfilled, the user interface was designed in a simple and intuitive way. Thus, students with little knowledge of airplane performance and even of aeronautics are able to carry out the analyses quickly and correctly. It is only necessary that the student enter the required data and select ”perform analysis”, and all graphs and values are calculated and presented. Also, a large number of explanatory remarks where included in the source code, for those with access to it. For users with only access to the executable code, the program was externally documented.

For the development of the toolbox, several scripts were implemented. Each one contained the equations necessary to perform one of the analyses, as well as the read commands for the input data and the write or plot commands for the output data. The set of scripts works interactively to allow the output data from one analysis to be used as input data from another.

A graphical interface was also created to allow the user to create a database designed to feed the scripts and to make their execution easy and intuitive. New databases can be created simply by filling in the blanks with the value of the variables. A previously created bank can be loaded via the ‘Load Database’ button. After inserting a data set, the user has the option of saving them to a text file using the ‘Save Database’ button. Once the database is available, the user must select the type of analysis he wants to perform by clicking on the corresponding option. By clicking on the button ‘Perform analysis’ the data entered in the interface at that moment will become available to the system. Before the code executes each script, it will verify that all fields have been filled in, and whether the indexed variables have the same number of elements. If not, a warning message will appear. After the execution of the script the numeric output data is presented in the Matlab command window and the graphs are displayed. The user can then save the graphic in one of several supported formats: png, jpg, bmp, pdf, etc.

At this point, the user can make other analysis, make changes to the database or perform the same analysis again. He can also load another database or simply clear the output data using the ‘Close graphs’ and ‘Clear Command Window’ buttons. If the user chooses to perform more than one analysis without closing the existing graphs, the graphs that relate to the same quantities will be displayed in the same window. This feature will allow the user to compare two different aircraft or to analyze the variation of the output parameters given the change in an input parameter.

6.2. Didactic experience

To access the usefulness of the tollbox, a survey was conducted among students who work and have worked in the performance area of the team, to learn their experiences in the analysis before and after using the software. Unfortunately, the success of the tool meant that just a few students had experience before and after its creation. Indeed, from 2015 onward, no student ever made the calculations by hand again.

The research was carried out through the Google Forms platform. According to the 8 forms received, students think that the toolbox considerably speeds up the laborious calculations and helps organize the output data, making it easy to test different aircraft configurations in the initial phase of the project. In addition, it is believed to help in the preparation of the final report and presentations. The tool was also credited with providing greater interaction between the areas, as it helps to organize the input and output data.

Another interesting feature mentioned by the students is that the toolbox can be used both in the forward and backward directions. In the forward direction, it is used for its original purpose, calculating performance characteristics of the UAV from the input data. In the reverse direction, it may be used to calculate input data from measured performance characteristics. An example of this is the possibility of calculating the value of the rolling resistance coefficient for a new set of tires, from the take-off distance measured experimentally.

The only disadvantage of APT detected in our diagnosis was that, as students understand that the code obtain the solutions “on its own”, they tend not to study the theory as they should. This leads to two problems: (i) students do not acquire the proper understanding about the subject; (ii) students spend a long time using the toolbox in trail and error tests, just to arrive at conclusions that were obvious from the theory (that they did not studied properly!).

Unfortunately, this situation seems to be commonplace in modern engineering courses, whenever a software is available to help in the calculations, and is rather worrying. In a research conducted more than one decade ago, [²²] describes characteristics and study habits of first-year students of Iowa State University’s electrical engineering course. Their findings, consistent with our own observations, identifies 15 common traits and behaviors of the 21^st century students. One of then was particularly associated with the problems just mentioned:

Students do not seek to find understanding, only answers. As a result, the students are resentful of the confusion that arises naturally from, and is necessary for, the learning process to be successful – that is, to foster understanding, and create lifelong learners…While the answer is important, the understanding of its discovery, or meaning, is of little value to them – that one has obtained an answer is sufficient.

Another observation closely related with what we consider a misuse of the computational tools is:

Students lack an understanding of the learning process. In effect, the students are attempting to gain theoretical knowledge…merely by working through numerous examples. Indeed, many students seem to believe that the process of learning consists only of working through such examples.

And last but not least:

Students lack an understanding of the meaning of hard work. Generally, unless they are captivated by a subject, they will not undertake it in earnest. Consequently, they are unwilling to perform the routine…necessary to develop their skills.

As for future versions of the toolbox, students showed interest in a module to calculate manual launch take-offs. The development of toolboxes for other areas of the project, with the capacity of interaction between them, was also mentioned and is already underway. Priority was requested in the future development of a module for a complete aerodynamic analysis toolbox, capable of generating the polar drag curve, to remove the need for the user to enter the data manually. An important, although obvious, detail about this request is that it depends on the availability of dedicated members of the Team, with adequate training in some modern programming language. Of course, being a member is not mandatory, but in our experience is highly recommended, as it generally guarantee the necessary theoretical basis. Indeed, the APT was developed by a second year electrical engineering undergraduate, member of the performance sector.

6.3. Conclusions about the APT

We believe that APT met all the requirements stipulated, including the didactic ones. In particular, the data output is made in a simple and clear way, allowing for a quick comparison between the different configurations studied. This is an important quality for such a software, once time is restricted in the kind of competition APT was designed to help with.

The evaluations of the tool returned by the (admittedly few) users that worked in the area before and after the implementation of APT, were considered rather positive. The only inconvenience, i.e., the tendency of the students to treat the software as a black box, is being currently addressed but it seems to be part of a much larger problem: the 21^st century students study habits. This is an issue that will certainly have to be dealt with, as computational tools like the APT each day occupy a larger share of the available tools for complex or laborious calculations. Our recommendation for colleagues who are thinking of creating and implementing similar tools, is to consider this question in as an early stage of the development of the tool as possible.

Regarding the results obtained with the aid of APT, the titles conquered by the Team in the competitions may be considered a good indicator of its success. Of course, the group relied on many other skills to obtain those titles, but the software was certainly one important factor in the process, as pointed out by most students involved.

7. Conclusions

This work presented an analysis of the take-off and landing dynamics of aircraft subjected to velocity dependent forces. It is an expansion of a previous work by Pellegrini and Rodrigues [⁶], including a number of new results: the generalization of the propeller-traction equation, the generalization of the motion equations, the inclusion of the average forces solution and the analysis of the $\kappa$ factor, the inclusion of the use of breaks with different braking profiles, the inclusion of an equation to calculate the static gliding wind velocity and, finally, the description of a toolbox created to perform the calculations.

The main results of the study are closed form analytical motion equations for the takeoff and landing of an airplane. These equations can be used to analytically calculate the takeoff and landing distances given the total weight of the airplane, or to numerically calculate the maximum weight given the available runaway length. They also allow the designers to investigate the behavior of the prototypes under different values of parameters like headwind velocity, maximum lift coefficient, wing area, braking profiles and others. The simplicity of the analytical results obtained allows them to be implemented on the optimization phase of the UAVs projects with small computational effort. This is particularly important in the AeroDesing competition, where the time schedule is often very tight.

Our equations were successfully applied to two commercial size airplanes and to four small scale UAVs. For most output variables, the model described the experimental data of the commercial size airplanes better than the UAVs’. The reason seems to be more related to the quality of the experimental data than to the accuracy of the theoretical model. Small UAVs like our prototypes are difficult to get instruments on board because of space and financial limitations. To partially circumvent this problem, a low cost traction measurement apparatus was designed and built in 2015, to generate the power unit curves. Efforts – a bit hindered by the COVID-19 pandemic – are under way to improve the traction measurements and to embark wireless telemetry systems in order to obtain reliable displacement and velocity data. A low cost system was also developed in 2019 to obtain more reliable values of the rolling friction coefficient.

The behavior of our solutions was further investigated varying a number of input parameters to evaluate their influence in the performance of the 2014 UAV prototype. Takeoff distances were calculated for different total weights of the UAV, head and tail winds and total drag coefficients. For the landing, the same analysis were made in addition to a thorough study of the influence of the braking profiles adopted and of the power unit thrust regime.

Acknowledgments

The authors would like to thank student Gabriel V. C. Pereira for his help debugging the MatLab codes.

References

• ^[1]
  J.D. Anderson, Aircraft Performance and Design (WCB/McGraw Hill, Boston, 1998).
• ^[2]
  R.S. Shevell, Fundamental of Flight (Prentice-Hall, Englewood Cliff, 1983).
• ^[3]
  C.D. Perkings and R.E. Hage, Airplane Performance Stability and Control (Wiley, New York, 1991).
• ^[4]
  J. Roskam and C.E. Lan, Airplane Aerodynamics and Perfomance (DARcorporation, Lawrence, 2000).
• ^[5]
  J.B. Russell, Performance and Stability of Aircraft (Butterworth-Heinemann, Burlington, 1996).
• ^[6]
  C.C. Pellegrini and M.S. Rodrigues, Rev. Bras. Ensino Fís. 37, 2307 (2015).
• ^[7]
  https://www2.anac.gov.br/anacpedia/glossario.htm, accessed in 08/11/2021.
  » https://www2.anac.gov.br/anacpedia/glossario.htm
• ^[8]
  C.C. Pellegrini and V.R. Alves, Rev. Bras. Ensino Fís. 42, e20190321 (2020).
• ^[9]
  C.C. Pellegrini, M.S. Rodrigues and E.D.O. Moreira, arXiv:2109.00355 (2021)
• ^[10]
  G.A. Williamson, B.D. McGranahan, B.A. Broughton, R.W. Deters, J.B. Brandt and M.S. Selig, Summary of low-speed airfoil data (University of Illinois, Champaign, 2012), v. 5.
• ^[11]
  S.K. Ojha, Flight Performance of Aircraft (American Institute of Aeronautics and Astronautics, Washington, 1995), p. 166.
• ^[12]
  C. Tjhai, Developing Stochastic Model of Thrust and Flight Dynamics for Small UAVs Masters Dissertation, University of Minnesota, Minneapolis (2013).
• ^[13]
  J. McIver, Cessna Skyhawk II/100 Performance Assessment (Temporal Images, Melbourne, 2003).
• ^[14]
  I.G. Goryacheva, Contact Mechanics in Tribology (Springer, Dordrecht, 1998).
• ^[15]
  A.V. Andrade-Neto, J.A. Cruz, M.S.R. Miltao and E.S. Ferreira, Rev. Bras. Ensino Fís. 35, 3704 (2013).
• ^[16]
  F.M. White, Fluid Mechanics (McGraw-Hill Education, New York, 2015), 8th ed.
• ^[17]
  Federal Aviation Administration. Federal Aviation Regulations – Airworthiness Standards: Normal, Utility, Acrobatic and Commuter Category Airplanes, Part 23, Sec. 23.51 (Federal Aviation Administration, Washington, 1989).
• ^[18]
  Cessna Aircraft Company, Pilot’s Operational Handbook: Skyhawk Cessna Model 172N (Cessna, Wichita, 1978).
• ^[19]
  D.O. Dommasch, S.S. Sherby and T.F. Connolly, Airplane aerodynamics (Pitman, New York, 1961).
• ^[20]
  J. Skorupski and H. Wierzbicki, Transp. Research Part C Emerging Technologies 80, 467 (2017).
• ^[21]
  Cessna Aircraft Company, Cessna Skyhawk 172S, Specification and Description (Cessna, Wichita, 2012).
• ^[22]
  M. Mina and R.M. Gerdes, European J. of Eng. Education 31, 509 (2006).

• 1
  According to a very traditional source [¹], Sir George Cayley was the first to realize that, in 1799. Nevertheless, all information coming form North American sources should be considered with caution, because the story tends to be told rather tendentiously.
• 2
  A thorough compilation of aeronautical terms can be found in the ANACpédia [⁷]. It is a data base for technical terms and abbreviations in the area, in English, Portuguese, Spanish and French. Hereafter, the aeronautical terms not defined in the text will be emphasized at their first appearance.
• 3
  The complete description may be found in the preprint deposited in the arXiv system [⁹].
• 4
  For this reason, in standard aircraft design, the wings are mounted to the fuselage nearly horizontally.
• 5
  In fact, they must depend on the propeller geometry: external diameter, hub diameter, airfoil used, aspect ratio, chord profile, twist profile, etc.
• 6
  The correct way to avoid the form $n^2<0$ which gives the impression of $n\in ℂ$, is to rewrite equation (5) as $T=\rho n|n|D^4{C}_T$. However this would complicate our notation quite a bit and thus we choose to keep the simpler form with $n^2$.
• 7
  Of course, in this sense, it is also useful to increase the value of the normal force. This is done by the spoilers in medium to large size aircraft. Mathematically, they act on equation (24) increasing all terms in the denominator of the integral.
• 8
  Data for “new” airplanes, as an Embraer 190 for example, are still protected by industrial secrecy.
• 9
  In aeronautical project, given the initial specifications of the aircraft and the flight mission, the project team has to make initial estimates for variables as the power unit traction, the propeller diameter, the weight of the aircraft or the length available for the takeoff and the total drag. Given these values, the other relevant quantities can be calculated. However, when the detailed project phase is reached, some of the preceding values may have to be changed, thus the importance in knowing the impact of such changes.
• 10
  As this story might seem a little difficult to believe, the following video shows one of our empty UAVs standing still in the air (02:30) at the year 2016 competition: https://youtu.be/QY8Dt6J0bGc. The 2016 rules required that the airplanes dumped its cargo on parachute.

Publication Dates

• Publication in this collection
  24 Jan 2022
• Date of issue
  2022

History

• Received
  23 Nov 2021
• Accepted
  13 Dec 2021