# Interior Point Algorithm for Multi-UAVs Formation Autonomous Reconfiguration

Authors: Wang Jian-hong; Rana Javed Masood
Journal: Journal of Control Science and Engineering
 (2016)

---

## Abstract

Here the problem of designing multi-UAVs formation autonomous reconfiguration is considered. Combined with three kinds of cost functions, nonlinear dynamic equations, and four inequality constraints, one nonlinear multiobjective optimization problem is constructed. After applying weighted sum method and separating all equality or inequality constraints, the former nonlinear multiobjective optimization problem can be converted into a standard nonlinear single objective optimization problem. Then the interior point algorithm is applied to solve it. Further some improvements are proposed to avoid rank deficiency of some matrices. The equivalence property between multiobjective optimization and single objective optimization through weighted sum method is proved. Finally the efficiency of the proposed strategy can be confirmed by the simulation example results.


---

## Body

## 1. Introduction

Multi-UAVs formation is the basis and prerequisite for the UAV mission. When formation task or battlefield environment changes with time, the entire multi-UAVs formation needs to be adjusted. This adjustment is called formation autonomous reconfiguration. In reconfiguration process, each UAV adjusts its own position as a new geometry and plans the flight trajectory from the original location to the new terminal location. Thus the planned flight path ensures each UAV is secure while considering the nonlinear dynamics of UAV, various formations cost function, and a variety of constraints. So multi-UAVs formation autonomous reconfiguration is designed as a mathematical optimization problem in which decision variables are the control input sequences. The cost functions include three kinds of functions such as investigation UAV cost function, interference drone missile costs function, and radar jamming UAV cost function. Moreover the constraints include radar threat constraint, missile threat constraint, artillery positions constraint, and formations anticollision constraint.An overview about the auxiliary role of the multi-UAVs formation trajectory design in a cooperative investigative process is given in [1]. In [2] an intelligent multiagents system is introduced and this multiagent system is represented as the communication structure between various UAVs. In [3] descriptions of many various methods of UAV path planning are summarized. In [4] the weights which exist in many intelligent algorithms are adjusted based on the gray system theory. In [5] the methods about how to fuse the multiple sensor information are analyzed to obtain state estimates from the perspective of information fusion. In [6] the particle filter algorithm is applied to track path under the non-Gaussian condition. In [7] one interacting multiple model algorithm based on particle filter is proposed. In [8] the ant colony algorithm from the multiobjective optimization is studied. In [9] one consensus genetic algorithm is used to solve the multiobjective optimization problem. In [10] the use of game theory can transform the original multiobjective optimization problem into a Nash bargaining process. In [11] the communication delay is considered and one information compensated method based on the information filter algorithm is presented.Here the interior point algorithm is applied in multi-UAVs formation autonomous reconfiguration to design the autonomous trajectory. In the optimization model, three different cost functions are established and the dynamic equations of each UAV are regarded as equality constraints. After combining radar threats, missile threat, artillery positions threat, and formations anticollision, four kinds of inequality constraints are obtained. Based on these cost functions and equality and inequality constraints, an multiobjective optimization problem with equality and inequality constraints is constructed. The multiobjective optimization problem is solved to obtain the control input sequences which are used to plan trajectory. The main contribution of this paper is to analyze this multiobjective optimization problem by using operations algorithm. For this multiobjective optimization problem, the weighted sum strategy is applied to transform that multiobjective optimization into a single objective optimization and the feasibility of the weighted sum strategy is proved theoretically. In order to rewrite the original optimization problem into a standard form of nonlinear optimization, the kinetic equation for each equation UAV is further expanded and some forms are used to describe the set of state variables and control inputs. For the standard nonlinear optimization form, the major steps of the proposed interior point algorithm are given here. After calculating the optimal input sequences of nonlinear optimization problem by using the interior point algorithm, the idea of nonlinear receding horizon predictive control can be introduced. This idea is to choose the first term in the obtained optimal control input sequence and discard the remaining terms.


## 2. Multi-UAVs Formation Autonomous Reconfiguration

Multi-UAVs formation autonomous reconfiguration means that at initial time instant, multiple UAVs fly in a formation pattern. When the battlefield environment changes with time, the UAV swarm adjusts each UAV’s former formation pattern and chooses one new flying pattern independently. In formation reconstruction processes, each UAV’s position in the new formation mode requires to be redesigned in order to produce a new flight path. During the process of designing new flight path, the dynamics of each UAV, flight cost function, and constraints are considered. Multi-UAVs formation autonomous reconfiguration process is shown in Figure1, where four different formation modes are described. It means that the groups of multiple UAVs sequentially select their own new formation modes according to the different surrounding environment.Figure 1
Multi-UAV formation autonomous reconfiguration process.Assume there aren UAVs flying in formation, and the length of time is N. Since each UAV is decoupled with each other, this assumption means dynamic coupling phenomenon does not exist. Set state vector of the i UAV as follows:(1)xi=xi1xi2⋯xiNT.The control input sequence is that(2)ui=ui1ui2⋯uiNT.The time invariant discrete time state space equation is given as(3)xik+1=fixik,uikxi∈Ξi,uik∈Θi.In (3), fi is a nonlinear mapping; this nonlinear mapping combines the state vector xik at time k and control input vector uik to predict one new state variable xik+1 at the next time instant k+1; Ξi and Θi are the feasible state set and control input constraint set with respect to the i UAV.The main contribution of this paper is to study how to transform the multiobjective optimization problem into a standard nonlinear optimization problem, and the classical interior point algorithm is applied to analyze the equivalence between the multiobjective optimization and single objective optimization by using weighted sum strategy.When UAV formation leaps across the battlefield environment, there are many methods to define the cost function. For example, here define the position of the virtual lead aircraft at timek as(4)plk=xlkylkT.The position of i UAV vi whose role is to carry out the reconnaissance mission is denoted as(5)pik=xikyikT.The cost function corresponding to UAV vi during the whole reconnaissance mission is defined as(6)Fi1xi,ui=∑k=1Nplk-pik2+uikRi2.In (6), Ri represents a positive definite weighting matrix; the second term is used to normalize the original optimization problem and ensure that the optimal solution will not depart the true solution. Similarly two other cost functions coming from the missile and radar interference are given:(7)Fi2xi,ui=∑k=1Npa1k-pik2+uikRi2(8)Fi3xi,ui=∑k=1Npa2k-pik2+uikRi2.In (7), the cost function Fi2xi,ui is constructed by the distance between the current position pa1 and ideal position pa1 of the missile interference that interferes the i UAV. pa1 is the ideal position of the missile interference. In (8), the cost function Fi3xi,ui is constructed by the distance between the current position pi and ideal position pa2 of the radar interference that interferes the i UAV. This cost function Fi3xi,ui can achieve the maximum protection under formation flight path.In total there aren UAVs lying in the formation geometry. Combining (6), (7), and (8), we consider the following cost function:(9)minuFi1xi,uiFi2xi,uiFi3xi,uii=1nu=u1u2⋯unT.In (9), the total number of elements in the minimization operation is 3n. The minimum solution u which guarantees all 3n elements can achieve their own minimization simultaneously does not exist. Fortunately one compromise solution can be found and this compromise solution is called the efficient solution in multiobjective optimization theory. Under not any constraints, the optimization problem (9) can be solved by the weighted sum strategy. But in multi-UAVs formation autonomous reconfiguration problem, four kinds of inequality constraints are considered; the interior point algorithm can only be used to get an efficient solution.All kinds of constraints include the radar threat constraint, missile threat constraint, artillery positions threat, and formations anticollision constraint. Their respective inequalities are defined as follows in turn.Setting the position and detection radius of the radarrj as prj=xrjyrjT and Rrjk, respectively, then the radar threat constraint corresponding to the UAV vi is defined as follows:(10)gi1xik,uik=Rrjk2-pik-prjk2≤0.Setting the position of missile mj as pmj=xmjymjT, the safety distance and safety angle cosine after disturbing missile mj are dmj and csmjk, respectively. Then missile threat constraint corresponding to the UAV vi is defined as follows:(11)gi2xik,uik=csmjk-cos⁡θ2≤0.Assuming that the radiation radius of artillery positions is not influenced by external disturbance, we set this radiation radius as a constant. The position and radiation radius of artillery nj are defined as pnj=xnjynjT and Rnjk, respectively; then artillery positions threat corresponding to the UAV vi is defined as follows:(12)gi3xik,uik=Rnjk2-pik-pnjk2≤0.The minimum safe distance among UAV formation is denoted by dmin, and then formations anticollision constraint corresponding to the UAV vi is defined as follows:(13)gi4xik,uik=dmin2-pik-pjk2≤0.Unifying the above multiobjective function (9), nonlinear dynamical equations (3), and four inequality constraints (10)–(13), multi-UAVs formation autonomous reconfiguration problem can be formulated as a nonlinear multiobjective optimization problem:(14)minuFi1xi,uiFi2xi,uiFi3xi,uii=1nsubjecttoxik+1=fixik,uikxi∈Ξi,uik∈Θii=1,2,…,n,k=1,2,…,Ngi1xik,uik≤0,gi2xik,uik≤0,gi3xik,uik≤0,gi4xik,uik≤0.Combining all the inequality constraints in (14) and vectoring them, we get(15)gix,u=gi1xik,uikgi2xik,uikgi3xik,uikgi4xik,uik,gx,u=g1x,u⋮gnx,u≤0,where the state variable and control input vector are given, respectively, as follows:(16)x=x1x2⋯xnT,u=u1u2⋯unT.Similarly combining all the nonlinear dynamics equations of each UAV, we get(17)xk+1=x1k+1x2k+1⋮xnk+1=f1x1k,u1kf2x2k,u2k⋮fnxnk,unk=fxk,uk.Using the vector form (15), (16), and (17), nonlinear multiobjective optimization problem can be simplified as(18)minuFi1xi,uiFi2xi,uiFi3xi,uii=1nsubjecttoxk+1=fxk,ukxi∈Ξi,uik∈Θii=1,2,…,n,k=1,2,…,Ngx,u≤0.Before solving (18) by the interior point algorithm, we transform the multiobjective optimization problem (18) into a nonlinear optimization problem with equality and inequality constraints.


## 3. Standardized Model

In (18), the cost function is a multiobjective vector:(19)minuFi1xi,uiFi2xi,uiFi3xi,uii=1n.Applying the weighting sum strategy and transforming to one single objective optimization problem,(20)minuFu=∑i=1nλ1iFi1xi,ui+λ2iFi2xi,ui+λ3iFi3xi,ui.In (20), vector λ1i,λ2i,λ3ii=1n indicates positive weighted scalar values and the following conditions are satisfied:(21)∑i=1nλ1i+λ2i+λ3i=1λ1i>0,λ2i>0,λ3i>0.The nonlinear dynamic equations of UAV (18) can be rewritten as follows:(22)x0-x0=0xk+1-fxk,uk=0,k=1,2⋯N-1.In (22), the first equation represents the initial state. For the set of feasible state Ξi and control input constraint set Θi, define the following constraints:(23)Ξi=xi∃ui∈Θi,withGisxi,ui=0,Hisxi,ui≤0Θi=ui∃xi∈Ξi,withGisxi,ui=0,Hisxi,ui≤0.The two above functions Gisxi,ui and Hisxi,ui do not depend on two variables xi and ui simultaneously.Combining all the equations and inequalities together and using (20), we obtain a standard nonlinear optimization problem:(24)minuFu=∑i=1nλ1iFi1xi,ui+λ2iFi2xi,ui+λ3iFi3xi,uiwithrespecttou=u1u2⋯unTsubjecttoGu=x0-x0xk+1-fxk,ukk=1N-1Gisxi,uii=1n=0,Hu=gx,uHisxi,uii=1n≤0.In order to solve the standard nonlinear optimization problem (20), the interior point algorithm is proposed.


## 4. Interior Point Algorithm

The purpose of the interior point algorithm is to generate an iterative sequenceuj; here the superscript symbol uj is different with control input uk at time instant k. This generated sequence will be included in the control input set. In the iteration process of each generating sequence, each element of the inequality constraint Hu is considered. After introducing the slack variable s, those inequality constraints can be converted into equality constraints [12]. The standard nonlinear optimization problem (24) can be rewritten as(25)minuFusubjecttoGu=0,Hu+s=0,s≥0.In (25), the slack variable s is chosen as a vector with appropriate dimension, and each of its elements is nonnegative. Then we construct one Lagrangian function corresponding to (25) as follows:(26)Lu,s,v,w=Fu-vTGu-wTHu+s.According to the necessary condition from the optimality theory, the equation Hu+s=0 holds at the minimum value.Applying the optimality Karush-Kuhn-Tucker sufficient and necessary condition [13], we obtain(27)∇Fu-∇GTuv-∇HTuw=0Sw-μe=0Gu=0,Hu+s=0.From (27), we have that(28)s≥0,w≥0,μ=0,e=11⋯1T.In (27), define the following two matrices:(29)S=diag⁡s,W=diag⁡w.In (27), the introduction of perturbation parameter μ is to guarantee that if we choose μ>0, the iteration sequence will be far away from the boundary of the control input set. The choice of the perturbation parameters μ satisfies that(30)μj⟶0asj⟶∞.In order to improve the performance of the optimal conditions (27), we use the interior point algorithm to solve (25) and introduce barrier function to eliminate the nonnegative condition s≥0:(31)minuFu-μ∑i=1nslog⁡sisubjecttoGu=0,Hu+s=0.The barrier function -μ∑i=1nslog⁡si is added here to prevent each element of slack variable s to approach 0 closely. Also applying the generalized optimality KKT sufficient and necessary condition [14] in (31) with one barrier function, we obtain(32)∇Fu-∇GTuv-∇HTuw=0μS-1e+w=0Gu=0,Hu+s=0.Using Newton increment steps in (32), we obtain the following system equation:(33)∇uu2Lu,s,v,w0-∇GTu-∇HTu0W0S∇Gu000∇Hu-I00ΔuΔsΔvΔw=-∇Fu-∇GTuv-∇HTuwSw-μeGuHu+s.After calculating the increment,(34)ΔuΔsΔvΔw.The new recursive iteration values are calculated as(35)uj+1=uj+αsmaxΔu,sj+1=sj+αsmaxΔsvj+1=vj+αwmaxΔv,wj+1=wj+αwmaxΔw,where two steps αsmax and αwmax are chosen as(36)αsmax=max⁡α∈0,1suchthats+αΔs≥1-τsαwmax=max⁡w∈0,1suchthatw+αΔw≥1-τw,where τ∈0,1, and take τ=0.995. To ensure the second block matrix in (33) has full row rank and this matrix is not singular, we rewrite (33) as follows:(37)∇uu2Lu,s,v,w0-∇GTu-∇HTu0∑0I-∇Gu000-∇HuI00ΔuΔsΔvΔw=-∇Fu-∇GTuv-∇HTuww-μS-1e-Gu-Hu-s,where ∑=S-1W, as in the iterative process, the normed function is a decreasing function. It means that the incremental vector ΔuΔsΔvΔw is a descent direction, so it is desirable to choose matrix(38)∇uu2Lu,s,v,w00∑as a positive definite matrix which is on the null space with respect to matrix (39)-∇Gu0-∇HuI.From the construction of ∑, it means ∑ is a positive definite matrix. But Hessian matrix ∇uu2Lu,s,v,w may be negative definite. To compensate for this defection, we use ∇uu2Lu,s,v,w+δI to replace Hessian matrix ∇uu2Lu,s,v,w. Scalar δ>0 is chosen sufficiently large to ensure the positive definiteness of the Hessian matrix. Additionally in the interior point algorithm, the rank deficiency of the gradient matrix ∇Gu is considered. So we make the following modifications to the primal-dual matrices:(40)∇uu2Lu,s,v,w+δI0-∇GTu-∇HTu0∑0I-∇Gu0γI0-∇HuI00ΔuΔsΔvΔw=-∇Fu-∇GTuv-∇HTuww-μS-1e-Gu-Hu-s.One normalized parameter γ>0 is added in the first matrix (40). As the iterative expressions (35) do not terminate in a limited period of time, one error criterion function can be applied to determine when to stop the entire iterative algorithm. This error criterion function may be that(41)Eu,s,v,w=max⁡∇Fu-∇GTuv-∇HTuw,Sw-μe,Gu,Hu+s.With this error criterion function as a judge to terminate the iterative algorithm, we summarize the basic steps of the interior point algorithm.Step 1.
Assume a pair of initial value is given asu0,s0, and let j=0.Step 2.
Calculate Lagrange multipliersv0 and w0 and define the parameters(42)μ0>0,σ,τ∈0,1.Step 3.
Verify whether the error criterion function holds:(43)Euj,sj,vj,wj≥ε.ε is a very small positive number; if above equality does not hold, then stop the iterative algorithm. Then the optimization variables uj,sj,vj,wj at this time can be regarded as the optimal solution to a nonlinear optimization problem.Step 4.
Solve (40) and calculate the search direction ΔuΔsΔvΔw.Step 5.
Use (36) to determine two steps αsmax and αwmax.Step 6.
Use the iterative expression (35) to solve the new iterative value.Step 7.
Letμj+1=μj,  k=k+1.Step 8.
Return to Step3.In above steps, we make some modifications to compensate the drawback of the usual interior point algorithm. These modifications can guarantee the optimal solution of the original nonlinear optimization problem [15].


## 5. Equivalence between the Two Optimization Problems

When solving the multiobjective optimization problem, firstly we apply the weighted sum strategy to convert it into one single objective optimization problem. Under the conditions of the positive weighted scalars (21), the equivalence between these two optimization problems holds. As there are 3n elements in (19), rewrite equality (19) as follows:(44)minuF1uF2u⋯F3nu.Similarly rewrite equality (20) as follows:(45)minuFu=∑i=13nλiFiu.Using (44) and (45), the conditions of equivalence between (44) and (45) can be summarized as the following proposition.Proposition 1.
Letλi>0,i=1,2,…,3n, and ∑i=13nλi=1. If u^ is an optimal solution of (45), then u^ is an effective solution of (44).Proof.
Letu^ be an optimal solution of (45), and it satisfies Fu′≤Fu^. Applying the use of positive weights Fiu′≤Fiu^,i=1,2⋯3n, we obtain the following equation:(46)∑i=13nλiFiu′<∑i=13nλiFiu^.Choose a sufficiently large positive number M such that the negative result about optimal solution u^ holds. This positive number M is chosen as(47)M=3n-1maxi,j⁡λjλi.Next we use the contradiction method to prove this proposition. Let u^ be not an effective solution of (44), then there exists one i∈1,2⋯3n and u such that(48)Fiu≤Fiu^Fiu^-Fiu>MFiu-Fiu^.For all j∈1,2⋯3n, we have(49)Fju≥Fju^.Continuing to formulate, we get(50)Fiu^-Fiu>3n-1λiλjFju-Fju^,j≠i.Multiplying both sides of (50) by λi/(3n-1) and taking sum operations, we have(51)λiFiu^-Fiu>∑j≠iλjFju-Fju^⟺λiFiu^-λiFiu>∑j≠iλjFju-∑j≠iλjFju^⟺λiFiu^+∑j≠iλjFju^>λiFiu+∑j≠iλjFju⟺∑i=13nλiFiu<∑i=13nλiFiu^.In (51), the last inequality contradicts the assumption that u^ is an optimal solution of (45). So it means that u^ is an effective solution of (44).


## 6. Numerical Example

To verify the interior point algorithm in multi-UAVs formation autonomous reconfiguration, the formation includes three UAVs: one radar jamming UAV, one missile jamming UAV, and one investigation UAV. The initial positions of three UAVs are all located at starting point coordinates (0,0) and the terminate positions are concentrated on coordinates (700,700). The vector of the flight maximum speed, minimum speed, and the speed deviation is that(52)8015±5.Surrounding battlefield environment contains radar threat, missile threat, and antiaircraft positions threat. The deployment coordinates of the radar threat are (300,300), and the deployment coordinates of the missile threat are (250,200). The region of the artillery positions threat is a rectangle range with a height of 300 and a width of 300. This rectangle range belongs to the no-fly zone.Applying the interior point algorithm in multi-UAVs formation autonomous reconfiguration, the weighted factor for each UAV is(53)Qi=diag⁡10,10.The weighted matrix in the cost function is (54)Ri=diag⁡1,10.The discrete time sampling period is taken as(55)ΔT=0.05s.The length of time is N=500s, the number of UAVs is n=3, and the initial value during the interior point algorithm is chosen as follows:(56)u0,s0=0.010.010.010.01T.Perturbation parameter is μ0=0.05, a positive number is ε=0.01, the normalized parameter is γ=0.5, and the scalar is δ=1.5.The simulation trajectory is shown in Figure2. The coordinates of threat in surrounding battlefield environment are set as (296,346) and (229,173); the effect radius of the threat is 100 m. Assume the threat of the threat zero is infinite. In order to increase the flight performance, the density of sampling points in the vicinity of the threat can be higher than the area of no threat. In Figure 2, ○ denotes the global optimal navigation point and ※ denotes the position of the optimal navigation point. Two annular regions represent a range of threats.Figure 2
Multi-UAV formation autonomous reconfiguration trajectory.From Figure2, three UAVs formations experience a total of three formation autonomous reconfiguration processes. The first autonomous reconfiguration occurs in coordinates (200,100); the second occurs in coordinates (280,280). After these two formation autonomous reconfigurations, three UAVs can bypass these three threats and fly from initial position to terminate position.Figure3 shows the iterative convergence curves of cost functions with respect to three UAVs. Here each cost function is written as formula (20) with the weighted scalar value being λ=1/9. From Figure 3, with the interior point algorithm runs, after we substitute the optimal control input solution u into the corresponding weighted sum cost function, the cost function will approach to zero closely with increasing number of iteration steps.Figure 3
Convergence curves of each UAV cost function.
(a)
UAV1 (b)
UAV2 (c)
UAV3


## 7. Conclusion

After establishing one nonlinear multiobjective optimization problem in multi-UAVs formation autonomous reconfiguration, we use weighted sum strategy and combine all the equations or inequalities to derive a standard single objective nonlinear optimization problem. Furthermore the interior point algorithm is proposed to solve the optimization problem and some improvements are made to ensure the optimal solution of the original nonlinear optimization problem. But in this paper, the asymptotic property and sensitivity analysis are not studied, so these two are the next research subjects.



---

*Source: 9095372-2016-08-09.xml*