An Efficient Non-Linear Application Algorithm Predictive Model for a Multi Aircraft Landing Dynamic System AIRLADYS R2019A+

Abstract

The aim of this paper is to set up an efficient nonlinear application algorithm predictive model for a multi aircraft landing dynamic system called “Aircraft Landing Dynamic System, Release 2019A+ version “AIRLADYS R2019A+”. This programming software combines dynamic programming technic for mathematical computing and optimisation run under AMPL and KNITRO Solver. It uses also a descriptive programming technic for software design. The user interfaces designed in Glade are saved as XML, and by using the GtkBuilder GTK+ object these can be loaded by applications dynamically as needed. By using GtkBuilder, Glade XML files can be used in numerous programming languages including C, C++, C#, Java, Perl, Python, AMPL, etc. Glade is Free Software released under the GNU GPL License. By these tools, the solved problem is a mathematical modelization problem as a non-convex optimal control governed by ordinary non-linear differential equations. The dynamic programming technic is applied because it is a sufficiently high order and it does not require computation of the partial derivatives of the aircraft dynamic. This application will be coded with Linux system on 64 bit operating system, but it can also be run on the windows system. High running performances are obtained with results giving feasible trajectories with a robust optimizing of the objective function.

Share and Cite:

Nahayo, F. and Khardi, S. (2020) An Efficient Non-Linear Application Algorithm Predictive Model for a Multi Aircraft Landing Dynamic System AIRLADYS R2019A+. Open Journal of Optimization, 9, 15-26. doi: 10.4236/ojop.2020.92002.

1. Introduction

In this work, an efficient nonlinear application algorithm predictive model for a multi aircraft landing dynamic system is developed while maintaining a reliable evolution of the flight procedures of aircraft dynamic system on approach. The aircraft is landing successively on one runway while maintaining the separation constraints [1].

This programming software combines dynamic programming technic for mathematical computing and optimisation run under AMPL and KNITRO Solver [2] [3]. It uses also a descriptive programming technic for software design. The user interfaces designed in Glade are saved as XML, and by using the GtkBuilder GTK+ object these can be loaded by applications dynamically as needed. By using GtkBuilder, Glade XML files can be used in numerous programming languages including C, C++, C#, Java, Perl, Python, AMPL, etc. Glade is Free Software released under the GNU GPL License.

This application will be coded with Linux system on 64 bit operating system. High running performances are obtained with results giving feasible trajectories with a robust optimizing of the objective function. The user interfaces are designed in Glade by using the GtkBuilder GTK+ object and this can be loaded by applications dynamically as needed.

2. Mathematical Modeling of the Aircraft Dynamic System

2.1. Aircraft Dynamic Equations

The equations of 3D-motion of each aircraft A i , i { 1,2 } read [4] [5] [6] [7]:

{ r ˙ i = E A C E 2 { r i q i ( B C ) + E p i q i + 1 2 ρ i S l V a i 2 C l i + j = 1 2 F j [ y M i j b cos β m i j sin α m i j z M i j b sin β m i j ] } , + A A C E 2 { p i q i ( A B ) E r i q i + 1 2 ρ i S l V a i 2 C n i + j = 1 2 F j [ x M i j b sin β m i j y M i j b cos β m i j cos α m i j ] } , X ˙ G i = V a i cos γ a i cos χ a i + u w , Y ˙ G i = V a i cos γ a i sin χ a i + v w , Z ˙ G i = V a i sin γ a i + w w , ϕ ˙ i = p i + q i sin ϕ i tan θ i + r i cos ϕ i tan θ i , θ ˙ i = q i cos ϕ i r i sin ϕ i , ψ ˙ i = sin ϕ i cos θ i q i + cos ϕ i cos θ i r i , (1)

where j { 1,2 } stands for the first and second engine of each aircraft i, the expressions A = I x x , B = I y y , C = I z z , E = I x z are the inertia moments of the aircraft, ρ i is the air density, S is the aircraft reference area, l is the aircraft reference length, g is the acceleration due to gravity, C D i = C D 0 + k C L i 2 is the

drag coefficient, C y i = C y β β + C y p p l V + C y r r l V + C Y δ l δ l i + C Y δ n δ n i is the lateral forces coefficient, C L i = C L α ( α a α a 0 ) + C L δ m δ m i + C L M M i + C L q q a b l V is the lift coefficient, C l i = C l β β + C l p p l V + C l r r l V + C l δ l δ l i + C l δ n δ n i is the rolling moment coefficient, C m i = C m 0 + C m α ( α α 0 ) + C m δ m δ m i is the pitching moment coefficient, C n i = C n β β + C n p p l V + C n r r l V + C n δ l δ l i + C n δ n δ n i is the yawing moment

coefficient, ( x M i j b , x M i j b , x M i j b ) is the position of the engine in the body frame, P 0 is the full thrust, ρ 0 is the atmospheric density at the ground, F = ( F x i , F y i , F z i ) is the propulsive force, V a i = ( u i , v i , w i ) is the aerodynamic speed, ( Δ A u i , Δ A v i , Δ A w i ) is the complementary acceleration, ( u w , v w , w w ) is the wind velocity, β m i j is the yaw setting of the engine and α m i j is the pitch setting of the engine. The mass change is reflected in the aircraft fuel consumption as described by E. Torenbeek [7] where the specific consumption is

C S R i = 2.01 × 10 5 ( Φ μ i K i η c ) Θ 5 η n ( 1 + η t f i λ ) G i + 0.2 M i 2 η d i η t f i λ ( 1 λ ) M i with the generating function G i = ( Φ K i η c ) ( 1 1.01 η i ν 1 ν ( K i + μ i ) ( 1 K i Φ η c η t ) ) ,

K i = μ i ( ε c ν 1 ν 1 ) , μ i = 1 + ν 1 2 M i 2 . The nomenclature of engine performance

variables are given by G i the gas generator power function, G0 the gas generator power function (static, sea level), K the temperature function of compression process, M i the flight Mach number, T4 the turbine Entry total Temperature, T0 the ambient temperature at sea level, T the flight temperature, while the nomenclature of engines yields is η c = 0.85 the isentropic compressor efficiency,

η d i = 1 1.3 ( 0.05 R e 1 5 ) 2 ( 0.5 M i ) 2 L D , the isentropic fan intake duct efficiency, L the

duct length, D the inlet diameter, Re the Reynolds number at the entrance of the nozzle, η f i = 0.86 3.13 × 10 2 M i the isentropic fan efficiency,

η i = 1 + η d i γ 1 2 M i 2 1 + γ 1 2 M i 2 the gas generator intake stagnation pressure ratio, η n = 0.97

the isentropic efficiency of expansion process in nozzle, η t = 0.88 the isentropic turbine efficiency η t f i = η t η f i , ε c the overall pressure ratio (compressor), ν the ratio of specific heats ν = 1.4 , λ the bypass ratio, μ i the ratio of stagnation to static temperature of ambient air, Φ the nondimensional turbine entry

temperature Φ = T 4 T and Θ the relative ambient temperature Θ = T T 0 . The

expressions α a i ( t ) , β a i ( t ) , θ i ( t ) , ψ i ( t ) , ϕ i ( t ) , V a i ( t ) , X G i ( t ) , Y G i ( t ) , Z G i ( t ) , p i ( t ) , q i ( t ) , r i ( t ) , m i ( t ) are respectively the attack angle, the aerodynamic sideslip angle, the inclination angle, the cup, the roll angle, the airspeed, the position vectors, the roll velocity of the aircraft relative to the earth, the pitch velocity of the aircraft relative to the earth, the yaw velocity of the aircraft relative to the earth and the aircraft mass. The system (2) could be written in a simplified form

d y i ( t ) d t = f i ( y i ( t ) , u i ( t ) ) , y i ( t ) = ( α a i ( t ) , β a i ( t ) , θ a i ( t ) , ψ a i ( t ) , ϕ i ( t ) , V a i ( t ) , X G i ( t ) , Y G i ( t ) , Z G i ( t ) , p i ( t ) , q i ( t ) , r i ( t ) , m i ( t ) ) u i ( t ) = ( δ l i ( t ) , δ m i ( t ) , δ n i ( t ) , δ x i ( t ) ) (2)

henceforth y i is called a state function and the expressions δ l i ( t ) , δ m i ( t ) , δ n i ( t ) , δ x i ( t ) are respectively the roll control, the pitch control, the yaw control and the thrust control. The dynamics relationship can be written as: y ˙ i ( t ) = f i ( y i , u i , t ) , t [ 0 , T ] , y i ( 0 ) = y i 0 [8].

2.2. The Optimal Objective Function Model

Let us define the quantity named the Sound Exposure Level “SEL” [4]:

S E L = 10 log [ t 10 0.1 L A 1 , d t ( t ) d t ]

where t is the noise event interval. [ t 10 , t 1 f ] and [ t 20 , t 2 f ] are the respective approach intervals for the first and second aircraft, the objective function is calculated as:

S E L G = 10 l o g { 1 t 2 f t 10 [ ( t 20 t 10 ) t 10 t 20 10 0.1 L A 1 ( t ) d t + ( t 1 f t 20 ) t 20 t 1 f 10 0.1 L A 1 ( t ) d t + ( t 1 f t 20 ) t 20 t 1 f 10 0.1 L A 2 ( t ) d t + ( t 2 f t 1 f ) t 1 f t 2 f 10 0.1 L A 2 ( t ) d t ] } , t [ t 10 , t 2 f ] (3)

where the cost function S E L G is the cumulated two-aircraft noise. Expressions L A 1 ( t ) , L A 2 ( t ) are equivalent and reflect the aircraft jet noise given by the formula:

L A 1 ( t ) = 141 + 10 log ( ρ 1 ρ i ) w + 10 log ( V e c ) 7.5 + 10 log s 1 + 3 log ( 2 s 1 π d 1 2 + 0.5 ) + 5 log τ 1 τ 2 + 10 log [ ( 1 v 2 v 1 ) m e + 1.2 ( 1 + s 2 v 2 2 s 1 v 1 2 ) 4 ( 1 + s 2 s 1 ) 3 ] 20 log R + Δ V + 10 log [ ( ρ i ρ I S A ) 2 ( c c I S A ) 4 ]

where v 1 is the jet speed at the entrance of the nozzle, v 2 the jet speed at the nozzle exit, τ 1 the inlet temperature of the nozzle, τ 2 the temperature at the nozzle exit, ρ i the density of air, ρ 1 the atmospheric density at the entrance of the nozzle, ρ I S A the atmospheric density at ground, s 1 the entrance area of the nozzle hydraulic engine, s 2 the emitting surface of the nozzle hydraulic engine, d 1 the inlet diameter of the nozzle hydraulic engine,

V e = v 1 [ 1 ( V / v 1 ) c o s ( α p ) ] 2 / 3 the effective speed ( α p is the angle between the axis of the motor and the axis of the aircraft), R the source observer distance, w

the exponent variable defined by w = 3 ( V e / c ) 3.5 0.6 + ( V e / c ) 3.5 1 , c the sound velocity

(m/s), m e the exhibiting variable depending on the type of aircraft:

m e = 1.1 s 2 s 1 , s 2 s 1 < 29.7 ; m e = 6.0 , s 2 s 1 29.7 , the term

Δ V = 15 log ( C D ( M c , θ ) ) 10 log ( 1 M cos θ ) , means the Doppler convection when C D ( M c , θ ) = [ ( 1 + M c cos θ ) 2 + 0.04 M c 2 ] , M the aircraft Mach Number,

M c the convection Mach Number: M c = 0.62 ( v 1 V cos ( α p ) ) / c , θ is the Beam angle. The objective formula above could be written in the following simplified form J G 12 ( y ( . ) , u ( . ) ) = t g ( y ( t ) , u ( t ) ) d t .

2.3. Constraints

The considered constraints concern aircraft flight speeds and altitudes, flight angles and control positions, energy constraint, aircraft separation, flight velocities of aircraft relative to the earth and the aircraft mass [6]. On the whole, the constraints come together under the relationship k 1 i ( y i , u i ) 0, k 2 i ( y i , u i ) 0 where

k 1 i ( y i , u i ) = ( α i ( t ) α i f , θ i ( t ) θ i f , ψ i ( t ) ψ i f , ϕ i ( t ) ϕ i f , V a i ( t ) V a i f , X G i ( t ) X G i f , Y G i ( t ) Y G i f , Z G i ( t ) Z G i f , p i ( t ) p i f , q i ( t ) q i f , r i ( t ) r i f , δ l i ( t ) δ l i f , δ m i ( t ) δ m i f , δ n i ( t ) δ n i f , δ x i ( t ) δ x i f , m i ( t ) m i f ) ,

k 2 i ( y i , u i ) = ( α i ( t ) α i 0 , θ i ( t ) θ i 0 , ψ i ( t ) ψ i 0 , ϕ i ( t ) ϕ i 0 , V a i ( t ) V a i 0 , X G i ( t ) X G i 0 , Y G i ( t ) Y G i 0 , Z G i ( t ) Z G i 0 , p i ( t ) p i 0 , q i ( t ) q i 0 , r i ( t ) r i 0 , δ l i ( t ) δ l i 0 , δ m i ( t ) δ m i 0 , δ n i ( t ) δ n i 0 , δ x i ( t ) δ x i 0 , m i ( t ) m i 0 ) .

2.4. The Aircraft Optimal Control Problem

The combination of the aircraft dynamic equation, the aircraft objective function and the aircraft flight constraints, the two-aircraft acoustic optimal control problem is given as follows:

{ min u U J G 12 ( y ( . ) , u ( . ) ) = t 10 t 1 f g 1 ( y 1 ( t ) , u 1 ( t ) , t ) d t + t 20 t 1 f g 12 ( y 1 ( t ) , u 1 ( t ) , y 2 ( t ) , u 2 ( t ) , t ) d t + t 20 t 2 f g 2 ( y 2 ( t ) , u 2 ( t ) , t ) d t + ϕ ( y ( t f ) ) y ( t ) = f ( u ( t ) , y ( t ) ) , u ( t ) = ( u 1 ( t ) , u 2 ( t ) ) , y ( t ) = ( y 1 ( t ) , y 2 ( t ) ) , k 1 i ( y i , u i ) 0 , k 2 i ( y i , u i ) 0 , t [ t 10 , t 2 f ] , t 10 = 0 , y ( 0 ) = y 0 , u ( 0 ) = u 0 (4)

where g 12 shows the aircraft coupling noise function and J G 12 is the SEL of the two A300-aircraft.

3. The Numerical Processing

The above system is an optimal control problem with mixed constraints on the state and control. In order to apply the formulation of Pontryagin, we rewrite directly this system as follows:

{ m i n J 12 ( y , u ) = t 0 t f g ( y ( t ) , u ( t ) , t ) d t y ˙ ( t ) = F ( y ( t ) , u ( t ) ) k 1 ( y , u , t ) 0, k 2 ( y , u , t ) 0, (5)

By applying optimization Bell man theory and the Runge-Kutta symplectic method, the following algorithm is developed [9] [10] [11].

Partitioned symplectic Algorithm of Runge-Kutta (SPRK) [4].

1) Subdividing the interval time [ t 0 , t f ] in N step h = t n + 1 t n = t f t 0 N , N is the maximum number of iteration.

2) For 0 n N ,

H u ( u k i , y k i , p k i ) = 0 y n + 1 = y n + h k = 1 s b k H p ( u n k , y n k ) , y ( t 0 ) = y 0 , y n k = y n + h j = 1 s a k j H p ( u n j , y n j ) , k = 1 , , s p n + 1 = p n h k = 1 s b ˜ k H y ( p n k , y n k ) , p N = ϕ ( t f ) ,

p n k = p n h j = 1 s a ˜ k j H y ( p n j , y n j ) , H ( y n k , p n k ) = H ( ψ ( y n k , p n k ) , y n k , p n k ) , a ˜ k j : = b j b j b k a j k , b j = b ˜ j , k 1 ( y n , u n k ) 0 , k 2 ( y n , u n k ) 0, display t n + 1 , y n + 1 , p n + 1 . (6)

3) Stop.

4. Graphic User Interface for an Efficient Nonlinear Application Algorithm Predictive Model for a Multi Aircraft Landing Dynamic System

The user interfaces is designed in Glade and by using the GtkBuilder GTK+ object, this can be loaded by applications dynamically as needed. By using GtkBuilder, Glade XML files can be used in numerous programming languages including C, C++, C#, Vala, Java, Perl, Python, AMPL and others. Glade is Free Software released under the GNU GPL License. Figure 1 shows AIRLADY SR2018A+ Graphic User Interface and all the menu toolbar functions programmed for its running and exploitation.

Figure 2(a) and Figure 2(b) show the aircraft trajectories and speeds characterized by a part of constant flight level followed by a continuous descent till the aircraft touch point. Constraints on speeds are considered, allowing a subsequent landing on the same runway. The maximum altitudes considered are 3500 m and 4100 m for the first and second aircraft. The approach duration is 600 s for the first aircraft and 690 s for the second. The aircraft speeds decrease from 200 m/s to 69 m/s. This shows the aircraft trajectory resulting from the two trajectories combination.

Figure 1. GNU general public license: graphic user interface.

(a)(b)

Figure 2. Aircraft altitudes and speeds. (a) First aircraft flight path and speed; (b) Second aircraft flight path and speed.

Figure 3(a) and Figure 3(b) show the two-Aircraft throttle and roll control versus time. The optimal standards procedures are confirmed as described in operational aircraft flights paths.

Figure 4(a) and Figure 4(b) show the two-Aircraft pitch and yaw control versus time. The optimal standards procedures are confirmed as described in operational aircraft flights paths.

Figure 5 shows also the two-aircraft flight-angles and throttles evolution versus time as recommended by ICAO during aircraft landing. As specified in this figure, the aircraft roll angles oscillate around zero. The flight-path angles are negative and bang-bang. They keep the recommended position for aircraft landing procedures. The attack angles stand between 2˚ and 20˚. Since the trajectory of the aircraft is aligned with the runway, the yaw angle are small as shown in Figure 5(b) and Figure 5(d).

Figure 6 shows the noise levels when the optimization is applied and the solutions obtained. The observation positions are (−20,000 m, −20,000 m, 0 m) for AONL1, (−19,800 m, −19,800 m, 0 m) for AONL2, ..., (−200 m, −200 m, 0 m) for AONL10. In this figure, AONL means Aircraft Optimal Noise Level. As specified, noise levels increase and are maximum when the observation point lies below the aircraft. Noise levels decrease gradually as the aircraft moves away from the observation point. By comparison, this result is also close to standard values of jet noise on approach as shown by Harvey [8] [12] [13] [14] [15] [16].

5. Conclusion

This paper develops mathematical solving methods of an optimal control dynamic system of two aircrafts landing successively in one run way. Numerical

(a)(b)

Figure 3. Aircraft throttle and roll control. (a) First aircraft throttle and roll control versus time; (b) Second aircraft throttle and roll control versus time.

(a)(b)

Figure 4. Aircraft pitch and yaw control. (a) First aircraft pitch and yaw ontrol versus time; (b) Second aircraft pitch and yaw control versus time.

results are found through a robust software application. The Numerical program had been coded with Ubuntu Linux system. The programming software combines dynamic programming technic for mathematical computing and optimisation run under AMPL and KNITRO Solver. It uses also a descriptive programming technic for software design. High running performances are obtained with results giving feasible trajectories with a robust optimizing of the objective function. The user interfaces had been designed in Glade by using the

(a) (b) (c) (d)

Figure 5. A1 aircraft flight-path angles. (a) A1 attack and flight-path angle; (b) A1 Roll and yaw angle versus time; (c) A2 attack and flight-path angle versus time; (d) A2 Roll and yaw angle versus time.

(a) (b) (c) (d)

Figure 6. Aircraft optimal noise levels. (a) Noise levels versus time; (b) Noise levels versus time; (c) Noise levels versus time; (d) Noise levels versus time.

GtkBuilder GTK+ object which can be loaded by applications dynamically as needed through many numerous programming languages.

6. Future Perspective of This Application

Challenges are so many in computer sciences tools and in applied mathematics.

Considering the operational procedures and the many types of aircraft in operating society, this application must grow up and include all those considerations. With this, this application is under construction...

Authors are ready to produce “Aircraft Landing Dynamic System, Release 2020 B+” rsion “AIRLADYS R2020B+” before the end of 2020 year.

Acknowledgements

This paper was developed in Institute of Applied Statistics, University of Burundi and supported by the Government of Burundi through the Ministry of Higher Education and Scientific Research and the French Institute of Science and Technology for Transport, Development and Networks—Transport and Environment Laboratory “IFSTTAR-LTE”, Lyon—France.

Funding

This work is supported by University of Burundi and the French Institute of Science and Technology for Transport, Development and Networks—Transport and Environment Laboratory “IFSTTAR-LTE”, Lyon—France.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1] Roux, E. (2006) Modèle de longueur de piste au décollage-atterrissage, Avions de Transport Civil. SUPAERO-ONERA, 345.
[2] Fourer, R., Gay, D.-M. and Kernigham, B.-W. (2003) A Modelling Language for Mathematical Programming. 2nd Edition, Thomson Brooks, New York.
http://www.ampl.com
[3] Waltz, R.-A. and Plantenga, T.-D. (2008) KNITRO User’s Manual, Version 5.2. University of Colorado, Boulder, CO.
http://www.ziena.com
[4] Nahayo, F. (2012) Modèle mathématique d’optimisation non-inéaire du bruit des avions commerciaux en approche sous contraintes energétiques. Thèse de Mathématiques Appliquées de l’UCBL I.
[5] Blin, K. (2000) Stochastic Conflict Detection for Air Traffic Management. Eurocontrol Experimental Centre Publications Office, France.
https://doi.org/10.2514/6.2000-4270
[6] Boiffier, J.-L. (1999) The Dynamics of Flight, The Equations, SUPAéRO (Ecole Nationale Supérieure de l’Aéronautique et de l’Espace) et ONERA-CERT. Toulouse 25 Janvier.
[7] Roux, E. (2005) Pour une approche analytique de la dynamique du vol, Thèse, SUPAERO-ONERA.
[8] Abdallah, L., Haddou, M. and Khardi, S. (2010) Optimization of Operational Aircraft Parameters Reducing Noise Emission. Applied Mathematical Sciences, 4, 515-535.
[9] Destunder, P. (2010) Méthodes numériques pour ingénieurs. Hermes sciences publications.
[10] Chryssoverghi, I., Colestos, J. and Kokkinis, B. (2007) Classical and Relaxed Optimization Methods for Optimal Control Problems. International Mathematical Forum, 2, 1477-1498.
https://doi.org/10.12988/imf.2007.07135
[11] Fortin, A. (2008) Analyse numérique pour ingénieurs. Troisième edition, Presses internationales polytechnique.
[12] Harris, M.-M. and Mary, E. How Do We Describe Aircraft Noise? NASA TM-82712, FICAN.
http://www.fican.org
[13] Harvey Hubbard, H. (1994) Aeroacoustics of Flight Vehicles, Theory and Practices Volume 1: Noise sources and Volume 2: Noise Control. NASA Langley Research Center, Hampton, VA.
[14] James Stone, R., Groesbeck, D.E. and Zola, C.L. (1991) An Improved Prediction Method for Noise Generated by Conventional Profile Coaxial Jets. NASA TM-82712, AIAA-81-1991.
[15] Khardi, S., Nahayo, F. and Haddou, M. (2011) The Trust Region Sequential Quadratic Programming Method Applied to Two-Aircraft Acoustic Optimal Control Problem. Applied Mathematical Sciences, 5, 1953-1976.
[16] Khardi, S. (2010) Mathematical Model for Advanced CDA and Takeoff Procedures Minimizing Aircraft Environmental Impact. International Mathematical Forum, 5, 1747-1774.

Copyright © 2024 by authors and Scientific Research Publishing Inc.

Creative Commons License

This work and the related PDF file are licensed under a Creative Commons Attribution-NonCommercial 4.0 International License.