Albatross-Like Utilization of Wind Gradient for Unpowered Flight of Fixed-Wing Aircraft

: The endurance of an aircraft can be considerably extended by its exploitation of the hidden energy of a wind gradient, as an albatross does. The process is referred to as dynamic soaring and there are two methods for its implementation, namely, sustainable climbing and the Rayleigh cycle. In this study, the criterion for sustainable climbing was determined, and a bio-inspired method for implementing the Rayleigh cycle in a shear wind was developed. The determined sustainable climbing criterion promises to facilitate the development of an unpowered aircraft and the choice of a more appropriate soaring environment, as was demonstrated in this study. The criterion consists of three factors, namely, the environment, aerodynamics, and wing loading. We develop an intuitive explanation of the Raleigh cycle and analyze the energy mechanics of utilizing a wind gradient in unpowered ﬂight. The energy harvest boundary and extreme power point were determined and used to design a simple bio-inspired guidance strategy for implementing the Rayleigh cycle. The proposed strategy, which involves the tuning of a single parameter, can be easily implemented in real-time applications. In the results and discussions, the effects of each factor on climbing performance are examined and the sensitivity of the aircraft factor is discussed using ﬁve examples. Experimental MATLAB simulations of the proposed strategy and the comparison of the results with those of Gauss Pseudospectral Optimization Software conﬁrm the feasibility of the proposed strategy.


Introduction
An aircraft usually consumes a large amount of energy during flight. There is, however, the question of whether an aircraft can fly indefinitely without carrying a power source. This might be realizable if a method could be developed for the aircraft to continuously extract energy from its environment. A wind field in which the wind gradient varies with the altitude is known to contain a large amount of energy [1][2][3]. In addition, albatrosses have been observed to soar for long periods without flapping their wings [4][5][6][7]. Lord Rayleigh was perhaps the first to discuss the problem relevant to albatrosses soaring in a horizontal but non-uniform wind field [8]. This type of flight involves the extraction of energy from a wind gradient and is referred to as dynamic soaring.
The phenomenon of dynamic soaring has attracted the interest of many scholars. The technique utilized by an albatross to obtain extra airspeed for soaring is referred to as the Rayleigh cycle [9] and is illustrated in Figure 1. Idrac [10] showed that wind speed increases with the altitude, which generally occurs over the sea where albatrosses fly. Different researchers have observed birds to understand their soaring mechanics [11][12][13], and theories regarding their utilization of wind gradient have been developed. Wood [14] used a two-dimensional model to investigate two types of complete soaring flight cycles in a logarithmic wind velocity profile. Hendriks [15]  Because the wind speed below the wave is nearly zero and the wind above the wave has a certain speed, the wind can be considered to comprise two layers with a wind shear between them. Suppose that the initial airspeed of the albatross is 17 m/s in the wind layer below the wave. If the wind speed above the wave is 7 m/s higher than that below the wave, when the albatross faces the wind and crosses the wind shear between the lower and upper wind layers, the airspeed would suddenly increase to 24 m/s. This sudden gain in speed is due to the lack of time for the albatross to change its ground speed.
With the development of computer science, the use of numerical algorithms has become a popular approach to determine the optimal trajectory under a wind gradient. Zhao [19,20] used NPSOL, a numerical optimization software, to minimize the dynamic soaring cycles of unpowered aircraft and powered aircraft. Sachs [21,22] also solved the optimal control problem of dynamic soaring and computed the minimum shear wind strength under certain conditions. Deittert et al. [23] utilized the differential flatness property to determine the optimal trajectory for cross-country flight. Kaushik et al. [24] used the GPOPS-II optimal control software for trajectory computation. There have been several other studies on the optimization of dynamic soaring [25][26][27][28][29].
Recent improvements were achieved in the practical application of dynamic soaring. Lawrance et al. [30,31] proposed an algorithm for determining the flight direction for the quickest energy harvest. Langelaan et al. [32] and Rodriguez et al. [33] proposed a method for estimating the wind field for trajectory design. Liu et al. [34] presented a nonlinear model of a predictive control strategy for trajectory determination for energy extraction. Barate et al. [35] controlled a simulated glider based on a set of Takagi-Sugeno-Kang fuzzy rules to exploit wind gradient energy. Although these contributions improved the sophistication of the concept of dynamic soaring, this was accompanied by a more complicated practical implementation process. In contrast, Gao et al. [36] and Yang [37] proposed a simple method for generating a guidance strategy for autonomous dynamic soaring. This involved the division of the Rayleigh cycle into four sections, the characteristics of which were determined by simple equations. Gao and Liu [38] also developed a dynamic soaring plan in which the Dubins path was used to improve the computational efficiency. They showed that dynamic soaring could be achieved using modern technology. However, their strategy has the shortcoming of involving numerous control parameters, which makes practical implementation cumbersome.
In the present study, we followed the simplification route of Gao et al. [36] in developing a bio-inspired Rayleigh cycle strategy that is more efficient in terms of computation time than the traditional optimization method. In addition, we developed a method for utilizing a wind gradient for sustainable climbing [39] in unpowered flight. Simulation experiments and analysis were used to Because the wind speed below the wave is nearly zero and the wind above the wave has a certain speed, the wind can be considered to comprise two layers with a wind shear between them. Suppose that the initial airspeed of the albatross is 17 m/s in the wind layer below the wave. If the wind speed above the wave is 7 m/s higher than that below the wave, when the albatross faces the wind and crosses the wind shear between the lower and upper wind layers, the airspeed would suddenly increase to 24 m/s. This sudden gain in speed is due to the lack of time for the albatross to change its ground speed.
With the development of computer science, the use of numerical algorithms has become a popular approach to determine the optimal trajectory under a wind gradient. Zhao [19,20] used NPSOL, a numerical optimization software, to minimize the dynamic soaring cycles of unpowered aircraft and powered aircraft. Sachs [21,22] also solved the optimal control problem of dynamic soaring and computed the minimum shear wind strength under certain conditions. Deittert et al. [23] utilized the differential flatness property to determine the optimal trajectory for cross-country flight. Kaushik et al. [24] used the GPOPS-II optimal control software for trajectory computation. There have been several other studies on the optimization of dynamic soaring [25][26][27][28][29].
Recent improvements were achieved in the practical application of dynamic soaring. Lawrance et al. [30,31] proposed an algorithm for determining the flight direction for the quickest energy harvest. Langelaan et al. [32] and Rodriguez et al. [33] proposed a method for estimating the wind field for trajectory design. Liu et al. [34] presented a nonlinear model of a predictive control strategy for trajectory determination for energy extraction. Barate et al. [35] controlled a simulated glider based on a set of Takagi-Sugeno-Kang fuzzy rules to exploit wind gradient energy. Although these contributions improved the sophistication of the concept of dynamic soaring, this was accompanied by a more complicated practical implementation process. In contrast, Gao et al. [36] and Yang [37] proposed a simple method for generating a guidance strategy for autonomous dynamic soaring. This involved the division of the Rayleigh cycle into four sections, the characteristics of which were determined by simple equations. Gao and Liu [38] also developed a dynamic soaring plan in which the Dubins path was used to improve the computational efficiency. They showed that dynamic soaring could be achieved using modern technology. However, their strategy has the shortcoming of involving numerous control parameters, which makes practical implementation cumbersome.
In the present study, we followed the simplification route of Gao et al. [36] in developing a bio-inspired Rayleigh cycle strategy that is more efficient in terms of computation time than the traditional optimization method. In addition, we developed a method for utilizing a wind gradient for sustainable climbing [39] in unpowered flight. Simulation experiments and analysis were used to verify the effectiveness of the proposed flight method. The results showed that the proposed method could be used to efficiently guide a micro fixed-wing aircraft through effortless soar in a shear wind.
The rest of this paper is organized as follows. In Section 2, the problem of dynamic soaring is defined and the development of a dynamic model of the concept is described. In Sections 3 and 4, the criterion for sustainable climbing is derived, the Rayleigh cycle is explained from an energy viewpoint, and a bio-inspired strategy for implementing the cycle is proposed. In Section 4, the effects of different factors in the criterion are subsequently analyzed and discussed. The results of experimental simulations involving the bio-inspired strategy for implementing the Rayleigh cycle are then presented and validated by the traditional optimized trajectory. Finally, the conclusions drawn from the study are presented. The innovative ideas of this paper are as follows: On the basis of equilibrium theory, a novel criterion for sustainable climbing that is easily implemented for engineers is derived, and this criterion decouples the effects of aerodynamics, environment and wing loading. Moreover, a simple bio-inspired strategy is proposed to generate a Rayleigh cycle, which consumes much less time than the traditional numerical optimization method.

Dynamic Model
A wind shear is generated by the wind boundary layer phenomenon, and its mathematical model can take different forms. The following example is presented in the Military Specification MIL-F-8785C [40]: where W denotes the mean wind velocity, W 6 denotes the wind velocity at an altitude of 6 m, z is the altitude, and z 0 is generally set to 0.5 m.
A wind shear can also be modeled by an exponential wind velocity profile [36]: where W R denotes the reference wind velocity at altitude Z R , and p denotes the slope of the boundary layer. Let us assume a horizontal wind, the magnitude of which varies with the altitude [36]. The wind gradient can be expressed as follows: There are two common approaches to the modelling of dynamic soaring, namely, those proposed by Zhao and Qi [20] and Sachs and Bussotis [21]. The latter utilizes Cartesian coordinates, while the former uses spherical coordinates. A non-inertial coordinate frame that moves with the wind at the aircraft was adopted in the present study (see Figure 2). The z-axis of the frame points skyward. If the wind velocity increases with the altitude, the x-axis points in the opposite direction of the wind; otherwise, it points in the direction of the wind. The y-axis is determined by the right-hand principle. In Figure 2, W is the wind velocity vector at the aircraft position, and O 0 -x 0 y 0 z 0 represents the origin at the initial time. The problem formulation and modeling procedures adopted in this study are similar to those of some relevant previous works [16,21,23,36]. The aircraft is under the action of its weight, aerodynamics, and inertial force. However, the model equations differ from ordinary equations in that the atmosphere-fixed frame is non-inertial. Further, the model transforms the wind gradient into an inertial force. Figure 1 illustrates an intuitive model that facilitates an understanding of the generation of the inertial force. Because the wind speed below the wave is nearly zero, and that above the wave has a certain magnitude, the wind can be considered to comprise two layers separated by a wind shear, around which a large gradient occurs. Let us assume that the airspeed relative to the albatross shown in Figure 1 is 17 m/s in the wind layer below the wave, and the wind speed above the wave is 7 m/s higher than that below the wave. If the albatross faces the wind and crosses the wind shear from the lower layer to the upper layer, the airspeed relative to it would suddenly increase to 24 m/s because there is no time for a gradual change in the ground speed. The wind gradient thus causes the albatross to gain extra airspeed, which is transformed into an inertial force, referred to as a dynamic soaring thrust. The dynamic model can be directly taken from Deittert's [27] or Gao's [36] work as follows: c o s c o s s i n s i n 2 and where v denotes the airspeed, θ denotes the pitch angle of the flight path, ϕ denotes the bank angle and ψ denotes the yaw angle measured anticlockwise from x-axis. For the parameters, g is the gravitational acceleration, ρ is the density of air, s the wing area, m is the aircraft mass, and CL and CD are respectively the lift and drag coefficients. According to the geometric relationship shown in Figure 2, we transform the equations above to the following equations in Cartesian coordinates by vx = vcosθcosψ, vy = vcosθsinψ and vz =vsinθ, where vx, vy, and vz denote the velocities of the aircraft along the respective axes.  The problem formulation and modeling procedures adopted in this study are similar to those of some relevant previous works [16,21,23,36]. The aircraft is under the action of its weight, aerodynamics, and inertial force. However, the model equations differ from ordinary equations in that the atmosphere-fixed frame is non-inertial. Further, the model transforms the wind gradient into an inertial force. Figure 1 illustrates an intuitive model that facilitates an understanding of the generation of the inertial force. Because the wind speed below the wave is nearly zero, and that above the wave has a certain magnitude, the wind can be considered to comprise two layers separated by a wind shear, around which a large gradient occurs. Let us assume that the airspeed relative to the albatross shown in Figure 1 is 17 m/s in the wind layer below the wave, and the wind speed above the wave is 7 m/s higher than that below the wave. If the albatross faces the wind and crosses the wind shear from the lower layer to the upper layer, the airspeed relative to it would suddenly increase to 24 m/s because there is no time for a gradual change in the ground speed. The wind gradient thus causes the albatross to gain extra airspeed, which is transformed into an inertial force, referred to as a dynamic soaring thrust. The dynamic model can be directly taken from Deittert's [27] or Gao's [36] work as follows: where v denotes the airspeed, θ denotes the pitch angle of the flight path, ϕ denotes the bank angle and ψ denotes the yaw angle measured anticlockwise from x-axis. For the parameters, g is the gravitational acceleration, ρ is the density of air, s the wing area, m is the aircraft mass, and C L and C D are respectively the lift and drag coefficients. According to the geometric relationship shown in Figure 2, we transform the equations above to the following equations in Cartesian coordinates by v x = vcosθcosψ, v y = vcosθsinψ and v z = vsinθ, where v x , v y , and v z denote the velocities of the aircraft along the respective axes. .
where P Lm and P Dm are the relevant aerodynamic parameters, given by where the term s/m in parentheses represents the reciprocal of the wing loading, and the lift-to-drag ratio L/D is given by

Equilibrium Points and Sustainable Climbing
An equilibrium point is an important concept of a nonlinear system . x = f (x). A point in the state space is an equilibrium point if, when the state of the system begins at the point, it remain there at all future times. Therefore, at an equilibrium point, .
x ≡ f (x) ≡ 0 must be satisfied. In the dynamic soaring system represented by Equations (4)-(6), the equilibrium points are defined by the constant velocities that can be maintained by the aircraft, namely, For a constant bank angle ϕ, Zhu et al. [39] proved the calculability of the equilibrium points of a dynamic soaring system. The equilibrium also indicates the force balance during flight through a gradient wind. If the value of v z at an equilibrium point is negative, the aircraft would continue to dive at the corresponding velocity. This type of equilibrium point is thus referred to as a descending equilibrium point. Conversely, if the value of v z at an equilibrium point is positive, the aircraft would continue climbing. This type of equilibrium point is referred to as an ascending equilibrium point. The existence of a descending equilibrium point has been proved to be trivial, while an ascending one cannot be guaranteed to exist. Should an ascending equilibrium point occur, the aircraft could utilize it for sustainable climbing. Zhu et al. [39] proved that the sufficient and necessary condition for the existence of an ascending equilibrium point was the satisfaction of the following inequality: where where R gam and η are abbreviated parameters defined as The above equations constitute the initial criterion for sustainable climbing. It is, however, difficult to apply this criterion because it requires the determination of the maximum of f 2 (v z ) at every time. In Section 3.2, we propose a simplified criterion based on Equation (12). The derivation of the proposed criterion is presented in Section 3.3.

Proposed Criterion for Unpowered Climbing
The simplified criterion for unpowered climbing is explicitly proposed in this subsection and derived in the next subsection. The criterion for an unpowered aircraft to continue climbing in a gradient wind field can be expressed as follows: where Π a and Π e are the aircraft and environmental fractions, respectively. A smaller Π a Π e increases the climbing probability of the aircraft. The environmental fraction Π e is given by which is only relevant to the environment. The aircraft fraction Π a can also be divided into two terms as follows: where Π A and Π S are the aerodynamic and wing loading fractions, respectively. Π A is defined as where h is a function of the lift-to-drag ratio, and can be approximated as The aerodynamic fraction Π A depends on only the aerodynamic characteristic of the aircraft, while the wing loading fraction Π S is given by Hence, the unpowered climbing criterion can be divided into three parts, namely, the environmental, aerodynamic, and wing loading fractions. The derivation procedure of the criterion is presented below.

Derivation of Proposed Criterion
We define Thus, the initial criterion in Equation (12) can be written as Before the derivation, v z is normalized as By substituting Equation (25) into Equation (13) and canceling the common factor, the function f 2 becomes where The range of the variant v zn is [0, 1]. Thus, the function f 2n of v zn was determined to depend on only the unique parameter η. Therefore, K 1 naturally becomes a function of η: It is difficult to obtain an analytical expression of Equation (28). Indeed, because the function of K 1 includes a maximization operation, an analytical expression is almost impossible to obtain, although the values of the equation can be obtained for discrete samples. Such values are, however, actually obtained by numerical calculations. For any value of η, the maximum of f 2n (v zn ) is searched for in the interval [0, 1]. A series of discrete samples are thus obtained. The numerically determined values of K 1 are plotted in Figure 3, which clearly shows a monotonic increase. Before the derivation, vz is normalized as By substituting Equation (25) into Equation (13) and canceling the common factor, the function f2 becomes where Ω η + +η The range of the variant vzn is [0, 1]. Thus, the function f2n of vzn was determined to depend on only the unique parameter η. Therefore, K1 naturally becomes a function of η: It is difficult to obtain an analytical expression of Equation (28). Indeed, because the function of K1 includes a maximization operation, an analytical expression is almost impossible to obtain, although the values of the equation can be obtained for discrete samples. Such values are, however, actually obtained by numerical calculations. For any value of η, the maximum of f2n(vzn) is searched for in the interval [0, 1]. A series of discrete samples are thus obtained. The numerically determined values of K1 are plotted in Figure 3, which clearly shows a monotonic increase. The value of η should be larger than 1; otherwise, the inequality in Equation (24) would not be satisfied because of the negative value that K1(η) would have. Hence, K1 guarantees its own inverse function K −1 . Based on Equation (24), the following inequality would also have to be satisfied to achieve the initial unpowered climbing criterion: The value of η should be larger than 1; otherwise, the inequality in Equation (24) would not be satisfied because of the negative value that K 1 (η) would have. Hence, K 1 guarantees its own inverse function K −1 . Based on Equation (24), the following inequality would also have to be satisfied to achieve the initial unpowered climbing criterion: Appl. Sci. 2017, 7, 1061 8 of 25 The expression of η can thus be rewritten as follows: By substituting Equation (30) into Equation (29), the latter becomes where Utilizing the numerical results of K 1 , we sought backward for the corresponding L/D values, and thus obtained discrete samples of h. The function h can be approximated by familiar functions such as Equation (21), the numerical results of which are plotted in Figure 3. The approximation is for L/D values of 0.3-60, because the lift-drag ratio for a fixed-wing aircraft is generally within 1-60 [41]. This range covers most fixed-wing aircraft in practice. The highest lift-to-drag ratio for an albatross is about 20-30 [22]. Hence, Equation (31) establishes the criterion for unpowered sustainable climbing.

Rayleigh Cycle Method
In the real world, the wind gradient is generally not sufficiently large to enable sustainable climbing by an aircraft. The Rayleigh cycle affords another flight strategy for utilizing a gradient wind field for unpowered flight. This strategy is based on the observation that albatrosses are capable of flying for long periods without flapping their wings. Idrac [10] observed that the wind layer where the albatrosses fly increases according to the altitude. A Rayleigh cycle can be divided into four phases, as shown in Figure 4.
Let us suppose that the wind over the ocean surface in Figure 4 is divided by the wind shear into two layers, with the wind speed above the wave being 7 m/s and that below it being 0 m/s. The first phase of a Rayleigh cycle is the upwind climb segment. At the beginning of this phase, both the airspeed and the ground speed of the albatross are 17 m/s upwind. After the albatross crosses the wind shear, the ground speed is maintained at 24 m/s, but the airspeed is changed to 24 m/s. The second phase is the high-altitude turn segment. In this phase, the albatross executes an approximately uniform circular motion in the air, with its heading changing from upwind to downwind. Thus, the airspeed is still 17 m/s, but the ground speed changes to 31 m/s. The third phase is the downwind dive segment. Similar to the upwind climb phase, in this phase, the albatross crosses the wind shear in the wind direction. The ground speed is maintained at 31 m/s, and the airspeed also changes to 31 m/s. The fourth phase is the low altitude turn segment, which completes the Raleigh cycle, linking back to the first phase. Similar to the high-altitude turn phase, in this fourth phase, the albatross executes an approximately uniform circular motion in the air, and then returns to the beginning states. The final airspeed and ground speed are both 31 m/s, indicating that the albatross almost gains double the upper layer wind speed during a Raleigh cycle. In other words, it harvests energy from the wind shear. Figure 4 is only an intuitive illustration of the Rayleigh cycle. In practice, the wind shear does not produce a sudden change in the wind speed, because it has a continuous wind gradient profile. The altitude fluctuation during a Rayleigh cycle is also considerable. Thus, there is always a transformation between gravitational potential energy and kinetic energy over the entire cycle, with the drag also Appl. Sci. 2017, 7, 1061 9 of 25 consuming energy. Therefore, we prove below that the upwind climb and downwind dive can be used to harvest extra energy. Figure 4 is only an intuitive illustration of the Rayleigh cycle. In practice, the wind shear does not produce a sudden change in the wind speed, because it has a continuous wind gradient profile. The altitude fluctuation during a Rayleigh cycle is also considerable. Thus, there is always a transformation between gravitational potential energy and kinetic energy over the entire cycle, with the drag also consuming energy. Therefore, we prove below that the upwind climb and downwind dive can be used to harvest extra energy.

Energy Perspective
Here, E denotes the mechanical energy, including the gravitational potential energy and kinetic energy, and the O-xyz coordinate system is as shown in Figure 2. According to the kinetic energy theorem, the change in the mechanical energy is equal to the work done against the drag and inertial force. Because the lift is perpendicular to the movement, it does no work. The following expression can thus be derived: (33) Here, dl is defined as dl = vdt (34) This implies that Hence, the normalized rate of the mechanical energy acceleration is equal to It should be noted that the temporal variation of the mechanical energy is independent of the control input ϕ. The coordinates can thus be transformed as follows:

Energy Perspective
Here, E denotes the mechanical energy, including the gravitational potential energy and kinetic energy, and the O-xyz coordinate system is as shown in Figure 2. According to the kinetic energy theorem, the change in the mechanical energy is equal to the work done against the drag and inertial force. Because the lift is perpendicular to the movement, it does no work. The following expression can thus be derived: Here, dl is defined as This implies that Hence, the normalized rate of the mechanical energy acceleration is equal to It should be noted that the temporal variation of the mechanical energy is independent of the control input ϕ. The coordinates can thus be transformed as follows: where r is the radial velocity, θ is the pitch angle, and ψ is the heading angle. Equation (36) can therefore be rewritten as follows: (E/m) = −P Dm r 3 + G w r 2 sin θ cos θ cos ψ (38) Figure 5 is a sectional drawing of a contour plot of the normalized rate of the mechanical energy acceleration d(E/m)/dt for P Dm = 0.002 and G w = 0.4. It clearly shows that Equation (38) is symmetrical about the vxOv z plane. The equation of the surface on which the energy is invariant (i.e., d(E/m)/dt = 0) is as follows: This corresponds to the green surface in Figure 5. It exists where v xvz > 0 and comprises two leaves joined at the origin, having a spindle-like shape. Energy is only harvested inside the spindle, which is thus called an energy-harvest spindle. In practice, the Rayleigh cycle is a loop across the two leaves of the spindle. The upwind climb (v x > 0, v z > 0) and downwind dive (v x < 0, v z < 0) obviously occur in the area where v xvz > 0. Energy would obviously be harvested if the aircraft was positioned within the two leaves. where r is the radial velocity, θ is the pitch angle, and ψ is the heading angle. Equation (36) can therefore be rewritten as follows: This corresponds to the green surface in Figure 5. It exists where vxvz > 0 and comprises two leaves joined at the origin, having a spindle-like shape. Energy is only harvested inside the spindle, which is thus called an energy-harvest spindle. In practice, the Rayleigh cycle is a loop across the two leaves of the spindle. The upwind climb (vx > 0, vz > 0) and downwind dive (vx < 0, vz < 0) obviously occur in the area where vxvz > 0. Energy would obviously be harvested if the aircraft was positioned within the two leaves. There is a maximum rate of the mechanical energy acceleration in the spindle. In addition, Equation (38) can be rewritten as There is a maximum rate of the mechanical energy acceleration in the spindle. In addition, Equation (38) can be rewritten as To maximize Equation (40), the term sin2θcosψ should be 1, because 1 is the maximum of sin2θcosψ when θ = π/4, ψ = 0 or θ = −π/4, ψ = π. This implies that the maximum point is on the straight line defined by v x = v z and v y = 0. The maximum point is obtained by setting the derivative of Equation (40) to zero as d dr from which the following is obtained: Equation (41) can be expressed in the Cartesian coordinate frame as The maximum power is thus given by This leads to the conclusion that a larger gradient, smaller drag coefficient, lower air density, and larger wing loading enhance energy harvesting for unpowered flight. This agrees with the results of the proposed criterion presented below in Section 5.1.

Optimization Method for Generating Rayleigh Cycle
The generation of a Rayleigh cycle is generally transformed into an optimal trajectory problem, which has been solved by many numerical optimization methods such as the Gauss pseudo-spectral method of the optimal software called Gauss Pseudospectral Optimization Software (GPOPS) [24][25][26][27][28]. The optimized trajectory should satisfy the dynamics Equations (4)- (8), and the performance function gives the total energy harvest in a cycle, expressed as where t f and t 0 are the final and initial times, respectively. The periodic boundary condition is as follows: The constraints on the lift coefficient C L and bank angle ϕ are Though this method for generating a Rayleigh cycle provides an optimized solution, the numerical computation is too time-consuming for real-time applications. In the present work, GPOPS [42], which is a MATLAB-based optimal process toolbox, was used to optimize the problem. The optimized results were then used to validate the proposed bio-inspired strategy for implementing the Rayleigh cycle, described in the next subsection.

Bio-Inspired Strategy
Although the optimization method affords an optimal means of generating a Rayleigh Cycle [18][19][20][21][22][23][24][25][26][27][28], the computation procedure, as noted above, is too time-consuming for real-time applications. Hence, we used a state machine to imitate the four phases of the Rayleigh cycle of an albatross, namely, the upwind climb, high-altitude turn, downwind dive, and low-altitude turn phases. The transfer logic among the phase states is illustrated in Figure 6. The switch variables are v z and the heading angle. The state machine used simple switch rules to generate the Rayleigh cycle, and no further optimization was thus required.

Bio-Inspired Strategy
Although the optimization method affords an optimal means of generating a Rayleigh Cycle [18][19][20][21][22][23][24][25][26][27][28], the computation procedure, as noted above, is too time-consuming for real-time applications. Hence, we used a state machine to imitate the four phases of the Rayleigh cycle of an albatross, namely, the upwind climb, high-altitude turn, downwind dive, and low-altitude turn phases. The transfer logic among the phase states is illustrated in Figure 6. The switch variables are vz and the heading angle. The state machine used simple switch rules to generate the Rayleigh cycle, and no further optimization was thus required. Figure 6. Transfer logic among the four phases in the bio-inspired strategy for implementing the Raleigh cycle.

Upwind Climb
When this state is active, the kinetic energy of the aircraft is transformed into gravitational potential energy. When the vertical speed vz falls below zero, as shown in Figure 6, the next state of high-altitude turn becomes active. The control input in this state is where (CL)opt is the optimum lift coefficient in Figure 7, the choice of which enables the determination of the best climbing performance, as described in Section 5.1.1.

High-Altitude Turn
In this state, the aircraft loses energy. To decrease the energy loss, the aircraft needs to move on to the state of downwind dive as soon as possible. Because the airspeed is quite low, the maximum bank angle and lift coefficient are given by The downwind dive state becomes active when the aircraft heading becomes the same as the wind direction.

Downwind Dive
Similar to the upwind climb state of energy harvest, the control input in the downwind dive state is also given by Equation (47). When the z acceleration of the aircraft begins to increase, The low-altitude turn state then becomes active.

Upwind Climb
When this state is active, the kinetic energy of the aircraft is transformed into gravitational potential energy. When the vertical speed v z falls below zero, as shown in Figure 6, the next state of high-altitude turn becomes active. The control input in this state is where (C L ) opt is the optimum lift coefficient in Figure 7, the choice of which enables the determination of the best climbing performance, as described in Section 5.1.1.

High-Altitude Turn
In this state, the aircraft loses energy. To decrease the energy loss, the aircraft needs to move on to the state of downwind dive as soon as possible. Because the airspeed is quite low, the maximum bank angle and lift coefficient are given by The downwind dive state becomes active when the aircraft heading becomes the same as the wind direction.

Downwind Dive
Similar to the upwind climb state of energy harvest, the control input in the downwind dive state is also given by Equation (47). When the z acceleration of the aircraft begins to increase, The low-altitude turn state then becomes active.

Low-Altitude Turn
Although energy is lost in this state, like in the high-altitude turn state, the airspeed is quite high. If a large lift coefficient is assigned, the corresponding large drag coefficient would produce a relatively large drag, which will result in the consumption of a large amount of energy. Consequently, a discount is set for the maximum lift coefficient to compromise for fast turning: where (C L ) min is the minimum lift coefficient, and λ is the discount rate, which is a tunable parameter requiring tuning in an experiment.

Initial State
The initial state can be defined from an energy perspective. If the initial conditions are v x > 0 and v z > 0, the upwind climb state would be active. If the initial velocity components v x and v z are both negative, the downwind dive state would be the initial active state. If v x > 0 and v z < 0, the high-altitude turn state would be the initial active state. If v x < 0 and v z > 0, the low-altitude turn state would be the initial active state. This logic is summarized in Table 1.

Optimum Lift Coefficient
The drag coefficient is generally expressed as [43] C where n and C D0 are parameters that depend on the airfoil itself. When an aircraft is in flight, its angle of attack can be adjusted by the elevators to vary the lift coefficient C L . The drag coefficient C D varies with C L according to Equation (51), and C L can thus be considered as a control input. Substituting Equations (51) and (11) into Equation (20) yields Equation (52) expresses Π A as a function of C L . Minimizing Π A optimizes C L and thus minimizes the climbing criterion in Equation (17). In the case of an albatross, C D0 and n are respectively 0.033 and 0.019 [22]. The function Π A (C L ) for this case is plotted in Figure 7. Figure 7 indicates that the optimal lift coefficient C L for an albatross is approximately 0.1. The lift coefficient corresponding to the maximum lift-to-drag ratio can be calculated using [21] as C D0 /n = 1.32 Appl As can be observed, these two C L values significantly differ. To understand the reason for the small optimal C L , we will consider a simple calculation. For C L = 1.32, we have L/D = 20, h = 1.08, and Thus, from Equation (20), Π A = 1.24.
Further, for C L = 0.1, we have L/D = 3, h = 1.8, and From Equation (20), Π A = 0.58. Although the maximum L/D minimizes h, the other item in Equation (20) has a much larger effect on the calculated value of Π A . In the general soaring problem without a wind gradient, the maximum lift-to-drag ratio is optimal for the soaring range [44]. The situation is, however, different in the presence of a wind gradient. An intuitive explanation of the difference can be presented from an energy viewpoint. In the absence of a wind gradient, maintaining the maximum lift-to-drag ratio enables energy saving, based on traditional flight mechanics. When there is a wind gradient, the aircraft can harvest extra energy from it, and a smaller lift coefficient would cause the velocity to increase, because of the higher velocity needed to generate the same lift to offset the weight. Thus, the aircraft crosses more wind layers during the same period, acquiring energy more efficiently, although it also saves less energy compared to that at the maximum lift-to-drag ratio. It is thus necessary to make a compromise between harvesting and saving energy. Therefore, the calculated optimal lift coefficient in Figure 4 is smaller than the corresponding maximum lift-to-drag ratio.  As can be observed, these two CL values significantly differ. To understand the reason for the small optimal CL, we will consider a simple calculation.

Wing Loading
According to Equation (22), a large wing loading produces a small wing loading fraction Π S . Therefore, a larger wing loading is preferred for effortless climbing. This may seem incredible considering that traditional aeronautics textbooks suppose that a small wing loading and high aspect ratio wing are favorable to climbing, with the former inducing a low speed, which reduces the drag and power consumption. However, the important consideration here is the energy harvested from the wind gradient for unpowered flight, rather than power conservation in powered flight. The speed increases due to the larger wing loading, which results in the aircraft crossing more wind layers and collecting more energy during a given period. This is intuitive from an energy viewpoint.
A large aircraft generally has a larger wing load compared to a small one. McCormick [44] proposed a model for estimating the wing loading based on a cubic scaling law. He described the relationship by defining the following upper and lower boundaries: Noth [45] validated this estimate using an aircraft database and the Great Flight diagram. Many other models have been proposed for different types of aerial vehicles [46]. Although a large wing loading is favorable for climbing, the applicable loading is limited by the structural strength, cruising speed, and mission configuration.

Sensitivity Analysis of Aircraft Fraction
To compare the sensitivity of the two aircraft fractions in Equation (19), we define the logarithm of the aircraft fractions as follows: The aircraft fraction can thus be obtained in the form of a summation from Equation (19), as follows: Smaller values of these three parameters are also preferable, considering the multiplication criterion in Equation (17). Five example aircrafts were used to analyze the sensitivity of the fractions. The relevant parameters are presented in Table 2. The parameters for four of the aircrafts were obtained from references, while those for the National University of Defense Technology (NUDT) solar aircraft were measured and estimated. Σ A was calculated with respect to C L for each aircraft using Equation (52). The results are plotted in Figure 8. As can be observed, the plots for Avatar P2, SBXC, and Sky Sailor are very close, as well as those for the albatross and NUDT Solar Aircraft.
It can also be seen from Figure 8 that all of the optimal lift coefficients are quite close to zero. The practical lift coefficient during flight is generally higher than the optimal value, and we therefore consider the operational Σ A to range between the minimum value and the value corresponding to the maximum lift coefficient. On this basis, we calculated the values of Σ A , Σ S , and Σ a for the different aircraft, as plotted in Figure 9.
corresponding to the maximum lift coefficient. On this basis, we calculated the values of ΣA, ΣS, and Σa for the different aircraft, as plotted in Figure 9.  In Figure 9, the wide and narrow bars respectively indicate the maximum and minimum values. The lift coefficient can be adjusted to assign an aircraft fraction between the minimum and maximum values, and the aerodynamic configuration and wing loading can also be redesigned to modify the aircraft fraction. From Figure 9, it can be seen that the variation of ΣA among the different aircrafts is as large as that of the wing loading fraction. This suggests that the aerodynamic fraction is generally as important as the wing loading fraction. Neither of these fractions can thus be ignored in a bid to enhance the sustainable climbing performance in a gradient wind.  corresponding to the maximum lift coefficient. On this basis, we calculated the values of ΣA, ΣS, and Σa for the different aircraft, as plotted in Figure 9.  In Figure 9, the wide and narrow bars respectively indicate the maximum and minimum values. The lift coefficient can be adjusted to assign an aircraft fraction between the minimum and maximum values, and the aerodynamic configuration and wing loading can also be redesigned to modify the aircraft fraction. From Figure 9, it can be seen that the variation of ΣA among the different aircrafts is as large as that of the wing loading fraction. This suggests that the aerodynamic fraction is generally as important as the wing loading fraction. Neither of these fractions can thus be ignored in a bid to enhance the sustainable climbing performance in a gradient wind. In Figure 9, the wide and narrow bars respectively indicate the maximum and minimum values. The lift coefficient can be adjusted to assign an aircraft fraction between the minimum and maximum values, and the aerodynamic configuration and wing loading can also be redesigned to modify the aircraft fraction. From Figure 9, it can be seen that the variation of Σ A among the different aircrafts is as large as that of the wing loading fraction. This suggests that the aerodynamic fraction is generally as important as the wing loading fraction. Neither of these fractions can thus be ignored in a bid to enhance the sustainable climbing performance in a gradient wind.

Minimum Gradient and Environmental Effects
According to the inequality in Equation (31), the minimum gradient for climbing (G w ) min can be readily determined as follows: This promises to facilitate the choice of appropriate parameters and locations for dynamic soaring, and indicates the dependency of the aircraft on the wind field.
The gravitational acceleration g and air density ρ vary little near the earth's surface. However, the situation is quite different near space and on other planets. The near space altitude is higher than that at which traditional aircrafts fly, but lower than that at which traditional astrovehicles usually orbit. It ranges between 20 and 100 km above sea level. The effects of gravity gradually weaken and the air thins with increasing altitude, as indicated in Figure 10. The data presented in the figure were obtained using the International Standard Atmosphere model [47]. The lower values of g and ρ lower the minimum gradient near space, relaxing the required wind field for climbing. An advantage of near space is that the wind field has a steady gradient, especially between 40 and 60 km above sea level ( Figure 11). The horizontal wind data in Changsha (28 • 11 49" N, 112 • 58 42" E) presented in Figure 11 were obtained by the MATLAB function atmoshwm07 using the Naval Research Laboratory Horizontal Window Model [48]. Gao et al. [49] attempted to use the wind shear at the bottom of near space to achieve high-altitude long-endurance flight.

Minimum Gradient and Environmental Effects
According to the inequality in Equation (31), the minimum gradient for climbing (Gw)min can be readily determined as follows: This promises to facilitate the choice of appropriate parameters and locations for dynamic soaring, and indicates the dependency of the aircraft on the wind field.
The gravitational acceleration g and air density ρ vary little near the earth's surface. However, the situation is quite different near space and on other planets. The near space altitude is higher than that at which traditional aircrafts fly, but lower than that at which traditional astrovehicles usually orbit. It ranges between 20 and 100 km above sea level. The effects of gravity gradually weaken and the air thins with increasing altitude, as indicated in Figure 10. The data presented in the figure were obtained using the International Standard Atmosphere model [47]. The lower values of g and ρ lower the minimum gradient near space, relaxing the required wind field for climbing. An advantage of near space is that the wind field has a steady gradient, especially between 40 and 60 km above sea level ( Figure 11). The horizontal wind data in Changsha (28°11′49″ N, 112°58′42″ E) presented in Figure 11 were obtained by the MATLAB function atmoshwm07 using the Naval Research Laboratory Horizontal Window Model [48]. Gao et al. [49] attempted to use the wind shear at the bottom of near space to achieve high-altitude long-endurance flight. Let us consider a flight at an altitude of 50 km above sea level. The air density at this altitude is very low at approximately 10 −4 kg/m 3 , compared to 1 kg/m 3 at sea level. There is, however, little change in the gravitational acceleration. From Equation (35), (G w ) min | 50 km above sea level (G w ) min | at sea level = (gρ) 50 km above sea level (gρ) at sea level ≈ 0.01 This indicates that the minimum gradient for climbing at an altitude of 50 km is 1% of that near the earth's surface. Climbing at the high altitude thus requires a relatively weak wind gradient, implying the potential for dynamic soaring over a large area.
NASA [50,51], ESA [52], and JAXA [53] are considering using a fixed-wing aircraft for the exploration of Mars. Early in the 1970s, NASA built and tested the Mini-Sniffer [54], which had a wingspan of 6.7 m. The Sky-Sailor [45,52] is an ultra-lightweight solar autonomous model aircraft designed by the Autonomous System Lab of EPFL under a contract with ESA, for scientific missions to Mars. The NUDT also developed a prototype solar aircraft for Mars exploration ( Figure 12). The exploitation of a wind shear on Mars may be a possible approach to extending the flight endurance of Mars exploration aircrafts. The gravitational acceleration on Mars is less than 4/9 of that on Earth, and the atmosphere is also much thinner ( Figure 13). The minimum gradient for This indicates that the minimum gradient for climbing at an altitude of 50 km is 1% of that near the earth's surface. Climbing at the high altitude thus requires a relatively weak wind gradient, implying the potential for dynamic soaring over a large area.
NASA [50,51], ESA [52], and JAXA [53] are considering using a fixed-wing aircraft for the exploration of Mars. Early in the 1970s, NASA built and tested the Mini-Sniffer [54], which had a wingspan of 6.7 m. The Sky-Sailor [45,52] is an ultra-lightweight solar autonomous model aircraft designed by the Autonomous System Lab of EPFL under a contract with ESA, for scientific missions to Mars. The NUDT also developed a prototype solar aircraft for Mars exploration (Figure 12). This indicates that the minimum gradient for climbing at an altitude of 50 km is 1% of that near the earth's surface. Climbing at the high altitude thus requires a relatively weak wind gradient, implying the potential for dynamic soaring over a large area.
NASA [50,51], ESA [52], and JAXA [53] are considering using a fixed-wing aircraft for the exploration of Mars. Early in the 1970s, NASA built and tested the Mini-Sniffer [54], which had a wingspan of 6.7 m. The Sky-Sailor [45,52] is an ultra-lightweight solar autonomous model aircraft designed by the Autonomous System Lab of EPFL under a contract with ESA, for scientific missions to Mars. The NUDT also developed a prototype solar aircraft for Mars exploration ( Figure 12). The exploitation of a wind shear on Mars may be a possible approach to extending the flight endurance of Mars exploration aircrafts. The gravitational acceleration on Mars is less than 4/9 of that on Earth, and the atmosphere is also much thinner ( Figure 13). The minimum gradient for The exploitation of a wind shear on Mars may be a possible approach to extending the flight endurance of Mars exploration aircrafts. The gravitational acceleration on Mars is less than 4/9 of that on Earth, and the atmosphere is also much thinner ( Figure 13). The minimum gradient for sustainable climbing on Mars is thus considered to be smaller than that on Earth. In addition, the wind field on Mars significantly varies with the altitude [55], which is also favorable to climbing. The data presented in Figure 13 were obtained using the standard Mars-GRAM 2005 [55]. , 7, 1061 sustainable climbing on Mars is thus considered to be smaller than that on Earth. In addition, the wind field on Mars significantly varies with the altitude [55], which is also favorable to climbing. The data presented in Figure 13 were obtained using the standard Mars-GRAM 2005 [55].

Rayleigh Cycle Simulation Results
The parameters of the albatross in Table 2 were used for the simulation of its Raleigh cycle. The wind shear model, which was in accordance with Equation (1), was obtained from MIL-F-8785C [40]. W6 was set to 15 m/s, which is a common value on Earth [48]. A dynamic model block was built and used to perform the simulation with the aid of the Simulink and Stateflow software (MathWorks, Natick, MA, USA). Figure 14 shows a snapshot of the Simulink model, while Figure 15 shows the Stateflow chart corresponding to the Mode Switch block in Figure 14.

Rayleigh Cycle Simulation Results
The parameters of the albatross in Table 2 were used for the simulation of its Raleigh cycle. The wind shear model, which was in accordance with Equation (1), was obtained from MIL-F-8785C [40]. W 6 was set to 15 m/s, which is a common value on Earth [48]. A dynamic model block was built and used to perform the simulation with the aid of the Simulink and Stateflow software (MathWorks, Natick, MA, USA). Figure 14 shows a snapshot of the Simulink model, while Figure 15 shows the Stateflow chart corresponding to the Mode Switch block in Figure 14.
Appl. Sci. 2017, 7, 1061 19 of 25 sustainable climbing on Mars is thus considered to be smaller than that on Earth. In addition, the wind field on Mars significantly varies with the altitude [55], which is also favorable to climbing. The data presented in Figure 13 were obtained using the standard Mars-GRAM 2005 [55].

Rayleigh Cycle Simulation Results
The parameters of the albatross in Table 2 were used for the simulation of its Raleigh cycle. The wind shear model, which was in accordance with Equation (1), was obtained from MIL-F-8785C [40]. W6 was set to 15 m/s, which is a common value on Earth [48]. A dynamic model block was built and used to perform the simulation with the aid of the Simulink and Stateflow software (MathWorks, Natick, MA, USA). Figure 14 shows a snapshot of the Simulink model, while Figure 15 shows the Stateflow chart corresponding to the Mode Switch block in Figure 14.   The Stateflow chart in Figure 15 realizes the switch logic in Figure 6, with the Controller block particularly realizing the algorithm described in Section 4.3. In Figure 15, the transition logic between the climbing state and the high-altitude turn is vz < 1, because 1 is a margin that ensures that the aircraft enters the high-altitude turn state before it begins to dive. To minimize the energy loss, the discount rate λ should be tuned. When the discount rate is set to 0, the low-altitude turn requires about 10 s, and the final velocity is so small that the aircraft is incapable of connecting to the next upwind climb phase. When the discount rate λ is set to 1, although the time required for the turn is significantly shortened, the energy consumption also increases. The bisection method was used to iteratively tune the discount rate in the present work, by which a value of 0.9 was found to be suitable.
For a given aircraft, the lift coefficient CL varies with the angle of attack, which can be approximated as the angle between the total velocity vector and the horizontal plane of the aircraft. The angle of attack is thus the actual control input. Owing to the availability of very limited data on the aerodynamics of albatrosses, the specific relationship between their angle of attack and lift coefficient is unclear. In the present study, an albatross was simply considered as an aircraft, and its lift coefficient was assumed to be proportional to the angle of attack for small values of the latter. Although this is consistent with regular aerodynamics theory [43], it was still difficult to determine the proportionality. Nevertheless, with the assumption, the lift coefficient could be directly controlled by the angle of attack in the simulation. Figures 16-18 compares the simulation results for the proposed bio-inspired strategy with those for the traditional GPOPS optimized trajectory. The figures cover the same period in the entire time series for the two methods, coinciding with a full cycle of the GPOPS method, and a little longer than a cycle generated by the bio-inspired strategy. The yellow background in Figure 16 exactly covers a cycle of the bio-inspired strategy. The figure shows the variations in the energy and altitude with time. Figure 17 shows the variations in the control inputs with time. As can be observed, the computed control inputs for the GPOPS method form smooth signals that also comprise four phases corresponding to those of the Rayleigh cycle, although the phase switch boundaries are not as obvious as those of the bio-inspired strategy. In the upwind climb and downwind dive segments, the bank angles for both methods are very small, sometimes approaching zero. Conversely, the bank angles during the high-altitude and low-altitude turns are large. In addition, the lift coefficients for both methods are maximum during the high-altitude turn, and minimum during the upwind climb and downwind dive. In the low-altitude turn segment, both methods visibly compromise on the assigned lift coefficient, opting for medium values. The Stateflow chart in Figure 15 realizes the switch logic in Figure 6, with the Controller block particularly realizing the algorithm described in Section 4.3. In Figure 15, the transition logic between the climbing state and the high-altitude turn is v z < 1, because 1 is a margin that ensures that the aircraft enters the high-altitude turn state before it begins to dive. To minimize the energy loss, the discount rate λ should be tuned. When the discount rate is set to 0, the low-altitude turn requires about 10 s, and the final velocity is so small that the aircraft is incapable of connecting to the next upwind climb phase. When the discount rate λ is set to 1, although the time required for the turn is significantly shortened, the energy consumption also increases. The bisection method was used to iteratively tune the discount rate in the present work, by which a value of 0.9 was found to be suitable.
For a given aircraft, the lift coefficient C L varies with the angle of attack, which can be approximated as the angle between the total velocity vector and the horizontal plane of the aircraft. The angle of attack is thus the actual control input. Owing to the availability of very limited data on the aerodynamics of albatrosses, the specific relationship between their angle of attack and lift coefficient is unclear. In the present study, an albatross was simply considered as an aircraft, and its lift coefficient was assumed to be proportional to the angle of attack for small values of the latter. Although this is consistent with regular aerodynamics theory [43], it was still difficult to determine the proportionality. Nevertheless, with the assumption, the lift coefficient could be directly controlled by the angle of attack in the simulation. Figures 16-18 compares the simulation results for the proposed bio-inspired strategy with those for the traditional GPOPS optimized trajectory. The figures cover the same period in the entire time series for the two methods, coinciding with a full cycle of the GPOPS method, and a little longer than a cycle generated by the bio-inspired strategy. The yellow background in Figure 16 exactly covers a cycle of the bio-inspired strategy. The figure shows the variations in the energy and altitude with time. Figure 17 shows the variations in the control inputs with time. As can be observed, the computed control inputs for the GPOPS method form smooth signals that also comprise four phases corresponding to those of the Rayleigh cycle, although the phase switch boundaries are not as obvious as those of the bio-inspired strategy. In the upwind climb and downwind dive segments, the bank angles for both methods are very small, sometimes approaching zero. Conversely, the bank angles during the high-altitude and low-altitude turns are large. In addition, the lift coefficients for both methods are maximum during the high-altitude turn, and minimum during the upwind climb and downwind dive. In the low-altitude turn segment, both methods visibly compromise on the assigned lift coefficient, opting for medium values.  In Figure 16, the indicated boundaries of the four phases (upwind climb, high-altitude turn, downwind dive, and low-altitude turn) are for the bio-inspired strategy. The phase switches of the continuous GPOPS optimized trajectory are not obvious. The plots for both methods in Figure 16 validate the Rayleigh cycle energy theory discussed in Section 4.1. Energy is harvested during the upwind climb and downwind dive, whereas energy is lost during the low-altitude and high-altitude turns. However, the total energy changes during a cycle for the two methods are positive. For the GPOPS method, the energy increase during a cycle is about 1300 J, accompanied by  In Figure 16, the indicated boundaries of the four phases (upwind climb, high-altitude turn, downwind dive, and low-altitude turn) are for the bio-inspired strategy. The phase switches of the continuous GPOPS optimized trajectory are not obvious. The plots for both methods in Figure 16 validate the Rayleigh cycle energy theory discussed in Section 4.1. Energy is harvested during the upwind climb and downwind dive, whereas energy is lost during the low-altitude and high-altitude turns. However, the total energy changes during a cycle for the two methods are positive. For the GPOPS method, the energy increase during a cycle is about 1300 J, accompanied by In Figure 16, the indicated boundaries of the four phases (upwind climb, high-altitude turn, downwind dive, and low-altitude turn) are for the bio-inspired strategy. The phase switches of the continuous GPOPS optimized trajectory are not obvious. The plots for both methods in Figure 16 validate the Rayleigh cycle energy theory discussed in Section 4.1. Energy is harvested during the upwind climb and downwind dive, whereas energy is lost during the low-altitude and high-altitude turns. However, the total energy changes during a cycle for the two methods are positive. For the GPOPS method, the energy increase during a cycle is about 1300 J, accompanied by an altitude increase of about 17 m. In the case of the bio-inspired strategy, the energy and altitude increases during a cycle are about 900 J and 12 m, respectively. The GPOPS method produces an optimized solution, hence its larger cyclic energy and altitude increases. Nevertheless, the proposed bio-inspired strategy is simpler and affords a significant saving in the computation time. The larger cyclic energy increase of the GPOPS method is afforded by its lower energy loss during the high-altitude turn. Figure 18 shows phase portraits of the GPOPS and proposed bio-inspired methods. The plot for each method covers a cycle in the velocity space. The altitude corresponding to each point is color-coded, and the four phases of the Rayleigh cycle of the bio-inspired strategy are indicated by numbers, where the 1-2, 2-3, 3-4, and 4-1 segments respectively correspond to the upwind climb, high-altitude turn, downwind dive, and low-altitude turn phases. , 7, 1061 an altitude increase of about 17 m. In the case of the bio-inspired strategy, the energy and altitude increases during a cycle are about 900 J and 12 m, respectively. The GPOPS method produces an optimized solution, hence its larger cyclic energy and altitude increases. Nevertheless, the proposed bio-inspired strategy is simpler and affords a significant saving in the computation time. The larger cyclic energy increase of the GPOPS method is afforded by its lower energy loss during the high-altitude turn. Figure 18 shows phase portraits of the GPOPS and proposed bio-inspired methods. The plot for each method covers a cycle in the velocity space. The altitude corresponding to each point is color-coded, and the four phases of the Rayleigh cycle of the bio-inspired strategy are indicated by numbers, where the 1-2, 2-3, 3-4, and 4-1 segments respectively correspond to the upwind climb, high-altitude turn, downwind dive, and low-altitude turn phases.  Table 3 summarizes the computation time for GPOPS and the bio-inspired strategy. GPOPS generally consumes about 15 s to compute a dynamic soaring trajectory in a personal computer with an Intel ® Core TM i3-2100 CPU running at 3.10 GHz, while the bio-inspired strategy consumes 2 ms because of the simple judging logic (the computation time can be checked in the Simulink report interface). This comparison shows that the bio-inspired strategy is more efficient at real-time applications. Table 3. Computation times for GPOPS and the bio-inspired strategy.

Method
GPOPS Bio-Inspired Strategy Computation time 15 s 2 ms

Conclusions
A criterion for sustainable climbing of an unpowered aircraft in a gradient wind was proposed and verified. The criterion is based on the product of the aircraft and environmental fractions, where the aircraft fraction is composed of the aerodynamic and wing loading fractions. The criterion was used to develop a method for calculating the optimal lift coefficient from the aerodynamic fraction. An expression of the minimum wind gradient required for unpowered flight was also derived from  Table 3 summarizes the computation time for GPOPS and the bio-inspired strategy. GPOPS generally consumes about 15 s to compute a dynamic soaring trajectory in a personal computer with an Intel ® Core TM i3-2100 CPU running at 3.10 GHz, while the bio-inspired strategy consumes 2 ms because of the simple judging logic (the computation time can be checked in the Simulink report interface). This comparison shows that the bio-inspired strategy is more efficient at real-time applications. Table 3. Computation times for GPOPS and the bio-inspired strategy.

GPOPS Bio-Inspired Strategy
Computation time 15 s 2 ms

Conclusions
A criterion for sustainable climbing of an unpowered aircraft in a gradient wind was proposed and verified. The criterion is based on the product of the aircraft and environmental fractions, where the aircraft fraction is composed of the aerodynamic and wing loading fractions. The criterion was used to develop a method for calculating the optimal lift coefficient from the aerodynamic fraction. An expression of the minimum wind gradient required for unpowered flight was also derived from the environmental fraction, and the feasibility of implementing the flight in near space and on Mars was analyzed. A sensitivity analysis of the aircraft fraction indicated that the aerodynamic fraction was generally as important as the wing loading fraction. Further, the energy harvesting boundaries of the Rayleigh cycle for unpowered flight was investigated in the velocity space, and the maximum power was determined. On the basis of the Rayleigh cycle energy theory, a bio-inspired Rayleigh cycle guidance strategy for unpowered flight was proposed, involving the tuning of a single parameter. A comparison of the simulation results for the proposed strategy and the traditional optimization method confirmed the effectiveness of the proposed strategy for real-time applications.