- Information
- AI Chat
Dynamic control allocation using constrained quadratic programming
Newcastle University
Preview text
Dynamic Control Allocation Using
Constrained Quadratic Programming
Ola H ̈arkeg ̊ard∗
Link ̈oping University, 581 83 Link ̈oping, Sweden
Abstract
Control allocation deals with the problem of distributing a given control demand among an available set of actuators. Most existing methods are static in the sense that the resulting control distribution depends only on the current control demand. In this paper we propose a method for dynamic control allocation, in which the resulting control distribution also depends on the distribution in the previous sampling instant. The method extends regular quadratic programming control allocation by also penalizing the actuator rates. This leads to a frequency dependent control distribution which can be designed to, e., account for different actuator bandwidths. The control allocation problem is posed as a constrained quadratic program which provides automatic redistribution of the control effort when one actuator saturates in position or in rate. When no saturations occur, the resulting control distribution coincides with the control demand fed through a linear filter.
∗Research engineer, Div. of Automatic Control, Dept. of Electrical Engineering, Link ̈oping University, 581 83 Link ̈oping, Sweden,ola@isy.liu,control.isy.liu/∼ola
Introduction
In recent years, nonlinear flight control design methods, like dynamic inversion1–3and backstepping,4,5 have gained increased attention. These methods result in control laws specifying the moments to be produced in pitch, roll, and yaw, rather than which particular control surface deflections to produce. How to transform these virtual, or generalized, control commands into actual control commands is known as the control allocation problem. Figure 1 illustrates the resulting control configuration. r v u y
x
Control Control law allocation Aircraft
Figure 1. Control configuration when control allocation is used.
With a redundant actuator suite there are several combinations of actuator positions which all produce the same virtual control, and hence give the same overall system behavior. This design freedom is often used to optimize some static performance index, like minimum control, or to prioritize among the actuators. This can be thought of as affecting the distribu- tion of control effect in magnitude among the actuators. Regardless of method(optimization based allocation,6–10daisy chain allocation,11–13direct allocation,10,14,15 etc.), the resulting mapping from the virtual control command,v(t), to true control input,u(t), can be written as a static relationship u(t)=h(v(t)) (1) A possibility that has been little explored is to also affect the distribution of thecontrol effect in the frequency domain, and use the redundancy to have different actuators operate in
Hence, there are practical cases where dynamic control allocation is desirable. In this paper, a new systematic method for dynamic control allocation is proposed. Themethod is an extension of regular quadratic programming control allocation. The key ideais to add an extra term to the optimization criterion to also penalize actuator rates. When no saturations occur, the control allocation mapping becomes a linear filter of the form
u(t)=Fu(t−T)+Gv(t)(3)
from the virtual control commandvto the actuator commands u. The frequency char- acteristics of this filter are decided by weighting matrices selected by the control designer. Thus, unlike most previous methods, no filters are to be explicitly constructed by the control designer. Two design examples are included to illustrate the potential benefits of using the proposed scheme for dynamic control allocation.
Control Allocation Problem Formulation
As stated in the introduction, an important application of control allocationis nonlinear flight control. Consider a general nonlinear dynamical model of an aircraft givenby
x ̇=f(x, δ) (4a) δ ̇=g(δ, u)(4b)
wherex = aircraft state vector, δ = actuator positions, andu = commanded actuator positions. To incorporate the actuator position and rate constraints we impose that
δmin≤δ≤δmax, ∣∣δ ̇∣∣≤δrate (5)
whereδmin andδmaxare the lower and upper position constraints, andδratespecifies the maximal individual actuator rates. Even in the case whenfandgare linear, it is nontrivial to design a control law which gives the desired closed loop dynamics while assuring that the actuator constraints are met. A common approach is therefore to split the design task into two subtasks. Neglecting the typically fast actuator dynamics, i., assumingδ=u, and viewing the actuators as pure moment generators yields the approximate model
x ̇=fM(x, M(x, u)) (6)
where M(x, u) is the mapping from the commanded actuator positions to the resulting aerodynamic moment acting on the aircraft andfMdescribes how the aerodynamic moment affects the aircraft dynamics. The control design can now be performed in two steps. First, design a control law in terms of the moment to be produced,
M(x, u)=k(r, x)(7)
that yields some desired closed loop dynamics, wherer= pilot command. Second, determine u, constrained by (5) (withδ=u), that satisfies (7). The latter step is the control allocation step. Since modern aircraft use digital control systems, it is reasonable to merge the constraints (5) into an overall time varying position constraint given by u(t)≤u(t)≤u(t)(8)
whereu∈Rmis the true control input,us∈Rmis the desired steady state control input, v∈Rk is the virtual control command,B∈Rk×m is the control effectiveness matrix, and W 1 ,W 2 ,andWvare square matrices of the proper dimensions. Bis assumed to have full row rankk.∥∥·∥∥denotes the Euclidean 2-norm defined by∥∥u∥∥=√uTu. Equation (12) should be interpreted as follows: Given Ω, the set of feasible control inputs (with respect to position and rate constraints) that minimize the virtual control error Bu(t)−v(t) (weighted byWv), pick the control input that minimizes the cost function in (12a). Hence, satisfying the virtual control demand (11) has the highest priority. When this is not possible due to the actuator constraints, (12b) corresponds to solving (11) in the least squares sense. The design matrixWv can then be used to affect the way that command limiting is performed by weighting the virtual control errors differently to prioritize certain components ofv. When there are several control inputs that give the same virtual control error (not neces- sarily zero), i., when Ω does not contain only a single point,uis made unique by minimizing the criterion in (12a). This criterion is a mix of 1) keeping the control input close to the desired steady state valueus, and 2) minimizing the change in the control input compared to the previous sampling instant. The trade-off between these two requirements isgoverned by the weighting matricesW 1 andW 2. A large diagonal entry inW 1 will make the correspond- ing actuator converge quickly to its desired position while a largeW 2 entry will prevent the actuator from moving too quickly. Note however that these weighting matrices only affect the control input ifuis not uniquely determined by (12b). The following assumption certifies that the overall control allocation problem (12) has a unique optimal solution.
Assumption 1 Assume that the weighting matrices W 1 andW 2 are symmetric and such
that W=(W 12 +W 22 ) 1 / 2 (13)
is nonsingular.
The symmetry assumption is no restriction since if, e.,W 1 is not symmetric, it can be replaced by the symmetric matrix square root (W 1 TW 1 ) 1 / 2 without affecting the solution. Equation (12) specifies which solution to the control allocation problem that issought but not how to find it. To actually solve the optimization problem, the two terms in (12a) can first be merged into one term without affecting the solution. Then, any QP solver suitable for real-time implementation6,7,9,20can be used to find the solution. Since the optimization problem (12) is to be solved at each sampling instant, no variables need to be constant. This means that the control efficiency matrixB can be updated continuously, which allows for reconfiguration after an actuator failure, and that different weighting matrices can be used for different flight cases. By including the previous control input in the optimization problem (12), the resulting control distribution will clearly be a mapping of the form
u(t)=h(v(t),u(t−T)) (14)
It is important to point out that despite the control allocator now being a dynamical system, no extra lag is introduced into the control loop since minimizing the virtual control error has top priority in (12). Let us now investigate some characteristics of the discrete time dynamicalsystem (14). In the following section we consider the nonsaturated case in whichhcan be found analytically, and investigate the issues of stability and steady state distribution.
Now, adding the linear constraint (15b), whereB has full row rank, gives the weighted, shifted pseudoinverse solution
u(t)=(I−GB)u 0 (t)+Gv(t) G=W− 1 (BW− 1 )†
from which it follows that
u(t)=(︸I−GB︷︷)W− 2 W 1 ︸ 2 E
us(t)+(︸I−GB︷︷)W− 2 W 22 ︸ F
u(t−T)+Gv(t)
which completes the proof. A more detailed proof can be found in Ref. 21.
The†symbol denotes the pseudoinverse operator defined as 22 B†=BT(BBT)− 1 for a k×mmatrixBwith full row rankk. The theorem shows that the optimal solution to the control allocation problem (15) is given by the linear filter (16). The properties of this filter will be investigatedin the two following sections.
Dynamic Properties
Let us first study the dynamic properties of the filter (16). Note that the optimization criterion in (15) does not consider future values ofu(t). It is therefore not obvious that the resulting filter (16) is stable. The poles of the filter, which can be found as the eigenvalues of the matrixF, are characterized by the following theorem.
Theorem 2 LetFbe defined as in Theorem 1 and let Assumption 1 hold. Then the eigen- values ofF,λ(F), satisfy 0 ≤λ(F)≤ 1 (20)
IfW 1 is nonsingular, the upper eigenvalue limit becomes strict, i.,
0 ≤λ(F)< 1 (21)
Proof: We wish to characterize the eigenvalues of
F=(I−GB)W− 2 W 22 =(I−W− 1 (BW− 1 )†B)W− 2 W 22 =W− 1 (I−(BW− 1 )†BW− 1 )W− 1 W 22
(22)
Let the singular value decomposition ofBW− 1 be given by
BW− 1 =UΣVT=U
[
Σr 0
]
V
rT V 0 T
=UΣrVrT
whereU andV are orthogonal matrices and Σris ak×k diagonal matrix with strictly positive diagonal entries (sinceBW− 1 has rankk). This yields
I−(BW− 1 )†BW− 1 =I−VrΣ−r 1 UTUΣrVrT=I−VrVrT=V 0 V 0 T
sinceVVT=VrVrT+V 0 V 0 T=I. Inserting this into (22) gives us
F=W− 1 V 0 V 0 TW− 1 W 22
Now use the fact 23 that the nonzero eigenvalues of a matrix productAB,λnz(AB), satisfy λnz(AB)=λnz(BA)toget
λnz(F)=λnz(V 0 TW− 1 W 22 W− 1 V 0 )
stability can be guaranteed (although asymptotic stability may hold).
- The fact that the poles lie on the positive real axis implies that the actuator responses to a step in the virtual control input are not oscillatory.
Steady State Properties
In the previous section we showed that the control allocation filter (16) is asymptotically stable under mild assumptions. Let us therefore investigate the steady state solution for a constant virtual control input.
Theorem 3 Letussatisfy Bus=v 0 (23)
wherev(t)=v 0 is the desired virtual control input. Then, ifW 1 is nonsingular, the steady state control distribution of (16) is given by
tlim→∞u(t)=us (24)
Proof: IfW 1 is nonsingular, the linear filter (16) is asymptotically stable according to Theorem 2. This means that in the limit,u(t)=u(t−T) holds. Then (15) reduces to
minu ∥∥W 1 (u−us)∥∥ 2 subject to Bu=v 0
(25)
IfussatisfiesBus=v 0 ,thenu=usis obviously one optimal solution to (25). Further, if W 1 is nonsingular,u=usis the unique optimal solution.
Sinceus can be time varying in (12), the theorem condition (23) can be fulfilled by selectingus(t)=Sv(t)whereBS = I. For example, selecting S = B† minimizes the
control input norm∥∥u∥∥at steady state. Ifusdoes not satisfy (23), the steady state control distribution will also depend onW 1. This is undesirable since it makes the role of the design parameterW 1 unclear.
Design Examples
Let us now apply the proposed method to two different design examples to see what dynamic control allocation can offer and how to select the tuning variables.
Actuator Dynamics
One application of dynamic control allocation is to account for actuator dynamics. Actuator dynamics can be an obstacle to performing control allocation since most allocation schemes – including the one proposed in this paper – assume a static relationship between the actuator commands and the resulting total control effort, see (11). Disregarding thesedynamics in cases when one or several of the actuators has a low bandwidth may deteriorate the overall system behavior, and possibly even lead to instability. A previously proposed strategy is to modify the natural actuator dynamics, using feed- back or feedforward compensation, or a combination of the two, to effectively increase the actuator bandwidth. This has proven to work well in several applications–26 However, there are situations when this is not practically feasible solution. In this section weconsider one such example and show that combining this type of compensation with dynamic control allocation can give better results. Consider the system depicted in Figure 2, with two actuators whose outputs,ua 1 andua 2 , produce a total control effort of
va=2ua 1 +ua 2 =Bua,B=
[
21
]
(26)
the sampling time isT=0 s. Both actuator outputs satisfy the position constraints and the produced control effortva(bottom) responds as expected to the command. However, in a practical situation there may be additional constraints which makes this an undesirable solution. For example, if the actuators are electrical motors, the large position error in the first actuator response may lead to an input voltage that is infeasible. To account for the difference in bandwidth in a more suitable way, let us design a dy- namic control allocator which uses both actuators to produce the low frequency part of the total control demand, but only the fast actuator for the high frequency part. This can be accomplished by selecting the tuning parameters in (12) as
us(t)=B†v(t),W 1 =diag(1,1),W 2 =diag(12,0),Wv= 1 (27)
The precompensation of the first actuator command is still necessary in order to smoothly merge the two actuator responses. The overall discrete time transfer functions from the virtual control commandvto the actuator commandsuare shown in Figure 5. Figure 4 (middle) shows the resulting step response. Initially, the fast actuator is used to produce most of the control effort, but after about 3 s the actuator commands have converged to the desired static distribution which is the same as before,u(t)=B†v(t). Thus, without affecting the static control distribution between the actuators, the transient distribution has been designed to better account for the difference in actuator bandwidth. Figure 6 shows the response whenv= 3 is commanded. This choice makes the desired steady state distributionusinfeasible. As seen from the figure, the algorithm responds by utilizing the second actuator more to compensate for the saturated first actuator.
00 0 1 1 2 2 3 3 4
1 Static allocation u
uua 1 a 2
00 0 1 1 2 2 3 3 4
1 Dynamic allocation u
uua 1 a 2
00 0 1 1 2 2 3 3 4
1
Virtual control
v
Time (s)
va
Figure 4. Simulation results for static and dynamic allocation. Dashed: commanded values. Solid: actual values.
10 −2 10 −1 100 101 102 103
−
−
−
−
0
Magnitude (dB)uu 11 (w/o precomp.) u 2
Control distribution
Frequency (rad/sec) Figure 5. Transfer functions fromvtouwith dynamic allocation.
Control reallocation
r ADMIRE FCS uAdm M(x, u) Lineariza-tion ADMIRE
Control allocation Eq. (12)
v u
x Figure 7. Overview of the closed loop system used for simulation.
vector. In the ADMIRE model the sampling time isT =0 s. The constrained least squares problem (12) is solved at each sampling instant using the weighted least squares active set solver from Ref. 20. The control input consists of the commanded deflections for the canard wings (u 1 ), the right elevons (u 2 ), the left elevons (u 3 ), and for the rudder (u 4 ). The actuator position and rate constraints in (5) are given by
δmin=
(
− 55 − 30 − 30 − 30
)T
· 180 π rad (28) δmax=
(
25 30 30 30
)T
· 180 π rad (29) δrate=
(
50 150 150 100
)T
· 180 π rad/s (30)
At trimmed flight at Mach 0, 1000 m, the control effectiveness matrix, containing the partial derivatives of the aerodynamic moment coefficients in roll (Cl), pitch (Cm), and yaw
(Cn) with respect to the control inputs, is given by
B=10− 2 ×
0 − 9. 09. 02. 7
19. 7 − 22. 4 − 22. 40
0 − 3. 33. 3 − 8. 0
rad− 1 (31)
from which it can be seen, for example, that the elevons are the most effectiveactuators for producing rolling moment while the rudder provides good yaw control, as expected. Thisis theB-matrix used in the design and analysis of the control allocation filter below. Let us now consider the requirements regarding the control distribution. At trimmed flight, it is beneficial not to deflect the canards at all to achieve low drag. We therefore select the steady state distributionusas the solution to
minus ∥∥us∥∥ subject to Bus=v us, 1 =
(32)
which yields
us(t)=Sv(t),S=
00 0
− 5. 0 − 2. 2 − 1. 7
5. 0 − 2. 21. 7
4. 10 − 11. 2
(33)
During maneuvering, corresponding to higher frequencies ofv, we seek a distribution that splits the pitch command between the canards and the elevons. Further, since the elevons have a higher rate limit than the rudder we put a higher rate penalty on the rudder.
Dynamic control allocation using constrained quadratic programming
University: Newcastle University
