New Trends on Modeling, Design, and Control of Chaotic SystemsView this Special Issue
Sliding Mode Control of Discrete Chaotic System Based on Multimodal Function Series Coupling
A new sliding mode control model of discrete chaotic systems based on multimodal function series coupling is proposed to overcome the shortcomings of the standard PSO algorithm in multimodal function optimization. Firstly, a series coupled PSO algorithm (PP algorithm) based on multimodal function is constructed, which is optimized by multipeak solution on the basis of the standard PSO algorithm. Secondly, the improved PSO algorithm is applied to search all the extreme points in the feasible domain. Thirdly, the Powell method is used to perform the local optimization of the search results, and the newly generated extreme points are added to the extreme point database according to the same peak judgment operator. Finally, the long training time of PP algorithm can be overcome by the characteristics of fast convergence rate of the cloud mutation model. And also, both the population size and the redundancy can be reduced. Then, the clonal selection algorithm is used to keep the diversity of the population effectively. Simulation results of the sliding mode control of discrete chaotic systems show that the improved PSO algorithm obviously improves the response speed, overshoot, and so on.
Chaotic system can be divided into integer order chaotic system and fractional chaotic system. Fractional chaotic system is a very complex nonlinear system . Sliding mode control has become an important branch of modern control theory because of its good robustness and antijamming characteristics . Fractional chaotic system sliding mode has advantages of fractional chaotic system and sliding mode control, which can further improve the control performance of the system on the basis of traditional sliding mode control and become an important research area of modern nonlinear control.
Before introducing the sliding mode control of fractional order chaotic system, we mainly studied the control and synchronization of the fractional chaotic systems. At present, the control and synchronization of the fractional chaotic system are divided into two categories. One is the classical system stability theory and control method based on traditional integer order system. On the basis of the theory and control method, the stability analysis of the fractional order is carried out to design the integer order controller . The other is to utilize the characteristics of the fractional order system to put forward the fractional order controller design method; then, the control of the fractional order system can be realized . The control methods can be divided into the following areas:(1)Feedback control based on pole placement of the fractional order system: Matignon  discussed the relationship between internal stability and external stability of the fractional differential equations. The system is linearized near the equilibrium point. And the author also described the characteristic equation of the system and the eigenvector based on the fractional calculus theory and the traditional algebraic polynomial, analyzed the relationship between the fractional system eigenvalues and the system matrix pole, and presented the sufficient conditions for the stability of the fractional order system.(2)The control of Gronwall inequality based on the fractional order: the Gronwall inequality is extensively applied to the integer order differential equation, and the solution of the system is expressed analytically. The Gronwall inequality is taken to analyze the system solution to obtain the design of the controller, which is a common control design method for integer order systems. Lazarević and Spasić  applied the Gronwall inequality to the fractional time-delay system, analyzed the solution of the fractional order system, and realized the finite-time stabilization control of the fractional order. N’Doye and others  generalized the Bellman-Gronwall inequality, analyzed a class of nonlinear affine fractional order systems, and achieved the asymptotic stability of the fractional systems. Ye et al.  extended the integer order Gronwall inequality to the fractional systems, proposed the generalized Gronwall inequality, and successfully implemented the asymptotic stability control for the fractional time-delay systems.(3)Fractional order PID control: in the classical control theory and applications, PID controller has been developed very well and is widely applied to the actual system. In accordance with the PID control, the integral controller and differential controller of the PID controller can be partially extended to the fractional order; that is, the fractional integral and the fractional differential are introduced into the control design. Vasundhara Devi et al.  first proposed the fractional order PID controller and proved that the integer order PID controller belongs to a subset of the fractional PID controller. By adjusting the fractional order, the control effect of the fractional integral and fractional differential controller can be changed.(4)Feedback control based on linear matrix inequality: Sabatier and others  gave a sufficient and necessary condition for the stability of the fractional order time-invariant interval uncertain systems. Because of the state boundedness characteristic of chaotic systems, this condition has been applied to the control and synchronization of the fractional order chaotic systems. Chen and others  presented a sufficient and necessary condition for the fractional order systems with parametric perturbations and applied it to the fractional order chaotic systems to achieve robust synchronization of the fractional order chaotic systems. Dadras and Momeni  introduced the passive theory into the fractional order system and designed the fractional integral sliding surface by the linear matrix inequality to realize the passive sliding mode control of the time-varying uncertain fractional chaotic systems.(5)Sliding mode control of the fractional order system: the application of the integer order sliding mode to the fractional systems has been proved to be effective. The integral sliding mode surface and approach law are designed, respectively. And the effectiveness of the sliding surface is proved. The Lyapunov function is designed to analyze the stability of the fractional order system, and then, the sliding mode control law is obtained . On the other hand, by combining the design idea of the integral sliding mode surface with the characteristics of the fractional system, we can design the fractional sliding mode surface, and the integer integral of the sliding surface is changed to the fractional integral . Similarly, the approximation of the integer order can be extended to the fractional systems; that is, the differential equation of integer order approximation law description can be generalized into the fractional order equation, which can realize the fractional order sliding mode approximation law .
In this paper, a PSO multiobjective optimization algorithm is introduced. Since the algorithm is easy to fall into the local optimization, a sliding mode control of discrete chaotic system based on multimodal function series coupling is proposed combining the characteristics of multimodal function optimization problems. Simulations show the effectiveness of this method.
2. Multimodal Function Optimization Problem Modeling
In many engineering optimization, such as complex system parameters and structural optimization, neural network weights, and structural optimization, we not only need to find the global optimal solution in the feasible region but also need to search several global optimal solutions and other valued local optimal solutions. Thus, multichoices and multi-information can be provided for decision makers. The above can be classified as multimodal function optimization problem or multipeak function optimization (MFO).
The multimodal optimization problem, as shown in (1), consists of two parts: the objective function and the search range.
In (1), the search range that can meet all the variables is called the feasible region, and the solution in feasible region is called the feasible solution. In the feasible region, the solution with the smallest objective function is the optimal solution. The optimal solution can be divided into global optimal solution and local optimal solution. In the multimodal function optimization, not only all the global optimal solutions of the optimization problem need to be found, but also the local optimal solutions need to be found as many as possible.
For the MFO problem of the objective function minimum value solution, the fitness function is defined as
In (2), is the function value of , and is the fitness value.
For the MFO problem of the function maximum value solution, the limit construction method is used to solve the minimization problem. The fitness function is defined as
In (3), is the maximum estimate of .
3. Solution Algorithm of Multimodal Function Problem Based on Series Coupling
Based on PSO algorithm and Powell method, this paper presents PP algorithm based on series coupling. The PP algorithm has the maximum repetition search () to avoid the situations in which not all the extreme values can be searched in an evolutionary iteration. In each iteration, firstly, the improved PSO algorithm is used to carry out the global search of all the extreme points in the feasible region. Secondly, Powell method is used to perform the local search towards the optimal solutions found by the improved PSO to improve the accuracy of the solution. Finally, the peak judgement operator adds the newly generated extreme points to the extreme point database and rerandomizes the initialization of the particles.
3.1. PSO Algorithm Improvement
In the case of the particle swarm with the numbers of particles in the D-dimensional search space, denote the flying speed of particle as . The current position is , the searched optimal position of particle is , and the searched optimal position of the swarm is . In order to reduce the possibility of the particles leaving the search space during the process of evolutionary iterations, the speed is limited as . The standard PSO is shown in
In (4), is the current iteration number, represents the th dimensional component of the particle velocity vector , and . represents the cognitive learning factor of the particle individual and represents the social learning factor. is the random number generated by the uniform distribution on .
In the multimodal optimization, in order to make each local optimal value become the “only way” of some particles in PSO algorithm, this paper improves the standard PSO and obtains the improved PSO.(1)Ideally, the self-cognition of each individual in the particle swarm represents each local optimization of the multipeak optimization. Thus, let ; we can obtain not just the global optimal values, but all the local optimal values.(2)Consider that the random function in (4) mainly increases the randomness of the particle motion so as to make the local optimal point jump out. In order to let each particle of the particle swarm converge to each local optimal value as soon as possible, in (4) should be removed.(3)The linear decreasing strategy of inertia weight is simple and intuitive and has better searching ability, which makes the algorithm balance between the global search and the local search. Therefore, uses the above strategy, which is shown in
In (5), represents the current PSO iteration algebra, represents PSO maximum iteration number, and and represent the maximum and the minimum value of .
Based on the above improvements, (4) can be changed to
The position change is performed according to (7). In order to avoid the particles leaving from the search space during the process of evolution iteration, the following condition should be satisfied: .
3.2. Nonlinear Powell Direct Search Method
The nonlinear Powell direct search method is a nonlinear direct local search method for solving unconstrained optimization problems without using derivatives. It is considered to be a relatively effective method in direct search. Powell search method is known as the direction acceleration method as well. Using this method, the calculation of the derivative is not required and only the function value is used to carry out the one-dimensional search from one point towards two directions. Then, the minimal value can be obtained. Nonlinearity means that the search direction of the initial point is not fixed, the acceleration direction will change with the change of the initial point position, and the search is not along a straight line.
The steps of the Powell search method are as follows:
Step 1. The particle of the last search generation of PSO is taken as the initial point . Set the search precision of the Powell method as . numbers of initial linearly independent search directions are given, and . Set .
Step 2. According to (9) below, one-dimensional search is performed accurately from point to each search direction in turn . Thus, is obtained.In (8), and represent the step sizes and is obtained by the exact linear search; that is, is the solution of the one-dimensional optimization problem above. can be negative, which indicates that the exact linear search takes place over the entire real axis.
Step 3. Let . If , and the solution is obtained, calculation ends. Otherwise, the exact linear search is performed from point to the direction of to get .
Step 4. According to (10), the indicator determined by the maximal decrease is calculated.
Step 5. When (11) holds, it shows that are still linearly independent, and the search directions of the next round are still . Then, . Return to Step2.
Step 6. When (11) does not hold, it shows that are linearly dependent. Thus, set to ensure linearly independence of the new generated search direction, , . Return to Step2.
3.3. The Same Peak Judgement of the Extreme Points
The extreme points found by the PP algorithm are compared with all the extreme points in the extreme point database in turn and make the same peak judgement. If it is judged as the same peak and greater than the fitness value of the existing extreme point, the extreme point is discarded. If it is judged as the same peak but less than the fitness value of the extreme point, the extreme point is replaced. And if it is judged to be different peak, it is added to set .
The same peak judgement operator uses the peak exploration method to perform the same peak judgement towards the global extreme points. In this paper, the peak exploration method is applied in the same peak judgement of all the global and local extreme points. The method is as follows.
Set as an extreme point in and as the extreme point found by the PP algorithm. points are inserted between and , and is divided into equal parts. The same peak operators of parts are calculated, which is shown as
In (12), and are two adjacent points between and , . When , is point , and when , is point . If value of changes from positive to negative, it means that there is a mountain between points and . Then, and are different peak extreme points. Otherwise, they are the same peak extreme points.
3.4. Performance Analysis of PP Algorithm in Multimodal Function Optimization Problem
Let the numbers of particles be 20105. Set the acceleration factor . Set the inertia weight: , , and , respectively. The maximal repeated search number is = 210, and . Select as the test function; the comparison of the results of the PP algorithm, multipopulation coevolution strategy ES, and the standard PSO algorithm is shown in Table 1.
In Table 1, represents the mean value of the optimal solution location at the time of the test. is the average of the global solutions of the 50 tests.
From Table 1, we can see that, compared with ES and the standard PSO, the PP algorithm can search all four global extreme points of the function in 50 tests, and also there are no redundant or missing extreme points. The solution quality is significantly improved. Therefore, the PP algorithm has obvious advantages in search precision.
For function , the standard PSO has no mean computing time, and we use “—” to express time. Time of the ES is 11.2 seconds, and time of the PP algorithm is 2.5 seconds. Thus, the PP algorithm has obvious advantages in the speed of convergence.
Then, the global convergence of PP algorithm is analyzed.
Suppose. () Definition domain of the multimodal function optimization problem is a bounded closed region of . () The objective function is a continuous function on the region .
Theorem 1. Set ( is the cumulative iteration algebra) as the population sequence generated by the PP algorithm, where is the optimal point in the th generation of the population; that is, . If the objective function and the feasible region of problem satisfy hypothesis1, holds, where ; that is, the population sequence converges to the global optimal solution with probability 1.
Demonstration. Define the cumulative iteration algebra as , where is the repeated search algebra (), is the maximum iteration algebra, and is the iterative algebra.
The change of the displacement of the two generations before and after the PP algorithm is . Thus, when the population particles converge, the social cognitive part tends to 0 and less than 1, the particle velocity will drop rapidly to 0, and the particle swarm will stop moving. Therefore, the population particles may converge on the local solution and not on the global solution. For this reason, when the population particles fall into the evolutionary stagnation state, the particles are subjected to the mutation operation through the basic normal cloud generator. And the basic normal cloud generator based on the cloud model is subject to normal distribution.
Due to the mutation operation of the basic normal cloud generator with normal distribution, the displacement of particle is subject to normal distribution; that is, , where represents the normal distribution with average value 0 and variance .
Therefore, we can conclude that Theorem 1 holds.
4. Optimization of PP Algorithm Based on Cloud Mutation Clonal Selection
In order to solve the problems of the multimodal function optimization, based on the above improved methods, this paper utilizes the characteristics of fast convergent speed of the cloud mutation model to compensate the shortcomings of fast training time of the PP algorithm. Meanwhile, it can reduce the population size and redundancy. And the clonal selection algorithm can effectively maintain the diversity of population. Hence, a particle swarm optimization algorithm based on cloud mutation clonal selection is proposed (WCPP).
4.1. Optimization of Particle Space Distribution Based on Cloud Mutation
The cloud model has the characteristics of randomness and stability tendency. The randomness can keep the individual diversity so as to search for more local extreme points, and the stable tendency can protect the better individual to perform adaptive localization of the global optimal values. Thus, in order to improve the accuracy of the PP algorithm and expand the search range to find other extreme points, a cloud mutation operator based on cloud model is introduced on the basis of the PP algorithm to update the particles.
Thresholds and are given. When a particle satisfies (14), that is, the particle degenerates and degradation amplitude is smaller than , it is considered that the particle has a generation in the evolutionary stagnation state.
In (15), represents the fitness value of in the th generation. is a constant. When , it indicates the evolution of particle . And when , it indicates the degradation of particle .
The cloud mutation condition can be defined as the degradation of a particle of successive generations and the degradation amplitude is less than . When a particle satisfies the condition of cloud mutation, mutation operation is performed to the particle by a basic normal cloud generator. is set as 5 here. Since is set too small, the mutation occurs too frequently and it is easy to miss the extreme point. On the other hand, the convergent speed of the algorithm can slow down.
One-dimensional normal cloud operator can be extended to dimensions. And the cloud mutation operator for each dimension of is as follows: Begin For ; ; ; ; ; A normal random number is generated by the basic normal cloud generator, at the same time a random number is also generated. When , is updated by . End For End
represents the horizontal width of the cloud, and the greater the horizontal width, the greater the scope of the particle search. When is set as , the search range of the initial evolution is very large, which is helpful to find more extreme points. With the dynamic narrowing of the search range of , it is favorable to the fine search of particles.
“Dimensional catastrophe problem” exists in dimension function. As the dimension increases, it is difficult to find the optimal solution in each dimension through the general evolutionary algorithms. The abnormal values on individual dimension lead to poor quality of the final solution or an inability to find the optimal solution to the optimization problem. Thus, in order to prevent the particles from falling into the local extreme points and failing to find the global optimal solution in high-dimensional multimodal function, the search particles are reinitialized with a certain probability to find better solutions. When in the search state, the particles are rerandomly distributed in the search space as shown in
After the initialization, particles are updated according to (16).
4.2. Population Diversity Optimization Based on Clonal Selection Algorithm
In order to ensure the diversity of the population, optimization is performed based on the clonal selection algorithm. The steps are divided into clonal amplification, adaptive wavelet mutation, and immune selection.
(1) Clonal Amplification. In the PSO immune system, the problem to be solved is the antigen. Each particle in the particle swarm has the potential to represent each local optimization value in the multipeak optimization. Therefore, particles searched in the last iteration of the PP algorithm are all used as antibodies. The affinity function of antibodies and antigens uses the fitness function of the PP algorithm; that is,
In (17), is the fitness value of particle and is the affinity value of the antibody (the particle ). From (17), we can see that the smaller the affinity value of the antibody, the better the affinity.
The antibodies in the antibody population were sorted in ascending order according to the affinity value before the clonal amplification. The clonal amplification operator for the sorted th antibody is calculated according to
In (18), is a constant greater than 1, is a round function, and can ensure that each antibody has a constant larger than 1 with a certain number of clones. From (18), we can see that the higher the affinity, the greater the clonal multiple. In the clonal selection algorithm, when the current iteration number is , the th cloned amplified antibody is changed to after times of clone.
(2) Adaptive Wavelet Mutation. Morlet wavelet has better probability distribution of upper and lower wings. Thus, Morlet wavelet is applied to immune mutation operation to make CSA have a wider range of searching, higher precision, and good regulating performance. Morlet wavelet mother function is shown in
With the progress of , the affinity of the antibody increases constantly. When the adaptive dynamic of the antibody mutation probability reduces, the convergent speed of the algorithm can be improved. is shown as
In (20), is the maximum iteration number. When , the initial value is set as 0.7. mutation uses adaptive wavelet mutation operator based on the time-varying scale. mutation formula is shown as follows: where is an increase constant.
In (22), is the time-based mutation scale, is the random number between , and is the system parameter, which determines the dependency degree of the random number disturbance on .
In (23), is the wavelet mutation function, and is a random number between . Function is shown as
In (24), is the shape parameter of and is set to 1. is the upper limit of and is set to 10000. Thus, the value of increases with the increase of between 1 and 10000, while the range of shrinks with the increase of .
From (21), we can see that, at the early stage of the evolution, for smaller , . When the value of is large, the mutation space is also large, so that the probability of searching other extreme points can be increased. And in the middle and late stages of the evolution, adaptive of is smaller. When is approaching , , and the value of is relatively small. Thus, by improving the local fine-tuning ability, the precision of the algorithm can be improved effectively.
(3) Immune Selection. Let , and . represents an antibody obtained by clonal amplification and mutation of the th iteration. Compression of population antibody is achieved by local selection, and the optimal solution of population is not deteriorated.
4.3. Performance Analysis of WCPP Algorithm in Multimodal Function Optimization
Select as test function, the result comparison of the PP algorithm and the WCPP algorithm is shown in Table 2.
In Table 2, represents the mean value of the optimal solution location at the time of the test. represents the average of the global solutions of the 50 tests.
The mean value of the local solution of the PP algorithm is larger than that of the WCPP algorithm, which shows that the local solution of the PP algorithm has higher accuracy than the WCPP algorithm in the case that 16 local extreme points can be searched for objective function .
5. Application of WCPP Algorithm to Sliding Mode Control of Discrete Chaotic Systems
5.1. Fractional Order Chaotic System Sliding Mode Control Model
Suppose a control system as , , , . Switching function needs to be determined, , . Then, the control function can be solved as follows:
In (26), . Then,(1)the sliding mode exists: that is, (26) holds;(2)the accessibility condition can be satisfied, and the movement points outside the switching surface will reach the switching surface with a limited time;(3)the stability of the sliding mode movement is ensured;(4)the dynamic quality requirement of the control system is satisfied.
If the first three basic problems can be satisfied, we can call the control as sliding mode control.
Henon fractional chaotic system is taken as an example; the sliding mode control model of fractional chaotic system is constructed. The Henon chaotic system with control is shown aswhere , .
Define the switching function as , where . Then, , , where is the sampling time. And we can obtain (28) as follows:
When the fractional order chaotic system enters the ideal sliding mode, can satisfy the following:
That is, where . Define ; then, we can have .
From the above, we can obtain the sliding mode control equivalent control part of the fractional order chaotic system , which is shown in
5.2. Sliding Mode Control Flow of the Discrete Chaotic System Based on WCPP Algorithm
(1) Evaluation Function Selection. The selection of the evaluation function not only needs to consider the stability, fastness, and accuracy of the system but also needs to consider the control energy problem, the form of which is shown as
In (32), is the system error, is the control variable, is the initial step time, and , , are the weights. The smaller the value of the evaluation function is, the closer the corresponding function is to the global optimal solution. The evaluation function value decreases as the algorithm runs.
(2) Sliding Mode Control Flow of the Discrete Chaotic System Based on WCPP Algorithm. By using the WCPP algorithm, the parameters of the sliding mode controller are globally optimized to construct the sliding mode control method. The flow is as follows:
Step 1. Initialize the sliding mode controller, other parameters, and a group of particles with size . The initialization parameters include maximum repetition search number of hybrid algorithm and the WCPP maximum number of iterations .
Step 2. The evaluation function value of each particle is obtained by control of the discrete chaotic system according to the position of each particle and the optimal position of individual particle and the optimal position of all particles are recorded.
Step 3. If the repeated search algebra is larger than , the iteration is stopped, and the algorithm is finished. The output of the sliding mode controller controls the discrete chaotic system. Otherwise, turn to Step4.
Step 4. If the WCPP algorithm iteration number , If the condition of the cloud mutation is satisfied, the cloud mutation is performed on the particles. Else Speed and positions are updated on the particles. End If Update particles and , and . End If
Step 5. The particles in the last search generation of the WCPP algorithm are cloned and amplified.
Adaptive wavelet mutation operator performance is carried out.
Immune selection performance is carried out.
update particles and ,
Step 6. , the positions and the speed are reinitialized.
Step 7. The searched optimal fitness particles are assigned to the initialized particles and turn to Step2.
In order to improve the fastness of the sliding mode control method for the discrete chaotic systems based on the WCPP algorithm, the values of the algorithm parameters and are small.
Figure 2 is the discrete chaotic system sliding mode control structure based on the WCPP algorithm.
5.3. Algorithm Simulations
A discrete chaotic system is taken as the simulation object; the model of the controlled object is shown as
For the above discrete chaotic system, simulations of step input signal are implemented based on the WCPP and other two algorithms, which is shown in Figure 3. In this simulation, the disturbance signal is set as . And at this time the WCPP optimized parameters are .
From Figure 3, we can see that the WCPP algorithm has very strong search ability. It can search the global optimal value of the sliding mode controller parameters without overshoot, and the error, the response speed, and adjustment time all achieved the optimal performance. Therefore, the control effect of the WCPP algorithm is superior to the control effect of the PSO and the PP algorithm.
For the above discrete chaotic system, the simulations of the square wave input signal are implemented by the WCPP algorithm and other two algorithms, respectively, which is shown in Figure 4.
Another discrete chaotic system is also taken as the simulation object; the model of the controlled object is shown as
For the above discrete chaotic system, simulations of step input signal are implemented based on the WCPP and other two algorithms, which is shown in Figure 5. The WCPP optimized parameters are .
From Figure 5, we can find that the WCPP algorithm is basically no overshoot and is able to achieve the desired value in the second iteration. The performances of the WCPP algorithm in error, response speed, and adjustment time are significantly better than those of PSO and PP algorithms. Therefore, WCPP algorithm has strong searching ability and good control quality.
For the above discrete chaotic system, the simulations of the square wave input signal are implemented by the WCPP algorithm and other two algorithms, respectively, which is shown in Figure 6.
Chaos and chaos systems are widely used in physics, engineering, biology, finance, and so forth. The modeling of complex systems based on the theory of calculus can more accurately reflect the dynamic characteristics of the systems. Discrete chaotic systems can be widely used in the field of chaotic secure communications. Thus, it is of great theoretical and applicable value in the applications of control and synchronization of the chaotic dynamic systems. In order to overcome the shortcoming of the standard PSO algorithm in multimodal function optimization problems, a sliding mode control model of discrete chaotic systems based on the coupled multimode functions is proposed. The simulation results show that the improved algorithm proposed in this paper has significantly improved the response speed and overshoot performance compared with the standard PSO algorithm.
Conflicts of Interest
The author declares that there are no conflicts of interest regarding the publication of this paper.
This work was supported by the National Natural Science Foundation of China (Grant no. 51675490).
Y.-H. Chang, C.-I. Wu, H.-W. Lin, H.-C. Chen, and C.-W. Chang, “Fractional order integral sliding-mode flux observer for direct field-oriented induction machines,” International Journal of Innovative Computing, Information and Control, vol. 8, no. 7, pp. 4851–4868, 2012.View at: Google Scholar
H. Liu, L. Guo, and Y. Zhang, “Composite attitude control for flexible spacecraft with simultaneous disturbance attenuation and rejection performance,” Proceedings of the Institution of Mechanical Engineers. Part I: Journal of Systems and Control Engineering, vol. 226, no. 2, pp. 154–161, 2012.View at: Publisher Site | Google Scholar
I. N'Doye, M. Zasadzinski, M. Darouach, and N.-E. Radhy, “Observer-based control for fractional-order continuous-time systems,” in Proceedings of the 48th IEEE Conference on Decision and Control held jointly with 28th Chinese Control Conference (CDC/CCC '09), pp. 1932–1937, December 2009.View at: Publisher Site | Google Scholar