Computation of robustly stabilizing PID controllers for interval systems

The paper is focused on the computation of all possible robustly stabilizing Proportional-Integral-Derivative (PID) controllers for plants with interval uncertainty. The main idea of the proposed method is based on Tan’s (et al.) technique for calculation of (nominally) stabilizing PI and PID controllers or robustly stabilizing PI controllers by means of plotting the stability boundary locus in either P-I plane or P-I-D space. Refinement of the existing method by consideration of 16 segment plants instead of 16 Kharitonov plants provides an elegant and efficient tool for finding all robustly stabilizing PID controllers for an interval system. The validity and relatively effortless application of presented theoretical concepts are demonstrated through a computation and simulation example in which the uncertain mathematical model of an experimental oblique wing aircraft is robustly stabilized.

investigation of robust stability for systems with parametric uncertainty-see e.g. (Barmish 1994;Bhattacharyya et al. 1995Bhattacharyya et al. , 2009. Typical problem of practical PI(D) controller design is to ensure, that the calculated controller will guarantee stability not only for one assumed nominal controlled system but also for the whole family of systems described by a model with parametric uncertainty. Such closed-loop control system is called as "robustly stable" and the controller itself is then robustly stabilizing one.
An array of techniques for calculation of (nominally) stabilizing PI and PID controllers have been already published, such as rules presented in (Söylemez et al. 2003), the Tan's method described in (Tan and Kaya 2003;Tan et al. 2006) or the Kronecker summation method from (Fang et al. 2009). Furthermore, these methods have been also extended for robust stabilization of interval plants by their combination with the sixteen plant theorem (Barmish et al. 1992;Barmish 1994). Nevertheless, this extension works only for PI but not for PID controllers.
The main aim of this paper is to present a method for computation of all possible robustly stabilizing PID controllers for interval plants and to demonstrate its serviceability by robust stabilization of an oblique wing aircraft model. More specifically, the goal is to refine the elegant and effective Tan's method (Tan and Kaya 2003;Tan et al. 2006) by the ideas from (Ho et al. 1998(Ho et al. , 2001, i.e. to use 16 segment plants instead of 16 Kharitonov plants, and to make the existing method applicable for computation of robustly stabilizing PID controllers. Previously, the computation of all (nominally) stabilizing PI or PID controllers, robustly stabilizing PI controllers and consequent choice of the specific controller with desired performance on the basis of the desired model method (formerly known as dynamics inversion method) (Vítečková 2000) is shown in (Matušů 2011). Then, the application of Kronecker summation method (Fang et al. 2009) to robust stabilization of a chemical reactor or robust stabilization of a third order nonlinear electronic model is given in  or (Matušů et al. 2010a), respectively. The robust stabilization of the same nonlinear electronic plant using the Tan's method (Tan and Kaya 2003;Tan et al. 2006) is presented e.g. in (Matušů et al. 2010b).
The paper is organized as follows. In "Computation of (nominal) stability regions for PI controllers" section, a graphical method for computation of (nominally) stabilizing PI controllers is recalled. "Computation of (nominal) stability regions for PID controllers" section has the same purpose but for PID controllers. Next, the computation of robustly stabilizing PI controllers for interval plants is presented in "Robust stabilization using PI controllers" section. The key "Robust stabilization using PID controllers" section extends the existing Tan's (et al.) method, combines it with the segment plants concept and makes it applicable for calculation of robustly stabilizing PID controllers. Further, the extensive "Illustrative example: Robust stabilization of oblique wing aircraft" section confirms the obtained results by means of the simulation example with an experimental oblique wing aircraft model. And finally, "Conclusion" section offers some conclusion remarks.

Computation of (nominal) stability regions for PI controllers
First, the fundamentals related to computation of (nominal) stability regions for PI controllers are going to be summarized.
Suppose the classical closed-loop control system according to Fig. 1, where C(s) represents a controller, G(s) stands for a controlled system, and signals w(t), e(t), u(t) and y(t) denote a reference value, tracking (control) error, actuating (control) signal and output (controlled) variable, respectively.
The controller is assumed in the well-known PI form: where k P , k I represent the proportional and integral gain, respectively. The principal task is to determine the parameters k P , k I which guarantee stabilization of the controlled plant: Several effective methods for the computation of stabilizing PI controllers have been already published, e.g. (Söylemez et al. 2003;Tan and Kaya 2003;Tan et al. 2006;Fang et al. 2009). Here, the Tan's method from (Tan and Kaya 2003;Tan et al. 2006) will be revisited and extended. This graphical approach is based on plotting the stability boundary locus. The substitution of s for jω in the plant transfer function (2) and subsequent decomposition of the numerator and denominator into their even and odd parts result in: Further, expressing the closed-loop characteristic polynomial and equating both real and imaginary parts to zero lead to the relations for the proportional and integral gains k P , k I : k P (ω) = P 5 (ω)P 4 (ω) − P 6 (ω)P 2 (ω) P 1 (ω)P 4 (ω) − P 2 (ω)P 3 (ω) k I (ω) = P 6 (ω)P 1 (ω) − P 5 (ω)P 3 (ω) P 1 (ω)P 4 (ω) − P 2 (ω)P 3 (ω) Simultaneous calculations of the Eq. (4) for a suitable range of ω and plotting the obtained values into the (k P , k I ) plane determine the stability boundary locus. The obtained curve together with the line k I = 0 split the (k P , k I ) plane into the stable and unstable regions. The decision if the respective region represents stabilizing or unstabilizing area can be done simply using a test point within each region. Nonetheless, the appropriate frequency gridding could represent a potential problem. Thus, the Nyquist plot based technique from (Söylemez et al. 2003) can be used for the improvement of the method. In this improvement, the frequency ω can be separated into several intervals within which the stability or instability can not change. The borders of such intervals are defined by the real values of ω which fulfill the equation: The obtained intervals could be helpful for the proper frequency scaling.

Computation of (nominal) stability regions for PID controllers
Now, the issue of (nominal) feedback stabilization will be elaborated again, but for the case of ideal PID controller given by: The principal idea for obtaining the relevant stability regions is to fix one controller parameter to a certain value and calculate the stability boundary locus using two remaining parameters analogously to the procedure presented in the previous "Computation of (nominal) stability regions for PI controllers" section.
The expression for the stability boundary locus in the (k P , k I ) plane for a fixed value of k D leads to a bit modified equations for proportional and integral gains: where Note that the last two terms in (9) depend on derivative constant k D . From the viewpoint of practical computation, k D is considered to be chosen and the corresponding set of boundary parameters k P , k I is consequently calculated while this process is repeated for several selected values of k D . Thus, the final stability regions are successively plotted through the "(k P , k I ) sections" in the (k P , k I , k D ) space.
Alternatively, the stability boundary locus in the (k P , k D ) plane for a fixed value of k I can be computed. This scenario would change the Eqs. (8) and (9) to, respectively: where Obviously, the final stability regions are given by the "(k P , k D ) sections" in the (k P , k I , k D ) space.
Nonetheless, the third option of obtaining the stability boundary, which consists in fixing k P and calculating the curves in (k I , k D ) plane, is not so straightforward as the previous two alternatives, because for this case it holds true: However, the stability region in the (k I , k D ) plane for a fixed k P can be acquired using the stability region in the (k P , k I ) plane and (k P , k D ) plane together as it has been presented in (Tan et al. 2006). In accordance with a linear programming based approach from (Ho et al. 1997), the stability region in the (k I , k D ) plane under fixed k P is a convex polygon which can be sometimes advantageous for easier plotting.

Robust stabilization using PI controllers
So far, the previous Sections were focused only on nominal stabilization of controlled plants with fixed parameters. The following parts will deal with robust stabilization. It means that the plant whose coefficients can vary within given intervals (interval system) is considered to represent a controlled object and the aim is to find all controllers which assure stabilization of all possible members of the interval system family. The works (Tan and Kaya 2003;Tan et al. 2006;Fang et al. 2009) have improved a feedback stabilization technique using PI controllers also for interval plants simply by using its combination with the sixteen plant theorem (Barmish et al. 1992;Barmish 1994;Ho et al. 1998). The sixteen plant theorem itself says that a first order controller robustly stabilizes an interval plant: k P (ω, k I ) = P 5 (ω)P 4 (ω) − P 6 (ω)P 2 (ω) where i, j ∊ {1, 2, 3, 4}; and B 1 (s) to B 4 (s) and A 1 (s) to A 4 (s) are the Kharitonov polynomials for the numerator and denominator of the interval plant (13).
Recall that the construction of Kharitonov polynomials e.g. for the numerator interval polynomial: is based on the use of the lower and upper bounds of interval parameters in compliance with the principle (Kharitonov 1978): Consequently, the robust stabilization of an interval plant directly follows from the simultaneous stabilization of all 16 fixed Kharitonov plants. Hence, the final area of stability for original interval plant is given by the intersection of all 16 related partial areas obtained individually using the techniques from the "Computation of (nominal) stability regions for PI controllers" section.

Robust stabilization using PID controllers
Unfortunately, the sixteen plant theorem is not applicable for robust stabilization of interval systems by PID controllers as it is not valid anymore (Pujara and Roy 2001). However, the suitable method based on the generalized Kharitonov theorem and linear programming techniques has been presented e.g. in (Ho et al. 1998(Ho et al. , 2001. This paper adopts the idea of Kharitonov segments used in the generalized Kharitonov theorem (Chapellat and Bhattacharyya 1989) and similar thirty-two edge theorem (Barmish 1994; Chapellat and Bhattacharyya 1989) and combines it with the stability boundary locus technique (Tan and Kaya 2003;Tan et al. 2006).
Consider an interval plant: and a PID controller (7). The family of interval systems (17) is stabilized by a fixed PID controller if and only if each of sixteen segment plants related to the interval family is stabilized by the same PID controller (Ho et al. 2001).
The mentioned 16 segment plants are defined as: where i, j ∊ {1, 2, 3, 4}; A 1 (s) to A 4 (s) are the Kharitonov polynomials for the denominator of the interval plant (17); and B S1 (s, λ) to B S4 (s, λ) are four Kharitonov segments (Chapellat and Bhattacharyya 1989;Ho et al. 2001;Barmish 1994) which can be written as: where λ ∊ 〈0, 1〉 and B 1 (s) to B 4 (s) are the Kharitonov polynomials for the numerator of the interval plant (17). The computation of robustly stabilizing PID controllers can be performed as follows: First, a certain value of controller parameter k D is chosen and fixed (alternatively, also the parameter k I or k P can be fixed according to "Computation of (nominal) stability regions for PID controllers" section, but the fixed k D is supposed here). Then, the stability boundary for one of segment plants (18) is calculated for several sampled values of λ ∊ 〈0, 1〉 using the Eqs. (8) and (9). The intersection of the obtained areas in (k P , k I ) plane gives the stability boundary locus for this specific segment plant. The calculations are repeated for all the remaining segment plants and the robust stability region for the original interval plant and chosen value of k D is determined by the intersection of areas for all 16 segment plants. From the practical viewpoint, the curves for all sampled λ ∊ 〈0, 1〉 and all 16 segment plants can be plotted in one figure and intersection can be found at a time. Anyway, the whole process should be repeated for the other selected values of k D and the very final robust stability region can be visualized by the simultaneous plotting of the "(k P , k I ) sections" into one graph in (k P , k I , k D ) space.

Illustrative example: Robust stabilization of oblique wing aircraft
This Section is intended to practically demonstrate the theoretical results from the previous parts by means of the illustrative example.
The controlled plant is supposed to be given by the uncertain mathematical model of an experimental oblique wing aircraft from (Barmish 1994;Dorf 1974): and the aim is to find all robustly stabilizing PI and PID controllers.

PI controller
The first of the Kharitonov plants constructed according to (14) is:  (4) and (5) for (22) and plotting the obtained results together with the line k I = 0 into the (k P , k I ) plane lead to the stability region which is visualized in Fig. 2. A decision on which part represents the area of stability can be simply done using an arbitrary pair (k P , k I ) from this region, calculating the corresponding closed-loop characteristic polynomial and verifying its stability. In this case, the stabilizing area lies inside the depicted shape.
The similar graphs can be plotted for remaining Kharitonov plants. The curves for all 16 Kharitonov plants are shown in Fig. 3.
The final robust stability region is determined by the intersection of regions for all 16 Kharitonov plants. It is zoomed and highlighted in Fig. 4.
Three PI controllers with specific position in relation to stability region from Fig. 4 have been chosen, namely: Their location is shown in Fig. 5.  The control results obtained by means of the selected controllers visually confirm the validity of robust stability region. The Fig. 6 shows the control responses of the loop with the controller C 1 (23) and 729 "representative" systems from the interval family (20). Each interval parameter has been divided into 2 subintervals and thus these 3 values and 6 parameters have resulted in 3 6 = 729 systems for simulation. Moreover, the red curve represents the output variable for the system with average values of uncertain parameters from (20): As can be clearly seen, some members of the interval family (20) are not stabilized by the controller C 1 (23) and thus the closed-loop system is really robustly unstable. The analogical simulations for controller C 2 (24) and 729 + 1 "representative" systems lead to the set of control responses from Fig. 7. In this case, the closed-loop system is on the stability border for some member of the interval family (20) (the worst case) which concurs with the position of the controller C 2 in Fig. 5. Finally, the same set of control responses is plotted in Fig. 8 for controller C 3 (25). Now, the control loop is obviously robustly stable.

PID controller
In this part, all robustly stabilizing PID controllers are going to be found for the same oblique wing aircraft model (20).
Initially, the derivative constant is chosen and fixed as k D = 1. The first of the segment plants (18) is constructed using: where B 1 (s), B 3 (s) and A 1 (s) are relevant Kharitonov polynomials and λ ∊ 〈0, 1〉. More specifically: (26) G A (s) = 64s + 128 s 4 + 3.7s 3 + 65.6s 2 + 32s The stability boundary locus for the segment plant (28) is given by the intersection of the stability areas for several sampled values of λ ∊ 〈0, 1〉. The Fig. 9 shows 11 curves for the range λ = 0:0.1:1 and the highlighted area represents the intersection.
The same process can be analogously repeated for the remaining 15 segment plants and then the intersection of all 16 partial intersections would lead to the stability boundary locus for the original interval family (20) under the assumption of k D = 1.
The Fig. 10 presents the curves for all 16 segment plants and sampled λ = 0:0.1:1 in a single plot. The zoomed and highlighted intersection representing the robust stability region for closed loop containing PID controller with k D = 1 and original interval system (20) is depicted in Fig. 11. The whole previous process was repeated for the other values of derivative constant k D , more specifically for k D = 0:0.5:5. The very final robust stability region visualized by means of corresponding eleven "(k P , k I ) sections" in (k P , k I , k D ) space is shown in Fig. 12.
The example of robustly stabilizing PID controller, which obviously lies inside the robust stability region from Fig. 12, can be chosen as: The Fig. 13 demonstrates the set of control responses of the loop with the controller C 4 (29) and 729 + 1 "representative" systems from the interval family (20) obtained analogously as in the previous PI control cases. As can be seen, the control loop is really robustly stable.
(29) C 4 (s) = k P + k I s + k D s = 1 + 0.5 s + 0.5s  The intentional choice of a controller outside the robust stability region from Fig. 12 evidently leads to the robustly unstable closed control loop. The representative of such robustly non-stabilizing controller is e.g.: Then, the corresponding set of control responses obtained under the same conditions as in the previous case is depicted in Fig. 14.

Conclusion
The main aim of paper has been to present the improved method for computation of stabilizing controllers with the conventional structure on the basis of plotting the stability boundary locus in either P-I plane or P-I-D space. Now, thanks to the combination of the original method with stabilization of so-called segment plants, the modified technique can be conveniently used for determination of all possible robustly stabilizing PID controllers for interval plants. In the illustrative example, the model of an experimental oblique wing aircraft is considered as a controlled object. Two final robust stability regions have been computed and visualized, one for PI and the other for PID controller, and selected representatives from stable or even intentionally unstable areas have been chosen and used for supporting control simulations.