LMI–based robust controller design approach in aircraft multidisciplinary design optimization problem

Abstract

This article proposes a linear matrix inequality–based robust controller design approach to implement the synchronous design of aircraft control discipline and other disciplines, in which the variation in design parameters is treated as equivalent perturbations. Considering the complicated mapping relationships between the coefficient arrays of aircraft motion model and the aircraft design parameters, the robust controller designed is directly based on the variation in these coefficient arrays so conservative that the multidisciplinary design optimization problem would be too difficult to solve, or even if there is a solution, the robustness of design result is generally poor. Therefore, this article derives the uncertainty model of disciplinary design parameters based on response surface approximation, converts the design problem of the robust controller into a problem of solving a standard linear matrix inequality, and theoretically gives a less conservative design method of the robust controller which is based on the variation in design parameters. Furthermore, the concurrent subspace approach is applied to the multidisciplinary system with this kind of robust controller in the design loop. A multidisciplinary design optimization of a tailless aircraft as example is shown that control discipline can be synchronous optimal design with other discipline, especially this method will greatly reduce the calculated amount of multidisciplinary design optimization and make multidisciplinary design optimization results more robustness of flight performance.

Introduction

Control discipline is one of the most important disciplines of multidisciplinary design optimization (MDO). Since control discipline design is always based on certain objects and its design parameters, such as aircraft, which results in the control discipline could be started to design only after completed the design of other disciplines such as aerodynamics, structure, propulsion. Vehicle MDO with control discipline, therefore, is often the only single discipline optimization design under the multidisciplinary design environment, and it is difficult to synchronously optimize the control discipline and other disciplines.^1–3 Although Dr Zheng Anbo et al. ⁴ put control disciplinary into the MDO design loop, the control disciplinary design and trajectory simulation all put the top system level and the control parameters fixed after complete design of the other disciplinaries, thus control disciplinary does not still synchronously optimize with other disciplinaries. Recently, Fan Wen-feng et al. ⁵ only take index of stability and maneuverability as the quantitative evaluation of MDO by the simple analysis of control disciplinary. With the continuous advance of MDO technology, more and more attentions are paid to optimize the robustness of multidisciplinary integrated design with the exception of the control discipline in order to improve this shortage. Researchers mostly focus on seeking the optimization of a design space with “flat peak.” Through improving MDO approach, it is expected to find a design space in which response values of the system performance index can be near the peak value and constitute “a flat response surface,”^6,7 when the design parameters vary in a super sphere whose radius is $\overrightarrow{r}$ . Although the explorations of these approaches have significant theoretical and engineering meaning, it is a passive approach because of the lack of consideration about the role of control discipline, and it may even be insoluble for a specific aircraft MDO problem.

In order to improve the robustness of the MDO design, control discipline must be put into the design loop to carry out the synchronous design between control discipline and various disciplines of controlled objects. This is an active approach which essentially levels off the performance response surface generated in the original design space. Actually, there are a set of complete theoretical method in robust control system design for solving problems of modeling errors and runtime ambient interference. ⁸ In the perspective of multidisciplinary design and synthesis, the variation in design parameters excepting the control discipline can be treated as “modeling errors or uncertainty of model parameters,” whereby a set of robust control optimization propositions can be constructed, in which multidisciplinary design parameters and state feedback control parameters are regarded as optimization variables, functions, and performance constraints of aircraft. With this approach, synchronous design and optimization of the control disciplinary with other disciplines become feasible.

Since the aircraft geometric parameters, engine design parameters, and so on are time-invariant parameters, the performance index, stability, and control in the conceptual phase are always calculated by the simplified theoretical formulas, which convert the related design variables into time-invariant variables. So, the MDO problem can be described as constrained optimization problem of nonlinear static parameters. However, if aircraft control discipline is brought into the design loop, control structures and parameters would be closely related to time-variant state variables of aircraft space motion. For example, the control of pitch channel is always related to time-variant parameters such as pitching angle, pitching angular speed, and height. In addition, for a complex multidisciplinary system, the optimal design only based on the design parameters may be inexistent or just a sharp “trough (peak) point,” and the robustness of the design is very poor.^2,9 So, it is necessary to introduce control discipline to widen the envelopes of optimized results. From the viewpoint of the characteristics of aircraft design, this article intends to describe the parameter variations in linear model with the direct utilization of the perturbations of various disciplinary design parameters, to propose a less conservative design method for robust controller. This can unify the optimization of controller in control discipline with the optimization of design parameters in other disciplines, so that the synchronous design optimization of control discipline and other disciplines is achieved.

The remainder of this article is organized as follows. In section “Description of MDO problem with control discipline in loop,” a basic description of MDO problem is presented. Then some basic concepts such as the design parameters, coefficient array, and state variable are described. In section “Design method of robust controller,” the reason why the new uncertain description method should be explored is discussed. And then this article proposes parameter uncertain description based on response surface approximation (RSM) and corresponding design method of linear matrix inequality (LMI) robust controller. In section “MDO approach with robust controller in design loop,” the LMI robust controller is applied to MDO design problem, and an MDO design of a tailless aircraft as example is illustrated to validate this method. Finally, conclusions are made in section “Conclusion.”

Description of MDO problem with control discipline in loop

This article considers the multidisciplinary analysis model of aircraft with control discipline in loop transformable to the following form

$$\{\begin{array}{c}\hfill \stackrel{·}{x}(t)=F(x(t),u(t),t,p_1,p_2,\dots ,p_i,\dots )\hfill \\ \multicolumn{1}{c}{\begin{array}{c}\hfill y(t)=G(x(t),u(t),t,p_1,p_2,\dots ,p_i,\dots )\hfill \\ \multicolumn{1}{c}{p_i=S_i({\mathrm{ds}}_1,{\mathrm{ds}}_2,\dots ,{\mathrm{ds}}_m)i=1,2,\dots ,n}\\ \multicolumn{1}{c}{u(t)=S_c(x(t),t,p_1,p_2,\dots ,p_i,\dots )}\\ \multicolumn{1}{c}{\left\{d\right\}=\{{\mathrm{ds}}_1,{\mathrm{ds}}_2,\dots ,{\mathrm{ds}}_m,u(t)\}}\end{array}}\end{array}$$

(1)

where $\{F,G\}$ represents the aircraft motion model; $\{S_i\}$ represents the analysis model of ith sub-discipline; $\{S_c\}$ represents the analysis model of control discipline; x, u, and y are the state, input, and output variables, respectively; $p_i$ denotes the coefficient array which is the design result of ith sub-discipline; ${\mathrm{ds}}_j$ is the static design parameters of jth component or subsystem, later taking the compact forms $\{\mathrm{ds}\}$ to represent the set of static design parameters $\{{\mathrm{ds}}_1,{\mathrm{ds}}_2,\dots ,{\mathrm{ds}}_m\}$ ; $\{d\}$ denotes all design parameters which contain static design variables $\{\mathrm{ds}\}$ and time-variable design variables $\{u(t)\}$ ; here, there are n + 1 sub-discipline in total such as aerodynamics, structure, and propulsion, while there are m component or subsystem, for example, geometry, structure, engine, and control system.

Suppose the constraint condition and objective function are represented by the following forms:

• The equality constraint: $g_1(x(t_f),t_f,d,p_i)=0$

• The inequality constraint: $g_2(x(t_f),t_f,d,p_i)\le 0$

• The objective function of static parameters: $J(d)=\Psi (d,p_0,p_s)$

• The objective function of control discipline: ${J(u)|}_{\mathrm{sy}=F({\mathrm{sd}}_0)}=\Vert G_{\mathrm{zw}}{\Vert }_2\mathrm{or}\Vert G_{\mathrm{zw}}{\Vert }_{\infty }$

The aircraft MDO is essentially to determine a set of multidisciplinary design schemes $\{d\}$ under constraints of the given equations and inequalities. Then the general description of aircraft MDO problem with control discipline (MDO-C for short) is given: MDO-C is to confirm an optimized design scheme $\left\{d^{*}\right\}=\{d_{s1}^{*},d_{s2}^{*},\dots ,d_{s(n-1)}^{*},d_{\mathrm{ss}}^{*},u^{*}(t)\}$ , which can make the states in the system meet the various limiting conditions, and make the performance indexes $J\left(d\right)$ and $J\left(u\right)$ with extremum.

System (1) includes not only the static design parameters $\{\mathrm{ds}\}$ which describe the system features but also the dynamic design parameters $\{u^{*}(t)\}$ which reflect the effects of system control, so that the system can get optimal value under the prescribed performance index (or objective function). That is to say, MDO is to find a set of static design parameters $\{{\mathrm{ds}}_1,{\mathrm{ds}}_2,\dots ,{\mathrm{ds}}_m\}$ , which contain aircraft geometry parameters (including wing, fuselage, control surface, etc.), quality parameters, performance parameters, engine parameters, and so on; to convert dynamic system from the initial state to a specific terminal state, or make some states of the dynamic system vary following the required laws, and the required performance index is guaranteed to achieve the optimal value, under the effect of corresponding controller $\{u^{*}(t)\}$ . So the aircraft multidisciplinary model could be divided into two categories, one is the static parameter model (SM) in which design variables are time-invariant, such as discipline model of aerodynamics, structural, and propulsion; the other is dynamic parameter model in which design variables are time-variant, or related to the time-variant state variables such as aircraft motion model and control system model.

In the conceptual design phase, aiming to reduce the MDO complexity, the model described by $\{F,G\}$ can be transformed into a linear model, namely, in the premise of small perturbations, aircraft motion model can be decomposed into longitudinal and lateral motion models through the linearization at each characteristic point. Taking the integration of quadratic functions of state variables and control variables as the functions of performance index, the optimization problem of controller can be transformed into the $H_2$ optimal control problem, and the optimal solution $\{u^{*}(t)\}$ can be expressed with a unified analytic expression. Thus, a simple linear feedback control law would be obtained in the form of $u=\mathrm{Kx}$ .^10,11 The solving of $\{u^{*}(t)\}$ essentially means the confirmation of feedback matrix K. When the terminal time $t_f$ is infinitely large, the feedback matrix K becomes a constant matrix in a linear time-invariant system. By this way, both the optimization of static parameters and dynamic parameters can be unified into searching the optimal solution of static parameters.

Design method of robust controller

Robust control is a perfect theoretical system which can solve the problems of modeling errors and disturbances of operating environment. However, in the perspective of conceptual design, the variations in static design parameters $\{\mathrm{ds}\}$ within specified ranges will inevitably lead to the variations in various coefficient arrays $\{p_i\}$ and result in the parameter uncertainty in linear aircraft model. In fact, these variations can be virtually treated as a kind of “modeling errors.” Robust controller can handle various uncertainties in the precise frame of a method, and the designed controller is able to cover the certain variable intervals of aircraft design schemes $\{\mathrm{ds}\}$ . Consequently, it provides possibilities for concurrent design and optimization of control discipline and other disciplines. In general, the boundary of various disciplinary design schemes $\{\mathrm{ds}\}$ is much larger than the perturbation boundary which is allowed by the controlled object in a specific design scheme. Without the consideration of the physical nature of uncertainties, system uncertainty area is determined directly by ranges of the coefficient arrays $\{p_i\}$ in linear model, thus the designed robust controller always tends to be very conservative and sometimes even insoluble.^12,13 Therefore, an approach of uncertainty description with less conservative should be sought during the aircraft MDO research which directly employs variations in various design parameters to map the perturbations of coefficient arrays in linear model. On the basis of these, MDO approach with control discipline in design loop can be researched.

Parameter uncertainty description based on response surface approximation

Supposing that $a_1,a_2\in \{p_i\}$ are two perturbation parameters in linear motion model of aircraft, which correspond to the coefficient arrays $\{p_i\}$ of sub-disciplinary design results which depend on static design parameters $\{\mathrm{ds}\}$ . In general, the following procedures can determine the uncertainty of $a_1$ and $a_2$ : first, define the value of design variable $\{\mathrm{ds}\}$ and construct a design space by experimental design methods; then, confirm the number of simulation calculation experiments and the number of variable combinations in each experiment; and finally, calculate the response of $a_1$ and $a_2$ at each test point and get the minimum values $a_{1min},a_{2min}$ and the maximum values $a_{1max},a_{2max}$ based on calculations, therefore the uncertainty of $a_1,a_2$ can be determined. When the variation in design parameters $\{\mathrm{ds}\}$ becomes larger, this method will show obvious disadvantages. As indicated by the curves in Figure 1, because $a_1$ and $a_2$ are not independent of each other, when the value of $a_1$ varies at the interval of $[a_{1min},a_{1max}]$ , the value of $a_2$ will only be a specific value at the interval of $[a_{2min},a_{2max}]$ (or vary at a smaller specific interval) instead of varying randomly at the interval of $[a_{2min},a_{2max}]$ . In the proposed method, various uncertainties have no strict physical meaning, and the described uncertainty ranges are always much larger than the possible uncertain ranges of actual object, so the designed control system is very conservative and even insoluble, and it is difficult to balance between robustness and performance index.

Figure 1.

Schematic diagram of parameter uncertainty range.

Response surface method (RSM) has been widely used in MDO. It can conveniently establish the mathematical model of the problems influenced by multiple factors and approximate the complicated functions in relatively simple forms.^14–16

Since the various design parameters generally vary in small ranges in a particular aircraft MDO scheme, it is able to divide $\{\mathrm{ds}\}$ into multiple intervals, and then construct multiple central points of $\{\mathrm{ds}\}$ and the corresponding perturbation boundaries. With the application of RSM at each interval, we can establish the response surface model of $\mathrm{ds}\to p_i$ and describe system uncertainty while regarding $\{\mathrm{ds}\}$ as perturbation parameters. Thus, the uncertainty description will be less conservative and can meet the requirements of MDO. In general, RSM model contains the major procedures of central composition design, multi-variable regression, response surface generation, regression equation validation, nonlinear optimization, verification of model’s forecasted values, and so on. ¹⁴ Now, it is supposed that the response surface model of $\mathrm{ds}\to p_i$ has been obtained as follows

$$a_i={\beta }_0+\sum _{i=1}^n{\beta }_i{\mathrm{ds}}_i+\sum _{i=1}^n{\beta }_{\mathrm{ii}}{\mathrm{ds}}_i^2+\sum _{i=1}^{n-1}\sum _{j=i+1}^n{\beta }_{\mathrm{ij}}{\mathrm{ds}}_i{\mathrm{ds}}_j$$

(2)

where $a_i\in \{p_j\}$ indicates ith coefficient array of jth sub-discipline; $\mathrm{Ds}=[{\mathrm{ds}}_1,{\mathrm{ds}}_2,\dots ,{\mathrm{ds}}_n{]}^T$ indicates various static design parameters; the coefficients ${\beta }_0,{\beta }_i,{\beta }_{\mathrm{ii}}$ , and ${\beta }_{\mathrm{ij}}$ can be determined by the central composition design experiment and the regression analysis approach.

Converting formula (2) into matrix form

$$a_i={\beta }_0+b^T·\mathrm{Ds}+{\mathrm{Ds}}^T·\mathrm{B}·\mathrm{Ds}$$

(3)

where

$$\begin{array}{cc}\multicolumn{1}{c}{}& \mathrm{Ds}={\left[\begin{array}{cccc}\hfill d_1\hfill & \hfill d_2\hfill & \hfill \cdots \hfill & \hfill d_n\hfill \end{array}\right]}^T\hfill \\ \multicolumn{1}{c}{}& b={\left[\begin{array}{cccc}\hfill {\beta }_1\hfill & \hfill {\beta }_2\hfill & \hfill \cdots \hfill & \hfill {\beta }_n\hfill \end{array}\right]}^T\hfill \\ \multicolumn{1}{c}{}& \mathrm{B}=\left[\begin{array}{cccc}\hfill {\beta }_{11}\hfill & \hfill {\beta }_{12}/2\hfill & \hfill \cdots \hfill & \hfill {\beta }_{1n}/2\hfill \\ \multicolumn{1}{c}{{\beta }_{12}/2}& \hfill {\beta }_{12}\hfill & \hfill \hfill & \hfill {\beta }_{2n}/2\hfill \\ \multicolumn{1}{c}{⋮}& \hfill \hfill & \hfill \ddots \hfill & \hfill ⋮\hfill \\ \multicolumn{1}{c}{{\beta }_{1n}/2}& \hfill {\beta }_{2n}/2\hfill & \hfill \cdots \hfill & \hfill {\beta }_{\mathrm{nn}}\hfill \end{array}\right]\hfill \end{array}$$

Assuming that central point of design parameter is ${\mathrm{DS}}_0$ , and the perturbation parameters $\mathrm{E}\in [{\mathrm{Ds}}_{0min},{\mathrm{Ds}}_{0max}]$ , then

$$\mathrm{Ds}={\mathrm{Ds}}_0+\mathrm{E}$$

(4)

where $\mathrm{E}$ indicates perturbation vector of parameters: $\mathrm{E}=[\begin{array}{cccc}\hfill {\epsilon }_1\hfill & \hfill {\epsilon }_2\hfill & \hfill \dots \hfill & \hfill {\epsilon }_n\hfill \end{array}{]}^T$ .

Substitute formula (4) into formula (3) to get

$$a_i={\beta }_0+b^T({\mathrm{Ds}}_0+\mathrm{E})+({\mathrm{Ds}}_0+\mathrm{E}{)}^T\mathrm{B}({\mathrm{Ds}}_0+\mathrm{E})$$

(5)

Normalizing $\mathrm{E}$ at the interval of $[{\mathrm{Ds}}_{0min},{\mathrm{Ds}}_{0max}]$ , the perturbation vector $\Delta \in [-1,+1]$ can be obtained, and

$$\mathrm{E}=(max({\mathrm{Ds}}_{0max})-min({\mathrm{Ds}}_{0min}))\Delta =k\Delta$$

(6)

Let $I=[\begin{array}{cccc}\hfill 1\hfill & \hfill 1\hfill & \hfill \dots \hfill & \hfill 1\hfill \end{array}{]}^T$ , then $\Delta \le I$ .

Simplifying and deducing formula (5) result in

$$\begin{array}{cc}\multicolumn{1}{c}{a_i}& =a_{i0}+b^T\mathrm{E}+2{\mathrm{Ds}}_0^T\mathrm{B}\mathrm{E}+{\mathrm{E}}^T\mathrm{B}\mathrm{E}\hfill \\ \multicolumn{1}{c}{}& =a_{i0}+k(b^T+2{\mathrm{Ds}}_0^T\mathrm{B}+k{\Delta }^T\mathrm{B})\Delta \hfill \\ \multicolumn{1}{c}{}& \le a_{i0}+k(b^T+2{\mathrm{Ds}}_0^T\mathrm{B}+{\mathrm{kI}}^T\mathrm{B})\Delta \hfill \\ \multicolumn{1}{c}{}& =a_{i0}+{\Lambda }_i\Delta \hfill \end{array}$$

(7)

where $a_{i0}$ indicates the value of the ith perturbation parameter in the linear model, which is at the central points of design parameters ${\mathrm{Ds}}_0$ , while ${\Lambda }_i$ indicates the uncertainty weight, and ${\Lambda }_i=k(b^T+2{\mathrm{Ds}}_0^T\mathrm{B}+{\mathrm{kI}}^T\mathrm{B})$ is a row vector.

The uncertainty expressions of the other coefficient in linear model are solved in the similar way as above, so the uncertainty models based on the variation in design parameters can be obtained which are expressed with the normalized perturbation vector $\Delta$ of design parameters. $\Lambda \Delta$ reflects the variations in actual physical parameters so that this model has the characteristic of less conservative. And $\Delta$ has the same dimension as that of design parameters, excluding redundant variables. Therefore, compared with the traditional design approach of robust controller, the order of uncertainty matrix in this approach is much lower

$$p_i={\mathrm{Ds}}_0+\Lambda \Delta ={\left[\begin{array}{c}\hfill a_1\hfill \\ \multicolumn{1}{c}{a_2}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{a_m}\end{array}\right]}_0+\left[\begin{array}{c}\hfill {\Lambda }_1\hfill \\ \multicolumn{1}{c}{{\Lambda }_2}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{{\Lambda }_m}\end{array}\right]\Delta$$

(8)

Using formula (8), it is possible to describe the uncertainty of the gain blocks in the longitudinal or lateral linear model of aircraft, and the aircraft linear fraction transformation (LFT) model can be deduced.

Let $\begin{array}{cc}\hfill {\overline{\Lambda }}_i\hfill & \hfill =b^T+2{\mathrm{Ds}}_0^T\mathrm{B}+{\mathrm{kI}}^T\mathrm{B}\hfill \\ \multicolumn{1}{c}{}& \hfill =[\begin{array}{cccc}\hfill {\lambda }_1\hfill & \hfill {\lambda }_2\hfill & \hfill \cdots \hfill & \hfill {\lambda }_n\hfill \end{array}]\hfill \end{array}$ , formula (7) can be transformed into the following form

$$\begin{array}{cc}\multicolumn{1}{c}{a_i}& \le a_{i0}+{\overline{\Lambda }}_{i}k\Delta \hfill \\ \multicolumn{1}{c}{}& =a_{i0}+{\overline{\Lambda }}_i\mathrm{E}\hfill \end{array}$$

(9)

According to formula (9), the conclusion can be shown: in the aircraft multidisciplinary design environment, if the mapping relationship between various static design parameters $\mathrm{Ds}$ and various disciplinary coefficient arrays $\{p_i\}$ can be approximated by the quadratic response surface model as described by formula (3), all the perturbation parameters in linear aircraft motion model can be linearly expressed by the perturbation vector $\mathrm{E}$ of design parameters. If this property is applied to the uncertainty description of the state equation interval model, ¹⁷ the order of the matrix that describes uncertainty will decline substantially, which result in the great reduction in the calculated amount.

Design method of robust controller based on design parameters

Supposing that the state-space model $\left(\begin{array}{ccc}\hfill A\hfill & \hfill B\hfill & \hfill \begin{array}{cc}\hfill C\hfill & \hfill D\hfill \end{array}\hfill \end{array}\right)$ can be used to describe the linear aircraft motion model on aircraft motion model $\{F,G\}$ , then $\left(\begin{array}{ccc}\hfill A\hfill & \hfill B\hfill & \hfill \begin{array}{cc}\hfill C\hfill & \hfill D\hfill \end{array}\hfill \end{array}\right)$ can be expressed as a function space of the various design parameters $\{\mathrm{ds}\}$ , and it will be a constant matrix when $\{\mathrm{ds}\}$ is given. Thus, the linear system of uncertainty parameters is obtained

$$\{\begin{array}{c}\hfill \stackrel{·}{x}=A(\mathrm{ds})x+B(\mathrm{ds})u\hfill \\ \multicolumn{1}{c}{y=C(\mathrm{ds})x+D(\mathrm{ds})u}\end{array}$$

(10)

According to formula (9), $A(·),B(·),C(·)$ , and $D(·)$ indicate the matrix functions of the perturbation vector $\mathrm{E}=[\begin{array}{cccc}\hfill {\epsilon }_1\hfill & \hfill {\epsilon }_2\hfill & \hfill \cdots \hfill & \hfill {\epsilon }_n\hfill \end{array}{]}^T$ , namely

$$\{\begin{array}{c}\hfill A(·)=A+{\epsilon }_1{A}_1+\cdots +{\epsilon }_n{A}_n\hfill \\ \multicolumn{1}{c}{B(·)=B+{\epsilon }_1{B}_1+\cdots +{\epsilon }_n{B}_n}\\ \multicolumn{1}{c}{C(·)=C+{\epsilon }_1{C}_1+\cdots +{\epsilon }_n{C}_n}\\ \multicolumn{1}{c}{D(·)=D+{\epsilon }_1{D}_1+\cdots +{\epsilon }_n{D}_n}\end{array}$$

(11)

Let $\begin{array}{cc}\multicolumn{1}{c}{}& A(·)=A+E_a{\Delta }_a{F}_a\hfill \\ \multicolumn{1}{c}{}& B(·)=B+E_b{\Delta }_b{F}_b\hfill \\ \multicolumn{1}{c}{}& C(·)=C+E_c{\Delta }_c{F}_c\hfill \\ \multicolumn{1}{c}{}& D(·)=D+E_d{\Delta }_d{F}_d\hfill \end{array}$

where $A,B,C$ , and $D$ are real constant matrices with appropriate dimensions, which is the nominal model, namely, the models when the system uncertainties are neglected; ${\Delta }_a,{\Delta }_b,{\Delta }_c$ , and ${\Delta }_d$ are the matrices of uncertainty parameters which reflect the uncertainties of parameters in the system model; $E_a,E_b,E_c,E_d,F_a,F_b,F_c,\mathrm{and}{F}_d$ are constant matrices with appropriate dimensions, which reflect how the uncertain parameters affect the system model, namely, denote the uncertain structure of the model. Although the uncertainty parameters in the model are unknown, they vary in a bounded range and the span of the range directly affects the determination of system performance. Multiplying the correlated coefficient matrices by the appropriate scale matrices, to normalize the uncertainty matrices, namely, $\Vert {\Delta }_i\Vert \le 1,i=a,b,c,d$ . Where the norm takes the maximum singular value of the matrix.

$A(·)$ is taken as an example to illustrate the construction method of matrices $E_a,F_a$ , and ${\Delta }_a$ . Assuming the matrices $A_i$ is $m\times m$ dimensional matrices, and let

$$A_i={\left[\begin{array}{cccc}\hfill a_{11}^i\hfill & \hfill a_{12}^i\hfill & \hfill \cdots \hfill & \hfill a_{1m}^i\hfill \\ \multicolumn{1}{c}{a_{21}^i}& \hfill a_{22}^i\hfill & \hfill \cdots \hfill & \hfill a_{2m}^i\hfill \\ \multicolumn{1}{c}{⋮}& \hfill ⋮\hfill & \hfill \ddots \hfill & \hfill ⋮\hfill \\ \multicolumn{1}{c}{a_{m1}^i}& \hfill a_{m2}^i\hfill & \hfill \cdots \hfill & \hfill a_{\mathrm{mm}}^i\hfill \end{array}\right]}_{m\times m}$$

(12)

Then, the uncertainty portion of $A(·)$ can be decomposed into the form as shown by formula (13)

$$\begin{array}{cc}\multicolumn{1}{c}{}& \sum _{i=1}^n{\epsilon }_i{A}_i=\left[\begin{array}{ccc}\hfill \sum _{i=1}^n{\epsilon }_i{a}_{11}^i\hfill & \hfill \cdots \hfill & \hfill \sum _{i=1}^n{\epsilon }_i{a}_{1n}^i\hfill \\ \multicolumn{1}{c}{⋮}& \hfill \ddots \hfill & \hfill ⋮\hfill \\ \multicolumn{1}{c}{\sum _{i=1}^n{\epsilon }_i{a}_{n1}^i}& \hfill \cdots \hfill & \hfill \sum _{i=1}^n{\epsilon }_i{a}_{\mathrm{nn}}^i\hfill \end{array}\right]\hfill \\ \multicolumn{1}{c}{}& {=\stackrel{E_a}{\overbrace{\left[\stackrel{m\times n}{\overbrace{\begin{array}{cccc}\hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill 1\hfill & \hfill \cdots \hfill & \hfill 1\hfill \\ \multicolumn{1}{c}{0}& \hfill \cdots \hfill & \hfill 0\hfill \end{array}}}\hfill & \hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill 0\hfill & \hfill \cdots \hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{1}& \hfill \cdots \hfill & \hfill 1\hfill \end{array}}}\hfill & \hfill \cdots \hfill & \hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill 0\hfill & \hfill \cdots \hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill \cdots \hfill & \hfill 0\hfill \end{array}}}\hfill \\ \multicolumn{1}{c}{}& \hfill \hfill & \hfill \ddots \hfill & \hfill \hfill \\ \multicolumn{1}{c}{0}& \hfill 0\hfill & \hfill \cdots \hfill & \hfill \begin{array}{ccc}\hfill 1\hfill & \hfill \cdots \hfill & \hfill 1\hfill \end{array}\hfill \end{array}}}\right]}}}_{m\times \mathrm{mn}}\hfill \\ \multicolumn{1}{c}{}& \begin{array}{cc}\hfill \begin{array}{cc}\hfill \hfill & \hfill \hfill \end{array}\hfill & \hfill \hfill \end{array}\times {\stackrel{{\Delta }_a}{\overbrace{\left[\stackrel{m\times n}{\overbrace{\begin{array}{ccc}\hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill {\epsilon }_1\hfill & \hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \ddots \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \hfill & \hfill {\epsilon }_n\hfill \end{array}}}\hfill & \hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \ddots \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \hfill & \hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill {\epsilon }_1\hfill & \hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \ddots \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \hfill & \hfill {\epsilon }_n\hfill \end{array}}}\hfill \end{array}}}\right]}}}_{\mathrm{mn}\times \mathrm{mn}}\hfill \\ \multicolumn{1}{c}{}& \times \stackrel{F_a}{\overbrace{{\left[\stackrel{m\times n}{\overbrace{\begin{array}{ccc}\hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill a_{11}^1\hfill & \hfill \cdots \hfill & \hfill a_{11}^n\hfill \end{array}}}\hfill & \hfill \cdots \hfill & \hfill \stackrel{n}{\overbrace{\begin{array}{ccc}\hfill a_{m1}^1\hfill & \hfill \cdots \hfill & \hfill a_{m1}^n\hfill \end{array}}}\hfill \\ \multicolumn{1}{c}{}& \hfill \ddots \hfill & \hfill \hfill \\ \multicolumn{1}{c}{\begin{array}{ccc}\hfill a_{1m}^1\hfill & \hfill \cdots \hfill & \hfill a_{1m}^n\hfill \end{array}}& \hfill \cdots \hfill & \hfill \begin{array}{ccc}\hfill a_{\mathrm{mm}}^1\hfill & \hfill \cdots \hfill & \hfill a_{\mathrm{mm}}^n\hfill \end{array}\hfill \end{array}}}\right]}_{m\times \mathrm{mn}}^T}}\hfill \end{array}$$

(13)

The decomposition matrices of $B(·),C(·)$ , and $D(·)$ can be constructed similarly. Finally, appropriately process the output matrixes (13) and make $D=0$ . Thus, the uncertainty model of system can be expressed as

$$\begin{array}{cc}\multicolumn{1}{c}{\stackrel{·}{x}(t)}& =\mathrm{Ax}(t)+B_{1}w(t)+B_{2}u(t)\hfill \\ \multicolumn{1}{c}{z(t)}& =C_{1}x(t)+D_{12}u(t)\hfill \\ \multicolumn{1}{c}{y(t)}& =C_{2}x(t)+D_{21}w(t)\hfill \\ \multicolumn{1}{c}{w(t)}& =\Delta q(t)\hfill \end{array}$$

(14)

where $\begin{array}{cc}\multicolumn{1}{c}{B_1}& =[\begin{array}{cccc}\hfill E_a\hfill & \hfill E_b\hfill & \hfill 0\hfill & \hfill 0\hfill \end{array}]\hfill \\ \multicolumn{1}{c}{D_{21}}& =[\begin{array}{cccc}\hfill 0\hfill & \hfill 0\hfill & \hfill E_c\hfill & \hfill E_d\hfill \end{array}]\hfill \\ \multicolumn{1}{c}{C_1^T}& =[\begin{array}{cccc}\hfill F_a\hfill & \hfill 0\hfill & \hfill F_c\hfill & \hfill 0\hfill \end{array}]\hfill \\ \multicolumn{1}{c}{D_{12}^T}& =[\begin{array}{cccc}\hfill 0\hfill & \hfill F_b\hfill & \hfill 0\hfill & \hfill F_d\hfill \end{array}]\hfill \\ \multicolumn{1}{c}{\Delta }& =\mathrm{diag}\{{\Delta }_a,{\Delta }_b,{\Delta }_c,{\Delta }_d\},\Vert \Delta \Vert \le 1\hfill \end{array}$

For the mentioned uncertainty system, try to find an $n_c$ -order ( $n_c>0$ ) linear time-invariant dynamic controller

$$\begin{array}{cc}\multicolumn{1}{c}{{\stackrel{·}{x}}_c}& =A_c{x}_c+B_{c}y\hfill \\ \multicolumn{1}{c}{u}& =C_c{x}_c+D_{c}y\hfill \end{array}$$

(15)

where $x_c\in R^{n_c}$ represents the state variables of the controller, and y represents the deferent values of controlled object.

System would inevitably cause external disturbances $\omega (t)$ under the perturbations of parameters. The system output $z(t)$ is required to keep minimal, that is to search the optimal control $\{u^{*}(t)\}$ when $\Vert G_{\mathrm{wz}}(s){\Vert }_{\infty }$ takes the minimum. This problem can be converted into the solving of a linear matrix inequality,^18,19 namely

$$\begin{array}{cc}\multicolumn{1}{c}{}& \begin{array}{cc}\hfill min\hfill & \hfill \gamma \hfill \end{array}\hfill \\ \multicolumn{1}{c}{}& \begin{array}{cc}\hfill s.t.\hfill & \hfill \{\begin{array}{c}\hfill {\left[\begin{array}{cc}\hfill W_1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{W_2}& \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill I\hfill \end{array}\right]}^T\left[\begin{array}{ccc}\hfill \mathrm{AX}+{\mathrm{XA}}^T\hfill & \hfill {\mathrm{XC}}_1^T\hfill & \hfill B_1\hfill \\ \multicolumn{1}{c}{C_{1}X}& \hfill -I\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{B_1^T}& \hfill 0\hfill & \hfill -I\hfill \end{array}\right]\left[\begin{array}{cc}\hfill W_1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{W_2}& \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill I\hfill \end{array}\right]<0\hfill \\ \multicolumn{1}{c}{{\left[\begin{array}{cc}\hfill V_1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{V_2}& \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill I\hfill \end{array}\right]}^T\left[\begin{array}{ccc}\hfill A^{T}Y+\mathrm{YA}\hfill & \hfill {\mathrm{YB}}_1\hfill & \hfill C_1^T\hfill \\ \multicolumn{1}{c}{B_1^{T}Y}& \hfill -I\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{C_1}& \hfill 0\hfill & \hfill -I\hfill \end{array}\right]\left[\begin{array}{cc}\hfill V_1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{V_2}& \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill I\hfill \end{array}\right]<0}\\ \multicolumn{1}{c}{\left[\begin{array}{cc}\hfill X\hfill & \hfill I\hfill \\ \multicolumn{1}{c}{I}& \hfill Y\hfill \end{array}\right]\ge 0}\end{array}\hfill \end{array}\hfill \end{array}$$

(16)

This convex optimization problem can be solved by LMI toolbox in MATLAB.^20,21 Once $(X,Y)$ is solved, subsequently the parameter matrices $A_c,B_c,C_c$ , and $D_c$ of the H_∞ controller will be solved.

According to these, the control disciplinary analysis model of the system can be obtained which is described in formula (10): assuming a set of the central points ${\mathrm{ds}}_0$ and the corresponding perturbation range $[{\mathrm{ds}}_0-\Delta ,{\mathrm{ds}}_0+\Delta ]$ of the static design parameters in various disciplines have been given, Sys is short for the uncertainty model $\left(\begin{array}{ccccccc}\hfill A\hfill & \hfill B_1\hfill & \hfill B_2\hfill & \hfill C_1\hfill & \hfill C_2\hfill & \hfill D_{12}\hfill & \hfill D_{21}\hfill \end{array}\right)$ which is described in formula (14). The construction parameters of $\mathrm{ds}\to p_i$ response surface approximation model and central points ${\mathrm{ds}}_0$ of design parameters are used in the solving of Sys. With the application of formula (16), the optimal performance index $\gamma$ of robustness in control discipline is solved. The analysis model can be simply expressed as

$$\gamma =F_{\mathrm{cm}}(d_s,p_0,\Delta ,\mathrm{Sys})$$

(17)

MDO approach with robust controller in design loop

MDO problem description based on robust control

The form of MDO description based on robust control is exactly the same as the form of MDO problem description with control discipline in design loop described previously, only the parts related to control discipline are replaced by robust controller. In order to simplify the discussion, it is assumed that the system comprises only two disciplinary models. One is SM, while the other is control disciplinary model (CM for short). Multiple disciplinary models can be processed in a similar method. Thus, the description of MDO problem is given below:

• The design variables: $\{d_s,d_c\}$ , where $\{d_s\}$ and $\{d_c\}$ are the design variables of SM and CM, respectively.

• The objective function: $J_{\mathrm{tot}}=w_1\gamma +w_2{J}_s$ , where $\gamma$ indicates the performance index of H_∞, $J_s$ indicates the objective function of SM, and $w_1$ and $w_2$ indicate the weight coefficients.

• System analysis: $p_s=S(d_s,p_{s0})$ for SM, while the model as formula (14) for CM.

The limiting conditions

$$\begin{array}{cc}\multicolumn{1}{c}{}& g_c^{min}\le g_c(d_s,d_c)\le g_c^{max}\hfill \\ \multicolumn{1}{c}{}& g_s^{min}\le g_s(d_s)\le g_s^{max}\hfill \end{array}$$

where $g_c(d_s,d_c)$ and $g_s(d_s)$ are the limiting conditions of CM and SM, respectively.

Concurrent subspace optimization approach of MDO-C problem

Now the algorithms of concurrent subspace optimization (CSO) and collaboration optimization are the two effective MDO algorithms which are widely used in the MDO. When the design parameters $\{d_s\}$ of SM perturb at a specific sub-interval, CM conducts the solving computation on robust controller only once. In overall perspective, the optimal solution of robust controller is discrete in design space, especially in the linear system, the response surface approximation is necessary in the description of relationship between model parameters and design parameters. CSO algorithm is more appropriate for the solving of MDO-C problem. In the CSO algorithm, the model of response surface approximate is introduced to simplify the analysis of sub-disciplines, and the results of sub-disciplinary optimization are treated as the design points which can be applied in the further construction of response surface. In the iterative process, the accuracy of corresponding response surface keeps increasing until the design variables converge in the system coordination. Considering the characteristics of uncertainties of various disciplinary design parameters $\{d_s\}$ , this article makes proper improvements in the description of system uncertainty and further reduces the amount of calculation in the system analysis of control discipline. The solving process of the problem is illustrated in Figure 2.

Figure 2.

Concurrent subspace optimization algorithm diagram in MDO-C.

Example

Now taking a tailless aircraft design as sample, the article introduces robust controller for the control of pitch angle, brings control discipline into the design loop, and applies CSO algorithm to realize the concurrent design and optimization of aerodynamic discipline and control discipline. MDO description of the problem is given below.

The objective function: $J=w_1\gamma +w_{2}f+w_3{S}_{\mathrm{ele}\mathrm{va}\mathrm{tor}}$ .

where $f=\sqrt{[max(l_{\mathrm{axial}}){]}^2+[max(l_{\mathrm{span}}){]}^2}$ indicates the minimum equivalent size of aircraft; $S_{\mathrm{ele}\mathrm{va}\mathrm{tor}}$ indicates elevator area; $w_1,w_2$ , and $w_3$ indicate weight coefficients, which can transform all the variables in the objective function into the same order of magnitude.

The aerodynamic disciplinary design variables $\{d_s\}$ include wing span, wing root-chord length, wingtip-chord length, sweepback angle of wing leading edge, wing positive dihedral angle, wing installation angle, distance from leading edge to head, fuselage length, fuselage maximum diameter, position of fuselage maximum cross section, diameter at the end of fuselage, position of centroid, flat-tail span, flat-tail root-chord length, flat-tail wingtip-chord length, vertical tail span, vertical tail root-chord length, vertical tail wingtip-chord length, elevator area, and rudder area (20 design variables in total).

The control disciplinary design variables $\{d_c\}$ include $V_1,V_2,W_1,W_2,X$ and $Y$ are variables of robust controller as formula (16) shows.

The constraint conditions: apart from the constraints of basic performance, the following restrictions should be added:

• Hang time: 30 min

• Taper ratio: $\eta \in [\begin{array}{cc}\hfill 0.2\hfill & \hfill 1\hfill \end{array}]$

• Aspect ratio: $\lambda \in [\begin{array}{cc}\hfill 2\hfill & \hfill 8\hfill \end{array}]$

• Mass of whole aircraft: 100 g

• Engine power: 11.2 W

The ranges of design parameters are shown in Table 1.

Table 1.

Ranges of aircraft design parameters and perturbation parameters.

	Design variable	Boundary condition		Design variable	Boundary condition
		Lower bound	Upper bound		Lower bound	Upper bound
Aerodynamic disciplinary parameter	Wing span (m)	0.05	0.20	Diameter at the end of fuselage (m)	30%−90% maximum diameter
	Wing root-chord length (m)	15%−50% span		Position of centroid (m)	20%−40% fuselage length
	Wingtip-chord length (m)	≤Root-chord length		Tailplane span (m)	40%−80% wing span
	Sweepback angle of wing leading edge (°)	0	45	Tailplane root-chord length (m)	15%−50% tailplane span
	Positive dihedral angle of wing (°)	0	10	Multidisciplinary
	Wing installation angle (°)	0	5	Vertical tail span (m)	20%−40% wing span
	Fuselage length (m)	60%−120% span		Vertical tail root-chord length (m)	≤Tailplane root-chord length
	Distance from leading edge to head (m)	Conventional: 5%−30% fuselage length tailless: 0%−40% fuselage length canard: 40%−60% fuselage length		Vertical tail wingtip-chord length (m)	≤Vertical tail root-chord length
	Fuselage maximum diameter (m)	10%−30% fuselage length		Elevator area (m2)	20%−40% flat-tail area
	Position of fuselage maximum cross section (m)	20%−50% fuselage length		Rudder area (m2)	20%−40% area of vertical tail

In this article, CSO approach is adopted in the optimization research, and the iterative process of optimization is given in Figure 3. As Figure 3 shows, the value of objective function converges after 294 times of iteration, and it has fallen into the optimizing interval after iterating almost 156 times. There is no obvious improvement for optimized values after extra 138 times of iteration. These illustrate that CSO can approach the vicinity of the optimal solution at a faster speed. Since the controllers adopted in control discipline have robustness for parameter perturbations, CSO approach needs to simply divide the variation ranges of these design parameters into several intervals and execute the system analysis in finite times, so that the calculated amount of MDO will be greatly reduced.

Figure 3.

Iterative process of aerodynamic and control disciplinary concurrent optimization.

Finally, the optimized result of the controller is obtained as equation (18). The designed optimal controller is a two-order linear time-invariant dynamic controller, and the controller manipulates the pitch angle by elevator. In this way, the impact of external disturbances on system output can be reduced to minimum level, as $\Vert G_{\mathrm{wz}}(s){\Vert }_{\infty }$ least

$$\begin{array}{cc}\multicolumn{1}{c}{}& {\mathbf{A}}_{\mathrm{c}}=\left[\begin{array}{cc}\hfill -180.82\hfill & \hfill 325.5\hfill \\ \multicolumn{1}{c}{-0.62}& \hfill 1.10\hfill \end{array}\right],{\mathbf{B}}_{\mathrm{c}}=\left[\begin{array}{c}\hfill -1.54\hfill \\ \multicolumn{1}{c}{0.094}\end{array}\right]\hfill \\ \multicolumn{1}{c}{}& {\mathbf{C}}_{\mathrm{c}}={\left[\begin{array}{c}\hfill -545.28\hfill \\ \multicolumn{1}{c}{983.36}\end{array}\right]}^T,{\mathbf{D}}_{\mathrm{c}}=0\hfill \end{array}$$

(18)

The assessment size of equivalent to be 167 mm, Compared with the MDO optimizing result which is obtained without considering control discipline in the design loop, this value increases a little. However, the design parameters in the area of performance index have robustness, and the results are more rational and effective.

The time response of total speed, flight-path angle, and pitch angular rate on ideal system, the nominal closed-loop system, and perturbed system with the variation in system parameters is shown in Figure 4. Assuming that unmanned air vehicle (UAV) straight and level fly at initial value of 11.28 m/s, the step command of total speed command signal with amplitude 12.19 m/s at t = 5 s lasts for 5 s; then at t = 10 s, the flight-path angle command signal with 5 m/s amplitude and duration of 4 s is given. The response curve is shown in Figure 4. Among them, the ideal system’s response is the step curve to satisfy the performance index requirements, and the nominal closed-loop system’s response is the real-time commands’ tracking curves with the zero disturbances of design parameters, and the disturbance system refers to the response curve with maximum uncertainty disturbance of design parameters. As shown in the figure, whether it is nominal response or disturbance response, the response curve is pretty close; they all can track speed and angle control command very well, which shows the good robust performance of the system.

Figure 4.

A tailless aircraft command step response curve on different system model.

Conclusion

Based on MDO approach research with control discipline in the design loop, this article proposes an effective method to design robust controller of aircraft and takes an example to validate this method. The results show that (1) considering the complexity of the mapping relationship between the coefficients in the aircraft motion model and the various design parameters, the uncertainty description based on RSM of design parameters under aircraft multidisciplinary design environment is a effective approach with less conservative; (2) the LMI-based robust controller design approach can simplify the solve complication of control law and improve the convergence rate of robust controller optimization; (3) with the employment of CSO approach, the synchronous design of aircraft MDO can be realized, and the robustness of flight performance index is strengthened. Therefore, this approach can be widely used in MDO design with control discipline in other fields.