Brogan, W.L., Lee, G.K.F., Sage, A.P., Kuo, B.C., Phillips, C.L., Harbor, R.D., Jacquot, R.G., McInroy, J.E., Atherton, D.P., Bay, J.S., Baumann, W.T., Chow, M-Y. “Control Systems” The Electrical Engineering Handbook Ed. Richard C. Dorf Boca Raton: CRC Press LLC, 2000 100 Control Systems 100.1 Models Classes of Systems to Be Modeled ? Two Major Approaches to Modeling ? Forms of the Model ? Nonuniqueness ? Approximation of Continuous Systems by Discrete Models 100.2 Dynamic Response Computing the Dynamic System Response ? Measures of the Dynamic System Response 100.3 Frequency Response Methods: Bode Diagram Approach Frequency Response Analysis Using the Bode Diagram ? Bode Diagram Design-Series Equalizers ? Composite Equalizers ? Minor-Loop Design 100.4 Root Locus Root Locus Properties ? Root Loci of Digital Control Systems ? Design with Root Locus 100.5 Compensation Control System Specifications ? Design ? Modern Control Design ? Other Modern Design Procedures 100.6 Digital Control Systems A Simple Example ? Single-Loop Linear Control Laws ? Proportional Control ? PID Control Algorithm ? The Closed- Loop System ? A Linear Control Example 100.7 Nonlinear Control Systems The Describing Function Method ? The Sinusoidal Describing Function ? Evaluation of the Describing Function ? Limit Cycles and Stability ? Stability and Accuracy ? Compensator Design ? Closed-Loop Frequency Response ? The Phase Plane Method ? Piecewise Linear Characteristics ? Discussion 100.8 Optimal Control and Estimation Linear Quadratic Regulators ? Optimal Estimation: The Kalman Filter ? Linear-Quadratic-Gaussian (LQG) Control ? H ∞ Control ? Example ? Other Approaches 100.9 Neural Control Brief Introduction to Artificial Neural Networks ? Neural Observer ? Neural Control ? HVAC Illustration ? Conclusion 100.1 Models William L. Brogan A naive trial-and-error approach to the design of a control system might consist of constructing a controller, installing it into the system to be controlled, performing tests, and then modifying the controller until satisfactory performance is achieved. This approach could be dangerous and uneconomical, if not impossible. A more rational approach to control system design uses mathematical models. A model is a mathematical description of system behavior, as influenced by input variables or initial conditions. The model is a stand-in for the actual system during the control system design stage. It is used to predict performance; to carry out stability, sensitivity, and trade-off William L. Brogan University of Nevada, Las Vegas Gordon K. F. Lee North Carolina State University Andrew P. Sage George Mason University Benjamin C. Kuo University of Illinois (Urbana- Champaign) Charles L. Phillips Auburn University Royce D. Harbor University of West Florida Raymond G. Jacquot University of Wyoming John E. McInroy University of Wyoming Derek P. Atherton University of Sussex John S. Bay Virginia Polytechnic Institute and State University William T. Baumann Virginia Polytechnic Institute and State University Mo-Yuen Chow North Carolina State University ? 2000 by CRC Press LLC ? 2000 by CRC Press LLC C ONTROL M ECHANISM FOR R OCKET A PPARATUS Robert H. Goddard Patented April 2, 1946 #2,397,657 A n excerpt from Robert Goddard’s patent application: This invention relates to rockets and rocket craft which are propelled by combustion apparatus using liquid fuel and a liquid to support combustion, such as liquid oxygen. Such combustion apparatus is disclosed in my prior application Serial No. 327,257 filed April 1, 1940. It is the general object of my present invention to provide control mechanism by which the necessary operative steps and adjustments for such mechanism will be affected automatically and in predetermined and orderly sequence. To the attainment of this object, I provide control mechanism which will automatically discontinue flight in a safe and orderly manner. Dr. Goddard was instrumental in developing rocket propulsion in this country, both solid-fuel rocket engines and later liquid-fuel rocket motors used in missile and spaceflight applications. Goddard died in 1945, before this pivotal patent (filed June 23, 1941) on automatic control of liquid-fuel rockets was granted. He assigned half the rights to the Guggenheim Foundation in New York. (Copyright ? 1995, Dewray Products, Inc. Used with permission.) studies; and answer various “what-if” questions in a safe and efficient manner. Of course, the validation of the model, and all conclusions derived from it, must ultimately be based upon test results with the physical hardware. The final form of the mathematical model depends upon the type of physical system, the method used to develop the model, and mathematical manipulations applied to it. These issues are discussed next. Classes of Systems to Be Modeled Most control problems are multidisciplinary. The system may consist of electrical, mechanical, thermal, optical, fluidic, or other physical components, as well as economic, biological, or ecological systems. Analogies exist between these various disciplines, based upon the similarity of the equations that describe the phenomena. The discussion of models in this section will be given in mathematical terms and therefore will apply to several disciplines. Figure 100.1 [Brogan, 1991] shows the classes of systems that might be encountered in control systems modeling. Several branches of this tree diagram are terminated with a dashed line indicating that additional branches have been omitted, similar to those at the same level on other paths. Distributed parameter systems have variables that are functions of both space and time (such as the voltage along a transmission line or the deflection of a point on an elastic structure). They are described by partial differential equations. These are often approximately modeled as a set of lumped parameter systems (described by ordinary differential or difference equations) by using modal expansions, finite element methods, or other approx- imations [Brogan, 1968]. The lumped parameter continuous-time and discrete-time families are stressed here. Two Major Approaches to Modeling In principle, models of a given physical system can be developed by two distinct approaches. Figure 100.2 shows the steps involved in analytical modeling. The real-world system is represented by an interconnection of idealized elements. Table 100.1 [Dorf, 1989] shows model elements from several disciplines and their elemental equations. An electrical circuit diagram is a typical result of this physical modeling step (box 3 of Fig. 100.2). Application of the appropriate physical laws (Kirchhoff, Newton, etc.) to the idealized physical model (consisting of point masses, ideal springs, lumped resistors, etc.) leads to a set of mathematical equations. For a circuit these will FIGURE 100.1 Major classes of system equations. (Source: W.L. Brogan, Modern Control Theory, 3rd ed., Englewood Cliffs, N.J.: Prentice-Hall, 1991, p. 13. With permission.) ? 2000 by CRC Press LLC be mesh or node equations in terms of elemental currents and voltages. Box 6 of Fig. 100.2 suggests a generalization to other disciplines, in terms of continuity and compatibility laws, using through variables (generalization of current that flows through an element) and across variables (generalization of voltage, which has a differential value across an element) [Shearer et al., 1967; Dorf, 1989]. Experimental or empirical modeling typically assumes an a priori form for the model equations and then uses available measurements to estimate the coefficient values that cause the assumed form to best fit the data. The assumed form could be based upon physical knowledge or it could be just a credible assumption. Time- series models include autoregressive (AR) models, moving average (MA) models, and the combination, called ARMA models. All are difference equations relating the input variables to the output variables at the discrete measurement times, of the form y(k + 1) = a 0 y(k) + a 1 y(k – 1) + a 2 y(k – 2) + … + a n y(k – n) + b 0 u(k + 1) + b 1 u(k) + … + b p u(k + 1 – p) + v(k) (100.1) where v(k) is a random noise term. The z-transform transfer function relating u to y is (100.2) FIGURE 100.2 Modeling considerations. (Source: W.L. Brogan, Modern Control Theory, 3rd ed., Englewood Cliffs, N.J.: Prentice-Hall, 1991, p. 5. With permission.) yz uz bbz bz az a z Hz p p n n () () () () –– – – = +++ ?++ = ? 01 1 0 1 1 1 L L ? 2000 by CRC Press LLC In the MA model all a i = 0. This is alternatively called an all-zero model or a finite impulse response (FIR) model. In the AR model all b j terms are zero except b 0 . This is called an all-pole model or an infinite impulse response (IIR) model. The ARMA model has both poles and zeros and also is an IIR model [Makhoul, 1975]. Adaptive and learning control systems have an experimental modeling aspect. The data fitting is carried out on-line, in real time, as part of the system operation. The modeling described above is normally done off-line [Astrom and Wittenmark, 1989]. Forms of the Model Regardless of whether a model is developed from knowledge of the physics of the process or from empirical data fitting, it can be further manipulated into several different but equivalent forms. This manipulation is box 7 in Fig. 100.2. The class that is most widely used in control studies is the deterministic lumped-parameter continuous-time constant-coefficient system. A simple example has one input u and one output y. This might be a circuit composed of one ideal source and an interconnection of ideal resistors, capacitors, and inductors. TABLE 100.1 Summary of Describing Differential Equations for Ideal Elements Type of Physical Describing Energy E or Element Element Equation Power P Symbol Electrical inductance Inductive Translational storage spring Rotational spring Fluid inertia Electrical capacitance Translational mass Capacitive Rotational storage mass Fluid capacitance Thermal capacitance Electrical resistance Translational damper Energy Rotational dissipators damper Fluid resistance Thermal resistance v 21 v 2 v 1 L i L= di dt v 21 = dF dt 1 K 1 K 1 R w 21 = dT dt P 21 I= dQ dt iC= dv 21 dt dP 21 dt dv 2 dt dw 2 dt dt 2 dt FM= TJ= QC ? = qC t = i= F v 21 ?v 21 ?w 21 = T= QP 21 = 1 R ? q t 21 = 1 R t ELi 2 = E= F 2 K 1 2 1 2 T 2 K 1 2 1 2 1 R E= EIQ 2 1 2 Cv 2 = E= E= E= E= EC t t 2 = = ?v 21 ?w 21 = = P 21 = 1 R ? t 21 = 1 R t 21 1 2 Mv 2 2 v 2 21 1 2 J w 2 2 1 2 C ? P 2 21 2 2 2 v 2 v 1 C i v 2 v 2 v 1 v 1 R i ? P 2 F P 1 I Q v 2 v 1 K F w 2 w 1 K T v 2 v 1 = M constant T w 2 w 1 = J constant constant QP 1 P 2 C ? F q C t 2 2 = P 2 P 1 R ? Q 21 R t q w 2 w 1 ? T ? 2000 by CRC Press LLC The equations for this system might consist of a set of mesh or node equations. These could be reduced to a single nth-order linear ordinary differential equation by eliminating extraneous variables. (100.3) This nth-order equation can be replaced by an input-output transfer function (100.4) The inverse Laplace transform L –1 {H(s)} = h(t) is the system impulse response function. Alternatively, by selecting a set of n internal state variables, Eq.(100.3) can be written as a coupled set of first-order differential equations plus an algebraic equation relating the states to the original output y. These equations are called state equations, and one possible choice for this example is, assuming m = n, and y(t) = [100…0]x(t) + b n u(t) (100.5) In matrix notation these are written more succinctly as · x = Ax + Bu and y = Cx + Du (100.6) Any one of these six possible model forms, or others, might constitute the result of box 8 in Fig. 100.2. Discrete- time system models have similar choices of form, including an nth-order difference equation as given in Eq. (100.1) or a z-transform input-output transfer function as given in Eq. (100.2). A set of n first-order difference equations (state equations) analogous to Eq. (100.5) or (100.6) also can be written. Extensions to systems with r inputs and m outputs lead to a set of m coupled equations similar to Eq. (100.3), one for each output y i . These higher-order equations can be reduced to n first-order state differential equations and m algebraic output equations as in Eq. (100.5) or (100.6). The A matrix is again of dimension n ′ n, but B is now n ′ r, C is m ′ n, and D is m ′ r. In all previous discussions, the number of state variables, n, is the order of the model. In transfer function form, an m ′ r matrix H(s) of transfer functions will describe the input-output behavior Y(s) = H(s)U(s) (100.7) Other transfer function forms are also applicable, including the left and right forms of the matrix fraction description (MFD) of the transfer functions [Kailath, 1980] dy dt a dy dt a dy dt ay bub du dt b du dt n n n n n m m m ++++=+++ - - - 1 1 1 1001 LL Ys Us Hs bs bs bsb sas asa m m m m n n n () () ()== ++++ ++++ - - - - 1 1 10 1 1 10 L L ˙ () – – – – () – – – – xxt a a a a t bab bab bab bab n n nnn nnn n n = é ? ê ê ê ê ê ê ù ? ú ú ú ú ú ú + é ? ê ê ê ê ê ê ù - - -- -- 1 2 1 0 11 22 11 00 100 0 011 0 000 1 000 0 L L MMMMMM L L M ? ú ú ú ú ú ú ut() ? 2000 by CRC Press LLC H(s) = P(s) –1 N(s)orH(s) = N(s)P(s) –1 (100.8) Both P and N are matrices whose elements are polynomials in s. Very similar model forms apply to continuous- time and discrete-time systems, with the major difference being whether Laplace transform or z-transform transfer functions are involved. When time-variable systems are encountered, the option of using high-order differential or difference equations versus sets of first-order state equations is still open. The system coefficients a i (t), b j (t) and/or the matrices A(t), B(t), C(t), and D(t) will now be time-varying. Transfer function approaches lose most of their utility in time-varying cases and are seldom used. With nonlinear systems all the options relating to the order and number of differential or difference equation still apply. The form of the nonlinear state equations is · x = f(x, u, t) y = h(x, u, t) (100.9) where the nonlinear vector-valued functions f(x, u, t) and h(x, u, t) replace the right-hand sides of Eq. (100.6). The transfer function forms are of no value in nonlinear cases. Stochastic systems [Maybeck, 1979] are modeled in similar forms, except the coefficients of the model and/or the inputs are described in probabilistic terms. Nonuniqueness There is not a unique correct model of a given system for several reasons. The selection of idealized elements to represent the system requires judgment based upon the intended purpose. For example, a satellite might be modeled as a point mass in a study of its gross motion through space. A detailed flexible structure model might be required if the goal is to control vibration of a crucial on-board sensor. In empirical modeling, the assumed starting form, Eq. (100.1), can vary. There is a trade-off between the complexity of the model form and the fidelity with which it will match the data set. For example, a pth-degree polynomial can exactly fit to p + 1 data points, but a straight line might be a better model of the underlying physics. Deviations from the line might be caused by extraneous measure- ment noise. Issues such as these are addressed in Astrom [1980]. The preceding paragraph addresses nonuniqueness in determining an input-output system description. In addition, state models developed from input-output descriptions are not unique. Suppose the transfer function of a single-input, single-output linear system is known exactly. The state variable model of this system is not unique for at least two reasons. An arbitrarily high-order state variable model can be found that will have this same transfer function. There is, however, a unique minimal or irreducible order n min from among all state models that have the specified transfer function. A state model of this order will have the desirable properties of controllability and observability. It is interesting to point out that the minimal order may be less than the actual order of the physical system. The second aspect of the nonuniqueness issue relates not to order, i.e., the number of state variables, but to choice of internal variables (state variables). Mathematical and physical methods of selecting state variables are available [Brogan, 1991]. An infinite number of choices exist, and each leads to a different set {A, B, C, D}, called a realization. Some state variable model forms are more convenient for revealing key system properties such as stability, controllability, observability, stabilizability, and detectability. Common forms include the controllable canonical form, the observable canonical form, the Jordan canonical form, and the Kalman canonical form. The reverse process is unique in that every valid realization leads to the same model transfer function H(s) = C{sI – A} –1 B + D (100.10) ? 2000 by CRC Press LLC Approximation of Continuous Systems by Discrete Models Modern control systems often are implemented digitally, and many modern sensors provide digital output, as shown in Fig. 100.3. In designing or analyzing such systems discrete-time approximate models of continuous- time systems are frequently needed. There are several general ways of proceeding, as shown in Fig. 100.4. Many choices exist for each path on the figure. Alternative choices of states or of approximation methods, such as forward or backward differences, lead to an infinite number of valid models. Defining Terms Controllability: A property that in the linear system case depends upon the A,B matrix pair which ensures the existence of some control input that will drive any arbitrary initial state to zero in finite time. Detectability: A system is detectable if all its unstable modes are observable. Observability:A property that in the linear system case depends upon the A,C matrix pair which ensures the ability to determine the initial values of all states by observing the system outputs for some finite time interval. Stabilizable:A system is stabilizable if all its unstable modes are controllable. State variables:A set of variables that completely summarize the system’s status in the following sense. If all states x i are known at time t 0 , then the values of all states and outputs can be determined uniquely for any time t 1 > t 0 , provided the inputs are known from t 0 onward. State variables are components in the state vector. State space is a vector space containing the state vectors. FIGURE 100.3Digital output provided by modern sensor. FIGURE 100.4State variable modeling paradigm. Computer Zero order hold Continuous system u j (t k ) y (t k )y (t) Digital controller Sampling sensor Continuous-time model: differential equations [e.g., Eq. (93.3) or Transfer functions e.g., Eq. (93.4)] Approximation or sampling Approximation or sampling Select states Select states Discrete-time model: difference equations [e.g., Eq. (93.1) or Transfer functions e.g., Eq. (93.2)] Continuous-time state equations: x = A(t)x + B(t)u(t) y(t) = Cx(t) + D(t)u(t) Discrete-time state equations: y(k + 1) = A(k)x(k) + B(k)u(k) y(k) = C(k)x(k) + D(k)u(k) ? 2000 by CRC Press LLC Related Topic 6.1 Definitions and Properties References K.J. Astrom, “Maximum likelihood and prediction error methods,” Automatica, vol. 16, pp. 551–574, 1980. K.J. Astrom and B. Wittenmark, Adaptive Control, Reading, Mass.: Addison-Wesley, 1989. W.L. Brogan, “Optimal control theory applied to systems described by partial differential equations,” in Advances in Control Systems, vol. 6, C. T. Leondes (ed.), New York: Academic Press, 1968, chap. 4. W.L. Brogan, Modern Control Theory, 3rd ed., Englewood Cliffs, N.J.: Prentice-Hall, 1991. R.C. Dorf, Modern Control Systems, 5th ed., Reading, Mass.: Addison-Wesley, 1989. T. Kailath, Linear Systems, Englewood Cliffs, N.J.: Prentice-Hall, 1980. J. Makhoul, “Linear prediction: A tutorial review,” Proc. IEEE, vol. 63, no. 4, pp. 561–580, 1975. P.S. Maybeck, Stochastic Models, Estimation and Control, vol. 1, New York: Academic Press, 1979. J.L. Shearer, A.T. Murphy, and H.H. Richardson, Introduction to Dynamic Systems, Reading, Mass.: Addison- Wesley, 1967. Further Information The monthly IEEE Control Systems Magazine frequently contains application articles involving models of interesting physical systems. The monthly IEEE Transactions on Automatic Control is concerned with theoretical aspects of systems. Models as discussed here are often the starting point for these investigations. Automatica is the source of many related articles. In particular an extended survey on system identification is given by Astrom and Eykhoff in vol. 7, pp. 123–162, 1971. Early developments of the state variable approach are given by R. E. Kalman in “Mathematical description of linear dynamical systems,” SIAM J. Control Ser., vol. A1, no. 2, pp. 152–192, 1963. 100.2 Dynamic Response Gordon K. F. Lee Computing the Dynamic System Response Consider a linear time-invariant dynamic system represented by a differential equation form (100.11) where y(t) and f(t) represent the output and input, respectively, of the system. Let p k (·) D = (d k /dt k )(·) define the differential operator so that (100.11) becomes (p n + a n–1 p n–1 + … + a 1 p + a 0 )y(t) = (b m p m + … + b 1 p + b 0 )f (t) (100.12) The solution to (100.11) is given by y(t) = y S (t) + y I (t) (100.13) dyt dt a dyt dt a dy t dt ayt b dft dt b df t dt bft n n n n n m m m () () () () () () () ++++ =+++ - - - 1 1 1 10 10 L L ? 2000 by CRC Press LLC where y S (t) is the zero-input response, or that part of the response due to the initial conditions (or states) only, and y I (t) is the zero-state response, or that part of the response due to the input f(t) only. Zero-Input Response: y S (t) Here f(t) = 0, and thus (100.11) becomes (p n + a n–1 p n–1 + … + a 1 p + a 0 )y(t) = 0 (100.14) That is, D(p)y(t) = 0 The roots of D(p) = 0 can be categorized as either distinct or multiple. That is, in general, where there are r distinct roots and q sets of multiple roots (each set has multiplicity k i ). Note that r + s = n, where s D = ? q i=1 k i . Each distinct root contributes a term to y S (t) of the form c i e lit , where c i is a constant, while each set of multiple roots contributes a set of terms to y S (t) of the form ? j=0 ki–1 c i,j t j e lit , where c i,j is some constant. Thus, the zero-input response is given by (100.15) The coefficients c i,j and c s+i are selected to satisfy the initial conditions. Special Case If all the roots of D(p) = 0 are distinct and the initial conditions for (100.11) are given by then the coefficients of (100.15) are given by the solution of (100.16) Zero-State Response: y I (t) Here the initial conditions are made identically zero. Observing (100.11), let Dp p p i k i q qi i r i () ( ) ( ) =- - = + = ?? ll 11 yt cte ce Sij jt j k i q i t i r i i i () , =+ = - = + = ??? + l s l s 0 1 11 y dy dt dy dt n n (), () , () 0 00 1 1 . . . , - - ì í ? ? ? ü y ? t ? y dy dt dy dt c c c n n n nn n n n () () () ()() () 0 0 0 11 1 1 1 12 1 1 2 11 1 2 M L L MMMM L M - - -- - é ? ê ê ê ê ê ê ê ù ? ú ú ú ú ú ú ú = é ? ê ê ê ê ê ù ? ú ú ú ú ú é ? ê ê ê ê ê ù ? ú ú ú ú ú ll l ll l Hp bp bpb pap apa m m n n n ()= +++ ++++ - - L L 10 1 1 10 ? 2000 by CRC Press LLC denote a rational function in the p operator. Consider using partial-fraction expansion on H(p) as (100.17) when the first term corresponds to the sets of multiple roots and the second term corresponds to the distinct roots. Note the constant residuals are computed as g s+i = [(p – l q+i )H(p)] p = lq+i and Then (100.18) is the impulse response of the system (100.11). Then the zero-state response is given by (100.19) that is, y I (t) is the time convolution between input f(t) and impulse response h(t). In some instances, it may be easier to find y s (t) and y I (t) using Laplace Transform methods. Measures of the Dynamic System Response Several measures may be employed to investigate dynamic response performance. These include: 1.Speed of the response—how quickly does the system reach its final value 2.Accuracy—how close is the final response to the desired response 3.Relative stability—how stable is the system or how close is the system to instability 4.Sensitivity—what happens to the system response if the system parameters change Objectives 3 and 4 may be analyzed by frequency domain methods (Section 100.3). Time-domain measures classically analyze the dynamic response by partitioning the total response into its steady-state (objective 2) and transient (objective 1) components. The steady-state response is that part of the response which remains as time approaches infinity; the transient response is that part of the response which vanishes as time approaches infinity. Measures of the Steady-State Response In the steady state, the accuracy of the time response is an indication of how well the dynamic response follows a desired time trajectory. Usually a test signal (reference signal) is selected to measure accuracy. Consider Fig. 100.5. In this configuration, the objective is to force y(t) to track a reference signal r(t) as close as possible. Hp g p g p ij i j j k i q i qii ri () () , = - + - == + += ??? l l s 111 g kj d dp pHp ij i kj kj i k p i i i i , () () ()! ()()= - - {} - - = 1 l l ht g j te ge ij j t j k i q i i r t i i i () ()! , = - + - == + = ??? + 1 1 111 l s l s yt fht d I t () ()( )=- ò ttt 0 ? 2000 by CRC Press LLC The steady-state error is a measure of the accuracy of the output y(t) in tracking the reference input r(t). Other configurations with different performance measures would result in other definitions of the steady-state error between two signals. From Fig. 100.5, the error e(t) is e(t) = r(t) – y(t) (100.20) and the steady-state error is (100.21) assuming the limits exists, where E(s) is the Laplace transform of e(t), and s is the Laplacian operator. With G(s) the transfer function of the system and H(s) the transfer function of the controller, the transfer function between y(t) and r(t) is found to be (100.22) with (100.23) Direct application of the steady-state error for various inputs yields Table 100.2. Note u(t) is the unit step function. This table can be extended to an mth-order input in a straightforward manner. Note that for e SS (t) to go to zero with a reference signal Ct m u(t), the term G(s)H(s) must have at least m poles at the origin (a type m system). FIGURE 100.5A tracking controller configuration. TABLE 100.2Steady-State Error Constants r(t) Error Test Signal eSS(t) Constant Rut Rtut R tut () (): (): :step function ramp function parabolic function 2 2 R K R K R K p v a 1+ KGsHs KsGsHs KsGsHs p s v s a s lim()() lim ()() lim ()() 0 0 0 2 = = = ? ? ? e t et sEs SS ts () () ()== ?¥ ?¥ lim lim Ts GsHs GsHs () ()() ()() = +1 Es Rs GsHs () () ()() = +1 ? 2000 by CRC Press LLC Measures of the Transient Response In general, analysis of the transient response of a dynamic system to a reference input is difficult. Hence formulating a standard measure of performance becomes complicated. In many cases, the response is dominated by a pair of poles and thus acts like a second-order system. Consider a reference unit step input to a dynamic system (Fig. 100.6). Critical parameters that measure transient response include: 1.M: maximum overshoot 2.% overshoot = M/A ′ 100%, where A is the final value of the time response 3.t d : delay time—the time required to reach 50% of A 4.t r : rise time—the time required to go from 10% of A to 90% of A 5.t s : settling time—the time required for the response to reach and stay within 5% of A To calculate these measures, consider a second-order system (100.24) where x is the damping coefficient and w n is the natural frequency of oscillation. For the range 0 < x < 1, the system response is underdamped, resulting in a damped oscillatory output. For a unit step input, the response is given by (100.25) The eigenvalues (poles) of the system [roots of the denominator of T(s)] provide some measure of the time constants of the system. For the system under study, the eigenvalues are at FIGURE 100.6Step response. Ts ss n nn ()= ++ w xw w 2 22 2 yt e t n t n () – ( )=+ - -- - ? è ? ? ? ? ÷ ÷ << - 1 1 1 1 01 2 21 2 xw x wx x x xsin tan – –xw w x nn jj±- =-11 2 where D ? 2000 by CRC Press LLC From the expression of y(t), one sees that the term xw n affects the rise time and exponential decay time. The effects of the damping coefficient on the transient response are seen in Fig. 100.7. The effects of the natural frequency of oscillation w n of the transient response can be seen in Fig. 100.8. As w n increases, the frequency of oscillation increases. For the case when 0 < x < 1, the underdamped case, one can analyze the critical transient response parameters. FIGURE 100.7Effect of the damping coefficient on the dynamic response. FIGURE 100.8Effect of the natural frequency of oscillation on the dynamic response. ? 2000 by CRC Press LLC To measure the peaks of Fig. 100.6, one finds (100.26) occurring at (100.27) Hence (100.28) occurring at (100.29) With these parameters, one finds and Note that increasing x decreases the % overshoot and decreases the settling time but increases t d and t r . When x = 1, the system has a double pole at –w n , resulting in a critically damped response. This is the point when the response just changes from oscillatory to exponential in form. For a unit step input, the response is given by y(t) = 1 – e –wnt (1 + w n t)(x = 1) (100.30) For the range x > 1, the system is overdamped due to two real system poles. For a unit step input, the response is given by yt n n n peak exp – . . .( ) ( ) , ,=+- - = - 11 1 01 1 2 px x t n n n n = - p wx1 2 : odd (overshoot) even (undershoot): y max exp=+ - 1 1 2 –px x t n max = - p wx1 2 t t t d n r n s n ? + ? ++ ? 107 111 14 3 2 . .. x w xx w xw yt ccc e c e ct c t nn ( ) ( )=+ - - ? è ? ? ? ÷ >1 11 1 1 121 2 12 ww x ? 2000 by CRC Press LLC (100.31) Finally, when x = 0, the response is purely sinusoidal. For a unit step, the response is given by y(t) = 1 – cos w n t (x = 0) (100.32) Defining Terms Impulse response:The response of a system when the input is an impulse function. Steady-state error: The difference between the desired reference signal and the actual signal in steady-state, i.e., when time approaches infinity. Steady-state response:That part of the response which remains as time approaches infinity. Transient response:That part of the response which vanishes as time approaches infinity. Zero-input response:That part of the response due to the initial condition only. Zero-state response:That part of the response due to the input only. Related Topics 6.1 Definitions and Properties?7.1 Introduction?112.2 A Brief History of CACSD Further Information J.J. D’Azzo and C.H. Harpis, Linear Control System Analysis and Design, New York: McGraw-Hill, 1981. R.C. Dorf, Modern Control Systems, 5th ed., Reading, Mass.: Addison-Wesley, 1989. M.E. El-Hawary, Control Systems Engineering, Reston, Va.: Reston, 1984. G.H. Hostetter, C. J. Savant, Jr., and R. T. Stefani, Design of Feedback Control Systems, Philadelphia: Saunders, 1989. B.C. Kuo, Automatic Control Systems, Englewood Cliffs, N.J.: Prentice-Hall, 1987. K. Ogata, Modern Control Engineering, Englewood Cliffs, N.J.: Prentice-Hall, 1970. N.K. Sinha, Control Systems, New York: Holt, 1986. 100.3 Frequency Response Methods: Bode Diagram Approach Andrew P. Sage Our efforts in this section are concerned with analysis and design of linear control systems by frequency response methods. Design generally involves trial-and-error repetition of analysis until a set of design specifications has been met. Thus, analysis methods are most useful in the design process, which is one phase of the systems engineering life cycle [Sage, 1992]. We will discuss one design method based on Bode diagrams. We will discuss the use of both simple series equalizers and composite equalizers as well as the use of minor-loop feedback in systems design. Figure 100.9 presents a flowchart of the frequency response method design process and indicates the key role of analysis in linear systems control design. The flowchart of Fig. 100.9 is applicable to control system design methods in general. There are several iterative loops, generally calling for trial-and-error efforts, that comprise the suggested design process. An experienced designer will often be able, primarily due to successful prior experience, to select a system structure and generic components such that the design specifications can be met with no or perhaps a very few iterations through the iterative loop involving adjustment of equalizer or compensation parameters to best meet specifications. If the parameter optimization, or parameter refinement such as to lead to maximum phase margin, approach shows the specifications cannot be met, we are then assured that no equalizer of the specific form selected will meet specifications. The next design step, if needed, would consist of modification of the equalizer form or cc 1 2 2 2 11=-+ - =-- -xx xx ? 2000 by CRC Press LLC structure and repetition of the analysis process to determine equalizer parameter values to best meet specifications. If specifications still cannot be met, we will usually next modify generic fixed components used in the system. This iterative design and analysis process is again repeated. If no reasonable fixed components can be obtained to meet specifications, then structural changes in the proposed system are next contemplated. If no structure can be found that allows satisfaction of specifications, either the client must be requested to relax the frequency response specifications or the project may be rejected as infeasible using present technology. As we might suspect, economics will play a dominant role in this design process. Changes made due to iteration in the inner loops of Fig. 100.9 normally involve little additional costs, whereas those made due to iterations in the outer loops will often involve major cost changes. Frequency Response Analysis Using the Bode Diagram The steady-state response of a stable linear constant-coefficient system has particular significance, as we know from an elementary study of electrical networks and circuits and of dynamics. We consider a stable linear system with input-output transfer function We assume a sinusoidal input u(t) = cos wt so that we have for the Laplace transform of the system output We expand this ratio of polynomials using the partial-fraction approach and obtain In this expression, F(s) contains all the poles of H(s). All of these lie in the left half plane since the system, represented by H(s), is assumed to be stable. The coefficients a 1 and a 2 are easily determined as FIGURE 100.9 System design life cycle for frequency-response-based design. Hs Zs Us () () () = Zs sH s s () () = + 22 w Zs Fs a sj a sj () ()=+ + + - 12 ww ? 2000 by CRC Press LLC We can represent the complex transfer function H(jw) in either of two forms, H(jw) = B(w) + jC(w) H(–jw) = B(w) – jC(w) The inverse Laplace transform of the system transfer function will result in a transient term due to the inverse transform of F(s), which will decay to zero as time progresses. A steady-state component will remain, and this is, from the inverse transform of the system equation, given by z(t) = a 1 e –jwt + a 2 jwt We combine several of these relations and obtain the result This result becomes, using the Euler identity, 1 z(t) = B(w) coswt – C(w) sinwt = [B 2 (w) + C 2 (w)] 1/2 cos(w + b) = *H(jw)* cos(wt + b) where tan b(w) = C(w)/B(w). As we see from this last result, there is a very direct relationship between the transfer function of a linear constant-coefficient system, the time response of a system to any known input, and the sinusoidal steady-state response of the system. We can always determine any of these if we are given any one of them. This is a very important result. This important conclusion justifies a design procedure for linear systems that is based only on sinusoidal steady-state response, as it is possible to determine transient responses, or responses to any given system input, from a knowledge of steady-state sinusoidal responses, at least in theory. In practice, this might be rather difficult computationally without some form of automated assistance. Bode Diagram Design-Series Equalizers In this subsection we consider three types of series equalization: 1.Gain adjustment, normally attenuation by a constant at all frequencies 2.Increasing the phase lead, or reducing the phase lag, at the crossover frequency by use of a phase lead network 3.Attenuation of the gain at middle and high frequencies such that the crossover frequency will be decreased to a lower value where the phase lag is less, by use of a lag network 1 The Euler identity is e jwt = cos wt + jsinwt. a Hj a Hj 1 2 2 2 = - = () () w w zt B ee C ee j jt jt jt jt () () ()= + ? è ? ? ? ÷ - - ? è ? ? ? ÷ -- ww ww ww 22 ? 2000 by CRC Press LLC In the subsection that follows this, we will first consider use of a composite or lag-lead network near crossover to attenuate gain only to reduce the crossover frequency to a value where the phase lag is less. Then we will consider more complex composite equalizers and state some general guidelines for Bode diagram design. Here, we will use Bode diagram frequency domain design techniques to develop a design procedure for each of three elementary types of series equalization. Gain Reduction Many linear control systems can be made sufficiently stable merely by reduction of the open-loop system gain to a sufficiently low value. This approach ignores all performance specifications, however, except that of phase margin (PM) and is, therefore, usually not a satisfactory approach. It is a very simple one, however, and serves to illustrate the approach to be taken in more complex cases. The following steps constitute an appropriate Bode diagram design procedure for compensation by gain adjustment: 1.Determine the required PM and the corresponding phase shift b c = –p + PM. 2.Determine the frequency w c at which the phase shift is such as to yield the phase shift at crossover required to give the desired PM. 3.Adjust the gain such that the actual crossover frequency occurs at the value computed in step 2. Phase-Lead Compensation In compensation using a phase-lead network, we increase the phase lead at the crossover frequency such that we meet a performance specification concerning phase shift. A phase-lead-compensating network transfer function is Figure 100.10 illustrates the gain versus frequency and phase versus frequency curves for a simple lead network with the transfer function of the foregoing equation. The maximum phase lead obtainable from a phase-lead network depends upon the ratio w 2 /w 1 that is used in designing the network. From the expression for the phase shift of the transfer function for this system, which is given by FIGURE 100.10Phase shift and gain curves for a simple lead network. Gs ss c () =+ ? è ? ? ? ÷ + ? è ? ? ? ÷ <11 12 12 ww ww ? 2000 by CRC Press LLC we see that the maximum amount of phase lead occurs at the point where the first derivative with respect to frequency is zero, or or at the frequency where w m = (w 1 w 2 ) 0.5 This frequency is easily seen to be at the center of the two break frequencies for the lead network on a Bode log asymptotic gain plot. It is interesting to note that this is exactly the same frequency that we would obtain using an arctangent approximation 2 with the assumption that w 1 < w < w 2 . There are many ways of realizing a simple phase-lead network. All methods require the use of an active element since the gain of the lead network at high frequencies is greater than 1. A simple electrical network realization is shown in Fig. 100.11. We now consider a simple design example. Suppose that we have an open-loop system with transfer function It turns out that this is often called a type-two system due to the presence of the double integration. This system will have a steady-state error of zero for a constant acceleration input. The crossover frequency, that is to say 2 The arctangent approximation is tan –1 (w/a) = w/a for w < a and tan –1 (w/a) = p/2 – a/w for w > a. This approximation is rather easily obtained through use of a Taylor series approximation. FIGURE 100.11A simple electrical lead network. b w w w w =- -- tan tan 1 1 1 2 d d m b w ww= =0 Gs s f ()= 10 4 2 ? 2000 by CRC Press LLC the frequency where the magnitude of the open-loop gain is 1, is 100 rad/s. The PM for the system without equalization is zero. We will design a simple lead network compensation for this zero PM system. If uncom- pensated, the closed-loop transfer function will be such that the system is unstable and any disturbance at all will result in a sinusoidally oscillating output. The asymptotic gain diagram for this example is easily obtained from the open-loop transfer function and we wish to select the break frequencies w 1 and w 2 such that the phase shift at crossover is maximum. Further, we want this maximum phase shift to be such that we obtain the specified PM. We use the procedure suggested in Fig. 100.12. Since the crossover frequency is such that w 1 < w c < w 2 , we have for the arctangent approximation to the phase shift in the vicinity of the crossover frequency In order to maximize the phase shift at crossover, we set and obtain as a result w m = (w 1 w 2 ) 0.5 FIGURE 100.12Life cycle of frequency domain design incorporating lead network compensation. GsGs Ks ss fc ()() () () = + + 1 1 1 2 2 / / w w bw p w w w w pw w w w () tan tan =-+ - ? - -- --1 1 1 2 1 2 2 d d m b w ww= =0 ? 2000 by CRC Press LLC We see that the crossover frequency obtained is halfway between the two break frequencies w 1 and w 2 on a logarithmic frequency coordinate. The phase shift at this optimum value of crossover frequency becomes For a PM of –3p/4, for example, we have –3p/4 = –p/2 – 2(w 1 /w 2 ) 0.5 , and we obtain w 1 /w 2 = 0.1542 as the ratio of frequencies. We see that we have need for a lead network with a gain of w 1 /w 2 = 6.485. The gain at the crossover frequency is 1, and from the asymptotic gain approximation that is valid for w 1 < w < w 2 , we have the expressions *G(jw)* = K/ww 1 and *G(jw c )* = 1 = K/w c w 1 which for a known K can be solved for w c and w 1 . Now that we have illustrated the design computation with a very simple example, we are in a position to state some general results. In the direct approach to design for a specified PM we assume a single lead network equalizer such that the open-loop system to transfer function results. This approach to design results in the following steps that are applicable for Bode diagram design to achieve maximum PM within an experientially determined control system structure that comprises a fixed plant and a compensation network with adjustable parameters: 1.We find an equation for the gain at the crossover frequency in terms of the compensated open-loop system break frequency. 2.We find an equation of the phase shift at crossover. 3.We find the relationship between equalizer parameters and crossover frequency such that the phase shift at crossover is the maximum possible and a minimum of additional gain is needed. 4.We determine all parameter specifications to meet the PM specifications. 5.We check to see that all design specifications have been met. If they have not, we iterate the design process. Figure 100.12 illustrates the steps involved in implementing this frequency domain design approach. Phase-Lag Compensation In the phase-lag-compensation frequency domain design approach, we reduce the gain at low frequencies such that crossover, the frequency where the gain magnitude is 1, occurs before the phase lag has had a chance to become intolerably large. A simple single-stage phase-lag-compensating network transfer function is Figure 100.13 illustrates the gain and phase versus frequency curves for a simple lag network with this transfer function. The maximum phase lag obtainable from a phase-lag network depends upon the ratio w 2 /w 1 that is used in designing the network. From the expression for the phase shift of this transfer function, we see that maximum phase lag occurs at that frequency w c where db/dw = 0. We obtain for this value w m = (w 1 w 2 ) 0.5 which is at the center of the two break frequencies for the lag network when the frequency response diagram is illustrated on a Bode diagram log-log asymptotic gain plot. The maximum value of the phase lag obtained at w = w m is bbw pw w cc == - - ? è ? ? ? ÷ () . 2 2 1 2 05 Gs s s c () = + + < 1 1 2 1 12 / / w w ww b w w w w =- -- tan tan 1 2 1 1 ? 2000 by CRC Press LLC which can be approximated in a more usable form, using the arctangent approximation, as The attenuation of the lag network at the frequency of minimum phase shift, or maximum phase lag, is obtained from the asymptotic approximation as Figure 100.13 presents a curve of attenuation magnitude obtainable at the frequency of maximum phase lag and the amount of the phase lag for various ratios w 2 /w 1 for this simple lag network. There are many ways to physically realize a lag network transfer function. Since the network only attenuates at some frequencies, as it never has a gain greater than 1 at any frequency, it can be realized with passive components only. Figure 100.14 presents an electrical realization of the simple lag network. Figure 100.15 presents a flowchart illustrating the design procedure envisioned here for lag network design. This is concep- tually very similar to that for a lead network and makes use of the five-step parameter optimization procedure suggested earlier. The object of lag network design is to reduce the gain at frequencies lower than the original crossover frequency in order to reduce the open-loop gain to unity before the phase shift becomes so excessive that the system PM is too small. A disadvantage of lag network compensation is that the attenuation introduced reduces FIGURE 100.13Phase shift and gain curves for a simple lag network. bw pw w pw w mm ( ) tan tan . . =- ? è ? ? ? ÷ =- ? è ? ? ? ÷ - - 2 2 2 2 1 2 1 05 1 1 2 05 bw pw w mm ()? 2 2 1 **G cm () . w w w = ? è ? ? ? ÷ 1 2 05 ? 2000 by CRC Press LLC the crossover frequency and makes the system slower in terms of its transient response. Of course, this would be advantageous if high-frequency noise is present and we wish to reduce its effect. The lag network is an entirely passive device and thus is more economical to instrument than the lead network. In lead network compensation we actually insert phase lead in the vicinity of the crossover frequency to increase the PM. Thus we realize a specified PM without lowering the medium-frequency system gain. We see that the disadvantages of the lag network are the advantages of the lead network and the advantages of the lag network are the disadvantages of the lead network. We can attempt to combine the lag network with the lead network into an all-passive structure called a lag- lead network. Generally we obtain better results than we can achieve using either a lead or a lag network. We will consider design using lag-lead networks in our next subsection as well as more complex composite equalization networks. Composite Equalizers In the previous subsection we examined the simplest forms of series equalization: gain adjustment, lead network compensation, and lag network compensation. In this subsection we will consider more complex design examples in which composite equalizers will be used for series compensation. The same design principles used earlier in this section will be used here as well. FIGURE 100.14 A simple electrical lag network. FIGURE 100.15 Life cycle of frequency domain design incorporating lag network compensation. ? 2000 by CRC Press LLC Lag-Lead Network Design The prime purpose of a lead network is to add phase lead near the crossover frequency to increase the PM. Accompanied with this phase lead is a gain increase that will increase the crossover frequency. This will sometimes cause difficulties if there is much phase lag in the uncompensated system at high frequencies. There may be situations where use of a phase-lead network to achieve a given PM is not possible due to too many high-frequency poles. The basic idea behind lag network design is to reduce the gain at “middle” frequencies such as to reduce the crossover frequency to a lower value than for the uncompensated system. If the phase lag is less at this lower frequency, then the PM will be increased by use of the lag network. We have seen that is not possible to use a lag network in situations in which there is not a frequency where an acceptable PM would exist if this frequency were the crossover frequency. Even if use of a lag network is possible, the significantly reduced crossover frequency resulting from its use may make the system so slow and sluggish in response to an input that system performance is unacceptable even though the relative stability of the system is acceptable. Examination of these characteristics or attributes of lead network and lag network compensation suggests that it might be possible to combine the two approaches to achieve the desirable features of each approach. Thus we will attempt to provide attenuation below the crossover frequency to decrease the phase lag at crossover and phase lead closer to the crossover frequency in order to increase the phase lead of the uncompensated system at the crossover frequency. The transfer function of the basic lag-lead network is where w 4 > w 3 > w 2 > w 1 . Often it is desirable that w 2 w 3 = w 1 w 4 such that the high-frequency gain of the equalizer is unity. It is generally not desirable that w 1 w 4 > w 2 w 3 as this indicates a high-frequency gain greater than 1, and this will require an active network, or gain, and a passive equalizer. It is a fact that we should always be able to realize a linear minimum phase network using passive components only if the network has a rational transfer function with a gain magnitude that is no greater than 1 at any real frequency. Figure 100.16 illustrates the gain magnitude and phase shift curves for a single-stage lag-lead network equalizer or compensator transfer function. Figure 100.17 illustrates an electrical network realization of a passive lag-lead network equalizer. Parameter matching can be used to determine the electrical network parameters that yield a specified transfer function. Because the relationships between the break frequencies and the equalizer component values are complex, it may be desirable, particularly in preliminary instrumentation of the control system, to use analog or digital computer programming techniques to construct the equalizer. Traditionally, there has been much analog computer simulation of control systems. The more contemporary approach suggests use of digital computer approaches that require numerical approximation of continuous-time physical systems. Figure 100.18 presents a flowchart that we may use for lag-lead network design. We see that this flowchart has much in common with the charts and design procedures for lead network and lag network design and that each of these approaches first involves determining or obtaining a set of desired specifications for the control system. Next, the form of a trial compensating network and the number of break frequencies in the network are selected. We must then obtain a number of equations, equal to the number of network break frequencies plus 1. One of these equations shows that the gain magnitude is 1 at the crossover frequency. The second equation will be an equation for the phase shift at crossover. It is generally desirable that there be at least two unspecified compensating network break frequencies such that we may use a third equation, the optimality of the phase shift at crossover equation, in which we set db/dwU w=wc = 0. If other equations are needed to represent the design situation, we obtain these from the design specifications themselves. General Bode Diagram Design Figure 100.19 presents a flowchart of a general design procedure for Bode diagram design. As we will see in the next subsection, a minor modification of this flowchart can be used to accomplish design using minor-loop feedback Gs ss ss c () ()() ()() = ++ ++ 11 11 23 14 // // ww ww ? 2000 by CRC Press LLC or a combination of minor-loop and series equations. These detailed flowcharts for Bode diagram design are, of course, part of the overall design procedure of Fig. 100.9. Much experience leads to the conclusion that satisfactory linear systems control design using frequency response approaches is such that the crossover frequency occurs on a gain magnitude curve which has a –1 slope at the crossover frequency. In the vicinity of crossover we may approximate any minimum phase transfer function, with crossover on a –1 slope, by FIGURE 100.16Phase shift and gain curves for a simple lag-lead network. FIGURE 100.17Simple electrical lag-lead network. Gs GsGs s ss fc c nn nm c () ()() () () == + + >> -- - ww w w www 1 1 1 1 2 1 12 1 1 / / for ? 2000 by CRC Press LLC Here w 1 is the break frequency just prior to crossover and w 2 is the break frequency just after crossover. It is easy to verify that we have *G(jw c )* = 1 if w 1 > w c > w 2 . Figure 100.20 illustrates this rather general approximation to a compensated system Bode diagram in the vicinity of the crossover frequency. We will conclude this subsection by determining some general design requirements for a system with this transfer function and the associated Bode asymptotic gain magnitude diagram of Fig. 100.20. There are three unknown frequencies in the foregoing equation. Thus we need three requirements or equations to determine design parameters. We will use the same three requirements used thus far in all our efforts in this section, namely: FIGURE 100.18Life cycle of frequency domain design incorporating lag-lead network compensation. FIGURE 100.19Life cycle of frequency domain design incorporating general Bode diagram compensation procedure. ? 2000 by CRC Press LLC 1. The gain at crossover is 1. 2. The PM is some specified value. 3. The PM at crossover is the maximum possible for a given w 2 /w 1 ratio. We see that the first requirement, that the gain is 1 at the crossover frequency, is satisfied by the foregoing equation if the crossover frequency occurs on the –1 slope portion of the gain curve as assumed in Fig. 100.20. We use the arctangent approximation to obtain the phase shift in the vicinity of crossover as To satisfy requirement 3 we set and obtain as the optimum setting for the crossover frequency. Substitution of the “optimum” frequency given by the foregoing into the phase shift equation results in We desire a specific PM here, and so the equalizer break frequency locations are specified. There is a single parameter here that is unspecified, and an additional equation must be found in any specific application. Alternately, we could simply assume a nominal crossover frequency of unity or simply normalize frequencies w 1 and w 2 in terms of the crossover frequency by use of the normalized frequencies w 1 = W 1 w c and w 2 = W 2 w c . FIGURE 100.20 Illustration of generic gain magnitude in the vicinity of crossover. bw p pw w w w () ( ) ( )=- + - - ? è ? ? ? ÷ -- n nm 2 1 2 1 1 2 d d n m c c bw w w w w ww () ( ) = == - - - 0 1 1 1 2 2 www c n m 2 12 1 1 = - - bw pw w () ( )( ) c mn= - --- 2 211 1 2 ? 2000 by CRC Press LLC It is a relatively simple matter to show that for a specified PM expressed in radians, we obtain for the normalized break frequencies It is relatively easy to implement this suggested Bode diagram design procedure which is based upon considering only these break frequencies immediately above and below crossover and which approximate all others. Break frequencies far below crossover are approximated by integrations or differentiations, that is, poles or zeros at s = 0, and break frequencies far above the crossover frequency are ignored. Minor-Loop Design In our efforts thus far in this section we have assumed that compensating networks would be placed in series with the fixed plant and then a unity feedback ratio loop closed around these elements to yield the closed-loop system. In many applications it may be physically convenient, perhaps due to instrumentation considerations, to use one or more minor loops to obtain a desired compensation of a fixed plant transfer function. For a single-input–single-output linear system there are no theoretical advantages whatever to any minor- loop compensation to series compensation as the same closed-loop transfer function can be realized by all procedures. However, when there are multiple inputs or outputs, then there may be considerable advantages to minor-loop design as contrasted to series compensation design. Multiple inputs often occur when there is a single-signal input and one or more noise or disturbance inputs present and a task of the system is to pass the signal inputs and reject the noise inputs. Also there may be saturation-type nonlinearities present, and we may be concerned not only with the primary system output but also with keeping the output at the saturation point within bounds such that the system remains linear. Thus there are reasons why minor-loop design may be preferable to series equalization. We have discussed block diagrams elsewhere in this handbook. It is desirable here to review some concepts that will be of value for our discussion of minor-loop design. Figure 100.21 illustrates a relatively general linear control system with a single minor loop. This block diagram could represent many simple control systems. G 1 (s) could represent a discriminator and series compensation and G 2 (s) could represent an amplifier and that part of a motor transfer function excluding the final integration to convert velocity to position. G 3 (s) might then represent an integrator. G 4 (s) would then represent a minor-loop compensation transfer function, such as that of a tachometer. The closed-loop transfer function for this system is given by FIGURE 100.21Feedback control system with a single minor loop and output disturbance. W n W m c c 1 1 2 2 21 21 == - == - w w w w PM PM () () Zs Us Hs GsGsGs GsGsGsGsGs () () () ()()() ()() ()()() == ++ 123 24 123 1 ? 2000 by CRC Press LLC It is convenient to define several other transfer functions that are based on the block diagram in Fig. 100.21. First there is the minor-loop gain G m (s) = G 2 (s)G 4 (s) which is just the loop gain of the minor loop only. The minor loop has the transfer function There will usually be a range or ranges of frequency for which the minor-loop gain magnitude is much less than 1, and we then have There will also generally be ranges of frequency for which the minor-loop gain magnitude is much greater than 1, and we then have We may use these two relations to considerably simplify our approach to the minor-loop design problem. For illustrative purposes, we will use two major-loop gain functions. First we will consider the major-loop gain with the minor-loop-compensating network removed such that G 4 (s) = 0. This represents the standard situation we have examined in the last subsection. This uncompensated major-loop transfer function is G Mu (s) = G 1 (s)G 2 (s)G 3 (s) With the minor-loop compensation inserted, the major-loop gain, the input-output transfer function with the unity ratio feedback open, is We may express the complete closed-loop transfer function in the form A particularly useful relationship may be obtained by combining the last three equations into one equation of the form Zs Us Hs Gs GsGs Gs Gs m m m m () () () () () () () () == + = + 2 24 2 11 Zs Us Hs Gs G m m mm () () () () ( )=? < 2 1**w Zs Us Hs Gs G m m mm () () () () ()=? > 1 1 4 **w Gs GsG sGs Gs Mc m () () () () () = + 12 3 1 Zs Us Hs Gs Gs Mc Mc () () () () () == +1 Gs Gs Gs Mc Mu m () () () = +1 ? 2000 by CRC Press LLC We may give this latter equation a particularly simple interpretation. For frequencies where the minor-loop gain G m (s) is low, the minor-loop–closed major-loop transfer function G Mc (s) is approximately that of the minor-loop–open major-loop transfer function G Mu in that G Mc (s) ? G Mu (s) *G m (w)* << 1 For frequencies where the minor-loop gain G m (s) is high, the minor-loop–closed major-loop transfer function is just This has an especially simple interpretation on the logarithmic frequency plots we use for Bode diagrams for we may simply subtract the minor-loop gain G m (s) from the minor-loop–open major-loop gain G Mu (s) to obtain the compensated system gain as the transfer function G Mc (s). The last several equations are the key relations for minor-loop design using this frequency response approach. These relations indicate that some forms of series compensation yield a given major-loop transfer function G Mc (s) which will not be appropriate for realization by minor-loop compensation. In particular, a lead network series compensation cannot be realized by means of equivalent minor-loop compensation. The gain of the fixed plant G Mu (s) is raised at high frequencies due to the use of a lead network compensation. Also, we see that G Mc (s) can only be lowered by use of a minor-loop gain G m (s). A lag network used for series compensation will result in a reduction in the fixed plant gain *G Mu (w)* at all high frequencies. This can only be achieved if the minor-loop transfer gain G m (s) is constant for high frequencies. In some cases this may be achievable but often will not be. It is possible to realize the equivalent of lag network series equalization by means of a minor-loop equalization for systems where the low- and high-frequency behavior of G Mu (s), or G f (s), and G Mc (s) are the same and where the gain magnitude of the compensated system *G Mc (s)* is at no frequency any greater than is the gain magnitude of the fixed plant *G f (s)* or the minor-loop–open major-loop transfer function *G Mu (s)*. Thus we see that lag-lead network series equalization is an ideal type of equalization to realize by means of equivalent minor-loop equalization. Figures 100.9 and 100.19 represent flowcharts of a suggested general design procedure for minor-loop compensator design as well as for the series equalization approaches we examined previously. In our work thus far we have assumed that parameters were constant and known. Such is seldom the case, and we must naturally be concerned with the effects of parameter variations, disturbances, and nonlinearities upon system performance. Suppose, for example, that we design a system with a certain gain assumed as K 1 . If the system operates open loop and the gain K 1 is in cascade or series with the other input-output components, then the overall transfer function changes by precisely the same factor as K 1 changes. If we have an amplifier with unity ratio feedback around a gain K 1 , the situation is much different. The closed-loop gain would nominally be K 1 /(1 + K 1 ), and a change to 2K 1 would give a closed-loop gain 2K 1 /(1 + 2K 1 ). If K 1 is large, say 10 3 , then the new gain is 0.99950025, which is a percentage change of less than 0.05% for a change in gain of 100%. Another advantage of minor-loop feedback occurs when there are output disturbances such as those due to wind gusts on an antenna. We consider the system illustrated in Fig. 100.21. The response due to D(s) alone is When we use the relation for the minor-loop gain G m (s) = G 2 (s)G 4 (s) Gs Gs Gs G Mc Mu m m () () () ()?>>**w 1 Zs Ds GsGsGsGs () () ()() ()() = ++ 1 1 24 12 ? 2000 by CRC Press LLC and the major-loop gain we can rewrite the response due to D(s) as Over the range of frequency where *G Mc (jw)* >> 1, such that the corrected loop gain is large, the attenuation of a load disturbance is proportional to the uncorrected loop gain. This is generally larger over a wider frequency range than the corrected loop gain magnitude *G Mc (jw)*, which is what the attenuation would be if series compensation were used. Over the range of frequencies where the minor-loop gain is large but where the corrected loop gain is small, that is, where *G m (jw)* > 1 and *G Mc (jw)* < 1, we obtain for the approximate response due to the disturbance and the output disturbance is therefore seen to be attenuated by the minor-loop gain rather than unattenuated as would be the case if series compensation had been used. This is, of course, highly desirable. At frequencies where both the minor-loop gain transfer and the major-loop gain are small we have Z(s)/D(s) ? 1, and over this range of frequencies neither minor-loop compensation nor series equalization is useful in reducing the effect of a load disturbance. Thus, we have shown here that there are quite a number of advantages to minor-loop compensation as compared to series equalization. Of course, there are limitations as well. Summary In this section, we have examined the subject of linear system compensation by means of the frequency response method of Bode diagrams. Our approach has been entirely in the frequency domain. We have discussed a variety of compensation networks, including: 1. Gain attenuation 2. Lead networks 3. Lag networks 4. Lag-lead networks and composite equalizers 5. Minor-loop feedback Despite its age, the frequency domain design approach represents a most useful approach for the design of linear control systems. It has been tested and proven in a great many practical design situations. Defining Terms Bode diagram: A graph of the gain magnitude and frequency response of a linear circuit or system, generally plotted on log-log coordinates. A major advantage of Bode diagrams is that the gain magnitude plot will look like straight lines or be asymptotic to straight lines. H.W. Bode, a well-known Bell Telephone Laboratories researcher, published Network Analysis and Feedback Amplifier Design in 1945. The approach, first described there, has been refined by a number of other workers over the past half-century. Crossover frequency: The frequency where the magnitude of the open-loop gain is 1. Gs GsG s GsGs Mc () () () () () = + 12 24 1 Zs Ds G s G s mMc () () [ ()] () = + 1 1 Zs Ds Gs m () () ()? ? 2000 by CRC Press LLC Equalizer:A network inserted into a system that has a transfer function or frequency response designed to compensate for undesired amplitude, phase, and frequency characteristics of the initial system. Filter and equalizer are generally synonymous terms. Lag network: In a simple phase-lag network, the phase angle associated with the input-output transfer function is always negative, or lagging. Figures 100.13 and 100.14 illustrate the essential characteristics of a lag network. Lag-lead network:The phase shift versus frequency curve in a phase lag-lead network is negative, or lagging, for low frequencies and positive, or leading, for high frequencies. The phase angle associated with the input-output transfer function is always positive, or leading. Figures 100.16 and 100.17 illustrate the essential characteristics of a lag-lead network, or composite equalizer. Lead network: In a simple phase-lead network, the phase angle associated with the input-output transfer function is always positive, or leading. Figures 100.10 and 100.11 illustrate the essential characteristics of a lead network. Series equalizer: In a single-loop feedback system, a series equalizer is placed in the single loop, generally at a point along the forward path from input to output where the equalizer itself consumes only a small amount of energy. In Fig. 100.21, G 1 (s) could represent a series equalizer. G 1 (s) could also be a series equalizer if G 4 (s) = 0. Specification: A statement of the design or development requirements to be satisfied by a system or product. Systems engineering: An approach to the overall life cycle evolution of a product or system. Generally, the systems engineering process comprises a number of phases. There are three essential phases in any systems engineering life cycle: formulation of requirements and specifications, design and development of the system or product, and deployment of the system. Each of these three basic phases may be further expanded into a larger number. For example, deployment generally comprises operational test and evaluation, maintenance over an extended operational life of the system, and modification and retrofit (or replacement) to meet new and evolving user needs. Related Topic 11.1 Introduction References J.L. Bower and P.M. Schultheiss, Introduction to the Design of Servomechanisms, New York: Wiley, 1958. A.P. Sage, Linear Systems Control, Champaign, Ill.: Matrix Press, 1978. A.P. Sage, Systems Engineering, New York: Wiley, 1992. M.G. Singh, Ed., Systems and Control Encyclopedia, Oxford: Pergamon, 1987. Further Information Many of the practical design situations used to test the frequency domain design approach are described in the excellent classic text by Bower and Schultheiss [1958]. A rather detailed discussion of the approach may also be found in Sage [1978] on which this discussion is, in part, based. A great variety of control systems design approaches, including frequency domain design approaches, are discussed in a recent definitive control systems encyclopedia [Singh, 1987], and there are a plethora of new introductory control systems textbooks that discuss it as well. As noted earlier, frequency domain design, in particular, and control systems design, in general, constitute one facet of systems engineering effort, such as described in Sage [1992]. 100.4 Root Locus Benjamin C. Kuo Root locus represents a trajectory of the roots of an algebraic equation with constant coefficients when a parameter varies. The technique is used extensively for the analysis and design of linear time-invariant control ? 2000 by CRC Press LLC systems. For linear time-invariant control systems the roots of the characteristic equation determine the stability of the system. For a stable continuous-data system the roots must all lie in the left half of the s plane. For a digital control system to be stable, the roots of the characteristic equation must all lie inside the unit circle *z* = 1 in the z plane. Thus, in the s plane, the imaginary axis is the stability boundary, whereas in the z plane the stability boundary is the unit circle. The location of the characteristic equation roots with respect to the stability boundary also determine the relative stability, i.e., the degree of stability, of the system. For a linear time-invariant system with continuous data, the characteristic equation can be written as F(s) = P(s) + KQ(s) = 0 (100.33) where P(s) is an Nth-order polynomial of s, P(s) = s N + a 1 s N–1 + … + a N–1 s + a N (100.34) and Q(s) is the Mth-order polynomial of s, Q(s) = s M + b 1 s M–1 + … + b M–1 s + b M (100.35) where N and M are positive integers. The real constant K can vary from –¥ to +¥. The coefficients a 1 , a 2 , …, a N , b 1 , b 2 , …, b M are real. As K is varied from –¥ to +¥, the roots of Eq. (100.33) trace out continuous trajectories in the s plane called the root loci. The above development can be extended to digital control systems by replacing s with z in Eqs. (100.33) through (100.35). Root Locus Properties The root locus problem can be formulated from Eq. (100.33) by dividing both sides of the equation by the terms that do not contain the variable parameter K. The result is (100.36) For a closed-loop control system with the loop transfer function KG(s)H(s), where the gain factor K has been factored out, the characteristic equation is known to be the zeros of the rational function 1 + KG(s)H(s) = 0 (100.37) Since Eqs. (100.36) and (100.37) have the same form, the general root locus problem can be formulated using Eq. (100.36). Equation (100.37) is written (100.38) To satisfy the last equation, the following conditions must be met simultaneously: (100.39) 10+= KQ s Ps () () GsHs K () ()=- 1 Condition on magnitude: ** ** GsHs K () () = 1 ? 2000 by CRC Press LLC where K varies between –¥ and +¥. Conditions on angles: DG(s)H(s) = (2k + l)p K 3 0 = odd multiples of p rad (100.40) DG(s)H(s) = 2kp K £ 0 = even multiples of p rad (100.41) where k = 0, ±1, ±2, …, ± any integer. In general, the conditions on angles in Eqs. (100.40) and (100.41) are used for the construction of the root loci, whereas the condition on magnitude in Eq. (100.39) is used to find the value of K on the loci once the loci are drawn. Let KG(s)H(s) be of the form (100.42) Applying Eqs. (100.40) and (100.41) to the last equation, the angles conditions are K 3 0: DG(s)H(s) = D(s + a) – D(s + b) – D(s + c) = (2k + 1)p (100.43) K £ 0: DG(s)H(s) = D(s + a) – D(s + b) – D(s + c) = 2kp (100.44) where k = 0, ±1, ±2, … . The graphical interpretation of the last two equations is shown in Fig. 100.22. For the point s 1 to be a point on the root locus, the angles of the phasors drawn from the poles and zeros of G(s)H(s) to s 1 must satisfy Eq. (100.43) or (100.44) depending on the sign of K. Applying the magnitude condition of Eq. (100.39) to (100.42), the magnitude of K is expressed as (100.45) FIGURE 100.22Graphical interpretation of magnitude and angle conditions of root loci. KGsHs Ksa sbsc ()() () ()() = + ++ ** ** * ** K sbsc sa BC A = ++ + = × ? 2000 by CRC Press LLC where A, B, and C are the lengths of the phasors drawn from the poles and zeros of G(s)H(s) to the point s 1 . The following properties of the root loci are useful for sketching the root loci based on the pole-zero configuration of G(s)H(s). Many computer programs, such as the ROOTLOCI in the ACSP software package [Kuo, 1991b], are available for computing and plotting the root loci. The proofs and derivations of these properties can be carried out from Eqs. (100.39), (100.40), and (100.41) [Kuo, 1991a]. Starting Points (K 5 0 Points).The points at which K = 0 on the root loci are at the poles of G(s)H(s). Ending Points (K 56` Points).The points at which K = ±¥ on the root loci are at the zeros of G(s)H(s). The poles and zeros referred to above include those at s = ¥. Number of Root Loci.The total number of root loci of Eq. (100.37) equals the higher of the number of poles and zeros of G(s)H(s). Symmetry of Root Loci. The root loci are symmetrical with respect to the axes of symmetry of the pole-zero configuration of G(s)H(s). In general, the root loci are symmetrical at least to the real axis of the complex s plane. Asymptotes of the Root Loci.Asymptotes of the root loci refer to the behavior of the root loci at *s* = ¥ when the number of poles and zeros of G(s)H(s) is not equal. Let N denote the number of finite poles of G(s)H(s) and M be the number of finite zeros of G(s)H(s). In general, 2*N – M* of the loci will approach infinity in the s plane. The properties of the root loci at *s* = ¥ are described by the angles and the intersects of the asymptotes. When N 1 M, the angles of the asymptotes are given by where k = 0, 1, 2, . . ., *N – M* – 1. The asymptotes intersect on the real axis at (100.48) Root Loci on the Real Axis.The entire real axis of the s plane is occupied by the root loci. When K > 0, root loci are found in sections of the real axis to the right of which the total number of poles and zeros of G(s)H(s) is odd. When K < 0, root loci are found in sections to the right of which the total number of poles and zeros of G(s)H(s) is even. As a summary of the root locus properties discussed above, the properties of the root loci of the following equation are displayed in Fig. 100.23. s 3 + 2s 2 + 2s + K(s + 3) = 0 (100.49) Dividing both sides of the last equation by the terms that do not contain K we get (100.50) Thus, the poles of G(s)H(s) are at s = 0, s = –1 + j, and s = –1 – j. The zero of G(s)H(s) is at z = –3. As shown in Fig. 100.23, the K = 0 points on the root loci are at the poles of G(s)H(s), and the K = ±¥ points are at the zeros of G(s)H(s). Since G(s)H(s) has two zeros at s = ¥, two of the three root loci approach f p p k k NM K k NM K = + - 3 - £ ì í ? ? ? ? ? () (.) (.) 21 0 10046 2 0 10047 ** ** s= - - ?? finite poles of finite zeros of GsHs GsHs NM ()() ()() 11 3 22 2 +=+ + ++ KGsHs Ks ss s ()() () ? 2000 by CRC Press LLC infinity in the s plane. The root loci are symmetrical to the real axis of the s plane, since the pole-zero configuration of G(s)H(s) is symmetrical to the axis. The asymptotes of the two root loci that approach infinity are characterized by Eqs. (100.46) through (100.48). The angles of the asymptotes are: (100.51) (100.52) Thus, for K 3 0, the angles of the asymptotes are f 0 = 90° and f 1 = 270°. For K £ 0, f 0 = 0° and f 1 = 180°. The intersect of the asymptotes is at (100.53) The root loci on the real axis are as indicated in Fig. 100.23. Angles of Departure and Arrival. The slope of the root locus in the vicinity of a pole of G(s)H(s) is measured at the angle of departure and that in the vicinity of a zero of G(s)H(s) is measured at the angle of arrival. The angle of departure (arrival) of a root locus at a pole (zero) of G(s)H(s) is determined by assigning a point s 1 to the root locus that is very close to the pole (zero) and applying the angle conditions of Eqs. (100.40) or (100.41). Figure 100.24 illustrates the calculation of the angles of arrival and departure of the root locus at the pole s = –1 + j. We assign a point s 1 that is on the locus for K > 0 near the pole and draw phasors from all the poles and the zero of G(s)H(s) to this point. The angles made by the phasors with respect to the real axis must satisfy the angle condition in Eq. (100.46). Let the angle of the phasor drawn from –1 + j to s 1 be designated as q, which is the angle of departure; the angles drawn from the other poles and zero can be approximated by regarding s 1 as being very close to –1 + j. Thus, Eq. (100.46) leads to DG(s 1 )H(s 1 ) = –q – 135° – 90° + 26.6° = –180° (100.54) or q = –18.4°. For the angle of arrival of the root locus at the pole s = –1 + j, we assign a point s 1 on the root loci for K < 0 near the pole. Drawing phasors from all the poles and the zero of G(s)H(s) to s 1 and applying the angle condition in Eq. (100.47), we have FIGURE 100.23Some properties of the root loci of G(s)H(s) = K(s + 3)/[s(s 2 + 2s + 2)]. K k k K k k k k 3= + - = £= - = 0 21 31 01 0 2 31 01 : () , : , f p f p s= -+---- - = 113 31 1 2 jj() ? 2000 by CRC Press LLC DG(s 1 )H(s 1 ) = –q – 135° – 90° + 26.6° = 0° (100.55) Thus, the angle of arrival of the root locus for K < 0 is q = 198.4°. Similarly, we can show that the angles of arrival and departure of the root locus at s = –3 are 180° and 0°, respectively. Intersection of the Root Loci with the Imaginary Axis.The points where the root loci intersect the imaginary axis (if there is any) in the s plane, and the corresponding values of K, may be determined by means of the Routh-Hurwitz stability criterion [Kuo, 1991a]. The root locus program can also be used on a computer to give the intersects. The complete root loci in Fig. 100.25 show that the root loci intersect the imaginary axis at s = ±j2.45, and the value of K is 4. The system is stable for 0 £ K < 4. Breakaway Points of the Root Loci. Breakaway points on the root loci correspond to multiple-order roots of the equation. At a breakaway point several root loci converge and then break away in different directions. The breakaway point can be real or complex. The latter case must be in complex conjugate pairs. FIGURE 100.24Angle of arrival and departure calculations. FIGURE 100.25The complete root loci of G(s)H(s) = K(s + 3)/[s(s 2 + 2s + 2)]. ? 2000 by CRC Press LLC The breakaway points of the root loci of Eq. (100.37) must satisfy the following condition: (100.56) On the other hand, not all solutions of Eq. (100.56) are breakaway points. To satisfy as a breakaway point, the point must also lie on the root loci, or satisfy Eq. (100.37). Applying Eq. (100.56) to the function G(s)H(s) given in Eq. (100.50), we have the equation that the breakaway points must satisfy, 2s 3 + 11s 2 + 12s + 6 = 0 (100.57) The roots of the last equation are s = –4.256, –0.622 + j0.564 and –0.622 – j0.564. As shown in Fig. 100.25, only the solution s = –4.256 is a breakaway point on the root loci. Root Loci of Digital Control Systems The root locus analysis presented in the preceding subsections can be applied to digital control systems without modifying the basic principles. For a linear time-invariant digital control system, the transfer functions are expressed in terms of the z-transform rather than the Laplace transform. The relationship between the z- transform variable z and the Laplace transform variable s is z = e Ts (100.58) where T is the sampling period in seconds. Typically, the characteristic equation roots are solutions of the equation 1 + KGH(z) = 0 (100.59) where K is the variable gain parameter. The root loci for a digital control system are constructed in the complex z plane. All the properties of the root loci in the s plane apply readily to the z plane, with the exception that the stability boundary is now the unit circle *z * = 1. That is, the system is stable if all the characteristic equation roots lie inside the unit circle. As an illustration, the open-loop transfer function of a digital control system is given as (100.60) The characteristic equation of the closed-loop system is z(z – 1) + K(z + 0.1) = 0 The root loci of the system are shown in Fig. 100.26. Notice that the system is stable for 0 £ K < 2.22. When K = 2.22, one root is at z = –1, which is on the stability boundary. Design with Root Locus The root locus diagram of the characteristic equation of a closed-loop control system can be used for design purposes. The roots of the characteristic equation can be positioned in the s plane (or the z plane for digital control systems) to realize a certain desired relative stability or damping of the system. It should be kept in mind that the zeros of the closed-loop transfer function also affect the relative stability of the system, although the absolute stability is strictly governed by the characteristic equation roots. dG s H s ds () () = 0 Gz Kz zz () (.) () = + - 01 1 ? 2000 by CRC Press LLC As an illustrative example, the constant-damping ratio line for z = 0.5 is shown in Fig. 100.25. The intersect of the z = 0.5 line and the root locus corresponds to K = 0.38. Let us assume that we want to keep the relative damping at approximately 0.5 but the gain K should be increased tenfold. The following cascade controller is applied to the system [Evans, 1948]: (100.61) The open-loop transfer function of the compensated system is now (100.62) Figure 100.27 shows the root locus diagram of the compensated system for K 3 0. The shape of the complex root loci is not appreciably affected by the controller, but the value of K that corresponds to a relative damping ratio of 0.5 is now approximately 3.9. In a similar manner the root loci of digital control systems can be reshaped in the z plane for design. The constant-damping ratio locus in the z plane is shown in Fig. 100.26. Defining Terms Angles of departure and arrival:The slope of the root locus in the vicinity of a pole of G(s)H(s) is measured as the angle of departure, and that in the vicinity of a zero of G(s)H(s) is measured as the angle of arrival. Asymptotes of root loci:The behavior of the root loci at *s* = ¥ when the number of poles and zeros of G(s)H(s) is not equal. Breakaway points of the root loci: Breakaway points on the root loci correspond to multiple-order roots of the equation. Root locus: The trajectory of the roots of an algebraic equation with constant coefficient when a parameter varies. Related Topics 6.1 Definitions and Properties?12.1 Introduction FIGURE 100.26Root loci in the z plane for a digital control system. Gs s s c ()= + + 15 150 GsGsHs Ks s ss s s c ()()() .()(.) ( .)( ) = ++ +++ 01 3 02 002 2 2 2 ? 2000 by CRC Press LLC References and Further Information R. C. Dorf, Modern Control Systems, 5th ed., Reading, Mass.: Addison-Wesley, 1989. W. R. Evans, “Graphical analysis of control systems,” Trans. AIEE, vol. 67, pp. 547–551, 1948. B. C. Kuo, Automatic Control Systems, 6th ed., Englewood Cliffs, N.J.: Prentice-Hall, 1991a. B. C. Kuo, ACSP Software and Manual, Englewood Cliffs, N.J.: Prentice-Hall, 1991b. B. C. Kuo, Digital Control Systems, 2nd ed., New York: Holt, 1992a. B. C. Kuo, DCSP Software and Manual, Champaign, Ill.: SRL, Inc., 1992b. 100.5 Compensation Charles L. Phillips and Royce D. Harbor Compensation is the process of modifying a closed-loop control system (usually by adding a compensator or controller) in such a way that the compensated system satisfies a given set of design specifications. This section presents the fundamentals of compensator design; actual techniques are available in the references. A single-loop control system is shown in Fig. 100.28. This system has the transfer function from input R(s) to output C(s) (100.63) and the characteristic equation is 1 + G c (s)G p (s)H(s) = 0 (100.64) FIGURE 100.27Root loci of Eq. (100.62). Ts Cs Rs GsGs GsGsHs cp cp () () () ()() ()()() == +1 ? 2000 by CRC Press LLC where G c (s) is the compensator transfer function, G p (s) is the plant transfer function, and H(s) is the sensor transfer function. The transfer function from the disturbance input D(s) to the output is G d (s)/[1 + G c (s)G p (s)H(s)]. The function G c (s)G p (s)H(s) is called the open-loop function. Control System Specifications The compensator transfer function G c (s) is designed to give the closed-loop system certain specified character- istics, which are realized through achieving one or more of the following: 1. Improving the transient response. Increasing the speed of response is generally accomplished by increas- ing the open-loop gain G c (jw)G p (jw)H(jw) at higher frequencies such that the system bandwidth is increased. Reducing overshoot (ringing) in the response generally involves increasing the phase margin f m of the system, which tends to remove any resonances in the system. The phase margin f m occurs at the frequency w 1 and is defined by the relationship *G c (jw 1 )G p (jw 1 )H(jw 1 )* = 1 with the angle of G c (jw 1 )G p (jw 1 )H(jw 1 ) equal to (180° + f m ). 2. Reducing the steady-state errors. Steady-state errors are decreased by increasing the open-loop gain G c (jw)G p (jw)H(jw) in the frequency range of the errors. Low-frequency errors are reduced by increasing the low-frequency open-loop gain and by increasing the type number of the system [the number of poles at the origin in the open-loop function G c (s)G p (s)H(s)]. 3. Reducing the sensitivity to plant parameters. Increasing the open-loop gain G c (jw)G p (jw) H(jw) tends to reduce the variations in the system characteristics due to variations in the parameters of the plant. 4. Rejecting disturbances. Increasing the open-loop gain G c (jw)G p (jw)H(jw) tends to reduce the effects of disturbances [D(s) in Fig. 100.28] on the system output, provided that the increase in gain does not appear in the direct path from disturbance inputs to the system output. 5. Increasing the relative stability. Increasing the open-loop gain tends to reduce phase and gain margins, which generally increases the overshoot in the system response. Hence, a trade-off exists between the beneficial effects of increasing the open-loop gain and the resulting detrimental effects of reducing the stability margins. Design Design procedures for compensators are categorized as either classical methods or modern methods. Classical methods discussed are: FIGURE 100.28 A closed-loop control system. H(s) G d (s) Plant G c (s) Compensation Sensor R(s) M(s) C(s) D(s) ++ + – G p (s) ? 2000 by CRC Press LLC ?Ph ase-lag frequency response ?Phase-lead frequency response ?Phase-lag root locus ?Phase-lead root locus Modern methods discussed are: ?Pole placement ?State estimation ?Optimal Frequency Response Design Classical design procedures are normally based on the open-loop function of the uncompensated system, G p (s)H(s). Two compensators are used in classical design; the first is called a phase-lag compensator, and the second is called a phase-lead compensator. The general characteristics of phase-lag-compensated systems are as follows: 1.The low-frequency behavior of a system is improved. This improvement appears as reduced errors at low frequencies, improved rejection of low-frequency disturbances, and reduced sensitivity to plant parameters in the low-frequency region. 2.The system bandwidth is reduced, resulting in a slower system time response and better rejection of high-frequency noise in the sensor output signal. The general characteristics of phase-lead-compensated systems are as follows: 1.The high-frequency behavior of a system is improved. This improvement appears as faster responses to inputs, improved rejection of high-frequency disturbances, and reduced sensitivity to changes in the plant parameters. 2.The system bandwidth is increased, which can increase the response to high-frequency noise in the sensor output signal. The transfer function of a first-order compensator can be expressed as (100.65) where –w 0 is the compensator zero, –w p is its pole, and K c is its dc gain. If w p < w 0 , the compensator is phase- lag. The Bode diagram of a phase-lag compensator is given in Fig. 100.29 for K c = 1. It is seen from Fig. 100.29 that the phase-lag compensator reduces the high-frequency gain of the open-loop function relative to the low-frequency gain. This effect allows a higher low-frequency gain, with the advantages listed above. The pole and zero of the compensator must be placed at very low frequencies relative to the compensated-system bandwidth so that the destabilizing effects of the negative phase of the compensator are negligible. If w p > w 0 the compensator is phase-lead. The Bode diagram of a phase-lead compensator is given in Fig. 100.30 for K c = 1. It is seen from Fig. 100.30 that the phase-lead compensator increases the high-frequency gain of the open- loop function relative to its low-frequency gain. Hence, the system has a larger bandwidth, with the advantages listed above. The pole and zero of the compensator are generally difficult to place, since the increased gain of the open-loop function tends to destabilize the system, while the phase lead of the compensator tends to stabilize the system. The pole-zero placement for the phase-lead compensator is much more critical than that of the phase-lag compensator. A typical Nyquist diagram of an uncompensated system is given in Fig. 100.31. The pole and the zero of a phase-lag compensator are placed in the frequency band labeled A. This placement negates the destabilizing Gs Ks s c c p () (/ ) / = + + w w 0 1 1 ? 2000 by CRC Press LLC effect of the negative phase of the compensator. The pole and zero of a phase-lead compensator are placed in the frequency band labeled B. This placement utilizes the stabilizing effect of the positive phase of the compensator. PID Controllers Proportional-plus-integral-plus-derivative (PID) compensators are probably the most utilized form for com- pensators. These compensators are essentially equivalent to a phase-lag compensator cascaded with a phase- lead compensator. The transfer function of this compensator is given by (100.66) FIGURE 100.29 Bode diagram for a phase-lag compensator. FIGURE 100.30 Bode diagram for a phase-lead compensator. Gs K K s Ks cP I D () =++ ? 2000 by CRC Press LLC A block diagram portrayal of the compensator is shown in Fig. 100.32. The integrator in this compensator increases the system type by one, resulting in an improved low-frequency response. The Bode diagram of a PID compensator is given in Fig. 100.33. FIGURE 100.31A typical Nyquist diagram for G p (s)H(s). FIGURE 100.32Block diagram of a PID compensator. FIGURE 100.33Bode diagram of a PID compensator. ? 2000 by CRC Press LLC With K D = 0, the compensator is phase-lag, with the pole in (100.65) moved to w p = 0. As a result the compensator is type one. The zero of the compensator is placed in the low-frequency range to correspond to the zero of the phase-lag compensator discussed above. With K I = 0, the compensator is phase-lead, with a single zero and the pole moved to infinity. Hence, the gain continues to increase with increasing frequency. If high-frequency noise is a problem, it may be necessary to add one or more poles to the PD or PID compensators. These poles must be placed at high frequencies relative to the phase-margin frequency such that the phase margin (stability characteristics) of the system is not degraded. PD compensators realized using rate sensors minimize noise problems [Phillips and Harbor, 1991]. Root Locus Design Root locus design procedures generally result in the placement of the two dominant poles of the closed-loop system transfer function. A system has two dominant poles if its behavior approximates that of a second-order system. The differences in root locus designs and frequency response designs appear only in the interpretation of the control-system specifications. A root locus design that improves the low-frequency characteristics of the system will result in a phase-lag controller; a phase-lead compensator results if the design improves the high- frequency response of the system. If a root locus design is performed, the frequency response characteristics of the system should be investigated. Also, if a frequency response design is performed, the poles of the closed- loop transfer function should be calculated. Modern Control Design The classical design procedures above are based on a transfer-function model of a system. Modern design procedures are based on a state-variable model of the plant. The plant transfer function is given by (100.67) where we use u(t) for the plant input and y(t) for the plant output. If the system model is nth order, the denominator of G p (s) is an nth-order polynomial. The state-variable model, or state model, for a single-input–single-output plant is given by (100.68) y(t) = Cx(t) where x(t) is the n ′ 1 state vector, u(t) is the plant input, y(t) is the plant output, A is the n ′ n system matrix, B is the n ′ 1 input matrix, and C is the 1 ′ n output matrix. The transfer function of (100.67) is an input- output model; the state model of (100.68) yields the same input-output model and in addition includes an internal model of the system. The state model of (100.68) is readily adaptable to a multiple-input–multiple- output system (a multivariable system); for that case, u(t) and y(t) are vectors. We will consider only single- input–single-output systems. The plant transfer function of (100.67) is related to the state model of (100.68) by G p (s) = C(sI – A) –1 B (100.69) The state model is not unique; many combinations of the matrices A, B, and C can be found to satisfy (100.69) for a given transfer function G p (s). Ys Us Gs p () () ()= dt dt tut yt t x Ax B Cx () () () () () =+ = ? 2000 by CRC Press LLC Classical compensator design procedures are based on the open-loop function G p (s)H(s) of Fig. 100.28. It is common to present modern design procedures as being based on only the plant model of (100.68). However, the models of the sensors that measure the signals for feedback must be included in the state model. This problem will become more evident as the modern procedures are presented. Pole Placement Probably the simplest modern design procedure is pole placement. Recall that root locus design was presented as placing the two dominant poles of the closed-loop transfer function at desired locations. The pole-placement procedure places all poles of the closed-loop transfer function, or equivalently, all roots of the closed-loop system characteristic equation, at desirable locations. The system design specifications are used to generate the desired closed-loop system characteristic equation a c (s), where a c (s) = s n + a n–1 s n–1 + … + a 1 s + a 0 = 0 (100.70) for an nth-order plant. This characteristic equation is realized by requiring the plant input to be a linear combination of the plant states, that is, u(t) = –K 1 x 1 (t) – K 2 x 2 (t) – … – K n x n (t) = –Kx(t) (100.71) where K is the 1 ′ n feedback-gain matrix. Hence all states must be measured and fed back. This operation is depicted in Fig. 100.34. The feedback-gain matrix K is determined from the desired characteristic equation for the closed-loop system of (100.70): a c (s) = *sI – A + BK* = 0 (100.72) The state feedback gain matrix K which yields the specified closed-loop characteristic equation a c (s) is K = [0 0 … 0 1][B AB … A n–1 B] –1 a c (A) (100.73) where a c (A) is (100.70) with the scalar s replaced with the matrix A. A plant is said to be controllable if the inverse matrix in (100.73) exists. Calculation of K completes the design process. A simple computer algorithm is available for solving (100.73) for K. FIGURE 100.34Implementation of pole-placement design. ? 2000 by CRC Press LLC State Estimation In general, modern design procedures require that the state vector x(t) be fed back, as in (100.71). The measurement of all state variables is difficult to implement for high-order systems. The usual procedure is to estimate the states of the system from the measurement of the output y(t), with the estimated states then fed back. Let the estimated state vector be x?. One procedure for estimating the system states is by an observer, which is a dynamic system realized by the equations (100.74) with the feedback equation of (100.71) now realized by (100.75) The matrix G in (100.74) is calculated by assuming an nth-order characteristic equation for the observer of the form a e (s) = *sI – A + GC* = 0 (100.76) The estimator gain matrix G which yields the specified estimator characteristic equation a e (s) is G = a e (A)[C CA … CA n–1 ] –T [0 0 … 0 1] T (100.77) where [·] T denotes the matrix transpose. A plant is said to be observable if the inverse matrix in (100.77) exists. An implementation of the closed-loop system is shown in Fig. 100.35. The observer is usually implemented on a digital computer. The plant and the observer in Fig. 100.35 are both nth-order; hence, the closed-loop system is of order 2n. The observer-pole-placement system of Fig. 100.35 is equivalent to the system of Fig. 100.36, which is of the form of closed-loop systems designed by classical procedures. The transfer function of the controller-estimator (equivalent compensator) of Fig. 100.36 is given by G ec (s) = K[sI – A + GC + BK] –1 G (100.78) FIGURE 100.35Implementation of observer-pole-placement design. dt dt tutyt ? () () ? () () () x AGCx B G=- + + ut t() ? ()=-Kx ? 2000 by CRC Press LLC This compensator is nth-order for an nth-order plant; hence, the total system is of order 2n. The characteristic equation for the compensated system is given by *s I – A + BK ** s I – A + GC* = a c (s)a e (s) = 0 (100.79) The roots of this equation are the roots of the pole-placement design plus those of the observer design. For this reason, the roots of the characteristic equation for the observer are usually chosen to be faster than those of the pole-placement design. Linear Quadratic Optimal Control We define an optimal control system as one for which some mathematical function is minimized. The function to be minimized is called the cost function. For steady-state linear quadratic optimal control the cost function is given by (100.80) where Q and R are chosen to satisfy the design criteria. In general, the choices are not straightforward. Minimization of (100.80) requires that the plant input be given by u(t) = –R –1 B T M ¥ x(t) (100.81) where the n ′ n matrix M ¥ is the solution to the algebraic Riccati equation M ¥ A + A T M ¥ – M ¥ BR –1 B T M ¥ + Q = 0 (100.82) The existence of a solution for this equation is involved [Friedland, 1986] and is not presented here. Optimal control systems can be designed for cost functions other than that of (100.80). Other Modern Design Procedures Other modern design procedures exist; for example, self-tuning control systems continually estimate certain plant parameters and adjust the compensator based on this estimation. These control systems are a type of adaptive control systems and usually require that the control algorithms be implemented using a digital computer. These control systems are beyond the scope of this book (see, for example, Astrom and Wittenmark, 1984). Defining Term Compensation: The process of physically altering a closed-loop system such that the system has specified characteristics. This alteration is achieved either by changing certain parameters in the system or by adding a physical system to the closed-loop system; in some cases both methods are used. FIGURE 100.36 Equivalent system for pole-placement design. VRud T t ¥ ¥ =+ ò [()() ()]xQxtt tt 2 ? 2000 by CRC Press LLC Related Topics 100.3 Frequency Response Methods: Bode Diagram Approach?100.4 Root Locus References K. J. Astrom and B. Wittenmark, Computer Controlled Systems, Englewood Cliffs, N.J.: Prentice-Hall, 1984. W. L. Brogan, Modern Control Theory, Englewood Cliffs, N.J.: Prentice-Hall, 1985. R. C. Dorf, Modern Control Systems, 7th ed., Reading, Mass.: Addison-Wesley, 1995. G. F. Franklin, J. D. Powell, and A. Emami-Naeini, Feedback Control of Dynamic Systems, Reading, Mass.: Addison-Wesley, 1986. B. Friedland, Control System Design, New York: McGraw-Hill, 1986. B. C. Kuo, Automatic Control Systems, Englewood Cliffs, N.J.: Prentice-Hall, 1987. C. L. Phillips and R. D. Harbor, Feedback Control Systems, 2nd ed., Englewood Cliffs, N.J.: Prentice-Hall, 1991. 100.6 Digital Control Systems Raymond G. Jacquot and John E. McInroy The use of the digital computer to control physical processes has been a topic of discussion in the technical literature for over four decades, but the actual use of a digital computer for control of industrial processes was reserved only for massive and slowly varying processes such that the high cost and slow computing speed of available computers could be tolerated. The invention of the integrated circuit microprocessor in the early 1970s radically changed all that; now microprocessors are used in control tasks in automobiles and household appliances, applications where high cost is not justifiable. When the term digital control is used, it usually refers to the process of employing a digital computer to control some process that is characterized by continuous-in-time dynamics. The control can be of the open- loop variety where the control strategy output by the digital computer is dictated without regard to the status of the process variables. An alternative technique is to supply the digital computer with digital data about the process variables to be controlled, and thus the control strategy output by the computer depends on the process variables that are to be controlled. This latter strategy is a feedback control strategy wherein the computer, the process, and interface hardware form a closed loop of information flow. Examples of dynamic systems that are controlled in such a closed-loop digital fashion are flight control of civilian and military aircraft, control of process variables in chemical processing plants, and position and force control in industrial robot manipulators. The simplest form of feedback control strategy provides an on-off control to the controlling variables based on measured values of the process variables. This strategy will be illustrated by a simple example in a following subsection. In the past decade and a half many excellent textbooks on the subject of digital control systems have been written, and most of them are in their second edition. The texts in the References provide in-depth development of the theory by which such systems are analyzed and designed. A Simple Example Such a closed-loop or feedback control situation is illustrated in Fig. 100.37, which illustrates the feedback control of the temperature in a simple environmental chamber that is to be kept at a constant temperature somewhat above room temperature. Heat is provided by turning on a relay that supplies power to a heater coil. The on-off signal to the relay can be supplied by 1 bit of an output port of the microprocessor (typically the port would be 8 bits wide). A second bit of the port can be used to turn a fan on and off to supply cooling air to the chamber. An analog-to-digital (A/D) converter is employed to convert the amplified thermocouple signal to a digital word that is then supplied to the input port of the microprocessor. The program being executed by the microprocessor reads the temper- ature data supplied to the input port and compares the binary number representing the temperature to a binary ? 2000 by CRC Press LLC version of the desired temperature and makes a decision whether or not to turn on the heater or the fan or to do nothing. The program being executed runs in a continuous loop, repeating the operations discussed above. This simple on-off control strategy is often not the best when extremely precise control of the process variables is required. A more precise control may be obtained if the controlling variable levels can be adjusted to be somewhat larger if the deviation of the process variable from the desired value is larger. Single-Loop Linear Control Laws Consider the case where a single variable of the process is to be controlled, as illustrated in Fig. 100.38. The output of the plant y(t) is to be sampled every T seconds by an A/D converter, and this sequence of numbers will be denoted as y(kT), k = 0, 1, 2, . . . . The goal is to make the sequence y(kT) follow some desired known sequence [the reference sequence r(kT)]. Consequently, the sequence y(kT) is subtracted from r(kT) to obtain the so-called error sequence e(kT). The control computer then acts on the error sequence, using some control algorithms, to produce the control effort sequence u(kT) that is supplied to the digital-to-analog (D/A) converter which then drives the actuating hardware with a signal proportional to u(kT). The output of the D/A converter is then held constant on the current time interval, and the control computer waits for the next sample of the variable to be controlled, the arrival of which repeats the sequence. The most commonly employed control algorithm or control law is a linear difference equation of the form u(kT) = a n e(kT) + a n–1 e((k – 1)T) + … + a 0 e((k – n)T) + b n–1 u((k – 1)T) + … + b 0 u((k – n)T) (100.83) The question remains as to how to select the coefficients a 0 , . . ., a n and b 0 , . . ., b n–1 in expression (100.83) to give an acceptable degree of control of the plant. FIGURE 100.37Microprocessor control of temperature in a simple environmental chamber. FIGURE 100.38Closed-loop control of a single process variable. ? 2000 by CRC Press LLC Proportional Control This is the simplest possible control algorithm for the digital processor wherein the most current control effort is proportional to the current error or using only the first term of relation (100.83) u(kT) = a n e(kT) (100.84) This algorithm has the advantage that it is simple to program, while, on the other hand, its disadvantage lies in the fact that it has poor disturbance rejection properties in that if a n is made large enough for good disturbance rejection, the closed-loop system can be unstable (i.e., have transient responses which increase with time). Since the object is to regulate the system output in a known way, these unbounded responses preclude this regulation. PID Control Algorithm A common technique employed for decades in chemical process control loops is that of proportional-plus- integral-plus-derivative (PID) control wherein a continuous-time control law would be given by (100.85) This would have to be implemented by an analog filter. To implement the design in digital form the proportional term can be carried forward as in relation (100.84); however, the integral can be replaced by trapezoidal integration using the error sequence, while the derivative can be replaced with the backward difference resulting in a computer control law of the form [Jacquot, 1995] (100.86) where T is the duration of the sampling interval. The selection of the coefficients in this algorithm (K i , K d , and K p ) is best accomplished by the Ziegler-Nichols tuning process [Franklin et al., 1990]. The Closed-Loop System When the plant process is linear or may be linearized about an operating point and the control law is linear as in expressions (100.83), (100.84), or (100.86), then an appropriate representation of the complete closed- loop system is by the so-called z-transform. The z-transform plays the role for linear, constant-coefficient difference equations that the Laplace transform plays for linear, constant-coefficient differential equations. This z-domain representation allows the system designer to investigate system time response, frequency response, and stability in a single analytical framework. If the plant can be represented by an s-domain transfer function G(s), then the discrete-time (z-domain) transfer function of the plant, the analog-to-digital converter, and the driving digital-to-analog converter is (100.87) ut K et K e d K de dt pi d t () () ( )=+ + ò tt 0 ukT u k T K KT K T ekT KT K K T ek T K T ek T p id i p dd () ( )) () (( ) ) (( ) ) =-+++ ? è ? ? ? ÷ +-- ? è ? ? ? ÷ -+ - 1 2 2 2 12 Gz z z ZL Gs s () () – = -? è ? ? ? ÷ é ? ê ù ? ú ì í ? ü y t 1 1 ? 2000 by CRC Press LLC where Z(·) is the z-transform and L –1 (·) is the inverse Laplace transform. The transfer function of the control law of (100.83) is (100.88) For the closed-loop system of Fig. 100.38 the closed-loop z-domain transfer function is (100.89) where G(z) and D(z) are specified above. The characteristic equation of the closed-loop system is 1 + G(z)D(z) = 0 (100.90) The dynamics and stability of the system can be assessed by the locations of the zeros of (100.90) (the closed- loop poles) in the complex z plane. For stability the zeros of (100.90) above must be restricted to the unit circle of the complex z plane. A Linear Control Example Consider the temperature control of a chemical mixing tank shown in Fig. 100.39. From a transient power balance the differential equation relating the rate of heat added q(t) to the deviation in temperature from the ambient q(t) is given as (100.91) where t is the time constant of the process and mc is the heat capacity of the tank. The transfer function of the tank is (100.92) FIGURE 100.39A computer-controlled thermal mixing tank. Dz Uz Ez az az a zbz b n n n n n n n () () () == +++ --- - - - - 1 1 0 1 1 0 L L Mz Yz Rz GzDz GzDz () () () ()() ()() == +1 d dt mc qt q t q+= 11 () Q() () () s Qs Gs mc s == + 1 1 / /t ? 2000 by CRC Press LLC The heater is driven by a D/A converter, and the temperature measurement is sampled with an A/D converter. The data converters are assumed to operate synchronously, so the discrete-time transfer function of the tank and the two data converters is from expression (100.87): (100.93) If a proportional control law is chosen, the transfer function associated with the control law is the gain a n = K or D(z) = K (100.94) The closed-loop characteristic equation is from (100.90): (100.95) If a common denominator is found, the resulting numerator is (100.96) The root of this equation is (100.97) If this root location is investigated as the gain parameter K is varied upward from zero, it is seen that the root starts at z = e –T/t for K = 0 and moves to the left along the real axis as K increases. Initially it is seen that the system becomes faster, but at some point the responses become damped and oscillatory, and as K is further increased the oscillatory tendency becomes less damped, and finally a value of K is reached where the oscillations are sustained at constant amplitude. A further increase in K will yield oscillations that increase with time. Typical unit step responses for r(k) = 1 and T/t = 0.2 are shown in Fig. 100.40. It is easy to observe this tendency toward oscillation as K increases, but a problem that is clear from Fig. 100.40 is that in the steady state there is a persistent error between the response and the reference [r(k) = 1]. Increasing the gain K will make this error smaller at the expense of more oscillations. As a remedy for this steady-state error problem and control of the dynamics, a control law transfer function D(z) will be sought that inserts integrator action into the loop while simultaneously canceling the pole of the plant. This dictates that the controller have a transfer function of the form (100.98) Typical unit step responses are illustrated in Fig. 100.41 for several values of the gain parameter. The control law that must be programmed in the digital processor is u(kT) = u((k – 1)T) + K[e(kT) – e –T/t e((k – 1)T)] (100.99) Gz z Qz mc e ze T T () () () / / == - - - - Qt t t 1 1 1 0+ - - = - - K mc e ze T T t t t / / ze K mc e TT -+-= -// () tt t 10 ze K mc e T T =+ - - - / / () t t t 1 Dz Uz Ez Kze z T () () () () / == - - -t 1 ? 2000 by CRC Press LLC The additional effort to program this over that required to program the proportional control law of (100.94) is easily justified since K and e –T/t are simply constants. Defining Terms Digital computer: A collection of digital devices including an arithmetic logic unit (ALU), read-only memory (ROM), random-access memory (RAM), and control and interface hardware. Feedback control: The regulation of a response variable of a system in a desired manner using measurements of that variable in the generation of the strategy of manipulation of the controlling variables. Related Topics 8.1 Introduction ? 112.1 Introduction ? 112.3 The State of the Art in CACSD References K.J. Astrom and B. Wittenmark, Computer Controlled Systems: Theory and Design, Englewood Cliffs, N.J.: Prentice-Hall, 1984. FIGURE 100.40 Step responses of proportionally controlled thermal mixing tank. FIGURE 100.41 Step responses of the compensated thermal mixing tank. ? 2000 by CRC Press LLC G.F. Franklin, J.D. Powell, and M.L. Workman, Digital Control of Dynamic Systems, 2nd ed., Reading, Mass.: Addison-Wesley, 1990. C.H. Houpis and G.B. Lamont, Digital Control Systems: Theory, Hardware, Software, 2nd ed., New York: McGraw-Hill, 1992. R.G. Jacquot, Modern Digital Control Systems, 2nd ed., New York: Marcel Dekker, 1995. B.C. Kuo, Digital Control Systems, 2nd ed., Orlando, Fla.: Saunders, 1992. C.L. Phillips and H.T. Nagle, Digital Control System Analysis and Design, 3rd ed., Englewood Cliffs, N.J.: Prentice- Hall, 1995. R. J. Vaccaro, Digital: A State-Space Approach, New York: McGraw-Hill, 1995. Further Information The IEEE Control Systems Magazine is a useful information source on control systems in general and digital control in particular. Highly technical articles on the state of the art in digital control may be found in the IEEE Transactions on Automatic Control, the IEEE Transactions on Control Systems Technology, and the ASME Journal of Dynamic Systems, Measurement and Control. 100.7 Nonlinear Control Systems 3 Derek P. Atherton The Describing Function Method The describing function method, abbreviated as DF, was developed in several countries in the 1940s [Atherton, 1982], to answer the question: “What are the necessary and sufficient conditions for the nonlinear feedback system of Fig. 100.42 to be stable?” The problem still remains unanswered for a system with static nonlinearity, n(x), and linear plant G(s). All of the original investigators found limit cycles in control systems and observed that, in many instances with structures such as Fig. 100.42, the wave form of the oscillation at the input to the nonlinearity was almost sinusoidal. If, for example, the nonlinearity in Fig. 100.42 is an ideal relay, that is has an on-off characteristic, so that an odd symmetrical input wave form will produce a square wave at its output, the output of G(s) will be almost sinusoidal when G(s) is a low pass filter which attenuates the higher harmonics in the square wave much more than the fundamental. It was, therefore, proposed that the nonlinearity should be represented by its gain to a sinusoid and that the conditions for sustaining a sinusoidal limit cycle be evaluated to assess the stability of the feedback loop. Because of the nonlinearity, this gain in response to a sinusoid is a function of the amplitude of the sinusoid and is known as the describing function. Because describing function methods can be used other than for a single sinusoidal input, the technique is referred to as the single sinusoidal DF or sinusoidal DF. The Sinusoidal Describing Function For the reasons explained above, if we assume in Fig. 100.42 that x(t) = a cos q, where q = wt and n(x) is a symmetrical odd nonlinearity, then the output y(t) will be given by the Fourier series, (100.100) (100.101) 3 The material in this section was previously published by CRC Press in The Control Handbook, William S. Levine, Ed., 1996. yanbn nn n qqq () =+ = ¥ ? cos sin , 0 where a 0 0= , ? 2000 by CRC Press LLC (100.102) and (100.103) The fundamental output from the nonlinearity is a 1 cos q + b 1 sin q, so that the describing function, DF, defined as the fundamental output divided by the input amplitude, is complex and given by (100.104) which may be written (100.105) where (100.106) Alternatively, in polar coordinates, (100.107) where and (100.108) If n(x) is single valued, then b 1 = 0 and (100.109) FIGURE 100.42 Block diagram of a nonlinear system. ayd 1 0 2 1= ( ) ( ) ò pqqq p cos , byd 1 0 2 1= ( ) ( ) ò pqqq p sin . Na a jb a ( ) =- ( ) 11 N a N a jN a pq ( ) = ( ) + ( ) Na aa Na ba pq ( ) = ( ) =- 11 and . Na Mae ja ( ) = ( ) y() Ma a b a ( ) =+ ( ) 1 2 1 2 12 Y aba ( ) =- ( ) - tan . 1 11 ayd 1 0 2 4= ( ) ( ) ò pqq p cos ? 2000 by CRC Press LLC giving (100.110) Although Eqs. (100.102) and (100.103) are an obvious approach to evaluating the fundamental output of a nonlinearity, they are indirect, because one must first determine the output wave form y(q) from the known nonlinear characteristic and sinusoidal input wave form. This is avoided if the substitution q =cos –1 (x/a) is made. After some simple manipulations, (100.111) and (100.112) The function p(x) is the amplitude probability density function of the input sinusoidal signal given by (100.113) The nonlinear characteristics n p (x) and n q (x), called the inphase and quadrature nonlinearities, are defined by (100.114) and (100.115) where n 1 (x) and n 2 (x) are the portions of a double-valued characteristic traversed by the input for x · > 0 and x · < 0, respectively, When the nonlinear characteristic is single-valued, n 1 (x) = n 2 (x), so n p (x) = n(x) and n q (x) = 0. Integrating Eq. (100.111) by parts yields (100.116) where n¢(x) = dn(x)/dx and n(0 + ) = lim e?¥ n(e), a useful alternative expression for evaluating a 1 . An additional advantage of using Eqs. (100.111) and (100.112) is that they yield proofs of some properties of the DF for symmetrical odd nonlinearities. These include the following: 1. For a double-valued nonlinearity, the quadrature component N q (a) is proportional to the area of the nonlinearity loop, that is, (100.117) Na a a a y d ( ) == ( ) ( ) ò 1 0 2 4 pqq p cos a a xn x p x dx p a 1 0 4= ( ) ( ) ( ) ò banxdx q a 1 0 4= ( ) ( ) ò p . px a x ( ) = ( ) - ( ) - 1 22 12 p . nx nx nx p ( ) = ( ) + ( ) [] 12 2 nx nx nx q ( ) = ( ) - ( ) [] 21 2 an anxaxdx a 1 22 12 0 40 4= ( ) ( ) + ( ) ¢ ( ) - ( ) + ò pp Na a q ( ) =- ( )( ) 1 2 p area of nonlinearity loop ? 2000 by CRC Press LLC 2.For two single-valued nonlinearities n a (x) and n b (x), with n a (x) < n b (x) for all 0 < x < b, N a (a) < N b (a) for input amplitudes less than b. 3.For a single-valued nonlinearity with k 1 x < n(x) < k 2 x for all 0 < x < b, k 1 < N(a) < k 2 for input amplitudes less than b. This is the sector property of the DF; a similar result can be obtained for a double-valued nonlinearity [Cook, 1973]. When the nonlinearity is single valued, from the properties of Fourier series, the DF, N(a), may also be defined as: 1.the variable gain, K, having the same sinusoidal input as the nonlinearity, which minimizes the mean squared value of the error between the output from the nonlinearity and that from the variable gain, and 2.the covariance of the input sinusoid and the nonlinearity output divided by the variance of the input. Evaluation of the Describing Function To illustrate the evaluation of the DF two simple examples are considered. Saturation Nonlinearity To calculate the DF, the input can alternatively be taken as a sin q. For an ideal saturation characteristic, the nonlinearity output wave form y(q) is as shown in Fig. 100.43. Because of the symmetry of the nonlinearity, the fundamental of the output can be evaluated from the integral over a quarter period so that which, for a > d, gives where a = sin –1 d/a. Evaluation of the integrals gives which, on substituting for d, give the result (100.118) FIGURE 100.43Saturation nonlinearity. Na a yd () = () ò 4 0 2 p qqq p sin , Na a ma d m d () =+ é ? ê ù ? ú òò 4 2 0 2 p qq d qq a a p sin sin Na m () = ( ) -+ é ? ê ù ? ú 4 2 2 4 p aa da sin cos Na m () = ( ) + ( ) pa a22sin . ? 2000 by CRC Press LLC Because, for a < d, the characteristic is linear giving N(a) = m, the DF for ideal saturation is mN s (d/a) where (100.119) where a = sin –1 d/a. Alternatively, one can evaluate N(a) from Eq. (100.116), yielding Using the substitution x = a sin q, as before. Relay with Dead Zone and Hysteresis The characteristic is shown in Fig. 100.44 together with the corresponding input, assumed equal to a cos q, and the corresponding output wave form. Using Eqs. (100.102) and (100.103) over the interval -p/2 to p/2 and assuming that the input amplitude a is greater than d + D, FIGURE 100.44Relay with dead zone and hysteresis. Na a a s d d pa a d ( ) = < ( ) + [] ì í ? ? ? 1 22 ,, sin , for and 1 for >, Na aa a ma x dx () == ( ) - ( ) ò 1 222 12 0 4 p d . Na m d m () = ( ) = ( ) + ( ) ò 422 2 0 pq paa a cos sin ahd h 1 2 2 = ? è ? ? ? ÷ = ( ) + ( ) - ò pq pba a b cos sin sin , ? 2000 by CRC Press LLC Thus (100.120) For the alternative approach, one must first obtain the in-phase and quadrature nonlinearities shown in Fig. 100.45. Using Eqs. (100.111) and (100.112), and as before. The DF of two nonlinearities in parallel equals the sum of their individual DFs, a result very useful for determining DFs, particularly of linear segmented characteristics with multiple break points. Several procedures [Altherton, 1982] are available for approximating the DF of a given nonlinearity either by numerical integration or by evaluating the DF of an approximating nonlinear characteristic defined, for example, by a quantized characteristic, linear segmented characteristic, or Fourier series. Table 100.3 gives a list of DFs for some commonly used approximations of nonlinear elements. Several of the results are in terms of the DF for an ideal saturation characteristic of unit slope, N s (d/a), defined in Eq. (100.119). FIGURE 100.45Function n p (x) and n q (x) for the relay of Figure 100.44. where and andad bd pq p d d p a b =- ( ) [] =+ ( ) [] = ( ) =- ( ) + ( ) - - ? è ? ? ? ÷ = -- - ò cos cos , sin . 11 1 2 24 DD D D D aa bhd h aa ha Na h a aa jh a () =-+ ( ) é ? ê ù ? ú +-- ( ) é ? ê ù ? ú ì í ? ü y t - 2 2 2 12 2 2 12 2 p dd p DD D . a a xh pxdx xhpxdx h a aa a 1 2 2 12 2 2 12 42 2 = ( ) ( )() + () =-+ ( ) é ? ê ù ? ú +-- ( ) é ? ê ù ? ú ì í ? ü y t - + + òò d d d p dd D D D DD , , bahdxha a d 1 424= ( ) ( ) = = ( ) - + ò pp p d D D D Area of nonlinearity loop ? 2000 by CRC Press LLC TABLE 100.3 DFs of Single-Valued Nonlinearities a< d 1 N p = 0 d M+1 > a > d M a < d (2M + 1)d > a > (2M – 1)d N p = 0 n = (2m – 1)/2 a < d N p = 0 a > d N p = 4h(a 2 – d 2 ) 1/2 /a 2 p N p = 4h/ap N p = (4h/ap) + m a < d 1 N p = (4h/ap) + m 1 d M+1 > a > d M N p = (4h/ap) + m M+1 N P = mN s (d/a) N p = m[1 – N s (d/a) Naha mm m M p = () - () = ? 4 222 12 1 pd Nha an m M p = () - () = ? 4 222 12 1 pd +- ()()+ = ? mmN a jjsj i M 1 1 d ? 2000 by CRC Press LLC TABLE 100.3 (continued) DFs of Single-Valued Nonlinearities N p = (m 1 – m 2 )N s (d/a) + m 2 N p = m[N s (d 2 /a) – N s (d 1 /a)] N p = –m 1 N s (d 1 /a) + (m 1 – m 2 )N s (d 2 /a) + m 2 a < d N p = 0 a > d N p = 4h(a 2 – d 2 ) 1/2 /a 2 p + m – mN s (d/a) a < d N p = m 1 a > d N p = (m 1 – m 2 )N s (d/a) + m 2 + 4h(a 2 – d 2 ) 1/2 /a 2 p a < d N p = 4h/ap a > d N p = 4h/[a – (a 2 – d 2 ) 1/2 ]/a 2 p N p = (m 1 + m 2 )N s (d/a) – m 2 N s [(m 1 + m 2 )d/m 2 a] a < d N p = m 1 a > d N p = m 1 N s (d/a) – 4m 1 d(a 2 – d 2 ) 1/2 /a 2 p m > –2 G is the gamma function N ma mm ma m p m m m = + () + () [] + () [] = + () [] + () [] - - - G GG G G 1 23 21 2 2 22 32 1 1 1 p ? 2000 by CRC Press LLC Limit Cycles and Stability To investigate the possibility of limit cycles in the autonomous closed loop system of Fig. 100.42, the input to the nonlinearity n(x) is assumed to be a sinusoid so that it can be replaced by the amplitude-dependent DF gain N(a). The open loop gain to a sinusoid is thus N(a)G(jw) and, therefore, a limit cycle exists if (100.121) where G(jw) = G c (jw)G 1 (jw). As in general, G(jw) is a complex function of w and N(a) is a complex function of a, solving Eq. (100.121) will yield both the frequency w and amplitude a of a possible limit cycle. A common procedure to examine solutions of Eq. (100.120) is to use a Nyquist diagram, where the G(jw) and C(a) = –1/N(a) loci are plotted as in Fig. 100.46, where they are shown intersecting for a = a 0 and w = w 0 . The DF method indicates therefore that the system has a limit cycle with the input sinusoid to the nonlinearity, x, equal to a 0 sin (w 0 t + f), where f depends on the initial conditions. When the G(jw) and C(a) loci do not intersect, the DF method predicts that no limit cycle will exist if the Nyquist stability criterion is satisfied for G(jw) with respect to any point on the C(a) locus. Obviously, if the nonlinearity has unit gain for small inputs, the point (–1, j0) will lie on C(a) and may be used as the critical point, analogous to a linear system. For a stable case, it is possible to use the gain and phase margin to judge the relative stability of the system. However, a gain and phase margin can be found for every amplitude a on the C(a) locus, so it is usually appropriate to use the minimum values of the quantities [Atherton, 1982]. When the nonlinear block includes dynamics so that its response is both amplitude and frequency dependent, that is N(a, w), then a limit cycle will exist if (100.122) To check for possible solutions of this equation, a family of C(a, w) loci, usually as functions of a for fixed values of w, is drawn on the Nyquist diagram. An additional point of interest is whether when a solution to Eq. (100.120) exists the predicted limit cycle is stable. When there is only one intersection point, the stability of the limit cycle can be found using the Loeb criterion which states that if the Nyquist stability criterion indicates instability (stability) for the point on C(a) with a < a 0 and stability (instability) for the point on C(a) with a > a 0 the limit cycle is stable (unstable). FIGURE 100.46Nyquist plot showing solution for a limit cycle. NaGj ()( ) =-w 1 Gj Na Cawww ( ) =- ( ) = ( ) 1, ,. ? 2000 by CRC Press LLC When multiple solutions exist, the situation is more complicated and the criterion above is a necessary but not sufficient result for the stability of the limit cycle [Choudhury and Atherton, 1974]. Normally in these cases, the stability of the limit cycle can be ascertained by examining the roots of the characteristic equation (100.123) where N iy (a) is known as the incremental describing function (IDF). N iy (a) for a single valued nonlinearity can be evaluated from (100.124) where n¢(x) and p(x) are as previously defined. N ig (a) is related to N(a) by the equation (100.125) Thus, for example, for an ideal relay, making d = D = 0 in Eq. (100.120) gives N(a) = 4h/ap, also found directly from Eq. (100.116), and, substituting this value in Eq. (100.125) yields N ig (a) = 2h/ap. Some examples of feedback system analysis using the DF follow. Autotuning in Process Control In 1943 Ziegler and Nichols [1943] suggested a technique for tuning the parameters of a PID controller. Their method was based on testing the plant in a closed loop with the PID controller in the proportional mode. The proportional gain was increased until the loop started to oscillate and then the value of gain and the oscillation frequency were measured. Formulae were given for setting the controller parameters based on the gain named the critical gain, K c , and the frequency called the critical frequency, w c . Assuming that the plant has a linear transfer function G 1 (s), then K c is its gain margin and w c the frequency at which its phase shift is 180°. Performing this test in practice may prove difficult. If the plant has a linear transfer function and the gain is adjusted too quickly, a large amplitude oscillation may start to build up. In 1984 Astrom and Hagglund [1984] suggested replacing the proportional control by a relay element to control the amplitude of the oscillation. Consider therefore the feedback loop of Fig. 100.42 with n(x) an ideal relay, G c (s) = 1, and the plant with a transfer function G 1 (s) = 10/(s + 1) 3 . The C(a) locus, –1/N(a) = –ap/4h, and the Nyquist locus G(jw) in Fig. 100.47 intersect. The values of a and w at the intersection can be calculated from (100.126) which can be written (100.127) (100.128) 10+ ()() =NaGs ig Na nxpxdx i a a g () = ¢ ()() - ò NaNaadNada ig () = () + ( ) () 2. -= + ( ) ah j p w 4 10 1 3 Arg and 10 1 180 3 + ( ) ? è ? ? ? ? ÷ ÷ =° jw , a h p w 4 10 1 2 32 = + ( ) . ? 2000 by CRC Press LLC The solution for w c from Eq. (100.127) is tan –1 w c = 60°, giving w c = . Because the DF solution is approximate, the actual measured frequency of oscillation will differ from this value by an amount which will be smaller the closer the oscillation is to a sinusoid. The exact frequency of oscillation in this case will be 1.708 rads/sec in error by a relatively small amount. For a square wave input to the plant at this frequency, the plant output signal will be distorted by a small percentage. The distortion, d, is defined by (100.129) Solving Eq. (100.128) gives the amplitude of oscillation a as 5h /p. The gain through the relay is N(a) equal to the critical gain K c . In the practical situation where a is measured, K c equal to 4h /ap, should be close to the known value of 0.8 for this transfer function. If the relay has an hysteresis of D, then with d = 0 in Eq. (100.120) gives from which Thus on the Nyquist plot, C(a) is a line parallel to the real axis at a distance pD/4h below it, as shown in Fig. 100.47 for D = 1 and h = p/4 giving C(a) = –(a 2 – 1) 1/2 – j. If the same transfer function is used for the plant, then the limit cycle solution is given by (100.130) FIGURE 100.47 Nyquist plot 10/(s + 1) 3 and C(a) loci for D = 0 and 4h/p. 3 d = é ? ê ù ? ú M.S. value of signal — M.S. value of fundamental harmonic M.S. value of fundamental harmonic 12 Na ha a j h a ( ) = - ( ) - 4 4 22 12 D D pp Ca Na h aj ( ) = - ( ) = - - ( ) + é ? ê ù ? ú 1 4 22 12 p DD. -- ( ) -= + ( ) aj j 2 12 3 1 10 1 w ? 2000 by CRC Press LLC where w = 1.266, which compares with an exact solution value of 1.254, and a = 1.91. For the oscillation with the ideal relay, Eq. (100.123) with N ig (a) = 2h/ap shows that the limit cycle is stable. This agrees with the perturbation approach which also shows that the limit cycle is stable when the relay has hysteresis. Feedback Loop with a Relay with Dead Zone For this example the feedback loop of Fig. 100.42 is considered with n(x) a relay with dead zone and G(s) = 2/s(s + 1) 2 . From Equation 19.22 with D = 0, the DF for this relay, given by (100.131) is real because the nonlinearity is single valued. A graph of N(a) against a is in Fig. 100.48, and shows that N(a) starts at zero, when a = d, increases to a maximum, with a value of 2h/pd at a = d , and then decreases toward zero for larger inputs. The C(a) locus, shown in Fig. 100.49, lies on the negative real axis starting at –¥ and returning there after reaching a maximum value of –pd/2h. The given transfer function G(jw) crosses the negative real axis, as shown in Fig. 100.49, at a frequency of tan –1 w = 45°, that, is w = 1 rad/sec and, therefore, cuts the C(a) locus twice. The two possible limit cycle amplitudes at this frequency can be found by solving FIGURE 100.48N(a) for ideal relay with dead zone. FIGURE 100.49Two limit cycles: a 1 , unstable; a 2 , stable. Na ha a a () =- ( ) >4 22 12 2 dpdfor . 2 a ha 2 22 12 4 1 p d- ( ) = ? 2000 by CRC Press LLC which gives a = 1.04 and 3.86 for d = 1 and h = p. Using the perturbation method or the IDF criterion, the smallest amplitude limit cycle is unstable and the larger one is stable. If a condition similar to the lower amplitude limit cycle is excited in the system, an oscillation will build up and stabilize at the higher amplitude limit cycle. Other techniques show that the exact frequencies of the limit cycles for the smaller and larger amplitudes are 0.709 and 0.989, respectively. Although the transfer function is a good low pass filter, the frequency of the smallest amplitude limit cycle is not predicted accurately because the output from the relay, a wave form with narrow pulses, is highly distorted. If the transfer function of G(s) is K/s(s + 1) 2 , then no limit cycle will exist in the feedback loop, and it will be stable if that is, K < pd/h. If d = 1 and h = p, K < 1 which may be compared with the exact result for stability of K < 0.96. Stability and Accuracy Because the DF method is an approximate procedure, it is desirable to judge its accuracy. Predicting that a system will be stable, when in practice it is not, may have unfortunate consequences. Many attempts have been made to solve this problem, but those obtained are difficult to apply or produce too conservative results [Mess and Bergen, 1975]. The problem is illustrated by the system of Fig. 100.42 with a symmetrical odd single-valued nonlinearity confined to a sector between lines of slope k 1 and k 2 , that is, k 1 x < n(x) < k 2 x for x > 0. For absolute stability, the circle criterion requires satisfying the Nyquist criterion for the locus G(jw) for all points within a circle having its diameter on the negative real axis of the Nyquist diagram between the points (–1/k 1 , 0) and (–1/k 2 , 0), as shown in Fig. 100.50. On the other hand, because the DF for this nonlinearity lies within the diameter of the circle, the DF method requires satisfying the Nyquist criterion for G(jw) for all points on the circle diameter, if the autonomous system is to be stable. Therefore, for a limit cycle in the system of Fig. 100.42, errors in the DF method relate to its inability to predict a phase shift, which the fundamental harmonic may experience in passing through the nonlinearity, rather than an incorrect magnitude of the gain. When the input to a single-valued nonlinearity is a sinusoid together with some of its harmonics, the fundamental output is not necessarily in phase with the fundamental FIGURE 100.50Circle criterion and stability. Kd h ww p w 1 2 2 1 + ( ) < = , ? 2000 by CRC Press LLC input, that is, the fundamental gain has a phase shift. The actual phase shift varies with the harmonic content of the input signal in a complex manner, because the phase shift depends on the amplitudes and phases of the individual input components. From an engineering viewpoint one can judge the accuracy of DF results by estimating the distortion, d, in the input to the nonlinearity. This is straightforward when a limit-cycle solution is given by the DF method; the loop may be considered opened at the nonlinearity input, the sinusoidal signal corresponding to the DF solution can be applied to the nonlinearity, and the harmonic content of the signal fed back to the nonlinearity input can be calculated. Experience indicates that the percentage accuracy of the DF method in predicting the fundamental amplitude and frequency of the limit cycle is less than the percentage distortion in the fedback signal. As mentioned previously, the DF method may incorrectly predict stability. To investigate this problem, the procedure above can be used again, by taking, as the nonlinearity input, a sinusoid with amplitude and frequency corresponding to values of those parameters where the phase margin is small. If the calculated fedback distortion is high, say greater than 2% per degree of phase margin, the DF result should not be relied on. The limit-cycle amplitude predicted by the DF is an approximation to the fundamental harmonic. The accuracy of this prediction cannot be assessed by using the peak value of the limit cycle to estimate an equivalent sinusoid. It is possible to estimate the limit cycle more accurately by balancing more harmonics, as mentioned earlier. Although this is difficult algebraically other than with loops whose nonlinearity is mathematically simply described, for example a cubic, software is available for this purpose [McNamara and Atherton, 1987]. The procedure involves solving sets of nonlinear algebraic equations but good starting guesses can usually be obtained for the magnitudes and phases of the other harmonic components from the wave form fedback to the nonlinearity, assuming its input is the DF solution. Compensator Design Although the design specifications for a control system are often in terms of step-response behavior, frequency domain design methods rely on the premise that the correlation between the frequency and a step response yields a less oscillatory step response if the gain and phase margins are increased. Therefore the design of a suitable linear compensator for the system of Fig. 100.42 using the DF method, is usually done by selecting for example a lead network to provide adequate gain and phase margins for all amplitudes. This approach may be used in example 2 of the previous section where a phase lead network could be added to stabilize the system, say for a gain of 1.5, for which it is unstable without compensation. Other approaches are the use of additional feedback signals or modification of the nonlinearity n(x) directly or indirectly [Atherton, 1982; Gelb and van der Velde, 1968]. When the plant is nonlinear, its frequency response also depends on the input sinusoidal amplitude repre- sented as G (jw, a). In recent years several approaches [Nanka-Bruce and Atherton, 1990; Taylor and Strobel, 1984] use the DF method to design a nonlinear compensator for the plant, with the objective of closed-loop performance independent of the input amplitude. Closed-Loop Frequency Response When the closed-loop system of Fig. 100.42 has a sinusoidal input r(t) = R sin(wt + q), it is possible to evaluate the closed-loop frequency response using the DF. If the feedback loop has no limit cycle when r(t) = 0 and, in addition, the sinusoidal input r(t) does not induce a limit cycle, then, provided that G c (s)G 1 (s) gives good filtering, x(t), the nonlinearity input, almost equals the sinusoid a sin wt. Balancing the components of frequency w around the loop, (100.132) gR t aggMa tat ccc c sin _ sin sin wqf wffy w + ( ) - ( ) +++ ( ) [] = 1 1 ? 2000 by CRC Press LLC where G c (jw) = g c e jf c and G 1 (jw) = g 1 e jf 1 . In principle Eq. (100.132), which can be written as two nonlinear algebraic equations, can be solved for the two unknowns a and q and the fundamental output signal can then be found from (100.133) to obtain the closed-loop frequency for R and w. Various graphical procedures have been proposed for solving the two nonlinear algebraic equations resulting from Eq. (100.132) [Levinson, 1953; Singh, 1965; West and Douce, 1954]. If the system is lightly damped, the nonlinear equations may have more than one solution, indicating that the frequency response of the system has a jump resonance. This phenomenon of a nonlinear system has been studied by many authors, both theoretically and practically [Lamba and Kavanagh, 1971; West et al., 1954]. The Phase Plane Method The phase plane method was the first method used by control engineers for studying the effects of nonlinearity in feedback systems. The technique which can only be used for systems with second order models was examined and further developed for control engineering purposes for several major reasons, 1. The phase plane approach has been used for several studies of second order nonlinear differential equations arising in fields such as planetary motion, nonlinear mechanics and oscillations in vacuum tube circuits. 2. Many of the control systems of interest, such as servomechanisms, could be approximated by second order nonlinear differential equations. 3. The phase plane was particularly appropriate for dealing with nonlinearities with linear segmented characteristics which were good approximations for the nonlinear phenomena encountered in control systems. The next section considers the basic aspects of the phase plane approach but later concentration is focused on control engineering applications where the nonlinear effects are approximated by linear segmented nonlin- earities. Background Early analytical work [Andronov et al., 1966], on second order models assumed the equations (100.134) for two first-order nonlinear differential equations. Equilibrium, or singular points, occur when and the slope of any solution curve, or trajectory, in the x 1 – x 2 state plane is (100.135) A second order nonlinear differential equation representing a control system can be written ct aMag t a ( ) = ( ) + ( ) + [] 11 sin wy f ˙ , ˙ , xPxx xQxx 112 212 = ( ) = ( ) ˙˙ xx 12 0== dx dx x x Qx x Px x 2 1 2 1 12 12 == ( ) ( ) ˙ ˙ , , ? 2000 by CRC Press LLC (100.136) If this is rearranged as two first-order equations, choosing the phase variables as the state variables, that is x 1 = x, x 2 = x · ,then Eq. (100.136) can be written as (100.137) which is a special case of Eq. (100.135). A variety of procedures has been proposed for sketching state [phase] plane trajectories for Eqs. (100.135) and (100.137). A complete plot showing trajectory motions throughout the entire state (phase) plane is known as a state (phase) portrait. Knowledge of these methods, despite the improvements in computation since they were originally proposed, can be particularly helpful for obtaining an appreciation of the system behavior. When simulation studies are undertaken, phase plane graphs are easily obtained and they are often more helpful for understanding the system behavior than displays of the variables x 1 and x 2 against time. Many investigations using the phase plane technique were concerned with the possibility of limit cycles in the nonlinear differential equations When a limit cycle exists, this results in a closed trajectory in the phase plane. Typical of such investigations was the work of Van der Pol, who considered the equation (100.138) where m is a positive constant. The phase plane form of this equation can be written as (100.139) The slope of a trajectory in the phase plane is (100.140) and this is only singular (that is, at an equilibrium point), when the right hand side of Eq. (100.140) is 0/0, that is x 1 = x 2 = 0. The form of this singular point which is obtained from linearization of the equation at the origin depends upon m, being an unstable focus for m < 2 and an unstable node for m > 2. All phase plane trajectories have a slope of r when they intersect the curve (100.141) One way of sketching phase plane behavior is to draw a set of curves given for various values of r by Eq. (100.141) and marking the trajectory slope r on the curves. This procedure is known as the method of isoclines and has been used to obtain the limit cycles shown in Fig. 100.51 for the Van der Pol equation with m = 0.2 and 4. Piecewise Linear Characteristics When the nonlinear elements occurring in a second order model can be approximated by linear segmented characteristics then the phase plane approach is usually easy to use because the nonlinearities divide the phase ˙˙ , ˙ xfxx+ ( ) =0 ˙˙˙ ,xxxfxx 122 12 ==- ( ) ˙˙ ˙ xxxx-m- ( ) +=10 2 ˙ ˙ , xx xfxx xxx 12 212 1 2 21 1 = =- ( ) =m- ( ) - dx dx x x xx x x 2 1 2 1 1 2 21 2 1 == m- ( ) - ˙ ˙ rx xx x 21 2 21 1=m- ( ) - ? 2000 by CRC Press LLC plane into various regions within which the motion may be described by different linear second-order equations [Atherton, 1982]. The procedure is illustrated by the simple relay system in Fig. 100.52. The block diagram represents an “ideal” relay position control system with velocity feedback. The plant is a double integrator, ignoring viscous (linear) friction, hysteresis in the relay, or backlash in the gearing. If the system output is denoted by x 1 and its derivative by x 2 , then the relay switches when –x 1 –x 2 = ±1; the equations of the dotted lines are marked switching lines on Fig. 100.53. Because the relay output provides constant values of ±2 and 0 to the double integrator plant, if we denote the constant value by h, then the state equations for the motion are FIGURE 100.51Phase portraits of the Van der Pol equation for different values, of m. ? 2000 by CRC Press LLC (100.142) which can be solved to give the phase plane equation (100.143) which is a parabola for h finite and the straight line x 2 = x 20 for h = 0, where x 20 and x 10 are the initial values of x 2 and x 1 . Similarly, more complex equations can be derived for other second-order transfer functions. Using Eq. (100.143) with the appropriate values of h for the three regions in the phase plane, the step response for an input of 4.6 units can be obtained as shown in Fig. 100.53. In the step response, when the trajectory meets the switching line x 1 + x 2 = –1 for the second time, trajectory motions at both sides of the line are directed towards it, resulting in a sliding motion down the switching line. Completing the phase portrait by drawing responses from other initial conditions shows that the autonomous system is stable and also that all responses will finally slide down a switching line to equilibrium at x 1 = ±1. An advantage of the phase plane method is that it can be used for systems with more than one nonlinearity and for those situations where parameters change as functions of the phase variables. For example, Fig. 100.54 shows the block diagram of an approximate model of a servomechanism with nonlinear effects due to torque saturation and Coulomb friction. FIGURE 100.52Relay system. FIGURE 100.53Phase plane for relay system. ˙ ˙ xx xh 12 2 = = xxhxx 2 2 20 2 110 2-=- ( ) ? 2000 by CRC Press LLC The differential equation of motion in phase variable form is (100.144) where f s denotes the saturation nonlinearity and sgn the signum function, which is +1 for x 2 > 0 and –1 for x 2 < 0. There are six linear differential equations describing the motion in different regions of the phase plane. For x 2 positive, Eq. (100.144) can be written so that for (a) x 2 +ve, x 1 < –2, we have x · 1 = x 2 , x · 2 = 3/2, a parabola in the phase plane. (b) x 2 +ve|x 1 | < 2, we have x · 1 = x 2 , x · 2 + x 1 + 1/2 = 0. (c) x 2 +ve, x 1 > 2, we have x · 1 = x 2 , x · 2 = –5/2, a parabola in the phase plane. Similarly for x 2 negative, (d) x 2 –ve, x 1 – 2, we have x · 1 = x 2 , x · 2 = –5/2, a parabola in the phase plane. (e) x 2 –ve, |x 2 | < 2, we have x · 1 = x 2 , x · 2 + x 1 – 1/2 = 0, a circle in the phase plane. (f) x 2 –ve,x 1 > 2, we have x · 1 = x 2 , x · 2 = –3/2, a parabola in the phase plane. Because all the phase plane trajectories are described by simple mathematical expressions, it is straightforward to calculate specific phase plane trajectories. Discussion The phase plane approach is useful for understanding the effects of nonlinearity in second order systems, particularly if it may be approximated by a linear segmented characteristic. Solutions for the trajectories with other nonlinear characteristics may not be possible analytically so that approximate sketching techniques were used in early work on nonlinear control. These approaches are described in many books, for example, [Blaquiere, 1966; Cosgriff, 1958; Cunningham, 1958; Gibson, 1963; Graham and McRuer, 1961; Hayashi, 1964, Thaler and Pastel, 1962; West, 1960]. Although the trajectories are now easily obtained with modern simulation techniques, knowledge of the topological aspects of the phase plane are still useful for interpreting the responses in different regions of the phase plane and appreciating the system behavior. Related Topics 5.2 Limiters ? 12.1 Introduction ? 12.3 Lyapunov Stability Theory References A.A. Andronov, A.A. Vitt, and S.E. Khaikin, Theory of Oscillators, Reading, Mass.: Addison-Wesley, 1966. (First edition published in Russia in 1937.) FIGURE 100.54 Block diagram of servomechanism. ˙ xfx x s21 2 12=- ( ) - ( ) sgn ˙ xfx s11 12 0+ ( ) += ? 2000 by CRC Press LLC K.J. Astrom, and T. Haggland, Automatic tuning of single regulators, Budapest: Proc IFAC Congress, Vol. 4, 267–272, 1984. D.P. Atherton, Nonlinear Control Engineering, Describing Function Analysis and Design, London: Van Nostrand Reinhold, 1975. D.P. Atherton, Non Linear Control Engineering, Student Ed., New York: Van Nostrand Reinhold, 1982. A. Blaquiere, Nonlinear Systems Analysis, New York: Academic Press, 1966. S.K. Choudhury, and D.P. Atherton, “Limit cycles in high order nonlinear systems,” Proc. Inst. Electr. Eng., 121, 717–724, 1974. P.A. Cook, “Describing function for a sector nonlinearity,” Proc. Inst. Electr. Eng., 120, 143–144, 1973. R. Cosgriff, Nonlinear Control Systems, New York: McGraw-Hill, 1958. W.J. Cunningham, Introduction to Nonlinear Analysis, New York: McGraw-Hill, 1958. A. Gelb and W.E. van der Velde, Multiple Input Describing Functions and Nonlinear Systems Design, New York: McGraw-Hill, 1968. J.E. Gibson, Nonlinear Automatic Control, New York: McGraw-Hill, 1963. D. Graham and D. McRuer, Analysis of Nonlinear Control Systems, New York: John Wiley & Sons, 1961. C. Hayashi, Nonlinear Oscillations in Physical Systems, New York, McGraw-Hill, 1964. S.S. Lamba and R.J. Kavanagh, “The phenomenon of isolated jump resonance and its application,” Proc. Inst. Electr. Eng., 118, 1047–1050, 1971. E. Levinson, “Some saturation phenomena in servomechanims with emphasis on the techometer stabilised system,” Trans, Am. Inst. Electr. Eng., Part 2, 72, 1–9, 1953. O.P. McNamara, and D.P. Atherton, “Limit cycle prediction in free structured nonlinear systems,” IFAC Congress, Munich, 8, 23–28, July 1987. A.I. Mees and A.R. Bergen, “Describing function revisited,” IEEE Trans. Autom. Control, 20, 473–478, 1975. O. Nanka-Bruce and D.P. Atherton, “Design of nonlinear controllers for nonlinear plants,” IFAC Congress, Tallinn, 6, 75–80, 1990. T.P. Singh, “Graphical method for finding the closed loop frequency response of nonlinear feedback control systems,” Proc. Inst. Electr. Eng., 112, 2167–2170, 1965. J.H. Taylor and K.L. Strobel, “Applications of a nonlinear controller design approach based on the quasilinear system models,” Prof ACC, San Diego, 817–824, 1984. G.J. Thaler and M.P. Pastel, Analysis and Design of Nonlinear Feedback Control Systems, New York: McGraw- Hill, 1962. J.C. West, Analytical Techniques of Nonlinear Control Systems, London: E.U.P., 1960. J.C. West and J.L. Douce, “The frequency response of a certain class of nonlinear feedback systems,” Br. J. Appl. Phys., 5, 201–210, 1954. J.C. West, B.W. Jayawant, and D.P. Rea, “Transition characteristics of the jump phenomenon in nonlinear resonant circuits,” Proc. Inst. Electr. Eng., 114, 381–392, 1967. J.G. Ziegler and N.B. Nichols, “Optimal setting for automatic controllers,” Trans. ASME, 65, 433–444, 1943. Further Information Many control engineering text books contain material on nonlinear systems where the describing function is discussed. The coverage, however, is usually restricted to the basic sinusoidal DF for determining limit cycles in feedback systems. The basic DF method, which is one of quasilinearisation, can be extended to cover other signals, such as random signals, and also to cover multiple input signals to nonlinearities and feedback system analysis. The two books with the most comprehensive coverage of this are Gelb and Van der Velde [1968] and Atherton [1975]. More specialized books on nonlinear feedback systems usually cover the phase plane method and the DF, together with other topics such as absolute stability, exact linearization, etc. ? 2000 by CRC Press LLC 100.8 Optimal Control and Estimation John S. Bay and William T. Baumann Consider the closed-loop feedback control of linear time-invariant, multi-input/multi-output (MIMO) state- space systems of the form: (100.145) In this form, the vector x ∈ H5228 n represents the internal state, u ∈ H5228 p represents the input, and y ∈ H5228 q represents the measured outputs. It is well known that if the system in Eq. (100.145) is stabilizable (i.e., its unstable part is controllable), then it can be asymptotically stabilized with static state feedback. If it is detectable (i.e., its unstable part is observable), then a state estimator can be found, whose state variables asymptotically approach the true state variables in Eq. (100.145). However, merely determining the state feedback gain or observer gain leaves considerable design freedom for satisfying criteria other than stabilization and asymptotic observation. In this chapter section, we will provide results of some basic techniques of optimal control and estimation, which provide a mechanism to find optimal feedback and observer (estimator) gains according to selected optimality criteria. Linear Quadratic Regulators The linear quadratic regulator (LQR) problem is to find an optimal control input u*(t) that minimizes the performance criterion (100.146) where S and Q are symmetric, positive-semidefinite weighting matrices; and R is a symmetric, positive-definite weighting matrix. In this criterion, the term ?x T (t f )Sx(t f ) represents a penalty for the state at the final time t f being different from zero. The term inside the integral, x T (t)Qx(t), represents a penalty on the transient response of the state vector. The term u T (t)Ru(t) represents a penalty on the size of the control input u(t). We allow S and Q to be positive-semidefinite because we can generally tolerate unbounded state variables, provided they are not observed at the output. However, by forcing R to be positive-definite, we can guarantee that the process of minimizing Eq. (100.146) gives a bounded input. Minimization of this control energy is one of the primary reasons for using optimal control. The optimal control u*(t) can be found via a number of techniques, including dynamic programming [Bay, 1999] and variational techniques [Kirk, 1970]. The result of any of these methods is that the optimal control u*(t) is a linear function of the state vector (linear state feedback) of the form: (100.147) where the n × n matrix function P(t) satisfies the following equation: (100.148) Equation (100.148) is known as the differential matrix Riccati equation, and it is solved in backward time, with the end-time condition P(t f ) = S. ˙ xt t t tt t Ax Bu yCxDu ( ) ( ) ( ) ( ) ( ) ( ) =+ =+ J x u x t sx t x t Qx t u t Ru t dt ff t t o f , ( ) = ( ) ( ) + ( ) ( ) + ( ) ( ) [] ∫ 1 2 1 2 TT ut RBPtxt* T ( ) =? ( ) ( ) ?1 ˙ P PBR B P Q PA A P=??? ?1T T ? 2000 by CRC Press LLC It may be noted that a reasonable optimization criterion may have no finite final time t f . Instead, it may be desired that the controller be continually active, implying t f → ∞ and eliminating the possibility of the final state term in Eq. (100.146). In this case, the optimization criterion is more properly written as (100.149) Fortunately, this criterion simplifies the optimal control solution to the steady-state solution of the finite-time problem. That is, for the criterion of Eq. (100.149), the optimal control is (100.150) where in this case P is the matrix solution to the following algebraic Riccati equation: (100.151) Such a steady-state optimal solution exists whenever the system (A, B) is stabilizable. Furthermore, this constant P is the unique positive-definite solution of Eq. (100.151) if and only if the pair (A, T) is detectable, where T is defined as the square-root of Q, Q = T T T. Stabilizability of (A, B) ensures convergence of the cost criterion integral Eq. (100.149), and detectability of (A, T) guarantees that no unstable part of x(t) escapes integration as part of the integrand of Eq. (100.149). We should point out also that for the corresponding discrete-time system: (100.152) minimization of the cost criterion: (100.153) over all inputs u(k) results in the optimal control (100.154) where S is the solution to the corresponding discrete-time algebraic Riccati equation: (100.155) Note that in both the continuous- and the discrete-time infinite-horizon cases, the optimal control is actually static state feedback of the form u = Kx. J x u x t Qx t u t Ru t dt t o , ( ) = ( ) ( ) + ( ) ( ) [] ∞ ∫ 1 2 TT u t R B Px t* ( ) =? ( ) ?1T 0 1 =??? ? PBR B P Q PA A P TT xk Axk Buk yk Cxk Duk + ( ) = ( ) + ( ) ( ) = ( ) + ( ) 1 J x k Qx k u k Ru k kk = ( ) ( ) + ( ) ( ) [] = ∞ ∑ 1 2 0 TT u k R B SB B SAx k* ( ) =? + [] ( ) ? TT 1 S A SA A SB R B SB B SA Q=? + [] + ? TT T T 1 ? 2000 by CRC Press LLC Optimal Estimation: The Kalman Filter It was noted above that the optimal controller for the linear quadratic cost criterion takes the form of full-state feedback. However, it is often the case that the full state is not physically available for feedback. Rather, it is usually the output that is measurable, so that we prefer a technique that uses y(t) (and possibly u(t)) to construct the control signal instead of the state x(t). The simple solution to this problem is to design an observer (or estimator), which produces an estimated state vector ? x(t). If this estimate asymptotically approaches the true state x(t), then we can simply combine the observer and the state feedback to produce a feedback control u(t) = K ? x(t). That we can simply substitute the observed value ? x(t) for x(t) in the feedback function is a fortunate property called the separation principle, which ensures that our controller calculations and observer calculations do not interfere with one another. Just as we have improved static state feedback by introducing the optimal LQR above, we can take the principles of observer design and extend them with some guarantees of optimality. Such an optimal estimator is the Kalman filter. The Kalman filter is derived assuming the system model: (100.156) where v(t) and w(t) are two white, Gaussian, zero-mean, mutually uncorrelated noise signals with E[w(t)v T (τ)] = 0 and (100.157) and (100.158) where δ(t) is the Dirac delta. Noise v(t) is called the plant noise, and w(t) is the measurement noise, often representing sensor noise. (By assuming these signals are first passed through auxiliary filters with specified dynamics, these noises can be made to resemble harmonic or narrow-band disturbances.) The goal is now to design a system that produces an estimated state ? x(t) for Eq. (100.156) while rejecting the influence of the signals v(t) and w(t). To do this, we need the further assumptions that the system’s initial state is guessed to be x(t 0 ) = x 0 and that this guess is uncorrelated with the plant and measurement noise: (100.159) The covariance of the initial guess is defined as (100.160) The Kalman filter is the system that performs this estimation, rejecting the influence of the noise. The filter itself is given by the equation: (100.161) ˙ x t Ax t Bu t Gv t y t Cx t w t ( ) = ( ) + ( ) + ( ) ( ) = ( ) + ( ) Evt Evtv V t ( ) [] = ( ) ( ) [] =? ( ) 0, T τδτ Ewt Ewtw W t ( ) [] = ( ) ( ) [] =? ( ) 0, T τδτ Exv Exw 00 00 TT ττ ( ) [] = ( ) [] = Ex Ex x Ex P 0000 0 ? ( ) [] ? ( ) [] ? ? ? ? ? ? = T ? ? ˙ ?? xt Axt But Lt yt Cxt ( ) = ( ) + ( ) + ( ) ( ) ? ( ) [] ? 2000 by CRC Press LLC which can be seen to resemble the standard full-order observer [Bay, 1999]. However, rather than choosing an observer gain L in Eq. (100.161) to simply stabilize the error dynamics of the observer, the Kalman gain L(t) is (100.162) where P(t) is the solution to the following differential Riccati equation: (100.163) whose initial condition is P(t 0 ) = P 0 . For the discrete-time system (100.164) with assumptions analogous to Eq. (100.157) through (100.158) for plant noise v(k) and measurement noise w(t), and with E[(x 0 – E(x 0 )(x 0 – E(x 0 )) T ] ? S 0 , the Kalman filter is given by the two-stage estimator (100.165) (100.166) where the Kalman gain L(k + 1) is computed from the equation (100.167) In Eq. (100.167), the term S(k) is determined from the Riccati equation: (100.168) with initial condition S(k 0 ) = S 0 . (Note that these equations can be combined and rearranged to produce a number of alternate formulations.) Equation (100.165) is often referred to as the time update equation. It represents the estimate of the state vector that results from knowledge of the system dynamics. Equation (100.166) is sometimes called the mea- surement update equation because it revises the time update with a so-called innovations term L(y – y) that adjusts this time update according to the error between the output expected from the time update, y(k + 1), and the actual, measured output, y(k + 1). The matrices P(t) in the continuous-time filter, and S(k) in the discrete-time filter are the error covariance matrices for the state estimate. That is, (100.169) Lt PtC W ( ) = ( ) ?T1 ˙ P t AP t P t A P t C W CP t GVG ( ) = ( ) + ( ) ? ( ) ( ) + ?TT T1 x k Ax k Bu k Gv k y k Cx k w k + ( ) = ( ) + ( ) + ( ) ( ) = ( ) + ( ) 1 xk Axk Buk+ ( ) = ( ) + ( ) 1 ? ? xk xk Lk yk Cxk+ ( ) =+ ( ) ++ ( ) + ( ) ?+ ( ) [] 1111 1 L k AS k A GVG C C AS k A GVG C W T + ( ) = ( ) + [] ( ) + [] + {} ? 1 1 TTT TT S k I L k C AS k A GVG I L k C L k WL k + ( ) =? + ( ) [] ( ) + [] ?+ ( ) [] ++ ( ) + ( ) 11 1 11 TT T T Sk eke k t ete tEPE ( ) = ( ) ( ) [] ( ) = ( ) ( ) [] ?? TT and ? 2000 by CRC Press LLC where e(k) ? x(k) – ? x(k) and e(t) ? x(t) – ? x(t). Thus, the size of these matrices is an indicator of the error variance in various components in the estimates, and it can be seen in Eq. (100.161) and (100.167) that as these covariances decrease, the estimators rely less and less on the innovations term and more on the system dynamics for an accurate estimate. Early in the estimation process, the situation is usually reversed, with the innovations having the larger effect. Linear-Quadratic-Gaussian (LQG) Control It can be shown that the Kalman filter is the optimal estimator for the state of system Eq. (100.156) or (100.164) in the sense that it minimizes the squared error due to the noise input terms. Therefore, it becomes a likely candidate for combination with the LQR of the previous section. Together, the combination of LQR and Kalman filter is known as the LQG (linear-quadratic-Gaussian) controller, and is a useful controller in many situations. This controller is optimal in the sense that it minimizes the expected root-mean-square value of the optimization criterion (100.170) when the noise inputs are unit variance white noise. In the frequency domain, this is equivalent to minimizing the H 2 -norm (100.171) of the transfer function G from the white noise inputs to the output z = where Q = T T T, as before. We should point out that the so-called H 2 controller is equivalent to the LQG controller but provides a unified framework for control and observation that explains the striking similarity between the LQ controller equations and the Kalman filter equations (for example, in the Riccati equations of (100.148) and (100.63)). See Zhou and Doyle [1998] for further information. H ∞ Control The standard H ∞ control problem considers a system of the form (100.172) where w is a deterministic disturbance input vector, y is the measured output vector, and z is a vector of variables to be controlled. The objective of H ∞ control is to minimize the H ∞ norm (100.173) (where σ denotes the maximum singular value of a matrix, and sup denotes “supremum,” or least upper bound) of the transfer function G from the disturbance input w to the output z. In the time domain, the square of this lim T E T x t Qx t u t Ru t dt →∞ ( ) ( ) + ( ) ( ) ? ? ? ? ? ? ? ? ? ? ∫ 1 0 12 TT T G trace G j G j d 2 12 1 2 = π ( ) ( )( ) ? ? ? ? ? ? ? ? ? ? ?∞ ∞ ∫ * ωωω v w[] Tx Ru 12/ ? ? ? ? ? ? ˙ xAxBwBu zCxDu yCxDw =+ + =+ =+ 12 112 221 GGj ∞ = ( ) [] sup ω σω ? 2000 by CRC Press LLC norm corresponds to the maximum possible energy magnification between the input and the output, where energy is defined as the integral of the square of a signal; for example, ∫ 0 ∞ w T (t)w(t)dt. One of the major reasons for the development of H ∞ control is that many performance and robustness criteria for MIMO systems involve the maximum singular value of certain system transfer functions. To optimize the performance or robustness of these systems requires the minimization of the maximum singular value of these transfer functions, which is exactly the objective of H ∞ control. From a disturbance rejection point of view, the H ∞ controller can be used to minimize the root-mean-square value of the controlled variable z due to the worst-case unit-energy disturbance w. This is to be contrasted with LQG (or H 2 ) controller, which minimizes the average response to a unit-variance random disturbance. In practice, it is common to solve the suboptimal H ∞ problem where it is desired to find an output feedback controller such that H14067GH14067 ∞ < γ, where γ is specified by the designer. For large values of γ, there will always be a solution to the problem. In fact, as γ approaches infinity in the equations below, the central H ∞ controller will approach the LQG controller. By decreasing the value of γ until just before a solution to the problem ceases to exist, the designer can get as close to the optimal H ∞ controller as desired. To ensure that a solution for some value of γ exists, the following standard assumptions on the system are made [Green and Limebeer, 1995]: 1. The pair (A, B 2 ) is stabilizable and the pair (A, C 2 ) is detectable 2. D T 12 D 12 = I and D 21 D T 21 =I 3. Rank = n + m for all real ω 4. Rank = n + q for all real ω where n is the dimension of x, m is the dimension of u, and q is the dimension of y. Under these assumptions, it can be shown that there exists a stabilizing, measurement feedback solution to the suboptimal H ∞ control problem if and only if the following three conditions are met. 1. The algebraic Riccati equation X ∞ ? + ? T X ∞ + ? C T ? C – X ∞ (B 2 B T 2 – γ –2 B 1 B T 1 )X ∞ = 0 has a positive semi- definite solution such that ? – (B 2 B T 2 – γ –2 B 1 B T 1 )X ∞ is stable, where ? = A – B 2 D T 12 C 1 and ? C T ? C = C T 1 (I – D 12 D 12 T )C 1 . 2. The algebraic Riccati equation AY ∞ + Y ∞ A T + BB T – Y ∞ (C T 2 C 2 – γ –2 C T 1 C 1 )Y ∞ = 0 has a positive semi- definite solution such that A – Y ∞ (C T 2 C 2 – γ –2 C T 1 C 1 ) is stable, where A = A – B 1 D T 12 C 2 and BB T = B 1 (I – D T 21 D 12 )B T 1 . 3. ρ(X ∞ Y ∞ ) < γ 2 , where ρ(·) denotes the maximum of the absolute values of the matrix’s eigenvalues. The so-called central controller that solves the suboptimal H ∞ problem can be written in a form that closely resembles the state-estimate feedback form of the LQG controller: (100.174) where C 2z = C 2 + γ –2 D 21 B T 1 X ∞ , F ∞ = D T 12 C 1 + B T 2 X ∞ , and Z ∞ = Y ∞ (I – γ –2 X ∞ Y ∞ ) –1 . The dynamic part of the above compensator can be interpreted as an estimate of the state assuming that the worst-case disturbance w* is present. The control signal is a linear feedback of the estimated state, just as in the LQG case. Although the controller formulas above look more complicated than in the LQG case, this is largely due to the fact that the controlled variable z has a more general form in the H ∞ problem statement above. It should be noted, however, that unlike the LQG case, the solution of the Riccati equations is coupled in the H ∞ case due to condition 3 above. AjI B CD ? ? ? ? ? ? ? ω 2 112 AjI B CD ? ? ? ? ? ? ? ω 1 221 ? ˙ ?? * ?? * ? ? * ? xAxBw BuBD ZC yCxDw uFx wBXx z =+ ++ + [] ?? ( ) =? = ∞ ∞ ? ∞ 12121 2 221 2 1 TT T γ ? 2000 by CRC Press LLC Example Consider the following state-space system, which represents a plant with two lightly damped modes at approx- imately ω 1 ≈ 1 and ω 2 ≈ 3.2: (100.175) Here, the term d(t) = [d 1 (t) d 2 (t)] T is a vector whose first term represents the plant disturbance, and whose second term represents the measurement disturbance (deterministic). To pose the LQG control problem, we can propose minimizing the cost function (100.176) where (100.177) and T = [.5 0 –1 0]. However, Eq. (100.175) and (100.177) are also in the form of Eq. (100.172), facilitating an H ∞ design that minimizes the H ∞ norm of the transfer function from the disturbance d to the controlled variable z. We can compare this design to an H 2 controller design that minimizes the H 2 norm (Eq. (100.176)). The two controllers will therefore minimize different norms of the same transfer function. The results of this comparison are shown in the curves of Fig. 100.55. The distinguishing feature of this comparison is the flattening effect of the H ∞ controller. Although this is a plot of a frequency response magnitude and not the maximum singular value, it is apparent that the H ∞ controller is reducing the peak response, while the H 2 controller is reducing the average response and providing faster roll-off in the frequency domain. Other Approaches Although the LQG and H ∞ design methodologies are probably the most commonly used techniques for linear MIMO controller design, there are many other optimization-based techniques available. A well-developed theory exists for L 1 control, which minimizes the maximum magnification of the peak of the input signal using a linear programming algorithm [Dahleh and Diaz-Bobillo, 1995]. This approach departs from traditional controller design methodologies that have attempted to arrive at a set of formulas for the optimal controller. But with the advent of powerful computers, it makes sense to consider a problem solved if it can be reduced to a problem that can be efficiently solved using a computer algorithm. This is also the premise underlying the design methodology of linear matrix inequalities, in which the control design problem is reduced to a convex programming problem [Boyd et al., 1994]. ˙ . .. .. xt xt dt ut yt xt ( ) = ?? ?? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ( ) + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ( ) + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ( ) ( ) = []( ) + 0100 110 0 0001 0 0 10 1 00 10 00 20 0 1 0 1 10 50 0 0101 []( ) dt x T Tx ru dt z z dt TT T + ( ) = ∞∞ ∫∫ 2 00 z T x r u? = ? ? ? ? ? ? + ? ? ? ? ? ? 0 0 ? 2000 by CRC Press LLC Defining Terms Controllability: A linear system is said to be controllable if a control input exists that will drive a system with an arbitrary initial condition to a desired final state in a finite time. Stabilizability: A linear system is said to be stabilizable if its unstable part is controllable. Observability: A linear system is said to be observable if its state vector can be reconstructed from finite- length observations of its input and output. Detectability: A linear system is said to be detectable if its unstable part is observable. Positive-(semi)definite: A positive- (semi)definite matrix is a symmetric matrix A such that for any nonzero vector x, the quadratic form x T Ax is positive (non-negative). Observer (or estimator): A linear system whose state output approximates the state vector of a different system, rejecting noise and disturbances in the process. Singular value: Singular values are non-negative real numbers that measure the magnification effect of an operator in the different basis directions of the operator’s space. References B. D. O. Anderson and J. B. Moore, Optimal Filtering, Prentice-Hall, 1979. J. S. Bay, Fundamentals of Linear State Space Systems, WCB/McGraw-Hill, 1999. S. Boyd et al., Linear Matrix Inequalities in System and Control Theory, Society for Industrial and Applied Mathematics, 1994. A. E. Bryson, Jr. and Y. C. Ho, Applied Optimal Control, Hemisphere Publishing, 1975. M. Dahleh and I. J. Diaz-Bobillo, Control of Uncertain Systems: A Linear Programming Approach, Prentice-Hall, 1995. J. C. Doyle, K. Glover, P. P. Khargonekar, and B. A. Francis, State-space solutions to standard H 2 and H ∞ control problems, IEEE Trans. on Automatic Control, AC-34, 831–847, 1988. FIGURE 100.55 Frequency response of the closed-loop transfer function from d 1 to Tx, comparing the H 2 (LQG) and H ∞ control designs, using r = 0.1. ? 2000 by CRC Press LLC B. A. Francis, A Course in H ∞ Control Theory, Lecture Notes in Control and Information Sciences, Vol. 88, Springer-Verlag, 1987. M. Green and D. J. N. Limebeer, Linear Robust Control, Prentice-Hall, 1995. R. E. Kalman, A new approach to linear filtering and prediction problems, Transactions of the ASME, Journal of Basic Engineering, 82, 35–45, 1960. D. E. Kirk, Optimal Control Theory, Prentice-Hall, 1970. K. Zhou and J. C. Doyle, Essentials of Robust Control, Prentice-Hall, 1998. Further Information Optimal control and estimation is an actively growing field of control systems research. Classic texts in the area include [Kirk, 1970] and [Bryson and Ho, 1975] for optimal control systems, and [Anderson and Moore, 1979] for Kalman filtering, with the original source being [Kalman, 1960]. Further information on H 2 and H ∞ theory can be found in [Doyle, et al., 1989], [Zhou and Doyle, 1998], [Green and Limebeer, 1995], and [Francis, 1987]. 100.9 Neural Control Mo-Yuen Chow Artificial intelligence had strong ties with automatic control during its early development stages several decades ago. Typical examples of these ties are the development of cybernetics, robotics, and early learning systems. Recent efforts to incorporate aspects of artificial intelligence into the design and operation of automatic control systems have focused on using techniques such as artificial neural networks, fuzzy logic, and expert systems. The application of one or more of these techniques in the design of control systems has come to be known as intelligent control [Antsaklis et al., 1994], a term questioned by some for its implications. Whether or not such systems should be classified as intelligent, they represent significant contributions to the field of automatic control, as evidenced by the rapidly growing wealth of literature devoted to the successful application of such systems to complex control problems [Chow and Menozzi, 1994a; 1994b; Chow and Teeter, 1997; Chow and Yee, 1991; Hunt et al., 1992; Miller et al., 1990; Nguyen and Widrow, 1990; Psaltis et al., 1988; Werbos, 1990]. The nonlinear functional mapping properties of Artificial Neural Networks (ANNs) are central to their use in system identification and control [Narendra and Parthasarathy, 1990]. Although a number of key theoretical problems remain, results pertaining to the approximation capabilities of neural networks demonstrate that they have great promise in the modeling and control of nonlinear systems. The artificial neural network technology has become increasingly popular as a tool for performing tasks such as automatic control, system identification, pattern recognition, and time series prediction. Most of the conventional methods, such as PI control, are based on mathematical and statistical procedures for the modeling of the system and the estimation of the optimal controller parameters. In practice, the plant to be controlled is often highly nonlinear and a mathematical model may be difficult to derive. In such cases, conventional techniques may prove to be suboptimal and may lack robustness in the face of modeling error, because they are only as accurate as the model that was used to design them. With the advancement of technology, however, sophisticated control using artificial neural network techniques has been developed and used successfully to improve the control of systems that cannot be easily handled by conventional control, thus giving rise to terminology such as neural control and intelligent control. Usually, a human operator is responsible for adjusting the controller’s parameters in order to use his/her own idea of good performance. Indirectly, the operator is performing a minimization of a cost function based on his/her knowledge. More generally, when confronted with a problem situation, humans execute a mapping between a set of events and the set of corresponding appropriate actions. The appropriateness of these actions is due to some basic acquired knowledge, or even instinct, that guides the initial stages of the mapping. Then, through experience and a set of implicit guidelines, the human learns to perform a better mapping. This is an ongoing process throughout a person’s life. Similarly, an ANN, if given initial guidance, can learn to improve its performance through a set of guidelines, (e.g., minimize a cost function). In fact, a properly structured ANN can learn any arbitrarily complicated mapping [Cybenko, 1989; Rumelhart and McClelland, 1986; Werbos, 1974]. ? 2000 by CRC Press LLC Successful adaptation of online neural controllers in many cases requires careful and substantial design effort. One of the common practices is to pre-train a neural controller offline, based on some simplified design methods, in the same way that a PI controller is designed based on approximated system models that can provide reasonable control performance, before putting the neural controller for online fine-tuning [Chow and Menozzi, 1993; Chow and Teeter, 1995; Teeter et al., 1996]. This approach, as shown in Fig. 100.56, can speed up the neural controller online adaptation process and increase the closed-loop online adaptation stability because the initial weights of the neural controller are much closer to the optimal final weights (if they exist) after the pre-training process. By learning online, the ANN controller can adapt to changing operating environments. This chapter section briefly describes the feedforward net paradigm to facilitate the discussion of the neural observer and neural controller concepts in later sections. An example of using neural control for an HVAC system will then be provided to demonstrate its effectiveness for solving real-world problems. Brief Introduction to Artificial Neural Networks Increasing interest in studying the mechanisms and structure of the brain has led to the development of new computational models for solving problems such as pattern recognition, fast information processing, and adaptation. In the 1940s, McCulloch and Pitts studied the potential and capabilities of the interconnection of components based on a model of a biological neuron. Since then, many different models and architectures have been developed and analyzed for a variety of applications [Zurada, 1992]. One of the most common neuron models is shown in Fig. 100.57. The inputs x to the neuron model are scaled by connection weights w and summed. An additional input b, often referred to as a bias, is added to this sum and the result becomes the input to a function f(·), called the activation function, which computes the output of the neuron. The bias can be considered as a connection weight for a constant input of +1. The terms neuron model and neuron are used interchangeably in this chapter section. The individual neurons are not very powerful in terms of computation or representation, but the intercon- nection of neurons to form an artificial neural network (ANN) can provide a means of encoding complex relationships between input and output variables. Of the many ANN architectures that have been proposed, the multilayer feedforward artificial neural network (MFANN) shown in Fig. 100.58 is one of the most popular. Bias terms have been omitted in Fig. 100.58 for simplicity. The layers between the input and output layers of an MFANN are usually referred to as hidden layers because their inputs and outputs are not measurable at the inputs or outputs of the network. It has been shown that an MFANN with a single hidden layer can approximate any continuous function to an arbitrary degree of accuracy [Cybenko, 1989; Werbos, 1974]. The process of adjusting ANN connection weights in an effort to obtain a desired input/output mapping is usually referred to as training. Training represents an optimization problem for which the solution is a set of weights that minimizes some measure of approximation error. The choices of activation functions, number of neurons, error measures, and optimization methods can significantly affect training results. FIGURE 100.56 A two-step training process for a neural controller. ? 2000 by CRC Press LLC When the desired mapping is described by a set of input/output data, weights are usually modified after the presentation of each input/output pattern. This method is referred to as pattern update and represents an approximation of true gradient descent that is generally valid when a sufficiently small stepsize is used. Batch update, for which weight changes are accumulated over one sweep of the set of training patterns before being applied, is sometimes used in an effort to more closely mimic true gradient descent. Variants of back-propa- gation and other training methods can be found in [Zurada, 1992]. Neural Observer The nonlinear functional mapping properties of neural networks are central to their use in identification and control [Chow and Teeter, 1995; Hunt et al., 1992; Poggio and Girosi, 1990; Teeter and Chow, 1998; Teeter et al. 1994]. Although a number of key theoretical problems remain, results pertaining to the approximation capa- bilities of neural networks demonstrate that they have great promise in the modeling of nonlinear systems. An important question in system identification is whether a system under study can be adequately represented within a given model structure [Hunt et al., 1992]. In the absence of such concrete theoretical results for neural networks, it is usually assumed that the system under consideration belongs to the class of systems that the chosen network is able to represent. Two system identification techniques are now introduced: forward modeling and inverse modeling. The procedure of training a neural network to represent the forward dynamics of a system is often referred to as the forward system identification approach [Hunt et al., 1992]. A schematic diagram of this process is shown in Fig. 100.59. The neural network is placed in parallel with the system, and the error, e, between the system outputs, y, and network outputs, ? y, is used to train the network. This represents a classical supervised learning problem for which the teacher (i.e., the system) provides target values (i.e., system outputs) directly in the output coordinate system of the learner (i.e., the network model) [Jordan and Rumelhart, 1991]. In an inverse system identification approach, a network is trained in an effort to model the inverse of the plant mapping [Hunt et al., 1992]. One of the simplest approaches, known as direct inverse system identification, is shown schematically in Fig. 100.60. A synthetic training signal, s, is introduced to the system, and the system output, y, is used as the input to the network. The network output is compared to the training signal and this error is used to train the network. FIGURE 100.57 Computational model of a biological neuron. FIGURE 100.58 A multilayer feedforward ANN. ? 2000 by CRC Press LLC The inverse modeling structure shown in Fig. 100.60 tends to force the network to represent the inverse of the plant, but there are potential drawbacks to this approach. The training signal must be chosen to sample over a wide range of system inputs, and the actual operational inputs may be difficult to define a priori [Jordan and Rumelhart, 1991]. This point is strongly related to the concept of persistent excitation discussed in adaptive control literature. A second drawback is that an incorrect inverse model can be obtained if the nonlinear system mapping is not one-to-one. An approach called specialized inverse modeling has been proposed in an effort to overcome these problems. The details of this approach can be found in [Psaltis et al., 1988]. The neural network identification models can be used in the adaptive control of unknown nonlinear plants. Neural Control A method of direct adaptive control is depicted in Fig. 100.61. Methods for directly adjusting control parameters based on the output error e c are generally not available. This is because the unknown nonlinear plant in Fig. 100.61 lies between the controller and the output error. Until such methods are developed, adaptive control of nonlinear plants must be performed using indirect methods [Narendra and Parthasarathy, 1990]. Figure 100.62 depicts a general method of indirect adaptive control using artificial neural networks. Tapped delay lines (TDL) provide delayed inputs and outputs of the plant to the neural controller and neural observer. Error e i is used to adapt the neural observer, while the parameters of the neural observer along with FIGURE 100.59 Forward system identification approach. FIGURE 100.60 Direct inverse system identification approach. ? 2000 by CRC Press LLC error e c are used to adapt the neural controller. The model reference approach depicted in Fig. 100.62 is commonly referred to as implicit model-following in adaptive control literature [Narendra and Annaswamy, 1989]. This type of model-following is performed when the dynamics of the closed-loop system are forced to asymptotically match the dynamics of the reference model. The method of explicit model-following is depicted in Fig. 100.63. In this case, the reference model acts as a prefilter, and the value of the system output is forced to asymptotically track the value of the reference model output [Narendra and Annaswamy, 1989]. A neural control and identification scheme can be described as the following: at each time step, the plant states are measured and a neural network controller computes the plant inputs. Controller parameters are then FIGURE 100.61 Direct adaptive control. FIGURE 100.62 A method of indirect adaptive control using neural networks. FIGURE 100.63 Explicit model-following. ? 2000 by CRC Press LLC adjusted between samples so that the system approaches optimal performance with respect to a given cost index. The general form of the cost index used to adapt the neural controller is: (100.178) where L is the cost function as a function of system state x and control u. Thus, at each time step, the goal is to minimize J k subject to system dynamics and control constraints, with k denoting the current time step and N the prediction horizon. The back-propagation training algorithm can be used to adapt neural networks for the identification and control of nonlinear plants [Teeter and Chow, 1998; Werbos, 1990]. For system identification, a network can be trained offline using plant input/output data obtained from simulation of a mathematical model or from observation of the physical system. When the network is used for adaptive identification, training can be performed using batch update with a window of sampled data, or with the pattern update method in which training patterns consist of past inputs and outputs measured at each sample time. In order to adaptively train a neural controller using gradient descent, the partial derivatives of a cost index, J k , with respect to the network weights, w, must be obtained [Chow and Yee, 1991]. Let J k have the form J k = L(y(k+1),u(k+1))+L(y(k+2),u(k+2)) + … + L(y(k+n),u(k+n)) where k is the current sample. For simplicity of notation, let L(y(k+i),u(k+i)) be denoted by L(k+i). Application of the chain rule yields (100.179) The ?u(k)/?w term can be calculated using the backpropagation approach since the controller is a neural network. The ?y(k+i)/?u(k) and ?u(k+i)/?u(k) terms are obtained using the neural controller and identifiers. First, future inputs and outputs of the plant are predicted. The partial derivatives are then obtained by recursively computing the input/output sensitivities of the plant and controller through i samples. This approach is often referred to as back-propagation through time [Chow and Yee, 1991; Werbos, 1990]. The training algorithm resembles methods used by Nguyen and Widrow [1990] and others for training a neural controller to achieve an end goal. In this case, however, the output trajectory is of interest and the training is performed in realtime (i.e., output values must be repeatedly predicted rather than observed over several trials). A flowchart of the control methodology is shown in Fig. 100.64. After controller outputs are computed, the weights of the controller are adjusted N times before the next sample time. The value of N can be selected based on time constraints or convergence properties of the neural controller and observers. If N is large, the neural observers are inaccurate; and if a large prediction horizon is used, the adaptation of controller parameters may cause performance to deteriorate. HVAC Illustration In order to demonstrate the ability of the neural identification and control schemes to handle disturbances and changes in nonlinear plant dynamics, a Heat, Ventilation, and Air-Conditioning (HVAC) system where a thermal load is added and the actual output vs. the commanded output of each actuator is modified. The neural observers are adapted at each time step in one simulation, while adaptation of the observers is stopped at time step k = 200 in another simulation. Both simulations are performed for 1000 time steps. The reference trajectory is plotted in Fig. 100.65, along with the output of the system that uses nonadaptive identifiers after time step k = 200. Tracking errors for both simulations are plotted in Fig. 100.66, where T 3 is the temperature to be controlled and r is the reference signal to be tracked. The performance of the system with nonadaptive observers deteriorates due to the disturbance and the changes in plant dynamics. In this case, adapting the neural observers enables them to more accurately predict JLii k ik N = ( ) ( )( ) =+ ∑ xu, 1 ?+ ( ) ? = ?+ ( ) ? ( ) ? ( ) ? = ?+ ( ) ?+ ( ) ?+ ( ) ? ( ) + ?+ ( ) ?+ ( ) ?+ ( ) ? ( ) ? ? ? ? ? ? ? ? ? ( ) ? Lk i Lk i k kLki ki ki k Lk i ki ki k k w u u w y y uu u u u w ? 2000 by CRC Press LLC future states and inputs, and results in better tracking. It is important to select appropriate adaptation stepsizes for the observers because only one training pattern, consisting of the most recent set of plant inputs and states, is available for training. In order to compare the neural identification and control methodology with another potential control methodology, a PI-type controller has been designed for the HVAC system. Typical responses of the system with PI-type controller are shown in Fig. 100.67. The simple PI-type controller satisfies the conservative performance specifications for the cases tested, but does not always use its resources efficiently. For example, if the outside air temperature of the HVAC system is close to the steady-state reference temperature, it may be more efficient to increase the room volumetric flow rate for a period of time in order to reduce the amount of heating or cooling performed by the heat exchanger. The neural and PI-type control schemes are tested using initial conditions T 3 (0) = 15°C, a constant outside temperature of 22°C, and a steady-state reference value of 20°C. The tracking errors and heat exchanger outputs for both methods are shown in Figs. 100.68 and 100.69, respectively. FIGURE 100.64 Flowchart of the neural control training methodology. FIGURE 100.65 Reference and output trajectories using nonadaptive neural identifiers. ? 2000 by CRC Press LLC Conclusion The use of neural networks for system identification and control provides a means of adapting a controller on- line in an effort to minimize a given cost index. The cost index includes typical measures associated with system performance and can be modified without significantly increasing the computational complexity of the adap- tation process. For nonlinear systems, the identification networks demonstrate the capacity to learn changes in the plant dynamics. The performance of the neural control and identification methodology compares favorably with many types of conventional approaches. References Antsaklis, P. J., Albus, S., Lemmon, M. D., Mystel, A., Passino, K. M., Saridis, G. N., and Werbos, P. (1994). Defining intelligent control, IEEE Control Systems, 4, 5, 58–66. Chow, M.-Y. and Menozzi, A. (1993). Design Methodology of an Intelligent Controller Using Artificial Neural Networks, IECON’93, Maui, Hawaii. Chow, M.-Y. and Menozzi, A. (1994). A Self-Organized CMAC Controller for Robot Arm Movements, IEEE International Conference on Industrial Technology, Guangzhou, China. Chow, M.-Y. and Menozzi, A. (1994). Using a Cerebellar Model for FES Control of the Upper Limb, 16th Annual International IEEE Engineering in Medicine and Biology Society Conference, Baltimore, MD. FIGURE 100.66 Tracking errors for the system using adaptive and nonadaptive neural identifiers. FIGURE 100.67 Typical responses of the system with PI-type controller. FIGURE 100.68 Tracking errors of the PI-type and neural control systems. FIGURE 100.69 Heat exchanger outputs for the PI-type and neural control systems. ? 2000 by CRC Press LLC Chow, M.-Y. and Teeter, J. (1995). A knowledge-based approach for improved neural network control of a servomotor system with nonlinear friction characteristics, Mechatronics, 5(8), 949–962. Chow, M.-Y. and Teeter, J. (1997). Reduced-Order Functional Link Neural Network for HVAC Thermal System Identification and Modeling, 1997 International Conference on Neural Networks, Houston, TX. Chow, M.-Y. and Yee, S. O. (1991). An adaptive backpropagation through time training algorithm for a neural controller, 1991 IEEE International Symposium on Intelligent Control, Arlington, VA, 170–175. Cybenko, G. (1989). Approximation by superpositions of a sigmoidal function, Mathematics of Control, Signals, and Systems, 2, 303–314. Hunt, K. J., Sbarbaro, D., Zbikowski, R., and Gawthrop, P. J. (1992). Neural networks for control systems — a survey, Automatica, 28(6), 1083–1112. Jordan, M. I. and Rumelhart, D. E. (1991). Forward models: supervised learning with a distal teacher. Occasional Paper No. 40, Center for Cognitive Science, MIT. Miller, III, W. Thomas, Sutton, Richard S., Werbos, Paul J. (1990). Neural Networks for Control, The MIT Press, Cambridge, MA. Narendra, K. S. and Annaswamy, A. M. (1989). Stable Adaptive Systems, Prentice-Hall, Englewood Cliffs, NJ. Narendra, K. S. and Parthasarathy, K. (1990). Identification and control of dynamical systems using neural networks, IEEE Transactions on Neural Networks, 1(1), 4–27. Nguyen, D. and Widrow, B. (1990). The truck backer-upper: an example of self-learning in neural networks, Neural Networks for Control, The MIT Press, Cambridge, MA, 287–299. Poggio, T. and Girosi, F. (1990). Networks for approximation and learning, Proceedings of the IEEE, 78(9), 1481–1497. Psaltis, D., Sideris, A., and Yamamura, A. A. (1988). A multilayered neural network controller, IEEE Control Systems Magazine, 17–21. Rumelhart, D. E. and McClelland, J. L. (1986). Parallel Distributed Processing: Explorations in the Microstructure of Cognition, The MIT Press, Cambridge, MA. Teeter, J. and Chow, M.-Y. (1998). Application of functional link neural network to HVAC thermal dynamic system identification, IEEE Transactions on Industrial Electronics, 45(1), 170–176. Teeter, J., Chow, M.-Y., and Brickley, J. J. Jr. Use of a Fuzzy Gain Tuner for Improved Control of a DC Motor System with Nonlinearities, IEEE International Conference on Industrial Technology, Guangzhou, China. Teeter, J. T., Chow, M.-Y., and Brickley, J. J. Jr.(1996). A novel fuzzy friction compensation approach to improve the performance of a dc motor control system, IEEE Transactions on Industrial Electronics, 43(1), 113–120. Werbos, P. J. (1974). Beyond Regression: New Tools for Prediction and Analysis in Behavioral Science, Harvard University Press, Cambridge, MA. Werbos, P. J. (1990). Backpropagation through time: What it does and how to do it, Proceedings of the IEEE, 78(10), 1550–1560. Zurada, J. M. (1992). Introduction to Artificial Neural Systems, West Publishing Company, St. Paul, MN. ? 2000 by CRC Press LLC