UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Modeling of linear systems with parameter variations : applications in hard disk and ball screw drives Sepasi, Daniel 2011-12-31

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
[if-you-see-this-DO-NOT-CLICK]
ubc_2011_fall_sepasi_daniel.pdf [ 6.95MB ]
[if-you-see-this-DO-NOT-CLICK]
Metadata
JSON: 1.0080734.json
JSON-LD: 1.0080734+ld.json
RDF/XML (Pretty): 1.0080734.xml
RDF/JSON: 1.0080734+rdf.json
Turtle: 1.0080734+rdf-turtle.txt
N-Triples: 1.0080734+rdf-ntriples.txt
Original Record: 1.0080734 +original-record.json
Full Text
1.0080734.txt
Citation
1.0080734.ris

Full Text

Modeling of Linear Systems with Parameter Variations: Applications in Hard Disk and Ball Screw Drives by Daniel Sepasi B.Sc., Mechanical Engineering, Sharif University of Technology, 2005 M.A.Sc., Mechanical Engineering, University of British Columbia, 2007 a thesis submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in the faculty of graduate studies (Mechanical Engineering) The University Of British Columbia (Vancouver) August 2011 c© Daniel Sepasi, 2011Abstract This thesis considers variations in the parameters of the dynamics of linear systems, and tackles modeling of Linear Time-Invariant (LTI) and Linear Parameter Varying (LPV) plants. The variations in the dynamics make the controller design challenging, and to successfully overcome this challenge, two methods are proposed in this thesis. One method generates a connected model set. The idea of the multi- dimensional principal curves methodology is employed to detect the nonlinear correlations between parameters of the given set of system dynamics. The connected model set is simple and tight, leading to both nonconservatism and reduced computational complexity in subsequent controller design, and hence, to improve the controller performance. The other method is developed to derive a family of discrete model sets for a given set of system response data. A relaxed version of the normalized cut methodology is developed and used in an algorithm to divide a given set of system responses into the smallest possible number of partitions in such a way that a desired performance objective is satisfied for all partitions by designing one controller for each partition. Using the proposed method, a tight uncertainty model is derived for Hard Disk Drive (HDD) systems, and an H∞ controller is synthesized. The dy- namics of HDDs is studied from a controller design point of view. Especially, the variations in the dynamics due to the change in temperatures and limited precision in the production line are examined. iiAlso, the variations in the dynamics of Ball Screw Drive (BSD) systems due to the structural flexibility, runout, and workpiece mass variation are studied. These three factors are explicitly incorporated in LPV models. To build the LPV models, it is determined how the system parameters are af- fected by two variables, namely, the measurable table position and the un- certain mass of the table. We design robust gain scheduling controllers which are scheduled by the table position and are robust over the table mass. iiiPreface This thesis entitled “Modeling of Linear Systems with Parameter Variations: Applications in Hard Disk and Ball Screw Drives” presents the research per- formed by Daniel Sepasi1. The research conducted in this thesis was super- vised by Dr. Ryozo Nagamune and co-supervised by Dr. Farrokh Sassani. In this section, we briefly explain the contents of the papers that are published or submitted for publications from this thesis [3, 4, 96–98]. We also clarify the relative contributions of co-authors in the papers. • M. Sepasi, F. Sassani, and R. Nagamune, “Parameter uncer- tainty modeling using the multi-dimensional principal curves”, Journal of Dynamic Systems, Measurement and Control, 2010, vol. 132, Issue 5, pp. 054501-054507. This paper proposes a tech- nique to model parametric uncertainties associated with linear time- invariant systems. The method is based on nonconvex optimization, involving a linear matrix inequality, a local optimization technique, and multi-dimensional principal curves. The proposed technique is ex- plained in Chapter 2 of this thesis. The author of this thesis was the principal researcher of this publication. Drs. Ryozo Nagamune and Farrokh Sassani assisted with formulating the problem and writing the paper. • E. Azadi Yazdi, M. Sepasi, F. Sassani and R. Nagamune, 1The author’s given name was changed from Mohammad to Daniel during his Ph.D. program. iv“Automated multiple robust track-following control system design in hard disk drives”, IEEE Transactions on Control System Technology, DOI: 10.1109/TCST.2010.2053541.2 This paper proposes a new design procedure for track-following control sys- tems in hard disk drives. The procedure is automated, in the sense that, for given experimental frequency response data of the suspension arm dynamics and a model structure, it automatically constructs a model set with parametric uncertainties. Subsequently, for the transfer function set it automatically designs a partition of the uncertainties and corresponding multiple robust controllers. The first step of the procedure, i.e. model set construction, is developed by the author of the thesis, and the second step, i.e. multiple robust controllers design, is developed by Dr. E. Azadi Yazdi. Experiments on actual hard disk drives demonstrate the usefulness and efficiency of the proposed proce- dure. The experiments are performed by Dr. E. Azadi Yazdi and the author of the thesis. The results of this paper is partly presented in Chapter 4 of the thesis. Drs. Ryozo Nagamune and Farrokh Sassani provided practical insight to the problem, and contributed significantly to the writing of this paper. • D. Sepasi, R. Nagamune, and F. Sassani, “Tracking control of flexible ball screw drives with runout effect and mass varia- tion”, Accepted in IEEE Transactions on Industrial Electron- ics, DOI: 10-1720-TIE.R2.3 In this paper, tracking controllers for a ball screw drive are designed, which consider flexibility and runout, as well as workpiece mass variation. The flexibility, runout, and mass vari- 2A brief version is also published in: E. Azadi Yazdi, M. Sepasi, F. Sassani and R. Nagamune, “Automated multiple robust track-following control system design in hard disk drives”, 2010 ASME Dynamic Systems and Control Conference, Boston, MA. 3A brief version is also published in: M. Sepasi, F. Sassani, and R. Nagamune, “Tracking control of flexible ball screw drives with runout effect compensation”, 2010 ASME Dynamic Systems and Control Conference, Boston, MA. vation are explicitly incorporated in a linear parameter-varying (LPV) model. To build an LPV model, it is determined through the principal curve method how the system parameters are affected by two time- varying variables, namely, the measurable position and the uncertain mass. For the LPV model, we design robust gain scheduling controllers which are scheduled by the measurable position and are robust over the uncertain mass. The designed controllers are implemented on a ball screw drive system. Drs. Ryozo Nagamune and Farrokh Sassani supervised the research and assisted with conducting the experiment and writing the paper. viTable of Contents Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii Preface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iv Table of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . vii List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xi List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xii Notation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xv Glossary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xvi Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . xvii Dedication . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xix 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Problem definition and methodology . . . . . . . . . . . . . . 3 1.2.1 Connected model set . . . . . . . . . . . . . . . . . . . 5 1.2.2 Discrete model set . . . . . . . . . . . . . . . . . . . . 10 1.3 Objectives of the thesis . . . . . . . . . . . . . . . . . . . . . . 15 1.4 Literature review . . . . . . . . . . . . . . . . . . . . . . . . . 16 1.4.1 Connected model set . . . . . . . . . . . . . . . . . . . 16 vii1.4.2 Discrete model set . . . . . . . . . . . . . . . . . . . . 20 1.5 The layout of the thesis . . . . . . . . . . . . . . . . . . . . . 22 2 Connected Model Set for a Set of System Response Data . 23 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 2.2 Connected set derivation problem . . . . . . . . . . . . . . . . 25 2.3 Connected set synthesis algorithm . . . . . . . . . . . . . . . . 26 2.3.1 General structure of the parameterizing function . . . . 26 2.3.2 Special structure of the parameterizing function . . . . 29 2.4 Examples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 2.4.1 Connected model set for a numerical example . . . . . 32 2.4.2 Uncertainty modeling of a dual-input dual-output Ball Screw Drive (BSD) . . . . . . . . . . . . . . . . . . . . 33 2.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3 Family of Discrete Model Sets for a Set of System Response Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39 3.2 Problem of deriving a family of discrete model sets . . . . . . 40 3.3 Synthesis of a family of discrete model sets . . . . . . . . . . . 41 3.3.1 Background material . . . . . . . . . . . . . . . . . . . 44 3.3.2 Parameter set partitioning and derivation of the best partition sets . . . . . . . . . . . . . . . . . . . . . . . 45 3.3.3 Optimum partition set selection . . . . . . . . . . . . . 51 3.3.4 Optimum partition set derivation algorithm . . . . . . 52 3.4 Numerical examples . . . . . . . . . . . . . . . . . . . . . . . . 53 3.4.1 Illustrative example . . . . . . . . . . . . . . . . . . . . 53 3.4.2 Closed-loop performance comparison of connected sets 55 3.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 viii4 Modeling and Robust Track-Following Controller Design for Hard Disk Drives . . . . . . . . . . . . . . . . . . . . . . . 59 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 4.2 HDD experimental setup . . . . . . . . . . . . . . . . . . . . . 61 4.3 Dynamics of HDDs . . . . . . . . . . . . . . . . . . . . . . . . 62 4.3.1 Variations in HDD dynamics due to temperature . . . 63 4.3.2 Variations in HDD dynamics due to the manufacturing limits . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63 4.4 Modeling of HDDs . . . . . . . . . . . . . . . . . . . . . . . . 64 4.5 Robust controller design for HDDs . . . . . . . . . . . . . . . 66 4.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 5 Modeling and Robust Tracking Controller Design for Flexi- ble Ball Screw Drives with Runout Effect and Mass Variation 71 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71 5.2 BSD experimental setup . . . . . . . . . . . . . . . . . . . . . 73 5.3 Variations in the dynamics of BSDs . . . . . . . . . . . . . . . 74 5.3.1 Position-dependent variations . . . . . . . . . . . . . . 75 5.3.2 Mass-dependent variations . . . . . . . . . . . . . . . . 76 5.4 Linear parameter varying model of BSDs . . . . . . . . . . . . 77 5.4.1 Linear time-invariant system identification based on frequency response . . . . . . . . . . . . . . . . . . . . 78 5.4.2 Uncertain LPV modeling . . . . . . . . . . . . . . . . . 81 5.5 Controller design for the BSD . . . . . . . . . . . . . . . . . . 86 5.5.1 Non-robust gain scheduling controller design . . . . . . 89 5.5.2 Robust gain scheduling controller design . . . . . . . . 92 5.5.3 Disturbance observer design . . . . . . . . . . . . . . . 94 5.6 Controller results . . . . . . . . . . . . . . . . . . . . . . . . . 96 5.6.1 Single controller for the BSD without mass variations . 98 5.6.2 Performance sensitivity of the BSD without mass vari- ations . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 ix5.6.3 Single controller for the BSD with mass variations . . . 100 5.6.4 Multiple controllers for the BSD with mass variations . 102 5.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 6 Conclusions, Contributions and Future Research Directions 104 6.1 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104 6.2 Summary of contributions . . . . . . . . . . . . . . . . . . . . 106 6.3 Future research directions . . . . . . . . . . . . . . . . . . . . 108 6.3.1 Uncertainty modeling for stochastic robust controller . 108 6.3.2 Performance oriented connected model set derivation . 108 6.3.3 Advanced performance oriented family of discrete model sets derivation . . . . . . . . . . . . . . . . . . . . . . . 109 6.3.4 Switching controllers for BSDs . . . . . . . . . . . . . . 109 6.3.5 BSD table mass estimation in real time . . . . . . . . . 110 Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111 Appendix A: Relaxed Form of the Optimization (3.9) . . . . . 124 Appendix B: Standard SDP Form of the Optimization (3.10) . 126 xList of Tables Table 1.1 Three special cases in this thesis. . . . . . . . . . . . . . . . 4 Table 2.1 Errors for different models . . . . . . . . . . . . . . . . . . 33 Table 2.2 Ball screw uncertainty model normalized errors for different values of nλ and m . . . . . . . . . . . . . . . . . . . . . . . 36 Table 3.1 The achieved closed-loop performance γ for different ap- proaches. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 Table 4.1 Numerical values of the coefficients of polynomial functions in (4.3). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 Table 5.1 Estimated parameters of the polynomial and sinusoidal func- tions. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 Table 5.2 Estimated parameters of the polynomial and sinusoidal func- tions for partition µ(1). . . . . . . . . . . . . . . . . . . . . 86 Table 5.3 Estimated parameters of the polynomial and sinusoidal func- tions for partition µ(2). . . . . . . . . . . . . . . . . . . . . 86 Table 5.4 Controllers tracking error results (µm). . . . . . . . . . . . 99 Table 5.5 Controllers tracking error (µm), and the percentage in- crease of errors in comparison with Krunout tracking error in Table 5.4. . . . . . . . . . . . . . . . . . . . . . . . . . . 100 Table 5.6 Tracking error results of the controllers in the third scenario.101 Table 5.7 Tracking errors and performance improvement (%) calcu- lated by 100(MAEsingle −MAEmultiple)/MAEsingle. . . . . 102 xiList of Figures Figure 1.1 Five different dynamical systems with same structure . . . 7 Figure 1.2 Model parameters (black dots) and approximating mani- folds. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 Figure 1.3 Parameter set of (1.7). . . . . . . . . . . . . . . . . . . . . 11 Figure 1.4 Mass-spring-damper closed-loop block diagram. . . . . . . 11 Figure 1.5 Violation or satisfaction of the disturbance rejection re- quirement. The dashed line shows the magnitude plot of W−1 and solid lines show sensitivity functions of the systems 13 Figure 1.6 Three ways of partitioning. . . . . . . . . . . . . . . . . . 14 Figure 2.1 Connected model set derivation. . . . . . . . . . . . . . . . 24 Figure 2.2 Visual explanations of distance notation. . . . . . . . . . . 27 Figure 2.3 A closed-loop configuration for the special case . . . . . . 30 Figure 2.4 Frequency responses of the ball screw machine (dotted line), and the estimated transfer function (solid line). . . . 35 Figure 2.5 The estimated transfer function (solid line), and the model set calculated by applying linear PCA (dashed line). The results for nλ ≥ 2 almost overlap the solid lines. . . . . . . 37 Figure 3.1 Family of discrete model set derivation. . . . . . . . . . . . 40 Figure 3.2 Main idea of the procedure to seek the optimum partition set. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 xiiFigure 3.3 Two ways of partitioning a graph. Vertices and edges are shown by dots and solid lines, respectively. The dashed lines are partition boundaries. . . . . . . . . . . . . . . . . 45 Figure 3.4 A short-form flowchart of the procedure explained in Sec- tion 3.3.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 Figure 3.5 The parameter set of the example in Section 1.2.2. . . . . 54 Figure 3.6 The normalized values of the costs of the different partition sets. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 Figure 3.7 Three cluster sets with least cost functions. . . . . . . . . 55 Figure 3.8 Two intuitive ways of clustering. . . . . . . . . . . . . . . 57 Figure 4.1 A schematic diagram of an HDD. . . . . . . . . . . . . . . 60 Figure 4.2 The solid model, which includes the arms, suspensions, and heads of an HDD. . . . . . . . . . . . . . . . . . . . . 61 Figure 4.3 HDD experimental setup. . . . . . . . . . . . . . . . . . . 62 Figure 4.4 Variations due to the changes in the temperature. . . . . . 64 Figure 4.5 10 FRF data (2 FRF are taken for each HDD). . . . . . . 65 Figure 4.6 10 random samples from the connected model set. . . . . . 67 Figure 4.7 Frequency-domain response of closed-loop systems. . . . . 69 Figure 4.8 Time-domain response of open-loop and closed-loop systems. 70 Figure 5.1 The ball screw experimental setup. . . . . . . . . . . . . . 74 Figure 5.2 A schematic diagram of a ball screw assembly. . . . . . . . 75 Figure 5.3 Run out effect in the ball screw shaft. . . . . . . . . . . . . 76 Figure 5.4 Frequency responses for 37 different positions of the table along the shaft when no mass is added to the table. . . . . 76 Figure 5.5 Frequency responses for four different masses added to the table when the position is at 0.25 m. . . . . . . . . . . . . 77 Figure 5.6 Estimated transfer function parameters. . . . . . . . . . . 80 xiiiFigure 5.7 Transfer function parameters (solid lines), values in the estimated single model (dash lines), boundary of partitions (vertical dotted lines). . . . . . . . . . . . . . . . . . . . . 84 Figure 5.8 Transfer function parameters (solid lines) and values in the estimated multi model (dash lines). . . . . . . . . . . . . . 87 Figure 5.9 The closed-loop block diagram. . . . . . . . . . . . . . . . 89 Figure 5.10 Synthesis closed-loop configurations for the case that the uncertainty is ignored in the plant G. . . . . . . . . . . . . 89 Figure 5.11 Synthesis closed-loop configurations for the uncertain plant G. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 Figure 5.12 Controller with disturbance observer scheme. . . . . . . . . 96 Figure 5.13 Frequency responses of perturbed T (s). . . . . . . . . . . . 97 Figure 5.14 Tracking errors and control inputs. . . . . . . . . . . . . . 99 Figure 5.15 Tracking errors of robust and non-robust Krunout controller for different added masses to the table. . . . . . . . . . . . 101 xivNotation Symbol Description 1n Vector with n elements of one 0n Vector with n elements of zero R Set of all real numbers Rn Set of all n-dimensional real vectors Rn×m Set of all n×m-dimensional real matrices C Set of all complex numbers Tr(A) The trace of a square matrix A A  0 Positive semidefiniteness A ≥ 0 Element-wise non-negativity Sn Set of all n-dimensional symmetric matrices xvGlossary BSD Ball Screw Drive CNC Computer Numerically Controlled DM Diffusion Map DOB Disturbance Observer FFT Fast Fourier Transform FRF Frequency Response Function HDD Hard Disk Drive LDV Laser Doppler Vibrometer LFT Linear Fractional Transformation LLE Local Linear Embedding LMI Linear Matrix Inequality LPV Linear Parameter Varying LTI Linear Time-Invariant LTV Linear Time-Varying MAE Mean Absolute Error MDS Multidimensional Scaling MSD Mass Spring Damper NCut Normalized Cut NLPCA Nonlinear Principal Component Analysis PES Position Error Signal SDP Semidefinite Programming VCM Voice Coil Motor xviAcknowledgments I would like to express my deepest appreciation to my advisors, Drs. Ryozo Nagamune and Farrokh Sassani, for their inspiration, friendship, encourage- ment, patience and unconditional support. I am extremely fortunate for having the opportunity to work with them and having the freedom to ex- plore science in their laboratory. I am really thankful for their support over the last years. I am indebted to Drs. Yusuf Altintas and Mu Chiao for sharing his lab equipment. I would like to convey my gratitude to our kind graduate secretary of Mechanical Engineering Department, Ms. Yuki Matsumura. In addition, I would like to extend my sincere thanks to Mr. Glenn Jolly, for his essential technical assistance during this research. Additionally, I would like to thank Mr. Perry Yabuno for timely acquisition of electro-mechanical devices required for the fabrication of the set-up. Additionally, I am grateful to my dear colleagues in the Control Engi- neering Laboratory: Dr. Ehsan Azadi, Mr. Marius Postma, and Mr. Massih Hanifzadegan. Their friendship and support made the lab a motivating envi- ronment during the period of my work. I am so lucky to have such supportive friends: Vahid Bazargan, Pirooz Darabi, Saghar Mohajeri, Sina Radmard, Amir Rasuli, Behnam Razavi, and Hamidreza Yamini. Lastly, I wish to express my genuine gratitude to my wonderful family for their never-ending love, support, and guidance. My fathers perfectionism and seek for the truth and my mothers energy, encouragement and moral support xviihave always been a great motivation for me. I also thank my brother, my durable source of inspiration, for his support in every step of my life. xviiiDedication To my family, an insufficient token of my appreciation of their unwavering love and faithfulness xixChapter 1 Introduction 1.1 Motivation To achieve satisfactory performance for control systems during their opera- tion, controllers must be designed for any conceivable situation. In different operating conditions, the dynamics of the system varies, and these varia- tions need to be taken into account in high performance controller design. In many servo systems, such as Hard Disk Drives (HDD) and Ball Screw Drives (BSD), it is often the case that the governing dynamic equations do not vary, but there can be variations associated with parameters in the equations. In this thesis, we assume that uncertain and scheduling variables cause variations in the parameters of the plant dynamics with fixed structures. In other words, the parameters are functions of these two types of variables. An uncertain variable is a parameter whose value is unknown and unmeasurable during the operations, but its variation is known to be bounded within a specific range. A varying parameter is called a scheduling variable if its value is available in real time while the dynamics varies. Such value can be either measured or estimated, and it can be used for scheduling controllers. Since parameter variations affect the system dynamics, they need to be taken into consideration in controller design. A well-established theoretical 1tool to deal with the variations in plant dynamics is robust control [120]. Using the robust control theory, one can design controllers that guarantee robust stability and performance for a model set with uncertainties. To avoid unnecessary conservatism which is inherent in robust control system design, it is of interest to derive a tight and accurate model set [13, 53, 69]. Moreover, in many applications, the variations of the system dynamics depend on the scheduling variables. The scheduling variables can be used to improve the closed-loop performance by adjusting controller parameters. To effectively address these issues in controller design, adaptive controllers [43, 108] have been utilized. We focus on an adaptive control method using the gain scheduling approach [91]. To successfully design a gain scheduling controller, it is essential to estimate the relations between the parameters of the system dynamics and the scheduling variables. One servo system that is studied in this thesis is an HDD. One of the most important characteristics determining the quality of an HDD is the areal storage density. It is essential to reduce the tracking error in order to increase the areal density of HDDs. To achieve a desired performance we require a precise model of the system dynamics, based on which the position control of the read/write head is designed [15, 33, 45, 70]. The dynamics can vary due to many factors such as variations in the fabrication environment, temperature change, and mechanical imperfections due to the elapse of time. Such variations are studied in detail in Chapter 4. Another servo system investigated in this thesis is a BSD. BSDs are mostly used for high precision motion applications, such as in CNC ma- chines and wire bonding. In these applications, the objective is to accurately position the workpiece relative to the tool. The quality of the machining product depends greatly on the tracking performance of the machine over a desired trajectory for the workpiece position. In order to achieve a satisfac- tory tracking performance during the operations, servo controllers must be designed to take into account any possible situation during the operations. 2Therefore, it is critical to derive a model that precisely presents the BSD dynamics, which has variations. These variations occur due to many un- avoidable factors such as changes in the BSD table position and mass during operations. Such variations are studied in detail in Chapter 5. This thesis develops two modeling methodologies to tackle the challenges caused by the variations in the plant dynamics. The effectiveness of the proposed methodologies are investigated using HDD and BSD plants. The dynamics of these two servo systems are also studied in detail. 1.2 Problem definition and methodology The following is a modeling problem from a controller design point of view. Here, the problem is stated in a very general form. Problem 1.2.1. For a given system response data set with the same govern- ing dynamic equation but possibly different parameter values, find a model set, which represents the given data set accurately, in such a way that a satisfac- tory closed-loop performance can be achieved by designing the corresponding controller set. Several special cases of the above problem have been addressed in the literature (See Section 1.4). Here are a number of special cases regarding • system response data: D1. the given system response data is in the frequency domain. D2. the given system response data is in the time domain. • closed-loop performance objective: P1. there is no given desired controller objective. However, the model is derived to be tight to avoid unnecessary conservatism inher- ent to any controller, which leads to the closed-loop performance enhancement. 3Table 1.1: Three special cases in this thesis. P1 P2 M1 Chapter 2 Future work M2 − Chapter 3 P2. a desired controller objective is given, and the model set is derived in such a way that this objective is satisfied for all the given system response data by designing the corresponding controller set. • the model set: M1. the derived model set is connected, i.e., not only systems in the given set but also “intermediate” systems are considered. M2. the derived model set is discrete, i.e., only systems in the given set are considered. In this thesis, the special cases are addressed as shown in Table 1.1. Two methodologies are developed to address two special cases. One method [97] considers the special cases M1, P1, D1, and D2. This method generates a tight connected model set to avoid unnecessary conservatism. The other method considers the special cases M2, P2, D1, and D2. This method derives a family of discrete model sets, in such a way that a given satisfactory closed- loop performance can be achieved by designing a family of robust controllers, each of which is in charge of one discrete model set. Both of these methods can deal with the time and frequency responses of systems. This section briefly describes the methodology, provides a simple exam- ple, and explains problem formulation for each of these methods deriving connected and discrete model sets. Here, to prevent confusing notations and excessive introductory material, the modeling problems are stated in an in- formal form. The problems are revisited and reformulated precisely in later 4chapters of the thesis. 1.2.1 Connected model set To understand the variations in the dynamics of a large number of systems with a common structure, we often take a number of samples, study their dynamics, and deduce a model representing not only the samples but also intermediate systems between the sampled ones. Here are two scenarios to clarify these applications, when the variations are caused by scheduling and uncertain variables. 1. Assume that the dynamics of a plant varies by changing the tempera- ture during operations. For the controller design, we need to derive a model representing the plant dynamics for the entire considered range of temperature. To this end, we can take system responses for some temperature samples, and identify the system based on these samples. However, the temperature changes continuously in reality, and there- fore, we need to drive a connected model set which approximately in- terpolates the dynamics of the sampled models. In the case where the temperature is measurable, it is considered as a scheduling variable, and the system is Linear Parameter Varying (LPV). 2. Suppose that we want to design a controller for products of a production line. We take system responses of some product samples to study the dynamics. These plants are considered as Linear Time-Invariant (LTI). Due to some factors, e.g., limited tolerances in the production line, there are variations associated with the dynamics of these samples. It is desired to derive a connected model set to cover these variations and represent the sampled products as well as the unsampled ones. Since we often can not measure nor estimate the sources of variations, we deal with an uncertainty modeling problem. 5In both scenarios, we need to model parameter variations. Such modeling is important for designing robust and gain scheduling controllers. A smaller model set generally leads to controllers with less conservative performances. Therefore, it is of interest to derive a set of models with the following char- acteristics: • it precisely represents each plant of the given sample set. • it is a connected set since it is also required to represent the interme- diate plants between the sampled ones. • a small number of independent variables are used to parameterize the set to simplify the model set expression, leading to both nonconser- vatism and reduced computational complexity in the subsequent con- troller design. A simple numerical example is studied to illustrate the advantages of vari- ations modeling in controller design. Assume there are five samples whose dynamics are governed by a common transfer function in the frequency do- main as [G(θ∗)](s) = b0(s+ b1)(s+ b2) , θ ∗ := [b0, b1, b2]T , (1.1) and the parameters are θ∗ =  0.022 3  +   0.01 00 1.5 2 0   [ λ1 λ2 ] , (1.2) where λ1 and λ2 are two independent varying terms,[ λ1 λ2 ] = {[ 1 1 ] , [ 1 −1 ] , [ −1 1 ] , [ −1 −1 ] , [ 0.1 0 ]} . (1.3) Figure 1.1 shows the Frequency Response Function (FRF) data of these sys- tems. 610−1 100 101 102 −180 −90 0 Frequency (rad/sec) Phase (deg ) −150 −100 −50 0 Magnitude (dB ) Figure 1.1: Five different dynamical systems with same structure Now, assume that the FRF data of the samples, depicted in Figure 1.1, are given without any knowledge about the information in (1.1), (1.2), and (1.3). The objective is to derive a connected model set to estimate all samples as well as the unsampled systems. First, we estimate the order of a proper rational transfer function to retrieve the structure in (1.1). The order can be selected by systematic methods such as the Akaike’s information criterion [1], inspection of the FRF data, or by trial and error. In this example, the order of the transfer function can be easily estimated by inspection, [G(θ)](s) = a0 s2 + a1s+ a2 , θ := [ a0, a1, a2 ]T . (1.4) Then, using the least-squares method, such as commands invfreqs.m or fitfrd.m in the Matlab software, the set Θ consisting of all estimated pa- rameter vectors is computed as Θ := {θ`}5`=1 =     0.038.5 17.5   ,  0.035.5 2.5   ,  0.014.5 3.5   ,  0.011.5 0.5   ,  0.0215.2 6.4      . (1.5) 7Figure 1.2(a) illustrates the parameter set, in which each vector in Θ is shown by a dot. In this figure, the 3-dimensional grey box lies in between the minimum and the maximum values of all the parameters. Such a box covers all the parameter variations and therefore, a controller can be designed for the same region. On the other hand, a manifold with fewer dimensions, one or two here, that can approximate all vectors in Θ may exist. Subsequently, a controller designed for this manifold can achieve a better performance in comparison with the case that the controller covers the entire grey box. We have applied Principal Component Analysis (PCA) to find a lin- ear manifold approximating the parameter vectors. The resultant plane is shown in Figure 1.2(b). For this simple example, an exact but unknown 2-dimensional manifold is obtained from (1.1), (1.2) , and (1.4) as a0a1 a2   =  0.025 6  +   0.012 4  λ1 +   01.5 4.5  λ2 +   00 3  λ1λ2. (1.6) Notice that the manifold in (1.6) is nonlinear in variable terms, and therefore, linear methods, such as PCA, can not estimate the parameter set as precisely as nonlinear ones, such as Nonlinear PCA (NLPCA). The goal is to detect a manifold, which interpolates the samples to cover the variations associated with not only the samples but also the interme- diate plant dynamics. Detecting such a manifold is essential to reduce the conservativeness of the controllers, and hence, to improve the performance. By generalizing the simple example above, we now formulate a problem of constructing a connected model set. Problem 1.2.2. Assume a set of transfer functions with a common struc- ture, but possibly different parameter values, is given. The goal is to derive a connected set of transfer functions, which represents each element of the given set accurately. The connected model set should be tight1, leading to 1For instance, “tight model set” refers to a model set with few dimensions in the 8(a) The gray box covers all the variations in a conservative way. (b) The gray plane shows a manifold approxi- mating the parameter set as calculated by linear PCA. Figure 1.2: Model parameters (black dots) and approximating mani- folds. 9nonconservatism in subsequent controller design. 1.2.2 Discrete model set When large parameter variations exist, it may be infeasible to design a sin- gle robust controller to satisfy a desired performance objective for the entire range of variations. One way to overcome this infeasibility is to divide the parameter set into a finite number of partitions, and design one controller for each partition. Since each controller deals with smaller parameter varia- tions, the performance can be enhanced. Generally, the overall performance strongly depends on the characteristics of the partitions, i.e., their size and the way in which partitions are separated. Therefore, it is of interest to di- vide the parameter set optimally, and derive a family of discrete model sets with the following characteristics: • it represents each plant of the given sample set by one element. • a given performance objective can be satisfied for all partitions by de- signing one controller for each partition. As a motivating example, let us consider a control problem for a set G including L = 12 mass-spring-damper (MSD) systems with force input and displacement output. The dynamics of such systems can be represented in a normalized form in the frequency domain as G`(s) := ω 2 n,` s2 + 2ζ`ωn,`s+ ω2n,` , ` = 1, . . . , L, (1.7) where ωn and ζ are modal parameters. The corresponding parameter set Θ := {(ωn,`, ζ`)}L`=1 is shown in Figure 1.3. The control objective is to reject the displacement disturbance for all the systems in G using feedback control. The disturbance rejection requirement previous example. 101 1.5 2 0.1 0.15 0.2 ω n ζ Figure 1.3: Parameter set of (1.7). Figure 1.4: Mass-spring-damper closed-loop block diagram. for a system G and a controller K can be expressed in terms of the sensitivity function as ΓMSD(G,K) := ∥∥∥∥ W 1 +GK ∥∥∥∥ ∞ < 1, (1.8) for a given weighting function W . In this example, we assume the weighting function is2 W (s) = 0.2s+ 391.9 s+ 0.3919 . (1.9) Figure 1.4 shows the closed-loop block diagram, where d is the displacement disturbance signal. One straightforward approach to controller design for multiple systems is to design one stabilizing feedback controller K` for each system G`. However, this approach can be demanding, especially for a large number of systems, because it leads to time-consuming design of a large number of controllers. Therefore, our goal is to design a small number of controllers, in such a way 2The weighting function is a design specification which is given based on the perfor- mance requirement [120]. 11that each system is controlled by one of the designed controllers. In this example, we want to design the smallest possible number of con- trollers. First, we try to design a single controller by solving the following optimization problem, γ := min K∈K max `=1,...,L ΓMSD(G`, K), (1.10) where K shows a set of controllers, which robustly stabilize all G`. Here, we follow an output-feedback controller design procedure proposed in [50], which is based on the Linear Matrix Inequality (LMI) technique. The optimum γ that is achieved by applying only one controller is 2.172. Since γ > 1, the disturbance rejection requirement (1.8) is violated for some plants, see Figure 1.5(a)3. One way to possibly fulfill this requirement is to divide the set G into two partitions G(1) and G(2), where G = G(1) ∪ G(2) and G(1) ∩ G(2) = Ø. Then, solving the optimization problem (1.10) for each partition, the optimum γ1 and γ2 are achieved. The value of γopt := max{γ1, γ2} determines if the control objective is satisfied. Since the dynamics of systems are governed by a common transfer function, partitions in parameter domain {Θ(q)}2q=1 are equivalent to {G(q)}2q=1, where Θ(q) includes parameter vectors of systems in G(q). Two intuitive ways of partitioning Θ, among many others, are shown in Figures 1.6(a) and 1.6(b). Here, elements of different partitions are shown by different types of markers. The optimum values of γopt are 1.243 and 1.045 respectively for partition sets shown in Figures 1.6(a) and 1.6(b), and thus, these two partition sets violate the requirement. On the other hand, a parti- tion set shown in Figure 1.6(c) provides γopt = 0.759 < 1, see Figure 1.5(b). In fact, this value of γopt is the global minimum for the case of two parti- tions. This optimal partition set (Figure 1.6(c)) is derived by “full search”, 3Based on the requirement (1.8), the magnitude of W−1 is compared with those of the sensitivity functions to examine the violation 1210−1 101 103 105 −150 −100 −50 0 50 Frequency (rad/sec) Magnitude (dB ) (a) One controller is designed for all sys- tems, and the performance objective is not satisfied for some systems, γopt > 1 10−1 101 103 105 −150 −100 −50 0 50 Frequency (rad/sec) Magnitude (dB ) (b) Two controllers are designed for an optimum partition set, and the perfor- mance objective is satisfied for all systems, γopt < 1 Figure 1.5: Violation or satisfaction of the disturbance rejection re- quirement. The dashed line shows the magnitude plot of W−1 and solid lines show sensitivity functions of the systems i.e., assessing the closed-loop performance of all 2509 possible partition sets, combinations of twelve points in two partitions. In many practical applications, the parameter set can be large in terms of dimension (more complex transfer function) and the number of elements (more systems to be controlled). In such cases, the “full search” might be- 131 1.5 2 0.1 0.15 0.2 ω n ζ (a) A partition set lead- ing to γopt = 1.242 1 1.5 2 0.1 0.15 0.2 ω n ζ (b) A partition set lead- ing to γopt = 1.045 1 1.5 2 0.1 0.15 0.2 ω n ζ (c) The global optimum partition set leading to γopt = 0.759 Figure 1.6: Three ways of partitioning. come impractical even for a small number of partitions. The goal is to obtain an optimum partition set in a systematic manner to satisfy the closed-loop performance objective. By generalizing the simple example above, we now formulate a problem of constructing a family of discrete model set. Problem 1.2.3. Assume a set of transfer functions and a desired closed- loop performance are given. The transfer functions have a common structure with possibly different parameter values. The goal is to divide the given set into the smallest possible number of partitions in such a way that the desired performance objective is satisfied for all partitions by designing one controller for each partition. 141.3 Objectives of the thesis Both of the connected and discrete model set identification problems have been previously addressed in the literature. Nevertheless, the solutions pro- posed previously often suffer from a lack of generality and/or conservatism. Therefore, the solutions to these problem are still being developed by re- searchers. Moreover, the dynamics of HDDs and BSDs have been thoroughly studied and reported in the literature. However, there is still a lack of attention to the dynamics of these systems from a controller design point of view. The objectives of this thesis are to O1. develop a modeling technique to derive connected model sets for Linear Parameter Varying (LPV) and Linear Time-Invariant (LTI) systems, O2. develop a systematic algorithm to derive a family of discrete model sets, in such a way that a desired closed-loop performance can be achieved by designing the corresponding controller set, O3. study the dynamics of HDDs from a controller design point of view, O4. model BSDs and design a robust gain scheduling controller to precisely control the position of the machine table. In this thesis, general linear plant dynamics are considered. The sys- tems are general in a sense that they can be single-input single-output or multi-input multi-output, and discrete-time or continuous-time. The pro- posed methods will be compared to existing methods available and reported in the literature. The developed methods are also numerically and/or ex- perimentally validated, using typical benchmark problems from the control literature, such as mass-spring-damper systems, as well as practical control problems, namely, the hard disk drive and ball screw drive servo systems. 151.4 Literature review 1.4.1 Connected model set As discussed earlier, we may need to derive a connected model set to cover all the parameter variations. There have been some investigations into pa- rameter variations modeling in various systems to generate connected model sets [11, 14, 52]. Most of the proposed methods assume that how each in- dependent variable affects the parameters in the transfer functions is known [74, 76]. This assumption is reasonable when we know the physical laws governing the system dynamics, e.g., the structures of (1.1) and (1.2) in the motivating example. However, this assumption is not valid in general. The parameter variation modeling problem that we consider is closely related to the parameter reduction. Hence, we review the literature in the area of parameter reduction, and utilize the related best methodology as the basis for modeling formulation. Parameter dimension reduction is the transformation of a high-dimensional parameter set into a meaningful representation of reduced dimensionality. This representation has a dimensionality that corresponds to the “intrinsic dimensionality” of the parameter set. The intrinsic dimensionality of the pa- rameter set is the minimum number of parameters needed to account for the observed properties of the model [34] without losing any critical information about the plant. To this end, the possible correlation of original parameters must be identified. In some applications, parameters with little contribution to the input- output behavior can be ignored in the estimation process [64, 116]. One application of this method of parameter reduction can be found in [24] where PCA and sensitivity analysis were used to reduce the number of parameters in a model representing a complex metabolic network. Sun and Hahn [101] applied three techniques to reduce the parameter set of fundamental models. They extended the described methods to nonlinear systems. However, in 16many cases including our application, i.e., high precision controller design, it is generally undesirable to neglect any parameter of the system model. Dontchev et al. [27] proposed a numerical method for reduction of a priori bounds on the values of uncertain plant parameters. The method successively eliminates parts of the parameter domain which are inconsistent with the plant measurement. The goal is to achieve the smallest possible area which contains all the information about the uncertainty of the parameter set. This method is computationally expensive. A large number of nonlinear techniques for dimensionality reduction have been proposed [7, 26, 41, 61, 105]. These techniques have the ability to deal with complex nonlinear data. One can divide such techniques into supervised and unsupervised learning methods. Most commonly, supervised learning generates a global model that maps inputs to desired outputs. This means that the data set consists of pairs of input objects and desired outputs. Many well-known techniques for su- pervised learning methodology are based on the linear discriminant analysis [31]. Since, in our applications, the data set is a parameter set and not a set of input-output data, the supervised learning method is not suitable here. On the other hand, the modeling problem of the parameter variations can be solved by implementing an unsupervised learning technique such as PCA [42]. Nonlinear techniques for parameter reduction can also be divided into two main categories: • Local techniques. These techniques attempts to preserve local prop- erties of the original data set. Four most common methodologies in this category are, Local Linear Embedding (LLE) [90], Laplacian eigenmaps [7], Hessian LLE [26], and local tangent space analysis [119]. There are two main disadvantages associated with this type: (a) Local properties do not necessarily follow the global structure, as noted in [9, 89], specially in the presence of noise. (b) Since the distribution of 17the original data set is not necessarily uniform, the neighbor selection must be adaptive; otherwise the performance of parameter reduction would be degraded. • Global techniques. These techniques attempt to preserve global properties of the original data set. Four well-known global nonlinear techniques are, Multidimensional Scaling (MDS) [20], Isomap [104], Diffusion Map (DM) [58], and NLPCA [55]. – MDS maps the original data to a low-dimensional representation while preserving the distance between the data points in a pairwise fashion. The quality of the MDS is expressed in the stress function ∑ i,j (‖xi − xj‖ − ‖yi − yj‖)2, (1.11) where ‖xi − xj‖ and ‖yi − yj‖ are the Euclidean distances be- tween the high-dimensional and low-dimensional data points re- spectively. By modifying its cost function, one may put more em- phasis on preserving distances which were originally small, such as Sammon cost function [92]. Although it is well known that MDS has been successful in many applications, it does not take into consideration the distribution of the neighboring points, and it is based on Euclidean distance. Therefore, MDS for data near to a curved manifold, such as “Swiss roll” data set, does not perform satisfactorily. – Isomap is a methodology that resolves the drawback of MDS by attempting to retain a curvilinear distance between data points in- stead of a Euclidean distance. Curvilinear distance is the distance between points measured over the manifold. – DM first constructs a graph of original data, x ∈ X. The weights of the edges in the graph are computed using a Gaussian kernel 18function, leading to a matrix with entries, pij. This matrix is used to calculate the diffusion distance: D(xi, xj) := ∑ k (pik − pjk)2 ψ(xk) , (1.12) where ψ(xk) is a term attributing more weight to part of the graph with high density. The main idea is that the DM is based on many paths through the graph. This makes the DM more robust than, e.g., the methods based on the curvilinear distance, such as the Isomap. – NLPCA is an extension of PCA. Traditional linear PCA is a data analysis technique identifying patterns and expressing the data with independent variables in lower dimensional space. In other words, if vector θ ∈ Rnθ represents the observations, there is a transformation matrix P to produce independent variables stored in λ ∈ Rnλ , nλ ≤ nθ, λ = P T θ. (1.13) There are some nonlinear extension of PCA that have been pro- posed over the past two decades. Studies on NLPCA can be di- vided into the utilization of neural network [54], principal curves [38], and kernel approaches [94, 95]. Although Isomap and DM show a capability of outperforming some other techniques to reduce the dimension of data, in uncertainty modeling aspect they are not so powerful as NLPCA. Therefore, we propose an optimization problem based on the concept of PCA. Here, we prefer to adopt the principal curves method over neural network and kernel PCA for the following reasons: - The computational cost of neural network formulations increases consid- erably if the number of observations (here, the number of elements in 19the parameter set) rises. On the other hand, based on some literature such as [25], principal curve computation is not normally subject to this computational concern. - Kernel PCA performs PCA in a feature space of arbitrarily large (pos- sibly infinite) dimensionality. The size of the kernel matrix increases quadratically with an increase in the number of samples. Also, com- pared to principal curves, kernel PCA is harder to interpret in the input space [57]. Principal curves presented by Hastie and Stuetzle [37, 38] are smooth one-dimensional curves that pass through the middle4 of multi-dimensional distributions or data sets. The shape of the principal curve is determined by the structure of the data set, and it provides a nonlinear summary of the data. The principal curve is formally defined to be smooth self-consistent curve for a data set. In other words, any point on the curve is the average of all of the data which project onto that point. Although, in the original definition, the principal curve is defined as a one-dimensional curve, in this study the concept of the principal curve is extended to develop a multi-dimensional curve. 1.4.2 Discrete model set As discussed above, the variations in the system dynamics increase the size of the model set, and consequently degrade the closed-loop performance, and thus need to be taken into consideration. In many applications, achieving the desired performance of the controller can be demanding due to the size of the model set [120]. There are two main approaches to overcome this challenging problem. The more common approach solves the problem by modifying the con- troller synthesis methodology for a given model set. Depending on the essence 4The interpretation of “middle” is given in Section 2.3.1. 20of the variations in the model set, different controller synthesis methodologies are available, such as the switching ([72, 115]), multiple robust ([18, 117]), adaptive ([32, 88]), and model predictive controllers [29]. The common idea of these methods is to design a number of controllers, such that, each of which covers one part of the variations. However, these methods lead to a more conservative closed-loop system in comparison with the case where a multiple-model set is derived based on the closed-loop performance. The other approach attempts to derive a multiple-model set. This ap- proach has a history of two generations [63, 65]. The first generation was rep- resented by Magill [73] and Lainiotis [60]. Blom and Bar-Shaloms pioneered an interacting multiple-model algorithm [8] and introduced the second gener- ation. The interacting multiple-model has earned an enviable reputation for multiple-model estimation via a significant number of successful applications (see, e.g., [6]). Since then, different aspects of multiple-model set derivation have been investigated [66, 87]. All the above literature in the multiple-model approach considers switch- ing between local models. To the best of our knowledge, there is no literature on derivation of a family of discrete model sets based on the desired closed- loop performance, as defined in the Problem 3.2.1. In cases where the system variations are only due to the parameter vari- ations, the parameter set can be divided into a finite number of partitions. Then, one controller is designed for each partition, and since it deals with less parametric variations the performance can be enhanced. Generally, the overall performance strongly depends on the characteristics of the partitions, i.e., their size and the way that they are formed. Therefore, it is of interest to divide the parameter set optimally. In general, dividing a given parameter set into partitions is based on the “similarity” among the parameters. A comprehensive introduction and survey in clustering can be found in [47]. Classical algorithms optimize sim- ple objectives, for instance K-means minimizes the spread over centroids. 21However, due to their simplifying assumptions about partition structure and intuitive interpretations, they may provide poor clustering solutions. One well-recognized method for clustering is the Normalized Cut (NCut) [99]. However, it involves an optimization with a nonlinear objective and a combinatorial nature of the feasible set, which leads to an NP-hard problem. Spectral relaxation is a legitimate approach with significant practical suc- cesses [99] for relaxing the NCut optimization. Nevertheless, Guattery and Miller discussed some valid drawbacks of spectral relaxation [35]. Semidefi- nite programming (SDP) has been powerful in approximating similar difficult clustering optimization [51, 113], which results in a tighter relaxation of NCut in comparison with the spectral relaxation [122]. In this thesis, a relaxed version of NCut is applied in an algorithm to achieve an optimum set of partitions in such a way that a given performance objective is satisfied for all partitions by designing one controller for each partition. Further details about this approach is provided in Chapter 3. 1.5 The layout of the thesis This thesis is organized as follows. Chapter 2 is devoted to Objective (O1) of the thesis, presented in Section 1.3. To be more specific, this chapter proposes a method to derive connected model sets. In Chapter 3 an algorithm is developed to achieve Objective (O2). The dynamics of HDD read/write head positioning systems is investigated in Chapter 4, where Objective (O3) is accomplished. The proposed method in Chapter 2 is then used to model the dynamics of HDD systems to design track-following controllers. The variations in the dynamics of BSDs are studied in Chapter 5 to achieve Objective (O4). Three main sources of variations are considered in tracking controller design. Finally, Chapter 6 concludes the thesis, summarizes its major contributions, and provides directions for future research. 22Chapter 2 Connected Model Set for a Set of System Response Data 1 2.1 Introduction This chapter proposes a systematic method to obtain a connected model set for a given set of system response data. It is assumed that the govern- ing equations of these system are the same but associated with parameter variations. This method achieves Objective (O1) in Section 1.3. The proposed practical procedure of this method is briefly illustrated in Figure 2.1, and is described with mathematical expressions in this chapter. This algorithm develops a mapping from experimental system response data of a batch of sampled plants to a connected model set. The resultant set is a set of real rational transfer functions, of which parameters are parameterized by a small number of uncorrelated variables that capture the differences in the dynamics of the sampled plants. The transfer function structure can be selected by systematic methods such as the Akaike’s information criterion [1], inspection of the FRF data, 1Most of this chapter is based on the following publication: M. Sepasi, F. Sassani, and R. Nagamune, “Parameter uncertainty modeling using the multi-dimensional principal curves”, Journal of Dynamic Systems, Measurement, and Control, 2010, Vol. 132, No. 5, pp. 054501-054508. 23Figure 2.1: Connected model set derivation. or by trial and error. The corresponding parameter set can be estimated through many powerful tools, such as the System Identification [68] and Signal Processing [67] Toolboxes in the Matlab software. Therefore, the pro- cedure of obtaining the parameter set from the system responses is ignored in this chapter, and in the problem formulation, it is assumed that the pa- rameter set is given. After deriving a connected parameter set, we generate the corresponding connected model set, which can be simply done since the transfer function structure is available. Therefore, in the problem formulation the goal is to obtain a connected parameter set. An example is explained in Section 2.4.2, in which the connected model set is derived for a general case where we deal with both uncertain and scheduling variables. This chapter is organized as follows. Section 2.2 formulates a connected model set derivation problem. In Section 2.3, a method is proposed to solve the formulated problem. Section 2.4 provides numerical and practical exam- ples to validate the proposed method. 242.2 Connected set derivation problem In Section 1.2.1, we have formulated a problem of generating a connected model set in a general term; See Problem 1.2.2. Here, using mathematical notation, we will reformulate it more rigorously. Problem 2.2.1. Assume a set of L parameter vectors of transfer functions with a common structure is given as Θ := {θ` ∈ Rnθ , ` = 1, . . . , L}. (2.1) The goal is to derive a tight connected parameter set ̂Θ(f) := {f(λ) ∈ Rnθ , λ ∈ Rnλ : nλ ≤ nθ}, (2.2) which represents each element of the given set Θ. The function f is the parameterizing operator, and λ represents uncorrelated variables. Remark. There are three remarks regarding the above problem; 1. The parameter λ can include uncertain and scheduling independent variables. 2. A small number of uncorrelated variables λ are used to parameterize the set Θ to obtain a tight model set, leading to both nonconservatism and reduced computational complexity in subsequent controller design. 3. The number of independent variables λ as well as the structure of the function f are assumed to be known. However, if they are unknown, they can be obtained by trial and error. See the example explained in Section 2.4.1. 252.3 Connected set synthesis algorithm The goal of this algorithm is to derive the set ̂Θ in (2.2), which interpo- lates L members of Θ in (2.1), and represent the intermediate plants. Also, it is desired that the resultant set ̂Θ has the smallest possible size2. To this end, possible correlations among parameter vectors are determined by parametrization with a small number of uncorrelated variables. In the following, first we explain the synthesis algorithm for a general function f . Then, in Section 2.3.2, a special case of this function is studied. 2.3.1 General structure of the parameterizing function Parameters λ represent the sources of variations, which can be scheduling variables (e.g. temperature), denoted by λs Λs := {λs ∈ Rns}, (2.3) or uncertain variables (e.g. product differences), denoted by λu which are normalized3, Λu := {λu ∈ Rnu : λu ∈ [−1, 1]}, (2.4) or a combination of both. Here, we develop the methodology based on the principal curves for the most general case, which takes into account both scheduling and uncertain 2Our intention of “the size of the parameter set” was clarified through an example in Section 1.2.1. In different problems, the word “size” can refer to the dimension, area, etc. 3It is common in robust controller design methods to normalize the uncertain variables. Therefore, we normalize the uncertain variables λu, contrary to the scheduling variables λs. 26Figure 2.2: Visual explanations of distance notation. variables Λ := Λs ⊕ Λu = {λ ∈ Rnλ : nλ = ns + nu}. (2.5) The distance between ̂Θ and a sample parameter vector θ` in Θ is defined by d(̂Θ(f), θ`) := min λ ‖f(λ) − θ`‖, (2.6) with a proper norm ‖ · ‖. See the visual explanations of the notation in Figure 2.2. Such distance is desired to be minimized for all L samples. The function f which yields the minimum distance is obtained by a minimax optimization problem min f∈F max `=1,...,L d(̂Θ(f), θ`), (2.7) where F represents the given class of functions. Remark. Here, the accuracy of the model is equally important for all samples. For example, instead of optimization (2.7), if we solve min f L∑ `=1 d(̂Θ(f), θ`) (2.8) to derive the function f , we may face a situation that the optimum f ∗ pro- vides a model that has a large error for one sample θ`. Since the error for 27other samples are small, the cost function of the optimization problem (2.8) is minimized. However, this result may not be acceptable from system model- ing point of view. Therefore, we prefer the minimax optimization formulation rather than minimization of other cost functions, such as mean absolute error and mean square error. Since the optimization problem (2.7) is nonconvex in general, we propose an iterative procedure to solve it as follows. Algorithm 2.3.1. Inputs: The set Θ = {θ` ∈ Rnθ ; ` = 1, . . . , L}, the structure of the func- tion f , and values of {λs`}L`=1, nu and ns. 1. Initialize the function f . 2. Discretize the set Λu by selecting N values for λu. Then, solve the optimization (2.6) for discretized λun min λun ‖f( [ λun λs` ] ) − θ`‖, for n = 1, . . . , N, (2.9) for each ` and denote the solution by λu∗` . 3. Solve the optimization (2.7) min f max `=1,...,L ‖f( [ λu∗` λs` ] ) − θ`‖ (2.10) to update the function f . 4. Iterate Steps 2 and Step 3 until change in optimal values between iter- ations is small or a satisfactory error is achieved. The initialization of the function f in the first step of the above algorithm can be done by any suitable algorithm. We suggest to use linear PCA, which provides the principal directions of correlation between the parameters. 28Depending on the structure of the function f , the third step can be a convex or nonconvex optimization problem. 2.3.2 Special structure of the parameterizing function In some cases, the nonlinear function f can be chosen such that it can be decomposed into two multiplicative matrices of the form f(λ) = CV (λ). Matrices C and V are the coefficient and the variable matrices, respectively. The details of this decomposition will be explained in a subsequent section through mathematical expressions. There are two main advantages of selecting a function f such that the varying terms can be extracted by the above decomposition. • Varying terms can be extracted from the system and the model can simply be transformed into a Linear Fractional Transformation (LFT) form shown in Figure 2.3, where K is the controller, G is the linear varying system, Λ is a matrix representing the variations, and P is the corresponding linear invariant system. For the LFT form, there are many robust control techniques available. • The third step of the Algorithm 2.3.1 becomes a convex optimization problem min C max `=1,...,L ‖CV ( [ λu∗` λs` ] ) − θ`‖, (2.11) and can be solved by well-developed techniques, such as LMI. As an example, assume that the function f : Rnλ 7−→ Rnθ is in the form of f(λ) =  0.025 6  +   0.012 4  λ1 +   01.5 4.5  λ2 +   00 3  λ1λ2, 29Figure 2.3: A closed-loop configuration for the special case where nλ = 2 and nθ = 3. Then, the matrices C ∈ R3×6 and V ∈ R6×1 can be decomposed as C =  0.02 0.01 0 0 0 05 2 1.5 0 0 0 6 4 4.5 0 0 3   , V (λ) =   1 λ1 λ2 λ21 λ22 λ1λ2   . In general, the variable matrix V (λ) is in the form of V (λ) :=   v0(λ) v1(λ) . . . vm(λ)   , (2.12) where m is the order of the polynomial function f , and vi is the Veronese map of degree i [36]. A veronese map vi : Rnλ 7−→ RJi(nλ) is a map from [λ1, · · · , λnλ ] to all Ji(nλ) monomials of degree i in λ1, . . . , λnλ . In the above 30example, v0 = 1, J0(2) = 1, v1(λ) = [ λ1 λ2 ] , J1(2) = 2, v2(λ) =   λ 2 1 λ22 λ1λ2   , J2(2) = 3. The size of the general matrix C is nθ × p, where p := ∑m i=0 Ji(nλ). The special case, m = 0, corresponds to a single transfer function. Pa- rameter m expresses the order of the uncertainty model. One may achieve a more accurate model by increasing the polynomial order, m. Notice that the parameter m assigns the size of the V matrix and consequently the size of the Λ block in Figure 2.3. Therefore, increasing m leads to a higher computa- tional cost in controller design and possibly to an unacceptable conservative performance. We would like to obtain a coefficient matrix C such that the connected set ̂Θ satisfactorily interpolates all the elements of Θ. For this purpose, we follow the Algorithm 2.3.1 with the special convex form in the third step shown in (2.11). 2.4 Examples In this section, the developed algorithm is validated through numerical and practical examples. First, we study the results of our method applied to the motivating example explained in Section 1.2.1. Then, the uncertainty asso- ciated with a multi-input multi-output BSD is modeled using the proposed method. 312.4.1 Connected model set for a numerical example We show the performance of the proposed algorithm by applying to the example explained in Section 1.2.1. Let us assume that the variations of the FRF data (Figure 1.1) are caused by only uncertain variables. The goal is to model these uncertainties. For this purpose, we employ a polynomial function f in the formulation (2.2). First, we apply linear PCA in order to model the uncertainties. Linear PCA is one possible method in dealing with parameter variations for robust control design [56, 75, 76]. To evaluate the final resultant connected set ̂Θ, we study how precisely it approximates the original parameter vectors. To this end, we obtain the nearest point in this connected set to each of the original parameter vector. The resultant parameter vectors are ̂ΘPCA = {ˆθ`}L`=1 =      0.038.701 17.061   ,  0.0285.523 2.222   ,  0.0184.350 4.185   ,  0.0092.166 0.168   ,  0.0215.215 6.388      . By comparing this resultant parameter set with the original one shown in Section 1.2.1, Equation (1.5), it can be seen that there are some deviations which might not be acceptable for many applications. Therefore, one may want to derive a nonlinear uncertainty model rather than a linear one in order to cover the variations of the parameter set Θ more accurately. By applying the proposed method, we derive different nonlinear uncer- tainty models for different values of nλ, number of uncertain terms, and m, the order of the polynomial function. To compare the results quantitatively, the error is defined as  := L∑ `=1 ‖θ` −̂θ`‖, (2.13) where θ` is the original parameter vector from (1.5), and ˆθ` is the nearest point in the connected set ̂Θ to θ`. 32Table 2.1: Errors for different models nλ = 1 nλ = 2 m = 1 1 0.12 m = 2 0.9 9 × 10−9 Error values  for different models are shown in Table 2.1, where the error of the linear PCA is normalized to 1. From (1.6), we know that the actual values are nλ = 2 and m = 2, which concur with the results shown in the table. 2.4.2 Uncertainty modeling of a dual-input dual-output Ball Screw Drive (BSD) As a practical example, we consider to model the uncertainties associated with a machine tool with a ball screw drive. A picture of this machine is illustrated in Figure 5.1. This machine is a dual-input dual-output system. A motor controls the position of a table along a shaft. For different locations of the table, the mass distribution of the system varies, and as a result, the characteristics of the plant change during the operation. It is assumed that this dual-input dual-output plant is governed by the following structure: G(s) = [ M1(s) M2(s) M2(s) M3(s) ] . To control the position of the table for all locations along the shaft, one way is to build a continuous uncertainty model such that it covers all the variations. For this purpose, we take some sample FRF data illustrating the system characteristics for different locations of the table, {(ωw,̂Miw`) : ωw ∈ R, ̂Miw` ∈ C, i = 1, 2, 3, w = 1, . . . ,W, ` = 1, . . . , L}, where ̂Miw` represents the system response of Mi(s) at the `-th location at 33the frequency of ωw, and W is the number of frequencies. Figure 2.4 shows the dynamic responses (dotted lines) of the plant for three different table locations. For the given FRF data, no a priori knowl- edge about the system transfer functions M is assumed except being time- invariant and linear. Therefore, we estimate the general rational form of the transfer functions whose orders are assumed to be identical for all the sam- ples, {G`(s)}L`=1. This implies that the variations shown in Figure 2.4 are due to parameter variations only. The first step is to obtain the order of the numerator and denominator polynomials of the transfer functions M . The results are as follows: M1(s) has 2 zeros and 4 poles, M2(s) has 1 zeros and 5 poles, and M3(s) has 3 zeros and 5 poles, which are selected by inspecting the shapes of the FRF data. Now, the least squares optimization is used to estimate the parameter sets Θi by minimizing the error as follows for i = 1, 2, 3. Θi := {θ` = arg min θ W∑ w=1 |̂Miw` − [Mi(θ)](jωw)|2, ` = 1, . . . , L}, (2.14) where W is the number of frequencies for the FRF data. The responses of the transfer functions corresponding to these parameter sets are shown with solid lines in Figure 2.4. Then, we define the parameter set Θ := {Θ1 ⊕ Θ2 ⊕ Θ3}, and model the uncertainty of this parameter set. First, we apply linear PCA (m = 1) with one uncertain variable (nλ = nu = 1). The nearest responses in the connected set to the original ones are shown by dashed lines in Figure 2.5. As can be seen in the figure, this model captures the main features but it is not capable to capture the entire ranges of variations. This is concluded because the dashed lines lie in between the solid lines. By increasing either nλ or m, the model captures the variations more satisfactorily. The numerical comparison is made in Table 2.2. Similar to Table 2.1, to evaluate the resultant uncertainty models, we compare the 34100 200 300 −180 −135 −90 −45 0 Frequency (rad/sec) Phase (deg ) −180 −160 −140 −120 Magnitude (dB ) (a) M1 100 200 300 −400 −200 Frequency (rad/sec) Phase (deg ) −170 −150 −130 Magnitude (dB ) (b) M2 100 200 300 −300 −200 −100 0 Frequency (rad/sec) Phase (deg ) −170 −150 Magnitude (dB ) (c) M3 Figure 2.4: Frequency responses of the ball screw machine (dotted line), and the estimated transfer function (solid line). 35Table 2.2: Ball screw uncertainty model normalized errors for different values of nλ and m nλ = 1 nλ = 2 nλ = 3 m = 1 1 2.3 × 10−9 7.6 × 10−9 m = 2 0.15 2.6 × 10−11 2 × 10−15 m = 3 0.12 1 × 10−15 1.5 × 10−15 errors defined in (2.13) between the parameter sets. The error of the linear PCA is normalized to unity. Depending on the complexity and accuracy requirements in a particular application, the values of nλ and m should be selected accordingly. As can be seen in Table 2.2, if the linear PCA is applied, m = 1, we are not able to achieve the normalized error less than the order of 10−9, and by increasing the polynomial order to m = 2, a significant improvement might be accomplished. For instance, in the case of nλ = 3, the error decreases by a factor of 106. The resultant frequency responses for nλ ≥ 2 almost overlap the estimated dynamics responses. 2.5 Conclusions This chapter proposed a systematic method to derive a connected set for a given set of system response data associated with parameter variations. This method achieved Objective (O1) in Section 1.3. The algorithm illustrated in Figure 1.1 was described with mathematical expressions in this chapter. This algorithm developed a mapping from experimental system response data of a batch of sampled plants to a connected set. The resultant set was a set of real rational transfer functions, of which parameters were parameterized by a small number of uncorrelated variables that captured the differences in the dynamics among the sampled plants. The method was based on a minimax optimization problem for the multi- dimensional principal curves. We applied our algorithm to an illustrative and 36100 200 300 −180 −135 −90 −45 0 Frequency (rad/sec) Phase (deg ) −180 −160 −140 −120 Magnitude (dB ) (a) M1 100 200 300 −400 −200 Frequency (rad/sec) Phase (deg ) −170 −150 −130 Magnitude (dB ) (b) M2 100 200 300 −300 −200 −100 0 Frequency (rad/sec) Phase (deg ) −170 −150 Magnitude (dB ) (c) M3 Figure 2.5: The estimated transfer function (solid line), and the model set calculated by applying linear PCA (dashed line). The results for nλ ≥ 2 almost overlap the solid lines. 37a practical example, and obtained accurate models in both cases, and hence circumvented the possibility of conservative performance of the closed-loop control systems. 38Chapter 3 Family of Discrete Model Sets for a Set of System Response Data 1 3.1 Introduction This chapter addresses Objective (O3) in Section 1.3. The algorithm that we propose to achieve this objective is briefly illustrated in Figure 3.1, and is described in detail with mathematical expressions in this chapter. This algorithm develops a mapping from the system response data to a family of discrete model sets. It is assumed that the governing equations of the given system responses are the same but associated with parameter variations. The goal is to divide the given set of system responses into the smallest possible number of partitions, and generate a family of discrete model sets, in such a way that the given performance objective is satisfied for all partitions by designing one controller for each partition. The common transfer function structure for the system responses can be selected by systematic methods such as the Akaike’s information criterion 1This chapter is based on the following articles which is under preparation: D. Sepasi, R. Nagamune and F. Sassani, “Performance-oriented multiple model set estimation using normalized cut” 39Figure 3.1: Family of discrete model set derivation. [1], inspection of the system responses, or by trial and error. The corre- sponding parameter set can be estimated through many powerful tools, such as the System Identification [68] and Signal Processing [67] Toolboxes in the Matlab software. Therefore, these two procedures, i.e., deriving the transfer function set from the system responses and obtaining the parameter set, are not discussed in this chapter. This chapter is organized as follows. Section 3.2 formulates a problem of deriving a family of discrete model sets. In Section 3.3, a method is proposed to solve the formulated problem. Section 3.4 provides numerical examples to validate the proposed method as well as the effectiveness of combining the method developed in this chapter with that in Chapter 2. 3.2 Problem of deriving a family of discrete model sets In Section 1.2.2, we have formulated a problem of generating a connected model set in a general term; See Problem 1.2.3. Here, using mathematical notation, we will reformulate it more rigorously. 40Problem 3.2.1. Assume a set of L transfer functions and a desired closed- loop performance are given. The transfer functions have a common structure with possibly different parameter vectors, G := {[G(θ)](s), θ ∈ Θ}, (3.1) where Θ := {θ` ∈ Rnθ , ` = 1, . . . , L} (3.2) is the related parameter set. The goal is to divide the set (3.1) into the smallest possible number of partitions, Q, and derive a partition set, ̂G := {G(q)}Qq=1, G(q) := {[G(θ)](s) : θ ∈ Θ(q)}, (3.3) such that a certain closed-loop performance is satisfied by designing one con- troller K(q) for each partition G(q). Remark. There are two remarks regarding the above problem. 1. For simplicity, the problem is formulated in a continuous-time setting. However, the discrete-time case can be treated analogously. 2. The resultant pairs of {G(q), K(q)} may be used for switching systems. However, we do not investigate potential issues related to switching between the controllers. 3.3 Synthesis of a family of discrete model sets In this section, we explain the main idea of how to tackle Problem 3.2.1 and derive the optimum partition set2. The proposed algorithm is summarized as a flowchart in Figure 3.2. Each step is briefly addressed first, and then, detailed explanations are provided in the following sections. 2Although this partition set is called optimum, it may not be the global optimum. 41Figure 3.2: Main idea of the procedure to seek the optimum partition set. (1,2) In the flowchart, there are two distinct steps labeled as (1) and (2). However, the developed algorithm performs these two steps simulta- neously as follows. For the given transfer function set, we employ a procedure, which is explained in Section 3.3.2, to obtain a family of P 42“best”3 partition sets, {̂Gp}Pp=1 = {̂G1, . . . ,̂GP }. (3.4) Note that each of these sets of partitions, ̂Gp, consists of Q partitions as shown in (3.3), ̂Gp = {G(1)p , . . . ,G(Q)p }. (3) The closed-loop performances of these P partition sets are assessed as explained in Section 3.3.3. Therefore, the partition set, which provides the best performance, is obtained. (4) If the best performance in Step (3) is satisfactory, we select the corre- sponding partition set as the optimum one. Otherwise, we go back to step (1) with Q = Q+ 1. Since the dynamics of all systems are governed by a common transfer function, partitions in the model set i.e., {G(q)}Qq=1, are equivalent to the partitions in the parameter domain, i.e., {Θ(q)}Qq=1, where Θ(q) includes pa- rameter vectors of the systems in G(q). Hence, we divide the parameter set Θ into Q partitions {Θ(q)}Qq=1 such that there exists a controller K(q) ∈ K to satisfy the performance objective of the systems in partition q, where K represents a class of stabilizing feedback controllers. In this section, first, we provide some background material, which is nec- essary before explaining the main algorithm. A procedure, which provides the family of P best partition sets, is explained in Section 3.3.2. Obtaining the optimum partition set, which leads to the best closed-loop performance is described in Section 3.3.3. The entire procedure is summarized as an algorithm in Section 3.3.4. 3Our meaning of the word “best”, and how to select the value of P will be explained in detail in Section 3.3.2. 433.3.1 Background material Notations used in this section are standard. If each element of the vector d is shown by a nonnegative number di, the vector d1/2 is defined as a vector with elements of d1/2i . Definition: The nonempty sets Θ(1), . . . ,Θ(Q) form partitions of the data points Θ if Θ(i) ∩ Θ(j) = Ø and Θ(1) ∪ . . . ∪ Θ(Q) = Θ. Given a set of data points Θ, the similarity graph S models the data set in the form of a graph for which the vertices represent the data points. If the similarity sij between two data points, θi and θj, is greater than a certain threshold, the corresponding two vertices are connected in the graph by an edge with the weight of sij. For a visual explanation see Figure 3.3, where data points and edges are shown by dots and lines, respectively. For notational brevity, only weight s12 is labeled. The similarity sij represents the local neighborhood relationships between θi and θj. Generally, the higher value sij holds, the more similar the corre- sponding points are. The goal of clustering is to find partitions of the graph S such that the edges within the partitions have high weights while the edges between the partitions have small weights. Three mostly used similarity graphs are neighborhood [48], k-nearest neighbor [30], and fully connected [80] graphs. A comprehensive survey on these graphs is given in [110], based on which, the fully connected similarity graph with Gaussian similarity weights is employed in our algorithm. Definition: The affinity matrix A = (aij)i,j=1,...,L is defined as aij = sij := e − ‖θi−θj‖2 2σ2 if i 6= j, aii = 0, (3.5) with a proper norm ‖ · ‖. The constant σ governs all the similarities throughout the entire data set, and regulates the width of the neighborhoods in the Gaussian similarity 44(a) (b) Figure 3.3: Two ways of partitioning a graph. Vertices and edges are shown by dots and solid lines, respectively. The dashed lines are partition boundaries. function. The parameter σ is in the order of the average distance of a point to its r-th nearest neighbors, where r is chosen to be the nearest integer to log(L) + 1 [10]. Definition: The degree of a vertex θi is defined as di := ∑L j=1 aij , and the degree matrix D is a diagonal L× L matrix with {d1, . . . , dL} on the diagonal. 3.3.2 Parameter set partitioning and derivation of the best partition sets This section describes Steps (1) and (2) in Figure 3.2 in detail. These steps attempt to minimize a cost function by dividing the estimated parameter set Θ into Q partitions. The outputs of this step are P partition sets, which 45provide the smallest cost function values as will be explained later. An intuitive objective of data clustering is to partition the data in such a way that the sum of the pairwise similarities of the points from different partitions is minimized. However, such a naive objective may yield a trivial solution. For example, some partitions may include only one point as shown in Figure 3.3(a). The Normalized Cut (NCut) methodology normalizes the previous naive cost function, which was the sum of the pairwise similarities of the points from different partitions, with the total weighted degree of the points in each partition. Consequently, partitions such as the one illustrated in Fig- ure 3.3(b) have higher chance to be derived in comparison with the case that a naive objective was considered. Definition: The indicator matrix is defined as ˆY ∈ RL×Q with columns of yˆj := [yˆ1j, yˆ2j, . . . , yˆLj]T , where yˆij = { 1 : if i ∈ Θj 0 : if i /∈ Θj . (3.6) For example, for the partition sets shown in Figure 3.3, we have For (a): ˆY =   1 0 1 0 1 0 0 1 1 0 1 0   , For (b): ˆY =   1 0 1 0 1 0 0 1 0 1 0 1   . (3.7) The weight of the NCut can be expressed as [113] WNCut( ˆY ) := ∑ j yˆTj (D − A)yˆj yˆTj Dyˆj = Tr[( ˆY TD ˆY )−1( ˆY T (D − A) ˆY )]. (3.8) 46Consequently, the clustering problem can be defined as the following opti- mization problem min ˆY WNCut( ˆY ). (3.9) It is an optimization problem with a nonlinear objective function, shown in (3.8), and a combinatorial nature of the feasible set ˆY , which leads to an NP-hard problem. One way to relax the NCut optimization problem (3.9) is to define an indicator matrix Y ∈ RL×Q which is the “normalized” ˆY by an unknown fac- tor, such that Y = D1/2 ˆY ( ˆY TD ˆY )−1/2 [113]. If Z := Y Y T , the optimization (3.9) can be expressed in the relaxed version as (The proof is provided in Appendix A) max Z Tr(WZ) (3.10) Zd1/2 = d1/2,Tr(Z) = Q, Z ≥ 0, I  Z  0, where I is the identity matrix with the right size, and W := D−1/2AD−1/2, (3.11) where D and A are the degree and affinity matrices, respectively. The con- straint Zd1/2 = d1/2 indicates ΣUTd1/2 = UTd1/2 where Σ and U are derived from eigenvalue decomposition, Z = UΣUT . Hence, the i-th eigenvalue Σii equals to one when (UTd1/2)i 6= 0, and can be arbitrary when (UTd1/2)i = 0. Since at the optimal solution Z∗ = Y Y T , the eigenvalues of Z∗ corresponding to Y are 1, and the indicator matrix Y is the Q eigenvectors in U correspond- ing to the Q largest eigenvalues in Σ. In order to introduce a standard Semidefinite Programming (SDP), which 47is equivalent to optimization (3.10), we define ˜Z := [ Z 1 − Z ] . (3.12) As discussed above, to obtain the indicator matrix Y , the eigenvectors of Z∗ should be extracted. The eigenvalue decomposition of ˜Z∗ is ˜Z∗ = [ Z∗ 1 − Z∗ ] = [ UΣUT U(1 − Σ)UT ] = [ U U ] [ Σ (1 − Σ) ] [ UT UT ] . (3.13) Therefore, eigenvectors of Z∗, i.e., U , can be selected from eigenvectors of ˜Z∗. The optimization problem (3.10) can be written in the standard SDP form as described in Appendix B. A standard dual pair of the SDP optimization shown in Appendix B can be modeled as a Linear Matrix Inequality (LMI), min x tTx (3.14) F (x) := Fx− [ W 0 ]  0, where x ∈ RM : M = 1.5(L2 + L) + 1, t = [ 12L+1 01.5L2−0.5L ] , (3.15) 48and F : RM → S2L is a linear operator as Fx = L∑ i=1 [ ( ˜Di + ( ˜Di)T )/2 0 ] xi + [ I/Q 0 ] xL+1 (3.16) + L∑ i=1 [ Bi Bi ] xi+L+1 + L−1∑ m=1 L∑ n=m+1 [ Cmn Cmn ] x n+(m+1)L− m2+m2 +1 + L∑ i=1 L∑ j=1 [ H ij (H ij)T ] xj+iL+L2+L2 +1 , where elements of matrices ˜Di, Bi, Cmn and H ij are zero except that • the i-th row of ˜Di is d1/2/d1/2i , • Biii = 1, • Cmnmn = Cmnnm = 1, and • H ijij = 1, for i, j = 1, . . . , L and 1 ≤ m < n ≤ L. Now, the optimal solution x∗ is used to achieve the eigenvectors of ˜Z∗. In the above primal-dual pair of SDP, optimizations (B.1) and (3.14), the strong duality is satisfied in the sense that for any optimal primal solution ˜Z∗ and any optimal dual x∗ we have4 ˜Z∗F (x∗) = 0. (3.17) Therefore, the large eigenvalues of ˜Z∗ correspond to the small eigenvalues of F (x∗). Consequently, eigenvectors of Z∗ can be selected from eigenvectors of F (x∗) by extracting its L eigenvectors, whose largest entries concentrate at 4In general, for any standard primal-dual pair of SDP, such a relationship (3.17) be- tween the optimal solutions holds [86]. 49the first half (see the decomposition in (3.13)) without knowing the primal solution ˜Z∗. Since the Q smallest eigenvalues of F (x∗) correspond to the Q largest eigenvalues of Z∗, Y can be simply selected from U by conducting eigenvalue decomposition of F (x∗). Then, the un-normalized indicator matrix ˆY should be approximated. To this end, we generate the matrix D−1/2Y ∈ RL×Q, and consider rows of this matrix as L points in RQ. These points represent the original points {θ`}L`=1 in a different domain. By applying the K-means technique, these points are divided into Q partitions. The K-means method attempts to minimize the cost υ, which is the sum of the point-to-centroid distances in each partition, summed over all Q partitions. The K-means method is well developed in the Statistics Toolbox of the Matlab software [49]. One of the outputs of the K-means algorithm is the centroid for each partition. One can put an accurate interpretation on the distances of each point, rows of the matrix D−1/2Y , to every centroid in order to rank the points of each partition. Therefore, we can determine how close each point is to the other members of its partition and to those of other partitions. As a result, the partition sets, indicated by ˆY , can be assorted by the cost υ, {{ ˆY1, υ1}, { ˆY2, υ2}, . . .} : υi < υj if i < j. (3.18) We ignore the partition sets with high costs, and keep P partition sets with the least costs. One may choose one of the following ways to select a value for P . 1. If there is a jump between υi and υi+1, we select P = i. 2. The value of P should be large enough to possibly cover all the satis- factory partition sets. After selecting P , the partition sets {̂Gp}Pp=1 corresponding to the indicator matrices ˆYp are entitled as the output of step (2) in Figure 3.2. The procedure developed in this section is illustrated in Figure 3.4. 50Figure 3.4: A short-form flowchart of the procedure explained in Sec- tion 3.3.2. 3.3.3 Optimum partition set selection This section describes Step (3) in Figure 3.2 in detail. This step attempts to obtain the optimum partition set, which leads to the best closed-loop performance, among the given P partition sets from Step (2). We suggest the following two approaches. Note that one controller K(q) satisfies the performance objective of all the systems in partition G(q). 1. Synthesize a controller set for each of the P partition sets, and the system-controller pairs are derived as {{G(q)p , K(q)p }Qq=1}Pp=1. The set with the best closed-loop performance is chosen as the optimum parti- tion set. Assume that the closed-loop performances for system-controller pairs can be denoted by {{γ(q)p }Qq=1}Pp=1, and that the better perfor- mance is achieved for the smaller γ. Therefore, the best performance can be selected as γopt := min p max q γ(q)p . (3.19) The partition set corresponding to the pair with the best closed-loop performance, γopt, is chosen as the optimum partition set. 512. Synthesize a controller set for the partition set with the minimum cost, i.e, υ1, and generate {G(q)1 , K(q)1 }Qq=1 and {γ(q)1 }Qq=1. Then, we define γ∗1 := maxq γ (q) 1 . (3.20) The performance of this controller set {K(q)1 }Qq=1 is analyzed on the other P − 1 partition sets. Therefore, we obtain {γ∗p}Pp=2 for the pairs of {{G(q)p , K(q)1 }Qq=1}Pp=2. Similarly, the best performance is selected as γopt := min p=1,...,P γ∗p . (3.21) In the first approach, P controller sets are synthesized while in the second approach, one controller set is synthesized, and analyzed for the remaining P − 1 partition sets. The user chooses one approach based on the advantages of each for a specific application. In general, controller analyzing is more straightforward than synthesizing. However, the first approach may achieve a partition set with a better closed-loop performance since controllers are synthesized for each partition set separately. 3.3.4 Optimum partition set derivation algorithm To improve the performance of the algorithm, the parameter set Θ should be normalized such that all parameters are in the range of [−1, 1]. Also, we may need to scale the normalized Θ to address the controller design limit. We multiply the i-th component of θ by a weight ρi. For example, if the parameter vector includes the natural frequency and the damping ratio of a system, and the variations in the natural frequency is more problematic in controller design, the weight for the natural frequency is higher than that for the damping ratio. The weight ρ ∈ Rnθ is either given as a priori knowledge or obtained by trial and error. Here, we summarize the proposed methods in Sections 3.3.2 and 3.3.3 as 52an algorithm. Algorithm 3.3.1. Inputs: Parameter set Θ = {θ` ∈ Rnθ ; ` = 1, . . . , L}, the number of partitions Q, and the scale factor ρ ∈ Rnθ . 1. Normalize Θ such that all parameters are in the range of [−1, 1]. 2. Scale Θ by ρ to address the controller design limit. 3. Generate constant matrices for optimization problem (3.14). 4. Follow the flowchart in Figure 3.4 to obtain P best partition sets. 5. Follow the procedure in Section 3.3.3 to derive the optimum partition set. 3.4 Numerical examples In this section, the developed method is validated through numerical ex- amples. First, we discuss the results of applying Algorithm 3.3.1 to the motivating example explained in Section 1.2.2. Then, we apply the method developed in Chapter 2 to the discrete sets to derive connected sets. Five different approaches are compared to highlight the effectiveness of the devel- oped methods. 3.4.1 Illustrative example It has been shown in Section 1.2.2 that, to achieve the desired closed-loop performance of systems (1.7), the parameter set shown in Figure 3.5 needs to be divided into at least two subsets. We apply Algorithm 3.3.1 step by step as follows to derive the optimum partition set. 1. Parameters ωn and ζ are normalized to the range of [−1, 1]. 531 1.5 2 0.1 0.15 0.2 ω n ζ Figure 3.5: The parameter set of the example in Section 1.2.2. 2. The weight ρ is chosen such that the parameter ωn is multiplied by 2, because variations in ωn deteriorate the controller performance more in comparison with that in ζ in this example. 3, 4. The flowchart in Figure 3.4 is followed, and the sorted pairs of cost- partition are derived as explained in (3.18). The normalized values of the costs υ for P = 12 partition sets are shown in Figure 3.6. Figure 3.7 shows three partition sets with the least costs. 5. As it can be seen, the global optimum partition set shown in Figure1.6(c), which obtained by checking the closed-loop performances of all the pos- sible partition sets, is among three partition sets derived in the previous step. Obviously, after checking the closed-loop performance, the global optimum partition set is chosen. The proposed algorithm takes 254 s of completion time using a computer with a 2.93 GHz processor and 2 GB of RAM to extract the optimum par- tition set, while examining the closed-loop performances of the entire 2509 possible partition sets takes 10246 s. It shows the efficiency of the proposed method. In many practical applications, the parameter set is large in terms of the dimension (more complex transfer function) and the number of ele- ments (more systems to be controlled). In such cases, the “full search” might become impractical even for a small number of partitions. 541 2 3 4 5 6 7 8 9 10 11 12 0 0.5 1 Partition set Normalized cos t Figure 3.6: The normalized values of the costs of the different partition sets. 1 1.5 2 0.1 0.15 0.2 ω n ζ (a) 1 1.5 2 0.1 0.15 0.2 ω n ζ (b) 1 1.5 2 0.1 0.15 0.2 ω n ζ (c) Figure 3.7: Three cluster sets with least cost functions. 3.4.2 Closed-loop performance comparison of connected sets In this section, the effectiveness of combining the algorithm in this chapter with that presented in Chapter 2 is studied. We consider five different ap- proaches to divide the given parameter set in the previous example into two partitions, where each of which is a connected set as follows. 551. The entire parameter variation region, {ωn ∈ [1, 2], ζ ∈ [0.1, 0.2]} shown by a rectangle in Figure 3.5, is divided into two subregions in ωn direc- tion, i.e., ωn1 ∈ [1, 1.5] and ωn2 ∈ [1.5, 2], and a controller is synthesized for each subregion. These two controllers cover the entire region, which leads to a conservative solution. 2. The entire parameter variation region is divided into two subregions in ζ direction, i.e., ζ1 ∈ [0.1, 0.15] and ζ2 ∈ [0.15, 0.2], and a controller is designed for each subregion. Similar to the previous approach, this approach is conservative. 3. Algorithm 2.3.1 is applied to each of two partitions, which are generated intuitively and shown in Figure 3.8(a), to derive connected sets. The user inputs to the algorithm are ns = 0, nu = 1, N = 50, and the functions f are chosen as second order polynomials for both partitions. 4. Similar to the previous approach, the connected sets are derived for the partitions, which are generated intuitively and shown in Figure 3.8(b). 5. Algorithm 2.3.1 is applied to the optimum partition set, Figure 3.7(b). For both partitions, ns = 0, nu = 1, N = 50, while the functions f are chosen as second and third order polynomials for the partitions shown by dots and circles, respectively. Robust controllers are designed for the resultant connected sets using the Matlab Robust Control Toolbox software [5]. The resultant closed-loop per- formances γ for above approaches are shown in Table 3.1 for both partitions for each approach. Note that the controller performances deteriorate in com- parison with the performances shown in the example in Section 1.2.2. The main reason is that the controllers guarantee the performances γ for infinite number of systems (connected sets) here. On the other hand, in the previous case (discrete sets), the performances γ are guaranteed for a finite number of systems. 56According to Table 3.1, the best performance is achieved for the forth approach. However, in this approach, the partitions are derived intuitively, which might be impractical for applications with the higher dimension of the parameter domain. The second best performance is obtained by the pro- posed method, which shows the effectiveness of the developed method. For a bigger parameter domain, which can be the case in practical applications, naive approaches (similar to the first and second ones above) lead to more conservative controllers, and intuitive clustering (similar to the third and forth approaches) becomes more difficult. 1 1.5 2 0.1 0.15 0.2 ω n ζ (a) 1 1.5 2 0.1 0.15 0.2 ω n ζ (b) Figure 3.8: Two intuitive ways of clustering. 3.5 Conclusions We proposed a technique to derive a family of discrete model sets, in the form of a partition set, for a given set of system response data. It was assumed that the given response data set was governed by a common transfer func- tion with variations in parameters. A systematic algorithm was developed 57Table 3.1: The achieved closed-loop performance γ for different ap- proaches. Approach γ1 γ2 max{γ1, γ2} 1 3.3168 2.8941 3.3168 2 3.8516 3.8526 3.8526 3 3.3124 3.1582 3.3124 4 3.2249 3.2231 3.2249 5 3.2999 2.0135 3.2999 to estimate a family of discrete model sets such that a certain closed-loop performance objective is fulfilled for all the given systems by designing a cor- responding controller set. A relaxed version of Normalized Cut was applied in an algorithm to divide the system set into a few partitions. The effective- ness of the proposed method was verified through an illustrative example. Also, the effectiveness of combining this method and the one developed in Chapter 2 was shown. 58Chapter 4 Modeling and Robust Track-Following Controller Design for Hard Disk Drives 1 4.1 Introduction Hard Disk Drives (HDDs) have been used widely in many consumer elec- tronics, such as commercial computer systems, digital music players, and video-cameras for more than 50 years. They have been continuously evolv- ing to achieve higher storage capacity and miniaturized sizes. In magnetic disks, data is stored on a recording medium (in industry commonly referred to as the media), which is responsive to the presence of strong magnetic fields, but stable in their absence. Figure 4.1 shows the schematic diagram of a typical single-stage HDD. Main components are a spindle motor, one or more disks with data written on their surfaces, suspensions, heads/sliders, and a Voice-Coil Motor (VCM), 1This chapter is based on the following articles: E. Azadi Yazdi, M. Sepasi, F. Sassani and R. Nagamune, “Automated multiple robust track-following control system design in hard disk drives”, 2010 ASME Dynamic Systems and Control Conference, Boston, MA, and E. Azadi Yazdi, M. Sepasi, F. Sassani and R. Nagamune, “Automated multiple robust track-following control system design in hard disk drives”, to appear in IEEE Transactions on Control System Technology. 59Figure 4.1: A schematic diagram of an HDD. which rotates the arms around a pivot. Figure 4.2 shows a solid model, which includes the arms, suspensions, and heads. During operations, the disk may spin at speeds as high as 10,000 RPM by the spindle motor and generates high velocity airflow between the disk surface and the head. This high speed airflow has the effects of air bearing. A dynamic balance keeps the slider at a flying height of several nanometers over the disk surface. The VCM positions the head at the right data track, and thereby data can be read from or written to the disk. The main HDD characteristic, which is the focus of most literature in this area, is areal storage density. It is essential to decrease the tracking error of the read/write head in order to increase the areal density of HDDs. A practical approach to achieve small tracking error is to add a secondary MEMS actuator to the servo assembly, and manufacture dual-stage HDDs (see, e.g., [106], [102], and [12]). Moreover, to achieve a good performance, we require a precise positioning control of the read/write head, such as the designs explained in [70], [33], [15], and [45]. To be able to design a high performance controller, the dynamics of HDDs should be examined thoroughly. Their dynamics can vary due to many factors such as variations in the fabrication environment, the tempera- ture change during the operations, and mechanical deteriorations due to the elapse of time. This chapter addresses the Objective (O3) in Section 1.3. Here, we study the variations in the dynamics of HDDs, and derive a math- 60ematical model for tracking controller design based on the FRF data. This chapter is organized as follows. The experimental setup is briefly described in Section 4.2. The dynamics of this setup is explained in Sec- tion 4.3. The variations in the dynamics due to the manufacturing limits and the temperature change are discussed in detail. Modeling of the sys- tem is presented in Section 4.4. Section 4.5 explains the design of a robust controller, and demonstrates its track-following performance. Figure 4.2: The solid model, which includes the arms, suspensions, and heads of an HDD. 4.2 HDD experimental setup As a prelude to the demonstration and verification of modeling and controller design methods in the following sections, we will first describe an HDD ex- perimental setup at the University of British Columbia, shown in Figure 4.3. The equipment, listed below, is quite standard for HDD servo experiments: • laser doppler vibrometer (LDV) OFV-5000 and OFV-551 (Polytec); • anti-vibration table RS3000 (Newport); • amplifier TA105-A14 (TRUST Automation Inc.); • FFT dynamic signal analyzer 35670A (Agilent Technologies); • controller board DS1103 (dSPACE Inc.); • five hard disk drives N256 (Maxtor); • A blower heater to change the HDD temperature. 61The input to the system is the voltage to the VCM while the output is the position of the head tip, which is measured by the LDV. Figure 4.3: HDD experimental setup. 4.3 Dynamics of HDDs It has been well-recognized that not only performance of track-following con- trollers but also their robustness are of great importance in HDDs [19, 77]. We study two types of factors, which introduce variations in dynamics. One type causes time-invariant variations, such as the product variability. The transfer function between the VCM voltage and the read/write head position may have differences between products due to the limited precision in the manufacturing line. The other type causes time-varying variations, such as the ones due to the temperature change during the HDD operations. The 62temperature change occurs owing to the heat caused by the spinning spindles and the cooling air generated by the fan to avoid overheating. To design a precise robust control system the knowledge about the plant dynamics as well as its variations is essential. However, such knowledge is frequently unavailable a priori. Therefore, some experiments must be carried out with the system in order to estimate the lacking information. In this section, the variations in the system dynamics are studied. 4.3.1 Variations in HDD dynamics due to temperature It is well known that the temperature influences the dynamics of the sys- tems by affecting such factors as the geometry and material properties. In all probability, the temperature of the arm varies during the operation, and therefore, it is essential to consider the influence of the temperature on the HDD dynamics. Temperature effects are studied for an HDD micro-actuator in [46]. However, the temperature effects on HDD dynamics are not suffi- ciently addressed in the literature from a controller design point of view (for one of the few examples see [82]). We derived the FRF data for three distinct temperatures of the arm, 35, 45, and 55 degree Celsius. This range of temperature change is realistic for a standard HDD during operations [45]. In the experiments, the temperatures are measured at the position of the pivot. However, one can assume that the temperature of the arm is fairly uniform [79]. Figure 4.4 shows the FRF data for different temperatures. As it can be seen, the frequency and damping ratio of each mode change slightly over temperature. 4.3.2 Variations in HDD dynamics due to the manufacturing limits To investigate the variations between the dynamics of different HDDs from one production line, the frequency responses of five sample HDDs are studied. 63103 104 −140 −120 −100 −80 Gain (dB )   103 104 −400 −300 −200 −100 Frequency (rad/sec) Phase (deg )   T = 55 T = 45 T = 35 (a) FRF data for HDD at different temperatures. 104 −100 −95 −90 −85 Gain (dB )   104 −300 −250 −200 −150 Frequency (rad/sec) Phase (deg )   T = 55 T = 45 T = 35 (b) Zoomed FRF data for HDD at different temperatures. Figure 4.4: Variations due to the changes in the temperature. Two sets of FRF data for each of the five HDDs are derived, and thus we have ten FRF data sets in total. As can be seen in Figure 4.5, all ten FRF data sets have similar gain and phase curves, but manufacturing variations obviously exist. 4.4 Modeling of HDDs Based on the results shown in Section 4.3, there are variations in the HDD dynamics, which need to be taken into consideration. By comparing the results shown in Figures 4.4 and 4.5, it can be concluded that the variations in dynamics due to the manufacturing limits are more significant than those due to the temperature. Therefore, in the modeling and controller design we 64Figure 4.5: 10 FRF data (2 FRF are taken for each HDD). ignore the effect of the time-varying temperature on the system dynamics. Since only variations due to the time-invariant source, i.e. limited toler- ance in the production line, are considered, the HDD dynamics is modeled as an LTI system. The sources of this type of variations may neither be measured nor estimated. Therefore, we deal with an uncertain LTI plant. It is desired to derive a connected model set to cover these variations and to represent the sampled products as well as the unsampled ones with a tight model. By inspection of the FRF data, we have selected the model structure as [G(θ)](s) := a s2 4∏ i=1 s2 + 2ζniωnis+ ω2ni s2 + 2ζdiωdis+ ω2di (4.1) with θ ∈ R17 consisting of a and {ζni, ωni, ζdi, ωdi}4i=1. Now, it is of interest to generate a connected set of transfer functions G(̂Θ) := {[G(θ)](s) : θ ∈ ̂Θ} (4.2) or equivalently, a connected parameter set ̂Θ, such that we have the following: • each FRF data in Figure 4.5 is represented properly by one element in G(̂Θ); • the members of the set ̂Θ are parameterized with a small number of 65uncorrelated parameters. The method developed in Chapter 2 is employed to estimate the con- nected parameter set ̂Θ. The inputs to the Algorithm 2.3.1 are • ns = 0, • nu = 1 and N = 100. The parameterizing functions f is chosen by trial and error as2 f(λ) := f0 + f1λ+ f2λ2 (4.3) where λ ∈ Λ ⊂ R and fj ∈ R17, j = 0, 1, 2. The uncertainty set Λ is Λ := {λ ∈ R : λ ∈ [−1, 1]}. (4.4) Since the function f is chosen as a polynomial, we applied the special case of the algorithm explained in Section 2.3.2 to estimate the coefficient vectors f0, f1, and f2. The resulting coefficient vectors are provided in Ta- ble 4.1. Figure 4.6 shows the Bode plots of ten samples from the set G(̂Θ(f)). The figure illustrates that the obtained model set captures well the major characteristics of the 10 FRF data sets. 4.5 Robust controller design for HDDs To design robust track-following controllers, various methodologies have been proposed. These methodologies include, for example, H∞ control [28], adap- tive robust control [103], and sliding mode control [44]. In pursuing to im- prove robustness, it is inevitable to worsen the tracking performance, due to the trade-off relationship between robustness and performance inherent in 2In order to obtain a more accurate model set, one may increase the order of the polynomial, but the subsequent controller design would be numerically more demanding. The second order parametrization provides sufficient resolution in this study. 66Table 4.1: Numerical values of the coefficients of polynomial functions in (4.3). parameter f0 f1 f2 a 1360.4 -389.1 181.1 ζn1 0.612 0.0101 -0.6 ωn1 13345 127 -719 ζd1 0.2104 0.0544 -0.0825 ωd1 11077 -1 15 ζn2 0.1746 0.2306 0.5195 ωn2 14276 573 -144 ζd2 0.0329 0.013 0.0346 ωd2 15650 148 -604 ζn3 0.0196 -0.01 -0.007 ωn3 17746 -469 1300 ζd3 0.059 0.011 -0.044 ωd3 18432 -420 695 ζn4 0.466 0.0705 -0.1871 ωn4 31419 3366 2889 ζd4 0.055 -0.0129 -0.0251 ωd4 30740 1926 -1272 Figure 4.6: 10 random samples from the connected model set. controller design. Therefore, in HDD applications, it is essential to compro- mise this conflict between high performance requirements and robust perfor- mance limitations. To meet such requirements by overcoming the limitations, a robust controller design technique has been proposed in [117]. 67We employ the proposed method in [117] for the derived connected model set (4.2) to meet the design specifications by designing an H∞ robust con- troller. A typical performance specification for HDD track-following control is as follows. Robust sensitivity shaping: The sensitivity function for a plant G(s) and a controller K(s) is defined as S(s) := F(G,K) := 1 1 +G(s)K(s) . (4.5) For HDD servo control, the function S represents the transfer function from the track reference signal to the Position Error Signal (PES), and from the output disturbance to the read/write head tip position. By interpreting the shape of the Bode plot of S, the tracking performance can be determined. In particular, the low frequency gain and the bandwidth of S indicate degree of tracking accuracy and tracking speed, respectively. Therefore, constraints on the FRF S(jω) is one of the key specifications for track-following. Let us define a class of fixed-structure controllers robustly stabilizing G(̂Θ(f)) as K. How to shape S can be expressed as a weighted H∞ problem: design a controller K ∈ K, which satisfies the inequality constraint max G∈G( ̂Θ(f)) ‖F(G,K)‖∞ < γ. (4.6) Here, γ is a given positive scalar, the function F contains weighting functions, and ‖F‖∞ denotes the H∞ norm of F , i.e. the maximum singular value of F over all frequencies. Since the transfer function G(̂Θ(f)) is parameterized by the variable λ ∈ Λ, we can rewrite the left-hand side of (4.6) in terms of λ as max λ∈Λ ‖F(G(λ), K)‖∞ < γ. (4.7) where Λ is defined in (4.4). For the connected transfer function set G(̂Θ(f)), the controller is synthe- 68Figure 4.7: Frequency-domain response of closed-loop systems. sized to fulfill the above explained track-following specification. The values of the given pre-specified parameters, e.g., γ, and tuning parameters, e.g., weighting functions, are provided in [4]3. The resultant sensitivity function of the closed-loop system is shown in Figure 4.7. As can be seen, the exper- iments show generally • low gain peaks, which leads to large stability and less oscillatory time- domain responses, • low gain at low frequencies, which leads to disturbance attenuation in these frequencies, and • high bandwidth, which leads to high speed tracking. For one of the five HDDs, the open-loop and the closed-loop responses are shown in Figure 4.8. From this time-domain signal, it is evident that the controller not only suppresses the head vibration but also eliminates the drift of the head position caused by effects, such as friction and air-flow turbulence. Note that there is an offset of around 11 µm. 3The design of the robust controller for HDDs is not within the scope of this thesis. 69Figure 4.8: Time-domain response of open-loop and closed-loop sys- tems. 4.6 Conclusions In this Chapter, we studied the dynamics of HDD systems, especially, the variations in the dynamics due to the change in the temperatures and limited tolerances in the production line. It was shown that the influence of the temperature on the system dynamics was not significant, and hence was ignored in the modeling and controller design. A tight connected model set was derived based on a set of experimental FRF data. Then, the controller synthesis was described. The experimental results were demonstrated in the frequency and time domains. 70Chapter 5 Modeling and Robust Tracking Controller Design for Flexible Ball Screw Drives with Runout Effect and Mass Variation1 5.1 Introduction Most machine tools rely on precision Ball Screw Drives (BSD) to accurately position the workpiece relative to the tool. The quality of the machining outcome depends significantly on the tracking performance of the workpiece position over the desired trajectory. In order to minimize the tracking errors at all times during machining processes, feedback servo controllers must be designed carefully for any conceivable condition [16]. To achieve small tracking errors in various conditions, the servo controllers should compensate for the variations in the dynamics of the systems. Such variations occur because of nonlinearities and uncertainties inherent to real 1This chapter is based on the following articles: M. Sepasi, F. Sassani and R. Nagamune, “Tracking Control of Flexible Ball Screw Drives with Runout Effect Compensation”, 2010 ASME Dynamic Systems and Control Conference, Boston, MA, and M. Sepasi, F. Sassani and R. Nagamune, “Tracking Control of Flexible Ball Screw Drives with Runout Effect and Mass Variation”, to appear in IEEE Transactions on Industrial Electronics. 71plants. Two of the most common sources of nonlinearities in BSD systems are the structural flexibility [123] and runout [78]. The former is a natural characteristic of physical materials, and the latter results from the variability in ball screw manufacturing and assembly. In addition, during operations, the table mass, which refers to the combined equivalent mass of the ball screw table and the mass of the workpiece attached to it, is normally varying and not measurable, and hence, classified as uncertain. Such uncertainty in the mass leads to the system parametric uncertainty, which along with the nonlinearities in the dynamics, makes the servo controller design a challenging process. Extensive research has been conducted on servo control methods applied to precision motion mechanisms. Classical controllers are found to be widely used [17, 111] because of their high adaptability, simplicity, and ease of un- derstanding, designing and tuning. Despite the popularity of the classical controllers, their performances are limited due to the uncertainties and non- linearities in the systems. To effectively address these issues in controller design, sliding mode controllers [84, 112] and adaptive controllers [43] have been utilized. This chapter focuses on adaptive control, and uses the gain scheduling approach. Since the variations due to the nonlinearities of the BSDs depend on the table position [85], and the speed of the table is bounded in reality, the tracking performance can be improved by adjusting controller parameters using the table position information. Therefore, one may consider the table position as a scheduling variable and design a gain scheduling controller [91] accordingly. Furthermore, the intended control system must have a good tracking performance in the presence of uncertainty in the system dynamics. The table mass uncertainty influences the transfer function parameters, and thus, makes them uncertain. It is critical to detect the correlations between the parameters and the table mass in order to develop a tight uncertain set for the controller synthesis. It is a non-trivial task to detect such correlations 72through physical laws alone due to the complex coupling effects in the BSDs [83]. In such cases, experimental and numerical analysis methods can be used in a complementary manner. This chapter addresses the Objective (O4) in Section 1.3. Here, we study the variations in the dynamics of BSDs and develop a systematic method to make a mathematical model and design tracking controllers for BSDs based on their FRF data. We consider structural flexibility and runout of the shaft, as well as the table mass variation. Although one specific experimental setup is used for investigation of the plant dynamics and validating the proposed method, the discussion and the methodology are general enough to be applied to other BSDs. This chapter is organized as follows. The experimental setup is described in Section 5.2. The dynamics of the BSD systems is explained in Section 5.3. The position-dependent and mass-dependent variations in the dynamics are discussed in detail. Modeling of the system is presented in Section 5.4. Sec- tion 5.5 explains the design of a number of controllers, and demonstrates their track-following performances in the presence of flexibility, runout, and mass variations. 5.2 BSD experimental setup The investigation of the variations in the dynamics and the demonstration and verification of modeling and controller design methods are carried out on an experimental BSD system at the University of British Columbia. The components of the setup, which is depicted in Figure 5.1, are listed below. • A brushless DC motor. • A linear encoder with a resolution of 50 nm. • A ball screw with 20 mm shaft diameter, 20 mm screw pitch, and 360 mm stroke. • A 20 kg table sliding on roller bearing guideways. 73Figure 5.1: The ball screw experimental setup. • A controller board. • An amplifier. • An FFT dynamic signal analyzer. (not shown in the picture) All experiments reported in the following sections are carried out within the range of [0.12, 0.3] m along the ball screw shaft. 5.3 Variations in the dynamics of BSDs A schematic diagram of the mechanical structure of the BSD is shown in Figure 5.2. The objective is to control the position of the table ` by applying the motor torque τ , while the disturbance d is applied. The dynamics of the BSD varies by changing the position and the mass of the table due to some factors, such as different mass distribution of the plant. From the controller design viewpoint, it is essential to know the manner in which these factors influence the dynamics of the system, and consequently, the positioning of the table. To this end, we take FRF data for different table positions and masses. 745.3.1 Position-dependent variations The actuation torque τ , transmitted to the table, passes through the active length of the ball screw shaft, i.e. the part of the shaft between the motor and the table. The equivalent stiffness of the ball screw within its active length mainly depends on the bearing, the shaft itself, and the ball screw- nut interface. The position of the table, `, affects the active length of the shaft and the corresponding stiffness, and hence, the dynamics of the system. On the other hand, the runout phenomenon mainly occurs due to the tolerances in the bearings and in the shaft manufacturing. Figure 5.3 visually explains this phenomenon. Figure 5.3(a) shows the ideal position of the ball screw shaft which is aligned with the motor shaft, while Figure 5.3(b) shows one possible configuration of the ball screw shaft in reality. Because of this phenomenon, the table positioning dynamics depends on the rotational angle of the shaft. The cyclic dependency is experimentally revealed later. The position of the table is measurable with a bounded rate of change within the plant limits ˙`(t) ∈ [−0.1, 0.1] m/s. (5.1) Hence, table position ` is employed as a scheduling variable. Consequently, we design a gain scheduling controller, which is adapted online by the value of this variable. To achieve a high performance controller, we need to extract accurately how this scheduling variable influences the dynamics. Figure 5.2: A schematic diagram of a ball screw assembly. 75(a) Ideal configuration. (b) Exaggerated configura- tion of a case which leads to runout. Figure 5.3: Run out effect in the ball screw shaft. 600 1000 1400 1800 −80 −70 −60 −50 −40 Frequency (rad/s) Gain (dB ) Figure 5.4: Frequency responses for 37 different positions of the table along the shaft when no mass is added to the table. To obtain an accurate representation of the dynamics needed for designing controllers, the FRF data of the system is examined for several positions of the table along the ball screw shaft. In these experiments, the input to the open-loop system is the motor voltage while the output is the table position. The FRF data for the entire range of table motion, [0.12, 0.3] m, is derived for every 5 mm, and thus, for 37 frequency points in total. Figure 5.4 shows the resultant 37 frequency responses when there is no added mass to the table. Here, the system response varies by changing the table location mostly due to the nonlinearities. 5.3.2 Mass-dependent variations The mass of the workpiece changes significantly as material is removed during the cutting operations. The situation is further complicated by the fact that 76600 1000 1400 1800 −80 −60 −40 Frequency (rad/s) Gain (dB )    0 kg  3 kg  6 kg  9 kg Figure 5.5: Frequency responses for four different masses added to the table when the position is at 0.25 m. the rate of mass change cannot be easily predicted because it is different for each tool, workpiece, machine tool and cutting process combination. The mass of the table is considered uncertain but expected to remain within a bounded range during the machining operation. However, the rate of change is not bounded. Therefore, the tracking gain scheduling controller must be robust enough to handle this time-varying uncertainty. To achieve a high performance robust controller, we need to extract accurately how this uncertain variable affects the parameters of the model. To obtain an accurate representation of the dynamics needed for designing controllers, the FRF data of the system is examined for different masses attached to the table. The FRF data for four different masses (0, 3, 6, 9 kg) is taken. Figure 5.5 shows the FRF data when the table position is at 0.25 m and different masses are placed on the table. Likewise, the frequency responses vary according to table mass changes. 5.4 Linear parameter varying model of BSDs As explained earlier, we have samples of FRF data for 37 different table po- sitions and 4 different masses, totally D = 148 samples. In this section, first we identify an LTI system for each sampled FRF data. These systems have 77the same structures but different parameters. Then, we apply the method developed in Chapter 2 to estimate the correlations between the system pa- rameters and the variables, which are scheduling and uncertain variables. Next, an uncertain Linear Parameter Varying (LPV) model is derived to cover the entire variations over the table position and mass. Moreover, we apply the method developed in Chapter 3 to divide the mass variations into a few subsets, and consequently, a multiple model is derived. 5.4.1 Linear time-invariant system identification based on frequency response The structure of the transfer function for LTI models may be selected by inspection of the FRF data, application of systematic methods such as the Akaike information criterion [1], or by trial and error. A model, which cap- tures the first mode dynamics, has been developed in [109], where the ball screw system was represented by a uniform beam and some rigid bodies con- nected through springs. Then, the beam model was transformed into a two degree of freedom system. Here, suppose that the model structure is given and fixed as G(s) := G1(s)G2(s), (5.2) G1(s) := 1Js2 + β s, G2(s) := kω2n s2 + 2ζωn s+ ω2n , where J and β are the equivalent mass and viscous damping factor respec- tively, while k, ζ and ωn are modal parameters. The transfer function G1(s) models the low frequency characteristics while G2(s) captures the first struc- tural mode of the ball screw. All parameters are estimated using a nonlinear least squares optimization formulation [76]. For the sake of convenience in the formulation, instead of J , parameter J∗ := 1/J is used. The estimated parameters show that the viscous damping β remains almost constant at the 78value of β = 604.01 over different masses m and different table positions `, m ∈ [20, 29] kg, ` ∈ [0.12, 0.3] m. (5.3) On the other hand, other parameters, k, ζ, ωn and J∗, vary depending on ` and m. Hence, the parameter set Θ is introduced as Θ := {θd := [kd, ζd, ωnd, J∗d ]T , d = 1, . . . , D}. (5.4) Figure 5.6(a) shows the values of estimated parameters at the table posi- tions of every 20 mm when no mass is added to the table. Since the sampled table position interval is equal to the ball screw pitch, the runout effect is not observable here. Hence, variations in estimated parameters are mostly due to the structural flexibility. Further, the FRF data is studied for more table positions, i.e. for every 5 mm. In this case, quasi-sinusoidal variations due to the runout can be observed, as shown in Figure 5.6(b). A normalized uncertain term µ is defined as√ 1 m ∈ [x, y] = { x+ y 2 + x− y 2 µ } : x := √ 1 max(m) , y := √ 1 min(m) , µ ∈ [−1, 1]. (5.5) Value of µ can be explicitly derived from µ = [√ 1 m − x+ y 2 ] 2 x− y . (5.6) The variations in the parameters of the system dynamics due to µ are shown in Figure 5.6(c). By this transformation, m = 20 and m = 29 correspond to µ = −1 and µ = 1, respectively. The figure demonstrates that the effect of table mass on the dynamics is as significant as the effect of the table position 790.15 0.2 0.25 0.3 3.5 6 x 10 4 k 0.15 0.2 0.25 0.3 0.01 0.04 0.07 ζ 0.15 0.2 0.25 0.3 1250 1400 ω n 0.15 0.2 0.25 0.3 0.01 0.015 Position (m) J* (a) 10 different table positions when no mass is added. 0.15 0.2 0.25 0.3 3.5 6 x 10 4 k 0.15 0.2 0.25 0.3 0.01 0.04 0.07 ζ 0.15 0.2 0.25 0.3 1250 1400 ω n 0.15 0.2 0.25 0.3 0.01 0.015 Position (m) J* (b) 37 different table positions when no mass is added. −1 −0.5 0 0.5 1 0.01 0.04 0.07 ζ −1 −0.5 0 0.5 1 1200 1400 ω n −1 −0.5 0 0.5 1 0.01 0.015J* µ −1 −0.5 0 0.5 1 3.5 6 x 10 4 k (c) 4 different masses are added while the table position is at 0.25 m. Figure 5.6: Estimated transfer function parameters. 80`. It is worthwhile to note that, contrary to parameters J∗, ζ and ωn, the parameter k is theoretically independent of the mass. However, due to some coupling inherent to the physical systems and unmodeled dynamics [83], k is influenced by the mass of the table, and hence, some variations can be observed in Figure 5.6(c). The way how the table position and mass affect the varying parameters is expressed by a function ˜f , ̂Θ( ˜f) := { ˜f(`,m) ∈ R4, ` ∈ [0.12, 0.3] m, m ∈ [20, 29] kg}, (5.7) where ̂Θ estimates the parameter set Θ. The goal is to obtain an accurate and yet simple function ˜f . The uncertainty parameter µ is introduced such that it is proportional to √ 1/m, see (5.5). This is inspired since, based on the physical laws, the parameters ζ and ωn are proportional to √ 1/m and the parameter J∗ is pro- portional to 1/m, and hence, they are proportional to µ and µ2, respectively. Consequently, the set ̂Θ in (5.7) can be reformulated as ̂Θ(f) := {f(`, µ) ∈ R4, ` ∈ [0.12, 0.3] m, µ ∈ [−1, 1]}, (5.8) where the function f is readily parameterized as second order polynomials of µ. 5.4.2 Uncertain LPV modeling In general, parameters of an LPV model vary with respect to independent variables called scheduling variables. As discussed above, the parameters of the linear model (5.2) are varying in time due to the scheduling variable. Therefore, we construct an LPV model to represent the dynamics of the system. In order to derive the LPV model, we interpolate the local LTI models. This interpolation is performed through the function f in (5.8). Remark. It is not always possible to directly interpolate local models to 81derive the overall LPV model [107]. However, according to [71], such deriva- tion is implementable if the parameters of the identified local models show a smooth variation over the scheduling parameter (position ` here), and there is no sign change due to the non-uniqueness of the balancing transformation [62]. Here, the system is qualified for the direct interpolation (see Figure 5.6). Moreover, our LPV system is uncertain. Therefore, it is essential to model the parametric uncertainty associated with the table mass. Such a model should be as tight as possible to reduce the unnecessary conservatism inherent to robust controllers. Based on the knowledge about the BSD dynamics and the variations shown in Figure 5.6 the function f is expressed as f(`, µ) = θ0 + P`(`) + Pµ(µ) + α sin(ω`+ φ), (5.9) where P` and Pµ are polynomials. The polynomials P`(`) approximate the variations due to the table position change shown in Figure 5.6(a), while the sinusoidal term approximates the quasi-sinusoidal variations shown in Figure 5.6(b). The polynomials Pµ(µ) model the variations due to the mass change shown in Figure 5.6(c). The parameters φ and ω are assumed to be the same for k, ζ, ωn and J∗. The phase φ is manually selected as φ = 1 rad by trial and error, and ω is calculated, ω = 2pi/(screw pitch) = 100pi rad/m. The estimation of the other parameters in the function f , i.e., θ0, coefficients in polynomials P` and Pµ, and α, will be explained later. For the controller design procedure proposed is Section 5.5, we need to derive the LPV model in the state space form. Therefore, the system (5.2) with uncertain time-varying parameters k, ζ, ωn and J∗ and constant β is expressed in a quadruple of the state space data (AG(`, µ), BG(`, µ), CG, DG). 82The stae-space matrices are obtained for an observable canonical form as AG(`, µ) =   0 0 0 J∗(`, µ) 1 − βJ∗(`, µ) 0 0 0 0 0 − ωn(`, µ)2 0 0 1 − 2ζ(`, µ)ωn(`, µ)  , BG(`, µ) = [ 0 0 k(`, µ)ωn(`, µ)2 0 ]T , CG = [ 0 1 0 0 ] , DG = 0. (5.10) Single uncertain LPV model First, we derive a single uncertain LPV model, which estimates the system dynamics for the entire range of table position and mass variation. The method developed in Chapter 2 is employed to estimate parameters of the functions f in (5.9). The inputs to the Algorithm 2.3.1 are • ns = 1 and λs = [0.12 : 0.005 : 0.3], • nu = 1 and N = 20. The resultant estimated parameters are given in Table 5.1, and the corre- sponding values of k, ζ, ωn and J∗ are shown in Figure 5.7 with dash lines over the entire range of the table positions for all the values of the added mass. Some model errors inherent with the modeling can be observed. These errors diminish by utilizing a more detailed (complex) model, e.g. higher or- der polynomials P` and Pµ, at the cost of higher complexity of controller design and implementations. Multiple uncertain LPV model In order to improve the closed-loop performance, we employed the method developed in Chapter 3 to divide the mass variation region m ∈ [20, 29] kg 830.15 0.2 0.25 0.3 3.5 6 x 10 4 k 0.15 0.2 0.25 0.3 0.01 0.04 0.07 ζ 0.15 0.2 0.25 0.3 1200 1400 ω n 0.15 0.2 0.25 0.3 0.01 0.015 Position (m) J* (a) Different table positions when no mass is added, µ = −1. −1 −0.5 0 0.5 1 0.01 0.04 0.07 ζ −1 −0.5 0 0.5 1 1200 1400 ω n −1 −0.5 0 0.5 1 0.01 0.015J* µ −1 −0.5 0 0.5 1 3.5 6 x 10 4 k (b) Different masses are added while the table position is at 0.25 m. Figure 5.7: Transfer function parameters (solid lines), values in the es- timated single model (dash lines), boundary of partitions (ver- tical dotted lines). 84Table 5.1: Estimated parameters of the polynomial and sinusoidal functions. P` coefficient Pµ coefficients θ0 α ` µ µ2 k 58849.8 2367.4 −77822.8 1885.3 4815.0 ζ 0.0567 0.0084 −0.132 0.0053 0.01841 ωn 1299.5 −28.64 −256 −81.13 −2.28 J∗ 0.0104 0.0013 0.0081 −0.00116 −3.742 × 10−5 into two subregions. By choosing ρi = 1 for i = 1, ..., 4 in the Algorithm 3.3.1, the optimum partition set, {Θ(q)}2q=1, is derived, which consists of the follow- ing portions of mass variations m(1) ∈ [20, 23] kg, m(2) ∈ (23, 29] kg. (5.11) Then, normalized parameters µ(1) and µ(2) can be introduced using (5.5) for m(1) and m(2), respectively. Vertical dotted lines in Figure 5.7 show the boundary of the partitions. The parametrization function f (q) can be expressed as f (q)(`, µ) = θ0 + P`(`) + P (q)µ (µ) + α sin(ω`+ φ). (5.12) The polynomials P (q)µ (µ) model the variations due to the mass change within partition q shown in Figure 5.7b, or equivalent Figures 5.8(b) and 5.8(c). According to these figures, the order of the polynomials P`, P (1)µ and P (2)µ are selected as 1, 1 and 2, respectively. Subsequently, the Algorithm 2.3.1 is applied independently to each partition with the similar inputs as for single model derivation. The estimated parameters for function f (q) are given in Tables 5.2 and 5.3 for partitions µ(1) and µ(2), respectively. The resultant estimated transfer 85Table 5.2: Estimated parameters of the polynomial and sinusoidal functions for partition µ(1). P` coeff. Pµ coeff. θ0 α ` µ k 59613.2 2367.8 −77823.2 −1665.7 ζ 0.0577 0.0081 −0.133 −0.0073 ωn 1353 −28.1 −255.4 −25.1 J∗ 0.0113 0.0011 0.00809 −2.05 × 10−4 Table 5.3: Estimated parameters of the polynomial and sinusoidal functions for partition µ(2). P` coefficient Pµ coefficients θ0 α ` µ µ2 k 60458.1 2367.5 −77798.3 3450.2 1150.8 ζ 0.0754 0.009 −0.131 0.01213 −0.01211 ωn 1291.17 −28.4 −256.3 −44 −14.67 J∗ 0.0111 0.0014 0.00811 −0.001 −2.22 × 10−4 function parameters are presented in Figure 5.8 with dash lines. As it can be seen, the multiple model leads to less errors in the model in comparison with the single model, see Figure 5.7. 5.5 Controller design for the BSD A gain scheduling output-feedback controller proposed by Apkarian and Adams [2] is applied to the LPV models derived earlier. The goal is to enforce the stability and achieve a good tracking performance of the closed- loop system. The design of the controller is an iterative procedure. It does not give a global solution to the problem, but it has been demonstrated in practice to result in acceptable solutions [22, 39, 118]. 860.15 0.2 0.25 0.3 3.5 6 x 10 4 k 0.15 0.2 0.25 0.3 0.01 0.04 0.07 ζ 0.15 0.2 0.25 0.3 1200 1400 ω n 0.15 0.2 0.25 0.3 0.01 0.015 Position (m) J* (a) Different table positions when no mass is added, i.e., µ(1) = −1. −1 −0.5 0 0.5 1 0.01 0.04 0.07 ζ −1 −0.5 0 0.5 1 1200 1400 ω n −1 −0.5 0 0.5 1 0.01 0.015J* µ(1) −1 −0.5 0 0.5 1 3.5 6 x 10 4 k (b) Different table masses in µ(1) region. −1 −0.5 0 0.5 1 0.01 0.04 0.07 ζ −1 −0.5 0 0.5 1 1200 1400 ω n −1 −0.5 0 0.5 1 0.01 0.015J* µ(2) −1 −0.5 0 0.5 1 3.5 6 x 10 4 k (c) Different table masses in µ(2) region. Figure 5.8: Transfer function parameters (solid lines) and values in the estimated multi model (dash lines). 87To examine the effectiveness of the robustness over the mass variation, controllers with and without robustness are designed, and referred to as the robust controllers and non-robust controllers, respectively. The non-robust controllers are designed for the case where there is no added mass to the table. While, the robust controllers take into consideration the mass variations. Figure 5.9 shows the closed-loop block diagram used for both robust and non-robust controller synthesis schemes. In this configuration, a con- stant weight Wu is chosen based on the physical plant specifications to limit the control input u. The weighting function We in the state space form (Ae, Be, Ce, De) is included to shape the frequency response, and is tuned to reach the desired performance of the controller K, which attempts to track the reference signal r in the presence of the disturbance d. The aim of tuning parameters in We is to shape the closed-loop sensitivity function to meet the following criteria: 1. Minimize the gain in low frequencies to achieve a good disturbance rejection capability and tracking performance. 2. Maximize the cross-over frequency to obtain a fast response and good tracking performance. 3. Minimize the gain peak to provide a large stability margin and less oscillatory time-domain response. For each controller synthesis, a distinct second order weighting function We is tuned to obtain the best achievable performance. However, the constant Wu is selected as Wu = 2.5 (5.13) to limit the control input in all controller designs. 88Figure 5.9: The closed-loop block diagram. Figure 5.10: Synthesis closed-loop configurations for the case that the uncertainty is ignored in the plant G. 5.5.1 Non-robust gain scheduling controller design The LPV model in (5.10) is a general model in which variations in table position and mass are considered. To address a special case where the mass uncertainty is ignored, the parametrization (5.9) is written in the form f(`) = θ0 + P`(`) + α sin(ω`+ φ). (5.14) 89Accordingly, the quadruple of the state space data for system G in the block diagram in Figure 5.9 is AG(`) =   0 0 0 J∗(`) 1 − βJ∗(`) 0 0 0 0 0 − ωn(`)2 0 0 1 − 2ζ(`)ωn(`)  , BG(`) = [ 0 0 k(`)ωn(`)2 0 ]T , CG = [ 0 1 0 0 ] , DG = 0. (5.15) The dynamics of this LPV system can be written in the form x˙ = AG(`(t))x+BG(`(t))u, y = CGx+DGu. (5.16) Figure 5.10 is a Linear Fractional Transformation (LFT) form [121] of the configuration shown in Figure 5.9. The model PLPV (`) in Figure 5.10 is generated by combining the system (5.16) with the weighting functions. The mathematical expression of the configuration shown in this figure is  [ x˙ x˙e ] [ ee eu ] [ e ]   =   AL B1L B2LC1L D11L D12L C2L D21L D22L   ︸ ︷︷ ︸ PLPV (`)   [ x xe ] [ d r ] [ u ]   , (5.17) where x and xe are the state vectors of the plant (5.16) and the weighting 90function We in the state space form (Ae, Be, Ce, De), respectively, and AL = [ AG(`) 0 −BeCG Ae ] , B1L = [ 0 0 −Be Be ] , B2L = [ BG(`) −BeDG ]T , C1L = [ −DeCG Ce 0 0 ] , C2L = [ −CG 0 ] , D11L = [ −De De 0 0 ] , D12L = [ −DeDG Wu ]T , D21L = [ −1 1 ] , D22L = −DG(`). The resultant controller is in the state space form (AK(`), BK(`), CK(`), DK(`)) for which the order is six, equal to the sum of the orders of G and We. The matrices are AK =N−1[ ˆAK −X(AL −B2L ˆDKC2L)Y0 − ˆBKC2LY0 −XB2L ˆCK ], BK =N−1( ˆBK −XB2L ˆDK), CK = ˆCK − ˆDKC2LY0, DK = ˆDK , (5.18) where N = I − XY0, and matrices Y0 and X are Lyapunov variables [93, 91Figure 5.11: Synthesis closed-loop configurations for the uncertain plant G. page 902]. Matrices ˆAK , ˆBK , ˆCK , ˆDK , and X are in an affine fashion as ˆAK(`) = ˆAK0 + ` ˆAK1 + sin(ω`+ φ) ˆAK2 , ˆBK(`) = ˆBK0 + ` ˆBK1 + sin(ω`+ φ) ˆBK2 , ˆCK(`) = ˆCK0 + ` ˆCK1 + sin(ω`+ φ) ˆCK2 , ˆDK(`) = ˆDK0 + ` ˆDK1 + sin(ω`+ φ) ˆDK2 , X(`) = X0 + `X1 + sin(ω`+ φ)X2, (5.19) which can be written in a generic form as O = O0 + `O1 + sin(ω`+ φ)O2, (5.20) where O = { ˆAK(`), ˆBK(`), ˆCK(`), ˆDK(`), X(`)} An iterative procedure is explained in [2] to obtain the optimum values for Y0, O0, O1, and O2 for a non-robust gain scheduling controller. 5.5.2 Robust gain scheduling controller design Figure 5.11 is also a Linear Fractional Transformation (LFT) form [121] of the configuration Figure 5.9 when the system G is uncertain. To derive such a form, the uncertainty of the plant G is extracted into an upper LFT form2. 2Extracting the LFT form is well developed in the Matlab software [5]. 92The dynamics of the certain part, which is a function of `, can be written in the form x˙ = A(`(t))x+B1(`(t))w +B2(`(t))u, z = C1(`(t))x+D11(`(t))w +D12(`(t))u, y = C2(`(t))x+D21(`(t))w +D22(`(t))u, (5.21) and the uncertainty block is Λ(t) := µ(t)I, (5.22) where the size of I is the same as the order of the polynomial Pµ(µ). Then, the system (5.21) combines with the weighting functions, generating the model PLPV (`) in Figure 5.11. The mathematical expression of the configuration shown in Figure 5.11 is  [ x˙ x˙e ]   zee eu   [ e ]   =   AL B1L B2LC1L D11L D12L C2L D21L D22L   ︸ ︷︷ ︸ PLPV (`)   [ x xe ]   wd r   [ u ]   , (5.23) 93where x and xe are the state vectors of the plant (5.10) and the weighting function We in the state space form (Ae, Be, Ce, De), respectively, and AL = [ A(`) 0 −BeC2(`) Ae ] , B1L = [ B1(`) 0 0 −BeD21(`) −Be Be ] , B2L = [ B2(`) −BeD22(`) ]T , C1L =   C1(`) 0−DeC2(`) Ce 0 0   , C2L = [ −C2(`) 0 ] , D11L =   D11(`) 0 0−DeD21(`) −De De 0 0 0   , D12L = [ D12(`) −DeD22(`) Wu ]T , D21L = [ −D21(`) −1 1 ] , D22L = −D22(`). Similar to Section 5.5.1, the resultant controller is in the state space form (AK(`), BK(`), CK(`), DK(`)) for which the order is six. The optimum values for the matrices are obtained by following a procedure given in [2] for a robust gain scheduling controller. 5.5.3 Disturbance observer design By tuning the weighting function We during controller design, we try to achieve the best possible tracking performance for each controller. Moreover, in the implementations of all the designed controllers, a disturbance observer 94(DOB) has been employed in order to further reduce the tracking error by attenuating the low frequency disturbances [23, 59, 81]. The addition of the DOB does not change the open-loop dynamics as well as the transfer function between r (the reference signal) and ` (the output signal). Therefore, the dynamics of the DOB does not need to be considered in the controller design procedure, which is explained in Section5.5.2. Figure 5.12 shows the block diagram of the implemented motion con- trol system. The following transfer functions can be obtained for this block diagram `(s) =Q(s)d(s) +G(s)T (s)u(s), Q(s) := G0(s)(1 − F (s))G0(s) −G0(s)F (s) +G(s)F (s) , T (s) := [ 1 − F (s) + F (s)G0(s)G(s) ] −1 . (5.24) The transfer function G0 includes the low frequency dynamics of the nominal plant. Here, G0 has the structure of G1(s) in (5.2) with nominal values of J and β. The block F represents a stable low-pass filter with the following characteristics. Its bandwidth is limited by the frequency where the unmod- eled dynamics is significant, it has a DC gain of unity, and its relative degree is greater or equal to that of G0. Here, a filter with all these characteristics is selected as F (s) = a0(τs)2 + a1τs+ a0 , (5.25) where τ = 1/(20pi), a0 = 120, and a1 = 200. The structure of the filter is chosen from [100], and the parameters are tuned to fulfill above characteristic requirements. The low frequency disturbance d is attenuated since at low frequencies F (s) ≈ 1, which means that Q(s) vanishes. Moreover, the DOB has small 95Figure 5.12: Controller with disturbance observer scheme. influence on the plant dynamics, i.e. open-loop transfer function between u and `, because T (s) ≈ 1 for the entire frequency range. The proximity of T (s) to unity is due to the facts that the gain of F (s) attenuates at high frequencies, and that G(s) ≈ G0(s) at low frequencies. This phenomenon is proven by simulation in Figure 5.13. This figure shows various frequency responses of T (s) for different values of parameters k, ζ, ωn and J∗. Based on the magnitude and the phase of the transfer function T (s), the effect of T (s) on the system dynamics can be neglected over the entire frequency range. Remark. The performance of the DOB depends strongly on the accuracy of the estimation of low frequency dynamics, e.g. the nominal values of J and β. Since the table mass is uncertain, the DOB is designed for the conservative case. In other words, the cutoff frequency of the low pass filter F (s) is chosen based on the plant with the smallest crossover frequency. 5.6 Controller results In order to compare the tracking performance of the controllers, a trajectory for the machine table position is generated. Considering the test bed limits, the following values are chosen in the trajectory generation: stroke of 0.18 m, velocity of 0.09 m/s, acceleration of 0.7 m/s2, and jerk of 50 m/s3. The trajectory moves the table along the shaft from one end, farthest from the motor, to the other and then returns to the starting position. Four sets of experiments are carried out in this section. 96−2 0 2 4 Magnitude (dB ) 102 103 104 −20 0 20 Phase (deg ) Frequency  (rad/sec) Figure 5.13: Frequency responses of perturbed T (s). 1. Three non-robust controllers are designed for the BSD without mass variations. Hereby, we study the importance of including the position dependent nonlinearities in the modeling and controller design steps. 2. We perturb the system parameters, and design new controllers to study the performance sensitivity to the estimated parameters provided in Table 5.1. 3. A robust gain scheduling controller for the BSD with mass variations is synthesized to study the effectiveness of the robustness in presence of mass variations. 4. A controller set for the derived multiple model is designed to study the effect of dividing the range of the mass variation on the closed-loop performance. All the above is discussed in detail in the following sections. 975.6.1 Single controller for the BSD without mass variations Considering the parametrization (5.20), the resultant controller, denoted by Krunout, compensates for the flexibility and runout effects. To synthesize a controller, denoted by Kno−runout which ignores the runout effect, the sinu- soidal terms in the parametrization (5.20) are excluded, and only the poly- nomials are considered O = O0 + `O1. (5.26) The tracking performances of two controllers, non-robust Krunout and Kno−runout, as well as a PID controller are investigated with no added mass to the table. In the design of Krunout and Kno−runout, the state space form (Ae, Be, Ce, De) of the tuned We is(   −10 −10 8 0  ,   64 0  , [ 2.3 82.3 ], [ 0.25 ] ) . (5.27) The experimental results are shown in Figure 5.14. The control input signal of the PID controller is slightly higher than those of the other two controllers. The spikes in the errors at the beginning and in the middle, when the direction of motion reverses, occur mostly due to the momentary static friction and a backlash-like effects in the nut assembly [21]. A numerical comparison based on the Mean Absolute Error (MAE) is given in Table 5.4. Tracking performance of the PID controller is the worst due to the plant nonlinearities, specially the structural flexibility. Recall that the flexibility is taken into account in both Kno−runout and Krunout, but the runout effect is considered in Krunout only. Two cases are reported in this table. “Full run” addresses the data for the entire time span of the run, i.e. t ∈ [0, 4.4] s, while “Tracking” refers to the time span of t ∈ [0.25, 2.2] ⋃ [2.52, 4.4] s, which ignores the error spikes. The track- 98Table 5.4: Controllers tracking error results (µm). PID Kno−runout Krunout Full run 32.6 21.2 17.1 Tracking 21.0 11.6 9.3 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −200 0 200 400 Tracking Error (µ m ) −2 −1 0 1 2 Control Input(V ) Time(s)   Kpoly Ksin PID Figure 5.14: Tracking errors and control inputs. ing performance of Krunout is 55.7 % and 19.8 % better than that of the PID and Kno−runout controllers, respectively. It justifies the use of Krunout in high precision application even though the design procedure is more involved. 5.6.2 Performance sensitivity of the BSD without mass variations Now, we study the sensitivity of the closed-loop performance to the accu- racy of the estimated model. The estimated nominal values of parameters k, ζ, ωn, J∗ (in Table 5.1) are perturbed one at a time by ±5 %. Different non-robust Krunout controllers, with common weighting functions (5.13) and (5.27), are designed based on the perturbed models. Tracking errors are pro- 99Table 5.5: Controllers tracking error (µm), and the percentage increase of errors in comparison with Krunout tracking error in Table 5.4. k ζ ωn J∗ −5 % 9.9 (6.4%) 9.4 (1.1%) 11.4 (22.6%) 9.2 (-1.1%) +5 % 11.3 (21.5%) 9.5 (2.1 %) 10.7 (15.1%) 11.1 (19.3%) vided in Table 5.5, when nominal values increase and decrease by 5 %. Also, a percentage showing how much the tracking error of Krunout increases by perturbing each parameter is shown in parentheses. Small tracking error changes can be neglected since they are related to the experimental setup repeatability. By perturbing some parameters, the per- formance degrades more than 20%. Therefore, it is critical to either estimate the model accurately, or if not possible, to design robust controllers. 5.6.3 Single controller for the BSD with mass variations In this section, we focus on the role and importance of the robustness of the controllers. A new robust controller Krunout is designed, where the state space form (Ae, Be, Ce, De) of the tuned We is selected as(   −3 −1 1 0  ,   16 0  , [ 0.3 21.7 ], [ 0.2 ] ) . (5.28) Figure 5.15 shows the tracking performances of the closed-loop plant when the runout effect is taken into consideration, and different masses are added to the table. A comparison between the tracking errors is given in Table 5.6. Recall that the non-robust Krunout is synthesized for the plant without added mass while the robust Krunout is robust over the table mass variation. Be- cause of not considering robustness, the non-robust Krunout yields the best 1000 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −200 0 200 400 Time(s) Tracking Error (µ m )   0 kg 3 kg 6 kg (a) Non-robust controller 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −200 0 200 400 Time(s) Tracking Error (µ m )   0 kg 3 kg 6 kg 9 kg (b) Robust controller Figure 5.15: Tracking errors of robust and non-robust Krunout con- troller for different added masses to the table. Table 5.6: Tracking error results of the controllers in the third scenario. Added mass (kg) Robust Krunout (µm) Non-robust Krunout (µm) 0 20.2 9.3 3 19.1 29.4 6 19.0 54.3 9 20.9 Unstable performance with an MAE of 9.3 µm. However, its performance deteriorates by adding some mass to the table, and eventually, the closed-loop system becomes unstable when a mass of 9 kg is added. On the other hand, the per- formance of the robust Krunout remains uniformly acceptable over different added masses. 101Table 5.7: Tracking errors and performance improvement (%) calcu- lated by 100(MAEsingle −MAEmultiple)/MAEsingle. Added mass (kg) 0 3 6 9 Single (µm) 20.2 19.1 19.0 20.9 multiple (µm) 15.1 14.5 15.9 16.6 Improvement (%) 25.2 24.1 16.3 20.6 5.6.4 Multiple controllers for the BSD with mass variations In this section, we study the effect of dividing the range of the mass variation on the closed-loop performance. A robust controller is designed for each LPV model in the derived multiple model. In controller design for both partitions, an identical We is selected as(   −5 −2.5 4 0  ,   32 0  , [ 0.44 31.23 ], [ 0.21 ] ) . (5.29) The results for the robust multiple controller show similar patterns to those shown in Figure 5.15(b) but with less tracking error. A numerical comparison is made in Table 5.7. These results imply that the multiple controller had potential for improving the track-following performance of the BSD by more than 16%. 5.7 Conclusions The variations in the dynamics of BSD systems due to the structural flex- ibility, runout, and workpiece mass variation were studied in Chapter 5. Tracking controllers for a BSD were designed, which considered flexibility and runout, as well as mass change. These three factors were explicitly in- corporated in LPV models. To build the LPV models, it was determined 102how the system parameters were affected by two scheduling and uncertain variables, namely, the measurable table position and the uncertain mass of the table. For the LPV models, we designed controllers which were scheduled by the table position and were robust over the table mass. The performances of the designed controllers were examined on the BSD experimental setup. It was experimentally demonstrated that the tracking performance im- proved significantly by taking into account the runout effect in modeling and controller design. Also, it was shown that more than 20% performance degradation occurs by perturbing some parameters by 5%. Therefore, it is critical to either estimate the model accurately, or if not possible, to de- sign robust controllers. In addition, it was verified that the consideration of robustness against mass variation in the design stage was necessary for maintaining the stability and a uniform tracking performance. Also, multi- ple model derivation and the performance of the corresponding controller set were demonstrated and discussed thoroughly. It was shown that the multiple controller had potential for improving the track-following performance of the BSD by more than 16%. 103Chapter 6 Conclusions, Contributions and Future Research Directions 6.1 Conclusions This thesis considered variations in the dynamics of linear systems, and tack- led modeling of Linear Time-Invariant (LTI) and Linear Parameter Varying (LPV) plants. These variations were assumed to be parametric, and caused by two types of variables, uncertain and scheduling. The variations in the dy- namics make the controller design challenging, and, to successfully overcome this challenge, two methods were proposed in this thesis. The method developed in Chapter 2 generated a connected model set based on a given set of system response data. This method interpolated the given system dynamics to cover the variations associated with not only these systems but also the intermediate plant dynamics. The connected model set was constructed to become simple and tight, leading to both noncon- servatism and reduced computational complexity in subsequent controller design, and hence, to improve the performance. We applied our algorithm to an illustrative and a practical example, and obtained accurate model in both cases. 104In Chapter 3, a method was developed to derive a family of discrete model sets for a given set of system response data. The idea was to divide the given set into the smallest possible number of partitions in such a way that a desired closed-loop performance was satisfied for all partitions by designing one controller for each partition. The effectiveness of the proposed method was verified through an illustrative example. Also, the effectiveness of combining this method and the one developed in Chapter 2 was shown through an example. In Chapter 4, we studied the dynamics of Hard Disk Drive (HDD) sys- tems, especially, the variations in the dynamics due to the change in tempera- tures and limited precision in the production line. A tight uncertainty model was derived based on a set of experimental frequency response data, and an H∞ controller was synthesized. The experimental results were demonstrated and discussed in the frequency and time domains. The variations in the dynamics of Ball Screw Drive (BSD) systems due to the structural flexibility, runout, and workpiece mass variation were studied in Chapter 5. Tracking controllers for a BSD were designed, which considered flexibility and runout, as well as mass change. These three factors were explicitly incorporated in LPV models. To build the LPV models, it was determined how the system parameters were affected by two scheduling and uncertain variables, namely, the measurable table position and the uncertain mass of the table. We designed controllers which were scheduled by the table position and were robust over the table mass. The performances of the designed controllers were examined on the BSD experimental setup. It was experimentally demonstrated that the tracking performance im- proved significantly by taking into account the runout effect in modeling and controller design. Also, it was shown that more than 20% performance degradation occurs by perturbing some parameters by 5%. Therefore, it is critical to either estimate the model accurately, or if not possible, to de- sign robust controllers. In addition, it was verified that the consideration 105of robustness against mass variation in the design stage was necessary for maintaining the stability and a uniform tracking performance. Also, multi- ple model derivation and the performance of the corresponding controller set were demonstrated and discussed thoroughly. It was shown that the multiple controller had potential for improving the track-following performance of the BSD by more than 16%. 6.2 Summary of contributions The contributions of this thesis are as follows. • The connected model set derivation method – The idea of the principal curves methodology in a multi-dimensional fashion is employed to detect the nonlinear correlations between parameters of the system dynamics. Therefore, the model can be parameterized by the minimum number of independent variables. The number of independent variables can be detected readily by trial and error using the developed method. – This method does not need any information about the way that uncertain and scheduling variables affect the physical parameters of the system, such as natural frequency and damping ratio, in contrast to most of the literature in this field [74, 76]. Therefore, the method is applicable to any form of the transfer functions, e.g., the general form. – The developed method is applicable to LTI and LPV systems. Such applications are demonstrated through examples in this the- sis. • The family of discrete model sets derivation method 106– A relaxed version of the normalized cut methodology is developed and used in an algorithm to divide a given set of LTI system re- sponses into the smallest possible number of partitions in such a way that a desired performance objective is satisfied for all parti- tions by designing one controller for each partition. To the best of our knowledge, there is no literature on derivation of a family of discrete model sets based on the desired closed-loop performance as described in this thesis. • Controller design for flexible BSDs with runout effect and mass varia- tion – The dynamics of BSDs is studied, where the position and mass de- pendent variations are examined in detail from a controller design point of view. LPV models are derived to represent the dynam- ics of the BSDs, which happens to be time-varying and uncer- tain. The modeling results prove the effectiveness of the proposed method. To the best of our knowledge, there is no literature con- sidering the effects of the structural flexibility, runout, and mass variations in the BSD systems simultaneously. – Tracking controller design method is proposed for BSDs, which consists of a disturbance observer and a robust gain scheduling controller. Tracking performances of a number of different con- trollers are compared, and it is demonstrated the importance of including flexibility, runout effect, and mass variations in modeling and controller design. Also, experimental results show more than 16% improvement in the tracking performance by implementing a multiple controller. 1076.3 Future research directions This section recommends a number of potential future research directions. 6.3.1 Uncertainty modeling for stochastic robust controller The method, which is developed to derive connected model sets, is success- ful in uncertainty modeling for HDD and BSD systems. One assumption is inherently made that these systems are deterministic, and not stochastic. Re- cently, stochastic robust control has received an increasing attention. Many results about conventional robust control are extended in stochastic setting since Hinrichsen [40] proposed the stochastic H∞ control. For instance, Xu and Chen [114] proposed the sufficient condition for the solvability of robust H∞ control problem for uncertain stochastic delay systems. One topic of research which emanates from this thesis is to extend the developed method to model the uncertainty for the robust controller design, when stochastic perturbations exist. 6.3.2 Performance oriented connected model set derivation The algorithm, which is explained in Chapter 2, attempts to circumvent the possibility of conservative performance of the closed-loop control systems. We do not take into account if there is a desired closed-loop performance. Refer to the special case in Table 1.1 which is denoted by “Chapter 2”. Since the closed-loop performance depends on the characteristics of the model set, considering the performance in deriving the connected model sets can be beneficial. Therefore, One research topic is to consider the special case denoted by “Future work” in Table 1.1. By considering the controller performance, we may need to divide the connected model set, and derive a family of connected model sets. The 108algorithm, which is explained in Chapter 2, derives a connected model set with the minimum size, and does not take into consideration the possibility of dividing the set, and deriving a family of connected model sets. On possible research direction is to extend the developed method in such a way that the model set is divided into a number of subsets if necessary based on the closed-loop performance. 6.3.3 Advanced performance oriented family of discrete model sets derivation In the algorithm explained in Chapter 3, first we generate a number of par- tition sets based on a clustering cost function, and then, the best partition set is chosen based on the closed-loop performance. The performance of this method can be improved by combining these two steps, and consequently, an optimization can be formulated to provide the best partition set in one step. The new method may be faster and computationally simpler and results in a better family of discrete model sets. 6.3.4 Switching controllers for BSDs In Chapter 5, we derive a multiple closed-loop system by dividing the mass variation to improve the tracking performance. However, deriving the mul- tiple system should agree with the machining operation. In this study, we divide the range of mass variations into two partitions, and assume that the table mass stays within one of these partitions during machining. However, in the machining operations, we may need to violate this boundary. For ex- ample, the machining starts with a mass within the range of (23, 29] kg, and during cutting operations the mass reduces to the range of the other parti- tion [20, 23] kg1. In this case, we need to consider switching between two controllers. Therefore, it is essential to design switching controllers, which 1These numbers are from the results in Chapter 5. 109guarantee stability and robust tracking, when switching occurs. 6.3.5 BSD table mass estimation in real time We assume that the value of the BSD table mass is not available in real time during the operations. This lack of information causes loss of performance due to the mass variations. One extension of the method proposed in this thesis is to estimate the table mass in real time. Therefore, the tracking performance can be improved by updating the mass-dependent parameters of the controllers as a function of the table mass, similar to the position- dependent gain scheduling controllers described in Chapter 5. 110Bibliography [1] H. Akaike, A new look at the statistical model identification, Automatic Control, IEEE Transactions on 19 (2002), no. 6, 716–723. → pages 7, 23, 40, 78 [2] P. Apkarian and R.J. Adams, Advanced gain-scheduling techniques for uncertain systems, IEEE Transactions on Control Systems Technology 6 (1998), no. 1, 21–32. → pages 86, 92, 94 [3] E. Azadi Yazdi, M. Sepasi, F. Sassani, and R. Nagamune, Automated multiple robust track-following control system design in hard disk drives, To appear in IEEE Transactions on Control System Technology. → pages iv [4] , Automated multiple robust track-following control system de- sign in hard disk drives, ASME Dynamic Systems and Control Con- ference (Boston, MA), September 2010, pp. 4163(1)–4163(6). → pages iv, 69 [5] G. Balas, R. Chiang, A. Packard, and M. Safonov, Robust control tool- box, Matlab Users Guide 3 (2005). → pages 56, 92 [6] Y. Bar-Shalom, X.R. Li, T. Kirubarajan, and J. Wiley, Estimation with applications to tracking and navigation, Wiley Online Library, 2001. → pages 21 [7] M. Belkin and P. Niyogi, Laplacian eigenmaps and spectral techniques for embedding and clustering, Advances in Neural Information Process- ing System 1 (2002), 585–592. → pages 17 111[8] H.A.P. Blom and Y. Bar-Shalom, The interacting multiple model al- gorithm for systems with Markovian switching coefficients, Automatic Control, IEEE Transactions on 33 (1988), no. 8, 780–783. → pages 21 [9] M. Brand, Charting a manifold, Advances in Neural Information Pro- cessing System (2003), 985–992. → pages 17 [10] MR Brito, EL Chavez, AJ Quiroz, and JE Yukich, Connectivity of the mutual k-nearest-neighbor graph in clustering and outlier detection, Statistics & Probability Letters 35 (1997), no. 1, 33–42. → pages 45 [11] G.D. Buckner, H. Choi, and N.S. Gibson, Estimating model uncer- tainty using confidence interval networks: Applications to robust con- trol, Journal of Dynamic Systems, Measurement, and Control 128 (2006), 626–635. → pages 16 [12] B. Chen, J. Miao, and F.E.H. Tay, Fabrication and characterization of DRIE-micromachined electrostatic microactuators for hard disk drives, Microsystem Technologies 13 (2007), no. 1, 11–19. → pages 60 [13] J. Chen and G. Gu, Control-oriented System Identification: An H∞ Approach, John Wiley & Sons, 2000. → pages 2 [14] L.K. Chen and A.G. Ulsoy, Identification of a driver steering model, and model mncertainty, mrom driving simulator data, Journal of Dy- namic Systems, Measurement, and Control 123 (2001), 623–629. → pages 16 [15] T.L. Chen and R. Horowitz, Design, fabrication and dynamic analysis of a PZT-actuated siliconsuspension, Proceedings of American Control Conference, vol. 2, 2001. → pages 2, 60 [16] G. Cheng and K. Peng, Robust composite nonlinear feedback control with application to a servo positioning system, IEEE Transactions on Industrial Electronics 54 (2007), no. 2, 1132–1140. → pages 71 [17] J.U. Cho, Q.N. Le, and J.W. Jeon, An FPGA-based multiple-axis motion control chip, IEEE Transactions on Industrial Electronics 56 (2009), no. 3, 856–870. → pages 72 112[18] J. Choi, R. Nagamune, and R. Horowitz, Multiple robust controller design based on parameter dependent lyapunov functions, Proceedings of 17th International Symposium on Mathematical Theory of Networks and Systems, 2006. → pages 21 [19] R. Conway and R. Horowitz, A µ-synthesis approach to guaranteed cost control in track-following servos, The International Federation of Automatic Control, 2008, pp. 833–838. → pages 62 [20] T.F. Cox and M.A.A. Cox, Multidimensional scaling, CRC Press, 2001. → pages 18 [21] JF Cuttino, TA Dow, and BF Knight, Analytical and experimental identification of nonlinearities in a single-nut, preloaded ball screw, Journal of Mechanical Design 119 (1997), 15–19. → pages 98 [22] F.A. Cuzzola, A multivariable and multi-objective approach for the con- trol of hot-strip mills, Journal of Dynamic Systems, Measurement, and Control 128 (2006), 856–868. → pages 86 [23] M. Defoort and T. Murakami, Sliding-mode control scheme for an intel- ligent bicycle, IEEE Transactions on Industrial Electronics 56 (2009), no. 9, 3357–3368. → pages 95 [24] D. Degenring, C. Froemel, G. Dikta, and R. Takors, Sensitivity analysis for the reduction of complex metabolism models, Journal of Process Control 14 (2004), no. 7, 729–745. → pages 16 [25] D. Dong and TJ McAvoy, Nonlinear principal component analysisBased on principal curves and neural networks, Computers and Chemical En- gineering 20 (1996), no. 1, 65–78. → pages 20 [26] D.L. Donoho and C. Grimes, Hessian eigenmaps: Locally linear embed- ding techniques for high-dimensional data, Proceedings of the National Academy of Sciences 100 (2003), no. 10, 5591–5596. → pages 17 [27] AL Dontchev, RP Polis, and VM Veliov, A dual method for param- eter identification under deterministicuncertainty, IEEE Transactions on Automatic Control 45 (2000), no. 7, 1341–1346. → pages 17 113[28] C. Du, SS Ge, and FL Lewis, H∞ compensation of external vibration impact on servo performance of hard disk drives in mobile applications, International Journal of Adaptive Control and Signal Processing 22 (2008), no. 4, 374–387. → pages 66 [29] J. Du, C. Song, and P. Li, Application of gap metric to model bank de- termination in multilinear model approach, Journal of Process Control 19 (2009), no. 2, 231–240. → pages 21 [30] S.A. Dudani, The distance-weighted k-nearest-neighbor rule, IEEE Transactions on Systems, Man and Cybernetics (1976), no. 4, 325–327. → pages 44 [31] R.A. Fisher, The use of multiple measurements in taxonomic problems, Annals of Eugenics 7 (1936), no. 2, 179–188. → pages 17 [32] Y. Fu and T. Chai, Nonlinear multivariable adaptive control using mul- tiple models and neural networks, Automatica 43 (2007), no. 6, 1101– 1110. → pages 21 [33] H. Fujita, K. Suzuki, M. Ataka, and S. Nakamura, A microactuator for head positioning system of hard disk drives, IEEE Transactions on Magnetics 35 (1999), no. 2 Part 1, 1006–1010. → pages 2, 60 [34] K. Fukunaga, Introduction to statistical pattern recognition, Academic Press, 1990. → pages 16 [35] S. Guattery and G.L. Miller, On the quality of spectral separators, SIAM Journal on Matrix Analysis and Applications 19 (1998), no. 3, 701–719. → pages 22 [36] J. Harris, Algebraic Geometry: A First Course, Springer, 1992. → pages 30 [37] T. Hastie, Principal curves and surfaces, Tech. report, 1984. → pages 20 [38] T. Hastie and W. Stuetzle, Principal curves., Journal of the American Statistical Association 84 (1989), no. 406, 502–516. → pages 19, 20 114[39] R. Hibino, M. Osawa, K. Kono, and K. Yoshizawa, Robust and sim- plified design of slip control system for torque converter lock-up clutch, Journal of Dynamic Systems, Measurement, and Control 131 (2009), 011008–011017. → pages 86 [40] D. Hinrichsen and A.J. Pritchard, Stochastic H∞, SIAM Journal on Control and Optimization 36 (1998), no. 5, 1504–1538. → pages 108 [41] GE Hinton and RR Salakhutdinov, Reducing the dimensionality of data with neural networks, Science 313 (2006), no. 5786, 504–507. → pages 17 [42] H. Hotelling, Analysis of a complex of statistical variables into principal components, Journal of Educational Psychology 24 (1933), no. 6, 417– 441. → pages 17 [43] C. Hu, B. Yao, and Q. Wang, Coordinated adaptive robust contour- ing control of an industrial biaxial precision gantry with cogging force compensations, IEEE Transactions on Industrial Electronics 57 (2010), no. 5, 1746–1754. → pages 2, 72 [44] Q. Hu, C. Du, L. Xie, and Y. Wang, Discrete-time sliding mode control with time-varying surface for hard disk drives, Control Systems Tech- nology, IEEE Transactions on 17 (2009), no. 1, 175–183. → pages 66 [45] T. Imamura, M. Katayama, Y. Ikegawa, T. Ohwe, R. Koishi, and T. Koshikawa, MEMS-based integrated head/actuator/slider for hard disk drives, IEEE/ASME Transactions on Mechatronics 3 (1998), no. 3, 166–174. → pages 2, 60, 63 [46] , MEMS-based integrated head/actuator/slider for hard disk drives, Mechatronics, IEEE/ASME Transactions on 3 (2002), no. 3, 166–174. → pages 63 [47] A.K. Jain, M.N. Murty, and P.J. Flynn, Data clustering: a review, ACM computing surveys (CSUR) 31 (1999), no. 3, 264–323. → pages 21 115[48] J.W. Jaromczyk and G.T. Toussaint, Relative neighborhood graphs and their relatives, Proceedings of the IEEE 80 (2002), no. 9, 1502–1517. → pages 44 [49] B. Jones, Matlab Statistics Toolbox, The MathWorks, Inc. Natick, MA, USA (1997). → pages 50 [50] S. Kanev, C. Scherer, M. Verhaegen, and B. De Schutter, Robust output-feedback controller design via local BMI optimization, Automat- ica 40 (2004), no. 7, 1115–1127. → pages 12 [51] S.E. Karisch and F. Rendl, Semidefinite programming and graph equipartition, Topics in Semidefinite and Interior-Point Methods (1998), 77–95. → pages 22 [52] P.V. Kokotovic, L. Menini, T. Nicosia, L. Zaccarian, and C.T. Abdal- lah, Current Trends in Nonlinear Systems and Control, Springer, 2006. → pages 16 [53] RL Kosut, MK Lau, SP Boyd, I.S. Inc, and CA Santa Clara, Set- membership identification of systems with parametric and nonparamet- ric uncertainty, IEEE Transactions on Automatic Control 37 (1992), no. 7, 929–941. → pages 2 [54] M.A. Kramer, Nonlinear principal component analysis using autoasso- ciative neural networks, AIChE Journal 37 (1991), no. 2, 233–243. → pages 19 [55] U. Kruger, J. Zhang, and L. Xie, Developments and applications of nonlinear principal component analysis-a review, Lecture Notes in Con- putational Science and Engineering 58 (2007), 1. → pages 18 [56] A. Kwiatkowski and H. Werner, Parameter reduction for LPV sys- tems via principle components analysis, Proceedings of the Interna- tional Federation of Automatic Control, 2005. → pages 32 [57] J.T.Y. Kwok and I.W.H. Tsang, The pre-image problem in kernel meth- ods, IEEE Transactions on Neural Networks 15 (2004), no. 6, 1517– 1525. → pages 20 116[58] S. Lafon and A.B. Lee, Diffusion maps and coarse-graining: A uni- fied framework for dimensionality reduction, graph partitioning, and data set parameterization, IEEE Transactions on Patern Analysis and Machine Intelligence (2006), 1393–1403. → pages 18 [59] C.Y. Lai, F.L. Lewis, V. Venkataramanan, X. Ren, S.S. Ge, and T. Liew, Disturbance and friction compensations in hard disk drives using neural networks, IEEE Transactions on Industrial Electronics 57 (2010), no. 2, 784–792. → pages 95 [60] D. Lainiotis, Optimal adaptive estimation: Structure and parameter adaption, Automatic Control, IEEE Transactions on 16 (1971), no. 2, 160–170. → pages 21 [61] G.R.G. Lanckriet, N. Cristianini, P. Bartlett, L. El Ghaoui, and M.I. Jordan, Learning the kernel matrix with semidefinite programming, The Journal of Machine Learning Research 5 (2004), 27–72. → pages 17 [62] A. Laub, M. Heath, C. Paige, and R. Ward, Computation of system balancing transformations and other applications of simultaneous diag- onalization algorithms, IEEE Transactions on Automatic Control 32 (2002), no. 2, 115–122. → pages 82 [63] R. Li et al., Hybrid estimation techniques, Control and Dynamic Sys- tems 76 (1996), 213–287. → pages 21 [64] R. Li, MA Henson, and MJ Kurtz, Selection of model parameters for off-line parameter estimation, IEEE Transactions on Control Systems Technology 12 (2004), no. 3, 402–412. → pages 16 [65] R. Li and V.P. Jilkov, Survey of maneuvering target tracking. Part V. Multiple-model methods, Aerospace and Electronic Systems, IEEE Transactions on 41 (2005), no. 4, 1255–1321. → pages 21 [66] X.R. Li, Engineer’s guide to variable-structure multiple-model estima- tion for tracking, Multitarget-multisensor tracking: Applications and advances. 3 (2000), 499–567. → pages 21 [67] J.B. Little, Mass.) MathWorks (Natick, and L. Shure, Signal processing toolbox for use with MATLAB: user’s guide, MathWorks, 1988. → pages 24, 40 117[68] L. Ljung, System identification toolbox, The Matlab users guide. → pages 24, 40 [69] , System Identification: Theory for the User., Prentice-Hall, 1999. → pages 2 [70] Y. Lou, P. Gao, B. Qin, G. Guo, E.H. Ong, A. Takada, and K. Okada, Dual-Stage servo with on-slider PZT microactuator for hard disk drives, IEEE Transactions on Magnetics 38 (2002), no. 5, 2183. → pages 2, 60 [71] M. Lovera and G. Mercere, Identification for gain-scheduling: a bal- anced subspace approach, American Control Conference, 2007, pp. 858– 863. → pages 82 [72] B. Lu and F. Wu, Switching LPV control designs using multi- ple parameter-dependent Lyapunov functions, Automatica 40 (2004), no. 11, 1973–1980. → pages 21 [73] D. Magill, Optimal adaptive estimation of sampled stochastic processes, Automatic Control, IEEE Transactions on 10 (1965), no. 4, 434–439. → pages 21 [74] R. Nagamune and J. Choi, Parameter reduction of nonlinear least- squares estimates via nonconvex optimization, 2008 American Control Conference, Seattle, Washington, USA, 1298–1303. → pages 16, 106 [75] , Parameter reduction of nonlinear least-squares estimates via the singular value decomposition, Proceedings of the International Fed- eration of Automatic Control, 12383–12388. → pages 32 [76] , Parameter Reduction in Estimated Model Sets for Robust Con- trol, Journal of Dynamic Systems, Measurement, and Control 132 (2010), 021002.1–021002.10. → pages 16, 32, 78, 106 [77] R. Nagamune, X. Huang, and R. Horowitz, Robust control synthesis techniques for multirate and multi-sensing track-following servo sys- tems in HDDs, ASME Journal of Dynamic Systems, Measurement and Control. → pages 62 118[78] K. Nagaoka, A. Matsubara, T. Fujita, and T. Sato, Analysis method of motion accuracy using nc system with synchronized measurement of tool-tip position, International Journal of Automation Technology 3 (2009), no. 4, 394–400. → pages 72 [79] P. Naphon and S. Maharchon, Temperature distribution of read/write head soldering with ribbon cable of HDD, International Communica- tions in Heat and Mass Transfer 37 (2010), no. 4, 379–384. → pages 63 [80] A. Ng, M. Jordan, and Y. Weiss, On spectral clustering: Analysis and an algorithm, Advances in Neural Information Processing Systems, 2001, pp. 849–856. → pages 44 [81] K. Ohnishi, N. Matsui, and Y. Hori, Estimation, identification, and sensorless control in motion control system, Proceedings of the IEEE 82 (1994), no. 8, 1253–1265. → pages 95 [82] K. Ohno, Y. Abe, and T. Maruyama, Robust following control design for hard disk drives, Control Applications, 2001.(CCA’01). Proceedings of the 2001 IEEE International Conference on, IEEE, 2002, pp. 930– 935. → pages 63 [83] C. Okwudire, Modeling and control of high speed machine tool feed drives, Ph.D. Thesis, University of British Columbia (2009). → pages 73, 81 [84] C. Okwudire and Y. Altintas, Minimum tracking error control of flex- ible ball screw drives using a discrete-time sliding mode controller, Journal of Dynamic Systems, Measurement, and Control 131 (2009), 051006–1–051006–12. → pages 72 [85] C.E. Okwudire and Y. Altintas, Hybrid modeling of ball screw drives with coupled axial, torsional, and lateral dynamics, Journal of Mechan- ical Design 131 (2009), 071002. → pages 72 [86] M.V. Ramana, L. Tunc¸el, and H. Wolkowicz, Strong duality for semidefinite programming, SIAM Journal on Optimization 7 (1997), no. 3, 641–662. → pages 49 119[87] Li Rong et al., Multiple-model estimation with variable structure-Part VI: Expected-mode augmentation, IEEE transactions on aerospace and electronic systems 41 (2005), no. 3, 853–867. → pages 21 [88] P. Rosa, C. Silvestre, J.S. Shamma, and M. Athans, Multiple-model adaptive control with set-valued observers, Proceedings of the 48th IEEE Conference on Decision and Control., IEEE, 2010, pp. 2441– 2447. → pages 21 [89] S. Roweis, L.K. Saul, and G.E. Hinton, Global Coordination of Local Linear Models, Proceedings of Advances in Neural Information Pro- cessing Systems, MIT Press, 2002. → pages 17 [90] S.T. Roweis and L.K. Saul, Nonlinear dimensionality reduction by lo- cally linear embedding, Science 290 (2000), no. 5500, 2323–2326. → pages 17 [91] WJ Rugh and J.S. Shamma, A survey of research on gain-scheduling, Automatica 36 (2000), no. 10, 1401–1425. → pages 2, 72 [92] J.W. Sammon, A nonlinear mapping for data structure analysis, IEEE Transactions on Computers 18 (1969), no. 5, 401–409. → pages 18 [93] C. Scherer, P. Gahinet, and M. Chilali, Multiobjective output-feedback control via LMI optimization, IEEE Transactions on Automatic Con- trol 42 (1997), no. 7, 896–911. → pages 91 [94] B. Scholkopf, A. Smola, and K.R. Muller, Nonlinear component anal- ysis as a kernel eigenvalue problem, Neural Computation 10 (1998), no. 5, 1299–1319. → pages 19 [95] B. Scholkopf and A.J. Smola, Learning with kernels: support vector machines, regularization, optimization, and beyond, MIT Press, 2002. → pages 19 [96] D. Sepasi, F. Sassani, and R. Nagamune, Tracking control of flexible ball screw drives with runout effect and mass variation, To appear in IEEE Transactions on Industrial Electronics. → pages iv [97] M. Sepasi, F. Sassani, and R. Nagamune, Parameter uncertainty mod- eling using the multi-dimensional principal curves, Journal of Dynamic 120Systems, Measurement and Control 132 (2010), 054501–054507. → pages 4 [98] , Tracking control of flexible ball screw drives with runout ef- fect compensation, ASME Dynamic Systems and Control Conference (Boston, MA), September 2010, pp. 4039(1)–4039(6). → pages iv [99] J. Shi and J. Malik, Normalized cuts and image segmentation, IEEE Transactions on Pattern Analysis and Machine Intelligence 22 (2002), no. 8, 888–905. → pages 22 [100] H. Shim and Y.J. Joo, State space analysis of disturbance observer and a robust stability condition, Decision and Control, 2007 46th IEEE Conference on, IEEE, 2007, pp. 2193–2198. → pages 95 [101] C. Sun and J. Hahn, Parameter reduction for stable dynamical systems based on Hankel singular values and sensitivity analysis, Chemical En- gineering Science 61 (2006), no. 16, 5393–5403. → pages 16 [102] N. Tagawa, K.I. Kitamura, and A. Mori, Design and fabrication of MEMS-based active slider using double-layered composite PZT thin film in hard disk drives, IEEE Transactions on Magnetics 39 (2003), no. 2, 926–931. → pages 60 [103] H.D. Taghirad and E. Jamei, Robust performance verification of adap- tive robust controller for hard disk drives, Industrial Electronics, IEEE Transactions on 55 (2008), no. 1, 448–456. → pages 66 [104] J.B. Tenenbaum, Mapping a manifold of perceptual observations, Ad- vances in Neural Information Processing System (1998), 682–688. → pages 18 [105] J.B. Tenenbaum, V. Silva, and J.C. Langford, A global geometric framework for nonlinear dimensionality reduction, 2000, pp. 2319– 2323. → pages 17 [106] H. Toshiyoshi, M. Mita, and H. Fujita, A MEMS piggyback actuator for hard-disk drives, Journal of Microelectromechanical Systems 11 (2002), no. 6, 648–654. → pages 60 121[107] R. To´th, F. Felici, P.S.C. Heuberger, and P.M.J. Van den Hof, Discrete time lpv i/o and state space representations, differences of behavior and pitfalls of interpolation, Proceedings of the European Control Confer- ence, 2007, pp. 5418–5425. → pages 82 [108] C.C. Tsai, H.C. Huang, and S.C. Lin, Adaptive neural network control of a self-balancing two-wheeled scooter, IEEE Transactions on Indus- trial Electronics 57 (2010), no. 4, 1420–1428. → pages 2 [109] K.K. Varanasi and S.A. Nayfeh, The dynamics of lead-screw drives: Low-order modeling and experiments, Journal of Dynamic Systems, Measurement, and Control 126 (2004), 388–396. → pages 78 [110] U. Von Luxburg, A tutorial on spectral clustering, Statistics and Com- puting 17 (2007), no. 4, 395–416. → pages 44 [111] R. Wai, J. Lee, and K. Chuang, Real-Time PID control strategy for Maglev transportation system via particle swarm optimization, IEEE Transactions on Industrial Electronics (2010), no. 99, 1. → pages 72 [112] Y. Xia, M. Fu, H. Yang, and G.P. Liu, Robust sliding-mode control for uncertain time-delay systems based on delta operator, IEEE Transac- tions on Industrial Electronics 56 (2009), no. 9, 3646–3655. → pages 72 [113] E. Xing, E.P. Xing, M. Jordan, and M.I. Jordan, On semidefinite re- laxations for normalized k-cut and connections to spectral clustering, Tech. report, 2003. → pages 22, 46, 47, 125 [114] S. Xu and T. Chen, Robust H∞ control for uncertain stochastic systems with state delay, IEEE Transactions on Automatic Control 47 (2002), no. 12, 2089–2094. → pages 108 [115] P. Yan and H. Ozbay, On switching controllers for a class of linear parameter varying systems, Systems & Control Letters 56 (2007), no. 7- 8, 504–511. → pages 21 [116] K.Z. Yao, B.M. Shaw, B. Kou, K.B. McAuley, and DW Bacon, Model- ing ethylene/butene copolymerization with multi-site catalysts: parame- ter estimability and experimental design, Polymer Reaction Engineering 11 (2003), no. 3, 563–588. → pages 16 122[117] E.A. Yazdi and R. Nagamune, Multiple robust H∞ controller design using the nonsmooth optimization method, International Journal of Ro- bust and Nonlinear Control 20 (2010), no. 11, 1197–1212. → pages 21, 67, 68 [118] F. Zhang, K.M. Grigoriadis, M.A. Franchek, and I.H. Makki, Linear parameter-varying lean burn air-fuel ratio control for a spark ignition engine, Journal of Dynamic Systems, Measurement, and Control 129 (2007), 404–414. → pages 86 [119] Z. Zhang and H. Zha, Principal manifolds and nonlinear dimensionality reduction via tangent space alignment, Journal of Shanghai University (English Edition) 8 (2004), no. 4, 406–424. → pages 17 [120] K. Zhou and J.C. Doyle, Essentials of robust control, vol. 104, Prentice Hall New Jersey, 1998. → pages 2, 11, 20 [121] K. Zhou, J.C. Doyle, and K. Glover, Robust and optimal control, Pren- tice Hall Englewood Cliffs, NJ, 1996. → pages 90, 92 [122] T. Zhou and D. Tao, Fast gradient clustering, NIPS 2009 Workshop on Discrete Optimization in Machine Learning: Submodularity, Sparsity & Polyhedra (NIPS: DISCML), 2009, pp. 1–6. → pages 22 [123] O. Zirn, C. Jaeger, and T. Scholler, Model based control of machine tool manipulators, IEEE International Symposium on Industrial Elec- tronics, 2008, pp. 1267–1274. → pages 72 123Appendix A: Relaxed Form of the Optimization (3.9) It is desired to show that the combinatorial optimization problem min ˆY Tr[( ˆY TD ˆY )−1( ˆY T (D − A) ˆY )] ˆY T ˆY = IQ, ˆYij ∈ {1, 0} (A.1) can be relaxed to a non-combinatorial one as max Z Tr(WZ) Zd1/2 = d1/2,Tr(Z) = Q, Z ≥ 0, I  Z  0, (A.2) where W := D−1/2AD−1/2, A and D are L × L diagonal matrices, and d is a vector whose entries are diagonal elements of D. The cost function of the optimization (A.1) can be expanded and rewritten as Tr[( ˆY TD ˆY )−1( ˆY TD ˆY ) − ( ˆY TD ˆY )−1( ˆY TA ˆY )]. (A.3) 124Therefore, the optimization (A.1) is equivalent to max ˆY Tr[( ˆY TD ˆY )−1( ˆY TA ˆY )], ˆY ∈ RL×Q, ˆYij ∈ {1, 0}. (A.4) The trace operator is invariant under cyclic permutations, i.e., Tr(KPMN) = Tr(MNKP ), for matrices K, P , M , and N with proper sizes. Hence, the known and unknown parts of the cost function in (A.4) can be separated as Tr(A · ˆY ( ˆY TD ˆY )−1 ˆY T ). Then, by normalizing the unknown part, the cost can be written as Tr(D−1/2AD−1/2 ·D1/2 ˆY ( ˆY TD ˆY )−1 ˆY TD1/2). Using the definition of W , the above can be written as Tr(WZ), where Z := D1/2 ˆY ( ˆY TD ˆY )−1 ˆY TD1/2. The matrix Z as defined above has the following properties [113]: Zd1/2 = d1/2,Tr(Z) = Q, Z ≥ 0, I  Z  0. (A.5) By deriving the optimizer for the relaxed optimization problem (A.2), we approximate the optimizer of the original optimization problem (A.1). 125Appendix B: Standard SDP Form of the Optimization (3.10) The optimization problem (3.10) can be written in the standard SDP form as max ˜Z Tr ([ W 0 ] ˜Z ) , (B.1) E ˜Z = t, ˜Z  0, where W is defined in (3.11), E : S2L → RM is a linear operator, where M = 1.5(L2 + L) + 1, and t = [ 12L+1 01.5L2−0.5L ] , (B.2) The set S2L represents a set of symmetric matrices with the order of 2L. The detailed version of the first constraint in the optimization (B.1), E ˜Z = t, is 126Tr ([ ( ˜Di + ( ˜Di)T )/2 0 ] ˜Z ) = 1, i = 1, . . . , L, (B.3) Tr ([ I/Q 0 ] ˜Z ) = 1, (B.4) Tr ([ Bi Bi ] ˜Z ) = 1, i = 1, . . . , L, (B.5) Tr ([ Cmn Cmn ] ˜Z ) = 0, 1 ≤ m < n ≤ L, (B.6) Tr [ H ij (H ij)T ]T ˜Z  = 0, i, j = 1, . . . , L, (B.7) where the matrices ˜Di, Bi, Cmn and H ij are similar to the ones defined for the optimization in (3.14). In the constrains, (B.3) and (B.4) represent Zd1/2 = d1/2 and Tr(Z) = Q in (3.10), respectively. The constraintes (B.5) and (B.6) guarantee that the diagonal blocks of ˜Z are Z and (1 −Z), and (B.7) keeps the other entries in ˜Z zeros. 127

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Country Views Downloads
China 23 11
United States 11 0
Canada 7 0
Iran 6 0
Tunisia 4 0
Russia 2 0
France 2 0
Bangladesh 1 0
United Kingdom 1 0
Indonesia 1 0
Taiwan 1 1
Malaysia 1 0
City Views Downloads
Unknown 14 12
Shenzhen 12 7
Vancouver 7 0
Ashburn 6 0
Beijing 4 0
Saint Petersburg 2 0
Nanjing 2 4
Xi'an 2 0
Lanzhou 2 0
Mountain View 2 0
Isfahan 2 0
Elk Grove Village 1 0
London 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.24.1-0080734/manifest

Comment

Related Items