1. Introduction

The cell mapping methods were created by Hsu in 1980s for global analysis of nonlinear dynamical systems that can have multiple steady-state responses including equilibrium states, periodic motions, chaotic attractors as well as domains of attraction of these steadystate responses. The cell mapping methods have been applied to global analysis and control design of deterministic, stochastic, and fuzzy dynamical systems. There have been several survey articles that present thorough discussions of the literature on the cell mapping methods and their applications (Sun & Luo 2012, Hong & Sun 2006d, Xu et al . 2013). A comprehensive review of the global analysis with the cell mapping method by Sun and Luo (2012) provides rich content on engineering applications and algorithm development. A thorough review of the progress in global analysis of nonlinear dynamics and its influence on the analysis, control, and design of mechanical and structural systems is presented by Rega and Lenci (2015). This paper presents a discussion of the literature of some control applications and algorithm developments of the cell mapping methods.

Two important extensions of the cell mapping methods have been developed to improve the accuracy of the solutions obtained in the cell state space. The first is the interpolated cell mapping which uses the cell mappings as a foundation to calculate point-wise solutions without further numerical integrations of differential equations. The second is the subdivision technique of the set-oriented method for improving the accuracy of the invariant solutions obtained with the cell mapping methods. For a long time, the cell mapping methods have been applied to dynamical systems with low dimension until now. With the advent of inexpensive computer memories and massively parallel computing technologies such as the graphical processing units (GPUs), global analysis of moderate- to high-dimensional nonlinear dynamical systems becomes feasible.

The cell mapping methods propose to discretize the continuum state space and the time. The discrete space consists of a finite collection of cells. The dynamical systems that originally obey ordinary or partial differential equations are now represented by the mappings in the cell state space, called the cell-to-cell mapping, or cell mapping for short. The cell mappings describe the system evolution over a short time in a finite region of interest in the cell state space. More importantly, long-term system responses such as periodic motion, equilibrium points, limit cycle, chaotic motion, domains of attraction, and stable and unstable manifolds of saddle points can be all obtained from the cell mappings.

2. Cell mapping methods

There are two versions of the cell mapping methods: the simple cell mapping (SCM) and the generalized cell mapping (GCM). This section presents a brief introduction of both SCM and GCM methods.

2.1 Simple cell mapping

We assume that the nonlinear dynamical system is described by a point mapping given by

where k is the iteration step, x_k is the n-dimensional state vector at the k-th step. Consider a finite region U ⊂ Rⁿ where a sufficiently rich dynamics of the system resides. We discretize U into a collection of small, finite size boxes, known as the cells. Since the region U is finite, there will be a finite number of cells in the discretized region. Hence, each cell in the collection can be numbered by one integer, denoted as z.

The SCM accepts only one image cell for a given pre-image cell, or domain cell. In other words, in SCM, the dynamics of the system starting from one cell with a small but finite volume is represented by that starting from a point in the cell, usually the center of the cell, leading to an integer-valued mapping

where C(·) is a symbolical notation of the integer mapping, and z_k is an integer representing the cell where the system resides at the k-th step. Usually, C(·) has to be constructed numerically. The region out of the domain U is called the sink cell. If the image of a cell is out of the domain of interest, we say that it is mapped to the sink cell. The sink cell always maps to itself.

Properties of SCM

Because there are only a finite number of cells in U, the integer mapping in Eq. (2) eventually will revisit the cells in the path. The revisited cells hence form closed groups called periodic groups. The minimum period of these groups is one, while the maximum possible period is equal to the total number of cells in U. For the group with period one, we have

for the cell z in the group.

The simple cell mappings z_k+1=C(z_k) are stored in a single array of length N_t where N_t is the total number of cells in U. For example, let C(i) denote the image array. If C(i)=j, then cell z=j is the image of cell z=i. This array can be viewed as the storage of a sparse matrix representing the simple cell mappings over where denotes the set of integers indexing the cells in U including the sink cell. The sparse matrix reads

The array C(i) contains the forward dynamics of the system in time. The stable steadystate responses of the system including equilibrium points, periodic orbits and chaotic motion form periodic groups in C(i).

We can also store the pre-image information of an image cell in an array, denoted as C^-1(j). That is, i=C^-1(j). In terms of the matrix p_ji, the sparse matrix of the backward dynamics is simply the transpose of the forward dynamics matrix.

The backward dynamics provides an important role in the global analysis of nonlinear dynamical systems. In the backward dynamics, the unstable responses appear to be stable. Consider a search starting from the stable steady-state responses, i.e., the identified periodic groups. If we search along the backward dynamics using C^-1(j), we would identify the domains of attraction of the stable responses. The backward search, on the other hand, reveals the boundaries of the domains of attraction, which are usually outlined by the unstable manifolds in the saddles.

2.2 Generalized cell mapping

The GCM accepts multiple images for a pre-image cell. This is consistent with the fact that the cell with a finite volume will evolve to cover multiple cells under the system dynamics over a finite time. For deterministic and stochastic systems, the GCM leads to a Markov chain representation of the dynamical system with the transition of probabilities given by

or in the component form

where p(k)={p_i(k)} denotes the probability that the system resides in the i-th cell at the k-th step, and P (k)={P_ij(k)} is the one step transition probability from the j-th cell to i-th cell at the k-th step. N_t is the total number of cells in the computational domain. When the matrix P (k) is independent of k, the Markov chain is said to be stationary. Otherwise, it is non-stationary.

The rich literature on Markov chains and later the graph theory has provided us highly effective algorithms for analyzing the GCM (Hsu 1982, 1995). The analysis of the GCM leads to the discovery of invariant sets, stable and unstable manifolds of saddle-like equilibrium states, domains of attraction and their boundaries. The invariant sets represent stable equilibrium states, periodic or chaotic motions. The invariant sets are called the persistent groups in the Markov chain literature (Chung 1967).

The stable and unstable manifolds of saddle-like equilibrium states, domains of attraction and their boundaries are represented by the so-called transient cells.

The stationary transition probability matrix P can be partitioned into the following canonical form, also known as the normal form

where P_i is a square matrix representing the transition probability matrix among the cells in the i-th persistent group, Q_i is associated with the i-th open communicating group. The cells in the group Q_i are transient. T_ij and R_ij represent the evolution paths from transient cells to stable and unstable attractors, respectively. Q_i often contains the saddle like attractors, unstable equilibrium points and unstable periodic orbits.

The ability of the GCM method to conduct global analysis of nonlinear dynamics is fully illustrated by the topological structure of the transition probability matrix P in the normal form. We can use the GCM method to discover invariant sets, stable and unstable manifolds of saddles, unstable solutions and domains of attraction of invariant sets of nonlinear dynamical systems. The unstable solutions as well as stable manifolds of saddles can be found by the backward cell mapping (Sun & Luo 2012, Dellnitz & Hohmann 1997, Dellnitz & Junge 2002).

Consider the transition probability sub-matrix P _k associated with the k-th persistent group. If this group of cells represents the period-K motion of the system, it can be partitioned into the following form

where P_{k, j} (1≤j≤K) is a sub-matrix of a certain dimension. When P_k has at least one non-zero diagonal element, the period of the persistent group is one. It is called an aperiodic group.

3. Set-oriented method

As an extension of cell mapping, Dellnitz and colleagues introduced the set-oriented method with the sub-division technique that is capable of obtaining the invariant sets of nonlinear dynamical systems with high accuracy (Dellnitz & Hohmann 1997). The setoriented method starts with relatively large cells and removes the cells that do not contain part of the invariant set by sampling a number of initial conditions from each cell. The sub-division is then applied to the retained cells. This is how the set-oriented method gains computational efficiency. The set-oriented method is an effectively mixed application of SCM and GCM with the sub-division applied to the covering set of invariant solutions.

There have been many studies of the set-oriented method. An adaptive sub-division algorithm was developed (Dellnitz & Junge 1998) that allows the existence of multiple different cell sizes to cover the solution. A study of non-smooth mechanical system was carried out by Neumann et al . (2007) with the set-oriented method to find global attractors. The algorithm for extracting unstable manifold and saddle solutions was introduced by Dellnitz and Junge (2002). The set-oriented method is also a robust tool for designing optimal controls (Junge & Osinga 2004, Grune & Junge 2005), especially for multi-objective optimal controls (Schütze et al . 2013, Blesken et al . 2009).

The set-oriented method with sub-division has not been applied to investigate the transient dynamics of the system such as the domains of attraction and basin boundary. On the other hand, the cell mapping methods were developed for comprehensive global analyses of nonlinear dynamical systems including the discovery of invariant sets and transient dynamics.

4. Interpolated cell mapping

The sub-division technique leads to smaller and smaller cells, and therefore, the accuracy of the solution for invariant sets improves. At some point, the sub-division has to stop. This is where another important extension of the cell mapping methods comes in: the interpolated cell mapping (ICM) method (Tongue 1987; Tongue & Gu 1988b, 1988c). The ICM uses the simple cell mappings to interpolate the image of a point without integrating the differential equation using this point as an initial condition. The simple cell mappings in the refined cells provide a database for interpolation. The ICM method is able to construct very fine solutions of invariant sets, which assumes that the simple cell mappings are on a sufficiently small grid and that the underlying dynamics of the system is smooth enough for interpolation.

The local interpolation error of ICM is of order O(h²) with the linear interpolation, where h is the cell size, whereas the accuracy of SCM is of order O(h) (Lee & Hsu 1994). More adjacent cells around the point of interest can be used to construct high order interpolations to further improve the accuracy (Tongue & Gu 1988a). A modified ICM by introducing the sampling idea of GCM was proposed to further increase the capability of ICM to capture the boundaries of domains of attraction (Ge & Lee 1997). Several typical nonlinear systems have been studied with the ICM method including the Lorenz system (White & Tongue 1995), a forced beam with cubic nonlinearity (Lee & Ghang 1994) and a spring-pendulum system (Lee & Hsu 1994). The nonlinear system identification approach using the method of interpolated cell mapping is proposed by Bursal and Tongue (1992).

Previous studies of the ICM have dealt with low dimensional state spaces. Extension of the interpolation scheme to much higher dimensional state space is non-trivial. This is because the number of adjacent cells is given by 3ⁿ -1 where n is the dimension of the state space. If all the adjacent cells are used for interpolation of the mapping inside the cell of interest, the number of evaluations will grow exponentially. Xiong et al . (2015) have developed a scheme using only 2n adjacent cells plus the cell under consideration for interpolation in high dimensional state space with the accuracy of order O(h²). The fine structure of the strange attractor of a six-dimensional Lorenz system has been depicted in their study, as shown in Fig. 1.

Figure  1.  Two-dimensional projections of a six-dimensional strange attractor of the Lorenz system. Blue dots are the centers of the cells in the invariant set. Red dots showing the fine structure of the attractor are generated with interpolation

DownLoad: Full-Size Img PowerPoint

If we put the set-oriented method with sub-division and the ICM method in the framework of the cell mapping methods, it becomes apparent that the ICM method represents a post-processing step to extract point mappings from the cell mappings on a refined partition of the cell state space.

A note

Both the set-oriented method and ICM method represent efforts to increase the computational efficiency for finding invariant sets of nonlinear dynamical systems with much improved accuracy. The accuracy of the solutions obtained by the cell mapping method was compared with the accuracy of the point-wise solutions obtained by numerical integration. Such a comparison and pursuit of point-wise accuracy are beyond the original purpose of the cell mapping methods.

As discussed by Hsu on several occasions, the goal of the cell mapping methods is to quickly discover the general global structure of the responses of nonlinear dynamical systems with a reasonable accuracy. In other words, the cell mapping methods can answer these questions with high efficiency: How many stable and unstable responses of a nonlinear dynamical system exist in a certain region of the state space? Where are they and how are they connected?

The cell mapping methods cannot deliver the fine structure of fractal dimensional objects such as basin boundaries and strange attractors, because highly accurate numerical integrations must be done to find the fine structure of fractal dimensional geometry. Nevertheless, the GCM can tell where in the state space fractal dimensional objects may exist and can outline their shape.

5. Global analysis of complex systems

The cell mapping methods were originally created for the global analysis of nonlinear dynamical systems. These methods also caught the attention of the controls community early on. It turns out that there are more challenging issues of nonlinear control systems that require global analyses.

5.1 Nonlinear control systems

The global behavior including equilibria, periodic motions, and their domain of attractions of a fuzzy dynamical system is studied with the cell mapping method with the min-max operation by Chen and Tsao (1988, 1989). The periodic solutions of a fuzzy-knowledge-based mobile robot motion control system and the domain of attraction of the periodic solution are found with the cell mapping method by Fei and Isik (1990). The concept of cell state mapping is combined with the synthesis techniques of maximum entropy self-organizing net in the fuzzy model identification to design fuzzy controls (Lin & Isik 1997). The cell mapping method is utilized to generate the rules of the fuzzy controller by systematically generating near-optimal trajectories for all possible initial states in the parking lot maneuvering a car (Leu & Kim 1998). A continued propagation cell mapping (CPCM) is developed to design a TSK-type fuzzy logic controller and a sliding mode-type controller for the uncertain nonlinear system (Rizk & Smith 1997). The cell mapping method with variable time steps is used to design fuzzy logic controls (Smith & Comer 1990, 1991). A systematic method for designing and evaluating a fuzzy logic controller based on the cell-to-cell mapping technique and genetic optimization is presented by Tai et al .(1995, 1996). The resulting controller exhibits near time-optimal control and improved system performance. A fuzzy cell-mapping control algorithm, formerly proposed for motor control, is extended to the application of satellite attitude manoeuvring/stabilizing control (Yen & Tarng 1996). The cell-mapping method treats the complex satellite attitude dynamics, and the fuzzy interpretation helps achieve a suitable control effort. The result is a smooth transition from manoeuvring to stabilization without overshoot while maintaining an optimal control performance.

The cell mapping method is used to study the effect of quantization in digital control systems (Wang et al . 2007). The periodic solutions and their domains of attraction for flexible systems under nonlinear feedback control are studied by Borre and Flashner (2011, 2012a, 2012b). The switching surface of the relay type control is chosen as the Poincaré section for the construction of the cell mappings. The simple cell mapping method is used for the global analysis of the dynamical behavior of nonlinear gear system with different degrees of wear faults (Wang et al . 2011, 2012).

5.2 Robotics

A number of interesting studies of global analysis of biped robots are available in the literature. The cell mapping method for a Poincaré return map is used to identify the initial states, i.e., the domains of attraction, of the stable gaits that match the structure parameter values (Liu & Tian 2009). The method is combined with the control strategy and is also used to search for the tracking target point, when the slope of the ground or the structure parameters of the underactuated biped robots have changed (Liu et al . 2009). The basins of attraction of stable passive walking models including the straight leg model and the model with knees are discovered with the cell mapping method (Liu et al . 2007). The cell mapping method is used to find stable limit cycles of a compass-gait biped with gait asymmetry as the parameters are varied (Moon & Spong 2010). The analysis shows that passive dynamic walking has multiple attractors, and marginally stable limit cycles exhibit not only period doubling, but also period reemerging, disconnecting, and disappearing. The global stability of passive biped robot was analyzed by a gradual point mapping-cell mapping (GPCM) method (Zang et al . 2013). The cell mapping method is combined with Newton-Raphson iteration to obtain the limit cycle of the periodic gait in the model, and the track stability of the limit cycle is analyzed by Zhao et al . (2011).

5.3 Path planning

Path planning is another application of the cell mapping method. An efficient cell space approach to represent obstacles in robotic navigation. A novel and efficient search algorithm for neighboring cells with the fast generation of all the neighboring cell addresses is proposed by Kim and Lee (2009). A cell-mapping method is introduced for planning global trajectories of robotic manipulators in cases where the cell space is composed of combination pairs of plane cells. Optimal trajectory problems in the free field and in the obstacle-constrained field are studied by Zhu and Leu (1990). The cell mapping method is applied to robot path planning problem with obstacle avoidance (Xue et al . 2014). A multi-objective optimization approach based on the cellular automata and cell mapping is proposed for robot path planning in complex terrain under radar surveillance (Naranjani & Sun 2015). The landscape analysis is carried out with the cell mapping methods by Kerschke et al . (2014) and Hernández et al . (2014).

5.4 Optimal control

The cell mapping methods have been applied to optimal control problems of deterministic and stochastic dynamic systems (Hsu 1985, Bursal & Hsu 1989, Crespo & Sun 2003b). Other interesting applications of the cell mapping methods include optimal space craft momentum unloading (Flashner & Burns 1990), single and multiple manipulators of robots (Zhu & Leu 1990), optimum trajectory planning in robotic systems by Wang and Lever (1994), tracking control of the read-write head of computer hard disks (Yen 1992), and airfoil flutter analysis Ding et al . (2005). Sun and his group studied the fixed final state optimal control problems with the simple cell mapping method (Crespo & Sun 2000a, 2000b), and applied the cell mapping methods to the optimal control of deterministic systems described by Bellman’s principle of optimality (Crespo & Sun 2003a). The cell mapping method is used to generate the general optimum trajectories for driving intelligent vehicles equipped with digital maps (Li & Wang 2002, 2003a, 2003b). Multiple objectives are considered including the minimum time, energy, and jerk trajectories. A summary of control studies of nonlinear dynamic systems using the cell mapping method is presented (Sun 2013).

The optimal control studies with the cell mapping methods represent the state space design approach. In most of the early studies, the optimal controls were obtained as a function of the state. Another control design approach is done in the parameter space. This will be discussed in Sect. 6.2.

5.5 Global analysis

Global analysis of nonlinear dynamical systems remains a popular application area of the cell mapping methods. The chaotic boundary crisis in the Duffing van der Pol vibroimpact oscillator is studied with the generalized cell mapping method (Feng & Xu 2011). The results suggest that the boundary crisis is associated with the tangency of the stable and unstable manifolds of the saddle. The evolution of the global structure of the coupled neural oscillators into the chaotic itinerancy is investigated by using an extended point mapping under cell reference (PMUCR) method (Jiang & Guo 2011). This method aims to retain the accuracy of point mapping while enhancing its computational efficiency. Global domains of attraction of multiple lock-in attractors of the small aerodynamic asymmetric rolling flying system are determined by numerical method of PMUCR (Sun et al . 2015).

The stochastic response of nonlinear oscillators under periodic and Gaussian white noise excitations is studied with the generalized cell mapping based on short-time Gaussian approximation (GCM/STGA) method. Both the transient and steady-state probability density functions (PDFs) of a smooth and discontinuous (SD) oscillator are computed to illustrate the application of the method. The effect of a chaotic saddle on the stochastic response is also studied. The stochastic P-bifurcation in terms of the steady state PDFs occurs with the decrease of the smoothness parameter, which corresponds to the deterministic pitchfork bifurcation (Han et al . 2016). The fixed interval smoothing is numerically approximated by recursive computation of the conditional density as a piecewise constant function, which is a coarse-grained representation of the system dynamics as an approximate aggregate Markov chain in discretized state space or cell space (Ungarala 2012).

The cell mapping method has been adopted for analysis and control design of fuzzy dynamical systems as early as in 1980s. However, the global analysis of fuzzy nonlinear dynamical systems remains a challenge. The response of fuzzy nonlinear dynamical systems is naturally global in the sense that it must be described by both a membership distribution and its geometry. Hong and Sun (2006a, 2006b, 2006c) have been leading the effort in developing fuzzy generalized cell mapping (FGCM) method for global analysis of nonlinear dynamical systems. They have studied various global fuzzy dynamics by using the FGCM method including bifurcations and blue sky catastrophes and crises in chaotic systems, as well as transient fuzzy responses with evolutionary membership distributions (Hong et al . 2015a, 2015b).

5.6 Parallel computing

For a long time, the cell mapping methods have been applied to dynamical systems with low dimension until now. With the advent of cheap dynamic memory and massively parallel computing technologies, such as the GPUs, global analysis of moderate-to high-dimensional nonlinear dynamical systems becomes feasible. Recent application of parallel computing with cell mapping technique has been reported by Eason and Dick (2014) where multi-core CPU architecture is used to speed up global analysis of nonlinear systems. In another study (Xiong et al . 2015), the SCM and GCM are implemented in a hybrid manner combined with the subdivision technique to enhance the accuracy of the steady-state responses. The ICM is used as a post-processing step to generate the point-wise approximation of the solutions without additional numerical integrations of differential equations. The cell mapping methods are applied to a nonlinear dynamical system with six-dimensional state space in this work.

6. Multi-objective optimization

Engineering systems such as controls are often designed to meet multiple and often conflicting objectives. To design these systems to meet the conflicting objectives in an optimal manner leads to multi-objective optimization problems (MOPs). An MOP can be stated as follows

where f_i : Q → R¹, F : Q → R^k. f_i are objective functions, k ∈ Q is a q-dimensional vector of design parameters. The domain Q⊂ R^q is the design space, and can in general be expressed in terms of inequality and equality constraints

Next, we define optimal solutions of the MOP by using the concept of dominance (Pareto 1971).

Definition 1

(a) Let v, w ∈ R^k. The vector v is said to be less than w (in short: v < _p w), if v_i < w_i for all i ∈ {1, 2, …, k}. The relation≤_p is defined analogously.

(b) A vector v ∈ Q is called dominated by another vector w ∈ Q (w ≺ v) with respect to MOP (10) if F (w)≤_p F (v) and F (w) ≠ F (v), otherwise v is called non-dominated by w.

If a vector w dominates a vector v, then w can be considered as a “better” solution of the MOP. The definition of optimality or the “best” solution of the MOP is now straightforward.

Definition 2

(a) A point w ∈ Q is called Pareto optimal or a Pareto point of MOP (10) if there is no v ∈ Q which dominates w.

(b) The set of all Pareto optimal solutions is called the Pareto set denoted as

(c) The image F() of is called the Pareto front.

The Pareto front typically forms (k -1)-dimensional manifolds under certain mild assumptions on the MOP (Hillermeier 2001).

6.1 MOP search algorithm as cell mapping

The SCM method divides the continuous parameter space Q into a collection of finite size cells. The number of cells in a finite domain Q is finite. Hence all the cells in Q can be sequentially indexed with one integer. A cell is represented by its central point in the search for the Pareto set.

There are two ways to construct the cell-to-cell mapping for MOP application with gradient-based and gradient-free search algorithms. Gradient-based methods generate point mappings, which can be converted to cell mappings. Gradient-free methods directly build cell mappings by using the objective functions evaluated at the centers of cells in Q.

Let v denote the searching direction such that

where γ is a step length for the search in the parameter space Q. γ is selected to meet the following dominance condition

If the dominance condition is not satisfied, the pre-selected step length will be iteratively cut by half. Two outcome can occur after a few iterations: (1) The dominance condition is satisfied for a k_n+1 which lies in a cell denoted as z_n+1. In this case, a cell mapping can be constructed and denoted as z_n+1=C [z_n]. (2) γ is so small that kn+1 lies in the same cell as k _n. In this case, the cell k_n is considered to be a candidate for the Pareto set and the cell mapping is set to be z_n=C [z_n].

The cyclic cells of the SCMs represent the solutions of the MOP. The cell mapping methods with sub-division have been applied to MOPs (Dellnitz & Schutze 2005, Jahn 2006, Schütze et al . 2009). Several interesting benchmark mathematical problems are studied. The Pareto solutions are compared with the available exact solutions (Naranjani et al . 2014, 2016). A parallel computing algorithm of SCM for MOPs is first proposed in (Fernández et al . 2016). This work lays a foundation for the cell mapping method to attack high dimensional problems. A hybrid method that combines the popular evolutionary algorithms such as genetic algorithm with the cell mapping method is developed by Naranjani et al . (2016). The new hybrid method takes advantages of both the evolutionary algorithms and the cell mapping method. It first implements an evolutionary algorithm to generate a set of random Pareto solutions in the design space. A set of cells in the discretized design space is then identified that contains all the random Pareto solutions. This is known as the covering set. The SCM is then applied only to the covering set to search and recover the global solution of the Pareto set. This is a highly efficient and effective way to solve for the high dimensional MOPs. The cell mapping method can also be implemented to find nearly-optimal solutions of the Pareto set with pre-selected tolerance (Hernández et al . 2013b).

6.2 Multi-objective optimal control design

Full state feedback control is an important part of the modern control theory. Since feedback controls are often designed to meet multiple and possibly conflicting performance goals, comprehensive studies are usually carried out to tune control gains in order to achieve the best overall performance (Cominos & Munro 2002, Wang et al . 1999). Designing feedback controls to meet multi-objectives naturally leads to an MOP defined in the parameter space. Since the solutions to an MOP form a Pareto set in the parameter space, such a design approach provides a wide range of choice representing various compromises of the conflicting objectives.

Multi-objective optimal control design can be carried out in time domain or frequency domain. Time domain approach uses the time domain specifications of the closed-loop response as the objective functions such as overshoot, peak time, settling time and tracking error (Kumar & Nair 2011). On the other hand, frequency domain design uses phase and gain margins as the objectives, and can consider robust issues such as model uncertainty, load disturbance and measurement noise. Multi-objective optimization with robustness often involves the optimization among several norms. Vroemen and De Jager (1997) reviewed the multi-objective design of robust controls for linear systems. They examined different combinations of , and L₂ norms to formulate the robust control synthesis problems. A more recent overview by Gambier and Badreddin (2007) summarized most available methods for multi-objective optimal control design in both time and frequency domain. They stated that despite the significant development of multi-objective optimization in control engineering, on-line design methods with multi-objective optimization are still at the beginning phase.

Even though there have been many studies of multi-objective optimization control design for linear systems, only a handful references are available for nonlinear systems, and are scattered in different disciplines. Since the concept of frequency domain in nonlinear system is not as well studied as in linear systems, the control design for nonlinear systems is usually done in time domain. A nonlinear fuzzy controller based on Pareto rule-base design is carried out by examining the temporal response by Zhao and Tsu (2003). A variable complexity modelling technique with multi-objective optimization design was studied by Silva et al . (2006) to tune the multivariable proportional integral (PI) control of a nonlinear thermodynamic model in gas turbine. A more theoretical research of multi-objective nonlinear control is presented by Azhmyakov (2008) where the multi-objective optimization algorithm is combined with the classical variational method.

Many algorithms for obtaining the Pareto set and Pareto front of MOPs have been developed. There are biologically inspired optimization algorithms such as genetic algorithm (Panda 2011), ant colony optimization (Chiha et al . 2012), immune algorithm (Khoie et al . 2011), and particle swarm optimization (Solihin et al . 2011). All these methods have been successfully applied to feedback control design including proportional-integral derivative (PID) controls to meet multiple objectives. Fliege and Saviter (2000) have developed several gradient-based algorithms by converting MOP to single-objective optimization problem (SOP) for point-wise evolution and step length determination of the steepest descend search for MOP solutions. Bosman (2012) expands the concept of gradient by introducing novel geometric transformations and combines it with the genetic algorithm for MOPs. A gradient-free approach is introduced by Zhong et al . (2010) to address MOPs with undifferentiable objective functions. In the work of Custodio et al . (2011), methods for pattern searching are adopted to direct gradient-free search.

Another approach to find the Pareto set is to use the set oriented methods with subdivision techniques (Dellnitz & Schütze 2005, Jahn 2006, Schütze et al . 2009). The advantage of the set oriented methods is that they generate an approximation of the global Pareto set in one single run of the algorithm. The SCM method can discover the global Pareto fronts with fine structures in a quite effective manner for low and moderate dimensional problems (Hernández et al . 2013a, Naranjani et al . 2013). Sun and his colleagues studied the multiobjective optimal control design using both SCM and GCM (Naranjani et al . 2013, Xiong et al . 2014b).

The multi-objective optimal designs developed by Xiong et al . (2014b) are for nonlinear sliding mode controls. The MOP is moderately high dimensional with five objective functions and six design parameters. The SCM is implemented in parallel computing on a PC with GPUs. The multi-objective optimal controls are implemented on an underactuated flexible mechanical link mounted on a rotating platform. The design parameters k=[α_a, λ_a, α_u, λ_u, η, ϕ] are the gains and parameters of the sliding mode control. The objectives t _{p, x1}, M_{p, x1}, max |x₂|, IAE_x1, and IAE_x2 consist of transient response characteristics and integrated tracking performance measures.

The Pareto front of the multi-objective optimal design of the sliding mode control is shown in Fig. 2. Figure 3 shows an example of the experimental results of tracking control of the flexible link following the command of the square wave. The optimal designs are first validated in the numerical simulations and then are used in the experiments. The ability of the optimally designed sliding mode control is clearly demonstrated in the experiments.

Figure  2.  The 5-dimensional Pareto front of the objective functions, projected on to 2-dimensional sub-spaces of the objective space. The color code indicates the level of M_{p, x1}, IAE_x2, max |x₂|, and t_{p, x1} in subplots from upper left in counterclockwise order. The conflicting nature among these objectives can be observed clearly. For example, the overshoot and peak time are conflicting for tracking control, which can be seen in the upper left plot of x₁

DownLoad: Full-Size Img PowerPoint

Figure  3.  Experimental square wave tracking response of the rotary flexible joint under a Pareto optimal sliding mode control. x₁ is the base angle and x₂ is the angle of the flexible link

DownLoad: Full-Size Img PowerPoint

7. Zero finding of nonlinear equations

Finding zeros of multi-variable nonlinear functions is a common problem existing in many scientific and engineering fields. In the area of dynamics, finding equilibrium states of nonlinear systems, bifurcation and stability analysis of the system all lead to zero finding of nonlinear functions. In control systems, the stability region in the controller parameter space can also be transformed to a zero finding problem. General zero finding problems can be expressed as f(x)=0 with f : R^m → PRⁿ and x ∈ U ⊂ PR^m where U is in a bounded region in R^m.

Either gradient based or gradient free point-to-point iterative search algorithm for finding zeros forms a dynamical system that evolves to the potential solutions. Hence, finding function zeros can be equivalently treated as finding global invariant sets of such iterative dynamical systems. Both the cell mapping and set-oriented methods were originally developed for finding global invariant sets. The set-oriented method has shown great performance with the capability of locating all solutions of nonlinear algebraic equations in both real and complex domains (Dellnitz et al . 2002a, 2002b).

Since analytical solutions for zeros of nonlinear functions are in general difficult to obtain, there have been many studies of numerical methods for zero finding. Classical Newton’s method with gradient information has been successfully applied to various problems for a long time (Madsen 1973). A number of novel variations of Newton’s method are popular choices for many applications (Carniel 1994, Bhaya & Kaszkurewicz 2004). Other algorithms are focused on the non-smooth or complex functions where the derivatives are not needed (Chandrupatla 1997).

To address the problem of finding global solutions in certain domain, intensive studies have been carried out to take both gradient based or gradient free algorithms as underlying dynamics and study their long term evolutionary status in parameter space. The homo-topy continuation method (Liu 1990), cell mapping (Carniel 1994), and set-oriented method (Dellnitz et al . 2002a, 2002b) have been applied by many scholars to attack the problem of global searching. The homotopy continuation method is performed in continuous pointwise parameter space while the latter two methods are performed in discrete cellular space. For problems with moderate to high dimensions, the point-wise methods become less feasible due to the increasing need of computational efforts. The set-oriented method and its predecessor, the cell mapping method, are computationally more effective.

An algorithm using the simple cell mapping and generalized cell mapping that can find zeros of multi-variable nonlinear functions in an efficient manner is presented by Xiong et al . (2014a). The SCM with sub-division is applied to find zeros of nonlinear algebraic equations and stability boundaries of control systems in the parameter space (Xiong et al . 2016).

Stability boundary as zero finding problem

We now show an example of finding the stability boundary of a linear time varying system subject to delayed feedback control. The system is known as the Mathieu equation.

By setting ε=1, δ=4, the uncontrolled system is parametrically unstable. We consider a proportional and derivative delayed feedback control and look for the stability boundary in the gain space (k_p, k_d).

The semi-discretization method by Sun and Song (2012) is applied to find a mapping of the extended state vector y_j+1=Φ_yj over a period. The stability boundary is determined when the maximum absolute value max |λ| of eigenvalues of mapping Φ is equal to unity. Hence, the problem of finding stability boundary becomes a zero finding problem of the implicit function defined below

The hybrid GCM-SCM algorithm is applied to find the stability boundary. We choose the time delay τ=π/4 here. A coarse cell space partition 10 × 10 of the domain [k_p, k_d] ∈ [-8, -5] × [-8, 1] is taken initially. After three subdivisions, we find the stability boundary shown in Fig. 4.

Figure  4.  Stability boundary of the Mathieu system in the feedback gain space (k_p, k_d)

DownLoad: Full-Size Img PowerPoint

8. Concluding remarks

In the past three decades, the cell mapping methods have received continuous attention from research communities all over the world. New applications and new algorithm developments have occurred and will continue. It is anticipated that more applications of the cell mapping method will emerge in the decades to come as scientific and engineering research is becoming more and more data-driven and computationally intensive. Such a trend is further supported by the ever more powerful and inexpensive new supercomputer technologies.