Table of Contents Author Guidelines Submit a Manuscript
Mathematical Problems in Engineering
Volume 2008, Article ID 817063, 26 pages
Research Article

Suboptimal Regulation of a Class of Bilinear Interconnected Systems with Finite-Time Sliding Planning Horizons

1Department of Electricity and Electronics, Institute of Research and Development of Processes (IIDP), Faculty of Science and Technology, University of the Basque Country, Leioa (Bizkaia), P.O. Box 644, 48080 Bilbao, Spain
2Department of Automatic Control and Systems Engineering, College of Industrial Technical Engineering (EUITI) Bilbao, University of the Basque Country, Bilbao (Bizkaia), Plaza de la Casilla 3, 48012 Bilbao, Spain
3Department of Applied Mathematics, College of Industrial Technical Engineering (EUITI) Bilbao, University of the Basque Country, Bilbao (Bizkaia), Plaza de la Casilla 3, 48012 Bilbao, Spain

Received 3 November 2006; Revised 20 July 2007; Accepted 19 December 2007

Academic Editor: Angelo Luongo

Copyright © 2008 M. de la Sen et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


This paper focuses on the suboptimization of a class of multivariable discrete-time bilinear systems consisting of interconnected bilinear subsystems with respect to a linear quadratic optimal regulation criterion which involves the use of state weighting terms only. Conditions which ensure the controllability of the overall system are given as a previous requirement for optimization. Three transformations of variables are made on the system equations in order to implement the scheme on an equivalent linear system. This leads to an equivalent representation of the used quadratic performance index that involves the appearance of quadratic weighting terms related to both transformed input and state variables. In this way, a Riccati-matrix sequence, allowing the synthesis of a standard feedback control law, is obtained. Finally, the proposed control scheme is tested on realistic examples.

1. Introduction

Dynamic bilinear systems have received great attention by researchers in the last decades, from the classical works of Anderson and Moore [1], Feldbaum [2], Tarn [3], or Tarn et al. [4], to the recent works of Al-Baiyat [5], Kotta et al. [6] or Garrido et al. [7]. They form a transitional class between the linear and the general nonlinear systems. The importance of such systems is well-known in diverse areas in which bilinear systems appear naturally for a variety of dynamic processes such as those taking place in electrical systems [8], industrial processes [9], power plants [10], chemical processes [11], nuclear fusion [12], biomedical applications [13], information management [14], mechanical systems [15], and aerospace and avionics [16]. which can not be satisfactorily modeled under the assumption of linearity. On the other hand, the design of controllers for bilinear systems has been an area of major research during the recent years. This growing interest in practical applications requires the development of suitable algorithms for control problems associated with bilinear systems. Many of the results obtained rely on the optimization theory with an appropriate performance index, as it is the case of the proposed control scheme.

Well-known works as Goka et al. [17] and Tarn et al. [18] stated the controllability conditions for standard classes of discrete bilinear systems. In the first paper, an equivalent system description was derived with the equivalent input being dependent on the products of the state and the original input. Since the equivalent feed-forward loop is of a linear nature, the analysis becomes greatly simplified. In the second paper, the uniform controllability of such systems using bounded input was studied. On the other hand, it is well known that the systems are often interconnected and, in many cases, several dynamic subsystems can be distinguished for analysis purposes. This sometimes takes place in computer communication, transportation networks, control and power stations, and so on. There exists important literature dealing with both the associate multivariable and the large-scale problems [1925]. For instance, standard adaptive control techniques were applied in [26] to compensate the undesirable deviations of the process parameters from their nominal values, being the overall process modeled as a bilinear system prior to the application of the adaptive scheme.

This paper reports suboptimal optimization techniques which are applied to bilinear models. Such models can be considered as direct extensions to the linear continuous interconnected systems stated by the work of Ramakrishna and Viswanadham [27]. A class of invariant discrete-time multivariable bilinear systems with interconnection subsystems is studied. First, the system is shown to be equivalently described by the linear feed-forward multivariable structure with multiplicative inputs including a deterministic disturbance vector. Also, other interpretations are stated. As a previous requirement for optimization, controllability results for the overall systems performance are given by extending those ones given in [17, 18, 2731]. This is achieved through the above equivalent linear system with the use of centralized control methods. A central coordinator is supposed to be available for the local controllers to supply information to each control station from the remaining ones. That control technique is applied to the optimization of the aforementioned system class. The used performance index involves state-weighting terms only. An explicit solution of Riccati type (which is unusual for such an optimization criterion) is found out as the suboptimal solution by using manipulations on the input/state variables of the problem statement [1, 32]. The importance of this strategy arises from problems where constrains on the input rather than input weights are introduced in the optimization criterion, what leads to a feedback-type control law.

The paper is organized as follows. In Section 2, the problem statement for the class of the interconnected invariant discrete-time bilinear systems is given. Also, controllability necessary conditions obtained by the decomposition techniques are studied. Section 3 is devoted to the development of suboptimal strategies for the regulation criterion, which involves quadratic terms of the state variables. Three transformations of variables are used to turn the loss function into another one, which involves state and input quadratic terms leading to a Riccati-type suboptimal solution. The control strategies may be considered as extensions of those ones used in [1, 29] or [25], for the linear continuous problem. The formulation is accomplished from a centralized viewpoint. In Section 4, two simulated examples are presented. The first one consists of a simulated numerical example; and the second one consists of a real regulation problem of two amplidynes acting on a DC-motor. In the second example, it is not necessary to consider constant current or voltage as it is usual in DC-motor control, with their product—namely, the electric par—a bilinear term describing the resulting dynamical system. Finally, conclusions end the paper. Appendix A is related to the meaning explanation of certain equivalent interpretations of the overall system structure. Appendix B concerns with the development of the transformations of variables involved in the redefinition of the performance index. By convenience, equations and results from the appendices are sometimes invoked in the main body of the paper order not to repeat mathematical material.

2. On the System Structure and Its Controllability Properties

2.1. Controllability of a Class of Bilinear Systems

In [17], the following typical class of bilinear systems was considered: 𝑥(𝑘+1)=𝐴+𝑢(𝑘)𝐵𝑥(𝑘)+𝐶𝑢(𝑘);𝑘=0,1,2,,𝑥()𝑅𝑛;𝑢()𝑅.(2.1) Equation (2.1) can equivalently be described, if rank (𝐵,𝐶)=1 and then rank (𝐵)=1 if 𝐶 and are both nonzero, with 𝐵=𝐶𝑇 (𝐶, being unique 𝑛-vectors within a scalar factor), by the linear feed-forward system with multiplicative input: 𝑥(𝑘+1)=𝐴𝑥(𝑘)+𝐶𝜈(𝑘);𝜈(𝑘)𝑢(𝑘)𝑇𝑥(𝑘);if𝐶=0,𝑢(𝑘)[𝑇𝑥(𝑘)+1];if𝐶0(2.2) with 𝑘=0,1,2,; and where 𝑥(𝑘) and 𝑢(𝑘) are the state and the control at time 𝑘, respectively; and 𝐴 and 𝐵 are real constant coefficient matrices of compatible dimensions.

Note that for arbitrary ranks of 𝐵 or (𝐵,𝐶) in (2.1), the decomposition 𝐵=𝑝𝑖=1𝑒𝑖𝑏𝑇𝑖 can be used, with 𝑒𝑖 and 𝑏𝑇𝑖   being, respectively, the 𝑖th unity vector (i.e., the 𝑖th component is unity and the remaining ones are zero) and the 𝑖th row vector of 𝐵. The controllability of the bilinear system was related in that paper to that of the linear system (2.2) as follows.

Theorem 2.1. (summary of results in [17]). For the equivalent systems (2.1)-(2.2), the following controllability results hold.
(i)If the bilinear system (2.1) is controllable, then the pair (𝐴,𝐶) in (2.2) is controllable.(ii)If (𝐴𝑇,) is not controllable, then the homogeneous bilinear system is not controllable.(iii)If (𝐴,𝐶) and (𝐴𝑇,) are both controllable pairs and if 𝑇𝐶, 𝑇𝐴1𝐶0, then the homogenous bilinear system is controllable in 𝑅𝑛 outside the origin.(iv)The inhomogeneous bilinear system is not controllable in 𝑅𝑛 if(a)rank [controllability matrix of (𝐴𝑇,)] = rank [controllability matrix of (𝐴𝑇,);1]; (1𝑛-vector with every element being unity),(b)𝑛𝑖=1𝑎𝑖=1, with a() being the elements of the last row of the 𝐴-matrix in the canonical state-space representation of system (2.2).(v)The inhomogeneous bilinear system is controllable in 𝑅𝑛, if (𝐴,𝐶) is a controllable pair, 𝑇𝐶, 𝑇𝐴1𝐶0, and at least one of the conditions (iv-a) or (iv-b) is not fulfilled.
Although these results are restricted to a special class of discrete bilinear systems, many natural discrete bilinear systems [10, 15, 3335] do satisfy these assumptions. More general controllability conditions for bilinear systems can be found in [18, 36].

Remark 2.2. More general results were derived in [30] for the case in which 𝐴1 does not exist. Such a paper studies the controllability conditions for three cases of violation of the conditions given in [37] (namely, (𝐴,𝐶) and (𝐴𝑇,) are controllable pairs and 𝑘0=g.c.d{𝑖/𝑇𝐴𝑖1𝐶0,0<𝑖𝑛2}=1) (g.c.d denotes the greatest common divisor). Such cases include the following:
(a)rank (𝐴,𝐶)=𝑛, rank Ω(𝑇,𝐴)<𝑛;(b)rank (𝐴,𝐶)<𝑛, rank Ω(𝑇,𝐴)=𝑛; and(c)rank (𝐴;𝐶)<𝑛, rank Ω(𝑇,𝐴)<𝑛, with (𝐴,𝐶)[𝐶,𝐴𝐶,,𝐴𝑛1𝐶] and Ω(𝑇,𝐴)[,𝐴𝑇,

2.2. Models for a Class of Bilinear Systems with Interconnected Subsystems

In this paper, the following multivariable invariant discrete-time structure of bilinear systems (,𝐴(𝑛1)𝑇]𝑇.) of homogeneous type is considered 𝑆𝑖 with the interconnection bilinear subsystems (𝑆𝑖𝑥𝑖(𝑘+1)=𝐴𝑖+𝑢𝑖(𝑘)𝐵𝑖𝑥𝑖(𝑘)+𝐷𝑖𝑧𝑖(𝑘),(2.3)𝑦𝑖(𝑘)=𝐶𝑖𝑥𝑖(𝑘);𝑖𝐼;𝐼=1,2,,𝑝;𝑘=0,1,2,(2.4)), which include linear and bilinear coupling terms, given by 𝐻i where 𝐻𝑖𝑧𝑖(𝑘+1)=𝑀𝑖+𝑢𝑖(𝑘)𝑃𝑖𝑧𝑖(𝑘)+𝑝𝑗=1𝐿𝑖𝑗+𝑢𝑖(𝑘)𝑀𝑖𝑗𝑦𝑗(𝑘);𝑖𝐼;𝑘=0,1,2,,(2.5) are, respectively, the state vector of the 𝑥𝑖(𝑘)𝑅𝑛𝑖,𝑧𝑖(𝑘)𝑅𝑎𝑖,𝑦𝑖(𝑘)𝑅𝑟𝑖,𝑢𝑖(𝑘)𝑅-subsystems, the state vectors of the interaction 𝑆𝑖-subsystems, the corresponding output and input vectors at time 𝐻𝑖, and the different coefficient matrices of appropriate dimensions (see Figure 1).

Figure 1: The class of multivariable discrete invariant-time bilinear subsystems with the interaction bilinear subsystems.

The linear continuous models of [27] for Figure 1 are more general than those ones above since additional measurements from the interconnection subsystems are used to generate the inputs to the plant. Structures of this type arise from practical systems such as countercurrent heat exchangers having weak nonlinearities modeled as bilinear.

Now, the extended state vector of the overall system is defined as 𝑘 all 𝑥(𝑘)𝑥𝑇1(𝑘),𝑥𝑇2(𝑘),,𝑥𝑇𝑝(𝑘)𝑇𝑅𝑛+𝑎,where𝑥𝑖(𝑘)𝑥𝑇𝑖(𝑘),𝑧𝑖𝑇(𝑘)𝑇𝑅𝑛𝑖+𝑎𝑖,(2.6) with 𝑖𝐼 and the extended output vector is defined as 𝑛𝑝𝑗=1𝑛𝑗,𝑎𝑝𝑗=1𝑎𝑗, with 𝑦(𝑘)[𝑦𝑇1(𝑘),𝑦𝑇2(𝑘),,𝑦𝑇𝑝(𝑘)]𝑇𝑅𝑟, Thus, the following multivariable structure, equivalent to (2.3)–(2.5), is deduced: 𝑟𝑝𝑗=1𝑟𝑗. where the matrices of parameters have the following structures: 𝑥(𝑘+1)=𝐴𝑥(𝑘)+𝑝𝑖=1𝑢𝑖(𝑘)𝐵𝑖𝑥(𝑘);𝑘=0,1,2,3,,(2.7)𝑦(𝑘)=𝐶𝑥(𝑘);𝑘=0,1,2,3,,(2.8) for all 𝐴=𝐴𝑇1,𝐴𝑇2,,𝐴𝑇𝑝𝑇;𝐴𝑖=𝐴𝑖1,𝐴𝑖2,,𝐴𝑖𝑖,,𝐴𝑖𝑝𝐴𝑖𝑖=[𝐴𝑖𝐷𝑖𝐿𝑖𝑖𝐶𝑖𝑀𝑖]𝐴𝑖𝑗=[00𝐿𝑖𝑗𝐶𝑗0](for𝑗𝑖)𝐵𝑖=0𝑖1,,0,𝐵𝑇𝑖,0𝑝1,,0𝑇;𝐵𝑖=𝐵𝑖1,𝐵𝑖2,,𝐵𝑖𝑝𝐵𝑖𝑖=[𝐵𝑖0𝑀𝑖𝑖𝐶𝑖𝑃𝑖]𝐵𝑖𝑗=[00𝑀𝑖𝑗𝐶𝑗0](for𝑗𝑖)𝐶=𝐶𝑇1,𝐶𝑇2,,𝐶𝑇𝑝;𝐶𝑖=0𝑖1,,0,𝐶𝑖,0,0𝑝1,,0(2.9)

Remark 2.3. In general, the number of subsystems 𝑗,𝑖𝐼. and the number of interaction subsystems 𝑆() are not equal. Thus, assume that 𝐻() and 𝑝1 are, respectively, the number of subsystems 𝑝2(𝑝1𝑝2) and the number of subsystems 𝑆𝑡,𝑡𝐼1={1,2,3,,𝑝1}, Then, an extended state vector 𝐻𝑠,𝑠𝐼2={1,2,,𝑝2}. can be obtained similarly, if 𝜉() and 𝑝=max(𝑝1𝑝2)
Taking into account the dimensions and the structural zeros of the 𝐼=𝐼1𝑈𝐼2.-matrices, they can be decomposed into a sum of dyads as 𝐵𝑖 where 𝐵𝑖=𝑛𝑖+𝑎𝑖𝑗=1𝑒𝑗+𝑙𝑖1𝑏𝑇𝑗+𝑙𝑖1(𝑖);forall𝑖𝐼,(2.10) is the unity 𝑒𝑗 vector, 𝑗th is the 𝑏𝑇𝑗(𝑖) vector of the 𝑗th-matrix, 𝐵𝑖, and the lower subscripts 𝑖𝐼 are related to the first nonstructural zero rows of the 1𝑗,all𝑗𝐼,-matrix. If they have additional (nonstructural) rows being zero, the corresponding 𝐵𝑖-vectors can be neglected in (2.10). From (2.7), it is clear that the rows of the 𝑒()-matrices have the following structures: 𝐵i𝑏𝑇𝑗(𝑖)[(𝑖1)block0,,0,𝑒𝑇𝑗𝑙𝑖+1𝐵𝑖,(𝑝𝑖)blocks0,,0];if1𝑖𝑗<1𝑖+𝑛𝑖,[(𝑖1)blocks𝑒𝑇𝑗𝑙𝑖𝑛𝑖+1𝑀𝑖1𝐶1,0,,𝑒𝑇𝑗𝑙𝑖.𝑛𝑖+1𝑀𝑖(𝑖1)𝐶𝑖1,0,𝑒𝑇𝑗𝑙𝑖𝑛𝑖+1𝑃𝑖,(𝑝1)blocks𝑒𝑇𝑗𝑙𝑖𝑛𝑖+1𝑀𝑖(𝑖+1)𝐶𝑖+1,0,,𝑒𝑇𝑗𝑙𝑖𝑛𝑖+1𝑀𝑖𝑝𝐶𝑖,0];if𝑙𝑖+𝑛𝑖𝑖𝑗<1𝑖+𝑛𝑖+𝑎𝑖,(2.11) Then one deduces from (2.11) that forall𝑖𝐼.𝑏𝑇𝑗+1𝑖1(𝑖)𝑥𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘);if1𝑗𝑛𝑖,𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖+1𝑦𝑙(𝑘);if𝑛𝑖+1𝑗𝑛𝑖+𝑎𝑖,0;if𝑛𝑖+𝑎𝑖+1𝑗𝑛+𝑎,(2.12), when with𝑘=0,1,2,3,, 𝑏𝑇𝑗(𝑖), 𝑝𝑇𝑗(𝑖) denote, respectively, the 𝑚𝑇𝑗(𝑖,1) vector row of the matrices 𝑗th From (2.7) through (2.12) one obtains 𝐵𝑖,𝑃𝑖,and𝑀𝑖,forall𝐼,1𝐼. where the sequences of scalars 𝑥(𝑘+1)=𝐴𝑥(𝑘)+𝑝𝑖=1𝑛𝑖+𝑎𝑒𝑗=1𝑒𝑗+𝑙𝑖1𝜈𝑗𝑖(𝑘);𝑘=0,1,2,,(2.13) are defined by 𝜈𝑗𝑖(𝑘)𝜈𝑗𝑖(𝑘)𝑢𝑖(𝑘)𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘);if1𝑗𝑛𝑖,𝑢𝑖(𝑘)𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘);if𝑛𝑖+1𝑗𝑛𝑖+𝑎𝑖,(2.14) are considered as equivalent (primary) inputs for the equivalent linear systems (2.13) to be determined prior to the generation of the forall𝑖𝐼,𝑘=0,1,2,,-sequence in the design implementation.

2.3. Obtaining the 𝑢()()-Sequences from the 𝑢𝑖(𝑘)-Sequences

Equations (2.13) and (2.14) constitute a feed-forward linear equivalent system with multiplicative inputs. But the 𝛾𝑗𝑖(𝑘)𝑝-inputs cannot, in general, be determined from the 𝑢𝑖(𝑘) equivalent (𝑛+𝑎)>𝑝-inputs. Then, it is necessary to take a set of p variables between the last ones, which will be optimized independently. Let us call 𝜈𝑗𝑖(𝑘) Thus, let us define the sets

𝜈𝑗𝑖(𝑘),all𝑖𝐼. such that for the 𝑁𝑖=𝑁𝑖𝑈𝑁+𝑖,𝑁𝑖=1,2,,𝑛𝑖,𝑁+𝑖=𝑛𝑖+1,𝑛𝑖+2,,𝑛𝑖+𝑎𝑖,all𝑖𝐼,(2.15) current sampling instant, 𝑘th and 𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)0,if𝑥𝑖(𝑘)0,if𝑗𝑖𝑁𝑖,𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦𝑖(𝑘)0,if𝑗𝑖𝑁+𝑖,(2.16)

Then from (2.14), one has the auxiliary inputs 𝑘=0,1,2,. with the real inputs for the nonsaturated case being𝜈𝑗𝑖(𝑘)𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)𝑖𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘)𝜈𝑗𝑖(𝑘);𝑗𝑁𝑖;𝑖𝐼,𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)1×𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)𝜈𝑗𝑖(𝑘),𝑗𝑁+𝑖;𝑖𝐼;𝑘=0,1,2,3,,(2.17)

Remark 2.4. Note that the minimum energy control (i.e., 𝑢𝑖(𝑘)𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)1𝜈𝑗𝑖(𝑘);if𝑗𝑖𝑁𝑖,ifthereexistssome𝑗𝑖𝑁𝑖suchthat𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)0,𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦𝑖(𝑘)1𝜈𝑗𝑖(𝑘);if𝑗𝑖𝑁𝑖,ifthereexistssome𝑗𝑖𝑁+𝑖suchthat𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)0,0;if𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)=0forall𝑗𝑖𝑁𝑖andalso𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)=0forall𝑗𝑖𝑁+𝑖,all𝑖𝐼;𝑘=0,1,2,3,.(2.18) some 𝑢𝑖()=0,) has been chosen when the system becomes insensitive to the corresponding control component. In fact, any arbitrary control could be applied.

Remark 2.5. Since there are 𝑖𝐼 equivalent inputs, deals with an over determined problem when the (𝑛+𝑎) inputs are solved from them. Their solution cannot generally be satisfied exactly, and there exist many possible ways of defining the “best” approximate solution if there are 𝑝<𝑛+𝑎 equivalent nonzero inputs. Also, note that the decomposition of the 𝑝<𝑠𝑛+𝑎-matrices (here, into a sum of dyads, as shown in (2.10)) is not unique. Then, the least-squares technique (which is known to give the best approximation) can be considered by taking into account the overall control term in (2.13) instead of the control itself. Namely, 𝐵𝑖 where 𝑢𝑖(𝑘)𝜈𝑇𝑖(𝑘)𝐵𝑖𝑥(𝑘)𝑥𝑇(𝑘)𝐵𝑇𝑖𝐵𝑖𝑥(𝑘);if𝐵𝑖𝑥(𝑘)0,0,otherwise,all𝑖𝐼;𝑘=0,1,2,3,,(2.19)
Note that in general, the (𝜈𝑖(𝑘)[𝑖1𝑗=1(𝑛𝑗+𝑎𝑗)zeros0,,0,𝜈𝑙𝑖(𝑘),,𝜈𝑛𝑖,𝑖(𝑘),𝜈𝑛𝑖+1,𝑖(𝑘),,𝜈𝑛𝑖+𝑎𝑖,𝑖(𝑘),𝑝𝑗=𝑖+1(𝑛𝑗+𝑎𝑗)zeros0,,0]𝑇𝑅𝑛+𝑎,all𝑖𝐼;𝑘=0,1,2,3,.(2.20)) vector has only nonzero components. Equation (2.19) can be rewritten as 𝑘𝑢𝑖(𝑘)[𝑛𝑖+𝑎𝑖𝑗=1{𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘)𝛿𝑗(𝑖)+𝑝𝑇𝑗𝑛𝑖0(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)𝛿𝑗(𝑖)}2]1×[𝑛𝑖+𝑎𝑖𝑗=1𝜈𝑗𝑖(𝑘){𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘)𝛿𝑗(𝑖)+𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)𝛿𝑗(𝑖)}],(2.21)

2.4. Two Alternative Interpretation Schemes of the Equivalent Linear System

Equation (2.18) leads to the following equivalences to the linear system (2.13): where𝛿(𝑖)1𝛿𝑗(𝑖),with𝛿𝑗(𝑖)=1if1𝑗𝑛and0,if𝑛𝑖+1𝑗𝑛𝑖+𝑎𝑖all𝑖𝐼;𝑘=0,1,2,. where the 𝑥(𝑘+1)=𝐴𝑥(𝑘)+𝑝𝑖=1𝑐𝑗+𝑙𝑖1(𝑘)𝜈𝑗𝑖(𝑘)(2.22)=𝐴𝑥(𝑘)+𝑝𝑖=1𝑒𝑗+𝑙𝑖1𝜈𝑗𝑖(𝑘)+𝜉(𝑘),(2.23)-vectors (𝑛+𝑎) and 𝑐()() are defined in Appendix A.

Equations (2.22)-(2.23) provide two alternative interpretations of (2.13). Representation (2.22) involves an extended state-dependent control vector, while (2.23) includes a known deterministic input-dependent vector. These equivalent representations will lead to two different suboptimal control designs in the next section.

2.5. Controllability of the Bilinear System with Interconnections

The controllability outside the origin of the bilinear system (2.3)–(2.5) can be studied by splitting the problem as follows:

(i)controllability of the feed-forward linear system,(ii)capability of the 𝜉()-sequences of actually acting upon the equivalent 𝑢()()-sequences.The following result on controllability conditions is useful as a previous requirement for optimization.

Theorem 2.6. (controllability of the bilinear system involving interconnections). Controllability of the bilinear system (2.3)–(2.5) is ensured on an interval 𝜈()(), [𝑘,𝑘]𝑅+, if the following set of conditions holds.
(i)The controllability grammian of the overall system is positive definite in 𝑘=𝑘+𝑚𝑎𝑥1𝑖𝑝(𝑛𝑖+𝑎𝑖)(ii)[𝑘,𝑘]. with max[rank(𝑝𝑙=1𝑀𝑠1𝐶1)]max[rank(𝐿𝑖𝑗𝐶𝑗)](iii)All the pairs of block diagram, 1𝑠𝑝,1𝑖𝑝and1𝑗𝑝., are nonsingular.(iv)[(𝐴𝑇𝑖,𝑀𝑇𝑖),𝐵𝑖],𝑖𝐼, (𝑏𝑇𝑗𝑖(𝑖)𝐶𝑖(𝑘)0 diagram 𝑏𝑇𝑗𝑖(𝑖), (𝐴1𝑖,𝑀1𝑖)𝐶𝑖(𝑘)0, all𝑖𝐼, where all𝑘𝑘𝑘 are, respectively, the subvectors of the matrices 𝐶()and𝑏() in (2.3)-(2.4) associated with each subsystem 𝐶()and𝐵().

Sketch of the Proof
Since the pair 𝑆() is controllable for all arbitrary matrix [𝐴+𝐵()𝐺(),𝐵()] if the pair 𝐺() is controllable [29, 38], then the controllability of the pair [𝐴,𝐵()] implies the controllability of the pair [𝐴𝑑,𝐵()] if there exists a [𝐴𝑑+𝐴,𝐵()],-matrix such that 𝐺() with 𝐴=𝐵()𝐺(), 𝐴=𝐴𝐴𝑑  block diagram 𝐴𝑑=. For the equality above relating ((𝐴1,𝑀1),(𝐴2,𝑀2),,(𝐴𝑝,𝑀𝑝)) and 𝐴,𝐵(), to hold, the Rouché-Fröbenius theorem implies that condition (ii) must be fulfilled according to (2.9) and (2.11). Also, since controllability grammian of 𝐺() on [𝐴𝑑,𝐵()] block diagram controllability grammian [[block diagram [𝑘,𝑘]],[block diagram (𝐴1,𝑀1),𝐵1], (𝐴2,𝑀2),𝐵2, [block diagram , which is positive if proposition (i) holds. This proves sufficiency of (i)-(ii) versus controllability of the equivalent feed-forward system. Propositions (iii)-(iv) ensure from Theorem 2.1 that controllability of the equivalent linear system implies controllability of the system (2.3)-(2.4) by allowing the determination of the scalars (𝐴𝑝,𝑀𝑝),𝐵𝑝]] from 𝑢𝑖(), 𝜈𝑗𝑖(), in order to drive each subsystem 𝑖𝐼 from an arbitrary initial point to a predefined final state.

Remark 2.7. (2.5.1) If condition (iv) is such that the matrices 𝑆𝑖 and 𝐴1𝑖, some 𝑀1𝑖, the same extensions referenced in Theorem 2.1(iii) are applicable (namely, see Remark 2.2).
(2.5.2) Theorem 2.6(iv) ensures the existence of an admissible control sequence which leads to the system outside the hyperplane of 𝑖𝐼 of insensitivity of the subsystem 𝑏𝑇𝑙𝑖+𝑗𝑖1(𝑖)𝑥𝑖(𝑘)0 to the 𝑆𝑖-control [17, 28]. A weak additive signal can be added to the right-hand side of (2.3) at the 𝑢𝑖() sample in case of the failure of condition (iv) in Theorem 2.6 at the next step.

3. Optimization Techniques

In this section, a regulation criterion which involves quadratic weighting terms in the state variables is given [1, 39]. The loss function used is 𝑘th where “𝐽𝑁=12𝑁𝑘=0𝑥𝑇(𝑘)𝑄(𝑘)𝑥(𝑘);𝑄()=𝑄𝑇()0;𝑁<,(3.1)” denotes positive semidefiniteness. The reason of choosing a finite-time horizon N is then seen from a practical context viewpoint.

In the sequel, the various inverse matrices which take place are assumed to exist. Because of their structures, this hypothesis is not strong.

3.1. Centralized Control

The overall system comprises subsystems, which are interconnected, and the implementation design must deal with the interactions that exist. Centralized control techniques are used for optimization. In the centralized approach, contrarily to the decentralized case [19, 27, 40], a central coordinator exists to take into account the interactions; namely, the coordinator is supposed to be available to the local controller to supply information to each control station from the remaining ones.

The philosophy involved is to define transformed state and input variables, which allows the redefinition of the loss function including quadratic terms of the redefined state and input variables. The idea of using transformations of variables was first pointed out in [1] for the linear continuous case. Its major interest appears in the case of constrained input sequences because, despite this fact, a recursive Riccati expression can be found leading to the optimal feedback solution. Thus, three transformations of variables are made on the equivalent feed-forward linear system of (2.22) (subsequently (2.23) will be invoked) to redefine the state vector equation of standard type (see (B.2) in Appendix B); that is, 0 On the other hand, such transformations are applied to the loss function (3.1) to obtain a standard performance index on a given planning horizon 𝑓(𝑘+1)=𝐴(3)(𝑘)𝑓(𝑘)+𝐶(3)(𝑘)𝜈(3)(𝑘);𝑘=0,1,2,.(3.2)

The necessary manipulation on (2.22)–(3.1) to derive (3.2)-(3.3) is detailed in Appendix B.

Now, since (3.3) is a quadratic criterion of standard type (i.e., it involves quadratic terms in both the input and the state variables), a Riccati-matrix sequence may be found, which leads to the optimal solution. The main optimization results are now summarized.

Theorem 3.1. (optimization of the system equivalent inputs). Assume that Theorem 2.6 holds. Then, the optimal equivalent input sequence for the equivalent feed-forward linear system with respect to the loss function (3.1) is given by the following expression: 𝐽𝑁=12𝑁𝑘=0𝑓𝑇(𝑘)𝑄(3)(𝑘)𝑓(𝑘)+𝜈(3)𝑇(𝑘)𝑅(3)(𝑘)𝜈(3)(𝑘);𝑁<.(3.3) where the optimal equivalent redefined input sequence is 𝜈𝑗𝑖(𝑘)=𝜈(3)𝑗𝑖(𝑘+1)𝐶𝑇𝑖(𝑘)𝑄(𝑘+1)𝐶𝑖(𝑘)1𝐶𝑇𝑖(𝑘)𝑄(𝑘+1)𝐴𝑥(𝑘)𝑒𝑇𝑗𝑖𝑅(3)1(𝑘+1)𝑊(2)𝑇(𝑘+1)𝐴𝑥(𝑘);𝑖𝐼;𝑘=0,1,2,,𝑁,(3.4) with 𝜈(3)𝑗𝑖(𝑘+1)=𝑒𝑇𝑗𝑖𝑅(3)(𝑘+1)+𝐶(3)𝑇(𝑘+1)𝑃(𝑘+2)𝐶(3)(𝑘+1)1×𝐶(3)𝑇(𝑘+1)𝑃(𝑘+2)𝐴(3)(𝑘+1)𝐴𝑥(𝑘);𝑘=0,1,2,,𝑁(3.5) being the recursive Riccati matrices; that is, 𝑃()

Outline of the Proof
It is obvious by direct substitution while dealing with the aforementioned optimization strategies [1, 29, 38] applied to the feed-forward linear equivalent system of (2.22).

Corollary 3.2 (generation of the optimal input sequences). Under the same hypotheses as in Theorem 3.1, consider the sequences of scalars 𝑃(𝑘)=𝑄(3)(𝑘)+𝐴(3)(𝑘)𝑃(𝑘+1)×{𝐼𝐶(3)(𝑘)𝑅(3)(𝑘)+𝐶(3)(𝑘)𝑃(𝑘+1)𝐶(3)(𝑘)1𝐶(3)(𝑘)𝑃(𝑘+1)}𝐴(3)(𝑘);𝑃(𝑁)=𝑄(3)(𝑁);𝑘=0,1,,𝑁.(3.6) with 𝑢𝑖(𝑘)𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)1𝜈𝑗𝑖(𝑘);all𝑗𝑁𝑖,𝑖𝑓𝑗𝑖𝑁𝑖,𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)1𝜈𝑗𝑖(𝑘);all𝑗𝑁𝑖,𝑖𝑓𝑗𝑖𝑁+𝑖,0;𝑖𝑓𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)=0,𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)=0𝑗𝑁𝑖;𝑘=0,1,2,,𝑁(3.7) being the 𝜈𝑗𝑖(𝑘),𝑘=0,1,2, component of the optimal equivalent inputs obtained from Theorem 3.1.
Then, the following propositions hold.
(a) The optimal input sequence with respect to the loss function (3.1) is 𝑗𝑖th if the interval 𝑢𝑖(𝑘)=𝑢𝑖(𝑘),𝑘=0,1,2, is allowed as input definition domain, (,). (b)all𝑖𝐼  If|𝑢𝑖(𝑘)|𝑀<, then the optimal input sequence with respect to the loss function (3.1) is all𝑖𝐼,𝑢𝑖=sat[𝑢𝑖(𝑘)], where the saturation function is defined by all𝑖𝐼,𝑘=0,1,2,;

Proof. It follows from Theorem 3.1 and the fact that, from (3.7), propositions (a)-(b) in Corollary 3.2 imply that sat𝑔(𝑥)𝑀sgn𝑔(𝑥);𝑖𝑓|||𝑔(𝑥)|||>𝑀,𝑔(𝑥);𝑜𝑡𝑒𝑟𝑤𝑖𝑠𝑒.(3.8) belongs to the boundary of its definition domain (and vice versa), and the fact that the optimal Hamiltonian function associated with the loss function (3.1) is a strictly convex function of the equivalent input sequence 𝑢()().

Remark 3.3. (3.1.1) Note that the equivalent and the true input sequences are obtained in a centralized way because of the coupling between the different 𝜈()() through the states of the remaining 𝑆𝑖
(3.1.2) Also, note that (2.7) can be rewritten as𝑆𝑗,𝑖,𝑗𝐼. where 𝑥(𝑘+1)=𝐴𝑥(𝑘)+𝑝𝑗=1𝐻𝑗(𝑘)𝑢𝑗(𝑘);𝑘=0,1,2,,(3.9) By applying optimization techniques to (3.9), one could obtain, at first, an optimal 𝐻𝑖(𝑘)𝐵𝑖𝑥(𝑘);all𝑖𝐼;𝑘=0,1,2,.(3.10)-sequence. But this would be nonsense because the Riccati-matrix sequence would be dependent on the future measurements of the state vector.
(3.1.3) Due to the structure of the 𝑢()()-vectors which are state-dependent (see (A.1)) and so unknown in advance, both the optimal equivalent and the true inputs cannot be determined through Riccati-matrix sequence from the results in Theorem 3.1 and Corollary 3.2.
Suboptimal schemes are now presented which are alternative to the optimal one of Theorem 3.1.

3.2. Centralized Suboptimal Schemes
3.2.1. Centralized Suboptimal Modified Scheme 1

It consists of the following variations on the optimal scheme from Theorem 3.1 and Corollary 3.2:

(i)to use finite (sufficient small for the problem coherence) time-sliding optimization horizons. Namely, the time horizon for optimization is one-step advanced when the input is determined at each step. Only the first input associated with each optimization horizon is in fact applied;(ii)to approximate the 𝑐()-vectors by their values by using each current value at the first sample of each optimization horizon; that is, 𝐶()(), 𝐶(𝑡)=𝐶(𝑘)

3.2.2. Centralized Suboptimal Modified Scheme 2

It is based on the application of (2.23) and (A.2) in Appendix A, by implementing the optimization scheme while neglecting the dependence of all𝑘𝑡𝑘+𝑁;𝑘=0,1,2,. on the 𝜉()-optimal equivalent input sequence. This leads to a suboptimal Riccati-matrix sequence, which does not depend on the state vector, and thus, implementable. The presence of the deterministic disturbance 𝜈(3)()() obliges to use an additional vector (denoted by 𝜉()) in the optimization procedure in Theorem 3.1 in order to maintain a Riccati-type solution. The costate associated with the Hamiltonian of the loss function (3.1) verifies the so-called “modified Riccati transformation” [26]; namely, 𝜌() Thus, by taking into account that now 𝜂(𝑡)(costate)𝑃(𝑡)𝑓(𝑡)𝜌(𝑡). one finds out the following result instead of Theorem 3.1. See Appendix B for particular mathematical details.

Theorem 3.4. (optimization of the equivalent inputs for the system representation including a deterministic disturbance). Under the same assumptions as in Theorem 3.1, it follows that the optimal equivalent input sequence for the equivalent feed-forward linear representation involving a deterministic disturbance with respect to the loss function (3.1) is given by 𝑓(𝑘)𝐴𝑥(𝑘1)+𝜉(𝑘1), where the optimal equivalent redefined input sequence is 𝜈𝑗𝑖(𝑘)=𝜈(3)𝑗𝑖(𝑘+1)𝑒𝑇𝑗𝑖𝑄(𝑘+1)𝑒𝑗𝑖1𝑒𝑇𝑗𝑖𝑄(𝑘+1)𝑒𝑇𝑗𝑖𝑅(3)1(𝑘+1)𝑊(2)𝑇(𝑘+1)𝐴𝑥(𝑘)+𝜉(𝑘);𝑎𝑙𝑙𝑖𝐼;𝑘=0,1,2,,(3.11) The recursive Riccati-matrix sequence is now given by 𝜈(3)𝑗𝑖(𝑘)=𝑒𝑇𝑗𝑖𝑅(3)(𝑘)+𝐶(3)𝑇(𝑘)𝑃(𝑘+1)𝐶(3)(𝑘)1𝐶(3)𝑇(𝑘)×𝜌(𝑘+1)𝑃(𝑘+1)𝐴(3)(𝑘)𝐴𝑥(𝑘1)+𝜉(𝑘1)+𝜉(𝑘);𝑎𝑙𝑙𝑖𝐼;𝑘=1,2,3,.(3.12)

Proof. It follows immediately from applying the three transformations of variables of Appendix B to (2.23) and (A.2), and from the fact that 𝑃(𝑘)=𝑄(3)(𝑘)+𝐴(3)(𝑘)×{𝑃(𝑘+1)𝑃(𝑘+1)𝐶(3)(𝑘)𝑅(3)(𝑘)+𝐶(3)𝑇(𝑘)𝑃(𝑘+1)𝐶(3)(𝑘)1𝐶(3)(𝑘)𝑃(𝑘+1)}𝐴(3)(𝑘);𝑃(𝑁)=𝑄(3)(𝑁);𝑘=0,1,2,,𝑁,(3.13)𝜌(𝑘+1)={𝐼+𝑄(3)(𝑘)𝑃(𝑘)𝐴(3)1(𝑘)𝑅(3)1(𝑘)𝐶(3)(𝑘)𝐶(3)𝑇(𝑘)𝐴(3)𝑇(𝑘)}𝐴(3)𝑇(𝑘){𝜌(𝑘)+𝑃(𝑘)𝑄(3)(𝑘)𝐴(3)1(𝑘)𝜉(𝑘)};𝜌(𝑁)=0;𝑘=0,1,2,,𝑁.(3.14) or equivalently 𝑓(𝑘)𝐴𝑥(𝑘1)+𝜉(𝑘1), and the fact that the 𝑓(𝑘+1)=𝐴(3)(𝑘)𝑓(𝑘)+𝐶(3)(𝑘)𝜈(3)(𝑘)+𝜉(𝑘);𝑘=0,1,2,,(3.15) vectors are now 𝐶(3)()

Remark 3.5. (3.2.1) Note that Corollary 3.2 also applies here to obtain the optimal system input sequence under constraints (see (3.7) and (3.8)).
(3.2.2) The 𝐶(3)(𝑘)=𝐶(3)𝐴𝑒𝑗1+𝑙11,𝐴𝑒𝑗2+𝑙21,,𝐴𝑒𝑗𝑝+𝑙𝑝1with𝑘=0,1,2,3,.(3.16)-vector (related to the deterministic disturbance in the modified Riccati transformation is dependent on the 𝜌()-sequence. To solve such a dependence, the 𝜈(3)()-vector of (3.14) can be decomposed into 𝜌(𝑘+1) being 𝜌(𝑘+1)=𝑚(𝑘)𝜌(𝑘)+𝑀(𝑘)𝜉(𝑘);𝑘=0,1,2,,(3.17), 𝑀() matrices defined by𝑀()𝑅(𝑛+𝑎)𝑥(𝑛+𝑎) Thus, the following result yields.

Corollary 3.6 (useful implementability result for obtaining the optimal equivalent input). If the same assumptions from Theorem 3.4 are fulfilled, then the optimal equivalent input for the system interpretation including a deterministic disturbance with respect to the loss function (3.1) can be rewritten as follows (see (3.11)): 𝑀(𝑘)=𝐼+𝑄(3)(𝑘)𝑃(𝑘)𝐴(3)1(𝑘)𝑅(3)1(𝑘)𝐶(3)(𝑘)𝐶(3)𝑇(𝑘)𝐴(3)𝑇(𝑘)1;𝑘=0,1,2,3,,𝑀(𝑘)=𝐼+𝑄(3)(𝑘)𝑃(𝑘)𝐴(3)1(𝑘)𝑅(3)1(𝑘)𝐶(3)(𝑘)𝐶(3)𝑇(𝑘)𝐴(3)𝑇(𝑘)1×𝐴(3)𝑇(𝑘)𝑃(𝑘)𝑄(3)(𝑘)𝐴(3)1(𝑘);𝑘=0,1,2,3,.(3.18)where the matrix 𝜈𝑗𝑖(𝑘)={1+𝑒𝑇𝑗𝑖{𝑅(3)(𝑘+1)+𝐶(3)𝑇(𝑘+1)𝑃(𝑘+2)𝐶(3)(𝑘+1)1𝐶(3)𝑇(𝑘+1)𝑃(𝑘+2)𝐴(3)(𝑘+1)+𝑅(3)1(𝑘+1)𝑊(2)(𝑘+1)𝐹(𝑘+1)+𝑒𝑇𝑗𝑖𝑄(𝑘+1)𝑒𝑗𝑖1}𝑒𝑗𝑖𝑄(𝑘+1)𝐹(𝑘+1)}1×{𝑒𝑇𝑗𝑖𝑅(3)(𝑘+1)+𝐶(3)𝑇(𝑘+1)𝑃(𝑘+2)𝐶(3)(𝑘+1)1𝐶(3)𝑇(𝑘+1)×𝑀(𝑘+1)𝑃(𝑘+2)𝜉(𝑘+1)+𝑀(𝑘+1)𝜌(𝑘+1)𝑃(𝑘+2)𝐴(3)(𝑘+1)𝐴𝑥(𝑘)𝑒𝑇𝑗𝑖𝑄(𝑘+1)𝑒𝑗𝑖1𝑒𝑇𝑗𝑖𝑄(𝑘+1)𝐴𝑥(𝑘)𝑒𝑇𝑗𝑖𝑅(3)1(𝑘+1)𝑊(2)𝑇(𝑘)𝐴𝑥(𝑘)};all𝑖𝐼;𝑘=0,1,2,3,,(3.19) is partitioned as 𝐹()𝑅(𝑛+𝑎)×𝑝 So that each vector 𝐹(𝑘)𝐹1(𝑘),𝐹2(𝑘),,𝐹𝑝(𝑘);𝑘=0,1,2,3,.(3.20)𝐹𝑖(𝑘)𝑅𝑛+𝑎

Proof (outline). Equation (3.19) can be obtained after direct calculations if (3.17) and (3.18) are substituted into (3.11) while taking Remark 3.5 into account.

3.3. Considerations about Implementation

The following remarks must be pointed out.

(1)The system implementation scheme, which has been reported in the paper, has appeared to be suboptimal twice. First of all, the optimization scheme is dependent on the state/input vectors (some estimates must be made) according to the two alternative interpretations of the feed-forward linear system (2.13). This fact is due to the nature of the decomposition methods which have been applied and also due to the over determination problem one must deal with when the system inputs are generated from the associate equivalent ones. The performance degradation (i.e., the optimality losses being inherent to the applied suboptimization procedures) could be studied through direct calculations. However, the hypotheses that have been taken into account allow the system implementation.(2)In summary, the steps that the designer ought to follow in the implementation environment are as follows:(a)to set up a finite-time sliding optimization horizon, 𝐹𝑖(𝑘)𝑝𝑖=1𝑛𝑖+𝑎𝑖𝑗=1,𝑗𝑖{𝑏𝑇𝑗𝑖𝑛𝑖(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖}1×{𝑏𝑇𝑗𝑛𝑖(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖};all𝑖𝐼,𝑘=0,1,2,.(3.21), [𝑘,𝑘+𝑁], integers,(b)to test the system controllability according to Theorem 2.6 and/or Corollary3.2,(c)to estimate the state vector (suboptimal modified scheme 1 of Section 3.2.1) and the equivalent deterministic disturbance (suboptimal modified scheme 2 of Section 3.2.2) on the optimization horizon 1𝑁< To do that, the estimated system is chosen from (3.2) as [𝑘,𝑘+𝑁]. and the equivalent input estimates become ̂𝑓𝑘+1=𝐴(3)1,𝑘+1̂𝑓𝑘+𝐶(3)1,𝑘𝜈(3)𝑘1,̂𝑓𝑘+𝑗=𝐴(3)1,𝑘+𝑗̂𝑓𝑘+𝑗1+𝐶(3)𝑘+𝑗1𝜈(3)𝑘+𝑗1;𝑗2,𝑁,(3.22)(d)to compute the Riccati matrix from (3.6) or (3.13), respectively, and the system matrices from (B.6) and (B.8) of Appendix B, or (B.8) and (3.16), respectively.(e)to generate the suboptimal input sequence ̂𝜈𝑘(3)𝑘=̂𝜈𝑘1(3)𝑘1=𝜈(3)𝑘1,(3.23)̂𝜈𝑘+𝑗(3)𝑘=̂𝜈𝑘+𝑗(3)𝑘1;allinteger𝑘𝑗,𝑗1,𝑁,(3.24)̂𝜈𝑘+𝑗(3)𝑘=0;if𝑘=𝑗,(3.25) from (3.4) or (3.11), respectively.(f)to obtain the suboptimal input sequence 𝜈𝑗𝑖(𝑘),alli𝐼,𝑘=0,1,2, by applying Corollary 3.2. One must have only input controls as inputs to the system; namely, p inputs.

Remark 3.7. (3.3.1) Note that in (3.22)–(3.25), the following system representation is assumed to be 𝑢𝑖(𝑘),alli𝐼,𝑘=0,1,2,3,,𝑁, instead of (3.2), and the subscript “𝑓𝑘+1=𝐴(3)𝑘𝑓𝑘+𝐶(3)𝑘𝜈(3)𝑘;𝑘=0,1,2,(3.26)” in the matrices is related to the suboptimal modified scheme 1. Similar expressions can be derived for the alternative interpretation scheme (3.15) of Section 3.2.2.
(3.3.2) Also, note that such “fictitious” states must be estimated in order to later reupdate the input sequence. In fact, only the 1 equivalent input is applied in each optimization horizon.
(3.3.3) The usual asymptotic stability tests via Lyapunov's theorem (which is applied to linear time-invariant optimal systems with respect to quadratic criteria) are not applicable here because the system under study is neither invariant nor linear. Besides, standard stability proofs for optimal regulators cannot be applied to the optimal scheme, which has been obtained because of its nonimplementability. The suboptimal modified schemes can lead to instability (or at least stability can not be proved) when 𝜈(3)(𝑘1) . But this fact does not affect the system stability in the current context of this paper since a finite-time sliding for the optimization horizon is used.
In particular, these circumstances make the implementation of a decentralized design quite difficult because some links of the equivalent feed-forward system (2.13) are cut off. However, stability can be ensured if 𝑁 is a Hurwitz's matrix and 𝐴 is a bounded sequence, if a saturation type rule is applied when computing the optimal redefined equivalent inputs 𝑢() in the various schemes.

4. Numerical Simulation

4.1. Example 1

This first example implemented deals with the following ninth-order discrete model: 𝜈(3)()()𝑆1𝑥1(𝑘+1)𝐴1(𝑇)+𝑢1(𝑘)𝐵1(𝑇)𝑥1(𝑘)+𝐷1(𝑇)𝑧1(𝑘);𝑢1𝑅,𝑦1(𝑘)𝑥1(𝑘);𝑥1()𝑅2𝑥1,𝐻1𝑧1(𝑘+1)𝑀1(𝑇)+𝑢1(𝑘)𝑃1(𝑇)𝑧1(𝑘)+𝐿11(𝑇)𝑥1(𝑘);𝑧1()𝑅2𝑥1,𝐻2𝑧2(𝑘+1)𝑀2(𝑇)𝑧2(𝑘)+𝐿21(𝑇)𝑥1(𝑘);𝑧2()𝑅3𝑥1,𝐻3𝑧3(𝑘+1)𝑀3(𝑇)𝑧3(𝑘)+𝐿31(𝑇)𝑥1(𝑘);𝑧3()𝑅2𝑥1,(4.1) Table 1 displays the entries of the constant discrete matrices related to the composite homogeneous bilinear structures all𝑘0,𝑘𝑧., 𝑆1, 𝐻𝑖, in (4.1). The simulations have been performed with the centralized suboptimal modified scheme 1, including free control type (namely, without saturation) only according to Corollary 3.2.

Table 1: Matrix entries of the simulated system structure.

The steady state of the signal trajectory under no control leads to the following parameters: 𝐿2100(6.9022+44.5212𝑇)𝑇0, 𝑦1,st=21.8, 𝑦2,st=238.63,  𝑢st=400 samples. The system performance is studied as a function of the optimization horizon, 𝑁st=57 Also, the plant settings are: Initial conditions, 𝑁., 𝑥1(0)=[20.0,250.0]𝑇, 𝑧1(0)=[0.0,0.0]𝑇, and 𝑧2(0)=[0.0,0.0,0.0]𝑇; sampling period 𝑧3(0)=[0.0,0.0]𝑇 seconds and working horizon 𝑇=0.03 samples.

The obtained results are shown in Table 2 where the deviations (percent) are related to the normal operation signals (i.e., WH=150, 𝑦1,st and 𝑦2,st) and the CCT is the average time the computer needs to perform a computation iteration of the optimization horizon.

Table 2: Typical results of the discretely controlled plant without saturation as a function of the optimization horizon.

Exhaustive simulations about variations of the optimization horizon have shown that the discrete control improves the system response from a better transient characteristic point of view. However, there exist deviations of the steady-state signal tracking related to the no-control operation mode signals. Also, the experiments become time-consuming and memory-storage expensive as the optimization horizon does increase. For instance, some of the simulated examples have been plotted in Figure 2, consisting of the two components of the plant output and the waveform of the control input to the system (see Table 2 for the definitions of the involved signals related to the optimization horizon sizes).

Figure 2
4.2. Example 2

A discrete bilinear system which corresponds to the discretization of a bilinear continuous one for small sampling period and zero-order hold (ZOH) is given as defined in Section 2.2, by where 𝑥(𝑘+1)=𝐴(𝑇)+𝑢(𝑘)𝑁(𝑇)𝑥(𝑘)+𝐸(𝑇)𝑢(𝑘)+𝑤(𝑘),𝑦(𝑘)=𝐶𝑥(𝑘),(4.2), 𝑥(𝑘)[𝑥𝑇1(𝑘),𝑥𝑇2(𝑘),𝑥𝑇3(𝑘)]𝑇 and there exists a deterministic disturbance 𝑦(𝑘)[𝑦𝑇1(𝑘),0,0]𝑇 so that 𝑤(𝑘)[𝑤𝑇1(𝑘),0,0]𝑇 being 𝑤1(𝑘)𝑆1(𝑇)𝑢2(𝑘)+[𝑆2(𝑇)+𝑢(𝑘)𝑆3(𝑘)]𝑞(𝑘), the sampling period, 𝑇 denotes transpose and the sampling instant is defined as ()𝑇

In (4.2), the coefficients of the block matrices are composed of the following structures: 𝑡𝑘+1=𝑡𝑘+𝑇. In particular, (4.2) and (4.3) correspond to the linearized discretization of an electromechanical system consisting of two amplidynes acting on the DC-motor shown in Figure 3, with a regulation mechanism incorporated as introduced in Section 3. This regulation scheme is not restricted to considering either the relevant control current or voltage to be constant, which emphasizes the interest of the bilinear modelling.

Figure 3: Electromechanical system with interconnected subsystems.

Figure 4 plots a typical tracking response output of the system as a function of the sampling period, for a fixed optimization horizon 𝐴(𝑇)=𝐴1(𝑇)0𝑀1(𝑇)𝑃2(𝑇)𝐴2(𝑇)𝑅2(𝑇)𝑃3(𝑇)𝑅3(𝑇)𝐴3(𝑇);𝑁(𝑇)=𝑁1(𝑇)0𝑆2(𝑇)000000;𝐸(𝑇)=𝐸1(𝑇)𝐸2(𝑇)𝐸3(𝑇);𝐶=𝐼00.(4.3) samples and 𝑁=50 being unity in the quadratic regulation performance criterion (see (3.1)).

Figure 4: Motor angular speed output as function of the sampling period with respect to a step reference signal of 277 rad/seg.

The control action implemented can be used as an alternative to the traditional strategies of linear control on electrical DC-machines. In general, it has been noticed that the discrete control implies better transients related to the linear approach because those transient responses are faster, that is, present lower settling times with greater shooting parameters. Besides, it has also been observed that the discrete control performance improves, in general, as the optimization horizon is increased for different experiments, with respect to that of the nominal plant when no control action is reached.

5. Conclusions

A multivariable invariant discrete-time bilinear system being composed of interconnected subsystems has been studied. An equivalent feed-forward linear system with equivalent inputs, which are derived from products state-input, has been given. Then, the system has been suboptimized with respect to a quadratic finite-time optimization horizon in order to drive each subsystem from any arbitrary initial point to a predefined final state.

The suboptimization has been made by neglecting either the time-dependence of the control vectors on the state vector (modified suboptimal scheme 1) or the dependence of a deterministic disturbance vector on the equivalent input sequence (modified suboptimal scheme 2). These approaches have effect only on the system implementability rather than on its stability. Besides, stability proofs lead to drawbacks when the optimization horizon is infinite because of the suboptimal real implementation. This is also translated into drawbacks when implementing decentralized expected results. The proposed suboptimal schemes have been proven by means of realistic examples.


A. Vector Variation for the Two Alternative Interpretation Schemes of the Feed-Forward Linear Systems

Section 2.4 deals with two alternative representations of the feed-forward linear system (2.13). This appendix is devoted to the derivation of the vectors 𝑄 and 𝐶()() involved in this section.

Substituting (2.18) into (2.13) and grouping terms, one obtains

𝜉() or alternatively 𝑐𝑗+𝑙𝑖1(𝑘)𝑒𝑗𝑖+𝑙𝑖1+𝑛𝑖+𝑎𝑖𝑗=1/𝑗𝑗𝑖{{𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖}1×{𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖}𝑒𝑗+𝑙𝑖1},(A.1) where 𝜉(𝑘)𝑝𝑖=1𝑛𝑖+𝑎𝑖𝑗=1/𝑗𝑗𝑖{𝑏𝑇𝑗𝑖(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑖𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑖𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖}1×{𝑏𝑇𝑗(𝑖)𝑥𝑖(𝑘)𝛿𝑗𝑖+[𝑝𝑇𝑗𝑛𝑖(𝑖)𝑧𝑖(𝑘)+𝑝𝑙=1𝑚𝑇𝑗𝑛𝑖(𝑖,𝑙)𝑦1(𝑘)]𝛿𝑗𝑖}𝜈𝑗𝑖(𝑘)𝑒𝑗+𝑙𝑖1};ifatleastanadmissible𝑗𝑖𝑁𝑖existsatthecurrentsamplinginstant,0;otherwise,(A.2)

B. Derivation of the Optimization Equations

In the sequel, in order not to repeat tedious notation, the 𝛿𝑗𝑖1𝛿𝑗𝑖,with𝛿𝑗𝑖=1,if𝑗𝑖𝑁𝑖,and𝛿𝑗𝑖=0,if𝑗𝑖𝑁+𝑖,all𝑖𝐼,𝑘=0,1,2,.-index (related to the system inputs and associate equations of the feed-forward linear system 2.13)) will be denoted by the subscript 𝑗𝑖 Also, three modified auxiliary inputs 𝑖,𝑖𝐼., 𝜈(1)(𝑘) and 𝜈(2)(𝑘), are calculated from 𝜈(3)(𝑘) and introduced due to the bilinear terms.

B.1. First Transformation

Let us define new variables as 𝜈(𝑘), with 𝐶(1)𝑗(𝑘)𝐴𝐶𝑗(𝑘1);𝜈(1)𝑗(𝑘)𝜈𝑗(𝑘1),𝑓(𝑘)𝐴𝑥(𝑘1)=𝑥(𝑘)𝑝𝑖=1𝐶𝑖(𝑘1)𝜈𝑖(𝑘1)(B.1) and 𝑗𝐼.

Substituting (B.1) into (2.22), one obtains 𝑘1 Also, the loss function (3.1) can be equivalently rewritten as 𝑓(𝑘+1)=𝐴𝑓(𝑘)+𝑝𝑖=1𝐶(1)𝑖(𝑘)𝜈(1)𝑖(𝑘);𝑘=0,1,2,3,.(B.2) where 𝐽𝑁=12{𝑁𝑘=0𝑓𝑇(𝑘)𝑄(𝑘)𝑝𝑖=1𝑟1𝑖(𝑘)𝑖(𝑘)𝑇𝑖(𝑘)𝑓(𝑘)+𝑝𝑖=1{𝑟𝑖(𝑘)𝜈(1)𝑖(𝑘)+𝑟1𝑖(𝑘)𝑇𝑖(𝑘)𝑓(𝑘)2+2𝑝𝑗=1/𝑗>𝑖𝐶(1)𝑇𝑖(𝑘)𝐴𝑇𝑗(𝑘)𝜈(1)𝑖(𝑘)𝜈(1)𝑗(𝑘)}};𝑁<,(B.3)𝑟𝑗(𝑘)𝐶(1)𝑇𝑗(𝑘)𝐴𝑇𝑄(𝑘)𝐴1𝐶(1)𝑗(𝑘)=𝐶𝑇𝑗(𝑘1)𝑄(𝑘)𝐶𝑗(𝑘1),𝑗(𝑘)𝑄(𝑘)𝐴1𝐶(1)𝑗(𝑘)=𝑄(𝑘)𝐶𝑗(𝑘1),(B.4).

B.2. Second Transformation

The following new variables are defined as follows: all𝑗𝐼;𝑘=1,2,3,𝐶(2)𝑗(𝑘)𝐶(1)𝑗(𝑘),𝜈(2)𝑗(𝑘)𝜈(1)𝑗(𝑘)+𝑟1𝑗(𝑘)𝑇𝑗(𝑘)𝑓(𝑘)=𝜈𝑗(𝑘1)+𝑥𝑇(𝑘1)𝐴𝑇𝑄(𝑘)𝐶𝑗(𝑘1)𝐶𝑇𝑗(𝑘1)𝑄(𝑘)𝐶𝑗(𝑘1)1,𝐴(2)(𝑘)𝐴𝑝𝑖=1𝑟1𝑖(𝑘)𝐶(1)𝑖(𝑘)𝑇𝑖(𝑘)=𝐴𝑝𝑖=1𝐶𝑇𝑖(𝑘1)𝑄(𝑘)𝐶𝑖(𝑘1)1𝐴𝐶𝑖(𝑘1)𝐶𝑇𝑖(𝑘1)𝑄(𝑘),(B.5).

From (B.5), one has equivalently to (B.3);all𝑗𝐼;𝑘=1,2,3, where 𝑓(𝑘+1)=𝐴(2)(𝑘)𝑓(𝑘)+𝐶(2)(𝑘)𝜈(2)(𝑘);𝑘=0,1,2,3,,(B.6) Also, taking into account the mentioned (B.5), (B.3) becomes as follows: 𝐶(2)()[𝐶(2)1(),𝐶(2)2(),,𝐶(2)𝑝()]𝑅(𝑛+𝑎)𝑥(𝑛+𝑎),𝜈(2)()[𝜈(2)1(),𝜈(2)2(),,𝜈(2)𝑝()]𝑅𝑝.(B.7) with 𝐽𝑁=12𝑁𝑘=0{𝑓𝑇(𝑘)𝑄(2)(𝑘)𝑓(𝑘)+𝑝𝑖=1𝑟𝑖(𝑘)𝜈(2)2𝑖(𝑘)+2𝑝𝑗=1/𝑗>𝑖𝑝𝑖𝑗(𝑘)𝜈2𝑖(𝑘)𝜈2𝑗(𝑘)2𝑝𝑖=1𝑝𝑗=1/𝑗>𝑖𝑓𝑇(𝑘)𝑝𝑖𝑗(𝑘)𝑡𝑇𝑗(𝑘)𝜈2𝑖(𝑘)+𝑡𝑇𝑖(𝑘)𝜈2𝑗(𝑘)};𝑁<(B.8)𝑝𝑖𝑗(𝑘)=𝑝𝑗𝑖(𝑘)𝐶(2)𝑇𝑖𝐴𝑇𝑗(𝑘)=𝐶𝑇𝑖(𝑘1)𝑄(𝑘)𝐶𝑗(𝑘1);𝑗𝑖,𝑡𝑖(𝑘)𝑟1𝑖(𝑘)𝑇𝑖(𝑘)=𝐶𝑇𝑖(𝑘1)𝑄(𝑘)𝐶𝑖(𝑘1)1𝐶𝑇𝑖(𝑘1)𝑄(𝑘),(B.9)

Equation (B.8) can be compactly rewritten as forall𝑖,𝑗𝐼;𝑘=0,1,2,. where the matrices involved are 𝐽𝑛=12𝑁𝑘=0𝑓𝑇(𝑘)𝑄(2)(𝑘)𝑓(𝑘)+𝜈(2)𝑇(𝑘)𝑅(2)(𝑘)𝜈(2)(𝑘)+2𝑓𝑇(𝑘)𝑊(2)(𝑘)𝜈(2)(𝑘);𝑁<,(B.10) with 𝑄(2)(𝑘)𝑄(𝑘)𝑝𝑖=1{𝐶𝑇𝑖(𝑘1)𝑄(𝑘)𝐶𝑖(𝑘1)1𝑄(𝑘)𝐶𝑖(𝑘1)𝐶𝑇𝑖(𝑘1)𝑄(𝑘)2𝑝𝑗=1/𝑗>𝑖{𝐶𝑇𝑗(𝑘1)𝑄(𝑘)𝐶𝑖(𝑘1)1𝐶𝑇𝑗(𝑘1)𝑄(𝑘)𝐶𝑗(𝑘1)1×𝐶𝑖(𝑘1)𝑄(𝑘)𝐶𝑗(𝑘1)𝑄(𝑘)𝐶𝑖(𝑘1)𝐶𝑇𝑗(𝑘1)𝑄(𝑘)}}(B.11)𝑅(2)(𝑘)=𝑅(2)𝑖𝑗(𝑘);𝑅(2)𝑖𝑗(𝑘)=𝑅(2)𝑗𝑖(𝑘)𝑟𝑖(𝑘),if𝑗=𝑖,𝑝𝑖𝑗(𝑘),if𝑗𝑖,(B.12)𝑊(2)(𝑘)𝑊1(𝑘)𝑊2(𝑘)(B.13)𝑊1(𝑘)𝑡𝑇1(𝑘),𝑡𝑇2(𝑘),,𝑡𝑇𝑝(𝑘),𝑊2(𝑘)𝑊𝑖𝑗(𝑘);𝑊𝑖𝑗(𝑘)=𝑊𝑗𝑖(𝑘)0if𝑗=𝑖𝑝𝑖𝑗(𝑘),if𝑗𝑖(B.14), and 𝑟𝑖(),𝑝𝑖𝑗(), in (B.12) to (B.14); 𝑡𝑖();𝑗𝑖, are given, respectively, in (B.4) and (B.9).

Note that the factorization of matrix 𝑖,𝑗𝐼in (B.13) is possible from (B.8).

B.3. Third Transformation

Although it is not necessary, this transformation of variables provides immediately a loss function (redefinition of (3.1)), which includes weighting terms associated with both the state vector and the transformed input. Let us define them as follows: 𝑊(2)() Substitution of (B.15) and (B.16) into (B.6) yields directly (3.2). Also, substituting (B.15) into (B.10), one obtains (3.3) if the weighting matrices are defined by 𝐶(3)(𝑘)𝐶(2)(𝑘),𝜈(3)(𝑘)𝜈(2)(𝑘)+𝑅(2)1(𝑘)𝑊(2)𝑇(𝑘)𝑓(𝑘),(B.15)𝐴(3)(𝑘)𝐴(2)(𝑘)𝐶(2)(𝑘)𝑅(2)1(𝑘)𝑊(2)𝑇(𝑘)𝑘=0,1,2,3,.(B.16) with 𝑄(3)(𝑘)𝑄(2)(𝑘)𝑊(2)(𝑘)𝑅(2)1(𝑘)𝑊(2)𝑇(𝑘),𝑅(3)(𝑘)𝑅(2)(𝑘)(B.17) and 𝑄(3)()0


The authors are very grateful to the Spanish MEC for supporting this work through the research Projects DPI2006-00714 and DPI2006-01677. They are also grateful to the UPV/EHU and the Basque Government MEC for its partial support through Projects EHU 06/88 and S-PE06UN10, respectively. They would like also to thank the comments and suggestions of the anonymous reviewers who have helped to improve the previous version of this manuscript.


  1. B. D. O. Anderson and J. B. Moore, Linear Optimal Control, Prentice-Hall, Englewood Cliffs, NJ, USA, 1971. View at Zentralblatt MATH · View at MathSciNet
  2. A. A. Feldbaum, Optimal Control Systems, vol. 22 of Mathematics in Science and Engineering, Academic Press, New York, NY, USA, 1965. View at MathSciNet
  3. T. J. Tarn, “Singular control of bilinear discrete systems,” Information and Computation, vol. 21, pp. 211–234, 1972. View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  4. T. J. Tarn, S. K. Rao, and J. Zaborszky, “Singular control of linear-discrete systems,” IEEE Transactions on Automatic Control, vol. 16, pp. 401–410, 1971. View at Publisher · View at Google Scholar · View at MathSciNet
  5. S. A. Al-Baiyat, “Model reduction of bilinear systems described by input-output difference equation,” International Journal of Systems Science, vol. 35, no. 9, pp. 503–510, 2004. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  6. Ü Kotta, S. Nõmm, and A. Zinober, “Classical state space realizability of input-output bilinear models,” International Journal of Control, vol. 76, no. 12, pp. 1224–1232, 2003. View at Publisher · View at Google Scholar · View at Zentralblatt MATH · View at MathSciNet
  7. A. J. Garrido, O. Barambones, I. Garrido, P. Alkorta, and F. Artaza, “Linear models for plasma current control,” in Proceedings of the 10th International Conference on Intelligent Systems and Control, pp. 240–245, Cambridge, Mass, USA. View at Google Scholar
  8. O. Barambones, A. J. Garrido, and F. J. Maseda, “Integral sliding-mode controller for induction motor based on field-oriented control theory,” IET Control Theory and Applications, vol. 1, no. 3, pp. 786–794, 2007. View at Publisher · View at Google Scholar